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