ENH: Add get/set_zooms(units='norm') to manipulate zooms in mm/s units

nibabel/analyze.py

@@ -399,7 +399,7 @@ def from_header(klass, header=None, check=True):
                                   f"but output header {klass} does not support it")
         obj.set_data_dtype(header.get_data_dtype())
         obj.set_data_shape(header.get_data_shape())
-        obj.set_zooms(header.get_zooms())
+        obj.set_zooms(header.get_zooms(units='raw'), units='raw')
         if check:
             obj.check_fix()
         return obj
@@ -661,13 +661,22 @@ def get_base_affine(self):
 
     get_best_affine = get_base_affine
 
-    def get_zooms(self):
-        """ Get zooms from header
+    def get_zooms(self, units='norm', raise_unknown=False):
+        """ Get zooms (spacing between voxels along each axis) from header
+
+        Parameters
+        ----------
+        units : {'norm', 'raw'}, optional
+            Return zooms in normalized units of mm/sec for spatial/temporal or
+            as raw values stored in header.
+        raise_unkown : bool, optional
+            If normalized units are requested and the units are ambiguous, raise
+            a ``ValueError``
 
         Returns
         -------
-        z : tuple
-           tuple of header zoom values
+        zooms : tuple
+            tuple of header zoom values
 
         Examples
         --------
@@ -681,6 +690,8 @@ def get_zooms(self):
         >>> hdr.get_zooms()
         (3.0, 4.0)
         """
+        if units not in ('norm', 'raw'):
+            raise ValueError("`units` parameter must be 'norm' or 'raw'")
         hdr = self._structarr
         dims = hdr['dim']
         ndim = dims[0]
@@ -689,11 +700,22 @@ def get_zooms(self):
         pixdims = hdr['pixdim']
         return tuple(pixdims[1:ndim + 1])
 
-    def set_zooms(self, zooms):
+    def set_zooms(self, zooms, units='norm'):
         """ Set zooms into header fields
 
         See docstring for ``get_zooms`` for examples
+
+        Parameters
+        ----------
+        zooms : sequence of floats
+            Zoom values to set in header
+        units : {'norm', 'raw'}, optional
+            Zooms are specified in normalized units of mm/sec for
+            spatial/temporal or as raw values to be interpreted according to
+            format specification.
         """
+        if units not in ('norm', 'raw'):
+            raise ValueError("`units` parameter must be 'norm' or 'raw'")
         hdr = self._structarr
         dims = hdr['dim']
         ndim = dims[0]

nibabel/ecat.py

@@ -578,8 +578,10 @@ def get_frame_affine(self, frame=0):
                                                               z_off])
         return aff
 
-    def get_zooms(self, frame=0):
+    def get_zooms(self, frame=0, units='norm', raise_unknown=False):
         """returns zooms  ...pixdims"""
+        if units not in ('norm', 'raw'):
+            raise ValueError("`units` parameter must be 'norm' or 'raw'")
         subhdr = self.subheaders[frame]
         x_zoom = subhdr['x_pixel_size'] * 10
         y_zoom = subhdr['y_pixel_size'] * 10

nibabel/freesurfer/mghformat.py

@@ -237,40 +237,75 @@ def _ndims(self):
         """
         return 3 + (self._structarr['dims'][3] > 1)
 
-    def get_zooms(self):
+    def get_zooms(self, units='norm', raise_unknown=False):
         """ Get zooms from header
 
         Returns the spacing of voxels in the x, y, and z dimensions.
         For four-dimensional files, a fourth zoom is included, equal to the
-        repetition time (TR) in ms (see `The MGH/MGZ Volume Format
-        <mghformat>`_).
+        repetition time (TR).
+        TR is stored in milliseconds (see `The MGH/MGZ Volume Format
+        <mghformat>`_),
+        so if ``units == 'raw'``, the fourth zoom will be in ms.
+        If ``units == 'norm'`` (default), the fourth zoom will be in
+        seconds.
 
         To access only the spatial zooms, use `hdr['delta']`.
 
+        Parameters
+        ----------
+        units : {'norm', 'raw'}, optional
+            Return zooms in normalized units of mm/sec for spatial/temporal or
+            as raw values stored in header.
+        raise_unkown : bool, optional
+            If normalized units are requested and the units are ambiguous, raise
+            a ``ValueError``
+
         Returns
         -------
-        z : tuple
-           tuple of header zoom values
+        zooms : tuple
+            tuple of header zoom values
 
         .. _mghformat: https://surfer.nmr.mgh.harvard.edu/fswiki/FsTutorial/MghFormat#line-82
         """
