OGS
EarClippingTriangulation.cpp
Go to the documentation of this file.
1
16
17#include <numeric>
18
19#include "BaseLib/Algorithm.h"
21#include "Point.h"
22#include "Polygon.h"
23#include "Triangle.h"
24
25namespace GeoLib
26{
28 GeoLib::Polygon const& polygon, std::list<GeoLib::Triangle>& triangles,
29 bool rot)
30{
31 copyPolygonPoints(polygon);
32
33 if (rot)
34 {
37 }
38
40 initLists();
41 clipEars();
43
44 std::vector<GeoLib::Point*> const& ref_pnts_vec(polygon.getPointsVec());
46 {
47 for (auto const& t : _triangles)
48 {
49 const std::size_t i0(polygon.getPointID(t[0]));
50 const std::size_t i1(polygon.getPointID(t[1]));
51 const std::size_t i2(polygon.getPointID(t[2]));
52 triangles.emplace_back(ref_pnts_vec, i0, i1, i2);
53 }
54 }
55 else
56 {
57 std::size_t n_pnts(_pnts.size() - 1);
58 for (auto const& t : _triangles)
59 {
60 const std::size_t i0(polygon.getPointID(n_pnts - t[0]));
61 const std::size_t i1(polygon.getPointID(n_pnts - t[1]));
62 const std::size_t i2(polygon.getPointID(n_pnts - t[2]));
63 triangles.emplace_back(ref_pnts_vec, i0, i1, i2);
64 }
65 }
66}
67
69{
70 for (auto p : _pnts)
71 {
72 delete p;
73 }
74}
75
77{
78 // copy points - last point is identical to the first
79 std::size_t n_pnts(polygon.getNumberOfPoints() - 1);
80 for (std::size_t k(0); k < n_pnts; k++)
81 {
82 _pnts.push_back(new GeoLib::Point(*(polygon.getPoint(k))));
83 }
84}
85
87{
88 std::size_t n_pnts(_pnts.size());
89 // get the left most upper point
90 std::size_t min_x_max_y_idx(0); // for orientation check
91 for (std::size_t k(0); k < n_pnts; k++)
92 {
93 if ((*(_pnts[k]))[0] <= (*(_pnts[min_x_max_y_idx]))[0])
94 {
95 if ((*(_pnts[k]))[0] < (*(_pnts[min_x_max_y_idx]))[0])
96 {
97 min_x_max_y_idx = k;
98 }
99 else
100 {
101 if ((*(_pnts[k]))[1] > (*(_pnts[min_x_max_y_idx]))[1])
102 {
103 min_x_max_y_idx = k;
104 }
105 }
106 }
107 }
108 // determine orientation
109 if (0 < min_x_max_y_idx && min_x_max_y_idx < n_pnts - 1)
110 {
112 GeoLib::getOrientation(*_pnts[min_x_max_y_idx - 1],
113 *_pnts[min_x_max_y_idx],
114 *_pnts[min_x_max_y_idx + 1]);
115 }
116 else
117 {
118 if (0 == min_x_max_y_idx)
119 {
121 *_pnts[n_pnts - 1], *_pnts[0], *_pnts[1]);
122 }
123 else
124 {
126 *_pnts[n_pnts - 2], *_pnts[n_pnts - 1], *_pnts[0]);
127 }
128 }
130 {
131 // switch orientation
132 for (std::size_t k(0); k < n_pnts / 2; k++)
133 {
134 std::swap(_pnts[k], _pnts[n_pnts - 1 - k]);
135 }
136 }
137}
138
139bool EarClippingTriangulation::isEar(std::size_t v0, std::size_t v1,
140 std::size_t v2) const
141{
142 for (auto const& v : _vertex_list)
143 {
144 if (v != v0 && v != v1 && v != v2)
145 {
146 if (MathLib::isPointInTriangle(*_pnts[v], *_pnts[v0], *_pnts[v1],
147 *_pnts[v2]))
148 {
149 return false;
150 }
151 }
152 }
153 return true;
154}
155
157{
158 _vertex_list.resize(_pnts.size());
159 std::iota(begin(_vertex_list), end(_vertex_list), 0);
160}
161
163{
164 // go through points checking ccw, cw or collinear order and identifying
165 // ears
166 auto it = _vertex_list.begin();
167 auto prev = _vertex_list.end();
168 --prev;
169 auto next = it;
170 ++next;
171 GeoLib::Orientation orientation;
172 bool first_run(
173 true); // saves special handling of the last case with identical code
174 while (_vertex_list.size() >= 3 && first_run)
175 {
176 if (next == _vertex_list.end())
177 {
178 first_run = false;
179 next = _vertex_list.begin();
180 }
181 orientation = getOrientation(*_pnts[*prev], *_pnts[*it], *_pnts[*next]);
182 if (orientation == GeoLib::COLLINEAR)
183 {
184 WARN(
185 "EarClippingTriangulation::initLists(): collinear points "
186 "({:f}, {:f}, {:f}), ({:f}, {:f}, {:f}), ({:f}, {:f}, {:f})",
187 (*_pnts[*prev])[0], (*_pnts[*prev])[1], (*_pnts[*prev])[2],
188 (*_pnts[*it])[0], (*_pnts[*it])[1], (*_pnts[*it])[2],
189 (*_pnts[*next])[0], (*_pnts[*next])[1], (*_pnts[*next])[2]);
190 it = _vertex_list.erase(it);
191 ++next;
192 }
193 else
194 {
195 if (orientation == GeoLib::CW)
196 {
197 _convex_vertex_list.push_back(*it);
198 if (isEar(*prev, *it, *next))
199 _ear_list.push_back(*it);
200 }
201 prev = it;
202 it = next;
203 ++next;
204 }
205 }
206}
207
209{
210 // *** clip an ear
211 while (_vertex_list.size() > 3)
212 {
213 // pop ear from list
214 std::size_t ear = _ear_list.front();
215 _ear_list.pop_front();
216 // remove ear tip from _convex_vertex_list
217 _convex_vertex_list.remove(ear);
218
219 // remove ear from vertex_list, apply changes to _ear_list,
220 // _convex_vertex_list
221 bool nfound(true);
222 auto it = _vertex_list.begin();
223 auto prev = _vertex_list.end();
224 --prev;
225 while (nfound && it != _vertex_list.end())
226 {
227 if (*it == ear)
228 {
229 nfound = false;
230 it = _vertex_list.erase(it); // remove ear tip
231 auto next = it;
232 if (next == _vertex_list.end())
233 {
234 next = _vertex_list.begin();
235 prev = _vertex_list.end();
236 --prev;
237 }
238 // add triangle
239 _triangles.emplace_back(_pnts, *prev, *next, ear);
240
241 // check the orientation of prevprev, prev, next
242 auto prevprev = prev;
243 if (prev == _vertex_list.begin())
244 {
245 prevprev = _vertex_list.end();
246 }
247 --prevprev;
248
249 // apply changes to _convex_vertex_list and _ear_list looking
250 // "backward"
252 *_pnts[*prevprev], *_pnts[*prev], *_pnts[*next]);
253 if (orientation == GeoLib::CW)
254 {
256 // prev is convex
257 if (isEar(*prevprev, *prev, *next))
258 {
259 // prev is an ear tip
261 }
262 else
263 {
264 // if necessary remove prev
265 _ear_list.remove(*prev);
266 }
267 }
268 else
269 {
270 // prev is not convex => reflex or collinear
271 _convex_vertex_list.remove(*prev);
272 _ear_list.remove(*prev);
273 if (orientation == GeoLib::COLLINEAR)
274 {
275 prev = _vertex_list.erase(prev);
276 if (prev == _vertex_list.begin())
277 {
278 prev = _vertex_list.end();
279 --prev;
280 }
281 else
282 {
283 --prev;
284 }
285 }
286 }
287
288 // check the orientation of prev, next, nextnext
289 auto nextnext = next;
290 ++nextnext;
291 auto help_it = _vertex_list.end();
292 --help_it;
293 if (next == help_it)
294 {
295 nextnext = _vertex_list.begin();
296 }
297
298 // apply changes to _convex_vertex_list and _ear_list looking
299 // "forward"
300 orientation = getOrientation(*_pnts[*prev], *_pnts[*next],
301 *_pnts[*nextnext]);
302 if (orientation == GeoLib::CW)
303 {
305 // next is convex
306 if (isEar(*prev, *next, *nextnext))
307 {
308 // next is an ear tip
310 }
311 else
312 {
313 // if necessary remove *next
314 _ear_list.remove(*next);
315 }
316 }
317 else
318 {
319 // next is not convex => reflex or collinear
320 _convex_vertex_list.remove(*next);
321 _ear_list.remove(*next);
322 if (orientation == GeoLib::COLLINEAR)
323 {
324 _vertex_list.erase(next);
325 }
326 }
327 }
328 else
329 {
330 prev = it;
331 ++it;
332 }
333 }
334 }
335}
336
338{
339 auto next = _vertex_list.begin();
340 auto prev = next;
341 ++next;
342 if (next == _vertex_list.end())
343 {
344 return;
345 }
346 auto it = next;
347 ++next;
348 if (next == _vertex_list.end())
349 {
350 return;
351 }
352
353 if (getOrientation(*_pnts[*prev], *_pnts[*it], *_pnts[*next]) ==
355 {
356 _triangles.emplace_back(_pnts, *prev, *it, *next);
357 }
358 else
359 {
360 _triangles.emplace_back(_pnts, *prev, *next, *it);
361 }
362}
363} // end namespace GeoLib
Definition of the EarClippingTriangulation class.
Definition of the Point class.
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:40
Definition of the Polygon class.
bool isEar(std::size_t v0, std::size_t v1, std::size_t v2) const
std::vector< GeoLib::Point * > _pnts
EarClippingTriangulation(GeoLib::Polygon const &polygon, std::list< GeoLib::Triangle > &triangles, bool rot=true)
std::list< std::size_t > _convex_vertex_list
void copyPolygonPoints(GeoLib::Polygon const &polygon)
std::list< GeoLib::Triangle > _triangles
std::size_t getPointID(std::size_t const i) const
Definition Polyline.cpp:160
std::size_t getNumberOfPoints() const
Definition Polyline.cpp:109
const Point * getPoint(std::size_t i) const
returns the i-th point contained in the polyline
Definition Polyline.cpp:179
std::vector< Point * > const & getPointsVec() const
Definition Polyline.cpp:185
void uniquePushBack(Container &container, typename Container::value_type const &element)
Definition Algorithm.h:212
Eigen::Matrix3d rotatePointsToXY(InputIterator1 p_pnts_begin, InputIterator1 p_pnts_end, InputIterator2 r_pnts_begin, InputIterator2 r_pnts_end)
Orientation getOrientation(MathLib::Point3d const &p0, MathLib::Point3d const &p1, MathLib::Point3d const &p2)
bool isPointInTriangle(MathLib::Point3d const &p, MathLib::Point3d const &a, MathLib::Point3d const &b, MathLib::Point3d const &c, double eps_pnt_out_of_plane, double eps_pnt_out_of_tri, MathLib::TriangleTest algorithm)