move explicit bezier to vg

diff --git a/vg_m.h b/vg_m.h
index dae258c9f8a8b4722bd23ae73952df78bf9da571..5ab6a22d8620d459498e15f99194a3d0ebdc33fc 100644 (file)
--- a/vg_m.h
+++ b/vg_m.h
@@ -2303,6 +2303,56 @@ static void eval_bezier3( v3f p0, v3f p1, v3f p2, f32 t, v3f p )
    v3_muladds( p, p2, t*t, p );
 }
 
+static f32 explicit_bezier( f32 A[2], f32 B[2], f32 C[2], f32 D[2], f32 x )
+{
+   f32 dAxDx = D[0]-A[0],
+       unitBx = (B[0] - A[0]) / dAxDx,
+       unitCx = (C[0] - A[0]) / dAxDx,
+
+       /* cubic coefficients */
+       a =  3.0f*unitBx - 3.0f*unitCx + 1.0f,
+       b = -6.0f*unitBx + 3.0f*unitCx,
+       c =  3.0f*unitBx,
+       d = -(x - A[0]) / dAxDx,
+
+        t0 = 0.0f,
+       Ft0 = d,
+        t1 = 1.0f,
+       Ft1 = a+b+c+d,
+        tc, Ftcx;
+
+   /* Illinois method to find root */
+   for( u32 j=0; j<8; j ++ )
+   {
+      tc = t1 - Ft1*(t1-t0)/(Ft1-Ft0);
+      Ftcx = tc*tc*tc*a + tc*tc*b + tc*c + d;
+
+      if( fabsf(Ftcx) < 0.00001f )
+         break;
+      
+      if( Ft1*Ftcx < 0.0f )
+      {
+         t0 = t1;
+         Ft0 = Ft1;
+      }
+      else
+         Ft0 *= 0.5f;
+
+      t1 = tc;
+      Ft1 = Ftcx;
+   }
+
+   /* Evaluate parametric bezier */
+   f32 t2 = tc*tc,
+       t3 = tc*tc*tc;
+   
+   return  D[1] * t3
+         + C[1] * (-3.0f*t3 + 3.0f*t2)
+         + B[1] * ( 3.0f*t3 - 6.0f*t2 + 3.0f*tc)
+         + A[1] * (-1.0f*t3 + 3.0f*t2 - 3.0f*tc + 1.0f);
+}
+
+
 /*
  * -----------------------------------------------------------------------------
  * Section 5.f                      Volumes