Skip to content

Commit

Permalink
Improve mesh deletion by introducing InsideNotIntersected and InsideA…
Browse files Browse the repository at this point in the history
…ndIntersected options (#238 | GRIDEDIT-622)
  • Loading branch information
lucacarniato authored Oct 10, 2023
1 parent 225f4b4 commit 4bd6fec
Show file tree
Hide file tree
Showing 21 changed files with 884 additions and 561 deletions.
2 changes: 2 additions & 0 deletions libs/MeshKernel/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,7 @@ set(
${SRC_DIR}/Mesh.cpp
${SRC_DIR}/Mesh1D.cpp
${SRC_DIR}/Mesh2D.cpp
${SRC_DIR}/Mesh2DIntersections.cpp
${SRC_DIR}/MeshRefinement.cpp
${SRC_DIR}/Network1D.cpp
${SRC_DIR}/Operations.cpp
Expand Down Expand Up @@ -110,6 +111,7 @@ set(
${DOMAIN_INC_DIR}/Mesh.hpp
${DOMAIN_INC_DIR}/Mesh1D.hpp
${DOMAIN_INC_DIR}/Mesh2D.hpp
${DOMAIN_INC_DIR}/Mesh2DIntersections.hpp
${DOMAIN_INC_DIR}/MeshInterpolation.hpp
${DOMAIN_INC_DIR}/MeshRefinement.hpp
${DOMAIN_INC_DIR}/Network1D.hpp
Expand Down
1 change: 0 additions & 1 deletion libs/MeshKernel/include/MeshKernel/LandBoundary.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,6 @@
#include <vector>

#include "MeshKernel/BoundingBox.hpp"
#include "MeshKernel/Entities.hpp"
#include "MeshKernel/Polygons.hpp"

namespace meshkernel
Expand Down
21 changes: 0 additions & 21 deletions libs/MeshKernel/include/MeshKernel/Mesh.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -109,27 +109,6 @@ namespace meshkernel
{Location::Edges, "Edges"},
{Location::Unknown, "Unknown"}};

/// edge-segment intersection
struct EdgeMeshPolylineIntersection
{
int polylineSegmentIndex{constants::missing::intValue}; ///< The intersected segment index (a polyline can formed by several segments)
double polylineDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as distance from the polyline start
double adimensionalPolylineSegmentDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as an adimensional distance from the segment start
UInt edgeIndex{constants::missing::uintValue}; ///< The edge index
UInt edgeFirstNode{constants::missing::uintValue}; ///< The first node of the edge is on the left (the virtual node)
UInt edgeSecondNode{constants::missing::uintValue}; ///< The second node of the edge is on the right (the inner node)
double edgeDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as an adimensional distance from the edge start
};

/// face-segment intersection
struct FaceMeshPolylineIntersection
{
double polylineDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as an adimensional distance from the polyline start
UInt faceIndex{constants::missing::uintValue}; ///< The face index
std::vector<UInt> edgeIndexses; ///< The indexes of crossed edges
std::vector<UInt> edgeNodes; ///< The indexes of the nodes defining the crossed edges
};

/// @brief Default constructor
Mesh() = default;

Expand Down
12 changes: 2 additions & 10 deletions libs/MeshKernel/include/MeshKernel/Mesh2D.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -55,9 +55,8 @@ namespace meshkernel
/// Enumerator describing the different options to delete a mesh
enum DeleteMeshOptions
{
AllNodesInside = 0,
FacesWithIncludedCircumcenters = 1,
FacesCompletelyIncluded = 2
InsideNotIntersected = 0,
InsideAndIntersected = 1
};

/// Enumerator describing the different node types
Expand Down Expand Up @@ -271,13 +270,6 @@ namespace meshkernel
/// @return A tuple with the intersectedFace face index and intersected edge index
[[nodiscard]] std::tuple<UInt, UInt> IsSegmentCrossingABoundaryEdge(const Point& firstPoint, const Point& secondPoint) const;

/// @brief Gets the edges and faces intersected by a polyline, with additional information on the intersections
/// @param[in] polyLine An input polyline, defined as a series of points
/// @return A tuple containing a vector of EdgeMeshPolylineIntersections and FaceMeshPolylineIntersections
[[nodiscard]] std::tuple<std::vector<EdgeMeshPolylineIntersection>,
std::vector<FaceMeshPolylineIntersection>>
GetPolylineIntersections(const std::vector<Point>& polyLine);

/// @brief Masks the edges of all faces entirely included in all polygons
/// @param[in] polygons The selection polygon
/// @param[in] invertSelection Invert selection
Expand Down
170 changes: 170 additions & 0 deletions libs/MeshKernel/include/MeshKernel/Mesh2DIntersections.hpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,170 @@
//---- 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: [email protected]
// 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"

namespace meshkernel
{
/// An intersection with a mesh edge
struct EdgeMeshPolyLineIntersection
{
int polylineSegmentIndex{constants::missing::intValue}; ///< The intersected segment index (a polyline can formed by several segments)
double polylineDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as distance from the polyline start
double adimensionalPolylineSegmentDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as an adimensional distance from the segment start
UInt edgeIndex{constants::missing::uintValue}; ///< The edge index
UInt edgeFirstNode{constants::missing::uintValue}; ///< The first node of the edge is on the left (the virtual node)
UInt edgeSecondNode{constants::missing::uintValue}; ///< The second node of the edge is on the right (the inner node)
double edgeDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as an adimensional distance from the edge start
};

/// An intersection with a mesh face
struct FaceMeshPolyLineIntersection
{
double polylineDistance{constants::missing::doubleValue}; ///< The location of the intersection expressed as an adimensional distance from the polyline start
UInt faceIndex{constants::missing::uintValue}; ///< The face index
std::vector<UInt> edgeIndices; ///< The indexes of crossed edges
std::vector<UInt> edgeNodes; ///< The indexes of the nodes defining the crossed edges
};

/// @brief Compute the intersections of polygon inner and outer perimeters
///
/// @note Uses a breadth first algorithm to reduce runtime complexity
class Mesh2DIntersections final
{
public:
/// @brief Constructor
explicit Mesh2DIntersections(Mesh2D& mesh);

/// @brief Compute intersection with a polygon, possibly containing multiple polylines
/// @param[in] polygon An input polygon
void Compute(const Polygons& polygon);

/// @brief Compute intersection with a single polyline
/// @param[in] polyLine An input polyline
void Compute(const std::vector<Point>& polyLine);

/// @brief Gets the edge intersections
/// @returns The edges intersections
[[nodiscard]] const auto& EdgeIntersections() const { return m_edgesIntersections; }

/// @brief Gets the face intersections
/// @returns The faces intersections
[[nodiscard]] const auto& FaceIntersections() const { return m_faceIntersections; }

/// @brief Sort intersections by polyline distance and erase entries with no intersections
/// @tparam T An intersection type, \ref EdgeMeshPolyLineIntersection or \ref FaceMeshPolyLineIntersection
/// @param intersections a vector containing the intersections
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:
/// @brief Enumeration defining the directions for searching the next segment index
enum class Direction
{
Forward,
Backward
};

/// @brief Gets one edge intersection
/// @returns The intersection seed
std::tuple<UInt, UInt> GetIntersectionSeed(const Mesh2D& mesh,
const std::vector<Point>& polyLine,
const std::vector<bool>& vistedEdges) const;

/// @brief Gets the next edge intersection
/// @returns The intersection seed
std::tuple<bool, UInt, UInt, double, double, double> GetNextEdgeIntersection(const std::vector<Point>& polyLine,
UInt edgeIndex,
UInt firstIndex,
UInt secondIndex,
Direction direction) const;

/// @brief Gets the next edge intersection
/// @returns The intersection seed
void IntersectFaceEdges(const std::vector<Point>& polyLine,
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);

/// @brief Update edge intersections
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;
}

/// @brief Update face intersections
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; ///< The mesh where the edges should be found
std::vector<EdgeMeshPolyLineIntersection> m_edgesIntersectionsCache; ///< A cache for saving the edge intersections of one inner or outer
std::vector<FaceMeshPolyLineIntersection> m_facesIntersectionsCache; ///< A cache for saving the local face intersections of one inner or outer
std::vector<EdgeMeshPolyLineIntersection> m_edgesIntersections; ///< A vector collecting all edge intersection results
std::vector<FaceMeshPolyLineIntersection> m_faceIntersections; ///< A vector collecting all face intersection results
static constexpr UInt maxSearchSegments = 10; ///< mex number of steps in polyline intersection algorithm
};

} // namespace meshkernel
1 change: 0 additions & 1 deletion libs/MeshKernel/include/MeshKernel/MeshRefinement.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,6 @@
#pragma once

#include <MeshKernel/AveragingInterpolation.hpp>
#include <MeshKernel/Entities.hpp>
#include <MeshKernel/Parameters.hpp>
#include <MeshKernel/Polygons.hpp>
#include <MeshKernel/Utilities/RTree.hpp>
Expand Down
27 changes: 12 additions & 15 deletions libs/MeshKernel/include/MeshKernel/Operations.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -451,21 +451,18 @@ namespace meshkernel
/// @param[in] secondSegmentSecondPoint The second point of the second segment
/// @param[in] adimensionalCrossProduct Whether to compute the dimensionless cross product
/// @param[in] projection The coordinate system projection
/// @param[out] intersectionPoint The intersection point
/// @param[out] crossProduct The cross product of the intersection
/// @param[out] ratioFirstSegment The distance of the intersection from the first node of the first segment, expressed as a ratio of the segment length
/// @param[out] ratioSecondSegment The distance of the intersection from the first node of the second segment, expressed as a ratio of the segment length
/// @return If the two segments are crossing
[[nodiscard]] bool AreSegmentsCrossing(const Point& firstSegmentFirstPoint,
const Point& firstSegmentSecondPoint,
const Point& secondSegmentFirstPoint,
const Point& secondSegmentSecondPoint,
bool adimensionalCrossProduct,
const Projection& projection,
Point& intersectionPoint,
double& crossProduct,
double& ratioFirstSegment,
double& ratioSecondSegment);
/// @return A tuple with:
/// If the two segments are crossing
/// The intersection point
/// The cross product of the intersection
/// The distance of the intersection from the first node of the first segment, expressed as a ratio of the segment length
/// The distance of the intersection from the first node of the second segment, expressed as a ratio of the segment length
[[nodiscard]] std::tuple<bool, Point, double, double, double> AreSegmentsCrossing(const Point& firstSegmentFirstPoint,
const Point& firstSegmentSecondPoint,
const Point& secondSegmentFirstPoint,
const Point& secondSegmentSecondPoint,
bool adimensionalCrossProduct,
const Projection& projection);

/// @brief Computes the coordinate of a point on a spline, given the dimensionless distance from the first corner point (splint)
/// @param[in] coordinates The spline node coordinates
Expand Down
2 changes: 0 additions & 2 deletions libs/MeshKernel/include/MeshKernel/Polygons.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -27,11 +27,9 @@

#pragma once

#include <unordered_map>
#include <vector>

#include <MeshKernel/BoundingBox.hpp>
#include <MeshKernel/Entities.hpp>
#include <MeshKernel/Exceptions.hpp>
#include <MeshKernel/Polygon.hpp>
#include <MeshKernel/PolygonalEnclosure.hpp>
Expand Down
26 changes: 10 additions & 16 deletions libs/MeshKernel/src/ConnectMeshes.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -562,22 +562,16 @@ void meshkernel::ConnectMeshes::FreeHangingNodes(Mesh2D& mesh,
UInt startNode = mesh.m_edges[oppositeEdgeId].first;
UInt endNode = mesh.m_edges[oppositeEdgeId].second;

double crossProduct;
double normalisedPolylineSegmentDistance;
double normalisedEdgeDistance;
Point intersectionPoint;

const bool segmentsCross = AreSegmentsCrossing(mesh.m_nodes[boundaryEdge.first],
mesh.m_nodes[startNode],
mesh.m_nodes[boundaryEdge.second],
mesh.m_nodes[endNode],
false,
mesh.m_projection,
intersectionPoint,
crossProduct,
normalisedPolylineSegmentDistance,
normalisedEdgeDistance);

const auto [segmentsCross,
intersectionPoint,
crossProduct,
normalisedPolylineSegmentDistance,
normalisedEdgeDistance] = AreSegmentsCrossing(mesh.m_nodes[boundaryEdge.first],
mesh.m_nodes[startNode],
mesh.m_nodes[boundaryEdge.second],
mesh.m_nodes[endNode],
false,
mesh.m_projection);
if (segmentsCross)
{
std::swap(startNode, endNode);
Expand Down
Loading

0 comments on commit 4bd6fec

Please sign in to comment.