Triangulation.h

Go to the documentation of this file.

/* This Source Code Form is subject to the terms of the Mozilla Public
 * License, v. 2.0. If a copy of the MPL was not distributed with this
 * file, You can obtain one at https://mozilla.org/MPL/2.0/. */
 
#ifndef CDT_vW1vZ0lO8rS4gY4uI4fB
#define CDT_vW1vZ0lO8rS4gY4uI4fB
 
#include "CDTUtils.h"
#include "LocatorKDTree.h"
 
#include <algorithm>
#include <cstdlib>
#include <iterator>
#include <stack>
#include <stdexcept>
#include <utility>
#include <vector>
 
namespace CDT
{
 
 
struct CDT_EXPORT VertexInsertionOrder
{
    enum Enum
    {
        Auto,
        AsProvided,
    };
};
 
struct CDT_EXPORT SuperGeometryType
{
    enum Enum
    {
        SuperTriangle, 
        Custom,        
    };
};
 
struct CDT_EXPORT IntersectingConstraintEdges
{
    enum Enum
    {
        Ignore,  
        Resolve, 
    };
};
 
typedef unsigned short LayerDepth;
typedef LayerDepth BoundaryOverlapCount;
 
 
template <typename T, typename TNearPointLocator = LocatorKDTree<T> >
class CDT_EXPORT Triangulation
{
public:
    typedef std::vector<V2d<T> > V2dVec; 
    V2dVec vertices;                     
    TriangleVec triangles;               
    EdgeUSet fixedEdges; 
 
    unordered_map<Edge, BoundaryOverlapCount> overlapCount;
 
    unordered_map<Edge, EdgeVec> pieceToOriginals;
 
    /*____ API _____*/
    Triangulation();
    explicit Triangulation(VertexInsertionOrder::Enum vertexInsertionOrder);
    Triangulation(
        VertexInsertionOrder::Enum vertexInsertionOrder,
        IntersectingConstraintEdges::Enum intersectingEdgesStrategy,
        T minDistToConstraintEdge);
    Triangulation(
        VertexInsertionOrder::Enum vertexInsertionOrder,
        const TNearPointLocator& nearPtLocator,
        IntersectingConstraintEdges::Enum intersectingEdgesStrategy,
        T minDistToConstraintEdge);
    template <
        typename TVertexIter,
        typename TGetVertexCoordX,
        typename TGetVertexCoordY>
    void insertVertices(
        TVertexIter first,
        TVertexIter last,
        TGetVertexCoordX getX,
        TGetVertexCoordY getY);
    void insertVertices(const std::vector<V2d<T> >& vertices);
    template <
        typename TEdgeIter,
        typename TGetEdgeVertexStart,
        typename TGetEdgeVertexEnd>
    void insertEdges(
        TEdgeIter first,
        TEdgeIter last,
        TGetEdgeVertexStart getStart,
        TGetEdgeVertexEnd getEnd);
    void insertEdges(const std::vector<Edge>& edges);
    template <
        typename TEdgeIter,
        typename TGetEdgeVertexStart,
        typename TGetEdgeVertexEnd>
    void conformToEdges(
        TEdgeIter first,
        TEdgeIter last,
        TGetEdgeVertexStart getStart,
        TGetEdgeVertexEnd getEnd);
    void conformToEdges(const std::vector<Edge>& edges);
    void eraseSuperTriangle();
    void eraseOuterTriangles();
    void eraseOuterTrianglesAndHoles();
    void initializedWithCustomSuperGeometry();
 
    bool isFinalized() const;
 
    std::vector<LayerDepth> calculateTriangleDepths() const;
 
 
    void flipEdge(TriInd iT, TriInd iTopo);
 
    void flipEdge(
        TriInd iT,
        TriInd iTopo,
        VertInd v1,
        VertInd v2,
        VertInd v3,
        VertInd v4,
        TriInd n1,
        TriInd n2,
        TriInd n3,
        TriInd n4);
 
    void removeTriangles(const TriIndUSet& removedTriangles);
 
    TriIndVec& VertTrisInternal();
 
private:
    /*____ Detail __*/
    void addSuperTriangle(const Box2d<T>& box);
    void addNewVertex(const V2d<T>& pos, TriInd iT);
    void insertVertex(VertInd iVert);
    void insertVertex(VertInd iVert, VertInd walkStart);
    void ensureDelaunayByEdgeFlips(
        const V2d<T>& v1,
        VertInd iV1,
        std::stack<TriInd>& triStack);
    std::vector<Edge> insertVertex_FlipFixedEdges(VertInd iV1);
 
    typedef tuple<IndexSizeType, IndexSizeType, TriInd, TriInd, Index>
        TriangulatePseudopolygonTask;
 
    void insertEdge(
        Edge edge,
        Edge originalEdge,
        EdgeVec& remaining,
        std::vector<TriangulatePseudopolygonTask>& tppIterations);
 
