# How does rioxarray reproject match handle nan values?

I have a couple GeoTiffs where the "no data" value is set to 32767. I want to set those values to NaN and then use nearest neighbor interpolation to resample the data and match the projection and resolution of the other GeoTiff. To do this, I use `rio.reproject_match`. The problem is that when I reproject the layer, the sum of all items returns `0` which I don't think is correct. Can someone find what I'm doing wrong here?

``````import numpy as np
import rioxarray
import xarray
from rasterio.enums import Resampling

# Source Raster
xds = rioxarray.open_rasterio(LS)
# Set no data value 32767 to NaN
xds = xds.where(xds != 32767, np.nan)
xds.rio.write_nodata(np.nan, inplace=True)

print(xds)
# <xarray.DataArray (band: 1, y: 1006, x: 950)>
# array([[[nan, nan, nan, ..., nan, nan, nan],
#         [nan, nan, nan, ..., nan, nan, nan],
#         [nan, nan, nan, ..., nan, nan, nan],
#         ...,
#         [nan, nan, nan, ..., nan, nan, nan],
#         [nan, nan, nan, ..., nan, nan, nan],
#         [nan, nan, nan, ..., nan, nan, nan]]])
# Coordinates:
#   * band         (band) int64 1
#   * x            (x) float64 -4.56e+03 -4.53e+03 ... 2.388e+04 2.391e+04
#   * y            (y) float64 2.029e+06 2.029e+06 ... 1.999e+06 1.998e+06
#     spatial_ref  int64 0
# Attributes: (12/13)
#     AREA_OR_POINT:           Area
#     DataType:                Generic
#     RepresentationType:      ATHEMATIC
#     STATISTICS_COVARIANCES:  20807216.58662666
#     STATISTICS_MAXIMUM:      11111
#     STATISTICS_MEAN:         2561.3518011412
#     ...                      ...
#     STATISTICS_SKIPFACTORX:  1
#     STATISTICS_SKIPFACTORY:  1
#     STATISTICS_STDDEV:       4561.4928024307
#     _FillValue:              nan
#     scale_factor:            1.0
#     add_offset:              0.0

print_raster(xds)
# Source Raster:
# ----------------
# shape: (1006, 950)
# resolution: (30.0, -30.0)
# bounds: (-4575.0, 1998465.0, 23925.0, 2028645.0)
# sum: 2284131573.0
# CRS: EPSG:9822

# Target Raster
xds_match = rioxarray.open_rasterio(precip)
xds_match = xds_match.where(xds_match != 32767, np.nan)
xds_match.rio.write_nodata(np.nan, inplace=True)

print_raster(xds_match)
# Raster to Match:
# ----------------
# shape: (2898, 2772)
# resolution: (10.0, -10.0)
# bounds: (747720.0, 4543380.0, 775440.0, 4572360.0)
# sum: 7517694.093363524
# CRS: EPSG:32614

# Resample and reproject the source raster to match the target raster
xds_repr_match = xds.rio.reproject_match(xds_match, resampling=Resampling.nearest, nodata=np.nan)

print_raster(xds_repr_match)
# Reprojected Raster:
# -------------------
# shape: (2898, 2772)
# resolution: (10.0, -10.0)
# bounds: (747720.0, 4543380.0, 775440.0, 4572360.0)
# sum: 0.0
# CRS: EPSG:32614

``````

I don't expect the sum to match the target raster sum, but I do expect it to be more than 0.

I've read this thread, but they're talking about `reproject` instead of `reproject_match` mostly. How to deal with Nan value when using rioxarray rio.reproject()?. Nonetheless, I built a little test case using `reproject` with some dummy data just to see if I understood how `sum` works in conjunction with NaN values and it does it properly (doesn't return 0):

``````import xarray as xr
import rioxarray as rxr
import numpy as np
from rasterio.enums import Resampling

data = [[1, 2, 3, 4], [5, 32767, 7, 8], [9, 10, 11, 12], [13, 14, 15, 16]]
x = np.arange(100, 116, 4)
y = np.arange(36, 20, -4)

foo = xr.DataArray(data, coords=[y, x], dims=['y', 'x'])
foo = foo.rio.write_crs('EPSG:4326')
foo = foo.where(foo != 32767, np.nan)
foo.rio.write_nodata(np.nan, inplace=True)
foo_resampled = foo.rio.reproject(
foo.rio.crs,
resolution=8,
resampling=Resampling.nearest,
nodata=np.nan
)

print(foo.sum().item())
# 130.0

print(foo_resampled.sum().item())
# 38.0
``````

Why doesn't `reproject_match` work in the way I expect with `np.nan` values and how can I fix it?

## 1 Answer

This turned out to a be a problem with the initial projection of the source raster (EPSG: 9822 is supposed to be Albers Equal Area Conic but it's not a CRS that is valid in QGIS and looking up the info on epsg.org, it seems to point to something related to France...). I defined the original projection of the source raster to the correct projection of EPSG 6350 and after running the script again, it appears to have fixed the issue.