Merge pull request #1440 from BruegelN/fix-remesh_along_isoline

solve fragmentation issue with remesh_along_isoline (cleaned up)

include/igl/remesh_along_isoline.cpp

@@ -7,6 +7,7 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "remesh_along_isoline.h"
 #include "list_to_matrix.h"
+#include <unordered_map>
 
 template <
   typename DerivedV,
@@ -71,6 +72,7 @@ template <
   }
 
   // Loop over each face
+  std::unordered_map<int, std::unordered_map<int, int>> edgeToBirthVert;
   for(int f = 0;f<F.rows();f++)
   {
     bool Psign[2];
@@ -116,18 +118,23 @@ template <
         // Create two new vertices
         for(int i = 0;i<2;i++)
         {
-          const double bci = (isoval - S(F(f,(P[i]+1)%3)))/
-            (S(F(f,P[i]))-S(F(f,(P[i]+1)%3)));
-          vBC.emplace_back(Ucount,F(f,P[i]),bci);
-          vBC.emplace_back(Ucount,F(f,(P[i]+1)%3),1.0-bci);
-          Ucount++;
+          if ((edgeToBirthVert.find(F(f, P[i])) == edgeToBirthVert.end()) || (edgeToBirthVert.at(F(f, P[i])).find(F(f, (P[i] + 1) % 3)) == edgeToBirthVert.at(F(f, P[i])).end()))
+          {
+            const double bci = (isoval - S(F(f,(P[i]+1)%3)))/
+              (S(F(f,P[i]))-S(F(f,(P[i]+1)%3)));
+            vBC.emplace_back(Ucount,F(f,P[i]),bci);
+            vBC.emplace_back(Ucount,F(f,(P[i]+1)%3),1.0-bci);
+            edgeToBirthVert[F(f, P[i])][F(f, (P[i] + 1) % 3)] = Ucount;
+            edgeToBirthVert[F(f, (P[i] + 1) % 3)][F(f, P[i])] = Ucount;
+            Ucount++;
+          }
         }
         const int v0 = F(f,P[0]);
-        const int v01 = Ucount-2;
         assert(((P[0]+1)%3) == P[1]);
         const int v1 = F(f,P[1]);
-        const int v12 = Ucount-1;
         const int v2 = F(f,(P[1]+1)%3);
+        const int v01 = edgeToBirthVert[v0][v1];
+        const int v12 = edgeToBirthVert[v1][v2];
         // v0
         // |  \
         // |   \
@@ -156,5 +163,5 @@ template <
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
-template void igl::remesh_along_isoline<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, double, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, 1, 0, -1, 1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> > const&, Eigen::Matrix<double, -1, 1, 0, -1, 1>::Scalar, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::SparseMatrix<double, 0, int>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
+template void igl::remesh_along_isoline<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, double, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> > const&, Eigen::Matrix<double, -1, 1, 0, -1, 1>::Scalar, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::SparseMatrix<double, 0, int>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
 #endif

tests/include/igl/remesh_along_isoline.cpp

@@ -0,0 +1,32 @@
+#include <test_common.h>
+#include <igl/remesh_along_isoline.h>
+#include <igl/facet_components.h>
+#include <Eigen/Sparse>
+
+TEST_CASE("remesh_along_isoline: triangle_mesh", "[igl]")
+{
+  Eigen::MatrixXd V;
+  Eigen::MatrixXi F;
+  //This is a cube of dimensions 1.0x1.0x1.0
+  igl::read_triangle_mesh(test_common::data_path("cube.obj"), V, F);
+  const double mean_z = V.col(2).mean();
+
+  Eigen::VectorXi C;
+  igl::facet_components(F, C);
+  const int fc_count_before = C.maxCoeff();
+
+  Eigen::VectorXd S = V.col(2);
+  Eigen::MatrixXd U;
+  Eigen::MatrixXi G;
+  Eigen::VectorXd SU;
+  Eigen::VectorXi J;
+  Eigen::SparseMatrix<double> BC;
+  Eigen::VectorXi L;
+  igl::remesh_along_isoline(V,F,S,mean_z,U,G,SU,J,BC,L);
+
+  igl::facet_components(G, C);
+  const int fc_count_after = C.maxCoeff();
+
+  // number of face connected components should not change
+  REQUIRE(fc_count_before == fc_count_after);
+}