diff --git a/App/src/IsGISAXSTools.cpp b/App/src/IsGISAXSTools.cpp
index 5d3fa1d84b4592028e593854b6c56026aa4a54a0..ecc8d25960ce12e4a64a41bfb17c2a65bd81cb32 100644
--- a/App/src/IsGISAXSTools.cpp
+++ b/App/src/IsGISAXSTools.cpp
@@ -81,7 +81,7 @@ void IsGISAXSTools::drawOutputDataInPad(const OutputData<double>& output,
     }
 
     TH1 *hist(0);
-    if(output.getNdimensions() == 2 ) {
+    if(output.getRank() == 2 ) {
         hist = IsGISAXSTools::getOutputDataTH2D(output, "p_hist1D");
     } else {
         hist = IsGISAXSTools::getOutputDataTH123D(output, "p_hist1D");
@@ -120,15 +120,15 @@ TH2D *IsGISAXSTools::getOutputDataTH2D(const OutputData<double>& output,
                                        const std::string& histo_name)
 {
     assert(&output);
-    if (output.getNdimensions() !=2)
+    if (output.getRank() !=2)
         throw( "IsGISAXSTools::getOutputDataTH2D() -> "
                "Warning! Expected number of dimensiobs is 2.");
 
     std::vector<AxisStructure > haxises;
-    haxises.resize(output.getNdimensions());
+    haxises.resize(output.getRank());
 
     // we assume variable bin size and prepare [nbins+1] array of left edges of each bin plus right edge of the last bin
-    for(size_t i_axis=0; i_axis<output.getNdimensions(); ++i_axis) {
+    for(size_t i_axis=0; i_axis<output.getRank(); ++i_axis) {
         const IAxis *axis = output.getAxis(i_axis);
         if( !axis ) throw("IsGISAXSTools::getOutputDataTH123D() -> Error! Can't cast axis");
         double dx(0);
@@ -192,16 +192,16 @@ TH2D *IsGISAXSTools::getOutputDataTH2D(const OutputData<double>& output,
 TH1 *IsGISAXSTools::getOutputDataTH123D(const OutputData<double>& output, const std::string &histo_name)
 {
     assert(&output);
-    if (output.getNdimensions() >3) {
+    if (output.getRank() >3) {
         std::cout << "IsGISAXSTools::getOutputDataTH123D() -> Warning! Expected number of dimensions should be not more than 3" << std::endl;
         return 0;
     }
 
     std::vector<AxisStructure > haxises;
-    haxises.resize(output.getNdimensions());
+    haxises.resize(output.getRank());
 
     // we assume variable bin size and prepare [nbins+1] array of left edges of each bin plus right edge of the last bin
-    for(size_t i_axis=0; i_axis<output.getNdimensions(); ++i_axis) {
+    for(size_t i_axis=0; i_axis<output.getRank(); ++i_axis) {
         const IAxis *axis = output.getAxis(i_axis);
         if( !axis ) throw("IsGISAXSTools::getOutputDataTH123D() -> Error! Can't cast axis");
         double dx(0);
@@ -240,16 +240,16 @@ TH1 *IsGISAXSTools::getOutputDataTH123D(const OutputData<double>& output, const
     TH1D *hist1(0);
     TH2D *hist2(0);
     TH3D *hist3(0);
-    if(output.getNdimensions() == 1) {
+    if(output.getRank() == 1) {
         hist1 = new TH1D(histo_name.c_str(), histo_name.c_str(), (int)haxises[0].nbins, &haxises[0].xbins[0]);
         hist1->GetXaxis()->SetTitle( haxises[0].name.c_str() );
         hist = hist1;
-    } else if(output.getNdimensions() == 2) {
+    } else if(output.getRank() == 2) {
         hist2 = new TH2D(histo_name.c_str(), histo_name.c_str(), (int)haxises[0].nbins, &haxises[0].xbins[0], (int)haxises[1].nbins, &haxises[1].xbins[0]);
         hist2->GetXaxis()->SetTitle( haxises[0].name.c_str() );
         hist2->GetYaxis()->SetTitle( haxises[1].name.c_str() );
         hist = hist2;
-    } else if(output.getNdimensions() == 3) {
+    } else if(output.getRank() == 3) {
         hist3 = new TH3D(histo_name.c_str(), histo_name.c_str(), (int)haxises[0].nbins, &haxises[0].xbins[0], (int)haxises[1].nbins, &haxises[1].xbins[0], (int)haxises[1].nbins, &haxises[1].xbins[0]);
         hist3->GetXaxis()->SetTitle( haxises[0].name.c_str() );
         hist3->GetYaxis()->SetTitle( haxises[1].name.c_str() );
@@ -478,7 +478,7 @@ TLine *IsGISAXSTools::getOutputDataScanLine(const OutputData<double> &data)
 {
     assert(&data);
 
-    if(data.getNdimensions() != 2) throw LogicErrorException("IsGISAXSTools::getOutputDataScanLine() -> Error! Number of dimensions should be 2");
+    if(data.getRank() != 2) throw LogicErrorException("IsGISAXSTools::getOutputDataScanLine() -> Error! Number of dimensions should be 2");
     double x1(0), x2(0), y1(0), y2(0);
     if( data.getAxis("alpha_f") && data.getAxis("alpha_f")->getSize() == 1) {
         // horizontal line
@@ -509,7 +509,7 @@ TH1D *IsGISAXSTools::getOutputDataScanHist(const OutputData<double> &data, const
 {
     assert(&data);
 
-    if(data.getNdimensions() != 2) throw LogicErrorException("IsGISAXSTools::getOutputDataScanHist() -> Error! Number of dimensions should be 2");
+    if(data.getRank() != 2) throw LogicErrorException("IsGISAXSTools::getOutputDataScanHist() -> Error! Number of dimensions should be 2");
     // one of axis should have dimension 1
     if( (data.getAxis("alpha_f") && data.getAxis("alpha_f")->getSize() != 1) && (data.getAxis("phi_f") && data.getAxis("phi_f")->getSize() != 1))
     {
diff --git a/App/src/TestIsGISAXS12.cpp b/App/src/TestIsGISAXS12.cpp
index 9996bc02cff8f15c53d59420b03e6b672792eb98..2f7341bf9f7107e1a05aa1373f4c35b00f22cc50 100644
--- a/App/src/TestIsGISAXS12.cpp
+++ b/App/src/TestIsGISAXS12.cpp
@@ -199,7 +199,6 @@ void TestIsGISAXS12::plot_isgisaxs_fit_results()
     c1->cd(6); leg3->Draw();
 
     c1->Update();
-
 }
 
 
@@ -532,7 +531,7 @@ void TestIsGISAXS12::print_axes(IsGISAXSData::DataSet_t &data)
 {
     for(size_t i_set=0; i_set<data.size(); ++i_set) {
         std::cout << "scan #" << i_set << "  ";
-        for(size_t i_axis=0; i_axis<data[(int)i_set]->getNdimensions(); ++i_axis) {
+        for(size_t i_axis=0; i_axis<data[(int)i_set]->getRank(); ++i_axis) {
             const IAxis *axis = data[(int)i_set]->getAxis(i_axis);
             std::cout << "( " << axis->getName() << ", " << axis->getSize() << ", " << axis->getMin() << ", " << axis->getMax() << " )   ";
         }
diff --git a/Core/Algorithms/inc/DWBASimulation.h b/Core/Algorithms/inc/DWBASimulation.h
index e3d04f808b423a96b112904522cda8bb62dece24..23524adc78cc14fd15d591164518960d21e29c39 100644
--- a/Core/Algorithms/inc/DWBASimulation.h
+++ b/Core/Algorithms/inc/DWBASimulation.h
@@ -38,7 +38,7 @@ class DWBASimulation : public ISimulation
     virtual void setThreadInfo(const ThreadInfo& thread_info)
     { m_thread_info = thread_info; }
 
-    //! return output data containing calculated intensity
+    //! Return output data containing calculated intensity
     const OutputData<double>& getDWBAIntensity() const
     { return m_dwba_intensity; }
 
@@ -46,7 +46,6 @@ class DWBASimulation : public ISimulation
     void addDWBAIntensity(const OutputData<double>& data_to_add)
     { m_dwba_intensity += data_to_add; }
 
-    //! clone DWBA simulation
     virtual DWBASimulation *clone() const;
 
     // ---------------------------------
@@ -57,16 +56,16 @@ class DWBASimulation : public ISimulation
     typedef OutputDataIterator<const double, const OutputData<double> >
         const_iterator;
 
-    //! return a read/write iterator that points to the first element
+    //! Return a read/write iterator that points to the first element
     iterator begin();
 
-    //! return a read-only iterator that points to the first element
+    //! Return a read-only iterator that points to the first element
     const_iterator begin() const;
 
-    //! return a read/write iterator that points to the one past last element
+    //! Return a read/write iterator that points to the one past last element
     const iterator end() { return m_dwba_intensity.end(); }
 
-    //! return a read-only iterator that points to the one past last element
+    //! Return a read-only iterator that points to the one past last element
     const const_iterator end() const { return m_dwba_intensity.end(); }
 
  protected:
diff --git a/Core/Algorithms/inc/DiffuseDWBASimulation.h b/Core/Algorithms/inc/DiffuseDWBASimulation.h
index e6da8b69fd24ec6d6c70c2c0c56f0e5040d7ad27..4fa7b781d20d737fd691e9a15d846379d1427285 100644
--- a/Core/Algorithms/inc/DiffuseDWBASimulation.h
+++ b/Core/Algorithms/inc/DiffuseDWBASimulation.h
@@ -24,15 +24,26 @@ class DiffuseDWBASimulation: public LayerDWBASimulation
  public:
     DiffuseDWBASimulation()
         : m_refractive_index(1., 0.), m_surface_density(1.) {}
+
     virtual ~DiffuseDWBASimulation() {}
 
     virtual void run();
+
     size_t getSize() const { return m_np_infos.size(); }
-    void addParticleInfo(DiffuseParticleInfo *p_info);
+
+    void addParticleInfo(DiffuseParticleInfo *p_info)
+    { m_np_infos.push_back(p_info); }
+
     void setRefractiveIndex(complex_t n) { m_refractive_index = n; }
+
     void setSurfaceDensity(double surface_density)
     { m_surface_density = surface_density; }
-    void rescaleAbundances(double factor);
+
+    void rescaleAbundances(double factor)
+    {
+        for (size_t np_index=0; np_index<m_np_infos.size(); ++np_index)
+            m_np_infos[np_index]->scaleAbundance(factor);
+    }
 
  protected:
     complex_t m_refractive_index;
@@ -41,7 +52,11 @@ class DiffuseDWBASimulation: public LayerDWBASimulation
     struct DiffuseFormFactorTerm {
         DiffuseFormFactorTerm()
             : m_form_factors(), m_probabilities(), m_factor(0.) {}
-        ~DiffuseFormFactorTerm();
+        ~DiffuseFormFactorTerm()
+        {
+            for (size_t i=0; i<m_form_factors.size(); ++i)
+                delete m_form_factors[i];
+        }
         std::vector<IFormFactor*> m_form_factors;
         std::vector<double> m_probabilities;
         double m_factor;
diff --git a/Core/Algorithms/inc/Mask.h b/Core/Algorithms/inc/Mask.h
index dff929b3651a998cefb34df48c1a6350d35629b0..f30db9a997e83e615a1394990d6715c151be80c3 100644
--- a/Core/Algorithms/inc/Mask.h
+++ b/Core/Algorithms/inc/Mask.h
@@ -18,7 +18,7 @@
 
 #include "MaskCoordinateFunction.h"
 
-//! Base class for masking OutputData elements
+//! Base class for masking OutputData elements.
 
 class Mask : public ICloneable
 {
@@ -46,7 +46,7 @@ class Mask : public ICloneable
     size_t nextSubIndex(size_t total_index);
 };
 
-//! Mask based on the index modulo a given number
+//! Mask based on the index modulo a given number.
 
 class MaskIndexModulus : public Mask
 {
@@ -64,7 +64,7 @@ class MaskIndexModulus : public Mask
     size_t m_remainder;
 };
 
-//! Mask based on the coordinates
+//! Mask based on the coordinates.
 
 class MaskCoordinates : public Mask
 {
diff --git a/Core/Algorithms/inc/ThreadInfo.h b/Core/Algorithms/inc/ThreadInfo.h
index bf2ac3ca62fb0eb5212ce51dccbeff7e564265e9..287fd03c3e04368197565d935fe42ad309af72e9 100644
--- a/Core/Algorithms/inc/ThreadInfo.h
+++ b/Core/Algorithms/inc/ThreadInfo.h
@@ -20,15 +20,9 @@
 
 struct ThreadInfo
 {
-    ThreadInfo();
+    ThreadInfo() : n_threads(1), i_thread(0) {}
     size_t n_threads;
     size_t i_thread;
 };
 
-inline ThreadInfo::ThreadInfo()
-: n_threads(1)
-, i_thread(0)
-{
-}
-
 #endif /* THREADINFO_H_ */
diff --git a/Core/Algorithms/src/DWBASimulation.cpp b/Core/Algorithms/src/DWBASimulation.cpp
index 043f8eadb7d61e254820ca0f78b3dd38eedc535c..769efef08b4057017dad10f6e3e54cbd82ff6601 100644
--- a/Core/Algorithms/src/DWBASimulation.cpp
+++ b/Core/Algorithms/src/DWBASimulation.cpp
@@ -47,7 +47,6 @@ DWBASimulation *DWBASimulation::clone() const
     p_result->m_thread_info = m_thread_info;
     if (mp_simulation)
         p_result->mp_simulation = mp_simulation->clone();
-
     return p_result;
 }
 
diff --git a/Core/Algorithms/src/DiffuseDWBASimulation.cpp b/Core/Algorithms/src/DiffuseDWBASimulation.cpp
index 806fabe168cbedf80587c6752fec849995f6e266..95967d9e42d4192dbcd5ba5c9af3fbf0b87fc583 100644
--- a/Core/Algorithms/src/DiffuseDWBASimulation.cpp
+++ b/Core/Algorithms/src/DiffuseDWBASimulation.cpp
@@ -16,30 +16,32 @@
 #include "DiffuseDWBASimulation.h"
 #include "FormFactorDWBAConstZ.h"
 
-// Run a simulation.
+#define SQR(x) ((x)*(x))
+
+//! Carry out one simulation thread.
 
 void DiffuseDWBASimulation::run()
 {
+    // Set diffuse terms.
     std::vector<DiffuseFormFactorTerm*> diffuse_terms;
     size_t nbr_heights = 50;
     size_t samples_per_particle = 9;
     initDiffuseFormFactorTerms(
         diffuse_terms, nbr_heights, samples_per_particle);
-    double wavevector_scattering_factor = M_PI/getWaveLength()/getWaveLength();
+
+    double wavevector_scattering_factor = M_PI/SQR(getWaveLength());
     cvector_t k_ij = m_ki;
     k_ij.setZ(-mp_kz_function->evaluate(-m_alpha_i));
 
-    DWBASimulation::iterator it_intensity = begin();
-    while ( it_intensity != end() ) {
+    for (DWBASimulation::iterator it_intensity =
+             begin(); it_intensity != end(); ++it_intensity ) {
         Bin1D phi_bin = getDWBAIntensity().getBinOfAxis(
             "phi_f", it_intensity.getIndex());
         Bin1D alpha_bin = getDWBAIntensity().getBinOfAxis(
             "alpha_f", it_intensity.getIndex());
         double alpha_f = alpha_bin.getMidPoint();
-        if (alpha_f<0) {
-            ++it_intensity;
+        if (alpha_f<0)
             continue;
-        }
         Bin1DCVector k_f_bin = getKfBin(getWaveLength(), alpha_bin, phi_bin);
 
         double total_intensity = 0;
@@ -55,29 +57,23 @@ void DiffuseDWBASimulation::run()
                 intensity += p_diffuse_term->m_probabilities[j]*std::norm(amp);
             }
             total_intensity +=
-                p_diffuse_term->m_factor*(intensity - std::norm(amplitude));
+                p_diffuse_term->m_factor * (intensity - std::norm(amplitude));
         }
-        *it_intensity = total_intensity*
-            wavevector_scattering_factor*wavevector_scattering_factor;
-        ++it_intensity;
+        *it_intensity = total_intensity * SQR(wavevector_scattering_factor);
     }
 
     for (size_t i=0; i<diffuse_terms.size(); ++i)
         delete diffuse_terms[i];
 }
 
-void DiffuseDWBASimulation::addParticleInfo(DiffuseParticleInfo *p_info)
-{
-    m_np_infos.push_back(p_info);
-}
-
-void DiffuseDWBASimulation::rescaleAbundances(double factor)
-{
-    for (size_t np_index=0; np_index<m_np_infos.size(); ++np_index) {
-        m_np_infos[np_index]->scaleAbundance(factor);
-    }
-}
+//! Initialize vector<DiffuseFormFactorTerm*> term.
 
+//! Called near beginning of this->run().
+//! Collect one entry (p_diffuse_term) per particle type and layer[?].
+//! For each entry, set
+//! - m_factor
+//! - m_form_factors <- vector<FormFactorDWBAConstZ*>
+//!
 void DiffuseDWBASimulation::initDiffuseFormFactorTerms(
         std::vector<DiffuseFormFactorTerm*>& terms,
         size_t nbr_heights,
@@ -93,7 +89,7 @@ void DiffuseDWBASimulation::initDiffuseFormFactorTerms(
             DiffuseFormFactorTerm *p_diffuse_term = new DiffuseFormFactorTerm;
             p_diffuse_term->m_factor = density_per_height_per_particle;
             double depth = p_diff_info->getDepth() -
-                (double)j*p_diff_info->getHeightRange()/(nbr_heights-1.0);
+                j*p_diff_info->getHeightRange()/(nbr_heights-1.0);
             std::vector<IFormFactor*> form_factors;
             p_diff_info->getParticle()->getSimpleFormFactor()->createDistributedFormFactors(
                 form_factors, p_diffuse_term->m_probabilities,
@@ -113,9 +109,3 @@ void DiffuseDWBASimulation::initDiffuseFormFactorTerms(
     }
 }
 
-DiffuseDWBASimulation::DiffuseFormFactorTerm::~DiffuseFormFactorTerm()
-{
-    for (size_t i=0; i<m_form_factors.size(); ++i) {
-        delete m_form_factors[i];
-    }
-}
diff --git a/Core/Algorithms/src/Instrument.cpp b/Core/Algorithms/src/Instrument.cpp
index 2ef9107b4d2ea8b4c91574844b5d7809d574fa02..572b4712ca7321f3bbac419fce6fb297411dd452 100644
--- a/Core/Algorithms/src/Instrument.cpp
+++ b/Core/Algorithms/src/Instrument.cpp
@@ -35,7 +35,7 @@ Instrument::Instrument(const Instrument& other)
 void Instrument::matchDetectorParameters(const OutputData<double>& output_data)
 {
     m_detector.clear();
-    for(size_t i_axis=0; i_axis<output_data.getNdimensions(); ++i_axis) {
+    for(size_t i_axis=0; i_axis<output_data.getRank(); ++i_axis) {
         const IAxis *axis = output_data.getAxis(i_axis);
         m_detector.addAxis(*axis);
     }
diff --git a/Core/Fitting/src/MinimizerScan.cpp b/Core/Fitting/src/MinimizerScan.cpp
index d421337255dff9a253c264555695cb9e6c4d4479..4fb41264ea279296202995068dc0306efd299940 100644
--- a/Core/Fitting/src/MinimizerScan.cpp
+++ b/Core/Fitting/src/MinimizerScan.cpp
@@ -25,7 +25,7 @@ void MinimizerScan::minimize()
 
     // scanning values of fit parameters
     for(OutputData<double>::iterator it = m_fcnvalues_map->begin(); it!=m_fcnvalues_map->end(); ++it) {
-        for(size_t i_axis=0; i_axis<m_fcnvalues_map->getNdimensions(); ++i_axis) {
+        for(size_t i_axis=0; i_axis<m_fcnvalues_map->getRank(); ++i_axis) {
             size_t xbin = m_fcnvalues_map->toCoordinate(it.getIndex(), i_axis);
             double value = (*m_fcnvalues_map->getAxis(i_axis))[xbin];
             std::string parname = m_fcnvalues_map->getAxis(i_axis)->getName();
@@ -53,7 +53,7 @@ void MinimizerScan::construct_fcnvalues_map()
             m_fcnvalues_map->addAxis(axis);
         }
     }
-    if( !m_fcnvalues_map->getNdimensions() ) {
+    if( !m_fcnvalues_map->getRank() ) {
         throw LogicErrorException("MinimizerScan::construct_parameter_map() -> Error! No parameters with TAttLimit::limited(left,right) attribute were found.");
     }
     m_fcnvalues_map->setAllTo(0.0);
@@ -64,7 +64,7 @@ void MinimizerScan::set_parvalues_to_minimum()
 {
     assert(m_fcnvalues_map);
     OutputData<double>::iterator it = std::min_element(m_fcnvalues_map->begin(), m_fcnvalues_map->end());
-    for(size_t i_axis=0; i_axis<m_fcnvalues_map->getNdimensions(); ++i_axis) {
+    for(size_t i_axis=0; i_axis<m_fcnvalues_map->getRank(); ++i_axis) {
         size_t xbin = m_fcnvalues_map->toCoordinate(it.getIndex(), i_axis);
         double value = (*m_fcnvalues_map->getAxis(i_axis))[xbin];
         std::string parname = m_fcnvalues_map->getAxis(i_axis)->getName();
diff --git a/Core/FormFactors/inc/IFormFactor.h b/Core/FormFactors/inc/IFormFactor.h
index 6937e523f805a438148b9ed8652e5b041d665cc8..2657214d65192f0803f57078dcbab366776a1e33 100644
--- a/Core/FormFactors/inc/IFormFactor.h
+++ b/Core/FormFactors/inc/IFormFactor.h
@@ -3,7 +3,7 @@
 //  BornAgain: simulate and fit scattering at grazing incidence
 //
 //! @file      FormFactors/inc/IFormFactor.h
-//! @brief     Defines class IFormFactor.
+//! @brief     Defines pure virtual interface class IFormFactor.
 //!
 //! @homepage  http://apps.jcns.fz-juelich.de/BornAgain
 //! @license   GNU General Public License v3 or higher (see COPYING)
diff --git a/Core/FormFactors/inc/IFormFactorBornSeparable.h b/Core/FormFactors/inc/IFormFactorBornSeparable.h
index 865a060f4d77160308d98ff5b8c32edab5966be4..cf54933d5b4045628a28d7adc306ac7792336f74 100644
--- a/Core/FormFactors/inc/IFormFactorBornSeparable.h
+++ b/Core/FormFactors/inc/IFormFactorBornSeparable.h
@@ -3,7 +3,7 @@
 //  BornAgain: simulate and fit scattering at grazing incidence
 //
 //! @file      FormFactors/inc/IFormFactorBornSeparable.h
-//! @brief     Defines class IFormFactorBornSeparable.
+//! @brief     Defines pure virtual interface class IFormFactorBornSeparable.
 //!
 //! @homepage  http://apps.jcns.fz-juelich.de/BornAgain
 //! @license   GNU General Public License v3 or higher (see COPYING)
diff --git a/Core/FormFactors/inc/IFormFactorDecorator.h b/Core/FormFactors/inc/IFormFactorDecorator.h
index 66bb359bb34ebd9fe112865a27f8b05fdce53f08..15eac8d098eb80527cc555cad0708ac07ab67277 100644
--- a/Core/FormFactors/inc/IFormFactorDecorator.h
+++ b/Core/FormFactors/inc/IFormFactorDecorator.h
@@ -3,7 +3,7 @@
 //  BornAgain: simulate and fit scattering at grazing incidence
 //
 //! @file      FormFactors/inc/IFormFactorDecorator.h
-//! @brief     Defines class IFormFactorDecorator.
+//! @brief     Defines pure vurtual interface class IFormFactorDecorator.
 //!
 //! @homepage  http://apps.jcns.fz-juelich.de/BornAgain
 //! @license   GNU General Public License v3 or higher (see COPYING)
diff --git a/Core/PythonAPI/src/PythonOutputData.cpp b/Core/PythonAPI/src/PythonOutputData.cpp
index 205afd51a97b37ef16046e1c6468da9ec2d9b52e..d83f6eae85022e56b58a7a2fd4de5b70a5385250 100644
--- a/Core/PythonAPI/src/PythonOutputData.cpp
+++ b/Core/PythonAPI/src/PythonOutputData.cpp
@@ -38,7 +38,7 @@ PyObject *ExportOutputDataAxis(const OutputData<double > &output, int naxis);
 int GetOutputDataNdimensions(const Simulation &simulation)
 {
     const OutputData<double > *data = simulation.getOutputData();
-    int ndims = data->getNdimensions();
+    int ndims = data->getRank();
     return ndims;
 }
 
@@ -66,7 +66,7 @@ PyObject *ExportOutputData(const OutputData<double > &output_data)
 {
     // getting size of dimensions from output_data
     std::vector<int > dimensions;
-    for(size_t i=0; i<output_data.getNdimensions(); i++) {
+    for(size_t i=0; i<output_data.getRank(); i++) {
         //const AxisDouble *axis = output_data.getAxis(i);
         const IAxis *axis = output_data.getAxis(i);
         dimensions.push_back( axis->getSize() );
@@ -106,7 +106,7 @@ PyObject *ExportOutputDataAxis(const OutputData<double > &output_data, int naxis
 {
     // getting size of dimensions from output_data
     std::vector<int > dimensions;
-    for(size_t i=0; i<output_data.getNdimensions(); i++) {
+    for(size_t i=0; i<output_data.getRank(); i++) {
         const IAxis *axis = output_data.getAxis(i);
         dimensions.push_back( axis->getSize() );
     }
diff --git a/Core/Tools/inc/ICloneable.h b/Core/Tools/inc/ICloneable.h
index 9c3880bb6acf26aa30ae6a8a8c51ade1f1da42c4..944b724408c62399a25a5f9da2c6e71ca41bc009 100644
--- a/Core/Tools/inc/ICloneable.h
+++ b/Core/Tools/inc/ICloneable.h
@@ -18,7 +18,7 @@
 
 #include "Exceptions.h"
 
-//! Disable copy constructor and = operator; child classes must provide clone().
+//! Disable copy constructor and operator=; child classes should provide clone().
 
 class ICloneable
 {
@@ -28,8 +28,10 @@ class ICloneable
     virtual ICloneable *clone() const = 0;
 
  private:
-    ICloneable(const ICloneable &) { throw NotImplementedException("ICloneable(const ICloneable &) -> Error: not implemented."); }
-    ICloneable &operator=(const ICloneable &) { throw NotImplementedException("ICloneable &operator=(const ICloneable &) -> Error: not implemented."); }
+    ICloneable(const ICloneable &) { throw NotImplementedException(
+        "ICloneable(const ICloneable &) -> Error: not implemented."); }
+    ICloneable &operator=(const ICloneable &) { throw NotImplementedException(
+        "ICloneable &operator=(const ICloneable &) -> Error: not implemented."); }
 };
 
 #endif // ICLONEABLE_H
diff --git a/Core/Tools/inc/IFactory.h b/Core/Tools/inc/IFactory.h
index 1f586cf9305de39d7ec89bc24b5d02f8b86b9232..a19f21f9e01ceb287668488b7191a4f5822f220b 100644
--- a/Core/Tools/inc/IFactory.h
+++ b/Core/Tools/inc/IFactory.h
@@ -3,7 +3,7 @@
 //  BornAgain: simulate and fit scattering at grazing incidence
 //
 //! @file      Tools/inc/IFactory.h
-//! @brief     Defines class IFactory.
+//! @brief     Defines interface class IFactory.
 //!
 //! @homepage  http://apps.jcns.fz-juelich.de/BornAgain
 //! @license   GNU General Public License v3 or higher (see COPYING)
diff --git a/Core/Tools/inc/INamed.h b/Core/Tools/inc/INamed.h
index a2ff2cb0cf666b68f58deca51f87f032e442a2da..98e6bd64e4af1b6487789c9c8666e224a5092583 100644
--- a/Core/Tools/inc/INamed.h
+++ b/Core/Tools/inc/INamed.h
@@ -3,7 +3,7 @@
 //  BornAgain: simulate and fit scattering at grazing incidence
 //
 //! @file      Tools/inc/INamed.h
-//! @brief     Defines class INamed.
+//! @brief     Defines the standard mix-in INamed.
 //!
 //! @homepage  http://apps.jcns.fz-juelich.de/BornAgain
 //! @license   GNU General Public License v3 or higher (see COPYING)
diff --git a/Core/Tools/inc/OutputData.h b/Core/Tools/inc/OutputData.h
index cca0cdbbb85b16a49a2bd8dd65348e949acd0e40..e7a4558c8533336956d2861fd28e7f238b631051 100644
--- a/Core/Tools/inc/OutputData.h
+++ b/Core/Tools/inc/OutputData.h
@@ -17,23 +17,20 @@
 #define OUTPUTDATA_H
 
 #include "AxisDouble.h"
-//#include "Exceptions.h"
 #include "Types.h"
 #include "LLData.h"
 #include "OutputDataIterator.h"
 #include "SafePointerVector.h"
 
-//#include <string>
 #include <sstream>
 
-//! Store data of any type in multi-dimensional space
+//! Store data of any type in multi-dimensional space.
 
 template <class T> class OutputData
 {
  public:
     OutputData();
     ~OutputData();
-    //! make object clone
     OutputData* clone() const;
 
     void copyFrom(const OutputData<T> &x);
@@ -49,29 +46,19 @@ template <class T> class OutputData
     // retrieve basic info
     // ---------------------------------
 
-    //! return number of dimensions
+    //! Return number of dimensions.
     size_t getRank() const { return m_value_axes.size(); }
 
-    //! return number of dimensions (same as above)
-    size_t getNdimensions() const { return getRank(); }
-
-    //! return total size of data buffer (product of bin number in every dimension)
+   //! Return total size of data buffer (product of bin number in every dimension).
     size_t getAllocatedSize() const {
         if (mp_ll_data) return mp_ll_data->getTotalSize();
         return 0;
     }
 
-    //! return total size of data buffer (product of bin number in every dimension)
-    //!
-    //! copy of getAllocatedSize() for use in OutputDataIterator
-    size_t size() const {
-        return getAllocatedSize();
-    }
-
-    //! return all sizes of its axes
+    //! Return all sizes of its axes
     std::vector<size_t> getAllSizes() const;
 
-    //! return copy of raw data vector
+    //! Return copy of raw data vector
     std::vector<T> getRawDataVector() const;
 
     //! fill raw array with data
@@ -87,22 +74,22 @@ template <class T> class OutputData
     friend class OutputDataIterator<T, OutputData<T> >;
     friend class OutputDataIterator<const T, const OutputData<T> >;
 
-    //! read/write iterator type
+    //! Read/write iterator type
     typedef OutputDataIterator<T, OutputData<T> > iterator;
 
-    //! read-only iterator type
+    //! Read-only iterator type
     typedef OutputDataIterator<const T, const OutputData<T> > const_iterator;
 
-    //! return a read/write iterator that points to the first element
+    //! Return a read/write iterator that points to the first element
     iterator begin();
 
-    //! return a read-only iterator that points to the first element
+    //! Return a read-only iterator that points to the first element
     const_iterator begin() const;
 
-    //! return a read/write iterator that points to the one past last element
+    //! Return a read/write iterator that points to the one past last element
     iterator end() { return iterator(this, getAllocatedSize()); }
 
-    //! return a read-only iterator that points to the one past last element
+    //! Return a read-only iterator that points to the one past last element
     const_iterator end() const  { return const_iterator(this, getAllocatedSize()); }
 
     //! get mask that will be used by iterators
@@ -114,29 +101,29 @@ template <class T> class OutputData
     //! add mask that will be used by iterators
     void addMask(const Mask &mask);
 
-    //! remove all masks
+    //! Remove all masks
     void removeAllMasks();
 
     // ---------------------------------
     // coordinate and index functions
     // ---------------------------------
 
-    //! return vector of coordinates for given index
+    //! Return vector of coordinates for given index
     std::vector<int> toCoordinates(size_t index) const;
 
-    //! return coordinate for given index and axis number
+    //! Return coordinate for given index and axis number
     int toCoordinate(size_t index, size_t i_selected_axis) const;
 
-    //! return index for specified coordinates
+    //! Return index for specified coordinates
     size_t toIndex(std::vector<int> coordinates) const;
 
-    //! return index of axis with given name for given total index
+    //! Return index of axis with given name for given total index
     size_t getIndexOfAxis(const std::string &axis_name, size_t total_index) const;
 
-    //! return value of axis with given name at given index
+    //! Return value of axis with given name at given index
     double getValueOfAxis(const std::string &axis_name, size_t index) const;
 
-    //! return bin of axis with given name and index
+    //! Return bin of axis with given name and index
     Bin1D getBinOfAxis(const std::string &axis_name, size_t index) const;
 
     // ---------
@@ -173,7 +160,7 @@ template <class T> class OutputData
     //! multiplication-assignment operator for two output data
     const OutputData<T> &operator*=(const OutputData<T> &right);
 
-   //! indexed accessor
+    //! indexed accessor
     T &operator[](size_t index) {
         if (mp_ll_data) return (*mp_ll_data)[index];
         throw ClassInitializationException("Low-level data object was not yet initialized");
@@ -189,10 +176,10 @@ template <class T> class OutputData
     // helpers
     // --------
 
-    //! return true if object have same dimensions
+    //! Return true if object have same dimensions
     bool hasSameDimensions(const OutputData<T> &right) const;
 
-    //! return true if object have same dimensions and shape of axises
+    //! Return true if object have same dimensions and shape of axises
     bool hasSameShape(const OutputData<T> &right) const;
  private:
     //! disabled copy constructor and assignment operators
@@ -256,9 +243,12 @@ template <class T> void OutputData<T>::copyFrom(const OutputData<T> &other)
 
 template <class T> void OutputData<T>::addAxis(const IAxis &new_axis)
 {
-    if( getAxis(new_axis.getName()) ) throw LogicErrorException("OutputData<T>::addAxis(const IAxis &new_axis) -> Error! Attempt to add axis with already existing name '"+ new_axis.getName()+std::string("'"));
-    if (new_axis.getSize()>0)
-    {
+    if( getAxis(new_axis.getName()) )
+        throw LogicErrorException(
+            "OutputData<T>::addAxis(const IAxis &new_axis) -> "
+            "Error! Attempt to add axis with already existing name '" +
+            new_axis.getName() + "'");
+    if (new_axis.getSize()>0) {
         m_value_axes.push_back(new_axis.clone());
         allocate();
     }
@@ -279,10 +269,8 @@ template <class T> const IAxis *OutputData<T>::getAxis(size_t index) const
 
 template <class T> const IAxis *OutputData<T>::getAxis(const std::string &label) const
 {
-    for (size_t i = 0; i < m_value_axes.size(); ++i)
-    {
-        if (m_value_axes[i]->getName() == label)
-        {
+    for (size_t i = 0; i < m_value_axes.size(); ++i) {
+        if (m_value_axes[i]->getName() == label) {
             return m_value_axes[i];
         }
     }
@@ -359,7 +347,9 @@ template <class T> void OutputData<T>::setMask(const Mask &mask)
 template <class T> void OutputData<T>::addMask(const Mask &mask)
 {
     if (mask.mp_submask) {
-        throw RuntimeErrorException("OutputData<T>::addMask() -> Error! One can only add single masks to OutputDataIterator at a time");
+        throw RuntimeErrorException(
+                    "OutputData<T>::addMask() -> "
+                    "Error! One can only add single masks to OutputDataIterator at a time");
     }
     Mask *p_old_mask = getMask();
     mp_mask = mask.clone();
@@ -380,13 +370,15 @@ template<class T> std::vector<int> OutputData<T>::toCoordinates(size_t index) co
     result.resize(mp_ll_data->getRank());
     for (size_t i=0; i<mp_ll_data->getRank(); ++i)
     {
-        result[mp_ll_data->getRank()-1-i] = (int)(remainder % m_value_axes[mp_ll_data->getRank()-1-i]->getSize());
+        result[mp_ll_data->getRank()-1-i] =
+                (int)(remainder % m_value_axes[mp_ll_data->getRank()-1-i]->getSize());
         remainder /= m_value_axes[mp_ll_data->getRank()-1-i]->getSize();
     }
     return result;
 }
 
-template<class T> int OutputData<T>::toCoordinate(size_t index, size_t i_selected_axis) const
+template<class T> int OutputData<T>::toCoordinate(
+        size_t index, size_t i_selected_axis) const
 {
     size_t remainder(index);
     for (size_t i=0; i<mp_ll_data->getRank(); ++i)
@@ -396,15 +388,16 @@ template<class T> int OutputData<T>::toCoordinate(size_t index, size_t i_selecte
         if(i_selected_axis == i_axis ) return result;
         remainder /= m_value_axes[i_axis]->getSize();
     }
-    throw LogicErrorException("OutputData<T>::toCoordinate() -> Error! No axis with given number");
+    throw LogicErrorException("OutputData<T>::toCoordinate() -> "
+                              "Error! No axis with given number");
 }
 
-
 template <class T> size_t OutputData<T>::toIndex(std::vector<int> coordinates) const
 {
-    if (coordinates.size() != mp_ll_data->getRank()) {
-        throw LogicErrorException("size_t OutputData<T>::toIndex() -> Error! Number of coordinates must match rank of data structure");
-    }
+    if (coordinates.size() != mp_ll_data->getRank())
+        throw LogicErrorException(
+                    "size_t OutputData<T>::toIndex() -> "
+                    "Error! Number of coordinates must match rank of data structure");
     size_t result = 0;
     int step_size = 1;
     for (size_t i=mp_ll_data->getRank(); i>0; --i)
@@ -415,7 +408,8 @@ template <class T> size_t OutputData<T>::toIndex(std::vector<int> coordinates) c
     return result;
 }
 
-template <class T> size_t OutputData<T>::getIndexOfAxis(const std::string &axis_name, size_t total_index) const
+template <class T> size_t OutputData<T>::getIndexOfAxis(
+        const std::string &axis_name, size_t total_index) const
 {
     std::vector<int> coordinates = toCoordinates(total_index);
     for (size_t i=0; i<m_value_axes.size(); ++i) {
@@ -423,10 +417,11 @@ template <class T> size_t OutputData<T>::getIndexOfAxis(const std::string &axis_
             return coordinates[i];
         }
     }
-    throw LogicErrorException("OutputData<T>::getIndexOfAxis() -> Error! Axis with given name not found '"+std::string(axis_name)+std::string("'"));
+    throw LogicErrorException(
+                "OutputData<T>::getIndexOfAxis() -> "
+                "Error! Axis with given name not found '" + axis_name+"'");
 }
 
-
 template <class T>
 double OutputData<T>::getValueOfAxis(const std::string &axis_name, size_t index) const
 {
@@ -436,7 +431,9 @@ double OutputData<T>::getValueOfAxis(const std::string &axis_name, size_t index)
             return (*m_value_axes[i])[axis_index];
         }
     }
-    throw LogicErrorException("OutputData<T>::getValueOfAxis() -> Error! Axis with given name not found '"+std::string(axis_name)+std::string("'"));
+    throw LogicErrorException(
+                "OutputData<T>::getValueOfAxis() -> "
+                "Error! Axis with given name not found '" + axis_name + "'");
 }
 
 template <class T>
@@ -448,7 +445,9 @@ Bin1D OutputData<T>::getBinOfAxis(const std::string &axis_name, size_t index) co
             return m_value_axes[i]->getBin(axis_index);
         }
     }
-    throw LogicErrorException("OutputData<T>::getBinOfAxis() -> Error! Axis with given name not found '"+std::string(axis_name)+std::string("'"));
+    throw LogicErrorException(
+                "OutputData<T>::getBinOfAxis() -> "
+                "Error! Axis with given name not found '" + axis_name + "'");
 }
 
 template<class T>
@@ -543,21 +542,25 @@ template<class T> inline void OutputData<T>::setRawDataArray(const T *source)
 }
 
 
-//! return true if object have same dimensions
-template<class T> inline bool OutputData<T>::hasSameDimensions(const OutputData<T> &right) const
+//! Return true if object have same dimensions
+template<class T> inline bool OutputData<T>::hasSameDimensions(
+    const OutputData<T> &right) const
 {
     if(!mp_ll_data || !right.mp_ll_data ) return false;
     return HaveSameDimensions(*mp_ll_data, *right.mp_ll_data);
 }
 
-//! return true if object have same dimensions and shape of axis
+//! Return true if object have same dimensions and shape of axis
 template<class T>
 bool OutputData<T>::hasSameShape(const OutputData<T> &right) const
 {
     if(!hasSameDimensions(right)) return false;
 
-    if( (mp_ll_data->getRank() != m_value_axes.size()) || (right.mp_ll_data->getRank() != right.m_value_axes.size()) ) {
-        throw LogicErrorException("OutputData<T>::hasSameShape() -> Panic! Inconsistent dimensions in LLData and axes");
+    if( (mp_ll_data->getRank() != m_value_axes.size()) ||
+        (right.mp_ll_data->getRank() != right.m_value_axes.size()) ) {
+        throw LogicErrorException(
+             "OutputData<T>::hasSameShape() -> "
+             "Panic! Inconsistent dimensions in LLData and axes");
     }
     for (size_t i=0; i<m_value_axes.size(); ++i) {
         const IAxis *p_axis_left = m_value_axes[i];
diff --git a/Core/Tools/inc/OutputDataIterator.h b/Core/Tools/inc/OutputDataIterator.h
index 9318fb4a0daad0ef1032c90ed8234cf40a85674c..66fe55e104a27544f29c56ffc1bed3dfd00a6c6c 100644
--- a/Core/Tools/inc/OutputDataIterator.h
+++ b/Core/Tools/inc/OutputDataIterator.h
@@ -173,7 +173,7 @@ template<class TValue, class TContainer> OutputDataIterator<TValue, TContainer>
         m_current_index = mp_mask->getNextIndex(m_current_index);
     }
     else {
-        if (m_current_index<mp_output_data->size()) {
+        if (m_current_index<mp_output_data->getAllocatedSize()) {
             ++m_current_index;
         }
     }
@@ -203,7 +203,7 @@ template<class TValue, class TContainer> void OutputDataIterator<TValue, TContai
     if (mp_mask != &mask) {
         delete mp_mask;
         mp_mask = mask.clone();
-        mp_mask->setMaxIndex(mp_output_data->size());
+        mp_mask->setMaxIndex(mp_output_data->getAllocatedSize());
     }
     m_current_index = mp_mask->getFirstValidIndex();
 }
@@ -216,7 +216,7 @@ template<class TValue, class TContainer> void OutputDataIterator<TValue, TContai
     Mask *p_old_mask = getMask();
     mp_mask = mask.clone();
     mp_mask->mp_submask = p_old_mask;
-    mp_mask->setMaxIndex(mp_output_data->size());
+    mp_mask->setMaxIndex(mp_output_data->getAllocatedSize());
     m_current_index = mp_mask->getFirstValidIndex();
 }
 
diff --git a/Core/Tools/inc/Units.h b/Core/Tools/inc/Units.h
index 0d699601024a211b0d4b65066422edc9183d5402..ace56137b9489784a7657fda3dd21cab6eeb7c6c 100644
--- a/Core/Tools/inc/Units.h
+++ b/Core/Tools/inc/Units.h
@@ -53,5 +53,4 @@ static const double PI_SQR      = 9.86960440108935861883449099987615114f;
 
 }
 
-
 #endif // UNITS_H
diff --git a/Core/Tools/src/OutputDataFunctions.cpp b/Core/Tools/src/OutputDataFunctions.cpp
index 3a1618201028b99b60b477f5b59fbb0b35b5ca48..bf06f0873d6ef29bae5667f865313e8db9bafa70 100644
--- a/Core/Tools/src/OutputDataFunctions.cpp
+++ b/Core/Tools/src/OutputDataFunctions.cpp
@@ -212,7 +212,7 @@ OutputData<double>* OutputDataFunctions::sliceAccrossOneAxis(
     const OutputData<double > &data, const std::string &fixed_axis_name,
     double fixed_axis_value)
 {
-    if (data.getNdimensions() != 2) {
+    if (data.getRank() != 2) {
         throw LogicErrorException("OutputDataFunctions::sliceAccrossOneAxis() -> Error! It was checked only with number of dimensions equal 2.");
     }
     if( !data.getAxis(fixed_axis_name) ) {
@@ -223,7 +223,7 @@ OutputData<double>* OutputDataFunctions::sliceAccrossOneAxis(
 
     const IAxis *fixed_axis(0);
     int fixed_axis_index(-1);
-    for(size_t i_axis=0; i_axis<data.getNdimensions(); i_axis++) {
+    for(size_t i_axis=0; i_axis<data.getRank(); i_axis++) {
         const IAxis *axis = data.getAxis(i_axis);
         if( axis->getName() != fixed_axis_name ) {
             sliced_data->addAxis(*axis);
@@ -261,7 +261,7 @@ OutputData<double>* OutputDataFunctions::selectRangeOnOneAxis(
     const OutputData<double > &data, const std::string &selected_axis_name,
     double axis_value1,  double axis_value2)
 {
-    if (data.getNdimensions() != 2) {
+    if (data.getRank() != 2) {
         throw LogicErrorException("OutputDataFunctions::selectRangeOnOneAxis() -> Error! It was checked only with number of dimensions equal 2.");
     }
 
@@ -282,7 +282,7 @@ OutputData<double>* OutputDataFunctions::selectRangeOnOneAxis(
 
     // preparing new data with modified axes
     OutputData<double > *new_data = new OutputData<double >;
-    for(size_t i_axis=0; i_axis<data.getNdimensions(); i_axis++) {
+    for(size_t i_axis=0; i_axis<data.getRank(); i_axis++) {
         const IAxis *axis = data.getAxis(i_axis);
         if( axis->getName() != selected_axis_name ) {
             new_data->addAxis(*axis);
diff --git a/Tests/UnitTests/TestCore/SimulationTest.h b/Tests/UnitTests/TestCore/SimulationTest.h
index b9c8ab1a68715fd8a370e8b01cf188d5078c895e..765bda383d044abcce0d72eef8e58828a389fd39 100644
--- a/Tests/UnitTests/TestCore/SimulationTest.h
+++ b/Tests/UnitTests/TestCore/SimulationTest.h
@@ -36,12 +36,10 @@ SimulationTest::SimulationTest()
     test_data.addAxis("phi_f", 10, 0., 10.);
     test_data.addAxis("theta_f", 20, 0., 20.);
     test_data.setAllTo(2.0);
-
 }
 
 SimulationTest::~SimulationTest()
 {
-
 }
 
 
@@ -49,8 +47,8 @@ TEST_F(SimulationTest, SimulationInitialState)
 {
     EXPECT_EQ( NULL, emptySimulation.getSample());
     EXPECT_EQ( size_t(1), emptySimulation.getOutputData()->getAllocatedSize());
-    EXPECT_EQ( size_t(0), emptySimulation.getOutputData()->getNdimensions());
-    EXPECT_TRUE(emptySimulation.getOutputData()->getNdimensions() == emptySimulation.getInstrument().getDetectorDimension() );
+    EXPECT_EQ( size_t(0), emptySimulation.getOutputData()->getRank());
+    EXPECT_TRUE(emptySimulation.getOutputData()->getRank() == emptySimulation.getInstrument().getDetectorDimension() );
 }
 
 
@@ -83,8 +81,8 @@ TEST_F(SimulationTest, SimulationInitialStateOfClone)
     Simulation *emptyClonedSimulation = emptySimulation.clone();
     EXPECT_EQ( NULL, emptyClonedSimulation->getSample());
     EXPECT_EQ( size_t(1), emptyClonedSimulation->getOutputData()->getAllocatedSize());
-    EXPECT_EQ( size_t(0), emptyClonedSimulation->getOutputData()->getNdimensions());
-    EXPECT_TRUE(emptyClonedSimulation->getOutputData()->getNdimensions() == emptyClonedSimulation->getInstrument().getDetector().getDimension() );
+    EXPECT_EQ( size_t(0), emptyClonedSimulation->getOutputData()->getRank());
+    EXPECT_TRUE(emptyClonedSimulation->getOutputData()->getRank() == emptyClonedSimulation->getInstrument().getDetector().getDimension() );
     EXPECT_EQ( double(1), emptyClonedSimulation->getInstrument().getIntensity());
     delete emptyClonedSimulation;
 }
@@ -107,5 +105,4 @@ TEST_F(SimulationTest, SimulationClone)
     delete clonedSimulation;
 }
 
-
 #endif // SIMULATIONTEST_H