+        if units == 'norm':
+            tfactor = 0.001
+        elif units == 'raw':
+            tfactor = 1
+        else:
+            raise ValueError("`units` parameter must be 'norm' or 'raw'")
+
         # Do not return time zoom (TR) if 3D image
-        tzoom = (self['tr'],) if self._ndims() > 3 else ()
+        tzoom = (self['tr'] * tfactor,) if self._ndims() > 3 else ()
         return tuple(self._structarr['delta']) + tzoom
 
-    def set_zooms(self, zooms):
+    def set_zooms(self, zooms, units='norm'):
         """ Set zooms into header fields
 
         Sets the spacing of voxels in the x, y, and z dimensions.
-        For four-dimensional files, a temporal zoom (repetition time, or TR, in
-        ms) may be provided as a fourth sequence element.
+        For four-dimensional files, a temporal zoom (repetition time, or TR)
+        may be provided as a fourth sequence element.
+
+        TR is stored in milliseconds (see `The MGH/MGZ Volume Format
+        <mghformat>`_),
+        so if ``units == 'raw'``, the fourth zoom will be interpreted as ms
+        and stored unmodified.
+        If ``units == 'norm'`` (default), the fourth zoom will be interpreted
+        as seconds, and converted to ms before storing in the header.
 
         Parameters
         ----------
         zooms : sequence
             sequence of floats specifying spatial and (optionally) temporal
             zooms
+        units : {'norm', 'raw'}, optional
+            Zooms are specified in normalized units of mm/sec for
+            spatial/temporal dimensions or as raw values to be stored in
+            header.
+
+        .. _mghformat: https://surfer.nmr.mgh.harvard.edu/fswiki/FsTutorial/MghFormat#line-82
         """
+        if units not in ('norm', 'raw'):
+            raise ValueError("`units` parameter must be 'norm' or 'raw'")
         hdr = self._structarr
         zooms = np.asarray(zooms)
         ndims = self._ndims()
@@ -283,7 +318,13 @@ def set_zooms(self, zooms):
         if len(zooms) == 4:
             if zooms[3] < 0:
                 raise HeaderDataError(f'TR must be non-negative; got {zooms[3]}')
