Adding PSSK implmentation (based on persistent heat maps).

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

Former-commit-id: 2b5a6fc2fa34d098cc53b53e6cee0d3dd65193db

9 files changed, 583 insertions, 2 deletions

diff --git a/src/Gudhi_stat/example/CMakeLists.txt b/src/Gudhi_stat/example/CMakeLists.txt
index e7573cc6..b0c3b15f 100644
--- a/src/Gudhi_stat/example/CMakeLists.txt
+++ b/src/Gudhi_stat/example/CMakeLists.txt
@@ -102,3 +102,20 @@ target_link_libraries(utilities/persistence_vectors/compute_scalar_product_of_pe
 add_executable ( utilities/persistence_vectors/plot_persistence_vectors utilities/persistence_vectors/plot_persistence_vectors.cpp )
 target_link_libraries(utilities/persistence_vectors/plot_persistence_vectors ${Boost_SYSTEM_LIBRARY})
 
+
+#PSSK
+add_executable ( utilities/PSSK/average_pssk utilities/PSSK/average_pssk.cpp )
+target_link_libraries(utilities/PSSK/average_pssk ${Boost_SYSTEM_LIBRARY})
+
+add_executable ( utilities/PSSK/create_pssk utilities/PSSK/create_pssk.cpp )
+target_link_libraries(utilities/PSSK/create_pssk ${Boost_SYSTEM_LIBRARY})
+
+add_executable ( utilities/PSSK/plot_pssk utilities/PSSK/plot_pssk.cpp )
+target_link_libraries(utilities/PSSK/plot_pssk ${Boost_SYSTEM_LIBRARY})
+
+add_executable ( utilities/PSSK/compute_distance_of_pssk utilities/PSSK/compute_distance_of_pssk.cpp )
+target_link_libraries(utilities/PSSK/compute_distance_of_pssk ${Boost_SYSTEM_LIBRARY})
+
+add_executable ( utilities/PSSK/compute_scalar_product_of_pssk utilities/PSSK/compute_scalar_product_of_pssk.cpp )
+target_link_libraries(utilities/PSSK/compute_scalar_product_of_pssk ${Boost_SYSTEM_LIBRARY})
+
diff --git a/src/Gudhi_stat/example/utilities/PSSK/average_pssk.cpp b/src/Gudhi_stat/example/utilities/PSSK/average_pssk.cpp
new file mode 100644
index 00000000..08ea60d3
--- /dev/null
+++ b/src/Gudhi_stat/example/utilities/PSSK/average_pssk.cpp
@@ -0,0 +1,73 @@
+/*    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 <http://www.gnu.org/licenses/>.
+ */
+
+
+
+//gudhi include
+#include <gudhi/concretizations/PSSK.h>
+
+
+
+using namespace Gudhi;
+using namespace Gudhi::Gudhi_stat;
+
+#include <iostream>
+
+
+int main( int argc , char** argv )
+{
+	std::cout << "This program computes average persistence landscape of persistence landscapes created based on persistence diagrams provided as an input. Please call this program with the names of files with persistence diagrams \n";
+	std::vector< const char* > filenames;
+	
+	if ( argc == 1 )
+	{
+		std::cout << "No input files given, the program will now terminate \n";
+		return 1;
+	}
+	
+	for ( int i = 1 ; i < argc ; ++i )
+	{
+		filenames.push_back( argv[i] );
+	}
+	
+	std::cout << "Creating persistence landscapes...\n";
+	std::vector< Abs_Topological_data_with_averages* > maps;
+	for ( size_t i = 0 ; i != filenames.size() ; ++i )
+	{
+		PSSK* l = new PSSK;
+		l->load_from_file( filenames[i] );
+		maps.push_back( (Abs_Topological_data_with_averages*)l );
+	}
+	
+	PSSK av;
+	av.compute_average( maps );
+	
+	av.print_to_file( "average.pssk" );
+	
+	for ( size_t i = 0 ; i != filenames.size() ; ++i )
+	{
+		delete maps[i];
+	}
+	
+	std::cout << "Done \n";
+	return 0;
+}
diff --git a/src/Gudhi_stat/example/utilities/PSSK/compute_distance_of_pssk.cpp b/src/Gudhi_stat/example/utilities/PSSK/compute_distance_of_pssk.cpp
new file mode 100644
index 00000000..8eacbd70
--- /dev/null
+++ b/src/Gudhi_stat/example/utilities/PSSK/compute_distance_of_pssk.cpp
@@ -0,0 +1,107 @@
+/*    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 <http://www.gnu.org/licenses/>.
+ */
+
+
+
+//gudhi include
+#include <gudhi/concretizations/PSSK.h>
+
+
+
+using namespace Gudhi;
+using namespace Gudhi::Gudhi_stat;
+
+#include <iostream>
+#include <sstream>
+
+
+int main( int argc , char** argv )
+{
+	std::cout << "This program compute dsitance of persistence landscapes stored in a file (the file needs to be created beforehand). \n";	
+	std::cout << "The first parameter of a program is an interger p. The program compute L^p distance of the two landscapes. For L^infty distance choose p = -1. \n";
+	std::cout << "The remaining parameters of this programs are names of files with persistence landscapes.\n";
+	
+	if ( argc < 3 )
+	{
+		std::cout << "Wrong number of parameters, the program will now terminate \n";
+		return 1;
+	}
+	
+	int p = atoi( argv[1] );
+
+	std::vector< const char* > filenames;
+	for ( int i = 2 ; i < argc ; ++i )
+	{
+		filenames.push_back( argv[i] );
+	}
+	std::vector< Abs_Topological_data_with_distances* > maps;
+	maps.reserve( filenames.size() );
+	for ( size_t file_no = 0 ; file_no != filenames.size() ; ++file_no )
+	{
+		PSSK* l = new PSSK;
+		l->load_from_file( filenames[file_no] );
+		maps.push_back( (Abs_Topological_data_with_distances*)l );
+	}
+		
+	//and now we will compute the scalar product of landscapes.
+	
+	//first we prepare an array:
+	std::vector< std::vector< double > > distance( filenames.size() );
+	for ( size_t i = 0 ; i != filenames.size() ; ++i )
+	{
+		std::vector< double > v( filenames.size() , 0 );
+		distance[i] = v;
+	}
+	
+	//and now we can compute the distnaces:
+	for ( size_t i = 0 ; i != filenames.size() ; ++i )
+	{
+		for ( size_t j = i ; j != filenames.size() ; ++j )
+		{
+			distance[i][j] = distance[j][i] = ((Persistence_heat_maps*)maps[i])->distance( maps[j] , p ) ;
+		}
+	}
+	
+	//and now output the result to the screen and a file:
+	ofstream out;
+	out.open( "distance" );
+	for ( size_t i = 0 ; i != distance.size() ; ++i )
+	{
+		for ( size_t j = 0 ; j != distance.size() ; ++j )
+		{
+			cout << distance[i][j] << " ";
+			out << distance[i][j] << " ";
+		}
+		cout << endl;
+		out << endl;
+	}
+	out.close();
+		
+	return 0;
+}
+
+
+
+
+
+
+	
diff --git a/src/Gudhi_stat/example/utilities/PSSK/compute_scalar_product_of_pssk.cpp b/src/Gudhi_stat/example/utilities/PSSK/compute_scalar_product_of_pssk.cpp
new file mode 100644
index 00000000..048b581c
--- /dev/null
+++ b/src/Gudhi_stat/example/utilities/PSSK/compute_scalar_product_of_pssk.cpp
@@ -0,0 +1,90 @@
+/*    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 <http://www.gnu.org/licenses/>.
+ */
+
+
+//gudhi include
+#include <gudhi/concretizations/PSSK.h>
+
+
+
+using namespace Gudhi;
+using namespace Gudhi::Gudhi_stat;
+
+#include <iostream>
+#include <sstream>
+
+
+int main( int argc , char** argv )
+{
+	std::cout << "This program compute scalar product of persistence landscapes stored in a file (the file needs to be created beforehand). \n";	
+	std::cout << "The parameters of this programs are names of files with persistence landscapes.\n";
+
+	std::vector< const char* > filenames;
+	for ( int i = 1 ; i < argc ; ++i )
+	{
+		filenames.push_back( argv[i] );
+	}
+	std::vector< PSSK* > maps;
+	maps.reserve( filenames.size() );
+	for ( size_t file_no = 0 ; file_no != filenames.size() ; ++file_no )
+	{
+		PSSK* l = new PSSK;
+		l->load_from_file( filenames[file_no] );
+		maps.push_back( l );
+	}
+		
+	//and now we will compute the scalar product of landscapes.
+	
+	//first we prepare an array:
+	std::vector< std::vector< double > > scalar_product( filenames.size() );
+	for ( size_t i = 0 ; i != filenames.size() ; ++i )
+	{
+		std::vector< double > v( filenames.size() , 0 );
+		scalar_product[i] = v;
+	}
+	
+	//and now we can compute the scalar product:
+	for ( size_t i = 0 ; i != maps.size() ; ++i )
+	{
+		for ( size_t j = i ; j != maps.size() ; ++j )
+		{
+			scalar_product[i][j] = scalar_product[j][i] = ((Persistence_heat_maps*)maps[i])->compute_scalar_product( maps[j] ) ;
+		}
+	}
+	
+	//and now output the result to the screen and a file:
+	ofstream out;
+	out.open( "scalar_product" );
+	for ( size_t i = 0 ; i != scalar_product.size() ; ++i )
+	{
+		for ( size_t j = 0 ; j != scalar_product.size() ; ++j )
+		{
+			cout << scalar_product[i][j] << " ";
+			out << scalar_product[i][j] << " ";
+		}
+		cout << endl;
+		out << endl;
+	}
+	out.close();
+		
+	return 0;
+}
diff --git a/src/Gudhi_stat/example/utilities/PSSK/create_pssk.cpp b/src/Gudhi_stat/example/utilities/PSSK/create_pssk.cpp
new file mode 100644
index 00000000..ccd884ec
--- /dev/null
+++ b/src/Gudhi_stat/example/utilities/PSSK/create_pssk.cpp
@@ -0,0 +1,77 @@
+/*    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 <http://www.gnu.org/licenses/>.
+ */
+
+
+
+#include <gudhi/abstract_classes/Abs_Topological_data.h>
+#include <gudhi/concretizations/PSSK.h>
+
+
+
+using namespace Gudhi;
+using namespace Gudhi::Gudhi_stat;
+
+#include <iostream>
+#include <sstream>
+
+
+
+int main( int argc , char** argv )
+{
+	std::cout << "This program creates PSSK of diagrams provided as an input.\n";
+	std::cout << "The first parameter of a program is an integer, a size of a grid.\n";
+	std::cout << "The second and third parameters are min and max of the grid. If you want those numbers to be computed based on the data, set them both to -1 \n";
+	std::cerr << "The fourth parameter is an integer, the standard deviation of a gaussian kernel expressed in a number of pixels \n";
+	std::cout << "The remaining parameters are the names of files with persistence diagrams. \n";
+	
+	if ( argc < 5 )
+	{
+		std::cout << "Wrong parameter list, the program will now terminate \n";
+		return 1;
+	}
+	
+	size_t size_of_grid = (size_t)atoi( argv[1] );
+	double min_ = atof( argv[2] );
+	double max_ = atof( argv[3] );
+	size_t stdiv = atof( argv[4] );
+	
+	std::vector< const char* > filenames;
+	for ( int i = 5 ; i < argc ; ++i )
+	{
+		filenames.push_back( argv[i] );
+	}	
+	
+	std::cout << "Creating persistence heat maps...\n";	
+	std::vector< std::vector<double> > filter = create_Gaussian_filter(stdiv,1);	
+	
+	for ( size_t i = 0 ; i != filenames.size() ; ++i )
+	{		
+		std::cout << "Creating a heat map based on a file : " << filenames[i] << std::endl;
+		PSSK l( filenames[i] , filter , size_of_grid , min_ , max_ );		
+		
+		std::stringstream ss;
+		ss << filenames[i] << ".pssk";
+		l.print_to_file( ss.str().c_str() );
+	}
+	std::cout << "Done \n";
+	return 0;
+}
diff --git a/src/Gudhi_stat/example/utilities/PSSK/plot_pssk.cpp b/src/Gudhi_stat/example/utilities/PSSK/plot_pssk.cpp
new file mode 100644
index 00000000..141024a5
--- /dev/null
+++ b/src/Gudhi_stat/example/utilities/PSSK/plot_pssk.cpp
@@ -0,0 +1,44 @@
+/*    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 <http://www.gnu.org/licenses/>.
+ */
+
+
+
+//gudhi include
+#include <gudhi/concretizations/PSSK.h>
+
+
+
+using namespace Gudhi;
+using namespace Gudhi::Gudhi_stat;
+
+#include <iostream>
+#include <sstream>
+
+
+int main( int argc , char** argv )
+{
+	std::cout << "This program plot persistence landscape stored in a file (the file needs to be created beforehand). Please call the code with the name of a landsape file \n";	
+	PSSK l;
+	l.load_from_file( argv[1] );		
+	l.plot( argv[1] );	
+	return 0;
+}
diff --git a/src/Gudhi_stat/example/utilities/persistence_heat_maps/create_persistence_heat_maps.cpp b/src/Gudhi_stat/example/utilities/persistence_heat_maps/create_persistence_heat_maps.cpp
index 14a25918..3d8dc2e8 100644
--- a/src/Gudhi_stat/example/utilities/persistence_heat_maps/create_persistence_heat_maps.cpp
+++ b/src/Gudhi_stat/example/utilities/persistence_heat_maps/create_persistence_heat_maps.cpp
@@ -37,7 +37,7 @@ using namespace Gudhi::Gudhi_stat;
 
 int main( int argc , char** argv )
 {
-	std::cout << "This program creates persistence landscape on grid of diagrams provided as an input.\n";
+	std::cout << "This program creates persistence heat map of diagrams provided as an input.\n";
 	std::cout << "The first parameter of a program is an integer, a size of a grid.\n";
 	std::cout << "The second and third parameters are min and max of the grid. If you want those numbers to be computed based on the data, set them both to -1 \n";
 	std::cerr << "The fourth parameter is an integer, the standard deviation of a gaussian kernel expressed in a number of pixels \n";
diff --git a/src/Gudhi_stat/include/gudhi/concretizations/PSSK.h b/src/Gudhi_stat/include/gudhi/concretizations/PSSK.h
new file mode 100644
index 00000000..61597ee2
--- /dev/null
+++ b/src/Gudhi_stat/include/gudhi/concretizations/PSSK.h
@@ -0,0 +1,173 @@
+/*    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 <http://www.gnu.org/licenses/>.
+ */
+
+#pragma once
+#ifndef PSSK_H
+#define PSSK_H
+
+//gudhi include
+#include <gudhi/concretizations/Persistence_heat_maps.h>
+
+
+namespace Gudhi
+{
+namespace Gudhi_stat
+{
+
+/**
+* This is a version of a representation presented in https://arxiv.org/abs/1412.6821
+* In that paper the authors are using the representation just to compute kernel. Over here, we extend the usability by far.
+* Note that the version presented here is not exact, since we are discretizing the kernel.
+* The only difference with respect to the original class is the method of creation. We have full (square) image, and for every point (p,q), we add a kernel at (p,q) and the negative kernel
+* at (q,p)
+**/
+
+class PSSK : public Persistence_heat_maps
+{
+public:
+	PSSK():Persistence_heat_maps(){}
+	
+	PSSK(const std::vector< std::pair< double,double > >  & interval , std::vector< std::vector<double> > filter = create_Gaussian_filter(5,1) , size_t number_of_pixels = 1000 , double min_ = -1 , double max_ = -1 )
+	:Persistence_heat_maps()
+	{
+		this->construct( interval  , filter , number_of_pixels , min_ , max_ );
+	}
+	
+	 
+	 PSSK( const char* filename , std::vector< std::vector<double> > filter = create_Gaussian_filter(5,1) , size_t number_of_pixels = 1000 , double min_ = -1 , double max_ = -1  ):
+	    Persistence_heat_maps()
+	    {
+			std::vector< std::pair< double , double > > intervals_ = read_standard_file( filename );         
+			this->construct( intervals_ , filter , number_of_pixels , min_ , max_ );
+		}
+
+protected:
+	void construct( const std::vector< std::pair<double,double> >& intervals_  , 
+					std::vector< std::vector<double> > filter = create_Gaussian_filter(5,1),
+                    size_t number_of_pixels = 1000 , double min_ = -1 , double max_ = -1 );
+};
+
+//if min_ == max_, then the program is requested to set up the values itself based on persistence intervals
+void PSSK::construct( const std::vector< std::pair<double,double> >& intervals_  ,  
+									   std::vector< std::vector<double> > filter,
+									   size_t number_of_pixels , double min_ , double max_ )
+{	
+    bool dbg = false;       
+    if ( dbg ){std::cerr << "Entering construct procedure \n";getchar();}
+
+    if ( min_ == max_ )
+    {
+        //in this case, we want the program to set up the min_ and max_ values by itself.
+        min_ = std::numeric_limits<int>::max();
+        max_ = -std::numeric_limits<int>::max();
+        
+        
+        for ( size_t i = 0 ; i != intervals_.size() ; ++i )
+        {
+			if ( intervals_[i].first < min_ )min_ = intervals_[i].first;
+			if ( intervals_[i].second > max_ )max_ = intervals_[i].second;
+		}
+        //now we have the structure filled in, and moreover we know min_ and max_ values of the interval, so we know the range.
+
+        //add some more space:
+        min_ -= fabs(max_ - min_)/100;
+        max_ += fabs(max_ - min_)/100;
+    }
+
+	if ( dbg )
+	{
+		cerr << "min_ : " << min_ << endl;
+		cerr << "max_ : " << max_ << endl;
+		cerr << "number_of_pixels : " << number_of_pixels << endl;
+		getchar();
+	}
+
+    this->min_ = min_;
+    this->max_ = max_;    
+    this->scalling_function_with_respect_to_distance_from_diagonal = scalling_function_with_respect_to_distance_from_diagonal;
+
+
+    //initialization of the structure heat_map
+    std::vector< std::vector<double> > heat_map_;
+    for ( size_t i = 0 ; i != number_of_pixels ; ++i )
+    {       
+		std::vector<double> v( number_of_pixels , 0 );
+        heat_map_.push_back( v );
+    }
+    this->heat_map = heat_map_;
+
+    if (dbg)cerr << "Done creating of the heat map, now we will fill in the structure \n";
+
+	for ( size_t pt_nr = 0 ; pt_nr != intervals_.size() ; ++pt_nr )
+	{
+		//compute the value of intervals_[pt_nr] in the grid:
+		int x_grid = (int)((intervals_[pt_nr].first - this->min_)/( this->max_-this->min_ )*number_of_pixels);
+		int y_grid = (int)((intervals_[pt_nr].second - this->min_)/( this->max_-this->min_ )*number_of_pixels);
+		
+		if ( dbg )
+		{
+			std::cerr << "point : " << intervals_[pt_nr].first << " , " << intervals_[pt_nr].second << endl;
+			std::cerr << "x_grid : " << x_grid << endl;
+			std::cerr << "y_grid : " << y_grid << endl;
+		}
+		
+		//x_grid and y_grid gives a center of the kernel. We want to have its lower left cordner. To get this, we need to shift x_grid and y_grid by a grid diameter.		
+		x_grid -= filter.size()/2;
+		y_grid -= filter.size()/2;
+		//note that the numbers x_grid and y_grid may be negative. 
+		
+		if ( dbg )
+		{
+			std::cerr << "After shift : \n";;
+			std::cerr << "x_grid : " << x_grid << endl;
+			std::cerr << "y_grid : " << y_grid << endl;
+			std::cerr << "filter.size() : " << filter.size() << std::endl;
+			getchar();
+		}					
+
+		
+		for ( size_t i = 0 ; i != filter.size() ; ++i )
+		{
+			for ( size_t j = 0 ; j != filter.size() ; ++j )
+			{			
+				//if the point (x_grid+i,y_grid+j) is the correct point in the grid.						
+				if ( 
+					  ((x_grid+i)>=0) && (x_grid+i<this->heat_map.size()) 
+					  &&
+					  ((y_grid+j)>=0) && (y_grid+j<this->heat_map.size()) 
+				   )
+				{
+					if ( dbg ){std::cerr << y_grid+j << " " <<  x_grid+i << std::endl;}
+					this->heat_map[ y_grid+j ][ x_grid+i ] += filter[i][j];
+					this->heat_map[ x_grid+i ][ y_grid+j ] += -filter[i][j];
+				}
+			}
+		}
+		
+	}
+}//construct
+
+
+#endif
+
+}//namespace Gudhi_stat
+}//namespace Gudhi
diff --git a/src/Gudhi_stat/include/gudhi/concretizations/Persistence_heat_maps.h b/src/Gudhi_stat/include/gudhi/concretizations/Persistence_heat_maps.h
index 6d06ac8e..9d00a221 100644
--- a/src/Gudhi_stat/include/gudhi/concretizations/Persistence_heat_maps.h
+++ b/src/Gudhi_stat/include/gudhi/concretizations/Persistence_heat_maps.h
@@ -314,7 +314,7 @@ public:
 	void compute_average( std::vector< Abs_Topological_data_with_averages* > to_average );	
 	double compute_scalar_product( const Abs_Topological_data_with_scalar_product* second_ );
 	
-private:
+protected:
 	//private methods
 	std::vector< std::vector<double> > check_and_initialize_maps( const std::vector<Persistence_heat_maps*>& maps );
     void construct( const std::vector< std::pair<double,double> >& intervals_  ,