diff --git a/Fit/Minimizer/BasicMinimizer.cpp b/Fit/Minimizer/BasicMinimizer.cpp index 62e06ecd28c645488e76295d15791df5d7750f0b..4582a8a6faf383fd7b2aa5cd2184226ef909bb27 100644 --- a/Fit/Minimizer/BasicMinimizer.cpp +++ b/Fit/Minimizer/BasicMinimizer.cpp @@ -14,6 +14,11 @@ // ************************************************************************** // #include "BasicMinimizer.h" +#include "ROOTMinimizerFunction.h" +#include "Math/Minimizer.h" +#include "FitParameter.h" +#include "FitSuiteParameters.h" + BasicMinimizer::BasicMinimizer(const std::string &minimizerName, const std::string &algorithmName) : m_minimizerName(minimizerName) @@ -22,6 +27,17 @@ BasicMinimizer::BasicMinimizer(const std::string &minimizerName, const std::stri } +BasicMinimizer::~BasicMinimizer() +{ + +} + +void BasicMinimizer::minimize() +{ + propagateOptions(); + rootMinimizer()->Minimize(); +} + std::string BasicMinimizer::minimizerName() const { return m_minimizerName; @@ -37,8 +53,116 @@ void BasicMinimizer::setAlgorithmName(const std::string &algorithmName) m_algorithmName = algorithmName; } -void BasicMinimizer::propagateOptions() +void BasicMinimizer::setParameter(size_t index, const FitParameter *par) +{ + bool success; + if (par->isFixed()) { + success = rootMinimizer()->SetFixedVariable((int)index, par->getName().c_str(), + par->getValue()); + + } + + else if (par->hasLowerAndUpperLimits()) { + success = rootMinimizer()->SetLimitedVariable((int)index, par->getName().c_str(), + par->getValue(), par->getStep(), + par->getLowerLimit(), par->getUpperLimit()); + } + + else if (par->hasLowerLimit() && !par->hasUpperLimit()) { + success = rootMinimizer()->SetLowerLimitedVariable((int)index, par->getName().c_str(), + par->getValue(), par->getStep(), + par->getLowerLimit()); + } + + else if (par->hasUpperLimit() && !par->hasLowerLimit()) { + success = rootMinimizer()->SetUpperLimitedVariable((int)index, par->getName().c_str(), + par->getValue(), par->getStep(), + par->getUpperLimit()); + } + + else if (!par->hasUpperLimit() && !par->hasLowerLimit() && !par->isFixed()) { + success = rootMinimizer()->SetVariable((int)index, par->getName().c_str(), par->getValue(), + par->getStep()); + } + + else { + throw std::runtime_error("BasicMinimizer::setParameter() -> Error! Unexpected parameter."); + } + + if( !success ) { + std::ostringstream ostr; + ostr << "BasicMinimizer::setParameter() -> Error! Can't set minimizer's fit parameter"; + ostr << "Index:" << index << " name '" << par->getName() << "'"; + throw std::runtime_error(ostr.str()); + } +} + +void BasicMinimizer::setParameters(const FitSuiteParameters ¶meters) { + size_t index(0); + for (auto par: parameters) + setParameter(index++, par ); + + if( (int)parameters.size() != fitParameterCount()) { + std::ostringstream ostr; + ostr << "BasicMinimizer::setParameters() -> Error! Unconsistency in fit parameter number: "; + ostr << "fitParameterCount = " << fitParameterCount(); + ostr << "parameters.size = " << parameters.size(); + throw std::runtime_error(ostr.str()); + } +} +void BasicMinimizer::setChiSquaredFunction(IMinimizer::function_chi2_t fun_chi2, size_t nparameters) +{ + m_chi2_func.reset(new ROOTMinimizerChiSquaredFunction(fun_chi2, (int)nparameters)); + if( !isGradientBasedAgorithm() ) rootMinimizer()->SetFunction(*m_chi2_func); +} + +void BasicMinimizer::setGradientFunction(IMinimizer::function_gradient_t fun_gradient, size_t nparameters, size_t ndatasize) +{ + m_gradient_func.reset(new ROOTMinimizerGradientFunction(fun_gradient, nparameters, ndatasize)); + if( isGradientBasedAgorithm() ) rootMinimizer()->SetFunction(*m_gradient_func); +} + +std::vector<double> BasicMinimizer::getValueOfVariablesAtMinimum() const +{ + std::vector<double > result; + result.resize(fitParameterCount(), 0.0); + std::copy(rootMinimizer()->X(), rootMinimizer()->X()+fitParameterCount(), result.begin()); + return result; +} + +std::vector<double> BasicMinimizer::getErrorOfVariables() const +{ + std::vector<double > result; + result.resize(fitParameterCount(), 0.0); + if(rootMinimizer()->Errors() != 0 ) { + std::copy(rootMinimizer()->Errors(), rootMinimizer()->Errors()+fitParameterCount(), result.begin()); + } + return result; +} + +void BasicMinimizer::printResults() const +{ + std::cout << toResultString() << std::endl; +} + +std::string BasicMinimizer::toResultString() const +{ +// std::cout << "DONE" << std::endl; +// std::cout << toOptionString() << std::endl; + return toOptionString(); +} + +//! Returns number of fit parameters defined (i.e. dimension of the function to be minimized). + +int BasicMinimizer::fitParameterCount() const +{ + return rootMinimizer()->NDim(); +} + +BasicMinimizer::root_minimizer_t *BasicMinimizer::rootMinimizer() +{ + return const_cast<root_minimizer_t *>(static_cast<const BasicMinimizer*>(this)->rootMinimizer()); } diff --git a/Fit/Minimizer/BasicMinimizer.h b/Fit/Minimizer/BasicMinimizer.h index c91dc226f6bb57679549a57699c973cdcb79103e..bc20658344e9e17e7a09f1b25c77c07831370abe 100644 --- a/Fit/Minimizer/BasicMinimizer.h +++ b/Fit/Minimizer/BasicMinimizer.h @@ -18,6 +18,15 @@ #include "IMinimizer.h" #include <string> +#include <memory> + +class ROOTMinimizerChiSquaredFunction; +class ROOTMinimizerGradientFunction; +namespace BA_ROOT { +namespace Math { +class Minimizer; +} +} //! @class BasicMinimizer //! @ingroup fitting_internal @@ -26,9 +35,14 @@ class BA_CORE_API_ BasicMinimizer : public IMinimizer { public: + typedef BA_ROOT::Math::Minimizer root_minimizer_t; + explicit BasicMinimizer(const std::string &minimizerName, const std::string &algorithmName = std::string()); - virtual ~BasicMinimizer(){} + virtual ~BasicMinimizer(); + + virtual void minimize(); + //! Returns name of the minimizer. std::string minimizerName() const; @@ -39,12 +53,43 @@ public: //! Sets minimization algorithm. void setAlgorithmName(const std::string &algorithmName); + + virtual void setParameter(size_t index, const FitParameter *par); + + virtual void setParameters(const FitSuiteParameters& parameters); + + + //! Sets chi squared function to minimize + virtual void setChiSquaredFunction(function_chi2_t fun_chi2, size_t nparameters); + + //! Sets gradient function to minimize + virtual void setGradientFunction( + function_gradient_t fun_gradient, size_t nparameters, size_t ndatasize); + + + virtual bool isGradientBasedAgorithm() { return false;} + + + virtual std::vector<double> getValueOfVariablesAtMinimum() const; + virtual std::vector<double> getErrorOfVariables() const; + + void printResults() const; + + std::string toResultString() const; + + protected: - virtual void propagateOptions(); + int fitParameterCount() const; + virtual void propagateOptions(){} + virtual const root_minimizer_t* rootMinimizer() const = 0; + root_minimizer_t* rootMinimizer(); private: std::string m_minimizerName; std::string m_algorithmName; + + std::unique_ptr<ROOTMinimizerChiSquaredFunction> m_chi2_func; + std::unique_ptr<ROOTMinimizerGradientFunction> m_gradient_func; }; #endif diff --git a/Fit/Minimizer/Configurable.cpp b/Fit/Minimizer/Configurable.cpp index d50eaeebaf58f244f9fdddc6755641869d0e415c..8844eee1e44a136aec9d5ffde3902db32baad478 100644 --- a/Fit/Minimizer/Configurable.cpp +++ b/Fit/Minimizer/Configurable.cpp @@ -14,6 +14,8 @@ // ************************************************************************** // #include "Configurable.h" +#include <sstream> +#include <variant_io.hpp> //! Returns true if option with such name already exists. Configurable::Configurable(const Configurable &other) @@ -33,7 +35,6 @@ Configurable &Configurable::operator=(const Configurable &other) Configurable::option_t Configurable::option(const std::string &optionName) { -// const_cast<option_t >(static_cast<const Configurable*>(this)->option(optionName)); for(auto option: m_options) { if(option->name() == optionName) return option; @@ -55,6 +56,15 @@ const Configurable::option_t Configurable::option(const std::string &optionName) } +std::string Configurable::toOptionString(const std::string &delimeter) const +{ + std::ostringstream result; + for(auto option: m_options) { + result << option->name() << std::string("=") << option->value() << delimeter; + } + return result.str(); +} + bool Configurable::exists(const std::string &name) { for(auto option: m_options) { diff --git a/Fit/Minimizer/Configurable.h b/Fit/Minimizer/Configurable.h index 8d7fa6ce2ba4aa4d4a92501fc169a8fe0ea0da61..a45778ed8ecdb90d9df9a0220966a9fd03d6fcc3 100644 --- a/Fit/Minimizer/Configurable.h +++ b/Fit/Minimizer/Configurable.h @@ -50,6 +50,9 @@ public: template<class T> void setOptionValue(const std::string& optionName, T value); + //! Returns string with all options using given delimeter + std::string toOptionString(const std::string &delimeter=";") const; + private: bool exists(const std::string &name); void swapContent(Configurable& other); diff --git a/Fit/Minimizer/Minuit2Minimizer.cpp b/Fit/Minimizer/Minuit2Minimizer.cpp index 9735999e65fb2ba38c6586c9ce89908b7ce855c4..c0c7fc1334c8d1e2684a04062a590a32316f78bc 100644 --- a/Fit/Minimizer/Minuit2Minimizer.cpp +++ b/Fit/Minimizer/Minuit2Minimizer.cpp @@ -107,3 +107,8 @@ void Minuit2Minimizer::propagateOptions() m_minuit2_minimizer->SetPrintLevel(printLevel()); } +const BasicMinimizer::root_minimizer_t *Minuit2Minimizer::rootMinimizer() const +{ + return m_minuit2_minimizer.get(); +} + diff --git a/Fit/Minimizer/Minuit2Minimizer.h b/Fit/Minimizer/Minuit2Minimizer.h index 711daca9c70f236c93c2476849dc17f0bb64f15f..17895ca3688a666302ea75796df7f7e3e8b8fc25 100644 --- a/Fit/Minimizer/Minuit2Minimizer.h +++ b/Fit/Minimizer/Minuit2Minimizer.h @@ -68,6 +68,7 @@ public: protected: void propagateOptions(); + const root_minimizer_t* rootMinimizer() const; private: std::unique_ptr<BA_ROOT::Minuit2::Minuit2Minimizer> m_minuit2_minimizer; diff --git a/Fit/Minimizer/variant_io.hpp b/Fit/Minimizer/variant_io.hpp new file mode 100644 index 0000000000000000000000000000000000000000..eb313b37ba9f5b13e6080263d140d97fb5bf9326 --- /dev/null +++ b/Fit/Minimizer/variant_io.hpp @@ -0,0 +1,45 @@ +#ifndef MAPBOX_UTIL_VARIANT_IO_HPP +#define MAPBOX_UTIL_VARIANT_IO_HPP + +#include <iosfwd> + +#include <variant.hpp> + +namespace mapbox { +namespace util { + +namespace detail { +// operator<< helper +template <typename Out> +class printer +{ +public: + explicit printer(Out& out) + : out_(out) {} + printer& operator=(printer const&) = delete; + + // visitor + template <typename T> + void operator()(T const& operand) const + { + out_ << operand; + } + +private: + Out& out_; +}; +} + +// operator<< +template <typename CharT, typename Traits, typename... Types> +VARIANT_INLINE std::basic_ostream<CharT, Traits>& +operator<<(std::basic_ostream<CharT, Traits>& out, variant<Types...> const& rhs) +{ + detail::printer<std::basic_ostream<CharT, Traits>> visitor(out); + apply_visitor(visitor, rhs); + return out; +} +} // namespace util +} // namespace mapbox + +#endif // MAPBOX_UTIL_VARIANT_IO_HPP diff --git a/Tests/Functional/Fit/CMakeLists.txt b/Tests/Functional/Fit/CMakeLists.txt index 215e1f161547db63f41832e39b3c6031146aad3c..1987d45e8f2d03d6e57fca38f160eb6f5d3d620b 100644 --- a/Tests/Functional/Fit/CMakeLists.txt +++ b/Tests/Functional/Fit/CMakeLists.txt @@ -27,6 +27,7 @@ include_directories( ${EIGEN3_INCLUDE_DIR} ${GSL_INCLUDE_DIR} ${CMAKE_CURRENT_SOURCE_DIR}/../TestMachinery + ${RootMinimizers_INCLUDE_DIRS} ) #file(GLOB source_files "*.cpp") diff --git a/Tests/Functional/Fit/ExperimentalFitTest.cpp b/Tests/Functional/Fit/ExperimentalFitTest.cpp index 0bdccc17f8f62e7666032814492a6551f06d5c90..02f74db2de6943221c5f190f3423bda1fe46a491 100644 --- a/Tests/Functional/Fit/ExperimentalFitTest.cpp +++ b/Tests/Functional/Fit/ExperimentalFitTest.cpp @@ -14,13 +14,27 @@ // ************************************************************************** // #include "ExperimentalFitTest.h" -#include "GISASSimulation.h" -#include "Histogram2D.h" -#include "Rectangle.h" -#include "RectangularDetector.h" -#include "Units.h" +#include "Minuit2Minimizer.h" +#include "FitSuite.h" ExperimentalFitTest::ExperimentalFitTest() : IMinimizerTest("Minuit2", "Migrad") { } + +std::unique_ptr<FitSuite> ExperimentalFitTest::createFitSuite() +{ + std::unique_ptr<FitSuite> result(new FitSuite()); + result->initPrint(10); + IMinimizer* minimizer = new Minuit2Minimizer(); + + result->setMinimizer(minimizer); + + for (size_t i = 0; i < m_parameters.size(); ++i) { + result->addFitParameter( + m_parameters[i].m_name, m_parameters[i].m_start_value, + AttLimits::lowerLimited(0.01), m_parameters[i].m_start_value / 100.); + } + return result; + +} diff --git a/Tests/Functional/Fit/ExperimentalFitTest.h b/Tests/Functional/Fit/ExperimentalFitTest.h index 735d209d637bb6cfa19d9d706a60cac8440ec5e0..5df4e31ee004697f0a9d9e7c87cc1d948f5f8da5 100644 --- a/Tests/Functional/Fit/ExperimentalFitTest.h +++ b/Tests/Functional/Fit/ExperimentalFitTest.h @@ -27,6 +27,9 @@ class ExperimentalFitTest : public IMinimizerTest public: ExperimentalFitTest(); +protected: + virtual std::unique_ptr<FitSuite> createFitSuite(); + }; #endif // EXPERIMENTALFITTEST_H diff --git a/Tests/Functional/Fit/IMinimizerTest.cpp b/Tests/Functional/Fit/IMinimizerTest.cpp index f792ee3dd29042996e01b89ece1c8e2bb138678e..84146e1ffa65e8a6270eb9fa417d2b563d07492c 100644 --- a/Tests/Functional/Fit/IMinimizerTest.cpp +++ b/Tests/Functional/Fit/IMinimizerTest.cpp @@ -64,8 +64,9 @@ void IMinimizerTest::runTest() // run fit fitSuite->runFit(); + std::vector<double> valuesAtMinimum = fitSuite->getMinimizer()->getValueOfVariablesAtMinimum(); for (size_t i = 0; i < m_parameters.size(); ++i) - m_parameters[i].m_found_value = fitSuite->getMinimizer()->getValueOfVariableAtMinimum(i); + m_parameters[i].m_found_value = valuesAtMinimum[i]; // analyze results for (size_t i = 0; i < m_parameters.size(); ++i) {