World Geodetic System (1984)#

The World Geodetic System 1984 (WGS 84) [198414] describes the best geodetic reference system for the Earth available for the practical applications of mapping, charting, geopositioning, and navigation, using data, techniques and technology available through 2013 by the United States of America’s National Geospatial-Intelligence Agency (NGA.)

The WGS 84 Coordinate System is a Conventional Terrestrial Reference System (CTRS), that follows the criteria outlined in the International Earth Rotation and Reference Systems Service (IERS):

  • It is geocentric, meaning the center of mass is defined for the whole Earth including oceans and atmosphere.

  • Its scale is that of the local Earth frame, in the meaning of a relativistic theory of gravitation.

  • Its orientation was initially given by the Bureau International de l’Heure (BIH) at the epoch 1984.0

  • It is defined as a right-handed, orthogonal, and Earth-fixed coordinate system, where the Z-axis serves as the rotational axis of the ellipsoid revolution.

WGS 84 (G1762) is the sixth and latest update (16 Oct 2013) to the realization of the WGS 84 Reference Frame, and is updated to incorporate international conventions and alignements to the International Terrestrial Reference Frame 2008 (ITRF2008.)

A simple and practical approach is a mathematically manageable reference surface (an ellipsoid) approximating the broad features of the figure and of the gravity field of the Earth, which is used to describe the Earth surface and the forces acting on and above it.

A second, and more precise, solution considers an equipotential surface of the gravity field of the Earth, called the geoid, which is used to define the Earth Gravitational Model (EGM.) The implementation of this second model falls out of the scope of this module.

The WGS 84 ellipsoid is defined by the semi-major axis (\(a\)), the reciprocal flattening (\(1/f\)) of an oblate ellipsoid of revolution, the Geocentric Gravitational Constant (\(GM\)) and the angular velocity (\(\omega\))

Symbol

Definition

Value

Unit

\(a\)

Semi-major Axis

6378137.0

m

\(1/f\)

Flattening Factor of the Earth

298.257223563

\(GM\)

Geocentric Gravitational Constant

3.986004418 x 10^14

m^3/s^2

\(\omega\)

Earth’s Nominal Mean Angular Velocity

7.292115 x 10^-5

rad/s

The first two parameters (\(a\), \(1/f\)) define the geometry of the rotational ellipsoid, while the other two parameters (\(GM\), \(\omega\)) permit the unique determination of its associated normal gravity field.

Having these 4 elements defined, it is possible to estimate most of the WGS84 parameters directly.

The class WGS already sets these values by default (and estimates the semi-minor axis b), but with slightly different notation:

>>> wgs = ahrs.utils.WGS()
>>> wgs.a       # Semi-major Axis
6378137.0
>>> 1/wgs.f     # Flattening Factor of the Earth
298.257223563
>>> wgs.gm      # Geocentric Gravitational Constant
398600441800000.0
>>> wgs.w       # Earth's Nominal Mean Angular Velocity
7.292115e-05
>>> wgs.b       # Semi-minor Axis
6356752.314245179

Furthermore, there are three unitless values used in the computation of the Moments of Inertia:

Symbol

Definition

Value

C 2,0 dyn[2008]

Dynamic Second Degree Zonal Harmonics

-4.84165143790815 x 10^-4

C 2,2 dyn[2008]

Dynamic Second Degree Sectorial Harmonics

2.43938357328313 x 10^-6

H

Dynamic Ellipticity

3.2737949 x 10^-3

The dynamic harmonics are recovered from the data obtained empirically in EGM2008, and are NOT derived from the ellipsoid parameters. The subscript dyn[2008] denote their origin.

The dynamic ellipticity is a factor in the theoretical value of the rate of precession of the equinoxes, which is known from observation too.

The WGS 84 Ellipsoid is identified as a geocentric, equipotential ellipsoid of revolution, i.e., an ellipsoid with a surface on which the value of the gravity potential is the same everywhere.

Earth’s Gravity Field#

