geomodels package¶

Geographic data models.

Geoid model¶

class geomodels.GeoidModel¶

Bases: object

Looking up the height of the geoid above the ellipsoid.

This class evaluates the height of one of the standard geoids, EGM84, EGM96, or EGM2008 by bi-linear or cubic interpolation into a rectangular grid of data.

See https://geographiclib.sourceforge.io/html/geoid.html.

The geoids are defined in terms of spherical harmonics. However in order to provide a quick and flexible method of evaluating the geoid heights, this class evaluates the height by interpolation into a grid of pre-computed values.

The height of the geoid above the ellipsoid, N, is sometimes called the geoid undulation. It can be used to convert a height above the ellipsoid, h, to the corresponding height above the geoid (the orthometric height, roughly the height above mean sea level), H, using the relations:

h = N + H
H = -N + h

This class is typically not thread safe in that a single instantiation cannot be safely used by multiple threads because of the way the object reads the data set and because it maintains a single-cell cache. If multiple threads need to calculate geoid heights they should all construct thread-local instantiations. Alternatively, set the optional thread-safe parameter to true in the constructor. This causes the constructor to read all the data into memory and to turn off the single-cell caching which results in a Geoid object which is thread safe.

See https://geographiclib.sourceforge.io/html/geoid.html.

Parameters:

name – the name of the geoid
path – (optional) directory for data file
cubic – (optional) interpolation method; false means bilinear, true (the default) means cubic.
threadsafe – (optional), if true, construct a thread safe object. The default is false

Raises:

GeographicErr – if the data file cannot be found, is unreadable, or is corrupt
GeographicErr – if threadsafe is True but the memory necessary for caching the data can’t be allocated.

The data file is formed by appending “.pgm” to the name. If path is specified (and is non-empty), then the file is loaded from directory, path. Otherwise the path is given by default_geoid_path(). If the threadsafe parameter is True, the data set is read into memory, the data file is closed, and single-cell caching is turned off; this results in a Geoid object which is thread safe.

Example¶

# Example of using the `geomodels.GeoidModel` class.
# This requires that the `egm96-5` geoid model be installed; see
# https://geographiclib.sourceforge.io/C++/doc/geoid.html#geoidinst

from geomodels import GeoidModel, EHeightConvDir

GEOIDTOELLIPSOID = EHeightConvDir.GEOIDTOELLIPSOID

egm96 = Geoid("egm96-5")

# Convert height above egm96 to height above the ellipsoid
lat = 42
lon = -75
height_above_geoid = 20

geoid_height = egm96(lat, lon)
height_above_ellipsoid = (
    height_above_geoid + GEOIDTOELLIPSOID * geoid_height
)

print(f"Height above ellipsoid: {height_above_ellipsoid}")

cache() → bool¶: Return true if a data cache is active.

cache_all()¶

Cache all the data.

Raises:

GeographicErr – if the memory necessary for caching the data can’t be allocated (in this case, you will have no cache and can try again with a smaller area).
GeographicErr – if there’s a problem reading the data.
GeographicErr – if this is called on a thread-safe Geoid.

On most computers, this is fast for data sets with grid resolution of 5’ or coarser. For a 1’ grid, the required RAM is 450MB; a 2.5’ grid needs 72MB; and a 5’ grid needs 18MB.

cache_area(south: float, west: float, north: float, east: float)¶

Set up a cache.

Parameters:

south (float) – south latitude (degrees) of the south edge of the cached area.
west – west longitude (degrees) of the west edge of the cached area.
north – north latitude (degrees) of the north edge of the cached area.
east – east longitude (degrees) of the east edge of the cached area.

Raises:

GeographicErr – if the memory necessary for caching the data can’t be allocated (in this case, you will have no cache and can try again with a smaller area).
GeographicErr – if there’s a problem reading the data.
GeographicErr – if this is called on a thread-safe Geoid.

Cache the data for the specified “rectangular” area bounded by the parallels south and north and the meridians west and east. east is always interpreted as being east of west, if necessary by adding 360 deg; to its value. south and north should be in the range [-90 deg, +90 deg].

