GammaLib  1.7.0.dev
 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-2018 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 //#include "GXmlElement.hpp"
37 
38 /* __ Forward declarations _______________________________________________ */
39 class GRan;
40 class GTime;
41 class GXmlElement;
42 
43 
44 /***********************************************************************//**
45  * @class GModelSpectralSuperExpPlaw
46  *
47  * @brief Super exponential cut off power law spectral class
48  *
49  * This class implements a power law spectrum with super exponential cut off.
50  * The model is defined by
51  *
52  * \f[
53  * S_{\rm E}(E | t) = {\tt m\_norm}
54  * \left( \frac{E}{\tt m\_pivot} \right)^{\tt m\_index1}
55  * \exp \left( \frac{-E}{\tt m\_ecut}^{\tt m\_index2} \right)
56  * \f]
57  *
58  * where
59  *
60  * \f${\tt m\_norm}\f$ is the normalization or prefactor,
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\_pivot}\f$ is the pivot energy.
64  * \f${\tt m\_index2}\f$ is index defining the cutoff shape,
65  ***************************************************************************/
67 
68 public:
69  // Constructors and destructors
71  GModelSpectralSuperExpPlaw(const std::string& type,
72  const std::string& prefactor,
73  const std::string& index1,
74  const std::string& pivot,
75  const std::string& cutoff,
76  const std::string& index2);
77  GModelSpectralSuperExpPlaw(const double& prefactor,
78  const double& index1,
79  const GEnergy& pivot,
80  const GEnergy& cutoff,
81  const double& index2);
82  explicit GModelSpectralSuperExpPlaw(const GXmlElement& xml);
84  virtual ~GModelSpectralSuperExpPlaw(void);
85 
86  // Operators
88 
89  // Implemented pure virtual methods
90  virtual void clear(void);
91  virtual GModelSpectralSuperExpPlaw* clone(void) const;
92  virtual std::string classname(void) const;
93  virtual std::string type(void) const;
94  virtual double eval(const GEnergy& srcEng,
95  const GTime& srcTime = GTime(),
96  const bool& gradients = false) const;
97  virtual double flux(const GEnergy& emin,
98  const GEnergy& emax) const;
99  virtual double eflux(const GEnergy& emin,
100  const GEnergy& emax) const;
101  virtual GEnergy mc(const GEnergy& emin,
102  const GEnergy& emax,
103  const GTime& time,
104  GRan& ran) const;
105  virtual void read(const GXmlElement& xml);
106  virtual void write(GXmlElement& xml) const;
107  virtual std::string print(const GChatter& chatter = NORMAL) const;
108 
109  // Other methods
110  void type(const std::string& type);
111  double prefactor(void) const;
112  void prefactor(const double& prefactor);
113  double index1(void) const;
114  void index1(const double& index1);
115  double index2(void) const;
116  void index2(const double& index2);
117  GEnergy cutoff(void) const;
118  void cutoff(const GEnergy& cutoff);
119  GEnergy pivot(void) const;
120  void pivot(const GEnergy& pivot);
121 
122 protected:
123  // Protected methods
124  void init_members(void);
125  void copy_members(const GModelSpectralSuperExpPlaw& model);
126  void free_members(void);
127  void update_eval_cache(const GEnergy& energy) const;
128  void update_mc_cache(const GEnergy& emin, const GEnergy& emax) const;
129 
130  // Photon flux integration kernel
131  class flux_kernel : public GFunction {
132  public:
133  flux_kernel(const double& norm,
134  const double& index1,
135  const double& pivot,
136  const double& ecut,
137  const double& index2) :
138  m_norm(norm),
139  m_index1(index1),
140  m_inv_pivot(1.0/pivot),
141  m_inv_ecut(1.0/ecut),
142  m_index2(index2) {}
143  double eval(const double& eng);
144  protected:
145  double m_norm; //!< Normalization
146  double m_index1; //!< Index1
147  double m_inv_pivot; //!< 1 / Pivot energy
148  double m_inv_ecut; //!< 1 / Cut off energy
149  double m_index2; //!< Index2
150  };
151 
152  // Energy flux integration kernel
153  class eflux_kernel : public GFunction {
154  public:
155  eflux_kernel(const double& norm,
156  const double& index1,
157  const double& pivot,
158  const double& ecut,
159  const double& index2) :
160  m_norm(norm),
161  m_index1(index1),
162  m_inv_pivot(1.0/pivot),
163  m_inv_ecut(1.0/ecut),
164  m_index2(index2) {}
165  double eval(const double& eng);
166  protected:
167  double m_norm; //!< Normalization
168  double m_index1; //!< Index1
169  double m_inv_pivot; //!< 1 / Pivot energy
170  double m_inv_ecut; //!< 1 / Cut off energy
171  double m_index2; //!< Index2
172  };
173 
174  // Protected members
175  std::string m_type; //!< Model type
176  GModelPar m_norm; //!< Normalization factor
177  GModelPar m_index1; //!< Spectral index
178  GModelPar m_index2; //!< Index of cutoff
179  GModelPar m_ecut; //!< Exponential cut off energy
180  GModelPar m_pivot; //!< Pivot energy
181 
182  // Cached members used for pre-computations
183  mutable GEnergy m_last_energy; //!< Last energy value
184  mutable double m_last_index1; //!< Last index1 parameter
185  mutable double m_last_ecut; //!< Last energy cut-off parameter
186  mutable double m_last_pivot; //!< Last pivot parameter
187  mutable double m_last_index2; //!< Last index2 parameter
188  mutable double m_last_e_norm; //!< Last E/Epivot value
189  mutable double m_last_e_cut; //!< Last E/Ecut value
190  mutable double m_last_exponent; //!< last pow(E/Ecut,index2) value
191  mutable double m_last_power; //!< Last power value
192  mutable double m_mc_emin; //!< Minimum energy
193  mutable double m_mc_emax; //!< Maximum energy
194  mutable double m_mc_exponent; //!< Exponent (index+1)
195  mutable double m_mc_pow_emin; //!< Power of minimum energy
196  mutable double m_mc_pow_ewidth; //!< Power of energy width
197 };
198 
199 
200 /***********************************************************************//**
201  * @brief Return class name
202  *
203  * @return String containing the class name ("GModelSpectralSuperExpPlaw").
204  ***************************************************************************/
205 inline
207 {
208  return ("GModelSpectralSuperExpPlaw");
209 }
210 
211 
212 /***********************************************************************//**
213  * @brief Return model type
214  *
215  * @return "SuperExpCutoff".
216  *
217  * Returns the type of the super exponentially cut off power law model.
218  ***************************************************************************/
219 inline
220 std::string GModelSpectralSuperExpPlaw::type(void) const
221 {
222  return (m_type);
223 }
224 
225 
226 /***********************************************************************//**
227  * @brief Set model type
228  *
229  * @param[in] type Model type.
230  *
231  * Set the type of the super exponentially cut off power law model.
232  ***************************************************************************/
233 inline
234 void GModelSpectralSuperExpPlaw::type(const std::string& type)
235 {
236  m_type = type;
237  return;
238 }
239 
240 
241 /***********************************************************************//**
242  * @brief Return pre factor
243  *
244  * @return Pre factor (ph/cm2/s/MeV).
245  *
246  * Returns the pre factor.
247  ***************************************************************************/
248 inline
250 {
251  return (m_norm.value());
252 }
253 
254 
255 /***********************************************************************//**
256  * @brief Set pre factor
257  *
258  * @param[in] prefactor Pre factor (ph/cm2/s/MeV).
259  *
260  * Sets the pre factor.
261  ***************************************************************************/
262 inline
263 void GModelSpectralSuperExpPlaw::prefactor(const double& prefactor)
264 {
265  m_norm.value(prefactor);
266  return;
267 }
268 
269 
270 /***********************************************************************//**
271  * @brief Return power law index
272  *
273  * @return Power law index.
274  *
275  * Returns the power law index.
276  ***************************************************************************/
277 inline
279 {
280  return (m_index1.value());
281 }
282 
283 
284 /***********************************************************************//**
285  * @brief Set power law index
286  *
287  * @param[in] index1 Power law index.
288  *
289  * Sets the power law index.
290  ***************************************************************************/
291 inline
292 void GModelSpectralSuperExpPlaw::index1(const double& index1)
293 {
294  m_index1.value(index1);
295  return;
296 }
297 
298 /***********************************************************************//**
299  * @brief Return cut off index
300  *
301  * @return cut off index.
302  *
303  * Returns the cut off index.
304  ***************************************************************************/
305 inline
307 {
308  return (m_index2.value());
309 }
310 
311 
312 /***********************************************************************//**
313  * @brief Set cut off index
314  *
315  * @param[in] index2 Cut off index.
316  *
317  * Sets the cut off index.
318  ***************************************************************************/
319 inline
320 void GModelSpectralSuperExpPlaw::index2(const double& index2)
321 {
322  m_index2.value(index2);
323  return;
324 }
325 
326 
327 /***********************************************************************//**
328  * @brief Return pivot energy
329  *
330  * @return Pivot energy.
331  *
332  * Returns the pivot energy.
333  ***************************************************************************/
334 inline
336 {
337  GEnergy energy;
338  energy.MeV(m_pivot.value());
339  return energy;
340 }
341 
342 
343 /***********************************************************************//**
344  * @brief Set pivot energy
345  *
346  * @param[in] pivot Pivot energy.
347  *
348  * Sets the pivot energy.
349  ***************************************************************************/
350 inline
352 {
353  m_pivot.value(pivot.MeV());
354  return;
355 }
356 
357 
358 /***********************************************************************//**
359  * @brief Return exponential cut-off energy
360  *
361  * @return Exponential cut-off energy.
362  *
363  * Returns the exponential cut-off energy.
364  ***************************************************************************/
365 inline
367 {
368  GEnergy energy;
369  energy.MeV(m_ecut.value());
370  return energy;
371 }
372 
373 
374 /***********************************************************************//**
375  * @brief Set exponential cut-off energy
376  *
377  * @param[in] cutoff Exponential cut-off energy.
378  *
379  * Sets the exponential cut-off energy.
380  ***************************************************************************/
381 inline
383 {
384  m_ecut.value(cutoff.MeV());
385  return;
386 }
387 
388 #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:821
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:47
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:54
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