pyvista / pyvista-support

[moved] Please head over to the Discussions tab of the PyVista repository
https://github.com/pyvista/pyvista/discussions
59 stars 4 forks source link

Delaunay_2d on the outermost layer #213

Open YezheWang opened 4 years ago

YezheWang commented 4 years ago

Description

I want to perform a triangulation to turn some discrete points (around the skull) into a connected surface, just on the outer most layer. But now , inside too! I haven't seen any parameters that can control the connection distance of points.

import pyvista as pv
import numpy as np
import SimpleITK as sitk

path = 'cmask21.nii.gz'
itk_img = sitk.ReadImage(path)
data = sitk.GetArrayFromImage(itk_img)

kline = np.where(data==1)
points = np.array(list(map(tuple, zip(*kline))))

cloud = pv.PolyData(points)
cloud.plot(point_size=5)

surf = cloud.delaunay_2d()
surf.plot(show_edges=True)

image

Example Data

cmask21.nii.gz

banesullivan commented 4 years ago

Apologies for such a delay on this. I had initially taken a look and had a tough time making any progress. Coming back to this, I thought it might help to at least share what I have:

surf = cloud.delaunay_3d(alpha=6.0, progress_bar=True).extract_surface().triangulate().clean()

p = pv.Plotter(notebook=False)
p.add_mesh(cloud, color='red')
p.add_mesh(surf, color=True, opacity=0.75)
p.enable_depth_peeling()
p.show()
Screen Shot 2020-08-24 at 7 33 21 PM

I also tried projecting the data to a common plane and doing the delaunay there (because this filter does the triangulation along a best fitting plane) but this wasn't totally successful (note that it might error on vtk<9)

import vtk

def shell_delaunay_2d(mesh):
    copy = mesh.copy()
    copy.points[:,1] = 0.0 # Smash shell down to a plane
    # Exagerate x and z
    copy.points[:,0] *= 10
    copy.points[:,2] *= 10
    tri = copy.delaunay_2d()
    # Put connectivity on plane to the original mesh
    return pv.PolyData(mesh.points, tri.faces)

shell = shell_delaunay_2d(cloud)

p = pv.Plotter(notebook=False)
p.add_mesh(cloud, color='red')
p.add_mesh(shell, color=True, opacity=0.75)
p.add_axes()
p.enable_depth_peeling()
p.show()
Screen Shot 2020-08-24 at 7 38 02 PM