diff --git a/App/inc/TestIsGISAXS9.h b/App/inc/TestIsGISAXS9.h
index 70699483b98a53e0ba76640e2b61cbbba75f8c0d..86260bc4665d8eb91c8d2dc6759415c219c42479 100644
--- a/App/inc/TestIsGISAXS9.h
+++ b/App/inc/TestIsGISAXS9.h
@@ -17,6 +17,7 @@
#include "IFunctionalTest.h"
#include "OutputData.h"
#include "ISample.h"
+#include "SafePointerVector.h"
//- -------------------------------------------------------------------
@@ -34,7 +35,7 @@ public:
void clear();
private:
- std::vector<OutputData<double> *> m_results;
+ SafePointerVector<OutputData<double> > m_results;
};
diff --git a/App/src/IsGISAXSTools.cpp b/App/src/IsGISAXSTools.cpp
index 3a979480e9bf842dbf6138d1045b8cd3a675444b..4310e8c220d8a3a9d7954369c17a5863218a97e4 100644
--- a/App/src/IsGISAXSTools.cpp
+++ b/App/src/IsGISAXSTools.cpp
@@ -90,7 +90,7 @@ TH2D *IsGISAXSTools::getOutputDataTH2D(const OutputData<double>& output, const s
// 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) {
- const AxisDouble *axis = output.getAxis(i_axis);
+ const IAxis *axis = output.getAxis(i_axis);
if( !axis ) throw("IsGISAXSTools::getOutputDataTH123D() -> Error! Can't cast axis");
double dx(0);
haxises[i_axis].nbins = axis->getSize();
@@ -154,7 +154,7 @@ TH1 *IsGISAXSTools::getOutputDataTH123D(const OutputData<double>& output, const
// 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) {
- const AxisDouble *axis = output.getAxis(i_axis);
+ const IAxis *axis = output.getAxis(i_axis);
if( !axis ) throw("IsGISAXSTools::getOutputDataTH123D() -> Error! Can't cast axis");
double dx(0);
haxises[i_axis].nbins = axis->getSize();
@@ -465,8 +465,8 @@ void IsGISAXSTools::exportOutputDataInVectors2D(const OutputData<double> &output
{
if (output_data.getRank() != 2) return;
- const AxisDouble *p_axis0 = output_data.getAxis(0);
- const AxisDouble *p_axis1 = output_data.getAxis(1);
+ const IAxis *p_axis0 = output_data.getAxis(0);
+ const IAxis *p_axis1 = output_data.getAxis(1);
//std::string axis0_name = p_axis0->getName();
//std::string axis1_name = p_axis1->getName();
size_t axis0_size = p_axis0->getSize();
diff --git a/App/src/TestFormFactor.cpp b/App/src/TestFormFactor.cpp
index 45df3456288963c44fe0d7b24c97814b9e64126f..acd7f3844d91fa8e9a9445f4a45762c1ee8412b3 100644
--- a/App/src/TestFormFactor.cpp
+++ b/App/src/TestFormFactor.cpp
@@ -27,8 +27,8 @@ TestFormFactor::~TestFormFactor()
void TestFormFactor::execute()
{
OutputData<double>::iterator it = mp_intensity_output->begin();
- const AxisDouble *p_y_axis = mp_intensity_output->getAxis("detector y-axis");
- const AxisDouble *p_z_axis = mp_intensity_output->getAxis("detector z-axis");
+ const IAxis *p_y_axis = mp_intensity_output->getAxis("detector y-axis");
+ const IAxis *p_z_axis = mp_intensity_output->getAxis("detector z-axis");
double lambda = 1.0;
double alpha_i = 0.2*M_PI/180.0;
cvector_t k_i;
@@ -53,8 +53,8 @@ void TestFormFactor::draw()
TCanvas *c1 = new TCanvas("c1_test_formfactor", "Cylinder Formfactor", 0, 0, 1024, 768);
(void)c1;
- const AxisDouble *p_y_axis = mp_intensity_output->getAxis("detector y-axis");
- const AxisDouble *p_z_axis = mp_intensity_output->getAxis("detector z-axis");
+ const IAxis *p_y_axis = mp_intensity_output->getAxis("detector y-axis");
+ const IAxis *p_z_axis = mp_intensity_output->getAxis("detector z-axis");
size_t y_size = p_y_axis->getSize();
size_t z_size = p_z_axis->getSize();
double y_start = (*p_y_axis)[0];
diff --git a/App/src/TestFresnelCoeff.cpp b/App/src/TestFresnelCoeff.cpp
index 2476675cb306a0f28bf5c2f08a026f1b9508d3a9..671108d6276cfd03db2e05aa99d5b8a77b0f0d4f 100644
--- a/App/src/TestFresnelCoeff.cpp
+++ b/App/src/TestFresnelCoeff.cpp
@@ -268,7 +268,7 @@ void TestFresnelCoeff::draw_roughness_set()
size_t ncoeffs = 2;
enum key_coeffs { kCoeffR, kCoeffT};
- const AxisDouble *p_rough = mp_coeffs->getAxis("roughness");
+ const IAxis *p_rough = mp_coeffs->getAxis("roughness");
size_t nroughness = p_rough->getSize();
diff --git a/App/src/TestFumiliLMA.cpp b/App/src/TestFumiliLMA.cpp
index 6ab407368921a163737af8057203d5d2eaa4232b..dc5cf8a1b3a51e52c9cb1face8c977d9e02f1922 100644
--- a/App/src/TestFumiliLMA.cpp
+++ b/App/src/TestFumiliLMA.cpp
@@ -165,8 +165,8 @@ void TestFumiliLMA::FillOutputDataFromFunction(OutputData<double> &data, TF2 *fu
data.setAllTo(0.0);
OutputData<double>::iterator it=data.begin();
- const AxisDouble *xaxis = data.getAxis(0);
- const AxisDouble *yaxis = data.getAxis(1);
+ const IAxis *xaxis = data.getAxis(0);
+ const IAxis *yaxis = data.getAxis(1);
while( it!= data.end() )
{
size_t index_x = data.toCoordinates(it.getIndex())[0];
@@ -185,8 +185,8 @@ double MyChi2Function::DataElement(const double *pars, unsigned int index, doubl
std::cout << " DataElement() -> " << g << " " << index;
for(int ipar=0; ipar<m_test->m_func->GetNpar(); ++ipar) std::cout << " " << pars[ipar];
std::cout << std::endl;
- const AxisDouble *xaxis = m_test->m_real_data->getAxis(0);
- const AxisDouble *yaxis = m_test->m_real_data->getAxis(1);
+ const IAxis *xaxis = m_test->m_real_data->getAxis(0);
+ const IAxis *yaxis = m_test->m_real_data->getAxis(1);
size_t index_x = m_test->m_real_data->toCoordinates(index)[0];
size_t index_y = m_test->m_real_data->toCoordinates(index)[1];
double x = (*xaxis)[index_x];
diff --git a/App/src/TestIsGISAXS12.cpp b/App/src/TestIsGISAXS12.cpp
index 8157ba801ca60aadbe5a729486b9128d97ca78b8..21ba23fe8eb9e1e1300690debf4c3a7fd5d189d3 100644
--- a/App/src/TestIsGISAXS12.cpp
+++ b/App/src/TestIsGISAXS12.cpp
@@ -210,7 +210,7 @@ void TestIsGISAXS12::run_isgisaxs_fit()
leg1->SetBorderSize(1);
leg1->SetFillStyle(0);
for(size_t i_set=0; i_set<m_fitSuite->getFitObjects()->size(); ++i_set) {
- c2->cd(i_set+1);
+ c2->cd((int)i_set+1);
const FitObject *obj = m_fitSuite->getFitObjects()->getObject(i_set);
TH1D *hreal = IsGISAXSTools::getOutputDataScanHist(*obj->getChiSquaredModule()->getRealData(),"gisasfw_real");
TH1D *hsimul = IsGISAXSTools::getOutputDataScanHist(*obj->getChiSquaredModule()->getSimulationData(),"gisasfw_simul");
@@ -234,7 +234,7 @@ void TestIsGISAXS12::run_isgisaxs_fit()
OutputData<double > *real = obj->getChiSquaredModule()->getRealData()->clone();
OutputData<double > *simul = obj->getChiSquaredModule()->getSimulationData()->clone();
- c2->cd(i_set+3);
+ c2->cd((int)(i_set+3));
*simul /= *real;
TH1D *hratio = IsGISAXSTools::getOutputDataScanHist(*simul,"gisasfw_real_simul_ratio");
hratio->DrawCopy();
@@ -271,7 +271,7 @@ void TestIsGISAXS12::run_test_chimodule()
OutputDataNormalizerScaleAndShift normalizer(1.31159E+05, -8.10009E-02);
double max_intensity(0);
- for(size_t i=0; i<isgi_results.size(); ++i) {
+ for(int i=0; i<(int)isgi_results.size(); ++i) {
OutputData<double>::const_iterator cit = std::max_element(isgi_results[i]->begin(), isgi_results[i]->end());
max_intensity = std::max(max_intensity, *cit);
}
@@ -281,10 +281,10 @@ void TestIsGISAXS12::run_test_chimodule()
chiModule.setOutputDataNormalizer( normalizer );
double chi_sum(0);
- for(size_t i=0; i<isgi_scans.size(); ++i) {
+ for(int i=0; i<(int)isgi_scans.size(); ++i) {
chiModule.setRealAndSimulatedData(*isgi_scans[i], *isgi_results[i]);
std::cout << " AAA " << isgi_scans.size()*isgi_results[i]->getAllocatedSize() - 12 << std::endl;
- chiModule.setNdegreeOfFreedom(isgi_scans.size()*isgi_results[i]->getAllocatedSize() - 12);
+ chiModule.setNdegreeOfFreedom((int)(isgi_scans.size()*isgi_results[i]->getAllocatedSize()) - 12);
double chi = 0.5*0.5*chiModule.calculateChiSquared();
chi_sum += chi;
std::cout << "chi : " << chi << " chi_sum:" << chi_sum << std::endl;
@@ -322,7 +322,7 @@ void TestIsGISAXS12::plot_isgisaxs_data()
leg1->SetBorderSize(1);
leg1->SetFillStyle(0);
- for(size_t i_set=0; i_set<isgi_scans.size(); ++i_set) {
+ for(int i_set=0; i_set<(int)isgi_scans.size(); ++i_set) {
c1->cd(1+i_set); gPad->SetLogy();
TH1D *hdata = IsGISAXSTools::getOutputDataScanHist(*isgi_scans[i_set], "data");
hdata->DrawCopy();
@@ -339,7 +339,7 @@ void TestIsGISAXS12::plot_isgisaxs_data()
TLegend *leg2 = new TLegend(0.5,0.6,0.85,0.85);
leg2->SetBorderSize(1);
leg2->SetFillStyle(0);
- for(size_t i_set=0; i_set<isgi_scans.size(); ++i_set) {
+ for(int i_set=0; i_set<(int)isgi_scans.size(); ++i_set) {
c1->cd(3+i_set); gPad->SetLogy();
TH1D *hdata = IsGISAXSTools::getOutputDataScanHist(*isgi_results[i_set], "data");
hdata->SetLineColor(kRed);
@@ -357,7 +357,7 @@ void TestIsGISAXS12::plot_isgisaxs_data()
TLegend *leg3 = new TLegend(0.5,0.6,0.85,0.85);
leg3->SetBorderSize(1);
leg3->SetFillStyle(0);
- for(size_t i_set=0; i_set<isgi_scans.size(); ++i_set) {
+ for(int i_set=0; i_set<(int)isgi_scans.size(); ++i_set) {
c1->cd(5+i_set);
*isgi_results[i_set] /= *isgi_scans_smoothed[i_set];
TH1D *hdata = IsGISAXSTools::getOutputDataScanHist(*isgi_results[i_set], "isgisaxs_results_data");
@@ -426,7 +426,7 @@ void TestIsGISAXS12::run_test_minimizer()
TLegend *leg1 = new TLegend(0.5,0.6,0.85,0.85);
leg1->SetBorderSize(1);
leg1->SetFillStyle(0);
- for(size_t i_set=0; i_set<isgi_results.size(); ++i_set) {
+ for(int i_set=0; i_set<(int)isgi_results.size(); ++i_set) {
c1->cd(1+i_set);
gPad->SetLogy();
TH1D *hdata = IsGISAXSTools::getOutputDataScanHist(*isgi_results[i_set], "data");
@@ -446,7 +446,7 @@ void TestIsGISAXS12::run_test_minimizer()
TLegend *leg2 = new TLegend(0.5,0.6,0.85,0.85);
leg2->SetBorderSize(1);
leg2->SetFillStyle(0);
- for(size_t i_set=0; i_set<isgi_results.size(); ++i_set) {
+ for(int i_set=0; i_set<(int)isgi_results.size(); ++i_set) {
c1->cd(3+i_set);
OutputData<double > *data = m_fitSuite->getFitObjects()->getObject(i_set)->getChiSquaredModule()->getSimulationData()->clone();
*data /= *isgi_results[i_set];
@@ -764,8 +764,8 @@ void TestIsGISAXS12::print_axes(DataSet &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[i_set]->getNdimensions(); ++i_axis) {
- const AxisDouble *axis = data[i_set]->getAxis(i_axis);
+ for(size_t i_axis=0; i_axis<data[(int)i_set]->getNdimensions(); ++i_axis) {
+ const IAxis *axis = data[(int)i_set]->getAxis(i_axis);
std::cout << "( " << axis->getName() << ", " << axis->getSize() << ", " << axis->getMin() << ", " << axis->getMax() << " ) ";
}
std::cout << std::endl;
diff --git a/App/src/TestIsGISAXS8.cpp b/App/src/TestIsGISAXS8.cpp
index afd55fc599c36efcd4d0b4254274f8a1710c6886..bf731c5e799e68e6adb3cada1fc932b7245a3357 100644
--- a/App/src/TestIsGISAXS8.cpp
+++ b/App/src/TestIsGISAXS8.cpp
@@ -17,8 +17,6 @@ TestIsGISAXS8::TestIsGISAXS8() : IFunctionalTest("TestIsGISAXS8")
void TestIsGISAXS8::execute()
{
- gsl_set_error_handler_off();
-
GISASExperiment experiment(mp_options);
experiment.setDetectorParameters(100, 0.0*Units::degree, 2.0*Units::degree, 100, 0.0*Units::degree, 2.0*Units::degree, true);
experiment.setBeamParameters(1.0*Units::angstrom, -0.2*Units::degree, 0.0*Units::degree);
diff --git a/App/src/TestIsGISAXS9.cpp b/App/src/TestIsGISAXS9.cpp
index ad7c418235761e42d18af6022af3976a66cc63e7..a5c933d00dd1718991b6d5783a779f9ce8f328ad 100644
--- a/App/src/TestIsGISAXS9.cpp
+++ b/App/src/TestIsGISAXS9.cpp
@@ -22,9 +22,6 @@ TestIsGISAXS9::~TestIsGISAXS9()
void TestIsGISAXS9::clear()
{
- for(std::vector<OutputData<double> *>::iterator it=m_results.begin(); it!=m_results.end(); ++it) {
- delete (*it);
- }
m_results.clear();
}
@@ -38,22 +35,25 @@ void TestIsGISAXS9::execute()
experiment.setDetectorParameters(100, 0.0*Units::degree, 2.0*Units::degree, 100, 0.0*Units::degree, 2.0*Units::degree, true);
experiment.setBeamParameters(1.0*Units::angstrom, -0.2*Units::degree, 0.0*Units::degree);
+ MultiLayer *p_sample;
+
// pyramid
- MultiLayer *sample = dynamic_cast<MultiLayer *>(SampleFactory::instance().createItem("IsGISAXS9_Pyramid"));
- experiment.setSample(*sample);
+ p_sample = dynamic_cast<MultiLayer *>(SampleFactory::instance().createItem("IsGISAXS9_Pyramid"));
+ experiment.setSample(*p_sample);
experiment.runSimulation();
OutputData<double> *data = experiment.getOutputDataClone();
m_results.push_back( data );
IsGISAXSTools::writeOutputDataToFile(*m_results.back(), Utils::FileSystem::GetHomePath()+"./Examples/IsGISAXS_examples/ex-9/this_pyramid_Z0.ima");
+ delete p_sample;
// rotated pyramid
- sample = dynamic_cast<MultiLayer *>(SampleFactory::instance().createItem("IsGISAXS9_RotatedPyramid"));
- experiment.setSample(*sample);
+ p_sample = dynamic_cast<MultiLayer *>(SampleFactory::instance().createItem("IsGISAXS9_RotatedPyramid"));
+ experiment.setSample(*p_sample);
experiment.runSimulation();
data = experiment.getOutputDataClone();
m_results.push_back( data );
IsGISAXSTools::writeOutputDataToFile(*m_results.back(), Utils::FileSystem::GetHomePath()+"./Examples/IsGISAXS_examples/ex-9/this_pyramid_Z45.ima");
-
+ delete p_sample;
}
diff --git a/App/src/TestMesoCrystal2.cpp b/App/src/TestMesoCrystal2.cpp
index 083b28231f99e8edd37f36b2417d12d08c1197c6..c6e2a5e20d63cec5633d430ab81214c398da09e4 100644
--- a/App/src/TestMesoCrystal2.cpp
+++ b/App/src/TestMesoCrystal2.cpp
@@ -226,8 +226,8 @@ void TestMesoCrystal2::initializeExperiment(const OutputData<double> *output_dat
mp_experiment->setDetectorParameters(200, 0.3*Units::degree, 0.073, 200, -0.4*Units::degree, 0.066);
} else {
// if there is output_data as input parameter, build detector using output_data axises
- const AxisDouble *axis0 = output_data->getAxis(0);
- const AxisDouble *axis1 = output_data->getAxis(1);
+ const IAxis *axis0 = output_data->getAxis(0);
+ const IAxis *axis1 = output_data->getAxis(1);
//std::cout << axis0->getSize() << " " << axis0->getMin() << " " << axis0->getMax() << " " << axis1->getSize() << " " << axis1->getMin() << " " << axis1->getMax() << std::endl;
mp_experiment->setDetectorParameters(axis0->getSize(), axis0->getMin(), axis0->getMax(), axis1->getSize(), axis1->getMin(), axis1->getMax());
}
diff --git a/App/src/TestMiscellaneous.cpp b/App/src/TestMiscellaneous.cpp
index c2ddc317ccc0612eeaf31f0d2c843858c1afd3bc..1add96973a7b1933b19e0cb9d0046ac9b077eef7 100644
--- a/App/src/TestMiscellaneous.cpp
+++ b/App/src/TestMiscellaneous.cpp
@@ -300,7 +300,7 @@ void TestMiscellaneous::test_DoubleToComplexInterpolatingFunction()
OpticalFresnel fresnelCalculator;
- const AxisDouble *p_alpha_axis = data_alpha->getAxis("alpha_f");
+ const IAxis *p_alpha_axis = data_alpha->getAxis("alpha_f");
std::map<double, OpticalFresnel::MultiLayerCoeff_t> fresnel_coeff_map;
for (size_t i=0; i<p_alpha_axis->getSize(); ++i) {
double angle = (*p_alpha_axis)[i];
diff --git a/App/src/TestRootTree.cpp b/App/src/TestRootTree.cpp
index a2e50f9aea6a01e5b730d536a383673d9a061f62..ffbd26dfa4e67d8361b0b9cac58aefbf752f740c 100644
--- a/App/src/TestRootTree.cpp
+++ b/App/src/TestRootTree.cpp
@@ -259,8 +259,8 @@ void TestRootTree::simple_write()
mp_data = mp_experiment->getOutputDataClone();
// accessing to scattering data
- const AxisDouble *axis0 = mp_data->getAxis(0);
- const AxisDouble *axis1 = mp_data->getAxis(1);
+ const IAxis *axis0 = mp_data->getAxis(0);
+ const IAxis *axis1 = mp_data->getAxis(1);
std::string axis0_name = axis0->getName();
std::string axis1_name = axis1->getName();
diff --git a/Core/Algorithms/inc/Detector.h b/Core/Algorithms/inc/Detector.h
index 26ad22825b35bf4064d4dd8d6d5395d95149be79..8bc3ad098c2a8b40c8b42574981a9d9d428365a2 100644
--- a/Core/Algorithms/inc/Detector.h
+++ b/Core/Algorithms/inc/Detector.h
@@ -15,7 +15,9 @@
//! @date Jun 21, 2012
#include "IDetectorResolution.h"
+#include "DetectorParameters.h"
#include "IParameterized.h"
+#include "SafePointerVector.h"
#include <vector>
@@ -33,8 +35,9 @@ public:
virtual ~Detector();
- void addAxis(const AxisDouble &axis);
- AxisDouble getAxis(size_t index) const;
+ void addAxis(const IAxis &axis);
+ void addAxis(const AxisParameters &axis_params);
+ const IAxis &getAxis(size_t index) const;
size_t getDimension() const { return m_axes.size(); }
void clear();
void setDetectorResolution(IDetectorResolution *p_detector_resolution) { delete mp_detector_resolution; mp_detector_resolution = p_detector_resolution; }
@@ -52,7 +55,10 @@ private:
//! swap function
void swapContent(Detector &other);
- std::vector<AxisDouble > m_axes;
+ //! initialize axis the way IsGISAXS does
+ void initializeAnglesIsgisaxs(AxisDouble *p_axis, const TSampledRange<double> &axis_range);
+
+ SafePointerVector<IAxis> m_axes;
IDetectorResolution *mp_detector_resolution;
};
diff --git a/Core/Algorithms/inc/DetectorParameters.h b/Core/Algorithms/inc/DetectorParameters.h
new file mode 100644
index 0000000000000000000000000000000000000000..25e9bbd2d3f7058812c0e35434d51e0f2bdaebb5
--- /dev/null
+++ b/Core/Algorithms/inc/DetectorParameters.h
@@ -0,0 +1,54 @@
+#ifndef DETECTORPARAMETERS_H_
+#define DETECTORPARAMETERS_H_
+// ********************************************************************
+// * The BornAgain project *
+// * Simulation of neutron and x-ray scattering at grazing incidence *
+// * *
+// * LICENSE AND DISCLAIMER *
+// * Lorem ipsum dolor sit amet, consectetur adipiscing elit. Mauris *
+// * eget quam orci. Quisque porta varius dui, quis posuere nibh *
+// * mollis quis. Mauris commodo rhoncus porttitor. *
+// ********************************************************************
+//! @file DetectorParameters.h
+//! @brief Definition of DetectorParameters class
+//! @author Scientific Computing Group at FRM II
+//! @date Dec 12, 2012
+
+#include "TRange.h"
+
+//- -------------------------------------------------------------------
+//! @class AxisParameters
+//! @brief Definition of AxisParameters class to store parameters
+//! for specifying an axis' data
+//- -------------------------------------------------------------------
+struct AxisParameters
+{
+ AxisParameters();
+ std::string m_name;
+ TSampledRange<double> m_range;
+ enum {
+ E_DEFAULT,
+ E_ISGISAXS,
+ } m_sample_method;
+};
+
+inline AxisParameters::AxisParameters()
+: m_name("")
+, m_range(0, 0.0, 0.0)
+, m_sample_method(AxisParameters::E_DEFAULT)
+{
+}
+
+//- -------------------------------------------------------------------
+//! @class DetectorParameters
+//! @brief Definition of DetectorParameters class to store parameters
+//! for specifying a 2D detector
+//- -------------------------------------------------------------------
+struct DetectorParameters
+{
+ AxisParameters m_phi_params;
+ AxisParameters m_alpha_params;
+};
+
+#endif /* DETECTORPARAMETERS_H_ */
+
diff --git a/Core/Algorithms/inc/GISASExperiment.h b/Core/Algorithms/inc/GISASExperiment.h
index 18deb506bfbace1ea2c18a6465beb7ff637fca4e..c2a95bc5f60ad2e557c4e4496ef53dff9dafa15f 100644
--- a/Core/Algorithms/inc/GISASExperiment.h
+++ b/Core/Algorithms/inc/GISASExperiment.h
@@ -16,6 +16,7 @@
#include "Experiment.h"
+#include "DetectorParameters.h"
#include "IResolutionFunction2D.h"
@@ -40,12 +41,16 @@ public:
void setDetectorParameters(size_t n_phi, double phi_f_min, double phi_f_max,
size_t n_alpha, double alpha_f_min, double alpha_f_max, bool isgisaxs_style=false);
- void setDetectorResolutionFunction(IResolutionFunction2D *p_resolution_function);
+ void setDetectorParameters(const DetectorParameters ¶ms);
+
+ void setDetectorResolutionFunction(IResolutionFunction2D *p_resolution_function);
void smearIntensityFromZAxisTilting();
virtual GISASExperiment *clone() const;
+ static const std::string PHI_AXIS_NAME;
+ static const std::string ALPHA_AXIS_NAME;
protected:
// hiding copy constructor and disabling assignment operator
GISASExperiment(const GISASExperiment &other);
@@ -55,13 +60,11 @@ private:
//! initialize pool parameters, i.e. register some of class members for later access via parameter pool
virtual void init_parameters();
- void initializeAnglesIsgisaxs(AxisDouble *p_axis, double start, double end, size_t size);
double getSolidAngle(size_t index) const;
double deltaAlpha(double alpha, double zeta) const;
double deltaPhi(double alpha, double phi, double zeta);
void createZetaAndProbVectors(std::vector<double> &zetas, std::vector<double> &probs, size_t nbr_zetas, double zeta_sigma);
void addToIntensityMap(double alpha, double phi, double value);
- int findClosestIndex(const AxisDouble *p_axis, double value);
};
#endif /* GISASEXPERIMENT_H_ */
diff --git a/Core/Algorithms/src/ConvolutionDetectorResolution.cpp b/Core/Algorithms/src/ConvolutionDetectorResolution.cpp
index a77ee5a546764083e17759b78477ce1bf2cb7655..ecaac74237405286f944caca8e2cce51423c1659 100644
--- a/Core/Algorithms/src/ConvolutionDetectorResolution.cpp
+++ b/Core/Algorithms/src/ConvolutionDetectorResolution.cpp
@@ -98,7 +98,7 @@ void ConvolutionDetectorResolution::apply1dConvolution(OutputData<double>* p_int
if (p_intensity_map->getRank() != 1) {
throw LogicErrorException("ConvolutionDetectorResolution::apply1dConvolution() -> Error! Number of axes for intensity map does not correspond to the dimension of the map.");
}
- const AxisDouble *p_axis = p_intensity_map->getAxis(0);
+ const IAxis *p_axis = p_intensity_map->getAxis(0);
// Construct source vector from original intensity map
std::vector<double> source_vector = p_intensity_map->getRawDataVector();
size_t data_size = source_vector.size();
@@ -131,8 +131,8 @@ void ConvolutionDetectorResolution::apply2dConvolution(OutputData<double>* p_int
if (p_intensity_map->getRank() != 2) {
throw LogicErrorException("ConvolutionDetectorResolution::apply2dConvolution() -> Error! Number of axes for intensity map does not correspond to the dimension of the map.");
}
- const AxisDouble *p_axis_1 = p_intensity_map->getAxis(0);
- const AxisDouble *p_axis_2 = p_intensity_map->getAxis(1);
+ const IAxis *p_axis_1 = p_intensity_map->getAxis(0);
+ const IAxis *p_axis_2 = p_intensity_map->getAxis(1);
size_t axis_size_1 = p_axis_1->getSize();
size_t axis_size_2 = p_axis_2->getSize();
if (axis_size_1 < 2 || axis_size_2 < 2) {
diff --git a/Core/Algorithms/src/Detector.cpp b/Core/Algorithms/src/Detector.cpp
index 133d45c58a56bf4a4e1c836127757e0c2d45f003..df8d67a39bf96182c5ee66e1037e3de4e9883322 100644
--- a/Core/Algorithms/src/Detector.cpp
+++ b/Core/Algorithms/src/Detector.cpp
@@ -1,4 +1,6 @@
#include "Detector.h"
+#include "AxisBin.h"
+#include "AxisDouble.h"
#include "Exceptions.h"
@@ -47,15 +49,39 @@ void Detector::swapContent(Detector &other)
/* ************************************************************************* */
// other methods
/* ************************************************************************* */
-void Detector::addAxis(const AxisDouble &axis)
+void Detector::addAxis(const IAxis &axis)
{
- m_axes.push_back(axis);
+ m_axes.push_back(axis.clone());
}
-AxisDouble Detector::getAxis(size_t index) const
+void Detector::addAxis(const AxisParameters &axis_params)
+{
+ IAxis *p_new_axis(0);
+ switch (axis_params.m_sample_method)
+ {
+ case AxisParameters::E_DEFAULT:
+ {
+ p_new_axis = new AxisBin(axis_params.m_name, axis_params.m_range.getNSamples(),
+ axis_params.m_range.getMin(), axis_params.m_range.getMax());
+ break;
+ }
+ case AxisParameters::E_ISGISAXS:
+ {
+ AxisDouble *p_axis = new AxisDouble(axis_params.m_name);
+ initializeAnglesIsgisaxs(p_axis, axis_params.m_range);
+ p_new_axis = p_axis;
+ break;
+ }
+ default:
+ throw RuntimeErrorException("Invalid sample method for axis.");
+ }
+ if (p_new_axis) m_axes.push_back(p_new_axis);
+}
+
+const IAxis &Detector::getAxis(size_t index) const
{
if (isCorrectAxisIndex(index)) {
- return m_axes[index];
+ return *m_axes[index];
}
throw OutOfBoundsException("Not so many axes in this detector.");
}
@@ -94,3 +120,13 @@ std::string Detector::addParametersToExternalPool(std::string path,
void Detector::init_parameters()
{
}
+
+void Detector::initializeAnglesIsgisaxs(AxisDouble* p_axis, const TSampledRange<double>& axis_range)
+{
+ double start_sin = std::sin(axis_range.getMin());
+ double end_sin = std::sin(axis_range.getMax());
+ double step = (end_sin-start_sin)/(axis_range.getNSamples()-1);
+ for(size_t i=0; i<axis_range.getNSamples(); ++i) {
+ p_axis->push_back(std::asin(start_sin + step*i));
+ }
+}
diff --git a/Core/Algorithms/src/Experiment.cpp b/Core/Algorithms/src/Experiment.cpp
index acd279d1e8d8e3f678efce2a97549d55b05b6e89..2367c8a1ac34525cb0b3404356a08d889d802878 100644
--- a/Core/Algorithms/src/Experiment.cpp
+++ b/Core/Algorithms/src/Experiment.cpp
@@ -87,11 +87,7 @@ void Experiment::normalize()
{
double incident_intensity = m_beam.getIntensity();
if (!m_is_normalized && incident_intensity!=1.0) {
- OutputData<double>::iterator it = m_intensity_map.begin();
- while (it != m_intensity_map.end()) {
- *it *= incident_intensity;
- ++it;
- }
+ m_intensity_map.scaleAll(incident_intensity);
m_is_normalized = true;
}
}
@@ -191,7 +187,7 @@ void Experiment::setDetectorParameters(const OutputData<double > &output_data)
//std::cout << "Experiment::setDetectorParameters() -> Info. Adjusting detector to have shape as in given output data" << std::endl;
m_detector.clear();
for(size_t i_axis=0; i_axis<output_data.getNdimensions(); ++i_axis) {
- const AxisDouble *axis = output_data.getAxis(i_axis);
+ const IAxis *axis = output_data.getAxis(i_axis);
m_detector.addAxis(*axis);
}
updateIntensityMapAxes();
diff --git a/Core/Algorithms/src/GISASExperiment.cpp b/Core/Algorithms/src/GISASExperiment.cpp
index d2add40bcf771fb8a2deb4eb0b693ed86e5a9acc..70235f2dfd23025b4f1e7c3178506317efaf1531 100644
--- a/Core/Algorithms/src/GISASExperiment.cpp
+++ b/Core/Algorithms/src/GISASExperiment.cpp
@@ -35,6 +35,8 @@ GISASExperiment::GISASExperiment(const GISASExperiment &other) : Experiment(othe
}
+const std::string GISASExperiment::PHI_AXIS_NAME = "phi_f";
+const std::string GISASExperiment::ALPHA_AXIS_NAME = "alpha_f";
/* ************************************************************************* */
// clone method
@@ -147,25 +149,33 @@ void GISASExperiment::normalize()
void GISASExperiment::setDetectorParameters(size_t n_phi, double phi_f_min, double phi_f_max,
size_t n_alpha, double alpha_f_min, double alpha_f_max, bool isgisaxs_style)
{
- const std::string s_phi_f("phi_f");
- const std::string s_alpha_f("alpha_f");
-
- m_detector.clear();
- AxisDouble phi_axis(s_phi_f);
- AxisDouble alpha_axis(s_alpha_f);
+ AxisParameters phi_params;
+ phi_params.m_name = PHI_AXIS_NAME;
+ phi_params.m_range = TSampledRange<double>(n_phi, phi_f_min, phi_f_max);
+ AxisParameters alpha_params;
+ alpha_params.m_name = ALPHA_AXIS_NAME;
+ alpha_params.m_range = TSampledRange<double>(n_alpha, alpha_f_min, alpha_f_max);
if (isgisaxs_style) {
- initializeAnglesIsgisaxs(&phi_axis, phi_f_min, phi_f_max, n_phi);
- initializeAnglesIsgisaxs(&alpha_axis, alpha_f_min, alpha_f_max, n_alpha);
+ phi_params.m_sample_method = AxisParameters::E_ISGISAXS;
+ alpha_params.m_sample_method = AxisParameters::E_ISGISAXS;
}
else {
- phi_axis.initElements(n_phi, phi_f_min, phi_f_max);
- alpha_axis.initElements(n_alpha, alpha_f_min, alpha_f_max);
+ phi_params.m_sample_method = AxisParameters::E_DEFAULT;
+ alpha_params.m_sample_method = AxisParameters::E_DEFAULT;
}
- m_detector.addAxis(phi_axis);
- m_detector.addAxis(alpha_axis);
- updateIntensityMapAxes();
+ DetectorParameters detector_params = { phi_params, alpha_params };
+ setDetectorParameters(detector_params);
}
+void GISASExperiment::setDetectorParameters(const DetectorParameters ¶ms)
+{
+ m_detector.clear();
+
+ m_detector.addAxis(params.m_phi_params);
+ m_detector.addAxis(params.m_alpha_params);
+
+ updateIntensityMapAxes();
+}
void GISASExperiment::setDetectorResolutionFunction(IResolutionFunction2D *p_resolution_function)
{
@@ -174,9 +184,6 @@ void GISASExperiment::setDetectorResolutionFunction(IResolutionFunction2D *p_res
void GISASExperiment::smearIntensityFromZAxisTilting()
{
- const std::string s_phi_f("phi_f");
- const std::string s_alpha_f("alpha_f");
-
size_t nbr_zetas = 5;
double zeta_sigma = 45.0*Units::degree;
std::vector<double> zetas;
@@ -187,8 +194,8 @@ void GISASExperiment::smearIntensityFromZAxisTilting()
m_intensity_map.setAllTo(0.0);
OutputData<double>::const_iterator it_clone = p_clone->begin();
while (it_clone != p_clone->end()) {
- double old_phi = p_clone->getValueOfAxis(s_phi_f, it_clone.getIndex());
- double old_alpha = p_clone->getValueOfAxis(s_alpha_f, it_clone.getIndex());
+ double old_phi = p_clone->getValueOfAxis(PHI_AXIS_NAME, it_clone.getIndex());
+ double old_alpha = p_clone->getValueOfAxis(ALPHA_AXIS_NAME, it_clone.getIndex());
for (size_t zeta_index=0; zeta_index<zetas.size(); ++zeta_index) {
double newphi = old_phi + deltaPhi(old_alpha, old_phi, zetas[zeta_index]);
double newalpha = old_alpha + deltaAlpha(old_alpha, zetas[zeta_index]);
@@ -202,27 +209,13 @@ void GISASExperiment::init_parameters()
{
}
-void GISASExperiment::initializeAnglesIsgisaxs(AxisDouble *p_axis, double start, double end, size_t size) {
- double start_sin = std::sin(start);
- double end_sin = std::sin(end);
- double step = (end_sin-start_sin)/(size-1);
- for(size_t i=0; i<size; ++i) {
- p_axis->push_back(std::asin(start_sin + step*i));
- }
- return;
-}
-
-
double GISASExperiment::getSolidAngle(size_t index) const
{
- const std::string s_phi_f("phi_f");
- const std::string s_alpha_f("alpha_f");
-
- const AxisDouble *p_alpha_axis = m_intensity_map.getAxis(s_alpha_f);
- const AxisDouble *p_phi_axis = m_intensity_map.getAxis(s_phi_f);
- size_t alpha_index = m_intensity_map.getIndexOfAxis(s_alpha_f, index);
+ const IAxis *p_alpha_axis = m_intensity_map.getAxis(ALPHA_AXIS_NAME);
+ const IAxis *p_phi_axis = m_intensity_map.getAxis(PHI_AXIS_NAME);
+ size_t alpha_index = m_intensity_map.getIndexOfAxis(ALPHA_AXIS_NAME, index);
size_t alpha_size = p_alpha_axis->getSize();
- size_t phi_index = m_intensity_map.getIndexOfAxis(s_phi_f, index);
+ size_t phi_index = m_intensity_map.getIndexOfAxis(PHI_AXIS_NAME, index);
size_t phi_size = p_phi_axis->getSize();
if (alpha_size<2 && phi_size<2) {
// Cannot determine detector cell size!
@@ -230,18 +223,18 @@ double GISASExperiment::getSolidAngle(size_t index) const
}
double dalpha(0), dphi(0);
- double alpha_f = m_intensity_map.getValueOfAxis(s_alpha_f, index);
+ double alpha_f = m_intensity_map.getValueOfAxis(ALPHA_AXIS_NAME, index);
double cos_alpha_f = std::cos(alpha_f);
if(alpha_size>1) {
if (alpha_index==0) {
- dalpha = p_alpha_axis->operator[](1) - p_alpha_axis->operator[](0);
+ dalpha = (*p_alpha_axis)[1] - (*p_alpha_axis)[0];
}
else if (alpha_index==alpha_size-1) {
- dalpha = p_alpha_axis->operator[](alpha_size-1) - p_alpha_axis->operator[](alpha_size-2);
+ dalpha = (*p_alpha_axis)[alpha_size-1] - (*p_alpha_axis)[alpha_size-2];
}
else {
- dalpha = (p_alpha_axis->operator[](alpha_index+1) - p_alpha_axis->operator[](alpha_index-1))/2.0;
+ dalpha = ((*p_alpha_axis)[alpha_index+1] - (*p_alpha_axis)[alpha_index-1])/2.0;
}
dalpha = std::abs(dalpha);
} else {
@@ -249,13 +242,13 @@ double GISASExperiment::getSolidAngle(size_t index) const
}
if(phi_size > 1) {
if (phi_index==0) {
- dphi = p_phi_axis->operator[](1) - p_phi_axis->operator[](0);
+ dphi = (*p_phi_axis)[1] - (*p_phi_axis)[0];
}
else if (phi_index==phi_size-1) {
- dphi = p_phi_axis->operator[](phi_size-1) - p_phi_axis->operator[](phi_size-2);
+ dphi = (*p_phi_axis)[phi_size-1] - (*p_phi_axis)[phi_size-2];
}
else {
- dphi = (p_phi_axis->operator[](phi_index+1) - p_phi_axis->operator[](phi_index-1))/2.0;
+ dphi = ((*p_phi_axis)[phi_index+1] - (*p_phi_axis)[phi_index-1])/2.0;
}
dphi = std::abs(dphi);
} else {
@@ -307,23 +300,23 @@ void GISASExperiment::createZetaAndProbVectors(std::vector<double>& zetas,
void GISASExperiment::addToIntensityMap(double alpha, double phi, double value)
{
- const AxisDouble *p_alpha_axis = m_intensity_map.getAxis("alpha_f");
- const AxisDouble *p_phi_axis = m_intensity_map.getAxis("phi_f");
+ const IAxis *p_alpha_axis = m_intensity_map.getAxis(ALPHA_AXIS_NAME);
+ const IAxis *p_phi_axis = m_intensity_map.getAxis(PHI_AXIS_NAME);
std::vector<int> coordinates;
- coordinates.push_back(findClosestIndex(p_alpha_axis, alpha));
- coordinates.push_back(findClosestIndex(p_phi_axis, phi));
+ coordinates.push_back((int)p_alpha_axis->findClosestIndex(alpha));
+ coordinates.push_back((int)p_phi_axis->findClosestIndex(phi));
m_intensity_map[m_intensity_map.toIndex(coordinates)] += value;
}
-int GISASExperiment::findClosestIndex(const AxisDouble *p_axis, double value)
-{
- int result = 0;
- double smallest_diff = std::abs(value-(*p_axis)[0]);
- for (size_t i=1; i<p_axis->getSize(); ++i) {
- double new_diff = std::abs(value-(*p_axis)[i]);
- if (new_diff > smallest_diff) break;
- result = (int)i;
- smallest_diff = new_diff;
- }
- return result;
-}
+//int GISASExperiment::findClosestIndex(const AxisDouble *p_axis, double value)
+//{
+// int result = 0;
+// double smallest_diff = std::abs(value-(*p_axis)[0]);
+// for (size_t i=1; i<p_axis->getSize(); ++i) {
+// double new_diff = std::abs(value-(*p_axis)[i]);
+// if (new_diff > smallest_diff) break;
+// result = (int)i;
+// smallest_diff = new_diff;
+// }
+// return result;
+//}
diff --git a/Core/Algorithms/src/MultiLayerDWBASimulation.cpp b/Core/Algorithms/src/MultiLayerDWBASimulation.cpp
index 5e9eaaea6a8dcb54428892c23266cdabfc9ff4a1..b53ef1e5ce99460c5f7d712339dfbd36e7b6ec4a 100644
--- a/Core/Algorithms/src/MultiLayerDWBASimulation.cpp
+++ b/Core/Algorithms/src/MultiLayerDWBASimulation.cpp
@@ -60,7 +60,7 @@ void MultiLayerDWBASimulation::run()
kvector_t m_ki_real(m_ki.x().real(), m_ki.y().real(), m_ki.z().real());
m_dwba_intensity.setAllTo(0.0);
- const AxisDouble *p_alpha_axis = getDWBAIntensity().getAxis("alpha_f");
+ const IAxis *p_alpha_axis = getDWBAIntensity().getAxis("alpha_f");
double lambda = 2.0*M_PI/m_ki_real.mag();
typedef std::pair<double, OpticalFresnel::MultiLayerCoeff_t > doublefresnel_t;
diff --git a/Core/Core.pro b/Core/Core.pro
index 143396d92f5d86783a82db5d875f27a368c75d7b..9f838bb2a619e5f2076c41e8824507fb83db940f 100644
--- a/Core/Core.pro
+++ b/Core/Core.pro
@@ -88,7 +88,8 @@ SOURCES += \
Samples/src/ParticleDecoration.cpp \
Samples/src/ParticleInfo.cpp \
\
- Tools/src/Axisdouble.cpp \
+ Tools/src/AxisBin.cpp \
+ Tools/src/AxisDouble.cpp \
Tools/src/Convolve.cpp \
Tools/src/CoreOptionsDescription.cpp \
Tools/src/DoubleToComplexInterpolatingFunction.cpp \
@@ -134,6 +135,7 @@ HEADERS += \
Algorithms/inc/ConvolutionDetectorResolution.h \
Algorithms/inc/DecouplingApproximationStrategy.h \
Algorithms/inc/Detector.h \
+ Algorithms/inc/DetectorParameters.h \
Algorithms/inc/DiffuseDWBASimulation.h \
Algorithms/inc/DWBADiffuseReflection.h \
Algorithms/inc/DWBASimulation.h \
@@ -223,8 +225,10 @@ HEADERS += \
Samples/inc/ParticleDecoration.h \
Samples/inc/ParticleInfo.h \
\
+ Tools/inc/AxisBin.h \
+ Tools/inc/AxisDouble.h \
Tools/inc/AttLimits.h \
- Tools/inc/Axisdouble.h \
+ Tools/inc/Bin.h \
Tools/inc/Convolve.h \
Tools/inc/Coordinate3D.h \
Tools/inc/DoubleToComplexInterpolatingFunction.h \
@@ -237,6 +241,7 @@ HEADERS += \
Tools/inc/FitSuiteObjects.h \
Tools/inc/FitSuiteParameters.h \
Tools/inc/FitSuiteStrategy.h \
+ Tools/inc/IAxis.h \
Tools/inc/IDoubleToComplexFunction.h \
Tools/inc/IFactory.h \
Tools/inc/IMinimizer.h \
diff --git a/Core/FormFactors/src/FormFactorPyramid.cpp b/Core/FormFactors/src/FormFactorPyramid.cpp
index 1914cd6034d5076b7dbb77fd5a2e27d812685e50..dffcb051a8bf65c2bfeb63c45980cbc66247f8cd 100644
--- a/Core/FormFactors/src/FormFactorPyramid.cpp
+++ b/Core/FormFactors/src/FormFactorPyramid.cpp
@@ -84,7 +84,7 @@ complex_t FormFactorPyramid::evaluate_for_q(const cvector_t &q) const
K2 = -MathFunctions::Sinc(q1)*std::exp(im*q1)*im + MathFunctions::Sinc(q2)*std::exp(-im*q2)*im;
K3 = MathFunctions::Sinc(q3)*std::exp(im*q3) + MathFunctions::Sinc(q4)*std::exp(-im*q4);
K4 = -MathFunctions::Sinc(q3)*std::exp(im*q3)*im + MathFunctions::Sinc(q4)*std::exp(-im*q4)*im;
- F = K1*cos( (qx-qy)*R ) + K2*sin( (qx-qy)*R ) - K3*cos( (qx+qy)*R ) - K4*sin( (qx+qy)*R );
+ F = K1*std::cos( (qx-qy)*R ) + K2*std::sin( (qx-qy)*R ) - K3*std::cos( (qx+qy)*R ) - K4*std::sin( (qx+qy)*R );
F = F*H/(qx*qy);
} else if(std::abs(qx) <= Numeric::double_epsilon && std::abs(qy) <= Numeric::double_epsilon) {
F=(4.*im*(-2./tga/tga - 2.*im*qz*R/tga + qz*qz*R*R - std::exp(im*H*qz) * ((-1.+im + H*qz)/tga - qz*R)*((1.+im + H*qz)/tga - qz*R) ))/std::pow(qz,3);
@@ -95,13 +95,13 @@ complex_t FormFactorPyramid::evaluate_for_q(const cvector_t &q) const
} else{
qxy=qy;
}
- F=(4.*(qxy*tga*(-(qxy*qxy*R) + qz*tga*(complex_t(0,-2) + qz*R*tga))*std::cos(qxy*R) -
- std::exp(im*H*qz)*qxy*(H*std::pow(qxy,2) - qxy*qxy*R*tga - qz*(complex_t(0,2) + H*qz)*std::pow(tga,2) +
+ F=(4.*(qxy*tga*(-(qxy*qxy*R) + qz*tga*(complex_t(0.0,-2.0) + qz*R*tga))*std::cos(qxy*R) -
+ std::exp(im*H*qz)*qxy*(H*std::pow(qxy,2) - qxy*qxy*R*tga - qz*(complex_t(0.0,2.0) + H*qz)*std::pow(tga,2) +
std::pow(qz,2)*R*std::pow(tga,3))*std::cos(qxy*(R - H/tga)) +
- tga*(std::pow(qxy,2)*(1. - complex_t(0,1)*qz*R*tga) + std::pow(qz,2)*std::pow(tga,2)*(1. + complex_t(0,1)*qz*R*tga))*
- std::sin(qxy*R) + complex_t(0,1)*std::exp(im*H*qz)*tga*
- (std::pow(qz,2)*std::pow(tga,2)*(complex_t(0,1) + H*qz - qz*R*tga) +
- std::pow(qxy,2)*(complex_t(0,1) - H*qz + qz*R*tga))*std::sin(qxy*(R - H/tga))))/
+ tga*(std::pow(qxy,2)*(1. - complex_t(0.0,1.0)*qz*R*tga) + std::pow(qz,2)*std::pow(tga,2)*(1. + complex_t(0.0,1.0)*qz*R*tga))*
+ std::sin(qxy*R) + complex_t(0.0,1.0)*std::exp(im*H*qz)*tga*
+ (std::pow(qz,2)*std::pow(tga,2)*(complex_t(0.0,1.0) + H*qz - qz*R*tga) +
+ std::pow(qxy,2)*(complex_t(0.0,1.0) - H*qz + qz*R*tga))*std::sin(qxy*(R - H/tga))))/
(qxy*std::pow(qxy - qz*tga,2)*std::pow(qxy + qz*tga,2));
}
return F;
diff --git a/Core/Tools/inc/AxisBin.h b/Core/Tools/inc/AxisBin.h
new file mode 100644
index 0000000000000000000000000000000000000000..8e2f0e2614a78bdb0bd4afd8d2a018b95b036fed
--- /dev/null
+++ b/Core/Tools/inc/AxisBin.h
@@ -0,0 +1,69 @@
+#ifndef AXISBIN_H_
+#define AXISBIN_H_
+// ********************************************************************
+// * The BornAgain project *
+// * Simulation of neutron and x-ray scattering at grazing incidence *
+// * *
+// * LICENSE AND DISCLAIMER *
+// * Lorem ipsum dolor sit amet, consectetur adipiscing elit. Mauris *
+// * eget quam orci. Quisque porta varius dui, quis posuere nibh *
+// * mollis quis. Mauris commodo rhoncus porttitor. *
+// ********************************************************************
+//! @file BinAxis.h
+//! @brief Definition of BinAxis class
+//! @author Scientific Computing Group at FRM II
+//! @date Dec 6, 2012
+
+#include "IAxis.h"
+
+#include <string>
+#include <vector>
+
+//- -------------------------------------------------------------------
+//! @class BinAxis
+//! @brief Definition of BinAxis class that stores the bins of an axis
+//- -------------------------------------------------------------------
+class AxisBin : public IAxis
+{
+public:
+ //! constructors
+ AxisBin(std::string name);
+ AxisBin(std::string name, size_t nbr_bins, double start, double end);
+
+ //! clone function
+ virtual AxisBin *clone() const;
+
+ //! create a new axis with half the number of bins
+ virtual AxisBin *createDoubleBinSize() const;
+
+ //! destructor
+ virtual ~AxisBin() {}
+
+ //! add new bin limit to the end
+ void push_back(double limit) { m_value_vector.push_back(limit); }
+
+ virtual size_t getSize() const;
+
+ //! indexed accessor retrieves midpoint of given bin
+ virtual double operator[](size_t index) const;
+
+ virtual Bin1D getBin(size_t index) const;
+
+ virtual double getMin() const { return m_value_vector.front(); }
+
+ virtual double getMax() const { return m_value_vector.back(); }
+
+ //! initialize axis bins
+ void initBins(size_t nbr_bins, double start, double end);
+
+ virtual size_t findClosestIndex(double value) const;
+
+// //! find the bin that contains the given value
+// Bin1D findMatchingBin(double value) const;
+protected:
+ virtual bool equals(const IAxis &other) const;
+private:
+ std::vector<double> m_value_vector; //!< vector containing the bin limits
+};
+
+#endif /* AXISBIN_H_ */
diff --git a/Core/Tools/inc/AxisDouble.h b/Core/Tools/inc/AxisDouble.h
index f4ee265c1619f87fa3e17beaeb46d37cd10bde14..6f49aadfe9b2f7e4b9c56980840dfda3644608e2 100644
--- a/Core/Tools/inc/AxisDouble.h
+++ b/Core/Tools/inc/AxisDouble.h
@@ -14,58 +14,56 @@
//! @author Scientific Computing Group at FRM II
//! @date Dec 4, 2012
+#include "IAxis.h"
+
#include <string>
#include <vector>
+// Forward declaration of BinAxis class, as needed for conversion constructor
+class AxisBin;
+
//- -------------------------------------------------------------------
//! @class AxisDouble
//! @brief Definition of AxisDouble class that stores the points of an axis
//- -------------------------------------------------------------------
-class AxisDouble
+class AxisDouble : public IAxis
{
public:
//! constructors
AxisDouble(std::string name);
AxisDouble(std::string name, size_t size, double start, double end);
- //! clone function
- AxisDouble* clone() const;
+ //! explicit conversion from BinAxis
+ //TODO: make explicit
+ AxisDouble(const AxisBin &source);
- //! create a new axis with half the number of bins
- AxisDouble createDoubleBinSize() const;
-
- //! destructor
- ~AxisDouble() {}
+ virtual AxisDouble *clone() const;
- //! retrieve the number of bins
- inline size_t getSize() const { return m_value_vector.size(); }
+ virtual AxisDouble *createDoubleBinSize() const;
- //! retrieve the label of the axis
- inline std::string getName() const { return m_name; }
+ //! destructor
+ virtual ~AxisDouble() {}
- //! set the axis label
- inline void setName(std::string name) { m_name = name; }
+ virtual size_t getSize() const { return m_sample_vector.size(); }
//! add point to the end of the axis
- inline void push_back(double element) { m_value_vector.push_back(element); }
+ void push_back(double element) { m_sample_vector.push_back(element); }
- //! indexed accessor
- inline double& operator[](size_t index) { return m_value_vector.at(index); }
+ virtual double operator[](size_t index) const { return m_sample_vector[index]; }
- //! const indexed accessor
- inline const double& operator[](size_t index) const { return m_value_vector.at(index); }
+ virtual Bin1D getBin(size_t index) const;
//! get value of first point of axis
- inline double getMin() const { return m_value_vector.front(); }
+ virtual double getMin() const { return m_sample_vector.front(); }
//! get value of last point of axis
- inline double getMax() const { return m_value_vector.back(); }
+ virtual double getMax() const { return m_sample_vector.back(); }
//! initialize axis points
void initElements(size_t size, double start, double end);
//! find number of bin which is closest to given value
- size_t findClosestIndex(double value) const;
+ virtual size_t findClosestIndex(double value) const;
//! find the index that corresponds to the given lower bound (index is inclusive)
size_t getLowerBoundIndex(double value) const;
@@ -73,13 +71,9 @@ public:
//! find the index that corresponds to the given upper bound (index is inclusive)
size_t getUpperBoundIndex(double value) const;
private:
- std::string m_name; //!< axis label
- std::vector<double> m_value_vector; //!< vector containing the axis points
+ virtual bool equals(const IAxis &other) const;
+ std::vector<double> m_sample_vector; //!< vector containing the axis points
+ double m_bin_size;
};
-//! global helper function for comparison of axes
-bool HaveSameNameAndShape(const AxisDouble &left, const AxisDouble &right);
-
-
-
#endif /* AXISDOUBLE_H_ */
diff --git a/Core/Tools/inc/Bin.h b/Core/Tools/inc/Bin.h
new file mode 100644
index 0000000000000000000000000000000000000000..655d6c4bd1a1e39bb98b269306d592cc23cc07c2
--- /dev/null
+++ b/Core/Tools/inc/Bin.h
@@ -0,0 +1,46 @@
+#ifndef BIN_H_
+#define BIN_H_
+// ********************************************************************
+// * The BornAgain project *
+// * Simulation of neutron and x-ray scattering at grazing incidence *
+// * *
+// * LICENSE AND DISCLAIMER *
+// * Lorem ipsum dolor sit amet, consectetur adipiscing elit. Mauris *
+// * eget quam orci. Quisque porta varius dui, quis posuere nibh *
+// * mollis quis. Mauris commodo rhoncus porttitor. *
+// ********************************************************************
+//! @file Bin.h
+//! @brief Definition of Bin classes
+//! @author Scientific Computing Group at FRM II
+//! @date Dec 10, 2012
+
+#include "Numeric.h"
+
+#include <cmath>
+
+//- -------------------------------------------------------------------
+//! @class Bin1D
+//! @brief Definition of Bin1D class that stores the bounds of a
+//! one-dimensional bin
+//- -------------------------------------------------------------------
+struct Bin1D
+{
+ double m_lower; //!< lower bound of the bin
+ double m_upper; //!< upper bound of the bin
+ double getMidPoint() const { return (m_lower + m_upper)/2.0; }
+};
+
+//! equality operator for bins
+inline bool operator==(const Bin1D &left, const Bin1D &right)
+{
+ if (std::abs(left.m_lower - right.m_lower) > Numeric::double_epsilon) return false;
+ if (std::abs(left.m_upper - right.m_upper) > Numeric::double_epsilon) return false;
+ return true;
+}
+
+//! inequality operator for bins
+inline bool operator!=(const Bin1D &left, const Bin1D &right) {
+ return !(left==right);
+}
+
+#endif /* BIN_H_ */
diff --git a/Core/Tools/inc/IAxis.h b/Core/Tools/inc/IAxis.h
new file mode 100644
index 0000000000000000000000000000000000000000..ced92e9c13e7c5a32852a0c5bbd2e89ffa6ab37e
--- /dev/null
+++ b/Core/Tools/inc/IAxis.h
@@ -0,0 +1,96 @@
+#ifndef IAXIS_H_
+#define IAXIS_H_
+// ********************************************************************
+// * The BornAgain project *
+// * Simulation of neutron and x-ray scattering at grazing incidence *
+// * *
+// * LICENSE AND DISCLAIMER *
+// * Lorem ipsum dolor sit amet, consectetur adipiscing elit. Mauris *
+// * eget quam orci. Quisque porta varius dui, quis posuere nibh *
+// * mollis quis. Mauris commodo rhoncus porttitor. *
+// ********************************************************************
+//! @file IAxis.h
+//! @brief Definition of IAxis interface
+//! @author Scientific Computing Group at FRM II
+//! @date Dec 10, 2012
+
+#include "Bin.h"
+
+#include <string>
+
+//- -------------------------------------------------------------------
+//! @class IAxis
+//! @brief Definition of IAxis interface for 1d axes
+//- -------------------------------------------------------------------
+
+class IAxis
+{
+public:
+ //! constructors
+ IAxis(std::string name) : m_name(name) {}
+
+ //! clone function
+ virtual IAxis *clone() const=0;
+
+ //! create a new axis with half the number of bins
+ virtual IAxis *createDoubleBinSize() const=0;
+
+ //! destructor
+ virtual ~IAxis() {}
+
+ //! retrieve the number of bins
+ virtual size_t getSize() const=0;
+
+ //! retrieve the label of the axis
+ std::string getName() const { return m_name; }
+
+ //! set the axis label
+ void setName(std::string name) { m_name = name; }
+
+ //! indexed accessor retrieves a sample
+ virtual double operator[](size_t index) const=0;
+
+ //! retrieve a 1d bin for the given index
+ virtual Bin1D getBin(size_t index) const=0;
+
+ //! get value of first point of axis
+ virtual double getMin() const=0;
+
+ //! get value of last point of axis
+ virtual double getMax() const=0;
+
+ //! find bin index which is best match for given value
+ virtual size_t findClosestIndex(double value) const=0;
+
+ //! test for equality
+ friend bool operator==(const IAxis &left, const IAxis &right) {
+ return left.equals(right);
+ }
+
+// //! find the index that corresponds to the given lower bound (index is inclusive)
+// virtual size_t getLowerBoundIndex(double value) const=0;
+//
+// //! find the index that corresponds to the given upper bound (index is inclusive)
+// virtual size_t getUpperBoundIndex(double value) const=0;
+protected:
+ virtual bool equals(const IAxis &other) const;
+ std::string m_name; //!< axis label
+};
+
+inline bool IAxis::equals(const IAxis& other) const
+{
+ return getName()==other.getName();
+}
+
+//! test for inequality
+inline bool operator!=(const IAxis &left, const IAxis &right) {
+ return !(left == right);
+}
+
+//! global helper function for comparison of axes
+inline bool HaveSameNameAndShape(const IAxis &left, const IAxis &right)
+{
+ return left == right;
+}
+
+#endif /* IAXIS_H_ */
diff --git a/Core/Tools/inc/IMinimizer.h b/Core/Tools/inc/IMinimizer.h
index b25f3cabaa354a1834a5bffee965f4b2a342c85a..7df17f8678ff098e121528e3af44c6b68a56e952 100644
--- a/Core/Tools/inc/IMinimizer.h
+++ b/Core/Tools/inc/IMinimizer.h
@@ -110,7 +110,7 @@ public:
//! return pointer to the parameters values at the minimum
virtual double getValueOfVariableAtMinimum(size_t i) const
{
- std::map<int, double >::const_iterator pos = m_values.find(i);
+ std::map<int, double >::const_iterator pos = m_values.find((int)i);
if(pos != m_values.end()){
return pos->second;
} else {
diff --git a/Core/Tools/inc/OutputData.h b/Core/Tools/inc/OutputData.h
index 37fdf160bc9272d93f1d2c568453fb247f719bf9..e942ff832590c535cb71c8c3bf5a21aa001016ca 100644
--- a/Core/Tools/inc/OutputData.h
+++ b/Core/Tools/inc/OutputData.h
@@ -19,6 +19,7 @@
#include "Types.h"
#include "LLData.h"
#include "OutputDataIterator.h"
+#include "SafePointerVector.h"
#include <string>
#include <sstream>
@@ -32,17 +33,17 @@ template <class T> class OutputData
{
public:
OutputData();
- virtual ~OutputData();
+ ~OutputData();
//! make object clone
OutputData* clone() const;
void copyFrom(const OutputData<T> &x);
- void addAxis(const AxisDouble &new_axis);
+ void addAxis(const IAxis &new_axis);
void addAxis(std::string name, size_t size, double start, double end);
- const AxisDouble *getAxis(size_t index) const;
- const AxisDouble *getAxis(std::string label) const;
+ const IAxis *getAxis(size_t index) const;
+ const IAxis *getAxis(std::string label) const;
size_t getAxisIndex(const std::string &label) const;
// ---------------------------------
@@ -199,7 +200,7 @@ private:
//! memory allocation for current dimensions configuration
void allocate();
- std::vector<AxisDouble> m_value_axes;
+ SafePointerVector<IAxis> m_value_axes;
LLData<T> *mp_ll_data;
Mask *mp_mask;
};
@@ -250,12 +251,12 @@ template <class T> void OutputData<T>::copyFrom(const OutputData<T> &other)
-template <class T> void OutputData<T>::addAxis(const AxisDouble &new_axis)
+template <class T> void OutputData<T>::addAxis(const IAxis &new_axis)
{
- if( getAxis(new_axis.getName()) ) throw LogicErrorException("OutputData<T>::addAxis(AxisDouble *) -> Error! Attempt to add axis with already existing name '"+ new_axis.getName()+std::string("'"));
+ 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)
{
- m_value_axes.push_back(new_axis);
+ m_value_axes.push_back(new_axis.clone());
allocate();
}
}
@@ -268,18 +269,18 @@ void OutputData<T>::addAxis(std::string name, size_t size, double start, double
addAxis(new_axis);
}
-template <class T> const AxisDouble *OutputData<T>::getAxis(size_t index) const
+template <class T> const IAxis *OutputData<T>::getAxis(size_t index) const
{
- return &m_value_axes.at(index);
+ return m_value_axes[index];
}
-template <class T> const AxisDouble *OutputData<T>::getAxis(std::string label) const
+template <class T> const IAxis *OutputData<T>::getAxis(std::string label) const
{
- for (std::vector<AxisDouble>::const_iterator it = m_value_axes.begin(); it != m_value_axes.end(); ++it)
+ for (size_t i = 0; i < m_value_axes.size(); ++i)
{
- if (it->getName() == label)
+ if (m_value_axes[i]->getName() == label)
{
- return &(*it);
+ return m_value_axes[i];
}
}
return 0;
@@ -288,10 +289,9 @@ template <class T> const AxisDouble *OutputData<T>::getAxis(std::string label) c
// return index of axis
template <class T> size_t OutputData<T>::getAxisIndex(const std::string &label) const
{
- size_t index(0);
- for (std::vector<AxisDouble>::const_iterator it = m_value_axes.begin(); it != m_value_axes.end(); ++it, ++index)
+ for (size_t i = 0; i < m_value_axes.size(); ++i)
{
- if (it->getName() == label) return index;
+ if (m_value_axes[i]->getName() == label) return i;
}
throw LogicErrorException("OutputData<T>::getIndexOfAxis() -> Error! Axis with given name not found '"+label+std::string("'"));
}
@@ -377,8 +377,8 @@ 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());
- 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;
}
@@ -389,9 +389,9 @@ template<class T> int OutputData<T>::toCoordinate(size_t index, size_t i_selecte
for (size_t i=0; i<mp_ll_data->getRank(); ++i)
{
size_t i_axis = mp_ll_data->getRank()-1-i;
- int result = (int)(remainder % m_value_axes[i_axis].getSize());
+ int result = (int)(remainder % m_value_axes[i_axis]->getSize());
if(i_selected_axis == i_axis ) return result;
- remainder /= m_value_axes[i_axis].getSize();
+ remainder /= m_value_axes[i_axis]->getSize();
}
throw LogicErrorException("OutputData<T>::toCoordinate() -> Error! No axis with given number");
}
@@ -407,7 +407,7 @@ template <class T> size_t OutputData<T>::toIndex(std::vector<int> coordinates) c
for (size_t i=mp_ll_data->getRank(); i>0; --i)
{
result += coordinates[i-1]*step_size;
- step_size *= m_value_axes[i-1].getSize();
+ step_size *= m_value_axes[i-1]->getSize();
}
return result;
}
@@ -416,34 +416,21 @@ template <class T> size_t OutputData<T>::getIndexOfAxis(const std::string &axis_
{
std::vector<int> coordinates = toCoordinates(total_index);
for (size_t i=0; i<m_value_axes.size(); ++i) {
- if (m_value_axes[i].getName() == axis_name) {
+ if (m_value_axes[i]->getName() == axis_name) {
return coordinates[i];
}
}
throw LogicErrorException("OutputData<T>::getIndexOfAxis() -> Error! Axis with given name not found '"+std::string(axis_name)+std::string("'"));
}
-//template <class T>
-//double OutputData<T>::getValueOfAxis(const char *axis_name, size_t index) const
-////double OutputData<T>::getValueOfAxis(std::string axis_name, size_t index) const
-//{
-// std::vector<int> coordinates = toCoordinates(index);
-// for (size_t i=0; i<m_value_axes.size(); ++i) {
-// if (m_value_axes[i].getName() == axis_name) {
-// return m_value_axes[i][coordinates[i]];
-// }
-// }
-// throw LogicErrorException("OutputData<T>::getValueOfAxis() -> Error! Axis with given name not found '"+std::string(axis_name)+std::string("'"));
-//}
-
template <class T>
double OutputData<T>::getValueOfAxis(const std::string &axis_name, size_t index) const
{
for (size_t i=0; i<m_value_axes.size(); ++i) {
- if (m_value_axes[i].getName() == axis_name) {
+ if (m_value_axes[i]->getName() == axis_name) {
int axis_index = toCoordinate(index, i);
- return m_value_axes[i][axis_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("'"));
@@ -559,10 +546,10 @@ bool OutputData<T>::hasSameShape(const OutputData<T> &right) const
throw LogicErrorException("OutputData<T>::hasSameShape() -> Panic! Inconsistent dimensions in LLData and axes");
}
for (size_t i=0; i<m_value_axes.size(); ++i) {
- const AxisDouble axis_left = m_value_axes[i];
- const AxisDouble axis_right = right.m_value_axes[i];
+ const IAxis *p_axis_left = m_value_axes[i];
+ const IAxis *p_axis_right = right.m_value_axes[i];
- if( !HaveSameNameAndShape(axis_left, axis_right)) return false;
+ if( !HaveSameNameAndShape(*p_axis_left, *p_axis_right)) return false;
}
return true;
}
diff --git a/Core/Tools/inc/SafePointerVector.h b/Core/Tools/inc/SafePointerVector.h
new file mode 100644
index 0000000000000000000000000000000000000000..bfdc59626dc1d7d86fca410128100ac4e9e2e7a8
--- /dev/null
+++ b/Core/Tools/inc/SafePointerVector.h
@@ -0,0 +1,111 @@
+#ifndef SAFEPOINTERVECTOR_H_
+#define SAFEPOINTERVECTOR_H_
+// ********************************************************************
+// * The BornAgain project *
+// * Simulation of neutron and x-ray scattering at grazing incidence *
+// * *
+// * LICENSE AND DISCLAIMER *
+// * Lorem ipsum dolor sit amet, consectetur adipiscing elit. Mauris *
+// * eget quam orci. Quisque porta varius dui, quis posuere nibh *
+// * mollis quis. Mauris commodo rhoncus porttitor. *
+// ********************************************************************
+//! @file SafePointerVector.h
+//! @brief Definition of SafePointerVector template
+//! @author Scientific Computing Group at FRM II
+//! @date Dec 10, 2012
+
+#include <vector>
+
+//- -------------------------------------------------------------------
+//! @class SafePointerVector
+//! @brief Definition of SafePointerVector<T> template for safe handling
+//! of vectors of pointers that are owned by the vector
+//- -------------------------------------------------------------------
+template <class T> class SafePointerVector
+{
+public:
+ typedef typename std::vector<T *>::iterator iterator;
+ typedef typename std::vector<T *>::const_iterator const_iterator;
+ SafePointerVector();
+ SafePointerVector(const SafePointerVector &other);
+ ~SafePointerVector();
+
+ SafePointerVector &operator=(const SafePointerVector &right);
+ size_t size() const;
+ void push_back(T *pointer);
+ T *operator[](size_t index);
+ const T *operator[](size_t index) const;
+ iterator begin() { return m_pointers.begin(); }
+ const_iterator begin() const { return m_pointers.begin(); }
+ iterator end() { return m_pointers.end(); }
+ const_iterator end() const { return m_pointers.end(); }
+
+ T *back() { return m_pointers.back(); }
+ T *const back() const { return m_pointers.back(); }
+ void clear();
+private:
+ std::vector<T *> m_pointers;
+};
+
+template<class T> SafePointerVector<T>::SafePointerVector()
+{
+}
+
+template<class T> SafePointerVector<T>::SafePointerVector(const SafePointerVector<T> &other)
+{
+ for (const_iterator it = other.begin();
+ it != other.end(); ++it) {
+ m_pointers.push_back((*it)->clone());
+ }
+}
+
+template<class T> SafePointerVector<T>::~SafePointerVector()
+{
+ clear();
+}
+
+template<class T> SafePointerVector<T> &SafePointerVector<T>::operator=(
+ const SafePointerVector<T> &right)
+{
+ if (this == &right) return *this;
+ clear();
+ for (const_iterator it = right.begin();
+ it != right.end(); ++it) {
+ m_pointers.push_back((*it)->clone());
+ }
+ return *this;
+}
+
+template<class T>
+inline size_t SafePointerVector<T>::size() const
+{
+ return m_pointers.size();
+}
+
+template<class T>
+inline void SafePointerVector<T>::push_back(T* pointer)
+{
+ m_pointers.push_back(pointer);
+}
+
+template<class T>
+inline T* SafePointerVector<T>::operator[](size_t index)
+{
+ return m_pointers[index];
+}
+
+template<class T>
+inline const T* SafePointerVector<T>::operator[](size_t index) const
+{
+ return m_pointers[index];
+}
+
+template<class T> void SafePointerVector<T>::clear()
+{
+ for (iterator it = begin(); it != end(); ++it) {
+ delete (*it);
+ }
+ m_pointers.clear();
+}
+
+#endif /* SAFEPOINTERVECTOR_H_ */
diff --git a/Core/Tools/inc/TRange.h b/Core/Tools/inc/TRange.h
index 684c24fd6ffb633cde815a002738525766949519..6ab79366235bcc54a1ca13cb2f9acadf43188191 100644
--- a/Core/Tools/inc/TRange.h
+++ b/Core/Tools/inc/TRange.h
@@ -28,4 +28,14 @@ private:
T m_min, m_max;
};
+template <class T> class TSampledRange : public TRange<T>
+{
+public:
+ TSampledRange(size_t n_samples, T min, T max) : TRange<T>(min, max), m_n_samples(n_samples) {}
+
+ size_t getNSamples() const { return m_n_samples; }
+private:
+ size_t m_n_samples;
+};
+
#endif /* TRANGE_H_ */
diff --git a/Core/Tools/src/AxisBin.cpp b/Core/Tools/src/AxisBin.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..e9fb04e868ae434ac6cab0aec7ef93ec6ab471e2
--- /dev/null
+++ b/Core/Tools/src/AxisBin.cpp
@@ -0,0 +1,97 @@
+#include "AxisBin.h"
+#include "Numeric.h"
+#include "Exceptions.h"
+
+AxisBin::AxisBin(std::string name)
+: IAxis(name)
+{
+}
+
+AxisBin::AxisBin(std::string name, size_t nbr_bins, double start, double end)
+: IAxis(name)
+{
+ initBins(nbr_bins, start, end);
+}
+
+AxisBin *AxisBin::clone() const
+{
+ AxisBin *p_clone = new AxisBin(getName());
+ p_clone->m_value_vector = this->m_value_vector;
+ return p_clone;
+}
+
+AxisBin *AxisBin::createDoubleBinSize() const
+{
+ if (getSize() < 2) {
+ return clone();
+ }
+ AxisBin *p_result = new AxisBin(getName());
+ for (size_t source_index=0; source_index<getSize(); source_index+=2)
+ {
+ p_result->push_back(m_value_vector[source_index]);
+ }
+ p_result->push_back(m_value_vector[getSize()]);
+ return p_result;
+}
+
+size_t AxisBin::getSize() const
+{
+ if (m_value_vector.size()<2) {
+ return 0;
+ }
+ return m_value_vector.size() - 1;
+}
+
+double AxisBin::operator[](size_t index) const
+{
+ return getBin(index).getMidPoint();
+}
+
+Bin1D AxisBin::getBin(size_t index) const
+{
+ Bin1D result = { m_value_vector[index], m_value_vector[index+1] };
+ return result;
+}
+
+void AxisBin::initBins(size_t nbr_bins, double start, double end)
+{
+ double step = (end - start)/(nbr_bins);
+ for (size_t i=0; i<nbr_bins+1; ++i)
+ {
+ push_back(start + step*(int)i);
+ }
+}
+
+size_t AxisBin::findClosestIndex(double value) const
+{
+ if(m_value_vector.size()<2) {
+ throw ClassInitializationException("AxisBin not (yet) correctly initialized");
+ }
+ if (value < getMin() || value > getMax()) {
+ throw OutOfBoundsException("Given value not in any bin");
+ }
+ std::vector<double>::const_iterator top_limit = std::lower_bound(m_value_vector.begin(), m_value_vector.end(), value);
+ if(top_limit != m_value_vector.begin() ) --top_limit;
+ size_t nbin = top_limit - m_value_vector.begin();
+ return nbin;
+}
+
+//Bin1D AxisBin::findMatchingBin(double value) const
+//{
+// return (*this)[findMatchingBinIndex(value)];
+//}
+//
+bool AxisBin::equals(const IAxis& other) const
+{
+ if (!IAxis::equals(other)) return false;
+ if (const AxisBin *p_other_cast = dynamic_cast<const AxisBin *>(&other)) {
+ if (getSize() != p_other_cast->getSize()) return false;
+ for(size_t i=0; i<m_value_vector.size(); ++i) {
+ if( std::abs(m_value_vector[i] - p_other_cast->m_value_vector[i]) > Numeric::double_epsilon ) {
+ return false;
+ }
+ }
+ return true;
+ }
+ return false;
+}
diff --git a/Core/Tools/src/AxisDouble.cpp b/Core/Tools/src/AxisDouble.cpp
index 16b12f965ce4a007d6a989af1ca943e5f7601127..729fa2963c5aa668a84bbb49f0ba2bc2dece6cf1 100644
--- a/Core/Tools/src/AxisDouble.cpp
+++ b/Core/Tools/src/AxisDouble.cpp
@@ -1,45 +1,62 @@
#include "AxisDouble.h"
+#include "AxisBin.h"
#include "Numeric.h"
#include "Exceptions.h"
-#include <cmath>
-#include <algorithm>
-
AxisDouble::AxisDouble(std::string name)
-: m_name(name)
+: IAxis(name)
+, m_bin_size(0)
{
}
AxisDouble::AxisDouble(std::string name, size_t size, double start, double end)
-: m_name(name)
+: IAxis(name)
+, m_bin_size(0)
{
initElements(size, start, end);
}
-AxisDouble* AxisDouble::clone() const
+AxisDouble::AxisDouble(const AxisBin& source)
+: IAxis(source.getName())
+, m_bin_size(0)
+{
+ for (size_t i=0; i<source.getSize(); ++i) {
+ m_sample_vector.push_back(source[i]);
+ }
+}
+
+AxisDouble *AxisDouble::clone() const
{
AxisDouble *p_clone = new AxisDouble(getName());
- p_clone->m_value_vector = this->m_value_vector;
+ p_clone->m_sample_vector = this->m_sample_vector;
+ p_clone->m_bin_size = this->m_bin_size;
return p_clone;
}
-AxisDouble AxisDouble::createDoubleBinSize() const
+AxisDouble *AxisDouble::createDoubleBinSize() const
{
if (getSize() < 2) {
- return *this;
+ return clone();
}
- AxisDouble result(getName());
+ AxisDouble *p_result = new AxisDouble(getName());
for (size_t source_index=0; source_index<getSize(); source_index+=2)
{
double value;
if (source_index==getSize()-1) {
- value = (3.0*m_value_vector.at(source_index) - m_value_vector.at(source_index-1))/2.0;
+ value = (3.0*m_sample_vector.at(source_index) - m_sample_vector.at(source_index-1))/2.0;
}
else {
- value = (m_value_vector.at(source_index) + m_value_vector.at(source_index+1))/2.0;
+ value = (m_sample_vector.at(source_index) + m_sample_vector.at(source_index+1))/2.0;
}
- result.push_back(value);
+ p_result->push_back(value);
}
+ p_result->m_bin_size = m_bin_size*2.0;
+ return p_result;
+}
+
+Bin1D AxisDouble::getBin(size_t index) const
+{
+ Bin1D result = { m_sample_vector[index] - m_bin_size/2.0, m_sample_vector[index] + m_bin_size/2.0 };
return result;
}
@@ -54,38 +71,45 @@ void AxisDouble::initElements(size_t size, double start, double end)
size_t AxisDouble::findClosestIndex(double value) const
{
- if(m_value_vector.size()<2) return 0;
- std::vector<double>::const_iterator before = std::lower_bound(m_value_vector.begin(), m_value_vector.end(), value);
- if(before == m_value_vector.end() ) --before;
- if(before == m_value_vector.begin() ) ++before;
+ if(m_sample_vector.size()<2) return 0;
+ std::vector<double>::const_iterator before = std::lower_bound(m_sample_vector.begin(), m_sample_vector.end(), value);
+ if(before == m_sample_vector.end() ) --before;
+ if(before == m_sample_vector.begin() ) ++before;
std::vector<double>::const_iterator after = before;
--before;
size_t nbin(0);
- ( *after-value) < (value - *before) ? nbin = std::distance(m_value_vector.begin(), after) : nbin = std::distance(m_value_vector.begin(), before);
+ ( *after-value) < (value - *before) ? nbin = std::distance(m_sample_vector.begin(), after) : nbin = std::distance(m_sample_vector.begin(), before);
return nbin;
}
size_t AxisDouble::getLowerBoundIndex(double value) const
{
- if(m_value_vector.size()<2) return 0;
- std::vector<double>::const_iterator lbound_it = std::lower_bound(m_value_vector.begin(), m_value_vector.end(), value);
- if(lbound_it == m_value_vector.end() ) {
+ if(m_sample_vector.size()<2) return 0;
+ std::vector<double>::const_iterator lbound_it = std::lower_bound(m_sample_vector.begin(), m_sample_vector.end(), value);
+ if(lbound_it == m_sample_vector.end() ) {
throw RuntimeErrorException("Given lower bound higher than highest present value");
}
- return lbound_it - m_value_vector.begin();
+ return lbound_it - m_sample_vector.begin();
}
size_t AxisDouble::getUpperBoundIndex(double value) const
{
- if(m_value_vector.size()<2) return 0;
- std::vector<double>::const_iterator lbound_it = std::upper_bound(m_value_vector.begin(), m_value_vector.end(), value);
- return (lbound_it - m_value_vector.begin()) - 1;
+ if(m_sample_vector.size()<2) return 0;
+ std::vector<double>::const_iterator lbound_it = std::upper_bound(m_sample_vector.begin(), m_sample_vector.end(), value);
+ return (lbound_it - m_sample_vector.begin()) - 1;
}
-bool HaveSameNameAndShape(const AxisDouble& left, const AxisDouble& right)
+bool AxisDouble::equals(const IAxis& other) const
{
- if(left.getSize() != right.getSize()) return false;
- if(left.getName() != right.getName()) return false;
- for(size_t i=0; i<left.getSize(); ++i) if( std::abs(left[i] - right[i]) > Numeric::double_epsilon) return false;
- return true;
+ if (!IAxis::equals(other)) return false;
+ if (const AxisDouble *p_other_cast = dynamic_cast<const AxisDouble *>(&other)) {
+ if (getSize() != p_other_cast->getSize()) return false;
+ for(size_t i=0; i<m_sample_vector.size(); ++i) {
+ if( std::abs(m_sample_vector[i] - p_other_cast->m_sample_vector[i]) > Numeric::double_epsilon ) {
+ return false;
+ }
+ }
+ return true;
+ }
+ return false;
}
diff --git a/Core/Tools/src/FitSuite.cpp b/Core/Tools/src/FitSuite.cpp
index 6ef0468a6d8ba0a8554e7e730c6ab9016cf07ab8..350278e42c05b223b068ad26c8bff6c90e5f2ef2 100644
--- a/Core/Tools/src/FitSuite.cpp
+++ b/Core/Tools/src/FitSuite.cpp
@@ -168,7 +168,7 @@ double FitSuite::functionToMinimize(const double *pars_current_values)
m_fit_objects.runSimulation();
// caclulate chi2 value
- double chi_squared = m_fit_objects.getChiSquaredValue(m_fit_parameters.getNfreeParameters());
+ double chi_squared = m_fit_objects.getChiSquaredValue((int)m_fit_parameters.getNfreeParameters());
notifyObservers();
m_n_call++;
diff --git a/Core/Tools/src/FitSuiteObjects.cpp b/Core/Tools/src/FitSuiteObjects.cpp
index 90a01a1f9b92b191de042d59027e5d7c69ce7864..393e7043f22997f4e00a6ec070059c916b2fd92b 100644
--- a/Core/Tools/src/FitSuiteObjects.cpp
+++ b/Core/Tools/src/FitSuiteObjects.cpp
@@ -52,7 +52,7 @@ double FitSuiteObjects::getChiSquaredValue(int n_free_fit_parameters)
for(FitObjects_t::iterator it = m_fit_objects.begin(); it!= m_fit_objects.end(); ++it) {
IChiSquaredModule *chi = (*it)->getChiSquaredModule();
- chi->setNdegreeOfFreedom( m_fit_objects.size() * (*it)->getRealData()->getAllocatedSize() - n_free_fit_parameters);
+ chi->setNdegreeOfFreedom( (int)(m_fit_objects.size() * (*it)->getRealData()->getAllocatedSize() - n_free_fit_parameters) );
// normalizing datasets to the maximum intensity over all fit objects defined
OutputDataNormalizerScaleAndShift *data_normalizer = dynamic_cast<OutputDataNormalizerScaleAndShift *>(chi->getOutputDataNormalizer());
if( data_normalizer) data_normalizer->setMaximumIntensity( max_intensity );
diff --git a/Core/Tools/src/OutputDataFunctions.cpp b/Core/Tools/src/OutputDataFunctions.cpp
index 6e2e36819b9f3d0b34745b19c22a97243d8b7bb6..06dff54e45d29bb1b24155055e6ca3aa133999b8 100644
--- a/Core/Tools/src/OutputDataFunctions.cpp
+++ b/Core/Tools/src/OutputDataFunctions.cpp
@@ -23,8 +23,10 @@ OutputData<double> *OutputDataFunctions::doubleBinSize(const OutputData<double>
// create new axes
for (size_t i=0; i<dimension; ++i) {
needs_resizing.push_back(source_sizes[i] > 1);
- const AxisDouble *source_axis = source.getAxis(i);
- p_result->addAxis(source_axis->createDoubleBinSize());
+ const IAxis *source_axis = source.getAxis(i);
+ IAxis *p_new_axis = source_axis->createDoubleBinSize();
+ p_result->addAxis(*p_new_axis);
+ delete p_new_axis;
}
// calculate new data content
OutputData<double>::const_iterator it_source = source.begin();
@@ -210,10 +212,10 @@ OutputData<double> *OutputDataFunctions::sliceAccrossOneAxis(const OutputData<do
OutputData<double > *sliced_data = new OutputData<double >;
- const AxisDouble *fixed_axis(0);
+ const IAxis *fixed_axis(0);
int fixed_axis_index(-1);
for(size_t i_axis=0; i_axis<data.getNdimensions(); i_axis++) {
- const AxisDouble *axis = data.getAxis(i_axis);
+ const IAxis *axis = data.getAxis(i_axis);
if( axis->getName() != fixed_axis_name ) {
sliced_data->addAxis(*axis);
} else {
@@ -252,7 +254,7 @@ OutputData<double> *OutputDataFunctions::selectRangeOnOneAxis(const OutputData<d
throw LogicErrorException("OutputDataFunctions::selectRangeOnOneAxis() -> Error! It was checked only with number of dimensions equal 2.");
}
- const AxisDouble *selected_axis = data.getAxis(selected_axis_name);
+ const IAxis *selected_axis = data.getAxis(selected_axis_name);
if( !selected_axis ) {
throw LogicErrorException("OutputDataFunctions::selectRangeOnOneAxis() -> Error! No axis with name "+selected_axis_name);
}
@@ -270,7 +272,7 @@ OutputData<double> *OutputDataFunctions::selectRangeOnOneAxis(const OutputData<d
// 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++) {
- const AxisDouble *axis = data.getAxis(i_axis);
+ const IAxis *axis = data.getAxis(i_axis);
if( axis->getName() != selected_axis_name ) {
new_data->addAxis(*axis);
} else {
@@ -343,9 +345,9 @@ Mask* OutputDataFunctions::CreateRectangularMask(const OutputData<double>& data,
int *maxima_i = new int[rank];
int *dims_i = new int[rank];
for (size_t i=0; i<rank; ++i) {
- const AxisDouble *p_axis = data.getAxis(i);
- minima_i[i] = (int)p_axis->getLowerBoundIndex(minima[i]);
- maxima_i[i] = (int)p_axis->getUpperBoundIndex(maxima[i]);
+ const IAxis *p_axis = data.getAxis(i);
+ minima_i[i] = (int)p_axis->findClosestIndex(minima[i]);
+ maxima_i[i] = (int)p_axis->findClosestIndex(maxima[i]);
dims_i[i] = (int)p_axis->getSize();
}
MaskCoordinateRectangleFunction *p_rectangle_function = new MaskCoordinateRectangleFunction(rank, minima_i, maxima_i);
@@ -366,9 +368,9 @@ Mask* OutputDataFunctions::CreateEllipticMask(const OutputData<double>& data,
int *radii_i = new int[rank];
int *dims_i = new int[rank];
for (size_t i=0; i<rank; ++i) {
- const AxisDouble *p_axis = data.getAxis(i);
- center_i[i] = (int)p_axis->getLowerBoundIndex(center[i]);
- int lower_index = (int)p_axis->getLowerBoundIndex((*p_axis)[center_i[i]] - radii[i]);
+ const IAxis *p_axis = data.getAxis(i);
+ center_i[i] = (int)p_axis->findClosestIndex(center[i]);
+ int lower_index = (int)p_axis->findClosestIndex((*p_axis)[center_i[i]] - radii[i]);
radii_i[i] = center_i[i] - lower_index;
dims_i[i] = (int)p_axis->getSize();
}
diff --git a/UnitTests/TestCore/DetectorTest.h b/UnitTests/TestCore/DetectorTest.h
index 89be1d35b359d6b96280d46b18c7ea77742b40a8..b6d5def59ce449be9d248e596289b66b28f925b3 100644
--- a/UnitTests/TestCore/DetectorTest.h
+++ b/UnitTests/TestCore/DetectorTest.h
@@ -57,9 +57,10 @@ TEST_F(DetectorTest, DetectorConstruction)
AxisDouble axis1("axis1", 20, 0.0, 20.0);
constructedDetector.addAxis(axis0);
constructedDetector.addAxis(axis1);
+
EXPECT_EQ((size_t)2, constructedDetector.getDimension());
- AxisDouble axis0copy = constructedDetector.getAxis(0);
- AxisDouble axis1copy = constructedDetector.getAxis(1);
+ const IAxis &axis0copy = constructedDetector.getAxis(0);
+ const IAxis &axis1copy = constructedDetector.getAxis(1);
ASSERT_TRUE(axis0.getMin() == axis0copy.getMin() );
ASSERT_TRUE(axis0.getMax() == axis0copy.getMax() );
ASSERT_TRUE(axis0.getName() == axis0copy.getName() );
diff --git a/XCode_GISASFW.xcodeproj/project.pbxproj b/XCode_GISASFW.xcodeproj/project.pbxproj
index f6af84fb7eb7ba043463d43d540c2d6785e3ae36..2c9fc394500a5d38f3a4e4a1d7a155a547c7a154 100644
--- a/XCode_GISASFW.xcodeproj/project.pbxproj
+++ b/XCode_GISASFW.xcodeproj/project.pbxproj
@@ -7,6 +7,7 @@
objects = {
/* Begin PBXBuildFile section */
+ 62027ED516775A7700A0D41C /* AxisBin.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 62027ED416775A7700A0D41C /* AxisBin.cpp */; };
6204B7991667A408003BF358 /* IOutputDataNormalizer.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 6204B7981667A407003BF358 /* IOutputDataNormalizer.cpp */; };
6218B466161B2577007FFA5C /* FormFactorParallelepiped.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 6218B465161B2577007FFA5C /* FormFactorParallelepiped.cpp */; };
6218B469161B25BC007FFA5C /* TestIsGISAXS11.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 6218B468161B25BC007FFA5C /* TestIsGISAXS11.cpp */; };
@@ -235,6 +236,11 @@
/* End PBXCopyFilesBuildPhase section */
/* Begin PBXFileReference section */
+ 62027ECF16775A6B00A0D41C /* AxisBin.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = AxisBin.h; sourceTree = "<group>"; };
+ 62027ED016775A6B00A0D41C /* Bin.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = Bin.h; sourceTree = "<group>"; };
+ 62027ED116775A6B00A0D41C /* IAxis.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = IAxis.h; sourceTree = "<group>"; };
+ 62027ED216775A6B00A0D41C /* SafePointerVector.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = SafePointerVector.h; sourceTree = "<group>"; };
+ 62027ED416775A7700A0D41C /* AxisBin.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = AxisBin.cpp; sourceTree = "<group>"; };
6204B7961667A3FC003BF358 /* IIntensityFunction.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = IIntensityFunction.h; sourceTree = "<group>"; };
6204B7971667A3FC003BF358 /* IOutputDataNormalizer.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = IOutputDataNormalizer.h; sourceTree = "<group>"; };
6204B7981667A407003BF358 /* IOutputDataNormalizer.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = IOutputDataNormalizer.cpp; sourceTree = "<group>"; };
@@ -1675,6 +1681,10 @@
622222F2160CB745008205AC /* inc */ = {
isa = PBXGroup;
children = (
+ 62027ECF16775A6B00A0D41C /* AxisBin.h */,
+ 62027ED016775A6B00A0D41C /* Bin.h */,
+ 62027ED116775A6B00A0D41C /* IAxis.h */,
+ 62027ED216775A6B00A0D41C /* SafePointerVector.h */,
62D7A934166E3F81009A45CE /* AxisDouble.h */,
62CC097C166504D600A5B720 /* AttLimits.h */,
62CC097D166504D600A5B720 /* FitObject.h */,
@@ -1725,6 +1735,7 @@
6222230A160CB745008205AC /* src */ = {
isa = PBXGroup;
children = (
+ 62027ED416775A7700A0D41C /* AxisBin.cpp */,
62D7A935166E3F8D009A45CE /* AxisDouble.cpp */,
62CC0986166504E700A5B720 /* FitObject.cpp */,
62CC0987166504E700A5B720 /* FitParameterLinked.cpp */,
@@ -68683,6 +68694,7 @@
6254C2651666652E0098EE7E /* IFormFactorBorn.cpp in Sources */,
6204B7991667A408003BF358 /* IOutputDataNormalizer.cpp in Sources */,
62D7A936166E3F8D009A45CE /* AxisDouble.cpp in Sources */,
+ 62027ED516775A7700A0D41C /* AxisBin.cpp in Sources */,
);
runOnlyForDeploymentPostprocessing = 0;
};
diff --git a/shared.pri b/shared.pri
index 8a97ab6dd4c5ed36cb6417b46122f7a933ede75d..305feea9f38f82a1ba591b7de6b398d68f830c6b 100644
--- a/shared.pri
+++ b/shared.pri
@@ -105,8 +105,7 @@ CONFIG(JCNS) {
QMAKE_CXXFLAGS_DEBUG += -fdiagnostics-show-option # to find out in gcc which option control warning
#QMAKE_CXXFLAGS_RELEASE += -O3 -ffast-math -msse3
QMAKE_CXXFLAGS_RELEASE -= -O2
-QMAKE_CXXFLAGS_RELEASE += -g -O3
-#QMAKE_CXXFLAGS_RELEASE += -g -O3
+QMAKE_CXXFLAGS_RELEASE += -g -O3 # -ffast-math removed because of problems with NaNs
# uncommenting line below produces non-stripped (very large) libraries
#QMAKE_STRIP=: