GammaLib  2.1.0.dev
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GResponse Class Referenceabstract

Abstract instrument response base class. More...

#include <GResponse.hpp>

Inheritance diagram for GResponse:
GBase GCOMResponse GCTAResponse GLATResponse GMWLResponse GSPIResponse GCTAResponseCube GCTAResponseIrf

Classes

class  edisp_kerns
 
class  irf_elliptical_kern_phi
 
class  irf_elliptical_kern_theta
 
class  irf_radial_kern_phi
 
class  irf_radial_kern_theta
 

Public Member Functions

 GResponse (void)
 Void constructor. More...
 
 GResponse (const GResponse &rsp)
 Copy constructor. More...
 
virtual ~GResponse (void)
 Destructor. More...
 
virtual GResponseoperator= (const GResponse &rsp)
 Assignment operator. More...
 
virtual void clear (void)=0
 Clear object. More...
 
virtual GResponseclone (void) const =0
 Clones object. More...
 
virtual std::string classname (void) const =0
 Return class name. More...
 
virtual bool use_edisp (void) const =0
 
virtual bool use_tdisp (void) const =0
 
virtual double irf (const GEvent &event, const GPhoton &photon, const GObservation &obs) const =0
 
virtual double nroi (const GModelSky &model, const GEnergy &obsEng, const GTime &obsTime, const GObservation &obs) const =0
 
virtual GEbounds ebounds (const GEnergy &obsEng) const =0
 
virtual std::string print (const GChatter &chatter=NORMAL) const =0
 Print content of object. More...
 
virtual double convolve (const GModelSky &model, const GEvent &event, const GObservation &obs, const bool &grad=true) const
 Convolve sky model with the instrument response. More...
 
virtual GVector convolve (const GModelSky &model, const GObservation &obs, GMatrixSparse *gradients=NULL) const
 Convolve sky model with the instrument response. More...
 
virtual double irf_spatial (const GEvent &event, const GSource &source, const GObservation &obs) const
 Return instrument response integrated over the spatial model. More...
 
virtual GVector irf_spatial (const GModelSky &model, const GObservation &obs, GMatrix *gradients=NULL) const
 Return instrument response vector integrated over the spatial model. More...
 
virtual void remove_response_cache (const std::string &name)
 Remove response cache for model. More...
 
- Public Member Functions inherited from GBase
virtual ~GBase (void)
 Destructor. More...
 

Protected Member Functions

void init_members (void)
 Initialise class members. More...
 
void copy_members (const GResponse &rsp)
 Copy class members. More...
 
void free_members (void)
 Delete class members. More...
 
double eval_prob (const GModelSky &model, const GEvent &event, const GEnergy &srcEng, const GTime &srcTime, const GObservation &obs, const bool &grad) const
 Convolve sky model with the instrument response. More...
 
GVector eval_probs (const GModelSky &model, const GObservation &obs, GMatrixSparse *gradients=NULL) const
 Convolve sky model with the instrument response. More...
 
int size_edisp_vector (const GModelSky &model, const GObservation &obs, const bool &grad) const
 Return size of vector for energy dispersion computation. More...
 
GEbounds ebounds_model (const GModelSky &model) const
 Return true energy intervals for sky model. More...
 
virtual double irf_ptsrc (const GEvent &event, const GSource &source, const GObservation &obs) const
 Return instrument response to point source. More...
 
virtual double irf_radial (const GEvent &event, const GSource &source, const GObservation &obs) const
 Return instrument response to radial source. More...
 
virtual double irf_elliptical (const GEvent &event, const GSource &source, const GObservation &obs) const
 Return instrument response to elliptical source. More...
 
virtual double irf_diffuse (const GEvent &event, const GSource &source, const GObservation &obs) const
 Return instrument response to diffuse source. More...
 
virtual double irf_composite (const GEvent &event, const GSource &source, const GObservation &obs) const
 Return instrument response to composite source. More...
 
virtual GVector irf_ptsrc (const GModelSky &model, const GObservation &obs, GMatrix *gradients=NULL) const
 Return instrument response to point source sky model. More...
 
