\( \newcommand{\E}{\mathrm{E}} \) \( \newcommand{\A}{\mathrm{A}} \) \( \newcommand{\R}{\mathrm{R}} \) \( \newcommand{\N}{\mathrm{N}} \) \( \newcommand{\Q}{\mathrm{Q}} \) \( \newcommand{\Z}{\mathrm{Z}} \) \( \def\ccSum #1#2#3{ \sum_{#1}^{#2}{#3} } \def\ccProd #1#2#3{ \sum_{#1}^{#2}{#3} }\)
CGAL 4.13 - Surface Mesh
User Manual

Author
Mario Botsch, Daniel Sieger, Philipp Moeller and Andreas Fabri
clown_fish.jpg

The class Surface_mesh is an implementation of a halfedge data structure and can be used to represent a polyhedral surface. It is an alternative to the CGAL packages Halfedge Data Structures and 3D Polyhedral Surface. The main difference is that it is indexed based and not pointer based. Additionally, the mechanism for adding information to vertices, halfedges, edges, and faces is much simpler and is done at runtime and not at compile time.

Because the data structure uses integer indices as descriptors for vertices, halfedges, edges and faces it has a lower memory footprint than a 64-bit pointer based version. As the indices are contiguous, they can be used as index into vectors which store properties.

When elements are removed, they are only marked as removed, and a garbage collection function must be called to really remove them.

The class Surface_mesh can be used through its class member functions as well as through the BGL API as described in the package CGAL and the Boost Graph Library, as it is a model of the concepts MutableFaceGraph and FaceListGraph. Therefore it is possible to apply the algorithms of the packages Triangulated Surface Mesh Simplification, Triangulated Surface Mesh Segmentation, and Triangulated Surface Mesh Deformation on a surface mesh.

Usage

The main class Surface_mesh provides four nested classes that represent the basic elements of the halfedge data structure:

These types are just wrappers for an integer and their main purpose is to guarantee type safety. They are default constructible, which yields an invalid element. New elements can be added and removed to the Surface_mesh through a set of low-level functions which do not maintain connectivity. One exception is Surface_mesh::add_face(), which tries to add a new face to the mesh (defined by a sequence of vertices), and fails if the operation is not topologically valid.

typedef Surface_mesh<Point> Mesh;
Mesh m;
Mesh::Vertex_index u = m.add_vertex(Point(0,1,0));
Mesh::Vertex_index v = m.add_vertex(Point(0,0,0));
Mesh::Vertex_index w = m.add_vertex(Point(1,0,0));
m.add_face(u, v, w);

As Surface_mesh is index-based Vertex_index, Halfedge_index, Edge_index, and Face_index do not have member functions to access connectivity or properties. The functions of the Surface_mesh instance from which they were created must be used to obtain this information.

Connectivity

A surface mesh is an edge-centered data structure capable of maintaining incidence information of vertices, edges, and faces. Each edge is represented by two halfedges with opposite orientation. Each halfedge stores a reference to an incident face and to an incident vertex. Additionally, it stores a reference to the next and previous halfedge incident to its incident face. For each face and each vertex an incident halfedge is stored. Halfedges do not store the index of the opposite halfedge, as Surface_mesh stores opposite halfedges consecutively in memory.

The following figure illustrates the functions which allow to navigate in a surface mesh: Surface_mesh::opposite(), Surface_mesh::next(), Surface_mesh::prev(), Surface_mesh::target(), and Surface_mesh::face(). Additionally, the functions Surface_mesh::halfedge() allows to obtain the halfedge associated to a vertex and to a face. Alternatively, one may use the free functions with the same names defined in the package CGAL and the Boost Graph Library.

connectivity.svg
Figure 27.1 Connectivity of halfedges and vertices in a surface mesh seen from outside.

The halfedges incident to a face form a cycle. Depending on from which side we look at the surface, the sequence of halfedges appears to be oriented clockwise or counterclockwise. When in this manual we speak about the orientation of a traversal then we look at the surface such that the halfedges around a face are oriented counterclockwise, as illustrated in Figure 27.1

The connectivity does not allow to represent faces with holes.

Ranges and Iterators

Surface_mesh provides iterator ranges to enumerate all vertices, halfedges, edges, and faces. It provides member functions returning ranges of elements which are compatible with the Boost.Range library.

Example

The following example shows how to obtain the iterator type from a range, alternatives for obtaining the begin and end iterator, and alternatives for range-based loops.


File Surface_mesh/sm_iterators.cpp

#include <vector>
#include <boost/foreach.hpp>
#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
typedef Mesh::Vertex_index vertex_descriptor;
typedef Mesh::Face_index face_descriptor;
int main()
{
Mesh m;
// u x
// +------------+
// | |
// | |
// | f |
// | |
// | |
// +------------+
// v w
// Add the points as vertices
vertex_descriptor u = m.add_vertex(K::Point_3(0,1,0));
vertex_descriptor v = m.add_vertex(K::Point_3(0,0,0));
vertex_descriptor w = m.add_vertex(K::Point_3(1,0,0));
vertex_descriptor x = m.add_vertex(K::Point_3(1,1,0));
/* face_descriptor f = */ m.add_face(u,v,w,x);
{
std::cout << "all vertices " << std::endl;
// The vertex iterator type is a nested type of the Vertex_range
Mesh::Vertex_range::iterator vb, ve;
Mesh::Vertex_range r = m.vertices();
// The iterators can be accessed through the C++ range API
vb = r.begin();
ve = r.end();
// or the boost Range API
vb = boost::begin(r);
ve = boost::end(r);
// or with boost::tie, as the CGAL range derives from std::pair
for(boost::tie(vb, ve) = m.vertices(); vb != ve; ++vb){
std::cout << *vb << std::endl;
}
// Instead of the classical for loop one can use
// the boost macro for a range
BOOST_FOREACH(vertex_descriptor vd, m.vertices()){
std::cout << vd << std::endl;
}
// or the C++11 for loop. Note that there is a ':' and not a ',' as in BOOST_FOREACH
#ifndef CGAL_CFG_NO_CPP0X_RANGE_BASED_FOR
for(vertex_descriptor vd : m.vertices()){
std::cout << vd << std::endl;
}
#endif
}
return 0;
}

Circulators

Circulators around faces and around vertices are provided as class templates in the package CGAL and the Boost Graph Library.

Circulators around faces basically call Surface_mesh::next() in order to go from halfedge to halfedge counterclockwise around the face, and when dereferenced return the halfedge or the incident vertex or the opposite face.

Circulators around the target vertex of an edge basically call Surface_mesh::opposite(Surface_mesh::next()) in order to go from halfedge to halfedge clockwise around the same target vertex.

All circulators model BidirectionalCirculator. In addition to that they also support a conversion to bool for more convenient checking of emptiness.

Example

The following example shows how to enumerate the vertices around the target of a given halfedge. The second loop shows that each of these circulator types comes with an equivalent iterator and a free function to create an iterator range.


File Surface_mesh/sm_circulators.cpp

#include <vector>
#include <boost/foreach.hpp>
#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
typedef Mesh::Vertex_index vertex_descriptor;
typedef Mesh::Face_index face_descriptor;
int main()
{
Mesh m;
// u x
// +------------+
// | |
// | |
// | f |
// | |
// | |
// +------------+
// v w
// Add the points as vertices
vertex_descriptor u = m.add_vertex(K::Point_3(0,1,0));
vertex_descriptor v = m.add_vertex(K::Point_3(0,0,0));
vertex_descriptor w = m.add_vertex(K::Point_3(1,0,0));
vertex_descriptor x = m.add_vertex(K::Point_3(1,1,0));
face_descriptor f = m.add_face(u,v,w,x);
{
std::cout << "vertices around vertex " << v << std::endl;
CGAL::Vertex_around_target_circulator<Mesh> vbegin(m.halfedge(v),m), done(vbegin);
do {
std::cout << *vbegin++ << std::endl;
} while(vbegin != done);
}
{
std::cout << "vertices around face " << f << std::endl;
for(boost::tie(vbegin, vend) = vertices_around_face(m.halfedge(f), m);
vbegin != vend;
++vbegin){
std::cout << *vbegin << std::endl;
}
}
// or the same again, but directly with a range based loop
BOOST_FOREACH(vertex_descriptor vd,vertices_around_face(m.halfedge(f), m)){
std::cout << vd << std::endl;
}
return 0;
}

Properties

Surface_mesh provides a mechanism to specify new properties for vertices, halfedges, edges, and faces at run-time. Each property is identified by a string and its key type. All the values of a given property are stored as consecutive blocks of memory. References to properties are invalidated whenever new elements of the key type are added to the data-structure or when the function Surface_mesh::collect_garbage() is performed. Properties of an element will continue to exist after the element has been deleted. Trying to access a property through an invalidated element will result in undefined behavior.

One property is maintained by default, namely "v:point". The value of this property has to be supplied when adding a new point to the data structure via Surface_mesh::add_vertex(). The property can be directly accessed using Surface_mesh::points() or Surface_mesh::point(Surface_mesh::Vertex_index v).

When an element is removed, it is only marked as removed, and it gets really removed when Surface_mesh::collect_garbage() is called. Garbage collection will also really remove the properties of these elements.

The connectivity is also stored in properties, namely the properties named "v:connectivity", "h:connectivity", and "f:connectivity". It is quite similar for the marker of deleted element, where we have "v:removed", "e:removed", and "f:removed".

Example

This example shows how to use the most common features of the property system.


File Surface_mesh/sm_properties.cpp

#include <string>
#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
#include <boost/foreach.hpp>
typedef Mesh::Vertex_index vertex_descriptor;
typedef Mesh::Face_index face_descriptor;
int main()
{
Mesh m;
vertex_descriptor v0 = m.add_vertex(K::Point_3(0,2,0));
vertex_descriptor v1 = m.add_vertex(K::Point_3(2,2,0));
vertex_descriptor v2 = m.add_vertex(K::Point_3(0,0,0));
vertex_descriptor v3 = m.add_vertex(K::Point_3(2,0,0));
vertex_descriptor v4 = m.add_vertex(K::Point_3(1,1,0));
m.add_face(v3, v1, v4);
m.add_face(v0, v4, v1);
m.add_face(v0, v2, v4);
m.add_face(v2, v3, v4);
// give each vertex a name, the default is empty
Mesh::Property_map<vertex_descriptor,std::string> name;
bool created;
boost::tie(name, created) = m.add_property_map<vertex_descriptor,std::string>("v:name","");
assert(created);
// add some names to the vertices
name[v0] = "hello";
name[v2] = "world";
{
// You get an existing property, and created will be false
Mesh::Property_map<vertex_descriptor,std::string> name;
bool created;
boost::tie(name, created) = m.add_property_map<vertex_descriptor,std::string>("v:name", "");
assert(! created);
}
// You can't get a property that does not exist
Mesh::Property_map<face_descriptor,std::string> gnus;
bool found;
boost::tie(gnus, found) = m.property_map<face_descriptor,std::string>("v:gnus");
assert(! found);
// retrieve the point property for which exists a convenience function
Mesh::Property_map<vertex_descriptor, K::Point_3> location = m.points();
BOOST_FOREACH( vertex_descriptor vd, m.vertices()) {
std::cout << name[vd] << " @ " << location[vd] << std::endl;
}
std::vector<std::string> props = m.properties<vertex_descriptor>();
BOOST_FOREACH(std::string p, props){
std::cout << p << std::endl;
}
// delete the string property again
m.remove_property_map(name);
return 0;
}

Borders

A halfedge stores a reference to a face, its incident face. A halfedge h is on the border, if it has no incident face, that is if sm.face(h) == Surface_mesh::null_face(). An edge is on the border, if any of its halfedges is on the border. A vertex is on the border, if any of its incident halfedges is on the border.

A vertex has only one associated halfedge. If the user takes care that the associated halfedge is a border halfedge, in case the vertex is on the border, there is no need to look at all incident halfedges in the is_border() function for vertices. In order to only check if the associated halfedge is on the border the function Surface_mesh::is_border(Vertex_index v, bool check_all_incident_halfedges = true) must be called with check_all_incident_halfedges = false.

The user is in charge to correctly set the halfedge associated to a vertex after having applied an operation that might invalidate this property. The functions Surface_mesh::set_vertex_halfedge_to_border_halfedge(Vertex_index v), Surface_mesh::set_vertex_halfedge_to_border_halfedge(Halfedge_index h), and Surface_mesh::set_vertex_halfedge_to_border_halfedge() enable to set the border halfedge for a single vertex v, for all vertices on the boundary of the face of h, and for all vertices of the surface mesh, respectively.

Surface Mesh and the BGL API

The class Surface_mesh is a model of the concept IncidenceGraph defined in the Boost Graph Library. This enables to apply algorithms such as Dijkstra shortest path, or Kruskal minimum spanning tree directly on a surface mesh.

The types and free functions of the BGL API have each a similar type or member function, for example

BGL Surface_mesh Remark
boost::graph_traits<G>::vertex_descriptor Surface_mesh::Vertex_index
boost::graph_traits<G>::edge_descriptor Surface_mesh::Edge_index
vertices(const G& g) sm.vertices()
edges(const G& g) sm.edges()
vd = source(ed,g) vd = sm.source(ed)
na n = sm.number_of_vertices() counts non-deleted vertices and has no BGL equivalent
n = num_vertices(g) n = sm.number_of_vertices() + sm.number_of_removed_vertices() counts used and deleted vertices in order to have an upper bound on the largest vertex index used

It would be nicer to return the number of vertices without taking removed vertices into account, but this would interact badly with the underlying vertex/edge index mappings. The index mapping would no longer fall in the range [0,num_vertices(g)) which is assumed in many of the algorithms.

The class Surface_mesh is also a model of the concept MutableFaceGraph defined in the package CGAL and the Boost Graph Library. This and similar concepts like HalfedgeGraph refine the graph concepts of the BGL by introducing the notion of halfedges and faces, as well as cycles of halfedges around faces and around vertices. Again, there are similar types and functions, for example:

BGL Surface_mesh
boost::graph_traits<G>::halfedge_descriptor Surface_mesh::Halfedge_index
boost::graph_traits<G>::face_descriptor Surface_mesh::Face_index
halfedges(const G& g) sm.halfedges()
faces(const G& g) sm.faces()
hd = next(hd, g) hd = sm.next(hd)
hd = prev(hd, g) hd = sm.prev(hd)
hd = opposite(hd,g) hd = sm.opposite(hd)
hd = halfedge(vd,g) hd = sm.halfedge(vd)
etc.

The BGL API described in the package CGAL and the Boost Graph Library enables us to write geometric algorithms operating on surface meshes, that work for any model of FaceGraph, or MutableFaceGraph. That is surface mesh simplification, deformation, or segmentation algorithms work for Surface_mesh and Polyhedron_3.

BGL algorithms use property maps in order to associate information to vertices and edges. One important property is the index, an integer between 0 and num_vertices(g) for the vertices of a graph g. This allows algorithms to create a vector of the approriate size in order to store per vertex information. For example a Boolean for storing if a vertex has already been visited during a graph traversal.

The BGL way of retrieving the vertex index property map of a graph g is vipm = get(boost::vertex_index, g), and get(vipm, vd) in order then to retrieve the index for a vertex descriptor vd, and it is get(vertex_index, g, vd) to obtain the vertex index directly.

Example

The first example shows that we can apply Kruskal's minimum spanning tree algorithm directly on a surface mesh.


File Surface_mesh/sm_kruskal.cpp

#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
#include <iostream>
#include <fstream>
#include <list>
#include <boost/graph/kruskal_min_spanning_tree.hpp>
typedef Kernel::Point_3 Point;
typedef boost::graph_traits<Mesh>::vertex_descriptor vertex_descriptor;
typedef boost::graph_traits<Mesh>::vertex_iterator vertex_iterator;
typedef boost::graph_traits<Mesh>::edge_descriptor edge_descriptor;
void
kruskal(const Mesh& sm)
{
// We use the default edge weight which is the squared length of the edge
std::list<edge_descriptor> mst;
boost::kruskal_minimum_spanning_tree(sm,
std::back_inserter(mst));
std::cout << "#VRML V2.0 utf8\n"
"Shape {\n"
" appearance Appearance {\n"
" material Material { emissiveColor 1 0 0}}\n"
" geometry\n"
" IndexedLineSet {\n"
" coord Coordinate {\n"
" point [ \n";
vertex_iterator vb,ve;
for(boost::tie(vb, ve) = vertices(sm); vb!=ve; ++vb){
std::cout << " " << sm.point(*vb) << "\n";
}
std::cout << " ]\n"
" }\n"
" coordIndex [\n";
for(std::list<edge_descriptor>::iterator it = mst.begin(); it != mst.end(); ++it)
{
edge_descriptor e = *it ;
vertex_descriptor s = source(e,sm);
vertex_descriptor t = target(e,sm);
std::cout << " " << s << ", " << t << ", -1\n";
}
std::cout << "]\n"
" }#IndexedLineSet\n"
"}# Shape\n";
}
int main(int,char** argv) {
Mesh sm;
std::ifstream input(argv[1]);
input >> sm;
kruskal(sm);
return 0;
}

The second example shows how we can use property maps for algorithms such as Prim's minimum spanning tree. The algorithm internally also uses a vertex index property map calling get(boost::vertex_index_t,sm). For the class Surface_mesh this boils down to an identity function as vertices are indices.


File Surface_mesh/sm_bgl.cpp

#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
#include <iostream>
#include <fstream>
// workaround a bug in Boost-1.54
#include <CGAL/boost/graph/dijkstra_shortest_paths.h>
#include <boost/graph/prim_minimum_spanning_tree.hpp>
#include <boost/foreach.hpp>
typedef Kernel::Point_3 Point;
typedef boost::graph_traits<Mesh>::vertex_descriptor vertex_descriptor;
int main(int /* argc */, char* argv[])
{
Mesh sm;
std::ifstream in(argv[1]);
in >> sm;
Mesh::Property_map<vertex_descriptor,vertex_descriptor> predecessor;
predecessor = sm.add_property_map<vertex_descriptor,vertex_descriptor>("v:predecessor").first;
boost::prim_minimum_spanning_tree(sm, predecessor, boost::root_vertex(*vertices(sm).first));
std::cout << "#VRML V2.0 utf8\n"
"DirectionalLight {\n"
"direction 0 -1 0\n"
"}\n"
"Shape {\n"
" appearance Appearance {\n"
" material Material { emissiveColor 1 0 0}}\n"
" geometry\n"
" IndexedLineSet {\n"
" coord Coordinate {\n"
" point [ \n";
BOOST_FOREACH(vertex_descriptor vd, vertices(sm)){
std::cout << " " << sm.point(vd) << "\n";
}
std::cout << " ]\n"
" }\n"
" coordIndex [\n";
BOOST_FOREACH(vertex_descriptor vd, vertices(sm)){
if(predecessor[vd]!=vd){
std::cout << " " << std::size_t(vd) << ", " << std::size_t(predecessor[vd]) << ", -1\n";
}
}
std::cout << "]\n"
" }#IndexedLineSet\n"
"}# Shape\n";
sm.remove_property_map(predecessor);
return 0;
}

Memory Managment

Memory management is semi-automatic. Memory grows as more elements are added to the structure but does not shrink when elements are removed.

When you add elements and the capacity of the underlying vector is exhausted, the vector reallocates memory. As descriptors are basically indices, they refer to the same element after a reallocation.

When you remove an element it is only marked as removed. Internally it is put in a free list, and when you add elements to the surface mesh, they are taken from the free list in case it is not empty.

For all elements we offer a function to obtain the number of used elements, as well as the number of used and removed elements. For vertices the functions are Surface_mesh::number_of_vertices() and Surface_mesh::number_of_removed_vertices(), respectively. The first function is slightly different from the free function num_vertices(const G&) of the BGL package. As BGL style algorithms use the indices of elements to access data in temporary vectors of size num_vertices() this function must return a number larger than the largest index of the elements.

Iterators such as Surface_mesh::Vertex_iterator only enumerate elements that are not marked as deleted.

To really shrink the used memory, Surface_mesh::collect_garbage() must be called. Garbage collection also compacts the properties associated with the surface mesh.

Note however that by garbage collecting elements get new indices. In case you keep vertex descriptors they are most probably no longer refering to the right vertices.

Example


File Surface_mesh/sm_memory.cpp

#include <iostream>
#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
typedef Mesh::Vertex_index vertex_descriptor;
int main()
{
Mesh m;
Mesh::Vertex_index u;
for(unsigned int i=0; i < 5; ++i){
Mesh::Vertex_index v = m.add_vertex(K::Point_3(0,0,i+1));
if(i==2) u=v;
}
m.remove_vertex(u);
std::cout << "After insertion of 5 vertices and removal of the 3. vertex\n"
<< "# vertices / # vertices + # removed vertices = "
<< m.number_of_vertices()
<< " / " << m.number_of_vertices() + m.number_of_removed_vertices() << std::endl;
std::cout << "Iterate over vertices\n";
{
BOOST_FOREACH(vertex_descriptor vd, m.vertices()){
std::cout << m.point(vd) << std::endl;
}
}
// The status of being used or removed is stored in a property map
Mesh::Property_map<Mesh::Vertex_index,bool> removed
= m.property_map<Mesh::Vertex_index,bool>("v:removed").first;
std::cout << "\nIterate over vertices and deleted vertices\n"
<< "# vertices / # vertices + # removed vertices = "
<< m.number_of_vertices()
<< " / " << m.number_of_vertices() + m.number_of_removed_vertices() << std::endl;
{
unsigned int i = 0, end = m.number_of_vertices() + m.number_of_removed_vertices();
for( ; i < end; ++i) {
vertex_descriptor vh(i);
assert(m.is_removed(vh) == removed[vh]);
std::cout << m.point(vh) << ((m.is_removed(vh)) ? " R\n" : "\n");
}
}
m.collect_garbage();
std::cout << "\nAfter garbage collection\n"
<< "# vertices / # vertices + # removed vertices = "
<< m.number_of_vertices()
<< " / " << m.number_of_vertices() + m.number_of_removed_vertices() << std::endl;
{
unsigned int i = 0, end = m.number_of_vertices() + m.number_of_removed_vertices();
for( ; i < end; ++i) {
vertex_descriptor vh(i);
std::cout << m.point(vh) << ((m.is_removed(vh)) ? " R\n" : "\n");
}
}
return 0;
}

Draw a Surface Mesh

A surface mesh can be visualized by calling the CGAL::draw() function as shown in the following example. This function opens a new window showing the given surface mesh. The function is blocking, that is the program continues as soon as the user closes the window.


File Surface_mesh/draw_surface_mesh.cpp

#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/draw_surface_mesh.h>
#include <fstream>
typedef Kernel::Point_3 Point;
int main(int argc, char* argv[])
{
Mesh sm1;
std::ifstream in1((argc>1)?argv[1]:"data/triangle.off");
in1 >> sm1;
CGAL::draw(sm1);
return EXIT_SUCCESS;
}

This function requires CGAL_Qt5, and is only available if the flag CGAL_USE_BASIC_VIEWER is defined at compile time.

draw_surface_mesh.png
Figure 27.2 Result of the run of the draw_surface_mesh program. A window shows the surface mesh and allows to navigate through the 3D scene.

Implementation Details

As integer type for the indices we have chosen boost::uint32_t. On 64 bit operating systems they take only half the size of a pointer. They still allow to have meshes with 2 billion elements.

We use std::vector for storing properties. So by accessing the adress of the 0th element of a property map you can access the underlying raw array. This may be useful, for example for passing an array of points to OpenGL.

We use a freelist for removed elements. This mean when a vertex gets removed and later add_vertex is called, the memory of the removed element is reused. This especially means that the n'th inserted element has not necessarily the index n-1, and when iterating over elements they will not be enumerated in the insertion order.

Implementation History

This package is derived from an early version of Daniel Sieger and Mario Botsch package Surface_mesh [1], which is inspired from the design of OpenMesh and the CGAL package 3D Polyhedral Surface.

Philipp Moeller and Andreas Fabri worked on the code so that iterators fulfill the requirements of the STL iterator concepts, and changed the API so that it becomes a model of the concepts MutableFaceGraph and FaceListGraph of the package CGAL and the Boost Graph Library.