conic_arc_2.h

很多二维 三维几何计算算法 C++ 类库

// Copyright (c) 2005  Tel-Aviv University (Israel).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org); you may redistribute it under
// the terms of the Q Public License version 1.0.
// See the file LICENSE.QPL distributed with CGAL.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $URL: svn+ssh://scm.gforge.inria.fr/svn/cgal/branches/CGAL-3.3-branch/Arrangement_2/include/CGAL/Arr_traits_2/Conic_arc_2.h $
// $Id: Conic_arc_2.h 36449 2007-02-19 09:39:17Z ophirset $
// 
//
// Author(s)     : Ron Wein <[e-mail]>

#ifndef CGAL_CONIC_ARC_2_H
#define CGAL_CONIC_ARC_2_H

/*! \file
 * Header file for the _Conic_arc_2<Int_kernel, Alg_kernel, Nt_traits> class.
 */

#include <CGAL/Arr_traits_2/Conic_point_2.h>
#include <CGAL/Arr_traits_2/Conic_intersections_2.h>
#include <CGAL/Bbox_2.h>

#include <ostream>

CGAL_BEGIN_NAMESPACE


/*!
 * Representation of a conic arc -- a bounded segment that lies on a conic
 * curve, the loci of all points satisfying the equation:
 *   r*x^2 + s*y^2 + t*xy + u*x + v*y +w = 0
 *
 * The class is templated with three parameters: 
 * Rat_kernel A kernel that provides the input objects or coefficients.
 *            Rat_kernel::FT should be an integral or a rational type.
 * Alg_kernel A geometric kernel, where Alg_kernel::FT is the number type
 *            for the coordinates of arrangement vertices, which are algebraic
 *            numbers of degree up to 4 (preferably it is CORE::Expr).
 * Nt_traits A traits class for performing various operations on the integer,
 *           rational and algebraic types. 
 */

template <class Rat_kernel_, class Alg_kernel_, class Nt_traits_>
class _Conic_arc_2
{
public:

  typedef Rat_kernel_                                      Rat_kernel;
  typedef Alg_kernel_                                      Alg_kernel;
  typedef Nt_traits_                                       Nt_traits;

  typedef _Conic_arc_2<Rat_kernel, Alg_kernel, Nt_traits>  Self;
 
  typedef typename Rat_kernel::FT                          Rational;
  typedef typename Rat_kernel::Point_2                     Rat_point_2;
  typedef typename Rat_kernel::Segment_2                   Rat_segment_2;
  typedef typename Rat_kernel::Circle_2                    Rat_circle_2;

  typedef typename Nt_traits::Integer                      Integer;

  typedef typename Alg_kernel::FT                          Algebraic;
  typedef typename Alg_kernel::Point_2                     Point_2;
  typedef _Conic_point_2<Alg_kernel>                       Conic_point_2;

protected:

  Integer        _r;       //
  Integer        _s;       // The coefficients of the supporting conic curve:
  Integer        _t;       //
  Integer        _u;       //
  Integer        _v;       //   r*x^2 + s*y^2 + t*xy + u*x + v*y +w = 0 .
  Integer        _w;       //

  Orientation    _orient;  // The orientation of the conic.

  // Bit masks for the _info field.
  enum
  {
    IS_VALID = 1,
    IS_FULL_CONIC = 2
  };

  int            _info;    // Does the arc represent a full conic curve.
  Conic_point_2  _source;  // The source of the arc (if not a full curve).
  Conic_point_2  _target;  // The target of the arc (if not a full curve).

  /*! \struct
   * For arcs whose base is a hyperbola we store the axis (a*x + b*y + c = 0)
   * which separates the two bracnes of the hyperbola. We also store the side
   * (NEGATIVE or POSITIVE) that the arc occupies.
   * In case of line segments connecting two algebraic endpoints, we use this
   * structure two store the coefficients of the line supporting this segment.
   * In this case we set the side field to be ZERO.
   */   
  struct Extra_data
  {
    Algebraic     a;
    Algebraic     b;
    Algebraic     c;
    Sign          side;
  };

  Extra_data    *_extra_data_P;  // The extra data stored with the arc
                                 // (may be NULL).

public:

  /// \name Construction and destruction fucntions.
  //@{

