40#ifndef GEOGRAM_DELAUNAY_CDT_2D
41#define GEOGRAM_DELAUNAY_CDT_2D
60 struct CDT2d_ConstraintWalker;
111 delaunay_ = delaunay;
156 return find_3(T_.data()+3*t, v);
169 return Tadj_[3*t+le];
182 return find_3(Tadj_.data()+3*t1, t2);
200 virtual void save(
const std::string& filename)
const = 0;
274 return ((Tflags_[t] & (1u << flag)) != 0);
291 T_MARKED_FLAG = DLIST_NB
327 cdt_(cdt), list_id_(list_id),
341 cdt_(cdt), list_id_(
index_t(-1)),
359 return (list_id_ !=
index_t(-1));
371 (back_==index_t(-1))==(front_==index_t(-1))
373 return (back_==index_t(-1));
376 bool contains(index_t t)
const {
378 return cdt_.Tflag_is_set(t, list_id_);
391 index_t next(index_t t)
const {
394 return cdt_.Tnext_[t];
397 index_t prev(index_t t)
const {
400 return cdt_.Tprev_[t];
404 for(index_t t=front_; t!=
index_t(-1); t = cdt_.Tnext_[t]) {
405 cdt_.Treset_flag(t,list_id_);
414 for(index_t t=front(); t!=
index_t(-1); t = next(t)) {
420 void push_back(index_t t) {
423 cdt_.Tset_flag(t,list_id_);
431 cdt_.Tnext_[back_] = t;
432 cdt_.Tprev_[t] = back_;
441 back_ = cdt_.Tprev_[back_];
446 cdt_.Tnext_[back_] =
index_t(-1);
449 cdt_.Treset_flag(t,list_id_);
453 void push_front(index_t t) {
456 cdt_.Tset_flag(t,list_id_);
464 cdt_.Tprev_[front_] = t;
465 cdt_.Tnext_[t] = front_;
474 front_ = cdt_.Tnext_[front_];
479 cdt_.Tprev_[front_] =
index_t(-1);
482 cdt_.Treset_flag(t,list_id_);
486 void remove(index_t t) {
490 }
else if(t == back_) {
494 index_t t_prev = cdt_.Tprev_[t];
495 index_t t_next = cdt_.Tnext_[t];
496 cdt_.Tprev_[t_next] = t_prev;
497 cdt_.Tnext_[t_prev] = t_next;
498 cdt_.Treset_flag(t,list_id_);
502 void show(std::ostream& out = std::cerr)
const {
514 out <<
"<uninitialized list>";
517 out <<
"<unknown list id:" << list_id_ <<
">";
521 for(index_t t=front(); t!=
index_t(-1); t = next(t)) {
551 insert_vertex_in_edge(v,t,le,S);
660 Tecnstr_[3*t] = e1cnstr;
661 Tecnstr_[3*t+1] = e2cnstr;
662 Tecnstr_[3*t+2] = e3cnstr;
685 Tadj_[i], Tadj_[j], Tadj_[k],
686 Tecnstr_[i], Tecnstr_[j], Tecnstr_[k]
740 index_t le2 = Tadj_find(t2,prev_t2_adj_e2);
742 Tset_edge_cnstr(t1,le1,Tedge_cnstr(t2,le2));
754 Tecnstr_.resize(nc,
index_t(-1));
755 Tflags_.resize(t+1,0);
756 Tnext_.resize(t+1,
index_t(-1));
757 Tprev_.resize(t+1,
index_t(-1));
771 return Tecnstr_[3*t+le];
785 Tecnstr_[3*t+le] = cnstr_id;
799 Tset_edge_cnstr(t, le, cnstr_id);
803 Tset_edge_cnstr(t2,le2,cnstr_id);
815 return (Tedge_cnstr(t,le) !=
index_t(-1));
837 t = Tadj(t, (lv+1)%3);
838 }
while(t != vT(v) && t !=
index_t(-1));
850 t = Tadj(t, (lv+2)%3);
856 t = Tadj(t, (lv+2)%3);
1004 for(
index_t t=0; t<nT(); ++t) {
1016 check_combinatorics();
1025 if(delaunay_ && exact_incircle_) {
1026 for(
index_t t=0; t<nT(); ++t) {
1027 for(
index_t le=0; le<3; ++le) {
1052 check_combinatorics();
1063 debug_check_combinatorics();
1064 debug_check_geometry();
1078 if(orient2d(u1,u2,v1)*orient2d(u1,u2,v2) > 0) {
1081 return (orient2d(v1,v2,u1)*orient2d(v1,v2,u2) < 0);
1097 index_t u1 = Tv(t,(le + 1)%3);
1098 index_t u2 = Tv(t,(le + 2)%3);
1099 return segment_segment_intersect(u1,u2,v1,v2);
1113 typedef std::pair<index_t, index_t> Edge;
1124 for_each_T_around_v(
1126 if(Tv(t, (lv+1)%3) == v2) {
1127 if(Tv(t, (lv+2)%3) != v1) {
1132 }
else if(Tv(t, (lv+1)%3) == v1) {
1133 if(Tv(t, (lv+2)%3) != v2) {
1144 (Tv(result,1) == v1 && Tv(result,2) == v2) ||
1145 (Tv(result,1) == v2 && Tv(result,2) == v1)
1275 double x1,
double y1,
double x2,
double y2
1277 create_enclosing_quad(
1294 debug_check_consistency();
1295 point_.push_back(p);
1296 index_t v = CDTBase2d::insert(point_.size()-1, hint);
1299 if(point_.size() > nv()) {
1302 debug_check_consistency();
1331 index_t nb_points,
const double* points,
1333 bool remove_unreferenced_vertices =
false
1339 void save(
const std::string& filename)
const override;
#define geo_assert(x)
Verifies that a condition is met.
#define geo_debug_assert(x)
Verifies that a condition is met.
Common include file, providing basic definitions. Should be included before anything else by all head...
Constrained Delaunay triangulation.
index_t create_intersection(index_t E1, index_t i, index_t j, index_t E2, index_t k, index_t l) override
Given two segments that have an intersection, create the intersection.
void create_enclosing_quad(const vec2 &p1, const vec2 &p2, const vec2 &p3, const vec2 &p4)
Creates a first large enclosing quad.
index_t insert(const vec2 &p, index_t hint=index_t(-1))
Inserts a point.
Sign orient2d(index_t i, index_t j, index_t k) const override
Sign incircle(index_t i, index_t j, index_t k, index_t l) const override
Incircle predicate.
void insert(index_t nb_points, const double *points, index_t *indices=nullptr, bool remove_unreferenced_vertices=false)
Batch-inserts a set of point.
const vec2 point(index_t v) const
Gets a point by index.
void clear() override
Removes everything from this triangulation.
void create_enclosing_triangle(const vec2 &p1, const vec2 &p2, const vec2 &p3)
Creates a first large enclosing triangle.
void create_enclosing_rectangle(double x1, double y1, double x2, double y2)
Creates a first large enclosing rectangle.
void save(const std::string &filename) const override
Saves this CDT to a geogram mesh file.
Base class for constrained Delaunay triangulation.
index_t insert(index_t v, index_t hint=index_t(-1))
Inserts a new point.
bool Tis_in_list(index_t t) const
Tests whether a triangle is in a DList.
index_t locate(index_t v, index_t hint=index_t(-1), Sign *orient=nullptr) const
Locates a vertex.
index_t nv() const
Gets the number of vertices.
void Trot(index_t t, index_t lv)
Rotates indices in triangle t in such a way that a given vertex becomes vertex 0.
void debug_check_consistency() const
Checks both combinatorics and geometry in debug mode, ignored in release mode, aborts on unconsistenc...
CDTBase2d()
CDTBase2d constructor.
void Delaunayize_vertex_neighbors(index_t v, DList &S)
Restores Delaunay condition starting from the triangles incident to a given vertex.
void check_combinatorics() const
Consistency combinatorial check for all the triangles.
static index_t find_3(const index_t *T, index_t v)
Finds the index of an integer in an array of three integers.
index_t find_intersected_edges(index_t i, index_t j, DList &Q)
Finds the edges intersected by a constraint.
void Delaunayize_new_edges(DList &N)
Restores Delaunay condition for a set of edges after inserting a constrained edge.
vector< uint8_t > Tflags_
void insert_vertex_in_edge(index_t v, index_t t, index_t le, DList &S)
Inserts a vertex in an edge.
bool Tedge_is_constrained(index_t t, index_t le) const
Tests whether an edge is constrained.
void insert_vertex_in_triangle(index_t v, index_t t, DList &S)
Inserts a vertex in a triangle.
virtual ~CDTBase2d()
CDTBase2d destructor.
index_t Tv(index_t t, index_t lv) const
Gets a vertex of a triangle.
index_t vT(index_t v) const
Gets a triangle incident to a given vertex.
void create_enclosing_triangle(index_t v1, index_t v2, index_t v3)
Creates the combinatorics for a first large enclosing triangle.
void remove_marked_triangles()
Removes all the triangles that have the flag T_MARKED_FLAG set.
bool segment_edge_intersect(index_t v1, index_t v2, index_t t, index_t le) const
Tests whether an edge triangle and a segment have a frank intersection.
virtual index_t create_intersection(index_t E1, index_t i, index_t j, index_t E2, index_t k, index_t l)=0
Given two segments that have an intersection, create the intersection.
void Tset_edge_cnstr_with_neighbor(index_t t, index_t le, index_t cnstr_id)
Sets an edge as constrained in triangle and in neighbor.
virtual Sign orient2d(index_t i, index_t j, index_t k) const =0
Orientation predicate.
bool Tflag_is_set(index_t t, index_t flag)
Tests a triangle flag.
index_t Tnew()
Creates a new triangle.
virtual Sign incircle(index_t i, index_t j, index_t k, index_t l) const =0
Incircle predicate.
void walk_constraint_t(CDT2d_ConstraintWalker &W, DList &Q)
Used by find_intersected_edges()
void Tset_flag(index_t t, index_t flag)
Sets a triangle flag.
void remove_external_triangles()
Recursively removes all the triangles adjacent to the border, and keeps what's surrounded by constrai...
void debug_check_geometry() const
Consistency geometrical check for all the triangles in debug mode, ignored in release mode.
void check_consistency() const
Checks both combinatorics and geometry, aborts on unconsistency.
void Tadj_back_connect(index_t t1, index_t le1, index_t prev_t2_adj_e2)
After having changed connections from triangle to a neighbor, creates connections from neighbor to tr...
virtual void save(const std::string &filename) const =0
Saves this CDT to a geogram mesh file.
void swap_edge(index_t t1, bool swap_t1_t2=false)
Swaps an edge.
void create_enclosing_quad(index_t v1, index_t v2, index_t v3, index_t v4)
Creates the combinatorics for a first large enclosing quad.
virtual void clear()
Removes everything from this triangulation.
index_t nT() const
Gets the number of triangles.
void debug_Tcheck(index_t t) const
Consistency check for a triangle in debug mode, ignored in release mode.
void walk_constraint_v(CDT2d_ConstraintWalker &W)
Used by find_intersected_edges()
index_t eT(Edge E)
Gets a triangle incident a a given edge.
void insert_vertex_in_edge(index_t v, index_t t, index_t le)
Inserts a vertex in an edge.
index_t Tedge_cnstr(index_t t, index_t le) const
Gets the constraint associated with an edge.
void for_each_T_around_v(index_t v, std::function< bool(index_t t, index_t lv)> doit)
Calls a user-defined function for each triangle around a vertex.
void Tadj_set(index_t t, index_t le, index_t adj)
Sets a triangle adjacency relation.
index_t Tadj_find(index_t t1, index_t t2) const
Finds the edge accross which a triangle is adjacent to another one.
void debug_check_combinatorics() const
Consistency combinatorial check for all the triangles in debug mode, ignored in release mode.
void Tset_edge_cnstr(index_t t, index_t le, index_t cnstr_id)
Sets an edge as constrained.
index_t Tv_find(index_t t, index_t v) const
Finds the local index of a vertex in a triangle.
vector< index_t > Tecnstr_
void constrain_edges_naive(index_t i, index_t j, DList &Q, vector< Edge > &N)
Simpler version of constrain_edges() kept for reference.
void Tset(index_t t, index_t v1, index_t v2, index_t v3, index_t adj1, index_t adj2, index_t adj3, index_t e1cnstr=index_t(-1), index_t e2cnstr=index_t(-1), index_t e3cnstr=index_t(-1))
Sets all the combinatorial information of a triangle and edge flags.
void constrain_edges(index_t i, index_t j, DList &Q, DList &N)
Constrains an edge by iteratively flipping the intersected edges.
void Tcheck(index_t t) const
Consistency check for a triangle.
index_t ncnstr() const
Gets the number of constraints.
void Treset_flag(index_t t, index_t flag)
Resets a triangle flag.
index_t locate_naive(index_t v, index_t hint=index_t(-1), Sign *orient=nullptr) const
Simpler version of locate() kept for reference.
bool is_convex_quad(index_t t) const
Tests whether triange t and its neighbor accross edge 0 form a strictly convex quad.
void insert_constraint(index_t i, index_t j)
Inserts a constraint.
void check_geometry() const
Consistency geometrical check for all the triangles.
bool segment_segment_intersect(index_t u1, index_t u2, index_t v1, index_t v2) const
Tests whether two segments have a frank intersection.
void Delaunayize_new_edges_naive(vector< Edge > &N)
Simpler version of Delaunayize_new_edges() that uses a vector instead of a DList, kept for reference.
index_t Tadj(index_t t, index_t le) const
Gets a triangle adjacent to a triangle.
void set_delaunay(bool delaunay)
Specifies whether a constrained Delaunay triangulation should be constructed, or just a plain constra...
void check_edge_intersections(index_t v1, index_t v2, const DList &Q)
Checks that the edges stored in a DList exactly correspond to all edge intersections between a segmen...
bool Tedge_is_Delaunay(index_t t, index_t le) const
Tests whether a triangle edge is Delaunay.
Vector with aligned memory allocation.
Geometric functions in 2d and 3d.
void clear(void *addr, size_t size)
Clears a memory block.
Global Vorpaline namespace.
void geo_argused(const T &)
Suppresses compiler warnings about unused parameters.
geo_index_t index_t
The type for storing and manipulating indices.
Sign
Integer constants that represent the sign of a value.
Filtered exact predicates for restricted Voronoi diagrams.
Doubly connected triangle list.
bool initialized() const
Tests whether a DList is initialized.
void initialize(index_t list_id)
Initializes a list.
DList(CDTBase2d &cdt, index_t list_id)
Constructs an empty DList.
DList(CDTBase2d &cdt)
Creates an uninitialized DList.