[Bf-blender-cvs] [19a0df25e36] master: Cleanup: move plane array intersection into a function

Campbell Barton noreply at git.blender.org
Sat Nov 7 11:27:35 CET 2020


Commit: 19a0df25e36d9e75d486e8c51ccb874095a93e91
Author: Campbell Barton
Date:   Sat Nov 7 20:14:40 2020 +1100
Branches: master
https://developer.blender.org/rB19a0df25e36d9e75d486e8c51ccb874095a93e91

Cleanup: move plane array intersection into a function

Also add check to ensure a point isn't occluded by it's own plane,
which could happen if a small epsilon values are passed in.

===================================================================

M	source/blender/blenlib/BLI_math_geom.h
M	source/blender/blenlib/intern/math_geom.c
M	source/blender/python/mathutils/mathutils_geometry.c

===================================================================

diff --git a/source/blender/blenlib/BLI_math_geom.h b/source/blender/blenlib/BLI_math_geom.h
index a9a55b10a9e..c0a9ea91e75 100644
--- a/source/blender/blenlib/BLI_math_geom.h
+++ b/source/blender/blenlib/BLI_math_geom.h
@@ -358,6 +358,14 @@ bool isect_plane_plane_v3(const float plane_a[4],
                           float r_isect_co[3],
                           float r_isect_no[3]) ATTR_WARN_UNUSED_RESULT;
 
