CGAL users discussion list

[Martin Baeker <>]

Sorry, I forgot about the headers (they are pretty standard anyway).

The full program follows after the line:


---------------------------------------------------------------



/// Trying to understand facet planes


#include <CGAL/Gmpq.h>
#include <CGAL/Cartesian.h>
#include <CGAL/Extended_cartesian.h>

#include <CGAL/Polyhedron_3.h>
#include <CGAL/intersections.h>
#include <CGAL/IO/Polyhedron_iostream.h>
#include <CGAL/Nef_polyhedron_3.h>
#include <CGAL/IO/Nef_polyhedron_iostream_3.h>
#include <iostream>
#include <fstream>
#include <CGAL/Cartesian_converter.h>
#include <stdarg.h>
#include <stdlib.h>


typedef CGAL::Cartesian<CGAL::Gmpq>   Kernel;
typedef CGAL::Polyhedron_3<Kernel>                           Polyhedron;
typedef Kernel::Plane_3  Plane_3;
typedef Kernel::Point_3    Point;
typedef Kernel::Segment_3    Segment;
typedef Kernel::Vector_3  Vector_3;
typedef Polyhedron::Halfedge_handle                          Halfedge_handle;
typedef Polyhedron::Facet_handle                          Facet_handle;
typedef Polyhedron::Facet                          Facet;


#define pretty(p) " " << ( #p ) << ": " << CGAL::to_double( p.x() ) << " " << CGAL::to_double( 
p.y() ) << " " << CGAL::to_double( p.z() )
#define OUT(variable) " " << ( #variable ) << ": " << ( variable )
#define WRITE(variable) std::cout << OUT(variable) << std::endl

const double INVALID=-1.e20;
const double TOLERANCE=1.e-5;



typedef Polyhedron::Facet_iterator                   Facet_iterator;
typedef Polyhedron::Halfedge_around_facet_circulator 
Halfedge_facet_circulator;
typedef Polyhedron::Halfedge_iterator Halfedge_iterator;

template <class Poly>
typename Poly::Halfedge_handle make_cube_3( Poly& P) {
    // appends a cube of size [0,1]^3 to the polyhedron P.
    CGAL_precondition( P.is_valid());
    typedef typename Poly::Point_3         Point;
    typedef typename Poly::Plane_3         Plane;
    typedef typename Poly::Halfedge_handle Halfedge_handle;
    Halfedge_handle h = P.make_tetrahedron( Point( 1, 0, 0),
                                            Point( 0, 0, 1),
                                            Point( 0, 0, 0),
                                            Point( 0, 1, 0));
    Halfedge_handle g = h->next()->opposite()->next();
    P.split_edge( h->next());
    P.split_edge( g->next());
    P.split_edge( g);
    h->next()->vertex()->point()     = Point( 1, 0, 1);
    g->next()->vertex()->point()     = Point( 0, 1, 1);
    g->opposite()->vertex()->point() = Point( 1, 1, 0);
    Halfedge_handle f = P.split_facet( g->next(), g->next()->next()->next());
    Halfedge_handle e = P.split_edge( f);
    e->vertex()->point() = Point( 1, 1, 1);
    P.split_facet( e, f->next()->next());
    CGAL_postcondition( P.is_valid());
    g = h;
    g->facet()->plane() = Plane( g->vertex()->point(),
                                 g->next()->vertex()->point(),
                                 g->next()->next()->vertex()->point());
    g = h->opposite();
    g->facet()->plane() = Plane( g->vertex()->point(),
                                 g->next()->vertex()->point(),
                                 g->next()->next()->vertex()->point());
    g = h->next()->opposite();
    g->facet()->plane() = Plane( g->vertex()->point(),
                                 g->next()->vertex()->point(),
                                 g->next()->next()->vertex()->point());
    g = h->next()->next()->opposite();
    g->facet()->plane() = Plane( g->vertex()->point(),
                                 g->next()->vertex()->point(),
                                 g->next()->next()->vertex()->point());
    g = h->next()->next()->next()->opposite();
    g->facet()->plane() = Plane( g->vertex()->point(),
                                 g->next()->vertex()->point(),
                                 g->next()->next()->vertex()->point());
    g = g->next()->next()->opposite();
    g->facet()->plane() = Plane( g->vertex()->point(),
                                 g->next()->vertex()->point(),
                                 g->next()->next()->vertex()->point());
    return h;
}

