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 50 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)
 
virtual 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 27 of file Polyline.cpp.

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

Referenced by constructPolylineFromSegments().

◆ Polyline() [2/2]

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

Copy constructor

Parameters
plyPolyline

Definition at line 29 of file Polyline.cpp.

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

◆ ~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 34 of file Polyline.cpp.

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

References _ply_pnt_ids, and _ply_pnts.

Referenced by GeoLib::PolygonWithSegmentMarker::addPoint(), GeoLib::PolylineWithSegmentMarker::addPoint(), closePolyline(), createPolyline(), insertPoint(), mergeGeometries(), FileIO::SwmmInterface::readLinksAsPolylines(), FileIO::Legacy::readPolylinePointVector(), GeoLib::rotatePolygonToXY(), 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 299 of file Polyline.cpp.

300 {
301  if (getNumberOfPoints() < 2)
302  {
303  ERR("Polyline::closePolyline(): Input polyline needs to be composed of "
304  "at least three points.");
305  }
306  if (!isClosed())
307  {
308  addPoint(getPointID(0));
309  }
310 }
void ERR(char const *fmt, Args const &... args)
Definition: Logging.h:42
std::size_t getPointID(std::size_t i) const
Definition: Polyline.cpp:150
std::size_t getNumberOfPoints() const
Definition: Polyline.cpp:99
bool isClosed() const
Definition: Polyline.cpp:109
virtual bool addPoint(std::size_t pnt_id)
Definition: Polyline.cpp:34

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

Referenced by GEOModels::connectPolylineSegments(), and main().

◆ 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 187 of file Polyline.cpp.

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

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 353 of file Polyline.cpp.

355 {
356  double dist(-1.0);
357  double lambda;
358  bool found = false;
359  double act_length_of_ply = 0.0;
360  // loop over all line segments of the polyline
361  for (std::size_t k = 0; k < getNumberOfSegments(); k++)
362  {
363  auto const a =
364  Eigen::Map<Eigen::Vector3d const>(getPoint(k)->getCoords());
365  auto const b =
366  Eigen::Map<Eigen::Vector3d const>(getPoint(k + 1)->getCoords());
367  double const seg_length((b - a).norm());
368  act_length_of_ply += seg_length;
369  // is the orthogonal projection of the j-th node to the
370  // line g(lambda) = _ply->getPoint(k) + lambda * (_ply->getPoint(k+1) -
371  // _ply->getPoint(k)) at the k-th line segment of the polyline, i.e. 0
372  // <= lambda <= 1?
374  *getPoint(k + 1), lambda,
375  dist) <= epsilon_radius)
376  {
377  double const lower_lambda(-epsilon_radius / seg_length);
378  double const upper_lambda(1 + epsilon_radius / seg_length);
379 
380  if (lower_lambda <= lambda && lambda <= upper_lambda)
381  {
382  found = true;
383  dist = act_length_of_ply + dist;
384  break;
385  } // end if lambda
386  }
387  } // end line segment loop
388 
389  if (!found)
390  {
391  dist = -1.0;
392  }
393  return dist;
394 }
const Point * getPoint(std::size_t i) const
returns the i-th point contained in the polyline
Definition: Polyline.cpp:176
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 107 of file Polyline.h.

107 { 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 312 of file Polyline.cpp.

314 {
315  assert(k < _ply_pnt_ids.size() - 1);
316 
317  GeoLib::Point const& source(*(_ply_pnts[_ply_pnt_ids[k]]));
318  GeoLib::Point const& dest(*(_ply_pnts[_ply_pnt_ids[k + 1]]));
319  long double a[2] = {dest[0] - source[0], dest[1] - source[1]}; // vector
320  long double b[2] = {pnt[0] - source[0], pnt[1] - source[1]}; // vector
321 
322  long double det_2x2(a[0] * b[1] - a[1] * b[0]);
323 
324  if (det_2x2 > std::numeric_limits<double>::epsilon())
325  {
326  return Location::LEFT;
327  }
328  if (std::numeric_limits<double>::epsilon() < std::abs(det_2x2))
329  {
330  return Location::RIGHT;
331  }
332  if (a[0] * b[0] < 0.0 || a[1] * b[1] < 0.0)
333  {
334  return Location::BEHIND;
335  }
336  if (a[0] * a[0] + a[1] * a[1] < b[0] * b[0] + b[1] * b[1])
337  {
338  return Location::BEYOND;
339  }
340  if (MathLib::sqrDist(pnt, *_ply_pnts[_ply_pnt_ids[k]]) <
341  pow(std::numeric_limits<double>::epsilon(), 2))
342  {
343  return Location::SOURCE;
344  }
345  if (MathLib::sqrDist(pnt, *_ply_pnts[_ply_pnt_ids[k + 1]]) <
346  std::sqrt(std::numeric_limits<double>::epsilon()))
347  {
348  return Location::DESTINATION;
349  }
350  return Location::BETWEEN;
351 }
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 104 of file Polyline.cpp.

105 {
106  return _ply_pnt_ids.empty() ? 0 : _ply_pnt_ids.size() - 1;
107 }

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 163 of file Polyline.cpp.

164 {
165  assert(i < getNumberOfSegments());
166  return LineSegment(_ply_pnts[_ply_pnt_ids[i]],
167  _ply_pnts[_ply_pnt_ids[i + 1]], false);
168 }

References _ply_pnt_ids, _ply_pnts, and getNumberOfSegments().

◆ getSegment() [2/2]

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

Definition at line 156 of file Polyline.cpp.

157 {
158  assert(i < getNumberOfSegments());
159  return LineSegment(_ply_pnts[_ply_pnt_ids[i]],
160  _ply_pnts[_ply_pnt_ids[i + 1]], false);
161 }

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 51 of file Polyline.cpp.

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

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 109 of file Polyline.cpp.

110 {
111  if (_ply_pnt_ids.size() < 3)
112  {
113  return false;
114  }
115 
116  return _ply_pnt_ids.front() == _ply_pnt_ids.back();
117 }

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 119 of file Polyline.cpp.

120 {
121  std::size_t const n_points(_ply_pnt_ids.size());
122  if (n_points < 4)
123  {
124  return true;
125  }
126 
127  GeoLib::Point const& p0(*this->getPoint(0));
128  GeoLib::Point const& p1(*this->getPoint(1));
129  GeoLib::Point const& p2(*this->getPoint(2));
130  for (std::size_t i = 3; i < n_points; ++i)
131  {
132  if (!MathLib::isCoplanar(p0, p1, p2, *this->getPoint(i)))
133  {
134  DBUG(
135  "Point {:d} is not coplanar to the first three points of the "
136  "line.",
137  i);
138  return false;
139  }
140  }
141  return true;
142 }
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().

Referenced by main().

◆ 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 144 of file Polyline.cpp.

145 {
146  return std::find(_ply_pnt_ids.begin(), _ply_pnt_ids.end(), pnt_id) !=
147  _ply_pnt_ids.end();
148 }

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)
virtual

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 87 of file Polyline.cpp.

88 {
89  if (pos >= _ply_pnt_ids.size())
90  {
91  return;
92  }
93 
94  auto const pos_dt(
95  static_cast<std::vector<std::size_t>::difference_type>(pos));
96  _ply_pnt_ids.erase(_ply_pnt_ids.begin() + pos_dt);
97 }

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 170 of file Polyline.cpp.

171 {
172  assert(idx < _ply_pnt_ids.size());
173  _ply_pnt_ids[idx] = id;
174 }

References _ply_pnt_ids.

Friends And Related Function Documentation

◆ Polygon

friend class Polygon
friend

Definition at line 92 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 225 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: