TTK
Searching...
No Matches
Geometry.cpp
Go to the documentation of this file.
1#include <Geometry.h>
2
3#include <algorithm>
4
5using namespace std;
6using namespace ttk;
7
8static const double PREC_DBL{Geometry::pow(10.0, -DBL_DIG)};
9static const float PREC_FLT{powf(10.0F, -FLT_DIG)};
10static const float PREC_FLT_1{powf(10.0F, -FLT_DIG + 1)};
11
12template <typename T>
13T Geometry::angle(const T *vA0, const T *vA1, const T *vB0, const T *vB1) {
14 return M_PI
15 - acos(dotProduct(vA0, vA1, vB0, vB1)
16 / (magnitude(vA0, vA1) * magnitude(vB0, vB1)));
17}
18
19template <typename T>
20T Geometry::angle2D(const T *vA0, const T *vA1, const T *vB0, const T *vB1) {
21 double vecA[2] = {vA1[0] - vA0[0], vA1[1] - vA0[1]};
22 double vecB[2] = {vB1[0] - vB0[0], vB1[1] - vB0[1]};
23 return atan2(vecB[1], vecB[0]) - atan2(vecA[1], vecA[0]);
24}
25
26template <typename T>
28 const T *vA1,
29 const T *vB0,
30 const T *vB1,
31 std::array<T, 3> *coefficients,
32 const T *tolerance) {
33
34 int aNullComponents = 0, bNullComponents = 0;
35 std::array<T, 3> a{}, b{};
36 for(int i = 0; i < 3; i++) {
37 a[i] = vA1[i] - vA0[i];
38 if(fabs(a[i]) < PREC_FLT) {
39 aNullComponents++;
40 }
41 b[i] = vB1[i] - vB0[i];
42 if(fabs(b[i]) < PREC_FLT) {
43 bNullComponents++;
44 }
45 }
46
47 if((aNullComponents == 3) || (bNullComponents == 3)) {
48 return true;
49 }
50
51 // check for axis aligned vectors
52 if((aNullComponents > 1) || (bNullComponents > 1)) {
53 if(aNullComponents == bNullComponents) {
54 // only one non-null component for both vectors
55 return true;
56 }
57 }
58
59 bool useDenominatorA = false;
60 T sumA = 0, sumB = 0;
61 for(int i = 0; i < 3; i++) {
62 sumA += fabs(a[i]);
63 sumB += fabs(b[i]);
64 }
65 if(sumA > sumB) {
66 useDenominatorA = true;
67 }
68
69 std::array<T, 3> k{};
70
71 T maxDenominator = 0;
72 int isNan = -1, maximizer = 0;
73 for(int i = 0; i < 3; i++) {
74 if(useDenominatorA) {
75 if(fabs(a[i]) > PREC_FLT) {
76 k[i] = b[i] / a[i];
77 } else {
78 isNan = i;
79 }
80 } else {
81 if(fabs(b[i]) > PREC_FLT) {
82 k[i] = a[i] / b[i];
83 } else {
84 isNan = i;
85 }
86 }
87
88 if(!i) {
89 maxDenominator = fabs(k[i]);
90 maximizer = i;
91 } else {
92 if(fabs(k[i]) > maxDenominator) {
93 maxDenominator = fabs(k[i]);
94 maximizer = i;
95 }
96 }
97 }
98
99 T colinearityThreshold;
100
101 colinearityThreshold = PREC_FLT;
102 if(tolerance) {
103 colinearityThreshold = *tolerance;
104 }
105
106 if(coefficients) {
107 (*coefficients) = k;
108 }
109
110 if(isNan == -1) {
111
112 if((fabs(1 - fabs(k[(maximizer + 1) % 3] / k[maximizer]))
113 < colinearityThreshold)
114 && (fabs(1 - fabs(k[(maximizer + 2) % 3] / k[maximizer]))
115 < colinearityThreshold)) {
116 return true;
117 }
118 } else {
119 if(fabs(1 - fabs(k[(isNan + 1) % 3] / k[(isNan + 2) % 3]))
120 < colinearityThreshold) {
121 return true;
122 }
123 }
124 k[0] = k[1] = k[2] = 0;
125
126 return false;
127}
128
129template <typename T>
131 const T *pptB,
132 const T *pptC,
133 const T tolerance) {
134 double ptA[2] = {pptA[0], pptA[1]}, ptB[2] = {pptB[0], pptB[1]},
135 ptC[2] = {pptC[0], pptC[1]};
136 return fabs(ptA[0] * (ptB[1] - ptC[1]) + ptB[0] * (ptC[1] - ptA[1])
137 + ptC[0] * (ptA[1] - ptB[1]))
138 <= tolerance;
139}
140
141template <typename T>
143 const T *p1,
144 const T *p,
145 std::array<T, 2> &baryCentrics,
146 const int &dimension) {
147
148 if(dimension > 3) {
149 return -1;
150 }
151
152 int bestI = 0;
153 T maxDenominator = 0;
154
155 for(int i = 0; i < dimension; i++) {
156
157 T denominator = fabs(p0[i] - p1[i]);
158 if(!i) {
159 maxDenominator = denominator;
160 bestI = i;
161 } else {
162 if(denominator > maxDenominator) {
163 maxDenominator = denominator;
164 bestI = i;
165 }
166 }
167 }
168
169 baryCentrics[0] = p0[bestI] - p1[bestI];
170 baryCentrics[0] = (p[bestI] - p1[bestI]) / baryCentrics[0];
171
172 baryCentrics[1] = 1 - baryCentrics[0];
173
174 // check if the point lies in the edge
175 std::array<T, 3> test{};
176 for(int i = 0; i < dimension; i++) {
177 test[i] = baryCentrics[0] * p0[i] + baryCentrics[1] * p1[i];
178 }
179
180 if((!((fabs(test[0] - p[0]) < PREC_FLT_1)
181 && (fabs(test[1] - p[1]) < PREC_FLT_1)))) {
182 for(int i = 0; i < 2; i++) {
183 baryCentrics[i] = -baryCentrics[i];
184 }
185 }
186
187 return 0;
188}
189template <typename T>
191 const T *p1,
192 const T *p2,
193 const T *p,
194 std::array<T, 3> &baryCentrics) {
195
196 // find the pair of coordinates that maximize the sum of the denominators
197 // (more stable computations)
198 int bestI = 0, bestJ = 1;
199 T maxDenominator = 0;
200
201 for(int i = 0; i < 2; i++) {
202 for(int j = i + 1; j < 3; j++) {
203
204 baryCentrics[0]
205 = (p1[j] - p2[j]) * (p0[i] - p2[i]) + (p2[i] - p1[i]) * (p0[j] - p2[j]);
206 baryCentrics[1]
207 = (p1[j] - p2[j]) * (p0[i] - p2[i]) + (p2[i] - p1[i]) * (p0[j] - p2[j]);
208
209 T denominator = fabs(baryCentrics[0]);
210
211 if(fabs(baryCentrics[1]) < denominator) {
212 denominator = fabs(baryCentrics[1]);
213 }
214
215 if((i == 0) && (j == 1)) {
216 maxDenominator = denominator;
217 } else {
218 if(denominator > maxDenominator) {
219 maxDenominator = denominator;
220 bestI = i;
221 bestJ = j;
222 }
223 }
224 }
225 }
226
227 baryCentrics[0] = (p1[bestJ] - p2[bestJ]) * (p0[bestI] - p2[bestI])
228 + (p2[bestI] - p1[bestI]) * (p0[bestJ] - p2[bestJ]);
229 // (y1 - y2)*(x - x2) + (x2 - x1)*(y - y2)
230 baryCentrics[0] = ((p1[bestJ] - p2[bestJ]) * (p[bestI] - p2[bestI])
231 + (p2[bestI] - p1[bestI]) * (p[bestJ] - p2[bestJ]))
232 / baryCentrics[0];
233
234 // (y1 - y2)*(x0 - x2) + (x2 - x1)*(y0 - y2)
235 baryCentrics[1] = (p1[bestJ] - p2[bestJ]) * (p0[bestI] - p2[bestI])
236 + (p2[bestI] - p1[bestI]) * (p0[bestJ] - p2[bestJ]);
237 // (y2 - y0)*(x - x2) + (x0 - x2)*(y - y2)
238 baryCentrics[1] = ((p2[bestJ] - p0[bestJ]) * (p[bestI] - p2[bestI])
239 + (p0[bestI] - p2[bestI]) * (p[bestJ] - p2[bestJ]))
240 / baryCentrics[1];
241
242 baryCentrics[2] = 1 - baryCentrics[0] - baryCentrics[1];
243
244 return 0;
245}
246
247template <typename T>
249 const T &yA,
250 const T &xB,
251 const T &yB,
252 const T &xC,
253 const T &yC,
254 const T &xD,
255 const T &yD,
256 T &x,
257 T &y) {
258
259 T d = (xA - xB) * (yC - yD) - (yA - yB) * (xC - xD);
260
261 if(fabs(d) < PREC_DBL) {
262 return false;
263 }
264
265 x = ((xC - xD) * (xA * yB - yA * xB) - (xA - xB) * (xC * yD - yC * xD)) / d;
266
267 y = ((yC - yD) * (xA * yB - yA * xB) - (yA - yB) * (xC * yD - yC * xD)) / d;
268
269 if((x < std::min(xA, xB) - PREC_FLT) || (x > std::max(xA, xB) + PREC_FLT)) {
270 return false;
271 }
272
273 if((x < std::min(xC, xD) - PREC_FLT) || (x > std::max(xC, xD) + PREC_FLT)) {
274 return false;
275 }
276
277 return true;
278}
279
280template <typename T>
282 const T *p1,
283 const T *p2,
284 T &area) {
285
286 std::array<T, 3> cross{};
287
288 crossProduct(p0, p1, p1, p2, cross);
289
290 area = 0.5 * magnitude(cross.data());
291
292 return 0;
293}
294
295template <typename T>
297 const T s1,
298 const T s2,
299 T &area) {
300
301 double const s = (s0 + s1 + s2) / 2.0;
302 area = std::sqrt(s * (s - s0) * (s - s1) * (s - s2));
303
304 return 0;
305}
306
307template <typename T>
309 const T *p1,
310 const T *p2,
311 std::array<T, 3> &angles) {
312
313 angles[0] = angle(p0, p1, p1, p2);
314 angles[1] = angle(p1, p2, p2, p0);
315 angles[2] = angle(p2, p0, p0, p1);
316
317 return 0;
318}
319
320template <typename T>
322 const T s1,
323 const T s2,
324 T &angle) {
325
326 angle = std::acos((s0 * s0 + s1 * s1 - s2 * s2) / (2.0 * s0 * s1));
327
328 return 0;
329}
330
331template <typename T>
332int Geometry::crossProduct(const T *vA0,
333 const T *vA1,
334 const T *vB0,
335 const T *vB1,
336 std::array<T, 3> &crossProduct) {
337
338 std::array<T, 3> a{}, b{};
339
340 for(int i = 0; i < 3; i++) {
341 a[i] = vA1[i] - vA0[i];
342 b[i] = vB1[i] - vB0[i];
343 }
344
345 for(int i = 0; i < 3; i++) {
346 crossProduct[i]
347 = a[(i + 1) % 3] * b[(i + 2) % 3] - a[(i + 2) % 3] * b[(i + 1) % 3];
348 }
349
350 return 0;
351}
352
353template <typename T>
354int Geometry::crossProduct(const T *vA, const T *vB, T *crossProduct) {
355 crossProduct[0] = vA[1] * vB[2] - vA[2] * vB[1];
356 crossProduct[1] = vA[2] * vB[0] - vA[0] * vB[2];
357 crossProduct[2] = vA[0] * vB[1] - vA[1] * vB[0];
358 return 0;
359}
360
361template <typename T>
362T Geometry::distance(const T *p0, const T *p1, const int &dimension) {
363
364 T distance = 0;
365
366 for(int i = 0; i < dimension; i++) {
367 distance += (p0[i] - p1[i]) * (p0[i] - p1[i]);
368 }
369
370 return sqrt(distance);
371}
372
373template <typename T>
374T Geometry::distance(const std::vector<T> &p0, const std::vector<T> &p1) {
375 return distance(p0.data(), p1.data(), p0.size());
376}
377
378template <typename T>
379T Geometry::distanceFlatten(const std::vector<std::vector<T>> &p0,
380 const std::vector<std::vector<T>> &p1) {
381 std::vector<T> p0_flatten, p1_flatten;
382 flattenMultiDimensionalVector(p0, p0_flatten);
383 flattenMultiDimensionalVector(p1, p1_flatten);
384 return distance(p0_flatten, p1_flatten);
385}
386
387template <typename T>
388T Geometry::dotProduct(const T *vA0, const T *vA1, const T *vB0, const T *vB1) {
389
390 T dotProduct = 0;
391 for(int i = 0; i < 3; i++) {
392 dotProduct += (vA1[i] - vA0[i]) * (vB1[i] - vB0[i]);
393 }
394
395 return dotProduct;
396}
397
398template <typename T>
399T Geometry::dotProduct(const T *vA, const T *vB, const int &dimension) {
400 T dotProduct = 0;
401 for(int i = 0; i < dimension; ++i)
402 dotProduct += vA[i] * vB[i];
403 return dotProduct;
404}
405
406template <typename T>
407T Geometry::dotProduct(const std::vector<T> &vA, const std::vector<T> &vB) {
408 return dotProduct(vA.data(), vB.data(), vA.size());
409}
410
411template <typename T>
412T Geometry::dotProductFlatten(const std::vector<std::vector<T>> &vA,
413 const std::vector<std::vector<T>> &vB) {
414 std::vector<T> vA_flatten, vB_flatten;
415 flattenMultiDimensionalVector(vA, vA_flatten);
416 flattenMultiDimensionalVector(vB, vB_flatten);
417 return dotProduct(vA_flatten, vB_flatten);
418}
419
420template <typename T>
422 const T *p1,
423 const T *p2,
424 const T *p) {
425
426 std::array<T, 3> barycentrics{};
427
428 Geometry::computeBarycentricCoordinates(p0, p1, p2, p, barycentrics);
429
430 for(int i = 0; i < static_cast<int>(barycentrics.size()); i++) {
431 if(barycentrics[i] < -PREC_DBL) {
432 return false;
433 }
434 if(barycentrics[i] > 1 + PREC_DBL) {
435 return false;
436 }
437 }
438
439 return true;
440}
441
442template <typename T>
444 const T &x, const T &y, const T &xA, const T &yA, const T &xB, const T &yB) {
445
446 std::array<T, 2> pA{xA, yA}, pB{xB, yB}, p{x, y};
447 return Geometry::isPointOnSegment(p.data(), pA.data(), pB.data(), 2);
448}
449
450template <typename T>
452 const T *pA,
453 const T *pB,
454 const int &dimension) {
455
456 std::array<T, 2> baryCentrics{};
457
458 Geometry::computeBarycentricCoordinates(pA, pB, p, baryCentrics, dimension);
459
460 return (
461 ((baryCentrics[0] > -PREC_DBL) && (baryCentrics[0] < 1 + PREC_DBL))
462 && ((baryCentrics[1] > -PREC_DBL) && (baryCentrics[1] < 1 + PREC_DBL)));
463}
464
465template <typename T>
467 const T *p1,
468 const T *p2,
469 const T *tolerance) {
470
471 bool maxDecision = false;
472 T maxCoefficient = 0;
473 std::array<T, 3> coefficients{};
474
475 bool decision = areVectorsColinear(p0, p1, p1, p2, &coefficients, tolerance);
476 maxDecision = decision;
477 for(int i = 0; i < 3; i++) {
478 if(!i) {
479 maxCoefficient = fabs(coefficients[i]);
480 maxDecision = decision;
481 } else {
482 if(fabs(coefficients[i]) > maxCoefficient) {
483 maxCoefficient = fabs(coefficients[i]);
484 maxDecision = decision;
485 }
486 }
487 }
488
489 decision = areVectorsColinear(p0, p2, p2, p1, &coefficients, tolerance);
490 for(int i = 0; i < 3; i++) {
491 if(fabs(coefficients[i]) > maxCoefficient) {
492 maxCoefficient = fabs(coefficients[i]);
493 maxDecision = decision;
494 }
495 }
496
497 decision = areVectorsColinear(p1, p0, p0, p2, &coefficients, tolerance);
498 for(int i = 0; i < 3; i++) {
499 if(fabs(coefficients[i]) > maxCoefficient) {
500 maxCoefficient = fabs(coefficients[i]);
501 maxDecision = decision;
502 }
503 }
504
505 return maxDecision;
506}
507
508template <typename T>
509T Geometry::magnitude(const T *v, const int &dimension) {
510 return sqrt(dotProduct(v, v, dimension));
511}
512
513template <typename T>
514T Geometry::magnitude(const std::vector<T> &v) {
515 return magnitude(v.data(), v.size());
516}
517
518template <typename T>
519T Geometry::magnitudeFlatten(const std::vector<std::vector<T>> &v) {
520 std::vector<T> v_flatten;
521 flattenMultiDimensionalVector(v, v_flatten);
522 return magnitude(v_flatten);
523}
524
525template <typename T>
526T Geometry::magnitude(const T *o, const T *d) {
527
528 T mag = 0;
529
530 for(int i = 0; i < 3; i++) {
531 mag += (o[i] - d[i]) * (o[i] - d[i]);
532 }
533
534 return sqrt(mag);
535}
536
537template <typename T>
539 const T *b,
540 T *out,
541 const int &dimension) {
542 for(int i = 0; i < dimension; ++i)
543 out[i] = b[i] - a[i];
544 return 0;
545}
546
547template <typename T>
548int Geometry::subtractVectors(const std::vector<T> &a,
549 const std::vector<T> &b,
550 std::vector<T> &out) {
551 out.resize(a.size());
552 return subtractVectors(a.data(), b.data(), out.data(), a.size());
553}
554
555template <typename T>
556int Geometry::addVectors(const T *a, const T *b, T *out, const int &dimension) {
557 for(int i = 0; i < dimension; ++i)
558 out[i] = b[i] + a[i];
559 return 0;
560}
561
562template <typename T>
564 const std::vector<T> &b,
565 std::vector<T> &out) {
566 out.resize(a.size());
567 return addVectors(a.data(), b.data(), out.data(), a.size());
568}
569
570template <typename T>
572 const std::vector<std::vector<T>> &b,
573 std::vector<std::vector<T>> &out) {
574 out.resize(a.size());
575 for(unsigned int i = 0; i < a.size(); ++i)
577 return 0;
578}
579
580template <typename T>
582 const std::vector<std::vector<std::vector<T>>> &a,
583 const std::vector<std::vector<std::vector<T>>> &b,
584 std::vector<std::vector<T>> &out) {
585 std::vector<std::vector<T>> a_flatten, b_flatten;
589 return 0;
590}
591
592template <typename T>
594 const T factor,
595 T *out,
596 const int &dimension) {
597 for(int i = 0; i < dimension; ++i)
598 out[i] = a[i] * factor;
599 return 0;
600}
601
602template <typename T>
603int Geometry::scaleVector(const std::vector<T> &a,
604 const T factor,
605 std::vector<T> &out) {
606 out.resize(a.size());
607 return scaleVector(a.data(), factor, out.data(), a.size());
608}
609
610template <typename T>
612 const T *b,
613 T *out,
614 const int &dimension) {
615 T dotProdBB = dotProduct(b, b, dimension);
616 T dotProdAB;
617 if(dotProdBB > PREC_DBL) {
618 dotProdAB = dotProduct(a, b, dimension);
619 dotProdAB /= dotProdBB;
620 } else
621 dotProdAB = 0; // Gram-Schmidt convention
622 for(int i = 0; i < dimension; ++i)
623 out[i] = b[i] * dotProdAB;
624 return 0;
625}
626
627template <typename T>
628int Geometry::vectorProjection(const std::vector<T> &a,
629 const std::vector<T> &b,
630 std::vector<T> &out) {
631 out.resize(a.size(), 0.0);
632 return vectorProjection(a.data(), b.data(), out.data(), a.size());
633}
634
635template <typename T>
637 const std::vector<T> &b,
638 std::vector<T> &a_out,
639 std::vector<T> &b_out) {
640 std::vector<T> sumV;
642 vectorProjection(a, sumV, a_out);
643 vectorProjection(b, sumV, b_out);
644}
645
646template <typename T>
647void Geometry::gramSchmidt(const std::vector<std::vector<T>> &a,
648 std::vector<std::vector<T>> &out) {
649 out.resize(a.size());
650 out[0] = a[0];
651 for(unsigned int i = 1; i < a.size(); ++i) {
652 std::vector<T> projecSum;
653 vectorProjection(a[i], out[0], projecSum);
654 for(unsigned int j = 1; j < i; ++j) {
655 std::vector<T> projecTemp, projecSumTemp;
656 vectorProjection(a[i], out[j], projecTemp);
658 projecSum = projecSumTemp;
659 }
660 subtractVectors(projecSum, a[i], out[i]);
661 }
662}
663
664template <typename T>
665bool Geometry::isVectorUniform(const std::vector<T> &a) {
666 for(unsigned int i = 0; i < a.size() - 1; ++i)
667 if(not(std::abs(a[i] - a[i + 1]) < PREC_DBL))
668 return false;
669 return true;
670}
671
672template <typename T>
673bool Geometry::isVectorNull(const std::vector<T> &a) {
674 for(unsigned int i = 0; i < a.size(); ++i)
675 if(not(std::abs(a[i]) < PREC_DBL))
676 return false;
677 return true;
678}
679
680template <typename T>
681bool Geometry::isVectorNullFlatten(const std::vector<std::vector<T>> &a) {
682 std::vector<T> a_flatten;
683 flattenMultiDimensionalVector(a, a_flatten);
684 return isVectorNull(a_flatten);
685}
686
687template <typename T>
689 const std::vector<std::vector<T>> &a, std::vector<T> &out) {
690 out.resize(a.size() * a[0].size());
691 for(unsigned int i = 0; i < a.size(); ++i)
692 for(unsigned int j = 0; j < a[0].size(); ++j)
693 out[i * a[0].size() + j] = a[i][j];
694 return 0;
695}
696
697template <typename T>
699 const std::vector<std::vector<std::vector<T>>> &a,
700 std::vector<std::vector<T>> &out) {
701 out.resize(a.size());
702 for(unsigned int i = 0; i < a.size(); ++i)
703 flattenMultiDimensionalVector(a[i], out[i]);
704 return 0;
705}
706
707template <typename T>
709 std::vector<std::vector<T>> &out,
710 const int &no_columns) {
711 if(a.size() % no_columns != 0)
712 return -1;
713 out.resize(a.size() / no_columns);
714 for(unsigned int i = 0; i < out.size(); ++i) {
715 out[i].resize(no_columns);
716 for(unsigned int j = 0; j < out[i].size(); ++j)
717 out[i][j] = a[i * no_columns + j];
718 }
719 return 0;
720}
721
722template <typename T>
723void Geometry::matrixMultiplication(const std::vector<std::vector<T>> &a,
724 const std::vector<std::vector<T>> &b,
725 std::vector<std::vector<T>> &out) {
726 out.resize(a.size(), std::vector<T>(b[0].size(), 0.0));
727 for(unsigned int i = 0; i < out.size(); ++i)
728 for(unsigned int j = 0; j < out[i].size(); ++j)
729 for(unsigned int k = 0; k < a[i].size(); ++k)
730 out[i][j] += a[i][k] * b[k][j];
731}
732
733template <typename T>
734void Geometry::subtractMatrices(const std::vector<std::vector<T>> &a,
735 const std::vector<std::vector<T>> &b,
736 std::vector<std::vector<T>> &out) {
737 out.resize(a.size(), std::vector<T>(a[0].size()));
738 for(unsigned int i = 0; i < out.size(); ++i)
739 for(unsigned int j = 0; j < out[0].size(); ++j)
740 out[i][j] = b[i][j] - a[i][j];
741}
742
743template <typename T>
745 const std::vector<std::vector<T>> &b,
746 std::vector<std::vector<T>> &out) {
747 out.resize(a.size(), std::vector<T>(a[0].size()));
748 for(unsigned int i = 0; i < a.size(); ++i)
749 for(unsigned int j = 0; j < a[0].size(); ++j)
750 out[i][j] = a[i][j] + b[i][j];
751}
752
753template <typename T>
754void Geometry::scaleMatrix(const std::vector<std::vector<T>> &a,
755 const T factor,
756 std::vector<std::vector<T>> &out) {
757 out.resize(a.size(), std::vector<T>(a[0].size()));
758 for(unsigned int i = 0; i < out.size(); ++i)
759 for(unsigned int j = 0; j < out[i].size(); ++j)
760 out[i][j] = a[i][j] * factor;
761}
762
763template <typename T>
764void Geometry::transposeMatrix(const std::vector<std::vector<T>> &a,
765 std::vector<std::vector<T>> &out) {
766 out.resize(a[0].size(), std::vector<T>(a.size()));
767 for(unsigned int i = 0; i < a.size(); ++i)
768 for(unsigned int j = 0; j < a[0].size(); ++j)
769 out[j][i] = a[i][j];
770}
771
772#define GEOMETRY_SPECIALIZE(TYPE) \
773 template TYPE Geometry::angle<TYPE>( \
774 TYPE const *, TYPE const *, TYPE const *, TYPE const *); \
775 template TYPE Geometry::angle2D<TYPE>( \
776 TYPE const *, TYPE const *, TYPE const *, TYPE const *); \
777 template bool Geometry::areVectorsColinear<TYPE>( \
778 TYPE const *, TYPE const *, TYPE const *, TYPE const *, \
779 std::array<TYPE, 3> *, TYPE const *); \
780 template bool Geometry::isTriangleColinear2D<TYPE>( \
781 TYPE const *, TYPE const *, TYPE const *, TYPE const); \
782 template int Geometry::computeBarycentricCoordinates<TYPE>( \
783 TYPE const *, TYPE const *, TYPE const *, std::array<TYPE, 2> &, \
784 int const &); \
785 template int Geometry::computeBarycentricCoordinates<TYPE>( \
786 TYPE const *, TYPE const *, TYPE const *, TYPE const *, \
787 std::array<TYPE, 3> &); \
788 template bool Geometry::computeSegmentIntersection<TYPE>( \
789 TYPE const &, TYPE const &, TYPE const &, TYPE const &, TYPE const &, \
790 TYPE const &, TYPE const &, TYPE const &, TYPE &, TYPE &); \
791 template int Geometry::computeTriangleAngles<TYPE>( \
792 TYPE const *, TYPE const *, TYPE const *, std::array<TYPE, 3> &); \
793 template int Geometry::computeTriangleAngleFromSides<TYPE>( \
794 TYPE const, TYPE const, TYPE const, TYPE &); \
795 template int Geometry::computeTriangleArea<TYPE>( \
796 TYPE const *, TYPE const *, TYPE const *, TYPE &); \
797 template int Geometry::computeTriangleAreaFromSides<TYPE>( \
798 TYPE const, TYPE const, TYPE const, TYPE &); \
799 template int Geometry::crossProduct<TYPE>(TYPE const *, TYPE const *, \
800 TYPE const *, TYPE const *, \
801 std::array<TYPE, 3> &); \
802 template int Geometry::crossProduct<TYPE>( \
803 TYPE const *, TYPE const *, TYPE *); \
804 template TYPE Geometry::distance<TYPE>( \
805 TYPE const *, TYPE const *, int const &); \
806 template TYPE Geometry::distance<TYPE>( \
807 std::vector<TYPE> const &, std::vector<TYPE> const &); \
808 template TYPE Geometry::distanceFlatten<TYPE>( \
809 std::vector<std::vector<TYPE>> const &, \
810 std::vector<std::vector<TYPE>> const &); \
811 template TYPE Geometry::dotProduct<TYPE>( \
812 TYPE const *, TYPE const *, TYPE const *, TYPE const *); \
813 template TYPE Geometry::dotProduct<TYPE>( \
814 TYPE const *, TYPE const *, int const &); \
815 template TYPE Geometry::dotProduct<TYPE>( \
816 std::vector<TYPE> const &, std::vector<TYPE> const &); \
817 template TYPE Geometry::dotProductFlatten<TYPE>( \
818 std::vector<std::vector<TYPE>> const &, \
819 std::vector<std::vector<TYPE>> const &); \
820 template bool Geometry::isPointInTriangle<TYPE>( \
821 TYPE const *, TYPE const *, TYPE const *, TYPE const *); \
822 template bool Geometry::isPointOnSegment<TYPE>(TYPE const &, TYPE const &, \
823 TYPE const &, TYPE const &, \
824 TYPE const &, TYPE const &); \
825 template bool Geometry::isPointOnSegment<TYPE>( \
826 TYPE const *, TYPE const *, TYPE const *, int const &); \
827 template bool Geometry::isTriangleColinear<TYPE>( \
828 TYPE const *, TYPE const *, TYPE const *, TYPE const *); \
829 template TYPE Geometry::magnitude<TYPE>(TYPE const *, int const &); \
830 template TYPE Geometry::magnitude<TYPE>(std::vector<TYPE> const &); \
831 template TYPE Geometry::magnitudeFlatten<TYPE>( \
832 std::vector<std::vector<TYPE>> const &); \
833 template TYPE Geometry::magnitude<TYPE>(TYPE const *, TYPE const *); \
834 template int Geometry::subtractVectors<TYPE>( \
835 TYPE const *, TYPE const *, TYPE *, int const &); \
836 template int Geometry::subtractVectors<TYPE>(std::vector<TYPE> const &, \
837 std::vector<TYPE> const &, \
838 std::vector<TYPE> &); \
840 TYPE const *, TYPE const *, TYPE *, int const &); \
841 template int Geometry::addVectors<TYPE>(std::vector<TYPE> const &, \
842 std::vector<TYPE> const &, \
843 std::vector<TYPE> &); \
845 std::vector<std::vector<TYPE>> const &, \
846 std::vector<std::vector<TYPE>> const &, std::vector<std::vector<TYPE>> &); \
848 std::vector<std::vector<std::vector<TYPE>>> const &, \
849 std::vector<std::vector<std::vector<TYPE>>> const &, \
850 std::vector<std::vector<TYPE>> &); \
851 template int Geometry::scaleVector<TYPE>( \
852 TYPE const *, TYPE const, TYPE *, int const &); \
853 template int Geometry::scaleVector<TYPE>( \
854 std::vector<TYPE> const &, TYPE const, std::vector<TYPE> &); \
855 template int Geometry::vectorProjection<TYPE>( \
856 TYPE const *, TYPE const *, TYPE *, int const &); \
857 template int Geometry::vectorProjection<TYPE>(std::vector<TYPE> const &, \
858 std::vector<TYPE> const &, \
859 std::vector<TYPE> &); \
861 std::vector<TYPE> const &, std::vector<TYPE> const &, std::vector<TYPE> &, \
862 std::vector<TYPE> &); \
863 template void Geometry::gramSchmidt<TYPE>( \
864 std::vector<std::vector<TYPE>> const &, std::vector<std::vector<TYPE>> &); \
865 template bool Geometry::isVectorUniform<TYPE>(std::vector<TYPE> const &); \
866 template bool Geometry::isVectorNull<TYPE>(std::vector<TYPE> const &); \
867 template bool Geometry::isVectorNullFlatten<TYPE>( \
868 std::vector<std::vector<TYPE>> const &); \
869 template int Geometry::flattenMultiDimensionalVector<TYPE>( \
870 std::vector<std::vector<TYPE>> const &, std::vector<TYPE> &); \
871 template int Geometry::multiFlattenMultiDimensionalVector<TYPE>( \
872 std::vector<std::vector<std::vector<TYPE>>> const &, \
873 std::vector<std::vector<TYPE>> &); \
874 template int Geometry::unflattenMultiDimensionalVector<TYPE>( \
875 std::vector<TYPE> const &, std::vector<std::vector<TYPE>> &, int const &); \
876 template void Geometry::matrixMultiplication<TYPE>( \
877 std::vector<std::vector<TYPE>> const &, \
878 std::vector<std::vector<TYPE>> const &, std::vector<std::vector<TYPE>> &); \
879 template void Geometry::subtractMatrices<TYPE>( \
880 std::vector<std::vector<TYPE>> const &, \
881 std::vector<std::vector<TYPE>> const &, std::vector<std::vector<TYPE>> &); \
883 std::vector<std::vector<TYPE>> const &, \
884 std::vector<std::vector<TYPE>> const &, std::vector<std::vector<TYPE>> &); \
885 template void Geometry::scaleMatrix<TYPE>( \
886 std::vector<std::vector<TYPE>> const &, TYPE const, \
887 std::vector<std::vector<TYPE>> &); \
888 template void Geometry::transposeMatrix<TYPE>( \
889 std::vector<std::vector<TYPE>> const &, std::vector<std::vector<TYPE>> &);
890
891// explicit specializations for float and double
#define GEOMETRY_SPECIALIZE(TYPE)
Definition Geometry.cpp:772
#define M_PI
Definition Os.h:50
int scaleVector(const T *a, const T factor, T *out, const int &dimension=3)
Definition Geometry.cpp:593
bool areVectorsColinear(const T *vA0, const T *vA1, const T *vB0, const T *vB1, std::array< T, 3 > *coefficients=nullptr, const T *tolerance=NULL)
Definition Geometry.cpp:27
int computeTriangleArea(const T *p0, const T *p1, const T *p2, T &area)
Definition Geometry.cpp:281
int addVectors(const T *a, const T *b, T *out, const int &dimension=3)
Definition Geometry.cpp:556
T dotProductFlatten(const std::vector< std::vector< T > > &vA, const std::vector< std::vector< T > > &vB)
Definition Geometry.cpp:412
int multiFlattenMultiDimensionalVector(const std::vector< std::vector< std::vector< T > > > &a, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:698
void subtractMatrices(const std::vector< std::vector< T > > &a, const std::vector< std::vector< T > > &b, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:734
bool isTriangleColinear2D(const T *pptA, const T *pptB, const T *pptC, const T tolerance)
Definition Geometry.cpp:130
T dotProduct(const T *vA0, const T *vA1, const T *vB0, const T *vB1)
Definition Geometry.cpp:388
void matrixMultiplication(const std::vector< std::vector< T > > &a, const std::vector< std::vector< T > > &b, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:723
bool isVectorNullFlatten(const std::vector< std::vector< T > > &a)
Definition Geometry.cpp:681
bool isPointOnSegment(const T &x, const T &y, const T &xA, const T &yA, const T &xB, const T &yB)
Definition Geometry.cpp:443
bool computeSegmentIntersection(const T &xA, const T &yA, const T &xB, const T &yB, const T &xC, const T &yC, const T &xD, const T &yD, T &x, T &y)
Definition Geometry.cpp:248
int computeBarycentricCoordinates(const T *p0, const T *p1, const T *p, std::array< T, 2 > &baryCentrics, const int &dimension=3)
Definition Geometry.cpp:142
void addMatrices(const std::vector< std::vector< T > > &a, const std::vector< std::vector< T > > &b, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:744
void gramSchmidt(const std::vector< std::vector< T > > &a, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:647
void scaleMatrix(const std::vector< std::vector< T > > &a, const T factor, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:754
int multiAddVectors(const std::vector< std::vector< T > > &a, const std::vector< std::vector< T > > &b, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:571
int flattenMultiDimensionalVector(const std::vector< std::vector< T > > &a, std::vector< T > &out)
Definition Geometry.cpp:688
bool isVectorUniform(const std::vector< T > &a)
Definition Geometry.cpp:665
int subtractVectors(const T *a, const T *b, T *out, const int &dimension=3)
Definition Geometry.cpp:538
int computeTriangleAreaFromSides(const T s0, const T s1, const T s2, T &area)
Definition Geometry.cpp:296
T magnitudeFlatten(const std::vector< std::vector< T > > &v)
Definition Geometry.cpp:519
bool isVectorNull(const std::vector< T > &a)
Definition Geometry.cpp:673
T1 pow(const T1 val, const T2 n)
Definition Geometry.h:454
bool isTriangleColinear(const T *p0, const T *p1, const T *p2, const T *tolerance=nullptr)
Definition Geometry.cpp:466
T angle(const T *vA0, const T *vA1, const T *vB0, const T *vB1)
Definition Geometry.cpp:13
void addVectorsProjection(const std::vector< T > &a, const std::vector< T > &b, std::vector< T > &a_out, std::vector< T > &b_out)
Definition Geometry.cpp:636
int computeTriangleAngles(const T *p0, const T *p1, const T *p2, std::array< T, 3 > &angles)
Definition Geometry.cpp:308
int computeTriangleAngleFromSides(const T s0, const T s1, const T s2, T &angle)
Definition Geometry.cpp:321
int vectorProjection(const T *a, const T *b, T *out, const int &dimension=3)
Definition Geometry.cpp:611
int crossProduct(const T *vA0, const T *vA1, const T *vB0, const T *vB1, std::array< T, 3 > &crossProduct)
Definition Geometry.cpp:332
void transposeMatrix(const std::vector< std::vector< T > > &a, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:764
int unflattenMultiDimensionalVector(const std::vector< T > &a, std::vector< std::vector< T > > &out, const int &no_columns=2)
Definition Geometry.cpp:708
T magnitude(const T *v, const int &dimension=3)
Definition Geometry.cpp:509
T distance(const T *p0, const T *p1, const int &dimension=3)
Definition Geometry.cpp:362
T angle2D(const T *vA0, const T *vA1, const T *vB0, const T *vB1)
Definition Geometry.cpp:20
T distanceFlatten(const std::vector< std::vector< T > > &p0, const std::vector< std::vector< T > > &p1)
Definition Geometry.cpp:379
int multiAddVectorsFlatten(const std::vector< std::vector< std::vector< T > > > &a, const std::vector< std::vector< std::vector< T > > > &b, std::vector< std::vector< T > > &out)
Definition Geometry.cpp:581
bool isPointInTriangle(const T *p0, const T *p1, const T *p2, const T *p)
Definition Geometry.cpp:421
The Topology ToolKit.