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 "SkIntersections.h" 8 #include "SkPathOpsCubic.h" 9 #include "SkPathOpsLine.h" 10 11 /* 12 Find the interection of a line and cubic by solving for valid t values. 13 14 Analogous to line-quadratic intersection, solve line-cubic intersection by 15 representing the cubic as: 16 x = a(1-t)^3 + 2b(1-t)^2t + c(1-t)t^2 + dt^3 17 y = e(1-t)^3 + 2f(1-t)^2t + g(1-t)t^2 + ht^3 18 and the line as: 19 y = i*x + j (if the line is more horizontal) 20 or: 21 x = i*y + j (if the line is more vertical) 22 23 Then using Mathematica, solve for the values of t where the cubic intersects the 24 line: 25 26 (in) Resultant[ 27 a*(1 - t)^3 + 3*b*(1 - t)^2*t + 3*c*(1 - t)*t^2 + d*t^3 - x, 28 e*(1 - t)^3 + 3*f*(1 - t)^2*t + 3*g*(1 - t)*t^2 + h*t^3 - i*x - j, x] 29 (out) -e + j + 30 3 e t - 3 f t - 31 3 e t^2 + 6 f t^2 - 3 g t^2 + 32 e t^3 - 3 f t^3 + 3 g t^3 - h t^3 + 33 i ( a - 34 3 a t + 3 b t + 35 3 a t^2 - 6 b t^2 + 3 c t^2 - 36 a t^3 + 3 b t^3 - 3 c t^3 + d t^3 ) 37 38 if i goes to infinity, we can rewrite the line in terms of x. Mathematica: 39 40 (in) Resultant[ 41 a*(1 - t)^3 + 3*b*(1 - t)^2*t + 3*c*(1 - t)*t^2 + d*t^3 - i*y - j, 42 e*(1 - t)^3 + 3*f*(1 - t)^2*t + 3*g*(1 - t)*t^2 + h*t^3 - y, y] 43 (out) a - j - 44 3 a t + 3 b t + 45 3 a t^2 - 6 b t^2 + 3 c t^2 - 46 a t^3 + 3 b t^3 - 3 c t^3 + d t^3 - 47 i ( e - 48 3 e t + 3 f t + 49 3 e t^2 - 6 f t^2 + 3 g t^2 - 50 e t^3 + 3 f t^3 - 3 g t^3 + h t^3 ) 51 52 Solving this with Mathematica produces an expression with hundreds of terms; 53 instead, use Numeric Solutions recipe to solve the cubic. 54 55 The near-horizontal case, in terms of: Ax^3 + Bx^2 + Cx + D == 0 56 A = (-(-e + 3*f - 3*g + h) + i*(-a + 3*b - 3*c + d) ) 57 B = 3*(-( e - 2*f + g ) + i*( a - 2*b + c ) ) 58 C = 3*(-(-e + f ) + i*(-a + b ) ) 59 D = (-( e ) + i*( a ) + j ) 60 61 The near-vertical case, in terms of: Ax^3 + Bx^2 + Cx + D == 0 62 A = ( (-a + 3*b - 3*c + d) - i*(-e + 3*f - 3*g + h) ) 63 B = 3*( ( a - 2*b + c ) - i*( e - 2*f + g ) ) 64 C = 3*( (-a + b ) - i*(-e + f ) ) 65 D = ( ( a ) - i*( e ) - j ) 66 67 For horizontal lines: 68 (in) Resultant[ 69 a*(1 - t)^3 + 3*b*(1 - t)^2*t + 3*c*(1 - t)*t^2 + d*t^3 - j, 70 e*(1 - t)^3 + 3*f*(1 - t)^2*t + 3*g*(1 - t)*t^2 + h*t^3 - y, y] 71 (out) e - j - 72 3 e t + 3 f t + 73 3 e t^2 - 6 f t^2 + 3 g t^2 - 74 e t^3 + 3 f t^3 - 3 g t^3 + h t^3 75 */ 76 77 class LineCubicIntersections { 78 public: 79 enum PinTPoint { 80 kPointUninitialized, 81 kPointInitialized 82 }; 83 84 LineCubicIntersections(const SkDCubic& c, const SkDLine& l, SkIntersections* i) 85 : fCubic(c) 86 , fLine(l) 87 , fIntersections(i) 88 , fAllowNear(true) { 89 } 90 91 void allowNear(bool allow) { 92 fAllowNear = allow; 93 } 94 95 // see parallel routine in line quadratic intersections 96 int intersectRay(double roots[3]) { 97 double adj = fLine[1].fX - fLine[0].fX; 98 double opp = fLine[1].fY - fLine[0].fY; 99 SkDCubic r; 100 for (int n = 0; n < 4; ++n) { 101 r[n].fX = (fCubic[n].fY - fLine[0].fY) * adj - (fCubic[n].fX - fLine[0].fX) * opp; 102 } 103 double A, B, C, D; 104 SkDCubic::Coefficients(&r[0].fX, &A, &B, &C, &D); 105 return SkDCubic::RootsValidT(A, B, C, D, roots); 106 } 107 108 int intersect() { 109 addExactEndPoints(); 110 double rootVals[3]; 111 int roots = intersectRay(rootVals); 112 for (int index = 0; index < roots; ++index) { 113 double cubicT = rootVals[index]; 114 double lineT = findLineT(cubicT); 115 SkDPoint pt; 116 if (pinTs(&cubicT, &lineT, &pt, kPointUninitialized)) { 117 #if ONE_OFF_DEBUG 118 SkDPoint cPt = fCubic.ptAtT(cubicT); 119 SkDebugf("%s pt=(%1.9g,%1.9g) cPt=(%1.9g,%1.9g)\n", __FUNCTION__, pt.fX, pt.fY, 120 cPt.fX, cPt.fY); 121 #endif 122 fIntersections->insert(cubicT, lineT, pt); 123 } 124 } 125 if (fAllowNear) { 126 addNearEndPoints(); 127 } 128 return fIntersections->used(); 129 } 130 131 int horizontalIntersect(double axisIntercept, double roots[3]) { 132 double A, B, C, D; 133 SkDCubic::Coefficients(&fCubic[0].fY, &A, &B, &C, &D); 134 D -= axisIntercept; 135 return SkDCubic::RootsValidT(A, B, C, D, roots); 136 } 137 138 int horizontalIntersect(double axisIntercept, double left, double right, bool flipped) { 139 addExactHorizontalEndPoints(left, right, axisIntercept); 140 double rootVals[3]; 141 int roots = horizontalIntersect(axisIntercept, rootVals); 142 for (int index = 0; index < roots; ++index) { 143 double cubicT = rootVals[index]; 144 SkDPoint pt = fCubic.ptAtT(cubicT); 145 double lineT = (pt.fX - left) / (right - left); 146 if (pinTs(&cubicT, &lineT, &pt, kPointInitialized)) { 147 fIntersections->insert(cubicT, lineT, pt); 148 } 149 } 150 if (fAllowNear) { 151 addNearHorizontalEndPoints(left, right, axisIntercept); 152 } 153 if (flipped) { 154 fIntersections->flip(); 155 } 156 return fIntersections->used(); 157 } 158 159 int verticalIntersect(double axisIntercept, double roots[3]) { 160 double A, B, C, D; 161 SkDCubic::Coefficients(&fCubic[0].fX, &A, &B, &C, &D); 162 D -= axisIntercept; 163 return SkDCubic::RootsValidT(A, B, C, D, roots); 164 } 165 166 int verticalIntersect(double axisIntercept, double top, double bottom, bool flipped) { 167 addExactVerticalEndPoints(top, bottom, axisIntercept); 168 double rootVals[3]; 169 int roots = verticalIntersect(axisIntercept, rootVals); 170 for (int index = 0; index < roots; ++index) { 171 double cubicT = rootVals[index]; 172 SkDPoint pt = fCubic.ptAtT(cubicT); 173 double lineT = (pt.fY - top) / (bottom - top); 174 if (pinTs(&cubicT, &lineT, &pt, kPointInitialized)) { 175 fIntersections->insert(cubicT, lineT, pt); 176 } 177 } 178 if (fAllowNear) { 179 addNearVerticalEndPoints(top, bottom, axisIntercept); 180 } 181 if (flipped) { 182 fIntersections->flip(); 183 } 184 return fIntersections->used(); 185 } 186 187 protected: 188 189 void addExactEndPoints() { 190 for (int cIndex = 0; cIndex < 4; cIndex += 3) { 191 double lineT = fLine.exactPoint(fCubic[cIndex]); 192 if (lineT < 0) { 193 continue; 194 } 195 double cubicT = (double) (cIndex >> 1); 196 fIntersections->insert(cubicT, lineT, fCubic[cIndex]); 197 } 198 } 199 200 void addNearEndPoints() { 201 for (int cIndex = 0; cIndex < 4; cIndex += 3) { 202 double cubicT = (double) (cIndex >> 1); 203 if (fIntersections->hasT(cubicT)) { 204 continue; 205 } 206 double lineT = fLine.nearPoint(fCubic[cIndex]); 207 if (lineT < 0) { 208 continue; 209 } 210 fIntersections->insert(cubicT, lineT, fCubic[cIndex]); 211 } 212 } 213 214 void addExactHorizontalEndPoints(double left, double right, double y) { 215 for (int cIndex = 0; cIndex < 4; cIndex += 3) { 216 double lineT = SkDLine::ExactPointH(fCubic[cIndex], left, right, y); 217 if (lineT < 0) { 218 continue; 219 } 220 double cubicT = (double) (cIndex >> 1); 221 fIntersections->insert(cubicT, lineT, fCubic[cIndex]); 222 } 223 } 224 225 void addNearHorizontalEndPoints(double left, double right, double y) { 226 for (int cIndex = 0; cIndex < 4; cIndex += 3) { 227 double cubicT = (double) (cIndex >> 1); 228 if (fIntersections->hasT(cubicT)) { 229 continue; 230 } 231 double lineT = SkDLine::NearPointH(fCubic[cIndex], left, right, y); 232 if (lineT < 0) { 233 continue; 234 } 235 fIntersections->insert(cubicT, lineT, fCubic[cIndex]); 236 } 237 // FIXME: see if line end is nearly on cubic 238 } 239 240 void addExactVerticalEndPoints(double top, double bottom, double x) { 241 for (int cIndex = 0; cIndex < 4; cIndex += 3) { 242 double lineT = SkDLine::ExactPointV(fCubic[cIndex], top, bottom, x); 243 if (lineT < 0) { 244 continue; 245 } 246 double cubicT = (double) (cIndex >> 1); 247 fIntersections->insert(cubicT, lineT, fCubic[cIndex]); 248 } 249 } 250 251 void addNearVerticalEndPoints(double top, double bottom, double x) { 252 for (int cIndex = 0; cIndex < 4; cIndex += 3) { 253 double cubicT = (double) (cIndex >> 1); 254 if (fIntersections->hasT(cubicT)) { 255 continue; 256 } 257 double lineT = SkDLine::NearPointV(fCubic[cIndex], top, bottom, x); 258 if (lineT < 0) { 259 continue; 260 } 261 fIntersections->insert(cubicT, lineT, fCubic[cIndex]); 262 } 263 // FIXME: see if line end is nearly on cubic 264 } 265 266 double findLineT(double t) { 267 SkDPoint xy = fCubic.ptAtT(t); 268 double dx = fLine[1].fX - fLine[0].fX; 269 double dy = fLine[1].fY - fLine[0].fY; 270 if (fabs(dx) > fabs(dy)) { 271 return (xy.fX - fLine[0].fX) / dx; 272 } 273 return (xy.fY - fLine[0].fY) / dy; 274 } 275 276 bool pinTs(double* cubicT, double* lineT, SkDPoint* pt, PinTPoint ptSet) { 277 if (!approximately_one_or_less(*lineT)) { 278 return false; 279 } 280 if (!approximately_zero_or_more(*lineT)) { 281 return false; 282 } 283 double cT = *cubicT = SkPinT(*cubicT); 284 double lT = *lineT = SkPinT(*lineT); 285 if (lT == 0 || lT == 1 || (ptSet == kPointUninitialized && cT != 0 && cT != 1)) { 286 *pt = fLine.ptAtT(lT); 287 } else if (ptSet == kPointUninitialized) { 288 *pt = fCubic.ptAtT(cT); 289 } 290 return true; 291 } 292 293 private: 294 const SkDCubic& fCubic; 295 const SkDLine& fLine; 296 SkIntersections* fIntersections; 297 bool fAllowNear; 298 }; 299 300 int SkIntersections::horizontal(const SkDCubic& cubic, double left, double right, double y, 301 bool flipped) { 302 SkDLine line = {{{ left, y }, { right, y }}}; 303 LineCubicIntersections c(cubic, line, this); 304 return c.horizontalIntersect(y, left, right, flipped); 305 } 306 307 int SkIntersections::vertical(const SkDCubic& cubic, double top, double bottom, double x, 308 bool flipped) { 309 SkDLine line = {{{ x, top }, { x, bottom }}}; 310 LineCubicIntersections c(cubic, line, this); 311 return c.verticalIntersect(x, top, bottom, flipped); 312 } 313 314 int SkIntersections::intersect(const SkDCubic& cubic, const SkDLine& line) { 315 LineCubicIntersections c(cubic, line, this); 316 c.allowNear(fAllowNear); 317 return c.intersect(); 318 } 319 320 int SkIntersections::intersectRay(const SkDCubic& cubic, const SkDLine& line) { 321 LineCubicIntersections c(cubic, line, this); 322 fUsed = c.intersectRay(fT[0]); 323 for (int index = 0; index < fUsed; ++index) { 324 fPt[index] = cubic.ptAtT(fT[0][index]); 325 } 326 return fUsed; 327 } 328
