GammaLib  2.0.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GModelSpatialRadialProfileGauss.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  * GModelSpatialRadialProfileGauss.cpp - Gaussian radial profile class *
3  * ----------------------------------------------------------------------- *
4  * copyright (C) 2016-2021 by Juergen Knoedlseder *
5  * ----------------------------------------------------------------------- *
6  * *
7  * This program is free software: you can redistribute it and/or modify *
8  * it under the terms of the GNU General Public License as published by *
9  * the Free Software Foundation, either version 3 of the License, or *
10  * (at your option) any later version. *
11  * *
12  * This program is distributed in the hope that it will be useful, *
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15  * GNU General Public License for more details. *
16  * *
17  * You should have received a copy of the GNU General Public License *
18  * along with this program. If not, see <http://www.gnu.org/licenses/>. *
19  * *
20  ***************************************************************************/
21 /**
22  * @file GModelSpatialRadialProfileGauss.cpp
23  * @brief Radial Gaussian profile model class implementation
24  * @author Juergen Knoedlseder
25  */
26 
27 /* __ Includes ___________________________________________________________ */
28 #ifdef HAVE_CONFIG_H
29 #include <config.h>
30 #endif
31 #include "GException.hpp"
32 #include "GTools.hpp"
33 #include "GMath.hpp"
34 #include "GXmlElement.hpp"
37 
38 /* __ Constants __________________________________________________________ */
39 
40 /* __ Globals ____________________________________________________________ */
42 const GModelSpatialRegistry g_radial_disk_registry(&g_radial_disk_seed);
43 
44 /* __ Method name definitions ____________________________________________ */
45 #define G_READ "GModelSpatialRadialProfileGauss::read(GXmlElement&)"
46 #define G_WRITE "GModelSpatialRadialProfileGauss::write(GXmlElement&)"
47 
48 /* __ Macros _____________________________________________________________ */
49 
50 /* __ Coding definitions _________________________________________________ */
51 
52 /* __ Debug definitions __________________________________________________ */
53 
54 
55 /*==========================================================================
56  = =
57  = Constructors/destructors =
58  = =
59  ==========================================================================*/
60 
61 /***********************************************************************//**
62  * @brief Void constructor
63  *
64  * Constructs empty radial Gaussian profile
65  ***************************************************************************/
68 {
69  // Initialise members
70  init_members();
71 
72  // Return
73  return;
74 }
75 
76 
77 /***********************************************************************//**
78  * @brief XML constructor
79  *
80  * @param[in] xml XML element.
81  *
82  * Constructs radial Gaussian profile model by extracting information from
83  * an XML element. See the read() method for more information about the
84  * expected structure of the XML element.
85  ***************************************************************************/
88 {
89  // Initialise members
90  init_members();
91 
92  // Read information from XML element
93  read(xml);
94 
95  // Return
96  return;
97 }
98 
99 
100 /***********************************************************************//**
101  * @brief Constructor
102  *
103  * @param[in] dir Sky position of Gaussian.
104  * @param[in] sigma Width of Gaussian (in degrees).
105  *
106  * Constructs a Gaussian spatial model using a sky direction (@p dir) and
107  * a Gaussian width parameter @p sigma in degrees.
108  ***************************************************************************/
110  const double& sigma) :
112 {
113  // Initialise members
114  init_members();
115 
116  // Assign direction and sigma
117  this->dir(dir);
118  this->sigma(sigma);
119 
120  // Return
121  return;
122 }
123 
124 
125 /***********************************************************************//**
126  * @brief Copy constructor
127  *
128  * @param[in] model Radial Gaussian profile model.
129  *
130  * Copies radial Gaussian profile model from another radial profile model.
131  ***************************************************************************/
134 {
135  // Initialise members
136  init_members();
137 
138  // Copy members
139  copy_members(model);
140 
141  // Return
142  return;
143 }
144 
145 
146 /***********************************************************************//**
147  * @brief Destructor
148  *
149  * Destructs radial Gaussian profile model.
150  ***************************************************************************/
152 {
153  // Free members
154  free_members();
155 
156  // Return
157  return;
158 }
159 
160 
161 /*==========================================================================
162  = =
163  = Operators =
164  = =
165  ==========================================================================*/
166 
167 /***********************************************************************//**
168  * @brief Assignment operator
169  *
170  * @param[in] model Radial Gaussian profile model.
171  * @return Radial Gaussian profile model.
172  *
173  * Assigns radial Gaussian profile model.
174  ***************************************************************************/
176 {
177  // Execute only if object is not identical
178  if (this != &model) {
179 
180  // Copy base class members
182 
183  // Free members
184  free_members();
185 
186  // Initialise members
187  init_members();
188 
189  // Copy members
190  copy_members(model);
191 
192  } // endif: object was not identical
193 
194  // Return
195  return *this;
196 }
197 
198 
199 /*==========================================================================
200  = =
201  = Public methods =
202  = =
203  ==========================================================================*/
204 
205 /***********************************************************************//**
206  * @brief Clear radial Gaussian profile model
207  *
208  * Clears radial Gaussian profile model.
209  ***************************************************************************/
211 {
212  // Free class members
213  free_members();
217 
218  // Initialise members
222  init_members();
223 
224  // Return
225  return;
226 }
227 
228 
229 /***********************************************************************//**
230  * @brief Clone radial Gaussian profile model
231  *
232  * @return Pointer to deep copy of radial Gaussian profile model.
233  *
234  * Returns a deep copy of the radial Gaussian profile model.
235  ***************************************************************************/
237 {
238  // Clone radial disk model
239  return new GModelSpatialRadialProfileGauss(*this);
240 }
241 
242 /***********************************************************************//**
243  * @brief Return maximum model radius (in radians)
244  *
245  * @return Maximum model radius (in radians).
246  ***************************************************************************/
248 {
249  // Return value
250  return 0.0 ;
251 }
252 
253 /***********************************************************************//**
254  * @brief Return maximum model radius (in radians)
255  *
256  * @return Maximum model radius (in radians).
257  ***************************************************************************/
259 {
260  // Return value
261  return (m_sigma.value() * gammalib::deg2rad * 5.0);
262 }
263 
264 
265 /***********************************************************************//**
266  * @brief Read model from XML element
267  *
268  * @param[in] xml XML element.
269  *
270  * Reads the Gaussian radial profile model information from an XML element.
271  * The XML element shall have either the format
272  *
273  * <spatialModel type="GaussianProfile">
274  * <parameter name="RA" scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
275  * <parameter name="DEC" scale="1.0" value="22.0145" min="-90" max="90" free="1"/>
276  * <parameter name="Sigma" scale="1.0" value="0.45" min="0.01" max="10" free="1"/>
277  * </spatialModel>
278  *
279  * or
280  *
281  * <spatialModel type="GaussianProfile">
282  * <parameter name="GLON" scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
283  * <parameter name="GLAT" scale="1.0" value="22.0145" min="-90" max="90" free="1"/>
284  * <parameter name="Sigma" scale="1.0" value="0.45" min="0.01" max="10" free="1"/>
285  * </spatialModel>
286  ***************************************************************************/
288 {
289  // Verify number of model parameters
291 
292  // Read Gaussian location
294 
295  // Read Sigma parameter
296  const GXmlElement* par = gammalib::xml_get_par(G_READ, xml, m_sigma.name());
297  m_sigma.read(*par);
298 
299  // Return
300  return;
301 }
302 
303 
304 /***********************************************************************//**
305  * @brief Write model into XML element
306  *
307  * @param[in] xml XML element into which model information is written.
308  *
309  * Writes the Gaussian radial profile model information into an XML element.
310  * The XML element will have the format
311  *
312  * <spatialModel type="GaussianProfile">
313  * <parameter name="RA" scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
314  * <parameter name="DEC" scale="1.0" value="22.0145" min="-90" max="90" free="1"/>
315  * <parameter name="Sigma" scale="1.0" value="0.45" min="0.01" max="10" free="1"/>
316  * </spatialModel>
317  ***************************************************************************/
319 {
320  // Verify model type
322 
323  // Write Gaussian location
325 
326  // Write Sigma parameter
328  m_sigma.write(*par);
329 
330  // Return
331  return;
332 }
333 
334 
335 /***********************************************************************//**
336  * @brief Print information
337  *
338  * @param[in] chatter Chattiness (defaults to NORMAL).
339  * @return String containing model information.
340  ***************************************************************************/
341 std::string GModelSpatialRadialProfileGauss::print(const GChatter& chatter) const
342 {
343  // Initialise result string
344  std::string result;
345 
346  // Continue only if chatter is not silent
347  if (chatter != SILENT) {
348 
349  // Append header
350  result.append("=== GModelSpatialRadialProfileGauss ===");
351 
352  // Append parameters
353  result.append("\n"+gammalib::parformat("Number of parameters"));
354  result.append(gammalib::str(size()));
355  for (int i = 0; i < size(); ++i) {
356  result.append("\n"+m_pars[i]->print(chatter));
357  }
358 
359  } // endif: chatter was not silent
360 
361  // Return result
362  return result;
363 }
364 
365 
366 /*==========================================================================
367  = =
368  = Private methods =
369  = =
370  ==========================================================================*/
371 
372 /***********************************************************************//**
373  * @brief Initialise class members
374  ***************************************************************************/
376 {
377  // Initialise model type
378  m_type = "GaussianProfile";
379 
380  // Initialise Gaussian sigma
381  m_sigma.clear();
382  m_sigma.name("Sigma");
383  m_sigma.unit("deg");
384  m_sigma.value(2.778e-4); // 1 arcsec
385  m_sigma.min(2.778e-4); // 1 arcsec
386  m_sigma.free();
387  m_sigma.scale(1.0);
388  m_sigma.gradient(0.0);
389  m_sigma.has_grad(false); // Radial components never have gradients
390 
391  // Set parameter pointer(s)
392  m_pars.push_back(&m_sigma);
393 
394  // Profile is independent of spatial coordinates of model center
395  m_coord_indep = true;
396 
397  // Return
398  return;
399 }
400 
401 
402 /***********************************************************************//**
403  * @brief Copy class members
404  *
405  * @param[in] model Radial Gaussian model.
406  *
407  * Copies class members from another radial profile model.
408  ***************************************************************************/
410 {
411  // Copy members. We do not have to push back the members on the parameter
412  // stack as this should have been done by init_members() that was called
413  // before. Otherwise we would have sigma twice on the stack.
414  m_type = model.m_type; // Needed to conserve model type
415  m_sigma = model.m_sigma;
416 
417  // Return
418  return;
419 }
420 
421 
422 /***********************************************************************//**
423  * @brief Delete class members
424  ***************************************************************************/
426 {
427  // Return
428  return;
429 }
430 
431 
432 /***********************************************************************//**
433  * @brief Radial profile
434  *
435  * @param[in] theta Angular distance from Gaussian centre (radians).
436  * @return Profile value.
437  ***************************************************************************/
438 double GModelSpatialRadialProfileGauss::profile_value(const double& theta) const
439 {
440  // Compute value
441  double sigma_rad = m_sigma.value() * gammalib::deg2rad;
442  double sigma2 = sigma_rad * sigma_rad;
443  double theta2 = theta * theta;
444  double value = std::exp(-0.5 * theta2 / sigma2) /
445  (gammalib::twopi * sigma2);
446 
447  // Return value
448  return value;
449 }
virtual GModelSpatialRadialProfileGauss & operator=(const GModelSpatialRadialProfileGauss &model)
Assignment operator.
void free_members(void)
Delete class members.
const std::string & name(void) const
Return parameter name.
XML element node class interface definition.
double gradient(void) const
Return parameter gradient.
virtual void write(GXmlElement &xml) const
Write model into XML element.
int size(void) const
Return number of parameters.
void write(GXmlElement &xml) const
Set or update parameter attributes in XML element.
Definition: GModelPar.cpp:351
XML element node class.
Definition: GXmlElement.hpp:48
const GModelSpatialRadialDisk g_radial_disk_seed
Gammalib tools definition.
Interface definition for the spatial model registry class.
double min(void) const
Return parameter minimum boundary.
virtual double theta_min(void) const
Return maximum model radius (in radians)
Radial Gaussian profile model class interface definition.
void init_members(void)
Initialise class members.
void init_members(void)
Initialise class members.
virtual double theta_max(void) const
Return maximum model radius (in radians)
const double & scale(void) const
Return parameter scale.
const double deg2rad
Definition: GMath.hpp:43
bool has_grad(void) const
Signal if parameter gradient is computed analytically.
const GSkyDir & dir(void) const
Return position of radial spatial model.
void free_members(void)
Delete class members.
void free_members(void)
Delete class members.
void free(void)
Free a parameter.
std::vector< GModelPar * > m_pars
Parameter pointers.
GXmlElement * xml_need_par(const std::string &origin, GXmlElement &xml, const std::string &name)
Return pointer to parameter with given name in XML element.
Definition: GTools.cpp:1637
virtual void read(const GXmlElement &xml)
Read model from XML element.
void clear(void)
Clear parameter.
Spatial model registry class definition.
GChatter
Definition: GTypemaps.hpp:33
void init_members(void)
Initialise class members.
std::string type(void) const
Return model type.
virtual void clear(void)
Clear radial Gaussian profile model.
void xml_check_parnum(const std::string &origin, const GXmlElement &xml, const int &number)
Checks number of parameters.
Definition: GTools.cpp:1777
virtual ~GModelSpatialRadialProfileGauss(void)
Destructor.
virtual void write(GXmlElement &xml) const
Write model into XML element.
double sigma(void) const
Return Gaussian sigma.
void free_members(void)
Delete class members.
std::string m_type
Spatial model type.
virtual GModelSpatialRadialProfileGauss * clone(void) const
Clone radial Gaussian profile model.
void read(const GXmlElement &xml)
Extract parameter attributes from XML element.
Definition: GModelPar.cpp:229
double value(void) const
Return parameter value.
virtual GModelSpatialRadialProfile & operator=(const GModelSpatialRadialProfile &model)
Assignment operator.
const std::string & unit(void) const
Return parameter unit.
Exception handler interface definition.
virtual double profile_value(const double &theta) const
Radial profile.
GVector exp(const GVector &vector)
Computes exponential of vector elements.
Definition: GVector.cpp:1232
void init_members(void)
Initialise class members.
const double twopi
Definition: GMath.hpp:36
void copy_members(const GModelSpatialRadialProfileGauss &model)
Copy class members.
std::string parformat(const std::string &s, const int &indent=0)
Convert string in parameter format.
Definition: GTools.cpp:1143
Sky direction class.
Definition: GSkyDir.hpp:62
virtual void read(const GXmlElement &xml)
Read model from XML element.
const GXmlElement * xml_get_par(const std::string &origin, const GXmlElement &xml, const std::string &name)
Return pointer to parameter with given name in XML element.
Definition: GTools.cpp:1689
bool m_coord_indep
True if model independent of sky coordinates.
Mathematical function definitions.
virtual std::string print(const GChatter &chatter=NORMAL) const
Print information.
std::string str(const unsigned short int &value)
Convert unsigned short integer value into string.
Definition: GTools.cpp:489
void xml_check_type(const std::string &origin, GXmlElement &xml, const std::string &type)
Checks the model type.
Definition: GTools.cpp:1819