cache_clear()¶

Clear the cache.

This never throws an error. This does nothing with a thread safe Geoid.

cache_east() → float¶

East edge of the cached area.

The cache excludes this edge.

cache_north() → float¶

North edge of the cached area.

The cache includes this edge.

cache_south() → float¶

South edge of the cached area.

The cache excludes this edge unless it’s the south pole.

cache_west() → float¶

West edge of the cached area.

The cache includes this edge.

convert_height(lat, lon, h, direction: EHeightConvDir)¶

Convert a height above the geoid to a height above the ellipsoid and vice versa.

Parameters:

lat – latitude of the point (degrees).
lon – longitude of the point (degrees).
h – height of the point (meters).
direction – a EHeightConvDir specifying the direction of the conversion; EHeightConvDir.GEOIDTOELLIPSOID means convert a height above the geoid to a height above the ellipsoid; Geoid.ELLIPSOIDTOGEOID means convert a height above the ellipsoid to a height above the geoid.

Raises:

GeographicErr – if there’s a problem reading the data; this never happens if (lat, lon) is within a successfully cached area.

Returns:

converted height (meters).

datetime() → str¶

Return date of the data file.

If absent, return “UNKNOWN”.

static default_geoid_name() → str¶

The default name for the geoid.

This is the value of the environment variable GEOGRAPHICLIB_GEOID_NAME, if set; otherwise, it is “egm96-5”. The GeoidModel class does not use this function; it is just provided as a convenience for a calling program when constructing a GeoidModel object.

static default_geoid_path() → str¶

The default path for geoid data files.

This is the value of the environment variable GEOGRAPHICLIB_GEOID_PATH, if set; otherwise, it is ${GEOGRAPHICLIB_DATA}/geoids if the environment variable GEOGRAPHICLIB_DATA is set; otherwise, it is a compile-time default.

description() → str¶

Return geoid description, if available, in the data file.

If the geoid description is absent in the data file return “NONE”.

equatorial_radius() → float¶

The equatorial radius of the ellipsoid (meters).

(The WGS84 value is returned because the supported geoid models are all based on this ellipsoid.)

flattening() → float¶

The flattening of the ellipsoid.

(The WGS84 value is returned because the supported geoid models are all based on this ellipsoid.

geoid_directory() → str¶: Return the directory used to load the geoid data.

geoid_file() → str¶: Return full file name used to load the geoid data.

geoid_name() → str¶

Return the “name” used to load the geoid data.

“name” is the first argument of the constructor.

interpolation() → str¶: Return interpolation method (cubic or bilinear).

max_error() → float¶

Return an estimate of the maximum interpolation and quantization error [m].

This relies on the value being stored in the data file. If the value is absent, return -1.

offset() → float¶

Offset (meters).

This in used in converting from the pixel values in the data file to geoid heights.

rms_error() → float¶

Return an estimate of the RMS interpolation and quantization error [m].

This relies on the value being stored in the data file. If the value is absent, return -1.

scale() → float¶

Scale (meters).

This in used in converting from the pixel values in the data file to geoid heights.

threadsafe() → bool¶: Return true if the object is constructed to be thread safe.

class geomodels.EHeightConvDir(value, names=<not given>, *values, module=None, qualname=None, type=None, start=1, boundary=None)¶

Bases: IntEnum

ELLIPSOIDTOGEOID = -1¶

GEOIDTOELLIPSOID = 1¶

NONE = 0¶

Gravity model¶

class geomodels.GravityModel¶

Bases: object

Model of the earth’s gravity field.

Evaluate the earth’s gravity field according to a model. The supported models treat only the gravitational field exterior to the mass of the earth. When computing the field at points near (but above) the surface of the earth a small correction can be applied to account for the mass of the atmosphere above the point in question; see The effect of the mass of the atmosphere. Determining the height of the geoid above the ellipsoid entails correcting for the mass of the earth above the geoid. The egm96 and egm2008 include separate correction terms to account for this mass.

Definitions and terminology (from Heiskanen and Moritz, Sec 2-13):

$V$ = gravitational potential
$\Phi$ = rotational potential
$W = V + \Phi = T + U$ = total potential
$V_{0}$ = normal gravitation potential
$U = V_{0} + \Phi$ = total normal potential
$T = W - U = V - V_{0}$ = anomalous or disturbing potential
$g = \nabla W = \gamma + \delta$
$f = \nabla \Phi$
$\Gamma = \nabla V_{0}$
$\gamma = \nabla U$
$\delta = \nabla T$ = gravity disturbance vector = $g_{P} - \gamma_{P}$
$\delta g$ = gravity disturbance = $g_{P} - \gamma_{P}$
$\Delta g$ = gravity anomaly vector = $g_{P} - \gamma_{Q}$ here the line $PQ$ is perpendicular to ellipsoid and the potential at $P$ equals the normal potential at $Q$
$\Delta g$ = gravity anomaly = $g{P} - \gamma_{Q}$
$(\xi, \eta)$ deflection of the vertical, the difference in directions of $g_{P}$ and $\gamma_{Q}, \xi = NS, \eta = EW$
$X$, $Y$, $Z$, geocentric coordinates
$x$, $y$, $z$, local cartesian coordinates used to denote the east, north and up directions.

References:

W. A. Heiskanen and H. Moritz, Physical Geodesy (Freeman, San Francisco, 1967).

See https://geographiclib.sourceforge.io/html/gravity.html.

Parameters:

name (str) – the name of the model
path – (optional) directory for data file
max_degree – (optional) if non-negative, truncate the degree of the model this value
max_order – (optional) if non-negative, truncate the order of the model this value

Raises:

GeographicErr – if the data file cannot be found, is unreadable, or is corrupt
MemoryError – if the memory necessary for storing the model can’t be allocated

A filename is formed by appending “.egm” (World Gravity Model) to the name. If path is specified (and is non-empty), then the file is loaded from directory, path. Otherwise the path is given by default_gravity_path().

This file contains the metadata which specifies the properties of the model. The coefficients for the spherical harmonic sums are obtained from a file obtained by appending “.cof” to metadata file (so the filename ends in “.egm.cof”).

Example¶

# Example of using the `geomodels.GravityModel` class.
# This requires that the `egm96` gravity model be installed; see
# https://geographiclib.sourceforge.io/C++/doc/gravity.html#gravityinst

from geomodels import GravityModel

grav = GravityModel("egm96")

# Mt Everest
lat = 27.99
lon = 86.93
h = 8820

w, gx, gy, gz = grav.gravity(lat, lon, h)

print(f"sum of the gravitational and centrifugal potentials (m**2/s**2): {w}")
print(f"easterly component of the acceleration (m/s**2): {gx}")
print(f"northerly component of the acceleration (m/s**2): {gy}")
print(f"upward component of the acceleration (m/s**2): {gz} (usually negative")

angular_velocity() → float¶

The angular velocity of the model and the reference ellipsoid.

$[\omega] = rad/s$

datetime() → str¶

Date of the gravity model.

If absent, return “UNKNOWN”.

static default_gravity_name() → str¶

The default name for the gravity model.

This is the value of the environment variable GEOGRAPHICLIB_GRAVITY_NAME, if set; otherwise, it is “egm96”. The GravityModel class does not use this function; it is just provided as a convenience for a calling program when constructing a GravityModel object.

static default_gravity_path() → str¶

the default path for gravity model data files.

This is the value of the environment variable GEOGRAPHICLIB_GRAVITY_PATH, if set; otherwise, it is $GEOGRAPHICLIB_DATA/gravity if the environment variable GEOGRAPHICLIB_DATA is set; otherwise, it is a compile-time default.

degree() → int¶: The maximum degree of the components of the model (Nmax).

description() → str¶

The description of the gravity model.

If the gravity model description is absent in the data file return “NONE”.

disturbance(lat, lon, h)¶

Evaluate the gravity disturbance vector at an arbitrary point above (or below) the ellipsoid.

Parameters:

lat – the geographic latitude (degrees)
lon – the geographic longitude (degrees)
h – the height above the ellipsoid (meters)

Returns:

T: the corresponding disturbing potential ([m**2 / s**2])
delta_x: the easterly component of the disturbance vector ([m / s**2])
delta_y: the northerly component of the disturbance vector ([m / s**2])
delta_z: the upward component of the disturbance vector ([m / s**2])

equatorial_radius() → float¶: The equatorial radius of the ellipsoid (meters).

flattening() → float¶: The flattening of the ellipsoid.

geoid_height(lat, lon)¶

Evaluate the geoid height.

Parameters:

lat – the geographic latitude (degrees)
lon – the geographic longitude (degrees)

Returns:

N the height of the geoid above the reference ellipsoid (meters)

Some approximations are made in computing the geoid height so that the results of the NGA codes are reproduced accurately.

gravity(lat, lon, h)¶

Evaluate the gravity at an arbitrary point above (or below) the ellipsoid.

The function includes the effects of the earth’s rotation.

Parameters:

lat – the geographic latitude (degrees)
lon – the geographic longitude (degrees)
h – the height above the ellipsoid (meters)

Returns:

W: the sum of the gravitational and centrifugal potentials ([m**2 / s**2])
gx: the easterly component of the acceleration ([m / s**2])
gy: the northerly component of the acceleration ([m / s**2])
gz: the upward component of the acceleration ([m / s**2]); this is usually negative

gravity_file() → str¶: Full file name used to load the gravity model.

gravity_model_directory() → str¶: Return the directory used to load the gravity model.

gravity_model_name() → str¶

Return the “name” used to load the gravity model.

“name” is the first argument of the constructor, but this may be overridden by the model file).

mass_constant() → float¶

The mass constant of the model ([GM] = m**3 /s**2).

It is the product of G the gravitational constant and M the mass of the earth (usually including the mass of the earth’s atmosphere).

order() → int¶: The maximum order of the components of the model (Mmax).

phi(x, y)¶

Evaluate the centrifugal acceleration in geocentric coordinates.

Parameters:

x – geocentric coordinate of point (meters)
y – geocentric coordinate of point (meters)

Returns:

Phi: the centrifugal potential ([m**2 / s**2])
fx: the x component of the centrifugal acceleration (m / s**2])
fy: the y component of the centrifugal acceleration (m / s**2])

reference_mass_constant() → float¶

The mass constant of the reference ellipsoid.

[GM] = m**3 / s**2.

spherical_anomaly(lat, lon, h)¶

Evaluate the components of the gravity anomaly vector using the spherical approximation.

Parameters:

lat – the geographic latitude (degrees)
lon – the geographic longitude (degrees)
h – the height above the ellipsoid (meters)

Returns:

Dg01: the gravity anomaly ([m / s**2])
xi: the northerly component of the deflection of the vertical (degrees)
eta: the easterly component of the deflection of the vertical (degrees).

The spherical approximation (see Heiskanen and Moritz, Sec 2-14) is used so that the results of the NGA codes are reproduced accurately.

t(x, y, z)¶

Evaluate disturbing potential in geocentric coordinates.

Parameters:

x – geocentric coordinate of point (meters)
y – geocentric coordinate of point (meters)
z – geocentric coordinate of point (meters)

Returns:

T = W - U, the disturbing potential (also called the anomalous potential) ([m**2 / s**2])

t_components(x, y, z)¶

Evaluate the components of the gravity disturbance in geocentric coordinates.

Parameters:

x – geocentric coordinate of point (meters)
y – geocentric coordinate of point (meters)
z – geocentric coordinate of point (meters)

Returns:

T = W - U: the disturbing potential (also called the anomalous potential) ([m**2 / s**2])
delta_x: the x component of the gravity disturbance ([m / s**2])
delta_y: the Y component of the gravity disturbance ([m / s**2])
delta_z: the Z component of the gravity disturbance ([m / s**2])

u(x, y, z)¶

Evaluate the components of the acceleration due to normal gravity and the centrifugal acceleration in geocentric coordinates.

Parameters:

