The class Delaunay_triangulation_2<Traits,Tds> is designed to represent the Delaunay triangulation of a set of points in a plane. A Delaunay triangulation of a set of points is a triangulation of the sets of points that fulfills the following empty circle property (also called Delaunay property): the circumscribing circle of any facet of the triangulation contains no point of the set in its interior. For a point set with no case of co-circularity of more than three points, the Delaunay triangulation is unique, it is the dual of the Voronoi diagram of the points.
#include <CGAL/Delaunay_triangulation_2.h>
The geometric traits Traits is to be instantiated with a model of DelaunayTriangulationTraits_2. The concept DelaunayTriangulationTraits_2 refines the concept TriangulationTraits_2, providing a predicate type to check the empty circle property.
Changing this predicate type allows to build Delaunay triangulations for different metrics such that L1 or L∞ or any metric defined by a convex object. However, the user of an exotic metric must be careful that the constructed triangulation has to be a triangulation of the convex hull which means that convex hull edges have to be Delaunay edges. This is granted for any smooth convex metric (like L2) and can be ensured for other metrics (like L∞) by the addition to the point set of well chosen sentinel points The concept of DelaunayTriangulationTraits_2 is described .
When dealing with a large triangulations, the user is advised to encapsulate the Delaunay triangulation class into a triangulation hierarchy, which means to use the class Triangulation_hierarchy_2<Tr> with the template parameter instantiated with Delaunay_triangulation_2<Traits,Tds> . The triangulation hierarchy will then offer the same functionalities but be much more for efficient for locations and insertions.
Delaunay_triangulation_2<Traits,Tds> dt ( Traits gt = Traits()); | |
default constructor.
| |
Delaunay_triangulation_2<Traits,Tds> dt ( tr); | |
copy constructor. All the vertices and faces are duplicated.
|
The following insertion and removal functions overwrite the functions inherited from the class Triangulation_2<Traits,Tds> to maintain the Delaunay property.
In the degenerate case when there are co-circular 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, Face_handle f=Face_handle()) | |||
inserts point p. If point p coincides with an already existing vertex, this vertex is returned and the triangulation is not updated. Optional parameter f is used to initialize the location of p. | ||||
Vertex_handle | dt.insert ( Point p, Locate_type& lt, Face_handle loc, int li) | |||
inserts a point p, the location of which is supposed to be given by (lt,loc,li), see the description of member function locate in class Triangulation_2<Traits,Tds>. | ||||
Vertex_handle | dt.push_back ( Point p) | equivalent to insert(p). | ||
template < class PointInputIterator > | ||||
std::ptrdiff_t | dt.insert ( PointInputIterator first, PointInputIterator last) | |||
inserts the points in the 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.
| ||||
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.
| ||||
void | dt.remove ( Vertex_handle v) | removes the vertex from the triangulation. |
Note that the other modifier functions of Triangulation_2<Traits,Tds> are not overwritten. Thus a call to insert_in_face insert_in_edge, insert_outside_convex_hull, insert_outside_affine_hull or flip on a valid Delaunay triangulation might lead to a triangulation which is no longer a Delaunay triangulation.
Vertex_handle | dt.nearest_vertex ( Point p, Face_handle f=Face_handle()) | |||
returns any nearest vertex of p. The implemented function begins with a location step and f may be used to initialize the location. | ||||
template <class OutputItFaces, class OutputItBoundaryEdges> | ||||
std::pair<OutputItFaces,OutputItBoundaryEdges> | ||||
| ||||
OutputItFaces is an output iterator with Face_handle as value type.
OutputItBoundaryEdges stands for an output iterator with Edge as value type.
This members function outputs in the container pointed to by fit
the faces which are in conflict with point p
i. e. the faces whose circumcircle contains p.
It outputs in the container pointed to by eit the
the boundary of the zone in conflict with p.
The boundary edges
of the conflict zone are output in counter-clockwise order
and each edge is described through its incident face
which is not in conflict with p.
The function returns in a std::pair the resulting output iterators.
| ||||
template <class OutputItFaces> | ||||
OutputItFaces | dt.get_conflicts ( Point p, OutputItFaces fit, Face_handle start) const | |||
same as above except that only the faces in conflict with p
are output. The function returns the resulting output iterator.
| ||||
template <class OutputItBoundaryEdges> | ||||
OutputItBoundaryEdges | dt.get_boundary_of_conflicts ( Point p, OutputItBoundaryEdges eit, Face_handle start) const | |||
OutputItBoundaryEdges stands for an output iterator with Edge as value type. This function outputs in the container pointed to by eit, the boundary of the zone in conflict with p. The boundary edges of the conflict zone are output in counterclockwise order and each edge is described through the incident face which is not in conflict with p. The function returns the resulting output iterator. |
Point | dt.dual ( Face_handle f) const |
Returns the center of the circle circumscribed to face f.
| ||
Object | dt.dual ( Edge e) const | returns a segment, a ray or a line supported by the bisector of the endpoints of e. If faces incident to e are both finite, a segment whose endpoints are the duals of each incident face is returned. If only one incident face is finite, a ray whose endpoint is the dual of the finite incident face is returned. Otherwise both incident faces are infinite and the bisector line is returned. | ||
Object | dt.dual ( Edge_circulator ec) const | |||
Idem | ||||
Object | dt.dual ( Edge_iterator ei) const | |||
Idem | ||||
template < class Stream> | ||||
Stream& | dt.draw_dual ( Stream & ps) | output the dual Voronoi diagram to stream ps. |
Oriented_side | dt.side_of_oriented_circle ( Face_handle f, Point p) const | |
Returns the side of p with respect to the circle circumscribing the triangle associated with f |
The checking function is_valid() is also overwritten to additionally test the empty circle property.
bool | dt.is_valid ( bool verbose = false, int level = 0) const | |
tests the validity of the triangulation as a Triangulation_2 and additionally tests the Delaunay property. This method is mainly useful for debugging Delaunay triangulation algorithms designed by the user. |
CGAL::Triangulation_2<Traits,Tds>,
TriangulationDataStructure_2,
DelaunayTriangulationTraits_2,
Triangulation_hierarchy_2<Tr>.
Insertion is implemented by inserting in the triangulation, then performing a sequence of Delaunay flips. The number of flips is O(d) if the new vertex is of degree d in the new triangulation. For points distributed uniformly at random, insertion takes time O(1) on average.
Removal calls the removal in the triangulation and then re-triangulates the hole in such a way that the Delaunay criterion is satisfied. Removal of a vertex of degree d takes time O(d2). The degree d is O(1) for a random vertex in the triangulation.
After a point location step, the nearest neighbor is found in time O(n) in the worst case, but in time O(1) for vertices distributed uniformly at random and any query point.