From 88a269c282aa44e7ada8fa611801eb153c3b0005 Mon Sep 17 00:00:00 2001 From: Stefan van der Walt Date: Mon, 29 Sep 2014 01:14:57 +0200 Subject: [PATCH] Add polygon perimiter drawing --- skimage/_shared/_geometry.py | 54 +++++++++++++++++ skimage/_shared/tests/test_geometry.py | 75 ++++++++++++++++++++++++ skimage/draw/__init__.py | 2 +- skimage/draw/draw.py | 81 ++++++++++++++++++++++++++ skimage/draw/tests/test_draw.py | 33 ++++++++++- 5 files changed, 243 insertions(+), 2 deletions(-) create mode 100644 skimage/_shared/_geometry.py create mode 100644 skimage/_shared/tests/test_geometry.py diff --git a/skimage/_shared/_geometry.py b/skimage/_shared/_geometry.py new file mode 100644 index 00000000..602692cb --- /dev/null +++ b/skimage/_shared/_geometry.py @@ -0,0 +1,54 @@ +__all__ = ['polygon_clip', 'polygon_area'] + +import numpy as np + + +def polygon_clip(rr, cc, r0, c0, r1, c1): + """Clip a polygon to the given bounding box. + + Parameters + ---------- + yp, xp : (N,) ndarray of double + Row and column coordinates of the polygon. + ytop, xleft, ybottom, xright : double + Coordinates of the bounding box. Note that the following + must hold true: ``x_left < x_right`` and ``y_top < y_bottom``. + + Returns + ------- + y_clipped, x_clipped : (M,) ndarray of double + Coordinates of clipped polygon. + + Notes + ----- + The algorithm is a translation of the Pascal code found in [1]_ and + includes fixes from Anti-Grain Geometry v2.4. + + References + ---------- + .. [1] You-Dong Liang and Brian A. Barsky, + An Analysis and Algorithm for Polygon Clipping, + Communications of the ACM, Vol 26, No 11, November 1983. + + """ + from matplotlib import path, transforms + + poly = path.Path(np.vstack((rr, cc)).T, closed=True) + clip_rect = transforms.Bbox([[r0, c0], [r1, c1]]) + + poly_clipped = poly.clip_to_bbox(clip_rect).to_polygons()[0] + + return poly_clipped[:, 0], poly_clipped[:, 1] + + +def polygon_area(py, px): + """Compute the area of a polygon. + + Parameters + ---------- + py, px : (N,) array of float + Polygon coordinates. + """ + py = np.asarray(py) + px = np.asarray(px) + return 0.5 * np.abs(np.sum((px[:-1] * py[1:]) - (px[1:] * py[:-1]))) diff --git a/skimage/_shared/tests/test_geometry.py b/skimage/_shared/tests/test_geometry.py new file mode 100644 index 00000000..7466d653 --- /dev/null +++ b/skimage/_shared/tests/test_geometry.py @@ -0,0 +1,75 @@ +from skimage._shared._geometry import polygon_clip, polygon_area + +import numpy as np +from numpy.testing import assert_equal, assert_almost_equal + + +hand = np.array( + [[ 1.64516129, 1.16145833 ], + [ 1.64516129, 1.59375 ], + [ 1.35080645, 1.921875 ], + [ 1.375 , 2.18229167 ], + [ 1.68548387, 1.9375 ], + [ 1.60887097, 2.55208333 ], + [ 1.68548387, 2.69791667 ], + [ 1.76209677, 2.56770833 ], + [ 1.83064516, 1.97395833 ], + [ 1.89516129, 2.75 ], + [ 1.9516129 , 2.84895833 ], + [ 2.01209677, 2.76041667 ], + [ 1.99193548, 1.99479167 ], + [ 2.11290323, 2.63020833 ], + [ 2.2016129 , 2.734375 ], + [ 2.25403226, 2.60416667 ], + [ 2.14919355, 1.953125 ], + [ 2.30645161, 2.36979167 ], + [ 2.39112903, 2.36979167 ], + [ 2.41532258, 2.1875 ], + [ 2.1733871 , 1.703125 ], + [ 2.07782258, 1.16666667 ]]) + + +def test_polygon_area(): + x = [0, 0, 1, 1] + y = [0, 1, 1, 0] + + assert_almost_equal(polygon_area(y, x), 1) + + x = [0, 0, 1] + y = [0, 1, 1] + + assert_almost_equal(polygon_area(y, x), 0.5) + + x = [0, 0, 0.5, 1, 1, 0.5] + y = [0, 1, 0.5, 1, 0, 0.5] + + assert_almost_equal(polygon_area(y, x), 0.5) + + +def test_poly_clip(): + x = [0, 1, 2, 1] + y = [0, -1, 0, 1] + + yc, xc = polygon_clip(y, x, 0, 0, 1, 1) + assert_equal(polygon_area(yc, xc), 0.5) + + x = [-1, 1.5, 1.5, -1] + y = [.5, 0.5, 1.5, 1.5] + yc, xc = polygon_clip(y, x, 0, 0, 1, 1) + assert_equal(polygon_area(yc, xc), 0.5) + + +def test_hand_clip(): + (r0, c0, r1, c1) = (1.0, 1.5, 2.1, 2.5) + clip_r, clip_c = polygon_clip(hand[:, 1], hand[:, 0], r0, c0, r1, c1) + assert_equal(clip_r.size, 19) + assert_equal(clip_r[0], clip_r[-1]) + assert_equal(clip_c[0], clip_c[-1]) + + (r0, c0, r1, c1) = (1.0, 1.5, 1.7, 2.5) + clip_r, clip_c = polygon_clip(hand[:, 1], hand[:, 0], r0, c0, r1, c1) + assert_equal(clip_r.size, 6) + + (r0, c0, r1, c1) = (1.0, 1.5, 1.5, 2.5) + clip_r, clip_c = polygon_clip(hand[:, 1], hand[:, 0], r0, c0, r1, c1) + assert_equal(clip_r.size, 5) diff --git a/skimage/draw/__init__.py b/skimage/draw/__init__.py index 5f788ea7..01dfc5a7 100644 --- a/skimage/draw/__init__.py +++ b/skimage/draw/__init__.py @@ -1,4 +1,4 @@ -from .draw import circle, ellipse, set_color +from .draw import circle, ellipse, polygon_perimiter, set_color from .draw3d import ellipsoid, ellipsoid_stats from ._draw import (line, line_aa, polygon, ellipse_perimeter, circle_perimeter, circle_perimeter_aa, diff --git a/skimage/draw/draw.py b/skimage/draw/draw.py index 1624b115..f9ee214b 100644 --- a/skimage/draw/draw.py +++ b/skimage/draw/draw.py @@ -1,7 +1,11 @@ # coding: utf-8 import numpy as np + from ._draw import _coords_inside_image +from .._shared._geometry import polygon_clip +from ._draw import line + def _ellipse_in_shape(shape, center, radiuses): """Generate coordinates of points within ellipse bounded by shape.""" @@ -125,6 +129,83 @@ def circle(r, c, radius, shape=None): return ellipse(r, c, radius, radius, shape) +def polygon_perimiter(cy, cx, shape=None, clip=False): + """Generate polygon perimiter coordinates. + + Parameters + ---------- + y : (N,) ndarray + Y-coordinates of vertices of polygon. + x : (N,) ndarray + X-coordinates of vertices of polygon. + shape : tuple, optional + Image shape which is used to determine maximum extents of output pixel + coordinates. This is useful for polygons which exceed the image size. + By default the full extents of the polygon are used. + clip : bool, optional + Whether to clip the polygon to the provided shape. If this is set + to True, the drawn figure will always be a closed polygon with all + edges visible. + + Returns + ------- + rr, cc : ndarray of int + Pixel coordinates of polygon. + May be used to directly index into an array, e.g. + ``img[rr, cc] = 1``. + + Examples + -------- + >>> from skimage.draw import polygon_perimiter + >>> img = np.zeros((10, 10), dtype=np.uint8) + >>> rr, cc = polygon_perimiter([5, -1, 5, 10], + ... [-1, 5, 11, 5], + ... shape=img.shape, clip=True) + >>> img[rr, cc] = 1 + >>> img + array([[0, 0, 0, 0, 1, 1, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0, 0, 1, 0, 0], + [0, 0, 1, 0, 0, 0, 0, 0, 1, 0], + [0, 1, 0, 0, 0, 0, 0, 0, 0, 1], + [1, 0, 0, 0, 0, 0, 0, 0, 0, 1], + [1, 0, 0, 0, 0, 0, 0, 0, 0, 1], + [1, 0, 0, 0, 0, 0, 0, 0, 0, 1], + [0, 1, 1, 0, 0, 0, 0, 0, 0, 1], + [0, 0, 0, 1, 0, 0, 0, 1, 1, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0]], dtype=uint8) + + """ + if clip: + if shape is None: + raise ValueError("Must specify clipping shape") + clip_box = np.array([0, 0, shape[0] - 1, shape[1] - 1]) + else: + clip_box = np.array([np.min(cy), np.min(cx), + np.max(cy), np.max(cx)]) + + # Do the clipping irrespective of whether clip is set. This + # ensures that the returned polygon is closed and is an array. + cy, cx = polygon_clip(cy, cx, *clip_box) + + cy = np.round(cy).astype(int) + cx = np.round(cx).astype(int) + + # Construct line segments + rr, cc = [], [] + for i in range(len(cy) - 1): + line_r, line_c = line(cy[i], cx[i], cy[i + 1], cx[i + 1]) + rr.extend(line_r) + cc.extend(line_c) + + rr = np.asarray(rr) + cc = np.asarray(cc) + + if shape is None: + return rr, cc + else: + return _coords_inside_image(rr, cc, shape) + + def set_color(img, coords, color, alpha=1): """Set pixel color in the image at the given coordinates. diff --git a/skimage/draw/tests/test_draw.py b/skimage/draw/tests/test_draw.py index 214f09b0..c1cbe4b8 100644 --- a/skimage/draw/tests/test_draw.py +++ b/skimage/draw/tests/test_draw.py @@ -2,7 +2,7 @@ from numpy.testing import assert_array_equal, assert_equal, assert_raises import numpy as np from skimage._shared.testing import test_parallel -from skimage.draw import (set_color, line, line_aa, polygon, +from skimage.draw import (set_color, line, line_aa, polygon, polygon_perimiter, circle, circle_perimeter, circle_perimeter_aa, ellipse, ellipse_perimeter, _bezier_segment, bezier_curve) @@ -809,6 +809,37 @@ def test_bezier_curve_shape(): assert_array_equal(img, img_[shift:-shift, :]) +def test_polygon_perimiter(): + expected = np.array( + [[1, 1, 1, 1], + [1, 0, 0, 1], + [1, 1, 1, 1]] + ) + out = np.zeros_like(expected) + + rr, cc = polygon_perimiter([0, 2, 2, 0], + [0, 0, 3, 3]) + + out[rr, cc] = 1 + assert_array_equal(out, expected) + + out = np.zeros_like(expected) + rr, cc = polygon_perimiter([-1, -1, 3, 3], + [-1, 4, 4, -1], + shape=out.shape, clip=True) + out[rr, cc] = 1 + assert_array_equal(out, expected) + + assert_raises(ValueError, polygon_perimiter, [0], [1], clip=True) + + +def test_polygon_perimiter_outside_image(): + rr, cc = polygon_perimiter([-1, -1, 3, 3], + [-1, 4, 4, -1], shape=(3, 4)) + assert_equal(len(rr), 0) + assert_equal(len(cc), 0) + + if __name__ == "__main__": from numpy.testing import run_module_suite run_module_suite()