From 781bbaf5bceb19100835e1fd52441d71e2be9ee0 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 22 Jul 2013 22:17:14 +0530 Subject: [PATCH 01/77] Added partial implementation of censure for mode=DoB --- skimage/feature/censure.py | 47 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 47 insertions(+) create mode 100644 skimage/feature/censure.py diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py new file mode 100644 index 00000000..402aa4e6 --- /dev/null +++ b/skimage/feature/censure.py @@ -0,0 +1,47 @@ +import numpy as np +from scipy.ndimage.filters import maximum_filter, minimum_filter +from skimage.transform import integral_image +import time + +def _get_filtered_image(image, n, mode='DoB'): + if mode == 'DoB': + inner_wt = (1.0 / (2*n + 1)**2) + outer_wt = (1.0 / (12*n**2 + 4*n)) + integral_img = integral_image(image) + filtered_image = np.zeros(image.shape) + # TODO : Outsource to Cython + start = time.time() + for i in range(2 * n, image.shape[0] - 2 * n): + for j in range(2 * n, image.shape[1] - 2 * n): + inner = integral_img[i + n, j + n] + integral_img[i - n, j - n] - integral_img[i + n, j - n] - integral_img[i - n, j + n] + outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n, j - 2 * n] - integral_img[i + 2 * n, j - 2 * n] - integral_img[i - 2 * n, j + 2 * n] + filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner + print time.time() - start + return filtered_image + + +def censure_keypoints(image, mode='DoB', threshold=0.1): + # TODO : Decide mode for convolve function + # TODO : Decide number of scales. Image-size dependent? + image = np.squeeze(image) + if image.ndim != 2: + raise ValueError("Only 2-D gray-scale images supported.") + # Generating all the scales + scale1 = _get_filtered_image(image, 1, mode) + scale2 = _get_filtered_image(image, 2, mode) + scale3 = _get_filtered_image(image, 3, mode) + scale4 = _get_filtered_image(image, 4, mode) + scale5 = _get_filtered_image(image, 5, mode) + scale6 = _get_filtered_image(image, 6, mode) + scale7 = _get_filtered_image(image, 7, mode) + # Stacking all the scales in the 3rd dimension + scales = np.dstack((scale1, scale2, scale3, scale4, scale5, scale6, scale7)) + # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 + # neighbourhood to zero + minimas = (minimum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales + maximas = (maximum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales + # Suppressing minimas and maximas weaker than threshold + minimas[np.abs(minimas) < threshold] = 0 + maximas[np.abs(maximas) < threshold] = 0 + response = maximas + np.abs(minimas) + return response From c40ab969afe1d40876f0342adc1a1d0c458e5f52 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 25 Jul 2013 00:03:50 +0530 Subject: [PATCH 02/77] Line Suppression of the response --- skimage/feature/__init__.py | 4 +++- skimage/feature/censure.py | 22 ++++++++++++++++++---- 2 files changed, 21 insertions(+), 5 deletions(-) diff --git a/skimage/feature/__init__.py b/skimage/feature/__init__.py index d0c51fb0..53c21d31 100644 --- a/skimage/feature/__init__.py +++ b/skimage/feature/__init__.py @@ -9,6 +9,7 @@ from .corner_cy import corner_moravec from .template import match_template from ._brief import brief, match_keypoints_brief from .util import pairwise_hamming_distance +from .censure import censure_keypoints __all__ = ['daisy', 'hog', @@ -26,4 +27,5 @@ __all__ = ['daisy', 'match_template', 'brief', 'pairwise_hamming_distance', - 'match_keypoints_brief'] + 'match_keypoints_brief', + 'censure_keypoints'] diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 402aa4e6..e4a862d5 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -1,6 +1,7 @@ import numpy as np from scipy.ndimage.filters import maximum_filter, minimum_filter from skimage.transform import integral_image +from skimage.feature import _compute_auto_correlation import time def _get_filtered_image(image, n, mode='DoB'): @@ -13,15 +14,23 @@ def _get_filtered_image(image, n, mode='DoB'): start = time.time() for i in range(2 * n, image.shape[0] - 2 * n): for j in range(2 * n, image.shape[1] - 2 * n): - inner = integral_img[i + n, j + n] + integral_img[i - n, j - n] - integral_img[i + n, j - n] - integral_img[i - n, j + n] - outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n, j - 2 * n] - integral_img[i + 2 * n, j - 2 * n] - integral_img[i - 2 * n, j + 2 * n] + inner = integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n] + outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n - 1, j - 2 * n - 1] - integral_img[i + 2 * n, j - 2 * n - 1] - integral_img[i - 2 * n - 1, j + 2 * n] filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner print time.time() - start return filtered_image -def censure_keypoints(image, mode='DoB', threshold=0.1): - # TODO : Decide mode for convolve function +def _suppress_line(response, sigma): + Axx, Axy, Ayy = _compute_auto_correlation(response, sigma) + detA = Axx * Ayy - Axy**2 + traceA = Axx + Ayy + # ratio of principal curvatures + rpc = traceA / detA + rpc[rpc > 10] = 0 + return rpc + +def censure_keypoints(image, mode='DoB', threshold=0.03): # TODO : Decide number of scales. Image-size dependent? image = np.squeeze(image) if image.ndim != 2: @@ -44,4 +53,9 @@ def censure_keypoints(image, mode='DoB', threshold=0.1): minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 response = maximas + np.abs(minimas) + response[:, :, 1] = _suppress_line(response[:, :, 1], 1.33) + response[:, :, 2] = _suppress_line(response[:, :, 2], 1.33) + response[:, :, 3] = _suppress_line(response[:, :, 3], 1.33) + response[:, :, 4] = _suppress_line(response[:, :, 4], 1.33) + response[:, :, 5] = _suppress_line(response[:, :, 5], 1.33) return response From c599b342361e7c4d8a01a624a3f342205afa2976 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 25 Jul 2013 21:06:17 +0530 Subject: [PATCH 03/77] Correcting import for a private function --- skimage/feature/censure.py | 22 ++++++++++++---------- 1 file changed, 12 insertions(+), 10 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index e4a862d5..f7981c98 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -1,9 +1,10 @@ import numpy as np from scipy.ndimage.filters import maximum_filter, minimum_filter from skimage.transform import integral_image -from skimage.feature import _compute_auto_correlation +from skimage.feature.corner import _compute_auto_correlation import time + def _get_filtered_image(image, n, mode='DoB'): if mode == 'DoB': inner_wt = (1.0 / (2*n + 1)**2) @@ -21,14 +22,15 @@ def _get_filtered_image(image, n, mode='DoB'): return filtered_image -def _suppress_line(response, sigma): +def _suppress_line(response, sigma, rpc_threshold): Axx, Axy, Ayy = _compute_auto_correlation(response, sigma) detA = Axx * Ayy - Axy**2 traceA = Axx + Ayy # ratio of principal curvatures - rpc = traceA / detA - rpc[rpc > 10] = 0 - return rpc + rpc = traceA / (detA + 0.001) + response[rpc > rpc_threshold] = 0 + return response + def censure_keypoints(image, mode='DoB', threshold=0.03): # TODO : Decide number of scales. Image-size dependent? @@ -53,9 +55,9 @@ def censure_keypoints(image, mode='DoB', threshold=0.03): minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 response = maximas + np.abs(minimas) - response[:, :, 1] = _suppress_line(response[:, :, 1], 1.33) - response[:, :, 2] = _suppress_line(response[:, :, 2], 1.33) - response[:, :, 3] = _suppress_line(response[:, :, 3], 1.33) - response[:, :, 4] = _suppress_line(response[:, :, 4], 1.33) - response[:, :, 5] = _suppress_line(response[:, :, 5], 1.33) + response[:, :, 1] = _suppress_line(response[:, :, 1], 1.33, 10) + response[:, :, 2] = _suppress_line(response[:, :, 2], 1.33, 10) + response[:, :, 3] = _suppress_line(response[:, :, 3], 1.33, 10) + response[:, :, 4] = _suppress_line(response[:, :, 4], 1.33, 10) + response[:, :, 5] = _suppress_line(response[:, :, 5], 1.33, 10) return response From 97990fbea562bfc441b6b7ae0862f04198982d52 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 25 Jul 2013 21:28:48 +0530 Subject: [PATCH 04/77] Adding TODO's --- skimage/feature/censure.py | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index f7981c98..2aa2057e 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -6,6 +6,7 @@ import time def _get_filtered_image(image, n, mode='DoB'): + # TODO : Implement the STAR and Octagon mode if mode == 'DoB': inner_wt = (1.0 / (2*n + 1)**2) outer_wt = (1.0 / (12*n**2 + 4*n)) @@ -32,7 +33,7 @@ def _suppress_line(response, sigma, rpc_threshold): return response -def censure_keypoints(image, mode='DoB', threshold=0.03): +def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): # TODO : Decide number of scales. Image-size dependent? image = np.squeeze(image) if image.ndim != 2: @@ -55,9 +56,12 @@ def censure_keypoints(image, mode='DoB', threshold=0.03): minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 response = maximas + np.abs(minimas) - response[:, :, 1] = _suppress_line(response[:, :, 1], 1.33, 10) - response[:, :, 2] = _suppress_line(response[:, :, 2], 1.33, 10) - response[:, :, 3] = _suppress_line(response[:, :, 3], 1.33, 10) - response[:, :, 4] = _suppress_line(response[:, :, 4], 1.33, 10) - response[:, :, 5] = _suppress_line(response[:, :, 5], 1.33, 10) + # TODO : Decide the rpc_threshold and sigma for all the scales. The paper only discusses + # values for scale2 i.e. response[:, :, 1] + response[:, :, 1] = _suppress_line(response[:, :, 1], 1.33, rpc_threshold) + response[:, :, 2] = _suppress_line(response[:, :, 2], 1.33, rpc_threshold) + response[:, :, 3] = _suppress_line(response[:, :, 3], 1.33, rpc_threshold) + response[:, :, 4] = _suppress_line(response[:, :, 4], 1.33, rpc_threshold) + response[:, :, 5] = _suppress_line(response[:, :, 5], 1.33, rpc_threshold) + # TODO : Return key-points from all the scales? return response From 3d9534b8448dae075e0717ffbba23a2bcf9ae4fc Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sat, 27 Jul 2013 01:58:36 +0530 Subject: [PATCH 05/77] Added mode Octagon using np.convolve --- skimage/feature/censure.py | 66 +++++++++++++++++++++++++++++++++----- 1 file changed, 58 insertions(+), 8 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 2aa2057e..6a707970 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -1,10 +1,12 @@ import numpy as np -from scipy.ndimage.filters import maximum_filter, minimum_filter +from scipy.ndimage.filters import maximum_filter, minimum_filter, convolve from skimage.transform import integral_image from skimage.feature.corner import _compute_auto_correlation +from skimage.morphology import convex_hull_image import time +""" def _get_filtered_image(image, n, mode='DoB'): # TODO : Implement the STAR and Octagon mode if mode == 'DoB': @@ -21,6 +23,52 @@ def _get_filtered_image(image, n, mode='DoB'): filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner print time.time() - start return filtered_image + elif mode == 'Octagon': + outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] + inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] +""" + + +def _oct(m, n): + f = np.zeros((m + 2*n, m + 2*n)) + f[0, n] = 1 + f[n, 0] = 1 + f[0, m + n -1] = 1 + f[m + n - 1, 0] = 1 + f[-1, n] = 1 + f[n, -1] = 1 + f[-1, m + n - 1] = 1 + f[m + n - 1, -1] = 1 + return convex_hull_image(f).astype(int) + + +def _octagon_filter(mo, no, mi, ni): + outer = (mo + 2 * no)**2 - 2 * no * (no + 1) + inner = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) + outer_wt = 1.0 / (outer - inner) + inner_wt = 1.0 / inner + c = ((mo + 2 * no) - (mi + 2 * ni)) / 2 + outer_oct = _oct(mo, no) + inner_oct = np.zeros((mo + 2 * no, mo + 2 * no)) + inner_oct[c:-c, c:-c] = _oct(mi, ni) + bfilter = outer_wt * outer_oct - (outer_wt + inner_wt) * inner_oct + return bfilter + + +def _filter_using_convolve(image, n, mode='DoB'): + + if mode == 'DoB': + inner_wt = (1.0 / (2*n + 1)**2) + outer_wt = (1.0 / (12*n**2 + 4*n)) + dob_filter = np.zeros((4 * n + 1, 4 * n + 1)) + dob_filter[:] = outer_wt + dob_filter[n : 3 * n + 1, n : 3 * n + 1] = - inner_wt + return convolve(image, dob_filter) + + elif mode == 'Octagon': + outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] + inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] + return convolve(image, _octagon_filter(outer_shape[n - 1][0], outer_shape[n - 1][1], inner_shape[n - 1][0], inner_shape[n - 1][1])) def _suppress_line(response, sigma, rpc_threshold): @@ -39,13 +87,15 @@ def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") # Generating all the scales - scale1 = _get_filtered_image(image, 1, mode) - scale2 = _get_filtered_image(image, 2, mode) - scale3 = _get_filtered_image(image, 3, mode) - scale4 = _get_filtered_image(image, 4, mode) - scale5 = _get_filtered_image(image, 5, mode) - scale6 = _get_filtered_image(image, 6, mode) - scale7 = _get_filtered_image(image, 7, mode) + start = time.time() + scale1 = _filter_using_convolve(image, 1, mode) + scale2 = _filter_using_convolve(image, 2, mode) + scale3 = _filter_using_convolve(image, 3, mode) + scale4 = _filter_using_convolve(image, 4, mode) + scale5 = _filter_using_convolve(image, 5, mode) + scale6 = _filter_using_convolve(image, 6, mode) + scale7 = _filter_using_convolve(image, 7, mode) + print time.time - start # Stacking all the scales in the 3rd dimension scales = np.dstack((scale1, scale2, scale3, scale4, scale5, scale6, scale7)) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 From 1484b9a72ad5ca9fb83f2a2a03ff26d527ae8521 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sat, 27 Jul 2013 02:07:19 +0530 Subject: [PATCH 06/77] Correcting typo --- skimage/feature/censure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 6a707970..f682735c 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -95,7 +95,7 @@ def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): scale5 = _filter_using_convolve(image, 5, mode) scale6 = _filter_using_convolve(image, 6, mode) scale7 = _filter_using_convolve(image, 7, mode) - print time.time - start + print time.time() - start # Stacking all the scales in the 3rd dimension scales = np.dstack((scale1, scale2, scale3, scale4, scale5, scale6, scale7)) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 From 1bbeeb45d6b29e756c3e7442d9506f918d47a2c5 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sat, 27 Jul 2013 16:34:28 +0530 Subject: [PATCH 07/77] For loop for mode=DoB in Cython --- bento.info | 3 +++ skimage/feature/censure.py | 41 ++++++++++++++++------------------ skimage/feature/censure_cy.pyx | 21 +++++++++++++++++ skimage/feature/setup.py | 3 +++ 4 files changed, 46 insertions(+), 22 deletions(-) create mode 100644 skimage/feature/censure_cy.pyx diff --git a/bento.info b/bento.info index 648445d6..acdb155d 100644 --- a/bento.info +++ b/bento.info @@ -90,6 +90,9 @@ Library: Extension: skimage.morphology._greyreconstruct Sources: skimage/morphology/_greyreconstruct.pyx + Extension: skimage.feature.censure_cy + Sources: + skimage/feature/censure_cy.pyx Extension: skimage.feature._brief_cy Sources: skimage/feature/_brief_cy.pyx diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index f682735c..b6a1d31f 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -1,32 +1,27 @@ import numpy as np from scipy.ndimage.filters import maximum_filter, minimum_filter, convolve -from skimage.transform import integral_image -from skimage.feature.corner import _compute_auto_correlation -from skimage.morphology import convex_hull_image + +from ..transform import integral_image +from ..feature.corner import _compute_auto_correlation +from ..morphology import convex_hull_image +from ..util import img_as_float + +from .censure_cy import _censure_dob_loop import time -""" def _get_filtered_image(image, n, mode='DoB'): # TODO : Implement the STAR and Octagon mode if mode == 'DoB': - inner_wt = (1.0 / (2*n + 1)**2) - outer_wt = (1.0 / (12*n**2 + 4*n)) + inner_wt = (1.0 / (2 * n + 1)**2) + outer_wt = (1.0 / (12 * n**2 + 4 * n)) integral_img = integral_image(image) filtered_image = np.zeros(image.shape) # TODO : Outsource to Cython start = time.time() - for i in range(2 * n, image.shape[0] - 2 * n): - for j in range(2 * n, image.shape[1] - 2 * n): - inner = integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n] - outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n - 1, j - 2 * n - 1] - integral_img[i + 2 * n, j - 2 * n - 1] - integral_img[i - 2 * n - 1, j + 2 * n] - filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner + _censure_dob_loop(image, n, integral_img, filtered_image, inner_wt, outer_wt) print time.time() - start return filtered_image - elif mode == 'Octagon': - outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] - inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] -""" def _oct(m, n): @@ -86,15 +81,17 @@ def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): image = np.squeeze(image) if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") + image = img_as_float(image) # Generating all the scales + image = np.ascontiguousarray(image) start = time.time() - scale1 = _filter_using_convolve(image, 1, mode) - scale2 = _filter_using_convolve(image, 2, mode) - scale3 = _filter_using_convolve(image, 3, mode) - scale4 = _filter_using_convolve(image, 4, mode) - scale5 = _filter_using_convolve(image, 5, mode) - scale6 = _filter_using_convolve(image, 6, mode) - scale7 = _filter_using_convolve(image, 7, mode) + scale1 = _get_filtered_image(image, 1, mode) + scale2 = _get_filtered_image(image, 2, mode) + scale3 = _get_filtered_image(image, 3, mode) + scale4 = _get_filtered_image(image, 4, mode) + scale5 = _get_filtered_image(image, 5, mode) + scale6 = _get_filtered_image(image, 6, mode) + scale7 = _get_filtered_image(image, 7, mode) print time.time() - start # Stacking all the scales in the 3rd dimension scales = np.dstack((scale1, scale2, scale3, scale4, scale5, scale6, scale7)) diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx new file mode 100644 index 00000000..b28bc5e1 --- /dev/null +++ b/skimage/feature/censure_cy.pyx @@ -0,0 +1,21 @@ +#cython: cdivision=True +#cython: boundscheck=False +#cython: nonecheck=False +#cython: wraparound=False + +cimport numpy as cnp + + +def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, + double[:, ::1] integral_img, + double[:, ::1] filtered_image, + cnp.float_t inner_wt, cnp.float_t outer_wt): + + cdef Py_ssize_t i, j + cdef double inner, outer + + for i in range(2 * n, image.shape[0] - 2 * n): + for j in range(2 * n, image.shape[1] - 2 * n): + inner = integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n] + outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n - 1, j - 2 * n - 1] - integral_img[i + 2 * n, j - 2 * n - 1] - integral_img[i - 2 * n - 1, j + 2 * n] + filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner diff --git a/skimage/feature/setup.py b/skimage/feature/setup.py index f4c62957..7df64c32 100644 --- a/skimage/feature/setup.py +++ b/skimage/feature/setup.py @@ -13,12 +13,15 @@ def configuration(parent_package='', top_path=None): config.add_data_dir('tests') cython(['corner_cy.pyx'], working_path=base_path) + cython(['censure_cy.pyx'], working_path=base_path) cython(['_brief_cy.pyx'], working_path=base_path) cython(['_texture.pyx'], working_path=base_path) cython(['_template.pyx'], working_path=base_path) config.add_extension('corner_cy', sources=['corner_cy.c'], include_dirs=[get_numpy_include_dirs()]) + config.add_extension('censure_cy', sources=['censure_cy.c'], + include_dirs=[get_numpy_include_dirs()]) config.add_extension('_brief_cy', sources=['_brief_cy.c'], include_dirs=[get_numpy_include_dirs()]) config.add_extension('_texture', sources=['_texture.c'], From a0266c3834e7d5c9d635f19249166d3e758222b5 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sun, 28 Jul 2013 00:55:41 +0530 Subject: [PATCH 08/77] Commenting out _filter_using_convolve --- skimage/feature/censure.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index b6a1d31f..cacf8e2f 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -50,6 +50,7 @@ def _octagon_filter(mo, no, mi, ni): return bfilter +""" def _filter_using_convolve(image, n, mode='DoB'): if mode == 'DoB': @@ -64,7 +65,7 @@ def _filter_using_convolve(image, n, mode='DoB'): outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] return convolve(image, _octagon_filter(outer_shape[n - 1][0], outer_shape[n - 1][1], inner_shape[n - 1][0], inner_shape[n - 1][1])) - +""" def _suppress_line(response, sigma, rpc_threshold): Axx, Axy, Ayy = _compute_auto_correlation(response, sigma) From ab729a00cac21ace8a118dc045c115778490e299 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sun, 28 Jul 2013 13:10:44 +0530 Subject: [PATCH 09/77] Removing time() statements and correcting types in censure_cy --- skimage/feature/censure.py | 12 +++++------- skimage/feature/censure_cy.pyx | 2 +- 2 files changed, 6 insertions(+), 8 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index cacf8e2f..ca361cf7 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -7,7 +7,6 @@ from ..morphology import convex_hull_image from ..util import img_as_float from .censure_cy import _censure_dob_loop -import time def _get_filtered_image(image, n, mode='DoB'): @@ -17,10 +16,8 @@ def _get_filtered_image(image, n, mode='DoB'): outer_wt = (1.0 / (12 * n**2 + 4 * n)) integral_img = integral_image(image) filtered_image = np.zeros(image.shape) - # TODO : Outsource to Cython - start = time.time() _censure_dob_loop(image, n, integral_img, filtered_image, inner_wt, outer_wt) - print time.time() - start + return filtered_image @@ -83,9 +80,10 @@ def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") image = img_as_float(image) - # Generating all the scales + image = np.ascontiguousarray(image) - start = time.time() + + # Generating all the scales scale1 = _get_filtered_image(image, 1, mode) scale2 = _get_filtered_image(image, 2, mode) scale3 = _get_filtered_image(image, 3, mode) @@ -93,7 +91,7 @@ def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): scale5 = _get_filtered_image(image, 5, mode) scale6 = _get_filtered_image(image, 6, mode) scale7 = _get_filtered_image(image, 7, mode) - print time.time() - start + # Stacking all the scales in the 3rd dimension scales = np.dstack((scale1, scale2, scale3, scale4, scale5, scale6, scale7)) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index b28bc5e1..d5e8bf02 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -9,7 +9,7 @@ cimport numpy as cnp def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, double[:, ::1] integral_img, double[:, ::1] filtered_image, - cnp.float_t inner_wt, cnp.float_t outer_wt): + double inner_wt, double outer_wt): cdef Py_ssize_t i, j cdef double inner, outer From 92864fdfb0501ab0ebaf22b6939c8bc5ea8e5f8c Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sun, 28 Jul 2013 13:30:04 +0530 Subject: [PATCH 10/77] Removing hard-coding in censure_keypoints() --- skimage/feature/censure.py | 25 ++++++++----------------- 1 file changed, 8 insertions(+), 17 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index ca361cf7..89114a07 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -74,7 +74,7 @@ def _suppress_line(response, sigma, rpc_threshold): return response -def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): +def censure_keypoints(image, mode='DoB', no_of_scales=7, threshold=0.03, rpc_threshold=10): # TODO : Decide number of scales. Image-size dependent? image = np.squeeze(image) if image.ndim != 2: @@ -84,16 +84,10 @@ def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): image = np.ascontiguousarray(image) # Generating all the scales - scale1 = _get_filtered_image(image, 1, mode) - scale2 = _get_filtered_image(image, 2, mode) - scale3 = _get_filtered_image(image, 3, mode) - scale4 = _get_filtered_image(image, 4, mode) - scale5 = _get_filtered_image(image, 5, mode) - scale6 = _get_filtered_image(image, 6, mode) - scale7 = _get_filtered_image(image, 7, mode) + scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) + for i in xrange(no_of_scales): + scales[:, :, i] = _get_filtered_image(image, i + 1, mode) - # Stacking all the scales in the 3rd dimension - scales = np.dstack((scale1, scale2, scale3, scale4, scale5, scale6, scale7)) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero minimas = (minimum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales @@ -102,12 +96,9 @@ def censure_keypoints(image, mode='DoB', threshold=0.03, rpc_threshold=10): minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 response = maximas + np.abs(minimas) - # TODO : Decide the rpc_threshold and sigma for all the scales. The paper only discusses - # values for scale2 i.e. response[:, :, 1] - response[:, :, 1] = _suppress_line(response[:, :, 1], 1.33, rpc_threshold) - response[:, :, 2] = _suppress_line(response[:, :, 2], 1.33, rpc_threshold) - response[:, :, 3] = _suppress_line(response[:, :, 3], 1.33, rpc_threshold) - response[:, :, 4] = _suppress_line(response[:, :, 4], 1.33, rpc_threshold) - response[:, :, 5] = _suppress_line(response[:, :, 5], 1.33, rpc_threshold) + + for i in xrange(1, no_of_scales - 1): + response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) + # TODO : Return key-points from all the scales? return response From 1e4ae4609b1b9a98983c59bebf433e94bcb5f3ed Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sun, 28 Jul 2013 13:43:33 +0530 Subject: [PATCH 11/77] Returning keypoints with the scale info --- skimage/feature/censure.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 89114a07..cd34e010 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -100,5 +100,6 @@ def censure_keypoints(image, mode='DoB', no_of_scales=7, threshold=0.03, rpc_thr for i in xrange(1, no_of_scales - 1): response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) - # TODO : Return key-points from all the scales? - return response + # Returning keypoints with its scale + keypoints = np.transpose(np.nonzero(response[:, :, 1:no_of_scales])) + [0, 0, 1] + return keypoints From 25f89eeaa74cf7d87daec27ecab271a312a907d0 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sun, 28 Jul 2013 13:52:07 +0530 Subject: [PATCH 12/77] Correcting the expression of response in censure_keypoints() --- skimage/feature/censure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index cd34e010..a38edf04 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -95,7 +95,7 @@ def censure_keypoints(image, mode='DoB', no_of_scales=7, threshold=0.03, rpc_thr # Suppressing minimas and maximas weaker than threshold minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 - response = maximas + np.abs(minimas) + response = maximas + minimas for i in xrange(1, no_of_scales - 1): response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) From fdbdbfe108813b6978e7f5319d029f4d1be80a65 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sun, 28 Jul 2013 14:13:09 +0530 Subject: [PATCH 13/77] Correcting the definition of Ratio of Principal Curvatures --- skimage/feature/censure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index a38edf04..c5d81373 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -69,7 +69,7 @@ def _suppress_line(response, sigma, rpc_threshold): detA = Axx * Ayy - Axy**2 traceA = Axx + Ayy # ratio of principal curvatures - rpc = traceA / (detA + 0.001) + rpc = traceA**2 / (detA + 0.001) response[rpc > rpc_threshold] = 0 return response From abb7bb65d426a7cdb4723c81dc7f0208a209e05f Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 29 Jul 2013 09:20:56 +0530 Subject: [PATCH 14/77] Adding _slanted_integral_image. Removed _oct and _filter_using_convolve --- skimage/feature/censure.py | 139 +++++++++++++++++++++++++------------ 1 file changed, 93 insertions(+), 46 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index c5d81373..0b8d2b1b 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -1,12 +1,12 @@ import numpy as np -from scipy.ndimage.filters import maximum_filter, minimum_filter, convolve +from scipy.ndimage.filters import maximum_filter, minimum_filter from ..transform import integral_image from ..feature.corner import _compute_auto_correlation from ..morphology import convex_hull_image from ..util import img_as_float -from .censure_cy import _censure_dob_loop +from .censure_cy import _censure_dob_loop, _slanted_integral_image def _get_filtered_image(image, n, mode='DoB'): @@ -17,52 +17,99 @@ def _get_filtered_image(image, n, mode='DoB'): integral_img = integral_image(image) filtered_image = np.zeros(image.shape) _censure_dob_loop(image, n, integral_img, filtered_image, inner_wt, outer_wt) - return filtered_image - - -def _oct(m, n): - f = np.zeros((m + 2*n, m + 2*n)) - f[0, n] = 1 - f[n, 0] = 1 - f[0, m + n -1] = 1 - f[m + n - 1, 0] = 1 - f[-1, n] = 1 - f[n, -1] = 1 - f[-1, m + n - 1] = 1 - f[m + n - 1, -1] = 1 - return convex_hull_image(f).astype(int) - - -def _octagon_filter(mo, no, mi, ni): - outer = (mo + 2 * no)**2 - 2 * no * (no + 1) - inner = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) - outer_wt = 1.0 / (outer - inner) - inner_wt = 1.0 / inner - c = ((mo + 2 * no) - (mi + 2 * ni)) / 2 - outer_oct = _oct(mo, no) - inner_oct = np.zeros((mo + 2 * no, mo + 2 * no)) - inner_oct[c:-c, c:-c] = _oct(mi, ni) - bfilter = outer_wt * outer_oct - (outer_wt + inner_wt) * inner_oct - return bfilter - - -""" -def _filter_using_convolve(image, n, mode='DoB'): - - if mode == 'DoB': - inner_wt = (1.0 / (2*n + 1)**2) - outer_wt = (1.0 / (12*n**2 + 4*n)) - dob_filter = np.zeros((4 * n + 1, 4 * n + 1)) - dob_filter[:] = outer_wt - dob_filter[n : 3 * n + 1, n : 3 * n + 1] = - inner_wt - return convolve(image, dob_filter) - elif mode == 'Octagon': outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] - return convolve(image, _octagon_filter(outer_shape[n - 1][0], outer_shape[n - 1][1], inner_shape[n - 1][0], inner_shape[n - 1][1])) -""" + # Take these out of the loop. No need to compute again for different scales. + integral_img = integral_image(image) + integral_img1 = _slanted_integral_image_modes(image, 1) + integral_img2 = _slanted_integral_image_modes(image, 2) + integral_img3 = _slanted_integral_image_modes(image, 3) + integral_img4 = _slanted_integral_image_modes(image, 4) + filtered_image = np.zeros(image.shape) + mo = outer_shape[n - 1][0] + no = outer_shape[n - 1][1] + mi = inner_shape[n - 1][0] + ni = inner_shape[n - 1][1] + outer_pixels = (mo + 2 * no)**2 - 2 * no * (no + 1) + inner_pixels = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) + outer_wt = 1.0 / (outer_pixels - inner_pixels) + inner_wt = 1.0 / inner_pixels + o_m = (mo - 1) / 2 + i_m = (mi - 1) / 2 + o_set = o_m + no + i_set = i_m + ni + # Outsource to Cython + for i in range(o_set + 1, image.shape[0] - o_set - 1): + for j in range(o_set + 1, image.shape[1] - o_set - 1): + outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m, j + o_set] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m, j - o_m] + outer += integral_img[i + (mo - 3) / 2, j + (mo - 3) / 2] - integral_img[i - o_m, j + (mo - 3) / 2] - integral_img[i + (mo - 3) / 2, j - o_m] + integral_img[i - o_m, j - o_m] + outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set, j - o_m] - integral_img[i - o_m, j - (mo + 1) / 2] + integral_img[i - o_m, j - o_set - 1] + outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m, j - o_set] - integral_img[i - (mo + 1) / 2, -1] - integral_img[i - o_set - 1, j + (mo - 3) / 2] + integral_img[i - (mo + 1) / 2, j + (mo - 3) / 2] + integral_img[i - o_set - 1, -1] + outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set, j + o_m] - integral_img[-1, j + o_set + 1] - integral_img[i + (mo - 3) / 2, j + o_m] + integral_img[-1, j + o_m] + integral_img[i + (mo - 3) / 2, j + o_set + 1] + + inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m, j + i_set] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m, j - i_m] + inner += integral_img[i + (mi - 3) / 2, j + (mi - 3) / 2] - integral_img[i - i_m, j + (mi - 3) / 2] - integral_img[i + (mi - 3) / 2, j - i_m] + integral_img[i - i_m, j - i_m] + inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set, j - i_m] - integral_img[i - i_m, j - (mi + 1) / 2] + integral_img[i - i_m, j - i_set - 1] + inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m, j - i_set] - integral_img[i - (mi + 1) / 2, -1] - integral_img[i - i_set - 1, j + (mi - 3) / 2] + integral_img[i - (mi + 1) / 2, j + (mi - 3) / 2] + integral_img[i - i_set - 1, -1] + inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set, j + i_m] - integral_img[-1, j + i_set + 1] - integral_img[i + (mi - 3) / 2, j + i_m] + integral_img[-1, j + i_m] + integral_img[i + (mi - 3) / 2, j + i_set + 1] + + filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner + return filtered_image +def _censure_octagon_loop(): + + +# Outsource to Cython + +def _slanted_integral_image(image): + flipped_lr = np.fliplr(image) + left_sum = np.zeros(image.shape[0]) + for i in range(image.shape[1] - image.shape[0], image.shape[1]): + left_sum[image.shape[1] - 1 - i] = np.sum(flipped_lr.diagonal(i)) + left_sum = left_sum.cumsum(0) + right_sum = np.sum(image, 1).cumsum(0) + image[:, 0] = left_sum + image[:, -1] = right_sum + integral_img = np.zeros((image.shape[0] + 1, image.shape[1])) + integral_img[1:, :] = image + for i in range(1, integral_img.shape[0]): + for j in range(1, integral_img.shape[1] - 1): + integral_img[i, j] += integral_img[i, j - 1] + integral_img[i - 1, j + 1] - integral_img[i - 1, j] + return integral_img[1:, :integral_img.shape[1]] + + +def _slanted_integral_image_modes(img, mode=1): + if mode == 1: + image = np.copy(img) + mode1 = _slanted_integral_image(image) + return mode1 + + elif mode == 2: + image = np.copy(img) + image = np.fliplr(image) + image = np.flipud(image) + mode2 = _slanted_integral_image(image) + mode2 = np.fliplr(mode2) + mode2 = np.flipud(mode2) + return mode2 + + elif mode == 3: + image = np.copy(img) + image = np.flipud(image) + image = image.T + mode3 = _slanted_integral_image(image) + mode3 = np.flipud(mode3.T) + return mode3 + + else: + image = np.copy(img) + image = np.fliplr(image) + image = image.T + mode4 = _slanted_integral_image(image) + mode4 = np.fliplr(mode4.T) + return mode4 + def _suppress_line(response, sigma, rpc_threshold): Axx, Axy, Ayy = _compute_auto_correlation(response, sigma) @@ -85,7 +132,7 @@ def censure_keypoints(image, mode='DoB', no_of_scales=7, threshold=0.03, rpc_thr # Generating all the scales scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) - for i in xrange(no_of_scales): + for i in range(no_of_scales): scales[:, :, i] = _get_filtered_image(image, i + 1, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 @@ -97,7 +144,7 @@ def censure_keypoints(image, mode='DoB', no_of_scales=7, threshold=0.03, rpc_thr maximas[np.abs(maximas) < threshold] = 0 response = maximas + minimas - for i in xrange(1, no_of_scales - 1): + for i in range(1, no_of_scales - 1): response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) # Returning keypoints with its scale From 474c7382ddcf42c070f5c6663fa9168df07f31a1 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 29 Jul 2013 09:24:38 +0530 Subject: [PATCH 15/77] Removing unused stuff --- skimage/feature/censure.py | 11 ++++------- 1 file changed, 4 insertions(+), 7 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 0b8d2b1b..c30da272 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -3,7 +3,6 @@ from scipy.ndimage.filters import maximum_filter, minimum_filter from ..transform import integral_image from ..feature.corner import _compute_auto_correlation -from ..morphology import convex_hull_image from ..util import img_as_float from .censure_cy import _censure_dob_loop, _slanted_integral_image @@ -57,11 +56,9 @@ def _get_filtered_image(image, n, mode='DoB'): filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner return filtered_image -def _censure_octagon_loop(): # Outsource to Cython - def _slanted_integral_image(image): flipped_lr = np.fliplr(image) left_sum = np.zeros(image.shape[0]) @@ -82,14 +79,14 @@ def _slanted_integral_image(image): def _slanted_integral_image_modes(img, mode=1): if mode == 1: image = np.copy(img) - mode1 = _slanted_integral_image(image) + mode1 = _slanted_integral_image(image, 1) return mode1 elif mode == 2: image = np.copy(img) image = np.fliplr(image) image = np.flipud(image) - mode2 = _slanted_integral_image(image) + mode2 = _slanted_integral_image(image, 2) mode2 = np.fliplr(mode2) mode2 = np.flipud(mode2) return mode2 @@ -98,7 +95,7 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img) image = np.flipud(image) image = image.T - mode3 = _slanted_integral_image(image) + mode3 = _slanted_integral_image(image, 3) mode3 = np.flipud(mode3.T) return mode3 @@ -106,7 +103,7 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img) image = np.fliplr(image) image = image.T - mode4 = _slanted_integral_image(image) + mode4 = _slanted_integral_image(image, 4) mode4 = np.fliplr(mode4.T) return mode4 From d15741c58c031ffb2171c96a70e46de606bfd5b1 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 29 Jul 2013 15:11:04 +0530 Subject: [PATCH 16/77] Cleaning up the non-Cython version for mode=Octagon --- skimage/feature/censure.py | 95 ++++++++++++++++++++------------------ 1 file changed, 51 insertions(+), 44 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index c30da272..6e1b87d1 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -5,57 +5,64 @@ from ..transform import integral_image from ..feature.corner import _compute_auto_correlation from ..util import img_as_float -from .censure_cy import _censure_dob_loop, _slanted_integral_image +from .censure_cy import _censure_dob_loop -def _get_filtered_image(image, n, mode='DoB'): +def _get_filtered_image(image, no_of_scales, mode): # TODO : Implement the STAR and Octagon mode if mode == 'DoB': - inner_wt = (1.0 / (2 * n + 1)**2) - outer_wt = (1.0 / (12 * n**2 + 4 * n)) - integral_img = integral_image(image) - filtered_image = np.zeros(image.shape) - _censure_dob_loop(image, n, integral_img, filtered_image, inner_wt, outer_wt) - return filtered_image + scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) + for i in range(no_of_scales): + n = i + 1 + inner_wt = (1.0 / (2 * n + 1)**2) + outer_wt = (1.0 / (12 * n**2 + 4 * n)) + integral_img = integral_image(image) + filtered_image = np.zeros(image.shape) + _censure_dob_loop(image, n, integral_img, filtered_image, inner_wt, outer_wt) + scales[:, :, i] = filtered_image + return scales elif mode == 'Octagon': outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] - # Take these out of the loop. No need to compute again for different scales. + scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) integral_img = integral_image(image) integral_img1 = _slanted_integral_image_modes(image, 1) integral_img2 = _slanted_integral_image_modes(image, 2) integral_img3 = _slanted_integral_image_modes(image, 3) integral_img4 = _slanted_integral_image_modes(image, 4) - filtered_image = np.zeros(image.shape) - mo = outer_shape[n - 1][0] - no = outer_shape[n - 1][1] - mi = inner_shape[n - 1][0] - ni = inner_shape[n - 1][1] - outer_pixels = (mo + 2 * no)**2 - 2 * no * (no + 1) - inner_pixels = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) - outer_wt = 1.0 / (outer_pixels - inner_pixels) - inner_wt = 1.0 / inner_pixels - o_m = (mo - 1) / 2 - i_m = (mi - 1) / 2 - o_set = o_m + no - i_set = i_m + ni - # Outsource to Cython - for i in range(o_set + 1, image.shape[0] - o_set - 1): - for j in range(o_set + 1, image.shape[1] - o_set - 1): - outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m, j + o_set] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m, j - o_m] - outer += integral_img[i + (mo - 3) / 2, j + (mo - 3) / 2] - integral_img[i - o_m, j + (mo - 3) / 2] - integral_img[i + (mo - 3) / 2, j - o_m] + integral_img[i - o_m, j - o_m] - outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set, j - o_m] - integral_img[i - o_m, j - (mo + 1) / 2] + integral_img[i - o_m, j - o_set - 1] - outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m, j - o_set] - integral_img[i - (mo + 1) / 2, -1] - integral_img[i - o_set - 1, j + (mo - 3) / 2] + integral_img[i - (mo + 1) / 2, j + (mo - 3) / 2] + integral_img[i - o_set - 1, -1] - outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set, j + o_m] - integral_img[-1, j + o_set + 1] - integral_img[i + (mo - 3) / 2, j + o_m] + integral_img[-1, j + o_m] + integral_img[i + (mo - 3) / 2, j + o_set + 1] + for k in range(no_of_scales): + n = k + 1 + filtered_image = np.zeros(image.shape) + mo = outer_shape[n - 1][0] + no = outer_shape[n - 1][1] + mi = inner_shape[n - 1][0] + ni = inner_shape[n - 1][1] + outer_pixels = (mo + 2 * no)**2 - 2 * no * (no + 1) + inner_pixels = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) + outer_wt = 1.0 / (outer_pixels - inner_pixels) + inner_wt = 1.0 / inner_pixels + o_m = (mo - 1) / 2 + i_m = (mi - 1) / 2 + o_set = o_m + no + i_set = i_m + ni + # Outsource to Cython + for i in range(o_set + 1, image.shape[0] - o_set - 1): + for j in range(o_set + 1, image.shape[1] - o_set - 1): + outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m, j + o_set] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m, j - o_m] + outer += integral_img[i + (mo - 3) / 2, j + (mo - 3) / 2] - integral_img[i - o_m, j + (mo - 3) / 2] - integral_img[i + (mo - 3) / 2, j - o_m] + integral_img[i - o_m, j - o_m] + outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set, j - o_m] - integral_img[i - o_m, j - (mo + 1) / 2] + integral_img[i - o_m, j - o_set - 1] + outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m, j - o_set] - integral_img[i - (mo + 1) / 2, -1] - integral_img[i - o_set - 1, j + (mo - 3) / 2] + integral_img[i - (mo + 1) / 2, j + (mo - 3) / 2] + integral_img[i - o_set - 1, -1] + outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set, j + o_m] - integral_img[-1, j + o_set + 1] - integral_img[i + (mo - 3) / 2, j + o_m] + integral_img[-1, j + o_m] + integral_img[i + (mo - 3) / 2, j + o_set + 1] - inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m, j + i_set] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m, j - i_m] - inner += integral_img[i + (mi - 3) / 2, j + (mi - 3) / 2] - integral_img[i - i_m, j + (mi - 3) / 2] - integral_img[i + (mi - 3) / 2, j - i_m] + integral_img[i - i_m, j - i_m] - inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set, j - i_m] - integral_img[i - i_m, j - (mi + 1) / 2] + integral_img[i - i_m, j - i_set - 1] - inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m, j - i_set] - integral_img[i - (mi + 1) / 2, -1] - integral_img[i - i_set - 1, j + (mi - 3) / 2] + integral_img[i - (mi + 1) / 2, j + (mi - 3) / 2] + integral_img[i - i_set - 1, -1] - inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set, j + i_m] - integral_img[-1, j + i_set + 1] - integral_img[i + (mi - 3) / 2, j + i_m] + integral_img[-1, j + i_m] + integral_img[i + (mi - 3) / 2, j + i_set + 1] + inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m, j + i_set] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m, j - i_m] + inner += integral_img[i + (mi - 3) / 2, j + (mi - 3) / 2] - integral_img[i - i_m, j + (mi - 3) / 2] - integral_img[i + (mi - 3) / 2, j - i_m] + integral_img[i - i_m, j - i_m] + inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set, j - i_m] - integral_img[i - i_m, j - (mi + 1) / 2] + integral_img[i - i_m, j - i_set - 1] + inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m, j - i_set] - integral_img[i - (mi + 1) / 2, -1] - integral_img[i - i_set - 1, j + (mi - 3) / 2] + integral_img[i - (mi + 1) / 2, j + (mi - 3) / 2] + integral_img[i - i_set - 1, -1] + inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set, j + i_m] - integral_img[-1, j + i_set + 1] - integral_img[i + (mi - 3) / 2, j + i_m] + integral_img[-1, j + i_m] + integral_img[i + (mi - 3) / 2, j + i_set + 1] - filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner - return filtered_image + filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner + scales[:, :, k] = filtered_image + return scales # Outsource to Cython @@ -79,14 +86,14 @@ def _slanted_integral_image(image): def _slanted_integral_image_modes(img, mode=1): if mode == 1: image = np.copy(img) - mode1 = _slanted_integral_image(image, 1) + mode1 = _slanted_integral_image(image) return mode1 elif mode == 2: image = np.copy(img) image = np.fliplr(image) image = np.flipud(image) - mode2 = _slanted_integral_image(image, 2) + mode2 = _slanted_integral_image(image) mode2 = np.fliplr(mode2) mode2 = np.flipud(mode2) return mode2 @@ -95,7 +102,7 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img) image = np.flipud(image) image = image.T - mode3 = _slanted_integral_image(image, 3) + mode3 = _slanted_integral_image(image) mode3 = np.flipud(mode3.T) return mode3 @@ -103,7 +110,7 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img) image = np.fliplr(image) image = image.T - mode4 = _slanted_integral_image(image, 4) + mode4 = _slanted_integral_image(image) mode4 = np.fliplr(mode4.T) return mode4 @@ -118,7 +125,7 @@ def _suppress_line(response, sigma, rpc_threshold): return response -def censure_keypoints(image, mode='DoB', no_of_scales=7, threshold=0.03, rpc_threshold=10): +def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_threshold=10): # TODO : Decide number of scales. Image-size dependent? image = np.squeeze(image) if image.ndim != 2: @@ -129,13 +136,13 @@ def censure_keypoints(image, mode='DoB', no_of_scales=7, threshold=0.03, rpc_thr # Generating all the scales scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) - for i in range(no_of_scales): - scales[:, :, i] = _get_filtered_image(image, i + 1, mode) + scales = _get_filtered_image(image, no_of_scales, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero minimas = (minimum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales maximas = (maximum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales + # Suppressing minimas and maximas weaker than threshold minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 From 7c32dc0c15ad6a85c09e68a6babed1d386984042 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 29 Jul 2013 17:16:15 +0530 Subject: [PATCH 17/77] Cythonizing the for loops --- skimage/feature/censure.py | 52 ++++++----------------------- skimage/feature/censure_cy.pyx | 60 ++++++++++++++++++++++++++++++++++ 2 files changed, 70 insertions(+), 42 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 6e1b87d1..47269a15 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -5,11 +5,11 @@ from ..transform import integral_image from ..feature.corner import _compute_auto_correlation from ..util import img_as_float -from .censure_cy import _censure_dob_loop - +from .censure_cy import _censure_dob_loop, _slanted_integral_image, _censure_octagon_loop +from time import time def _get_filtered_image(image, no_of_scales, mode): - # TODO : Implement the STAR and Octagon mode + # TODO : Implement the STAR mode if mode == 'DoB': scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) for i in range(no_of_scales): @@ -22,6 +22,7 @@ def _get_filtered_image(image, no_of_scales, mode): scales[:, :, i] = filtered_image return scales elif mode == 'Octagon': + # TODO : Decide the shapes of Octagon filters for scales > 7 outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) @@ -41,48 +42,13 @@ def _get_filtered_image(image, no_of_scales, mode): inner_pixels = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) outer_wt = 1.0 / (outer_pixels - inner_pixels) inner_wt = 1.0 / inner_pixels - o_m = (mo - 1) / 2 - i_m = (mi - 1) / 2 - o_set = o_m + no - i_set = i_m + ni - # Outsource to Cython - for i in range(o_set + 1, image.shape[0] - o_set - 1): - for j in range(o_set + 1, image.shape[1] - o_set - 1): - outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m, j + o_set] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m, j - o_m] - outer += integral_img[i + (mo - 3) / 2, j + (mo - 3) / 2] - integral_img[i - o_m, j + (mo - 3) / 2] - integral_img[i + (mo - 3) / 2, j - o_m] + integral_img[i - o_m, j - o_m] - outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set, j - o_m] - integral_img[i - o_m, j - (mo + 1) / 2] + integral_img[i - o_m, j - o_set - 1] - outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m, j - o_set] - integral_img[i - (mo + 1) / 2, -1] - integral_img[i - o_set - 1, j + (mo - 3) / 2] + integral_img[i - (mo + 1) / 2, j + (mo - 3) / 2] + integral_img[i - o_set - 1, -1] - outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set, j + o_m] - integral_img[-1, j + o_set + 1] - integral_img[i + (mo - 3) / 2, j + o_m] + integral_img[-1, j + o_m] + integral_img[i + (mo - 3) / 2, j + o_set + 1] - inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m, j + i_set] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m, j - i_m] - inner += integral_img[i + (mi - 3) / 2, j + (mi - 3) / 2] - integral_img[i - i_m, j + (mi - 3) / 2] - integral_img[i + (mi - 3) / 2, j - i_m] + integral_img[i - i_m, j - i_m] - inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set, j - i_m] - integral_img[i - i_m, j - (mi + 1) / 2] + integral_img[i - i_m, j - i_set - 1] - inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m, j - i_set] - integral_img[i - (mi + 1) / 2, -1] - integral_img[i - i_set - 1, j + (mi - 3) / 2] + integral_img[i - (mi + 1) / 2, j + (mi - 3) / 2] + integral_img[i - i_set - 1, -1] - inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set, j + i_m] - integral_img[-1, j + i_set + 1] - integral_img[i + (mi - 3) / 2, j + i_m] + integral_img[-1, j + i_m] + integral_img[i + (mi - 3) / 2, j + i_set + 1] + _censure_octagon_loop(image, integral_img, integral_img1, integral_img2, integral_img3, integral_img4, filtered_image, outer_wt, inner_wt, mo, no, mi, ni) - filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner scales[:, :, k] = filtered_image return scales -# Outsource to Cython -def _slanted_integral_image(image): - flipped_lr = np.fliplr(image) - left_sum = np.zeros(image.shape[0]) - for i in range(image.shape[1] - image.shape[0], image.shape[1]): - left_sum[image.shape[1] - 1 - i] = np.sum(flipped_lr.diagonal(i)) - left_sum = left_sum.cumsum(0) - right_sum = np.sum(image, 1).cumsum(0) - image[:, 0] = left_sum - image[:, -1] = right_sum - integral_img = np.zeros((image.shape[0] + 1, image.shape[1])) - integral_img[1:, :] = image - for i in range(1, integral_img.shape[0]): - for j in range(1, integral_img.shape[1] - 1): - integral_img[i, j] += integral_img[i, j - 1] + integral_img[i - 1, j + 1] - integral_img[i - 1, j] - return integral_img[1:, :integral_img.shape[1]] - - def _slanted_integral_image_modes(img, mode=1): if mode == 1: image = np.copy(img) @@ -135,22 +101,24 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr image = np.ascontiguousarray(image) # Generating all the scales + start = time() scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) scales = _get_filtered_image(image, no_of_scales, mode) + print time() - start # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero minimas = (minimum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales maximas = (maximum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales - + print time() - start # Suppressing minimas and maximas weaker than threshold minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 response = maximas + minimas - + print time() - start for i in range(1, no_of_scales - 1): response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) - + print time() - start # Returning keypoints with its scale keypoints = np.transpose(np.nonzero(response[:, :, 1:no_of_scales])) + [0, 0, 1] return keypoints diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index d5e8bf02..104d3348 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -4,6 +4,7 @@ #cython: wraparound=False cimport numpy as cnp +import numpy as np def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, @@ -19,3 +20,62 @@ def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, inner = integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n] outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n - 1, j - 2 * n - 1] - integral_img[i + 2 * n, j - 2 * n - 1] - integral_img[i - 2 * n - 1, j + 2 * n] filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner + + +def _slanted_integral_image(double[:, :] image): + + cdef Py_ssize_t i, j + cdef double[:, :] flipped_lr = np.fliplr(image) + cdef double[:] left_sum = np.zeros(image.shape[0], dtype=np.float) + cdef double[:] right_sum + cdef double[:, :] integral_img = np.zeros((image.shape[0] + 1, image.shape[1]), dtype=np.float) + + flipped_lr = np.fliplr(image) + for i in range(image.shape[1] - image.shape[0], image.shape[1]): + left_sum[image.shape[1] - 1 - i] = np.sum(flipped_lr.diagonal(i)) + left_sum = left_sum.cumsum(0) + + right_sum = np.sum(image, 1).cumsum(0) + + image[:, 0] = left_sum + image[:, -1] = right_sum + + integral_img[1:, :] = image + + for i in range(1, integral_img.shape[0]): + for j in range(1, integral_img.shape[1] - 1): + integral_img[i, j] += integral_img[i, j - 1] + integral_img[i - 1, j + 1] - integral_img[i - 1, j] + return integral_img[1:, :integral_img.shape[1]] + + +def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, + double[:, ::1] integral_img1, + double[:, ::1] integral_img2, + double[:, ::1] integral_img3, + double[:, ::1] integral_img4, + double[:, ::1] filtered_image, + double outer_wt, double inner_wt, + int mo, int no, int mi, int ni): + + cdef Py_ssize_t i, j, o_m, i_m, o_set, i_set + + o_m = (mo - 1) / 2 + i_m = (mi - 1) / 2 + o_set = o_m + no + i_set = i_m + ni + + for i in range(o_set + 1, image.shape[0] - o_set - 1): + for j in range(o_set + 1, image.shape[1] - o_set - 1): + outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m, j + o_set] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m, j - o_m] + outer += integral_img[i + o_m - 1, j + o_m - 1] - integral_img[i - o_m, j + o_m - 1] - integral_img[i + o_m - 1, j - o_m] + integral_img[i - o_m, j - o_m] + outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set, j - o_m] - integral_img[i - o_m, j - o_m + 1] + integral_img[i - o_m, j - o_set - 1] + outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m, j - o_set] - integral_img[i - o_m + 1, -1] - integral_img[i - o_set - 1, j + o_m - 1] + integral_img[i - o_m + 1, j + o_m - 1] + integral_img[i - o_set - 1, -1] + outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set, j + o_m] - integral_img[-1, j + o_set + 1] - integral_img[i + o_m - 1, j + o_m] + integral_img[-1, j + o_m] + integral_img[i + o_m - 1, j + o_set + 1] + + inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m, j + i_set] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m, j - i_m] + inner += integral_img[i + i_m - 1, j + i_m - 1] - integral_img[i - i_m, j + i_m - 1] - integral_img[i + i_m - 1, j - i_m] + integral_img[i - i_m, j - i_m] + inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set, j - i_m] - integral_img[i - i_m, j - i_m + 1] + integral_img[i - i_m, j - i_set - 1] + inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m, j - i_set] - integral_img[i - i_m + 1, -1] - integral_img[i - i_set - 1, j + i_m - 1] + integral_img[i - i_m + 1, j + i_m - 1] + integral_img[i - i_set - 1, -1] + inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set, j + i_m] - integral_img[-1, j + i_set + 1] - integral_img[i + i_m - 1, j + i_m] + integral_img[-1, j + i_m] + integral_img[i + i_m - 1, j + i_set + 1] + + filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner From 5ecbb57e829b51e1aa8cf1de68be67f8de712420 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 29 Jul 2013 22:57:13 +0530 Subject: [PATCH 18/77] Adding the docs for censure_keypoints --- skimage/feature/censure.py | 45 +++++++++++++++++++++++++++++++++++++- 1 file changed, 44 insertions(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 47269a15..618a29b9 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -92,7 +92,50 @@ def _suppress_line(response, sigma, rpc_threshold): def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_threshold=10): - # TODO : Decide number of scales. Image-size dependent? + """ + Extracts Censure keypoints along with the corresponding scale using + either Difference of Boxes, Octagon or STAR bilevel filter. + + Parameters + ---------- + image : 2D ndarray + Input image. + + no_of_scales : positive integer + Number of scales to extract keypoints from. Default is 7. + + mode : 'DoB' + Type of bilevel filter used to get the scales of input image. Possible + values are 'DoB', 'Octagon' and 'STAR'. Default is 'DoB'. + + threshold : + Threshold value used to suppress maximas and minimas with a weak + magnitude response obtained after Non-Maximal Suppression. Default + is 0.03. + + rpc_threshold : + Threshold for rejecting interest points which have ratio of principal + curvatures greater than this value. Default is 10. + + Returns + ------- + keypoints : (N, 3) array + Location of extracted keypoints along with the corresponding scale. + + References + ---------- + .. [1] Motilal Agrawal, Kurt Konolige and Morten Rufus Blas + "CenSurE: Center Surround Extremas for Realtime Feature + Detection and Matching", + http://link.springer.com/content/pdf/10.1007%2F978-3-540-88693-8_8.pdf + + .. [2] Adam Schmidt, Marek Kraft, Michal Fularz and Zuzanna Domagala + "Comparative Assessment of Point Feature Detectors and + Descriptors in the Context of Robot Navigation" + http://www.jamris.org/01_2013/saveas.php?QUEST=JAMRIS_No01_2013_P_11-20.pdf + + """ + image = np.squeeze(image) if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") From 6243d4bc0944a034345c63f84376cf29e0c4aac4 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 29 Jul 2013 22:58:30 +0530 Subject: [PATCH 19/77] Trying to resolve Cython TypeErrors --- skimage/feature/censure_cy.pyx | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 104d3348..4491363d 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -22,23 +22,23 @@ def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner -def _slanted_integral_image(double[:, :] image): +def _slanted_integral_image(double[:, ::1] image): cdef Py_ssize_t i, j - cdef double[:, :] flipped_lr = np.fliplr(image) cdef double[:] left_sum = np.zeros(image.shape[0], dtype=np.float) - cdef double[:] right_sum cdef double[:, :] integral_img = np.zeros((image.shape[0] + 1, image.shape[1]), dtype=np.float) - flipped_lr = np.fliplr(image) + flipped_lr = np.asarray(image[:, ::-1]) for i in range(image.shape[1] - image.shape[0], image.shape[1]): left_sum[image.shape[1] - 1 - i] = np.sum(flipped_lr.diagonal(i)) - left_sum = left_sum.cumsum(0) + left_sum_np = np.asarray(left_sum) + #image = np.asarray(image) + left_sum_np = left_sum_np.cumsum(0) + #right_sum_np = np.asarray() + right_sum_np = np.sum(image, 1).cumsum(0) - right_sum = np.sum(image, 1).cumsum(0) - - image[:, 0] = left_sum - image[:, -1] = right_sum + image[:, 0] = left_sum_np + image[:, -1] = right_sum_np integral_img[1:, :] = image From 4a19b60721f3be1fb2d35ac638e1deab8598eb05 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 30 Jul 2013 02:32:24 +0530 Subject: [PATCH 20/77] Trying to Debug the SegFault --- skimage/feature/censure.py | 27 ++++++++++++++++----------- skimage/feature/censure_cy.pyx | 26 +++++++++++++++++--------- 2 files changed, 33 insertions(+), 20 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 618a29b9..f7d75680 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -6,7 +6,7 @@ from ..feature.corner import _compute_auto_correlation from ..util import img_as_float from .censure_cy import _censure_dob_loop, _slanted_integral_image, _censure_octagon_loop -from time import time + def _get_filtered_image(image, no_of_scales, mode): # TODO : Implement the STAR mode @@ -52,14 +52,17 @@ def _get_filtered_image(image, no_of_scales, mode): def _slanted_integral_image_modes(img, mode=1): if mode == 1: image = np.copy(img) - mode1 = _slanted_integral_image(image) - return mode1 + mode1 = np.zeros((image.shape[0] + 1, image.shape[1])) + _slanted_integral_image(image, mode1) + return mode1[1:, :mode1.shape[1]] elif mode == 2: image = np.copy(img) image = np.fliplr(image) image = np.flipud(image) - mode2 = _slanted_integral_image(image) + mode2 = np.zeros((image.shape[0] + 1, image.shape[1])) + _slanted_integral_image(image, mode2) + mode2 = mode2[1:, :mode2.shape[1]] mode2 = np.fliplr(mode2) mode2 = np.flipud(mode2) return mode2 @@ -68,7 +71,9 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img) image = np.flipud(image) image = image.T - mode3 = _slanted_integral_image(image) + mode3 = np.zeros((image.shape[0] + 1, image.shape[1])) + _slanted_integral_image(image, mode3) + mode3 = mode3[1:, :mode3.shape[1]] mode3 = np.flipud(mode3.T) return mode3 @@ -76,7 +81,9 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img) image = np.fliplr(image) image = image.T - mode4 = _slanted_integral_image(image) + mode4 = np.zeros((image.shape[0] + 1, image.shape[1])) + _slanted_integral_image(image, mode4) + mode4 = mode4[1:, :mode4.shape[1]] mode4 = np.fliplr(mode4.T) return mode4 @@ -144,24 +151,22 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr image = np.ascontiguousarray(image) # Generating all the scales - start = time() scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) scales = _get_filtered_image(image, no_of_scales, mode) - print time() - start # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero minimas = (minimum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales maximas = (maximum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales - print time() - start + # Suppressing minimas and maximas weaker than threshold minimas[np.abs(minimas) < threshold] = 0 maximas[np.abs(maximas) < threshold] = 0 response = maximas + minimas - print time() - start + for i in range(1, no_of_scales - 1): response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) - print time() - start + # Returning keypoints with its scale keypoints = np.transpose(np.nonzero(response[:, :, 1:no_of_scales])) + [0, 0, 1] return keypoints diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 4491363d..9783d201 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -22,30 +22,37 @@ def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner -def _slanted_integral_image(double[:, ::1] image): +def _slanted_integral_image(double[:, ::1] image, + double[:, ::1] integral_img): cdef Py_ssize_t i, j cdef double[:] left_sum = np.zeros(image.shape[0], dtype=np.float) - cdef double[:, :] integral_img = np.zeros((image.shape[0] + 1, image.shape[1]), dtype=np.float) flipped_lr = np.asarray(image[:, ::-1]) for i in range(image.shape[1] - image.shape[0], image.shape[1]): left_sum[image.shape[1] - 1 - i] = np.sum(flipped_lr.diagonal(i)) left_sum_np = np.asarray(left_sum) - #image = np.asarray(image) + left_sum_np = left_sum_np.cumsum(0) - #right_sum_np = np.asarray() + right_sum_np = np.sum(image, 1).cumsum(0) - image[:, 0] = left_sum_np - image[:, -1] = right_sum_np + print '1' + for i in range(image.shape[0]): + image[i, 0] = left_sum_np[i] + image[i, -1] = right_sum_np[i] - integral_img[1:, :] = image + print '2' + for i in range(1, integral_img.shape[0]): + for j in range(integral_img.shape[1]): + integral_img[i, j] = image[i - 1, j] + print '3' for i in range(1, integral_img.shape[0]): for j in range(1, integral_img.shape[1] - 1): integral_img[i, j] += integral_img[i, j - 1] + integral_img[i - 1, j + 1] - integral_img[i - 1, j] - return integral_img[1:, :integral_img.shape[1]] + print '4' + def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, @@ -63,7 +70,7 @@ def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, i_m = (mi - 1) / 2 o_set = o_m + no i_set = i_m + ni - + print '5' for i in range(o_set + 1, image.shape[0] - o_set - 1): for j in range(o_set + 1, image.shape[1] - o_set - 1): outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m, j + o_set] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m, j - o_m] @@ -79,3 +86,4 @@ def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set, j + i_m] - integral_img[-1, j + i_set + 1] - integral_img[i + i_m - 1, j + i_m] + integral_img[-1, j + i_m] + integral_img[i + i_m - 1, j + i_set + 1] filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner + print '6' From f75fbb986be6bb7ef5a56064d39555e4ce4caaff Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 30 Jul 2013 02:54:50 +0530 Subject: [PATCH 21/77] Trying to Debug SegFault : 2 --- skimage/feature/censure.py | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index f7d75680..86bafc4b 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -28,9 +28,13 @@ def _get_filtered_image(image, no_of_scales, mode): scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) integral_img = integral_image(image) integral_img1 = _slanted_integral_image_modes(image, 1) + print '8' integral_img2 = _slanted_integral_image_modes(image, 2) + print '9' integral_img3 = _slanted_integral_image_modes(image, 3) + print '10' integral_img4 = _slanted_integral_image_modes(image, 4) + print '11' for k in range(no_of_scales): n = k + 1 filtered_image = np.zeros(image.shape) @@ -54,7 +58,8 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img) mode1 = np.zeros((image.shape[0] + 1, image.shape[1])) _slanted_integral_image(image, mode1) - return mode1[1:, :mode1.shape[1]] + print '7' + return mode1[1:, :] elif mode == 2: image = np.copy(img) @@ -62,7 +67,8 @@ def _slanted_integral_image_modes(img, mode=1): image = np.flipud(image) mode2 = np.zeros((image.shape[0] + 1, image.shape[1])) _slanted_integral_image(image, mode2) - mode2 = mode2[1:, :mode2.shape[1]] + print '7' + mode2 = mode2[1:, :] mode2 = np.fliplr(mode2) mode2 = np.flipud(mode2) return mode2 @@ -73,7 +79,8 @@ def _slanted_integral_image_modes(img, mode=1): image = image.T mode3 = np.zeros((image.shape[0] + 1, image.shape[1])) _slanted_integral_image(image, mode3) - mode3 = mode3[1:, :mode3.shape[1]] + print '7' + mode3 = mode3[1:, :] mode3 = np.flipud(mode3.T) return mode3 @@ -83,7 +90,8 @@ def _slanted_integral_image_modes(img, mode=1): image = image.T mode4 = np.zeros((image.shape[0] + 1, image.shape[1])) _slanted_integral_image(image, mode4) - mode4 = mode4[1:, :mode4.shape[1]] + print '7' + mode4 = mode4[1:, :] mode4 = np.fliplr(mode4.T) return mode4 From 345c1690cca2213398e262e280132becce6e7873 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 30 Jul 2013 22:14:23 +0530 Subject: [PATCH 22/77] Trying to debug the Segfault-3 : Making the arrays C-contiguous --- skimage/feature/censure.py | 16 ++++++++-------- skimage/feature/censure_cy.pyx | 6 +++--- 2 files changed, 11 insertions(+), 11 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 86bafc4b..36f57836 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -55,17 +55,17 @@ def _get_filtered_image(image, no_of_scales, mode): def _slanted_integral_image_modes(img, mode=1): if mode == 1: - image = np.copy(img) - mode1 = np.zeros((image.shape[0] + 1, image.shape[1])) + image = np.copy(img, order='C') + mode1 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode1) print '7' return mode1[1:, :] elif mode == 2: - image = np.copy(img) + image = np.copy(img, order='C') image = np.fliplr(image) image = np.flipud(image) - mode2 = np.zeros((image.shape[0] + 1, image.shape[1])) + mode2 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode2) print '7' mode2 = mode2[1:, :] @@ -74,10 +74,10 @@ def _slanted_integral_image_modes(img, mode=1): return mode2 elif mode == 3: - image = np.copy(img) + image = np.copy(img, order='C') image = np.flipud(image) image = image.T - mode3 = np.zeros((image.shape[0] + 1, image.shape[1])) + mode3 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode3) print '7' mode3 = mode3[1:, :] @@ -85,10 +85,10 @@ def _slanted_integral_image_modes(img, mode=1): return mode3 else: - image = np.copy(img) + image = np.copy(img, order='C') image = np.fliplr(image) image = image.T - mode4 = np.zeros((image.shape[0] + 1, image.shape[1])) + mode4 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode4) print '7' mode4 = mode4[1:, :] diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 9783d201..d8cc80b5 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -1,7 +1,7 @@ #cython: cdivision=True -#cython: boundscheck=False -#cython: nonecheck=False -#cython: wraparound=False +#cython: boundscheck=True +#cython: nonecheck=True +#cython: wraparound=True cimport numpy as cnp import numpy as np From efba964dde645619aaacf7e8efcbe4044dda44c5 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 05:34:03 +0530 Subject: [PATCH 23/77] Fixing many bugs; Working Censure for mode=Octagon --- skimage/feature/censure.py | 16 +++++++--------- skimage/feature/censure_cy.pyx | 24 ++++++++++-------------- 2 files changed, 17 insertions(+), 23 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 36f57836..9997ad50 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -28,13 +28,12 @@ def _get_filtered_image(image, no_of_scales, mode): scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) integral_img = integral_image(image) integral_img1 = _slanted_integral_image_modes(image, 1) - print '8' integral_img2 = _slanted_integral_image_modes(image, 2) - print '9' + integral_img2 = np.ascontiguousarray(integral_img2) integral_img3 = _slanted_integral_image_modes(image, 3) - print '10' + integral_img3 = np.ascontiguousarray(integral_img3) integral_img4 = _slanted_integral_image_modes(image, 4) - print '11' + integral_img4 = np.ascontiguousarray(integral_img4) for k in range(no_of_scales): n = k + 1 filtered_image = np.zeros(image.shape) @@ -58,16 +57,15 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img, order='C') mode1 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode1) - print '7' return mode1[1:, :] elif mode == 2: image = np.copy(img, order='C') image = np.fliplr(image) image = np.flipud(image) + image = np.ascontiguousarray(image) mode2 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode2) - print '7' mode2 = mode2[1:, :] mode2 = np.fliplr(mode2) mode2 = np.flipud(mode2) @@ -77,9 +75,9 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img, order='C') image = np.flipud(image) image = image.T + image = np.ascontiguousarray(image) mode3 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode3) - print '7' mode3 = mode3[1:, :] mode3 = np.flipud(mode3.T) return mode3 @@ -88,9 +86,9 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img, order='C') image = np.fliplr(image) image = image.T + image = np.ascontiguousarray(image) mode4 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode4) - print '7' mode4 = mode4[1:, :] mode4 = np.fliplr(mode4.T) return mode4 @@ -176,5 +174,5 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) # Returning keypoints with its scale - keypoints = np.transpose(np.nonzero(response[:, :, 1:no_of_scales])) + [0, 0, 1] + keypoints = np.transpose(np.nonzero(response[:, :, 1:no_of_scales - 1])) + [0, 0, 2] return keypoints diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index d8cc80b5..29c91378 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -37,21 +37,18 @@ def _slanted_integral_image(double[:, ::1] image, right_sum_np = np.sum(image, 1).cumsum(0) - print '1' for i in range(image.shape[0]): image[i, 0] = left_sum_np[i] image[i, -1] = right_sum_np[i] - print '2' for i in range(1, integral_img.shape[0]): for j in range(integral_img.shape[1]): integral_img[i, j] = image[i - 1, j] - print '3' for i in range(1, integral_img.shape[0]): for j in range(1, integral_img.shape[1] - 1): integral_img[i, j] += integral_img[i, j - 1] + integral_img[i - 1, j + 1] - integral_img[i - 1, j] - print '4' + @@ -70,20 +67,19 @@ def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, i_m = (mi - 1) / 2 o_set = o_m + no i_set = i_m + ni - print '5' + for i in range(o_set + 1, image.shape[0] - o_set - 1): for j in range(o_set + 1, image.shape[1] - o_set - 1): - outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m, j + o_set] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m, j - o_m] + outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m - 1, j + o_set + 1] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m - 1, j - o_m] outer += integral_img[i + o_m - 1, j + o_m - 1] - integral_img[i - o_m, j + o_m - 1] - integral_img[i + o_m - 1, j - o_m] + integral_img[i - o_m, j - o_m] - outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set, j - o_m] - integral_img[i - o_m, j - o_m + 1] + integral_img[i - o_m, j - o_set - 1] - outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m, j - o_set] - integral_img[i - o_m + 1, -1] - integral_img[i - o_set - 1, j + o_m - 1] + integral_img[i - o_m + 1, j + o_m - 1] + integral_img[i - o_set - 1, -1] - outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set, j + o_m] - integral_img[-1, j + o_set + 1] - integral_img[i + o_m - 1, j + o_m] + integral_img[-1, j + o_m] + integral_img[i + o_m - 1, j + o_set + 1] + outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set + 1, j - o_m + 1] - integral_img[i - o_m, j - o_m] + integral_img[i - o_m, j - o_set - 1] + outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m + 1, j - o_set - 1] - integral_img[i - o_m, -1] - integral_img[i - o_set - 1, j + o_m - 1] + integral_img[i - o_m, j + o_m - 1] + integral_img[i - o_set - 1, -1] + outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set - 1, j + o_m - 1] - integral_img[-1, j + o_set] - integral_img[i + o_m - 1, j + o_m - 1] + integral_img[-1, j + o_m - 1] + integral_img[i + o_m - 1, j + o_set] - inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m, j + i_set] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m, j - i_m] + inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m - 1, j + i_set + 1] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m - 1, j - i_m] inner += integral_img[i + i_m - 1, j + i_m - 1] - integral_img[i - i_m, j + i_m - 1] - integral_img[i + i_m - 1, j - i_m] + integral_img[i - i_m, j - i_m] - inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set, j - i_m] - integral_img[i - i_m, j - i_m + 1] + integral_img[i - i_m, j - i_set - 1] - inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m, j - i_set] - integral_img[i - i_m + 1, -1] - integral_img[i - i_set - 1, j + i_m - 1] + integral_img[i - i_m + 1, j + i_m - 1] + integral_img[i - i_set - 1, -1] - inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set, j + i_m] - integral_img[-1, j + i_set + 1] - integral_img[i + i_m - 1, j + i_m] + integral_img[-1, j + i_m] + integral_img[i + i_m - 1, j + i_set + 1] + inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set + 1, j - i_m + 1] - integral_img[i - i_m, j - i_m] + integral_img[i - i_m, j - i_set - 1] + inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m + 1, j - i_set - 1] - integral_img[i - i_m, -1] - integral_img[i - i_set - 1, j + i_m - 1] + integral_img[i - i_m, j + i_m - 1] + integral_img[i - i_set - 1, -1] + inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set - 1, j + i_m - 1] - integral_img[-1, j + i_set] - integral_img[i + i_m - 1, j + i_m - 1] + integral_img[-1, j + i_m - 1] + integral_img[i + i_m - 1, j + i_set] filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner - print '6' From 50a18383df275d3f131b275baccd8e8f63c58fac Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 12:40:49 +0530 Subject: [PATCH 24/77] Adding doc for the different modes of _slanted_integral_inage_modes --- skimage/feature/censure.py | 42 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 42 insertions(+) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 9997ad50..cea2b4c3 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -54,12 +54,34 @@ def _get_filtered_image(image, no_of_scales, mode): def _slanted_integral_image_modes(img, mode=1): if mode == 1: + """ + The following figures describe area that is summed up to calculate + the value at point @ in slanted integral image. + _________________ + |********/ | + |*******/ | + |******/ | + |-----@ | + | | + | | + |_________________| + """ image = np.copy(img, order='C') mode1 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode1) return mode1[1:, :] elif mode == 2: + """ + _________________ + | | + | | + | | + | @_____| + | /******| + | /*******| + |________/________| + """ image = np.copy(img, order='C') image = np.fliplr(image) image = np.flipud(image) @@ -72,6 +94,16 @@ def _slanted_integral_image_modes(img, mode=1): return mode2 elif mode == 3: + """ + _________________ + | | + |\\ | + |*\\ | + |**\\ | + |***@ | + |***| | + |___|_____________| + """ image = np.copy(img, order='C') image = np.flipud(image) image = image.T @@ -83,6 +115,16 @@ def _slanted_integral_image_modes(img, mode=1): return mode3 else: + """ + ________________ + | |****| + | |****| + | @****| + | \\**| + | \\*| + | \\| + |________________| + """ image = np.copy(img, order='C') image = np.fliplr(image) image = image.T From 8c38b5c031e9ff765e594dc9fd12f9b2f82b3120 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 12:50:52 +0530 Subject: [PATCH 25/77] Correcting the docs and making them explicit --- skimage/feature/censure.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index cea2b4c3..9a77dfcd 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -157,20 +157,20 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr Input image. no_of_scales : positive integer - Number of scales to extract keypoints from. Default is 7. + Number of scales to extract keypoints from. The keypoints will be + extracted from all the scales except the first and the last. - mode : 'DoB' + mode : {'DoB', 'Octagon', 'STAR'} Type of bilevel filter used to get the scales of input image. Possible - values are 'DoB', 'Octagon' and 'STAR'. Default is 'DoB'. + values are 'DoB', 'Octagon' and 'STAR'. - threshold : + threshold : float Threshold value used to suppress maximas and minimas with a weak - magnitude response obtained after Non-Maximal Suppression. Default - is 0.03. + magnitude response obtained after Non-Maximal Suppression. - rpc_threshold : + rpc_threshold : float Threshold for rejecting interest points which have ratio of principal - curvatures greater than this value. Default is 10. + curvatures greater than this value. Returns ------- From 3b9ee0094d59f2310d66dec79eadcdc973e9bfb3 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 12:55:04 +0530 Subject: [PATCH 26/77] Commenting the sigma passed to _suppress_line --- skimage/feature/censure.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 9a77dfcd..01316ba9 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -213,6 +213,9 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr response = maximas + minimas for i in range(1, no_of_scales - 1): + # sigma = (window_size - 1) / 6.0 + # window_size = 7 + 2 * i + # Hence sigma = 1 + i / 3.0 response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) # Returning keypoints with its scale From 7df7f3a1f0c969ec9099d44277d6f0133a85a57a Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 12:56:07 +0530 Subject: [PATCH 27/77] Removing unnecessary data-type conversion --- skimage/feature/censure.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 01316ba9..efb0df81 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -204,8 +204,8 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero - minimas = (minimum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales - maximas = (maximum_filter(scales, (3, 3, 3)) == scales).astype(int) * scales + minimas = (minimum_filter(scales, (3, 3, 3)) == scales) * scales + maximas = (maximum_filter(scales, (3, 3, 3)) == scales) * scales # Suppressing minimas and maximas weaker than threshold minimas[np.abs(minimas) < threshold] = 0 From 6c1424732b9384f518736575117843101be77d6c Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 12:58:07 +0530 Subject: [PATCH 28/77] Removing unnecessary array pre-allocation --- skimage/feature/censure.py | 1 - 1 file changed, 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index efb0df81..0b99b529 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -199,7 +199,6 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr image = np.ascontiguousarray(image) # Generating all the scales - scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) scales = _get_filtered_image(image, no_of_scales, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 From ec7ea0007c97f602fed00e5661e63f6293abea01 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 12:59:35 +0530 Subject: [PATCH 29/77] Changing no_of_scales to n_scales --- skimage/feature/censure.py | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 0b99b529..b63c5207 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -8,11 +8,11 @@ from ..util import img_as_float from .censure_cy import _censure_dob_loop, _slanted_integral_image, _censure_octagon_loop -def _get_filtered_image(image, no_of_scales, mode): +def _get_filtered_image(image, n_scales, mode): # TODO : Implement the STAR mode if mode == 'DoB': - scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) - for i in range(no_of_scales): + scales = np.zeros((image.shape[0], image.shape[1], n_scales)) + for i in range(n_scales): n = i + 1 inner_wt = (1.0 / (2 * n + 1)**2) outer_wt = (1.0 / (12 * n**2 + 4 * n)) @@ -25,7 +25,7 @@ def _get_filtered_image(image, no_of_scales, mode): # TODO : Decide the shapes of Octagon filters for scales > 7 outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] - scales = np.zeros((image.shape[0], image.shape[1], no_of_scales)) + scales = np.zeros((image.shape[0], image.shape[1], n_scales)) integral_img = integral_image(image) integral_img1 = _slanted_integral_image_modes(image, 1) integral_img2 = _slanted_integral_image_modes(image, 2) @@ -34,7 +34,7 @@ def _get_filtered_image(image, no_of_scales, mode): integral_img3 = np.ascontiguousarray(integral_img3) integral_img4 = _slanted_integral_image_modes(image, 4) integral_img4 = np.ascontiguousarray(integral_img4) - for k in range(no_of_scales): + for k in range(n_scales): n = k + 1 filtered_image = np.zeros(image.shape) mo = outer_shape[n - 1][0] @@ -146,7 +146,7 @@ def _suppress_line(response, sigma, rpc_threshold): return response -def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_threshold=10): +def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, rpc_threshold=10): """ Extracts Censure keypoints along with the corresponding scale using either Difference of Boxes, Octagon or STAR bilevel filter. @@ -156,7 +156,7 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr image : 2D ndarray Input image. - no_of_scales : positive integer + n_scales : positive integer Number of scales to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last. @@ -199,7 +199,7 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr image = np.ascontiguousarray(image) # Generating all the scales - scales = _get_filtered_image(image, no_of_scales, mode) + scales = _get_filtered_image(image, n_scales, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero @@ -211,12 +211,12 @@ def censure_keypoints(image, no_of_scales=7, mode='DoB', threshold=0.03, rpc_thr maximas[np.abs(maximas) < threshold] = 0 response = maximas + minimas - for i in range(1, no_of_scales - 1): + for i in range(1, n_scales - 1): # sigma = (window_size - 1) / 6.0 # window_size = 7 + 2 * i # Hence sigma = 1 + i / 3.0 response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) # Returning keypoints with its scale - keypoints = np.transpose(np.nonzero(response[:, :, 1:no_of_scales - 1])) + [0, 0, 2] + keypoints = np.transpose(np.nonzero(response[:, :, 1:n_scales - 1])) + [0, 0, 2] return keypoints From 555da870fe82bb0714d4effda7c6d200879bcfac Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 13:18:04 +0530 Subject: [PATCH 30/77] Renaming *_wt to *_weight; correcting the object type of n in _censure_dob_loop --- skimage/feature/censure.py | 15 +++++++++------ skimage/feature/censure_cy.pyx | 10 +++++----- 2 files changed, 14 insertions(+), 11 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index b63c5207..1175235f 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -14,11 +14,14 @@ def _get_filtered_image(image, n_scales, mode): scales = np.zeros((image.shape[0], image.shape[1], n_scales)) for i in range(n_scales): n = i + 1 - inner_wt = (1.0 / (2 * n + 1)**2) - outer_wt = (1.0 / (12 * n**2 + 4 * n)) + # Constant multipliers for the outer region and the inner region + # of the bilevel filters with the constraint of keeping the DC bias + # 0. + inner_weight = (1.0 / (2 * n + 1)**2) + outer_weight = (1.0 / (12 * n**2 + 4 * n)) integral_img = integral_image(image) filtered_image = np.zeros(image.shape) - _censure_dob_loop(image, n, integral_img, filtered_image, inner_wt, outer_wt) + _censure_dob_loop(image, n, integral_img, filtered_image, inner_weight, outer_weight) scales[:, :, i] = filtered_image return scales elif mode == 'Octagon': @@ -43,10 +46,10 @@ def _get_filtered_image(image, n_scales, mode): ni = inner_shape[n - 1][1] outer_pixels = (mo + 2 * no)**2 - 2 * no * (no + 1) inner_pixels = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) - outer_wt = 1.0 / (outer_pixels - inner_pixels) - inner_wt = 1.0 / inner_pixels + outer_weight = 1.0 / (outer_pixels - inner_pixels) + inner_weight = 1.0 / inner_pixels - _censure_octagon_loop(image, integral_img, integral_img1, integral_img2, integral_img3, integral_img4, filtered_image, outer_wt, inner_wt, mo, no, mi, ni) + _censure_octagon_loop(image, integral_img, integral_img1, integral_img2, integral_img3, integral_img4, filtered_image, outer_weight, inner_weight, mo, no, mi, ni) scales[:, :, k] = filtered_image return scales diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 29c91378..75f58278 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -7,10 +7,10 @@ cimport numpy as cnp import numpy as np -def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, +def _censure_dob_loop(double[:, ::1] image, Py_ssize_t n, double[:, ::1] integral_img, double[:, ::1] filtered_image, - double inner_wt, double outer_wt): + double inner_weight, double outer_weight): cdef Py_ssize_t i, j cdef double inner, outer @@ -19,7 +19,7 @@ def _censure_dob_loop(double[:, ::1] image, cnp.int16_t n, for j in range(2 * n, image.shape[1] - 2 * n): inner = integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n] outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n - 1, j - 2 * n - 1] - integral_img[i + 2 * n, j - 2 * n - 1] - integral_img[i - 2 * n - 1, j + 2 * n] - filtered_image[i, j] = outer_wt * outer - (inner_wt + outer_wt) * inner + filtered_image[i, j] = outer_weight * outer - (inner_weight + outer_weight) * inner def _slanted_integral_image(double[:, ::1] image, @@ -58,7 +58,7 @@ def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, double[:, ::1] integral_img3, double[:, ::1] integral_img4, double[:, ::1] filtered_image, - double outer_wt, double inner_wt, + double outer_weight, double inner_weight, int mo, int no, int mi, int ni): cdef Py_ssize_t i, j, o_m, i_m, o_set, i_set @@ -82,4 +82,4 @@ def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m + 1, j - i_set - 1] - integral_img[i - i_m, -1] - integral_img[i - i_set - 1, j + i_m - 1] + integral_img[i - i_m, j + i_m - 1] + integral_img[i - i_set - 1, -1] inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set - 1, j + i_m - 1] - integral_img[-1, j + i_set] - integral_img[i + i_m - 1, j + i_m - 1] + integral_img[-1, j + i_m - 1] + integral_img[i + i_m - 1, j + i_set] - filtered_image[i, j] = outer_wt * outer - (outer_wt + inner_wt) * inner + filtered_image[i, j] = outer_weight * outer - (outer_weight + inner_weight) * inner From 4d9baceb8a389bd9d47245455d61482cc3eecb44 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 17:05:25 +0530 Subject: [PATCH 31/77] Documenting the code; Removing ascontiguousarray statements --- skimage/feature/censure.py | 26 +++++++++++++++------ skimage/feature/censure_cy.pyx | 41 +++++++++++++++++++++++++--------- 2 files changed, 50 insertions(+), 17 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 1175235f..c6ab7f38 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -24,6 +24,7 @@ def _get_filtered_image(image, n_scales, mode): _censure_dob_loop(image, n, integral_img, filtered_image, inner_weight, outer_weight) scales[:, :, i] = filtered_image return scales + elif mode == 'Octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] @@ -32,11 +33,9 @@ def _get_filtered_image(image, n_scales, mode): integral_img = integral_image(image) integral_img1 = _slanted_integral_image_modes(image, 1) integral_img2 = _slanted_integral_image_modes(image, 2) - integral_img2 = np.ascontiguousarray(integral_img2) integral_img3 = _slanted_integral_image_modes(image, 3) - integral_img3 = np.ascontiguousarray(integral_img3) integral_img4 = _slanted_integral_image_modes(image, 4) - integral_img4 = np.ascontiguousarray(integral_img4) + for k in range(n_scales): n = k + 1 filtered_image = np.zeros(image.shape) @@ -59,7 +58,11 @@ def _slanted_integral_image_modes(img, mode=1): if mode == 1: """ The following figures describe area that is summed up to calculate - the value at point @ in slanted integral image. + the value at point @ in slanted integral image. The subtended at @ is + 135 degrees. + + censure_cy._slanted_integral_image performs the mode1 + _slanted_integral_image _________________ |********/ | |*******/ | @@ -70,12 +73,17 @@ def _slanted_integral_image_modes(img, mode=1): |_________________| """ image = np.copy(img, order='C') + mode1 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode1) return mode1[1:, :] elif mode == 2: """ + For mode2, the image can be first flipped left-right and then up-down. + Then we can use censure_cy._slanted_integral_image and the returned + result can be flipped left-right and then up-down to get the following + mode. _________________ | | | | @@ -88,9 +96,10 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img, order='C') image = np.fliplr(image) image = np.flipud(image) - image = np.ascontiguousarray(image) + mode2 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode2) + mode2 = mode2[1:, :] mode2 = np.fliplr(mode2) mode2 = np.flipud(mode2) @@ -110,9 +119,10 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img, order='C') image = np.flipud(image) image = image.T - image = np.ascontiguousarray(image) + mode3 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode3) + mode3 = mode3[1:, :] mode3 = np.flipud(mode3.T) return mode3 @@ -131,9 +141,10 @@ def _slanted_integral_image_modes(img, mode=1): image = np.copy(img, order='C') image = np.fliplr(image) image = image.T - image = np.ascontiguousarray(image) + mode4 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') _slanted_integral_image(image, mode4) + mode4 = mode4[1:, :] mode4 = np.fliplr(mode4.T) return mode4 @@ -143,6 +154,7 @@ def _suppress_line(response, sigma, rpc_threshold): Axx, Axy, Ayy = _compute_auto_correlation(response, sigma) detA = Axx * Ayy - Axy**2 traceA = Axx + Ayy + # ratio of principal curvatures rpc = traceA**2 / (detA + 0.001) response[rpc > rpc_threshold] = 0 diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 75f58278..b65b0d96 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -22,8 +22,8 @@ def _censure_dob_loop(double[:, ::1] image, Py_ssize_t n, filtered_image[i, j] = outer_weight * outer - (inner_weight + outer_weight) * inner -def _slanted_integral_image(double[:, ::1] image, - double[:, ::1] integral_img): +def _slanted_integral_image(double[:, :] image, + double[:, :] integral_img): cdef Py_ssize_t i, j cdef double[:] left_sum = np.zeros(image.shape[0], dtype=np.float) @@ -33,8 +33,10 @@ def _slanted_integral_image(double[:, ::1] image, left_sum[image.shape[1] - 1 - i] = np.sum(flipped_lr.diagonal(i)) left_sum_np = np.asarray(left_sum) + # Initializing the leftmost column of the slanted integral image left_sum_np = left_sum_np.cumsum(0) + # Initializing the rightmost column of the slanted integral image right_sum_np = np.sum(image, 1).cumsum(0) for i in range(image.shape[0]): @@ -50,32 +52,51 @@ def _slanted_integral_image(double[:, ::1] image, integral_img[i, j] += integral_img[i, j - 1] + integral_img[i - 1, j + 1] - integral_img[i - 1, j] - - -def _censure_octagon_loop(double[:, ::1] image, double[:, ::1] integral_img, - double[:, ::1] integral_img1, - double[:, ::1] integral_img2, - double[:, ::1] integral_img3, - double[:, ::1] integral_img4, - double[:, ::1] filtered_image, +def _censure_octagon_loop(double[:, :] image, double[:, :] integral_img, + double[:, :] integral_img1, + double[:, :] integral_img2, + double[:, :] integral_img3, + double[:, :] integral_img4, + double[:, :] filtered_image, double outer_weight, double inner_weight, int mo, int no, int mi, int ni): cdef Py_ssize_t i, j, o_m, i_m, o_set, i_set + """ + For a (5, 2) octagon, i.e. mo = 5 and no = 2, + + |---o_set---| + [0, 0, 1, 1, 1, 1, 1, 0, 0] + [0, 1, 1, 1, 1, 1, 1, 1, 0] + [1, 1, 1, 1, 1, 1, 1, 1, 1] + [1, 1, 1, 1, 1, 1, 1, 1, 1] + [1, 1, 1, 1, 1, 1, 1, 1, 1] + [1, 1, 1, 1, 1, 1, 1, 1, 1] + [1, 1, 1, 1, 1, 1, 1, 1, 1] + [0, 1, 1, 1, 1, 1, 1, 1, 0] + [0, 0, 1, 1, 1, 1, 1, 0, 0] + |-o_m-| + """ o_m = (mo - 1) / 2 i_m = (mi - 1) / 2 + + # o_set and i_set are the distances of the center of the octagon + # from the horizontal or vertical sides of the octagon, + # for outer and inner octagon respectively o_set = o_m + no i_set = i_m + ni for i in range(o_set + 1, image.shape[0] - o_set - 1): for j in range(o_set + 1, image.shape[1] - o_set - 1): + # Calculating the sum of pixels in the outer octagon outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m - 1, j + o_set + 1] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m - 1, j - o_m] outer += integral_img[i + o_m - 1, j + o_m - 1] - integral_img[i - o_m, j + o_m - 1] - integral_img[i + o_m - 1, j - o_m] + integral_img[i - o_m, j - o_m] outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set + 1, j - o_m + 1] - integral_img[i - o_m, j - o_m] + integral_img[i - o_m, j - o_set - 1] outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m + 1, j - o_set - 1] - integral_img[i - o_m, -1] - integral_img[i - o_set - 1, j + o_m - 1] + integral_img[i - o_m, j + o_m - 1] + integral_img[i - o_set - 1, -1] outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set - 1, j + o_m - 1] - integral_img[-1, j + o_set] - integral_img[i + o_m - 1, j + o_m - 1] + integral_img[-1, j + o_m - 1] + integral_img[i + o_m - 1, j + o_set] + # Calculating the sum of pixels in the inner octagon inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m - 1, j + i_set + 1] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m - 1, j - i_m] inner += integral_img[i + i_m - 1, j + i_m - 1] - integral_img[i - i_m, j + i_m - 1] - integral_img[i + i_m - 1, j - i_m] + integral_img[i - i_m, j - i_m] inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set + 1, j - i_m + 1] - integral_img[i - i_m, j - i_m] + integral_img[i - i_m, j - i_set - 1] From 7478f2c79656d61eb1b6d886334b1c861bd57b5c Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 21:40:25 +0530 Subject: [PATCH 32/77] Correcting indentation --- skimage/feature/censure_cy.pyx | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index b65b0d96..92d5a142 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -8,9 +8,9 @@ import numpy as np def _censure_dob_loop(double[:, ::1] image, Py_ssize_t n, - double[:, ::1] integral_img, - double[:, ::1] filtered_image, - double inner_weight, double outer_weight): + double[:, ::1] integral_img, + double[:, ::1] filtered_image, + double inner_weight, double outer_weight): cdef Py_ssize_t i, j cdef double inner, outer From 3318886ff30e2c2a564e5eb277061ad6da9c157c Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 22:40:12 +0530 Subject: [PATCH 33/77] Using convolve in _get_filtered_image for mode=Octagon --- skimage/feature/censure.py | 56 ++++++++++++++++++++++++++++++++++---- 1 file changed, 51 insertions(+), 5 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index c6ab7f38..7f2ae816 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -1,17 +1,18 @@ import numpy as np -from scipy.ndimage.filters import maximum_filter, minimum_filter +from scipy.ndimage.filters import maximum_filter, minimum_filter, convolve from ..transform import integral_image from ..feature.corner import _compute_auto_correlation from ..util import img_as_float +from ..morphology import convex_hull_image from .censure_cy import _censure_dob_loop, _slanted_integral_image, _censure_octagon_loop def _get_filtered_image(image, n_scales, mode): # TODO : Implement the STAR mode + scales = np.zeros((image.shape[0], image.shape[1], n_scales), dtype=np.double) if mode == 'DoB': - scales = np.zeros((image.shape[0], image.shape[1], n_scales)) for i in range(n_scales): n = i + 1 # Constant multipliers for the outer region and the inner region @@ -23,13 +24,15 @@ def _get_filtered_image(image, n_scales, mode): filtered_image = np.zeros(image.shape) _censure_dob_loop(image, n, integral_img, filtered_image, inner_weight, outer_weight) scales[:, :, i] = filtered_image - return scales + elif mode == 'Octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] - scales = np.zeros((image.shape[0], image.shape[1], n_scales)) + for i in range(n_scales): + scales[:, :, i] = convolve(image, _octagon_filter(outer_shape[i][0], outer_shape[i][1], inner_shape[i][0], inner_shape[i][1])) + """ integral_img = integral_image(image) integral_img1 = _slanted_integral_image_modes(image, 1) integral_img2 = _slanted_integral_image_modes(image, 2) @@ -51,7 +54,50 @@ def _get_filtered_image(image, n_scales, mode): _censure_octagon_loop(image, integral_img, integral_img1, integral_img2, integral_img3, integral_img4, filtered_image, outer_weight, inner_weight, mo, no, mi, ni) scales[:, :, k] = filtered_image - return scales + """ + return scales + + +def _oct(m, n): + f = np.zeros((m + 2*n, m + 2*n)) + f[0, n] = 1 + f[n, 0] = 1 + f[0, m + n -1] = 1 + f[m + n - 1, 0] = 1 + f[-1, n] = 1 + f[n, -1] = 1 + f[-1, m + n - 1] = 1 + f[m + n - 1, -1] = 1 + return convex_hull_image(f).astype(int) + + +def _octagon_filter(mo, no, mi, ni): + outer = (mo + 2 * no)**2 - 2 * no * (no + 1) + inner = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) + outer_wt = 1.0 / (outer - inner) + inner_wt = 1.0 / inner + c = ((mo + 2 * no) - (mi + 2 * ni)) / 2 + outer_oct = _oct(mo, no) + inner_oct = np.zeros((mo + 2 * no, mo + 2 * no)) + inner_oct[c:-c, c:-c] = _oct(mi, ni) + bfilter = outer_wt * outer_oct - (outer_wt + inner_wt) * inner_oct + return bfilter + + +def _filter_using_convolve(image, n, mode='DoB'): + + if mode == 'DoB': + inner_wt = (1.0 / (2*n + 1)**2) + outer_wt = (1.0 / (12*n**2 + 4*n)) + dob_filter = np.zeros((4 * n + 1, 4 * n + 1)) + dob_filter[:] = outer_wt + dob_filter[n : 3 * n + 1, n : 3 * n + 1] = - inner_wt + return convolve(image, dob_filter) + + elif mode == 'Octagon': + outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] + inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] + return convolve(image, _octagon_filter(outer_shape[n - 1][0], outer_shape[n - 1][1], inner_shape[n - 1][0], inner_shape[n - 1][1])) def _slanted_integral_image_modes(img, mode=1): From 344374e290ee820d698e67dd68a66e0179b6f973 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 23:15:45 +0530 Subject: [PATCH 34/77] Removing slanted_integral_image Cython functions for mode=Octagon --- skimage/feature/censure_cy.pyx | 84 ---------------------------------- 1 file changed, 84 deletions(-) diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 92d5a142..2d71f0dd 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -20,87 +20,3 @@ def _censure_dob_loop(double[:, ::1] image, Py_ssize_t n, inner = integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n] outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n - 1, j - 2 * n - 1] - integral_img[i + 2 * n, j - 2 * n - 1] - integral_img[i - 2 * n - 1, j + 2 * n] filtered_image[i, j] = outer_weight * outer - (inner_weight + outer_weight) * inner - - -def _slanted_integral_image(double[:, :] image, - double[:, :] integral_img): - - cdef Py_ssize_t i, j - cdef double[:] left_sum = np.zeros(image.shape[0], dtype=np.float) - - flipped_lr = np.asarray(image[:, ::-1]) - for i in range(image.shape[1] - image.shape[0], image.shape[1]): - left_sum[image.shape[1] - 1 - i] = np.sum(flipped_lr.diagonal(i)) - left_sum_np = np.asarray(left_sum) - - # Initializing the leftmost column of the slanted integral image - left_sum_np = left_sum_np.cumsum(0) - - # Initializing the rightmost column of the slanted integral image - right_sum_np = np.sum(image, 1).cumsum(0) - - for i in range(image.shape[0]): - image[i, 0] = left_sum_np[i] - image[i, -1] = right_sum_np[i] - - for i in range(1, integral_img.shape[0]): - for j in range(integral_img.shape[1]): - integral_img[i, j] = image[i - 1, j] - - for i in range(1, integral_img.shape[0]): - for j in range(1, integral_img.shape[1] - 1): - integral_img[i, j] += integral_img[i, j - 1] + integral_img[i - 1, j + 1] - integral_img[i - 1, j] - - -def _censure_octagon_loop(double[:, :] image, double[:, :] integral_img, - double[:, :] integral_img1, - double[:, :] integral_img2, - double[:, :] integral_img3, - double[:, :] integral_img4, - double[:, :] filtered_image, - double outer_weight, double inner_weight, - int mo, int no, int mi, int ni): - - cdef Py_ssize_t i, j, o_m, i_m, o_set, i_set - - """ - For a (5, 2) octagon, i.e. mo = 5 and no = 2, - - |---o_set---| - [0, 0, 1, 1, 1, 1, 1, 0, 0] - [0, 1, 1, 1, 1, 1, 1, 1, 0] - [1, 1, 1, 1, 1, 1, 1, 1, 1] - [1, 1, 1, 1, 1, 1, 1, 1, 1] - [1, 1, 1, 1, 1, 1, 1, 1, 1] - [1, 1, 1, 1, 1, 1, 1, 1, 1] - [1, 1, 1, 1, 1, 1, 1, 1, 1] - [0, 1, 1, 1, 1, 1, 1, 1, 0] - [0, 0, 1, 1, 1, 1, 1, 0, 0] - |-o_m-| - """ - o_m = (mo - 1) / 2 - i_m = (mi - 1) / 2 - - # o_set and i_set are the distances of the center of the octagon - # from the horizontal or vertical sides of the octagon, - # for outer and inner octagon respectively - o_set = o_m + no - i_set = i_m + ni - - for i in range(o_set + 1, image.shape[0] - o_set - 1): - for j in range(o_set + 1, image.shape[1] - o_set - 1): - # Calculating the sum of pixels in the outer octagon - outer = integral_img1[i + o_set, j + o_m] - integral_img1[i + o_m - 1, j + o_set + 1] - integral_img[i + o_set, j - o_m] + integral_img[i + o_m - 1, j - o_m] - outer += integral_img[i + o_m - 1, j + o_m - 1] - integral_img[i - o_m, j + o_m - 1] - integral_img[i + o_m - 1, j - o_m] + integral_img[i - o_m, j - o_m] - outer += integral_img4[i + o_m, j - o_set] - integral_img4[i + o_set + 1, j - o_m + 1] - integral_img[i - o_m, j - o_m] + integral_img[i - o_m, j - o_set - 1] - outer += integral_img2[i - o_set, j - o_m] - integral_img2[i - o_m + 1, j - o_set - 1] - integral_img[i - o_m, -1] - integral_img[i - o_set - 1, j + o_m - 1] + integral_img[i - o_m, j + o_m - 1] + integral_img[i - o_set - 1, -1] - outer += integral_img3[i - o_m, j + o_set] - integral_img3[i - o_set - 1, j + o_m - 1] - integral_img[-1, j + o_set] - integral_img[i + o_m - 1, j + o_m - 1] + integral_img[-1, j + o_m - 1] + integral_img[i + o_m - 1, j + o_set] - - # Calculating the sum of pixels in the inner octagon - inner = integral_img1[i + i_set, j + i_m] - integral_img1[i + i_m - 1, j + i_set + 1] - integral_img[i + i_set, j - i_m] + integral_img[i + i_m - 1, j - i_m] - inner += integral_img[i + i_m - 1, j + i_m - 1] - integral_img[i - i_m, j + i_m - 1] - integral_img[i + i_m - 1, j - i_m] + integral_img[i - i_m, j - i_m] - inner += integral_img4[i + i_m, j - i_set] - integral_img4[i + i_set + 1, j - i_m + 1] - integral_img[i - i_m, j - i_m] + integral_img[i - i_m, j - i_set - 1] - inner += integral_img2[i - i_set, j - i_m] - integral_img2[i - i_m + 1, j - i_set - 1] - integral_img[i - i_m, -1] - integral_img[i - i_set - 1, j + i_m - 1] + integral_img[i - i_m, j + i_m - 1] + integral_img[i - i_set - 1, -1] - inner += integral_img3[i - i_m, j + i_set] - integral_img3[i - i_set - 1, j + i_m - 1] - integral_img[-1, j + i_set] - integral_img[i + i_m - 1, j + i_m - 1] + integral_img[-1, j + i_m - 1] + integral_img[i + i_m - 1, j + i_set] - - filtered_image[i, j] = outer_weight * outer - (outer_weight + inner_weight) * inner From 5ee71cba9b0bd97fa7ee15c76936f1f0a2d3c49c Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 23:16:47 +0530 Subject: [PATCH 35/77] Reverting back to non-debugging mode in Cython --- skimage/feature/censure_cy.pyx | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 2d71f0dd..1c007609 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -1,7 +1,7 @@ #cython: cdivision=True -#cython: boundscheck=True -#cython: nonecheck=True -#cython: wraparound=True +#cython: boundscheck=False +#cython: nonecheck=False +#cython: wraparound=False cimport numpy as cnp import numpy as np From 599ce8f7bd606d7b884c81b926c63b865b722391 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 23:20:59 +0530 Subject: [PATCH 36/77] Remove all the slanted integral image functions and the vivid graphics from censure.py --- skimage/feature/censure.py | 138 +------------------------------------ 1 file changed, 2 insertions(+), 136 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 7f2ae816..9acbf7cd 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -6,7 +6,7 @@ from ..feature.corner import _compute_auto_correlation from ..util import img_as_float from ..morphology import convex_hull_image -from .censure_cy import _censure_dob_loop, _slanted_integral_image, _censure_octagon_loop +from .censure_cy import _censure_dob_loop def _get_filtered_image(image, n_scales, mode): @@ -25,39 +25,17 @@ def _get_filtered_image(image, n_scales, mode): _censure_dob_loop(image, n, integral_img, filtered_image, inner_weight, outer_weight) scales[:, :, i] = filtered_image - elif mode == 'Octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] for i in range(n_scales): scales[:, :, i] = convolve(image, _octagon_filter(outer_shape[i][0], outer_shape[i][1], inner_shape[i][0], inner_shape[i][1])) - """ - integral_img = integral_image(image) - integral_img1 = _slanted_integral_image_modes(image, 1) - integral_img2 = _slanted_integral_image_modes(image, 2) - integral_img3 = _slanted_integral_image_modes(image, 3) - integral_img4 = _slanted_integral_image_modes(image, 4) - for k in range(n_scales): - n = k + 1 - filtered_image = np.zeros(image.shape) - mo = outer_shape[n - 1][0] - no = outer_shape[n - 1][1] - mi = inner_shape[n - 1][0] - ni = inner_shape[n - 1][1] - outer_pixels = (mo + 2 * no)**2 - 2 * no * (no + 1) - inner_pixels = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) - outer_weight = 1.0 / (outer_pixels - inner_pixels) - inner_weight = 1.0 / inner_pixels - - _censure_octagon_loop(image, integral_img, integral_img1, integral_img2, integral_img3, integral_img4, filtered_image, outer_weight, inner_weight, mo, no, mi, ni) - - scales[:, :, k] = filtered_image - """ return scales +# TODO : Import from selem after getting #669 merged. def _oct(m, n): f = np.zeros((m + 2*n, m + 2*n)) f[0, n] = 1 @@ -84,118 +62,6 @@ def _octagon_filter(mo, no, mi, ni): return bfilter -def _filter_using_convolve(image, n, mode='DoB'): - - if mode == 'DoB': - inner_wt = (1.0 / (2*n + 1)**2) - outer_wt = (1.0 / (12*n**2 + 4*n)) - dob_filter = np.zeros((4 * n + 1, 4 * n + 1)) - dob_filter[:] = outer_wt - dob_filter[n : 3 * n + 1, n : 3 * n + 1] = - inner_wt - return convolve(image, dob_filter) - - elif mode == 'Octagon': - outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] - inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] - return convolve(image, _octagon_filter(outer_shape[n - 1][0], outer_shape[n - 1][1], inner_shape[n - 1][0], inner_shape[n - 1][1])) - - -def _slanted_integral_image_modes(img, mode=1): - if mode == 1: - """ - The following figures describe area that is summed up to calculate - the value at point @ in slanted integral image. The subtended at @ is - 135 degrees. - - censure_cy._slanted_integral_image performs the mode1 - _slanted_integral_image - _________________ - |********/ | - |*******/ | - |******/ | - |-----@ | - | | - | | - |_________________| - """ - image = np.copy(img, order='C') - - mode1 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') - _slanted_integral_image(image, mode1) - return mode1[1:, :] - - elif mode == 2: - """ - For mode2, the image can be first flipped left-right and then up-down. - Then we can use censure_cy._slanted_integral_image and the returned - result can be flipped left-right and then up-down to get the following - mode. - _________________ - | | - | | - | | - | @_____| - | /******| - | /*******| - |________/________| - """ - image = np.copy(img, order='C') - image = np.fliplr(image) - image = np.flipud(image) - - mode2 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') - _slanted_integral_image(image, mode2) - - mode2 = mode2[1:, :] - mode2 = np.fliplr(mode2) - mode2 = np.flipud(mode2) - return mode2 - - elif mode == 3: - """ - _________________ - | | - |\\ | - |*\\ | - |**\\ | - |***@ | - |***| | - |___|_____________| - """ - image = np.copy(img, order='C') - image = np.flipud(image) - image = image.T - - mode3 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') - _slanted_integral_image(image, mode3) - - mode3 = mode3[1:, :] - mode3 = np.flipud(mode3.T) - return mode3 - - else: - """ - ________________ - | |****| - | |****| - | @****| - | \\**| - | \\*| - | \\| - |________________| - """ - image = np.copy(img, order='C') - image = np.fliplr(image) - image = image.T - - mode4 = np.zeros((image.shape[0] + 1, image.shape[1]), order='C') - _slanted_integral_image(image, mode4) - - mode4 = mode4[1:, :] - mode4 = np.fliplr(mode4.T) - return mode4 - - def _suppress_line(response, sigma, rpc_threshold): Axx, Axy, Ayy = _compute_auto_correlation(response, sigma) detA = Axx * Ayy - Axy**2 From 6472dd10766d36a6e54ba492c9b6441c45504400 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 31 Jul 2013 23:40:23 +0530 Subject: [PATCH 37/77] Wrapping the long lines --- skimage/feature/censure.py | 33 ++++++++++++++++++++++++--------- skimage/feature/censure_cy.pyx | 15 ++++++++++++--- 2 files changed, 36 insertions(+), 12 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 9acbf7cd..8b8e1ae4 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -11,26 +11,38 @@ from .censure_cy import _censure_dob_loop def _get_filtered_image(image, n_scales, mode): # TODO : Implement the STAR mode - scales = np.zeros((image.shape[0], image.shape[1], n_scales), dtype=np.double) + scales = np.zeros((image.shape[0], image.shape[1], n_scales), + dtype=np.double) + if mode == 'DoB': for i in range(n_scales): n = i + 1 + # Constant multipliers for the outer region and the inner region - # of the bilevel filters with the constraint of keeping the DC bias - # 0. + # of the bilevel filters with the constraint of keeping the + # DC bias 0. inner_weight = (1.0 / (2 * n + 1)**2) outer_weight = (1.0 / (12 * n**2 + 4 * n)) + integral_img = integral_image(image) + filtered_image = np.zeros(image.shape) - _censure_dob_loop(image, n, integral_img, filtered_image, inner_weight, outer_weight) + _censure_dob_loop(image, n, integral_img, filtered_image, + inner_weight, outer_weight) + scales[:, :, i] = filtered_image elif mode == 'Octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 - outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] + outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), + (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] + for i in range(n_scales): - scales[:, :, i] = convolve(image, _octagon_filter(outer_shape[i][0], outer_shape[i][1], inner_shape[i][0], inner_shape[i][1])) + scales[:, :, i] = convolve(image, + _octagon_filter(outer_shape[i][0], + outer_shape[i][1], inner_shape[i][0], + inner_shape[i][1])) return scales @@ -73,7 +85,8 @@ def _suppress_line(response, sigma, rpc_threshold): return response -def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, rpc_threshold=10): +def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, + rpc_threshold=10): """ Extracts Censure keypoints along with the corresponding scale using either Difference of Boxes, Octagon or STAR bilevel filter. @@ -142,8 +155,10 @@ def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, rpc_thresho # sigma = (window_size - 1) / 6.0 # window_size = 7 + 2 * i # Hence sigma = 1 + i / 3.0 - response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), rpc_threshold) + response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), + rpc_threshold) # Returning keypoints with its scale - keypoints = np.transpose(np.nonzero(response[:, :, 1:n_scales - 1])) + [0, 0, 2] + keypoints = (np.transpose(np.nonzero(response[:, :, 1:n_scales - 1])) + + [0, 0, 2]) return keypoints diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 1c007609..93c7c142 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -17,6 +17,15 @@ def _censure_dob_loop(double[:, ::1] image, Py_ssize_t n, for i in range(2 * n, image.shape[0] - 2 * n): for j in range(2 * n, image.shape[1] - 2 * n): - inner = integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n] - outer = integral_img[i + 2 * n, j + 2 * n] + integral_img[i - 2 * n - 1, j - 2 * n - 1] - integral_img[i + 2 * n, j - 2 * n - 1] - integral_img[i - 2 * n - 1, j + 2 * n] - filtered_image[i, j] = outer_weight * outer - (inner_weight + outer_weight) * inner + inner = (integral_img[i + n, j + n] + + integral_img[i - n - 1, j - n - 1] + - integral_img[i + n, j - n - 1] + - integral_img[i - n - 1, j + n]) + + outer = (integral_img[i + 2 * n, j + 2 * n] + + integral_img[i - 2 * n - 1, j - 2 * n - 1] + - integral_img[i + 2 * n, j - 2 * n - 1] + - integral_img[i - 2 * n - 1, j + 2 * n]) + + filtered_image[i, j] = (outer_weight * outer + - (inner_weight + outer_weight) * inner) From 361652cec209e8efb975af26b77bf87cbc893dda Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 1 Aug 2013 01:11:26 +0530 Subject: [PATCH 38/77] Added a NOTE explaining the preference of convolve over slanted integral image --- skimage/feature/censure.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 8b8e1ae4..d0ac94a9 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -32,12 +32,18 @@ def _get_filtered_image(image, n_scales, mode): scales[:, :, i] = filtered_image + # NOTE : For the Octagon shaped filter, we implemented and evaluated the + # slanted integral image based image filtering but the performance was + # more or less equal to image filtering using + # scipy.ndimage.filters.convolve(). Hence we have decided to use the + # later for a much cleaner implementation. elif mode == 'Octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] + # for i in range(n_scales): scales[:, :, i] = convolve(image, _octagon_filter(outer_shape[i][0], From 212874fae340f34b1ff14ccbc7bcfffd9476e0e4 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Fri, 2 Aug 2013 13:35:29 +0530 Subject: [PATCH 39/77] Adding the STAR mode --- skimage/feature/censure.py | 52 ++++++++++++++++++++++++++++++++++---- 1 file changed, 47 insertions(+), 5 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index d0ac94a9..179eac70 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -49,6 +49,15 @@ def _get_filtered_image(image, n_scales, mode): _octagon_filter(outer_shape[i][0], outer_shape[i][1], inner_shape[i][0], inner_shape[i][1])) + else: + shape = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, + 128] + filter_shape = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), + (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] + for i in range(n_scales): + scales[:, :, i] = convolve(image, + _star_filter(shape[filter_shape[i][0]], + shape[filter_shape[i][1]])) return scales @@ -70,13 +79,46 @@ def _oct(m, n): def _octagon_filter(mo, no, mi, ni): outer = (mo + 2 * no)**2 - 2 * no * (no + 1) inner = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) - outer_wt = 1.0 / (outer - inner) - inner_wt = 1.0 / inner + outer_weight = 1.0 / (outer - inner) + inner_weight = 1.0 / inner c = ((mo + 2 * no) - (mi + 2 * ni)) / 2 outer_oct = _oct(mo, no) inner_oct = np.zeros((mo + 2 * no, mo + 2 * no)) - inner_oct[c:-c, c:-c] = _oct(mi, ni) - bfilter = outer_wt * outer_oct - (outer_wt + inner_wt) * inner_oct + inner_oct[c: -c, c: -c] = _oct(mi, ni) + bfilter = (outer_weight * outer_oct - + (outer_weight + inner_weight) * inner_oct) + return bfilter + + +def _star(a): + if a == 1: + bfilter = np.zeros((3, 3)) + bfilter[:] = 1 + return bfilter + m = 2 * a + 1 + n = a / 2 + selem_square = np.zeros((m + 2 * n, m + 2 * n), dtype=np.uint8) + selem_square[n: m + n, n: m + n] = 1 + selem_triangle = np.zeros((m + 2 * n, m + 2 * n), dtype=np.uint8) + selem_triangle[(m + 2 * n - 1) / 2, 0] = 1 + selem_triangle[(m + 1) / 2, n - 1] = 1 + selem_triangle[(m + 4 * n - 3) / 2, n - 1] = 1 + selem_triangle = convex_hull_image(selem_triangle).astype(int) + selem_triangle += selem_triangle[:, ::-1] + selem_triangle.T + selem_triangle.T[::-1, :] + return selem_square + selem_triangle + + +def _star_filter(m, n): + outer = 4 * m**2 + 4 * m + 1 + 4 * (m / 2)**2 + inner = 4 * n**2 + 4 * n + 1 + 4 * (n / 2)**2 + outer_weight = 1.0 / (outer - inner) + inner_weight = 1.0 / inner + c = m + m / 2 - n - n / 2 + outer_star = _star(m) + inner_star = np.zeros((outer_star.shape)) + inner_star[c: -c, c: -c] = _star(n) + bfilter = (outer_weight * outer_star - + (outer_weight + inner_weight) * inner_star) return bfilter @@ -106,7 +148,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, Number of scales to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last. - mode : {'DoB', 'Octagon', 'STAR'} + mode : ('DoB', 'Octagon', 'STAR') Type of bilevel filter used to get the scales of input image. Possible values are 'DoB', 'Octagon' and 'STAR'. From 8ef8ba59d0c8792fe5618da0fe4ab415e355c855 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Fri, 2 Aug 2013 19:02:16 +0530 Subject: [PATCH 40/77] Making the weight calculation statements more readable --- skimage/feature/censure.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 179eac70..1f3506bb 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -109,14 +109,12 @@ def _star(a): def _star_filter(m, n): - outer = 4 * m**2 + 4 * m + 1 + 4 * (m / 2)**2 - inner = 4 * n**2 + 4 * n + 1 + 4 * (n / 2)**2 - outer_weight = 1.0 / (outer - inner) - inner_weight = 1.0 / inner c = m + m / 2 - n - n / 2 outer_star = _star(m) inner_star = np.zeros((outer_star.shape)) inner_star[c: -c, c: -c] = _star(n) + outer_weight = 1.0 / (np.sum(outer_star - inner_star)) + inner_weight = 1.0 / np.sum(inner_star) bfilter = (outer_weight * outer_star - (outer_weight + inner_weight) * inner_star) return bfilter From 973ab0a5480ad4c2fa6d997aa88f01d6dbb5e36f Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Sat, 3 Aug 2013 14:14:35 +0530 Subject: [PATCH 41/77] PEP8 corrections --- skimage/feature/censure.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 1f3506bb..6ca6fef4 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -51,7 +51,7 @@ def _get_filtered_image(image, n_scales, mode): inner_shape[i][1])) else: shape = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, - 128] + 128] filter_shape = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] for i in range(n_scales): @@ -104,7 +104,8 @@ def _star(a): selem_triangle[(m + 1) / 2, n - 1] = 1 selem_triangle[(m + 4 * n - 3) / 2, n - 1] = 1 selem_triangle = convex_hull_image(selem_triangle).astype(int) - selem_triangle += selem_triangle[:, ::-1] + selem_triangle.T + selem_triangle.T[::-1, :] + selem_triangle += (selem_triangle[:, ::-1] + selem_triangle.T + + selem_triangle.T[::-1, :]) return selem_square + selem_triangle From 87ef133e809e6399ed65638ae7cb2ada4446121d Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 5 Aug 2013 13:26:59 +0530 Subject: [PATCH 42/77] Replacing threshold by nms_threshold --- skimage/feature/censure.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 6ca6fef4..0ad69fca 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -132,7 +132,7 @@ def _suppress_line(response, sigma, rpc_threshold): return response -def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, +def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.03, rpc_threshold=10): """ Extracts Censure keypoints along with the corresponding scale using @@ -151,7 +151,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, Type of bilevel filter used to get the scales of input image. Possible values are 'DoB', 'Octagon' and 'STAR'. - threshold : float + nms_threshold : float Threshold value used to suppress maximas and minimas with a weak magnitude response obtained after Non-Maximal Suppression. @@ -193,9 +193,9 @@ def censure_keypoints(image, n_scales=7, mode='DoB', threshold=0.03, minimas = (minimum_filter(scales, (3, 3, 3)) == scales) * scales maximas = (maximum_filter(scales, (3, 3, 3)) == scales) * scales - # Suppressing minimas and maximas weaker than threshold - minimas[np.abs(minimas) < threshold] = 0 - maximas[np.abs(maximas) < threshold] = 0 + # Suppressing minimas and maximas weaker than nms_threshold + minimas[np.abs(minimas) < nms_threshold] = 0 + maximas[np.abs(maximas) < nms_threshold] = 0 response = maximas + minimas for i in range(1, n_scales - 1): From 54f9b06e464741631e6dbd59fe6f2791834a9ebb Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 5 Aug 2013 23:02:23 +0530 Subject: [PATCH 43/77] Returning the keypoints and scales separately --- skimage/feature/censure.py | 19 +++++++++++++------ 1 file changed, 13 insertions(+), 6 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 0ad69fca..fa7bd133 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -10,7 +10,7 @@ from .censure_cy import _censure_dob_loop def _get_filtered_image(image, n_scales, mode): - # TODO : Implement the STAR mode + scales = np.zeros((image.shape[0], image.shape[1], n_scales), dtype=np.double) @@ -132,7 +132,7 @@ def _suppress_line(response, sigma, rpc_threshold): return response -def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.03, +def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.15, rpc_threshold=10): """ Extracts Censure keypoints along with the corresponding scale using @@ -161,8 +161,11 @@ def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.03, Returns ------- - keypoints : (N, 3) array - Location of extracted keypoints along with the corresponding scale. + keypoints : (N, 2) array + Location of extracted keypoints. + + scale : (N, 1) array + The corresponding scale of the N extracted keypoints. References ---------- @@ -206,6 +209,10 @@ def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.03, rpc_threshold) # Returning keypoints with its scale - keypoints = (np.transpose(np.nonzero(response[:, :, 1:n_scales - 1])) + keypoints_with_scale = (np.transpose(np.nonzero(response[:, :, 1:n_scales - 1])) + [0, 0, 2]) - return keypoints + + keypoints = keypoints_with_scale[:, :2] + scale = keypoints_with_scale[:, -1] + + return keypoints, scale From ccbca1349b18eea9c1eb00b9197427223972290b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Johannes=20Sch=C3=B6nberger?= Date: Mon, 5 Aug 2013 23:34:16 +0200 Subject: [PATCH 44/77] Fix bugs in censure keypoint detector and improve code --- skimage/feature/censure.py | 67 ++++++++++++++++---------------------- 1 file changed, 28 insertions(+), 39 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index fa7bd133..966921f2 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -1,12 +1,12 @@ import numpy as np from scipy.ndimage.filters import maximum_filter, minimum_filter, convolve -from ..transform import integral_image -from ..feature.corner import _compute_auto_correlation -from ..util import img_as_float -from ..morphology import convex_hull_image +from skimage.transform import integral_image +from skimage.feature.corner import _compute_auto_correlation +from skimage.util import img_as_float +from skimage.morphology import convex_hull_image -from .censure_cy import _censure_dob_loop +from skimage.feature.censure_cy import _censure_dob_loop def _get_filtered_image(image, n_scales, mode): @@ -43,15 +43,13 @@ def _get_filtered_image(image, n_scales, mode): (15, 10)] inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] - # for i in range(n_scales): scales[:, :, i] = convolve(image, _octagon_filter(outer_shape[i][0], outer_shape[i][1], inner_shape[i][0], inner_shape[i][1])) else: - shape = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, - 128] + shape = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] filter_shape = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] for i in range(n_scales): @@ -69,7 +67,7 @@ def _oct(m, n): f[n, 0] = 1 f[0, m + n -1] = 1 f[m + n - 1, 0] = 1 - f[-1, n] = 1 + f[-1, n] = 1 f[n, -1] = 1 f[-1, m + n - 1] = 1 f[m + n - 1, -1] = 1 @@ -121,19 +119,15 @@ def _star_filter(m, n): return bfilter -def _suppress_line(response, sigma, rpc_threshold): - Axx, Axy, Ayy = _compute_auto_correlation(response, sigma) - detA = Axx * Ayy - Axy**2 - traceA = Axx + Ayy - - # ratio of principal curvatures - rpc = traceA**2 / (detA + 0.001) - response[rpc > rpc_threshold] = 0 - return response +def _suppress_lines(feature_mask, image, sigma, line_threshold): + Axx, Axy, Ayy = _compute_auto_correlation(image, sigma) + feature_mask[(Axx + Ayy) * (Axx + Ayy) + < line_threshold * (Axx * Ayy - Axy * Axy)] = 0 + return feature_mask -def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.15, - rpc_threshold=10): +def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, + line_threshold=10): """ Extracts Censure keypoints along with the corresponding scale using either Difference of Boxes, Octagon or STAR bilevel filter. @@ -151,11 +145,11 @@ def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.15, Type of bilevel filter used to get the scales of input image. Possible values are 'DoB', 'Octagon' and 'STAR'. - nms_threshold : float + non_max_threshold : float Threshold value used to suppress maximas and minimas with a weak magnitude response obtained after Non-Maximal Suppression. - rpc_threshold : float + line_threshold : float Threshold for rejecting interest points which have ratio of principal curvatures greater than this value. @@ -176,7 +170,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.15, .. [2] Adam Schmidt, Marek Kraft, Michal Fularz and Zuzanna Domagala "Comparative Assessment of Point Feature Detectors and - Descriptors in the Context of Robot Navigation" + Descriptors in the Context of Robot Navigation" http://www.jamris.org/01_2013/saveas.php?QUEST=JAMRIS_No01_2013_P_11-20.pdf """ @@ -189,30 +183,25 @@ def censure_keypoints(image, n_scales=7, mode='DoB', nms_threshold=0.15, image = np.ascontiguousarray(image) # Generating all the scales - scales = _get_filtered_image(image, n_scales, mode) + filter_response = _get_filtered_image(image, n_scales, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero - minimas = (minimum_filter(scales, (3, 3, 3)) == scales) * scales - maximas = (maximum_filter(scales, (3, 3, 3)) == scales) * scales + minimas = minimum_filter(filter_response, (3, 3, 3)) == filter_response + maximas = maximum_filter(filter_response, (3, 3, 3)) == filter_response - # Suppressing minimas and maximas weaker than nms_threshold - minimas[np.abs(minimas) < nms_threshold] = 0 - maximas[np.abs(maximas) < nms_threshold] = 0 - response = maximas + minimas + feature_mask = minimas | maximas + feature_mask[filter_response < non_max_threshold] = False for i in range(1, n_scales - 1): # sigma = (window_size - 1) / 6.0 # window_size = 7 + 2 * i # Hence sigma = 1 + i / 3.0 - response[:, :, i] = _suppress_line(response[:, :, i], (1 + i / 3.0), - rpc_threshold) + feature_mask[:, :, i] = _suppress_lines(feature_mask[:, :, i], image, + (1 + i / 3.0), line_threshold) - # Returning keypoints with its scale - keypoints_with_scale = (np.transpose(np.nonzero(response[:, :, 1:n_scales - 1])) - + [0, 0, 2]) + rows, cols, scales = np.nonzero(feature_mask[..., 1:n_scales - 1]) + keypoints = np.column_stack([rows, cols]) + scales = scales + 2 - keypoints = keypoints_with_scale[:, :2] - scale = keypoints_with_scale[:, -1] - - return keypoints, scale + return keypoints, scales From bea9aa4414d923431cc1b3c3cc14d35a4321e55f Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 6 Aug 2013 05:12:25 +0530 Subject: [PATCH 45/77] Changing the variable names of constant objects --- skimage/feature/censure.py | 28 +++++++++++++++------------- 1 file changed, 15 insertions(+), 13 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 966921f2..c04f7d24 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -14,6 +14,13 @@ def _get_filtered_image(image, n_scales, mode): scales = np.zeros((image.shape[0], image.shape[1], n_scales), dtype=np.double) + OCTAGON_OUTER_SHAPE = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), + (15, 10)] + OCTAGON_INNER_SHAPE = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] + + STAR_SHAPE = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] + STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), + (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] if mode == 'DoB': for i in range(n_scales): n = i + 1 @@ -39,23 +46,18 @@ def _get_filtered_image(image, n_scales, mode): # later for a much cleaner implementation. elif mode == 'Octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 - outer_shape = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), - (15, 10)] - inner_shape = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] for i in range(n_scales): scales[:, :, i] = convolve(image, - _octagon_filter(outer_shape[i][0], - outer_shape[i][1], inner_shape[i][0], - inner_shape[i][1])) + _octagon_filter_kernel(OCTAGON_OUTER_SHAPE[i][0], + OCTAGON_OUTER_SHAPE[i][1], OCTAGON_INNER_SHAPE[i][0], + OCTAGON_INNER_SHAPE[i][1])) else: - shape = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] - filter_shape = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), - (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] + for i in range(n_scales): scales[:, :, i] = convolve(image, - _star_filter(shape[filter_shape[i][0]], - shape[filter_shape[i][1]])) + _star_filter_kernel(STAR_SHAPE[STAR_FILTER_SHAPE[i][0]], + STAR_SHAPE[STAR_FILTER_SHAPE[i][1]])) return scales @@ -74,7 +76,7 @@ def _oct(m, n): return convex_hull_image(f).astype(int) -def _octagon_filter(mo, no, mi, ni): +def _octagon_filter_kernel(mo, no, mi, ni): outer = (mo + 2 * no)**2 - 2 * no * (no + 1) inner = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) outer_weight = 1.0 / (outer - inner) @@ -107,7 +109,7 @@ def _star(a): return selem_square + selem_triangle -def _star_filter(m, n): +def _star_filter_kernel(m, n): c = m + m / 2 - n - n / 2 outer_star = _star(m) inner_star = np.zeros((outer_star.shape)) From 5f46fd01be76ac8c7252ab35ba58609e81ad0ca0 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 6 Aug 2013 19:45:44 +0530 Subject: [PATCH 46/77] Correcting a bug in Line Suppression --- skimage/feature/censure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index c04f7d24..726af80d 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -124,7 +124,7 @@ def _star_filter_kernel(m, n): def _suppress_lines(feature_mask, image, sigma, line_threshold): Axx, Axy, Ayy = _compute_auto_correlation(image, sigma) feature_mask[(Axx + Ayy) * (Axx + Ayy) - < line_threshold * (Axx * Ayy - Axy * Axy)] = 0 + > line_threshold * (Axx * Ayy - Axy * Axy)] = 0 return feature_mask From 54247e3a2667c757b3d358f21ccb20917313ab5a Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 6 Aug 2013 20:29:49 +0530 Subject: [PATCH 47/77] Putting constant variables outside the definition --- skimage/feature/censure.py | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 726af80d..6178335c 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -9,18 +9,20 @@ from skimage.morphology import convex_hull_image from skimage.feature.censure_cy import _censure_dob_loop +OCTAGON_OUTER_SHAPE = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), + (15, 10)] +OCTAGON_INNER_SHAPE = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] + +STAR_SHAPE = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] +STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), + (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] + + def _get_filtered_image(image, n_scales, mode): scales = np.zeros((image.shape[0], image.shape[1], n_scales), dtype=np.double) - OCTAGON_OUTER_SHAPE = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), - (15, 10)] - OCTAGON_INNER_SHAPE = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] - - STAR_SHAPE = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] - STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), - (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] if mode == 'DoB': for i in range(n_scales): n = i + 1 From 91bb4baaf8da6b288d0028662b65ed69fceec2bd Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 6 Aug 2013 20:36:05 +0530 Subject: [PATCH 48/77] Making the docs more explicit --- skimage/feature/_brief.py | 6 +++--- skimage/feature/censure.py | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/skimage/feature/_brief.py b/skimage/feature/_brief.py index 22b8e75a..35ee04bc 100644 --- a/skimage/feature/_brief.py +++ b/skimage/feature/_brief.py @@ -16,7 +16,7 @@ def brief(image, keypoints, descriptor_size=256, mode='normal', patch_size=49, image : 2D ndarray Input image. keypoints : (P, 2) ndarray - Array of keypoint locations. + Array of keypoint locations in the format (row, col). descriptor_size : int Size of BRIEF descriptor about each keypoint. Sizes 128, 256 and 512 preferred by the authors. Default is 256. @@ -44,8 +44,8 @@ def brief(image, keypoints, descriptor_size=256, mode='normal', patch_size=49, (i, j) either being True or False representing the outcome of Intensity comparison about ith keypoint on jth decision pixel-pair. keypoints : (Q, 2) ndarray - Keypoints after removing out those that are near border. - Returned only if return_keypoints is True. + Location i.e. (row, col) of keypoints after removing out those that + are near border. References ---------- diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 6178335c..6e403e49 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -160,7 +160,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, Returns ------- keypoints : (N, 2) array - Location of extracted keypoints. + Location of the extracted keypoints in the (row, col) format. scale : (N, 1) array The corresponding scale of the N extracted keypoints. From 14411590bf869ed89503db99fe187d8dd4b2b6a7 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 6 Aug 2013 20:38:08 +0530 Subject: [PATCH 49/77] Correcting the imports --- skimage/feature/tests/test_brief.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/skimage/feature/tests/test_brief.py b/skimage/feature/tests/test_brief.py index b78270f7..1d26cbbd 100644 --- a/skimage/feature/tests/test_brief.py +++ b/skimage/feature/tests/test_brief.py @@ -2,9 +2,9 @@ import numpy as np from numpy.testing import assert_array_equal, assert_raises from skimage import data from skimage import transform as tf -from skimage.feature.corner import corner_peaks, corner_harris from skimage.color import rgb2gray -from skimage.feature import brief, match_keypoints_brief +from skimage.feature import (brief, match_keypoints_brief, corner_peaks, + corner_harris) def test_brief_color_image_unsupported_error(): From c3e18b0aeeb945cf1a73422ebc25d2f09e72793c Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 6 Aug 2013 20:46:41 +0530 Subject: [PATCH 50/77] Assigning boolean to the feature mask in line suppression --- skimage/feature/censure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 6e403e49..774250bd 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -126,7 +126,7 @@ def _star_filter_kernel(m, n): def _suppress_lines(feature_mask, image, sigma, line_threshold): Axx, Axy, Ayy = _compute_auto_correlation(image, sigma) feature_mask[(Axx + Ayy) * (Axx + Ayy) - > line_threshold * (Axx * Ayy - Axy * Axy)] = 0 + > line_threshold * (Axx * Ayy - Axy * Axy)] = False return feature_mask From 4b0b342eadfcf17820de758524efd8fbf35073c1 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 7 Aug 2013 14:55:32 +0530 Subject: [PATCH 51/77] Actively filtering border keypoints on all scales --- skimage/feature/censure.py | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 774250bd..79765ccf 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -204,6 +204,22 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, feature_mask[:, :, i] = _suppress_lines(feature_mask[:, :, i], image, (1 + i / 3.0), line_threshold) + if mode == 'Octagon': + for i in range(1, n_scales - 1): + c = (OCTAGON_OUTER_SHAPE[i][0] - 1) / 2 + OCTAGON_OUTER_SHAPE[i][1] + feature_mask[:c, :, i] = False + feature_mask[:, :c, i] = False + feature_mask[-c:, :, i] = False + feature_mask[:, -c:, i] = False + + elif mode == 'STAR': + for i in range(1, n_scales - 1): + c = STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] / 2 + feature_mask[:c, :, i] = False + feature_mask[:, :c, i] = False + feature_mask[-c:, :, i] = False + feature_mask[:, -c:, i] = False + rows, cols, scales = np.nonzero(feature_mask[..., 1:n_scales - 1]) keypoints = np.column_stack([rows, cols]) scales = scales + 2 From 3752f19cf92c3afa6e988d0df3e1f9e84447539a Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 8 Aug 2013 05:57:55 +0530 Subject: [PATCH 52/77] Adding tests for all the modes in censure --- skimage/feature/tests/test_censure.py | 72 +++++++++++++++++++++++++++ 1 file changed, 72 insertions(+) create mode 100644 skimage/feature/tests/test_censure.py diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py new file mode 100644 index 00000000..2a582a9a --- /dev/null +++ b/skimage/feature/tests/test_censure.py @@ -0,0 +1,72 @@ +import numpy as np +from numpy.testing import assert_array_equal, assert_raises +from skimage.data import moon +from skimage.feature import censure_keypoints + + +def test_censure_keypoints_color_image_unsupported_error(): + """Censure keypoints can be extracted from gray-scale images only.""" + img = np.zeros((20, 20, 3)) + assert_raises(ValueError, censure_keypoints, img) + + +def test_censure_keypoints_moon_image_DoB(): + """Verify the actual Censure keypoints and their corresponding scale with + the expected values for DoB filter.""" + img = moon() + actual_kp_DoB, actual_scale = censure_keypoints(img, 7, 'DoB', 0.15) + expected_kp_DoB = np.array([[ 4, 507], + [ 8, 503], + [ 12, 499], + [ 21, 497], + [ 36, 46], + [119, 350], + [185, 177], + [287, 250], + [357, 239], + [463, 116], + [464, 132], + [467, 260]]) + expected_scale = np.array([2, 4, 6, 3, 4, 4, 2, 2, 3, 2, 2, 2]) + assert_array_equal(expected_kp_DoB, actual_kp_DoB) + assert_array_equal(expected_scale, actual_scale) + + +def test_censure_keypoints_moon_image_Octagon(): + """Verify the actual Censure keypoints and their corresponding scale with + the expected values for Octagon filter.""" + img = moon() + actual_kp_Octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', 0.15) + expected_kp_Octagon = np.array([[ 21, 496], + [ 35, 46], + [287, 250], + [356, 239], + [463, 116]]) + expected_scale = np.array([3, 4, 2, 2, 2]) + assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) + assert_array_equal(expected_scale, actual_scale) + + +def test_censure_keypoints_moon_image_STAR(): + """Verify the actual Censure keypoints and their corresponding scale with + the expected values for STAR filter.""" + img = moon() + actual_kp_STAR, actual_scale = censure_keypoints(img, 7, 'STAR', 0.15) + expected_kp_STAR = np.array([[ 21, 497], + [ 36, 46], + [117, 356], + [185, 177], + [260, 227], + [287, 250], + [357, 239], + [451, 281], + [463, 116], + [467, 260]]) + expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2]) + assert_array_equal(expected_kp_STAR, actual_kp_STAR) + assert_array_equal(expected_scale, actual_scale) + + +if __name__ == '__main__': + from numpy import testing + testing.run_module_suite() From c9a93a1ce39fe2a5fdd915b46fb9658e454e807f Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 8 Aug 2013 07:51:03 +0530 Subject: [PATCH 53/77] Debugging Travis error --- skimage/feature/tests/test_censure.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 2a582a9a..3d4e1a1c 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -28,6 +28,8 @@ def test_censure_keypoints_moon_image_DoB(): [464, 132], [467, 260]]) expected_scale = np.array([2, 4, 6, 3, 4, 4, 2, 2, 3, 2, 2, 2]) + print actual_kp_DoB + print actual_scale assert_array_equal(expected_kp_DoB, actual_kp_DoB) assert_array_equal(expected_scale, actual_scale) @@ -43,6 +45,8 @@ def test_censure_keypoints_moon_image_Octagon(): [356, 239], [463, 116]]) expected_scale = np.array([3, 4, 2, 2, 2]) + print actual_kp_Octagon + print actual_scale assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) assert_array_equal(expected_scale, actual_scale) @@ -63,6 +67,8 @@ def test_censure_keypoints_moon_image_STAR(): [463, 116], [467, 260]]) expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2]) + print actual_kp_STAR + print actual_scale assert_array_equal(expected_kp_STAR, actual_kp_STAR) assert_array_equal(expected_scale, actual_scale) From fb4050ed37168bb283410ba34a896d96d62c60cf Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 8 Aug 2013 09:57:01 +0530 Subject: [PATCH 54/77] Reverting back to last but one commit --- skimage/feature/tests/test_censure.py | 9 +++------ 1 file changed, 3 insertions(+), 6 deletions(-) diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 3d4e1a1c..deebc9c1 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -28,8 +28,7 @@ def test_censure_keypoints_moon_image_DoB(): [464, 132], [467, 260]]) expected_scale = np.array([2, 4, 6, 3, 4, 4, 2, 2, 3, 2, 2, 2]) - print actual_kp_DoB - print actual_scale + assert_array_equal(expected_kp_DoB, actual_kp_DoB) assert_array_equal(expected_scale, actual_scale) @@ -45,8 +44,7 @@ def test_censure_keypoints_moon_image_Octagon(): [356, 239], [463, 116]]) expected_scale = np.array([3, 4, 2, 2, 2]) - print actual_kp_Octagon - print actual_scale + assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) assert_array_equal(expected_scale, actual_scale) @@ -67,8 +65,7 @@ def test_censure_keypoints_moon_image_STAR(): [463, 116], [467, 260]]) expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2]) - print actual_kp_STAR - print actual_scale + assert_array_equal(expected_kp_STAR, actual_kp_STAR) assert_array_equal(expected_scale, actual_scale) From 963f681e48d1722c67cbda25892e4506e711e444 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 8 Aug 2013 11:26:53 +0530 Subject: [PATCH 55/77] Making the code compatible with Python 3 --- skimage/feature/censure.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 79765ccf..1a2f495c 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -83,7 +83,7 @@ def _octagon_filter_kernel(mo, no, mi, ni): inner = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) outer_weight = 1.0 / (outer - inner) inner_weight = 1.0 / inner - c = ((mo + 2 * no) - (mi + 2 * ni)) / 2 + c = ((mo + 2 * no) - (mi + 2 * ni)) // 2 outer_oct = _oct(mo, no) inner_oct = np.zeros((mo + 2 * no, mo + 2 * no)) inner_oct[c: -c, c: -c] = _oct(mi, ni) @@ -98,13 +98,13 @@ def _star(a): bfilter[:] = 1 return bfilter m = 2 * a + 1 - n = a / 2 + n = a // 2 selem_square = np.zeros((m + 2 * n, m + 2 * n), dtype=np.uint8) selem_square[n: m + n, n: m + n] = 1 selem_triangle = np.zeros((m + 2 * n, m + 2 * n), dtype=np.uint8) - selem_triangle[(m + 2 * n - 1) / 2, 0] = 1 - selem_triangle[(m + 1) / 2, n - 1] = 1 - selem_triangle[(m + 4 * n - 3) / 2, n - 1] = 1 + selem_triangle[(m + 2 * n - 1) // 2, 0] = 1 + selem_triangle[(m + 1) // 2, n - 1] = 1 + selem_triangle[(m + 4 * n - 3) // 2, n - 1] = 1 selem_triangle = convex_hull_image(selem_triangle).astype(int) selem_triangle += (selem_triangle[:, ::-1] + selem_triangle.T + selem_triangle.T[::-1, :]) @@ -112,7 +112,7 @@ def _star(a): def _star_filter_kernel(m, n): - c = m + m / 2 - n - n / 2 + c = m + m // 2 - n - n // 2 outer_star = _star(m) inner_star = np.zeros((outer_star.shape)) inner_star[c: -c, c: -c] = _star(n) @@ -206,7 +206,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, if mode == 'Octagon': for i in range(1, n_scales - 1): - c = (OCTAGON_OUTER_SHAPE[i][0] - 1) / 2 + OCTAGON_OUTER_SHAPE[i][1] + c = (OCTAGON_OUTER_SHAPE[i][0] - 1) // 2 + OCTAGON_OUTER_SHAPE[i][1] feature_mask[:c, :, i] = False feature_mask[:, :c, i] = False feature_mask[-c:, :, i] = False @@ -214,7 +214,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, elif mode == 'STAR': for i in range(1, n_scales - 1): - c = STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] / 2 + c = STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] // 2 feature_mask[:c, :, i] = False feature_mask[:, :c, i] = False feature_mask[-c:, :, i] = False From 1227507c9e42f44e5f8de08664315f7d4641063c Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 8 Aug 2013 18:08:49 +0530 Subject: [PATCH 56/77] Actively filtering border keypoints for all the scales part 2 --- skimage/feature/censure.py | 39 +++++++++++++++------------ skimage/feature/tests/test_censure.py | 34 ++++++++++++----------- skimage/feature/util.py | 13 +++++---- 3 files changed, 49 insertions(+), 37 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 1a2f495c..7500106d 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -5,6 +5,7 @@ from skimage.transform import integral_image from skimage.feature.corner import _compute_auto_correlation from skimage.util import img_as_float from skimage.morphology import convex_hull_image +from skimage.feature.util import _remove_border_keypoints from skimage.feature.censure_cy import _censure_dob_loop @@ -204,24 +205,28 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, feature_mask[:, :, i] = _suppress_lines(feature_mask[:, :, i], image, (1 + i / 3.0), line_threshold) - if mode == 'Octagon': - for i in range(1, n_scales - 1): - c = (OCTAGON_OUTER_SHAPE[i][0] - 1) // 2 + OCTAGON_OUTER_SHAPE[i][1] - feature_mask[:c, :, i] = False - feature_mask[:, :c, i] = False - feature_mask[-c:, :, i] = False - feature_mask[:, -c:, i] = False - - elif mode == 'STAR': - for i in range(1, n_scales - 1): - c = STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] // 2 - feature_mask[:c, :, i] = False - feature_mask[:, :c, i] = False - feature_mask[-c:, :, i] = False - feature_mask[:, -c:, i] = False - rows, cols, scales = np.nonzero(feature_mask[..., 1:n_scales - 1]) keypoints = np.column_stack([rows, cols]) scales = scales + 2 - return keypoints, scales + if mode == 'DoB': + return keypoints, scales + + filtered_keypoints = np.empty((0, 2), dtype=np.int32) + filtered_scales = np.empty((0), dtype=np.int32) + + if mode == 'Octagon': + for i in range(2, n_scales): + c = (OCTAGON_OUTER_SHAPE[i - 1][0] - 1) // 2 + OCTAGON_OUTER_SHAPE[i - 1][1] + filtered_keypoints_for_scale = _remove_border_keypoints(image, keypoints[scales == i], c) + filtered_keypoints = np.vstack((filtered_keypoints, filtered_keypoints_for_scale)) + filtered_scales = np.hstack((filtered_scales, np.asarray(len(filtered_keypoints_for_scale) * [i], dtype=np.int32))) + + elif mode == 'STAR': + for i in range(2, n_scales): + c = STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] // 2 + filtered_keypoints_for_scale = _remove_border_keypoints(image, keypoints[scales == i], c) + filtered_keypoints = np.vstack((filtered_keypoints, filtered_keypoints_for_scale)) + filtered_scales = np.hstack((filtered_scales, np.asarray(len(filtered_keypoints_for_scale) * [i], dtype=np.int32))) + + return filtered_keypoints, filtered_scales diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index deebc9c1..4f1ddcfb 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -28,7 +28,8 @@ def test_censure_keypoints_moon_image_DoB(): [464, 132], [467, 260]]) expected_scale = np.array([2, 4, 6, 3, 4, 4, 2, 2, 3, 2, 2, 2]) - + print actual_scale + print actual_kp_DoB assert_array_equal(expected_kp_DoB, actual_kp_DoB) assert_array_equal(expected_scale, actual_scale) @@ -38,13 +39,15 @@ def test_censure_keypoints_moon_image_Octagon(): the expected values for Octagon filter.""" img = moon() actual_kp_Octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', 0.15) - expected_kp_Octagon = np.array([[ 21, 496], - [ 35, 46], - [287, 250], + expected_kp_Octagon = np.array([[287, 250], [356, 239], - [463, 116]]) - expected_scale = np.array([3, 4, 2, 2, 2]) + [463, 116], + [ 21, 496], + [ 35, 46]]) + expected_scale = np.array([2, 2, 2, 3, 4], dtype=np.int32) + print actual_scale + print actual_kp_Octagon assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) assert_array_equal(expected_scale, actual_scale) @@ -54,18 +57,19 @@ def test_censure_keypoints_moon_image_STAR(): the expected values for STAR filter.""" img = moon() actual_kp_STAR, actual_scale = censure_keypoints(img, 7, 'STAR', 0.15) - expected_kp_STAR = np.array([[ 21, 497], - [ 36, 46], - [117, 356], - [185, 177], - [260, 227], + expected_kp_STAR = np.array([[185, 177], [287, 250], + [463, 116], + [467, 260], + [ 21, 497], + [ 36, 46], + [260, 227], [357, 239], [451, 281], - [463, 116], - [467, 260]]) - expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2]) - + [117, 356]]) + expected_scale = np.array([2, 2, 2, 2, 3, 3, 3, 3, 5, 6], dtype=np.int32) + print actual_scale + print actual_kp_STAR assert_array_equal(expected_kp_STAR, actual_kp_STAR) assert_array_equal(expected_scale, actual_scale) diff --git a/skimage/feature/util.py b/skimage/feature/util.py index aec4dfc8..f327cff5 100644 --- a/skimage/feature/util.py +++ b/skimage/feature/util.py @@ -1,3 +1,4 @@ +import numpy as np def _remove_border_keypoints(image, keypoints, dist): @@ -5,11 +6,13 @@ def _remove_border_keypoints(image, keypoints, dist): width = image.shape[0] height = image.shape[1] - keypoints = keypoints[(dist - 1 < keypoints[:, 0]) - & (keypoints[:, 0] < width - dist + 1) - & (dist - 1 < keypoints[:, 1]) - & (keypoints[:, 1] < height - dist + 1)] - + try: + keypoints = keypoints[(dist - 1 < keypoints[:, 0]) + & (keypoints[:, 0] < width - dist + 1) + & (dist - 1 < keypoints[:, 1]) + & (keypoints[:, 1] < height - dist + 1)] + except IndexError: + return np.empty((0, 2), dtype=np.int32) return keypoints From 6fc84dc7e6cbcd69fbb16559f974ca8ae248a8f4 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 8 Aug 2013 23:50:29 +0530 Subject: [PATCH 57/77] Removing print statements --- skimage/feature/tests/test_censure.py | 9 +++------ 1 file changed, 3 insertions(+), 6 deletions(-) diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 4f1ddcfb..1822ec0a 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -28,8 +28,7 @@ def test_censure_keypoints_moon_image_DoB(): [464, 132], [467, 260]]) expected_scale = np.array([2, 4, 6, 3, 4, 4, 2, 2, 3, 2, 2, 2]) - print actual_scale - print actual_kp_DoB + assert_array_equal(expected_kp_DoB, actual_kp_DoB) assert_array_equal(expected_scale, actual_scale) @@ -46,8 +45,7 @@ def test_censure_keypoints_moon_image_Octagon(): [ 35, 46]]) expected_scale = np.array([2, 2, 2, 3, 4], dtype=np.int32) - print actual_scale - print actual_kp_Octagon + assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) assert_array_equal(expected_scale, actual_scale) @@ -68,8 +66,7 @@ def test_censure_keypoints_moon_image_STAR(): [451, 281], [117, 356]]) expected_scale = np.array([2, 2, 2, 2, 3, 3, 3, 3, 5, 6], dtype=np.int32) - print actual_scale - print actual_kp_STAR + assert_array_equal(expected_kp_STAR, actual_kp_STAR) assert_array_equal(expected_scale, actual_scale) From 4f109d18e2542d8308366cc82678d53b57bbe7fb Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 8 Aug 2013 23:55:44 +0530 Subject: [PATCH 58/77] Returning mask in _remove_border_keypoints --- skimage/feature/util.py | 15 ++++++--------- 1 file changed, 6 insertions(+), 9 deletions(-) diff --git a/skimage/feature/util.py b/skimage/feature/util.py index f327cff5..c644af50 100644 --- a/skimage/feature/util.py +++ b/skimage/feature/util.py @@ -1,4 +1,3 @@ -import numpy as np def _remove_border_keypoints(image, keypoints, dist): @@ -6,14 +5,12 @@ def _remove_border_keypoints(image, keypoints, dist): width = image.shape[0] height = image.shape[1] - try: - keypoints = keypoints[(dist - 1 < keypoints[:, 0]) - & (keypoints[:, 0] < width - dist + 1) - & (dist - 1 < keypoints[:, 1]) - & (keypoints[:, 1] < height - dist + 1)] - except IndexError: - return np.empty((0, 2), dtype=np.int32) - return keypoints + keypoints_filtering_mask = (dist - 1 < keypoints[:, 0] + & keypoints[:, 0] < width - dist + 1 + & dist - 1 < keypoints[:, 1] + & keypoints[:, 1] < height - dist + 1) + + return keypoints_filtering_mask def pairwise_hamming_distance(array1, array2): From 62eb4ef998fbdf830b9abfb8454ef68a26fd8507 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Fri, 9 Aug 2013 00:42:58 +0530 Subject: [PATCH 59/77] Filtering out border keypoints using masking --- skimage/feature/_brief.py | 2 +- skimage/feature/censure.py | 13 ++++--------- skimage/feature/tests/test_censure.py | 25 +++++++++++++------------ skimage/feature/util.py | 8 ++++---- 4 files changed, 22 insertions(+), 26 deletions(-) diff --git a/skimage/feature/_brief.py b/skimage/feature/_brief.py index 35ee04bc..73230188 100644 --- a/skimage/feature/_brief.py +++ b/skimage/feature/_brief.py @@ -142,7 +142,7 @@ def brief(image, keypoints, descriptor_size=256, mode='normal', patch_size=49, # Removing keypoints that are within (patch_size / 2) distance from the # image border - keypoints = _remove_border_keypoints(image, keypoints, patch_size // 2) + keypoints = keypoints[_remove_border_keypoints(image, keypoints, patch_size // 2)] keypoints = np.ascontiguousarray(keypoints) descriptors = np.zeros((keypoints.shape[0], descriptor_size), dtype=bool, diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 7500106d..62f2976c 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -212,21 +212,16 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, if mode == 'DoB': return keypoints, scales - filtered_keypoints = np.empty((0, 2), dtype=np.int32) - filtered_scales = np.empty((0), dtype=np.int32) + cumulative_mask = np.zeros(keypoints.shape[0], dtype=np.bool) if mode == 'Octagon': for i in range(2, n_scales): c = (OCTAGON_OUTER_SHAPE[i - 1][0] - 1) // 2 + OCTAGON_OUTER_SHAPE[i - 1][1] - filtered_keypoints_for_scale = _remove_border_keypoints(image, keypoints[scales == i], c) - filtered_keypoints = np.vstack((filtered_keypoints, filtered_keypoints_for_scale)) - filtered_scales = np.hstack((filtered_scales, np.asarray(len(filtered_keypoints_for_scale) * [i], dtype=np.int32))) + cumulative_mask = cumulative_mask | (_remove_border_keypoints(image, keypoints, c) & (scales == i)) elif mode == 'STAR': for i in range(2, n_scales): c = STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] // 2 - filtered_keypoints_for_scale = _remove_border_keypoints(image, keypoints[scales == i], c) - filtered_keypoints = np.vstack((filtered_keypoints, filtered_keypoints_for_scale)) - filtered_scales = np.hstack((filtered_scales, np.asarray(len(filtered_keypoints_for_scale) * [i], dtype=np.int32))) + cumulative_mask = cumulative_mask | (_remove_border_keypoints(image, keypoints, c) & (scales == i)) - return filtered_keypoints, filtered_scales + return keypoints[cumulative_mask], scales[cumulative_mask] diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 1822ec0a..1f8a01dd 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -38,13 +38,13 @@ def test_censure_keypoints_moon_image_Octagon(): the expected values for Octagon filter.""" img = moon() actual_kp_Octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', 0.15) - expected_kp_Octagon = np.array([[287, 250], + expected_kp_Octagon = np.array([[ 21, 496], + [ 35, 46], + [287, 250], [356, 239], - [463, 116], - [ 21, 496], - [ 35, 46]]) + [463, 116]]) - expected_scale = np.array([2, 2, 2, 3, 4], dtype=np.int32) + expected_scale = np.array([3, 4, 2, 2, 2], dtype=np.int32) assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) assert_array_equal(expected_scale, actual_scale) @@ -55,17 +55,18 @@ def test_censure_keypoints_moon_image_STAR(): the expected values for STAR filter.""" img = moon() actual_kp_STAR, actual_scale = censure_keypoints(img, 7, 'STAR', 0.15) - expected_kp_STAR = np.array([[185, 177], - [287, 250], - [463, 116], - [467, 260], - [ 21, 497], + expected_kp_STAR = np.array([[ 21, 497], [ 36, 46], + [117, 356], + [185, 177], [260, 227], + [287, 250], [357, 239], [451, 281], - [117, 356]]) - expected_scale = np.array([2, 2, 2, 2, 3, 3, 3, 3, 5, 6], dtype=np.int32) + [463, 116], + [467, 260]]) + + expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2], dtype=np.intp) assert_array_equal(expected_kp_STAR, actual_kp_STAR) assert_array_equal(expected_scale, actual_scale) diff --git a/skimage/feature/util.py b/skimage/feature/util.py index c644af50..513fa05d 100644 --- a/skimage/feature/util.py +++ b/skimage/feature/util.py @@ -5,10 +5,10 @@ def _remove_border_keypoints(image, keypoints, dist): width = image.shape[0] height = image.shape[1] - keypoints_filtering_mask = (dist - 1 < keypoints[:, 0] - & keypoints[:, 0] < width - dist + 1 - & dist - 1 < keypoints[:, 1] - & keypoints[:, 1] < height - dist + 1) + keypoints_filtering_mask = ((dist - 1 < keypoints[:, 0]) & + (keypoints[:, 0] < width - dist + 1) & + (dist - 1 < keypoints[:, 1]) & + (keypoints[:, 1] < height - dist + 1)) return keypoints_filtering_mask From 1bd261c23fccd8a3deb0cfe4c3840a9ef7f8e67d Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Fri, 9 Aug 2013 00:44:18 +0530 Subject: [PATCH 60/77] Replacing int32 by intp --- skimage/feature/tests/test_censure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 1f8a01dd..41d1886e 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -44,7 +44,7 @@ def test_censure_keypoints_moon_image_Octagon(): [356, 239], [463, 116]]) - expected_scale = np.array([3, 4, 2, 2, 2], dtype=np.int32) + expected_scale = np.array([3, 4, 2, 2, 2], dtype=np.intp) assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) assert_array_equal(expected_scale, actual_scale) From 157c22bbc9bf592b003f1562ed28bd1050d5234e Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Fri, 9 Aug 2013 14:54:41 +0530 Subject: [PATCH 61/77] Changing _remove_border_keypoints to _mask_border_keypoints --- skimage/feature/_brief.py | 4 ++-- skimage/feature/censure.py | 6 +++--- skimage/feature/util.py | 2 +- 3 files changed, 6 insertions(+), 6 deletions(-) diff --git a/skimage/feature/_brief.py b/skimage/feature/_brief.py index 73230188..27a8d085 100644 --- a/skimage/feature/_brief.py +++ b/skimage/feature/_brief.py @@ -2,7 +2,7 @@ import numpy as np from scipy.ndimage.filters import gaussian_filter from ..util import img_as_float -from .util import _remove_border_keypoints, pairwise_hamming_distance +from .util import _mask_border_keypoints, pairwise_hamming_distance from ._brief_cy import _brief_loop @@ -142,7 +142,7 @@ def brief(image, keypoints, descriptor_size=256, mode='normal', patch_size=49, # Removing keypoints that are within (patch_size / 2) distance from the # image border - keypoints = keypoints[_remove_border_keypoints(image, keypoints, patch_size // 2)] + keypoints = keypoints[_mask_border_keypoints(image, keypoints, patch_size // 2)] keypoints = np.ascontiguousarray(keypoints) descriptors = np.zeros((keypoints.shape[0], descriptor_size), dtype=bool, diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 62f2976c..e3264436 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -5,7 +5,7 @@ from skimage.transform import integral_image from skimage.feature.corner import _compute_auto_correlation from skimage.util import img_as_float from skimage.morphology import convex_hull_image -from skimage.feature.util import _remove_border_keypoints +from skimage.feature.util import _mask_border_keypoints from skimage.feature.censure_cy import _censure_dob_loop @@ -217,11 +217,11 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, if mode == 'Octagon': for i in range(2, n_scales): c = (OCTAGON_OUTER_SHAPE[i - 1][0] - 1) // 2 + OCTAGON_OUTER_SHAPE[i - 1][1] - cumulative_mask = cumulative_mask | (_remove_border_keypoints(image, keypoints, c) & (scales == i)) + cumulative_mask = cumulative_mask | (_mask_border_keypoints(image, keypoints, c) & (scales == i)) elif mode == 'STAR': for i in range(2, n_scales): c = STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] // 2 - cumulative_mask = cumulative_mask | (_remove_border_keypoints(image, keypoints, c) & (scales == i)) + cumulative_mask = cumulative_mask | (_mask_border_keypoints(image, keypoints, c) & (scales == i)) return keypoints[cumulative_mask], scales[cumulative_mask] diff --git a/skimage/feature/util.py b/skimage/feature/util.py index 513fa05d..eb3817e8 100644 --- a/skimage/feature/util.py +++ b/skimage/feature/util.py @@ -1,6 +1,6 @@ -def _remove_border_keypoints(image, keypoints, dist): +def _mask_border_keypoints(image, keypoints, dist): """Removes keypoints that are within dist pixels from the image border.""" width = image.shape[0] height = image.shape[1] From 0b37611df60bf1da7542324f77e9d52093aab824 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Johannes=20Sch=C3=B6nberger?= Date: Sun, 11 Aug 2013 11:02:24 +0200 Subject: [PATCH 62/77] Fix several bugs in DoB method and improve overall code quality --- skimage/feature/censure.py | 76 +++++++++++++++------------ skimage/feature/censure_cy.pyx | 21 ++++---- skimage/feature/tests/test_censure.py | 15 +++--- 3 files changed, 59 insertions(+), 53 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index e3264436..583e3d0a 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -24,7 +24,16 @@ def _get_filtered_image(image, n_scales, mode): scales = np.zeros((image.shape[0], image.shape[1], n_scales), dtype=np.double) - if mode == 'DoB': + if mode == 'dob': + + # make scales[:, :, i] contiguous memory block + item_size = scales.itemsize + scales.strides = (item_size * scales.shape[0], + item_size, + item_size * scales.shape[0] * scales.shape[1]) + + integral_img = integral_image(image) + for i in range(n_scales): n = i + 1 @@ -34,33 +43,29 @@ def _get_filtered_image(image, n_scales, mode): inner_weight = (1.0 / (2 * n + 1)**2) outer_weight = (1.0 / (12 * n**2 + 4 * n)) - integral_img = integral_image(image) - - filtered_image = np.zeros(image.shape) - _censure_dob_loop(image, n, integral_img, filtered_image, + _censure_dob_loop(n, integral_img, scales[:, :, i], inner_weight, outer_weight) - scales[:, :, i] = filtered_image - # NOTE : For the Octagon shaped filter, we implemented and evaluated the # slanted integral image based image filtering but the performance was # more or less equal to image filtering using # scipy.ndimage.filters.convolve(). Hence we have decided to use the # later for a much cleaner implementation. - elif mode == 'Octagon': + elif mode == 'octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 for i in range(n_scales): + mo, no = OCTAGON_OUTER_SHAPE[i] + mi, ni = OCTAGON_INNER_SHAPE[i] scales[:, :, i] = convolve(image, - _octagon_filter_kernel(OCTAGON_OUTER_SHAPE[i][0], - OCTAGON_OUTER_SHAPE[i][1], OCTAGON_INNER_SHAPE[i][0], - OCTAGON_INNER_SHAPE[i][1])) - else: + _octagon_filter_kernel(mo, no, mi, ni)) + elif mode == 'star': for i in range(n_scales): + m = STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] + n = STAR_SHAPE[STAR_FILTER_SHAPE[i][1]] scales[:, :, i] = convolve(image, - _star_filter_kernel(STAR_SHAPE[STAR_FILTER_SHAPE[i][0]], - STAR_SHAPE[STAR_FILTER_SHAPE[i][1]])) + _star_filter_kernel(m, n)) return scales @@ -115,7 +120,7 @@ def _star(a): def _star_filter_kernel(m, n): c = m + m // 2 - n - n // 2 outer_star = _star(m) - inner_star = np.zeros((outer_star.shape)) + inner_star = np.zeros_like(outer_star) inner_star[c: -c, c: -c] = _star(n) outer_weight = 1.0 / (np.sum(outer_star - inner_star)) inner_weight = 1.0 / np.sum(inner_star) @@ -128,7 +133,6 @@ def _suppress_lines(feature_mask, image, sigma, line_threshold): Axx, Axy, Ayy = _compute_auto_correlation(image, sigma) feature_mask[(Axx + Ayy) * (Axx + Ayy) > line_threshold * (Axx * Ayy - Axy * Axy)] = False - return feature_mask def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, @@ -141,19 +145,15 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, ---------- image : 2D ndarray Input image. - n_scales : positive integer Number of scales to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last. - mode : ('DoB', 'Octagon', 'STAR') Type of bilevel filter used to get the scales of input image. Possible values are 'DoB', 'Octagon' and 'STAR'. - non_max_threshold : float Threshold value used to suppress maximas and minimas with a weak magnitude response obtained after Non-Maximal Suppression. - line_threshold : float Threshold for rejecting interest points which have ratio of principal curvatures greater than this value. @@ -162,8 +162,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, ------- keypoints : (N, 2) array Location of the extracted keypoints in the (row, col) format. - - scale : (N, 1) array + scales : (N, 1) array The corresponding scale of the N extracted keypoints. References @@ -183,10 +182,15 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, image = np.squeeze(image) if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") - image = img_as_float(image) + image = img_as_float(image) image = np.ascontiguousarray(image) + mode = mode.lower() + + if mode not in ('dob', 'octagon', 'star'): + raise ValueError('Mode must be one of "DoB", "Octagon", "STAR".') + # Generating all the scales filter_response = _get_filtered_image(image, n_scales, mode) @@ -199,29 +203,33 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, feature_mask[filter_response < non_max_threshold] = False for i in range(1, n_scales - 1): - # sigma = (window_size - 1) / 6.0 + # sigma = (window_size - 1) / 6.0, so the window covers > 99% of the + # kernel's distribution # window_size = 7 + 2 * i # Hence sigma = 1 + i / 3.0 - feature_mask[:, :, i] = _suppress_lines(feature_mask[:, :, i], image, - (1 + i / 3.0), line_threshold) + _suppress_lines(feature_mask[:, :, i], image, + (1 + i / 3.0), line_threshold) rows, cols, scales = np.nonzero(feature_mask[..., 1:n_scales - 1]) keypoints = np.column_stack([rows, cols]) scales = scales + 2 - if mode == 'DoB': + if mode == 'dob': return keypoints, scales cumulative_mask = np.zeros(keypoints.shape[0], dtype=np.bool) - if mode == 'Octagon': + if mode == 'octagon': for i in range(2, n_scales): - c = (OCTAGON_OUTER_SHAPE[i - 1][0] - 1) // 2 + OCTAGON_OUTER_SHAPE[i - 1][1] - cumulative_mask = cumulative_mask | (_mask_border_keypoints(image, keypoints, c) & (scales == i)) - - elif mode == 'STAR': + c = (OCTAGON_OUTER_SHAPE[i - 1][0] - 1) // 2 \ + + OCTAGON_OUTER_SHAPE[i - 1][1] + cumulative_mask |= _mask_border_keypoints(image, keypoints, c) \ + & (scales == i) + elif mode == 'star': for i in range(2, n_scales): - c = STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] + STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] // 2 - cumulative_mask = cumulative_mask | (_mask_border_keypoints(image, keypoints, c) & (scales == i)) + c = STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] \ + + STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] // 2 + cumulative_mask |= _mask_border_keypoints(image, keypoints, c) \ + & (scales == i) return keypoints[cumulative_mask], scales[cumulative_mask] diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index 93c7c142..cfd1260f 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -3,29 +3,28 @@ #cython: nonecheck=False #cython: wraparound=False -cimport numpy as cnp -import numpy as np - -def _censure_dob_loop(double[:, ::1] image, Py_ssize_t n, +def _censure_dob_loop(Py_ssize_t n, double[:, ::1] integral_img, double[:, ::1] filtered_image, double inner_weight, double outer_weight): cdef Py_ssize_t i, j cdef double inner, outer + cdef Py_ssize_t n2 = 2 * n + cdef double total_weight = inner_weight + outer_weight - for i in range(2 * n, image.shape[0] - 2 * n): - for j in range(2 * n, image.shape[1] - 2 * n): + for i in range(n2, integral_img.shape[0] - n2): + for j in range(n2, integral_img.shape[1] - n2): inner = (integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] - integral_img[i - n - 1, j + n]) - outer = (integral_img[i + 2 * n, j + 2 * n] - + integral_img[i - 2 * n - 1, j - 2 * n - 1] - - integral_img[i + 2 * n, j - 2 * n - 1] - - integral_img[i - 2 * n - 1, j + 2 * n]) + outer = (integral_img[i + n2, j + n2] + + integral_img[i - n2 - 1, j - n2 - 1] + - integral_img[i + n2, j - n2 - 1] + - integral_img[i - n2 - 1, j + n2]) filtered_image[i, j] = (outer_weight * outer - - (inner_weight + outer_weight) * inner) + - total_weight * inner) diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 41d1886e..14e97a39 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -1,6 +1,7 @@ import numpy as np from numpy.testing import assert_array_equal, assert_raises from skimage.data import moon +from skimage.util import img_as_ubyte from skimage.feature import censure_keypoints @@ -15,10 +16,7 @@ def test_censure_keypoints_moon_image_DoB(): the expected values for DoB filter.""" img = moon() actual_kp_DoB, actual_scale = censure_keypoints(img, 7, 'DoB', 0.15) - expected_kp_DoB = np.array([[ 4, 507], - [ 8, 503], - [ 12, 499], - [ 21, 497], + expected_kp_DoB = np.array([[ 21, 497], [ 36, 46], [119, 350], [185, 177], @@ -27,7 +25,7 @@ def test_censure_keypoints_moon_image_DoB(): [463, 116], [464, 132], [467, 260]]) - expected_scale = np.array([2, 4, 6, 3, 4, 4, 2, 2, 3, 2, 2, 2]) + expected_scale = np.array([3, 4, 4, 2, 2, 3, 2, 2, 2]) assert_array_equal(expected_kp_DoB, actual_kp_DoB) assert_array_equal(expected_scale, actual_scale) @@ -37,14 +35,15 @@ def test_censure_keypoints_moon_image_Octagon(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for Octagon filter.""" img = moon() - actual_kp_Octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', 0.15) + actual_kp_Octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', + 0.15) expected_kp_Octagon = np.array([[ 21, 496], [ 35, 46], [287, 250], [356, 239], [463, 116]]) - expected_scale = np.array([3, 4, 2, 2, 2], dtype=np.intp) + expected_scale = np.array([3, 4, 2, 2, 2]) assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) assert_array_equal(expected_scale, actual_scale) @@ -66,7 +65,7 @@ def test_censure_keypoints_moon_image_STAR(): [463, 116], [467, 260]]) - expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2], dtype=np.intp) + expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2]) assert_array_equal(expected_kp_STAR, actual_kp_STAR) assert_array_equal(expected_scale, actual_scale) From 7665281c4a76f2d0799bf154b4999853f401fe8b Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 12 Aug 2013 16:27:35 +0530 Subject: [PATCH 63/77] Lowercasing variables; removing unused import; replacing _remove by _mask --- skimage/feature/tests/test_censure.py | 21 ++++++++++----------- 1 file changed, 10 insertions(+), 11 deletions(-) diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 14e97a39..9acdd2e9 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -1,7 +1,6 @@ import numpy as np from numpy.testing import assert_array_equal, assert_raises from skimage.data import moon -from skimage.util import img_as_ubyte from skimage.feature import censure_keypoints @@ -11,12 +10,12 @@ def test_censure_keypoints_color_image_unsupported_error(): assert_raises(ValueError, censure_keypoints, img) -def test_censure_keypoints_moon_image_DoB(): +def test_censure_keypoints_moon_image_dob(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for DoB filter.""" img = moon() - actual_kp_DoB, actual_scale = censure_keypoints(img, 7, 'DoB', 0.15) - expected_kp_DoB = np.array([[ 21, 497], + actual_kp_dob, actual_scale = censure_keypoints(img, 7, 'DoB', 0.15) + expected_kp_dob = np.array([[ 21, 497], [ 36, 46], [119, 350], [185, 177], @@ -27,7 +26,7 @@ def test_censure_keypoints_moon_image_DoB(): [467, 260]]) expected_scale = np.array([3, 4, 4, 2, 2, 3, 2, 2, 2]) - assert_array_equal(expected_kp_DoB, actual_kp_DoB) + assert_array_equal(expected_kp_dob, actual_kp_dob) assert_array_equal(expected_scale, actual_scale) @@ -35,9 +34,9 @@ def test_censure_keypoints_moon_image_Octagon(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for Octagon filter.""" img = moon() - actual_kp_Octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', + actual_kp_octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', 0.15) - expected_kp_Octagon = np.array([[ 21, 496], + expected_kp_octagon = np.array([[ 21, 496], [ 35, 46], [287, 250], [356, 239], @@ -45,7 +44,7 @@ def test_censure_keypoints_moon_image_Octagon(): expected_scale = np.array([3, 4, 2, 2, 2]) - assert_array_equal(expected_kp_Octagon, actual_kp_Octagon) + assert_array_equal(expected_kp_octagon, actual_kp_octagon) assert_array_equal(expected_scale, actual_scale) @@ -53,8 +52,8 @@ def test_censure_keypoints_moon_image_STAR(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for STAR filter.""" img = moon() - actual_kp_STAR, actual_scale = censure_keypoints(img, 7, 'STAR', 0.15) - expected_kp_STAR = np.array([[ 21, 497], + actual_kp_star, actual_scale = censure_keypoints(img, 7, 'STAR', 0.15) + expected_kp_star = np.array([[ 21, 497], [ 36, 46], [117, 356], [185, 177], @@ -67,7 +66,7 @@ def test_censure_keypoints_moon_image_STAR(): expected_scale = np.array([3, 3, 6, 2, 3, 2, 3, 5, 2, 2]) - assert_array_equal(expected_kp_STAR, actual_kp_STAR) + assert_array_equal(expected_kp_star, actual_kp_star) assert_array_equal(expected_scale, actual_scale) From 229910efe739b2a53841af3faaaa9bb69f83d177 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 12 Aug 2013 16:37:09 +0530 Subject: [PATCH 64/77] Replacing censure_keypoints by keypoints_censure, n_scales by max_scale --- skimage/feature/__init__.py | 4 ++-- skimage/feature/censure.py | 24 ++++++++++++------------ skimage/feature/tests/test_censure.py | 18 +++++++++--------- 3 files changed, 23 insertions(+), 23 deletions(-) diff --git a/skimage/feature/__init__.py b/skimage/feature/__init__.py index 53c21d31..bde81e59 100644 --- a/skimage/feature/__init__.py +++ b/skimage/feature/__init__.py @@ -9,7 +9,7 @@ from .corner_cy import corner_moravec from .template import match_template from ._brief import brief, match_keypoints_brief from .util import pairwise_hamming_distance -from .censure import censure_keypoints +from .censure import keypoints_censure __all__ = ['daisy', 'hog', @@ -28,4 +28,4 @@ __all__ = ['daisy', 'brief', 'pairwise_hamming_distance', 'match_keypoints_brief', - 'censure_keypoints'] + 'keypoints_censure'] diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 583e3d0a..3d11177b 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -19,9 +19,9 @@ STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] -def _get_filtered_image(image, n_scales, mode): +def _get_filtered_image(image, max_scale, mode): - scales = np.zeros((image.shape[0], image.shape[1], n_scales), + scales = np.zeros((image.shape[0], image.shape[1], max_scale), dtype=np.double) if mode == 'dob': @@ -34,7 +34,7 @@ def _get_filtered_image(image, n_scales, mode): integral_img = integral_image(image) - for i in range(n_scales): + for i in range(max_scale): n = i + 1 # Constant multipliers for the outer region and the inner region @@ -54,14 +54,14 @@ def _get_filtered_image(image, n_scales, mode): elif mode == 'octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 - for i in range(n_scales): + for i in range(max_scale): mo, no = OCTAGON_OUTER_SHAPE[i] mi, ni = OCTAGON_INNER_SHAPE[i] scales[:, :, i] = convolve(image, _octagon_filter_kernel(mo, no, mi, ni)) elif mode == 'star': - for i in range(n_scales): + for i in range(max_scale): m = STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] n = STAR_SHAPE[STAR_FILTER_SHAPE[i][1]] scales[:, :, i] = convolve(image, @@ -135,7 +135,7 @@ def _suppress_lines(feature_mask, image, sigma, line_threshold): > line_threshold * (Axx * Ayy - Axy * Axy)] = False -def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, +def keypoints_censure(image, max_scale=7, mode='DoB', non_max_threshold=0.15, line_threshold=10): """ Extracts Censure keypoints along with the corresponding scale using @@ -145,7 +145,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, ---------- image : 2D ndarray Input image. - n_scales : positive integer + max_scale : positive integer Number of scales to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last. mode : ('DoB', 'Octagon', 'STAR') @@ -192,7 +192,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, raise ValueError('Mode must be one of "DoB", "Octagon", "STAR".') # Generating all the scales - filter_response = _get_filtered_image(image, n_scales, mode) + filter_response = _get_filtered_image(image, max_scale, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero @@ -202,7 +202,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, feature_mask = minimas | maximas feature_mask[filter_response < non_max_threshold] = False - for i in range(1, n_scales - 1): + for i in range(1, max_scale - 1): # sigma = (window_size - 1) / 6.0, so the window covers > 99% of the # kernel's distribution # window_size = 7 + 2 * i @@ -210,7 +210,7 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, _suppress_lines(feature_mask[:, :, i], image, (1 + i / 3.0), line_threshold) - rows, cols, scales = np.nonzero(feature_mask[..., 1:n_scales - 1]) + rows, cols, scales = np.nonzero(feature_mask[..., 1:max_scale - 1]) keypoints = np.column_stack([rows, cols]) scales = scales + 2 @@ -220,13 +220,13 @@ def censure_keypoints(image, n_scales=7, mode='DoB', non_max_threshold=0.15, cumulative_mask = np.zeros(keypoints.shape[0], dtype=np.bool) if mode == 'octagon': - for i in range(2, n_scales): + for i in range(2, max_scale): c = (OCTAGON_OUTER_SHAPE[i - 1][0] - 1) // 2 \ + OCTAGON_OUTER_SHAPE[i - 1][1] cumulative_mask |= _mask_border_keypoints(image, keypoints, c) \ & (scales == i) elif mode == 'star': - for i in range(2, n_scales): + for i in range(2, max_scale): c = STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] \ + STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] // 2 cumulative_mask |= _mask_border_keypoints(image, keypoints, c) \ diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 9acdd2e9..5401eed1 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -1,20 +1,20 @@ import numpy as np from numpy.testing import assert_array_equal, assert_raises from skimage.data import moon -from skimage.feature import censure_keypoints +from skimage.feature import keypoints_censure -def test_censure_keypoints_color_image_unsupported_error(): +def test_keypoints_censure_color_image_unsupported_error(): """Censure keypoints can be extracted from gray-scale images only.""" img = np.zeros((20, 20, 3)) - assert_raises(ValueError, censure_keypoints, img) + assert_raises(ValueError, keypoints_censure, img) -def test_censure_keypoints_moon_image_dob(): +def test_keypoints_censure_moon_image_dob(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for DoB filter.""" img = moon() - actual_kp_dob, actual_scale = censure_keypoints(img, 7, 'DoB', 0.15) + actual_kp_dob, actual_scale = keypoints_censure(img, 7, 'DoB', 0.15) expected_kp_dob = np.array([[ 21, 497], [ 36, 46], [119, 350], @@ -30,11 +30,11 @@ def test_censure_keypoints_moon_image_dob(): assert_array_equal(expected_scale, actual_scale) -def test_censure_keypoints_moon_image_Octagon(): +def test_keypoints_censure_moon_image_Octagon(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for Octagon filter.""" img = moon() - actual_kp_octagon, actual_scale = censure_keypoints(img, 7, 'Octagon', + actual_kp_octagon, actual_scale = keypoints_censure(img, 7, 'Octagon', 0.15) expected_kp_octagon = np.array([[ 21, 496], [ 35, 46], @@ -48,11 +48,11 @@ def test_censure_keypoints_moon_image_Octagon(): assert_array_equal(expected_scale, actual_scale) -def test_censure_keypoints_moon_image_STAR(): +def test_keypoints_censure_moon_image_STAR(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for STAR filter.""" img = moon() - actual_kp_star, actual_scale = censure_keypoints(img, 7, 'STAR', 0.15) + actual_kp_star, actual_scale = keypoints_censure(img, 7, 'STAR', 0.15) expected_kp_star = np.array([[ 21, 497], [ 36, 46], [117, 356], From 6d1be9aee40cf2b51205c1be491b927d829cb500 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 12 Aug 2013 16:59:09 +0530 Subject: [PATCH 65/77] Documenting the code --- skimage/feature/censure.py | 10 ++++++++-- skimage/feature/censure_cy.pyx | 4 ++++ 2 files changed, 12 insertions(+), 2 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 3d11177b..18313eaf 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -149,8 +149,14 @@ def keypoints_censure(image, max_scale=7, mode='DoB', non_max_threshold=0.15, Number of scales to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last. mode : ('DoB', 'Octagon', 'STAR') - Type of bilevel filter used to get the scales of input image. Possible - values are 'DoB', 'Octagon' and 'STAR'. + Type of bilevel filter used to get the scales of the input image. + Possible values are 'DoB', 'Octagon' and 'STAR'. The three modes + represent the shape of the bilevel filters i.e. box(square), octagon + and star respectively. For instance, a bilevel octagon filter consists + of a smaller inner octagon and a larger outer octagon with the filter + weights being uniformly negative in both the inner octagon while + uniformly positive in the difference region. Use STAR and Octagon for + better features and DoB for better performance. non_max_threshold : float Threshold value used to suppress maximas and minimas with a weak magnitude response obtained after Non-Maximal Suppression. diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index cfd1260f..a9980071 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -8,6 +8,10 @@ def _censure_dob_loop(Py_ssize_t n, double[:, ::1] integral_img, double[:, ::1] filtered_image, double inner_weight, double outer_weight): + # This function calculates the value in the DoB filtered image using + # integral images. If r = right. l = left, u = up, d = down, the sum of + # pixel values in the rectangle formed by (u, l), (u, r), (d, r), (d, l) + # is calculated as I(d, r) + I(u - 1, l - 1) - I(u - 1, r) - I(d, l - 1). cdef Py_ssize_t i, j cdef double inner, outer From c3b6ce2c8270abd5801afcd1633f7600f3b2f071 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Mon, 12 Aug 2013 23:51:19 +0530 Subject: [PATCH 66/77] Including min_scale as another input parameter --- skimage/feature/censure.py | 55 +++++++++++++++------------ skimage/feature/censure_cy.pyx | 42 +++++++++++++++++++- skimage/feature/tests/test_censure.py | 10 ++--- 3 files changed, 76 insertions(+), 31 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 18313eaf..996f7e16 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -19,10 +19,10 @@ STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] -def _get_filtered_image(image, max_scale, mode): +def _get_filtered_image(image, min_scale, max_scale, mode): - scales = np.zeros((image.shape[0], image.shape[1], max_scale), - dtype=np.double) + scales = np.zeros((image.shape[0], image.shape[1], + max_scale - min_scale +1), dtype=np.double) if mode == 'dob': @@ -34,8 +34,8 @@ def _get_filtered_image(image, max_scale, mode): integral_img = integral_image(image) - for i in range(max_scale): - n = i + 1 + for i in range(max_scale - min_scale + 1): + n = min_scale + i # Constant multipliers for the outer region and the inner region # of the bilevel filters with the constraint of keeping the @@ -54,16 +54,16 @@ def _get_filtered_image(image, max_scale, mode): elif mode == 'octagon': # TODO : Decide the shapes of Octagon filters for scales > 7 - for i in range(max_scale): - mo, no = OCTAGON_OUTER_SHAPE[i] - mi, ni = OCTAGON_INNER_SHAPE[i] + for i in range(max_scale - min_scale + 1): + mo, no = OCTAGON_OUTER_SHAPE[min_scale + i - 1] + mi, ni = OCTAGON_INNER_SHAPE[min_scale + i - 1] scales[:, :, i] = convolve(image, _octagon_filter_kernel(mo, no, mi, ni)) elif mode == 'star': - for i in range(max_scale): - m = STAR_SHAPE[STAR_FILTER_SHAPE[i][0]] - n = STAR_SHAPE[STAR_FILTER_SHAPE[i][1]] + for i in range(max_scale - min_scale + 1): + m = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][0]] + n = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][1]] scales[:, :, i] = convolve(image, _star_filter_kernel(m, n)) @@ -135,8 +135,8 @@ def _suppress_lines(feature_mask, image, sigma, line_threshold): > line_threshold * (Axx * Ayy - Axy * Axy)] = False -def keypoints_censure(image, max_scale=7, mode='DoB', non_max_threshold=0.15, - line_threshold=10): +def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', + non_max_threshold=0.15, line_threshold=10): """ Extracts Censure keypoints along with the corresponding scale using either Difference of Boxes, Octagon or STAR bilevel filter. @@ -145,9 +145,12 @@ def keypoints_censure(image, max_scale=7, mode='DoB', non_max_threshold=0.15, ---------- image : 2D ndarray Input image. + min_scale : positive integer + Minimum scale to extract keypoints from. max_scale : positive integer - Number of scales to extract keypoints from. The keypoints will be - extracted from all the scales except the first and the last. + Maximum scale to extract keypoints from. The keypoints will be + extracted from all the scales except the first and the last i.e. + from the scales in the range [min_scale + 1, max_scale - 1]. mode : ('DoB', 'Octagon', 'STAR') Type of bilevel filter used to get the scales of the input image. Possible values are 'DoB', 'Octagon' and 'STAR'. The three modes @@ -197,8 +200,12 @@ def keypoints_censure(image, max_scale=7, mode='DoB', non_max_threshold=0.15, if mode not in ('dob', 'octagon', 'star'): raise ValueError('Mode must be one of "DoB", "Octagon", "STAR".') + if max_scale - min_scale < 2: + raise ValueError('The number of scales should be greater than or' + 'equal to 3.') + # Generating all the scales - filter_response = _get_filtered_image(image, max_scale, mode) + filter_response = _get_filtered_image(image, min_scale, max_scale, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero @@ -208,17 +215,17 @@ def keypoints_censure(image, max_scale=7, mode='DoB', non_max_threshold=0.15, feature_mask = minimas | maximas feature_mask[filter_response < non_max_threshold] = False - for i in range(1, max_scale - 1): + for i in range(1, max_scale - min_scale): # sigma = (window_size - 1) / 6.0, so the window covers > 99% of the # kernel's distribution - # window_size = 7 + 2 * i - # Hence sigma = 1 + i / 3.0 + # window_size = 7 + 2 * (min_scale - 1 + i) + # Hence sigma = 1 + (min_scale - 1 + i)/ 3.0 _suppress_lines(feature_mask[:, :, i], image, - (1 + i / 3.0), line_threshold) + (1 + (min_scale + i - 1) / 3.0), line_threshold) - rows, cols, scales = np.nonzero(feature_mask[..., 1:max_scale - 1]) + rows, cols, scales = np.nonzero(feature_mask[..., 1:max_scale - min_scale]) keypoints = np.column_stack([rows, cols]) - scales = scales + 2 + scales = scales + min_scale + 1 if mode == 'dob': return keypoints, scales @@ -226,13 +233,13 @@ def keypoints_censure(image, max_scale=7, mode='DoB', non_max_threshold=0.15, cumulative_mask = np.zeros(keypoints.shape[0], dtype=np.bool) if mode == 'octagon': - for i in range(2, max_scale): + for i in range(min_scale + 1, max_scale): c = (OCTAGON_OUTER_SHAPE[i - 1][0] - 1) // 2 \ + OCTAGON_OUTER_SHAPE[i - 1][1] cumulative_mask |= _mask_border_keypoints(image, keypoints, c) \ & (scales == i) elif mode == 'star': - for i in range(2, max_scale): + for i in range(min_scale + 1, max_scale): c = STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] \ + STAR_SHAPE[STAR_FILTER_SHAPE[i - 1][0]] // 2 cumulative_mask |= _mask_border_keypoints(image, keypoints, c) \ diff --git a/skimage/feature/censure_cy.pyx b/skimage/feature/censure_cy.pyx index a9980071..1c352fde 100644 --- a/skimage/feature/censure_cy.pyx +++ b/skimage/feature/censure_cy.pyx @@ -18,8 +18,46 @@ def _censure_dob_loop(Py_ssize_t n, cdef Py_ssize_t n2 = 2 * n cdef double total_weight = inner_weight + outer_weight - for i in range(n2, integral_img.shape[0] - n2): - for j in range(n2, integral_img.shape[1] - n2): + # top-left pixel + inner = (integral_img[n2 + n, n2 + n] + + integral_img[n2 - n - 1, n2 - n - 1] + - integral_img[n2 + n, n2 - n - 1] + - integral_img[n2 - n - 1, n2 + n]) + + outer = integral_img[2 * n2, 2 * n2] + + filtered_image[n2, n2] = (outer_weight * outer + - total_weight * inner) + + # left column + for i in range(n2 + 1, integral_img.shape[0] - n2): + inner = (integral_img[i + n, n2 + n] + + integral_img[i - n - 1, n2 - n - 1] + - integral_img[i + n, n2 - n - 1] + - integral_img[i - n - 1, n2 + n]) + + outer = (integral_img[i + n2, 2 * n2] + - integral_img[i - n2 - 1, 2 * n2]) + + filtered_image[i, n2] = (outer_weight * outer + - total_weight * inner) + + # top row + for j in range(n2 + 1, integral_img.shape[1] - n2): + inner = (integral_img[n2 + n, j + n] + + integral_img[n2 - n - 1, j - n - 1] + - integral_img[n2 + n, j - n - 1] + - integral_img[n2 - n - 1, j + n]) + + outer = (integral_img[2 * n2, j + n2] + - integral_img[2 * n2, j - n2 - 1]) + + filtered_image[n2, j] = (outer_weight * outer + - total_weight * inner) + + # remaining block + for i in range(n2 + 1, integral_img.shape[0] - n2): + for j in range(n2 + 1, integral_img.shape[1] - n2): inner = (integral_img[i + n, j + n] + integral_img[i - n - 1, j - n - 1] - integral_img[i + n, j - n - 1] diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 5401eed1..26555697 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -14,7 +14,7 @@ def test_keypoints_censure_moon_image_dob(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for DoB filter.""" img = moon() - actual_kp_dob, actual_scale = keypoints_censure(img, 7, 'DoB', 0.15) + actual_kp_dob, actual_scale = keypoints_censure(img, 1, 7, 'DoB', 0.15) expected_kp_dob = np.array([[ 21, 497], [ 36, 46], [119, 350], @@ -30,11 +30,11 @@ def test_keypoints_censure_moon_image_dob(): assert_array_equal(expected_scale, actual_scale) -def test_keypoints_censure_moon_image_Octagon(): +def test_keypoints_censure_moon_image_octagon(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for Octagon filter.""" img = moon() - actual_kp_octagon, actual_scale = keypoints_censure(img, 7, 'Octagon', + actual_kp_octagon, actual_scale = keypoints_censure(img, 1, 7, 'Octagon', 0.15) expected_kp_octagon = np.array([[ 21, 496], [ 35, 46], @@ -48,11 +48,11 @@ def test_keypoints_censure_moon_image_Octagon(): assert_array_equal(expected_scale, actual_scale) -def test_keypoints_censure_moon_image_STAR(): +def test_keypoints_censure_moon_image_star(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for STAR filter.""" img = moon() - actual_kp_star, actual_scale = keypoints_censure(img, 7, 'STAR', 0.15) + actual_kp_star, actual_scale = keypoints_censure(img, 1, 7, 'STAR', 0.15) expected_kp_star = np.array([[ 21, 497], [ 36, 46], [117, 356], From 6c27d278d0f6ae107ea4a24b6ca1664cfb5c3fb5 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 13 Aug 2013 00:21:25 +0530 Subject: [PATCH 67/77] Adding remaining tests --- skimage/feature/censure.py | 7 +++---- skimage/feature/tests/test_censure.py | 14 ++++++++++++++ 2 files changed, 17 insertions(+), 4 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 996f7e16..d52583c2 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -192,11 +192,7 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") - image = img_as_float(image) - image = np.ascontiguousarray(image) - mode = mode.lower() - if mode not in ('dob', 'octagon', 'star'): raise ValueError('Mode must be one of "DoB", "Octagon", "STAR".') @@ -204,6 +200,9 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', raise ValueError('The number of scales should be greater than or' 'equal to 3.') + image = img_as_float(image) + image = np.ascontiguousarray(image) + # Generating all the scales filter_response = _get_filtered_image(image, min_scale, max_scale, mode) diff --git a/skimage/feature/tests/test_censure.py b/skimage/feature/tests/test_censure.py index 26555697..4cd2ad68 100644 --- a/skimage/feature/tests/test_censure.py +++ b/skimage/feature/tests/test_censure.py @@ -10,6 +10,20 @@ def test_keypoints_censure_color_image_unsupported_error(): assert_raises(ValueError, keypoints_censure, img) +def test_keypoints_censure_mode_validity_error(): + """Mode argument in keypoints_censure can be either DoB, Octagon or + STAR.""" + img = np.zeros((20, 20)) + assert_raises(ValueError, keypoints_censure, img, mode='dummy') + + +def test_keypoints_censure_scale_range_error(): + """Difference between the the max_scale and min_scale parameters in + keypoints_censure should be greater than or equal to two.""" + img = np.zeros((20, 20)) + assert_raises(ValueError, keypoints_censure, img, min_scale=1, max_scale=2) + + def test_keypoints_censure_moon_image_dob(): """Verify the actual Censure keypoints and their corresponding scale with the expected values for DoB filter.""" From f2632f26bb4ec2586b4150250a54c471b818a667 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Tue, 13 Aug 2013 20:11:37 +0530 Subject: [PATCH 68/77] Adding gallery example for plotting censure keypoints --- doc/examples/plot_censure_keypoints.py | 46 ++++++++++++++++++++++++++ 1 file changed, 46 insertions(+) create mode 100644 doc/examples/plot_censure_keypoints.py diff --git a/doc/examples/plot_censure_keypoints.py b/doc/examples/plot_censure_keypoints.py new file mode 100644 index 00000000..c2119b46 --- /dev/null +++ b/doc/examples/plot_censure_keypoints.py @@ -0,0 +1,46 @@ +""" +=========================== +Censure Keypoints Detection +=========================== + +In this example, we detect and plot the Censure Keypoints at various scales +using Difference of Boxes, Octagon and Star shaped bi-level filters. + +""" +from skimage.feature import keypoints_censure +from skimage.data import lena +from skimage.color import rgb2gray +import matplotlib.pyplot as plt + +# Initializing the parameters for Censure keypoints +img = lena() +gray_img = rgb2gray(img) +min_scale = 1 +max_scale = 7 +nms_threshold = 0.15 +rpc_threshold = 10 + + +# Plotting features for the following modes +for mode in ['dob', 'octagon', 'star']: + + kp_censure, scale = keypoints_censure(gray_img, min_scale, max_scale, + mode, nms_threshold,rpc_threshold) + f, axarr = plt.subplots((max_scale - min_scale + 1) // 3, 3) + + # Plotting Censure features at all the scales + for i in range(max_scale - min_scale - 1): + keypoints = kp_censure[scale == i + min_scale + 1] + num = len(keypoints) + x = keypoints[:, 1] + y = keypoints[:, 0] + s = 5 * 2**(i + min_scale + 1) + axarr[i // 3, i - (i // 3) * 3].imshow(img) + axarr[i // 3, i - (i // 3) * 3].scatter(x, y, s, facecolors='none', + edgecolors='g') + axarr[i // 3, i - (i // 3) * 3].set_title(' %s %s-Censure features at ' + 'scale %d' % (num, mode, i + + min_scale + 1)) + + plt.suptitle('NMS threshold = %f, RPC threshold = %d' % (nms_threshold, rpc_threshold)) +plt.show() From b61b3303559da8c8ca4018c09228f1c32b933d2b Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 14 Aug 2013 14:25:15 +0530 Subject: [PATCH 69/77] Adding a short description of the algorithm --- skimage/feature/censure.py | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index d52583c2..272218f6 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -187,7 +187,17 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', http://www.jamris.org/01_2013/saveas.php?QUEST=JAMRIS_No01_2013_P_11-20.pdf """ - + # (1) First we generate the required scales on the input grayscale image + # using a bilevel filter and stack them up in `filter_response`. + # (2) We then perform Non-Maximal suppression in 3 x 3 x 3 window on the + # filter_response to suppress points that are neither minima or maxima in + # 3 x 3 x 3 neighbourhood. We obtain a boolean ndarray `feature_mask` + # containing all the minimas and maximas in `filter_response` as True. + # (3) Then we suppress all the points in the `feature_mask` for which the + # corresponding point in the image at a particular scale has the ratio of + # principal curvatures greater than `line_threshold`. + # (4) Finally, we remove the border keypoints and return the keypoints + # along with its corresponding scale. image = np.squeeze(image) if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") From 3ee183fdf8714ca5d9f7aabc3e0581b2c13c0833 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Wed, 14 Aug 2013 14:36:56 +0530 Subject: [PATCH 70/77] PEP8 corrections and stylistic changes --- doc/examples/plot_censure_keypoints.py | 10 +++++----- skimage/feature/censure.py | 9 +++++---- 2 files changed, 10 insertions(+), 9 deletions(-) diff --git a/doc/examples/plot_censure_keypoints.py b/doc/examples/plot_censure_keypoints.py index c2119b46..7f08d409 100644 --- a/doc/examples/plot_censure_keypoints.py +++ b/doc/examples/plot_censure_keypoints.py @@ -20,17 +20,16 @@ max_scale = 7 nms_threshold = 0.15 rpc_threshold = 10 - -# Plotting features for the following modes +# Detecting Censure keypoints for the following filters for mode in ['dob', 'octagon', 'star']: kp_censure, scale = keypoints_censure(gray_img, min_scale, max_scale, - mode, nms_threshold,rpc_threshold) + mode, nms_threshold, rpc_threshold) f, axarr = plt.subplots((max_scale - min_scale + 1) // 3, 3) # Plotting Censure features at all the scales for i in range(max_scale - min_scale - 1): - keypoints = kp_censure[scale == i + min_scale + 1] + keypoints = kp_censure[scale == i + min_scale + 1] num = len(keypoints) x = keypoints[:, 1] y = keypoints[:, 0] @@ -42,5 +41,6 @@ for mode in ['dob', 'octagon', 'star']: 'scale %d' % (num, mode, i + min_scale + 1)) - plt.suptitle('NMS threshold = %f, RPC threshold = %d' % (nms_threshold, rpc_threshold)) + plt.suptitle('NMS threshold = %f, RPC threshold = %d' + % (nms_threshold, rpc_threshold)) plt.show() diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 272218f6..d90b36cf 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -16,13 +16,13 @@ OCTAGON_INNER_SHAPE = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] STAR_SHAPE = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), - (9, 6),(11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] + (9, 6), (11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] def _get_filtered_image(image, min_scale, max_scale, mode): scales = np.zeros((image.shape[0], image.shape[1], - max_scale - min_scale +1), dtype=np.double) + max_scale - min_scale + 1), dtype=np.double) if mode == 'dob': @@ -75,7 +75,7 @@ def _oct(m, n): f = np.zeros((m + 2*n, m + 2*n)) f[0, n] = 1 f[n, 0] = 1 - f[0, m + n -1] = 1 + f[0, m + n - 1] = 1 f[m + n - 1, 0] = 1 f[-1, n] = 1 f[n, -1] = 1 @@ -198,6 +198,7 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', # principal curvatures greater than `line_threshold`. # (4) Finally, we remove the border keypoints and return the keypoints # along with its corresponding scale. + image = np.squeeze(image) if image.ndim != 2: raise ValueError("Only 2-D gray-scale images supported.") @@ -230,7 +231,7 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', # window_size = 7 + 2 * (min_scale - 1 + i) # Hence sigma = 1 + (min_scale - 1 + i)/ 3.0 _suppress_lines(feature_mask[:, :, i], image, - (1 + (min_scale + i - 1) / 3.0), line_threshold) + (1 + (min_scale + i - 1) / 3.0), line_threshold) rows, cols, scales = np.nonzero(feature_mask[..., 1:max_scale - min_scale]) keypoints = np.column_stack([rows, cols]) From 7c4152b62bae274228272b3966cb400681cc0d4a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Johannes=20Sch=C3=B6nberger?= Date: Thu, 15 Aug 2013 08:36:09 +0200 Subject: [PATCH 71/77] Fix censure example and fix some minor issues --- doc/examples/plot_censure_keypoints.py | 63 +++++++++++++++----------- skimage/feature/censure.py | 58 ++++++++++++------------ 2 files changed, 65 insertions(+), 56 deletions(-) diff --git a/doc/examples/plot_censure_keypoints.py b/doc/examples/plot_censure_keypoints.py index 7f08d409..95fa0e58 100644 --- a/doc/examples/plot_censure_keypoints.py +++ b/doc/examples/plot_censure_keypoints.py @@ -1,12 +1,14 @@ """ -=========================== -Censure Keypoints Detection -=========================== +========================= +CenSurE Feature Detection +========================= -In this example, we detect and plot the Censure Keypoints at various scales -using Difference of Boxes, Octagon and Star shaped bi-level filters. +In this example we detect and plot the CenSurE (Center Surround Extrema) +features at various scales using Difference of Boxes, Octagon and Star shaped +bi-level filters. """ + from skimage.feature import keypoints_censure from skimage.data import lena from skimage.color import rgb2gray @@ -15,32 +17,39 @@ import matplotlib.pyplot as plt # Initializing the parameters for Censure keypoints img = lena() gray_img = rgb2gray(img) -min_scale = 1 -max_scale = 7 -nms_threshold = 0.15 -rpc_threshold = 10 +min_scale = 2 +max_scale = 6 +non_max_threshold = 0.15 +line_threshold = 10 + + +f, ax = plt.subplots(nrows=(max_scale - min_scale - 1), ncols=3, + figsize=(6, 6)) +plt.subplots_adjust(wspace=0.02, hspace=0.02, top=0.94, + bottom=0.02, left=0.06, right=0.98) # Detecting Censure keypoints for the following filters -for mode in ['dob', 'octagon', 'star']: +for col, mode in enumerate(['dob', 'octagon', 'star']): - kp_censure, scale = keypoints_censure(gray_img, min_scale, max_scale, - mode, nms_threshold, rpc_threshold) - f, axarr = plt.subplots((max_scale - min_scale + 1) // 3, 3) + ax[0, col].set_title(mode.upper(), fontsize=12) + + keypoints, scales = keypoints_censure(gray_img, min_scale, max_scale, + mode, non_max_threshold, + line_threshold) # Plotting Censure features at all the scales - for i in range(max_scale - min_scale - 1): - keypoints = kp_censure[scale == i + min_scale + 1] - num = len(keypoints) - x = keypoints[:, 1] - y = keypoints[:, 0] - s = 5 * 2**(i + min_scale + 1) - axarr[i // 3, i - (i // 3) * 3].imshow(img) - axarr[i // 3, i - (i // 3) * 3].scatter(x, y, s, facecolors='none', - edgecolors='g') - axarr[i // 3, i - (i // 3) * 3].set_title(' %s %s-Censure features at ' - 'scale %d' % (num, mode, i + - min_scale + 1)) + for row, scale in enumerate(range(min_scale + 1, max_scale)): + keypoints_i = keypoints[scales == scale] + num = len(keypoints_i) + x = keypoints_i[:, 1] + y = keypoints_i[:, 0] + s = 0.5 * 2 ** (scale + min_scale + 1) + ax[row, col].imshow(img) + ax[row, col].scatter(x, y, s, facecolors='none', edgecolors='b') + ax[row, col].set_xticks([]) + ax[row, col].set_yticks([]) + ax[row, col].axis((0, img.shape[1], img.shape[0], 0)) + if col == 0: + ax[row, col].set_ylabel('Scale %d' % scale, fontsize=12) - plt.suptitle('NMS threshold = %f, RPC threshold = %d' - % (nms_threshold, rpc_threshold)) plt.show() diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index d90b36cf..2dda9306 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -19,18 +19,17 @@ STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6), (11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] -def _get_filtered_image(image, min_scale, max_scale, mode): +def _filter_image(image, min_scale, max_scale, mode): - scales = np.zeros((image.shape[0], image.shape[1], - max_scale - min_scale + 1), dtype=np.double) + response = np.zeros((image.shape[0], image.shape[1], + max_scale - min_scale + 1), dtype=np.double) if mode == 'dob': - # make scales[:, :, i] contiguous memory block - item_size = scales.itemsize - scales.strides = (item_size * scales.shape[0], - item_size, - item_size * scales.shape[0] * scales.shape[1]) + # make response[:, :, i] contiguous memory block + item_size = response.itemsize + response.strides = (item_size * response.shape[0], item_size, + item_size * response.shape[0] * response.shape[1]) integral_img = integral_image(image) @@ -38,12 +37,12 @@ def _get_filtered_image(image, min_scale, max_scale, mode): n = min_scale + i # Constant multipliers for the outer region and the inner region - # of the bilevel filters with the constraint of keeping the + # of the bi-level filters with the constraint of keeping the # DC bias 0. inner_weight = (1.0 / (2 * n + 1)**2) outer_weight = (1.0 / (12 * n**2 + 4 * n)) - _censure_dob_loop(n, integral_img, scales[:, :, i], + _censure_dob_loop(n, integral_img, response[:, :, i], inner_weight, outer_weight) # NOTE : For the Octagon shaped filter, we implemented and evaluated the @@ -57,17 +56,17 @@ def _get_filtered_image(image, min_scale, max_scale, mode): for i in range(max_scale - min_scale + 1): mo, no = OCTAGON_OUTER_SHAPE[min_scale + i - 1] mi, ni = OCTAGON_INNER_SHAPE[min_scale + i - 1] - scales[:, :, i] = convolve(image, - _octagon_filter_kernel(mo, no, mi, ni)) + response[:, :, i] = convolve(image, + _octagon_filter_kernel(mo, no, mi, ni)) + elif mode == 'star': for i in range(max_scale - min_scale + 1): m = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][0]] n = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][1]] - scales[:, :, i] = convolve(image, - _star_filter_kernel(m, n)) + response[:, :, i] = convolve(image, _star_filter_kernel(m, n)) - return scales + return response # TODO : Import from selem after getting #669 merged. @@ -138,24 +137,24 @@ def _suppress_lines(feature_mask, image, sigma, line_threshold): def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', non_max_threshold=0.15, line_threshold=10): """ - Extracts Censure keypoints along with the corresponding scale using - either Difference of Boxes, Octagon or STAR bilevel filter. + Extracts CenSurE keypoints along with the corresponding scale using + either Difference of Boxes, Octagon or STAR bi-level filter. Parameters ---------- image : 2D ndarray Input image. - min_scale : positive integer + min_scale : int Minimum scale to extract keypoints from. - max_scale : positive integer + max_scale : int Maximum scale to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last i.e. from the scales in the range [min_scale + 1, max_scale - 1]. - mode : ('DoB', 'Octagon', 'STAR') - Type of bilevel filter used to get the scales of the input image. + mode : {'DoB', 'Octagon', 'STAR'} + Type of bi-level filter used to get the scales of the input image. Possible values are 'DoB', 'Octagon' and 'STAR'. The three modes - represent the shape of the bilevel filters i.e. box(square), octagon - and star respectively. For instance, a bilevel octagon filter consists + represent the shape of the bi-level filters i.e. box(square), octagon + and star respectively. For instance, a bi-level octagon filter consists of a smaller inner octagon and a larger outer octagon with the filter weights being uniformly negative in both the inner octagon while uniformly positive in the difference region. Use STAR and Octagon for @@ -170,7 +169,7 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', Returns ------- keypoints : (N, 2) array - Location of the extracted keypoints in the (row, col) format. + Location of the extracted keypoints in the ``(row, col)`` format. scales : (N, 1) array The corresponding scale of the N extracted keypoints. @@ -187,8 +186,9 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', http://www.jamris.org/01_2013/saveas.php?QUEST=JAMRIS_No01_2013_P_11-20.pdf """ + # (1) First we generate the required scales on the input grayscale image - # using a bilevel filter and stack them up in `filter_response`. + # using a bi-level filter and stack them up in `filter_response`. # (2) We then perform Non-Maximal suppression in 3 x 3 x 3 window on the # filter_response to suppress points that are neither minima or maxima in # 3 x 3 x 3 neighbourhood. We obtain a boolean ndarray `feature_mask` @@ -207,15 +207,15 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', if mode not in ('dob', 'octagon', 'star'): raise ValueError('Mode must be one of "DoB", "Octagon", "STAR".') - if max_scale - min_scale < 2: - raise ValueError('The number of scales should be greater than or' - 'equal to 3.') + if min_scale < 1 or max_scale < 1 or max_scale - min_scale < 2: + raise ValueError('The scales must be >= 1 and the number of scales ' + 'should be >= 3.') image = img_as_float(image) image = np.ascontiguousarray(image) # Generating all the scales - filter_response = _get_filtered_image(image, min_scale, max_scale, mode) + filter_response = _filter_image(image, min_scale, max_scale, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero From 8d6d5f7fa783ad4cd51daef42b6129f067d3bf0a Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 15 Aug 2013 14:47:10 +0530 Subject: [PATCH 72/77] Adding Octagon shapes for higher scales --- skimage/feature/censure.py | 63 +++++++++++++++++++------------------- 1 file changed, 32 insertions(+), 31 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 2dda9306..199cfbaa 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -11,25 +11,27 @@ from skimage.feature.censure_cy import _censure_dob_loop OCTAGON_OUTER_SHAPE = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), - (15, 10)] -OCTAGON_INNER_SHAPE = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5)] + (15, 10), (15, 11), (15, 12), (17, 13), (17, 14)] +OCTAGON_INNER_SHAPE = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5), + (7, 5), (7, 6), (9, 6), (9, 7)] STAR_SHAPE = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6), (11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] -def _filter_image(image, min_scale, max_scale, mode): +def _get_filtered_image(image, min_scale, max_scale, mode): - response = np.zeros((image.shape[0], image.shape[1], - max_scale - min_scale + 1), dtype=np.double) + scales = np.zeros((image.shape[0], image.shape[1], + max_scale - min_scale + 1), dtype=np.double) if mode == 'dob': - # make response[:, :, i] contiguous memory block - item_size = response.itemsize - response.strides = (item_size * response.shape[0], item_size, - item_size * response.shape[0] * response.shape[1]) + # make scales[:, :, i] contiguous memory block + item_size = scales.itemsize + scales.strides = (item_size * scales.shape[0], + item_size, + item_size * scales.shape[0] * scales.shape[1]) integral_img = integral_image(image) @@ -37,12 +39,12 @@ def _filter_image(image, min_scale, max_scale, mode): n = min_scale + i # Constant multipliers for the outer region and the inner region - # of the bi-level filters with the constraint of keeping the + # of the bilevel filters with the constraint of keeping the # DC bias 0. inner_weight = (1.0 / (2 * n + 1)**2) outer_weight = (1.0 / (12 * n**2 + 4 * n)) - _censure_dob_loop(n, integral_img, response[:, :, i], + _censure_dob_loop(n, integral_img, scales[:, :, i], inner_weight, outer_weight) # NOTE : For the Octagon shaped filter, we implemented and evaluated the @@ -56,17 +58,17 @@ def _filter_image(image, min_scale, max_scale, mode): for i in range(max_scale - min_scale + 1): mo, no = OCTAGON_OUTER_SHAPE[min_scale + i - 1] mi, ni = OCTAGON_INNER_SHAPE[min_scale + i - 1] - response[:, :, i] = convolve(image, - _octagon_filter_kernel(mo, no, mi, ni)) - + scales[:, :, i] = convolve(image, + _octagon_filter_kernel(mo, no, mi, ni)) elif mode == 'star': for i in range(max_scale - min_scale + 1): m = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][0]] n = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][1]] - response[:, :, i] = convolve(image, _star_filter_kernel(m, n)) + scales[:, :, i] = convolve(image, + _star_filter_kernel(m, n)) - return response + return scales # TODO : Import from selem after getting #669 merged. @@ -137,24 +139,24 @@ def _suppress_lines(feature_mask, image, sigma, line_threshold): def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', non_max_threshold=0.15, line_threshold=10): """ - Extracts CenSurE keypoints along with the corresponding scale using - either Difference of Boxes, Octagon or STAR bi-level filter. + Extracts Censure keypoints along with the corresponding scale using + either Difference of Boxes, Octagon or STAR bilevel filter. Parameters ---------- image : 2D ndarray Input image. - min_scale : int + min_scale : positive integer Minimum scale to extract keypoints from. - max_scale : int + max_scale : positive integer Maximum scale to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last i.e. from the scales in the range [min_scale + 1, max_scale - 1]. - mode : {'DoB', 'Octagon', 'STAR'} - Type of bi-level filter used to get the scales of the input image. + mode : ('DoB', 'Octagon', 'STAR') + Type of bilevel filter used to get the scales of the input image. Possible values are 'DoB', 'Octagon' and 'STAR'. The three modes - represent the shape of the bi-level filters i.e. box(square), octagon - and star respectively. For instance, a bi-level octagon filter consists + represent the shape of the bilevel filters i.e. box(square), octagon + and star respectively. For instance, a bilevel octagon filter consists of a smaller inner octagon and a larger outer octagon with the filter weights being uniformly negative in both the inner octagon while uniformly positive in the difference region. Use STAR and Octagon for @@ -169,7 +171,7 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', Returns ------- keypoints : (N, 2) array - Location of the extracted keypoints in the ``(row, col)`` format. + Location of the extracted keypoints in the (row, col) format. scales : (N, 1) array The corresponding scale of the N extracted keypoints. @@ -186,9 +188,8 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', http://www.jamris.org/01_2013/saveas.php?QUEST=JAMRIS_No01_2013_P_11-20.pdf """ - # (1) First we generate the required scales on the input grayscale image - # using a bi-level filter and stack them up in `filter_response`. + # using a bilevel filter and stack them up in `filter_response`. # (2) We then perform Non-Maximal suppression in 3 x 3 x 3 window on the # filter_response to suppress points that are neither minima or maxima in # 3 x 3 x 3 neighbourhood. We obtain a boolean ndarray `feature_mask` @@ -207,15 +208,15 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', if mode not in ('dob', 'octagon', 'star'): raise ValueError('Mode must be one of "DoB", "Octagon", "STAR".') - if min_scale < 1 or max_scale < 1 or max_scale - min_scale < 2: - raise ValueError('The scales must be >= 1 and the number of scales ' - 'should be >= 3.') + if max_scale - min_scale < 2: + raise ValueError('The number of scales should be greater than or' + 'equal to 3.') image = img_as_float(image) image = np.ascontiguousarray(image) # Generating all the scales - filter_response = _filter_image(image, min_scale, max_scale, mode) + filter_response = _get_filtered_image(image, min_scale, max_scale, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero From cd5316c36de726fa8b7b9a2a46153fe94f6ea8a9 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 15 Aug 2013 15:30:08 +0530 Subject: [PATCH 73/77] Importing star and octagon from morphplogy.selem --- skimage/feature/censure.py | 43 +++++--------------------------------- 1 file changed, 5 insertions(+), 38 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index 199cfbaa..a15195a2 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -4,7 +4,7 @@ from scipy.ndimage.filters import maximum_filter, minimum_filter, convolve from skimage.transform import integral_image from skimage.feature.corner import _compute_auto_correlation from skimage.util import img_as_float -from skimage.morphology import convex_hull_image +from skimage.morphology import convex_hull_image, octagon, star from skimage.feature.util import _mask_border_keypoints from skimage.feature.censure_cy import _censure_dob_loop @@ -71,58 +71,25 @@ def _get_filtered_image(image, min_scale, max_scale, mode): return scales -# TODO : Import from selem after getting #669 merged. -def _oct(m, n): - f = np.zeros((m + 2*n, m + 2*n)) - f[0, n] = 1 - f[n, 0] = 1 - f[0, m + n - 1] = 1 - f[m + n - 1, 0] = 1 - f[-1, n] = 1 - f[n, -1] = 1 - f[-1, m + n - 1] = 1 - f[m + n - 1, -1] = 1 - return convex_hull_image(f).astype(int) - - def _octagon_filter_kernel(mo, no, mi, ni): outer = (mo + 2 * no)**2 - 2 * no * (no + 1) inner = (mi + 2 * ni)**2 - 2 * ni * (ni + 1) outer_weight = 1.0 / (outer - inner) inner_weight = 1.0 / inner c = ((mo + 2 * no) - (mi + 2 * ni)) // 2 - outer_oct = _oct(mo, no) + outer_oct = octagon(mo, no) inner_oct = np.zeros((mo + 2 * no, mo + 2 * no)) - inner_oct[c: -c, c: -c] = _oct(mi, ni) + inner_oct[c: -c, c: -c] = octagon(mi, ni) bfilter = (outer_weight * outer_oct - (outer_weight + inner_weight) * inner_oct) return bfilter -def _star(a): - if a == 1: - bfilter = np.zeros((3, 3)) - bfilter[:] = 1 - return bfilter - m = 2 * a + 1 - n = a // 2 - selem_square = np.zeros((m + 2 * n, m + 2 * n), dtype=np.uint8) - selem_square[n: m + n, n: m + n] = 1 - selem_triangle = np.zeros((m + 2 * n, m + 2 * n), dtype=np.uint8) - selem_triangle[(m + 2 * n - 1) // 2, 0] = 1 - selem_triangle[(m + 1) // 2, n - 1] = 1 - selem_triangle[(m + 4 * n - 3) // 2, n - 1] = 1 - selem_triangle = convex_hull_image(selem_triangle).astype(int) - selem_triangle += (selem_triangle[:, ::-1] + selem_triangle.T + - selem_triangle.T[::-1, :]) - return selem_square + selem_triangle - - def _star_filter_kernel(m, n): c = m + m // 2 - n - n // 2 - outer_star = _star(m) + outer_star = star(m) inner_star = np.zeros_like(outer_star) - inner_star[c: -c, c: -c] = _star(n) + inner_star[c: -c, c: -c] = star(n) outer_weight = 1.0 / (np.sum(outer_star - inner_star)) inner_weight = 1.0 / np.sum(inner_star) bfilter = (outer_weight * outer_star - From ffbeeaee27207b924f5e966e45cc202d4b9a2d70 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 15 Aug 2013 16:37:07 +0530 Subject: [PATCH 74/77] Correcting selem.star for Python 3 --- skimage/morphology/selem.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/skimage/morphology/selem.py b/skimage/morphology/selem.py index f85b3a95..4df52e94 100644 --- a/skimage/morphology/selem.py +++ b/skimage/morphology/selem.py @@ -280,10 +280,10 @@ def star(a, dtype=np.uint8): bfilter[:] = 1 return bfilter m = 2 * a + 1 - n = a / 2 + n = a // 2 selem_square = np.zeros((m + 2 * n, m + 2 * n)) selem_square[n: m + n, n: m + n] = 1 - c = (m + 2 * n - 1) / 2 + c = (m + 2 * n - 1) // 2 selem_rotated = np.zeros((m + 2 * n, m + 2 * n)) selem_rotated[0, c] = selem_rotated[-1, c] = selem_rotated[c, 0] = selem_rotated[c, -1] = 1 selem_rotated = convex_hull_image(selem_rotated).astype(int) From 56686fe809aca16ce4a02d7238fd68b773543d04 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 15 Aug 2013 17:28:08 +0530 Subject: [PATCH 75/77] Reverting back to changes made by Johannes --- skimage/feature/censure.py | 60 +++++++++++++++++++------------------- 1 file changed, 30 insertions(+), 30 deletions(-) diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index a15195a2..dee2c81c 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -4,7 +4,7 @@ from scipy.ndimage.filters import maximum_filter, minimum_filter, convolve from skimage.transform import integral_image from skimage.feature.corner import _compute_auto_correlation from skimage.util import img_as_float -from skimage.morphology import convex_hull_image, octagon, star +from skimage.morphology import octagon, star from skimage.feature.util import _mask_border_keypoints from skimage.feature.censure_cy import _censure_dob_loop @@ -20,18 +20,17 @@ STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6), (11, 8), (13, 10), (14, 11), (15, 12), (16, 14)] -def _get_filtered_image(image, min_scale, max_scale, mode): +def _filter_image(image, min_scale, max_scale, mode): - scales = np.zeros((image.shape[0], image.shape[1], - max_scale - min_scale + 1), dtype=np.double) + response = np.zeros((image.shape[0], image.shape[1], + max_scale - min_scale + 1), dtype=np.double) if mode == 'dob': - # make scales[:, :, i] contiguous memory block - item_size = scales.itemsize - scales.strides = (item_size * scales.shape[0], - item_size, - item_size * scales.shape[0] * scales.shape[1]) + # make response[:, :, i] contiguous memory block + item_size = response.itemsize + response.strides = (item_size * response.shape[0], item_size, + item_size * response.shape[0] * response.shape[1]) integral_img = integral_image(image) @@ -39,12 +38,12 @@ def _get_filtered_image(image, min_scale, max_scale, mode): n = min_scale + i # Constant multipliers for the outer region and the inner region - # of the bilevel filters with the constraint of keeping the + # of the bi-level filters with the constraint of keeping the # DC bias 0. inner_weight = (1.0 / (2 * n + 1)**2) outer_weight = (1.0 / (12 * n**2 + 4 * n)) - _censure_dob_loop(n, integral_img, scales[:, :, i], + _censure_dob_loop(n, integral_img, response[:, :, i], inner_weight, outer_weight) # NOTE : For the Octagon shaped filter, we implemented and evaluated the @@ -58,17 +57,17 @@ def _get_filtered_image(image, min_scale, max_scale, mode): for i in range(max_scale - min_scale + 1): mo, no = OCTAGON_OUTER_SHAPE[min_scale + i - 1] mi, ni = OCTAGON_INNER_SHAPE[min_scale + i - 1] - scales[:, :, i] = convolve(image, - _octagon_filter_kernel(mo, no, mi, ni)) + response[:, :, i] = convolve(image, + _octagon_filter_kernel(mo, no, mi, ni)) + elif mode == 'star': for i in range(max_scale - min_scale + 1): m = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][0]] n = STAR_SHAPE[STAR_FILTER_SHAPE[min_scale + i - 1][1]] - scales[:, :, i] = convolve(image, - _star_filter_kernel(m, n)) + response[:, :, i] = convolve(image, _star_filter_kernel(m, n)) - return scales + return response def _octagon_filter_kernel(mo, no, mi, ni): @@ -106,24 +105,24 @@ def _suppress_lines(feature_mask, image, sigma, line_threshold): def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', non_max_threshold=0.15, line_threshold=10): """ - Extracts Censure keypoints along with the corresponding scale using - either Difference of Boxes, Octagon or STAR bilevel filter. + Extracts CenSurE keypoints along with the corresponding scale using + either Difference of Boxes, Octagon or STAR bi-level filter. Parameters ---------- image : 2D ndarray Input image. - min_scale : positive integer + min_scale : int Minimum scale to extract keypoints from. - max_scale : positive integer + max_scale : int Maximum scale to extract keypoints from. The keypoints will be extracted from all the scales except the first and the last i.e. from the scales in the range [min_scale + 1, max_scale - 1]. - mode : ('DoB', 'Octagon', 'STAR') - Type of bilevel filter used to get the scales of the input image. + mode : {'DoB', 'Octagon', 'STAR'} + Type of bi-level filter used to get the scales of the input image. Possible values are 'DoB', 'Octagon' and 'STAR'. The three modes - represent the shape of the bilevel filters i.e. box(square), octagon - and star respectively. For instance, a bilevel octagon filter consists + represent the shape of the bi-level filters i.e. box(square), octagon + and star respectively. For instance, a bi-level octagon filter consists of a smaller inner octagon and a larger outer octagon with the filter weights being uniformly negative in both the inner octagon while uniformly positive in the difference region. Use STAR and Octagon for @@ -138,7 +137,7 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', Returns ------- keypoints : (N, 2) array - Location of the extracted keypoints in the (row, col) format. + Location of the extracted keypoints in the ``(row, col)`` format. scales : (N, 1) array The corresponding scale of the N extracted keypoints. @@ -155,8 +154,9 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', http://www.jamris.org/01_2013/saveas.php?QUEST=JAMRIS_No01_2013_P_11-20.pdf """ + # (1) First we generate the required scales on the input grayscale image - # using a bilevel filter and stack them up in `filter_response`. + # using a bi-level filter and stack them up in `filter_response`. # (2) We then perform Non-Maximal suppression in 3 x 3 x 3 window on the # filter_response to suppress points that are neither minima or maxima in # 3 x 3 x 3 neighbourhood. We obtain a boolean ndarray `feature_mask` @@ -175,15 +175,15 @@ def keypoints_censure(image, min_scale=1, max_scale=7, mode='DoB', if mode not in ('dob', 'octagon', 'star'): raise ValueError('Mode must be one of "DoB", "Octagon", "STAR".') - if max_scale - min_scale < 2: - raise ValueError('The number of scales should be greater than or' - 'equal to 3.') + if min_scale < 1 or max_scale < 1 or max_scale - min_scale < 2: + raise ValueError('The scales must be >= 1 and the number of scales ' + 'should be >= 3.') image = img_as_float(image) image = np.ascontiguousarray(image) # Generating all the scales - filter_response = _get_filtered_image(image, min_scale, max_scale, mode) + filter_response = _filter_image(image, min_scale, max_scale, mode) # Suppressing points that are neither minima or maxima in their 3 x 3 x 3 # neighbourhood to zero From 550dfce134083ada14d23dd9ccf3baf87a1b75c5 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 15 Aug 2013 18:20:03 +0530 Subject: [PATCH 76/77] Final changes in censure example --- doc/examples/plot_censure_keypoints.py | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/doc/examples/plot_censure_keypoints.py b/doc/examples/plot_censure_keypoints.py index 95fa0e58..eb6f2a16 100644 --- a/doc/examples/plot_censure_keypoints.py +++ b/doc/examples/plot_censure_keypoints.py @@ -23,7 +23,7 @@ non_max_threshold = 0.15 line_threshold = 10 -f, ax = plt.subplots(nrows=(max_scale - min_scale - 1), ncols=3, +_, ax = plt.subplots(nrows=(max_scale - min_scale - 1), ncols=3, figsize=(6, 6)) plt.subplots_adjust(wspace=0.02, hspace=0.02, top=0.94, bottom=0.02, left=0.06, right=0.98) @@ -39,10 +39,9 @@ for col, mode in enumerate(['dob', 'octagon', 'star']): # Plotting Censure features at all the scales for row, scale in enumerate(range(min_scale + 1, max_scale)): - keypoints_i = keypoints[scales == scale] - num = len(keypoints_i) - x = keypoints_i[:, 1] - y = keypoints_i[:, 0] + mask = scales == scale + x = keypoints[mask][:, 1] + y = keypoints[mask][:, 0] s = 0.5 * 2 ** (scale + min_scale + 1) ax[row, col].imshow(img) ax[row, col].scatter(x, y, s, facecolors='none', edgecolors='b') From 9ad7c78c8e8cea70a8c27c66e5292e281f0d2d16 Mon Sep 17 00:00:00 2001 From: Ankit Agrawal Date: Thu, 15 Aug 2013 19:14:47 +0530 Subject: [PATCH 77/77] Documenting the choice of filter sizes --- doc/examples/plot_censure_keypoints.py | 4 ++-- skimage/feature/censure.py | 7 +++++++ 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/doc/examples/plot_censure_keypoints.py b/doc/examples/plot_censure_keypoints.py index eb6f2a16..fbba4e1b 100644 --- a/doc/examples/plot_censure_keypoints.py +++ b/doc/examples/plot_censure_keypoints.py @@ -40,8 +40,8 @@ for col, mode in enumerate(['dob', 'octagon', 'star']): # Plotting Censure features at all the scales for row, scale in enumerate(range(min_scale + 1, max_scale)): mask = scales == scale - x = keypoints[mask][:, 1] - y = keypoints[mask][:, 0] + x = keypoints[mask, 1] + y = keypoints[mask, 0] s = 0.5 * 2 ** (scale + min_scale + 1) ax[row, col].imshow(img) ax[row, col].scatter(x, y, s, facecolors='none', edgecolors='b') diff --git a/skimage/feature/censure.py b/skimage/feature/censure.py index dee2c81c..00be16a8 100644 --- a/skimage/feature/censure.py +++ b/skimage/feature/censure.py @@ -10,11 +10,18 @@ from skimage.feature.util import _mask_border_keypoints from skimage.feature.censure_cy import _censure_dob_loop +# The paper(Reference [1]) mentions the sizes of the Octagon shaped filter +# kernel for the first seven scales only. The sizes of the later scales +# have been extrapolated based on the following statement in the paper. +# "These octagons scale linearly and were experimentally chosen to correspond +# to the seven DOBs described in the previous section." OCTAGON_OUTER_SHAPE = [(5, 2), (5, 3), (7, 3), (9, 4), (9, 7), (13, 7), (15, 10), (15, 11), (15, 12), (17, 13), (17, 14)] OCTAGON_INNER_SHAPE = [(3, 0), (3, 1), (3, 2), (5, 2), (5, 3), (5, 4), (5, 5), (7, 5), (7, 6), (9, 6), (9, 7)] +# The sizes for the STAR shaped filter kernel for different scales have been +# taken from the OpenCV implementation. STAR_SHAPE = [1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128] STAR_FILTER_SHAPE = [(1, 0), (3, 1), (4, 2), (5, 3), (7, 4), (8, 5), (9, 6), (11, 8), (13, 10), (14, 11), (15, 12), (16, 14)]