Answers to Marc's comments.

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

Former-commit-id: 481841d609bae68df85a93527e6d3cc20f1b06e2

1 files changed, 196 insertions, 159 deletions

diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h
index 22b703a9..807be335 100644
--- a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h
+++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h
@@ -23,13 +23,14 @@
 #pragma once
 
 #include <iostream>
+#include <numeric>
 #include <string>
 #include <vector>
 #include <string>
 #include <fstream>
 #include <algorithm>
 #include <iterator>
-#include <limits>

+#include <limits>
 #include <ctime>
 #include "Bitmap_cubical_complex/counter.h"
 
@@ -43,7 +44,11 @@ namespace Cubical_complex
 {
 
 
-
+/**
+ *@class Bitmap_cubical_complex_base
+ *@brief Cubical complex represented as a bitmap, class with basic implementation.
+ *@ingroup cubical_complex
+ */
 /**
  * This is a class implementing a basic bitmap data structure to store cubical complexes.
  * It implements only the most basic subroutines.
@@ -64,13 +69,13 @@ namespace Cubical_complex
 template <typename T>
 class Bitmap_cubical_complex_base
 {
-public:

-    typedef T filtration_type;

-    /**

-    *Default constructor

-    **/

-    Bitmap_cubical_complex_base()

-    {

+public:
+    typedef T filtration_type;
+    /**
+    *Default constructor
+    **/
+    Bitmap_cubical_complex_base()
+    {
     }
     /**
     * There are a few constructors of a Bitmap_cubical_complex_base class.
@@ -93,6 +98,11 @@ public:
     Bitmap_cubical_complex_base( const std::vector<unsigned>& dimensions , const std::vector<T>& top_dimensional_cells );
 
     /**
+    * Destructor of the Bitmap_cubical_complex_base class.
+    **/
+    virtual ~Bitmap_cubical_complex_base(){}
+
+    /**
     * The functions get_boundary_of_a_cell, get_coboundary_of_a_cell, get_dimension_of_a_cell
     * and get_cell_data are the basic
      * functions that compute boundary / coboundary / dimension and the filtration
@@ -152,18 +162,27 @@ public:
     * Writing to stream operator.
     **/
     template <typename K>
-    friend ostream& operator << ( ostream & os , const Bitmap_cubical_complex_base<K>& b );

-

-

+    friend ostream& operator << ( ostream & os , const Bitmap_cubical_complex_base<K>& b );
+
+
     /**
-    * Functions that put the input data to bins.
+    * Function that put the input data to bins. Sometimes if most of the cells have different birth-death times, the performance of the algorithms to compute persistence gets
+    * worst. When dealing with this type of data, one may want to put different values on cells to some number of bins. The function put_data_toBins( size_t number_of_bins )
+    * ais designed for that purpose. The parameter of the function is the number of bins (distinct values) we want to have in the cubical complex.
     **/
     void put_data_toBins( size_t number_of_bins );
-    void put_data_toBins( T diameter_of_bin );

-

+
+    /**
+    * Function that put the input data to bins. Sometimes if most of the cells have different birth-death times, the performance of the algorithms to compute persistence gets
+    * worst. When dealing with this type of data, one may want to put different values on cells to some number of bins. The function put_data_toBins( T diameter_of_bin ) is
+    * designed for that purpose. The parameter of it is the diameter of each bin. Note that the bottleneck distance between the persistence diagram of the cubical complex
+    * before and after using such a function will be bounded by the parameter diameter_of_bin.
+    **/
+    void put_data_toBins( T diameter_of_bin );
+
     /**
     * Functions to find min and max values of filtration.
-    **/

+    **/
     std::pair< T ,T > min_max_filtration();
 
     //ITERATORS
@@ -172,42 +191,60 @@ public:
     * Iterator through all cells in the complex (in order they appear in the structure -- i.e.
     * in lexicographical order).
     **/
-    typedef typename std::vector< T >::iterator all_cells_iterator;

-

+    typedef typename std::vector< T >::iterator all_cells_iterator;
+
+
     /**
-    * Function returning an iterator to the first cell of the bitmap.
+    * Constant iterator through all cells in the complex (in order they appear in the structure -- i.e.
+    * in lexicographical order).
     **/
-    all_cells_iterator all_cells_begin()
+    typedef typename std::vector< T >::const_iterator all_cells_const_iterator;
+
+    /**
+    * Function returning a constant iterator to the first cell of the bitmap.
+    **/
+    all_cells_const_iterator all_cells_const_begin()const
     {
         return this->data.begin();
-    }

-

-     /**
-    * Function returning an iterator to the last cell of the bitmap.
+    }
+
+
+    /**
+    * Function returning a constant iterator to the last cell of the bitmap.
     **/
-    all_cells_iterator all_cells_end()const
+    all_cells_const_iterator all_cells_const_end()const
     {
         return this->data.end();
     }
 
     /**
-    * Constant iterator through all cells in the complex (in order they appear in the structure -- i.e.
-    * in lexicographical order).
+    * Function returning an iterator to the first cell of the bitmap.
     **/
-    typedef typename std::vector< T >::const_iterator all_cells_const_iterator;

-

+    all_cells_iterator all_cells_begin()
+    {
+        return this->data.begin();
+    }
+
     /**
     * Function returning a constant iterator to the first cell of the bitmap.
     **/
-    all_cells_const_iterator all_cells_const_begin()const
+    all_cells_const_iterator all_cells_begin() const
     {
         return this->data.begin();
-    }

-

+    }
+
+    /**
+    * Function returning an iterator to the last cell of the bitmap.
+    **/
+    all_cells_iterator all_cells_end()
+    {
+        return this->data.end();
+    }
+
     /**
     * Function returning a constant iterator to the last cell of the bitmap.
     **/
-    all_cells_const_iterator all_cells_const_end()const
+    all_cells_const_iterator all_cells_end() const
     {
         return this->data.end();
     }
@@ -303,8 +340,8 @@ public:
         protected:
             std::vector< size_t > counter;
             Bitmap_cubical_complex_base& b;
-    };

