diff --git a/src/eddymotion/data/filtering.py b/src/eddymotion/data/filtering.py
index d4fe4a33..e9a3e6e8 100644
--- a/src/eddymotion/data/filtering.py
+++ b/src/eddymotion/data/filtering.py
@@ -196,9 +196,9 @@ def downsample(
             smooth = datashape[:3] / shape[:3]
         data = gaussian_filter(data, smooth)
 
-    extents = (
-        apply_affine(imnii.affine, datashape - 0.5)
-        - apply_affine(imnii.affine, (-0.5, -0.5, -0.5))
+    extents = np.abs(
+        apply_affine(imnii.affine, datashape - 1)
+        - apply_affine(imnii.affine, (0.0, 0.0, 0.0))
     )
     newzooms = extents / shape
 
@@ -231,8 +231,7 @@ def downsample(
         data,
         locations.T,
         order=order,
-        mode="constant",
-        cval=0,
+        mode="mirror",
         prefilter=True,
     ).reshape(shape)
 
@@ -252,7 +251,6 @@ def decimate(
     in_file: str,
     factor: int | tuple[int, int, int],
     smooth: bool | tuple[int, int, int] = True,
-    order: int = 3,
     nonnegative: bool = True,
 ) -> Nifti1Image:
     """
@@ -278,9 +276,6 @@ def decimate(
         Alternatively, a tuple of three integers can be provided to specify
         different smoothing kernel sizes for each spatial dimension. Setting to
         False disables smoothing.
-    order : :obj:`int`, optional (default=3)
-        The order of the spline interpolation used for downsampling. Higher
-        orders provide smoother results but are computationally more expensive.
     nonnegative : :obj:`bool`, optional (default=``True``)
         If True, negative values in the downsampled data are set to zero.
 
diff --git a/test/test_filtering.py b/test/test_filtering.py
index 9c57f0e6..1923c182 100644
--- a/test/test_filtering.py
+++ b/test/test_filtering.py
@@ -38,23 +38,19 @@
 )
 @pytest.mark.parametrize(
     ("zoom_x", ),
-    # [(1.0, ), (-1.0, ), (2.0, ), (-2.0, )],
-    [(2.0,)],
+    [(1.0, ), (-1.0, ), (2.0, ), (-2.0, )],
 )
 @pytest.mark.parametrize(
     ("zoom_y", ),
-    # [(1.0, ), (-1.0, ), (2.0, ), (-2.0, )],
-    [(-2.0,)],
+    [(1.0, ), (-1.0, ), (2.0, ), (-2.0, )],
 )
 @pytest.mark.parametrize(
     ("zoom_z", ),
-    # [(1.0, ), (-1.0, ), (2.0, ), (-2.0, )],
-    [(-2.0,)],
+    [(1.0, ), (-1.0, ), (2.0, ), (-2.0, )],
 )
 @pytest.mark.parametrize(
     ("angle_x", ),
-    # [(0.0, ), (0.2, ), (-0.05, )],
-    [(-0.05,)]
+    [(0.0, ), (0.2, ), (-0.05, )],
 )
 @pytest.mark.parametrize(
     ("angle_y", ),
@@ -82,6 +78,7 @@ def test_decimation(
     angle_y,
     angle_z,
     offsets,
+    outdir,
 ):
     """Exercise decimation."""
 
@@ -120,8 +117,30 @@ def test_decimation(
     test_image.to_filename(fname)
 
     # Need to define test oracle. For now, just see if it doesn't smoke.
-    out = decimate(fname, factor=2, smooth=False, order=1)
-    out.to_filename(tmp_path / "decimated.nii.gz")
+    out = decimate(fname, factor=2, smooth=False)
+
+    out = downsample(fname, shape=(10, 10, 10), smooth=False, order=0)
+
+    if outdir:
+        from niworkflows.interfaces.reportlets.registration import (
+            SimpleBeforeAfterRPT as SimpleBeforeAfter,
+        )
+
+        out.to_filename(tmp_path / "decimated.nii.gz")
+
+        SimpleBeforeAfter(
+            after_label="Decimated",
+            before_label="Original",
+            after=str(tmp_path / "decimated.nii.gz"),
+            before=str(fname),
+            out_report=str(outdir / f'decimated-{tmp_path.name}.svg'),
+        ).run()
 
-    out = downsample(fname, shape=(10, 10, 10), smooth=False, order=1)
-    out.to_filename(tmp_path / "downsampled.nii.gz")
+        out.to_filename(tmp_path / "downsampled.nii.gz")
+        SimpleBeforeAfter(
+            after_label="Downsampled",
+            before_label="Original",
+            after=str(tmp_path / "downsampled.nii.gz"),
+            before=str(fname),
+            out_report=str(outdir / f'downsampled-{tmp_path.name}.svg'),
+        ).run()