Skip to content
Snippets Groups Projects
Commit 714007cc authored by Pospelov, Gennady's avatar Pospelov, Gennady
Browse files

Refactoring in FitObjects

parent cd26e767
No related branches found
No related tags found
No related merge requests found
...@@ -17,6 +17,7 @@ ...@@ -17,6 +17,7 @@
#define FITELEMENT_H #define FITELEMENT_H
#include "WinDllMacros.h" #include "WinDllMacros.h"
#include <cstddef>
//! @class FitElement //! @class FitElement
//! @ingroup fitting_internal //! @ingroup fitting_internal
...@@ -27,10 +28,12 @@ class BA_CORE_API_ FitElement ...@@ -27,10 +28,12 @@ class BA_CORE_API_ FitElement
{ {
public: public:
FitElement(); FitElement();
FitElement(double simul_value, double real_value); FitElement(size_t index, double simul_value, double real_value);
FitElement(const FitElement &other); FitElement(const FitElement &other);
FitElement &operator=(const FitElement &other); FitElement &operator=(const FitElement &other);
size_t getIndex() const;
double getSimulValue() const; double getSimulValue() const;
double getRealValue() const; double getRealValue() const;
...@@ -45,7 +48,7 @@ public: ...@@ -45,7 +48,7 @@ public:
private: private:
void swapContent(FitElement &other); void swapContent(FitElement &other);
size_t m_index;
double m_simul_value; double m_simul_value;
double m_real_value; double m_real_value;
double m_weight; double m_weight;
...@@ -54,6 +57,11 @@ private: ...@@ -54,6 +57,11 @@ private:
}; };
inline size_t FitElement::getIndex() const
{
return m_index;
}
inline double FitElement::getSimulValue() const inline double FitElement::getSimulValue() const
{ {
return m_simul_value; return m_simul_value;
......
...@@ -17,7 +17,8 @@ ...@@ -17,7 +17,8 @@
#include <algorithm> #include <algorithm>
FitElement::FitElement() FitElement::FitElement()
: m_simul_value(0.0) : m_index(0)
, m_simul_value(0.0)
, m_real_value(0.0) , m_real_value(0.0)
, m_weight(1.0) , m_weight(1.0)
, m_squared_difference(0.0) , m_squared_difference(0.0)
...@@ -25,8 +26,9 @@ FitElement::FitElement() ...@@ -25,8 +26,9 @@ FitElement::FitElement()
{ {
} }
FitElement::FitElement(double simul_value, double real_value) FitElement::FitElement(size_t index, double simul_value, double real_value)
: m_simul_value(simul_value) : m_index(index)
, m_simul_value(simul_value)
, m_real_value(real_value) , m_real_value(real_value)
, m_weight(1.0) , m_weight(1.0)
, m_squared_difference(0.0) , m_squared_difference(0.0)
...@@ -36,7 +38,8 @@ FitElement::FitElement(double simul_value, double real_value) ...@@ -36,7 +38,8 @@ FitElement::FitElement(double simul_value, double real_value)
FitElement::FitElement(const FitElement &other) FitElement::FitElement(const FitElement &other)
: m_simul_value(other.m_simul_value) : m_index(other.m_index)
, m_simul_value(other.m_simul_value)
, m_real_value(other.m_real_value) , m_real_value(other.m_real_value)
, m_weight(other.m_weight) , m_weight(other.m_weight)
, m_squared_difference(other.m_squared_difference) , m_squared_difference(other.m_squared_difference)
...@@ -57,6 +60,7 @@ FitElement &FitElement::operator=(const FitElement &other) ...@@ -57,6 +60,7 @@ FitElement &FitElement::operator=(const FitElement &other)
void FitElement::swapContent(FitElement &other) void FitElement::swapContent(FitElement &other)
{ {
std::swap(this->m_index, other.m_index);
std::swap(this->m_simul_value, other.m_simul_value); std::swap(this->m_simul_value, other.m_simul_value);
std::swap(this->m_real_value, other.m_real_value); std::swap(this->m_real_value, other.m_real_value);
std::swap(this->m_weight, other.m_weight); std::swap(this->m_weight, other.m_weight);
......
...@@ -101,7 +101,7 @@ class DrawObserver(IObserver): ...@@ -101,7 +101,7 @@ class DrawObserver(IObserver):
pylab.colorbar(im) pylab.colorbar(im)
pylab.title('Simulated data') pylab.title('Simulated data')
# plotting difference map # plotting difference map
#diff_map = (real_data - simulated_data)/(real_data + 1) # diff_map = (real_data - simulated_data)/(real_data + 1)
diff_map= fit_suite.getFitObjects().getChiSquaredMap().getArray() diff_map= fit_suite.getFitObjects().getChiSquaredMap().getArray()
pylab.subplot(2, 2, 3) pylab.subplot(2, 2, 3)
im = pylab.imshow(diff_map, norm=matplotlib.colors.LogNorm(), extent=[-1.0, 1.0, 0, 2.0], vmin = 0.001, vmax = 1.0) im = pylab.imshow(diff_map, norm=matplotlib.colors.LogNorm(), extent=[-1.0, 1.0, 0, 2.0], vmin = 0.001, vmax = 1.0)
...@@ -133,7 +133,7 @@ def run_fitting(): ...@@ -133,7 +133,7 @@ def run_fitting():
fit_suite = FitSuite() fit_suite = FitSuite()
# fit_suite.setMinimizer("GSLLMA") # fit_suite.setMinimizer("GSLLMA")
fit_suite.setMinimizer("Minuit2", "Fumili") # fit_suite.setMinimizer("Minuit2", "Fumili")
# Here we are settings masks to the detector plane to simulate image which looks like a Pac-Man # Here we are settings masks to the detector plane to simulate image which looks like a Pac-Man
# from ancient arcade game # from ancient arcade game
......
...@@ -52,7 +52,7 @@ class BA_CORE_API_ FitObject : public IParameterized ...@@ -52,7 +52,7 @@ class BA_CORE_API_ FitObject : public IParameterized
//! Returns simulated data //! Returns simulated data
const OutputData<double> *getSimulationData() const { const OutputData<double> *getSimulationData() const {
return m_chi2_module->getSimulationData(); } return m_simulation_data.get(); }
//! Returns chi2 module //! Returns chi2 module
...@@ -79,7 +79,11 @@ class BA_CORE_API_ FitObject : public IParameterized ...@@ -79,7 +79,11 @@ class BA_CORE_API_ FitObject : public IParameterized
//! Returns size of data //! Returns size of data
size_t getSizeOfData() const; size_t getSizeOfData() const;
std::vector<FitElement> calculateFitElements(); void calculateFitElements(std::vector<FitElement> &fit_elements);
OutputData<double> *getChiSquaredMap(std::vector<FitElement>::const_iterator first,
std::vector<FitElement>::const_iterator last) const;
protected: protected:
//! Registers some class members for later access via parameter pool //! Registers some class members for later access via parameter pool
...@@ -92,6 +96,8 @@ class BA_CORE_API_ FitObject : public IParameterized ...@@ -92,6 +96,8 @@ class BA_CORE_API_ FitObject : public IParameterized
GISASSimulation* m_simulation; //!< external simulation (not owned by this) GISASSimulation* m_simulation; //!< external simulation (not owned by this)
OutputData<double>* m_real_data; //!< real data OutputData<double>* m_real_data; //!< real data
IChiSquaredModule* m_chi2_module; //!< chi2 module IChiSquaredModule* m_chi2_module; //!< chi2 module
std::auto_ptr<OutputData<double> > m_simulation_data;
double m_weight; //!< weight of data set in chi2 calculations double m_weight; //!< weight of data set in chi2 calculations
}; };
......
...@@ -118,6 +118,14 @@ class BA_CORE_API_ FitSuite : public IObservable ...@@ -118,6 +118,14 @@ class BA_CORE_API_ FitSuite : public IObservable
//! Returns total wall time in seconds which was spend for run fit //! Returns total wall time in seconds which was spend for run fit
double getRunTime() const; double getRunTime() const;
const OutputData<double> * getRealData(size_t i_item = 0) const { return m_fit_objects.getRealData(i_item); }
const OutputData<double> * getSimulationData(size_t i_item = 0) const { return m_fit_objects.getSimulationData(i_item); }
const OutputData<double> * getChiSquaredMap(size_t i_item = 0) const { return m_fit_objects.getChiSquaredMap(i_item); }
private: private:
//! disabled copy constructor and assignment operator //! disabled copy constructor and assignment operator
FitSuite& operator=(const FitSuite& ); FitSuite& operator=(const FitSuite& );
......
...@@ -35,7 +35,7 @@ class BA_CORE_API_ FitSuiteObjects : public IParameterized ...@@ -35,7 +35,7 @@ class BA_CORE_API_ FitSuiteObjects : public IParameterized
typedef SafePointerVector<FitObject > FitObjects_t; typedef SafePointerVector<FitObject > FitObjects_t;
FitSuiteObjects(); FitSuiteObjects();
virtual ~FitSuiteObjects(){} virtual ~FitSuiteObjects();
//! clear all data //! clear all data
void clear(); void clear();
...@@ -105,8 +105,9 @@ class BA_CORE_API_ FitSuiteObjects : public IParameterized ...@@ -105,8 +105,9 @@ class BA_CORE_API_ FitSuiteObjects : public IParameterized
m_nfree_parameters = nfree_parameters; } m_nfree_parameters = nfree_parameters; }
//! Returns chi-squared map //! Returns chi-squared map
OutputData<double> * getChiSquaredMap(size_t i_item = 0); const OutputData<double> * getChiSquaredMap(size_t i_item = 0) const;
void setChiSquaredModule(const IChiSquaredModule &chi2_module);
protected: protected:
//! Registers some class members for later access via parameter pool //! Registers some class members for later access via parameter pool
...@@ -145,6 +146,8 @@ class BA_CORE_API_ FitSuiteObjects : public IParameterized ...@@ -145,6 +146,8 @@ class BA_CORE_API_ FitSuiteObjects : public IParameterized
double m_chi_squared_value; double m_chi_squared_value;
std::vector<FitElement> m_fit_elements; std::vector<FitElement> m_fit_elements;
std::auto_ptr<IChiSquaredModule> m_chi2_module;
}; };
#endif // FITSUITEKIT_H #endif // FITSUITEKIT_H
......
...@@ -107,24 +107,33 @@ size_t FitObject::getSizeOfData() const ...@@ -107,24 +107,33 @@ size_t FitObject::getSizeOfData() const
return result; return result;
} }
std::vector<FitElement> FitObject::calculateFitElements() void FitObject::calculateFitElements(std::vector<FitElement> &fit_elements)
{ {
std::vector<FitElement> result;
m_simulation->runSimulation(); m_simulation->runSimulation();
boost::scoped_ptr<OutputData<double> > sim_data(m_simulation->getIntensityData()); m_simulation_data.reset(m_simulation->getIntensityData());
const OutputData<bool> *masks(0); const OutputData<bool> *masks(0);
if(m_simulation->getInstrument().getDetector()->hasMasks()) { if(m_simulation->getInstrument().getDetector()->hasMasks()) {
masks = m_simulation->getInstrument().getDetector()->getDetectorMask()->getMaskData(); masks = m_simulation->getInstrument().getDetector()->getDetectorMask()->getMaskData();
} }
for(size_t index=0; index<sim_data->getAllocatedSize(); ++index) { for(size_t index=0; index<m_simulation_data->getAllocatedSize(); ++index) {
if(masks && (*masks)[index]) continue; if(masks && (*masks)[index]) continue;
FitElement element((*sim_data)[index], (*m_real_data)[index]); FitElement element(index, (*m_simulation_data)[index], (*m_real_data)[index]);
result.push_back(element); fit_elements.push_back(element);
}
}
OutputData<double> *FitObject::getChiSquaredMap(std::vector<FitElement>::const_iterator first,
std::vector<FitElement>::const_iterator last) const
{
OutputData<double> *result = new OutputData<double>;
result->copyShapeFrom(*m_simulation_data);
for(std::vector<FitElement>::const_iterator it=first; it!=last; ++it) {
(*result)[it->getIndex()] = it->getSquaredDifference();
} }
m_chi2_module->processFitElements(result.begin(), result.end());
return result; return result;
} }
...@@ -142,7 +142,6 @@ void FitSuite::minimize() ...@@ -142,7 +142,6 @@ void FitSuite::minimize()
IMinimizer::function_chi2_t fun_chi2 = boost::bind(&FitSuiteChiSquaredFunction::evaluate, &m_function_chi2, _1); IMinimizer::function_chi2_t fun_chi2 = boost::bind(&FitSuiteChiSquaredFunction::evaluate, &m_function_chi2, _1);
m_minimizer->setChiSquaredFunction( fun_chi2, m_fit_parameters.size()); m_minimizer->setChiSquaredFunction( fun_chi2, m_fit_parameters.size());
std::cout << "m_fit_objects.getSizeOfDataSet() " << m_fit_objects.getSizeOfDataSet() << std::endl;
IMinimizer::function_gradient_t fun_gradient = boost::bind(&FitSuiteGradientFunction::evaluate, &m_function_gradient, _1, _2, _3); IMinimizer::function_gradient_t fun_gradient = boost::bind(&FitSuiteGradientFunction::evaluate, &m_function_gradient, _1, _2, _3);
m_minimizer->setGradientFunction( fun_gradient, m_fit_parameters.size(), m_fit_objects.getSizeOfDataSet() ); m_minimizer->setGradientFunction( fun_gradient, m_fit_parameters.size(), m_fit_objects.getSizeOfDataSet() );
......
...@@ -21,11 +21,17 @@ FitSuiteObjects::FitSuiteObjects() ...@@ -21,11 +21,17 @@ FitSuiteObjects::FitSuiteObjects()
// , m_simulation_normalize(false), // , m_simulation_normalize(false),
,m_nfree_parameters(0) ,m_nfree_parameters(0)
, m_chi_squared_value(0) , m_chi_squared_value(0)
, m_chi2_module(new ChiSquaredModule())
{ {
setName("FitSuiteObjects"); setName("FitSuiteObjects");
init_parameters(); init_parameters();
} }
FitSuiteObjects::~FitSuiteObjects()
{
}
//! clear all data //! clear all data
void FitSuiteObjects::clear() void FitSuiteObjects::clear()
{ {
...@@ -50,19 +56,29 @@ void FitSuiteObjects::add( ...@@ -50,19 +56,29 @@ void FitSuiteObjects::add(
//! loop through all defined simulations and run them //! loop through all defined simulations and run them
void FitSuiteObjects::runSimulations() void FitSuiteObjects::runSimulations()
{ {
// for(FitObjects_t::iterator it =
// m_fit_objects.begin(); it!= m_fit_objects.end(); ++it) {
// (*it)->getSimulation()->runSimulation();
//// if(m_simulation_normalize) {
//// (*it)->getSimulation()->normalize();
//// }
// }
// m_chi_squared_value = calculateChiSquaredValue();
m_fit_elements.clear();
m_fit_elements.reserve(getSizeOfDataSet());
for(FitObjects_t::iterator it = for(FitObjects_t::iterator it =
m_fit_objects.begin(); it!= m_fit_objects.end(); ++it) { m_fit_objects.begin(); it!= m_fit_objects.end(); ++it) {
(*it)->getSimulation()->runSimulation(); (*it)->calculateFitElements(m_fit_elements);
// if(m_simulation_normalize) {
// (*it)->getSimulation()->normalize();
// }
} }
m_chi_squared_value = calculateChiSquaredValue();
std::cout << "aaa " << m_chi_squared_value;
for(FitObjects_t::iterator it = if(m_fit_elements.size() != getSizeOfDataSet()) {
m_fit_objects.begin(); it!= m_fit_objects.end(); ++it) { std::ostringstream message;
m_fit_elements = (*it)->calculateFitElements(); message << "FitSuiteObjects::runSimulations() -> Error. Dataset size mismatch. "
<< " m_fit_elements.size():" << m_fit_elements.size()
<< " getSizeOfDataset():" << getSizeOfDataSet() << std::endl;
throw LogicErrorException(message.str());
} }
m_chi_squared_value = calculateChiSquaredValueNew(); m_chi_squared_value = calculateChiSquaredValueNew();
...@@ -78,7 +94,6 @@ size_t FitSuiteObjects::getSizeOfDataSet() const ...@@ -78,7 +94,6 @@ size_t FitSuiteObjects::getSizeOfDataSet() const
m_fit_objects.begin(); it!= m_fit_objects.end(); ++it) { m_fit_objects.begin(); it!= m_fit_objects.end(); ++it) {
result += (*it)->getSizeOfData(); result += (*it)->getSizeOfData();
} }
std::cout << "XXX " << result << std::endl;
if(result == 0) { if(result == 0) {
throw LogicErrorException("Panic, zero dataset size"); throw LogicErrorException("Panic, zero dataset size");
} }
...@@ -131,6 +146,8 @@ double FitSuiteObjects::calculateChiSquaredValue() ...@@ -131,6 +146,8 @@ double FitSuiteObjects::calculateChiSquaredValue()
double FitSuiteObjects::calculateChiSquaredValueNew() double FitSuiteObjects::calculateChiSquaredValueNew()
{ {
m_chi2_module->processFitElements(m_fit_elements.begin(), m_fit_elements.end());
double result(0); double result(0);
for(std::vector<FitElement>::iterator it=m_fit_elements.begin(); it!=m_fit_elements.end(); ++it) { for(std::vector<FitElement>::iterator it=m_fit_elements.begin(); it!=m_fit_elements.end(); ++it) {
result += it->getSquaredDifference(); result += it->getSquaredDifference();
...@@ -146,8 +163,8 @@ double FitSuiteObjects::calculateChiSquaredValueNew() ...@@ -146,8 +163,8 @@ double FitSuiteObjects::calculateChiSquaredValueNew()
double FitSuiteObjects::getResidualValue(size_t global_index) double FitSuiteObjects::getResidualValue(size_t global_index)
{ {
if(global_index >= getSizeOfDataSet()) { if(global_index >= m_fit_elements.size()) {
throw LogicErrorException(" FitSuiteObjects::getResidualValue() -> Error. Wrong size of dataset"); throw LogicErrorException(" FitSuiteObjects::getResidualValue() -> Error. Index exceeds size of dataset.");
} }
// size_t index(0); // size_t index(0);
// const FitObject *fitObject = // const FitObject *fitObject =
...@@ -195,7 +212,21 @@ std::string FitSuiteObjects::addParametersToExternalPool( ...@@ -195,7 +212,21 @@ std::string FitSuiteObjects::addParametersToExternalPool(
return new_path; return new_path;
} }
OutputData<double> *FitSuiteObjects::getChiSquaredMap(size_t i_item) const OutputData<double> *FitSuiteObjects::getChiSquaredMap(size_t i_item) const
{
check_index(i_item);
size_t istart(0);
for(size_t i=0; i<i_item; ++i) {
std::cout << "FIXME HERE" << std::endl;
assert(0);
istart += m_fit_objects[i]->getSizeOfData();
}
std::vector<FitElement>::const_iterator start = m_fit_elements.begin() + istart;
std::vector<FitElement>::const_iterator end = start + m_fit_objects[i_item]->getSizeOfData();
return m_fit_objects[check_index(i_item)]->getChiSquaredMap(start, end);
}
void FitSuiteObjects::setChiSquaredModule(const IChiSquaredModule &chi2_module)
{ {
return m_fit_objects[check_index(i_item)]->getChiSquaredModule()->createChi2DifferenceMap(); m_chi2_module.reset(chi2_module.clone());
} }
...@@ -129,6 +129,7 @@ TEST_F(DetectorMaskTest, AssignmentOperator) ...@@ -129,6 +129,7 @@ TEST_F(DetectorMaskTest, AssignmentOperator)
EXPECT_FALSE(mask.getMask(index)); EXPECT_FALSE(mask.getMask(index));
} }
} }
EXPECT_EQ(mask.getNumberOfMaskedChannels(), 32);
} }
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment