Fixed the changing WL matrix bug in witness complex. Added is_witness_complex test

git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/witness@615 636b058d-ea47-450e-bf9e-a15bfbe3eedb

Former-commit-id: 60917e75766b53d08240b510c0508e7d781b56d6

1 files changed, 68 insertions, 25 deletions

diff --git a/src/Witness_complex/include/gudhi/Witness_complex.h b/src/Witness_complex/include/gudhi/Witness_complex.h
index 04fcc98f..201d6525 100644
--- a/src/Witness_complex/include/gudhi/Witness_complex.h
+++ b/src/Witness_complex/include/gudhi/Witness_complex.h
@@ -140,7 +140,6 @@ namespace Gudhi {
     //void witness_complex(std::vector< std::vector< Vertex_handle > > & knn)
     {
       std::cout << "**Start the procedure witness_complex" << std::endl;
-      int k=2; /* current dimension in iterative construction */
       //Construction of the active witness list
       int nbW = knn.size();
       //int nbL = knn.at(0).size();
@@ -157,13 +156,14 @@ namespace Gudhi {
         // by doing it we don't assume that landmarks are necessarily witnesses themselves anymore
         counter++;
         vv = {i};
-        /* TODO Filtration */
         returnValue = insert_simplex(vv, Filtration_value(0.0));
         /* TODO Error if not inserted : normally no need here though*/
       }
+      int k=1; /* current dimension in iterative construction */
       //std::cout << "Successfully added landmarks" << std::endl;
       // PRINT2
       //print_sc(root()); std::cout << std::endl;
+      /*
       int u,v;     // two extremities of an edge
       int count = 0;
       if (nbL > 1) // if the supposed dimension of the complex is >0
@@ -202,39 +202,37 @@ namespace Gudhi {
               active_w.push_back(i);
             }
         }
-      std::cout << "k=1, active witnesses: " << active_w.size() << std::endl;
+      */
+      for (int i=0; i != nbW; ++i)
+        active_w.push_back(i);
+      std::cout << "k=0, active witnesses: " << active_w.size() << std::endl;
       //std::cout << "Successfully added edges" << std::endl;
-      count_good = {0,0};
-      count_bad = {0,0};
+      count_good = {0};
+      count_bad = {0};
       int D = knn[0].size();
       while (!active_w.empty() && k < D )
         {
 	  count_good.push_back(0);
 	  count_bad.push_back(0);
-          count++;
           //std::cout << "Started the step k=" << k << std::endl;
           typename ActiveWitnessList::iterator it = active_w.begin();
           while (it != active_w.end())
             {
               typeVectorVertex simplex_vector;
               /* THE INSERTION: Checking if all the subfaces are in the simplex tree*/
-              // First sort the first k landmarks
-              VertexHandle inserted_vertex = knn[*it][k];
-              bool ok = all_faces_in(knn, *it, k, inserted_vertex);
+              bool ok = all_faces_in(knn, *it, k);
               if (ok)
                 {
                   for (int i = 0; i != k+1; ++i)
                     simplex_vector.push_back(knn[*it][i]);
                   returnValue = insert_simplex(simplex_vector,0.0);
-                  if (returnValue.second)
-                    count++;
                   it++;
                 }
               else
                 active_w.erase(it++); //First increase the iterator and then erase the previous element
             }
 	  std::cout << "k=" << k << ", active witnesses: " << active_w.size() << std::endl;
-          std::cout << "** k=" << k << ", num_simplices: " <<count << std::endl;
+          //std::cout << "** k=" << k << ", num_simplices: " <<count << std::endl;
           k++;
         }
       //print_sc(root()); std::cout << std::endl;
@@ -324,7 +322,7 @@ private:
      *  inserted_vertex is the handle of the (k+1)-th vertex witnessed by witness_id
      */
     template <typename KNearestNeighbours>
-    bool all_faces_in(KNearestNeighbours &knn, int witness_id, int k, VertexHandle inserted_vertex)
+    bool all_faces_in(KNearestNeighbours &knn, int witness_id, int k)
     {
       //std::cout << "All face in with the landmark " << inserted_vertex << std::endl;
       std::vector< VertexHandle > facet;
@@ -332,20 +330,17 @@ private:
       // CHECK ALL THE FACETS
       for (int i = 0; i != k+1; ++i)
         {
-          if (knn[witness_id][i] != inserted_vertex)
+          facet = {};
+          for (int j = 0; j != k+1; ++j)
             {
-              facet = {};
-              for (int j = 0; j != k+1; ++j)
+              if (j != i)
                 {
-                  if (j != i)
-                    {
-                      facet.push_back(knn[witness_id][j]);
-                    }
-                }//endfor
-              if (find(facet) == null_simplex())
-                return false;
-              //std::cout << "++++ finished loop safely\n";
-            }//endif
+                  facet.push_back(knn[witness_id][j]);
+                }
+            }//endfor
+          if (find(facet) == null_simplex())
+            return false;
+          //std::cout << "++++ finished loop safely\n";
         } //endfor
       return true;
     }
@@ -1059,6 +1054,54 @@ bool complex_is_pseudomanifold(int dimension)
         }
     }
 
+  public:
+    /** \brief Verification if every simplex in the complex is witnessed
+     */
+    template< class KNearestNeighbors >
+    bool is_witness_complex(KNearestNeighbors WL)
+    {
+      //bool final_result = true;
+      for (Simplex_handle sh: complex_simplex_range())
+        {
+          bool is_witnessed = false;
+          typeVectorVertex simplex;
+          int nbV = 0; //number of verticed in the simplex
+          for (int v: simplex_vertex_range(sh))
+            simplex.push_back(v);
+          nbV = simplex.size();
+          for (typeVectorVertex w: WL)
+            {
+              bool has_vertices = true;  
+              for (int v: simplex)
+                if (std::find(w.begin(), w.begin()+nbV, v) == w.begin()+nbV)
+                  {
+                    has_vertices = false;
+                    //break;
+                  }
+              if (has_vertices)
+                {
+                  is_witnessed = true;
+                  std::cout << "The simplex ";
+                  print_vector(simplex);
+                  std::cout << " is witnessed by the witness ";
+                  print_vector(w);
+                  std::cout << std::endl;
+                  break;
+                }
+            }
+          if (!is_witnessed)
+            {
+              std::cout << "The following simplex is not witnessed ";
+              print_vector(simplex);
+              std::cout << std::endl;
+              assert(is_witnessed);
+              return false;
+            }
+        }
+      return true; // Arrive here if the not_witnessed check failed all the time
+    }
+
+    
 }; //class Witness_complex