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

Radial Polynom CTA model class. More...

#include <GCTAModelRadialPolynom.hpp>

Inheritance diagram for GCTAModelRadialPolynom:
GCTAModelRadial GCTAModelSpatial GBase

Classes

class  integrand
 

Public Member Functions

 GCTAModelRadialPolynom (void)
 Void constructor. More...
 
 GCTAModelRadialPolynom (const std::vector< double > &coeffs)
 Constructor. More...
 
 GCTAModelRadialPolynom (const GXmlElement &xml)
 XML constructor. More...
 
 GCTAModelRadialPolynom (const GCTAModelRadialPolynom &model)
 Copy constructor. More...
 
virtual ~GCTAModelRadialPolynom (void)
 Destructor. More...
 
virtual GCTAModelRadialPolynomoperator= (const GCTAModelRadialPolynom &model)
 Assignment operator. More...
 
virtual void clear (void)
 Clear instance. More...
 
virtual GCTAModelRadialPolynomclone (void) const
 Clone instance. More...
 
virtual std::string classname (void) const
 Return class name. More...
 
virtual std::string type (void) const
 Return model type. More...
 
virtual double eval (const double &offset, const bool &gradients=false) const
 Evaluate function. More...
 
virtual GCTAInstDir mc (GRan &ran) const
 Returns MC instrument direction. More...
 
virtual double mc_max_value (const GCTAObservation &obs) const
 Return maximum function value for Monte Carlo simulations. More...
 
virtual double omega (void) const
 Returns integral over radial model (in steradians) More...
 
virtual void read (const GXmlElement &xml)
 Read model from XML element. More...
 
virtual void write (GXmlElement &xml) const
 Write model into XML element. More...
 
virtual std::string print (const GChatter &chatter=NORMAL) const
 Print point source information. More...
 
- Public Member Functions inherited from GCTAModelRadial
 GCTAModelRadial (void)
 Void constructor. More...
 
 GCTAModelRadial (const GCTAModelRadial &model)
 Copy constructor. More...
 
virtual ~GCTAModelRadial (void)
 Destructor. More...
 
virtual GCTAModelRadialoperator= (const GCTAModelRadial &model)
 Assignment operator. More...
 
virtual double eval (const GCTAInstDir &dir, const GEnergy &energy, const GTime &time, const bool &gradients=false) const
 Evaluate function. More...
 
virtual GCTAInstDir mc (const GEnergy &energy, const GTime &time, const GCTAObservation &obs, GRan &ran) const
 Returns MC instrument direction. More...
 
- Public Member Functions inherited from GCTAModelSpatial
 GCTAModelSpatial (void)
 Void constructor. More...
 
 GCTAModelSpatial (const GCTAModelSpatial &model)
 Copy constructor. More...
 
virtual ~GCTAModelSpatial (void)
 Destructor. More...
 
virtual GCTAModelSpatialoperator= (const GCTAModelSpatial &model)
 Assignment operator. More...
 
virtual GModelParoperator[] (const int &index)
 Returns model parameter. More...
 
virtual const GModelParoperator[] (const int &index) const
 Returns model parameter (const version) More...
 
virtual GModelParoperator[] (const std::string &name)
 Returns reference to model parameter. More...
 
virtual const GModelParoperator[] (const std::string &name) const
 Returns reference to model parameter (const version) More...
 
int size (void) const
 Return number of model parameters. More...
 
virtual double npred (const GEnergy &energy, const GTime &time, const GObservation &obs) const
 Return integral of spatial model component. 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 GCTAModelRadialPolynom &model)
 Copy class members. More...
 
void free_members (void)
 Delete class members. More...
 
void update_pars (void)
 Update parameter mapping. More...
 
- Protected Member Functions inherited from GCTAModelRadial
void init_members (void)
 Initialise class members. More...
 
void copy_members (const GCTAModelRadial &model)
 Copy class members. More...
 
void free_members (void)
 Delete class members. More...
 
- Protected Member Functions inherited from GCTAModelSpatial
void init_members (void)
 Initialise class members. More...
 
void copy_members (const GCTAModelSpatial &model)
 Copy class members. More...
 
void free_members (void)
 Delete class members. More...
 

Protected Attributes

std::vector< GModelParm_coeffs
 Coefficients. More...
 
- Protected Attributes inherited from GCTAModelSpatial
std::vector< GModelPar * > m_pars
 Parameter pointers. More...
 

Detailed Description

Radial Polynom CTA model class.

This class implements the radial function

\[f(\theta) = \sum_{i=0}^m c_i \theta^i\]

where \(\theta\) is the offset angle (in degrees), and \(c_i\) are the polynomial coefficients.

This function represents a Polynom in \(\theta\).

Definition at line 59 of file GCTAModelRadialPolynom.hpp.

Constructor & Destructor Documentation

GCTAModelRadialPolynom::GCTAModelRadialPolynom ( void  )

Void constructor.