virtual GVector irf_radial (const GModelSky &model, const GObservation &obs, GMatrix *gradients=NULL) const
 Return instrument response to radial source sky model. More...
 
virtual GVector irf_elliptical (const GModelSky &model, const GObservation &obs, GMatrix *gradients=NULL) const
 Return instrument response to ellipitical source sky model. More...
 
virtual GVector irf_diffuse (const GModelSky &model, const GObservation &obs, GMatrix *gradients=NULL) const
 Return instrument response to diffuse source sky model. More...
 
virtual GVector irf_composite (const GModelSky &model, const GObservation &obs, GMatrix *gradients=NULL) const
 Return instrument response to composite source sky model. More...
 

Protected Attributes

bool m_use_irf_cache
 Control usage of irf cache. More...
 
bool m_use_nroi_cache
 Control usage of nroi cache. More...
 
int m_irf_radial_iter_theta
 Radial model integration theta iterations. More...
 
int m_irf_radial_iter_phi
 Radial model integration phi iterations. More...
 
int m_irf_elliptical_iter_theta
 Elliptical model integration theta iterations. More...
 
int m_irf_elliptical_iter_phi
 Elliptical model integration phi iterations. More...
 
double m_irf_diffuse_resolution
 Angular resolution for diffuse model. More...
 
GResponseCache m_irf_cache
 
GResponseCache m_nroi_cache
 
GResponseVectorCache m_irf_vector_cache
 

Detailed Description

Abstract instrument response base class.

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 photon. An alternative method exists that returns the response for a specific source.

The nroi method returns the spatial integral of the instrument response function times the sky model over the region of interest. This method is only required for unbinned analysis.

The ebounds method returns the true energy boundaries for a specified measured event energy. This method is used for computing the energy dispersion.

Definition at line 77 of file GResponse.hpp.

Constructor & Destructor Documentation

GResponse::GResponse ( void  )

Void constructor.

Definition at line 86 of file GResponse.cpp.

References init_members().

GResponse::GResponse ( const GResponse rsp)

Copy constructor.

Parameters
[in]rspResponse.

Definition at line 101 of file GResponse.cpp.

References copy_members(), and init_members().

GResponse::~GResponse ( void  )
virtual

Destructor.

Definition at line 117 of file GResponse.cpp.

References free_members().

Member Function Documentation

virtual std::string GResponse::classname ( void  ) const
pure virtual

Return class name.

Returns
String containing the class name.

Returns the class name for non-abstract classes in a human readable way.

Implements GBase.

Implemented in GCTAResponseCube, GCTAResponseIrf, GSPIResponse, GCOMResponse, GLATResponse, GCTAResponse, and GMWLResponse.

virtual void GResponse::clear ( void  )
pure virtual

Clear object.

Sets the object to a clean initial state. After calling the method the object will be in the same state as it were if an empty instance of the object would have been created.

Implements GBase.

Implemented in GCTAResponseCube, GCTAResponseIrf, GSPIResponse, GCOMResponse, GLATResponse, GCTAResponse, and GMWLResponse.

virtual GResponse* GResponse::clone ( void  ) const
pure virtual

Clones object.

Returns
Pointer to deep copy of object.

Creates a deep copy of the object and returns a pointer to the object.

Implements GBase.

Implemented in GCTAResponseCube, GCTAResponseIrf, GSPIResponse, GCOMResponse, GLATResponse, GCTAResponse, and GMWLResponse.

double GResponse::convolve ( const GModelSky model,
const GEvent event,
const GObservation obs,
const bool &  grad = true 
) const
virtual

Convolve sky model with the instrument response.

Parameters
[in]modelSky model.
[in]eventEvent.
[in]obsObservation.
[in]gradShould model gradients be computed?
Returns
Event probability.

Computes the event probability

\[ P(p',E',t') = \int \int \int S(p,E,t) \times R(p',E',t'|p,E,t) \, dp \, dE \, dt \]

without taking into account any time dispersion. Energy dispersion is correctly handled by this method. If time dispersion is indeed needed, an instrument specific method needs to be provided.

Definition at line 186 of file GResponse.cpp.

