Creating

There are several ways to create a RasterIndex, depending on the source of your data and the metadata conventions it uses.

import numpy as np
import xarray as xr
import rasterix

Using assign_index

The easiest way to create a RasterIndex is using assign_index(). This function automatically detects the affine transform from various metadata conventions, including:

See Heuristics for the full priority order and detection logic.

import pyproj

# Example: dataset with GDAL-style GeoTransform metadata
ds = xr.Dataset(
    {"temperature": (("y", "x"), np.random.rand(100, 100))},
    coords={
        "spatial_ref": (
            (),
            0,
            pyproj.CRS.from_epsg(32610).to_cf() | {"GeoTransform": "400000.0 10.0 0.0 5000000.0 0.0 -10.0"},
        )
    },
    attrs={"grid_mapping": "spatial_ref"},
)

# assign_index auto-detects the convention and creates the RasterIndex
ds = rasterix.assign_index(ds)
ds

<xarray.Dataset> Size: 82kB
Dimensions:      (y: 100, x: 100)
Coordinates:
  * y            (y) float64 800B 5e+06 5e+06 5e+06 ... 4.999e+06 4.999e+06
  * x            (x) float64 800B 4e+05 4e+05 4e+05 ... 4.01e+05 4.01e+05
    spatial_ref  int64 8B 0
Data variables:
    temperature  (y, x) float64 80kB 0.2836 0.3261 0.7299 ... 0.9677 0.01361
Indexes:
  ┌ x        RasterIndex (crs=None)
  └ y
Attributes:
    grid_mapping:  spatial_ref

Zarr spatial convention

For data following the Zarr Spatial Convention, spatial: attributes are detected on the variable (array-level) or the Dataset (group-level), with array-level attributes taking precedence. The convention must be registered in the zarr_conventions attribute. When present, spatial:dimensions is used to auto-detect the spatial dimension names, so non-standard names work without passing x_dim/y_dim; the CRS is detected from the companion proj: convention if present:

ds = xr.Dataset(
    {"temperature": (("Y", "X"), np.random.rand(100, 100))},
    attrs={
        "zarr_conventions": [{"name": "spatial:"}, {"name": "proj:"}],
        "spatial:dimensions": ["Y", "X"],
        "spatial:transform": [10.0, 0.0, 400000.0, 0.0, -10.0, 5000000.0],
        "proj:code": "EPSG:32610",
    },
)

ds = rasterix.assign_index(ds)
ds

<xarray.Dataset> Size: 82kB
Dimensions:      (Y: 100, X: 100)
Coordinates:
  * Y            (Y) float64 800B 5e+06 5e+06 5e+06 ... 4.999e+06 4.999e+06
  * X            (X) float64 800B 4e+05 4e+05 4e+05 ... 4.01e+05 4.01e+05
Data variables:
    temperature  (Y, X) float64 80kB 0.2417 0.1806 0.1136 ... 0.5166 0.5543
Indexes:
  ┌ X        RasterIndex (crs=EPSG:32610)
  └ Y
Attributes:
    zarr_conventions:    [{'name': 'spatial:'}, {'name': 'proj:'}]
    spatial:dimensions:  ['Y', 'X']
    proj:code:           EPSG:32610

Direct construction with class methods

For more control, you can create a RasterIndex directly using class methods. This is useful when you have the transform parameters available directly, or when working with data that doesn’t have embedded metadata.

From an Affine transform

Use RasterIndex.from_transform() to create from an Affine transform directly (corresponds to GDAL GeoTransform convention):

from affine import Affine

transform = Affine(10.0, 0.0, 400000.0, 0.0, -10.0, 5000000.0)
index = rasterix.RasterIndex.from_transform(
    transform,
    width=100,
    height=100,
    crs="EPSG:32610",
)
index

RasterIndex(crs=EPSG:32610)
    AxisAffineTransformIndex(AxisAffineTransform(a=10, b=0, c=4e+05, d=0, e=-10, f=5e+06, axis=X, dim='x'))
    AxisAffineTransformIndex(AxisAffineTransform(a=10, b=0, c=4e+05, d=0, e=-10, f=5e+06, axis=Y, dim='y'))

# Assign to data
ds_manual = xr.Dataset(
    {"data": (("y", "x"), np.random.rand(100, 100))},
    coords=xr.Coordinates.from_xindex(index),
)
ds_manual

<xarray.Dataset> Size: 82kB
Dimensions:  (y: 100, x: 100)
Coordinates:
  * y        (y) float64 800B 5e+06 5e+06 5e+06 ... 4.999e+06 4.999e+06
  * x        (x) float64 800B 4e+05 4e+05 4e+05 ... 4.01e+05 4.01e+05 4.01e+05
