OGS

Detailed Description

Class Polyline consists mainly of a reference to a point vector and a vector that stores the indices in the point vector. A polyline consists of at least one line segment. The polyline is specified by the points of the line segments. The class Polyline stores ids of pointers to the points in the _ply_pnt_ids vector.

Definition at line 52 of file Polyline.h.

#include <Polyline.h>

Inheritance diagram for GeoLib::Polyline:
[legend]
Collaboration diagram for GeoLib::Polyline:
[legend]

Classes

class  SegmentIterator
 

Public Member Functions

 Polyline (const std::vector< Point * > &pnt_vec)
 
 Polyline (const Polyline &ply)
 
Polylineoperator= (Polyline const &other)=delete
 
 ~Polyline () override=default
 
GEOTYPE getGeoType () const override
 return a geometry type More...
 
virtual bool addPoint (std::size_t pnt_id)
 
virtual bool insertPoint (std::size_t pos, std::size_t pnt_id)
 
void removePoint (std::size_t pos)
 
void closePolyline ()
 
std::size_t getNumberOfPoints () const
 
std::size_t getNumberOfSegments () const
 
bool isClosed () const
 
bool isCoplanar () const
 
bool isPointIDInPolyline (std::size_t pnt_id) const
 
std::size_t getPointID (std::size_t i) const
 
void setPointID (std::size_t idx, std::size_t id)
 
const PointgetPoint (std::size_t i) const
 returns the i-th point contained in the polyline More...
 
SegmentIterator begin () const
 
SegmentIterator end () const
 
std::vector< Point * > const & getPointsVec () const
 
double getDistanceAlongPolyline (const MathLib::Point3d &pnt, const double epsilon_radius) const
 
- Public Member Functions inherited from GeoLib::GeoObject
virtual ~GeoObject ()=default
 

Static Public Member Functions

static PolylineconstructPolylineFromSegments (const std::vector< Polyline * > &ply_vec, double prox=0.0)
 

Protected Member Functions

Location getLocationOfPoint (std::size_t k, MathLib::Point3d const &pnt) const
 

Protected Attributes

const std::vector< Point * > & _ply_pnts
 
std::vector< std::size_t > _ply_pnt_ids
 

Private Member Functions

LineSegment getSegment (std::size_t i) const
 
LineSegment getSegment (std::size_t i)
 

Friends

class Polygon
 

Constructor & Destructor Documentation

◆ Polyline() [1/2]

GeoLib::Polyline::Polyline ( const std::vector< Point * > &  pnt_vec)
explicit

constructor

Parameters
pnt_veca reference to the point vector

Definition at line 28 of file Polyline.cpp.

28 : _ply_pnts(pnt_vec) {}
const std::vector< Point * > & _ply_pnts
Definition: Polyline.h:227

Referenced by constructPolylineFromSegments().

◆ Polyline() [2/2]

GeoLib::Polyline::Polyline ( const Polyline ply)

Copy constructor

Parameters
plyPolyline

Definition at line 30 of file Polyline.cpp.

31  : _ply_pnts(ply._ply_pnts), _ply_pnt_ids(ply._ply_pnt_ids)
32 {
33 }
std::vector< std::size_t > _ply_pnt_ids
Definition: Polyline.h:229

◆ ~Polyline()

GeoLib::Polyline::~Polyline ( )
overridedefault

Member Function Documentation

◆ addPoint()

bool GeoLib::Polyline::addPoint ( std::size_t  pnt_id)
virtual

Adds an id of a point at the end of the polyline if and only if the resulting segment won't be empty. The id have to be inside the (internal) _ply_pnts vector the polyline is based on.

Returns
If the point could be added the return value is true. If the addition of the point would result in empty line segment false is returned.

Reimplemented in GeoLib::PolylineWithSegmentMarker, and GeoLib::PolygonWithSegmentMarker.

Definition at line 35 of file Polyline.cpp.

