From f456245d57aa13caeba4356b0f195e1a6be864c4 Mon Sep 17 00:00:00 2001
From: Vighnesh Birodkar <vighneshbirodkar@gmail.com>
Date: Fri, 14 Mar 2014 20:57:17 +0530
Subject: [PATCH] added cython file

---
 skimage/feature/_hessian_det_appx.pyx | 69 +++++++++++++++++++++++++++
 skimage/feature/setup.py              |  3 ++
 2 files changed, 72 insertions(+)
 create mode 100644 skimage/feature/_hessian_det_appx.pyx

diff --git a/skimage/feature/_hessian_det_appx.pyx b/skimage/feature/_hessian_det_appx.pyx
new file mode 100644
index 00000000..fa5f4cef
--- /dev/null
+++ b/skimage/feature/_hessian_det_appx.pyx
@@ -0,0 +1,69 @@
+import numpy as np
+cimport numpy as np
+
+from skimage.transform import integral_image, integrate
+from skimage import util
+
+
+cdef inline int  clip(np.int_t x, np.int_t low, np.int_t high):
+    if(x > high):
+        return high
+    if(x < low):
+        return low
+    return x
+
+
+cdef inline int integ(np.int_t[:, :] img, np.int_t r1, np.int_t c1, np.int_t rl, np.int_t cl):
+
+    r1 = clip(r1, 0, img.shape[0] - 1)
+    c1 = clip(c1, 0, img.shape[1] - 1)
+
+    r2 = clip(r1 + rl, 0, img.shape[0] - 1)
+    c2 = clip(c1 + cl, 0, img.shape[1] - 1)
+
+    cdef np.int_t r = img[r2, c2] + img[r1, c1] - img[r1, c2] - img[r2, c1]
+
+    if (r < 0):
+        return 0
+    return r
+
+
+def hessian_det_appx(np.ndarray[np.int_t, ndim=2] image, float sigma):
+
+    cdef np.int_t[:, :] img = image
+    cdef int size = int(3 * sigma)
+    cdef np.ndarray[np.float_t, ndim = 2] out = np.zeros_like(img).astype(np.float)
+
+    cdef int height = img.shape[0]
+    cdef int width = img.shape[1]
+    cdef int r, c
+    cdef int s2 = (size - 1) / 2
+    cdef int s3 = size / 3
+    cdef int l = size / 3
+    cdef int w = size
+
+    cdef float dxx, dyy, dxy
+
+    if not size % 2:
+        size += 1
+
+    for r in range(height):
+        for c in range(width):
+        
+            dxy =  integ(img, r - s3, c + 1, s3, s3) + \
+                   integ(img, r + 1, c - s3, s3, s3) - \
+                   integ(img, r - s3, c - s3, s3, s3) - \
+                   integ(img, r + 1, c + 1, s3, s3)
+            dxy = -dxy / w / w
+
+            dxx = integ(img, r - s3 + 1, c - s2, 2 * s3 - 1,w) - \
+                  integ(img, r - s3 + 1, c - s3 / 2, 2 * s3 - 1, s3) * 3
+            dxx = -dxx / w / w
+
+            dyy = integ(img, r - s2, c - s2 + 1, w, 2 * s3 - 1) - \
+                  integ(img, r - s3 / 2, c - s3 + 1, s3, 2 * s3 - 1) * 3
+            dyy = -dyy / w / w
+
+            out[r, c] = (dxx * dyy - 0.81 * (dxy * dxy))
+
+    return out
diff --git a/skimage/feature/setup.py b/skimage/feature/setup.py
index 9d8b4fbb..d7be3dcf 100644
--- a/skimage/feature/setup.py
+++ b/skimage/feature/setup.py
@@ -17,6 +17,7 @@ def configuration(parent_package='', top_path=None):
     cython(['orb_cy.pyx'], working_path=base_path)
     cython(['brief_cy.pyx'], working_path=base_path)
     cython(['_texture.pyx'], working_path=base_path)
+    cython(['_hessian_det_appx.pyx'], working_path=base_path)
 
     config.add_extension('corner_cy', sources=['corner_cy.c'],
                          include_dirs=[get_numpy_include_dirs()])
@@ -28,6 +29,8 @@ def configuration(parent_package='', top_path=None):
                          include_dirs=[get_numpy_include_dirs()])
     config.add_extension('_texture', sources=['_texture.c'],
                          include_dirs=[get_numpy_include_dirs(), '../_shared'])
+    config.add_extension('_hessian_det_appx', sources=['_hessian_det_appx.c'],
+                         include_dirs=[get_numpy_include_dirs()])
 
     return config