GeographicLib
1.21
|
00001 /** 00002 * \file SphericalHarmonic2.hpp 00003 * \brief Header for GeographicLib::SphericalHarmonic2 class 00004 * 00005 * Copyright (c) Charles Karney (2011, 2012) <charles@karney.com> and licensed 00006 * under the MIT/X11 License. For more information, see 00007 * http://geographiclib.sourceforge.net/ 00008 **********************************************************************/ 00009 00010 #if !defined(GEOGRAPHICLIB_SPHERICALHARMONIC2_HPP) 00011 #define GEOGRAPHICLIB_SPHERICALHARMONIC2_HPP \ 00012 "$Id: ce4cda614c1966dea65610bc73bc4db562677fa8 $" 00013 00014 #include <vector> 00015 #include <GeographicLib/Constants.hpp> 00016 #include <GeographicLib/SphericalEngine.hpp> 00017 #include <GeographicLib/CircularEngine.hpp> 00018 00019 namespace GeographicLib { 00020 00021 /** 00022 * \brief Spherical Harmonic series with two corrections to the coefficients. 00023 * 00024 * This classes is similar to SphericalHarmonic, except that the coefficients 00025 * \e C<sub>\e nm</sub> are replaced by \e C<sub>\e nm</sub> + \e tau' 00026 * C'<sub>\e nm</sub> + \e tau'' C''<sub>\e nm</sub> (and similarly for \e 00027 * S<sub>\e nm</sub>). 00028 * 00029 * Example of use: 00030 * \include example-SphericalHarmonic2.cpp 00031 **********************************************************************/ 00032 00033 class GEOGRAPHIC_EXPORT SphericalHarmonic2 { 00034 public: 00035 /** 00036 * Supported normalizations for associate Legendre polynomials. 00037 **********************************************************************/ 00038 enum normalization { 00039 /** 00040 * Fully normalized associated Legendre polynomials. See 00041 * SphericalHarmonic::FULL for documentation. 00042 * 00043 * @hideinitializer 00044 **********************************************************************/ 00045 FULL = SphericalEngine::FULL, 00046 /** 00047 * Schmidt semi-normalized associated Legendre polynomials. See 00048 * SphericalHarmonic::SCHMIDT for documentation. 00049 * 00050 * @hideinitializer 00051 **********************************************************************/ 00052 SCHMIDT = SphericalEngine::SCHMIDT, 00053 /// \cond SKIP 00054 // These are deprecated... 00055 full = FULL, 00056 schmidt = SCHMIDT, 00057 /// \endcond 00058 }; 00059 00060 private: 00061 typedef Math::real real; 00062 SphericalEngine::coeff _c[3]; 00063 real _a; 00064 unsigned _norm; 00065 00066 public: 00067 /** 00068 * Constructor with a full set of coefficients specified. 00069 * 00070 * @param[in] C the coefficients \e C<sub>\e nm</sub>. 00071 * @param[in] S the coefficients \e S<sub>\e nm</sub>. 00072 * @param[in] N the maximum degree and order of the sum 00073 * @param[in] C1 the coefficients \e C'<sub>\e nm</sub>. 00074 * @param[in] S1 the coefficients \e S'<sub>\e nm</sub>. 00075 * @param[in] N1 the maximum degree and order of the first correction 00076 * coefficients \e C'<sub>\e nm</sub> and \e S'<sub>\e nm</sub>. 00077 * @param[in] C2 the coefficients \e C''<sub>\e nm</sub>. 00078 * @param[in] S2 the coefficients \e S''<sub>\e nm</sub>. 00079 * @param[in] N2 the maximum degree and order of the second correction 00080 * coefficients \e C'<sub>\e nm</sub> and \e S'<sub>\e nm</sub>. 00081 * @param[in] a the reference radius appearing in the definition of the 00082 * sum. 00083 * @param[in] norm the normalization for the associated Legendre 00084 * polynomials, either SphericalHarmonic2::FULL (the default) or 00085 * SphericalHarmonic2::SCHMIDT. 00086 * 00087 * See SphericalHarmonic for the way the coefficients should be stored. \e 00088 * N1 and \e N2 should satisfy \e N1 <= \e N and \e N2 <= \e N. 00089 * 00090 * The class stores <i>pointers</i> to the first elements of \e C, \e S, \e 00091 * C', \e S', \e C'', and \e S''. These arrays should not be altered or 00092 * destroyed during the lifetime of a SphericalHarmonic object. 00093 **********************************************************************/ 00094 SphericalHarmonic2(const std::vector<real>& C, 00095 const std::vector<real>& S, 00096 int N, 00097 const std::vector<real>& C1, 00098 const std::vector<real>& S1, 00099 int N1, 00100 const std::vector<real>& C2, 00101 const std::vector<real>& S2, 00102 int N2, 00103 real a, unsigned norm = FULL) 00104 : _a(a) 00105 , _norm(norm) { 00106 if (!(N1 <= N && N2 <= N)) 00107 throw GeographicErr("N1 and N2 cannot be larger that N"); 00108 _c[0] = SphericalEngine::coeff(C, S, N); 00109 _c[1] = SphericalEngine::coeff(C1, S1, N1); 00110 _c[2] = SphericalEngine::coeff(C2, S2, N2); 00111 } 00112 00113 /** 00114 * Constructor with a subset of coefficients specified. 00115 * 00116 * @param[in] C the coefficients \e C<sub>\e nm</sub>. 00117 * @param[in] S the coefficients \e S<sub>\e nm</sub>. 00118 * @param[in] N the degree used to determine the layout of \e C and \e S. 00119 * @param[in] nmx the maximum degree used in the sum. The sum over \e n is 00120 * from 0 thru \e nmx. 00121 * @param[in] mmx the maximum order used in the sum. The sum over \e m is 00122 * from 0 thru min(\e n, \e mmx). 00123 * @param[in] C1 the coefficients \e C'<sub>\e nm</sub>. 00124 * @param[in] S1 the coefficients \e S'<sub>\e nm</sub>. 00125 * @param[in] N1 the degree used to determine the layout of \e C' and \e 00126 * S'. 00127 * @param[in] nmx1 the maximum degree used for \e C' and \e S'. 00128 * @param[in] mmx1 the maximum order used for \e C' and \e S'. 00129 * @param[in] C2 the coefficients \e C''<sub>\e nm</sub>. 00130 * @param[in] S2 the coefficients \e S''<sub>\e nm</sub>. 00131 * @param[in] N2 the degree used to determine the layout of \e C'' and \e 00132 * S''. 00133 * @param[in] nmx2 the maximum degree used for \e C'' and \e S''. 00134 * @param[in] mmx2 the maximum order used for \e C'' and \e S''. 00135 * @param[in] a the reference radius appearing in the definition of the 00136 * sum. 00137 * @param[in] norm the normalization for the associated Legendre 00138 * polynomials, either SphericalHarmonic2::FULL (the default) or 00139 * SphericalHarmonic2::SCHMIDT. 00140 * 00141 * The class stores <i>pointers</i> to the first elements of \e C, \e S, \e 00142 * C', \e S', \e C'', and \e S''. These arrays should not be altered or 00143 * destroyed during the lifetime of a SphericalHarmonic object. 00144 **********************************************************************/ 00145 SphericalHarmonic2(const std::vector<real>& C, 00146 const std::vector<real>& S, 00147 int N, int nmx, int mmx, 00148 const std::vector<real>& C1, 00149 const std::vector<real>& S1, 00150 int N1, int nmx1, int mmx1, 00151 const std::vector<real>& C2, 00152 const std::vector<real>& S2, 00153 int N2, int nmx2, int mmx2, 00154 real a, unsigned norm = FULL) 00155 : _a(a) 00156 , _norm(norm) { 00157 if (!(nmx1 <= nmx && nmx2 <= nmx)) 00158 throw GeographicErr("nmx1 and nmx2 cannot be larger that nmx"); 00159 if (!(mmx1 <= mmx && mmx2 <= mmx)) 00160 throw GeographicErr("mmx1 and mmx2cannot be larger that mmx"); 00161 _c[0] = SphericalEngine::coeff(C, S, N, nmx, mmx); 00162 _c[1] = SphericalEngine::coeff(C1, S1, N1, nmx1, mmx1); 00163 _c[2] = SphericalEngine::coeff(C2, S2, N2, nmx2, mmx2); 00164 } 00165 00166 /** 00167 * A default constructor so that the object can be created when the 00168 * constructor for another object is initialized. This default object can 00169 * then be reset with the default copy assignment operator. 00170 **********************************************************************/ 00171 SphericalHarmonic2() {} 00172 00173 /** 00174 * Compute a spherical harmonic sum with two correction terms. 00175 * 00176 * @param[in] tau1 multiplier for correction coefficients \e C' and \e S'. 00177 * @param[in] tau2 multiplier for correction coefficients \e C'' and \e S''. 00178 * @param[in] x cartesian coordinate. 00179 * @param[in] y cartesian coordinate. 00180 * @param[in] z cartesian coordinate. 00181 * @return \e V the spherical harmonic sum. 00182 * 00183 * This routine requires constant memory and thus never throws an 00184 * exception. 00185 **********************************************************************/ 00186 Math::real operator()(real tau1, real tau2, real x, real y, real z) 00187 const throw() { 00188 real f[] = {1, tau1, tau2}; 00189 real v = 0; 00190 real dummy; 00191 switch (_norm) { 00192 case FULL: 00193 v = SphericalEngine::Value<false, SphericalEngine::FULL, 3> 00194 (_c, f, x, y, z, _a, dummy, dummy, dummy); 00195 break; 00196 case SCHMIDT: 00197 v = SphericalEngine::Value<false, SphericalEngine::SCHMIDT, 3> 00198 (_c, f, x, y, z, _a, dummy, dummy, dummy); 00199 break; 00200 } 00201 return v; 00202 } 00203 00204 /** 00205 * Compute a spherical harmonic sum with two correction terms and its 00206 * gradient. 00207 * 00208 * @param[in] tau1 multiplier for correction coefficients \e C' and \e S'. 00209 * @param[in] tau2 multiplier for correction coefficients \e C'' and \e S''. 00210 * @param[in] x cartesian coordinate. 00211 * @param[in] y cartesian coordinate. 00212 * @param[in] z cartesian coordinate. 00213 * @param[out] gradx \e x component of the gradient 00214 * @param[out] grady \e y component of the gradient 00215 * @param[out] gradz \e z component of the gradient 00216 * @return \e V the spherical harmonic sum. 00217 * 00218 * This is the same as the previous function, except that the components of 00219 * the gradients of the sum in the \e x, \e y, and \e z directions are 00220 * computed. This routine requires constant memory and thus never throws 00221 * an exception. 00222 **********************************************************************/ 00223 Math::real operator()(real tau1, real tau2, real x, real y, real z, 00224 real& gradx, real& grady, real& gradz) const throw() { 00225 real f[] = {1, tau1, tau2}; 00226 real v = 0; 00227 switch (_norm) { 00228 case FULL: 00229 v = SphericalEngine::Value<true, SphericalEngine::FULL, 3> 00230 (_c, f, x, y, z, _a, gradx, grady, gradz); 00231 break; 00232 case SCHMIDT: 00233 v = SphericalEngine::Value<true, SphericalEngine::SCHMIDT, 3> 00234 (_c, f, x, y, z, _a, gradx, grady, gradz); 00235 break; 00236 } 00237 return v; 00238 } 00239 00240 /** 00241 * Create a CircularEngine to allow the efficient evaluation of several 00242 * points on a circle of latitude at fixed values of \e tau1 and \e tau2. 00243 * 00244 * @param[in] tau1 multiplier for correction coefficients \e C' and \e S'. 00245 * @param[in] tau2 multiplier for correction coefficients \e C'' and \e S''. 00246 * @param[in] p the radius of the circle. 00247 * @param[in] z the height of the circle above the equatorial plane. 00248 * @param[in] gradp if true the returned object will be able to compute the 00249 * gradient of the sum. 00250 * @return the CircularEngine object. 00251 * 00252 * SphericalHarmonic2::operator()() exchanges the order of the sums in the 00253 * definition, i.e., sum(n = 0..N)[sum(m = 0..n)[...]] becomes sum(m = 00254 * 0..N)[sum(n = m..N)[...]]. SphericalHarmonic2::Circle performs the 00255 * inner sum over degree \e n (which entails about <i>N</i><sup>2</sup> 00256 * operations). Calling CircularEngine::operator()() on the returned 00257 * object performs the outer sum over the order \e m (about \e N 00258 * operations). This routine may throw a bad_alloc exception in the 00259 * CircularEngine constructor. 00260 * 00261 * See SphericalHarmonic::Circle for an example of its use. 00262 **********************************************************************/ 00263 CircularEngine Circle(real tau1, real tau2, real p, real z, bool gradp) 00264 const { 00265 real f[] = {1, tau1, tau2}; 00266 switch (_norm) { 00267 case FULL: 00268 return gradp ? 00269 SphericalEngine::Circle<true, SphericalEngine::FULL, 3> 00270 (_c, f, p, z, _a) : 00271 SphericalEngine::Circle<false, SphericalEngine::FULL, 3> 00272 (_c, f, p, z, _a); 00273 break; 00274 case SCHMIDT: 00275 default: // To avoid compiler warnings 00276 return gradp ? 00277 SphericalEngine::Circle<true, SphericalEngine::SCHMIDT, 3> 00278 (_c, f, p, z, _a) : 00279 SphericalEngine::Circle<false, SphericalEngine::SCHMIDT, 3> 00280 (_c, f, p, z, _a); 00281 break; 00282 } 00283 } 00284 00285 /** 00286 * @return the zeroth SphericalEngine::coeff object. 00287 **********************************************************************/ 00288 const SphericalEngine::coeff& Coefficients() const throw() 00289 { return _c[0]; } 00290 /** 00291 * @return the first SphericalEngine::coeff object. 00292 **********************************************************************/ 00293 const SphericalEngine::coeff& Coefficients1() const throw() 00294 { return _c[1]; } 00295 /** 00296 * @return the second SphericalEngine::coeff object. 00297 **********************************************************************/ 00298 const SphericalEngine::coeff& Coefficients2() const throw() 00299 { return _c[2]; } 00300 }; 00301 00302 } // namespace GeographicLib 00303 00304 #endif // GEOGRAPHICLIB_SPHERICALHARMONIC2_HPP