-

+    };
+
     /**
     * Function returning a Top_dimensional_cells_iterator to the first top dimensional cell cell of the bitmap.
     **/
@@ -312,8 +349,8 @@ public:
     {
         Top_dimensional_cells_iterator a(*this);
         return a;
-    }

-

+    }
+
     /**
     * Function returning a Top_dimensional_cells_iterator to the last top dimensional cell cell of the bitmap.
     **/
@@ -333,11 +370,11 @@ public:
 //****************************************************************************************************************//
 //****************************************************************************************************************//
 //****************************************************************************************************************//
-

-

-inline size_t number_cells()const

-{

-    return this->total_number_of_cells;

+
+
+inline size_t number_cells()const
+{
+    return this->total_number_of_cells;
 }
 
 //****************************************************************************************************************//
@@ -375,7 +412,7 @@ protected:
 
      std::vector<unsigned> compute_counter_for_given_cell( size_t cell )const
     {
-        std::vector<unsigned> counter;

+        std::vector<unsigned> counter;
         counter.reserve( this->sizes.size() );
         for ( size_t dim = this->sizes.size() ; dim != 0 ; --dim )
         {
@@ -384,59 +421,59 @@ protected:
         }
         std::reverse( counter.begin() , counter.end() );
         return counter;
-    }

-    void read_perseus_style_file( const char* perseus_style_file );

-    void setup_bitmap_based_on_top_dimensional_cells_list(const std::vector<unsigned>& sizes_in_following_directions , const std::vector<T>& top_dimensional_cells);

-    Bitmap_cubical_complex_base( const char* perseus_style_file , std::vector<bool> directions );

-    Bitmap_cubical_complex_base( const std::vector<unsigned>& sizes , std::vector<bool> directions );

+    }
+    void read_perseus_style_file( const char* perseus_style_file );
+    void setup_bitmap_based_on_top_dimensional_cells_list(const std::vector<unsigned>& sizes_in_following_directions , const std::vector<T>& top_dimensional_cells);
+    Bitmap_cubical_complex_base( const char* perseus_style_file , std::vector<bool> directions );
+    Bitmap_cubical_complex_base( const std::vector<unsigned>& sizes , std::vector<bool> directions );
     Bitmap_cubical_complex_base( const std::vector<unsigned>& dimensions , const std::vector<T>& top_dimensional_cells , std::vector<bool> directions );
 };
 
 
-template <typename T>

-void Bitmap_cubical_complex_base<T>::put_data_toBins( size_t number_of_bins )

