diff --git a/skimage/feature/__init__.py b/skimage/feature/__init__.py index 7d1d9add..6a933056 100644 --- a/skimage/feature/__init__.py +++ b/skimage/feature/__init__.py @@ -6,7 +6,7 @@ from .corner import (corner_kitchen_rosenfeld, corner_harris, corner_shi_tomasi, corner_foerstner, corner_subpix, corner_peaks, corner_fast, structure_tensor, structure_tensor_eigvals, hessian_matrix, - hessian_matrix_eigvals) + hessian_matrix_eigvals, hessian_matrix_det) from .corner_cy import corner_moravec, corner_orientations from .template import match_template from .brief import BRIEF @@ -15,7 +15,6 @@ from .orb import ORB from .match import match_descriptors from .util import plot_matches from .blob import blob_dog, blob_log, blob_doh -from ._hessian_det_appx import hessian_det_appx __all__ = ['daisy', @@ -27,6 +26,7 @@ __all__ = ['daisy', 'structure_tensor', 'structure_tensor_eigvals', 'hessian_matrix', + 'hessian_matrx_det', 'hessian_matrix_eigvals', 'corner_kitchen_rosenfeld', 'corner_harris', @@ -45,5 +45,5 @@ __all__ = ['daisy', 'plot_matches', 'blob_dog', 'blob_doh', - 'hessian_det_appx', - 'blob_log'] + 'blob_log', + 'hessian_matrix_det'] diff --git a/skimage/feature/_hessian_det_appx.pyx b/skimage/feature/_hessian_det_appx.pyx index 5d3f7511..e66c0d53 100644 --- a/skimage/feature/_hessian_det_appx.pyx +++ b/skimage/feature/_hessian_det_appx.pyx @@ -78,7 +78,7 @@ cdef inline cnp.double_t _integ( return ans -def hessian_det_appx(cnp.double_t[:, ::1] img, float sigma): +def _hessian_matrix_det(cnp.double_t[:, ::1] img, float sigma): """Computes the approximate Hessian Determinant over an image. This method uses box filters over integral images to compute the diff --git a/skimage/feature/blob.py b/skimage/feature/blob.py index 57818a98..9c35d2a5 100644 --- a/skimage/feature/blob.py +++ b/skimage/feature/blob.py @@ -1,3 +1,4 @@ + import numpy as np from scipy.ndimage.filters import gaussian_filter, gaussian_laplace import itertools as itt @@ -6,7 +7,7 @@ from math import sqrt, hypot, log from numpy import arccos from skimage.util import img_as_float from .peak import peak_local_max -from ._hessian_det_appx import hessian_det_appx +from ._hessian_det_appx import _hessian_matrix_det from skimage.transform import integral_image @@ -398,7 +399,7 @@ def blob_doh(image, min_sigma=1, max_sigma=30, num_sigma=10, threshold=0.01, else: sigma_list = np.linspace(min_sigma, max_sigma, num_sigma) - hessian_images = [hessian_det_appx(image, s) for s in sigma_list] + hessian_images = [_hessian_matrix_det(image, s) for s in sigma_list] image_cube = np.dstack(hessian_images) local_maxima = peak_local_max(image_cube, threshold_abs=threshold, diff --git a/skimage/feature/corner.py b/skimage/feature/corner.py index f4894aa5..67bed0a1 100644 --- a/skimage/feature/corner.py +++ b/skimage/feature/corner.py @@ -7,6 +7,8 @@ from skimage.util import img_as_float, pad from skimage.feature import peak_local_max from skimage.feature.util import _prepare_grayscale_input_2D from skimage.feature.corner_cy import _corner_fast +from ._hessian_det_appx import _hessian_matrix_det +from ..transform import integral_image def _compute_derivatives(image, mode='constant', cval=0): @@ -29,7 +31,7 @@ def _compute_derivatives(image, mode='constant', cval=0): imy : ndarray Derivative in y-direction. - """ +v """ imy = ndimage.sobel(image, axis=0, mode=mode, cval=cval) imx = ndimage.sobel(image, axis=1, mode=mode, cval=cval) @@ -170,6 +172,52 @@ def hessian_matrix(image, sigma=1, mode='constant', cval=0): return Hxx, Hxy, Hyy +def hessian_matrix_det(image, sigma, integral=True): + """Computes the approximate Hessian Determinant over an image. + + This method uses box filters over integral images to compute the + approximate Hessian Determinant as described in [1]_. + + Parameters + ---------- + image : array + The integral image over which to compute Hessian Determinant. + sigma : float + Standard deviation used for the Gaussian kernel, used for the Hessian + matrix + integral : bool + If `False`, `image` is assumed to be integral and intergral image is + not computed. If `True` the integral image is computed for `image` + and used for finding the Hessian Determinant. + + + Returns + ------- + out : array + The array of the Determinant of Hessians. + + References + ---------- + .. [1] Herbert Bay, Andreas Ess, Tinne Tuytelaars, Luc Van Gool, + "SURF: Speeded Up Robust Features" + ftp://ftp.vision.ee.ethz.ch/publications/articles/eth_biwi_00517.pdf + + Notes + ----- + The running time of this method only depends on size of the image. It is + independent of `sigma` as one would expect. The downside is that the + result for `sigma` less than `3` is not accurate, i.e., not similar to + the result obtained if someone computed the Hessian and took it's + determinant. + """ + + image = img_as_float(image) + if(integral): + image = integral_image(image) + + return np.array(_hessian_matrix_det(image, sigma)) + + def _image_orthogonal_matrix22_eigvals(M00, M01, M11): l1 = (M00 + M11) / 2 + np.sqrt(4 * M01 ** 2 + (M00 - M11) ** 2) / 2 l2 = (M00 + M11) / 2 - np.sqrt(4 * M01 ** 2 + (M00 - M11) ** 2) / 2 diff --git a/skimage/feature/tests/test_corner.py b/skimage/feature/tests/test_corner.py index 21070779..dd9f7703 100644 --- a/skimage/feature/tests/test_corner.py +++ b/skimage/feature/tests/test_corner.py @@ -12,7 +12,8 @@ from skimage.feature import (corner_moravec, corner_harris, corner_shi_tomasi, corner_kitchen_rosenfeld, corner_foerstner, corner_fast, corner_orientations, structure_tensor, structure_tensor_eigvals, - hessian_matrix, hessian_matrix_eigvals) + hessian_matrix, hessian_matrix_eigvals, + hessian_matrix_det) def test_structure_tensor(): @@ -91,6 +92,12 @@ def test_hessian_matrix_eigvals(): [0, 0, 0, 0, 0]])) +def test_hessian_matrix_det(): + image = np.ones((5, 5)) + det = hessian_matrix_det(image, 3, False) + assert_array_equal(det, 0) + + def test_square_image(): im = np.zeros((50, 50)).astype(float) im[:25, :25] = 1.