function_objects.h

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

	{
	  pV = phy*qhw;
	  qV = qhy*phw;
	}
      return CGAL_NTS compare(pV, qV);
    }
  };

  template <typename K>
  class Compare_yx_2
  {
    typedef typename K::Point_2            Point_2;
  public:
    typedef typename K::Comparison_result  result_type;
    typedef Arity_tag< 2 >                 Arity;

    result_type
    operator()( const Point_2& p, const Point_2& q) const
    {
      typedef typename K::RT RT;

      const RT& phx = p.hx();
      const RT& phy = p.hy();
      const RT& phw = p.hw();
      const RT& qhx = q.hx();
      const RT& qhy = q.hy();
      const RT& qhw = q.hw();

      RT pV = phy*qhw;
      RT qV = qhy*phw;
      if ( pV == qV )
	{
	  pV = phx*qhw;
	  qV = qhx*phw;
	}
      return CGAL_NTS compare(pV, qV);
    }
  };

  template <typename K>
  class Compare_xy_3
  {
    typedef typename K::Point_3            Point_3;
  public:
    typedef typename K::Comparison_result  result_type;
    typedef Arity_tag< 2 >                 Arity;

    result_type
    operator()( const Point_3& p, const Point_3& q) const
    {
      typedef typename K::RT RT;
      RT pV = p.hx()*q.hw();
      RT qV = q.hx()*p.hw();
      if ( pV < qV )
	{
	  return SMALLER;
	}
      if ( qV < pV )    //   ( pV > qV )
	{
	  return LARGER;
	}
      // same x
      pV = p.hy()*q.hw();
      qV = q.hy()*p.hw();
      return CGAL_NTS compare(pV, qV);
    }
  };

  template <typename K>
  class Compare_x_2
  {
    typedef typename K::Point_2            Point_2;
    typedef typename K::Line_2             Line_2;
  public:
    typedef typename K::Comparison_result  result_type;
    typedef Arity_tag< 2 >                 Arity;

    result_type
    operator()( const Point_2& p, const Point_2& q) const
    {
      return CGAL_NTS compare(p.hx()*q.hw(), q.hx()*p.hw());
    }

    result_type
    operator()( const Point_2& p, const Line_2& l1, const Line_2& l2) const
    {
      Point_2 ip = gp_linear_intersection( l1, l2 );
      return this->operator()(p, ip);
    } // FIXME

    result_type
    operator()( const Line_2& l, const Line_2& h1, const Line_2& h2) const
    {
      return this->operator()(l, h1, l, h2);
    } // FIXME

    result_type
    operator()( const Line_2& l1, const Line_2& l2,
	        const Line_2& h1, const Line_2& h2) const
    {
      Point_2 lip = gp_linear_intersection( l1, l2 );
      Point_2 hip = gp_linear_intersection( h1, h2 );
      return this->operator()(lip, hip);
    } // FIXME
  };

  template <typename K>
  class Compare_x_3
  {
    typedef typename K::Point_3            Point_3;
  public:
    typedef typename K::Comparison_result  result_type;
    typedef Arity_tag< 2 >                 Arity;

    result_type
    operator()( const Point_3& p, const Point_3& q) const
    { return CGAL_NTS compare(p.hx() * q.hw(), q.hx() * p.hw() ); }
  };

  template <typename K>
  class Compare_y_at_x_2
  {
    typedef typename K::Point_2            Point_2;
    typedef typename K::Line_2             Line_2;
    typedef typename K::Segment_2          Segment_2;
  public:
    typedef typename K::Comparison_result  result_type;
    typedef Arity_tag< 3 >                 Arity;

    result_type
    operator()( const Point_2& p, const Line_2& h) const
    {
      typedef typename K::RT RT;
      CGAL_kernel_precondition( ! h.is_vertical() );
      typename K::Oriented_side ors = h.oriented_side( p );
      if ( h.b() < 0 )
	  ors = opposite( ors );
      return enum_cast<Comparison_result>(ors);
    } // FIXME

    result_type
    operator()( const Point_2& p, const Line_2& h1, const Line_2& h2) const
    { return CGAL_NTS compare(h1.y_at_x( p.x() ), h2.y_at_x( p.x() )); }
    // FIXME

    result_type
    operator()( const Line_2& l1, const Line_2& l2, const Line_2& h) const
    { return compare_y_at_x( gp_linear_intersection( l1, l2 ), h); }
    // FIXME

    result_type
    operator()( const Line_2& l1, const Line_2& l2,
	        const Line_2& h1, const Line_2& h2) const
    { return compare_y_at_x( gp_linear_intersection( l1, l2 ), h1, h2 ); }
    // FIXME

    result_type
    operator()( const Point_2& p, const Segment_2& s) const
    {
      // compares the y-coordinates of p and the vertical projection of p on s.
      // Precondition : p is in the x-range of s.

      if (compare_x(s.source(), s.target()) == SMALLER) {
        CGAL_kernel_precondition(compare_x(s.source(), p) != LARGER
				 && compare_x(p, s.target()) != LARGER);
        return (Comparison_result) orientation(p, s.source(), s.target());
      }
      else if (compare_x(s.source(), s.target()) == LARGER) {
        CGAL_kernel_precondition(compare_x(s.target(), p) != LARGER
				 && compare_x(p, s.source()) != LARGER);
        return (Comparison_result) orientation(p, s.target(), s.source());
      }
      else {
        CGAL_kernel_precondition(compare_x(s.target(), p) == EQUAL);
	if (compare_y(p, s.source()) == SMALLER &&
	    compare_y(p, s.target()) == SMALLER)
	  return SMALLER;
	if (compare_y(p, s.source()) == LARGER &&
	    compare_y(p, s.target()) == LARGER)
	  return LARGER;
	return EQUAL;
      }
    } // FIXME

    result_type
    operator()( const Point_2& p,
	        const Segment_2& s1, const Segment_2& s2) const
    {
      // compares the y-coordinates of the vertical projections
      //   of p on s1 and s2
      // Precondition : p is in the x-range of s1 and s2.
      // - if one or two segments are vertical :
      //   - if the segments intersect, return EQUAL
      //   - if not, return the obvious SMALLER/LARGER.

      typedef typename K::FT FT;
      FT px = p.x();
      FT s1sx = s1.source().x();
      FT s1sy = s1.source().y();
      FT s1tx = s1.target().x();
      FT s1ty = s1.target().y();
      FT s2sx = s2.source().x();
      FT s2sy = s2.source().y();
      FT s2tx = s2.target().x();
      FT s2ty = s2.target().y();

      CGAL_kernel_precondition(px >= (CGAL::min)(s1sx, s1tx) &&
	                       px <= (CGAL::max)(s1sx, s1tx));
      CGAL_kernel_precondition(px >= (CGAL::min)(s2sx, s2tx) &&
	                       px <= (CGAL::max)(s2sx, s2tx));

      if (s1sx != s1tx && s2sx != s2tx) {
	FT s1stx = s1sx-s1tx;
	FT s2stx = s2sx-s2tx;

	return enum_cast<Comparison_result>(
	 CGAL_NTS compare(s1sx, s1tx) *
	 CGAL_NTS compare(s2sx, s2tx) *
	 CGAL_NTS compare(-(s1sx-px)*(s1sy-s1ty)*s2stx,
			  (s2sy-s1sy)*s2stx*s1stx
			  -(s2sx-px)*(s2sy-s2ty)*s1stx ));
      }
      else {
	if (s1sx == s1tx) { // s1 is vertical
	  typename K::Comparison_result c1, c2;
	  c1 = compare_y_at_x(s1.source(), s2);
	  c2 = compare_y_at_x(s1.target(), s2);
	  if (c1 == c2)
	    return c1;
	  return EQUAL;
	}
	// s2 is vertical
	typename K::Comparison_result c3, c4;
	c3 = compare_y_at_x(s2.source(), s1);
	c4 = compare_y_at_x(s2.target(), s1);
	if (c3 == c4)
	  return opposite(c3);
	return EQUAL;
      }
    } // FIXME
  };

  template <typename K>
  class Compare_y_2
  {
    typedef typename K::Point_2            Point_2;
    typedef typename K::Line_2             Line_2;
  public:
    typedef typename K::Comparison_result  result_type;
    typedef Arity_tag< 2 >                 Arity;

    result_type
    operator()( const Point_2& p, const Point_2& q) const
    {
      typedef typename K::RT RT;

      const RT& phy = p.hy();
      const RT& phw = p.hw();
      const RT& qhy = q.hy();
      const RT& qhw = q.hw();
      return CGAL_NTS compare(phy * qhw, qhy * phw);
    }

    result_type
    operator()( const Point_2& p, const Line_2& l1, const Line_2& l2) const
    {
      Point_2 ip = gp_linear_intersection( l1, l2 );
      return compare_y( p, ip );
    } // FIXME

    result_type
    operator()( const Line_2& l, const Line_2& h1, const Line_2& h2) const
    {
      return this->operator()(l, h1, l, h2);
    }

    result_type
    operator()( const Line_2& l1, const Line_2& l2,
	        const Line_2& h1, const Line_2& h2) const
    {
      Point_2 lip = gp_linear_intersection( l1, l2 );
      Point_2 hip = gp_linear_intersection( h1, h2 );
      return this->operator()( lip, hip );
    } // FIXME
  };

  template <typename K>
  class Compare_y_3
  {
    typedef typename K::Point_3            Point_3;
  public:
    typedef typename K::Comparison_result  result_type;
    typedef Arity_tag< 2 >                 Arity;

    result_type
    operator()( const Point_3& p, const Point_3& q) const
    { return CGAL_NTS compare(p.hy() * q.hw(), q.hy() * p.hw() ); }
  };

  template <typename K>
  class Compare_z_3
  {
    typedef typename K::Point_3             Point_3;
  public:
    typedef typename K::Comparison_result   result_type;
    typedef Arity_tag< 2 >                  Arity;

    result_type
    operator()( const Point_3& p, const Point_3& q) const
    { return CGAL_NTS compare(p.hz() * q.hw(), q.hz() * p.hw() ); }
  };

  template <typename K>
  class Compute_area_2
  {
    typedef typename K::RT                RT;
    typedef typename K::FT                FT;
    typedef typename K::Iso_rectangle_2   Iso_rectangle_2;
    typedef typename K::Triangle_2        Triangle_2;
    typedef typename K::Point_2           Point_2;
    typedef typename K::Vector_2          Vector_2;
    typedef typename K::Construct_vector_2 Construct_vector_2;
    Construct_vector_2 co;
  public:
    typedef FT               result_type;
    typedef Arity_tag< 1 >   Arity;

    FT
    operator()( const Point_2& p, const Point_2& q, const Point_2& r ) const
    {
      Vector_2 v1 = co(p, q);
      Vector_2 v2 = co(p, r);

      RT num = v1.hx()*v2.hy() - v2.hx()*v1.hy();
      RT den = RT(2) * v1.hw() * v2.hw();
      return FT(num)/FT(den);
    }

    FT
    operator()( const Iso_rectangle_2& r ) const
    { return (r.xmax()-r.xmin()) * (r.ymax()-r.ymin()); }

    FT
    operator()( const Triangle_2& t ) const
    { return t.area(); }
  };

  template <typename K>
  class Compute_determinant_2
  {
    typedef typename K::FT                FT;
    typedef typename K::Vector_2          Vector_2;
  public:
    typedef FT               result_type;
    typedef Arity_tag< 2 >   Arity;

    result_type
    operator()(const Vector_2& v, const Vector_2& w) const
    {
	return det2x2_by_formula(v.hx(), v.hy(),
                                 w.hx(), w.hy()) / FT(v.hw() * w.hw());
    }
  };

  template <typename K>
  class Compute_determinant_3
  {
    typedef typename K::FT                FT;
    typedef typename K::Vector_3          Vector_3;
  public:
    typedef FT               result_type;
    typedef Arity_tag< 3 >   Arity;

    result_type
    operator()(const Vector_3& v, const Vector_3& w, const Vector_3& t) const
    {
	return det3x3_by_formula(v.hx(), v.hy(), v.hz(),
                                 w.hx(), w.hy(), w.hz(),
                                 t.hx(), t.hy(), t.hz())
                              / FT(v.hw() * w.hw() * t.hw());
    }
  };

  template <typename K>
  class Compute_scalar_product_2
  {
    typedef typename K::RT                RT;
    typedef typename K::FT                FT;
    typedef typename K::Vector_2          Vector_2;
  public:
    typedef FT               result_type;
    typedef Arity_tag< 2 >   Arity;

    FT
    operator()(const Vector_2& v, const Vector_2& w) const
    {
        return FT( RT(v.hx()*w.hx() + v.hy()*w.hy()) ) /
               FT( RT(v.hw()*w.hw() ) );
    }
  };

  template <typename K>
  class Compute_scalar_product_3
  {
    typedef typename K::RT                RT;
    typedef typename K::FT                FT;
    typedef typename K::Vector_3          Vector_3;
  public:
    typedef FT               result_type;
    typedef Arity_tag< 2 >   Arity;

    FT
    operator()(const Vector_3& v, const Vector_3& w) const
    {
        return FT( RT(v.hx()*w.hx() + v.hy()*w.hy()) + v.hz()*w.hz() ) /
               FT( RT(v.hw()*w.hw() ) );
    }
  };

  // TODO ...
  template <typename K>
  class Compute_squared_radius_2
  {
    typedef typename K::FT          FT;
    typedef typename K::Point_2     Point_2;
    typedef typename K::Circle_2    Circle_2;
  public:
    typedef FT               result_type;
    typedef Arity_tag< 1 >   Arity;

    FT
    operator()( const Circle_2& c) const
    { return c.rep().squared_radius(); }

    FT
    operator()( const Point_2& p, const Point_2& q) const
    {
      typedef typename K::FT FT;
      return squared_distance(p, q)/FT(4);
    }  // FIXME

    FT
    operator()( const Point_2& p, const Point_2& q, const Point_2& r) const
    { return squared_distance(p, circumcenter(p, q, r)); }
    // FIXME
  };

  template <typename K>
  class Compute_squared_radius_3
  {
    typedef typename K::FT          FT;
    typedef typename K::Point_3     Point_3;
    typedef typename K::Sphere_3    Sphere_3;
  public:
    typedef FT               result_type;
    typedef Arity_tag< 1 >   Arity;

    FT
    operator()( const Sphere_3& s) const
    { return s.rep().squared_radius(); }

    FT