Mesh2DIntersections.hpp

//---- GPL ---------------------------------------------------------------------
//
// Copyright (C)  Stichting Deltares, 2011-2023.
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation version 3.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program.  If not, see <http://www.gnu.org/licenses/>.
//
// contact: delft3d.support@deltares.nl
// Stichting Deltares
// P.O. Box 177
// 2600 MH Delft, The Netherlands
//
// All indications and logos of, and references to, "Delft3D" and "Deltares"
// are registered trademarks of Stichting Deltares, and remain the property of
// Stichting Deltares. All rights reserved.
//
//------------------------------------------------------------------------------
 
#pragma once
#include <algorithm>
#include <vector>
 
#include "MeshKernel/Definitions.hpp"
#include "MeshKernel/Mesh2D.hpp"
#include "MeshKernel/Polygons.hpp"
 
#include <queue>
 
namespace meshkernel
{
    struct EdgeMeshPolyLineIntersection
    {
        int polylineSegmentIndex{constants::missing::intValue};                      
        double polylineDistance{constants::missing::doubleValue};                    
        double adimensionalPolylineSegmentDistance{constants::missing::doubleValue}; 
        UInt edgeIndex{constants::missing::uintValue};                               
        UInt edgeFirstNode{constants::missing::uintValue};                           
        UInt edgeSecondNode{constants::missing::uintValue};                          
        double edgeDistance{constants::missing::doubleValue};                        
    };
 
    struct FaceMeshPolyLineIntersection
    {
        double polylineDistance{constants::missing::doubleValue}; 
        UInt faceIndex{constants::missing::uintValue};            
        std::vector<UInt> edgeIndices;                            
        std::vector<UInt> edgeNodes;                              
    };
 
    class Mesh2DIntersections final
    {
    public:
        explicit Mesh2DIntersections(Mesh2D& mesh);
 
        void Compute(const Polygons& polygon);
 
        void Compute(const std::vector<Point>& polyLine);
 
        [[nodiscard]] const auto& EdgeIntersections() const { return m_edgesIntersections; }
 
        [[nodiscard]] const auto& FaceIntersections() const { return m_faceIntersections; }
 
        template <typename T>
        static void sortAndEraseIntersections(std::vector<T>& intersections)
        {
            std::ranges::sort(intersections,
                              [](const T& first, const T& second)
                              { return first.polylineDistance < second.polylineDistance; });
 
            std::erase_if(intersections, [](const T& v)
                          { return v.polylineDistance < 0; });
        }
 
    private:
        enum class Direction
        {
            Forward,
            Backward
        };
 
        std::tuple<UInt, UInt> GetIntersectionSeed(const Mesh2D& mesh,
                                                   const std::vector<Point>& polyLine,
                                                   const std::vector<BoundingBox>& polyLineBoundingBoxes,
                                                   const std::vector<bool>& vistedEdges) const;
 
        std::tuple<bool, UInt, UInt, double, double, double> GetNextEdgeIntersection(const std::vector<Point>& polyLine,
                                                                                     const std::vector<BoundingBox>& polyLineBoundingBoxes,
                                                                                     UInt edgeIndex,
                                                                                     UInt firstIndex,
                                                                                     UInt secondIndex,
                                                                                     Direction direction) const;
 
        void IntersectFaceEdges(const std::vector<Point>& polyLine,
                                const std::vector<BoundingBox>& polyLineBoundingBoxes,
                                const std::vector<double>& cumulativeLength,
                                UInt currentCrossingEdge,
                                UInt currentFaceIndex,
                                UInt segmentIndex,
                                std::vector<bool>& vistedEdges,
                                std::vector<bool>& vistedFace,
                                std::queue<std::array<UInt, 2>>& crossingEdges);
 
        static void updateEdgeIntersections(const UInt segmentIndex,
                                            const UInt edgeIndex,
                                            const Edge edge,
                                            const std::vector<double>& cumulativeLength,
                                            const double crossProductValue,
                                            const double adimensionalEdgeDistance,
                                            const double adimensionalPolylineSegmentDistance,
                                            std::vector<EdgeMeshPolyLineIntersection>& intersections)
        {
            const auto [edgeFirstNode, edgeSecondNode] = edge;
 
            intersections[edgeIndex].polylineSegmentIndex = static_cast<int>(segmentIndex);
            intersections[edgeIndex].polylineDistance = cumulativeLength[segmentIndex] +
                                                        adimensionalPolylineSegmentDistance * (cumulativeLength[segmentIndex + 1] - cumulativeLength[segmentIndex]);
            intersections[edgeIndex].adimensionalPolylineSegmentDistance = adimensionalPolylineSegmentDistance;
            intersections[edgeIndex].edgeFirstNode = crossProductValue < 0 ? edgeSecondNode : edgeFirstNode;
            intersections[edgeIndex].edgeSecondNode = crossProductValue < 0 ? edgeFirstNode : edgeSecondNode;
            intersections[edgeIndex].edgeDistance = adimensionalEdgeDistance;
            intersections[edgeIndex].edgeIndex = edgeIndex;
        }
 
        static void updateFaceIntersections(const UInt faceIndex,
                                            const UInt edgeIndex,
                                            std::vector<FaceMeshPolyLineIntersection>& intersections)
        {
            intersections[faceIndex].faceIndex = faceIndex;
            intersections[faceIndex].edgeIndices.emplace_back(edgeIndex);
        }
 
        Mesh2D& m_mesh;                                                      
        std::vector<EdgeMeshPolyLineIntersection> m_edgesIntersectionsCache; 
        std::vector<FaceMeshPolyLineIntersection> m_facesIntersectionsCache; 
        std::vector<EdgeMeshPolyLineIntersection> m_edgesIntersections;      
        std::vector<FaceMeshPolyLineIntersection> m_faceIntersections;       
        BoundingBox m_meshBoundingBox;                                       
        std::vector<BoundingBox> m_meshEdgesBoundingBoxes;                   
        static constexpr UInt maxSearchSegments = 1000;                      
    };
 
} // namespace meshkernel