diff --git a/skimage/morphology/__init__.py b/skimage/morphology/__init__.py index 4579794f..d1d35ee5 100644 --- a/skimage/morphology/__init__.py +++ b/skimage/morphology/__init__.py @@ -6,7 +6,7 @@ from .selem import (square, rectangle, diamond, disk, cube, octahedron, ball, octagon, star) from .watershed import watershed from ._skeletonize import skeletonize, medial_axis -from ._skeletonize_3d import skeletonize_3d as compute_thin_image +from ._skeletonize_3d import skeletonize_3d from .convex_hull import convex_hull_image, convex_hull_object from .greyreconstruct import reconstruction from .misc import remove_small_objects, remove_small_holes diff --git a/skimage/morphology/_skeletonize_3d.py b/skimage/morphology/_skeletonize_3d.py index cd63c372..8e2460da 100644 --- a/skimage/morphology/_skeletonize_3d.py +++ b/skimage/morphology/_skeletonize_3d.py @@ -1,10 +1,11 @@ from __future__ import division, print_function, absolute_import import numpy as np +from ..util import img_as_ubyte from ._skeletonize_3d_cy import _compute_thin_image -def skeletonize_3d(img_in): +def skeletonize_3d(img): """Compute the skeleton of a binary image. Thinning is used to reduce each connected component in a binary image @@ -12,7 +13,7 @@ def skeletonize_3d(img_in): Parameters ---------- - image : ndarray, 2D or 3D + img : ndarray, 2D or 3D A binary image containing the objects to be skeletonized. Zeros represent background, nonzero values are foreground. @@ -33,17 +34,19 @@ def skeletonize_3d(img_in): """ # make sure the image is 3D or 2D (if it is, temporarily upcast to 3D) - if img_in.ndim < 2 or img_in.ndim > 3: - raise ValueError('expect 2D, got ndim = %s' % img_in.ndim) + if img.ndim < 2 or img.ndim > 3: + raise ValueError('expect 2D, got ndim = %s' % img.ndim) - img = img_in.copy() + img = img_as_ubyte(img) + img = np.ascontiguousarray(img) + + img = img.copy() if img.ndim == 2: img = img[None, ...] # normalize to binary maxval = img.max() img[img != 0] = 1 - img = img.astype(np.uint8) # pad w/ zeros to simplify dealing w/ boundaries img_o = np.pad(img, pad_width=1, mode='constant') diff --git a/skimage/morphology/tests/test_skeletonize_3d.py b/skimage/morphology/tests/test_skeletonize_3d.py index a3832dde..a3f90b9c 100644 --- a/skimage/morphology/tests/test_skeletonize_3d.py +++ b/skimage/morphology/tests/test_skeletonize_3d.py @@ -1,14 +1,111 @@ from __future__ import division, print_function, absolute_import import os +import warnings import numpy as np -from numpy.testing import assert_equal, run_module_suite +from numpy.testing import (assert_equal, run_module_suite, assert_raises, + assert_) + +import scipy.ndimage as ndi import skimage -from skimage import io +from skimage import io, draw, data_dir +#from skimage import draw +from skimage.util import img_as_ubyte -from skimage.morphology import compute_thin_image +from skimage.morphology import skeletonize_3d + + +# basic behavior tests (mostly copied over from 2D skeletonize) + +def test_skeletonize_wrong_dim(): + im = np.zeros(5, dtype=np.uint8) + assert_raises(ValueError, skeletonize_3d, im) + + im = np.zeros((5, 5, 5, 5), dtype=np.uint8) + assert_raises(ValueError, skeletonize_3d, im) + + +def test_skeletonize_no_foreground(): + im = np.zeros((5, 5), dtype=np.uint8) + result = skeletonize_3d(im) + assert_equal(result, im) + + +def test_skeletonize_all_foreground(): + im = np.ones((3, 4), dtype=np.uint8) + assert_equal(skeletonize_3d(im), + np.array([[0, 0, 0, 0], + [1, 1, 1, 1], + [0, 0, 0, 0]], dtype=np.uint8)) + + +def test_skeletonize_single_point(): + im = np.zeros((5, 5), dtype=np.uint8) + im[3, 3] = 1 + result = skeletonize_3d(im) + assert_equal(result, im) + + +def test_skeletonize_already_thinned(): + im = np.zeros((5, 5), dtype=np.uint8) + im[3, 1:-1] = 1 + im[2, -1] = 1 + im[4, 0] = 1 + result = skeletonize_3d(im) + assert_equal(result, im) + + +def test_dtype_conv(): + # check that the operation does the right thing with floats etc + # also check non-contiguous input + img = np.random.random((16, 16))[::2, ::2] + img[img < 0.5] = 0 + + orig = img.copy() + + with warnings.catch_warnings(): + # UserWarning for possible precision loss, expected + warnings.simplefilter('ignore', UserWarning) + res = skeletonize_3d(img) + + assert_equal(res.dtype, np.uint8) + assert_equal(img, orig) # operation does not clobber the original + assert_equal(res.max(), + img_as_ubyte(img).max()) # the intensity range is preserved + + +def test_skeletonize_num_neighbours(): + # an empty image + image = np.zeros((300, 300)) + + # foreground object 1 + image[10:-10, 10:100] = 1 + image[-100:-10, 10:-10] = 1 + image[10:-10, -100:-10] = 1 + + # foreground object 2 + rs, cs = draw.line(250, 150, 10, 280) + for i in range(10): + image[rs + i, cs] = 1 + rs, cs = draw.line(10, 150, 250, 280) + for i in range(20): + image[rs + i, cs] = 1 + + # foreground object 3 + ir, ic = np.indices(image.shape) + circle1 = (ic - 135)**2 + (ir - 150)**2 < 30**2 + circle2 = (ic - 135)**2 + (ir - 150)**2 < 20**2 + image[circle1] = 1 + image[circle2] = 0 + result = skeletonize_3d(image) + + # there should never be a 2x2 block of foreground pixels in a skeleton + mask = np.array([[1, 1], + [1, 1]], np.uint8) + blocks = ndi.correlate(result, mask, mode='constant') + assert_(not np.any(blocks == 4)) # nose test generators: @@ -44,7 +141,7 @@ def check_skel(fname): dtype=np.uint8) # compute - img1_2d = compute_thin_image(img) + img1_2d = skeletonize_3d(img) # and compare to FIJI img_f = np.loadtxt(os.path.join(get_data_path(), fname + '_fiji.txt'), @@ -57,7 +154,7 @@ def check_skel_3d(fname): img = io.imread(os.path.join(get_data_path(), fname + '.tif')) img_f = io.imread(os.path.join(get_data_path(), fname + '_fiji.tif')) - img_s = compute_thin_image(img) + img_s = skeletonize_3d(img) assert_equal(img_s, img_f)