CGAL users discussion list

[Laurent Rineau <>]

On Thursday 16 October 2008 05:11:56 Sebastian Schoellhammer wrote:
> Hello Dear List,
>
> I'm a total beginner with CGAL (and I have to say, a beginning programmer
> as well) but I would love to extract triangles from this remeshing example.
> Here is what I am doing:
>
> // meshing surface, with the "manifold without boundary" algorithm
>   CGAL::make_surface_mesh(c2t3, surface, criteria, CGAL::Manifold_tag());

Hi Sebastien,

I am the author of the Surface mesh generation package (and I now work at 
GeometryFactory).

The examples in CGAL-3.3.1 are too short, and they will be improved in next 
release. Here is an extended examples/Surface_mesher/mesh_a_3d_gray_image.cpp 
(see attachments), where I have added several lines of code that show how to 
get the facets after the meshing. Output facets are not oriented. If you use 
OpenGL to display the facets, you will need to make OpenGL ignore the 
orientation. The example also show the use of a non-documented function of 
CGAL-3.3.1, which output the 2D complex in an OFF file. OFF files can be 
viewed by the tool geomview (among other software).

-- 
Laurent Rineau, PhD
Engineer at GeometryFactory
http://www.geometryfactory.com/
#include <CGAL/Surface_mesh_default_triangulation_3.h>
#include <CGAL/Complex_2_in_triangulation_3.h>
#include <CGAL/make_surface_mesh.h>
#include <CGAL/Gray_level_image_3.h>
#include <CGAL/Implicit_surface_3.h>

#include <fstream>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.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 CGAL::Gray_level_image_3<GT::FT, GT::Point_3> Gray_level_image;
typedef CGAL::Implicit_surface_3<GT, Gray_level_image> Surface_3;

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

  // the 'function' is a 3D gray level image
  Gray_level_image image("data/skull_2.9.inr", 2.9);

  // Carefully choosen bounding sphere: the center must be inside the
  // surface defined by 'image' and the radius must be high enough so that
  // the sphere actually bounds the whole image.
  GT::Point_3 bounding_sphere_center(122., 102., 117.);
  GT::FT bounding_sphere_squared_radius = 200.*200.*2.;
  GT::Sphere_3 bounding_sphere(bounding_sphere_center,
                                   bounding_sphere_squared_radius);

  // definition of the surface, with 10^-5 as relative precision
  Surface_3 surface(image, bounding_sphere, 1e-5);

  // defining meshing criteria
  CGAL::Surface_mesh_default_criteria_3<Tr> criteria(30.,
                                                     5.,
                                                     5.);

  // meshing surface, with the "manifold without boundary" algorithm
  CGAL::make_surface_mesh(c2t3, surface, criteria, CGAL::Manifold_tag());

  std::vector<double> vertices;

  // for all facets of the 2D complex...
  for(C2t3::Facet_iterator fit = c2t3.facets_begin(), end = c2t3.facets_end();
      fit != end; ++fit)
  {
    const Tr::Cell_handle& facet_cell = fit->first;
    const int& facet_index = fit->second;

    // ... get all vertices of the facet...
    for(int i = 1; i <= 3; ++i)
    {
      const GT::Point_3& p = facet_cell->vertex((facet_index+i)&3)->point();

      // ...and store their coordinates in an array of double.
      vertices.push_back(p.x());
      vertices.push_back(p.y());
      vertices.push_back(p.z());
    }
  }
  // At this point, &vertices[0] is a pointer to an array of double than
  // can be passed to OpenGL. Facets are not oriented.

  // Output the 2D complex to an OFF file.
  std::ofstream out("out.off");
  CGAL::output_surface_facets_to_off (out, c2t3);

  std::cout << "Final number of points: " << tr.number_of_vertices() << "\n";
}