Geodetic Coordinates¶

This tutorial is about the satkit.itrfcoord data type — how to describe a specific point on or near the Earth's surface and convert between geodetic and Cartesian representations of that point.

For the abstract definition of the ITRF reference frame, and for rotating vectors between ITRF and other frames (GCRF, TEME, ...), see the Coordinate Frames tutorial.

All points here live in the International Terrestrial Reference Frame (ITRF) — the Earth-fixed, Earth-rotating frame used in geodesy and satellite operations. satkit.itrfcoord represents a position in this frame and provides conversions between:

• Geodetic coordinates: latitude, longitude, and altitude above the WGS-84 ellipsoid (ellipsoid-normal, not geocentric latitude)
• Cartesian coordinates: a 3-element (x, y, z) vector in meters with origin at Earth's center of mass
• Local tangent planes: East-North-Up (ENU) and North-East-Down (NED) for range/bearing-style computations near a reference point
• Geodesic distances: great-circle-equivalent distances along the WGS-84 ellipsoid between two itrfcoord values

This tutorial demonstrates constructing coordinates from geodetic parameters, converting to and from Cartesian vectors, and computing distances between locations.

In [1]:

import satkit as sk

# Specify a coordinate via Geodetic latitude & longitude (if altitude is left out, default is 0
boston = sk.itrfcoord(latitude_deg=42.14, longitude_deg=-71.15)

# Denver, CO is at a higher altitude; altitude is specified in meters
denver = sk.itrfcoord(latitude_deg=32.7392, longitude_deg=-104.99, altitude=1600)
print(f'Denver = {denver}')

# Get the Cartesian vector for Denver
print(f'Denver vector coordinate is {denver.vector}')

# Create a duplicate of the Denver coordinate, using vector as input
denver2 = sk.itrfcoord(denver.vector)
print(f'Denver created from Cartesian is {denver2}')

import satkit as sk # Specify a coordinate via Geodetic latitude & longitude (if altitude is left out, default is 0 boston = sk.itrfcoord(latitude_deg=42.14, longitude_deg=-71.15) # Denver, CO is at a higher altitude; altitude is specified in meters denver = sk.itrfcoord(latitude_deg=32.7392, longitude_deg=-104.99, altitude=1600) print(f'Denver = {denver}') # Get the Cartesian vector for Denver print(f'Denver vector coordinate is {denver.vector}') # Create a duplicate of the Denver coordinate, using vector as input denver2 = sk.itrfcoord(denver.vector) print(f'Denver created from Cartesian is {denver2}')

Denver = ITRFCoord(lat:  32.7392 deg, lon: -104.9900 deg, hae:  1.60 km)
Denver vector coordinate is [-1389345.60167272 -5188730.62268842  3430531.07678629]
Denver created from Cartesian is ITRFCoord(lat:  32.7392 deg, lon: -104.9900 deg, hae:  1.60 km)

Creating Coordinates¶

Coordinates can be created from geodetic parameters (latitude, longitude, altitude) or from a Cartesian position vector. The altitude defaults to 0 (on the ellipsoid surface) if not specified.

Geodesic Distance¶

The geodesic_distance method computes the shortest path between two points along the surface of the WGS-84 ellipsoid, returning the distance in meters along with the initial and final headings of the path.

In [2]:

import satkit as sk

newyork = sk.itrfcoord(latitude_deg=40.7128, longitude_deg=-74.0060)
destinations = [
    {'name': 'London', 'latitude_deg': 51.5072, 'longitude_deg': -0.1276},
    {'name': 'Paris', 'latitude_deg': 48.8566, 'longitude_deg': 2.3522},
    {'name': 'Toronto', 'latitude_deg': 43.6532, 'longitude_deg': -79.3832},
    {'name': 'Tokyo', 'latitude_deg': 35.6763, 'longitude_deg': 139.65},
    {'name': 'Sau Paulo', 'latitude_deg': -23.5558, 'longitude_deg': -46.6396}
]
for dest in destinations:
    destcoord = sk.itrfcoord(latitude_deg=dest['latitude_deg'], longitude_deg=dest['longitude_deg'])
    # Distance, start, and end heading along great circle to destination coordinate
    dist, _heading_start, _heading_end = newyork.geodesic_distance(destcoord)
    print(f'New York to {dest['name']} distance is {dist/1e3:.0f} km')

import satkit as sk newyork = sk.itrfcoord(latitude_deg=40.7128, longitude_deg=-74.0060) destinations = [ {'name': 'London', 'latitude_deg': 51.5072, 'longitude_deg': -0.1276}, {'name': 'Paris', 'latitude_deg': 48.8566, 'longitude_deg': 2.3522}, {'name': 'Toronto', 'latitude_deg': 43.6532, 'longitude_deg': -79.3832}, {'name': 'Tokyo', 'latitude_deg': 35.6763, 'longitude_deg': 139.65}, {'name': 'Sau Paulo', 'latitude_deg': -23.5558, 'longitude_deg': -46.6396} ] for dest in destinations: destcoord = sk.itrfcoord(latitude_deg=dest['latitude_deg'], longitude_deg=dest['longitude_deg']) # Distance, start, and end heading along great circle to destination coordinate dist, _heading_start, _heading_end = newyork.geodesic_distance(destcoord) print(f'New York to {dest['name']} distance is {dist/1e3:.0f} km')

