diff --git a/tests/test_dem.py b/tests/test_dem.py
index c44aa352..e5b0ebdc 100644
--- a/tests/test_dem.py
+++ b/tests/test_dem.py
@@ -321,7 +321,7 @@ def test_to_vcrs__grids(self, grid_shifts: dict[str, Any]) -> None:
 
         # Using an arbitrary elevation of 100 m (no influence on the transformation)
         dem = DEM.from_array(
-            data=np.array([[100]]),
+            data=np.array([[100, 100]]),
             transform=rio.transform.from_bounds(
                 grid_shifts["lon"], grid_shifts["lat"], grid_shifts["lon"] + 0.01, grid_shifts["lat"] + 0.01, 0.01, 0.01
             ),
diff --git a/tests/test_spatialstats.py b/tests/test_spatialstats.py
index b067faa9..a6e3f170 100644
--- a/tests/test_spatialstats.py
+++ b/tests/test_spatialstats.py
@@ -853,6 +853,8 @@ def test_check_params_variogram_model(self) -> None:
     def test_estimate_model_spatial_correlation_and_infer_from_stable(self) -> None:
         """Test consistency of outputs and errors in wrapper functions for estimation of spatial correlation"""
 
+        warnings.filterwarnings("ignore", category=RuntimeWarning, message="Mean of empty slice")
+
         # Keep only data on stable
         diff_on_stable = self.diff.copy()
         diff_on_stable.set_mask(self.mask)
diff --git a/tests/test_terrain.py b/tests/test_terrain.py
index b08898e0..e542ac7a 100644
--- a/tests/test_terrain.py
+++ b/tests/test_terrain.py
@@ -21,6 +21,7 @@ def run_gdaldem(filepath: str, processing: str, options: str | None = None) -> M
     """Run GDAL's DEMProcessing and return the read numpy array."""
     # Rasterio strongly recommends against importing gdal along rio, so this is done here instead.
     from osgeo import gdal
+    gdal.UseExceptions()
 
     # Converting string into gdal processing options here to avoid import gdal outside this function:
     # Riley or Wilson for Terrain Ruggedness, and Zevenberg or Horn for slope, aspect and hillshade
@@ -123,7 +124,10 @@ def test_attribute_functions_against_gdaldem(self, attribute: str) -> None:
 
         # For hillshade, we round into an integer to match GDAL's output
         if attribute in ["hillshade_Horn", "hillshade_Zevenberg"]:
-            attr_xdem = attr_xdem.astype("int").astype("float32")
+            with warnings.catch_warnings():
+                # Normal that a warning would be raised here, so we catch it
+                warnings.filterwarnings("ignore", message="invalid value encountered in cast", category=RuntimeWarning)
+                attr_xdem = attr_xdem.astype("int").astype("float32")
 
         # We compute the difference and keep only valid values
         diff = (attr_xdem - attr_gdal).filled(np.nan)
@@ -171,9 +175,6 @@ def test_attribute_functions_against_gdaldem(self, attribute: str) -> None:
         # Validate that this doesn't raise weird warnings after introducing nans.
         functions[attribute](dem)
 
-    @pytest.mark.skip(
-        "richdem wheels don't build on latest GDAL versions, " "need to circumvent that problem..."
-    )  # type: ignore
     @pytest.mark.parametrize(
         "attribute",
         ["slope_Horn", "aspect_Horn", "hillshade_Horn", "curvature", "profile_curvature", "planform_curvature"],