  /*!
   * Default constructor.
   */
  _Conic_arc_2 () :
    _r(0), _s(0), _t(0), _u(0), _v(0), _w(0),
    _orient (COLLINEAR),
    _info (0),
    _extra_data_P (NULL)
  {}

  /*!
   * Copy constructor.
   * \param arc The copied arc.
   */
  _Conic_arc_2 (const Self& arc) :
    _r(arc._r), _s(arc._s), _t(arc._t), _u(arc._u), _v(arc._v), _w(arc._w),
    _orient (arc._orient),
    _info (arc._info),
    _source(arc._source),
    _target(arc._target)
  {
    if (arc._extra_data_P != NULL)
      _extra_data_P = new Extra_data (*(arc._extra_data_P));
    else
      _extra_data_P = NULL;
  }

  /*! 
   * Construct a conic arc which is the full conic:
   *   C: r*x^2 + s*y^2 + t*xy + u*x + v*y + w = 0
   * \pre The conic C must be an ellipse (so 4rs - t^2 > 0).
   */
  _Conic_arc_2 (const Rational& r, const Rational& s, const Rational& t,
                const Rational& u, const Rational& v, const Rational& w) :
    _extra_data_P (NULL)
  {
    // Make sure the given curve is an ellipse (4rs - t^2 should be positive).
    CGAL_precondition (CGAL::sign (4*r*s - t*t) == POSITIVE);

    // Set the arc to be the full conic (and compute the orientation).
    Rational    rat_coeffs [6];

    rat_coeffs[0] = r;
    rat_coeffs[1] = s;
    rat_coeffs[2] = t;
    rat_coeffs[3] = u;
    rat_coeffs[4] = v;
    rat_coeffs[5] = w;

    _set_full (rat_coeffs, true);
  }

  /*! 
   * Construct a conic arc which lies on the conic:
   *   C: r*x^2 + s*y^2 + t*xy + u*x + v*y + w = 0
   * \param orient The orientation of the arc (clockwise or counterclockwise).
   * \param source The source point.
   * \param target The target point.
   * \pre The source and the target must be on the conic boundary and must
   * not be the same.
   */
  _Conic_arc_2 (const Rational& r, const Rational& s, const Rational& t,
                const Rational& u, const Rational& v, const Rational& w,
                const Orientation& orient,
                const Point_2& source, const Point_2& target) :
    _orient (orient),
    _source (source),
    _target (target),
    _extra_data_P (NULL)
  {
    // Make sure that the source and the taget are not the same.
    CGAL_precondition (Alg_kernel().compare_xy_2_object() (source,
                                                           target) != EQUAL);

    // Set the arc properties (no need to compute the orientation).
    Rational    rat_coeffs [6];

    rat_coeffs[0] = r;
    rat_coeffs[1] = s;
    rat_coeffs[2] = t;
    rat_coeffs[3] = u;
    rat_coeffs[4] = v;
    rat_coeffs[5] = w;

    _set (rat_coeffs);
  }

  /*!
   * Construct a conic arc from the given line segment.
   * \param seg The line segment with rational endpoints.
   */
  _Conic_arc_2 (const Rat_segment_2& seg) :
    _orient (COLLINEAR),
    _extra_data_P (NULL)
  {
    // Set the source and target.
    Rat_kernel        ker;
    Rat_point_2       source = ker.construct_vertex_2_object() (seg, 0); 
    Rat_point_2       target = ker.construct_vertex_2_object() (seg, 1); 
    Rational          x1 = source.x();
    Rational          y1 = source.y();
    Rational          x2 = target.x();
    Rational          y2 = target.y();
    Nt_traits         nt_traits;

    _source = Point_2 (nt_traits.convert (x1), nt_traits.convert (y1));
    _target = Point_2 (nt_traits.convert (x2), nt_traits.convert (y2));
 
    // Make sure that the source and the taget are not the same.
    CGAL_precondition (Alg_kernel().compare_xy_2_object() (_source,
                                                           _target) != EQUAL);

    // The supporting conic is r=s=t=0, and u*x + v*y + w = 0 should hold
    // for both the source (x1,y1) and the target (x2, y2).
    const Rational    _zero (0);
    const Rational    _one (1);
    Rational          rat_coeffs [6];

    rat_coeffs[0] = _zero;
    rat_coeffs[1] = _zero;
    rat_coeffs[2] = _zero;

    if (CGAL::compare (x1, x2) == EQUAL)
    {
      // The supporting conic is a vertical line, of the form x = CONST.
      rat_coeffs[3] = _one;
      rat_coeffs[4] = _zero;
      rat_coeffs[5] = -x1;
    }
    else
    {
      // The supporting line is A*x + B*y + C = 0, where:
      //
      //  A = y2 - y1,    B = x1 - x2,    C = x2*y1 - x1*y2 
      //
      rat_coeffs[3] = y2 - y1;
      rat_coeffs[4] = x1 - x2;
      rat_coeffs[5] = x2*y1 - x1*y2;
    }

    // Set the arc properties (no need to compute the orientation).
    _set (rat_coeffs);
  }

