GammaLib  2.0.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GModelSpectralSuperExpPlaw.hpp
Go to the documentation of this file.
1 /***************************************************************************
2  * GModelSpectralSuperExpPlaw.hpp - Super exp. cut off power law model *
3  * ----------------------------------------------------------------------- *
4  * copyright (C) 2014-2021 by Michael Mayer *
5  * ----------------------------------------------------------------------- *
6  * *
7  * This program is free software: you can redistribute it and/or modify *
8  * it under the terms of the GNU General Public License as published by *
9  * the Free Software Foundation, either version 3 of the License, or *
10  * (at your option) any later version. *
11  * *
12  * This program is distributed in the hope that it will be useful, *
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15  * GNU General Public License for more details. *
16  * *
17  * You should have received a copy of the GNU General Public License *
18  * along with this program. If not, see <http://www.gnu.org/licenses/>. *
19  * *
20  ***************************************************************************/
21 /**
22  * @file GModelSpectralSuperExpPlaw.hpp
23  * @brief Super exponential cut off power law spectral class interface definition
24  * @author Michael Mayer
25  */
26 
27 #ifndef GMODELSPECTRALSUPEREXPPLAW_HPP
28 #define GMODELSPECTRALSUPEREXPPLAW_HPP
29 
30 /* __ Includes ___________________________________________________________ */
31 #include <string>
32 #include "GModelSpectral.hpp"
33 #include "GFunction.hpp"
34 #include "GEnergy.hpp"
35 #include "GModelPar.hpp"
36 
37 /* __ Forward declarations _______________________________________________ */
38 class GRan;
39 class GTime;
40 class GXmlElement;
41 
42 
43 /***********************************************************************//**
44  * @class GModelSpectralSuperExpPlaw
45  *
46  * @brief Super exponential cut off power law spectral class
47  *
48  * This class implements a power law spectrum with super exponential cut off.
49  * The model is defined by
50  *
51  * \f[
52  * S_{\rm E}(E | t) = {\tt m\_norm}
53  * \left( \frac{E}{\tt m\_pivot} \right)^{\tt m\_index1}
54  * \exp \left( - \left( \frac{E}{\tt m\_ecut} \right)^{\tt m\_index2}
55  * \right)
56  * \f]
57  *
58  * where
59  * - \f${\tt m\_norm}\f$ is the normalization or prefactor,
60  * - \f${\tt m\_pivot}\f$ is the pivot energy,
61  * - \f${\tt m\_index1}\f$ is the spectral index
62  * - \f${\tt m\_ecut}\f$ is the cut off energy, and
63  * - \f${\tt m\_index2}\f$ is the index of the cut off.
64  ***************************************************************************/
66 
67 public:
68  // Constructors and destructors
70  GModelSpectralSuperExpPlaw(const std::string& type,
71  const std::string& prefactor,
72  const std::string& index1,
73  const std::string& pivot,
74  const std::string& cutoff,
75  const std::string& index2);
76  GModelSpectralSuperExpPlaw(const double& prefactor,
77  const double& index1,
78  const GEnergy& pivot,
79  const GEnergy& cutoff,
80  const double& index2);
81  explicit GModelSpectralSuperExpPlaw(const GXmlElement& xml);
83  virtual ~GModelSpectralSuperExpPlaw(void);
84 
85  // Operators
87 
88  // Implemented pure virtual methods
89  virtual void clear(void);
90  virtual GModelSpectralSuperExpPlaw* clone(void) const;
91  virtual std::string classname(void) const;
92  virtual std::string type(void) const;
93  virtual double eval(const GEnergy& srcEng,
94  const GTime& srcTime = GTime(),
95  const bool& gradients = false) const;
96  virtual double flux(const GEnergy& emin,
97  const GEnergy& emax) const;
98  virtual double eflux(const GEnergy& emin,
99  const GEnergy& emax) const;
100  virtual GEnergy mc(const GEnergy& emin,
101  const GEnergy& emax,
102  const GTime& time,
103  GRan& ran) const;
104  virtual void read(const GXmlElement& xml);
105  virtual void write(GXmlElement& xml) const;
106  virtual std::string print(const GChatter& chatter = NORMAL) const;
107 
108  // Other methods
109  void type(const std::string& type);
110  double prefactor(void) const;
111  void prefactor(const double& prefactor);
112  double index1(void) const;
113  void index1(const double& index1);
114  double index2(void) const;
115  void index2(const double& index2);
116  GEnergy cutoff(void) const;
117  void cutoff(const GEnergy& cutoff);
118  GEnergy pivot(void) const;
119  void pivot(const GEnergy& pivot);
120 
121 protected:
122  // Protected methods
123  void init_members(void);
124  void copy_members(const GModelSpectralSuperExpPlaw& model);
125  void free_members(void);
126  void update_eval_cache(const GEnergy& energy) const;
127  void update_mc_cache(const GEnergy& emin, const GEnergy& emax) const;
128 
129  // Photon flux integration kernel
130  class flux_kernel : public GFunction {
131  public:
132  flux_kernel(const double& norm,
133  const double& index1,
134  const double& pivot,
135  const double& ecut,
136  const double& index2) :
137  m_norm(norm),
138  m_index1(index1),
139  m_inv_pivot(1.0/pivot),
140  m_inv_ecut(1.0/ecut),
141  m_index2(index2) {}
142  double eval(const double& eng);
143  protected:
144  double m_norm; //!< Normalization
145  double m_index1; //!< Index1
146  double m_inv_pivot; //!< 1 / Pivot energy
147  double m_inv_ecut; //!< 1 / Cut off energy
148  double m_index2; //!< Index2
149  };
150 
151  // Energy flux integration kernel
152  class eflux_kernel : public GFunction {
153  public:
154  eflux_kernel(const double& norm,
155  const double& index1,
156  const double& pivot,
157  const double& ecut,
158  const double& index2) :
159  m_norm(norm),
160  m_index1(index1),
161  m_inv_pivot(1.0/pivot),
162  m_inv_ecut(1.0/ecut),
163  m_index2(index2) {}
164  double eval(const double& eng);
165  protected:
166  double m_norm; //!< Normalization
167  double m_index1; //!< Index1
168  double m_inv_pivot; //!< 1 / Pivot energy
169  double m_inv_ecut; //!< 1 / Cut off energy
170  double m_index2; //!< Index2
171  };
172 
173  // Protected members
174  std::string m_type; //!< Model type
175  GModelPar m_norm; //!< Normalization factor
176  GModelPar m_index1; //!< Spectral index
177  GModelPar m_index2; //!< Index of cutoff
178  GModelPar m_ecut; //!< Exponential cut off energy
179  GModelPar m_pivot; //!< Pivot energy
180 
181  // Cached members used for pre-computations
182  mutable GEnergy m_last_energy; //!< Last energy value
183  mutable double m_last_index1; //!< Last index1 parameter
184  mutable double m_last_ecut; //!< Last energy cut-off parameter
185  mutable double m_last_pivot; //!< Last pivot parameter
186  mutable double m_last_index2; //!< Last index2 parameter
187  mutable double m_last_e_norm; //!< Last E/Epivot value
188  mutable double m_last_e_cut; //!< Last E/Ecut value
189  mutable double m_last_exponent; //!< last pow(E/Ecut,index2) value
190  mutable double m_last_power; //!< Last power value
191  mutable double m_mc_emin; //!< Minimum energy
192  mutable double m_mc_emax; //!< Maximum energy
193  mutable double m_mc_exponent; //!< Exponent (index+1)
194  mutable double m_mc_pow_emin; //!< Power of minimum energy
195  mutable double m_mc_pow_ewidth; //!< Power of energy width
196 };
197 
198 
199 /***********************************************************************//**
200  * @brief Return class name
201  *
202  * @return String containing the class name ("GModelSpectralSuperExpPlaw").
203  ***************************************************************************/
204 inline
206 {
207  return ("GModelSpectralSuperExpPlaw");
208 }
209 
210 
211 /***********************************************************************//**
212  * @brief Return model type
213  *
214  * @return "SuperExpCutoff".
215  *
216  * Returns the type of the super exponentially cut off power law model.
217  ***************************************************************************/
218 inline
219 std::string GModelSpectralSuperExpPlaw::type(void) const
220 {
221  return (m_type);
222 }
223 
224 
225 /***********************************************************************//**
226  * @brief Set model type
227  *
228  * @param[in] type Model type.
229  *
230  * Set the type of the super exponentially cut off power law model.
231  ***************************************************************************/
232 inline
233 void GModelSpectralSuperExpPlaw::type(const std::string& type)
234 {
235  m_type = type;
236  return;
237 }
238 
239 
240 /***********************************************************************//**
241  * @brief Return pre factor
242  *
243  * @return Pre factor (ph/cm2/s/MeV).
244  *
245  * Returns the pre factor.
246  ***************************************************************************/
247 inline
249 {
250  return (m_norm.value());
251 }
252 
253 
254 /***********************************************************************//**
255  * @brief Set pre factor
256  *
257  * @param[in] prefactor Pre factor (ph/cm2/s/MeV).
258  *
259  * Sets the pre factor.
260  ***************************************************************************/
261 inline
262 void GModelSpectralSuperExpPlaw::prefactor(const double& prefactor)
263 {
264  m_norm.value(prefactor);
265  return;
266 }
267 
268 
269 /***********************************************************************//**
270  * @brief Return power law index
271  *
272  * @return Power law index.
273  *
274  * Returns the power law index.
275  ***************************************************************************/
276 inline
278 {
279  return (m_index1.value());
280 }
281 
282 
283 /***********************************************************************//**
284  * @brief Set power law index
285  *
286  * @param[in] index1 Power law index.
287  *
288  * Sets the power law index.
289  ***************************************************************************/
290 inline
291 void GModelSpectralSuperExpPlaw::index1(const double& index1)
292 {
293  m_index1.value(index1);
294  return;
295 }
296 
297 /***********************************************************************//**
298  * @brief Return cut off index
299  *
300  * @return cut off index.
301  *
302  * Returns the cut off index.
303  ***************************************************************************/
304 inline
306 {
307  return (m_index2.value());
308 }
309 
310 
311 /***********************************************************************//**
312  * @brief Set cut off index
313  *
314  * @param[in] index2 Cut off index.
315  *
316  * Sets the cut off index.
317  ***************************************************************************/
318 inline
319 void GModelSpectralSuperExpPlaw::index2(const double& index2)
320 {
321  m_index2.value(index2);
322  return;
323 }
324 
325 
326 /***********************************************************************//**
327  * @brief Return pivot energy
328  *
329  * @return Pivot energy.
330  *
331  * Returns the pivot energy.
332  ***************************************************************************/
333 inline
335 {
336  GEnergy energy;
337  energy.MeV(m_pivot.value());
338  return energy;
339 }
340 
341 
342 /***********************************************************************//**
343  * @brief Set pivot energy
344  *
345  * @param[in] pivot Pivot energy.
346  *
347  * Sets the pivot energy.
348  ***************************************************************************/
349 inline
351 {
352  m_pivot.value(pivot.MeV());
353  return;
354 }
355 
356 
357 /***********************************************************************//**
358  * @brief Return exponential cut-off energy
359  *
360  * @return Exponential cut-off energy.
361  *
362  * Returns the exponential cut-off energy.
363  ***************************************************************************/
364 inline
366 {
367  GEnergy energy;
368  energy.MeV(m_ecut.value());
369  return energy;
370 }
371 
372 
373 /***********************************************************************//**
374  * @brief Set exponential cut-off energy
375  *
376  * @param[in] cutoff Exponential cut-off energy.
377  *
378  * Sets the exponential cut-off energy.
379  ***************************************************************************/
380 inline
382 {
383  m_ecut.value(cutoff.MeV());
384  return;
385 }
386 
387 #endif /* GMODELSPECTRALSUPEREXPPLAW_HPP */
void free_members(void)
Delete class members.
GModelPar m_index2
Index of cutoff.
double norm(const GVector &vector)
Computes vector norm.
Definition: GVector.cpp:864
double eval(const double &eng)
Kernel for energy flux integration.
GModelPar m_ecut
Exponential cut off energy.
Energy value class definition.
virtual void clear(void)
Clear exponentially cut off power law model.
Abstract spectral model base class.
void update_mc_cache(const GEnergy &emin, const GEnergy &emax) const
Update Monte Carlo pre computation cache.
double m_last_exponent
last pow(E/Ecut,index2) value
GEnergy pivot(void) const
Return pivot energy.
virtual double flux(const GEnergy &emin, const GEnergy &emax) const
Returns model photon flux between emin, emax
double m_last_index2
Last index2 parameter.
void update_eval_cache(const GEnergy &energy) const
Update eval precomputation cache.
double m_last_e_norm
Last E/Epivot value.
virtual std::string type(void) const
Return model type.
double m_mc_pow_ewidth
Power of energy width.
XML element node class.
Definition: GXmlElement.hpp:48
Random number generator class.
Definition: GRan.hpp:44
virtual double eval(const GEnergy &srcEng, const GTime &srcTime=GTime(), const bool &gradients=false) const
Evaluate function.
double MeV(void) const
Return energy in MeV.
Definition: GEnergy.cpp:321
Time class.
Definition: GTime.hpp:55
void init_members(void)
Initialise class members.
GEnergy cutoff(void) const
Return exponential cut-off energy.
void copy_members(const GModelSpectralSuperExpPlaw &model)
Copy class members.
virtual void read(const GXmlElement &xml)
Read model from XML element.
Model parameter class interface definition.
virtual GModelSpectralSuperExpPlaw * clone(void) const
Clone super exponentially cut off power law model.
virtual ~GModelSpectralSuperExpPlaw(void)
Destructor.
Model parameter class.
Definition: GModelPar.hpp:87
GEnergy m_last_energy
Last energy value.
flux_kernel(const double &norm, const double &index1, const double &pivot, const double &ecut, const double &index2)
virtual std::string print(const GChatter &chatter=NORMAL) const
Print model information.
virtual void write(GXmlElement &xml) const
Write model into XML element.
double m_mc_pow_emin
Power of minimum energy.
Single parameter function abstract base class definition.
GModelPar m_norm
Normalization factor.
double m_last_pivot
Last pivot parameter.
double prefactor(void) const
Return pre factor.
double eval(const double &eng)
Kernel for photon flux integration.
GChatter
Definition: GTypemaps.hpp:33
GModelSpectralSuperExpPlaw(void)
Void constructor.
Abstract spectral model base class interface definition.
virtual std::string classname(void) const
Return class name.
double value(void) const
Return parameter value.
Single parameter function abstract base class.
Definition: GFunction.hpp:44
virtual GModelSpectralSuperExpPlaw & operator=(const GModelSpectralSuperExpPlaw &model)
Assignment operator.
virtual double eflux(const GEnergy &emin, const GEnergy &emax) const
Returns model energy flux between emin, emax
double m_last_e_cut
Last E/Ecut value.
double m_last_ecut
Last energy cut-off parameter.
Super exponential cut off power law spectral class.
virtual GEnergy mc(const GEnergy &emin, const GEnergy &emax, const GTime &time, GRan &ran) const
Returns MC energy between [emin, emax].
double m_mc_exponent
Exponent (index+1)
double index2(void) const
Return cut off index.
double m_last_index1
Last index1 parameter.
double index1(void) const
Return power law index.
eflux_kernel(const double &norm, const double &index1, const double &pivot, const double &ecut, const double &index2)
Class that handles energies in a unit independent way.
Definition: GEnergy.hpp:48