adaptations for cylindical beams

doc/gh_DFMergeAssemblies.rst

@@ -0,0 +1,8 @@
+.. image:: ../src/gh/components/DF_merge_assemblies/icon.png
+    :align: left
+    :width: 40px
+
+``DFMergeAssemblies`` component
+===============================
+
+.. ghcomponent_to_rst:: ../src/gh/components/DF_merge_assemblies

doc/gh_components.rst

@@ -81,6 +81,10 @@ DF has a Grasshopper_ plugin with a set of components that allows the user to in
      - .. image:: ../src/gh/components/DF_remove_statistical_outliers/icon.png
      - `gh_DFRemoveStatisticalOutliers <gh_DFRemoveStatisticalOutliers.html>`_
 
+   * - .. image:: ../src/gh/components/DF_merge_assemblies/icon.png
+     - `gh_DFMergeAssemblies <gh_DFMergeAssemblies.html>`_
+     - 
+     - 
 
 
 .. toctree::
@@ -114,4 +118,5 @@ DF has a Grasshopper_ plugin with a set of components that allows the user to in
    gh_DFRemoveBeam
    gh_DFColorizeCloud
    gh_DFBrepToCloud
-   gh_DFRemoveStatisticalOutliers
+   gh_DFRemoveStatisticalOutliers
+   gh_DFMergeAssemblies

src/diffCheck/geometry/DFMesh.cc

@@ -112,31 +112,67 @@ namespace diffCheck::geometry
             Eigen::Vector3d v1 = this->Vertices[triangle[1]];
             Eigen::Vector3d v2 = this->Vertices[triangle[2]];
             Eigen::Vector3d n = (v1 - v0).cross(v2 - v0);
-            double normOfNormal = n.norm();
             n.normalize();
 
-            Eigen::Vector3d projectedPoint = point - n * (n.dot(point - v0)) ;
+            // Project the point onto the plane of the triangle
+            Eigen::Vector3d projectedPoint = point - n * (n.dot(point - v0));
 
-            double referenceTriangleArea = normOfNormal*0.5;
-            Eigen::Vector3d n1 = (v1 - v0).cross(projectedPoint - v0);
-            double area1 = n1.norm()*0.5;
-            Eigen::Vector3d n2 = (v2 - v1).cross(projectedPoint - v1);
-            double area2 = n2.norm()*0.5;
-            Eigen::Vector3d n3 = (v0 - v2).cross(projectedPoint - v2);
-            double area3 = n3.norm()*0.5;
-            double res = (area1 + area2 + area3 - referenceTriangleArea) / referenceTriangleArea;
+            // Compute vectors
+            Eigen::Vector3d v0v1 = v1 - v0;
+            Eigen::Vector3d v0v2 = v2 - v0;
+            Eigen::Vector3d v0p = projectedPoint - v0;
 
-            // arbitrary value to avoid false positives (points that, when projected on the triangle, are in it, but that are actually located too far from the mesh to actually belong to it)
-            double maxProjectionDistance = std::min({(v1 - v0).norm(), (v2 - v1).norm(), (v0 - v2).norm()}) / 2;
+            // Compute dot products
+            double dot00 = v0v2.dot(v0v2);
+            double dot01 = v0v2.dot(v0v1);
+            double dot02 = v0v2.dot(v0p);
+            double dot11 = v0v1.dot(v0v1);
+            double dot12 = v0v1.dot(v0p);
 
-            if (std::abs(res) < associationThreshold && (projectedPoint - point).norm() < maxProjectionDistance)
+            // Compute barycentric coordinates
+            double invDenom = 1.0 / (dot00 * dot11 - dot01 * dot01);
+            double u = (dot11 * dot02 - dot01 * dot12) * invDenom;
+            double v = (dot00 * dot12 - dot01 * dot02) * invDenom;
+
+            // Check if point is in triangle
+            if ((u >= -associationThreshold) && (v >= -associationThreshold) && (u + v <= 1 + associationThreshold))
             {
-                return true;
+                // Check if the point is close enough to the face
+                double maxProjectionDistance = std::min({(v1 - v0).norm(), (v2 - v1).norm(), (v0 - v2).norm()}) ;
+                if ((projectedPoint - point).norm() < maxProjectionDistance)
+                {
+                    return true;
+                }
             }
         }
         return false;
     }
 
+    std::tuple<Eigen::Vector3d, Eigen::Vector3d> DFMesh::ComputeOBBCenterAndAxis()
+    {
+        Eigen::Vector3d center = Eigen::Vector3d::Zero();
+        Eigen::Vector3d axis = Eigen::Vector3d::Zero();
+
+        std::vector<Eigen::Vector3d> tightBoundingBox = this->GetTightBoundingBox();
+
+        Eigen::Vector3d deltaFirstDir = tightBoundingBox[1] - tightBoundingBox[0];
+        Eigen::Vector3d deltaSecondDir = tightBoundingBox[2] - tightBoundingBox[0];
+        Eigen::Vector3d deltaThirdDir = tightBoundingBox[3] - tightBoundingBox[0];
+
+        for (Eigen::Vector3d direction : {deltaFirstDir, deltaSecondDir, deltaThirdDir})
+        {
+            if (direction.norm() > axis.norm())
+            {
+                axis = direction;
+            }
+        }
+        axis.normalize();
+
+        center = tightBoundingBox[0] + deltaFirstDir/2 + deltaSecondDir/2 + deltaThirdDir/2;
+
+        return std::make_tuple(center, axis);
+    }
+
     void DFMesh::LoadFromPLY(const std::string &path)
     {
         std::shared_ptr<diffCheck::geometry::DFMesh> tempMesh_ptr = diffCheck::io::ReadPLYMeshFromFile(path);

src/diffCheck/geometry/DFMesh.hh

@@ -86,6 +86,13 @@ namespace diffCheck::geometry
          */
         bool IsPointOnFace(Eigen::Vector3d point, double associationThreshold = 0.1);
 
+        /**
+         * @brief Get the center and main axis of oriented boundung box of the mesh. It was developped for the cylinder case, but can be used for other shapes.
+         * 
+         * @return std::tuple<Eigen::Vector3d, Eigen::Vector3d> the first element is the center of the obb of the mesh, the second element is the main axis of the obb of the mesh
+         */
+        std::tuple<Eigen::Vector3d, Eigen::Vector3d> ComputeOBBCenterAndAxis();
+
     public:  ///< I/O loader
         /**
          * @brief Read a mesh from a file