40#ifndef GEOGRAM_DELAUNAY_CDT_2D
41#define GEOGRAM_DELAUNAY_CDT_2D
62 struct CDT2d_ConstraintWalker;
114 bool remove_internal_holes=
false
125 delaunay_ = delaunay;
170 return find_3(T_.data()+3*t, v);
183 return Tadj_[3*t+le];
196 return find_3(Tadj_.data()+3*t1, t2);
240 return Tecnstr_first_[3*t+le];
251 return ecnstr_next_[ecit];
262 return ecnstr_val_[ecit];
279 index_t ecit = Tedge_cnstr_first(t,le);
281 ecit = edge_cnstr_next(ecit)
293 virtual void save(
const std::string& filename)
const = 0;
302 virtual void begin_insert_transaction();
303 virtual void commit_insert_transaction();
304 virtual void rollback_insert_transaction();
371 return ((Tflags_[t] & (1u << flag)) != 0);
388 T_MARKED_FLAG = DLIST_NB,
389 T_VISITED_FLAG = DLIST_NB+1
425 cdt_(cdt), list_id_(list_id),
439 cdt_(cdt), list_id_(
index_t(-1)),
457 return (list_id_ !=
index_t(-1));
469 (back_==index_t(-1))==(front_==index_t(-1))
471 return (back_==index_t(-1));
474 bool contains(index_t t)
const {
476 return cdt_.Tflag_is_set(t, list_id_);
489 index_t next(index_t t)
const {
492 return cdt_.Tnext_[t];
495 index_t prev(index_t t)
const {
498 return cdt_.Tprev_[t];
502 for(index_t t=front_; t!=
index_t(-1); t = cdt_.Tnext_[t]) {
503 cdt_.Treset_flag(t,list_id_);
512 for(index_t t=front(); t!=
index_t(-1); t = next(t)) {
518 void push_back(index_t t) {
521 cdt_.Tset_flag(t,list_id_);
529 cdt_.Tnext_[back_] = t;
530 cdt_.Tprev_[t] = back_;
539 back_ = cdt_.Tprev_[back_];
544 cdt_.Tnext_[back_] =
index_t(-1);
547 cdt_.Treset_flag(t,list_id_);
551 void push_front(index_t t) {
554 cdt_.Tset_flag(t,list_id_);
562 cdt_.Tprev_[front_] = t;
563 cdt_.Tnext_[t] = front_;
572 front_ = cdt_.Tnext_[front_];
577 cdt_.Tprev_[front_] =
index_t(-1);
580 cdt_.Treset_flag(t,list_id_);
584 void remove(index_t t) {
588 }
else if(t == back_) {
592 index_t t_prev = cdt_.Tprev_[t];
593 index_t t_next = cdt_.Tnext_[t];
594 cdt_.Tprev_[t_next] = t_prev;
595 cdt_.Tnext_[t_prev] = t_next;
596 cdt_.Treset_flag(t,list_id_);
600 void show(std::ostream& out = std::cerr)
const {
612 out <<
"<uninitialized list>";
615 out <<
"<unknown list id:" << list_id_ <<
">";
619 for(index_t t=front(); t!=
index_t(-1); t = next(t)) {
649 insert_vertex_in_edge(v,t,le,S);
768 Tecnstr_first_[3*t] = e1cnstr;
769 Tecnstr_first_[3*t+1] = e2cnstr;
770 Tecnstr_first_[3*t+2] = e3cnstr;
793 Tadj_[i], Tadj_[j], Tadj_[k],
794 Tecnstr_first_[i], Tecnstr_first_[j], Tecnstr_first_[k]
861 index_t le2 = Tadj_find(t2,prev_t2_adj_e2);
863 Tset_edge_cnstr_first(t1,le1,Tedge_cnstr_first(t2,le2));
875 Tecnstr_first_.resize(nc,
index_t(-1));
876 Tflags_.resize(t+1,0);
877 Tnext_.resize(t+1,
index_t(-1));
878 Tprev_.resize(t+1,
index_t(-1));
894 Tecnstr_first_[3*t+le] = ecit;
916 index_t ecit = Tedge_cnstr_first(t,le);
918 ecit = edge_cnstr_next(ecit)
920 if(edge_cnstr(ecit) == cnstr_id) {
924 ecnstr_val_.push_back(cnstr_id);
925 ecnstr_next_.push_back(Tedge_cnstr_first(t,le));
926 Tset_edge_cnstr_first(t,le, ecnstr_val_.size()-1);
942 index_t t_e_cnstr_first = Tedge_cnstr_first(t,le);
944 Tadd_edge_cnstr(t, le, cnstr_id);
951 Tset_edge_cnstr_first(t2,le2,Tedge_cnstr_first(t,le));
963 return (Tedge_cnstr_first(t,le) !=
index_t(-1));
985 t = Tadj(t, (lv+1)%3);
986 }
while(t != vT(v) && t !=
index_t(-1));
998 t = Tadj(t, (lv+2)%3);
1004 t = Tadj(t, (lv+2)%3);
1123 if(Tadj(t,e) ==
index_t(-1)) {
1152 for(
index_t t=0; t<nT(); ++t) {
1164 check_combinatorics();
1192 check_combinatorics();
1203 debug_check_combinatorics();
1204 debug_check_geometry();
1218 if(orient2d(u1,u2,v1)*orient2d(u1,u2,v2) > 0) {
1221 return (orient2d(v1,v2,u1)*orient2d(v1,v2,u2) < 0);
1237 index_t u1 = Tv(t,(le + 1)%3);
1238 index_t u2 = Tv(t,(le + 2)%3);
1239 return segment_segment_intersect(u1,u2,v1,v2);
1253 typedef std::pair<index_t, index_t> Edge;
1264 for_each_T_around_v(
1266 if(Tv(t, (lv+1)%3) == v2) {
1267 if(Tv(t, (lv+2)%3) != v1) {
1272 }
else if(Tv(t, (lv+1)%3) == v1) {
1273 if(Tv(t, (lv+2)%3) != v2) {
1284 (Tv(result,1) == v1 && Tv(result,2) == v2) ||
1285 (Tv(result,1) == v2 && Tv(result,2) == v1)
1422 double x1,
double y1,
double x2,
double y2
1424 create_enclosing_quad(
1441 debug_check_consistency();
1442 point_.push_back(p);
1443 index_t v = CDTBase2d::insert(point_.size()-1, hint);
1446 if(point_.size() > nv()) {
1449 debug_check_consistency();
1478 index_t nb_points,
const double* points,
1480 bool remove_unreferenced_vertices =
false
1486 void save(
const std::string& filename)
const override;
1581 constraints_.push_back(
bindex(v1,v2,bindex::KEEP_ORDER));
1582 cnstr_operand_bits_.push_back(operand_bits);
1583 CDTBase2d::insert_constraint(v1,v2);
1605 double x1,
double y1,
double x2,
double y2
1607 create_enclosing_quad(
1660 const std::string& boolean_expression,
bool mark_only=
false
1666 void save(
const std::string& filename)
const override;
1670 void begin_insert_transaction()
override;
1671 void commit_insert_transaction()
override;
1672 void rollback_insert_transaction()
override;
1694#ifndef GEOGRAM_USE_EXACT_NT
1700 mutable std::map<trindex, Sign> pred_cache_;
1701 bool use_pred_cache_insert_buffer_;
1702 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.
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.
vector< index_t > ecnstr_val_
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.
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.
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.
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.
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.
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.
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.
index_t insert(const ExactPoint &p, index_t id=0, index_t hint=index_t(-1))
Inserts a point.
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.
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.