Move math and vector double routines into blenlib from delaunay code

7 files changed, 127 insertions, 94 deletions

diff --git a/source/blender/blenlib/BLI_math_base.h b/source/blender/blenlib/BLI_math_base.h
index deb29605552..7b770f2dd3e 100644
--- a/source/blender/blenlib/BLI_math_base.h
+++ b/source/blender/blenlib/BLI_math_base.h
@@ -121,6 +121,9 @@ MINLINE float max_fff(float a, float b, float c);
 MINLINE float min_ffff(float a, float b, float c, float d);
 MINLINE float max_ffff(float a, float b, float c, float d);
 
+MINLINE double min_dd(double a, double b);
+MINLINE double max_dd(double a, double b);
+
 MINLINE int min_ii(int a, int b);
 MINLINE int max_ii(int a, int b);
 MINLINE int min_iii(int a, int b, int c);
diff --git a/source/blender/blenlib/BLI_math_geom.h b/source/blender/blenlib/BLI_math_geom.h
index 39b1b96d009..caecc0aebc2 100644
--- a/source/blender/blenlib/BLI_math_geom.h
+++ b/source/blender/blenlib/BLI_math_geom.h
@@ -188,6 +188,10 @@ float dist_squared_to_projected_aabb_simple(const float projmat[4][4],
                                             const float bbmax[3]);
 
 float closest_to_line_v2(float r_close[2], const float p[2], const float l1[2], const float l2[2]);