-{

-    bool bdg = false;

-

-    std::pair< T ,T > min_max = this->min_max_filtration();

-    T dx = (min_max.second-min_max.first)/(T)number_of_bins;

-

-    //now put the data into the appropriate bins:

-    for ( size_t i = 0 ; i != this->data.size() ; ++i )

-    {

-        if ( bdg ){cerr << "Before binning : " << this->data[i] << endl;}

-        this->data[i] = min_max.first + dx*(this->data[i]-min_max.first)/number_of_bins;

-        if ( bdg ){cerr << "After binning : " << this->data[i] << endl;getchar();}

-    }

-}

-

 template <typename T>
-void Bitmap_cubical_complex_base<T>::put_data_toBins( T diameter_of_bin )

-{

-    bool bdg = false;

-    std::pair< T ,T > min_max = this->min_max_filtration();

-

-    size_t number_of_bins = (min_max.second - min_max.first)/diameter_of_bin;

-    //now put the data into the appropriate bins:

-    for ( size_t i = 0 ; i != this->data.size() ; ++i )

-    {

-        if ( bdg ){cerr << "Before binning : " << this->data[i] << endl;}

-        this->data[i] = min_max.first + diameter_of_bin*(this->data[i]-min_max.first)/number_of_bins;

-        if ( bdg ){cerr << "After binning : " << this->data[i] << endl;getchar();}

-    }

-}

-

-template <typename T>

-std::pair< T ,T > Bitmap_cubical_complex_base<T>::min_max_filtration()

-{

-    std::pair< T ,T > min_max( std::numeric_limits<T>::max() , std::numeric_limits<T>::min() );

-    for ( size_t i = 0 ; i != this->data.size() ; ++i )

-    {

-        if ( this->data[i] < min_max.first )min_max.first = this->data[i];

-        if ( this->data[i] > min_max.second )min_max.second = this->data[i];

-    }

-    return min_max;

-}

+void Bitmap_cubical_complex_base<T>::put_data_toBins( size_t number_of_bins )
+{
+    bool bdg = false;
+
+    std::pair< T ,T > min_max = this->min_max_filtration();
+    T dx = (min_max.second-min_max.first)/(T)number_of_bins;
+
+    //now put the data into the appropriate bins:
+    for ( size_t i = 0 ; i != this->data.size() ; ++i )
+    {
+        if ( bdg ){cerr << "Before binning : " << this->data[i] << endl;}
+        this->data[i] = min_max.first + dx*(this->data[i]-min_max.first)/number_of_bins;
+        if ( bdg ){cerr << "After binning : " << this->data[i] << endl;getchar();}
+    }
+}
+
+template <typename T>
+void Bitmap_cubical_complex_base<T>::put_data_toBins( T diameter_of_bin )
+{
+    bool bdg = false;
+    std::pair< T ,T > min_max = this->min_max_filtration();
+
+    size_t number_of_bins = (min_max.second - min_max.first)/diameter_of_bin;
+    //now put the data into the appropriate bins:
+    for ( size_t i = 0 ; i != this->data.size() ; ++i )
+    {
+        if ( bdg ){cerr << "Before binning : " << this->data[i] << endl;}
+        this->data[i] = min_max.first + diameter_of_bin*(this->data[i]-min_max.first)/number_of_bins;
+        if ( bdg ){cerr << "After binning : " << this->data[i] << endl;getchar();}
+    }
+}
+
+template <typename T>
+std::pair< T ,T > Bitmap_cubical_complex_base<T>::min_max_filtration()
+{
+    std::pair< T ,T > min_max( std::numeric_limits<T>::max() , std::numeric_limits<T>::min() );
+    for ( size_t i = 0 ; i != this->data.size() ; ++i )
+    {
+        if ( this->data[i] < min_max.first )min_max.first = this->data[i];
+        if ( this->data[i] > min_max.second )min_max.second = this->data[i];
+    }
+    return min_max;
+}
 
 
 template <typename K>
@@ -457,10 +494,10 @@ Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base
 {
     this->set_up_containers( sizes );
 }
-

-template <typename T>

-void Bitmap_cubical_complex_base<T>::setup_bitmap_based_on_top_dimensional_cells_list(const std::vector<unsigned>& sizes_in_following_directions , const std::vector<T>& top_dimensional_cells)

