Commit 4c752bf0c21f4dd29fb03c50380a9bda75fad3fc
1 parent
259defff
GFunction::eval method now takes a const double precision reference as argument.
The GFunction::eval interface has been revised and all affected GammaLib classes have been adapted.
Showing
27 changed files
with
310 additions
and
180 deletions
NEWS
| ... | ... | @@ -39,6 +39,10 @@ GFitsTable |
| 39 | 39 | The append_column() and insert_column() methods have been |
| 40 | 40 | renamed to append() and insert(). |
| 41 | 41 | |
| 42 | +GFunction | |
| 43 | + The eval() method now takes a const reference to a double precisio | |
| 44 | + value as argument. | |
| 45 | + | |
| 42 | 46 | GMatrixBase |
| 43 | 47 | GMatrix |
| 44 | 48 | GMatrixSparse (formerly GSparseMatrix) | ... | ... |
include/GDerivative.hpp
| ... | ... | @@ -36,7 +36,7 @@ |
| 36 | 36 | /***********************************************************************//** |
| 37 | 37 | * @class GDerivative |
| 38 | 38 | * |
| 39 | - * @brief Numerical derivatives class interface defintion | |
| 39 | + * @brief Numerical derivatives class | |
| 40 | 40 | * |
| 41 | 41 | * This class allows to compute numerical derivatives using various methods. |
| 42 | 42 | * The function to be derived is implemented by the abstract GFunction |
| ... | ... | @@ -53,23 +53,23 @@ public: |
| 53 | 53 | virtual ~GDerivative(void); |
| 54 | 54 | |
| 55 | 55 | // Operators |
| 56 | - GDerivative& operator= (const GDerivative& dx); | |
| 56 | + GDerivative& operator=(const GDerivative& dx); | |
| 57 | 57 | |
| 58 | 58 | // Methods |
| 59 | 59 | void clear(void); |
| 60 | 60 | GDerivative* clone(void) const; |
| 61 | - void max_iter(const int& max_iter) { m_max_iter=max_iter; } | |
| 62 | - void eps(const double& eps) { m_eps=eps; } | |
| 63 | - void step_frac(const double& f) { m_step_frac=f; } | |
| 64 | - void silent(const bool& silent) { m_silent=silent; } | |
| 65 | - const int& iter(void) const { return m_iter; } | |
| 66 | - const int& max_iter(void) const { return m_max_iter; } | |
| 67 | - const double& eps(void) const { return m_eps; } | |
| 68 | - const double& step_frac(void) const { return m_step_frac; } | |
| 69 | - const bool& silent(void) const { return m_silent; } | |
| 70 | - void function(GFunction* func) { m_func = func; } | |
| 71 | - const GFunction* function(void) const { return m_func; } | |
| 72 | - double value(const double& x, double step = 0.0); | |
| 61 | + void max_iter(const int& max_iter); | |
| 62 | + void eps(const double& eps); | |
| 63 | + void step_frac(const double& fraction); | |
| 64 | + void silent(const bool& silent); | |
| 65 | + const int& iter(void) const; | |
| 66 | + const int& max_iter(void) const; | |
| 67 | + const double& eps(void) const; | |
| 68 | + const double& step_frac(void) const; | |
| 69 | + const bool& silent(void) const; | |
| 70 | + void function(GFunction* func); | |
| 71 | + const GFunction* function(void) const; | |
| 72 | + double value(const double& x, const double& step = 0.0); | |
| 73 | 73 | double ridder(const double& x, const double& h, double& err); |
| 74 | 74 | double minuit2(const double& x, double& err); |
| 75 | 75 | double difference(const double& x, const double& h); |
| ... | ... | @@ -103,4 +103,143 @@ protected: |
| 103 | 103 | bool m_silent; //!< Suppress warnings |
| 104 | 104 | }; |
| 105 | 105 | |
| 106 | + | |
| 107 | +/***********************************************************************//** | |
| 108 | + * @brief Return number of iterations | |
| 109 | + * | |
| 110 | + * @return Number of iterations. | |
| 111 | + ***************************************************************************/ | |
| 112 | +inline | |
| 113 | +const int& GDerivative::iter(void) const | |
| 114 | +{ | |
| 115 | + return m_iter; | |
| 116 | +} | |
| 117 | + | |
| 118 | + | |
| 119 | +/***********************************************************************//** | |
| 120 | + * @brief Set maximum number of iterations | |
| 121 | + * | |
| 122 | + * @param[in] max_iter Maximum number of iterations. | |
| 123 | + ***************************************************************************/ | |
| 124 | +inline | |
| 125 | +void GDerivative::max_iter(const int& max_iter) | |
| 126 | +{ | |
| 127 | + m_max_iter = max_iter; | |
| 128 | + return; | |
| 129 | +} | |
| 130 | + | |
| 131 | + | |
| 132 | +/***********************************************************************//** | |
| 133 | + * @brief Return maximum number of iterations | |
| 134 | + * | |
| 135 | + * @return Maximum number of iterations. | |
| 136 | + ***************************************************************************/ | |
| 137 | +inline | |
| 138 | +const int& GDerivative::max_iter(void) const | |
| 139 | +{ | |
| 140 | + return m_max_iter; | |
| 141 | +} | |
| 142 | + | |
| 143 | + | |
| 144 | +/***********************************************************************//** | |
| 145 | + * @brief Set precision | |
| 146 | + * | |
| 147 | + * @param[in] eps Precision. | |
| 148 | + ***************************************************************************/ | |
| 149 | +inline | |
| 150 | +void GDerivative::eps(const double& eps) | |
| 151 | +{ | |
| 152 | + m_eps = eps; | |
| 153 | + return; | |
| 154 | +} | |
| 155 | + | |
| 156 | + | |
| 157 | +/***********************************************************************//** | |
| 158 | + * @brief Get precision | |
| 159 | + * | |
| 160 | + * @return Precision. | |
| 161 | + ***************************************************************************/ | |
| 162 | +inline | |
| 163 | +const double& GDerivative::eps(void) const | |
| 164 | +{ | |
| 165 | + return m_eps; | |
| 166 | +} | |
| 167 | + | |
| 168 | + | |
| 169 | +/***********************************************************************//** | |
| 170 | + * @brief Set step fraction | |
| 171 | + * | |
| 172 | + * @param[in] fraction Step fraction. | |
| 173 | + ***************************************************************************/ | |
| 174 | +inline | |
| 175 | +void GDerivative::step_frac(const double& fraction) | |
| 176 | +{ | |
| 177 | + m_step_frac = fraction; | |
| 178 | + return; | |
| 179 | +} | |
| 180 | + | |
| 181 | + | |
| 182 | +/***********************************************************************//** | |
| 183 | + * @brief Get step fraction | |
| 184 | + * | |
| 185 | + * @return Step fraction. | |
| 186 | + ***************************************************************************/ | |
| 187 | +inline | |
| 188 | +const double& GDerivative::step_frac(void) const | |
| 189 | +{ | |
| 190 | + return m_step_frac; | |
| 191 | +} | |
| 192 | + | |
| 193 | + | |
| 194 | +/***********************************************************************//** | |
| 195 | + * @brief Set silence flag | |
| 196 | + * | |
| 197 | + * @param[in] silent Silence flag. | |
| 198 | + ***************************************************************************/ | |
| 199 | +inline | |
| 200 | +void GDerivative::silent(const bool& silent) | |
| 201 | +{ | |
| 202 | + m_silent = silent; | |
| 203 | + return; | |
| 204 | +} | |
| 205 | + | |
| 206 | + | |
| 207 | +/***********************************************************************//** | |
| 208 | + * @brief Get silence flag | |
| 209 | + * | |
| 210 | + * @return True is class is silent, false otherwise. | |
| 211 | + ***************************************************************************/ | |
| 212 | +inline | |
| 213 | +const bool& GDerivative::silent(void) const | |
| 214 | +{ | |
| 215 | + return m_silent; | |
| 216 | +} | |
| 217 | + | |
| 218 | + | |
| 219 | +/***********************************************************************//** | |
| 220 | + * @brief Set function | |
| 221 | + * | |
| 222 | + * @param[in] function Function. | |
| 223 | + * | |
| 224 | + * Sets the function for which the derivate should be determined. | |
| 225 | + ***************************************************************************/ | |
| 226 | +inline | |
| 227 | +void GDerivative::function(GFunction* function) | |
| 228 | +{ | |
| 229 | + m_func = function; | |
| 230 | + return; | |
| 231 | +} | |
| 232 | + | |
| 233 | + | |
| 234 | +/***********************************************************************//** | |
| 235 | + * @brief Get function | |
| 236 | + * | |
| 237 | + * @return Function. | |
| 238 | + ***************************************************************************/ | |
| 239 | +inline | |
| 240 | +const GFunction* GDerivative::function(void) const | |
| 241 | +{ | |
| 242 | + return m_func; | |
| 243 | +} | |
| 244 | + | |
| 106 | 245 | #endif /* GDERIVATIVE_HPP */ | ... | ... |
include/GFunction.hpp
| 1 | 1 | /*************************************************************************** |
| 2 | - * GFunction.hpp - Single parameter function abstract base class * | |
| 2 | + * GFunction.hpp - Single parameter function abstract base class * | |
| 3 | 3 | * ----------------------------------------------------------------------- * |
| 4 | - * copyright (C) 2011 by Jurgen Knodlseder * | |
| 4 | + * copyright (C) 2011-2013 by Juergen Knoedlseder * | |
| 5 | 5 | * ----------------------------------------------------------------------- * |
| 6 | 6 | * * |
| 7 | 7 | * This program is free software: you can redistribute it and/or modify * |
| ... | ... | @@ -20,8 +20,8 @@ |
| 20 | 20 | ***************************************************************************/ |
| 21 | 21 | /** |
| 22 | 22 | * @file GFunction.hpp |
| 23 | - * @brief GFunction abstract virtual base class interface definition. | |
| 24 | - * @author J. Knodlseder | |
| 23 | + * @brief Single parameter function abstract base class definition | |
| 24 | + * @author Juergen Knoedlseder | |
| 25 | 25 | */ |
| 26 | 26 | |
| 27 | 27 | #ifndef GFUNCTION_HPP |
| ... | ... | @@ -33,7 +33,7 @@ |
| 33 | 33 | /***********************************************************************//** |
| 34 | 34 | * @class GFunction |
| 35 | 35 | * |
| 36 | - * @brief GFunction class interface defintion. | |
| 36 | + * @brief Single parameter function abstract base class | |
| 37 | 37 | * |
| 38 | 38 | * This class implements the abstract interface for a one parameter function. |
| 39 | 39 | * This function is for example used for integration or numerical computation |
| ... | ... | @@ -47,19 +47,19 @@ public: |
| 47 | 47 | |
| 48 | 48 | // Constructors and destructors |
| 49 | 49 | GFunction(void); |
| 50 | - GFunction(const GFunction& func); | |
| 50 | + GFunction(const GFunction& function); | |
| 51 | 51 | virtual ~GFunction(void); |
| 52 | 52 | |
| 53 | 53 | // Operators |
| 54 | - GFunction& operator= (const GFunction& func); | |
| 54 | + GFunction& operator=(const GFunction& function); | |
| 55 | 55 | |
| 56 | 56 | // Methods |
| 57 | - virtual double eval(double x) = 0; | |
| 57 | + virtual double eval(const double& x) = 0; | |
| 58 | 58 | |
| 59 | 59 | protected: |
| 60 | 60 | // Protected methods |
| 61 | 61 | void init_members(void); |
| 62 | - void copy_members(const GFunction& func); | |
| 62 | + void copy_members(const GFunction& function); | |
| 63 | 63 | void free_members(void); |
| 64 | 64 | }; |
| 65 | 65 | ... | ... |
include/GModelSpectralExpPlaw.hpp
| ... | ... | @@ -121,7 +121,7 @@ protected: |
| 121 | 121 | m_index(index), |
| 122 | 122 | m_inv_pivot(1.0/pivot), |
| 123 | 123 | m_inv_ecut(1.0/ecut) {} |
| 124 | - double eval(double eng); | |
| 124 | + double eval(const double& eng); | |
| 125 | 125 | protected: |
| 126 | 126 | const double& m_norm; //!< Normalization |
| 127 | 127 | const double& m_index; //!< Index |
| ... | ... | @@ -140,7 +140,7 @@ protected: |
| 140 | 140 | m_index(index), |
| 141 | 141 | m_inv_pivot(1.0/pivot), |
| 142 | 142 | m_inv_ecut(1.0/ecut) {} |
| 143 | - double eval(double eng); | |
| 143 | + double eval(const double& eng); | |
| 144 | 144 | protected: |
| 145 | 145 | const double& m_norm; //!< Normalization |
| 146 | 146 | const double& m_index; //!< Index | ... | ... |
include/GModelSpectralLogParabola.hpp
| ... | ... | @@ -127,7 +127,7 @@ protected: |
| 127 | 127 | m_pivot(pivot) {}; |
| 128 | 128 | |
| 129 | 129 | // Method |
| 130 | - double eval(double x) { | |
| 130 | + double eval(const double& x) { | |
| 131 | 131 | double xrel = x/m_pivot.MeV(); |
| 132 | 132 | return m_norm*std::pow(xrel, m_index + m_curvature * std::log(xrel)); |
| 133 | 133 | } |
| ... | ... | @@ -151,7 +151,7 @@ protected: |
| 151 | 151 | flux_kern(norm, index, curvature, pivot) {}; |
| 152 | 152 | |
| 153 | 153 | // Method |
| 154 | - double eval(double x) { | |
| 154 | + double eval(const double& x) { | |
| 155 | 155 | return x * flux_kern::eval(x); |
| 156 | 156 | } |
| 157 | 157 | }; | ... | ... |
include/GObservation.hpp
| ... | ... | @@ -83,9 +83,11 @@ public: |
| 83 | 83 | virtual std::string print(const GChatter& chatter = NORMAL) const = 0; |
| 84 | 84 | |
| 85 | 85 | // Virtual methods |
| 86 | - virtual double model(const GModels& models, const GEvent& event, | |
| 86 | + virtual double model(const GModels& models, | |
| 87 | + const GEvent& event, | |
| 88 | + GVector* gradient = NULL) const; | |
| 89 | + virtual double npred(const GModels& models, | |
| 87 | 90 | GVector* gradient = NULL) const; |
| 88 | - virtual double npred(const GModels& models, GVector* gradient = NULL) const; | |
| 89 | 91 | |
| 90 | 92 | // Implemented methods |
| 91 | 93 | void name(const std::string& name); |
| ... | ... | @@ -98,8 +100,11 @@ public: |
| 98 | 100 | const std::string& statistics(void) const { return m_statistics; } |
| 99 | 101 | |
| 100 | 102 | // Other methods |
| 101 | - virtual double model_grad(const GModel& model, const GEvent& event, int ipar) const; | |
| 102 | - virtual double npred_grad(const GModel& model, int ipar) const; | |
| 103 | + virtual double model_grad(const GModel& model, | |
| 104 | + const GEvent& event, | |
| 105 | + const int& ipar) const; | |
| 106 | + virtual double npred_grad(const GModel& model, | |
| 107 | + const int& ipar) const; | |
| 103 | 108 | |
| 104 | 109 | protected: |
| 105 | 110 | // Protected methods |
| ... | ... | @@ -119,7 +124,7 @@ protected: |
| 119 | 124 | m_model(&model), |
| 120 | 125 | m_event(&event), |
| 121 | 126 | m_ipar(ipar) { } |
| 122 | - double eval(double x); | |
| 127 | + double eval(const double& x); | |
| 123 | 128 | protected: |
| 124 | 129 | const GObservation* m_parent; //!< Pointer to parent |
| 125 | 130 | const GModel* m_model; //!< Pointer to model |
| ... | ... | @@ -138,7 +143,7 @@ protected: |
| 138 | 143 | const GModel* model) : |
| 139 | 144 | m_parent(parent), |
| 140 | 145 | m_model(model) { } |
| 141 | - double eval(double x); | |
| 146 | + double eval(const double& x); | |
| 142 | 147 | protected: |
| 143 | 148 | const GObservation* m_parent; //!< Pointer to parent |
| 144 | 149 | const GModel* m_model; //!< Pointer to source model |
| ... | ... | @@ -152,7 +157,7 @@ protected: |
| 152 | 157 | m_parent(parent), |
| 153 | 158 | m_model(model), |
| 154 | 159 | m_time(obsTime) { } |
| 155 | - double eval(double x); | |
| 160 | + double eval(const double& x); | |
| 156 | 161 | protected: |
| 157 | 162 | const GObservation* m_parent; //!< Pointer to parent |
| 158 | 163 | const GModel* m_model; //!< Pointer to source model |
| ... | ... | @@ -168,7 +173,7 @@ protected: |
| 168 | 173 | m_parent(parent), |
| 169 | 174 | m_model(&model), |
| 170 | 175 | m_ipar(ipar) { return; } |
| 171 | - double eval(double x); | |
| 176 | + double eval(const double& x); | |
| 172 | 177 | protected: |
| 173 | 178 | const GObservation* m_parent; //!< Pointer to parent |
| 174 | 179 | const GModel* m_model; //!< Pointer to model | ... | ... |
include/GResponse.hpp
| ... | ... | @@ -129,7 +129,7 @@ protected: |
| 129 | 129 | m_srcTime(srcTime), |
| 130 | 130 | m_obs(obs), |
| 131 | 131 | m_rot(rot) { } |
| 132 | - double eval(double theta); | |
| 132 | + double eval(const double& theta); | |
| 133 | 133 | protected: |
| 134 | 134 | const GResponse& m_rsp; //!< Response |
| 135 | 135 | const GModelSpatialRadial& m_spatial; //!< Spatial model |
| ... | ... | @@ -157,7 +157,7 @@ protected: |
| 157 | 157 | m_theta(theta), |
| 158 | 158 | m_cos_theta(std::cos(theta)), |
| 159 | 159 | m_sin_theta(sin_theta) { } |
| 160 | - double eval(double phi); | |
| 160 | + double eval(const double& phi); | |
| 161 | 161 | protected: |
| 162 | 162 | const GResponse& m_rsp; //!< Response |
| 163 | 163 | const GEnergy& m_srcEng; //!< True photon energy |
| ... | ... | @@ -184,7 +184,7 @@ protected: |
| 184 | 184 | m_srcTime(srcTime), |
| 185 | 185 | m_obs(obs), |
| 186 | 186 | m_rot(rot) { } |
| 187 | - double eval(double theta); | |
| 187 | + double eval(const double& theta); | |
| 188 | 188 | protected: |
| 189 | 189 | const GResponse& m_rsp; //!< Response |
| 190 | 190 | const GModelSpatialElliptical& m_spatial; //!< Spatial model |
| ... | ... | @@ -214,7 +214,7 @@ protected: |
| 214 | 214 | m_theta(theta), |
| 215 | 215 | m_cos_theta(std::cos(theta)), |
| 216 | 216 | m_sin_theta(sin_theta) { } |
| 217 | - double eval(double phi); | |
| 217 | + double eval(const double& phi); | |
| 218 | 218 | protected: |
| 219 | 219 | const GResponse& m_rsp; //!< Response |
| 220 | 220 | const GModelSpatialElliptical& m_spatial; //!< Spatial model | ... | ... |
inst/cta/include/GCTAModelBackground.hpp
| ... | ... | @@ -115,7 +115,7 @@ protected: |
| 115 | 115 | m_cosdist(std::cos(dist)), |
| 116 | 116 | m_sindist(std::sin(dist)), |
| 117 | 117 | m_omega0(omega0) { } |
| 118 | - double eval(double theta); | |
| 118 | + double eval(const double& theta); | |
| 119 | 119 | protected: |
| 120 | 120 | const GModelSpatial* m_model; //!< Spatial model |
| 121 | 121 | const GEnergy& m_obsEng; //!< True photon energy |
| ... | ... | @@ -144,7 +144,7 @@ protected: |
| 144 | 144 | m_theta(theta), |
| 145 | 145 | m_cos_theta(std::cos(theta)), |
| 146 | 146 | m_sin_theta(sin_theta) { } |
| 147 | - double eval(double phi); | |
| 147 | + double eval(const double& phi); | |
| 148 | 148 | protected: |
| 149 | 149 | const GModelSpatial* m_model; //!< Spatial model |
| 150 | 150 | const GEnergy& m_obsEng; //!< True photon energy | ... | ... |
inst/cta/include/GCTAModelRadialAcceptance.hpp
| ... | ... | @@ -105,7 +105,7 @@ protected: |
| 105 | 105 | m_dist(dist), |
| 106 | 106 | m_cosdist(std::cos(dist)), |
| 107 | 107 | m_sindist(std::sin(dist)) { } |
| 108 | - double eval(double r); | |
| 108 | + double eval(const double& r); | |
| 109 | 109 | protected: |
| 110 | 110 | const GCTAModelRadial* m_parent; //!< Pointer to radial model |
| 111 | 111 | double m_roi; //!< ROI radius in radians | ... | ... |
inst/cta/include/GCTAModelRadialGauss.hpp
| ... | ... | @@ -91,7 +91,7 @@ protected: |
| 91 | 91 | class integrand : public GFunction { |
| 92 | 92 | public: |
| 93 | 93 | integrand(double sigma) : m_sigma(sigma) { } |
| 94 | - double eval(double x) { | |
| 94 | + double eval(const double& x) { | |
| 95 | 95 | double arg = x * x / m_sigma; |
| 96 | 96 | double arg2 = arg * arg; |
| 97 | 97 | double f = std::exp(-0.5 * arg2); | ... | ... |
inst/cta/include/GCTAModelRadialPolynom.hpp
| ... | ... | @@ -94,7 +94,7 @@ protected: |
| 94 | 94 | class integrand : public GFunction { |
| 95 | 95 | public: |
| 96 | 96 | integrand(const GCTAModelRadialPolynom* model) : m_model(model) { } |
| 97 | - double eval(double x) { | |
| 97 | + double eval(const double& x) { | |
| 98 | 98 | return (std::sin(x)*m_model->eval(x*gammalib::rad2deg)); |
| 99 | 99 | } |
| 100 | 100 | private: | ... | ... |
inst/cta/include/GCTAModelRadialProfile.hpp
| ... | ... | @@ -97,7 +97,7 @@ protected: |
| 97 | 97 | class integrand : public GFunction { |
| 98 | 98 | public: |
| 99 | 99 | integrand(const GCTAModelRadialProfile* model) : m_model(model) { } |
| 100 | - double eval(double x) { | |
| 100 | + double eval(const double& x) { | |
| 101 | 101 | return (std::sin(x)*m_model->eval(x*gammalib::rad2deg)); |
| 102 | 102 | } |
| 103 | 103 | private: | ... | ... |
inst/cta/src/GCTAModelBackground.cpp
| ... | ... | @@ -1212,7 +1212,7 @@ GModelTemporal* GCTAModelBackground::xml_temporal(const GXmlElement& temporal) c |
| 1212 | 1212 | * is limited to an arc around the vector connecting the model centre to |
| 1213 | 1213 | * the ROI centre. This limitation assures proper converges properly. |
| 1214 | 1214 | ***************************************************************************/ |
| 1215 | -double GCTAModelBackground::npred_roi_kern_theta::eval(double theta) | |
| 1215 | +double GCTAModelBackground::npred_roi_kern_theta::eval(const double& theta) | |
| 1216 | 1216 | { |
| 1217 | 1217 | // Initialise value |
| 1218 | 1218 | double value = 0.0; |
| ... | ... | @@ -1289,7 +1289,7 @@ double GCTAModelBackground::npred_roi_kern_theta::eval(double theta) |
| 1289 | 1289 | * multiplication is definitely not the fastest way to do that |
| 1290 | 1290 | * computation). |
| 1291 | 1291 | ***************************************************************************/ |
| 1292 | -double GCTAModelBackground::npred_roi_kern_phi::eval(double phi) | |
| 1292 | +double GCTAModelBackground::npred_roi_kern_phi::eval(const double& phi) | |
| 1293 | 1293 | { |
| 1294 | 1294 | // Initialise value |
| 1295 | 1295 | double value = 0.0; | ... | ... |
inst/cta/src/GCTAModelRadialAcceptance.cpp
| ... | ... | @@ -1041,7 +1041,7 @@ GModelTemporal* GCTAModelRadialAcceptance::xml_temporal(const GXmlElement& tempo |
| 1041 | 1041 | * the arc of the azimuth angle that lies within the ROI. This integration |
| 1042 | 1042 | * is done analytically using the "cta_roi_arclength" support function. |
| 1043 | 1043 | ***************************************************************************/ |
| 1044 | -double GCTAModelRadialAcceptance::roi_kern::eval(double offset) | |
| 1044 | +double GCTAModelRadialAcceptance::roi_kern::eval(const double& offset) | |
| 1045 | 1045 | { |
| 1046 | 1046 | // Circumvent const correctness |
| 1047 | 1047 | GCTAModelRadial* radial = const_cast<GCTAModelRadial*>(m_parent); |
| ... | ... | @@ -1071,10 +1071,3 @@ double GCTAModelRadialAcceptance::roi_kern::eval(double offset) |
| 1071 | 1071 | // Return |
| 1072 | 1072 | return value; |
| 1073 | 1073 | } |
| 1074 | - | |
| 1075 | - | |
| 1076 | -/*========================================================================== | |
| 1077 | - = = | |
| 1078 | - = Friends = | |
| 1079 | - = = | |
| 1080 | - ==========================================================================*/ | ... | ... |
inst/cta/src/GCTAModelRadialGauss.cpp
| ... | ... | @@ -618,10 +618,3 @@ void GCTAModelRadialGauss::free_members(void) |
| 618 | 618 | // Return |
| 619 | 619 | return; |
| 620 | 620 | } |
| 621 | - | |
| 622 | - | |
| 623 | -/*========================================================================== | |
| 624 | - = = | |
| 625 | - = Friends = | |
| 626 | - = = | |
| 627 | - ==========================================================================*/ | ... | ... |
inst/cta/src/GCTAResponse_helpers.cpp
| ... | ... | @@ -72,7 +72,7 @@ |
| 72 | 72 | * symmetric, hence it just multiplies the PSF value by the arclength times |
| 73 | 73 | * the sinus of the offset angle. |
| 74 | 74 | ***************************************************************************/ |
| 75 | -double cta_npsf_kern_rad_azsym::eval(double delta) | |
| 75 | +double cta_npsf_kern_rad_azsym::eval(const double& delta) | |
| 76 | 76 | { |
| 77 | 77 | // Initialise PSF value |
| 78 | 78 | double value = 0.0; |
| ... | ... | @@ -132,7 +132,7 @@ double cta_npsf_kern_rad_azsym::eval(double delta) |
| 132 | 132 | * |
| 133 | 133 | * for the zenith angle integration of radial models. |
| 134 | 134 | ***************************************************************************/ |
| 135 | -double cta_irf_radial_kern_rho::eval(double rho) | |
| 135 | +double cta_irf_radial_kern_rho::eval(const double& rho) | |
| 136 | 136 | { |
| 137 | 137 | // Compute half length of arc that lies within PSF validity circle |
| 138 | 138 | // (in radians) |
| ... | ... | @@ -243,7 +243,7 @@ double cta_irf_radial_kern_rho::eval(double rho) |
| 243 | 243 | * The angle \f$\theta\f$ is used in the computation of the IRF (no |
| 244 | 244 | * azimuthal dependence is so far implemented for the IRF computation). |
| 245 | 245 | ***************************************************************************/ |
| 246 | -double cta_irf_radial_kern_omega::eval(double omega) | |
| 246 | +double cta_irf_radial_kern_omega::eval(const double& omega) | |
| 247 | 247 | { |
| 248 | 248 | // Compute PSF offset angle [radians] |
| 249 | 249 | double delta = std::acos(m_cos_psf + m_sin_psf * std::cos(omega)); |
| ... | ... | @@ -301,7 +301,7 @@ double cta_irf_radial_kern_omega::eval(double omega) |
| 301 | 301 | * the ROI centre. This limitation assures that the integration converges |
| 302 | 302 | * properly. |
| 303 | 303 | ***************************************************************************/ |
| 304 | -double cta_npred_radial_kern_rho::eval(double rho) | |
| 304 | +double cta_npred_radial_kern_rho::eval(const double& rho) | |
| 305 | 305 | { |
| 306 | 306 | // Initialise Npred value |
| 307 | 307 | double npred = 0.0; |
| ... | ... | @@ -371,7 +371,7 @@ double cta_npred_radial_kern_rho::eval(double rho) |
| 371 | 371 | * multiplication is definitely not the fastest way to do that |
| 372 | 372 | * computation). |
| 373 | 373 | ***************************************************************************/ |
| 374 | -double cta_npred_radial_kern_omega::eval(double omega) | |
| 374 | +double cta_npred_radial_kern_omega::eval(const double& omega) | |
| 375 | 375 | { |
| 376 | 376 | // Compute sky direction vector in native coordinates |
| 377 | 377 | double cos_omega = std::cos(omega); |
| ... | ... | @@ -429,7 +429,7 @@ double cta_npred_radial_kern_omega::eval(double omega) |
| 429 | 429 | * d\omega |
| 430 | 430 | * \f] |
| 431 | 431 | ***************************************************************************/ |
| 432 | -double cta_irf_elliptical_kern_rho::eval(double rho) | |
| 432 | +double cta_irf_elliptical_kern_rho::eval(const double& rho) | |
| 433 | 433 | { |
| 434 | 434 | // Compute half length of arc that lies within PSF validity circle |
| 435 | 435 | // (in radians) |
| ... | ... | @@ -536,7 +536,7 @@ double cta_irf_elliptical_kern_rho::eval(double rho) |
| 536 | 536 | * The angle \f$\theta\f$ is used in the computation of the IRF (no |
| 537 | 537 | * azimuthal dependence is so far implemented for the IRF computation). |
| 538 | 538 | ***************************************************************************/ |
| 539 | -double cta_irf_elliptical_kern_omega::eval(double omega) | |
| 539 | +double cta_irf_elliptical_kern_omega::eval(const double& omega) | |
| 540 | 540 | { |
| 541 | 541 | // Initialise IRF value |
| 542 | 542 | double irf = 0.0; |
| ... | ... | @@ -607,7 +607,7 @@ double cta_irf_elliptical_kern_omega::eval(double omega) |
| 607 | 607 | * is limited to an arc around the vector connecting the model centre to |
| 608 | 608 | * the ROI centre. This limitation assures proper converges properly. |
| 609 | 609 | ***************************************************************************/ |
| 610 | -double cta_npred_elliptical_kern_rho::eval(double rho) | |
| 610 | +double cta_npred_elliptical_kern_rho::eval(const double& rho) | |
| 611 | 611 | { |
| 612 | 612 | // Initialise Npred value |
| 613 | 613 | double npred = 0.0; |
| ... | ... | @@ -677,7 +677,7 @@ double cta_npred_elliptical_kern_rho::eval(double rho) |
| 677 | 677 | * S_{\rm p}(\rho,\omega | E, t) \, N_{\rm pred}(\rho,\omega) |
| 678 | 678 | * \f] |
| 679 | 679 | ***************************************************************************/ |
| 680 | -double cta_npred_elliptical_kern_omega::eval(double omega) | |
| 680 | +double cta_npred_elliptical_kern_omega::eval(const double& omega) | |
| 681 | 681 | { |
| 682 | 682 | // Initialise Npred value |
| 683 | 683 | double npred = 0.0; |
| ... | ... | @@ -755,7 +755,7 @@ double cta_npred_elliptical_kern_omega::eval(double omega) |
| 755 | 755 | * Note that the integration is only performed for \f$\theta>0\f$. Otherwise |
| 756 | 756 | * zero is returned. |
| 757 | 757 | ***************************************************************************/ |
| 758 | -double cta_irf_diffuse_kern_theta::eval(double theta) | |
| 758 | +double cta_irf_diffuse_kern_theta::eval(const double& theta) | |
| 759 | 759 | { |
| 760 | 760 | // Initialise result |
| 761 | 761 | double irf = 0.0; |
| ... | ... | @@ -854,7 +854,7 @@ double cta_irf_diffuse_kern_theta::eval(double theta) |
| 854 | 854 | * @todo Implement azimuth angle computation of true photon in camera |
| 855 | 855 | * @todo Replace (theta,phi) by (delta,alpha) |
| 856 | 856 | ***************************************************************************/ |
| 857 | -double cta_irf_diffuse_kern_phi::eval(double phi) | |
| 857 | +double cta_irf_diffuse_kern_phi::eval(const double& phi) | |
| 858 | 858 | { |
| 859 | 859 | // Initialise result |
| 860 | 860 | double irf = 0.0; |
| ... | ... | @@ -943,7 +943,7 @@ double cta_irf_diffuse_kern_phi::eval(double phi) |
| 943 | 943 | * computation time and computation precision. A value of 1e-4 was judged |
| 944 | 944 | * appropriate. |
| 945 | 945 | ***************************************************************************/ |
| 946 | -double cta_npred_diffuse_kern_theta::eval(double theta) | |
| 946 | +double cta_npred_diffuse_kern_theta::eval(const double& theta) | |
| 947 | 947 | { |
| 948 | 948 | // Initialise Npred value |
| 949 | 949 | double npred = 0.0; |
| ... | ... | @@ -1003,7 +1003,7 @@ double cta_npred_diffuse_kern_theta::eval(double theta) |
| 1003 | 1003 | * multiplication is definitely not the fastest way to do that |
| 1004 | 1004 | * computation). |
| 1005 | 1005 | ***************************************************************************/ |
| 1006 | -double cta_npred_diffuse_kern_phi::eval(double phi) | |
| 1006 | +double cta_npred_diffuse_kern_phi::eval(const double& phi) | |
| 1007 | 1007 | { |
| 1008 | 1008 | // Initialise Npred value |
| 1009 | 1009 | double npred = 0.0; | ... | ... |
inst/cta/src/GCTAResponse_helpers.hpp
| ... | ... | @@ -78,7 +78,7 @@ public: |
| 78 | 78 | m_phi(phi), |
| 79 | 79 | m_zenith(zenith), |
| 80 | 80 | m_azimuth(azimuth) { } |
| 81 | - double eval(double delta); | |
| 81 | + double eval(const double& delta); | |
| 82 | 82 | protected: |
| 83 | 83 | const GCTAResponse& m_rsp; //!< CTA response function |
| 84 | 84 | const double& m_roi; //!< ROI radius in radians |
| ... | ... | @@ -153,7 +153,7 @@ public: |
| 153 | 153 | m_omega0(omega0), |
| 154 | 154 | m_delta_max(delta_max), |
| 155 | 155 | m_cos_delta_max(std::cos(delta_max)) { } |
| 156 | - double eval(double rho); | |
| 156 | + double eval(const double& rho); | |
| 157 | 157 | protected: |
| 158 | 158 | const GCTAResponse& m_rsp; //!< CTA response |
| 159 | 159 | const GModelSpatialRadial& m_model; //!< Radial spatial model |
| ... | ... | @@ -221,7 +221,7 @@ public: |
| 221 | 221 | m_sin_psf(sin_psf), |
| 222 | 222 | m_cos_ph(cos_ph), |
| 223 | 223 | m_sin_ph(sin_ph) { } |
| 224 | - double eval(double omega); | |
| 224 | + double eval(const double& omega); | |
| 225 | 225 | protected: |
| 226 | 226 | const GCTAResponse& m_rsp; //!< CTA response |
| 227 | 227 | const double& m_zenith; //!< Zenith angle |
| ... | ... | @@ -292,7 +292,7 @@ public: |
| 292 | 292 | m_radius(radius), |
| 293 | 293 | m_cos_radius(std::cos(radius)), |
| 294 | 294 | m_omega0(omega0) { } |
| 295 | - double eval(double rho); | |
| 295 | + double eval(const double& rho); | |
| 296 | 296 | protected: |
| 297 | 297 | const GCTAResponse& m_rsp; //!< CTA response |
| 298 | 298 | const GModelSpatialRadial& m_model; //!< Radial spatial model |
| ... | ... | @@ -345,7 +345,7 @@ public: |
| 345 | 345 | m_rot(rot), |
| 346 | 346 | m_sin_rho(sin_rho), |
| 347 | 347 | m_cos_rho(cos_rho) { } |
| 348 | - double eval(double omega); | |
| 348 | + double eval(const double& omega); | |
| 349 | 349 | protected: |
| 350 | 350 | const GCTAResponse& m_rsp; //!< CTA response |
| 351 | 351 | const GEnergy& m_srcEng; //!< True photon energy |
| ... | ... | @@ -419,7 +419,7 @@ public: |
| 419 | 419 | m_omega0(omega0), |
| 420 | 420 | m_delta_max(delta_max), |
| 421 | 421 | m_cos_delta_max(std::cos(delta_max)) { } |
| 422 | - double eval(double rho); | |
| 422 | + double eval(const double& rho); | |
| 423 | 423 | public: |
| 424 | 424 | const GCTAResponse& m_rsp; //!< CTA response |
| 425 | 425 | const GModelSpatialElliptical& m_model; //!< Spatial model |
| ... | ... | @@ -493,7 +493,7 @@ public: |
| 493 | 493 | m_sin_psf(sin_psf), |
| 494 | 494 | m_cos_ph(cos_ph), |
| 495 | 495 | m_sin_ph(sin_ph) { } |
| 496 | - double eval(double omega); | |
| 496 | + double eval(const double& omega); | |
| 497 | 497 | public: |
| 498 | 498 | const GCTAResponse& m_rsp; //!< CTA response |
| 499 | 499 | const GModelSpatialElliptical& m_model; //!< Spatial model |
| ... | ... | @@ -567,7 +567,7 @@ public: |
| 567 | 567 | m_radius(radius), |
| 568 | 568 | m_cos_radius(std::cos(radius)), |
| 569 | 569 | m_omega0(omega0) { } |
| 570 | - double eval(double rho); | |
| 570 | + double eval(const double& rho); | |
| 571 | 571 | protected: |
| 572 | 572 | const GCTAResponse& m_rsp; //!< CTA response |
| 573 | 573 | const GModelSpatialElliptical& m_model; //!< Spatial model |
| ... | ... | @@ -624,7 +624,7 @@ public: |
| 624 | 624 | m_rot(rot), |
| 625 | 625 | m_sin_rho(sin_rho), |
| 626 | 626 | m_cos_rho(cos_rho) { } |
| 627 | - double eval(double omega); | |
| 627 | + double eval(const double& omega); | |
| 628 | 628 | protected: |
| 629 | 629 | const GCTAResponse& m_rsp; //!< CTA response |
| 630 | 630 | const GModelSpatialElliptical& m_model; //!< Model |
| ... | ... | @@ -694,7 +694,7 @@ public: |
| 694 | 694 | m_rot(rot), |
| 695 | 695 | m_sin_eta(std::sin(eta)), |
| 696 | 696 | m_cos_eta(std::cos(eta)) { } |
| 697 | - double eval(double theta); | |
| 697 | + double eval(const double& theta); | |
| 698 | 698 | protected: |
| 699 | 699 | const GCTAResponse& m_rsp; //!< CTA response |
| 700 | 700 | const GModelSpatial& m_model; //!< Spatial model |
| ... | ... | @@ -764,7 +764,7 @@ public: |
| 764 | 764 | m_cos_theta(cos_theta), |
| 765 | 765 | m_sin_ph(sin_ph), |
| 766 | 766 | m_cos_ph(cos_ph) { } |
| 767 | - double eval(double phi); | |
| 767 | + double eval(const double& phi); | |
| 768 | 768 | protected: |
| 769 | 769 | const GCTAResponse& m_rsp; //!< CTA response |
| 770 | 770 | const GModelSpatial& m_model; //!< Spatial model |
| ... | ... | @@ -826,7 +826,7 @@ public: |
| 826 | 826 | m_srcTime(srcTime), |
| 827 | 827 | m_obs(obs), |
| 828 | 828 | m_rot(rot) { } |
| 829 | - double eval(double theta); | |
| 829 | + double eval(const double& theta); | |
| 830 | 830 | protected: |
| 831 | 831 | const GCTAResponse& m_rsp; //!< CTA response |
| 832 | 832 | const GModelSpatial& m_model; //!< Spatial model |
| ... | ... | @@ -876,7 +876,7 @@ public: |
| 876 | 876 | m_theta(theta), |
| 877 | 877 | m_cos_theta(std::cos(theta)), |
| 878 | 878 | m_sin_theta(sin_theta) { } |
| 879 | - double eval(double phi); | |
| 879 | + double eval(const double& phi); | |
| 880 | 880 | protected: |
| 881 | 881 | const GCTAResponse& m_rsp; //!< CTA response |
| 882 | 882 | const GModelSpatial& m_model; //!< Spatial model | ... | ... |
inst/lat/include/GLATPsfV1.hpp
| ... | ... | @@ -82,7 +82,7 @@ private: |
| 82 | 82 | double gcore, double gtail) : |
| 83 | 83 | m_ncore(ncore), m_ntail(ntail), m_sigma(sigma), |
| 84 | 84 | m_gcore(gcore), m_gtail(gtail) { } |
| 85 | - double eval(double x) { | |
| 85 | + double eval(const double& x) { | |
| 86 | 86 | double r = x / m_sigma; |
| 87 | 87 | double u = 0.5 * r * r; |
| 88 | 88 | double f = m_ncore * base_fct(u, m_gcore) + | ... | ... |
inst/lat/include/GLATPsfV3.hpp
| ... | ... | @@ -92,7 +92,7 @@ private: |
| 92 | 92 | m_ncore(ncore), m_ntail(ntail), |
| 93 | 93 | m_score(score), m_stail(stail), |
| 94 | 94 | m_gcore(gcore), m_gtail(gtail) { } |
| 95 | - double eval(double x) { | |
| 95 | + double eval(const double& x) { | |
| 96 | 96 | double rc = x / m_score; |
| 97 | 97 | double uc = 0.5 * rc * rc; |
| 98 | 98 | double rt = x / m_stail; | ... | ... |
pyext/GDerivative.i
| ... | ... | @@ -33,7 +33,7 @@ |
| 33 | 33 | /***********************************************************************//** |
| 34 | 34 | * @class GDerivative |
| 35 | 35 | * |
| 36 | - * @brief Numerical derivatives class interface defintion | |
| 36 | + * @brief Numerical derivatives class | |
| 37 | 37 | ***************************************************************************/ |
| 38 | 38 | class GDerivative : public GBase { |
| 39 | 39 | public: |
| ... | ... | @@ -48,7 +48,7 @@ public: |
| 48 | 48 | GDerivative* clone(void) const; |
| 49 | 49 | void max_iter(const int& max_iter); |
| 50 | 50 | void eps(const double& eps); |
| 51 | - void step_frac(const double& f); | |
| 51 | + void step_frac(const double& fraction); | |
| 52 | 52 | void silent(const bool& silent); |
| 53 | 53 | const int& iter(void) const; |
| 54 | 54 | const int& max_iter(void) const; |
| ... | ... | @@ -57,7 +57,7 @@ public: |
| 57 | 57 | const bool& silent(void) const; |
| 58 | 58 | void function(GFunction* func); |
| 59 | 59 | const GFunction* function(void) const; |
| 60 | - double value(const double& x, double step = 0.0); | |
| 60 | + double value(const double& x, const double& step = 0.0); | |
| 61 | 61 | double ridder(const double& x, const double& h, double& err); |
| 62 | 62 | double minuit2(const double& x, double& err); |
| 63 | 63 | double difference(const double& x, const double& h); | ... | ... |
pyext/GFunction.i
| 1 | 1 | /*************************************************************************** |
| 2 | 2 | * GFunction.i - Single parameter function abstract base class * |
| 3 | 3 | * ----------------------------------------------------------------------- * |
| 4 | - * copyright (C) 2011-2012 by Juergen Knoedlseder * | |
| 4 | + * copyright (C) 2011-2013 by Juergen Knoedlseder * | |
| 5 | 5 | * ----------------------------------------------------------------------- * |
| 6 | 6 | * * |
| 7 | 7 | * This program is free software: you can redistribute it and/or modify * |
| ... | ... | @@ -38,9 +38,9 @@ class GFunction { |
| 38 | 38 | public: |
| 39 | 39 | // Constructors and destructors |
| 40 | 40 | GFunction(void); |
| 41 | - GFunction(const GFunction& func); | |
| 41 | + GFunction(const GFunction& function); | |
| 42 | 42 | virtual ~GFunction(void); |
| 43 | 43 | |
| 44 | 44 | // Methods |
| 45 | - virtual double eval(double x) = 0; | |
| 45 | + virtual double eval(const double& x) = 0; | |
| 46 | 46 | }; | ... | ... |
src/model/GModelSpectralExpPlaw.cpp
| ... | ... | @@ -982,7 +982,7 @@ void GModelSpectralExpPlaw::update_mc_cache(const GEnergy& emin, |
| 982 | 982 | * |
| 983 | 983 | * @param[in] energy Energy (MeV). |
| 984 | 984 | ***************************************************************************/ |
| 985 | -double GModelSpectralExpPlaw::flux_kernel::eval(double energy) | |
| 985 | +double GModelSpectralExpPlaw::flux_kernel::eval(const double& energy) | |
| 986 | 986 | { |
| 987 | 987 | // Evaluate function value |
| 988 | 988 | double e_norm = energy * m_inv_pivot; |
| ... | ... | @@ -1000,7 +1000,7 @@ double GModelSpectralExpPlaw::flux_kernel::eval(double energy) |
| 1000 | 1000 | * |
| 1001 | 1001 | * @param[in] energy Energy (MeV). |
| 1002 | 1002 | ***************************************************************************/ |
| 1003 | -double GModelSpectralExpPlaw::eflux_kernel::eval(double energy) | |
| 1003 | +double GModelSpectralExpPlaw::eflux_kernel::eval(const double& energy) | |
| 1004 | 1004 | { |
| 1005 | 1005 | // Evaluate function value |
| 1006 | 1006 | double e_norm = energy * m_inv_pivot; | ... | ... |
src/numerics/GDerivative.cpp
| ... | ... | @@ -33,7 +33,9 @@ |
| 33 | 33 | #include "GTools.hpp" |
| 34 | 34 | |
| 35 | 35 | /* __ Method name definitions ____________________________________________ */ |
| 36 | +#define G_VALUE "GDerivative::value(double&, double&)" | |
| 36 | 37 | #define G_RIDDER "GDerivative::ridder(double&, double&, double&)" |
| 38 | +#define G_MINUIT2 "GDerivative::minuit2(double&, double&)" | |
| 37 | 39 | |
| 38 | 40 | /* __ Macros _____________________________________________________________ */ |
| 39 | 41 | |
| ... | ... | @@ -127,8 +129,9 @@ GDerivative::~GDerivative(void) |
| 127 | 129 | * @brief Assignment operator |
| 128 | 130 | * |
| 129 | 131 | * @param[in] dx Derivative. |
| 132 | + * @return Derivative. | |
| 130 | 133 | ***************************************************************************/ |
| 131 | -GDerivative& GDerivative::operator= (const GDerivative& dx) | |
| 134 | +GDerivative& GDerivative::operator=(const GDerivative& dx) | |
| 132 | 135 | { |
| 133 | 136 | // Execute only if object is not identical |
| 134 | 137 | if (this != &dx) { |
| ... | ... | @@ -156,7 +159,7 @@ GDerivative& GDerivative::operator= (const GDerivative& dx) |
| 156 | 159 | ==========================================================================*/ |
| 157 | 160 | |
| 158 | 161 | /***********************************************************************//** |
| 159 | - * @brief Clear instance | |
| 162 | + * @brief Clear derivative | |
| 160 | 163 | ***************************************************************************/ |
| 161 | 164 | void GDerivative::clear(void) |
| 162 | 165 | { |
| ... | ... | @@ -172,7 +175,9 @@ void GDerivative::clear(void) |
| 172 | 175 | |
| 173 | 176 | |
| 174 | 177 | /***********************************************************************//** |
| 175 | - * @brief Clone instance | |
| 178 | + * @brief Clone derivative | |
| 179 | + * | |
| 180 | + * @return Pointer to deep copy of derivative. | |
| 176 | 181 | ***************************************************************************/ |
| 177 | 182 | GDerivative* GDerivative::clone(void) const |
| 178 | 183 | { |
| ... | ... | @@ -195,7 +200,7 @@ GDerivative* GDerivative::clone(void) const |
| 195 | 200 | * The maximum number of iterations is controlled by the max_iter() |
| 196 | 201 | * method (by default, the maximum is set to 5). |
| 197 | 202 | ***************************************************************************/ |
| 198 | -double GDerivative::value(const double& x, double step) | |
| 203 | +double GDerivative::value(const double& x, const double& step) | |
| 199 | 204 | { |
| 200 | 205 | // Set constants |
| 201 | 206 | const double min_h = 1.0e-6; |
| ... | ... | @@ -248,14 +253,13 @@ double GDerivative::value(const double& x, double step) |
| 248 | 253 | // Compile option: signal if we exceed the tolerance |
| 249 | 254 | if (!silent()) { |
| 250 | 255 | if (err >= m_eps) { |
| 251 | - std::cout << "WARNING: GDerivative::value("; | |
| 252 | - std::cout << "iter=" << m_iter; | |
| 253 | - std::cout << ", x=" << x; | |
| 254 | - std::cout << ", dy/dx=" << result; | |
| 255 | - std::cout << ", h=" << h; | |
| 256 | - std::cout << ", err=" << err; | |
| 257 | - std::cout << "): error exceeds tolerance of "; | |
| 258 | - std::cout << m_eps << std::endl; | |
| 256 | + std::string msg = "Derivative uncertainty "+gammalib::str(err)+ | |
| 257 | + " exceeds tolerance of "+gammalib::str(m_eps)+ | |
| 258 | + " at function value "+gammalib::str(x)+ | |
| 259 | + " (df/dx="+gammalib::str(result)+ | |
| 260 | + ", iter="+gammalib::str(m_iter)+ | |
| 261 | + ", h="+gammalib::str(h)+")."; | |
| 262 | + gammalib::warning(G_VALUE, msg); | |
| 259 | 263 | } |
| 260 | 264 | } |
| 261 | 265 | |
| ... | ... | @@ -263,7 +267,6 @@ double GDerivative::value(const double& x, double step) |
| 263 | 267 | return result; |
| 264 | 268 | } |
| 265 | 269 | |
| 266 | - | |
| 267 | 270 | /***********************************************************************//** |
| 268 | 271 | * @brief Returns derivative by Ridders' method |
| 269 | 272 | * |
| ... | ... | @@ -480,15 +483,13 @@ double GDerivative::minuit2(const double& x, double& err) |
| 480 | 483 | // Compile option: signal if we exceed the tolerance |
| 481 | 484 | if (!silent()) { |
| 482 | 485 | if (err >= grad_tol) { |
| 483 | - std::cout << "WARNING: GDerivative::minuit("; | |
| 484 | - std::cout << "iter=" << m_iter; | |
| 485 | - std::cout << ", x=" << x; | |
| 486 | - std::cout << ", grd=" << grd; | |
| 487 | - std::cout << ", step=" << step; | |
| 488 | - std::cout << ", step_change=" << step_change; | |
| 489 | - std::cout << ", err=" << err; | |
| 490 | - std::cout << "): error exceeds tolerance of "; | |
| 491 | - std::cout << grad_tol << std::endl; | |
| 486 | + std::string msg = "Derivative uncertainty "+gammalib::str(err)+ | |
| 487 | + " exceeds tolerance of "+gammalib::str(grad_tol)+ | |
| 488 | + " at function value "+gammalib::str(x)+ | |
| 489 | + " (df/dx="+gammalib::str(grd)+ | |
| 490 | + ", step="+gammalib::str(step)+ | |
| 491 | + ", step_change="+gammalib::str(step_change)+")."; | |
| 492 | + gammalib::warning(G_MINUIT2, msg); | |
| 492 | 493 | } |
| 493 | 494 | } |
| 494 | 495 | ... | ... |
src/numerics/GFunction.cpp
| 1 | 1 | /*************************************************************************** |
| 2 | - * GFunction.cpp - Single parameter function abstract base class * | |
| 2 | + * GFunction.hpp - Single parameter function abstract base class * | |
| 3 | 3 | * ----------------------------------------------------------------------- * |
| 4 | - * copyright (C) 2011 by Jurgen Knodlseder * | |
| 4 | + * copyright (C) 2011-2013 by Juergen Knoedlseder * | |
| 5 | 5 | * ----------------------------------------------------------------------- * |
| 6 | 6 | * * |
| 7 | 7 | * This program is free software: you can redistribute it and/or modify * |
| ... | ... | @@ -20,8 +20,8 @@ |
| 20 | 20 | ***************************************************************************/ |
| 21 | 21 | /** |
| 22 | 22 | * @file GFunction.cpp |
| 23 | - * @brief GFunction abstract virtual base class implementation. | |
| 24 | - * @author J. Knodlseder | |
| 23 | + * @brief Single parameter function abstract base class implementation | |
| 24 | + * @author Juergen Knoedlseder | |
| 25 | 25 | */ |
| 26 | 26 | |
| 27 | 27 | /* __ Includes ___________________________________________________________ */ |
| ... | ... | @@ -58,15 +58,15 @@ GFunction::GFunction(void) |
| 58 | 58 | /***********************************************************************//** |
| 59 | 59 | * @brief Copy constructor |
| 60 | 60 | * |
| 61 | - * @param[in] func Function. | |
| 61 | + * @param[in] function Function. | |
| 62 | 62 | ***************************************************************************/ |
| 63 | -GFunction::GFunction(const GFunction& func) | |
| 63 | +GFunction::GFunction(const GFunction& function) | |
| 64 | 64 | { |
| 65 | 65 | // Initialise members |
| 66 | 66 | init_members(); |
| 67 | 67 | |
| 68 | 68 | // Copy members |
| 69 | - copy_members(func); | |
| 69 | + copy_members(function); | |
| 70 | 70 | |
| 71 | 71 | // Return |
| 72 | 72 | return; |
| ... | ... | @@ -95,12 +95,13 @@ GFunction::~GFunction(void) |
| 95 | 95 | /***********************************************************************//** |
| 96 | 96 | * @brief Assignment operator |
| 97 | 97 | * |
| 98 | - * @param[in] func Function. | |
| 98 | + * @param[in] function Function. | |
| 99 | + * @return Functions. | |
| 99 | 100 | ***************************************************************************/ |
| 100 | -GFunction& GFunction::operator= (const GFunction& func) | |
| 101 | +GFunction& GFunction::operator=(const GFunction& function) | |
| 101 | 102 | { |
| 102 | 103 | // Execute only if object is not identical |
| 103 | - if (this != &func) { | |
| 104 | + if (this != &function) { | |
| 104 | 105 | |
| 105 | 106 | // Free members |
| 106 | 107 | free_members(); |
| ... | ... | @@ -109,7 +110,7 @@ GFunction& GFunction::operator= (const GFunction& func) |
| 109 | 110 | init_members(); |
| 110 | 111 | |
| 111 | 112 | // Copy members |
| 112 | - copy_members(func); | |
| 113 | + copy_members(function); | |
| 113 | 114 | |
| 114 | 115 | } // endif: object was not identical |
| 115 | 116 | |
| ... | ... | @@ -143,9 +144,9 @@ void GFunction::init_members(void) |
| 143 | 144 | /***********************************************************************//** |
| 144 | 145 | * @brief Copy class members |
| 145 | 146 | * |
| 146 | - * @param[in] func Function. | |
| 147 | + * @param[in] function Function. | |
| 147 | 148 | ***************************************************************************/ |
| 148 | -void GFunction::copy_members(const GFunction& func) | |
| 149 | +void GFunction::copy_members(const GFunction& function) | |
| 149 | 150 | { |
| 150 | 151 | // Return |
| 151 | 152 | return; |
| ... | ... | @@ -160,10 +161,3 @@ void GFunction::free_members(void) |
| 160 | 161 | // Return |
| 161 | 162 | return; |
| 162 | 163 | } |
| 163 | - | |
| 164 | - | |
| 165 | -/*========================================================================== | |
| 166 | - = = | |
| 167 | - = Friends = | |
| 168 | - = = | |
| 169 | - ==========================================================================*/ | ... | ... |
src/obs/GObservation.cpp
| ... | ... | @@ -503,8 +503,9 @@ void GObservation::free_members(void) |
| 503 | 503 | * Eventually, the higher level method should avoid going in a |
| 504 | 504 | * parameter domain that is not defined. |
| 505 | 505 | ***************************************************************************/ |
| 506 | -double GObservation::model_grad(const GModel& model, const GEvent& event, | |
| 507 | - int ipar) const | |
| 506 | +double GObservation::model_grad(const GModel& model, | |
| 507 | + const GEvent& event, | |
| 508 | + const int& ipar) const | |
| 508 | 509 | { |
| 509 | 510 | // Initialise gradient |
| 510 | 511 | double grad = 0.0; |
| ... | ... | @@ -600,7 +601,7 @@ double GObservation::model_grad(const GModel& model, const GEvent& event, |
| 600 | 601 | * |
| 601 | 602 | * @param[in] x Function value. |
| 602 | 603 | ***************************************************************************/ |
| 603 | -double GObservation::model_func::eval(double x) | |
| 604 | +double GObservation::model_func::eval(const double& x) | |
| 604 | 605 | { |
| 605 | 606 | // Get non-const model pointer (circumvent const correctness) |
| 606 | 607 | GModel* model = const_cast<GModel*>(m_model); |
| ... | ... | @@ -670,7 +671,7 @@ double GObservation::model_func::eval(double x) |
| 670 | 671 | * Eventually, the higher level method should avoid going in a |
| 671 | 672 | * parameter domain that is not defined. |
| 672 | 673 | ***************************************************************************/ |
| 673 | -double GObservation::npred_grad(const GModel& model, int ipar) const | |
| 674 | +double GObservation::npred_grad(const GModel& model, const int& ipar) const | |
| 674 | 675 | { |
| 675 | 676 | // Initialise result |
| 676 | 677 | double grad = 0.0; |
| ... | ... | @@ -752,7 +753,7 @@ double GObservation::npred_grad(const GModel& model, int ipar) const |
| 752 | 753 | * |
| 753 | 754 | * @param[in] x Function value. |
| 754 | 755 | ***************************************************************************/ |
| 755 | -double GObservation::npred_func::eval(double x) | |
| 756 | +double GObservation::npred_func::eval(const double& x) | |
| 756 | 757 | { |
| 757 | 758 | // Get non-const model pointer (circumvent const correctness) |
| 758 | 759 | GModel* model = const_cast<GModel*>(m_model); |
| ... | ... | @@ -860,7 +861,7 @@ double GObservation::npred_temp(const GModel& model) const |
| 860 | 861 | * |
| 861 | 862 | * @param[in] x Function value. |
| 862 | 863 | ***************************************************************************/ |
| 863 | -double GObservation::npred_temp_kern::eval(double x) | |
| 864 | +double GObservation::npred_temp_kern::eval(const double& x) | |
| 864 | 865 | { |
| 865 | 866 | // Convert argument in native reference in seconds |
| 866 | 867 | GTime time; |
| ... | ... | @@ -952,42 +953,43 @@ double GObservation::npred_spec(const GModel& model, |
| 952 | 953 | * |
| 953 | 954 | * This method implements the integration kernel needed for the npred_spec() |
| 954 | 955 | * method. If G_LN_ENERGY_INT is defined the energy integration is done |
| 955 | - * logarithmically. | |
| 956 | + * logarithmically, i.e. @p x is given in ln(energy) instead of energy. | |
| 956 | 957 | ***************************************************************************/ |
| 957 | -double GObservation::npred_spec_kern::eval(double x) | |
| 958 | +double GObservation::npred_spec_kern::eval(const double& x) | |
| 958 | 959 | { |
| 959 | - #if defined(G_LN_ENERGY_INT) | |
| 960 | - // Variable substitution | |
| 961 | - #if defined(G_NAN_CHECK) | |
| 962 | - double x_in = x; | |
| 963 | - #endif | |
| 964 | - x = exp(x); | |
| 965 | - #endif | |
| 966 | - | |
| 967 | - // Set energy in MeV | |
| 960 | + // Set energy | |
| 968 | 961 | GEnergy eng; |
| 962 | + #if defined(G_LN_ENERGY_INT) | |
| 963 | + double expx = std::exp(x); | |
| 964 | + eng.MeV(expx); | |
| 965 | + #else | |
| 969 | 966 | eng.MeV(x); |
| 967 | + #endif | |
| 970 | 968 | |
| 971 | 969 | // Get function value |
| 972 | 970 | double value = m_model->npred(eng, *m_time, *m_parent); |
| 973 | 971 | |
| 974 | - #if defined(G_LN_ENERGY_INT) | |
| 975 | - // Correct for variable substitution | |
| 972 | + // Save value if needed | |
| 976 | 973 | #if defined(G_NAN_CHECK) |
| 977 | 974 | double value_out = value; |
| 978 | 975 | #endif |
| 979 | - value *= x; | |
| 976 | + | |
| 977 | + // Correct for variable substitution | |
| 978 | + #if defined(G_LN_ENERGY_INT) | |
| 979 | + value *= expx; | |
| 980 | 980 | #endif |
| 981 | 981 | |
| 982 | 982 | // Compile option: Check for NaN |
| 983 | 983 | #if defined(G_NAN_CHECK) |
| 984 | 984 | if (gammalib::isnotanumber(value) || gammalib::isinfinite(value)) { |
| 985 | 985 | std::cout << "*** ERROR: GObservation::npred_spec_kern::eval"; |
| 986 | - std::cout << "(x=" << x_in << "): "; | |
| 986 | + std::cout << "(x=" << x << "): "; | |
| 987 | 987 | std::cout << " NaN/Inf encountered"; |
| 988 | 988 | std::cout << " (value=" << value; |
| 989 | 989 | std::cout << " (value_out=" << value_out; |
| 990 | - std::cout << " x=" << x; | |
| 990 | + #if defined(G_LN_ENERGY_INT) | |
| 991 | + std::cout << " exp(x)=" << expx; | |
| 992 | + #endif | |
| 991 | 993 | std::cout << ")" << std::endl; |
| 992 | 994 | } |
| 993 | 995 | #endif | ... | ... |
src/obs/GResponse.cpp
| ... | ... | @@ -603,7 +603,7 @@ void GResponse::free_members(void) |
| 603 | 603 | * |
| 604 | 604 | * @param[in] theta Radial model offset angle (radians). |
| 605 | 605 | ***************************************************************************/ |
| 606 | -double GResponse::npred_radial_kern_theta::eval(double theta) | |
| 606 | +double GResponse::npred_radial_kern_theta::eval(const double& theta) | |
| 607 | 607 | { |
| 608 | 608 | // Initialise Npred value |
| 609 | 609 | double npred = 0.0; |
| ... | ... | @@ -650,7 +650,7 @@ double GResponse::npred_radial_kern_theta::eval(double theta) |
| 650 | 650 | * |
| 651 | 651 | * @param[in] phi Azimuth angle (radians). |
| 652 | 652 | ***************************************************************************/ |
| 653 | -double GResponse::npred_radial_kern_phi::eval(double phi) | |
| 653 | +double GResponse::npred_radial_kern_phi::eval(const double& phi) | |
| 654 | 654 | { |
| 655 | 655 | // Compute sky direction vector in native coordinates |
| 656 | 656 | double cos_phi = std::cos(phi); |
| ... | ... | @@ -693,7 +693,7 @@ double GResponse::npred_radial_kern_phi::eval(double phi) |
| 693 | 693 | * |
| 694 | 694 | * @param[in] theta Radial model offset angle (radians). |
| 695 | 695 | ***************************************************************************/ |
| 696 | -double GResponse::npred_elliptical_kern_theta::eval(double theta) | |
| 696 | +double GResponse::npred_elliptical_kern_theta::eval(const double& theta) | |
| 697 | 697 | { |
| 698 | 698 | // Initialise Npred value |
| 699 | 699 | double npred = 0.0; |
| ... | ... | @@ -737,7 +737,7 @@ double GResponse::npred_elliptical_kern_theta::eval(double theta) |
| 737 | 737 | * |
| 738 | 738 | * @param[in] phi Azimuth angle (radians). |
| 739 | 739 | ***************************************************************************/ |
| 740 | -double GResponse::npred_elliptical_kern_phi::eval(double phi) | |
| 740 | +double GResponse::npred_elliptical_kern_phi::eval(const double& phi) | |
| 741 | 741 | { |
| 742 | 742 | // Compute sky direction vector in native coordinates |
| 743 | 743 | double cos_phi = std::cos(phi); | ... | ... |
test/test_GNumerics.hpp
| 1 | 1 | /*************************************************************************** |
| 2 | - * test_GNumerics.hpp - test numerics modules * | |
| 2 | + * test_GNumerics.hpp - test numerics modules * | |
| 3 | 3 | * ----------------------------------------------------------------------- * |
| 4 | 4 | * copyright (C) 2010-2013 by Juergen Knoedlseder * |
| 5 | 5 | * ----------------------------------------------------------------------- * |
| ... | ... | @@ -42,7 +42,7 @@ class Gauss : public GFunction { |
| 42 | 42 | public: |
| 43 | 43 | Gauss(const double& sigma) : m_sigma(sigma) { return; } |
| 44 | 44 | virtual ~Gauss(void) { return; } |
| 45 | - double eval(double x) { | |
| 45 | + double eval(const double& x) { | |
| 46 | 46 | double arg = -0.5*x*x/m_sigma/m_sigma; |
| 47 | 47 | double val = 1.0/std::sqrt(gammalib::twopi)/m_sigma * std::exp(arg); |
| 48 | 48 | return val; |
| ... | ... | @@ -53,20 +53,19 @@ protected: |
| 53 | 53 | |
| 54 | 54 | class TestGNumerics : public GTestSuite |
| 55 | 55 | { |
| 56 | - public: | |
| 57 | - // Constructors and destructors | |
| 58 | - TestGNumerics(void) : GTestSuite(), m_sigma(1.0) { return; } | |
| 59 | - virtual ~TestGNumerics(void){ return; } | |
| 60 | - | |
| 61 | - // Methods | |
| 62 | - virtual void set(void); | |
| 63 | - void test_integral(void); | |
| 64 | - void test_romberg_integration(void); | |
| 56 | +public: | |
| 57 | + // Constructors and destructors | |
| 58 | + TestGNumerics(void) : GTestSuite(), m_sigma(1.0) { } | |
| 59 | + virtual ~TestGNumerics(void){ } | |
| 65 | 60 | |
| 66 | - // Private attributes | |
| 67 | - private: | |
| 68 | - double m_sigma; | |
| 61 | + // Methods | |
| 62 | + virtual void set(void); | |
| 63 | + void test_integral(void); | |
| 64 | + void test_romberg_integration(void); | |
| 69 | 65 | |
| 66 | +private: | |
| 67 | + // Private members | |
| 68 | + double m_sigma; | |
| 70 | 69 | }; |
| 71 | 70 | |
| 72 | 71 | #endif /* TEST_NUMERICS_HPP */ | ... | ... |