GCC Code Coverage Report

Directory: ./
File: mesh/mesh_surface_intersection_internal.cpp
Date: 2026-09-27 03:22:43
Exec Total Coverage
Lines: 541 665 81.4%
Functions: 27 32 84.4%
Branches: 547 1354 40.4%
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 4334 times.
✓ Branch 4 taken 1350 times.
 5684 if(result.w.sign() == ZERO) {
83 4334 return false;
84 }
85
86
1/2
✓ Branch 1 taken 1350 times.
✗ Branch 2 not taken.
 2700 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 1350 times.
✗ Branch 2 not taken.
 1350 );
91
92
1/2
✓ Branch 1 taken 1350 times.
✗ Branch 2 not taken.
 2700 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 1350 times.
✗ Branch 2 not taken.
 1350 );
97
98
1/2
✓ Branch 1 taken 1350 times.
✗ Branch 2 not taken.
 2700 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 1350 times.
✗ Branch 2 not taken.
 1350 );
103
104 1350 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 67564 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 67564 times.
✗ Branch 2 not taken.
 67564 exact::vec3 D = make_vec3<exact::vec3>(q1,q2);
123
1/2
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
 67564 exact::vec3 E1 = make_vec3<exact::vec3>(p1,p2);
124
1/2
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
 67564 exact::vec3 E2 = make_vec3<exact::vec3>(p1,p3);
125
1/2
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
 67564 exact::vec3 AO = make_vec3<exact::vec3>(p1,q1);
126
1/2
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
 67564 exact::vec3 N = cross(E1,E2);
127
2/4
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 67564 times.
✗ Branch 5 not taken.
 67564 exact::scalar d = -dot(D,N);
128
2/8
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 67564 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 67564 geo_debug_assert(d.sign() != ZERO);
129
2/4
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 67564 times.
✗ Branch 5 not taken.
 67564 exact::rational t(dot(AO,N),d);
130
1/2
✓ Branch 1 taken 67564 times.
✗ Branch 2 not taken.
 135128 return mix(t,q1,q2);
