Table Of Contents

Previous topic

Observation handling

Next topic

Sky maps, sky coordinates and sky regions

Model handling

Overview

The following figure presents an overview over the C++ classes of the model module and their relations.

../../_images/uml_model.png

Overview over the model module

The central C++ class of the model module is the abstract base class GModel which defines a model component. Model components are combined using the GModels container C++ class.

Evaluating models for observations

Models are evaluated by observations using the GObservation::model method. The GObservation::npred provides the number of predicted events by all models within the relevant dataspace for a given observation.

Call tree for model evaluation shows the call tree of the GObservation::model (left) and GObservation::npred (right) methods.

For model evaluation, the :doxy:GObservation::model method loops over all models in the GModels container and calls the GModel::eval_gradients() method of each model. For GModelSky models this method invokes the GResponse::convolve() method that calls among others the GResponse::irf() method. That method needs to be implement for each instruments.

The GObservation::npred method performs an integration over observed time, observed energy and observed arrival direction. The integrate over observed time is for the moment a simple multiplication with ontime. The integral over observed energy is performed by the GObservation::npred_spec method that calls the GModel::eval_gradients() method. For GModelSky models this method calls the GResponse::nroi() method that needs to be implement for each instruments.

../../_images/calltree_model.png

Call tree for model evaluation

Sky models

Sky models describe the spatial, spectral and temporal properties of a gamma-ray source using the following factorization:

\[M(\alpha,\delta,E,t) = M_{\rm spatial}(\alpha,\delta | E,t) \times M_{\rm spectral}(E | t) \times M_{\rm temporal}(t)\]

The spatial model component \(M_{\rm spatial}(\alpha,\delta|E,t)\) is defined by the abstract GModelSpatial class, the spectral model component \(M_{\rm spectral}(E|t)\) is defined by the abstract GModelSpectral class and the temporal component \(M_{\rm temporal}(t)\) is defined by the abstract GModelTemporal class.

The spatial model component describes the energy and time dependent morphology of the source. It satisfies

\[\int_{\Omega} M_{\rm spatial}(\alpha,\delta|E,t) d\Omega = 1\]

for all \(E\) and \(t\), hence the spatial component does not impact the spatially integrated spectral and temporal properties of the source (the integration is done here over the spatial parameters \(\alpha\) and \(\delta\) in a spherical coordinate system).

The spectral model component describes the spatially integrated time dependent spectral distribution of the source. It satisfies

\[\int_{E} M_{\rm spectral}(E | t) dE = \Phi\]

for all \(t\), where \(\Phi\) is the spatially and spectrally integrated total source flux. The spectral component does not impact the temporal properties of the integrated flux \(\Phi\).

The temporal model component describes the relative variation of the source flux with respect to the mean value given by the spectral model component.

Spatial components

Point source

<source name="Crab" type="PointSource">
  <spatialModel type="SkyDirFunction">
    <parameter name="RA"  scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
    <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="1"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Radial disk

<source name="Crab" type="ExtendedSource">
  <spatialModel type="DiskFunction">
    <parameter name="RA"     scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
    <parameter name="DEC"    scale="1.0" value="22.0145" min="-90"  max="90"  free="1"/>
    <parameter name="Radius" scale="1.0" value="0.20"    min="0.01" max="10"  free="1"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Radial Gaussian

<source name="Crab" type="ExtendedSource">
  <spatialModel type="GaussFunction">
    <parameter name="RA"    scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
    <parameter name="DEC"   scale="1.0" value="22.0145" min="-90"  max="90"  free="1"/>
    <parameter name="Sigma" scale="1.0" value="0.20"    min="0.01" max="10"  free="1"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Radial shell

<source name="Crab" type="ExtendedSource">
  <spatialModel type="ShellFunction">
    <parameter name="RA"     scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
    <parameter name="DEC"    scale="1.0" value="22.0145" min="-90"  max="90"  free="1"/>
    <parameter name="Radius" scale="1.0" value="0.30"    min="0.01" max="10"  free="1"/>
    <parameter name="Width"  scale="1.0" value="0.10"    min="0.01" max="10"  free="1"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Elliptical disk

