upsample, limit faces to subset of indices, reorder vertex and face list so face'd indices are first, all pairs distances

Former-commit-id: 1e9a117fb46ae79d0b467845558bd0f458f951dd

adjacency_list.h

@@ -1,5 +1,5 @@
-#ifndef IGL_ADJACENCY_MATRIX_H
-#define IGL_ADJACENCY_MATRIX_H
+#ifndef IGL_ADJACENCY_LIST_H
+#define IGL_ADJACENCY_LIST_H
 
 #define EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
 #include <Eigen/Dense>
@@ -58,7 +58,7 @@ inline void igl::adjacency_list(
   }
   
   // Remove duplicates
-  for(int i=0; i<A.size();++i)
+  for(int i=0; i<(int)A.size();++i)
   {
     std::sort(A[i].begin(), A[i].end());
     A[i].erase(std::unique(A[i].begin(), A[i].end()), A[i].end());
@@ -91,7 +91,7 @@ inline void igl::adjacency_list(
       }
 	  }
 	  
-    for(int v=0; v<SR.size();++v)
+    for(int v=0; v<(int)SR.size();++v)
     {
       std::vector<int>& vv = A.at(v);
       std::vector<std::vector<int> >& sr = SR[v];
@@ -99,21 +99,21 @@ inline void igl::adjacency_list(
       std::vector<std::vector<int> > pn = sr;
       
       // Compute previous/next for every element in sr
-      for(int i=0;i<sr.size();++i)
+      for(int i=0;i<(int)sr.size();++i)
       {
         int a = sr[i][0];
         int b = sr[i][1];
         
         // search for previous
         int p = -1;
-        for(int j=0;j<sr.size();++j)
+        for(int j=0;j<(int)sr.size();++j)
           if(sr[j][1] == a)
             p = j;
         pn[i][0] = p;
         
         // search for next
         int n = -1;
-        for(int j=0;j<sr.size();++j)
+        for(int j=0;j<(int)sr.size();++j)
           if(sr[j][0] == b)
             n = j;
         pn[i][1] = n;
@@ -122,14 +122,14 @@ inline void igl::adjacency_list(
       
       // assume manifoldness (look for beginning of a single chain)
       int c = 0;
-      for(int j=0; j<=sr.size();++j)
+      for(int j=0; j<=(int)sr.size();++j)
         if (pn[c][0] != -1)
           c = pn[c][0];
       
       if (pn[c][0] == -1) // border case
       {
         // finally produce the new vv relation
-        for(int j=0; j<sr.size();++j)
+        for(int j=0; j<(int)sr.size();++j)
         {
           vv[j] = sr[c][0];
           if (pn[c][1] != -1)
@@ -140,7 +140,7 @@ inline void igl::adjacency_list(
       else
       {
         // finally produce the new vv relation
-        for(int j=0; j<sr.size();++j)
+        for(int j=0; j<(int)sr.size();++j)
         {
           vv[j] = sr[c][0];
           

all_pairs_distances.h

@@ -0,0 +1,55 @@
+#ifndef IGL_ALL_PAIRS_DISTANCES_H
+#define IGL_ALL_PAIRS_DISTANCES_H
+
+namespace igl
+{
+  // ALL_PAIRS_DISTANCES compute distances between each point i in V and point j
+  // in U
+  // 
+  // D = all_pairs_distances(V,U)
+  // 
+  // Templates:
+  //   Mat  matrix class like MatrixXd
+  // Inputs:
+  //   V  #V by dim list of points
+  //   U  #U by dim list of points
+  //   squared  whether to return squared distances
+  // Outputs:
+  //   D  #V by #U matrix of distances, where D(i,j) gives the distance or
+  //     squareed distance between V(i,:) and U(j,:)
+  // 
+  template <typename Mat>
+  void all_pairs_distances(
+    const Mat & V,
+    const Mat & U,
+    const bool squared, 
+    Mat & D);
+}
+
+// Implementation
+
+template <typename Mat>
+void igl::all_pairs_distances(
+  const Mat & V,
+  const Mat & U,
+  const bool squared,
+  Mat & D)
+{
+  // dimension should be the same
+  assert(V.cols() == U.cols());
+  // resize output
+  D.resize(V.rows(),U.rows());
+  for(int i = 0;i<V.rows();i++)
+  {
+    for(int j=0;j<U.rows();j++)
+    {
+      D(i,j) = (V.row(i)-U.row(j)).array().pow(2).sum();
+      if(!squared)
+      {
+        D(i,j) = sqrt(D(i,j));
+      }
+    }
+  }
+}
+
+#endif

edges.h

@@ -16,18 +16,14 @@ namespace igl
   //   E #E by 2 list of edges in no particular order
   //
   // See also: adjacency_matrix
-  inline void edges(
-    const Eigen::MatrixXi& F, 
-    Eigen::MatrixXi& E);
+  inline void edges( const Eigen::MatrixXi& F, Eigen::MatrixXi& E);
 }
 
 // Implementation
 #include <map>
-#include "adjacency_matrix.h"
+#include <adjacency_matrix.h>
 
-inline void igl::edges(
-  const Eigen::MatrixXi& F, 
-  Eigen::MatrixXi& E)
+inline void igl::edges( const Eigen::MatrixXi& F, Eigen::MatrixXi& E)
 {
   // build adjacency matrix
   Eigen::SparseMatrix<int> A;

faces_first.h

@@ -0,0 +1,110 @@
+#ifndef IGL_FACES_FIRST_H
+#define IGL_FACES_FIRST_H
+namespace igl
+{
+  // FACES_FIRST Reorder vertices so that vertices in face list come before
+  // vertices that don't appear in the face list. This is especially useful if
+  // the face list contains only surface faces and you want surface vertices
+  // listed before internal vertices
+  //
+  // [RV,RT,RF,IM] = faces_first(V,T,F);
+  //
+  // Templates:
+  //   MatV  matrix for vertex positions, e.g. MatrixXd
+  //   MatF  matrix for vertex positions, e.g. MatrixXi
+  //   VecI  matrix for vertex positions, e.g. VectorXi
+  // Input:
+  //  V  # vertices by 3 vertex positions
+  //  F  # faces by 3 list of face indices
+  // Output: 
+  //  RV  # vertices by 3 vertex positions, order such that if the jth vertex is
+  //    some face in F, and the kth vertex is not then j comes before k
+  //  RF  # faces by 3 list of face indices, reindexed to use RV
+  //  IM  # faces by 1 list of indices such that: RF = IM(F) and RT = IM(T)
+  //    and RV(IM,:) = V
+  //
+  template <typename MatV, typename MatF, typename VecI>
+  void faces_first(
+    const MatV & V, 
+    const MatF & F, 
+    MatV & RV, 
+    MatF & RF, 
+    VecI & IM);
+}
+
+// Implementation
+#include <vector>
+#include <Eigen/Dense>
+
+template <typename MatV, typename MatF, typename VecI>
+void igl::faces_first(
+  const MatV & V, 
+  const MatF & F, 
+  MatV & RV, 
+  MatF & RF, 
+  VecI & IM)
+{
+  using namespace std;
+  using namespace Eigen;
+  vector<bool> in_face(V.rows());
+  for(int i = 0; i<F.rows(); i++)
+  {
+    for(int j = 0; j<F.cols(); j++)
+    {
+      in_face[F(i,j)] = true;
+    }
+  }
+  // count number of vertices not in faces
+  int num_in_F = 0;
+  for(int i = 0;i<V.rows();i++)
+  {
+    num_in_F += (in_face[i]?1:0);
+  }
+  // list of unique vertices that occur in F
+  VectorXi U(num_in_F);
+  // list of unique vertices that do not occur in F
+  VectorXi NU(V.rows()-num_in_F);
+  int Ui = 0;
+  int NUi = 0;
+  // loop over vertices
+  for(int i = 0;i<V.rows();i++)
+  {
+    if(in_face[i])
+    {
+      U(Ui) = i;
+      Ui++;
+    }else
+    {
+      NU(NUi) = i;
+      NUi++;
+    }
+  }
+  IM.resize(V.rows());
+  // reindex vertices that occur in faces to be first
+  for(int i = 0;i<U.size();i++)
+  {
+    IM(U(i)) = i;
+  }
+  // reindex vertices that do not occur in faces to come after those that do
+  for(int i = 0;i<NU.size();i++)
+  {
+    IM(NU(i)) = i+U.size();
+  }
+  RF.resize(F.rows(),F.cols());
+  // Reindex faces
+  for(int i = 0; i<F.rows(); i++)
+  {
+    for(int j = 0; j<F.cols(); j++)
+    {
+      RF(i,j) = IM(F(i,j));
+    }
+  }
+  RV.resize(V.rows(),V.cols());
+  // Reorder vertices
+  for(int i = 0;i<V.rows();i++)
+  {
+    RV.row(IM(i)) = V.row(i);
+  }
+}
+
+#endif

ismanifold.h

@@ -32,7 +32,7 @@ namespace igl
       }
     std::sort(TTT.begin(),TTT.end());
     
-    for(int i=2;i<TTT.size();++i)
+    for(int i=2;i<(int)TTT.size();++i)
     {
       std::vector<int>& r1 = TTT[i-2];
       std::vector<int>& r2 = TTT[i-1];

limit_faces.h

@@ -0,0 +1,82 @@
+#ifndef IGL_LIMIT_FACES
+#define IGL_LIMIT_FACES
+namespace igl
+{
+  // LIMIT_FACES limit given faces F to those which contain (only) indices found
+  // in L.
+  //
+  // [LF] = limit_faces(F,L,exclusive);
+  // [LF,in] = limit_faces(F,L,exclusive);
+  //
+  // Templates:
+  //   MatF matrix type of faces, matrixXi
+  //   VecL  matrix type of vertex indices, VectorXi
+  // Inputs:
+  //   F  #F by 3 list of face indices
+  //   L  #L by 1 list of allowed indices
+  //   exclusive  flag specifying whether a face is included only if all its
+  //     indices are in L, default is false
+  // Outputs:
+  //   LF  #LF by 3 list of remaining faces after limiting
+  //   in  #F list of whether given face was included
+  //
+  template <typename MatF, typename VecL>
+  void limit_faces(
+    const MatF & F, 
+    const VecL & L, 
+    const bool exclusive,
+    MatF & LF);
+}
+
+// Implementation
+#include <vector>
+
+template <typename MatF, typename VecL>
+void igl::limit_faces(
+  const MatF & F, 
+  const VecL & L, 
+  const bool exclusive,
+  MatF & LF)
+{
+  using namespace std;
+  using namespace Eigen;
+  vector<bool> in(F.rows(),false);
+  int num_in = 0;
+  // loop over faces
+  for(int i = 0;i<F.rows();i++)
+  {
+    bool all = true;
+    bool any = false;
+    for(int j = 0;j<F.cols();j++)
+    {
+      bool found = false;
+      // loop over L
+      for(int l = 0;l<L.size();l++)
+      {
+        if(F(i,j) == L(l))
+        {
+          found = true;
+          break;
+        }
+      }
+      any |= found;
+      all &= found;
+    }
+    in[i] = (exclusive?all:any);
+    num_in += (in[i]?1:0);
+  }
+
+  LF.resize(num_in,F.cols());
+  // loop over faces
+  int lfi = 0;
+  for(int i = 0;i<F.rows();i++)
+  {
+    if(in[i])
+    {
+      LF.row(lfi) = F.row(i);
+      lfi++;
+    }
+  }
+}
+
+#endif

material_colors.h

@@ -11,4 +11,5 @@ namespace igl
   const float CYAN_AMBIENT[4] =   {  59.0/255.0, 68.0/255.0,255.0/255.0,1.0f };
   const float CYAN_DIFFUSE[4] =   {  94.0/255.0,185.0/255.0,238.0/255.0,1.0f };
   const float CYAN_SPECULAR[4] =   { 163.0/255.0,221.0/255.0,255.0/255.0,1.0f };
+  const float DENIS_PURPLE_DIFFUSE[4] =   { 80.0/255.0,64.0/255.0,255.0/255.0,1.0f };
 }

tt.h

@@ -39,7 +39,7 @@ namespace igl
   {
     TT = Eigen::MatrixXi::Constant((int)(F.rows()),3,-1);
     
-    for(int i=1;i<TTT.size();++i)
+    for(int i=1;i<(int)TTT.size();++i)
     {
       std::vector<int>& r1 = TTT[i-1];
       std::vector<int>& r2 = TTT[i];
@@ -56,7 +56,7 @@ namespace igl
   {
     TTi = Eigen::MatrixXi::Constant((int)(F.rows()),3,-1);
     
-    for(int i=1;i<TTT.size();++i)
+    for(int i=1;i<(int)TTT.size();++i)
     {
       std::vector<int>& r1 = TTT[i-1];
       std::vector<int>& r2 = TTT[i];

upsample.h

@@ -0,0 +1,115 @@
+#ifndef IGL_UPSAMPLE_H
+#define IGL_UPSAMPLE_H
+namespace igl
+{
+  // Subdivide a mesh without moving vertices: loop subdivision but odd
+  // vertices stay put and even vertices are just edge midpoints
+  // 
+  // Templates:
+  //   MatV  matrix for vertex positions, e.g. MatrixXd
+  //   MatF  matrix for vertex positions, e.g. MatrixXi
+  // Inputs:
+  //   V  #V by dim  mesh vertices
+  //   F  #F by 3  mesh triangles
+  // Outputs:
+  //   NV new vertex positions, V is guaranteed to be at top
+  //   NF new list of face indices
+  //
+  // NOTE: V should not be the same as NV,
+  // NOTE: F should not be the same as NF, use other proto
+  template <typename MatV, typename MatF>
+  void upsample( const MatV & V, const MatF & F, MatV & NV, MatF & NF);
+  // Virtually in place wrapper
+  template <typename MatV, typename MatF>
+  void upsample( MatV & V,MatF & F);
+}
+
+// Implementation
+#include <tt.h>
+#include <adjacency_list.h>
+#include <Eigen/Dense>
+
+template <typename MatV, typename MatF>
+void igl::upsample( const MatV & V, const MatF & F, MatV & NV, MatF & NF)
+{
+  // Use "in place" wrapper instead
+  assert(&V != &NV);
+  assert(&F != &NF);
+  using namespace igl;
+  using namespace std;
+  using namespace Eigen;
+
+  MatF FF, FFi;
+  tt<double>(V,F,FF,FFi);
+
+  // TODO: Cache optimization missing from here, it is a mess
+  
+  // Compute the number and positions of the vertices to insert (on edges)
+  MatF NI = MatF::Constant(FF.rows(),FF.cols(),-1);
+  int counter = 0;
+  
+  for(int i=0;i<FF.rows();++i)
+  {
+    for(int j=0;j<3;++j)
+    {
+      if(NI(i,j) == -1)
+      {
+        NI(i,j) = counter;
+        if (FF(i,j) != -1) // If it is not a border
+          NI(FF(i,j),FFi(i,j)) = counter;
+        ++counter;
+      }
+    }
+  }
+
+  int n_odd = V.rows();
+  int n_even = counter;
+
+  Eigen::DynamicSparseMatrix<double> SUBD(V.rows()+n_even,V.rows());
+  SUBD.reserve(15 * (V.rows()+n_even));
+  
+  // Preallocate NV and NF
+  NV = MatV(V.rows()+n_even,V.cols());
+  NF = MatF(F.rows()*4,3);
+  
+  // Fill the odd vertices position
+  NV.block(0,0,V.rows(),V.cols()) = V;
+
+  // Fill the even vertices position
+  for(int i=0;i<FF.rows();++i)
+  {
+    for(int j=0;j<3;++j)
+    {
+      NV.row(NI(i,j) + n_odd) = 0.5 * V.row(F(i,j)) + 0.5 * V.row(F(i,(j+1)%3));
+    }
+  }
+
+  // Build the new topology (Every face is replaced by four)
+  for(int i=0; i<F.rows();++i)
+  {
+    VectorXi VI(6);
+    VI << F(i,0), F(i,1), F(i,2), NI(i,0) + n_odd, NI(i,1) + n_odd, NI(i,2) + n_odd;
+    
+    VectorXi f0(3), f1(3), f2(3), f3(3);
+    f0 << VI(0), VI(3), VI(5);
+    f1 << VI(1), VI(4), VI(3);
+    f2 << VI(3), VI(4), VI(5);
+    f3 << VI(4), VI(2), VI(5);
+    
+    NF.row((i*4)+0) = f0;
+    NF.row((i*4)+1) = f1;
+    NF.row((i*4)+2) = f2;
+    NF.row((i*4)+3) = f3;
+  }
+  
+}
+
+template <typename MatV, typename MatF>
+void igl::upsample( MatV & V,MatF & F)
+{
+  const MatV V_copy = V;
+  const MatF F_copy = F;
+  return upsample(V_copy,F_copy,V,F);
+}
+
+#endif

vf.h

@@ -1,9 +1,8 @@
 #ifndef IGL_VF_H
 #define IGL_VF_H
 
-#define EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
 #include <Eigen/Dense>
-#include <Eigen/Sparse>
+#include <vector>
 
 namespace igl 
 {
@@ -20,7 +19,7 @@ namespace igl
   inline void vf( 
     const Eigen::Matrix<S, Eigen::Dynamic, Eigen::Dynamic> & V, 
     const Eigen::MatrixXi & F, 
-    vector<vector<T> >& VF, vector<vector<T> >& VFi);
+    std::vector<std::vector<T> >& VF, std::vector<std::vector<T> >& VFi);
 }
 
 // Implementation
@@ -30,7 +29,7 @@ template <typename T, typename S>
 inline void igl::vf(
   const Eigen::Matrix<S, Eigen::Dynamic, Eigen::Dynamic> & V, 
   const Eigen::MatrixXi & F, 
-  vector<vector<T> >& VF, vector<vector<T> >& VFi)
+  std::vector<std::vector<T> >& VF, std::vector<std::vector<T> >& VFi)
 {
   VF.clear();
   VFi.clear();