Merge pull request #844 from odebeir/add_rank_sum

Add sum filter to rank filters (minor increment)
This commit is contained in:
Johannes Schönberger
2013-12-18 01:11:01 -08:00
8 changed files with 276 additions and 5 deletions
+6 -3
View File
@@ -1,11 +1,11 @@
from .generic import (autolevel, bottomhat, equalize, gradient, maximum, mean,
subtract_mean, median, minimum, modal, enhance_contrast,
pop, threshold, tophat, noise_filter, entropy, otsu)
pop, threshold, tophat, noise_filter, entropy, otsu, sum)
from ._percentile import (autolevel_percentile, gradient_percentile,
mean_percentile, subtract_mean_percentile,
enhance_contrast_percentile, percentile,
pop_percentile, threshold_percentile)
from .bilateral import mean_bilateral, pop_bilateral
pop_percentile, sum_percentile, threshold_percentile)
from .bilateral import mean_bilateral, pop_bilateral, sum_bilateral
from skimage._shared.utils import deprecated
@@ -51,6 +51,9 @@ __all__ = ['autolevel',
'pop',
'pop_percentile',
'pop_bilateral',
'sum',
'sum_bilateral',
'sum_percentile',
'threshold',
'threshold_percentile',
'tophat',
+36
View File
@@ -310,6 +310,42 @@ def pop_percentile(image, selem, out=None, mask=None, shift_x=False,
image, selem, out=out, mask=mask, shift_x=shift_x,
shift_y=shift_y, p0=p0, p1=p1)
def sum_percentile(image, selem, out=None, mask=None, shift_x=False,
shift_y=False, p0=0, p1=1):
"""Return greyscale local sum of an image.
sum is computed on the given structuring element. Only levels between
percentiles [p0, p1] are used. Result is truncated (8bit or 16bit).
Parameters
----------
image : ndarray (uint8, uint16)
Image array.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray (same dtype as input)
If None, a new array will be allocated.
mask : ndarray
Mask array that defines (>0) area of the image included in the local
neighborhood. If None, the complete image is used (default).
shift_x, shift_y : int
Offset added to the structuring element center point. Shift is bounded
to the structuring element sizes (center must be inside the given
structuring element).
p0, p1 : float in [0, ..., 1]
Define the [p0, p1] percentile interval to be considered for computing
the value.
Returns
-------
out : ndarray (same dtype as input image)
Output image.
"""
return _apply(percentile_cy._sum,
image, selem, out=out, mask=mask, shift_x=shift_x,
shift_y=shift_y, p0=p0, p1=p1)
def threshold_percentile(image, selem, out=None, mask=None, shift_x=False,
shift_y=False, p0=0):
+60 -1
View File
@@ -30,7 +30,7 @@ from . import bilateral_cy
from .generic import _handle_input
__all__ = ['mean_bilateral', 'pop_bilateral']
__all__ = ['mean_bilateral', 'pop_bilateral', 'sum_bilateral']
def _apply(func, image, selem, out, mask, shift_x, shift_y, s0, s1,
@@ -155,3 +155,62 @@ def pop_bilateral(image, selem, out=None, mask=None, shift_x=False,
return _apply(bilateral_cy._pop, image, selem, out=out,
mask=mask, shift_x=shift_x, shift_y=shift_y, s0=s0, s1=s1)
def sum_bilateral(image, selem, out=None, mask=None, shift_x=False,
shift_y=False, s0=10, s1=10):
"""Apply a flat kernel bilateral filter.
This is an edge-preserving and noise reducing denoising filter. It averages
pixels based on their spatial closeness and radiometric similarity.
Spatial closeness is measured by considering only the local pixel
neighborhood given by a structuring element (selem).
Radiometric similarity is defined by the greylevel interval [g-s0, g+s1]
where g is the current pixel greylevel. Only pixels belonging to the
structuring element AND having a greylevel inside this interval are
summed. Return greyscale local bilateral sum of an image.
Result is truncated (8bit or 16bit).
Parameters
----------
image : ndarray (uint8, uint16)
Image array.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray (same dtype as input)
If None, a new array will be allocated.
mask : ndarray
Mask array that defines (>0) area of the image included in the local
neighborhood. If None, the complete image is used (default).
shift_x, shift_y : int
Offset added to the structuring element center point. Shift is bounded
to the structuring element sizes (center must be inside the given
structuring element).
s0, s1 : int
Define the [s0, s1] interval around the greyvalue of the center pixel
to be considered for computing the value.
Returns
-------
out : ndarray (same dtype as input image)
Output image.
See also
--------
skimage.filter.denoise_bilateral for a gaussian bilateral filter.
Examples
--------
>>> from skimage import data
>>> from skimage.morphology import disk
>>> from skimage.filter.rank import sum_bilateral
>>> # Load test image / cast to uint16
>>> ima = data.camera().astype(np.uint16)
>>> # bilateral filtering of cameraman image using a flat kernel
>>> bilat_ima = sum_bilateral(ima, disk(20), s0=10,s1=10)
"""
return _apply(bilateral_cy._sum, image, selem, out=out,
mask=mask, shift_x=shift_x, shift_y=shift_y, s0=s0, s1=s1)
+31
View File
@@ -47,6 +47,27 @@ cdef inline double _kernel_pop(Py_ssize_t* histo, double pop, dtype_t g,
else:
return 0
cdef inline double _kernel_sum(Py_ssize_t* histo, double pop, dtype_t g,
Py_ssize_t max_bin, Py_ssize_t mid_bin,
double p0, double p1,
Py_ssize_t s0, Py_ssize_t s1):
cdef Py_ssize_t i
cdef Py_ssize_t bilat_pop = 0
cdef Py_ssize_t sum = 0
if pop:
for i in range(max_bin):
if (g > (i - s0)) and (g < (i + s1)):
bilat_pop += histo[i]
sum += histo[i] * i
if bilat_pop:
return sum
else:
return 0
else:
return 0
def _mean(dtype_t[:, ::1] image,
char[:, ::1] selem,
@@ -68,3 +89,13 @@ def _pop(dtype_t[:, ::1] image,
_core(_kernel_pop[dtype_t], image, selem, mask, out,
shift_x, shift_y, 0, 0, s0, s1, max_bin)
def _sum(dtype_t[:, ::1] image,
char[:, ::1] selem,
char[:, ::1] mask,
dtype_t_out[:, ::1] out,
char shift_x, char shift_y, Py_ssize_t s0, Py_ssize_t s1,
Py_ssize_t max_bin):
_core(_kernel_sum[dtype_t], image, selem, mask, out,
shift_x, shift_y, 0, 0, s0, s1, max_bin)
+45 -1
View File
@@ -508,7 +508,6 @@ def pop(image, selem, out=None, mask=None, shift_x=False, shift_y=False):
Examples
--------
>>> # Local mean
>>> from skimage.morphology import square
>>> import skimage.filter.rank as rank
>>> ima = 255 * np.array([[0, 0, 0, 0, 0],
@@ -528,6 +527,51 @@ def pop(image, selem, out=None, mask=None, shift_x=False, shift_y=False):
return _apply(generic_cy._pop, image, selem, out=out,
mask=mask, shift_x=shift_x, shift_y=shift_y)
def sum(image, selem, out=None, mask=None, shift_x=False, shift_y=False):
"""Return the sum of pixels inside the neighborhood (truncated to uint8 or uint16).
Parameters
----------
image : ndarray (uint8, uint16)
Image array.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray (same dtype as input)
If None, a new array will be allocated.
mask : ndarray
Mask array that defines (>0) area of the image included in the local
neighborhood. If None, the complete image is used (default).
shift_x, shift_y : int
Offset added to the structuring element center point. Shift is bounded
to the structuring element sizes (center must be inside the given
structuring element).
Returns
-------
out : ndarray (same dtype as input image)
Output image.
Examples
--------
>>> from skimage.morphology import square
>>> import skimage.filter.rank as rank
>>> ima = np.array([[0, 0, 0, 0, 0],
... [0, 1, 1, 1, 0],
... [0, 1, 1, 1, 0],
... [0, 1, 1, 1, 0],
... [0, 0, 0, 0, 0]], dtype=np.uint8)
>>> rank.sum(ima, square(3))
array([[1, 2, 3, 2, 1],
[2, 4, 6, 4, 2],
[3, 6, 9, 6, 3],
[2, 4, 6, 4, 2],
[1, 2, 3, 2, 1]], dtype=uint8)
"""
return _apply(generic_cy._sum, image, selem, out=out,
mask=mask, shift_x=shift_x, shift_y=shift_y)
def threshold(image, selem, out=None, mask=None, shift_x=False, shift_y=False):
"""Return greyscale local threshold of an image.
+24
View File
@@ -221,6 +221,21 @@ cdef inline double _kernel_pop(Py_ssize_t* histo, double pop, dtype_t g,
return pop
cdef inline double _kernel_sum(Py_ssize_t* histo, double pop,dtype_t g,
Py_ssize_t max_bin, Py_ssize_t mid_bin,
double p0, double p1,
Py_ssize_t s0, Py_ssize_t s1):
cdef Py_ssize_t i
cdef Py_ssize_t sum = 0
if pop:
for i in range(max_bin):
sum += histo[i] * i
return sum
else:
return 0
cdef inline double _kernel_threshold(Py_ssize_t* histo, double pop, dtype_t g,
Py_ssize_t max_bin, Py_ssize_t mid_bin,
@@ -455,6 +470,15 @@ def _pop(dtype_t[:, ::1] image,
_core(_kernel_pop[dtype_t], image, selem, mask, out,
shift_x, shift_y, 0, 0, 0, 0, max_bin)
def _sum(dtype_t[:, ::1] image,
char[:, ::1] selem,
char[:, ::1] mask,
dtype_t_out[:, ::1] out,
char shift_x, char shift_y, Py_ssize_t max_bin):
_core(_kernel_sum[dtype_t], image, selem, mask,
out, shift_x, shift_y, 0, 0, 0, 0, max_bin)
def _threshold(dtype_t[:, ::1] image,
char[:, ::1] selem,
+32
View File
@@ -90,6 +90,29 @@ cdef inline double _kernel_mean(Py_ssize_t* histo, double pop, dtype_t g,
else:
return 0
cdef inline double _kernel_sum(Py_ssize_t* histo, double pop, dtype_t g,
Py_ssize_t max_bin, Py_ssize_t mid_bin,
double p0, double p1,
Py_ssize_t s0, Py_ssize_t s1):
cdef Py_ssize_t i, sum, sum_g, n
if pop:
sum = 0
sum_g = 0
n = 0
for i in range(max_bin):
sum += histo[i]
if (sum >= p0 * pop) and (sum <= p1 * pop):
n += histo[i]
sum_g += histo[i] * i
if n > 0:
return sum_g
else:
return 0
else:
return 0
cdef inline double _kernel_subtract_mean(Py_ssize_t* histo, double pop,
dtype_t g,
@@ -245,6 +268,15 @@ def _mean(dtype_t[:, ::1] image,
_core(_kernel_mean[dtype_t], image, selem, mask, out,
shift_x, shift_y, p0, p1, 0, 0, max_bin)
def _sum(dtype_t[:, ::1] image,
char[:, ::1] selem,
char[:, ::1] mask,
dtype_t_out[:, ::1] out,
char shift_x, char shift_y, double p0, double p1,
Py_ssize_t max_bin):
_core(_kernel_sum[dtype_t], image, selem, mask, out,
shift_x, shift_y, p0, p1, 0, 0, max_bin)
def _subtract_mean(dtype_t[:, ::1] image,
char[:, ::1] selem,
+42
View File
@@ -498,6 +498,48 @@ def test_percentile_median():
img_max = rank.median(img16, selem=selem)
assert_array_equal(img_p0, img_max)
def test_sum():
# check the number of valid pixels in the neighborhood
image8 = np.array([[0, 0, 0, 0, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 0, 0, 0, 0]], dtype=np.uint8)
image16 = 400*np.array([[0, 0, 0, 0, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 0, 0, 0, 0]], dtype=np.uint16)
elem = np.ones((3, 3), dtype=np.uint8)
out8 = np.empty_like(image8)
out16 = np.empty_like(image16)
mask = np.ones(image8.shape, dtype=np.uint8)
r = np.array([[1, 2, 3, 2, 1],
[2, 4, 6, 4, 2],
[3, 6, 9, 6, 3],
[2, 4, 6, 4, 2],
[1, 2, 3, 2, 1]], dtype=np.uint8)
rank.sum(image=image8, selem=elem, out=out8, mask=mask)
assert_array_equal(r, out8)
rank.sum_percentile(image=image8, selem=elem, out=out8, mask=mask,p0=.0,p1=1.)
assert_array_equal(r, out8)
rank.sum_bilateral(image=image8, selem=elem, out=out8, mask=mask,s0=255,s1=255)
assert_array_equal(r, out8)
r = 400* np.array([[1, 2, 3, 2, 1],
[2, 4, 6, 4, 2],
[3, 6, 9, 6, 3],
[2, 4, 6, 4, 2],
[1, 2, 3, 2, 1]], dtype=np.uint16)
rank.sum(image=image16, selem=elem, out=out16, mask=mask)
assert_array_equal(r, out16)
rank.sum_percentile(image=image16, selem=elem, out=out16, mask=mask,p0=.0,p1=1.)
assert_array_equal(r, out16)
rank.sum_bilateral(image=image16, selem=elem, out=out16, mask=mask,s0=1000,s1=1000)
assert_array_equal(r, out16)
if __name__ == "__main__":
run_module_suite()