Geogram Version 1.8.5
A programming library of geometric algorithms
Loading...
Searching...
No Matches
CDT_2d.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2000-2022 Inria
3 * All rights reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 *
8 * * Redistributions of source code must retain the above copyright notice,
9 * this list of conditions and the following disclaimer.
10 * * Redistributions in binary form must reproduce the above copyright notice,
11 * this list of conditions and the following disclaimer in the documentation
12 * and/or other materials provided with the distribution.
13 * * Neither the name of the ALICE Project-Team nor the names of its
14 * contributors may be used to endorse or promote products derived from this
15 * software without specific prior written permission.
16 *
17 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
18 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
21 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27 * POSSIBILITY OF SUCH DAMAGE.
28 *
29 * Contact: Bruno Levy
30 *
31 * https://www.inria.fr/fr/bruno-levy
32 *
33 * Inria,
34 * Domaine de Voluceau,
35 * 78150 Le Chesnay - Rocquencourt
36 * FRANCE
37 *
38 */
39
40#ifndef GEOGRAM_DELAUNAY_CDT_2D
41#define GEOGRAM_DELAUNAY_CDT_2D
42
46#include <functional>
47
54namespace GEO {
55
60 struct CDT2d_ConstraintWalker;
61
73 class GEOGRAM_API CDTBase2d {
74 public:
79
83 virtual ~CDTBase2d();
84
88 virtual void clear();
89
96
102
110 void set_delaunay(bool delaunay) {
111 delaunay_ = delaunay;
112 }
113
117 index_t nT() const {
118 return T_.size()/3;
119 }
120
124 index_t nv() const {
125 return nv_;
126 }
127
131 index_t ncnstr() const {
132 return ncnstr_;
133 }
134
141 index_t Tv(index_t t, index_t lv) const {
142 geo_debug_assert(t<nT());
143 geo_debug_assert(lv<3);
144 return T_[3*t+lv];
145 }
146
154 geo_debug_assert(t<nT());
155 geo_debug_assert(v<nv());
156 return find_3(T_.data()+3*t, v);
157 }
158
167 geo_debug_assert(t<nT());
168 geo_debug_assert(le<3);
169 return Tadj_[3*t+le];
170 }
171
180 geo_debug_assert(t1<nT());
181 geo_debug_assert(t2<nT());
182 return find_3(Tadj_.data()+3*t1, t2);
183 }
184
191 index_t vT(index_t v) const {
192 geo_debug_assert(v < nv());
193 return v2T_[v];
194 }
195
200 virtual void save(const std::string& filename) const = 0;
201
207
208 protected:
219
229
239 index_t v1, index_t v2, index_t v3, index_t v4
240 );
241
247 void Tset_flag(index_t t, index_t flag) {
248 geo_debug_assert(t < nT());
249 geo_debug_assert(flag < 8);
250 Tflags_[t] |= Numeric::uint8(1u << flag);
251 }
252
259 geo_debug_assert(t < nT());
260 geo_debug_assert(flag < 8);
261 Tflags_[t] &= Numeric::uint8(~(1u << flag));
262 }
263
272 geo_debug_assert(t < nT());
273 geo_debug_assert(flag < 8);
274 return ((Tflags_[t] & (1u << flag)) != 0);
275 }
276
280 enum {
281 DLIST_S_ID=0,
282 DLIST_Q_ID=1,
283 DLIST_N_ID=2,
284 DLIST_NB=3
285 };
286
290 enum {
291 T_MARKED_FLAG = DLIST_NB
292 };
293
300 bool Tis_in_list(index_t t) const {
301 return (
302 (Tflags_[t] &
303 Numeric::uint8((1 << DLIST_NB)-1)
304 ) != 0
305 );
306 }
307
320 struct DList {
326 DList(CDTBase2d& cdt, index_t list_id) :
327 cdt_(cdt), list_id_(list_id),
328 back_(index_t(-1)), front_(index_t(-1)) {
329 geo_debug_assert(list_id < DLIST_NB);
330 }
331
341 cdt_(cdt), list_id_(index_t(-1)),
342 back_(index_t(-1)), front_(index_t(-1)) {
343 }
344
349 void initialize(index_t list_id) {
350 geo_debug_assert(!initialized());
351 geo_debug_assert(list_id < DLIST_NB);
352 list_id_ = list_id;
353 }
354
358 bool initialized() const {
359 return (list_id_ != index_t(-1));
360 }
361
362 ~DList() {
363 if(initialized()) {
364 clear();
365 }
366 }
367
368 bool empty() const {
369 geo_debug_assert(initialized());
371 (back_==index_t(-1))==(front_==index_t(-1))
372 );
373 return (back_==index_t(-1));
374 }
375
376 bool contains(index_t t) const {
377 geo_debug_assert(initialized());
378 return cdt_.Tflag_is_set(t, list_id_);
379 }
380
381 index_t front() const {
382 geo_debug_assert(initialized());
383 return front_;
384 }
385
386 index_t back() const {
387 geo_debug_assert(initialized());
388 return back_;
389 }
390
391 index_t next(index_t t) const {
392 geo_debug_assert(initialized());
393 geo_debug_assert(contains(t));
394 return cdt_.Tnext_[t];
395 }
396
397 index_t prev(index_t t) const {
398 geo_debug_assert(initialized());
399 geo_debug_assert(contains(t));
400 return cdt_.Tprev_[t];
401 }
402
403 void clear() {
404 for(index_t t=front_; t!=index_t(-1); t = cdt_.Tnext_[t]) {
405 cdt_.Treset_flag(t,list_id_);
406 }
407 back_ = index_t(-1);
408 front_ = index_t(-1);
409 }
410
411 index_t size() const {
412 geo_debug_assert(initialized());
413 index_t result = 0;
414 for(index_t t=front(); t!=index_t(-1); t = next(t)) {
415 ++result;
416 }
417 return result;
418 }
419
420 void push_back(index_t t) {
421 geo_debug_assert(initialized());
422 geo_debug_assert(!cdt_.Tis_in_list(t));
423 cdt_.Tset_flag(t,list_id_);
424 if(empty()) {
425 back_ = t;
426 front_ = t;
427 cdt_.Tnext_[t] = index_t(-1);
428 cdt_.Tprev_[t] = index_t(-1);
429 } else {
430 cdt_.Tnext_[t] = index_t(-1);
431 cdt_.Tnext_[back_] = t;
432 cdt_.Tprev_[t] = back_;
433 back_ = t;
434 }
435 }
436
437 index_t pop_back() {
438 geo_debug_assert(initialized());
439 geo_debug_assert(!empty());
440 index_t t = back_;
441 back_ = cdt_.Tprev_[back_];
442 if(back_ == index_t(-1)) {
443 geo_debug_assert(front_ == t);
444 front_ = index_t(-1);
445 } else {
446 cdt_.Tnext_[back_] = index_t(-1);
447 }
448 geo_debug_assert(contains(t));
449 cdt_.Treset_flag(t,list_id_);
450 return t;
451 }
452
453 void push_front(index_t t) {
454 geo_debug_assert(initialized());
455 geo_debug_assert(!cdt_.Tis_in_list(t));
456 cdt_.Tset_flag(t,list_id_);
457 if(empty()) {
458 back_ = t;
459 front_ = t;
460 cdt_.Tnext_[t] = index_t(-1);
461 cdt_.Tprev_[t] = index_t(-1);
462 } else {
463 cdt_.Tprev_[t] = index_t(-1);
464 cdt_.Tprev_[front_] = t;
465 cdt_.Tnext_[t] = front_;
466 front_ = t;
467 }
468 }
469
470 index_t pop_front() {
471 geo_debug_assert(initialized());
472 geo_debug_assert(!empty());
473 index_t t = front_;
474 front_ = cdt_.Tnext_[front_];
475 if(front_ == index_t(-1)) {
476 geo_debug_assert(back_ == t);
477 back_ = index_t(-1);
478 } else {
479 cdt_.Tprev_[front_] = index_t(-1);
480 }
481 geo_debug_assert(contains(t));
482 cdt_.Treset_flag(t,list_id_);
483 return t;
484 }
485
486 void remove(index_t t) {
487 geo_debug_assert(initialized());
488 if(t == front_) {
489 pop_front();
490 } else if(t == back_) {
491 pop_back();
492 } else {
493 geo_debug_assert(contains(t));
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_);
499 }
500 }
501
502 void show(std::ostream& out = std::cerr) const {
503 switch(list_id_) {
504 case DLIST_S_ID:
505 out << "S";
506 break;
507 case DLIST_Q_ID:
508 out << "Q";
509 break;
510 case DLIST_N_ID:
511 out << "N";
512 break;
513 case index_t(-1):
514 out << "<uninitialized list>";
515 break;
516 default:
517 out << "<unknown list id:" << list_id_ << ">";
518 break;
519 }
520 out << "=";
521 for(index_t t=front(); t!=index_t(-1); t = next(t)) {
522 out << t << ";";
523 }
524 out << std::endl;
525 }
526
527 private:
528 CDTBase2d& cdt_;
529 index_t list_id_;
530 index_t back_;
531 index_t front_;
532 };
533
542
550 DList S(*this);
551 insert_vertex_in_edge(v,t,le,S);
552 }
553
561
577
581 void walk_constraint_v(CDT2d_ConstraintWalker& W);
582
586 void walk_constraint_t(CDT2d_ConstraintWalker& W, DList& Q);
587
599
612
620
621
632 void Tset(
633 index_t t,
634 index_t v1, index_t v2, index_t v3,
635 index_t adj1, index_t adj2, index_t adj3,
636 index_t e1cnstr = index_t(-1),
637 index_t e2cnstr = index_t(-1),
638 index_t e3cnstr = index_t(-1)
639 ) {
640 geo_debug_assert(t < nT());
641 geo_debug_assert(v1 < nv());
642 geo_debug_assert(v2 < nv());
643 geo_debug_assert(v3 < nv());
644 geo_debug_assert(adj1 < nT() || adj1 == index_t(-1));
645 geo_debug_assert(adj2 < nT() || adj2 == index_t(-1));
646 geo_debug_assert(adj3 < nT() || adj3 == index_t(-1));
647 geo_debug_assert(v1 != v2);
648 geo_debug_assert(v2 != v3);
649 geo_debug_assert(v3 != v1);
650 geo_debug_assert(adj1 != adj2 || adj1 == index_t(-1));
651 geo_debug_assert(adj2 != adj3 || adj2 == index_t(-1));
652 geo_debug_assert(adj3 != adj1 || adj3 == index_t(-1));
653 geo_debug_assert(orient2d(v1,v2,v3) != ZERO);
654 T_[3*t ] = v1;
655 T_[3*t+1] = v2;
656 T_[3*t+2] = v3;
657 Tadj_[3*t ] = adj1;
658 Tadj_[3*t+1] = adj2;
659 Tadj_[3*t+2] = adj3;
660 Tecnstr_[3*t] = e1cnstr;
661 Tecnstr_[3*t+1] = e2cnstr;
662 Tecnstr_[3*t+2] = e3cnstr;
663 v2T_[v1] = t;
664 v2T_[v2] = t;
665 v2T_[v3] = t;
666 }
667
675 void Trot(index_t t, index_t lv) {
676 geo_debug_assert(t < nT());
677 geo_debug_assert(lv < 3);
678 if(lv != 0) {
679 index_t i = 3*t+lv;
680 index_t j = 3*t+((lv+1)%3);
681 index_t k = 3*t+((lv+2)%3);
682 Tset(
683 t,
684 T_[i], T_[j], T_[k],
685 Tadj_[i], Tadj_[j], Tadj_[k],
686 Tecnstr_[i], Tecnstr_[j], Tecnstr_[k]
687 );
688 }
689 }
690
703 void swap_edge(index_t t1, bool swap_t1_t2=false);
704
712 void Tadj_set(index_t t, index_t le, index_t adj) {
713 geo_debug_assert(t < nT());
714 geo_debug_assert(adj < nT());
715 geo_debug_assert(le < 3);
716 Tadj_[3*t+le] = adj;
717 }
718
732 index_t t1, index_t le1, index_t prev_t2_adj_e2
733 ) {
734 geo_debug_assert(t1 < nT());
735 geo_debug_assert(le1 < 3);
736 index_t t2 = Tadj(t1,le1);
737 if(t2 == index_t(-1)) {
738 return;
739 }
740 index_t le2 = Tadj_find(t2,prev_t2_adj_e2);
741 Tadj_set(t2,le2,t1);
742 Tset_edge_cnstr(t1,le1,Tedge_cnstr(t2,le2));
743 }
744
750 index_t t = nT();
751 index_t nc = (t+1)*3; // new number of corners
752 T_.resize(nc, index_t(-1));
753 Tadj_.resize(nc, index_t(-1));
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));
758 return t;
759 }
760
769 geo_debug_assert(t < nT());
770 geo_debug_assert(le < 3);
771 return Tecnstr_[3*t+le];
772 }
773
781 index_t t, index_t le, index_t cnstr_id
782 ) {
783 geo_debug_assert(t < nT());
784 geo_debug_assert(le < 3);
785 Tecnstr_[3*t+le] = cnstr_id;
786 }
787
795 index_t t, index_t le, index_t cnstr_id
796 ) {
797 geo_debug_assert(t < nT());
798 geo_debug_assert(le < 3);
799 Tset_edge_cnstr(t, le, cnstr_id);
800 index_t t2 = Tadj(t,le);
801 if(t2 != index_t(-1)) {
802 index_t le2 = Tadj_find(t2,t);
803 Tset_edge_cnstr(t2,le2,cnstr_id);
804 }
805 }
806
815 return (Tedge_cnstr(t,le) != index_t(-1));
816 }
817
828 index_t v, std::function<bool(index_t t, index_t lv)> doit
829 ) {
830 index_t t = vT(v);
831 index_t lv = index_t(-1);
832 do {
833 lv = Tv_find(t,v);
834 if(doit(t,lv)) {
835 return;
836 }
837 t = Tadj(t, (lv+1)%3);
838 } while(t != vT(v) && t != index_t(-1));
839
840 // We are done, this was an interior vertex
841 if(t != index_t(-1)) {
842 return;
843 }
844
845 // It was a vertex on the border, so we need
846 // to traverse the triangle fan in the other
847 // direction until we reach the border again
848 t = vT(v);
849 lv = Tv_find(t,v);
850 t = Tadj(t, (lv+2)%3);
851 while(t != index_t(-1)) {
852 lv = Tv_find(t,v);
853 if(doit(t,lv)) {
854 return;
855 }
856 t = Tadj(t, (lv+2)%3);
857 }
858 }
859
860
872 index_t v, index_t hint = index_t(-1), Sign* orient = nullptr
873 ) const;
874
882 bool is_convex_quad(index_t t) const;
883
889 virtual Sign orient2d(index_t i,index_t j,index_t k) const=0;
890
902 virtual Sign incircle(index_t i,index_t j,index_t k,index_t l) const=0;
903
924 index_t E1, index_t i, index_t j,
925 index_t E2, index_t k, index_t l
926 ) = 0;
927
936 static inline index_t find_3(const index_t* T, index_t v) {
937 // The following expression is 10% faster than using
938 // if() statements. This uses the C++ norm, that
939 // ensures that the 'true' boolean value converted to
940 // an int is always 1. With most compilers, this avoids
941 // generating branching instructions.
942 // Thank to Laurent Alonso for this idea.
943 index_t result = index_t( (T[1] == v) | ((T[2] == v) * 2) );
944 // Sanity check, important if it was T[0], not explicitly
945 // tested (detects input that does not meet the precondition).
946 geo_debug_assert(T[result] == v);
947 return result;
948 }
949
950 /*******************************************************************/
951
961
962 /******************** Debugging ************************************/
963
969 void Tcheck(index_t t) const {
970 if(t == index_t(-1)) {
971 return;
972 }
973 for(index_t e=0; e<3; ++e) {
974 geo_assert(Tv(t,e) != Tv(t,(e+1)%3));
975 if(Tadj(t,e) == index_t(-1)) {
976 continue;
977 }
978 geo_assert(Tadj(t,e) != Tadj(t,(e+1)%3));
979 index_t t2 = Tadj(t,e);
980 index_t e2 = Tadj_find(t2,t);
981 geo_assert(Tadj(t2,e2) == t);
982 }
983 }
984
991 void debug_Tcheck(index_t t) const {
992#ifdef GEO_DEBUG
993 Tcheck(t);
994#else
995 geo_argused(t);
996#endif
997 }
998
1003 void check_combinatorics() const {
1004 for(index_t t=0; t<nT(); ++t) {
1005 Tcheck(t);
1006 }
1007 }
1008
1015#ifdef GEO_DEBUG
1016 check_combinatorics();
1017#endif
1018 }
1019
1024 void check_geometry() const {
1025 if(delaunay_ && exact_incircle_) {
1026 for(index_t t=0; t<nT(); ++t) {
1027 for(index_t le=0; le<3; ++le) {
1028 geo_assert(Tedge_is_Delaunay(t,le));
1029 }
1030 }
1031 }
1032 }
1033
1040#ifdef GEO_DEBUG
1041 check_geometry();
1042#endif
1043 }
1044
1045
1046 public:
1051 void check_consistency() const {
1052 check_combinatorics();
1053 check_geometry();
1054 }
1055
1056 protected:
1063 debug_check_combinatorics();
1064 debug_check_geometry();
1065 }
1066
1076 index_t u1, index_t u2, index_t v1, index_t v2
1077 ) const {
1078 if(orient2d(u1,u2,v1)*orient2d(u1,u2,v2) > 0) {
1079 return false;
1080 }
1081 return (orient2d(v1,v2,u1)*orient2d(v1,v2,u2) < 0);
1082 }
1083
1095 index_t v1, index_t v2, index_t t, index_t le
1096 ) const {
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);
1100 }
1101
1110 index_t v1, index_t v2, const DList& Q
1111 );
1112
1113 typedef std::pair<index_t, index_t> Edge;
1114
1120 index_t eT(Edge E) {
1121 index_t v1 = E.first;
1122 index_t v2 = E.second;
1123 index_t result = index_t(-1);
1124 for_each_T_around_v(
1125 v1, [&](index_t t, index_t lv)->bool {
1126 if(Tv(t, (lv+1)%3) == v2) {
1127 if(Tv(t, (lv+2)%3) != v1) {
1128 Trot(t, (lv+2)%3);
1129 }
1130 result = t;
1131 return true;
1132 } else if(Tv(t, (lv+1)%3) == v1) {
1133 if(Tv(t, (lv+2)%3) != v2) {
1134 Trot(t, (lv+2)%3);
1135 }
1136 result = t;
1137 return true;
1138 }
1139 return false;
1140 }
1141 );
1142 geo_debug_assert(result != index_t(-1));
1144 (Tv(result,1) == v1 && Tv(result,2) == v2) ||
1145 (Tv(result,1) == v2 && Tv(result,2) == v1)
1146 );
1147 return result;
1148 }
1149
1155 index_t v, index_t hint = index_t(-1), Sign* orient = nullptr
1156 ) const;
1157
1163 index_t i, index_t j, DList& Q, vector<Edge>& N
1164 );
1165
1171
1172 protected:
1173 index_t nv_;
1174 index_t ncnstr_;
1185 };
1186
1187 /*****************************************************************/
1188
1234 class GEOGRAM_API CDT2d: public CDTBase2d {
1235 public:
1236
1237 CDT2d();
1238
1239 ~CDT2d() override;
1240
1244 void clear() override;
1245
1253 const vec2& p1, const vec2& p2, const vec2& p3
1254 );
1255
1264 const vec2& p1, const vec2& p2, const vec2& p3, const vec2& p4
1265 );
1266
1267
1275 double x1, double y1, double x2, double y2
1276 ) {
1277 create_enclosing_quad(
1278 vec2(x1,y1),
1279 vec2(x2,y1),
1280 vec2(x2,y2),
1281 vec2(x1,y2)
1282 );
1283 }
1284
1293 index_t insert(const vec2& p, index_t hint = index_t(-1)) {
1294 debug_check_consistency();
1295 point_.push_back(p);
1296 index_t v = CDTBase2d::insert(point_.size()-1, hint);
1297 // If inserted point already existed in
1298 // triangulation, then nv() did not increase
1299 if(point_.size() > nv()) {
1300 point_.pop_back();
1301 }
1302 debug_check_consistency();
1303 return v;
1304 }
1305
1331 index_t nb_points, const double* points,
1332 index_t* indices = nullptr,
1333 bool remove_unreferenced_vertices = false
1334 );
1335
1339 void save(const std::string& filename) const override;
1340
1346 const vec2 point(index_t v) const {
1347 geo_debug_assert(v < nv());
1348 return point_[v];
1349 }
1350
1351 protected:
1355 Sign orient2d(index_t i, index_t j, index_t k) const override;
1356
1360 Sign incircle(index_t i,index_t j,index_t k,index_t l) const override;
1361
1366 index_t E1, index_t i, index_t j,
1367 index_t E2, index_t k, index_t l
1368 ) override;
1369
1370 protected:
1371 vector<vec2> point_;
1372 };
1373
1374 /*****************************************************************/
1375}
1376
1377#endif
#define geo_assert(x)
Verifies that a condition is met.
Definition assert.h:149
#define geo_debug_assert(x)
Verifies that a condition is met.
Definition assert.h:195
Common include file, providing basic definitions. Should be included before anything else by all head...
Constrained Delaunay triangulation.
Definition CDT_2d.h:1234
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.
Definition CDT_2d.h:1293
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.
Definition CDT_2d.h:1346
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.
Definition CDT_2d.h:1274
void save(const std::string &filename) const override
Saves this CDT to a geogram mesh file.
Base class for constrained Delaunay triangulation.
Definition CDT_2d.h:73
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.
Definition CDT_2d.h:300
bool exact_incircle_
Definition CDT_2d.h:1184
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.
Definition CDT_2d.h:124
void Trot(index_t t, index_t lv)
Rotates indices in triangle t in such a way that a given vertex becomes vertex 0.
Definition CDT_2d.h:675
void debug_check_consistency() const
Checks both combinatorics and geometry in debug mode, ignored in release mode, aborts on unconsistenc...
Definition CDT_2d.h:1062
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.
Definition CDT_2d.h:1003
static index_t find_3(const index_t *T, index_t v)
Finds the index of an integer in an array of three integers.
Definition CDT_2d.h:936
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_
Definition CDT_2d.h:1178
vector< index_t > Tadj_
Definition CDT_2d.h:1176
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.
Definition CDT_2d.h:814
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.
Definition CDT_2d.h:141
index_t vT(index_t v) const
Gets a triangle incident to a given vertex.
Definition CDT_2d.h:191
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.
Definition CDT_2d.h:1094
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.
Definition CDT_2d.h:794
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.
Definition CDT_2d.h:271
index_t Tnew()
Creates a new triangle.
Definition CDT_2d.h:749
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()
vector< index_t > Tnext_
Definition CDT_2d.h:1180
void Tset_flag(index_t t, index_t flag)
Sets a triangle flag.
Definition CDT_2d.h:247
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.
Definition CDT_2d.h:1039
void check_consistency() const
Checks both combinatorics and geometry, aborts on unconsistency.
Definition CDT_2d.h:1051
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...
Definition CDT_2d.h:731
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.
Definition CDT_2d.h:117
void debug_Tcheck(index_t t) const
Consistency check for a triangle in debug mode, ignored in release mode.
Definition CDT_2d.h:991
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.
Definition CDT_2d.h:1120
Sign orient_012_
Definition CDT_2d.h:1183
void insert_vertex_in_edge(index_t v, index_t t, index_t le)
Inserts a vertex in an edge.
Definition CDT_2d.h:549
index_t Tedge_cnstr(index_t t, index_t le) const
Gets the constraint associated with an edge.
Definition CDT_2d.h:768
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.
Definition CDT_2d.h:827
void Tadj_set(index_t t, index_t le, index_t adj)
Sets a triangle adjacency relation.
Definition CDT_2d.h:712
index_t Tadj_find(index_t t1, index_t t2) const
Finds the edge accross which a triangle is adjacent to another one.
Definition CDT_2d.h:179
void debug_check_combinatorics() const
Consistency combinatorial check for all the triangles in debug mode, ignored in release mode.
Definition CDT_2d.h:1014
void Tset_edge_cnstr(index_t t, index_t le, index_t cnstr_id)
Sets an edge as constrained.
Definition CDT_2d.h:780
index_t Tv_find(index_t t, index_t v) const
Finds the local index of a vertex in a triangle.
Definition CDT_2d.h:153
vector< index_t > T_
Definition CDT_2d.h:1175
vector< index_t > v2T_
Definition CDT_2d.h:1177
vector< index_t > Tecnstr_
Definition CDT_2d.h:1179
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.
Definition CDT_2d.h:632
void constrain_edges(index_t i, index_t j, DList &Q, DList &N)
Constrains an edge by iteratively flipping the intersected edges.
vector< index_t > Tprev_
Definition CDT_2d.h:1181
void Tcheck(index_t t) const
Consistency check for a triangle.
Definition CDT_2d.h:969
index_t ncnstr() const
Gets the number of constraints.
Definition CDT_2d.h:131
void Treset_flag(index_t t, index_t flag)
Resets a triangle flag.
Definition CDT_2d.h:258
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.
Definition CDT_2d.h:1024
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.
Definition CDT_2d.h:1075
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.
Definition CDT_2d.h:166
void set_delaunay(bool delaunay)
Specifies whether a constrained Delaunay triangulation should be constructed, or just a plain constra...
Definition CDT_2d.h:110
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.
Definition memory.h:623
Geometric functions in 2d and 3d.
void clear(void *addr, size_t size)
Clears a memory block.
Definition memory.h:104
uint8_t uint8
Definition numeric.h:90
Global Vorpaline namespace.
Definition algorithm.h:64
void geo_argused(const T &)
Suppresses compiler warnings about unused parameters.
Definition argused.h:60
geo_index_t index_t
The type for storing and manipulating indices.
Definition numeric.h:287
Sign
Integer constants that represent the sign of a value.
Definition numeric.h:224
@ ZERO
Definition numeric.h:228
Filtered exact predicates for restricted Voronoi diagrams.
Doubly connected triangle list.
Definition CDT_2d.h:320
bool initialized() const
Tests whether a DList is initialized.
Definition CDT_2d.h:358
void initialize(index_t list_id)
Initializes a list.
Definition CDT_2d.h:349
DList(CDTBase2d &cdt, index_t list_id)
Constructs an empty DList.
Definition CDT_2d.h:326
DList(CDTBase2d &cdt)
Creates an uninitialized DList.
Definition CDT_2d.h:340