diff --git a/rioxarray/raster_array.py b/rioxarray/raster_array.py index 73c795a3..a0c703de 100644 --- a/rioxarray/raster_array.py +++ b/rioxarray/raster_array.py @@ -9,6 +9,7 @@ - https://github.com/opendatacube/datacube-core/blob/1d345f08a10a13c316f81100936b0ad8b1a374eb/LICENSE # noqa: E501 """ + import copy import os from collections.abc import Hashable, Iterable, Mapping @@ -128,7 +129,13 @@ def _make_dst_affine( **kwargs, ): """Determine the affine of the new projected `xarray.DataArray`""" - src_bounds = () if "gcps" in kwargs else src_data_array.rio.bounds() + src_bounds = () + if ( + "gcps" not in kwargs + and "rpcs" not in kwargs + and "src_geoloc_array" not in kwargs + ): + src_bounds = src_data_array.rio.bounds() src_height, src_width = src_data_array.rio.shape dst_height, dst_width = dst_shape if dst_shape is not None else (None, None) # pylint: disable=isinstance-second-argument-not-valid-type @@ -453,8 +460,12 @@ def reproject( if gcps: kwargs.setdefault("gcps", gcps) - gcps_or_rpcs = "gcps" in kwargs or "rpcs" in kwargs - src_affine = None if gcps_or_rpcs else self.transform(recalc=True) + use_affine = ( + "gcps" not in kwargs + and "rpcs" not in kwargs + and "src_geoloc_array" not in kwargs + ) + src_affine = None if not use_affine else self.transform(recalc=True) if transform is None: dst_affine, dst_width, dst_height = _make_dst_affine( src_data_array=self._obj, @@ -474,7 +485,6 @@ def reproject( dst_data = self._create_dst_data(dst_height=dst_height, dst_width=dst_width) dst_nodata = self._get_dst_nodata(nodata) - rasterio.warp.reproject( source=self._obj.values, destination=dst_data, @@ -506,7 +516,7 @@ def reproject( dst_affine=dst_affine, dst_width=dst_width, dst_height=dst_height, - force_generate=gcps_or_rpcs, + force_generate=not use_affine, ), dims=tuple(dst_dims), attrs=new_attrs, diff --git a/test/conftest.py b/test/conftest.py index 46b8ad47..59f3d88d 100644 --- a/test/conftest.py +++ b/test/conftest.py @@ -19,6 +19,7 @@ GDAL_GE_36 = version.parse(rasterio.__gdal_version__) >= version.parse("3.6.0") GDAL_GE_361 = version.parse(rasterio.__gdal_version__) >= version.parse("3.6.1") GDAL_GE_364 = version.parse(rasterio.__gdal_version__) >= version.parse("3.6.4") +RASTERIO_GTE_1_4 = version.parse(rasterio.__version__) >= version.parse("1.4a1") # xarray.testing.assert_equal(input_xarray, compare_xarray) diff --git a/test/integration/test_integration_rioxarray.py b/test/integration/test_integration_rioxarray.py index 95bd8daa..113089a6 100644 --- a/test/integration/test_integration_rioxarray.py +++ b/test/integration/test_integration_rioxarray.py @@ -3184,3 +3184,29 @@ def test_bounds__ordered__dataarray(): def test_bounds__ordered__dataset(): xds = xarray.Dataset(None, coords={"x": range(5), "y": range(5)}) assert xds.rio.bounds() == (-0.5, -0.5, 4.5, 4.5) + + +@pytest.mark.skipif(not RASTERIO_GE_14, reason="Requires rasterio 1.4+") +def test_non_rectilinear__reproject(open_rasterio): + test_file = os.path.join(TEST_INPUT_DATA_DIR, "2d_test.tif") + with open_rasterio(test_file) as xds: + xds_1d = xds.rio.reproject( + "EPSG:4326", + src_geoloc_array=( + xds.coords["xc"].values, + xds.coords["yc"].values, + ), + georeferencing_convention="PIXEL_CENTER", + ) + assert xds_1d.coords["x"].shape == (14,) + assert xds_1d.coords["y"].shape == (10,) + xds_1d.rio.transform().almost_equals( + Affine( + 216.8587081056465, + 0.0, + 115698.25, + 0.0, + -216.8587081056465, + 2818720.0, + ) + )