From 27b5dbca48dc85da2f4cbd658009a840e5685b03 Mon Sep 17 00:00:00 2001 From: Evgeni Burovski Date: Tue, 2 Feb 2016 23:47:09 +0000 Subject: [PATCH] MAINT: skel3d: address review comments Tidy up the python wrapper a little. --- skimage/morphology/_skeletonize_3d.py | 71 +++++++++++-------- .../morphology/tests/test_skeletonize_3d.py | 66 ++++------------- 2 files changed, 55 insertions(+), 82 deletions(-) diff --git a/skimage/morphology/_skeletonize_3d.py b/skimage/morphology/_skeletonize_3d.py index f92590fb..fe0bf86a 100644 --- a/skimage/morphology/_skeletonize_3d.py +++ b/skimage/morphology/_skeletonize_3d.py @@ -1,50 +1,59 @@ from __future__ import division, print_function, absolute_import import numpy as np - from ._skeletonize_3d_cy import _compute_thin_image -def _prepare_image(img_in): - """Convert to a binary image, pad the it w/ zeros, and ensure it's 3D. +def skeletonize_3d(img_in): + """Compute the skeleton of a binary image. + + Thinning is used to reduce each connected component in a binary image + to a single-pixel wide skeleton. + + Parameters + ---------- + image : ndarray, 2D or 3D + A binary image containing the objects to be skeletonized. Zeros + represent background, nonzero values are foreground. + + Returns + ------- + skeleton : ndarray + The thinned image. + + See also + -------- + skeletonize, medial_axis + + References + ---------- + .. [Lee94] Lee et al, Building skeleton models via 3-D medial surface/axis + thinning algorithms. Computer Vision, Graphics, and Image Processing, + 56(6):462–478, 1994. + """ + # 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) img = img_in.copy() if img.ndim == 2: - img = img.reshape((1,) + img.shape) - - intensity = img.max() + 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/ neighborhood of a pixel - img_o = np.zeros(tuple(s + 2 for s in img.shape), - dtype=np.uint8) - img_o[1:-1, 1:-1, 1:-1] = img.astype(np.uint8) - return img_o, intensity + # pad w/ zeros to simplify dealing w/ boundaries + img_o = np.pad(img, pad_width=1, mode='constant') + # do the computation + img_o = np.asarray(_compute_thin_image(img_o)) -def _postprocess_image(img_o, intensity): - """Clip the image (padding is an implementation detail), convert to b/w. - If the original was 2D, convert back to 2D. - """ - img_oo = img_o[1:-1, 1:-1, 1:-1] - img_oo = img_oo.squeeze() - img_oo *= intensity - return img_oo + # clip it back and restore the original intensity range + img_o = img_o[1:-1, 1:-1, 1:-1] + img_o = img_o.squeeze() + img_o *= maxval - -def skeletonize_3d(img_in): - """Compute the thin image. - """ - img, intensity = _prepare_image(img_in) - img = np.asarray(_compute_thin_image(img)) - img = _postprocess_image(img, intensity) - return img - - -if __name__ == "__main__": - pass + return img_o diff --git a/skimage/morphology/tests/test_skeletonize_3d.py b/skimage/morphology/tests/test_skeletonize_3d.py index 985a5bde..a3832dde 100644 --- a/skimage/morphology/tests/test_skeletonize_3d.py +++ b/skimage/morphology/tests/test_skeletonize_3d.py @@ -3,31 +3,30 @@ from __future__ import division, print_function, absolute_import import os import numpy as np -from numpy.testing import assert_equal +from numpy.testing import assert_equal, run_module_suite import skimage from skimage import io from skimage.morphology import compute_thin_image -import matplotlib.pyplot as plt -import matplotlib.ticker as ticker # nose test generators: - # 2D images def test_simple_2d_images(): for fname in ("strip", "loop", "cross", "two-hole"): yield check_skel, fname + # trivial 3D images def test_simple_3d(): for fname in ['3/stack', '4/stack']: yield check_skel_3d, fname + # 'slow' test: Bat Cochlea from FIJI collections. def test_large(): - for fname in ['bat/bat-cochlea-volume',]: + for fname in ['bat/bat-cochlea-volume']: yield check_skel_3d, fname @@ -39,63 +38,28 @@ def get_data_path(): 'data') -def check_skel(fname, viz=False): +def check_skel(fname): # compute the thin image and compare the result to that of ImageJ - - img = np.loadtxt(os.path.join(get_data_path(), fname+'.txt'), dtype=np.uint8) - - if viz: - ax = _viz(img, **dict(marker='s', color='b', s=99, alpha=0.2)) + img = np.loadtxt(os.path.join(get_data_path(), fname + '.txt'), + dtype=np.uint8) # compute img1_2d = compute_thin_image(img) - if viz: - ax = _viz(img1_2d, ax, **dict(marker='o', color='r', - s=80, alpha=0.7, label='us')) + # and compare to FIJI + img_f = np.loadtxt(os.path.join(get_data_path(), fname + '_fiji.txt'), + dtype=np.uint8) - # compare to FIJI - img_f = np.loadtxt(os.path.join(get_data_path(), fname+'_fiji.txt'), dtype=np.uint8) - - if not viz: - # actually compare images - assert_equal(img1_2d, img_f) - else: - ax = _viz(img_f, ax, **dict(marker='o', color='g', s=45, label='fiji')) - - ax.legend() - ax.grid(True) - - def yformatter(val, pos): - return int(img.shape[1] - val + 1) - def xformatter(val, pos): - return int(val + 1) - ax.xaxis.set_major_formatter(ticker.FuncFormatter(xformatter)) - ax.yaxis.set_major_formatter(ticker.FuncFormatter(yformatter)) - - plt.show() - - -def _viz(img, ax=None, **kwds): - if ax is None: - import matplotlib.pyplot as plt - fix, ax = plt.subplots() - - x, y = np.nonzero(img) - ax.scatter(y, img.shape[1] - x, **kwds) - return ax + assert_equal(img1_2d, img_f) 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 = 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) assert_equal(img_s, img_f) -if __name__ == "__main__": - import sys - if len(sys.argv) < 2: - sys.exit("Expect an image name from the data/ directory.") - check_skel(sys.argv[1], True) +if __name__ == '__main__': + run_module_suite()