diff --git a/App/App.pro b/App/App.pro index e0ca2ee5acfa01092bb13a861038427c90702976..62ab0ae08a72193843cea5c4aa867f8e2bc4bec2 100644 --- a/App/App.pro +++ b/App/App.pro @@ -15,7 +15,6 @@ SOURCES += \ src/AppOptionsDescription.cpp \ src/CommandLine.cpp \ src/DrawHelper.cpp \ - src/FitSuiteHelper.cpp \ src/FunctionalTestFactory.cpp \ src/IFunctionalTest.cpp \ src/IsGISAXSData.cpp \ @@ -57,7 +56,10 @@ SOURCES += \ src/TestToyExperiment.cpp \ src/TreeEventStructure.cpp \ src/ROOTGSLNLSMinimizer.cpp \ - $${FUNCTIONAL_TESTS}/IsGISAXS01/IsGISAXS01.cpp + $${FUNCTIONAL_TESTS}/IsGISAXS01/IsGISAXS01.cpp \ + src/ROOTGSLSimAnMinimizer.cpp \ + src/FitSuiteObserverFactory.cpp \ + src/MinimizerFactory.cpp HEADERS += \ inc/App.h \ @@ -65,7 +67,6 @@ HEADERS += \ inc/AppOptionsDescription.h \ inc/CommandLine.h \ inc/DrawHelper.h \ - inc/FitSuiteHelper.h \ inc/FunctionalTestFactory.h \ inc/IFunctionalTest.h \ inc/IsGISAXSData.h \ @@ -107,7 +108,10 @@ HEADERS += \ inc/TreeEventStructure.h \ inc/ROOTMinimizerFunction.h \ inc/ROOTGSLNLSMinimizer.h \ - $${FUNCTIONAL_TESTS}/IsGISAXS01/IsGISAXS01.h + $${FUNCTIONAL_TESTS}/IsGISAXS01/IsGISAXS01.h \ + inc/ROOTGSLSimAnMinimizer.h \ + inc/FitSuiteObserverFactory.h \ + inc/MinimizerFactory.h # additional locations LOCATIONS = ./inc $${FUNCTIONAL_TESTS}/IsGISAXS01 diff --git a/App/inc/FitSuiteHelper.h b/App/inc/FitSuiteObserverFactory.h similarity index 82% rename from App/inc/FitSuiteHelper.h rename to App/inc/FitSuiteObserverFactory.h index f26ec6b98599ddfe7db585126fa6f5a909391f3b..f711ffbc7b72fcd3e1b768078a0aca71bd0aa296 100644 --- a/App/inc/FitSuiteHelper.h +++ b/App/inc/FitSuiteObserverFactory.h @@ -23,6 +23,7 @@ #include <map> #include <time.h> #include <sys/time.h> +#include <boost/shared_ptr.hpp> #include "TH1.h" class TPaveText; @@ -87,6 +88,26 @@ private: }; +class FitSuiteObserverFactory +{ +public: + typedef boost::shared_ptr<FitSuiteObserverPrint > observer_print_t; + typedef boost::shared_ptr<FitSuiteObserverDraw > observer_draw_t; + typedef boost::shared_ptr<FitSuiteObserverWriteTree > observer_tree_t; + + static observer_print_t createPrintObserver(int run_every_nth=20) { + return observer_print_t( new FitSuiteObserverPrint(run_every_nth) ); + } + + static observer_draw_t createDrawObserver(int run_every_nth=20) { + return observer_draw_t( new FitSuiteObserverDraw(run_every_nth) ); + } + + static observer_tree_t createTreeObserver() { + return observer_tree_t( new FitSuiteObserverWriteTree() ); + } +}; + #endif // FITSUITEHELPER_H diff --git a/App/inc/MinimizerFactory.h b/App/inc/MinimizerFactory.h new file mode 100644 index 0000000000000000000000000000000000000000..35146d3d6a7c4367c8fa792cd99c37943c5076ab --- /dev/null +++ b/App/inc/MinimizerFactory.h @@ -0,0 +1,53 @@ +#ifndef MINIMIZERFACTORY_H +#define MINIMIZERFACTORY_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 MinimizerFactory.h +//! @brief Definition of SampleFactory class +//! @author Scientific Computing Group at FRM II +//! @date 23.01.2013 + +#include "IMinimizer.h" +#include <string> +#include <vector> +#include <map> + + +//- ------------------------------------------------------------------- +//! @class MinimizerFactory +//! @brief Factory to create minimizers +//- ------------------------------------------------------------------- +class MinimizerFactory +{ +public: + static IMinimizer *createMinimizer(const std::string &minimizer, const std::string &algorithm = std::string(), const std::string &options=std::string() ); + static void print_catalogue(); + +private: + + //! @class map of minimizer names holding list of defined algorithms for every minimizer + class Catalogue { + public: + typedef std::map<std::string, std::vector<std::string > > catalogue_t; + typedef catalogue_t::const_iterator const_iterator; + Catalogue(); + const_iterator begin() const { return m_data.begin(); } + const_iterator end() const { return m_data.end(); } + bool isValid(const std::string &minimizer, const std::string &algorithm) const; + friend std::ostream &operator<<(std::ostream &ostr, const Catalogue &m) { m.print(ostr); return ostr; } + private: + void print(std::ostream &ostr) const; + catalogue_t m_data; + }; + + static Catalogue m_catalogue; +}; + +#endif // MINIMIZERFACTORY_H diff --git a/App/inc/ROOTGSLSimAnMinimizer.h b/App/inc/ROOTGSLSimAnMinimizer.h new file mode 100644 index 0000000000000000000000000000000000000000..9b7734d2813c3b2384571dd4f84da9c9730c057a --- /dev/null +++ b/App/inc/ROOTGSLSimAnMinimizer.h @@ -0,0 +1,219 @@ +#ifndef ROOTGSLSIMANMINIMIZER_H +#define ROOTGSLSIMANMINIMIZER_H + +// @(#)root/mathmore:$Id: GSLSimAnMinimizer.h 32583 2010-03-12 09:57:42Z moneta $ +// Author: L. Moneta Wed Dec 20 17:16:32 2006 + +/********************************************************************** + * * + * Copyright (c) 2006 LCG ROOT Math Team, CERN/PH-SFT * + * * + * This library is free software; you can redistribute it and/or * + * modify it under the terms of the GNU General Public License * + * as published by the Free Software Foundation; either version 2 * + * of the License, or (at your option) any later version. * + * * + * This library is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * + * General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this library (see file COPYING); if not, write * + * to the Free Software Foundation, Inc., 59 Temple Place, Suite * + * 330, Boston, MA 02111-1307 USA, or contact the author. * + * * + **********************************************************************/ + +// Header file for class GSLSimAnMinimizer + + +#ifndef ROOT_Math_Minimizer +#include "Math/Minimizer.h" +#endif + + +#ifndef ROOT_Math_IFunctionfwd +#include "Math/IFunctionfwd.h" +#endif + +#ifndef ROOT_Math_IParamFunctionfwd +#include "Math/IParamFunctionfwd.h" +#endif + + + +#ifndef ROOT_Math_GSLSimAnnealing +#include "Math/GSLSimAnnealing.h" +#endif + +#ifndef ROOT_Math_MinimizerVariable +#include "Math/MinimizerVariable.h" +#endif + +#include <vector> +#include <map> + +#include "Fit/Fitter.h" +#include "Math/MinimTransformFunction.h" +#include "Math/MultiNumGradFunction.h" +#include "Math/MinimizerVariable.h" + + +namespace ROOT { + + namespace Patch { + //typedef ROOT::Math::MinimTransformFunction MinimTransformFunction; + + using ROOT::Math::kDefault; + using ROOT::Math::kLowBound; + using ROOT::Math::kUpBound; + using ROOT::Math::kBounds; + using ROOT::Math::kFix; + using ROOT::Math::MultiNumGradFunction; + using ROOT::Math::MinimTransformFunction; + +// class MinimTransformFunction; + + +//_____________________________________________________________________________________ +/** + GSLSimAnMinimizer class for minimization using simulated annealing + using the algorithm from + <A HREF="http://www.gnu.org/software/gsl/manual/html_node/Simulated-Annealing.html"> + GSL</A>. + It implements the ROOT::Minimizer interface and + a plug-in (name "GSLSimAn") exists to instantiate this class via the plug-in manager + + @ingroup MultiMin +*/ +class GSLSimAnMinimizer : public ROOT::Math::Minimizer { + +public: + + /** + Default constructor + */ + GSLSimAnMinimizer (int type = 0); + + /** + Destructor (no operations) + */ + ~GSLSimAnMinimizer (); + +private: + // usually copying is non trivial, so we make this unaccessible + + /** + Copy constructor + */ + GSLSimAnMinimizer(const GSLSimAnMinimizer &) : ROOT::Math::Minimizer() {} + + /** + Assignment operator + */ + GSLSimAnMinimizer & operator = (const GSLSimAnMinimizer & rhs) { + if (this == &rhs) return *this; // time saving self-test + return *this; + } + +public: + + /// set the function to minimize + virtual void SetFunction(const ROOT::Math::IMultiGenFunction & func); + + /// set gradient the function to minimize + virtual void SetFunction(const ROOT::Math::IMultiGradFunction & func); + + /// set free variable + virtual bool SetVariable(unsigned int ivar, const std::string & name, double val, double step); + + /// set fixed variable (override if minimizer supports them ) + virtual bool SetFixedVariable(unsigned int /* ivar */, const std::string & /* name */, double /* val */); + + /// set lower limit variable (override if minimizer supports them ) + virtual bool SetLowerLimitedVariable(unsigned int ivar , const std::string & name , double val , double step , double lower ); + /// set upper limit variable (override if minimizer supports them ) + virtual bool SetUpperLimitedVariable(unsigned int ivar , const std::string & name , double val , double step , double upper ); + /// set upper/lower limited variable (override if minimizer supports them ) + virtual bool SetLimitedVariable(unsigned int ivar , const std::string & name , double val , double step , double /* lower */, double /* upper */); + + /// set the value of an existing variable + virtual bool SetVariableValue(unsigned int ivar, double val ); + /// set the values of all existing variables (array must be dimensioned to the size of the existing parameters) + virtual bool SetVariableValues(const double * x); + + + /// method to perform the minimization + virtual bool Minimize(); + + /// return minimum function value + virtual double MinValue() const { return fMinVal; } + + /// return expected distance reached from the minimum + virtual double Edm() const { return 0; } // not impl. } + + /// return pointer to X values at the minimum + virtual const double * X() const { return &fValues.front(); } + + /// return pointer to gradient values at the minimum + virtual const double * MinGradient() const { return 0; } // not impl. + + /// number of function calls to reach the minimum + virtual unsigned int NCalls() const { return 0; } // not yet ipl. + + /// this is <= Function().NDim() which is the total + /// number of variables (free+ constrained ones) + virtual unsigned int NDim() const { return fDim; } + + /// number of free variables (real dimension of the problem) + /// this is <= Function().NDim() which is the total + virtual unsigned int NFree() const { return fDim; } + + /// minimizer provides error and error matrix + virtual bool ProvidesError() const { return false; } + + /// return errors at the minimum + virtual const double * Errors() const { return 0; } + + /** return covariance matrices elements + if the variable is fixed the matrix is zero + The ordering of the variables is the same as in errors + */ + virtual double CovMatrix(unsigned int , unsigned int ) const { return 0; } + + + /// return reference to the objective function + ///virtual const ROOT::Math::IGenFunction & Function() const; + + ROOT::Math::GSLSimAnnealing &getSolver() { return fSolver; } + + //void SetAnnealingParameters(const GSLSimAnParams &newpars) { fSolver.Params() + +protected: + +private: + + unsigned int fDim; // dimension of the function to be minimized + bool fOwnFunc; // flag to indicate if objective function is managed + + ROOT::Math::GSLSimAnnealing fSolver; + const ROOT::Math::IMultiGenFunction * fObjFunc; + + double fMinVal; // minimum values + + mutable std::vector<double> fValues; + + std::vector<double> fSteps; + std::vector<std::string> fNames; + std::vector<ROOT::Math::EMinimVariableType> fVarTypes; // vector specifyng the type of variables + std::map< unsigned int, std::pair<double, double> > fBounds; // map specifying the bound using as key the parameter index + +}; + + } // end namespace Math + +} // end namespace ROOT + + +#endif // ROOTGSLSIMANMINIMIZER_H diff --git a/App/inc/ROOTMinimizer.h b/App/inc/ROOTMinimizer.h index c610b715fe0780384cb61e612a93635318728c15..05066e0d89b748a85f08571fdc4927253b69a68b 100644 --- a/App/inc/ROOTMinimizer.h +++ b/App/inc/ROOTMinimizer.h @@ -64,9 +64,6 @@ public: //! return created minimizer ROOT::Math::Minimizer *getROOTMinimizer() { return m_root_minimizer; } - //! checking validity of the combination minimizer_name and algo_type - bool isValidNames(const std::string &minimizer_name, const std::string &algo_type); - //! check if type of algorithm is Levenberg-Marquardt or similar bool isGradientBasedAgorithm(); diff --git a/App/inc/TestToyExperiment.h b/App/inc/TestToyExperiment.h index e2ceb4172b5a4a30659f1664f2dadbcfde5cce83..b2a4511e35f7e2649d3795fba4a72b9b740fdd90 100644 --- a/App/inc/TestToyExperiment.h +++ b/App/inc/TestToyExperiment.h @@ -34,7 +34,6 @@ public: ToyExperiment(TF2 *func) : m_func(func) { pars.resize(func->GetNpar(), 0.0); setName("ToyExperiment"); init_parameters(); } virtual ~ToyExperiment() {} virtual void runSimulation(); - virtual void runSimulationElement(size_t index); virtual ToyExperiment *clone() const { return new ToyExperiment(*this); } void setParameter(size_t i, double value) { pars[i] = value; } protected: diff --git a/App/src/FitSuiteHelper.cpp b/App/src/FitSuiteObserverFactory.cpp similarity index 99% rename from App/src/FitSuiteHelper.cpp rename to App/src/FitSuiteObserverFactory.cpp index efb64cd7cef1067c693d624180d22f2194074c35..935658883d2603aefce167e0df01e5d1c937b861 100644 --- a/App/src/FitSuiteHelper.cpp +++ b/App/src/FitSuiteObserverFactory.cpp @@ -1,8 +1,7 @@ #include "FitSuite.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "TreeEventStructure.h" #include "IsGISAXSTools.h" -#include "ROOTMinimizer.h" #include "FitSuiteParameters.h" #include "ExperimentConstants.h" diff --git a/App/src/MinimizerFactory.cpp b/App/src/MinimizerFactory.cpp new file mode 100644 index 0000000000000000000000000000000000000000..765986e91a59de1f940f218cb54042472015f197 --- /dev/null +++ b/App/src/MinimizerFactory.cpp @@ -0,0 +1,88 @@ +#include "MinimizerFactory.h" +#include "ROOTGSLNLSMinimizer.h" +#include "ROOTGSLSimAnMinimizer.h" +#include "MinimizerTest.h" +#include "MinimizerScan.h" +#include "ROOTMinimizer.h" + +#include <boost/assign/list_of.hpp> +#include <iomanip> + +MinimizerFactory::Catalogue MinimizerFactory::m_catalogue = MinimizerFactory::Catalogue(); + + +MinimizerFactory::Catalogue::Catalogue() +{ + // our minimizers + m_data["Test"] = boost::assign::list_of(""); + m_data["Scan"] = boost::assign::list_of(""); + // ROOT minimizers + m_data["Minuit"] = boost::assign::list_of("Migrad")("Simplex")("Combined")("Scan"); + m_data["Minuit2"] = boost::assign::list_of("Migrad")("Simplex")("Combined")("Scan")("Fumili"); + m_data["Fumili"] = boost::assign::list_of(""); + m_data["GSLMultiMin"] = boost::assign::list_of("ConjugateFR")("ConjugatePR")("BFGS")("BFGS2")("SteepestDescent"); + m_data["GSLMultiFit"] = boost::assign::list_of(""); + m_data["GSLSimAn"] = boost::assign::list_of(""); + m_data["Genetic"] = boost::assign::list_of(""); +} + + +void MinimizerFactory::Catalogue::print(std::ostream &ostr) const +{ + for(MinimizerFactory::Catalogue::const_iterator it=m_data.begin(); it!=m_data.end(); ++it) { + ostr << std::setw(20) << std::left<< it->first << " : "; + for(size_t i=0; i<it->second.size(); ++i ) { + ostr << it->second[i] << " "; + } + ostr << std::endl; + } + ostr << std::endl; +} + + +bool MinimizerFactory::Catalogue::isValid(const std::string &minimizer, const std::string &algorithm) const +{ + // check minimizers names + MinimizerFactory::Catalogue::const_iterator it = m_data.find(minimizer); + if(it != m_data.end() ) { + // check minimizer's algorithm type + for(size_t i=0; i<it->second.size(); ++i ) if(it->second[i] == algorithm ) return true; + } + return false; +} + +void MinimizerFactory::print_catalogue() +{ + std::cout << m_catalogue; +} + + +IMinimizer *MinimizerFactory::createMinimizer(const std::string &minimizer, const std::string &algorithm, const std::string &options) +{ + if( !m_catalogue.isValid(minimizer, algorithm) ) { + std::ostringstream ostr; + ostr << "MinimizerFactory::MinimizerFactory() -> Error! Wrong minimizer name. Possible names are:" << std::endl; + ostr << m_catalogue; + throw LogicErrorException(ostr.str()); + } + + IMinimizer *result(0); + if( minimizer == "Test" ) { + result = new MinimizerTest(); + } else if( minimizer == "Scan" ) { + result = new MinimizerScan(); + } else { + result = new ROOTMinimizer(minimizer, algorithm); + } + if( !options.empty() ) { + try { + result->setOptions(options); + } catch (NotImplementedException &e) { + std::cout << "MinimizerFactory::createMinimizer() -> Warning! Minimizer doesn't have method implemented" << e.what() << std::endl; + } + } + + return result; +} + + diff --git a/App/src/ROOTGSLSimAnMinimizer.cpp b/App/src/ROOTGSLSimAnMinimizer.cpp new file mode 100644 index 0000000000000000000000000000000000000000..2dc06326113e94afcb48754fa2a5daadff9bf695 --- /dev/null +++ b/App/src/ROOTGSLSimAnMinimizer.cpp @@ -0,0 +1,219 @@ +#include "ROOTGSLSimAnMinimizer.h" + + +#include "Math/GSLSimAnMinimizer.h" +#include "Math/WrappedParamFunction.h" +#include "Math/Error.h" + +#include "Math/MinimTransformFunction.h" +#include "Math/MultiNumGradFunction.h" // needed to use transformation function + +#include <iostream> +#include <cassert> + +namespace ROOT { + + namespace Patch { + + +// GSLSimAnMinimizer implementation + +GSLSimAnMinimizer::GSLSimAnMinimizer( int /* ROOT::Math::EGSLSimAnMinimizerType type */ ) : + fDim(0), + fOwnFunc(false), + fObjFunc(0), + fMinVal(0) +{ + // Constructor implementation : create GSLMultiFit wrapper object + + fValues.reserve(10); + fNames.reserve(10); + fSteps.reserve(10); + + SetMaxIterations(100); + SetPrintLevel(0); +} + +GSLSimAnMinimizer::~GSLSimAnMinimizer () { + if ( fOwnFunc && fObjFunc) delete fObjFunc; +} + +bool GSLSimAnMinimizer::SetVariable(unsigned int ivar, const std::string & name, double val, double step) { + // set variable in minimizer - support only free variables + // no transformation implemented - so far + if (ivar > fValues.size() ) return false; + if (ivar == fValues.size() ) { + fValues.push_back(val); + fNames.push_back(name); + // step is the simmulated annealing scale + fSteps.push_back( step ); + fVarTypes.push_back(kDefault); + } + else { + fValues[ivar] = val; + fNames[ivar] = name; + fSteps[ivar] = step; + fVarTypes[ivar] = kDefault; + + // remove bounds if needed + std::map<unsigned int, std::pair<double, double> >::iterator iter = fBounds.find(ivar); + if ( iter != fBounds.end() ) fBounds.erase (iter); + } + return true; + +} + +bool GSLSimAnMinimizer::SetLowerLimitedVariable(unsigned int ivar, const std::string & name, double val, double step, double lower ) { + bool ret = SetVariable(ivar, name, val, step); + if (!ret) return false; + fBounds[ivar] = std::make_pair( lower, lower); + fVarTypes[ivar] = kLowBound; + return true; +} +bool GSLSimAnMinimizer::SetUpperLimitedVariable(unsigned int ivar, const std::string & name, double val, double step, double upper) { + bool ret = SetVariable(ivar, name, val, step); + if (!ret) return false; + fBounds[ivar] = std::make_pair( upper, upper); + fVarTypes[ivar] = kUpBound; + return true; +} +bool GSLSimAnMinimizer::SetLimitedVariable(unsigned int ivar, const std::string & name, double val, double step, double lower, double upper ) { + bool ret = SetVariable(ivar, name, val, step); + if (!ret) return false; + fBounds[ivar] = std::make_pair( lower, upper); + fVarTypes[ivar] = kBounds; + return true; +} + + + +bool GSLSimAnMinimizer::SetFixedVariable(unsigned int ivar, const std::string & name, double val) { + /// set fixed variable (override if minimizer supports them ) + // use zero step size + bool ret = SetVariable(ivar, name, val, 0.); + if (!ret) return false; + fVarTypes[ivar] = kFix; + return true; +} + +bool GSLSimAnMinimizer::SetVariableValue(unsigned int ivar, double val) { + // set variable value in minimizer + // no transformation implemented - so far + if (ivar > fValues.size() ) return false; + fValues[ivar] = val; + return true; +} + +bool GSLSimAnMinimizer::SetVariableValues( const double * x) { + // set all variable values in minimizer + if (x == 0) return false; + std::copy(x,x+fValues.size(), fValues.begin() ); + return true; +} + + +void GSLSimAnMinimizer::SetFunction(const ROOT::Math::IMultiGenFunction & func) { + // set the function to minimize + + // keep pointers to the chi2 function + fObjFunc = &func; + fDim = func.NDim(); +} + +void GSLSimAnMinimizer::SetFunction(const ROOT::Math::IMultiGradFunction & func ) { + // set the function to minimize + // use the other methods + SetFunction( static_cast<const ROOT::Math::IMultiGenFunction &>(func) ); +} + + +bool GSLSimAnMinimizer::Minimize() { + // set initial parameters of the minimizer + int debugLevel = PrintLevel(); + + if (debugLevel >=1 ) std::cout <<"Minimize using GSLSimAnMinimizer " << std::endl; + + + // adapt the steps (use largers) + for (unsigned int i = 0; i < fSteps.size() ; ++i) + fSteps[i] *= 10; + + // vector of internal values (copied by default) + std::vector<double> xvar (fValues ); + std::vector<double> steps (fSteps); + + + // check if a transformation is needed + bool doTransform = (fBounds.size() > 0); + unsigned int ivar = 0; + while (!doTransform && ivar < fVarTypes.size() ) { + doTransform = (fVarTypes[ivar++] != kDefault ); + } + + // if needed do transformation and wrap objective function in a new transformation function + // and transform from external variables (and steps) to internals one + MinimTransformFunction * trFunc = 0; + if (doTransform) { + // since objective function is gradient build a gradient function for the transformation + // although gradient is not needed + + trFunc = new MinimTransformFunction ( new MultiNumGradFunction( *fObjFunc), fVarTypes, fValues, fBounds ); + + trFunc->InvTransformation(&fValues.front(), &xvar[0]); + + // need to transform also the steps + trFunc->InvStepTransformation(&fValues.front(), &fSteps.front(), &steps[0]); + + xvar.resize( trFunc->NDim() ); + steps.resize( trFunc->NDim() ); + + fObjFunc = trFunc; + fOwnFunc = true; // flag to indicate we need to delete the function + } + + assert (xvar.size() == steps.size() ); + + // output vector + std::vector<double> xmin(xvar.size() ); + + int iret = fSolver.Solve(*fObjFunc, &xvar.front(), &steps.front(), &xmin[0], (debugLevel > 1) ); + + fMinVal = (*fObjFunc)(&xmin.front() ); + + // get the result (transform if needed) + if (trFunc != 0) { + const double * xtrans = trFunc->Transformation(&xmin.front()); + assert(fValues.size() == trFunc->NTot() ); + assert( trFunc->NTot() == NDim() ); + std::copy(xtrans, xtrans + trFunc->NTot(), fValues.begin() ); + } + else { + // case of no transformation applied + assert( fValues.size() == xmin.size() ); + std::copy(xmin.begin(), xmin.end(), fValues.begin() ); + } + + + + if (debugLevel >=1 ) { + if (iret == 0) + std::cout << "GSLSimAnMinimizer: Minimum Found" << std::endl; + else + std::cout << "GSLSimAnMinimizer: Error in solving" << std::endl; + + int pr = std::cout.precision(18); + std::cout << "FVAL = " << fMinVal << std::endl; + std::cout.precision(pr); + for (unsigned int i = 0; i < fDim; ++i) + std::cout << fNames[i] << "\t = " << fValues[i] << std::endl; + } + + + return ( iret == 0) ? true : false; +} + + + + } // end namespace Patch + +} // end namespace ROOT diff --git a/App/src/ROOTMinimizer.cpp b/App/src/ROOTMinimizer.cpp index 45f7efc9bfe5a46a4390c80eaf5e00e2846c17b2..9d8ddd820064ace0ced970075506a671a3c9f6dc 100644 --- a/App/src/ROOTMinimizer.cpp +++ b/App/src/ROOTMinimizer.cpp @@ -7,6 +7,7 @@ #include <sstream> #include <boost/assign/list_of.hpp> #include "ROOTGSLNLSMinimizer.h" +#include "ROOTGSLSimAnMinimizer.h" /* ************************************************************************* */ @@ -22,15 +23,15 @@ ROOTMinimizer::ROOTMinimizer(const std::string &minimizer_name, const std::strin , m_chi2_func(0) , m_gradient_func(0) { - if( !isValidNames(m_minimizer_name, m_algo_type) ) throw LogicErrorException("ROOTMinimizer::ROOTMinimizer() -> Error! Wrong minimizer initialization parameters."); - if( m_minimizer_name == "GSLMultiFit") { // hacked version of ROOT's GSL Levenberg-Marquardt minimizer m_root_minimizer = new ROOT::Patch::GSLNLSMinimizer(2); + }else if( m_minimizer_name == "GSLSimAn") { + // hacked version of ROOT's GSL Simulated annealing minimizer + m_root_minimizer = new ROOT::Patch::GSLSimAnMinimizer(); } else { m_root_minimizer = ROOT::Math::Factory::CreateMinimizer(minimizer_name, algo_type ); } - if( !m_root_minimizer ) throw NullPointerException("ROOTMinimizer::ROOTMinimizer() -> Error! Can't create minimizer."); } @@ -42,45 +43,6 @@ ROOTMinimizer::~ROOTMinimizer() } -/* ************************************************************************* */ -// checking validity of the combination minimizer_name and algo_type -/* ************************************************************************* */ -bool ROOTMinimizer::isValidNames(const std::string &minimizer_name, const std::string &algo_type) -{ - // valid minimizer names and algo types - typedef std::map<std::string, std::vector<std::string > > algotypes_t; - algotypes_t algoTypes; - algoTypes["Minuit"] = boost::assign::list_of("Migrad")("Simplex")("Combined")("Scan"); - algoTypes["Minuit2"] = boost::assign::list_of("Migrad")("Simplex")("Combined")("Scan")("Fumili"); - algoTypes["Fumili"] = boost::assign::list_of(""); - algoTypes["GSLMultiMin"] = boost::assign::list_of("ConjugateFR")("ConjugatePR")("BFGS")("BFGS2")("SteepestDescent"); - algoTypes["GSLMultiFit"] = boost::assign::list_of(""); - algoTypes["GSLSimAn"] = boost::assign::list_of(""); - algoTypes["Genetic"] = boost::assign::list_of(""); - - // check minimizers names - algotypes_t::iterator it = algoTypes.find(minimizer_name); - if(it != algoTypes.end() ) { - // check minimizer's algorithm type - for(size_t i=0; i<it->second.size(); ++i ) if(it->second[i] == algo_type ) return true; - } - - // if not, print complaining and return false - std::cout << "ROOTMinimizer::isValidNames() -> Warning! Wrong minimizer name '" << minimizer_name << "' and/or algorithm type '" << algo_type << "'." << std::endl; - std::cout << "List of allowed minimizers (and algos)" << std::endl; - for(it = algoTypes.begin(); it!= algoTypes.end(); ++it) { - std::cout << std::setw(20) << std::left<< it->first << " : "; - for(size_t i=0; i<it->second.size(); ++i ) { - std::cout << it->second[i] << " "; - } - std::cout << std::endl; - } - std::cout << std::endl; - - return false; -} - - /* ************************************************************************* */ // check if type of algorithm is Levenberg-Marquardt or similar // (that means that it requires manual gradient calculations) diff --git a/App/src/SampleFactory.cpp b/App/src/SampleFactory.cpp index 6134abfcfaa3dc119d013555461468160f81727f..7b2f98c26d1f34f558486741b8900b8296b383df 100644 --- a/App/src/SampleFactory.cpp +++ b/App/src/SampleFactory.cpp @@ -8,8 +8,7 @@ SampleFactory::SampleFactory() { - // Experiment will take care about samples - setOwnObjects(true); + setOwnObjects(true); //factory will take care about samples // samples used for fresnel coefficients validation registerItem("AirOnSubstrate", StandardSamples::AirOnSubstrate); diff --git a/App/src/TestFittingModule1.cpp b/App/src/TestFittingModule1.cpp index 4c7e799ac8f1bf2c45adc1cf7c6d35af03e3d395..b193f95668aacb5197efa3cd417b8393df015ff5 100644 --- a/App/src/TestFittingModule1.cpp +++ b/App/src/TestFittingModule1.cpp @@ -13,7 +13,7 @@ #include "FormFactors.h" #include "Exceptions.h" #include "DrawHelper.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "ResolutionFunction2DSimple.h" #include "AttLimits.h" #include "ISquaredFunction.h" @@ -21,7 +21,7 @@ #include "IObserver.h" #include "FitSuite.h" -#include "ROOTMinimizer.h" +#include "MinimizerFactory.h" #include "TROOT.h" #include "TCanvas.h" @@ -62,7 +62,7 @@ void TestFittingModule1::execute() //m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Migrad") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("GSLMultiFit", "") ); - m_fitSuite->setMinimizer( new ROOTMinimizer("Fumili", "") ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("Minuit2", "Migrad") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Fumili") ); @@ -84,9 +84,9 @@ void TestFittingModule1::execute() //minim->SetMaxIterations(50); // for GSL //minim->SetPrintLevel(4); - m_fitSuite->attachObserver( new FitSuiteObserverPrint() ); - m_fitSuite->attachObserver( new FitSuiteObserverDraw() ); - //fitSuite->attachObserver( new FitSuiteObserverWriteTree() ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver() ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createDrawObserver() ); + //fitSuite->attachObserver( ObserverFactory::createTreeObserver() ); m_fitSuite->runFit(); } diff --git a/App/src/TestFittingModule2.cpp b/App/src/TestFittingModule2.cpp index 627a9d88578c7d7b1e60176865b9780ba211dc81..42825482bd2665fd1c4acf72aa252242a5a8de58 100644 --- a/App/src/TestFittingModule2.cpp +++ b/App/src/TestFittingModule2.cpp @@ -13,15 +13,15 @@ #include "FormFactors.h" #include "Exceptions.h" #include "DrawHelper.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "ResolutionFunction2DSimple.h" #include "AttLimits.h" #include "IIntensityFunction.h" #include "IObserver.h" #include "FitSuite.h" -#include "ROOTMinimizer.h" #include "FitSuiteStrategy.h" +#include "MinimizerFactory.h" #include "TROOT.h" #include "TCanvas.h" @@ -41,9 +41,9 @@ TestFittingModule2::TestFittingModule2() // setting up fitSuite m_fitSuite = new FitSuite(); - m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Migrad") ); - m_fitSuite->attachObserver( new FitSuiteObserverPrint() ); - m_fitSuite->attachObserver( new FitSuiteObserverDraw() ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("Minuit2", "Migrad") ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver() ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createDrawObserver() ); } @@ -163,7 +163,7 @@ void TestFittingModule2::fit_example_strategies() m_fitSuite->addExperimentAndRealData(*mp_experiment, *mp_real_data); - m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Migrad") ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("Minuit2", "Migrad") ); m_fitSuite->runFit(); diff --git a/App/src/TestFittingModule3.cpp b/App/src/TestFittingModule3.cpp index bf89f2182f950d80743d474133faacc17847d5b0..9b506b371752eb0757b3f2c2cd2248ae0b8bdc6f 100644 --- a/App/src/TestFittingModule3.cpp +++ b/App/src/TestFittingModule3.cpp @@ -13,7 +13,7 @@ #include "FormFactors.h" #include "Exceptions.h" #include "DrawHelper.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "ResolutionFunction2DSimple.h" #include "AttLimits.h" #include "OutputDataFunctions.h" @@ -22,7 +22,7 @@ #include "IObserver.h" #include "FitSuite.h" -#include "ROOTMinimizer.h" +#include "MinimizerFactory.h" #include "TROOT.h" #include "TCanvas.h" @@ -79,9 +79,9 @@ void TestFittingModule3::execute() m_fitSuite->addExperimentAndRealData(*m_experiment, *(*it)); } - m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Migrad") ); - m_fitSuite->attachObserver( new FitSuiteObserverPrint() ); - m_fitSuite->attachObserver( new FitSuiteObserverDraw(1) ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("Minuit2", "Migrad") ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver() ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createDrawObserver(1) ); m_fitSuite->runFit(); diff --git a/App/src/TestIsGISAXS12.cpp b/App/src/TestIsGISAXS12.cpp index 30c6edbfa931f4bf0ea75915ceaa9453ff9c6800..bf38847d3ad66c3d95618b746e99dcb3282bb791 100644 --- a/App/src/TestIsGISAXS12.cpp +++ b/App/src/TestIsGISAXS12.cpp @@ -19,14 +19,14 @@ #include "GISASExperiment.h" #include "DrawHelper.h" #include "FitSuite.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "ResolutionFunction2DSimple.h" #include "MathFunctions.h" -#include "ROOTMinimizer.h" #include "MinimizerTest.h" #include "OutputDataFunctions.h" #include "ExperimentConstants.h" #include "OutputDataIOFactory.h" +#include "MinimizerFactory.h" #include <iostream> #include <fstream> @@ -203,12 +203,14 @@ void TestIsGISAXS12::run_isgisaxs_fit() // creating fit suite m_fitSuite = new FitSuite(); - m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Migrad") ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("Minuit2", "Migrad") ); + //m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Fumili") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("GSLMultiFit", "") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("Fumili", "") ); - m_fitSuite->attachObserver( new FitSuiteObserverPrint(10) ); - m_fitSuite->attachObserver( new FitSuiteObserverDraw(10) ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver(10) ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createDrawObserver(10) ); + // ROOT::Math::Minimizer *minim = (dynamic_cast<ROOTMinimizer *>(m_fitSuite->getMinimizer()))->getROOTMinimizer(); // minim->SetPrintLevel(4); diff --git a/App/src/TestIsGISAXS13.cpp b/App/src/TestIsGISAXS13.cpp index 7919f3c86bb81814559517ba5e76e99b83b99d39..a468b9a997d9c88a5927a46a5221c37f55223b6b 100644 --- a/App/src/TestIsGISAXS13.cpp +++ b/App/src/TestIsGISAXS13.cpp @@ -20,15 +20,16 @@ #include "GISASExperiment.h" #include "DrawHelper.h" #include "FitSuite.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "ResolutionFunction2DSimple.h" #include "MathFunctions.h" -#include "ROOTMinimizer.h" #include "OutputDataFunctions.h" #include "ExperimentConstants.h" #include "OutputDataIOFactory.h" #include "IsGISAXSData.h" #include "MinimizerScan.h" +#include "ROOTGSLSimAnMinimizer.h" +#include "MinimizerFactory.h" #include <iostream> #include <fstream> @@ -81,7 +82,9 @@ void TestIsGISAXS13::run_isgisaxs_fit() //m_fitSuite->setMinimizer( new ROOTMinimizer("GSLMultiFit", "") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("Fumili", "") ); -// m_fitSuite->setMinimizer( new ROOTMinimizer("GSLSimAn") ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("GSLSimAn") ); + //m_fitSuite->setMinimizer( new ROOTMinimizer("Genetic") ); + // ROOT::Math::Minimizer *minim = (dynamic_cast<ROOTMinimizer *>(m_fitSuite->getMinimizer()))->getROOTMinimizer(); // minim->SetPrintLevel(4); // minim->SetMaxIterations(400); @@ -92,11 +95,25 @@ void TestIsGISAXS13::run_isgisaxs_fit() // minim->SetPrecision(1.); - m_fitSuite->setMinimizer( new MinimizerScan(20) ); +// ROOT::Math::Minimizer *minim = (dynamic_cast<ROOTMinimizer *>(m_fitSuite->getMinimizer()))->getROOTMinimizer(); +// if( !minim ) throw NullPointerException("TestToyExperiment::execute() -> Can't cast minimizer #1"); +// ROOT::Patch::GSLSimAnMinimizer *siman = dynamic_cast<ROOT::Patch::GSLSimAnMinimizer *>(minim); +// if( !siman ) throw NullPointerException("TestToyExperiment::execute() -> Can't cast minimizer #2"); +// ROOT::Math::GSLSimAnParams &pars = siman->getSolver().Params(); +// siman->SetPrintLevel(2); +// pars.n_tries = 100; +// pars.iters_fixed_T = 10; +// pars.step_size = 0.5; +// pars.k = 1; +// pars.t_initial = 50; +// pars.mu = 1.05; +// pars.t_min = 0.5; + +// m_fitSuite->setMinimizer( new MinimizerScan(20) ); -// m_fitSuite->attachObserver( new FitSuiteObserverPrint(10) ); -// m_fitSuite->attachObserver( new FitSuiteObserverDraw(10) ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver(10) ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createDrawObserver(50) ); // Migrad @@ -128,14 +145,22 @@ void TestIsGISAXS13::run_isgisaxs_fit() // m_fitSuite->addFitParameter("*Normalizer/scale", 1e5, 1e5, AttLimits::limited(1e4, 2e5)); // m_fitSuite->addFitParameter("*Normalizer/shift", 10, 10, AttLimits::limited(0., 20.)); - m_fitSuite->addFitParameter("*Normalizer/scale", 1.333150e+05, 1e5, AttLimits::fixed()); - m_fitSuite->addFitParameter("*Normalizer/shift", 1.212894e+00, 10, AttLimits::fixed()); - m_fitSuite->addFitParameter("*SampleBuilder/particle_radius", 4.0e+00*Units::nanometer, 0.01*Units::nanometer, AttLimits::limited(2.0, 8.0) ); - m_fitSuite->addFitParameter("*SampleBuilder/dispersion_radius", 1.991729e-01, 0.2, AttLimits::fixed() ); + + m_fitSuite->addFitParameter("*Normalizer/scale", 1e5, 1e3, AttLimits::limited(1e4, 2e5)); + m_fitSuite->addFitParameter("*Normalizer/shift", 10, 0.1, AttLimits::limited(0., 20.)); + +// m_fitSuite->addFitParameter("*Normalizer/scale", 1.333150e+05, 1e5, AttLimits::fixed()); +// m_fitSuite->addFitParameter("*Normalizer/shift", 1.212894e+00, 10, AttLimits::fixed()); + m_fitSuite->addFitParameter("*SampleBuilder/particle_radius", 4.0e+00*Units::nanometer, 0.04*Units::nanometer, AttLimits::limited(2.0, 8.0) ); +// m_fitSuite->addFitParameter("*SampleBuilder/dispersion_radius", 1.991729e-01, 0.2, AttLimits::fixed() ); + m_fitSuite->addFitParameter("*SampleBuilder/dispersion_radius", 0.2, 0.002, AttLimits::limited(0.01, 1.0) ); //m_fitSuite->addFitParameter("*SampleBuilder/height_aspect_ratio", 9.937262e-01, 1.0, AttLimits::fixed() ); - m_fitSuite->addFitParameter("*SampleBuilder/height_aspect_ratio", 9.937262e-01, 1.0, AttLimits::limited(0.5, 1.5) ); - m_fitSuite->addFitParameter("*SampleBuilder/interf_distance", 1.497518e+01*Units::nanometer, 15.*Units::nanometer, AttLimits::fixed() ); - m_fitSuite->addFitParameter("*SampleBuilder/interf_width", 3.023607e+00*Units::nanometer, 3.0*Units::nanometer, AttLimits::fixed() ); + m_fitSuite->addFitParameter("*SampleBuilder/height_aspect_ratio", 0.8, 0.08, AttLimits::limited(0.5, 1.5) ); + +// m_fitSuite->addFitParameter("*SampleBuilder/interf_distance", 1.497518e+01*Units::nanometer, 15.*Units::nanometer, AttLimits::fixed() ); +// m_fitSuite->addFitParameter("*SampleBuilder/interf_width", 3.023607e+00*Units::nanometer, 3.0*Units::nanometer, AttLimits::fixed() ); + m_fitSuite->addFitParameter("*SampleBuilder/interf_distance", 15*Units::nanometer, 0.015*Units::nanometer, AttLimits::limited(0.01, 50.0) ); + m_fitSuite->addFitParameter("*SampleBuilder/interf_width", 3*Units::nanometer, 0.03*Units::nanometer, AttLimits::limited(0.01, 10.) ); // setting up fitSuite @@ -143,7 +168,7 @@ void TestIsGISAXS13::run_isgisaxs_fit() chiModule.setChiSquaredFunction( SquaredFunctionWithSystematicError(0.08) ); //chiModule.setChiSquaredFunction(SquaredFunctionDefault());// isgisaxs uses epsion=0, which correspond to our SquaredFunctionDefault chiModule.setOutputDataNormalizer( OutputDataNormalizerScaleAndShift() ); - //chiModule.setIntensityFunction( IntensityFunctionSqrt() ); + //chiModule.setIntensityFunction( IntensityFunctionLog() ); for(IsGISAXSData::DataSet_t::iterator it=isgi_scans.begin(); it!= isgi_scans.end(); ++it) { m_fitSuite->addExperimentAndRealData(*m_experiment, *(*it), chiModule); @@ -151,73 +176,73 @@ void TestIsGISAXS13::run_isgisaxs_fit() m_fitSuite->runFit(); - TCanvas *c1 = new TCanvas("c1","c1",1024,768); - c1->cd(); +// TCanvas *c1 = new TCanvas("c1","c1",1024,768); +// c1->cd(); - MinimizerScan *minim = dynamic_cast<MinimizerScan *>(m_fitSuite->getMinimizer()); - std::cout << "Results " << minim->getMinValue() << std::endl; - for(size_t i=0; i<minim->getNumberOfVariables(); ++i) { - std::cout << i << " " << minim->getValueOfVariableAtMinimum(i) << std::endl; - } +// MinimizerScan *minim = dynamic_cast<MinimizerScan *>(m_fitSuite->getMinimizer()); +// std::cout << "Results " << minim->getMinValue() << std::endl; +// for(size_t i=0; i<minim->getNumberOfVariables(); ++i) { +// std::cout << i << " " << minim->getValueOfVariableAtMinimum(i) << std::endl; +// } - const OutputData<double> *data = minim->getOutputData(); - if(data->getNdimensions()==1) { - TH1 * h1 = IsGISAXSTools::getOutputDataTH123D(*data,"hist"); - h1->Draw(); - gPad->SetLogy(); - } else { - TH2D *hist = IsGISAXSTools::getOutputDataTH2D( *data, "hist"); - hist->Draw("colz"); - gPad->SetLogz(); - } +// const OutputData<double> *data = minim->getOutputData(); +// if(data->getNdimensions()==1) { +// TH1 * h1 = IsGISAXSTools::getOutputDataTH123D(*data,"hist"); +// h1->Draw(); +// gPad->SetLogy(); +// } else { +// TH2D *hist = IsGISAXSTools::getOutputDataTH2D( *data, "hist"); +// hist->Draw("colz"); +// gPad->SetLogz(); +// } -// // drawing results -// TCanvas *c2 = new TCanvas("c2","GISASFW fit results",800,600); -// c2->Divide(2,2); -// 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<m_fitSuite->getFitObjects()->size(); ++i_set) { -// 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"); -// hreal->SetLineColor(kBlue); -// gPad->SetLogy(); -// hreal->DrawCopy(); -// hsimul->DrawCopy("same"); -// if(i_set==0) leg1->AddEntry(hreal,"GISASFW data","lp"); -// if(i_set==0) leg1->AddEntry(hsimul,"GISASFW simul","lp"); -// } -// c2->cd(1); leg1->Draw(); -// c2->cd(2); leg1->Draw(); - -// // drawing ratio -// 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<m_fitSuite->getFitObjects()->size(); ++i_set) { -// c2->cd(3+1); -// const FitObject *obj = m_fitSuite->getFitObjects()->getObject(i_set); -// OutputData<double > *real = obj->getChiSquaredModule()->getRealData()->clone(); -// OutputData<double > *simul = obj->getChiSquaredModule()->getSimulationData()->clone(); - -// c2->cd((int)(i_set+3)); -// *simul /= *real; -// TH1D *hratio = IsGISAXSTools::getOutputDataScanHist(*simul,"gisasfw_real_simul_ratio"); -// hratio->DrawCopy(); -// if(i_set==0) { -// leg2->AddEntry(hratio,"GISASFW simul/real","lp"); -// } -// delete real; -// delete simul; -// } -// c2->cd(3); leg2->Draw(); -// c2->cd(4); leg2->Draw(); + // drawing results + TCanvas *c2 = new TCanvas("c2","GISASFW fit results",800,600); + c2->Divide(2,2); + 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<m_fitSuite->getFitObjects()->size(); ++i_set) { + 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"); + hreal->SetLineColor(kBlue); + gPad->SetLogy(); + hreal->DrawCopy(); + hsimul->DrawCopy("same"); + if(i_set==0) leg1->AddEntry(hreal,"GISASFW data","lp"); + if(i_set==0) leg1->AddEntry(hsimul,"GISASFW simul","lp"); + } + c2->cd(1); leg1->Draw(); + c2->cd(2); leg1->Draw(); + + // drawing ratio + 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<m_fitSuite->getFitObjects()->size(); ++i_set) { + c2->cd(3+1); + const FitObject *obj = m_fitSuite->getFitObjects()->getObject(i_set); + OutputData<double > *real = obj->getChiSquaredModule()->getRealData()->clone(); + OutputData<double > *simul = obj->getChiSquaredModule()->getSimulationData()->clone(); + + c2->cd((int)(i_set+3)); + *simul /= *real; + TH1D *hratio = IsGISAXSTools::getOutputDataScanHist(*simul,"gisasfw_real_simul_ratio"); + hratio->DrawCopy(); + if(i_set==0) { + leg2->AddEntry(hratio,"GISASFW simul/real","lp"); + } + delete real; + delete simul; + } + c2->cd(3); leg2->Draw(); + c2->cd(4); leg2->Draw(); -// c2->Update(); + c2->Update(); } diff --git a/App/src/TestIsGISAXS5.cpp b/App/src/TestIsGISAXS5.cpp index 8bcdeafb1ffb547665adf4f658992f6805d55bc1..34402849d89306c30a6219b03cfd25df9068fce6 100644 --- a/App/src/TestIsGISAXS5.cpp +++ b/App/src/TestIsGISAXS5.cpp @@ -19,10 +19,10 @@ #include "GISASExperiment.h" #include "DrawHelper.h" #include "FitSuite.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "ResolutionFunction2DSimple.h" #include "MathFunctions.h" -#include "ROOTMinimizer.h" +#include "MinimizerFactory.h" #include "OutputDataFunctions.h" #include "ExperimentConstants.h" #include "OutputDataIOFactory.h" @@ -156,12 +156,14 @@ void TestIsGISAXS5::run_isgisaxs_fit() // creating fit suite m_fitSuite = new FitSuite(); - m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Migrad") ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("Minuit2", "Migrad") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Fumili") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("GSLMultiFit", "") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("Fumili", "") ); - m_fitSuite->attachObserver( new FitSuiteObserverPrint(10) ); - m_fitSuite->attachObserver( new FitSuiteObserverDraw(10) ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver(10) ); + m_fitSuite->attachObserver( FitSuiteObserverFactory::createDrawObserver(10) ); + + // ROOT::Math::Minimizer *minim = (dynamic_cast<ROOTMinimizer *>(m_fitSuite->getMinimizer()))->getROOTMinimizer(); // minim->SetPrintLevel(4); // minim->SetMaxIterations(400); diff --git a/App/src/TestMesoCrystal2.cpp b/App/src/TestMesoCrystal2.cpp index db8207d23db65709e8589acb32f9f6112d70a6cd..2f3f912047161bccfda0ee30aa3fb90368595758 100644 --- a/App/src/TestMesoCrystal2.cpp +++ b/App/src/TestMesoCrystal2.cpp @@ -21,12 +21,12 @@ #include "OutputDataReader.h" #include "OutputDataIOFactory.h" #include "FitSuite.h" -#include "ROOTMinimizer.h" #include "SampleFactory.h" #include "TRange.h" -#include "FitSuiteHelper.h" +#include "FitSuiteObserverFactory.h" #include "AttLimits.h" #include "ProgramOptions.h" +#include "MinimizerFactory.h" #include "TCanvas.h" #include "TH2D.h" @@ -90,9 +90,11 @@ void TestMesoCrystal2::execute() fitSuite->addExperimentAndRealData(*mp_experiment, *real_data); // fitSuite->addExperimentAndRealData(mp_experiment, real_data_quarter); - fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Combined") ); - ROOT::Math::Minimizer *minim = (dynamic_cast<ROOTMinimizer *>(fitSuite->getMinimizer()))->getROOTMinimizer(); - minim->SetStrategy(2); // 0- not accurate, 1 - normal, 2 - acurate (maximum FCN calls) + fitSuite->setMinimizer( MinimizerFactory::createMinimizer("Minuit2", "Combined") ); +// ROOT::Math::Minimizer *minim = (dynamic_cast<ROOTMinimizer *>(fitSuite->getMinimizer()))->getROOTMinimizer(); +// minim->SetStrategy(2); // 0- not accurate, 1 - normal, 2 - acurate (maximum FCN calls) + + fitSuite->addFitParameter("*/lattice_length_a", 6.2*Units::nanometer, 1.0*Units::nanometer, AttLimits::limited(4.0*Units::nanometer, 8.0*Units::nanometer) ); @@ -197,9 +199,9 @@ void TestMesoCrystal2::execute() // FitSuiteObserverDraw *drawObserver = new FitSuiteObserverDraw(canvas_name); // fitSuite->attachObserver(drawObserver); - fitSuite->attachObserver( new FitSuiteObserverPrint() ); - fitSuite->attachObserver( new FitSuiteObserverDraw() ); - fitSuite->attachObserver( new FitSuiteObserverWriteTree() ); + fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver() ); + fitSuite->attachObserver( FitSuiteObserverFactory::createDrawObserver() ); + fitSuite->attachObserver( FitSuiteObserverFactory::createTreeObserver() ); fitSuite->runFit(); diff --git a/App/src/TestToyExperiment.cpp b/App/src/TestToyExperiment.cpp index a70c19e148e2c4a89b44c18aff211e0ef6c528c9..dccf77e6d049425433a200b8a42f14a7144d7770 100644 --- a/App/src/TestToyExperiment.cpp +++ b/App/src/TestToyExperiment.cpp @@ -2,12 +2,12 @@ #include "Exceptions.h" #include "IsGISAXSTools.h" #include "FitSuite.h" -#include "FitSuiteHelper.h" -#include "ROOTMinimizer.h" +#include "FitSuiteObserverFactory.h" #include "ExperimentConstants.h" +#include "ROOTGSLSimAnMinimizer.h" #include <iostream> - +#include "MinimizerFactory.h" /* ************************************************************************* */ // @@ -30,22 +30,6 @@ void ToyExperiment::runSimulation() } } -void ToyExperiment::runSimulationElement(size_t index) -{ - (void)index; -// if( !m_func ) throw NullPointerException("ToyExperiment::runSimulation() -> Error! No function is defined."); - -// m_func->SetParameters(&pars[0]); -// const std::string s_phi_f("phi_f"); -// const std::string s_alpha_f("alpha_f"); -// double phi_f = m_intensity_map.getValueOfAxis(s_phi_f, index); -// double alpha_f = m_intensity_map.getValueOfAxis(s_alpha_f, index); -// double value = m_func->Eval(phi_f, alpha_f); -// m_intensity_map[index] = value; - - throw NotImplementedException("ToyExperiment::runSimulationElement"); -} - void ToyExperiment::init_parameters() { @@ -71,13 +55,13 @@ TestToyExperiment::TestToyExperiment() , m_fitSuite(0) { - -// return; m_sigma_noise = 0.01; m_func_object = new SincXSincYFunctionObject(); // m_func = new TF2("sincxy", m_func_object, -5.,5., -5.,5., 3, "SincXSincYFunctionObject"); m_func = new TF2("sincxy", m_func_object, -10.,10., -10.,10., 3, "SincXSincYFunctionObject"); //m_func->SetParameters(1.0, 1.0, 0.5); // parameters we have to find + + } @@ -95,21 +79,49 @@ void TestToyExperiment::execute() { std::cout << "TestToyExperiment()::execute() -> Hello World!" << std::endl; + MinimizerFactory::print_catalogue(); + return; + + + initializeExperimentAndRealData(); // setting up fitSuite FitSuite *m_fitSuite = new FitSuite(); //m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Migrad") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("Minuit2", "Fumili") ); - m_fitSuite->setMinimizer( new ROOTMinimizer("Fumili") ); + //m_fitSuite->setMinimizer( new ROOTMinimizer("Fumili") ); //m_fitSuite->setMinimizer( new ROOTMinimizer("GSLMultiFit") ); + m_fitSuite->setMinimizer( MinimizerFactory::createMinimizer("GSLSimAn") ); + +// ROOT::Math::Minimizer *minim = (dynamic_cast<ROOTMinimizer *>(m_fitSuite->getMinimizer()))->getROOTMinimizer(); +// if( !minim ) throw NullPointerException("TestToyExperiment::execute() -> Can't cast minimizer #1"); +// ROOT::Patch::GSLSimAnMinimizer *siman = dynamic_cast<ROOT::Patch::GSLSimAnMinimizer *>(minim); +// if( !siman ) throw NullPointerException("TestToyExperiment::execute() -> Can't cast minimizer #2"); +// ROOT::Math::GSLSimAnParams &pars = siman->getSolver().Params(); +// siman->SetPrintLevel(2); +// pars.n_tries = 100; +// pars.iters_fixed_T = 10; +// pars.step_size = 1.0; +// pars.k = 1; +// pars.t_initial = 50; +// pars.mu = 1.05; +// pars.t_min = 0.1; + + + + //m_fitSuite->setMinimizer( new ROOTMinimizer("Genetic") ); + + + + + m_fitSuite->attachObserver( FitSuiteObserverFactory::createPrintObserver() ); +// m_fitSuite->attachObserver( ObserverFactory::createDrawObserver() ); -// m_fitSuite->attachObserver( new FitSuiteObserverPrint() ); -// m_fitSuite->attachObserver( new FitSuiteObserverDraw() ); - m_fitSuite->addFitParameter("*/par0", 1.0, 0.01); - m_fitSuite->addFitParameter("*/par1", 0.0, 0.01); - m_fitSuite->addFitParameter("*/par2", 0.0, 0.01); + m_fitSuite->addFitParameter("*/par0", 1.0, 0.01, AttLimits::limited(0.5, 1.5)); + m_fitSuite->addFitParameter("*/par1", 0.0, 0.01, AttLimits::limited(0.0, 3.0)); + m_fitSuite->addFitParameter("*/par2", 0.0, 0.01, AttLimits::limited(0.0, 3.0)); //m_fitSuite->addFitParameter("*/par2", -2.5, 0.01, AttLimits::fixed()); ChiSquaredModule chi_module; @@ -138,8 +150,8 @@ void TestToyExperiment::initializeExperimentAndRealData() // generating real data delete m_real_data; m_experiment->setParameter(0, 1.0); - m_experiment->setParameter(1, -2.0); - m_experiment->setParameter(2, -2.5); + m_experiment->setParameter(1, 2.0); + m_experiment->setParameter(2, 2.5); // m_experiment->setParameter(0, 2.0); // m_experiment->setParameter(1, 1.0); // m_experiment->setParameter(2, 0.5); diff --git a/Core/Algorithms/inc/Experiment.h b/Core/Algorithms/inc/Experiment.h index b67d9ea11f7a9d9119111541fabb6bbc78da05b6..43067de10c8dd91082682b40f0dd60c8ab4c2400 100644 --- a/Core/Algorithms/inc/Experiment.h +++ b/Core/Algorithms/inc/Experiment.h @@ -37,9 +37,6 @@ public: //! run a simulation with the current parameter settings virtual void runSimulation(); - //! run a simulation with the current parameter settings for single OutputData element (temp) - virtual void runSimulationElement(size_t index); - //! normalize the detector counts virtual void normalize(); diff --git a/Core/Algorithms/src/Experiment.cpp b/Core/Algorithms/src/Experiment.cpp index cb8015f0e22b109e12a945601bb37401c3ef59f4..be7c586a8143bdb1080217ca348cfdb42d605fd3 100644 --- a/Core/Algorithms/src/Experiment.cpp +++ b/Core/Algorithms/src/Experiment.cpp @@ -100,12 +100,6 @@ void Experiment::runSimulation() } -void Experiment::runSimulationElement(size_t /* index */) -{ - throw NotImplementedException("Experiment::runSimulationElement() -> Error! Not implemented."); -} - - void Experiment::normalize() { double incident_intensity = m_beam.getIntensity(); diff --git a/Core/Fitting/inc/IMinimizer.h b/Core/Fitting/inc/IMinimizer.h index 4b151e72ff312cc938131bdde3f526e87d042b8a..5f15887f7d6db3c654bc4f8d3f2b2dd4d0b0a53b 100644 --- a/Core/Fitting/inc/IMinimizer.h +++ b/Core/Fitting/inc/IMinimizer.h @@ -69,6 +69,9 @@ public: //! print fit results virtual void printResults() const; + + //! set minimizer command + virtual void setOptions(const std::string &option); }; @@ -132,4 +135,10 @@ inline void IMinimizer::printResults() const throw NotImplementedException("IMinimizer::printResults() -> Not implemented."); } +inline void IMinimizer::setOptions(const std::string &options) +{ + (void)options; + throw NotImplementedException("IMinimizer::setOptions() -> Not implemented."); +} + #endif // IMINIMIZER_H diff --git a/Core/Tools/inc/IObserver.h b/Core/Tools/inc/IObserver.h index 5d4abff69a503a7015fb3bd8c58ae263bfa90ebf..a5efb3b84124ba97e94e23635f2c1e6faf82d316 100644 --- a/Core/Tools/inc/IObserver.h +++ b/Core/Tools/inc/IObserver.h @@ -16,6 +16,8 @@ #include "Exceptions.h" #include <list> +#include <boost/shared_ptr.hpp> + class IObservable; //- ------------------------------------------------------------------- @@ -25,16 +27,16 @@ class IObservable; //- ------------------------------------------------------------------- class IObserver { public: - //! destructor detach observer from observed subject - virtual ~IObserver(); +// IObserver() : m_observed_subject(0) {} + virtual ~IObserver() {} + //! method which is used by observable subject to notify change in status virtual void update (IObservable *subject) = 0; - //! set pointer to observed subject - virtual void setObservedSubject(IObservable *subject); -protected: - IObserver() : m_observed_subject(0) {} -private: - IObservable *m_observed_subject; + +// //! set pointer to observed subject +// virtual void setObservedSubject(IObservable *subject); +//private: +// IObservable *m_observed_subject; }; @@ -45,22 +47,19 @@ private: //- ------------------------------------------------------------------- class IObservable { public: - typedef std::list<IObserver *> observers_t; + typedef boost::shared_ptr<IObserver > observer_t; + typedef std::list<observer_t > observerlist_t; + virtual ~IObservable(){} //! attach observer to the list of observers - virtual void attachObserver(IObserver *obj); - - //! detach observer from observers list - virtual void detachObserver(IObserver *obj); + virtual void attachObserver(observer_t obj); //! notify observers about change in status virtual void notifyObservers(); -protected: - IObservable(){} private: - observers_t m_observers; + observerlist_t m_observers; }; diff --git a/Core/Tools/src/IObserver.cpp b/Core/Tools/src/IObserver.cpp index 42c874f96cda01f81728499c9c775d96b48cc544..baf35dcc3329dcd31b301f1aa7a82e6ccb8b6e4b 100644 --- a/Core/Tools/src/IObserver.cpp +++ b/Core/Tools/src/IObserver.cpp @@ -1,43 +1,29 @@ #include "IObserver.h" - -/* ************************************************************************* */ -// Observer -/* ************************************************************************* */ -IObserver::~IObserver() -{ - if(m_observed_subject) m_observed_subject->detachObserver(this); -} - - -void IObserver::setObservedSubject(IObservable *subject) -{ - m_observed_subject = subject; -} +//void IObserver::setObservedSubject(IObservable *subject) +//{ +// m_observed_subject = subject; +//} - -/* ************************************************************************* */ -// Observable -/* ************************************************************************* */ -void IObservable::attachObserver(IObserver *obj) +void IObservable::attachObserver(observer_t obj) { - obj->setObservedSubject(this); +// obj->setObservedSubject(this); m_observers.push_back(obj); } -void IObservable::detachObserver(IObserver *obj) -{ - m_observers.remove(obj); -} - - void IObservable::notifyObservers() { - for(observers_t::iterator it = m_observers.begin(); it!=m_observers.end(); ++it) { + for(observerlist_t::iterator it = m_observers.begin(); it!=m_observers.end(); ++it) { (*it)->update(this); } } + +//void IObservable::detachObserver(IObserver *obj) +//{ +// m_observers.remove(obj); +//} + diff --git a/Tests/UnitTests/TestCore/TestCore b/Tests/UnitTests/TestCore/TestCore index 24bf8fa1f1d558f78c018fd07a67d30b5a67689b..ad8396d70203ec984177a19e12131242b3a6cbc3 100755 Binary files a/Tests/UnitTests/TestCore/TestCore and b/Tests/UnitTests/TestCore/TestCore differ diff --git a/shared.pri b/shared.pri index 2ae4058fad70efbe8a00db615ede67b750eb8f30..be0875e4c9583b1fca6768aac9f633afb55da259 100644 --- a/shared.pri +++ b/shared.pri @@ -98,8 +98,8 @@ DEPENDPATH += $${LOCATIONS} # optimization flag used in release builds (the -O2 is the default used by qmake) 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 # -ffast-math removed because of problems with NaNs +#QMAKE_CXXFLAGS_RELEASE -= -O2 +#QMAKE_CXXFLAGS_RELEASE += -g -O3 # -ffast-math removed because of problems with NaNs #QMAKE_STRIP=: # produces non-stripped (very large) libraries # to compile with GPERFTOOLS support for code profiling