Definition at line 72 of file GCTAModelRadialPolynom.cpp.

References init_members().

Referenced by clone().

GCTAModelRadialPolynom::GCTAModelRadialPolynom ( const std::vector< double > &  coeffs)
explicit
GCTAModelRadialPolynom::GCTAModelRadialPolynom ( const GXmlElement xml)
explicit

XML constructor.

Parameters
[in]xmlXML element.

Creates instance of a radial Polynom model by extracting information from an XML element. See GCTAModelRadialPolynom::read() for more information about the expected structure of the XML element.

Definition at line 133 of file GCTAModelRadialPolynom.cpp.

References init_members(), and read().

GCTAModelRadialPolynom::GCTAModelRadialPolynom ( const GCTAModelRadialPolynom model)

Copy constructor.

Parameters
[in]modelRadial Polynom model.

Definition at line 152 of file GCTAModelRadialPolynom.cpp.

References copy_members(), and init_members().

GCTAModelRadialPolynom::~GCTAModelRadialPolynom ( void  )
virtual

Destructor.

Definition at line 169 of file GCTAModelRadialPolynom.cpp.

References free_members().

Member Function Documentation

std::string GCTAModelRadialPolynom::classname ( void  ) const
inlinevirtual

Return class name.

Returns
String containing the class name ("GCTAModelRadialPolynom").

Implements GCTAModelRadial.

Definition at line 120 of file GCTAModelRadialPolynom.hpp.

void GCTAModelRadialPolynom::clear ( void  )
virtual
GCTAModelRadialPolynom * GCTAModelRadialPolynom::clone ( void  ) const
virtual

Clone instance.

Implements GCTAModelRadial.

Definition at line 241 of file GCTAModelRadialPolynom.cpp.

References GCTAModelRadialPolynom().

void GCTAModelRadialPolynom::copy_members ( const GCTAModelRadialPolynom model)
protected

Copy class members.

Parameters
[in]modelRadial Polynomian model.

Definition at line 727 of file GCTAModelRadialPolynom.cpp.

References m_coeffs, and update_pars().

Referenced by GCTAModelRadialPolynom(), and operator=().

double GCTAModelRadialPolynom::eval ( const double &  offset,
const bool &  gradients = false 
) const
virtual

Evaluate function.

Parameters
[in]offsetOffset angle [degrees].
[in]gradientsCompute gradients?
Returns
Function value

Evaluates the radial polynomial model for a given offset. The model is defined as

\[f(\theta) = \sum_{i=0}^m c_i \theta^i\]

where \(\theta\) is the offset angle (in degrees), and \(c_i\) are the polynomial coefficients.

If the gradients flag is true the method will also compute the partial derivatives of the parameters. The partial derivative of the Polynom are given by

\[\frac{df}{d{c_i}_v} = {c_i}_s \theta^i\]

where \({c_i}_v\) is the value part, \({c_i}_s\) is the scaling part, and \(c_i = {c_i}_v {c_i}_s\).

Implements GCTAModelRadial.

Definition at line 270 of file GCTAModelRadialPolynom.cpp.

References gammalib::is_infinite(), gammalib::is_notanumber(), and m_coeffs.

Referenced by GCTAModelRadialPolynom::integrand::eval(), and mc().

void GCTAModelRadialPolynom::free_members ( void  )
protected

Delete class members.

Definition at line 743 of file GCTAModelRadialPolynom.cpp.

Referenced by clear(), operator=(), and ~GCTAModelRadialPolynom().

void GCTAModelRadialPolynom::init_members ( void  )
protected

Initialise class members.

Definition at line 709 of file GCTAModelRadialPolynom.cpp.

References m_coeffs, and update_pars().

Referenced by clear(), GCTAModelRadialPolynom(), and operator=().

GCTAInstDir GCTAModelRadialPolynom::mc ( GRan ran) const
virtual

Returns MC instrument direction.

Parameters
[in,out]ranRandom number generator.
Returns
CTA instrument direction.

Draws an arbitrary CTA instrument position from

\[f(\theta) = \sin(\theta) \left( \sum_{i=0}^m c_i \theta^i \right)\]

where \(\theta\) is the offset angle (in degrees), and \(c_i\) are the polynomial coefficients, using the rejection method.

Note that the maximum offset angle is fixed by the constant g_cta_radial_polynom_offset_max. This needs eventually adjusting on real data. The main reason for the tight limit is to avoid divergence at large offset angles.

Todo:
This method actually assumes that the polynom is always < 1, which may not be necessarily the case. Ideally, the method should determine the maximum of the polynomial to throw events. This is a severe limitation and should rapidly be corrected.

Implements GCTAModelRadial.

Definition at line 353 of file GCTAModelRadialPolynom.cpp.

References cos(), gammalib::deg2rad, eval(), g_cta_radial_polynom_offset_max, sin(), and GRan::uniform().

