diff --git a/libs/MeshKernel/CMakeLists.txt b/libs/MeshKernel/CMakeLists.txt index c39bdc675..a06794096 100644 --- a/libs/MeshKernel/CMakeLists.txt +++ b/libs/MeshKernel/CMakeLists.txt @@ -116,6 +116,7 @@ set( ${DOMAIN_INC_DIR}/Mesh2DIntersections.hpp ${DOMAIN_INC_DIR}/MeshInterpolation.hpp ${DOMAIN_INC_DIR}/MeshRefinement.hpp + ${DOMAIN_INC_DIR}/MeshTransformation.hpp ${DOMAIN_INC_DIR}/Network1D.hpp ${DOMAIN_INC_DIR}/Operations.hpp ${DOMAIN_INC_DIR}/OrthogonalizationAndSmoothing.hpp diff --git a/libs/MeshKernel/include/MeshKernel/MeshTransformation.hpp b/libs/MeshKernel/include/MeshKernel/MeshTransformation.hpp new file mode 100644 index 000000000..9d97caae0 --- /dev/null +++ b/libs/MeshKernel/include/MeshKernel/MeshTransformation.hpp @@ -0,0 +1,287 @@ +//---- 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 . +// +// 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 +#include + +#include "MeshKernel/Definitions.hpp" +#include "MeshKernel/Exceptions.hpp" +#include "MeshKernel/Mesh.hpp" +#include "MeshKernel/Point.hpp" +#include "MeshKernel/Vector.hpp" + +namespace meshkernel +{ + + /// @brief Ensure any instantiation of the MeshTransformation Compute function is with the correct operation + template + concept HasTransformationFunction = requires(Function op, Point p) {{ op(p)} -> std::same_as; }; + + /// @brief Ensure any instantiation of the MeshTransformation Compute function is able to determine the projection. + template + concept HasTransformationProjection = requires(Function op) {{ op.TransformationProjection()} -> std::same_as; }; + + /// @brief To ensure the MeshTransformation Compute template parameter has a valid interface. + template + concept TransformationFunction = HasTransformationFunction && HasTransformationProjection; + + /// @brief Apply a translation transformation to a point or a vector. + class Translation + { + public: + /// @brief Default constructor, default is no translation + Translation() = default; + + /// @brief Construct with user defined translation + explicit Translation(const Vector& trans) : m_translation(trans) {} + + /// @brief Get the projection required for the translation + [[nodiscard]] Projection TransformationProjection() const + { + return Projection::cartesian; + } + + /// @brief Reset translation to identity translation (i.e. no translation) + void identity() + { + m_translation = {0.0, 0.0}; + } + + /// @brief Reset the translation to a new translation quantity + void reset(const Vector& trans) + { + m_translation = trans; + } + + /// @brief Get the current defined translation vector. + const Vector& vector() const + { + return m_translation; + } + + /// @brief Compose two translation objects. + Translation compose(const Translation& trans) const + { + return Translation(m_translation + trans.m_translation); + } + + /// @brief Apply the translation to a point in Cartesian coordinate system + Point operator()(const Point& pnt) const + { + return pnt + m_translation; + } + + /// @brief Apply the translation to a vector in Cartesian coordinate system + Vector operator()(const Vector& vec) const + { + return vec + m_translation; + } + + private: + /// @brief The translation values + Vector m_translation{0.0, 0.0}; + }; + + /// @brief Apply a rotation transformation to a point or a vector. + class Rotation + { + public: + /// @brief Default constructor, default is theta = 0. + Rotation() = default; + + /// @brief Construct with user defined rotation angle, in degrees. + explicit Rotation(const double angle) + { + reset(angle); + } + + /// @brief Get the projection required for the rotation + [[nodiscard]] Projection TransformationProjection() const + { + return Projection::cartesian; + } + + /// @brief Reset rotation to identity translation (i.e. no rotation, theta = 0) + void identity() + { + m_theta = 0.0; + m_cosTheta = 1.0; + m_sinTheta = 0.0; + } + + /// @brief Reset the rotation to a new rotation angle + /// + /// The angle in degrees. + void reset(const double angle) + { + m_theta = angle; + m_cosTheta = std::cos(m_theta * constants::conversion::degToRad); + m_sinTheta = std::sin(m_theta * constants::conversion::degToRad); + } + + /// @brief Get the current defined rotation angle in degrees + double angle() const + { + return m_theta; + } + + /// @brief Compose two rotation objects. + Rotation compose(const Rotation& rot) const + { + return Rotation(m_theta + rot.m_theta); + } + + /// @brief Apply the rotation to a point in Cartesian coordinate system + Point operator()(const Point& pnt) const + { + Point result({m_cosTheta * pnt.x - m_sinTheta * pnt.y, + m_sinTheta * pnt.x + m_cosTheta * pnt.y}); + return result; + } + + /// @brief Apply the rotation to a vector in Cartesian coordinate system + Vector operator()(const Vector& vec) const + { + Vector result({m_cosTheta * vec.x() - m_sinTheta * vec.y(), + m_sinTheta * vec.x() + m_cosTheta * vec.y()}); + return result; + } + + private: + /// @brief The rotation angle, theta + double m_theta = 0.0; + + /// @brief cosine of the rotation angle + double m_cosTheta = 1.0; + + /// @brief sine of the rotation angle + double m_sinTheta = 0.0; + }; + + /// @brief A composition of translation and rotation transformations + class RigidBodyTransformation + { + public: + /// @brief Default constructor, default is no transformation + RigidBodyTransformation() = default; + + /// @brief Get the projection required for the transformation + [[nodiscard]] Projection TransformationProjection() const + { + return Projection::cartesian; + } + + /// @brief Reset transformation to identity transformation (i.e. no transformation) + void identity() + { + m_rotation.identity(); + m_translation.identity(); + } + + /// @brief Compose rotation and transformation object. + /// + /// Will be applied: rot (transformation). + void compose(const Rotation& rot) + { + m_rotation = rot.compose(m_rotation); + m_translation.reset(rot(m_translation.vector())); + } + + /// @brief Compose translation and transformation object. + /// + /// Will be applied: translation (transformation). + void compose(const Translation& trans) + { + m_translation = trans.compose(m_translation); + } + + /// @brief Get the current rotation. + const Rotation& rotation() const + { + return m_rotation; + } + + /// @brief Get the current translation. + const Translation& translation() const + { + return m_translation; + } + + /// @brief Apply the transformation to a point in Cartesian coordinate system + Point operator()(const Point& pnt) const + { + Point result = m_rotation(pnt); + result = m_translation(result); + return result; + } + + /// @brief Apply the transformation to a vector in Cartesian coordinate system + Vector operator()(const Vector& vec) const + { + Vector result = m_rotation(vec); + result = m_translation(result); + return result; + } + + private: + /// @brief The rotation part of the transformation + Rotation m_rotation; + + /// @brief The translation part of the transformation + Translation m_translation; + }; + + /// @brief Apply a transformation to a mesh + class MeshTransformation + { + public: + /// @brief Apply a transformation to a mesh with a Cartesian projection + template + static void Compute(Mesh& mesh, Transformation transformation) + { + if (mesh.m_projection != transformation.TransformationProjection()) + { + throw MeshKernelError("Incorrect mesh coordinate system, expecting '{}', found '{}'", + ToString(transformation.TransformationProjection()), ToString(mesh.m_projection)); + } + +#pragma omp parallel for + for (int i = 0; i < static_cast(mesh.GetNumNodes()); ++i) + { + if (mesh.m_nodes[i].IsValid()) + { + mesh.m_nodes[i] = transformation(mesh.m_nodes[i]); + } + } + + mesh.Administrate(); + } + }; + +} // namespace meshkernel diff --git a/libs/MeshKernel/tests/CMakeLists.txt b/libs/MeshKernel/tests/CMakeLists.txt index 1af27935c..a91e878f6 100644 --- a/libs/MeshKernel/tests/CMakeLists.txt +++ b/libs/MeshKernel/tests/CMakeLists.txt @@ -36,6 +36,7 @@ set( ${SRC_DIR}/Mesh2DTest.cpp ${SRC_DIR}/MeshRefinementTests.cpp ${SRC_DIR}/MeshTests.cpp + ${SRC_DIR}/MeshTransformationTest.cpp ${SRC_DIR}/OrthogonalizationTests.cpp ${SRC_DIR}/ParametersTests.cpp ${SRC_DIR}/PointTests.cpp diff --git a/libs/MeshKernel/tests/src/MeshTransformationTest.cpp b/libs/MeshKernel/tests/src/MeshTransformationTest.cpp new file mode 100644 index 000000000..a638167ff --- /dev/null +++ b/libs/MeshKernel/tests/src/MeshTransformationTest.cpp @@ -0,0 +1,245 @@ +#include +#include +#include + +#include "MeshKernel/Constants.hpp" +#include "MeshKernel/Entities.hpp" +#include "MeshKernel/Mesh2D.hpp" +#include "MeshKernel/MeshTransformation.hpp" +#include "MeshKernel/Polygons.hpp" + +#include "TestUtils/Definitions.hpp" +#include "TestUtils/MakeMeshes.hpp" + +namespace mk = meshkernel; + +TEST(MeshTransformationTest, BasicTranslationTest) +{ + // Test basic functionality of the translation class + double translationX = 1.0; + double translationY = 2.0; + + mk::Vector vec(translationX, translationY); + mk::Translation translation(vec); + + EXPECT_EQ(translationX, translation.vector().x()); + EXPECT_EQ(translationY, translation.vector().y()); + + translation.identity(); + EXPECT_EQ(0.0, translation.vector().x()); + EXPECT_EQ(0.0, translation.vector().y()); + + std::swap(vec.x(), vec.y()); + vec.x() = -vec.x(); + vec.y() = -vec.y(); + translation.reset(vec); + + EXPECT_EQ(-translationY, translation.vector().x()); + EXPECT_EQ(-translationX, translation.vector().y()); +} + +TEST(MeshTransformationTest, BasicRotationTest) +{ + // Test basic functionality of the rotation class + double theta = 45.0; + + mk::Rotation rotation(theta); + + EXPECT_EQ(theta, rotation.angle()); + + theta *= 1.5; + rotation.reset(theta); + EXPECT_EQ(theta, rotation.angle()); + + rotation.identity(); + EXPECT_EQ(0.0, rotation.angle()); +} + +TEST(MeshTransformationTest, PointTranslationTest) +{ + // Test a simple translation of a point in Cartesian coordinates. + double originX = 12.0; + double originY = -19.0; + double translationX = 1.0; + double translationY = 2.0; + + mk::Vector vec(translationX, translationY); + mk::Translation translation(vec); + mk::Point pnt(originX, originY); + pnt = translation(pnt); + + EXPECT_EQ(originX + translationX, pnt.x); + EXPECT_EQ(originY + translationY, pnt.y); +} + +TEST(MeshTransformationTest, BasicRigidBodyTransformationTest) +{ + // Test the rigid body transformation + // The test will rotate a series of points about a point that is not the origin + // + // Step 1: create composition of all simple transformations + // 1. translate the rotation point (-1,1) to origin (0,0) + // 2. rotate by 90 degrees + // 3. translate back to original rotation point (-1,1) + // + // Step 2 apply composite transformation to a series of points + constexpr double tolerance = 1.0e-10; + + double theta = -90.0; + mk::Point rotationPoint(-1.0, 1.0); + + mk::RigidBodyTransformation transformation; + + // Expect initial state to be the identity + EXPECT_EQ(0.0, transformation.rotation().angle()); + EXPECT_EQ(0.0, transformation.translation().vector().x()); + EXPECT_EQ(0.0, transformation.translation().vector().y()); + + // First translate + mk::Vector vec(-rotationPoint.x, -rotationPoint.y); + mk::Translation translation(vec); + transformation.compose(translation); + + EXPECT_EQ(0.0, transformation.rotation().angle()); + EXPECT_EQ(-rotationPoint.x, transformation.translation().vector().x()); + EXPECT_EQ(-rotationPoint.y, transformation.translation().vector().y()); + + // Second rotate + mk::Rotation rotation(theta); + transformation.compose(rotation); + + mk::Vector rotatedVector = rotation(vec); + + EXPECT_EQ(theta, transformation.rotation().angle()); + EXPECT_EQ(rotatedVector.x(), transformation.translation().vector().x()); + EXPECT_EQ(rotatedVector.y(), transformation.translation().vector().y()); + + // Third translate back + vec = {rotationPoint.x, rotationPoint.y}; + translation.reset(vec); + transformation.compose(translation); + + EXPECT_EQ(theta, transformation.rotation().angle()); + + double expectedTranslationX = -2.0; + double expectedTranslationY = 0.0; + + EXPECT_NEAR(expectedTranslationX, transformation.translation().vector().x(), tolerance); + EXPECT_NEAR(expectedTranslationY, transformation.translation().vector().y(), tolerance); + + //------------------------------------------------------------ + // Apply the composite transformation to a series of points + + mk::Point pnt1{-3.0, 4.0}; + mk::Point pnt2{-2.0, 2.0}; + mk::Point pnt3{-3.0, 1.0}; + + mk::Point expected1{2.0, 3.0}; + mk::Point expected2{0.0, 2.0}; + mk::Point expected3{-1.0, 3.0}; + + pnt1 = transformation(pnt1); + EXPECT_NEAR(expected1.x, pnt1.x, tolerance); + EXPECT_NEAR(expected1.y, pnt1.y, tolerance); + + pnt2 = transformation(pnt2); + EXPECT_NEAR(expected2.x, pnt2.x, tolerance); + EXPECT_NEAR(expected2.y, pnt2.y, tolerance); + + pnt3 = transformation(pnt3); + EXPECT_NEAR(expected3.x, pnt3.x, tolerance); + EXPECT_NEAR(expected3.y, pnt3.y, tolerance); +} + +TEST(MeshTransformationTest, PointRotationTest) +{ + // Test a simple rotation of a point in Cartesian coordinates. + double originX = 12.0; + double originY = -19.0; + + double theta = 45.0; + double cosTheta = std::cos(theta * M_PI / 180.0); + double sinTheta = std::sin(theta * M_PI / 180.0); + + mk::Rotation rotation(theta); + mk::Point pnt(originX, originY); + pnt = rotation(pnt); + + EXPECT_EQ(originX * cosTheta - originY * sinTheta, pnt.x); + EXPECT_EQ(originX * sinTheta + originY * cosTheta, pnt.y); +} + +TEST(MeshTransformationTest, MeshTranslationTest) +{ + // Test the translation on the entire mesh. + + mk::UInt nx = 11; + mk::UInt ny = 11; + + double delta = 10.0; + + // Make two copies of the same mesh so that they can be compared in the test. + std::shared_ptr originalMesh = MakeRectangularMeshForTesting(nx, ny, delta, mk::Projection::cartesian); + std::shared_ptr mesh = MakeRectangularMeshForTesting(nx, ny, delta, mk::Projection::cartesian); + + mk::Vector vec(1.0, 2.0); + + mk::Translation translation(vec); + + mk::MeshTransformation::Compute(*mesh, translation); + + for (mk::UInt i = 0; i < mesh->GetNumNodes(); ++i) + { + EXPECT_EQ(originalMesh->m_nodes[i].x + vec.x(), mesh->m_nodes[i].x); + EXPECT_EQ(originalMesh->m_nodes[i].y + vec.y(), mesh->m_nodes[i].y); + } +} + +TEST(MeshTransformationTest, MeshRotationTest) +{ + // Test the rotation on the entire mesh. + + mk::UInt nx = 11; + mk::UInt ny = 11; + + double delta = 10.0; + + // Make two copies of the same mesh so that they can be compared in the test. + std::shared_ptr originalMesh = MakeRectangularMeshForTesting(nx, ny, delta, mk::Projection::cartesian); + std::shared_ptr mesh = MakeRectangularMeshForTesting(nx, ny, delta, mk::Projection::cartesian); + + double theta = 45.0; + + mk::Rotation rotation(theta); + double cosTheta = std::cos(theta * M_PI / 180.0); + double sinTheta = std::sin(theta * M_PI / 180.0); + + mk::MeshTransformation::Compute(*mesh, rotation); + + for (mk::UInt i = 0; i < mesh->GetNumNodes(); ++i) + { + mk::Point expected{originalMesh->m_nodes[i].x * cosTheta - originalMesh->m_nodes[i].y * sinTheta, + originalMesh->m_nodes[i].x * sinTheta + originalMesh->m_nodes[i].y * cosTheta}; + EXPECT_EQ(expected.x, mesh->m_nodes[i].x); + EXPECT_EQ(expected.y, mesh->m_nodes[i].y); + } +} + +TEST(MeshTransformationTest, IncorrectProjectionTest) +{ + // Test correct failure with non Cartesian projection. + + // Generate mesh in spherical coordinate system + std::shared_ptr mesh = MakeRectangularMeshForTesting(11, 11, 10.0, mk::Projection::spherical); + + mk::Rotation rotation(45.0); + + // Should throw an exception with spherical coordinate system + EXPECT_THROW(mk::MeshTransformation::Compute(*mesh, rotation), mk::MeshKernelError); + + // Change projection to Projection::sphericalAccurate + mesh->m_projection = mk::Projection::sphericalAccurate; + + // Should throw an exception with spherical-accurate coordinate system + EXPECT_THROW(mk::MeshTransformation::Compute(*mesh, rotation), mk::MeshKernelError); +} diff --git a/libs/MeshKernelApi/include/MeshKernelApi/MeshKernel.hpp b/libs/MeshKernelApi/include/MeshKernelApi/MeshKernel.hpp index c1dfcd4dc..6eb48b58f 100644 --- a/libs/MeshKernelApi/include/MeshKernelApi/MeshKernel.hpp +++ b/libs/MeshKernelApi/include/MeshKernelApi/MeshKernel.hpp @@ -710,6 +710,15 @@ namespace meshkernelapi /// @returns Error code MKERNEL_API int mkernel_mesh2d_compute_inner_ortogonalization_iteration(int meshKernelId); + /// @brief Rotate a mesh2d about a point. + /// + /// @param[in] meshKernelId The id of the mesh state + /// @param[in] centreX X-coordinate of the centre of rotation + /// @param[in] centreY Y-coordinate of the centre of rotation + /// @param[in] theta Angle of rotation + /// @returns Error code + MKERNEL_API int mkernel_mesh2d_rotate(int meshKernelId, double centreX, double centreY, double theta); + /// @brief Snaps the spline (or splines) to the land boundary /// /// @param[in] meshKernelId The id of the mesh state @@ -724,6 +733,14 @@ namespace meshkernelapi int startSplineIndex, int endSplineIndex); + /// @brief Translate a mesh2d + /// + /// @param[in] meshKernelId The id of the mesh state + /// @param[in] translationX X-component of the translation vector + /// @param[in] translationY Y-component of the translation vector + /// @returns Error code + MKERNEL_API int mkernel_mesh2d_translate(int meshKernelId, double translationX, double translationY); + /// The function modifies the mesh for achieving orthogonality between the edges and the segments connecting the face circumcenters. /// The amount of orthogonality is traded against the mesh smoothing (in this case the equality of face areas). /// The parameter to regulate the amount of orthogonalization is contained in \ref meshkernel::OrthogonalizationParameters::orthogonalization_to_smoothing_factor diff --git a/libs/MeshKernelApi/src/MeshKernel.cpp b/libs/MeshKernelApi/src/MeshKernel.cpp index 2fcc78685..63ee6d1a8 100644 --- a/libs/MeshKernelApi/src/MeshKernel.cpp +++ b/libs/MeshKernelApi/src/MeshKernel.cpp @@ -55,6 +55,7 @@ #include #include #include +#include #include #include #include @@ -1787,6 +1788,57 @@ namespace meshkernelapi return lastExitCode; } + MKERNEL_API int mkernel_mesh2d_rotate(int meshKernelId, double centreX, double centreY, double theta) + { + lastExitCode = meshkernel::ExitCode::Success; + + try + { + if (!meshKernelState.contains(meshKernelId)) + { + throw meshkernel::MeshKernelError("The selected mesh kernel id does not exist."); + } + + meshkernel::RigidBodyTransformation transformation; + meshkernel::Translation translation; + + // Translate centre of rotation to origin + transformation.compose(meshkernel::Translation(meshkernel::Vector(-centreX, -centreY))); + // Add rotation + transformation.compose(meshkernel::Rotation(theta)); + // Translate origin back to centre of rotation + transformation.compose(meshkernel::Translation(meshkernel::Vector(centreX, centreY))); + + meshkernel::MeshTransformation::Compute(*meshKernelState[meshKernelId].m_mesh2d, transformation); + } + catch (...) + { + lastExitCode = HandleException(); + } + return lastExitCode; + } + + MKERNEL_API int mkernel_mesh2d_translate(int meshKernelId, double translationX, double translationY) + { + lastExitCode = meshkernel::ExitCode::Success; + + try + { + if (!meshKernelState.contains(meshKernelId)) + { + throw meshkernel::MeshKernelError("The selected mesh kernel id does not exist."); + } + + meshkernel::Translation translation(meshkernel::Vector(translationX, translationY)); + meshkernel::MeshTransformation::Compute(*meshKernelState[meshKernelId].m_mesh2d, translation); + } + catch (...) + { + lastExitCode = HandleException(); + } + return lastExitCode; + } + MKERNEL_API int mkernel_splines_snap_to_landboundary(int meshKernelId, const GeometryList& land, GeometryList& splines, diff --git a/libs/MeshKernelApi/tests/src/ApiTest.cpp b/libs/MeshKernelApi/tests/src/ApiTest.cpp index 06bf5c62c..a2bc37487 100644 --- a/libs/MeshKernelApi/tests/src/ApiTest.cpp +++ b/libs/MeshKernelApi/tests/src/ApiTest.cpp @@ -1,6 +1,7 @@ #include #include +#include #include #include @@ -2487,3 +2488,183 @@ TEST_F(CartesianApiTestFixture, GenerateTriangularGridThroughApi_OnClockWisePoly ASSERT_EQ(5, mesh2d.num_nodes); ASSERT_EQ(8, mesh2d.num_edges); } + +TEST_F(CartesianApiTestFixture, TranslateMesh) +{ + // Prepare + + MakeMesh(); + auto const meshKernelId = GetMeshKernelId(); + + meshkernelapi::Mesh2D mesh2d{}; + [[maybe_unused]] int errorCode = mkernel_mesh2d_get_dimensions(meshKernelId, mesh2d); + + // Get copy of original mesh + // Will be used later to compare with translated mesh + std::vector edge_faces(mesh2d.num_edges * 2); + std::vector edge_nodes(mesh2d.num_edges * 2); + std::vector face_nodes(mesh2d.num_face_nodes); + std::vector face_edges(mesh2d.num_face_nodes); + std::vector nodes_per_face(mesh2d.num_faces); + std::vector node_x(mesh2d.num_nodes); + std::vector node_y(mesh2d.num_nodes); + std::vector edge_x(mesh2d.num_edges); + std::vector edge_y(mesh2d.num_edges); + std::vector face_x(mesh2d.num_faces); + std::vector face_y(mesh2d.num_faces); + + mesh2d.edge_faces = edge_faces.data(); + mesh2d.edge_nodes = edge_nodes.data(); + mesh2d.face_nodes = face_nodes.data(); + mesh2d.face_edges = face_edges.data(); + mesh2d.nodes_per_face = nodes_per_face.data(); + mesh2d.node_x = node_x.data(); + mesh2d.node_y = node_y.data(); + mesh2d.edge_x = edge_x.data(); + mesh2d.edge_y = edge_y.data(); + mesh2d.face_x = face_x.data(); + mesh2d.face_y = face_y.data(); + errorCode = mkernel_mesh2d_get_data(meshKernelId, mesh2d); + + double translationX = 10.0; + double translationY = 5.0; + + meshkernelapi::mkernel_mesh2d_translate(meshKernelId, translationX, translationY); + + meshkernelapi::Mesh2D mesh2dTranslated{}; + errorCode = mkernel_mesh2d_get_dimensions(meshKernelId, mesh2dTranslated); + + // Data for translated mesh + std::vector edge_faces_t(mesh2dTranslated.num_edges * 2); + std::vector edge_nodes_t(mesh2dTranslated.num_edges * 2); + std::vector face_nodes_t(mesh2dTranslated.num_face_nodes); + std::vector face_edges_t(mesh2dTranslated.num_face_nodes); + std::vector nodes_per_face_t(mesh2dTranslated.num_faces); + std::vector node_x_t(mesh2dTranslated.num_nodes); + std::vector node_y_t(mesh2dTranslated.num_nodes); + std::vector edge_x_t(mesh2dTranslated.num_edges); + std::vector edge_y_t(mesh2dTranslated.num_edges); + std::vector face_x_t(mesh2dTranslated.num_faces); + std::vector face_y_t(mesh2dTranslated.num_faces); + + mesh2dTranslated.edge_faces = edge_faces_t.data(); + mesh2dTranslated.edge_nodes = edge_nodes_t.data(); + mesh2dTranslated.face_nodes = face_nodes_t.data(); + mesh2dTranslated.face_edges = face_edges_t.data(); + mesh2dTranslated.nodes_per_face = nodes_per_face_t.data(); + mesh2dTranslated.node_x = node_x_t.data(); + mesh2dTranslated.node_y = node_y_t.data(); + mesh2dTranslated.edge_x = edge_x_t.data(); + mesh2dTranslated.edge_y = edge_y_t.data(); + mesh2dTranslated.face_x = face_x_t.data(); + mesh2dTranslated.face_y = face_y_t.data(); + + errorCode = mkernel_mesh2d_get_data(meshKernelId, mesh2dTranslated); + + const double tolerance = 1e-12; + + for (int i = 0; i < mesh2d.num_nodes; ++i) + { + EXPECT_NEAR(mesh2d.node_x[i] + translationX, mesh2dTranslated.node_x[i], tolerance); + EXPECT_NEAR(mesh2d.node_y[i] + translationY, mesh2dTranslated.node_y[i], tolerance); + } +} + +TEST_F(CartesianApiTestFixture, RotateMesh) +{ + // Prepare + + MakeMesh(); + auto const meshKernelId = GetMeshKernelId(); + + meshkernelapi::Mesh2D mesh2d{}; + [[maybe_unused]] int errorCode = mkernel_mesh2d_get_dimensions(meshKernelId, mesh2d); + + // Get copy of original mesh + // Will be used later to compare with translated mesh + std::vector edge_faces(mesh2d.num_edges * 2); + std::vector edge_nodes(mesh2d.num_edges * 2); + std::vector face_nodes(mesh2d.num_face_nodes); + std::vector face_edges(mesh2d.num_face_nodes); + std::vector nodes_per_face(mesh2d.num_faces); + std::vector node_x(mesh2d.num_nodes); + std::vector node_y(mesh2d.num_nodes); + std::vector edge_x(mesh2d.num_edges); + std::vector edge_y(mesh2d.num_edges); + std::vector face_x(mesh2d.num_faces); + std::vector face_y(mesh2d.num_faces); + + mesh2d.edge_faces = edge_faces.data(); + mesh2d.edge_nodes = edge_nodes.data(); + mesh2d.face_nodes = face_nodes.data(); + mesh2d.face_edges = face_edges.data(); + mesh2d.nodes_per_face = nodes_per_face.data(); + mesh2d.node_x = node_x.data(); + mesh2d.node_y = node_y.data(); + mesh2d.edge_x = edge_x.data(); + mesh2d.edge_y = edge_y.data(); + mesh2d.face_x = face_x.data(); + mesh2d.face_y = face_y.data(); + errorCode = mkernel_mesh2d_get_data(meshKernelId, mesh2d); + + double centreX = 10.0; + double centreY = 5.0; + double angle = 45.0; + + meshkernelapi::mkernel_mesh2d_rotate(meshKernelId, centreX, centreY, angle); + + meshkernelapi::Mesh2D mesh2dTranslated{}; + errorCode = mkernel_mesh2d_get_dimensions(meshKernelId, mesh2dTranslated); + + // Data for translated mesh + std::vector edge_faces_t(mesh2dTranslated.num_edges * 2); + std::vector edge_nodes_t(mesh2dTranslated.num_edges * 2); + std::vector face_nodes_t(mesh2dTranslated.num_face_nodes); + std::vector face_edges_t(mesh2dTranslated.num_face_nodes); + std::vector nodes_per_face_t(mesh2dTranslated.num_faces); + std::vector node_x_t(mesh2dTranslated.num_nodes); + std::vector node_y_t(mesh2dTranslated.num_nodes); + std::vector edge_x_t(mesh2dTranslated.num_edges); + std::vector edge_y_t(mesh2dTranslated.num_edges); + std::vector face_x_t(mesh2dTranslated.num_faces); + std::vector face_y_t(mesh2dTranslated.num_faces); + + mesh2dTranslated.edge_faces = edge_faces_t.data(); + mesh2dTranslated.edge_nodes = edge_nodes_t.data(); + mesh2dTranslated.face_nodes = face_nodes_t.data(); + mesh2dTranslated.face_edges = face_edges_t.data(); + mesh2dTranslated.nodes_per_face = nodes_per_face_t.data(); + mesh2dTranslated.node_x = node_x_t.data(); + mesh2dTranslated.node_y = node_y_t.data(); + mesh2dTranslated.edge_x = edge_x_t.data(); + mesh2dTranslated.edge_y = edge_y_t.data(); + mesh2dTranslated.face_x = face_x_t.data(); + mesh2dTranslated.face_y = face_y_t.data(); + + errorCode = mkernel_mesh2d_get_data(meshKernelId, mesh2dTranslated); + + const double tolerance = 1e-12; + + meshkernel::RigidBodyTransformation transformation; + + // First translate + meshkernel::Translation translation(meshkernel::Vector(-centreX, -centreY)); + transformation.compose(translation); + + // Second rotate + transformation.compose(meshkernel::Rotation(angle)); + + // Third translate back + translation.reset(meshkernel::Vector(centreX, centreY)); + transformation.compose(translation); + + for (int i = 0; i < mesh2d.num_nodes; ++i) + { + // Apply expected transformed to the original mesh node + meshkernel::Point expected = transformation(meshkernel::Point(mesh2d.node_x[i], mesh2d.node_y[i])); + + // The compare the expected node with the transformed mesh node + EXPECT_NEAR(expected.x, mesh2dTranslated.node_x[i], tolerance); + EXPECT_NEAR(expected.y, mesh2dTranslated.node_y[i], tolerance); + } +}