fix: clarify what is a 3D transform chain and a 4D transform

3D transform chains resulting of composing several transformations
(e.g., affine and deformation fields in spatial normalization) should
not be split into its components.

This is in contrast to lists of 3D transforms such as head-motion
correcting affines, where each applies to one timepoint.
These should be considered 4D and in some future they may integrate
slice timing correction in them.

nitransforms/resampling.py

@@ -98,9 +98,10 @@ def apply(
     # Avoid opening the data array just yet
     input_dtype = get_obj_dtype(spatialimage.dataobj)
 
-    # Number of transformations
+    # Number of data volumes
     data_nvols = 1 if spatialimage.ndim < 4 else spatialimage.shape[-1]
-    xfm_nvols = len(transform)
+    # Number of transforms: transforms chains (e.g., affine + field, are a single transform)
+    xfm_nvols = 1 if transform.ndim < 4 else len(transform)
 
     if data_nvols != xfm_nvols and min(data_nvols, xfm_nvols) > 1:
         raise ValueError(