Unable to determine bounds from coordinates after a projection+combo

[Scolymus]

I want to calculate an index map that is the substraction/addition of two different bands (I'll use Sentinel-2 2A images). I found in stackoverflow that this can be archieved by doing:

import rasterio as rio
import rioxarray as riox

CH1 = riox.open_rasterio(filename_CH1).sel(band=1).astype('float32')
CH2 = riox.open_rasterio(filename_CH2).sel(band=1).astype('float32')

combo = ((CH1-CH2)/(CH1+CH2)).where((CH1+CH2)!=0., 0)
combo.attrs = CH1.attrs
combo = combo.rio.write_crs(CH1.rio.crs)
combo.rio.to_raster(filename_out)

CH1.close()
CH2.close()
combo.close()

This code works perfectly. However, I need to construct an index map where one band is 10m resolution and the other is 20m. Since both images refer to the same area, I introduced:

import rasterio as rio
from rasterio.enums import Resampling
import rioxarray as riox

CH1 = riox.open_rasterio(filename_CH1).sel(band=1).astype('float32')
CH2 = riox.open_rasterio(filename_CH2).sel(band=1).astype('float32')

# Resize if different image size
# Notice: This will work because we refer always to the same field!
if CH1.shape > CH2.shape:
    CH1 = CH1.rio.reproject(
            CH1.rio.crs,
            shape=(CH2.shape[0], CH2.shape[1]),
            resampling=Resampling.bilinear,
    )
elif CH1.shape < CH2.shape:
    CH2 = CH2.rio.reproject(
            CH2.rio.crs,
            shape=(CH1.shape[0], CH1.shape[1]),
            resampling=Resampling.bilinear,
    )

combo = ((CH1-CH2)/(CH1+CH2)).where((CH1+CH2)!=0., 0)
combo.attrs = CH1.attrs
combo = combo.rio.write_crs(CH1.rio.crs)
combo.rio.to_raster(filename_out)

CH1.close()
CH2.close()
combo.close()

However, this code is producing me a warning and an error.

• First, if I try to save the reprojected image I have this warning: /opt/homebrew/lib/python3.9/site-packages/rioxarray/raster_writer.py:110: UserWarning: The nodata value (3.402823466e+38) has been automatically changed to (3.4028234663852886e+38) to match the dtype of the data.
  warnings.warn( I guess this is normal because of the bilinear calculation, and I shouldn't worry about.
• Second, I cannot save combo. I get the following error, which I couldn't fix it:
  Traceback (most recent call last):
  File "/Users/scolymus/Documents/code/tests.py", line 87, in
  combo.rio.to_raster("/Users/scolymus/image_test.tif")
  File "/opt/homebrew/lib/python3.9/site-packages/rioxarray/raster_array.py", line 1112, in to_raster
  transform=self.transform(recalc=recalc_transform),
  File "/opt/homebrew/lib/python3.9/site-packages/rioxarray/rioxarray.py", line 656, in transform
  src_left, _, _, src_top = self._unordered_bounds(recalc=recalc)
  File "/opt/homebrew/lib/python3.9/site-packages/rioxarray/rioxarray.py", line 1028, in _unordered_bounds
  resolution_x, resolution_y = self.resolution(recalc=recalc)
  File "/opt/homebrew/lib/python3.9/site-packages/rioxarray/rioxarray.py", line 995, in resolution
  left, bottom, right, top = self._internal_bounds()
  File "/opt/homebrew/lib/python3.9/site-packages/rioxarray/rioxarray.py", line 963, in _internal_bounds
  raise NoDataInBounds(
  rioxarray.exceptions.NoDataInBounds: Unable to determine bounds from coordinates.

How can I fix it? Thanks!
PS. Btw, is it correct to substitue 0 for NaN in the where operation or should I assign something different/play with masks?

[snowman2]

You may be interested in rio.reproject_match: https://corteva.github.io/rioxarray/stable/examples/reproject_match.html

Also, you may want to write your metadata nodata value as NaN if that is what you expect it to be: https://corteva.github.io/rioxarray/stable/getting_started/nodata_management.html

xds.rio.write_nodata(numpy.nan, inplace=True)

[Scolymus]

ok! nice! It's working!