131 67564 }
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 286659 MeshInTriangle::ExactPoint MeshInTriangle::Vertex::compute_geometry() {
153 // Case 1: f1 vertex
154
2/2
✓ Branch 1 taken 46518 times.
✓ Branch 2 taken 240141 times.
 286659 if(region_dim(sym.R1) == 0) {
155 46518 index_t lv = index_t(sym.R1) - index_t(T1_RGN_P0);
156
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 46518 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 46518 geo_debug_assert(lv < 3);
157
2/4
✓ Branch 1 taken 46518 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 46518 times.
✗ Branch 5 not taken.
 46518 mesh_vertex_index = mesh().facets.vertex(sym.f1,lv);
158
1/2
✓ Branch 1 taken 46518 times.
✗ Branch 2 not taken.
 46518 vec3 p = mit->mesh_vertex(mesh_vertex_index);
159
1/2
✓ Branch 1 taken 46518 times.
✗ Branch 2 not taken.
 46518 return ExactPoint(p);
160 }
161
162
2/8
✓ Branch 0 taken 240141 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 240141 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 240141 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 218486 times.
 240141 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 218486 if(
175
5/6
✓ Branch 1 taken 191064 times.
✓ Branch 2 taken 27422 times.
✓ Branch 4 taken 191064 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 191064 times.
✓ Branch 7 taken 27422 times.
 436972 (region_dim(sym.R1) == 2 || region_dim(sym.R1) == 1) &&
176
2/2
✓ Branch 1 taken 191064 times.
✓ Branch 2 taken 27422 times.
 218486 region_dim(sym.R2) == 1
177 ) {
178
1/2
✓ Branch 1 taken 191064 times.
✗ Branch 2 not taken.
 191064 vec3 p1 = mit->mesh_facet_vertex(sym.f1, 0);
179
1/2
✓ Branch 1 taken 191064 times.
✗ Branch 2 not taken.
 191064 vec3 p2 = mit->mesh_facet_vertex(sym.f1, 1);
180
1/2
✓ Branch 1 taken 191064 times.
✗ Branch 2 not taken.
 191064 vec3 p3 = mit->mesh_facet_vertex(sym.f1, 2);
181 191064 index_t e = index_t(sym.R2)-index_t(T2_RGN_E0);
182
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 191064 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 191064 geo_debug_assert(e<3);
183
1/2
✓ Branch 1 taken 191064 times.
✗ Branch 2 not taken.
 191064 vec3 q1 = mit->mesh_facet_vertex(sym.f2, (e+1)%3);
184
1/2
✓ Branch 1 taken 191064 times.
✗ Branch 2 not taken.
 191064 vec3 q2 = mit->mesh_facet_vertex(sym.f2, (e+2)%3);
185
186 bool seg_seg_two_D = (
187
1/2
✓ Branch 1 taken 191064 times.
✗ Branch 2 not taken.
 191064 region_dim(sym.R1) == 1 &&
188
5/6
✓ Branch 0 taken 163642 times.
✓ Branch 1 taken 27422 times.
✓ Branch 3 taken 163642 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 150982 times.
✓ Branch 6 taken 12660 times.
 342046 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 40082 times.
✓ Branch 1 taken 150982 times.
 191064 if(!seg_seg_two_D) {
192
1/2
✓ Branch 1 taken 40082 times.
✗ Branch 2 not taken.
 40082 return plane_line_intersection(p1,p2,p3,q1,q2);
193 }
194 }
195
196 // case 4: f1 edge /\ f2
197
5/6
✓ Branch 1 taken 178404 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 27422 times.
✓ Branch 5 taken 150982 times.
✓ Branch 6 taken 27422 times.
✓ Branch 7 taken 150982 times.
 178404 if(region_dim(sym.R1) == 1 && region_dim(sym.R2) == 2) {
198 27422 index_t e = index_t(sym.R1)-index_t(T1_RGN_E0);
199
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 27422 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 27422 geo_debug_assert(e<3);
200
1/2
✓ Branch 1 taken 27422 times.
✗ Branch 2 not taken.
 27422 vec3 p1 = mit->mesh_facet_vertex(sym.f2,0);
201
1/2
✓ Branch 1 taken 27422 times.
✗ Branch 2 not taken.
 27422 vec3 p2 = mit->mesh_facet_vertex(sym.f2,1);
202
1/2
✓ Branch 1 taken 27422 times.
✗ Branch 2 not taken.
 27422 vec3 p3 = mit->mesh_facet_vertex(sym.f2,2);
203
1/2
✓ Branch 1 taken 27422 times.
✗ Branch 2 not taken.
 27422 vec3 q1 = mit->mesh_facet_vertex(sym.f1, (e+1)%3);
204
1/2
✓ Branch 1 taken 27422 times.
✗ Branch 2 not taken.
 27422 vec3 q2 = mit->mesh_facet_vertex(sym.f1, (e+2)%3);
205
1/2
✓ Branch 1 taken 27422 times.
✗ Branch 2 not taken.
 27422 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 292343 void MeshInTriangle::Vertex::init_geometry(const ExactPoint& P) {
235 292343 point_exact = P;
236 292343 Numeric::optimize_number_representation(point_exact);
237 #ifndef GEOGRAM_USE_EXACT_NT
238 292343 l = ExactCDT2d::squared_length(P[mit->u_], P[mit->v_], P.w);
239 #endif
240 292343 }
241
242 228 MeshInTriangle::MeshInTriangle(MeshSurfaceIntersection& EM) :
243 228 exact_mesh_(EM),
244 456 mesh_(EM.readonly_mesh()),
245 228 f1_(NO_INDEX),
246 228 dry_run_(false),
247 228 use_pred_cache_insert_buffer_(false)
248 {
249 #ifdef GEOGRAM_USE_EXACT_NT
250 CDTBase2d::exact_incircle_ = true;
251 #else
252 // Since incircle() with expansions computes approximated
253 // lifted coordinate, we need to activate additional
254 // checks for Delaunayization.
255 228 CDTBase2d::exact_incircle_ = false;
256 #endif
257 228 }
258
259 15506 void MeshInTriangle::clear() {
260 15506 vertex_.resize(0);
261 15506 edges_.resize(0);
262 15506 f1_ = NO_INDEX;
263 15506 pred_cache_.clear();
264 15506 pred_cache_insert_buffer_.resize(0);
265 15506 use_pred_cache_insert_buffer_ = false;
266 15506 CDTBase2d::clear();
267 15506 }
268
269 15506 void MeshInTriangle::begin_facet(index_t f) {
270 15506 f1_ = f;
271
272
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 vec3 p1 = mesh_facet_vertex(f,0);
273
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 vec3 p2 = mesh_facet_vertex(f,1);
274
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 vec3 p3 = mesh_facet_vertex(f,2);
275
276
2/8
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 15506 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 15506 geo_debug_assert(!PCK::aligned_3d(p1,p2,p3));
277
278
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 f1_normal_axis_ = PCK::triangle_normal_axis(
279 p1,p2,p3
280 );
281
282 15506 u_ = coord_index_t((f1_normal_axis_ + 1) % 3);
283 15506 v_ = coord_index_t((f1_normal_axis_ + 2) % 3);
284
2/2
✓ Branch 0 taken 46518 times.
✓ Branch 1 taken 15506 times.
 62024 for(index_t lv=0; lv<3; ++lv) {
285
1/2
✓ Branch 1 taken 46518 times.
✗ Branch 2 not taken.
 46518 vertex_.emplace_back(this, f, lv);
286 }
287
288
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 CDTBase2d::create_enclosing_triangle(0,1,2);
289
290
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 edges_.emplace_back(1,2);
291
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 edges_.emplace_back(2,0);
292
1/2
✓ Branch 1 taken 15506 times.
✗ Branch 2 not taken.
 15506 edges_.emplace_back(0,1);
293 15506 }
294
295 262926 index_t MeshInTriangle::add_vertex(
296 index_t f2, TriangleRegion R1, TriangleRegion R2
297 ) {
298
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 262926 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 262926 geo_debug_assert(f1_ != NO_INDEX);
299
300 // If vertex is a macro-vertex, return it directly.
301
2/2
✓ Branch 1 taken 22785 times.
✓ Branch 2 taken 240141 times.
 262926 if(region_dim(R1) == 0) {
302 22785 return index_t(R1);
303 }
304
305 // Create the vertex
306
1/2
✓ Branch 1 taken 240141 times.
✗ Branch 2 not taken.
 240141 vertex_.emplace_back(this, f1_, f2, R1, R2);
307
308 // Insert it into the triangulation
309 240141 index_t v = CDTBase2d::insert(vertex_.size()-1);
310
311 // If it was an existing vertex, return the existing vertex
312
2/2
✓ Branch 2 taken 151848 times.
✓ Branch 3 taken 88293 times.
 240141 if(vertex_.size() > CDTBase2d::nv()) {
313 151848 vertex_.pop_back();
314 }
315 240141 return v;
316 }
317
318 127160 void MeshInTriangle::add_edge(
319 index_t f2,
320 TriangleRegion AR1, TriangleRegion AR2,
321 TriangleRegion BR1, TriangleRegion BR2
322 ) {
323
1/2
✓ Branch 1 taken 127160 times.
✗ Branch 2 not taken.
 127160 index_t v1 = add_vertex(f2, AR1, AR2);
324
1/2
✓ Branch 1 taken 127160 times.
✗ Branch 2 not taken.
 127160 index_t v2 = add_vertex(f2, BR1, BR2);
325
326 // If both extremities are on the same edge of f1,
327 // we do not add the edge, because it will be generated
328 // when remeshing the edge of f1
329
4/6
✓ Branch 1 taken 127160 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 127160 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 49028 times.
✓ Branch 7 taken 78132 times.
 127160 if(region_dim(regions_convex_hull(AR1,BR1)) == 1) {
330 49028 return;
331 }
332
333 // Generate also the combinatorial information of the edge,
334 // that indicates whether both extremities are on the same
335 // edge of f2 (useful later to compute the intersections)
336
2/4
✓ Branch 1 taken 78132 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 78132 times.
✗ Branch 5 not taken.
 78132 edges_.emplace_back(v1,v2,f2,regions_convex_hull(AR2,BR2));
337
338 // Constraints will be added to the triangulation during commit()
339 }
340
341 15506 void MeshInTriangle::commit() {
342
343
2/2
✓ Branch 2 taken 124650 times.
✓ Branch 3 taken 15506 times.
 155662 for(const Edge& E: edges_) {
344
1/2
✓ Branch 1 taken 124650 times.
✗ Branch 2 not taken.
 124650 CDTBase2d::insert_constraint(E.v1, E.v2);
345 }
346
347
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 15506 times.
 15506 if(dry_run_) {
348 ✗ return;
349 }
350
351 // Protect global mesh from concurrent accesses
352 15506 exact_mesh_.lock();
353
354 // Create vertices and facets in target mesh
355
2/2
✓ Branch 1 taken 140495 times.
✓ Branch 2 taken 15506 times.
 156001 for(index_t i=0; i<vertex_.size(); ++i) {
356 // Vertex already exists in this MeshInTriangle
357
2/2
✓ Branch 1 taken 49254 times.
✓ Branch 2 taken 91241 times.
 140495 if(vertex_[i].mesh_vertex_index != NO_INDEX) {
358 49254 continue;
359 }
360 91241 vertex_[i].mesh_vertex_index =
361 91241 exact_mesh_.find_or_create_exact_vertex(
362 91241 vertex_[i].point_exact
363 );
364 }
365
366 // Create facets in target mesh
367
2/2
✓ Branch 1 taken 131392 times.
✓ Branch 2 taken 15506 times.
 146898 for(index_t t=0; t<CDTBase2d::nT(); ++t) {
368 131392 index_t i = CDTBase2d::Tv(t,0);
369 131392 index_t j = CDTBase2d::Tv(t,1);
370 131392 index_t k = CDTBase2d::Tv(t,2);
371 131392 i = vertex_[i].mesh_vertex_index;
372 131392 j = vertex_[j].mesh_vertex_index;
373 131392 k = vertex_[k].mesh_vertex_index;
374 131392 index_t new_t = target_mesh().facets.create_triangle(i,j,k);
375 // Copy all attributes from initial facet
376 131392 target_mesh().facets.attributes().copy_item(new_t, f1_);
377 }
378
379 // We are done with modification in the mesh
380 15506 exact_mesh_.unlock();
381 }
382
383 ✗ void MeshInTriangle::get_constraints(Mesh& M, bool with_edges) const {
384 ✗ if(M.vertices.nb() == 0) {
385 ✗ M.vertices.set_dimension(2);
386 ✗ for(index_t v=0; v<vertex_.size(); ++v) {
387 ✗ vec2 p = vertex_[v].get_UV_approx();
388 ✗ M.vertices.create_vertex(p.data());
389 }
390 }
391 ✗ if(with_edges && M.edges.nb() == 0) {
392 ✗ for(const Edge& E: edges_) {
393 ✗ M.edges.create_edge(E.v1, E.v2);
394 }
395 }
396 ✗ }
397
398 /**
399 * \brief Tests the parity of the permutation of a list of
400 * three distinct indices with respect to the canonical order.
401 */
402 5734310 static bool odd_order(index_t i, index_t j, index_t k) {
403 // Implementation: sort the elements (bubble sort is OK for
404 // such a small number), and invert parity each time
405 // two elements are swapped.
406 5734310 index_t tab[3] = { i, j, k};
407 5734310 const int N = 3;
408 5734310 bool result = false;
409
2/2
✓ Branch 0 taken 11468620 times.
✓ Branch 1 taken 5734310 times.
 17202930 for (int I = 0; I < N - 1; ++I) {
410
2/2
✓ Branch 0 taken 17202930 times.
✓ Branch 1 taken 11468620 times.
 28671550 for (int J = 0; J < N - I - 1; ++J) {
411
2/2
✓ Branch 0 taken 10353758 times.
✓ Branch 1 taken 6849172 times.
 17202930 if (tab[J] > tab[J + 1]) {
412 10353758 std::swap(tab[J], tab[J + 1]);
413 10353758 result = !result;
414 }
415 }
416 }
417 5734310 return result;
418 }
419
420 245825 void MeshInTriangle::begin_insert_transaction() {
421 245825 use_pred_cache_insert_buffer_ = true;
422 245825 }
423
424 88293 void MeshInTriangle::commit_insert_transaction() {
425
2/2
✓ Branch 2 taken 891614 times.
✓ Branch 3 taken 88293 times.
 1068200 for(const auto& it: pred_cache_insert_buffer_) {
426
1/2
✓ Branch 1 taken 891614 times.
✗ Branch 2 not taken.
 891614 pred_cache_[it.first] = it.second;
427 }
428 88293 pred_cache_insert_buffer_.resize(0);
429 88293 use_pred_cache_insert_buffer_ = false;
430 88293 }
431
432 157532 void MeshInTriangle::rollback_insert_transaction() {
433 157532 pred_cache_insert_buffer_.resize(0);
434 157532 use_pred_cache_insert_buffer_ = false;
435 157532 }
436
437 5734310 Sign MeshInTriangle::orient2d(index_t vx1,index_t vx2,index_t vx3) const {
438
439
1/2
✓ Branch 1 taken 5734310 times.
✗ Branch 2 not taken.
 5734310 trindex K(vx1, vx2, vx3);
440
441
2/2
✓ Branch 0 taken 3632851 times.
✓ Branch 1 taken 2101459 times.
 5734310 if(use_pred_cache_insert_buffer_) {
442 10898553 Sign result = PCK::orient_2d_projected(
443
1/2
✓ Branch 1 taken 3632851 times.
✗ Branch 2 not taken.
 3632851 vertex_[K.indices[0]].point_exact,
444
1/2
✓ Branch 1 taken 3632851 times.
✗ Branch 2 not taken.
 3632851 vertex_[K.indices[1]].point_exact,
445
1/2
✓ Branch 1 taken 3632851 times.
✗ Branch 2 not taken.
 3632851 vertex_[K.indices[2]].point_exact,
446
1/2
✓ Branch 1 taken 3632851 times.
✗ Branch 2 not taken.
 3632851 f1_normal_axis_
447 3632851 );
448
2/4
✓ Branch 1 taken 3632851 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3632851 times.
✗ Branch 5 not taken.
 3632851 pred_cache_insert_buffer_.push_back(std::make_pair(K, result));
449
2/2
✓ Branch 1 taken 1783521 times.
✓ Branch 2 taken 1849330 times.
 3632851 if(odd_order(vx1,vx2,vx3)) {
450 1783521 result = Sign(-result);
451 }
452 3632851 return result;
453 }
454
455 bool inserted;
456 2101459 std::map<trindex, Sign>::iterator it;
457
2/4
✓ Branch 1 taken 2101459 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2101459 times.
✗ Branch 5 not taken.
 2101459 std::tie(it,inserted) = pred_cache_.insert(std::make_pair(K,ZERO));
458 Sign result;
459
460
2/2
✓ Branch 0 taken 691388 times.
✓ Branch 1 taken 1410071 times.
 2101459 if(inserted) {
461 2074164 result = PCK::orient_2d_projected(
462
1/2
✓ Branch 1 taken 691388 times.
✗ Branch 2 not taken.
 691388 vertex_[K.indices[0]].point_exact,
463
1/2
✓ Branch 1 taken 691388 times.
✗ Branch 2 not taken.
 691388 vertex_[K.indices[1]].point_exact,
464
1/2
✓ Branch 1 taken 691388 times.
✗ Branch 2 not taken.
 691388 vertex_[K.indices[2]].point_exact,
465
1/2
✓ Branch 1 taken 691388 times.
✗ Branch 2 not taken.
 691388 f1_normal_axis_
466 );
467 691388 it->second = result;
468 } else {
469 1410071 result = it->second;
470 }
471
472
2/2
✓ Branch 1 taken 963237 times.
✓ Branch 2 taken 1138222 times.
 2101459 if(odd_order(vx1,vx2,vx3)) {
473 963237 result = Sign(-result);
474 }
475
476 2101459 return result;
477 }
478
479 187312 Sign MeshInTriangle::incircle(
480 index_t v1,index_t v2,index_t v3,index_t v4
481 ) const {
482 exact::vec2h p1(
483 374624 vertex_[v1].point_exact[u_],
484 374624 vertex_[v1].point_exact[v_],
485
2/4
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
 187312 vertex_[v1].point_exact.w
486
4/8
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 187312 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 187312 times.
✗ Branch 11 not taken.
 187312 );
487 exact::vec2h p2(
488 374624 vertex_[v2].point_exact[u_],
489 374624 vertex_[v2].point_exact[v_],
490
2/4
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
 187312 vertex_[v2].point_exact.w
491
4/8
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 187312 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 187312 times.
✗ Branch 11 not taken.
 187312 );
492 exact::vec2h p3(
493 374624 vertex_[v3].point_exact[u_],
494 374624 vertex_[v3].point_exact[v_],
495
2/4
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
 187312 vertex_[v3].point_exact.w
496
4/8
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 187312 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 187312 times.
✗ Branch 11 not taken.
 187312 );
497 exact::vec2h p4(
498 374624 vertex_[v4].point_exact[u_],
499 374624 vertex_[v4].point_exact[v_],
500
2/4
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
 187312 vertex_[v4].point_exact.w
501
4/8
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 187312 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 187312 times.
✗ Branch 11 not taken.
 187312 );
502 #ifdef GEOGRAM_USE_EXACT_NT
503 return PCK::incircle_2d_SOS(p1,p2,p3,p4);
504 #else
505 374624 return PCK::incircle_2d_SOS_with_lengths(
506 p1,p2,p3,p4,
507
1/2
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
 187312 vertex_[v1].l,
508
1/2
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
 187312 vertex_[v2].l,
509
1/2
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
 187312 vertex_[v3].l,
510
2/4
✓ Branch 1 taken 187312 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 187312 times.
✗ Branch 5 not taken.
 187312 vertex_[v4].l
511 374624 );
512 #endif
513 187312 }
514
515 5684 index_t MeshInTriangle::create_intersection(
516 index_t e1, index_t i, index_t j,
517 index_t e2, index_t k, index_t l
518 ) {
519 5684 geo_argused(i);
520 5684 geo_argused(j);
521 5684 geo_argused(k);
522 5684 geo_argused(l);
523
524 5684 ExactPoint I;
525
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 get_edge_edge_intersection(e1,e2,I);
526
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 vertex_.emplace_back(this,I);
527 5684 index_t x = vertex_.size()-1;
528
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 CDTBase2d::v2T_.push_back(NO_INDEX);
529
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_);
530 5684 ++CDTBase2d::nv_;
531 5684 return x;
532 5684 }
533
534 5684 void MeshInTriangle::get_edge_edge_intersection(
535 index_t e1, index_t e2, ExactPoint& I
536 ) const {
537 5684 index_t f1 = f1_;
538
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 index_t f2 = edges_[e1].sym.f2;
539
1/2
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
 5684 index_t f3 = edges_[e2].sym.f2;
540
541
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);
542
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);
543
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);
544
545 vec3 P[9] = {
546 11368 mesh_facet_vertex(f1,0), mesh_facet_vertex(f1,1),
547 5684 mesh_facet_vertex(f1,2),
548 11368 mesh_facet_vertex(f2,0), mesh_facet_vertex(f2,1),
549 5684 mesh_facet_vertex(f2,2),
550 11368 mesh_facet_vertex(f3,0), mesh_facet_vertex(f3,1),
551 5684 mesh_facet_vertex(f3,2)
552
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 };
553
554
3/4
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4334 times.
✓ Branch 4 taken 1350 times.
 5684 if(!get_three_planes_intersection(
555 I,
556 P[0], P[1], P[2],
557 P[3], P[4], P[5],
558 P[6], P[7], P[8]
559 )) {
560
1/2
✓ Branch 1 taken 4334 times.
✗ Branch 2 not taken.
 4334 get_edge_edge_intersection_2D(e1,e2,I);
561 4334 return;
562 }
563 }
564
565 4334 void MeshInTriangle::get_edge_edge_intersection_2D(
566 index_t e1, index_t e2, ExactPoint& I
567 ) const {
568 4334 const Edge& E1 = edges_[e1];
569 4334 const Edge& E2 = edges_[e2];
570
571 4334 if(
572
4/4
✓ Branch 1 taken 4288 times.
✓ Branch 2 taken 46 times.
✓ Branch 3 taken 4274 times.
✓ Branch 4 taken 60 times.
 8622 region_dim(E1.sym.R2) == 1 &&
573
2/2
✓ Branch 1 taken 4274 times.
✓ Branch 2 taken 14 times.
 4288 region_dim(E2.sym.R2) == 1
574 ) {
575 4274 index_t le1 = index_t(E1.sym.R2)-index_t(T2_RGN_E0);
576 4274 index_t le2 = index_t(E2.sym.R2)-index_t(T2_RGN_E0);
577
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 4274 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 4274 geo_debug_assert(le1 < 3);
578
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 4274 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 4274 geo_debug_assert(le2 < 3);
579
580
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 vec2 p1_uv = mesh_facet_vertex_UV(E1.sym.f2, (le1+1)%3);
581
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 vec2 p2_uv = mesh_facet_vertex_UV(E1.sym.f2, (le1+2)%3);
582
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 vec2 q1_uv = mesh_facet_vertex_UV(E2.sym.f2, (le2+1)%3);
583
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 vec2 q2_uv = mesh_facet_vertex_UV(E2.sym.f2, (le2+2)%3);
584
585
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 exact::vec2 C1 = make_vec2<exact::vec2>(p1_uv, p2_uv);
586
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 exact::vec2 C2 = make_vec2<exact::vec2>(q2_uv, q1_uv);
587
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 exact::vec2 B = make_vec2<exact::vec2>(p1_uv, q1_uv);
588
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 exact::scalar d = det(C1,C2);
589
2/8
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 4274 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 4274 geo_debug_assert(d.sign() != ZERO);
590
2/4
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4274 times.
✗ Branch 5 not taken.
 4274 exact::rational t(det(B,C2),d);
591
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 8548 I = mix(
592 t,
593
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 mesh_facet_vertex(E1.sym.f2,(le1+1)%3),
594
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 mesh_facet_vertex(E1.sym.f2,(le1+2)%3)
595
1/2
✓ Branch 1 taken 4274 times.
✗ Branch 2 not taken.
 4274 );
596
597 4274 } else {
598
6/14
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 46 times.
✓ Branch 4 taken 14 times.
✓ Branch 6 taken 46 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 46 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 60 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
 60 geo_debug_assert(
599 region_dim(E1.sym.R2) == 1 || region_dim(E2.sym.R2) == 1
600 );
601 60 index_t f1 = E1.sym.f2;
602 60 TriangleRegion R1 = E1.sym.R2;
603 60 index_t f2 = E2.sym.f2;
604 60 TriangleRegion R2 = E2.sym.R2;
605
3/4
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 14 times.
✓ Branch 4 taken 46 times.
 60 if(region_dim(R1) == 1) {
606 14 std::swap(f1,f2);
607 14 std::swap(R1,R2);
608 }
609
610 60 index_t e = index_t(R2) - index_t(T2_RGN_E0);
611
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 60 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 60 geo_debug_assert(e < 3);
612
613
1/2
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
 120 I = plane_line_intersection(
614
1/2
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
 60 mesh_facet_vertex(f1,0),
615
1/2
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
 60 mesh_facet_vertex(f1,1),
616
1/2
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
 60 mesh_facet_vertex(f1,2),
617
1/2
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
 60 mesh_facet_vertex(f2,(e+1)%3),
618
1/2
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
 120 mesh_facet_vertex(f2,(e+2)%3)
619
1/2
✓ Branch 1 taken 60 times.
✗ Branch 2 not taken.
 60 );
620 }
621 4334 }
622
623 ✗ void MeshInTriangle::save(const std::string& filename) const {
624 ✗ Mesh M;
625 ✗ M.vertices.set_dimension(2);
626 ✗ for(index_t v=0; v<CDTBase2d::nv(); ++v) {
627 ✗ vec2 p = vertex_[v].get_UV_approx();
628 ✗ M.vertices.create_vertex(p.data());
629 }
630 ✗ for(index_t t=0; t<CDTBase2d::nT(); ++t) {
631 ✗ M.facets.create_triangle(
632 CDTBase2d::Tv(t,0),
633 CDTBase2d::Tv(t,1),
634 CDTBase2d::Tv(t,2)
635 );
636 }
637
638 ✗ Attribute<double> tex_coord;
639 ✗ tex_coord.create_vector_attribute(
640 ✗ M.facet_corners.attributes(), "tex_coord", 2
641 );
642 static double triangle_tex[3][2] = {
643 {0.0, 0.0},
644 {1.0, 0.0},
645 {0.0, 1.0}
646 };
647 ✗ for(index_t c: M.facet_corners) {
648 ✗ tex_coord[2*c] = triangle_tex[c%3][0];
649 ✗ tex_coord[2*c+1] = triangle_tex[c%3][1];
650 }
651 ✗ mesh_save(M, filename);
652 ✗ }
653
654 /**************************************************************************/
655
656 285 CoplanarFacets::CoplanarFacets(
657 MeshSurfaceIntersection& I, bool clear_attributes,
658 double angle_tolerance
659 285 ) :
660 285 I_(I),
661 285 mesh_(I.target_mesh()),
662 285 mesh_copy_(I.readonly_mesh()),
663 285 angle_tolerance_(angle_tolerance),
664
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"),
665
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"),
666
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 570 c_is_coplanar_(
667
1/2
✓ Branch 1 taken 285 times.
✗ Branch 2 not taken.
 570 I.target_mesh().facet_corners.attributes(),"is_coplanar"
668 ),
669
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"),
670
1/2
✓ Branch 1 taken 285 times.
✗ Branch 2 not taken.
 285 halfedges_(*this),
