00001 /** 00002 * \file GravityCircle.hpp 00003 * \brief Header for GeographicLib::GravityCircle class 00004 * 00005 * Copyright (c) Charles Karney (2011) <charles@karney.com> and licensed under 00006 * the MIT/X11 License. For more information, see 00007 * http://geographiclib.sourceforge.net/ 00008 **********************************************************************/ 00009 00010 #if !defined(GEOGRAPHICLIB_GRAVITYCIRCLE_HPP) 00011 #define GEOGRAPHICLIB_GRAVITYCIRCLE_HPP 1 00012 00013 #include <vector> 00014 #include <GeographicLib/Constants.hpp> 00015 #include <GeographicLib/CircularEngine.hpp> 00016 #include <GeographicLib/GravityModel.hpp> 00017 00018 namespace GeographicLib { 00019 00020 /** 00021 * \brief Gravity on a circle of latitude 00022 * 00023 * Evaluate the earth's gravity field on a circle of constant height and 00024 * latitude. This uses a CircularEngine to pre-evaluate the inner sum of the 00025 * spherical harmonic sum, allowing the values of the field at several 00026 * different longitudes to be evaluated rapidly. 00027 * 00028 * Use GravityModel::Circle to create a GravityCircle object. (The 00029 * constructor for this class is private.) 00030 * 00031 * See \ref gravityparallel for an example of using GravityCircle (together 00032 * with OpenMP) to speed up the computation of geoid heights. 00033 * 00034 * Example of use: 00035 * \include example-GravityCircle.cpp 00036 * 00037 * <a href="Gravity.1.html">Gravity</a> is a command-line utility providing 00038 * access to the functionality of GravityModel and GravityCircle. 00039 **********************************************************************/ 00040 00041 class GEOGRAPHICLIB_EXPORT GravityCircle { 00042 private: 00043 typedef Math::real real; 00044 enum mask { 00045 NONE = GravityModel::NONE, 00046 GRAVITY = GravityModel::GRAVITY, 00047 DISTURBANCE = GravityModel::DISTURBANCE, 00048 DISTURBING_POTENTIAL = GravityModel::DISTURBING_POTENTIAL, 00049 GEOID_HEIGHT = GravityModel::GEOID_HEIGHT, 00050 SPHERICAL_ANOMALY = GravityModel::SPHERICAL_ANOMALY, 00051 ALL = GravityModel::ALL, 00052 }; 00053 00054 unsigned _caps; 00055 real _a, _f, _lat, _h, _Z, _Px, _invR, _cpsi, _spsi, 00056 _cphi, _sphi, _amodel, _GMmodel, _dzonal0, 00057 _corrmult, _gamma0, _gamma, _frot; 00058 CircularEngine _gravitational, _disturbing, _correction; 00059 00060 GravityCircle(mask caps, real a, real f, real lat, real h, 00061 real Z, real P, real cphi, real sphi, 00062 real amodel, real GMmodel, real dzonal0, real corrmult, 00063 real gamma0, real gamma, real frot, 00064 const CircularEngine& gravitational, 00065 const CircularEngine& disturbing, 00066 const CircularEngine& correction) 00067 : _caps(caps) 00068 , _a(a) 00069 , _f(f) 00070 , _lat(lat) 00071 , _h(h) 00072 , _Z(Z) 00073 , _Px(P) 00074 , _invR(1 / Math::hypot(_Px, _Z)) 00075 , _cpsi(_Px * _invR) 00076 , _spsi(_Z * _invR) 00077 , _cphi(cphi) 00078 , _sphi(sphi) 00079 , _amodel(amodel) 00080 , _GMmodel(GMmodel) 00081 , _dzonal0(dzonal0) 00082 , _corrmult(corrmult) 00083 , _gamma0(gamma0) 00084 , _gamma(gamma) 00085 , _frot(frot) 00086 , _gravitational(gravitational) 00087 , _disturbing(disturbing) 00088 , _correction(correction) 00089 {} 00090 00091 friend class GravityModel; // GravityModel calls the private constructor 00092 Math::real W(real clam, real slam, 00093 real& gX, real& gY, real& gZ) const; 00094 Math::real V(real clam, real slam, 00095 real& gX, real& gY, real& gZ) const; 00096 Math::real InternalT(real clam, real slam, 00097 real& deltaX, real& deltaY, real& deltaZ, 00098 bool gradp, bool correct) const; 00099 public: 00100 /** 00101 * A default constructor for the normal gravity. This sets up an 00102 * uninitialized object which can be later replaced by the 00103 * GravityModel::Circle. 00104 **********************************************************************/ 00105 GravityCircle() : _a(-1) {} 00106 00107 /** \name Compute the gravitational field 00108 **********************************************************************/ 00109 ///@{ 00110 /** 00111 * Evaluate the gravity. 00112 * 00113 * @param[in] lon the geographic longitude (degrees). 00114 * @param[out] gx the easterly component of the acceleration 00115 * (m s<sup>−2</sup>). 00116 * @param[out] gy the northerly component of the acceleration 00117 * (m s<sup>−2</sup>). 00118 * @param[out] gz the upward component of the acceleration 00119 * (m s<sup>−2</sup>); this is usually negative. 00120 * @return \e W the sum of the gravitational and centrifugal potentials. 00121 * 00122 * The function includes the effects of the earth's rotation. 00123 **********************************************************************/ 00124 Math::real Gravity(real lon, real& gx, real& gy, real& gz) const; 00125 00126 /** 00127 * Evaluate the gravity disturbance vector. 00128 * 00129 * @param[in] lon the geographic longitude (degrees). 00130 * @param[out] deltax the easterly component of the disturbance vector 00131 * (m s<sup>−2</sup>). 00132 * @param[out] deltay the northerly component of the disturbance vector 00133 * (m s<sup>−2</sup>). 00134 * @param[out] deltaz the upward component of the disturbance vector 00135 * (m s<sup>−2</sup>). 00136 * @return \e T the corresponding disturbing potential. 00137 **********************************************************************/ 00138 Math::real Disturbance(real lon, real& deltax, real& deltay, real& deltaz) 00139 const; 00140 00141 /** 00142 * Evaluate the geoid height. 00143 * 00144 * @param[in] lon the geographic longitude (degrees). 00145 * @return \e N the height of the geoid above the reference ellipsoid 00146 * (meters). 00147 * 00148 * Some approximations are made in computing the geoid height so that the 00149 * results of the NGA codes are reproduced accurately. Details are given 00150 * in \ref gravitygeoid. 00151 **********************************************************************/ 00152 Math::real GeoidHeight(real lon) const; 00153 00154 /** 00155 * Evaluate the components of the gravity anomaly vector using the 00156 * spherical approximation. 00157 * 00158 * @param[in] lon the geographic longitude (degrees). 00159 * @param[out] Dg01 the gravity anomaly (m s<sup>−2</sup>). 00160 * @param[out] xi the northerly component of the deflection of the vertical 00161 * (degrees). 00162 * @param[out] eta the easterly component of the deflection of the vertical 00163 * (degrees). 00164 * 00165 * The spherical approximation (see Heiskanen and Moritz, Sec 2-14) is used 00166 * so that the results of the NGA codes are reproduced accurately. 00167 * approximations used here. Details are given in \ref gravitygeoid. 00168 **********************************************************************/ 00169 void SphericalAnomaly(real lon, real& Dg01, real& xi, real& eta) 00170 const; 00171 00172 /** 00173 * Evaluate the components of the acceleration due to gravity and the 00174 * centrifugal acceleration in geocentric coordinates. 00175 * 00176 * @param[in] lon the geographic longitude (degrees). 00177 * @param[out] gX the \e X component of the acceleration 00178 * (m s<sup>−2</sup>). 00179 * @param[out] gY the \e Y component of the acceleration 00180 * (m s<sup>−2</sup>). 00181 * @param[out] gZ the \e Z component of the acceleration 00182 * (m s<sup>−2</sup>). 00183 * @return \e W = \e V + Φ the sum of the gravitational and 00184 * centrifugal potentials (m<sup>2</sup> s<sup>−2</sup>). 00185 **********************************************************************/ 00186 Math::real W(real lon, real& gX, real& gY, real& gZ) const { 00187 real clam, slam; 00188 CircularEngine::cossin(lon, clam, slam); 00189 return W(clam, slam, gX, gY, gZ); 00190 } 00191 00192 /** 00193 * Evaluate the components of the acceleration due to gravity in geocentric 00194 * coordinates. 00195 * 00196 * @param[in] lon the geographic longitude (degrees). 00197 * @param[out] GX the \e X component of the acceleration 00198 * (m s<sup>−2</sup>). 00199 * @param[out] GY the \e Y component of the acceleration 00200 * (m s<sup>−2</sup>). 00201 * @param[out] GZ the \e Z component of the acceleration 00202 * (m s<sup>−2</sup>). 00203 * @return \e V = \e W - Φ the gravitational potential 00204 * (m<sup>2</sup> s<sup>−2</sup>). 00205 **********************************************************************/ 00206 Math::real V(real lon, real& GX, real& GY, real& GZ) const { 00207 real clam, slam; 00208 CircularEngine::cossin(lon, clam, slam); 00209 return V(clam, slam, GX, GY, GZ); 00210 } 00211 00212 /** 00213 * Evaluate the components of the gravity disturbance in geocentric 00214 * coordinates. 00215 * 00216 * @param[in] lon the geographic longitude (degrees). 00217 * @param[out] deltaX the \e X component of the gravity disturbance 00218 * (m s<sup>−2</sup>). 00219 * @param[out] deltaY the \e Y component of the gravity disturbance 00220 * (m s<sup>−2</sup>). 00221 * @param[out] deltaZ the \e Z component of the gravity disturbance 00222 * (m s<sup>−2</sup>). 00223 * @return \e T = \e W - \e U the disturbing potential (also called the 00224 * anomalous potential) (m<sup>2</sup> s<sup>−2</sup>). 00225 **********************************************************************/ 00226 Math::real T(real lon, real& deltaX, real& deltaY, real& deltaZ) 00227 const { 00228 real clam, slam; 00229 CircularEngine::cossin(lon, clam, slam); 00230 return InternalT(clam, slam, deltaX, deltaY, deltaZ, true, true); 00231 } 00232 00233 /** 00234 * Evaluate disturbing potential in geocentric coordinates. 00235 * 00236 * @param[in] lon the geographic longitude (degrees). 00237 * @return \e T = \e W - \e U the disturbing potential (also called the 00238 * anomalous potential) (m<sup>2</sup> s<sup>−2</sup>). 00239 **********************************************************************/ 00240 Math::real T(real lon) const { 00241 real clam, slam, dummy; 00242 CircularEngine::cossin(lon, clam, slam); 00243 return InternalT(clam, slam, dummy, dummy, dummy, false, true); 00244 } 00245 00246 ///@} 00247 00248 /** \name Inspector functions 00249 **********************************************************************/ 00250 ///@{ 00251 /** 00252 * @return true if the object has been initialized. 00253 **********************************************************************/ 00254 bool Init() const { return _a > 0; } 00255 00256 /** 00257 * @return \e a the equatorial radius of the ellipsoid (meters). This is 00258 * the value inherited from the GravityModel object used in the 00259 * constructor. 00260 **********************************************************************/ 00261 Math::real MajorRadius() const 00262 { return Init() ? _a : Math::NaN(); } 00263 00264 /** 00265 * @return \e f the flattening of the ellipsoid. This is the value 00266 * inherited from the GravityModel object used in the constructor. 00267 **********************************************************************/ 00268 Math::real Flattening() const 00269 { return Init() ? _f : Math::NaN(); } 00270 00271 /** 00272 * @return the latitude of the circle (degrees). 00273 **********************************************************************/ 00274 Math::real Latitude() const 00275 { return Init() ? _lat : Math::NaN(); } 00276 00277 /** 00278 * @return the height of the circle (meters). 00279 **********************************************************************/ 00280 Math::real Height() const 00281 { return Init() ? _h : Math::NaN(); } 00282 00283 /** 00284 * @return \e caps the computational capabilities that this object was 00285 * constructed with. 00286 **********************************************************************/ 00287 unsigned Capabilities() const { return _caps; } 00288 00289 /** 00290 * @param[in] testcaps a set of bitor'ed GeodesicLine::mask values. 00291 * @return true if the GeodesicLine object has all these capabilities. 00292 **********************************************************************/ 00293 bool Capabilities(unsigned testcaps) const { 00294 return (_caps & testcaps) == testcaps; 00295 } 00296 ///@} 00297 }; 00298 00299 } // namespace GeographicLib 00300 00301 #endif // GEOGRAPHICLIB_GRAVITYCIRCLE_HPP