GeographicLib  1.43
CircularEngine.hpp
Go to the documentation of this file.
1 /**
2  * \file CircularEngine.hpp
3  * \brief Header for GeographicLib::CircularEngine class
4  *
5  * Copyright (c) Charles Karney (2011) <charles@karney.com> and licensed under
6  * the MIT/X11 License. For more information, see
7  * http://geographiclib.sourceforge.net/
8  **********************************************************************/
9 
10 #if !defined(GEOGRAPHICLIB_CIRCULARENGINE_HPP)
11 #define GEOGRAPHICLIB_CIRCULARENGINE_HPP 1
12 
13 #include <vector>
16 
17 #if defined(_MSC_VER)
18 // Squelch warnings about dll vs vector
19 # pragma warning (push)
20 # pragma warning (disable: 4251)
21 #endif
22 
23 namespace GeographicLib {
24 
25  /**
26  * \brief Spherical harmonic sums for a circle
27  *
28  * The class is a companion to SphericalEngine. If the results of a
29  * spherical harmonic sum are needed for several points on a circle of
30  * constant latitude \e lat and height \e h, then SphericalEngine::Circle can
31  * compute the inner sum, which is independent of longitude \e lon, and
32  * produce a CircularEngine object. CircularEngine::operator()() can
33  * then be used to perform the outer sum for particular vales of \e lon.
34  * This can lead to substantial improvements in computational speed for high
35  * degree sum (approximately by a factor of \e N / 2 where \e N is the
36  * maximum degree).
37  *
38  * CircularEngine is tightly linked to the internals of SphericalEngine. For
39  * that reason, the constructor for this class is private. Use
40  * SphericalHarmonic::Circle, SphericalHarmonic1::Circle, and
41  * SphericalHarmonic2::Circle to create instances of this class.
42  *
43  * CircularEngine stores the coefficients needed to allow the summation over
44  * order to be performed in 2 or 6 vectors of length \e M + 1 (depending on
45  * whether gradients are to be calculated). For this reason the constructor
46  * may throw a std::bad_alloc exception.
47  *
48  * Example of use:
49  * \include example-CircularEngine.cpp
50  **********************************************************************/
51 
53  private:
54  typedef Math::real real;
55  enum normalization {
56  FULL = SphericalEngine::FULL,
57  SCHMIDT = SphericalEngine::SCHMIDT,
58  };
59  int _M;
60  bool _gradp;
61  unsigned _norm;
62  real _a, _r, _u, _t;
63  std::vector<real> _wc, _ws, _wrc, _wrs, _wtc, _wts;
64  real _q, _uq, _uq2;
65 
66  Math::real Value(bool gradp, real cl, real sl,
67  real& gradx, real& grady, real& gradz) const;
68 
69  static inline void cossin(real x, real& cosx, real& sinx) {
70  using std::abs; using std::cos; using std::sin;
71  x = x >= 180 ? x - 360 : (x < -180 ? x + 360 : x);
72  real xi = x * Math::degree();
73  cosx = abs(x) == 90 ? 0 : cos(xi);
74  sinx = x == -180 ? 0 : sin(xi);
75  }
76 
77  friend class SphericalEngine;
78  friend class GravityCircle; // Access to cossin
79  friend class MagneticCircle; // Access to cossin
80  CircularEngine(int M, bool gradp, unsigned norm,
81  real a, real r, real u, real t)
82  : _M(M)
83  , _gradp(gradp)
84  , _norm(norm)
85  , _a(a)
86  , _r(r)
87  , _u(u)
88  , _t(t)
89  , _wc(std::vector<real>(_M + 1, 0))
90  , _ws(std::vector<real>(_M + 1, 0))
91  , _wrc(std::vector<real>(_gradp ? _M + 1 : 0, 0))
92  , _wrs(std::vector<real>(_gradp ? _M + 1 : 0, 0))
93  , _wtc(std::vector<real>(_gradp ? _M + 1 : 0, 0))
94  , _wts(std::vector<real>(_gradp ? _M + 1 : 0, 0))
95  {
96  _q = _a / _r;
97  _uq = _u * _q;
98  _uq2 = Math::sq(_uq);
99  }
100 
101  void SetCoeff(int m, real wc, real ws)
102  { _wc[m] = wc; _ws[m] = ws; }
103 
104  void SetCoeff(int m, real wc, real ws,
105  real wrc, real wrs, real wtc, real wts) {
106  _wc[m] = wc; _ws[m] = ws;
107  if (_gradp) {
108  _wrc[m] = wrc; _wrs[m] = wrs;
109  _wtc[m] = wtc; _wts[m] = wts;
110  }
111  }
112 
113  public:
114 
115  /**
116  * A default constructor. CircularEngine::operator()() on the resulting
117  * object returns zero. The resulting object can be assigned to the result
118  * of SphericalHarmonic::Circle.
119  **********************************************************************/
121  : _M(-1)
122  , _gradp(true)
123  , _u(0)
124  , _t(1)
125  {}
126 
127  /**
128  * Evaluate the sum for a particular longitude given in terms of its
129  * cosine and sine.
130  *
131  * @param[in] coslon the cosine of the longitude.
132  * @param[in] sinlon the sine of the longitude.
133  * @return \e V the value of the sum.
134  *
135  * The arguments must satisfy <i>coslon</i><sup>2</sup> +
136  * <i>sinlon</i><sup>2</sup> = 1.
137  **********************************************************************/
138  Math::real operator()(real coslon, real sinlon) const {
139  real dummy;
140  return Value(false, coslon, sinlon, dummy, dummy, dummy);
141  }
142 
143  /**
144  * Evaluate the sum for a particular longitude.
145  *
146  * @param[in] lon the longitude (degrees).
147  * @return \e V the value of the sum.
148  **********************************************************************/
149  Math::real operator()(real lon) const {
150  real coslon, sinlon;
151  cossin(lon, coslon, sinlon);
152  return (*this)(coslon, sinlon);
153  }
154 
155  /**
156  * Evaluate the sum and its gradient for a particular longitude given in
157  * terms of its cosine and sine.
158  *
159  * @param[in] coslon the cosine of the longitude.
160  * @param[in] sinlon the sine of the longitude.
161  * @param[out] gradx \e x component of the gradient.
162  * @param[out] grady \e y component of the gradient.
163  * @param[out] gradz \e z component of the gradient.
164  * @return \e V the value of the sum.
165  *
166  * The gradients will only be computed if the CircularEngine object was
167  * created with this capability (e.g., via \e gradp = true in
168  * SphericalHarmonic::Circle). If not, \e gradx, etc., will not be
169  * touched. The arguments must satisfy <i>coslon</i><sup>2</sup> +
170  * <i>sinlon</i><sup>2</sup> = 1.
171  **********************************************************************/
172  Math::real operator()(real coslon, real sinlon,
173  real& gradx, real& grady, real& gradz) const {
174  return Value(true, coslon, sinlon, gradx, grady, gradz);
175  }
176 
177  /**
178  * Evaluate the sum and its gradient for a particular longitude.
179  *
180  * @param[in] lon the longitude (degrees).
181  * @param[out] gradx \e x component of the gradient.
182  * @param[out] grady \e y component of the gradient.
183  * @param[out] gradz \e z component of the gradient.
184  * @return \e V the value of the sum.
185  *
186  * The gradients will only be computed if the CircularEngine object was
187  * created with this capability (e.g., via \e gradp = true in
188  * SphericalHarmonic::Circle). If not, \e gradx, etc., will not be
189  * touched.
190  **********************************************************************/
192  real& gradx, real& grady, real& gradz) const {
193  real coslon, sinlon;
194  cossin(lon, coslon, sinlon);
195  return (*this)(coslon, sinlon, gradx, grady, gradz);
196  }
197  };
198 
199 } // namespace GeographicLib
200 
201 #if defined(_MSC_VER)
202 # pragma warning (pop)
203 #endif
204 
205 #endif // GEOGRAPHICLIB_CIRCULARENGINE_HPP
#define GEOGRAPHICLIB_EXPORT
Definition: Constants.hpp:90
The evaluation engine for SphericalHarmonic.
GeographicLib::Math::real real
Definition: GeodSolve.cpp:32
Geomagnetic field on a circle of latitude.
static T sq(T x)
Definition: Math.hpp:244
Namespace for GeographicLib.
Definition: Accumulator.cpp:12
static T degree()
Definition: Math.hpp:228
Math::real operator()(real coslon, real sinlon) const
Spherical harmonic sums for a circle.
Header for GeographicLib::Constants class.
Math::real operator()(real lon, real &gradx, real &grady, real &gradz) const
Math::real operator()(real coslon, real sinlon, real &gradx, real &grady, real &gradz) const
Header for GeographicLib::SphericalEngine class.
Math::real operator()(real lon) const
Gravity on a circle of latitude.