CDT  1.4.1
C++ library for constrained Delaunay triangulation
Loading...
Searching...
No Matches
CDTUtils.hpp
Go to the documentation of this file.
1/* This Source Code Form is subject to the terms of the Mozilla Public
2 * License, v. 2.0. If a copy of the MPL was not distributed with this
3 * file, You can obtain one at https://mozilla.org/MPL/2.0/. */
4
10#include "CDTUtils.h"
11
12#include "predicates.h" // robust predicates: orient, in-circle
13
14#include <stdexcept>
15
16namespace CDT
17{
18
19//*****************************************************************************
20// V2d
21//*****************************************************************************
22template <typename T>
23V2d<T> V2d<T>::make(const T x, const T y)
24{
25 V2d<T> out = {x, y};
26 return out;
27}
28
29//*****************************************************************************
30// Box2d
31//*****************************************************************************
32template <typename T>
33Box2d<T> envelopBox(const std::vector<V2d<T> >& vertices)
34{
35 return envelopBox<T>(
36 vertices.begin(), vertices.end(), getX_V2d<T>, getY_V2d<T>);
37}
38
39//*****************************************************************************
40// Edge
41//*****************************************************************************
43 : m_vertices(
44 iV1 < iV2 ? std::make_pair(iV1, iV2) : std::make_pair(iV2, iV1))
45{}
46
48{
49 return m_vertices == other.m_vertices;
50}
51
53{
54 return !(this->operator==(other));
55}
56
58{
59 return m_vertices.first;
60}
61
63{
64 return m_vertices.second;
65}
66
67CDT_INLINE_IF_HEADER_ONLY const std::pair<VertInd, VertInd>& Edge::verts() const
68{
69 return m_vertices;
70}
71
72//*****************************************************************************
73// Utility functions
74//*****************************************************************************
76{
77 return Index((i + 1) % 3);
78}
79
81{
82 return Index((i + 2) % 3);
83}
84
86{
87 return location == PtTriLocation::OnEdge1 ||
88 location == PtTriLocation::OnEdge2 ||
89 location == PtTriLocation::OnEdge3;
90}
91
93{
94 assert(isOnEdge(location));
95 return static_cast<Index>(location - PtTriLocation::OnEdge1);
96}
97
98template <typename T>
99T orient2D(const V2d<T>& p, const V2d<T>& v1, const V2d<T>& v2)
100{
101 return predicates::adaptive::orient2d(v1.x, v1.y, v2.x, v2.y, p.x, p.y);
102}
103
104template <typename T>
106 const V2d<T>& p,
107 const V2d<T>& v1,
108 const V2d<T>& v2,
109 const T orientationTolerance)
110{
111 return classifyOrientation(orient2D(p, v1, v2), orientationTolerance);
112}
113
114template <typename T>
116classifyOrientation(const T orientation, const T orientationTolerance)
117{
118 if(orientation < -orientationTolerance)
119 return PtLineLocation::Right;
120 if(orientation > orientationTolerance)
121 return PtLineLocation::Left;
122 return PtLineLocation::OnLine;
123}
124
125template <typename T>
127 const V2d<T>& p,
128 const V2d<T>& v1,
129 const V2d<T>& v2,
130 const V2d<T>& v3)
131{
132 using namespace predicates::adaptive;
133 PtTriLocation::Enum result = PtTriLocation::Inside;
134 PtLineLocation::Enum edgeCheck = locatePointLine(p, v1, v2);
135 if(edgeCheck == PtLineLocation::Right)
136 return PtTriLocation::Outside;
137 if(edgeCheck == PtLineLocation::OnLine)
138 result = PtTriLocation::OnEdge1;
139 edgeCheck = locatePointLine(p, v2, v3);
140 if(edgeCheck == PtLineLocation::Right)
141 return PtTriLocation::Outside;
142 if(edgeCheck == PtLineLocation::OnLine)
143 {
144 result = (result == PtTriLocation::Inside) ? PtTriLocation::OnEdge2
145 : PtTriLocation::OnVertex;
146 }
147 edgeCheck = locatePointLine(p, v3, v1);
148 if(edgeCheck == PtLineLocation::Right)
149 return PtTriLocation::Outside;
150 if(edgeCheck == PtLineLocation::OnLine)
151 {
152 result = (result == PtTriLocation::Inside) ? PtTriLocation::OnEdge3
153 : PtTriLocation::OnVertex;
154 }
155 return result;
156}
157
159{
160 if(vertIndex == Index(0))
161 return Index(1);
162 if(vertIndex == Index(1))
163 return Index(2);
164 if(vertIndex == Index(2))
165 return Index(0);
166 assert(false && "Invalid vertex index");
167 throw std::runtime_error("Invalid vertex index");
168}
169
171{
172 if(neighborIndex == Index(0))
173 return Index(2);
174 if(neighborIndex == Index(1))
175 return Index(0);
176 if(neighborIndex == Index(2))
177 return Index(1);
178 assert(false && "Invalid neighbor index");
179 throw std::runtime_error("Invalid neighbor index");
180}
181
184{
185 assert(vv[0] == iVert || vv[1] == iVert || vv[2] == iVert);
186 if(vv[0] == iVert)
187 return Index(1);
188 if(vv[1] == iVert)
189 return Index(2);
190 return Index(0);
191}
192
194 const VerticesArr3& vv,
195 const VertInd iVedge1,
196 const VertInd iVedge2)
197{
198 assert(vv[0] == iVedge1 || vv[1] == iVedge1 || vv[2] == iVedge1);
199 assert(vv[0] == iVedge2 || vv[1] == iVedge2 || vv[2] == iVedge2);
200 assert(
201 (vv[0] != iVedge1 && vv[0] != iVedge2) ||
202 (vv[1] != iVedge1 && vv[1] != iVedge2) ||
203 (vv[2] != iVedge1 && vv[2] != iVedge2));
204 /*
205 * vv[2]
206 * /\
207 * n[2]/ \n[1]
208 * /____\
209 * vv[0] n[0] vv[1]
210 */
211 if(vv[0] == iVedge1)
212 {
213 if(vv[1] == iVedge2)
214 return Index(0);
215 return Index(2);
216 }
217 if(vv[0] == iVedge2)
218 {
219 if(vv[1] == iVedge1)
220 return Index(0);
221 return Index(2);
222 }
223 return Index(1);
224}
225
227opposedVertexInd(const NeighborsArr3& nn, const TriInd iTopo)
228{
229 assert(nn[0] == iTopo || nn[1] == iTopo || nn[2] == iTopo);
230 if(nn[0] == iTopo)
231 return Index(2);
232 if(nn[1] == iTopo)
233 return Index(0);
234 return Index(1);
235}
236
238vertexInd(const VerticesArr3& vv, const VertInd iV)
239{
240 assert(vv[0] == iV || vv[1] == iV || vv[2] == iV);
241 if(vv[0] == iV)
242 return Index(0);
243 if(vv[1] == iV)
244 return Index(1);
245 return Index(2);
246}
247
249opposedTriangle(const Triangle& tri, const VertInd iVert)
250{
251 return tri.neighbors[opposedTriangleInd(tri.vertices, iVert)];
252}
253
255opposedVertex(const Triangle& tri, const TriInd iTopo)
256{
257 return tri.vertices[opposedVertexInd(tri.neighbors, iTopo)];
258}
259
262edgeNeighbor(const Triangle& tri, VertInd iVedge1, VertInd iVedge2)
263{
264 return tri.neighbors[edgeNeighborInd(tri.vertices, iVedge1, iVedge2)];
265}
266
267template <typename T>
269 const V2d<T>& p,
270 const V2d<T>& v1,
271 const V2d<T>& v2,
272 const V2d<T>& v3)
273{
274 using namespace predicates::adaptive;
275 return incircle(v1.x, v1.y, v2.x, v2.y, v3.x, v3.y, p.x, p.y) > T(0);
276}
277
279bool verticesShareEdge(const TriIndVec& aTris, const TriIndVec& bTris)
280{
281 for(TriIndVec::const_iterator it = aTris.begin(); it != aTris.end(); ++it)
282 if(std::find(bTris.begin(), bTris.end(), *it) != bTris.end())
283 return true;
284 return false;
285}
286
287template <typename T>
288T distanceSquared(const T ax, const T ay, const T bx, const T by)
289{
290 const T dx = bx - ax;
291 const T dy = by - ay;
292 return dx * dx + dy * dy;
293}
294
295template <typename T>
296T distance(const T ax, const T ay, const T bx, const T by)
297{
298 return std::sqrt(distanceSquared(ax, ay, bx, by));
299}
300
301template <typename T>
302T distance(const V2d<T>& a, const V2d<T>& b)
303{
304 return distance(a.x, a.y, b.x, b.y);
305}
306
307template <typename T>
308T distanceSquared(const V2d<T>& a, const V2d<T>& b)
309{
310 return distanceSquared(a.x, a.y, b.x, b.y);
311}
312
314{
315 return t.vertices[0] < 3 || t.vertices[1] < 3 || t.vertices[2] < 3;
316}
317
318} // namespace CDT
Utilities and helpers.
#define CDT_INLINE_IF_HEADER_ONLY
Macro for inlining non-template functions when in header-only mode to avoid multiple declaration erro...
Definition CDTUtils.h:38
Namespace containing triangulation functionality.
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY VertInd opposedVertex(const Triangle &tri, TriInd iTopo)
Given two triangles, return vertex of first triangle opposed to the second.
Definition CDTUtils.hpp:255
std::vector< TriInd > TriIndVec
Vector of triangle indices.
Definition CDTUtils.h:155
Box2d< T > envelopBox(TVertexIter first, TVertexIter last, TGetVertexCoordX getX, TGetVertexCoordY getY)
Bounding box of a collection of custom 2D points given coordinate getters.
Definition CDTUtils.h:187
CDT_INLINE_IF_HEADER_ONLY Index edgeNeighborInd(const VerticesArr3 &vv, VertInd iVedge1, VertInd iVedge2)
Index of triangle's neighbor opposed to an edge.
Definition CDTUtils.hpp:193
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opposedTriangleInd(const VerticesArr3 &vv, VertInd iVert)
Index of triangle's neighbor opposed to a vertex.
Definition CDTUtils.hpp:183
CDT_EXPORT PtLineLocation::Enum classifyOrientation(T orientation, T orientationTolerance=T(0))
Classify value of orient2d predicate.
Definition CDTUtils.hpp:116
array< TriInd, 3 > NeighborsArr3
array of three neighbors
Definition CDTUtils.h:157
CDT_EXPORT T distance(const V2d< T > &a, const V2d< T > &b)
Distance between two 2D points.
Definition CDTUtils.hpp:302
IndexSizeType VertInd
Vertex index.
Definition CDTUtils.h:142
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY bool verticesShareEdge(const TriIndVec &aTris, const TriIndVec &bTris)
Test if two vertices share at least one common triangle.
Definition CDTUtils.hpp:279
CDT_EXPORT Index cw(Index i)
Advance vertex or neighbor index clockwise.
Definition CDTUtils.hpp:80
array< VertInd, 3 > VerticesArr3
array of three vertex indices
Definition CDTUtils.h:156
CDT_INLINE_IF_HEADER_ONLY bool touchesSuperTriangle(const Triangle &t)
Check if any of triangle's vertices belongs to a super-triangle.
Definition CDTUtils.hpp:313
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opoVrt(Index neighborIndex)
Opposed vertex index from neighbor index.
Definition CDTUtils.hpp:170
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opoNbr(Index vertIndex)
Opposed neighbor index from vertex index.
Definition CDTUtils.hpp:158
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index vertexInd(const VerticesArr3 &vv, VertInd iV)
If triangle has a given vertex return vertex-index.
Definition CDTUtils.hpp:238
CDT_EXPORT T orient2D(const V2d< T > &p, const V2d< T > &v1, const V2d< T > &v2)
Orient p against line v1-v2 2D: robust geometric predicate.
Definition CDTUtils.hpp:99
CDT_EXPORT PtLineLocation::Enum locatePointLine(const V2d< T > &p, const V2d< T > &v1, const V2d< T > &v2, T orientationTolerance=T(0))
Check if point lies to the left of, to the right of, or on a line.
Definition CDTUtils.hpp:105
CDT_EXPORT T distanceSquared(const V2d< T > &a, const V2d< T > &b)
Squared distance between two 2D points.
Definition CDTUtils.hpp:308
CDT_EXPORT Index ccw(Index i)
Advance vertex or neighbor index counter-clockwise.
Definition CDTUtils.hpp:75
CDT_EXPORT Index edgeNeighbor(PtTriLocation::Enum location)
Neighbor index from a on-edge location.
Definition CDTUtils.hpp:92
const T & getX_V2d(const V2d< T > &v)
X- coordinate getter for V2d.
Definition CDTUtils.h:106
unsigned char Index
Index in triangle.
Definition CDTUtils.h:140
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY TriInd opposedTriangle(const Triangle &tri, VertInd iVert)
Given triangle and a vertex find opposed triangle.
Definition CDTUtils.hpp:249
CDT_EXPORT bool isInCircumcircle(const V2d< T > &p, const V2d< T > &v1, const V2d< T > &v2, const V2d< T > &v3)
Test if point lies in a circumscribed circle of a triangle.
Definition CDTUtils.hpp:268
CDT_EXPORT bool isOnEdge(PtTriLocation::Enum location)
Check if location is classified as on any of three edges.
Definition CDTUtils.hpp:85
IndexSizeType TriInd
Triangle index.
Definition CDTUtils.h:144
CDT_EXPORT PtTriLocation::Enum locatePointTriangle(const V2d< T > &p, const V2d< T > &v1, const V2d< T > &v2, const V2d< T > &v3)
Check if point a lies inside of, outside of, or on an edge of a triangle.
Definition CDTUtils.hpp:126
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opposedVertexInd(const NeighborsArr3 &nn, TriInd iTopo)
Index of triangle's vertex opposed to a triangle.
Definition CDTUtils.hpp:227
const T & getY_V2d(const V2d< T > &v)
Y-coordinate getter for V2d.
Definition CDTUtils.h:113
2D bounding box
Definition CDTUtils.h:162
Edge connecting two vertices: vertex with smaller index is always first.
Definition CDTUtils.h:209
bool operator==(const Edge &other) const
Equals operator.
Definition CDTUtils.hpp:47
bool operator!=(const Edge &other) const
Not-equals operator.
Definition CDTUtils.hpp:52
const std::pair< VertInd, VertInd > & verts() const
Edges' vertices.
Definition CDTUtils.hpp:67
Edge(VertInd iV1, VertInd iV2)
Constructor.
Definition CDTUtils.hpp:42
VertInd v2() const
V2 getter.
Definition CDTUtils.hpp:62
VertInd v1() const
V1 getter.
Definition CDTUtils.hpp:57
Triangulation triangle (counter-clockwise winding)
Definition CDTUtils.h:259
VerticesArr3 vertices
triangle's three vertices
Definition CDTUtils.h:260
NeighborsArr3 neighbors
triangle's three neighbors
Definition CDTUtils.h:261
2D vector
Definition CDTUtils.h:96
T y
Y-coordinate.
Definition CDTUtils.h:98
static V2d make(T x, T y)
Create vector from X and Y coordinates.
Definition CDTUtils.hpp:23
T x
X-coordinate.
Definition CDTUtils.h:97