diff --git a/Core/PythonAPI/inc/vdouble2d_t.pypp.h b/Core/PythonAPI/inc/vdouble2d_t.pypp.h
deleted file mode 100644
index bcf3743f8517a7623d5f8cfbf70e57f74b6f1f3d..0000000000000000000000000000000000000000
--- a/Core/PythonAPI/inc/vdouble2d_t.pypp.h
+++ /dev/null
@@ -1,23 +0,0 @@
-// This file has been generated by Py++.
-
-// ************************************************************************** //
-//
-//  BornAgain: simulate and fit scattering at grazing incidence
-//
-//! @file      Automatically generated boost::python code for BornAgain Python bindings
-//! @brief     Automatically generated boost::python code for BornAgain Python bindings
-//!
-//! @homepage  http://bornagainproject.org
-//! @license   GNU General Public License v3 or higher (see COPYING)
-//! @copyright Forschungszentrum Juelich GmbH 2015
-//! @authors   Scientific Computing Group at MLZ Garching
-//! @authors   C. Durniak, M. Ganeva, G. Pospelov, W. Van Herck, J. Wuttke
-//
-// ************************************************************************** //
-
-#ifndef vdouble2d_t_hpp__pyplusplus_wrapper
-#define vdouble2d_t_hpp__pyplusplus_wrapper
-
-void register_vdouble2d_t_class();
-
-#endif//vdouble2d_t_hpp__pyplusplus_wrapper
diff --git a/Core/PythonAPI/src/ConstKBinAxis.pypp.cpp b/Core/PythonAPI/src/ConstKBinAxis.pypp.cpp
index 3d0dfb22515c10b652b60564d8b026c71e617e39..fb65d930ee732fc18434d68137608683b48a0037 100644
--- a/Core/PythonAPI/src/ConstKBinAxis.pypp.cpp
+++ b/Core/PythonAPI/src/ConstKBinAxis.pypp.cpp
@@ -119,6 +119,18 @@ struct ConstKBinAxis_wrapper : ConstKBinAxis, bp::wrapper< ConstKBinAxis > {
         return VariableBinAxis::getBinBoundaries( );
     }
 
+    virtual double getBinCenter( ::std::size_t index ) const  {
+        if( bp::override func_getBinCenter = this->get_override( "getBinCenter" ) )
+            return func_getBinCenter( index );
+        else{
+            return this->VariableBinAxis::getBinCenter( index );
+        }
+    }
+    
+    double default_getBinCenter( ::std::size_t index ) const  {
+        return VariableBinAxis::getBinCenter( index );
+    }
+
     virtual ::std::vector< double > getBinCenters(  ) const  {
         if( bp::override func_getBinCenters = this->get_override( "getBinCenters" ) )
             return func_getBinCenters(  );
@@ -260,6 +272,18 @@ void register_ConstKBinAxis_class(){
                 , getBinBoundaries_function_type(&::VariableBinAxis::getBinBoundaries)
                 , default_getBinBoundaries_function_type(&ConstKBinAxis_wrapper::default_getBinBoundaries) );
         
+        }
+        { //::VariableBinAxis::getBinCenter
+        
+            typedef double ( ::VariableBinAxis::*getBinCenter_function_type)( ::std::size_t ) const;
+            typedef double ( ConstKBinAxis_wrapper::*default_getBinCenter_function_type)( ::std::size_t ) const;
+            
+            ConstKBinAxis_exposer.def( 
+                "getBinCenter"
+                , getBinCenter_function_type(&::VariableBinAxis::getBinCenter)
+                , default_getBinCenter_function_type(&ConstKBinAxis_wrapper::default_getBinCenter)
+                , ( bp::arg("index") ) );
+        
         }
         { //::VariableBinAxis::getBinCenters
         
diff --git a/Core/PythonAPI/src/CustomBinAxis.pypp.cpp b/Core/PythonAPI/src/CustomBinAxis.pypp.cpp
index 378004f1ece05a023bb8b1f890c9ecd700dd7ac5..61099c7d5c9dab68d92e6c9a66e4886040423eae 100644
--- a/Core/PythonAPI/src/CustomBinAxis.pypp.cpp
+++ b/Core/PythonAPI/src/CustomBinAxis.pypp.cpp
@@ -131,6 +131,18 @@ struct CustomBinAxis_wrapper : CustomBinAxis, bp::wrapper< CustomBinAxis > {
         return VariableBinAxis::getBinBoundaries( );
     }
 
+    virtual double getBinCenter( ::std::size_t index ) const  {
+        if( bp::override func_getBinCenter = this->get_override( "getBinCenter" ) )
+            return func_getBinCenter( index );
+        else{
+            return this->VariableBinAxis::getBinCenter( index );
+        }
+    }
+    
+    double default_getBinCenter( ::std::size_t index ) const  {
+        return VariableBinAxis::getBinCenter( index );
+    }
+
     virtual double getMax(  ) const  {
         if( bp::override func_getMax = this->get_override( "getMax" ) )
             return func_getMax(  );
@@ -271,6 +283,18 @@ void register_CustomBinAxis_class(){
                 , getBinBoundaries_function_type(&::VariableBinAxis::getBinBoundaries)
                 , default_getBinBoundaries_function_type(&CustomBinAxis_wrapper::default_getBinBoundaries) );
         
+        }
+        { //::VariableBinAxis::getBinCenter
+        
+            typedef double ( ::VariableBinAxis::*getBinCenter_function_type)( ::std::size_t ) const;
+            typedef double ( CustomBinAxis_wrapper::*default_getBinCenter_function_type)( ::std::size_t ) const;
+            
+            CustomBinAxis_exposer.def( 
+                "getBinCenter"
+                , getBinCenter_function_type(&::VariableBinAxis::getBinCenter)
+                , default_getBinCenter_function_type(&CustomBinAxis_wrapper::default_getBinCenter)
+                , ( bp::arg("index") ) );
+        
         }
         { //::VariableBinAxis::getMax
         
diff --git a/Core/PythonAPI/src/FixedBinAxis.pypp.cpp b/Core/PythonAPI/src/FixedBinAxis.pypp.cpp
index b31f6a96e8c3dd790b935e0fa6b684ddd4f42ee1..2875892522ed44450709e0fe3f86bbb88e8c1c34 100644
--- a/Core/PythonAPI/src/FixedBinAxis.pypp.cpp
+++ b/Core/PythonAPI/src/FixedBinAxis.pypp.cpp
@@ -95,6 +95,18 @@ struct FixedBinAxis_wrapper : FixedBinAxis, bp::wrapper< FixedBinAxis > {
         return FixedBinAxis::getBinBoundaries( );
     }
 
+    virtual double getBinCenter( ::std::size_t index ) const  {
+        if( bp::override func_getBinCenter = this->get_override( "getBinCenter" ) )
+            return func_getBinCenter( index );
+        else{
+            return this->FixedBinAxis::getBinCenter( index );
+        }
+    }
+    
+    double default_getBinCenter( ::std::size_t index ) const  {
+        return FixedBinAxis::getBinCenter( index );
+    }
+
     virtual ::std::vector< double > getBinCenters(  ) const  {
         if( bp::override func_getBinCenters = this->get_override( "getBinCenters" ) )
             return func_getBinCenters(  );
@@ -248,6 +260,18 @@ void register_FixedBinAxis_class(){
                 , getBinBoundaries_function_type(&::FixedBinAxis::getBinBoundaries)
                 , default_getBinBoundaries_function_type(&FixedBinAxis_wrapper::default_getBinBoundaries) );
         
+        }
+        { //::FixedBinAxis::getBinCenter
+        
+            typedef double ( ::FixedBinAxis::*getBinCenter_function_type)( ::std::size_t ) const;
+            typedef double ( FixedBinAxis_wrapper::*default_getBinCenter_function_type)( ::std::size_t ) const;
+            
+            FixedBinAxis_exposer.def( 
+                "getBinCenter"
+                , getBinCenter_function_type(&::FixedBinAxis::getBinCenter)
+                , default_getBinCenter_function_type(&FixedBinAxis_wrapper::default_getBinCenter)
+                , ( bp::arg("index") ) );
+        
         }
         { //::FixedBinAxis::getBinCenters
         
diff --git a/Core/PythonAPI/src/Histogram1D.pypp.cpp b/Core/PythonAPI/src/Histogram1D.pypp.cpp
index 2f27936e4e0ca1f260aa1ae37588c14f2a9b0050..0f0275bf015505e5d940977351479009aaa730fd 100644
--- a/Core/PythonAPI/src/Histogram1D.pypp.cpp
+++ b/Core/PythonAPI/src/Histogram1D.pypp.cpp
@@ -92,16 +92,16 @@ struct Histogram1D_wrapper : Histogram1D, bp::wrapper< Histogram1D > {
         return IHistogram::getXaxis( );
     }
 
-    virtual double getXaxisValue( ::std::size_t binGlobalIndex ) {
+    virtual double getXaxisValue( ::std::size_t globalbin ) {
         if( bp::override func_getXaxisValue = this->get_override( "getXaxisValue" ) )
-            return func_getXaxisValue( binGlobalIndex );
+            return func_getXaxisValue( globalbin );
         else{
-            return this->IHistogram::getXaxisValue( binGlobalIndex );
+            return this->IHistogram::getXaxisValue( globalbin );
         }
     }
     
-    double default_getXaxisValue( ::std::size_t binGlobalIndex ) {
-        return IHistogram::getXaxisValue( binGlobalIndex );
+    double default_getXaxisValue( ::std::size_t globalbin ) {
+        return IHistogram::getXaxisValue( globalbin );
     }
 
     virtual double getXmax(  ) const  {
@@ -140,16 +140,16 @@ struct Histogram1D_wrapper : Histogram1D, bp::wrapper< Histogram1D > {
         return IHistogram::getYaxis( );
     }
 
-    virtual double getYaxisValue( ::std::size_t binGlobalIndex ) {
+    virtual double getYaxisValue( ::std::size_t globalbin ) {
         if( bp::override func_getYaxisValue = this->get_override( "getYaxisValue" ) )
-            return func_getYaxisValue( binGlobalIndex );
+            return func_getYaxisValue( globalbin );
         else{
-            return this->IHistogram::getYaxisValue( binGlobalIndex );
+            return this->IHistogram::getYaxisValue( globalbin );
         }
     }
     
-    double default_getYaxisValue( ::std::size_t binGlobalIndex ) {
-        return IHistogram::getYaxisValue( binGlobalIndex );
+    double default_getYaxisValue( ::std::size_t globalbin ) {
+        return IHistogram::getYaxisValue( globalbin );
     }
 
     virtual double getYmax(  ) const  {
@@ -219,6 +219,15 @@ void register_Histogram1D_class(){
                 , getBinCenters_function_type( &::Histogram1D::getBinCenters )
                 , "Increment bin with abscissa x with a weight." );
         
+        }
+        { //::Histogram1D::getBinErrors
+        
+            typedef ::std::vector< double > ( ::Histogram1D::*getBinErrors_function_type)(  ) const;
+            
+            Histogram1D_exposer.def( 
+                "getBinErrors"
+                , getBinErrors_function_type( &::Histogram1D::getBinErrors ) );
+        
         }
         { //::Histogram1D::getBinValues
         
@@ -272,7 +281,7 @@ void register_Histogram1D_class(){
                 "getXaxisValue"
                 , getXaxisValue_function_type(&::IHistogram::getXaxisValue)
                 , default_getXaxisValue_function_type(&Histogram1D_wrapper::default_getXaxisValue)
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") ) );
         
         }
         { //::IHistogram::getXmax
@@ -318,7 +327,7 @@ void register_Histogram1D_class(){
                 "getYaxisValue"
                 , getYaxisValue_function_type(&::IHistogram::getYaxisValue)
                 , default_getYaxisValue_function_type(&Histogram1D_wrapper::default_getYaxisValue)
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") ) );
         
         }
         { //::IHistogram::getYmax
diff --git a/Core/PythonAPI/src/Histogram2D.pypp.cpp b/Core/PythonAPI/src/Histogram2D.pypp.cpp
index 202a93948849940308fe46132b0a7abd39a07b20..6474db301cb9ff9de65e75bfd6d5c12a754df8dc 100644
--- a/Core/PythonAPI/src/Histogram2D.pypp.cpp
+++ b/Core/PythonAPI/src/Histogram2D.pypp.cpp
@@ -92,16 +92,16 @@ struct Histogram2D_wrapper : Histogram2D, bp::wrapper< Histogram2D > {
         return IHistogram::getXaxis( );
     }
 
-    virtual double getXaxisValue( ::std::size_t binGlobalIndex ) {
+    virtual double getXaxisValue( ::std::size_t globalbin ) {
         if( bp::override func_getXaxisValue = this->get_override( "getXaxisValue" ) )
-            return func_getXaxisValue( binGlobalIndex );
+            return func_getXaxisValue( globalbin );
         else{
-            return this->IHistogram::getXaxisValue( binGlobalIndex );
+            return this->IHistogram::getXaxisValue( globalbin );
         }
     }
     
-    double default_getXaxisValue( ::std::size_t binGlobalIndex ) {
-        return IHistogram::getXaxisValue( binGlobalIndex );
+    double default_getXaxisValue( ::std::size_t globalbin ) {
+        return IHistogram::getXaxisValue( globalbin );
     }
 
     virtual double getXmax(  ) const  {
@@ -140,16 +140,16 @@ struct Histogram2D_wrapper : Histogram2D, bp::wrapper< Histogram2D > {
         return IHistogram::getYaxis( );
     }
 
-    virtual double getYaxisValue( ::std::size_t binGlobalIndex ) {
+    virtual double getYaxisValue( ::std::size_t globalbin ) {
         if( bp::override func_getYaxisValue = this->get_override( "getYaxisValue" ) )
-            return func_getYaxisValue( binGlobalIndex );
+            return func_getYaxisValue( globalbin );
         else{
-            return this->IHistogram::getYaxisValue( binGlobalIndex );
+            return this->IHistogram::getYaxisValue( globalbin );
         }
     }
     
-    double default_getYaxisValue( ::std::size_t binGlobalIndex ) {
-        return IHistogram::getYaxisValue( binGlobalIndex );
+    double default_getYaxisValue( ::std::size_t globalbin ) {
+        return IHistogram::getYaxisValue( globalbin );
     }
 
     virtual double getYmax(  ) const  {
@@ -199,6 +199,18 @@ void register_Histogram2D_class(){
         Histogram2D_exposer.def( bp::init< int, std::vector< double > const &, int, std::vector< double > const & >(( bp::arg("nbinsx"), bp::arg("xbins"), bp::arg("nbinsy"), bp::arg("ybins") ), "Constructor for variable bin size histograms. @param nbinsx number of bins on X-axis @param xbins Array of size nbins+1 containing low-edges for each bin and upper edge of last bin. @param nbinsy number of bins on Y-axis @param ybins Array of size nbins+1 containing low-edges for each bin and upper edge of last bin. \n\n:Parameters:\n  - 'nbinsx' - number of bins on X-axis\n  - 'xbins' - Array of size nbins+1 containing low-edges for each\n  - 'nbinsy' - number of bins on Y-axis\n  - 'ybins' - Array of size nbins+1 containing low-edges for each\n") );
         Histogram2D_exposer.def( bp::init< IAxis const &, IAxis const & >(( bp::arg("axis_x"), bp::arg("axis_y") ), "Constructor for 2D histogram with custom axes.") );
         Histogram2D_exposer.def( bp::init< OutputData< double > const & >(( bp::arg("data") ), "Constructor for 2D histograms from basic OutputData object.") );
+        { //::Histogram2D::crop
+        
+            typedef ::Histogram2D * ( ::Histogram2D::*crop_function_type)( double,double,double,double ) ;
+            
+            Histogram2D_exposer.def( 
+                "crop"
+                , crop_function_type( &::Histogram2D::crop )
+                , ( bp::arg("xmin"), bp::arg("ymin"), bp::arg("xmax"), bp::arg("ymax") )
+                , bp::return_value_policy< bp::reference_existing_object >()
+                , "Create new histogram by applying rectangular clip." );
+        
+        }
         { //::Histogram2D::fill
         
             typedef int ( ::Histogram2D::*fill_function_type)( double,double,double ) ;
@@ -209,15 +221,6 @@ void register_Histogram2D_class(){
                 , ( bp::arg("x"), bp::arg("y"), bp::arg("weight")=1.0e+0 )
                 , "Increment bin with abscissa x and ordinate y with a weight." );
         
-        }
-        { //::Histogram2D::getData
-        
-            typedef ::vdouble2d_t ( ::Histogram2D::*getData_function_type)(  ) const;
-            
-            Histogram2D_exposer.def( 
-                "getData"
-                , getData_function_type( &::Histogram2D::getData ) );
-        
         }
         { //::Histogram2D::getRank
         
@@ -334,7 +337,7 @@ void register_Histogram2D_class(){
                 "getXaxisValue"
                 , getXaxisValue_function_type(&::IHistogram::getXaxisValue)
                 , default_getXaxisValue_function_type(&Histogram2D_wrapper::default_getXaxisValue)
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") ) );
         
         }
         { //::IHistogram::getXmax
@@ -380,7 +383,7 @@ void register_Histogram2D_class(){
                 "getYaxisValue"
                 , getYaxisValue_function_type(&::IHistogram::getYaxisValue)
                 , default_getYaxisValue_function_type(&Histogram2D_wrapper::default_getYaxisValue)
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") ) );
         
         }
         { //::IHistogram::getYmax
diff --git a/Core/PythonAPI/src/IAxis.pypp.cpp b/Core/PythonAPI/src/IAxis.pypp.cpp
index dbe88a03159b7ad2117f7a801a8533e4596f671e..1796240a109b5ca1db77f62c7a107c5e52bb3074 100644
--- a/Core/PythonAPI/src/IAxis.pypp.cpp
+++ b/Core/PythonAPI/src/IAxis.pypp.cpp
@@ -98,6 +98,11 @@ struct IAxis_wrapper : IAxis, bp::wrapper< IAxis > {
         return IAxis::getBinBoundaries( );
     }
 
+    virtual double getBinCenter( ::std::size_t index ) const {
+        bp::override func_getBinCenter = this->get_override( "getBinCenter" );
+        return func_getBinCenter( index );
+    }
+
     virtual ::std::vector< double > getBinCenters(  ) const  {
         if( bp::override func_getBinCenters = this->get_override( "getBinCenters" ) )
             return func_getBinCenters(  );
@@ -225,6 +230,17 @@ void register_IAxis_class(){
                 , getBinBoundaries_function_type(&::IAxis::getBinBoundaries)
                 , default_getBinBoundaries_function_type(&IAxis_wrapper::default_getBinBoundaries) );
         
+        }
+        { //::IAxis::getBinCenter
+        
+            typedef double ( ::IAxis::*getBinCenter_function_type)( ::std::size_t ) const;
+            
+            IAxis_exposer.def( 
+                "getBinCenter"
+                , bp::pure_virtual( getBinCenter_function_type(&::IAxis::getBinCenter) )
+                , ( bp::arg("index") )
+                , "Returns value of last point of axis." );
+        
         }
         { //::IAxis::getBinCenters
         
diff --git a/Core/PythonAPI/src/IHistogram.pypp.cpp b/Core/PythonAPI/src/IHistogram.pypp.cpp
index ae9721951cb47420f69d4db414db1c3c3546ab3a..779870c105a0d310125c38e3c3d8cf1869426abd 100644
--- a/Core/PythonAPI/src/IHistogram.pypp.cpp
+++ b/Core/PythonAPI/src/IHistogram.pypp.cpp
@@ -92,16 +92,16 @@ struct IHistogram_wrapper : IHistogram, bp::wrapper< IHistogram > {
         return IHistogram::getXaxis( );
     }
 
-    virtual double getXaxisValue( ::std::size_t binGlobalIndex ) {
+    virtual double getXaxisValue( ::std::size_t globalbin ) {
         if( bp::override func_getXaxisValue = this->get_override( "getXaxisValue" ) )
-            return func_getXaxisValue( binGlobalIndex );
+            return func_getXaxisValue( globalbin );
         else{
-            return this->IHistogram::getXaxisValue( binGlobalIndex );
+            return this->IHistogram::getXaxisValue( globalbin );
         }
     }
     
-    double default_getXaxisValue( ::std::size_t binGlobalIndex ) {
-        return IHistogram::getXaxisValue( binGlobalIndex );
+    double default_getXaxisValue( ::std::size_t globalbin ) {
+        return IHistogram::getXaxisValue( globalbin );
     }
 
     virtual double getXmax(  ) const  {
@@ -140,16 +140,16 @@ struct IHistogram_wrapper : IHistogram, bp::wrapper< IHistogram > {
         return IHistogram::getYaxis( );
     }
 
-    virtual double getYaxisValue( ::std::size_t binGlobalIndex ) {
+    virtual double getYaxisValue( ::std::size_t globalbin ) {
         if( bp::override func_getYaxisValue = this->get_override( "getYaxisValue" ) )
-            return func_getYaxisValue( binGlobalIndex );
+            return func_getYaxisValue( globalbin );
         else{
-            return this->IHistogram::getYaxisValue( binGlobalIndex );
+            return this->IHistogram::getYaxisValue( globalbin );
         }
     }
     
-    double default_getYaxisValue( ::std::size_t binGlobalIndex ) {
-        return IHistogram::getYaxisValue( binGlobalIndex );
+    double default_getYaxisValue( ::std::size_t globalbin ) {
+        return IHistogram::getYaxisValue( globalbin );
     }
 
     virtual double getYmax(  ) const  {
@@ -213,15 +213,81 @@ void register_IHistogram_class(){
                 , getArray_function_type( &::IHistogram::getArray ) );
         
         }
-        { //::IHistogram::getBinValue
+        { //::IHistogram::getBinContent
         
-            typedef double ( ::IHistogram::*getBinValue_function_type)( ::std::size_t ) const;
+            typedef double ( ::IHistogram::*getBinContent_function_type)( int ) const;
             
             IHistogram_exposer.def( 
-                "getBinValue"
-                , getBinValue_function_type( &::IHistogram::getBinValue )
-                , ( bp::arg("binGlobalIndex") )
-                , "Returns total number of histogram bins. For 2D histograms the result will be the product of bin numbers along X and Y axes. " );
+                "getBinContent"
+                , getBinContent_function_type( &::IHistogram::getBinContent )
+                , ( bp::arg("bin") )
+                , "Returns content of the bin with given index. For 1D histograms bin index is related to x-axis. For 2D histograms bin index is global bin index. @param bin Bin index @return The content of the bin (which is normally the value accumulated by the bin) \n\n:Parameters:\n  - 'bin' - Bin index\n" );
+        
+        }
+        { //::IHistogram::getBinContent
+        
+            typedef double ( ::IHistogram::*getBinContent_function_type)( int,int ) const;
+            
+            IHistogram_exposer.def( 
+                "getBinContent"
+                , getBinContent_function_type( &::IHistogram::getBinContent )
+                , ( bp::arg("binx"), bp::arg("biny") )
+                , "Returns content of the bin of 2D histogram with given axes indices. @param binx X-axis bin index @param biny Y-axis bin index @return The content of the bin (which is normally the value accumulated by the bin) \n\n:Parameters:\n  - 'binx' - X-axis bin index\n  - 'biny' - Y-axis bin index\n" );
+        
+        }
+        { //::IHistogram::getBinError
+        
+            typedef double ( ::IHistogram::*getBinError_function_type)( int ) const;
+            
+            IHistogram_exposer.def( 
+                "getBinError"
+                , getBinError_function_type( &::IHistogram::getBinError )
+                , ( bp::arg("bin") )
+                , "Returns error of the bin with given index." );
+        
+        }
+        { //::IHistogram::getBinError
+        
+            typedef double ( ::IHistogram::*getBinError_function_type)( int,int ) const;
+            
+            IHistogram_exposer.def( 
+                "getBinError"
+                , getBinError_function_type( &::IHistogram::getBinError )
+                , ( bp::arg("binx"), bp::arg("biny") )
+                , "Returns error of the bin of 2D histogram with given axes indices." );
+        
+        }
+        { //::IHistogram::getBinNumberOfEntries
+        
+            typedef int ( ::IHistogram::*getBinNumberOfEntries_function_type)( int ) const;
+            
+            IHistogram_exposer.def( 
+                "getBinNumberOfEntries"
+                , getBinNumberOfEntries_function_type( &::IHistogram::getBinNumberOfEntries )
+                , ( bp::arg("bin") )
+                , "Returns number of entries in the bin with given index." );
+        
+        }
+        { //::IHistogram::getBinNumberOfEntries
+        
+            typedef int ( ::IHistogram::*getBinNumberOfEntries_function_type)( int,int ) const;
+            
+            IHistogram_exposer.def( 
+                "getBinNumberOfEntries"
+                , getBinNumberOfEntries_function_type( &::IHistogram::getBinNumberOfEntries )
+                , ( bp::arg("binx"), bp::arg("biny") )
+                , "Returns number of entries in the bin of 2D histogram with given axes indices." );
+        
+        }
+        { //::IHistogram::getGlobalBin
+        
+            typedef int ( ::IHistogram::*getGlobalBin_function_type)( int,int ) const;
+            
+            IHistogram_exposer.def( 
+                "getGlobalBin"
+                , getGlobalBin_function_type( &::IHistogram::getGlobalBin )
+                , ( bp::arg("binx"), bp::arg("biny")=(int)(0) )
+                , "Returns global bin index for given axes indices. For 1D histogram the global bin index coinside with axis index. @param binx X-axis bin index @param biny Y-axis bin index @return The global bin index \n\n:Parameters:\n  - 'binx' - X-axis bin index\n  - 'biny' - Y-axis bin index\n" );
         
         }
         { //::IHistogram::getRank
@@ -265,7 +331,8 @@ void register_IHistogram_class(){
             IHistogram_exposer.def( 
                 "getXaxisIndex"
                 , getXaxisIndex_function_type( &::IHistogram::getXaxisIndex )
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") )
+                , "Returns x-axis bin index for given globalbin. For 1D histograms returned value conicide with globalbin value " );
         
         }
         { //::IHistogram::getXaxisValue
@@ -277,7 +344,7 @@ void register_IHistogram_class(){
                 "getXaxisValue"
                 , getXaxisValue_function_type(&::IHistogram::getXaxisValue)
                 , default_getXaxisValue_function_type(&IHistogram_wrapper::default_getXaxisValue)
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") ) );
         
         }
         { //::IHistogram::getXmax
@@ -321,7 +388,8 @@ void register_IHistogram_class(){
             IHistogram_exposer.def( 
                 "getYaxisIndex"
                 , getYaxisIndex_function_type( &::IHistogram::getYaxisIndex )
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") )
+                , "Returns x-axis bin index for given globalbin. For 1D histograms returned value conicide with globalbin value " );
         
         }
         { //::IHistogram::getYaxisValue
@@ -333,7 +401,7 @@ void register_IHistogram_class(){
                 "getYaxisValue"
                 , getYaxisValue_function_type(&::IHistogram::getYaxisValue)
                 , default_getYaxisValue_function_type(&IHistogram_wrapper::default_getYaxisValue)
-                , ( bp::arg("binGlobalIndex") ) );
+                , ( bp::arg("globalbin") ) );
         
         }
         { //::IHistogram::getYmax
diff --git a/Core/PythonAPI/src/PythonModule.cpp b/Core/PythonAPI/src/PythonModule.cpp
index c13181b567e0d60b1860786a620596d12f7b1e81..37a0f89c86b30662b6f25e8f725a33171079d74d 100644
--- a/Core/PythonAPI/src/PythonModule.cpp
+++ b/Core/PythonAPI/src/PythonModule.cpp
@@ -24,7 +24,6 @@ GCC_DIAG_ON(missing-field-initializers)
 #include "InterferenceFunctionNone.pypp.h"
 #include "FTDistribution2DGate.pypp.h"
 #include "FormFactorPrism6.pypp.h"
-#include "vdouble2d_t.pypp.h"
 #include "FormFactorTruncatedCube.pypp.h"
 #include "IFTDistribution2D.pypp.h"
 #include "FormFactorLorentz.pypp.h"
@@ -157,7 +156,6 @@ BOOST_PYTHON_MODULE(libBornAgainCore){
     boost::python::docstring_options doc_options(true, true, false);
 
     register_vector_longinteger_t_class();
-    register_vdouble2d_t_class();
     register_vector_string_t_class();
     register_vector_complex_t_class();
     register_vector_integer_t_class();
diff --git a/Core/PythonAPI/src/VariableBinAxis.pypp.cpp b/Core/PythonAPI/src/VariableBinAxis.pypp.cpp
index 42fe70314fa03a06ff3939d7da5f0acdb04b0ae4..2c516db899d2e2e2f2f4693c80cbe0f51781d86c 100644
--- a/Core/PythonAPI/src/VariableBinAxis.pypp.cpp
+++ b/Core/PythonAPI/src/VariableBinAxis.pypp.cpp
@@ -95,6 +95,18 @@ struct VariableBinAxis_wrapper : VariableBinAxis, bp::wrapper< VariableBinAxis >
         return VariableBinAxis::getBinBoundaries( );
     }
 
+    virtual double getBinCenter( ::std::size_t index ) const  {
+        if( bp::override func_getBinCenter = this->get_override( "getBinCenter" ) )
+            return func_getBinCenter( index );
+        else{
+            return this->VariableBinAxis::getBinCenter( index );
+        }
+    }
+    
+    double default_getBinCenter( ::std::size_t index ) const  {
+        return VariableBinAxis::getBinCenter( index );
+    }
+
     virtual ::std::vector< double > getBinCenters(  ) const  {
         if( bp::override func_getBinCenters = this->get_override( "getBinCenters" ) )
             return func_getBinCenters(  );
@@ -248,6 +260,18 @@ void register_VariableBinAxis_class(){
                 , getBinBoundaries_function_type(&::VariableBinAxis::getBinBoundaries)
                 , default_getBinBoundaries_function_type(&VariableBinAxis_wrapper::default_getBinBoundaries) );
         
+        }
+        { //::VariableBinAxis::getBinCenter
+        
+            typedef double ( ::VariableBinAxis::*getBinCenter_function_type)( ::std::size_t ) const;
+            typedef double ( VariableBinAxis_wrapper::*default_getBinCenter_function_type)( ::std::size_t ) const;
+            
+            VariableBinAxis_exposer.def( 
+                "getBinCenter"
+                , getBinCenter_function_type(&::VariableBinAxis::getBinCenter)
+                , default_getBinCenter_function_type(&VariableBinAxis_wrapper::default_getBinCenter)
+                , ( bp::arg("index") ) );
+        
         }
         { //::VariableBinAxis::getBinCenters
         
diff --git a/Core/PythonAPI/src/vdouble2d_t.pypp.cpp b/Core/PythonAPI/src/vdouble2d_t.pypp.cpp
deleted file mode 100644
index e400d5a1cf26f14e7579bc9db8a6fcdec87c9889..0000000000000000000000000000000000000000
--- a/Core/PythonAPI/src/vdouble2d_t.pypp.cpp
+++ /dev/null
@@ -1,40 +0,0 @@
-// This file has been generated by Py++.
-
-// ************************************************************************** //
-//
-//  BornAgain: simulate and fit scattering at grazing incidence
-//
-//! @file      Automatically generated boost::python code for BornAgain Python bindings
-//! @brief     Automatically generated boost::python code for BornAgain Python bindings
-//!
-//! @homepage  http://bornagainproject.org
-//! @license   GNU General Public License v3 or higher (see COPYING)
-//! @copyright Forschungszentrum Juelich GmbH 2015
-//! @authors   Scientific Computing Group at MLZ Garching
-//! @authors   C. Durniak, M. Ganeva, G. Pospelov, W. Van Herck, J. Wuttke
-//
-// ************************************************************************** //
-
-#include "Macros.h"
-GCC_DIAG_OFF(unused-parameter)
-GCC_DIAG_OFF(missing-field-initializers)
-#include "boost/python.hpp"
-#include "boost/python/suite/indexing/vector_indexing_suite.hpp"
-GCC_DIAG_ON(unused-parameter)
-GCC_DIAG_ON(missing-field-initializers)
-#include "PythonCoreList.h"
-#include "vdouble2d_t.pypp.h"
-
-namespace bp = boost::python;
-
-void register_vdouble2d_t_class(){
-
-    { //::std::vector< std::vector< double > >
-        typedef bp::class_< std::vector< std::vector< double > > > vdouble2d_t_exposer_t;
-        vdouble2d_t_exposer_t vdouble2d_t_exposer = vdouble2d_t_exposer_t( "vdouble2d_t" );
-        bp::scope vdouble2d_t_scope( vdouble2d_t_exposer );
-        //WARNING: the next line of code will not compile, because "::std::vector<double, std::allocator<double> >" does not have operator== !
-        vdouble2d_t_exposer.def( bp::vector_indexing_suite< ::std::vector< std::vector< double > > >() );
-    }
-
-}
diff --git a/Core/Tools/inc/FixedBinAxis.h b/Core/Tools/inc/FixedBinAxis.h
index f7ed817702778b2957ac18e5f00ae3b57587b4d6..7b356ff40ff39b7b9e1dfea2b32687631a47ce04 100644
--- a/Core/Tools/inc/FixedBinAxis.h
+++ b/Core/Tools/inc/FixedBinAxis.h
@@ -46,6 +46,8 @@ public:
 
     double getMax() const;
 
+    double getBinCenter(size_t index) const;
+
     size_t findClosestIndex(double value) const;
 
     std::vector<double > getBinCenters() const;
diff --git a/Core/Tools/inc/Histogram1D.h b/Core/Tools/inc/Histogram1D.h
index 6985e22a99bd434f22fc5636032e82b7f37b04c0..5ddd8a358d6faf90021b8dd0ce821b6a0fad0e0c 100644
--- a/Core/Tools/inc/Histogram1D.h
+++ b/Core/Tools/inc/Histogram1D.h
@@ -51,6 +51,9 @@ public:
     std::vector<double> getBinCenters() const;
 
     std::vector<double> getBinValues() const;
+
+    std::vector<double> getBinErrors() const;
+
 };
 
 
diff --git a/Core/Tools/inc/Histogram2D.h b/Core/Tools/inc/Histogram2D.h
index 33d533fe69ae4896cb72d9feb2035fcdb39124b1..724c3c4ecd72fbaed702e5590462b60f5af22ba4 100644
--- a/Core/Tools/inc/Histogram2D.h
+++ b/Core/Tools/inc/Histogram2D.h
@@ -102,30 +102,13 @@ public:
     Histogram1D *projectionY(double xlow, double xup, ProjectionType projectionType = INTEGRAL);
 
 
+    //! Create new histogram by applying rectangular clip.
+    Histogram2D *crop(double xmin, double ymin, double xmax, double ymax);
 
-    vdouble2d_t getData() const;
 
-//    //! Project a 2D histogram into 1D histogram along Y.
-//    Histogram1D *projectionY(ProjectionType = INTEGRAL);
-
-//    //! Project a 2D histogram into 1D histogram along X, the projection
-//    //! is made from the channels along the Y axis at 'yvalue'.
-//    Histogram1D *profileX(double yvalue);
+//    vdouble2d_t getData() const;
 
-//    //! Project a 2D histogram into 1D histogram along X, the projection
-//    //! is made from the channels along the Y axis ranging from ymin to ymax.
-//    Histogram1D *profileX(double ymin, double ymax, ProjectionType = INTEGRAL);
 
-//    //! Project a 2D histogram into 1D histogram along X, the projection
-//    //! is made from the channels along the Y axis at 'yvalue'.
-//    Histogram1D *profileY(double yvalue);
-
-//    //! Project a 2D histogram into 1D histogram along X, the projection
-//    //! is made from the channels along the Y axis ranging from ymin to ymax.
-//    Histogram1D *profileX(double ymin, double ymax, ProjectionType projectionType = INTEGRAL);
-
-//    //! Create new intensity data by applying rectangular clip.
-//    Histogram1D *createClipped(double xmin, double ymin, double xmax, double ymax);
 
 //    //!
 //    //! The function returns the corresponding global bin number which has its content
diff --git a/Core/Tools/inc/IAxis.h b/Core/Tools/inc/IAxis.h
index 7701e546e156ad55952012ed78c869afc8ed5a77..f3fa56d68c415e23f4868c273b6a40c2deece7a8 100644
--- a/Core/Tools/inc/IAxis.h
+++ b/Core/Tools/inc/IAxis.h
@@ -58,6 +58,8 @@ public:
     //! Returns value of last point of axis
     virtual double getMax() const=0;
 
+    virtual double getBinCenter(size_t index) const=0;
+
     //! find bin index which is best match for given value
     virtual size_t findClosestIndex(double value) const=0;
 
diff --git a/Core/Tools/inc/IHistogram.h b/Core/Tools/inc/IHistogram.h
index 8fabc62bdf05c552d4b8ce5e131b2126573eaf06..99383bcdd969f49cffd429cdc9d71a1e4ecbdf6c 100644
--- a/Core/Tools/inc/IHistogram.h
+++ b/Core/Tools/inc/IHistogram.h
@@ -48,40 +48,85 @@ public:
 
     //! Returns total number of histogram bins. For 2D histograms the result will be the product
     //! of bin numbers along X and Y axes.
-    virtual size_t getTotalNumberOfBins() const;
-
-    double getBinValue(size_t binGlobalIndex) const;
+    size_t getTotalNumberOfBins() const;
 
     //! returns x-axis
-    virtual const IAxis *getXaxis() const;
+    const IAxis *getXaxis() const;
+
+    //! returns y-axis for 2D histograms
+    const IAxis *getYaxis() const;
+
+    //! Returns x-axis min (lower edge of first bin).
+    double getXmin() const;
+
+    //! Returns x-axis max (upper edge of last bin).
+    double getXmax() const;
+
+    //! Returns y-axis min (lower edge of first bin) for 2D histograms.
+    double getYmin() const;
 
-    //! returns y-axis (throws an exception for 1D histograms)
-    virtual const IAxis *getYaxis() const;
+    //! Returns y-axis max (upper edge of last bin) for 2D histograms.
+    double getYmax() const;
 
+    //! Returns global bin index for given axes indices. For 1D histogram the global bin
+    //! index coinside with axis index.
+    //! @param binx X-axis bin index
+    //! @param biny Y-axis bin index
+    //! @return The global bin index
+    int getGlobalBin(int binx, int biny = 0) const;
 
-    virtual double getXmin() const;
-    virtual double getXmax() const;
-    virtual double getYmin() const;
-    virtual double getYmax() const;
+    //! Returns x-axis bin index for given globalbin. For 1D histograms returned value conicide
+    //! with globalbin value
+    int getXaxisIndex(size_t globalbin) const;
 
-    int getXaxisIndex(size_t binGlobalIndex) const;
-    int getYaxisIndex(size_t binGlobalIndex) const;
+
+    //! Returns x-axis bin index for given globalbin. For 1D histograms returned value conicide
+    //! with globalbin value
+    int getYaxisIndex(size_t globalbin) const;
 
     //! Returns the value on x-axis corresponding to the global bin index.
     //! @param binGlobalIndex The global bin index
-    virtual double getXaxisValue(size_t binGlobalIndex);
+    //! @return The center of axis's corresponding bin
+    double getXaxisValue(size_t globalbin);
 
-    //! Returns the value on y-axis corresponding to the global bin index.
+    //! Returns the value on y-axis corresponding to the global bin index (for 2D histograms).
     //! @param globalbin The global bin index
-    virtual double getYaxisValue(size_t binGlobalIndex);
+    //! @return The center of axis's corresponding bin
+    double getYaxisValue(size_t globalbin);
+
+
+    //! Returns content of the bin with given index. For 1D histograms bin index is related
+    //! to x-axis. For 2D histograms bin index is global bin index.
+    //! @param bin Bin index
+    //! @return The content of the bin (which is normally the value accumulated by the bin)
+    double getBinContent(int bin) const;
+
+    //! Returns content of the bin of 2D histogram with given axes indices.
+    //! @param binx X-axis bin index
+    //! @param biny Y-axis bin index
+    //! @return The content of the bin (which is normally the value accumulated by the bin)
+    double getBinContent(int binx, int biny) const;
+
+    //! Returns error of the bin with given index.
+    double getBinError(int bin) const;
+
+    //! Returns error of the bin of 2D histogram with given axes indices.
+    double getBinError(int binx, int biny) const;
+
+    //! Returns number of entries in the bin with given index.
+    int getBinNumberOfEntries(int bin) const;
+
+    //! Returns number of entries in the bin of 2D histogram with given axes indices.
+    int getBinNumberOfEntries(int binx, int biny) const;
 
-    //! Reset histogram content (axes remains)
-    virtual void reset();
 
 #ifdef BORNAGAIN_PYTHON
     PyObject *getArray() const;
 #endif
 
+    //! Reset histogram content (axes remains)
+    virtual void reset();
+
 protected:
     void check_x_axis() const;
     void check_y_axis() const;
diff --git a/Core/Tools/inc/VariableBinAxis.h b/Core/Tools/inc/VariableBinAxis.h
index 6548c05958789dd80ecea6263963bd8e040e0828..722a4334f2e61901412c04c7942b697dd580942b 100644
--- a/Core/Tools/inc/VariableBinAxis.h
+++ b/Core/Tools/inc/VariableBinAxis.h
@@ -46,6 +46,8 @@ public:
 
     double getMax() const;
 
+    double getBinCenter(size_t index) const;
+
     size_t findClosestIndex(double value) const;
 
     std::vector<double > getBinCenters() const;
diff --git a/Core/Tools/src/FixedBinAxis.cpp b/Core/Tools/src/FixedBinAxis.cpp
index fe45fd0fc66027c996304db315c3c1957c9f09fd..eea317cc5757a43c603447b9684f43877293d0a3 100644
--- a/Core/Tools/src/FixedBinAxis.cpp
+++ b/Core/Tools/src/FixedBinAxis.cpp
@@ -73,6 +73,11 @@ double FixedBinAxis::getMax() const
     return m_end;
 }
 
+double FixedBinAxis::getBinCenter(size_t index) const
+{
+    return (*this)[index];
+}
+
 
 size_t FixedBinAxis::findClosestIndex(double value) const
 {
diff --git a/Core/Tools/src/Histogram1D.cpp b/Core/Tools/src/Histogram1D.cpp
index 0c47fd54293838f7aebcbb8ebc6139fdd6e2464a..9235923a65fc0b1590228e27bd36fea1cea8cfac 100644
--- a/Core/Tools/src/Histogram1D.cpp
+++ b/Core/Tools/src/Histogram1D.cpp
@@ -66,9 +66,19 @@ std::vector<double> Histogram1D::getBinCenters() const
 std::vector<double> Histogram1D::getBinValues() const
 {
     std::vector<double> result;
-    result.resize(getXaxis()->getSize(), 0.0);
-    for(size_t index = 0; index<m_data.getAllocatedSize(); ++index) {
-        result[index] = m_data[index].getValue();
+    result.resize(getTotalNumberOfBins(), 0.0);
+    for(size_t index=0; index<getTotalNumberOfBins(); ++index) {
+        result[index] = getBinContent(index);
+    }
+    return result;
+}
+
+std::vector<double> Histogram1D::getBinErrors() const
+{
+    std::vector<double> result;
+    result.resize(getTotalNumberOfBins(), 0.0);
+    for(size_t index=0; index<getTotalNumberOfBins(); ++index) {
+        result[index] = getBinError(index);
     }
     return result;
 }
diff --git a/Core/Tools/src/Histogram2D.cpp b/Core/Tools/src/Histogram2D.cpp
index d3fc98c321375dc47ee5b5146f5570686858ff37..8ef5b7321d5851cbafdfe29621e37a4aede3da4e 100644
--- a/Core/Tools/src/Histogram2D.cpp
+++ b/Core/Tools/src/Histogram2D.cpp
@@ -17,7 +17,7 @@
 #include "Histogram1D.h"
 #include "FixedBinAxis.h"
 #include "VariableBinAxis.h"
-
+#include <boost/scoped_ptr.hpp>
 
 Histogram2D::Histogram2D(int nbinsx, double xlow, double xup, int nbinsy, double ylow, double yup)
 {
@@ -99,30 +99,52 @@ Histogram1D *Histogram2D::projectionY(double xlow, double xup, IHistogram::Proje
     return create_projectionY(xbinlow, xbinup, projectionType);
 }
 
-
-vdouble2d_t Histogram2D::getData() const
+Histogram2D *Histogram2D::crop(double xmin, double ymin, double xmax, double ymax)
 {
-    vdouble2d_t result;
-    result.resize(getXaxis()->getSize());
-
-    vdouble1d_t column;
-    column.resize(getYaxis()->getSize());
-
-    for(size_t i=0; i< result.size(); ++i) {
-        result[i] = column;
-    }
-
-    size_t index = 0;
-    for (size_t ix=0; ix<getXaxis()->getSize(); ++ix) {
-        for(size_t iy=0; iy<getYaxis()->getSize(); ++iy) {
-            result[ix][iy] = m_data[index].getValue();
-            ++index;
+    boost::scoped_ptr<IAxis > xaxis(getXaxis()->createClippedAxis(xmin, xmax));
+    boost::scoped_ptr<IAxis > yaxis(getYaxis()->createClippedAxis(ymin, ymax));
+
+    Histogram2D *result = new Histogram2D(*xaxis, *yaxis);
+    OutputData<CumulativeValue>::const_iterator it_origin = m_data.begin();
+    OutputData<CumulativeValue>::iterator it_result = result->m_data.begin();
+    while (it_origin != m_data.end())
+    {
+        double x = m_data.getAxisValue(it_origin.getIndex(), 0);
+        double y = m_data.getAxisValue(it_origin.getIndex(), 1);
+        if(result->getXaxis()->contains(x) && result->getYaxis()->contains(y)) {
+            *it_result = *it_origin;
+            ++it_result;
         }
+        ++it_origin;
     }
-
     return result;
 }
 
+
+
+//vdouble2d_t Histogram2D::getData() const
+//{
+//    vdouble2d_t result;
+//    result.resize(getXaxis()->getSize());
+
+//    vdouble1d_t column;
+//    column.resize(getYaxis()->getSize());
+
+//    for(size_t i=0; i< result.size(); ++i) {
+//        result[i] = column;
+//    }
+
+//    size_t index = 0;
+//    for (size_t ix=0; ix<getXaxis()->getSize(); ++ix) {
+//        for(size_t iy=0; iy<getYaxis()->getSize(); ++iy) {
+//            result[ix][iy] = m_data[index].getValue();
+//            ++index;
+//        }
+//    }
+
+//    return result;
+//}
+
 Histogram1D *Histogram2D::create_projectionX(int ybinlow, int ybinup, IHistogram::ProjectionType projectionType)
 {
     Histogram1D *result = new Histogram1D(*this->getXaxis());
@@ -132,7 +154,7 @@ Histogram1D *Histogram2D::create_projectionX(int ybinlow, int ybinup, IHistogram
         int ybin = getYaxisIndex(index);
 
         if(ybin >= ybinlow && ybin <= ybinup) {
-            result->fill(getXaxisValue(index), getBinValue(index));
+            result->fill(getXaxisValue(index), getBinContent(index));
         }
     }
     return result;
@@ -147,7 +169,7 @@ Histogram1D *Histogram2D::create_projectionY(int xbinlow, int xbinup, IHistogram
         int xbin = getXaxisIndex(index);
 
         if(xbin >= xbinlow && xbin <= xbinup) {
-            result->fill(getYaxisValue(index), getBinValue(index));
+            result->fill(getYaxisValue(index), getBinContent(index));
         }
     }
     return result;
diff --git a/Core/Tools/src/IHistogram.cpp b/Core/Tools/src/IHistogram.cpp
index 8abf98d22e8318ff828963d772f1c77e1daee097..cf0f4d0b4640ce32f8683e408a85ae6ebb95778e 100644
--- a/Core/Tools/src/IHistogram.cpp
+++ b/Core/Tools/src/IHistogram.cpp
@@ -18,6 +18,7 @@
 #include "VariableBinAxis.h"
 #include "Exceptions.h"
 #include <sstream>
+#include <boost/assign/list_of.hpp>
 
 IHistogram::IHistogram(const IAxis &axis_x)
 {
@@ -46,10 +47,10 @@ size_t IHistogram::getTotalNumberOfBins() const
     return m_data.getAllocatedSize();
 }
 
-double IHistogram::getBinValue(size_t binGlobalIndex) const
-{
-    return m_data[binGlobalIndex].getValue();
-}
+//double IHistogram::getBinValue(size_t binGlobalIndex) const
+//{
+//    return m_data[binGlobalIndex].getValue();
+//}
 
 const IAxis *IHistogram::getXaxis() const
 {
@@ -83,31 +84,64 @@ double IHistogram::getYmax() const
     return getYaxis()->getMax();
 }
 
-int IHistogram::getXaxisIndex(size_t binGlobalIndex) const
+int IHistogram::getGlobalBin(int binx, int biny) const
 {
-    return m_data.getAxisBinIndex(binGlobalIndex, 0);
+    std::vector<int > axes_indices;
+    axes_indices.push_back(binx);
+    if(getRank() == 2) axes_indices.push_back(biny);
+    return m_data.toGlobalIndex(axes_indices);
 }
 
-int IHistogram::getYaxisIndex(size_t binGlobalIndex) const
+int IHistogram::getXaxisIndex(size_t globalbin) const
 {
-    return m_data.getAxisBinIndex(binGlobalIndex, 1);
+    return m_data.getAxisBinIndex(globalbin, 0);
 }
 
-double IHistogram::getXaxisValue(size_t binGlobalIndex)
+int IHistogram::getYaxisIndex(size_t globalbin) const
+{
+    return m_data.getAxisBinIndex(globalbin, 1);
+}
+
+double IHistogram::getXaxisValue(size_t globalbin)
 {
     check_x_axis();
-    return m_data.getAxisValue(binGlobalIndex, 0);
+    return m_data.getAxisValue(globalbin, 0);
 }
 
-double IHistogram::getYaxisValue(size_t binGlobalIndex)
+double IHistogram::getYaxisValue(size_t globalbin)
 {
     check_y_axis();
-    return m_data.getAxisValue(binGlobalIndex, 1);
+    return m_data.getAxisValue(globalbin, 1);
 }
 
-void IHistogram::reset()
+double IHistogram::getBinContent(int bin) const
 {
-    m_data.setAllTo(CumulativeValue());
+    return m_data[bin].getValue();
+}
+
+double IHistogram::getBinContent(int binx, int biny) const
+{
+    return getBinContent(getGlobalBin(binx, biny));
+}
+
+double IHistogram::getBinError(int bin) const
+{
+    return m_data[bin].getRMS();
+}
+
+double IHistogram::getBinError(int binx, int biny) const
+{
+    return getBinError(getGlobalBin(binx, biny));
+}
+
+int IHistogram::getBinNumberOfEntries(int bin) const
+{
+    return m_data[bin].getNumberOfEntries();
+}
+
+int IHistogram::getBinNumberOfEntries(int binx, int biny) const
+{
+    return getBinNumberOfEntries(getGlobalBin(binx, biny));
 }
 
 PyObject *IHistogram::getArray() const
@@ -120,6 +154,12 @@ PyObject *IHistogram::getArray() const
     return array.getArray();
 }
 
+void IHistogram::reset()
+{
+    m_data.setAllTo(CumulativeValue());
+}
+
+
 void IHistogram::check_x_axis() const
 {
     if(getRank() <1) {
diff --git a/Core/Tools/src/OutputData.cpp b/Core/Tools/src/OutputData.cpp
index b1bc0ee11cf1278eadcfe59db44dc8699a9bcbb3..e37c158b874d2faa424a8d02fb017ccf3e16e828 100644
--- a/Core/Tools/src/OutputData.cpp
+++ b/Core/Tools/src/OutputData.cpp
@@ -57,8 +57,24 @@ PyObject *OutputData<double>::getArray() const
         *array_buffer++ = (*this)[index];
     }
 
+//    for(size_t index=0; index<getAllocatedSize(); ++index) {
+
+//        std::vector<int> indices = getAxesBinIndices(index);
+//        int offset(0);
+//        int step_size = 1;
+//        for(size_t i_axis=0; i_axis<getRank(); ++i_axis ) {
+//            offset += indices[i_axis]*step_size;
+//            step_size *= dimensions[i_axis];
+//        }
+//        std::cout << index << " offset" << offset << std::endl;
+
+//        *(array_buffer+offset) = (*this)[index];
+//    }
+
     return pyarray;
 }
 
+
+
 #endif
 
diff --git a/Core/Tools/src/VariableBinAxis.cpp b/Core/Tools/src/VariableBinAxis.cpp
index dd5add00a1fe2318b08da2a1ffc8faafeb8e082d..1942c84f1f7737f39cd4abc10b3f9993c60585e0 100644
--- a/Core/Tools/src/VariableBinAxis.cpp
+++ b/Core/Tools/src/VariableBinAxis.cpp
@@ -25,7 +25,8 @@ VariableBinAxis::VariableBinAxis(const std::string &name, size_t nbins, const st
     , m_nbins(nbins)
 {
     if(m_nbins != bin_boundaries.size()-1)
-        throw Exceptions::LogicErrorException("VariableBinAxis::VariableBinAxis() -> Error! The size of value_vector should be of size [nbins+1].");
+        throw Exceptions::LogicErrorException("VariableBinAxis::VariableBinAxis() -> Error! "
+            "The size of bin_boundaries should be of size [nbins+1].");
 
     setBinBoundaries(bin_boundaries);
 }
@@ -79,6 +80,11 @@ double VariableBinAxis::getMax() const
     return m_bin_boundaries.back();
 }
 
+double VariableBinAxis::getBinCenter(size_t index) const
+{
+    return getBin(index).getMidPoint();
+}
+
 
 size_t VariableBinAxis::findClosestIndex(double value) const
 {
diff --git a/Tests/UnitTests/TestCore/CumulativeValueTest.h b/Tests/UnitTests/TestCore/CumulativeValueTest.h
index f3390ee7d26e5457ad8341b797038aa699f91d42..39ca91ec42754fc871af13ff2b0e4224e025818a 100644
--- a/Tests/UnitTests/TestCore/CumulativeValueTest.h
+++ b/Tests/UnitTests/TestCore/CumulativeValueTest.h
@@ -73,7 +73,6 @@ TEST_F(CumulativeValueTest, AddValuesWithWeights)
     EXPECT_DOUBLE_EQ(10.0, cv1.getValue());
     EXPECT_DOUBLE_EQ(2.0, cv1.getAverage());
     EXPECT_FLOAT_EQ(1.0, cv1.getRMS());
-
 }
 
 
diff --git a/Tests/UnitTests/TestCore/FixedBinAxisTest.h b/Tests/UnitTests/TestCore/FixedBinAxisTest.h
index 0464a4f20ddf62c8c5d972b7c43c2d53331b4754..1a017bdf63ba66bd6ab68b84e46e173e5eb46aa3 100644
--- a/Tests/UnitTests/TestCore/FixedBinAxisTest.h
+++ b/Tests/UnitTests/TestCore/FixedBinAxisTest.h
@@ -140,6 +140,11 @@ TEST_F(FixedBinAxisTest, BinCenters)
     EXPECT_DOUBLE_EQ(-1.0, centers[0]);
     EXPECT_DOUBLE_EQ(0.0, centers[1]);
     EXPECT_DOUBLE_EQ(1.0, centers[2]);
+
+    EXPECT_DOUBLE_EQ(axis.getBinCenter(0), centers[0]);
+    EXPECT_DOUBLE_EQ(axis.getBinCenter(1), centers[1]);
+    EXPECT_DOUBLE_EQ(axis.getBinCenter(2), centers[2]);
+
 }
 
 TEST_F(FixedBinAxisTest, BinBoundaries)
diff --git a/Tests/UnitTests/TestCore/Histogram1DTest.h b/Tests/UnitTests/TestCore/Histogram1DTest.h
new file mode 100644
index 0000000000000000000000000000000000000000..1191962beba08218494c9d1099c0fcf128ed9aac
--- /dev/null
+++ b/Tests/UnitTests/TestCore/Histogram1DTest.h
@@ -0,0 +1,114 @@
+#ifndef HISTOGRAM1DTEST_H
+#define HISTOGRAM1DTEST_H
+
+#include "Histogram1D.h"
+#include "Exceptions.h"
+#include <boost/assign/list_of.hpp>
+#include "gtest/gtest.h"
+
+class Histogram1DTest : public ::testing::Test
+{
+ protected:
+    Histogram1DTest(){}
+    virtual ~Histogram1DTest(){}
+};
+
+TEST_F(Histogram1DTest, FixedBinConstructor)
+{
+    Histogram1D hist(5, 0.0, 5.0);
+
+    EXPECT_EQ(1, hist.getRank());
+    EXPECT_EQ(5, hist.getTotalNumberOfBins());
+    EXPECT_EQ(0.0, hist.getXmin());
+    EXPECT_EQ(5.0, hist.getXmax());
+    EXPECT_EQ(std::string("x-axis"), hist.getXaxis()->getName());
+    EXPECT_THROW(hist.getYaxis(), LogicErrorException);
+    for(size_t index=0; index<hist.getTotalNumberOfBins(); ++index) {
+        EXPECT_EQ(index, hist.getGlobalBin(index));
+        EXPECT_EQ(index, hist.getXaxisIndex(index));
+    }
+}
+
+TEST_F(Histogram1DTest, FixedBinDefaultContent)
+{
+    Histogram1D hist(5, 0.0, 5.0);
+
+    // bin centers
+    std::vector<double> bin_centers = boost::assign::list_of(0.5)(1.5)(2.5)(3.5)(4.5);
+    std::vector<double> centers = hist.getBinCenters();
+    for(size_t index=0; index < bin_centers.size(); ++index) {
+        EXPECT_EQ(centers[index], bin_centers[index]);
+        EXPECT_EQ(hist.getXaxisValue(index), bin_centers[index]);
+        EXPECT_EQ(hist.getXaxis()->getBinCenter(index), bin_centers[index]);
+    }
+
+    // default bin values
+    std::vector<double> values = hist.getBinValues();
+    for(size_t index=0; index < bin_centers.size(); ++index) {
+        EXPECT_EQ(hist.getBinContent(index), 0.0);
+        EXPECT_EQ(values[index], 0.0);
+    }
+
+    // default bin errors
+    std::vector<double> errors = hist.getBinErrors();
+    for(size_t index=0; index < bin_centers.size(); ++index) {
+        EXPECT_EQ(hist.getBinError(index), 0.0);
+        EXPECT_EQ(errors[index], 0.0);
+    }
+
+    // default bin entries
+    for(size_t index=0; index < bin_centers.size(); ++index) {
+        EXPECT_EQ(hist.getBinNumberOfEntries(index), 0);
+    }
+
+}
+
+TEST_F(Histogram1DTest, FixedBinFill)
+{
+    Histogram1D hist(5, 0.0, 5.0);
+
+    // filling two different bins
+
+    hist.fill(0.5, 88.0);
+    hist.fill(4.5, 99.0);
+    EXPECT_EQ(hist.getBinContent(0), 88.0);
+    EXPECT_EQ(hist.getBinNumberOfEntries(0), 1);
+    EXPECT_EQ(hist.getBinError(0), 0.0);
+
+    EXPECT_EQ(hist.getBinContent(4), 99.0);
+    EXPECT_EQ(hist.getBinNumberOfEntries(4), 1);
+    EXPECT_EQ(hist.getBinError(4), 0.0);
+
+    std::vector<double> values = boost::assign::list_of(88.0)(0.0)(0.0)(0.0)(99.0);
+    for(size_t index=0; index<hist.getTotalNumberOfBins(); ++index) {
+        EXPECT_EQ(hist.getBinValues()[index], values[index]);
+        EXPECT_EQ(hist.getBinErrors()[index], 0.0);
+    }
+
+    // resetting histograms
+    hist.reset();
+    EXPECT_EQ(hist.getBinContent(0), 0.0);
+    EXPECT_EQ(hist.getBinNumberOfEntries(0), 0);
+    EXPECT_EQ(hist.getBinError(0), 0.0);
+    EXPECT_EQ(hist.getBinContent(4), 0.0);
+    EXPECT_EQ(hist.getBinNumberOfEntries(4), 0);
+    EXPECT_EQ(hist.getBinError(4), 0.0);
+
+    // another fill
+    const double xvalue(1.5);
+    const int xbin = 1;
+
+    hist.fill(xvalue);
+    hist.fill(xvalue, 3.0);
+    EXPECT_EQ(2, hist.getBinNumberOfEntries(xbin));
+    EXPECT_EQ(4.0, hist.getBinContent(xbin));
+    EXPECT_EQ(1.0, hist.getBinError(xbin));
+
+}
+
+
+
+
+
+#endif
+
diff --git a/Tests/UnitTests/TestCore/Histogram2DTest.h b/Tests/UnitTests/TestCore/Histogram2DTest.h
new file mode 100644
index 0000000000000000000000000000000000000000..64bafd30a1d3d9015bc5719c3d5f1851c9a71e4c
--- /dev/null
+++ b/Tests/UnitTests/TestCore/Histogram2DTest.h
@@ -0,0 +1,139 @@
+#ifndef HISTOGRAM2DTEST_H
+#define HISTOGRAM2DTEST_H
+
+#include "Histogram2D.h"
+#include <boost/assign/list_of.hpp>
+#include "gtest/gtest.h"
+
+class Histogram2DTest : public ::testing::Test
+{
+ protected:
+    Histogram2DTest();
+    virtual ~Histogram2DTest();
+
+    Histogram2D *hist;
+};
+
+
+// y
+// 4.0   -----------------------------------
+//       |     |          |     |          |
+//       |  2  |    5     |  8  |    11    |
+// 2.0   -----------------------------------
+//       |  1  |    4     |  7  |    10    |
+// 1.0   -----------------------------------
+//       |  0  |    3     |  6  |    9     |
+// 0.0   -----------------------------------
+//     -1.0  -0.5        0.5   1.0        2.0  X
+
+Histogram2DTest::Histogram2DTest()
+{
+    std::vector<double> xbin_edges = boost::assign::list_of(-1.0)(-0.5)(0.5)(1.0)(2.0);
+    std::vector<double> ybin_edges = boost::assign::list_of(0.0)(1.0)(2.0)(4.0);
+    hist = new Histogram2D(4, xbin_edges, 3, ybin_edges);
+}
+
+Histogram2DTest::~Histogram2DTest()
+{
+    delete hist;
+}
+
+TEST_F(Histogram2DTest, VariableHist)
+{
+    hist->reset();
+
+    // basic axes check
+    EXPECT_EQ(12, hist->getTotalNumberOfBins());
+    EXPECT_EQ(hist->getRank(), 2);
+    EXPECT_EQ(hist->getXaxis()->getName(), std::string("x-axis"));
+    EXPECT_EQ(hist->getXaxis()->getSize(), 4);
+    EXPECT_EQ(hist->getXmin(), -1.0);
+    EXPECT_EQ(hist->getXmax(), 2.0);
+    EXPECT_EQ(hist->getYaxis()->getName(), std::string("y-axis"));
+    EXPECT_EQ(hist->getYaxis()->getSize(), 3);
+    EXPECT_EQ(hist->getYmin(), 0.0);
+    EXPECT_EQ(hist->getYmax(), 4.0);
+
+    // globalbin -> axes indices
+    EXPECT_EQ(hist->getXaxisIndex(0), 0);
+    EXPECT_EQ(hist->getXaxisIndex(1), 0);
+    EXPECT_EQ(hist->getXaxisIndex(2), 0);
+    EXPECT_EQ(hist->getXaxisIndex(3), 1);
+    EXPECT_EQ(hist->getXaxisIndex(4), 1);
+    EXPECT_EQ(hist->getXaxisIndex(5), 1);
+    EXPECT_EQ(hist->getXaxisIndex(9), 3);
+    EXPECT_EQ(hist->getXaxisIndex(10), 3);
+    EXPECT_EQ(hist->getXaxisIndex(11), 3);
+
+    EXPECT_EQ(hist->getYaxisIndex(0), 0);
+    EXPECT_EQ(hist->getYaxisIndex(1), 1);
+    EXPECT_EQ(hist->getYaxisIndex(2), 2);
+    EXPECT_EQ(hist->getYaxisIndex(3), 0);
+    EXPECT_EQ(hist->getYaxisIndex(4), 1);
+    EXPECT_EQ(hist->getYaxisIndex(5), 2);
+    EXPECT_EQ(hist->getYaxisIndex(9), 0);
+    EXPECT_EQ(hist->getYaxisIndex(10), 1);
+    EXPECT_EQ(hist->getYaxisIndex(11), 2);
+
+    // axes indices -> global bin
+    EXPECT_EQ(hist->getGlobalBin(0,0), 0);
+    EXPECT_EQ(hist->getGlobalBin(0,2), 2);
+    EXPECT_EQ(hist->getGlobalBin(1,1), 4);
+    EXPECT_EQ(hist->getGlobalBin(3,2), 11);
+
+    // bin centers
+    EXPECT_EQ(hist->getXaxisValue(0), -0.75);
+    EXPECT_EQ(hist->getXaxisValue(2), -0.75);
+    EXPECT_EQ(hist->getXaxisValue(4), 0.0);
+    EXPECT_EQ(hist->getXaxisValue(10), 1.5);
+    EXPECT_EQ(hist->getXaxisValue(11), 1.5);
+
+    EXPECT_EQ(hist->getYaxisValue(0), 0.5);
+    EXPECT_EQ(hist->getYaxisValue(2), 3.0);
+    EXPECT_EQ(hist->getYaxisValue(4), 1.5);
+    EXPECT_EQ(hist->getYaxisValue(10), 1.5);
+    EXPECT_EQ(hist->getYaxisValue(11), 3.0);
+
+}
+
+// y
+// 4.0   -----------------------------------
+//       |     |          |     |          |
+//       |  2  |    5     |  8  |    11    |
+// 2.0   -----------------------------------
+//       |  1  |    4     |  7  |    10    |
+// 1.0   -----------------------------------
+//       |  0  |    3     |  6  |    9     |
+// 0.0   -----------------------------------
+//     -1.0  -0.5        0.5   1.0        2.0  X
+
+TEST_F(Histogram2DTest, VariableHistFill)
+{
+    hist->reset();
+
+    // values to fill all histogram
+    std::vector<double> xvalues = boost::assign::list_of(-0.75)(-0.75)(-0.75)(0.0)(0.0)(0.0)(0.75)(0.75)(0.75)(1.5)(1.5)(1.5);
+    std::vector<double> yvalues = boost::assign::list_of(0.5)(1.5)(3.0)(0.5)(1.5)(3.0)(0.5)(1.5)(3.0)(0.5)(1.5)(3.0);
+
+    for(size_t i=0; i<xvalues.size(); ++i) {
+        hist->fill(xvalues[i], yvalues[i], i*10.0);
+    }
+
+    for(size_t globalbin=0; globalbin<hist->getTotalNumberOfBins(); ++globalbin) {
+        EXPECT_EQ(globalbin*10.0, hist->getBinContent(globalbin));
+        EXPECT_EQ(1.0, hist->getBinNumberOfEntries(globalbin));
+
+    }
+
+    for(size_t binx=0; binx<hist->getXaxis()->getSize(); ++binx){
+        for(size_t biny=0; biny<hist->getYaxis()->getSize(); ++biny){
+            int globalbin = hist->getGlobalBin(binx, biny);
+            EXPECT_EQ(globalbin*10.0, hist->getBinContent(binx, biny));
+            EXPECT_EQ(1.0, hist->getBinNumberOfEntries(binx, biny));
+        }
+    }
+
+}
+
+
+#endif
diff --git a/Tests/UnitTests/TestCore/main.cpp b/Tests/UnitTests/TestCore/main.cpp
index b5659bdc08b33221573d036df6bb4055af5edba1..bc5dd1e39c05d05e4c992d8edbe6fb8b91b66090 100644
--- a/Tests/UnitTests/TestCore/main.cpp
+++ b/Tests/UnitTests/TestCore/main.cpp
@@ -51,6 +51,9 @@
 #include "ParameterDistributionTest.h"
 #include "UtilsTest.h"
 #include "CumulativeValueTest.h"
+#include "Histogram1DTest.h"
+#include "Histogram2DTest.h"
+
 
 
 struct ErrorStreamRedirect {