    void insertEdgeIteration(
        Edge edge,
        Edge originalEdge,
        EdgeVec& remaining,
        std::vector<TriangulatePseudopolygonTask>& tppIterations);
 
    IndexSizeType previousEdgeEncounter(
        const VertInd v,
        const std::vector<VertInd>& poly) const;
 
    typedef tuple<Edge, EdgeVec, BoundaryOverlapCount> ConformToEdgeTask;
 
    void conformToEdge(
        Edge edge,
        EdgeVec originals,
        BoundaryOverlapCount overlaps,
        std::vector<ConformToEdgeTask>& remaining);
 
    void conformToEdgeIteration(
        Edge edge,
        const EdgeVec& originals,
        BoundaryOverlapCount overlaps,
        std::vector<ConformToEdgeTask>& remaining);
 
    tuple<TriInd, VertInd, VertInd> intersectedTriangle(
        VertInd iA,
        const V2d<T>& a,
        const V2d<T>& b,
        T orientationTolerance = T(0)) const;
    std::stack<TriInd> insertVertexInsideTriangle(VertInd v, TriInd iT);
    std::stack<TriInd> insertVertexOnEdge(VertInd v, TriInd iT1, TriInd iT2);
    array<TriInd, 2> trianglesAt(const V2d<T>& pos) const;
    array<TriInd, 2>
    walkingSearchTrianglesAt(const V2d<T>& pos, VertInd startVertex) const;
    TriInd walkTriangles(VertInd startVertex, const V2d<T>& pos) const;
    void edgeFlipInfo(
        TriInd iT,
        VertInd iV1,
        TriInd& iTopo,
        VertInd& iV2,
        VertInd& iV3,
        VertInd& iV4,
        TriInd& n1,
        TriInd& n2,
        TriInd& n3,
        TriInd& n4);
    bool isFlipNeeded(
        const V2d<T>& v,
        VertInd iV1,
        VertInd iV2,
        VertInd iV3,
        VertInd iV4) const;
    void changeNeighbor(TriInd iT, TriInd oldNeighbor, TriInd newNeighbor);
    void changeNeighbor(
        TriInd iT,
        VertInd iVedge1,
        VertInd iVedge2,
        TriInd newNeighbor);
    void triangulatePseudopolygon(
        const std::vector<VertInd>& poly,
        const std::vector<TriInd>& outerTris,
        TriInd iT,
        TriInd iN,
        std::vector<TriangulatePseudopolygonTask>& iterations);
    void triangulatePseudopolygonIteration(
        const std::vector<VertInd>& poly,
        const std::vector<TriInd>& outerTris,
        std::vector<TriangulatePseudopolygonTask>& iterations);
    IndexSizeType findDelaunayPoint(
        const std::vector<VertInd>& poly,
        IndexSizeType iA,
        IndexSizeType iB) const;
    TriInd addTriangle(const Triangle& t); // note: invalidates iterators!
    TriInd addTriangle(); // note: invalidates triangle iterators!
    void finalizeTriangulation(const TriIndUSet& removedTriangles);
    TriIndUSet growToBoundary(std::stack<TriInd> seeds) const;
    void fixEdge(const Edge& edge, BoundaryOverlapCount overlaps);
    void fixEdge(const Edge& edge);
    void fixEdge(const Edge& edge, const Edge& originalEdge);
    void makeDummy(TriInd iT);
    void eraseDummies();
    unordered_map<TriInd, LayerDepth> peelLayer(
        std::stack<TriInd> seeds,
        LayerDepth layerDepth,
        std::vector<LayerDepth>& triDepths) const;
 
    void insertVertices_AsProvided(VertInd superGeomVertCount);
    void insertVertices_Randomized(VertInd superGeomVertCount);
    void insertVertices_KDTreeBFS(
        VertInd superGeomVertCount,
        V2d<T> boxMin,
        V2d<T> boxMax);
    bool hasEdge(VertInd a, VertInd b) const;
    void setAdjacentTriangle(const VertInd v, const TriInd t);
    void pivotVertexTriangleCW(VertInd v);
    void removeAdjacentTriangle(VertInd v);
    void tryAddVertexToLocator(const VertInd v);
    void tryInitNearestPointLocator();
 
