Geogram Version 1.8.5
A programming library of geometric algorithms
Loading...
Searching...
No Matches
mesh_surface_intersection.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2000-2022 Inria
3 * All rights reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 *
8 * * Redistributions of source code must retain the above copyright notice,
9 * this list of conditions and the following disclaimer.
10 * * Redistributions in binary form must reproduce the above copyright notice,
11 * this list of conditions and the following disclaimer in the documentation
12 * and/or other materials provided with the distribution.
13 * * Neither the name of the ALICE Project-Team nor the names of its
14 * contributors may be used to endorse or promote products derived from this
15 * software without specific prior written permission.
16 *
17 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
18 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
21 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27 * POSSIBILITY OF SUCH DAMAGE.
28 *
29 * Contact: Bruno Levy
30 *
31 * https://www.inria.fr/fr/bruno-levy
32 *
33 * Inria,
34 * Domaine de Voluceau,
35 * 78150 Le Chesnay - Rocquencourt
36 * FRANCE
37 *
38 */
39
40#ifndef GEOGRAM_MESH_MESH_SURFACE_INTERSECTION
41#define GEOGRAM_MESH_MESH_SURFACE_INTERSECTION
42
44#include <geogram/mesh/mesh.h>
48#include <functional>
49
56namespace GEO {
57
58 /********************************************************************/
59
64 class GEOGRAM_API MeshSurfaceIntersection {
65 public:
68
77 void intersect();
78
79
80 void remove_external_shell();
81
82 void remove_internal_shells();
83
88 void set_verbose(bool x) {
89 verbose_ = x;
90 }
91
98 void set_delaunay(bool x) {
99 delaunay_ = x;
100 }
101
106 void set_approx_incircle(bool x) {
107 approx_incircle_ = x;
108 }
109
115 radial_sort_.set_approx_predicates(x);
116 }
117
124 detect_intersecting_neighbors_ = x;
125 }
126
134 void set_radial_sort(bool x) {
135 use_radial_sort_ = x;
136 }
137
145 void set_normalize(bool x) {
146 normalize_ = x;
147 }
148
149 void save_exact(const std::string& filename);
150
151 protected:
152
161 void lock() {
162 Process::acquire_spinlock(lock_);
163 }
164
168 void unlock() {
169 Process::release_spinlock(lock_);
170 }
171
182
194
201 return mesh_;
202 }
203
209 const Mesh& target_mesh() const {
210 return mesh_;
211 }
212
225 const Mesh& readonly_mesh() const {
226 return mesh_copy_;
227 }
228
238 );
239
240 void build_Weiler_model();
241
245 void mark_external_shell(vector<index_t>& on_external_shell);
246
247 index_t halfedge_vertex(index_t h, index_t dlv) const {
248 index_t f = h/3;
249 index_t lv = (h+dlv)%3;
250 return mesh_.facets.vertex(f,lv);
251 }
252
253 index_t alpha2(index_t h) const {
254 index_t t1 = h/3;
255 index_t t2 = mesh_.facet_corners.adjacent_facet(h);
256 if(t2 == index_t(-1)) {
257 return index_t(-1);
258 }
259 for(index_t h2: mesh_.facets.corners(t2)) {
260 if(mesh_.facet_corners.adjacent_facet(h2) == t1) {
261 return h2;
262 }
263 }
265 }
266
267 void sew2(index_t h1, index_t h2) {
269 halfedge_vertex(h1,0) == halfedge_vertex(h2,1)
270 );
272 halfedge_vertex(h2,0) == halfedge_vertex(h1,1)
273 );
274 index_t t1 = h1/3;
275 index_t t2 = h2/3;
276 mesh_.facet_corners.set_adjacent_facet(h1,t2);
277 mesh_.facet_corners.set_adjacent_facet(h2,t1);
278 }
279
280 index_t alpha3(index_t h) const {
281 return facet_corner_alpha3_[h];
282 }
283
284 index_t alpha3_facet(index_t f) const {
285 return alpha3(3*f)/3;
286 }
287
288 void sew3(index_t h1, index_t h2) {
290 halfedge_vertex(h1,0) == halfedge_vertex(h2,1)
291 );
293 halfedge_vertex(h2,0) == halfedge_vertex(h1,1)
294 );
295 facet_corner_alpha3_[h1] = h2;
296 facet_corner_alpha3_[h2] = h1;
297 }
298
299 void save_triangle(const std::string& name, index_t h) {
300 std::ofstream out(name + ".obj");
301 index_t v1 = halfedge_vertex(h,0);
302 index_t v2 = halfedge_vertex(h,1);
303 index_t v3 = halfedge_vertex(h,2);
304 out << "v " << vec3(mesh_.vertices.point_ptr(v1)) << std::endl;
305 out << "v " << vec3(mesh_.vertices.point_ptr(v2)) << std::endl;
306 out << "v " << vec3(mesh_.vertices.point_ptr(v3)) << std::endl;
307 out << "f 1 2 3" << std::endl;
308 }
309
310 void save_radial(
311 const std::string& name,
312 vector<index_t>::iterator b, vector<index_t>::iterator e
313 ) {
314 std::ofstream out(name + ".obj");
315 index_t v_ofs = 0;
316 for(vector<index_t>::iterator i=b; i!=e; ++i) {
317 index_t t = (*i/3);
318 index_t v1 = mesh_.facets.vertex(t,0);
319 index_t v2 = mesh_.facets.vertex(t,1);
320 index_t v3 = mesh_.facets.vertex(t,2);
321 out << "v "
322 << vec3(mesh_.vertices.point_ptr(v1)) << std::endl;
323 out << "v "
324 << vec3(mesh_.vertices.point_ptr(v2)) << std::endl;
325 out << "v "
326 << vec3(mesh_.vertices.point_ptr(v3)) << std::endl;
327 out << "f "
328 << v_ofs+1 << " " << v_ofs+2 << " " << v_ofs+3
329 << std::endl;
330 v_ofs += 3;
331 }
332 }
333
334 protected:
335
339 class GEOGRAM_API RadialSort {
340 public:
346 const MeshSurfaceIntersection& mesh
347 ) : mesh_(mesh),
348 approx_predicates_(false),
349 h_ref_(index_t(-1)),
350 degenerate_(false)
351 {
352 }
353
360 approx_predicates_ = x;
361 }
362
367 void init(index_t h_ref);
368
375 bool operator()(index_t h1, index_t h2) const;
376
383 bool degenerate() const {
384 return degenerate_;
385 }
386
394 static vec3E exact_direction(const vec3HE& p1, const vec3HE& p2);
395
396 protected:
397
407
415
416
417 public:
418 void test(index_t h1, index_t h2) {
419 (*this)(h1,h2);
420 Sign o_ref1 = h_orient(h_ref_, h1);
421 Sign o_ref2 = h_orient(h_ref_, h2);
422 Sign oN_ref1 = h_refNorient(h1);
423 Sign oN_ref2 = h_refNorient(h2);
424 Sign o_12 = h_orient(h1,h2);
425 std::cerr
426 << " o_ref1=" << int(o_ref1) << " o_ref2=" << int(o_ref2)
427 << " oN_ref1=" << int(oN_ref1) << " oN_ref2=" << int(oN_ref2)
428 << " o_12=" << int(o_12)
429 << std::endl;
430 }
431
432 private:
433 const MeshSurfaceIntersection& mesh_;
434 bool approx_predicates_;
435 index_t h_ref_; // ---reference halfedge
436 vec3E U_ref_; // -.
437 vec3E V_ref_; // +-reference basis
438 vec3E N_ref_; // -'
439 mutable vector< std::pair<index_t, Sign> > refNorient_cache_;
440 mutable bool degenerate_;
441 };
442
443
444 protected:
445 Process::spinlock lock_;
446 Mesh& mesh_;
447 Mesh mesh_copy_;
448 Attribute<const vec3HE*> vertex_to_exact_point_;
449 Attribute<index_t> facet_corner_alpha3_;
450 Attribute<bool> facet_corner_degenerate_;
451 std::map<vec3HE,index_t,vec3HELexicoCompare> exact_point_to_vertex_;
452 RadialSort radial_sort_;
453 bool verbose_;
454 bool delaunay_;
455 bool detect_intersecting_neighbors_;
456 bool approx_incircle_;
457 bool use_radial_sort_;
458 bool normalize_;
459 vec3 normalize_center_;
460 double normalize_radius_;
461 friend class MeshInTriangle;
462 };
463
464 /********************************************************************/
465
482 Mesh& M, std::function<bool(index_t)> eqn,
483 const std::string& attribute="", bool reorder=true
484 );
485
511 Mesh& M, const std::string& expr,
512 const std::string& attribute="", bool reorder=true
513 );
514
515 /*************************************************************************/
516
527 void GEOGRAM_API mesh_boolean_operation(
528 Mesh& result, Mesh& A, Mesh& B, const std::string& operation
529 );
530
540 void GEOGRAM_API mesh_union(Mesh& result, Mesh& A, Mesh& B);
541
551 void GEOGRAM_API mesh_intersection(Mesh& result, Mesh& A, Mesh& B);
552
562 void GEOGRAM_API mesh_difference(Mesh& result, Mesh& A, Mesh& B);
563
568 void GEOGRAM_API mesh_remove_intersections(Mesh& M, index_t max_iter = 3);
569
580 Mesh& M, index_t f1, index_t f2
581 );
582
583 /**************************************************************************/
584}
585
586#endif
587
#define geo_assert_not_reached
Sets a non reachable point in the program.
Definition assert.h:177
#define geo_debug_assert(x)
Verifies that a condition is met.
Definition assert.h:195
Generic mechanism for attributes.
Common include file, providing basic definitions. Should be included before anything else by all head...
Manages an attribute attached to a set of object.
Meshes a single triangle with the constraints that come from the intersections with the other triangl...
Sign h_orient(index_t h1, index_t h2) const
Computes the relative orientations of two halfedges.
Sign h_refNorient(index_t h2) const
Computes the normal orientation of a halfedge relative to h_ref.
void init(index_t h_ref)
Initializes radial sorting around a given halfedge.
void set_approx_predicates(bool x)
Specifies whether approximate predicates shoud be used.
bool degenerate() const
Tests if a degeneracy was encountered.
RadialSort(const MeshSurfaceIntersection &mesh)
RadialSort constructor.
static vec3E exact_direction(const vec3HE &p1, const vec3HE &p2)
Computes a vector of arbitrary length with its direction given by two points.
bool operator()(index_t h1, index_t h2) const
Compares two halfedges.
Computes surface intersections.
void set_verbose(bool x)
Display information while computing the intersection. Default is unset.
vec3HE exact_vertex(index_t v) const
Gets the exact point associated with a vertex.
void unlock()
Releases the lock associated with this mesh.
void set_detect_intersecting_neighbors(bool x)
detect and compute intersections between facets that share a facet or an edge. Set to false if input ...
void set_normalize(bool x)
Specifies whether coordinates should be normalized during computation. If set, original coordinates a...
void set_approx_incircle(bool x)
If set, then Delaunay mode uses approximated incircle predicate (else it uses exact arithmetics)....
const Mesh & readonly_mesh() const
Gets a copy of the initial mesh passed to the constructor.
const Mesh & target_mesh() const
Gets the target mesh.
void set_radial_sort(bool x)
Specifies whether surfaces should be duplicated and radial edges sorted in order to create the volume...
index_t find_or_create_exact_vertex(const vec3HE &p)
Finds or creates a vertex in the mesh, by exact coordinates.
void lock()
Acquires a lock on this mesh.
bool radial_sort(vector< index_t >::iterator b, vector< index_t >::iterator e)
Sorts a range of halfedges in radial order.
void mark_external_shell(vector< index_t > &on_external_shell)
Marks all the facets that are on the external shell.
Mesh & target_mesh()
Gets the target mesh.
void set_delaunay(bool x)
If set, compute constrained Delaunay triangulation in the intersected triangles. If there are interse...
void set_approx_radial_sort(bool x)
If set, do not use exact geometry for ordering triangles around radial edge. Default is unset.
Represents a mesh.
Definition mesh.h:2648
Vector with aligned memory allocation.
Definition memory.h:623
Exact predicates and constructs.
The class that represents a mesh.
int spinlock
Definition psm.h:125
Global Vorpaline namespace.
Definition algorithm.h:64
bool mesh_facets_have_intersection(Mesh &M, index_t f1, index_t f2)
Tests whether two mesh facets have a non-degenerate intersection.
void mesh_boolean_operation(Mesh &result, Mesh &A, Mesh &B, const std::string &operation)
Computes the union of two surface meshes.
void mesh_intersection(Mesh &result, Mesh &A, Mesh &B)
Computes the intersection of two surface meshes.
vecng< 3, Numeric::float64 > vec3
Represents points and vectors in 3d.
Definition geometry.h:65
void mesh_classify_intersections(Mesh &M, std::function< bool(index_t)> eqn, const std::string &attribute="", bool reorder=true)
Classifies the facets of the result of mesh_intersect_surface() based on a boolean function.
geo_index_t index_t
The type for storing and manipulating indices.
Definition numeric.h:287
Sign
Integer constants that represent the sign of a value.
Definition numeric.h:224
void mesh_union(Mesh &result, Mesh &A, Mesh &B)
Computes the union of two surface meshes.
void mesh_difference(Mesh &result, Mesh &A, Mesh &B)
Computes the difference of two surface meshes.
void mesh_remove_intersections(Mesh &M, index_t max_iter=3)
Attempts to make a surface mesh conformal by removing intersecting facets and re-triangulating the ho...
Function and classes for process manipulation.
3D vector in homogeneous coordinates with coordinates as arithmetic expansions.