Implement Conics in SkCurveMeasure

src/utils/SkCurveMeasure.cpp

 /*
  * Copyright 2016 Google Inc.
  *
  * Use of this source code is governed by a BSD-style license that can be
  * found in the LICENSE file.
  */
 
 #include "SkCurveMeasure.h"
 #include "SkGeometry.h"
 
 // for abs
 #include <cmath>
 
 #define UNIMPLEMENTED SkDEBUGF(("%s:%d unimplemented\n", __FILE__, __LINE__))

[caryclark, 2016/08/15 19:13:04]

Make this an assert.

[Harry Stern, 2016/08/15 20:55:16]

Done, depends on https://codereview.chromium.org/2244253003/ to build.

 
 /// Used inside SkCurveMeasure::getTime's Newton's iteration
 static inline SkPoint evaluate(const SkPoint pts[4], SkSegType segType,
                                SkScalar t) {
     SkPoint pos;
     switch (segType) {
         case kQuad_SegType:
             pos = SkEvalQuadAt(pts, t);
             break;
         case kLine_SegType:
@@ skipping 37 matching lines @@
         default:
             UNIMPLEMENTED;
     }
 
     return tan;
 }
 /// Used in ArcLengthIntegrator::computeLength
 static inline Sk8f evaluateDerivativeLength(const Sk8f& ts,
                                             const Sk8f (&xCoeff)[3],
                                             const Sk8f (&yCoeff)[3],
+                                            const SkScalar w,
                                             const SkSegType segType) {
     Sk8f x;
     Sk8f y;
     switch (segType) {
         case kQuad_SegType:
             x = xCoeff[0]*ts + xCoeff[1];
             y = yCoeff[0]*ts + yCoeff[1];
             break;
         case kCubic_SegType:
             x = (xCoeff[0]*ts + xCoeff[1])*ts + xCoeff[2];
             y = (yCoeff[0]*ts + yCoeff[1])*ts + yCoeff[2];
             break;
-        case kConic_SegType:
-            UNIMPLEMENTED;
+        case kConic_SegType: {
+            Sk8f xnum = (xCoeff[0]*ts + xCoeff[1])*ts + xCoeff[2];
+            Sk8f ynum = (yCoeff[0]*ts + yCoeff[1])*ts + yCoeff[2];
+
+            Sk8f denom =  ts*(-4 + 4*w + ts*(8 - 12*w + 4*w*w + ts*(-8 + 16*w -
+                            8*w*w + ts*(4 - 8*w + 4*w*w)))) + 1;

[caryclark, 2016/08/15 19:13:04]

Where did this come from? That is, please document your derivation.

C++ compilers will not do a good job of refactoring this. I see trivial common
terms:
-8 + 16*w - 8*w*w
 4 -  8*w + 4*w*w

Have you looked at the compiled output? Are you getting the SIMD you expect?

[Harry Stern, 2016/08/15 20:55:16]

I took out those terms. I honestly didn't notice them when I looked because I
was too fixated on the 12.

I am getting SIMD from this function in the output. `objdump -S -w -d
out/Release/obj/src/utils/utils.SkCurveMeasure.o` is showing me that the
function was inlined and is using e.g. mulps, addps, sqrtps

[Harry Stern, 2016/08/17 02:48:46]

Also, I have added the mathematica commands used to generate the coefficients.
Is there something better I could do?

+
+            x = xnum / denom;
+            y = ynum / denom;

[caryclark, 2016/08/15 19:13:04]

Is there a unit test for this?

[Harry Stern, 2016/08/15 20:55:16]

It's on my todo by tomorrow.

+        }
             break;
         default:
             UNIMPLEMENTED;
     }
 
     x = x * x;
     y = y * y;
 
     return (x + y).sqrt();
 }
@@ skipping 31 matching lines @@
             // precompute coefficients for derivative
             xCoeff[0] = Sk8f(3.0f*(-Ax + 3.0f*(Bx - Cx) + Dx));
             xCoeff[1] = Sk8f(3.0f*(2.0f*(Ax - 2.0f*Bx + Cx)));
             xCoeff[2] = Sk8f(3.0f*(-Ax + Bx));
 
             yCoeff[0] = Sk8f(3.0f*(-Ay + 3.0f*(By - Cy) + Dy));
             yCoeff[1] = Sk8f(3.0f*(2.0f*(Ay - 2.0f*By + Cy)));
             yCoeff[2] = Sk8f(3.0f*(-Ay + By));
         }
             break;
-        case kConic_SegType:
-            UNIMPLEMENTED;
+        case kConic_SegType: {
+            float Ax = pts[0].x();
+            float Bx = pts[1].x();
+            float Cx = pts[2].x();
+            float Ay = pts[0].y();
+            float By = pts[1].y();
+            float Cy = pts[2].y();
+
+            fConicW = pts[3].x();

[Harry Stern, 2016/08/15 18:01:21]

fConicW gets set here

+
+            SkScalar w = fConicW;
+
+            // precompute coefficients for derivative
+            xCoeff[0] = Sk8f(2*(Ax - Cx - Ax*w + Cx*w));
+            xCoeff[1] = Sk8f(-2*(Ax + Cx + 2*(Ax*w - Bx*w)));
+            xCoeff[2] = Sk8f(-2*(Ax*w + Bx*w));
+
+            yCoeff[0] = Sk8f(2*(Ay - Cy - Ay*w + Cy*w));
+            yCoeff[1] = Sk8f(-2*(Ay + Cy + 2*(Ay*w - By*w)));
+            yCoeff[2] = Sk8f(-2*(Ay*w + By*w));
+        }
             break;
         default:
             UNIMPLEMENTED;
     }
 }
 
 // We use Gaussian quadrature
 // (https://en.wikipedia.org/wiki/Gaussian_quadrature)
 // to approximate the arc length integral here, because it is amenable to SIMD.
 SkScalar ArcLengthIntegrator::computeLength(SkScalar t) {
     SkScalar length = 0.0f;
 
-    Sk8f lengths = evaluateDerivativeLength(absc*t, xCoeff, yCoeff, fSegType);
+    Sk8f lengths = evaluateDerivativeLength(absc*t, xCoeff, yCoeff, fConicW, fSegType);
     lengths = weights*lengths;
     // is it faster or more accurate to sum and then multiply or vice versa?
     lengths = lengths*(t*0.5f);
 
     // Why does SkNx index with ints? does negative index mean something?
     for (int i = 0; i < 8; i++) {
         length += lengths[i];
     }
     return length;
 }
@@ skipping 139 matching lines @@
     if (time) {
         *time = t;
     }
     if (pos) {
         *pos = evaluate(fPts, fSegType, t);
     }
     if (tan) {
         *tan = evaluateDerivative(fPts, fSegType, t);
     }
 }