In a rectangular system a point \(\mathbf{p}\) is located with the coordinates \(\begin{pmatrix}x & y & z\end{pmatrix}\), but in ellipsoidal coordinates, this point is found with \(\begin{pmatrix}u & \phi & \lambda\end{pmatrix}\), where \(u\) is the semi-minor axis of the ellipsoid, \(\phi\) is the angle between the plumb line and the equatorial plane called the geographical latitude, and \(\lambda\) is the angle between \(\mathbf{p}\) and the meridian plane of Greenwich called the geographical longitude.

We assume Earth is an ellipsoid of revolution which is an equipotential surface of a normal gravity field. Although the Earth is not a perfect ellipsoid, its gravity field is easier to handle assuming it is one. The deviations of the field are so small that they are considered linear for this model.

The definition of the potential of the normal gravity field \(U\) [HWM06] is:

\[U = V + \Phi\]

where \(V\) is the potential of gravitational force defined as [1]:

\[V = \frac{GM}{E}\arctan\Big(\frac{E}{u}\Big) + \frac{1}{2}\omega^2a^2\frac{q}{q_0}\Big(\sin^2\beta-\frac{1}{3}\Big)\]

and \(\Phi\) is the potential of centrifugal force:

\[\Phi = \frac{1}{2} \omega^2(u^2+E^2)\cos^2\beta\]

Here, the ellipsoid’s linear eccentricity is simply \(E = \sqrt{a^2-b^2}\), and we set the helper constants:

\[\begin{split}\begin{array}{ll} q &= \frac{1}{2} \Big[\Big(1+\frac{3u^2}{E^2}\Big)\arctan\frac{E}{u}-\frac{3u}{E}\Big] \\ q_0 &= \frac{1}{2} \Big[\Big(1+\frac{3b^2}{E^2}\Big)\arctan\frac{E}{b}-\frac{3b}{E}\Big] \\ q_0' &= 3\Big[\Big(1+\frac{1}{e'^2}\Big)\Big(1-\frac{1}{e'}\arctan e'\Big)\Big] - 1 \end{array}\end{split}\]

The normal gravity vector \(\vec{\mathbf{g}}\) is obtained from its potential \(U\) :

\[\begin{split}\vec{\mathbf{g}} = \nabla U = \begin{bmatrix}\mathrm{g}_u \\ \mathrm{g}_\beta \\ \mathrm{g}_\lambda \end{bmatrix} = \begin{bmatrix} \frac{1}{w}\frac{\partial U}{\partial u} \\ \frac{1}{w\sqrt{u^2+E^2}}\frac{\partial U}{\partial \beta} \\ 0 \end{bmatrix}\end{split}\]

We see that \(\mathrm{g}_\lambda=0\), because the effects of the potential along the longitudinal direction are neglected.

Normal Gravity on the Surface#

The magnitude of the normal gravity vector, simply called normal gravity, is:

\[\mathrm{g} = \|\vec{\mathbf{g}}\| = \sqrt{\mathrm{g}_u^2+\mathrm{g}_\beta^2}\]

where \(u=b\) at the surface of the ellipsoid making \(\mathrm{g}_{\beta,0}=0\). So, the total gravity on its surface, \(\mathrm{g}_0\), is just:

\[\mathrm{g}_0 = \|\vec{\mathbf{g}_0}\| = \frac{GM}{a\sqrt{a^2\sin^2\beta+b^2\cos^2\beta}} \Big[\Big(1+\frac{me'q_0'}{3q_0}\Big)\sin^2\beta + \Big(1-m-\frac{me'q_0'}{6q_0}\Big)\cos^2\beta\Big]\]

where:

\[m = \frac{\omega^2a^2b}{GM}\]

