Geogram Version 1.8.5
A programming library of geometric algorithms
Loading...
Searching...
No Matches
exact_geometry.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2000-2023 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_NUMERICS_EXACT_GEOMETRY
41#define GEOGRAM_NUMERICS_EXACT_GEOMETRY
42
46
58namespace GEO {
59
65
71
78 struct vec2HE {
79
84 x(expansion_nt::UNINITIALIZED),
85 y(expansion_nt::UNINITIALIZED),
86 w(expansion_nt::UNINITIALIZED)
87 {
88 }
89
90 vec2HE(
91 const expansion_nt& x_in,
92 const expansion_nt& y_in,
93 const expansion_nt& w_in
94 ) : x(x_in), y(y_in), w(w_in) {
95 }
96
97 vec2HE(
98 expansion_nt&& x_in,
99 expansion_nt&& y_in,
100 expansion_nt&& w_in
101 ) : x(x_in), y(y_in), w(w_in) {
102 }
103
104 vec2HE(const vec2HE& rhs) :
105 x(rhs.x), y(rhs.y), w(rhs.w) {
106 }
107
108 vec2HE(vec2HE&& rhs) :
109 x(rhs.x), y(rhs.y), w(rhs.w) {
110 }
111
112 explicit vec2HE(const vec2& rhs) :
113 x(rhs.x), y(rhs.y), w(1.0) {
114 }
115
116 vec2HE& operator=(const vec2HE& rhs) {
117 if(&rhs != this) {
118 x=rhs.x;
119 y=rhs.y;
120 w=rhs.w;
121 }
122 return *this;
123 }
124
125 vec2HE& operator=(vec2HE&& rhs) {
126 if(&rhs != this) {
127 x=rhs.x;
128 y=rhs.y;
129 w=rhs.w;
130 }
131 return *this;
132 }
133
134 expansion_nt* data() {
135 return &x;
136 }
137
138 const expansion_nt* data() const {
139 return &x;
140 }
141
142 expansion_nt& operator[](coord_index_t i) {
143 geo_debug_assert(i < 2);
144 return data()[i];
145 }
146
147 const expansion_nt& operator[](coord_index_t i) const {
148 geo_debug_assert(i < 2);
149 return data()[i];
150 }
151
156 void optimize() {
157 x.optimize();
158 y.optimize();
159 w.optimize();
160 }
161
162 expansion_nt x;
163 expansion_nt y;
164 expansion_nt w;
165 };
166
173 struct vec3HE {
174
179 x(expansion_nt::UNINITIALIZED),
180 y(expansion_nt::UNINITIALIZED),
181 z(expansion_nt::UNINITIALIZED),
182 w(expansion_nt::UNINITIALIZED)
183 {
184 }
185
186 vec3HE(
187 const expansion_nt& x_in,
188 const expansion_nt& y_in,
189 const expansion_nt& z_in,
190 const expansion_nt& w_in
191 ) : x(x_in), y(y_in), z(z_in), w(w_in) {
192 }
193
194 vec3HE(
195 expansion_nt&& x_in,
196 expansion_nt&& y_in,
197 expansion_nt&& z_in,
198 expansion_nt&& w_in
199 ) : x(x_in), y(y_in), z(z_in), w(w_in) {
200 }
201
202 vec3HE(
203 double x_in,
204 double y_in,
205 double z_in,
206 double w_in
207 ) : x(x_in), y(y_in), z(z_in), w(w_in) {
208 }
209
210 vec3HE(const vec3HE& rhs) :
211 x(rhs.x), y(rhs.y), z(rhs.z), w(rhs.w) {
212 }
213
214 vec3HE(vec3HE&& rhs) :
215 x(rhs.x), y(rhs.y), z(rhs.z), w(rhs.w) {
216 }
217
218 explicit vec3HE(const vec3& rhs) :
219 x(rhs.x), y(rhs.y), z(rhs.z), w(1.0) {
220 }
221
222 vec3HE& operator=(const vec3HE& rhs) {
223 if(&rhs != this) {
224 x=rhs.x;
225 y=rhs.y;
226 z=rhs.z;
227 w=rhs.w;
228 }
229 return *this;
230 }
231
232 vec3HE& operator=(vec3HE&& rhs) {
233 if(&rhs != this) {
234 x=rhs.x;
235 y=rhs.y;
236 z=rhs.z;
237 w=rhs.w;
238 }
239 return *this;
240 }
241
242 expansion_nt* data() {
243 return &x;
244 }
245
246 const expansion_nt* data() const {
247 return &x;
248 }
249
250 expansion_nt& operator[](coord_index_t i) {
251 geo_debug_assert(i < 3);
252 return data()[i];
253 }
254
255 const expansion_nt& operator[](coord_index_t i) const {
256 geo_debug_assert(i < 3);
257 return data()[i];
258 }
259
264 void optimize() {
265 x.optimize();
266 y.optimize();
267 z.optimize();
268 w.optimize();
269 }
270
271 expansion_nt x;
272 expansion_nt y;
273 expansion_nt z;
274 expansion_nt w;
275 };
276
277 vec2HE GEOGRAM_API operator-(const vec2HE& p1, const vec2HE& p2);
278
279 vec3HE GEOGRAM_API operator-(const vec3HE& p1, const vec3HE& p2);
280
285 class GEOGRAM_API vec2HELexicoCompare {
286 public:
293 bool operator()(const vec2HE& v1, const vec2HE& v2) const;
294 };
295
296
301 class GEOGRAM_API vec3HELexicoCompare {
302 public:
309 bool operator()(const vec3HE& v1, const vec3HE& v2) const;
310 };
311
316 public:
318 u = coord_index_t((coord+1)%3);
319 v = coord_index_t((coord+2)%3);
320 }
327 bool operator()(const vec3HE& v1, const vec3HE& v2) const;
330 };
331
332 vec3HE GEOGRAM_API mix(
333 const rational_nt& t, const vec3& p1, const vec3& p2
334 );
335
336 vec2HE GEOGRAM_API mix(
337 const rational_nt& t, const vec2& p1, const vec2& p2
338 );
339
340 vec2HE GEOGRAM_API mix(
341 const rational_nt& t, const vec2HE& p1, const vec2HE& p2
342 );
343
344 vec3HE GEOGRAM_API mix(
345 const rational_nt& t, const vec3HE& p1, const vec3HE& p2
346 );
347
351 template<> expansion_nt GEOGRAM_API det(const vec2E& v1, const vec2E& v2);
352
356 template<> expansion_nt GEOGRAM_API dot(const vec2E& v1, const vec2E& v2);
357
361 template<> expansion_nt GEOGRAM_API dot(const vec3E& v1, const vec3E& v2);
362
363 namespace PCK {
364
365 template <class T> inline bool same_point(
366 const vecng<3,T>& v1, const vecng<3,T>& v2
367 ) {
368 // operator== is well optimized for expansion_nt and rational_nt
369 return (v1.x == v2.x) && (v1.y == v2.y) && (v2.z == v1.z);
370 }
371
372 template <class T> inline bool same_point(
373 const vecng<2,T>& v1, const vecng<2,T>& v2
374 ) {
375 // operator== is well optimized for expansion_nt and rational_nt
376 return (v1.x == v2.x) && (v1.y == v2.y);
377 }
378
379 bool GEOGRAM_API same_point(const vec2HE& v1, const vec2HE& v2);
380 bool GEOGRAM_API same_point(const vec3HE& v1, const vec3HE& v2);
381
382 Sign GEOGRAM_API orient_2d(
383 const vec2HE& p0, const vec2HE& p1, const vec2HE& p2
384 );
385
386 Sign GEOGRAM_API orient_2d_projected(
387 const vec3HE& p0, const vec3HE& p1, const vec3HE& p2,
388 coord_index_t axis
389 );
390
391 // TODO: check, it seems that orient_3d(vec3,...) and
392 // orient_3d(vec3HE,...) have opposite orientations !
393 Sign GEOGRAM_API orient_3d(
394 const vec3HE& p0, const vec3HE& p1,
395 const vec3HE& p2, const vec3HE& p3
396 );
397
398 Sign GEOGRAM_API dot_2d(
399 const vec2HE& p0, const vec2HE& p1, const vec2HE& p2
400 );
401
421 const vec2HE& p0, const vec2HE& p1,
422 const vec2HE& p2, const vec2HE& p3,
423 double h0, double h1, double h2, double h3
424 );
425
426 Sign GEOGRAM_API orient_2dlifted_SOS_projected(
427 const vec3HE& p0, const vec3HE& p1,
428 const vec3HE& p2, const vec3HE& p3,
429 double h0, double h1, double h2, double h3,
430 coord_index_t axis
431 );
432
433 }
434
442 template <class VEC3>
443 inline VEC3 make_vec3(const vec3& p) {
444 typedef typename VEC3::value_type value_type;
445 return VEC3(value_type(p.x),value_type(p.y),value_type(p.z));
446 }
447
455 template <class VEC3 = vec3>
456 inline VEC3 make_vec3(const vec3& p1, const vec3& p2) {
457 typedef typename VEC3::value_type value_type;
458 return VEC3(
459 value_type(p2.x) - value_type(p1.x),
460 value_type(p2.y) - value_type(p1.y),
461 value_type(p2.z) - value_type(p1.z)
462 );
463 }
464
468 template <>
469 inline vec3E make_vec3<vec3E>(const vec3& p1, const vec3& p2) {
470 return vec3E(
471 expansion_nt(expansion_nt::DIFF, p2.x, p1.x),
472 expansion_nt(expansion_nt::DIFF, p2.y, p1.y),
473 expansion_nt(expansion_nt::DIFF, p2.z, p1.z)
474 );
475 }
476
484 template <class VEC2>
485 inline VEC2 make_vec2(
486 const vec2& p1, const vec2& p2
487 ) {
488 typedef typename VEC2::value_type value_type;
489 return VEC2(
490 value_type(p2.x) - value_type(p1.x),
491 value_type(p2.y) - value_type(p1.y)
492 );
493 }
494
498 template <>
499 inline vec2E make_vec2<vec2E>(const vec2& p1, const vec2& p2) {
500 return vec2E(
501 expansion_nt(expansion_nt::DIFF, p2.x, p1.x),
502 expansion_nt(expansion_nt::DIFF, p2.y, p1.y)
503 );
504 }
505
514 template <class VEC3>
515 inline VEC3 triangle_normal(
516 const vec3& p1, const vec3& p2, const vec3& p3
517 ) {
518 return cross(
519 make_vec3<VEC3>(p1,p2),
520 make_vec3<VEC3>(p1,p3)
521 );
522 }
523
527 template <> GEOGRAM_API vec3E triangle_normal<vec3E>(
528 const vec3& p1, const vec3& p2, const vec3& p3
529 );
530
543 vec3HE& result,
544 const vec3& p1, const vec3& p2, const vec3& p3,
545 const vec3& q1, const vec3& q2, const vec3& q3,
546 const vec3& r1, const vec3& r2, const vec3& r3
547 );
548
558 const vec3& p1, const vec3& p2, const vec3& p3,
559 const vec3& q1, const vec3& q2
560 );
561
562 /************************************************************************/
563
573 const vec3& p1, const vec3& p2, const vec3& p3
574 );
575}
576
577#endif
#define geo_debug_assert(x)
Verifies that a condition is met.
Definition assert.h:195
Common include file, providing basic definitions. Should be included before anything else by all head...
Expansion_nt (expansion Number Type) is used to compute the sign of polynoms exactly.
void optimize()
Optimizes the internal representation without changing the represented value.
Rational_nt (rational Number Type) is used to compute the sign of rational fractions exactly.
Comparator class for vec3HE \detail Used to create maps indexed by vec3HE.
bool operator()(const vec2HE &v1, const vec2HE &v2) const
Compares two vec3HE.
Comparator class for vec3HE \detail Used to create maps indexed by vec3HE.
bool operator()(const vec3HE &v1, const vec3HE &v2) const
Compares two vec3HE.
Comparator class for projected vec3HE.
bool operator()(const vec3HE &v1, const vec3HE &v2) const
Compares two vec3HE.
Generic maths vector.
Definition vecg.h:69
High-level interface to multi-precision arithmetics.
Geometric functions in 2d and 3d.
Sign orient_2dlifted_SOS(const vec2HE &p0, const vec2HE &p1, const vec2HE &p2, const vec2HE &p3, double h0, double h1, double h2, double h3)
Computes the 3d orientation test with lifted points.
Global Vorpaline namespace.
Definition algorithm.h:64
expansion_nt det(const vec2E &v1, const vec2E &v2)
Specialization optimized using low-level API.
Quaternion operator-(const Quaternion &a, const Quaternion &b)
Computes the difference between two Quaternion.
Definition quaternion.h:252
T dot(const vecng< 3, T > &v1, const vecng< 3, T > &v2)
Computes the dot product of 2 vectors. vecng
Definition vecg.h:869
vecng< 2, expansion_nt > vec2E
vec2 with coordinates as expansions
vecng< 3, Numeric::float64 > vec3
Represents points and vectors in 3d.
Definition geometry.h:65
vec3E make_vec3< vec3E >(const vec3 &p1, const vec3 &p2)
Specialization for vec3E.
VEC3 make_vec3(const vec3 &p)
Converts a 3d vector with double coordinates into a 3d vector with coordinates of arbitrary type.
vec2E make_vec2< vec2E >(const vec2 &p1, const vec2 &p2)
Specialization for vec2E.
Sign
Integer constants that represent the sign of a value.
Definition numeric.h:224
bool get_three_planes_intersection(vec3HE &result, const vec3 &p1, const vec3 &p2, const vec3 &p3, const vec3 &q1, const vec3 &q2, const vec3 &q3, const vec3 &r1, const vec3 &r2, const vec3 &r3)
Computes the exact intersection between the support planes of three triangles.
VEC2 make_vec2(const vec2 &p1, const vec2 &p2)
Creates a vector with coordinates of arbitrary type from two points with double coordinates.
vec3E triangle_normal< vec3E >(const vec3 &p1, const vec3 &p2, const vec3 &p3)
Specialization for vec3E.
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
vec3HE plane_line_intersection(const vec3 &p1, const vec3 &p2, const vec3 &p3, const vec3 &q1, const vec3 &q2)
Computes the exact intersection between the support plane of a triangle and the support line of a seg...
VEC3 triangle_normal(const vec3 &p1, const vec3 &p2, const vec3 &p3)
Computes the normal to a triangle from its three vertices.
vecng< 3, expansion_nt > vec3E
vec3 with coordinates as expansions
coord_index_t triangle_normal_axis_exact(const vec3 &p1, const vec3 &p2, const vec3 &p3)
Finds an axis along which a triangle can be projected without degeneracy.
2D vector in homogeneous coordinates with coordinates as arithmetic expansions
vec2HE()
Creates an uninitialized vec2HE.
void optimize()
Optimizes the internal storage of the expansions used to store the coordinates.
3D vector in homogeneous coordinates with coordinates as arithmetic expansions.
void optimize()
Optimizes the internal storage of the expansions used to store the coordinates.
vec3HE()
Creates an uninitialized vec3HE.