wassname/scikit-image
1
0
Fork 0
mirror of https://github.com/wassname/scikit-image.git synced 2026-07-11 21:25:42 +08:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #1803 from JDWarner/mesh_orientation_cleanup

Browse Source 
Guarantee correct mesh orientation from marching cubes
This commit is contained in:
Stefan van der Walt
2015-12-07 20:13:58 -08:00
parent 4eb4dcd52b 443f4cc0eb
commit c133d5fd59
2 changed files with 117 additions and 17 deletions
Download Patch File Download Diff File Expand all files Collapse all files
skimage/measure/_marching_cubes.py
+99 -16
View File
@@ -1,8 +1,10 @@
import numpy as np
import scipy.ndimage as ndi
from . import _marching_cubes_cy
def marching_cubes(volume, level, spacing=(1., 1., 1.)):
def marching_cubes(volume, level, spacing=(1., 1., 1.),
gradient_direction='descent'):
"""
Marching cubes algorithm to find iso-valued surfaces in 3d volumetric data
@@ -15,6 +17,12 @@ def marching_cubes(volume, level, spacing=(1., 1., 1.)):
spacing : length-3 tuple of floats
Voxel spacing in spatial dimensions corresponding to numpy array
indexing dimensions (M, N, P) as in `volume`.
gradient_direction : string
Controls if the mesh was generated from an isosurface with gradient
descent toward objects of interest (the default), or the opposite.
The two options are:
* descent : Object was greater than exterior
* ascent : Exterior was greater than object
Returns
-------
@@ -68,14 +76,10 @@ def marching_cubes(volume, level, spacing=(1., 1., 1.)):
lexicographical order) coordinate in the contour. This is a side-effect
of how the input array is traversed, but can be relied upon.
The generated mesh does not guarantee coherent orientation because of how
symmetry is used in the algorithm. If this is required, e.g. due to a
particular visualization package or for generating 3D printing STL files,
the utility ``skimage.measure.correct_mesh_orientation`` is available to
fix this in post-processing.
The generated mesh guarantees coherent orientation as of version 0.12.
To quantify the area of an isosurface generated by this algorithm, pass
the outputs directly into `skimage.measure.mesh_surface_area`.
outputs directly into `skimage.measure.mesh_surface_area`.
Regarding visualization of algorithm output, the ``mayavi`` package
is recommended. To contour a volume named `myvolume` about the level 0.0::
@@ -122,8 +126,18 @@ def marching_cubes(volume, level, spacing=(1., 1., 1.)):
# Returns a true mesh with no degenerate faces.
verts, faces = _marching_cubes_cy.unpack_unique_verts(raw_faces)
verts = np.asarray(verts)
faces = np.asarray(faces)
# Calculate gradient of `volume`, then interpolate to vertices in `verts`
grad_x, grad_y, grad_z = np.gradient(volume)
# Fancy indexing to define two vector arrays from triangle vertices
faces = _correct_mesh_orientation(volume, verts[faces], faces, spacing,
gradient_direction)
# Adjust for non-isotropic spacing in `verts` at time of return
return np.asarray(verts) * np.r_[spacing], np.asarray(faces)
return verts * np.r_[spacing], faces
def mesh_surface_area(verts, faces):
@@ -187,7 +201,7 @@ def correct_mesh_orientation(volume, verts, faces, spacing=(1., 1., 1.),
indexing dimensions (M, N, P) as in `volume`.
gradient_direction : string
Controls if the mesh was generated from an isosurface with gradient
ascent toward objects of interest (the default), or the opposite.
descent toward objects of interest (the default), or the opposite.
The two options are:
* descent : Object was greater than exterior
* ascent : Exterior was greater than object
@@ -225,17 +239,86 @@ def correct_mesh_orientation(volume, verts, faces, spacing=(1., 1., 1.),
skimage.measure.mesh_surface_area
"""
import scipy.ndimage as ndi
import warnings
warnings.warn(
DeprecationWarning("`correct_mesh_orientation` is deprecated for "
"removal as `marching_cubes` now guarantess "
"correct mesh orientation."))
# Calculate gradient of `volume`, then interpolate to vertices in `verts`
grad_x, grad_y, grad_z = np.gradient(volume)
verts = verts.copy()
verts[:, 0] /= spacing[0]
verts[:, 1] /= spacing[1]
verts[:, 2] /= spacing[2]
# Fancy indexing to define two vector arrays from triangle vertices
actual_verts = verts[faces]
actual_verts[:, 0] /= spacing[0]
actual_verts[:, 1] /= spacing[1]
actual_verts[:, 2] /= spacing[2]
return _correct_mesh_orientation(volume, actual_verts, faces, spacing,
gradient_direction)
def _correct_mesh_orientation(volume, actual_verts, faces,
spacing=(1., 1., 1.),
gradient_direction='descent'):
"""
Correct orientations of mesh faces.
Parameters
----------
volume : (M, N, P) array of doubles
Input data volume to find isosurfaces. Will be cast to `np.float64`.
actual_verts : (F, 3, 3) array of floats
Array with (face, vertex, coords) index coordinates.
faces : (F, 3) array of ints
List of length-3 lists of integers, referencing vertex coordinates as
provided in `verts`.
spacing : length-3 tuple of floats
Voxel spacing in spatial dimensions corresponding to numpy array
indexing dimensions (M, N, P) as in `volume`.
gradient_direction : string
Controls if the mesh was generated from an isosurface with gradient
descent toward objects of interest (the default), or the opposite.
The two options are:
* descent : Object was greater than exterior
* ascent : Exterior was greater than object
Returns
-------
faces_corrected (F, 3) array of ints
Corrected list of faces referencing vertex coordinates in `verts`.
Notes
-----
Certain applications and mesh processing algorithms require all faces
to be oriented in a consistent way. Generally, this means a normal vector
points "out" of the meshed shapes. This algorithm corrects the output from
`skimage.measure.marching_cubes` by flipping the orientation of
mis-oriented faces.
Because marching cubes could be used to find isosurfaces either on
gradient descent (where the desired object has greater values than the
exterior) or ascent (where the desired object has lower values than the
exterior), the ``gradient_direction`` kwarg allows the user to inform this
algorithm which is correct. If the resulting mesh appears to be oriented
completely incorrectly, try changing this option.
The arguments expected by this function are the exact outputs from
`skimage.measure.marching_cubes` except `actual_verts`, which is an
uncorrected version of the fancy indexing operation `verts[faces]`.
Only `faces` is corrected and returned as the vertices do not change,
only the order in which they are referenced.
This algorithm assumes ``faces`` provided are exclusively triangles.
See Also
--------
skimage.measure.marching_cubes
skimage.measure.mesh_surface_area
"""
# Calculate gradient of `volume`, then interpolate to vertices in `verts`
grad_x, grad_y, grad_z = np.gradient(volume)
a = actual_verts[:, 0, :] - actual_verts[:, 1, :]
b = actual_verts[:, 0, :] - actual_verts[:, 2, :]
skimage/measure/tests/test_marching_cubes.py
+18 -1
View File
@@ -39,7 +39,24 @@ def test_invalid_input():
def test_correct_mesh_orientation():
sphere_small = ellipsoid(1, 1, 1, levelset=True)
verts, faces = marching_cubes(sphere_small, 0.)
# Mesh with incorrectly oriented faces which was previously returned from
# `marching_cubes`, before it guaranteed correct mesh orientation
verts = np.array([[1., 2., 2.],
[2., 2., 1.],
[2., 1., 2.],
[2., 2., 3.],
[2., 3., 2.],
[3., 2., 2.]])
faces = np.array([[0, 1, 2],
[2, 0, 3],
[1, 0, 4],
[4, 0, 3],
[1, 2, 5],
[2, 3, 5],
[1, 4, 5],
[5, 4, 3]])
# Correct mesh orientation - descent
corrected_faces1 = correct_mesh_orientation(sphere_small, verts, faces,