-{

+
+template <typename T>
+void Bitmap_cubical_complex_base<T>::setup_bitmap_based_on_top_dimensional_cells_list(const std::vector<unsigned>& sizes_in_following_directions , const std::vector<T>& top_dimensional_cells)
+{
     this->set_up_containers( sizes_in_following_directions );
 
     size_t number_of_top_dimensional_elements = 1;
@@ -489,19 +526,19 @@ void Bitmap_cubical_complex_base<T>::setup_bitmap_based_on_top_dimensional_cells
         (*it) = top_dimensional_cells[index];
         ++index;
     }
-    this->impose_lower_star_filtration();

-}

+    this->impose_lower_star_filtration();
+}
 
 template <typename T>
 Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base
 ( const std::vector<unsigned>& sizes_in_following_directions , const std::vector<T>& top_dimensional_cells )
 {
     this->setup_bitmap_based_on_top_dimensional_cells_list( sizes_in_following_directions , top_dimensional_cells );
-}

-

-template <typename T>

-void Bitmap_cubical_complex_base<T>::read_perseus_style_file( const char* perseus_style_file )

-{

+}
+
+template <typename T>
+void Bitmap_cubical_complex_base<T>::read_perseus_style_file( const char* perseus_style_file )
+{
     bool dbg = false;
     ifstream inFiltration, inIds;
     inFiltration.open( perseus_style_file );
@@ -510,7 +547,7 @@ void Bitmap_cubical_complex_base<T>::read_perseus_style_file( const char* perseu
 
     if (dbg){cerr << "dimensionOfData : " << dimensionOfData << endl;getchar();}
 
-    std::vector<unsigned> sizes;

+    std::vector<unsigned> sizes;
     sizes.reserve( dimensionOfData );
     for ( size_t i = 0 ; i != dimensionOfData ; ++i )
     {
@@ -541,31 +578,31 @@ void Bitmap_cubical_complex_base<T>::read_perseus_style_file( const char* perseu
         ++it;
     }
     inFiltration.close();
-    this->impose_lower_star_filtration();

-}

-

-template <typename T>

-Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base( const char* perseus_style_file , std::vector<bool> directions )

-{

-    //this constructor is here just for compatibility with a class that creates cubical complexes with periodic bundary conditions.

-    //It ignores the last parameter of the function.

-    this->read_perseus_style_file( perseus_style_file );

-}

-

-template <typename T>

-Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base( const std::vector<unsigned>& sizes , std::vector<bool> directions )

-{

-    //this constructor is here just for compatibility with a class that creates cubical complexes with periodic bundary conditions.

-    //It ignores the last parameter of the function.

-    this->set_up_containers( sizes );

-}

-

-template <typename T>

+    this->impose_lower_star_filtration();
+}
+
+template <typename T>
+Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base( const char* perseus_style_file , std::vector<bool> directions )
+{
+    //this constructor is here just for compatibility with a class that creates cubical complexes with periodic bundary conditions.
+    //It ignores the last parameter of the function.
+    this->read_perseus_style_file( perseus_style_file );
+}
+
+template <typename T>
+Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base( const std::vector<unsigned>& sizes , std::vector<bool> directions )
+{
+    //this constructor is here just for compatibility with a class that creates cubical complexes with periodic bundary conditions.
+    //It ignores the last parameter of the function.
+    this->set_up_containers( sizes );
+}
+
+template <typename T>
 Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base( const std::vector<unsigned>& dimensions , const std::vector<T>& top_dimensional_cells , std::vector<bool> directions )
-{

-    //this constructor is here just for compatibility with a class that creates cubical complexes with periodic bundary conditions.

-    //It ignores the last parameter of the function.

-    this->setup_bitmap_based_on_top_dimensional_cells_list( dimensions , top_dimensional_cells );

+{
+    //this constructor is here just for compatibility with a class that creates cubical complexes with periodic bundary conditions.
+    //It ignores the last parameter of the function.
+    this->setup_bitmap_based_on_top_dimensional_cells_list( dimensions , top_dimensional_cells );
 }
 
 template <typename T>
@@ -578,29 +615,29 @@ Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base( const char* perseus
 template <typename T>
 std::vector< size_t > Bitmap_cubical_complex_base<T>::get_boundary_of_a_cell( size_t cell )const
 {
-    std::vector< size_t > boundary_elements;

-

-    //Speed traded of for memory. Check if it is better in practice.

-    boundary_elements.reserve( this->dimension()*2 );

+    std::vector< size_t > boundary_elements;
+
+    //Speed traded of for memory. Check if it is better in practice.
+    boundary_elements.reserve( this->dimension()*2 );
 
     size_t cell1 = cell;
     for ( size_t i = this->multipliers.size() ; i != 0 ; --i )
     {
         unsigned position = cell1/this->multipliers[i-1];
         if ( position%2 == 1 )
-        {

-            boundary_elements.push_back( cell - this->multipliers[ i-1 ] );

+        {
+            boundary_elements.push_back( cell - this->multipliers[ i-1 ] );
             boundary_elements.push_back( cell + this->multipliers[ i-1 ] );
         }
         cell1 = cell1%this->multipliers[i-1];
     }
     return boundary_elements;
-}

+}
+
+
+
 
 
-

-

-

 template <typename T>
 std::vector< size_t > Bitmap_cubical_complex_base<T>::get_coboundary_of_a_cell( size_t cell )const
 {
@@ -625,8 +662,8 @@ std::vector< size_t > Bitmap_cubical_complex_base<T>::get_coboundary_of_a_cell(
         cell1 = cell1%this->multipliers[i-1];
     }
     return coboundary_elements;
-}

-

+}
+
 
 
 
@@ -675,14 +712,14 @@ void Bitmap_cubical_complex_base<T>::impose_lower_star_filtration()
     //this vector will be used to check which elements have already been taken care of
     //in imposing lower star filtration:
     std::vector<bool> is_this_cell_considered( this->data.size() , false );
-

-    size_t size_to_reserve = 1;

-    for ( size_t i = 0 ; i != this->multipliers.size() ; ++i )

-    {

-        size_to_reserve *= (size_t)((this->multipliers[i]-1)/2);

-    }

-
-    std::vector<size_t> indices_to_consider;

+
+    size_t size_to_reserve = 1;
+    for ( size_t i = 0 ; i != this->multipliers.size() ; ++i )
+    {
+        size_to_reserve *= (size_t)((this->multipliers[i]-1)/2);
+    }
+
+    std::vector<size_t> indices_to_consider;
     indices_to_consider.reserve( size_to_reserve );
     //we assume here that we already have a filtration on the top dimensional cells and
     //we have to extend it to lower ones.
@@ -708,20 +745,20 @@ void Bitmap_cubical_complex_base<T>::impose_lower_star_filtration()
         {
             std::vector<size_t> bd = this->get_boundary_of_a_cell( indices_to_consider[i] );
             for ( size_t boundaryIt = 0 ; boundaryIt != bd.size() ; ++boundaryIt )
-            {

-                if ( dbg )

-                {

-                    cerr << "filtration of a cell : " << bd[boundaryIt] << " is : " << this->data[ bd[boundaryIt] ] << " while of a cell: " << indices_to_consider[i] << " is: " << this->data[ indices_to_consider[i] ] << endl;

-                    getchar();

-

+            {
+                if ( dbg )
+                {
+                    cerr << "filtration of a cell : " << bd[boundaryIt] << " is : " << this->data[ bd[boundaryIt] ] << " while of a cell: " << indices_to_consider[i] << " is: " << this->data[ indices_to_consider[i] ] << endl;
+                    getchar();
+
                 }
                 if ( this->data[ bd[boundaryIt] ] > this->data[ indices_to_consider[i] ] )
                 {
-                    this->data[ bd[boundaryIt] ] = this->data[ indices_to_consider[i] ];

-                    if ( dbg )

-                    {

-                        cerr << "Setting the value of a cell : " << bd[boundaryIt] << " to : " << this->data[ indices_to_consider[i] ] << endl;

-                        getchar();

+                    this->data[ bd[boundaryIt] ] = this->data[ indices_to_consider[i] ];
+                    if ( dbg )
+                    {
+                        cerr << "Setting the value of a cell : " << bd[boundaryIt] << " to : " << this->data[ indices_to_consider[i] ] << endl;
+                        getchar();
                     }
                 }
                 if ( is_this_cell_considered[ bd[boundaryIt] ] == false )
@@ -760,4 +797,4 @@ bool compareFirstElementsOfTuples( const std::pair< std::pair< T , size_t > , ch
 
 }
 
-}

+}