Merge pull request #695 from ahojnnes/memoryviews

Use typed memoryviews for draw, feature, filter.
This commit is contained in:
Stefan van der Walt
2013-08-17 15:40:40 -07:00
7 changed files with 68 additions and 100 deletions
+5 -7
View File
@@ -47,8 +47,6 @@ def line(Py_ssize_t y, Py_ssize_t x, Py_ssize_t y2, Py_ssize_t x2):
"""
cdef cnp.ndarray[cnp.intp_t, ndim=1, mode="c"] rr, cc
cdef char steep = 0
cdef Py_ssize_t dx = abs(x2 - x)
cdef Py_ssize_t dy = abs(y2 - y)
@@ -69,8 +67,8 @@ def line(Py_ssize_t y, Py_ssize_t x, Py_ssize_t y2, Py_ssize_t x2):
sx, sy = sy, sx
d = (2 * dy) - dx
rr = np.zeros(int(dx) + 1, dtype=np.intp)
cc = np.zeros(int(dx) + 1, dtype=np.intp)
cdef Py_ssize_t[:] rr = np.zeros(int(dx) + 1, dtype=np.intp)
cdef Py_ssize_t[:] cc = np.zeros(int(dx) + 1, dtype=np.intp)
for i in range(dx):
if steep:
@@ -88,7 +86,7 @@ def line(Py_ssize_t y, Py_ssize_t x, Py_ssize_t y2, Py_ssize_t x2):
rr[dx] = y2
cc[dx] = x2
return rr, cc
return np.asarray(rr), np.asarray(cc)
def polygon(y, x, shape=None):
@@ -149,8 +147,8 @@ def polygon(y, x, shape=None):
# make contigous arrays for r, c coordinates
cdef cnp.ndarray contiguous_rdata, contiguous_cdata
contiguous_rdata = np.ascontiguousarray(y, 'double')
contiguous_cdata = np.ascontiguousarray(x, 'double')
contiguous_rdata = np.ascontiguousarray(y, dtype=np.double)
contiguous_cdata = np.ascontiguousarray(x, dtype=np.double)
cdef cnp.double_t* rptr = <cnp.double_t*>contiguous_rdata.data
cdef cnp.double_t* cptr = <cnp.double_t*>contiguous_cdata.data
+4 -4
View File
@@ -50,9 +50,9 @@ from skimage.transform import integral
def match_template(cnp.ndarray[float, ndim=2, mode="c"] image,
cnp.ndarray[float, ndim=2, mode="c"] template):
cdef cnp.ndarray[float, ndim=2, mode="c"] corr
cdef cnp.ndarray[float, ndim=2, mode="c"] image_sat
cdef cnp.ndarray[float, ndim=2, mode="c"] image_sqr_sat
cdef float[:, ::1] corr
cdef float[:, ::1] image_sat
cdef float[:, ::1] image_sqr_sat
cdef float template_mean = np.mean(template)
cdef float template_ssd
cdef float inv_area
@@ -94,4 +94,4 @@ def match_template(cnp.ndarray[float, ndim=2, mode="c"] image,
den = sqrt((window_sqr_sum - window_mean_sqr) * template_ssd)
corr[r, c] /= den
return corr
return np.asarray(corr)
+20 -24
View File
@@ -8,15 +8,9 @@ from libc.math cimport sin, cos, abs
from skimage._shared.interpolation cimport bilinear_interpolation
def _glcm_loop(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
negative_indices=False, mode='c'] image,
cnp.ndarray[dtype=cnp.float64_t, ndim=1,
negative_indices=False, mode='c'] distances,
cnp.ndarray[dtype=cnp.float64_t, ndim=1,
negative_indices=False, mode='c'] angles,
int levels,
cnp.ndarray[dtype=cnp.uint32_t, ndim=4,
negative_indices=False, mode='c'] out):
def _glcm_loop(cnp.uint8_t[:, ::1] image, double[:] distances,
double[:] angles, Py_ssize_t levels,
cnp.uint32_t[:, :, :, ::1] out):
"""Perform co-occurrence matrix accumulation.
Parameters
@@ -81,7 +75,7 @@ cdef inline int _bit_rotate_right(int value, int length):
return (value >> 1) | ((value & 1) << (length - 1))
def _local_binary_pattern(cnp.ndarray[double, ndim=2] image,
def _local_binary_pattern(double[:, ::1] image,
int P, float R, char method='D'):
"""Gray scale and rotation invariant LBP (Local Binary Patterns).
@@ -92,8 +86,8 @@ def _local_binary_pattern(cnp.ndarray[double, ndim=2] image,
image : (N, M) double array
Graylevel image.
P : int
Number of circularly symmetric neighbour set points (quantization of the
angular space).
Number of circularly symmetric neighbour set points (quantization of
the angular space).
R : float
Radius of circle (spatial resolution of the operator).
method : {'D', 'R', 'U', 'V'}
@@ -111,19 +105,20 @@ def _local_binary_pattern(cnp.ndarray[double, ndim=2] image,
"""
# texture weights
cdef cnp.ndarray[int, ndim=1] weights = 2 ** np.arange(P, dtype=np.int32)
cdef int[:] weights = 2 ** np.arange(P, dtype=np.int32)
# local position of texture elements
rp = - R * np.sin(2 * np.pi * np.arange(P, dtype=np.double) / P)
cp = R * np.cos(2 * np.pi * np.arange(P, dtype=np.double) / P)
cdef cnp.ndarray[double, ndim=2] coords = np.round(np.vstack([rp, cp]).T, 5)
rr = - R * np.sin(2 * np.pi * np.arange(P, dtype=np.double) / P)
cc = R * np.cos(2 * np.pi * np.arange(P, dtype=np.double) / P)
cdef double[:] rp = np.round(rr, 5)
cdef double[:] cp = np.round(cc, 5)
# pre allocate arrays for computation
cdef cnp.ndarray[double, ndim=1] texture = np.zeros(P, np.double)
cdef cnp.ndarray[char, ndim=1] signed_texture = np.zeros(P, np.int8)
cdef cnp.ndarray[int, ndim=1] rotation_chain = np.zeros(P, np.int32)
# pre-allocate arrays for computation
cdef double[:] texture = np.zeros(P, dtype=np.double)
cdef char[:] signed_texture = np.zeros(P, dtype=np.int8)
cdef int[:] rotation_chain = np.zeros(P, dtype=np.int32)
output_shape = (image.shape[0], image.shape[1])
cdef cnp.ndarray[double, ndim=2] output = np.zeros(output_shape, np.double)
cdef double[:, ::1] output = np.zeros(output_shape, dtype=np.double)
cdef Py_ssize_t rows = image.shape[0]
cdef Py_ssize_t cols = image.shape[1]
@@ -133,8 +128,9 @@ def _local_binary_pattern(cnp.ndarray[double, ndim=2] image,
for r in range(image.shape[0]):
for c in range(image.shape[1]):
for i in range(P):
texture[i] = bilinear_interpolation(<double*>image.data,
rows, cols, r + coords[i, 0], c + coords[i, 1], 'C', 0)
texture[i] = bilinear_interpolation(&image[0, 0], rows, cols,
r + rp[i], c + cp[i],
'C', 0)
# signed / thresholded texture
for i in range(P):
if texture[i] - image[r, c] >= 0:
@@ -181,4 +177,4 @@ def _local_binary_pattern(cnp.ndarray[double, ndim=2] image,
output[r, c] = lbp
return output
return np.asarray(output)
+6 -15
View File
@@ -59,16 +59,8 @@ def corner_moravec(image, Py_ssize_t window_size=1):
cdef Py_ssize_t rows = image.shape[0]
cdef Py_ssize_t cols = image.shape[1]
cdef cnp.ndarray[dtype=cnp.double_t, ndim=2, mode='c'] cimage, out
if image.ndim == 3:
cimage = rgb2grey(image)
cimage = np.ascontiguousarray(img_as_float(image))
out = np.zeros(image.shape, dtype=np.double)
cdef double* image_data = <double*>cimage.data
cdef double* out_data = <double*>out.data
cdef double[:, ::1] cimage = np.ascontiguousarray(img_as_float(image))
cdef double[:, ::1] out = np.zeros(image.shape, dtype=np.double)
cdef double msum, min_msum
cdef Py_ssize_t r, c, br, bc, mr, mc, a, b
@@ -81,11 +73,10 @@ def corner_moravec(image, Py_ssize_t window_size=1):
msum = 0
for mr in range(- window_size, window_size + 1):
for mc in range(- window_size, window_size + 1):
a = (r + mr) * cols + c + mc
b = (br + mr) * cols + bc + mc
msum += (image_data[a] - image_data[b]) ** 2
msum += (cimage[r + mr, c + mc]
- cimage[br + mr, bc + mc]) ** 2
min_msum = min(msum, min_msum)
out_data[r * cols + c] = min_msum
out[r, c] = min_msum
return out
return np.asarray(out)
+1 -1
View File
@@ -271,6 +271,6 @@ def local_binary_pattern(image, P, R, method='default'):
'uniform': ord('U'),
'var': ord('V')
}
image = np.array(image, dtype='double', copy=True)
image = np.ascontiguousarray(image, dtype=np.double)
output = _local_binary_pattern(image, P, R, methods[method.lower()])
return output
+7 -14
View File
@@ -731,13 +731,8 @@ cdef int c_median_filter(Py_ssize_t rows,
return 0
def median_filter(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
negative_indices=False, mode='c'] data,
cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
negative_indices=False, mode='c'] mask,
cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
negative_indices=False, mode='c'] output,
int radius,
def median_filter(cnp.uint8_t[:, ::1] data, cnp.uint8_t[:, ::1] mask,
cnp.uint8_t[:, ::1] output, int radius,
cnp.int32_t percent):
"""Median filter with octagon shape and masking.
@@ -773,12 +768,10 @@ def median_filter(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
raise ValueError('Data shape (%d, %d) is not output shape (%d, %d)' %
(data.shape[0], data.shape[1],
output.shape[0], output.shape[1]))
if c_median_filter(<cnp.int32_t>data.shape[0],
<cnp.int32_t>data.shape[1],
<cnp.int32_t>data.strides[0],
<cnp.int32_t>data.strides[1],
if c_median_filter(data.shape[0], data.shape[1],
data.strides[0], data.strides[1],
radius, percent,
<cnp.uint8_t*>data.data,
<cnp.uint8_t*>mask.data,
<cnp.uint8_t*>output.data):
&data[0, 0],
&mask[0, 0],
&output[0, 0]):
raise MemoryError('Failed to allocate scratchpad memory')
+25 -35
View File
@@ -113,11 +113,8 @@ def denoise_bilateral(image, Py_ssize_t win_size=5, sigma_range=None,
double max_value
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] cimage
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] out
double* image_data
double* out_data
double[:, :, ::1] cimage
double[:, :, ::1] out
double* color_lut
double* range_lut
@@ -143,8 +140,6 @@ def denoise_bilateral(image, Py_ssize_t win_size=5, sigma_range=None,
cimage = np.ascontiguousarray(image)
out = np.zeros((rows, cols, dims), dtype=np.double)
image_data = <double*>cimage.data
out_data = <double*>out.data
color_lut = _compute_color_lut(bins, csigma_range, max_value)
range_lut = _compute_range_lut(win_size, sigma_spatial)
dist_scale = bins / dims / max_value
@@ -159,11 +154,10 @@ def denoise_bilateral(image, Py_ssize_t win_size=5, sigma_range=None,
for r in range(rows):
for c in range(cols):
pixel_addr = r * cols * dims + c * dims
total_weight = 0
for d in range(dims):
total_values[d] = 0
centres[d] = image_data[pixel_addr + d]
centres[d] = cimage[r, c, d]
for wr in range(-window_ext, window_ext + 1):
rr = wr + r
kr = wr + window_ext
@@ -175,7 +169,7 @@ def denoise_bilateral(image, Py_ssize_t win_size=5, sigma_range=None,
# distance between centre stack and current position
dist = 0
for d in range(dims):
value = get_pixel3d(image_data, rows, cols, dims,
value = get_pixel3d(&cimage[0, 0, 0], rows, cols, dims,
rr, cc, d, cmode, cval)
values[d] = value
dist += (centres[d] - value)**2
@@ -189,7 +183,7 @@ def denoise_bilateral(image, Py_ssize_t win_size=5, sigma_range=None,
total_values[d] += values[d] * weight
total_weight += weight
for d in range(dims):
out_data[pixel_addr + d] = total_values[d] / total_weight
out[r, c, d] = total_values[d] / total_weight
free(color_lut)
free(range_lut)
@@ -197,11 +191,11 @@ def denoise_bilateral(image, Py_ssize_t win_size=5, sigma_range=None,
free(centres)
free(total_values)
return np.squeeze(out)
return np.squeeze(np.asarray(out))
def denoise_tv_bregman(image, double weight, int max_iter=100, double eps=1e-3,
isotropic=True):
char isotropic=True):
"""Perform total-variation denoising using split-Bregman optimization.
Total-variation denoising (also know as total-variation regularization)
@@ -258,21 +252,17 @@ def denoise_tv_bregman(image, double weight, int max_iter=100, double eps=1e-3,
Py_ssize_t total = rows * cols * dims
shape_ext = (rows2, cols2, dims)
shape_ext = (rows2, cols2, dims)
u = np.zeros(shape_ext, dtype=np.double)
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] cimage = \
np.ascontiguousarray(image)
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] u = \
np.zeros(shape_ext, dtype=np.double)
cdef:
double[:, :, ::1] cimage = np.ascontiguousarray(image)
double[:, :, ::1] cu = u
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] dx = \
np.zeros(shape_ext, dtype=np.double)
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] dy = \
np.zeros(shape_ext, dtype=np.double)
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] bx = \
np.zeros(shape_ext, dtype=np.double)
cnp.ndarray[dtype=cnp.double_t, ndim=3, mode='c'] by = \
np.zeros(shape_ext, dtype=np.double)
double[:, :, ::1] dx = np.zeros(shape_ext, dtype=np.double)
double[:, :, ::1] dy = np.zeros(shape_ext, dtype=np.double)
double[:, :, ::1] bx = np.zeros(shape_ext, dtype=np.double)
double[:, :, ::1] by = np.zeros(shape_ext, dtype=np.double)
double ux, uy, uprev, unew, bxx, byy, dxx, dyy, s
int i = 0
@@ -296,19 +286,19 @@ def denoise_tv_bregman(image, double weight, int max_iter=100, double eps=1e-3,
for r in range(1, rows + 1):
for c in range(1, cols + 1):
uprev = u[r, c, k]
uprev = cu[r, c, k]
# forward derivatives
ux = u[r, c + 1, k] - uprev
uy = u[r + 1, c, k] - uprev
ux = cu[r, c + 1, k] - uprev
uy = cu[r + 1, c, k] - uprev
# Gauss-Seidel method
unew = (
lam * (
+ u[r + 1, c, k]
+ u[r - 1, c, k]
+ u[r, c + 1, k]
+ u[r, c - 1, k]
+ cu[r + 1, c, k]
+ cu[r - 1, c, k]
+ cu[r, c + 1, k]
+ cu[r, c - 1, k]
+ dx[r, c - 1, k]
- dx[r, c, k]
@@ -321,7 +311,7 @@ def denoise_tv_bregman(image, double weight, int max_iter=100, double eps=1e-3,
+ by[r, c, k]
) + weight * cimage[r - 1, c - 1, k]
) / norm
u[r, c, k] = unew
cu[r, c, k] = unew
# update root mean square error
rmse += (unew - uprev)**2
@@ -360,4 +350,4 @@ def denoise_tv_bregman(image, double weight, int max_iter=100, double eps=1e-3,
rmse = sqrt(rmse / total)
i += 1
return np.squeeze(u[1:-1, 1:-1])
return np.squeeze(np.asarray(u[1:-1, 1:-1]))