<source name="Crab" type="ExtendedSource">
  <spatialModel type="EllipticalDisk">
    <parameter name="RA"          scale="1.0" value="83.6331" min="-360"  max="360" free="1"/>
    <parameter name="DEC"         scale="1.0" value="22.0145" min="-90"   max="90"  free="1"/>
    <parameter name="PA"          scale="1.0" value="45.0"    min="-360"  max="360" free="1"/>
    <parameter name="MinorRadius" scale="1.0" value="0.5"     min="0.001" max="10"  free="1"/>
    <parameter name="MajorRadius" scale="1.0" value="2.0"     min="0.001" max="10"  free="1"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Elliptical Gaussian

<source name="Crab" type="ExtendedSource">
  <spatialModel type="EllipticalGauss">
    <parameter name="RA"          scale="1.0" value="83.6331" min="-360"  max="360" free="1"/>
    <parameter name="DEC"         scale="1.0" value="22.0145" min="-90"   max="90"  free="1"/>
    <parameter name="PA"          scale="1.0" value="45.0"    min="-360"  max="360" free="1"/>
    <parameter name="MinorRadius" scale="1.0" value="0.5"     min="0.001" max="10"  free="1"/>
    <parameter name="MajorRadius" scale="1.0" value="2.0"     min="0.001" max="10"  free="1"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Isotropic source

<source name="Crab" type="DiffuseSource">
  <spatialModel type="ConstantValue">
     <parameter name="Value" scale="1" value="1" min="1"  max="1" free="0"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Diffuse map

<source name="Crab" type="DiffuseSource">
  <spatialModel type="SpatialMap" file="map.fits">
     <parameter name="Prefactor" scale="1" value="1" min="0.001" max="1000.0" free="0"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Diffuse map cube

<source name="Crab" type="DiffuseSource">
  <spatialModel type="MapCubeFunction" file="map_cube.fits">
    <parameter name="Normalization" scale="1" value="1" min="0.001" max="1000.0" free="0"/>
  </spatialModel>
  <spectrum type="...">
    ...
  </spectrum>
</source>

Spectral components

Constant

The GModelSpectralConst class implements the constant function

\[\frac{dN}{dE} = N_0\]

where the parameters in the XML definition have the following mappings:

  • \(N_0\) = Normalization

The XML format for specifying a constant is:

<spectrum type="ConstantValue">
 <parameter name="Normalization" scale="1e-16" value="5.7"  min="1e-07" max="1000.0" free="1"/>
</spectrum>

Node function

The generalisation of the broken power law is the node function, which is defined by a set of energy and intensity values, the so called nodes, which are connected by power laws.

The XML format for specifying a node function is:

<spectrum type="NodeFunction">
 <node>
   <parameter scale="1.0"   name="Energy"    min="0.1"   max="1.0e20" value="1.0"  free="0"/>
   <parameter scale="1e-07" name="Intensity" min="1e-07" max="1000.0" value="1.0"  free="1"/>
 </node>
 <node>
   <parameter scale="1.0"   name="Energy"    min="0.1"   max="1.0e20" value="10.0" free="0"/>
   <parameter scale="1e-07" name="Intensity" min="1e-07" max="1000.0" value="0.1"  free="1"/>
 </node>
</spectrum>

(in this example there are two nodes; the number of nodes in a node function is arbitrary).

File function

A function defined using an input ASCII file with columns of energy and differential flux values. The energy values are assumed to be in units of MeV, the flux values are normally assumed to be in units of \({\rm cm}^{-2} {\rm s}^{-1} {\rm MeV}^{-1}\). The only exception to this rule is the isotropic diffuse model GModelSpatialDiffuseConst for which the flux values are given in units of \({\rm cm}^{-2} {\rm s}^{-1} {\rm MeV}^{-1} {\rm sr}^{-1}\).

The only parameter of the model is a multiplicative normalization:

\[\frac{dN}{dE} = N_0 \left. \frac{dN}{dE} \right\rvert_{\rm file}\]

where the parameters in the XML definition have the following mappings:

  • \(N_0\) = Normalization

The XML format for specifying a file function is:

<spectrum type="FileFunction" file="data/filefunction.txt">
 <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
</spectrum>

If the file is given as relative path, the path is relative to the working directory of the executable. Alternatively, an absolute path may be specified. Any environment variable present in the path name will be expanded.

Power law

The GModelSpectralPlaw class implements the power law function

\[\frac{dN}{dE} = k_0 \left( \frac{E}{E_0} \right)^{\gamma}\]

