From c471d475eb5f62b930bf10ec8e5da08682a33751 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Johannes=20Sch=C3=B6nberger?= Date: Tue, 21 Aug 2012 15:43:26 +0200 Subject: [PATCH] Refactor pnpoly function in Cython and add to shared package --- skimage/_shared/geometry.pxd | 7 ++++ skimage/_shared/geometry.pyx | 27 +++++++++++++ skimage/_shared/setup.py | 2 + skimage/draw/_draw.pyx | 8 +--- skimage/draw/setup.py | 2 +- skimage/morphology/_pnpoly.h | 72 ---------------------------------- skimage/morphology/_pnpoly.pyx | 17 +++----- skimage/morphology/setup.py | 2 +- 8 files changed, 45 insertions(+), 92 deletions(-) create mode 100644 skimage/_shared/geometry.pxd create mode 100644 skimage/_shared/geometry.pyx delete mode 100644 skimage/morphology/_pnpoly.h diff --git a/skimage/_shared/geometry.pxd b/skimage/_shared/geometry.pxd new file mode 100644 index 00000000..b228b23a --- /dev/null +++ b/skimage/_shared/geometry.pxd @@ -0,0 +1,7 @@ + +cdef inline unsigned char point_in_polygon(int nr_verts, double *xp, double *yp, + double x, double y) + +cdef void points_in_polygon(int nr_verts, double *xp, double *yp, + int nr_points, double *x, double *y, + unsigned char *result) diff --git a/skimage/_shared/geometry.pyx b/skimage/_shared/geometry.pyx new file mode 100644 index 00000000..de1ed4f4 --- /dev/null +++ b/skimage/_shared/geometry.pyx @@ -0,0 +1,27 @@ +#cython: cdivison=True +#cython: boundscheck=False +#cython: nonecheck=False +#cython: wraparound=False + + +cdef inline unsigned char point_in_polygon(int nr_verts, double *xp, double *yp, + double x, double y): + cdef int i + cdef unsigned char c = 0 + cdef int j = nr_verts - 1 + for i in range(nr_verts): + if ( + (((yp[i] <= y) and (y < yp[j])) or + ((yp[j] <= y) and (y < yp[i]))) + and (x < (xp[j] - xp[i]) * (y - yp[i]) / (yp[j] - yp[i]) + xp[i]) + ): + c = not c + j = i + return c + +cdef void points_in_polygon(int nr_verts, double *xp, double *yp, + int nr_points, double *x, double *y, + unsigned char *result): + cdef int n + for n in range(nr_points): + result[n] = point_in_polygon(nr_verts, xp, yp, x[n], y[n]) diff --git a/skimage/_shared/setup.py b/skimage/_shared/setup.py index ed48916e..81113656 100644 --- a/skimage/_shared/setup.py +++ b/skimage/_shared/setup.py @@ -13,9 +13,11 @@ def configuration(parent_package='', top_path=None): config = Configuration('_shared', parent_package, top_path) config.add_data_dir('tests') + cython(['geometry.pyx'], working_path=base_path) cython(['interpolation.pyx'], working_path=base_path) cython(['transform.pyx'], working_path=base_path) + config.add_extension('geometry', sources=['geometry.c']) config.add_extension('interpolation', sources=['interpolation.c'], include_dirs=[get_numpy_include_dirs()]) config.add_extension('transform', sources=['transform.c'], diff --git a/skimage/draw/_draw.pyx b/skimage/draw/_draw.pyx index 9dc1d307..fbd525fa 100644 --- a/skimage/draw/_draw.pyx +++ b/skimage/draw/_draw.pyx @@ -3,11 +3,7 @@ import math from libc.math cimport sqrt cimport numpy as np cimport cython - - -cdef extern from "../morphology/_pnpoly.h": - int pnpoly(int nr_verts, double *xp, double *yp, - double x, double y) +from skimage._shared.geometry cimport point_in_polygon @cython.boundscheck(False) @@ -119,7 +115,7 @@ def polygon(y, x, shape=None): for r in range(minr, maxr+1): for c in range(minc, maxc+1): - if pnpoly(nr_verts, cptr, rptr, c, r): + if point_in_polygon(nr_verts, cptr, rptr, c, r): rr.append(r) cc.append(c) diff --git a/skimage/draw/setup.py b/skimage/draw/setup.py index 4503d664..59067f9a 100644 --- a/skimage/draw/setup.py +++ b/skimage/draw/setup.py @@ -15,7 +15,7 @@ def configuration(parent_package='', top_path=None): cython(['_draw.pyx'], working_path=base_path) config.add_extension('_draw', sources=['_draw.c'], - include_dirs=[get_numpy_include_dirs()]) + include_dirs=[get_numpy_include_dirs(), '../shared']) return config diff --git a/skimage/morphology/_pnpoly.h b/skimage/morphology/_pnpoly.h deleted file mode 100644 index 95c89bcb..00000000 --- a/skimage/morphology/_pnpoly.h +++ /dev/null @@ -1,72 +0,0 @@ -/* `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) + out[m, n] = point_in_polygon(V, vx.data, vy.data, m, n) return out.view(bool) - + def points_inside_poly(points, verts): """Test whether points lie inside a polygon. @@ -84,8 +77,8 @@ def points_inside_poly(points, verts): 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, + + points_in_polygon(vx.shape[0], vx.data, vy.data, x.shape[0], x.data, y.data, out.data) diff --git a/skimage/morphology/setup.py b/skimage/morphology/setup.py index 2dd4a378..bbd176e4 100644 --- a/skimage/morphology/setup.py +++ b/skimage/morphology/setup.py @@ -28,7 +28,7 @@ def configuration(parent_package='', top_path=None): config.add_extension('_skeletonize_cy', sources=['_skeletonize_cy.c'], include_dirs=[get_numpy_include_dirs()]) config.add_extension('_pnpoly', sources=['_pnpoly.c'], - include_dirs=[get_numpy_include_dirs()]) + include_dirs=[get_numpy_include_dirs(), '../shared']) config.add_extension('_convex_hull', sources=['_convex_hull.c'], include_dirs=[get_numpy_include_dirs()])