Sergio Hernández Cadena / gammalib

Commit 5af896d0966cf6035fcbdce7619fe9d51a799a61

Authored by knodlseder
1 parent 190926e2

declare eval() and eval_gradients() const

Inline Side-by-side
Showing 17 changed files with 455 additions and 401 deletions
include/GModelSpectralConst.hpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralConst.hpp
15   - * @brief GModelSpectralConst class interface definition.
  15 + * @brief Constant spectral model class interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -30,15 +30,15 @@
30 30 /***********************************************************************//**
31 31 * @class GModelSpectralConst
32 32 *
33   - * @brief Spectral constant interface definition
  33 + * @brief Constant spectral model class
34 34 *
35 35 * This class implements a constant as the spectral component of the
36   - * gamma-ray data model. The function is defined as
  36 + * gamma-ray sky model. The function is defined as
37 37 * \f[I(E)=norm\f]
38 38 * where
39 39 * \f$norm\f$ is the normalization constant.
40 40 ***************************************************************************/
41   -class GModelSpectralConst : public GModelSpectral {
  41 +class GModelSpectralConst : public GModelSpectral {
42 42  
43 43 public:
44 44 // Constructors and destructors
... ... @@ -48,20 +48,19 @@ public:
48 48 virtual ~GModelSpectralConst(void);
49 49  
50 50 // Operators
51   - GModelSpectralConst& operator= (const GModelSpectralConst& model);
  51 + virtual GModelSpectralConst& operator=(const GModelSpectralConst& model);
52 52  
53 53 // Implemented pure virtual methods
54   - void clear(void);
55   - GModelSpectralConst* clone(void) const;
56   - int size(void) const { return m_npars; }
57   - std::string type(void) const { return "ConstantValue"; }
58   - double eval(const GEnergy& srcEng);
59   - double eval_gradients(const GEnergy& srcEng);
60   - double flux(const GEnergy& emin, const GEnergy& emax) const;
61   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
62   - void read(const GXmlElement& xml);
63   - void write(GXmlElement& xml) const;
64   - std::string print(void) const;
  54 + virtual void clear(void);
  55 + virtual GModelSpectralConst* clone(void) const;
  56 + virtual std::string type(void) const { return "ConstantValue"; }
  57 + virtual double eval(const GEnergy& srcEng) const;
  58 + virtual double eval_gradients(const GEnergy& srcEng) const;
  59 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  60 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  61 + virtual void read(const GXmlElement& xml);
  62 + virtual void write(GXmlElement& xml) const;
  63 + virtual std::string print(void) const;
65 64  
66 65 // Other methods
67 66 double norm(void) const { return m_norm.real_value(); }
... ... @@ -72,13 +71,8 @@ protected:
72 71 void copy_members(const GModelSpectralConst& model);
73 72 void free_members(void);
74 73  
75   - // Implemented pure virtual methods
76   - GModelPar** par(void) { return m_par; }
77   -
78 74 // Protected members
79   - int m_npars; //!< Number of parameters
80   - GModelPar* m_par[1]; //!< Pointers to parameters
81   - GModelPar m_norm; //!< Normalization factor
  75 + GModelPar m_norm; //!< Normalization factor
82 76 };
83 77  
84 78 #endif /* GMODELSPECTRALCONST_HPP */
... ...
include/GModelSpectralExpPlaw.hpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralExpPlaw.hpp
15   - * @brief GModelSpectralExpPlaw class interface definition.
  15 + * @brief Exponential cut off power law spectral class interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -20,6 +20,7 @@
20 20 #define GMODELSPECTRALEXPPLAW_HPP
21 21  
22 22 /* __ Includes ___________________________________________________________ */
  23 +#include <string>
23 24 #include "GModelPar.hpp"
24 25 #include "GModelSpectral.hpp"
25 26 #include "GEnergy.hpp"
... ... @@ -29,10 +30,10 @@
29 30 /***********************************************************************//**
30 31 * @class GModelSpectralExpPlaw
31 32 *
32   - * @brief Spectral power law interface definition.
  33 + * @brief Exponential cut off power law spectral class
33 34 *
34 35 * This class implements a power law as the spectral component of the
35   - * gamma-ray data model. The power law is defined as
  36 + * gamma-ray sky model. The power law is defined as
36 37 * \f[I(E)=norm (E/pivot)^{index} \exp(-E/ecut)\f]
37 38 * where
38 39 * \f$norm\f$ is the normalization or prefactor,
... ... @@ -40,7 +41,7 @@
40 41 * \f$index\f$ is the spectral index, and
41 42 * \f$ecut\f$ is the cut off energy.
42 43 ***************************************************************************/
43   -class GModelSpectralExpPlaw : public GModelSpectral {
  44 +class GModelSpectralExpPlaw : public GModelSpectral {
44 45  
45 46 public:
46 47 // Constructors and destructors
... ... @@ -51,20 +52,19 @@ public:
51 52 virtual ~GModelSpectralExpPlaw(void);
52 53  
53 54 // Operators
54   - GModelSpectralExpPlaw& operator= (const GModelSpectralExpPlaw& model);
  55 + virtual GModelSpectralExpPlaw& operator= (const GModelSpectralExpPlaw& model);
55 56  
56 57 // Implemented pure virtual methods
57   - void clear(void);
58   - GModelSpectralExpPlaw* clone(void) const;
59   - int size(void) const { return m_npars; }
60   - std::string type(void) const { return "ExpCutoff"; }
61   - double eval(const GEnergy& srcEng);
62   - double eval_gradients(const GEnergy& srcEng);
63   - double flux(const GEnergy& emin, const GEnergy& emax) const;
64   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
65   - void read(const GXmlElement& xml);
66   - void write(GXmlElement& xml) const;
67   - std::string print(void) const;
  58 + virtual void clear(void);
  59 + virtual GModelSpectralExpPlaw* clone(void) const;
  60 + virtual std::string type(void) const { return "ExpCutoff"; }
  61 + virtual double eval(const GEnergy& srcEng) const;
  62 + virtual double eval_gradients(const GEnergy& srcEng) const;
  63 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  64 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  65 + virtual void read(const GXmlElement& xml);
  66 + virtual void write(GXmlElement& xml) const;
  67 + virtual std::string print(void) const;
68 68  
69 69 // Other methods
70 70 void autoscale(void);
... ... @@ -79,16 +79,11 @@ protected:
79 79 void copy_members(const GModelSpectralExpPlaw& model);
80 80 void free_members(void);
81 81  
82   - // Implemented pure virtual methods
83   - GModelPar** par(void) { return m_par; }
84   -
85 82 // Protected members
86   - int m_npars; //!< Number of parameters
87   - GModelPar* m_par[4]; //!< Pointers to parameters
88   - GModelPar m_norm; //!< Normalization factor
89   - GModelPar m_index; //!< Spectral index
90   - GModelPar m_ecut; //!< Exponential cut off energy
91   - GModelPar m_pivot; //!< Pivot energy
  83 + GModelPar m_norm; //!< Normalization factor
  84 + GModelPar m_index; //!< Spectral index
  85 + GModelPar m_ecut; //!< Exponential cut off energy
  86 + GModelPar m_pivot; //!< Pivot energy
92 87 };
93 88  
94 89 #endif /* GMODELSPECTRALEXPPLAW_HPP */
... ...
include/GModelSpectralFunc.hpp
  View file @5af896d
... ... @@ -32,7 +32,7 @@
32 32 /***********************************************************************//**
33 33 * @class GModelSpectralFunc
34 34 *
35   - * @brief Spectral function interface definition
  35 + * @brief Spectral function model class
36 36 *
37 37 * This class implements an arbitrary function as the spectral component of
38 38 * the model. The function is defined as
... ... @@ -40,7 +40,7 @@
40 40 * where
41 41 * \f$norm\f$ is the normalization of the function.
42 42 ***************************************************************************/
43   -class GModelSpectralFunc : public GModelSpectral {
  43 +class GModelSpectralFunc : public GModelSpectral {
44 44  
45 45 public:
46 46 // Constructors and destructors
... ... @@ -51,20 +51,19 @@ public:
51 51 virtual ~GModelSpectralFunc(void);
52 52  
53 53 // Operators
54   - GModelSpectralFunc& operator= (const GModelSpectralFunc& model);
  54 + virtual GModelSpectralFunc& operator=(const GModelSpectralFunc& model);
55 55  
56 56 // Implemented pure virtual methods
57   - void clear(void);
58   - GModelSpectralFunc* clone(void) const;
59   - int size(void) const { return m_npars; }
60   - std::string type(void) const { return "FileFunction"; }
61   - double eval(const GEnergy& srcEng);
62   - double eval_gradients(const GEnergy& srcEng);
63   - double flux(const GEnergy& emin, const GEnergy& emax) const;
64   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
65   - void read(const GXmlElement& xml);
66   - void write(GXmlElement& xml) const;
67   - std::string print(void) const;
  57 + virtual void clear(void);
  58 + virtual GModelSpectralFunc* clone(void) const;
  59 + virtual std::string type(void) const { return "FileFunction"; }
  60 + virtual double eval(const GEnergy& srcEng) const;
  61 + virtual double eval_gradients(const GEnergy& srcEng) const;
  62 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  63 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  64 + virtual void read(const GXmlElement& xml);
  65 + virtual void write(GXmlElement& xml) const;
  66 + virtual std::string print(void) const;
68 67  
69 68 // Other methods
70 69 double norm(void) const { return m_norm.real_value(); }
... ... @@ -76,12 +75,7 @@ protected:
76 75 void free_members(void);
77 76 void load_nodes(const std::string& filename);
78 77  
79   - // Implemented pure virtual methods
80   - GModelPar** par(void) { return m_par; }
81   -
82 78 // Protected members
83   - int m_npars; //!< Number of parameters
84   - GModelPar* m_par[1]; //!< Pointers to parameters
85 79 GModelPar m_norm; //!< Normalization factor
86 80 GNodeArray m_nodes; //!< Nodes of function
87 81 std::vector<double> m_values; //!< Function values at nodes
... ...
include/GModelSpectralPlaw.hpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralPlaw.hpp
15   - * @brief GModelSpectralPlaw class interface definition.
  15 + * @brief Power law spectral model class interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -20,6 +20,7 @@
20 20 #define GMODELSPECTRALPLAW_HPP
21 21  
22 22 /* __ Includes ___________________________________________________________ */
  23 +#include <string>
23 24 #include "GModelPar.hpp"
24 25 #include "GModelSpectral.hpp"
25 26 #include "GEnergy.hpp"
... ... @@ -29,17 +30,17 @@
29 30 /***********************************************************************//**
30 31 * @class GModelSpectralPlaw
31 32 *
32   - * @brief Spectral power law interface definition.
  33 + * @brief Power law spectral model class
33 34 *
34 35 * This class implements a power law as the spectral component of the
35   - * gamma-ray data model. The power law is defined as
  36 + * gamma-ray sky model. The power law is defined as
36 37 * \f[I(E)=norm (E/pivot)^{index}\f]
37 38 * where
38 39 * \f$norm\f$ is the normalization or prefactor,
39 40 * \f$pivot\f$ is the pivot energy, and
40 41 * \f$index\f$ is the spectral index.
41 42 ***************************************************************************/
42   -class GModelSpectralPlaw : public GModelSpectral {
  43 +class GModelSpectralPlaw : public GModelSpectral {
43 44  
44 45 public:
45 46 // Constructors and destructors
... ... @@ -50,20 +51,19 @@ public:
50 51 virtual ~GModelSpectralPlaw(void);
51 52  
52 53 // Operators
53   - GModelSpectralPlaw& operator= (const GModelSpectralPlaw& model);
  54 + virtual GModelSpectralPlaw& operator=(const GModelSpectralPlaw& model);
54 55  
55 56 // Implemented pure virtual methods
56   - void clear(void);
57   - GModelSpectralPlaw* clone(void) const;
58   - int size(void) const { return m_npars; }
59   - std::string type(void) const { return "PowerLaw"; }
60   - double eval(const GEnergy& srcEng);
61   - double eval_gradients(const GEnergy& srcEng);
62   - double flux(const GEnergy& emin, const GEnergy& emax) const;
63   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
64   - void read(const GXmlElement& xml);
65   - void write(GXmlElement& xml) const;
66   - std::string print(void) const;
  57 + virtual void clear(void);
  58 + virtual GModelSpectralPlaw* clone(void) const;
  59 + virtual std::string type(void) const { return "PowerLaw"; }
  60 + virtual double eval(const GEnergy& srcEng) const;
  61 + virtual double eval_gradients(const GEnergy& srcEng) const;
  62 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  63 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  64 + virtual void read(const GXmlElement& xml);
  65 + virtual void write(GXmlElement& xml) const;
  66 + virtual std::string print(void) const;
67 67  
68 68 // Other methods
69 69 void autoscale(void);
... ... @@ -77,15 +77,10 @@ protected:
77 77 void copy_members(const GModelSpectralPlaw& model);
78 78 void free_members(void);
79 79  
80   - // Implemented pure virtual methods
81   - GModelPar** par(void) { return m_par; }
82   -
83 80 // Protected members
84   - int m_npars; //!< Number of parameters
85   - GModelPar* m_par[3]; //!< Pointers to parameters
86   - GModelPar m_norm; //!< Normalization factor
87   - GModelPar m_index; //!< Spectral index
88   - GModelPar m_pivot; //!< Pivot energy
  81 + GModelPar m_norm; //!< Normalization factor
  82 + GModelPar m_index; //!< Spectral index
  83 + GModelPar m_pivot; //!< Pivot energy
89 84 };
90 85  
91 86 #endif /* GMODELSPECTRALPLAW_HPP */
... ...
include/GModelSpectralPlaw2.hpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralPlaw2.hpp
15   - * @brief GModelSpectralPlaw2 class interface definition.
  15 + * @brief Flux normalized power law spectral model class interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -20,6 +20,7 @@
20 20 #define GMODELSPECTRALPLAW2_HPP
21 21  
22 22 /* __ Includes ___________________________________________________________ */
  23 +#include <string>
23 24 #include "GModelPar.hpp"
24 25 #include "GModelSpectral.hpp"
25 26 #include "GEnergy.hpp"
... ... @@ -29,10 +30,10 @@
29 30 /***********************************************************************//**
30 31 * @class GModelSpectralPlaw2
31 32 *
32   - * @brief Spectral power law interface definition.
  33 + * @brief Flux normalized power law spectral model class
33 34 *
34 35 * This class implements a power law as the spectral component of the
35   - * gamma-ray data model.
  36 + * gamma-ray sky model.
36 37 * The power law function is defined as
37 38 * \f[I(E)=integral (index+1)/(emax^{index+1}-emin^{index+1}) E^{index}\f]
38 39 * for \f$index \ne -1\f$ and
... ... @@ -43,7 +44,7 @@
43 44 * \f$emin\f$ is the lower energy limit, and
44 45 * \f$emax\f$ is the upper energy limit.
45 46 ***************************************************************************/
46   -class GModelSpectralPlaw2 : public GModelSpectral {
  47 +class GModelSpectralPlaw2 : public GModelSpectral {
47 48  
48 49 public:
49 50 // Constructors and destructors
... ... @@ -54,20 +55,19 @@ public:
54 55 virtual ~GModelSpectralPlaw2(void);
55 56  
56 57 // Operators
57   - GModelSpectralPlaw2& operator= (const GModelSpectralPlaw2& model);
  58 + virtual GModelSpectralPlaw2& operator= (const GModelSpectralPlaw2& model);
58 59  
59 60 // Implemented pure virtual methods
60   - void clear(void);
61   - GModelSpectralPlaw2* clone(void) const;
62   - int size(void) const { return m_npars; }
63   - std::string type(void) const { return "PowerLaw2"; }
64   - double eval(const GEnergy& srcEng);
65   - double eval_gradients(const GEnergy& srcEng);
66   - double flux(const GEnergy& emin, const GEnergy& emax) const;
67   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
68   - void read(const GXmlElement& xml);
69   - void write(GXmlElement& xml) const;
70   - std::string print(void) const;
  61 + virtual void clear(void);
  62 + virtual GModelSpectralPlaw2* clone(void) const;
  63 + virtual std::string type(void) const { return "PowerLaw2"; }
  64 + virtual double eval(const GEnergy& srcEng) const;
  65 + virtual double eval_gradients(const GEnergy& srcEng) const;
  66 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  67 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  68 + virtual void read(const GXmlElement& xml);
  69 + virtual void write(GXmlElement& xml) const;
  70 + virtual std::string print(void) const;
71 71  
72 72 // Other methods
73 73 double integral(void) const { return m_integral.real_value(); }
... ... @@ -82,31 +82,26 @@ protected:
82 82 void free_members(void);
83 83 void update(const GEnergy& srcEng);
84 84  
85   - // Implemented pure virtual methods
86   - GModelPar** par(void) { return m_par; }
87   -
88 85 // Protected members
89   - int m_npars; //!< Number of parameters
90   - GModelPar* m_par[4]; //!< Pointers to parameters
91   - GModelPar m_integral; //!< Integral flux
92   - GModelPar m_index; //!< Spectral index
93   - GModelPar m_emin; //!< Lower energy limit (MeV)
94   - GModelPar m_emax; //!< Upper energy limit (MeV)
95   - double m_log_emin; //!< Log(emin)
96   - double m_log_emax; //!< Log(emax)
97   - double m_pow_emin; //!< emin^(index+1)
98   - double m_pow_emax; //!< emax^(index+1)
99   - double m_norm; //!< Power-law normalization (for pivot energy 1 MeV)
100   - double m_g_norm; //!< Power-law normalization gradient
101   - double m_power; //!< Power-law factor
102   - double m_last_integral; //!< Last integral flux
103   - double m_last_index; //!< Last spectral index (MeV)
104   - double m_last_emin; //!< Last lower energy limit (MeV)
105   - double m_last_emax; //!< Last upper energy limit (MeV)
106   - GEnergy m_last_energy; //!< Last source energy
107   - double m_last_value; //!< Last function value
108   - double m_last_g_integral; //!< Last integral flux gradient
109   - double m_last_g_index; //!< Last spectral index gradient
  86 + GModelPar m_integral; //!< Integral flux
  87 + GModelPar m_index; //!< Spectral index
  88 + GModelPar m_emin; //!< Lower energy limit (MeV)
  89 + GModelPar m_emax; //!< Upper energy limit (MeV)
  90 + double m_log_emin; //!< Log(emin)
  91 + double m_log_emax; //!< Log(emax)
  92 + double m_pow_emin; //!< emin^(index+1)
  93 + double m_pow_emax; //!< emax^(index+1)
  94 + double m_norm; //!< Power-law normalization (for pivot energy 1 MeV)
  95 + double m_g_norm; //!< Power-law normalization gradient
  96 + double m_power; //!< Power-law factor
  97 + double m_last_integral; //!< Last integral flux
  98 + double m_last_index; //!< Last spectral index (MeV)
  99 + double m_last_emin; //!< Last lower energy limit (MeV)
  100 + double m_last_emax; //!< Last upper energy limit (MeV)
  101 + GEnergy m_last_energy; //!< Last source energy
  102 + double m_last_value; //!< Last function value
  103 + double m_last_g_integral; //!< Last integral flux gradient
  104 + double m_last_g_index; //!< Last spectral index gradient
110 105 };
111 106  
112 107 #endif /* GMODELSPECTRALPLAW2_HPP */
... ...
include/GModelTemporalConst.hpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelTemporalConst.hpp
15   - * @brief GModelTemporalConst class interface definition.
  15 + * @brief Constant temporal model class interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -20,22 +20,24 @@
20 20 #define GMODELTEMPORALCONST_HPP
21 21  
22 22 /* __ Includes ___________________________________________________________ */
  23 +#include <string>
23 24 #include "GModelPar.hpp"
24 25 #include "GModelTemporal.hpp"
25 26 #include "GTime.hpp"
  27 +#include "GXmlElement.hpp"
26 28  
27 29  
28 30 /***********************************************************************//**
29 31 * @class GModelTemporalConst
30 32 *
31   - * @brief Interface definition for the constant model class
  33 + * @brief Constant temporal model class
32 34 *
33 35 * This class implements the temporal component of the factorised model for
34 36 * a model that is constant in time. It has a single parameter norm that is
35 37 * used for scaling, and that normally is set to 1. The parameter is not
36 38 * supposed to be fitted.
37 39 ***************************************************************************/
38   -class GModelTemporalConst : public GModelTemporal {
  40 +class GModelTemporalConst : public GModelTemporal {
39 41  
40 42 public:
41 43 // Constructors and destructors
... ... @@ -44,21 +46,19 @@ public:
44 46 virtual ~GModelTemporalConst(void);
45 47  
46 48 // Operators
47   - GModelTemporalConst& operator= (const GModelTemporalConst& model);
  49 + virtual GModelTemporalConst& operator= (const GModelTemporalConst& model);
48 50  
49 51 // Implemented virtual methods
50   - void clear(void);
51   - GModelTemporalConst* clone(void) const;
52   - int size(void) const { return m_npars; }
53   - std::string type(void) const { return "Constant"; }
54   - double eval(const GTime& srcTime);
55   - double eval_gradients(const GTime& srcTime);
56   - GTimes mc(const double& rate,
57   - const GTime& tmin, const GTime& tmax,
58   - GRan& ran);
59   - void read(const GXmlElement& xml);
60   - void write(GXmlElement& xml) const;
61   - std::string print(void) const;
  52 + virtual void clear(void);
  53 + virtual GModelTemporalConst* clone(void) const;
  54 + virtual std::string type(void) const { return "Constant"; }
  55 + virtual double eval(const GTime& srcTime) const;
  56 + virtual double eval_gradients(const GTime& srcTime) const;
  57 + virtual GTimes mc(const double& rate, const GTime& tmin,
  58 + const GTime& tmax, GRan& ran) const;
  59 + virtual void read(const GXmlElement& xml);
  60 + virtual void write(GXmlElement& xml) const;
  61 + virtual std::string print(void) const;
62 62  
63 63 // Other methods
64 64 double norm(void) const { return m_norm.real_value(); }
... ... @@ -69,13 +69,8 @@ protected:
69 69 void copy_members(const GModelTemporalConst& model);
70 70 void free_members(void);
71 71  
72   - // Implemented pure virtual methods
73   - GModelPar** par(void) { return m_par; }
74   -
75 72 // Protected members
76   - int m_npars; //!< Number of parameters
77   - GModelPar* m_par[1]; //!< Pointers to parameters
78   - GModelPar m_norm; //!< Constant
  73 + GModelPar m_norm; //!< Constant
79 74 };
80 75  
81 76 #endif /* GMODELTEMPORALCONST_HPP */
... ...
pyext/GModelSpectralExpPlaw.i
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralExpPlaw.i
15   - * @brief GModelSpectralExpPlaw class python interface
  15 + * @brief Exponential cut off power law spectral class Python interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18 %{
... ... @@ -25,9 +25,9 @@
25 25 /***********************************************************************//**
26 26 * @class GModelSpectralExpPlaw
27 27 *
28   - * @brief Implements an exponentially cut off power law.
  28 + * @brief Exponential cut off power law spectral class
29 29 ***************************************************************************/
30   -class GModelSpectralExpPlaw : public GModelSpectral {
  30 +class GModelSpectralExpPlaw : public GModelSpectral {
31 31 public:
32 32 // Constructors and destructors
33 33 GModelSpectralExpPlaw(void);
... ... @@ -38,23 +38,22 @@ public:
38 38 virtual ~GModelSpectralExpPlaw(void);
39 39  
40 40 // Implemented pure virtual methods
41   - void clear(void);
42   - GModelSpectralExpPlaw* clone(void) const;
43   - int size(void) const { return m_npars; }
44   - std::string type(void) const { return "ExpCutoff"; }
45   - double eval(const GEnergy& srcEng);
46   - double eval_gradients(const GEnergy& srcEng);
47   - double flux(const GEnergy& emin, const GEnergy& emax) const;
48   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
49   - void read(const GXmlElement& xml);
50   - void write(GXmlElement& xml) const;
  41 + virtual void clear(void);
  42 + virtual GModelSpectralExpPlaw* clone(void) const;
  43 + virtual std::string type(void) const;
  44 + virtual double eval(const GEnergy& srcEng) const;
  45 + virtual double eval_gradients(const GEnergy& srcEng) const;
  46 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  47 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  48 + virtual void read(const GXmlElement& xml);
  49 + virtual void write(GXmlElement& xml) const;
51 50  
52 51 // Other methods
53 52 void autoscale(void);
54   - double norm(void) const { return m_norm.real_value(); }
55   - double index(void) const { return m_index.real_value(); }
56   - double ecut(void) const { return m_ecut.real_value(); }
57   - double pivot(void) const { return m_pivot.real_value(); }
  53 + double norm(void) const;
  54 + double index(void) const;
  55 + double ecut(void) const;
  56 + double pivot(void) const;
58 57 };
59 58  
60 59  
... ...
pyext/GModelSpectralPlaw.i
  View file @5af896d
1 1 /***************************************************************************
2   - * GModelSpectralPlaw.i - Spectral power law model class python I/F *
  2 + * GModelSpectralPlaw.i - Spectral power law model class *
3 3 * ----------------------------------------------------------------------- *
4 4 * copyright (C) 2009-2011 by Jurgen Knodlseder *
5 5 * ----------------------------------------------------------------------- *
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralPlaw.i
15   - * @brief GModelSpectralPlaw class python interface
  15 + * @brief Power law spectral model class Python interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18 %{
... ... @@ -25,9 +25,9 @@
25 25 /***********************************************************************//**
26 26 * @class GModelSpectralPlaw
27 27 *
28   - * @brief Powerlaw python interface definition
  28 + * @brief Power law spectral model class
29 29 ***************************************************************************/
30   -class GModelSpectralPlaw : public GModelSpectral {
  30 +class GModelSpectralPlaw : public GModelSpectral {
31 31 public:
32 32 // Constructors and destructors
33 33 GModelSpectralPlaw(void);
... ... @@ -37,22 +37,21 @@ public:
37 37 virtual ~GModelSpectralPlaw(void);
38 38  
39 39 // Implemented pure virtual methods
40   - void clear(void);
41   - GModelSpectralPlaw* clone(void) const;
42   - int size(void) const { return m_npars; }
43   - std::string type(void) const { return "PowerLaw"; }
44   - double eval(const GEnergy& srcEng);
45   - double eval_gradients(const GEnergy& srcEng);
46   - double flux(const GEnergy& emin, const GEnergy& emax) const;
47   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
48   - void read(const GXmlElement& xml);
49   - void write(GXmlElement& xml) const;
  40 + virtual void clear(void);
  41 + virtual GModelSpectralPlaw* clone(void) const;
  42 + virtual std::string type(void) const;
  43 + virtual double eval(const GEnergy& srcEng) const;
  44 + virtual double eval_gradients(const GEnergy& srcEng) const;
  45 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  46 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  47 + virtual void read(const GXmlElement& xml);
  48 + virtual void write(GXmlElement& xml) const;
50 49  
51 50 // Other methods
52 51 void autoscale(void);
53   - double norm(void) const { return m_norm.real_value(); }
54   - double index(void) const { return m_index.real_value(); }
55   - double pivot(void) const { return m_pivot.real_value(); }
  52 + double norm(void) const;
  53 + double index(void) const;
  54 + double pivot(void) const;
56 55 };
57 56  
58 57  
... ...
pyext/GModelSpectralPlaw2.i
  View file @5af896d
1 1 /***************************************************************************
2   - * GModelSpectralPlaw2.i - Spectral power law model class python I/F *
  2 + * GModelSpectralPlaw2.i - Spectral power law model class *
3 3 * ----------------------------------------------------------------------- *
4 4 * copyright (C) 2011 by Jurgen Knodlseder *
5 5 * ----------------------------------------------------------------------- *
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralPlaw2.i
15   - * @brief GModelSpectralPlaw2 class python interface
  15 + * @brief Flux normalized power law spectral model class Python interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18 %{
... ... @@ -25,9 +25,9 @@
25 25 /***********************************************************************//**
26 26 * @class GModelSpectralPlaw2
27 27 *
28   - * @brief Powerlaw python interface definition
  28 + * @brief Flux normalized power law spectral model class
29 29 ***************************************************************************/
30   -class GModelSpectralPlaw2 : public GModelSpectral {
  30 +class GModelSpectralPlaw2 : public GModelSpectral {
31 31 public:
32 32 // Constructors and destructors
33 33 GModelSpectralPlaw2(void);
... ... @@ -37,16 +37,15 @@ public:
37 37 virtual ~GModelSpectralPlaw2(void);
38 38  
39 39 // Implemented pure virtual methods
40   - void clear(void);
41   - GModelSpectralPlaw2* clone(void) const;
42   - int size(void) const;
43   - std::string type(void) const;
44   - double eval(const GEnergy& srcEng);
45   - double eval_gradients(const GEnergy& srcEng);
46   - double flux(const GEnergy& emin, const GEnergy& emax) const;
47   - GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
48   - void read(const GXmlElement& xml);
49   - void write(GXmlElement& xml) const;
  40 + virtual void clear(void);
  41 + virtual GModelSpectralPlaw2* clone(void) const;
  42 + virtual std::string type(void) const;
  43 + virtual double eval(const GEnergy& srcEng) const;
  44 + virtual double eval_gradients(const GEnergy& srcEng) const;
  45 + virtual double flux(const GEnergy& emin, const GEnergy& emax) const;
  46 + virtual GEnergy mc(const GEnergy& emin, const GEnergy& emax, GRan& ran) const;
  47 + virtual void read(const GXmlElement& xml);
  48 + virtual void write(GXmlElement& xml) const;
50 49  
51 50 // Other methods
52 51 double integral(void) const;
... ...
pyext/GModelTemporalConst.i
  View file @5af896d
1 1 /***************************************************************************
2   - * GModelTemporalConst.i - Temporal constant model class python I/F *
  2 + * GModelTemporalConst.i - Temporal constant model class *
3 3 * ----------------------------------------------------------------------- *
4 4 * copyright (C) 2009-2011 by Jurgen Knodlseder *
5 5 * ----------------------------------------------------------------------- *
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelTemporalConst.i
15   - * @brief GModelTemporalConst class python interface
  15 + * @brief Constant temporal model class Python interface definition
16 16 * @author J. Knodlseder
17 17 */
18 18 %{
... ... @@ -25,7 +25,7 @@
25 25 /***********************************************************************//**
26 26 * @class GModelTemporalConst
27 27 *
28   - * @brief Python interface definition for the constant model class.
  28 + * @brief Constant temporal model class
29 29 ***************************************************************************/
30 30 class GModelTemporalConst : public GModelTemporal {
31 31 public:
... ... @@ -35,20 +35,18 @@ public:
35 35 virtual ~GModelTemporalConst(void);
36 36  
37 37 // Implemented virtual methods
38   - void clear(void);
39   - GModelTemporalConst* clone(void) const;
40   - int size(void) const { return m_npars; }
41   - std::string type(void) const { return "Constant"; }
42   - double eval(const GTime& srcTime);
43   - double eval_gradients(const GTime& srcTime);
44   - GTimes mc(const double& rate,
45   - const GTime& tmin, const GTime& tmax,
46   - GRan& ran);
47   - void read(const GXmlElement& xml);
48   - void write(GXmlElement& xml) const;
  38 + virtual void clear(void);
  39 + virtual GModelTemporalConst* clone(void) const;
  40 + virtual std::string type(void) const;
  41 + virtual double eval(const GTime& srcTime) const;
  42 + virtual double eval_gradients(const GTime& srcTime) const;
  43 + virtual GTimes mc(const double& rate, const GTime& tmin,
  44 + const GTime& tmax, GRan& ran) const;
  45 + virtual void read(const GXmlElement& xml);
  46 + virtual void write(GXmlElement& xml) const;
49 47  
50 48 // Other methods
51   - double norm(void) const { return m_norm.real_value(); }
  49 + double norm(void) const;
52 50 };
53 51  
54 52  
... ...
src/model/GModelSpatialCube.cpp
  View file @5af896d
... ... @@ -435,9 +435,6 @@ void GModelSpatialCube::free_members(void)
435 435 *
436 436 * @param[in] filename File name.
437 437 *
438   - * @exception GException::feature_not_implemented
439   - * Map cube loading not yet implemented.
440   - *
441 438 * @todo Implement method.
442 439 ***************************************************************************/
443 440 void GModelSpatialCube::load_cube(const std::string& filename)
... ... @@ -445,9 +442,6 @@ void GModelSpatialCube::load_cube(const std::string&amp; filename)
445 442 // Set filename
446 443 m_filename = filename;
447 444  
448   - // Dump warning that method is not yet implemented
449   - throw GException::feature_not_implemented(G_LOAD_CUBE);
450   -
451 445 // Return
452 446 return;
453 447 }
... ...
src/model/GModelSpectralConst.cpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralConst.cpp
15   - * @brief GModelSpectralConst class implementation.
  15 + * @brief Constant spectral model class implementation
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -54,7 +54,7 @@ const GModelSpectralRegistry g_spectral_const_registry(&amp;g_spectral_const_seed);
54 54 ***************************************************************************/
55 55 GModelSpectralConst::GModelSpectralConst(void) : GModelSpectral()
56 56 {
57   - // Initialise private members for clean destruction
  57 + // Initialise members
58 58 init_members();
59 59  
60 60 // Return
... ... @@ -63,16 +63,18 @@ GModelSpectralConst::GModelSpectralConst(void) : GModelSpectral()
63 63  
64 64  
65 65 /***********************************************************************//**
66   - * @brief Constructor
  66 + * @brief XML constructor
67 67 *
68 68 * @param[in] xml XML element.
69 69 *
70   - * Construct an arbitrary spectral function from a XML element.
  70 + * Creates instance of a constant spectral model by extracting information
  71 + * from an XML element. See GModelSpectralConst::read() for more information
  72 + * about the expected structure of the XML element.
71 73 ***************************************************************************/
72 74 GModelSpectralConst::GModelSpectralConst(const GXmlElement& xml) :
73 75 GModelSpectral()
74 76 {
75   - // Initialise private members for clean destruction
  77 + // Initialise members
76 78 init_members();
77 79  
78 80 // Read information from XML element
... ... @@ -86,12 +88,12 @@ GModelSpectralConst::GModelSpectralConst(const GXmlElement&amp; xml) :
86 88 /***********************************************************************//**
87 89 * @brief Copy constructor
88 90 *
89   - * @param[in] model Spectral component.
  91 + * @param[in] model Spectral constant model.
90 92 ***************************************************************************/
91 93 GModelSpectralConst::GModelSpectralConst(const GModelSpectralConst& model) :
92 94 GModelSpectral(model)
93 95 {
94   - // Initialise private members for clean destruction
  96 + // Initialise members
95 97 init_members();
96 98  
97 99 // Copy members
... ... @@ -124,7 +126,7 @@ GModelSpectralConst::~GModelSpectralConst(void)
124 126 /***********************************************************************//**
125 127 * @brief Assignment operator
126 128 *
127   - * @param[in] model Spectral component.
  129 + * @param[in] model Spectral constant model.
128 130 ***************************************************************************/
129 131 GModelSpectralConst& GModelSpectralConst::operator= (const GModelSpectralConst& model)
130 132 {
... ... @@ -137,7 +139,7 @@ GModelSpectralConst&amp; GModelSpectralConst::operator= (const GModelSpectralConst&amp;
137 139 // Free members
138 140 free_members();
139 141  
140   - // Initialise private members for clean destruction
  142 + // Initialise members
141 143 init_members();
142 144  
143 145 // Copy members
... ... @@ -193,7 +195,7 @@ GModelSpectralConst* GModelSpectralConst::clone(void) const
193 195 * where
194 196 * \f$norm\f$ is the normalization of the function.
195 197 ***************************************************************************/
196   -double GModelSpectralConst::eval(const GEnergy& srcEng)
  198 +double GModelSpectralConst::eval(const GEnergy& srcEng) const
197 199 {
198 200 // Compute function value
199 201 double value = norm();
... ... @@ -219,16 +221,16 @@ double GModelSpectralConst::eval(const GEnergy&amp; srcEng)
219 221 * The partial derivative of the normalization value is given by
220 222 * \f[dI/dn_v=n_s\f]
221 223 ***************************************************************************/
222   -double GModelSpectralConst::eval_gradients(const GEnergy& srcEng)
  224 +double GModelSpectralConst::eval_gradients(const GEnergy& srcEng) const
223 225 {
224 226 // Compute function value
225 227 double value = norm();
226 228  
227 229 // Compute partial derivatives of the parameter values
228   - double g_norm = (m_norm.isfree()) ? m_norm.scale() : 0.0;
  230 + double g_norm = (m_norm.isfree()) ? m_norm.scale() : 0.0;
229 231  
230   - // Set gradients
231   - m_norm.gradient(g_norm);
  232 + // Set gradients (circumvent const correctness)
  233 + ((GModelSpectralConst*)this)->m_norm.gradient(g_norm);
232 234  
233 235 // Return
234 236 return value;
... ... @@ -387,7 +389,7 @@ std::string GModelSpectralConst::print(void) const
387 389 result.append("=== GModelSpectralConst ===\n");
388 390 result.append(parformat("Number of parameters")+str(size()));
389 391 for (int i = 0; i < size(); ++i)
390   - result.append("\n"+m_par[i]->print());
  392 + result.append("\n"+m_pars[i]->print());
391 393  
392 394 // Return result
393 395 return result;
... ... @@ -405,17 +407,19 @@ std::string GModelSpectralConst::print(void) const
405 407 ***************************************************************************/
406 408 void GModelSpectralConst::init_members(void)
407 409 {
408   - // Initialise parameters
409   - m_npars = 1;
410   - m_par[0] = &m_norm;
411   -
412 410 // Initialise powerlaw normalisation
413   - m_norm = GModelPar();
  411 + m_norm.clear();
414 412 m_norm.name("Value");
415 413 m_norm.value(1.0);
416 414 m_norm.scale(1.0);
417 415 m_norm.range(0.0, 1000.0);
418 416 m_norm.free();
  417 + m_norm.gradient(0.0);
  418 + m_norm.hasgrad(true);
  419 +
  420 + // Set parameter pointer(s)
  421 + m_pars.clear();
  422 + m_pars.push_back(&m_norm);
419 423  
420 424 // Return
421 425 return;
... ... @@ -425,13 +429,17 @@ void GModelSpectralConst::init_members(void)
425 429 /***********************************************************************//**
426 430 * @brief Copy class members
427 431 *
428   - * @param[in] model Spectral model component.
  432 + * @param[in] model Spectral constant model.
429 433 ***************************************************************************/
430 434 void GModelSpectralConst::copy_members(const GModelSpectralConst& model)
431 435 {
432   - // Copy model parameters
  436 + // Copy members
433 437 m_norm = model.m_norm;
434 438  
  439 + // Set parameter pointer(s)
  440 + m_pars.clear();
  441 + m_pars.push_back(&m_norm);
  442 +
435 443 // Return
436 444 return;
437 445 }
... ...
src/model/GModelSpectralExpPlaw.cpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralExpPlaw.cpp
15   - * @brief GModelSpectralExpPlaw class implementation.
  15 + * @brief Exponential cut off power law spectral class implementation
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -56,7 +56,7 @@ const GModelSpectralRegistry g_spectral_expplaw_registry(&amp;g_spectral_expplaw_see
56 56 ***************************************************************************/
57 57 GModelSpectralExpPlaw::GModelSpectralExpPlaw(void) : GModelSpectral()
58 58 {
59   - // Initialise private members for clean destruction
  59 + // Initialise members
60 60 init_members();
61 61  
62 62 // Return
... ... @@ -69,7 +69,7 @@ GModelSpectralExpPlaw::GModelSpectralExpPlaw(void) : GModelSpectral()
69 69 *
70 70 * @param[in] norm Power law normalization.
71 71 * @param[in] index Power law index.
72   - * @param[in] ecut Cut off energy.
  72 + * @param[in] ecut Cut off energy (in MeV).
73 73 *
74 74 * Construct an exponential cut off power law from a normalization value,
75 75 * a spectral index and a cut off energy.
... ... @@ -78,7 +78,7 @@ GModelSpectralExpPlaw::GModelSpectralExpPlaw(const double&amp; norm,
78 78 const double& index,
79 79 const double& ecut) : GModelSpectral()
80 80 {
81   - // Initialise private members for clean destruction
  81 + // Initialise members
82 82 init_members();
83 83  
84 84 // Set parameters
... ... @@ -92,15 +92,18 @@ GModelSpectralExpPlaw::GModelSpectralExpPlaw(const double&amp; norm,
92 92  
93 93  
94 94 /***********************************************************************//**
95   - * @brief Constructor
  95 + * @brief XML constructor
96 96 *
97   - * @param[in] xml XML element containing position information.
  97 + * @param[in] xml XML element.
98 98 *
99   - * Construct a spectral power law from a XML element.
  99 + * Creates instance of an exponential cut off power law spectral model by
  100 + * extracting information from an XML element.
  101 + * See GModelSpectralExpPlaw::read() for more information about the expected
  102 + * structure of the XML element.
100 103 ***************************************************************************/
101 104 GModelSpectralExpPlaw::GModelSpectralExpPlaw(const GXmlElement& xml) : GModelSpectral()
102 105 {
103   - // Initialise private members for clean destruction
  106 + // Initialise members
104 107 init_members();
105 108  
106 109 // Read information from XML element
... ... @@ -114,12 +117,12 @@ GModelSpectralExpPlaw::GModelSpectralExpPlaw(const GXmlElement&amp; xml) : GModelSpe
114 117 /***********************************************************************//**
115 118 * @brief Copy constructor
116 119 *
117   - * @param[in] model Model from which the instance should be built.
  120 + * @param[in] model Exponential cut off power law model.
118 121 ***************************************************************************/
119 122 GModelSpectralExpPlaw::GModelSpectralExpPlaw(const GModelSpectralExpPlaw& model) :
120 123 GModelSpectral(model)
121 124 {
122   - // Initialise private members for clean destruction
  125 + // Initialise members
123 126 init_members();
124 127  
125 128 // Copy members
... ... @@ -152,7 +155,7 @@ GModelSpectralExpPlaw::~GModelSpectralExpPlaw(void)
152 155 /***********************************************************************//**
153 156 * @brief Assignment operator
154 157 *
155   - * @param[in] model Model which should be assigned.
  158 + * @param[in] model Exponential cut off power law model.
156 159 ***************************************************************************/
157 160 GModelSpectralExpPlaw& GModelSpectralExpPlaw::operator= (const GModelSpectralExpPlaw& model)
158 161 {
... ... @@ -165,7 +168,7 @@ GModelSpectralExpPlaw&amp; GModelSpectralExpPlaw::operator= (const GModelSpectralExp
165 168 // Free members
166 169 free_members();
167 170  
168   - // Initialise private members for clean destruction
  171 + // Initialise members
169 172 init_members();
170 173  
171 174 // Copy members
... ... @@ -229,7 +232,7 @@ GModelSpectralExpPlaw* GModelSpectralExpPlaw::clone(void) const
229 232 * in the futur. Furthermore, the method expects that pivot!=0 and ecut!=0.
230 233 * Otherwise Inf or NaN may result.
231 234 ***************************************************************************/
232   -double GModelSpectralExpPlaw::eval(const GEnergy& srcEng)
  235 +double GModelSpectralExpPlaw::eval(const GEnergy& srcEng) const
233 236 {
234 237 // Compute function value
235 238 double energy = srcEng.MeV();
... ... @@ -266,11 +269,11 @@ double GModelSpectralExpPlaw::eval(const GEnergy&amp; srcEng)
266 269 * \f[dI/dp_v=norm (E/pivot)^{index} \exp(-E/ecut) (-index) / p_v\f]
267 270 *
268 271 * @todo For the moment the pivot and the cut off energies are fixed to
269   - * units of MeV. This may not be ideal and should eventually be improved
270   - * in the futur. Furthermore, the method expects that pivot!=0 and ecut!=0.
271   - * Otherwise Inf or NaN may result.
  272 + * units of MeV. This may not be ideal and should eventually be
  273 + * improved in the futur. Furthermore, the method expects that
  274 + * pivot!=0 and ecut!=0. Otherwise Inf or NaN may result.
272 275 ***************************************************************************/
273   -double GModelSpectralExpPlaw::eval_gradients(const GEnergy& srcEng)
  276 +double GModelSpectralExpPlaw::eval_gradients(const GEnergy& srcEng) const
274 277 {
275 278 // Compute function value
276 279 double energy = srcEng.MeV();
... ... @@ -285,11 +288,11 @@ double GModelSpectralExpPlaw::eval_gradients(const GEnergy&amp; srcEng)
285 288 double g_ecut = (m_ecut.isfree()) ? value * e_cut/m_ecut.value() : 0.0;
286 289 double g_pivot = (m_pivot.isfree()) ? -value * index() / m_pivot.value() : 0.0;
287 290  
288   - // Set gradients
289   - m_norm.gradient(g_norm);
290   - m_index.gradient(g_index);
291   - m_ecut.gradient(g_ecut);
292   - m_pivot.gradient(g_pivot);
  291 + // Set gradients (circumvent const correctness)
  292 + ((GModelSpectralExpPlaw*)this)->m_norm.gradient(g_norm);
  293 + ((GModelSpectralExpPlaw*)this)->m_index.gradient(g_index);
  294 + ((GModelSpectralExpPlaw*)this)->m_ecut.gradient(g_ecut);
  295 + ((GModelSpectralExpPlaw*)this)->m_pivot.gradient(g_pivot);
293 296  
294 297 // Return
295 298 return value;
... ... @@ -302,6 +305,9 @@ double GModelSpectralExpPlaw::eval_gradients(const GEnergy&amp; srcEng)
302 305 * @param[in] emin Minimum photon energy.
303 306 * @param[in] emax Maximum photon energy.
304 307 *
  308 + * @exception GException::feature_not_implemented
  309 + * Method not yet implemented.
  310 + *
305 311 * Computes
306 312 * \f[\int_{E_{\rm min}}^{E_{\rm max}} I(E) dE\f]
307 313 * where
... ... @@ -309,7 +315,7 @@ double GModelSpectralExpPlaw::eval_gradients(const GEnergy&amp; srcEng)
309 315 * energy, respectively, and
310 316 * \f$I(E)\f$ is the spectral model (units: ph/cm2/s/MeV).
311 317 *
312   - * @todo Implement method
  318 + * @todo Implement method.
313 319 ***************************************************************************/
314 320 double GModelSpectralExpPlaw::flux(const GEnergy& emin, const GEnergy& emax) const
315 321 {
... ... @@ -328,7 +334,10 @@ double GModelSpectralExpPlaw::flux(const GEnergy&amp; emin, const GEnergy&amp; emax) con
328 334 * @param[in] emax Maximum photon energy.
329 335 * @param[in] ran Random number generator.
330 336 *
331   - * @todo To be implemented
  337 + * @exception GException::feature_not_implemented
  338 + * Method not yet implemented.
  339 + *
  340 + * @todo Implement method.
332 341 ***************************************************************************/
333 342 GEnergy GModelSpectralExpPlaw::mc(const GEnergy& emin, const GEnergy& emax,
334 343 GRan& ran) const
... ... @@ -347,8 +356,9 @@ GEnergy GModelSpectralExpPlaw::mc(const GEnergy&amp; emin, const GEnergy&amp; emax,
347 356 /***********************************************************************//**
348 357 * @brief Autoscale normalization
349 358 *
350   - * Set the internal scaling of the normalization parameter so that the
351   - * parameter value equals 1.
  359 + * Based on the actual value of the m_norm parameter, set the scale of m_norm
  360 + * so that the value will be 1. If minimum and/or maximum value boundaries
  361 + * exist, the boundaries are also modified accordingly.
352 362 ***************************************************************************/
353 363 void GModelSpectralExpPlaw::autoscale(void)
354 364 {
... ... @@ -537,7 +547,7 @@ std::string GModelSpectralExpPlaw::print(void) const
537 547 result.append("=== GModelSpectralExpPlaw ===\n");
538 548 result.append(parformat("Number of parameters")+str(size()));
539 549 for (int i = 0; i < size(); ++i)
540   - result.append("\n"+m_par[i]->print());
  550 + result.append("\n"+m_pars[i]->print());
541 551  
542 552 // Return result
543 553 return result;
... ... @@ -555,43 +565,51 @@ std::string GModelSpectralExpPlaw::print(void) const
555 565 ***************************************************************************/
556 566 void GModelSpectralExpPlaw::init_members(void)
557 567 {
558   - // Initialise parameters
559   - m_npars = 4;
560   - m_par[0] = &m_norm;
561   - m_par[1] = &m_index;
562   - m_par[2] = &m_ecut;
563   - m_par[3] = &m_pivot;
564   -
565 568 // Initialise powerlaw normalisation
566   - m_norm = GModelPar();
  569 + m_norm.clear();
567 570 m_norm.name("Prefactor");
568 571 m_norm.unit("ph/cm2/s/MeV");
569 572 m_norm.value(1.0);
570 573 m_norm.scale(1.0);
571 574 m_norm.free();
  575 + m_norm.gradient(0.0);
  576 + m_norm.hasgrad(true);
572 577  
573 578 // Initialise powerlaw index
574   - m_index = GModelPar();
  579 + m_index.clear();
575 580 m_index.name("Index");
576 581 m_index.value(-2.0);
577 582 m_index.scale(1.0);
578 583 m_index.free();
  584 + m_index.gradient(0.0);
  585 + m_index.hasgrad(true);
579 586  
580 587 // Initialise cut off energy
581   - m_ecut = GModelPar();
  588 + m_ecut.clear();
582 589 m_ecut.name("Cutoff");
583 590 m_ecut.unit("MeV");
584 591 m_ecut.value(1000.0);
585 592 m_ecut.scale(1.0);
586 593 m_ecut.free();
  594 + m_ecut.gradient(0.0);
  595 + m_ecut.hasgrad(true);
587 596  
588 597 // Initialise pivot energy
589   - m_pivot = GModelPar();
  598 + m_pivot.clear();
590 599 m_pivot.name("PivotEnergy");
591 600 m_pivot.unit("MeV");
592 601 m_pivot.value(100.0);
593 602 m_pivot.scale(1.0);
594 603 m_pivot.fix();
  604 + m_pivot.gradient(0.0);
  605 + m_pivot.hasgrad(true);
  606 +
  607 + // Set parameter pointer(s)
  608 + m_pars.clear();
  609 + m_pars.push_back(&m_norm);
  610 + m_pars.push_back(&m_index);
  611 + m_pars.push_back(&m_ecut);
  612 + m_pars.push_back(&m_pivot);
595 613  
596 614 // Return
597 615 return;
... ... @@ -601,17 +619,23 @@ void GModelSpectralExpPlaw::init_members(void)
601 619 /***********************************************************************//**
602 620 * @brief Copy class members
603 621 *
604   - * @param[in] model GModelSpectralExpPlaw members which should be copied.
  622 + * @param[in] model Exponential cut off power law model.
605 623 ***************************************************************************/
606 624 void GModelSpectralExpPlaw::copy_members(const GModelSpectralExpPlaw& model)
607 625 {
608   - // Copy model parameters (we do not need to copy the rest since it is
609   - // static)
  626 + // Copy members
610 627 m_norm = model.m_norm;
611 628 m_index = model.m_index;
612 629 m_ecut = model.m_ecut;
613 630 m_pivot = model.m_pivot;
614 631  
  632 + // Set parameter pointer(s)
  633 + m_pars.clear();
  634 + m_pars.push_back(&m_norm);
  635 + m_pars.push_back(&m_index);
  636 + m_pars.push_back(&m_ecut);
  637 + m_pars.push_back(&m_pivot);
  638 +
615 639 // Return
616 640 return;
617 641 }
... ...
src/model/GModelSpectralFunc.cpp
  View file @5af896d
... ... @@ -58,7 +58,7 @@ const GModelSpectralRegistry g_spectral_func_registry(&amp;g_spectral_func_seed);
58 58 ***************************************************************************/
59 59 GModelSpectralFunc::GModelSpectralFunc(void) : GModelSpectral()
60 60 {
61   - // Initialise private members for clean destruction
  61 + // Initialise members
62 62 init_members();
63 63  
64 64 // Return
... ... @@ -67,17 +67,18 @@ GModelSpectralFunc::GModelSpectralFunc(void) : GModelSpectral()
67 67  
68 68  
69 69 /***********************************************************************//**
70   - * @brief Constructor
  70 + * @brief File constructor
71 71 *
72 72 * @param[in] filename File name of nodes.
73 73 *
74   - * Construct an arbitrary spectral function from a list of nodes that is
75   - * found in the specified file.
  74 + * Creates instance of a spectral function model from a list of nodes that is
  75 + * found in the specified file. See GModelSpectralFunc::load_nodes() for more
  76 + * information about the expected structure of the file.
76 77 ***************************************************************************/
77 78 GModelSpectralFunc::GModelSpectralFunc(const std::string& filename) :
78 79 GModelSpectral()
79 80 {
80   - // Initialise private members for clean destruction
  81 + // Initialise members
81 82 init_members();
82 83  
83 84 // Load nodes
... ... @@ -89,16 +90,18 @@ GModelSpectralFunc::GModelSpectralFunc(const std::string&amp; filename) :
89 90  
90 91  
91 92 /***********************************************************************//**
92   - * @brief Constructor
  93 + * @brief XML constructor
93 94 *
94 95 * @param[in] xml XML element.
95 96 *
96   - * Construct an arbitrary spectral function from a XML element.
  97 + * Creates instance of a spectral function model by extracting information
  98 + * from an XML element. See GModelSpectralFunc::read() for more information
  99 + * about the expected structure of the XML element.
97 100 ***************************************************************************/
98 101 GModelSpectralFunc::GModelSpectralFunc(const GXmlElement& xml) :
99 102 GModelSpectral()
100 103 {
101   - // Initialise private members for clean destruction
  104 + // Initialise members
102 105 init_members();
103 106  
104 107 // Read information from XML element
... ... @@ -112,12 +115,12 @@ GModelSpectralFunc::GModelSpectralFunc(const GXmlElement&amp; xml) :
112 115 /***********************************************************************//**
113 116 * @brief Copy constructor
114 117 *
115   - * @param[in] model Spectral component.
  118 + * @param[in] model Spectral function model.
116 119 ***************************************************************************/
117 120 GModelSpectralFunc::GModelSpectralFunc(const GModelSpectralFunc& model) :
118 121 GModelSpectral(model)
119 122 {
120   - // Initialise private members for clean destruction
  123 + // Initialise members
121 124 init_members();
122 125  
123 126 // Copy members
... ... @@ -150,7 +153,7 @@ GModelSpectralFunc::~GModelSpectralFunc(void)
150 153 /***********************************************************************//**
151 154 * @brief Assignment operator
152 155 *
153   - * @param[in] model Spectral component.
  156 + * @param[in] model Spectral function model.
154 157 ***************************************************************************/
155 158 GModelSpectralFunc& GModelSpectralFunc::operator= (const GModelSpectralFunc& model)
156 159 {
... ... @@ -163,7 +166,7 @@ GModelSpectralFunc&amp; GModelSpectralFunc::operator= (const GModelSpectralFunc&amp; mod
163 166 // Free members
164 167 free_members();
165 168  
166   - // Initialise private members for clean destruction
  169 + // Initialise members
167 170 init_members();
168 171  
169 172 // Copy members
... ... @@ -221,7 +224,7 @@ GModelSpectralFunc* GModelSpectralFunc::clone(void) const
221 224 * Note that the node energies are stored as log10 of energy in units of
222 225 * MeV.
223 226 ***************************************************************************/
224   -double GModelSpectralFunc::eval(const GEnergy& srcEng)
  227 +double GModelSpectralFunc::eval(const GEnergy& srcEng) const
225 228 {
226 229 // Interpolate function
227 230 double func = m_nodes.interpolate(srcEng.log10MeV(), m_values);
... ... @@ -250,7 +253,7 @@ double GModelSpectralFunc::eval(const GEnergy&amp; srcEng)
250 253 * The partial derivative of the normalization value is given by
251 254 * \f[dI/dn_v=n_s f(E)\f]
252 255 ***************************************************************************/
253   -double GModelSpectralFunc::eval_gradients(const GEnergy& srcEng)
  256 +double GModelSpectralFunc::eval_gradients(const GEnergy& srcEng) const
254 257 {
255 258 // Interpolate function
256 259 double func = m_nodes.interpolate(srcEng.log10MeV(), m_values);
... ... @@ -261,8 +264,8 @@ double GModelSpectralFunc::eval_gradients(const GEnergy&amp; srcEng)
261 264 // Compute partial derivatives of the parameter values
262 265 double g_norm = (m_norm.isfree()) ? m_norm.scale() * func : 0.0;
263 266  
264   - // Set gradients
265   - m_norm.gradient(g_norm);
  267 + // Set gradients (circumvent const correctness)
  268 + ((GModelSpectralFunc*)this)->m_norm.gradient(g_norm);
266 269  
267 270 // Return
268 271 return value;
... ... @@ -275,6 +278,9 @@ double GModelSpectralFunc::eval_gradients(const GEnergy&amp; srcEng)
275 278 * @param[in] emin Minimum photon energy.
276 279 * @param[in] emax Maximum photon energy.
277 280 *
  281 + * @exception GException::feature_not_implemented
  282 + * Method not yet implemented
  283 + *
278 284 * Computes
279 285 * \f[\int_{E_{\rm min}}^{E_{\rm max}} I(E) dE\f]
280 286 * where
... ... @@ -301,6 +307,9 @@ double GModelSpectralFunc::flux(const GEnergy&amp; emin, const GEnergy&amp; emax) const
301 307 * @param[in] emax Maximum photon energy.
302 308 * @param[in] ran Random number generator.
303 309 *
  310 + * @exception GException::feature_not_implemented
  311 + * Method not yet implemented
  312 + *
304 313 * @todo To be implemented
305 314 ***************************************************************************/
306 315 GEnergy GModelSpectralFunc::mc(const GEnergy& emin, const GEnergy& emax,
... ... @@ -426,7 +435,7 @@ std::string GModelSpectralFunc::print(void) const
426 435 result.append(parformat("Number of nodes")+str(m_nodes.size()));
427 436 result.append(parformat("Number of parameters")+str(size()));
428 437 for (int i = 0; i < size(); ++i)
429   - result.append("\n"+m_par[i]->print());
  438 + result.append("\n"+m_pars[i]->print());
430 439  
431 440 // Return result
432 441 return result;
... ... @@ -444,18 +453,20 @@ std::string GModelSpectralFunc::print(void) const
444 453 ***************************************************************************/
445 454 void GModelSpectralFunc::init_members(void)
446 455 {
447   - // Initialise parameters
448   - m_npars = 1;
449   - m_par[0] = &m_norm;
450   -
451 456 // Initialise powerlaw normalisation
452   - m_norm = GModelPar();
  457 + m_norm.clear();
453 458 m_norm.name("Normalization");
454 459 m_norm.min(0.0);
455 460 m_norm.max(1000.0);
456 461 m_norm.value(1.0);
457 462 m_norm.scale(1.0);
458 463 m_norm.free();
  464 + m_norm.gradient(0.0);
  465 + m_norm.hasgrad(true);
  466 +
  467 + // Set parameter pointer(s)
  468 + m_pars.clear();
  469 + m_pars.push_back(&m_norm);
459 470  
460 471 // Initialise other members
461 472 m_nodes.clear();
... ... @@ -470,17 +481,20 @@ void GModelSpectralFunc::init_members(void)
470 481 /***********************************************************************//**
471 482 * @brief Copy class members
472 483 *
473   - * @param[in] model Spectral model component.
  484 + * @param[in] model Spectral function model.
474 485 ***************************************************************************/
475 486 void GModelSpectralFunc::copy_members(const GModelSpectralFunc& model)
476 487 {
477   - // Copy model parameters (we do not need to copy the rest since it is
478   - // static)
  488 + // Copy members
479 489 m_norm = model.m_norm;
480 490 m_nodes = model.m_nodes;
481 491 m_values = model.m_values;
482 492 m_filename = model.m_filename;
483 493  
  494 + // Set parameter pointer(s)
  495 + m_pars.clear();
  496 + m_pars.push_back(&m_norm);
  497 +
484 498 // Return
485 499 return;
486 500 }
... ...
src/model/GModelSpectralPlaw.cpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralPlaw.cpp
15   - * @brief GModelSpectralPlaw class implementation.
  15 + * @brief Power law spectral model class implementation
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -74,7 +74,7 @@ GModelSpectralPlaw::GModelSpectralPlaw(void) : GModelSpectral()
74 74 GModelSpectralPlaw::GModelSpectralPlaw(const double& norm,
75 75 const double& index) : GModelSpectral()
76 76 {
77   - // Initialise private members for clean destruction
  77 + // Initialise members
78 78 init_members();
79 79  
80 80 // Set parameters
... ... @@ -87,15 +87,17 @@ GModelSpectralPlaw::GModelSpectralPlaw(const double&amp; norm,
87 87  
88 88  
89 89 /***********************************************************************//**
90   - * @brief Constructor
  90 + * @brief XML constructor
91 91 *
92 92 * @param[in] xml XML element containing position information.
93 93 *
94   - * Construct a spectral power law from a XML element.
  94 + * Creates instance of a power law spectral model by extracting information
  95 + * from an XML element. See GModelSpectralPlaw::read() for more information
  96 + * about the expected structure of the XML element.
95 97 ***************************************************************************/
96 98 GModelSpectralPlaw::GModelSpectralPlaw(const GXmlElement& xml) : GModelSpectral()
97 99 {
98   - // Initialise private members for clean destruction
  100 + // Initialise members
99 101 init_members();
100 102  
101 103 // Read information from XML element
... ... @@ -109,12 +111,12 @@ GModelSpectralPlaw::GModelSpectralPlaw(const GXmlElement&amp; xml) : GModelSpectral(
109 111 /***********************************************************************//**
110 112 * @brief Copy constructor
111 113 *
112   - * @param[in] model Model from which the instance should be built.
  114 + * @param[in] model Spectral power law model.
113 115 ***************************************************************************/
114 116 GModelSpectralPlaw::GModelSpectralPlaw(const GModelSpectralPlaw& model) :
115 117 GModelSpectral(model)
116 118 {
117   - // Initialise private members for clean destruction
  119 + // Initialise members
118 120 init_members();
119 121  
120 122 // Copy members
... ... @@ -147,9 +149,9 @@ GModelSpectralPlaw::~GModelSpectralPlaw(void)
147 149 /***********************************************************************//**
148 150 * @brief Assignment operator
149 151 *
150   - * @param[in] model Model which should be assigned.
  152 + * @param[in] model Spectral power law model.
151 153 ***************************************************************************/
152   -GModelSpectralPlaw& GModelSpectralPlaw::operator= (const GModelSpectralPlaw& model)
  154 +GModelSpectralPlaw& GModelSpectralPlaw::operator=(const GModelSpectralPlaw& model)
153 155 {
154 156 // Execute only if object is not identical
155 157 if (this != &model) {
... ... @@ -160,7 +162,7 @@ GModelSpectralPlaw&amp; GModelSpectralPlaw::operator= (const GModelSpectralPlaw&amp; mod
160 162 // Free members
161 163 free_members();
162 164  
163   - // Initialise private members for clean destruction
  165 + // Initialise members
164 166 init_members();
165 167  
166 168 // Copy members
... ... @@ -223,7 +225,7 @@ GModelSpectralPlaw* GModelSpectralPlaw::clone(void) const
223 225 * Furthermore, the method expects that energy!=0. Otherwise Inf or NaN
224 226 * may result.
225 227 ***************************************************************************/
226   -double GModelSpectralPlaw::eval(const GEnergy& srcEng)
  228 +double GModelSpectralPlaw::eval(const GEnergy& srcEng) const
227 229 {
228 230 // Compute function value
229 231 double energy = srcEng.MeV() / pivot();
... ... @@ -260,7 +262,7 @@ double GModelSpectralPlaw::eval(const GEnergy&amp; srcEng)
260 262 * Furthermore, the method expects that energy!=0. Otherwise Inf or NaN
261 263 * may result.
262 264 ***************************************************************************/
263   -double GModelSpectralPlaw::eval_gradients(const GEnergy& srcEng)
  265 +double GModelSpectralPlaw::eval_gradients(const GEnergy& srcEng) const
264 266 {
265 267 // Compute function value
266 268 double energy = srcEng.MeV() / pivot();
... ... @@ -272,10 +274,10 @@ double GModelSpectralPlaw::eval_gradients(const GEnergy&amp; srcEng)
272 274 double g_index = (m_index.isfree()) ? value * m_index.scale() * std::log(energy) : 0.0;
273 275 double g_pivot = (m_pivot.isfree()) ? -value * index() / m_pivot.value() : 0.0;
274 276  
275   - // Set gradients
276   - m_norm.gradient(g_norm);
277   - m_index.gradient(g_index);
278   - m_pivot.gradient(g_pivot);
  277 + // Set gradients (circumvent const correctness)
  278 + ((GModelSpectralPlaw*)this)->m_norm.gradient(g_norm);
  279 + ((GModelSpectralPlaw*)this)->m_index.gradient(g_index);
  280 + ((GModelSpectralPlaw*)this)->m_pivot.gradient(g_pivot);
279 281  
280 282 // Return
281 283 return value;
... ... @@ -301,7 +303,8 @@ double GModelSpectralPlaw::flux(const GEnergy&amp; emin, const GEnergy&amp; emax) const
301 303 double flux = norm() * pow(pivot(), -index());
302 304 if (index() != -1.0) {
303 305 double exponent = index() + 1.0;
304   - flux *= (std::pow(emax.MeV(), exponent) - std::pow(emin.MeV(), exponent)) / exponent;
  306 + flux *= (std::pow(emax.MeV(), exponent)-std::pow(emin.MeV(), exponent)) /
  307 + exponent;
305 308 }
306 309 else
307 310 flux *= (std::log(emax.MeV()) - std::log(emin.MeV()));
... ... @@ -355,8 +358,9 @@ GEnergy GModelSpectralPlaw::mc(const GEnergy&amp; emin, const GEnergy&amp; emax,
355 358 /***********************************************************************//**
356 359 * @brief Autoscale normalization
357 360 *
358   - * Set the internal scaling of the normalization parameter so that the
359   - * parameter value equals 1.
  361 + * Based on the actual value of the m_norm parameter, set the scale of m_norm
  362 + * so that the value will be 1. If minimum and/or maximum value boundaries
  363 + * exist, the boundaries are also modified accordingly.
360 364 ***************************************************************************/
361 365 void GModelSpectralPlaw::autoscale(void)
362 366 {
... ... @@ -530,7 +534,7 @@ std::string GModelSpectralPlaw::print(void) const
530 534 result.append("=== GModelSpectralPlaw ===\n");
531 535 result.append(parformat("Number of parameters")+str(size()));
532 536 for (int i = 0; i < size(); ++i)
533   - result.append("\n"+m_par[i]->print());
  537 + result.append("\n"+m_pars[i]->print());
534 538  
535 539 // Return result
536 540 return result;
... ... @@ -548,34 +552,40 @@ std::string GModelSpectralPlaw::print(void) const
548 552 ***************************************************************************/
549 553 void GModelSpectralPlaw::init_members(void)
550 554 {
551   - // Initialise parameters
552   - m_npars = 3;
553   - m_par[0] = &m_norm;
554   - m_par[1] = &m_index;
555   - m_par[2] = &m_pivot;
556   -
557 555 // Initialise powerlaw normalisation
558   - m_norm = GModelPar();
  556 + m_norm.clear();
559 557 m_norm.name("Prefactor");
560 558 m_norm.unit("ph/cm2/s/MeV");
561 559 m_norm.value(1.0);
562 560 m_norm.scale(1.0);
563 561 m_norm.free();
  562 + m_norm.gradient(0.0);
  563 + m_norm.hasgrad(true);
564 564  
565 565 // Initialise powerlaw index
566   - m_index = GModelPar();
  566 + m_index.clear();
567 567 m_index.name("Index");
568 568 m_index.value(-2.0);
569 569 m_index.scale(1.0);
570 570 m_index.free();
  571 + m_index.gradient(0.0);
  572 + m_index.hasgrad(true);
571 573  
572 574 // Initialise pivot energy
573   - m_pivot = GModelPar();
  575 + m_pivot.clear();
574 576 m_pivot.name("PivotEnergy");
575 577 m_pivot.unit("MeV");
576 578 m_pivot.value(100.0);
577 579 m_pivot.scale(1.0);
578 580 m_pivot.fix();
  581 + m_pivot.gradient(0.0);
  582 + m_pivot.hasgrad(true);
  583 +
  584 + // Set parameter pointer(s)
  585 + m_pars.clear();
  586 + m_pars.push_back(&m_norm);
  587 + m_pars.push_back(&m_index);
  588 + m_pars.push_back(&m_pivot);
579 589  
580 590 // Return
581 591 return;
... ... @@ -589,12 +599,17 @@ void GModelSpectralPlaw::init_members(void)
589 599 ***************************************************************************/
590 600 void GModelSpectralPlaw::copy_members(const GModelSpectralPlaw& model)
591 601 {
592   - // Copy model parameters (we do not need to copy the rest since it is
593   - // static)
  602 + // Copy members
594 603 m_norm = model.m_norm;
595 604 m_index = model.m_index;
596 605 m_pivot = model.m_pivot;
597 606  
  607 + // Set parameter pointer(s)
  608 + m_pars.clear();
  609 + m_pars.push_back(&m_norm);
  610 + m_pars.push_back(&m_index);
  611 + m_pars.push_back(&m_pivot);
  612 +
598 613 // Return
599 614 return;
600 615 }
... ...
src/model/GModelSpectralPlaw2.cpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelSpectralPlaw2.cpp
15   - * @brief GModelSpectralPlaw2 class implementation.
  15 + * @brief Flux normalized power law spectral model class implementation
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -55,7 +55,7 @@ const GModelSpectralRegistry g_spectral_plaw2_registry(&amp;g_spectral_plaw2_seed);
55 55 ***************************************************************************/
56 56 GModelSpectralPlaw2::GModelSpectralPlaw2(void) : GModelSpectral()
57 57 {
58   - // Initialise private members for clean destruction
  58 + // Initialise members
59 59 init_members();
60 60  
61 61 // Return
... ... @@ -75,7 +75,7 @@ GModelSpectralPlaw2::GModelSpectralPlaw2(const double&amp; integral,
75 75 const double& index) :
76 76 GModelSpectral()
77 77 {
78   - // Initialise private members for clean destruction
  78 + // Initialise members
79 79 init_members();
80 80  
81 81 // Set parameters
... ... @@ -88,16 +88,18 @@ GModelSpectralPlaw2::GModelSpectralPlaw2(const double&amp; integral,
88 88  
89 89  
90 90 /***********************************************************************//**
91   - * @brief Constructor
  91 + * @brief XML constructor
92 92 *
93 93 * @param[in] xml XML element.
94 94 *
95   - * Construct a spectral power law from a XML element.
  95 + * Creates instance of a flux normalized power law spectral model by
  96 + * extracting information from an XML element. See GModelSpectralPlaw2::read()
  97 + * for more information about the expected structure of the XML element.
96 98 ***************************************************************************/
97 99 GModelSpectralPlaw2::GModelSpectralPlaw2(const GXmlElement& xml) :
98 100 GModelSpectral()
99 101 {
100   - // Initialise private members for clean destruction
  102 + // Initialise members
101 103 init_members();
102 104  
103 105 // Read information from XML element
... ... @@ -116,7 +118,7 @@ GModelSpectralPlaw2::GModelSpectralPlaw2(const GXmlElement&amp; xml) :
116 118 GModelSpectralPlaw2::GModelSpectralPlaw2(const GModelSpectralPlaw2& model) :
117 119 GModelSpectral(model)
118 120 {
119   - // Initialise private members for clean destruction
  121 + // Initialise members
120 122 init_members();
121 123  
122 124 // Copy members
... ... @@ -162,7 +164,7 @@ GModelSpectralPlaw2&amp; GModelSpectralPlaw2::operator= (const GModelSpectralPlaw2&amp;
162 164 // Free members
163 165 free_members();
164 166  
165   - // Initialise private members for clean destruction
  167 + // Initialise members
166 168 init_members();
167 169  
168 170 // Copy members
... ... @@ -223,10 +225,10 @@ GModelSpectralPlaw2* GModelSpectralPlaw2::clone(void) const
223 225 * \f$emin\f$ is the lower energy limit, and
224 226 * \f$emax\f$ is the upper energy limit.
225 227 ***************************************************************************/
226   -double GModelSpectralPlaw2::eval(const GEnergy& srcEng)
  228 +double GModelSpectralPlaw2::eval(const GEnergy& srcEng) const
227 229 {
228 230 // Update precomputed values
229   - update(srcEng);
  231 + ((GModelSpectralPlaw2*)this)->update(srcEng);
230 232  
231 233 // Compute function value
232 234 double value = integral() * m_norm * m_power;
... ... @@ -265,17 +267,17 @@ double GModelSpectralPlaw2::eval(const GEnergy&amp; srcEng)
265 267 * \f[dI/df_v=I(E) / f_v\f]
266 268 * \f[dI/di_v=I(E) \log(E) i_s\f]
267 269 * for \f$index = -1\f$.
268   - * The partial derivatives for the energy boundaries are fixed to 0 as it is
269   - * not expected the fit these parameters.
  270 + *
  271 + * No partial derivatives are supported for the energy boundaries.
270 272 ***************************************************************************/
271   -double GModelSpectralPlaw2::eval_gradients(const GEnergy& srcEng)
  273 +double GModelSpectralPlaw2::eval_gradients(const GEnergy& srcEng) const
272 274 {
273 275 // Initialise gradients
274 276 double g_integral = 0.0;
275 277 double g_index = 0.0;
276 278  
277 279 // Update precomputed values
278   - update(srcEng);
  280 + ((GModelSpectralPlaw2*)this)->update(srcEng);
279 281  
280 282 // Compute function value
281 283 double value = integral() * m_norm * m_power;
... ... @@ -289,10 +291,8 @@ double GModelSpectralPlaw2::eval_gradients(const GEnergy&amp; srcEng)
289 291 g_index = value * (m_g_norm + ln10*srcEng.log10MeV()) * m_index.scale();
290 292  
291 293 // Set gradients
292   - m_integral.gradient(g_integral);
293   - m_index.gradient(g_index);
294   - m_emin.gradient(0.0);
295   - m_emax.gradient(0.0);
  294 + ((GModelSpectralPlaw2*)this)->m_integral.gradient(g_integral);
  295 + ((GModelSpectralPlaw2*)this)->m_index.gradient(g_index);
296 296  
297 297 // Return
298 298 return value;
... ... @@ -526,7 +526,7 @@ std::string GModelSpectralPlaw2::print(void) const
526 526 result.append("=== GModelSpectralPlaw2 ===\n");
527 527 result.append(parformat("Number of parameters")+str(size()));
528 528 for (int i = 0; i < size(); ++i)
529   - result.append("\n"+m_par[i]->print());
  529 + result.append("\n"+m_pars[i]->print());
530 530  
531 531 // Return result
532 532 return result;
... ... @@ -544,47 +544,55 @@ std::string GModelSpectralPlaw2::print(void) const
544 544 ***************************************************************************/
545 545 void GModelSpectralPlaw2::init_members(void)
546 546 {
547   - // Initialise parameters
548   - m_npars = 4;
549   - m_par[0] = &m_integral;
550   - m_par[1] = &m_index;
551   - m_par[2] = &m_emin;
552   - m_par[3] = &m_emax;
553   -
554 547 // Initialise integral flux
555   - m_integral = GModelPar();
  548 + m_integral.clear();
556 549 m_integral.name("Integral");
557 550 m_integral.unit("ph/cm2/s");
558 551 m_integral.value(1.0);
559 552 m_integral.scale(1.0);
560 553 m_integral.range(0.0, 10.0);
561 554 m_integral.free();
  555 + m_integral.gradient(0.0);
  556 + m_integral.hasgrad(true);
562 557  
563 558 // Initialise powerlaw index
564   - m_index = GModelPar();
  559 + m_index.clear();
565 560 m_index.name("Index");
566 561 m_index.value(-2.0);
567 562 m_index.scale(1.0);
568 563 m_index.range(-5.0, -0.1);
569 564 m_index.free();
  565 + m_index.gradient(0.0);
  566 + m_index.hasgrad(true);
570 567  
571 568 // Initialise lower limit
572   - m_emin = GModelPar();
  569 + m_emin.clear();
573 570 m_emin.name("LowerLimit");
574 571 m_emin.unit("MeV");
575 572 m_emin.value(100.0);
576 573 m_emin.scale(1.0);
577 574 m_emin.range(0.001, 1.0e15);
578 575 m_emin.fix();
  576 + m_emin.gradient(0.0);
  577 + m_emin.hasgrad(false);
579 578  
580 579 // Initialise upper limit
581   - m_emax = GModelPar();
  580 + m_emax.clear();
582 581 m_emax.name("UpperLimit");
583 582 m_emax.unit("MeV");
584 583 m_emax.value(500000.0);
585 584 m_emax.scale(1.0);
586 585 m_emax.range(0.001, 1.0e15);
587 586 m_emax.fix();
  587 + m_emax.gradient(0.0);
  588 + m_emax.hasgrad(false);
  589 +
  590 + // Set parameter pointer(s)
  591 + m_pars.clear();
  592 + m_pars.push_back(&m_integral);
  593 + m_pars.push_back(&m_index);
  594 + m_pars.push_back(&m_emin);
  595 + m_pars.push_back(&m_emax);
588 596  
589 597 // Initialise last parameters (for fast computation)
590 598 m_log_emin = 0.0;
... ... @@ -614,12 +622,19 @@ void GModelSpectralPlaw2::init_members(void)
614 622 ***************************************************************************/
615 623 void GModelSpectralPlaw2::copy_members(const GModelSpectralPlaw2& model)
616 624 {
617   - // Copy model parameters (no need to copy others since they are static)
  625 + // Copy members
618 626 m_integral = model.m_integral;
619 627 m_index = model.m_index;
620 628 m_emin = model.m_emin;
621 629 m_emax = model.m_emax;
622 630  
  631 + // Set parameter pointer(s)
  632 + m_pars.clear();
  633 + m_pars.push_back(&m_integral);
  634 + m_pars.push_back(&m_index);
  635 + m_pars.push_back(&m_emin);
  636 + m_pars.push_back(&m_emax);
  637 +
623 638 // Copy bookkeeping information
624 639 m_log_emin = model.m_log_emin;
625 640 m_log_emax = model.m_log_emax;
... ...
src/model/GModelTemporalConst.cpp
  View file @5af896d
... ... @@ -12,7 +12,7 @@
12 12 ***************************************************************************/
13 13 /**
14 14 * @file GModelTemporalConst.cpp
15   - * @brief GModelTemporalConst class implementation.
  15 + * @brief Constant temporal model class interface implementation
16 16 * @author J. Knodlseder
17 17 */
18 18  
... ... @@ -51,7 +51,7 @@ const GModelTemporalRegistry g_temporal_const_registry(&amp;g_temporal_const_seed);
51 51 ***************************************************************************/
52 52 GModelTemporalConst::GModelTemporalConst(void) : GModelTemporal()
53 53 {
54   - // Initialise private members for clean destruction
  54 + // Initialise members
55 55 init_members();
56 56  
57 57 // Return
... ... @@ -62,12 +62,12 @@ GModelTemporalConst::GModelTemporalConst(void) : GModelTemporal()
62 62 /***********************************************************************//**
63 63 * @brief Copy constructor
64 64 *
65   - * @param[in] model Model from which the instance should be built.
  65 + * @param[in] model Constant temporal model
66 66 ***************************************************************************/
67 67 GModelTemporalConst::GModelTemporalConst(const GModelTemporalConst& model) :
68 68 GModelTemporal(model)
69 69 {
70   - // Initialise private members for clean destruction
  70 + // Initialise members
71 71 init_members();
72 72  
73 73 // Copy members
... ... @@ -100,7 +100,7 @@ GModelTemporalConst::~GModelTemporalConst(void)
100 100 /***********************************************************************//**
101 101 * @brief Assignment operator
102 102 *
103   - * @param[in] model Model which should be assigned.
  103 + * @param[in] model Constant temporal model
104 104 ***************************************************************************/
105 105 GModelTemporalConst& GModelTemporalConst::operator= (const GModelTemporalConst& model)
106 106 {
... ... @@ -113,7 +113,7 @@ GModelTemporalConst&amp; GModelTemporalConst::operator= (const GModelTemporalConst&amp;
113 113 // Free members
114 114 free_members();
115 115  
116   - // Initialise private members for clean destruction
  116 + // Initialise members
117 117 init_members();
118 118  
119 119 // Copy members
... ... @@ -167,10 +167,13 @@ GModelTemporalConst* GModelTemporalConst::clone(void) const
167 167 * This method implements the temporal component of a constant model. It
168 168 * simply returns the normalization constant (which typically is set to 1).
169 169 ***************************************************************************/
170   -double GModelTemporalConst::eval(const GTime& srcTime)
  170 +double GModelTemporalConst::eval(const GTime& srcTime) const
171 171 {
  172 + // Compute function value
  173 + double value = norm();
  174 +
172 175 // Return
173   - return (norm());
  176 + return value;
174 177 }
175 178  
176 179  
... ... @@ -179,17 +182,30 @@ double GModelTemporalConst::eval(const GTime&amp; srcTime)
179 182 *
180 183 * @param[in] srcTime True photon arrival time (not used).
181 184 *
182   - * This method implements the temporal component of a constant model. It
183   - * simply returns the normalization constant (which typically is set to 1)
184   - * and a gradient of 0 (as it is not expected to fit this parameter).
  185 + * The temporal model is defined as
  186 + * \f[I(E)=norm\f]
  187 + * where
  188 + * \f$norm=n_s n_v\f$ is the normalization of the function.
  189 + * Note that the normalization is factorised into a scaling factor and a
  190 + * value and that the method is expected to return the gradient with respect
  191 + * to the parameter value \f$n_v\f$.
  192 + *
  193 + * The partial derivative of the normalization value is given by
  194 + * \f[dI/dn_v=n_s\f]
185 195 ***************************************************************************/
186   -double GModelTemporalConst::eval_gradients(const GTime& srcTime)
  196 +double GModelTemporalConst::eval_gradients(const GTime& srcTime) const
187 197 {
188   - // Set gradient to 0
189   - m_norm.gradient(0.0);
  198 + // Compute function value
  199 + double value = norm();
  200 +
  201 + // Compute partial derivatives of the parameter values
  202 + double g_norm = (m_norm.isfree()) ? m_norm.scale() : 0.0;
  203 +
  204 + // Set gradients (circumvent const correctness)
  205 + ((GModelTemporalConst*)this)->m_norm.gradient(g_norm);
190 206  
191 207 // Return
192   - return (norm());
  208 + return value;
193 209 }
194 210  
195 211  
... ... @@ -205,7 +221,7 @@ double GModelTemporalConst::eval_gradients(const GTime&amp; srcTime)
205 221 * event rate that is specified by the rate parameter.
206 222 ***************************************************************************/
207 223 GTimes GModelTemporalConst::mc(const double& rate, const GTime& tmin,
208   - const GTime& tmax, GRan& ran)
  224 + const GTime& tmax, GRan& ran) const
209 225 {
210 226 // Allocates empty vector of times
211 227 GTimes times;
... ... @@ -240,9 +256,9 @@ GTimes GModelTemporalConst::mc(const double&amp; rate, const GTime&amp; tmin,
240 256 /***********************************************************************//**
241 257 * @brief Read model from XML element
242 258 *
243   - * @param[in] xml XML element containing power law model information.
  259 + * @param[in] xml XML element.
244 260 *
245   - * @todo Not yet implemented.
  261 + * @todo To be implemented.
246 262 ***************************************************************************/
247 263 void GModelTemporalConst::read(const GXmlElement& xml)
248 264 {
... ... @@ -254,9 +270,9 @@ void GModelTemporalConst::read(const GXmlElement&amp; xml)
254 270 /***********************************************************************//**
255 271 * @brief Write model into XML element
256 272 *
257   - * @param[in] xml XML element into which model information is written.
  273 + * @param[in] xml XML element.
258 274 *
259   - * @todo Not yet implemented.
  275 + * @todo To be implemented.
260 276 ***************************************************************************/
261 277 void GModelTemporalConst::write(GXmlElement& xml) const
262 278 {
... ... @@ -277,7 +293,7 @@ std::string GModelTemporalConst::print(void) const
277 293 result.append("=== GModelTemporalConst ===\n");
278 294 result.append(parformat("Number of parameters")+str(size()));
279 295 for (int i = 0; i < size(); ++i)
280   - result.append("\n"+m_par[i]->print());
  296 + result.append("\n"+m_pars[i]->print());
281 297  
282 298 // Return result
283 299 return result;
... ... @@ -295,16 +311,18 @@ std::string GModelTemporalConst::print(void) const
295 311 ***************************************************************************/
296 312 void GModelTemporalConst::init_members(void)
297 313 {
298   - // Initialise parameters
299   - m_npars = 1;
300   - m_par[0] = &m_norm;
301   -
302 314 // Initialise normalisation parameter
303   - m_norm = GModelPar();
  315 + m_norm.clear();
304 316 m_norm.name("Constant");
305 317 m_norm.unit("(relative value)");
306 318 m_norm.value(1.0);
307 319 m_norm.fix();
  320 + m_norm.gradient(0.0);
  321 + m_norm.hasgrad(true);
  322 +
  323 + // Set parameter pointer(s)
  324 + m_pars.clear();
  325 + m_pars.push_back(&m_norm);
308 326  
309 327 // Return
310 328 return;
... ... @@ -314,13 +332,16 @@ void GModelTemporalConst::init_members(void)
314 332 /***********************************************************************//**
315 333 * @brief Copy class members
316 334 *
317   - * @param[in] model GModelTemporalConst members which should be copied.
  335 + * @param[in] model Constant temporal model
318 336 ***************************************************************************/
319 337 void GModelTemporalConst::copy_members(const GModelTemporalConst& model)
320 338 {
321   - // Copy model parameters (we do not need to copy the rest since it is
322   - // static)
323   - m_norm = model.m_norm;
  339 + // Copy members
  340 + m_norm = model.m_norm;
  341 +
  342 + // Set parameter pointer(s)
  343 + m_pars.clear();
  344 + m_pars.push_back(&m_norm);
324 345  
325 346 // Return
326 347 return;
... ...