typhon.geodesy

Functions for handling geographical coordinate systems and reference ellipsoids.

Unless otherwise stated functions are ported from atmlab-2-3-181.

class typhon.geodesy.ellipsoidmodels[source]

Bases: object

Provide data for different reference ellipsoids.

The following models are covered:

  • SphericalEarth (radius set as constants.earth_radius)
  • WGS84
  • SphericalVenus (radius same as used in ARTS)
  • SphericalMars (radius same as used in ARTS)
  • EllipsoidMars
  • SphericalJupiter (radius same as used in ARTS)

Examples

>>> e = ellipsoidmodels()
get(model='WGS84')[source]

Return data for different reference ellipsoids.

Parameters:model (str) – Model ellipsoid.
Returns:Equatorial radius (r), eccentricity (e)
Return type:tuple

Examples

>>> e['WGS84']
(6378137, 0.0818191908426)
>>> e.get('WGS84')
(6378137, 0.0818191908426)
>>> ellipsoidmodel()['WGS84']
(6378137, 0.0818191908426)
models

List of available models.

Examples

>>> e.models
['SphericalVenus',
 'SphericalMars',
 'WGS84',
 'SphericalJupiter',
 'EllipsoidMars',
 'SphericalEarth']
typhon.geodesy.ellipsoid_r_geodetic(ellipsoid, lat)[source]

Geodetic radius of a reference ellipsoid.

Gives the distance from the Earth’s centre and the reference ellipsoid as a function of geodetic latitude.

Note

To obtain the radii for geocentric latitude, use ellipsoid_r_geocentric().

Parameters:
Returns:

Radii.

Return type:

float or ndarray

typhon.geodesy.ellipsoid_r_geocentric(ellipsoid, lat)[source]

Geocentric radius of a reference ellipsoid.

Gives the distance from the Earth’s centre and the reference ellipsoid as a function of geocentric latitude.

Note

To obtain the radii for geodetic latitude, use ellipsoid_r_geodetic().

Parameters:
Returns:

Radii.

Return type:

float or ndarray

typhon.geodesy.ellipsoid2d(ellipsoid, orbitinc)[source]

Approximate ellipsoid for 2D calculations.

Determines an approximate reference ellipsoid following an orbit track. The new ellipsoid is determined simply, by determining the radius at the maximum latitude and from this value calculate a new eccentricity. The orbit is specified by giving the orbit inclination, that is normally a value around 100 deg for polar sun-synchronous orbits.

Parameters:
Returns:

Modified ellipsoid vector.

Return type:

tuple

typhon.geodesy.ellipsoidcurvradius(ellipsoid, lat_gd, azimuth)[source]

Sets ellispoid to local curvature radius

Calculates the curvature radius for the given latitude and azimuth angle, and uses this to set a spherical reference ellipsoid suitable for 1D calculations. The curvature radius is a better local approximation than using the local ellipsoid radius.

For exact result the geodetic latitude shall be used.

Parameters:
  • lat_gd – Geodetic latitude.
  • azimuth – Azimuthal angle (angle from NS plane). If given curvature radii are returned, see above.
Returns:

Modified ellipsoid.

Return type:

tuple

typhon.geodesy.sind(x)[source]

Sine of argument in degrees.

typhon.geodesy.cosd(x)[source]

Cosine of argument in degrees.

typhon.geodesy.tand(x)[source]

Tangent of argument in degrees.

typhon.geodesy.asind(x)[source]

Inverse sine in degrees.

typhon.geodesy.inrange(x, minx, maxx, exclude='none', text=None)[source]

Test if x is within given bounds.

Parameters:
  • x – Variable to test.
  • minx – Lower boundary.
  • maxx – Upper boundary.
  • exclude (str) – Exclude boundaries. Possible values are: ‘none’, ‘lower’, ‘upper’ and ‘both’
  • text (str) – Addiitional warning text.
Raises:

Exception – If value is out of bounds.

typhon.geodesy.cart2geocentric(x, y, z, lat0=None, lon0=None, za0=None, aa0=None)[source]

Convert cartesian position to spherical coordinates.

The geocentric Cartesian coordinate system is fixed with respect to the Earth, with its origin at the center of the ellipsoid and its X-, Y-, and Z-axes intersecting the surface at the following points:

  • X-axis: Equator at the Prime Meridian (0°, 0°)
  • Y-axis: Equator at 90-degrees East (0°, 90°)
  • Z-axis: North Pole (90°, 0°)

A common synonym is Earth-Centered, Earth-Fixed coordinates, or ECEF.

If the optional arguments are given, it is ensured that latitude and longitude are kept constant for zenith or nadir cases, and the longitude for N-S cases. The optional input shall be interpreted as the [x,y,z] is obtained by moving from [lat0,lon0] in the direction of [za0,aa0].

Parameters:
  • x – Coordinate in x dimension.
  • y – Coordinate in y dimension.
  • z – Coordinate in z dimension.
  • lat0 – Original latitude.
  • lon0 – Original longitude.
  • za0 – Orignal zenith angle.
  • aa0 – Orignal azimuth angle.
Returns:

Radius, Latitude, Longitude

Return type:

tuple

typhon.geodesy.geocentric2cart(r, lat, lon)[source]

Convert from spherical coordinate to a cartesian position.

See cart2geocentric() for a defintion of the geocentric coordinate system.

Parameters:
r: Radius. lat: Latitude in degree. lon Longitude in degree.
Returns:
tuple: Coordinate in x, y, z dimension.
typhon.geodesy.cart2geodetic(x, y, z, ellipsoid=None)[source]

Converts from cartesian to geodetic coordinates.

The geodetic coordinates refer to the reference ellipsoid specified by input ellipsoid. See module docstring for a defintion of the geocentric coordinate system.

