0

I would like to load a raster from a GeoTIFF file, and then load a time series of data from the Open Data Cube into the same CRS and resolution so I can use the raster as a mask. I can use the raster's .geobox to do this nicely when my GeoTIFF is in projected coordinates:

import xarray as xr
import datacube
dc = datacube.Datacube()

# Load raster and remove the redundant "band" dimension
raster = xr.open_rasterio("raster.tif").squeeze("band")

# Load data from datacube into raster geobox
ds = dc.load(
    product="s2b_ard_granule",
    measurements=["nbart_red"],
    time=("2018-01", "2018-02"),
    like=raster.geobox
)

# Use raster data as a mask
ds.where(raster > 500)

However, if my raster is in geographic coordinates (degrees), datacube will load with dimensions named longitude and latitude, not x and y. This means the dimensions no longer match the dimensions in raster, which prevents me from using it as a mask in the final step.

I've tried simply re-naming the dimensions in ds using ds.rename({'longitude': 'x', 'latitude': 'y'}). This fixes the immediate problem and allows me to use the raster data as a mask, but ds still has a .geobox with dimensions named longitude and latitude after the rename which now don't match the xr.DataArray's dimensions. This causes further issues downstream where I want to use the .geobox for other steps.

Is there:

  1. Any way to load data from datacube to match a raster .geobox, but not switch to the longitude and latitude dimension names even if the raster has geographic coords?
  2. If not, is there any way to completely rename the coords so that they are renamed in both the xr.DataArray and the .geobox?

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.