+double closest_to_line_v2_db(double r_close[2],
+                             const double p[2],
+                             const double l1[2],
+                             const double l2[2]);
 float closest_to_line_v3(float r_close[3], const float p[3], const float l1[3], const float l2[3]);
 void closest_to_line_segment_v2(float r_close[2],
                                 const float p[2],
@@ -267,7 +271,12 @@ bool isect_seg_seg_v2_simple(const float v1[2],
                              const float v2[2],
                              const float v3[2],
                              const float v4[2]);
-
+int isect_seg_seg_v2_lambda_mu_db(const double v1[2],
+                                  const double v2[2],
+                                  const double v3[2],
+                                  const double v4[2],
+                                  double *r_lambda,
+                                  double *r_mu);
 int isect_line_sphere_v3(const float l1[3],
                          const float l2[3],
                          const float sp[3],
diff --git a/source/blender/blenlib/BLI_math_vector.h b/source/blender/blenlib/BLI_math_vector.h
index ccb42683d5a..31d725af5ad 100644
--- a/source/blender/blenlib/BLI_math_vector.h
+++ b/source/blender/blenlib/BLI_math_vector.h
@@ -103,6 +103,7 @@ MINLINE void add_v2_fl(float r[2], float f);
 MINLINE void add_v3_fl(float r[3], float f);
 MINLINE void add_v4_fl(float r[4], float f);
 MINLINE void add_v2_v2(float r[2], const float a[2]);
+MINLINE void add_v2_v2_db(double r[2], const double a[2]);
 MINLINE void add_v2_v2v2(float r[2], const float a[2], const float b[2]);
 MINLINE void add_v2_v2v2_int(int r[2], const int a[2], const int b[2]);
 MINLINE void add_v3_v3(float r[3], const float a[3]);
@@ -116,6 +117,7 @@ MINLINE void add_v3fl_v3fl_v3s(float r[3], const float a[3], const short b[3]);
 
 MINLINE void sub_v2_v2(float r[2], const float a[2]);
 MINLINE void sub_v2_v2v2(float r[2], const float a[2], const float b[2]);
+MINLINE void sub_v2_v2v2_db(double r[2], const double a[2], const double b[2]);
 MINLINE void sub_v2_v2v2_int(int r[2], const int a[2], const int b[2]);
 MINLINE void sub_v3_v3(float r[3], const float a[3]);
 MINLINE void sub_v3_v3v3(float r[3], const float a[3], const float b[3]);
@@ -184,6 +186,7 @@ MINLINE void abs_v4(float r[4]);
 MINLINE void abs_v4_v4(float r[4], const float a[4]);
 
 MINLINE float dot_v2v2(const float a[2], const float b[2]) ATTR_WARN_UNUSED_RESULT;
+MINLINE double dot_v2v2_db(const double a[2], const double b[2]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float dot_v3v3(const float a[3], const float b[3]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float dot_v3v3v3(const float p[3],
                          const float a[3],
@@ -212,8 +215,10 @@ MINLINE int len_manhattan_v2_int(const int v[2]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float len_manhattan_v3(const float v[3]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float len_v2(const float a[2]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float len_v2v2(const float a[2], const float b[2]) ATTR_WARN_UNUSED_RESULT;
+MINLINE double len_v2v2_db(const double a[2], const double b[2]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float len_v2v2_int(const int v1[2], const int v2[2]);
 MINLINE float len_squared_v2v2(const float a[2], const float b[2]) ATTR_WARN_UNUSED_RESULT;
+MINLINE double len_squared_v2v2_db(const double a[2], const double b[2]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float len_squared_v3v3(const float a[3], const float b[3]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float len_squared_v4v4(const float a[4], const float b[4]) ATTR_WARN_UNUSED_RESULT;
 MINLINE float len_manhattan_v2v2(const float a[2], const float b[2]) ATTR_WARN_UNUSED_RESULT;
diff --git a/source/blender/blenlib/intern/delaunay_2d.c b/source/blender/blenlib/intern/delaunay_2d.c
index 23f560c5463..102ea8eb234 100644
--- a/source/blender/blenlib/intern/delaunay_2d.c
+++ b/source/blender/blenlib/intern/delaunay_2d.c
@@ -115,99 +115,6 @@ static void validate_face_centroid(SymEdge *se);
 static void validate_cdt(CDT_state *cdt, bool check_all_tris);
 #endif
 
-/* TODO: move these to BLI_vector... and BLI_math... */
-static double max_dd(const double a, const double b)
-{
-  return (a > b) ? a : b;
-}
-
-static double len_v2v2_db(const double a[2], const double b[2])
-{
-  return sqrt((b[0] - a[0]) * (b[0] - a[0]) + (b[1] - a[1]) * (b[1] - a[1]));
-}
-
-static double len_squared_v2v2_db(const double a[2], const double b[2])
-{
-  return (b[0] - a[0]) * (b[0] - a[0]) + (b[1] - a[1]) * (b[1] - a[1]);
-}
-
-static void add_v2_v2_db(double a[2], const double b[2])
-{
-  a[0] += b[0];
-  a[1] += b[1];
-}
-
-static void sub_v2_v2v2_db(double *a, const double *b, const double *c)
-{
-  a[0] = b[0] - c[0];
-  a[1] = b[1] - c[1];
-}
-
-static double dot_v2v2_db(const double *a, const double *b)
-{
-  return a[0] * b[0] + a[1] * b[1];
-}
-
-static double closest_to_line_v2_db(double r_close[2],
-                                    const double p[2],
-                                    const double l1[2],
-                                    const double l2[2])
-{
-  double h[2], u[2], lambda, denom;
-  sub_v2_v2v2_db(u, l2, l1);
-  sub_v2_v2v2_db(h, p, l1);
-  denom = dot_v2v2_db(u, u);
-  if (denom < DBL_EPSILON) {
-    r_close[0] = l1[0];
-    r_close[1] = l1[1];
-    return 0.0;
-  }
-  lambda = dot_v2v2_db(u, h) / dot_v2v2_db(u, u);
-  r_close[0] = l1[0] + u[0] * lambda;
-  r_close[1] = l1[1] + u[1] * lambda;
-  return lambda;
-}
-
-/**
- * If intersection == ISECT_LINE_LINE_CROSS or ISECT_LINE_LINE_NONE:
- * <pre>
- * pt = v1 + lamba * (v2 - v1) = v3 + mu * (v4 - v3)
- * </pre>
- */
-static int isect_seg_seg_v2_lambda_mu_db(const double v1[2],
-                                         const double v2[2],
-                                         const double v3[2],
-                                         const double v4[2],
-                                         double *r_lambda,
-                                         double *r_mu)
-{
-  double div, lambda, mu;
-
-  div = (v2[0] - v1[0]) * (v4[1] - v3[1]) - (v2[1] - v1[1]) * (v4[0] - v3[0]);
-  if (fabs(div) < DBL_EPSILON) {
-    return ISECT_LINE_LINE_COLINEAR;
-  }
-
-  lambda = ((v1[1] - v3[1]) * (v4[0] - v3[0]) - (v1[0] - v3[0]) * (v4[1] - v3[1])) / div;
-
-  mu = ((v1[1] - v3[1]) * (v2[0] - v1[0]) - (v1[0] - v3[0]) * (v2[1] - v1[1])) / div;
-
-  if (r_lambda) {
-    *r_lambda = lambda;
-  }
-  if (r_mu) {
-    *r_mu = mu;
-  }
-
-  if (lambda >= 0.0 && lambda <= 1.0 && mu >= 0.0 && mu <= 1.0) {
-    if (lambda == 0.0 || lambda == 1.0 || mu == 0.0 || mu == 1.0) {
-      return ISECT_LINE_LINE_EXACT;
-    }
-    return ISECT_LINE_LINE_CROSS;
-  }
-  return ISECT_LINE_LINE_NONE;
-}
-
 /** return 1 if a,b,c forms CCW angle, -1 if a CW angle, 0 if straight  */
 static int CCW_test(const double a[2], const double b[2], const double c[2])
 {
diff --git a/source/blender/blenlib/intern/math_base_inline.c b/source/blender/blenlib/intern/math_base_inline.c
index 320f8a0f1ab..47327a878d4 100644
--- a/source/blender/blenlib/intern/math_base_inline.c
+++ b/source/blender/blenlib/intern/math_base_inline.c
@@ -353,6 +353,15 @@ MINLINE float max_ff(float a, float b)
   return (a > b) ? a : b;
 }
 
+MINLINE double min_dd(double a, double b)
+{
+  return (a < b) ? a : b;
+}
+MINLINE double max_dd(double a, double b)
+{
+  return (a > b) ? a : b;
+}
+
 MINLINE int min_ii(int a, int b)
 {
   return (a < b) ? a : b;
diff --git a/source/blender/blenlib/intern/math_geom.c b/source/blender/blenlib/intern/math_geom.c
index c4bdc73e0e3..0162f4552f0 100644
--- a/source/blender/blenlib/intern/math_geom.c
+++ b/source/blender/blenlib/intern/math_geom.c
@@ -1356,6 +1356,52 @@ bool isect_seg_seg_v2_simple(const float v1[2],
 }
 
 /**
+ * If intersection == ISECT_LINE_LINE_CROSS or ISECT_LINE_LINE_NONE:
+ * <pre>
+ * pt = v1 + lamba * (v2 - v1) = v3 + mu * (v4 - v3)
+ * </pre>
+ * \returns intersection type:
+ * - ISECT_LINE_LINE_COLINEAR: collinear.
+ * - ISECT_LINE_LINE_EXACT: intersection at an endpoint of either.
+ * - ISECT_LINE_LINE_CROSS: interection, not at an endpoint.
+ * - ISECT_LINE_LINE_NONE: no intersection.
+ * Also returns lambda and mu in r_lambda and r_mu.
+ */
+int isect_seg_seg_v2_lambda_mu_db(const double v1[2],
+                                  const double v2[2],
+                                  const double v3[2],
+                                  const double v4[2],
+                                  double *r_lambda,
+                                  double *r_mu)
+{
+  double div, lambda, mu;
+
+  div = (v2[0] - v1[0]) * (v4[1] - v3[1]) - (v2[1] - v1[1]) * (v4[0] - v3[0]);
+  if (fabs(div) < DBL_EPSILON) {
+    return ISECT_LINE_LINE_COLINEAR;
+  }
+
+  lambda = ((v1[1] - v3[1]) * (v4[0] - v3[0]) - (v1[0] - v3[0]) * (v4[1] - v3[1])) / div;
+
+  mu = ((v1[1] - v3[1]) * (v2[0] - v1[0]) - (v1[0] - v3[0]) * (v2[1] - v1[1])) / div;
+
+  if (r_lambda) {
+    *r_lambda = lambda;
+  }
+  if (r_mu) {
+    *r_mu = mu;
+  }
+
+  if (lambda >= 0.0 && lambda <= 1.0 && mu >= 0.0 && mu <= 1.0) {
+    if (lambda == 0.0 || lambda == 1.0 || mu == 0.0 || mu == 1.0) {
+      return ISECT_LINE_LINE_EXACT;
+    }
+    return ISECT_LINE_LINE_CROSS;
+  }
+  return ISECT_LINE_LINE_NONE;
+}
+
+/**
  * \param l1, l2: Coordinates (point of line).
  * \param sp, r:  Coordinate and radius (sphere).
  * \return r_p1, r_p2: Intersection coordinates.
@@ -3187,6 +3233,26 @@ float closest_to_line_v2(float r_close[2], const float p[2], const float l1[2],
   return lambda;
 }
 
+double closest_to_line_v2_db(double r_close[2],
+                             const double p[2],
+                             const double l1[2],
+                             const double l2[2])
+{
+  double h[2], u[2], lambda, denom;
+  sub_v2_v2v2_db(u, l2, l1);
+  sub_v2_v2v2_db(h, p, l1);
+  denom = dot_v2v2_db(u, u);
+  if (denom < DBL_EPSILON) {
+    r_close[0] = l1[0];
+    r_close[1] = l1[1];
+    return 0.0;
+  }
+  lambda = dot_v2v2_db(u, h) / dot_v2v2_db(u, u);
+  r_close[0] = l1[0] + u[0] * lambda;
+  r_close[1] = l1[1] + u[1] * lambda;
+  return lambda;
+}
+
 float ray_point_factor_v3_ex(const float p[3],
                              const float ray_origin[3],
                              const float ray_direction[3],
diff --git a/source/blender/blenlib/intern/math_vector_inline.c b/source/blender/blenlib/intern/math_vector_inline.c
index 8cb618ae14e..b4b53a1dd58 100644
--- a/source/blender/blenlib/intern/math_vector_inline.c
+++ b/source/blender/blenlib/intern/math_vector_inline.c
@@ -355,6 +355,12 @@ MINLINE void add_v2_v2(float r[2], const float a[2])
   r[1] += a[1];
 }
 
+MINLINE void add_v2_v2_db(double r[2], const double a[2])
+{
+  r[0] += a[0];
+  r[1] += a[1];
+}
+
 MINLINE void add_v2_v2v2(float r[2], const float a[2], const float b[2])
 {
   r[0] = a[0] + b[0];
@@ -430,6 +436,12 @@ MINLINE void sub_v2_v2v2(float r[2], const float a[2], const float b[2])
   r[1] = a[1] - b[1];
 }
 
+MINLINE void sub_v2_v2v2_db(double r[2], const double a[2], const double b[2])
+{
+  r[0] = a[0] - b[0];
+  r[1] = a[1] - b[1];
+}
+
 MINLINE void sub_v2_v2v2_int(int r[2], const int a[2], const int b[2])
 {
   r[0] = a[0] - b[0];
@@ -830,6 +842,11 @@ MINLINE float dot_v2v2(const float a[2], const float b[2])
   return a[0] * b[0] + a[1] * b[1];
 }
 
+MINLINE double dot_v2v2_db(const double a[2], const double b[2])
+{
+  return a[0] * b[0] + a[1] * b[1];
+}
+
 MINLINE float dot_v3v3(const float a[3], const float b[3])
 {
   return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
@@ -956,6 +973,15 @@ MINLINE float len_v2v2(const float v1[2], const float v2[2])
   return sqrtf(x * x + y * y);
 }
 
+MINLINE double len_v2v2_db(const double v1[2], const double v2[2])
+{
+  double x, y;
+
+  x = v1[0] - v2[0];
+  y = v1[1] - v2[1];
+  return sqrt(x * x + y * y);
+}
+
 MINLINE float len_v2v2_int(const int v1[2], const int v2[2])
 {
   float x, y;
@@ -978,6 +1004,14 @@ MINLINE float len_squared_v2v2(const float a[2], const float b[2])
   return dot_v2v2(d, d);
 }
 
+MINLINE double len_squared_v2v2_db(const double a[2], const double b[2])
+{
+  double d[2];
+
+  sub_v2_v2v2_db(d, b, a);
+  return dot_v2v2_db(d, d);
+}
+
 MINLINE float len_squared_v3v3(const float a[3], const float b[3])
 {
   float d[3];