From 9f843fa6195f0efe60e46cc404e1cb79d13075ac Mon Sep 17 00:00:00 2001 From: Almar Klein Date: Wed, 20 Apr 2016 02:26:04 +0200 Subject: [PATCH] marching cubes tests --- skimage/measure/tests/test_marching_cubes.py | 84 +++++++++++++++++++- 1 file changed, 80 insertions(+), 4 deletions(-) diff --git a/skimage/measure/tests/test_marching_cubes.py b/skimage/measure/tests/test_marching_cubes.py index 0e3d77ec..491fedfc 100644 --- a/skimage/measure/tests/test_marching_cubes.py +++ b/skimage/measure/tests/test_marching_cubes.py @@ -2,16 +2,23 @@ import numpy as np from numpy.testing import assert_raises from skimage.draw import ellipsoid, ellipsoid_stats -from skimage.measure import (marching_cubes, mesh_surface_area, - correct_mesh_orientation) +from skimage.measure import (marching_cubes, marching_cubes_lewiner, + mesh_surface_area, correct_mesh_orientation) def test_marching_cubes_isotropic(): ellipsoid_isotropic = ellipsoid(6, 10, 16, levelset=True) _, surf = ellipsoid_stats(6, 10, 16) + + # Classic verts, faces = marching_cubes(ellipsoid_isotropic, 0.) surf_calc = mesh_surface_area(verts, faces) - + # Test within 1% tolerance for isotropic. Will always underestimate. + assert surf > surf_calc and surf_calc > surf * 0.99 + + # Lewiner + verts, faces, *_ = marching_cubes_lewiner(ellipsoid_isotropic, 0.) + surf_calc = mesh_surface_area(verts, faces) # Test within 1% tolerance for isotropic. Will always underestimate. assert surf > surf_calc and surf_calc > surf * 0.99 @@ -21,20 +28,36 @@ def test_marching_cubes_anisotropic(): ellipsoid_anisotropic = ellipsoid(6, 10, 16, spacing=spacing, levelset=True) _, surf = ellipsoid_stats(6, 10, 16) + + # Classic verts, faces = marching_cubes(ellipsoid_anisotropic, 0., spacing=spacing) surf_calc = mesh_surface_area(verts, faces) - + # Test within 1.5% tolerance for anisotropic. Will always underestimate. + assert surf > surf_calc and surf_calc > surf * 0.985 + + # Lewiner + verts, faces, *_ = marching_cubes_lewiner(ellipsoid_anisotropic, 0., + spacing=spacing) + surf_calc = mesh_surface_area(verts, faces) # Test within 1.5% tolerance for anisotropic. Will always underestimate. assert surf > surf_calc and surf_calc > surf * 0.985 def test_invalid_input(): + # Classic assert_raises(ValueError, marching_cubes, np.zeros((2, 2, 1)), 0) assert_raises(ValueError, marching_cubes, np.zeros((2, 2, 1)), 1) assert_raises(ValueError, marching_cubes, np.ones((3, 3, 3)), 1, spacing=(1, 2)) assert_raises(ValueError, marching_cubes, np.zeros((20, 20)), 0) + + # Lewiner + assert_raises(ValueError, marching_cubes_lewiner, np.zeros((2, 2, 1)), 0) + assert_raises(ValueError, marching_cubes_lewiner, np.zeros((2, 2, 1)), 1) + assert_raises(ValueError, marching_cubes_lewiner, np.ones((3, 3, 3)), 1, + spacing=(1, 2)) + assert_raises(ValueError, marching_cubes_lewiner, np.zeros((20, 20)), 0) def test_correct_mesh_orientation(): @@ -74,5 +97,58 @@ def test_correct_mesh_orientation(): np.testing.assert_array_equal(expected, corrected_faces1) +def test_both_algs_same_result_ellipse(): + # Performing this test on data that does not have ambiguities + + sphere_small = ellipsoid(1, 1, 1, levelset=True) + + vertices1, faces1, *_ = marching_cubes(sphere_small, 0) + vertices2, faces2, *_ = marching_cubes_lewiner(sphere_small, 0, allow_degenerate=False) + vertices3, faces3, *_ = marching_cubes_lewiner(sphere_small, 0, allow_degenerate=False, use_classic=True) + + # Order id different, best we can do is test equal shape and same vertices present + assert _same_mesh(vertices1, faces1, vertices2, faces2) + assert _same_mesh(vertices1, faces1, vertices3, faces3) + + +def _same_mesh(vertices1, faces1, vertices2, faces2): + faces1 = [sorted(f) for f in faces1] + faces2 = [sorted(f) for f in faces2] + triangles1 = vertices1[np.array(faces1)] + triangles2 = vertices2[np.array(faces2)] + triang1 = set([tuple(t.flat) for t in triangles1]) + triang2 = set([tuple(t.flat) for t in triangles2]) + return triang1 == triang2 + + +def test_both_algs_same_result_donut(): + # Performing this test on data that does not have ambiguities + + n = 48 + a, b = 2.5/n, -1.25 + isovalue = 0.0 + # + vol = np.empty((n,n,n), 'float32') + for iz in range(vol.shape[0]): + for iy in range(vol.shape[1]): + for ix in range(vol.shape[2]): + z, y, x = float(iz)*a+b, float(iy)*a+b, float(ix)*a+b + vol[iz,iy,ix] = ( ( + (8*x)**2 + (8*y-2)**2 + (8*z)**2 + 16 - 1.85*1.85 ) * ( (8*x)**2 + + (8*y-2)**2 + (8*z)**2 + 16 - 1.85*1.85 ) - 64 * ( (8*x)**2 + (8*y-2)**2 ) + ) * ( ( (8*x)**2 + ((8*y-2)+4)*((8*y-2)+4) + (8*z)**2 + 16 - 1.85*1.85 ) + * ( (8*x)**2 + ((8*y-2)+4)*((8*y-2)+4) + (8*z)**2 + 16 - 1.85*1.85 ) - + 64 * ( ((8*y-2)+4)*((8*y-2)+4) + (8*z)**2 + ) ) + 1025 + + vertices1, faces1, *_ = marching_cubes(vol, 0) + vertices2, faces2, *_ = marching_cubes_lewiner(vol, 0) + vertices3, faces3, *_ = marching_cubes_lewiner(vol, 0, use_classic=True) + + assert not _same_mesh(vertices1, faces1, vertices2, faces2) + #assert _same_mesh(vertices1, faces1, vertices3, faces3) # would have been nice + assert not _same_mesh(vertices2, faces2, vertices3, faces3) + + if __name__ == '__main__': np.testing.run_module_suite()