void do_stuff(Polyhedron & pcube)
{
  for (Halfedge_iterator h=pcube.halfedges_begin(); h != 
pcube.halfedges_end(); ++h)
    {
      Facet face=*(h->facet());
      std::cout << face.plane().orthogonal_vector() << std::endl;
    }
}


int main() 
{
  Polyhedron pcube;
  make_cube_3(pcube);
  std::ofstream os1("poly2.off" ) ;
  os1 << pcube ;
  os1.close();
  do_stuff(pcube);
  std::cout << "\n";

  Polyhedron pcube2;
  std::ifstream is1("poly2.off" ) ;
  is1 >> pcube2;
  do_stuff(pcube2);



  return 0;
}



> Martin Baeker wrote:
> > Dear all,
> >
> > I'm getting a bit desperate while trying to work with Polyhedra-stuff.
> >
> > The manual says that I can access the plane that belongs to a facet with
> > Plane_3&         f.plane () if Supports_facet_plane  == CGAL::Tag_true
> >
> > However, I do not understand how this really works.
> >
> > When I create a Polyhedron, e.g., with the make_cube_3-routine that is
> > provided in the examples, I can access the plane corresponding to each
> > facet without problems.
> > However, when I write this Polyhedron to an off-file and read it in
> > again, I cannot access its planes anymore.
> >
> > Here is a simple test program to show this:
> >
> > void do_stuff(Polyhedron & pcube)
> > {
> >   for (Halfedge_iterator h=pcube.halfedges_begin(); h != 
> > pcube.halfedges_end(); ++h)
> >     {
> >       Facet face=*(h->facet());
> >       std::cout << face.plane().orthogonal_vector() << std::endl;
> >     }
> > }
> >
> > int main() {
> >   Polyhedron pcube;
> >   make_cube_3(pcube);
> >   std::ofstream os1("poly2.off" ) ;
> >   os1 << pcube ;
> >   os1.close();
> >   do_stuff(pcube);
> >   std::cout << "\n";
> >   Polyhedron pcube2;
> >   std::ifstream is1("poly2.off" ) ;
> >   is1 >> pcube2;
> >   do_stuff(pcube2);
> > }
> >
> > The output for pcube is correct; for pcube2, all orthogonal vectors are 
> > output as 0/1 0/1 0/1.
> >
> >
> > Is there a way to make sure that I can always access the plane
> > belonging to a facet? If not, why not? Could I do it by just taking
> > three points from the facet and create a plane with them?
> >
> >
> > Any help will be greatly appreciated,
> >
> > Martin.
> >
> >
> >
> >
> >
> >                    Priv.-Doz. Dr. Martin Bäker
> >                    Institut für Werkstoffe
> >                    Technische Universität Braunschweig
> >                    Langer Kamp 8
> >                    38106 Braunschweig
> >                    Germany
> >                    Tel.: 00-49-531-391-3073
> >                    Fax   00-49-531-391-3058
> >                    e-mail 
> > <>
>
> Can you please post the typedefs as well, so as to make your example 
> compilable.
>
>
> S.
>
>
>
> --
> You are currently subscribed to cgal-discuss.
> To unsubscribe or access the archives, go to
> https://lists-sop.inria.fr/wws/info/cgal-discuss

                   Priv.-Doz. Dr. Martin Bäker
                   Institut für Werkstoffe
                   Technische Universität Braunschweig
                   Langer Kamp 8
                   38106 Braunschweig
                   Germany
                   Tel.: 00-49-531-391-3073
                   Fax   00-49-531-391-3058
                   e-mail 
<>