Data variables:
    data     (y, x) float64 80kB 0.8056 0.3171 0.5731 ... 0.9262 0.1775 0.8366
Indexes:
  ┌ x        RasterIndex (crs=EPSG:32610)
  └ y

From a GeoTransform

Use RasterIndex.from_geotransform() to create from a GDAL GeoTransform, either as a sequence or a space-separated string:

# From a tuple
geotransform = (400000.0, 10.0, 0.0, 5000000.0, 0.0, -10.0)
index = rasterix.RasterIndex.from_geotransform(
    geotransform,
    width=100,
    height=100,
    crs="EPSG:32610",
)
index

RasterIndex(crs=EPSG:32610)
    AxisAffineTransformIndex(AxisAffineTransform(a=10, b=0, c=4e+05, d=0, e=-10, f=5e+06, axis=X, dim='x'))
    AxisAffineTransformIndex(AxisAffineTransform(a=10, b=0, c=4e+05, d=0, e=-10, f=5e+06, axis=Y, dim='y'))

# From a string (as commonly stored in netCDF attributes)
geotransform = "400000.0 10.0 0.0 5000000.0 0.0 -10.0"
index = rasterix.RasterIndex.from_geotransform(
    geotransform,
    width=100,
    height=100,
    crs="EPSG:32610",
)
index

RasterIndex(crs=EPSG:32610)
    AxisAffineTransformIndex(AxisAffineTransform(a=10, b=0, c=4e+05, d=0, e=-10, f=5e+06, axis=X, dim='x'))
    AxisAffineTransformIndex(AxisAffineTransform(a=10, b=0, c=4e+05, d=0, e=-10, f=5e+06, axis=Y, dim='y'))

From GeoTIFF tiepoint and scale

Use RasterIndex.from_tiepoint_and_scale() to create from GeoTIFF-style tiepoint and pixel scale attributes:

index = rasterix.RasterIndex.from_tiepoint_and_scale(
    tiepoint=[0.0, 0.0, 0.0, 323400.0, 4265400.0, 0.0],
    scale=[30.0, 30.0, 0.0],
    width=100,
    height=100,
    crs="EPSG:32610",
)
index

RasterIndex(crs=EPSG:32610)
    AxisAffineTransformIndex(AxisAffineTransform(a=30, b=0, c=3.234e+05, d=0, e=30, f=4.265e+06, axis=X, dim='x'))
    AxisAffineTransformIndex(AxisAffineTransform(a=30, b=0, c=3.234e+05, d=0, e=30, f=4.265e+06, axis=Y, dim='y'))

From STAC projection metadata

Use RasterIndex.from_stac_proj_metadata() to create from STAC projection extension metadata:

metadata = {"proj:transform": [30.0, 0.0, 323400.0, 0.0, -30.0, 4268400.0]}
index = rasterix.RasterIndex.from_stac_proj_metadata(
    metadata,
    width=100,
    height=100,
    crs="EPSG:32610",
)
index

RasterIndex(crs=EPSG:32610)
    AxisAffineTransformIndex(AxisAffineTransform(a=30, b=0, c=3.234e+05, d=0, e=-30, f=4.268e+06, axis=X, dim='x'))
    AxisAffineTransformIndex(AxisAffineTransform(a=30, b=0, c=3.234e+05, d=0, e=-30, f=4.268e+06, axis=Y, dim='y'))

Accessing transform information

Once created, you can access the affine transform using methods on RasterIndex:

ds = rasterix.assign_index(
    xr.Dataset(
        {"data": (("y", "x"), np.random.rand(100, 100))},
        coords={
            "spatial_ref": ((), 0, {"GeoTransform": "400000.0 10.0 0.0 5000000.0 0.0 -10.0"}),
        },
    )
)

# Top-left corner transform (GDAL convention) via transform()
ds.xindexes["x"].transform()

Affine(10.0, 0.0, 400000.0,
       0.0, -10.0, 5000000.0)

# Pixel center transform via center_transform()
ds.xindexes["x"].center_transform()

Affine(10.0, 0.0, 400005.0,
       0.0, -10.0, 4999995.0)

# Bounding box via the bbox property
ds.xindexes["x"].bbox

BoundingBox(left=400000.0, bottom=4999000.0, right=401000.0, top=5000000.0)

# As GeoTransform string (for saving) via as_geotransform()
ds.xindexes["x"].as_geotransform()

'400000.0 10.0 0.0 5000000.0 0.0 -10.0'