From 5624e7c0481d59fec8345895df18c4124c7a2ebd Mon Sep 17 00:00:00 2001 From: knodlseder <> Date: Sun, 26 Dec 2010 20:06:08 +0000 Subject: [PATCH] remove non-generic response methods --- include/GResponse.hpp | 40 ++++++++++++---------------------------- inst/cta/include/GCTAResponse.hpp | 33 ++++++++++++++++++++------------- inst/cta/pyext/GCTAResponse.i | 21 +++++++++++++-------- inst/cta/src/GCTAResponse.cpp | 292 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++--------------------------------------------------------------------------------------------- inst/lat/include/GLATResponse.hpp | 68 ++++++++++++++++++++++++++++++-------------------------------------- inst/lat/pyext/GLATResponse.i | 39 +++++++++++++++++++-------------------- inst/lat/src/GLATResponse.cpp | 85 +++++++++++++++++++++++++++++++++++++++++++++++++++++-------------------------------- inst/mwl/include/GMWLResponse.hpp | 9 +++------ inst/mwl/src/GMWLResponse.cpp | 4 ++-- pyext/GResponse.i | 22 ++++++---------------- src/obs/GResponse.cpp | 40 ++-------------------------------------- 11 files changed, 359 insertions(+), 294 deletions(-) diff --git a/include/GResponse.hpp b/include/GResponse.hpp index b7fbdc8..2ec9875 100644 --- a/include/GResponse.hpp +++ b/include/GResponse.hpp @@ -12,7 +12,7 @@ ***************************************************************************/ /** * @file GResponse.hpp - * @brief Abstract response base class definition. + * @brief Abstract response base class definition * @author J. Knodlseder */ @@ -24,12 +24,6 @@ #include #include "GLog.hpp" #include "GEvent.hpp" -#include "GPointing.hpp" -#include "GInstDir.hpp" -#include "GRoi.hpp" -#include "GGti.hpp" -#include "GEbounds.hpp" -#include "GSkyDir.hpp" #include "GEnergy.hpp" #include "GTime.hpp" @@ -43,10 +37,14 @@ class GObservation; * * @brief Abstract interface for the instrument response function * - * The response function provides conversion methods between physical - * parameters (such as source position, flux, ...) and the measured - * instrumental parameters (such as measured energy, photon interaction, - * ...). + * The response function provides conversion between physical parameters + * (such as source position, flux, ...) and the measured instrumental + * parameters (such as measured energy, photon interaction, ...). + * For a given observation, the irf method returns the instrument response + * for a given event and source model as function of the true photon energy + * and the true photon arrival time. + * The npred method returns the integral of the instrument response function + * over the dataspace. This method is only required for unbinned analysis. ***************************************************************************/ class GResponse { @@ -66,35 +64,21 @@ public: // Pure virtual methods virtual void clear(void) = 0; virtual GResponse* clone(void) const = 0; - virtual void load(const std::string& irfname) = 0; virtual bool hasedisp(void) const = 0; virtual bool hastdisp(void) const = 0; - virtual double irf(const GInstDir& obsDir, const GEnergy& obsEng, - const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const = 0; virtual double irf(const GEvent& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const = 0; - virtual double nirf(const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const = 0; + virtual double npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const = 0; virtual std::string print(void) const = 0; - // Other methods - virtual void caldb(const std::string& caldb); - virtual std::string caldb(void) const { return m_caldb; } - protected: // Protected methods void init_members(void); void copy_members(const GResponse& rsp); void free_members(void); - - // Protected data area - std::string m_caldb; //!< Name of or path to the calibration database - std::string m_rspname; //!< Name of the instrument response }; #endif /* GRESPONSE_HPP */ diff --git a/inst/cta/include/GCTAResponse.hpp b/inst/cta/include/GCTAResponse.hpp index c42f880..9b8dd21 100644 --- a/inst/cta/include/GCTAResponse.hpp +++ b/inst/cta/include/GCTAResponse.hpp @@ -58,26 +58,33 @@ public: // Implement pure virtual base class methods void clear(void); GCTAResponse* clone(void) const; - void load(const std::string& rspname); bool hasedisp(void) const { return false; } bool hastdisp(void) const { return false; } - double irf(const GInstDir& obsDir, const GEnergy& obsEng, const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; double irf(const GEvent& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; - double nirf(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; + double npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const; std::string print(void) const; // Other Methods + void caldb(const std::string& caldb); + std::string caldb(void) const { return m_caldb; } + void load(const std::string& rspname); + + // Other response methods + double irf(const GInstDir& obsDir, const GEnergy& obsEng, const GTime& obsTime, + const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; double irf(const GCTAEventAtom& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; double irf(const GCTAEventBin& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; + double npred(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; double live(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, const GPointing& pnt) const; double aeff(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, @@ -97,8 +104,6 @@ public: const GPointing& pnt, const GEbounds& ebds) const; double ntdisp(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, const GPointing& pnt, const GGti& gti) const; - - // Other Methods double psf(const double& theta, const double& sigma) const; double psf_sigma(const GEnergy& srcEng) const; double npsf(const double& psf, const double& radroi, const double& sigma) const; @@ -135,11 +140,13 @@ private: }; // Private data members - GNodeArray m_nodes; //!< log(E) nodes for interpolation - std::vector m_logE; //!< log(E) = log10(E/TeV) - bin centre - std::vector m_aeff; //!< Effective area in square metres after all cuts - std::vector m_r68; //!< 68% containment radius of PSF post cuts in degrees - std::vector m_r80; //!< 80% containment radius of PSF post cuts in degrees + std::string m_caldb; //!< Name of or path to the calibration database + std::string m_rspname; //!< Name of the instrument response + GNodeArray m_nodes; //!< log(E) nodes for interpolation + std::vector m_logE; //!< log(E) = log10(E/TeV) - bin centre + std::vector m_aeff; //!< Effective area in square metres after all cuts + std::vector m_r68; //!< 68% containment radius of PSF post cuts in degrees + std::vector m_r80; //!< 80% containment radius of PSF post cuts in degrees }; diff --git a/inst/cta/pyext/GCTAResponse.i b/inst/cta/pyext/GCTAResponse.i index 68fb19f..9379538 100644 --- a/inst/cta/pyext/GCTAResponse.i +++ b/inst/cta/pyext/GCTAResponse.i @@ -37,25 +37,32 @@ public: // Implement pure virtual base class methods void clear(void); GCTAResponse* clone(void) const; - void load(const std::string& rspname); bool hasedisp(void) const { return false; } bool hastdisp(void) const { return false; } - double irf(const GInstDir& obsDir, const GEnergy& obsEng, const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; double irf(const GEvent& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; - double nirf(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; + double npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const; // Other Methods + void caldb(const std::string& caldb); + std::string caldb(void) const { return m_caldb; } + void load(const std::string& rspname); + + // Other response methods + double irf(const GInstDir& obsDir, const GEnergy& obsEng, const GTime& obsTime, + const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; double irf(const GCTAEventAtom& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; double irf(const GCTAEventBin& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; + double npred(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; double live(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, const GPointing& pnt) const; double aeff(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, @@ -75,8 +82,6 @@ public: const GPointing& pnt, const GEbounds& ebds) const; double ntdisp(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, const GPointing& pnt, const GGti& gti) const; - - // Other Methods double psf(const double& theta, const double& sigma) const; double psf_sigma(const GEnergy& srcEng) const; double npsf(const double& psf, const double& radroi, const double& sigma) const; diff --git a/inst/cta/src/GCTAResponse.cpp b/inst/cta/src/GCTAResponse.cpp index f513f47..aa2773b 100644 --- a/inst/cta/src/GCTAResponse.cpp +++ b/inst/cta/src/GCTAResponse.cpp @@ -33,6 +33,13 @@ #include "GIntegrand.hpp" /* __ Method name definitions ____________________________________________ */ +#define G_CALDB "GCTAResponse::caldb(std::string&)" +#define G_IRF_ATOM "GCTAResponse::irf(GCTAEventAtom&, GModel&, GEnergy&,"\ + " GTime&, GObservation&)" +#define G_IRF_BIN "GCTAResponse::irf(GCTAEventBin&, GModel&, GEnergy&,"\ + " GTime&, GObservation&)" +#define G_NPRED "GCTAResponse::npred(GModel&, GEnergy&, GTime&,"\ + " GObservation&)" #define G_READ "GCTAResponse::read_performance_table(std::string&)" /* __ Macros _____________________________________________________________ */ @@ -51,7 +58,7 @@ ==========================================================================*/ /***********************************************************************//** - * @brief Constructor + * @brief Void constructor ***************************************************************************/ GCTAResponse::GCTAResponse(void) : GResponse() { @@ -66,7 +73,7 @@ GCTAResponse::GCTAResponse(void) : GResponse() /***********************************************************************//** * @brief Copy constructor * - * @param[in] rsp Response to be copied + * @param[in] rsp Response. ***************************************************************************/ GCTAResponse::GCTAResponse(const GCTAResponse& rsp) : GResponse(rsp) { @@ -103,7 +110,7 @@ GCTAResponse::~GCTAResponse(void) /***********************************************************************//** * @brief Assignment operator * - * @param[in] rsp Response to be assigned + * @param[in] rsp Response. ***************************************************************************/ GCTAResponse& GCTAResponse::operator= (const GCTAResponse& rsp) { @@ -136,36 +143,93 @@ GCTAResponse& GCTAResponse::operator= (const GCTAResponse& rsp) ==========================================================================*/ /***********************************************************************//** - * @brief Return value of point source instrument response function. + * @brief Clear instance +***************************************************************************/ +void GCTAResponse::clear(void) +{ + // Free class members (base and derived classes, derived class first) + free_members(); + this->GResponse::free_members(); + + // Initialise members + this->GResponse::init_members(); + init_members(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Clone instance +***************************************************************************/ +GCTAResponse* GCTAResponse::clone(void) const +{ + return new GCTAResponse(*this); +} + + +/***********************************************************************//** + * @brief Set path to the calibration database * - * @param[in] obsDir Observed photon direction. - * @param[in] obsEng Observed energy of photon. - * @param[in] obsTime Observed photon arrival time. - * @param[in] srcDir True photon direction. - * @param[in] srcEng True energy of photon. - * @param[in] srcTime True photon arrival time. - * @param[in] obs Observations. + * @param[in] caldb Path to calibration database + * + * @exception GException::caldb_not_found + * Calibration database repository not found. * - * This method implements the default and complete instrument response - * function (IRF). + * This default method simply checks if the calibration database directory + * exists. If the directory exists, the path will be stored. No checking is + * implemented that checks for the consistency of the calibration database. + * + * @todo Implement a GCalDB class that handles any calibration database + * issues. GCalDB may be an abstract class for which instrument + * specific methods are implement to handle any instrument specific + * IRF database issues. ***************************************************************************/ -double GCTAResponse::irf(const GInstDir& obsDir, const GEnergy& obsEng, - const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, const GObservation& obs) const +void GCTAResponse::caldb(const std::string& caldb) { - // Get pointing - const GPointing *pnt = obs.pointing(srcTime); + // Check if calibration database directory is accessible + if (access(caldb.c_str(), R_OK) != 0) + throw GException::caldb_not_found(G_CALDB, caldb); + + // Store the path to the calibration database + m_caldb = caldb; - // Get point source IRF components - double rsp = live(srcDir, srcEng, srcTime, *pnt); - rsp *= aeff(srcDir, srcEng, srcTime, *pnt); - rsp *= psf(obsDir, srcDir, srcEng, srcTime, *pnt); - rsp *= edisp(obsEng, srcDir, srcEng, srcTime, *pnt); - rsp *= tdisp(obsTime, srcDir, srcEng, srcTime, *pnt); + // Return + return; +} - // Return IRF value - return rsp; + +/***********************************************************************//** + * @brief Load CTA response. + * + * @param[in] irfname Name of CTA response (without any file extension). + * + * The actually dummy version of the CTA response loads a CTA performance + * table given in ASCII format into memory. + ***************************************************************************/ +void GCTAResponse::load(const std::string& irfname) +{ + // Save calibration database name + std::string caldb = m_caldb; + + // Clear instance + clear(); + + // Restore calibration database name + m_caldb = caldb; + + // Build filename + std::string filename = m_caldb + "/" + irfname + ".dat"; + + // Read performance table + read_performance_table(filename); + + // Store response name + m_rspname = irfname; + + // Return + return; } @@ -204,6 +268,9 @@ double GCTAResponse::irf(const GEvent& event, const GModel& model, * @param[in] srcEng True energy of photon. * @param[in] srcTime True photon arrival time. * @param[in] obs Observation. + * + * @exception GException::feature_not_implemented + * Diffuse IRF not yet implemented. ***************************************************************************/ double GCTAResponse::irf(const GCTAEventAtom& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, @@ -224,6 +291,15 @@ double GCTAResponse::irf(const GCTAEventAtom& event, const GModel& model, } // endif: model was point source + // ... otherwise return diffuse IRF + else { + + // Feature not yet implemented + throw GException::feature_not_implemented(G_IRF_ATOM, + "Diffuse IRF not yet implemented."); + + } // endelse: model was not a point source + // Return IRF value return rsp; } @@ -237,6 +313,9 @@ double GCTAResponse::irf(const GCTAEventAtom& event, const GModel& model, * @param[in] srcEng True energy of photon. * @param[in] srcTime True photon arrival time. * @param[in] obs Observation. + * + * @exception GException::feature_not_implemented + * Diffuse IRF not yet implemented. ***************************************************************************/ double GCTAResponse::irf(const GCTAEventBin& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, @@ -257,28 +336,106 @@ double GCTAResponse::irf(const GCTAEventBin& event, const GModel& model, } // endif: model was point source + // ... otherwise return diffuse IRF + else { + + // Feature not yet implemented + throw GException::feature_not_implemented(G_IRF_BIN, + "Diffuse IRF not yet implemented."); + + } // endelse: model was not a point source + // Return IRF value return rsp; } /***********************************************************************//** - * @brief Return integral of instrument response function. + * @brief Return value of point source instrument response function. * + * @param[in] obsDir Observed photon direction. + * @param[in] obsEng Observed energy of photon. + * @param[in] obsTime Observed photon arrival time. * @param[in] srcDir True photon direction. * @param[in] srcEng True energy of photon. * @param[in] srcTime True photon arrival time. + * @param[in] obs Observations. + ***************************************************************************/ +double GCTAResponse::irf(const GInstDir& obsDir, const GEnergy& obsEng, + const GTime& obsTime, + const GSkyDir& srcDir, const GEnergy& srcEng, + const GTime& srcTime, const GObservation& obs) const +{ + // Get pointing + const GPointing *pnt = obs.pointing(srcTime); + + // Get point source IRF components + double rsp = live(srcDir, srcEng, srcTime, *pnt); + rsp *= aeff(srcDir, srcEng, srcTime, *pnt); + rsp *= psf(obsDir, srcDir, srcEng, srcTime, *pnt); + rsp *= edisp(obsEng, srcDir, srcEng, srcTime, *pnt); + rsp *= tdisp(obsTime, srcDir, srcEng, srcTime, *pnt); + + // Return IRF value + return rsp; +} + + +/***********************************************************************//** + * @brief Return integral of instrument response function. + * + * @param[in] model Source model. + * @param[in] srcEng True energy of photon. + * @param[in] srcTime True photon arrival time. * @param[in] obs Observation. * @param[in] roi Region of interest of data selection. * @param[in] ebds Energy boundaries of data selection. * @param[in] gti Good Time Intervals of data selection. + ***************************************************************************/ +double GCTAResponse::npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, const GObservation& obs) const +{ + // Initialise response value + double rsp = 0.0; + + // If model is a point source then return the point source IRF + if (model.spatial()->isptsource()) { + + // Get point source location + GSkyDir srcDir = static_cast(model.spatial())->dir(); + + // Compute IRF + rsp = npred(srcDir, srcEng, srcTime, obs); + + } // endif: model was point source + + // ... otherwise return diffuse IRF + else { + + // Feature not yet implemented + throw GException::feature_not_implemented(G_NPRED, + "Diffuse IRF not yet implemented."); + + } // endelse: model was not a point source + + // Return response value + return rsp; +} + + +/***********************************************************************//** + * @brief Return integral of instrument response function. * - * This method implements the default and complete integral of the instrument - * response function (IRF). It may be overwritted by a specific method in the - * derived class that drops response terms that are not used. + * @param[in] srcDir True photon direction. + * @param[in] srcEng True energy of photon. + * @param[in] srcTime True photon arrival time. + * @param[in] obs Observation. + * @param[in] roi Region of interest of data selection. + * @param[in] ebds Energy boundaries of data selection. + * @param[in] gti Good Time Intervals of data selection. ***************************************************************************/ -double GCTAResponse::nirf(const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, const GObservation& obs) const +double GCTAResponse::npred(const GSkyDir& srcDir, const GEnergy& srcEng, + const GTime& srcTime, const GObservation& obs) const { // Get pointers const GPointing *pnt = obs.pointing(srcTime); @@ -509,61 +666,6 @@ double GCTAResponse::ntdisp(const GSkyDir& srcDir, const GEnergy& srcEng, /***********************************************************************//** - * @brief Load CTA response. - * - * @param[in] irfname Name of CTA response (without any file extension). - * - * The actually dummy version of the CTA response loads a CTA performance - * table given in ASCII format into memory. - ***************************************************************************/ -void GCTAResponse::load(const std::string& irfname) -{ - // Initialise response members - free_members(); - init_members(); - - // Build filename - std::string filename = m_caldb + "/" + irfname + ".dat"; - - // Read performance table - read_performance_table(filename); - - // Store response name - m_rspname = irfname; - - // Return - return; -} - - -/***********************************************************************//** - * @brief Clear instance -***************************************************************************/ -void GCTAResponse::clear(void) -{ - // Free class members (base and derived classes, derived class first) - free_members(); - this->GResponse::free_members(); - - // Initialise members - this->GResponse::init_members(); - init_members(); - - // Return - return; -} - - -/***********************************************************************//** - * @brief Clone instance -***************************************************************************/ -GCTAResponse* GCTAResponse::clone(void) const -{ - return new GCTAResponse(*this); -} - - -/***********************************************************************//** * @brief Print CTA response information ***************************************************************************/ std::string GCTAResponse::print(void) const @@ -786,6 +888,8 @@ double GCTAResponse::npsf_kern_azsym(const double& rad, void GCTAResponse::init_members(void) { // Initialise members + m_caldb.clear(); + m_rspname.clear(); m_logE.clear(); m_aeff.clear(); m_r68.clear(); @@ -804,11 +908,13 @@ void GCTAResponse::init_members(void) void GCTAResponse::copy_members(const GCTAResponse& rsp) { // Copy attributes - m_nodes = rsp.m_nodes; - m_logE = rsp.m_logE; - m_aeff = rsp.m_aeff; - m_r68 = rsp.m_r68; - m_r80 = rsp.m_r80; + m_caldb = rsp.m_caldb; + m_rspname = rsp.m_rspname; + m_nodes = rsp.m_nodes; + m_logE = rsp.m_logE; + m_aeff = rsp.m_aeff; + m_r68 = rsp.m_r68; + m_r80 = rsp.m_r80; // Return return; diff --git a/inst/lat/include/GLATResponse.hpp b/inst/lat/include/GLATResponse.hpp index 9587021..75f3ffe 100644 --- a/inst/lat/include/GLATResponse.hpp +++ b/inst/lat/include/GLATResponse.hpp @@ -1,5 +1,5 @@ /*************************************************************************** - * GLATResponse.hpp - GLAST LAT Response class * + * GLATResponse.hpp - Fermi LAT Response class * * ----------------------------------------------------------------------- * * copyright (C) 2008-2010 by Jurgen Knodlseder * * ----------------------------------------------------------------------- * @@ -12,7 +12,7 @@ ***************************************************************************/ /** * @file GLATResponse.hpp - * @brief GLATResponse class interface definition. + * @brief Fermi LAT Response class definition * @author J. Knodlseder */ @@ -21,9 +21,9 @@ /* __ Includes ___________________________________________________________ */ #include +#include #include "GLATEventAtom.hpp" #include "GLATEventBin.hpp" -#include "GLATResponseTable.hpp" #include "GLATAeff.hpp" #include "GLATPsf.hpp" #include "GLATEdisp.hpp" @@ -32,17 +32,8 @@ #include "GModel.hpp" #include "GObservation.hpp" #include "GResponse.hpp" -#include "GPointing.hpp" -#include "GInstDir.hpp" -#include "GRoi.hpp" -#include "GGti.hpp" -#include "GEbounds.hpp" -#include "GSkyDir.hpp" #include "GEnergy.hpp" #include "GTime.hpp" -#include "GVector.hpp" -#include "GFits.hpp" -#include "GFitsTable.hpp" /***********************************************************************//** @@ -64,34 +55,33 @@ public: // Implement pure virtual base class methods void clear(void); GLATResponse* clone(void) const; - void load(const std::string& rspname); bool hasedisp(void) const { return false; } bool hastdisp(void) const { return false; } - double irf(const GInstDir& obsDir, const GEnergy& obsEng, - const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const; double irf(const GEvent& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; - double nirf(const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const; + double npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const; std::string print(void) const; // Other Methods - int size(void) const { return m_aeff.size(); } - GLATAeff* aeff(const int& index) const; - GLATPsf* psf(const int& index) const; - GLATEdisp* edisp(const int& index) const; - void save(const std::string& rspname) const; - double irf(const GLATEventAtom& event, const GModel& model, - const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; - double irf(const GLATEventBin& event, const GModel& model, - const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; + void caldb(const std::string& caldb); + std::string caldb(void) const { return m_caldb; } + void load(const std::string& rspname); + int size(void) const { return m_aeff.size(); } + GLATAeff* aeff(const int& index) const; + GLATPsf* psf(const int& index) const; + GLATEdisp* edisp(const int& index) const; + void save(const std::string& rspname) const; + + // Reponse methods + double irf(const GLATEventAtom& event, const GModel& model, + const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; + double irf(const GLATEventBin& event, const GModel& model, + const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; private: // Private methods @@ -100,12 +90,14 @@ private: void free_members(void); // Private members - bool m_hasfront; //!< Front IRF loaded - bool m_hasback; //!< Back IRF loaded - std::vector m_aeff; //!< Effective areas - std::vector m_psf; //!< Point spread functions - std::vector m_edisp; //!< Energy dispersions - std::vector m_ptsrc; //!< Mean PSFs for point sources + std::string m_caldb; //!< Name of or path to the calibration database + std::string m_rspname; //!< Name of the instrument response + bool m_hasfront; //!< Front IRF loaded? + bool m_hasback; //!< Back IRF loaded? + std::vector m_aeff; //!< Effective areas + std::vector m_psf; //!< Point spread functions + std::vector m_edisp; //!< Energy dispersions + std::vector m_ptsrc; //!< Mean PSFs for point sources }; #endif /* GLATRESPONSE_HPP */ diff --git a/inst/lat/pyext/GLATResponse.i b/inst/lat/pyext/GLATResponse.i index 3e0ab50..59bdcfa 100644 --- a/inst/lat/pyext/GLATResponse.i +++ b/inst/lat/pyext/GLATResponse.i @@ -36,33 +36,32 @@ public: // Implement pure virtual base class methods void clear(void); GLATResponse* clone(void) const; - void load(const std::string& rspname); bool hasedisp(void) const { return false; } bool hastdisp(void) const { return false; } - double irf(const GInstDir& obsDir, const GEnergy& obsEng, - const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const; double irf(const GEvent& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const; - double nirf(const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const; + double npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const; // Other Methods - int size(void) const { return m_aeff.size(); } - GLATAeff* aeff(const int& index) const; - GLATPsf* psf(const int& index) const; - GLATEdisp* edisp(const int& index) const; - void save(const std::string& rspname) const; - double irf(const GLATEventAtom& event, const GModel& model, - const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; - double irf(const GLATEventBin& event, const GModel& model, - const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const; + void caldb(const std::string& caldb); + std::string caldb(void) const { return m_caldb; } + void load(const std::string& rspname); + int size(void) const { return m_aeff.size(); } + GLATAeff* aeff(const int& index) const; + GLATPsf* psf(const int& index) const; + GLATEdisp* edisp(const int& index) const; + void save(const std::string& rspname) const; + + // Reponse methods + double irf(const GLATEventAtom& event, const GModel& model, + const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; + double irf(const GLATEventBin& event, const GModel& model, + const GEnergy& srcEng, const GTime& srcTime, + const GObservation& obs) const; }; diff --git a/inst/lat/src/GLATResponse.cpp b/inst/lat/src/GLATResponse.cpp index 3790807..66e3183 100644 --- a/inst/lat/src/GLATResponse.cpp +++ b/inst/lat/src/GLATResponse.cpp @@ -12,7 +12,7 @@ ***************************************************************************/ /** * @file GLATResponse.cpp - * @brief GLATResponse class implementation. + * @brief Fermi LAT Response class implementation. * @author J. Knodlseder */ @@ -27,19 +27,23 @@ #include "GLATEventBin.hpp" #include "GLATEventCube.hpp" #include "GLATException.hpp" +#include "GException.hpp" #include "GModelSpatialPtsrc.hpp" -#include "GVector.hpp" #include "GFits.hpp" -#include "GFitsTable.hpp" #include "GTools.hpp" /* __ Method name definitions ____________________________________________ */ -#define G_LOAD "load(std::string&)" +#define G_CALDB "GLATResponse::caldb(std::string&)" +#define G_LOAD "GLATResponse::load(std::string&)" #define G_AEFF "GLATResponse::aeff(int&)" #define G_PSF "GLATResponse::psf(int&)" #define G_EDISP "GLATResponse::edisp(int&)" +#define G_IRF_ATOM "GLATResponse::irf(GLATEventAtom&, GModel&, GEnergy&,"\ + "GTime&, GObservation&)" #define G_IRF_BIN "GLATResponse::irf(GLATEventBin&, GModel&, GEnergy&,"\ "GTime&, GObservation&)" +#define G_NPRED "GLATResponse::npred(GModel&, GEnergy&, GTime&,"\ + " GObservation&)" /* __ Macros _____________________________________________________________ */ @@ -170,28 +174,33 @@ GLATResponse* GLATResponse::clone(void) const /***********************************************************************//** - * @brief Return value of point source instrument response function + * @brief Set path to the calibration database * - * @param[in] obsDir Observed photon direction. - * @param[in] obsEng Observed energy of photon. - * @param[in] obsTime Observed photon arrival time. - * @param[in] srcDir True photon direction. - * @param[in] srcEng True energy of photon. - * @param[in] srcTime True photon arrival time. - * @param[in] obs Observation. + * @param[in] caldb Path to calibration database * - * @todo Not yet implemented. + * @exception GException::caldb_not_found + * Calibration database repository not found. + * + * This default method simply checks if the calibration database directory + * exists. If the directory exists, the path will be stored. No checking is + * implemented that checks for the consistency of the calibration database. + * + * @todo Implement a GCalDB class that handles any calibration database + * issues. GCalDB may be an abstract class for which instrument + * specific methods are implement to handle any instrument specific + * IRF database issues. ***************************************************************************/ -double GLATResponse::irf(const GInstDir& obsDir, const GEnergy& obsEng, - const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, const GObservation& obs) const +void GLATResponse::caldb(const std::string& caldb) { - // Initialise - double rsp = 0.0; + // Check if calibration database directory is accessible + if (access(caldb.c_str(), R_OK) != 0) + throw GException::caldb_not_found(G_CALDB, caldb); - // Return IRF value - return rsp; + // Store the path to the calibration database + m_caldb = caldb; + + // Return + return; } @@ -243,8 +252,8 @@ double GLATResponse::irf(const GLATEventAtom& event, const GModel& model, // Initialise IRF with "no response" double irf = 0.0; - // Determine index for diffuse response - //int inx = event.diffinx(model.name()); + // Dump warning that integration is not yet implemented + throw GException::feature_not_implemented(G_IRF_ATOM); // Return IRF value return irf; @@ -385,22 +394,25 @@ double GLATResponse::irf(const GLATEventBin& event, const GModel& model, /***********************************************************************//** * @brief Return integral of instrument response function. * - * @param[in] srcDir True photon direction. + * @param[in] model Source model. * @param[in] srcEng True energy of photon. * @param[in] srcTime True photon arrival time. * @param[in] obs Observation. * * @todo Not yet implemented. ***************************************************************************/ -double GLATResponse::nirf(const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const +double GLATResponse::npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const { // Initialise - double nirf = 0.0; + double npred = 0.0; + + // Notify that method is not yet implemented + throw GException::feature_not_implemented(G_NPRED); // Return integrated IRF value - return nirf; + return npred; } @@ -422,9 +434,14 @@ double GLATResponse::nirf(const GSkyDir& srcDir, const GEnergy& srcEng, ***************************************************************************/ void GLATResponse::load(const std::string& rspname) { - // Free members (cannot call clear() here since this will overwrite the - // calibration database path m_caldb) - free_members(); + // Save calibration database name + std::string caldb = m_caldb; + + // Clear instance + clear(); + + // Restore calibration database name + m_caldb = caldb; // Determine response types to be loaded std::vector array = split(rspname, "::"); @@ -615,6 +632,8 @@ std::string GLATResponse::print(void) const void GLATResponse::init_members(void) { // Initialise members + m_caldb.clear(); + m_rspname.clear(); m_hasfront = false; m_hasback = false; m_aeff.clear(); @@ -640,6 +659,8 @@ void GLATResponse::init_members(void) void GLATResponse::copy_members(const GLATResponse& rsp) { // Copy members + m_caldb = rsp.m_caldb; + m_rspname = rsp.m_rspname; m_hasfront = rsp.m_hasfront; m_hasback = rsp.m_hasback; m_aeff = rsp.m_aeff; diff --git a/inst/mwl/include/GMWLResponse.hpp b/inst/mwl/include/GMWLResponse.hpp index 005cd43..2168395 100644 --- a/inst/mwl/include/GMWLResponse.hpp +++ b/inst/mwl/include/GMWLResponse.hpp @@ -44,19 +44,16 @@ public: GMWLResponse& operator= (const GMWLResponse & rsp); // Reponse function computation methods - double irf(const GInstDir& obsDir, const GEnergy& obsEng, const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const { return 1.0; } double irf(const GEvent& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const { return 1.0; } - double nirf(const GSkyDir& srcDir, const GEnergy& srcEng, const GTime& srcTime, - const GObservation& obs) const { return 1.0; } + double npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const { return 1.0; } // Pure virtual base class methods void clear(void); GMWLResponse* clone(void) const; - void load(const std::string& irfname) { return; } bool hasedisp(void) const { return false; } bool hastdisp(void) const { return false; } std::string print(void) const; diff --git a/inst/mwl/src/GMWLResponse.cpp b/inst/mwl/src/GMWLResponse.cpp index ef1bc05..0b94cf8 100644 --- a/inst/mwl/src/GMWLResponse.cpp +++ b/inst/mwl/src/GMWLResponse.cpp @@ -160,8 +160,8 @@ std::string GMWLResponse::print(void) const // Append header result.append("=== GMWLResponse ===\n"); - result.append(parformat("Calibration database")+m_caldb+"\n"); - result.append(parformat("Response functions")+m_rspname); + //result.append(parformat("Calibration database")+m_caldb+"\n"); + //result.append(parformat("Response functions")+m_rspname); // Return result return result; diff --git a/pyext/GResponse.i b/pyext/GResponse.i index 2b40728..2875d3b 100644 --- a/pyext/GResponse.i +++ b/pyext/GResponse.i @@ -1,5 +1,5 @@ /*************************************************************************** - * GResponse.i - Response abstract base class SWIG interface * + * GResponse.i - Abstract response base class * * ----------------------------------------------------------------------- * * copyright (C) 2008-2010 by Jurgen Knodlseder * * ----------------------------------------------------------------------- * @@ -12,7 +12,7 @@ ***************************************************************************/ /** * @file GResponse.i - * @brief GResponse class SWIG file. + * @brief Abstract response base class python bindings * @author J. Knodlseder */ %{ @@ -24,7 +24,7 @@ /***********************************************************************//** * @class GResponse * - * @brief Abstract SWIG interface for the instrument response function. + * @brief Abstract response base class python bindings ***************************************************************************/ class GResponse { public: @@ -36,24 +36,14 @@ public: // Pure virtual methods virtual void clear(void) = 0; virtual GResponse* clone(void) const = 0; - virtual void load(const std::string& irfname) = 0; virtual bool hasedisp(void) const = 0; virtual bool hastdisp(void) const = 0; - virtual double irf(const GInstDir& obsDir, const GEnergy& obsEng, - const GTime& obsTime, - const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const = 0; virtual double irf(const GEvent& event, const GModel& model, const GEnergy& srcEng, const GTime& srcTime, const GObservation& obs) const = 0; - virtual double nirf(const GSkyDir& srcDir, const GEnergy& srcEng, - const GTime& srcTime, - const GObservation& obs) const = 0; - - // Other methods - virtual void caldb(const std::string& caldb); - virtual std::string caldb(void) const { return m_caldb; } + virtual double npred(const GModel& model, const GEnergy& srcEng, + const GTime& srcTime, + const GObservation& obs) const = 0; }; diff --git a/src/obs/GResponse.cpp b/src/obs/GResponse.cpp index 8ee5b7c..d368e60 100644 --- a/src/obs/GResponse.cpp +++ b/src/obs/GResponse.cpp @@ -21,7 +21,6 @@ #include "GException.hpp" /* __ Method name definitions ____________________________________________ */ -#define G_CALDB "GResponse::caldb(std::string&)" /* __ Macros _____________________________________________________________ */ @@ -53,7 +52,7 @@ GResponse::GResponse(void) /***********************************************************************//** * @brief Copy constructor * - * @param[in] rsp Response from which the instance should be built. + * @param[in] rsp Response. ***************************************************************************/ GResponse::GResponse(const GResponse& rsp) { @@ -90,7 +89,7 @@ GResponse::~GResponse(void) /***********************************************************************//** * @brief Assignment operator * - * @param[in] rsp Response which should be assigned. + * @param[in] rsp Response. ***************************************************************************/ GResponse& GResponse::operator= (const GResponse& rsp) { @@ -119,33 +118,6 @@ GResponse& GResponse::operator= (const GResponse& rsp) = = ==========================================================================*/ -/***********************************************************************//** - * @brief Set path to the calibration database - * - * @param[in] caldb Path to calibration database - * - * This default method simply checks if the calibration database directory - * exists. If the directory exists, the path will be stored. No checking is - * implemented that checks for the consistency of the calibration database. - * - * @todo Implement a GCalDB class that handles any calibration database - * issues. GCalDB may be an abstract class for which instrument - * specific methods are implement to handle any instrument specific - * IRF database issues. - ***************************************************************************/ -void GResponse::caldb(const std::string& caldb) -{ - // Check if calibration database directory is accessible - if (access(caldb.c_str(), R_OK) != 0) - throw GException::caldb_not_found(G_CALDB, caldb); - - // Store the path to the calibration database - m_caldb = caldb; - - // Return - return; -} - /*========================================================================== = = @@ -158,10 +130,6 @@ void GResponse::caldb(const std::string& caldb) ***************************************************************************/ void GResponse::init_members(void) { - // Initialise members - m_caldb.clear(); - m_rspname.clear(); - // Return return; } @@ -174,10 +142,6 @@ void GResponse::init_members(void) ***************************************************************************/ void GResponse::copy_members(const GResponse& rsp) { - // Copy attributes - m_caldb = rsp.m_caldb; - m_rspname = rsp.m_rspname; - // Return return; } -- libgit2 0.21.4