| 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_fill_holes.h> | ||
| 41 | #include <geogram/mesh/mesh.h> | ||
| 42 | #include <geogram/mesh/mesh_repair.h> | ||
| 43 | #include <geogram/mesh/mesh_halfedges.h> | ||
| 44 | #include <geogram/mesh/mesh_io.h> | ||
| 45 | #include <geogram/mesh/index.h> | ||
| 46 | #include <geogram/basic/command_line.h> | ||
| 47 | #include <geogram/basic/logger.h> | ||
| 48 | |||
| 49 | #undef geo_debug_assert | ||
| 50 | #define geo_debug_assert(x) geo_assert(x) | ||
| 51 | |||
| 52 | namespace { | ||
| 53 | |||
| 54 | using namespace GEO; | ||
| 55 | |||
| 56 | /** | ||
| 57 | * \brief Checks whether a halfedge is | ||
| 58 | * incident to a vertex. | ||
| 59 | * \details Checks whether the origin of \p H | ||
| 60 | * is adjacent to \p v. In other words, | ||
| 61 | * returns true in one of the following | ||
| 62 | * configurations (the origin of \p H is | ||
| 63 | * denoted by \p x). | ||
| 64 | * \code | ||
| 65 | * | ||
| 66 | * H | ||
| 67 | * v--->x--->* | ||
| 68 | * | ||
| 69 | * H | ||
| 70 | * *--->x--->v | ||
| 71 | * | ||
| 72 | * \endcode | ||
| 73 | * \param[in] MH the mesh, wrapped with halfedge accessors | ||
| 74 | * \param[in] H the halfedge | ||
| 75 | * \param[in] v the index of the vertex | ||
| 76 | */ | ||
| 77 |
1/2✓ Branch 0 taken 120834 times.
✗ Branch 1 not taken.
|
120834 | bool halfedge_has_neighbor( |
| 78 | const MeshHalfedges& MH, | ||
| 79 | const MeshHalfedges::Halfedge& H, | ||
| 80 | index_t v | ||
| 81 | ) { | ||
| 82 | |||
| 83 | ✗ | geo_debug_assert(MH.halfedge_is_valid(H)); | |
| 84 | |||
| 85 | const Mesh& M = MH.mesh(); | ||
| 86 | index_t f = H.facet; | ||
| 87 | index_t c = H.corner; | ||
| 88 | |||
| 89 | { | ||
| 90 | index_t cnext = M.facets.next_corner_around_facet(f, c); | ||
| 91 |
1/2✓ Branch 0 taken 120834 times.
✗ Branch 1 not taken.
|
120834 | if(M.facet_corners.vertex(cnext) == v) { |
| 92 | return true; | ||
| 93 | } | ||
| 94 | } | ||
| 95 | |||
| 96 | { | ||
| 97 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 120834 times.
|
120834 | index_t cprev = M.facets.prev_corner_around_facet(f, c); |
| 98 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 120834 times.
|
120834 | if(M.facet_corners.vertex(cprev) == v) { |
| 99 | ✗ | return true; | |
| 100 | } | ||
| 101 | } | ||
| 102 | |||
| 103 | return false; | ||
| 104 | } | ||
| 105 | |||
| 106 | /** | ||
| 107 | * \brief Checks whether an halfedge exists between | ||
| 108 | * the two origins of two halfedges. | ||
| 109 | * \param[in] MH the mesh, wrapped with halfedge accessors | ||
| 110 | * \param[in] h1 first halfedge | ||
| 111 | * \param[in] h2 second halfedge | ||
| 112 | * \return true if an halfedge exists between the origins | ||
| 113 | * of \p h1 and \p h2, false otherwise | ||
| 114 | */ | ||
| 115 | 48019 | bool halfedge_exists_between_vertices( | |
| 116 | const MeshHalfedges& MH, | ||
| 117 | const MeshHalfedges::Halfedge& h1, | ||
| 118 | const MeshHalfedges::Halfedge& h2 | ||
| 119 | ) { | ||
| 120 | 48019 | index_t v2 = MH.mesh().facet_corners.vertex(h2.corner); | |
| 121 | 48019 | MeshHalfedges::Halfedge H = h1; | |
| 122 | do { | ||
| 123 |
1/2✓ Branch 1 taken 120834 times.
✗ Branch 2 not taken.
|
120834 | if(halfedge_has_neighbor(MH, H, v2)) { |
| 124 | return true; | ||
| 125 | } | ||
| 126 |
2/2✓ Branch 1 taken 73091 times.
✓ Branch 2 taken 47743 times.
|
120834 | if(!MH.move_to_prev_around_vertex(H)) { |
| 127 | break; | ||
| 128 | } | ||
| 129 |
2/2✓ Branch 0 taken 72815 times.
✓ Branch 1 taken 276 times.
|
73091 | } while(H != h1); |
| 130 | return false; | ||
| 131 | } | ||
| 132 | |||
| 133 | /** | ||
| 134 | * \brief Internal representation of a Hole. | ||
| 135 | * \details A Hole is an ordered sequence of Halfedge. | ||
| 136 | */ | ||
| 137 | typedef vector<MeshHalfedges::Halfedge> Hole; | ||
| 138 | |||
| 139 | /** | ||
| 140 | * \brief Computes a vector orthogonal to the border of a surface and in the | ||
| 141 | * tangent plane of the surface. | ||
| 142 | * \param[in] MH the mesh, wrapped with halfedge accessors | ||
| 143 | * \param[in] H the halfedge | ||
| 144 | * \return a 3d vector orthogonal to \p H, in the plange of the | ||
| 145 | * surface triangle incident to \p H and pointing towards the | ||
| 146 | * exterior of the surface. | ||
| 147 | */ | ||
| 148 | ✗ | vec3 border_normal( | |
| 149 | const MeshHalfedges& MH, | ||
| 150 | const MeshHalfedges::Halfedge& H | ||
| 151 | ) { | ||
| 152 | const Mesh& M = MH.mesh(); | ||
| 153 | ✗ | index_t c1 = H.corner; | |
| 154 | ✗ | index_t f = H.facet; | |
| 155 | index_t c2 = M.facets.next_corner_around_facet(f, c1); | ||
| 156 | vec3 E = M.facet_corners.point(c2) - M.facet_corners.point(c1); | ||
| 157 | ✗ | vec3 N = Geom::mesh_facet_normal(M, f); | |
| 158 | ✗ | return cross(E, N); | |
| 159 | } | ||
| 160 | |||
| 161 | /** | ||
| 162 | * \brief Splits a hole into two. | ||
| 163 | * \details This function is used recursively to triangulate the holes. | ||
| 164 | * \param[in] MH the mesh, wrapped with halfedge accessors | ||
| 165 | * \param[in] hole the hole to be split | ||
| 166 | * \param[out] hole1 one of the computed halves | ||
| 167 | * \param[out] hole2 the other half | ||
| 168 | * \param[in] use_normals if set, then couples of vertices | ||
| 169 | * (v1,v2) that have their border normals that match | ||
| 170 | * the new segment [v1,v2], are connected in priority. | ||
| 171 | * This improves the result on 788_raw_hand.off (but | ||
| 172 | * may cause some triangles overlap) | ||
| 173 | * \retval true on success | ||
| 174 | * \retval false otherwise (for instance, if no valid edge | ||
| 175 | * could be found to split the hole) | ||
| 176 | */ | ||
| 177 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 314 times.
|
314 | bool split_hole( |
| 178 | const MeshHalfedges& MH, const Hole& hole, | ||
| 179 | Hole& hole1, Hole& hole2, bool use_normals | ||
| 180 | ) { | ||
| 181 | |||
| 182 | // Step 0: compute normals to border | ||
| 183 | vector<vec3> N; | ||
| 184 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 314 times.
|
314 | if(use_normals) { |
| 185 | ✗ | N.assign(hole.size(), vec3(0.0, 0.0, 0.0)); | |
| 186 | ✗ | for(index_t i = 0; i < hole.size(); i++) { | |
| 187 | ✗ | index_t j = i + 1; | |
| 188 | ✗ | if(j == hole.size()) { | |
| 189 | j = 0; | ||
| 190 | } | ||
| 191 | ✗ | vec3 n = border_normal(MH, hole[i]); | |
| 192 | N[i] += n; | ||
| 193 | N[j] += n; | ||
| 194 | } | ||
| 195 | ✗ | for(index_t i = 0; i < N.size(); i++) { | |
| 196 | ✗ | N[i] = normalize(N[i]); | |
| 197 | } | ||
| 198 | } | ||
| 199 | |||
| 200 | // Step 1: compute total curve length and curvilinear abscissa | ||
| 201 | vector<double> s(hole.size()); | ||
| 202 | double cur_s = 0.0; | ||
| 203 | 314 | s[0] = cur_s; | |
| 204 |
2/2✓ Branch 0 taken 2670 times.
✓ Branch 1 taken 314 times.
|
3298 | for(index_t i = 1; i < hole.size(); i++) { |
| 205 | 2670 | const vec3& p1 = Geom::halfedge_vertex_from(MH.mesh(), hole[i - 1]); | |
| 206 | const vec3& p2 = Geom::halfedge_vertex_from(MH.mesh(), hole[i]); | ||
| 207 | 2670 | cur_s += length(p2 - p1); | |
| 208 | 2670 | s[i] = cur_s; | |
| 209 | } | ||
| 210 | 314 | const vec3& p1 = Geom::halfedge_vertex_from( | |
| 211 | MH.mesh(), hole[hole.size() - 1] | ||
| 212 | ); | ||
| 213 | const vec3& p2 = Geom::halfedge_vertex_from(MH.mesh(), hole[0]); | ||
| 214 | 314 | double total_length = cur_s + length(p2 - p1); | |
| 215 | |||
| 216 | // Step 2: find best pair to connect | ||
| 217 | double best_rij = Numeric::max_float64(); | ||
| 218 | signed_index_t v1 = -1; | ||
| 219 | signed_index_t v2 = -1; | ||
| 220 |
2/2✓ Branch 0 taken 2984 times.
✓ Branch 1 taken 314 times.
|
6596 | for(index_t i = 0; i < hole.size(); i++) { |
| 221 |
2/2✓ Branch 0 taken 48333 times.
✓ Branch 1 taken 2984 times.
|
102634 | for(index_t j = i + 2; j < hole.size(); j++) { |
| 222 | |||
| 223 | // Do not split using vertices | ||
| 224 | // already connected by an edge. | ||
| 225 | 314 | if( | |
| 226 |
5/6✓ Branch 0 taken 2356 times.
✓ Branch 1 taken 45977 times.
✓ Branch 2 taken 2042 times.
✓ Branch 3 taken 314 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 48019 times.
|
96352 | (i == 0 && j == hole.size() - 1) || |
| 227 |
1/2✓ Branch 1 taken 48019 times.
✗ Branch 2 not taken.
|
48019 | halfedge_exists_between_vertices(MH, hole[i], hole[j]) |
| 228 | ) { | ||
| 229 | continue; | ||
| 230 | } | ||
| 231 | |||
| 232 | double dsij = std::min( | ||
| 233 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 48019 times.
|
48019 | s[j] - s[i], total_length - (s[j] - s[i]) |
| 234 | ); | ||
| 235 | const vec3& pi = Geom::halfedge_vertex_from(MH.mesh(), hole[i]); | ||
| 236 | const vec3& pj = Geom::halfedge_vertex_from(MH.mesh(), hole[j]); | ||
| 237 | double dxij = length(pj - pi); | ||
| 238 | |||
| 239 | dsij = std::max(dsij, 1e-6); | ||
| 240 | dxij = std::max(dxij, 1e-6); | ||
| 241 | 48019 | double rij = dxij / dsij; | |
| 242 | |||
| 243 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 48019 times.
|
48019 | if(use_normals) { |
| 244 | const vec3& Pi = Geom::halfedge_vertex_from( | ||
| 245 | MH.mesh(), hole[i] | ||
| 246 | ); | ||
| 247 | const vec3& Pj = Geom::halfedge_vertex_from( | ||
| 248 | MH.mesh(), hole[j] | ||
| 249 | ); | ||
| 250 | ✗ | vec3 Dij = normalize(Pj - Pi); | |
| 251 | |||
| 252 | // between -1 (worse) and 1 (best) | ||
| 253 | double angle_factor = | ||
| 254 | ✗ | 0.5 * (dot(Dij, N[i]) - dot(Dij, N[j])); | |
| 255 | |||
| 256 | // between 0 (best) and 1 (worse) | ||
| 257 | ✗ | angle_factor = 0.5 * (1.0 - angle_factor); | |
| 258 | |||
| 259 | ✗ | rij *= angle_factor; | |
| 260 | } | ||
| 261 | |||
| 262 |
2/2✓ Branch 0 taken 1869 times.
✓ Branch 1 taken 46150 times.
|
48019 | if(rij < best_rij) { |
| 263 | best_rij = rij; | ||
| 264 | 1869 | v1 = signed_index_t(i); | |
| 265 | 1869 | v2 = signed_index_t(j); | |
| 266 | } | ||
| 267 | } | ||
| 268 | } | ||
| 269 | |||
| 270 |
2/4✓ Branch 0 taken 314 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 314 times.
✗ Branch 3 not taken.
|
314 | if(v1 == -1 || v2 == -1) { |
| 271 | return false; | ||
| 272 | } | ||
| 273 | |||
| 274 | // Now I do not think this can happen | ||
| 275 | // (to be checked) | ||
| 276 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 314 times.
|
314 | if(v2 < v1) { |
| 277 | std::swap(v1, v2); | ||
| 278 | } | ||
| 279 | |||
| 280 | // Step 3: copy the two "sub-holes" | ||
| 281 | hole1.clear(); | ||
| 282 | hole2.clear(); | ||
| 283 |
2/2✓ Branch 0 taken 2984 times.
✓ Branch 1 taken 314 times.
|
6282 | for(signed_index_t i = 0; i < signed_index_t(hole.size()); i++) { |
| 284 |
4/4✓ Branch 0 taken 2313 times.
✓ Branch 1 taken 671 times.
✓ Branch 2 taken 897 times.
✓ Branch 3 taken 1416 times.
|
2984 | if(i <= v1 || i >= v2) { |
| 285 |
2/2✓ Branch 0 taken 552 times.
✓ Branch 1 taken 1016 times.
|
1568 | hole1.push_back(hole[i]); |
| 286 | } | ||
| 287 |
4/4✓ Branch 0 taken 2627 times.
✓ Branch 1 taken 357 times.
✓ Branch 2 taken 2044 times.
✓ Branch 3 taken 583 times.
|
2984 | if(i >= v1 && i <= v2) { |
| 288 |
2/2✓ Branch 0 taken 931 times.
✓ Branch 1 taken 1113 times.
|
2044 | hole2.push_back(hole[i]); |
| 289 | } | ||
| 290 | } | ||
| 291 | |||
| 292 | return true; | ||
| 293 | } | ||
| 294 | |||
| 295 | /** | ||
| 296 | * \brief Triangulates a hole. | ||
| 297 | * \param[in] MH the mesh, wrapped with halfedge accessors | ||
| 298 | * \param[in] hole the hole, represented by a vector of halfedges | ||
| 299 | * \param[out] triangles the generated triangles | ||
| 300 | * \param[in] use_normals if set, then couples of vertices | ||
| 301 | * (v1,v2) that have their border normals that match | ||
| 302 | * the new segment [v1,v2], are connected in priority. | ||
| 303 | * This improves the result on 788_raw_hand.off (but | ||
| 304 | * may cause some triangles overlap) | ||
| 305 | * \param[in] clear if set, \p triangles is cleared | ||
| 306 | * \retval true on success | ||
| 307 | * \retval false otherwise (for instance, if no valid edge | ||
| 308 | * could be found to split the hole) | ||
| 309 | */ | ||
| 310 | 358 | bool triangulate_hole_loop_splitting( | |
| 311 | const MeshHalfedges& MH, const Hole& hole, | ||
| 312 | vector<trindex>& triangles, bool use_normals, | ||
| 313 | bool clear = true | ||
| 314 | ) { | ||
| 315 | bool ok = true; | ||
| 316 |
2/2✓ Branch 0 taken 8 times.
✓ Branch 1 taken 350 times.
|
358 | if(clear) { |
| 317 | triangles.clear(); | ||
| 318 | } | ||
| 319 |
2/2✓ Branch 0 taken 183 times.
✓ Branch 1 taken 175 times.
|
358 | if(hole.size() <= 3) { |
| 320 |
1/2✓ Branch 0 taken 183 times.
✗ Branch 1 not taken.
|
183 | if(hole.size() == 3) { |
| 321 | trindex T( | ||
| 322 | 183 | MH.mesh().facet_corners.vertex(hole[0].corner), | |
| 323 |
2/2✓ Branch 0 taken 167 times.
✓ Branch 1 taken 16 times.
|
183 | MH.mesh().facet_corners.vertex(hole[1].corner), |
| 324 |
2/2✓ Branch 0 taken 167 times.
✓ Branch 1 taken 16 times.
|
183 | MH.mesh().facet_corners.vertex(hole[2].corner), |
| 325 | trindex::KEEP_ORDER | ||
| 326 | ); | ||
| 327 |
2/2✓ Branch 0 taken 167 times.
✓ Branch 1 taken 16 times.
|
183 | triangles.push_back(T); |
| 328 | } | ||
| 329 | } else { | ||
| 330 | Hole hole1, hole2; | ||
| 331 |
1/2✓ Branch 1 taken 175 times.
✗ Branch 2 not taken.
|
175 | ok = split_hole(MH, hole, hole1, hole2, use_normals); |
| 332 |
3/6✓ Branch 0 taken 175 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 175 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 175 times.
✗ Branch 6 not taken.
|
175 | ok = ok && triangulate_hole_loop_splitting( |
| 333 | MH, hole1, triangles, use_normals, false | ||
| 334 | ); | ||
| 335 |
2/4✓ Branch 1 taken 175 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 175 times.
|
175 | ok = ok && triangulate_hole_loop_splitting( |
| 336 | MH, hole2, triangles, use_normals, false | ||
| 337 | ); | ||
| 338 | } | ||
| 339 | 358 | return ok; | |
| 340 | } | ||
| 341 | |||
| 342 | /************************************************************************/ | ||
| 343 | |||
| 344 | /** | ||
| 345 | * \brief Computes the score obtained when generating a triangle | ||
| 346 | * in the ear cutting algorithm triangulate_hole_ear_cutting(). | ||
| 347 | * \param[in] M the mesh | ||
| 348 | * \param[in] T1 first triangle, encoded as a vertex indices triplet | ||
| 349 | * \param[in] T2 second triangle, encoded as a vertex indices triplet | ||
| 350 | * \return the score obtained when generating an ear from \p T1 to \p T2 | ||
| 351 | */ | ||
| 352 | ✗ | double ear_score( | |
| 353 | const Mesh& M, | ||
| 354 | const trindex& T1, | ||
| 355 | const trindex& T2 | ||
| 356 | ) { | ||
| 357 | ✗ | geo_debug_assert(T1.indices[1] == T2.indices[0]); | |
| 358 | ✗ | const vec3& p10 = M.vertices.point(T1.indices[0]); | |
| 359 | const vec3& p11 = M.vertices.point(T1.indices[1]); | ||
| 360 | ✗ | const vec3& p12 = M.vertices.point(T1.indices[2]); | |
| 361 | const vec3& p20 = M.vertices.point(T2.indices[0]); | ||
| 362 | ✗ | const vec3& p21 = M.vertices.point(T2.indices[1]); | |
| 363 | ✗ | const vec3& p22 = M.vertices.point(T2.indices[2]); | |
| 364 | ✗ | vec3 n = normalize( | |
| 365 | ✗ | Geom::triangle_normal(p10, p11, p12) + | |
| 366 | ✗ | Geom::triangle_normal(p20, p21, p22) | |
| 367 | ); | ||
| 368 | ✗ | vec3 a = normalize(p11 - p10); | |
| 369 | ✗ | vec3 b = normalize(p21 - p20); | |
| 370 | ✗ | return -::atan2(dot(n, cross(a, b)), dot(a, b)); | |
| 371 | } | ||
| 372 | |||
| 373 | /** | ||
| 374 | * \brief Triangulates a hole using the ear cutting algorithm. | ||
| 375 | * \param[in] MH the mesh, wrapped with halfedge accessors | ||
| 376 | * \param[in] hole_in the hole, represented by a vector of halfedges | ||
| 377 | * \param[out] triangles the generated triangles | ||
| 378 | * \param[out] clear if set, \p triangles is cleared | ||
| 379 | */ | ||
| 380 | ✗ | void triangulate_hole_ear_cutting( | |
| 381 | const MeshHalfedges& MH, const Hole& hole_in, | ||
| 382 | vector<trindex>& triangles, bool clear = true | ||
| 383 | ) { | ||
| 384 | ✗ | if(clear) { | |
| 385 | triangles.clear(); | ||
| 386 | } | ||
| 387 | ✗ | if(hole_in.size() <= 3) { | |
| 388 | ✗ | if(hole_in.size() == 3) { | |
| 389 | trindex T( | ||
| 390 | ✗ | MH.mesh().facet_corners.vertex(hole_in[0].corner), | |
| 391 | ✗ | MH.mesh().facet_corners.vertex(hole_in[1].corner), | |
| 392 | ✗ | MH.mesh().facet_corners.vertex(hole_in[2].corner), | |
| 393 | trindex::KEEP_ORDER | ||
| 394 | ); | ||
| 395 | ✗ | triangles.push_back(T); | |
| 396 | } | ||
| 397 | } else { | ||
| 398 | const Mesh& M = MH.mesh(); | ||
| 399 | |||
| 400 | // Step 1: convert hole into easier-to-manipulate representation. | ||
| 401 | vector<trindex> hole; | ||
| 402 | ✗ | hole.reserve(hole_in.size()); | |
| 403 | ✗ | for(index_t i = 0; i < hole_in.size(); i++) { | |
| 404 | const MeshHalfedges::Halfedge& H = hole_in[i]; | ||
| 405 | ✗ | geo_debug_assert(H.facet != MeshHalfedges::Halfedge::NO_FACET); | |
| 406 | ✗ | index_t c = H.corner; | |
| 407 | index_t v1 = M.facet_corners.vertex(c); | ||
| 408 | c = M.facets.next_corner_around_facet(H.facet, c); | ||
| 409 | index_t v2 = M.facet_corners.vertex(c); | ||
| 410 | c = M.facets.next_corner_around_facet(H.facet, c); | ||
| 411 | index_t v3 = M.facet_corners.vertex(c); | ||
| 412 | ✗ | hole.push_back(trindex(v1, v2, v3, trindex::KEEP_ORDER)); | |
| 413 | } | ||
| 414 | |||
| 415 | // Step 2: ear cutting | ||
| 416 | ✗ | while(hole.size() > 3) { | |
| 417 | signed_index_t best_i1 = -1; | ||
| 418 | double best_score = Numeric::min_float64(); | ||
| 419 | // TODO: take existing edges into account. | ||
| 420 | ✗ | for(index_t i1 = 0; i1 < hole.size(); i1++) { | |
| 421 | ✗ | index_t i2 = i1 + 1; | |
| 422 | ✗ | if(i2 == hole.size()) { | |
| 423 | i2 = 0; | ||
| 424 | } | ||
| 425 | ✗ | double score = ear_score(M, hole[i1], hole[i2]); | |
| 426 | ✗ | if(score > best_score) { | |
| 427 | ✗ | best_i1 = signed_index_t(i1); | |
| 428 | best_score = score; | ||
| 429 | } | ||
| 430 | } | ||
| 431 | ✗ | geo_assert(best_i1 != -1); | |
| 432 | ✗ | index_t best_i2 = index_t(best_i1) + 1; | |
| 433 | ✗ | if(best_i2 == hole.size()) { | |
| 434 | best_i2 = 0; | ||
| 435 | } | ||
| 436 | const trindex& T1 = hole[best_i1]; | ||
| 437 | const trindex& T2 = hole[best_i2]; | ||
| 438 | ✗ | geo_debug_assert(T1.indices[1] == T2.indices[0]); | |
| 439 | trindex T( | ||
| 440 | ✗ | T1.indices[0], T2.indices[1], T1.indices[1], | |
| 441 | trindex::KEEP_ORDER | ||
| 442 | ); | ||
| 443 | ✗ | hole[best_i1] = T; | |
| 444 | hole.erase(hole.begin() + std::ptrdiff_t(best_i2)); | ||
| 445 | ✗ | triangles.push_back(T); | |
| 446 | } | ||
| 447 | |||
| 448 | // Step 3: last triangle | ||
| 449 | ✗ | geo_assert(hole.size() == 3); | |
| 450 | trindex T( | ||
| 451 | hole[0].indices[0], | ||
| 452 | hole[1].indices[0], | ||
| 453 | hole[2].indices[0], | ||
| 454 | trindex::KEEP_ORDER | ||
| 455 | ✗ | ); | |
| 456 | ✗ | triangles.push_back(T); | |
| 457 | } | ||
| 458 | ✗ | } | |
| 459 | |||
| 460 | /************************************************************************/ | ||
| 461 | |||
| 462 | /** | ||
| 463 | * \brief Computes the area of a hole, i.e. the area of the generated | ||
| 464 | * triangles that will fill the hole. | ||
| 465 | * \param[in] M the mesh | ||
| 466 | * \param[in] triangles the triangles that will fill the hole | ||
| 467 | * \return the summed ares of the triangles in \p triangles | ||
| 468 | */ | ||
| 469 | 8 | double hole_area( | |
| 470 | const Mesh& M, const vector<trindex>& triangles | ||
| 471 | ) { | ||
| 472 | double result = 0.0; | ||
| 473 |
2/2✓ Branch 0 taken 183 times.
✓ Branch 1 taken 8 times.
|
199 | for(index_t t = 0; t < triangles.size(); t++) { |
| 474 | 183 | index_t i = triangles[t].indices[0]; | |
| 475 | 183 | index_t j = triangles[t].indices[1]; | |
| 476 | 183 | index_t k = triangles[t].indices[2]; | |
| 477 | const vec3& p1 = M.vertices.point(i); | ||
| 478 | const vec3& p2 = M.vertices.point(j); | ||
| 479 | const vec3& p3 = M.vertices.point(k); | ||
| 480 | 183 | result += Geom::triangle_area(p1, p2, p3); | |
| 481 | } | ||
| 482 | 8 | return result; | |
| 483 | } | ||
| 484 | |||
| 485 | /** | ||
| 486 | * \brief Strategy used to fill the holes. | ||
| 487 | */ | ||
| 488 | enum HoleFilling { | ||
| 489 | LOOP_SPLIT, /**< Splits loops by generating small segments */ | ||
| 490 | NLOOP_SPLIT, /**< Takes normals into account */ | ||
| 491 | EAR_CUT /**< Uses the "ear cutting" strategy */ | ||
| 492 | }; | ||
| 493 | |||
| 494 | /************************************************************************/ | ||
| 495 | |||
| 496 | /* // commented-out for now, see issue #72 | ||
| 497 | * \brief Removes all the facets of a mesh that are | ||
| 498 | * on a bridge. | ||
| 499 | * \details A facet is said to be on a bridge if it is | ||
| 500 | * incident to a border and if when turning around the | ||
| 501 | * border it is encountered more than once. | ||
| 502 | * | ||
| 503 | * | ||
| 504 | * | ||
| 505 | void remove_bridges(Mesh& M) { | ||
| 506 | MeshHalfedges MH(M); | ||
| 507 | vector<bool> corner_is_visited(M.facet_corners.nb(),false); | ||
| 508 | vector<index_t> f_status(M.facets.nb(),0); | ||
| 509 | index_t f_stamp=1; | ||
| 510 | const index_t BRIDGE = NO_INDEX; | ||
| 511 | |||
| 512 | for(index_t f: M.facets) { | ||
| 513 | for(index_t c: M.facets.corners(f)) { | ||
| 514 | if( | ||
| 515 | M.facet_corners.adjacent_facet(c) == NO_FACET && | ||
| 516 | !corner_is_visited[c] | ||
| 517 | ) { | ||
| 518 | MeshHalfedges::Halfedge first(f, c); | ||
| 519 | MeshHalfedges::Halfedge H(f, c); | ||
| 520 | do { | ||
| 521 | corner_is_visited[H.corner] = true; | ||
| 522 | MH.move_to_next_around_facet(H); | ||
| 523 | while(MH.move_to_next_around_vertex(H)) { | ||
| 524 | if(f_status[H.facet] == f_stamp) { | ||
| 525 | f_status[H.facet] = BRIDGE; | ||
| 526 | } else if( | ||
| 527 | f_status[H.facet] != BRIDGE && | ||
| 528 | f_status[H.facet] != f_stamp) { | ||
| 529 | f_status[H.facet] = f_stamp; | ||
| 530 | } | ||
| 531 | } | ||
| 532 | } while(H != first); | ||
| 533 | ++f_stamp; | ||
| 534 | } | ||
| 535 | } | ||
| 536 | } | ||
| 537 | index_t nb_bridges = 0; | ||
| 538 | for(index_t f: M.facets) { | ||
| 539 | if(f_status[f] == BRIDGE) { | ||
| 540 | ++nb_bridges; | ||
| 541 | } else { | ||
| 542 | f_status[f] = 0; | ||
| 543 | } | ||
| 544 | } | ||
| 545 | if(nb_bridges != 0) { | ||
| 546 | M.facets.delete_elements(f_status); | ||
| 547 | Logger::out("Bridges") | ||
| 548 | << "Removed " << nb_bridges << " bridge(s)" | ||
| 549 | << std::endl; | ||
| 550 | } | ||
| 551 | } | ||
| 552 | */ | ||
| 553 | } | ||
| 554 | |||
| 555 | /****************************************************************************/ | ||
| 556 | |||
| 557 | namespace GEO { | ||
| 558 | |||
| 559 | 8 | void fill_holes( | |
| 560 | Mesh& M, double max_area, index_t max_edges, bool repair | ||
| 561 | ) { | ||
| 562 | |||
| 563 |
1/2✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
|
8 | if(max_area == 0.0 || max_edges == 0) { |
| 564 | 5 | return; | |
| 565 | } | ||
| 566 | |||
| 567 | // remove_bridges(M); // commented out for now, see issue #72 | ||
| 568 | |||
| 569 | MeshHalfedges MH(M); | ||
| 570 | |||
| 571 | vector<Hole> holes; | ||
| 572 | |||
| 573 | index_t nb_filled_holes = 0; | ||
| 574 | index_t nb_skipped_by_edges = 0; | ||
| 575 | index_t nb_skipped_by_area = 0; | ||
| 576 | index_t nb_could_not_fill = 0; | ||
| 577 | |||
| 578 | // Step 1: detect holes | ||
| 579 | { | ||
| 580 | 8 | vector<bool> corner_is_visited(M.facet_corners.nb(), false); | |
| 581 |
2/2✓ Branch 0 taken 136340 times.
✓ Branch 1 taken 8 times.
|
136348 | for(index_t f: M.facets) { |
| 582 |
2/2✓ Branch 0 taken 409020 times.
✓ Branch 1 taken 136340 times.
|
545360 | for(index_t c: M.facets.corners(f)) { |
| 583 | if( | ||
| 584 |
4/4✓ Branch 0 taken 2256 times.
✓ Branch 1 taken 406764 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 2246 times.
|
409020 | M.facet_corners.adjacent_facet(c) == NO_FACET && |
| 585 | !corner_is_visited[c] | ||
| 586 | ) { | ||
| 587 | 10 | holes.push_back(Hole()); | |
| 588 | MeshHalfedges::Halfedge first(f, c); | ||
| 589 | MeshHalfedges::Halfedge H(f, c); | ||
| 590 | do { | ||
| 591 |
2/2✓ Branch 0 taken 2202 times.
✓ Branch 1 taken 54 times.
|
2256 | holes.rbegin()->push_back(H); |
| 592 |
1/2✓ Branch 1 taken 2256 times.
✗ Branch 2 not taken.
|
2256 | corner_is_visited[H.corner] = true; |
| 593 |
1/2✓ Branch 1 taken 2256 times.
✗ Branch 2 not taken.
|
2256 | MH.move_to_next_around_border(H); |
| 594 |
2/2✓ Branch 0 taken 2246 times.
✓ Branch 1 taken 10 times.
|
2256 | } while(H != first); |
| 595 | |||
| 596 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
|
10 | if(holes.rbegin()->size() > max_edges) { |
| 597 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2 | ++nb_skipped_by_edges; |
| 598 | holes.pop_back(); | ||
| 599 | } | ||
| 600 | } | ||
| 601 | } | ||
| 602 | } | ||
| 603 | } | ||
| 604 | |||
| 605 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 3 times.
|
8 | if(holes.size() == 0) { |
| 606 | return; | ||
| 607 | } | ||
| 608 | |||
| 609 |
2/4✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
|
3 | Logger::out("FillHoles") << "Found " << holes.size() |
| 610 | << " holes" << std::endl; | ||
| 611 | |||
| 612 | HoleFilling algo = LOOP_SPLIT; | ||
| 613 |
2/4✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
|
3 | std::string algo_name = CmdLine::get_arg("algo:hole_filling"); |
| 614 |
1/2✓ Branch 0 taken 3 times.
✗ Branch 1 not taken.
|
3 | if(algo_name == "loop_split") { |
| 615 | algo = LOOP_SPLIT; | ||
| 616 | ✗ | } else if(algo_name == "Nloop_split") { | |
| 617 | algo = NLOOP_SPLIT; | ||
| 618 | ✗ | } else if(algo_name == "ear_cut") { | |
| 619 | algo = EAR_CUT; | ||
| 620 | } else { | ||
| 621 | ✗ | Logger::warn("FillHoles") | |
| 622 | << algo_name << ": no such hole filling method" | ||
| 623 | << std::endl; | ||
| 624 | ✗ | Logger::warn("FillHoles") | |
| 625 | << "falling back to \'loop_split\'" | ||
| 626 | << std::endl; | ||
| 627 | } | ||
| 628 | |||
| 629 | |||
| 630 |
2/2✓ Branch 0 taken 8 times.
✓ Branch 1 taken 3 times.
|
22 | for(index_t i = 0; i < holes.size(); i++) { |
| 631 | vector<trindex> triangles; | ||
| 632 | bool ok = true; | ||
| 633 |
1/3✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
|
8 | switch(algo) { |
| 634 | case LOOP_SPLIT: | ||
| 635 |
1/2✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
|
8 | ok = triangulate_hole_loop_splitting( |
| 636 | MH, holes[i], triangles, false | ||
| 637 | ); | ||
| 638 | break; | ||
| 639 | case NLOOP_SPLIT: | ||
| 640 | ✗ | ok = triangulate_hole_loop_splitting( | |
| 641 | MH, holes[i], triangles, true | ||
| 642 | ); | ||
| 643 | break; | ||
| 644 | case EAR_CUT: | ||
| 645 | ✗ | triangulate_hole_ear_cutting(MH, holes[i], triangles); | |
| 646 | break; | ||
| 647 | } | ||
| 648 | |||
| 649 |
1/2✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
|
8 | if(ok) { |
| 650 |
2/2✓ Branch 1 taken 7 times.
✓ Branch 2 taken 1 times.
|
8 | if(hole_area(M, triangles) < max_area) { |
| 651 |
2/2✓ Branch 0 taken 7 times.
✓ Branch 1 taken 7 times.
|
21 | for(index_t j = 0; j < triangles.size(); j++) { |
| 652 |
1/2✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
|
7 | M.facets.create_triangle( |
| 653 | 7 | triangles[j].indices[2], | |
| 654 | 7 | triangles[j].indices[1], | |
| 655 |
1/2✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
|
7 | triangles[j].indices[0] |
| 656 | ); | ||
| 657 | } | ||
| 658 | 7 | ++nb_filled_holes; | |
| 659 | } else { | ||
| 660 | 1 | ++nb_skipped_by_area; | |
| 661 | } | ||
| 662 | } else { | ||
| 663 | ✗ | ++nb_could_not_fill; | |
| 664 | } | ||
| 665 | |||
| 666 | } | ||
| 667 | |||
| 668 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 2 times.
|
3 | if(nb_skipped_by_area != 0) { |
| 669 |
2/4✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
|
2 | Logger::out("FillHoles") |
| 670 | << "Skipped " << nb_skipped_by_area | ||
| 671 | << " holes (area too large)" << std::endl; | ||
| 672 | } | ||
| 673 | |||
| 674 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 1 times.
|
3 | if(nb_skipped_by_edges != 0) { |
| 675 |
2/4✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
|
4 | Logger::out("FillHoles") |
| 676 | << "Skipped " << nb_skipped_by_edges | ||
| 677 | << " holes (too many edges)" << std::endl; | ||
| 678 | } | ||
| 679 | |||
| 680 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 3 times.
|
3 | if(nb_could_not_fill != 0) { |
| 681 | ✗ | Logger::out("FillHoles") | |
| 682 | << "Skipped " << nb_could_not_fill | ||
| 683 | << " holes (could not fill)" << std::endl; | ||
| 684 | } | ||
| 685 | |||
| 686 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 1 times.
|
3 | if(nb_filled_holes != 0 && repair) { |
| 687 | // Needed because we may generate zero-length edges | ||
| 688 | // and zero-area facets that need to be eliminated. | ||
| 689 | // Note: this also reconstructs the connections between the facets. | ||
| 690 | MeshRepairMode mode = MESH_REPAIR_DEFAULT; | ||
| 691 |
1/2✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
|
2 | mesh_repair(M, mode); |
| 692 | } | ||
| 693 | } | ||
| 694 | |||
| 695 | /** | ||
| 696 | * \brief Tessellates a hole. | ||
| 697 | * \param[in] MH the MeshHalfHeddes the hole belongs to. | ||
| 698 | * \param[in] H the hole. | ||
| 699 | * \param[in] max_nb_vertices maximum number of vertices in | ||
| 700 | * the new facets to create. | ||
| 701 | * \param[in] copy_facet_attrib an optional facet index. If | ||
| 702 | * specified, all the attributes of this facet will be copied | ||
| 703 | * to the created facets. | ||
| 704 | */ | ||
| 705 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 417 times.
|
417 | static void tessellate_hole( |
| 706 | MeshHalfedges& MH, Hole& H, index_t max_nb_vertices, | ||
| 707 | index_t copy_facet_attrib = NO_INDEX | ||
| 708 | ) { | ||
| 709 | Mesh& M = MH.mesh(); | ||
| 710 |
2/2✓ Branch 0 taken 278 times.
✓ Branch 1 taken 139 times.
|
417 | if(H.size() <= max_nb_vertices) { |
| 711 | 278 | index_t f = M.facets.create_polygon(H.size()); | |
| 712 |
2/2✓ Branch 0 taken 834 times.
✓ Branch 1 taken 278 times.
|
2224 | FOR(i,H.size()) { |
| 713 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 834 times.
|
834 | index_t v = M.facet_corners.vertex(H[i].corner); |
| 714 | M.facets.set_vertex(f,i,v); | ||
| 715 | } | ||
| 716 |
1/2✓ Branch 0 taken 278 times.
✗ Branch 1 not taken.
|
278 | if(copy_facet_attrib != NO_INDEX) { |
| 717 | 278 | M.facets.attributes().copy_item(f, copy_facet_attrib); | |
| 718 | } | ||
| 719 | } else { | ||
| 720 | Hole H1,H2; | ||
| 721 |
1/2✓ Branch 1 taken 139 times.
✗ Branch 2 not taken.
|
139 | split_hole(MH,H,H1,H2,false); |
| 722 |
1/2✓ Branch 1 taken 139 times.
✗ Branch 2 not taken.
|
139 | tessellate_hole(MH,H1,max_nb_vertices,copy_facet_attrib); |
| 723 |
1/2✓ Branch 1 taken 139 times.
✗ Branch 2 not taken.
|
139 | tessellate_hole(MH,H2,max_nb_vertices,copy_facet_attrib); |
| 724 | } | ||
| 725 | 417 | } | |
| 726 | |||
| 727 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | void tessellate_facets( |
| 728 | Mesh& M, index_t max_nb_vertices | ||
| 729 | ) { | ||
| 730 | MeshHalfedges MH(M); | ||
| 731 |
1/4✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
|
6 | vector<index_t> delete_f(M.facets.nb(),0); |
| 732 |
2/2✓ Branch 0 taken 143 times.
✓ Branch 1 taken 6 times.
|
149 | for(index_t f: M.facets) { |
| 733 |
2/2✓ Branch 0 taken 139 times.
✓ Branch 1 taken 4 times.
|
143 | if(M.facets.nb_vertices(f) > max_nb_vertices) { |
| 734 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 139 times.
|
139 | delete_f[f] = 1; |
| 735 | Hole h; | ||
| 736 |
2/2✓ Branch 0 taken 556 times.
✓ Branch 1 taken 139 times.
|
695 | for(index_t c: M.facets.corners(f)) { |
| 737 |
0/2✗ Branch 0 not taken.
✗ Branch 1 not taken.
|
556 | h.push_back(MeshHalfedges::Halfedge(f,c)); |
| 738 | } | ||
| 739 |
1/2✓ Branch 1 taken 139 times.
✗ Branch 2 not taken.
|
139 | tessellate_hole(MH, h, max_nb_vertices, f); |
| 740 | } | ||
| 741 | } | ||
| 742 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | delete_f.resize(M.facets.nb()); |
| 743 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | M.facets.delete_elements(delete_f); |
| 744 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | M.facets.connect(); |
| 745 |
1/2✓ Branch 0 taken 6 times.
✗ Branch 1 not taken.
|
6 | if(max_nb_vertices == 3) { |
| 746 | M.facets.is_simplicial(); | ||
| 747 | } | ||
| 748 | 6 | } | |
| 749 | |||
| 750 | } | ||
| 751 |