pastax.geo

Geographic utilities: unit conversions and Haversine distance.

chordal_projection¶

pastax.geo.chordal_projection(x)

Project [lon, lat] points onto a sphere of radius EARTH_RADIUS.

The projection is defined as:

$$x = R \cos(\varphi) \cos(\lambda) \\ y = R \cos(\varphi) \sin(\lambda) \\ z = R \sin(\varphi)$$

(1)

where $R$ is the Earth’s radius, $\varphi$ is the latitude, and $\lambda$ is the longitude (in radians).

Parameters

x (Float[jaxlib._jax.Array, '... 2']) – Point(s) [lon, lat] in degrees, shape (..., 2).

Returns

Point(s) [x, y, z] in metres, shape (..., 3).

Return type

Float[jaxlib._jax.Array, ‘… 3’]

EARTH_RADIUS¶

pastax.geo.EARTH_RADIUS: float = 6371008.8

Mean Earth radius in metres (IUGG 2015 mean radius).

haversine¶

pastax.geo.haversine(x1, x2)

Great-circle distance between [lon, lat] points.

Uses the spherical haversine formula with EARTH_RADIUS as the sphere radius $R$:

$$a = \sin^2\!\left(\tfrac{\Delta\varphi}{2}\right) + \cos\varphi_1 \cos\varphi_2 \sin^2\!\left(\tfrac{\Delta\lambda}{2}\right)$$

(2)

$$d = 2 R \arcsin\!\sqrt{a}$$

(3)

where $\varphi$ is latitude and $\lambda$ is longitude (in radians). The last axis of each input must have size 2 (lon, lat); leading axes broadcast under standard NumPy/JAX rules.

Parameters

• x1 (Float[jaxlib._jax.Array, '... 2']) – First point(s) [lon, lat] in degrees, shape (..., 2).

• x2 (Float[jaxlib._jax.Array, '... 2']) – Second point(s) [lon, lat] in degrees, shape (..., 2).

Returns

Great-circle distance in metres, with shape matching the broadcast of the leading axes of x1 and x2.

Return type

Float[jaxlib._jax.Array, ‘…’]

meters_to_degrees¶

pastax.geo.meters_to_degrees(disp_m, lat_deg)

Convert an [east, north] displacement in metres to [dlon, dlat] in degrees.

Uses a flat-Earth approximation around lat_deg: the meridional component is converted via EARTH_RADIUS; the zonal component is additionally divided by $\cos(\mathrm{lat})$ to account for shrinking longitude circles toward the poles.

The division uses pastax.safe_divide(), so an exactly zero $\cos(\mathrm{lat})$ yields 0 (with a finite gradient) rather than inf/nan. Note this is a corner-case guard only: the flat-Earth longitude scaling genuinely diverges at the poles, and near (but not at) $\pm 90^\circ$ — where $\cos(\mathrm{lat})$ is tiny but nonzero — the zonal component is still legitimately huge. Treat results within a cell of the poles as unreliable.

Parameters

• disp_m (Float[jaxlib._jax.Array, '... 2']) – Displacement(s) [east, north] in metres. The last axis must have size 2; leading axes are passed through unchanged.

• lat_deg (Float[jaxlib._jax.Array, '']) – Reference latitude in degrees, used for the longitude scaling.

Returns

Same shape as disp_m, but expressed as [dlon, dlat] in degrees.

Return type

Float[jaxlib._jax.Array, ‘… 2’]

degrees_to_meters¶

pastax.geo.degrees_to_meters(disp_deg, lat_deg)

Convert a [dlon, dlat] displacement in degrees to [east, north] in metres.

Inverse of meters_to_degrees(). Uses a flat-Earth approximation around lat_deg.

Parameters

• disp_deg (Float[jaxlib._jax.Array, '... 2']) – Displacement(s) [dlon, dlat] in degrees. The last axis must have size 2; leading axes are passed through unchanged.

• lat_deg (Float[jaxlib._jax.Array, '']) – Reference latitude in degrees, used for the longitude scaling.

Returns

Same shape as disp_deg, but expressed as [east, north] in metres.

Return type

Float[jaxlib._jax.Array, ‘… 2’]

wrap_longitude¶

pastax.geo.wrap_longitude(lon, period=360.0, lower=-180.0)

Wrap longitude(s) into the half-open window [lower, lower + period).

The solver integrates position without normalising it, so a particle that drifts across the antimeridian accumulates an unbounded longitude (181°, 200°, …) — which is the correct continuous representation for a trajectory. This is a post-processing / display helper that folds such longitudes back into a canonical window. The defaults give the [-180, 180) convention; pass lower=0.0 for [0, 360).

Element-wise and shape-preserving; pure arithmetic, so it is safe under jit / vmap / grad. For a [lon, lat] trajectory of shape (..., 2), wrap only the longitude column, e.g.:

traj = traj.at[..., 0].set(wrap_longitude(traj[..., 0]))

Parameters

• lon (Float[jaxlib._jax.Array, '...']) – Longitude value(s) in degrees, any shape.

• period (float) – Longitude period in degrees (default 360.0).

• lower (float) – Lower edge of the target window (default -180.0); the window is [lower, lower + period).

Returns

Longitudes folded into [lower, lower + period), same shape as lon.

Return type

Float[jaxlib._jax.Array, ‘…’]