Does 'xr_reproject' handle rpcs ?

[remi-braun]

Hello,

After some testing, I found that xr_reproject doesn't seem to handle well rpcs.
I cannot give you the example as the data is proprietary but here are some comparisons:

rioxarray

        kwargs.update(
            {
                "RPC_DEM": dem_path,
                "RPC_DEM_MISSING_VALUE": 0,
                "OSR_USE_ETMERC": "YES",
                "BIGTIFF": "IF_NEEDED",
            }
        )

        # Reproject
        out_xda = src_xda.rio.reproject(
            dst_crs=crs,
            resolution=pixel_size,
            resampling=Resampling.bilinear,
            nodata=nodata,
            num_threads=utils.get_max_cores(),
            rpcs=rpcs,
            dtype=src_xda.dtype,
            **kwargs
        )

src_xda:

<xarray.DataArray 'DIM_SPOT7_MS_201602150257025_SEN_1671661101' (band: 4,
                                                                 y: 5557,
                                                                 x: 5239)> Size: 466MB
dask.array<where, shape=(4, 5557, 5239), dtype=float32, chunksize=(1, 5557, 5239), chunktype=numpy.ndarray>
Coordinates:
  * x            (x) float64 42kB 0.5 1.5 2.5 ... 5.236e+03 5.238e+03 5.238e+03
  * y            (y) float64 44kB 0.5 1.5 2.5 ... 5.554e+03 5.556e+03 5.556e+03
    spatial_ref  int64 8B 0
  * band         (band) int64 32B 1 2 3 4

out_xda:

<xarray.DataArray 'DIM_SPOT7_MS_201602150257025_SEN_1671661101' (band: 4,
                                                                 y: 8634,
                                                                 x: 6937)> Size: 958MB
Coordinates:
  * x            (x) float64 55kB 6.801e+05 6.801e+05 ... 7.217e+05 7.217e+05
  * y            (y) float64 69kB 9.342e+06 9.342e+06 ... 9.29e+06 9.29e+06
  * band         (band) int64 32B 1 2 3 4
    spatial_ref  int64 8B 0

xr_reproject

        kwargs.update(
            {
                "RPC_DEM": dem_path,
                "RPC_DEM_MISSING_VALUE": 0,
                "OSR_USE_ETMERC": "YES",
                "BIGTIFF": "IF_NEEDED",
            }
        )

        out_xda = src_xda.odc.reproject(
            how=crs,
            resolution=pixel_size,
            resampling=Resampling.bilinear,
            dst_nodata=nodata,
            num_threads=utils.get_max_cores(),
            rpcs=rpcs,
            dtype=src_xda.dtype,
            **kwargs
        )

out_xda:

<xarray.DataArray 'reproject-25d99544f75c4cdb99c0b8dd990b39ee' (band: 4,
                                                                y: 6661998,
                                                                x: 4525385)> Size: 482TB
dask.array<reproject, shape=(4, 6661998, 4525385), dtype=float32, chunksize=(1, 5557, 5239), chunktype=numpy.ndarray>
Coordinates:
  * band         (band) int64 32B 1 2 3 4
  * y            (y) float64 53MB 3e+07 3e+07 3e+07 ... -9.976e+06 -9.976e+06
  * x            (x) float64 36MB -1.255e+07 -1.255e+07 ... 1.46e+07 1.46e+07
    spatial_ref  int32 4B 32748

As you can see, the output shape is completely wrong. My feeling is that the RPCs are propagated to ever chunk.

[Kirill888]

@remi-braun no, RPC based GeoBoxes are not supported, GCPs are, but only as input not as output. If interested please create a new issue for that feature request.

[remi-braun]

@Kirill888 thanks for the answer! I'll create it then :)

Do you have any documentation on how to use GCPs with xr_reproject?