[gtk/wip/matthiasc/lottie-stroke: 11/17] Add some more conic math

[Matthias Clasen <matthiasc src gnome org>]

commit c7355e38f8d21aa96d9eb6b91818ff8c094a18a7
Author: Matthias Clasen <mclasen redhat com>
Date:   Wed Dec 2 21:44:54 2020 -0500

    Add some more conic math
    
    Add a function to get bounding boxes for
    rational quadratic Beziers, and a method
    to get the 'shoulder point' of a conic.

 gsk/gskpathstroke.c | 93 ++++++++++++++++++++++++++++++++++++++++++++++-------
 1 file changed, 82 insertions(+), 11 deletions(-)
---
diff --git a/gsk/gskpathstroke.c b/gsk/gskpathstroke.c
index 8bc892e2f1..615137e59c 100644
--- a/gsk/gskpathstroke.c
+++ b/gsk/gskpathstroke.c
@@ -293,6 +293,86 @@ split_conic (const graphene_point_t points[3],
   *rw = right[1].z / sqrt (right[0].z);
 }
 
+/* not sure this is useful for anything in particular */
+static void
+get_conic_shoulder_point (const graphene_point_t  p[3],
+                          float                   w,
+                          graphene_point_t       *q)
+{
+  graphene_point_t m;
+
+  graphene_point_interpolate (&p[0], &p[2], 0.5, &m);
+  graphene_point_interpolate (&m, &p[1], w / (1 + w), q);
+}
+
+static gboolean
+acceptable (float t)
+{
+  return 0 <= t && t <= 1;
+}
+
+/* Solve N = 0 where N is the numerator of derivative of P/Q, with
+ * P = (1-t)^a + 2t(1-t)wb + t^2c
+ * Q = (1-t)^2 + 2t(1-t)w + t^2
+ */
+static int
+get_conic_extrema (float a, float b, float c, float w, float t[10])
+{
+  float q, tt;
+  int n = 0;
+  float w2 = w*w;
+  float wac = (w - 1)*(a - c);
+
+  if (wac != 0)
+    {
+      q = - sqrt (a*a - 4*a*b*w2 + 4*a*c*w2 - 2*a*c + 4*b*b*w2 - 4*b*c*w2 + c*c);
+
+      tt = (- q + 2*a*w - a - 2*b*w + c)/(2*wac);
+
+      if (acceptable (tt))
+        t[n++] = tt;
+
+      tt = (q + 2*a*w - a - 2*b*w + c)/(2*wac);
+
+      if (acceptable (tt))
+        t[n++] = tt;
+    }
+
+  if (w * (b - c) != 0 && a == c)
+    t[n++] = 0.5;
+
+  if (w == 1 && a - 2*b + c != 0)
+    {
+      tt = (a - b) / (a - 2*b + c);
+      if (acceptable (tt))
+        t[n++] = tt;
+    }
+
+  return n;
+}
+
+static void
+get_conic_bounds (const graphene_point_t  p[3],
+                  float                   w,
+                  graphene_rect_t        *bounds)
+{
+  graphene_point_t q;
+  float t[10];
+  int i, n;
+
+  graphene_rect_init (bounds, p[0].x, p[0].y, 0, 0);
+  graphene_rect_expand (bounds, &p[2], bounds);
+
+  n = get_conic_extrema (p[0].x, p[1].x, p[2].x, w, t);
+  n += get_conic_extrema (p[0].y, p[1].y, p[2].y, w, &t[n]);
+
+  for (i = 0; i < n; i++)
+    {
+      get_conic (p, w, t[i], &q);
+      graphene_rect_expand (bounds, &q, bounds);
+    }
+}
+
 /* compute the angle between a, b and c in the range of [0, 360] */
 static float
 three_point_angle (const graphene_point_t *a,
@@ -330,12 +410,6 @@ cubic_is_simple (const graphene_point_t *pts)
   return TRUE;
 }
 
-static gboolean
-acceptable (float t)
-{
-  return 0 <= t && t <= 1;
-}
-
 static float
 cuberoot (float v)
 {
@@ -597,11 +671,8 @@ get_curve_bounds (const graphene_point_t  pts[4],
   float t[2];
   int n, i;
 
-  get_cubic (pts, 0, &p);
-  graphene_rect_init (bounds, p.x, p.y, 0, 0);
-
-  get_cubic (pts, 1, &p);
-  graphene_rect_expand (bounds, &p, bounds);
+  graphene_rect_init (bounds, pts[0].x, pts[0].y, 0, 0);
+  graphene_rect_expand (bounds, &pts[3], bounds);
 
   n = get_cubic_extrema (pts[0].x, pts[1].x, pts[2].x, pts[3].x, t);
   for (i = 0; i < n; i++)