Parameters:
  • x – Coordinates in x dimension.
  • y – Coordinates in y dimension.
  • z – Coordinates in z dimension.
  • ellipsoid – A tuple with the form (semimajor axis, eccentricity). Default is ‘WGS84’ from ellipsoidmodels.
Returns:

Geodetic height, latitude and longitude

Return type:

tuple

typhon.geodesy.geodetic2cart(h, lat, lon, ellipsoid=None)[source]

Converts from geodetic to geocentric cartesian coordinates.

The geodetic coordinates refer to the reference ellipsoid specified by input ellipsoid. See module docstring for a defintion of the geocentric coordinate system.

Parameters:
  • h – Geodetic height (height above the reference ellipsoid).
  • lat – Geodetic latitude.
  • lon – Geodetic longitude.
  • ellipsoid – A tuple with the form (semimajor axis, eccentricity). Default is ‘WGS84’ from ellipsoidmodels.
Returns:

x, y, z coordinates.

Return type:

tuple

typhon.geodesy.geodetic2geocentric(h, lat, lon, ellipsoid=None, **kwargs)[source]

Converts from geodetic to geocentric coordinates.

The geodetic coordinates refer to the reference ellipsoid specified by input ellipsoid. See module docstring for a defintion of the geocentric coordinate system.

Parameters:
  • h – Geodetic height (height above the reference ellipsoid).
  • lat – Geodetic latitude.
  • lon – Geodetic longitude.
  • kwargs – Additional keyword arguments for cart2geocentric().
  • ellipsoid – A tuple with the form (semimajor axis, eccentricity). Default is ‘WGS84’ from ellipsoidmodels.
Returns:

Radius, geocentric latiude, geocentric longitude

Return type:

tuple

typhon.geodesy.geocentric2geodetic(r, lat, lon, ellipsoid=None)[source]

Converts from geocentric to geodetic coordinates.

The geodetic coordinates refer to the reference ellipsoid specified by input ellipsoid. See module docstring for a defintion of the geocentric coordinate system.

Returns:

Geodetic height, latitude and longitude

Return type:

tuple

Parameters:
  • r – Radius:
  • lat – Geocentric latitude.
  • lon – Geocentric longitude.
  • ellipsoid – A tuple with the form (semimajor axis, eccentricity). Default is ‘WGS84’ from ellipsoidmodels.
typhon.geodesy.great_circle_distance(lat1, lon1, lat2, lon2, r=None)[source]

Calculate the distance between two geograpgical positions.

“As-the-crow-flies” distance between two points, specified by their latitude and longitude.

If the optional argument r is given, the distance in m is returned. Otherwise the angular distance in degrees is returned.

Parameters:
  • lat1 – Latitude of position 1.
  • lon1 – Longitude of position 1.
  • lat2 – Latitude of position 2.
  • lon2 – Longitude of position 2.
  • r (float) – The radius (common for both points).
Returns:

Distance, either in degress or m.

typhon.geodesy.geographic_mean(lat, lon, h=0, ellipsoid=None)[source]

Calculate mean position for set of coordinates.

Parameters:
  • lat (float or ndarray) – Latitudes in degrees.
  • lon (float or ndarray) – Longitudes in degrees.
  • h (float or ndarray) – Optiional altitude for each coordinate (default is).
  • ellipsoid – A tuple with the form (semimajor axis, eccentricity). Default is ‘WGS84’ from ellipsoidmodels.
Returns:

Mean latitudes and longitudes in degrees.

Return type:

tuple

typhon.geodesy.cartposlos2geocentric(x, y, z, dx, dy, dz, ppc=None, lat0=None, lon0=None, za0=None, aa0=None)[source]

Converts cartesian POS/LOS to spherical coordinates.

Position is given as (x,y,z), while line-of-sight is given as (dx,dy,dz). The corresponding quantities in polar coordinates are (r,lat,lon) and (za,aa), respectively.

See Contents for defintion of coordinate systems.

If the optional arguments are given, it is ensured that latitude and longitude are kept constant for zenith or nadir cases, and the longitude and azimuth angle for N-S cases. The optional input shall be interpreted as the [x,y,z] is obtained by moving from [r0,lat0,lon0] in the direction of [za0,aa0].

This version is different from the atmlab version by normalizing the los- vector and demanding all or nothing for the optional arguments to work.

Parameters:
  • x – Coordinate in x dimension.
  • y – Coordinate in y dimension.
  • z – Coordinate in z dimension.
  • dx – LOS component in x dimension.
  • dy – LOS component in y dimension.
  • dz – LOS component in z dimension.
  • ppc – Propagation path constant = r0*sin(za0).
  • lat0 – Original latitude.
  • lon0 – Original longitude.
  • za0 – Orignal zenith angle.
  • aa0 – Orignal azimuth angle.
Returns:

Radius, Latitude, Longitude,

Zenith angle, Azimuth angle

Return type:

tuple(ndarray)

typhon.geodesy.geocentricposlos2cart(r, lat, lon, za, aa)[source]

GEOCENTRICPOSLOS2CART converts from spherical POS/LOS to cartesian POS/LOS

See Contents for a defintion of the geocentric coordinate system. The local LOS angles are defined following the EAST-NORTH-UP system:

za aa

90 0 points towards north

90 90 points towards east

0 aa points up

FORMAT [x,y,z,dx,dy,dz]=geocentricposlos2cart(r,lat,lon,za,az)

OUT

x Coordinate in x dimension

y Coordinate in y dimension

z Coordinate in z dimension

dx LOS component in x dimension

dy LOS component in y dimension

dz LOS

IN

r Radius

lat Latitude

lon Longitude

za zenith angle

aa azimuth angle

Ported from atmlab. Original author: Bengt Rydberg 2011-10-31