CGAL users discussion list

[Laurent Rineau <>]

On Sunday 13 April 2008 16:21:50 Philipp Jenke wrote:
> Hi Pierre,
>
> Polygon soup is ok, but perfect would be a vertex list and three
> indices for each triangle (shared vertex for all adjacent triangles).
>
> Thanks,
> Philipp

Hi Philipp,

The result of the CGAL surface mesh generator is a data structure embedded in 
a 3D triangulation, named the "2D Complex in 3D Triangulation", and 
represented by the class 
CGAL::Surface_mesh_complex_2_in_triangulation_3<Triangulation>. The 
documentation is available at that place:
http://www.cgal.org/Manual/3.3/doc_html/cgal_manual/Surface_mesher_ref/Class_Surface_mesh_complex_2_in_triangulation_3.html

Actually that reference page only says that the class is model of a concept, 
so the actual list of types and methods available in that class is here:
http://www.cgal.org/Manual/3.3/doc_html/cgal_manual/Surface_mesher_ref/Concept_SurfaceMeshComplex_2InTriangulation_3.html

What you need is the section "Traversal of the complex", and particularly the 
Facet_iterator.

To use that class, you also need to know a little about CGAL triangulations 
(types and methods to deal with vertices, cells and facets).
http://www.cgal.org/Manual/3.3/doc_html/cgal_manual/Triangulation_3_ref/Class_Triangulation_3.html


I have attached a modified version of the example 
examples/Surface_mesher/mesh_an_implicit_function.cpp you quoted in an other 
thread. After the computation, the code outputs the vertices and facets 
lists. A second attached modified version outputs the surface mesh in an OFF 
file.

-- 
Laurent Rineau
INRIA - Sophia Antipolis
BP 93, 2004 Route des Lucioles
06902 Sophia Antipolis Cedex FRANCE
Tel: +33 4 92 38 78 62 (Fax: +33.4.97.15.53.95)
#include <CGAL/Surface_mesh_default_triangulation_3.h>
#include <CGAL/Complex_2_in_triangulation_3.h>
#include <CGAL/make_surface_mesh.h>
#include <CGAL/Implicit_surface_3.h>
#include <fstream>

// default triangulation for Surface_mesher
typedef CGAL::Surface_mesh_default_triangulation_3 Tr;

// c2t3
typedef CGAL::Complex_2_in_triangulation_3<Tr> C2t3;

typedef Tr::Geom_traits GT;
typedef Tr::Vertex_handle Vertex_handle;
typedef GT::Sphere_3 Sphere_3;
typedef GT::Point_3 Point_3;
typedef GT::FT FT;

typedef FT (*Function)(Point_3);

typedef CGAL::Implicit_surface_3<GT, Function> Surface_3;

FT sphere_function (Point_3 p) {
  const FT x2=p.x()*p.x(), y2=p.y()*p.y(), z2=p.z()*p.z();
  return x2+y2+z2-1;
}

