40#ifndef GEOGRAM_DELAUNAY_CDT_2D
41#define GEOGRAM_DELAUNAY_CDT_2D
62 struct CDT2d_ConstraintWalker;
123 delaunay_ = delaunay;
168 return find_3(T_.data()+3*t, v);
181 return Tadj_[3*t+le];
194 return find_3(Tadj_.data()+3*t1, t2);
238 return Tecnstr_first_[3*t+le];
249 return ecnstr_next_[ecit];
260 return ecnstr_val_[ecit];
277 index_t ecit = Tedge_cnstr_first(t,le);
279 ecit = edge_cnstr_next(ecit)
291 virtual void save(
const std::string& filename)
const = 0;
300 virtual void begin_insert_transaction();
301 virtual void commit_insert_transaction();
302 virtual void rollback_insert_transaction();
369 return ((Tflags_[t] & (1u << flag)) != 0);
386 T_MARKED_FLAG = DLIST_NB,
387 T_VISITED_FLAG = DLIST_NB+1
423 cdt_(cdt), list_id_(list_id),
424 back_(NO_INDEX), front_(NO_INDEX) {
437 cdt_(cdt), list_id_(NO_INDEX),
438 back_(NO_INDEX), front_(NO_INDEX) {
455 return (list_id_ != NO_INDEX);
467 (back_==NO_INDEX)==(front_==NO_INDEX)
469 return (back_==NO_INDEX);
472 bool contains(index_t t)
const {
474 return cdt_.Tflag_is_set(t, list_id_);
487 index_t next(index_t t)
const {
490 return cdt_.Tnext_[t];
493 index_t prev(index_t t)
const {
496 return cdt_.Tprev_[t];
500 for(index_t t=front_; t!=NO_INDEX; t = cdt_.Tnext_[t]) {
501 cdt_.Treset_flag(t,list_id_);
510 for(index_t t=front(); t!=NO_INDEX; t = next(t)) {
516 void push_back(index_t t) {
519 cdt_.Tset_flag(t,list_id_);
523 cdt_.Tnext_[t] = NO_INDEX;
524 cdt_.Tprev_[t] = NO_INDEX;
526 cdt_.Tnext_[t] = NO_INDEX;
527 cdt_.Tnext_[back_] = t;
528 cdt_.Tprev_[t] = back_;
537 back_ = cdt_.Tprev_[back_];
538 if(back_ == NO_INDEX) {
542 cdt_.Tnext_[back_] = NO_INDEX;
545 cdt_.Treset_flag(t,list_id_);
549 void push_front(index_t t) {
552 cdt_.Tset_flag(t,list_id_);
556 cdt_.Tnext_[t] = NO_INDEX;
557 cdt_.Tprev_[t] = NO_INDEX;
559 cdt_.Tprev_[t] = NO_INDEX;
560 cdt_.Tprev_[front_] = t;
561 cdt_.Tnext_[t] = front_;
570 front_ = cdt_.Tnext_[front_];
571 if(front_ == NO_INDEX) {
575 cdt_.Tprev_[front_] = NO_INDEX;
578 cdt_.Treset_flag(t,list_id_);
582 void remove(index_t t) {
586 }
else if(t == back_) {
590 index_t t_prev = cdt_.Tprev_[t];
591 index_t t_next = cdt_.Tnext_[t];
592 cdt_.Tprev_[t_next] = t_prev;
593 cdt_.Tnext_[t_prev] = t_next;
594 cdt_.Treset_flag(t,list_id_);
598 void show(std::ostream& out = std::cerr)
const {
610 out <<
"<uninitialized list>";
613 out <<
"<unknown list id:" << list_id_ <<
">";
617 for(index_t t=front(); t!=NO_INDEX; t = next(t)) {
647 insert_vertex_in_edge(v,t,le,S);
766 Tecnstr_first_[3*t] = e1cnstr;
767 Tecnstr_first_[3*t+1] = e2cnstr;
768 Tecnstr_first_[3*t+2] = e3cnstr;
791 Tadj_[i], Tadj_[j], Tadj_[k],
792 Tecnstr_first_[i], Tecnstr_first_[j], Tecnstr_first_[k]
859 index_t le2 = Tadj_find(t2,prev_t2_adj_e2);
861 Tset_edge_cnstr_first(t1,le1,Tedge_cnstr_first(t2,le2));
871 T_.resize(nc, NO_INDEX);
872 Tadj_.resize(nc, NO_INDEX);
873 Tecnstr_first_.resize(nc, NO_INDEX);
874 Tflags_.resize(t+1,0);
875 Tnext_.resize(t+1,NO_INDEX);
876 Tprev_.resize(t+1,NO_INDEX);
892 Tecnstr_first_[3*t+le] = ecit;
914 index_t ecit = Tedge_cnstr_first(t,le);
916 ecit = edge_cnstr_next(ecit)
918 if(edge_cnstr(ecit) == cnstr_id) {
922 ecnstr_val_.push_back(cnstr_id);
923 ecnstr_next_.push_back(Tedge_cnstr_first(t,le));
924 Tset_edge_cnstr_first(t,le, ecnstr_val_.size()-1);
940 index_t t_e_cnstr_first = Tedge_cnstr_first(t,le);
942 Tadd_edge_cnstr(t, le, cnstr_id);
949 Tset_edge_cnstr_first(t2,le2,Tedge_cnstr_first(t,le));
961 return (Tedge_cnstr_first(t,le) != NO_INDEX);
983 t = Tadj(t, (lv+1)%3);
984 }
while(t != vT(v) && t != NO_INDEX);
996 t = Tadj(t, (lv+2)%3);
997 while(t != NO_INDEX) {
1002 t = Tadj(t, (lv+2)%3);
1121 if(Tadj(t,e) == NO_INDEX) {
1150 for(
index_t t=0; t<nT(); ++t) {
1162 check_combinatorics();
1190 check_combinatorics();
1201 debug_check_combinatorics();
1202 debug_check_geometry();
1216 if(orient2d(u1,u2,v1)*orient2d(u1,u2,v2) > 0) {
1219 return (orient2d(v1,v2,u1)*orient2d(v1,v2,u2) < 0);
1235 index_t u1 = Tv(t,(le + 1)%3);
1236 index_t u2 = Tv(t,(le + 2)%3);
1237 return segment_segment_intersect(u1,u2,v1,v2);
1251 typedef std::pair<index_t, index_t> Edge;
1262 for_each_T_around_v(
1264 if(Tv(t, (lv+1)%3) == v2) {
1265 if(Tv(t, (lv+2)%3) != v1) {
1270 }
else if(Tv(t, (lv+1)%3) == v1) {
1271 if(Tv(t, (lv+2)%3) != v2) {
1282 (Tv(result,1) == v1 && Tv(result,2) == v2) ||
1283 (Tv(result,1) == v2 && Tv(result,2) == v1)
1420 double x1,
double y1,
double x2,
double y2
1422 create_enclosing_quad(
1439 debug_check_consistency();
1440 point_.push_back(p);
1441 index_t v = CDTBase2d::insert(point_.size()-1, hint);
1444 if(point_.size() > nv()) {
1447 debug_check_consistency();
1476 index_t nb_points,
const double* points,
1478 bool remove_unreferenced_vertices =
false
1484 void save(
const std::string& filename)
const override;
1579 constraints_.push_back(
bindex(v1,v2,bindex::KEEP_ORDER));
1580 cnstr_operand_bits_.push_back(operand_bits);
1581 CDTBase2d::insert_constraint(v1,v2);
1603 double x1,
double y1,
double x2,
double y2
1605 create_enclosing_quad(
1658 const std::string& boolean_expression,
bool mark_only=
false
1664 void save(
const std::string& filename)
const override;
1668 void begin_insert_transaction()
override;
1669 void commit_insert_transaction()
override;
1670 void rollback_insert_transaction()
override;
1692#ifndef GEOGRAM_USE_EXACT_NT
1698 mutable std::map<trindex, Sign> pred_cache_;
1699 bool use_pred_cache_insert_buffer_;
1700 mutable std::vector<std::pair<trindex, Sign>> pred_cache_insert_buffer_;
#define geo_assert(x)
Verifies that a condition is met.
#define geo_debug_assert(x)
Verifies that a condition is met.
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=NO_INDEX)
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.
vector< index_t > ecnstr_val_
bool Tis_in_list(index_t t) const
Tests whether a triangle is in a DList.
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.
index_t Tedge_cnstr_first(index_t t, index_t le) const
Gets the constraint associated with an edge.
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.
vector< index_t > Tecnstr_first_
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 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=NO_INDEX, index_t e2cnstr=NO_INDEX, index_t e3cnstr=NO_INDEX)
Sets all the combinatorial information of a triangle and edge flags.
index_t edge_cnstr_next(index_t ecit) const
Gets the successor of an edge constraint iterator.
void Delaunayize_vertex_neighbors(index_t from_v)
Restores Delaunay condition starting from the triangles incident to a given vertex.
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.
index_t Topp(index_t t, index_t e=0) const
Gets the neighboring triangle vertex opposite to a given vertex.
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 Tadd_edge_cnstr_with_neighbor(index_t t, index_t le, index_t cnstr_id)
Adds a constraint to a triangle edge and to the neighboring edge if it exists.
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.
bool exact_intersections_
void remove_external_triangles(bool remove_internal_holes=false)
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 locate_naive(index_t v, index_t hint=NO_INDEX, Sign *orient=nullptr) const
Simpler version of locate() kept for reference.
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.
virtual void check_geometry() const
Consistency geometrical check for all the triangles.
index_t edge_cnstr(index_t ecit) const
Gets an edge constraint from an edge constraint iterator.
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 Tadd_edge_cnstr(index_t t, index_t le, index_t cnstr_id)
Adds a constraint to a triangle edge.
void Tset_edge_cnstr_first(index_t t, index_t le, index_t ecit)
Sets the constraints list associated with an edge.
index_t Tv_find(index_t t, index_t v) const
Finds the local index of a vertex in a triangle.
void constrain_edges_naive(index_t i, index_t j, DList &Q, vector< Edge > &N)
Simpler version of constrain_edges() kept for reference.
index_t insert(index_t v, index_t hint=NO_INDEX)
Inserts a new point.
void constrain_edges(index_t i, index_t j, DList &Q, DList &N)
Constrains an edge by iteratively flipping the intersected edges.
index_t locate(index_t v, index_t hint=NO_INDEX, Sign *orient=nullptr) const
Locates a vertex.
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.
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.
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 Tedge_cnstr_nb(index_t t, index_t le) const
Gets the number of constraints associated with a triange edge.
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< index_t > ecnstr_next_
Constrained Delaunay Triangulation with vertices that are exact points. Can be used to implement 2D C...
~ExactCDT2d() override
ExactCDT2d destructor.
index_t insert(const ExactPoint &p, index_t id=0, index_t hint=NO_INDEX)
Inserts a point.
ExactCDT2d()
ExactCDT2d constructor.
void create_enclosing_quad(const ExactPoint &p1, const ExactPoint &p2, const ExactPoint &p3, const ExactPoint &p4)
Creates a first large enclosing quad.
void classify_triangles(const std::string &boolean_expression, bool mark_only=false)
Used by 2D CSG operations, discards triangles according to a boolean operation.
Sign orient2d(index_t i, index_t j, index_t k) const override
void set_vertex_id(index_t v, index_t id)
Sets a vertex id by vertex index.
const ExactPoint & point(index_t v) const
Gets a point by vertex index.
Sign incircle(index_t i, index_t j, index_t k, index_t l) const override
Incircle predicate.
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.
index_t vertex_id(index_t v) const
Gets a vertex id by vertex index.
void clear() override
Removes everything from this triangulation.
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
void insert_constraint(index_t v1, index_t v2, index_t operand_bits=0)
Inserts a constraint.
2d vector with homogeneous coordinates
Vector with aligned memory allocation.
Exact predicates and constructs.
Common include file, providing basic definitions. Should be included before anything else by all head...
Geometric functions in 2d and 3d.
Classes for managing tuples of indices.
basic_bindex< index_t > bindex
A basic_bindex made of 2 unsigned integers.
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.