geo
geo
Geospatial utilities for planetary image data.
Coordinate transforms between pixel, map-projected, and geographic coordinate systems, built on rasterio and pyproj. No GDAL required.
Works with rioxarray DataArrays (the standard way to open planetary images in planetarypy) or standalone with affine transforms and CRS. Supports IAU planetary CRS codes for Mars, Moon, and all solar system bodies.
Functions for direct use: pixel_to_xy, xy_to_pixel, pixel_to_lonlat, lonlat_to_pixel, xy_to_lonlat, lonlat_to_xy, is_within, image_azimuth, image_azimuth_cw_from_right, pixel_resolution
Projected-raster geometry (GDAL-native, ISIS-free, format-agnostic): is_projected, raster_footprint, footprints_to_gdf, overlaps
Anti-meridian (±180° longitude) handling: split_at_antimeridian, normalise_lon_bounds
Point class: Combines lon/lat, pixel, and projected coordinates in one object with CRS awareness. Handles transforms, bounds checking, azimuths, and Shapely interop — filling a gap that Shapely (no CRS/pixels), pyproj (no point object), and rasterio (no lon/lat) each leave open.
Examples
>>> from planetarypy.geo import Point, pixel_to_lonlat
>>> p = Point(lon=137.85, lat=-5.08, crs="IAU_2015:49900")
>>> p.is_within(da)
True
>>> p.to_shapely()
<POINT (137.85 -5.08)>Classes
| Name | Description |
|---|---|
| Point | A geographic point on a planetary body. |
Point
geo.Point(
lon=None,
lat=None,
x=None,
y=None,
sample=None,
line=None,
crs=None,
source=None,
)A geographic point on a planetary body.
Can be created from lon/lat, projected (x, y), or pixel coordinates. Converts lazily between coordinate systems as needed.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| lon | float | Geographic coordinates in degrees. | None |
| lat | float | Geographic coordinates in degrees. | None |
| x | float | Map-projected coordinates. | None |
| y | float | Map-projected coordinates. | None |
| sample | float | Pixel coordinates. | None |
| line | float | Pixel coordinates. | None |
| crs | str or pyproj.CRS | CRS for projected/geographic coordinates (e.g. “IAU_2015:49900”). | None |
| source | xarray.DataArray or rasterio.DatasetReader | Image source — provides CRS and transform for pixel conversions. | None |
Examples
>>> p = Point(lon=137.85, lat=-5.08, crs="IAU_2015:49900")
>>> p.to_xy("IAU_2015:49910")
(8167439.5, -301670.8)>>> p = Point(sample=250, line=250, source=da)
>>> p.lon, p.lat
(137.4, -4.6)Methods
| Name | Description |
|---|---|
| azimuth_to | Azimuth angle from this point to another, clockwise from north. |
| is_within | Check if this point falls within an image. |
| to_pixel | Convert to pixel coordinates in an image. |
| to_shapely | Convert to a Shapely Point (lon, lat). |
| to_xy | Convert to map-projected coordinates. |
azimuth_to
geo.Point.azimuth_to(other)Azimuth angle from this point to another, clockwise from north.
Both points must have pixel coordinates (via source).
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| other | Point | required |
Returns
| Name | Type | Description |
|---|---|---|
| float | Degrees, clockwise from north. |
is_within
geo.Point.is_within(source=None)Check if this point falls within an image.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| source | xarray.DataArray or rasterio.DatasetReader | Image to check against. Uses the source from init if not given. | None |
Returns
| Name | Type | Description |
|---|---|---|
| bool |
to_pixel
geo.Point.to_pixel(source=None)Convert to pixel coordinates in an image.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| source | xarray.DataArray or rasterio.DatasetReader | Image to project into. Uses the source from init if not given. | None |
Returns
| Name | Type | Description |
|---|---|---|
| sample, line : float |
to_shapely
geo.Point.to_shapely()Convert to a Shapely Point (lon, lat).
Returns
| Name | Type | Description |
|---|---|---|
| shapely.geometry.Point |
to_xy
geo.Point.to_xy(target_crs=None)Convert to map-projected coordinates.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| target_crs | str or CRS | Target projected CRS. If None, uses the point’s own CRS. | None |
Returns
| Name | Type | Description |
|---|---|---|
| x, y : float |
Functions
| Name | Description |
|---|---|
| footprints_to_gdf | Assemble valid-data footprints of many rasters into a GeoDataFrame. |
| image_azimuth | Calculate azimuth angle between two image points. |
| image_azimuth_cw_from_right | Calculate azimuth angle clockwise from right (3 o’clock). |
| is_projected | Return True if the raster’s CRS is projected (not geographic). |
| is_within | Check if geographic coordinates fall within an image. |
| lonlat_to_pixel | Convert geographic (lon, lat) to pixel coordinates. |
| lonlat_to_xy | Convert geographic (lon, lat) to map-projected coordinates. |
| normalise_lon_bounds | Normalise a longitude bbox to [-180, 180], picking the short arc. |
| overlaps | Pairwise overlaps (positive area) between footprints in gdf. |
| pixel_resolution | Return the pixel resolution in map units. |
| pixel_to_lonlat | Convert pixel coordinates to geographic (lon, lat). |
| pixel_to_xy | Convert pixel coordinates to map-projected (x, y) coordinates. |
| raster_footprint | Valid-data footprint of a raster as a shapely (Multi)Polygon. |
| split_at_antimeridian | Build (possibly split) shapely polygon(s) from lon/lat corners. |
| xy_to_lonlat | Convert map-projected coordinates to geographic (lon, lat). |
| xy_to_pixel | Convert map-projected (x, y) coordinates to pixel coordinates. |
footprints_to_gdf
geo.footprints_to_gdf(sources, *, id_fn=None, simplify=None)Assemble valid-data footprints of many rasters into a GeoDataFrame.
One row per source with an id and a geometry column. All sources must share a single CRS (raises otherwise); the GeoDataFrame carries it.
id_fn maps a source to its identifier. The default uses the file name with extension (e.g. ESP_075205_0930_RED5_0.cub) — generic, instrument-free, and lossless: keeping the extension distinguishes the same product in different formats (e.g. a .cub and its .tif conversion, whose footprints may differ). Note that neither the name nor its stem is the PDS product id (filenames carry processing suffixes); pass your own id_fn for a domain key, or id_fn=str for full paths.
Raises if the resulting ids aren’t unique — duplicate ids would silently corrupt downstream id-keyed operations like :func:overlaps.
Requires geopandas (the [isis] extra).
image_azimuth
geo.image_azimuth(sample1, line1, sample2, line2)Calculate azimuth angle between two image points.
Returns the clockwise angle from image-north (up) in degrees. This matches the convention for images with origin at upper-left (standard for planetary data).
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| sample1 | float | Origin point (column, row). | required |
| line1 | float | Origin point (column, row). | required |
| sample2 | float | Target point (column, row). | required |
| line2 | float | Target point (column, row). | required |
Returns
| Name | Type | Description |
|---|---|---|
| azimuth | float | Angle in degrees [0, 360), clockwise from up. |
image_azimuth_cw_from_right
geo.image_azimuth_cw_from_right(sample1, line1, sample2, line2)Calculate azimuth angle clockwise from right (3 o’clock).
This is the convention used by HiRISE, where angles are measured clockwise from the +sample direction (right/east in image space). The y-axis is flipped (line increases downward), so clockwise in image space matches the trigonometric direction.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| sample1 | float | Origin point (column, row). | required |
| line1 | float | Origin point (column, row). | required |
| sample2 | float | Target point (column, row). | required |
| line2 | float | Target point (column, row). | required |
Returns
| Name | Type | Description |
|---|---|---|
| azimuth | float | Angle in degrees [0, 360), clockwise from right. |
is_projected
geo.is_projected(source)Return True if the raster’s CRS is projected (not geographic).
Accepts a file path, an open rasterio dataset, or a rioxarray DataArray. Format-agnostic (ISIS .cub, GeoTIFF, …).
is_within
geo.is_within(source, lon, lat)Check if geographic coordinates fall within an image.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| source | xarray.DataArray or rasterio.DatasetReader | Must have a CRS and transform. | required |
| lon | float or array - like | Geographic coordinates in degrees. | required |
| lat | float or array - like | Geographic coordinates in degrees. | required |
Returns
| Name | Type | Description |
|---|---|---|
| bool or np.ndarray of bool | True if the point(s) are within the image bounds. |
lonlat_to_pixel
geo.lonlat_to_pixel(source, lon, lat)Convert geographic (lon, lat) to pixel coordinates.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| source | xarray.DataArray or rasterio.DatasetReader | Must have a CRS and transform. | required |
| lon | float or array - like | Geographic coordinates in degrees. | required |
| lat | float or array - like | Geographic coordinates in degrees. | required |
Returns
| Name | Type | Description |
|---|---|---|
| sample, line : float or np.ndarray | Pixel column(s) and row(s). |
lonlat_to_xy
geo.lonlat_to_xy(crs, lon, lat)Convert geographic (lon, lat) to map-projected coordinates.
Works for both Earth (EPSG) and planetary CRS (IAU codes, proj4 strings).
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| crs | rasterio.crs.CRS or pyproj.CRS | The target coordinate reference system. | required |
| lon | float or array - like | Geographic coordinates in degrees. | required |
| lat | float or array - like | Geographic coordinates in degrees. | required |
Returns
| Name | Type | Description |
|---|---|---|
| x, y : float or np.ndarray | Projected coordinates. |
normalise_lon_bounds
geo.normalise_lon_bounds(lon_min, lon_max)Normalise a longitude bbox to [-180, 180], picking the short arc.
Distinguishes two wrap geometries that a raw [0, 360] bbox conflates:
- Antimeridian wrap (e.g. raw
(170, 190)): raw width < 180°. The[-180, 180]form haslon_min > lon_max; the interior is(lon >= lon_min) OR (lon <= lon_max). - Prime-meridian wrap (e.g. raw
(8, 359)for a narrow feature near 0°): raw width > 180°, so the vertices encode the long way; the real (short) interior is the single interval(raw_max - 360, raw_min).
Returns
| Name | Type | Description |
|---|---|---|
| (lon_min, lon_max, crosses) : tuple[float, float, bool] | Bounds in [-180, 180] and whether the interior crosses ±180° (i.e. callers must use OR rather than AND on the longitude test). |
overlaps
geo.overlaps(gdf)Pairwise overlaps (positive area) between footprints in gdf.
Takes a GeoDataFrame with id + geometry (as built by :func:footprints_to_gdf, whose id values must be unique) and returns one with id_1, id_2, geometry (the intersection) and area (CRS units) per unordered intersecting pair. Edge-only touches (zero area) are dropped.
Requires geopandas (the [isis] extra).
pixel_resolution
geo.pixel_resolution(transform)Return the pixel resolution in map units.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| transform | affine.Affine | The affine geotransform. | required |
Returns
| Name | Type | Description |
|---|---|---|
| res_x, res_y : float | Pixel size in x and y (map units, typically meters). res_y is returned as a positive value. |
pixel_to_lonlat
geo.pixel_to_lonlat(source, sample, line)Convert pixel coordinates to geographic (lon, lat).
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| source | xarray.DataArray or rasterio.DatasetReader | Must have a CRS and transform (e.g. opened via rioxarray). | required |
| sample | float or array - like | Pixel column(s). | required |
| line | float or array - like | Pixel row(s). | required |
Returns
| Name | Type | Description |
|---|---|---|
| lon, lat : float or np.ndarray | Geographic coordinates in degrees. |
pixel_to_xy
geo.pixel_to_xy(transform, sample, line)Convert pixel coordinates to map-projected (x, y) coordinates.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| transform | affine.Affine or rasterio.transform.Affine | The affine geotransform (from rasterio dataset or rioxarray). | required |
| sample | float or array - like | Pixel column(s) (0-based). | required |
| line | float or array - like | Pixel row(s) (0-based). | required |
Returns
| Name | Type | Description |
|---|---|---|
| x, y : float or np.ndarray | Map-projected coordinates. |
raster_footprint
geo.raster_footprint(source, *, simplify=None)Valid-data footprint of a raster as a shapely (Multi)Polygon.
Built from the dataset mask, so nodata borders — common on map-projected products — are excluded: the result is the actual data outline, not the bounding box. The polygon is in the raster’s own CRS. Returns None for an all-nodata raster.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| source | path | rasterio dataset | required | |
| simplify | float | Douglas-Peucker tolerance (CRS units) to thin the outline. | None |
split_at_antimeridian
geo.split_at_antimeridian(corners)Build (possibly split) shapely polygon(s) from lon/lat corners.
Three cases:
- No crossing (lon-span < 180°): one polygon, constructed directly.
- Antimeridian crossing (span >= 180°, but < 180° after shifting negative-lon corners by +360°): shift, split at lon=180, shift the over-180° piece back — returns the two hemisphere halves.
- Pole-containing (still spans >= 180° after the shift): delegate to :func:
antimeridian.fix_polygon, forcing the pole by the mean-latitude sign — returns one polar-cap polygon.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| corners | list[tuple[float, float]] | (lon, lat) vertices in degrees. |
required |
Returns
| Name | Type | Description |
|---|---|---|
| list[shapely.Polygon] | One polygon (no-crossing or polar cap), two (antimeridian crossing), or empty for degenerate input (< 3 corners / invalid geometry). |
xy_to_lonlat
geo.xy_to_lonlat(crs, x, y)Convert map-projected coordinates to geographic (lon, lat).
Works for both Earth (EPSG) and planetary CRS (IAU codes, proj4 strings).
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| crs | rasterio.crs.CRS or pyproj.CRS | The coordinate reference system of the projected data. | required |
| x | float or array - like | Projected coordinates. | required |
| y | float or array - like | Projected coordinates. | required |
Returns
| Name | Type | Description |
|---|---|---|
| lon, lat : float or np.ndarray | Geographic coordinates in degrees. |
xy_to_pixel
geo.xy_to_pixel(transform, x, y)Convert map-projected (x, y) coordinates to pixel coordinates.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| transform | affine.Affine | The affine geotransform. | required |
| x | float or array - like | Map-projected coordinates. | required |
| y | float or array - like | Map-projected coordinates. | required |
Returns
| Name | Type | Description |
|---|---|---|
| sample, line : float or np.ndarray | Pixel column(s) and row(s) (0-based, fractional). |