CircularEngine.hpp 6.8 KB

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