Geogram Version 1.9.6-rc
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
47#include <geogram/mesh/index.h>
48#include <functional>
49
56namespace GEO {
57
62 struct CDT2d_ConstraintWalker;
63
75 class GEOGRAM_API CDTBase2d {
76 public:
81
85 virtual ~CDTBase2d();
86
90 virtual void clear();
91
98
113 void remove_external_triangles(bool remove_internal_holes=false);
114
122 void set_delaunay(bool delaunay) {
123 delaunay_ = delaunay;
124 }
125
129 index_t nT() const {
130 return T_.size()/3;
131 }
132
136 index_t nv() const {
137 return nv_;
138 }
139
143 index_t ncnstr() const {
144 return ncnstr_;
145 }
146
153 index_t Tv(index_t t, index_t lv) const {
154 geo_debug_assert(t<nT());
155 geo_debug_assert(lv<3);
156 return T_[3*t+lv];
157 }
158
166 geo_debug_assert(t<nT());
167 geo_debug_assert(v<nv());
168 return find_3(T_.data()+3*t, v);
169 }
170
179 geo_debug_assert(t<nT());
180 geo_debug_assert(le<3);
181 return Tadj_[3*t+le];
182 }
183
192 geo_debug_assert(t1<nT());
193 geo_debug_assert(t2<nT());
194 return find_3(Tadj_.data()+3*t1, t2);
195 }
196
203 index_t vT(index_t v) const {
204 geo_debug_assert(v < nv());
205 return v2T_[v];
206 }
207
208
236 geo_debug_assert(t < nT());
237 geo_debug_assert(le < 3);
238 return Tecnstr_first_[3*t+le];
239 }
240
249 return ecnstr_next_[ecit];
250 }
251
260 return ecnstr_val_[ecit];
261 }
262
275 index_t result = 0;
276 for(
277 index_t ecit = Tedge_cnstr_first(t,le);
278 ecit != NO_INDEX;
279 ecit = edge_cnstr_next(ecit)
280 ) {
281 ++result;
282 }
283 return result;
284 }
285
286
291 virtual void save(const std::string& filename) const = 0;
292
298
299 protected:
300 virtual void begin_insert_transaction();
301 virtual void commit_insert_transaction();
302 virtual void rollback_insert_transaction();
303
313 index_t insert(index_t v, index_t hint = NO_INDEX);
314
324
334 index_t v1, index_t v2, index_t v3, index_t v4
335 );
336
342 void Tset_flag(index_t t, index_t flag) {
343 geo_debug_assert(t < nT());
344 geo_debug_assert(flag < 8);
345 Tflags_[t] |= Numeric::uint8(1u << flag);
346 }
347
354 geo_debug_assert(t < nT());
355 geo_debug_assert(flag < 8);
356 Tflags_[t] &= Numeric::uint8(~(1u << flag));
357 }
358
367 geo_debug_assert(t < nT());
368 geo_debug_assert(flag < 8);
369 return ((Tflags_[t] & (1u << flag)) != 0);
370 }
371
375 enum {
376 DLIST_S_ID=0,
377 DLIST_Q_ID=1,
378 DLIST_N_ID=2,
379 DLIST_NB=3
380 };
381
385 enum {
386 T_MARKED_FLAG = DLIST_NB,
387 T_VISITED_FLAG = DLIST_NB+1
388 };
389
396 bool Tis_in_list(index_t t) const {
397 return (
398 (Tflags_[t] &
399 Numeric::uint8((1 << DLIST_NB)-1)
400 ) != 0
401 );
402 }
403
416 struct DList {
422 DList(CDTBase2d& cdt, index_t list_id) :
423 cdt_(cdt), list_id_(list_id),
424 back_(NO_INDEX), front_(NO_INDEX) {
425 geo_debug_assert(list_id < DLIST_NB);
426 }
427
437 cdt_(cdt), list_id_(NO_INDEX),
438 back_(NO_INDEX), front_(NO_INDEX) {
439 }
440
445 void initialize(index_t list_id) {
446 geo_debug_assert(!initialized());
447 geo_debug_assert(list_id < DLIST_NB);
448 list_id_ = list_id;
449 }
450
454 bool initialized() const {
455 return (list_id_ != NO_INDEX);
456 }
457
458 ~DList() {
459 if(initialized()) {
460 clear();
461 }
462 }
463
464 bool empty() const {
465 geo_debug_assert(initialized());
467 (back_==NO_INDEX)==(front_==NO_INDEX)
468 );
469 return (back_==NO_INDEX);
470 }
471
472 bool contains(index_t t) const {
473 geo_debug_assert(initialized());
474 return cdt_.Tflag_is_set(t, list_id_);
475 }
476
477 index_t front() const {
478 geo_debug_assert(initialized());
479 return front_;
480 }
481
482 index_t back() const {
483 geo_debug_assert(initialized());
484 return back_;
485 }
486
487 index_t next(index_t t) const {
488 geo_debug_assert(initialized());
489 geo_debug_assert(contains(t));
490 return cdt_.Tnext_[t];
491 }
492
493 index_t prev(index_t t) const {
494 geo_debug_assert(initialized());
495 geo_debug_assert(contains(t));
496 return cdt_.Tprev_[t];
497 }
498
499 void clear() {
500 for(index_t t=front_; t!=NO_INDEX; t = cdt_.Tnext_[t]) {
501 cdt_.Treset_flag(t,list_id_);
502 }
503 back_ = NO_INDEX;
504 front_ = NO_INDEX;
505 }
506
507 index_t size() const {
508 geo_debug_assert(initialized());
509 index_t result = 0;
510 for(index_t t=front(); t!=NO_INDEX; t = next(t)) {
511 ++result;
512 }
513 return result;
514 }
515
516 void push_back(index_t t) {
517 geo_debug_assert(initialized());
518 geo_debug_assert(!cdt_.Tis_in_list(t));
519 cdt_.Tset_flag(t,list_id_);
520 if(empty()) {
521 back_ = t;
522 front_ = t;
523 cdt_.Tnext_[t] = NO_INDEX;
524 cdt_.Tprev_[t] = NO_INDEX;
525 } else {
526 cdt_.Tnext_[t] = NO_INDEX;
527 cdt_.Tnext_[back_] = t;
528 cdt_.Tprev_[t] = back_;
529 back_ = t;
530 }
531 }
532
533 index_t pop_back() {
534 geo_debug_assert(initialized());
535 geo_debug_assert(!empty());
536 index_t t = back_;
537 back_ = cdt_.Tprev_[back_];
538 if(back_ == NO_INDEX) {
539 geo_debug_assert(front_ == t);
540 front_ = NO_INDEX;
541 } else {
542 cdt_.Tnext_[back_] = NO_INDEX;
543 }
544 geo_debug_assert(contains(t));
545 cdt_.Treset_flag(t,list_id_);
546 return t;
547 }
548
549 void push_front(index_t t) {
550 geo_debug_assert(initialized());
551 geo_debug_assert(!cdt_.Tis_in_list(t));
552 cdt_.Tset_flag(t,list_id_);
553 if(empty()) {
554 back_ = t;
555 front_ = t;
556 cdt_.Tnext_[t] = NO_INDEX;
557 cdt_.Tprev_[t] = NO_INDEX;
558 } else {
559 cdt_.Tprev_[t] = NO_INDEX;
560 cdt_.Tprev_[front_] = t;
561 cdt_.Tnext_[t] = front_;
562 front_ = t;
563 }
564 }
565
566 index_t pop_front() {
567 geo_debug_assert(initialized());
568 geo_debug_assert(!empty());
569 index_t t = front_;
570 front_ = cdt_.Tnext_[front_];
571 if(front_ == NO_INDEX) {
572 geo_debug_assert(back_ == t);
573 back_ = NO_INDEX;
574 } else {
575 cdt_.Tprev_[front_] = NO_INDEX;
576 }
577 geo_debug_assert(contains(t));
578 cdt_.Treset_flag(t,list_id_);
579 return t;
580 }
581
582 void remove(index_t t) {
583 geo_debug_assert(initialized());
584 if(t == front_) {
585 pop_front();
586 } else if(t == back_) {
587 pop_back();
588 } else {
589 geo_debug_assert(contains(t));
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_);
595 }
596 }
597
598 void show(std::ostream& out = std::cerr) const {
599 switch(list_id_) {
600 case DLIST_S_ID:
601 out << "S";
602 break;
603 case DLIST_Q_ID:
604 out << "Q";
605 break;
606 case DLIST_N_ID:
607 out << "N";
608 break;
609 case NO_INDEX:
610 out << "<uninitialized list>";
611 break;
612 default:
613 out << "<unknown list id:" << list_id_ << ">";
614 break;
615 }
616 out << "=";
617 for(index_t t=front(); t!=NO_INDEX; t = next(t)) {
618 out << t << ";";
619 }
620 out << std::endl;
621 }
622
623 private:
624 CDTBase2d& cdt_;
625 index_t list_id_;
626 index_t back_;
627 index_t front_;
628 };
629
638
646 DList S(*this);
647 insert_vertex_in_edge(v,t,le,S);
648 }
649
657
673
677 void walk_constraint_v(CDT2d_ConstraintWalker& W);
678
682 void walk_constraint_t(CDT2d_ConstraintWalker& W, DList& Q);
683
695
705
718
726
727
738 void Tset(
739 index_t t,
740 index_t v1, index_t v2, index_t v3,
741 index_t adj1, index_t adj2, index_t adj3,
742 index_t e1cnstr = NO_INDEX,
743 index_t e2cnstr = NO_INDEX,
744 index_t e3cnstr = NO_INDEX
745 ) {
746 geo_debug_assert(t < nT());
747 geo_debug_assert(v1 < nv());
748 geo_debug_assert(v2 < nv());
749 geo_debug_assert(v3 < nv());
750 geo_debug_assert(adj1 < nT() || adj1 == NO_INDEX);
751 geo_debug_assert(adj2 < nT() || adj2 == NO_INDEX);
752 geo_debug_assert(adj3 < nT() || adj3 == NO_INDEX);
753 geo_debug_assert(v1 != v2);
754 geo_debug_assert(v2 != v3);
755 geo_debug_assert(v3 != v1);
756 geo_debug_assert(adj1 != adj2 || adj1 == NO_INDEX);
757 geo_debug_assert(adj2 != adj3 || adj2 == NO_INDEX);
758 geo_debug_assert(adj3 != adj1 || adj3 == NO_INDEX);
759 geo_debug_assert(orient2d(v1,v2,v3) != ZERO);
760 T_[3*t ] = v1;
761 T_[3*t+1] = v2;
762 T_[3*t+2] = v3;
763 Tadj_[3*t ] = adj1;
764 Tadj_[3*t+1] = adj2;
765 Tadj_[3*t+2] = adj3;
766 Tecnstr_first_[3*t] = e1cnstr;
767 Tecnstr_first_[3*t+1] = e2cnstr;
768 Tecnstr_first_[3*t+2] = e3cnstr;
769 v2T_[v1] = t;
770 v2T_[v2] = t;
771 v2T_[v3] = t;
772 }
773
781 void Trot(index_t t, index_t lv) {
782 geo_debug_assert(t < nT());
783 geo_debug_assert(lv < 3);
784 if(lv != 0) {
785 index_t i = 3*t+lv;
786 index_t j = 3*t+((lv+1)%3);
787 index_t k = 3*t+((lv+2)%3);
788 Tset(
789 t,
790 T_[i], T_[j], T_[k],
791 Tadj_[i], Tadj_[j], Tadj_[k],
792 Tecnstr_first_[i], Tecnstr_first_[j], Tecnstr_first_[k]
793 );
794 }
795 }
796
809 void swap_edge(index_t t1, bool swap_t1_t2=false);
810
818 void Tadj_set(index_t t, index_t le, index_t adj) {
819 geo_debug_assert(t < nT());
820 geo_debug_assert(adj < nT());
821 geo_debug_assert(le < 3);
822 Tadj_[3*t+le] = adj;
823 }
824
829 index_t Topp(index_t t, index_t e=0) const {
830 index_t t2 = Tadj(t,e);
831 if(t2 == NO_INDEX) {
832 return NO_INDEX;
833 }
834 index_t e2 = Tadj_find(t2,t);
835 return Tv(t2,e2);
836 }
837
851 index_t t1, index_t le1, index_t prev_t2_adj_e2
852 ) {
853 geo_debug_assert(t1 < nT());
854 geo_debug_assert(le1 < 3);
855 index_t t2 = Tadj(t1,le1);
856 if(t2 == NO_INDEX) {
857 return;
858 }
859 index_t le2 = Tadj_find(t2,prev_t2_adj_e2);
860 Tadj_set(t2,le2,t1);
861 Tset_edge_cnstr_first(t1,le1,Tedge_cnstr_first(t2,le2));
862 }
863
869 index_t t = nT();
870 index_t nc = (t+1)*3; // new number of corners
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);
877 return t;
878 }
879
888 index_t t, index_t le, index_t ecit
889 ) {
890 geo_debug_assert(t < nT());
891 geo_debug_assert(le < 3);
892 Tecnstr_first_[3*t+le] = ecit;
893 }
894
902 index_t t, index_t le, index_t cnstr_id
903 ) {
904 geo_debug_assert(t < nT());
905 geo_debug_assert(le < 3);
906 // Check whether the edge is already constrained with the
907 // same constraint.
908 // TODO (if possible): understand how this can happen and
909 // remove this bloc of code that is not super elegant
910 // (it seems to be when we arrive at j and coming from a vertex
911 // traversed by the edge, both conditions make the constraint
912 // added to the traversed edge).
913 for(
914 index_t ecit = Tedge_cnstr_first(t,le);
915 ecit != NO_INDEX;
916 ecit = edge_cnstr_next(ecit)
917 ) {
918 if(edge_cnstr(ecit) == cnstr_id) {
919 return;
920 }
921 }
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);
925 }
926
935 index_t t, index_t le, index_t cnstr_id
936 ) {
937 geo_debug_assert(t < nT());
938 geo_debug_assert(le < 3);
939#ifdef GEO_DEBUG
940 index_t t_e_cnstr_first = Tedge_cnstr_first(t,le);
941#endif
942 Tadd_edge_cnstr(t, le, cnstr_id);
943 index_t t2 = Tadj(t,le);
944 if(t2 != NO_INDEX) {
945 index_t le2 = Tadj_find(t2,t);
946 // Sanity check: make sure the two edges always share the
947 // same constraint list.
948 geo_debug_assert(Tedge_cnstr_first(t2,le2) == t_e_cnstr_first);
949 Tset_edge_cnstr_first(t2,le2,Tedge_cnstr_first(t,le));
950 }
951 }
952
961 return (Tedge_cnstr_first(t,le) != NO_INDEX);
962 }
963
974 index_t v, std::function<bool(index_t t, index_t lv)> doit
975 ) {
976 index_t t = vT(v);
977 index_t lv = NO_INDEX;
978 do {
979 lv = Tv_find(t,v);
980 if(doit(t,lv)) {
981 return;
982 }
983 t = Tadj(t, (lv+1)%3);
984 } while(t != vT(v) && t != NO_INDEX);
985
986 // We are done, this was an interior vertex
987 if(t != NO_INDEX) {
988 return;
989 }
990
991 // It was a vertex on the border, so we need
992 // to traverse the triangle fan in the other
993 // direction until we reach the border again
994 t = vT(v);
995 lv = Tv_find(t,v);
996 t = Tadj(t, (lv+2)%3);
997 while(t != NO_INDEX) {
998 lv = Tv_find(t,v);
999 if(doit(t,lv)) {
1000 return;
1001 }
1002 t = Tadj(t, (lv+2)%3);
1003 }
1004 }
1005
1006
1018 index_t v, index_t hint = NO_INDEX, Sign* orient = nullptr
1019 ) const;
1020
1029
1035 virtual Sign orient2d(index_t i,index_t j,index_t k) const=0;
1036
1048 virtual Sign incircle(index_t i,index_t j,index_t k,index_t l) const=0;
1049
1070 index_t E1, index_t i, index_t j,
1071 index_t E2, index_t k, index_t l
1072 ) = 0;
1073
1082 static inline index_t find_3(const index_t* T, index_t v) {
1083 // The following expression is 10% faster than using
1084 // if() statements. This uses the C++ norm, that
1085 // ensures that the 'true' boolean value converted to
1086 // an int is always 1. With most compilers, this avoids
1087 // generating branching instructions.
1088 // Thank to Laurent Alonso for this idea.
1089 index_t result = index_t( (T[1] == v) | ((T[2] == v) * 2) );
1090 // Sanity check, important if it was T[0], not explicitly
1091 // tested (detects input that does not meet the precondition).
1092 geo_debug_assert(T[result] == v);
1093 return result;
1094 }
1095
1096 /*******************************************************************/
1097
1107
1108 /******************** Debugging ************************************/
1109
1115 void Tcheck(index_t t) const {
1116 if(t == NO_INDEX) {
1117 return;
1118 }
1119 for(index_t e=0; e<3; ++e) {
1120 geo_assert(Tv(t,e) != Tv(t,(e+1)%3));
1121 if(Tadj(t,e) == NO_INDEX) {
1122 continue;
1123 }
1124 geo_assert(Tadj(t,e) != Tadj(t,(e+1)%3));
1125 index_t t2 = Tadj(t,e);
1126 index_t e2 = Tadj_find(t2,t);
1127 geo_assert(Tadj(t2,e2) == t);
1128 }
1129 }
1130
1137 void debug_Tcheck(index_t t) const {
1138#ifdef GEO_DEBUG
1139 Tcheck(t);
1140#else
1141 geo_argused(t);
1142#endif
1143 }
1144
1149 void check_combinatorics() const {
1150 for(index_t t=0; t<nT(); ++t) {
1151 Tcheck(t);
1152 }
1153 }
1154
1161#ifdef GEO_DEBUG
1162 check_combinatorics();
1163#endif
1164 }
1165
1170 virtual void check_geometry() const;
1171
1178#ifdef GEO_DEBUG
1179 check_geometry();
1180#endif
1181 }
1182
1183
1184 public:
1189 void check_consistency() const {
1190 check_combinatorics();
1191 check_geometry();
1192 }
1193
1194 protected:
1201 debug_check_combinatorics();
1202 debug_check_geometry();
1203 }
1204
1214 index_t u1, index_t u2, index_t v1, index_t v2
1215 ) const {
1216 if(orient2d(u1,u2,v1)*orient2d(u1,u2,v2) > 0) {
1217 return false;
1218 }
1219 return (orient2d(v1,v2,u1)*orient2d(v1,v2,u2) < 0);
1220 }
1221
1233 index_t v1, index_t v2, index_t t, index_t le
1234 ) const {
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);
1238 }
1239
1248 index_t v1, index_t v2, const DList& Q
1249 );
1250
1251 typedef std::pair<index_t, index_t> Edge;
1252
1258 index_t eT(Edge E) {
1259 index_t v1 = E.first;
1260 index_t v2 = E.second;
1261 index_t result = NO_INDEX;
1262 for_each_T_around_v(
1263 v1, [&](index_t t, index_t lv)->bool {
1264 if(Tv(t, (lv+1)%3) == v2) {
1265 if(Tv(t, (lv+2)%3) != v1) {
1266 Trot(t, (lv+2)%3);
1267 }
1268 result = t;
1269 return true;
1270 } else if(Tv(t, (lv+1)%3) == v1) {
1271 if(Tv(t, (lv+2)%3) != v2) {
1272 Trot(t, (lv+2)%3);
1273 }
1274 result = t;
1275 return true;
1276 }
1277 return false;
1278 }
1279 );
1280 geo_debug_assert(result != NO_INDEX);
1282 (Tv(result,1) == v1 && Tv(result,2) == v2) ||
1283 (Tv(result,1) == v2 && Tv(result,2) == v1)
1284 );
1285 return result;
1286 }
1287
1293 index_t v, index_t hint = NO_INDEX, Sign* orient = nullptr
1294 ) const;
1295
1301 index_t i, index_t j, DList& Q, vector<Edge>& N
1302 );
1303
1309
1310 protected:
1311 index_t nv_;
1312 index_t ncnstr_;
1326 };
1327
1328 /*****************************************************************/
1329
1379 class GEOGRAM_API CDT2d: public CDTBase2d {
1380 public:
1381
1382 CDT2d();
1383
1384 ~CDT2d() override;
1385
1389 void clear() override;
1390
1398 const vec2& p1, const vec2& p2, const vec2& p3
1399 );
1400
1409 const vec2& p1, const vec2& p2, const vec2& p3, const vec2& p4
1410 );
1411
1412
1420 double x1, double y1, double x2, double y2
1421 ) {
1422 create_enclosing_quad(
1423 vec2(x1,y1),
1424 vec2(x2,y1),
1425 vec2(x2,y2),
1426 vec2(x1,y2)
1427 );
1428 }
1429
1438 index_t insert(const vec2& p, index_t hint = NO_INDEX) {
1439 debug_check_consistency();
1440 point_.push_back(p);
1441 index_t v = CDTBase2d::insert(point_.size()-1, hint);
1442 // If inserted point already existed in
1443 // triangulation, then nv() did not increase
1444 if(point_.size() > nv()) {
1445 point_.pop_back();
1446 }
1447 debug_check_consistency();
1448 return v;
1449 }
1450
1476 index_t nb_points, const double* points,
1477 index_t* indices = nullptr,
1478 bool remove_unreferenced_vertices = false
1479 );
1480
1484 void save(const std::string& filename) const override;
1485
1491 const vec2 point(index_t v) const {
1492 geo_debug_assert(v < nv());
1493 return point_[v];
1494 }
1495
1496 protected:
1500 Sign orient2d(index_t i, index_t j, index_t k) const override;
1501
1505 Sign incircle(index_t i,index_t j,index_t k,index_t l) const override;
1506
1511 index_t E1, index_t i, index_t j,
1512 index_t E2, index_t k, index_t l
1513 ) override;
1514
1515 protected:
1516 vector<vec2> point_;
1517 };
1518
1519 /*****************************************************************/
1520
1535 class GEOGRAM_API ExactCDT2d : public CDTBase2d {
1536 public:
1537 typedef exact::vec2h ExactPoint;
1538
1543
1547 ~ExactCDT2d() override;
1548
1552 void clear() override;
1553
1566 const ExactPoint& p, index_t id=0, index_t hint = NO_INDEX
1567 );
1568
1578 void insert_constraint(index_t v1, index_t v2, index_t operand_bits=0) {
1579 constraints_.push_back(bindex(v1,v2,bindex::KEEP_ORDER));
1580 cnstr_operand_bits_.push_back(operand_bits);
1581 CDTBase2d::insert_constraint(v1,v2);
1582 }
1583
1592 const ExactPoint& p1, const ExactPoint& p2,
1593 const ExactPoint& p3, const ExactPoint& p4
1594 );
1595
1603 double x1, double y1, double x2, double y2
1604 ) {
1605 create_enclosing_quad(
1606 ExactPoint(vec2(x1,y1)),
1607 ExactPoint(vec2(x2,y1)),
1608 ExactPoint(vec2(x2,y2)),
1609 ExactPoint(vec2(x1,y2))
1610 );
1611 }
1612
1618 const ExactPoint& point(index_t v) const {
1619 geo_debug_assert(v < nv());
1620 return point_[v];
1621 }
1622
1629 geo_debug_assert(v < nv());
1630 return id_[v];
1631 }
1632
1639 geo_debug_assert(v < nv());
1640 id_[v] = id;
1641 }
1642
1658 const std::string& boolean_expression, bool mark_only=false
1659 );
1660
1664 void save(const std::string& filename) const override;
1665
1666 protected:
1667 void add_point(const ExactPoint& p, index_t id = NO_INDEX);
1668 void begin_insert_transaction() override;
1669 void commit_insert_transaction() override;
1670 void rollback_insert_transaction() override;
1671
1675 Sign orient2d(index_t i, index_t j, index_t k) const override;
1676
1680 Sign incircle(index_t i,index_t j,index_t k,index_t l) const override;
1681
1686 index_t E1, index_t i, index_t j,
1687 index_t E2, index_t k, index_t l
1688 ) override;
1689
1690 protected:
1691 vector<ExactPoint> point_;
1692#ifndef GEOGRAM_USE_EXACT_NT
1693 vector<double> length_;
1694#endif
1695 vector<index_t> id_;
1696 vector<index_t> cnstr_operand_bits_;
1697 vector<index_t> facet_inclusion_bits_;
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_;
1701 vector<bindex> constraints_;
1702 };
1703
1704 /*****************************************************************/
1705
1706}
1707
1708#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:196
Constrained Delaunay triangulation.
Definition CDT_2d.h:1379
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.
Definition CDT_2d.h:1438
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:1491
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:1419
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:75
vector< index_t > ecnstr_val_
Definition CDT_2d.h:1318
bool Tis_in_list(index_t t) const
Tests whether a triangle is in a DList.
Definition CDT_2d.h:396
bool exact_incircle_
Definition CDT_2d.h:1324
index_t nv() const
Gets the number of vertices.
Definition CDT_2d.h:136
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:781
void debug_check_consistency() const
Checks both combinatorics and geometry in debug mode, ignored in release mode, aborts on unconsistenc...
Definition CDT_2d.h:1200
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:1149
index_t Tedge_cnstr_first(index_t t, index_t le) const
Gets the constraint associated with an edge.
Definition CDT_2d.h:235
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:1082
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:1316
vector< index_t > Tadj_
Definition CDT_2d.h:1314
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:960
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_
Definition CDT_2d.h:1317
index_t Tv(index_t t, index_t lv) const
Gets a vertex of a triangle.
Definition CDT_2d.h:153
index_t vT(index_t v) const
Gets a triangle incident to a given vertex.
Definition CDT_2d.h:203
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.
Definition CDT_2d.h:738
index_t edge_cnstr_next(index_t ecit) const
Gets the successor of an edge constraint iterator.
Definition CDT_2d.h:248
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.
Definition CDT_2d.h:1232
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.
Definition CDT_2d.h:829
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:366
index_t Tnew()
Creates a new triangle.
Definition CDT_2d.h:868
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.
Definition CDT_2d.h:934
void walk_constraint_t(CDT2d_ConstraintWalker &W, DList &Q)
Used by find_intersected_edges()
vector< index_t > Tnext_
Definition CDT_2d.h:1320
void Tset_flag(index_t t, index_t flag)
Sets a triangle flag.
Definition CDT_2d.h:342
bool exact_intersections_
Definition CDT_2d.h:1325
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.
Definition CDT_2d.h:1177
void check_consistency() const
Checks both combinatorics and geometry, aborts on unconsistency.
Definition CDT_2d.h:1189
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:850
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:129
void debug_Tcheck(index_t t) const
Consistency check for a triangle in debug mode, ignored in release mode.
Definition CDT_2d.h:1137
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:1258
Sign orient_012_
Definition CDT_2d.h:1323
void insert_vertex_in_edge(index_t v, index_t t, index_t le)
Inserts a vertex in an edge.
Definition CDT_2d.h:645
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.
Definition CDT_2d.h:973
void Tadj_set(index_t t, index_t le, index_t adj)
Sets a triangle adjacency relation.
Definition CDT_2d.h:818
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.
Definition CDT_2d.h:259
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:191
void debug_check_combinatorics() const
Consistency combinatorial check for all the triangles in debug mode, ignored in release mode.
Definition CDT_2d.h:1160
void Tadd_edge_cnstr(index_t t, index_t le, index_t cnstr_id)
Adds a constraint to a triangle edge.
Definition CDT_2d.h:901
void Tset_edge_cnstr_first(index_t t, index_t le, index_t ecit)
Sets the constraints list associated with an edge.
Definition CDT_2d.h:887
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:165
vector< index_t > T_
Definition CDT_2d.h:1313
vector< index_t > v2T_
Definition CDT_2d.h:1315
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.
vector< index_t > Tprev_
Definition CDT_2d.h:1321
void Tcheck(index_t t) const
Consistency check for a triangle.
Definition CDT_2d.h:1115
index_t ncnstr() const
Gets the number of constraints.
Definition CDT_2d.h:143
void Treset_flag(index_t t, index_t flag)
Resets a triangle flag.
Definition CDT_2d.h:353
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.
Definition CDT_2d.h:1213
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.
Definition CDT_2d.h:274
index_t Tadj(index_t t, index_t le) const
Gets a triangle adjacent to a triangle.
Definition CDT_2d.h:178
void set_delaunay(bool delaunay)
Specifies whether a constrained Delaunay triangulation should be constructed, or just a plain constra...
Definition CDT_2d.h:122
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_
Definition CDT_2d.h:1319
Constrained Delaunay Triangulation with vertices that are exact points. Can be used to implement 2D C...
Definition CDT_2d.h:1535
~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.
Definition CDT_2d.h:1638
const ExactPoint & point(index_t v) const
Gets a point by vertex index.
Definition CDT_2d.h:1618
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.
Definition CDT_2d.h:1628
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.
Definition CDT_2d.h:1602
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.
Definition CDT_2d.h:1578
2d vector with homogeneous coordinates
Definition vechg.h:60
Vector with aligned memory allocation.
Definition memory.h:660
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.
Definition index.h:205
void clear(void *addr, size_t size)
Clears a memory block.
Definition memory.h:116
uint8_t uint8
Definition numeric.h:135
Global Vorpaline namespace.
Definition basic.h:55
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:329
Sign
Integer constants that represent the sign of a value.
Definition numeric.h:68
@ ZERO
Definition numeric.h:72
Filtered exact predicates for restricted Voronoi diagrams.
Doubly connected triangle list.
Definition CDT_2d.h:416
bool initialized() const
Tests whether a DList is initialized.
Definition CDT_2d.h:454
void initialize(index_t list_id)
Initializes a list.
Definition CDT_2d.h:445
DList(CDTBase2d &cdt, index_t list_id)
Constructs an empty DList.
Definition CDT_2d.h:422
DList(CDTBase2d &cdt)
Creates an uninitialized DList.
Definition CDT_2d.h:436