1import numpy as np
2import scipy.ndimage as ndi
3
4from skimage import io, draw
5from skimage.data import binary_blobs
6from skimage.util import img_as_ubyte
7from skimage.morphology import skeletonize, skeletonize_3d
8
9from skimage._shared import testing
10from skimage._shared.testing import assert_equal, assert_, parametrize, fetch
11
12# basic behavior tests (mostly copied over from 2D skeletonize)
13
14def test_skeletonize_wrong_dim():
15    im = np.zeros(5, dtype=np.uint8)
16    with testing.raises(ValueError):
17        skeletonize(im, method='lee')
18
19    im = np.zeros((5, 5, 5, 5), dtype=np.uint8)
20    with testing.raises(ValueError):
21        skeletonize(im, method='lee')
22
23
24def test_skeletonize_1D_old_api():
25    # a corner case of an image of a shape(1, N)
26    im = np.ones((5, 1), dtype=np.uint8)
27    res = skeletonize_3d(im)
28    assert_equal(res, im)
29
30
31def test_skeletonize_1D():
32    # a corner case of an image of a shape(1, N)
33    im = np.ones((5, 1), dtype=np.uint8)
34    res = skeletonize(im, method='lee')
35    assert_equal(res, im)
36
37
38def test_skeletonize_no_foreground():
39    im = np.zeros((5, 5), dtype=np.uint8)
40    result = skeletonize(im, method='lee')
41    assert_equal(result, im)
42
43
44def test_skeletonize_all_foreground():
45    im = np.ones((3, 4), dtype=np.uint8)
46    assert_equal(skeletonize(im, method='lee'),
47                 np.array([[0, 0, 0, 0],
48                           [1, 1, 1, 1],
49                           [0, 0, 0, 0]], dtype=np.uint8))
50
51
52def test_skeletonize_single_point():
53    im = np.zeros((5, 5), dtype=np.uint8)
54    im[3, 3] = 1
55    result = skeletonize(im, method='lee')
56    assert_equal(result, im)
57
58
59def test_skeletonize_already_thinned():
60    im = np.zeros((5, 5), dtype=np.uint8)
61    im[3, 1:-1] = 1
62    im[2, -1] = 1
63    im[4, 0] = 1
64    result = skeletonize(im, method='lee')
65    assert_equal(result, im)
66
67
68def test_dtype_conv():
69    # check that the operation does the right thing with floats etc
70    # also check non-contiguous input
71    img = np.random.random((16, 16))[::2, ::2]
72    img[img < 0.5] = 0
73
74    orig = img.copy()
75    res = skeletonize(img, method='lee')
76    img_max = img_as_ubyte(img).max()
77
78    assert_equal(res.dtype, np.uint8)
79    assert_equal(img, orig)  # operation does not clobber the original
80    assert_equal(res.max(), img_max)    # the intensity range is preserved
81
82
83@parametrize("img", [
84    np.ones((8, 8), dtype=float), np.ones((4, 8, 8), dtype=float)
85])
86def test_input_with_warning(img):
87    # check that the input is not clobbered
88    # for 2D and 3D images of varying dtypes
89    check_input(img)
90
91
92@parametrize("img", [
93    np.ones((8, 8), dtype=np.uint8), np.ones((4, 8, 8), dtype=np.uint8),
94    np.ones((8, 8), dtype=bool), np.ones((4, 8, 8), dtype=bool)
95])
96def test_input_without_warning(img):
97    # check that the input is not clobbered
98    # for 2D and 3D images of varying dtypes
99    check_input(img)
100
101
102def check_input(img):
103    orig = img.copy()
104    skeletonize(img, method='lee')
105    assert_equal(img, orig)
106
107
108def test_skeletonize_num_neighbours():
109    # an empty image
110    image = np.zeros((300, 300))
111
112    # foreground object 1
113    image[10:-10, 10:100] = 1
114    image[-100:-10, 10:-10] = 1
115    image[10:-10, -100:-10] = 1
116
117    # foreground object 2
118    rs, cs = draw.line(250, 150, 10, 280)
119    for i in range(10):
120        image[rs + i, cs] = 1
121    rs, cs = draw.line(10, 150, 250, 280)
122    for i in range(20):
123        image[rs + i, cs] = 1
124
125    # foreground object 3
126    ir, ic = np.indices(image.shape)
127    circle1 = (ic - 135)**2 + (ir - 150)**2 < 30**2
128    circle2 = (ic - 135)**2 + (ir - 150)**2 < 20**2
129    image[circle1] = 1
130    image[circle2] = 0
131    result = skeletonize(image, method='lee')
132
133    # there should never be a 2x2 block of foreground pixels in a skeleton
134    mask = np.array([[1,  1],
135                     [1,  1]], np.uint8)
136    blocks = ndi.correlate(result, mask, mode='constant')
137    assert_(not np.any(blocks == 4))
138
139
140def test_two_hole_image():
141    # test a simple 2D image against FIJI
142    img_o = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0],
143                      [0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0],
144                      [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0],
145                      [0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
146                      [0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0],
147                      [0, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0],
148                      [0, 0, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0],
149                      [0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0],
150                      [0, 0, 0, 1, 1, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0],
151                      [0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
152                      [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
153                      [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
154                      [0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0],
155                      [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]],
156                     dtype=np.uint8)
157    img_f = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
158                      [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
159                      [0, 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 0, 0, 0, 0, 0],
160                      [0, 0, 0, 0, 1, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0],
161                      [0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0],
162                      [0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0],
163                      [0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0],
164                      [0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0],
165                      [0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0],
166                      [0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0],
167                      [0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 0, 1, 0, 0, 0],
168                      [0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 0, 0, 0, 0],
169                      [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
170                      [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]],
171                     dtype=np.uint8)
172    res = skeletonize(img_o, method='lee')
173    assert_equal(res, img_f)
174
175
176def test_3d_vs_fiji():
177    # generate an image with blobs and compate its skeleton to
178    # the skeleton generated by FIJI (Plugins>Skeleton->Skeletonize)
179    img = binary_blobs(32, 0.05, n_dim=3, seed=1234)
180    img = img[:-2, ...]
181    img = img.astype(np.uint8)*255
182
183    img_s = skeletonize(img)
184    img_f = io.imread(fetch("data/_blobs_3d_fiji_skeleton.tif"))
185    assert_equal(img_s, img_f)
186