Merge pull request #33 from amueller/median_filter_minor

Median filter minor cleanup.
This commit is contained in:
Stefan van der Walt
2011-09-26 14:23:16 -07:00
2 changed files with 73 additions and 59 deletions
+25 -20
View File
@@ -15,12 +15,13 @@ import numpy as np
from . import _ctmf
from rank_order import rank_order
def median_filter(data, mask=None, radius=1, percent=50):
def median_filter(image, mask=None, radius=2, percent=50):
'''Masked median filter with octagon shape.
Parameters
----------
data : (M,N) ndarray, dtype uint8
image : (M,N) ndarray, dtype uint8
Input image.
mask : (M,N) ndarray, dtype uint8, optional
A value of 1 indicates a significant pixel, 0
@@ -43,42 +44,46 @@ def median_filter(data, mask=None, radius=1, percent=50):
'''
if image.ndim != 2:
raise TypeError("The input 'image' must be a two dimensional array.")
if mask is None:
mask = np.ones(data.shape, dtype=np.bool)
mask = np.ones(image.shape, dtype=np.bool)
mask = np.ascontiguousarray(mask, dtype=np.bool)
if np.all(~ mask):
return data.copy()
return image.copy()
#
# Normalize the ranked data to 0-255
# Normalize the ranked image to 0-255
#
if (not np.issubdtype(data.dtype, np.int) or
np.min(data) < 0 or np.max(data) > 255):
ranked_data,translation = rank_order(data[mask])
max_ranked_data = np.max(ranked_data)
if max_ranked_data == 0:
return data
if max_ranked_data > 255:
ranked_data = ranked_data * 255 // max_ranked_data
if (not np.issubdtype(image.dtype, np.int) or
np.min(image) < 0 or np.max(image) > 255):
ranked_image, translation = rank_order(image[mask])
max_ranked_image = np.max(ranked_image)
if max_ranked_image == 0:
return image
if max_ranked_image > 255:
ranked_image = ranked_image * 255 // max_ranked_image
was_ranked = True
else:
ranked_data = data[mask]
ranked_image = image[mask]
was_ranked = False
input = np.zeros(data.shape, np.uint8 )
input[mask] = ranked_data
input = np.zeros(image.shape, np.uint8)
input[mask] = ranked_image
mask.dtype = np.uint8
output = np.zeros(data.shape, np.uint8)
output = np.zeros(image.shape, np.uint8)
_ctmf.median_filter(input, mask, output, radius, percent)
if was_ranked:
#
# The translation gives the original value at each ranking.
# We rescale the output to the original ranking and then
# use the translation to look up the original value in the data.
# use the translation to look up the original value in the image.
#
if max_ranked_data > 255:
result = translation[output.astype(np.uint32) * max_ranked_data // 255]
if max_ranked_image > 255:
result = translation[output.astype(np.uint32) *
max_ranked_image // 255]
else:
result = translation[output]
else:
+48 -39
View File
@@ -1,50 +1,53 @@
import os.path
import numpy as np
from numpy.testing import *
from scikits.image.filter import median_filter
def test_00_00_zeros():
'''The median filter on an array of all zeros should be zero'''
result = median_filter(np.zeros((10,10)), np.ones((10,10),bool), 3)
result = median_filter(np.zeros((10, 10)), np.ones((10, 10), bool), 3)
assert np.all(result == 0)
def test_00_01_all_masked():
'''Test a completely masked image
Regression test of IMG-1029'''
result = median_filter(np.zeros((10,10)), np.zeros((10,10), bool), 3)
result = median_filter(np.zeros((10, 10)), np.zeros((10, 10), bool), 3)
assert (np.all(result == 0))
def test_00_02_all_but_one_masked():
mask = np.zeros((10,10), bool)
mask[5,5] = True
result = median_filter(np.zeros((10,10)), mask, 3)
mask = np.zeros((10, 10), bool)
mask[5, 5] = True
median_filter(np.zeros((10, 10)), mask, 3)
def test_01_01_mask():
'''The median filter, masking a single value'''
img = np.zeros((10,10))
img[5,5] = 1
mask = np.ones((10,10),bool)
mask[5,5] = False
img = np.zeros((10, 10))
img[5, 5] = 1
mask = np.ones((10, 10), bool)
mask[5, 5] = False
result = median_filter(img, mask, 3)
assert (np.all(result[mask] == 0))
assert_equal(result[5,5], 1)
np.testing.assert_equal(result[5, 5], 1)
def test_02_01_median():
'''A median filter larger than the image = median of image'''
np.random.seed(0)
img = np.random.uniform(size=(9,9))
result = median_filter(img, np.ones((9,9),bool), 20)
assert_equal(result[0,0], np.median(img))
img = np.random.uniform(size=(9, 9))
result = median_filter(img, np.ones((9, 9), bool), 20)
np.testing.assert_equal(result[0, 0], np.median(img))
assert (np.all(result == np.median(img)))
def test_02_02_median_bigger():
'''Use an image of more than 255 values to test approximation'''
np.random.seed(0)
img = np.random.uniform(size=(20,20))
result = median_filter(img, np.ones((20,20),bool),40)
img = np.random.uniform(size=(20, 20))
result = median_filter(img, np.ones((20, 20), bool), 40)
sorted = np.ravel(img)
sorted.sort()
min_acceptable = sorted[198]
@@ -52,54 +55,60 @@ def test_02_02_median_bigger():
assert (np.all(result >= min_acceptable))
assert (np.all(result <= max_acceptable))
def test_03_01_shape():
'''Make sure the median filter is the expected octagonal shape'''
radius = 5
a_2 = int(radius / 2.414213)
i,j = np.mgrid[-10:11,-10:11]
octagon = np.ones((21,21), bool)
i, j = np.mgrid[-10:11, -10:11]
octagon = np.ones((21, 21), bool)
#
# constrain the octagon mask to be the points that are on
# the correct side of the 8 edges
#
octagon[i < -radius] = False
octagon[i > radius] = False
octagon[i > radius] = False
octagon[j < -radius] = False
octagon[j > radius] = False
octagon[i+j < -radius-a_2] = False
octagon[j-i > radius+a_2] = False
octagon[i+j > radius+a_2] = False
octagon[i-j > radius+a_2] = False
octagon[j > radius] = False
octagon[i + j < -radius - a_2] = False
octagon[j - i > radius + a_2] = False
octagon[i + j > radius + a_2] = False
octagon[i - j > radius + a_2] = False
np.random.seed(0)
img = np.random.uniform(size=(21,21))
result = median_filter(img, np.ones((21,21),bool), radius)
img = np.random.uniform(size=(21, 21))
result = median_filter(img, np.ones((21, 21), bool), radius)
sorted = img[octagon]
sorted.sort()
min_acceptable = sorted[len(sorted)/2-1]
max_acceptable = sorted[len(sorted)/2+1]
assert (result[10,10] >= min_acceptable)
assert (result[10,10] <= max_acceptable)
min_acceptable = sorted[len(sorted) / 2 - 1]
max_acceptable = sorted[len(sorted) / 2 + 1]
assert (result[10, 10] >= min_acceptable)
assert (result[10, 10] <= max_acceptable)
def test_04_01_half_masked():
'''Make sure that the median filter can handle large masked areas.'''
img = np.ones((20, 20))
mask = np.ones((20, 20),bool)
mask = np.ones((20, 20), bool)
mask[10:, :] = False
img[~ mask] = 2
img[1, 1] = 0 # to prevent short circuit for uniform data.
img[1, 1] = 0 # to prevent short circuit for uniform data.
result = median_filter(img, mask, 5)
# in partial coverage areas, the result should be only from the masked pixels
# in partial coverage areas, the result should be only
# from the masked pixels
assert (np.all(result[:14, :] == 1))
# in zero coverage areas, the result should be the lowest valud in the valid area
# in zero coverage areas, the result should be the lowest
# value in the valid area
assert (np.all(result[15:, :] == np.min(img[mask])))
def test_default_values():
img = (np.random.random((20, 20)) * 255).astype(np.uint8)
mask = np.ones((20, 20), dtype=np.uint8)
result1 = median_filter(img, mask, radius=1, percent=50)
result1 = median_filter(img, mask, radius=2, percent=50)
result2 = median_filter(img)
assert_array_equal(result1, result2)
np.testing.assert_array_equal(result1, result2)
if __name__ == "__main__":
run_module_suite()
np.testing.run_module_suite()