36 {
37  if (pnt_id >= _ply_pnts.size())
38  {
39  return false;
40  }
41  std::size_t const n_pnts(_ply_pnt_ids.size());
42 
43  // don't insert point if this would result in identical IDs for two adjacent
44  // points
45  if (n_pnts > 0 && _ply_pnt_ids.back() == pnt_id)
46  {
47  return false;
48  }
49 
50  _ply_pnt_ids.push_back(pnt_id);
51 
52  return true;
53 }

References _ply_pnt_ids, and _ply_pnts.

Referenced by GeoLib::PolygonWithSegmentMarker::addPoint(), GeoLib::PolylineWithSegmentMarker::addPoint(), closePolyline(), createPolyline(), insertPoint(), MeshGeoToolsLib::markNodesOutSideOfPolygon(), mergeGeometries(), FileIO::SwmmInterface::readLinksAsPolylines(), FileIO::Legacy::readPolylinePointVector(), and GeoLib::Polygon::splitPolygonAtPoint().

◆ begin()

◆ closePolyline()

void GeoLib::Polyline::closePolyline ( )

Closes a polyline by adding a line segment that connects start- and end-point.

Definition at line 309 of file Polyline.cpp.

310 {
311  if (getNumberOfPoints() < 2)
312  {
313  ERR("Polyline::closePolyline(): Input polyline needs to be composed of "
314  "at least three points.");
315  }
316  if (!isClosed())
317  {
318  addPoint(getPointID(0));
319  }
320 }
void ERR(char const *fmt, Args const &... args)
Definition: Logging.h:42
std::size_t getPointID(std::size_t i) const
Definition: Polyline.cpp:160
std::size_t getNumberOfPoints() const
Definition: Polyline.cpp:109
bool isClosed() const
Definition: Polyline.cpp:119
virtual bool addPoint(std::size_t pnt_id)
Definition: Polyline.cpp:35

References addPoint(), ERR(), getNumberOfPoints(), getPointID(), and isClosed().

Referenced by GEOModels::connectPolylineSegments().

◆ constructPolylineFromSegments()

Polyline * GeoLib::Polyline::constructPolylineFromSegments ( const std::vector< Polyline * > &  ply_vec,
double  prox = 0.0 
)
static

Constructs one polyline from a vector of connected polylines. All polylines in this vector need to reference the same point vector.

Definition at line 197 of file Polyline.cpp.

