orthodrome

path_contains_points(verts, points)[source]
float_array_broadcast(*args)[source]
class Loc(lat, lon)[source]

Simple location representation

Attrib lat:Latitude degree
Attrib lon:Longitude degree
clip(x, mi, ma)[source]

Clipping data array x

Parameters:
Returns:

Clipped data

Return type:

numpy.ndarray

wrap(x, mi, ma)[source]

Wrapping continuous data to fundamental phase values.

x_{\mathrm{wrapped}} = x_{\mathrm{cont},i} -             \frac{ x_{\mathrm{cont},i} - r_{\mathrm{min}} }                  { r_{\mathrm{max}} -  r_{\mathrm{min}}}             \cdot ( r_{\mathrm{max}} -  r_{\mathrm{min}}),\ \quad
x_{\mathrm{wrapped}}\; \in
    \;[ r_{\mathrm{min}},\, r_{\mathrm{max}}].

Parameters:
  • x (numpy.ndarray) – Continunous data to be wrapped
  • mi (float) – Minimum value of wrapped data
  • ma (float) – Maximum value of wrapped data
Returns:

Wrapped data

Return type:

numpy.ndarray

cosdelta(*args)[source]

Cosine of the angular distance between two points a and b on a sphere.

This function (find implementation below) returns the cosine of the distance angle ‘delta’ between two points a and b, coordinates of which are expected to be given in geographical coordinates and in degrees. For numerical stability a maximum of 1.0 is enforced.

A_{\mathrm{lat'}} = \frac{ \pi}{180} \cdot A_{lat}, \quad         A_{\mathrm{lon'}} = \frac{ \pi}{180} \cdot A_{lon}, \quad         B_{\mathrm{lat'}} = \frac{ \pi}{180} \cdot B_{lat}, \quad         B_{\mathrm{lon'}} = \frac{ \pi}{180} \cdot B_{lon}\\[0.5cm]