References ebounds(), ebounds_model(), GEbounds::emax(), GEbounds::emin(), GEvent::energy(), eval_prob(), GOptimizerPar::factor_gradient(), GIntegrals::fixed_iter(), GOptimizerPar::has_grad(), GModel::has_scales(), GObservation::instrument(), GOptimizerPar::is_free(), gammalib::is_infinite(), gammalib::is_notanumber(), GEnergy::MeV(), GOptimizerPar::name(), GIntegrals::romberg(), GModel::scale(), GModel::scales(), GModelTemporal::size(), GEbounds::size(), GModelSpectral::size(), size_edisp_vector(), GModelSky::spatial(), GModelSky::spectral(), GModelSky::temporal(), and use_edisp().

GVector GResponse::convolve ( const GModelSky model,
const GObservation obs,
GMatrixSparse gradients = NULL 
) const
virtual

Convolve sky model with the instrument response.

Parameters
[in]modelSky model.
[in]obsObservation.
[out]gradientsPointer to matrix of gradients.
Returns
Vector of event probabilities.

Computes the event probability

\[ P(p',E',t') = \int \int \int S(p,E,t) \times R(p',E',t'|p,E,t) \, dp \, dE \, dt \]

without taking into account any time dispersion. Energy dispersion is correctly handled by this method. If time dispersion is indeed needed, an instrument specific method needs to be provided.

Definition at line 346 of file GResponse.cpp.

References GMatrixSparse::column(), GMatrixBase::columns(), ebounds(), ebounds_model(), GEbounds::emax(), GEbounds::emin(), eval_probs(), GObservation::events(), GIntegrals::fixed_iter(), G_CONVOLVE, GOptimizerPar::has_grad(), GModel::has_scales(), GObservation::instrument(), GOptimizerPar::is_free(), gammalib::is_infinite(), gammalib::is_notanumber(), GEnergy::MeV(), GOptimizerPar::name(), GIntegrals::romberg(), GMatrixBase::rows(), GModel::scale(), GModel::scales(), GModelTemporal::size(), GEvents::size(), GEbounds::size(), GModelSpectral::size(), GModelSpatial::size(), GModel::size(), size_edisp_vector(), GModelSky::spatial(), GModelSky::spectral(), gammalib::str(), GModelSky::temporal(), and use_edisp().

void GResponse::copy_members ( const GResponse rsp)
protected
virtual GEbounds GResponse::ebounds ( const GEnergy obsEng) const
pure virtual
GEbounds GResponse::ebounds_model ( const GModelSky model) const
protected

Return true energy intervals for sky model.

Parameters
[in]modelSky model.
Returns
True energy intervals.

Returns the true energy intervals for a sky model. For all spectral models other than the GModelSpectralGauss model the method will return an empty energy boundaries structure. For the GModelSpectralGauss model the method will return the interval [mean - 5 * sigma, mean + 5 * sigma].

Definition at line 1965 of file GResponse.cpp.

References GEbounds::append(), ebounds(), GModelSpectralGauss::mean(), GEnergy::MeV(), GModelSpectralGauss::sigma(), and GModelSky::spectral().

Referenced by convolve().

double GResponse::eval_prob ( const GModelSky model,
const GEvent event,
const GEnergy srcEng,
const GTime srcTime,
const GObservation obs,
const bool &  grad 
) const
protected

Convolve sky model with the instrument response.

Parameters
[in]modelSky model.
[in]eventEvent.
[in]srcEngSource energy.
[in]srcTimeSource time.
[in]obsObservation.
[in]gradShould model gradients be computed? (default: true)
Returns
Event probability.

Computes the event probability

\[ P(p',E',t'|E,t) = \int S(p,E,t) \times R(p',E',t'|p,E,t) \, dp \]

for a given true energy \(E\) and time \(t\).

Definition at line 1536 of file GResponse.cpp.

References GModelTemporal::eval(), GModelSpectral::eval(), GOptimizerPar::factor_gradient(), GModel::has_scales(), GObservation::instrument(), irf(), irf_spatial(), GOptimizerPar::is_free(), gammalib::is_infinite(), gammalib::is_notanumber(), GModel::name(), GOptimizerPar::name(), GOptimizerPar::scale(), GModel::scale(), GModel::scales(), GModelTemporal::size(), GModelSpectral::size(), GModelSky::spatial(), GModelSky::spectral(), and GModelSky::temporal().

