OGS 6.1.0-1699-ge946d4c5f
AnalyticalGeometry.cpp
Go to the documentation of this file.
1 
15 #include "AnalyticalGeometry.h"
16 
17 #include <algorithm>
18 #include <cmath>
19 #include <limits>
20 
21 #include <logog/include/logog.hpp>
22 
23 #include <Eigen/Dense>
24 
25 #include "BaseLib/StringTools.h"
26 
27 #include "Polyline.h"
28 #include "PointVec.h"
29 
31 
32 extern double orient2d(double *, double *, double *);
33 
34 namespace ExactPredicates
35 {
37  MathLib::Point3d const& b, MathLib::Point3d const& c)
38 {
39  return orient2d(const_cast<double*>(a.getCoords()),
40  const_cast<double*>(b.getCoords()),
41  const_cast<double*>(c.getCoords()));
42 }
43 }
44 
45 namespace GeoLib
46 {
47 Orientation getOrientation(const double& p0_x, const double& p0_y, const double& p1_x,
48  const double& p1_y, const double& p2_x, const double& p2_y)
49 {
50  double h1((p1_x - p0_x) * (p2_y - p0_y));
51  double h2((p2_x - p0_x) * (p1_y - p0_y));
52 
53  double tol(std::numeric_limits<double>::epsilon());
54  if (fabs(h1 - h2) <= tol * std::max(fabs(h1), fabs(h2)))
55  return COLLINEAR;
56  if (h1 - h2 > 0.0)
57  return CCW;
58 
59  return CW;
60 }
61 
63  const GeoLib::Point* p2)
64 {
65  return getOrientation((*p0)[0], (*p0)[1], (*p1)[0], (*p1)[1], (*p2)[0], (*p2)[1]);
66 }
67 
69 {
70  const double eps(std::numeric_limits<double>::epsilon());
71 
72  // check degenerated cases
73  if (v.getLength() < eps)
74  return false;
75 
76  if (w.getLength() < eps)
77  return false;
78 
79  v.normalize();
80  w.normalize();
81 
82  bool parallel(true);
83  if (std::abs(v[0]-w[0]) > eps)
84  parallel = false;
85  if (std::abs(v[1]-w[1]) > eps)
86  parallel = false;
87  if (std::abs(v[2]-w[2]) > eps)
88  parallel = false;
89 
90  if (! parallel) {
91  parallel = true;
92  // change sense of direction of v_normalised
93  v *= -1.0;
94  // check again
95  if (std::abs(v[0]-w[0]) > eps)
96  parallel = false;
97  if (std::abs(v[1]-w[1]) > eps)
98  parallel = false;
99  if (std::abs(v[2]-w[2]) > eps)
100  parallel = false;
101  }
102 
103  return parallel;
104 }
105 
107  GeoLib::LineSegment const& s1,
108  GeoLib::Point& s)
109 {
110  GeoLib::Point const& a{s0.getBeginPoint()};
111  GeoLib::Point const& b{s0.getEndPoint()};
112  GeoLib::Point const& c{s1.getBeginPoint()};
113  GeoLib::Point const& d{s1.getEndPoint()};
114 
115  if (!isCoplanar(a, b, c, d))
116  return false;
117 
118  // handle special cases here to avoid computing intersection numerical
119  if (MathLib::sqrDist(a, c) < std::numeric_limits<double>::epsilon() ||
120  MathLib::sqrDist(a, d) < std::numeric_limits<double>::epsilon()) {
121  s = a;
122  return true;
123  }
124  if (MathLib::sqrDist(b, c) < std::numeric_limits<double>::epsilon() ||
125  MathLib::sqrDist(b, d) < std::numeric_limits<double>::epsilon()) {
126  s = b;
127  return true;
128  }
129 
130  MathLib::Vector3 const v(a, b);
131  MathLib::Vector3 const w(c, d);
132  MathLib::Vector3 const qp(a, c);
133  MathLib::Vector3 const pq(c, a);
134 
135  auto isLineSegmentIntersectingAB = [&v](MathLib::Vector3 const& ap,
136  std::size_t i)
137  {
138  // check if p is located at v=(a,b): (ap = t*v, t in [0,1])
139  return 0.0 <= ap[i] / v[i] && ap[i] / v[i] <= 1.0;
140  };
141 
142  if (parallel(v,w)) { // original line segments (a,b) and (c,d) are parallel
143  if (parallel(pq,v)) { // line segment (a,b) and (a,c) are also parallel
144  // Here it is already checked that the line segments (a,b) and (c,d)
145  // are parallel. At this point it is also known that the line
146  // segment (a,c) is also parallel to (a,b). In that case it is
147  // possible to express c as c(t) = a + t * (b-a) (analog for the
148  // point d). Since the evaluation of all three coordinate equations
149  // (x,y,z) have to lead to the same solution for the parameter t it
150  // is sufficient to evaluate t only once.
151 
152  // Search id of coordinate with largest absolute value which is will
153  // be used in the subsequent computations. This prevents division by
154  // zero in case the line segments are parallel to one of the
155  // coordinate axis.
156  std::size_t i_max(std::abs(v[0]) <= std::abs(v[1]) ? 1 : 0);
157  i_max = std::abs(v[i_max]) <= std::abs(v[2]) ? 2 : i_max;
158  if (isLineSegmentIntersectingAB(qp, i_max)) {
159  s = c;
160  return true;
161  }
162  MathLib::Vector3 const ad(a, d);
163  if (isLineSegmentIntersectingAB(ad, i_max)) {
164  s = d;
165  return true;
166  }
167  return false;
168  }
169  return false;
170  }
171 
172  // general case
173  const double sqr_len_v(v.getSqrLength());
174  const double sqr_len_w(w.getSqrLength());
175 
176  Eigen::Matrix2d mat;
177  mat(0,0) = sqr_len_v;
178  mat(0,1) = -1.0 * MathLib::scalarProduct(v,w);
179  mat(1,1) = sqr_len_w;
180  mat(1,0) = mat(0,1);
181 
182  Eigen::Vector2d rhs;
183  rhs << MathLib::scalarProduct(v, qp), MathLib::scalarProduct(w, pq);
184 
185  rhs = mat.partialPivLu().solve(rhs);
186 
187  // no theory for the following tolerances, determined by testing
188  // lower tolerance: little bit smaller than zero
189  const double l(-1.0*std::numeric_limits<float>::epsilon());
190  // upper tolerance a little bit greater than one
191  const double u(1.0+std::numeric_limits<float>::epsilon());
192  if (rhs[0] < l || u < rhs[0] || rhs[1] < l || u < rhs[1]) {
193  return false;
194  }
195 
196  // compute points along line segments with minimal distance
197  GeoLib::Point const p0(a[0]+rhs[0]*v[0], a[1]+rhs[0]*v[1], a[2]+rhs[0]*v[2]);
198  GeoLib::Point const p1(c[0]+rhs[1]*w[0], c[1]+rhs[1]*w[1], c[2]+rhs[1]*w[2]);
199 
200  double const min_dist(sqrt(MathLib::sqrDist(p0, p1)));
201  double const min_seg_len(std::min(sqrt(sqr_len_v), sqrt(sqr_len_w)));
202  if (min_dist < min_seg_len * 1e-6) {
203  s[0] = 0.5 * (p0[0] + p1[0]);
204  s[1] = 0.5 * (p0[1] + p1[1]);
205  s[2] = 0.5 * (p0[2] + p1[2]);
206  return true;
207  }
208 
209  return false;
210 }
211 
215  GeoLib::Point& intersection_pnt)
216 {
217  std::size_t const n_segs(ply->getNumberOfSegments());
218  // Neighbouring segments always intersects at a common vertex. The algorithm
219  // checks for intersections of non-neighbouring segments.
220  for (seg_it0 = ply->begin(); seg_it0 != ply->end() - 2; ++seg_it0)
221  {
222  seg_it1 = seg_it0+2;
223  std::size_t const seg_num_0 = seg_it0.getSegmentNumber();
224  for ( ; seg_it1 != ply->end(); ++seg_it1) {
225  // Do not check first and last segment, because they are
226  // neighboured.
227  if (!(seg_num_0 == 0 && seg_it1.getSegmentNumber() == n_segs - 1)) {
228  if (lineSegmentIntersect(*seg_it0, *seg_it1, intersection_pnt)) {
229  return true;
230  }
231  }
232  }
233  }
234  return false;
235 }
236 
238 {
239  // *** some frequently used terms ***
240  // n_1^2 + n_2^2
241  const double h0(plane_normal[0] * plane_normal[0] + plane_normal[1] * plane_normal[1]);
242  // 1 / sqrt (n_1^2 + n_2^2)
243  const double h1(1.0 / sqrt(h0));
244  // 1 / sqrt (n_1^2 + n_2^2 + n_3^2)
245  const double h2(1.0 / sqrt(h0 + plane_normal[2] * plane_normal[2]));
246 
247  // calc rotation matrix
248  rot_mat(0, 0) = plane_normal[1] * h1;
249  rot_mat(0, 1) = -plane_normal[0] * h1;
250  rot_mat(0, 2) = 0.0;
251  rot_mat(1, 0) = plane_normal[0] * h2;
252  rot_mat(1, 1) = plane_normal[1] * h2;
253  rot_mat(1, 2) = plane_normal[2] * h2;
254  rot_mat(2, 0) = plane_normal[0] * plane_normal[2] * h1 * h2;
255  rot_mat(2, 1) = plane_normal[1] * plane_normal[2] * h1 * h2;
256  rot_mat(2, 2) = -sqrt(h0) * h2;
257 }
258 
259 void rotatePoints(MathLib::DenseMatrix<double> const& rot_mat, std::vector<GeoLib::Point*> &pnts)
260 {
261  rotatePoints(rot_mat, pnts.begin(), pnts.end());
262 }
263 
264 MathLib::DenseMatrix<double> rotatePointsToXY(std::vector<GeoLib::Point*>& pnts)
265 {
266  return rotatePointsToXY(pnts.begin(), pnts.end(), pnts.begin(), pnts.end());
267 }
268 
269 void rotatePointsToXZ(std::vector<GeoLib::Point*> &pnts)
270 {
271  assert(pnts.size()>2);
272  // calculate supporting plane
273  MathLib::Vector3 plane_normal;
274  double d;
275  // compute the plane normal
276  GeoLib::getNewellPlane(pnts, plane_normal, d);
277 
278  const double tol (std::numeric_limits<double>::epsilon());
279  if (std::abs(plane_normal[0]) > tol || std::abs(plane_normal[1]) > tol) {
280  // rotate copied points into x-z-plane
281  MathLib::DenseMatrix<double> rot_mat(3, 3);
282  computeRotationMatrixToXZ(plane_normal, rot_mat);
283  rotatePoints(rot_mat, pnts);
284  }
285 
286  for (auto & pnt : pnts)
287  (*pnt)[1] = 0.0; // should be -= d but there are numerical errors
288 }
289 
290 std::unique_ptr<GeoLib::Point> triangleLineIntersection(
291  MathLib::Point3d const& a, MathLib::Point3d const& b,
292  MathLib::Point3d const& c, MathLib::Point3d const& p,
293  MathLib::Point3d const& q)
294 {
295  const MathLib::Vector3 pq(p, q);
296  const MathLib::Vector3 pa(p, a);
297  const MathLib::Vector3 pb(p, b);
298  const MathLib::Vector3 pc(p, c);
299 
300  double u (MathLib::scalarTriple(pq, pc, pb));
301  if (u<0) return nullptr;
302  double v (MathLib::scalarTriple(pq, pa, pc));
303  if (v<0) return nullptr;
304  double w (MathLib::scalarTriple(pq, pb, pa));
305  if (w<0) return nullptr;
306 
307  const double denom (1.0/(u+v+w));
308  u*=denom;
309  v*=denom;
310  w*=denom;
311  return std::make_unique<GeoLib::Point>(u * a[0] + v * b[0] + w * c[0],
312  u * a[1] + v * b[1] + w * c[1],
313  u * a[2] + v * b[2] + w * c[2]);
314 }
315 
317  std::vector<GeoLib::Polyline*> & plys)
318 {
319  auto computeSegmentIntersections = [&pnt_vec](GeoLib::Polyline& poly0,
320  GeoLib::Polyline& poly1)
321  {
322  for (auto seg0_it(poly0.begin()); seg0_it != poly0.end(); ++seg0_it)
323  {
324  for (auto seg1_it(poly1.begin()); seg1_it != poly1.end(); ++seg1_it)
325  {
326  GeoLib::Point s(0.0, 0.0, 0.0, pnt_vec.size());
327  if (lineSegmentIntersect(*seg0_it, *seg1_it, s))
328  {
329  std::size_t const id(
330  pnt_vec.push_back(new GeoLib::Point(s)));
331  poly0.insertPoint(seg0_it.getSegmentNumber() + 1, id);
332  poly1.insertPoint(seg1_it.getSegmentNumber() + 1, id);
333  }
334  }
335  }
336  };
337 
338  for (auto it0(plys.begin()); it0 != plys.end(); ++it0) {
339  auto it1(it0);
340  ++it1;
341  for (; it1 != plys.end(); ++it1) {
342  computeSegmentIntersections(*(*it0), *(*it1));
343  }
344  }
345 }
346 
348  MathLib::Vector3 & plane_normal)
349 {
350  // 1 copy all points
351  auto* polygon_pnts(new std::vector<GeoLib::Point*>);
352  for (std::size_t k(0); k < polygon_in.getNumberOfPoints(); k++)
353  polygon_pnts->push_back (new GeoLib::Point (*(polygon_in.getPoint(k))));
354 
355  // 2 rotate points
356  double d_polygon (0.0);
357  GeoLib::getNewellPlane (*polygon_pnts, plane_normal, d_polygon);
358  MathLib::DenseMatrix<double> rot_mat(3,3);
359  GeoLib::computeRotationMatrixToXY(plane_normal, rot_mat);
360  GeoLib::rotatePoints(rot_mat, *polygon_pnts);
361 
362  // 3 set z coord to zero
363  std::for_each(polygon_pnts->begin(), polygon_pnts->end(),
364  [] (GeoLib::Point* p) { (*p)[2] = 0.0; }
365  );
366 
367  // 4 create new polygon
368  GeoLib::Polyline rot_polyline(*polygon_pnts);
369  for (std::size_t k(0); k < polygon_in.getNumberOfPoints(); k++)
370  rot_polyline.addPoint(k);
371  rot_polyline.addPoint(0);
372  return GeoLib::Polygon(rot_polyline);
373 }
374 
375 std::vector<MathLib::Point3d> lineSegmentIntersect2d(
376  GeoLib::LineSegment const& ab, GeoLib::LineSegment const& cd)
377 {
378  GeoLib::Point const& a{ab.getBeginPoint()};
379  GeoLib::Point const& b{ab.getEndPoint()};
380  GeoLib::Point const& c{cd.getBeginPoint()};
381  GeoLib::Point const& d{cd.getEndPoint()};
382 
383  double const orient_abc(ExactPredicates::getOrientation2d(a, b, c));
384  double const orient_abd(ExactPredicates::getOrientation2d(a, b, d));
385 
386  // check if the segment (cd) lies on the left or on the right of (ab)
387  if ((orient_abc > 0 && orient_abd > 0) || (orient_abc < 0 && orient_abd < 0)) {
388  return std::vector<MathLib::Point3d>();
389  }
390 
391  // check: (cd) and (ab) are on the same line
392  if (orient_abc == 0.0 && orient_abd == 0.0) {
393  double const eps(std::numeric_limits<double>::epsilon());
394  if (MathLib::sqrDist2d(a,c) < eps && MathLib::sqrDist2d(b,d) < eps)
395  return {{ a, b }};
396  if (MathLib::sqrDist2d(a,d) < eps && MathLib::sqrDist2d(b,c) < eps)
397  return {{ a, b }};
398 
399  // Since orient_ab and orient_abd vanish, a, b, c, d are on the same
400  // line and for this reason it is enough to check the x-component.
401  auto isPointOnSegment = [](double q, double p0, double p1)
402  {
403  double const t((q - p0) / (p1 - p0));
404  return 0 <= t && t <= 1;
405  };
406 
407  // check if c in (ab)
408  if (isPointOnSegment(c[0], a[0], b[0])) {
409  // check if a in (cd)
410  if (isPointOnSegment(a[0], c[0], d[0])) {
411  return {{a, c}};
412  }
413  // check b == c
414  if (MathLib::sqrDist2d(b,c) < eps) {
415  return {{b}};
416  }
417  // check if b in (cd)
418  if (isPointOnSegment(b[0], c[0], d[0])) {
419  return {{b, c}};
420  }
421  // check d in (ab)
422  if (isPointOnSegment(d[0], a[0], b[0])) {
423  return {{c, d}};
424  }
425  std::stringstream err;
426  err.precision(std::numeric_limits<double>::digits10);
427  err << ab << " x " << cd;
428  OGS_FATAL(
429  "The case of parallel line segments (%s) is not handled yet. "
430  "Aborting.",
431  err.str().c_str());
432  }
433 
434  // check if d in (ab)
435  if (isPointOnSegment(d[0], a[0], b[0])) {
436  // check if a in (cd)
437  if (isPointOnSegment(a[0], c[0], d[0])) {
438  return {{a, d}};
439  }
440  // check if b==d
441  if (MathLib::sqrDist2d(b, d) < eps) {
442  return {{b}};
443  }
444  // check if b in (cd)
445  if (isPointOnSegment(b[0], c[0], d[0])) {
446  return {{b, d}};
447  }
448  // d in (ab), b not in (cd): check c in (ab)
449  if (isPointOnSegment(c[0], a[0], b[0])) {
450  return {{c, d}};
451  }
452 
453  std::stringstream err;
454  err.precision(std::numeric_limits<double>::digits10);
455  err << ab << " x " << cd;
456  OGS_FATAL(
457  "The case of parallel line segments (%s) "
458  "is not handled yet. Aborting.",
459  err.str().c_str());
460  }
461  return std::vector<MathLib::Point3d>();
462  }
463 
464  // precondition: points a, b, c are collinear
465  // the function checks if the point c is onto the line segment (a,b)
466  auto isCollinearPointOntoLineSegment = [](MathLib::Point3d const& a,
467  MathLib::Point3d const& b,
468  MathLib::Point3d const& c) {
469  if (b[0] - a[0] != 0)
470  {
471  double const t = (c[0] - a[0]) / (b[0] - a[0]);
472  return 0.0 <= t && t <= 1.0;
473  }
474  if (b[1] - a[1] != 0)
475  {
476  double const t = (c[1] - a[1]) / (b[1] - a[1]);
477  return 0.0 <= t && t <= 1.0;
478  }
479  if (b[2] - a[2] != 0)
480  {
481  double const t = (c[2] - a[2]) / (b[2] - a[2]);
482  return 0.0 <= t && t <= 1.0;
483  }
484  return false;
485  };
486 
487  if (orient_abc == 0.0) {
488  if (isCollinearPointOntoLineSegment(a,b,c))
489  return {{c}};
490  return std::vector<MathLib::Point3d>();
491  }
492 
493  if (orient_abd == 0.0) {
494  if (isCollinearPointOntoLineSegment(a,b,d))
495  return {{d}};
496  return std::vector<MathLib::Point3d>();
497  }
498 
499  // check if the segment (ab) lies on the left or on the right of (cd)
500  double const orient_cda(ExactPredicates::getOrientation2d(c, d, a));
501  double const orient_cdb(ExactPredicates::getOrientation2d(c, d, b));
502  if ((orient_cda > 0 && orient_cdb > 0) || (orient_cda < 0 && orient_cdb < 0)) {
503  return std::vector<MathLib::Point3d>();
504  }
505 
506  // at this point it is sure that there is an intersection and the system of
507  // linear equations will be invertible
508  // solve the two linear equations (b-a, c-d) (t, s)^T = (c-a) simultaneously
509  Eigen::Matrix2d mat;
510  mat(0,0) = b[0]-a[0];
511  mat(0,1) = c[0]-d[0];
512  mat(1,0) = b[1]-a[1];
513  mat(1,1) = c[1]-d[1];
514  Eigen::Vector2d rhs;
515  rhs << c[0] - a[0], c[1] - a[1];
516 
517  rhs = mat.partialPivLu().solve(rhs);
518  if (0 <= rhs[1] && rhs[1] <= 1.0) {
519  return { MathLib::Point3d{std::array<double,3>{{
520  c[0]+rhs[1]*(d[0]-c[0]), c[1]+rhs[1]*(d[1]-c[1]),
521  c[2]+rhs[1]*(d[2]-c[2])}} } };
522  }
523  return std::vector<MathLib::Point3d>(); // parameter s not in the valid
524  // range
525 }
526 
528  MathLib::Point3d const& seg_beg_pnt,
529  std::vector<GeoLib::LineSegment>& sub_segments)
530 {
531  double const eps(std::numeric_limits<double>::epsilon());
532 
533  auto findNextSegment = [&eps](
534  MathLib::Point3d const& seg_beg_pnt,
535  std::vector<GeoLib::LineSegment>& sub_segments,
536  std::vector<GeoLib::LineSegment>::iterator& sub_seg_it)
537  {
538  if (sub_seg_it == sub_segments.end())
539  return;
540  // find appropriate segment for the given segment begin point
541  auto act_beg_seg_it = std::find_if(
542  sub_seg_it, sub_segments.end(),
543  [&seg_beg_pnt, &eps](GeoLib::LineSegment const& seg)
544  {
545  return MathLib::sqrDist(seg_beg_pnt, seg.getBeginPoint()) < eps ||
546  MathLib::sqrDist(seg_beg_pnt, seg.getEndPoint()) < eps;
547  });
548  if (act_beg_seg_it == sub_segments.end())
549  return;
550  // if necessary correct orientation of segment, i.e. swap beg and end
551  if (MathLib::sqrDist(seg_beg_pnt, act_beg_seg_it->getEndPoint()) <
552  MathLib::sqrDist(seg_beg_pnt, act_beg_seg_it->getBeginPoint()))
553  std::swap(act_beg_seg_it->getBeginPoint(),
554  act_beg_seg_it->getEndPoint());
555  assert(sub_seg_it != sub_segments.end());
556  // exchange segments within the container
557  if (sub_seg_it != act_beg_seg_it)
558  std::swap(*sub_seg_it, *act_beg_seg_it);
559  };
560 
561  // find start segment
562  auto seg_it = sub_segments.begin();
563  findNextSegment(seg_beg_pnt, sub_segments, seg_it);
564 
565  while (seg_it != sub_segments.end())
566  {
567  MathLib::Point3d & new_seg_beg_pnt(seg_it->getEndPoint());
568  seg_it++;
569  if (seg_it != sub_segments.end())
570  findNextSegment(new_seg_beg_pnt, sub_segments, seg_it);
571  }
572 }
573 
574 } // end namespace GeoLib
double getOrientation2d(MathLib::Point3d const &a, MathLib::Point3d const &b, MathLib::Point3d const &c)
Definition of the PolyLine class.
void computeAndInsertAllIntersectionPoints(GeoLib::PointVec &pnt_vec, std::vector< GeoLib::Polyline *> &plys)
MathLib::DenseMatrix< double > rotatePointsToXY(std::vector< GeoLib::Point *> &pnts)
double sqrDist(const double *p0, const double *p1)
Definition: MathTools.h:77
void rotatePoints(MathLib::DenseMatrix< double > const &rot_mat, std::vector< GeoLib::Point *> &pnts)
Definition of string helper functions.
std::vector< MathLib::Point3d > lineSegmentIntersect2d(GeoLib::LineSegment const &ab, GeoLib::LineSegment const &cd)
GeoLib::Point const & getBeginPoint() const
Definition: LineSegment.cpp:62
std::size_t size() const
Definition: TemplateVec.h:102
std::size_t getNumberOfSegments() const
Definition: Polyline.cpp:185
void rotatePoints(MathLib::DenseMatrix< double > const &rot_mat, InputIterator pnts_begin, InputIterator pnts_end)
double orient2d(double *, double *, double *)
double scalarTriple(MathLib::Vector3 const &u, MathLib::Vector3 const &v, MathLib::Vector3 const &w)
Calculates the scalar triple (u x v) . w.
Definition: Vector3.cpp:16
double sqrDist2d(MathLib::Point3d const &p0, MathLib::Point3d const &p1)
Definition: Point3d.h:58
double getSqrLength(void) const
Returns the squared length.
Definition: Vector3.h:121
bool lineSegmentsIntersect(const GeoLib::Polyline *ply, GeoLib::Polyline::SegmentIterator &seg_it0, GeoLib::Polyline::SegmentIterator &seg_it1, GeoLib::Point &intersection_pnt)
const Point * getPoint(std::size_t i) const
returns the i-th point contained in the polyline
Definition: Polyline.cpp:249
Definition of analytical geometry functions.
bool lineSegmentIntersect(GeoLib::LineSegment const &s0, GeoLib::LineSegment const &s1, GeoLib::Point &s)
Orientation getOrientation(const GeoLib::Point *p0, const GeoLib::Point *p1, const GeoLib::Point *p2)
void computeRotationMatrixToXZ(MathLib::Vector3 const &plane_normal, MathLib::DenseMatrix< double > &rot_mat)
std::size_t getNumberOfPoints() const
Definition: Polyline.cpp:180
std::unique_ptr< GeoLib::Point > triangleLineIntersection(MathLib::Point3d const &a, MathLib::Point3d const &b, MathLib::Point3d const &c, MathLib::Point3d const &p, MathLib::Point3d const &q)
Definition of the GEOObjects class.
Definition: BaseItem.h:20
GeoLib::Polygon rotatePolygonToXY(GeoLib::Polygon const &polygon_in, MathLib::Vector3 &plane_normal)
void sortSegments(MathLib::Point3d const &seg_beg_pnt, std::vector< GeoLib::LineSegment > &sub_segments)
bool isCoplanar(const MathLib::Point3d &a, const MathLib::Point3d &b, const MathLib::Point3d &c, const MathLib::Point3d &d)
Checks if the four given points are located on a plane.
static const double q
T scalarProduct(T const *const v0, T const *const v1)
Definition: MathTools.h:28
GeoLib::Point const & getEndPoint() const
Definition: LineSegment.cpp:72
const T * getCoords() const
Definition: TemplatePoint.h:77
virtual bool insertPoint(std::size_t pos, std::size_t pnt_id)
Definition: Polyline.cpp:69
void computeRotationMatrixToXY(MathLib::Vector3 const &n, T_MATRIX &rot_mat)
void getNewellPlane(InputIterator pnts_begin, InputIterator pnts_end, MathLib::Vector3 &plane_normal, double &d)
Class Polyline consists mainly of a reference to a point vector and a vector that stores the indices ...
Definition: Polyline.h:50
SegmentIterator begin() const
Definition: Polyline.h:194
double getLength(void) const
Returns the length.
Definition: Vector3.h:127
This class manages pointers to Points in a std::vector along with a name. It also handles the deletin...
Definition: PointVec.h:39
#define OGS_FATAL(fmt,...)
Definition: Error.h:71
void rotatePointsToXZ(std::vector< GeoLib::Point *> &pnts)
std::size_t push_back(Point *pnt)
Definition: PointVec.cpp:124
std::size_t getSegmentNumber() const
Definition: Polyline.cpp:470
Definition of the PointVec class.
bool parallel(MathLib::Vector3 v, MathLib::Vector3 w)
SegmentIterator end() const
Definition: Polyline.h:199