src/third_party/skia/tests/PathOpsAngleIdeas.cpp

   1 /*
   2  * Copyright 2013 Google Inc.
   3  *
   4  * Use of this source code is governed by a BSD-style license that can be
   5  * found in the LICENSE file.
   6  */
   7 #include "PathOpsTestCommon.h"
   8 #include "SkIntersections.h"
   9 #include "SkOpSegment.h"
  10 #include "SkPathOpsTriangle.h"
  11 #include "SkRandom.h"
  12 #include "SkTArray.h"
  13 #include "SkTSort.h"
  14 #include "Test.h"
  15
  16 static bool gPathOpsAngleIdeasVerbose = false;
  17 static bool gPathOpsAngleIdeasEnableBruteCheck = false;
  18
  19 class PathOpsAngleTester {
  20 public:
  21     static int ConvexHullOverlaps(const SkOpAngle& lh, const SkOpAngle& rh) {
  22         return lh.convexHullOverlaps(rh);
  23     }
  24
  25     static int EndsIntersect(const SkOpAngle& lh, const SkOpAngle& rh) {
  26         return lh.endsIntersect(rh);
  27     }
  28 };
  29
  30 struct TRange {
  31     double tMin1;
  32     double tMin2;
  33     double t1;
  34     double t2;
  35     double tMin;
  36     double a1;
  37     double a2;
  38     bool ccw;
  39 };
  40
  41 static double testArc(skiatest::Reporter* reporter, const SkDQuad& quad, const SkDQuad& arcRef,
  42         int octant) {
  43     SkDQuad arc = arcRef;
  44     SkDVector offset = {quad[0].fX, quad[0].fY};
  45     arc[0] += offset;
  46     arc[1] += offset;
  47     arc[2] += offset;
  48     SkIntersections i;
  49     i.intersect(arc, quad);
  50     if (i.used() == 0) {
  51         return -1;
  52     }
  53     int smallest = -1;
  54     double t = 2;
  55     for (int idx = 0; idx < i.used(); ++idx) {
  56         if (i[0][idx] > 1 || i[0][idx] < 0) {
  57             i.reset();
  58             i.intersect(arc, quad);
  59         }
  60         if (i[1][idx] > 1 || i[1][idx] < 0) {
  61             i.reset();
  62             i.intersect(arc, quad);
  63         }
  64         if (t > i[1][idx]) {
  65             smallest = idx;
  66             t = i[1][idx];
  67         }
  68     }
  69     REPORTER_ASSERT(reporter, smallest >= 0);
  70     REPORTER_ASSERT(reporter, t >= 0 && t <= 1);
  71     return i[1][smallest];
  72 }
  73
  74 static void orderQuads(skiatest::Reporter* reporter, const SkDQuad& quad, double radius,
  75         SkTArray<double, false>* tArray) {
  76     double r = radius;
  77     double s = r * SK_ScalarTanPIOver8;
  78     double m = r * SK_ScalarRoot2Over2;
  79     // construct circle from quads
  80     const SkDQuad circle[8] = {{{{ r,  0}, { r, -s}, { m, -m}}},
  81                                 {{{ m, -m}, { s, -r}, { 0, -r}}},
  82                                 {{{ 0, -r}, {-s, -r}, {-m, -m}}},
  83                                 {{{-m, -m}, {-r, -s}, {-r,  0}}},
  84                                 {{{-r,  0}, {-r,  s}, {-m,  m}}},
  85                                 {{{-m,  m}, {-s,  r}, { 0,  r}}},
  86                                 {{{ 0,  r}, { s,  r}, { m,  m}}},
  87                                 {{{ m,  m}, { r,  s}, { r,  0}}}};
  88     for (int octant = 0; octant < 8; ++octant) {
  89         double t = testArc(reporter, quad, circle[octant], octant);
  90         if (t < 0) {
  91             continue;
  92         }
  93         for (int index = 0; index < tArray->count(); ++index) {
  94             double matchT = (*tArray)[index];
  95             if (approximately_equal(t, matchT)) {
  96                 goto next;
  97             }
  98         }
  99         tArray->push_back(t);
 100 next:   ;
 101     }
 102 }
 103
 104 static double quadAngle(skiatest::Reporter* reporter, const SkDQuad& quad, double t) {
 105     const SkDVector& pt = quad.ptAtT(t) - quad[0];
 106     double angle = (atan2(pt.fY, pt.fX) + SK_ScalarPI) * 8 / (SK_ScalarPI * 2);
 107     REPORTER_ASSERT(reporter, angle >= 0 && angle <= 8);
 108     return angle;
 109 }
 110
 111 static bool angleDirection(double a1, double a2) {
 112     double delta = a1 - a2;
 113     return (delta < 4 && delta > 0) || delta < -4;
 114 }
 115
 116 static void setQuadHullSweep(const SkDQuad& quad, SkDVector sweep[2]) {
 117     sweep[0] = quad[1] - quad[0];
 118     sweep[1] = quad[2] - quad[0];
 119 }
 120
 121 static double distEndRatio(double dist, const SkDQuad& quad) {
 122     SkDVector v[] = {quad[2] - quad[0], quad[1] - quad[0], quad[2] - quad[1]};
 123     double longest = SkTMax(v[0].length(), SkTMax(v[1].length(), v[2].length()));
 124     return longest / dist;
 125 }
 126
 127 static bool checkParallel(skiatest::Reporter* reporter, const SkDQuad& quad1, const SkDQuad& quad2) {
 128     SkDVector sweep[2], tweep[2];
 129     setQuadHullSweep(quad1, sweep);
 130     setQuadHullSweep(quad2, tweep);
 131     // if the ctrl tangents are not nearly parallel, use them
 132     // solve for opposite direction displacement scale factor == m
 133     // initial dir = v1.cross(v2) == v2.x * v1.y - v2.y * v1.x
 134     // displacement of q1[1] : dq1 = { -m * v1.y, m * v1.x } + q1[1]
 135     // straight angle when : v2.x * (dq1.y - q1[0].y) == v2.y * (dq1.x - q1[0].x)
 136     //                       v2.x * (m * v1.x + v1.y) == v2.y * (-m * v1.y + v1.x)
 137     // - m * (v2.x * v1.x + v2.y * v1.y) == v2.x * v1.y - v2.y * v1.x
 138     // m = (v2.y * v1.x - v2.x * v1.y) / (v2.x * v1.x + v2.y * v1.y)
 139     // m = v1.cross(v2) / v1.dot(v2)
 140     double s0dt0 = sweep[0].dot(tweep[0]);
 141     REPORTER_ASSERT(reporter, s0dt0 != 0);
 142     double s0xt0 = sweep[0].crossCheck(tweep[0]);
 143     double m = s0xt0 / s0dt0;
 144     double sDist = sweep[0].length() * m;
 145     double tDist = tweep[0].length() * m;
 146     bool useS = fabs(sDist) < fabs(tDist);
 147     double mFactor = fabs(useS ? distEndRatio(sDist, quad1) : distEndRatio(tDist, quad2));
 148     if (mFactor < 5000) {  // empirically found limit
 149         return s0xt0 < 0;
 150     }
 151     SkDVector m0 = quad1.ptAtT(0.5) - quad1[0];
 152     SkDVector m1 = quad2.ptAtT(0.5) - quad2[0];
 153     return m0.crossCheck(m1) < 0;
 154 }
 155
 156 /* returns
 157    -1 if overlaps
 158     0 if no overlap cw
 159     1 if no overlap ccw
 160 */
 161 static int quadHullsOverlap(skiatest::Reporter* reporter, const SkDQuad& quad1,
 162         const SkDQuad& quad2) {
 163     SkDVector sweep[2], tweep[2];
 164     setQuadHullSweep(quad1, sweep);
 165     setQuadHullSweep(quad2, tweep);
 166     double s0xs1 = sweep[0].crossCheck(sweep[1]);
 167     double s0xt0 = sweep[0].crossCheck(tweep[0]);
 168     double s1xt0 = sweep[1].crossCheck(tweep[0]);
 169     bool tBetweenS = s0xs1 > 0 ? s0xt0 > 0 && s1xt0 < 0 : s0xt0 < 0 && s1xt0 > 0;
 170     double s0xt1 = sweep[0].crossCheck(tweep[1]);
 171     double s1xt1 = sweep[1].crossCheck(tweep[1]);
 172     tBetweenS |= s0xs1 > 0 ? s0xt1 > 0 && s1xt1 < 0 : s0xt1 < 0 && s1xt1 > 0;
 173     double t0xt1 = tweep[0].crossCheck(tweep[1]);
 174     if (tBetweenS) {
 175         return -1;
 176     }
 177     if ((s0xt0 == 0 && s1xt1 == 0) || (s1xt0 == 0 && s0xt1 == 0)) {  // s0 to s1 equals t0 to t1
 178         return -1;
 179     }
 180     bool sBetweenT = t0xt1 > 0 ? s0xt0 < 0 && s0xt1 > 0 : s0xt0 > 0 && s0xt1 < 0;
 181     sBetweenT |= t0xt1 > 0 ? s1xt0 < 0 && s1xt1 > 0 : s1xt0 > 0 && s1xt1 < 0;
 182     if (sBetweenT) {
 183         return -1;
 184     }
 185     // if all of the sweeps are in the same half plane, then the order of any pair is enough
 186     if (s0xt0 >= 0 && s0xt1 >= 0 && s1xt0 >= 0 && s1xt1 >= 0) {
 187         return 0;
 188     }
 189     if (s0xt0 <= 0 && s0xt1 <= 0 && s1xt0 <= 0 && s1xt1 <= 0) {
 190         return 1;
 191     }
 192     // if the outside sweeps are greater than 180 degress:
 193         // first assume the inital tangents are the ordering
 194         // if the midpoint direction matches the inital order, that is enough
 195     SkDVector m0 = quad1.ptAtT(0.5) - quad1[0];
 196     SkDVector m1 = quad2.ptAtT(0.5) - quad2[0];
 197     double m0xm1 = m0.crossCheck(m1);
 198     if (s0xt0 > 0 && m0xm1 > 0) {
 199         return 0;
 200     }
 201     if (s0xt0 < 0 && m0xm1 < 0) {
 202         return 1;
 203     }
 204     REPORTER_ASSERT(reporter, s0xt0 != 0);
 205     return checkParallel(reporter, quad1, quad2);
 206 }
 207
 208 static double radianSweep(double start, double end) {
 209     double sweep = end - start;
 210     if (sweep > SK_ScalarPI) {
 211         sweep -= 2 * SK_ScalarPI;
 212     } else if (sweep < -SK_ScalarPI) {
 213         sweep += 2 * SK_ScalarPI;
 214     }
 215     return sweep;
 216 }
 217
 218 static bool radianBetween(double start, double test, double end) {
 219     double startToEnd = radianSweep(start, end);
 220     double startToTest = radianSweep(start, test);
 221     double testToEnd = radianSweep(test, end);
 222     return (startToTest <= 0 && testToEnd <= 0 && startToTest >= startToEnd) ||
 223         (startToTest >= 0 && testToEnd >= 0 && startToTest <= startToEnd);
 224 }
 225
 226 static bool orderTRange(skiatest::Reporter* reporter, const SkDQuad& quad1, const SkDQuad& quad2,
 227         double r, TRange* result) {
 228     SkTArray<double, false> t1Array, t2Array;
 229     orderQuads(reporter, quad1, r, &t1Array);
 230     orderQuads(reporter,quad2, r, &t2Array);
 231     if (!t1Array.count() || !t2Array.count()) {
 232         return false;
 233     }
 234     SkTQSort<double>(t1Array.begin(), t1Array.end() - 1);
 235     SkTQSort<double>(t2Array.begin(), t2Array.end() - 1);
 236     double t1 = result->tMin1 = t1Array[0];
 237     double t2 = result->tMin2 = t2Array[0];
 238     double a1 = quadAngle(reporter,quad1, t1);
 239     double a2 = quadAngle(reporter,quad2, t2);
 240     if (approximately_equal(a1, a2)) {
 241         return false;
 242     }
 243     bool refCCW = angleDirection(a1, a2);
 244     result->t1 = t1;
 245     result->t2 = t2;
 246     result->tMin = SkTMin(t1, t2);
 247     result->a1 = a1;
 248     result->a2 = a2;
 249     result->ccw = refCCW;
 250     return true;
 251 }
 252
 253 static bool equalPoints(const SkDPoint& pt1, const SkDPoint& pt2, double max) {
 254     return approximately_zero_when_compared_to(pt1.fX - pt2.fX, max)
 255             && approximately_zero_when_compared_to(pt1.fY - pt2.fY, max);
 256 }
 257
 258 static double maxDist(const SkDQuad& quad) {
 259     SkDRect bounds;
 260     bounds.setBounds(quad);
 261     SkDVector corner[4] = {
 262         { bounds.fLeft - quad[0].fX, bounds.fTop - quad[0].fY },
 263         { bounds.fRight - quad[0].fX, bounds.fTop - quad[0].fY },
 264         { bounds.fLeft - quad[0].fX, bounds.fBottom - quad[0].fY },
 265         { bounds.fRight - quad[0].fX, bounds.fBottom - quad[0].fY }
 266     };
 267     double max = 0;
 268     for (unsigned index = 0; index < SK_ARRAY_COUNT(corner); ++index) {
 269         max = SkTMax(max, corner[index].length());
 270     }
 271     return max;
 272 }
 273
 274 static double maxQuad(const SkDQuad& quad) {
 275     double max = 0;
 276     for (int index = 0; index < 2; ++index) {
 277         max = SkTMax(max, fabs(quad[index].fX));
 278         max = SkTMax(max, fabs(quad[index].fY));
 279     }
 280     return max;
 281 }
 282
 283 static bool bruteMinT(skiatest::Reporter* reporter, const SkDQuad& quad1, const SkDQuad& quad2,
 284         TRange* lowerRange, TRange* upperRange) {
 285     double maxRadius = SkTMin(maxDist(quad1), maxDist(quad2));
 286     double maxQuads = SkTMax(maxQuad(quad1), maxQuad(quad2));
 287     double r = maxRadius / 2;
 288     double rStep = r / 2;
 289     SkDPoint best1 = {SK_ScalarInfinity, SK_ScalarInfinity};
 290     SkDPoint best2 = {SK_ScalarInfinity, SK_ScalarInfinity};
 291     int bestCCW = -1;
 292     double bestR = maxRadius;
 293     upperRange->tMin = 0;
 294     lowerRange->tMin = 1;
 295     do {
 296         do {  // find upper bounds of single result
 297             TRange tRange;
 298             bool stepUp = orderTRange(reporter, quad1, quad2, r, &tRange);
 299             if (stepUp) {
 300                 SkDPoint pt1 = quad1.ptAtT(tRange.t1);
 301                 if (equalPoints(pt1, best1, maxQuads)) {
 302                     break;
 303                 }
 304                 best1 = pt1;
 305                 SkDPoint pt2 = quad2.ptAtT(tRange.t2);
 306                 if (equalPoints(pt2, best2, maxQuads)) {
 307                     break;
 308                 }
 309                 best2 = pt2;
 310                 if (gPathOpsAngleIdeasVerbose) {
 311                     SkDebugf("u bestCCW=%d ccw=%d bestMin=%1.9g:%1.9g r=%1.9g tMin=%1.9g\n",
 312                             bestCCW, tRange.ccw, lowerRange->tMin, upperRange->tMin, r,
 313                             tRange.tMin);
 314                 }
 315                 if (bestCCW >= 0 && bestCCW != (int) tRange.ccw) {
 316                     if (tRange.tMin < upperRange->tMin) {
 317                         upperRange->tMin = 0;
 318                     } else {
 319                         stepUp = false;
 320                     }
 321                 }
 322                 if (upperRange->tMin < tRange.tMin) {
 323                     bestCCW = tRange.ccw;
 324                     bestR = r;
 325                     *upperRange = tRange;
 326                 }
 327                 if (lowerRange->tMin > tRange.tMin) {
 328                     *lowerRange = tRange;
 329                 }
 330             }
 331             r += stepUp ? rStep : -rStep;
 332             rStep /= 2;
 333         } while (rStep > FLT_EPSILON);
 334         if (bestCCW < 0) {
 335             REPORTER_ASSERT(reporter, bestR < maxRadius);
 336             return false;
 337         }
 338         double lastHighR = bestR;
 339         r = bestR / 2;
 340         rStep = r / 2;
 341         do {  // find lower bounds of single result
 342             TRange tRange;
 343             bool success = orderTRange(reporter, quad1, quad2, r, &tRange);
 344             if (success) {
 345                 if (gPathOpsAngleIdeasVerbose) {
 346                     SkDebugf("l bestCCW=%d ccw=%d bestMin=%1.9g:%1.9g r=%1.9g tMin=%1.9g\n",
 347                             bestCCW, tRange.ccw, lowerRange->tMin, upperRange->tMin, r,
 348                             tRange.tMin);
 349                 }
 350                 if (bestCCW != (int) tRange.ccw || upperRange->tMin < tRange.tMin) {
 351                     bestCCW = tRange.ccw;
 352                     *upperRange = tRange;
 353                     bestR = lastHighR;
 354                     break;  // need to establish a new upper bounds
 355                 }
 356                 SkDPoint pt1 = quad1.ptAtT(tRange.t1);
 357                 SkDPoint pt2 = quad2.ptAtT(tRange.t2);
 358                 if (equalPoints(pt1, best1, maxQuads)) {
 359                     goto breakOut;
 360                 }
 361                 best1 = pt1;
 362                 if (equalPoints(pt2, best2, maxQuads)) {
 363                     goto breakOut;
 364                 }
 365                 best2 = pt2;
 366                 if (equalPoints(pt1, pt2, maxQuads)) {
 367                     success = false;
 368                 } else {
 369                     if (upperRange->tMin < tRange.tMin) {
 370                         *upperRange = tRange;
 371                     }
 372                     if (lowerRange->tMin > tRange.tMin) {
 373                         *lowerRange = tRange;
 374                     }
 375                 }
 376                 lastHighR = SkTMin(r, lastHighR);
 377             }
 378             r += success ? -rStep : rStep;
 379             rStep /= 2;
 380         } while (rStep > FLT_EPSILON);
 381     } while (rStep > FLT_EPSILON);
 382 breakOut:
 383     if (gPathOpsAngleIdeasVerbose) {
 384         SkDebugf("l a2-a1==%1.9g\n", lowerRange->a2 - lowerRange->a1);
 385     }
 386     return true;
 387 }
 388
 389 static void bruteForce(skiatest::Reporter* reporter, const SkDQuad& quad1, const SkDQuad& quad2,
 390         bool ccw) {
 391     if (!gPathOpsAngleIdeasEnableBruteCheck) {
 392         return;
 393     }
 394     TRange lowerRange, upperRange;
 395     bool result = bruteMinT(reporter, quad1, quad2, &lowerRange, &upperRange);
 396     REPORTER_ASSERT(reporter, result);
 397     double angle = fabs(lowerRange.a2 - lowerRange.a1);
 398     REPORTER_ASSERT(reporter, angle > 3.998 || ccw == upperRange.ccw);
 399 }
 400
 401 static bool bruteForceCheck(skiatest::Reporter* reporter, const SkDQuad& quad1,
 402         const SkDQuad& quad2, bool ccw) {
 403     TRange lowerRange, upperRange;
 404     bool result = bruteMinT(reporter, quad1, quad2, &lowerRange, &upperRange);
 405     REPORTER_ASSERT(reporter, result);
 406     return ccw == upperRange.ccw;
 407 }
 408
 409 class PathOpsSegmentTester {
 410 public:
 411     static void ConstructQuad(SkOpSegment* segment, SkPoint shortQuad[3]) {
 412         segment->debugConstructQuad(shortQuad);
 413     }
 414 };
 415
 416 static void makeSegment(const SkDQuad& quad, SkPoint shortQuad[3], SkOpSegment* result) {
 417     shortQuad[0] = quad[0].asSkPoint();
 418     shortQuad[1] = quad[1].asSkPoint();
 419     shortQuad[2] = quad[2].asSkPoint();
 420     PathOpsSegmentTester::ConstructQuad(result, shortQuad);
 421 }
 422
 423 static void testQuadAngles(skiatest::Reporter* reporter, const SkDQuad& quad1, const SkDQuad& quad2,
 424         int testNo) {
 425     SkPoint shortQuads[2][3];
 426     SkOpSegment seg[2];
 427     makeSegment(quad1, shortQuads[0], &seg[0]);
 428     makeSegment(quad2, shortQuads[1], &seg[1]);
 429     int realOverlap = PathOpsAngleTester::ConvexHullOverlaps(*seg[0].debugLastAngle(),
 430             *seg[1].debugLastAngle());
 431     const SkDPoint& origin = quad1[0];
 432     REPORTER_ASSERT(reporter, origin == quad2[0]);
 433     double a1s = atan2(origin.fY - quad1[1].fY, quad1[1].fX - origin.fX);
 434     double a1e = atan2(origin.fY - quad1[2].fY, quad1[2].fX - origin.fX);
 435     double a2s = atan2(origin.fY - quad2[1].fY, quad2[1].fX - origin.fX);
 436     double a2e = atan2(origin.fY - quad2[2].fY, quad2[2].fX - origin.fX);
 437     bool oldSchoolOverlap = radianBetween(a1s, a2s, a1e)
 438         || radianBetween(a1s, a2e, a1e) || radianBetween(a2s, a1s, a2e)
 439         || radianBetween(a2s, a1e, a2e);
 440     int overlap = quadHullsOverlap(reporter, quad1, quad2);
 441     bool realMatchesOverlap = realOverlap == overlap || SK_ScalarPI - fabs(a2s - a1s) < 0.002;
 442     if (realOverlap != overlap) {
 443         SkDebugf("\nSK_ScalarPI - fabs(a2s - a1s) = %1.9g\n", SK_ScalarPI - fabs(a2s - a1s));
 444     }
 445     if (!realMatchesOverlap) {
 446         DumpQ(quad1, quad2, testNo);
 447     }
 448     REPORTER_ASSERT(reporter, realMatchesOverlap);
 449     if (oldSchoolOverlap != (overlap < 0)) {
 450         overlap = quadHullsOverlap(reporter, quad1, quad2);  // set a breakpoint and debug if assert fires
 451         REPORTER_ASSERT(reporter, oldSchoolOverlap == (overlap < 0));
 452     }
 453     SkDVector v1s = quad1[1] - quad1[0];
 454     SkDVector v1e = quad1[2] - quad1[0];
 455     SkDVector v2s = quad2[1] - quad2[0];
 456     SkDVector v2e = quad2[2] - quad2[0];
 457     double vDir[2] = { v1s.cross(v1e), v2s.cross(v2e) };
 458     bool ray1In2 = v1s.cross(v2s) * vDir[1] <= 0 && v1s.cross(v2e) * vDir[1] >= 0;
 459     bool ray2In1 = v2s.cross(v1s) * vDir[0] <= 0 && v2s.cross(v1e) * vDir[0] >= 0;
 460     if (overlap >= 0) {
 461         // verify that hulls really don't overlap
 462         REPORTER_ASSERT(reporter, !ray1In2);
 463         REPORTER_ASSERT(reporter, !ray2In1);
 464         bool ctrl1In2 = v1e.cross(v2s) * vDir[1] <= 0 && v1e.cross(v2e) * vDir[1] >= 0;
 465         REPORTER_ASSERT(reporter, !ctrl1In2);
 466         bool ctrl2In1 = v2e.cross(v1s) * vDir[0] <= 0 && v2e.cross(v1e) * vDir[0] >= 0;
 467         REPORTER_ASSERT(reporter, !ctrl2In1);
 468         // check answer against reference
 469         bruteForce(reporter, quad1, quad2, overlap > 0);
 470     }
 471     // continue end point rays and see if they intersect the opposite curve
 472     SkDLine rays[] = {{{origin, quad2[2]}}, {{origin, quad1[2]}}};
 473     const SkDQuad* quads[] = {&quad1, &quad2};
 474     SkDVector midSpokes[2];
 475     SkIntersections intersect[2];
 476     double minX, minY, maxX, maxY;
 477     minX = minY = SK_ScalarInfinity;
 478     maxX = maxY = -SK_ScalarInfinity;
 479     double maxWidth = 0;
 480     bool useIntersect = false;
 481     double smallestTs[] = {1, 1};
 482     for (unsigned index = 0; index < SK_ARRAY_COUNT(quads); ++index) {
 483         const SkDQuad& q = *quads[index];
 484         midSpokes[index] = q.ptAtT(0.5) - origin;
 485         minX = SkTMin(SkTMin(SkTMin(minX, origin.fX), q[1].fX), q[2].fX);
 486         minY = SkTMin(SkTMin(SkTMin(minY, origin.fY), q[1].fY), q[2].fY);
 487         maxX = SkTMax(SkTMax(SkTMax(maxX, origin.fX), q[1].fX), q[2].fX);
 488         maxY = SkTMax(SkTMax(SkTMax(maxY, origin.fY), q[1].fY), q[2].fY);
 489         maxWidth = SkTMax(maxWidth, SkTMax(maxX - minX, maxY - minY));
 490         intersect[index].intersectRay(q, rays[index]);
 491         const SkIntersections& i = intersect[index];
 492         REPORTER_ASSERT(reporter, i.used() >= 1);
 493         bool foundZero = false;
 494         double smallT = 1;
 495         for (int idx2 = 0; idx2 < i.used(); ++idx2) {
 496             double t = i[0][idx2];
 497             if (t == 0) {
 498                 foundZero = true;
 499                 continue;
 500             }
 501             if (smallT > t) {
 502                 smallT = t;
 503             }
 504         }
 505         REPORTER_ASSERT(reporter, foundZero == true);
 506         if (smallT == 1) {
 507             continue;
 508         }
 509         SkDVector ray = q.ptAtT(smallT) - origin;
 510         SkDVector end = rays[index][1] - origin;
 511         if (ray.fX * end.fX < 0 || ray.fY * end.fY < 0) {
 512             continue;
 513         }
 514         double rayDist = ray.length();
 515         double endDist = end.length();
 516         double delta = fabs(rayDist - endDist) / maxWidth;
 517         if (delta > 1e-4) {
 518             useIntersect ^= true;
 519         }
 520         smallestTs[index] = smallT;
 521     }
 522     bool firstInside;
 523     if (useIntersect) {
 524         int sIndex = (int) (smallestTs[1] < 1);
 525         REPORTER_ASSERT(reporter, smallestTs[sIndex ^ 1] == 1);
 526         double t = smallestTs[sIndex];
 527         const SkDQuad& q = *quads[sIndex];
 528         SkDVector ray = q.ptAtT(t) - origin;
 529         SkDVector end = rays[sIndex][1] - origin;
 530         double rayDist = ray.length();
 531         double endDist = end.length();
 532         SkDVector mid = q.ptAtT(t / 2) - origin;
 533         double midXray = mid.crossCheck(ray);
 534         if (gPathOpsAngleIdeasVerbose) {
 535             SkDebugf("rayDist>endDist:%d sIndex==0:%d vDir[sIndex]<0:%d midXray<0:%d\n",
 536                     rayDist > endDist, sIndex == 0, vDir[sIndex] < 0, midXray < 0);
 537         }
 538         SkASSERT(SkScalarSignAsInt(SkDoubleToScalar(midXray))
 539             == SkScalarSignAsInt(SkDoubleToScalar(vDir[sIndex])));
 540         firstInside = (rayDist > endDist) ^ (sIndex == 0) ^ (vDir[sIndex] < 0);
 541     } else if (overlap >= 0) {
 542         return;  // answer has already been determined
 543     } else {
 544         firstInside = checkParallel(reporter, quad1, quad2);
 545     }
 546     if (overlap < 0) {
 547         SkDEBUGCODE(int realEnds =)
 548                 PathOpsAngleTester::EndsIntersect(*seg[0].debugLastAngle(),
 549                 *seg[1].debugLastAngle());
 550         SkASSERT(realEnds == (firstInside ? 1 : 0));
 551     }
 552     bruteForce(reporter, quad1, quad2, firstInside);
 553 }
 554
 555 DEF_TEST(PathOpsAngleOverlapHullsOne, reporter) {
 556 //    gPathOpsAngleIdeasVerbose = true;
 557     const SkDQuad quads[] = {
 558 {{{939.4808349609375, 914.355224609375}, {-357.7921142578125, 590.842529296875}, {736.8936767578125, -350.717529296875}}},
 559 {{{939.4808349609375, 914.355224609375}, {-182.85418701171875, 634.4552001953125}, {-509.62615966796875, 576.1182861328125}}}
 560     };
 561     for (int index = 0; index < (int) SK_ARRAY_COUNT(quads); index += 2) {
 562         testQuadAngles(reporter, quads[index], quads[index + 1], 0);
 563     }
 564 }
 565
 566 DEF_TEST(PathOpsAngleOverlapHulls, reporter) {
 567     if (!gPathOpsAngleIdeasVerbose) {  // takes a while to run -- so exclude it by default
 568         return;
 569     }
 570     SkRandom ran;
 571     for (int index = 0; index < 100000; ++index) {
 572         if (index % 1000 == 999) SkDebugf(".");
 573         SkDPoint origin = {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)};
 574         SkDQuad quad1 = {{origin, {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)},
 575             {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)}}};
 576         if (quad1[0] == quad1[2]) {
 577             continue;
 578         }
 579         SkDQuad quad2 = {{origin, {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)},
 580             {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)}}};
 581         if (quad2[0] == quad2[2]) {
 582             continue;
 583         }
 584         SkIntersections i;
 585         i.intersect(quad1, quad2);
 586         REPORTER_ASSERT(reporter, i.used() >= 1);
 587         if (i.used() > 1) {
 588             continue;
 589         }
 590         testQuadAngles(reporter, quad1, quad2, index);
 591     }
 592 }
 593
 594 DEF_TEST(PathOpsAngleBruteT, reporter) {
 595     if (!gPathOpsAngleIdeasVerbose) {  // takes a while to run -- so exclude it by default
 596         return;
 597     }
 598     SkRandom ran;
 599     double smaller = SK_Scalar1;
 600     SkDQuad small[2];
 601     SkDEBUGCODE(int smallIndex);
 602     for (int index = 0; index < 100000; ++index) {
 603         SkDPoint origin = {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)};
 604         SkDQuad quad1 = {{origin, {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)},
 605             {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)}}};
 606         if (quad1[0] == quad1[2]) {
 607             continue;
 608         }
 609         SkDQuad quad2 = {{origin, {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)},
 610             {ran.nextRangeF(-1000, 1000), ran.nextRangeF(-1000, 1000)}}};
 611         if (quad2[0] == quad2[2]) {
 612             continue;
 613         }
 614         SkIntersections i;
 615         i.intersect(quad1, quad2);
 616         REPORTER_ASSERT(reporter, i.used() >= 1);
 617         if (i.used() > 1) {
 618             continue;
 619         }
 620         TRange lowerRange, upperRange;
 621         bool result = bruteMinT(reporter, quad1, quad2, &lowerRange, &upperRange);
 622         REPORTER_ASSERT(reporter, result);
 623         double min = SkTMin(upperRange.t1, upperRange.t2);
 624         if (smaller > min) {
 625             small[0] = quad1;
 626             small[1] = quad2;
 627             SkDEBUGCODE(smallIndex = index);
 628             smaller = min;
 629         }
 630     }
 631 #ifdef SK_DEBUG
 632     DumpQ(small[0], small[1], smallIndex);
 633 #endif
 634 }
 635
 636 DEF_TEST(PathOpsAngleBruteTOne, reporter) {
 637 //    gPathOpsAngleIdeasVerbose = true;
 638     const SkDQuad quads[] = {
 639 {{{-770.8492431640625, 948.2369384765625}, {-853.37066650390625, 972.0301513671875}, {-200.62042236328125, -26.7174072265625}}},
 640 {{{-770.8492431640625, 948.2369384765625}, {513.602783203125, 578.8681640625}, {960.641357421875, -813.69757080078125}}},
 641 {{{563.8267822265625, -107.4566650390625}, {-44.67724609375, -136.57452392578125}, {492.3856201171875, -268.79644775390625}}},
 642 {{{563.8267822265625, -107.4566650390625}, {708.049072265625, -100.77789306640625}, {-48.88226318359375, 967.9022216796875}}},
 643 {{{598.857421875, 846.345458984375}, {-644.095703125, -316.12921142578125}, {-97.64599609375, 20.6158447265625}}},
 644 {{{598.857421875, 846.345458984375}, {715.7142333984375, 955.3599853515625}, {-919.9478759765625, 691.611328125}}},
 645     };
 646     TRange lowerRange, upperRange;
 647     bruteMinT(reporter, quads[0], quads[1], &lowerRange, &upperRange);
 648     bruteMinT(reporter, quads[2], quads[3], &lowerRange, &upperRange);
 649     bruteMinT(reporter, quads[4], quads[5], &lowerRange, &upperRange);
 650 }
 651
 652 /*
 653 The sorting problem happens when the inital tangents are not a true indicator of the curve direction
 654 Nearly always, the initial tangents do give the right answer,
 655 so the trick is to figure out when the initial tangent cannot be trusted.
 656 If the convex hulls of both curves are in the same half plane, and not overlapping, sorting the
 657 hulls is enough.
 658 If the hulls overlap, and have the same general direction, then intersect the shorter end point ray
 659 with the opposing curve, and see on which side of the shorter curve the opposing intersection lies.
 660 Otherwise, if the control vector is extremely short, likely the point on curve must be computed
 661 If moving the control point slightly can change the sign of the cross product, either answer could
 662 be "right".
 663 We need to determine how short is extremely short. Move the control point a set percentage of
 664 the largest length to determine how stable the curve is vis-a-vis the initial tangent.
 665 */
 666
 667 static const SkDQuad extremeTests[][2] = {
 668     {
 669         {{{-708.0077926931004,-154.61669472244046},
 670             {-707.9234268635319,-154.30459999551294},
 671             {505.58447265625,-504.9130859375}}},
 672         {{{-708.0077926931004,-154.61669472244046},
 673             {-711.127526325141,-163.9446090624656},
 674             {-32.39227294921875,-906.3277587890625}}},
 675     }, {
 676         {{{-708.0077926931004,-154.61669472244046},
 677             {-708.2875337527566,-154.36676458635623},
 678             {505.58447265625,-504.9130859375}}},
 679         {{{-708.0077926931004,-154.61669472244046},
 680             {-708.4111557216864,-154.5366642875255},
 681             {-32.39227294921875,-906.3277587890625}}},
 682     }, {
 683         {{{-609.0230951752058,-267.5435593490574},
 684             {-594.1120809906336,-136.08492475411555},
 685             {505.58447265625,-504.9130859375}}},
 686         {{{-609.0230951752058,-267.5435593490574},
 687             {-693.7467719138988,-341.3259237831895},
 688             {-32.39227294921875,-906.3277587890625}}}
 689     }, {
 690         {{{-708.0077926931004,-154.61669472244046},
 691             {-707.9994640658723,-154.58588461064852},
 692             {505.58447265625,-504.9130859375}}},
 693         {{{-708.0077926931004,-154.61669472244046},
 694             {-708.0239418990758,-154.6403553507124},
 695             {-32.39227294921875,-906.3277587890625}}}
 696     }, {
 697         {{{-708.0077926931004,-154.61669472244046},
 698             {-707.9993222215099,-154.55999389855003},
 699             {68.88981098017803,296.9273945411635}}},
 700         {{{-708.0077926931004,-154.61669472244046},
 701             {-708.0509091919608,-154.64675214697067},
 702             {-777.4871194247767,-995.1470120113145}}}
 703     }, {
 704         {{{-708.0077926931004,-154.61669472244046},
 705             {-708.0060491116379,-154.60889321524968},
 706             {229.97088707895057,-430.0569357467175}}},
 707         {{{-708.0077926931004,-154.61669472244046},
 708             {-708.013911296257,-154.6219143988058},
 709             {138.13162892614037,-573.3689311737394}}}
 710     }, {
 711         {{{-543.2570545751013,-237.29243831075053},
 712             {-452.4119186056987,-143.47223056267802},
 713             {229.97088707895057,-430.0569357467175}}},
 714         {{{-543.2570545751013,-237.29243831075053},
 715             {-660.5330371214436,-362.0016148388},
 716             {138.13162892614037,-573.3689311737394}}},
 717     },
 718 };
 719
 720 static double endCtrlRatio(const SkDQuad quad) {
 721     SkDVector longEdge = quad[2] - quad[0];
 722     double longLen = longEdge.length();
 723     SkDVector shortEdge = quad[1] - quad[0];
 724     double shortLen = shortEdge.length();
 725     return longLen / shortLen;
 726 }
 727
 728 static void computeMV(const SkDQuad& quad, const SkDVector& v, double m, SkDVector mV[2]) {
 729         SkDPoint mPta = {quad[1].fX - m * v.fY, quad[1].fY + m * v.fX};
 730         SkDPoint mPtb = {quad[1].fX + m * v.fY, quad[1].fY - m * v.fX};
 731         mV[0] = mPta - quad[0];
 732         mV[1] = mPtb - quad[0];
 733 }
 734
 735 static double mDistance(skiatest::Reporter* reporter, bool agrees, const SkDQuad& q1,
 736         const SkDQuad& q2) {
 737     if (1 && agrees) {
 738         return SK_ScalarMax;
 739     }
 740     // how close is the angle from inflecting in the opposite direction?
 741     SkDVector v1 = q1[1] - q1[0];
 742     SkDVector v2 = q2[1] - q2[0];
 743     double dir = v1.crossCheck(v2);
 744     REPORTER_ASSERT(reporter, dir != 0);
 745     // solve for opposite direction displacement scale factor == m
 746     // initial dir = v1.cross(v2) == v2.x * v1.y - v2.y * v1.x
 747     // displacement of q1[1] : dq1 = { -m * v1.y, m * v1.x } + q1[1]
 748     // straight angle when : v2.x * (dq1.y - q1[0].y) == v2.y * (dq1.x - q1[0].x)
 749     //                       v2.x * (m * v1.x + v1.y) == v2.y * (-m * v1.y + v1.x)
 750     // - m * (v2.x * v1.x + v2.y * v1.y) == v2.x * v1.y - v2.y * v1.x
 751     // m = (v2.y * v1.x - v2.x * v1.y) / (v2.x * v1.x + v2.y * v1.y)
 752     // m = v1.cross(v2) / v1.dot(v2)
 753     double div = v1.dot(v2);
 754     REPORTER_ASSERT(reporter, div != 0);
 755     double m = dir / div;
 756     SkDVector mV1[2], mV2[2];
 757     computeMV(q1, v1, m, mV1);
 758     computeMV(q2, v2, m, mV2);
 759     double dist1 = v1.length() * m;
 760     double dist2 = v2.length() * m;
 761     if (gPathOpsAngleIdeasVerbose) {
 762         SkDebugf("%c r1=%1.9g r2=%1.9g m=%1.9g dist1=%1.9g dist2=%1.9g "
 763                 " dir%c 1a=%1.9g 1b=%1.9g 2a=%1.9g 2b=%1.9g\n", agrees ? 'T' : 'F',
 764                 endCtrlRatio(q1), endCtrlRatio(q2), m, dist1, dist2, dir > 0 ? '+' : '-',
 765                 mV1[0].crossCheck(v2), mV1[1].crossCheck(v2),
 766                 mV2[0].crossCheck(v1), mV2[1].crossCheck(v1));
 767     }
 768     if (1) {
 769         bool use1 = fabs(dist1) < fabs(dist2);
 770         if (gPathOpsAngleIdeasVerbose) {
 771             SkDebugf("%c dist=%1.9g r=%1.9g\n", agrees ? 'T' : 'F', use1 ? dist1 : dist2,
 772                 use1 ? distEndRatio(dist1, q1) : distEndRatio(dist2, q2));
 773         }
 774         return fabs(use1 ? distEndRatio(dist1, q1) : distEndRatio(dist2, q2));
 775     }
 776     return SK_ScalarMax;
 777 }
 778
 779 static void midPointAgrees(skiatest::Reporter* reporter, const SkDQuad& q1, const SkDQuad& q2,
 780         bool ccw) {
 781     SkDPoint mid1 = q1.ptAtT(0.5);
 782     SkDVector m1 = mid1 - q1[0];
 783     SkDPoint mid2 = q2.ptAtT(0.5);
 784     SkDVector m2 = mid2 - q2[0];
 785     REPORTER_ASSERT(reporter, ccw ? m1.crossCheck(m2) < 0 : m1.crossCheck(m2) > 0);
 786 }
 787
 788 DEF_TEST(PathOpsAngleExtreme, reporter) {
 789     if (!gPathOpsAngleIdeasVerbose) {  // takes a while to run -- so exclude it by default
 790         return;
 791     }
 792     double maxR = SK_ScalarMax;
 793     for (int index = 0; index < (int) SK_ARRAY_COUNT(extremeTests); ++index) {
 794         const SkDQuad& quad1 = extremeTests[index][0];
 795         const SkDQuad& quad2 = extremeTests[index][1];
 796         if (gPathOpsAngleIdeasVerbose) {
 797             SkDebugf("%s %d\n", __FUNCTION__, index);
 798         }
 799         REPORTER_ASSERT(reporter, quad1[0] == quad2[0]);
 800         SkIntersections i;
 801         i.intersect(quad1, quad2);
 802         REPORTER_ASSERT(reporter, i.used() == 1);
 803         REPORTER_ASSERT(reporter, i.pt(0) == quad1[0]);
 804         int overlap = quadHullsOverlap(reporter, quad1, quad2);
 805         REPORTER_ASSERT(reporter, overlap >= 0);
 806         SkDVector sweep[2], tweep[2];
 807         setQuadHullSweep(quad1, sweep);
 808         setQuadHullSweep(quad2, tweep);
 809         double s0xt0 = sweep[0].crossCheck(tweep[0]);
 810         REPORTER_ASSERT(reporter, s0xt0 != 0);
 811         bool ccw = s0xt0 < 0;
 812         bool agrees = bruteForceCheck(reporter, quad1, quad2, ccw);
 813         maxR = SkTMin(maxR, mDistance(reporter, agrees, quad1, quad2));
 814         if (agrees) {
 815             continue;
 816         }
 817         midPointAgrees(reporter, quad1, quad2, !ccw);
 818         SkDQuad q1 = quad1;
 819         SkDQuad q2 = quad2;
 820         double loFail = 1;
 821         double hiPass = 2;
 822         // double vectors until t passes
 823         do {
 824             q1[1].fX = quad1[0].fX * (1 - hiPass) + quad1[1].fX * hiPass;
 825             q1[1].fY = quad1[0].fY * (1 - hiPass) + quad1[1].fY * hiPass;
 826             q2[1].fX = quad2[0].fX * (1 - hiPass) + quad2[1].fX * hiPass;
 827             q2[1].fY = quad2[0].fY * (1 - hiPass) + quad2[1].fY * hiPass;
 828             agrees = bruteForceCheck(reporter, q1, q2, ccw);
 829             maxR = SkTMin(maxR, mDistance(reporter, agrees, q1, q2));
 830             if (agrees) {
 831                 break;
 832             }
 833             midPointAgrees(reporter, quad1, quad2, !ccw);
 834             loFail = hiPass;
 835             hiPass *= 2;
 836         } while (true);
 837         // binary search to find minimum pass
 838         double midTest = (loFail + hiPass) / 2;
 839         double step = (hiPass - loFail) / 4;
 840         while (step > FLT_EPSILON) {
 841             q1[1].fX = quad1[0].fX * (1 - midTest) + quad1[1].fX * midTest;
 842             q1[1].fY = quad1[0].fY * (1 - midTest) + quad1[1].fY * midTest;
 843             q2[1].fX = quad2[0].fX * (1 - midTest) + quad2[1].fX * midTest;
 844             q2[1].fY = quad2[0].fY * (1 - midTest) + quad2[1].fY * midTest;
 845             agrees = bruteForceCheck(reporter, q1, q2, ccw);
 846             maxR = SkTMin(maxR, mDistance(reporter, agrees, q1, q2));
 847             if (!agrees) {
 848                 midPointAgrees(reporter, quad1, quad2, !ccw);
 849             }
 850             midTest += agrees ? -step : step;
 851             step /= 2;
 852         }
 853 #ifdef SK_DEBUG
 854 //        DumpQ(q1, q2, 999);
 855 #endif
 856     }
 857     if (gPathOpsAngleIdeasVerbose) {
 858         SkDebugf("maxR=%1.9g\n", maxR);
 859     }
 860 }