GeographicLib
1.38
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
include
GeographicLib
NormalGravity.hpp
Go to the documentation of this file.
1
/**
2
* \file NormalGravity.hpp
3
* \brief Header for GeographicLib::NormalGravity class
4
*
5
* Copyright (c) Charles Karney (2011-2014) <charles@karney.com> and licensed
6
* under the MIT/X11 License. For more information, see
7
* http://geographiclib.sourceforge.net/
8
**********************************************************************/
9
10
#if !defined(GEOGRAPHICLIB_NORMALGRAVITY_HPP)
11
#define GEOGRAPHICLIB_NORMALGRAVITY_HPP 1
12
13
#include <
GeographicLib/Constants.hpp
>
14
#include <
GeographicLib/Geocentric.hpp
>
15
16
namespace
GeographicLib {
17
18
/**
19
* \brief The normal gravity of the earth
20
*
21
* "Normal" gravity refers to an idealization of the earth which is modeled
22
* as an rotating ellipsoid. The eccentricity of the ellipsoid, the rotation
23
* speed, and the distribution of mass within the ellipsoid are such that the
24
* surface of the ellipsoid is a surface of constant potential (gravitational
25
* plus centrifugal). The acceleration due to gravity is therefore
26
* perpendicular to the surface of the ellipsoid.
27
*
28
* There is a closed solution to this problem which is implemented here.
29
* Series "approximations" are only used to evaluate certain combinations of
30
* elementary functions where use of the closed expression results in a loss
31
* of accuracy for small arguments due to cancellation of the two leading
32
* terms. However these series include sufficient terms to give full machine
33
* precision.
34
*
35
* Definitions:
36
* - <i>V</i><sub>0</sub>, the gravitational contribution to the normal
37
* potential;
38
* - Φ, the rotational contribution to the normal potential;
39
* - \e U = <i>V</i><sub>0</sub> + Φ, the total
40
* potential;
41
* - <b>Γ</b> = ∇<i>V</i><sub>0</sub>, the acceleration due to
42
* mass of the earth;
43
* - <b>f</b> = ∇Φ, the centrifugal acceleration;
44
* - <b>γ</b> = ∇\e U = <b>Γ</b> + <b>f</b>, the normal
45
* acceleration;
46
* - \e X, \e Y, \e Z, geocentric coordinates;
47
* - \e x, \e y, \e z, local cartesian coordinates used to denote the east,
48
* north and up directions.
49
*
50
* References:
51
* - W. A. Heiskanen and H. Moritz, Physical Geodesy (Freeman, San
52
* Francisco, 1967), Secs. 1-19, 2-7, 2-8 (2-9, 2-10), 6-2 (6-3).
53
* - H. Moritz, Geodetic Reference System 1980, J. Geodesy 54(3), 395-405
54
* (1980) http://dx.doi.org/10.1007/BF02521480
55
*
56
* Example of use:
57
* \include example-NormalGravity.cpp
58
**********************************************************************/
59
60
class
GEOGRAPHICLIB_EXPORT
NormalGravity
{
61
private
:
62
static
const
int
maxit_ = 10;
63
typedef
Math::real
real
;
64
friend
class
GravityModel
;
65
real _a, _GM, _omega, _f, _J2, _omega2, _aomega2;
66
real _e2, _ep2, _b, _E, _U0, _gammae, _gammap, _q0, _m, _k, _fstar;
67
Geocentric
_earth;
68
// (atan(y)-(y-y^3/3))/y^5 (y = sqrt(x)) = 1/5-y/7+y^2/9-y^3/11...
69
static
real atan5(real x);
70
// (atan(y)-(y-y^3/3+y^5/5))/y^7 (y = sqrt(x)) = -1/7+x/9-x^2/11+x^3/13...
71
static
real atan7(real x);
72
static
real qf(real ep2);
73
static
real dq(real ep2);
74
static
real qpf(real ep2);
75
real Jn(
int
n)
const
;
76
public
:
77
78
/** \name Setting up the normal gravity
79
**********************************************************************/
80
///@{
81
/**
82
* Constructor for the normal gravity.
83
*
84
* @param[in] a equatorial radius (meters).
85
* @param[in] GM mass constant of the ellipsoid
86
* (meters<sup>3</sup>/seconds<sup>2</sup>); this is the product of \e G
87
* the gravitational constant and \e M the mass of the earth (usually
88
* including the mass of the earth's atmosphere).
89
* @param[in] omega the angular velocity (rad s<sup>−1</sup>).
90
* @param[in] f the flattening of the ellipsoid.
91
* @param[in] J2 the dynamical form factor.
92
* @exception if \e a is not positive or the other constants are
93
* inconsistent (see below).
94
*
95
* If \e omega is non-zero, then exactly one of \e f and \e J2 should be
96
* positive and this will be used to define the ellipsoid. The shape of
97
* the ellipsoid can be given in one of two ways:
98
* - geometrically, the ellipsoid is defined by the flattening \e f = (\e a
99
* − \e b) / \e a, where \e a and \e b are the equatorial radius
100
* and the polar semi-axis.
101
* - physically, the ellipsoid is defined by the dynamical form factor
102
* <i>J</i><sub>2</sub> = (\e C − \e A) / <i>Ma</i><sup>2</sup>,
103
* where \e A and \e C are the equatorial and polar moments of inertia
104
* and \e M is the mass of the earth.
105
* .
106
* If \e omega, \e f, and \e J2 are all zero, then the ellipsoid becomes a
107
* sphere.
108
**********************************************************************/
109
NormalGravity
(real a, real GM, real omega, real f, real J2);
110
111
/**
112
* A default constructor for the normal gravity. This sets up an
113
* uninitialized object and is used by GravityModel which constructs this
114
* object before it has read in the parameters for the reference ellipsoid.
115
**********************************************************************/
116
NormalGravity
() : _a(-1) {}
117
///@}
118
119
/** \name Compute the gravity
120
**********************************************************************/
121
///@{
122
/**
123
* Evaluate the gravity on the surface of the ellipsoid.
124
*
125
* @param[in] lat the geographic latitude (degrees).
126
* @return γ the acceleration due to gravity, positive downwards
127
* (m s<sup>−2</sup>).
128
*
129
* Due to the axial symmetry of the ellipsoid, the result is independent of
130
* the value of the longitude. This acceleration is perpendicular to the
131
* surface of the ellipsoid. It includes the effects of the earth's
132
* rotation.
133
**********************************************************************/
134
Math::real
SurfaceGravity(
real
lat)
const
;
135
136
/**
137
* Evaluate the gravity at an arbitrary point above (or below) the
138
* ellipsoid.
139
*
140
* @param[in] lat the geographic latitude (degrees).
141
* @param[in] h the height above the ellipsoid (meters).
142
* @param[out] gammay the northerly component of the acceleration
143
* (m s<sup>−2</sup>).
144
* @param[out] gammaz the upward component of the acceleration
145
* (m s<sup>−2</sup>); this is usually negative.
146
* @return \e U the corresponding normal potential.
147
*
148
* Due to the axial symmetry of the ellipsoid, the result is independent of
149
* the value of the longitude and the easterly component of the
150
* acceleration vanishes, \e gammax = 0. The function includes the effects
151
* of the earth's rotation. When \e h = 0, this function gives \e gammay =
152
* 0 and the returned value matches that of NormalGravity::SurfaceGravity.
153
**********************************************************************/
154
Math::real
Gravity(
real
lat,
real
h,
real
& gammay,
real
& gammaz)
155
const
;
156
157
/**
158
* Evaluate the components of the acceleration due to gravity and the
159
* centrifugal acceleration in geocentric coordinates.
160
*
161
* @param[in] X geocentric coordinate of point (meters).
162
* @param[in] Y geocentric coordinate of point (meters).
163
* @param[in] Z geocentric coordinate of point (meters).
164
* @param[out] gammaX the \e X component of the acceleration
165
* (m s<sup>−2</sup>).
166
* @param[out] gammaY the \e Y component of the acceleration
167
* (m s<sup>−2</sup>).
168
* @param[out] gammaZ the \e Z component of the acceleration
169
* (m s<sup>−2</sup>).
170
* @return \e U = <i>V</i><sub>0</sub> + Φ the sum of the
171
* gravitational and centrifugal potentials
172
* (m<sup>2</sup> s<sup>−2</sup>).
173
*
174
* The acceleration given by <b>γ</b> = ∇\e U =
175
* ∇<i>V</i><sub>0</sub> + ∇Φ = <b>Γ</b> + <b>f</b>.
176
**********************************************************************/
177
Math::real
U(
real
X,
real
Y,
real
Z,
178
real
& gammaX,
real
& gammaY,
real
& gammaZ)
const
;
179
180
/**
181
* Evaluate the components of the acceleration due to gravity alone in
182
* geocentric coordinates.
183
*
184
* @param[in] X geocentric coordinate of point (meters).
185
* @param[in] Y geocentric coordinate of point (meters).
186
* @param[in] Z geocentric coordinate of point (meters).
187
* @param[out] GammaX the \e X component of the acceleration due to gravity
188
* (m s<sup>−2</sup>).
189
* @param[out] GammaY the \e Y component of the acceleration due to gravity
190
* (m s<sup>−2</sup>).
191
* @param[out] GammaZ the \e Z component of the acceleration due to gravity
192
* (m s<sup>−2</sup>).
193
* @return <i>V</i><sub>0</sub> the gravitational potential
194
* (m<sup>2</sup> s<sup>−2</sup>).
195
*
196
* This function excludes the centrifugal acceleration and is appropriate
197
* to use for space applications. In terrestrial applications, the
198
* function NormalGravity::U (which includes this effect) should usually be
199
* used.
200
**********************************************************************/
201
Math::real
V0(
real
X,
real
Y,
real
Z,
202
real
& GammaX,
real
& GammaY,
real
& GammaZ)
const
;
203
204
/**
205
* Evaluate the centrifugal acceleration in geocentric coordinates.
206
*
207
* @param[in] X geocentric coordinate of point (meters).
208
* @param[in] Y geocentric coordinate of point (meters).
209
* @param[out] fX the \e X component of the centrifugal acceleration
210
* (m s<sup>−2</sup>).
211
* @param[out] fY the \e Y component of the centrifugal acceleration
212
* (m s<sup>−2</sup>).
213
* @return Φ the centrifugal potential (m<sup>2</sup>
214
* s<sup>−2</sup>).
215
*
216
* Φ is independent of \e Z, thus \e fZ = 0. This function
217
* NormalGravity::U sums the results of NormalGravity::V0 and
218
* NormalGravity::Phi.
219
**********************************************************************/
220
Math::real
Phi(
real
X,
real
Y,
real
& fX,
real
& fY)
const
;
221
///@}
222
223
/** \name Inspector functions
224
**********************************************************************/
225
///@{
226
/**
227
* @return true if the object has been initialized.
228
**********************************************************************/
229
bool
Init
()
const
{
return
_a > 0; }
230
231
/**
232
* @return \e a the equatorial radius of the ellipsoid (meters). This is
233
* the value used in the constructor.
234
**********************************************************************/
235
Math::real
MajorRadius()
const
236
{
return
Init() ? _a :
Math::NaN
(); }
237
238
/**
239
* @return \e GM the mass constant of the ellipsoid
240
* (m<sup>3</sup> s<sup>−2</sup>). This is the value used in the
241
* constructor.
242
**********************************************************************/
243
Math::real
MassConstant()
const
244
{
return
Init() ? _GM :
Math::NaN
(); }
245
246
/**
247
* @return <i>J</i><sub><i>n</i></sub> the dynamical form factors of the
248
* ellipsoid.
249
*
250
* If \e n = 2 (the default), this is the value of <i>J</i><sub>2</sub>
251
* used in the constructor. Otherwise it is the zonal coefficient of the
252
* Legendre harmonic sum of the normal gravitational potential. Note that
253
* <i>J</i><sub><i>n</i></sub> = 0 if \e n is odd. In most gravity
254
* applications, fully normalized Legendre functions are used and the
255
* corresponding coefficient is <i>C</i><sub><i>n</i>0</sub> =
256
* −<i>J</i><sub><i>n</i></sub> / sqrt(2 \e n + 1).
257
**********************************************************************/
258
Math::real
DynamicalFormFactor(
int
n = 2)
const
259
{
return
Init() ? ( n == 2 ? _J2 : Jn(n)) :
Math::NaN
(); }
260
261
/**
262
* @return ω the angular velocity of the ellipsoid (rad
263
* s<sup>−1</sup>). This is the value used in the constructor.
264
**********************************************************************/
265
Math::real
AngularVelocity()
const
266
{
return
Init() ? _omega :
Math::NaN
(); }
267
268
/**
269
* @return <i>f</i> the flattening of the ellipsoid (\e a − \e b)/\e
270
* a.
271
**********************************************************************/
272
Math::real
Flattening()
const
273
{
return
Init() ? _f :
Math::NaN
(); }
274
275
/**
276
* @return γ<sub>e</sub> the normal gravity at equator (m
277
* s<sup>−2</sup>).
278
**********************************************************************/
279
Math::real
EquatorialGravity()
const
280
{
return
Init() ? _gammae :
Math::NaN
(); }
281
282
/**
283
* @return γ<sub>p</sub> the normal gravity at poles (m
284
* s<sup>−2</sup>).
285
**********************************************************************/
286
Math::real
PolarGravity()
const
287
{
return
Init() ? _gammap :
Math::NaN
(); }
288
289
/**
290
* @return <i>f*</i> the gravity flattening (γ<sub>p</sub> −
291
* γ<sub>e</sub>) / γ<sub>e</sub>.
292
**********************************************************************/
293
Math::real
GravityFlattening()
const
294
{
return
Init() ? _fstar :
Math::NaN
(); }
295
296
/**
297
* @return <i>U</i><sub>0</sub> the constant normal potential for the
298
* surface of the ellipsoid (m<sup>2</sup> s<sup>−2</sup>).
299
**********************************************************************/
300
Math::real
SurfacePotential()
const
301
{
return
Init() ? _U0 :
Math::NaN
(); }
302
303
/**
304
* @return the Geocentric object used by this instance.
305
**********************************************************************/
306
const
Geocentric
&
Earth
()
const
{
return
_earth; }
307
///@}
308
309
/**
310
* A global instantiation of NormalGravity for the WGS84 ellipsoid.
311
**********************************************************************/
312
static
const
NormalGravity
& WGS84();
313
314
/**
315
* A global instantiation of NormalGravity for the GRS80 ellipsoid.
316
**********************************************************************/
317
static
const
NormalGravity
& GRS80();
318
319
/**
320
* Compute the flattening from the dynamical form factor.
321
*
322
* @param[in] a equatorial radius (meters).
323
* @param[in] GM mass constant of the ellipsoid
324
* (meters<sup>3</sup>/seconds<sup>2</sup>); this is the product of \e G
325
* the gravitational constant and \e M the mass of the earth (usually
326
* including the mass of the earth's atmosphere).
327
* @param[in] omega the angular velocity (rad s<sup>−1</sup>).
328
* @param[in] J2 the dynamical form factor.
329
* @return \e f the flattening of the ellipsoid.
330
**********************************************************************/
331
static
Math::real
J2ToFlattening(
real
a,
real
GM,
real
omega,
real
J2);
332
333
/**
334
* Compute the dynamical form factor from the flattening.
335
*
336
* @param[in] a equatorial radius (meters).
337
* @param[in] GM mass constant of the ellipsoid
338
* (meters<sup>3</sup>/seconds<sup>2</sup>); this is the product of \e G
339
* the gravitational constant and \e M the mass of the earth (usually
340
* including the mass of the earth's atmosphere).
341
* @param[in] omega the angular velocity (rad s<sup>−1</sup>).
342
* @param[in] f the flattening of the ellipsoid.
343
* @return \e J2 the dynamical form factor.
344
**********************************************************************/
345
static
Math::real
FlatteningToJ2(
real
a,
real
GM,
real
omega,
real
f);
346
};
347
348
}
// namespace GeographicLib
349
350
#endif // GEOGRAPHICLIB_NORMALGRAVITY_HPP
Generated by
1.8.3.1