How to get the vertices and the faces of a polygon mesh with CGAL?

Question

With CGAL, I'm able to get a polyhedron mesh in OFF format. For example the program below builds two tetrahedra, computes their intersection, and returns the result in an OFF file. The OFF output provides the vertices coordinates and the faces given by the vertices indices.

But I would like to get the vertices and the faces as C++ variables (for example double vectors for the vertices and int vectors for the faces). Is it possible, and how?

Of course I could extract what I want from the OFF file but there's probably a better way.

#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
#include <CGAL/Polyhedron_3.h>
#include <CGAL/Nef_polyhedron_3.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/boost/graph/convert_nef_polyhedron_to_polygon_mesh.h>
#include <iostream>
#include <sstream>
#include <fstream>
typedef CGAL::Exact_predicates_exact_constructions_kernel Exact_kernel;
typedef CGAL::Polyhedron_3<Exact_kernel>                  Polyhedron;
typedef CGAL::Surface_mesh<Exact_kernel::Point_3>         Surface_mesh;
typedef CGAL::Nef_polyhedron_3<Exact_kernel>              Nef_polyhedron;
typedef Exact_kernel::Point_3                             Point_3;

int main() {
    double phi = (1+sqrt(5))/2;
    double a = 1/sqrt(3);
    double b = a/phi;
    double c = a*phi;

    Point_3 p1( 0, b, c);
    Point_3 q1( b, -c, 0);
    Point_3 r1( a, a, -a);
    Point_3 s1( -c, 0, -b);
    Polyhedron P1;
    P1.make_tetrahedron( p1, q1, r1, s1);

    Point_3 p2( a, -a, -a);
    Point_3 q2( a, a, a);
    Point_3 r2( -a, -a, a);
    Point_3 s2( -a, a, -a);
    Polyhedron P2;
    P2.make_tetrahedron( p2, q2, r2, s2);

    Nef_polyhedron nef1(P1);
    Nef_polyhedron nef2(P2);
    /* compute the intersection */
    Nef_polyhedron nef = nef1*nef2;

    Surface_mesh smesh;
    CGAL::convert_nef_polyhedron_to_polygon_mesh(nef, smesh);

    std::ofstream outfile;
    outfile.open("output.off");
    outfile << smesh;
    outfile.close();
}

Answer 1

Well, I've finally managed. I post an answer in case it could help somebody.

Vertices

#include <CGAL/Gmpq.h>

typedef Surface_mesh::Vertex_index vertex_descriptor;

double gmpq2double(CGAL::Gmpq r){
  return r.numerator().to_double()/r.denominator().to_double();
}

  printf("Number of vertices: %d\n", smesh.number_of_vertices());
  double** vertices = (double**)malloc(smesh.number_of_vertices() * sizeof(double*));
  std::cout << "Iterate over vertices\n";
  {
    unsigned i_vertex = 0;
    BOOST_FOREACH(vertex_descriptor vd, smesh.vertices()){
      std::cout << smesh.point(vd) << std::endl;
      vertices[i_vertex] = (double*)malloc(3 * sizeof(double));
      for(unsigned k=0; k < 3; k++){
        vertices[i_vertex][k] = gmpq2double(smesh.point(vd)[k].exact());
      }
      i_vertex++;
    } 
  }

Faces

This is not exactly my code and I have not tested this one. Anyway it gives the idea.

typedef Surface_mesh::Face_index face_descriptor;

  printf("Number of faces: %d\n", smesh.number_of_faces());
  std::vector<unsigned>* faces = new std::vector<unsigned>[smesh.number_of_faces()];
  std::cout << "Iterate over faces\n";
  {
    unsigned i_face = 0;
    BOOST_FOREACH(face_descriptor fd, smesh.faces()){
      BOOST_FOREACH(vertex_descriptor vd, vertices_around_face(smesh.halfedge(fd), smesh)){
        printf("vertex: %u\n", vd);
        faces[i_face].push_back(vd);
      }
      i_face++;
    }
  }

Bonus: edges

typedef Surface_mesh::Edge_index edge_descriptor;

  printf("Number of edges: %d\n", smesh.number_of_edges());
  unsigned** edges = (unsigned**)malloc(smesh.number_of_edges() * sizeof(unsigned*));
  std::cout << "Iterate over edges\n";
  {
    unsigned i_edge = 0;
    BOOST_FOREACH(edge_descriptor ed, smesh.edges()){
      edges[i_edge] = (unsigned*)malloc(2 * sizeof(unsigned));
      edges[i_edge][0] = source(ed,smesh);
      edges[i_edge][1] = target(ed,smesh);
      i_edge++;
    }
  }