  /*!
   * Set a circular arc that corresponds to a full circle.
   * \param circ The circle (with rational center and rational squared radius).
   */ 
  _Conic_arc_2 (const Rat_circle_2& circ) :
    _orient (CLOCKWISE),
    _extra_data_P (NULL)
  {
    // Get the circle properties.
    Rat_kernel        ker;
    Rat_point_2       center = ker.construct_center_2_object() (circ); 
    Rational          x0 = center.x();
    Rational          y0 = center.y();
    Rational          R_sqr = ker.compute_squared_radius_2_object() (circ); 

    // Produce the correponding conic: if the circle center is (x0,y0)
    // and its squared radius is R^2, that its equation is:
    //   x^2 + y^2 - 2*x0*x - 2*y0*y + (x0^2 + y0^2 - R^2) = 0
    // Note that this equation describes a curve with a negative (clockwise) 
    // orientation.
    const Rational    _zero (0);
    const Rational    _one (1);
    const Rational    _minus_two (-2);
    Rational          rat_coeffs [6];

    rat_coeffs[0] = _one;
    rat_coeffs[1] = _one;
    rat_coeffs[2] = _zero;
    rat_coeffs[3] = _minus_two*x0;
    rat_coeffs[4] = _minus_two*y0;
    rat_coeffs[5] = x0*x0 + y0*y0 - R_sqr;

    // Set the arc to be the full conic (no need to compute the orientation).
    _set_full (rat_coeffs, false);
  }

  /*!
   * Set a circular arc that lies on the given circle:
   *   C: (x - x0)^2 + (y - y0)^2 = R^2
   * \param orient The orientation of the circle.
   * \param source The source point.
   * \param target The target point.
   * \pre The source and the target must be on the conic boundary and must
   *      not be the same.
   */  
  _Conic_arc_2 (const Rat_circle_2& circ,
                const Orientation& orient,
                const Point_2& source, const Point_2& target) :
    _orient(orient),
    _source(source),
    _target(target),
    _extra_data_P (NULL)
  {
    // Make sure that the source and the taget are not the same.
    CGAL_precondition (Alg_kernel().compare_xy_2_object() (source,
                                                           target) != EQUAL);
    CGAL_precondition (orient != COLLINEAR);

    // Get the circle properties.
    Rat_kernel        ker;
    Rat_point_2       center = ker.construct_center_2_object() (circ);
    Rational          x0 = center.x();
    Rational          y0 = center.y();
    Rational          R_sqr = ker.compute_squared_radius_2_object() (circ); 

    // Produce the correponding conic: if the circle center is (x0,y0)
    // and it squared radius is R^2, that its equation is:
    //   x^2 + y^2 - 2*x0*x - 2*y0*y + (x0^2 + y0^2 - R^2) = 0
    // Since this equation describes a curve with a negative (clockwise) 
    // orientation, we multiply it by -1 if necessary to obtain a positive
    // (counterclockwise) orientation.
    const Rational    _zero (0);
    Rational          rat_coeffs[6];

    if (_orient == COUNTERCLOCKWISE)
    {
      const Rational  _minus_one (-1);
      const Rational  _two (2);

      rat_coeffs[0] = _minus_one;
      rat_coeffs[1] = _minus_one;
      rat_coeffs[2] = _zero;
      rat_coeffs[3] = _two*x0;
      rat_coeffs[4] = _two*y0;
      rat_coeffs[5] = R_sqr - x0*x0 - y0*y0;
    }
    else
    {
      const Rational    _one (1);
      const Rational    _minus_two (-2);

      rat_coeffs[0] = _one;
      rat_coeffs[1] = _one;
      rat_coeffs[2] = _zero;
      rat_coeffs[3] = _minus_two*x0;
      rat_coeffs[4] = _minus_two*y0;
      rat_coeffs[5] = x0*x0 + y0*y0 - R_sqr;
    }

    // Set the arc properties (no need to compute the orientation).
    _set (rat_coeffs);
  }

