I’d like to “split” cells on a PolyData surface, adding vertices and cells. For example, in the screenshot below, I am starting with the black mesh, and I’d like to add all the points where the pink lines intersect the mesh, and adjust the cells connectivity appropriately.
Is there an existing filter to do something like this, or a method that anybody could suggest?
If you just want to add know the point ids inside the square, you could use a vtkAbstractPointLocator. You can use any concrete subclass, but if you are in 2D, I’d recommend using vtkStaticPointLocator2D.
You can use the method vtkAbstractPointLocator::FindPointsWithinRadius, querying with the center of your rectangle with a radius equal to half its diagonal. This method outputs a list of point ids that are within the input radius. Then you can filter them out by checking if they fit inside the rectangle (threshold x and y coordinates separately).
I’m not super familiar with cell locators, but they exist, and you can take a look if you want to extract the cells as well. Look for vtkAbstractCellLocator.
The short answer is no: there is no such filter in VTK. So you will have to do the “hard work” manually by creating a new mesh (polydata), adding all the existing points and cells that are not intersected by the square (using a cell locator), then perform the intersection, adding the intersection points and the connectivity of the new cells.
VTK (in particular the vtkPolyData data structure) is clearly not the best tool to perform mesh manipulation but this is doable with a bit of work.
Thanks @Yohann_Bearzi and @Joachim_P!
I am interested in the general 3D case; the original figure I provided was not very clear.
More formally, the problem I need to solve is this - given a PolyData and a closed polyline which lies on that surface, build a new PolyData containing:
I am trying to avoid doing as much of the hard work as possible 
. I could be way off, but I feel like a lot of the machinery already exists in either vtkClipPolyData or vtkBooleanOperationPolyDataFilter.
I am also open to suggestions if if you are aware of other libraries which might be able to accomplish this.
My first idea was to use vtkCookieCutter but it doesn’t to work as expected - the other possibility is to cut using a loop of points:
import vtk
from vedo import *
mesh = Mesh('https://www.dropbox.com/s/bu9otn8mcbtbtll/dolfin_fine.vtk').lw(0.1)
pts = [
[ 8.53826183e-01, 1.90910461e-01, 0],
[ 8.55859667e-01, 8.72127538e-01, 0],
[ 7.50018833e-02, 8.68060570e-01, 0],
[ 1.01437173e-01, 6.07674966e-02, 0],
]
cookie = Line(pts, closed=True).lw(3).c('red')
# ..this doesn't quite work:
# ccut = vtk.vtkCookieCutter()
# ccut.SetInputData(mesh.polydata())
# ccut.SetLoopsData(cookie.triangulate().polydata())
# ccut.Update()
# cmesh = Mesh(ccut).lw(0.1).c('tomato')
# ..this uses vtkSelectPolyData with SetLoop() + vtkClipPolyData
cmesh = mesh.clone().cutWithPointLoop(cookie, invert=1).c('tomato')
show([[mesh, cookie], cmesh], N=2, axes=1)
the vtkCookieCutter filter (commented out part above) gives me:
maybe I’m not using it properly…
In 3D the point loop method should work, see an example here (search “loop”).
…sorry I just realized I skipped your last message… maybe this fits better what you are looking for:
import vtk
from vedo import *
mesh = Mesh('https://www.dropbox.com/s/bu9otn8mcbtbtll/dolfin_fine.vtk').lw(0.1)
pts = [
[ 8.53826183e-01, 1.90910461e-01, 0],
[ 8.55859667e-01, 8.72127538e-01, 0],
[ 7.50018833e-02, 8.68060570e-01, 0],
[ 1.01437173e-01, 6.07674966e-02, 0],
]
line = Line(pts, closed=True).lw(3).c('red')
invert=1
onBoundary=0
ippd = vtk.vtkImplicitSelectionLoop()
ippd.SetLoop(line.polydata().GetPoints())
ippd.AutomaticNormalGenerationOn()
clipper = vtk.vtkExtractPolyDataGeometry()
clipper.SetInputData(mesh.polydata())
clipper.SetImplicitFunction(ippd)
clipper.SetExtractInside(invert)
clipper.SetExtractBoundaryCells(onBoundary)
clipper.Update()
cmesh = Mesh(clipper).lw(1).c('t')
show([[mesh, line], [cmesh, line]], N=2, axes=1)
Thanks for the responses, and for introducing me to some new VTK filters and vedo!
I am a bit hesitant to use clipping with an implicit function due to the discussion in https://www.vtkjournal.org/browse/publication/797.
I did find a solution using vtkIntersectionPolyDataFilter with the following approach:
The example below (uses pyvista) adds the edges of a trapezoid projected onto a sphere.
import pyvista as pv
import vtk
surf = pv.Sphere(center=(0, 0, 0), radius=200)
polyline = pv.PolyLine([(10, 0, 0),
(10, 50, 0),
(10, 25, 100),
(10, 0, 100)], closed=True)
skirt = polyline.extrude((400, 0, 0)).triangulate().clean().extract_surface()
f = vtk.vtkIntersectionPolyDataFilter()
f.SetInputDataObject(0, surf)
f.SetInputDataObject(1, skirt)
f.SetComputeIntersectionPointArray(True)
f.Update()
intersection_polyline = pv.wrap(f.GetOutputDataObject(0)) # intersection of the two surfaces
surf_remeshed = pv.wrap(f.GetOutputDataObject(1)) # the surface which has been remeshed along the split
skirt_remeshed = pv.wrap(f.GetOutputDataObject(2)) # the skirt which has been remeshed along the split
p = pv.Plotter(shape=(1, 2))
p.subplot(0, 0)
p.add_mesh(surf, style='wireframe', color='g')
p.add_mesh(surf.copy(), opacity=0.5, color='g')
p.add_mesh(skirt, style='wireframe', color='b')
p.add_mesh(skirt.copy(), opacity=0.5, color='b')
p.subplot(0, 1)
p.add_mesh(surf_remeshed, style='wireframe')
p.add_mesh(surf_remeshed.copy(), opacity=1.0, color='g')
p.link_views()
p.show()
I had the same problem with vtkCookieCutter. There was a bug in an internal method dealing with line intersections which is now fixed, like 4 weeks ago. You could try using/building the latest vtk version from gitlab and see if that fixes things for you.
@whophil, If you’d like to try another alternative, might I suggest TriSurfaceCutter?
All you need to do is turn on Embedding and turn off Removal.
import tsc
...
cutter = tsc.tscTriSurfaceCutter()
cutter.SetInputConnection(triangulation.GetOutputPort())
cutter.SetLoopsConnection(1, loops.GetOutputPort())
cutter.RemoveOff()
cutter.EmbedOn()
...
It is not yet integrated in VTK, so you have to build it separately and use the optional plugin with ParaView or in a Python/C++ application.
Cons:
Pros:
Embed: On | Remove: Off | Invert: NA
Embed: On | Remove: On | Invert: On
Embed: On | Remove: On | Invert: Off
Embed: Off | Remove: On | Invert: On
Embed: Off | Remove: On | Invert: Off
4, 5 look that way since no cell-line intersection was computed (Embed: Off)
I think the first set of options is sufficient for you!