671
1/2
✓ Branch 8 taken 285 times.
✗ Branch 9 not taken.
 855 polylines_(*this)
672 {
673
2/2
✓ Branch 0 taken 57 times.
✓ Branch 1 taken 228 times.
 285 if(clear_attributes) {
674
2/2
✓ Branch 4 taken 123494 times.
✓ Branch 5 taken 57 times.
 123551 for(index_t f: mesh_.facets) {
675
1/2
✓ Branch 1 taken 123494 times.
✗ Branch 2 not taken.
 123494 facet_group_[f] = NO_INDEX;
676 }
677
2/2
✓ Branch 4 taken 61753 times.
✓ Branch 5 taken 57 times.
 61810 for(index_t v: mesh_.vertices) {
678
1/2
✓ Branch 2 taken 61753 times.
✗ Branch 3 not taken.
 61753 keep_vertex_[v] = false;
679 }
680
1/2
✓ Branch 1 taken 57 times.
✗ Branch 2 not taken.
 57 find_coplanar_facets();
681 }
682
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 f_visited_.assign(mesh_.facets.nb(),false);
683
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 h_visited_.assign(mesh_.facet_corners.nb(),false);
684
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 v_visited_.assign(mesh_.vertices.nb(),false);
685
1/2
✓ Branch 2 taken 285 times.
✗ Branch 3 not taken.
 285 v_idx_.assign(mesh_.vertices.nb(),NO_INDEX);
686 285 }
687
688 57 void CoplanarFacets::find_coplanar_facets() {
689
690 // Positioned by MeshSurfaceIntersection::build_Weiler_model()
691 Attribute<bool> corner_is_on_border(
692
1/2
✓ Branch 1 taken 57 times.
✗ Branch 2 not taken.
 114 mesh_.facet_corners.attributes(), "is_on_border"
693
1/2
✓ Branch 2 taken 57 times.
✗ Branch 3 not taken.
 57 );
694
695
2/2
✓ Branch 4 taken 370482 times.
✓ Branch 5 taken 57 times.
 370539 for(index_t c: mesh_.facet_corners) {
696
1/2
✓ Branch 2 taken 370482 times.
✗ Branch 3 not taken.
 370482 c_is_coplanar_[c] = false;
697 }
698
699 // TODO: when there is an angle tolerance, one should check instead
700 // angle deviation w.r.t. a single seed facet per facet group, because
701 // with the present algorithm, if a large number of tiny facets are
702 // connected (e.g. highly tessellated cylinder), one may group facets
703 // with large angle deviation (without seeing it because each facet has
704 // small angle deviation w.r.t. its neighbors).
705
706
1/2
✓ Branch 1 taken 57 times.
✗ Branch 2 not taken.
 57 parallel_for(
707 57 0, mesh_.facet_corners.nb(),
708 114 [&](index_t c1) {
709 370482 index_t f1 = (c1 / 3);
710 370482 index_t le1 = (c1 % 3);
711
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t f2 = mesh_.facet_corners.adjacent_facet(c1);
712
713
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 370482 times.
 370482 if(f2 == NO_INDEX) {
714 216697 return;
715 }
716
717 // do not traverse true borders
718
2/4
✓ Branch 2 taken 370482 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 370482 times.
 370482 if(corner_is_on_border[c1]) {
719 ✗ return;
720 }
721
722
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t v11 = mesh_.facets.vertex(f1,le1);
723
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t v12 = mesh_.facets.vertex(f1,(le1+1)%3);
724
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t le2 = mesh_.facets.find_edge(f2,v12,v11);
725
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t c2 = mesh_.facets.corner(f2,le2);
726
727 #ifdef GEO_DEBUG
728
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t v21 = mesh_.facets.vertex(f2,le2);
729
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t v22 = mesh_.facets.vertex(f2,(le2+1)%3);
730
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t v13 = mesh_.facets.vertex(f1,(le1+2)%3);
731
1/2
✓ Branch 1 taken 370482 times.
✗ Branch 2 not taken.
 370482 index_t v23 = mesh_.facets.vertex(f2,(le2+2)%3);
732
733
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 370482 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 370482 geo_debug_assert(v11 == v22);
734
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 370482 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 370482 geo_debug_assert(v12 == v21);
735
3/10
✓ Branch 0 taken 370482 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 370482 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 370482 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 370482 geo_debug_assert(v11!=v12 && v12!=v13 && v13!=v11);
736
3/10
✓ Branch 0 taken 370482 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 370482 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 370482 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 370482 geo_debug_assert(v21!=v22 && v22!=v23 && v23!=v21);
737 #endif
738
739
2/2
✓ Branch 0 taken 185241 times.
✓ Branch 1 taken 185241 times.
 370482 if(c1 > c2) {
740 185241 return;
741 }
742
743 // Small optimization: if both triangles come from same
744 // original facet then they are coplanar
745
4/6
✓ Branch 1 taken 185241 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 185241 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 31456 times.
✓ Branch 7 taken 153785 times.
 185241 if(I_.get_initial_facet(f1) == I_.get_initial_facet(f2)) {
746
1/2
✓ Branch 2 taken 31456 times.
✗ Branch 3 not taken.
 31456 c_is_coplanar_[c1] = true;
747
1/2
✓ Branch 2 taken 31456 times.
✗ Branch 3 not taken.
 31456 c_is_coplanar_[c2] = true;
748 31456 return;
749 }
750
751 // Use original triangles for co-planarity test
752
1/2
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
 153785 auto [p1, p2, p3] = I_.get_initial_facet_vertices(f1);
753
1/2
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
 153785 auto [q1, q2, q3] = I_.get_initial_facet_vertices(f2);
754
3/4
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 30326 times.
✓ Branch 4 taken 123459 times.
 153785 if(triangles_are_coplanar(p1,p2,p3,q1,q2,q3)) {
755
1/2
✓ Branch 2 taken 30326 times.
✗ Branch 3 not taken.
 30326 c_is_coplanar_[c1] = true;
756
1/2
✓ Branch 2 taken 30326 times.
✗ Branch 3 not taken.
 30326 c_is_coplanar_[c2] = true;
757 }
758 }
759 );
760 57 }
761
762 146172 void CoplanarFacets::get(index_t f, index_t group_id) {
763
764 146172 facets_.resize(0);
765 146172 vertices_.resize(0);
766 146172 halfedges_.initialize();
767 146172 polylines_.initialize();
768
769 // Get facets
770 {
771
1/2
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
 146172 std::stack<index_t> S;
772
1/2
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
 146172 facet_group_[f] = group_id;
773 146172 f_visited_[f] = true;
774
1/2
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
 146172 S.push(f);
775
1/2
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
 146172 facets_.push_back(f);
776
2/2
✓ Branch 1 taken 246988 times.
✓ Branch 2 taken 146172 times.
 393160 while(!S.empty()) {
777 246988 index_t f1 = S.top();
778 246988 S.pop();
779
2/2
✓ Branch 0 taken 740964 times.
✓ Branch 1 taken 246988 times.
 987952 for(index_t le1=0; le1<3; ++le1) {
780
1/2
✓ Branch 1 taken 740964 times.
✗ Branch 2 not taken.
 740964 index_t f2 = mesh_.facets.adjacent(f1,le1);
781 740964 if(
782
3/4
✓ Branch 0 taken 740964 times.
✗ Branch 1 not taken.
✓ Branch 4 taken 594652 times.
✓ Branch 5 taken 146312 times.
 1335616 f2 != NO_INDEX && !f_visited_[f2] &&
783
6/8
✓ Branch 1 taken 594652 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 594652 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 100816 times.
✓ Branch 8 taken 493836 times.
✓ Branch 9 taken 100816 times.
✓ Branch 10 taken 640148 times.
 1335616 c_is_coplanar_[mesh_.facets.corner(f1,le1)]
784 ) {
785
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];
786 100816 f_visited_[f2] = true;
787
1/2
✓ Branch 1 taken 100816 times.
✗ Branch 2 not taken.
 100816 S.push(f2);
788
1/2
✓ Branch 1 taken 100816 times.
✗ Branch 2 not taken.
 100816 facets_.push_back(f2);
789 }
790 }
791 }
792
2/2
✓ Branch 2 taken 246988 times.
✓ Branch 3 taken 146172 times.
 539332 for(index_t cur_f: facets_) {
793 246988 f_visited_[cur_f] = false;
794 }
795 146172 }
796
797 146172 group_id_ = group_id;
798
799 // Initialize projection coordinates, using original facet
800 {
801
2/4
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 146172 times.
✗ Branch 5 not taken.
 146172 auto [p1, p2, p3] = I_.get_initial_facet_vertices(facets_[0]);
802
1/2
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
 146172 coord_index_t projection_axis = PCK::triangle_normal_axis(p1,p2,p3);
803 146172 u_ = coord_index_t((projection_axis+1)%3);
804 146172 v_ = coord_index_t((projection_axis+2)%3);
805
1/2
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
 146172 Sign o = PCK::orient_2d(
806
6/12
✓ Branch 1 taken 146172 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 146172 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 146172 times.
✗ Branch 9 not taken.
✓ Branch 11 taken 146172 times.
✗ Branch 12 not taken.
✓ Branch 15 taken 146172 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 146172 times.
✗ Branch 19 not taken.
 146172 vec2(p1[u_],p1[v_]), vec2(p2[u_],p2[v_]), vec2(p3[u_],p3[v_])
807 );
808
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 146172 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 146172 geo_debug_assert(o != ZERO);
809
2/2
✓ Branch 0 taken 71090 times.
✓ Branch 1 taken 75082 times.
 146172 if(o < 0) {
810 71090 std::swap(u_,v_);
811 }
812 }
813
814 // Get vertices and halfedges
815 {
816
2/2
✓ Branch 2 taken 246988 times.
✓ Branch 3 taken 146172 times.
 539332 for(index_t f1: facets_) {
817
2/2
✓ Branch 0 taken 740964 times.
✓ Branch 1 taken 246988 times.
 987952 for(index_t le=0; le<3; ++le) {
818
1/2
✓ Branch 1 taken 740964 times.
✗ Branch 2 not taken.
 740964 index_t f2 = mesh_.facets.adjacent(f1,le);
819 740964 if(
820
1/2
✓ Branch 0 taken 740964 times.
✗ Branch 1 not taken.
 1481928 f2 == NO_INDEX ||
821
6/8
✓ Branch 1 taken 740964 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 740964 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 493836 times.
✓ Branch 8 taken 247128 times.
✓ Branch 9 taken 493836 times.
✓ Branch 10 taken 247128 times.
 1481928 !c_is_coplanar_[mesh_.facets.corner(f1,le)]
822 ) {
823
2/4
✓ Branch 1 taken 493836 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 493836 times.
✗ Branch 5 not taken.
 493836 halfedges_.add(mesh_.facets.corners_begin(f1)+le);
824
1/2
✓ Branch 1 taken 493836 times.
✗ Branch 2 not taken.
 493836 index_t v1 = mesh_.facets.vertex(f1,le);
825
1/2
✓ Branch 1 taken 493836 times.
✗ Branch 2 not taken.
 493836 index_t v2 = mesh_.facets.vertex(f1,(le+1)%3);
826
2/2
✓ Branch 2 taken 174788 times.
✓ Branch 3 taken 319048 times.
 493836 if(!v_visited_[v1]) {
827
1/2
✓ Branch 2 taken 174788 times.
✗ Branch 3 not taken.
 174788 v_idx_[v1] = vertices_.size();
828
1/2
✓ Branch 1 taken 174788 times.
✗ Branch 2 not taken.
 174788 vertices_.push_back(v1);
829 174788 v_visited_[v1] = true;
830 }
831
2/2
✓ Branch 2 taken 307370 times.
✓ Branch 3 taken 186466 times.
 493836 if(!v_visited_[v2]) {
832
1/2
✓ Branch 2 taken 307370 times.
✗ Branch 3 not taken.
 307370 v_idx_[v2] = vertices_.size();
833
1/2
✓ Branch 1 taken 307370 times.
✗ Branch 2 not taken.
 307370 vertices_.push_back(v2);
834 307370 v_visited_[v2] = true;
835 }
836 } else {
837 // This one for the particular case of a non-manifold
838 // vertex, such as a cone apex touching a facet
839 // (ThingiCSG/Basic/cube_cone_1.scad)
840
1/2
✓ Branch 1 taken 247128 times.
✗ Branch 2 not taken.
 247128 index_t v = mesh_.facets.vertex(f1,le);
841
7/8
✓ Branch 2 taken 247128 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 71303 times.
✓ Branch 5 taken 175825 times.
✓ Branch 8 taken 11674 times.
✓ Branch 9 taken 59629 times.
✓ Branch 10 taken 11674 times.
✓ Branch 11 taken 235454 times.
 247128 if(keep_vertex_[v] && !v_visited_[v]) {
842
1/2
✓ Branch 1 taken 11674 times.
✗ Branch 2 not taken.
 11674 vertices_.push_back(v);
843 11674 v_visited_[v] = true;
844 }
845 }
846 }
847 }
848
2/2
✓ Branch 2 taken 493832 times.
✓ Branch 3 taken 146172 times.
 786176 for(index_t v: vertices_) {
849 493832 v_visited_[v] = false;
850 }
851 }
852
853 // Get polylines
854 {
855 // Get all polylines starting from vertices with more than
856 // 2 incident halfedges.
857
2/2
✓ Branch 2 taken 493832 times.
✓ Branch 3 taken 146172 times.
 786176 for(index_t v: vertices_) {
858
3/4
✓ Branch 1 taken 493832 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 493828 times.
 493832 if(halfedges_.nb_halfedges_around_vertex(v) > 1) {
859 4 for(
860
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 index_t h=halfedges_.vertex_first_halfedge(v);
861
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)
862 ) {
863
1/2
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
 8 if(!h_visited_[h]) {
864 8 polylines_.begin_polyline();
865 8 index_t h2 = h;
866 do {
867
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]);
868 76 h_visited_[h2] = true;
869
1/2
✓ Branch 1 taken 76 times.
✗ Branch 2 not taken.
 76 polylines_.add_halfedge(h2);
870
1/2
✓ Branch 1 taken 76 times.
✗ Branch 2 not taken.
 76 h2 = halfedges_.next_along_polyline(h2);
871
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);
872
1/2
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
 8 polylines_.end_polyline();
873 }
874 }
875 }
876 }
877 // There can be also closed halfedge loops with no irregular vertex
878
2/2
✓ Branch 2 taken 493836 times.
✓ Branch 3 taken 146172 times.
 786180 for(index_t h: halfedges_) {
879
2/2
✓ Branch 2 taken 146370 times.
✓ Branch 3 taken 347466 times.
 493836 if(!h_visited_[h]) {
880 146370 polylines_.begin_polyline();
881 146370 index_t h2 = h;
882 do {
883
1/6
✗ Branch 2 not taken.
✓ Branch 3 taken 493760 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 493760 geo_debug_assert(!h_visited_[h2]);
884 493760 h_visited_[h2] = true;
885
1/2
✓ Branch 1 taken 493760 times.
✗ Branch 2 not taken.
 493760 polylines_.add_halfedge(h2);
886
1/2
✓ Branch 1 taken 493760 times.
✗ Branch 2 not taken.
 493760 h2 = halfedges_.next_along_polyline(h2);
887
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 493760 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 493760 geo_debug_assert(h2 != NO_INDEX);
888
2/2
✓ Branch 0 taken 347390 times.
✓ Branch 1 taken 146370 times.
 493760 } while(h2 != h);
889
1/2
✓ Branch 1 taken 146370 times.
✗ Branch 2 not taken.
 146370 polylines_.end_polyline();
890 }
891 }
892
2/2
✓ Branch 2 taken 493836 times.
✓ Branch 3 taken 146172 times.
 786180 for(index_t h: halfedges_) {
893 493836 h_visited_[h] = false;
894 }
895 }
896 146172 }
897
898 73086 void CoplanarFacets::mark_vertices_to_keep() {
899
3/4
✓ Branch 2 taken 73086 times.
✗ Branch 3 not taken.
✓ Branch 7 taken 73189 times.
✓ Branch 8 taken 73086 times.
 146275 for(index_t P: polylines_) {
900
1/2
✓ Branch 1 taken 73189 times.
✗ Branch 2 not taken.
 73189 index_t first_v = polylines_.first_vertex(P);
901
1/2
✓ Branch 1 taken 73189 times.
✗ Branch 2 not taken.
 73189 index_t last_v = polylines_.last_vertex(P);
902
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 73189 times.
 73189 if(first_v != last_v) {
903 ✗ keep_vertex_[first_v] = true;
904 ✗ keep_vertex_[last_v] = true;
905 }
906
1/2
✓ Branch 1 taken 73189 times.
✗ Branch 2 not taken.
 73189 index_t v1 = polylines_.prev_first_vertex(P);
907 73189 index_t v2 = NO_INDEX;
908 73189 index_t v3 = NO_INDEX;
909
4/6
✓ Branch 1 taken 73189 times.
✗ Branch 2 not taken.
✓ Branch 6 taken 246918 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 246918 times.
✓ Branch 11 taken 73189 times.
 320107 for(index_t h: polylines_.halfedges(P)) {
910
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 246918 times.
 246918 if(v1 == NO_INDEX) {
911 ✗ continue;
912 }
913
1/2
✓ Branch 1 taken 246918 times.
✗ Branch 2 not taken.
 246918 v2 = halfedges_.vertex(h,0);
914
1/2
✓ Branch 1 taken 246918 times.
✗ Branch 2 not taken.
 246918 v3 = halfedges_.vertex(h,1);
915
1/2
✓ Branch 1 taken 246918 times.
✗ Branch 2 not taken.
 246918 ExactPoint p1 = I_.exact_vertex(v1);
916
1/2
✓ Branch 1 taken 246918 times.
✗ Branch 2 not taken.
 246918 ExactPoint p2 = I_.exact_vertex(v2);
917
1/2
✓ Branch 1 taken 246918 times.
✗ Branch 2 not taken.
 246918 ExactPoint p3 = I_.exact_vertex(v3);
918
3/4
✓ Branch 1 taken 246918 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 238473 times.
✓ Branch 4 taken 8445 times.
 246918 if(!edges_are_colinear(p1,p2,p3)) {
919
1/2
✓ Branch 2 taken 238473 times.
✗ Branch 3 not taken.
 238473 keep_vertex_[v2] = true;
920 }
921 246918 v1 = v2;
922 246918 }
923 }
924 73086 }
925
926 ✗ void CoplanarFacets::save_borders(const std::string& filename) {
927 ✗ Mesh borders;
928 ✗ borders.vertices.set_dimension(2);
929 ✗ index_t cur_idx = 0;
930 ✗ for(index_t v: vertices_) {
931 ✗ vec3 p = mesh_.vertices.point(v);
932 ✗ vec2 q(p[u_],p[v_]);
933 ✗ borders.vertices.create_vertex(q.data());
934 ✗ v_idx_[v] = cur_idx;
935 ✗ ++cur_idx;
936 }
937
938 ✗ for(index_t h: halfedges_) {
939 ✗ index_t v1 = halfedges_.vertex(h,0);
940 ✗ index_t v2 = halfedges_.vertex(h,1);
941 ✗ v1 = v_idx_[v1];
942 ✗ v2 = v_idx_[v2];
943 ✗ geo_debug_assert(v1 != NO_INDEX);
944 ✗ geo_debug_assert(v2 != NO_INDEX);
945 ✗ borders.edges.create_edge(v1,v2);
946 }
947
948 ✗ Attribute<bool> selection(borders.vertices.attributes(), "selection");
949 ✗ for(index_t v: vertices_) {
950 ✗ geo_debug_assert(v_idx_[v] != NO_INDEX);
951 ✗ selection[v_idx_[v]] = keep_vertex_[v];
952 }
953 ✗ mesh_save(borders,filename);
954 ✗ }
955
956 ✗ void CoplanarFacets::save_facet_group(const std::string& filename) {
957 ✗ Mesh M;
958 ✗ Attribute<bool> keep_vertex(M.vertices.attributes(),"keep");
959 ✗ M.vertices.set_dimension(2);
960 ✗ for(index_t f: facets_) {
961 ✗ for(index_t lv=0; lv<3; ++lv) {
962 ✗ index_t v = mesh_.facets.vertex(f,lv);
963 ✗ v_idx_[v] = NO_INDEX;
964 }
965 }
966 ✗ for(index_t f: facets_) {
967 ✗ for(index_t lv=0; lv<3; ++lv) {
968 ✗ index_t v = mesh_.facets.vertex(f,lv);
969 ✗ if(v_idx_[v] == NO_INDEX) {
970 ✗ vec3 p = mesh_.vertices.point(v);
971 ✗ vec2 q(p[u_], p[v_]);
972 ✗ v_idx_[v] = M.vertices.create_vertex(q.data());
973 ✗ keep_vertex[v_idx_[v]] = keep_vertex_[v];
974 }
975 }
976 ✗ M.facets.create_triangle(
977 ✗ v_idx_[mesh_.facets.vertex(f,0)],
978 ✗ v_idx_[mesh_.facets.vertex(f,1)],
979 ✗ v_idx_[mesh_.facets.vertex(f,2)]
980 );
981 }
982
983 ✗ for(index_t f: facets_) {
984 ✗ for(index_t lv=0; lv<3; ++lv) {
985 ✗ index_t v = mesh_.facets.vertex(f,lv);
986 ✗ v_idx_[v] = NO_INDEX;
987 }
988 }
989
990 ✗ M.facets.connect();
991 ✗ mesh_save(M,filename);
992 ✗ }
993
994 8067 void CoplanarFacets::triangulate() {
995
996 // Compute 2D projected BBOX
997 8067 double umin = Numeric::max_float64();
998 8067 double vmin = Numeric::max_float64();
999 8067 double umax = -Numeric::max_float64();
1000 8067 double vmax = -Numeric::max_float64();
1001
2/2
✓ Branch 2 taken 58475 times.
✓ Branch 3 taken 8067 times.
 74609 for(index_t f: facets_) {
1002
2/2
✓ Branch 0 taken 175425 times.
✓ Branch 1 taken 58475 times.
 233900 for(index_t lv=0; lv<3; ++lv) {
1003
1/2
✓ Branch 1 taken 175425 times.
✗ Branch 2 not taken.
 175425 index_t vx = mesh_.facets.vertex(f,lv);
1004
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_];
1005
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_];
1006 175425 umin = std::min(umin, u);
1007 175425 umax = std::max(umax, u);
1008 175425 vmin = std::min(vmin, v);
1009 175425 vmax = std::max(vmax, v);
1010 }
1011 }
1012 8067 double d = std::max(umax-umin, vmax-vmin);
1013 8067 d *= 10.0;
1014 8067 d = std::max(d, 1.0);
1015 8067 umin-=d; vmin-=d; umax+=d; vmax+=d;
1016
1017 // Create CDT
1018
1/2
✓ Branch 1 taken 8067 times.
✗ Branch 2 not taken.
 8067 CDT.clear();
