[ENH] Update FreeSurfer to 7.3.1, dmri-amico to 1.5.4

.circleci/config.yml

@@ -2,7 +2,7 @@ version: 2.1
 
 .dockersetup: &dockersetup
   docker:
-    - image: pennbbl/qsiprep_build:23.2.4
+    - image: pennbbl/qsiprep_build:23.3.1
   working_directory: /src/qsiprep
 
 runinstall: &runinstall

Dockerfile

@@ -1,4 +1,4 @@
-FROM pennbbl/qsiprep_build:23.2.4
+FROM pennbbl/qsiprep_build:23.3.1
 
 # WORKDIR /root/
 # # Installing qsiprep

qsiprep/interfaces/converters.py

@@ -72,7 +72,7 @@ def _run_interface(self, runtime):
         if isdefined(mask_file):
             mask_img = nb.load(mask_file)
         else:
-            ampl_data = amplitudes_img.get_data()
+            ampl_data = amplitudes_img.get_fdata()
             ampl_mask = ampl_data.sum(3) > 1e-6
             mask_img = nb.Nifti1Image(ampl_mask.astype(np.float),
                                       amplitudes_img.affine)
@@ -179,7 +179,7 @@ def amplitudes_to_fibgz(amplitudes_img, odf_dirs, odf_faces, output_file,
     odf_dirs: np.ndarray
         N x 3 array containing the directions corresponding to the
         amplitudes in ``amplitudes_img``. The values in
-        ``amplitudes_img.get_data()[..., i]`` are for the
+        ``amplitudes_img.get_fdata()[..., i]`` are for the
         direction in ``odf_dirs[i]``.
     odf_faces: np.ndarray
         triangles connecting the vertices in ``odf_dirs``
@@ -208,8 +208,8 @@ def amplitudes_to_fibgz(amplitudes_img, odf_dirs, odf_faces, output_file,
         raise ValueError("Differing grid between mask and amplitudes")
 
     # Get the flat mask
-    ampl_data = amplitudes_img.get_data()
-    flat_mask = mask_img.get_data().flatten(order="F") > 0
+    ampl_data = amplitudes_img.get_fdata()
+    flat_mask = mask_img.get_fdata().flatten(order="F") > 0
     odf_array = ampl_data.reshape(-1, ampl_data.shape[3], order='F')
     del ampl_data
     masked_odfs = odf_array[flat_mask, :]
@@ -285,7 +285,7 @@ def amico_directions_to_fibgz(directions_img, od_img, icvf_img, isovf_img,
     odf_dirs: np.ndarray
         N x 3 array containing the directions corresponding to the
         amplitudes in ``amplitudes_img``. The values in
-        ``amplitudes_img.get_data()[..., i]`` are for the
+        ``amplitudes_img.get_fdata()[..., i]`` are for the
         direction in ``odf_dirs[i]``.
     odf_faces: np.ndarray
         triangles connecting the vertices in ``odf_dirs``
@@ -361,7 +361,7 @@ def amplitudes_to_sh_mif(amplitudes_img, odf_dirs, output_file, working_dir):
     odf_dirs: np.ndarray
         2*N x 3 array containing the directions corresponding to the
         amplitudes in ``amplitudes_img``. The values in
-        ``amplitudes_img.get_data()[..., i]`` are for the
+        ``amplitudes_img.get_fdata()[..., i]`` are for the
         direction in ``odf_dirs[i]``. Here the second half of the
         directions are the opposite of the fist and therefore have the
         same amplitudes.

qsiprep/interfaces/dipy.py

@@ -214,18 +214,18 @@ def _run_interface(self, runtime):
         in_file = self.inputs.in_file
 
         uneq_img = nb.load(in_file)
-        uneq_data = uneq_img.get_data()
+        uneq_data = uneq_img.get_fdata()
 
         mask = nb.load(self.inputs.mask_file)
-        bool_mask = mask.get_data() > 0
+        bool_mask = mask.get_fdata() > 0
         data_voxels = uneq_data[bool_mask]
 
         # Do a clip on 2 to 98th percentile
         bottom_2, top_98 = np.percentile(data_voxels, np.array([1, 99]),
                                          axis=None)
         clipped_b0 = np.clip(data_voxels, 0, top_98)
         eq_data = histeq(clipped_b0, num_bins=512)
-        output = np.zeros_like(mask.get_data())
+        output = np.zeros_like(mask.get_fdata())
         output[bool_mask] = eq_data
         eq_img = nb.Nifti1Image(output, uneq_img.affine, uneq_img.header)
         self._results['out_file'] = fname_presuffix(
@@ -279,7 +279,7 @@ def _get_gtab(self, external_bvals=None, external_bvecs=None):
 
     def _get_mask(self, amplitudes_img, gtab):
         if not isdefined(self.inputs.mask_file):
-            dwi_data = amplitudes_img.get_data()
+            dwi_data = amplitudes_img.get_fdata()
             LOGGER.warning("Creating an Otsu mask, check that the whole brain is covered.")
             _, mask_array = median_otsu(dwi_data,
                                         vol_idx=gtab.b0s_mask,
@@ -292,7 +292,7 @@ def _get_mask(self, amplitudes_img, gtab):
                                       amplitudes_img.header)
         else:
             mask_img = nb.load(self.inputs.mask_file)
-            mask_array = mask_img.get_data() > 0
+            mask_array = mask_img.get_fdata() > 0
         return mask_img, mask_array
 
     def _save_scalar(self, data, suffix, runtime, ref_img):
@@ -732,7 +732,7 @@ def _run_interface(self, runtime):
                                        newpath=runtime.cwd, use_ext=True)
             nb.Nifti1Image(data, dwi_img.affine).to_filename(out_name)
             self._results[metric] = out_name
-        
+
         # Get the kurtosis metrics
         for metric in ["mk", "ak", "rk", "mkt"]:
             data = np.nan_to_num(

qsiprep/interfaces/dsi_studio.py

@@ -19,8 +19,8 @@
 
 class DSIStudioCommandLineInputSpec(CommandLineInputSpec):
     num_threads = traits.Int(
-        1, 
-        usedefault=True, 
+        1,
+        usedefault=True,
         argstr="--thread_count=%d",
         nohash=True)
 
@@ -640,7 +640,7 @@ def _run_interface(self, runtime):
 
         # No matter what, still use the correct affine
         nb.Nifti1Image(
-            reoriented_img.get_data()[::-1, ::-1, :],
+            reoriented_img.get_fdata()[::-1, ::-1, :],
             correct_img.affine).to_filename(new_file)
         self._results['out_file'] = new_file
 

qsiprep/interfaces/fmap.py

@@ -354,14 +354,14 @@ def _unwrap(fmap_data, mag_file, mask=None):
 
     nb.Nifti1Image(fmap_data, magnii.affine).to_filename('fmap_rad.nii.gz')
     nb.Nifti1Image(mask, magnii.affine).to_filename('fmap_mask.nii.gz')
-    nb.Nifti1Image(magnii.get_data(), magnii.affine).to_filename('fmap_mag.nii.gz')
+    nb.Nifti1Image(magnii.get_fdata(), magnii.affine).to_filename('fmap_mag.nii.gz')
 
     # Run prelude
     res = PRELUDE(phase_file='fmap_rad.nii.gz',
                   magnitude_file='fmap_mag.nii.gz',
                   mask_file='fmap_mask.nii.gz').run()
 
-    unwrapped = nb.load(res.outputs.unwrapped_phase_file).get_data() * (fmapmax / pi)
+    unwrapped = nb.load(res.outputs.unwrapped_phase_file).get_fdata() * (fmapmax / pi)
     return unwrapped
 
 
@@ -539,7 +539,7 @@ def _torads(in_file, fmap_range=None, newpath=None):
 
     out_file = fname_presuffix(in_file, suffix='_rad', newpath=newpath)
     fmapnii = nb.load(in_file)
-    fmapdata = fmapnii.get_data()
+    fmapdata = fmapnii.get_fdata()
 
     if fmap_range is None:
         fmap_range = max(abs(fmapdata.min()), fmapdata.max())
@@ -558,7 +558,7 @@ def _tohz(in_file, range_hz, newpath=None):
 
     out_file = fname_presuffix(in_file, suffix='_hz', newpath=newpath)
     fmapnii = nb.load(in_file)
-    fmapdata = fmapnii.get_data()
+    fmapdata = fmapnii.get_fdata()
     fmapdata = fmapdata * (range_hz / pi)
     out_img = nb.Nifti1Image(fmapdata, fmapnii.affine, fmapnii.header)
     out_img.set_data_dtype('float32')
@@ -592,7 +592,7 @@ def phdiff2fmap(in_file, delta_te, newpath=None):
 
     out_file = fname_presuffix(in_file, suffix='_fmap', newpath=newpath)
     image = nb.load(in_file)
-    data = (image.get_data().astype(np.float32) / (2. * math.pi * delta_te))
+    data = (image.get_fdata().astype(np.float32) / (2. * math.pi * delta_te))
     nii = nb.Nifti1Image(data, image.affine, image.header)
     nii.set_data_dtype(np.float32)
     nii.to_filename(out_file)

qsiprep/interfaces/freesurfer.py

@@ -307,9 +307,9 @@ def _run_interface(self, runtime):
 
         anatnii = nb.load(self.inputs.in_anat)
         msknii = nb.Nifti1Image(
-            grow_mask(anatnii.get_data(),
-                      nb.load(self.inputs.in_aseg).get_data(),
-                      nb.load(self.inputs.in_ants).get_data()),
+            grow_mask(anatnii.get_fdata(),
+                      nb.load(self.inputs.in_aseg).get_fdata(),
+                      nb.load(self.inputs.in_ants).get_fdata()),
             anatnii.affine,
             anatnii.header
         )
@@ -354,9 +354,9 @@ def inject_skullstripped(subjects_dir, subject_id, skullstripped):
     if not op.exists(bm_auto):
         img = nb.load(t1)
         mask = nb.load(skullstripped)
-        bmask = new_img_like(mask, mask.get_data() > 0)
+        bmask = new_img_like(mask, mask.get_fdata() > 0)
         resampled_mask = resample_to_img(bmask, img, 'nearest')
-        masked_image = new_img_like(img, img.get_data() * resampled_mask.get_data())
+        masked_image = new_img_like(img, img.get_fdata() * resampled_mask.get_fdata())
         masked_image.to_filename(bm_auto)
 
     if not op.exists(bm):
@@ -584,7 +584,7 @@ class _SynthStripOutputSpec(TraitedSpec):
 class SynthStrip(FSCommandOpenMP):
     input_spec = _SynthStripInputSpec
     output_spec = _SynthStripOutputSpec
-    _cmd = os.getenv("FREESURFER_HOME") + "/mri_synthstrip"
+    _cmd = os.getenv("FREESURFER_HOME") + "/bin/mri_synthstrip"
 
     def _num_threads_update(self):
         if self.inputs.num_threads:

qsiprep/interfaces/gradients.py

@@ -108,7 +108,7 @@ def is_unique_sample(vec):
         np.savetxt(output_bval, unique_bvals, fmt='%d', newline=' ')
         unique_bvecs = bvecs[unique_indices]
         np.savetxt(output_bvec, unique_bvecs.T, fmt='%.8f')
-        unique_data = original_image.get_data()[..., unique_indices]
+        unique_data = original_image.get_fdata()[..., unique_indices]
         nb.Nifti1Image(unique_data, original_image.affine, original_image.header
                        ).to_filename(output_nii)
         self._results['bval_file'] = output_bval
@@ -147,7 +147,7 @@ def _run_interface(self, runtime):
 
         output_npz = os.path.join(runtime.cwd, "slice_stats.npz")
         mask_img = nb.load(self.inputs.mask_image)
-        mask = mask_img.get_data() > 0
+        mask = mask_img.get_fdata() > 0
         masked_slices = (mask * np.arange(mask_img.shape[2])[np.newaxis, np.newaxis, :]
                          ).astype(np.int)
         slice_nums, slice_counts = np.unique(masked_slices[mask], return_counts=True)
@@ -813,7 +813,7 @@ def create_tensor_image(mask_img, direction, prefix):
     for direction in ['xx', 'xy', 'xz', 'yy', 'yz', 'zz']:
         this_component = prefix + '_temp_dtiComp_%s.nii.gz' % direction
         LOGGER.info("writing %s", this_component)
-        nb.Nifti1Image(mask_img.get_data()*tensor[tensor_index[direction]], mask_img.affine,
+        nb.Nifti1Image(mask_img.get_fdata()*tensor[tensor_index[direction]], mask_img.affine,
                        mask_img.header).to_filename(this_component)
         temp_components.append(this_component)
 
@@ -841,9 +841,9 @@ def reorient_tensor_image(tensor_image, warp_file, mask_img, prefix, output_fnam
 
     # Load the reoriented tensor and get the principal directions out
     reoriented_dt_img = nb.load(reoriented_tensor_fname)
-    reoriented_tensor_data = reoriented_dt_img.get_data().squeeze()
+    reoriented_tensor_data = reoriented_dt_img.get_fdata().squeeze()
 
-    mask_data = mask_img.get_data() > 0
+    mask_data = mask_img.get_fdata() > 0
     output_data = np.zeros(mask_img.shape + (3,))
 
     reoriented_tensors = reoriented_tensor_data[mask_data]
@@ -882,7 +882,7 @@ def local_bvec_rotation(original_bvecs, warp_transforms, mask_image, runtime, ou
     """Create a vector in each voxel that accounts for nonlinear warps."""
     prefix = os.path.join(runtime.cwd, "local_bvec_")
     mask_img = nb.load(mask_image)
-    mask_data = mask_img.get_data()
+    mask_data = mask_img.get_fdata()
     b0_image = get_vector_nii(np.stack([np.zeros_like(mask_data)] * 3, -1), mask_img.affine,
                               mask_img.header)
     commands = []
@@ -907,7 +907,7 @@ def local_bvec_rotation(original_bvecs, warp_transforms, mask_image, runtime, ou
         commands.append(rotate_cmd)
         rotated_vec_files.append(out_fname)
     concatenated = np.stack(
-        [nb.load(img, mmap=False).get_data().astype("<f4") for img in rotated_vec_files], -1)
+        [nb.load(img, mmap=False).get_fdata().astype("<f4") for img in rotated_vec_files], -1)
     nb.Nifti1Image(concatenated, mask_img.affine, mask_img.header).to_filename(output_fname)
     for temp_file in rotated_vec_files:
         os.remove(temp_file)

qsiprep/interfaces/images.py

@@ -135,7 +135,7 @@ def _run_interface(self, runtime):
 
         if len(in_files) == 1:
             filenii = nb.load(in_files[0])
-            filedata = filenii.get_data()
+            filedata = filenii.get_fdata()
 
             # magnitude files can have an extra dimension empty
             if filedata.ndim == 5:
@@ -148,7 +148,7 @@ def _run_interface(self, runtime):
                     nb.Nifti1Image(sqdata, filenii.affine,
                                    filenii.header).to_filename(in_files[0])
 
-            if np.squeeze(nb.load(in_files[0]).get_data()).ndim < 4:
+            if np.squeeze(nb.load(in_files[0]).get_fdata()).ndim < 4:
                 self._results['out_file'] = in_files[0]
                 self._results['out_avg'] = in_files[0]
                 # TODO: generate identity out_mats and zero-filled out_movpar
@@ -165,7 +165,7 @@ def _run_interface(self, runtime):
         self._results['out_file'] = mcres.outputs.out_file
 
         hmcnii = nb.load(mcres.outputs.out_file)
-        hmcdat = hmcnii.get_data().mean(axis=3)
+        hmcdat = hmcnii.get_fdata().mean(axis=3)
         if self.inputs.zero_based_avg:
             hmcdat -= hmcdat.min()
 
@@ -264,7 +264,7 @@ def _run_interface(self, runtime):
                 offset = (reoriented.affine[:3, 3] * size_factor - reoriented.affine[:3, 3])
                 target_affine[:3, 3] = reoriented.affine[:3, 3] + offset.astype(int)
 
-            data = nli.resample_img(reoriented, target_affine, target_shape).get_data()
+            data = nli.resample_img(reoriented, target_affine, target_shape).get_fdata()
             conform_xfm = np.linalg.inv(reoriented.affine).dot(target_affine)
             reoriented = reoriented.__class__(data, target_affine, reoriented.header)
 

qsiprep/interfaces/nilearn.py

@@ -82,13 +82,13 @@ def _run_interface(self, runtime):
         )
 
         if self.inputs.closing:
-            closed = sim.binary_closing(masknii.get_data().astype(
+            closed = sim.binary_closing(masknii.get_fdata().astype(
                 np.uint8), sim.ball(1)).astype(np.uint8)
             masknii = masknii.__class__(closed, masknii.affine,
                                         masknii.header)
 
         if self.inputs.fill_holes:
-            filled = binary_fill_holes(masknii.get_data().astype(
+            filled = binary_fill_holes(masknii.get_fdata().astype(
                 np.uint8), sim.ball(6)).astype(np.uint8)
             masknii = masknii.__class__(filled, masknii.affine,
                                         masknii.header)
@@ -139,7 +139,7 @@ def _run_interface(self, runtime):
             new_nii.header.set_zooms(list(new_nii.header.get_zooms()[:3]) +
                                      [src_hdr.get_zooms()[3]])
         if self.inputs.is_dwi:
-            new_nii = nb.Nifti1Image(np.abs(new_nii.get_data()), new_nii.affine, new_nii.header)
+            new_nii = nb.Nifti1Image(np.abs(new_nii.get_fdata()), new_nii.affine, new_nii.header)
 
         new_nii.to_filename(self._results['out_file'])
 
@@ -317,7 +317,7 @@ def _enhance_t2_contrast(in_file, newpath=None, offset=0.5):
     out_file = fname_presuffix(in_file, suffix='_t1enh',
                                newpath=newpath)
     nii = nb.load(in_file)
-    data = nii.get_data()
+    data = nii.get_fdata()
     maxd = data.max()
     newdata = np.log(offset + data / maxd)
     newdata -= newdata.min()

qsiprep/interfaces/niworkflows.py

@@ -51,7 +51,7 @@ def __init__(self, sliceqc_file, mask_file, confounds, usecols=None, units=None,
             slice_scores = np.loadtxt(sliceqc_file, skiprows=1)
             # Get the slice counts
             mask_img = nb.load(mask_file)
-            mask = mask_img.get_data() > 0
+            mask = mask_img.get_fdata() > 0
             masked_slices = (mask * np.arange(mask_img.shape[2])[np.newaxis, np.newaxis, :]
                              ).astype(np.int)
             slice_nums, slice_counts = np.unique(masked_slices[mask], return_counts=True)

qsiprep/interfaces/qc.py

@@ -64,7 +64,7 @@ def mplfigcontour(data, outfile=None, as_bytes=False):
 
 def load_and_reorient(filename):
     img = nib.load(filename)
-    data, aff = img.get_data(), img.affine
+    data, aff = img.get_fdata(), img.affine
     data = reorient_array(data, aff)
     return data
 

qsiprep/interfaces/shoreline.py

@@ -175,7 +175,7 @@ def _run_interface(self, runtime):
         pred_val = self.inputs.bval_to_predict
         # Load the mask image:
         mask_img = nb.load(self.inputs.aligned_mask)
-        mask_array = mask_img.get_data() > 1e-6
+        mask_array = mask_img.get_fdata() > 1e-6
         all_images = self.inputs.aligned_dwis
         if isinstance(self.inputs.aligned_bvecs, np.ndarray):
             bvecs = self.inputs.aligned_bvecs
@@ -217,7 +217,7 @@ def _run_interface(self, runtime):
             shore_array = shore_fit._shore_coef[mask_array]
             output_data = np.zeros(mask_array.shape)
             output_data[mask_array] = np.dot(shore_array, prediction_dir)
-        
+
         elif self.inputs.model == 'tensor':
             dti_wls = dti.TensorModel(training_gtab)
             fit_wls = dti_wls.fit(training_data, mask=mask_array)
@@ -226,7 +226,7 @@ def _run_interface(self, runtime):
 
         else:
             raise NotImplementedError('Unsupported model: ' + self.inputs.model)
-        
+
         prediction_file = op.join(
             runtime.cwd,
             "predicted_b%d_%.2f_%.2f_%.2f.nii.gz" % (
@@ -257,7 +257,7 @@ def _run_interface(self, runtime):
         model_images = quick_load_images(self.inputs.predicted_images)
         observed_images = quick_load_images(self.inputs.hmc_warped_images)
         mask_image = nb.load(self.inputs.mask_image)
-        mask = mask_image.get_data() > 1e-6
+        mask = mask_image.get_fdata() > 1e-6
         signal_vals = model_images[mask]
         b0 = signal_vals[:, 0][:, np.newaxis]
         signal_vals = signal_vals / b0