CubicCurve2D backports - review
Dr Andrew John Hughes
ahughes at redhat.com
Tue Feb 1 07:30:36 PST 2011
On 18:32 Mon 31 Jan , Denis Lila wrote:
> Hi.
>
> I'd like to backport the three attached patches.
>
>
> hg diff:
>
> diff -r cd6310f10fab ChangeLog
> --- a/ChangeLog Sun Jan 30 00:45:18 2011 +0000
> +++ b/ChangeLog Mon Jan 31 13:28:19 2011 -0500
> @@ -1,3 +1,11 @@
> +2011-01-31 Denis Lila <dlila at redhat.com>
> +
> + * NEWS: Update with the 3 backports
> + * Makefile.am (ICEDTEA_PATCHES): Add the patches
> + * patches/openjdk/4493128-CubicCurve2D.patch: New file.
> + * patches/openjdk/4645692-CubicCurve2D.solveCubic.patch: Likewise.
> + * patches/openjdk/4724552-CubicCurve2D.patch: Likewise.
> +
> 2011-01-29 Andrew John Hughes <ahughes at redhat.com>
>
> * patches/gcc-suffix.patch,
> diff -r cd6310f10fab Makefile.am
> --- a/Makefile.am Sun Jan 30 00:45:18 2011 +0000
> +++ b/Makefile.am Mon Jan 31 13:28:19 2011 -0500
> @@ -275,7 +275,10 @@
> patches/jtreg-international-fonts-styles.patch \
> patches/openjdk/6736649-text_bearings.patch \
> patches/openjdk/6797139-jbutton_truncation.patch \
> - patches/openjdk/6883341-text_bearing_exception.patch
> + patches/openjdk/6883341-text_bearing_exception.patch \
> + patches/openjdk/4724552-CubicCurve2D.patch \
> + patches/openjdk/4493128-CubicCurve2D.patch \
> + patches/openjdk/4645692-CubicCurve2D.solveCubic.patch
>
> if !WITH_ALT_HSBUILD
> ICEDTEA_PATCHES += \
> diff -r cd6310f10fab NEWS
> --- a/NEWS Sun Jan 30 00:45:18 2011 +0000
> +++ b/NEWS Mon Jan 31 13:28:19 2011 -0500
> @@ -403,6 +403,9 @@
> - S6736649: test/closed/javax/swing/JMenuItem/6458123/ManualBug6458123.java fails on Linux
> - S6797139: JButton title is truncating for some strings irrespective of preferred size.
> - S6883341: SWAT: jdk7-b72 swat build(2009-09-17) threw exceptions when running Java2D demo by clicking Paint ta
> + - S4493128: CubicCurve2D intersects method fails
> + - S4724552: CubicCurve2D.contains(Rectangle2D) returns true when only partially contained.
> + - S4645692: solveCubic does not return all solutions.
> * Bug fixes
> - RH647157, RH582455: Update fontconfig files for rhel 6
> - RH661505: JPEGs with sRGB IEC61966-2.1 color profiles have wrong colors
>
> ...(and here should be the diff of the patches, but those are attached).
>
Would be easier to just attach the result of hg diff which would include the patches :-)
Approved.
> Thank you,
> Denis.
> --- openjdk/jdk/src/share/classes/java/awt/geom/CubicCurve2D.java.old 2011-01-31 12:38:24.340733697 -0500
> +++ openjdk/jdk/src/share/classes/java/awt/geom/CubicCurve2D.java 2011-01-31 12:54:27.462601773 -0500
> @@ -1387,203 +1387,13 @@
> return false;
> }
>
> - // Trivially accept if either endpoint is inside the rectangle
> - // (not on its border since it may end there and not go inside)
> - // Record where they lie with respect to the rectangle.
> - // -1 => left, 0 => inside, 1 => right
> - double x1 = getX1();
> - double y1 = getY1();
> - int x1tag = getTag(x1, x, x+w);
> - int y1tag = getTag(y1, y, y+h);
> - if (x1tag == INSIDE && y1tag == INSIDE) {
> - return true;
> - }
> - double x2 = getX2();
> - double y2 = getY2();
> - int x2tag = getTag(x2, x, x+w);
> - int y2tag = getTag(y2, y, y+h);
> - if (x2tag == INSIDE && y2tag == INSIDE) {
> - return true;
> - }
> -
> - double ctrlx1 = getCtrlX1();
> - double ctrly1 = getCtrlY1();
> - double ctrlx2 = getCtrlX2();
> - double ctrly2 = getCtrlY2();
> - int ctrlx1tag = getTag(ctrlx1, x, x+w);
> - int ctrly1tag = getTag(ctrly1, y, y+h);
> - int ctrlx2tag = getTag(ctrlx2, x, x+w);
> - int ctrly2tag = getTag(ctrly2, y, y+h);
> -
> - // Trivially reject if all points are entirely to one side of
> - // the rectangle.
> - if (x1tag < INSIDE && x2tag < INSIDE &&
> - ctrlx1tag < INSIDE && ctrlx2tag < INSIDE)
> - {
> - return false; // All points left
> - }
> - if (y1tag < INSIDE && y2tag < INSIDE &&
> - ctrly1tag < INSIDE && ctrly2tag < INSIDE)
> - {
> - return false; // All points above
> - }
> - if (x1tag > INSIDE && x2tag > INSIDE &&
> - ctrlx1tag > INSIDE && ctrlx2tag > INSIDE)
> - {
> - return false; // All points right
> - }
> - if (y1tag > INSIDE && y2tag > INSIDE &&
> - ctrly1tag > INSIDE && ctrly2tag > INSIDE)
> - {
> - return false; // All points below
> - }
> -
> - // Test for endpoints on the edge where either the segment
> - // or the curve is headed "inwards" from them
> - // Note: These tests are a superset of the fast endpoint tests
> - // above and thus repeat those tests, but take more time
> - // and cover more cases
> - if (inwards(x1tag, x2tag, ctrlx1tag) &&
> - inwards(y1tag, y2tag, ctrly1tag))
> - {
> - // First endpoint on border with either edge moving inside
> - return true;
> - }
> - if (inwards(x2tag, x1tag, ctrlx2tag) &&
> - inwards(y2tag, y1tag, ctrly2tag))
> - {
> - // Second endpoint on border with either edge moving inside
> - return true;
> - }
> -
> - // Trivially accept if endpoints span directly across the rectangle
> - boolean xoverlap = (x1tag * x2tag <= 0);
> - boolean yoverlap = (y1tag * y2tag <= 0);
> - if (x1tag == INSIDE && x2tag == INSIDE && yoverlap) {
> - return true;
> - }
> - if (y1tag == INSIDE && y2tag == INSIDE && xoverlap) {
> - return true;
> - }
> -
> - // We now know that both endpoints are outside the rectangle
> - // but the 4 points are not all on one side of the rectangle.
> - // Therefore the curve cannot be contained inside the rectangle,
> - // but the rectangle might be contained inside the curve, or
> - // the curve might intersect the boundary of the rectangle.
> -
> - double[] eqn = new double[4];
> - double[] res = new double[4];
> - if (!yoverlap) {
> - // Both y coordinates for the closing segment are above or
> - // below the rectangle which means that we can only intersect
> - // if the curve crosses the top (or bottom) of the rectangle
> - // in more than one place and if those crossing locations
> - // span the horizontal range of the rectangle.
> - fillEqn(eqn, (y1tag < INSIDE ? y : y+h), y1, ctrly1, ctrly2, y2);
> - int num = solveCubic(eqn, res);
> - num = evalCubic(res, num, true, true, null,
> - x1, ctrlx1, ctrlx2, x2);
> - // odd counts imply the crossing was out of [0,1] bounds
> - // otherwise there is no way for that part of the curve to
> - // "return" to meet its endpoint
> - return (num == 2 &&
> - getTag(res[0], x, x+w) * getTag(res[1], x, x+w) <= 0);
> - }
> -
> - // Y ranges overlap. Now we examine the X ranges
> - if (!xoverlap) {
> - // Both x coordinates for the closing segment are left of
> - // or right of the rectangle which means that we can only
> - // intersect if the curve crosses the left (or right) edge
> - // of the rectangle in more than one place and if those
> - // crossing locations span the vertical range of the rectangle.
> - fillEqn(eqn, (x1tag < INSIDE ? x : x+w), x1, ctrlx1, ctrlx2, x2);
> - int num = solveCubic(eqn, res);
> - num = evalCubic(res, num, true, true, null,
> - y1, ctrly1, ctrly2, y2);
> - // odd counts imply the crossing was out of [0,1] bounds
> - // otherwise there is no way for that part of the curve to
> - // "return" to meet its endpoint
> - return (num == 2 &&
> - getTag(res[0], y, y+h) * getTag(res[1], y, y+h) <= 0);
> - }
> -
> - // The X and Y ranges of the endpoints overlap the X and Y
> - // ranges of the rectangle, now find out how the endpoint
> - // line segment intersects the Y range of the rectangle
> - double dx = x2 - x1;
> - double dy = y2 - y1;
> - double k = y2 * x1 - x2 * y1;
> - int c1tag, c2tag;
> - if (y1tag == INSIDE) {
> - c1tag = x1tag;
> - } else {
> - c1tag = getTag((k + dx * (y1tag < INSIDE ? y : y+h)) / dy, x, x+w);
> - }
> - if (y2tag == INSIDE) {
> - c2tag = x2tag;
> - } else {
> - c2tag = getTag((k + dx * (y2tag < INSIDE ? y : y+h)) / dy, x, x+w);
> - }
> - // If the part of the line segment that intersects the Y range
> - // of the rectangle crosses it horizontally - trivially accept
> - if (c1tag * c2tag <= 0) {
> - return true;
> - }
> -
> - // Now we know that both the X and Y ranges intersect and that
> - // the endpoint line segment does not directly cross the rectangle.
> - //
> - // We can almost treat this case like one of the cases above
> - // where both endpoints are to one side, except that we may
> - // get one or three intersections of the curve with the vertical
> - // side of the rectangle. This is because the endpoint segment
> - // accounts for the other intersection in an even pairing. Thus,
> - // with the endpoint crossing we end up with 2 or 4 total crossings.
> - //
> - // (Remember there is overlap in both the X and Y ranges which
> - // means that the segment itself must cross at least one vertical
> - // edge of the rectangle - in particular, the "near vertical side"
> - // - leaving an odd number of intersections for the curve.)
> - //
> - // Now we calculate the y tags of all the intersections on the
> - // "near vertical side" of the rectangle. We will have one with
> - // the endpoint segment, and one or three with the curve. If
> - // any pair of those vertical intersections overlap the Y range
> - // of the rectangle, we have an intersection. Otherwise, we don't.
> -
> - // c1tag = vertical intersection class of the endpoint segment
> - //
> - // Choose the y tag of the endpoint that was not on the same
> - // side of the rectangle as the subsegment calculated above.
> - // Note that we can "steal" the existing Y tag of that endpoint
> - // since it will be provably the same as the vertical intersection.
> - c1tag = ((c1tag * x1tag <= 0) ? y1tag : y2tag);
> -
> - // Now we have to calculate an array of solutions of the curve
> - // with the "near vertical side" of the rectangle. Then we
> - // need to sort the tags and do a pairwise range test to see
> - // if either of the pairs of crossings spans the Y range of
> - // the rectangle.
> - //
> - // Note that the c2tag can still tell us which vertical edge
> - // to test against.
> - fillEqn(eqn, (c2tag < INSIDE ? x : x+w), x1, ctrlx1, ctrlx2, x2);
> - int num = solveCubic(eqn, res);
> - num = evalCubic(res, num, true, true, null, y1, ctrly1, ctrly2, y2);
> -
> - // Now put all of the tags into a bucket and sort them. There
> - // is an intersection iff one of the pairs of tags "spans" the
> - // Y range of the rectangle.
> - int tags[] = new int[num+1];
> - for (int i = 0; i < num; i++) {
> - tags[i] = getTag(res[i], y, y+h);
> - }
> - tags[num] = c1tag;
> - Arrays.sort(tags);
> - return ((num >= 1 && tags[0] * tags[1] <= 0) ||
> - (num >= 3 && tags[2] * tags[3] <= 0));
> + int numCrossings = rectCrossings(x, y, w, h);
> + // the intended return value is
> + // numCrossings != 0 || numCrossings == Curve.RECT_INTERSECTS
> + // but if (numCrossings != 0) numCrossings == INTERSECTS won't matter
> + // and if !(numCrossings != 0) then numCrossings == 0, so
> + // numCrossings != RECT_INTERSECT
> + return numCrossings != 0;
> }
>
> /**
> --- openjdk/jdk/src/share/classes/java/awt/geom/CubicCurve2D.java.old 2011-01-31 12:54:27.462601773 -0500
> +++ openjdk/jdk/src/share/classes/java/awt/geom/CubicCurve2D.java 2011-01-31 12:58:12.813545289 -0500
> @@ -31,6 +31,10 @@
> import java.io.Serializable;
> import sun.awt.geom.Curve;
>
> +import static java.lang.Math.abs;
> +import static java.lang.Math.max;
> +import static java.lang.Math.ulp;
> +
> /**
> * The <code>CubicCurve2D</code> class defines a cubic parametric curve
> * segment in {@code (x,y)} coordinate space.
> @@ -1083,95 +1087,286 @@
> * @since 1.3
> */
> public static int solveCubic(double eqn[], double res[]) {
> - // From Numerical Recipes, 5.6, Quadratic and Cubic Equations
> - double d = eqn[3];
> - if (d == 0.0) {
> - // The cubic has degenerated to quadratic (or line or ...).
> + // From Graphics Gems:
> + // http://tog.acm.org/resources/GraphicsGems/gems/Roots3And4.c
> + final double d = eqn[3];
> + if (d == 0) {
> return QuadCurve2D.solveQuadratic(eqn, res);
> }
> - double a = eqn[2] / d;
> - double b = eqn[1] / d;
> - double c = eqn[0] / d;
> - int roots = 0;
> - double Q = (a * a - 3.0 * b) / 9.0;
> - double R = (2.0 * a * a * a - 9.0 * a * b + 27.0 * c) / 54.0;
> - double R2 = R * R;
> - double Q3 = Q * Q * Q;
> - a = a / 3.0;
> - if (R2 < Q3) {
> - double theta = Math.acos(R / Math.sqrt(Q3));
> - Q = -2.0 * Math.sqrt(Q);
> +
> + /* normal form: x^3 + Ax^2 + Bx + C = 0 */
> + final double A = eqn[2] / d;
> + final double B = eqn[1] / d;
> + final double C = eqn[0] / d;
> +
> +
> + // substitute x = y - A/3 to eliminate quadratic term:
> + // x^3 +Px + Q = 0
> + //
> + // Since we actually need P/3 and Q/2 for all of the
> + // calculations that follow, we will calculate
> + // p = P/3
> + // q = Q/2
> + // instead and use those values for simplicity of the code.
> + double sq_A = A * A;
> + double p = 1.0/3 * (-1.0/3 * sq_A + B);
> + double q = 1.0/2 * (2.0/27 * A * sq_A - 1.0/3 * A * B + C);
> +
> + /* use Cardano's formula */
> +
> + double cb_p = p * p * p;
> + double D = q * q + cb_p;
> +
> + final double sub = 1.0/3 * A;
> +
> + int num;
> + if (D < 0) { /* Casus irreducibilis: three real solutions */
> + // see: http://en.wikipedia.org/wiki/Cubic_function#Trigonometric_.28and_hyperbolic.29_method
> + double phi = 1.0/3 * Math.acos(-q / Math.sqrt(-cb_p));
> + double t = 2 * Math.sqrt(-p);
> +
> if (res == eqn) {
> - // Copy the eqn so that we don't clobber it with the
> - // roots. This is needed so that fixRoots can do its
> - // work with the original equation.
> - eqn = new double[4];
> - System.arraycopy(res, 0, eqn, 0, 4);
> + eqn = Arrays.copyOf(eqn, 4);
> }
> - res[roots++] = Q * Math.cos(theta / 3.0) - a;
> - res[roots++] = Q * Math.cos((theta + Math.PI * 2.0)/ 3.0) - a;
> - res[roots++] = Q * Math.cos((theta - Math.PI * 2.0)/ 3.0) - a;
> - fixRoots(res, eqn);
> - } else {
> - boolean neg = (R < 0.0);
> - double S = Math.sqrt(R2 - Q3);
> - if (neg) {
> - R = -R;
> +
> + res[ 0 ] = ( t * Math.cos(phi));
> + res[ 1 ] = (-t * Math.cos(phi + Math.PI / 3));
> + res[ 2 ] = (-t * Math.cos(phi - Math.PI / 3));
> + num = 3;
> +
> + for (int i = 0; i < num; ++i) {
> + res[ i ] -= sub;
> }
> - double A = Math.pow(R + S, 1.0 / 3.0);
> - if (!neg) {
> - A = -A;
> +
> + } else {
> + // Please see the comment in fixRoots marked 'XXX' before changing
> + // any of the code in this case.
> + double sqrt_D = Math.sqrt(D);
> + double u = Math.cbrt(sqrt_D - q);
> + double v = - Math.cbrt(sqrt_D + q);
> + double uv = u+v;
> +
> + num = 1;
> +
> + double err = 1200000000*ulp(abs(uv) + abs(sub));
> + if (iszero(D, err) || within(u, v, err)) {
> + if (res == eqn) {
> + eqn = Arrays.copyOf(eqn, 4);
> + }
> + res[1] = -(uv / 2) - sub;
> + num = 2;
> }
> - double B = (A == 0.0) ? 0.0 : (Q / A);
> - res[roots++] = (A + B) - a;
> + // this must be done after the potential Arrays.copyOf
> + res[ 0 ] = uv - sub;
> }
> - return roots;
> +
> + if (num > 1) { // num == 3 || num == 2
> + num = fixRoots(eqn, res, num);
> + }
> + if (num > 2 && (res[2] == res[1] || res[2] == res[0])) {
> + num--;
> + }
> + if (num > 1 && res[1] == res[0]) {
> + res[1] = res[--num]; // Copies res[2] to res[1] if needed
> + }
> + return num;
> + }
> +
> + // preconditions: eqn != res && eqn[3] != 0 && num > 1
> + // This method tries to improve the accuracy of the roots of eqn (which
> + // should be in res). It also might eliminate roots in res if it decideds
> + // that they're not real roots. It will not check for roots that the
> + // computation of res might have missed, so this method should only be
> + // used when the roots in res have been computed using an algorithm
> + // that never underestimates the number of roots (such as solveCubic above)
> + private static int fixRoots(double[] eqn, double[] res, int num) {
> + double[] intervals = {eqn[1], 2*eqn[2], 3*eqn[3]};
> + int critCount = QuadCurve2D.solveQuadratic(intervals, intervals);
> + if (critCount == 2 && intervals[0] == intervals[1]) {
> + critCount--;
> + }
> + if (critCount == 2 && intervals[0] > intervals[1]) {
> + double tmp = intervals[0];
> + intervals[0] = intervals[1];
> + intervals[1] = tmp;
> + }
> +
> + // below we use critCount to possibly filter out roots that shouldn't
> + // have been computed. We require that eqn[3] != 0, so eqn is a proper
> + // cubic, which means that its limits at -/+inf are -/+inf or +/-inf.
> + // Therefore, if critCount==2, the curve is shaped like a sideways S,
> + // and it could have 1-3 roots. If critCount==0 it is monotonic, and
> + // if critCount==1 it is monotonic with a single point where it is
> + // flat. In the last 2 cases there can only be 1 root. So in cases
> + // where num > 1 but critCount < 2, we eliminate all roots in res
> + // except one.
> +
> + if (num == 3) {
> + double xe = getRootUpperBound(eqn);
> + double x0 = -xe;
> +
> + Arrays.sort(res, 0, num);
> + if (critCount == 2) {
> + // this just tries to improve the accuracy of the computed
> + // roots using Newton's method.
> + res[0] = refineRootWithHint(eqn, x0, intervals[0], res[0]);
> + res[1] = refineRootWithHint(eqn, intervals[0], intervals[1], res[1]);
> + res[2] = refineRootWithHint(eqn, intervals[1], xe, res[2]);
> + return 3;
> + } else if (critCount == 1) {
> + // we only need fx0 and fxe for the sign of the polynomial
> + // at -inf and +inf respectively, so we don't need to do
> + // fx0 = solveEqn(eqn, 3, x0); fxe = solveEqn(eqn, 3, xe)
> + double fxe = eqn[3];
> + double fx0 = -fxe;
> +
> + double x1 = intervals[0];
> + double fx1 = solveEqn(eqn, 3, x1);
> +
> + // if critCount == 1 or critCount == 0, but num == 3 then
> + // something has gone wrong. This branch and the one below
> + // would ideally never execute, but if they do we can't know
> + // which of the computed roots is closest to the real root;
> + // therefore, we can't use refineRootWithHint. But even if
> + // we did know, being here most likely means that the
> + // curve is very flat close to two of the computed roots
> + // (or maybe even all three). This might make Newton's method
> + // fail altogether, which would be a pain to detect and fix.
> + // This is why we use a very stable bisection method.
> + if (oppositeSigns(fx0, fx1)) {
> + res[0] = bisectRootWithHint(eqn, x0, x1, res[0]);
> + } else if (oppositeSigns(fx1, fxe)) {
> + res[0] = bisectRootWithHint(eqn, x1, xe, res[2]);
> + } else /* fx1 must be 0 */ {
> + res[0] = x1;
> + }
> + // return 1
> + } else if (critCount == 0) {
> + res[0] = bisectRootWithHint(eqn, x0, xe, res[1]);
> + // return 1
> + }
> + } else if (num == 2 && critCount == 2) {
> + // XXX: here we assume that res[0] has better accuracy than res[1].
> + // This is true because this method is only used from solveCubic
> + // which puts in res[0] the root that it would compute anyway even
> + // if num==1. If this method is ever used from any other method, or
> + // if the solveCubic implementation changes, this assumption should
> + // be reevaluated, and the choice of goodRoot might have to become
> + // goodRoot = (abs(eqn'(res[0])) > abs(eqn'(res[1]))) ? res[0] : res[1]
> + // where eqn' is the derivative of eqn.
> + double goodRoot = res[0];
> + double badRoot = res[1];
> + double x1 = intervals[0];
> + double x2 = intervals[1];
> + // If a cubic curve really has 2 roots, one of those roots must be
> + // at a critical point. That can't be goodRoot, so we compute x to
> + // be the farthest critical point from goodRoot. If there are two
> + // roots, x must be the second one, so we evaluate eqn at x, and if
> + // it is zero (or close enough) we put x in res[1] (or badRoot, if
> + // |solveEqn(eqn, 3, badRoot)| < |solveEqn(eqn, 3, x)| but this
> + // shouldn't happen often).
> + double x = abs(x1 - goodRoot) > abs(x2 - goodRoot) ? x1 : x2;
> + double fx = solveEqn(eqn, 3, x);
> +
> + if (iszero(fx, 10000000*ulp(x))) {
> + double badRootVal = solveEqn(eqn, 3, badRoot);
> + res[1] = abs(badRootVal) < abs(fx) ? badRoot : x;
> + return 2;
> + }
> + } // else there can only be one root - goodRoot, and it is already in res[0]
> +
> + return 1;
> }
>
> - /*
> - * This pruning step is necessary since solveCubic uses the
> - * cosine function to calculate the roots when there are 3
> - * of them. Since the cosine method can have an error of
> - * +/- 1E-14 we need to make sure that we don't make any
> - * bad decisions due to an error.
> - *
> - * If the root is not near one of the endpoints, then we will
> - * only have a slight inaccuracy in calculating the x intercept
> - * which will only cause a slightly wrong answer for some
> - * points very close to the curve. While the results in that
> - * case are not as accurate as they could be, they are not
> - * disastrously inaccurate either.
> - *
> - * On the other hand, if the error happens near one end of
> - * the curve, then our processing to reject values outside
> - * of the t=[0,1] range will fail and the results of that
> - * failure will be disastrous since for an entire horizontal
> - * range of test points, we will either overcount or undercount
> - * the crossings and get a wrong answer for all of them, even
> - * when they are clearly and obviously inside or outside the
> - * curve.
> - *
> - * To work around this problem, we try a couple of Newton-Raphson
> - * iterations to see if the true root is closer to the endpoint
> - * or further away. If it is further away, then we can stop
> - * since we know we are on the right side of the endpoint. If
> - * we change direction, then either we are now being dragged away
> - * from the endpoint in which case the first condition will cause
> - * us to stop, or we have passed the endpoint and are headed back.
> - * In the second case, we simply evaluate the slope at the
> - * endpoint itself and place ourselves on the appropriate side
> - * of it or on it depending on that result.
> - */
> - private static void fixRoots(double res[], double eqn[]) {
> - final double EPSILON = 1E-5;
> + // use newton's method.
> + private static double refineRootWithHint(double[] eqn, double min, double max, double t) {
> + if (!inInterval(t, min, max)) {
> + return t;
> + }
> + double[] deriv = {eqn[1], 2*eqn[2], 3*eqn[3]};
> + double origt = t;
> for (int i = 0; i < 3; i++) {
> - double t = res[i];
> - if (Math.abs(t) < EPSILON) {
> - res[i] = findZero(t, 0, eqn);
> - } else if (Math.abs(t - 1) < EPSILON) {
> - res[i] = findZero(t, 1, eqn);
> + double slope = solveEqn(deriv, 2, t);
> + double y = solveEqn(eqn, 3, t);
> + double delta = - (y / slope);
> + double newt = t + delta;
> +
> + if (slope == 0 || y == 0 || t == newt) {
> + break;
> }
> +
> + t = newt;
> + }
> + if (within(t, origt, 1000*ulp(origt)) && inInterval(t, min, max)) {
> + return t;
> }
> + return origt;
> + }
> +
> + private static double bisectRootWithHint(double[] eqn, double x0, double xe, double hint) {
> + double delta1 = Math.min(abs(hint - x0) / 64, 0.0625);
> + double delta2 = Math.min(abs(hint - xe) / 64, 0.0625);
> + double x02 = hint - delta1;
> + double xe2 = hint + delta2;
> + double fx02 = solveEqn(eqn, 3, x02);
> + double fxe2 = solveEqn(eqn, 3, xe2);
> + while (oppositeSigns(fx02, fxe2)) {
> + if (x02 >= xe2) {
> + return x02;
> + }
> + x0 = x02;
> + xe = xe2;
> + delta1 /= 64;
> + delta2 /= 64;
> + x02 = hint - delta1;
> + xe2 = hint + delta2;
> + fx02 = solveEqn(eqn, 3, x02);
> + fxe2 = solveEqn(eqn, 3, xe2);
> + }
> + if (fx02 == 0) {
> + return x02;
> + }
> + if (fxe2 == 0) {
> + return xe2;
> + }
> +
> + return bisectRoot(eqn, x0, xe);
> + }
> +
> + private static double bisectRoot(double[] eqn, double x0, double xe) {
> + double fx0 = solveEqn(eqn, 3, x0);
> + double m = x0 + (xe - x0) / 2;
> + while (m != x0 && m != xe) {
> + double fm = solveEqn(eqn, 3, m);
> + if (fm == 0) {
> + return m;
> + }
> + if (oppositeSigns(fx0, fm)) {
> + xe = m;
> + } else {
> + fx0 = fm;
> + x0 = m;
> + }
> + m = x0 + (xe-x0)/2;
> + }
> + return m;
> + }
> +
> + private static boolean inInterval(double t, double min, double max) {
> + return min <= t && t <= max;
> + }
> +
> + private static boolean within(double x, double y, double err) {
> + double d = y - x;
> + return (d <= err && d >= -err);
> + }
> +
> + private static boolean iszero(double x, double err) {
> + return within(x, 0, err);
> + }
> +
> + private static boolean oppositeSigns(double x1, double x2) {
> + return (x1 < 0 && x2 > 0) || (x1 > 0 && x2 < 0);
> }
>
> private static double solveEqn(double eqn[], int order, double t) {
> @@ -1182,60 +1377,26 @@
> return v;
> }
>
> - private static double findZero(double t, double target, double eqn[]) {
> - double slopeqn[] = {eqn[1], 2*eqn[2], 3*eqn[3]};
> - double slope;
> - double origdelta = 0;
> - double origt = t;
> - while (true) {
> - slope = solveEqn(slopeqn, 2, t);
> - if (slope == 0) {
> - // At a local minima - must return
> - return t;
> - }
> - double y = solveEqn(eqn, 3, t);
> - if (y == 0) {
> - // Found it! - return it
> - return t;
> - }
> - // assert(slope != 0 && y != 0);
> - double delta = - (y / slope);
> - // assert(delta != 0);
> - if (origdelta == 0) {
> - origdelta = delta;
> - }
> - if (t < target) {
> - if (delta < 0) return t;
> - } else if (t > target) {
> - if (delta > 0) return t;
> - } else { /* t == target */
> - return (delta > 0
> - ? (target + java.lang.Double.MIN_VALUE)
> - : (target - java.lang.Double.MIN_VALUE));
> - }
> - double newt = t + delta;
> - if (t == newt) {
> - // The deltas are so small that we aren't moving...
> - return t;
> - }
> - if (delta * origdelta < 0) {
> - // We have reversed our path.
> - int tag = (origt < t
> - ? getTag(target, origt, t)
> - : getTag(target, t, origt));
> - if (tag != INSIDE) {
> - // Local minima found away from target - return the middle
> - return (origt + t) / 2;
> - }
> - // Local minima somewhere near target - move to target
> - // and let the slope determine the resulting t.
> - t = target;
> - } else {
> - t = newt;
> - }
> - }
> + /*
> + * Computes M+1 where M is an upper bound for all the roots in of eqn.
> + * See: http://en.wikipedia.org/wiki/Sturm%27s_theorem#Applications.
> + * The above link doesn't contain a proof, but I [dlila] proved it myself
> + * so the result is reliable. The proof isn't difficult, but it's a bit
> + * long to include here.
> + * Precondition: eqn must represent a cubic polynomial
> + */
> + private static double getRootUpperBound(double[] eqn) {
> + double d = eqn[3];
> + double a = eqn[2];
> + double b = eqn[1];
> + double c = eqn[0];
> +
> + double M = 1 + max(max(abs(a), abs(b)), abs(c)) / abs(d);
> + M += ulp(M) + 1;
> + return M;
> }
>
> +
> /**
> * {@inheritDoc}
> * @since 1.2
> @@ -1273,110 +1434,6 @@
> return contains(p.getX(), p.getY());
> }
>
> - /*
> - * Fill an array with the coefficients of the parametric equation
> - * in t, ready for solving against val with solveCubic.
> - * We currently have:
> - * <pre>
> - * val = P(t) = C1(1-t)^3 + 3CP1 t(1-t)^2 + 3CP2 t^2(1-t) + C2 t^3
> - * = C1 - 3C1t + 3C1t^2 - C1t^3 +
> - * 3CP1t - 6CP1t^2 + 3CP1t^3 +
> - * 3CP2t^2 - 3CP2t^3 +
> - * C2t^3
> - * 0 = (C1 - val) +
> - * (3CP1 - 3C1) t +
> - * (3C1 - 6CP1 + 3CP2) t^2 +
> - * (C2 - 3CP2 + 3CP1 - C1) t^3
> - * 0 = C + Bt + At^2 + Dt^3
> - * C = C1 - val
> - * B = 3*CP1 - 3*C1
> - * A = 3*CP2 - 6*CP1 + 3*C1
> - * D = C2 - 3*CP2 + 3*CP1 - C1
> - * </pre>
> - */
> - private static void fillEqn(double eqn[], double val,
> - double c1, double cp1, double cp2, double c2) {
> - eqn[0] = c1 - val;
> - eqn[1] = (cp1 - c1) * 3.0;
> - eqn[2] = (cp2 - cp1 - cp1 + c1) * 3.0;
> - eqn[3] = c2 + (cp1 - cp2) * 3.0 - c1;
> - return;
> - }
> -
> - /*
> - * Evaluate the t values in the first num slots of the vals[] array
> - * and place the evaluated values back into the same array. Only
> - * evaluate t values that are within the range <0, 1>, including
> - * the 0 and 1 ends of the range iff the include0 or include1
> - * booleans are true. If an "inflection" equation is handed in,
> - * then any points which represent a point of inflection for that
> - * cubic equation are also ignored.
> - */
> - private static int evalCubic(double vals[], int num,
> - boolean include0,
> - boolean include1,
> - double inflect[],
> - double c1, double cp1,
> - double cp2, double c2) {
> - int j = 0;
> - for (int i = 0; i < num; i++) {
> - double t = vals[i];
> - if ((include0 ? t >= 0 : t > 0) &&
> - (include1 ? t <= 1 : t < 1) &&
> - (inflect == null ||
> - inflect[1] + (2*inflect[2] + 3*inflect[3]*t)*t != 0))
> - {
> - double u = 1 - t;
> - vals[j++] = c1*u*u*u + 3*cp1*t*u*u + 3*cp2*t*t*u + c2*t*t*t;
> - }
> - }
> - return j;
> - }
> -
> - private static final int BELOW = -2;
> - private static final int LOWEDGE = -1;
> - private static final int INSIDE = 0;
> - private static final int HIGHEDGE = 1;
> - private static final int ABOVE = 2;
> -
> - /*
> - * Determine where coord lies with respect to the range from
> - * low to high. It is assumed that low <= high. The return
> - * value is one of the 5 values BELOW, LOWEDGE, INSIDE, HIGHEDGE,
> - * or ABOVE.
> - */
> - private static int getTag(double coord, double low, double high) {
> - if (coord <= low) {
> - return (coord < low ? BELOW : LOWEDGE);
> - }
> - if (coord >= high) {
> - return (coord > high ? ABOVE : HIGHEDGE);
> - }
> - return INSIDE;
> - }
> -
> - /*
> - * Determine if the pttag represents a coordinate that is already
> - * in its test range, or is on the border with either of the two
> - * opttags representing another coordinate that is "towards the
> - * inside" of that test range. In other words, are either of the
> - * two "opt" points "drawing the pt inward"?
> - */
> - private static boolean inwards(int pttag, int opt1tag, int opt2tag) {
> - switch (pttag) {
> - case BELOW:
> - case ABOVE:
> - default:
> - return false;
> - case LOWEDGE:
> - return (opt1tag >= INSIDE || opt2tag >= INSIDE);
> - case INSIDE:
> - return true;
> - case HIGHEDGE:
> - return (opt1tag <= INSIDE || opt2tag <= INSIDE);
> - }
> - }
> -
> /**
> * {@inheritDoc}
> * @since 1.2
> --- openjdk/jdk/src/share/classes/java/awt/geom/CubicCurve2D.java.old 2011-01-20 18:54:14.000000000 -0500
> +++ openjdk/jdk/src/share/classes/java/awt/geom/CubicCurve2D.java 2011-01-31 12:38:24.340733697 -0500
> @@ -1602,20 +1602,32 @@
> if (w <= 0 || h <= 0) {
> return false;
> }
> - // Assertion: Cubic curves closed by connecting their
> - // endpoints form either one or two convex halves with
> - // the closing line segment as an edge of both sides.
> - if (!(contains(x, y) &&
> - contains(x + w, y) &&
> - contains(x + w, y + h) &&
> - contains(x, y + h))) {
> - return false;
> +
> + int numCrossings = rectCrossings(x, y, w, h);
> + return !(numCrossings == 0 || numCrossings == Curve.RECT_INTERSECTS);
> + }
> +
> + private int rectCrossings(double x, double y, double w, double h) {
> + int crossings = 0;
> + if (!(getX1() == getX2() && getY1() == getY2())) {
> + crossings = Curve.rectCrossingsForLine(crossings,
> + x, y,
> + x+w, y+h,
> + getX1(), getY1(),
> + getX2(), getY2());
> + if (crossings == Curve.RECT_INTERSECTS) {
> + return crossings;
> + }
> }
> - // Either the rectangle is entirely inside one of the convex
> - // halves or it crosses from one to the other, in which case
> - // it must intersect the closing line segment.
> - Rectangle2D rect = new Rectangle2D.Double(x, y, w, h);
> - return !rect.intersectsLine(getX1(), getY1(), getX2(), getY2());
> + // we call this with the curve's direction reversed, because we wanted
> + // to call rectCrossingsForLine first, because it's cheaper.
> + return Curve.rectCrossingsForCubic(crossings,
> + x, y,
> + x+w, y+h,
> + getX2(), getY2(),
> + getCtrlX2(), getCtrlY2(),
> + getCtrlX1(), getCtrlY1(),
> + getX1(), getY1(), 0);
> }
>
> /**
--
Andrew :)
Free Java Software Engineer
Red Hat, Inc. (http://www.redhat.com)
Support Free Java!
Contribute to GNU Classpath and IcedTea
http://www.gnu.org/software/classpath
http://icedtea.classpath.org
PGP Key: 94EFD9D8 (http://subkeys.pgp.net)
Fingerprint = F8EF F1EA 401E 2E60 15FA 7927 142C 2591 94EF D9D8
More information about the distro-pkg-dev
mailing list