From 62937147e40a7d2da7aa7a7a604808feeccaa75e Mon Sep 17 00:00:00 2001 From: vrouvrea Date: Fri, 25 Sep 2015 14:57:29 +0000 Subject: Add bitmap cubical complex git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/bitmap@794 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 50d9b8eb80e0fe99f871afa5bdbb853add97e25e --- .../include/gudhi/Bitmap_cubical_complex.h | 706 +++++++++++++++++++++ .../include/gudhi/Bitmap_cubical_complex_base.h | 577 +++++++++++++++++ src/Bitmap_cubical_complex/include/gudhi/counter.h | 136 ++++ 3 files changed, 1419 insertions(+) create mode 100644 src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h create mode 100644 src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h create mode 100644 src/Bitmap_cubical_complex/include/gudhi/counter.h (limited to 'src/Bitmap_cubical_complex/include/gudhi') diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h new file mode 100644 index 00000000..61ae8105 --- /dev/null +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h @@ -0,0 +1,706 @@ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): Pawel Dlotko + * + * Copyright (C) 2015 INRIA Sophia-Saclay (France) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +#ifndef BITMAP_CUBICAL_COMPLEX_H_ +#define BITMAP_CUBICAL_COMPLEX_H_ + +#include + +#include + +//global variable, was used just for debugging. +bool globalDbg = false; + +template +class Bitmap_cubical_complex : public Bitmap_cubical_complex_base { + public: + //*********************************************************************************************************************************// + //Typedefs and typenames + //*********************************************************************************************************************************// + friend class Simplex_handle; + typedef size_t Simplex_key; + typedef T Filtration_value; + + + //*********************************************************************************************************************************// + //Simplex handle class + //*********************************************************************************************************************************// + + /** + * Handle of a cell, required for compatibility with the function to compute persistence in Gudhi. Elements of this class are: the pointer to the bitmap B in which the considered cell is + * together with a position of this cell in B. Given this data, one can get all the information about the considered cell. + **/ + class Simplex_handle { + public: + + Simplex_handle() { + if (globalDbg) { + cerr << "Simplex_handle()\n"; + } + this->b = 0; + this->position = 0; + } + + Simplex_handle(Bitmap_cubical_complex* b) { + if (globalDbg) { + cerr << "Simplex_handle(Bitmap_cubical_complex* b)\n"; + } + this->b = b; + this->position = 0; + } + + Simplex_handle(const Simplex_handle& org) : b(org.b) { + if (globalDbg) { + cerr << "Simplex_handle( const Simplex_handle& org )\n"; + } + this->position = org.position; + } + + Simplex_handle& operator=(const Simplex_handle& rhs) { + if (globalDbg) { + cerr << "Simplex_handle operator = \n"; + } + this->position = rhs.position; + this->b = rhs.b; + return *this; + } + + Simplex_handle(Bitmap_cubical_complex* b, Simplex_key position) { + if (globalDbg) { + cerr << "Simplex_handle(Bitmap_cubical_complex* b , Simplex_key position)\n"; + cerr << "Position : " << position << endl; + } + this->b = b; + this->position = position; + } + friend class Bitmap_cubical_complex; + private: + Bitmap_cubical_complex* b; + Simplex_key position; //Assumption -- this field always keep the REAL position of simplex in the bitmap, no matter what keys have been. + //to deal with the keys, the class Bitmap_cubical_complex have extra vectors: keyAssociatedToSimplex and simplexAssociatedToKey + //that allow to move between actual cell and the key assigned to it. + }; + + + //*********************************************************************************************************************************// + //Constructors + //*********************************************************************************************************************************// + //Over here we need to definie various input types. I am proposing the following ones: + //Perseus style + //H5 files? TODO + //binary files with little endiangs / big endians? TODO + //constructor from a vector of elements of a type T. TODO + + /** + * Constructor form a Perseus-style file. + **/ + Bitmap_cubical_complex(char* perseusStyleFile) : Bitmap_cubical_complex_base(perseusStyleFile) { + if (globalDbg) { + cerr << "Bitmap_cubical_complex( char* perseusStyleFile )\n"; + } + std::vector< size_t > keyAssociatedToSimplex(this->totalNumberOfCells + 1); + std::vector< size_t > simplexAssociatedToKey(this->totalNumberOfCells + 1); + + for (size_t i = 0; i != this->totalNumberOfCells; ++i) { + keyAssociatedToSimplex[i] = simplexAssociatedToKey[i] = i; + } + this->keyAssociatedToSimplex = keyAssociatedToSimplex; + this->simplexAssociatedToKey = simplexAssociatedToKey; + //we initialize this only once, in each constructor, when the bitmap is constructed. If the user decide to change some elements of the bitmap, then this procedure need + //to be called again. + this->initializeElementsOrderedAccordingToFiltration(); + } + + /** + * Constructor that requires vector of elements of type unsigned, which gives number of top dimensional cells in the following directions and vector of element of a type T + * with filtration on top dimensional cells. + **/ + Bitmap_cubical_complex(std::vector dimensions, std::vector topDimensionalCells) : Bitmap_cubical_complex_base(dimensions, topDimensionalCells) { + std::vector< size_t > keyAssociatedToSimplex(this->totalNumberOfCells + 1); + std::vector< size_t > simplexAssociatedToKey(this->totalNumberOfCells + 1); + + for (size_t i = 0; i != this->totalNumberOfCells; ++i) { + keyAssociatedToSimplex[i] = simplexAssociatedToKey[i] = i; + } + this->keyAssociatedToSimplex = keyAssociatedToSimplex; + this->simplexAssociatedToKey = simplexAssociatedToKey; + //we initialize this only once, in each constructor, when the bitmap is constructed. If the user decide to change some elements of the bitmap, then this procedure need + //to be called again. + this->initializeElementsOrderedAccordingToFiltration(); + } + + //*********************************************************************************************************************************// + //Other 'easy' functions + //*********************************************************************************************************************************// + + /** + * Returns number of all cubes in the complex. + **/ + size_t num_simplices()const { + return this->totalNumberOfCells; + } + + /** + * Returns a Simplex_handle to a cube that do not exist in this complex. + **/ + Simplex_handle null_simplex() { + return Simplex_handle(this, this->data.size()); + } + + /** + * Returns dimension of the complex. + **/ + size_t dimension() { + return this->sizes.size(); + } + + /** + * Return dimension of a cell pointed by the Simplex_handle. + **/ + size_t dimension(const Simplex_handle& sh) { + if (globalDbg) { + cerr << "int dimension(const Simplex_handle& sh)\n"; + } + if (sh.position != this->data.size()) return sh.b->get_dimension_of_a_cell(sh.position); + return std::numeric_limits::max(); + } + + /** + * Return the filtration of a cell pointed by the Simplex_handle. + **/ + T filtration(const Simplex_handle& sh) { + if (globalDbg) { + cerr << "T filtration(const Simplex_handle& sh)\n"; + } + //Returns the filtration value of a simplex. + if (sh.position != this->data.size()) return sh.b->data[ sh.position ]; + return INT_MAX; + } + + /** + * Return a key which is not a key of any cube in the considered data structure. + **/ + Simplex_key null_key() { + if (globalDbg) { + cerr << "Simplex_key null_key()\n"; + } + return this->data.size(); + } + + /** + * Return the key of a cube pointed by the Simplex_handle. + **/ + Simplex_key key(const Simplex_handle& sh) { + if (globalDbg) { + cerr << "Simplex_key key(const Simplex_handle& sh)\n"; + } + return sh.b->keyAssociatedToSimplex[ sh.position ]; + } + + /** + * Return the Simplex_handle given the key of the cube. + **/ + Simplex_handle simplex(Simplex_key key) { + if (globalDbg) { + cerr << "Simplex_handle simplex(Simplex_key key)\n"; + } + return Simplex_handle(this, this->simplexAssociatedToKey[ key ]); + } + + /** + * Assign key to a cube pointed by the Simplex_handle + **/ + void assign_key(Simplex_handle& sh, Simplex_key key) { + if (globalDbg) { + cerr << "void assign_key(Simplex_handle& sh, Simplex_key key)\n"; + } + this->keyAssociatedToSimplex[sh.position] = key; + this->simplexAssociatedToKey[key] = sh.position; + } + + /** + * Function called from a constructor. It is needed for Filtration_simplex_iterator to work. + **/ + void initializeElementsOrderedAccordingToFiltration(); + + + + //*********************************************************************************************************************************// + //Iterators + //*********************************************************************************************************************************// + + /** + * Boundary_simplex_iterator class allows iteration on boundary of each cube. + **/ + class Boundary_simplex_range; + + class Boundary_simplex_iterator : std::iterator< std::input_iterator_tag, Simplex_handle > { + //Iterator on the simplices belonging to the boundary of a simplex. + //value_type must be 'Simplex_handle'. + public: + + Boundary_simplex_iterator(Simplex_handle& sh) : sh(sh) { + if (globalDbg) { + cerr << "Boundary_simplex_iterator( Simplex_handle& sh )\n"; + } + this->position = 0; + this->boundaryElements = this->sh.b->get_boundary_of_a_cell(this->sh.position); + } + + Boundary_simplex_iterator operator++() { + if (globalDbg) { + cerr << "Boundary_simplex_iterator operator++()\n"; + } + ++this->position; + return *this; + } + + Boundary_simplex_iterator operator++(int) { + Boundary_simplex_iterator result = *this; + ++(*this); + return result; + } + + Boundary_simplex_iterator operator=(const Boundary_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "Boundary_simplex_iterator operator =\n"; + } + this->sh = rhs.sh; + this->boundaryElements.clear(); + this->boundaryElementsinsert(this->boundaryElements.end(), rhs.boundaryElements.begin(), rhs.boundaryElements.end()); + } + + bool operator==(const Boundary_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "bool operator ==\n"; + } + if (this->position == rhs.position) { + if (this->boundaryElements.size() != rhs.boundaryElements.size())return false; + for (size_t i = 0; i != this->boundaryElements.size(); ++i) { + if (this->boundaryElements[i] != rhs.boundaryElements[i])return false; + } + return true; + } + return false; + } + + bool operator!=(const Boundary_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "bool operator != \n"; + } + return !(*this == rhs); + } + + Simplex_handle operator*() { + if (globalDbg) { + cerr << "Simplex_handle operator*\n"; + } + return Simplex_handle(this->sh.b, this->boundaryElements[this->position]); + } + + friend class Boundary_simplex_range; + private: + Simplex_handle sh; + std::vector< size_t > boundaryElements; + size_t position; + }; + + /** + * Boundary_simplex_range class provides ranges for boundary iterators. + **/ + class Boundary_simplex_range { + //Range giving access to the simplices in the boundary of a simplex. + //.begin() and .end() return type Boundary_simplex_iterator. + public: + + Boundary_simplex_range(const Simplex_handle& sh) : sh(sh) { }; + + Boundary_simplex_iterator begin() { + if (globalDbg) { + cerr << "Boundary_simplex_iterator begin\n"; + } + Boundary_simplex_iterator it(this->sh); + return it; + } + + Boundary_simplex_iterator end() { + if (globalDbg) { + cerr << "Boundary_simplex_iterator end()\n"; + } + Boundary_simplex_iterator it(this->sh); + it.position = it.boundaryElements.size(); + return it; + } + private: + Simplex_handle sh; + }; + + + /** + * Filtration_simplex_iterator class provides an iterator though the whole structure in the order of filtration. Secondary criteria for filtration are: + * (1) Dimension of a cube (lower dimensional comes first). + * (2) Position in the data structure (the ones that are earlies in the data structure comes first). + **/ + class Filtration_simplex_range; + + class Filtration_simplex_iterator : std::iterator< std::input_iterator_tag, Simplex_handle > { + //Iterator over all simplices of the complex in the order of the indexing scheme. + //'value_type' must be 'Simplex_handle'. + public: + + Filtration_simplex_iterator(Bitmap_cubical_complex* b) : b(b), position(0) { }; + + Filtration_simplex_iterator() : b(NULL) { }; + + Filtration_simplex_iterator operator++() { + if (globalDbg) { + cerr << "Filtration_simplex_iterator operator++\n"; + } + ++this->position; + return (*this); + } + + Filtration_simplex_iterator operator++(int) { + Filtration_simplex_iterator result = *this; + ++(*this); + return result; + } + + Filtration_simplex_iterator operator=(const Filtration_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "Filtration_simplex_iterator operator =\n"; + } + this->b = rhs.b; + this->position = rhs.position; + } + + bool operator==(const Filtration_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "bool operator == ( const Filtration_simplex_iterator& rhs )\n"; + } + if (this->position == rhs.position) { + return true; + } + return false; + } + + bool operator!=(const Filtration_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "bool operator != ( const Filtration_simplex_iterator& rhs )\n"; + } + return !(*this == rhs); + } + + Simplex_handle operator*() { + if (globalDbg) { + cerr << "Simplex_handle operator*()\n"; + } + return Simplex_handle(this->b, this->b->elementsOrderedAccordingToFiltration[ this->position ]); + } + + friend class Filtration_simplex_range; + private: + Bitmap_cubical_complex* b; + size_t position; + }; + + /** + * Filtration_simplex_range provides the ranges for Filtration_simplex_iterator. + **/ + class Filtration_simplex_range { + //Range over the simplices of the complex in the order of the filtration. + //.begin() and .end() return type Filtration_simplex_iterator. + public: + + Filtration_simplex_range(Bitmap_cubical_complex* b) : b(b) { }; + + Filtration_simplex_iterator begin() { + if (globalDbg) { + cerr << "Filtration_simplex_iterator begin() \n"; + } + return Filtration_simplex_iterator(this->b); + } + + Filtration_simplex_iterator end() { + if (globalDbg) { + cerr << "Filtration_simplex_iterator end()\n"; + } + Filtration_simplex_iterator it(this->b); + it.position = this->b->elementsOrderedAccordingToFiltration.size(); + return it; + } + private: + Bitmap_cubical_complex* b; + }; + + + + //*********************************************************************************************************************************// + //Methods to access iterators from the container: + + /** + * boundary_simplex_range creates an object of a Boundary_simplex_range class that provides ranges for the Boundary_simplex_iterator. + **/ + Boundary_simplex_range boundary_simplex_range(Simplex_handle& sh) { + if (globalDbg) { + cerr << "Boundary_simplex_range boundary_simplex_range(Simplex_handle& sh)\n"; + } + //Returns a range giving access to all simplices of the boundary of a simplex, i.e. the set of codimension 1 subsimplices of the Simplex. + return Boundary_simplex_range(sh); + } + + /** + * filtration_simplex_range creates an object of a Filtration_simplex_range class that provides ranges for the Filtration_simplex_iterator. + **/ + Filtration_simplex_range filtration_simplex_range() { + if (globalDbg) { + cerr << "Filtration_simplex_range filtration_simplex_range()\n"; + } + //Returns a range over the simplices of the complex in the order of the filtration + return Filtration_simplex_range(this); + } + //*********************************************************************************************************************************// + + + + //*********************************************************************************************************************************// + //Elements which are in Gudhi now, but I (and in all the cases I asked also Marc) do not understand why they are there. + //TODO -- the file IndexingTag.h in the Gudhi library contains an empty structure, so I understand that this is something that was planned (for simplicial maps?) + //but was never finished. The only idea I have here is to use the same empty structure from IndexingTag.h file, but only if the compiler needs it. If the compiler + //do not need it, then I would rather not add here elements which I do not understand. + //typedef Indexing_tag + + /** + * Function needed for compatibility with Gudhi. Not useful for other purposes. + **/ + std::pair endpoints(Simplex_handle sh) { + std::vector< size_t > bdry = this->get_boundary_of_a_cell(sh.position); + if (globalDbg) { + cerr << "std::pair endpoints( Simplex_handle sh )\n"; + cerr << "bdry.size() : " << bdry.size() << endl; + } + //this method returns two first elements from the boundary of sh. + if (bdry.size() < 2)throw ("Error in endpoints in Bitmap_cubical_complex class. The cell for which this method was called have less than two elements in the boundary."); + return std::make_pair(Simplex_handle(this, bdry[0]), Simplex_handle(this, bdry[1])); + } + + + /** + * Class needed for compatibility with Gudhi. Not useful for other purposes. + **/ + class Skeleton_simplex_range; + + class Skeleton_simplex_iterator : std::iterator< std::input_iterator_tag, Simplex_handle > { + //Iterator over all simplices of the complex in the order of the indexing scheme. + //'value_type' must be 'Simplex_handle'. + public: + + Skeleton_simplex_iterator(Bitmap_cubical_complex* b, size_t d) : b(b), dimension(d) { + if (globalDbg) { + cerr << "Skeleton_simplex_iterator ( Bitmap_cubical_complex* b , size_t d )\n"; + } + //find the position of the first simplex of a dimension d + this->position = 0; + while ((this->position != b->data.size()) && (this->b->get_dimension_of_a_cell(this->position) != this->dimension)) { + ++this->position; + } + }; + + Skeleton_simplex_iterator() : b(NULL), dimension(0) { }; + + Skeleton_simplex_iterator operator++() { + if (globalDbg) { + cerr << "Skeleton_simplex_iterator operator++()\n"; + } + //increment the position as long as you did not get to the next element of the dimension dimension. + ++this->position; + while ((this->position != this->b->data.size()) && (this->b->get_dimension_of_a_cell(this->position) != this->dimension)) { + ++this->position; + } + return (*this); + } + + Skeleton_simplex_iterator operator++(int) { + Skeleton_simplex_iterator result = *this; + ++(*this); + return result; + } + + Skeleton_simplex_iterator operator=(const Skeleton_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "Skeleton_simplex_iterator operator =\n"; + } + this->b = rhs.b; + this->position = rhs.position; + } + + bool operator==(const Skeleton_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "bool operator ==\n"; + } + if (this->position == rhs.position) { + return true; + } + return false; + } + + bool operator!=(const Skeleton_simplex_iterator& rhs) { + if (globalDbg) { + cerr << "bool operator != ( const Skeleton_simplex_iterator& rhs )\n"; + } + return !(*this == rhs); + } + + Simplex_handle operator*() { + if (globalDbg) { + cerr << "Simplex_handle operator*() \n"; + } + return Simplex_handle(this->b, this->position); + } + + friend class Skeleton_simplex_range; + private: + Bitmap_cubical_complex* b; + size_t position; + int dimension; + }; + + /** + * Class needed for compatibility with Gudhi. Not useful for other purposes. + **/ + class Skeleton_simplex_range { + //Range over the simplices of the complex in the order of the filtration. + //.begin() and .end() return type Filtration_simplex_iterator. + public: + + Skeleton_simplex_range(Bitmap_cubical_complex* b, int dimension) : b(b), dimension(dimension) { }; + + Skeleton_simplex_iterator begin() { + if (globalDbg) { + cerr << "Skeleton_simplex_iterator begin()\n"; + } + return Skeleton_simplex_iterator(this->b, this->dimension); + } + + Skeleton_simplex_iterator end() { + if (globalDbg) { + cerr << "Skeleton_simplex_iterator end()\n"; + } + Skeleton_simplex_iterator it(this->b, this->dimension); + it.position = this->b->data.size(); + return it; + } + private: + Bitmap_cubical_complex* b; + int dimension; + }; + + /** + * Function needed for compatibility with Gudhi. Not useful for other purposes. + **/ + Skeleton_simplex_range skeleton_simplex_range(int dimension) { + if (globalDbg) { + cerr << "Skeleton_simplex_range skeleton_simplex_range( int dimension )\n"; + } + return Skeleton_simplex_range(this, dimension); + } + + + + //*********************************************************************************************************************************// + //functions used for debugging: + + /** + * Function used for debugging purposes. + **/ + void printKeyAssociatedToSimplex() { + for (size_t i = 0; i != this->data.size(); ++i) { + cerr << i << " -> " << this->simplexAssociatedToKey[i] << endl; + } + } + + /** + * Function used for debugging purposes. + **/ + size_t printRealPosition(const Simplex_handle& sh) { + return sh.position; + } + + private: + std::vector< size_t > keyAssociatedToSimplex; + std::vector< size_t > simplexAssociatedToKey; + std::vector< size_t > elementsOrderedAccordingToFiltration; //needed by Filtration_simplex_iterator. If this iterator is not used, this field is not initialized. +}; //Bitmap_cubical_complex + +template +bool compareElementsForElementsOrderedAccordingToFiltration(const std::pair< size_t, std::pair< T, char > >& f, const std::pair< size_t, std::pair< T, char > >& s) { + if (globalDbg) { + cerr << "ompareElementsForElementsOrderedAccordingToFiltration\n"; + } + if (f.second.first < s.second.first) { + return true; + } else { + if (f.second.first > s.second.first) { + return false; + } else { + //in this case f.second.first == s.second.first, and we use dimension to compare: + if (f.second.second < s.second.second) { + return true; + } else { + if (f.second.second > s.second.second) { + return false; + } else { + //in this case, both the filtration value and the dimensions for those cells are the same. Since it may be nice to have a stable sorting procedure, in this case, we compare positions in the bitmap: + return ( f.first < s.first); + } + } + } + } +} + +template +void Bitmap_cubical_complex::initializeElementsOrderedAccordingToFiltration() { + if (globalDbg) { + cerr << "void Bitmap_cubical_complex::initializeElementsOrderedAccordingToFiltration() \n"; + } + //( position , (filtration , dimension) ) + std::vector< std::pair< size_t, std::pair< T, char > > > dataOfElementsFromBitmap(this->data.size()); + for (size_t i = 0; i != this->data.size(); ++i) { + //TODO -- this can be optimized by having a counter here. We do not need to re-compute the dimension for every cell from scratch + dataOfElementsFromBitmap[i] = std::make_pair(i, std::make_pair(this->data[i], this->get_dimension_of_a_cell(i))); + } + std::sort(dataOfElementsFromBitmap.begin(), dataOfElementsFromBitmap.end(), compareElementsForElementsOrderedAccordingToFiltration); + + std::vector< size_t > elementsOfBitmapOrderedAccordingToFiltrationThenAccordingToDimensionThenAccordingToPositionInBitmap(this->data.size()); + for (size_t i = 0; i != dataOfElementsFromBitmap.size(); ++i) { + elementsOfBitmapOrderedAccordingToFiltrationThenAccordingToDimensionThenAccordingToPositionInBitmap[i] = dataOfElementsFromBitmap[i].first; + } + this->elementsOrderedAccordingToFiltration = elementsOfBitmapOrderedAccordingToFiltrationThenAccordingToDimensionThenAccordingToPositionInBitmap; +} + + +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// + +#endif // BITMAP_CUBICAL_COMPLEX_H_ 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 new file mode 100644 index 00000000..26c97872 --- /dev/null +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h @@ -0,0 +1,577 @@ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): Pawel Dlotko + * + * Copyright (C) 2015 INRIA Sophia-Saclay (France) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +#ifndef BITMAP_CUBICAL_COMPLEX_BASE_H_ +#define BITMAP_CUBICAL_COMPLEX_BASE_H_ + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +using namespace std; + +/** + * This is a class implementing a basic bitmap data structure to store cubical complexes. It implements only the most basic subroutines. + * The idea of the bitmap is the following. Our aim is to have a memory efficient data structure to store d-dimensional cubical complex C being a cubical decomposition + * of a rectangular region of a space. This is achieved by storing C as a vector of bits (this is where the name 'bitmap' came from). Each cell is represented by a single + * bit (in case of black and white bitmaps, or by a single element of a type T (here T is a filtration type of a bitmap, typically a double). All the informations needed for homology and + * persistent homology computations (like dimension of a cell, boundary and coboundary elements of a cell, are then obtained from the position of the element in C. + */ +template +class Bitmap_cubical_complex_base { + public: + /** + * There are a few constructors of a Bitmap_cubical_complex_base class. First one, that takes vector, creates an empty bitmap of a dimension equal to the number of elements in the + * input vector and size in the i-th dimension equal the number in the position i-of the input vector. + */ + Bitmap_cubical_complex_base(std::vector sizes_); + /** + * The second constructor takes as a input a Perseus style file. For more details, please consult the documentations of Perseus software as well as examples attached to this + * implementation. + **/ + Bitmap_cubical_complex_base(char* perseusStyleFile_); + /** + * The last constructor of a Bitmap_cubical_complex_base class accepts vector of dimensions (as the first one) together with vector of filtration values of top dimensional cells. + **/ + Bitmap_cubical_complex_base(std::vector dimensions_, std::vector topDimensionalCells_); + + /** + * The functions get_boundary_of_a_cell, get_coboundary_of_a_cell and get_cell_data are the basic functions that compute boundary / coboundary / dimension and the filtration + * value form a position of a cell in the structure of a bitmap. The input parameter of all of those function is a non-negative integer, indicating a position of a cube in the data structure. + * In the case of functions that compute (co)boundary, the output is a vector if non-negative integers pointing to the positions of (co)boundary element of the input cell. + */ + inline std::vector< size_t > get_boundary_of_a_cell(size_t cell_); + /** + * 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 + * value form a position of a cell in the structure of a bitmap. The input parameter of all of those function is a non-negative integer, indicating a position of a cube in the data structure. + * In the case of functions that compute (co)boundary, the output is a vector if non-negative integers pointing to the positions of (co)boundary element of the input cell. + **/ + inline std::vector< size_t > get_coboundary_of_a_cell(size_t cell_); + /** + * In the case of get_dimension_of_a_cell function, the output is a non-negative integer indicating the dimension of a cell. + **/ + inline unsigned get_dimension_of_a_cell(size_t cell_); + /** + * In the case of get_cell_data, the output parameter is a reference to the value of a cube in a given position. + **/ + inline T& get_cell_data(size_t cell_); + + + /** + * Typical input used to construct a baseBitmap class is a filtration given at the top dimensional cells. Then, there are a few ways one can pick the filtration of lower dimensional + * cells. The most typical one is by so called lower star filtration. This function is always called by any constructor which takes the top dimensional cells. If you use such a constructor, + * then there is no need to call this function. Call it only if you are putting the filtration of the cells by your own (for instance by using topDimensionalCellsIterator). + **/ + void impose_lower_star_filtration(); //assume that top dimensional cells are already set. + + /** + * Returns dimension of a complex. + **/ + inline unsigned dimension() { + return sizes.size(); + } + + /** + * Returns number of all cubes in the data structure. + **/ + inline unsigned size_of_bitmap() { + return this->data.size(); + } + + /** + * Writing to stream operator. + **/ + template + friend ostream& operator<<(ostream & os_, const Bitmap_cubical_complex_base& b_); + + //ITERATORS + + /** + * 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; + + all_cells_iterator all_cells_begin()const { + return this->data.begin(); + } + + all_cells_iterator all_cells_end()const { + return this->data.end(); + } + + + typedef typename std::vector< T >::const_iterator all_cells_const_iterator; + + all_cells_const_iterator all_cells_const_begin()const { + return this->data.begin(); + } + + all_cells_const_iterator all_cells_const_end()const { + return this->data.end(); + } + + /** + * Iterator through top dimensional cells of the complex. The cells appear in order they are stored in the structure (i.e. in lexicographical order) + **/ + class Top_dimensional_cells_iterator : std::iterator< std::input_iterator_tag, double > { + public: + + Top_dimensional_cells_iterator(Bitmap_cubical_complex_base& b_) : b(b_) { + for (size_t i = 0; i != b_.dimension(); ++i) { + this->counter.push_back(0); + } + } + + Top_dimensional_cells_iterator operator++() { + //first find first element of the counter that can be increased: + size_t dim = 0; + while ((dim != this->b.dimension()) && (this->counter[dim] == this->b.sizes[dim] - 1))++dim; + + if (dim != this->b.dimension()) { + ++this->counter[dim]; + for (size_t i = 0; i != dim; ++i) { + this->counter[i] = 0; + } + } else { + ++this->counter[0]; + } + return *this; + } + + Top_dimensional_cells_iterator operator++(int) { + Top_dimensional_cells_iterator result = *this; + ++(*this); + return result; + } + + Top_dimensional_cells_iterator operator=(const Top_dimensional_cells_iterator& rhs_) { + this->counter = rhs_.counter; + this->b = rhs_.b; + return *this; + } + + bool operator==(const Top_dimensional_cells_iterator& rhs_) { + if (&this->b != &rhs_.b)return false; + if (this->counter.size() != rhs_.counter.size())return false; + for (size_t i = 0; i != this->counter.size(); ++i) { + if (this->counter[i] != rhs_.counter[i])return false; + } + return true; + } + + bool operator!=(const Top_dimensional_cells_iterator& rhs_) { + return !(*this == rhs_); + } + + T& operator*() { + //given the counter, compute the index in the array and return this element. + unsigned index = 0; + for (size_t i = 0; i != this->counter.size(); ++i) { + index += (2 * this->counter[i] + 1) * this->b.multipliers[i]; + } + return this->b.data[index]; + } + + size_t computeIndexInBitmap() { + size_t index = 0; + for (size_t i = 0; i != this->counter.size(); ++i) { + index += (2 * this->counter[i] + 1) * this->b.multipliers[i]; + } + return index; + } + + void printCounter() { + for (size_t i = 0; i != this->counter.size(); ++i) { + cout << this->counter[i] << " "; + } + } + friend class Bitmap_cubical_complex_base; + protected: + std::vector< unsigned > counter; + Bitmap_cubical_complex_base& b; + }; + + Top_dimensional_cells_iterator top_dimensional_cells_begin() { + Top_dimensional_cells_iterator a(*this); + return a; + } + + Top_dimensional_cells_iterator top_dimensional_cells_end() { + Top_dimensional_cells_iterator a(*this); + for (size_t i = 0; i != this->dimension(); ++i) { + a.counter[i] = this->sizes[i] - 1; + } + a.counter[0]++; + return a; + } + + + //****************************************************************************************************************// + //****************************************************************************************************************// + //****************************************************************************************************************// + //****************************************************************************************************************// + + + //****************************************************************************************************************// + //****************************************************************************************************************// + //****************************************************************************************************************// + //****************************************************************************************************************// + + protected: + std::vector sizes; + std::vector multipliers; + std::vector data; + size_t totalNumberOfCells; + + void set_up_containers(std::vector sizes_) { + unsigned multiplier = 1; + for (size_t i = 0; i != sizes_.size(); ++i) { + this->sizes.push_back(sizes_[i]); + this->multipliers.push_back(multiplier); + //multiplier *= 2*(sizes[i]+1)+1; + multiplier *= 2 * sizes_[i] + 1; + } + //std::reverse( this->sizes.begin() , this->sizes.end() ); + std::vector data(multiplier); + std::fill(data.begin(), data.end(), INT_MAX); + this->totalNumberOfCells = multiplier; + this->data = data; + } + + size_t compute_position_in_bitmap(std::vector< int > counter_) { + size_t position = 0; + for (size_t i = 0; i != this->multipliers.size(); ++i) { + position += this->multipliers[i] * counter_[i]; + } + return position; + } + + std::vector compute_counter_for_given_cell(size_t cell_) { + std::vector counter; + for (size_t dim = this->sizes.size(); dim != 0; --dim) { + counter.push_back(cell_ / this->multipliers[dim - 1]); + cell_ = cell_ % this->multipliers[dim - 1]; + } + std::reverse(counter.begin(), counter.end()); + return counter; + } + + std::vector< size_t > generate_vector_of_shifts_for_bitmaps_with_periodic_boundary_conditions(std::vector< bool > directionsForPeriodicBCond_); +}; + +template +ostream& operator<<(ostream & out_, const Bitmap_cubical_complex_base& b_) { + //for ( typename bitmap::all_cells_const_iterator it = b.all_cells_const_begin() ; it != b.all_cells_const_end() ; ++it ) + for (typename Bitmap_cubical_complex_base::all_cells_const_iterator it = b_.all_cells_const_begin(); it != b_.all_cells_const_end(); ++it) { + out_ << *it << " "; + } + return out_; +} + +template +Bitmap_cubical_complex_base::Bitmap_cubical_complex_base(std::vector sizes_) { + this->set_up_containers(sizes_); +} + +template +Bitmap_cubical_complex_base::Bitmap_cubical_complex_base(std::vector sizesInFollowingDirections_, std::vector topDimensionalCells_) { + this->set_up_containers(sizesInFollowingDirections_); + + size_t numberOfTopDimensionalElements = 1; + for (size_t i = 0; i != sizesInFollowingDirections_.size(); ++i) { + numberOfTopDimensionalElements *= sizesInFollowingDirections_[i]; + } + if (numberOfTopDimensionalElements != topDimensionalCells_.size()) { + cerr << "Error in constructor Bitmap_cubical_complex_base( std::vector sizesInFollowingDirections_ , std::vector topDimensionalCells_ ). Number of top dimensional elements that follow from sizesInFollowingDirections vector is different than the size of topDimensionalCells vector." << endl; + throw ("Error in constructor Bitmap_cubical_complex_base( std::vector sizesInFollowingDirections_ , std::vector topDimensionalCells_ ). Number of top dimensional elements that follow from sizesInFollowingDirections vector is different than the size of topDimensionalCells vector."); + } + + Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); + size_t index = 0; + for (it = this->top_dimensional_cells_begin(); it != this->top_dimensional_cells_end(); ++it) { + (*it) = topDimensionalCells_[index]; + ++index; + } + this->impose_lower_star_filtration(); +} + +template +Bitmap_cubical_complex_base::Bitmap_cubical_complex_base(char* perseusStyleFile_) { + bool dbg = false; + ifstream inFiltration, inIds; + inFiltration.open(perseusStyleFile_); + unsigned dimensionOfData; + inFiltration >> dimensionOfData; + + if (dbg) { + cerr << "dimensionOfData : " << dimensionOfData << endl; + } + + std::vector sizes; + for (size_t i = 0; i != dimensionOfData; ++i) { + int sizeInThisDimension; + inFiltration >> sizeInThisDimension; + sizeInThisDimension = abs(sizeInThisDimension); + sizes.push_back(sizeInThisDimension); + if (dbg) { + cerr << "sizeInThisDimension : " << sizeInThisDimension << endl; + } + } + this->set_up_containers(sizes); + + Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); + it = this->top_dimensional_cells_begin(); + + //TODO -- over here we also need to read id's of cell and put them to bitmapElement structure! + while (!inFiltration.eof()) { + double filtrationLevel; + inFiltration >> filtrationLevel; + if (dbg) { + cerr << "Cell of an index : " << it.computeIndexInBitmap() << " and dimension: " << this->get_dimension_of_a_cell(it.computeIndexInBitmap()) << " get the value : " << filtrationLevel << endl; + } + *it = filtrationLevel; + ++it; + } + inFiltration.close(); + this->impose_lower_star_filtration(); +} + +template +std::vector< size_t > Bitmap_cubical_complex_base::get_boundary_of_a_cell(size_t cell_) { + bool bdg = false; + //first of all, we need to take the list of coordinates in which the cell has nonzero length. We do it by using modified version to compute dimension of a cell: + std::vector< unsigned > dimensionsInWhichCellHasNonzeroLength; + unsigned dimension = 0; + size_t cell1 = cell_; + for (size_t i = this->multipliers.size(); i != 0; --i) { + unsigned position = cell1 / multipliers[i - 1]; + if (position % 2 == 1) { + dimensionsInWhichCellHasNonzeroLength.push_back(i - 1); + dimension++; + } + cell1 = cell1 % multipliers[i - 1]; + } + + if (bdg) { + cerr << "dimensionsInWhichCellHasNonzeroLength : \n"; + for (size_t i = 0; i != dimensionsInWhichCellHasNonzeroLength.size(); ++i) { + cerr << dimensionsInWhichCellHasNonzeroLength[i] << endl; + } + getchar(); + } + + std::vector< size_t > boundaryElements; + if (dimensionsInWhichCellHasNonzeroLength.size() == 0)return boundaryElements; + for (size_t i = 0; i != dimensionsInWhichCellHasNonzeroLength.size(); ++i) { + boundaryElements.push_back(cell_ - multipliers[ dimensionsInWhichCellHasNonzeroLength[i] ]); + boundaryElements.push_back(cell_ + multipliers[ dimensionsInWhichCellHasNonzeroLength[i] ]); + + if (bdg) cerr << "multipliers[dimensionsInWhichCellHasNonzeroLength[i]] : " << multipliers[dimensionsInWhichCellHasNonzeroLength[i]] << endl; + if (bdg) cerr << "cell_ - multipliers[dimensionsInWhichCellHasNonzeroLength[i]] : " << cell_ - multipliers[dimensionsInWhichCellHasNonzeroLength[i]] << endl; + if (bdg) cerr << "cell_ + multipliers[dimensionsInWhichCellHasNonzeroLength[i]] : " << cell_ + multipliers[dimensionsInWhichCellHasNonzeroLength[i]] << endl; + } + return boundaryElements; +} + +template +std::vector< size_t > Bitmap_cubical_complex_base::get_coboundary_of_a_cell(size_t cell_) { + bool bdg = false; + //first of all, we need to take the list of coordinates in which the cell has nonzero length. We do it by using modified version to compute dimension of a cell: + std::vector< unsigned > dimensionsInWhichCellHasZeroLength; + unsigned dimension = 0; + size_t cell1 = cell_; + for (size_t i = this->multipliers.size(); i != 0; --i) { + unsigned position = cell1 / multipliers[i - 1]; + if (position % 2 == 0) { + dimensionsInWhichCellHasZeroLength.push_back(i - 1); + dimension++; + } + cell1 = cell1 % multipliers[i - 1]; + } + + std::vector counter = this->compute_counter_for_given_cell(cell_); + //reverse(counter.begin() , counter.end()); + + if (bdg) { + cerr << "dimensionsInWhichCellHasZeroLength : \n"; + for (size_t i = 0; i != dimensionsInWhichCellHasZeroLength.size(); ++i) { + cerr << dimensionsInWhichCellHasZeroLength[i] << endl; + } + cerr << "\n counter : " << endl; + for (size_t i = 0; i != counter.size(); ++i) { + cerr << counter[i] << endl; + } + getchar(); + } + + std::vector< size_t > coboundaryElements; + if (dimensionsInWhichCellHasZeroLength.size() == 0)return coboundaryElements; + for (size_t i = 0; i != dimensionsInWhichCellHasZeroLength.size(); ++i) { + if (bdg) { + cerr << "Dimension : " << i << endl; + if (counter[dimensionsInWhichCellHasZeroLength[i]] == 0) { + cerr << "In dimension : " << i << " we cannot substract, since we will jump out of a Bitmap_cubical_complex_base \n"; + } + if (counter[dimensionsInWhichCellHasZeroLength[i]] == 2 * this->sizes[dimensionsInWhichCellHasZeroLength[i]]) { + cerr << "In dimension : " << i << " we cannot substract, since we will jump out of a Bitmap_cubical_complex_base \n"; + } + } + + + if ((cell_ > multipliers[dimensionsInWhichCellHasZeroLength[i]]) && (counter[dimensionsInWhichCellHasZeroLength[i]] != 0)) + //if ( counter[dimensionsInWhichCellHasZeroLength[i]] != 0 ) + { + if (bdg)cerr << "Subtracting : " << cell_ - multipliers[dimensionsInWhichCellHasZeroLength[i]] << endl; + coboundaryElements.push_back(cell_ - multipliers[dimensionsInWhichCellHasZeroLength[i]]); + } + if ((cell_ + multipliers[dimensionsInWhichCellHasZeroLength[i]] < this->data.size()) && (counter[dimensionsInWhichCellHasZeroLength[i]] != 2 * this->sizes[dimensionsInWhichCellHasZeroLength[i]])) + //if ( counter[dimensionsInWhichCellHasZeroLength[i]] != 2*this->sizes[dimensionsInWhichCellHasZeroLength[i]] ) + { + coboundaryElements.push_back(cell_ + multipliers[dimensionsInWhichCellHasZeroLength[i]]); + if (bdg)cerr << "Adding : " << cell_ + multipliers[dimensionsInWhichCellHasZeroLength[i]] << endl; + } + } + return coboundaryElements; +} + +template +unsigned Bitmap_cubical_complex_base::get_dimension_of_a_cell(size_t cell_) { + bool dbg = false; + if (dbg)cerr << "\n\n\n Computing position o a cell of an index : " << cell_ << endl; + unsigned dimension = 0; + for (size_t i = this->multipliers.size(); i != 0; --i) { + unsigned position = cell_ / multipliers[i - 1]; + + if (dbg)cerr << "i-1 :" << i - 1 << endl; + if (dbg)cerr << "cell_ : " << cell_ << endl; + if (dbg)cerr << "position : " << position << endl; + if (dbg)cerr << "multipliers[" << i - 1 << "] = " << multipliers[i - 1] << endl; + if (dbg)getchar(); + + if (position % 2 == 1) { + if (dbg)cerr << "Nonzero length in this direction \n"; + dimension++; + } + cell_ = cell_ % multipliers[i - 1]; + } + return dimension; +} + +template +T& Bitmap_cubical_complex_base::get_cell_data(size_t cell_) { + return this->data[cell_]; +} + +template +void Bitmap_cubical_complex_base::impose_lower_star_filtration() { + bool dbg = false; + + //this vector will be used to check which elements have already been taken care of in imposing lower star filtration: + std::vector isThisCellConsidered(this->data.size(), false); + + std::vector indicesToConsider; + //we assume here that we already have a filtration on the top dimensional cells and we have to extend it to lower ones. + typename Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); + for (it = this->top_dimensional_cells_begin(); it != this->top_dimensional_cells_end(); ++it) { + indicesToConsider.push_back(it.computeIndexInBitmap()); + } + + while (indicesToConsider.size()) { + if (dbg) { + cerr << "indicesToConsider in this iteration \n"; + for (size_t i = 0; i != indicesToConsider.size(); ++i) { + cout << indicesToConsider[i] << " "; + } + getchar(); + } + std::vector newIndicesToConsider; + for (size_t i = 0; i != indicesToConsider.size(); ++i) { + std::vector bd = this->get_boundary_of_a_cell(indicesToConsider[i]); + for (size_t boundaryIt = 0; boundaryIt != bd.size(); ++boundaryIt) { + if (this->data[ bd[boundaryIt] ] > this->data[ indicesToConsider[i] ]) { + this->data[ bd[boundaryIt] ] = this->data[ indicesToConsider[i] ]; + } + if (isThisCellConsidered[ bd[boundaryIt] ] == false) { + newIndicesToConsider.push_back(bd[boundaryIt]); + isThisCellConsidered[ bd[boundaryIt] ] = true; + } + } + } + indicesToConsider.swap(newIndicesToConsider); + } +} + +template +std::vector< size_t > Bitmap_cubical_complex_base::generate_vector_of_shifts_for_bitmaps_with_periodic_boundary_conditions(std::vector< bool > directionsForPeriodicBCond_) { + bool dbg = false; + if (this->sizes.size() != directionsForPeriodicBCond_.size())throw "directionsForPeriodicBCond_ vector size is different from the size of the bitmap. The program will now terminate \n"; + + std::vector sizes(this->sizes.size()); + for (size_t i = 0; i != this->sizes.size(); ++i)sizes[i] = 2 * this->sizes[i]; + + counter c(sizes); + + std::vector< size_t > result; + + for (size_t i = 0; i != this->data.size(); ++i) { + size_t position; + if (!c.isFinal()) { + position = i; + //result.push_back( i ); + } else { + std::vector< bool > finals = c.directionsOfFinals(); + bool jumpInPosition = false; + for (size_t dir = 0; dir != finals.size(); ++dir) { + if (finals[dir] == false)continue; + if (directionsForPeriodicBCond_[dir]) { + jumpInPosition = true; + } + } + if (jumpInPosition == true) { + //in this case this guy is final, so we need to find 'the opposite one' + position = compute_position_in_bitmap(c.findOpposite(directionsForPeriodicBCond_)); + } else { + position = i; + } + } + result.push_back(position); + if (dbg) { + cerr << " position : " << position << endl; + cerr << c << endl; + getchar(); + } + + c.increment(); + } + + return result; +} + +#endif // BITMAP_CUBICAL_COMPLEX_BASE_H_ diff --git a/src/Bitmap_cubical_complex/include/gudhi/counter.h b/src/Bitmap_cubical_complex/include/gudhi/counter.h new file mode 100644 index 00000000..9df819b2 --- /dev/null +++ b/src/Bitmap_cubical_complex/include/gudhi/counter.h @@ -0,0 +1,136 @@ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): Pawel Dlotko + * + * Copyright (C) 2015 INRIA Sophia-Saclay (France) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +#ifndef COUNTER_H_ +#define COUNTER_H_ + +#include +#include + +using namespace std; + +/** + * This is an implementation of a simple counter. It is needed for the implementation of a bitmapCubicalComplex. + **/ + +class counter { + public: + + /** + * Constructor of a counter class. It takes only the parameter which is the end value of the counter. The default beginning value is a vector of the same length as the endd, filled-in with zeros. + **/ + counter(std::vector< int > endd) { + for (size_t i = 0; i != endd.size(); ++i) { + this->current.push_back(0); + this->begin.push_back(0); + this->end.push_back(endd[i]); + } + } + + /** + * Constructor of a counter class. It takes as the input beginn and end vector. It assumes that begin vector is lexicographically below the end vector. + **/ + counter(std::vector< int > beginn, std::vector< int > endd) { + if (beginn.size() != endd.size())throw "In constructor of a counter, begin and end vectors do not have the same size. Program terminate"; + for (size_t i = 0; i != endd.size(); ++i) { + this->current.push_back(0); + this->begin.push_back(0); + this->end.push_back(endd[i]); + } + } + + /** + * Function to increment the counter. If the value returned by the function is true, then the incrementation process was successful. + * If the value of the function is false, that means, that the counter have reached its end-value. + **/ + bool increment() { + size_t i = 0; + while ((i != this->end.size()) && (this->current[i] == this->end[i])) { + ++i; + } + + if (i == this->end.size())return false; + ++this->current[i]; + for (size_t j = 0; j != i; ++j) { + this->current[j] = this->begin[j]; + } + return true; + } + + /** + * Function to check if we are at the end of counter. + **/ + bool isFinal() { + for (size_t i = 0; i != this->current.size(); ++i) { + if (this->current[i] == this->end[i])return true; + } + return false; + } + + /** + * Function required in the implementation of bitmapCubicalComplexWPeriodicBoundaryCondition. Its aim is to find an counter corresponding to the element the following + * boundary element is identified with when periodic boundary conditions are imposed. + **/ + std::vector< int > findOpposite(std::vector< bool > directionsForPeriodicBCond) { + std::vector< int > result; + for (size_t i = 0; i != this->current.size(); ++i) { + if ((this->current[i] == this->end[i]) && (directionsForPeriodicBCond[i] == true)) { + result.push_back(this->begin[i]); + } else { + result.push_back(this->current[i]); + } + } + return result; + } + + /** + * Function checking at which positions the current value of a counter is the final value of the counter. + **/ + std::vector< bool > directionsOfFinals() { + std::vector< bool > result; + for (size_t i = 0; i != this->current.size(); ++i) { + if (this->current[i] == this->end[i]) { + result.push_back(true); + } else { + result.push_back(false); + } + } + return result; + } + + /** + * Function to write counter to the stream. + **/ + friend std::ostream& operator<<(std::ostream& out, const counter& c) { + //cerr << "c.current.size() : " << c.current.size() << endl; + for (size_t i = 0; i != c.current.size(); ++i) { + out << c.current[i] << " "; + } + return out; + } + private: + std::vector< int > begin; + std::vector< int > end; + std::vector< int > current; +}; + +#endif // COUNTER_H_ -- cgit v1.2.3