Add fused type specialization for moments functions

This commit is contained in:
Johannes Schönberger
2014-11-26 17:03:50 -05:00
parent 802b22a62f
commit 71a5437794
4 changed files with 81 additions and 61 deletions
+2 -2
View File
@@ -46,9 +46,9 @@ Library:
Extension: skimage.measure._find_contours_cy
Sources:
skimage/measure/_find_contours_cy.pyx
Extension: skimage.measure._moments
Extension: skimage.measure._moments_cy
Sources:
skimage/measure/_moments.pyx
skimage/measure/_moments_cy.pyx
Extension: skimage.measure._marching_cubes_cy
Sources:
skimage/measure/_marching_cubes_cy.pyx
@@ -1,11 +1,8 @@
#cython: cdivision=True
#cython: boundscheck=False
#cython: nonecheck=False
#cython: wraparound=False
import numpy as np
from . import _moments_cy
cpdef moments(double[:, :] image, Py_ssize_t order=3):
def moments(image, order=3):
"""Calculate all raw image moments up to a certain order.
The following properties can be calculated from raw image moments:
@@ -17,7 +14,7 @@ cpdef moments(double[:, :] image, Py_ssize_t order=3):
Parameters
----------
image : 2D double array
image : 2D double or uint8 array
Rasterized shape as image.
order : int, optional
Maximum order of moments. Default is 3.
@@ -39,11 +36,10 @@ cpdef moments(double[:, :] image, Py_ssize_t order=3):
.. [4] http://en.wikipedia.org/wiki/Image_moment
"""
return moments_central(image, 0, 0, order)
return _moments_cy.moments_central(image, 0, 0, order)
cpdef moments_central(double[:, :] image, double cr, double cc,
Py_ssize_t order=3):
def moments_central(image, cr, cc, order=3):
"""Calculate all central image moments up to a certain order.
Note that central moments are translation invariant but not scale and
@@ -51,7 +47,7 @@ cpdef moments_central(double[:, :] image, double cr, double cc,
Parameters
----------
image : 2D double array
image : 2D double or uint8 array
Rasterized shape as image.
cr : double
Center row coordinate.
@@ -77,25 +73,11 @@ cpdef moments_central(double[:, :] image, double cr, double cc,
.. [4] http://en.wikipedia.org/wiki/Image_moment
"""
cdef Py_ssize_t p, q, r, c
cdef double[:, ::1] mu = np.zeros((order + 1, order + 1), dtype=np.double)
cdef double val, dr, dc, dcp, drq
for r in range(image.shape[0]):
dr = r - cr
for c in range(image.shape[1]):
dc = c - cc
val = image[r, c]
dcp = 1
for p in range(order + 1):
drq = 1
for q in range(order + 1):
mu[p, q] += val * drq * dcp
drq *= dr
dcp *= dc
return np.asarray(mu)
return _moments_cy.moments_central(image, cr, cc, order)
def moments_normalized(double[:, :] mu, Py_ssize_t order=3):
def moments_normalized(mu, order=3):
"""Calculate all normalized central image moments up to a certain order.
Note that normalized central moments are translation and scale invariant
@@ -125,18 +107,10 @@ def moments_normalized(double[:, :] mu, Py_ssize_t order=3):
.. [4] http://en.wikipedia.org/wiki/Image_moment
"""
cdef Py_ssize_t p, q
cdef double[:, ::1] nu = np.zeros((order + 1, order + 1), dtype=np.double)
for p in range(order + 1):
for q in range(order + 1):
if p + q >= 2:
nu[p, q] = mu[p, q] / mu[0, 0] ** (<double>(p + q) / 2 + 1)
else:
nu[p, q] = np.nan
return np.asarray(nu)
return _moments_cy.moments_normalized(mu.astype(np.double), order)
def moments_hu(double[:, :] nu):
def moments_hu(nu):
"""Calculate Hu's set of image moments.
Note that this set of moments is proofed to be translation, scale and
@@ -167,23 +141,4 @@ def moments_hu(double[:, :] nu):
"""
cdef double[::1] hu = np.zeros((7, ), dtype=np.double)
cdef double t0 = nu[3, 0] + nu[1, 2]
cdef double t1 = nu[2, 1] + nu[0, 3]
cdef double q0 = t0 * t0
cdef double q1 = t1 * t1
cdef double n4 = 4 * nu[1, 1]
cdef double s = nu[2, 0] + nu[0, 2]
cdef double d = nu[2, 0] - nu[0, 2]
hu[0] = s
hu[1] = d * d + n4 * nu[1, 1]
hu[3] = q0 + q1
hu[5] = d * (q0 - q1) + n4 * t0 * t1
t0 *= q0 - 3 * q1
t1 *= 3 * q0 - q1
q0 = nu[3, 0]- 3 * nu[1, 2]
q1 = 3 * nu[2, 1] - nu[0, 3]
hu[2] = q0 * q0 + q1 * q1
hu[4] = q0 * t0 + q1 * t1
hu[6] = q1 * t0 - q0 * t1
return np.asarray(hu)
return _moments_cy.moments_hu(nu.astype(np.double))
+65
View File
@@ -0,0 +1,65 @@
#cython: cdivision=True
#cython: boundscheck=False
#cython: nonecheck=False
#cython: wraparound=False
import numpy as np
import cython
ctypedef fused image_t:
cython.uchar[:, :]
cython.double[:, :]
cpdef moments_central(image_t image, double cr, double cc, Py_ssize_t order):
cdef Py_ssize_t p, q, r, c
cdef double[:, ::1] mu = np.zeros((order + 1, order + 1), dtype=np.double)
cdef double val, dr, dc, dcp, drq
for r in range(image.shape[0]):
dr = r - cr
for c in range(image.shape[1]):
dc = c - cc
val = image[r, c]
dcp = 1
for p in range(order + 1):
drq = 1
for q in range(order + 1):
mu[p, q] += val * drq * dcp
drq *= dr
dcp *= dc
return np.asarray(mu)
def moments_normalized(double[:, :] mu, Py_ssize_t order=3):
cdef Py_ssize_t p, q
cdef double[:, ::1] nu = np.zeros((order + 1, order + 1), dtype=np.double)
for p in range(order + 1):
for q in range(order + 1):
if p + q >= 2:
nu[p, q] = mu[p, q] / mu[0, 0] ** (<double>(p + q) / 2 + 1)
else:
nu[p, q] = np.nan
return np.asarray(nu)
def moments_hu(double[:, :] nu):
cdef double[::1] hu = np.zeros((7, ), dtype=np.double)
cdef double t0 = nu[3, 0] + nu[1, 2]
cdef double t1 = nu[2, 1] + nu[0, 3]
cdef double q0 = t0 * t0
cdef double q1 = t1 * t1
cdef double n4 = 4 * nu[1, 1]
cdef double s = nu[2, 0] + nu[0, 2]
cdef double d = nu[2, 0] - nu[0, 2]
hu[0] = s
hu[1] = d * d + n4 * nu[1, 1]
hu[3] = q0 + q1
hu[5] = d * (q0 - q1) + n4 * t0 * t1
t0 *= q0 - 3 * q1
t1 *= 3 * q0 - q1
q0 = nu[3, 0]- 3 * nu[1, 2]
q1 = 3 * nu[2, 1] - nu[0, 3]
hu[2] = q0 * q0 + q1 * q1
hu[4] = q0 * t0 + q1 * t1
hu[6] = q1 * t0 - q0 * t1
return np.asarray(hu)
+2 -2
View File
@@ -14,7 +14,7 @@ def configuration(parent_package='', top_path=None):
cython(['_ccomp.pyx'], working_path=base_path)
cython(['_find_contours_cy.pyx'], working_path=base_path)
cython(['_moments.pyx'], working_path=base_path)
cython(['_moments_cy.pyx'], working_path=base_path)
cython(['_marching_cubes_cy.pyx'], working_path=base_path)
cython(['_pnpoly.pyx'], working_path=base_path)
@@ -22,7 +22,7 @@ def configuration(parent_package='', top_path=None):
include_dirs=[get_numpy_include_dirs()])
config.add_extension('_find_contours_cy', sources=['_find_contours_cy.c'],
include_dirs=[get_numpy_include_dirs()])
config.add_extension('_moments', sources=['_moments.c'],
config.add_extension('_moments_cy', sources=['_moments_cy.c'],
include_dirs=[get_numpy_include_dirs()])
config.add_extension('_marching_cubes_cy',
sources=['_marching_cubes_cy.c'],