GeographicLib  1.35
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 throw();
68 
69  static inline void cossin(real x, real& cosx, real& sinx) {
70  x = x >= 180 ? x - 360 : (x < -180 ? x + 360 : x);
71  real xi = x * Math::degree<real>();
72  cosx = std::abs(x) == 90 ? 0 : cos(xi);
73  sinx = x == -180 ? 0 : sin(xi);
74  }
75 
76  friend class SphericalEngine;
77  friend class GravityCircle; // Access to cossin
78  friend class MagneticCircle; // Access to cossin
79  CircularEngine(int M, bool gradp, unsigned norm,
80  real a, real r, real u, real t)
81  : _M(M)
82  , _gradp(gradp)
83  , _norm(norm)
84  , _a(a)
85  , _r(r)
86  , _u(u)
87  , _t(t)
88  , _wc(std::vector<real>(_M + 1, 0))
89  , _ws(std::vector<real>(_M + 1, 0))
90  , _wrc(std::vector<real>(_gradp ? _M + 1 : 0, 0))
91  , _wrs(std::vector<real>(_gradp ? _M + 1 : 0, 0))
92  , _wtc(std::vector<real>(_gradp ? _M + 1 : 0, 0))
93  , _wts(std::vector<real>(_gradp ? _M + 1 : 0, 0))
94  {
95  _q = _a / _r;
96  _uq = _u * _q;
97  _uq2 = Math::sq(_uq);
98  }
99 
100  void SetCoeff(int m, real wc, real ws)
101  { _wc[m] = wc; _ws[m] = ws; }
102 
103  void SetCoeff(int m, real wc, real ws,
104  real wrc, real wrs, real wtc, real wts) {
105  _wc[m] = wc; _ws[m] = ws;
106  if (_gradp) {
107  _wrc[m] = wrc; _wrs[m] = wrs;
108  _wtc[m] = wtc; _wts[m] = wts;
109  }
110  }
111 
112  public:
113 
114  /**
115  * A default constructor. CircularEngine::operator()() on the resulting
116  * object returns zero. The resulting object can be assigned to the result
117  * of SphericalHarmonic::Circle.
118  **********************************************************************/
120  : _M(-1)
121  , _gradp(true)
122  , _u(0)
123  , _t(1)
124  {}
125 
126  /**
127  * Evaluate the sum for a particular longitude given in terms of its
128  * cosine and sine.
129  *
130  * @param[in] coslon the cosine of the longitude.
131  * @param[in] sinlon the sine of the longitude.
132  * @return \e V the value of the sum.
133  *
134  * The arguments must satisfy <i>coslon</i><sup>2</sup> +
135  * <i>sinlon</i><sup>2</sup> = 1.
136  **********************************************************************/
137  Math::real operator()(real coslon, real sinlon) const throw() {
138  real dummy;
139  return Value(false, coslon, sinlon, dummy, dummy, dummy);
140  }
141 
142  /**
143  * Evaluate the sum for a particular longitude.
144  *
145  * @param[in] lon the longitude (degrees).
146  * @return \e V the value of the sum.
147  **********************************************************************/
148  Math::real operator()(real lon) const throw() {
149  real coslon, sinlon;
150  cossin(lon, coslon, sinlon);
151  return (*this)(coslon, sinlon);
152  }
153 
154  /**
155  * Evaluate the sum and its gradient for a particular longitude given in
156  * terms of its cosine and sine.
157  *
158  * @param[in] coslon the cosine of the longitude.
159  * @param[in] sinlon the sine of the longitude.
160  * @param[out] gradx \e x component of the gradient.
161  * @param[out] grady \e y component of the gradient.
162  * @param[out] gradz \e z component of the gradient.
163  * @return \e V the value of the sum.
164  *
165  * The gradients will only be computed if the CircularEngine object was
166  * created with this capability (e.g., via \e gradp = true in
167  * SphericalHarmonic::Circle). If not, \e gradx, etc., will not be
168  * touched. The arguments must satisfy <i>coslon</i><sup>2</sup> +
169  * <i>sinlon</i><sup>2</sup> = 1.
170  **********************************************************************/
171  Math::real operator()(real coslon, real sinlon,
172  real& gradx, real& grady, real& gradz) const throw() {
173  return Value(true, coslon, sinlon, gradx, grady, gradz);
174  }
175 
176  /**
177  * Evaluate the sum and its gradient for a particular longitude.
178  *
179  * @param[in] lon the longitude (degrees).
180  * @param[out] gradx \e x component of the gradient.
181  * @param[out] grady \e y component of the gradient.
182  * @param[out] gradz \e z component of the gradient.
183  * @return \e V the value of the sum.
184  *
185  * The gradients will only be computed if the CircularEngine object was
186  * created with this capability (e.g., via \e gradp = true in
187  * SphericalHarmonic::Circle). If not, \e gradx, etc., will not be
188  * touched.
189  **********************************************************************/
191  real& gradx, real& grady, real& gradz) const throw() {
192  real coslon, sinlon;
193  cossin(lon, coslon, sinlon);
194  return (*this)(coslon, sinlon, gradx, grady, gradz);
195  }
196  };
197 
198 } // namespace GeographicLib
199 
200 #if defined(_MSC_VER)
201 # pragma warning (pop)
202 #endif
203 
204 #endif // GEOGRAPHICLIB_CIRCULARENGINE_HPP
#define GEOGRAPHICLIB_EXPORT
Definition: Constants.hpp:52
The evaluation engine for SphericalHarmonic.
Math::real operator()(real lon) const
GeographicLib::Math::real real
Definition: GeodSolve.cpp:40
Geomagnetic field on a circle of latitude.
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
static T sq(T x)
Definition: Math.hpp:153
Spherical harmonic sums for a circle.
Math::real operator()(real coslon, real sinlon) const
Header for GeographicLib::Constants class.
Header for GeographicLib::SphericalEngine class.
Gravity on a circle of latitude.