GCC Code Coverage Report

Directory: ./
File: lib/geogram/mesh/mesh_surface_intersection_internal.cpp
Date: 2026-09-07 02:37:58
Exec Total Coverage
Lines: 551 675 81.6%
Functions: 27 32 84.4%
Branches: 562 1380 40.7%
Line Branch Exec Source
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 #include <geogram/mesh/mesh_surface_intersection_internal.h>
41 #include <geogram/mesh/mesh_surface_intersection.h>
42 #include <geogram/basic/debug_stream.h>
43 #include <geogram/basic/boolean_expression.h>
44 #include <stack>
45
46 namespace {
47 using namespace GEO;
48
49 /**
50 * \brief Computes the exact intersection between the support
51 * planes of three triangles
52 * \param[in] p1 , p2 , p3 the three vertices of the first triangle
53 * \param[in] q1 , q2 , q3 the three vertices of the second triangle
54 * \param[in] r1 , r2 , r3 the three vertices of the third triangle
55 * \param[out] result the exact intersection between the three planes
56 * if it exsists
57 * \retval true if the planes have an intersection
58 * \retval false otherwise
59 */
60 5684 bool get_three_planes_intersection(
61 MeshSurfaceIntersection::ExactPoint& result,
62 const vec3& p1, const vec3& p2, const vec3& p3,
63 const vec3& q1, const vec3& q2, const vec3& q3,
64 const vec3& r1, const vec3& r2, const vec3& r3
65 ) {
66
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 exact::vec3 N1 = triangle_normal<exact::vec3>(p1,p2,p3);
67
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 exact::vec3 N2 = triangle_normal<exact::vec3>(q1,q2,q3);
68
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 exact::vec3 N3 = triangle_normal<exact::vec3>(r1,r2,r3);
69
70 exact::vec3 B(
71
2/4
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5684 times.
✗ Branch 5 not taken.
 11368 dot(N1,exact::vec3(p1)),
72
2/4
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5684 times.
✗ Branch 5 not taken.
 11368 dot(N2,exact::vec3(q1)),
73
2/4
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5684 times.
✗ Branch 5 not taken.
 11368 dot(N3,exact::vec3(r1))
74
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 );
75
76
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 11368 result.w = det3x3(
77 N1.x, N1.y, N1.z,
78 N2.x, N2.y, N2.z,
79 N3.x, N3.y, N3.z
80
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 );
81
82
3/4
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4343 times.
✓ Branch 4 taken 1341 times.
 5684 if(result.w.sign() == ZERO) {
83 4343 return false;
84 }
85
86
1/2
✓ Branch 1 taken 1341 times.
✗ Branch 2 not taken.
 2682 result.x = det3x3(
87 B.x, N1.y, N1.z,
88 B.y, N2.y, N2.z,
89 B.z, N3.y, N3.z
90
1/2
✓ Branch 1 taken 1341 times.
✗ Branch 2 not taken.
 1341 );
91
92
1/2
✓ Branch 1 taken 1341 times.
✗ Branch 2 not taken.
 2682 result.y = det3x3(
93 N1.x, B.x, N1.z,
94 N2.x, B.y, N2.z,
95 N3.x, B.z, N3.z
96
1/2
✓ Branch 1 taken 1341 times.
✗ Branch 2 not taken.
 1341 );
97
98
1/2
✓ Branch 1 taken 1341 times.
✗ Branch 2 not taken.
 2682 result.z = det3x3(
99 N1.x, N1.y, B.x,
100 N2.x, N2.y, B.y,
101 N3.x, N3.y, B.z
102
1/2
✓ Branch 1 taken 1341 times.
✗ Branch 2 not taken.
 1341 );
103
104 1341 return true;
105 5684 }
106
107 /**
108 * \brief Computes the exact intersection between the support plane
109 * of a triangle and the support line of a segment
110 * \pre The intersection exists
111 * \param[in] p1 , p2 , p3 the three vertices of the triangle
112 * \param[in] q1 , q2 the two vertices of the segment
113 * \return the exact intersection between the plane and the line
114 */
115 67556 MeshSurfaceIntersection::ExactPoint plane_line_intersection(
116 const vec3& p1, const vec3& p2, const vec3& p3,
117 const vec3& q1, const vec3& q2
118 ) {
119 // Moller & Trumbore's algorithm
120 // see: https://stackoverflow.com/questions/42740765/
121 // intersection-between-line-and-triangle-in-3d
122
1/2
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
 67556 exact::vec3 D = make_vec3<exact::vec3>(q1,q2);
123
1/2
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
 67556 exact::vec3 E1 = make_vec3<exact::vec3>(p1,p2);
124
1/2
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
 67556 exact::vec3 E2 = make_vec3<exact::vec3>(p1,p3);
125
1/2
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
 67556 exact::vec3 AO = make_vec3<exact::vec3>(p1,q1);
126
1/2
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
 67556 exact::vec3 N = cross(E1,E2);
127
2/4
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 67556 times.
✗ Branch 5 not taken.
 67556 exact::scalar d = -dot(D,N);
128
2/8
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 67556 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 67556 geo_debug_assert(d.sign() != ZERO);
129
2/4
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 67556 times.
✗ Branch 5 not taken.
 67556 exact::rational t(dot(AO,N),d);
130
1/2
✓ Branch 1 taken 67556 times.
✗ Branch 2 not taken.
 135112 return mix(t,q1,q2);
131 67556 }
132 }
133
134 namespace GEO {
135
136 ✗ void MeshInTriangle::Vertex::print(std::ostream& out) const {
137 ✗ if(sym.f1 != NO_INDEX) {
138 ✗ out << " ( ";
139 ✗ out << sym.f1;
140 ✗ out << region_to_string(sym.R1).substr(2);
141 }
142 ✗ if(sym.f2 != NO_INDEX) {
143 ✗ out << " /\\ ";
144 ✗ out << sym.f2;
145 ✗ out << region_to_string(sym.R2).substr(2);
146 }
147 ✗ if(sym.f1 != NO_INDEX) {
148 ✗ out << " ) ";
149 }
150 ✗ }
151
152 286652 MeshInTriangle::ExactPoint MeshInTriangle::Vertex::compute_geometry() {
153 // Case 1: f1 vertex
154
2/2
✓ Branch 1 taken 46515 times.
✓ Branch 2 taken 240137 times.
 286652 if(region_dim(sym.R1) == 0) {
155 46515 index_t lv = index_t(sym.R1) - index_t(T1_RGN_P0);
156
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 46515 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 46515 geo_debug_assert(lv < 3);
157
2/4
✓ Branch 1 taken 46515 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 46515 times.
✗ Branch 5 not taken.
 46515 mesh_vertex_index = mesh().facets.vertex(sym.f1,lv);
158
1/2
✓ Branch 1 taken 46515 times.
✗ Branch 2 not taken.
 46515 vec3 p = mit->mesh_vertex(mesh_vertex_index);
159
1/2
✓ Branch 1 taken 46515 times.
✗ Branch 2 not taken.
 46515 return ExactPoint(p);
160 }
161
162
2/8
✓ Branch 0 taken 240137 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 240137 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 240137 geo_debug_assert(sym.f1 != NO_INDEX && sym.f2 != NO_INDEX);
163
164 // Case 2: f2 vertex
165
2/2
✓ Branch 1 taken 21655 times.
✓ Branch 2 taken 218482 times.
 240137 if(region_dim(sym.R2) == 0) {
166 21655 index_t lv = index_t(sym.R2) - index_t(T2_RGN_P0);
167
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 21655 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 21655 geo_debug_assert(lv < 3);
168
2/4
✓ Branch 1 taken 21655 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 21655 times.
✗ Branch 5 not taken.
 21655 mesh_vertex_index = mesh().facets.vertex(sym.f2, lv);
169
1/2
✓ Branch 1 taken 21655 times.
✗ Branch 2 not taken.
 21655 vec3 p = mit->mesh_vertex(mesh_vertex_index);
170
1/2
✓ Branch 1 taken 21655 times.
✗ Branch 2 not taken.
 21655 return ExactPoint(p);
171 }
172
173 // case 3: f1 /\ f2 edge in 3D or f1 edge /\ f2 edge in 3D
174 218482 if(
175
5/6
✓ Branch 1 taken 191062 times.
✓ Branch 2 taken 27420 times.
✓ Branch 4 taken 191062 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 191062 times.
✓ Branch 7 taken 27420 times.
 436964 (region_dim(sym.R1) == 2 || region_dim(sym.R1) == 1) &&
176
2/2
✓ Branch 1 taken 191062 times.
✓ Branch 2 taken 27420 times.
 218482 region_dim(sym.R2) == 1
177 ) {
178
1/2
✓ Branch 1 taken 191062 times.
✗ Branch 2 not taken.
 191062 vec3 p1 = mit->mesh_facet_vertex(sym.f1, 0);
179
1/2
✓ Branch 1 taken 191062 times.
✗ Branch 2 not taken.
 191062 vec3 p2 = mit->mesh_facet_vertex(sym.f1, 1);
180
1/2
✓ Branch 1 taken 191062 times.
✗ Branch 2 not taken.
 191062 vec3 p3 = mit->mesh_facet_vertex(sym.f1, 2);
181 191062 index_t e = index_t(sym.R2)-index_t(T2_RGN_E0);
182
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 191062 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 191062 geo_debug_assert(e<3);
183
1/2
✓ Branch 1 taken 191062 times.
✗ Branch 2 not taken.
 191062 vec3 q1 = mit->mesh_facet_vertex(sym.f2, (e+1)%3);
184
1/2
✓ Branch 1 taken 191062 times.
✗ Branch 2 not taken.
 191062 vec3 q2 = mit->mesh_facet_vertex(sym.f2, (e+2)%3);
185
186 bool seg_seg_two_D = (
187
1/2
✓ Branch 1 taken 191062 times.
✗ Branch 2 not taken.
 191062 region_dim(sym.R1) == 1 &&
188
5/6
✓ Branch 0 taken 163642 times.
✓ Branch 1 taken 27420 times.
✓ Branch 3 taken 163642 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 150982 times.
✓ Branch 6 taken 12660 times.
 342044 PCK::orient_3d(p1,p2,p3,q1) == ZERO &&
189
2/4
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 150982 times.
✗ Branch 4 not taken.
 150982 PCK::orient_3d(p1,p2,p3,q2) == ZERO) ;
190
191
2/2
✓ Branch 0 taken 40080 times.
✓ Branch 1 taken 150982 times.
 191062 if(!seg_seg_two_D) {
192
1/2
✓ Branch 1 taken 40080 times.
✗ Branch 2 not taken.
 40080 return plane_line_intersection(p1,p2,p3,q1,q2);
193 }
194 }
195
196 // case 4: f1 edge /\ f2
197
5/6
✓ Branch 1 taken 178402 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 27420 times.
✓ Branch 5 taken 150982 times.
✓ Branch 6 taken 27420 times.
✓ Branch 7 taken 150982 times.
 178402 if(region_dim(sym.R1) == 1 && region_dim(sym.R2) == 2) {
198 27420 index_t e = index_t(sym.R1)-index_t(T1_RGN_E0);
199
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 27420 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 27420 geo_debug_assert(e<3);
200
1/2
✓ Branch 1 taken 27420 times.
✗ Branch 2 not taken.
 27420 vec3 p1 = mit->mesh_facet_vertex(sym.f2,0);
201
1/2
✓ Branch 1 taken 27420 times.
✗ Branch 2 not taken.
 27420 vec3 p2 = mit->mesh_facet_vertex(sym.f2,1);
202
1/2
✓ Branch 1 taken 27420 times.
✗ Branch 2 not taken.
 27420 vec3 p3 = mit->mesh_facet_vertex(sym.f2,2);
203
1/2
✓ Branch 1 taken 27420 times.
✗ Branch 2 not taken.
 27420 vec3 q1 = mit->mesh_facet_vertex(sym.f1, (e+1)%3);
204
1/2
✓ Branch 1 taken 27420 times.
✗ Branch 2 not taken.
 27420 vec3 q2 = mit->mesh_facet_vertex(sym.f1, (e+2)%3);
205
1/2
✓ Branch 1 taken 27420 times.
✗ Branch 2 not taken.
 27420 return plane_line_intersection(p1,p2,p3,q1,q2);
206 }
207
208 // case 5: f1 edge /\ f2 edge in 2D
209
3/6
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 150982 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 150982 times.
✗ Branch 7 not taken.
 150982 if(region_dim(sym.R1) == 1 && region_dim(sym.R2) == 1) {
210 150982 index_t e1 = index_t(sym.R1) - index_t(T1_RGN_E0);
211
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 150982 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 150982 geo_debug_assert(e1 < 3);
212 150982 index_t e2 = index_t(sym.R2) - index_t(T2_RGN_E0);
213
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 150982 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 150982 geo_debug_assert(e2 < 3);
214
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 vec2 p1 = mit->mesh_facet_vertex_UV(sym.f1, (e1+1)%3);
215
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 vec2 p2 = mit->mesh_facet_vertex_UV(sym.f1, (e1+2)%3);
216
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 vec2 q1 = mit->mesh_facet_vertex_UV(sym.f2, (e2+1)%3);
217
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 vec2 q2 = mit->mesh_facet_vertex_UV(sym.f2, (e2+2)%3);
218
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 vec3 P1 = mit->mesh_facet_vertex(sym.f1, (e1+1)%3);
219
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 vec3 P2 = mit->mesh_facet_vertex(sym.f1, (e1+2)%3);
220
221
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 exact::vec2 D1 = make_vec2<exact::vec2>(p1,p2);
222
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 exact::vec2 D2 = make_vec2<exact::vec2>(q1,q2);
223
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 exact::scalar d = det(D1,D2);
224
2/8
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 150982 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 150982 geo_debug_assert(d.sign() != ZERO);
225
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 exact::vec2 AO = make_vec2<exact::vec2>(p1,q1);
226
2/4
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 150982 times.
✗ Branch 5 not taken.
 150982 exact::rational t(det(AO,D2),d);
227
1/2
✓ Branch 1 taken 150982 times.
✗ Branch 2 not taken.
 150982 return mix(t,P1,P2);
228 150982 }
229
230 // Normally we enumerated all possible cases
231 ✗ geo_assert_not_reached;
232 }
233
234 292336 void MeshInTriangle::Vertex::init_geometry(const ExactPoint& P) {
235 292336 point_exact = P;
236 292336 Numeric::optimize_number_representation(point_exact);
237 #ifndef GEOGRAM_USE_EXACT_NT
238
8/16
✓ Branch 1 taken 292336 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 292336 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 292336 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 292336 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 292336 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 292336 times.
✗ Branch 17 not taken.
✓ Branch 19 taken 292336 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 292336 times.
✗ Branch 23 not taken.
 584672 l = (geo_sqr(P[mit->u_]) + geo_sqr(P[mit->v_])).estimate() /
239
3/6
✓ Branch 1 taken 292336 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 292336 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 292336 times.
✗ Branch 8 not taken.
 292336 geo_sqr(P.w).estimate() ;
240 #endif
241 292336 }
242
243 229 MeshInTriangle::MeshInTriangle(MeshSurfaceIntersection& EM) :
244 229 exact_mesh_(EM),
245 458 mesh_(EM.readonly_mesh()),
246 229 f1_(NO_INDEX),
247 229 dry_run_(false),
248 229 use_pred_cache_insert_buffer_(false)
249 {
250 #ifdef GEOGRAM_USE_EXACT_NT
251 CDTBase2d::exact_incircle_ = true;
252 #else
253 // Since incircle() with expansions computes approximated
254 // lifted coordinate, we need to activate additional
255 // checks for Delaunayization.
256 229 CDTBase2d::exact_incircle_ = false;
257 #endif
258 229 }
259
260 15505 void MeshInTriangle::clear() {
261 15505 vertex_.resize(0);
262 15505 edges_.resize(0);
263 15505 f1_ = NO_INDEX;
264 15505 pred_cache_.clear();
265 15505 pred_cache_insert_buffer_.resize(0);
266 15505 use_pred_cache_insert_buffer_ = false;
267 15505 CDTBase2d::clear();
268 15505 }
269
270 15505 void MeshInTriangle::begin_facet(index_t f) {
271 15505 f1_ = f;
272
273 15505 latest_f2_ = NO_INDEX;
274 15505 latest_f2_count_ = 0;
275
276
1/2
✓ Branch 1 taken 15505 times.
✗ Branch 2 not taken.
 15505 vec3 p1 = mesh_facet_vertex(f,0);
277
1/2
✓ Branch 1 taken 15505 times.
✗ Branch 2 not taken.
 15505 vec3 p2 = mesh_facet_vertex(f,1);
278
1/2
✓ Branch 1 taken 15505 times.
✗ Branch 2 not taken.
 15505 vec3 p3 = mesh_facet_vertex(f,2);
279
280
2/8
✓ Branch 1 taken 15505 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 15505 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 15505 geo_debug_assert(!PCK::aligned_3d(p1,p2,p3));
281
282
1/2
✓ Branch 1 taken 15505 times.
✗ Branch 2 not taken.
 15505 f1_normal_axis_ = PCK::triangle_normal_axis(
283 p1,p2,p3
284 );
285
286 15505 u_ = coord_index_t((f1_normal_axis_ + 1) % 3);
287 15505 v_ = coord_index_t((f1_normal_axis_ + 2) % 3);
288
2/2
✓ Branch 0 taken 46515 times.
✓ Branch 1 taken 15505 times.
 62020 for(index_t lv=0; lv<3; ++lv) {
289
2/4
✓ Branch 1 taken 46515 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 46515 times.
✗ Branch 5 not taken.
 46515 vertex_.push_back(Vertex(this, f, lv));
290 }
291
292
1/2
✓ Branch 1 taken 15505 times.
✗ Branch 2 not taken.
 15505 CDTBase2d::create_enclosing_triangle(0,1,2);
293
294
1/2
✓ Branch 2 taken 15505 times.
✗ Branch 3 not taken.
 15505 edges_.push_back(Edge(1,2));
295
1/2
✓ Branch 2 taken 15505 times.
✗ Branch 3 not taken.
 15505 edges_.push_back(Edge(2,0));
296
1/2
✓ Branch 2 taken 15505 times.
✗ Branch 3 not taken.
 15505 edges_.push_back(Edge(0,1));
297
298 15505 has_planar_isect_ = false;
299 15505 }
300
301 262922 index_t MeshInTriangle::add_vertex(
302 index_t f2, TriangleRegion R1, TriangleRegion R2
303 ) {
304
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 262922 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 262922 geo_debug_assert(f1_ != NO_INDEX);
305
306 // If the same f2 comes more than twice, then
307 // we got a planar facet /\ facet intersection
308 // (and it is good to know it, see get_constraints())
309
3/4
✓ Branch 0 taken 262922 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 153431 times.
✓ Branch 3 taken 109491 times.
 262922 if(f2 != NO_INDEX && f2 == latest_f2_) {
310 153431 ++latest_f2_count_;
311
2/2
✓ Branch 0 taken 36687 times.
✓ Branch 1 taken 116744 times.
 153431 if(latest_f2_count_ > 2) {
312 36687 has_planar_isect_ = true;
313 }
314 } else {
315 109491 latest_f2_ = f2;
316 109491 latest_f2_count_ = 0;
317 }
318
319 // If vertex is a macro-vertex, return it directly.
320
2/2
✓ Branch 1 taken 22785 times.
✓ Branch 2 taken 240137 times.
 262922 if(region_dim(R1) == 0) {
321 22785 return index_t(R1);
322 }
323
324 // Create the vertex
325
2/4
✓ Branch 1 taken 240137 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 240137 times.
✗ Branch 5 not taken.
 240137 vertex_.push_back(Vertex(this, f1_, f2, R1, R2));
326
327 // Insert it into the triangulation
328 240137 index_t v = CDTBase2d::insert(vertex_.size()-1);
329
330 // If it was an existing vertex, return the existing vertex
331
2/2
✓ Branch 2 taken 151847 times.
✓ Branch 3 taken 88290 times.
 240137 if(vertex_.size() > CDTBase2d::nv()) {
332 151847 vertex_.pop_back();
333 }
334 240137 return v;
335 }
336
337 127158 void MeshInTriangle::add_edge(
338 index_t f2,
339 TriangleRegion AR1, TriangleRegion AR2,
340 TriangleRegion BR1, TriangleRegion BR2
341 ) {
342 127158 index_t v1 = add_vertex(f2, AR1, AR2);
343 127158 index_t v2 = add_vertex(f2, BR1, BR2);
344
345 // If both extremities are on the same edge of f1,
346 // we do not add the edge, because it will be generated
347 // when remeshing the edge of f1
348
2/2
✓ Branch 2 taken 49028 times.
✓ Branch 3 taken 78130 times.
 127158 if(region_dim(regions_convex_hull(AR1,BR1)) == 1) {
349 49028 return;
350 }
351
352 // Generate also the combinatorial information of the edge,
353 // that indicates whether both extremities are on the same
354 // edge of f2 (useful later to compute the intersections)
355
2/4
✓ Branch 1 taken 78130 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 78130 times.
✗ Branch 6 not taken.
 78130 edges_.push_back(Edge(v1,v2,f2,regions_convex_hull(AR2,BR2)));
356
357 // Constraints will be added to the triangulation during commit()
358 }
359
360 15505 void MeshInTriangle::commit() {
361
362
2/2
✓ Branch 2 taken 124645 times.
✓ Branch 3 taken 15505 times.
 155655 for(const Edge& E: edges_) {
363
1/2
✓ Branch 1 taken 124645 times.
✗ Branch 2 not taken.
 124645 CDTBase2d::insert_constraint(E.v1, E.v2);
364 }
365
366
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 15505 times.
 15505 if(dry_run_) {
367 ✗ return;
368 }
369
370 // Protect global mesh from concurrent accesses
371 15505 exact_mesh_.lock();
372
373 // Create vertices and facets in target mesh
374
2/2
✓ Branch 1 taken 140489 times.
✓ Branch 2 taken 15505 times.
 155994 for(index_t i=0; i<vertex_.size(); ++i) {
375 // Vertex already exists in this MeshInTriangle
376
2/2
✓ Branch 1 taken 49252 times.
✓ Branch 2 taken 91237 times.
 140489 if(vertex_[i].mesh_vertex_index != NO_INDEX) {
377 49252 continue;
378 }
379 91237 vertex_[i].mesh_vertex_index =
380 91237 exact_mesh_.find_or_create_exact_vertex(
381 91237 vertex_[i].point_exact
382 );
383 }
384
385 // Create facets in target mesh
386
2/2
✓ Branch 1 taken 131387 times.
✓ Branch 2 taken 15505 times.
 146892 for(index_t t=0; t<CDTBase2d::nT(); ++t) {
387 131387 index_t i = CDTBase2d::Tv(t,0);
388 131387 index_t j = CDTBase2d::Tv(t,1);
389 131387 index_t k = CDTBase2d::Tv(t,2);
390 131387 i = vertex_[i].mesh_vertex_index;
391 131387 j = vertex_[j].mesh_vertex_index;
392 131387 k = vertex_[k].mesh_vertex_index;
393 131387 index_t new_t = target_mesh().facets.create_triangle(i,j,k);
394 // Copy all attributes from initial facet
395 131387 target_mesh().facets.attributes().copy_item(new_t, f1_);
396 }
397
398 // We are done with modification in the mesh
399 15505 exact_mesh_.unlock();
400 }
401
402 ✗ void MeshInTriangle::get_constraints(Mesh& M, bool with_edges) const {
403 ✗ if(M.vertices.nb() == 0) {
404 ✗ M.vertices.set_dimension(2);
405 ✗ for(index_t v=0; v<vertex_.size(); ++v) {
406 ✗ vec2 p = vertex_[v].get_UV_approx();
407 ✗ M.vertices.create_vertex(p.data());
408 }
409 }
410 ✗ if(with_edges && M.edges.nb() == 0) {
411 ✗ for(const Edge& E: edges_) {
412 ✗ M.edges.create_edge(E.v1, E.v2);
413 }
414 }
415 ✗ }
416
417 /**
418 * \brief Tests the parity of the permutation of a list of
419 * three distinct indices with respect to the canonical order.
420 */
421 5727897 static bool odd_order(index_t i, index_t j, index_t k) {
422 // Implementation: sort the elements (bubble sort is OK for
423 // such a small number), and invert parity each time
424 // two elements are swapped.
425 5727897 index_t tab[3] = { i, j, k};
426 5727897 const int N = 3;
427 5727897 bool result = false;
428
2/2
✓ Branch 0 taken 11455794 times.
✓ Branch 1 taken 5727897 times.
 17183691 for (int I = 0; I < N - 1; ++I) {
429
2/2
✓ Branch 0 taken 17183691 times.
✓ Branch 1 taken 11455794 times.
 28639485 for (int J = 0; J < N - I - 1; ++J) {
430
2/2
✓ Branch 0 taken 10341232 times.
✓ Branch 1 taken 6842459 times.
 17183691 if (tab[J] > tab[J + 1]) {
431 10341232 std::swap(tab[J], tab[J + 1]);
432 10341232 result = !result;
433 }
434 }
435 }
436 5727897 return result;
437 }
438
439 245821 void MeshInTriangle::begin_insert_transaction() {
440 245821 use_pred_cache_insert_buffer_ = true;
441 245821 }
442
443 88290 void MeshInTriangle::commit_insert_transaction() {
444
2/2
✓ Branch 2 taken 891974 times.
✓ Branch 3 taken 88290 times.
 1068554 for(const auto& it: pred_cache_insert_buffer_) {
445
1/2
✓ Branch 1 taken 891974 times.
✗ Branch 2 not taken.
 891974 pred_cache_[it.first] = it.second;
446 }
447 88290 pred_cache_insert_buffer_.resize(0);
448 88290 use_pred_cache_insert_buffer_ = false;
449 88290 }
450
451 157531 void MeshInTriangle::rollback_insert_transaction() {
452 157531 pred_cache_insert_buffer_.resize(0);
453 157531 use_pred_cache_insert_buffer_ = false;
454 157531 }
455
456 5727897 Sign MeshInTriangle::orient2d(index_t vx1,index_t vx2,index_t vx3) const {
457
458
1/2
✓ Branch 1 taken 5727897 times.
✗ Branch 2 not taken.
 5727897 trindex K(vx1, vx2, vx3);
459
460
2/2
✓ Branch 0 taken 3631024 times.
✓ Branch 1 taken 2096873 times.
 5727897 if(use_pred_cache_insert_buffer_) {
461 10893072 Sign result = PCK::orient_2d_projected(
462
1/2
✓ Branch 1 taken 3631024 times.
✗ Branch 2 not taken.
 3631024 vertex_[K.indices[0]].point_exact,
463
1/2
✓ Branch 1 taken 3631024 times.
✗ Branch 2 not taken.
 3631024 vertex_[K.indices[1]].point_exact,
464
1/2
✓ Branch 1 taken 3631024 times.
✗ Branch 2 not taken.
 3631024 vertex_[K.indices[2]].point_exact,
465
1/2
✓ Branch 1 taken 3631024 times.
✗ Branch 2 not taken.
 3631024 f1_normal_axis_
466 3631024 );
467
2/4
✓ Branch 1 taken 3631024 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3631024 times.
✗ Branch 5 not taken.
 3631024 pred_cache_insert_buffer_.push_back(std::make_pair(K, result));
468
2/2
✓ Branch 1 taken 1788150 times.
✓ Branch 2 taken 1842874 times.
 3631024 if(odd_order(vx1,vx2,vx3)) {
469 1788150 result = Sign(-result);
470 }
471 3631024 return result;
472 }
473
474 bool inserted;
475 2096873 std::map<trindex, Sign>::iterator it;
476
2/4
✓ Branch 1 taken 2096873 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2096873 times.
✗ Branch 5 not taken.
 2096873 std::tie(it,inserted) = pred_cache_.insert(std::make_pair(K,ZERO));
477 Sign result;
478
479
2/2
✓ Branch 0 taken 690418 times.
✓ Branch 1 taken 1406455 times.
 2096873 if(inserted) {
480 2071254 result = PCK::orient_2d_projected(
481
1/2
✓ Branch 1 taken 690418 times.
✗ Branch 2 not taken.
 690418 vertex_[K.indices[0]].point_exact,
482
1/2
✓ Branch 1 taken 690418 times.
✗ Branch 2 not taken.
 690418 vertex_[K.indices[1]].point_exact,
483
1/2
✓ Branch 1 taken 690418 times.
✗ Branch 2 not taken.
 690418 vertex_[K.indices[2]].point_exact,
484
1/2
✓ Branch 1 taken 690418 times.
✗ Branch 2 not taken.
 690418 f1_normal_axis_
485 );
486 690418 it->second = result;
487 } else {
488 1406455 result = it->second;
489 }
490
491
2/2
✓ Branch 1 taken 961110 times.
✓ Branch 2 taken 1135763 times.
 2096873 if(odd_order(vx1,vx2,vx3)) {
492 961110 result = Sign(-result);
493 }
494
495 2096873 return result;
496 }
497
498 186174 Sign MeshInTriangle::incircle(
499 index_t v1,index_t v2,index_t v3,index_t v4
500 ) const {
501 exact::vec2h p1(
502 372348 vertex_[v1].point_exact[u_],
503 372348 vertex_[v1].point_exact[v_],
504
2/4
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
 186174 vertex_[v1].point_exact.w
505
4/8
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 186174 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 186174 times.
✗ Branch 11 not taken.
 186174 );
506 exact::vec2h p2(
507 372348 vertex_[v2].point_exact[u_],
508 372348 vertex_[v2].point_exact[v_],
509
2/4
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
 186174 vertex_[v2].point_exact.w
510
4/8
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 186174 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 186174 times.
✗ Branch 11 not taken.
 186174 );
511 exact::vec2h p3(
512 372348 vertex_[v3].point_exact[u_],
513 372348 vertex_[v3].point_exact[v_],
514
2/4
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
 186174 vertex_[v3].point_exact.w
515
4/8
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 186174 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 186174 times.
✗ Branch 11 not taken.
 186174 );
516 exact::vec2h p4(
517 372348 vertex_[v4].point_exact[u_],
518 372348 vertex_[v4].point_exact[v_],
519
2/4
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
 186174 vertex_[v4].point_exact.w
520
4/8
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 186174 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 186174 times.
✗ Branch 11 not taken.
 186174 );
521 #ifdef GEOGRAM_USE_EXACT_NT
522 return PCK::incircle_2d_SOS(p1,p2,p3,p4);
523 #else
524 372348 return PCK::incircle_2d_SOS_with_lengths(
525 p1,p2,p3,p4,
526
1/2
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
 186174 vertex_[v1].l,
527
1/2
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
 186174 vertex_[v2].l,
528
1/2
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
 186174 vertex_[v3].l,
529
2/4
✓ Branch 1 taken 186174 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 186174 times.
✗ Branch 5 not taken.
 186174 vertex_[v4].l
530 372348 );
531 #endif
532 186174 }
533
534 5684 index_t MeshInTriangle::create_intersection(
535 index_t e1, index_t i, index_t j,
536 index_t e2, index_t k, index_t l
537 ) {
538 5684 geo_argused(i);
539 5684 geo_argused(j);
540 5684 geo_argused(k);
541 5684 geo_argused(l);
542
543 5684 ExactPoint I;
544
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 get_edge_edge_intersection(e1,e2,I);
545
2/4
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5684 times.
✗ Branch 5 not taken.
 5684 vertex_.push_back(Vertex(this,I));
546 5684 index_t x = vertex_.size()-1;
547
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 CDTBase2d::v2T_.push_back(NO_INDEX);
548
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 5684 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 5684 geo_debug_assert(x == CDTBase2d::nv_);
549 5684 ++CDTBase2d::nv_;
550 5684 return x;
551 5684 }
552
553 5684 void MeshInTriangle::get_edge_edge_intersection(
554 index_t e1, index_t e2, ExactPoint& I
555 ) const {
556 5684 index_t f1 = f1_;
557
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 index_t f2 = edges_[e1].sym.f2;
558
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 index_t f3 = edges_[e2].sym.f2;
559
560
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 5684 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 5684 geo_debug_assert(f1 != NO_INDEX);
561
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 5684 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 5684 geo_debug_assert(f2 != NO_INDEX);
562
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 5684 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 5684 geo_debug_assert(f3 != NO_INDEX);
563
564 vec3 P[9] = {
565 11368 mesh_facet_vertex(f1,0), mesh_facet_vertex(f1,1),
566 5684 mesh_facet_vertex(f1,2),
567 11368 mesh_facet_vertex(f2,0), mesh_facet_vertex(f2,1),
568 5684 mesh_facet_vertex(f2,2),
569 11368 mesh_facet_vertex(f3,0), mesh_facet_vertex(f3,1),
570 5684 mesh_facet_vertex(f3,2)
571
9/18
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5684 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 5684 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 5684 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 5684 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 5684 times.
✗ Branch 17 not taken.
✓ Branch 19 taken 5684 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 5684 times.
✗ Branch 23 not taken.
✓ Branch 25 taken 5684 times.
✗ Branch 26 not taken.
 5684 };
572
573
3/4
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4343 times.
✓ Branch 4 taken 1341 times.
 5684 if(!get_three_planes_intersection(
574 I,
575 P[0], P[1], P[2],
576 P[3], P[4], P[5],
577 P[6], P[7], P[8]
578 )) {
579
1/2
✓ Branch 1 taken 4343 times.
✗ Branch 2 not taken.
 4343 get_edge_edge_intersection_2D(e1,e2,I);
580 4343 return;
581 }
582 }
583
584 4343 void MeshInTriangle::get_edge_edge_intersection_2D(
585 index_t e1, index_t e2, ExactPoint& I
586 ) const {
587 4343 const Edge& E1 = edges_[e1];
588 4343 const Edge& E2 = edges_[e2];
589
590 4343 if(
591
4/4
✓ Branch 1 taken 4299 times.
✓ Branch 2 taken 44 times.
✓ Branch 3 taken 4287 times.
✓ Branch 4 taken 56 times.
 8642 region_dim(E1.sym.R2) == 1 &&
592
2/2
✓ Branch 1 taken 4287 times.
✓ Branch 2 taken 12 times.
 4299 region_dim(E2.sym.R2) == 1
593 ) {
594 4287 index_t le1 = index_t(E1.sym.R2)-index_t(T2_RGN_E0);
595 4287 index_t le2 = index_t(E2.sym.R2)-index_t(T2_RGN_E0);
596
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 4287 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 4287 geo_debug_assert(le1 < 3);
597
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 4287 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 4287 geo_debug_assert(le2 < 3);
598
599
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 vec2 p1_uv = mesh_facet_vertex_UV(E1.sym.f2, (le1+1)%3);
600
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 vec2 p2_uv = mesh_facet_vertex_UV(E1.sym.f2, (le1+2)%3);
601
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 vec2 q1_uv = mesh_facet_vertex_UV(E2.sym.f2, (le2+1)%3);
602
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 vec2 q2_uv = mesh_facet_vertex_UV(E2.sym.f2, (le2+2)%3);
603
604
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 exact::vec2 C1 = make_vec2<exact::vec2>(p1_uv, p2_uv);
605
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 exact::vec2 C2 = make_vec2<exact::vec2>(q2_uv, q1_uv);
606
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 exact::vec2 B = make_vec2<exact::vec2>(p1_uv, q1_uv);
607
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 exact::scalar d = det(C1,C2);
608
2/8
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 4287 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 4287 geo_debug_assert(d.sign() != ZERO);
609
2/4
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4287 times.
✗ Branch 5 not taken.
 4287 exact::rational t(det(B,C2),d);
610
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 8574 I = mix(
611 t,
612
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 mesh_facet_vertex(E1.sym.f2,(le1+1)%3),
613
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 mesh_facet_vertex(E1.sym.f2,(le1+2)%3)
614
1/2
✓ Branch 1 taken 4287 times.
✗ Branch 2 not taken.
 4287 );
615
616 4287 } else {
617
6/14
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 44 times.
✓ Branch 4 taken 12 times.
✓ Branch 6 taken 44 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 44 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 56 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
 56 geo_debug_assert(
618 region_dim(E1.sym.R2) == 1 || region_dim(E2.sym.R2) == 1
619 );
620 56 index_t f1 = E1.sym.f2;
621 56 TriangleRegion R1 = E1.sym.R2;
622 56 index_t f2 = E2.sym.f2;
623 56 TriangleRegion R2 = E2.sym.R2;
624
3/4
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 12 times.
✓ Branch 4 taken 44 times.
 56 if(region_dim(R1) == 1) {
625 12 std::swap(f1,f2);
626 12 std::swap(R1,R2);
627 }
628
629 56 index_t e = index_t(R2) - index_t(T2_RGN_E0);
630
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 56 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 56 geo_debug_assert(e < 3);
631
632
1/2
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
 112 I = plane_line_intersection(
633
1/2
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
 56 mesh_facet_vertex(f1,0),
634
1/2
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
 56 mesh_facet_vertex(f1,1),
635
1/2
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
 56 mesh_facet_vertex(f1,2),
636
1/2
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
 56 mesh_facet_vertex(f2,(e+1)%3),
637
1/2
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
 112 mesh_facet_vertex(f2,(e+2)%3)
638
1/2
✓ Branch 1 taken 56 times.
✗ Branch 2 not taken.
 56 );
639 }
640 4343 }
641
642 ✗ void MeshInTriangle::save(const std::string& filename) const {
643 ✗ Mesh M;
644 ✗ M.vertices.set_dimension(2);
645 ✗ for(index_t v=0; v<CDTBase2d::nv(); ++v) {
646 ✗ vec2 p = vertex_[v].get_UV_approx();
647 ✗ M.vertices.create_vertex(p.data());
648 }
649 ✗ for(index_t t=0; t<CDTBase2d::nT(); ++t) {
650 ✗ M.facets.create_triangle(
651 CDTBase2d::Tv(t,0),
652 CDTBase2d::Tv(t,1),
653 CDTBase2d::Tv(t,2)
654 );
655 }
656
657 ✗ Attribute<double> tex_coord;
658 ✗ tex_coord.create_vector_attribute(
659 ✗ M.facet_corners.attributes(), "tex_coord", 2
660 );
661 static double triangle_tex[3][2] = {
662 {0.0, 0.0},
663 {1.0, 0.0},
664 {0.0, 1.0}
665 };
666 ✗ for(index_t c: M.facet_corners) {
667 ✗ tex_coord[2*c] = triangle_tex[c%3][0];
668 ✗ tex_coord[2*c+1] = triangle_tex[c%3][1];
669 }
670 ✗ mesh_save(M, filename);
671 ✗ }
672
673 /**************************************************************************/
674
675 285 CoplanarFacets::CoplanarFacets(
676 MeshSurfaceIntersection& I, bool clear_attributes,
677 double angle_tolerance
678 285 ) :
679 285 I_(I),
680 285 mesh_(I.target_mesh()),
681 285 mesh_copy_(I.readonly_mesh()),
682 285 angle_tolerance_(angle_tolerance),
683
2/4
✓ Branch 1 taken 285 times.
✗ Branch 2 not taken.
✓ Branch 6 taken 285 times.
✗ Branch 7 not taken.
 570 facet_group_(I.target_mesh().facets.attributes(),"group"),
684
2/4
✓ Branch 1 taken 285 times.
✗ Branch 2 not taken.
✓ Branch 6 taken 285 times.
✗ Branch 7 not taken.
 570 keep_vertex_(I.target_mesh().vertices.attributes(),"keep"),
685
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 570 c_is_coplanar_(
686
1/2
✓ Branch 1 taken 285 times.
✗ Branch 2 not taken.
 570 I.target_mesh().facet_corners.attributes(),"is_coplanar"
687 ),
688
2/4
✓ Branch 1 taken 285 times.
✗ Branch 2 not taken.
✓ Branch 6 taken 285 times.
✗ Branch 7 not taken.
 570 f_is_flipped_(I.target_mesh().facets.attributes(),"flipped"),
689
1/2
✓ Branch 1 taken 285 times.
✗ Branch 2 not taken.
 285 halfedges_(*this),
690
1/2
✓ Branch 8 taken 285 times.
✗ Branch 9 not taken.
 855 polylines_(*this)
691 {
692
2/2
✓ Branch 0 taken 57 times.
✓ Branch 1 taken 228 times.
 285 if(clear_attributes) {
693
2/2
✓ Branch 4 taken 123492 times.
✓ Branch 5 taken 57 times.
 123549 for(index_t f: mesh_.facets) {
694
1/2
✓ Branch 1 taken 123492 times.
✗ Branch 2 not taken.
 123492 facet_group_[f] = NO_INDEX;
695 }
696
2/2
✓ Branch 4 taken 61752 times.
✓ Branch 5 taken 57 times.
 61809 for(index_t v: mesh_.vertices) {
697
1/2
✓ Branch 2 taken 61752 times.
✗ Branch 3 not taken.
 61752 keep_vertex_[v] = false;
698 }
699
1/2
✓ Branch 1 taken 57 times.
✗ Branch 2 not taken.
 57 find_coplanar_facets();
700 }
701
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 f_visited_.assign(mesh_.facets.nb(),false);
702
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 h_visited_.assign(mesh_.facet_corners.nb(),false);
703
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 v_visited_.assign(mesh_.vertices.nb(),false);
704
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 v_idx_.assign(mesh_.vertices.nb(),NO_INDEX);
705 285 }
706
707 57 void CoplanarFacets::find_coplanar_facets() {
708
709 // Positioned by MeshSurfaceIntersection::build_Weiler_model()
710 Attribute<bool> corner_is_on_border(
711
1/2
✓ Branch 1 taken 57 times.
✗ Branch 2 not taken.
 114 mesh_.facet_corners.attributes(), "is_on_border"
712
1/2
✓ Branch 2 taken 57 times.
✗ Branch 3 not taken.
 57 );
713
714
2/2
✓ Branch 4 taken 370476 times.
✓ Branch 5 taken 57 times.
 370533 for(index_t c: mesh_.facet_corners) {
715
1/2
✓ Branch 2 taken 370476 times.
✗ Branch 3 not taken.
 370476 c_is_coplanar_[c] = false;
716 }
717
718 // TODO: when there is an angle tolerance, one should check instead
719 // angle deviation w.r.t. a single seed facet per facet group, because
720 // with the present algorithm, if a large number of tiny facets are
721 // connected (e.g. highly tessellated cylinder), one may group facets
722 // with large angle deviation (without seeing it because each facet has
723 // small angle deviation w.r.t. its neighbors).
724
725
1/2
✓ Branch 1 taken 57 times.
✗ Branch 2 not taken.
 57 parallel_for(
726 57 0, mesh_.facet_corners.nb(),
727 114 [&](index_t c1) {
728 370476 index_t f1 = (c1 / 3);
729 370476 index_t le1 = (c1 % 3);
730
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t f2 = mesh_.facet_corners.adjacent_facet(c1);
731
732
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 370476 times.
 370476 if(f2 == NO_INDEX) {
733 216682 return;
734 }
735
736 // do not traverse true borders
737
2/4
✓ Branch 2 taken 370476 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 370476 times.
 370476 if(corner_is_on_border[c1]) {
738 ✗ return;
739 }
740
741
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t v11 = mesh_.facets.vertex(f1,le1);
742
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t v12 = mesh_.facets.vertex(f1,(le1+1)%3);
743
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t le2 = mesh_.facets.find_edge(f2,v12,v11);
744
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t c2 = mesh_.facets.corner(f2,le2);
745
746 #ifdef GEO_DEBUG
747
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t v21 = mesh_.facets.vertex(f2,le2);
748
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t v22 = mesh_.facets.vertex(f2,(le2+1)%3);
749
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t v13 = mesh_.facets.vertex(f1,(le1+2)%3);
750
1/2
✓ Branch 1 taken 370476 times.
✗ Branch 2 not taken.
 370476 index_t v23 = mesh_.facets.vertex(f2,(le2+2)%3);
751
752
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 370476 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 370476 geo_debug_assert(v11 == v22);
753
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 370476 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 370476 geo_debug_assert(v12 == v21);
754
3/10
✓ Branch 0 taken 370476 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 370476 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 370476 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 370476 geo_debug_assert(v11!=v12 && v12!=v13 && v13!=v11);
755
3/10
✓ Branch 0 taken 370476 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 370476 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 370476 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 370476 geo_debug_assert(v21!=v22 && v22!=v23 && v23!=v21);
756 #endif
757
758
2/2
✓ Branch 0 taken 185238 times.
✓ Branch 1 taken 185238 times.
 370476 if(c1 > c2) {
759 185238 return;
760 }
761
762 // Small optimization: if both triangles come from same
763 // original facet then they are coplanar
764
4/6
✓ Branch 1 taken 185238 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 185238 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 31444 times.
✓ Branch 7 taken 153794 times.
 185238 if(I_.get_initial_facet(f1) == I_.get_initial_facet(f2)) {
765
1/2
✓ Branch 2 taken 31444 times.
✗ Branch 3 not taken.
 31444 c_is_coplanar_[c1] = true;
766
1/2
✓ Branch 2 taken 31444 times.
✗ Branch 3 not taken.
 31444 c_is_coplanar_[c2] = true;
767 31444 return;
768 }
769
770 // Use original triangles for co-planarity test
771
1/2
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
 153794 auto [p1, p2, p3] = I_.get_initial_facet_vertices(f1);
772
1/2
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
 153794 auto [q1, q2, q3] = I_.get_initial_facet_vertices(f2);
773
3/4
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 30338 times.
✓ Branch 4 taken 123456 times.
 153794 if(triangles_are_coplanar(p1,p2,p3,q1,q2,q3)) {
774
1/2
✓ Branch 2 taken 30338 times.
✗ Branch 3 not taken.
 30338 c_is_coplanar_[c1] = true;
775
1/2
✓ Branch 2 taken 30338 times.
✗ Branch 3 not taken.
 30338 c_is_coplanar_[c2] = true;
776 }
777 }
778 );
779 57 }
780
781 146168 void CoplanarFacets::get(index_t f, index_t group_id) {
782
783 146168 facets_.resize(0);
784 146168 vertices_.resize(0);
785 146168 halfedges_.initialize();
786 146168 polylines_.initialize();
787
788 // Get facets
789 {
790
1/2
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
 146168 std::stack<index_t> S;
791
1/2
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
 146168 facet_group_[f] = group_id;
792 146168 f_visited_[f] = true;
793
1/2
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
 146168 S.push(f);
794
1/2
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
 146168 facets_.push_back(f);
795
2/2
✓ Branch 1 taken 246984 times.
✓ Branch 2 taken 146168 times.
 393152 while(!S.empty()) {
796 246984 index_t f1 = S.top();
797 246984 S.pop();
798
2/2
✓ Branch 0 taken 740952 times.
✓ Branch 1 taken 246984 times.
 987936 for(index_t le1=0; le1<3; ++le1) {
799
1/2
✓ Branch 1 taken 740952 times.
✗ Branch 2 not taken.
 740952 index_t f2 = mesh_.facets.adjacent(f1,le1);
800 740952 if(
801
3/4
✓ Branch 0 taken 740952 times.
✗ Branch 1 not taken.
✓ Branch 4 taken 594640 times.
✓ Branch 5 taken 146312 times.
 1335592 f2 != NO_INDEX && !f_visited_[f2] &&
802
6/8
✓ Branch 1 taken 594640 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 594640 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 100816 times.
✓ Branch 8 taken 493824 times.
✓ Branch 9 taken 100816 times.
✓ Branch 10 taken 640136 times.
 1335592 c_is_coplanar_[mesh_.facets.corner(f1,le1)]
803 ) {
804
2/4
✓ Branch 1 taken 100816 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 100816 times.
✗ Branch 5 not taken.
 100816 facet_group_[f2] = facet_group_[f1];
805 100816 f_visited_[f2] = true;
806
1/2
✓ Branch 1 taken 100816 times.
✗ Branch 2 not taken.
 100816 S.push(f2);
807
1/2
✓ Branch 1 taken 100816 times.
✗ Branch 2 not taken.
 100816 facets_.push_back(f2);
808 }
809 }
810 }
811
2/2
✓ Branch 2 taken 246984 times.
✓ Branch 3 taken 146168 times.
 539320 for(index_t cur_f: facets_) {
812 246984 f_visited_[cur_f] = false;
813 }
814 146168 }
815
816 146168 group_id_ = group_id;
817
818 // Initialize projection coordinates, using original facet
819 {
820
2/4
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 146168 times.
✗ Branch 5 not taken.
 146168 auto [p1, p2, p3] = I_.get_initial_facet_vertices(facets_[0]);
821
1/2
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
 146168 coord_index_t projection_axis = PCK::triangle_normal_axis(p1,p2,p3);
822 146168 u_ = coord_index_t((projection_axis+1)%3);
823 146168 v_ = coord_index_t((projection_axis+2)%3);
824
1/2
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
 146168 Sign o = PCK::orient_2d(
825
6/12
✓ Branch 1 taken 146168 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 146168 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 146168 times.
✗ Branch 9 not taken.
✓ Branch 11 taken 146168 times.
✗ Branch 12 not taken.
✓ Branch 15 taken 146168 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 146168 times.
✗ Branch 19 not taken.
 146168 vec2(p1[u_],p1[v_]), vec2(p2[u_],p2[v_]), vec2(p3[u_],p3[v_])
826 );
827
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 146168 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 146168 geo_debug_assert(o != ZERO);
828
2/2
✓ Branch 0 taken 71084 times.
✓ Branch 1 taken 75084 times.
 146168 if(o < 0) {
829 71084 std::swap(u_,v_);
830 }
831 }
832
833 // Get vertices and halfedges
834 {
835
2/2
✓ Branch 2 taken 246984 times.
✓ Branch 3 taken 146168 times.
 539320 for(index_t f1: facets_) {
836
2/2
✓ Branch 0 taken 740952 times.
✓ Branch 1 taken 246984 times.
 987936 for(index_t le=0; le<3; ++le) {
837
1/2
✓ Branch 1 taken 740952 times.
✗ Branch 2 not taken.
 740952 index_t f2 = mesh_.facets.adjacent(f1,le);
838 740952 if(
839
1/2
✓ Branch 0 taken 740952 times.
✗ Branch 1 not taken.
 1481904 f2 == NO_INDEX ||
840
6/8
✓ Branch 1 taken 740952 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 740952 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 493824 times.
✓ Branch 8 taken 247128 times.
✓ Branch 9 taken 493824 times.
✓ Branch 10 taken 247128 times.
 1481904 !c_is_coplanar_[mesh_.facets.corner(f1,le)]
841 ) {
842
2/4
✓ Branch 1 taken 493824 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 493824 times.
✗ Branch 5 not taken.
 493824 halfedges_.add(mesh_.facets.corners_begin(f1)+le);
843
1/2
✓ Branch 1 taken 493824 times.
✗ Branch 2 not taken.
 493824 index_t v1 = mesh_.facets.vertex(f1,le);
844
1/2
✓ Branch 1 taken 493824 times.
✗ Branch 2 not taken.
 493824 index_t v2 = mesh_.facets.vertex(f1,(le+1)%3);
845
2/2
✓ Branch 2 taken 175216 times.
✓ Branch 3 taken 318608 times.
 493824 if(!v_visited_[v1]) {
846
1/2
✓ Branch 2 taken 175216 times.
✗ Branch 3 not taken.
 175216 v_idx_[v1] = vertices_.size();
847
1/2
✓ Branch 1 taken 175216 times.
✗ Branch 2 not taken.
 175216 vertices_.push_back(v1);
848 175216 v_visited_[v1] = true;
849 }
850
2/2
✓ Branch 2 taken 307112 times.
✓ Branch 3 taken 186712 times.
 493824 if(!v_visited_[v2]) {
851
1/2
✓ Branch 2 taken 307112 times.
✗ Branch 3 not taken.
 307112 v_idx_[v2] = vertices_.size();
852
1/2
✓ Branch 1 taken 307112 times.
✗ Branch 2 not taken.
 307112 vertices_.push_back(v2);
853 307112 v_visited_[v2] = true;
854 }
855 } else {
856 // This one for the particular case of a non-manifold
857 // vertex, such as a cone apex touching a facet
858 // (ThingiCSG/Basic/cube_cone_1.scad)
859
1/2
✓ Branch 1 taken 247128 times.
✗ Branch 2 not taken.
 247128 index_t v = mesh_.facets.vertex(f1,le);
860
7/8
✓ Branch 2 taken 247128 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 71279 times.
✓ Branch 5 taken 175849 times.
✓ Branch 8 taken 11492 times.
✓ Branch 9 taken 59787 times.
✓ Branch 10 taken 11492 times.
✓ Branch 11 taken 235636 times.
 247128 if(keep_vertex_[v] && !v_visited_[v]) {
861
1/2
✓ Branch 1 taken 11492 times.
✗ Branch 2 not taken.
 11492 vertices_.push_back(v);
862 11492 v_visited_[v] = true;
863 }
864 }
865 }
866 }
867
2/2
✓ Branch 2 taken 493820 times.
✓ Branch 3 taken 146168 times.
 786156 for(index_t v: vertices_) {
868 493820 v_visited_[v] = false;
869 }
870 }
871
872 // Get polylines
873 {
874 // Get all polylines starting from vertices with more than
875 // 2 incident halfedges.
876
2/2
✓ Branch 2 taken 493820 times.
✓ Branch 3 taken 146168 times.
 786156 for(index_t v: vertices_) {
877
3/4
✓ Branch 1 taken 493820 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 493816 times.
 493820 if(halfedges_.nb_halfedges_around_vertex(v) > 1) {
878 4 for(
879
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 index_t h=halfedges_.vertex_first_halfedge(v);
880
3/4
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 4 times.
 12 h != NO_INDEX; h = halfedges_.next_around_vertex(h)
881 ) {
882
1/2
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
 8 if(!h_visited_[h]) {
883 8 polylines_.begin_polyline();
884 8 index_t h2 = h;
885 do {
886
1/6
✗ Branch 2 not taken.
✓ Branch 3 taken 76 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 76 geo_debug_assert(!h_visited_[h2]);
887 76 h_visited_[h2] = true;
888
1/2
✓ Branch 1 taken 76 times.
✗ Branch 2 not taken.
 76 polylines_.add_halfedge(h2);
889
1/2
✓ Branch 1 taken 76 times.
✗ Branch 2 not taken.
 76 h2 = halfedges_.next_along_polyline(h2);
890
3/4
✓ Branch 0 taken 68 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 68 times.
✗ Branch 3 not taken.
 76 } while(h2 != NO_INDEX && h2 != h);
891
1/2
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
 8 polylines_.end_polyline();
892 }
893 }
894 }
895 }
896 // There can be also closed halfedge loops with no irregular vertex
897
2/2
✓ Branch 2 taken 493824 times.
✓ Branch 3 taken 146168 times.
 786160 for(index_t h: halfedges_) {
898
2/2
✓ Branch 2 taken 146366 times.
✓ Branch 3 taken 347458 times.
 493824 if(!h_visited_[h]) {
899 146366 polylines_.begin_polyline();
900 146366 index_t h2 = h;
901 do {
902
1/6
✗ Branch 2 not taken.
✓ Branch 3 taken 493748 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 493748 geo_debug_assert(!h_visited_[h2]);
903 493748 h_visited_[h2] = true;
904
1/2
✓ Branch 1 taken 493748 times.
✗ Branch 2 not taken.
 493748 polylines_.add_halfedge(h2);
905
1/2
✓ Branch 1 taken 493748 times.
✗ Branch 2 not taken.
 493748 h2 = halfedges_.next_along_polyline(h2);
906
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 493748 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 493748 geo_debug_assert(h2 != NO_INDEX);
907
2/2
✓ Branch 0 taken 347382 times.
✓ Branch 1 taken 146366 times.
 493748 } while(h2 != h);
908
1/2
✓ Branch 1 taken 146366 times.
✗ Branch 2 not taken.
 146366 polylines_.end_polyline();
909 }
910 }
911
2/2
✓ Branch 2 taken 493824 times.
✓ Branch 3 taken 146168 times.
 786160 for(index_t h: halfedges_) {
912 493824 h_visited_[h] = false;
913 }
914 }
915 146168 }
916
917 73084 void CoplanarFacets::mark_vertices_to_keep() {
918
3/4
✓ Branch 2 taken 73084 times.
✗ Branch 3 not taken.
✓ Branch 7 taken 73187 times.
✓ Branch 8 taken 73084 times.
 146271 for(index_t P: polylines_) {
919
1/2
✓ Branch 1 taken 73187 times.
✗ Branch 2 not taken.
 73187 index_t first_v = polylines_.first_vertex(P);
920
1/2
✓ Branch 1 taken 73187 times.
✗ Branch 2 not taken.
 73187 index_t last_v = polylines_.last_vertex(P);
921
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 73187 times.
 73187 if(first_v != last_v) {
922 ✗ keep_vertex_[first_v] = true;
923 ✗ keep_vertex_[last_v] = true;
924 }
925
1/2
✓ Branch 1 taken 73187 times.
✗ Branch 2 not taken.
 73187 index_t v1 = polylines_.prev_first_vertex(P);
926 73187 index_t v2 = NO_INDEX;
927 73187 index_t v3 = NO_INDEX;
928
4/6
✓ Branch 1 taken 73187 times.
✗ Branch 2 not taken.
✓ Branch 6 taken 246912 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 246912 times.
✓ Branch 11 taken 73187 times.
 320099 for(index_t h: polylines_.halfedges(P)) {
929
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 246912 times.
 246912 if(v1 == NO_INDEX) {
930 ✗ continue;
931 }
932
1/2
✓ Branch 1 taken 246912 times.
✗ Branch 2 not taken.
 246912 v2 = halfedges_.vertex(h,0);
933
1/2
✓ Branch 1 taken 246912 times.
✗ Branch 2 not taken.
 246912 v3 = halfedges_.vertex(h,1);
934
1/2
✓ Branch 1 taken 246912 times.
✗ Branch 2 not taken.
 246912 ExactPoint p1 = I_.exact_vertex(v1);
935
1/2
✓ Branch 1 taken 246912 times.
✗ Branch 2 not taken.
 246912 ExactPoint p2 = I_.exact_vertex(v2);
936
1/2
✓ Branch 1 taken 246912 times.
✗ Branch 2 not taken.
 246912 ExactPoint p3 = I_.exact_vertex(v3);
937
3/4
✓ Branch 1 taken 246912 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 238467 times.
✓ Branch 4 taken 8445 times.
 246912 if(!edges_are_colinear(p1,p2,p3)) {
938
1/2
✓ Branch 2 taken 238467 times.
✗ Branch 3 not taken.
 238467 keep_vertex_[v2] = true;
939 }
940 246912 v1 = v2;
941 246912 }
942 }
943 73084 }
944
945 ✗ void CoplanarFacets::save_borders(const std::string& filename) {
946 ✗ Mesh borders;
947 ✗ borders.vertices.set_dimension(2);
948 ✗ index_t cur_idx = 0;
949 ✗ for(index_t v: vertices_) {
950 ✗ vec3 p = mesh_.vertices.point(v);
951 ✗ vec2 q(p[u_],p[v_]);
952 ✗ borders.vertices.create_vertex(q.data());
953 ✗ v_idx_[v] = cur_idx;
954 ✗ ++cur_idx;
955 }
956
957 ✗ for(index_t h: halfedges_) {
958 ✗ index_t v1 = halfedges_.vertex(h,0);
959 ✗ index_t v2 = halfedges_.vertex(h,1);
960 ✗ v1 = v_idx_[v1];
961 ✗ v2 = v_idx_[v2];
962 ✗ geo_debug_assert(v1 != NO_INDEX);
963 ✗ geo_debug_assert(v2 != NO_INDEX);
964 ✗ borders.edges.create_edge(v1,v2);
965 }
966
967 ✗ Attribute<bool> selection(borders.vertices.attributes(), "selection");
968 ✗ for(index_t v: vertices_) {
969 ✗ geo_debug_assert(v_idx_[v] != NO_INDEX);
970 ✗ selection[v_idx_[v]] = keep_vertex_[v];
971 }
972 ✗ mesh_save(borders,filename);
973 ✗ }
974
975 ✗ void CoplanarFacets::save_facet_group(const std::string& filename) {
976 ✗ Mesh M;
977 ✗ Attribute<bool> keep_vertex(M.vertices.attributes(),"keep");
978 ✗ M.vertices.set_dimension(2);
979 ✗ for(index_t f: facets_) {
980 ✗ for(index_t lv=0; lv<3; ++lv) {
981 ✗ index_t v = mesh_.facets.vertex(f,lv);
982 ✗ v_idx_[v] = NO_INDEX;
983 }
984 }
985 ✗ for(index_t f: facets_) {
986 ✗ for(index_t lv=0; lv<3; ++lv) {
987 ✗ index_t v = mesh_.facets.vertex(f,lv);
988 ✗ if(v_idx_[v] == NO_INDEX) {
989 ✗ vec3 p = mesh_.vertices.point(v);
990 ✗ vec2 q(p[u_], p[v_]);
991 ✗ v_idx_[v] = M.vertices.create_vertex(q.data());
992 ✗ keep_vertex[v_idx_[v]] = keep_vertex_[v];
993 }
994 }
995 ✗ M.facets.create_triangle(
996 ✗ v_idx_[mesh_.facets.vertex(f,0)],
997 ✗ v_idx_[mesh_.facets.vertex(f,1)],
998 ✗ v_idx_[mesh_.facets.vertex(f,2)]
999 );
1000 }
1001
1002 ✗ for(index_t f: facets_) {
1003 ✗ for(index_t lv=0; lv<3; ++lv) {
1004 ✗ index_t v = mesh_.facets.vertex(f,lv);
1005 ✗ v_idx_[v] = NO_INDEX;
1006 }
1007 }
1008
1009 ✗ M.facets.connect();
1010 ✗ mesh_save(M,filename);
1011 ✗ }
1012
1013 8067 void CoplanarFacets::triangulate() {
1014
1015 // Compute 2D projected BBOX
1016 8067 double umin = Numeric::max_float64();
1017 8067 double vmin = Numeric::max_float64();
1018 8067 double umax = -Numeric::max_float64();
1019 8067 double vmax = -Numeric::max_float64();
1020
2/2
✓ Branch 2 taken 58475 times.
✓ Branch 3 taken 8067 times.
 74609 for(index_t f: facets_) {
1021
2/2
✓ Branch 0 taken 175425 times.
✓ Branch 1 taken 58475 times.
 233900 for(index_t lv=0; lv<3; ++lv) {
1022
1/2
✓ Branch 1 taken 175425 times.
✗ Branch 2 not taken.
 175425 index_t vx = mesh_.facets.vertex(f,lv);
1023
2/4
✓ Branch 1 taken 175425 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 175425 times.
✗ Branch 5 not taken.
 175425 double u = mesh_.vertices.point(vx)[u_];
1024
2/4
✓ Branch 1 taken 175425 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 175425 times.
✗ Branch 5 not taken.
 175425 double v = mesh_.vertices.point(vx)[v_];
1025 175425 umin = std::min(umin, u);
1026 175425 umax = std::max(umax, u);
1027 175425 vmin = std::min(vmin, v);
1028 175425 vmax = std::max(vmax, v);
1029 }
1030 }
1031 8067 double d = std::max(umax-umin, vmax-vmin);
1032 8067 d *= 10.0;
1033 8067 d = std::max(d, 1.0);
1034 8067 umin-=d; vmin-=d; umax+=d; vmax+=d;
1035
1036 // Create CDT
1037
1/2
✓ Branch 1 taken 8067 times.
✗ Branch 2 not taken.
 8067 CDT.clear();
1038
1/2
✓ Branch 1 taken 8067 times.
✗ Branch 2 not taken.
 8067 CDT.create_enclosing_rectangle(umin, vmin, umax, vmax);
1039
1040
2/2
✓ Branch 2 taken 51859 times.
✓ Branch 3 taken 8067 times.
 67993 for(index_t v: vertices_) {
1041
3/4
✓ Branch 2 taken 51859 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 43473 times.
✓ Branch 5 taken 8386 times.
 51859 if(keep_vertex_[v]) {
1042
1/2
✓ Branch 1 taken 43473 times.
✗ Branch 2 not taken.
 43473 ExactPoint P = I_.exact_vertex(v);
1043
5/10
✓ Branch 1 taken 43473 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 43473 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 43473 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 43473 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 43473 times.
✗ Branch 14 not taken.
 43473 v_idx_[v] = CDT.insert(exact::vec2h(P[u_], P[v_], P.w), v);
1044 43473 } else {
1045
1/2
✓ Branch 1 taken 8386 times.
✗ Branch 2 not taken.
 8386 v_idx_[v] = NO_INDEX;
1046 }
1047 }
1048
1049 // Insert constraints
1050
3/4
✓ Branch 2 taken 8067 times.
✗ Branch 3 not taken.
✓ Branch 7 taken 8170 times.
✓ Branch 8 taken 8067 times.
 16237 for(index_t P: polylines_) {
1051 8170 vector<index_t> Pvertices;
1052
1/2
✓ Branch 1 taken 8170 times.
✗ Branch 2 not taken.
 8170 index_t v = polylines_.first_vertex(P);
1053
3/4
✓ Branch 2 taken 8170 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 7203 times.
✓ Branch 5 taken 967 times.
 8170 if(keep_vertex_[v]) {
1054
1/2
✓ Branch 1 taken 7203 times.
✗ Branch 2 not taken.
 7203 Pvertices.push_back(v);
1055 }
1056
4/6
✓ Branch 1 taken 8170 times.
✗ Branch 2 not taken.
✓ Branch 6 taken 51861 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 51861 times.
✓ Branch 11 taken 8170 times.
 60031 for(index_t h: polylines_.halfedges(P)) {
1057
1/2
✓ Branch 1 taken 51861 times.
✗ Branch 2 not taken.
 51861 v = halfedges_.vertex(h,1);
1058
3/4
✓ Branch 2 taken 51861 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 43475 times.
✓ Branch 5 taken 8386 times.
 51861 if(keep_vertex_[v]) {
1059
1/2
✓ Branch 1 taken 43475 times.
✗ Branch 2 not taken.
 43475 Pvertices.push_back(v);
1060 }
1061 }
1062 8170 if(
1063
5/10
✓ Branch 1 taken 8170 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 8170 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 8170 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 8170 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 8170 times.
✗ Branch 11 not taken.
 16340 polylines_.first_vertex(P) == polylines_.last_vertex(P) &&
1064 8170 Pvertices.size() != 0
1065 ) {
1066
2/4
✓ Branch 1 taken 8170 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 8170 times.
✗ Branch 5 not taken.
 8170 Pvertices.push_back(Pvertices[0]);
1067 }
1068
1069
2/2
✓ Branch 1 taken 50678 times.
✓ Branch 2 taken 8170 times.
 58848 for(index_t i=0; i+1<Pvertices.size(); ++i) {
1070
2/4
✓ Branch 1 taken 50678 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 50678 times.
✗ Branch 5 not taken.
 50678 index_t v1 = v_idx_[Pvertices[i]];
1071
2/4
✓ Branch 1 taken 50678 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 50678 times.
✗ Branch 5 not taken.
 50678 index_t v2 = v_idx_[Pvertices[i+1]];
1072
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 50678 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 50678 geo_debug_assert(v1 != NO_INDEX);
1073
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 50678 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 50678 geo_debug_assert(v2 != NO_INDEX);
1074
1/2
✓ Branch 1 taken 50678 times.
✗ Branch 2 not taken.
 50678 CDT.insert_constraint(v1,v2,NO_INDEX);
1075 }
1076 8170 }
1077
1078
1/2
✓ Branch 1 taken 8067 times.
✗ Branch 2 not taken.
 8067 CDT.remove_external_triangles(true);
1079 8067 }
1080
1081 153794 bool CoplanarFacets::triangles_are_coplanar(
1082 const vec3& p1, const vec3& p2, const vec3& p3,
1083 const vec3& q1, const vec3& q2, const vec3& q3
1084 ) const {
1085 exact::vec3 N1 = cross(
1086
2/4
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 153794 times.
✗ Branch 5 not taken.
 307588 make_vec3<exact::vec3>(p1,p2), make_vec3<exact::vec3>(p1,p3)
1087
1/2
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
 153794 );
1088 exact::vec3 N2 = cross(
1089
2/4
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 153794 times.
✗ Branch 5 not taken.
 307588 make_vec3<exact::vec3>(q1,q2), make_vec3<exact::vec3>(q1,q3)
1090
1/2
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
 153794 );
1091
1092
9/14
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 35693 times.
✓ Branch 4 taken 118101 times.
✓ Branch 6 taken 35693 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 30251 times.
✓ Branch 9 taken 5442 times.
✓ Branch 11 taken 30251 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 30251 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 153794 times.
 153794 if(N1.x.sign() == ZERO && N1.y.sign() == ZERO && N1.z.sign() == ZERO) {
1093 ✗ std::cerr << std::endl;
1094 ✗ std::cerr << "degenerate triangle" << std::endl;
1095 ✗ std::cerr << "aligned: " << PCK::aligned_3d(p1,p2,p3) << std::endl;
1096 ✗ return false;
1097 }
1098
1099
9/14
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 36472 times.
✓ Branch 4 taken 117322 times.
✓ Branch 6 taken 36472 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 31046 times.
✓ Branch 9 taken 5426 times.
✓ Branch 11 taken 31046 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 31046 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 153794 times.
 153794 if(N2.x.sign() == ZERO && N2.y.sign() == ZERO && N2.z.sign() == ZERO) {
1100 ✗ std::cerr << std::endl;
1101 ✗ std::cerr << "degenerate triangle" << std::endl;
1102 ✗ std::cerr << "aligned: " << PCK::aligned_3d(q1,q2,q3) << std::endl;
1103 ✗ return false;
1104 }
1105
1106 // Tolerance for co-planarity test
1107
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 153794 times.
 153794 if(angle_tolerance_ != 0.0) {
1108 ✗ double threshold = cos(angle_tolerance_ * M_PI / 180.0);
1109 ✗ exact::scalar left = geo_sqr(dot(N1,N2));
1110 exact::scalar right =
1111 ✗ exact::scalar(threshold*threshold)*length2(N1)*length2(N2);
1112 ✗ return left > right;
1113 ✗ }
1114
1115 // Exact version
1116
1/2
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
 153794 exact::vec3 N12 = cross(N1,N2);
1117
11/14
✓ Branch 1 taken 153794 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 38745 times.
✓ Branch 4 taken 115049 times.
✓ Branch 6 taken 38745 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 32828 times.
✓ Branch 9 taken 5917 times.
✓ Branch 11 taken 32828 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 2490 times.
✓ Branch 14 taken 30338 times.
✓ Branch 15 taken 123456 times.
✓ Branch 16 taken 30338 times.
 153794 if((N12.x.sign()!=ZERO) || (N12.y.sign()!=ZERO) ||(N12.z.sign()!=ZERO)) {
1118 123456 return false;
1119 }
1120
1121 30338 return true;
1122 153794 }
1123
1124 /**************************************************************************/
1125
1126 246912 bool CoplanarFacets::edges_are_colinear(
1127 const ExactPoint& P1, const ExactPoint& P2, const ExactPoint& P3
1128 ) const {
1129
1130
1/2
✓ Branch 0 taken 246912 times.
✗ Branch 1 not taken.
 246912 if(angle_tolerance_ == 0.0) {
1131
1/2
✓ Branch 1 taken 246912 times.
✗ Branch 2 not taken.
 246912 return PCK::on_segment_3d(P2,P1,P3);
1132 }
1133
1134 ✗ ExactPoint UU = P1-P2;
1135 ✗ exact::vec3 U(UU.x, UU.y, UU.z);
1136 ✗ if(UU.w.sign() == NEGATIVE) {
1137 ✗ U.x.negate(); U.y.negate(); U.z.negate();
1138 }
1139 ✗ ExactPoint VV = P3-P2;
1140 ✗ exact::vec3 V(VV.x, VV.y, VV.z);
1141 ✗ if(VV.w.sign() == NEGATIVE) {
1142 ✗ V.x.negate(); V.y.negate(); V.z.negate();
1143 }
1144
1145 ✗ double threshold = cos(angle_tolerance_ * M_PI / 180.0);
1146
1147 ✗ exact::scalar left = dot(U,V);
1148
1149 ✗ if(left.sign() == POSITIVE) {
1150 ✗ return false;
1151 }
1152
1153 ✗ left = geo_sqr(left);
1154 exact::scalar right =
1155 ✗ exact::scalar(threshold*threshold)*length2(U)*length2(V);
1156
1157 ✗ return left > right;
1158 ✗ }
1159
1160 /**************************************************************************/
1161
1162 }
1163