diff --git a/docs/src/whatsnew/latest.rst b/docs/src/whatsnew/latest.rst index 737af1347e..cc19c8bd8f 100644 --- a/docs/src/whatsnew/latest.rst +++ b/docs/src/whatsnew/latest.rst @@ -36,7 +36,10 @@ This document explains the changes made to Iris for this release 🐛 Bugs Fixed ============= -#. N/A +#. `@PraneelBhatia`_ fixed :func:`iris.experimental.raster.export_geotiff`, + which failed with a ``TypeError`` under numpy >= 2.4: the geotransform + values are now passed to GDAL as Python scalars instead of 1-element + arrays. (:issue:`7133`, :pull:`7140`) 💣 Incompatible Changes @@ -142,6 +145,7 @@ This document explains the changes made to Iris for this release core dev names are automatically included by the common_links.inc: .. _@hdyson: https://github.com/hdyson +.. _@PraneelBhatia: https://github.com/PraneelBhatia .. comment diff --git a/lib/iris/experimental/raster.py b/lib/iris/experimental/raster.py index 0b5057136c..ea10f6ab7b 100644 --- a/lib/iris/experimental/raster.py +++ b/lib/iris/experimental/raster.py @@ -173,7 +173,10 @@ def export_geotiff(cube, fname): ) if not coord.is_contiguous(): raise ValueError("Coordinate {!r} bounds must be contiguous.".format(name)) - xy_step.append(np.diff(coord.bounds[0])) + # Take a Python scalar step value : GDAL's SetGeoTransform requires + # plain numbers, and numpy >= 2.4 no longer implicitly converts the + # 1-element array which np.diff produces here. + xy_step.append(np.diff(coord.bounds[0]).item()) if not np.allclose(np.diff(coord.bounds), xy_step[-1]): msg = "Coordinate {!r} bounds must be regularly spaced.".format(name) raise ValueError(msg) @@ -205,6 +208,6 @@ def export_geotiff(cube, fname): x_bounds = x_bounds.copy() x_bounds[big_indices] -= 360 - x_min = np.min(x_bounds) - y_max = np.max(coord_y.bounds) + x_min = np.min(x_bounds).item() + y_max = np.max(coord_y.bounds).item() _gdal_write_array(x_min, x_step, y_max, y_step, coord_system, data, fname, "GTiff") diff --git a/lib/iris/tests/unit/experimental/raster/test_export_geotiff.py b/lib/iris/tests/unit/experimental/raster/test_export_geotiff.py index e1c76f50de..62e4276b63 100644 --- a/lib/iris/tests/unit/experimental/raster/test_export_geotiff.py +++ b/lib/iris/tests/unit/experimental/raster/test_export_geotiff.py @@ -162,7 +162,7 @@ def test_ellipsoid(self, tmp_path): @_shared_utils.skip_gdal class TestGeoTransform: - def test_(self, tmp_path): + def _cube(self): data = np.arange(12).reshape(3, 4).astype(np.uint8) cube = Cube(data, "air_pressure_anomaly") coord = DimCoord([30, 40, 50], "latitude", units="degrees") @@ -171,7 +171,23 @@ def test_(self, tmp_path): coord = DimCoord([-10, -5, 0, 5], "longitude", units="degrees") coord.guess_bounds() cube.add_dim_coord(coord, 1) + return cube + + def test_geotransform(self, tmp_path): + cube = self._cube() temp_filename = str(tmp_path / "tmp.tif") export_geotiff(cube, temp_filename) dataset = gdal.Open(temp_filename, gdal.GA_ReadOnly) assert dataset.GetGeoTransform() == (-12.5, 5, 0, 55, 0, -10) + + def test_scalar_transform_values(self, mocker, tmp_path): + # The geotransform values passed to GDAL must be Python scalars : + # numpy >= 2.4 no longer implicitly converts the 1-element arrays + # which were previously passed. See issue #7133. + write_array = mocker.patch("iris.experimental.raster._gdal_write_array") + export_geotiff(self._cube(), str(tmp_path / "tmp.tif")) + x_min, x_step, y_max, y_step = write_array.call_args.args[:4] + assert (x_min, x_step, y_max, y_step) == (-12.5, 5.0, 55.0, -10.0) + for value in (x_min, x_step, y_max, y_step): + assert isinstance(value, (int, float)) + assert not isinstance(value, np.generic)