Referenced by convolve().

GVector GResponse::eval_probs ( const GModelSky model,
const GObservation obs,
GMatrixSparse gradients = NULL 
) const
protected

Convolve sky model with the instrument response.

Parameters
[in]modelSky model.
[in]obsObservation.
[out]gradientsPointer to matrix of gradients.
Returns
Event probabilities.

Computes the event probability

\[ P(p',E',t'|E,t) = \int S(p,E,t) \times R(p',E',t'|p,E,t) \, dp \]

for a all events.

Definition at line 1674 of file GResponse.cpp.

References GMatrixSparse::column(), GModelTemporal::eval(), GModelSpectral::eval(), GObservation::events(), GModel::has_scales(), GObservation::instrument(), irf_spatial(), GOptimizerPar::is_free(), gammalib::is_infinite(), gammalib::is_notanumber(), GOptimizerPar::name(), GOptimizerPar::scale(), GModel::scale(), GModel::scales(), GEvents::size(), GModelTemporal::size(), GModelSpectral::size(), GModelSpatial::size(), GModel::size(), GModelSky::spatial(), GModelSky::spectral(), and GModelSky::temporal().

Referenced by convolve().

void GResponse::free_members ( void  )
protected
double GResponse::irf_composite ( const GEvent event,
const GSource source,
const GObservation obs 
) const
protectedvirtual

Return instrument response to composite source.

Parameters
[in]eventObserved event.
[in]sourceSource.
[in]obsObservation.
Returns
Instrument response to composite source.

Returns the instrument response to a specified composite source.

Definition at line 1255 of file GResponse.cpp.

References GModelSpatialComposite::component(), GModelSpatialComposite::components(), GSource::energy(), irf(), irf_spatial(), GSource::model(), GSource::name(), GModelSpatialComposite::scale(), sum(), GModelSpatialComposite::sum_of_scales(), and GSource::time().

Referenced by GCTAResponseCube::irf_spatial(), and irf_spatial().

GVector GResponse::irf_composite ( const GModelSky model,
const GObservation obs,
GMatrix gradients = NULL 
) const
protectedvirtual

Return instrument response to composite source sky model.

Parameters
[in]modelSky model.
[in]obsObservation.
[in]gradientsGradients matrix.
Returns
Instrument response to composite source sky model.

Returns the instrument response to a composite source sky model for all events.

Definition at line 1478 of file GResponse.cpp.

References GModelSpatialComposite::component(), GModelSpatialComposite::components(), GObservation::events(), irf_spatial(), GModelSpatialComposite::scale(), GEvents::size(), GModelSky::spatial(), sum(), and GModelSpatialComposite::sum_of_scales().

double GResponse::irf_diffuse ( const GEvent event,
const GSource source,
const GObservation obs 
) const
protectedvirtual

Return instrument response to diffuse source.

Parameters
[in]eventObserved event.
[in]sourceSource.
[in]obsObservation.
Returns
Instrument response to diffuse source.

Returns the instrument response to a diffuse source for a given event, source and observation. The method computes

\[ {\tt irf}(p', E', t') = \sum_p M_{\rm S}(p | E, t) \, \Delta(p) R(p', E', t' | p, E, t) \]

