CDT  1.4.2
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
9
10#include "CDTUtils.h"
11
12#include "predicates.h" // robust predicates: orient, in-circle
13
14#include <stdexcept>
15
16namespace CDT
17{
18
19CDT_INLINE_IF_HEADER_ONLY Index ccw(Index i)
20{
21 return Index((i + 1) % 3);
22}
23
24CDT_INLINE_IF_HEADER_ONLY Index cw(Index i)
25{
26 return Index((i + 2) % 3);
27}
28
29CDT_INLINE_IF_HEADER_ONLY bool isOnEdge(const PtTriLocation::Enum location)
30{
31 return location == PtTriLocation::OnEdge1 ||
32 location == PtTriLocation::OnEdge2 ||
33 location == PtTriLocation::OnEdge3;
34}
35
36CDT_INLINE_IF_HEADER_ONLY Index edgeNeighbor(const PtTriLocation::Enum location)
37{
38 assert(isOnEdge(location));
39 return static_cast<Index>(location - PtTriLocation::OnEdge1);
40}
41
42template <typename T>
43T orient2D(const V2d<T>& p, const V2d<T>& v1, const V2d<T>& v2)
44{
45 return predicates::adaptive::orient2d(v1.x, v1.y, v2.x, v2.y, p.x, p.y);
46}
47
48template <typename T>
50 const V2d<T>& p,
51 const V2d<T>& v1,
52 const V2d<T>& v2,
53 const T orientationTolerance)
54{
55 return classifyOrientation(orient2D(p, v1, v2), orientationTolerance);
56}
57
58template <typename T>
60classifyOrientation(const T orientation, const T orientationTolerance)
61{
62 if(orientation < -orientationTolerance)
63 return PtLineLocation::Right;
64 if(orientation > orientationTolerance)
65 return PtLineLocation::Left;
66 return PtLineLocation::OnLine;
67}
68
69template <typename T>
71 const V2d<T>& p,
72 const V2d<T>& v1,
73 const V2d<T>& v2,
74 const V2d<T>& v3)
75{
76 using namespace predicates::adaptive;
77 PtTriLocation::Enum result = PtTriLocation::Inside;
78 PtLineLocation::Enum edgeCheck = locatePointLine(p, v1, v2);
79 if(edgeCheck == PtLineLocation::Right)
80 return PtTriLocation::Outside;
81 if(edgeCheck == PtLineLocation::OnLine)
82 result = PtTriLocation::OnEdge1;
83 edgeCheck = locatePointLine(p, v2, v3);
84 if(edgeCheck == PtLineLocation::Right)
85 return PtTriLocation::Outside;
86 if(edgeCheck == PtLineLocation::OnLine)
87 {
88 result = (result == PtTriLocation::Inside) ? PtTriLocation::OnEdge2
89 : PtTriLocation::OnVertex;
90 }
91 edgeCheck = locatePointLine(p, v3, v1);
92 if(edgeCheck == PtLineLocation::Right)
93 return PtTriLocation::Outside;
94 if(edgeCheck == PtLineLocation::OnLine)
95 {
96 result = (result == PtTriLocation::Inside) ? PtTriLocation::OnEdge3
97 : PtTriLocation::OnVertex;
98 }
99 return result;
100}
101
102CDT_INLINE_IF_HEADER_ONLY Index opoNbr(const Index vertIndex)
103{
104 if(vertIndex == Index(0))
105 return Index(1);
106 if(vertIndex == Index(1))
107 return Index(2);
108 if(vertIndex == Index(2))
109 return Index(0);
110 assert(false && "Invalid vertex index");
111 handleException(std::runtime_error("Invalid vertex index"));
112 return invalidIndex;
113}
114
115CDT_INLINE_IF_HEADER_ONLY Index opoVrt(const Index neighborIndex)
116{
117 if(neighborIndex == Index(0))
118 return Index(2);
119 if(neighborIndex == Index(1))
120 return Index(0);
121 if(neighborIndex == Index(2))
122 return Index(1);
123 assert(false && "Invalid neighbor index");
124 handleException(std::runtime_error("Invalid neighbor index"));
125 return invalidIndex;
126}
127
128CDT_INLINE_IF_HEADER_ONLY Index
130{
131 assert(vv[0] == iVert || vv[1] == iVert || vv[2] == iVert);
132 if(vv[0] == iVert)
133 return Index(1);
134 if(vv[1] == iVert)
135 return Index(2);
136 return Index(0);
137}
138
139CDT_INLINE_IF_HEADER_ONLY Index edgeNeighborInd(
140 const VerticesArr3& vv,
141 const VertInd iVedge1,
142 const VertInd iVedge2)
143{
144 assert(vv[0] == iVedge1 || vv[1] == iVedge1 || vv[2] == iVedge1);
145 assert(vv[0] == iVedge2 || vv[1] == iVedge2 || vv[2] == iVedge2);
146 assert(
147 (vv[0] != iVedge1 && vv[0] != iVedge2) ||
148 (vv[1] != iVedge1 && vv[1] != iVedge2) ||
149 (vv[2] != iVedge1 && vv[2] != iVedge2));
150 /*
151 * vv[2]
152 * /\
153 * n[2]/ \n[1]
154 * /____\
155 * vv[0] n[0] vv[1]
156 */
157 if(vv[0] == iVedge1)
158 {
159 if(vv[1] == iVedge2)
160 return Index(0);
161 return Index(2);
162 }
163 if(vv[0] == iVedge2)
164 {
165 if(vv[1] == iVedge1)
166 return Index(0);
167 return Index(2);
168 }
169 return Index(1);
170}
171
172CDT_INLINE_IF_HEADER_ONLY Index
173opposedVertexInd(const NeighborsArr3& nn, const TriInd iTopo)
174{
175 assert(nn[0] == iTopo || nn[1] == iTopo || nn[2] == iTopo);
176 if(nn[0] == iTopo)
177 return Index(2);
178 if(nn[1] == iTopo)
179 return Index(0);
180 return Index(1);
181}
182
183CDT_INLINE_IF_HEADER_ONLY Index
184vertexInd(const VerticesArr3& vv, const VertInd iV)
185{
186 assert(vv[0] == iV || vv[1] == iV || vv[2] == iV);
187 if(vv[0] == iV)
188 return Index(0);
189 if(vv[1] == iV)
190 return Index(1);
191 return Index(2);
192}
193
194CDT_INLINE_IF_HEADER_ONLY TriInd
195opposedTriangle(const Triangle& tri, const VertInd iVert)
196{
197 return tri.neighbors[opposedTriangleInd(tri.vertices, iVert)];
198}
199
200CDT_INLINE_IF_HEADER_ONLY VertInd
201opposedVertex(const Triangle& tri, const TriInd iTopo)
202{
203 return tri.vertices[opposedVertexInd(tri.neighbors, iTopo)];
204}
205
207CDT_INLINE_IF_HEADER_ONLY TriInd
208edgeNeighbor(const Triangle& tri, VertInd iVedge1, VertInd iVedge2)
209{
210 return tri.neighbors[edgeNeighborInd(tri.vertices, iVedge1, iVedge2)];
211}
212
213template <typename T>
215 const V2d<T>& p,
216 const V2d<T>& v1,
217 const V2d<T>& v2,
218 const V2d<T>& v3)
219{
220 using namespace predicates::adaptive;
221 return incircle(v1.x, v1.y, v2.x, v2.y, v3.x, v3.y, p.x, p.y) > T(0);
222}
223
224CDT_INLINE_IF_HEADER_ONLY
225bool verticesShareEdge(const TriIndVec& aTris, const TriIndVec& bTris)
226{
227 for(TriIndVec::const_iterator it = aTris.begin(); it != aTris.end(); ++it)
228 if(std::find(bTris.begin(), bTris.end(), *it) != bTris.end())
229 return true;
230 return false;
231}
232
233template <typename T>
234T distanceSquared(const T ax, const T ay, const T bx, const T by)
235{
236 const T dx = bx - ax;
237 const T dy = by - ay;
238 return dx * dx + dy * dy;
239}
240
241template <typename T>
242T distance(const T ax, const T ay, const T bx, const T by)
243{
244 return std::sqrt(distanceSquared(ax, ay, bx, by));
245}
246
247template <typename T>
248T distance(const V2d<T>& a, const V2d<T>& b)
249{
250 return distance(a.x, a.y, b.x, b.y);
251}
252
253template <typename T>
254T distanceSquared(const V2d<T>& a, const V2d<T>& b)
255{
256 return distanceSquared(a.x, a.y, b.x, b.y);
257}
258
260{
261 return t.vertices[0] < 3 || t.vertices[1] < 3 || t.vertices[2] < 3;
262}
263
264} // namespace CDT
Utilities and helpers.
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:201
std::vector< TriInd > TriIndVec
Vector of triangle indices.
Definition CDTUtils.h:209
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:139
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:129
CDT_EXPORT PtLineLocation::Enum classifyOrientation(T orientation, T orientationTolerance=T(0))
Classify value of orient2d predicate.
Definition CDTUtils.hpp:60
array< TriInd, 3 > NeighborsArr3
array of three neighbors
Definition CDTUtils.h:211
CDT_EXPORT T distance(const V2d< T > &a, const V2d< T > &b)
Distance between two 2D points.
Definition CDTUtils.hpp:248
IndexSizeType VertInd
Vertex index.
Definition CDTUtils.h:191
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:225
CDT_EXPORT Index cw(Index i)
Advance vertex or neighbor index clockwise.
Definition CDTUtils.hpp:24
array< VertInd, 3 > VerticesArr3
array of three vertex indices
Definition CDTUtils.h:210
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:259
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opoVrt(Index neighborIndex)
Opposed vertex index from neighbor index.
Definition CDTUtils.hpp:115
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opoNbr(Index vertIndex)
Opposed neighbor index from vertex index.
Definition CDTUtils.hpp:102
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:184
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:43
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:49
CDT_EXPORT T distanceSquared(const V2d< T > &a, const V2d< T > &b)
Squared distance between two 2D points.
Definition CDTUtils.hpp:254
CDT_EXPORT Index ccw(Index i)
Advance vertex or neighbor index counter-clockwise.
Definition CDTUtils.hpp:19
CDT_EXPORT Index edgeNeighbor(PtTriLocation::Enum location)
Neighbor index from a on-edge location.
Definition CDTUtils.hpp:36
unsigned char Index
Index in triangle.
Definition CDTUtils.h:189
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:195
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:214
CDT_EXPORT bool isOnEdge(PtTriLocation::Enum location)
Check if location is classified as on any of three edges.
Definition CDTUtils.hpp:29
IndexSizeType TriInd
Triangle index.
Definition CDTUtils.h:193
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:70
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:173
Triangulation triangle (counter-clockwise winding)
Definition CDTUtils.h:344
VerticesArr3 vertices
triangle's three vertices
Definition CDTUtils.h:345
NeighborsArr3 neighbors
triangle's three neighbors
Definition CDTUtils.h:346
2D vector
Definition CDTUtils.h:136
T y
Y-coordinate.
Definition CDTUtils.h:138
T x
X-coordinate.
Definition CDTUtils.h:137