diff --git a/skimage/feature/corner.py b/skimage/feature/corner.py index 6081e1a2..fef9b046 100644 --- a/skimage/feature/corner.py +++ b/skimage/feature/corner.py @@ -610,7 +610,7 @@ def corner_fast(image, n=12, threshold=0.15): [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]) - >>> corner_peaks(corner_fast(square, 9), min_distance=1) + >>> corner_peaks(corner_fast(square, 9), min_distance=1, threshold_rel=0.1) array([[3, 3], [3, 8], [8, 3], @@ -799,7 +799,7 @@ def corner_subpix(image, corners, window_size=11, alpha=0.99): return corners_subpix -def corner_peaks(image, min_distance=10, threshold_abs=0, threshold_rel=0.1, +def corner_peaks(image, min_distance=1, threshold_abs=None, threshold_rel=None, exclude_border=True, indices=True, num_peaks=np.inf, footprint=None, labels=None): """Find corners in corner measure response image. @@ -823,18 +823,13 @@ def corner_peaks(image, min_distance=10, threshold_abs=0, threshold_rel=0.1, [ 0., 0., 1., 1., 0.], [ 0., 0., 1., 1., 0.], [ 0., 0., 0., 0., 0.]]) - >>> peak_local_max(response, exclude_border=False) + >>> peak_local_max(response) array([[2, 2], [2, 3], [3, 2], [3, 3]]) - >>> corner_peaks(response, exclude_border=False) + >>> corner_peaks(response) array([[2, 2]]) - >>> corner_peaks(response, exclude_border=False, min_distance=0) - array([[2, 2], - [2, 3], - [3, 2], - [3, 3]]) """ diff --git a/skimage/feature/peak.py b/skimage/feature/peak.py index 421abeec..0cc1b26d 100644 --- a/skimage/feature/peak.py +++ b/skimage/feature/peak.py @@ -3,46 +3,49 @@ import scipy.ndimage as ndi from ..filters import rank_order -def peak_local_max(image, min_distance=10, threshold_abs=0, threshold_rel=0.1, - exclude_border=True, indices=True, num_peaks=np.inf, - footprint=None, labels=None): - """ - Find peaks in an image, and return them as coordinates or a boolean array. +def peak_local_max(image, min_distance=1, threshold_abs=None, + threshold_rel=None, exclude_border=True, indices=True, + num_peaks=np.inf, footprint=None, labels=None): + """Find peaks in an image as coordinate list or boolean mask. Peaks are the local maxima in a region of `2 * min_distance + 1` (i.e. peaks are separated by at least `min_distance`). - NOTE: If peaks are flat (i.e. multiple adjacent pixels have identical + If peaks are flat (i.e. multiple adjacent pixels have identical intensities), the coordinates of all such pixels are returned. + If both `threshold_abs` and `threshold_rel` are provided, the maximum + of the two is chosen as the minimum intensity threshold of peaks. + Parameters ---------- - image : ndarray of floats + image : ndarray Input image. - min_distance : int + min_distance : int, optional Minimum number of pixels separating peaks in a region of `2 * min_distance + 1` (i.e. peaks are separated by at least `min_distance`). If `exclude_border` is True, this value also excludes a border `min_distance` from the image boundary. To find the maximum number of peaks, use `min_distance=1`. - threshold_abs : float - Minimum intensity of peaks. - threshold_rel : float - Minimum intensity of peaks calculated as `max(image) * threshold_rel`. - exclude_border : bool + threshold_abs : float, optional + Minimum intensity of peaks. By default, the absolute threshold is + the minimum intensity of the image. + threshold_rel : float, optional + Minimum intensity of peaks, calculated as `max(image) * threshold_rel`. + exclude_border : bool, optional If True, `min_distance` excludes peaks from the border of the image as well as from each other. - indices : bool - If True, the output will be an array representing peak coordinates. - If False, the output will be a boolean array shaped as `image.shape` - with peaks present at True elements. - num_peaks : int + indices : bool, optional + If True, the output will be an array representing peak + coordinates. If False, the output will be a boolean array shaped as + `image.shape` with peaks present at True elements. + num_peaks : int, optional Maximum number of peaks. When the number of peaks exceeds `num_peaks`, return `num_peaks` peaks based on highest peak intensity. footprint : ndarray of bools, optional If provided, `footprint == 1` represents the local region within which to search for peaks at every point in `image`. Overrides - `min_distance`, except for border exclusion if `exclude_border=True`. + `min_distance` (also for `exclude_border`). labels : ndarray of ints, optional If provided, each unique region `labels == value` represents a unique region to search for peaks. Zero is reserved for background. @@ -58,10 +61,10 @@ def peak_local_max(image, min_distance=10, threshold_abs=0, threshold_rel=0.1, Notes ----- The peak local maximum function returns the coordinates of local peaks - (maxima) in a image. A maximum filter is used for finding local maxima. - This operation dilates the original image. After comparison between - dilated and original image, peak_local_max function returns the - coordinates of peaks where dilated image = original. + (maxima) in an image. A maximum filter is used for finding local maxima. + This operation dilates the original image. After comparison of the dilated + and original image, this function returns the coordinates or a mask of the + peaks where the dilated image equals the original image. Examples -------- @@ -90,7 +93,9 @@ def peak_local_max(image, min_distance=10, threshold_abs=0, threshold_rel=0.1, array([[10, 10, 10]]) """ + out = np.zeros_like(image, dtype=np.bool) + # In the case of labels, recursively build and return an output # operating on each label separately if labels is not None: @@ -123,7 +128,6 @@ def peak_local_max(image, min_distance=10, threshold_abs=0, threshold_rel=0.1, else: return out - image = image.copy() # Non maximum filter if footprint is not None: image_max = ndi.maximum_filter(image, footprint=footprint, @@ -131,25 +135,33 @@ def peak_local_max(image, min_distance=10, threshold_abs=0, threshold_rel=0.1, else: size = 2 * min_distance + 1 image_max = ndi.maximum_filter(image, size=size, mode='constant') - mask = (image == image_max) - image *= mask + mask = image == image_max - if exclude_border: + if exclude_border and (footprint is not None or min_distance > 0): # zero out the image borders - for i in range(image.ndim): - image = image.swapaxes(0, i) - image[:min_distance] = 0 - image[-min_distance:] = 0 - image = image.swapaxes(0, i) + for i in range(mask.ndim): + mask = mask.swapaxes(0, i) + remove = (footprint.shape[i] if footprint is not None + else 2 * min_distance) + mask[:remove // 2] = mask[-remove // 2:] = False + mask = mask.swapaxes(0, i) # find top peak candidates above a threshold - peak_threshold = max(np.max(image.ravel()) * threshold_rel, threshold_abs) + thresholds = [] + if threshold_abs is None: + threshold_abs = image.min() + thresholds.append(threshold_abs) + if threshold_rel is not None: + thresholds.append(threshold_rel * image.max()) + if thresholds: + mask &= image > max(thresholds) # get coordinates of peaks - coordinates = np.argwhere(image > peak_threshold) + coordinates = np.transpose(mask.nonzero()) if coordinates.shape[0] > num_peaks: - intensities = image.flat[np.ravel_multi_index(coordinates.transpose(),image.shape)] + intensities = image.flat[np.ravel_multi_index(coordinates.transpose(), + image.shape)] idx_maxsort = np.argsort(intensities)[::-1] coordinates = coordinates[idx_maxsort][:num_peaks] diff --git a/skimage/feature/tests/test_brief.py b/skimage/feature/tests/test_brief.py index 610dc984..7f6ae9b7 100644 --- a/skimage/feature/tests/test_brief.py +++ b/skimage/feature/tests/test_brief.py @@ -18,7 +18,8 @@ def test_normal_mode(): """Verify the computed BRIEF descriptors with expected for normal mode.""" img = data.coins() - keypoints = corner_peaks(corner_harris(img), min_distance=5) + keypoints = corner_peaks(corner_harris(img), min_distance=5, + threshold_abs=0, threshold_rel=0.1) extractor = BRIEF(descriptor_size=8, sigma=2) @@ -40,7 +41,8 @@ def test_uniform_mode(): """Verify the computed BRIEF descriptors with expected for uniform mode.""" img = data.coins() - keypoints = corner_peaks(corner_harris(img), min_distance=5) + keypoints = corner_peaks(corner_harris(img), min_distance=5, + threshold_abs=0, threshold_rel=0.1) extractor = BRIEF(descriptor_size=8, sigma=2, mode='uniform') diff --git a/skimage/feature/tests/test_corner.py b/skimage/feature/tests/test_corner.py index 953e47e3..c3919fe8 100644 --- a/skimage/feature/tests/test_corner.py +++ b/skimage/feature/tests/test_corner.py @@ -107,21 +107,25 @@ def test_square_image(): im[:25, :25] = 1. # Moravec - results = peak_local_max(corner_moravec(im)) + results = peak_local_max(corner_moravec(im), + min_distance=10, threshold_rel=0) # interest points along edge assert len(results) == 57 # Harris - results = peak_local_max(corner_harris(im, method='k')) + results = peak_local_max(corner_harris(im, method='k'), + min_distance=10, threshold_rel=0) # interest at corner assert len(results) == 1 - results = peak_local_max(corner_harris(im, method='eps')) + results = peak_local_max(corner_harris(im, method='eps'), + min_distance=10, threshold_rel=0) # interest at corner assert len(results) == 1 # Shi-Tomasi - results = peak_local_max(corner_shi_tomasi(im)) + results = peak_local_max(corner_shi_tomasi(im), + min_distance=10, threshold_rel=0) # interest at corner assert len(results) == 1 @@ -133,18 +137,22 @@ def test_noisy_square_image(): im = im + np.random.uniform(size=im.shape) * .2 # Moravec - results = peak_local_max(corner_moravec(im)) + results = peak_local_max(corner_moravec(im), + min_distance=10, threshold_rel=0) # undefined number of interest points assert results.any() # Harris - results = peak_local_max(corner_harris(im, sigma=1.5, method='k')) + results = peak_local_max(corner_harris(im, method='k'), + min_distance=10, threshold_rel=0) assert len(results) == 1 - results = peak_local_max(corner_harris(im, sigma=1.5, method='eps')) + results = peak_local_max(corner_harris(im, method='eps'), + min_distance=10, threshold_rel=0) assert len(results) == 1 # Shi-Tomasi - results = peak_local_max(corner_shi_tomasi(im, sigma=1.5)) + results = peak_local_max(corner_shi_tomasi(im, sigma=1.5), + min_distance=10, threshold_rel=0) assert len(results) == 1 @@ -156,11 +164,13 @@ def test_squared_dot(): # Moravec fails # Harris - results = peak_local_max(corner_harris(im)) + results = peak_local_max(corner_harris(im), + min_distance=10, threshold_rel=0) assert (results == np.array([[6, 6]])).all() # Shi-Tomasi - results = peak_local_max(corner_shi_tomasi(im)) + results = peak_local_max(corner_shi_tomasi(im), + min_distance=10, threshold_rel=0) assert (results == np.array([[6, 6]])).all() @@ -173,20 +183,26 @@ def test_rotated_img(): im_rotated = im.T # Moravec - results = peak_local_max(corner_moravec(im)) - results_rotated = peak_local_max(corner_moravec(im_rotated)) + results = peak_local_max(corner_moravec(im), + min_distance=10, threshold_rel=0) + results_rotated = peak_local_max(corner_moravec(im_rotated), + min_distance=10, threshold_rel=0) assert (np.sort(results[:, 0]) == np.sort(results_rotated[:, 1])).all() assert (np.sort(results[:, 1]) == np.sort(results_rotated[:, 0])).all() # Harris - results = peak_local_max(corner_harris(im)) - results_rotated = peak_local_max(corner_harris(im_rotated)) + results = peak_local_max(corner_harris(im), + min_distance=10, threshold_rel=0) + results_rotated = peak_local_max(corner_harris(im_rotated), + min_distance=10, threshold_rel=0) assert (np.sort(results[:, 0]) == np.sort(results_rotated[:, 1])).all() assert (np.sort(results[:, 1]) == np.sort(results_rotated[:, 0])).all() # Shi-Tomasi - results = peak_local_max(corner_shi_tomasi(im)) - results_rotated = peak_local_max(corner_shi_tomasi(im_rotated)) + results = peak_local_max(corner_shi_tomasi(im), + min_distance=10, threshold_rel=0) + results_rotated = peak_local_max(corner_shi_tomasi(im_rotated), + min_distance=10, threshold_rel=0) assert (np.sort(results[:, 0]) == np.sort(results_rotated[:, 1])).all() assert (np.sort(results[:, 1]) == np.sort(results_rotated[:, 0])).all() @@ -195,7 +211,8 @@ def test_subpix_edge(): img = np.zeros((50, 50)) img[:25, :25] = 255 img[25:, 25:] = 255 - corner = peak_local_max(corner_harris(img), num_peaks=1) + corner = peak_local_max(corner_harris(img), + min_distance=10, threshold_rel=0, num_peaks=1) subpix = corner_subpix(img, corner) assert_array_equal(subpix[0], (24.5, 24.5)) @@ -203,7 +220,8 @@ def test_subpix_edge(): def test_subpix_dot(): img = np.zeros((50, 50)) img[25, 25] = 255 - corner = peak_local_max(corner_harris(img), num_peaks=1) + corner = peak_local_max(corner_harris(img), + min_distance=10, threshold_rel=0, num_peaks=1) subpix = corner_subpix(img, corner) assert_array_equal(subpix[0], (25, 25)) @@ -214,7 +232,8 @@ def test_subpix_no_class(): assert_array_equal(subpix[0], (np.nan, np.nan)) img[25, 25] = 1e-10 - corner = peak_local_max(corner_harris(img), num_peaks=1) + corner = peak_local_max(corner_harris(img), + min_distance=10, threshold_rel=0, num_peaks=1) subpix = corner_subpix(img, np.array([[25, 25]])) assert_array_equal(subpix[0], (np.nan, np.nan)) @@ -223,7 +242,7 @@ def test_subpix_border(): img = np.zeros((50, 50)) img[1:25,1:25] = 255 img[25:-1,25:-1] = 255 - corner = corner_peaks(corner_harris(img), min_distance=1) + corner = corner_peaks(corner_harris(img), threshold_rel=0) subpix = corner_subpix(img, corner, window_size=11) ref = np.array([[ 0.52040816, 0.52040816], [ 0.52040816, 24.47959184], @@ -244,21 +263,23 @@ def test_num_peaks(): for i in range(20): n = np.random.random_integers(20) - results = peak_local_max(img_corners, num_peaks=n) + results = peak_local_max(img_corners, + min_distance=10, threshold_rel=0, num_peaks=n) assert (results.shape[0] == n) def test_corner_peaks(): - response = np.zeros((5, 5)) - response[2:4, 2:4] = 1 + response = np.zeros((10, 10)) + response[2:5, 2:5] = 1 - corners = corner_peaks(response, exclude_border=False) + corners = corner_peaks(response, exclude_border=False, min_distance=10, + threshold_rel=0) assert len(corners) == 1 - corners = corner_peaks(response, exclude_border=False, min_distance=0) + corners = corner_peaks(response, exclude_border=False, min_distance=1) assert len(corners) == 4 - corners = corner_peaks(response, exclude_border=False, min_distance=0, + corners = corner_peaks(response, exclude_border=False, min_distance=1, indices=False) assert np.sum(corners) == 4 @@ -323,7 +344,8 @@ def test_corner_fast_lena(): [492, 139], [494, 169], [496, 266]]) - actual = corner_peaks(corner_fast(img, 12, 0.3)) + actual = corner_peaks(corner_fast(img, 12, 0.3), + min_distance=10, threshold_rel=0) assert_array_equal(actual, expected) @@ -342,7 +364,8 @@ def test_corner_orientations_even_shape_error(): @test_parallel() def test_corner_orientations_astronaut(): img = rgb2gray(data.astronaut()) - corners = corner_peaks(corner_fast(img, 11, 0.35)) + corners = corner_peaks(corner_fast(img, 11, 0.35), + min_distance=10, threshold_abs=0, threshold_rel=0.1) expected = np.array([-1.75220190e+00, 2.01197383e+00, -2.01162417e+00, -1.88247204e-01, 1.19134149e+00, -6.61151410e-01, -2.99143370e+00, 2.17103132e+00, -7.52950306e-04, @@ -355,7 +378,6 @@ def test_corner_orientations_astronaut(): -4.40598471e-01, 3.14918803e-01, -1.76069982e+00, 3.05330950e+00, 2.39291733e+00, -1.22091334e-01, -3.09279990e-01, 1.45931342e+00]) - actual = corner_orientations(img, corners, octagon(3, 2)) assert_almost_equal(actual, expected) @@ -363,7 +385,8 @@ def test_corner_orientations_astronaut(): def test_corner_orientations_square(): square = np.zeros((12, 12)) square[3:9, 3:9] = 1 - corners = corner_peaks(corner_fast(square, 9), min_distance=1) + corners = corner_peaks(corner_fast(square, 9), + min_distance=1, threshold_rel=0) actual_orientations = corner_orientations(square, corners, octagon(3, 2)) actual_orientations_degrees = np.rad2deg(actual_orientations) expected_orientations_degree = np.array([ 45., 135., -45., -135.]) diff --git a/skimage/feature/tests/test_match.py b/skimage/feature/tests/test_match.py index 31d20d36..770dc37b 100644 --- a/skimage/feature/tests/test_match.py +++ b/skimage/feature/tests/test_match.py @@ -36,11 +36,13 @@ def test_binary_descriptors_lena_rotation_crosscheck_false(): extractor = BRIEF(descriptor_size=512) - keypoints1 = corner_peaks(corner_harris(img), min_distance=5) + keypoints1 = corner_peaks(corner_harris(img), min_distance=5, + threshold_abs=0, threshold_rel=0.1) extractor.extract(img, keypoints1) descriptors1 = extractor.descriptors - keypoints2 = corner_peaks(corner_harris(rotated_img), min_distance=5) + keypoints2 = corner_peaks(corner_harris(rotated_img), min_distance=5, + threshold_abs=0, threshold_rel=0.1) extractor.extract(rotated_img, keypoints2) descriptors2 = extractor.descriptors @@ -69,11 +71,13 @@ def test_binary_descriptors_lena_rotation_crosscheck_true(): extractor = BRIEF(descriptor_size=512) - keypoints1 = corner_peaks(corner_harris(img), min_distance=5) + keypoints1 = corner_peaks(corner_harris(img), min_distance=5, + threshold_abs=0, threshold_rel=0.1) extractor.extract(img, keypoints1) descriptors1 = extractor.descriptors - keypoints2 = corner_peaks(corner_harris(rotated_img), min_distance=5) + keypoints2 = corner_peaks(corner_harris(rotated_img), min_distance=5, + threshold_abs=0, threshold_rel=0.1) extractor.extract(rotated_img, keypoints2) descriptors2 = extractor.descriptors diff --git a/skimage/feature/tests/test_peak.py b/skimage/feature/tests/test_peak.py index e6256208..40292cd9 100644 --- a/skimage/feature/tests/test_peak.py +++ b/skimage/feature/tests/test_peak.py @@ -1,6 +1,6 @@ import numpy as np from numpy.testing import (assert_array_almost_equal as assert_close, - assert_equal) + assert_equal, assert_raises) from scipy import ndimage as ndi from skimage.feature import peak @@ -70,12 +70,14 @@ def test_num_peaks(): image[1, 5] = 12 image[3, 5] = 8 image[5, 3] = 7 - assert len(peak.peak_local_max(image, min_distance=1)) == 5 - peaks_limited = peak.peak_local_max(image, min_distance=1, num_peaks=2) + assert len(peak.peak_local_max(image, min_distance=1, threshold_abs=0)) == 5 + peaks_limited = peak.peak_local_max( + image, min_distance=1, threshold_abs=0, num_peaks=2) assert len(peaks_limited) == 2 assert (1, 3) in peaks_limited assert (1, 5) in peaks_limited - peaks_limited = peak.peak_local_max(image, min_distance=1, num_peaks=4) + peaks_limited = peak.peak_local_max( + image, min_distance=1, threshold_abs=0, num_peaks=4) assert len(peaks_limited) == 4 assert (1, 3) in peaks_limited assert (1, 5) in peaks_limited @@ -270,9 +272,11 @@ def test_disk(): result = peak.peak_local_max(image, labels=np.ones((10, 20)), footprint=footprint, min_distance=1, threshold_rel=0, - indices=False, exclude_border=False) + threshold_abs=-1, indices=False, + exclude_border=False) assert np.all(result) - result = peak.peak_local_max(image, footprint=footprint) + result = peak.peak_local_max(image, footprint=footprint, threshold_abs=-1, + indices=False, exclude_border=False) assert np.all(result) @@ -280,11 +284,14 @@ def test_3D(): image = np.zeros((30, 30, 30)) image[15, 15, 15] = 1 image[5, 5, 5] = 1 - assert_equal(peak.peak_local_max(image), [[15, 15, 15]]) - assert_equal(peak.peak_local_max(image, min_distance=6), [[15, 15, 15]]) - assert_equal(peak.peak_local_max(image, exclude_border=False), + assert_equal(peak.peak_local_max(image, min_distance=10, threshold_rel=0), + [[15, 15, 15]]) + assert_equal(peak.peak_local_max(image, min_distance=6, threshold_rel=0), + [[15, 15, 15]]) + assert_equal(peak.peak_local_max(image, min_distance=10, threshold_rel=0, + exclude_border=False), [[5, 5, 5], [15, 15, 15]]) - assert_equal(peak.peak_local_max(image, min_distance=5), + assert_equal(peak.peak_local_max(image, min_distance=5, threshold_rel=0), [[5, 5, 5], [15, 15, 15]]) @@ -292,14 +299,30 @@ def test_4D(): image = np.zeros((30, 30, 30, 30)) image[15, 15, 15, 15] = 1 image[5, 5, 5, 5] = 1 - assert_equal(peak.peak_local_max(image), [[15, 15, 15, 15]]) - assert_equal(peak.peak_local_max(image, min_distance=6), [[15, 15, 15, 15]]) - assert_equal(peak.peak_local_max(image, exclude_border=False), + assert_equal(peak.peak_local_max(image, min_distance=10, threshold_rel=0), + [[15, 15, 15, 15]]) + assert_equal(peak.peak_local_max(image, min_distance=6, threshold_rel=0), + [[15, 15, 15, 15]]) + assert_equal(peak.peak_local_max(image, min_distance=10, threshold_rel=0, + exclude_border=False), [[5, 5, 5, 5], [15, 15, 15, 15]]) - assert_equal(peak.peak_local_max(image, min_distance=5), + assert_equal(peak.peak_local_max(image, min_distance=5, threshold_rel=0), [[5, 5, 5, 5], [15, 15, 15, 15]]) +def test_threshold_rel_default(): + image = np.ones((5, 5)) + + image[2, 2] = 1 + assert len(peak.peak_local_max(image)) == 0 + + image[2, 2] = 2 + assert_equal(peak.peak_local_max(image), [[2, 2]]) + + image[2, 2] = 0 + assert len(peak.peak_local_max(image, min_distance=0)) == image.size - 1 + + if __name__ == '__main__': from numpy import testing testing.run_module_suite()