where

  • \(M_{\rm S}(p | E, t)\) is the diffuse model component,
  • \(\Delta(p)\) is the solid angle of the map pixel
  • \(R(p', E', t' | p, E, t)\) is the Instrument Response Function (IRF),
  • \(p'\) is the measured instrument direction,
  • \(E'\) is the measured or reconstructed energy,
  • \(t'\) is the measured arrival time,
  • \(p\) is the true photon arrival direction,
  • \(E\) is the true photon energy, and
  • \(t\) is the true trigger time.

The method simply adds up all sky map pixels multiplied with the IRF for all diffuse model pixels. For models that do not contain sky map pixels (such as the GModelSpatialDiffuseConst model) the method allocates an all-sky sky map with an angular resolution that is specified by the m_irf_diffuse_resolution member.

Reimplemented in GCTAResponseIrf, GCTAResponseCube, and GLATResponse.

Definition at line 1118 of file GResponse.cpp.

References GModelSpatialDiffuseCube::cube(), GSource::energy(), GModelSpatialDiffuse::eval(), GModelSpatialDiffuseCube::eval(), GSkyMap::inx2dir(), irf(), m_irf_diffuse_resolution, GModelSpatialDiffuseMap::map(), GSource::model(), GSkyMap::npix(), GSkyMap::solidangle(), and GSource::time().

Referenced by irf_diffuse(), and irf_spatial().

GVector GResponse::irf_diffuse ( const GModelSky model,
const GObservation obs,
GMatrix gradients = NULL 
) const
protectedvirtual

Return instrument response to diffuse source sky model.

Parameters
[in]modelSky model.
[in]obsObservation.
[in]gradientsGradients matrix.
Returns
Instrument response to diffuse source sky model.

Returns the instrument response to a diffuse source sky model for all events.

Reimplemented in GCOMResponse.

Definition at line 1434 of file GResponse.cpp.

References GObservation::events(), irf_diffuse(), GModel::name(), GEvents::size(), and GModelSky::spatial().

double GResponse::irf_elliptical ( const GEvent event,
const GSource source,
const GObservation obs 
) const
protectedvirtual

Return instrument response to elliptical source.

Parameters
[in]eventObserved event.
[in]sourceSource.
[in]obsObservation.
Returns
Instrument response to elliptical source.

Returns the instrument response to an elliptical source for a given event, source and observation. The method computes

\[ {\tt irf}(p', E', t') = \int_0^{\theta_{\rm max}} \int_0^{2\pi} M_{\rm S}(\theta, \phi | E, t) \, R(p', E', t' | \theta, \phi, E, t) \, sin \, \theta \, d\,\phi \, d\,\theta \]

where

  • \(M_{\rm S}(\theta, \phi | E, t)\) is the elliptical model component,
  • \(R(p', E', t' | \theta, \phi, E, t)\) is the Instrument Response Function (IRF),
  • \(p'\) is the measured instrument direction,
  • \(E'\) is the measured or reconstructed energy,
  • \(t'\) is the measured arrival time,
  • \(\theta\) is the radial distance from the model centre,
  • \(\phi\) is the azimuthal angle around the model centre,
  • \(E\) is the true photon energy, and
  • \(t\) is the true trigger time.

The azimuth integration is done over \(2\pi\) using the irf_elliptical_kern_phi::eval() method. The radial integration is done out to a maximum angle \(\theta_{\rm max}\), given by GModelSpatialElliptical::theta_max(), by the irf_elliptical_kern_theta::eval() method.

Reimplemented in GCTAResponseIrf, and GCTAResponseCube.

Definition at line 1015 of file GResponse.cpp.

References GSkyDir::dec_deg(), GModelSpatialElliptical::dir(), GSource::energy(), GMatrix::eulery(), GMatrix::eulerz(), GIntegral::fixed_iter(), irf(), gammalib::is_infinite(), gammalib::is_notanumber(), m_irf_elliptical_iter_phi, m_irf_elliptical_iter_theta, GSource::model(), GSkyDir::ra_deg(), GIntegral::romberg(), GModelSpatialElliptical::theta_max(), and GSource::time().

Referenced by irf_elliptical(), and irf_spatial().

GVector GResponse::irf_elliptical ( const GModelSky model,
const GObservation obs,
GMatrix gradients = NULL 
) const
protectedvirtual

Return instrument response to ellipitical source sky model.

Parameters
[in]modelSky model.
[in]obsObservation.
[in]gradientsGradients matrix.
Returns
Instrument response to ellipitical source sky model.

Returns the instrument response to a ellipitical source sky model for all events.

Reimplemented in GCOMResponse.

Definition at line 1390 of file GResponse.cpp.

References GObservation::events(), irf_elliptical(), GModel::name(), GEvents::size(), and GModelSky::spatial().

double GResponse::irf_ptsrc ( const GEvent event,
const GSource source,
const GObservation obs 
) const
protectedvirtual

Return instrument response to point source.

Parameters
[in]eventObserved event.
[in]sourceSource.
[in]obsObservation.
Returns
Instrument response to point source.

Returns the instrument response to a point source.

Reimplemented in GCTAResponseIrf, GCTAResponseCube, and GLATResponse.

Definition at line 852 of file GResponse.cpp.

References GModelSpatialPointSource::dir(), GSource::energy(), irf(), GSource::model(), and GSource::time().

Referenced by irf_ptsrc(), and irf_spatial().

GVector GResponse::irf_ptsrc ( const GModelSky model,
const GObservation obs,
GMatrix gradients = NULL 
) const
protectedvirtual

Return instrument response to point source sky model.

Parameters
[in]modelSky model.
[in]obsObservation.
[in]gradientsGradients matrix.
Returns
Instrument response to point source sky model.

Returns the instrument response to a point source sky model for all events.

Reimplemented in GCOMResponse.

Definition at line 1302 of file GResponse.cpp.

References GObservation::events(), irf_ptsrc(), GModel::name(), GEvents::size(), and GModelSky::spatial().

double GResponse::irf_radial ( const GEvent event,
const GSource source,
const GObservation obs 
) const
protectedvirtual

Return instrument response to radial source.

Parameters
[in]eventObserved event.
[in]sourceSource.
[in]obsObservation.
Returns
Instrument response to radial source.

Returns the instrument response to a radial source for a given event, source and observation. The method computes

\[ {\tt irf}(p', E', t') = \int_0^{\theta_{\rm max}} M_{\rm S}(\theta | E, t) \, sin \, \theta \int_0^{2\pi} R(p', E', t' | \theta, \phi, E, t) \, d\,\phi \, d\,\theta \]

where

  • \(M_{\rm S}(\theta | E, t)\) is the radial model component,
  • \(R(p', E', t' | \theta, \phi, E, t)\) is the Instrument Response Function (IRF),
  • \(p'\) is the measured instrument direction,
  • \(E'\) is the measured or reconstructed energy,
  • \(t'\) is the measured arrival time,
  • \(\theta\) is the radial distance from the model centre,
  • \(\phi\) is the azimuthal angle around the model centre,
  • \(E\) is the true photon energy, and
  • \(t\) is the true trigger time.

The azimuth integration is done over \(2\pi\) using the irf_radial_kern_phi::eval() method. The radial integration is done out to a maximum angle \(\theta_{\rm max}\), given by GModelSpatialRadial::theta_max(), by the irf_radial_kern_theta::eval() method.

Reimplemented in GCTAResponseIrf, and GCTAResponseCube.

Definition at line 908 of file GResponse.cpp.

References GSkyDir::dec_deg(), GModelSpatialRadial::dir(), GSource::energy(), GMatrix::eulery(), GMatrix::eulerz(), GIntegral::fixed_iter(), irf(), gammalib::is_infinite(), gammalib::is_notanumber(), m_irf_radial_iter_phi, m_irf_radial_iter_theta, GSource::model(), GSkyDir::ra_deg(), GIntegral::romberg(), GModelSpatialRadial::theta_max(), and GSource::time().

Referenced by irf_radial(), and irf_spatial().

GVector GResponse::irf_radial ( const GModelSky model,
const GObservation obs,
GMatrix gradients = NULL 
) const
protectedvirtual

Return instrument response to radial source sky model.

Parameters
[in]modelSky model.
[in]obsObservation.
[in]gradientsGradients matrix.
Returns
Instrument response to radial source sky model.

Returns the instrument response to a radial source sky model for all events.

Reimplemented in GCTAResponseCube, and GCOMResponse.

Definition at line 1346 of file GResponse.cpp.

References GObservation::events(), irf_radial(), GModel::name(), GEvents::size(), and GModelSky::spatial().

double GResponse::irf_spatial ( const GEvent event,
const GSource source,
const GObservation obs 
) const
virtual

Return instrument response integrated over the spatial model.

Parameters
[in]eventEvent.
[in]sourceSource.
[in]obsObservation.
Returns
Instrument response to a spatial model.

Returns the instrument response for a given event, source and observation integrated over the spatial model component. The method computes

\[ {\tt irf}(p', E', t') = \int_p M_{\rm S}(p | E, t) \, R(p', E', t' | p, E, t) \, d\,p \]

where

  • \(M_{\rm S}(p | E, t)\) is the spatial model component,
  • \(R(p', E', t' | p, E, t)\) is the Instrument Response Function (IRF),
  • \(p'\) is the measured instrument direction,
  • \(E'\) is the measured or reconstructed energy,
  • \(t'\) is the measured arrival time,
  • \(p\) is the true photon arrival direction,
  • \(E\) is the true photon energy, and
  • \(t\) is the true trigger time.

The integration is done over all relevant true sky directions \(p\).

Depending on the type of the source model the method branches to the following methods to perform the actual computations

The method implements a caching mechanism for spatial models that have all parameters fixed. For those models the instrument response for a given event and observation is only computed once and then stored in an internal cache from which it is fetched back in case that the method is called again for the same event and observation.

Reimplemented in GCTAResponseCube, GLATResponse, and GMWLResponse.

Definition at line 596 of file GResponse.cpp.

References GModelSpatial::code(), GResponseCache::contains(), GSource::energy(), GMODEL_SPATIAL_COMPOSITE, GMODEL_SPATIAL_DIFFUSE, GMODEL_SPATIAL_ELLIPTICAL, GMODEL_SPATIAL_POINT_SOURCE, GMODEL_SPATIAL_RADIAL, GModelSpatial::has_free_pars(), GObservation::id(), irf(), irf_composite(), irf_diffuse(), irf_elliptical(), irf_ptsrc(), irf_radial(), m_irf_cache, m_use_irf_cache, GSource::model(), GSource::name(), and GResponseCache::set().

Referenced by GCTAOnOffObservation::compute_arf(), eval_prob(), eval_probs(), and irf_composite().

GVector GResponse::irf_spatial ( const GModelSky model,
const GObservation obs,
GMatrix gradients = NULL 
) const
virtual

Return instrument response vector integrated over the spatial model.

Parameters
[in]modelSky model.
[in]obsObservation.
[in]gradientsGradients matrix.
Returns
Instrument response vector to a spatial model.

Returns the instrument response to a sky model integrated over the spatial model component for all events in a given observation. The method computes

\[ {\tt irf}(p', E', t') = \int_p M_{\rm S}(p | E, t) \, R(p', E', t' | p, E, t) \, d\,p \]

where

  • \(M_{\rm S}(p | E, t)\) is the spatial model component,
  • \(R(p', E', t' | p, E, t)\) is the Instrument Response Function (IRF),
  • \(p'\) is the measured instrument direction,
  • \(E'\) is the measured or reconstructed energy,
  • \(t'\) is the measured arrival time,
  • \(p\) is the true photon arrival direction,
  • \(E\) is the true photon energy, and
  • \(t\) is the true trigger time.

The integration is done over all relevant true sky directions \(p\).

Depending on the type of the source model the method branches to the following methods to perform the actual computations

The method implements a caching mechanism for spatial models that have all parameters fixed. For those models the instrument response for a given event and observation is only computed once and then stored in an internal cache from which it is fetched back in case that the method is called again for the same event and observation.

Definition at line 696 of file GResponse.cpp.

References GModelSpatial::code(), GResponseVectorCache::contains(), GObservation::events(), GMODEL_SPATIAL_COMPOSITE, GMODEL_SPATIAL_DIFFUSE, GMODEL_SPATIAL_ELLIPTICAL, GMODEL_SPATIAL_POINT_SOURCE, GMODEL_SPATIAL_RADIAL, GModelSpatial::has_free_pars(), GObservation::id(), irf_composite(), irf_diffuse(), irf_elliptical(), irf_ptsrc(), irf_radial(), m_irf_vector_cache, m_use_irf_cache, GModel::name(), GResponseVectorCache::set(), GEvents::size(), GModel::size(), and GModelSky::spatial().

virtual double GResponse::nroi ( const GModelSky model,
const GEnergy obsEng,
const GTime obsTime,
const GObservation obs 
) const
pure virtual
GResponse & GResponse::operator= ( const GResponse rsp)
virtual

Assignment operator.

Parameters
[in]rspResponse.
Returns
Response.

Definition at line 139 of file GResponse.cpp.

References copy_members(), free_members(), and init_members().

Referenced by GCTAResponse::operator=(), GMWLResponse::operator=(), GLATResponse::operator=(), GCOMResponse::operator=(), and GSPIResponse::operator=().

virtual std::string GResponse::print ( const GChatter chatter = NORMAL) const
pure virtual

Print content of object.

Parameters
[in]chatterChattiness (defaults to NORMAL).
Returns
String containing the content of the object.

Formats the content in a standard way and puts this content in a C++ string that is returned.

Implements GBase.

Implemented in GCTAResponseCube, GCTAResponseIrf, GSPIResponse, GCOMResponse, GCTAResponse, GLATResponse, and GMWLResponse.

void GResponse::remove_response_cache ( const std::string &  name)
virtual

Remove response cache for model.

Parameters
[in]nameModel name.

Remove response cache for model name from response cache.

Definition at line 766 of file GResponse.cpp.

References m_irf_cache, m_irf_vector_cache, m_nroi_cache, GResponseVectorCache::remove(), and GResponseCache::remove().

int GResponse::size_edisp_vector ( const GModelSky model,
const GObservation obs,
const bool &  grad 
) const
protected

Return size of vector for energy dispersion computation.

Parameters
[in]modelSky model.
[in]obsObservation.
[out]gradSignals whether gradients should be computed.
Returns
Size of vector for energy dispersion computation.

Computes the size of the vector that will be computed for the computation of the energy dispersion.

Definition at line 1836 of file GResponse.cpp.

References GOptimizerPar::has_grad(), GModel::has_scales(), GObservation::instrument(), GOptimizerPar::is_free(), GOptimizerPar::name(), GModel::scale(), GModel::scales(), GModelTemporal::size(), GModelSpectral::size(), GModelSky::spectral(), and GModelSky::temporal().

Referenced by convolve().

virtual bool GResponse::use_edisp ( void  ) const
pure virtual
virtual bool GResponse::use_tdisp ( void  ) const
pure virtual

Member Data Documentation

double GResponse::m_irf_diffuse_resolution
protected

Angular resolution for diffuse model.

Definition at line 319 of file GResponse.hpp.

Referenced by copy_members(), init_members(), and irf_diffuse().

int GResponse::m_irf_elliptical_iter_phi
protected

Elliptical model integration phi iterations.

Definition at line 318 of file GResponse.hpp.

Referenced by copy_members(), init_members(), and irf_elliptical().

int GResponse::m_irf_elliptical_iter_theta
protected

Elliptical model integration theta iterations.

Definition at line 317 of file GResponse.hpp.

Referenced by copy_members(), init_members(), and irf_elliptical().

int GResponse::m_irf_radial_iter_phi
protected

Radial model integration phi iterations.

Definition at line 316 of file GResponse.hpp.

Referenced by copy_members(), init_members(), and irf_radial().

int GResponse::m_irf_radial_iter_theta
protected

Radial model integration theta iterations.

Definition at line 315 of file GResponse.hpp.

Referenced by copy_members(), init_members(), and irf_radial().

GResponseVectorCache GResponse::m_irf_vector_cache
mutableprotected
GResponseCache GResponse::m_nroi_cache
mutableprotected
bool GResponse::m_use_irf_cache
protected

Control usage of irf cache.

Definition at line 313 of file GResponse.hpp.

Referenced by copy_members(), init_members(), GCTAResponseCube::irf_spatial(), and irf_spatial().

bool GResponse::m_use_nroi_cache
protected

Control usage of nroi cache.

Definition at line 314 of file GResponse.hpp.

Referenced by copy_members(), init_members(), and GCTAResponseIrf::nroi().


The documentation for this class was generated from the following files: