1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/CGAL/halfspaces_intersection.h	Thu Nov 03 15:28:12 2011 +0100
     1.3 @@ -0,0 +1,97 @@
     1.4 +#ifndef CGAL_HALFSPACES_INTERSECTION_HPP
     1.5 +#define CGAL_HALFSPACES_INTERSECTION_HPP
     1.6 +
     1.7 +#include <CGAL/Convex_hull_traits_3.h>
     1.8 +#include <CGAL/convex_hull_3.h>
     1.9 +#include <CGAL/Polyhedron_incremental_builder_3.h>
    1.10 +
    1.11 +CGAL_BEGIN_NAMESPACE
    1.12 +namespace internal
    1.13 +{
    1.14 +  template <class Polyhedron>
    1.15 +  class Build_dual_polyhedron :
    1.16 +    public CGAL::Modifier_base<typename Polyhedron::HalfedgeDS>
    1.17 +  { 
    1.18 +    typedef typename Polyhedron::HalfedgeDS HDS;
    1.19 +    const Polyhedron &_primal;
    1.20 +      
    1.21 +  public:
    1.22 +    Build_dual_polyhedron (const Polyhedron & primal):
    1.23 +      _primal (primal)
    1.24 +    {}
    1.25 +      
    1.26 +    void operator () (HDS &hds)
    1.27 +    {
    1.28 +      typedef typename Polyhedron::Facet Facet;
    1.29 +      typedef typename Polyhedron::Facet_const_handle Facet_const_handle;
    1.30 +      typedef typename Polyhedron::Facet_const_iterator Facet_const_iterator;
    1.31 +      typedef typename Polyhedron::Vertex_const_iterator Vertex_const_iterator;
    1.32 +      typename CGAL::Polyhedron_incremental_builder_3<HDS> B( hds, true);
    1.33 +      
    1.34 +      B.begin_surface(_primal.size_of_facets(),
    1.35 +		      _primal.size_of_vertices(),
    1.36 +		      _primal.size_of_vertices());
    1.37 +      
    1.38 +      // compute coordinates of extreme vertices in the dual polyhedron
    1.39 +      // from primal faces
    1.40 +      std::map<Facet_const_handle, size_t> extreme_points;
    1.41 +      size_t n = 0;
    1.42 +
    1.43 +      for (Facet_const_iterator it = _primal.facets_begin();
    1.44 +	   it != _primal.facets_end(); ++it, ++n)
    1.45 +	{
    1.46 +	  typename Facet::Halfedge_const_handle h = it->halfedge();
    1.47 +	  typename Facet::Plane_3 p ( h->vertex()->point(),
    1.48 +				      h->next()->vertex()->point(),
    1.49 +				      h->next()->next()->vertex()->point());
    1.50 +	  B.add_vertex(CGAL::ORIGIN + p.orthogonal_vector () / (-p.d()));
    1.51 +	  extreme_points[it] = n;
    1.52 +	}
    1.53 +      
    1.54 +      // build faces
    1.55 +      for (Vertex_const_iterator it = _primal.vertices_begin();
    1.56 +	   it != _primal.vertices_end(); ++it)
    1.57 +	{
    1.58 +	  assert (it->is_bivalent() == false);
    1.59 +
    1.60 +	  typename Polyhedron::Halfedge_around_vertex_const_circulator
    1.61 +	    h0 = it->vertex_begin(), hf = h0;
    1.62 +
    1.63 +	  B.begin_facet();
    1.64 +	  do
    1.65 +	    {
    1.66 +	      B.add_vertex_to_facet(extreme_points[hf->facet()]);
    1.67 +	    } while (++hf != h0);
    1.68 +	  B.end_facet();
    1.69 +	}
    1.70 +
    1.71 +      B.end_surface();
    1.72 +    }
    1.73 +  };
    1.74 +
    1.75 +  template <class PlaneIterator, class Polyhedron, class K>
    1.76 +  void
    1.77 +  intersection(PlaneIterator pbegin, PlaneIterator pend,
    1.78 +	       Polyhedron &P, const K &k)
    1.79 +  {
    1.80 +    typedef typename CGAL::Convex_hull_traits_3<K> Traits;
    1.81 +    typedef typename Polyhedron::HalfedgeDS HalfedgeDS;
    1.82 +    typedef typename K::Point_3 Point;
    1.83 +    typedef typename CGAL::internal::Build_dual_polyhedron<Polyhedron> Builder;
    1.84 +
    1.85 +    // construct dual points to apply the convex hull
    1.86 +    std::vector<Point> dual_points;
    1.87 +    for (PlaneIterator p = pbegin; p != pend; ++p)
    1.88 +      dual_points.push_back(CGAL::ORIGIN + p->orthogonal_vector () / (-p->d()));
    1.89 +
    1.90 +    Polyhedron ch;
    1.91 +    CGAL::convex_hull_3(dual_points.begin(), dual_points.end(), ch);    
    1.92 +
    1.93 +    Builder build_dual (ch);
    1.94 +    P.delegate(build_dual);
    1.95 +  }
    1.96 +}
    1.97 +
    1.98 +CGAL_END_NAMESPACE
    1.99 +
   1.100 +#endif