GammaLib  1.7.0.dev
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GModelSpatialRadialShell.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  * GModelSpatialRadialShell.cpp - Radial shell source model class *
3  * ----------------------------------------------------------------------- *
4  * copyright (C) 2011-2018 by Christoph Deil *
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 GModelSpatialRadialShell.cpp
23  * @brief Radial shell model class implementation
24  * @author Christoph Deil
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"
36 
37 /* __ Constants __________________________________________________________ */
38 
39 /* __ Globals ____________________________________________________________ */
41 const GModelSpatialRegistry g_radial_shell_registry(&g_radial_shell_seed);
42 #if defined(G_LEGACY_XML_FORMAT)
43 const GModelSpatialRadialShell g_radial_shell_legacy_seed(true, "ShellFunction");
44 const GModelSpatialRegistry g_radial_shell_legacy_registry(&g_radial_shell_legacy_seed);
45 #endif
46 
47 /* __ Method name definitions ____________________________________________ */
48 #define G_MC "GModelSpatialRadialShell::mc(GEnergy&, GTime&, GRan&)"
49 #define G_READ "GModelSpatialRadialShell::read(GXmlElement&)"
50 #define G_WRITE "GModelSpatialRadialShell::write(GXmlElement&)"
51 
52 /* __ Macros _____________________________________________________________ */
53 
54 /* __ Coding definitions _________________________________________________ */
55 
56 /* __ Debug definitions __________________________________________________ */
57 //#define G_DEBUG_MC //!< Debug MC method
58 
59 
60 /*==========================================================================
61  = =
62  = Constructors/destructors =
63  = =
64  ==========================================================================*/
65 
66 /***********************************************************************//**
67  * @brief Void constructor
68  *
69  * Constructs empty radial shell model.
70  ***************************************************************************/
72 {
73  // Initialise members
74  init_members();
75 
76  // Return
77  return;
78 }
79 
80 
81 /***********************************************************************//**
82  * @brief Model type constructor
83  *
84  * @param[in] dummy Dummy flag.
85  * @param[in] type Model type.
86  *
87  * Constructs empty radial shell model by specifying a model @p type.
88  ***************************************************************************/
90  const std::string& type) :
92 {
93  // Initialise members
94  init_members();
95 
96  // Set model type
97  m_type = type;
98 
99  // Return
100  return;
101 }
102 
103 
104 /***********************************************************************//**
105  * @brief Shell constructor
106  *
107  * @param[in] dir Sky position of shell centre.
108  * @param[in] radius Inner shell radius (degrees).
109  * @param[in] width Shell width (degrees).
110  *
111  * Constructs the shell model from the shell centre (@p dir), the inner
112  * shell @p radius, and the shell @p width.
113  ***************************************************************************/
115  const double& radius,
116  const double& width) :
118 {
119  // Initialise members
120  init_members();
121 
122  // Assign parameters
123  this->dir(dir);
124  this->radius(radius);
125  this->width(width);
126 
127  // Return
128  return;
129 }
130 
131 
132 /***********************************************************************//**
133  * @brief XML constructor
134  *
135  * @param[in] xml XML element.
136  *
137  * Constructs shell spatial model by extracting information from an XML
138  * element. See the read() method for more information about the expected
139  * structure of the XML element.
140  ***************************************************************************/
143 {
144  // Initialise members
145  init_members();
146 
147  // Read information from XML element
148  read(xml);
149 
150  // Return
151  return;
152 }
153 
154 
155 /***********************************************************************//**
156  * @brief Copy constructor
157  *
158  * @param[in] model Radial shell source model.
159  ***************************************************************************/
161  GModelSpatialRadial(model)
162 {
163  // Initialise members
164  init_members();
165 
166  // Copy members
167  copy_members(model);
168 
169  // Return
170  return;
171 }
172 
173 
174 /***********************************************************************//**
175  * @brief Destructor
176  ***************************************************************************/
178 {
179  // Free members
180  free_members();
181 
182  // Return
183  return;
184 }
185 
186 
187 /*==========================================================================
188  = =
189  = Operators =
190  = =
191  ==========================================================================*/
192 
193 /***********************************************************************//**
194  * @brief Assignment operator
195  *
196  * @param[in] model Radial shell source model.
197  * @return Radial shell source model.
198  ***************************************************************************/
200 {
201  // Execute only if object is not identical
202  if (this != &model) {
203 
204  // Copy base class members
205  this->GModelSpatialRadial::operator=(model);
206 
207  // Free members
208  free_members();
209 
210  // Initialise members
211  init_members();
212 
213  // Copy members
214  copy_members(model);
215 
216  } // endif: object was not identical
217 
218  // Return
219  return *this;
220 }
221 
222 
223 /*==========================================================================
224  = =
225  = Public methods =
226  = =
227  ==========================================================================*/
228 
229 /***********************************************************************//**
230  * @brief Clear radial shell model
231  ***************************************************************************/
233 {
234  // Free class members (base and derived classes, derived class first)
235  free_members();
238 
239  // Initialise members
242  init_members();
243 
244  // Return
245  return;
246 }
247 
248 
249 /***********************************************************************//**
250  * @brief Clone radial shell model
251  *
252  * @return Pointer to deep copy of radial shell model.
253  ***************************************************************************/
255 {
256  // Clone radial shell model
257  return new GModelSpatialRadialShell(*this);
258 }
259 
260 
261 /***********************************************************************//**
262  * @brief Evaluate function (in units of sr^-1)
263  *
264  * @param[in] theta Angular distance from shell centre (radians).
265  * @param[in] energy Photon energy.
266  * @param[in] time Photon arrival time.
267  * @param[in] gradients Compute gradients?
268  * @return Model value.
269  *
270  * Evaluates the spatial part for a shell source model. The shell source
271  * model is a radial function \f$f(\theta)\f$, where \f$\theta\f$ is the
272  * angular separation between shell centre and the actual location.
273  *
274  * The function is given by
275  * \f[
276  * S_{\rm p}(\vec{p} | E, t) = {\tt m\_norm} \left \{
277  * \begin{array}{l l}
278  * \displaystyle
279  * \sqrt{ \sin^2 \theta_{\rm out} - \sin^2 \theta } -
280  * \sqrt{ \sin^2 \theta_{\rm in} - \sin^2 \theta }
281  * & \mbox{if $\theta \le \theta_{\rm in}$} \\
282  * \\
283  * \displaystyle
284  * \sqrt{ \sin^2 \theta_{\rm out} - \sin^2 \theta }
285  * & \mbox{if $\theta_{\rm in} < \theta \le \theta_{\rm out}$} \\
286  * \\
287  * \displaystyle
288  * 0 & \mbox{if $\theta > \theta_{\rm out}$}
289  * \end{array}
290  * \right .
291  * \f]
292  *
293  * Here, \f$\theta_{\rm in}\f$ and \f$\theta_{\rm out}\f$ are the shell inner
294  * and outer radius.
295  *
296  * The method will not compute analytical parameter gradients, even if the
297  * @p gradients argument is set to true. Radial disk parameter gradients
298  * need to be computed numerically.
299  ***************************************************************************/
300 double GModelSpatialRadialShell::eval(const double& theta,
301  const GEnergy& energy,
302  const GTime& time,
303  const bool& gradients) const
304 {
305  // Update precomputation cache
306  update();
307 
308  // Set x appropriately
309  double x = std::sin(theta);
310  x *= x;
311 
312  // Compute result
313  double result = 0.0;
314  if (x < m_x_out) {
315  result = std::sqrt(m_x_out - x);
316  if (x < m_x_in) {
317  result -= std::sqrt(m_x_in - x);
318  }
319  result *= m_norm;
320  }
321 
322  // Compile option: Check for NaN/Inf
323  #if defined(G_NAN_CHECK)
324  if (gammalib::is_notanumber(result) || gammalib::is_infinite(result)) {
325  std::cout << "*** ERROR: GModelSpatialRadialShell::eval";
326  std::cout << "(theta=" << theta << "): NaN/Inf encountered";
327  std::cout << " (result=" << result;
328  std::cout << ", m_norm=" << m_norm;
329  std::cout << ", x=" << x;
330  std::cout << ", m_x_out=" << m_x_out;
331  std::cout << ", m_x_in=" << m_x_in;
332  std::cout << ")" << std::endl;
333  }
334  #endif
335 
336  // Return normalised value
337  return result;
338 }
339 
340 
341 /***********************************************************************//**
342  * @brief Returns MC sky direction
343  *
344  * @param[in] energy Photon energy.
345  * @param[in] time Photon arrival time.
346  * @param[in,out] ran Random number generator.
347  * @return Sky direction.
348  *
349  * @exception GException::invalid_return_value
350  * Invalid u_max value
351  *
352  * Draws an arbitrary sky position from the 2D shell distribution.
353  ***************************************************************************/
355  const GTime& time,
356  GRan& ran) const
357 {
358  // Update precomputation cache
359  update();
360 
361  // Initialise simulation
362  #if defined(G_DEBUG_MC)
363  int n_samples = 0;
364  #endif
365  double u_max = m_norm * m_x_out;
366  double value = 0.0;
367  double u = 1.0;
368  double theta = 0.0;
369 
370  // Throw an exception if the maximum value is zero
371  if (u_max <= 0.0) {
372  std::string msg = "Non positive maximum radial shell model value.";
374  }
375 
376  // Simulate offset from photon arrival direction using a rejection method
377  do {
378  theta = ran.uniform() * m_theta_out;
379  value = eval(theta, energy, time) * std::sin(theta);
380  u = ran.uniform() * u_max;
381  #if defined(G_DEBUG_MC)
382  n_samples++;
383  #endif
384  } while (u > value);
385  #if defined(G_DEBUG_MC)
386  std::cout << "#=" << n_samples << " ";
387  #endif
388 
389  // Simulate azimuth angle
390  double phi = 360.0 * ran.uniform();
391 
392  // Rotate pointing direction by offset and azimuth angle
393  GSkyDir sky_dir = dir();
394  sky_dir.rotate_deg(phi, theta * gammalib::rad2deg);
395 
396  // Return sky direction
397  return sky_dir;
398 }
399 
400 
401 /***********************************************************************//**
402  * @brief Checks where model contains specified sky direction
403  *
404  * @param[in] dir Sky direction.
405  * @param[in] margin Margin to be added to sky direction (degrees)
406  * @return True if the model contains the sky direction.
407  *
408  * Signals whether a sky direction is contained in the radial shell model.
409  ***************************************************************************/
411  const double& margin) const
412 {
413  // Compute distance to centre (radians)
414  double distance = dir.dist(this->dir());
415 
416  // Return flag
417  return (distance <= theta_max() + margin*gammalib::deg2rad);
418 }
419 
420 
421 /***********************************************************************//**
422  * @brief Return maximum model radius (in radians)
423  *
424  * @return Maximum model radius (in radians).
425  ***************************************************************************/
427 {
428  // Return value
429  return ((radius()+width()) * gammalib::deg2rad);
430 }
431 
432 
433 /***********************************************************************//**
434  * @brief Read model from XML element
435  *
436  * @param[in] xml XML element.
437  *
438  * @exception GException::model_invalid_parnum
439  * Invalid number of model parameters found in XML element.
440  * @exception GException::model_invalid_parnames
441  * Invalid model parameter names found in XML element.
442  *
443  * Reads the radial shell model information from an XML element. The XML
444  * element shall have either the format
445  *
446  * <spatialModel type="RadialShell">
447  * <parameter name="RA" scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
448  * <parameter name="DEC" scale="1.0" value="22.0145" min="-90" max="90" free="1"/>
449  * <parameter name="Radius" scale="1.0" value="0.30" min="0.01" max="10" free="1"/>
450  * <parameter name="Width" scale="1.0" value="0.10" min="0.01" max="10" free="1"/>
451  * </spatialModel>
452  *
453  * or
454  *
455  * <spatialModel type="RadialShell">
456  * <parameter name="GLON" scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
457  * <parameter name="GLAT" scale="1.0" value="22.0145" min="-90" max="90" free="1"/>
458  * <parameter name="Radius" scale="1.0" value="0.30" min="0.01" max="10" free="1"/>
459  * <parameter name="Width" scale="1.0" value="0.10" min="0.01" max="10" free="1"/>
460  * </spatialModel>
461  *
462  * @todo Implement tests of the radius and radius boundary and the width and
463  * width boundary. The radius and radius boundary should be >=0, the
464  * width and width boundary should be >0.
465  ***************************************************************************/
467 {
468  // Determine number of parameter nodes in XML element
469  int npars = xml.elements("parameter");
470 
471  // Verify that XML element has exactly 4 parameters
472  if (xml.elements() != 4 || npars != 4) {
474  "Shell model requires exactly 4 parameters.");
475  }
476 
477  // Read shell location
479 
480  // Get parameters
483 
484  // Read parameters
485  m_radius.read(*radius);
486  m_width.read(*width);
487 
488  // Return
489  return;
490 }
491 
492 
493 /***********************************************************************//**
494  * @brief Write model into XML element
495  *
496  * @param[in] xml XML element into which model information is written.
497  *
498  * @exception GException::model_invalid_spatial
499  * Existing XML element is not of type 'GaussFunction'
500  * @exception GException::model_invalid_parnum
501  * Invalid number of model parameters found in XML element.
502  * @exception GException::model_invalid_parnames
503  * Invalid model parameter names found in XML element.
504  *
505  * Writes the radial shell model information into an XML element. The XML
506  * element will have the format
507  *
508  * <spatialModel type="RadialShell">
509  * <parameter name="RA" scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
510  * <parameter name="DEC" scale="1.0" value="22.0145" min="-90" max="90" free="1"/>
511  * <parameter name="Radius" scale="1.0" value="0.30" min="0.01" max="10" free="1"/>
512  * <parameter name="Width" scale="1.0" value="0.10" min="0.01" max="10" free="1"/>
513  * </spatialModel>
514  *
515  ***************************************************************************/
517 {
518  // Write shell location
520 
521  // Get or create parameters
524 
525  // Write parameters
526  m_radius.write(*radius);
527  m_width.write(*width);
528 
529  // Return
530  return;
531 }
532 
533 
534 /***********************************************************************//**
535  * @brief Print information
536  *
537  * @param[in] chatter Chattiness (defaults to NORMAL).
538  * @return String containing model information.
539  ***************************************************************************/
540 std::string GModelSpatialRadialShell::print(const GChatter& chatter) const
541 {
542  // Initialise result string
543  std::string result;
544 
545  // Continue only if chatter is not silent
546  if (chatter != SILENT) {
547 
548  // Append header
549  result.append("=== GModelSpatialRadialShell ===");
550 
551  // Append parameters
552  result.append("\n"+gammalib::parformat("Number of parameters"));
553  result.append(gammalib::str(size()));
554  for (int i = 0; i < size(); ++i) {
555  result.append("\n"+m_pars[i]->print(chatter));
556  }
557 
558  } // endif: chatter was not silent
559 
560  // Return result
561  return result;
562 }
563 
564 
565 /*==========================================================================
566  = =
567  = Private methods =
568  = =
569  ==========================================================================*/
570 
571 /***********************************************************************//**
572  * @brief Initialise class members
573  ***************************************************************************/
575 {
576  // Initialise model type
577  m_type = "RadialShell";
578 
579  // Initialise Radius
580  m_radius.clear();
581  m_radius.name("Radius");
582  m_radius.unit("deg");
583  m_radius.value(0.0);
584  m_radius.min(0.0);
585  m_radius.free();
586  m_radius.scale(1.0);
587  m_radius.gradient(0.0);
588  m_radius.has_grad(false); // Radial components never have gradients
589 
590  // Initialise Width
591  m_width.clear();
592  m_width.name("Width");
593  m_width.unit("deg");
594  m_width.value(2.778e-4);
595  m_width.min(2.778e-4); // 1 arcsec
596  m_width.free();
597  m_width.scale(1.0);
598  m_width.gradient(0.0);
599  m_width.has_grad(false); // Radial components never have gradients
600 
601  // Set parameter pointer(s)
602  m_pars.push_back(&m_radius);
603  m_pars.push_back(&m_width);
604 
605  // Initialise other members
606  m_region.clear();
607 
608  // Initialise precomputation cache. Note that zero values flag
609  // uninitialised as a zero radius and width shell is not meaningful
610  m_last_radius = 0.0;
611  m_last_width = 0.0;
612  m_theta_in = 0.0;
613  m_x_in = 0.0;
614  m_theta_out = 0.0;
615  m_x_out = 0.0;
616  m_norm = 0.0;
617 
618  // Return
619  return;
620 }
621 
622 
623 /***********************************************************************//**
624  * @brief Copy class members
625  *
626  * @param[in] model Radial shell source model.
627  *
628  * We do not have to push back the members on the parameter stack as this
629  * should have been done by init_members() that was called before. Otherwise
630  * we would have the radius and width twice on the stack.
631  ***************************************************************************/
633 {
634  // Copy members
635  m_type = model.m_type;
636  m_radius = model.m_radius;
637  m_width = model.m_width;
638  m_region = model.m_region;
639 
640  // Copy precomputation cache
642  m_last_width = model.m_last_width;
643  m_theta_in = model.m_theta_in;
644  m_x_in = model.m_x_in;
645  m_theta_out = model.m_theta_out;
646  m_x_out = model.m_x_out;
647  m_norm = model.m_norm;
648 
649  // Return
650  return;
651 }
652 
653 
654 /***********************************************************************//**
655  * @brief Delete class members
656  ***************************************************************************/
658 {
659  // Return
660  return;
661 }
662 
663 
664 /***********************************************************************//**
665  * @brief Update precomputation cache
666  *
667  * Performs precomputations that are valid for a given pair of radius and
668  * width values. The following members are set by this method:
669  * m_last_radius, m_last_width, m_theta_in, m_theta_out, m_x_in, m_x_out
670  * and m_norm.
671  *
672  * m_theta_in contains the inner shell radius in radians, while m_theta_out
673  * contains the outer shell radius in radians.
674  *
675  * \f${\tt m\_x\_in} = \sin^2 {\tt m\_theta\_in}\f$,
676  * \f${\tt m\_x\_out} = \sin^2 {\tt m\_theta\_out}\f$, and
677  * \f[{\tt m\_norm} = \frac{1}{2 \pi}
678  * \frac{\sqrt{1-\cos 2 {\tt m\_theta\_out}} -
679  * \sqrt{1-\cos 2 {\tt m\_theta\_in}}}{2 \sqrt{2}} +
680  * \frac{1+\cos 2 {\tt m\_theta\_out}}{4} \ln \left(
681  * \frac{\sqrt{2} \cos {\tt m\_theta\_out}}
682  * {\sqrt{2} + \sqrt{1 - \cos 2 {\tt m\_theta\_out}}} \right) -
683  * \frac{1+\cos 2 {\tt m\_theta\_in}}{4} \ln \left(
684  * \frac{\sqrt{2} \cos {\tt m\_theta\_in}}
685  * {\sqrt{2} + \sqrt{1 - \cos 2 {\tt m\_theta\_in}}} \right)\f]
686  ***************************************************************************/
688 {
689  // Set constants
690  //const double c1 = gammalib::twopi / 3.0;
691  const double c2 = 1.0 / (2.0 * gammalib::sqrt_two);
692 
693  // Update if radius or width have changed
694  if (m_last_radius != radius() || m_last_width != width()) {
695 
696  // Store last values
697  m_last_radius = radius();
698  m_last_width = width();
699 
700  // Perform precomputations
703  double sin_theta_in = std::sin(m_theta_in);
704  double sin_theta_out = std::sin(m_theta_out);
705  double term1 = (f1(m_theta_out) - f1(m_theta_in)) * c2;
706  double term2 = f2(m_theta_out);
707  double term3 = f2(m_theta_in);
708  double denom = gammalib::twopi * (term1 + term2 - term3);
709  m_norm = (denom > 0.0) ? 1.0 / denom : 0.0;
710  m_x_in = sin_theta_in * sin_theta_in;
711  m_x_out = sin_theta_out * sin_theta_out;
712 
713  } // endif: update required
714 
715  // Compile option: Check for NaN/Inf
716  #if defined(G_NAN_CHECK)
718  std::cout << "*** ERROR: GModelSpatialRadialShell::update:";
719  std::cout << " NaN/Inf encountered";
720  std::cout << " (m_norm=" << m_norm;
721  std::cout << ", radius=" << radius();
722  std::cout << ", width=" << width();
723  std::cout << ", m_theta_in=" << m_theta_in;
724  std::cout << ", m_theta_out=" << m_theta_out;
725  std::cout << ", term1=" << (f1(m_theta_out) - f1(m_theta_in)) * c2;
726  std::cout << ", term2=" << f2(m_theta_out);
727  std::cout << ", term3=" << f2(m_theta_in);
728  std::cout << ")" << std::endl;
729  }
730  #endif
731 
732  // Return
733  return;
734 }
735 
736 
737 /***********************************************************************//**
738  * @brief Return function 1 value needed for precomputation
739  *
740  * Computes \f$f1(x) = \sqrt{1 - \cos 2 x}\f$.
741  ***************************************************************************/
743 {
744  // Compute value
745  double f1 = std::sqrt(1.0 - std::cos(2.0 * x));
746 
747  // Return value
748  return f1;
749 }
750 
751 
752 /***********************************************************************//**
753  * @brief Return function 2 value needed for precomputation
754  *
755  * Compute
756  * \f[f2(x) = \frac{1+\cos 2x}{4}
757  * \ln \left( \frac{\sqrt{2} \cos x}{\sqrt{2} + \sqrt{ 1 - \cos 2 x}}
758  * \right)\f].
759  ***************************************************************************/
761 {
762  // Compute value
763  double t1 = (1.0 + std::cos(2.0*x)) / 4.0;
764  double t2 = gammalib::sqrt_two * std::cos(x);
765  double t3 = gammalib::sqrt_two + f1(x);
766  double f2 = t1 * std::log(t2 / t3);
767 
768  // Return value
769  return f2;
770 }
771 
772 
773 /***********************************************************************//**
774  * @brief Set boundary sky region
775  ***************************************************************************/
777 {
778  // Set sky region centre to Gaussian centre
780 
781  // Set sky region radius to sum of shell radius and width
783 
784  // Return
785  return;
786 }
double width(void) const
Return shell width.
GModelPar m_width
Shell thickness (deg)
virtual GModelSpatialRadialShell & operator=(const GModelSpatialRadialShell &model)
Assignment operator.
Radial shell model class interface definition.
const std::string & name(void) const
Return parameter name.
double gradient(void) const
Return parameter gradient.
int size(void) const
Return number of parameters.
double m_last_width
Last shell width (deg)
void init_members(void)
Initialise class members.
void write(GXmlElement &xml) const
Set or update parameter attributes in XML element.
Definition: GModelPar.cpp:353
GVector cos(const GVector &vector)
Computes cosine of vector elements.
Definition: GVector.cpp:1100
const GModelSpatialRadialShell g_radial_shell_seed
XML element node class.
Definition: GXmlElement.hpp:47
void set_region(void) const
Set boundary sky region.
virtual void read(const GXmlElement &xml)
Read model from XML element.
Random number generator class.
Definition: GRan.hpp:44
#define G_MC
virtual int elements(void) const
Return number of GXMLElement children of node.
Definition: GXmlNode.cpp:580
Time class.
Definition: GTime.hpp:54
Gammalib tools definition.
GModelPar m_dec
Declination (deg)
#define G_READ
const double sqrt_two
Definition: GMath.hpp:52
Interface definition for the spatial model registry class.
double min(void) const
Return parameter minimum boundary.
void init_members(void)
Initialise class members.
static double f2(double x)
Return function 2 value needed for precomputation.
bool is_notanumber(const double &x)
Signal if argument is not a number.
Definition: GTools.hpp:184
const double & radius(void) const
Return circular region radius (in degrees)
double m_x_in
sin(m_theta_in)^2
const double & scale(void) const
Return parameter scale.
bool is_infinite(const double &x)
Signal if argument is infinite.
Definition: GTools.hpp:167
double radius(void) const
Return shell radius.
const double deg2rad
Definition: GMath.hpp:43
void free_members(void)
Delete class members.
GVector sqrt(const GVector &vector)
Computes square root of vector elements.
Definition: GVector.cpp:1268
bool has_grad(void) const
Signal if parameter gradient is computed analytically.
void copy_members(const GModelSpatialRadialShell &model)
Copy class members.
void update(void) const
Update precomputation cache.
void free_members(void)
Delete class members.
void free(void)
Free a parameter.
std::vector< GModelPar * > m_pars
Parameter pointers.
void rotate_deg(const double &phi, const double &theta)
Rotate sky direction by zenith and azimuth angle.
Definition: GSkyDir.cpp:321
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:1513
virtual void read(const GXmlElement &xml)
Read model from XML element.
virtual std::string type(void) const
Return model type.
GVector log(const GVector &vector)
Computes natural logarithm of vector elements.
Definition: GVector.cpp:1184
double uniform(void)
Returns random double precision floating value in range 0 to 1.
Definition: GRan.cpp:242
void clear(void)
Clear parameter.
virtual double theta_max(void) const
Return maximum model radius (in radians)
virtual void write(GXmlElement &xml) const
Write model into XML element.
virtual ~GModelSpatialRadialShell(void)
Destructor.
Spatial model registry class definition.
virtual GModelSpatialRadialShell * clone(void) const
Clone radial shell model.
GChatter
Definition: GTypemaps.hpp:33
double m_last_radius
Last shell radius (deg)
virtual bool contains(const GSkyDir &dir, const double &margin=0.0) const
Checks where model contains specified sky direction.
virtual GSkyDir mc(const GEnergy &energy, const GTime &time, GRan &ran) const
Returns MC sky direction.
GModelPar m_radius
Inner shell radius (deg)
#define G_WRITE
GModelPar m_ra
Right Ascension (deg)
virtual void write(GXmlElement &xml) const
Write model into XML element.
GModelSpatialRadialShell(void)
Void constructor.
virtual void clear(void)
Clear radial shell model.
double m_theta_in
Inner shell radius (rad)
virtual GModelSpatialRadial & operator=(const GModelSpatialRadial &model)
Assignment operator.
const GSkyDir & centre(void) const
Return circular region centre.
void free_members(void)
Delete class members.
GSkyDir dir(void) const
Return position of radial spatial model.
virtual double eval(const double &theta, const GEnergy &energy, const GTime &time, const bool &gradients=false) const
Evaluate function (in units of sr^-1)
void read(const GXmlElement &xml)
Extract parameter attributes from XML element.
Definition: GModelPar.cpp:229
double value(void) const
Return parameter value.
const std::string & unit(void) const
Return parameter unit.
double m_theta_out
Outer shell radius (rad)
GVector sin(const GVector &vector)
Computes sine of vector elements.
Definition: GVector.cpp:1226
Exception handler interface definition.
double dist(const GSkyDir &dir) const
Compute angular distance between sky directions in radians.
Definition: GSkyDir.hpp:265
static double f1(double x)
Return function 1 value needed for precomputation.
Abstract radial spatial model base class.
void init_members(void)
Initialise class members.
const double twopi
Definition: GMath.hpp:36
double m_norm
Shell normalization.
const double rad2deg
Definition: GMath.hpp:44
double m_x_out
sin(m_theta_out)^2
std::string parformat(const std::string &s, const int &indent=0)
Convert string in parameter format.
Definition: GTools.cpp:1022
Sky direction class.
Definition: GSkyDir.hpp:62
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:1562
GSkyRegionCircle m_region
Bounding circle.
Mathematical function definitions.
void clear(void)
Clear instance.
Class that handles energies in a unit independent way.
Definition: GEnergy.hpp:48
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:413