diff --git a/NEWS b/NEWS index 929c5d6..a6b8670 100644 --- a/NEWS +++ b/NEWS @@ -1,6 +1,6 @@ New Features and Important Changes in GammaLib 1.7.0 -17 February 2020 +9 April 2020 1. Introduction @@ -12,7 +12,7 @@ made since the last release of GammaLib. 2. Public interface modifications --------------------------------- The following classes have been added: -- none +- GResponseCache The following classes have been removed: - GCTAResponseCache::ndirs() @@ -49,7 +49,7 @@ The following methods have been removed: - none The following methods have been renamed: -- none +- GInstDir::hash() - Returns a cached value The arguments for the following methods have been changed: - GTime::utc(int&) - Add optional precision value @@ -99,6 +99,12 @@ None 10. Observation module ---------------------- +The GResponseCache class was added as a copy of GCTAResponseCache to implement +a generic response cache for response computation. + +An GInstDir::hash() method was added that returns a unique double precision +value for every class instance that can be used for response caching. + Add GObservations::npred() methods that compute the number of predicted events for a model with a given name or for a given reference to a model. Optionally an instrument string can be passed so that the number of predicted events are diff --git a/include/GResponseCache.hpp b/include/GResponseCache.hpp new file mode 100644 index 0000000..f79da85 --- /dev/null +++ b/include/GResponseCache.hpp @@ -0,0 +1,149 @@ +/*************************************************************************** + * GResponseCache.hpp - Response cache class * + * ----------------------------------------------------------------------- * + * copyright (C) 2020 by Juergen Knoedlseder * + * ----------------------------------------------------------------------- * + * * + * This program 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 3 of the License, or * + * (at your option) any later version. * + * * + * This program 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 program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GResponseCache.hpp + * @brief Response cache class definition + * @author Juergen Knoedlseder + */ + +#ifndef GRESPONSECACHE_HPP +#define GRESPONSECACHE_HPP + +/* __ Includes ___________________________________________________________ */ +#include +#include +#include "GBase.hpp" +#include "GEnergy.hpp" + +/* __ Forward declarations _______________________________________________ */ +class GInstDir; + +/* __ Constants __________________________________________________________ */ + +/* __ Typedefs ___________________________________________________________ */ +typedef std::map GResponseCacheElement; +typedef std::map GResponseCacheName; + + +/***********************************************************************//** + * @class GResponseCache + * + * @brief Response cache class + * + * The class implements a cache for the Instrument Response Function values + * so that the values do not need to be recomputed each time but can be + * fetched from the cache. + ***************************************************************************/ +class GResponseCache : public GBase { + +public: + // Constructors and destructors + GResponseCache(void); + GResponseCache(const GResponseCache& cache); + virtual ~GResponseCache(void); + + // Operators + GResponseCache& operator=(const GResponseCache& cache); + + // Methods + void clear(void); + GResponseCache* clone(void) const; + std::string classname(void) const; + int size(void) const; + bool is_empty(void) const; + void set(const std::string& name, + const GEnergy& ereco, + const GEnergy& etrue, + const double& value); + void set(const std::string& name, + const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue, + const double& value); + void remove(const std::string& name); + bool contains(const std::string& name, + const GEnergy& ereco, + const GEnergy& etrue, + double* value = NULL) const; + bool contains(const std::string& name, + const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue, + double* value = NULL) const; + std::string print(const GChatter& chatter = NORMAL) const; + +protected: + // Protected methods + void init_members(void); + void copy_members(const GResponseCache& cache); + void free_members(void); + double encode(const GEnergy& ereco, + const GEnergy& etrue) const; + double encode(const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue) const; + + // Protected members + GResponseCacheName m_cache; +}; + + +/***********************************************************************//** + * @brief Return class name + * + * @return String containing the class name ("GResponseCache"). + ***************************************************************************/ +inline +std::string GResponseCache::classname(void) const +{ + return ("GResponseCache"); +} + + +/***********************************************************************//** + * @brief Checks whether the cache is empty + * + * @return True if cache is empty, false otherwise. + * + * Checks whether the response cache is empty. + ***************************************************************************/ +inline +bool GResponseCache::is_empty(void) const +{ + return (m_cache.empty()); +} + + +/***********************************************************************//** + * @brief Remove cache for source + * + * @param[in] name Source name. + * + * Remove cache for source with @p name. + ***************************************************************************/ +inline +void GResponseCache::remove(const std::string& name) +{ + m_cache.erase(name); + return; +} + +#endif /* GRESPONSECACHE_HPP */ diff --git a/include/GammaLib.hpp b/include/GammaLib.hpp index a422250..11c0343 100644 --- a/include/GammaLib.hpp +++ b/include/GammaLib.hpp @@ -1,7 +1,7 @@ /*************************************************************************** * GammaLib.hpp - Gamma-Ray Astronomy Library Header file * * ----------------------------------------------------------------------- * - * copyright (C) 2006-2018 by Juergen Knoedlseder * + * copyright (C) 2006-2020 by Juergen Knoedlseder * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * @@ -188,6 +188,7 @@ #include "GPhases.hpp" #include "GInstDir.hpp" #include "GResponse.hpp" +#include "GResponseCache.hpp" #include "GPhotons.hpp" #include "GPhoton.hpp" #include "GSource.hpp" diff --git a/include/Makefile.am b/include/Makefile.am index 3618f1f..6a848bc 100644 --- a/include/Makefile.am +++ b/include/Makefile.am @@ -1,7 +1,7 @@ ############################################################################# # Makefile.am for include module # # ------------------------------------------------------------------------- # -# Copyright (C) 2016-2018 Juergen Knoedlseder # +# Copyright (C) 2016-2020 Juergen Knoedlseder # # ------------------------------------------------------------------------- # # # # This program is free software: you can redistribute it and/or modify # @@ -148,6 +148,7 @@ pkginclude_HEADERS = GammaLib.hpp \ GPhases.hpp \ GInstDir.hpp \ GResponse.hpp \ + GResponseCache.hpp \ GPhotons.hpp \ GPhoton.hpp \ GSource.hpp \ diff --git a/pyext/GResponseCache.i b/pyext/GResponseCache.i new file mode 100644 index 0000000..c4707c8 --- /dev/null +++ b/pyext/GResponseCache.i @@ -0,0 +1,92 @@ +/*************************************************************************** + * GResponseCache.i - Response cache class * + * ----------------------------------------------------------------------- * + * copyright (C) 2020 by Juergen Knoedlseder * + * ----------------------------------------------------------------------- * + * * + * This program 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 3 of the License, or * + * (at your option) any later version. * + * * + * This program 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 program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GResponseCache.i + * @brief Response cache class definition + * @author Juergen Knoedlseder + */ +%{ +/* Put headers and other declarations here that are needed for compilation */ +#include "GResponseCache.hpp" +%} + +/* __ Typemaps ___________________________________________________________ */ +/* + * Typemap so that the contains() method returns two values, first the + * flag indicating whether the value exists or not, and second the + * value, i.e. + * + * flag, value = cache.contains('Crab', GEnergy(1.0,'TeV'), GEnergy(1.0,'TeV')) + * + */ +%include "typemaps.i" +%apply double *OUTPUT { double* value }; + + +/***********************************************************************//** + * @class GResponseCache + * + * @brief Response cache class + ***************************************************************************/ +class GResponseCache : public GBase { + +public: + // Constructors and destructors + GResponseCache(void); + GResponseCache(const GResponseCache& cache); + virtual ~GResponseCache(void); + + // Methods + void clear(void); + GResponseCache* clone(void) const; + std::string classname(void) const; + int size(void) const; + bool is_empty(void) const; + void set(const std::string& name, + const GEnergy& ereco, + const GEnergy& etrue, + const double& value); + void set(const std::string& name, + const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue, + const double& value); + void remove(const std::string& name); + bool contains(const std::string& name, + const GEnergy& ereco, + const GEnergy& etrue, + double* value) const; + bool contains(const std::string& name, + const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue, + double* value) const; +}; + + +/***********************************************************************//** + * @brief GResponseCache class extension + ***************************************************************************/ +%extend GResponseCache { + GResponseCache copy() { + return (*self); + } +}; diff --git a/pyext/gammalib/obs.i b/pyext/gammalib/obs.i index a61ce4b..daf57b9 100644 --- a/pyext/gammalib/obs.i +++ b/pyext/gammalib/obs.i @@ -1,7 +1,7 @@ /*************************************************************************** * obs.i - Observation module * * ----------------------------------------------------------------------- * - * copyright (C) 2011-2017 by Juergen Knoedlseder * + * copyright (C) 2011-2020 by Juergen Knoedlseder * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * @@ -196,6 +196,7 @@ %include "GRoi.i" %include "GPhases.i" %include "GResponse.i" +%include "GResponseCache.i" %include "GPhoton.i" %include "GPhotons.i" %include "GSource.i" diff --git a/src/obs/GResponseCache.cpp b/src/obs/GResponseCache.cpp new file mode 100644 index 0000000..40e8103 --- /dev/null +++ b/src/obs/GResponseCache.cpp @@ -0,0 +1,486 @@ +/*************************************************************************** + * GResponseCache.cpp - Response cache class * + * ----------------------------------------------------------------------- * + * copyright (C) 2020 by Juergen Knoedlseder * + * ----------------------------------------------------------------------- * + * * + * This program 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 3 of the License, or * + * (at your option) any later version. * + * * + * This program 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 program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GResponseCache.cpp + * @brief Response cache class implementation + * @author Juergen Knoedlseder + */ + +/* __ Includes ___________________________________________________________ */ +#ifdef HAVE_CONFIG_H +#include +#endif +#include "GTools.hpp" +#include "GResponseCache.hpp" +#include "GInstDir.hpp" + +/* __ Method name definitions ____________________________________________ */ + +/* __ Macros _____________________________________________________________ */ + +/* __ Coding definitions _________________________________________________ */ + +/* __ Debug definitions __________________________________________________ */ + + +/*========================================================================== + = = + = Constructors/destructors = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Void constructor + ***************************************************************************/ +GResponseCache::GResponseCache(void) +{ + // Initialise class members + init_members(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Copy constructor + * + * @param[in] cache Response cache. + ***************************************************************************/ +GResponseCache::GResponseCache(const GResponseCache& cache) +{ + // Initialise class members + init_members(); + + // Copy members + copy_members(cache); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Destructor + ***************************************************************************/ +GResponseCache::~GResponseCache(void) +{ + // Free members + free_members(); + + // Return + return; +} + + +/*========================================================================== + = = + = Operators = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Assignment operator + * + * @param[in] cache Response cache. + * @return Response cache. + ***************************************************************************/ +GResponseCache& GResponseCache::operator=(const GResponseCache& cache) +{ + // Execute only if object is not identical + if (this != &cache) { + + // Free members + free_members(); + + // Initialise private members + init_members(); + + // Copy members + copy_members(cache); + + } // endif: object was not identical + + // Return this object + return *this; +} + + +/*========================================================================== + = = + = Public methods = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Clear response cache + ***************************************************************************/ +void GResponseCache::clear(void) +{ + // Free members + free_members(); + + // Initialise private members + init_members(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Return number of elements in cache + * + * @return Number of elements in cache + * + * Returns the number of elements in the response cache. + ***************************************************************************/ +int GResponseCache::size(void) const +{ + // Initialize size + int size = 0; + + // Compute size + for (GResponseCacheName::const_iterator it_name = m_cache.begin(); + it_name != m_cache.end(); ++it_name) { + size += it_name->second.size(); + } + + // Return size + return size; +} + + +/***********************************************************************//** + * @brief Set cache value + * + * @param[in] name Cache name. + * @param[in] ereco Reconstructed energy. + * @param[in] etrue True energy. + * @param[in] value Cache value. + * + * Set cache value for a given @p name, reconstructed energy @p ereco, and + * true energy @p etrue. + ***************************************************************************/ +void GResponseCache::set(const std::string& name, + const GEnergy& ereco, + const GEnergy& etrue, + const double& value) +{ + // Get element identifier + double element = encode(ereco, etrue); + + // Set cache value + m_cache[name][element] = value; + + // Return + return; +} + + +/***********************************************************************//** + * @brief Set cache value + * + * @param[in] name Cache name. + * @param[in] dir Instrument direction. + * @param[in] ereco Reconstructed energy. + * @param[in] etrue True energy. + * @param[in] value Cache value. + * + * Set cache value for a given @p name, instrument direction @p dir, + * reconstructed energy @p ereco, and true energy @p etrue. + ***************************************************************************/ +void GResponseCache::set(const std::string& name, + const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue, + const double& value) +{ + // Get element identifier + double element = encode(dir, ereco, etrue); + + // Set cache value + m_cache[name][element] = value; + + // Return + return; +} + + +/***********************************************************************//** + * @brief Check if cache contains a value for specific parameters + * + * @param[in] name Cache name. + * @param[in] ereco Reconstructed energy. + * @param[in] etrue True energy. + * @param[out] value Pointer to cached value (only if found). + * @return True if cached value was found, otherwise false. + * + * Check if the cache contains a value for a given @p name, reconstructed + * energy @p ereco, and true energy @p etrue. + * + * If the @p value pointer argument is not NULL, the method will return the + * cached value through this argument in case that the value exists. + ***************************************************************************/ +bool GResponseCache::contains(const std::string& name, + const GEnergy& ereco, + const GEnergy& etrue, + double* value) const +{ + // Initialise containment flag + bool contains = false; + + // Get element identifier + double element = encode(ereco, etrue); + + // Search for name in cache + GResponseCacheName::const_iterator it_name = m_cache.find(name); + if (it_name != m_cache.end()) { + + // Search for element in cache + GResponseCacheElement::const_iterator it_element = + it_name->second.find(element); + if (it_element != it_name->second.end()) { + contains = true; + if (value != NULL) { + *value = it_element->second; + } + } // endif: element found + + } // endif: name found + + // Return containment flag + return contains; +} + + +/***********************************************************************//** + * @brief Check if cache contains a value for specific parameters + * + * @param[in] name Cache name. + * @param[in] dir Instrument direction. + * @param[in] ereco Reconstructed energy. + * @param[in] etrue True energy. + * @param[out] value Pointer to cached value (only if found). + * @return True if cached value was found, otherwise false. + * + * Check if the cache contains a value for a given @p name, instrument + * direction @p dir, reconstructed energy @p ereco, and true energy @p etrue. + * + * If the @p value pointer argument is not NULL, the method will return the + * cached value through this argument in case that the value exists. + ***************************************************************************/ +bool GResponseCache::contains(const std::string& name, + const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue, + double* value) const +{ + // Initialise containment flag + bool contains = false; + + // Get element identifier + double element = encode(dir, ereco, etrue); + + // Search for name in cache + GResponseCacheName::const_iterator it_name = m_cache.find(name); + if (it_name != m_cache.end()) { + + // Search for element in cache + GResponseCacheElement::const_iterator it_element = + it_name->second.find(element); + if (it_element != it_name->second.end()) { + contains = true; + if (value != NULL) { + *value = it_element->second; + } + } // endif: element found + + } // endif: name found + + // Return containment flag + return contains; +} + + +/***********************************************************************//** + * @brief Clone response cache + * + * @return Pointer to deep copy of response cache. + ***************************************************************************/ +GResponseCache* GResponseCache::clone(void) const +{ + return new GResponseCache(*this); +} + + +/***********************************************************************//** + * @brief Print response cache + * + * @param[in] chatter Chattiness. + * @return String containing response cache information. + ***************************************************************************/ +std::string GResponseCache::print(const GChatter& chatter) const +{ + // Initialise result string + std::string result; + + // Continue only if chatter is not silent + if (chatter != SILENT) { + + // Append header + result.append("=== GResponseCache ==="); + + // Append total cache size + result.append("\n"+gammalib::parformat("Number of cached values")); + result.append(gammalib::str(size())); + + // Append information + for (GResponseCacheName::const_iterator it_name = m_cache.begin(); + it_name != m_cache.end(); ++it_name) { + + // Append name + result.append("\n"+gammalib::parformat("Name")+it_name->first); + + // Compute number of elements + int nelement = it_name->second.size(); + + // Append number of elements + result.append("\n"+gammalib::parformat("Elements")); + result.append(gammalib::str(nelement)); + + } // endfor: looped over names + + } // endif: chatter was not silent + + // Return result + return result; +} + + +/*========================================================================== + = = + = Private methods = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Initialise class members + ***************************************************************************/ +void GResponseCache::init_members(void) +{ + // Initialise members + m_cache.clear(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Copy class members + * + * @param[in] cache Response cache. + ***************************************************************************/ +void GResponseCache::copy_members(const GResponseCache& cache) +{ + // Copy members + m_cache = cache.m_cache; + + // Return + return; +} + + +/***********************************************************************//** + * @brief Delete class members + ***************************************************************************/ +void GResponseCache::free_members(void) +{ + // Return + return; +} + + +/***********************************************************************//** + * @brief Encode reconstructued and true energy + * + * @param[in] ereco Reconstructed energy. + * @param[in] etrue True energy. + * @return Encoded reconstructued and true energy + * + * Encodes the reconstructued and true energy in a single double precision + * value. The encoding is done using the following formula: + * + * \f[ + * {\tt encode} = E_{\rm reco} + E_{\rm true} \times 10^2 + * \f] + * + * where + * - \f$E_{\rm reco}\f$ is the reconstructued energy in TeV, and + * - \f$E_{\rm true}\f$ is the true energy in TeV. + ***************************************************************************/ +double GResponseCache::encode(const GEnergy& ereco, + const GEnergy& etrue) const +{ + // Construct encoded instrument direction and reconstructued energy + double encoded = ereco.MeV() + etrue.MeV() * 1.0e2; + + // Return encoded value + return encoded; +} + + +/***********************************************************************//** + * @brief Encode instrument direction, reconstructued and true energy + * + * @param[in] dir Instrument direction. + * @param[in] ereco Reconstructed energy. + * @param[in] etrue True energy. + * @return Encoded instrument direction, reconstructued and true energy + * + * Encodes the instrument direction, reconstructued and true energy in a + * single double precision value. The encoding is done using the following + * formula: + * + * \f[ + * {\tt encode} = {\tt dir.hash()} + E_{\rm reco} \times 10^6 + + E_{\rm true} \times 10^8 + * \f] + * + * where + * - \f${\tt dir.encode()}\f$ is an encoding value returned by @p dir, + * - \f$\delta\f$ is the Declination of the instrument direction in degrees, + * - \f$E_{\rm reco}\f$ is the reconstructued energy in TeV, and + * - \f$E_{\rm true}\f$ is the true energy in TeV. + * + * @todo Use instrument direction encode method + ***************************************************************************/ +double GResponseCache::encode(const GInstDir& dir, + const GEnergy& ereco, + const GEnergy& etrue) const +{ + // Construct encoded instrument direction and reconstructued energy + double encoded = dir.hash() + ereco.MeV() * 1.0e6 + etrue.MeV() * 1.0e8; + + // Return encoded value + return encoded; +} diff --git a/src/obs/Makefile.am b/src/obs/Makefile.am index 0041332..0d5e78f 100644 --- a/src/obs/Makefile.am +++ b/src/obs/Makefile.am @@ -1,7 +1,7 @@ ############################################################################# # Makefile.am for obs module # # ------------------------------------------------------------------------- # -# Copyright (C) 2016-2017 Juergen Knoedlseder # +# Copyright (C) 2016-2020 Juergen Knoedlseder # # ------------------------------------------------------------------------- # # # # This program is free software: you can redistribute it and/or modify # @@ -46,6 +46,7 @@ sources = GEnergies.cpp \ GEbounds.cpp \ GPhases.cpp \ GResponse.cpp \ + GResponseCache.cpp \ GInstDir.cpp \ GPhotons.cpp \ GPhoton.cpp \