double GCTAModelRadialPolynom::mc_max_value ( const GCTAObservation obs) const
virtual

Return maximum function value for Monte Carlo simulations.

Parameters
[in]obsCTA Observation (not used).
Returns
Maximum function value for Monte Carlo simulations.

Implements GCTAModelRadial.

Definition at line 420 of file GCTAModelRadialPolynom.cpp.

References gammalib::deg2rad, g_cta_radial_polynom_offset_max, and sin().

double GCTAModelRadialPolynom::omega ( void  ) const
virtual

Returns integral over radial model (in steradians)

Computes

\[\Omega = 2 \pi \int_0^{\pi} \sin \theta f(\theta) d\theta\]

where

\[f(\theta) = \sum_{i=0}^m c_i \theta^i\]

\(\theta\) is the offset angle (in degrees), and \(c_i\) are the polynomial coefficients.

The integration is performed numerically, and the upper integration bound \(\pi\) is fixed to < g_cta_radial_polynom_offset_max. This needs eventually adjusting on real data. The main reason for the tight limit is to avoid divergence at large offset angles.

Implements GCTAModelRadial.

Definition at line 447 of file GCTAModelRadialPolynom.cpp.

References gammalib::deg2rad, g_cta_radial_polynom_offset_max, GIntegral::romberg(), and gammalib::twopi.

GCTAModelRadialPolynom & GCTAModelRadialPolynom::operator= ( const GCTAModelRadialPolynom model)
virtual

Assignment operator.

Parameters
[in]modelRadial Polynom model.

Definition at line 190 of file GCTAModelRadialPolynom.cpp.

References copy_members(), free_members(), init_members(), and GCTAModelRadial::operator=().

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

Print point source information.

Parameters
[in]chatterChattiness (defaults to NORMAL).
Returns
String containing point source information.

Implements GCTAModelRadial.

Definition at line 675 of file GCTAModelRadialPolynom.cpp.

References GCTAModelSpatial::m_pars, gammalib::parformat(), SILENT, GCTAModelSpatial::size(), and gammalib::str().

void GCTAModelRadialPolynom::read ( const GXmlElement xml)
virtual

Read model from XML element.

Parameters
[in]xmlXML element.
Exceptions
GException::invalid_valueInvalid XML format encountered.

Read the radial Polynom model information from an XML element in the following format

<radialModel type="...">
  <parameter name="Coeff0"  scale="1.0" value="1.5" min="0.1" max="10.0" free="1"/>
  <parameter name="Coeff1"  scale="1.0" value="3.0" min="1.0" max="10.0" free="1"/>
  <parameter name="Coeff2"  scale="1.0" value="5.0" min="1.0" max="10.0" free="1"/>
  ...
</radialModel>
Todo:
Implement a test of the coefficient boundaries.

Implements GCTAModelRadial.

Definition at line 487 of file GCTAModelRadialPolynom.cpp.

References GXmlElement::attribute(), GXmlNode::element(), GXmlNode::elements(), G_READ, m_coeffs, GOptimizerPar::name(), GModelPar::read(), gammalib::str(), gammalib::toint(), GOptimizerPar::unit(), and update_pars().

Referenced by GCTAModelRadialPolynom().

std::string GCTAModelRadialPolynom::type ( void  ) const
inlinevirtual

Return model type.

Returns
Model type "Polynom".

Implements GCTAModelRadial.

Definition at line 132 of file GCTAModelRadialPolynom.hpp.

Referenced by write().

void GCTAModelRadialPolynom::update_pars ( void  )
protected

Update parameter mapping.

Definition at line 753 of file GCTAModelRadialPolynom.cpp.

References m_coeffs, and GCTAModelSpatial::m_pars.

Referenced by copy_members(), GCTAModelRadialPolynom(), init_members(), and read().

void GCTAModelRadialPolynom::write ( GXmlElement xml) const
virtual

Write model into XML element.

Parameters
[in]xmlXML element.

Write the polynomial radial model information into an XML element in the following format

<radialModel type="...">
  <parameter name="Coeff0"  scale="1.0" value="1.5" min="0.1" max="10.0" free="1"/>
  <parameter name="Coeff1"  scale="1.0" value="3.0" min="1.0" max="10.0" free="1"/>
  <parameter name="Coeff2"  scale="1.0" value="5.0" min="1.0" max="10.0" free="1"/>
  ...
</radialModel>

Implements GCTAModelRadial.

Definition at line 642 of file GCTAModelRadialPolynom.cpp.

References G_WRITE, m_coeffs, gammalib::str(), type(), gammalib::xml_check_type(), and gammalib::xml_need_par().

Member Data Documentation

std::vector<GModelPar> GCTAModelRadialPolynom::m_coeffs
protected

Coefficients.

Definition at line 110 of file GCTAModelRadialPolynom.hpp.

Referenced by copy_members(), eval(), GCTAModelRadialPolynom(), init_members(), read(), update_pars(), and write().


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