199 {
200  std::size_t nLines = ply_vec.size();
201 
202  auto* new_ply = new Polyline(*ply_vec[0]);
203  std::vector<GeoLib::Point*> pnt_vec(new_ply->getPointsVec());
204 
205  std::vector<Polyline*> local_ply_vec;
206  for (std::size_t i = 1; i < nLines; i++)
207  {
208  local_ply_vec.push_back(ply_vec[i]);
209  }
210 
211  while (!local_ply_vec.empty())
212  {
213  bool ply_found(false);
214  prox *= prox; // square distance once to save time later
215  for (auto it = local_ply_vec.begin(); it != local_ply_vec.end(); ++it)
216  {
217  if (pnt_vec == (*it)->getPointsVec())
218  {
219  std::size_t nPoints((*it)->getNumberOfPoints());
220 
221  // if (new_ply->getPointID(0) == (*it)->getPointID(0))
222  if (pointsAreIdentical(pnt_vec, new_ply->getPointID(0),
223  (*it)->getPointID(0), prox))
224  {
225  auto* tmp = new Polyline((*it)->getPointsVec());
226  for (std::size_t k = 0; k < nPoints; k++)
227  {
228  tmp->addPoint((*it)->getPointID(nPoints - k - 1));
229  }
230 
231  std::size_t new_ply_size(new_ply->getNumberOfPoints());
232  for (std::size_t k = 1; k < new_ply_size; k++)
233  {
234  tmp->addPoint(new_ply->getPointID(k));
235  }
236  delete new_ply;
237  new_ply = tmp;
238  ply_found = true;
239  }
240  // else if (new_ply->getPointID(0) ==
241  // (*it)->getPointID(nPoints-1))
242  else if (pointsAreIdentical(pnt_vec, new_ply->getPointID(0),
243  (*it)->getPointID(nPoints - 1),
244  prox))
245  {
246  auto* tmp = new Polyline(**it);
247  std::size_t new_ply_size(new_ply->getNumberOfPoints());
248  for (std::size_t k = 1; k < new_ply_size; k++)
249  {
250  tmp->addPoint(new_ply->getPointID(k));
251  }
252  delete new_ply;
253  new_ply = tmp;
254  ply_found = true;
255  }
256  // else if (new_ply->getPointID(new_ply->getNumberOfPoints()-1)
257  // == (*it)->getPointID(0))
258  else if (pointsAreIdentical(
259  pnt_vec,
260  new_ply->getPointID(new_ply->getNumberOfPoints() -
261  1),
262  (*it)->getPointID(0), prox))
263  {
264  for (std::size_t k = 1; k < nPoints; k++)
265  {
266  new_ply->addPoint((*it)->getPointID(k));
267  }
268  ply_found = true;
269  }
270  // else if (new_ply->getPointID(new_ply->getNumberOfPoints()-1)
271  // == (*it)->getPointID(nPoints-1))
272  else if (pointsAreIdentical(
273  pnt_vec,
274  new_ply->getPointID(new_ply->getNumberOfPoints() -
275  1),
276  (*it)->getPointID(nPoints - 1), prox))
277  {
278  for (std::size_t k = 1; k < nPoints; k++)
279  {
280  new_ply->addPoint((*it)->getPointID(nPoints - k - 1));
281  }
282  ply_found = true;
283  }
284  if (ply_found)
285  {
286  local_ply_vec.erase(it);
287  break;
288  }
289  }
290  else
291  {
292  ERR("Error in Polyline::contructPolylineFromSegments() - Line "
293  "segments use different point vectors.");
294  }
295  }
296 
297  if (!ply_found)
298  {
299  ERR("Error in Polyline::contructPolylineFromSegments() - Not all "
300  "segments are connected.");
301  delete new_ply;
302  new_ply = nullptr;
303  break;
304  }
305  }
306  return new_ply;
307 }
Polyline(const std::vector< Point * > &pnt_vec)
Definition: Polyline.cpp:28
bool pointsAreIdentical(const std::vector< Point * > &pnt_vec, std::size_t i, std::size_t j, double prox)
Definition: Polyline.cpp:567

References Polyline(), ERR(), and GeoLib::pointsAreIdentical().

Referenced by GEOModels::connectPolylineSegments().

◆ end()

◆ getDistanceAlongPolyline()

double GeoLib::Polyline::getDistanceAlongPolyline ( const MathLib::Point3d pnt,
const double  epsilon_radius 
) const

returns the distance along the polyline from the beginning of the polyline

Parameters
pntthe point on the polyline
epsilon_radiusthe epsilon
Returns
the distance along the polyline between the given point and the beginning of the polyline. If the given point is not on the polyine, negative value is returned.

Definition at line 363 of file Polyline.cpp.

365 {
366  double dist(-1.0);
367  double lambda;
368  bool found = false;
369  double act_length_of_ply = 0.0;
370  // loop over all line segments of the polyline
371  for (std::size_t k = 0; k < getNumberOfSegments(); k++)
372  {
373  auto const a =
374  Eigen::Map<Eigen::Vector3d const>(getPoint(k)->getCoords());
375  auto const b =
376  Eigen::Map<Eigen::Vector3d const>(getPoint(k + 1)->getCoords());
377  double const seg_length((b - a).norm());
378  act_length_of_ply += seg_length;
379  // is the orthogonal projection of the j-th node to the
380  // line g(lambda) = _ply->getPoint(k) + lambda * (_ply->getPoint(k+1) -
381  // _ply->getPoint(k)) at the k-th line segment of the polyline, i.e. 0
382  // <= lambda <= 1?
384  *getPoint(k + 1), lambda,
385  dist) <= epsilon_radius)
386  {
387  double const lower_lambda(-epsilon_radius / seg_length);
388  double const upper_lambda(1 + epsilon_radius / seg_length);
389 
390  if (lower_lambda <= lambda && lambda <= upper_lambda)
391  {
392  found = true;
393  dist = act_length_of_ply + dist;
394  break;
395  } // end if lambda
396  }
397  } // end line segment loop
398 
399  if (!found)
400  {
401  dist = -1.0;
402  }
403  return dist;
404 }
const Point * getPoint(std::size_t i) const
returns the i-th point contained in the polyline
Definition: Polyline.cpp:186
double norm(MatrixOrVector const &x, MathLib::VecNormType type)
Definition: LinAlg.h:88
double calcProjPntToLineAndDists(Point3d const &pp, Point3d const &pa, Point3d const &pb, double &lambda, double &d0)
Definition: MathTools.cpp:19

