Add Moravec interest point detection

skimage/feature/_interest.pyx

@@ -0,0 +1,39 @@
+#cython: cdivision=True
+#cython: boundscheck=False
+#cython: nonecheck=False
+#cython: wraparound=False
+import numpy as np
+cimport numpy as cnp
+from libc.float cimport DBL_MAX
+
+
+def moravec(cnp.ndarray[dtype=cnp.double_t, ndim=2, mode='c'] image,
+             int block_size):
+    cdef int rows = image.shape[0]
+    cdef int cols = image.shape[1]
+
+    cdef cnp.ndarray[dtype=cnp.double_t, ndim=2, mode='c'] out = \
+         np.zeros_like(image)
+
+    cdef double* image_data = <double*>image.data
+    cdef double* out_data = <double*>out.data
+
+    cdef double msum, min_msum
+    cdef int r, c, br, bc, mr, mc, a, b
+    for r in range(2 * block_size, rows - 2 * block_size):
+        for c in range(2 * block_size, cols - 2 * block_size):
+            min_msum = DBL_MAX
+            for br in range(r - block_size, r + block_size + 1):
+                for bc in range(c - block_size, c + block_size + 1):
+                    if br != r and bc != c:
+                        msum = 0
+                        for mr in range(- block_size, block_size + 1):
+                            for mc in range(- block_size, block_size + 1):
+                                a = (r + mr) * cols + c + mc
+                                b = (br + mr) * cols + bc + mc
+                                msum += (image_data[a] - image_data[b]) ** 2
+                        min_msum = min(msum, min_msum)
+
+            out_data[r * cols + c] = min_msum
+
+    return out

skimage/feature/interest.py

@@ -2,7 +2,58 @@ import numpy as np
 from scipy import ndimage
 from skimage.color import rgb2grey
 from skimage.util import img_as_float
-from . import peak
+from . import peak, _interest
+
+
+def moravec(image, block_size=3, mode='constant', cval=0):
+    """Compute Moravec response image.
+
+    This interest operator is comparatively fast but not rotation invariant.
+
+    Parameters
+    ----------
+    image : ndarray
+        Input image.
+    block_size : int, optional
+        Block size for mean filtering the squared gradients.
+    mode : {'reflect', 'constant', 'nearest', 'mirror', 'wrap'}, optional
+        The mode parameter determines how the array borders are handled, where
+        cval is the value when mode is equal to 'constant'.
+    cval : double, optional
+        Constant value to use for constant mode.
+
+    Returns
+    -------
+    coordinates : (N, 2) array
+        `(row, column)` coordinates of interest points.
+
+    Examples
+    -------
+    >>> from skimage.feature import moravec, peak_local_max
+    >>> square = np.zeros([10, 10])
+    >>> square[1:9,1:9] = 1
+    >>> square
+    array([[ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
+           [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
+           [ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.],
+           [ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.],
+           [ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.],
+           [ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.],
+           [ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.],
+           [ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.],
+           [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
+           [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])
+    >>> moravec(square)), square.shape)
+    (2, 6)
+
+    """
+
+    if image.ndim == 3:
+        image =  rgb2grey(image)
+
+    image = np.ascontiguousarray(img_as_float(image))
+
+    return _corner._moravec(image, block_size)
 
 
 def harris(image, eps=1e-6, gaussian_deviation=1):

skimage/feature/setup.py

@@ -12,9 +12,12 @@ def configuration(parent_package='', top_path=None):
     config = Configuration('feature', parent_package, top_path)
     config.add_data_dir('tests')
 
+    cython(['_interest.pyx'], working_path=base_path)
     cython(['_texture.pyx'], working_path=base_path)
     cython(['_template.pyx'], working_path=base_path)
 
+    config.add_extension('_interest', sources=['_interest.c'],
+                         include_dirs=[get_numpy_include_dirs()])
     config.add_extension('_texture', sources=['_texture.c'],
                          include_dirs=[get_numpy_include_dirs(), '../_shared'])
     config.add_extension('_template', sources=['_template.c'],