x – geocentric coordinate of point (meters)
y – geocentric coordinate of point (meters)
z – geocentric coordinate of point (meters)

Returns:

U = V0 + Phi:` the sum of the normal gravitational and centrifugal potentials ([m**2 / s**2])
gamma_x: the x component of the normal acceleration ([m / s**2])
gamma_y: the y component of the normal acceleration ([m / s**2])
gamma_z: the z component of the normal acceleration ([m / s**2])

v(x, y, z)¶

Evaluate the components of the gravity disturbance in geocentric coordinates.

Parameters:

x – geocentric coordinate of point (meters)
y – geocentric coordinate of point (meters)
z – geocentric coordinate of point (meters)

Returns:

V = W - Phi: the gravitational potential ([m**2 / s**2])
gx: the x component of the acceleration ([m / s**2])
gy: the y component of the acceleration ([m / s**2])
gz: the z component of the acceleration ([m / s**2])

w(x, y, z)¶

Evaluate the components of the acceleration due to gravity and the centrifugal acceleration in geocentric coordinates.

Parameters:

x – geocentric coordinate of point (meters)
y – geocentric coordinate of point (meters)
z – geocentric coordinate of point (meters)

Returns:

W = V + Phi: the sum of the gravitational and centrifugal potentials ([m**2 / s**2])
gx: the x component of the acceleration ([m / s**2])
gy: the y component of the acceleration ([m / s**2])
gz: the z component of the acceleration ([m / s**2])

Magnetic filed model¶

class geomodels.MagneticFieldModel¶

Bases: object

Model of the earth’s magnetic field.

Evaluate the earth’s magnetic field according to a model. At present only internal magnetic fields are handled. These are due to the earth’s code and crust; these vary slowly (over many years). Excluded are the effects of currents in the ionosphere and magnetosphere which have daily and annual variations.

See https://geographiclib.sourceforge.io/html/magnetic.html for details of how to install the magnetic models and the data format.

Parameters:

name (str) – the name of the model
path – (optional) directory for data file

Raises:

GeographicErr – if the data file cannot be found, is unreadable, or is corrupt
MemoryError – if the memory necessary for storing the model can’t be allocated

A filename is formed by appending “.wmm” (World Magnetic Model) to the name. If path is specified (and is non-empty), then the file is loaded from directory, path. Otherwise the path is given by the MagneticFieldModel.default_magnetic_path().

This file contains the metadata which specifies the properties of the model. The coefficients for the spherical harmonic sums are obtained from a file obtained by appending “.cof” to metadata file (so the filename ends in “.wwm.cof”).

The model is not tied to a particular ellipsoidal model of the earth. The final earth argument to the constructor specifies an ellipsoid to allow geodetic coordinates to the transformed into the spherical coordinates used in the spherical harmonic sum.

Example¶

# Example of using the `geomodels.MagneticModel` class.
# This requires that the wmm2010 magnetic model be installed; see
# https://geographiclib.sourceforge.io/C++/doc/magnetic.html#magneticinst

from geomodels import MagneticModel

mag MagneticModel("wmm2010")

# Mt Everest
lat = 27.99
lon = 86.93
h = 8820
t = 2012

Bx, By, Bz = mag(t, lat,lon, h, Bx, By, Bz)
H, F, D, I = MagneticModel.compute_field_components(Bx, By, Bz)

print(f"horizontal magnetic field (nT): {H}")
print(f"total magnetic field (nT): {F}")
print(f"declination of the field (degrees east of north): {D}")
print(f"inclination of the field (degrees down from horizontal): {I}")

compute_with_rate(t, lat, lon, h)¶

Compute the magnetic field and its rate.

Evaluate the components of the geomagnetic field and their time derivatives

Parameters:

t – the time (years)
lat – latitude of the point (degrees)
lon – longitude of the point (degrees)
h – the height of the point above the ellipsoid (meters)

Returns:

Bx, By, Bz, Bxt, Byt, Bzt: where Bx, By, Bz are the easterly, northerly and vertical (up) components of the magnetic field (nT) respectively, and Bxt, Byt, Bzt are the rate (nT/yr) of change of Bx, By and Bz respectively.

datetime() → str¶

Return date of the model.

Return date of the model, if available, from the ReleaseDate field in the data file; if absent, return “UNKNOWN”.

static default_magnetic_name() → str¶

The default name for the magnetic model.

This is the value of the environment variable GEOGRAPHICLIB_MAGNETIC_NAME, if set; otherwise, it is “wmm2015”. The MagneticFieldModel class does not use this function; it is just provided as a convenience for a calling program when constructing a MagneticFieldModel object.

static default_magnetic_path() → str¶

Return the default path for magnetic model data files.

This is the value of the environment variable GEOGRAPHICLIB_MAGNETIC_PATH, if set; otherwise, it is $GEOGRAPHICLIB_DATA/magnetic if the environment variable GEOGRAPHICLIB_DATA is set; otherwise, it is a compile-time default.

degree() → int¶: The maximum degree of the components of the model (Nmax).

description() → str¶

Return the description of the magnetic model.

Return the description of the magnetic model if available, from the MagneticFieldModel.description() file in the data file; if absent, return “NONE”.

equatorial_radius() → float¶

The equatorial radius of the ellipsoid (meters).

This is the value of ‘a’ inherited from the Geocentric object used in the constructor.

static field_components(Bx, By, Bz)¶

Compute various quantities dependent on the magnetic field.

Parameters:

Bx (float) – the x (easterly) component of the magnetic field (nT)
By (float) – the y (northerly) component of the magnetic field (nT)
Bz (float) – the z (vertical, up positive) component of the magnetic field (nT)

Returns:

H the horizontal magnetic field (nT)
F the total magnetic field (nT)
D the declination of the field (degrees east of north)
I the inclination of the field (degrees down from horizontal)

flattening() → float¶

The flattening of the ellipsoid.

This is the value inherited from the Geocentric object used in the constructor.

magnetic_file() → str¶: Full file name used to load the magnetic model.

magnetic_model_directory() → str¶: Directory used to load the magnetic model.

magnetic_model_name() → str¶

Name used to load the magnetic model.

The ‘name’ used to load the magnetic model (from the first argument of the constructor, but this may be overridden by the model file).

max_height() → float¶

Maximum height.

The maximum height above the ellipsoid (in meters) for which this MagneticModel should be used.

Because the model will typically provide useful results slightly outside the range of allowed heights, no check of t argument is made by the MagneticFieldModel.__call__() operator.

max_time() → float¶

The maximum time (in years) for which this model should be used.

Because the model will typically provide useful results slightly outside the range of allowed times, no check of the t argument is made by the MagneticFieldModel.__call__() operator.

min_height() → float¶

Minimum height.

The minimum height above the ellipsoid (in meters) for which this MagneticModel should be used.

Because the model will typically provide useful results slightly outside the range of allowed heights, no check of t argument is made by the MagneticFieldModel.__call__() operator.

min_time() → float¶

The minimum time (in years) for which this model should be used.

Because the model will typically provide useful results slightly outside the range of allowed times, no check of t argument is made by the MagneticFieldModel.__call__() operator.

order() → int¶: The maximum order of the components of the model (Mmax).

Functions¶

geomodels.lib_version_info() → VersionInfo[int, int, int]¶: Return the (major, minor, patch) version of GeographicLib.

Note

the returned information refers to the version of the GeographicLib library used at build time. The version of the shared library actually used at runtime could be different.

geomodels.lib_version_str() → str¶: Return the version string of GeographicLib.

Note

the returned information refers to the version of the GeographicLib library used at build time. The version of the shared library actually used at runtime could be different.

geomodels.test(verbosity: int = 1, failfast: bool = False)[source]¶

Run the test suite for the geomodels package.

Parameters:

verbosity (int) – verbosity level (higher is more verbose). Default: 1.
failfast (bool) – stop the test run on the first error or failure. Default: False.

geomodels package¶

Geoid model¶

Example¶

Gravity model¶

Example¶

Magnetic filed model¶

Example¶

Functions¶

Sub-modules¶

GeoModels

Navigation

Related Topics