New York to London distance is 5585 km
New York to Paris distance is 5853 km
New York to Toronto distance is 551 km
New York to Tokyo distance is 10876 km
New York to Sau Paulo distance is 7658 km

Local Tangent Planes and Working with Many Coordinates¶

The to_enu(origin) and to_ned(origin) methods return the position of self relative to origin, expressed in the East-North-Up or North-East-Down tangent plane anchored at origin — the natural frame for look angles, ranges, and bearings from a ground station.

itrfcoord is deliberately a scalar type: one object describes one location, and batched (Nx3) inputs are intentionally not part of the API. For workflows that touch many points at once — trajectories, ground tracks, station networks — keep the data in NumPy arrays and hoist the frame machinery out of the loop. A per-point Python loop works but costs roughly a microsecond per point in interpreter and object overhead, while the underlying math is tens of nanoseconds. Two patterns cover the common cases:

1. ENU / NED for many points relative to one origin. The tangent plane is a single rotation anchored at the origin, so an entire array converts with one matrix multiply — about 15 ms per million points, versus roughly 1 s per million with a per-point loop.
2. Cartesian → geodetic. The WGS-84 conversion is genuinely per-point; construct coordinates in a comprehension (≈0.7 s per million points). If this is ever the bottleneck, subsample or coarsen first — for ground-track plotting, a few hundred points per orbit is plenty.

In [3]:

import numpy as np
import satkit as sk

# A ground station as the local-frame origin
origin = sk.itrfcoord(latitude_deg=42.466, longitude_deg=-71.1516)

# Scalar usage: position of Boston in the station's ENU frame
boston = sk.itrfcoord(latitude_deg=42.14, longitude_deg=-71.15)
print(f'Boston in station ENU frame (m): {boston.to_enu(origin)}')

# --- Many points relative to one origin ---
# Example data: 100,000 ITRF positions within ~100 km of the station
rng = np.random.default_rng(7)
enu_offsets = rng.uniform([-1e5, -1e5, 0.0], [1e5, 1e5, 1e5], (100_000, 3))

# The ENU frame is one rotation: qenu2itrf maps ENU -> ITRF, and its
# conjugate maps ITRF -> ENU. As a matrix, the whole array converts
# in a single multiply -- no per-point itrfcoord objects needed.
R = origin.qenu2itrf.conj.as_rotation_matrix()
vecs_itrf = origin.vector + enu_offsets @ R  # ENU offsets -> ITRF points
enu = (vecs_itrf - origin.vector) @ R.T      # and back again

# Element-wise identical to the scalar API
assert np.allclose(enu[0], sk.itrfcoord(vecs_itrf[0]).to_enu(origin))
# (for NED, use origin.qned2itrf in place of origin.qenu2itrf)

# Cartesian -> geodetic is per-point: use a comprehension
lat, lon, alt = np.array(
    [(c.latitude_deg, c.longitude_deg, c.altitude)
     for c in map(sk.itrfcoord, vecs_itrf)]).T
print(f'latitude range:  [{lat.min():8.3f}, {lat.max():8.3f}] deg')
print(f'altitude range:  [{alt.min()/1e3:8.2f}, {alt.max()/1e3:8.2f}] km')

import numpy as np import satkit as sk # A ground station as the local-frame origin origin = sk.itrfcoord(latitude_deg=42.466, longitude_deg=-71.1516) # Scalar usage: position of Boston in the station's ENU frame boston = sk.itrfcoord(latitude_deg=42.14, longitude_deg=-71.15) print(f'Boston in station ENU frame (m): {boston.to_enu(origin)}') # --- Many points relative to one origin --- # Example data: 100,000 ITRF positions within ~100 km of the station rng = np.random.default_rng(7) enu_offsets = rng.uniform([-1e5, -1e5, 0.0], [1e5, 1e5, 1e5], (100_000, 3)) # The ENU frame is one rotation: qenu2itrf maps ENU -> ITRF, and its # conjugate maps ITRF -> ENU. As a matrix, the whole array converts # in a single multiply -- no per-point itrfcoord objects needed. R = origin.qenu2itrf.conj.as_rotation_matrix() vecs_itrf = origin.vector + enu_offsets @ R # ENU offsets -> ITRF points enu = (vecs_itrf - origin.vector) @ R.T # and back again # Element-wise identical to the scalar API assert np.allclose(enu[0], sk.itrfcoord(vecs_itrf[0]).to_enu(origin)) # (for NED, use origin.qned2itrf in place of origin.qenu2itrf) # Cartesian -> geodetic is per-point: use a comprehension lat, lon, alt = np.array( [(c.latitude_deg, c.longitude_deg, c.altitude) for c in map(sk.itrfcoord, vecs_itrf)]).T print(f'latitude range: [{lat.min():8.3f}, {lat.max():8.3f}] deg') print(f'altitude range: [{alt.min()/1e3:8.2f}, {alt.max()/1e3:8.2f}] km')

Boston in station ENU frame (m): [   132.27025425 -36211.61311894   -103.01860009]

latitude range:  [  41.562,   43.364] deg
altitude range:  [    0.02,   101.24] km