\cos(\Delta) = \min( 1.0, \quad  \sin( A_{\mathrm{lat'}})                      \sin( B_{\mathrm{lat'}} ) +                      \cos(A_{\mathrm{lat'}})  \cos( B_{\mathrm{lat'}} )                      \cos( B_{\mathrm{lon'}} - A_{\mathrm{lon'}} )

Parameters:
Returns:

cosdelta

Return type:

float

cosdelta_numpy(a_lats, a_lons, b_lats, b_lons)[source]

Cosine of the angular distance between two points a and b on a sphere.

This function returns the cosines of the distance angles delta between two points a and b given as numpy.ndarray. The coordinates are expected to be given in geographical coordinates and in degrees. For numerical stability a maximum of 1.0 is enforced.

Please find the details of the implementation in the documentation of the function pyrocko.orthodrome.cosdelta() above.

Parameters:
Returns:

cosdelta

azimuth(*args)[source]

Azimuth calculation

This function (find implementation below) returns azimuth … between points a and b, coordinates of which are expected to be given in geographical coordinates and in degrees.

A_{\mathrm{lat'}} = \frac{ \pi}{180} \cdot A_{lat}, \quad         A_{\mathrm{lon'}} = \frac{ \pi}{180} \cdot A_{lon}, \quad         B_{\mathrm{lat'}} = \frac{ \pi}{180} \cdot B_{lat}, \quad         B_{\mathrm{lon'}} = \frac{ \pi}{180} \cdot B_{lon}\\

\varphi_{\mathrm{azi},AB} = \frac{180}{\pi} \arctan \left[ \frac{                \cos( A_{\mathrm{lat'}}) \cos( B_{\mathrm{lat'}} )                 \sin(B_{\mathrm{lon'}} - A_{\mathrm{lon'}} )}                 {\sin ( B_{\mathrm{lat'}} ) - \sin( A_{\mathrm{lat'}}                   cosdelta) } \right]

Parameters:
Returns:

Azimuth in degree

azimuth_numpy(a_lats, a_lons, b_lats, b_lons, _cosdelta=None)[source]

Calculation of the azimuth (track angle) from a location A towards B.

This function returns azimuths (track angles) from locations A towards B given in numpy.ndarray. Coordinates are expected to be given in geographical coordinates and in degrees.

Please find the details of the implementation in the documentation of the function pyrocko.orthodrome.azimuth().

Parameters:
Returns:

Azimuths in degrees

Return type:

numpy.ndarray, (N)

azibazi(*args, **kwargs)[source]
azibazi_numpy(a_lats, a_lons, b_lats, b_lons, implementation='c')[source]
azidist_numpy(*args)[source]

Calculation of the azimuth (track angle) and the distance from locations A towards B on a sphere.

The assisting functions used are pyrocko.orthodrome.cosdelta() and pyrocko.orthodrome.azimuth()

Parameters:
Returns:

Azimuths in degrees, distances in degrees

Return type:

numpy.ndarray, (2xN)

distance_accurate50m(*args, **kwargs)[source]

Accurate distance calculation based on a spheroid of rotation.

Function returns distance in meter between points A and B, coordinates of which must be given in geographical coordinates and in degrees. The returned distance should be accurate to 50,m using WGS84. Values for the Earth’s equator radius and the Earth’s oblateness (f_oblate) are defined in the pyrocko configuration file pyrocko.config.

From wikipedia (http://de.wikipedia.org/wiki/Orthodrome), based on:

Meeus, J.: Astronomical Algorithms, S 85, Willmann-Bell, Richmond 2000 (2nd ed., 2nd printing), ISBN 0-943396-61-1

F = \frac{\pi}{180}                             \frac{(A_{lat} + B_{lat})}{2}, \quad         G = \frac{\pi}{180}                              \frac{(A_{lat} - B_{lat})}{2}, \quad         l = \frac{\pi}{180}                             \frac{(A_{lon} - B_{lon})}{2} \quad         \\[0.5cm]
S = \sin^2(G) \cdot \cos^2(l) +               \cos^2(F) \cdot \sin^2(l), \quad \quad
C = \cos^2(G) \cdot \cos^2(l) +                   \sin^2(F) \cdot \sin^2(l)

w = \arctan \left( \sqrt{ \frac{S}{C}} \right) , \quad         r =  \sqrt{\frac{S}{C} }

The spherical-earth distance D between A and B, can be given with:

D_{sphere} = 2w \cdot R_{equator}

The oblateness of the Earth requires some correction with correction factors h1 and h2:

h_1 = \frac{3r - 1}{2C}, \quad         h_2 = \frac{3r +1 }{2S}\\[0.5cm]

D = D_{\mathrm{sphere}} \cdot [ 1 + h_1 \,f_{\mathrm{oblate}}                  \cdot \sin^2(F)                  \cos^2(G) - h_2\, f_{\mathrm{oblate}}                  \cdot \cos^2(F) \sin^2(G)]

Parameters:
Returns:

Distance in meter

Return type:

float

distance_accurate50m_numpy(a_lats, a_lons, b_lats, b_lons, implementation='c')[source]

Accurate distance calculation based on a spheroid of rotation.

Function returns distance in meter between points a and b, coordinates of which must be given in geographical coordinates and in degrees. The returned distance should be accurate to 50,m using WGS84. Values for the Earth’s equator radius and the Earth’s oblateness (f_oblate) are defined in the pyrocko configuration file pyrocko.config.

From wikipedia (http://de.wikipedia.org/wiki/Orthodrome), based on:

Meeus, J.: Astronomical Algorithms, S 85, Willmann-Bell, Richmond 2000 (2nd ed., 2nd printing), ISBN 0-943396-61-1

F_i = \frac{\pi}{180}                            \frac{(a_{lat,i} + a_{lat,i})}{2}, \quad         G_i = \frac{\pi}{180}                             \frac{(a_{lat,i} - b_{lat,i})}{2}, \quad         l_i= \frac{\pi}{180}                             \frac{(a_{lon,i} - b_{lon,i})}{2} \\[0.5cm]
S_i = \sin^2(G_i) \cdot \cos^2(l_i) +               \cos^2(F_i) \cdot \sin^2(l_i), \quad \quad
C_i = \cos^2(G_i) \cdot \cos^2(l_i) +                   \sin^2(F_i) \cdot \sin^2(l_i)

w_i = \arctan \left( \sqrt{\frac{S_i}{C_i}} \right), \quad         r_i =  \sqrt{\frac{S_i}{C_i} }

The spherical-earth distance D between a and b, can be given with:

D_{\mathrm{sphere},i} = 2w_i \cdot R_{\mathrm{equator}}

The oblateness of the Earth requires some correction with correction factors h1 and h2:

h_{1.i} = \frac{3r - 1}{2C_i}, \quad         h_{2,i} = \frac{3r +1 }{2S_i}\\[0.5cm]

D_{AB,i} = D_{\mathrm{sphere},i} \cdot [1 + h_{1,i}
    \,f_{\mathrm{oblate}}             \cdot \sin^2(F_i)             \cos^2(G_i) - h_{2,i}\, f_{\mathrm{oblate}}             \cdot \cos^2(F_i) \sin^2(G_i)]

Parameters:
Returns:

Distances in degrees

Return type:

numpy.ndarray, (N)

ne_to_latlon(lat0, lon0, north_m, east_m)[source]

Transform local cartesian coordinates to latitude and longitude.

From east and north coordinates (x and y coordinate numpy.ndarray) relative to a reference differences in longitude and latitude are calculated, which are effectively changes in azimuth and distance, respectively:

\text{distance change:}\; \Delta {\bf{a}} &= \sqrt{{\bf{y}}^2 +                                            {\bf{x}}^2 }/ \mathrm{R_E},

\text{azimuth change:}\; \Delta \bf{\gamma} &= \arctan( \bf{x}                                                           / \bf{y}).

The projection used preserves the azimuths of the input points.

Parameters:
  • lat0 (float) – Latitude origin of the cartesian coordinate system.
  • lon0 (float) – Longitude origin of the cartesian coordinate system.
  • north_m (numpy.ndarray, (N)) – Northing distances from origin in meters.
  • east_m (numpy.ndarray, (N)) – Easting distances from origin in meters.
Returns:

Array with latitudes and longitudes

Return type:

numpy.ndarray, (2xN)

azidist_to_latlon(lat0, lon0, azimuth_deg, distance_deg)[source]

(Durchreichen??).

azidist_to_latlon_rad(lat0, lon0, azimuth_rad, distance_rad)[source]

Absolute latitudes and longitudes are calculated from relative changes.

For numerical stability a range between of -1.0 and 1.0 is enforced for c and alpha.

\Delta {\bf a}_i \; \text{and} \; \Delta \gamma_i \;             \text{are relative distances and azimuths from lat0 and lon0 for                 \textit{i} source points of a finite source.}

\mathrm{b} &= \frac{\pi}{2} -\frac{\pi}{180} \;\mathrm{lat_0}\\
{\bf c}_i &=\arccos[\; \cos(\Delta {\bf{a}}_i) \cos(\mathrm{b})               + |\Delta \gamma_i| \,\sin(\Delta {\bf a}_i)
        \sin(\mathrm{b})\; ] \\
\mathrm{lat}_i &=  \frac{180}{\pi}                            \left(\frac{\pi}{2} - {\bf c}_i \right)

\alpha_i &= \arcsin \left[ \; \frac{ \sin(\Delta {\bf a}_i )                          \sin(|\Delta \gamma_i|)}{\sin({\bf c}_i)}\;
                 \right] \\
\alpha_i &= \begin{cases}
         \alpha_i, &\text{if}  \; \cos(\Delta {\bf a}_i) -                  \cos(\mathrm{b}) \cos({\bf{c}}_i) > 0, \;                  \text{else} \\
         \pi - \alpha_i, & \text{if} \; \alpha_i > 0,\;
         \text{else}\\
        -\pi - \alpha_i, & \text{if} \; \alpha_i < 0.
        \end{cases} \\
\mathrm{lon}_i &=   \mathrm{lon_0} +              \frac{180}{\pi} \, \frac{\Delta \gamma_i }{|\Delta \gamma_i|}                             \cdot \alpha_i                               \text{, with $\alpha_i \in [-\pi,\pi]$}

Parameters:
  • lat0 (float) – Latitude origin of the cartesian coordinate system.
  • lon0 (float) – Longitude origin of the cartesian coordinate system.
  • distance_rad (numpy.ndarray, (N)) – Distances from origin in radians.
  • azimuth_rad (numpy.ndarray, (N)) – Azimuth from radians.
Returns:

Array with latitudes and longitudes

Return type:

numpy.ndarray, (2xN)

ne_to_latlon_alternative_method(lat0, lon0, north_m, east_m)[source]

Transform local cartesian coordinates to latitude and longitude.

Like pyrocko.orthodrome.ne_to_latlon(), but this method (implementation below), although it should be numerically more stable, suffers problems at points which are across the pole as seen from the cartesian origin.

\text{distance change:}\; \Delta {{\bf a}_i} &=\sqrt{{\bf{y}}^2_i +
                                {\bf{x}}^2_i }/ \mathrm{R_E},\\
\text{azimuth change:}\; \Delta {\bf \gamma}_i &=
                                \arctan( {\bf x}_i {\bf y}_i). \\
\mathrm{b} &= \frac{\pi}{2} -\frac{\pi}{180} \;\mathrm{lat_0}\\

{{\bf z}_1}_i &= \cos{\left( \frac{\Delta {\bf a}_i -                             \mathrm{b}}{2} \right)}                             \cos {\left( \frac{|\gamma_i|}{2} \right) }\\
{{\bf n}_1}_i &= \cos{\left( \frac{\Delta {\bf a}_i +                             \mathrm{b}}{2} \right)}                            \sin {\left( \frac{|\gamma_i|}{2} \right) }\\
{{\bf z}_2}_i &= \sin{\left( \frac{\Delta {\bf a}_i -                             \mathrm{b}}{2} \right)}                            \cos {\left( \frac{|\gamma_i|}{2} \right) }\\
{{\bf n}_2}_i &= \sin{\left( \frac{\Delta {\bf a}_i +                             \mathrm{b}}{2} \right)}                            \sin {\left( \frac{|\gamma_i|}{2} \right) }\\
{{\bf t}_1}_i &= \arctan{\left( \frac{{{\bf z}_1}_i}                                     {{{\bf n}_1}_i} \right) }\\
{{\bf t}_2}_i &= \arctan{\left( \frac{{{\bf z}_2}_i}                                     {{{\bf n}_2}_i} \right) } \\[0.5cm]
c &= \begin{cases}
      2 \cdot \arccos \left( {{\bf z}_1}_i / \sin({{\bf t}_1}_i)                               \right),\; \text{if }                               |\sin({{\bf t}_1}_i)| >                                 |\sin({{\bf t}_2}_i)|,\; \text{else} \\
      2 \cdot \arcsin{\left( {{\bf z}_2}_i /                                  \sin({{\bf t}_2}_i) \right)}.
     \end{cases}\\

{\bf {lat}}_i  &= \frac{180}{ \pi } \left( \frac{\pi}{2}                                               - {\bf {c}}_i \right) \\
{\bf {lon}}_i &=  {\bf {lon}}_0 + \frac{180}{ \pi }                                       \frac{\gamma_i}{|\gamma_i|},                                      \text{ with}\; \gamma \in [-\pi,\pi]

Parameters:
  • lat0 (float) – Latitude origin of the cartesian coordinate system.
  • lon0 (float) – Longitude origin of the cartesian coordinate system.
  • north_m (numpy.ndarray, (N)) – Northing distances from origin in meters.
  • east_m (numpy.ndarray, (N)) – Easting distances from origin in meters.
Returns:

Array with latitudes and longitudes

Return type:

numpy.ndarray, (2xN)

latlon_to_ne(*args)[source]

Relative cartesian coordinates with respect to a reference location.

For two locations, a reference location A and another location B, given in geographical coordinates in degrees, the corresponding cartesian coordinates are calculated. Assisting functions are pyrocko.orthodrome.azimuth' and :py:func:`pyrocko.orthodrome.distance_accurate50m()

D_{AB} &= \mathrm{distance\_accurate50m(}A, B \mathrm{)}, \quad                         \varphi_{\mathrm{azi},AB} = \mathrm{azimuth(}A,B                             \mathrm{)}\\[0.3cm]

n &= D_{AB} \cdot \cos( \frac{\pi }{180}                                     \varphi_{\mathrm{azi},AB} )\\
e &= D_{AB} \cdot \sin( \frac{\pi }{180} \varphi_{\mathrm{azi},AB})

Parameters:
Returns:

Northing and easting from refloc to location

Return type:

tuple, float

latlon_to_ne_numpy(lat0, lon0, lat, lon)[source]

Relative cartesian coordinates with respect to a reference location.

For two locations, a reference location (lat0, lon0) and another location B, given in geographical coordinates in degrees, the corresponding cartesian coordinates are calculated. Assisting functions are azimuth() and distance_accurate50m().

Parameters:
  • lat0 – reference location latitude
  • lon0 – reference location longitude
  • lat – absolute location latitude
  • lon – absolute location longitude
Returns:

(n, e): relative north and east positions

Return type:

numpy.ndarray, (2xN)

Implemented formulations:

D_{AB} &= \mathrm{distance\_accurate50m(}A, B \mathrm{)}, \quad            \varphi_{\mathrm{azi},AB} = \mathrm{azimuth(}A,B                                                 \mathrm{)}\\[0.3cm]

n &= D_{AB} \cdot \cos( \frac{\pi }{180} \varphi_{                                                 \mathrm{azi},AB} )\\
e &= D_{AB} \cdot \sin( \frac{\pi }{180} \varphi_{                                                 \mathrm{azi},AB} )

get_wgs84()[source]
amap(n)[source]
ne_to_latlon2(*args)[source]
latlon_to_ne2(*args)[source]
distance_accurate15nm(*args)[source]
positive_region(region)[source]

Normalize parameterization of a rectangular geographical region.

Parameters:region(west, east, south, north)
Returns:(west, east, south, north), where west <= east and where west and east are in the range [-180., 180.+360.[
points_in_region(p, region)[source]

Check what points are contained in a rectangular geographical region.

Parameters:
  • p – NumPy array of shape (N, 2) where each row is a (lat, lon) pair [deg]
  • region(west, east, south, north) [deg]
Returns:

NumPy array of shape (N), type bool

point_in_region(p, region)[source]

Check if a point is contained in a rectangular geographical region.

Parameters:
  • p(lat, lon) [deg]
  • region(west, east, south, north) [deg]
Returns:

bool

radius_to_region(lat, lon, radius)[source]

Get a rectangular region which fully contains a given circular region.

Parameters:
  • lat,lon – center of circular region [deg]
  • radius – radius of circular region [m]
Returns:

rectangular region as (east, west, south, north) [deg]

geographic_midpoint(lats, lons, weights=None)[source]

Calculate geographic midpoints by finding the center of gravity.

This method suffers from instabilities if points are centered around the poles.

Parameters:
Returns:

Latitudes and longitudes of the modpoints

Return type:

numpy.ndarray, (2xN)

geographic_midpoint_locations(locations, weights=None)[source]
geodetic_to_ecef(lat, lon, alt)[source]

Convert geodetic coordinates to Earth-Centered, Earth-Fixed (ECEF) Cartesian coordinates. [3] [4]

Parameters:
  • lat (float) – Geodetic latitude in [deg].
  • lon (float) – Geodetic longitude in [deg].
  • alt (float) – Geodetic altitude (height) in [m] (positive for points outside the geoid).
Returns:

ECEF Cartesian coordinates (X, Y, Z) in [m].

Return type:

tuple, float

[3]https://en.wikipedia.org/wiki/ECEF
[4]https://en.wikipedia.org/wiki/Geographic_coordinate_conversion #From_geodetic_to_ECEF_coordinates
ecef_to_geodetic(X, Y, Z)[source]

Convert Earth-Centered, Earth-Fixed (ECEF) Cartesian coordinates to geodetic coordinates (Ferrari’s solution).

Parameters:Y, Z (X,) – Cartesian coordinates in ECEF system in [m].
Returns:Geodetic coordinates (lat, lon, alt). Latitude and longitude are in [deg] and altitude is in [m] (positive for points outside the geoid).
Return type:tuple, float

See also

https://en.wikipedia.org/wiki/Geographic_coordinate_conversion #The_application_of_Ferrari.27s_solution

exception Farside[source]
latlon_to_xyz(latlons)[source]
xyz_to_latlon(xyz)[source]
rot_to_00(lat, lon)[source]
distances3d(a, b)[source]
circulation(points2)[source]
stereographic(points)[source]
stereographic_poly(points)[source]
contains_point(polygon, point)[source]
contains_points(polygon, points)[source]