scratchfoundation/paper.js
1
0
Fork 0
mirror of https://github.com/scratchfoundation/paper.js.git synced 2025-01-04 03:45:58 -05:00
Code Issues Projects Releases Packages Wiki Activity Actions

Improve precision of Numerical.solveCubic() and fix issues in Curve.getCrossings().

Browse source 
Closes #202.
This commit is contained in:
Jürg Lehni 2013-04-20 19:14:19 -07:00
parent 461def5383
commit 14aa8e5dea
2 changed files with 71 additions and 41 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

65
src/path/Curve.js
View file

@ -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
});

47
src/util/Numerical.js
View file

@ -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;
}
};
};