Merge branch 'damian-morphology' of git://github.com/deads/scikits.image into damian

Conflicts:
	scikits/image/setup.py
This commit is contained in:
Stefan van der Walt
2009-11-20 11:27:09 +02:00
24 changed files with 4622 additions and 2 deletions
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
+2
View File
@@ -0,0 +1,2 @@
from grey import *
from selem import *
File diff suppressed because it is too large Load Diff
+103
View File
@@ -0,0 +1,103 @@
"""
:author: Damian Eads, 2009
:license: modified BSD
"""
from __future__ import division
import numpy as np
cimport numpy as np
cimport cython
STREL_DTYPE = np.uint8
ctypedef np.uint8_t STREL_DTYPE_t
IMAGE_DTYPE = np.uint8
ctypedef np.uint8_t IMAGE_DTYPE_t
cdef inline int int_max(int a, int b): return a if a >= b else b
cdef inline int int_min(int a, int b): return a if a <= b else b
@cython.boundscheck(False)
def dilate(np.ndarray[IMAGE_DTYPE_t, ndim=2] image not None,
np.ndarray[IMAGE_DTYPE_t, ndim=2] selem not None,
np.ndarray[IMAGE_DTYPE_t, ndim=2] out):
cdef int hw = selem.shape[0] / 2
cdef int hh = selem.shape[1] / 2
cdef int width = image.shape[0], height = image.shape[1]
if out is None:
out = np.zeros([width, height], dtype=IMAGE_DTYPE)
assert out.shape[0] == image.shape[0]
assert out.shape[1] == image.shape[1]
cdef int x, y, ix, iy, cx, cy
cdef IMAGE_DTYPE_t max_so_far
cdef int sw = selem.shape[0], sh = selem.shape[1]
cdef np.ndarray[np.int_t, ndim=2] xinc = np.zeros([sw, sh], dtype=np.int)
cdef np.ndarray[np.int_t, ndim=2] yinc = np.zeros([sw, sh], dtype=np.int)
for x in range(sw):
for y in range(sh):
xinc[x, y] = (x - hw)
yinc[x, y] = (y - hh)
for x in range(width):
for y in range(height):
max_so_far = 0
for cx in range(0, sw):
for cy in range(0, sh):
ix = x + xinc[cx,cy]
iy = y + yinc[cx,cy]
if ix>=0 and iy>=0 and ix < width and iy < height \
and selem[cx, cy] == 1 \
and image[ix,iy] > max_so_far:
max_so_far = image[ix,iy]
out[x,y] = max_so_far
return out
@cython.boundscheck(False)
def erode(np.ndarray[IMAGE_DTYPE_t, ndim=2] image not None,
np.ndarray[IMAGE_DTYPE_t, ndim=2] selem not None,
np.ndarray[IMAGE_DTYPE_t, ndim=2] out):
cdef int hw = selem.shape[0] / 2
cdef int hh = selem.shape[1] / 2
cdef int width = image.shape[0], height = image.shape[1]
if out is None:
out = np.zeros([width, height], dtype=IMAGE_DTYPE)
assert out.shape[0] == image.shape[0]
assert out.shape[1] == image.shape[1]
cdef int x, y, ix, iy, cx, cy
cdef IMAGE_DTYPE_t min_so_far
cdef int sw = selem.shape[0], sh = selem.shape[1]
cdef np.ndarray[np.int_t, ndim=2] xinc = np.zeros([sw, sh], dtype=np.int)
cdef np.ndarray[np.int_t, ndim=2] yinc = np.zeros([sw, sh], dtype=np.int)
for x in range(sw):
for y in range(sh):
xinc[x, y] = (x - hw)
yinc[x, y] = (y - hh)
for x in range(width):
for y in range(height):
min_so_far = 255
for cx in range(0, sw):
for cy in range(0, sh):
ix = x + xinc[cx,cy]
iy = y + yinc[cx,cy]
if ix>=0 and iy>=0 and ix < width \
and iy < height and selem[cx, cy] == 1 \
and image[ix,iy] < min_so_far:
min_so_far = image[ix,iy]
out[x,y] = min_so_far
return out
+187
View File
@@ -0,0 +1,187 @@
"""
:author: Damian Eads, 2009
:license: modified BSD
"""
__docformat__ = 'restructuredtext en'
import numpy as np
eps = np.finfo(float).eps
def greyscale_erode(image, selem, out=None):
"""
Performs a greyscale morphological erosion on an image given a
structuring element. The eroded pixel at (i,j) is the minimum
over all pixels in the neighborhood centered at (i,j).
Parameters
----------
image : ndarray
The image as an ndarray.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray
The array to store the result of the morphology. If None is
passed, a new array will be allocated.
Returns
-------
eroded : ndarray
The result of the morphological erosion.
"""
if image is out:
raise NotImplementedError("In-place erosion not supported!")
try:
import scikits.image.morphology.cmorph as cmorph
out = cmorph.erode(image, selem, out=out)
return out;
except ImportError:
raise ImportError("cmorph extension not available.")
def greyscale_dilate(image, selem, out=None):
"""
Performs a greyscale morphological dilation on an image given a
structuring element. The dilated pixel at (i,j) is the maximum
over all pixels in the neighborhood centered at (i,j).
Parameters
----------
image : ndarray
The image as an ndarray.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray
The array to store the result of the morphology. If None, is
passed, a new array will be allocated.
Returns
-------
dilated : ndarray
The result of the morphological dilation.
"""
if image is out:
raise NotImplementedError("In-place dilation not supported!")
try:
import cmorph
out = cmorph.dilate(image, selem, out=out)
return out;
except ImportError:
raise ImportError("cmorph extension not available.")
def greyscale_open(image, selem, out=None):
"""
Performs a greyscale morphological opening on an image given a
structuring element defined as a erosion followed by a dilation.
Parameters
----------
image : ndarray
The image as an ndarray.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray
The array to store the result of the morphology. If None
is passed, a new array will be allocated.
Returns
-------
opening : ndarray
The result of the morphological opening.
"""
eroded = greyscale_erode(image, selem)
out = greyscale_dilate(eroded, selem, out=out)
return out
def greyscale_close(image, selem, out=None):
"""
Performs a greyscale morphological closing on an image given a
structuring element defined as a dilation followed by an erosion.
Parameters
----------
image : ndarray
The image as an ndarray.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray
The array to store the result of the morphology. If None,
is passed, a new array will be allocated.
Returns
-------
opening : ndarray
The result of the morphological opening.
"""
dilated = greyscale_dilate(image, selem)
out = greyscale_erode(dilated, selem, out=out)
return out
def greyscale_white_top_hat(image, selem, out=None):
"""
Applies a white top hat on an image given a structuring element.
Parameters
----------
image : ndarray
The image as an ndarray.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray
The array to store the result of the morphology. If None
is passed, a new array will be allocated.
Returns
-------
opening : ndarray
The result of the morphological white top hat.
"""
if image is out:
raise NotImplementedError("Cannot perform white top hat in place.")
eroded = greyscale_erode(image, selem)
out = greyscale_dilate(eroded, selem, out=out)
out = image - out
return out
def greyscale_black_top_hat(image, selem, out=None):
"""
Applies a black top hat on an image given a structuring element.
Parameters
----------
image : ndarray
The image as an ndarray.
selem : ndarray
The neighborhood expressed as a 2-D array of 1's and 0's.
out : ndarray
The array to store the result of the morphology. If None
is passed, a new array will be allocated.
Returns
-------
opening : ndarray
The result of the black top filter.
"""
if image is out:
raise NotImplementedError("Cannot perform white top hat in place.")
dilated = greyscale_dilate(image, selem)
out = greyscale_erode(dilated, selem, out=out)
out = out - image
if image is out:
raise NotImplementedError("Cannot perform black top hat in place.")
return out
+113
View File
@@ -0,0 +1,113 @@
"""
:author: Damian Eads, 2009
:license: modified BSD
"""
import numpy as np
def square(width, dtype=np.uint8):
"""
Generates a flat, square-shaped structuring element. Every pixel
along the perimeter has a chessboard distance no greater than radius
(radius=floor(width/2)) pixels.
Parameters
----------
width : int
The width and height of the square
Additional Parameters
---------------------
dtype : data-type
The data type of the structuring element.
Returns
-------
selem : ndarray
A structuring element consisting only of ones, i.e. every
pixel belongs to the neighborhood.
"""
return np.ones((width, width), dtype=dtype)
def rectangle(width, height, dtype=np.uint8):
"""
Generates a flat, rectangular-shaped structuring element of a
given width and height. Every pixel in the rectangle belongs
to the neighboorhood.
Parameters
----------
width : int
The width of the rectangle
height : int
The height of the rectangle
Additional Parameters
---------------------
dtype : data-type
The data type of the structuring element.
Returns
-------
selem : ndarray
A structuring element consisting only of ones, i.e. every
pixel belongs to the neighborhood.
"""
return np.ones((width, height), dtype=dtype)
def diamond(radius, dtype=np.uint8):
"""
Generates a flat, diamond-shaped structuring element of a given
radius. A pixel is part of the neighborhood (i.e. labeled 1) iff
the city block/manhattan distance between it and the center of the
neighborhood is no greater than radius.
Parameters
----------
radius : string
The radius of the disk-shaped structuring element.
dtype : data-type
The data type of the structuring element.
Returns
-------
selem : ndarray
The structuring element where elements of the neighborhood
are 1 and 0 otherwise.
"""
half = radius
(I, J) = np.meshgrid(xrange(0, radius*2+1), xrange(0, radius*2+1))
s = np.abs(I-half)+np.abs(J-half)
return np.array(s <= radius, dtype=dtype)
def disk(radius, dtype=np.uint8):
"""
Generates a flat, disk-shaped structuring element of a given radius.
A pixel is within the neighborhood iff the euclidean distance between
it and the origin is no greater than a radius.
Parameters
----------
radius : string
The radius of the disk-shaped structuring element.
dtype : data-type
The data type of the structuring element.
Returns
-------
selem : ndarray
The structuring element where elements of the neighborhood
are 1 and 0 otherwise.
"""
L = np.linspace(-radius, radius, 2*radius+1)
(X, Y) = np.meshgrid(L, L)
s = X**2
s += Y**2
return np.array(s <= radius * radius, dtype=dtype)
+90
View File
@@ -0,0 +1,90 @@
#!/usr/bin/env python
import os
import shutil
import hashlib
base_path = os.path.dirname(__file__)
def same_cython(f0, f1):
'''Compare two Cython generated C-files, based on their md5-sum.
Returns True if the files are identical, False if not. The first
lines are skipped, due to the timestamp printed there.
'''
def md5sum(f):
m = hashlib.new('md5')
while True:
d = f.read(8096)
if not d:
break
m.update(d)
return m.hexdigest()
f0 = file(f0)
f0.readline()
f1 = file(f1)
f1.readline()
return md5sum(f0) == md5sum(f1)
def configuration(parent_package='', top_path=None):
from numpy.distutils.misc_util import Configuration, get_numpy_include_dirs
config = Configuration('morphology', parent_package, top_path)
config.add_data_dir('tests')
# since distutils/cython has problems, we'll check to see if cython is
# installed and use that to rebuild the .c files, if not, we'll just build
# directly from the included .c files
cython_files = ['cmorph.pyx']
try:
import Cython
for pyxfile in [os.path.join(base_path, f) for f in cython_files]:
# make a backup of the good c files
c_file = pyxfile.rstrip('pyx') + 'c'
c_file_new = c_file + '.new'
# run cython compiler
os.system('cython -o %s %s' % (c_file_new, pyxfile))
# if the resulting file is small, cython compilation failed
size = os.path.getsize(c_file_new)
if size < 100:
print "Cython compilation of %s failed. Using " \
"pre-generated file." % os.path.basename(pyxfile)
continue
# if the generated .c file differs from the one provided,
# use that one instead
if not same_cython(c_file_new, c_file):
shutil.copy(c_file_new, c_file)
except ImportError:
# if cython is not found, we just build from the included .c files
pass
for pyxfile in cython_files:
c_file = pyxfile.rstrip('pyx') + 'c'
config.add_extension(pyxfile.rstrip('.pyx'),
sources=[c_file],
include_dirs=[get_numpy_include_dirs()])
return config
if __name__ == '__main__':
from numpy.distutils.core import setup
setup(maintainer = 'Scikits.Image Developers',
author = 'Damian Eads',
maintainer_email = 'scikits-image@googlegroups.com',
description = 'Morphology Wrapper',
url = 'http://stefanv.github.com/scikits.image/',
license = 'SciPy License (BSD Style)',
**(configuration(top_path='').todict())
)
@@ -0,0 +1,52 @@
import os.path
import numpy as np
from numpy.testing import *
from scikits.image import data_dir
from scikits.image.io import imread
from scikits.image import data_dir
from scikits.image.morphology import *
lena = np.load(os.path.join(data_dir, 'lena_GRAY_U8.npy'))
class TestMorphology():
def morph_worker(self, img, fn, morph_func, strel_func):
matlab_results = np.load(os.path.join(data_dir, fn))
k = 0
for expected_result in matlab_results:
mask = strel_func(k)
actual_result = morph_func(lena, mask)
assert_equal(expected_result, actual_result)
k = k + 1
def test_erode_diamond(self):
self.morph_worker(lena, "diamond-erode-matlab-output.npy", greyscale_erode, diamond)
def test_dilate_diamond(self):
self.morph_worker(lena, "diamond-dilate-matlab-output.npy", greyscale_dilate, diamond)
def test_open_diamond(self):
self.morph_worker(lena, "diamond-open-matlab-output.npy", greyscale_open, diamond)
def test_close_diamond(self):
self.morph_worker(lena, "diamond-close-matlab-output.npy", greyscale_close, diamond)
def test_tophat_diamond(self):
self.morph_worker(lena, "diamond-tophat-matlab-output.npy", greyscale_white_top_hat, diamond)
def test_bothat_diamond(self):
self.morph_worker(lena, "diamond-bothat-matlab-output.npy", greyscale_black_top_hat, diamond)
def test_erode_disk(self):
self.morph_worker(lena, "disk-erode-matlab-output.npy", greyscale_erode, disk)
def test_dilate_disk(self):
self.morph_worker(lena, "disk-dilate-matlab-output.npy", greyscale_dilate, disk)
def test_open_disk(self):
self.morph_worker(lena, "disk-open-matlab-output.npy", greyscale_open, disk)
def test_close_disk(self):
self.morph_worker(lena, "disk-close-matlab-output.npy", greyscale_close, disk)
@@ -0,0 +1,43 @@
# Author: Damian Eads
import os.path
import numpy as np
from numpy.testing import *
from scikits.image import data_dir
from scikits.image.io import *
from scikits.image import data_dir
from scikits.image.morphology import *
class TestSElem():
def test_square_selem(self):
for k in xrange(0, 5):
actual_mask = selem.square(k)
expected_mask = np.ones((k, k), dtype='uint8')
assert_equal(expected_mask, actual_mask)
def test_rectangle_selem(self):
for i in xrange(0, 5):
for j in xrange(0, 5):
actual_mask = selem.rectangle(i, j)
expected_mask = np.ones((i, j), dtype='uint8')
assert_equal(expected_mask, actual_mask)
def strel_worker(self, fn, func):
matlab_masks = np.load(os.path.join(data_dir, fn))
k = 0
for expected_mask in matlab_masks:
actual_mask = func(k)
if (expected_mask.shape == (1,)):
expected_mask = expected_mask[:,np.newaxis]
assert_equal(expected_mask, actual_mask)
k = k + 1
def test_selem_disk(self):
self.strel_worker("disk-matlab-output.npy", selem.disk)
def test_selem_diamond(self):
self.strel_worker("diamond-matlab-output.npy", selem.diamond)
+1
View File
@@ -8,6 +8,7 @@ def configuration(parent_package='', top_path=None):
config.add_subpackage('opencv')
config.add_subpackage('graph')
config.add_subpackage('io')
config.add_subpackage('morphology')
def add_test_directories(arg, dirname, fnames):
if dirname.split(os.path.sep)[-1] == 'tests':
+2 -2
View File
@@ -1,2 +1,2 @@
version='unbuilt-dev'
# THIS FILE IS GENERATED FROM THE SCIKITS.IMAGE SETUP.PY
version='0.2dev'