diff --git a/skimage/morphology/selem.py b/skimage/morphology/selem.py index 41be0737..2cccae22 100644 --- a/skimage/morphology/selem.py +++ b/skimage/morphology/selem.py @@ -93,7 +93,7 @@ def disk(radius, dtype=np.uint8): """ Generates a flat, disk-shaped structuring element of a given radius. A pixel is within the neighborhood if the euclidean distance between - it and the origin is no greater than a radius. + it and the origin is no greater than radius. Parameters ---------- @@ -114,3 +114,90 @@ def disk(radius, dtype=np.uint8): s = X**2 s += Y**2 return np.array(s <= radius * radius, dtype=dtype) + + +def cube(width, dtype=np.uint8): + """ + Generates a cube-shaped structuring element (the 3D equivalent of + a square). Every pixel along the perimeter has a chessboard distance + no greater than radius (radius=floor(width/2)) pixels. + + Parameters + ---------- + width : int + The width, height and depth of the cube + + Other Parameters + ---------------- + dtype : data-type + The data type of the structuring element. + + Returns + ------- + selem : ndarray + A structuring element consisting only of ones, i.e. every + pixel belongs to the neighborhood. + + """ + return np.ones((width, width, width), dtype=dtype) + + +def octahedron(radius, dtype=np.uint8): + """ + Generates a octahedron-shaped structuring element of a given radius + (the 3D equivalent of a diamond). A pixel is part of the + neighborhood (i.e. labeled 1) if the city block/manhattan distance + between it and the center of the neighborhood is no greater than + radius. + + Parameters + ---------- + radius : int + The radius of the octahedron-shaped structuring element. + + dtype : data-type + The data type of the structuring element. + + Returns + ------- + + selem : ndarray + The structuring element where elements of the neighborhood + are 1 and 0 otherwise. + """ + # note that in contrast to diamond(), this method allows non-integer radii + n = 2 * radius + 1 + Z, Y, X = np.mgrid[ -radius:radius:n*1j, + -radius:radius:n*1j, + -radius:radius:n*1j] + s = np.abs(X) + np.abs(Y) + np.abs(Z) + return np.array(s <= radius, dtype=dtype) + + +def ball(radius, dtype=np.uint8): + """ + Generates a ball-shaped structuring element of a given radius (the + 3D equivalent of a disk). A pixel is within the neighborhood if the + euclidean distance between it and the origin is no greater than + radius. + + Parameters + ---------- + radius : int + The radius of the ball-shaped structuring element. + + dtype : data-type + The data type of the structuring element. + + Returns + ------- + selem : ndarray + The structuring element where elements of the neighborhood + are 1 and 0 otherwise. + """ + n = 2 * radius + 1 + Z, Y, X = np.mgrid[ -radius:radius:n*1j, + -radius:radius:n*1j, + -radius:radius:n*1j] + s = X**2 + Y**2 + Z**2 + return np.array(s <= radius * radius, dtype=dtype) diff --git a/skimage/morphology/tests/test_selem.py b/skimage/morphology/tests/test_selem.py index eec6cd5b..ccaf43a0 100644 --- a/skimage/morphology/tests/test_selem.py +++ b/skimage/morphology/tests/test_selem.py @@ -36,9 +36,32 @@ class TestSElem(): expected_mask = expected_mask[:, np.newaxis] assert_equal(expected_mask, actual_mask) k = k + 1 - + + def strel_worker_3d(self, fn, func): + matlab_masks = np.load(os.path.join(data_dir, fn)) + k = 0 + for arrname in sorted(matlab_masks): + expected_mask = matlab_masks[arrname] + actual_mask = func(k) + if (expected_mask.shape == (1,)): + expected_mask = expected_mask[:, np.newaxis] + # Test center slice for each dimension. This gives a good + # indication of validity without the need for a 3D reference + # mask. + c = int(expected_mask.shape[0]/2) + assert_equal(expected_mask, actual_mask[c,:,:]) + assert_equal(expected_mask, actual_mask[:,c,:]) + assert_equal(expected_mask, actual_mask[:,:,c]) + k = k + 1 + def test_selem_disk(self): self.strel_worker("disk-matlab-output.npz", selem.disk) def test_selem_diamond(self): self.strel_worker("diamond-matlab-output.npz", selem.diamond) + + def test_selem_ball(self): + self.strel_worker_3d("disk-matlab-output.npz", selem.ball) + + def test_selem_octahedron(self): + self.strel_worker_3d("diamond-matlab-output.npz", selem.octahedron)