-            hdr['tr'] = zooms[3]
+            if units == 'norm':
+                tfactor = 1000
+            elif units == 'raw':
+                tfactor = 1
+            else:
+                raise ValueError("`units` parameter must be 'norm' or 'raw'")
+            hdr['tr'] = zooms[3] * tfactor
 
     def get_data_shape(self):
         """ Get shape of data

nibabel/freesurfer/tests/test_mghformat.py

@@ -74,7 +74,8 @@ def test_read_mgh():
     assert_almost_equal(h['flip_angle'], 0.0)
     assert_almost_equal(h['te'], 0.0)
     assert_almost_equal(h['ti'], 0.0)
-    assert_array_almost_equal(h.get_zooms(), [1, 1, 1, 2])
+    assert_array_almost_equal(h.get_zooms(units='raw'), [1, 1, 1, 2])
+    assert_array_almost_equal(h.get_zooms(units='norm'), [1, 1, 1, 0.002])
     assert_array_almost_equal(h.get_vox2ras(), v2r)
     assert_array_almost_equal(h.get_vox2ras_tkr(), v2rtkr)
 
@@ -147,9 +148,14 @@ def test_write_noaffine_mgh():
 def test_set_zooms():
     mgz = load(MGZ_FNAME)
     h = mgz.header
-    assert_array_almost_equal(h.get_zooms(), [1, 1, 1, 2])
-    h.set_zooms([1, 1, 1, 3])
-    assert_array_almost_equal(h.get_zooms(), [1, 1, 1, 3])
+    assert_array_almost_equal(h.get_zooms(units='raw'), [1, 1, 1, 2])
+    assert_array_almost_equal(h.get_zooms(units='norm'), [1, 1, 1, 0.002])
+    h.set_zooms([1, 1, 1, 3], units='raw')
+    assert_array_almost_equal(h.get_zooms(units='raw'), [1, 1, 1, 3])
+    assert_array_almost_equal(h.get_zooms(units='norm'), [1, 1, 1, 0.003])
+    h.set_zooms([1, 1, 1, 3], units='norm')
+    assert_array_almost_equal(h.get_zooms(units='raw'), [1, 1, 1, 3000])
+    assert_array_almost_equal(h.get_zooms(units='norm'), [1, 1, 1, 3])
     for zooms in ((-1, 1, 1, 1),
                   (1, -1, 1, 1),
                   (1, 1, -1, 1),
@@ -351,6 +357,63 @@ def check_dtypes(self, expected, actual):
         # MGH requires the actual to be a big endian version of expected
         assert expected.newbyteorder('>') == actual
 
+    def test_zooms_edge_cases(self):
+        img_klass = self.image_class
+        aff = np.eye(4)
+        arr = np.arange(120, dtype=np.int16).reshape((2, 3, 4, 5))
+        img = img_klass(arr, aff)
+
+        assert_array_almost_equal(img.header.get_zooms(units='raw'),
+                                  (1, 1, 1, 0))
+        assert_array_almost_equal(img.header.get_zooms(units='norm'),
+                                  (1, 1, 1, 0))
+
+        img.header.set_zooms((1, 1, 1, 2000), units='raw')
+        assert_array_almost_equal(img.header.get_zooms(units='raw'),
+                                  (1, 1, 1, 2000))
+        assert_array_almost_equal(img.header.get_zooms(units='norm'),
+                                  (1, 1, 1, 2))
+        assert_array_almost_equal(img.header.get_zooms(), (1, 1, 1, 2))
+
+        img.header.set_zooms((2, 2, 2, 3), units='norm')
+        assert_array_almost_equal(img.header.get_zooms(units='raw'),
+                                  (2, 2, 2, 3000))
+        assert_array_almost_equal(img.header.get_zooms(units='norm'),
+                                  (2, 2, 2, 3))
+        assert_array_almost_equal(img.header.get_zooms(), (2, 2, 2, 3))
+
+        # It's legal to set zooms for spatial dimensions only
+        img.header.set_zooms((3, 3, 3), units='norm')
+        assert_array_almost_equal(img.header.get_zooms(units='raw'),
+                                  (3, 3, 3, 3000))
+        assert_array_almost_equal(img.header.get_zooms(units='norm'),
+                                  (3, 3, 3, 3))
+        assert_array_almost_equal(img.header.get_zooms(), (3, 3, 3, 3))
+
+        arr = np.arange(24, dtype=np.int16).reshape((2, 3, 4))
+        img = img_klass(arr, aff)
+
+        assert_array_almost_equal(img.header.get_zooms(units='raw'),
+                                                       (1, 1, 1))
+        assert_array_almost_equal(img.header.get_zooms(units='norm'),
+                                                       (1, 1, 1))
+
+        img.header.set_zooms((2, 2, 2), units='raw')
+        assert_array_almost_equal(img.header.get_zooms(units='raw'),
+                                                       (2, 2, 2))
+        assert_array_almost_equal(img.header.get_zooms(units='norm'),
+                                                       (2, 2, 2))
+
+        img.header.set_zooms((3, 3, 3), units='norm')
+        assert_array_almost_equal(img.header.get_zooms(units='raw'),
+                                                       (3, 3, 3))
+        assert_array_almost_equal(img.header.get_zooms(units='norm'),
+                                                       (3, 3, 3))
+
+        # Cannot set TR as zoom for 3D image
+        assert_raises(HeaderDataError, img.header.set_zooms, (4, 4, 4, 5), 'raw')
+        assert_raises(HeaderDataError, img.header.set_zooms, (4, 4, 4, 5), 'norm')
+
 
 class TestMGHHeader(tws._TestLabeledWrapStruct):
     header_class = MGHHeader

nibabel/minc1.py

@@ -90,8 +90,10 @@ def get_data_dtype(self):
     def get_data_shape(self):
         return self._image.data.shape
 
-    def get_zooms(self):
+    def get_zooms(self, units='norm', raise_unknown=False):
         """ Get real-world sizes of voxels """
+        if units not in ('norm', 'raw'):
+            raise ValueError("`units` parameter must be 'norm' or 'raw'")
         # zooms must be positive; but steps in MINC can be negative
         return tuple([abs(float(dim.step)) if hasattr(dim, 'step') else 1.0
                       for dim in self._dims])