Add direct bezier cubic support for GPU shaders

BUG=
R=bsalomon@google.com, jvanverth@google.com, robertphillips@google.com

Author: egdaniel@google.com

Review URL: https://chromiumcodereview.appspot.com/22900007

git-svn-id: http://skia.googlecode.com/svn/trunk@10814 2bbb7eff-a529-9590-31e7-b0007b416f81

3 files changed, 473 insertions, 0 deletions

diff --git a/src/gpu/GrAAHairLinePathRenderer.cpp b/src/gpu/GrAAHairLinePathRenderer.cpp
index 5ef981be62..2b5f3cc055 100644
--- a/src/gpu/GrAAHairLinePathRenderer.cpp
+++ b/src/gpu/GrAAHairLinePathRenderer.cpp
@@ -738,6 +738,105 @@ void add_line(const SkPoint p[2],
 }
 
 /**
+ * Shader is based off of "Resolution Independent Curve Rendering using
+ * Programmable Graphics Hardware" by Loop and Blinn.
+ * The output of this effect is a hairline edge for non rational cubics.
+ * Cubics are specified by implicit equation K^3 - LM.
+ * K, L, and M, are the first three values of the vertex attribute,
+ * the fourth value is not used. Distance is calculated using a
+ * first order approximation from the taylor series.
+ * Coverage is max(0, 1-distance).
+ */
+class HairCubicEdgeEffect : public GrEffect {
+public:
+    static GrEffectRef* Create() {
+        GR_CREATE_STATIC_EFFECT(gHairCubicEdgeEffect, HairCubicEdgeEffect, ());
+        gHairCubicEdgeEffect->ref();
+        return gHairCubicEdgeEffect;
+    }
+
+    virtual ~HairCubicEdgeEffect() {}
+
+    static const char* Name() { return "HairCubicEdge"; }
+
+    virtual void getConstantColorComponents(GrColor* color,
+                                            uint32_t* validFlags) const SK_OVERRIDE {
+        *validFlags = 0;
+    }
+
+    virtual const GrBackendEffectFactory& getFactory() const SK_OVERRIDE {
+        return GrTBackendEffectFactory<HairCubicEdgeEffect>::getInstance();
+    }
+
+    class GLEffect : public GrGLEffect {
+    public:
+        GLEffect(const GrBackendEffectFactory& factory, const GrDrawEffect&)
+            : INHERITED (factory) {}
+
+        virtual void emitCode(GrGLShaderBuilder* builder,
+                              const GrDrawEffect& drawEffect,
+                              EffectKey key,
+                              const char* outputColor,
+                              const char* inputColor,
+                              const TextureSamplerArray& samplers) SK_OVERRIDE {
+            const char *vsName, *fsName;
+
+            SkAssertResult(builder->enableFeature(
+                    GrGLShaderBuilder::kStandardDerivatives_GLSLFeature));
+            builder->addVarying(kVec4f_GrSLType, "CubicCoeffs",
+                                &vsName, &fsName);
+            const SkString* attr0Name =
+                builder->getEffectAttributeName(drawEffect.getVertexAttribIndices()[0]);
+            builder->vsCodeAppendf("\t%s = %s;\n", vsName, attr0Name->c_str());
+
+            builder->fsCodeAppend("\t\tfloat edgeAlpha;\n");
+
+            builder->fsCodeAppendf("\t\tvec3 dklmdx = dFdx(%s.xyz);\n", fsName);
+            builder->fsCodeAppendf("\t\tvec3 dklmdy = dFdy(%s.xyz);\n", fsName);
+            builder->fsCodeAppendf("\t\tfloat dfdx =\n"
+                                   "\t\t3.0*%s.x*%s.x*dklmdx.x - %s.y*dklmdx.z - %s.z*dklmdx.y;\n",
+                                   fsName, fsName, fsName, fsName);
+            builder->fsCodeAppendf("\t\tfloat dfdy =\n"
+                                   "\t\t3.0*%s.x*%s.x*dklmdy.x - %s.y*dklmdy.z - %s.z*dklmdy.y;\n",
+                                   fsName, fsName, fsName, fsName);
+            builder->fsCodeAppend("\t\tvec2 gF = vec2(dfdx, dfdy);\n");
+            builder->fsCodeAppend("\t\tfloat gFM = sqrt(dot(gF, gF));\n");
+            builder->fsCodeAppendf("\t\tfloat func = abs(%s.x*%s.x*%s.x - %s.y*%s.z);\n",
+                                   fsName, fsName, fsName, fsName, fsName);
+            builder->fsCodeAppend("\t\tedgeAlpha = func / gFM;\n");
+            builder->fsCodeAppend("\t\tedgeAlpha = max(1.0 - edgeAlpha, 0.0);\n");
+            // Add line below for smooth cubic ramp
+            // builder->fsCodeAppend("\t\tedgeAlpha = edgeAlpha*edgeAlpha*(3.0-2.0*edgeAlpha);\n");
+
+            SkString modulate;
+            GrGLSLModulatef<4>(&modulate, inputColor, "edgeAlpha");
+            builder->fsCodeAppendf("\t%s = %s;\n", outputColor, modulate.c_str());
+        }
+
+        static inline EffectKey GenKey(const GrDrawEffect& drawEffect, const GrGLCaps&) {
+            return 0x0;
+        }
+
+        virtual void setData(const GrGLUniformManager&, const GrDrawEffect&) SK_OVERRIDE {}
+
+    private:
+        typedef GrGLEffect INHERITED;
+    };
+private:
+    HairCubicEdgeEffect() {
+        this->addVertexAttrib(kVec4f_GrSLType);
+    }
+
+    virtual bool onIsEqual(const GrEffect& other) const SK_OVERRIDE {
+        return true;
+    }
+
+    GR_DECLARE_EFFECT_TEST;
+
+    typedef GrEffect INHERITED;
+};
+
+/**
  * Shader is based off of Loop-Blinn Quadratic GPU Rendering
  * The output of this effect is a hairline edge for conics.
  * Conics specified by implicit equation K^2 - LM.
diff --git a/src/gpu/GrPathUtils.cpp b/src/gpu/GrPathUtils.cpp
index 2d853883d7..ca878338d3 100644
--- a/src/gpu/GrPathUtils.cpp
+++ b/src/gpu/GrPathUtils.cpp
@@ -476,3 +476,332 @@ void GrPathUtils::convertCubicToQuads(const GrPoint p[4],
     }
 
 }
+
+////////////////////////////////////////////////////////////////////////////////
+
+enum CubicType {
+    kSerpentine_CubicType,
+    kCusp_CubicType,
+    kLoop_CubicType,
+    kQuadratic_CubicType,
+    kLine_CubicType,
+    kPoint_CubicType
+};
+
+// discr(I) = d0^2 * (3*d1^2 - 4*d0*d2)
+// Classification:
+// discr(I) > 0        Serpentine
+// discr(I) = 0        Cusp
+// discr(I) < 0        Loop
+// d0 = d1 = 0         Quadratic
+// d0 = d1 = d2 = 0    Line
+// p0 = p1 = p2 = p3   Point
+static CubicType classify_cubic(const SkPoint p[4], const SkScalar d[3]) {
+    if (p[0] == p[1] && p[0] == p[2] && p[0] == p[3]) {
+        return kPoint_CubicType;
+    }
+    const SkScalar discr = d[0] * d[0] * (3.f * d[1] * d[1] - 4.f * d[0] * d[2]);
+    if (discr > SK_ScalarNearlyZero) {
+        return kSerpentine_CubicType;
+    } else if (discr < -SK_ScalarNearlyZero) {
+        return kLoop_CubicType;
+    } else {
+        if (0.f == d[0] && 0.f == d[1]) {
+            return (0.f == d[2] ? kLine_CubicType : kQuadratic_CubicType);
+        } else {
+            return kCusp_CubicType;
+        }
+    }
+}
+
+// Assumes the third component of points is 1.
+// Calcs p0 . (p1 x p2)
+static SkScalar calc_dot_cross_cubic(const SkPoint& p0, const SkPoint& p1, const SkPoint& p2) {
+    const SkScalar xComp = p0.fX * (p1.fY - p2.fY);
+    const SkScalar yComp = p0.fY * (p2.fX - p1.fX);
+    const SkScalar wComp = p1.fX * p2.fY - p1.fY * p2.fX;
+    return (xComp + yComp + wComp);
+}
+
+// Solves linear system to extract klm
+// P.K = k (similarly for l, m)
+// Where P is matrix of control points
+// K is coefficients for the line K
+// k is vector of values of K evaluated at the control points
+// Solving for K, thus K = P^(-1) . k
+static void calc_cubic_klm(const SkPoint p[4], const SkScalar controlK[4],
+                           const SkScalar controlL[4], const SkScalar controlM[4],
+                           SkScalar k[3], SkScalar l[3], SkScalar m[3]) {
+    SkMatrix matrix;
+    matrix.setAll(p[0].fX, p[0].fY, 1.f,
+                  p[1].fX, p[1].fY, 1.f,
+                  p[2].fX, p[2].fY, 1.f);
+    SkMatrix inverse;
+    if (matrix.invert(&inverse)) {
+       inverse.mapHomogeneousPoints(k, controlK, 1);
+       inverse.mapHomogeneousPoints(l, controlL, 1);
+       inverse.mapHomogeneousPoints(m, controlM, 1);
+    }
+
+}
+
+static void set_serp_klm(const SkScalar d[3], SkScalar k[4], SkScalar l[4], SkScalar m[4]) {
+    SkScalar tempSqrt = SkScalarSqrt(9.f * d[1] * d[1] - 12.f * d[0] * d[2]);
+    SkScalar ls = 3.f * d[1] - tempSqrt;
+    SkScalar lt = 6.f * d[0];
+    SkScalar ms = 3.f * d[1] + tempSqrt;
+    SkScalar mt = 6.f * d[0];
+
+    k[0] = ls * ms;
+    k[1] = (3.f * ls * ms - ls * mt - lt * ms) / 3.f;
+    k[2] = (lt * (mt - 2.f * ms) + ls * (3.f * ms - 2.f * mt)) / 3.f;
+    k[3] = (lt - ls) * (mt - ms);
+
+    l[0] = ls * ls * ls;
+    const SkScalar lt_ls = lt - ls;
+    l[1] = ls * ls * lt_ls * -1.f;
+    l[2] = lt_ls * lt_ls * ls;
+    l[3] = -1.f * lt_ls * lt_ls * lt_ls;
+
+    m[0] = ms * ms * ms;
+    const SkScalar mt_ms = mt - ms;
+    m[1] = ms * ms * mt_ms * -1.f;
+    m[2] = mt_ms * mt_ms * ms;
+    m[3] = -1.f * mt_ms * mt_ms * mt_ms;
+
+    // If d0 < 0 we need to flip the orientation of our curve
+    // This is done by negating the k and l values
+    // We want negative distance values to be on the inside
+    if ( d[0] > 0) {
+        for (int i = 0; i < 4; ++i) {
+            k[i] = -k[i];
+            l[i] = -l[i];
+        }
+    }
+}
+
+static void set_loop_klm(const SkScalar d[3], SkScalar k[4], SkScalar l[4], SkScalar m[4]) {
+    SkScalar tempSqrt = SkScalarSqrt(4.f * d[0] * d[2] - 3.f * d[1] * d[1]);
+    SkScalar ls = d[1] - tempSqrt;
+    SkScalar lt = 2.f * d[0];
+    SkScalar ms = d[1] + tempSqrt;
+    SkScalar mt = 2.f * d[0];
+
+    k[0] = ls * ms;
+    k[1] = (3.f * ls*ms - ls * mt - lt * ms) / 3.f;
+    k[2] = (lt * (mt - 2.f * ms) + ls * (3.f * ms - 2.f * mt)) / 3.f;
+    k[3] = (lt - ls) * (mt - ms);
+
+    l[0] = ls * ls * ms;
+    l[1] = (ls * (ls * (mt - 3.f * ms) + 2.f * lt * ms))/-3.f;
+    l[2] = ((lt - ls) * (ls * (2.f * mt - 3.f * ms) + lt * ms))/3.f;
+    l[3] = -1.f * (lt - ls) * (lt - ls) * (mt - ms);
+
+    m[0] = ls * ms * ms;
+    m[1] = (ms * (ls * (2.f * mt - 3.f * ms) + lt * ms))/-3.f;
+    m[2] = ((mt - ms) * (ls * (mt - 3.f * ms) + 2.f * lt * ms))/3.f;
+    m[3] = -1.f * (lt - ls) * (mt - ms) * (mt - ms);
+
+
+    // If (d0 < 0 && sign(k1) > 0) || (d0 > 0 && sign(k1) < 0),
+    // we need to flip the orientation of our curve.
+    // This is done by negating the k and l values
+    if ( (d[0] < 0 && k[1] < 0) || (d[0] > 0 && k[1] > 0)) {
+        for (int i = 0; i < 4; ++i) {
+            k[i] = -k[i];
+            l[i] = -l[i];
+        }
+    }
+}
+
+static void set_cusp_klm(const SkScalar d[3], SkScalar k[4], SkScalar l[4], SkScalar m[4]) {
+    const SkScalar ls = d[2];
+    const SkScalar lt = 3.f * d[1];
+
+    k[0] = ls;
+    k[1] = ls - lt / 3.f;
+    k[2] = ls - 2.f * lt / 3.f;
+    k[3] = ls - lt;
+
+    l[0] = ls * ls * ls;
+    const SkScalar ls_lt = ls - lt;
+    l[1] = ls * ls * ls_lt;
+    l[2] = ls_lt * ls_lt * ls;
+    l[3] = ls_lt * ls_lt * ls_lt;
+
+    m[0] = 1.f;
+    m[1] = 1.f;
+    m[2] = 1.f;
+    m[3] = 1.f;
+}
+
+// For the case when a cubic is actually a quadratic
+// M =
+// 0     0     0
+// 1/3   0     1/3
+// 2/3   1/3   2/3
+// 1     1     1
+static void set_quadratic_klm(const SkScalar d[3], SkScalar k[4], SkScalar l[4], SkScalar m[4]) {
+    k[0] = 0.f;
+    k[1] = 1.f/3.f;
+    k[2] = 2.f/3.f;
+    k[3] = 1.f;
+
+    l[0] = 0.f;
+    l[1] = 0.f;
+    l[2] = 1.f/3.f;
+    l[3] = 1.f;
+
+    m[0] = 0.f;
+    m[1] = 1.f/3.f;
+    m[2] = 2.f/3.f;
+    m[3] = 1.f;
+
+    // If d2 < 0 we need to flip the orientation of our curve
+    // This is done by negating the k and l values
+    if ( d[2] > 0) {
+        for (int i = 0; i < 4; ++i) {
+            k[i] = -k[i];
+            l[i] = -l[i];
+        }
+    }
+}
+
+// Calc coefficients of I(s,t) where roots of I are inflection points of curve
+// I(s,t) = t*(3*d0*s^2 - 3*d1*s*t + d2*t^2)
+// d0 = a1 - 2*a2+3*a3
+// d1 = -a2 + 3*a3
+// d2 = 3*a3
+// a1 = p0 . (p3 x p2)
+// a2 = p1 . (p0 x p3)
+// a3 = p2 . (p1 x p0)
+// Places the values of d1, d2, d3 in array d passed in
+static void calc_cubic_inflection_func(const SkPoint p[4], SkScalar d[3]) {
+    SkScalar a1 = calc_dot_cross_cubic(p[0], p[3], p[2]);
+    SkScalar a2 = calc_dot_cross_cubic(p[1], p[0], p[3]);
+    SkScalar a3 = calc_dot_cross_cubic(p[2], p[1], p[0]);
+
+    // need to scale a's or values in later calculations will grow to high
+    SkScalar max = SkScalarAbs(a1);
+    max = SkMaxScalar(max, SkScalarAbs(a2));
+    max = SkMaxScalar(max, SkScalarAbs(a3));
+    max = 1.f/max;
+    a1 = a1 * max;
+    a2 = a2 * max;
+    a3 = a3 * max;
+
+    d[2] = 3.f * a3;
+    d[1] = d[2] - a2;
+    d[0] = d[1] - a2 + a1;
+}
+
+int GrPathUtils::chopCubicAtLoopIntersection(const SkPoint src[4], SkPoint dst[10], SkScalar klm[9],
+                                             SkScalar klm_rev[3]) {
+    // Variable to store the two parametric values at the loop double point
+    SkScalar smallS = 0.f;
+    SkScalar largeS = 0.f;
+
+    SkScalar d[3];
+    calc_cubic_inflection_func(src, d);
+
+    CubicType cType = classify_cubic(src, d);
+
+    int chop_count = 0;
+    if (kLoop_CubicType == cType) {
+        SkScalar tempSqrt = SkScalarSqrt(4.f * d[0] * d[2] - 3.f * d[1] * d[1]);
+        SkScalar ls = d[1] - tempSqrt;
+        SkScalar lt = 2.f * d[0];
+        SkScalar ms = d[1] + tempSqrt;
+        SkScalar mt = 2.f * d[0];
+        ls = ls / lt;
+        ms = ms / mt;
+        // need to have t values sorted since this is what is expected by SkChopCubicAt
+        if (ls <= ms) {
+            smallS = ls;
+            largeS = ms;
+        } else {
+            smallS = ms;
+            largeS = ls;
+        }
+
+        SkScalar chop_ts[2];
+        if (smallS > 0.f && smallS < 1.f) {
+            chop_ts[chop_count++] = smallS;
+        }
+        if (largeS > 0.f && largeS < 1.f) {
+            chop_ts[chop_count++] = largeS;
+        }
+        if(dst) {
+            SkChopCubicAt(src, dst, chop_ts, chop_count);
+        }
+    } else {
+        if (dst) {
+            memcpy(dst, src, sizeof(SkPoint) * 4);
+        }
+    }
+
+    if (klm && klm_rev) {
+        // Set klm_rev to to match the sub_section of cubic that needs to have its orientation
+        // flipped. This will always be the section that is the "loop"
+        if (2 == chop_count) {
+            klm_rev[0] = 1.f;
+            klm_rev[1] = -1.f;
+            klm_rev[2] = 1.f;
+        } else if (1 == chop_count) {
+            if (smallS < 0.f) {
+                klm_rev[0] = -1.f;
+                klm_rev[1] = 1.f;
+            } else {
+                klm_rev[0] = 1.f;
+                klm_rev[1] = -1.f;
+            }
+        } else {
+            if (smallS < 0.f && largeS > 1.f) {
+                klm_rev[0] = -1.f;
+            } else {
+                klm_rev[0] = 1.f;
+            }
+        }
+        SkScalar controlK[4];
+        SkScalar controlL[4];
+        SkScalar controlM[4];
+
+        if (kSerpentine_CubicType == cType || (kCusp_CubicType == cType && 0.f != d[0])) {
+            set_serp_klm(d, controlK, controlL, controlM);
+        } else if (kLoop_CubicType == cType) {
+            set_loop_klm(d, controlK, controlL, controlM);
+        } else if (kCusp_CubicType == cType) {
+            SkASSERT(0.f == d[0]);
+            set_cusp_klm(d, controlK, controlL, controlM);
+        } else if (kQuadratic_CubicType == cType) {
+            set_quadratic_klm(d, controlK, controlL, controlM);
+        }
+
+        calc_cubic_klm(src, controlK, controlL, controlM, klm, &klm[3], &klm[6]);
+    }
+    return chop_count + 1;
+}
+
+void GrPathUtils::getCubicKLM(const SkPoint p[4], SkScalar klm[9]) {
+    SkScalar d[3];
+    calc_cubic_inflection_func(p, d);
+
+    CubicType cType = classify_cubic(p, d);
+
+    SkScalar controlK[4];
+    SkScalar controlL[4];
+    SkScalar controlM[4];
+
+    if (kSerpentine_CubicType == cType || (kCusp_CubicType == cType && 0.f != d[0])) {
+        set_serp_klm(d, controlK, controlL, controlM);
+    } else if (kLoop_CubicType == cType) {
+        set_loop_klm(d, controlK, controlL, controlM);
+    } else if (kCusp_CubicType == cType) {
+        SkASSERT(0.f == d[0]);
+        set_cusp_klm(d, controlK, controlL, controlM);
+    } else if (kQuadratic_CubicType == cType) {
+        set_quadratic_klm(d, controlK, controlL, controlM);
+    }
+
+    calc_cubic_klm(p, controlK, controlL, controlM, klm, &klm[3], &klm[6]);
+}
diff --git a/src/gpu/GrPathUtils.h b/src/gpu/GrPathUtils.h
index fc319ec50d..dc4ecd9ee3 100644
--- a/src/gpu/GrPathUtils.h
+++ b/src/gpu/GrPathUtils.h
@@ -115,5 +115,50 @@ namespace GrPathUtils {
                              bool constrainWithinTangents,
                              SkPath::Direction dir,
                              SkTArray<SkPoint, true>* quads);
+
+    // Chops the cubic bezier passed in by src, at the double point (intersection point)
+    // if the curve is a cubic loop. If it is a loop, there will be two parametric values for
+    // the double point: ls and ms. We chop the cubic at these values if they are between 0 and 1.
+    // Return value:
+    // Value of 3: ls and ms are both between (0,1), and dst will contain the three cubics,
+    //             dst[0..3], dst[3..6], and dst[6..9] if dst is not NULL
+    // Value of 2: Only one of ls and ms are between (0,1), and dst will contain the two cubics,
+    //             dst[0..3] and dst[3..6] if dst is not NULL
+    // Value of 1: Neither ls or ms are between (0,1), and dst will contain the one original cubic,
+    //             dst[0..3] if dst is not NULL
+    //
+    // Optional KLM Calculation:
+    // The function can also return the KLM linear functionals for the chopped cubic implicit form
+    // of K^3 - LM.
+    // It will calculate a single set of KLM values that can be shared by all sub cubics, except
+    // for the subsection that is "the loop" the K and L values need to be negated.
+    // Output:
+    // klm:     Holds the values for the linear functionals as:
+    //          K = (klm[0], klm[1], klm[2])
+    //          L = (klm[3], klm[4], klm[5])
+    //          M = (klm[6], klm[7], klm[8])
+    // klm_rev: These values are flags for the corresponding sub cubic saying whether or not
+    //          the K and L values need to be flipped. A value of -1.f means flip K and L and
+    //          a value of 1.f means do nothing.
+    //          *****DO NOT FLIP M, JUST K AND L*****
+    //
+    // Notice that the klm lines are calculated in the same space as the input control points.
+    // If you transform the points the lines will also need to be transformed. This can be done
+    // by mapping the lines with the inverse-transpose of the matrix used to map the points.
+    int chopCubicAtLoopIntersection(const SkPoint src[4], SkPoint dst[10] = NULL,
+                                    SkScalar klm[9] = NULL, SkScalar klm_rev[3] = NULL);
+
+    // Input is p which holds the 4 control points of a non-rational cubic Bezier curve.
+    // Output is the coefficients of the three linear functionals K, L, & M which
+    // represent the implicit form of the cubic as f(x,y,w) = K^3 - LM. The w term
+    // will always be 1. The output is stored in the array klm, where the values are:
+    // K = (klm[0], klm[1], klm[2])
+    // L = (klm[3], klm[4], klm[5])
+    // M = (klm[6], klm[7], klm[8])
+    //
+    // Notice that the klm lines are calculated in the same space as the input control points.
+    // If you transform the points the lines will also need to be transformed. This can be done
+    // by mapping the lines with the inverse-transpose of the matrix used to map the points.
+    void getCubicKLM(const SkPoint p[4], SkScalar klm[9]);
 };
 #endif