From 73ca6243b31e225e9fd5b75a96cbc82d62557de6 Mon Sep 17 00:00:00 2001 From: "caryclark@google.com" Date: Thu, 17 Jan 2013 21:02:47 +0000 Subject: [PATCH] shape ops work in progress mostly working on cubic/cubic intersect git-svn-id: http://skia.googlecode.com/svn/trunk@7266 2bbb7eff-a529-9590-31e7-b0007b416f81 --- experimental/Intersection/CubicBezierClip.cpp | 4 +- experimental/Intersection/CubicIntersection.cpp | 130 ++ .../Intersection/CubicIntersection_Test.cpp | 209 ++ .../Intersection/CubicIntersection_TestData.cpp | 1 - .../Intersection/CubicIntersection_TestData.h | 5 +- experimental/Intersection/CubicReduceOrder.cpp | 19 +- experimental/Intersection/CubicToQuadratics.cpp | 128 +- .../Intersection/CubicToQuadratics_Test.cpp | 117 +- experimental/Intersection/CubicUtilities.cpp | 28 +- experimental/Intersection/CubicUtilities.h | 7 +- experimental/Intersection/CurveIntersection.h | 6 +- experimental/Intersection/DataTypes.cpp | 7 - experimental/Intersection/DataTypes.h | 26 +- experimental/Intersection/DataTypes_Test.h | 7 + experimental/Intersection/EdgeWalker.cpp | 2 +- experimental/Intersection/Intersection_Tests.cpp | 6 +- experimental/Intersection/Intersection_Tests.h | 6 +- experimental/Intersection/Intersections.cpp | 99 + experimental/Intersection/Intersections.h | 129 +- .../Intersection/LineCubicIntersection.cpp | 277 +-- .../Intersection/LineCubicIntersection_Test.cpp | 6 +- experimental/Intersection/LineIntersection.cpp | 20 +- experimental/Intersection/LineIntersection.h | 5 +- experimental/Intersection/LineParameters.h | 12 +- .../Intersection/LineParameteters_Test.cpp | 6 +- experimental/Intersection/QuadraticImplicit.cpp | 345 +++- .../Intersection/QuadraticIntersection_Test.cpp | 51 +- .../Intersection/QuadraticIntersection_TestData.h | 5 +- experimental/Intersection/QuadraticUtilities.h | 8 +- experimental/Intersection/QuarticRoot.cpp | 60 +- experimental/Intersection/QuarticRoot_Test.cpp | 15 +- experimental/Intersection/Simplify.cpp | 202 +- experimental/Intersection/SimplifyNew_Test.cpp | 44 +- experimental/Intersection/TSearch.h | 29 + experimental/Intersection/op.htm | 29 +- experimental/Intersection/qc.htm | 2079 ++++++++++++++++++++ 36 files changed, 3674 insertions(+), 455 deletions(-) create mode 100644 experimental/Intersection/DataTypes_Test.h create mode 100644 experimental/Intersection/qc.htm diff --git a/experimental/Intersection/CubicBezierClip.cpp b/experimental/Intersection/CubicBezierClip.cpp index d0141e9..825a7b9 100644 --- a/experimental/Intersection/CubicBezierClip.cpp +++ b/experimental/Intersection/CubicBezierClip.cpp @@ -26,7 +26,8 @@ bool bezier_clip(const Cubic& cubic1, const Cubic& cubic2, double& minT, double& } double distance[2]; - endLine.controlPtDistance(cubic1, distance); + distance[0] = endLine.controlPtDistance(cubic1, 1); + distance[1] = endLine.controlPtDistance(cubic1, 2); // find fat line double top = distance[0]; @@ -87,4 +88,3 @@ bool bezier_clip(const Cubic& cubic1, const Cubic& cubic2, double& minT, double& return minT < maxT; // returns false if distance shows no intersection } - diff --git a/experimental/Intersection/CubicIntersection.cpp b/experimental/Intersection/CubicIntersection.cpp index fcd4119..d5a5bde 100644 --- a/experimental/Intersection/CubicIntersection.cpp +++ b/experimental/Intersection/CubicIntersection.cpp @@ -159,3 +159,133 @@ bool intersect(const Cubic& c1, const Cubic& c2, Intersections& i) { return c.intersect(); } +#include "CubicUtilities.h" + +// this flavor approximates the cubics with quads to find the intersecting ts +// OPTIMIZE: if this strategy proves successful, the quad approximations, or the ts used +// to create the approximations, could be stored in the cubic segment +// fixme: this strategy needs to add short line segments on either end, or similarly extend the +// initial and final quadratics +bool intersect2(const Cubic& c1, const Cubic& c2, Intersections& i) { + SkTDArray ts1; + double precision1 = calcPrecision(c1); + cubic_to_quadratics(c1, precision1, ts1); + SkTDArray ts2; + double precision2 = calcPrecision(c2); + cubic_to_quadratics(c2, precision2, ts2); + double t1Start = 0; + int ts1Count = ts1.count(); + for (int i1 = 0; i1 <= ts1Count; ++i1) { + const double t1 = i1 < ts1Count ? ts1[i1] : 1; + Cubic part1; + sub_divide(c1, t1Start, t1, part1); + Quadratic q1; + demote_cubic_to_quad(part1, q1); + // start here; + // should reduceOrder be looser in this use case if quartic is going to blow up on an + // extremely shallow quadratic? + // maybe quadratics to lines need the same sort of recursive solution that I hope to find + // for cubics to quadratics ... + Quadratic s1; + int o1 = reduceOrder(q1, s1); + double t2Start = 0; + int ts2Count = ts2.count(); + for (int i2 = 0; i2 <= ts2Count; ++i2) { + const double t2 = i2 < ts2Count ? ts2[i2] : 1; + Cubic part2; + sub_divide(c2, t2Start, t2, part2); + Quadratic q2; + demote_cubic_to_quad(part2, q2); + Quadratic s2; + double o2 = reduceOrder(q2, s2); + Intersections locals; + if (o1 == 3 && o2 == 3) { + intersect2(q1, q2, locals); + } else if (o1 <= 2 && o2 <= 2) { + i.fUsed = intersect((const _Line&) s1, (const _Line&) s2, i.fT[0], i.fT[1]); + } else if (o1 == 3 && o2 <= 2) { + intersect(q1, (const _Line&) s2, i); + } else { + SkASSERT(o1 <= 2 && o2 == 3); + intersect(q2, (const _Line&) s1, i); + for (int s = 0; s < i.fUsed; ++s) { + SkTSwap(i.fT[0][s], i.fT[1][s]); + } + } + for (int tIdx = 0; tIdx < locals.used(); ++tIdx) { + double to1 = t1Start + (t1 - t1Start) * locals.fT[0][tIdx]; + double to2 = t2Start + (t2 - t2Start) * locals.fT[1][tIdx]; + i.insert(to1, to2); + } + t2Start = t2; + } + t1Start = t1; + } + return i.intersected(); +} + +int intersect(const Cubic& cubic, const Quadratic& quad, Intersections& i) { + SkTDArray ts; + double precision = calcPrecision(cubic); + cubic_to_quadratics(cubic, precision, ts); + double tStart = 0; + Cubic part; + int tsCount = ts.count(); + for (int idx = 0; idx <= tsCount; ++idx) { + double t = idx < tsCount ? ts[idx] : 1; + Quadratic q1; + sub_divide(cubic, tStart, t, part); + demote_cubic_to_quad(part, q1); + Intersections locals; + intersect2(q1, quad, locals); + for (int tIdx = 0; tIdx < locals.used(); ++tIdx) { + double globalT = tStart + (t - tStart) * locals.fT[0][tIdx]; + i.insertOne(globalT, 0); + globalT = locals.fT[1][tIdx]; + i.insertOne(globalT, 1); + } + tStart = t; + } + return i.used(); +} + +bool intersect(const Cubic& cubic, Intersections& i) { + SkTDArray ts; + double precision = calcPrecision(cubic); + cubic_to_quadratics(cubic, precision, ts); + int tsCount = ts.count(); + if (tsCount == 1) { + return false; + } + double t1Start = 0; + Cubic part; + for (int idx = 0; idx < tsCount; ++idx) { + double t1 = ts[idx]; + Quadratic q1; + sub_divide(cubic, t1Start, t1, part); + demote_cubic_to_quad(part, q1); + double t2Start = t1; + for (int i2 = idx + 1; i2 <= tsCount; ++i2) { + const double t2 = i2 < tsCount ? ts[i2] : 1; + Quadratic q2; + sub_divide(cubic, t2Start, t2, part); + demote_cubic_to_quad(part, q2); + Intersections locals; + intersect2(q1, q2, locals); + for (int tIdx = 0; tIdx < locals.used(); ++tIdx) { + // discard intersections at cusp? (maximum curvature) + double t1sect = locals.fT[0][tIdx]; + double t2sect = locals.fT[1][tIdx]; + if (idx + 1 == i2 && t1sect == 1 && t2sect == 0) { + continue; + } + double to1 = t1Start + (t1 - t1Start) * t1sect; + double to2 = t2Start + (t2 - t2Start) * t2sect; + i.insert(to1, to2); + } + t2Start = t2; + } + t1Start = t1; + } + return i.intersected(); +} diff --git a/experimental/Intersection/CubicIntersection_Test.cpp b/experimental/Intersection/CubicIntersection_Test.cpp index 256529c..7d27bf4 100644 --- a/experimental/Intersection/CubicIntersection_Test.cpp +++ b/experimental/Intersection/CubicIntersection_Test.cpp @@ -56,3 +56,212 @@ void CubicIntersection_Test() { } } } + +static void oneOff(const Cubic& cubic1, const Cubic& cubic2) { + SkTDArray quads1; + cubic_to_quadratics(cubic1, calcPrecision(cubic1), quads1); + for (int index = 0; index < quads1.count(); ++index) { + const Quadratic& q = quads1[index]; + SkDebugf("{{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", q[0].x, q[0].y, + q[1].x, q[1].y, q[2].x, q[2].y); + } + SkDebugf("\n"); + SkTDArray quads2; + cubic_to_quadratics(cubic2, calcPrecision(cubic2), quads2); + for (int index = 0; index < quads2.count(); ++index) { + const Quadratic& q = quads2[index]; + SkDebugf("{{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", q[0].x, q[0].y, + q[1].x, q[1].y, q[2].x, q[2].y); + } + SkDebugf("\n"); + Intersections intersections2; + intersect2(cubic1, cubic2, intersections2); + for (int pt = 0; pt < intersections2.used(); ++pt) { + double tt1 = intersections2.fT[0][pt]; + double tx1, ty1; + xy_at_t(cubic1, tt1, tx1, ty1); + int pt2 = intersections2.fFlip ? intersections2.used() - pt - 1 : pt; + double tt2 = intersections2.fT[1][pt2]; + double tx2, ty2; + xy_at_t(cubic2, tt2, tx2, ty2); + SkDebugf("%s t1=%1.9g (%1.9g, %1.9g) (%1.9g, %1.9g) t2=%1.9g\n", __FUNCTION__, + tt1, tx1, ty1, tx2, ty2, tt2); + } +} + +static const Cubic testSet[] = { +{{67.426548091427676, 37.993772624988935}, {23.483695892376684, 90.476863174921306}, {35.597065061143162, 79.872482633158796}, {75.38634169631932, 18.244890038969412}}, +{{61.336508189019057, 82.693132843213675}, {44.639380902349664, 54.074825790745592}, {16.815615499771951, 20.049704667203923}, {41.866884958868326, 56.735503699973002}}, + +{{67.4265481, 37.9937726}, {23.4836959, 90.4768632}, {35.5970651, 79.8724826}, {75.3863417, 18.24489}}, +{{61.3365082, 82.6931328}, {44.6393809, 54.0748258}, {16.8156155, 20.0497047}, {41.866885, 56.7355037}}, + +{{18.1312339, 31.6473732}, {95.5711034, 63.5350219}, {92.3283165, 62.0158945}, {18.5656052, 32.1268808}}, +{{97.402018, 35.7169972}, {33.1127443, 25.8935163}, {1.13970027, 54.9424981}, {56.4860195, 60.529264}}, +}; + +const size_t testSetCount = sizeof(testSet) / sizeof(testSet[0]); + +void CubicIntersection_OneOffTest() { + for (size_t outer = 0; outer < testSetCount - 1; ++outer) { + SkDebugf("%s quads1[%d]\n", __FUNCTION__, outer); + const Cubic& cubic1 = testSet[outer]; + for (size_t inner = outer + 1; inner < testSetCount; ++inner) { + SkDebugf("%s quads2[%d]\n", __FUNCTION__, inner); + const Cubic& cubic2 = testSet[inner]; + oneOff(cubic1, cubic2); + } + } +} + +#define DEBUG_CRASH 1 + +class CubicChopper { +public: + +// only finds one intersection +CubicChopper(const Cubic& c1, const Cubic& c2) + : cubic1(c1) + , cubic2(c2) + , depth(0) { +} + +bool intersect(double minT1, double maxT1, double minT2, double maxT2) { + Cubic sub1, sub2; + // FIXME: carry last subdivide and reduceOrder result with cubic + sub_divide(cubic1, minT1, maxT1, sub1); + sub_divide(cubic2, minT2, maxT2, sub2); + Intersections i; + intersect2(sub1, sub2, i); + if (i.used() == 0) { + return false; + } + double x1, y1, x2, y2; + t1 = minT1 + i.fT[0][0] * (maxT1 - minT1); + t2 = minT2 + i.fT[1][0] * (maxT2 - minT2); + xy_at_t(cubic1, t1, x1, y1); + xy_at_t(cubic2, t2, x2, y2); + if (AlmostEqualUlps(x1, x2) && AlmostEqualUlps(y1, y2)) { + return true; + } + double half1 = (minT1 + maxT1) / 2; + double half2 = (minT2 + maxT2) / 2; + ++depth; + bool result; + if (depth & 1) { + result = intersect(minT1, half1, minT2, maxT2) || intersect(half1, maxT1, minT2, maxT2) + || intersect(minT1, maxT1, minT2, half2) || intersect(minT1, maxT1, half2, maxT2); + } else { + result = intersect(minT1, maxT1, minT2, half2) || intersect(minT1, maxT1, half2, maxT2) + || intersect(minT1, half1, minT2, maxT2) || intersect(half1, maxT1, minT2, maxT2); + } + --depth; + return result; +} + +const Cubic& cubic1; +const Cubic& cubic2; +double t1; +double t2; +int depth; +}; + +#define TRY_OLD 0 // old way fails on test == 1 + +void CubicIntersection_RandTest() { + srand(0); + const int tests = 1000000; // 10000000; + double largestFactor = DBL_MAX; + for (int test = 0; test < tests; ++test) { + Cubic cubic1, cubic2; + for (int i = 0; i < 4; ++i) { + cubic1[i].x = (double) rand() / RAND_MAX * 100; + cubic1[i].y = (double) rand() / RAND_MAX * 100; + cubic2[i].x = (double) rand() / RAND_MAX * 100; + cubic2[i].y = (double) rand() / RAND_MAX * 100; + } + if (test == 2513) { // the pair crosses three times, but the quadratic approximation + continue; // only sees one -- should be OK to ignore the other two? + } + if (test == 12932) { // this exposes a weakness when one cubic touches the other but + continue; // does not touch the quad approximation. Captured in qc.htm as cubic15 + } + #if DEBUG_CRASH + char str[1024]; + sprintf(str, "{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n" + "{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n", + cubic1[0].x, cubic1[0].y, cubic1[1].x, cubic1[1].y, cubic1[2].x, cubic1[2].y, + cubic1[3].x, cubic1[3].y, + cubic2[0].x, cubic2[0].y, cubic2[1].x, cubic2[1].y, cubic2[2].x, cubic2[2].y, + cubic2[3].x, cubic2[3].y); + #endif + _Rect rect1, rect2; + rect1.setBounds(cubic1); + rect2.setBounds(cubic2); + bool boundsIntersect = rect1.left <= rect2.right && rect2.left <= rect2.right + && rect1.top <= rect2.bottom && rect2.top <= rect1.bottom; + Intersections i1, i2; + #if TRY_OLD + bool oldIntersects = intersect(cubic1, cubic2, i1); + #else + bool oldIntersects = false; + #endif + if (test == -1) { + SkDebugf("ready...\n"); + } + bool newIntersects = intersect2(cubic1, cubic2, i2); + if (!boundsIntersect && (oldIntersects || newIntersects)) { + SkDebugf("%s %d unexpected intersection boundsIntersect=%d oldIntersects=%d" + " newIntersects=%d\n%s %s\n", __FUNCTION__, test, boundsIntersect, + oldIntersects, newIntersects, __FUNCTION__, str); + assert(0); + } + if (oldIntersects && !newIntersects) { + SkDebugf("%s %d missing intersection oldIntersects=%d newIntersects=%d\n%s %s\n", + __FUNCTION__, test, oldIntersects, newIntersects, __FUNCTION__, str); + assert(0); + } + if (!oldIntersects && !newIntersects) { + continue; + } + if (i2.used() > 1) { + continue; + // just look at single intercepts for simplicity + } + Intersections self1, self2; // self-intersect checks + if (intersect(cubic1, self1)) { + continue; + } + if (intersect(cubic2, self2)) { + continue; + } + // binary search for range necessary to enclose real intersection + CubicChopper c(cubic1, cubic2); + bool result = c.intersect(0, 1, 0, 1); + if (!result) { + // FIXME: a failure here probably means that a core routine used by CubicChopper is failing + continue; + } + double delta1 = fabs(c.t1 - i2.fT[0][0]); + double delta2 = fabs(c.t2 - i2.fT[1][0]); + double calc1 = calcPrecision(cubic1); + double calc2 = calcPrecision(cubic2); + double factor1 = calc1 / delta1; + double factor2 = calc2 / delta2; + SkDebugf("%s %d calc1=%1.9g delta1=%1.9g factor1=%1.9g calc2=%1.9g delta2=%1.9g" + " factor2=%1.9g\n", __FUNCTION__, test, + calc1, delta1, factor1, calc2, delta2, factor2); + if (factor1 < largestFactor) { + SkDebugf("WE HAVE A WINNER! %1.9g\n", factor1); + SkDebugf("%s\n", str); + oneOff(cubic1, cubic2); + largestFactor = factor1; + } + if (factor2 < largestFactor) { + SkDebugf("WE HAVE A WINNER! %1.9g\n", factor2); + SkDebugf("%s\n", str); + oneOff(cubic1, cubic2); + largestFactor = factor2; + } + } +} diff --git a/experimental/Intersection/CubicIntersection_TestData.cpp b/experimental/Intersection/CubicIntersection_TestData.cpp index 645d7b0..c138a67 100644 --- a/experimental/Intersection/CubicIntersection_TestData.cpp +++ b/experimental/Intersection/CubicIntersection_TestData.cpp @@ -4,7 +4,6 @@ * Use of this source code is governed by a BSD-style license that can be * found in the LICENSE file. */ -#define IN_TEST 1 #include "CubicIntersection_TestData.h" #include diff --git a/experimental/Intersection/CubicIntersection_TestData.h b/experimental/Intersection/CubicIntersection_TestData.h index 27436f8..9893b3c 100644 --- a/experimental/Intersection/CubicIntersection_TestData.h +++ b/experimental/Intersection/CubicIntersection_TestData.h @@ -4,11 +4,8 @@ * Use of this source code is governed by a BSD-style license that can be * found in the LICENSE file. */ -#if !defined(IN_TEST) - #define IN_TEST 1 -#endif - #include "DataTypes.h" +#include "DataTypes_Test.h" extern const Cubic pointDegenerates[]; extern const Cubic notPointDegenerates[]; diff --git a/experimental/Intersection/CubicReduceOrder.cpp b/experimental/Intersection/CubicReduceOrder.cpp index d1775d8..fdc7265 100644 --- a/experimental/Intersection/CubicReduceOrder.cpp +++ b/experimental/Intersection/CubicReduceOrder.cpp @@ -149,22 +149,17 @@ static int check_linear(const Cubic& cubic, Cubic& reduction, bool isLinear(const Cubic& cubic, int startIndex, int endIndex) { LineParameters lineParameters; lineParameters.cubicEndPoints(cubic, startIndex, endIndex); - double normalSquared = lineParameters.normalSquared(); - double distance[2]; // distance is not normalized + // FIXME: maybe it's possible to avoid this and compare non-normalized + lineParameters.normalize(); int mask = other_two(startIndex, endIndex); int inner1 = startIndex ^ mask; int inner2 = endIndex ^ mask; - lineParameters.controlPtDistance(cubic, inner1, inner2, distance); - double limit = normalSquared; - int index; - for (index = 0; index < 2; ++index) { - double distSq = distance[index]; - distSq *= distSq; - if (approximately_greater(distSq, limit)) { - return false; - } + double distance = lineParameters.controlPtDistance(cubic, inner1); + if (!approximately_zero(distance)) { + return false; } - return true; + distance = lineParameters.controlPtDistance(cubic, inner2); + return approximately_zero(distance); } /* food for thought: diff --git a/experimental/Intersection/CubicToQuadratics.cpp b/experimental/Intersection/CubicToQuadratics.cpp index 424836c..ff6d7a2 100644 --- a/experimental/Intersection/CubicToQuadratics.cpp +++ b/experimental/Intersection/CubicToQuadratics.cpp @@ -40,11 +40,18 @@ Tdiv>=1 - the entire cubic can be approximated by the mid-point approximation confirmed by (maybe stolen from) http://www.caffeineowl.com/graphics/2d/vectorial/cubic2quad01.html +// maybe in turn derived from http://www.cccg.ca/proceedings/2004/36.pdf +// also stored at http://www.cis.usouthal.edu/~hain/general/Publications/Bezier/bezier%20cccg04%20paper.pdf */ #include "CubicUtilities.h" #include "CurveIntersection.h" +#include "LineIntersection.h" + +const bool AVERAGE_END_POINTS = true; // results in better fitting curves + +#define USE_CUBIC_END_POINTS 1 static double calcTDiv(const Cubic& cubic, double precision, double start) { const double adjust = sqrt(3) / 36; @@ -62,54 +69,105 @@ static double calcTDiv(const Cubic& cubic, double precision, double start) { double dy = c[3].y - 3 * (c[2].y - c[1].y) - c[0].y; double dist = sqrt(dx * dx + dy * dy); double tDiv3 = precision / (adjust * dist); - return cube_root(tDiv3); + double t = cube_root(tDiv3); + if (start > 0) { + t = start + (1 - start) * t; + } + return t; } -static void demote(const Cubic& cubic, Quadratic& quad) { +void demote_cubic_to_quad(const Cubic& cubic, Quadratic& quad) { quad[0] = cubic[0]; - quad[1].x = (cubic[3].x - 3 * (cubic[2].x - cubic[1].x) - cubic[0].x) / 4; - quad[1].y = (cubic[3].y - 3 * (cubic[2].y - cubic[1].y) - cubic[0].y) / 4; +if (AVERAGE_END_POINTS) { + const _Point fromC1 = { (3 * cubic[1].x - cubic[0].x) / 2, (3 * cubic[1].y - cubic[0].y) / 2 }; + const _Point fromC2 = { (3 * cubic[2].x - cubic[3].x) / 2, (3 * cubic[2].y - cubic[3].y) / 2 }; + quad[1].x = (fromC1.x + fromC2.x) / 2; + quad[1].y = (fromC1.y + fromC2.y) / 2; +} else { + lineIntersect((const _Line&) cubic[0], (const _Line&) cubic[2], quad[1]); +} quad[2] = cubic[3]; } int cubic_to_quadratics(const Cubic& cubic, double precision, SkTDArray& quadratics) { - quadratics.setCount(1); // FIXME: every place I have setCount(), I also want setReserve() + SkTDArray ts; + cubic_to_quadratics(cubic, precision, ts); + int tsCount = ts.count(); + double t1Start = 0; + int order = 0; + for (int idx = 0; idx <= tsCount; ++idx) { + double t1 = idx < tsCount ? ts[idx] : 1; + Cubic part; + sub_divide(cubic, t1Start, t1, part); + Quadratic q1; + demote_cubic_to_quad(part, q1); + Quadratic s1; + int o1 = reduceOrder(q1, s1); + if (order < o1) { + order = o1; + } + memcpy(quadratics.append(), o1 < 2 ? s1 : q1, sizeof(Quadratic)); + t1Start = t1; + } + return order; +} + +static bool addSimpleTs(const Cubic& cubic, double precision, SkTDArray& ts) { + double tDiv = calcTDiv(cubic, precision, 0); + if (tDiv >= 1) { + return true; + } + if (tDiv >= 0.5) { + *ts.append() = 0.5; + return true; + } + return false; +} + +static void addTs(const Cubic& cubic, double precision, double start, double end, + SkTDArray& ts) { + double tDiv = calcTDiv(cubic, precision, 0); + double parts = ceil(1.0 / tDiv); + for (double index = 0; index < parts; ++index) { + double newT = start + (index / parts) * (end - start); + if (newT > 0 && newT < 1) { + *ts.append() = newT; + } + } +} + +// flavor that returns T values only, deferring computing the quads until they are needed +void cubic_to_quadratics(const Cubic& cubic, double precision, SkTDArray& ts) { Cubic reduced; int order = reduceOrder(cubic, reduced, kReduceOrder_QuadraticsAllowed); if (order < 3) { - memcpy(quadratics[0], reduced, order * sizeof(_Point)); - return order; + return; } - double tDiv = calcTDiv(cubic, precision, 0); - if (approximately_greater_than_one(tDiv)) { - demote(cubic, quadratics[0]); - return 3; + double inflectT[2]; + int inflections = find_cubic_inflections(cubic, inflectT); + SkASSERT(inflections <= 2); + if (inflections == 0 && addSimpleTs(cubic, precision, ts)) { + return; } - if (tDiv >= 0.5) { + if (inflections == 1) { CubicPair pair; - chop_at(cubic, pair, 0.5); - quadratics.setCount(2); - demote(pair.first(), quadratics[0]); - demote(pair.second(), quadratics[1]); - return 3; + chop_at(cubic, pair, inflectT[0]); + addTs(pair.first(), precision, 0, inflectT[0], ts); + addTs(pair.second(), precision, inflectT[0], 1, ts); + return; } - double start = 0; - int index = -1; - // if we iterate but make little progress, check to see if the curve is flat - // or if the control points are tangent to each other - do { - Cubic part; - sub_divide(cubic, start, tDiv, part); - quadratics.append(); - demote(part, quadratics[++index]); - if (tDiv == 1) { - break; + if (inflections == 2) { + if (inflectT[0] > inflectT[1]) { + SkTSwap(inflectT[0], inflectT[1]); } - start = tDiv; - tDiv = calcTDiv(cubic, precision, start); - if (tDiv > 1) { - tDiv = 1; - } - } while (true); - return 3; + Cubic part; + sub_divide(cubic, 0, inflectT[0], part); + addTs(part, precision, 0, inflectT[0], ts); + sub_divide(cubic, inflectT[0], inflectT[1], part); + addTs(part, precision, inflectT[0], inflectT[1], ts); + sub_divide(cubic, inflectT[1], 1, part); + addTs(part, precision, inflectT[1], 1, ts); + return; + } + addTs(cubic, precision, 0, 1, ts); } diff --git a/experimental/Intersection/CubicToQuadratics_Test.cpp b/experimental/Intersection/CubicToQuadratics_Test.cpp index 24c0ffe..843833f 100644 --- a/experimental/Intersection/CubicToQuadratics_Test.cpp +++ b/experimental/Intersection/CubicToQuadratics_Test.cpp @@ -3,21 +3,14 @@ #include "Intersection_Tests.h" #include "QuadraticIntersection_TestData.h" #include "TestUtilities.h" +#include "SkGeometry.h" -static double calcPrecision(const Cubic& cubic) { - _Rect dRect; - dRect.setBounds(cubic); - double width = dRect.right - dRect.left; - double height = dRect.bottom - dRect.top; - return (width > height ? width : height) / 256; -} - -static void test(const Cubic(& cubics)[], const char* name, int firstTest, size_t testCount) { +static void test(const Cubic* cubics, const char* name, int firstTest, size_t testCount) { SkTDArray quads; for (size_t index = firstTest; index < testCount; ++index) { const Cubic& cubic = cubics[index]; double precision = calcPrecision(cubic); - cubic_to_quadratics(cubic, precision, quads); + (void) cubic_to_quadratics(cubic, precision, quads); if (quads.count() != 1) { printf("%s [%d] cubic to quadratics failed count=%d\n", name, (int) index, quads.count()); @@ -25,14 +18,14 @@ static void test(const Cubic(& cubics)[], const char* name, int firstTest, size_ } } -static void test(const Quadratic(& quadTests)[], const char* name, int firstTest, size_t testCount) { +static void test(const Quadratic* quadTests, const char* name, int firstTest, size_t testCount) { SkTDArray quads; for (size_t index = firstTest; index < testCount; ++index) { const Quadratic& quad = quadTests[index]; Cubic cubic; quad_to_cubic(quad, cubic); double precision = calcPrecision(cubic); - cubic_to_quadratics(cubic, precision, quads); + (void) cubic_to_quadratics(cubic, precision, quads); if (quads.count() != 1) { printf("%s [%d] cubic to quadratics failed count=%d\n", name, (int) index, quads.count()); @@ -40,7 +33,7 @@ static void test(const Quadratic(& quadTests)[], const char* name, int firstTest } } -static void testC(const Cubic(& cubics)[], const char* name, int firstTest, size_t testCount) { +static void testC(const Cubic* cubics, const char* name, int firstTest, size_t testCount) { SkTDArray quads; // test if computed line end points are valid for (size_t index = firstTest; index < testCount; ++index) { @@ -63,15 +56,15 @@ static void testC(const Cubic(& cubics)[], const char* name, int firstTest, size } } -static void testC(const Cubic(& cubics)[][2], const char* name, int firstTest, size_t testCount) { +static void testC(const Cubic(* cubics)[2], const char* name, int firstTest, size_t testCount) { SkTDArray quads; for (size_t index = firstTest; index < testCount; ++index) { for (int idx2 = 0; idx2 < 2; ++idx2) { const Cubic& cubic = cubics[index][idx2]; double precision = calcPrecision(cubic); int order = cubic_to_quadratics(cubic, precision, quads); - assert(order != 4); - if (order < 3) { + assert(order != 4); + if (order < 3) { continue; } if (!AlmostEqualUlps(cubic[0].x, quads[0][0].x) @@ -136,6 +129,8 @@ void CubicToQuadratics_Test() { } static Cubic locals[] = { + {{14.5975863, 41.632436}, {16.3518929, 26.2639684}, {18.5165519, 7.68775139}, {8.03767257, 89.1628526}}, + {{69.7292014, 38.6877352}, {24.7648688, 23.1501713}, {84.9283191, 90.2588441}, {80.392774, 61.3533852}}, {{ 60.776536520932126, 71.249307306133829 @@ -148,39 +143,113 @@ static Cubic locals[] = { }, { 45.261946574441133, 17.536076632112298 - }} + }}, }; static size_t localsCount = sizeof(locals) / sizeof(locals[0]); +#define DEBUG_CRASH 0 +#define TEST_AVERAGE_END_POINTS 0 // must take const off to test +extern const bool AVERAGE_END_POINTS; + void CubicsToQuadratics_RandTest() { for (size_t x = 0; x < localsCount; ++x) { const Cubic& cubic = locals[x]; + const SkPoint skcubic[4] = {{(float) cubic[0].x, (float) cubic[0].y}, + {(float) cubic[1].x, (float) cubic[1].y}, {(float) cubic[2].x, (float) cubic[2].y}, + {(float) cubic[3].x, (float) cubic[3].y}}; + SkScalar skinflect[2]; + int skin = SkFindCubicInflections(skcubic, skinflect); + SkDebugf("%s %d %1.9g\n", __FUNCTION__, skin, skinflect[0]); SkTDArray quads; double precision = calcPrecision(cubic); - cubic_to_quadratics(cubic, precision, quads); + (void) cubic_to_quadratics(cubic, precision, quads); + SkDebugf("%s quads=%d\n", __FUNCTION__, quads.count()); } srand(0); - const int arrayMax = 1000; - const int tests = 1000000; + const int arrayMax = 8; + const int sampleMax = 10; + const int tests = 1000000; // 10000000; int quadDist[arrayMax]; bzero(quadDist, sizeof(quadDist)); + Cubic samples[arrayMax][sampleMax]; + int sampleCount[arrayMax]; + bzero(sampleCount, sizeof(sampleCount)); for (int x = 0; x < tests; ++x) { Cubic cubic; for (int i = 0; i < 4; ++i) { cubic[i].x = (double) rand() / RAND_MAX * 100; cubic[i].y = (double) rand() / RAND_MAX * 100; } + #if DEBUG_CRASH + char str[1024]; + sprintf(str, "{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n", + cubic[0].x, cubic[0].y, cubic[1].x, cubic[1].y, cubic[2].x, cubic[2].y, + cubic[3].x, cubic[3].y); + #endif SkTDArray quads; double precision = calcPrecision(cubic); - cubic_to_quadratics(cubic, precision, quads); - assert(quads.count() < arrayMax); - quadDist[quads.count()]++; + (void) cubic_to_quadratics(cubic, precision, quads); + int count = quads.count(); + assert(count > 0); + assert(--count < arrayMax); + quadDist[count]++; + int sCount = sampleCount[count]; + if (sCount < sampleMax) { + memcpy(samples[count][sCount], cubic, sizeof(Cubic)); + sampleCount[count]++; + } } for (int x = 0; x < arrayMax; ++x) { if (!quadDist[x]) { continue; } - printf("%d 1.9%g%%\n", x, (double) quadDist[x] / tests * 100); + SkDebugf("%d %1.9g%%\n", x + 1, (double) quadDist[x] / tests * 100); + } + SkDebugf("\n"); + for (int x = 0; x < arrayMax; ++x) { + for (int y = 0; y < sampleCount[x]; ++y) { +#if TEST_AVERAGE_END_POINTS + for (int w = 0; w < 2; ++w) { + AVERAGE_END_POINTS = w; +#else + int w = 0; +#endif + SkDebugf("
\n", x + 1, y, w ? "x" : ""); + const Cubic& cubic = samples[x][y]; + SkDebugf("{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n", + cubic[0].x, cubic[0].y, cubic[1].x, cubic[1].y, cubic[2].x, cubic[2].y, + cubic[3].x, cubic[3].y); + SkTDArray quads; + double precision = calcPrecision(cubic); + (void) cubic_to_quadratics(cubic, precision, quads); + for (int z = 0; z < quads.count(); ++z) { + const Quadratic& quad = quads[z]; + SkDebugf("{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n", + quad[0].x, quad[0].y, quad[1].x, quad[1].y, quad[2].x, quad[2].y); + } + SkDebugf("
\n\n"); +#if TEST_AVERAGE_END_POINTS + } +#endif + } + } + SkDebugf("\n\n"); + SkDebugf(" + + + + + + -- 2.7.4