This is very similar to How to extract surface from polydata
The vtk.vtkBooleanOperationPolyDataFilter is simply not producing a result which is causing the segfault down your pipeline. @cory.quammen would likely have more details on how to get the vtk.vtkBooleanOperationPolyDataFilter working
Also, note a lot of this code can be simplified with vtki. Perhaps if we figure out how to implement this, it could become a common routine in vtki…
import vtki
import vtk
def get_mesh_intersection(mesh1, mesh2):
'''
Find the intersection volume between mesh vtk objects
and return the volume of intersection
'''
alg = vtk.vtkBooleanOperationPolyDataFilter()
alg.SetInputData(0, mesh2)
alg.SetInputData(1, mesh1)
alg.SetOperationToIntersection()
alg.Update()
intersection = vtki.wrap(alg.GetOutput())
return intersection
# Load meshes
mesh_i = vtki.read('probe0.vtk')
mesh_j = vtki.read('probe1.vtk')
print('VOLUME of probe1:\t', mesh_i.volume, '# of Polys:\t', mesh_i.n_faces)
print('VOLUME of probe2:\t', mesh_j.volume, '# of Polys:\t', mesh_j.n_faces)
intersection = get_mesh_intersection(mesh_i, mesh_j)
print('INTERSECTION VOLUME:\t', intersection.volume)
with output:
VOLUME of probe1: 8.360612727741916 # of Polys: 344
VOLUME of probe2: 8.440116557580266 # of Polys: 480
INTERSECTION VOLUME: 0.0
Also, plotting is easy:
vtki.MultiBlock([mesh_i, mesh_j]).plot(show_edges=1, multi_colors=1)
or
p = vtki.Plotter()
p.add_mesh(mesh_i, color='g', show_edges=1)
p.add_mesh(mesh_j, color='b', show_edges=1)
p.show()