From 82c53c6fa15eba93ee1967bab649afa6d901f71e Mon Sep 17 00:00:00 2001 From: Jürgen Knödlseder Date: Tue, 7 Apr 2020 22:39:31 +0200 Subject: [PATCH] Add and use GModel::has_eval_indices() and GModel::eval_indices() methods (#745) --- include/GModel.hpp | 80 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++------------------------ pyext/GModel.i | 44 +++++++++++++++++++++++--------------------- src/model/GModel.cpp | 34 +++++++++++++++++++--------------- src/obs/GObservation.cpp | 103 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++------------------------------- 4 files changed, 170 insertions(+), 91 deletions(-) diff --git a/include/GModel.hpp b/include/GModel.hpp index 29246df..4b73c42 100644 --- a/include/GModel.hpp +++ b/include/GModel.hpp @@ -1,7 +1,7 @@ /*************************************************************************** * GModel.hpp - Abstract virtual model base class * * ----------------------------------------------------------------------- * - * copyright (C) 2009-2019 by Juergen Knoedlseder * + * copyright (C) 2009-2020 by Juergen Knoedlseder * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * @@ -126,29 +126,31 @@ public: virtual std::string print(const GChatter& chatter = NORMAL) const = 0; // Implemented methods - int size(void) const; - int scales(void) const; - GModelPar& at(const int& index); - const GModelPar& at(const int& index) const; - bool has_par(const std::string& name) const; - bool has_scales(void) const; - const std::string& name(void) const; - void name(const std::string& name); - const double& ts(void) const; - void ts(const double& ts); - const bool& tscalc(void) const; - void tscalc(const bool& tscalc); - const bool& has_ts(void) const; - std::string instruments(void) const; - void instruments(const std::string& instruments); - GModelPar& scale(const int& index); - const GModelPar& scale(const int& index) const; - GModelPar scale(const std::string& instrument) const; - void scale(const GModelPar& par); - std::string ids(void) const; - void ids(const std::string& ids); - bool is_valid(const std::string& instruments, - const std::string& ids) const; + int size(void) const; + int scales(void) const; + GModelPar& at(const int& index); + const GModelPar& at(const int& index) const; + bool has_par(const std::string& name) const; + bool has_scales(void) const; + const std::string& name(void) const; + void name(const std::string& name); + const double& ts(void) const; + void ts(const double& ts); + const bool& tscalc(void) const; + void tscalc(const bool& tscalc); + const bool& has_ts(void) const; + std::string instruments(void) const; + void instruments(const std::string& instruments); + GModelPar& scale(const int& index); + const GModelPar& scale(const int& index) const; + GModelPar scale(const std::string& instrument) const; + void scale(const GModelPar& par); + std::string ids(void) const; + void ids(const std::string& ids); + bool is_valid(const std::string& instruments, + const std::string& ids) const; + const bool& has_eval_indices(void) const; + const std::vector& eval_indices(void) const; protected: // Protected methods @@ -171,6 +173,12 @@ protected: bool m_has_tscalc; //!< Signals if tscalc attribute is available bool m_tscalc; //!< Signals if TS should be computed double m_ts; //!< Test Statistic of the model + + // Protected members for efficient gradient access + mutable bool m_has_eval_inx; //!< Signal that parameter indices + //!< updated by eval() method exist + mutable std::vector m_eval_inx; //!< List of parameter indices updated + //!< by eval() method }; @@ -352,4 +360,28 @@ const bool& GModel::has_ts(void) const return m_has_ts; } + +/***********************************************************************//** + * @brief Signals that parameter indices updated by eval() method were set + * + * @return True if parameter indices updated by eval() method were set. + ***************************************************************************/ +inline +const bool& GModel::has_eval_indices(void) const +{ + return m_has_eval_inx; +} + + +/***********************************************************************//** + * @brief Return vector of parameter indices updated by eval() method + * + * @return Vector of parameter indices updated by eval() method. + ***************************************************************************/ +inline +const std::vector& GModel::eval_indices(void) const +{ + return m_eval_inx; +} + #endif /* GMODEL_HPP */ diff --git a/pyext/GModel.i b/pyext/GModel.i index 9c4d3ac..e712fe2 100644 --- a/pyext/GModel.i +++ b/pyext/GModel.i @@ -1,7 +1,7 @@ /*************************************************************************** * GModel.i - Abstract virtual model base class * * ----------------------------------------------------------------------- * - * copyright (C) 2009-2019 by Juergen Knoedlseder * + * copyright (C) 2009-2020 by Juergen Knoedlseder * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * @@ -59,26 +59,28 @@ public: virtual void write(GXmlElement& xml) const = 0; // Implemented methods - int size(void) const; - int scales(void) const; - bool has_par(const std::string& name) const; - bool has_scales(void) const; - const std::string& name(void) const; - void name(const std::string& name); - const double& ts(void) const; - void ts(const double& ts); - const bool& tscalc(void) const; - void tscalc(const bool& tscalc); - const bool& has_ts(void) const; - std::string instruments(void) const; - void instruments(const std::string& instruments); - GModelPar& scale(const int& index); - GModelPar scale(const std::string& instrument) const; - void scale(const GModelPar& par); - std::string ids(void) const; - void ids(const std::string& ids); - bool is_valid(const std::string& instruments, - const std::string& ids) const; + int size(void) const; + int scales(void) const; + bool has_par(const std::string& name) const; + bool has_scales(void) const; + const std::string& name(void) const; + void name(const std::string& name); + const double& ts(void) const; + void ts(const double& ts); + const bool& tscalc(void) const; + void tscalc(const bool& tscalc); + const bool& has_ts(void) const; + std::string instruments(void) const; + void instruments(const std::string& instruments); + GModelPar& scale(const int& index); + GModelPar scale(const std::string& instrument) const; + void scale(const GModelPar& par); + std::string ids(void) const; + void ids(const std::string& ids); + bool is_valid(const std::string& instruments, + const std::string& ids) const; + const bool& has_eval_indices(void) const; + const std::vector& eval_indices(void) const; }; diff --git a/src/model/GModel.cpp b/src/model/GModel.cpp index 3c476d2..ae88b3d 100644 --- a/src/model/GModel.cpp +++ b/src/model/GModel.cpp @@ -1,7 +1,7 @@ /*************************************************************************** * GModel.cpp - Abstract virtual model base class * * ----------------------------------------------------------------------- * - * copyright (C) 2009-2019 by Juergen Knoedlseder * + * copyright (C) 2009-2020 by Juergen Knoedlseder * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * @@ -617,11 +617,13 @@ void GModel::init_members(void) m_scales.clear(); m_ids.clear(); m_pars.clear(); - m_ts = 0.0; - m_has_ts = false; - m_has_tscalc = false; - m_tscalc = false; - + m_ts = 0.0; + m_has_ts = false; + m_has_tscalc = false; + m_tscalc = false; + m_has_eval_inx = false; + m_eval_inx.clear(); + // Return return; } @@ -635,15 +637,17 @@ void GModel::init_members(void) void GModel::copy_members(const GModel& model) { // Copy members - m_name = model.m_name; - m_instruments = model.m_instruments; - m_scales = model.m_scales; - m_ids = model.m_ids; - m_pars = model.m_pars; - m_ts = model.m_ts; - m_has_ts = model.m_has_ts; - m_has_tscalc = model.m_has_tscalc; - m_tscalc = model.m_tscalc; + m_name = model.m_name; + m_instruments = model.m_instruments; + m_scales = model.m_scales; + m_ids = model.m_ids; + m_pars = model.m_pars; + m_ts = model.m_ts; + m_has_ts = model.m_has_ts; + m_has_tscalc = model.m_has_tscalc; + m_tscalc = model.m_tscalc; + m_has_eval_inx = model.m_has_eval_inx; + m_eval_inx = model.m_eval_inx; // Return return; diff --git a/src/obs/GObservation.cpp b/src/obs/GObservation.cpp index fe7112a..dcc5f38 100644 --- a/src/obs/GObservation.cpp +++ b/src/obs/GObservation.cpp @@ -1,7 +1,7 @@ /*************************************************************************** * GObservation.cpp - Abstract observation base class * * ----------------------------------------------------------------------- * - * copyright (C) 2008-2019 by Juergen Knoedlseder * + * copyright (C) 2008-2020 by Juergen Knoedlseder * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * @@ -302,39 +302,80 @@ double GObservation::model(const GModels& models, // over the energy dispersion. Here it's simpler to just use // numerical gradients. if (gradient != NULL) { - for (int ipar = 0; ipar < mptr->size(); ++ipar) { - - // Get reference to model parameter - const GModelPar& par = (*mptr)[ipar]; - - // Make sure that we have a slot for the gradient - #if defined(G_RANGE_CHECK) - if (igrad >= grad_size) { - std::string msg = "Vector has not enough elements " - "to store the model parameter " - "gradients. "+ - gammalib::str(models.npars())+ - " elements requested while vector " - "only contains "+ - gammalib::str(gradient->size())+ - " elements."; - throw GException::invalid_value(G_MODEL, msg); - } - #endif - - if (par.is_free()) { - if (par.has_grad() && !response()->use_edisp()) { - (*gradient)[igrad+ipar] = par.factor_gradient(); + + // Make sure that we have a slot for the gradient + #if defined(G_RANGE_CHECK) + if (igrad+mptr->size() > grad_size) { + std::string msg = "Vector has not enough elements to " + "store the model parameter gradients. "+ + gammalib::str(models.npars())+" elements " + "requested while vector only contains "+ + gammalib::str(gradient->size())+" elements."; + throw GException::invalid_value(G_MODEL, msg); + } + #endif + + // If model provides the parameter indices that were + // updated by the eval() method then use them as this is + // less time consuming than looping over all indices + if (mptr->has_eval_indices()) { + + // Get number of relevant parameters + int npars = mptr->eval_indices().size(); + + // Loop over all relevant parameters + for (int index = 0; index < npars; ++index) { + + // Retrieve parameter index + int ipar = mptr->eval_indices()[index]; + + // Get reference to model parameter + const GModelPar& par = (*mptr)[ipar]; + + // Set gradient + if (par.is_free()) { + if (par.has_grad() && !response()->use_edisp()) { + (*gradient)[igrad+ipar] = par.factor_gradient(); + } + else { + (*gradient)[igrad+ipar] = model_grad(*mptr, par, event); + } } else { - (*gradient)[igrad+ipar] = model_grad(*mptr, par, event); + (*gradient)[igrad+ipar] = 0.0; } - } - else { - (*gradient)[igrad+ipar] = 0.0; - } - } - } + + } // endfor: looped over all parameter indices + + } // endif: parameter indices were available + + // ... otherwise loop over all parameter indices + else { + + // Loop over all parameters + for (int ipar = 0; ipar < mptr->size(); ++ipar) { + + // Get reference to model parameter + const GModelPar& par = (*mptr)[ipar]; + + // Set gradient + if (par.is_free()) { + if (par.has_grad() && !response()->use_edisp()) { + (*gradient)[igrad+ipar] = par.factor_gradient(); + } + else { + (*gradient)[igrad+ipar] = model_grad(*mptr, par, event); + } + } + else { + (*gradient)[igrad+ipar] = 0.0; + } + + } // endfor: looped over all parameter indices + + } // endelse: no parameter indices were available + + } // endif: use model gradients } // endif: model component was valid for instrument -- libgit2 0.21.4