Geogram Version 1.8.5
A programming library of geometric algorithms
Loading...
Searching...
No Matches
mesh_surface_intersection_internal.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_INTERNAL
41#define GEOGRAM_MESH_MESH_SURFACE_INTERSECTION_INTERNAL
42
50#include <geogram/mesh/mesh.h>
55
56namespace GEO {
57
58 class MeshSurfaceIntersection;
59
68 class MeshInTriangle : public CDTBase2d {
69 public:
70
71 /***************************************************************/
72
82 class Edge {
83 public:
84 Edge(
85 index_t v1_in = index_t(-1),
86 index_t v2_in = index_t(-1),
87 index_t f2 = index_t(-1),
88 TriangleRegion R2 = T2_RGN_T
89 ) : v1(v1_in),
90 v2(v2_in) {
91 sym.f2 = f2;
92 sym.R2 = R2;
93 }
94 index_t v1,v2; // The two extremities of the edge
95 struct { // Symbolic information: this edge = f1 /\ f2.R2
96 index_t f2;
98 } sym;
99 };
100
101 /***************************************************************/
102
110 class Vertex {
111 public:
112
113 enum Type {
114 UNINITIALIZED, MESH_VERTEX, PRIMARY_ISECT, SECONDARY_ISECT
115 };
116
125 ) {
126 geo_assert(f == M->f1_);
127 type = MESH_VERTEX;
128 mit = M;
129 init_sym(f, index_t(-1), TriangleRegion(lv), T2_RGN_T);
131 }
132
141 index_t f1, index_t f2,
143 ) {
144 geo_assert(f1 == M->f1_);
145 type = PRIMARY_ISECT;
146 mit = M;
147 init_sym(f1,f2,R1,R2);
149 }
150
157 MeshInTriangle* M, const vec3HE& point_exact_in
158 ) {
159 type = SECONDARY_ISECT;
160 mit = M;
161 init_sym(index_t(-1), index_t(-1), T1_RGN_T, T2_RGN_T);
162 init_geometry(point_exact_in);
163 }
164
169 type = UNINITIALIZED;
170 mit = nullptr;
171 init_sym(index_t(-1), index_t(-1), T1_RGN_T, T2_RGN_T);
172 mesh_vertex_index = index_t(-1);
173 }
174
179 const Mesh& mesh() const {
180 return mit->mesh();
181 }
182
188 void print(std::ostream& out=std::cerr) const;
189
195 std::string to_string() const {
196 std::ostringstream out;
197 print(out);
198 return out.str();
199 }
200
201 vec2 get_UV_approx() const {
202 double u = point_exact[mit->u_].estimate();
203 double v = point_exact[mit->v_].estimate();
204 double w = point_exact.w.estimate();
205 return vec2(u/w,v/w);
206 }
207
208 protected:
209
218 ) {
219 sym.f1 = f1;
220 sym.f2 = f2;
221 sym.R1 = R1;
222 sym.R2 = R2;
223 mesh_vertex_index = index_t(-1);
224 }
225
232
237 void init_geometry(const vec3HE& P);
238
239 public:
240 MeshInTriangle* mit;
241 vec3HE point_exact; // Exact homogeneous coords using expansions
242 double h_approx; // Lifting coordinate for incircle
243 Type type; // MESH_VERTEX, PRIMARY_ISECT or SECONDARY_ISECT
244 index_t mesh_vertex_index; // Global mesh vertex index once created
245 struct { // Symbolic information - tri-tri isect
246 index_t f1,f2; // global facet indices in mesh
247 TriangleRegion R1,R2; // triangle regions
248 } sym;
249 };
250
251 /***************************************************************/
252
254
259 const Mesh& mesh() const {
260 return mesh_;
261 }
262
268 return exact_mesh_.target_mesh();
269 }
270
275 void set_approx_incircle(bool x) {
276 approx_incircle_ = x;
277 }
278
284 void save_constraints(const std::string& filename) {
285 Mesh M;
287 mesh_save(M,filename);
288 }
289
290 void begin_facet(index_t f);
291
292 index_t add_vertex(index_t f2, TriangleRegion R1, TriangleRegion R2);
293
294 void add_edge(
295 index_t f2,
298 );
299
300 void end_facet() {
301 commit();
302 clear();
303 }
304
305 protected:
306
310 void commit();
311
315 void get_constraints(Mesh& M, bool with_edges=true) const;
316
317 vec3 mesh_vertex(index_t v) const {
318 return vec3(mesh().vertices.point_ptr(v));
319 }
320
321 vec3 mesh_facet_vertex(index_t f, index_t lv) const {
322 index_t v = mesh().facets.vertex(f,lv);
323 return mesh_vertex(v);
324 }
325
326 vec2 mesh_vertex_UV(index_t v) const {
327 const double* p = mesh().vertices.point_ptr(v);
328 return vec2(p[u_], p[v_]);
329 }
330
331 vec2 mesh_facet_vertex_UV(index_t f, index_t lv) const {
332 index_t v = mesh().facets.vertex(f,lv);
333 return mesh_vertex_UV(v);
334 }
335
336 void clear() override {
337 vertex_.resize(0);
338 edges_.resize(0);
339 f1_ = index_t(-1);
341 }
342
343 void log_err() const {
344 std::cerr << "Houston, we got a problem (while remeshing facet "
345 << f1_ << "):" << std::endl;
346 }
347
348 protected:
349
350 /********************** CDTBase2d overrides ***********************/
351
359 Sign orient2d(index_t v1,index_t v2,index_t v3) const override;
360
372 index_t v1,index_t v2,index_t v3,index_t v4
373 ) const override;
374
395 index_t e1, index_t i, index_t j,
396 index_t e2, index_t k, index_t l
397 ) override;
398
405 index_t e1, index_t e2, vec3HE& I
406 ) const;
407
415 index_t e1, index_t e2, vec3HE& I
416 ) const;
417
418 public:
419 void save(const std::string& filename) const override;
420
421 private:
422 MeshSurfaceIntersection& exact_mesh_;
423 const Mesh& mesh_;
424 index_t f1_;
425 index_t latest_f2_;
426 index_t latest_f2_count_;
427 coord_index_t f1_normal_axis_;
428 coord_index_t u_; // = (f1_normal_axis_ + 1)%3
429 coord_index_t v_; // = (f1_normal_axis_ + 2)%3
430 vector<Vertex> vertex_;
431 vector<Edge> edges_;
432 bool has_planar_isect_;
433 bool approx_incircle_;
434 };
435
436 /*************************************************************************/
437
444 struct IsectInfo {
445 public:
446
451 void flip() {
452 std::swap(f1,f2);
453 A_rgn_f1 = swap_T1_T2(A_rgn_f1);
454 A_rgn_f2 = swap_T1_T2(A_rgn_f2);
455 std::swap(A_rgn_f1, A_rgn_f2);
456 B_rgn_f1 = swap_T1_T2(B_rgn_f1);
457 B_rgn_f2 = swap_T1_T2(B_rgn_f2);
458 std::swap(B_rgn_f1, B_rgn_f2);
459 }
460
466 bool is_point() const {
467 return
468 A_rgn_f1 == B_rgn_f1 &&
469 A_rgn_f2 == B_rgn_f2 ;
470 }
471
472 index_t f1;
473 index_t f2;
474 TriangleRegion A_rgn_f1;
475 TriangleRegion A_rgn_f2;
476 TriangleRegion B_rgn_f1;
477 TriangleRegion B_rgn_f2;
478 };
479
480}
481
482#endif
Simple constained Delaunay triangulation in 2D.
#define geo_assert(x)
Verifies that a condition is met.
Definition assert.h:149
Common include file, providing basic definitions. Should be included before anything else by all head...
Base class for constrained Delaunay triangulation.
Definition CDT_2d.h:73
virtual void clear()
Removes everything from this triangulation.
index_t vertex(index_t f, index_t lv) const
Gets a vertex by facet and local vertex index.
Definition mesh.h:1043
void init_geometry(const vec3HE &P)
Optimizes exact numbers in generated points and computes approximate coordinates.
Vertex(MeshInTriangle *M, index_t f1, index_t f2, TriangleRegion R1, TriangleRegion R2)
Constructor for intersections with other facets.
void init_sym(index_t f1, index_t f2, TriangleRegion R1, TriangleRegion R2)
Initializes the symbolic information of this Vertex.
std::string to_string() const
Gets a string representation of this Vertex.
Vertex(MeshInTriangle *M, const vec3HE &point_exact_in)
Constructor for intersections between constraints.
vec3HE compute_geometry()
Gets the geometry of this vertex.
const Mesh & mesh() const
Gets the mesh.
void print(std::ostream &out=std::cerr) const
Prints this vertex.
Vertex(MeshInTriangle *M, index_t f, index_t lv)
Constructor for macro-triangle vertices.
Meshes a single triangle with the constraints that come from the intersections with the other triangl...
index_t create_intersection(index_t e1, index_t i, index_t j, index_t e2, index_t k, index_t l) override
Given two segments that have an intersection, create the intersection.
Sign incircle(index_t v1, index_t v2, index_t v3, index_t v4) const override
Tests the relative position of a point with respect to the circumscribed circle of a triangle.
Sign orient2d(index_t v1, index_t v2, index_t v3) const override
Tests the orientation of three vertices.
void clear() override
Removes everything from this triangulation.
void get_edge_edge_intersection(index_t e1, index_t e2, vec3HE &I) const
Computes the intersection between two edges.
const Mesh & mesh() const
Gets the readonly initial mesh.
void save_constraints(const std::string &filename)
For debugging, save constraints to a file.
void set_approx_incircle(bool x)
If Delaunay is set, use approximated incircle predicate (default: use exact incircle)
void get_edge_edge_intersection_2D(index_t e1, index_t e2, vec3HE &I) const
Auxilliary function used by get_edge_edge_intersection() for the special case when the two edges are ...
Mesh & target_mesh()
Gets the target mesh.
void save(const std::string &filename) const override
Saves this CDT to a geogram mesh file.
void commit()
Creates new vertices and new triangles in target mesh.
void get_constraints(Mesh &M, bool with_edges=true) const
For debugging, copies the constraints to a mesh.
Computes surface intersections.
Mesh & target_mesh()
Gets the target mesh.
const double * point_ptr(index_t v) const
Gets a point.
Definition mesh.h:444
Represents a mesh.
Definition mesh.h:2648
double estimate() const
Computes an approximation of the stored value in this expansion.
Vector with aligned memory allocation.
Definition memory.h:623
Exact predicates and constructs.
The class that represents a mesh.
Functions to load and save meshes.
Functions for computing intersections between surfacic meshes and boolean operations.
Global Vorpaline namespace.
Definition algorithm.h:64
TriangleRegion
Encodes the location of a point within a triangle.
TriangleRegion swap_T1_T2(TriangleRegion R)
Replaces T1 with T2 or T2 with T1 in a region code.
vecng< 3, Numeric::float64 > vec3
Represents points and vectors in 3d.
Definition geometry.h:65
bool mesh_save(const Mesh &M, const std::string &filename, const MeshIOFlags &ioflags=MeshIOFlags())
Saves a mesh to a file.
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
vecng< 2, Numeric::float64 > vec2
Represents points and vectors in 2d.
Definition geometry.h:59
geo_coord_index_t coord_index_t
The type for storing coordinate indices, and iterating on the coordinates of a point.
Definition numeric.h:321
Stores information about a triangle-triangle intersection.
bool is_point() const
Tests whether intersection is just a point.
3D vector in homogeneous coordinates with coordinates as arithmetic expansions.
Symbolic computation of triangle-triangle intersection.