wassname/scikit-image
1
0
Fork 0
mirror of https://github.com/wassname/scikit-image.git synced 2026-07-03 07:10:47 +08:00
Code Issues Packages Projects Releases Wiki Activity

Fix Gaussian kernels in hessian_matrix

Browse Source
This commit is contained in:
Johannes Schönberger
2016-01-20 22:11:23 +01:00
parent ce4bc59c11
commit 6fd0cecffd
2 changed files with 33 additions and 30 deletions
Download Patch File Download Diff File Expand all files Collapse all files
skimage/feature/corner.py
+8 -5
View File
@@ -149,22 +149,25 @@ def hessian_matrix(image, sigma=1, mode='constant', cval=0):
image = _prepare_grayscale_input_2D(image)
# window extent to the left and right, which covers > 99% of the normal
# distribution
# Window extent which covers > 99% of the normal distribution.
window_ext = max(1, np.ceil(3 * sigma))
ky, kx = np.mgrid[-window_ext:window_ext + 1, -window_ext:window_ext + 1]
# second derivative Gaussian kernels
# Second derivative Gaussian kernels.
gaussian_exp = np.exp(-(kx ** 2 + ky ** 2) / (2 * sigma ** 2))
kernel_xx = 1 / (2 * np.pi * sigma ** 4) * (kx ** 2 / sigma ** 2 - 1)
kernel_xx *= gaussian_exp
kernel_xx /= kernel_xx.sum()
kernel_xy = 1 / (2 * np.pi * sigma ** 6) * (kx * ky)
kernel_xy *= gaussian_exp
kernel_xy /= kernel_xx.sum()
kernel_yy = kernel_xx.transpose()
# Remove small kernel values.
eps = np.finfo(kernel_xx.dtype).eps
kernel_xx[np.abs(kernel_xx) < eps * np.abs(kernel_xx).max()] = 0
kernel_xy[np.abs(kernel_xy) < eps * np.abs(kernel_xy).max()] = 0
kernel_yy[np.abs(kernel_yy) < eps * np.abs(kernel_yy).max()] = 0
Hxx = ndi.convolve(image, kernel_xx, mode=mode, cval=cval)
Hxy = ndi.convolve(image, kernel_xy, mode=mode, cval=cval)
Hyy = ndi.convolve(image, kernel_yy, mode=mode, cval=cval)
skimage/feature/tests/test_corner.py
+25 -25
View File
@@ -42,21 +42,21 @@ def test_hessian_matrix():
square = np.zeros((5, 5))
square[2, 2] = 1
Hxx, Hxy, Hyy = hessian_matrix(square, sigma=0.1)
assert_array_equal(Hxx, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 1, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
assert_array_equal(Hxy, np.array([[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]]))
assert_array_equal(Hyy, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 1, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
assert_almost_equal(Hxx, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, -1591.549431, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
assert_almost_equal(Hxy, np.array([[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]]))
assert_almost_equal(Hyy, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, -1591.549431, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
def test_structure_tensor_eigvals():
@@ -81,16 +81,16 @@ def test_hessian_matrix_eigvals():
square[2, 2] = 1
Hxx, Hxy, Hyy = hessian_matrix(square, sigma=0.1)
l1, l2 = hessian_matrix_eigvals(Hxx, Hxy, Hyy)
assert_array_equal(l1, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 1, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
assert_array_equal(l2, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 1, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
assert_almost_equal(l1, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, -1591.549431, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
assert_almost_equal(l2, np.array([[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, -1591.549431, 0, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0]]))
@test_parallel()