where the parameters in the XML definition have the following mappings:

  • \(k_0\) = Prefactor
  • \(\gamma\) = Index
  • \(E_0\) = Scale

The XML format for specifying a power law is:

<spectrum type="PowerLaw">
 <parameter name="Prefactor" scale="1e-16" value="5.7"  min="1e-07" max="1000.0" free="1"/>
 <parameter name="Index"     scale="-1"    value="2.48" min="0.0"   max="+5.0"   free="1"/>
 <parameter name="Scale"     scale="1e6"   value="0.3"  min="0.01"  max="1000.0" free="0"/>
</spectrum>

An alternative power law function is defined by the GModelSpectralPlaw2 class that uses the integral flux as parameter rather than the Prefactor:

\[\frac{dN}{dE} = \frac{N(\gamma+1)E^{\gamma}} {E_{\rm max}^{\gamma+1} - E_{\rm min}^{\gamma+1}}\]

where the parameters in the XML definition have the following mappings:

  • \(N\) = Integral
  • \(\gamma\) = Index
  • \(E_{\rm min}\) = LowerLimit
  • \(E_{\rm max}\) = UpperLimit

The XML format for specifying a power law defined by the integral flux is:

<spectrum type="PowerLaw2">
 <parameter scale="1e-07" name="Integral"   min="1e-07" max="1000.0"    value="1.0" free="1"/>
 <parameter scale="1.0"   name="Index"      min="-5.0"  max="+5.0"      value="-2.0" free="1"/>
 <parameter scale="1.0"   name="LowerLimit" min="10.0"  max="1000000.0" value="100.0" free="0"/>
 <parameter scale="1.0"   name="UpperLimit" min="10.0"  max="1000000.0" value="500000.0" free="0"/>
</spectrum>

Note

The UpperLimit and LowerLimit parameters are always treated as fixed and, as should be apparent from this definition, the flux given by the Integral parameter is over the range (LowerLimit, UpperLimit). Use of this model allows the errors on the integrated flux to be evaluated directly by likelihood, obviating the need to propagate the errors if one is using the PowerLaw form.

Exponentially cut-off power law

The GModelSpectralExpPlaw class implements the exponentially cut-off power law function

\[\frac{dN}{dE} = k_0 \left( \frac{E}{E_0} \right)^{\gamma} \exp \left( \frac{-E}{E_{\rm cut}} \right)\]

where the parameters in the XML definition have the following mappings:

  • \(k_0\) = Prefactor
  • \(\gamma\) = Index
  • \(E_0\) = Scale
  • \(E_{\rm cut}\) = Cutoff

The XML format for specifying an exponentially cut-off power law is:

<spectrum type="ExpCutoff">
 <parameter name="Prefactor" scale="1e-16" value="5.7"  min="1e-07" max="1000.0" free="1"/>
 <parameter name="Index"     scale="-1"    value="2.48" min="0.0"   max="+5.0"   free="1"/>
 <parameter name="Cutoff"    scale="1e6"   value="1.0"  min="0.01"  max="1000.0" free="1"/>
 <parameter name="Scale"     scale="1e6"   value="0.3"  min="0.01"  max="1000.0" free="0"/>
</spectrum>

Super exponentially cut-off power law

The GModelSpectralSuperExpPlaw class implements the super exponentially cut-off power law function

\[\frac{dN}{dE} = k_0 \left( \frac{E}{E_0} \right)^{\gamma} \exp \left( -\left( \frac{E}{E_{\rm cut}} \right)^{\alpha} \right)\]

where the parameters in the XML definition have the following mappings:

  • \(k_0\) = Prefactor
  • \(\gamma\) = Index1
  • \(\alpha\) = Index2
  • \(E_0\) = Scale
  • \(E_{\rm cut}\) = Cutoff
<spectrum type="PLSuperExpCutoff">
 <parameter name="Prefactor" scale="1e-16" value="1.0" min="1e-07" max="1000.0" free="1"/>
 <parameter name="Index1"    scale="-1"    value="2.0" min="0.0"   max="+5.0"   free="1"/>
 <parameter name="Cutoff"    scale="1e6"   value="1.0" min="0.01"  max="1000.0" free="1"/>
 <parameter name="Scale"     scale="1e6"   value="1.0" min="0.01"  max="1000.0" free="0"/>
 <parameter name="Index2"    scale="1.0"   value="1.5" min="0.1"   max="5.0"    free="1"/>