+bool isect_planes_v3_fn(
+    const float planes[][4],
+    const int planes_len,
+    const float eps_coplanar,
+    const float eps_isect,
+    void (*callback_fn)(const float co[3], int i, int j, int k, void *user_data),
+    void *user_data);
+
 /* line/ray triangle */
 bool isect_line_segment_tri_v3(const float p1[3],
                                const float p2[3],
diff --git a/source/blender/blenlib/intern/math_geom.c b/source/blender/blenlib/intern/math_geom.c
index 1d2480f4d62..2b0018e7662 100644
--- a/source/blender/blenlib/intern/math_geom.c
+++ b/source/blender/blenlib/intern/math_geom.c
@@ -2297,6 +2297,81 @@ bool isect_plane_plane_v3(const float plane_a[4],
   return false;
 }
 
+/**
+ * Intersect all planes, calling `callback_fn` for each point that intersects
+ * 3 of the planes that isn't outside any of the other planes.
+ *
+ * This can be thought of as calculating a convex-hull from an array of planes.
+ *
+ * \param eps_coplanar: Epsilon for testing if two planes are aligned (co-planar).
+ * \param eps_isect: Epsilon for testing of a point is behind any of the planes.
+ *
+ * \warning As complexity is a little under `O(N^3)`, this is only suitable for small arrays.
+ *
+ * \note This function could be optimized by some spatial structure.
+ */
+bool isect_planes_v3_fn(
+    const float planes[][4],
+    const int planes_len,
+    const float eps_coplanar,
+    const float eps_isect,
+    void (*callback_fn)(const float co[3], int i, int j, int k, void *user_data),
+    void *user_data)
+{
+  bool found = false;
+
+  float n1n2[3], n2n3[3], n3n1[3];
+
+  for (int i = 0; i < planes_len; i++) {
+    const float *n1 = planes[i];
+    for (int j = i + 1; j < planes_len; j++) {
+      const float *n2 = planes[j];
+      cross_v3_v3v3(n1n2, n1, n2);
+      if (len_squared_v3(n1n2) <= eps_coplanar) {
+        continue;
+      }
+      for (int k = j + 1; k < planes_len; k++) {
+        const float *n3 = planes[k];
+        cross_v3_v3v3(n2n3, n2, n3);
+        if (len_squared_v3(n2n3) <= eps_coplanar) {
+          continue;
+        }
+
+        cross_v3_v3v3(n3n1, n3, n1);
+        if (len_squared_v3(n3n1) <= eps_coplanar) {
+          continue;
+        }
+        const float quotient = -dot_v3v3(n1, n2n3);
+        if (fabsf(quotient) < eps_coplanar) {
+          continue;
+        }
+        const float co_test[3] = {
+            ((n2n3[0] * n1[3]) + (n3n1[0] * n2[3]) + (n1n2[0] * n3[3])) / quotient,
+            ((n2n3[1] * n1[3]) + (n3n1[1] * n2[3]) + (n1n2[1] * n3[3])) / quotient,
+            ((n2n3[2] * n1[3]) + (n3n1[2] * n2[3]) + (n1n2[2] * n3[3])) / quotient,
+        };
+        int i_test;
+        for (i_test = 0; i_test < planes_len; i_test++) {
+          const float *np_test = planes[i_test];
+          if (((dot_v3v3(np_test, co_test) + np_test[3]) > eps_isect)) {
+            /* For low epsilon values the point could intersect it's own plane. */
+            if (!ELEM(i_test, i, j, k)) {
+              break;
+            }
+          }
+        }
+
+        if (i_test == planes_len) { /* ok */
+          callback_fn(co_test, i, j, k, user_data);
+          found = true;
+        }
+      }
+    }
+  }
+
+  return found;
+}
+
 /**
  * Intersect two triangles.
  *
diff --git a/source/blender/python/mathutils/mathutils_geometry.c b/source/blender/python/mathutils/mathutils_geometry.c
index ad5f1a486b4..77ced169dab 100644
--- a/source/blender/python/mathutils/mathutils_geometry.c
+++ b/source/blender/python/mathutils/mathutils_geometry.c
@@ -1062,6 +1062,20 @@ static PyObject *M_Geometry_barycentric_transform(PyObject *UNUSED(self), PyObje
   return Vector_CreatePyObject(pt_dst, 3, NULL);
 }
 
+struct PointsInPlanes_UserData {
+  PyObject *py_verts;
+  char *planes_used;
+};
+
+static void points_in_planes_fn(const float co[3], int i, int j, int k, void *user_data_p)
+{
+  struct PointsInPlanes_UserData *user_data = user_data_p;
+  PyList_APPEND(user_data->py_verts, Vector_CreatePyObject(co, 3, NULL));
+  user_data->planes_used[i] = true;
+  user_data->planes_used[j] = true;
+  user_data->planes_used[k] = true;
+}
+
 PyDoc_STRVAR(M_Geometry_points_in_planes_doc,
              ".. function:: points_in_planes(planes)\n"
              "\n"
@@ -1073,7 +1087,6 @@ PyDoc_STRVAR(M_Geometry_points_in_planes_doc,
              "   :return: two lists, once containing the vertices inside the planes, another "
              "containing the plane indices used\n"
              "   :rtype: pair of lists\n");
-/* note: this function could be optimized by some spatial structure */
 static PyObject *M_Geometry_points_in_planes(PyObject *UNUSED(self), PyObject *args)
 {
   PyObject *py_planes;
@@ -1090,81 +1103,37 @@ static PyObject *M_Geometry_points_in_planes(PyObject *UNUSED(self), PyObject *a
   }
 
   /* note, this could be refactored into plain C easy - py bits are noted */
-  const float eps = 0.0001f;
-  const uint len = (uint)planes_len;
-  uint i, j, k, l;
 
-  float n1n2[3], n2n3[3], n3n1[3];
-  float potentialVertex[3];
-  char *planes_used = PyMem_Malloc(sizeof(char) * len);
+  struct PointsInPlanes_UserData user_data = {
+      .py_verts = PyList_New(0),
+      .planes_used = PyMem_Malloc(sizeof(char) * planes_len),
+  };
 
   /* python */
-  PyObject *py_verts = PyList_New(0);
   PyObject *py_plane_index = PyList_New(0);
 
-  memset(planes_used, 0, sizeof(char) * len);
+  memset(user_data.planes_used, 0, sizeof(char) * planes_len);
 
-  for (i = 0; i < len; i++) {
-    const float *N1 = planes[i];
-    for (j = i + 1; j < len; j++) {
-      const float *N2 = planes[j];
-      cross_v3_v3v3(n1n2, N1, N2);
-      if (len_squared_v3(n1n2) > eps) {
-        for (k = j + 1; k < len; k++) {
-          const float *N3 = planes[k];
-          cross_v3_v3v3(n2n3, N2, N3);
-          if (len_squared_v3(n2n3) > eps) {
-            cross_v3_v3v3(n3n1, N3, N1);
-            if (len_squared_v3(n3n1) > eps) {
-              const float quotient = dot_v3v3(N1, n2n3);
-              if (fabsf(quotient) > eps) {
-                /**
-                 * <pre>
-                 * potentialVertex = (
-                 *     (n2n3 * N1[3] + n3n1 * N2[3] + n1n2 * N3[3]) *
-                 *     (-1.0 / quotient));
-                 * </pre>
-                 */
-                const float quotient_ninv = -1.0f / quotient;
-                potentialVertex[0] = ((n2n3[0] * N1[3]) + (n3n1[0] * N2[3]) + (n1n2[0] * N3[3])) *
-                                     quotient_ninv;
-                potentialVertex[1] = ((n2n3[1] * N1[3]) + (n3n1[1] * N2[3]) + (n1n2[1] * N3[3])) *
-                                     quotient_ninv;
-                potentialVertex[2] = ((n2n3[2] * N1[3]) + (n3n1[2] * N2[3]) + (n1n2[2] * N3[3])) *
-                                     quotient_ninv;
-                for (l = 0; l < len; l++) {
-                  const float *NP = planes[l];
-                  if ((dot_v3v3(NP, potentialVertex) + NP[3]) > 0.000001f) {
-                    break;
-                  }
-                }
-
-                if (l == len) { /* ok */
-                  /* python */
-                  PyList_APPEND(py_verts, Vector_CreatePyObject(potentialVertex, 3, NULL));
-                  planes_used[i] = planes_used[j] = planes_used[k] = true;
-                }
-              }
-            }
-          }
-        }
-      }
-    }
-  }
+  const float eps_coplanar = 1e-4f;
+  const float eps_isect = 1e-6f;
 
+  const bool has_isect = isect_planes_v3_fn(
+      planes, planes_len, eps_coplanar, eps_isect, points_in_planes_fn, &user_data);
   PyMem_Free(planes);
 
-  /* now make a list of used planes */
-  for (i = 0; i < len; i++) {
-    if (planes_used[i]) {
-      PyList_APPEND(py_plane_index, PyLong_FromLong(i));
+  /* Now make user_data list of used planes. */
+  if (has_isect) {
+    for (int i = 0; i < planes_len; i++) {
+      if (user_data.planes_used[i]) {
+        PyList_APPEND(py_plane_index, PyLong_FromLong(i));
+      }
     }
   }
-  PyMem_Free(planes_used);
+  PyMem_Free(user_data.planes_used);
 
   {
     PyObject *ret = PyTuple_New(2);
-    PyTuple_SET_ITEMS(ret, py_verts, py_plane_index);
+    PyTuple_SET_ITEMS(ret, user_data.py_verts, py_plane_index);
     return ret;
   }
 }



More information about the Bf-blender-cvs mailing list