Merge pull request #1898 from ahojnnes/anntzer-peak-local-max

ENH: Improved defaults for peak_local_max
This commit is contained in:
Josh Warner
2016-01-29 20:53:19 -06:00
6 changed files with 153 additions and 94 deletions
+4 -9
View File
@@ -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]])
"""
+47 -35
View File
@@ -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]
+4 -2
View File
@@ -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')
+53 -30
View File
@@ -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.])
+8 -4
View File
@@ -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
+37 -14
View File
@@ -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()