1019
1/2
✓ Branch 1 taken 8067 times.
✗ Branch 2 not taken.
 8067 CDT.create_enclosing_rectangle(umin, vmin, umax, vmax);
1020
1021
2/2
✓ Branch 2 taken 51859 times.
✓ Branch 3 taken 8067 times.
 67993 for(index_t v: vertices_) {
1022
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]) {
1023
1/2
✓ Branch 1 taken 43473 times.
✗ Branch 2 not taken.
 43473 ExactPoint P = I_.exact_vertex(v);
1024
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);
1025 43473 } else {
1026
1/2
✓ Branch 1 taken 8386 times.
✗ Branch 2 not taken.
 8386 v_idx_[v] = NO_INDEX;
1027 }
1028 }
1029
1030 // Insert constraints
1031
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_) {
1032 8170 vector<index_t> Pvertices;
1033
1/2
✓ Branch 1 taken 8170 times.
✗ Branch 2 not taken.
 8170 index_t v = polylines_.first_vertex(P);
1034
3/4
✓ Branch 2 taken 8170 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 7193 times.
✓ Branch 5 taken 977 times.
 8170 if(keep_vertex_[v]) {
1035
1/2
✓ Branch 1 taken 7193 times.
✗ Branch 2 not taken.
 7193 Pvertices.push_back(v);
1036 }
1037
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)) {
1038
1/2
✓ Branch 1 taken 51861 times.
✗ Branch 2 not taken.
 51861 v = halfedges_.vertex(h,1);
1039
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]) {
1040
1/2
✓ Branch 1 taken 43475 times.
✗ Branch 2 not taken.
 43475 Pvertices.push_back(v);
1041 }
1042 }
1043 8170 if(
1044
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) &&
1045 8170 Pvertices.size() != 0
1046 ) {
1047
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]);
1048 }
1049
1050
2/2
✓ Branch 1 taken 50668 times.
✓ Branch 2 taken 8170 times.
 58838 for(index_t i=0; i+1<Pvertices.size(); ++i) {
1051
2/4
✓ Branch 1 taken 50668 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 50668 times.
✗ Branch 5 not taken.
 50668 index_t v1 = v_idx_[Pvertices[i]];
1052
2/4
✓ Branch 1 taken 50668 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 50668 times.
✗ Branch 5 not taken.
 50668 index_t v2 = v_idx_[Pvertices[i+1]];
1053
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 50668 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 50668 geo_debug_assert(v1 != NO_INDEX);
1054
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 50668 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 50668 geo_debug_assert(v2 != NO_INDEX);
1055
1/2
✓ Branch 1 taken 50668 times.
✗ Branch 2 not taken.
 50668 CDT.insert_constraint(v1,v2,NO_INDEX);
1056 }
1057 8170 }
1058
1059
1/2
✓ Branch 1 taken 8067 times.
✗ Branch 2 not taken.
 8067 CDT.remove_external_triangles(true);
