add 2d convex hull utility function, BLI_convexhull_2d(), and python api mathutils...
authorCampbell Barton <ideasman42@gmail.com>
Tue, 10 Sep 2013 07:52:10 +0000 (07:52 +0000)
committerCampbell Barton <ideasman42@gmail.com>
Tue, 10 Sep 2013 07:52:10 +0000 (07:52 +0000)
uses Andrew's monotone chain 2D convex hull algorithm.

source/blender/blenlib/BLI_convexhull2d.h [new file with mode: 0644]
source/blender/blenlib/CMakeLists.txt
source/blender/blenlib/intern/convexhull2d.c [new file with mode: 0644]
source/blender/python/mathutils/mathutils_geometry.c

diff --git a/source/blender/blenlib/BLI_convexhull2d.h b/source/blender/blenlib/BLI_convexhull2d.h
new file mode 100644 (file)
index 0000000..30de786
--- /dev/null
@@ -0,0 +1,31 @@
+/*
+ * ***** BEGIN GPL LICENSE BLOCK *****
+ *
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ *
+ * ***** END GPL LICENSE BLOCK *****
+ */
+
+#ifndef __BLI_CONVEXHULL2D__
+#define __BLI_CONVEXHULL2D__
+
+/** \file BLI_convexhull2d.h
+ *  \ingroup bli
+ */
+
+int BLI_convexhull_2d_presorted(const float (*points)[2], const int n, int r_points[]);
+int BLI_convexhull_2d(const float (*points)[2], const int n, int r_points[]);
+
+#endif  /* __BLI_CONVEXHULL2D__ */
index 82f07268521f650c384c55dfb498f1a5969ad5b2..31125fd315326c0581a1728a68fcfe3966ef5f45 100644 (file)
@@ -54,6 +54,7 @@ set(SRC
        intern/boxpack2d.c
        intern/buffer.c
        intern/callbacks.c
+       intern/convexhull2d.c
        intern/cpu.c
        intern/dynlib.c
        intern/edgehash.c
@@ -109,6 +110,7 @@ set(SRC
        BLI_buffer.h
        BLI_callbacks.h
        BLI_compiler_attrs.h
+       BLI_convexhull2d.h
        BLI_cpu.h
        BLI_dlrbTree.h
        BLI_dynlib.h
diff --git a/source/blender/blenlib/intern/convexhull2d.c b/source/blender/blenlib/intern/convexhull2d.c
new file mode 100644 (file)
index 0000000..9909f47
--- /dev/null
@@ -0,0 +1,221 @@
+/*
+ * ***** BEGIN GPL LICENSE BLOCK *****
+ *
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ *
+ * ***** END GPL LICENSE BLOCK *****
+ */
+
+#include <stdlib.h>
+#include <string.h>
+
+#include "MEM_guardedalloc.h"
+
+#include "BLI_convexhull2d.h"
+
+#include "BLI_strict_flags.h"
+
+
+/* Copyright 2001, softSurfer (www.softsurfer.com)
+ * This code may be freely used and modified for any purpose
+ * providing that this copyright notice is included with it.
+ * SoftSurfer makes no warranty for this code, and cannot be held
+ * liable for any real or imagined damage resulting from its use.
+ * Users of this code must verify correctness for their application.
+ * http://softsurfer.com/Archive/algorithm_0203/algorithm_0203.htm
+ */
+
+/**
+ * tests if a point is Left|On|Right of an infinite line.
+ *    Input:  three points P0, P1, and P2
+ * \returns > 0.0 for P2 left of the line through P0 and P1.
+ *          = 0.0 for P2 on the line.
+ *          < 0.0 for P2 right of the line.
+ */
+static float is_left(const float p0[2], const float p1[2], const float p2[2])
+{
+       return (p1[0] - p0[0]) * (p2[1] - p0[1]) - (p2[0] - p0[0]) * (p1[1] - p0[1]);
+}
+
+/**
+ * A.M. Andrew's monotone chain 2D convex hull algorithm
+ *
+ * \param  points  An array of 2D points presorted by increasing x and y-coords.
+ * \param  n  The number of points in points.
+ * \param  r_points  An array of the convex hull vertex indices (max is n).
+ * \returns the number of points in r_points.
+ */
+int BLI_convexhull_2d_presorted(const float (*points)[2], const int n, int r_points[])
+{
+       /* the output array r_points[] will be used as the stack */
+       int bot = 0;
+       int top = -1;  /* indices for bottom and top of the stack */
+       int i;  /* array scan index */
+       int minmin, minmax;
+       int maxmin, maxmax;
+       float xmax;
+
+       /* Get the indices of points with min x-coord and min|max y-coord */
+       float xmin = points[0][0];
+       for (i = 1; i < n; i++) {
+               if (points[i][0] != xmin) {
+                       break;
+               }
+       }
+
+       minmin = 0;
+       minmax = i - 1;
+       if (minmax == n - 1) {  /* degenerate case: all x-coords == xmin */
+               r_points[++top] = minmin;
+               if (points[minmax][1] != points[minmin][1])  /* a nontrivial segment */
+                       r_points[++top] = minmax;
+               r_points[++top] = minmin;  /* add polygon endpoint */
+               return top + 1;
+       }
+
+       /* Get the indices of points with max x-coord and min|max y-coord */
+
+       maxmax = n - 1;
+       xmax = points[n - 1][0];
+       for (i = n - 2; i >= 0; i--) {
+               if (points[i][0] != xmax) {
+                       break;
+               }
+       }
+       maxmin = i + 1;
+
+       /* Compute the lower hull on the stack r_points */
+       r_points[++top] = minmin;  /* push minmin point onto stack */
+       i = minmax;
+       while (++i <= maxmin) {
+               /* the lower line joins points[minmin] with points[maxmin] */
+               if (is_left(points[minmin], points[maxmin], points[i]) >= 0 && i < maxmin) {
+                       continue;  /* ignore points[i] above or on the lower line */
+               }
+
+               while (top > 0) {  /* there are at least 2 points on the stack */
+                       /* test if points[i] is left of the line at the stack top */
+                       if (is_left(points[r_points[top - 1]], points[r_points[top]], points[i]) > 0.0f) {
+                               break;  /* points[i] is a new hull vertex */
+                       }
+                       else {
+                               top--;  /* pop top point off stack */
+                       }
+               }
+
+               r_points[++top] = i;  /* push points[i] onto stack */
+       }
+
+       /* Next, compute the upper hull on the stack r_points above the bottom hull */
+       if (maxmax != maxmin) {  /* if distinct xmax points */
+               r_points[++top] = maxmax;  /* push maxmax point onto stack */
+       }
+
+       bot = top;  /* the bottom point of the upper hull stack */
+       i = maxmin;
+       while (--i >= minmax) {
+               /* the upper line joins points[maxmax] with points[minmax] */
+               if (is_left(points[maxmax], points[minmax], points[i]) >= 0 && i > minmax) {
+                       continue;  /* ignore points[i] below or on the upper line */
+               }
+
+               while (top > bot) {  /* at least 2 points on the upper stack */
+                       /* test if points[i] is left of the line at the stack top */
+                       if (is_left(points[r_points[top - 1]], points[r_points[top]], points[i]) > 0.0f) {
+                               break;  /* points[i] is a new hull vertex */
+                       }
+                       else {
+                               top--;  /* pop top point off stack */
+                       }
+               }
+
+               if (points[i][0] == points[r_points[0]][0] && points[i][1] == points[r_points[0]][1]) {
+                       return top + 1;  /* special case (mgomes) */
+               }
+
+               r_points[++top] = i;  /* push points[i] onto stack */
+       }
+
+       if (minmax != minmin) {
+               r_points[++top] = minmin;  /* push joining endpoint onto stack */
+       }
+
+       return top + 1;
+}
+
+struct PointRef {
+       const float *pt;  /* 2d vector */
+};
+
+static int pointref_cmp_x(const void *a_, const void *b_)
+{
+       const struct PointRef *a = a_;
+       const struct PointRef *b = b_;
+       if      (a->pt[0] > b->pt[0]) return  1;
+       else if (a->pt[0] < b->pt[0]) return -1;
+       else                          return  0;
+}
+
+static int pointref_cmp_y(const void *a_, const void *b_)
+{
+       const struct PointRef *a = a_;
+       const struct PointRef *b = b_;
+       if      (a->pt[1] > b->pt[1]) return  1;
+       else if (a->pt[1] < b->pt[1]) return -1;
+       else                          return  0;
+}
+
+/**
+ * A.M. Andrew's monotone chain 2D convex hull algorithm
+ *
+ * \param  points  An array of 2D points.
+ * \param  n  The number of points in points.
+ * \param  r_points  An array of the convex hull vertex indices (max is n).
+ * \returns the number of points in r_points.
+ */
+int BLI_convexhull_2d(const float (*points)[2], const int n, int r_points[])
+{
+       struct PointRef *points_ref = MEM_mallocN(sizeof(*points_ref) * (size_t)n, __func__);
+       float (*points_sort)[2] = MEM_mallocN(sizeof(*points_sort) * (size_t)n, __func__);
+       int *points_map;
+       int tot, i;
+
+       for (i = 0; i < n; i++) {
+               points_ref[i].pt = points[i];
+       }
+
+       /* Sort the points by X, then by Y (required by the algorithm) */
+       qsort(points_ref, (size_t)n, sizeof(struct PointRef), pointref_cmp_x);
+       qsort(points_ref, (size_t)n, sizeof(struct PointRef), pointref_cmp_y);
+
+       for (i = 0; i < n; i++) {
+               memcpy(points_sort[i], points_ref[i].pt, sizeof(float[2]));
+       }
+
+       tot = BLI_convexhull_2d_presorted((const float (*)[2])points_sort, n, r_points);
+
+       /* map back to the original index values */
+       points_map = (int *)points_sort;  /* abuse float array for temp storage */
+       for (i = 0; i < tot; i++) {
+               points_map[i] = (int)((const float(*)[2])points_ref[r_points[i]].pt - points);
+       }
+
+       memcpy(r_points, points_map, (size_t)tot * sizeof(*points_map));
+
+       MEM_freeN(points_ref);
+       MEM_freeN(points_sort);
+
+       return tot;
+}
index b69e3eb23f8d36d6af7bd5251ed95702002f237a..0113efb8dc27c152822a940460ff3fe9ae672e32 100644 (file)
@@ -34,6 +34,7 @@
 #  include "MEM_guardedalloc.h"
 #  include "BLI_blenlib.h"
 #  include "BLI_boxpack2d.h"
+#  include "BLI_convexhull2d.h"
 #  include "BKE_displist.h"
 #  include "BKE_curve.h"
 #endif
@@ -1503,6 +1504,70 @@ static PyObject *M_Geometry_box_pack_2d(PyObject *UNUSED(self), PyObject *boxlis
        return ret;
 }
 
+
+PyDoc_STRVAR(M_Geometry_convex_hull_2d_doc,
+".. function:: convex_hull_2d(points)\n"
+"\n"
+"   Returns a list of indices into the list given\n"
+"\n"
+"   :arg points: list of 2d points.\n"
+"   :type points: list\n"
+"   :return: a list of indices\n"
+"   :rtype: list of ints\n"
+);
+static PyObject *M_Geometry_convex_hull_2d(PyObject *UNUSED(self), PyObject *pointlist)
+{
+       Py_ssize_t len;
+
+       PyObject *ret;
+
+       if (!PyList_Check(pointlist)) {
+               PyErr_SetString(PyExc_TypeError,
+                               "expected a list of Vectors");
+               return NULL;
+       }
+
+       len = PyList_GET_SIZE(pointlist);
+       if (len) {
+               float (*points)[2] = MEM_mallocN(sizeof(*points) * len, __func__);
+               int *index_map;
+               Py_ssize_t len_ret, i;
+               PyObject *list_item;
+               bool ok = true;
+
+               for (i = 0; i < len; i++) {
+                       list_item = PyList_GET_ITEM(pointlist, i);
+                       if (mathutils_array_parse(points[i], 2, 2, list_item, "convex_hull") == -1) {
+                               ok = false;
+                               break;
+                       }
+               }
+
+               if (ok == false) {
+                       MEM_freeN(points);
+                       return NULL;
+               }
+
+               index_map  = MEM_mallocN(sizeof(*index_map) * len, __func__);
+
+               /* Non Python function */
+               len_ret = BLI_convexhull_2d((const float (*)[2])points, len, index_map);
+               MEM_freeN(points);
+
+               ret = PyList_New(len_ret);
+               for (i = 0; i < len_ret; i++) {
+                       PyList_SET_ITEM(ret, i, PyLong_FromLong(index_map[i]));
+               }
+
+               MEM_freeN(index_map);
+       }
+       else {
+               ret = PyList_New(0);
+       }
+
+       return ret;
+}
+
 #endif /* MATH_STANDALONE */
 
 
@@ -1527,6 +1592,7 @@ static PyMethodDef M_Geometry_methods[] = {
 #ifndef MATH_STANDALONE
        {"interpolate_bezier", (PyCFunction) M_Geometry_interpolate_bezier, METH_VARARGS, M_Geometry_interpolate_bezier_doc},
        {"tessellate_polygon", (PyCFunction) M_Geometry_tessellate_polygon, METH_O, M_Geometry_tessellate_polygon_doc},
+       {"convex_hull_2d", (PyCFunction) M_Geometry_convex_hull_2d, METH_O, M_Geometry_convex_hull_2d_doc},
        {"box_pack_2d", (PyCFunction) M_Geometry_box_pack_2d, METH_O, M_Geometry_box_pack_2d_doc},
 #endif
        {NULL, NULL, 0, NULL}