CubicIntersection_Test.cpp revision 1304bb25aa3b0baa61fc2e2900fabcef88801b59
1/* 2 * Copyright 2012 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 "CurveIntersection.h" 8#include "CurveUtilities.h" 9#include "CubicIntersection_TestData.h" 10#include "Intersection_Tests.h" 11#include "Intersections.h" 12#include "TestUtilities.h" 13 14#define SHOW_ORIGINAL 1 15 16const int firstCubicIntersectionTest = 9; 17 18static void standardTestCases() { 19 for (size_t index = firstCubicIntersectionTest; index < tests_count; ++index) { 20 const Cubic& cubic1 = tests[index][0]; 21 const Cubic& cubic2 = tests[index][1]; 22 Cubic reduce1, reduce2; 23 int order1 = reduceOrder(cubic1, reduce1, kReduceOrder_NoQuadraticsAllowed, 24 kReduceOrder_TreatAsFill); 25 int order2 = reduceOrder(cubic2, reduce2, kReduceOrder_NoQuadraticsAllowed, 26 kReduceOrder_TreatAsFill); 27 if (order1 < 4) { 28 printf("%s [%d] cubic1 order=%d\n", __FUNCTION__, (int) index, order1); 29 continue; 30 } 31 if (order2 < 4) { 32 printf("%s [%d] cubic2 order=%d\n", __FUNCTION__, (int) index, order2); 33 continue; 34 } 35 if (implicit_matches(reduce1, reduce2)) { 36 printf("%s [%d] coincident\n", __FUNCTION__, (int) index); 37 continue; 38 } 39 Intersections tIntersections; 40 intersect(reduce1, reduce2, tIntersections); 41 if (!tIntersections.intersected()) { 42 printf("%s [%d] no intersection\n", __FUNCTION__, (int) index); 43 continue; 44 } 45 for (int pt = 0; pt < tIntersections.used(); ++pt) { 46 double tt1 = tIntersections.fT[0][pt]; 47 double tx1, ty1; 48 xy_at_t(cubic1, tt1, tx1, ty1); 49 double tt2 = tIntersections.fT[1][pt]; 50 double tx2, ty2; 51 xy_at_t(cubic2, tt2, tx2, ty2); 52 if (!AlmostEqualUlps(tx1, tx2)) { 53 printf("%s [%d,%d] x!= t1=%g (%g,%g) t2=%g (%g,%g)\n", 54 __FUNCTION__, (int)index, pt, tt1, tx1, ty1, tt2, tx2, ty2); 55 } 56 if (!AlmostEqualUlps(ty1, ty2)) { 57 printf("%s [%d,%d] y!= t1=%g (%g,%g) t2=%g (%g,%g)\n", 58 __FUNCTION__, (int)index, pt, tt1, tx1, ty1, tt2, tx2, ty2); 59 } 60 } 61 } 62} 63 64static const Cubic testSet[] = { 65{{0,1}, {4,5}, {1,0}, {5,3}}, 66{{0,1}, {3,5}, {1,0}, {5,4}}, 67 68{{0, 1}, {1, 6}, {1, 0}, {1, 0}}, 69{{0, 1}, {0, 1}, {1, 0}, {6, 1}}, 70 71{{0,1}, {3,4}, {1,0}, {5,1}}, 72{{0,1}, {1,5}, {1,0}, {4,3}}, 73 74{{0,1}, {1,2}, {1,0}, {6,1}}, 75{{0,1}, {1,6}, {1,0}, {2,1}}, 76 77{{0,1}, {0,5}, {1,0}, {4,0}}, 78{{0,1}, {0,4}, {1,0}, {5,0}}, 79 80{{0,1}, {3,4}, {1,0}, {3,0}}, 81{{0,1}, {0,3}, {1,0}, {4,3}}, 82 83{{0, 0}, {1, 2}, {3, 4}, {4, 4}}, 84{{0, 0}, {1, 2}, {3, 4}, {4, 4}}, 85{{4, 4}, {3, 4}, {1, 2}, {0, 0}}, 86 87{{0,1}, {2,3}, {1,0}, {1,0}}, 88{{0,1}, {0,1}, {1,0}, {3,2}}, 89 90{{0,2}, {0,1}, {1,0}, {1,0}}, 91{{0,1}, {0,1}, {2,0}, {1,0}}, 92 93{{0, 0}, {0, 1}, {1, 1}, {1, 0}}, 94{{1, 0}, {0, 0}, {0, 1}, {1, 1}}, 95 96{{0, 1}, {0, 2}, {1, 0}, {1, 0}}, 97{{0, 1}, {0, 1}, {1, 0}, {2, 0}}, 98 99{{0, 1}, {1, 6}, {1, 0}, {2, 0}}, 100{{0, 1}, {0, 2}, {1, 0}, {6, 1}}, 101 102{{0, 1}, {5, 6}, {1, 0}, {1, 0}}, 103{{0, 1}, {0, 1}, {1, 0}, {6, 5}}, 104 105{{95.837747722788592, 45.025976907939643}, {16.564570095652982, 0.72959763963222402}, {63.209855865319199, 68.047528419665767}, {57.640240647662544, 59.524565264361243}}, 106{{51.593891741518817, 38.53849970667553}, {62.34752929878772, 74.924924725166022}, {74.810149322641152, 34.17966562983564}, {29.368398119401373, 94.66719277886078}}, 107 108{{39.765160968417838, 33.060396198677083}, {5.1922921581157908, 66.854301452103215}, {31.619281802149157, 25.269248720849514}, {81.541621071073038, 70.025341524754353}}, 109{{46.078911165743556, 48.259962651999651}, {20.24450549867214, 49.403916182650214}, {0.26325131778756683, 24.46489805563581}, {15.915006546264051, 83.515023059917155}}, 110 111{{65.454505973241524, 93.881892270353575}, {45.867360264932437, 92.723972719499827}, {2.1464054482739447, 74.636369140183717}, {33.774068594804994, 40.770872887582925}}, 112{{72.963387832494163, 95.659300729473728}, {11.809496633619768, 82.209921247423594}, {13.456139067865974, 57.329313623406605}, {36.060621606214262, 70.867335643091849}}, 113 114{{32.484981432782945, 75.082940782924624}, {42.467313093350882, 48.131159948246157}, {3.5963115764764657, 43.208665839959245}, {79.442476890721579, 89.709102357602262}}, 115{{18.98573861410177, 93.308887208490106}, {40.405250173250792, 91.039661826118675}, {8.0467721950480584, 42.100282172719147}, {40.883324221187891, 26.030185504830527}}, 116 117{{7.5374809128872498, 82.441702896003477}, {22.444346930107265, 22.138854312775123}, {66.76091829629658, 50.753805856571446}, {78.193478508942519, 97.7932997968948}}, 118{{97.700573130371311, 53.53260215070685}, {87.72443481149358, 84.575876772671876}, {19.215031396232092, 47.032676472809484}, {11.989686410869325, 10.659507480757082}}, 119 120{{26.192053931854691, 9.8504326817814416}, {10.174241480498686, 98.476562741434464}, {21.177712558385782, 33.814968789841501}, {75.329030899018534, 55.02231980442177}}, 121{{56.222082700683771, 24.54395039218662}, {95.589995289030483, 81.050822735322086}, {28.180450866082897, 28.837706255185282}, {60.128952916771617, 87.311672180570511}}, 122 123{{42.449716172390481, 52.379709366885805}, {27.896043159019225, 48.797373636065686}, {92.770268299044233, 89.899302036454571}, {12.102066544863426, 99.43241951960718}}, 124{{45.77532924980639, 45.958701495993274}, {37.458701356062065, 68.393691335056758}, {37.569326692060258, 27.673713456687381}, {60.674866037757539, 62.47349659096146}}, 125 126{{67.426548091427676, 37.993772624988935}, {23.483695892376684, 90.476863174921306}, {35.597065061143162, 79.872482633158796}, {75.38634169631932, 18.244890038969412}}, 127{{61.336508189019057, 82.693132843213675}, {44.639380902349664, 54.074825790745592}, {16.815615499771951, 20.049704667203923}, {41.866884958868326, 56.735503699973002}}, 128 129{{67.4265481, 37.9937726}, {23.4836959, 90.4768632}, {35.5970651, 79.8724826}, {75.3863417, 18.24489}}, 130{{61.3365082, 82.6931328}, {44.6393809, 54.0748258}, {16.8156155, 20.0497047}, {41.866885, 56.7355037}}, 131 132{{18.1312339, 31.6473732}, {95.5711034, 63.5350219}, {92.3283165, 62.0158945}, {18.5656052, 32.1268808}}, 133{{97.402018, 35.7169972}, {33.1127443, 25.8935163}, {1.13970027, 54.9424981}, {56.4860195, 60.529264}}, 134}; 135 136const size_t testSetCount = sizeof(testSet) / sizeof(testSet[0]); 137 138static const Cubic newTestSet[] = { 139{{0,5}, {0,5}, {5,4}, {6,4}}, 140{{4,5}, {4,6}, {5,0}, {5,0}}, 141 142{{0,4}, {1,3}, {5,4}, {4,2}}, 143{{4,5}, {2,4}, {4,0}, {3,1}}, 144 145{{0,2}, {1,5}, {3,2}, {4,1}}, 146{{2,3}, {1,4}, {2,0}, {5,1}}, 147 148{{0,2}, {2,3}, {5,1}, {3,2}}, 149{{1,5}, {2,3}, {2,0}, {3,2}}, 150 151{{2,6}, {4,5}, {1,0}, {6,1}}, 152{{0,1}, {1,6}, {6,2}, {5,4}}, 153 154{{0,1}, {1,2}, {6,5}, {5,4}}, 155{{5,6}, {4,5}, {1,0}, {2,1}}, 156 157{{2.5119999999999996, 1.5710000000000002}, {2.6399999999999983, 1.6599999999999997}, {2.8000000000000007, 1.8000000000000003}, {3, 2}}, 158{{2.4181876227114887, 1.9849772580462195}, {2.8269904869227211, 2.009330650246834}, {3.2004679292461624, 1.9942047174679169}, {3.4986199496818058, 2.0035994597094731}}, 159 160{{2,3}, {1,4}, {1,0}, {6,0}}, 161{{0,1}, {0,6}, {3,2}, {4,1}}, 162 163{{0,2}, {1,5}, {1,0}, {6,1}}, 164{{0,1}, {1,6}, {2,0}, {5,1}}, 165 166{{0,1}, {1,5}, {2,1}, {4,0}}, 167{{1,2}, {0,4}, {1,0}, {5,1}}, 168 169{{0,1}, {3,5}, {2,1}, {3,1}}, 170{{1,2}, {1,3}, {1,0}, {5,3}}, 171 172{{0,1}, {2,5}, {6,0}, {5,3}}, 173{{0,6}, {3,5}, {1,0}, {5,2}}, 174 175{{0,1}, {3,6}, {1,0}, {5,2}}, 176{{0,1}, {2,5}, {1,0}, {6,3}}, 177 178{{1,2},{5,6},{1,0},{1,0}}, 179{{0,1},{0,1},{2,1},{6,5}}, 180 181{{0,6},{1,2},{1,0},{1,0}}, 182{{0,1},{0,1},{6,0},{2,1}}, 183 184{{0,2},{0,1},{3,0},{1,0}}, 185{{0,3},{0,1},{2,0},{1,0}}, 186}; 187 188const size_t newTestSetCount = sizeof(newTestSet) / sizeof(newTestSet[0]); 189 190#if 0 191static void oneOff(const Cubic& cubic1, const Cubic& cubic2) { 192 SkTDArray<Quadratic> quads1; 193 cubic_to_quadratics(cubic1, calcPrecision(cubic1), quads1); 194#if SHOW_ORIGINAL 195 SkDebugf("computed quadratics given\n"); 196 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}}, {%1.9g,%1.9g}},\n", 197 cubic1[0].x, cubic1[0].y, cubic1[1].x, cubic1[1].y, 198 cubic1[2].x, cubic1[2].y, cubic1[3].x, cubic1[3].y)); 199 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}}, {%1.9g,%1.9g}},\n", 200 cubic2[0].x, cubic2[0].y, cubic2[1].x, cubic2[1].y, 201 cubic2[2].x, cubic2[2].y, cubic2[3].x, cubic2[3].y)); 202#endif 203#if ONE_OFF_DEBUG 204 SkDebugf("computed quadratics set 1\n"); 205 for (int index = 0; index < quads1.count(); ++index) { 206 const Quadratic& q = quads1[index]; 207 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", q[0].x, q[0].y, 208 q[1].x, q[1].y, q[2].x, q[2].y); 209 } 210#endif 211 SkTDArray<Quadratic> quads2; 212 cubic_to_quadratics(cubic2, calcPrecision(cubic2), quads2); 213#if ONE_OFF_DEBUG 214 SkDebugf("computed quadratics set 2\n"); 215 for (int index = 0; index < quads2.count(); ++index) { 216 const Quadratic& q = quads2[index]; 217 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", q[0].x, q[0].y, 218 q[1].x, q[1].y, q[2].x, q[2].y); 219 } 220#endif 221 Intersections intersections2, intersections3; 222 intersect2(cubic1, cubic2, intersections2); 223 intersect3(cubic1, cubic2, intersections3); 224 int pt1, pt2, pt3; 225 bool found; 226 double tt1, tt2, last = -1; 227 _Point xy1, xy2; 228 for (pt1 = 0; pt1 < intersections2.used(); ++pt1) { 229 tt1 = intersections2.fT[0][pt1]; 230 SkASSERT(!approximately_equal(last, tt1)); 231 last = tt1; 232 xy_at_t(cubic1, tt1, xy1.x, xy1.y); 233 pt2 = intersections2.fFlip ? intersections2.used() - pt1 - 1 : pt1; 234 tt2 = intersections2.fT[1][pt2]; 235 xy_at_t(cubic2, tt2, xy2.x, xy2.y); 236#if ONE_OFF_DEBUG 237 SkDebugf("%s t1=%1.9g (%1.9g, %1.9g) (%1.9g, %1.9g) (%1.9g, %1.9g) t2=%1.9g\n", 238 __FUNCTION__, tt1, xy1.x, xy1.y, intersections2.fPt[pt1].x, 239 intersections2.fPt[pt1].y, xy2.x, xy2.y, tt2); 240#endif 241 SkASSERT(xy1.approximatelyEqual(xy2)); 242#if SK_DEBUG 243 found = false; 244 for (pt3 = 0; pt3 < intersections3.used(); ++pt3) { 245 if (roughly_equal(tt1, intersections3.fT[0][pt3])) { 246 found = true; 247 break; 248 } 249 } 250 SkASSERT(found); 251#endif 252 } 253 last = -1; 254 for (pt3 = 0; pt3 < intersections3.used(); ++pt3) { 255 found = false; 256 double tt3 = intersections3.fT[0][pt3]; 257 SkASSERT(!approximately_equal(last, tt3)); 258 last = tt3; 259 for (pt1 = 0; pt1 < intersections2.used(); ++pt1) { 260 if (approximately_equal(tt3, intersections2.fT[0][pt1])) { 261 found = true; 262 break; 263 } 264 } 265 if (!found) { 266 tt1 = intersections3.fT[0][pt3]; 267 xy_at_t(cubic1, tt1, xy1.x, xy1.y); 268 pt2 = intersections3.fFlip ? intersections3.used() - pt3 - 1 : pt3; 269 tt2 = intersections3.fT[1][pt2]; 270 xy_at_t(cubic2, tt2, xy2.x, xy2.y); 271 #if ONE_OFF_DEBUG 272 SkDebugf("%s t1=%1.9g (%1.9g, %1.9g) (%1.9g, %1.9g) (%1.9g, %1.9g) t2=%1.9g\n", 273 __FUNCTION__, tt1, xy1.x, xy1.y, intersections3.fPt[pt1].x, 274 intersections3.fPt[pt1].y, xy2.x, xy2.y, tt2); 275 #endif 276 SkASSERT(xy1.approximatelyEqual(xy2)); 277 SkDebugf("%s missing in intersect2\n", __FUNCTION__); 278 } 279 } 280} 281#endif 282 283static void oneOff3(const Cubic& cubic1, const Cubic& cubic2) { 284#if ONE_OFF_DEBUG 285 SkDebugf("computed quadratics given\n"); 286 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", 287 cubic1[0].x, cubic1[0].y, cubic1[1].x, cubic1[1].y, 288 cubic1[2].x, cubic1[2].y, cubic1[3].x, cubic1[3].y); 289 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", 290 cubic2[0].x, cubic2[0].y, cubic2[1].x, cubic2[1].y, 291 cubic2[2].x, cubic2[2].y, cubic2[3].x, cubic2[3].y); 292#endif 293 SkTDArray<Quadratic> quads1; 294 cubic_to_quadratics(cubic1, calcPrecision(cubic1), quads1); 295#if ONE_OFF_DEBUG 296 SkDebugf("computed quadratics set 1\n"); 297 for (int index = 0; index < quads1.count(); ++index) { 298 const Quadratic& q = quads1[index]; 299 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", q[0].x, q[0].y, 300 q[1].x, q[1].y, q[2].x, q[2].y); 301 } 302#endif 303 SkTDArray<Quadratic> quads2; 304 cubic_to_quadratics(cubic2, calcPrecision(cubic2), quads2); 305#if ONE_OFF_DEBUG 306 SkDebugf("computed quadratics set 2\n"); 307 for (int index = 0; index < quads2.count(); ++index) { 308 const Quadratic& q = quads2[index]; 309 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", q[0].x, q[0].y, 310 q[1].x, q[1].y, q[2].x, q[2].y); 311 } 312#endif 313 Intersections intersections3; 314 intersect3(cubic1, cubic2, intersections3); 315 int pt2, pt3; 316 double tt1, tt2, last = -1; 317 _Point xy1, xy2; 318 for (pt3 = 0; pt3 < intersections3.used(); ++pt3) { 319 double tt3 = intersections3.fT[0][pt3]; 320 // SkASSERT(!approximately_equal(last, tt3)); 321 last = tt3; 322 tt1 = intersections3.fT[0][pt3]; 323 xy_at_t(cubic1, tt1, xy1.x, xy1.y); 324 pt2 = intersections3.fFlip ? intersections3.used() - pt3 - 1 : pt3; 325 tt2 = intersections3.fT[1][pt2]; 326 xy_at_t(cubic2, tt2, xy2.x, xy2.y); 327#if ONE_OFF_DEBUG 328 SkDebugf("%s t1=%1.9g (%1.9g, %1.9g) (%1.9g, %1.9g) (%1.9g, %1.9g) t2=%1.9g\n", 329 __FUNCTION__, tt1, xy1.x, xy1.y, intersections3.fPt[pt3].x, 330 intersections3.fPt[pt3].y, xy2.x, xy2.y, tt2); 331#endif 332 SkASSERT(xy1.approximatelyEqual(xy2)); 333 } 334} 335 336#if 0 337static int fails[][2] = { {0, 23}, // fails in intersect2 recursing 338 {2, 7}, // answers differ, but neither is correct ('3' is closer) 339 {3, 26}, // fails in intersect2 recursing 340 {4, 9}, // fails in intersect2 recursing 341 {4, 10}, // fails in intersect2 recursing 342 {10, 17}, // fails in intersect2 recursing 343 {12, 14}, // loops indefinitely 344 {12, 21}, // fails in intersect2 recursing 345 {13, 21}, // fails in intersect2 recursing 346 {14, 21}, // fails in intersect2 recursing 347 {17, 25}, // fails in intersect2 recursing 348 {23, 25}, // fails in intersect2 recursing 349}; 350 351static int failCount = sizeof(fails) / sizeof(fails[0]); 352#endif 353 354static void oneOff(int outer, int inner) { 355 const Cubic& cubic1 = testSet[outer]; 356 const Cubic& cubic2 = testSet[inner]; 357#if 0 358 bool failing = false; 359 for (int i = 0; i < failCount; ++i) { 360 if ((fails[i][0] == outer && fails[i][1] == inner) 361 || (fails[i][1] == outer && fails[i][0] == inner)) { 362 failing = true; 363 break; 364 } 365 } 366 if (!failing) { 367 oneOff(cubic1, cubic2); 368 } else { 369#endif 370 oneOff3(cubic1, cubic2); 371// } 372} 373 374void CubicIntersection_OneOffTest() { 375 oneOff(12, 14); 376} 377 378static void newOneOff(int outer, int inner) { 379 const Cubic& cubic1 = newTestSet[outer]; 380 const Cubic& cubic2 = newTestSet[inner]; 381 oneOff3(cubic1, cubic2); 382} 383 384void CubicIntersection_NewOneOffTest() { 385 newOneOff(0, 1); 386} 387 388static void oneOffTests() { 389 for (size_t outer = 0; outer < testSetCount - 1; ++outer) { 390 for (size_t inner = outer + 1; inner < testSetCount; ++inner) { 391 oneOff(outer, inner); 392 } 393 } 394} 395 396void CubicIntersection_OneOffTests() { 397 oneOffTests(); 398} 399 400#define DEBUG_CRASH 0 401 402class CubicChopper { 403public: 404 405// only finds one intersection 406CubicChopper(const Cubic& c1, const Cubic& c2) 407 : cubic1(c1) 408 , cubic2(c2) 409 , depth(0) { 410} 411 412bool intersect(double minT1, double maxT1, double minT2, double maxT2) { 413 Cubic sub1, sub2; 414 // FIXME: carry last subdivide and reduceOrder result with cubic 415 sub_divide(cubic1, minT1, maxT1, sub1); 416 sub_divide(cubic2, minT2, maxT2, sub2); 417 Intersections i; 418 intersect3(sub1, sub2, i); 419 if (i.used() == 0) { 420 return false; 421 } 422 double x1, y1, x2, y2; 423 t1 = minT1 + i.fT[0][0] * (maxT1 - minT1); 424 t2 = minT2 + i.fT[1][0] * (maxT2 - minT2); 425 xy_at_t(cubic1, t1, x1, y1); 426 xy_at_t(cubic2, t2, x2, y2); 427 if (AlmostEqualUlps(x1, x2) && AlmostEqualUlps(y1, y2)) { 428 return true; 429 } 430 double half1 = (minT1 + maxT1) / 2; 431 double half2 = (minT2 + maxT2) / 2; 432 ++depth; 433 bool result; 434 if (depth & 1) { 435 result = intersect(minT1, half1, minT2, maxT2) || intersect(half1, maxT1, minT2, maxT2) 436 || intersect(minT1, maxT1, minT2, half2) || intersect(minT1, maxT1, half2, maxT2); 437 } else { 438 result = intersect(minT1, maxT1, minT2, half2) || intersect(minT1, maxT1, half2, maxT2) 439 || intersect(minT1, half1, minT2, maxT2) || intersect(half1, maxT1, minT2, maxT2); 440 } 441 --depth; 442 return result; 443} 444 445const Cubic& cubic1; 446const Cubic& cubic2; 447double t1; 448double t2; 449int depth; 450}; 451 452#define TRY_OLD 0 // old way fails on test == 1 453 454void CubicIntersection_RandTestOld() { 455 srand(0); 456 const int tests = 1000000; // 10000000; 457 double largestFactor = DBL_MAX; 458 for (int test = 0; test < tests; ++test) { 459 Cubic cubic1, cubic2; 460 for (int i = 0; i < 4; ++i) { 461 cubic1[i].x = (double) rand() / RAND_MAX * 100; 462 cubic1[i].y = (double) rand() / RAND_MAX * 100; 463 cubic2[i].x = (double) rand() / RAND_MAX * 100; 464 cubic2[i].y = (double) rand() / RAND_MAX * 100; 465 } 466 if (test == 2513) { // the pair crosses three times, but the quadratic approximation 467 continue; // only sees one -- should be OK to ignore the other two? 468 } 469 if (test == 12932) { // this exposes a weakness when one cubic touches the other but 470 continue; // does not touch the quad approximation. Captured in qc.htm as cubic15 471 } 472 #if DEBUG_CRASH 473 char str[1024]; 474 sprintf(str, "{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n" 475 "{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n", 476 cubic1[0].x, cubic1[0].y, cubic1[1].x, cubic1[1].y, cubic1[2].x, cubic1[2].y, 477 cubic1[3].x, cubic1[3].y, 478 cubic2[0].x, cubic2[0].y, cubic2[1].x, cubic2[1].y, cubic2[2].x, cubic2[2].y, 479 cubic2[3].x, cubic2[3].y); 480 #endif 481 _Rect rect1, rect2; 482 rect1.setBounds(cubic1); 483 rect2.setBounds(cubic2); 484 bool boundsIntersect = rect1.left <= rect2.right && rect2.left <= rect2.right 485 && rect1.top <= rect2.bottom && rect2.top <= rect1.bottom; 486 Intersections i1, i2; 487 #if TRY_OLD 488 bool oldIntersects = intersect(cubic1, cubic2, i1); 489 #else 490 bool oldIntersects = false; 491 #endif 492 if (test == -1) { 493 SkDebugf("ready...\n"); 494 } 495 bool newIntersects = intersect3(cubic1, cubic2, i2); 496 if (!boundsIntersect && (oldIntersects || newIntersects)) { 497 #if DEBUG_CRASH 498 SkDebugf("%s %d unexpected intersection boundsIntersect=%d oldIntersects=%d" 499 " newIntersects=%d\n%s %s\n", __FUNCTION__, test, boundsIntersect, 500 oldIntersects, newIntersects, __FUNCTION__, str); 501 #endif 502 SkASSERT(0); 503 } 504 if (oldIntersects && !newIntersects) { 505 #if DEBUG_CRASH 506 SkDebugf("%s %d missing intersection oldIntersects=%d newIntersects=%d\n%s %s\n", 507 __FUNCTION__, test, oldIntersects, newIntersects, __FUNCTION__, str); 508 #endif 509 SkASSERT(0); 510 } 511 if (!oldIntersects && !newIntersects) { 512 continue; 513 } 514 if (i2.used() > 1) { 515 continue; 516 // just look at single intercepts for simplicity 517 } 518 Intersections self1, self2; // self-intersect checks 519 if (intersect(cubic1, self1)) { 520 continue; 521 } 522 if (intersect(cubic2, self2)) { 523 continue; 524 } 525 // binary search for range necessary to enclose real intersection 526 CubicChopper c(cubic1, cubic2); 527 bool result = c.intersect(0, 1, 0, 1); 528 if (!result) { 529 // FIXME: a failure here probably means that a core routine used by CubicChopper is failing 530 continue; 531 } 532 double delta1 = fabs(c.t1 - i2.fT[0][0]); 533 double delta2 = fabs(c.t2 - i2.fT[1][0]); 534 double calc1 = calcPrecision(cubic1); 535 double calc2 = calcPrecision(cubic2); 536 double factor1 = calc1 / delta1; 537 double factor2 = calc2 / delta2; 538 SkDebugf("%s %d calc1=%1.9g delta1=%1.9g factor1=%1.9g calc2=%1.9g delta2=%1.9g" 539 " factor2=%1.9g\n", __FUNCTION__, test, 540 calc1, delta1, factor1, calc2, delta2, factor2); 541 if (factor1 < largestFactor) { 542 SkDebugf("WE HAVE A WINNER! %1.9g\n", factor1); 543 #if DEBUG_CRASH 544 SkDebugf("%s\n", str); 545 #endif 546 oneOff3(cubic1, cubic2); 547 largestFactor = factor1; 548 } 549 if (factor2 < largestFactor) { 550 SkDebugf("WE HAVE A WINNER! %1.9g\n", factor2); 551 #if DEBUG_CRASH 552 SkDebugf("%s\n", str); 553 #endif 554 oneOff3(cubic1, cubic2); 555 largestFactor = factor2; 556 } 557 } 558} 559 560void CubicIntersection_RandTest() { 561 srand(0); 562 const int tests = 10000000; 563 for (int test = 0; test < tests; ++test) { 564 Cubic cubic1, cubic2; 565 for (int i = 0; i < 4; ++i) { 566 cubic1[i].x = (double) rand() / RAND_MAX * 100; 567 cubic1[i].y = (double) rand() / RAND_MAX * 100; 568 cubic2[i].x = (double) rand() / RAND_MAX * 100; 569 cubic2[i].y = (double) rand() / RAND_MAX * 100; 570 } 571 #if DEBUG_CRASH 572 char str[1024]; 573 sprintf(str, "{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n" 574 "{{%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}, {%1.9g, %1.9g}},\n", 575 cubic1[0].x, cubic1[0].y, cubic1[1].x, cubic1[1].y, cubic1[2].x, cubic1[2].y, 576 cubic1[3].x, cubic1[3].y, 577 cubic2[0].x, cubic2[0].y, cubic2[1].x, cubic2[1].y, cubic2[2].x, cubic2[2].y, 578 cubic2[3].x, cubic2[3].y); 579 #endif 580 _Rect rect1, rect2; 581 rect1.setBounds(cubic1); 582 rect2.setBounds(cubic2); 583 bool boundsIntersect = rect1.left <= rect2.right && rect2.left <= rect2.right 584 && rect1.top <= rect2.bottom && rect2.top <= rect1.bottom; 585 if (test == -1) { 586 SkDebugf("ready...\n"); 587 } 588 Intersections intersections2; 589 bool newIntersects = intersect3(cubic1, cubic2, intersections2); 590 if (!boundsIntersect && newIntersects) { 591 #if DEBUG_CRASH 592 SkDebugf("%s %d unexpected intersection boundsIntersect=%d " 593 " newIntersects=%d\n%s %s\n", __FUNCTION__, test, boundsIntersect, 594 newIntersects, __FUNCTION__, str); 595 #endif 596 SkASSERT(0); 597 } 598 for (int pt = 0; pt < intersections2.used(); ++pt) { 599 double tt1 = intersections2.fT[0][pt]; 600 _Point xy1, xy2; 601 xy_at_t(cubic1, tt1, xy1.x, xy1.y); 602 int pt2 = intersections2.fFlip ? intersections2.used() - pt - 1 : pt; 603 double tt2 = intersections2.fT[1][pt2]; 604 xy_at_t(cubic2, tt2, xy2.x, xy2.y); 605 #if 0 606 SkDebugf("%s t1=%1.9g (%1.9g, %1.9g) (%1.9g, %1.9g) t2=%1.9g\n", __FUNCTION__, 607 tt1, xy1.x, xy1.y, xy2.x, xy2.y, tt2); 608 #endif 609 SkASSERT(xy1.approximatelyEqual(xy2)); 610 } 611 } 612} 613 614static void intersectionFinder(int index0, int index1, double t1Seed, double t2Seed, 615 double t1Step, double t2Step) { 616 const Cubic& cubic1 = newTestSet[index0]; 617 const Cubic& cubic2 = newTestSet[index1]; 618 _Point t1[3], t2[3]; 619 bool toggle = true; 620 do { 621 xy_at_t(cubic1, t1Seed - t1Step, t1[0].x, t1[0].y); 622 xy_at_t(cubic1, t1Seed, t1[1].x, t1[1].y); 623 xy_at_t(cubic1, t1Seed + t1Step, t1[2].x, t1[2].y); 624 xy_at_t(cubic2, t2Seed - t2Step, t2[0].x, t2[0].y); 625 xy_at_t(cubic2, t2Seed, t2[1].x, t2[1].y); 626 xy_at_t(cubic2, t2Seed + t2Step, t2[2].x, t2[2].y); 627 double dist[3][3]; 628 dist[1][1] = t1[1].distance(t2[1]); 629 int best_i = 1, best_j = 1; 630 for (int i = 0; i < 3; ++i) { 631 for (int j = 0; j < 3; ++j) { 632 if (i == 1 && j == 1) { 633 continue; 634 } 635 dist[i][j] = t1[i].distance(t2[j]); 636 if (dist[best_i][best_j] > dist[i][j]) { 637 best_i = i; 638 best_j = j; 639 } 640 } 641 } 642 if (best_i == 0) { 643 t1Seed -= t1Step; 644 } else if (best_i == 2) { 645 t1Seed += t1Step; 646 } 647 if (best_j == 0) { 648 t2Seed -= t2Step; 649 } else if (best_j == 2) { 650 t2Seed += t2Step; 651 } 652 if (best_i == 1 && best_j == 1) { 653 if ((toggle ^= true)) { 654 t1Step /= 2; 655 } else { 656 t2Step /= 2; 657 } 658 } 659 } while (!t1[1].approximatelyEqual(t2[1])); 660 t1Step = t2Step = 0.1; 661 double t10 = t1Seed - t1Step * 2; 662 double t12 = t1Seed + t1Step * 2; 663 double t20 = t2Seed - t2Step * 2; 664 double t22 = t2Seed + t2Step * 2; 665 _Point test; 666 while (!approximately_zero(t1Step)) { 667 xy_at_t(cubic1, t10, test.x, test.y); 668 t10 += t1[1].approximatelyEqual(test) ? -t1Step : t1Step; 669 t1Step /= 2; 670 } 671 t1Step = 0.1; 672 while (!approximately_zero(t1Step)) { 673 xy_at_t(cubic1, t12, test.x, test.y); 674 t12 -= t1[1].approximatelyEqual(test) ? -t1Step : t1Step; 675 t1Step /= 2; 676 } 677 while (!approximately_zero(t2Step)) { 678 xy_at_t(cubic2, t20, test.x, test.y); 679 t20 += t2[1].approximatelyEqual(test) ? -t2Step : t2Step; 680 t2Step /= 2; 681 } 682 t2Step = 0.1; 683 while (!approximately_zero(t2Step)) { 684 xy_at_t(cubic2, t22, test.x, test.y); 685 t22 -= t2[1].approximatelyEqual(test) ? -t2Step : t2Step; 686 t2Step /= 2; 687 } 688#if ONE_OFF_DEBUG 689 SkDebugf("%s t1=(%1.9g<%1.9g<%1.9g) t2=(%1.9g<%1.9g<%1.9g)\n", __FUNCTION__, 690 t10, t1Seed, t12, t20, t2Seed, t22); 691 _Point p10 = xy_at_t(cubic1, t10); 692 _Point p1Seed = xy_at_t(cubic1, t1Seed); 693 _Point p12 = xy_at_t(cubic1, t12); 694 SkDebugf("%s p1=(%1.9g,%1.9g)<(%1.9g,%1.9g)<(%1.9g,%1.9g)\n", __FUNCTION__, 695 p10.x, p10.y, p1Seed.x, p1Seed.y, p12.x, p12.y); 696 _Point p20 = xy_at_t(cubic2, t20); 697 _Point p2Seed = xy_at_t(cubic2, t2Seed); 698 _Point p22 = xy_at_t(cubic2, t22); 699 SkDebugf("%s p2=(%1.9g,%1.9g)<(%1.9g,%1.9g)<(%1.9g,%1.9g)\n", __FUNCTION__, 700 p20.x, p20.y, p2Seed.x, p2Seed.y, p22.x, p22.y); 701#endif 702} 703 704void CubicIntersection_IntersectionFinder() { 705 706 double t1Seed = 0.5; 707 double t2Seed = 0.3; 708 double t1Step = 0.1; 709 double t2Step = 0.1; 710 if (false) intersectionFinder(0, 1, t1Seed, t2Seed, t1Step, t2Step); 711 intersectionFinder(1, 0, .3, .5, t1Step, t2Step); 712} 713 714static void coincidentTest() { 715#if 0 716 Cubic cubic1 = {{0, 1}, {0, 2}, {1, 0}, {1, 0}}; 717 Cubic cubic2 = {{0, 1}, {0, 2}, {1, 0}, {6, 1}}; 718#endif 719} 720 721void CubicIntersection_SelfTest() { 722 const Cubic selfSet[] = { 723 {{0,2}, {2,3}, {5,1}, {3,2}}, 724 {{0,2}, {3,5}, {5,0}, {4,2}}, 725 {{3.34,8.98}, {1.95,10.27}, {3.76,7.65}, {4.96,10.64}}, 726 {{3.13,2.74}, {1.08,4.62}, {3.71,0.94}, {2.01,3.81}}, 727 {{6.71,3.14}, {7.99,2.75}, {8.27,1.96}, {6.35,3.57}}, 728 {{12.81,7.27}, {7.22,6.98}, {12.49,8.97}, {11.42,6.18}}, 729 }; 730 size_t selfSetCount = sizeof(selfSet) / sizeof(selfSet[0]); 731 for (size_t index = 0; index < selfSetCount; ++index) { 732 const Cubic& cubic = selfSet[index]; 733 #if ONE_OFF_DEBUG 734 int idx2; 735 double max[3]; 736 int ts = find_cubic_max_curvature(cubic, max); 737 for (idx2 = 0; idx2 < ts; ++idx2) { 738 SkDebugf("%s max[%d]=%1.9g (%1.9g, %1.9g)\n", __FUNCTION__, idx2, 739 max[idx2], xy_at_t(cubic, max[idx2]).x, xy_at_t(cubic, max[idx2]).y); 740 } 741 SkTDArray<double> ts1; 742 SkTDArray<Quadratic> quads1; 743 cubic_to_quadratics(cubic, calcPrecision(cubic), ts1); 744 for (idx2 = 0; idx2 < ts1.count(); ++idx2) { 745 SkDebugf("%s t[%d]=%1.9g\n", __FUNCTION__, idx2, ts1[idx2]); 746 } 747 cubic_to_quadratics(cubic, calcPrecision(cubic), quads1); 748 for (idx2 = 0; idx2 < quads1.count(); ++idx2) { 749 const Quadratic& q = quads1[idx2]; 750 SkDebugf(" {{%1.9g,%1.9g}, {%1.9g,%1.9g}, {%1.9g,%1.9g}},\n", 751 q[0].x, q[0].y, q[1].x, q[1].y, q[2].x, q[2].y); 752 } 753 SkDebugf("\n"); 754 #endif 755 Intersections i; 756 SkDEBUGCODE(int result = ) intersect(cubic, i); 757 SkASSERT(result == 1); 758 SkASSERT(i.used() == 1); 759 SkASSERT(!approximately_equal(i.fT[0][0], i.fT[1][0])); 760 _Point pt1 = xy_at_t(cubic, i.fT[0][0]); 761 _Point pt2 = xy_at_t(cubic, i.fT[1][0]); 762 SkASSERT(pt1.approximatelyEqual(pt2)); 763 } 764} 765 766void CubicIntersection_Test() { 767 oneOffTests(); 768 coincidentTest(); 769 standardTestCases(); 770} 771