multiple fix

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

Former-commit-id: 0bcb4f121da2714128c3ae8ab0da2000e815d339

2 files changed, 134 insertions, 89 deletions

diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree.h b/src/Simplex_tree/include/gudhi/Simplex_tree.h
index 241b692f..a355c334 100644
--- a/src/Simplex_tree/include/gudhi/Simplex_tree.h
+++ b/src/Simplex_tree/include/gudhi/Simplex_tree.h
@@ -156,6 +156,12 @@ class Simplex_tree {
   typedef Simplex_tree_simplex_vertex_iterator<Simplex_tree> Simplex_vertex_iterator;
   /** \brief Range over the vertices of a simplex. */
   typedef boost::iterator_range<Simplex_vertex_iterator> Simplex_vertex_range;
+  /** \brief Range over the cofaces of a simplex. */
+  /** \brief Iterator over the cofaces of a simplex.
+   *
+   * 'value_type' is Vertex_handle. */
+  typedef typename std::vector<Simplex_handle>::iterator Coface_simplex_iterator;
+  typedef boost::iterator_range<Coface_simplex_iterator> Coface_simplex_range;
   /** \brief Iterator over the simplices of the boundary of a simplex.
    *
    * 'value_type' is Simplex_handle. */
@@ -187,8 +193,7 @@ class Simplex_tree {
   /** \name Range and iterator methods
    * @{ */
 
-  /** \brief Returns a range over the vertices of the simplicial complex.
-   *
+  /** \brief Returns a range over the vertices of the simplicial complex.  *
    * The order is increasing according to < on Vertex_handles.*/
   Complex_vertex_range complex_vertex_range() {
     return Complex_vertex_range(
@@ -206,6 +211,7 @@ class Simplex_tree {
                                  Complex_simplex_iterator());
   }
 
+
   /** \brief Returns a range over the simplices of the dim-skeleton of the simplicial complex.
    *
    * The \f$d\f$-skeleton of a simplicial complex \f$\mathbf{K}\f$ is the simplicial complex containing the
@@ -252,10 +258,12 @@ class Simplex_tree {
    * equal to \f$(-1)^{\text{dim} \sigma}\f$ the canonical orientation on the simplex.
    */
   Simplex_vertex_range simplex_vertex_range(Simplex_handle sh) {
+//    assert (simplex != null_simplex()); // Empty simplex
     return Simplex_vertex_range(Simplex_vertex_iterator(this, sh),
                                 Simplex_vertex_iterator(this));
   }
 
+
   /** \brief Returns a range over the simplices of the boundary of a simplex.
    *
    * The boundary of a simplex is the set of codimension \f$1\f$ subsimplices of the simplex.
@@ -386,6 +394,12 @@ class Simplex_tree {
   bool has_children(Simplex_handle sh) {
     return (sh->second.children()->parent() == sh->first);
   }
+  
+  /** \brief Returns true iff the node in the simplex tree pointed by
+   * sh has children.*/
+  bool has_children(Dit_value_t sh) {
+    return (sh.second.children()->parent() == sh.first);
+  }
 
   /** \brief Given a range of Vertex_handles, returns the Simplex_handle
    * of the simplex in the simplicial complex containing the corresponding
@@ -607,55 +621,57 @@ private:
 	 *\param vertices contains a list of vertices, which represent the simplex.
 	 *\param curr_res contains a list of vertices, which represent the current path in the tree.
 	 *\param cofaces contains a list of Simplex_handle, representing all the cofaces asked.
-	 *\param length length of the vertices list
+	 *\param length length of the vertex list "vertices"
 	 * Prefix actions : When the bottom vertex matches with the current vertex in the tree, we remove the bottom vertex from vertices.
 	 * Infix actions : Then we call or not the recursion.
 	 * Postfix actions : Finally, we add back the removed vertex into vertices, and remove this vertex from curr_res so that we didn't change the parameters.
 	 * If the vertices list is empty, we need to check if the size of the curr_res list matches with the dimension of the cofaces asked.
      */
-    void rec_coface(std::vector<Vertex_handle> &vertices, Siblings *curr_sib, std::vector<Vertex_handle> *curr_res, std::vector<Simplex_handle>& cofaces, int length, long codimension)
+    void rec_coface(std::vector<Vertex_handle> &vertices, Siblings *curr_sib, std::vector<Vertex_handle> &curr_res, std::vector<Simplex_handle>& cofaces, int length, int nbVertices)
     {
-        for (auto sib = curr_sib->members().begin(); sib != curr_sib->members().end() && (vertices.empty() || sib->first <= vertices.back()); ++sib)
+		if (!(nbVertices == length || (int)curr_res.size() <= nbVertices)) // dimension of actual simplex <= nbVertices
+			return;
+//        for (auto simplex = curr_sib->members().begin(); simplex != curr_sib->members().end(); ++simplex)
+        for (auto& simplex : curr_sib->members()) 
         {
-			bool continueRecursion = (codimension == length || (int)curr_res->size() <= codimension); // dimension of actual simplex <= codimension
             if (vertices.empty())
+			{
+				//                if ((int)curr_res.size() >= length && continueRecursion)
+				// If we reached the end of the vertices, and the simplex has more vertices than the given simplex, we found a coface
+				curr_res.push_back(simplex.first);
+				bool equalDim = (nbVertices == length || (int)curr_res.size() == nbVertices); // dimension of actual simplex == nbVertices
+				if (equalDim)
+					cofaces.push_back(curr_sib->members().find(simplex.first));
+				if (has_children(simplex))
+					rec_coface(vertices, simplex.second.children(), curr_res, cofaces, length, nbVertices);
+				curr_res.pop_back();
+			}
+            else 
             {
-                if ((int)curr_res->size() >= length && continueRecursion)
-                    // If we reached the end of the vertices, and the simplex has more vertices than the given simplex, we found a coface
+                if (simplex.first == vertices.back()) // If curr_sib matches with the top vertex
                 {
-                    curr_res->push_back(sib->first);
-					bool egalDim = (codimension == length || (int)curr_res->size() == codimension); // dimension of actual simplex == codimension
-					if (egalDim)
- 	                   cofaces.push_back(curr_sib->members().find(sib->first));
-                    if (has_children(sib))
-                        rec_coface(vertices, sib->second.children(), curr_res, cofaces, length, codimension);
-                    curr_res->pop_back();
-                }
-            }
-            else if (continueRecursion)
-            {
-                if (sib->first == vertices.back()) // If curr_sib matches with the top vertex
-                {
-                    curr_res->push_back(sib->first);
-					bool equalDim = (codimension == length || (int)curr_res->size() == codimension); // dimension of actual simplex == codimension
-                    if (vertices.size() == 1 && (int)curr_res->size() > length && equalDim)
-                        cofaces.push_back(curr_sib->members().find(sib->first));
-                    if (has_children(sib))
+                    curr_res.push_back(simplex.first);
+					bool equalDim = (nbVertices == length || (int)curr_res.size() == nbVertices); // dimension of actual simplex == nbVertices
+                    if (vertices.size() == 1 && (int)curr_res.size() > length && equalDim)
+                        cofaces.push_back(curr_sib->members().find(simplex.first));
+                    if (has_children(simplex))
                     { // Rec call
                         Vertex_handle tmp = vertices.back();
                         vertices.pop_back();
-                        rec_coface(vertices, sib->second.children(), curr_res, cofaces, length, codimension);
+                        rec_coface(vertices, simplex.second.children(), curr_res, cofaces, length, nbVertices);
                         vertices.push_back(tmp);
                     }
-                    curr_res->pop_back();
+                    curr_res.pop_back();
                 }
-                else // (sib->first < vertices.back()
+				else if (simplex.first > vertices.back())
+					return;
+                else // (simplex.first < vertices.back()
                 {
-                    if (has_children(sib))
+                    if (has_children(simplex))
                     {
-						curr_res->push_back(sib->first);
-                        rec_coface(vertices, sib->second.children(), curr_res, cofaces, length, codimension);
-						curr_res->pop_back();
+						curr_res.push_back(simplex.first);
+                        rec_coface(vertices, simplex.second.children(), curr_res, cofaces, length, nbVertices);
+						curr_res.pop_back();
                     }
                 }
             }
@@ -668,8 +684,8 @@ public:
      * \return Vector of Simplex_handle, empty vector if no cofaces found.
      */
     
-    std::vector<Simplex_handle> star_simplex_range(const Simplex_handle simplex)	{
-        return cofaces_simplex_range(simplex, 0);
+    Coface_simplex_range star_simplex_range(const Simplex_handle simplex)	{
+        return coface_simplex_range(simplex, 0);
     }
     
     
@@ -678,27 +694,24 @@ public:
      * \param simplex represent the n-simplex of which we search the n+codimension cofaces
      * \param codimension The function returns the n+codimension-cofaces of the n-simplex. If codimension = 0, return all cofaces (equivalent of star function)
      * \return Vector of Simplex_handle, empty vector if no cofaces found.
-     * \warning n+codimension must be lower than Simplex_tree dimension, otherwise an an empty vector is returned.
      */
     
-    std::vector<Simplex_handle> cofaces_simplex_range(const Simplex_handle simplex, int codimension)	{
+    Coface_simplex_range coface_simplex_range(const Simplex_handle simplex, int codimension)	{
         std::vector<Simplex_handle> cofaces;
-		if (dimension_ == -1) // Empty simplex tree
-			return cofaces;
-        if (simplex == null_simplex()) // Empty simplex
-            return cofaces;
-        if (codimension < 0) // codimension must be positive or null integer
-            return cofaces; 
+        assert (simplex != null_simplex()); // Empty simplex
+        assert (codimension >= 0); // codimension must be positive or null integer
 		Simplex_vertex_range rg = simplex_vertex_range(simplex);
         std::vector<Vertex_handle> copy(rg.begin(), rg.end());
         if (codimension + (int)copy.size() > dimension_ + 1 || (codimension == 0 && (int)copy.size() > dimension_) ) // n+codimension greater than dimension_
-            return cofaces;
-        std::sort(copy.begin(), copy.end(), std::greater<Vertex_handle>());  // must be sorted in decreasing order
+			return Coface_simplex_range(cofaces.begin(), cofaces.end());
+		assert(std::is_sorted(copy.begin(), copy.end(), std::greater<Vertex_handle>()));  // must be sorted in decreasing order
         std::vector<Vertex_handle> res;
-        rec_coface(copy, &root_, &res, cofaces, (int)copy.size(), codimension + (int)copy.size());
-        return cofaces;
+        rec_coface(copy, &root_, res, cofaces, (int)copy.size(), codimension + (int)copy.size());
+		return Coface_simplex_range(cofaces.begin(), cofaces.end());
     }
-    
+
+
+
 
  private:
   /** \brief Returns true iff the list of vertices of sh1
diff --git a/src/Simplex_tree/test/simplex_tree_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_unit_test.cpp
index f8de1c24..344cb009 100644
--- a/src/Simplex_tree/test/simplex_tree_unit_test.cpp
+++ b/src/Simplex_tree/test/simplex_tree_unit_test.cpp
@@ -174,25 +174,20 @@ void set_and_test_simplex_tree_dim_fil(typeST& simplexTree, int vectorSize, cons
   BOOST_CHECK(simplexTree.num_simplices() == nb_simplices);
 }
 
-void test_cofaces(typeST& st, std::vector<Vertex_handle> v, int dim, std::vector<std::vector<Vertex_handle>> res)
-{
-	std::vector<typeST::Simplex_handle> cofaces;
+void test_cofaces(typeST& st, std::vector<Vertex_handle> v, int dim, std::vector<typeST::Simplex_handle> res) {
+	typeST::Coface_simplex_range cofaces;
 	if (dim == 0)
 		cofaces = st.star_simplex_range(st.find(v));
 	else
-		cofaces = st.cofaces_simplex_range(st.find(v), dim);
-	std::vector<Vertex_handle> currentVertices;
-	for (unsigned long i = 0; i < cofaces.size(); ++i)
+		cofaces = st.coface_simplex_range(st.find(v), dim);
+	for (auto simplex = cofaces.begin(); simplex != cofaces.end(); ++simplex)
 	{
-		typeST::Simplex_vertex_range rg = st.simplex_vertex_range(cofaces[i]);
-		currentVertices.clear();
-		for (auto j = rg.begin(); j != rg.end(); ++j)
-		{
-			std::cout << "(" << *j << ")";
-			currentVertices.push_back(*j);
+		typeST::Simplex_vertex_range rg = st.simplex_vertex_range(*simplex);
+		for (auto vertex = rg.begin(); vertex != rg.end(); ++vertex) {
+			std::cout << "(" << *vertex << ")" ;
 		}
-		BOOST_CHECK(std::find(res.begin(), res.end(), currentVertices)!=res.end());
 		std::cout << std::endl;
+		BOOST_CHECK(std::find(res.begin(), res.end(), *simplex)!=res.end());
 	}
 }
 
@@ -620,18 +615,32 @@ BOOST_AUTO_TEST_CASE( NSimplexAndSubfaces_tree_insertion )
     st.set_dimension(3);
     std::cout << "COFACE ALGORITHM" << std::endl;
     std::vector<Vertex_handle> v;
-	std::vector<std::vector<Vertex_handle>> result;
+	std::vector<Vertex_handle> simplex;
+	std::vector<typeST::Simplex_handle> result;
     v.push_back(3);
 	std::cout << "First test : " << std::endl;
     std::cout << "Star of (3):" << std::endl;
-	int firstCoface[] = {3, 0};
-	int secondCoface[] = {4, 3};
-	int thirdCoface[] = {5, 4, 3};
-	int fourthCoface[] = {5, 3};
-	result.push_back(std::vector<Vertex_handle>(firstCoface, firstCoface + 2));
-	result.push_back(std::vector<Vertex_handle>(secondCoface, secondCoface + 2));
-	result.push_back(std::vector<Vertex_handle>(thirdCoface, thirdCoface + 3));
-	result.push_back(std::vector<Vertex_handle>(fourthCoface, fourthCoface + 2));
+	simplex.push_back(3);
+	simplex.push_back(0);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
+	simplex.push_back(4);
+	simplex.push_back(3);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
+	simplex.push_back(5);
+	simplex.push_back(4);
+	simplex.push_back(3);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
+	simplex.push_back(5);
+	simplex.push_back(3);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
     test_cofaces(st, v, 0, result);
     v.clear();
 	result.clear();
@@ -640,34 +649,57 @@ BOOST_AUTO_TEST_CASE( NSimplexAndSubfaces_tree_insertion )
     v.push_back(7);
 	std::cout << "Second test : " << std::endl;
     std::cout << "Star of (1,7): " << std::endl;
-	int firstCoface2[] = {7, 6, 1, 0};
-	int secondCoface2[] = {7, 1, 0};
-	int thirdCoface2[] = {7 ,6, 1};
-	result.push_back(std::vector<Vertex_handle>(firstCoface2, firstCoface2 + 4));
-	result.push_back(std::vector<Vertex_handle>(secondCoface2, secondCoface2 + 3));
-	result.push_back(std::vector<Vertex_handle>(thirdCoface2, thirdCoface2 + 3));
+
+	simplex.push_back(7);
+	simplex.push_back(6);
+	simplex.push_back(1);
+	simplex.push_back(0);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
+	simplex.push_back(7);
+	simplex.push_back(1);
+	simplex.push_back(0);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
+	simplex.push_back(7);
+	simplex.push_back(6);
+	simplex.push_back(1);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
     test_cofaces(st, v, 0, result);
 	result.clear();
 
 	std::cout << "Third test : " << std::endl;
     std::cout << "2-dimension Cofaces of simplex(1,7) : " << std::endl;
-	int secondCoface3[] = {7, 1, 0};
-	int firstCoface3[] = {7, 6, 1};
-	result.push_back(std::vector<Vertex_handle>(firstCoface3, firstCoface3 + 3));
-	result.push_back(std::vector<Vertex_handle>(secondCoface3, secondCoface3 + 3));
+
+	simplex.push_back(7);
+	simplex.push_back(1);
+	simplex.push_back(0);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
+	simplex.push_back(7);
+	simplex.push_back(6);
+	simplex.push_back(1);
+	result.push_back(st.find(simplex));
+	simplex.clear();
+
     test_cofaces(st, v, 1, result);
 	result.clear();
 
-    std::cout << "Cofaces with a codimension too high (codimension + vetices > tree.dimension)" << std::endl;
+    std::cout << "Cofaces with a codimension too high (codimension + vetices > tree.dimension) :" << std::endl;
     test_cofaces(st, v, 5, result);
-    std::cout << "Cofaces with an empty codimension" << std::endl;
-    test_cofaces(st, v, -1, result);
-    std::cout << "Cofaces in an empty simplex tree" << std::endl;
-    typeST empty_tree;
-    test_cofaces(empty_tree, v, 1, result);
-    std::cout << "Cofaces of an empty simplex" << std::endl;
-    v.clear();
-    test_cofaces(st, v, 1, result);
+//    std::cout << "Cofaces with an empty codimension" << std::endl;
+//    test_cofaces(st, v, -1, result);
+//    std::cout << "Cofaces in an empty simplex tree" << std::endl;
+//   typeST empty_tree;
+//    test_cofaces(empty_tree, v, 1, result);
+//    std::cout << "Cofaces of an empty simplex" << std::endl;
+//    v.clear();
+//    test_cofaces(st, v, 1, result);
  
  /*
     // TEST Off read