我们从Python开源项目中，提取了以下12个代码示例，用于说明如何使用skimage.measure.marching_cubes()。
def plot_3D(img, threshold=-400): verts, faces = measure.marching_cubes(img, threshold) fig = plt.figure(figsize=(10, 10)) ax = fig.add_subplot(111, projection='3d') mesh = Poly3DCollection(verts[faces], alpha=0.1) face_color = [0.5, 0.5, 1] mesh.set_facecolor(face_color) ax.add_collection3d(mesh) ax.set_xlim(0, img.shape[0]) ax.set_ylim(0, img.shape[1]) ax.set_zlim(0, img.shape[2]) plt.show()
def generate_mesh(image, isovalue=0, channel=0): """ Creates and returns a Mesh object :param image: an AICSImage object :param isovalue: The value that is used to pick the isosurface returned by the marching cubes algorithm For more info: https://www.youtube.com/watch?v=5fNbCFjqWao @ 40:00 mins :param channel: The channel in the image that is used to extract the isosurface :return: A Mesh object """ if not isinstance(image, AICSImage): raise ValueError("Meshes can only be generated with AICSImage objects!") if channel >= image.size_c: raise IndexError("Channel provided for mesh generation is out of bounds for image data!") image_stack = image.get_image_data("ZYX", C=channel) # Use marching cubes to obtain the surface mesh of the membrane wall verts, faces, normals, values = measure.marching_cubes(image_stack, isovalue, allow_degenerate=False) return Mesh(verts, faces, normals, values)
def marching_cubes(field,iso=0.5): try: from skimage.measure import marching_cubes surface_points, surface_triangles = marching_cubes(density_field,iso) except ImportError: print "Please try to install SciKit-Image!" from mayavi import mlab from mayavi.mlab import contour3d mlab.clf() surface = mlab.contour3d(field,contours=[iso]) my_actor=surface.actor.actors[0] poly_data_object=my_actor.mapper.input surface_points = (np.array(poly_data_object.points) - np.array([abs(grid_points/2.),abs(grid_points/2.),abs(grid_points/2.)])[np.newaxis,:])*(grid_max/abs(grid_points/2.)) surface_triangles = poly_data_object.polys.data.to_array().reshape([-1,4]) surface_triangles = surface_triangles[:,1:] return surface_points, surface_triangles
def plot_3d_cubic(image): ''' plot the 3D cubic :param image: image saved as npy file path :return: ''' from skimage import measure, morphology from mpl_toolkits.mplot3d.art3d import Poly3DCollection image = np.load(image) verts, faces = measure.marching_cubes(image,0) fig = plt.figure(figsize=(40, 40)) ax = fig.add_subplot(111, projection='3d') # Fancy indexing: `verts[faces]` to generate a collection of triangles mesh = Poly3DCollection(verts[faces], alpha=0.1) face_color = [0.5, 0.5, 1] mesh.set_facecolor(face_color) ax.add_collection3d(mesh) ax.set_xlim(0, image.shape[0]) ax.set_ylim(0, image.shape[1]) ax.set_zlim(0, image.shape[2]) plt.show() # LUNA2016 data prepare ,first step: truncate HU to -1000 to 400
def plot_3d(image, threshold=-300): # Position the scan upright, # so the head of the patient would be at the top facing the camera p = image.transpose(2,1,0) verts, faces = measure.marching_cubes(p, threshold) fig = plt.figure(figsize=(10, 10)) ax = fig.add_subplot(111, projection='3d') # Fancy indexing: `verts[faces]` to generate a collection of triangles mesh = Poly3DCollection(verts[faces], alpha=0.70) face_color = [0.45, 0.45, 0.75] mesh.set_facecolor(face_color) ax.add_collection3d(mesh) ax.set_xlim(0, p.shape[0]) ax.set_ylim(0, p.shape[1]) ax.set_zlim(0, p.shape[2]) plt.show()
def three_d_print(self): """This will produce a 3d printable stl based on self.volume_data. It is to be used for the final "print" button, and needs to be fed high quality data.""" name = raw_input('What should the filename be?') + '.stl' verts, faces = measure.marching_cubes(self.volume_data, 0) #Marching Cubes algorithm solid = mesh.Mesh(np.zeros(faces.shape[0], dtype=mesh.Mesh.dtype)) for i, f in enumerate(faces): for j in range(3): solid.vectors[i][j] = verts[f[j],:] solid.save(name) #Here be UI
def plot_3d(image, threshold=-300): # Position the scan upright, # so the head of the patient would be at the top facing the camera p = image.transpose(2,1,0) #p = image verts, faces = measure.marching_cubes(p, threshold) fig = plt.figure(figsize=(10, 10)) ax = fig.add_subplot(111, projection='3d') # Fancy indexing: `verts[faces]` to generate a collection of triangles mesh = Poly3DCollection(verts[faces], alpha=0.70) face_color = [0.45, 0.45, 0.75] mesh.set_facecolor(face_color) ax.add_collection3d(mesh) ax.set_xlim(0, p.shape[0]) ax.set_ylim(0, p.shape[1]) ax.set_zlim(0, p.shape[2]) plt.show()
def estimate_surface_area(self): """ Estimate the surface area by summing the areas of a trianglation of the nodules surface in 3d. Returned units are mm^2. """ mask = self.get_boolean_mask() mask = np.pad(mask, [(1,1), (1,1), (1,1)], 'constant') # Cap the ends. dist = dtrans(mask) - dtrans(~mask) rxy = self.scan.pixel_spacing rz = self.scan.slice_thickness verts, faces, _, _ = marching_cubes(dist, 0, spacing=(rxy, rxy, rz)) return mesh_surface_area(verts, faces)
def getVFByMarchingCubes(voxels, threshold=0.5): v, f = sk.marching_cubes(voxels, level=threshold) return v, f
def implicit_surface(density_field,size,resolution,iso=0.5): import numpy as np from scipy.cluster.vq import kmeans, vq from openalea.container import array_dict from skimage.measure import marching_cubes surface_points, surface_triangles = marching_cubes(density_field,iso) surface_points = (np.array(surface_points))*(size*resolution/np.array(density_field.shape)) - size*resolution/2. points_ids = np.arange(len(surface_points)) points_to_delete = [] for p,point in enumerate(surface_points): matching_points = np.sort(np.where(vq(surface_points,np.array([point]))[1] == 0)[0]) if len(matching_points) > 1: points_to_fuse = matching_points[1:] for m_p in points_to_fuse: surface_triangles[np.where(surface_triangles==m_p)] = matching_points[0] points_to_delete.append(m_p) points_to_delete = np.unique(points_to_delete) print len(points_to_delete),"points deleted" surface_points = np.delete(surface_points,points_to_delete,0) points_ids = np.delete(points_ids,points_to_delete,0) surface_triangles = array_dict(np.arange(len(surface_points)),points_ids).values(surface_triangles) for p,point in enumerate(surface_points): matching_points = np.where(vq(surface_points,np.array([point]))[1] == 0)[0] if len(matching_points) > 1: print p,point raw_input() triangles_to_delete = [] for t,triangle in enumerate(surface_triangles): if len(np.unique(triangle)) < 3: triangles_to_delete.append(t) # elif triangle.max() >= len(surface_points): # triangles_to_delete.append(t) surface_triangles = np.delete(surface_triangles,triangles_to_delete,0) return surface_points, surface_triangles
