Rioxarray reading COG (from Earth Engine) with bottom > top

Cloud
1

I noticed this when I was running a pad_xy on my dataset expanding its dimensions to fit a bounding box which I’ll be filling with other data (I expect a dataset filling the coordinates, with NA), this looks like so:

bounds = [30.233265, 50.071171, 31.274053, 50.742212]
ds = rxr.open_rasterio(uri, masked=True, mask_and_scale=False, band_as_variable=True, parse_coordinates=True)
ds = ds.rio.pad_box(minx=bounds[0], miny=bounds[1], maxx=bounds[2], maxy=bounds[3], constant_values=None)

Resulting in:

  File ".../xarray_complete_zonal.py", line 139, in get_dataset
    ds_expanded_time = ds_expanded_time.rio.pad_box(minx=minx, miny=miny, maxx=maxx, maxy=maxy, constant_values=None)
  File ".../rioxarray/rioxarray/raster_dataset.py", line 239, in pad_box
    .rio.pad_box(
         ~~~~~~~^
        minx, miny, maxx, maxy, constant_values=constant_values
        ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
    )
    ^
  File ".../rioxarray/rioxarray/raster_array.py", line 743, in pad_box
    pd_array = self.pad_xy(
        minx=pad_minx,
    ...<3 lines>...
        constant_values=constant_values,
    )
  File "/Users/jww/epoch/rioxarray/rioxarray/raster_array.py", line 668, in pad_xy
    top = y_coord[0]
          ~~~~~~~^^^
IndexError: index 0 is out of bounds for axis 0 with size 0

This seemed rather strange, and lookin at the rioxarray code it’s:

if top - resolution_y < maxy:
    new_y_coord: numpy.ndarray = numpy.arange(bottom, maxy, -resolution_y)[::-1]
    y_before = len(new_y_coord) - len(y_coord)
    y_coord = new_y_coord
    top = y_coord[0]

The issue seems to be with the line new_y_coord: numpy.ndarray = numpy.arange(bottom, maxy, -resolution_y)[::-1] and by printing out my values for bottom, top and resolution_y etc.

I find that my bottom is > than my top, and that resolution_y is positive, such that
numpy.arange(bottom, maxy, -resolution_y)[::-1]
results in an empty array, due to resolution_y being expected to be a negative number

Using the rasterio get transform method, the gdal transform is
[30.601289816588324, 8.983152841195215e-05, 0.0, 50.06239213175364, 0.0, 8.983152841195215e-05] the 6th value being the y resolution (positive) and the imagery being in EPSG:4326

I understand that generally speaking resolution_y is expected to be a negative value from the GDAL docs " GT(5) n-s pixel resolution / pixel height (negative value for a north-up image)." so is this a problem I should fix with the projections of my datasets? Or is this an implementation error with rioxarray assuming negative values? I can also contact the EE team if this seems to be an issue

2

Resolved via:

with rasterio.open(uri) as src:
    with rasterio.vrt.WarpedVRT(src) as vrt:
        ds = rxr.open_rasterio(vrt, masked=True, mask_and_scale=False, band_as_variable=True, parse_coordinates=True)