11#include "portable_nth_element.hpp"
22typedef std::deque<TriInd> TriDeque;
30const std::size_t nTargetVerts = 0;
36const float minDistToConstraintEdge(0);
42template <
typename T,
typename TNearPo
intLocator>
44 : m_nTargetVerts(detail::defaults::nTargetVerts)
45 , m_superGeomType(detail::defaults::superGeomType)
46 , m_vertexInsertionOrder(detail::defaults::vertexInsertionOrder)
47 , m_intersectingEdgesStrategy(detail::defaults::intersectingEdgesStrategy)
48 , m_minDistToConstraintEdge(detail::defaults::minDistToConstraintEdge)
49#ifdef CDT_ENABLE_CALLBACK_HANDLER
50 , m_callbackHandler(NULL)
54template <
typename T,
typename TNearPo
intLocator>
57 : m_nTargetVerts(detail::defaults::nTargetVerts)
58 , m_superGeomType(detail::defaults::superGeomType)
59 , m_vertexInsertionOrder(vertexInsertionOrder)
60 , m_intersectingEdgesStrategy(detail::defaults::intersectingEdgesStrategy)
61 , m_minDistToConstraintEdge(detail::defaults::minDistToConstraintEdge)
62#ifdef CDT_ENABLE_CALLBACK_HANDLER
63 , m_callbackHandler(NULL)
67template <
typename T,
typename TNearPo
intLocator>
71 const T minDistToConstraintEdge)
72 : m_nTargetVerts(detail::defaults::nTargetVerts)
73 , m_superGeomType(detail::defaults::superGeomType)
74 , m_vertexInsertionOrder(vertexInsertionOrder)
75 , m_intersectingEdgesStrategy(intersectingEdgesStrategy)
76 , m_minDistToConstraintEdge(minDistToConstraintEdge)
77#ifdef CDT_ENABLE_CALLBACK_HANDLER
78 , m_callbackHandler(NULL)
82template <
typename T,
typename TNearPo
intLocator>
85 const TNearPointLocator& nearPtLocator,
87 const T minDistToConstraintEdge)
88 : m_nearPtLocator(nearPtLocator)
89 , m_nTargetVerts(detail::defaults::nTargetVerts)
90 , m_superGeomType(detail::defaults::superGeomType)
91 , m_vertexInsertionOrder(vertexInsertionOrder)
92 , m_intersectingEdgesStrategy(intersectingEdgesStrategy)
93 , m_minDistToConstraintEdge(minDistToConstraintEdge)
94#ifdef CDT_ENABLE_CALLBACK_HANDLER
95 , m_callbackHandler(NULL)
99template <
typename T,
typename TNearPo
intLocator>
111 finalizeTriangulation(toErase);
114template <
typename T,
typename TNearPo
intLocator>
117 assert(m_vertTris[0] != noNeighbor);
118 const std::stack<TriInd> seed(std::deque<TriInd>(1, m_vertTris[0]));
119 const TriIndUSet toErase = growToBoundary(seed);
120 finalizeTriangulation(toErase);
123template <
typename T,
typename TNearPo
intLocator>
129 for(std::size_t iT = 0; iT !=
triangles.size(); ++iT)
131 if(triDepths[iT] % 2 == 0)
132 toErase.insert(
static_cast<TriInd>(iT));
134 finalizeTriangulation(toErase);
143template <
typename T,
typename TNearPo
intLocator>
147 if(removedTriangles.empty())
153 if(removedTriangles.count(iT))
155 triIndMap[iT] = iTnew;
166 for(NeighborsArr3::iterator n = nn.begin(); n != nn.end(); ++n)
168 if(removedTriangles.count(*n))
172 else if(*n != noNeighbor)
180template <
typename T,
typename TNearPo
intLocator>
186template <
typename T,
typename TNearPo
intLocator>
192template <
typename T,
typename TNearPo
intLocator>
193void Triangulation<T, TNearPointLocator>::finalizeTriangulation(
201 vertices.begin(), vertices.begin() + nSuperTriangleVertices);
205 typedef CDT::EdgeUSet::const_iterator It;
206 for(It e = fixedEdges.begin(); e != fixedEdges.end(); ++e)
210 fixedEdges = updatedFixedEdges;
213 unordered_map<Edge, BoundaryOverlapCount> updatedOverlapCount;
214 typedef unordered_map<Edge, BoundaryOverlapCount>::const_iterator
216 for(It it = overlapCount.begin(); it != overlapCount.end(); ++it)
218 updatedOverlapCount.insert(
222 overlapCount = updatedOverlapCount;
225 unordered_map<Edge, EdgeVec> updatedPieceToOriginals;
226 typedef unordered_map<Edge, EdgeVec>::const_iterator It;
227 for(It it = pieceToOriginals.begin(); it != pieceToOriginals.end();
231 for(EdgeVec::iterator eeIt = ee.begin(); eeIt != ee.end();
236 updatedPieceToOriginals.insert(
239 pieceToOriginals = updatedPieceToOriginals;
243 removeTriangles(removedTriangles);
247 for(TriangleVec::iterator t = triangles.begin(); t != triangles.end();
251 for(VerticesArr3::iterator v = vv.begin(); v != vv.end(); ++v)
253 *v -= nSuperTriangleVertices;
259template <
typename T,
typename TNearPo
intLocator>
262 m_nearPtLocator.initialize(
vertices);
263 m_nTargetVerts =
static_cast<IndexSizeType
>(
vertices.size());
267template <
typename T,
typename TNearPo
intLocator>
268TriIndUSet Triangulation<T, TNearPointLocator>::growToBoundary(
269 std::stack<TriInd> seeds)
const
272 while(!seeds.empty())
274 const TriInd iT = seeds.top();
276 traversed.insert(iT);
281 if(fixedEdges.count(opEdge))
284 if(iN != noNeighbor && traversed.count(iN) == 0)
291template <
typename T,
typename TNearPo
intLocator>
292TriInd Triangulation<T, TNearPointLocator>::addTriangle(
const Triangle& t)
294 const TriInd iT(triangles.size());
295 triangles.push_back(t);
299template <
typename T,
typename TNearPo
intLocator>
300TriInd Triangulation<T, TNearPointLocator>::addTriangle()
305template <
typename T,
typename TNearPo
intLocator>
307 const std::vector<Edge>& edges)
312template <
typename T,
typename TNearPo
intLocator>
314 const std::vector<Edge>& edges)
319template <
typename T,
typename TNearPo
intLocator>
320void Triangulation<T, TNearPointLocator>::fixEdge(
const Edge& edge)
322 if(!fixedEdges.insert(edge).second)
324 ++overlapCount[edge];
332template <
typename T,
typename Allocator1>
333void insert_unique(std::vector<T, Allocator1>& to,
const T& elem)
335 if(std::find(to.begin(), to.end(), elem) == to.end())
342template <
typename T,
typename Allocator1,
typename Allocator2>
344 std::vector<T, Allocator1>& to,
345 const std::vector<T, Allocator2>& from)
347 typedef typename std::vector<T, Allocator2>::const_iterator Cit;
348 to.reserve(to.size() + from.size());
349 for(Cit cit = from.begin(); cit != from.end(); ++cit)
351 insert_unique(to, *cit);
357template <
typename T,
typename TNearPo
intLocator>
358void Triangulation<T, TNearPointLocator>::splitFixedEdge(
363 const Edge half1(edge.v1(), iSplitVert);
364 const Edge half2(iSplitVert, edge.v2());
366 fixedEdges.erase(edge);
370 typedef unordered_map<Edge, BoundaryOverlapCount>::const_iterator It;
371 const It overlapIt = overlapCount.find(edge);
372 if(overlapIt != overlapCount.end())
374 overlapCount[half1] += overlapIt->second;
375 overlapCount[half2] += overlapIt->second;
376 overlapCount.erase(overlapIt);
380 const unordered_map<Edge, EdgeVec>::const_iterator originalsIt =
381 pieceToOriginals.find(edge);
382 if(originalsIt != pieceToOriginals.end())
384 newOriginals = originalsIt->second;
385 pieceToOriginals.erase(originalsIt);
387 detail::insert_unique(pieceToOriginals[half1], newOriginals);
388 detail::insert_unique(pieceToOriginals[half2], newOriginals);
391template <
typename T,
typename TNearPo
intLocator>
392VertInd Triangulation<T, TNearPointLocator>::addSplitEdgeVertex(
399 addNewVertex(splitVert, noNeighbor);
401#ifdef CDT_ENABLE_CALLBACK_HANDLER
402 if(m_callbackHandler)
404 m_callbackHandler->onAddVertexStart(
409 std::stack<TriInd> triStack = insertVertexOnEdge(iSplitVert, iT, iTopo);
410 tryAddVertexToLocator(iSplitVert);
411 ensureDelaunayByEdgeFlips(iSplitVert, triStack);
415template <
typename T,
typename TNearPo
intLocator>
416VertInd Triangulation<T, TNearPointLocator>::splitFixedEdgeAt(
422 const VertInd iSplitVert = addSplitEdgeVertex(splitVert, iT, iTopo);
423 splitFixedEdge(edge, iSplitVert);
427template <
typename T,
typename TNearPo
intLocator>
428void Triangulation<T, TNearPointLocator>::fixEdge(
430 const Edge& originalEdge)
433 if(edge != originalEdge)
434 detail::insert_unique(pieceToOriginals[edge], originalEdge);
441T lerp(
const T& a,
const T& b,
const T t)
443 return (T(1) - t) * a + t * b;
448V2d<T> intersectionPosition(
454 using namespace predicates::adaptive;
458 const T a_cd = orient2d(c.x, c.y, d.x, d.y, a.x, a.y);
459 const T b_cd = orient2d(c.x, c.y, d.x, d.y, b.x, b.y);
460 const T t_ab = a_cd / (a_cd - b_cd);
462 const T c_ab = orient2d(a.x, a.y, b.x, b.y, c.x, c.y);
463 const T d_ab = orient2d(a.x, a.y, b.x, b.y, d.x, d.y);
464 const T t_cd = c_ab / (c_ab - d_ab);
467 std::fabs(a.x - b.x) < std::fabs(c.x - d.x) ? lerp(a.x, b.x, t_ab)
468 : lerp(c.x, d.x, t_cd),
469 std::fabs(a.y - b.y) < std::fabs(c.y - d.y) ? lerp(a.y, b.y, t_ab)
470 : lerp(c.y, d.y, t_cd));
475template <
typename T,
typename TNearPo
intLocator>
476void Triangulation<T, TNearPointLocator>::insertEdgeIteration(
478 const Edge originalEdge,
480 std::vector<TriangulatePseudoPolygonTask>& tppIterations)
489 fixEdge(edge, originalEdge);
493 const V2d<T>& a = vertices[iA];
494 const V2d<T>& b = vertices[iB];
495 const T distanceTolerance =
496 m_minDistToConstraintEdge == T(0)
498 : m_minDistToConstraintEdge *
distance(a, b);
503 tie(iT, iVL, iVR) = intersectedTriangle(iA, a, b, distanceTolerance);
507 const Edge edgePart(iA, iVL);
508 fixEdge(edgePart, originalEdge);
509 remaining.push_back(
Edge(iVL, iB));
513 std::vector<TriInd> intersected(1, iT);
514 std::vector<VertInd> polyL, polyR;
517 polyL.push_back(iVL);
520 polyR.push_back(iVR);
521 unordered_map<Edge, TriInd> outerTris;
526 while(!t.containsVertex(iB))
529 const Triangle& tOpo = triangles[iTopo];
532 switch(m_intersectingEdgesStrategy)
535 if(fixedEdges.count(
Edge(iVL, iVR)))
538 Edge e1 = originalEdge;
540 e2 = pieceToOriginals.count(e2)
541 ? pieceToOriginals.at(e2).front()
545 VertInd(e1.v1() - m_nTargetVerts),
546 VertInd(e1.v2() - m_nTargetVerts));
548 VertInd(e2.v1() - m_nTargetVerts),
549 VertInd(e2.v2() - m_nTargetVerts));
552 pieceToOriginals.count(e2) ? pieceToOriginals.at(e2).front()
554 CDT_SOURCE_LOCATION));
558 if(!fixedEdges.count(
Edge(iVL, iVR)))
561 const V2d<T> newV = detail::intersectionPosition(
562 vertices[iA], vertices[iB], vertices[iVL], vertices[iVR]);
564 splitFixedEdgeAt(
Edge(iVL, iVR), newV, iT, iTopo);
567 remaining.push_back(
Edge(iA, iNewVert));
568 remaining.push_back(
Edge(iNewVert, iB));
572 assert(!fixedEdges.count(
Edge(iVL, iVR)));
578 if(loc == PtLineLocation::Left)
580 const Edge e(polyL.back(), iVopo);
582 if(!outerTris.insert(std::make_pair(e, outer)).second)
583 outerTris.at(e) = noNeighbor;
584 polyL.push_back(iVopo);
588 else if(loc == PtLineLocation::Right)
590 const Edge e(polyR.back(), iVopo);
592 if(!outerTris.insert(std::make_pair(e, outer)).second)
593 outerTris.at(e) = noNeighbor;
594 polyR.push_back(iVopo);
601 intersected.push_back(iTopo);
610 assert(!intersected.empty());
613 if(m_vertTris[iA] == intersected.front())
614 pivotVertexTriangleCW(iA);
615 if(m_vertTris[iB] == intersected.back())
616 pivotVertexTriangleCW(iB);
619#ifdef CDT_ENABLE_CALLBACK_HANDLER
620 if(m_callbackHandler)
622 m_callbackHandler->onReTriangulatePolygon(intersected);
627 std::reverse(polyR.begin(), polyR.end());
632 assert(intersected.size() >= 2);
633 const TriInd iTL = intersected.back();
634 intersected.pop_back();
635 const TriInd iTR = intersected.back();
636 intersected.pop_back();
638 triangulatePseudoPolygon(
639 polyL, outerTris, iTL, iTR, intersected, tppIterations);
640 triangulatePseudoPolygon(
641 polyR, outerTris, iTR, iTL, intersected, tppIterations);
642 assert(intersected.empty());
648 const Edge edgePart(iA, iB);
649 fixEdge(edgePart, originalEdge);
650 remaining.push_back(
Edge(iB, edge.v2()));
655 fixEdge(edge, originalEdge);
659template <
typename T,
typename TNearPo
intLocator>
660void Triangulation<T, TNearPointLocator>::insertEdge(
662 const Edge originalEdge,
664 std::vector<TriangulatePseudoPolygonTask>& tppIterations)
666#ifdef CDT_ENABLE_CALLBACK_HANDLER
667 if(m_callbackHandler)
669 m_callbackHandler->onAddEdgeStart(edge);
675 remaining.push_back(edge);
676 while(!remaining.empty())
678 edge = remaining.back();
679 remaining.pop_back();
680 insertEdgeIteration(edge, originalEdge, remaining, tppIterations);
684template <
typename T,
typename TNearPo
intLocator>
685void Triangulation<T, TNearPointLocator>::conformToEdgeIteration(
688 BoundaryOverlapCount overlaps,
689 std::vector<ConformToEdgeTask>& remaining)
700 overlapCount[edge] = overlaps;
702 if(!originals.empty() && edge != originals.front())
704 detail::insert_unique(pieceToOriginals[edge], originals);
709 const V2d<T>& a = vertices[iA];
710 const V2d<T>& b = vertices[iB];
711 const T distanceTolerance =
712 m_minDistToConstraintEdge == T(0)
714 : m_minDistToConstraintEdge *
distance(a, b);
717 tie(iT, iVleft, iVright) = intersectedTriangle(iA, a, b, distanceTolerance);
721 const Edge edgePart(iA, iVleft);
724 overlapCount[edgePart] = overlaps;
725 detail::insert_unique(pieceToOriginals[edgePart], originals);
726#ifdef CDT_CXX11_IS_SUPPORTED
727 remaining.emplace_back(
Edge(iVleft, iB), originals, overlaps);
729 remaining.push_back(make_tuple(
Edge(iVleft, iB), originals, overlaps));
736 while(std::find(t.vertices.begin(), t.vertices.end(), iB) ==
740 const Triangle& tOpo = triangles[iTopo];
742 const V2d<T> vOpo = vertices[iVopo];
744 switch(m_intersectingEdgesStrategy)
747 if(fixedEdges.count(
Edge(iVleft, iVright)))
750 Edge e1 = pieceToOriginals.count(edge)
751 ? pieceToOriginals.at(edge).front()
753 Edge e2(iVleft, iVright);
754 e2 = pieceToOriginals.count(e2)
755 ? pieceToOriginals.at(e2).front()
759 VertInd(e1.v1() - m_nTargetVerts),
760 VertInd(e1.v2() - m_nTargetVerts));
762 VertInd(e2.v1() - m_nTargetVerts),
763 VertInd(e2.v2() - m_nTargetVerts));
769 if(!fixedEdges.count(
Edge(iVleft, iVright)))
772 const V2d<T> newV = detail::intersectionPosition(
778 splitFixedEdgeAt(
Edge(iVleft, iVright), newV, iT, iTopo);
779#ifdef CDT_CXX11_IS_SUPPORTED
780 remaining.emplace_back(
Edge(iNewVert, iB), originals, overlaps);
781 remaining.emplace_back(
Edge(iA, iNewVert), originals, overlaps);
784 make_tuple(
Edge(iNewVert, iB), originals, overlaps));
786 make_tuple(
Edge(iA, iNewVert), originals, overlaps));
791 assert(!fixedEdges.count(
Edge(iVleft, iVright)));
800 if(loc == PtLineLocation::Left)
805 else if(loc == PtLineLocation::Right)
817#ifdef CDT_CXX11_IS_SUPPORTED
818 remaining.emplace_back(
Edge(iB, edge.v2()), originals, overlaps);
821 make_tuple(
Edge(iB, edge.v2()), originals, overlaps));
827 const V2d<T>& start = vertices[iA];
828 const V2d<T>& end = vertices[iB];
830 V2d<T>((start.x + end.x) / T(2), (start.y + end.y) / T(2)), noNeighbor);
832#ifdef CDT_ENABLE_CALLBACK_HANDLER
833 if(m_callbackHandler)
835 m_callbackHandler->onAddVertexStart(
840 const std::vector<Edge> flippedFixedEdges =
841 insertVertex_FlipFixedEdges(iMid);
843#ifdef CDT_CXX11_IS_SUPPORTED
844 remaining.emplace_back(
Edge(iMid, iB), originals, overlaps);
845 remaining.emplace_back(
Edge(iA, iMid), originals, overlaps);
847 remaining.push_back(make_tuple(
Edge(iMid, iB), originals, overlaps));
848 remaining.push_back(make_tuple(
Edge(iA, iMid), originals, overlaps));
853 for(std::vector<Edge>::const_iterator it = flippedFixedEdges.begin();
854 it != flippedFixedEdges.end();
857 const Edge& flippedFixedEdge = *it;
858 fixedEdges.erase(flippedFixedEdge);
860 BoundaryOverlapCount prevOverlaps = 0;
861 const unordered_map<Edge, BoundaryOverlapCount>::const_iterator
862 overlapsIt = overlapCount.find(flippedFixedEdge);
863 if(overlapsIt != overlapCount.end())
865 prevOverlaps = overlapsIt->second;
866 overlapCount.erase(overlapsIt);
869 EdgeVec prevOriginals(1, flippedFixedEdge);
870 const unordered_map<Edge, EdgeVec>::const_iterator originalsIt =
871 pieceToOriginals.find(flippedFixedEdge);
872 if(originalsIt != pieceToOriginals.end())
874 prevOriginals = originalsIt->second;
876#ifdef CDT_CXX11_IS_SUPPORTED
877 remaining.emplace_back(flippedFixedEdge, prevOriginals, prevOverlaps);
880 make_tuple(flippedFixedEdge, prevOriginals, prevOverlaps));
885template <
typename T,
typename TNearPo
intLocator>
886void Triangulation<T, TNearPointLocator>::conformToEdge(
889 BoundaryOverlapCount overlaps,
890 std::vector<ConformToEdgeTask>& remaining)
892#ifdef CDT_ENABLE_CALLBACK_HANDLER
893 if(m_callbackHandler)
895 m_callbackHandler->onAddEdgeStart(edge);
901#ifdef CDT_CXX11_IS_SUPPORTED
902 remaining.emplace_back(edge, originals, overlaps);
904 remaining.push_back(make_tuple(edge, originals, overlaps));
906 while(!remaining.empty())
908 tie(edge, originals, overlaps) = remaining.back();
909 remaining.pop_back();
910 conformToEdgeIteration(edge, originals, overlaps, remaining);
924template <
typename T,
typename TNearPo
intLocator>
925tuple<TriInd, VertInd, VertInd>
926Triangulation<T, TNearPointLocator>::intersectedTriangle(
930 const T orientationTolerance)
const
932 const TriInd startTri = m_vertTris[iA];
939 const T orientP2 =
orient2D(vertices[iP2], a, b);
941 if(locP2 == PtLineLocation::Right)
944 const T orientP1 =
orient2D(vertices[iP1], a, b);
946 if(locP1 == PtLineLocation::OnLine)
948 return make_tuple(noNeighbor, iP1, iP1);
950 if(locP1 == PtLineLocation::Left)
952 if(orientationTolerance)
956 if(std::abs(orientP1) <= std::abs(orientP2))
958 closestOrient = orientP1;
963 closestOrient = orientP2;
967 closestOrient, orientationTolerance) ==
968 PtLineLocation::OnLine)
970 return make_tuple(noNeighbor, iClosestP, iClosestP);
973 return make_tuple(iT, iP1, iP2);
976 iT = t.next(iA).first;
977 }
while(iT != startTri);
979 handleException(
Error(
980 "Could not find vertex triangle intersected by an edge.",
981 CDT_SOURCE_LOCATION));
982 return make_tuple(noNeighbor, noVertex, noVertex);
985template <
typename T,
typename TNearPo
intLocator>
986void Triangulation<T, TNearPointLocator>::addSuperTriangle(
const Box2d<T>& box)
988 m_nTargetVerts = nSuperTriangleVertices;
992 (box.min.x + box.max.x) / T(2), (box.min.y + box.max.y) / T(2)};
993 const T w = box.max.
x - box.min.x;
994 const T h = box.max.y - box.min.y;
995 T r = std::max(w, h);
1000 r = std::max(T(2) * r, T(1));
1006 while(center.y <= center.y - r)
1010 const T R = T(2) * r;
1011 const T cos_30_deg = T(0.8660254037844386);
1012 const T shiftX = R * cos_30_deg;
1013 const V2d<T> posV1 = {center.x - shiftX, center.y - r};
1014 const V2d<T> posV2 = {center.x + shiftX, center.y - r};
1015 const V2d<T> posV3 = {center.x, center.y + R};
1016 addNewVertex(posV1,
TriInd(0));
1017 addNewVertex(posV2,
TriInd(0));
1018 addNewVertex(posV3,
TriInd(0));
1020#ifdef CDT_ENABLE_CALLBACK_HANDLER
1021 if(m_callbackHandler)
1023 m_callbackHandler->onAddSuperTriangle();
1032 m_nearPtLocator.initialize(vertices);
1036template <
typename T,
typename TNearPo
intLocator>
1037void Triangulation<T, TNearPointLocator>::addNewVertex(
1041 vertices.push_back(pos);
1042 m_vertTris.push_back(iT);
1045template <
typename T,
typename TNearPo
intLocator>
1047Triangulation<T, TNearPointLocator>::insertVertex_FlipFixedEdges(
1050 std::vector<Edge> flippedFixedEdges;
1052 const V2d<T>& v1 = vertices[iV1];
1053 const VertInd startVertex = m_nearPtLocator.nearPoint(v1, vertices);
1054 array<TriInd, 2> trisAt = walkingSearchTrianglesAt(iV1, startVertex);
1055 std::stack<TriInd> triStack =
1056 trisAt[1] == noNeighbor ? insertVertexInsideTriangle(iV1, trisAt[0])
1057 : insertVertexOnEdge(iV1, trisAt[0], trisAt[1]);
1059 TriInd iTopo, n1, n2, n3, n4;
1061 while(!triStack.empty())
1063 const TriInd iT = triStack.top();
1066 edgeFlipInfo(iT, iV1, iTopo, iV2, iV3, iV4, n1, n2, n3, n4);
1067 if(iTopo != noNeighbor && isFlipNeeded(iV1, iV2, iV3, iV4))
1070 const Edge flippedEdge(iV2, iV4);
1071 if(!fixedEdges.empty() &&
1072 fixedEdges.find(flippedEdge) != fixedEdges.end())
1074 flippedFixedEdges.push_back(flippedEdge);
1077 flipEdge(iT, iTopo, iV1, iV2, iV3, iV4, n1, n2, n3, n4);
1079 triStack.push(iTopo);
1083 tryAddVertexToLocator(iV1);
1084 return flippedFixedEdges;
1087template <
typename T,
typename TNearPo
intLocator>
1088void Triangulation<T, TNearPointLocator>::insertVertex(
1089 const VertInd iVert,
1090 const VertInd walkStart)
1092#ifdef CDT_ENABLE_CALLBACK_HANDLER
1093 if(m_callbackHandler)
1095 m_callbackHandler->onAddVertexStart(iVert, AddVertexType::UserInput);
1099 const array<TriInd, 2> trisAt = walkingSearchTrianglesAt(iVert, walkStart);
1100 std::stack<TriInd> triStack =
1101 trisAt[1] == noNeighbor
1102 ? insertVertexInsideTriangle(iVert, trisAt[0])
1103 : insertVertexOnEdge(iVert, trisAt[0], trisAt[1]);
1104 ensureDelaunayByEdgeFlips(iVert, triStack);
1107template <
typename T,
typename TNearPo
intLocator>
1108void Triangulation<T, TNearPointLocator>::insertVertex(
const VertInd iVert)
1110 const V2d<T>& v = vertices[iVert];
1111 const VertInd walkStart = m_nearPtLocator.nearPoint(v, vertices);
1112 insertVertex(iVert, walkStart);
1113 tryAddVertexToLocator(iVert);
1116template <
typename T,
typename TNearPo
intLocator>
1117void Triangulation<T, TNearPointLocator>::ensureDelaunayByEdgeFlips(
1119 std::stack<TriInd>& triStack)
1121 TriInd iTopo, n1, n2, n3, n4;
1123 while(!triStack.empty())
1125 const TriInd iT = triStack.top();
1128 edgeFlipInfo(iT, iV1, iTopo, iV2, iV3, iV4, n1, n2, n3, n4);
1129 if(iTopo != noNeighbor && isFlipNeeded(iV1, iV2, iV3, iV4))
1131 flipEdge(iT, iTopo, iV1, iV2, iV3, iV4, n1, n2, n3, n4);
1133 triStack.push(iTopo);
1151template <
typename T,
typename TNearPo
intLocator>
1152void Triangulation<T, TNearPointLocator>::edgeFlipInfo(
1169 const Triangle& t = triangles[iT];
1170 if(t.vertices[0] == iV1)
1172 iV2 = t.vertices[1];
1173 iV4 = t.vertices[2];
1174 n1 = t.neighbors[0];
1175 n3 = t.neighbors[2];
1176 iTopo = t.neighbors[1];
1178 else if(t.vertices[1] == iV1)
1180 iV2 = t.vertices[2];
1181 iV4 = t.vertices[0];
1182 n1 = t.neighbors[1];
1183 n3 = t.neighbors[0];
1184 iTopo = t.neighbors[2];
1188 iV2 = t.vertices[0];
1189 iV4 = t.vertices[1];
1190 n1 = t.neighbors[2];
1191 n3 = t.neighbors[1];
1192 iTopo = t.neighbors[0];
1194 if(iTopo == noNeighbor)
1196 const Triangle& tOpo = triangles[iTopo];
1197 if(tOpo.neighbors[0] == iT)
1199 iV3 = tOpo.vertices[2];
1200 n2 = tOpo.neighbors[1];
1201 n4 = tOpo.neighbors[2];
1203 else if(tOpo.neighbors[1] == iT)
1205 iV3 = tOpo.vertices[0];
1206 n2 = tOpo.neighbors[2];
1207 n4 = tOpo.neighbors[0];
1211 iV3 = tOpo.vertices[1];
1212 n2 = tOpo.neighbors[0];
1213 n4 = tOpo.neighbors[1];
1240template <
typename T,
typename TNearPo
intLocator>
1241bool Triangulation<T, TNearPointLocator>::isFlipNeeded(
1245 const VertInd iV4)
const
1247 if(fixedEdges.count(Edge(iV2, iV4)))
1249 const V2d<T>& v1 = vertices[iV1];
1250 const V2d<T>& v2 = vertices[iV2];
1251 const V2d<T>& v3 = vertices[iV3];
1252 const V2d<T>& v4 = vertices[iV4];
1253 if(m_superGeomType == SuperGeometryType::SuperTriangle)
1310template <
typename T,
typename TNearPo
intLocator>
1315#ifdef CDT_ENABLE_CALLBACK_HANDLER
1316 if(m_callbackHandler)
1318 m_callbackHandler->onFlipEdge(iT, iTopo);
1324 const array<TriInd, 3>& triNs = t.
neighbors;
1325 const array<TriInd, 3>& triOpoNs = tOpo.
neighbors;
1326 const array<VertInd, 3>& triVs = t.
vertices;
1327 const array<VertInd, 3>& triOpoVs = tOpo.
vertices;
1332 const TriInd n1 = triNs[i];
1335 const VertInd v3 = triOpoVs[i];
1337 const TriInd n4 = triOpoNs[i];
1343 changeNeighbor(n1, iT, iTopo);
1344 changeNeighbor(n4, iTopo, iT);
1350 setAdjacentTriangle(v4, iT);
1351 setAdjacentTriangle(v2, iTopo);
1371template <
typename T,
typename TNearPo
intLocator>
1384#ifdef CDT_ENABLE_CALLBACK_HANDLER
1385 if(m_callbackHandler)
1387 m_callbackHandler->onFlipEdge(iT, iTopo);
1395 changeNeighbor(n1, iT, iTopo);
1396 changeNeighbor(n4, iTopo, iT);
1402 setAdjacentTriangle(v4, iT);
1403 setAdjacentTriangle(v2, iTopo);
1423template <
typename T,
typename TNearPo
intLocator>
1425Triangulation<T, TNearPointLocator>::insertVertexInsideTriangle(
1429 const TriInd iNewT1 = addTriangle();
1430 const TriInd iNewT2 = addTriangle();
1432#ifdef CDT_ENABLE_CALLBACK_HANDLER
1433 if(m_callbackHandler)
1435 m_callbackHandler->onInsertVertexInsideTriangle(iT, iNewT1, iNewT2);
1439 Triangle& t = triangles[iT];
1440 const array<VertInd, 3> vv = t.vertices;
1441 const array<TriInd, 3> nn = t.neighbors;
1442 const VertInd v1 = vv[0], v2 = vv[1], v3 = vv[2];
1443 const TriInd n1 = nn[0], n2 = nn[1], n3 = nn[2];
1446 triangles[iNewT1] = Triangle(arr3(v2, v3, v), arr3(n2, iNewT2, iT));
1447 triangles[iNewT2] = Triangle(arr3(v3, v1, v), arr3(n3, iT, iNewT1));
1448 t = Triangle(arr3(v1, v2, v), arr3(n1, iNewT1, iNewT2));
1450 setAdjacentTriangle(v, iT);
1451 setAdjacentTriangle(v3, iNewT1);
1453 changeNeighbor(n2, iT, iNewT1);
1454 changeNeighbor(n3, iT, iNewT2);
1456 std::stack<TriInd> newTriangles;
1457 newTriangles.push(iT);
1458 newTriangles.push(iNewT1);
1459 newTriangles.push(iNewT2);
1460 return newTriangles;
1476template <
typename T,
typename TNearPo
intLocator>
1477std::stack<TriInd> Triangulation<T, TNearPointLocator>::insertVertexOnEdge(
1482 const TriInd iTnew1 = addTriangle();
1483 const TriInd iTnew2 = addTriangle();
1485#ifdef CDT_ENABLE_CALLBACK_HANDLER
1486 if(m_callbackHandler)
1488 m_callbackHandler->onInsertVertexOnEdge(iT1, iT2, iTnew1, iTnew2);
1492 Triangle& t1 = triangles[iT1];
1493 Triangle& t2 = triangles[iT2];
1495 const VertInd v1 = t1.vertices[i];
1497 const TriInd n1 = t1.neighbors[i];
1498 const TriInd n4 = t1.neighbors[
cw(i)];
1500 const VertInd v3 = t2.vertices[i];
1502 const TriInd n3 = t2.neighbors[i];
1503 const TriInd n2 = t2.neighbors[
cw(i)];
1505 t1 = Triangle(
arr3(v, v1, v2),
arr3(iTnew1, n1, iT2));
1506 t2 = Triangle(
arr3(v, v2, v3),
arr3(iT1, n2, iTnew2));
1507 triangles[iTnew1] = Triangle(
arr3(v, v4, v1),
arr3(iTnew2, n4, iT1));
1508 triangles[iTnew2] = Triangle(
arr3(v, v3, v4),
arr3(iT2, n3, iTnew1));
1510 setAdjacentTriangle(v, iT1);
1511 setAdjacentTriangle(v4, iTnew1);
1513 changeNeighbor(n4, iT1, iTnew1);
1514 changeNeighbor(n3, iT2, iTnew2);
1516 std::stack<TriInd> newTriangles;
1517 newTriangles.push(iT1);
1518 newTriangles.push(iTnew2);
1519 newTriangles.push(iT2);
1520 newTriangles.push(iTnew1);
1521 return newTriangles;
1524template <
typename T,
typename TNearPo
intLocator>
1526Triangulation<T, TNearPointLocator>::trianglesAt(
const V2d<T>& pos)
const
1528 array<TriInd, 2> out = {noNeighbor, noNeighbor};
1529 for(TriInd i =
TriInd(0); i <
TriInd(triangles.size()); ++i)
1531 const Triangle& t = triangles[i];
1532 const V2d<T>& v1 = vertices[t.vertices[0]];
1533 const V2d<T>& v2 = vertices[t.vertices[1]];
1534 const V2d<T>& v3 = vertices[t.vertices[2]];
1536 if(loc == PtTriLocation::Outside)
1544 Error(
"No triangle was found at position", CDT_SOURCE_LOCATION));
1548template <
typename T,
typename TNearPo
intLocator>
1549TriInd Triangulation<T, TNearPointLocator>::walkTriangles(
1550 const VertInd startVertex,
1551 const V2d<T>& pos)
const
1554 TriInd currTri = m_vertTris[startVertex];
1556 detail::SplitMix64RandGen prng;
1559 const Triangle& t = triangles[currTri];
1562 const Index offset(prng() % 3);
1563 for(Index i_(0); i_ <
Index(3); ++i_)
1565 const Index i((i_ + offset) % 3);
1566 const V2d<T>& vStart = vertices[t.vertices[i]];
1567 const V2d<T>& vEnd = vertices[t.vertices[
ccw(i)]];
1568 const PtLineLocation::Enum edgeCheck =
1570 const TriInd iN = t.neighbors[i];
1571 if(edgeCheck == PtLineLocation::Right && iN != noNeighbor)
1582template <
typename T,
typename TNearPo
intLocator>
1583array<TriInd, 2> Triangulation<T, TNearPointLocator>::walkingSearchTrianglesAt(
1585 const VertInd startVertex)
const
1587 const V2d<T> v = vertices[iV];
1588 array<TriInd, 2> out = {noNeighbor, noNeighbor};
1589 const TriInd iT = walkTriangles(startVertex, v);
1591 const Triangle& t = triangles[iT];
1592 const V2d<T>& v1 = vertices[t.vertices[0]];
1593 const V2d<T>& v2 = vertices[t.vertices[1]];
1594 const V2d<T>& v3 = vertices[t.vertices[2]];
1597 if(loc == PtTriLocation::Outside)
1600 Error(
"No triangle was found at position", CDT_SOURCE_LOCATION));
1602 if(loc == PtTriLocation::OnVertex)
1604 const VertInd iDupe = v1 == v ? t.vertices[0]
1605 : v2 == v ? t.vertices[1]
1607 handleException(DuplicateVertexError(
1609 VertInd(iDupe - m_nTargetVerts),
1610 CDT_SOURCE_LOCATION));
1619template <
typename T,
typename TNearPo
intLocator>
1620void Triangulation<T, TNearPointLocator>::changeNeighbor(
1622 const TriInd oldNeighbor,
1623 const TriInd newNeighbor)
1625 if(iT == noNeighbor)
1629 nn[0] == oldNeighbor || nn[1] == oldNeighbor || nn[2] == oldNeighbor);
1630 if(nn[0] == oldNeighbor)
1631 nn[0] = newNeighbor;
1632 else if(nn[1] == oldNeighbor)
1633 nn[1] = newNeighbor;
1635 nn[2] = newNeighbor;
1638template <
typename T,
typename TNearPo
intLocator>
1639void Triangulation<T, TNearPointLocator>::changeNeighbor(
1641 const VertInd iVedge1,
1642 const VertInd iVedge2,
1643 const TriInd newNeighbor)
1645 assert(iT != noNeighbor);
1646 Triangle& t = triangles[iT];
1647 t.neighbors[
edgeNeighborInd(t.vertices, iVedge1, iVedge2)] = newNeighbor;
1650template <
typename T,
typename TNearPo
intLocator>
1651void Triangulation<T, TNearPointLocator>::triangulatePseudoPolygon(
1652 const std::vector<VertInd>& poly,
1653 unordered_map<Edge, TriInd>& outerTris,
1656 std::vector<TriInd>& trianglesToReuse,
1657 std::vector<TriangulatePseudoPolygonTask>& iterations)
1659 assert(poly.size() > 2);
1662 iterations.push_back(make_tuple(
1664 static_cast<IndexSizeType
>(poly.size() - 1),
1668 while(!iterations.empty())
1670 triangulatePseudoPolygonIteration(
1671 poly, outerTris, trianglesToReuse, iterations);
1675template <
typename T,
typename TNearPo
intLocator>
1676void Triangulation<T, TNearPointLocator>::triangulatePseudoPolygonIteration(
1677 const std::vector<VertInd>& poly,
1678 unordered_map<Edge, TriInd>& outerTris,
1679 std::vector<TriInd>& trianglesToReuse,
1680 std::vector<TriangulatePseudoPolygonTask>& iterations)
1682 IndexSizeType iA, iB;
1685 assert(!iterations.empty());
1686 tie(iA, iB, iT, iParent, iInParent) = iterations.back();
1687 iterations.pop_back();
1688 assert(iB - iA > 1 && iT != noNeighbor && iParent != noNeighbor);
1689 Triangle& t = triangles[iT];
1691 const IndexSizeType iC = findDelaunayPoint(poly, iA, iB);
1703 assert(!trianglesToReuse.empty());
1704 const TriInd iNext = trianglesToReuse.back();
1705 trianglesToReuse.pop_back();
1706 iterations.push_back(make_tuple(iC, iB, iNext, iT,
Index(1)));
1710 const Edge outerEdge(b, c);
1711 const TriInd outerTri = outerTris.at(outerEdge);
1712 if(outerTri != noNeighbor)
1714 assert(outerTri != iT);
1715 t.neighbors[1] = outerTri;
1716 changeNeighbor(outerTri, c, b, iT);
1719 outerTris.at(outerEdge) = iT;
1724 assert(!trianglesToReuse.empty());
1725 const TriInd iNext = trianglesToReuse.back();
1726 trianglesToReuse.pop_back();
1727 iterations.push_back(make_tuple(iA, iC, iNext, iT,
Index(2)));
1731 const Edge outerEdge(c, a);
1732 const TriInd outerTri = outerTris.at(outerEdge);
1733 if(outerTri != noNeighbor)
1735 assert(outerTri != iT);
1736 t.neighbors[2] = outerTri;
1737 changeNeighbor(outerTri, c, a, iT);
1740 outerTris.at(outerEdge) = iT;
1745 triangles[iParent].neighbors[iInParent] = iT;
1746 t.neighbors[0] = iParent;
1747 t.vertices =
arr3(a, b, c);
1748 setAdjacentTriangle(c, iT);
1751template <
typename T,
typename TNearPo
intLocator>
1752IndexSizeType Triangulation<T, TNearPointLocator>::findDelaunayPoint(
1753 const std::vector<VertInd>& poly,
1754 const IndexSizeType iA,
1755 const IndexSizeType iB)
const
1757 assert(iB - iA > 1);
1758 const V2d<T>& a = vertices[poly[iA]];
1759 const V2d<T>& b = vertices[poly[iB]];
1760 IndexSizeType out = iA + 1;
1761 const V2d<T>* c = &vertices[poly[out]];
1762 for(IndexSizeType i = iA + 1; i < iB; ++i)
1764 const V2d<T>& v = vertices[poly[i]];
1771 assert(out > iA && out < iB);
1775template <
typename T,
typename TNearPo
intLocator>
1777 const std::vector<
V2d<T> >& newVertices)
1783template <
typename T,
typename TNearPo
intLocator>
1786 return m_vertTris.empty() && !
vertices.empty();
1789template <
typename T,
typename TNearPo
intLocator>
1790unordered_map<TriInd, LayerDepth>
1791Triangulation<T, TNearPointLocator>::peelLayer(
1792 std::stack<TriInd> seeds,
1794 std::vector<LayerDepth>& triDepths)
const
1796 unordered_map<TriInd, LayerDepth> behindBoundary;
1797 while(!seeds.empty())
1799 const TriInd iT = seeds.top();
1801 triDepths[iT] = std::min(triDepths[iT], layerDepth);
1802 behindBoundary.erase(iT);
1808 if(iN == noNeighbor || triDepths[iN] <= layerDepth)
1810 if(fixedEdges.count(opEdge))
1812 const unordered_map<Edge, LayerDepth>::const_iterator cit =
1813 overlapCount.find(opEdge);
1814 const LayerDepth triDepth = cit == overlapCount.end()
1816 : layerDepth + cit->second + 1;
1817 behindBoundary[iN] = triDepth;
1823 return behindBoundary;
1826template <
typename T,
typename TNearPo
intLocator>
1827std::vector<LayerDepth>
1830 std::vector<LayerDepth> triDepths(
1831 triangles.size(), std::numeric_limits<LayerDepth>::max());
1832 std::stack<TriInd> seeds(TriDeque(1, m_vertTris[0]));
1836 unordered_map<LayerDepth, TriIndUSet> seedsByDepth;
1839 const unordered_map<TriInd, LayerDepth>& newSeeds =
1840 peelLayer(seeds, layerDepth, triDepths);
1842 seedsByDepth.erase(layerDepth);
1843 typedef unordered_map<TriInd, LayerDepth>::const_iterator Iter;
1844 for(Iter it = newSeeds.begin(); it != newSeeds.end(); ++it)
1846 deepestSeedDepth = std::max(deepestSeedDepth, it->second);
1847 seedsByDepth[it->second].insert(it->first);
1849 const TriIndUSet& nextLayerSeeds = seedsByDepth[layerDepth + 1];
1850 seeds = std::stack<TriInd>(
1851 TriDeque(nextLayerSeeds.begin(), nextLayerSeeds.end()));
1853 }
while(!seeds.empty() || deepestSeedDepth > layerDepth);
1858#ifdef CDT_ENABLE_CALLBACK_HANDLER
1859template <
typename T,
typename TNearPo
intLocator>
1863 m_callbackHandler = callbackHandler;
1867template <
typename T,
typename TNearPo
intLocator>
1868void Triangulation<T, TNearPointLocator>::insertVertices_AsProvided(
1871 for(
VertInd iV = superGeomVertCount; iV < vertices.size(); ++iV)
1873#ifdef CDT_ENABLE_CALLBACK_HANDLER
1874 if(m_callbackHandler && m_callbackHandler->isAbortCalculation())
1883template <
typename T,
typename TNearPo
intLocator>
1884void Triangulation<T, TNearPointLocator>::insertVertices_Randomized(
1885 VertInd superGeomVertCount)
1887 std::size_t vertexCount = vertices.size() - superGeomVertCount;
1888 std::vector<VertInd> ii(vertexCount);
1889 detail::iota(ii.begin(), ii.end(), superGeomVertCount);
1890 detail::random_shuffle(ii.begin(), ii.end());
1891 for(std::vector<VertInd>::iterator it = ii.begin(); it != ii.end(); ++it)
1893#ifdef CDT_ENABLE_CALLBACK_HANDLER
1894 if(m_callbackHandler && m_callbackHandler->isAbortCalculation())
1907template <
typename T>
1908inline double log2_bc(T x)
1910#ifdef CDT_CXX11_IS_SUPPORTED
1911 return std::log2(x);
1913 static double log2_constant = std::log(2.0);
1914 return std::log(
static_cast<double>(x)) / log2_constant;
1924 const int filledLayerPow2 =
1925 static_cast<int>(std::floor(log2_bc(vertexCount)) - 1);
1926 const std::size_t nodesInFilledTree =
1927 static_cast<std::size_t
>(std::pow(2., filledLayerPow2 + 1) - 1);
1928 const std::size_t nodesInLastFilledLayer =
1929 static_cast<std::size_t
>(std::pow(2., filledLayerPow2));
1930 const std::size_t nodesInLastLayer = vertexCount - nodesInFilledTree;
1931 return nodesInLastLayer >= nodesInLastFilledLayer
1932 ? nodesInLastFilledLayer + nodesInLastLayer -
1933 nodesInLastFilledLayer
1934 : nodesInLastFilledLayer;
1937template <
typename T>
1938class FixedCapacityQueue
1941 FixedCapacityQueue(
const std::size_t capacity)
1943 , m_front(m_vec.begin())
1944 , m_back(m_vec.begin())
1951 const T& front()
const
1959 if(m_front == m_vec.end())
1960 m_front = m_vec.begin();
1963 void push(
const T& t)
1965 assert(m_size < m_vec.size());
1968 if(m_back == m_vec.end())
1969 m_back = m_vec.begin();
1972#ifdef CDT_CXX11_IS_SUPPORTED
1973 void push(
const T&& t)
1975 assert(m_size < m_vec.size());
1978 if(m_back == m_vec.end())
1979 m_back = m_vec.begin();
1984 std::vector<T> m_vec;
1985 typename std::vector<T>::iterator m_front;
1986 typename std::vector<T>::iterator m_back;
1990template <
typename T>
1993 const std::vector<V2d<T> >& m_vertices;
1996 less_than_x(
const std::vector<
V2d<T> >& vertices)
1997 : m_vertices(vertices)
2001 return m_vertices[a].x < m_vertices[b].x;
2005template <
typename T>
2008 const std::vector<V2d<T> >& m_vertices;
2011 less_than_y(
const std::vector<
V2d<T> >& vertices)
2012 : m_vertices(vertices)
2016 return m_vertices[a].y < m_vertices[b].y;
2022template <
typename T,
typename TNearPo
intLocator>
2023void Triangulation<T, TNearPointLocator>::insertVertices_KDTreeBFS(
2028 const VertInd vertexCount(vertices.size() - superGeomVertCount);
2031 std::vector<VertInd> ii(vertexCount);
2032 detail::iota(ii.begin(), ii.end(), superGeomVertCount);
2034 typedef std::vector<VertInd>::iterator It;
2037 queue.push(make_tuple(ii.begin(), ii.end(), box.
min, box.
max,
VertInd(0)));
2040 V2d<T> newBoxMin, newBoxMax;
2046 while(!queue.empty())
2048#ifdef CDT_ENABLE_CALLBACK_HANDLER
2049 if(m_callbackHandler && m_callbackHandler->isAbortCalculation())
2054 tie(first, last, box.
min, box.
max, parent) = queue.front();
2056 assert(first != last);
2058 const std::ptrdiff_t len = std::distance(first, last);
2061 insertVertex(*first, parent);
2064 const It midIt = first + len / 2;
2067 detail::portable_nth_element(first, midIt, last, cmpX);
2069 const T split = vertices[mid].x;
2070 newBoxMin.
x = split;
2071 newBoxMin.
y = box.
min.y;
2072 newBoxMax.
x = split;
2073 newBoxMax.
y = box.
max.y;
2077 detail::portable_nth_element(first, midIt, last, cmpY);
2079 const T split = vertices[mid].y;
2080 newBoxMin.
x = box.
min.x;
2081 newBoxMin.
y = split;
2082 newBoxMax.
x = box.
max.x;
2083 newBoxMax.
y = split;
2085 insertVertex(mid, parent);
2088 queue.push(make_tuple(first, midIt, box.
min, newBoxMax, mid));
2090 if(midIt + 1 != last)
2092 queue.push(make_tuple(midIt + 1, last, newBoxMin, box.
max, mid));
2097template <
typename T,
typename TNearPo
intLocator>
2098std::pair<TriInd, TriInd> Triangulation<T, TNearPointLocator>::edgeTriangles(
2100 const VertInd b)
const
2102 const TriInd triStart = m_vertTris[a];
2103 assert(triStart != noNeighbor);
2104 TriInd iT = triStart, iTNext = triStart;
2108 const Triangle& t = triangles[iT];
2109 tie(iTNext, iV) = t.
next(a);
2110 assert(iTNext != noNeighbor);
2113 return std::make_pair(iT, iTNext);
2116 }
while(iT != triStart);
2117 return std::make_pair(noNeighbor, noNeighbor);
2120template <
typename T,
typename TNearPo
intLocator>
2121bool Triangulation<T, TNearPointLocator>::hasEdge(
2123 const VertInd b)
const
2125 return edgeTriangles(a, b).first != invalidIndexSizeType;
2128template <
typename T,
typename TNearPo
intLocator>
2129void Triangulation<T, TNearPointLocator>::setAdjacentTriangle(
2133 assert(t != noNeighbor);
2136 triangles[t].vertices[0] == v || triangles[t].vertices[1] == v ||
2137 triangles[t].vertices[2] == v);
2140template <
typename T,
typename TNearPo
intLocator>
2141void Triangulation<T, TNearPointLocator>::pivotVertexTriangleCW(
const VertInd v)
2143 assert(m_vertTris[v] != noNeighbor);
2144 m_vertTris[v] = triangles[m_vertTris[v]].next(v).first;
2145 assert(m_vertTris[v] != noNeighbor);
2147 triangles[m_vertTris[v]].vertices[0] == v ||
2148 triangles[m_vertTris[v]].vertices[1] == v ||
2149 triangles[m_vertTris[v]].vertices[2] == v);
2152template <
typename T,
typename TNearPo
intLocator>
2153void Triangulation<T, TNearPointLocator>::tryAddVertexToLocator(
const VertInd v)
2155 if(!m_nearPtLocator.empty())
2156 m_nearPtLocator.addPoint(v, vertices);
2159template <
typename T,
typename TNearPo
intLocator>
2160void Triangulation<T, TNearPointLocator>::tryInitNearestPointLocator()
2162 if(!vertices.empty() && m_nearPtLocator.empty())
2164 m_nearPtLocator.initialize(vertices);
void iota(ForwardIt first, ForwardIt last, T value)
backport from c++11
std::size_t maxQueueLengthBFSKDTree(const std::size_t vertexCount)
Since KD-tree bulk load builds a balanced tree the maximum length of a queue can be pre-calculated: i...
Interface for the callback handler that user can derive from and inject into the triangulation to mon...
Error thrown when intersecting constraint edges are detected, but triangulation is not configured to ...
void eraseOuterTriangles()
Erase triangles outside of constrained boundary using growing.
void conformToEdges(TEdgeIter first, TEdgeIter last, TGetEdgeVertexStart getStart, TGetEdgeVertexEnd getEnd)
Insert constraint edges into triangulation for Conforming Delaunay Triangulation (for example see fig...
std::vector< LayerDepth > calculateTriangleDepths() const
Calculate depth of each triangle in constraint triangulation.
bool isFinalized() const
Check if the triangulation was finalized with erase... method and super-triangle was removed.
void eraseOuterTrianglesAndHoles()
Erase triangles outside of constrained boundary and auto-detected holes.
V2dVec vertices
triangulation's vertices
void insertEdges(TEdgeIter first, TEdgeIter last, TGetEdgeVertexStart getStart, TGetEdgeVertexEnd getEnd)
Insert constraint edges into triangulation for Constrained Delaunay Triangulation (for example see fi...
void insertVertices(TVertexIter first, TVertexIter last, TGetVertexCoordX getX, TGetVertexCoordY getY)
Insert custom point-types specified by iterator range and X/Y-getters.
void eraseSuperTriangle()
Erase triangles adjacent to super triangle.
TriangleVec triangles
triangulation's triangles
Triangulation()
Default constructor.
void setCallbackHandler(ICallbackHandler *callbackHandler)
Set user-provided callback handler.
void initializedWithCustomSuperGeometry()
Call this method after directly setting custom super-geometry via vertices and triangles members.
unsigned short LayerDepth
Type used for storing layer depths for triangles.
void flipEdge(TriInd iT, TriInd iTopo)
Flip an edge between two triangle.
TriIndVec & VertTrisInternal()
Access internal vertex adjacent triangles.
void removeTriangles(const TriIndUSet &removedTriangles)
Remove triangles with specified indices.
Namespace containing triangulation functionality.
std::vector< Edge > EdgeVec
Vector of edges.
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.
unordered_set< Edge > EdgeUSet
Hash table of edges.
VertInd edge_get_v2(const Edge &e)
Get edge second vertex.
std::vector< TriInd > TriIndVec
Vector of triangle indices.
CDT_INLINE_IF_HEADER_ONLY Index edgeNeighborInd(const VerticesArr3 &vv, VertInd iVedge1, VertInd iVedge2)
Index of triangle's neighbor opposed to an edge.
VertInd edge_get_v1(const Edge &e)
Get edge first vertex.
CDT_EXPORT PtLineLocation::Enum classifyOrientation(T orientation, T orientationTolerance=T(0))
Classify value of orient2d predicate.
array< TriInd, 3 > NeighborsArr3
array of three neighbors
CDT_EXPORT T distance(const V2d< T > &a, const V2d< T > &b)
Distance between two 2D points.
IndexSizeType VertInd
Vertex index.
CDT_EXPORT Index cw(Index i)
Advance vertex or neighbor index clockwise.
array< VertInd, 3 > VerticesArr3
array of three vertex indices
CDT_INLINE_IF_HEADER_ONLY bool touchesSuperTriangle(const Triangle &t)
Check if any of triangle's vertices belongs to a super-triangle.
array< T, 3 > arr3(const T &v0, const T &v1, const T &v2)
Needed for c++03 compatibility (no uniform initialization available)
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opoNbr(Index vertIndex)
Opposed neighbor index from vertex index.
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index vertexInd(const VerticesArr3 &vv, VertInd iV)
If triangle has a given vertex return vertex-index.
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.
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.
unordered_set< TriInd > TriIndUSet
Hash table of triangles.
CDT_EXPORT Index ccw(Index i)
Advance vertex or neighbor index counter-clockwise.
CDT_EXPORT Index edgeNeighbor(PtTriLocation::Enum location)
Neighbor index from a on-edge location.
unordered_map< TriInd, TriInd > TriIndUMap
Triangle hash map.
const T & getX_V2d(const V2d< T > &v)
X- coordinate getter for V2d.
unsigned char Index
Index in triangle.
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY TriInd opposedTriangle(const Triangle &tri, VertInd iVert)
Given triangle and a vertex find opposed triangle.
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.
Edge RemapNoSuperTriangle(const Edge &e)
Remap removing super-triangle: subtract 3 from vertices.
CDT_EXPORT bool isOnEdge(PtTriLocation::Enum location)
Check if location is classified as on any of three edges.
IndexSizeType TriInd
Triangle index.
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.
CDT_EXPORT CDT_INLINE_IF_HEADER_ONLY Index opposedVertexInd(const NeighborsArr3 &nn, TriInd iTopo)
Index of triangle's vertex opposed to a triangle.
const T & getY_V2d(const V2d< T > &v)
Y-coordinate getter for V2d.
@ FixedEdgeMidpoint
During conforming triangulation edge mid-point is added.
@ FixedEdgesIntersection
Resolving fixed/constraint edges' intersection.
V2d< T > max
max box corner
V2d< T > min
min box corner
Edge connecting two vertices: vertex with smaller index is always first.
VertInd v1() const
V1 getter.
VertInd v2() const
V2 getter.
@ TryResolve
attempt to resolve constraint edge intersections
@ NotAllowed
constraint edge intersections are not allowed
@ DontCheck
No checks: slightly faster but less safe.
@ SuperTriangle
conventional super-triangle
@ Custom
user-specified custom geometry (e.g., grid)
Triangulation triangle (counter-clockwise winding)
VerticesArr3 vertices
triangle's three vertices
NeighborsArr3 neighbors
triangle's three neighbors
std::pair< TriInd, VertInd > next(const VertInd i) const
Next triangle adjacent to a vertex (clockwise)
@ Auto
Automatic insertion order optimized for better performance.