duplicate vertex support in selfintersect

Former-commit-id: 9454ba09ccffcc0f2c6ab5473642fbbf0d8416cd

examples/ambient-occlusion-mex/compile.m

@@ -4,6 +4,6 @@ mex -v -o ambient_occlusion -DMEX -largeArrayDims ...
   CXXFLAGS="\$CXXFLAGS -m64 -msse4.2 -fopenmp" ...
   -I/usr/local/igl/libigl/external/embree/embree/ ...
   -I/usr/local/igl/libigl/external/embree/ ...
-  -L/usr/local/igl/libigl/external/embree/bin -lembree -lsys ... 
+  -L/usr/local/igl/libigl/external/embree/build -lembree -lsys ... 
   mexFunction.cpp parse_rhs.cpp;
 

include/igl/WindingNumberTree.h

@@ -322,6 +322,7 @@ inline double igl::WindingNumberTree<Point>::winding_number_boundary(const Point
   using namespace std;
 
   double w = 0;
+  // `cap` is already flipped inside out, so we don't need to flip sign of w
   igl::winding_number_3(
     V.data(),
     V.rows(),

include/igl/cgal/SelfIntersectMesh.h

@@ -85,6 +85,8 @@ namespace igl
     public:
       // Constructs (VV,FF) a new mesh with self-intersections of (V,F)
       // subdivided
+      //
+      // See also: selfintersect.h
       inline SelfIntersectMesh(
         const Eigen::MatrixXd & V,
         const Eigen::MatrixXi & F,
@@ -92,7 +94,9 @@ namespace igl
         Eigen::MatrixXd & VV,
         Eigen::MatrixXi & FF,
         Eigen::MatrixXi & IF,
-        Eigen::VectorXi & J);
+        Eigen::VectorXi & J,
+        Eigen::VectorXi & IM
+        );
     private:
       // Helper function to mark a face as offensive
       //
@@ -252,7 +256,8 @@ inline igl::SelfIntersectMesh<Kernel>::SelfIntersectMesh(
   Eigen::MatrixXd & VV,
   Eigen::MatrixXi & FF,
   Eigen::MatrixXi & IF,
-  Eigen::VectorXi & J):
+  Eigen::VectorXi & J,
+  Eigen::VectorXi & IM):
   V(V),
   F(F),
   count(0),
@@ -357,8 +362,8 @@ inline igl::SelfIntersectMesh<Kernel>::SelfIntersectMesh(
           //  P[o].to_3d(vit->point())<<endl;
 #ifndef NDEBUG
           // I want to be sure that the original corners really show up as the
-          // original corners of the CDT. I.e. I don't trust CGAL to maintain the
-          // order
+          // original corners of the CDT. I.e. I don't trust CGAL to maintain
+          // the order
           assert(T[f].vertex(i) == P[o].to_3d(vit->point()));
 #endif
           // For first three, use original index in F
@@ -496,6 +501,38 @@ inline igl::SelfIntersectMesh<Kernel>::SelfIntersectMesh(
       i++;
     }
   }
+  IM.resize(VV.rows(),1);
+  map<Point_3,int> vv2i;
+  // Safe to check for duplicates using double for original vertices: if
+  // incoming reps are different then the points are unique.
+  for(int v = 0;v<V.rows();v++)
+  {
+    const Point_3 p(V(v,0),V(v,1),V(v,2));
+    if(vv2i.count(p)==0)
+    {
+      vv2i[p] = v;
+    }
+    assert(vv2i.count(p) == 1);
+    IM(v) = vv2i[p];
+  }
+  // Must check for duplicates of new vertices using exact.
+  {
+    int v = V.rows();
+    for(
+      typename list<Point_3>::const_iterator nvit = NV.begin();
+      nvit != NV.end();
+      nvit++)
+    {
+      const Point_3 & p = *nvit;
+      if(vv2i.count(p)==0)
+      {
+        vv2i[p] = v;
+      }
+      assert(vv2i.count(p) == 1);
+      IM(v) = vv2i[p];
+      v++;
+    }
+  }
 
   // Q: Does this give the same result as TETGEN?
   // A: For the cow and beast, yes.

include/igl/cgal/selfintersect.cpp

@@ -17,7 +17,8 @@ IGL_INLINE void igl::selfintersect(
   Eigen::MatrixXd & VV,
   Eigen::MatrixXi & FF,
   Eigen::MatrixXi & IF,
-  Eigen::VectorXi & J)
+  Eigen::VectorXi & J,
+  Eigen::VectorXi & IM)
 {
   using namespace std;
   if(params.detect_only)
@@ -35,7 +36,7 @@ IGL_INLINE void igl::selfintersect(
 //#endif
 
     typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
-    SelfIntersectMesh<Kernel> SIM = SelfIntersectMesh<Kernel>(V,F,params,VV,FF,IF,J);
+    SelfIntersectMesh<Kernel> SIM = SelfIntersectMesh<Kernel>(V,F,params,VV,FF,IF,J,IM);
 
 //#ifdef __APPLE__
 //    signal(SIGFPE,SIG_DFL);
@@ -44,6 +45,6 @@ IGL_INLINE void igl::selfintersect(
   }else
   {
     typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
-    SelfIntersectMesh<Kernel> SIM = SelfIntersectMesh<Kernel>(V,F,params,VV,FF,IF,J);
+    SelfIntersectMesh<Kernel> SIM = SelfIntersectMesh<Kernel>(V,F,params,VV,FF,IF,J,IM);
   }
 }

include/igl/cgal/selfintersect.h

@@ -47,6 +47,12 @@ namespace igl
   //   IF  #intersecting face pairs by 2  list of intersecting face pairs,
   //     indexing F
   //   J  #FF list of indices into F denoting birth triangle
+  //   IM  #VV list of indices into VV of unique vertices.
+  //
+  // Known bugs: If an existing edge in (V,F) lies exactly on another face then
+  // any resulting additional vertices along that edge may not get properly
+  // connected so that the output mesh has the same global topology. This is
+  // because 
   IGL_INLINE void selfintersect(
     const Eigen::MatrixXd & V,
     const Eigen::MatrixXi & F,
@@ -54,7 +60,8 @@ namespace igl
     Eigen::MatrixXd & VV,
     Eigen::MatrixXi & FF,
     Eigen::MatrixXi & IF,
-    Eigen::VectorXi & J);
+    Eigen::VectorXi & J,
+    Eigen::VectorXi & IM);
 }
 
 #ifdef IGL_HEADER_ONLY