and \(e'=\frac{E}{b}\) is the second eccentricity.

From here we can find two basic variables that will help us to linearly estimate the normal gravity at any point on the ellipsoid.

Keeping it on the surface, we can estimate the normal gravity at the equator:

\[\mathrm{g}_e = \mathrm{g}_0(\beta=0°) = \frac{GM}{ab}\Big(1-m-\frac{me'q_0'}{6q_0}\Big)\]

Similarly, we estimate the normal gravity at the poles:

\[\mathrm{g}_p = \mathrm{g}_0(\beta=90°) = \frac{GM}{a^2}\Big(1+\frac{me'q_0'}{3q_0}\Big)\]

With these two basic values, we can linearly approximate the normal at any latitude on the surface of the ellipsoid, but we need to do it in geographical coordinates.

Using the property \(\tan\beta=\frac{b}{a}\tan\phi\), we define a closed form formula to find the normal gravity \(\mathrm{g}\) at any given latitude \(\phi\) [Som29]:

\[\mathrm{g}(\phi) = \frac{ag_e \cos^2\phi + bg_p\sin^2\phi}{\sqrt{a^2\cos^2\phi + b^2\sin^2\phi}}\]

For numerical computation, a more convenient form is:

\[\mathrm{g}(\phi) = \mathrm{g}_e\frac{1+k\sin^2\phi}{\sqrt{1-e^2\sin^2\phi}}\]

using the helper variable \(k = \frac{b\mathrm{g}_p}{a\mathrm{g}_e} - 1\).

The estimation of the normal gravity is already simplified with the class WGS requiring the latitude only.

>>> wgs.normal_gravity(0.0)     # Normal gravity at Equator (latitude = 0.0 °)
9.78032533590406
>>> wgs.normal_gravity(50.0)    # Normal gravity at latitude = 50.0 °
9.810702135603085

Normal Gravity above the Surface#

At small heights above the surface, the normal gravity can be approximated with a truncated Taylor Series with a positive direction downward along the geodetic normal to the ellipsoid:

\[\mathrm{g}(\phi, h) = \mathrm{g}(\phi) \Big(1 - \frac{2}{a}\big(1+f+m-2f\sin^2\phi\big)h + \frac{3}{a^2}h^2\Big)\]

where \(h\) is the height, in meters, above the ellipsoid’s surface.

>>> wgs.normal_gravity(50.0, 1000.0)    # Gravity at latitude = 50.0 °, 1000 m above surface
9.807617683884756

Advanced Gravitational Methods#

All methods above can be used for cartography and basic gravitational references. For more advanced and precise estimations it is recommended to use the EGM2008, which is defined, like the WMM, with spherical harmonics, but of degree 2190 and order 2159.

Because of its high complexity and demanding computation, the implementation of the EGM2008 is left out of this module in favour of the more convenient applications developed in Fortran by the NGA.

Footnotes#

class ahrs.utils.wgs84.WGS(a: float = 6378137.0, f: float = 0.0033528106647474805, GM: float = 398600441800000.0, w: float = 7.292115e-05)#

Bases: ReferenceEllipsoid

World Geodetic System 1984

This class is a subclass of the ReferenceEllipsoid and sets the default values for the WGS 84 model:

Symbol

Definition

Value

Unit

\(a\)

Semi-major Axis

6378137.0

m

\(1/f\)

Flattening Factor of the Earth

298.257223563

\(GM\)

Geocentric Gravitational Constant

3.986004418 x 10^14

m^3/s^2

\(\omega\)

Earth’s Nominal Mean Angular Velocity

7.292115 x 10^-5

rad/s
Parameters:

  • a (float, default: 6378137.0) – Ellipsoid’s Semi-major axis (Equatorial Radius), in meters. Defaults to Earth’s semi-major axis.

  • f (float, default: 0.0033528106647474805) – Ellipsoid’s flattening factor. Defaults to Earth’s flattening.

  • GM (float, default: 3.986004418e14) – Ellipsoid’s Standard Gravitational Constant in m^3/s^2

  • w (float, default: 0.00007292115) – Ellipsoid’s rotation rate in rad/s

Examples

>>> wgs = ahrs.utils.WGS()
>>> wgs.a
6378137.0
>>> wgs.f
0.0033528106647474805
>>> wgs.gm
398600441800000.0
>>> wgs.w
7.292115e-05

It inherits methods and properties from the ReferenceEllipsoid class, and adds some specific methods for the WGS 84 model.

>>> wgs.mass
5.972186390142457e+24
>>> wgs.sidereal_day
86164.09053083288
>>> wgs.normal_gravity(0.0)     # Normal gravity at Equator (latitude = 0.0 °)
9.78032533590406
>>> wgs.normal_gravity(50.0)    # Normal gravity at latitude = 50.0 °
9.810702135603085
property atmosphere_gravitational_constant: float#

Gravitational Constant of the Atmosphere \(GM_A\), computed as:

\[GM_A = G\ M_A\]

where \(M_A\) is the total mean mass of the atmosphere (with water vapor) equal to \(5.148\times 10^{18} \mathrm{kg}\).

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.atmosphere_gravitational_constant
343591934.4
property dynamic_inertial_moment_about_X: float#

Dynamic Moment of Inertia (\(A\)), with respect to the X-Axis of Rotation, computed as:

\[A_{dyn} = -\sqrt{5}Ma^2\Big[\Big(1-\frac{1}{H}\Big)C_{2,0dyn} - \frac{C_{2,2dyn}}{\sqrt{3}}\Big]\]

where \(H=3.2737949 \times 10^{-3}\) is the Dynamic Ellipticity.

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.dynamic_inertial_moment_about_X
8.007921777277886e+37
property dynamic_inertial_moment_about_Y: float#

Dynamic Moment of Inertia (\(B\)), with respect to the X-Axis of Rotation, computed as:

\[B_{dyn} = -\sqrt{5}Ma^2\Big[\Big(1-\frac{1}{H}\Big)C_{2,0dyn} + \frac{C_{2,2dyn}}{\sqrt{3}}\Big]\]

where \(H=3.2737949 \times 10^{-3}\) is the Dynamic Ellipticity.

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.dynamic_inertial_moment_about_Y
8.008074799852911e+37
property dynamic_inertial_moment_about_Z: float#

Dynamic Moment of Inertia (\(C\)), with respect to the Z-Axis of Rotation, computed as:

\[C_{dyn} = -\sqrt{5}Ma^2\frac{C_{2,0dyn}}{H}\]

where \(H=3.2737949 \times 10^{-3}\) is the Dynamic Ellipticity.

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.dynamic_inertial_moment_about_Z
8.03430094201443e+37
property geometric_dynamic_ellipticity: float#

Geometric Solution for Dynamic Ellipticity \(H\), computed as:

\[H_{geo} = \frac{C_{geo}-A_{geo}}{C_{geo}}\]

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.geometric_dynamic_ellipticity
0.003258100628533992
property geometric_inertial_moment: float#

Geometric Moment of Inertia (\(A\)), with respect to Any Axis in the Equatorial Plane, computed as:

\[A_{geo} = C_{geo} + \sqrt{5}Ma^2 C_{2,0geo}\]

where \(C_{2,0geo} = -4.84166774985\times 10^{-4}\) is Earth’s Geographic Second Degree Zonal Harmonic.

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.geometric_inertial_moment
8.046726628049449e+37
property geometric_inertial_moment_about_Z: float#

Geometric Moment of Inertia (\(C\)), with respect to the Z-Axis of Rotation, computed as:

\[C_{geo} = \frac{2}{3}Ma^2\Big(1-\frac{2}{5}\sqrt{\frac{5m}{2f}-1}\Big)\]

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.geometric_inertial_moment_about_Z
8.073029370114392e+37
property gravitational_constant_without_atmosphere: float#

Geocentric Gravitational Constant with Earth’s Atmosphere Excluded \(GM'\), computed as:

\[GM' = GM - GM_A\]

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.gravitational_constant_without_atmosphere
398600098208065.6
property mass: float#

The Mass \(M\) of the Earth, in kg, computed as:

\[M = \frac{GM}{G}\]

where \(G\) is the universal constant of gravitation equal to \(6.67428\times 10^{-11} \frac{\mathrm{m}^3}{\mathrm{kg} \mathrm{s}^2}\)

Note

The universal constant of gravitation for WGS84 is not the same as the one defined by CODATA [TMNT21a] (\(6.67430\times 10^{-11}\)), which is the default value in the ReferenceEllipsoid class.

Example

>>> wgs = ahrs.utils.WGS()
>>> wgs.mass
5.972186390142457e+24