diff --git a/src/path/Curve.js b/src/path/Curve.js index 96b48ca7..5646495c 100644 --- a/src/path/Curve.js +++ b/src/path/Curve.js @@ -264,27 +264,48 @@ var Curve = this.Curve = Base.extend(/** @lends Curve# */{ }, getCrossings: function(point, roots) { - // Implement the crossing number algorithm: + // Implementation of the crossing number algorithm: // http://en.wikipedia.org/wiki/Point_in_polygon // Solve the y-axis cubic polynomial for point.y and count all solutions // to the right of point.x as crossings. var vals = this.getValues(), count = Curve.solveCubic(vals, 1, point.y, roots), - crossings = 0; + crossings = 0, + tolerance = /*#=*/ Numerical.TOLERANCE; for (var i = 0; i < count; i++) { var t = roots[i]; - if (t >= 0 && t < 1 && Curve.evaluate(vals, t, true, 0).x > point.x) { - // If we're close to 0 and are not changing y-direction from the - // previous curve, do not count this root, as we're merely - // touching a tip. Passing 1 for Curve.evaluate()'s type means - // we're calculating tangents, and then check their y-slope for - // a change of direction: - if (t < /*#=*/ Numerical.TOLERANCE - && Curve.evaluate(this.getPrevious().getValues(), 1, true, 1).y - * Curve.evaluate(vals, t, true, 1).y - >= /*#=*/ Numerical.TOLERANCE) - continue; - crossings++; + if (t >= -tolerance && t < 1 - tolerance) { + var pt = Curve.evaluate(vals, t, true, 0); +/*#*/ if (options.debug) { + console.log(t, point.y, pt.y); + new Path.Circle({ + center: Curve.evaluate(vals, t, true, 0), + radius: 2, + strokeColor: 'red', + strokeWidth: 0.25 + }); +/*#*/ } + if (pt.x >= point.x - tolerance) { + // Passing 1 for Curve.evaluate()'s type calculates tangents. + var tangent = Curve.evaluate(vals, t, true, 1); + if ( + // Skip touching stationary points (tips), but if the + // actual point is on one, do not skip this solution! + Math.abs(pt.x - point.x) > tolerance + && ( + // Check derivate for stationary points + Math.abs(tangent.y) < tolerance + // If root is close to 0 and not changing vertical + // orientation from the previous curve, do not count + // this root, as it's touching a corner. + || t < tolerance + // Check the y-slope for a change of orientation + && tangent.y * Curve.evaluate( + this.getPrevious().getValues(), 1, true, 1).y + < tolerance)) + continue; + crossings++; + } } } return crossings; @@ -513,7 +534,7 @@ statics: { b = 3 * (c2 - c1) - c, a = p2 - p1 - c - b; return Numerical.solveCubic(a, b, c, p1 - val, roots, - /*#=*/ Numerical.TOLERANCE); + /*#=*/ Numerical.EPSILON); }, getParameterOf: function(v, x, y) { @@ -646,11 +667,13 @@ statics: { var bounds1 = Curve.getBounds(v1), bounds2 = Curve.getBounds(v2); /*#*/ if (options.debug) { - new Path.Rectangle(bounds1).set({ + new Path.Rectangle({ + rectangle: bounds1, strokeColor: 'green', strokeWidth: 0.1 }); - new Path.Rectangle(bounds2).set({ + new Path.Rectangle({ + rectangle: bounds2, strokeColor: 'red', strokeWidth: 0.1 }); @@ -660,11 +683,15 @@ statics: { if (Curve.isFlatEnough(v1, /*#=*/ Numerical.TOLERANCE) && Curve.isFlatEnough(v2, /*#=*/ Numerical.TOLERANCE)) { /*#*/ if (options.debug) { - new Path.Line(v1[0], v1[1], v1[6], v1[7]).set({ + new Path.Line({ + from: [v1[0], v1[1]], + to: [v1[6], v1[7]], strokeColor: 'green', strokeWidth: 0.1 }); - new Path.Line(v2[0], v2[1], v2[6], v2[7]).set({ + new Path.Line({ + from: [v2[0], v2[1]], + to: [v2[6], v2[7]], strokeColor: 'red', strokeWidth: 0.1 }); diff --git a/src/util/Numerical.js b/src/util/Numerical.js index a13f5f16..779f8b69 100644 --- a/src/util/Numerical.js +++ b/src/util/Numerical.js @@ -171,35 +171,38 @@ var Numerical = this.Numerical = new function() { d /= a; // Compute discriminants var bb = b * b, - p = 1 / 3 * (-1 / 3 * bb + c), - q = 1 / 2 * (2 / 27 * b * bb - 1 / 3 * b * c + d), + p = (bb - 3 * c) / 9, + q = (2 * bb * b - 9 * b * c + 27 * d) / 54, // Use Cardano's formula ppp = p * p * p, - D = q * q + ppp; + D = q * q - ppp; // Substitute x = y - b/3 to eliminate quadric term: x^3 +px + q = 0 b /= 3; if (abs(D) < tolerance) { - if (abs(q) < tolerance) { // One triple solution. - roots[0] = - b; - return 1; - } else { // One single and one double solution. - var u = cbrt(-q); - roots[0] = 2 * u - b; - roots[1] = - u - b; - return 2; - } + if (abs(q) < tolerance) { // One triple solution. + roots[0] = - b; + return 1; + } else { // One single and one double solution. + var sqp = sqrt(p), + snq = q < 0 ? -1 : 1; + roots[0] = -snq * 2 * sqp - b; + roots[1] = snq * sqp - b; + return 2; + } } else if (D < 0) { // Casus irreducibilis: three real solutions - var phi = 1 / 3 * Math.acos(-q / sqrt(-ppp)); - var t = 2 * sqrt(-p); - roots[0] = t * cos(phi) - b; - roots[1] = - t * cos(phi + PI / 3) - b; - roots[2] = - t * cos(phi - PI / 3) - b; - return 3; - } else { // One real solution - D = sqrt(D); - roots[0] = cbrt(D - q) - cbrt(D + q) - b; - return 1; + var sqp = sqrt(p), + phi = Math.acos(q / (sqp * sqp * sqp)) / 3, + o = 2 * PI / 3, + t = -2 * sqp; + roots[0] = t * cos(phi) - b; + roots[1] = t * cos(phi + o) - b; + roots[2] = t * cos(phi - o) - b; + return 3; } + // One real solution + var sqD = sqrt(D); + roots[0] = cbrt(sqD - q) - cbrt(sqD + q) - b; + return 1; } }; };