FIX: Segmentation plots aligned with cardinal axes

While nilearn/nilearn#25 is considered, this PR advances the generation
of mosaics after aligning the affine matrix of the data with the
cardinal axes. This should overcome the glitches when nilearn is
requested to plot contours and the mask for the ROI has voxels that only
partially intersect with the visualization plane.

Resolves: #542.
References: #281, #304.
Depends: #543.

Author: oesteban

niworkflows/utils/images.py

@@ -3,6 +3,16 @@
 import numpy as np
 
 
+def as_canonical(img):
+    """Drop rotation w.r.t. cardinal axes of input image."""
+    img = nb.as_closest_canonical(img)
+    zooms = list(img.header.get_zooms()[:3])
+    newaff = np.diag(zooms + [1])
+    rot = newaff[:3, :3].dot(np.linalg.inv(img.affine[:3, :3]))
+    newaff[:3, 3] = rot.dot(img.affine[:3, 3])
+    return nb.Nifti1Image(img.dataobj, newaff, img.header)
+
+
 def unsafe_write_nifti_header_and_data(fname, header, data):
     """Write header and data without any consistency checks or data munging
 

niworkflows/viz/utils.py

@@ -20,6 +20,7 @@
 
 from nipype.utils import filemanip
 from .. import NIWORKFLOWS_LOG
+from ..utils.images import as_canonical
 
 
 SVGNS = "http://www.w3.org/2000/svg"
@@ -222,19 +223,27 @@ def plot_segs(
     compress="auto",
     **plot_params
 ):
-    """ plot segmentation as contours over the image (e.g. anatomical).
+    """
+    Generate a static mosaic with ROIs represented by their delimiting contour.
+
+    Plot segmentation as contours over the image (e.g. anatomical).
     seg_niis should be a list of files. mask_nii helps determine the cut
     coordinates. plot_params will be passed on to nilearn plot_* functions. If
     seg_niis is a list of size one, it behaves as if it was plotting the mask.
     """
     plot_params = {} if plot_params is None else plot_params
 
-    image_nii = _3d_in_file(image_nii)
+    image_nii = as_canonical(_3d_in_file(image_nii))
+    seg_niis = [as_canonical(_3d_in_file(f)) for f in seg_niis]
     data = image_nii.get_fdata()
 
     plot_params = robust_set_limits(data, plot_params)
 
-    bbox_nii = nb.load(image_nii if bbox_nii is None else bbox_nii)
+    bbox_nii = (
+        image_nii if bbox_nii is None
+        else as_canonical(_3d_in_file(bbox_nii))
+    )
+
     if masked:
         bbox_nii = nlimage.threshold_img(bbox_nii, 1e-3)