cubiccurve2d.java

gcc的组建

   * i.e. the length of the red line.
   */
  public double getFlatness()
  {
    return Math.sqrt(getFlatnessSq(getX1(), getY1(), getCtrlX1(), getCtrlY1(),
                                   getCtrlX2(), getCtrlY2(), getX2(), getY2()));
  }

  /**
   * Subdivides this curve into two halves.
   *
   * <p><img src="doc-files/CubicCurve2D-3.png" width="700"
   * height="180" alt="A drawing that illustrates the effects of
   * subdividing a CubicCurve2D" />
   *
   * @param left a curve whose geometry will be set to the left half
   * of this curve, or <code>null</code> if the caller is not
   * interested in the left half.
   *
   * @param right a curve whose geometry will be set to the right half
   * of this curve, or <code>null</code> if the caller is not
   * interested in the right half.
   */
  public void subdivide(CubicCurve2D left, CubicCurve2D right)
  {
    // Use empty slots at end to share single array.
    double[] d = new double[]
                 {
                   getX1(), getY1(), getCtrlX1(), getCtrlY1(), getCtrlX2(),
                   getCtrlY2(), getX2(), getY2(), 0, 0, 0, 0, 0, 0
                 };
    subdivide(d, 0, d, 0, d, 6);
    if (left != null)
      left.setCurve(d, 0);
    if (right != null)
      right.setCurve(d, 6);
  }

  /**
   * Subdivides a cubic curve into two halves.
   *
   * <p><img src="doc-files/CubicCurve2D-3.png" width="700"
   * height="180" alt="A drawing that illustrates the effects of
   * subdividing a CubicCurve2D" />
   *
   * @param src the curve to be subdivided.
   *
   * @param left a curve whose geometry will be set to the left half
   * of <code>src</code>, or <code>null</code> if the caller is not
   * interested in the left half.
   *
   * @param right a curve whose geometry will be set to the right half
   * of <code>src</code>, or <code>null</code> if the caller is not
   * interested in the right half.
   */
  public static void subdivide(CubicCurve2D src, CubicCurve2D left,
                               CubicCurve2D right)
  {
    src.subdivide(left, right);
  }

  /**
   * Subdivides a cubic curve into two halves, passing all coordinates
   * in an array.
   *
   * <p><img src="doc-files/CubicCurve2D-3.png" width="700"
   * height="180" alt="A drawing that illustrates the effects of
   * subdividing a CubicCurve2D" />
   *
   * <p>The left end point and the right start point will always be
   * identical. Memory-concious programmers thus may want to pass the
   * same array for both <code>left</code> and <code>right</code>, and
   * set <code>rightOff</code> to <code>leftOff + 6</code>.
   *
   * @param src an array containing the coordinates of the curve to be
   * subdivided.  The <i>x</i> coordinate of the start point P1 is
   * located at <code>src[srcOff]</code>, its <i>y</i> at
   * <code>src[srcOff + 1]</code>.  The <i>x</i> coordinate of the
   * first control point C1 is located at <code>src[srcOff +
   * 2]</code>, its <i>y</i> at <code>src[srcOff + 3]</code>.  The
   * <i>x</i> coordinate of the second control point C2 is located at
   * <code>src[srcOff + 4]</code>, its <i>y</i> at <code>src[srcOff +
   * 5]</code>. The <i>x</i> coordinate of the end point is located at
   * <code>src[srcOff + 6]</code>, its <i>y</i> at <code>src[srcOff +
   * 7]</code>.
   *
   * @param srcOff an offset into <code>src</code>, specifying
   * the index of the start point&#x2019;s <i>x</i> coordinate.
   *
   * @param left an array that will receive the coordinates of the
   * left half of <code>src</code>. It is acceptable to pass
   * <code>src</code>. A caller who is not interested in the left half
   * can pass <code>null</code>.
   *
   * @param leftOff an offset into <code>left</code>, specifying the
   * index where the start point&#x2019;s <i>x</i> coordinate will be
   * stored.
   *
   * @param right an array that will receive the coordinates of the
   * right half of <code>src</code>. It is acceptable to pass
   * <code>src</code> or <code>left</code>. A caller who is not
   * interested in the right half can pass <code>null</code>.
   *
   * @param rightOff an offset into <code>right</code>, specifying the
   * index where the start point&#x2019;s <i>x</i> coordinate will be
   * stored.
   */
  public static void subdivide(double[] src, int srcOff, double[] left,
                               int leftOff, double[] right, int rightOff)
  {
    // To understand this code, please have a look at the image
    // "CubicCurve2D-3.png" in the sub-directory "doc-files".
    double src_C1_x;
    double src_C1_y;
    double src_C2_x;
    double src_C2_y;
    double left_P1_x;
    double left_P1_y;
    double left_C1_x;
    double left_C1_y;
    double left_C2_x;
    double left_C2_y;
    double right_C1_x;
    double right_C1_y;
    double right_C2_x;
    double right_C2_y;
    double right_P2_x;
    double right_P2_y;
    double Mid_x; // Mid = left.P2 = right.P1
    double Mid_y; // Mid = left.P2 = right.P1

    left_P1_x = src[srcOff];
    left_P1_y = src[srcOff + 1];
    src_C1_x = src[srcOff + 2];
    src_C1_y = src[srcOff + 3];
    src_C2_x = src[srcOff + 4];
    src_C2_y = src[srcOff + 5];
    right_P2_x = src[srcOff + 6];
    right_P2_y = src[srcOff + 7];

    left_C1_x = (left_P1_x + src_C1_x) / 2;
    left_C1_y = (left_P1_y + src_C1_y) / 2;
    right_C2_x = (right_P2_x + src_C2_x) / 2;
    right_C2_y = (right_P2_y + src_C2_y) / 2;
    Mid_x = (src_C1_x + src_C2_x) / 2;
    Mid_y = (src_C1_y + src_C2_y) / 2;
    left_C2_x = (left_C1_x + Mid_x) / 2;
    left_C2_y = (left_C1_y + Mid_y) / 2;
    right_C1_x = (Mid_x + right_C2_x) / 2;
    right_C1_y = (Mid_y + right_C2_y) / 2;
    Mid_x = (left_C2_x + right_C1_x) / 2;
    Mid_y = (left_C2_y + right_C1_y) / 2;

    if (left != null)
      {
	left[leftOff] = left_P1_x;
	left[leftOff + 1] = left_P1_y;
	left[leftOff + 2] = left_C1_x;
	left[leftOff + 3] = left_C1_y;
	left[leftOff + 4] = left_C2_x;
	left[leftOff + 5] = left_C2_y;
	left[leftOff + 6] = Mid_x;
	left[leftOff + 7] = Mid_y;
      }

    if (right != null)
      {
	right[rightOff] = Mid_x;
	right[rightOff + 1] = Mid_y;
	right[rightOff + 2] = right_C1_x;
	right[rightOff + 3] = right_C1_y;
	right[rightOff + 4] = right_C2_x;
	right[rightOff + 5] = right_C2_y;
	right[rightOff + 6] = right_P2_x;
	right[rightOff + 7] = right_P2_y;
      }
  }

  /**
   * Finds the non-complex roots of a cubic equation, placing the
   * results into the same array as the equation coefficients. The
   * following equation is being solved:
   *
   * <blockquote><code>eqn[3]</code> &#xb7; <i>x</i><sup>3</sup>
   * + <code>eqn[2]</code> &#xb7; <i>x</i><sup>2</sup>
   * + <code>eqn[1]</code> &#xb7; <i>x</i>
   * + <code>eqn[0]</code>
   * = 0
   * </blockquote>
   *
   * <p>For some background about solving cubic equations, see the
   * article <a
   * href="http://planetmath.org/encyclopedia/CubicFormula.html"
   * >&#x201c;Cubic Formula&#x201d;</a> in <a
   * href="http://planetmath.org/" >PlanetMath</a>.  For an extensive
   * library of numerical algorithms written in the C programming
   * language, see the <a href= "http://www.gnu.org/software/gsl/">GNU
   * Scientific Library</a>, from which this implementation was
   * adapted.
   *
   * @param eqn an array with the coefficients of the equation. When
   * this procedure has returned, <code>eqn</code> will contain the
   * non-complex solutions of the equation, in no particular order.
   *
   * @return the number of non-complex solutions. A result of 0
   * indicates that the equation has no non-complex solutions. A
   * result of -1 indicates that the equation is constant (i.e.,
   * always or never zero).
   *
   * @see #solveCubic(double[], double[])
   * @see QuadCurve2D#solveQuadratic(double[],double[])
   *
   * @author Brian Gough ([e-mail])
   * (original C implementation in the <a href=
   * "http://www.gnu.org/software/gsl/">GNU Scientific Library</a>)
   *
   * @author Sascha Brawer ([e-mail])
   * (adaptation to Java)
   */
  public static int solveCubic(double[] eqn)
  {
    return solveCubic(eqn, eqn);
  }

  /**
   * Finds the non-complex roots of a cubic equation. The following
   * equation is being solved:
   *
   * <blockquote><code>eqn[3]</code> &#xb7; <i>x</i><sup>3</sup>
   * + <code>eqn[2]</code> &#xb7; <i>x</i><sup>2</sup>
   * + <code>eqn[1]</code> &#xb7; <i>x</i>
   * + <code>eqn[0]</code>
   * = 0
   * </blockquote>
   *
   * <p>For some background about solving cubic equations, see the
   * article <a
   * href="http://planetmath.org/encyclopedia/CubicFormula.html"
   * >&#x201c;Cubic Formula&#x201d;</a> in <a
   * href="http://planetmath.org/" >PlanetMath</a>.  For an extensive
   * library of numerical algorithms written in the C programming
   * language, see the <a href= "http://www.gnu.org/software/gsl/">GNU
   * Scientific Library</a>, from which this implementation was
   * adapted.
   *
   * @see QuadCurve2D#solveQuadratic(double[],double[])
   *
   * @param eqn an array with the coefficients of the equation.
   *
   * @param res an array into which the non-complex roots will be
   * stored.  The results may be in an arbitrary order. It is safe to
   * pass the same array object reference for both <code>eqn</code>
   * and <code>res</code>.
   *
   * @return the number of non-complex solutions. A result of 0
   * indicates that the equation has no non-complex solutions. A
   * result of -1 indicates that the equation is constant (i.e.,
   * always or never zero).
   *
   * @author Brian Gough ([e-mail])
   * (original C implementation in the <a href=
   * "http://www.gnu.org/software/gsl/">GNU Scientific Library</a>)
   *
   * @author Sascha Brawer ([e-mail])
   * (adaptation to Java)
   */
  public static int solveCubic(double[] eqn, double[] res)
  {
    // Adapted from poly/solve_cubic.c in the GNU Scientific Library
    // (GSL), revision 1.7 of 2003-07-26. For the original source, see
    // http://www.gnu.org/software/gsl/
    //
    // Brian Gough, the author of that code, has granted the
    // permission to use it in GNU Classpath under the GNU Classpath
    // license, and has assigned the copyright to the Free Software
    // Foundation.
    //
    // The Java implementation is very similar to the GSL code, but
    // not a strict one-to-one copy. For example, GSL would sort the
    // result.
    
    double a;
    double b;
    double c;
    double q;
    double r;
    double Q;
    double R;
    double c3;
    double Q3;
    double R2;
    double CR2;
    double CQ3;

    // If the cubic coefficient is zero, we have a quadratic equation.
    c3 = eqn[3];
    if (c3 == 0)
      return QuadCurve2D.solveQuadratic(eqn, res);

    // Divide the equation by the cubic coefficient.
    c = eqn[0] / c3;
    b = eqn[1] / c3;
    a = eqn[2] / c3;

    // We now need to solve x^3 + ax^2 + bx + c = 0.
    q = a * a - 3 * b;
    r = 2 * a * a * a - 9 * a * b + 27 * c;

    Q = q / 9;
    R = r / 54;

    Q3 = Q * Q * Q;
    R2 = R * R;

    CR2 = 729 * r * r;
    CQ3 = 2916 * q * q * q;

    if (R == 0 && Q == 0)
      {
	// The GNU Scientific Library would return three identical
	// solutions in this case.
	res[0] = -a / 3;
	return 1;
      }

    if (CR2 == CQ3)
      {
	/* this test is actually R2 == Q3, written in a form suitable
	   for exact computation with integers */
	/* Due to finite precision some double roots may be missed, and
	   considered to be a pair of complex roots z = x +/- epsilon i
	   close to the real axis. */
	double sqrtQ = Math.sqrt(Q);

	if (R > 0)
	  {
	    res[0] = -2 * sqrtQ - a / 3;
	    res[1] = sqrtQ - a / 3;
	  }
	else
	  {
	    res[0] = -sqrtQ - a / 3;
	    res[1] = 2 * sqrtQ - a / 3;
	  }
	return 2;
      }

    if (CR2 < CQ3) /* equivalent to R2 < Q3 */
      {
	double sqrtQ = Math.sqrt(Q);
	double sqrtQ3 = sqrtQ * sqrtQ * sqrtQ;
	double theta = Math.acos(R / sqrtQ3);
	double norm = -2 * sqrtQ;
	res[0] = norm * Math.cos(theta / 3) - a / 3;
	res[1] = norm * Math.cos((theta + 2.0 * Math.PI) / 3) - a / 3;
	res[2] = norm * Math.cos((theta - 2.0 * Math.PI) / 3) - a / 3;

	// The GNU Scientific Library sorts the results. We don't.
	return 3;
      }

    double sgnR = (R >= 0 ? 1 : -1);
    double A = -sgnR * Math.pow(Math.abs(R) + Math.sqrt(R2 - Q3), 1.0 / 3.0);
    double B = Q / A;
    res[0] = A + B - a / 3;
    return 1;
  }

  /**
   * Determines whether a position lies inside the area bounded
   * by the curve and the straight line connecting its end points.
   *
   * <p><img src="doc-files/CubicCurve2D-5.png" width="350" height="180"
   * alt="A drawing of the area spanned by the curve" />
   *
   * <p>The above drawing illustrates in which area points are
   * considered &#x201c;inside&#x201d; a CubicCurve2D.
   */
  public boolean contains(double x, double y)
  {
    if (! getBounds2D().contains(x, y))
      return false;

    return ((getAxisIntersections(x, y, true, BIG_VALUE) & 1) != 0);
  }

  /**
   * Determines whether a point lies inside the area bounded
   * by the curve and the straight line connecting its end points.
   *
   * <p><img src="doc-files/CubicCurve2D-5.png" width="350" height="180"
   * alt="A drawing of the area spanned by the curve" />
   *
   * <p>The above drawing illustrates in which area points are
   * considered &#x201c;inside&#x201d; a CubicCurve2D.
   */
  public boolean contains(Point2D p)
  {
    return contains(p.getX(), p.getY());
  }

  /**
   * Determines whether any part of a rectangle is inside the area bounded
   * by the curve and the straight line connecting its end points.
   *
   * <p><img src="doc-files/CubicCurve2D-5.png" width="350" height="180"
   * alt="A drawing of the area spanned by the curve" />
   *
   * <p>The above drawing illustrates in which area points are
   * considered &#x201c;inside&#x201d; in a CubicCurve2D.
   * @see #contains(double, double)
   */
  public boolean intersects(double x, double y, double w, double h)
  {
    if (! getBounds2D().contains(x, y, w, h))
      return false;

    /* Does any edge intersect? */
    if (getAxisIntersections(x, y, true, w) != 0 /* top */
        || getAxisIntersections(x, y + h, true, w) != 0 /* bottom */
        || getAxisIntersections(x + w, y, false, h) != 0 /* right */
        || getAxisIntersections(x, y, false, h) != 0) /* left */
      return true;

    /* No intersections, is any point inside? */
    if ((getAxisIntersections(x, y, true, BIG_VALUE) & 1) != 0)
      return true;

    return false;
  }

  /**