1060 8067 }
1061
1062 153785 bool CoplanarFacets::triangles_are_coplanar(
1063 const vec3& p1, const vec3& p2, const vec3& p3,
1064 const vec3& q1, const vec3& q2, const vec3& q3
1065 ) const {
1066 exact::vec3 N1 = cross(
1067
2/4
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 153785 times.
✗ Branch 5 not taken.
 307570 make_vec3<exact::vec3>(p1,p2), make_vec3<exact::vec3>(p1,p3)
1068
1/2
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
 153785 );
1069 exact::vec3 N2 = cross(
1070
2/4
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 153785 times.
✗ Branch 5 not taken.
 307570 make_vec3<exact::vec3>(q1,q2), make_vec3<exact::vec3>(q1,q3)
1071
1/2
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
 153785 );
1072
1073
9/14
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 35680 times.
✓ Branch 4 taken 118105 times.
✓ Branch 6 taken 35680 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 30302 times.
✓ Branch 9 taken 5378 times.
✓ Branch 11 taken 30302 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 30302 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 153785 times.
 153785 if(N1.x.sign() == ZERO && N1.y.sign() == ZERO && N1.z.sign() == ZERO) {
1074 ✗ std::cerr << std::endl;
1075 ✗ std::cerr << "degenerate triangle" << std::endl;
1076 ✗ std::cerr << "aligned: " << PCK::aligned_3d(p1,p2,p3) << std::endl;
1077 ✗ return false;
1078 }
1079
1080
9/14
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 36467 times.
✓ Branch 4 taken 117318 times.
✓ Branch 6 taken 36467 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 30991 times.
✓ Branch 9 taken 5476 times.
✓ Branch 11 taken 30991 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 30991 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 153785 times.
 153785 if(N2.x.sign() == ZERO && N2.y.sign() == ZERO && N2.z.sign() == ZERO) {
1081 ✗ std::cerr << std::endl;
1082 ✗ std::cerr << "degenerate triangle" << std::endl;
1083 ✗ std::cerr << "aligned: " << PCK::aligned_3d(q1,q2,q3) << std::endl;
1084 ✗ return false;
1085 }
1086
1087 // Tolerance for co-planarity test
1088
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 153785 times.
 153785 if(angle_tolerance_ != 0.0) {
1089 ✗ double threshold = cos(angle_tolerance_ * M_PI / 180.0);
1090 ✗ exact::scalar left = geo_sqr(dot(N1,N2));
1091 exact::scalar right =
1092 ✗ exact::scalar(threshold*threshold)*length2(N1)*length2(N2);
1093 ✗ return left > right;
1094 ✗ }
1095
1096 // Exact version
1097
1/2
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
 153785 exact::vec3 N12 = cross(N1,N2);
1098
11/14
✓ Branch 1 taken 153785 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 38733 times.
✓ Branch 4 taken 115052 times.
✓ Branch 6 taken 38733 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 32817 times.
✓ Branch 9 taken 5916 times.
✓ Branch 11 taken 32817 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 2491 times.
✓ Branch 14 taken 30326 times.
✓ Branch 15 taken 123459 times.
✓ Branch 16 taken 30326 times.
 153785 if((N12.x.sign()!=ZERO) || (N12.y.sign()!=ZERO) ||(N12.z.sign()!=ZERO)) {
1099 123459 return false;
1100 }
1101
1102 30326 return true;
1103 153785 }
1104
1105 /**************************************************************************/
1106
1107 246918 bool CoplanarFacets::edges_are_colinear(
1108 const ExactPoint& P1, const ExactPoint& P2, const ExactPoint& P3
1109 ) const {
1110
1111
1/2
✓ Branch 0 taken 246918 times.
✗ Branch 1 not taken.
 246918 if(angle_tolerance_ == 0.0) {
1112
1/2
✓ Branch 1 taken 246918 times.
✗ Branch 2 not taken.
 246918 return PCK::on_segment_3d(P2,P1,P3);
1113 }
1114
1115 ✗ ExactPoint UU = P1-P2;
1116 ✗ exact::vec3 U(UU.x, UU.y, UU.z);
1117 ✗ if(UU.w.sign() == NEGATIVE) {
1118 ✗ U.x.negate(); U.y.negate(); U.z.negate();
1119 }
1120 ✗ ExactPoint VV = P3-P2;
1121 ✗ exact::vec3 V(VV.x, VV.y, VV.z);
1122 ✗ if(VV.w.sign() == NEGATIVE) {
1123 ✗ V.x.negate(); V.y.negate(); V.z.negate();
1124 }
1125
1126 ✗ double threshold = cos(angle_tolerance_ * M_PI / 180.0);
1127
1128 ✗ exact::scalar left = dot(U,V);
1129
1130 ✗ if(left.sign() == POSITIVE) {
1131 ✗ return false;
1132 }
1133
1134 ✗ left = geo_sqr(left);
1135 exact::scalar right =
1136 ✗ exact::scalar(threshold*threshold)*length2(U)*length2(V);
1137
1138 ✗ return left > right;
1139 ✗ }
1140
1141 /**************************************************************************/
1142
1143 }
1144