int main() {
  Tr tr;            // 3D-Delaunay triangulation
  C2t3 c2t3 (tr);   // 2D-complex in 3D-Delaunay triangulation

  // defining the surface
  Surface_3 surface(sphere_function,             // pointer to function
                    Sphere_3(CGAL::ORIGIN, 2.)); // bounding sphere

  // defining meshing criteria
  CGAL::Surface_mesh_default_criteria_3<Tr> criteria(30.,  // angular bound
                                                     0.1,  // radius bound
                                                     0.1); // distance bound
  // meshing surface
  CGAL::make_surface_mesh(c2t3, surface, criteria, CGAL::Non_manifold_tag());

  std::ofstream out_file("output.off");

  out_file.imbue(std::locale("POSIX")); // use a POSIX locale
  
  out_file << "OFF\n" 
           << tr.number_of_vertices() << " "
           << c2t3.number_of_facets() << " 0\n";

  std::map<Vertex_handle, int> vertices_map;
  int vertices_counter = 0;

  // traversal of the vertices list
  for(Tr::Finite_vertices_iterator 
        vit = tr.finite_vertices_begin(),
        end = tr.finite_vertices_end();
      vit != end; ++vit)
  {
    // note: a vertex iterator is implicitly convertable to a vertex handle
    vertices_map[vit] = vertices_counter;
    out_file << vit->point() << "\n";
    ++vertices_counter;
  }

  // traversal of the 2D complex
  for(C2t3::Facet_iterator 
        fit = c2t3.facets_begin(),
        end = c2t3.facets_end();
      fit != end; ++fit)
  {
    // a facet is a pair (cell_handle, index)
    const Tr::Cell_handle& cell = fit->first;
    const int index = fit->second;

    const Vertex_handle va = cell->vertex((index+1)&3);
    const Vertex_handle vb = cell->vertex((index+2)&3);
    const Vertex_handle vc = cell->vertex((index+3)&3);

    out_file << "3 " 
             << vertices_map[va] << " "
             << vertices_map[vb] << " "
             << vertices_map[vc] << "\n";
  }
}
#include <CGAL/Surface_mesh_default_triangulation_3.h>
#include <CGAL/Complex_2_in_triangulation_3.h>
#include <CGAL/make_surface_mesh.h>
#include <CGAL/Implicit_surface_3.h>

// default triangulation for Surface_mesher
typedef CGAL::Surface_mesh_default_triangulation_3 Tr;

// c2t3
typedef CGAL::Complex_2_in_triangulation_3<Tr> C2t3;

typedef Tr::Geom_traits GT;
typedef Tr::Vertex_handle Vertex_handle;
typedef GT::Sphere_3 Sphere_3;
typedef GT::Point_3 Point_3;
typedef GT::FT FT;

typedef FT (*Function)(Point_3);

typedef CGAL::Implicit_surface_3<GT, Function> Surface_3;

FT sphere_function (Point_3 p) {
  const FT x2=p.x()*p.x(), y2=p.y()*p.y(), z2=p.z()*p.z();
  return x2+y2+z2-1;
}

int main() {
  Tr tr;            // 3D-Delaunay triangulation
  C2t3 c2t3 (tr);   // 2D-complex in 3D-Delaunay triangulation

  // defining the surface
  Surface_3 surface(sphere_function,             // pointer to function
                    Sphere_3(CGAL::ORIGIN, 2.)); // bounding sphere

  // defining meshing criteria
  CGAL::Surface_mesh_default_criteria_3<Tr> criteria(30.,  // angular bound
                                                     0.1,  // radius bound
                                                     0.1); // distance bound
  // meshing surface
  CGAL::make_surface_mesh(c2t3, surface, criteria, CGAL::Non_manifold_tag());

  std::cout << "Vertices list... ("
            << tr.number_of_vertices() << " vertices)\n";
  
  std::map<Vertex_handle, int> vertices_map;
  int vertices_counter = 0;

  // traversal of the vertices list
  for(Tr::Finite_vertices_iterator 
        vit = tr.finite_vertices_begin(),
        end = tr.finite_vertices_end();
      vit != end; ++vit)
  {
    // note: a vertex iterator is implicitly convertable to a vertex handle
    vertices_map[vit] = vertices_counter;
    std::cout << "vertex #" << vertices_counter << " = ("
              << vit->point() << ")\n";
    ++vertices_counter;
  }

  std::cout << "Facets list... ("
            << c2t3.number_of_facets() << " facets)\n";

  // traversal of the 2D complex
  for(C2t3::Facet_iterator 
        fit = c2t3.facets_begin(),
        end = c2t3.facets_end();
      fit != end; ++fit)
  {
    // a facet is a pair (cell_handle, index)
    const Tr::Cell_handle& cell = fit->first;
    const int index = fit->second;

    const Vertex_handle va = cell->vertex((index+1)&3);
    const Vertex_handle vb = cell->vertex((index+2)&3);
    const Vertex_handle vc = cell->vertex((index+3)&3);

    std::cout << "facet " 
              << vertices_map[va] << " "
              << vertices_map[vb] << " "
              << vertices_map[vc] << "\n";
  }
}