docs/history.rst

@@ -1,6 +1,10 @@
 History
 =======
 
+0.15.5
+------
+- BUG:reproject: Allow rotated rasters (issue #746)
+
 0.15.4
 ------
 - BUG:reproject_match: Remove setting spatial dims on output resampled dataset (issue #768)

pyproject.toml

@@ -4,7 +4,7 @@ requires = ["setuptools", "wheel"]
 
 [project]
 name = "rioxarray"
-version = "0.15.4"
+version = "0.15.5"
 description = "geospatial xarray extension powered by rasterio"
 maintainers = [
     {name = "rioxarray Contributors"},

rioxarray/rioxarray.py

@@ -141,8 +141,10 @@ def _get_nonspatial_coords(
         src_data_array.rio.x_dim,
         src_data_array.rio.y_dim,
         DEFAULT_GRID_MAP,
+        "xc",
+        "yc",
     }:
-        if src_data_array[coord].dims:
+        if src_data_array[coord].ndim == 1:
             coords[coord] = xarray.IndexVariable(
                 src_data_array[coord].dims,
                 src_data_array[coord].values,
@@ -166,9 +168,13 @@ def _make_coords(
 ) -> dict[Hashable, Any]:
     """Generate the coordinates of the new projected `xarray.DataArray`"""
     coords = _get_nonspatial_coords(src_data_array)
-    if force_generate or (
-        src_data_array.rio.x_dim in src_data_array.coords
-        and src_data_array.rio.y_dim in src_data_array.coords
+    if (
+        force_generate
+        or (
+            src_data_array.rio.x_dim in src_data_array.coords
+            and src_data_array.rio.y_dim in src_data_array.coords
+        )
+        or ("xc" in src_data_array.coords and "yc" in src_data_array.coords)
     ):
         new_coords = _generate_spatial_coords(dst_affine, dst_width, dst_height)
         new_coords.update(coords)

test/integration/test_integration_rioxarray.py

@@ -28,7 +28,7 @@
     OneDimensionalRaster,
     RioXarrayError,
 )
-from rioxarray.rioxarray import _make_coords
+from rioxarray.rioxarray import _generate_spatial_coords, _make_coords
 from test.conftest import (
     GDAL_GE_361,
     TEST_COMPARE_DATA_DIR,
@@ -803,6 +803,60 @@ def test_reproject_match__geographic_dataset():
     assert_almost_equal(ds_match.y.values, ds.lat.values)
 
 
+def test_reproject_match__rotated():
+    rotated_affine = Affine(1, 0.2, 0, 0, 1, 0)
+    image = numpy.random.randint(0, 255, size=(100, 100), dtype=numpy.uint8)
+    rotated_dataset = (
+        xarray.DataArray(
+            image,
+            dims=("y", "x"),
+            coords=_generate_spatial_coords(
+                affine=rotated_affine,
+                width=image.shape[1],
+                height=image.shape[0],
+            ),
+        )
+        .astype("uint16")
+        .rio.write_nodata(0, inplace=True)
+        .rio.write_crs("EPSG:4326", inplace=True)
+        .rio.write_transform(rotated_affine, inplace=True)
+        .rio.write_coordinate_system(inplace=True)
+    )
+    squared_affine = Affine(1, 0, 0, 0, 1, 0)
+    squared_dataset = (
+        xarray.DataArray(
+            image,
+            dims=("y", "x"),
+            coords=_generate_spatial_coords(
+                affine=squared_affine,
+                width=image.shape[1],
+                height=image.shape[0],
+            ),
+        )
+        .astype("uint16")
+        .rio.write_nodata(0, inplace=True)
+        .rio.write_crs("EPSG:4326", inplace=True)
+        .rio.write_transform(squared_affine, inplace=True)
+        .rio.write_coordinate_system(inplace=True)
+    )
+    reprojected_to_squared_dataset = rotated_dataset.rio.reproject_match(
+        squared_dataset
+    )
+    reprojected_to_rotated_dataset = squared_dataset.rio.reproject_match(
+        rotated_dataset
+    )
+    assert (
+        "xc" in reprojected_to_rotated_dataset.coords
+        and "yc" in reprojected_to_rotated_dataset.coords
+    )
+    assert (
+        "x" in reprojected_to_squared_dataset.coords
+        and "y" in reprojected_to_squared_dataset.coords
+    )
+    assert reprojected_to_squared_dataset.rio.transform() == squared_affine
+    assert reprojected_to_rotated_dataset.rio.transform() == rotated_affine
+
+
 def test_interpolate_na__non_geospatial(dummy_dataset_non_geospatial):
     pytest.importorskip("scipy")
     ds = dummy_dataset_non_geospatial