References MathLib::calcProjPntToLineAndDists(), getNumberOfSegments(), getPoint(), and MathLib::LinAlg::norm().

Referenced by MeshGeoToolsLib::MeshNodesAlongPolyline::MeshNodesAlongPolyline().

◆ getGeoType()

GEOTYPE GeoLib::Polyline::getGeoType ( ) const
inlineoverridevirtual

return a geometry type

Implements GeoLib::GeoObject.

Definition at line 109 of file Polyline.h.

109 { return GEOTYPE::POLYLINE; }

References GeoLib::POLYLINE.

◆ getLocationOfPoint()

Location GeoLib::Polyline::getLocationOfPoint ( std::size_t  k,
MathLib::Point3d const &  pnt 
) const
protected

2D method - ignores z coordinate. It calculates the location of the point relative to the k-th line segment of the polyline. (literature reference: Computational Geometry and Computer Graphics in C++; Michael J. Laszlo)

Parameters
kthe number of line segment
pntthe point
Returns
a value of enum LOCATION

Definition at line 322 of file Polyline.cpp.

324 {
325  assert(k < _ply_pnt_ids.size() - 1);
326 
327  GeoLib::Point const& source(*(_ply_pnts[_ply_pnt_ids[k]]));
328  GeoLib::Point const& dest(*(_ply_pnts[_ply_pnt_ids[k + 1]]));
329  long double a[2] = {dest[0] - source[0], dest[1] - source[1]}; // vector
330  long double b[2] = {pnt[0] - source[0], pnt[1] - source[1]}; // vector
331 
332  long double det_2x2(a[0] * b[1] - a[1] * b[0]);
333 
334  if (det_2x2 > std::numeric_limits<double>::epsilon())
335  {
336  return Location::LEFT;
337  }
338  if (std::numeric_limits<double>::epsilon() < std::abs(det_2x2))
339  {
340  return Location::RIGHT;
341  }
342  if (a[0] * b[0] < 0.0 || a[1] * b[1] < 0.0)
343  {
344  return Location::BEHIND;
345  }
346  if (a[0] * a[0] + a[1] * a[1] < b[0] * b[0] + b[1] * b[1])
347  {
348  return Location::BEYOND;
349  }
350  if (MathLib::sqrDist(pnt, *_ply_pnts[_ply_pnt_ids[k]]) <
351  pow(std::numeric_limits<double>::epsilon(), 2))
352  {
353  return Location::SOURCE;
354  }
355  if (MathLib::sqrDist(pnt, *_ply_pnts[_ply_pnt_ids[k + 1]]) <
356  std::sqrt(std::numeric_limits<double>::epsilon()))
357  {
358  return Location::DESTINATION;
359  }
360  return Location::BETWEEN;
361 }
double sqrDist(MathLib::Point3d const &p0, MathLib::Point3d const &p1)
Definition: Point3d.h:48

References _ply_pnt_ids, _ply_pnts, GeoLib::BEHIND, GeoLib::BETWEEN, GeoLib::BEYOND, GeoLib::DESTINATION, GeoLib::LEFT, GeoLib::RIGHT, GeoLib::SOURCE, and MathLib::sqrDist().

Referenced by GeoLib::Polygon::getEdgeType().

◆ getNumberOfPoints()

◆ getNumberOfSegments()

std::size_t GeoLib::Polyline::getNumberOfSegments ( ) const

Definition at line 114 of file Polyline.cpp.

115 {
116  return _ply_pnt_ids.empty() ? 0 : _ply_pnt_ids.size() - 1;
117 }

References _ply_pnt_ids.

Referenced by end(), getDistanceAlongPolyline(), getSegment(), and GeoLib::lineSegmentsIntersect().

◆ getPoint()

◆ getPointID()

◆ getPointsVec()

std::vector< Point * > const & GeoLib::Polyline::getPointsVec ( ) const

◆ getSegment() [1/2]

LineSegment GeoLib::Polyline::getSegment ( std::size_t  i)
private

Definition at line 173 of file Polyline.cpp.

174 {
175  assert(i < getNumberOfSegments());
176  return LineSegment(_ply_pnts[_ply_pnt_ids[i]],
177  _ply_pnts[_ply_pnt_ids[i + 1]], false);
178 }

References _ply_pnt_ids, _ply_pnts, and getNumberOfSegments().

◆ getSegment() [2/2]

LineSegment GeoLib::Polyline::getSegment ( std::size_t  i) const
private

Definition at line 166 of file Polyline.cpp.

167 {
168  assert(i < getNumberOfSegments());
169  return LineSegment(_ply_pnts[_ply_pnt_ids[i]],
170  _ply_pnts[_ply_pnt_ids[i + 1]], false);
171 }

References _ply_pnt_ids, _ply_pnts, and getNumberOfSegments().

◆ insertPoint()

bool GeoLib::Polyline::insertPoint ( std::size_t  pos,
std::size_t  pnt_id 
)
virtual

Method inserts a new point (that have to be inside the _ply_pnts vector) at the given position in the polyline if and only if the resulting segments won't be empty.

Parameters
posthe position in the polyline, pos have to be a value into the interval [0, number of points)
pnt_idthe id of the new point in the vector of points the polyline is based on
Returns
true if the point could be inserted, else false (if empty line segments would be created).

Reimplemented in GeoLib::PolylineWithSegmentMarker, and GeoLib::PolygonWithSegmentMarker.

Definition at line 55 of file Polyline.cpp.

56 {
57  if (pnt_id >= _ply_pnts.size())
58  {
59  return false;
60  }
61  if (pos > _ply_pnt_ids.size())
62  {
63  return false;
64  }
65 
66  if (pos == _ply_pnt_ids.size())
67  {
68  return addPoint(pnt_id);
69  }
70 
71  // check if inserting pnt_id would result in two identical IDs for adjacent
72  // points
73  if (pos == 0 && pnt_id == _ply_pnt_ids[0])
74  {
75  return false;
76  }
77  if (pos != 0)
78  {
79  if (pos == (_ply_pnt_ids.size() - 1) && pnt_id == _ply_pnt_ids[pos])
80  {
81  return false;
82  }
83  if (pnt_id == _ply_pnt_ids[pos - 1] || pnt_id == _ply_pnt_ids[pos])
84  {
85  return false;
86  }
87  }
88 
89  auto const pos_dt(
90  static_cast<std::vector<std::size_t>::difference_type>(pos));
91  auto it(_ply_pnt_ids.begin() + pos_dt);
92  _ply_pnt_ids.insert(it, pnt_id);
93 
94  return true;
95 }

References _ply_pnt_ids, _ply_pnts, and addPoint().

Referenced by GeoLib::computeAndInsertAllIntersectionPoints(), GeoLib::PolygonWithSegmentMarker::insertPoint(), GeoLib::PolylineWithSegmentMarker::insertPoint(), FileIO::GMSH::GMSHPolygonTree::insertPolyline(), and MeshGeoToolsLib::insertSubSegments().

◆ isClosed()

bool GeoLib::Polyline::isClosed ( ) const

returns true if the polyline is closed

Definition at line 119 of file Polyline.cpp.

120 {
121  if (_ply_pnt_ids.size() < 3)
122  {
123  return false;
124  }
125 
126  return _ply_pnt_ids.front() == _ply_pnt_ids.back();
127 }

