diff --git a/CONTRIBUTORS.txt b/CONTRIBUTORS.txt
index 307379f8..bff86ccc 100644
--- a/CONTRIBUTORS.txt
+++ b/CONTRIBUTORS.txt
@@ -84,3 +84,5 @@
- Nelle Varoquaux
Renaming of the package to ``skimage``.
+- W. Randolph Franklin
+ Point in polygon test.
diff --git a/doc/examples/plot_convex_hull.py b/doc/examples/plot_convex_hull.py
new file mode 100644
index 00000000..8308bd90
--- /dev/null
+++ b/doc/examples/plot_convex_hull.py
@@ -0,0 +1,35 @@
+"""
+===========
+Convex Hull
+===========
+
+The convex hull of a binary image is the set of pixels included in the
+smallest convex polygon that surround all white pixels in the input.
+
+In this example, we show how the input pixels (white) get filled in by the
+convex hull (white and grey).
+
+A good overview of the algorithm is given on `Steve Eddin's blog
+`__.
+
+"""
+
+import numpy as np
+import matplotlib.pyplot as plt
+
+from skimage.morphology import convex_hull
+
+image = np.array(
+ [[0, 0, 0, 0, 0, 0, 0, 0, 0],
+ [0, 0, 0, 0, 1, 0, 0, 0, 0],
+ [0, 0, 0, 1, 0, 1, 0, 0, 0],
+ [0, 0, 1, 0, 0, 0, 1, 0, 0],
+ [0, 1, 0, 0, 0, 0, 0, 1, 0],
+ [0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=float)
+
+chull = convex_hull(image)
+image[chull] += 1.7
+image -= -1.7
+
+plt.imshow(image, cmap=plt.cm.gray, interpolation='nearest')
+plt.show()
diff --git a/skimage/morphology/__init__.py b/skimage/morphology/__init__.py
index efb92507..e3176579 100644
--- a/skimage/morphology/__init__.py
+++ b/skimage/morphology/__init__.py
@@ -3,3 +3,4 @@ from selem import *
from .ccomp import label
from watershed import watershed, is_local_maximum
from skeletonize import skeletonize, medial_axis
+from .convex_hull import convex_hull_image
diff --git a/skimage/morphology/_convex_hull.pyx b/skimage/morphology/_convex_hull.pyx
new file mode 100644
index 00000000..dc426900
--- /dev/null
+++ b/skimage/morphology/_convex_hull.pyx
@@ -0,0 +1,59 @@
+# -*- python -*-
+
+cimport numpy as np
+import numpy as np
+
+def possible_hull(np.ndarray[dtype=np.uint8_t, ndim=2, mode="c"] img):
+ """Return positions of pixels that possibly belong to the convex hull.
+
+ Parameters
+ ----------
+ img : ndarray of bool
+ Binary input image.
+
+ Returns
+ -------
+ coords : ndarray (N, 2)
+ The ``(row, column)`` coordinates of all pixels that possibly belong to
+ the convex hull.
+
+ """
+ cdef int i, j, k
+ cdef unsigned int M, N
+
+ M = img.shape[0]
+ N = img.shape[1]
+
+ # Output: M storage slots for left boundary pixels
+ # N storage slots for top boundary pixels
+ # M storage slots for right boundary pixels
+ # N storage slots for bottom boundary pixels
+ cdef np.ndarray[dtype=np.int_t, ndim=2] nonzero = \
+ np.ones((2 * (M + N), 2), dtype=np.int)
+ nonzero *= -1
+
+ k = 0
+ for i in range(M):
+ for j in range(N):
+ if img[i, j] != 0:
+ # Left check
+ if nonzero[i, 1] == -1:
+ nonzero[i, 0] = i
+ nonzero[i, 1] = j
+
+ # Right check
+ elif nonzero[M + N + i, 1] < j:
+ nonzero[M + N + i, 0] = i
+ nonzero[M + N + i, 1] = j
+
+ # Top check
+ if nonzero[M + j, 1] == -1:
+ nonzero[M + j, 0] = i
+ nonzero[M + j, 1] = j
+
+ # Bottom check
+ elif nonzero[2 * M + N + j, 0] < i:
+ nonzero[2 * M + N + j, 0] = i
+ nonzero[2 * M + N + j, 1] = j
+
+ return nonzero[nonzero[:, 0] != -1]
diff --git a/skimage/morphology/_pnpoly.h b/skimage/morphology/_pnpoly.h
new file mode 100644
index 00000000..95c89bcb
--- /dev/null
+++ b/skimage/morphology/_pnpoly.h
@@ -0,0 +1,72 @@
+/* `pnpoly` is from
+
+ http://www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/pnpoly.html
+
+ Copyright (c) 1970-2003, Wm. Randolph Franklin
+
+ Permission is hereby granted, free of charge, to any person
+ obtaining a copy of this software and associated documentation
+ files (the "Software"), to deal in the Software without
+ restriction, including without limitation the rights to use, copy,
+ modify, merge, publish, distribute, sublicense, and/or sell copies
+ of the Software, and to permit persons to whom the Software is
+ furnished to do so, subject to the following conditions:
+
+ 1. Redistributions of source code must retain the above copyright
+ notice, this list of conditions and the following disclaimers.
+ 2. Redistributions in binary form must reproduce the above
+ copyright notice in the documentation and/or other materials
+ provided with the distribution.
+ 3. The name of W. Randolph Franklin may not be used to endorse or
+ promote products derived from this Software without specific
+ prior written permission.
+
+ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
+ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
+ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
+ NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS
+ BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN
+ ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
+ CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+ SOFTWARE. */
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+unsigned char pnpoly(int nr_verts, double *xp, double *yp, double x, double y)
+{
+ int i, j;
+ unsigned char c = 0;
+ for (i = 0, j = nr_verts-1; i < nr_verts; j = i++) {
+ if ((((yp[i]<=y) && (yvx.data, vy.data, m, n)
+
+ return out.view(bool)
+
+
+def points_inside_poly(points, verts):
+ """Test whether points lie inside a polygon.
+
+ Parameters
+ ----------
+ points : (N, 2) array
+ Input points, ``(x, y)``.
+ verts : (M, 2) array
+ Vertices of the polygon, sorted either clockwise or anti-clockwise.
+ The first point may (but does not need to be) duplicated.
+
+ Returns
+ -------
+ mask : (N,) array of bool
+ True if corresponding point is inside the polygon.
+
+ """
+ cdef np.ndarray[np.double_t, ndim=1, mode="c"] x, y, vx, vy
+
+ points = np.asarray(points)
+ verts = np.asarray(verts)
+
+ x = points[:, 0].astype(np.double)
+ y = points[:, 1].astype(np.double)
+
+ vx = verts[:, 0].astype(np.double)
+ vy = verts[:, 1].astype(np.double)
+
+ cdef np.ndarray[np.uint8_t, ndim=1] out = \
+ np.zeros(x.shape[0], dtype=np.uint8)
+
+ npnpoly(vx.shape[0], vx.data, vy.data,
+ x.shape[0], x.data, y.data,
+ out.data)
+
+ return out.astype(bool)
+
diff --git a/skimage/morphology/convex_hull.py b/skimage/morphology/convex_hull.py
new file mode 100644
index 00000000..f461dc4d
--- /dev/null
+++ b/skimage/morphology/convex_hull.py
@@ -0,0 +1,65 @@
+__all__ = ['convex_hull_image']
+
+import numpy as np
+from ._pnpoly import points_inside_poly, grid_points_inside_poly
+from ._convex_hull import possible_hull
+
+def convex_hull_image(image):
+ """Compute the convex hull image of a binary image.
+
+ The convex hull is the set of pixels included in the smallest convex
+ polygon that surround all white pixels in the input image.
+
+ Parameters
+ ----------
+ image : ndarray
+ Binary input image. This array is cast to bool before processing.
+
+ Returns
+ -------
+ hull : ndarray of uint8
+ Binary image with pixels in convex hull set to 255.
+
+ References
+ ----------
+ .. [1] http://blogs.mathworks.com/steve/2011/10/04/binary-image-convex-hull-algorithm-notes/
+
+ """
+
+ image = image.astype(bool)
+
+ # Here we do an optimisation by choosing only pixels that are
+ # the starting or ending pixel of a row or column. This vastly
+ # limits the number of coordinates to examine for the virtual
+ # hull.
+ coords = possible_hull(image.astype(np.uint8))
+ N = len(coords)
+
+ # Add a vertex for the middle of each pixel edge
+ coords_corners = np.empty((N * 4, 2))
+ for i, (x_offset, y_offset) in enumerate(zip((0, 0, -0.5, 0.5),
+ (-0.5, 0.5, 0, 0))):
+ coords_corners[i * N:(i + 1) * N] = coords + [x_offset, y_offset]
+
+ coords = coords_corners
+
+ try:
+ from scipy.spatial import Delaunay
+ except ImportError:
+ raise ImportError('Could not import scipy.spatial, only available in '
+ 'scipy >= 0.9.')
+
+ # Find the convex hull
+ chull = Delaunay(coords).convex_hull
+ v = coords[np.unique(chull)]
+
+ # Sort vertices clock-wise
+ v_centred = v - v.mean(axis=0)
+ angles = np.arctan2(v_centred[:, 0], v_centred[:, 1])
+ v = v[np.argsort(angles)]
+
+ # For each pixel coordinate, check whether that pixel
+ # lies inside the convex hull
+ mask = grid_points_inside_poly(image.shape[:2], v)
+
+ return mask
diff --git a/skimage/morphology/setup.py b/skimage/morphology/setup.py
index 46010989..e28446bd 100644
--- a/skimage/morphology/setup.py
+++ b/skimage/morphology/setup.py
@@ -15,6 +15,8 @@ def configuration(parent_package='', top_path=None):
cython(['cmorph.pyx'], working_path=base_path)
cython(['_watershed.pyx'], working_path=base_path)
cython(['_skeletonize.pyx'], working_path=base_path)
+ cython(['_pnpoly.pyx'], working_path=base_path)
+ cython(['_convex_hull.pyx'], working_path=base_path)
config.add_extension('ccomp', sources=['ccomp.c'],
include_dirs=[get_numpy_include_dirs()])
@@ -24,7 +26,10 @@ def configuration(parent_package='', top_path=None):
include_dirs=[get_numpy_include_dirs()])
config.add_extension('_skeletonize', sources=['_skeletonize.c'],
include_dirs=[get_numpy_include_dirs()])
-
+ config.add_extension('_pnpoly', sources=['_pnpoly.c'],
+ include_dirs=[get_numpy_include_dirs()])
+ config.add_extension('_convex_hull', sources=['_convex_hull.c'],
+ include_dirs=[get_numpy_include_dirs()])
return config
diff --git a/skimage/morphology/tests/test_convex_hull.py b/skimage/morphology/tests/test_convex_hull.py
new file mode 100644
index 00000000..21f45cd7
--- /dev/null
+++ b/skimage/morphology/tests/test_convex_hull.py
@@ -0,0 +1,67 @@
+import numpy as np
+from numpy.testing import assert_array_equal
+from numpy.testing.decorators import skipif
+from skimage.morphology import convex_hull_image
+from skimage.morphology._convex_hull import possible_hull
+
+try:
+ import scipy.spatial
+ scipy_spatial = True
+except ImportError:
+ scipy_spatial = False
+
+@skipif(not scipy_spatial)
+def test_basic():
+ image = np.array(
+ [[0, 0, 0, 0, 0, 0, 0, 0, 0],
+ [0, 0, 0, 0, 1, 0, 0, 0, 0],
+ [0, 0, 0, 1, 0, 1, 0, 0, 0],
+ [0, 0, 1, 0, 0, 0, 1, 0, 0],
+ [0, 1, 0, 0, 0, 0, 0, 1, 0],
+ [0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=bool)
+
+ expected = np.array(
+ [[0, 0, 0, 0, 0, 0, 0, 0, 0],
+ [0, 0, 0, 0, 1, 0, 0, 0, 0],
+ [0, 0, 0, 1, 1, 1, 0, 0, 0],
+ [0, 0, 1, 1, 1, 1, 1, 0, 0],
+ [0, 1, 1, 1, 1, 1, 1, 1, 0],
+ [0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=bool)
+
+ assert_array_equal(convex_hull_image(image), expected)
+
+@skipif(not scipy_spatial)
+def test_possible_hull():
+ image = np.array(
+ [[0, 0, 0, 0, 0, 0, 0, 0, 0],
+ [0, 0, 0, 0, 1, 0, 0, 0, 0],
+ [0, 0, 0, 1, 0, 1, 0, 0, 0],
+ [0, 0, 1, 1, 1, 1, 1, 0, 0],
+ [0, 1, 1, 1, 1, 1, 1, 1, 0],
+ [0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=np.uint8)
+
+ expected = np.array([[1, 4],
+ [2, 3],
+ [3, 2],
+ [4, 1],
+ [4, 1],
+ [3, 2],
+ [2, 3],
+ [1, 4],
+ [2, 5],
+ [3, 6],
+ [4, 7],
+ [2, 5],
+ [3, 6],
+ [4, 7],
+ [4, 2],
+ [4, 3],
+ [4, 4],
+ [4, 5],
+ [4, 6]])
+
+ ph = possible_hull(image)
+ assert_array_equal(ph, expected)
+
+if __name__ == "__main__":
+ np.testing.run_module_suite()
diff --git a/skimage/morphology/tests/test_pnpoly.py b/skimage/morphology/tests/test_pnpoly.py
new file mode 100644
index 00000000..33efe15f
--- /dev/null
+++ b/skimage/morphology/tests/test_pnpoly.py
@@ -0,0 +1,38 @@
+import numpy as np
+from numpy.testing import assert_array_equal
+
+from skimage.morphology._pnpoly import points_inside_poly, \
+ grid_points_inside_poly
+
+class test_npnpoly():
+ def test_square(self):
+ v = np.array([[0, 0],
+ [0, 1],
+ [1, 1],
+ [1, 0]])
+ assert(points_inside_poly([[0.5, 0.5]], v)[0])
+ assert(not points_inside_poly([[-0.1, 0.1]], v)[0])
+
+ def test_triangle(self):
+ v = np.array([[0, 0],
+ [1, 0],
+ [0.5, 0.75]])
+ assert(points_inside_poly([[0.5, 0.7]], v)[0])
+ assert(not points_inside_poly([[0.5, 0.76]], v)[0])
+ assert(not points_inside_poly([[0.7, 0.5]], v)[0])
+
+ def test_type(self):
+ assert(points_inside_poly([[0, 0]], [[0, 0]]).dtype == np.bool)
+
+def test_grid_points_inside_poly():
+ v = np.array([[0, 0],
+ [5, 0],
+ [5, 5]])
+
+ expected = np.tril(np.ones((5, 5), dtype=bool))
+
+ assert_array_equal(grid_points_inside_poly((5, 5), v),
+ expected)
+
+if __name__ == "__main__":
+ np.testing.run_module_suite()