    std::vector<TriInd> m_dummyTris;
    TNearPointLocator m_nearPtLocator;
    std::size_t m_nTargetVerts;
    SuperGeometryType::Enum m_superGeomType;
    VertexInsertionOrder::Enum m_vertexInsertionOrder;
    IntersectingConstraintEdges::Enum m_intersectingEdgesStrategy;
    T m_minDistToConstraintEdge;
    TriIndVec m_vertTris; 
};
 
 
namespace detail
{
 
struct SplitMix64RandGen
{
    typedef unsigned long long uint64;
    uint64 m_state;
    explicit SplitMix64RandGen(uint64 state)
        : m_state(state)
    {}
    explicit SplitMix64RandGen()
        : m_state(0)
    {}
    uint64 operator()()
    {
        uint64 z = (m_state += 0x9e3779b97f4a7c15);
        z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9;
        z = (z ^ (z >> 27)) * 0x94d049bb133111eb;
        return z ^ (z >> 31);
    }
};
 
template <class RandomIt>
void random_shuffle(RandomIt first, RandomIt last)
{
    detail::SplitMix64RandGen prng;
    typename std::iterator_traits<RandomIt>::difference_type i, n;
    n = last - first;
    for(i = n - 1; i > 0; --i)
    {
        std::swap(first[i], first[prng() % (i + 1)]);
    }
}
 
// backport from c++11
template <class ForwardIt, class T>
void iota(ForwardIt first, ForwardIt last, T value)
{
    while(first != last)
    {
        *first++ = value;
        ++value;
    }
}
 
} // namespace detail
 
//-----------------------
// Triangulation methods
//-----------------------
template <typename T, typename TNearPointLocator>
template <
    typename TVertexIter,
    typename TGetVertexCoordX,
    typename TGetVertexCoordY>
void Triangulation<T, TNearPointLocator>::insertVertices(
    const TVertexIter first,
    const TVertexIter last,
    TGetVertexCoordX getX,
    TGetVertexCoordY getY)
{
    if(isFinalized())
    {
        throw std::runtime_error(
            "Triangulation was finalized with 'erase...' method. Inserting new "
            "vertices is not possible");
    }
 
    const bool isFirstTime = vertices.empty();
    Box2d<T> box;
    if(vertices.empty()) // called first time
    {
        box = envelopBox<T>(first, last, getX, getY);
        addSuperTriangle(box);
    }
    tryInitNearestPointLocator();
 
    const VertInd nExistingVerts(vertices.size());
    const VertInd nVerts(nExistingVerts + std::distance(first, last));
    // optimization, try to pre-allocate tris
    triangles.reserve(triangles.size() + 2 * nVerts);
    vertices.reserve(nVerts);
    m_vertTris.reserve(nVerts);
    for(TVertexIter it = first; it != last; ++it)
        addNewVertex(V2d<T>::make(getX(*it), getY(*it)), noNeighbor);
 
    switch(m_vertexInsertionOrder)
    {
    case VertexInsertionOrder::AsProvided:
        insertVertices_AsProvided(nExistingVerts);
        break;
    case VertexInsertionOrder::Auto:
        isFirstTime ? insertVertices_KDTreeBFS(nExistingVerts, box.min, box.max)
                    : insertVertices_Randomized(nExistingVerts);
        break;
    }
}
 
template <typename T, typename TNearPointLocator>
template <
    typename TEdgeIter,
    typename TGetEdgeVertexStart,
    typename TGetEdgeVertexEnd>
void Triangulation<T, TNearPointLocator>::insertEdges(
    TEdgeIter first,
    const TEdgeIter last,
    TGetEdgeVertexStart getStart,
    TGetEdgeVertexEnd getEnd)
{
    // state shared between different runs for performance gains
    std::vector<TriangulatePseudopolygonTask> tppIterations;
    EdgeVec remaining;
    if(isFinalized())
    {
        throw std::runtime_error(
            "Triangulation was finalized with 'erase...' method. Inserting new "
            "edges is not possible");
    }
    for(; first != last; ++first)
    {
        // +3 to account for super-triangle vertices
        const Edge edge(
            VertInd(getStart(*first) + m_nTargetVerts),
            VertInd(getEnd(*first) + m_nTargetVerts));
        insertEdge(edge, edge, remaining, tppIterations);
    }
    eraseDummies();
}
 
template <typename T, typename TNearPointLocator>
template <
    typename TEdgeIter,
    typename TGetEdgeVertexStart,
    typename TGetEdgeVertexEnd>
void Triangulation<T, TNearPointLocator>::conformToEdges(
    TEdgeIter first,
    const TEdgeIter last,
    TGetEdgeVertexStart getStart,
    TGetEdgeVertexEnd getEnd)
{
    if(isFinalized())
    {
        throw std::runtime_error(
            "Triangulation was finalized with 'erase...' method. Conforming to "
            "new edges is not possible");
    }
    tryInitNearestPointLocator();
    // state shared between different runs for performance gains
    std::vector<ConformToEdgeTask> remaining;
    for(; first != last; ++first)
    {
        // +3 to account for super-triangle vertices
        const Edge e(
            VertInd(getStart(*first) + m_nTargetVerts),
            VertInd(getEnd(*first) + m_nTargetVerts));
        conformToEdge(e, EdgeVec(1, e), 0, remaining);
    }
    eraseDummies();
}
 
} // namespace CDT
 
#ifndef CDT_USE_AS_COMPILED_LIBRARY
#include "Triangulation.hpp"
#endif
 
#endif // header-guard