SPAR/AIModule/BWTA/vendors/CGAL/CGAL/Arr_geometry_traits/Arr_plane_3.h

Go to the documentation of this file.

// Copyright (c) 2006  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.5-branch/Arrangement_on_surface_2/include/CGAL/Arr_geometry_traits/Arr_plane_3.h $
// $Id: Arr_plane_3.h 42535 2008-03-20 16:46:25Z spion $
// 
//
// Author(s)     : Efi Fogel          <efif@post.tau.ac.il>

#ifndef CGAL_ARR_PLANE_3_h
#define CGAL_ARR_PLANE_3_h

#include <CGAL/basic.h>
#include <iostream>

#include <CGAL/Kernel/solve.h>

CGAL_BEGIN_NAMESPACE

template <class Kernel_>
class Arr_plane_3 {
public:
  typedef Kernel_                       Kernel;
  typedef Arr_plane_3<Kernel>           Self;
  typedef typename Kernel::Point_2      Point_2;
  typedef typename Kernel::Point_3      Point_3;
  typedef typename Kernel::Vector_3     Vector_3;
  typedef typename Kernel::Direction_3  Direction_3;
  typedef typename Kernel::Line_3       Line_3;
  typedef typename Kernel::FT           FT;

private:
  FT m_a;

  FT m_b;

  FT m_c;
  
public:
  Arr_plane_3() : m_a(0), m_b(0), m_c(0) {}

  Arr_plane_3(int a, int b, int c) : m_a(a), m_b(b), m_c(c) {}

  Arr_plane_3(typename Kernel::Plane_3 p)
  {
    CGAL_precondition_code(Kernel kernel;);
    CGAL_precondition_code(typename Kernel::Point_3 orig = kernel.construct_point_3_object()(ORIGIN););
    CGAL_precondition(kernel.has_on_3_object()(p, orig));

    m_a = p.a(); m_b = p.b(); m_c = p.c() ;
  }

  Arr_plane_3(const Point_3 & p, const Point_3 & r)
  {
    FT prx = r.x() - p.x();
    FT pry = r.y() - p.y();
    FT prz = r.z() - p.z();
    m_a = r.y() * prz - pry * r.z();
    m_b = r.z() * prx - prz * r.x();
    m_c = r.x() * pry - prx * r.y();
  }

  const FT & a() const { return m_a; }

  const FT & b() const { return m_b; }

  const FT & c() const { return m_c; }
  
  FT operator[](unsigned int i) const
  {
    CGAL_assertion(i < 3);
    return (i == 0) ? m_a : ((i == 1) ? m_b : m_c);
  }

  bool equal(Self plane) const
  {
    return ((a() == plane.a()) &&
            (b() == plane.b()) &&
            (c() == plane.c()));
  }

  operator typename Kernel::Plane_3 () const
  {
    Kernel kernel;
    return kernel.construct_plane_3_object() (m_a, m_b, m_c, 0);
  }

  Point_2 to_2d(const Point_3 & p) const
  {
    Kernel kernel;
    typename Kernel::Plane_3 base_plane(m_a, m_b, m_c, 0);
    Vector_3 v1 = kernel.construct_base_vector_3_object()(base_plane, 1);    
    Vector_3 v2 = kernel.construct_base_vector_3_object()(base_plane, 2);    
    Vector_3 v3 = kernel.construct_orthogonal_vector_3_object()(base_plane);

    FT denom = v2[1]*v3[0]*v1[2] - v2[0]*v3[1]*v1[2] + v3[2]*v2[0]*v1[1] +
               v2[2]*v3[1]*v1[0] - v3[0]*v2[2]*v1[1] - v2[1]*v3[2]*v1[0];
    FT x = - (v2[1]*v3[2]*p[0] - v2[1]*v3[0]*p[2] + v3[0]*v2[2]*p[1] +
              v2[0]*v3[1]*p[2] - v3[2]*v2[0]*p[1] - v2[2]*v3[1]*p[0]) / denom;
    FT y = (v1[1]*v3[2]*p[0] - v1[1]*v3[0]*p[2] - v3[1]*p[0]*v1[2] +
            v3[1]*v1[0]*p[2] + p[1]*v3[0]*v1[2] - p[1]*v3[2]*v1[0]) / denom;
    
    return Point_2(x, y);
  }

  Point_3 to_3d(const Point_2 & p_2, unsigned int i) const
  {
    CGAL_assertion(i < 3);

#if 0
    std::cout << "(a, b, c): " << a() << "," << b() << "," << c()
              << std::endl;
#endif
    
    if (i == 0) {
      CGAL_assertion(m_a != 0);
      FT y = p_2.x();
      FT z = p_2.y();
      FT x = -(m_b * y + m_c * z) / m_a;
      Point_3 p_3(x, y, z);
      return p_3;
    }

    if (i == 1) {
      CGAL_assertion(m_b != 0);
      FT z = p_2.x();
      FT x = p_2.y();
      FT y = -(m_a * x + m_c * z) / m_b;
      Point_3 p_3(x, y, z);
      return p_3;
    }

    // if (i == 2) return Base::to_3d(p_2);
    CGAL_assertion(m_c != 0);
    FT x = p_2.x();
    FT y = p_2.y();
    FT z = -(m_a * x + m_b * y) / m_c;
    Point_3 p_3(x, y, z);
    return p_3;
  }

  Oriented_side oriented_side(const Point_3 & p) const
  {
    return CGAL_NTS sign(m_a*p.x() + m_b*p.y() + m_c*p.z());
  }
};

template <class Kernel>
CGAL::Object intersect(const Arr_plane_3<Kernel> & plane1,
                       const Arr_plane_3<Kernel> & plane2)
{
  typedef typename Kernel::Point_3      Point_3;
  typedef typename Kernel::Direction_3  Direction_3;
  typedef typename Kernel::Line_3       Line_3;
  typedef typename Kernel::FT           FT;

  // We know that the plane goes throgh the origin
  const FT & a1 = plane1.a();
  const FT & b1 = plane1.b();
  const FT & c1 = plane1.c();

  const FT & a2 = plane2.a();
  const FT & b2 = plane2.b();
  const FT & c2 = plane2.c();

  FT det = a1*b2 - a2*b1;
  if (det != 0) {
    Point_3 is_pt = Point_3(0, 0, 0, det);
    Direction_3 is_dir = Direction_3(b1*c2 - c1*b2, a2*c1 - a1*c2, det);
    return make_object(Line_3(is_pt, is_dir));
  }

  det = a1*c2 - a2*c1;
  if (det != 0) {
    Point_3 is_pt = Point_3(0, 0, 0, det);
    Direction_3 is_dir = Direction_3(c1*b2 - b1*c2, det, a2*b1 - a1*b2);
    return make_object(Line_3(is_pt, is_dir));
  }
  det = b1*c2 - c1*b2;
  if (det != 0) {
    Point_3 is_pt = Point_3(0, 0, 0, det);
    Direction_3 is_dir = Direction_3(det, c1*a2 - a1*c2, a1*b2 - b1*a2);
    return make_object(Line_3(is_pt, is_dir));
  }

  // degenerate case
  return make_object(plane1);
}

template <class Kernel>
typename Kernel::Point_2
construct_projected_xy_point(const Arr_plane_3<Kernel> & plane,
                             const typename Kernel::Point_3 & p)
{
  return plane.to_2d(p);
}

template <class Kernel>
inline std::ostream & operator<<(std::ostream & os,
                                 const Arr_plane_3<Kernel> & plane)
{
  os << plane[0] << ", " << plane[1] << ", " << plane[2];
  return os;
}  

CGAL_END_NAMESPACE

#endif