</spectrum>

Broken power law

The GModelSpectralBrokenPlaw class implements the broken power law function

\[\begin{split}\frac{dN}{dE} = k_0 \times \left \{ \begin{eqnarray} \left( \frac{E}{E_b} \right)^{\gamma_1} & {\rm if\,\,} E < E_b \\ \left( \frac{E}{E_b} \right)^{\gamma_2} & {\rm otherwise} \end{eqnarray} \right .\end{split}\]

where the parameters in the XML definition have the following mappings:

  • \(k_0\) = Prefactor
  • \(\gamma_1\) = Index1
  • \(\gamma_2\) = Index2
  • \(E_b\) = BreakValue

The XML format for specifying a broken power law is:

<spectrum type="BrokenPowerLaw">
 <parameter name="Prefactor"  scale="1e-16" value="5.7"  min="1e-07" max="1000.0" free="1"/>
 <parameter name="Index1"     scale="-1"    value="2.48" min="0.0"   max="+5.0"   free="1"/>
 <parameter name="BreakValue" scale="1e6"   value="0.3"  min="0.01"  max="1000.0" free="1"/>
 <parameter name="Index2"     scale="-1"    value="2.70" min="0.01"  max="1000.0" free="1"/>
</spectrum>

Gaussian

The GModelSpectralGauss class implements the gaussian function

\[\frac{dN}{dE} = \frac{N_0}{\sqrt{2\pi}\sigma} \exp \left( \frac{-(E-\bar{E})^2}{2 \sigma^2} \right)\]

where the parameters in the XML definition have the following mappings:

  • \(N_0\) = Normalization
  • \(\bar{E}\) = Mean
  • \(\sigma\) = Sigma

The XML format for specifying a Gaussian is:

<spectrum type="Gaussian">
 <parameter name="Normalization" scale="1e-10" value="1.0"  min="1e-07" max="1000.0" free="1"/>
 <parameter name="Mean"          scale="1e6"   value="5.0"  min="0.01"  max="100.0"  free="1"/>
 <parameter name="Sigma"         scale="1e6"   value="1.0"  min="0.01"  max="100.0"  free="1"/>
</spectrum>

Log parabola

The GModelSpectralLogParabola class implements the log parabola function

\[\frac{dN}{dE} = k_0 \left( \frac{E}{E_0} \right)^{\gamma+\eta \ln(E/E_0)}\]

where the parameters in the XML definition have the following mappings:

  • \(k_0\) = Prefactor
  • \(\gamma\) = Index
  • \(\eta\) = Curvature
  • \(E_0\) = Scale

The XML format for specifying a log parabola spectrum is:

<spectrum type="LogParabola">
 <parameter name="Prefactor" scale="1e-17" value="5.878"   min="1e-07" max="1000.0" free="1"/>
 <parameter name="Index"     scale="-1"    value="2.32473" min="0.0"   max="+5.0"   free="1"/>
 <parameter name="Curvature" scale="-1"    value="0.074"   min="-5.0"  max="+5.0"   free="1"/>
 <parameter name="Scale"     scale="1e6"   value="1.0"     min="0.01"  max="1000.0" free="0"/>
</spectrum>

An alternative XML format is supported for compatibility with the Fermi/LAT XML format:

<spectrum type="LogParabola">
 <parameter name="Prefactor" scale="1e-17" value="5.878"   min="1e-07" max="1000.0" free="1"/>
 <parameter name="alpha"     scale="1"     value="2.32473" min="0.0"   max="+5.0"   free="1"/>
 <parameter name="beta"      scale="1"     value="0.074"   min="-5.0"  max="+5.0"   free="1"/>
 <parameter name="Scale"     scale="1e6"   value="1.0"     min="0.01"  max="1000.0" free="0"/>
</spectrum>

where

  • alpha = -Index
  • beta = -Curvature

Temporal components

Constant

The only temporal model component that exists so far is the constant model that is implemented by the GModelTemporalConst class.

The XML format for specifying a constant temporal model is:

<temporalModel type="Constant">
  <parameter name="Normalization" scale="1.0" value="1.0" min="0.1" max="10.0" free="0"/>
</temporalModel>

Background models