I am trying to accurately cut loops out of a 3D PolyData with a loop which lies on the surface (or very close to it). Similar topics have been discussed a number of times before:
-
Punching a polygonal hole into a triangulated surface demonstrates how to do this for 2D using
vtkImprintFilter
-
Splitting cells on a PolyData Extrudes the loop into a skirt and uses
vtkIntersectionPolyDataFilter
scribe the loop onto the surface. However, with that approach I have not been able to “cut out” the surface patch which is inside the loop - I don’t know how to select only the points/cells which are inside. Also, I have encountered simple cases where the intersection fails entirely.
Is there any updated guidance on this? It feels like the new vtkImprintFilter
has the capacity to do this, but I have been unsuccessful getting it to work on 3D geometries. CC @will.schroeder who implemented this filter.
An example is shown below, using the surface.vtp
and boundaryline.vtp
files provided by @cobo on this post. This example produces the correct result when the target mesh and loop are co-planar, but not generally.
#!/usr/bin/env python
import vtk
import pyvista as pv
import numpy as np
target = pv.read("surface.vtp")
loop = pv.read("boundaryline.vtp")
# # make the data planar - the imprint works fine if these lines are un-commented!
# target.points[:, 2] = 0.
# loop.points[:, 2] = 0.
# ray-cast the loop points onto the target
origins = loop.points
direction = np.array([0, 0, 1])
length = 100.
intersections = []
for p in loop.points:
# uses `vtkOBBTree.IntersectWithLine`
intersection, _ = target.extract_surface().ray_trace(p + direction * length,
p - direction * length,
first_point=True)
intersections.append(intersection)
intersections = np.array(intersections)
# update point data
loop_on_target = loop.copy(deep=True)
loop_on_target.points = intersections
# extract contour loop from loop polydata
loopExtract = vtk.vtkContourLoopExtraction()
loopExtract.SetInputDataObject(loop_on_target)
loopExtract.Update()
# clean
cleanLoop = vtk.vtkCleanPolyData()
cleanLoop.SetInputConnection(loopExtract.GetOutputPort())
# imprint
imp = vtk.vtkImprintFilter()
imp.SetTargetData(target)
imp.SetImprintConnection(cleanLoop.GetOutputPort())
imp.SetTolerance(0.01)
imp.BoundaryEdgeInsertionOn()
imp.TriangulateOutputOn()
imp.SetOutputTypeToImprintedRegion()
imp.Update()
# get result of filter pipeline
res = pv.filters._get_output(imp)
# plot
p = pv.Plotter()
p.add_mesh(target, style='wireframe')
p.add_mesh(res, color='green')
p.add_mesh(res.copy(), color='w', style='wireframe')
p.add_mesh(loop_on_target, color='red', style='wireframe', line_width=5)
p.view_xy()
p.show()
Correct result when target mesh and loop are co-planar. The green shaded region is the imprinted region.
Wrong result when the target mesh is 3D and the loop is projected onto the mesh.