  /*!
   * Construct a circular arc from the given three points.
   * \param p1 The arc source.
   * \param p2 A point in the interior of the arc.
   * \param p3 The arc target.
   * \pre The three points must not be collinear.
   */
  _Conic_arc_2 (const Rat_point_2& p1,
                const Rat_point_2& p2,
                const Rat_point_2& p3):
    _extra_data_P (NULL)
  {
    // Set the source and target.
    Rational          x1 = p1.x();
    Rational          y1 = p1.y();
    Rational          x2 = p2.x();
    Rational          y2 = p2.y();
    Rational          x3 = p3.x();
    Rational          y3 = p3.y();
    Nt_traits         nt_traits;

    _source = Point_2 (nt_traits.convert (x1), nt_traits.convert (y1));
    _target = Point_2 (nt_traits.convert (x3), nt_traits.convert (y3));
 
    // Make sure that the source and the taget are not the same.
    CGAL_precondition (Alg_kernel().compare_xy_2_object() (_source,
                                                           _target) != EQUAL);

    // Compute the lines: A1*x + B1*y + C1 = 0,
    //               and: A2*x + B2*y + C2 = 0,
    // where:
    const Rational  _two  = 2;

    const Rational  A1 = _two*(x1 - x2);
    const Rational  B1 = _two*(y1 - y2);
    const Rational  C1 = y2*y2 - y1*y1 + x2*x2 - x1*x1;

    const Rational  A2 = _two*(x2 - x3);
    const Rational  B2 = _two*(y2 - y3);
    const Rational  C2 = y3*y3 - y2*y2 + x3*x3 - x2*x2;

    // Compute the coordinates of the intersection point between the
    // two lines, given by (Nx / D, Ny / D), where:
    const Rational  Nx = B1*C2 - B2*C1;
    const Rational  Ny = A2*C1 - A1*C2;
    const Rational  D = A1*B2 - A2*B1;

    // Make sure the three points are not collinear.
    const bool  points_collinear = (CGAL::sign (D) == ZERO);

    if (points_collinear)
    {
      _info = 0;           // Inavlid arc.
      return;
    }

    // The equation of the underlying circle is given by:
    Rational          rat_coeffs[6];
    
    rat_coeffs[0] = D*D;
    rat_coeffs[1] = D*D;
    rat_coeffs[2] = 0;
    rat_coeffs[3] = -_two*D*Nx;
    rat_coeffs[4] = -_two*D*Ny;
    rat_coeffs[5] = 
      Nx*Nx + Ny*Ny - ((D*x2 - Nx)*(D*x2 - Nx) + (D*y2 - Ny)*(D*y2 - Ny));

    // Determine the orientation: If the mid-point forms a left-turn with
    // the source and the target points, the orientation is positive (going
    // counterclockwise).
    // Otherwise, it is negative (going clockwise).
    Alg_kernel                         ker;
    typename Alg_kernel::Orientation_2 orient_f = ker.orientation_2_object();


    Point_2  p_mid = Point_2 (nt_traits.convert (x2), nt_traits.convert (y2));
 
    if (orient_f(_source, p_mid, _target) == LEFT_TURN)
      _orient = COUNTERCLOCKWISE;
    else
      _orient = CLOCKWISE;

    // Set the arc properties (no need to compute the orientation).
    _set (rat_coeffs);
  }

  /*!
   * Construct a conic arc from the given five points, specified by the
   * points p1, p2, p3, p4 and p5.
   * \param p1 The source point of the given arc.
   * \param p2,p3,p4 Points lying on the conic arc, between p1 and p5.
   * \param p5 The target point of the given arc.
   * \pre No three points are collinear.
   */
  _Conic_arc_2 (const Rat_point_2& p1,
		const Rat_point_2& p2,
		const Rat_point_2& p3,
		const Rat_point_2& p4,
		const Rat_point_2& p5) :
    _extra_data_P(NULL)
  {
    // Make sure that no three points are collinear.
    Rat_kernel                         ker;