References _ply_pnt_ids.

Referenced by closePolyline(), FileIO::createSurface(), FileIO::createSurfaceWithEarClipping(), GeoLib::Polygon::initialise(), main(), MeshGeoToolsLib::BoundaryElementsAlongPolyline::modifyEdgeNodeOrdering(), and FileIO::FEFLOWMeshInterface::setMaterialIDs().

◆ isCoplanar()

bool GeoLib::Polyline::isCoplanar ( ) const

returns true if the polyline is coplanar

Definition at line 129 of file Polyline.cpp.

130 {
131  std::size_t const n_points(_ply_pnt_ids.size());
132  if (n_points < 4)
133  {
134  return true;
135  }
136 
137  GeoLib::Point const& p0(*this->getPoint(0));
138  GeoLib::Point const& p1(*this->getPoint(1));
139  GeoLib::Point const& p2(*this->getPoint(2));
140  for (std::size_t i = 3; i < n_points; ++i)
141  {
142  if (!MathLib::isCoplanar(p0, p1, p2, *this->getPoint(i)))
143  {
144  DBUG(
145  "Point {:d} is not coplanar to the first three points of the "
146  "line.",
147  i);
148  return false;
149  }
150  }
151  return true;
152 }
void DBUG(char const *fmt, Args const &... args)
Definition: Logging.h:27
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.

References _ply_pnt_ids, DBUG(), getPoint(), and MathLib::isCoplanar().

◆ isPointIDInPolyline()

bool GeoLib::Polyline::isPointIDInPolyline ( std::size_t  pnt_id) const

Method tests if the given id of a point (within the vector of points the polyline is based on) is inside the polyline.

Parameters
pnt_idthe id of the point
Returns
true if the point is part of the polyline, else false

Definition at line 154 of file Polyline.cpp.

155 {
156  return std::find(_ply_pnt_ids.begin(), _ply_pnt_ids.end(), pnt_id) !=
157  _ply_pnt_ids.end();
158 }

References _ply_pnt_ids.

Referenced by FileIO::GMSH::GMSHPolygonTree::insertPolyline().

◆ operator=()

Polyline& GeoLib::Polyline::operator= ( Polyline const &  other)
delete

◆ removePoint()

void GeoLib::Polyline::removePoint ( std::size_t  pos)

Method removes a point from the polyline. The connecting line segments will be removed and the length of the polyline will be changed.

Parameters
posa valid position within the polyline

Definition at line 97 of file Polyline.cpp.

98 {
99  if (pos >= _ply_pnt_ids.size())
100  {
101  return;
102  }
103 
104  auto const pos_dt(
105  static_cast<std::vector<std::size_t>::difference_type>(pos));
106  _ply_pnt_ids.erase(_ply_pnt_ids.begin() + pos_dt);
107 }

References _ply_pnt_ids.

◆ setPointID()

void GeoLib::Polyline::setPointID ( std::size_t  idx,
std::size_t  id 
)

Changes a point index for one point in a line

Parameters
idxIndex of point in line
idID of point in PointVec object

Definition at line 180 of file Polyline.cpp.

181 {
182  assert(idx < _ply_pnt_ids.size());
183  _ply_pnt_ids[idx] = id;
184 }

References _ply_pnt_ids.

Friends And Related Function Documentation

◆ Polygon

friend class Polygon
friend

Definition at line 94 of file Polyline.h.

Referenced by GeoLib::Polygon::Polygon(), and GeoLib::Polygon::splitPolygonAtPoint().

Member Data Documentation

◆ _ply_pnt_ids

std::vector<std::size_t> GeoLib::Polyline::_ply_pnt_ids
protected

◆ _ply_pnts

const std::vector<Point*>& GeoLib::Polyline::_ply_pnts
protected

a reference to the vector of pointers to the geometric points

Definition at line 227 of file Polyline.h.

Referenced by addPoint(), getLocationOfPoint(), getPoint(), getPointsVec(), getSegment(), insertPoint(), and GeoLib::Polygon::splitPolygonAtIntersection().


The documentation for this class was generated from the following files: