Class

CGAL::Delaunay_triangulation_3<DelaunayTriangulationTraits_3,TriangulationDataStructure_3,LocationPolicy>

Definition

The class Delaunay_triangulation_3 represents a three-dimensional Delaunay triangulation.

#include <CGAL/Delaunay_triangulation_3.h>

Parameters

The first template argument must be a model of the DelaunayTriangulationTraits_3 concept.

The second template argument must be a model of the TriangulationDataStructure_3 concept. It has the default value Triangulation_data_structure_3<Triangulation_vertex_base_3<DelaunayTriangulationTraits_3>, Triangulation_cell_base_3<DelaunayTriangulationTraits_3> >.

The third template argument is a tag which must be a Location_policy<Tag> : either Fast_location or Compact_location. Fast_location offers faster (O(log n) time) point location, which can be beneficial when performing point locations or random point insertions (with no good location hint) in large data sets. It is currently implemented using an additional triangulation hierarchy data structure [Dev02]. The default is Compact_location, which saves memory (3-5%) by avoiding the need for this separate data structure, and point location is then performed roughly in O(n1/3) time. Note that this argument can also come in second position, which can be useful when the default value for the TriangulationDataStructure_3 parameter is satisfactory (this is using so-called deduced parameters). Note that this argument replaces the functionality provided before Cgal 3.6 by Triangulation_hierarchy_3. An example of use can be found in the user manual 39.5.5.

Inherits From

Triangulation_3<DelaunayTriangulationTraits_3,Delaunay_triangulation_3<DelaunayTriangulationTraits_3,TriangulationDataStructure_3,LocationPolicy>::Triangulation_data_structure >

Types

typedef LocationPolicy Location_policy;

In addition to those inherited, the following types are defined, for use by the construction of the Voronoi diagram:

typedef DelaunayTriangulationTraits_3::Line_3
Line;
typedef DelaunayTriangulationTraits_3::Ray_3
Ray;
typedef DelaunayTriangulationTraits_3::Plane_3
Plane;
typedef DelaunayTriangulationTraits_3::Object_3
Object;

Creation

Delaunay_triangulation_3<DelaunayTriangulationTraits_3,TriangulationDataStructure_3,LocationPolicy> dt ( DelaunayTriangulationTraits_3 traits = DelaunayTriangulationTraits_3());
Creates an empty Delaunay triangulation, possibly specifying a traits class traits.


Delaunay_triangulation_3<DelaunayTriangulationTraits_3,TriangulationDataStructure_3,LocationPolicy> dt ( Delaunay_triangulation_3 dt1);
Copy constructor.


template < class InputIterator >
Delaunay_triangulation_3<DelaunayTriangulationTraits_3,TriangulationDataStructure_3,LocationPolicy> dt ( InputIterator first,
InputIterator last,
DelaunayTriangulationTraits_3 traits = DelaunayTriangulationTraits_3());
Equivalent to constructing an empty triangulation with the optional traits class argument and calling insert(first,last).

Operations

Insertion

The following methods overload the corresponding methods of triangulations to ensure the empty sphere property of Delaunay triangulations.

In the degenerate case when there are co-spherical points, the Delaunay triangulation is known not to be uniquely defined. In this case, Cgal chooses a particular Delaunay triangulation using a symbolic perturbation scheme [DT03].

Vertex_handle dt.insert ( Point p, Cell_handle start = Cell_handle())
Inserts point p in the triangulation and returns the corresponding vertex. Similar to the insertion in a triangulation, but ensures in addition the empty sphere property of all the created faces. The optional argument start is used as a starting place for the search.

Vertex_handle dt.insert ( Point p, Vertex_handle hint)
Same as above but uses hint as a starting place for the search.

Vertex_handle dt.insert ( Point p, Locate_type lt, Cell_handle loc, int li, int lj)
Inserts point p in the triangulation and returns the corresponding vertex. Similar to the above insert() function, but takes as additional parameter the return values of a previous location query. See description of Triangulation_3::locate().

The following method allows one to insert several points. It returns the number of inserted points.

template < class PointInputIterator >
std::ptrdiff_t dt.insert ( PointInputIterator first, PointInputIterator last)
Inserts the points in the iterator range [.first, last.). Returns the number of inserted points. Note that this function is not guaranteed to insert the points following the order of PointInputIterator, as spatial_sort is used to improve efficiency.
Precondition: The value_type of first and last is Point.

template < class PointWithInfoInputIterator >
std::ptrdiff_t dt.insert ( PointWithInfoInputIterator first, PointWithInfoInputIterator last)
Inserts the points in the iterator range [.first, last.). Returns the number of inserted points. Note that this function is not guaranteed to insert the points following the order of PointWithInfoInputIterator, as spatial_sort is used to improve efficiency. Given a pair (p,i), the vertex v storing p also stores i, that is v.point() == p and v.info() == i. If several pairs have the same point, only one vertex is created, and one of the objects of type Vertex::Info will be stored in the vertex.
Precondition: Vertex must be model of the concept TriangulationVertexBaseWithInfo_3. The value_type of first and last is std::pair<Point,Vertex::Info>.

Displacement

Vertex_handle dt.move_if_no_collision ( Vertex_handle v, Point p)
if there is not already another vertex placed on p, the triangulation is modified such that the new position of vertex v is p, and v is returned. Otherwise, the triangulation is not modified and the vertex at point p is returned.
Precondition: Vertex v must be finite.

Vertex_handle dt.move ( Vertex_handle v, Point p)
same as above if there is no collision. Otherwise, v is deleted and the vertex placed on p is returned.
Precondition: Vertex v must be finite.

Removal

When a vertex v is removed from a triangulation, all the cells incident to v must be removed, and the polyhedral region consisting of all the tetrahedra that are incident to v must be re-triangulated. So, the problem reduces to triangulating a polyhedral region, while preserving its boundary, or to compute a constrained triangulation. This is known to be sometimes impossible: the Schönhardt polyhedron cannot be triangulated [She98a].

However, when dealing with Delaunay triangulations, the case of such polyhedra that cannot be re-triangulated cannot happen, so Cgal proposes a vertex removal.

void dt.remove ( Vertex_handle v) Removes the vertex v from the triangulation.
Precondition: v is a finite vertex of the triangulation.

template < typename InputIterator >
int dt.remove ( InputIterator first, InputIterator beyond)
Removes the vertices specified by the iterator range [first, beyond) of value type Vertex_handle. remove() is called over each element of the range. The number of vertices removed is returned.
Precondition: (i) all vertices of the range are finite vertices of the triangulation; and (ii) no vertices are repeated in the range.

template < typename InputIterator >
int dt.remove_cluster ( InputIterator first, InputIterator beyond)
This function has exactly the same result and the same preconditions as remove(first, beyond). The difference is in the implementation and efficiency. This version does not re-triangulate the hole after each point removal but only after removing all vertices. This is more efficient if (and only if) the removed points are organized in a small number of connected components of the Delaunay triangulation.

If, due to some point removals, the size of the Delaunay triangulation decreases drastically, it might be interesting to defragment the CGAL::Compact_container (used by the Triangulation_data_structure_3).

Queries

Bounded_side dt.side_of_sphere ( Cell_handle c, Point p) const
Returns a value indicating on which side of the circumscribed sphere of c the point p lies. More precisely, it returns:
- ON_BOUNDED_SIDE if p is inside the sphere. For an infinite cell this means that p lies strictly either in the half space limited by its finite facet and not containing any other point of the triangulation, or in the interior of the disk circumscribing the finite facet.
- ON_BOUNDARY if p on the boundary of the sphere. For an infinite cell this means that p lies on the circle circumscribing the finite facet.
- ON_UNBOUNDED_SIDE if p lies outside the sphere. For an infinite cell this means that p does not satisfy either of the two previous conditions.
Precondition: dt.dimension() =3.

Bounded_side dt.side_of_circle ( Facet f, Point p) const
Returns a value indicating on which side of the circumscribed circle of f the point p lies. More precisely, it returns:
- in dimension 3:
- For a finite facet, ON_BOUNDARY if p lies on the circle, ON_UNBOUNDED_SIDE when it lies in the exterior of the disk, ON_BOUNDED_SIDE when it lies in its interior.
- For an infinite facet, it considers the plane defined by the finite facet of the same cell, and does the same as in dimension 2 in this plane.
- in dimension 2:
- For a finite facet, ON_BOUNDARY if p lies on the circle, ON_UNBOUNDED_SIDE when it lies in the exterior of the disk, ON_BOUNDED_SIDE when it lies in its interior.
- For an infinite facet, ON_BOUNDARY if the point lies on the finite edge of f (endpoints included), ON_BOUNDED_SIDE for a point in the open half plane defined by f and not containing any other point of the triangulation, ON_UNBOUNDED_SIDE elsewhere.
Precondition: dt.dimension() 2 and in dimension 3, p is coplanar with f.

Bounded_side dt.side_of_circle ( Cell_handle c, int i, Point p)
Same as the previous method for facet i of cell c.

Vertex_handle dt.nearest_vertex ( Point p, Cell_handle c = Cell_handle())
Returns any nearest vertex to the point p, or the default constructed handle if the triangulation is empty. The optional argument c is a hint specifying where to start the search.
Precondition: c is a cell of dt.

Vertex_handle dt.nearest_vertex_in_cell ( Point p, Cell_handle c)
Returns the vertex of the cell c that is nearest to p.

A point p is said to be in conflict with a cell c in dimension 3 (resp. a facet f in dimension 2) iff dt.side_of_sphere(c, p) (resp. dt.side_of_circle(f, p)) returns ON_BOUNDED_SIDE. The set of cells (resp. facets in dimension 2) which are in conflict with p is connected, and it forms a hole.

template <class OutputIteratorBoundaryFacets, class OutputIteratorCells>
std::pair<OutputIteratorBoundaryFacets, OutputIteratorCells>
dt.find_conflicts ( Point p,
Cell_handle c,
OutputIteratorBoundaryFacets bfit,
OutputIteratorCells cit)
Computes the conflict hole induced by p. The starting cell (resp. facet) c must be in conflict. Then this function returns respectively in the output iterators:
- cit: the cells (resp. facets) in conflict.
- bfit: the facets (resp. edges) on the boundary, that is, the facets (resp. edges) (t, i) where the cell (resp. facet) t is in conflict, but t->neighbor(i) is not.
This function can be used in conjunction with insert_in_hole() in order to decide the insertion of a point after seeing which elements of the triangulation are affected. Returns the pair composed of the resulting output iterators.
Precondition: dt.dimension() 2, and c is in conflict with p.

template <class OutputIteratorBoundaryFacets, class OutputIteratorCells, class OutputIteratorInternalFacets>
Triple<OutputIteratorBoundaryFacets, OutputIteratorCells, OutputIteratorInternalFacets>
dt.find_conflicts ( Point p,
Cell_handle c,
OutputIteratorBoundaryFacets bfit,
OutputIteratorCells cit,
OutputIteratorInternalFacets ifit)
Same as the other find_conflicts() function, except that it also computes the internal facets, i.e. the facets common to two cells which are in conflict with p. Then this function returns respectively in the output iterators:
- cit: the cells (resp. facets) in conflict.
- bfit: the facets (resp. edges) on the boundary, that is, the facets (resp. edges) (t, i) where the cell (resp. facet) t is in conflict, but t->neighbor(i) is not.
- ifit: the facets (resp. edges) inside the hole, that is, delimiting two cells (resp facets) in conflict.
Returns the Triple composed of the resulting output iterators.
Precondition: dt.dimension() 2, and c is in conflict with p.

template <class OutputIterator>
OutputIterator dt.vertices_in_conflict ( Point p, Cell_handle c, OutputIterator res)
This function is renamed vertices_on_conflict_zone_boundary since CGAL-3.8.
template <class OutputIterator>
OutputIterator dt.vertices_on_conflict_zone_boundary ( Point p, Cell_handle c, OutputIterator res)
Similar to find_conflicts(), but reports the vertices which are on the boundary of the conflict hole of p, in the output iterator res. Returns the resulting output iterator.
Precondition: dt.dimension() 2, and c is in conflict with p.

A face (cell, facet or edge) is said to be a Gabriel face iff its smallest circumscribing sphere do not enclose any vertex of the triangulation. Any Gabriel face belongs to the Delaunay triangulation, but the reciprocal is not true. The following member functions test the Gabriel property of Delaunay faces.

bool dt.is_Gabriel ( Cell_handle c, int i)
bool dt.is_Gabriel ( Cell_handle c, int i, int j)
bool dt.is_Gabriel ( Facet f)
bool dt.is_Gabriel ( Edge e)

Voronoi diagram

Cgal offers several functionalities to display the Voronoi diagram of a set of points in 3D.

Note that the user should use a kernel with exact constructions in order to guarantee the computation of the Voronoi diagram (as opposed to computing the triangulation only, which requires only exact predicates).

Point dt.dual ( Cell_handle c) const Returns the circumcenter of the four vertices of c.
Precondition: dt.dimension()=3 and c is not infinite.

Object dt.dual ( Facet f) const Returns the dual of facet f, which is
in dimension 3: either a segment, if the two cells incident to f are finite, or a ray, if one of them is infinite;
in dimension 2: a point.
Precondition: dt.dimension() 2 and f is not infinite.

Object dt.dual ( Cell_handle c, int i) const
same as the previous method for facet (c,i).

Line dt.dual_support ( Cell_handle c, int i) const
returns the line supporting the dual of the facet.
Precondition: dt.dimension() 2 and f is not infinite.

template <class Stream>
Stream & dt.draw_dual ( Stream & os) Sends the set of duals to all the facets of dt into os.

Checking

bool dt.is_valid ( bool verbose = false) const
Checks the combinatorial validity of the triangulation and the validity of its geometric embedding (see Section 39.1). Also checks that all the circumscribing spheres (resp. circles in dimension 2) of cells (resp. facets in dimension 2) are empty.
When verbose is set to true, messages describing the first invalidity encountered are printed.

bool dt.is_valid ( Cell_handle c, bool verbose = false) const
Checks the combinatorial and geometric validity of the cell (see Section 39.1). Also checks that the circumscribing sphere (resp. circle in dimension 2) of cells (resp. facet in dimension 2) is empty.
When verbose is set to true, messages are printed to give a precise indication of the kind of invalidity encountered.

These methods are mainly a debugging help for the users of advanced features.

See Also

CGAL::Regular_triangulation_3