40 #ifndef GEOGRAM_DELAUNAY_DELAUNAY_2D
41 #define GEOGRAM_DELAUNAY_DELAUNAY_2D
127 index_t nb_vertices,
const double* vertices
145 return has_empty_cells_;
156 abort_if_empty_cell_ = x;
204 const double* p,
index_t hint = NO_TRIANGLE,
205 bool thread_safe =
false,
206 Sign* orient =
nullptr
321 return cell_to_v_store_.size() / 3;
364 static const index_t END_OF_LIST = ~(NOT_IN_LIST_BIT);
386 return (cell_next_[t] & NOT_IN_LIST_BIT) == 0;
402 return cell_next_[t];
418 if(last == END_OF_LIST) {
421 cell_next_[t] = END_OF_LIST;
423 cell_next_[t] = first;
440 cell_next_[t] = NOT_IN_LIST;
463 cell_to_v_store_[3 * t] >= 0 &&
464 cell_to_v_store_[3 * t + 1] >= 0 &&
465 cell_to_v_store_[3 * t + 2] >= 0 ;
480 return !triangle_is_free(t) && triangle_is_finite(t);
494 !triangle_is_free(t) && (
495 cell_to_v_store_[3 * t] == VERTEX_AT_INFINITY ||
496 cell_to_v_store_[3 * t + 1] == VERTEX_AT_INFINITY ||
497 cell_to_v_store_[3 * t + 2] == VERTEX_AT_INFINITY
512 return triangle_is_in_list(t);
524 if(first_free_ == END_OF_LIST) {
525 cell_to_v_store_.resize(cell_to_v_store_.size() + 3, -1);
526 cell_to_cell_store_.resize(cell_to_cell_store_.size() + 3, -1);
530 cell_next_.push_back(
index_t(NOT_IN_LIST));
531 result = max_t() - 1;
533 result = first_free_;
534 first_free_ = triangle_next(first_free_);
535 remove_triangle_from_list(result);
538 cell_to_cell_store_[3 * result] = -1;
539 cell_to_cell_store_[3 * result + 1] = -1;
540 cell_to_cell_store_[3 * result + 2] = -1;
560 index_t result = new_triangle();
561 cell_to_v_store_[3 * result] = v1;
562 cell_to_v_store_[3 * result + 1] = v2;
563 cell_to_v_store_[3 * result + 2] = v3;
575 cur_stamp_ = (stamp | NOT_IN_LIST_BIT);
594 return cell_next_[t] == cur_stamp_;
611 cell_next_[t] = cur_stamp_;
635 return index_t(triangle_edge_vertex_[e][v]);
648 return cell_to_v_store_[3 * t + lv];
677 return index_t(cell_to_v_store_[3 * t + lv]);
689 cell_to_v_store_[3 * t + lv] = v;
768 cell_to_v_store_[3 * t] = v0;
769 cell_to_v_store_[3 * t + 1] = v1;
770 cell_to_v_store_[3 * t + 2] = v2;
797 pv[i] = (v == -1) ?
nullptr : vertex_ptr(
index_t(v));
802 for(
index_t le = 0; le < 3; ++le) {
804 if(pv[le] ==
nullptr) {
830 if(triangle_is_in_list(t2)) {
835 if(triangle_is_marked(t2)) {
839 return triangle_is_conflict(t2, p);
848 double h0 = heights_[finite_triangle_vertex(t, 0)];
849 double h1 = heights_[finite_triangle_vertex(t, 1)];
850 double h2 = heights_[finite_triangle_vertex(t, 2)];
852 (p - vertex_ptr(0)) /
int(vertex_stride_)
854 double h = heights_[pindex];
856 pv[0],pv[1],pv[2],p,h0,h1,h2,h
914 index_t first,
const std::string& list_name
945 bool verbose_debug_mode_;
950 bool benchmark_mode_;
960 static char triangle_edge_vertex_[3][2];
965 std::stack<index_t> S_;
970 bool has_empty_cells_;
976 bool abort_if_empty_cell_;
#define geo_debug_assert(x)
Verifies that a condition is met.
Implementation of Delaunay in 2d.
void find_conflict_zone_iterative(const double *p, index_t t, index_t &t_bndry, index_t &e_bndry, index_t &first, index_t &last)
This function is used to implement find_conflict_zone.
bool triangle_is_finite(index_t t) const
Tests whether a given triangle is a finite one.
static index_t triangle_edge_vertex(index_t e, index_t v)
Returns the local index of a vertex by edge and by local vertex index in the edge.
bool triangle_is_real(index_t t) const
Tests whether a triangle is a real one.
void set_triangle_vertex(index_t t, index_t lv, signed_index_t v)
Sets a triangle-to-vertex adjacency.
void set_vertices(index_t nb_vertices, const double *vertices) override
Sets the vertices of this Delaunay, and recomputes the cells.
bool triangle_is_in_list(index_t t) const
Tests whether a triangle belongs to a linked list.
index_t triangle_next(index_t t) const
Gets the index of a successor of a triangle.
signed_index_t triangle_vertex(index_t t, index_t lv) const
Gets the index of a vertex of a triangle.
void set_tet(index_t t, signed_index_t v0, signed_index_t v1, signed_index_t v2, index_t a0, index_t a1, index_t a2)
Sets the vertices and adjacent triangles of a triangle.
void set_triangle_mark_stamp(index_t stamp)
Generates a unique stamp for marking tets.
index_t new_triangle(signed_index_t v1, signed_index_t v2, signed_index_t v3)
Creates a new triangle.
signed_index_t triangle_adjacent(index_t t, index_t le) const
Gets the index of a triangle adjacent to another one.
void show_list(index_t first, const std::string &list_name) const
For debugging purposes, displays a triangle.
index_t locate_inexact(const double *p, index_t hint, index_t max_iter) const
Finds the triangle that (approximately) contains a point using inexact predicates.
index_t max_t() const
Maximum valid index for a triangle.
Delaunay2d(coord_index_t dimension=2)
Constructs a new Delaunay2d.
~Delaunay2d() override
Delaunay2d destructor.
void show_triangle(index_t t) const
For debugging purposes, displays a triangle.
void set_triangle_adjacent(index_t t1, index_t le1, index_t t2)
Sets a triangle-to-triangle adjacency.
void check_geometry(bool verbose=false) const
For debugging purposes, test some geometrical properties.
bool triangle_is_virtual(index_t t) const
Tests whether a triangle is a virtual one.
void abort_if_empty_cell(bool x)
Specifies behavior if an empty cell is detected.
index_t stellate_conflict_zone(index_t v, index_t t_bndry, index_t e_bndry)
Creates a star of triangles filling the conflict zone.
bool create_first_triangle(index_t &iv0, index_t &iv1, index_t &iv2)
Finds in the pointset a set of three non-colinear points.
index_t find_triangle_vertex(index_t t, signed_index_t v) const
Finds the index of the vertex in a triangle.
index_t insert(index_t v, index_t hint=NO_TRIANGLE)
Inserts a point in the triangulation.
void find_conflict_zone(index_t v, index_t t, const Sign *orient, index_t &t_bndry, index_t &e_bndry, index_t &first, index_t &last)
Determines the list of triangles in conflict with a given point.
index_t finite_triangle_vertex(index_t t, index_t lv) const
Gets the index of a vertex of a triangle.
index_t find_triangle_adjacent(index_t t1, index_t t2_in) const
Finds the index of the edge accros which t1 is adjacent to t2_in.
void remove_triangle_from_list(index_t t)
Removes a triangle from the linked list it belongs to.
void mark_triangle(index_t t)
Marks a triangle.
index_t locate(const double *p, index_t hint=NO_TRIANGLE, bool thread_safe=false, Sign *orient=nullptr) const
Finds the triangle that contains a point.
bool has_empty_cells() const
Tests whether the Laguerre diagram has empty cells.
bool triangle_is_marked(index_t t) const
Tests whether a triangle is marked.
void add_triangle_to_list(index_t t, index_t &first, index_t &last)
Adds a triangle to a linked list.
bool triangle_is_free(index_t t) const
Tests whether a triangle is in the free list.
static index_t find_3(const signed_index_t *T, signed_index_t v)
Finds the index of an integer in an array of three integers.
index_t nearest_vertex(const double *p) const override
Computes the nearest vertex from a query point.
void check_combinatorics(bool verbose=false) const
For debugging purposes, tests some combinatorial properties.
index_t new_triangle()
Creates a new triangle.
void show_triangle_adjacent(index_t t, index_t le) const
For debugging purposes, displays a triangle adjacency.
bool triangle_is_conflict(index_t t, const double *p) const
Tests whether a given triangle is in conflict with a given 3d point.
Abstract interface for Delaunay triangulation in Nd.
Regular Delaunay triangulation of weighted points.
~RegularWeightedDelaunay2d() override
RegularWeightedDelaunay2d destructor.
RegularWeightedDelaunay2d(coord_index_t dimension=3)
Constructs a new Regular Delaunay2d triangulation.
Abstract interface for Delaunay.
Common include file, providing basic definitions. Should be included before anything else by all head...
Geometric functions in 2d and 3d.
Sign in_circle_2d_SOS(const double *p0, const double *p1, const double *p2, const double *p3)
Tests whether a 2d point is inside the circumscribed circle of a 3d triangle.
Sign orient_2dlifted_SOS(const double *p0, const double *p1, const double *p2, const double *p3, double h0, double h1, double h2, double h3)
Computes the 3d orientation test with lifted points.
Sign orient_2d(const vec2HE &p0, const vec2HE &p1, const vec2HE &p2)
Computes the orientation predicate in 2d.
Global Vorpaline namespace.
geo_signed_index_t signed_index_t
The type for storing and manipulating indices differences.
geo_index_t index_t
The type for storing and manipulating indices.
Sign
Integer constants that represent the sign of a value.
geo_coord_index_t coord_index_t
The type for storing coordinate indices, and iterating on the coordinates of a point.
Filtered exact predicates for restricted Voronoi diagrams.