Merge pull request #418 from dPys/development

Development

Author: dPys

docs/_static/omnetome.png

@@ -403,14 +403,14 @@ def get_parser():
     parser.add_argument(
         "-em",
         metavar="Error margin",
-        default=12,
+        default=10,
         nargs="+",
         help="(Hyperparameter): Distance (in the units of the streamlines, "
              "usually mm). If any coordinate in the streamline is within this "
              "distance from the center of any voxel in the ROI, the filtering "
              "criterion is set to True for this streamline, otherwise False. "
              "Defaults to the distance between the center of each voxel and "
-             "the corner of the voxel.\n",
+             "the corner of the voxel. Default is 10.\n",
     )
     parser.add_argument(
         "-dg",
@@ -3114,9 +3114,6 @@ def wf_multi_subject(
             config.update_config(cfg_v)
             config.enable_debug_mode()
             config.enable_resource_monitor()
-
-            import logging
-
             callback_log_path = f"{wf_multi.base_dir}/run_stats.log"
             logger = logging.getLogger("callback")
             logger.setLevel(logging.DEBUG)
@@ -3171,7 +3168,7 @@ def wf_multi_subject(
                     shutil.rmtree(cnfnd_tmp_dir)
                 shutil.rmtree(f"{dir}/reg_fmri", ignore_errors=True)
                 for file_ in [i for i in glob.glob(
-                        f"{dir}/*/*") if os.path.isfile(i)]:
+                        f"{dir}/func/*") if os.path.isfile(i)]:
                     if ("reor-RAS" in file_) or ("res-" in file_):
                         try:
                             os.remove(file_)
@@ -3181,7 +3178,7 @@ def wf_multi_subject(
                 shutil.rmtree(f"{dir}/dmri_tmp", ignore_errors=True)
                 shutil.rmtree(f"{dir}/reg_dmri", ignore_errors=True)
                 for file_ in [i for i in glob.glob(
-                        f"{dir}/*/*") if os.path.isfile(i)]:
+                        f"{dir}/dwi/*") if os.path.isfile(i)]:
                     if ("reor-RAS" in file_) or ("res-" in file_):
                         try:
                             os.remove(file_)
@@ -3319,9 +3316,6 @@ def wf_multi_subject(
             config.update_config(cfg_v)
             config.enable_debug_mode()
             config.enable_resource_monitor()
-
-            import logging
-
             callback_log_path = f"{wf.base_dir}/run_stats.log"
             logger = logging.getLogger("callback")
             logger.setLevel(logging.DEBUG)
@@ -3372,15 +3366,15 @@ def wf_multi_subject(
         print("Cleaning up...")
         if func_file:
             for file_ in [i for i in glob.glob(
-                    f"{subj_dir}/*/*") if os.path.isfile(i)]:
+                    f"{subj_dir}/func/*") if os.path.isfile(i)]:
                 if ("reor-RAS" in file_) or ("res-" in file_):
                     try:
                         os.remove(file_)
                     except BaseException:
                         continue
         if dwi_file:
             for file_ in [i for i in glob.glob(
-                    f"{subj_dir}/*/*") if os.path.isfile(i)]:
+                    f"{subj_dir}/dwi/*") if os.path.isfile(i)]:
                 if ("reor-RAS" in file_) or ("res-" in file_) or \
                    ("_bvecs_reor.bvec" in file_):
                     try:

pynets/cli/pynets_run.py

@@ -498,11 +498,13 @@ def _run_interface(self, runtime):
         from pynets.fmri import clustools
         from pynets.registration.reg_utils import check_orient_and_dims
         from joblib import Parallel, delayed
+        from joblib.externals.loky.backend import resource_tracker
         from pynets.registration import reg_utils as regutils
         from pynets.core.utils import decompress_nifti
         import pkg_resources
         import shutil
         import tempfile
+        resource_tracker.warnings = None
 
         template = pkg_resources.resource_filename(
             "pynets", f"templates/{self.inputs.template_name}_brain_"
@@ -625,6 +627,7 @@ def _run_interface(self, runtime):
             if float(c_boot) > 1:
                 import random
                 from joblib import Memory
+                from joblib.externals.loky import get_reusable_executor
                 print(
                     f"Performing circular block bootstrapping with {c_boot}"
                     f" iterations..."
@@ -701,6 +704,8 @@ def run_bs_iteration(i, ts_data, work_dir, local_corr,
                                                    iter_bootedparcels if
                                                    i is not None])
                         counter = len(boot_parcellations)
+                        del iter_bootedparcels
+                        gc.collect()
 
                 print('Bootstrapped samples complete:')
                 print(boot_parcellations)
@@ -713,6 +718,10 @@ def run_bs_iteration(i, ts_data, work_dir, local_corr,
                 nib.save(consensus_parcellation, nip.uatlas)
                 memory.clear(warn=False)
                 shutil.rmtree(cache_dir, ignore_errors=True)
+                del parallel, memory, cache_dir
+                get_reusable_executor().shutdown(wait=True)
+                gc.collect()
+
                 for i in boot_parcellations:
                     if os.path.isfile(i):
                         os.remove(i)
@@ -1186,7 +1195,8 @@ def _run_interface(self, runtime):
                 copy=True,
                 use_hardlink=False)
         else:
-            if len(anat_mask_existing) > 0 and self.inputs.mask is None:
+            if len(anat_mask_existing) > 0 and \
+                 self.inputs.mask is None and op.isfile(anat_mask_existing[0]):
                 mask_tmp_path = fname_presuffix(
                     anat_mask_existing[0], suffix="_tmp", newpath=runtime.cwd
                 )
@@ -1195,9 +1205,6 @@ def _run_interface(self, runtime):
                     mask_tmp_path,
                     copy=True,
                     use_hardlink=False)
-                mask_tmp_path = regutils.check_orient_and_dims(
-                    mask_tmp_path, runtime.cwd,
-                    self.inputs.vox_size)
             else:
                 mask_tmp_path = None
 
@@ -1241,12 +1248,21 @@ def _run_interface(self, runtime):
         else:
             csf_mask = None
 
+        anat_file_tmp_path = fname_presuffix(
+            self.inputs.anat_file, suffix="_tmp", newpath=runtime.cwd
+        )
+        copyfile(
+            self.inputs.anat_file,
+            anat_file_tmp_path,
+            copy=True,
+            use_hardlink=False)
+
         reg = register.DmriReg(
             basedir_path=runtime.cwd,
             fa_path=fa_tmp_path,
             ap_path=ap_tmp_path,
             B0_mask=B0_mask_tmp_path,
-            anat_file=self.inputs.anat_file,
+            anat_file=anat_file_tmp_path,
             vox_size=self.inputs.vox_size,
             template_name=self.inputs.template_name,
             simple=self.inputs.simple,
@@ -1284,7 +1300,7 @@ def _run_interface(self, runtime):
         self._results["gm_in_dwi"] = reg.gm_in_dwi
         self._results["vent_csf_in_dwi"] = reg.vent_csf_in_dwi
         self._results["csf_mask_dwi"] = reg.csf_mask_dwi
-        self._results["anat_file"] = self.inputs.anat_file
+        self._results["anat_file"] = anat_file_tmp_path
         self._results["t1w2dwi"] = reg.t1w2dwi
         self._results["B0_mask"] = B0_mask_tmp_path
         self._results["ap_path"] = ap_tmp_path
@@ -1894,6 +1910,7 @@ class RegisterFunc(SimpleInterface):
     def _run_interface(self, runtime):
         import gc
         import glob
+        import os.path as op
         from pynets.registration import register
         from pynets.registration import reg_utils as regutils
         from nipype.utils.filemanip import fname_presuffix, copyfile
@@ -1917,7 +1934,8 @@ def _run_interface(self, runtime):
                 copy=True,
                 use_hardlink=False)
         else:
-            if len(anat_mask_existing) > 0 and self.inputs.mask is None:
+            if len(anat_mask_existing) > 0 and \
+                 self.inputs.mask is None and op.isfile(anat_mask_existing[0]):
                 mask_tmp_path = fname_presuffix(
                     anat_mask_existing[0], suffix="_tmp",
                     newpath=runtime.cwd
@@ -1927,9 +1945,6 @@ def _run_interface(self, runtime):
                     mask_tmp_path,
                     copy=True,
                     use_hardlink=False)
-                mask_tmp_path = regutils.check_orient_and_dims(
-                    mask_tmp_path, self.inputs.in_dir,
-                    self.inputs.vox_size)
             else:
                 mask_tmp_path = None
 
@@ -1961,9 +1976,18 @@ def _run_interface(self, runtime):
         else:
             wm_mask = None
 
+        anat_file_tmp_path = fname_presuffix(
+            self.inputs.anat_file, suffix="_tmp", newpath=runtime.cwd
+        )
+        copyfile(
+            self.inputs.anat_file,
+            anat_file_tmp_path,
+            copy=True,
+            use_hardlink=False)
+
         reg = register.FmriReg(
             basedir_path=runtime.cwd,
-            anat_file=self.inputs.anat_file,
+            anat_file=anat_file_tmp_path,
             vox_size=self.inputs.vox_size,
             template_name=self.inputs.template_name,
             simple=self.inputs.simple,
@@ -3063,7 +3087,7 @@ def _run_interface(self, runtime):
         from dipy.io import save_pickle
         from dipy.io import read_bvals_bvecs
         from dipy.core.gradients import gradient_table
-
+        from nipype.utils.filemanip import copyfile, fname_presuffix
         # from dipy.segment.mask import median_otsu
         from pynets.registration.reg_utils import median
         from pynets.dmri.dmri_utils import normalize_gradients, extract_b0
@@ -3075,8 +3099,26 @@ def _run_interface(self, runtime):
         gtab_file = f"{runtime.cwd}/gtab.pkl"
         all_b0s_file = f"{runtime.cwd}/all_b0s.nii.gz"
 
+        fbval_tmp_path = fname_presuffix(
+            self.inputs.fbval, suffix="_tmp", newpath=runtime.cwd
+        )
+        copyfile(
+            self.inputs.fbval,
+            fbval_tmp_path,
+            copy=True,
+            use_hardlink=False)
+
+        fbvec_tmp_path = fname_presuffix(
+            self.inputs.fbvec, suffix="_tmp", newpath=runtime.cwd
+        )
+        copyfile(
+            self.inputs.fbvec,
+            fbvec_tmp_path,
+            copy=True,
+            use_hardlink=False)
+
         # loading bvecs/bvals
-        bvals, bvecs = read_bvals_bvecs(self.inputs.fbval, self.inputs.fbvec)
+        bvals, bvecs = read_bvals_bvecs(fbval_tmp_path, fbvec_tmp_path)
         bvecs_norm, bvals_norm = normalize_gradients(
             bvecs, bvals, b0_threshold=self.inputs.b0_thr
         )
@@ -3103,10 +3145,19 @@ def _run_interface(self, runtime):
         # Save gradient table to pickle
         save_pickle(gtab_file, gtab)
 
+        dwi_file_tmp_path = fname_presuffix(
+            self.inputs.dwi_file, suffix="_tmp", newpath=runtime.cwd
+        )
+        copyfile(
+            self.inputs.dwi_file,
+            dwi_file_tmp_path,
+            copy=True,
+            use_hardlink=False)
+
         # Extract and Combine all b0s collected, make mean b0
         print("Extracting b0's...")
         all_b0s_file = extract_b0(
-            self.inputs.dwi_file,
+            dwi_file_tmp_path,
             b0_thr_ixs,
             all_b0s_file)
         med_b0_file = median(all_b0s_file)
@@ -3130,11 +3181,16 @@ def _run_interface(self, runtime):
         # Get mean B0 brain mask
         cmd = f"bet {med_b0_file} {B0_bet} -m -f 0.2"
         os.system(cmd)
-        time.sleep(1)
+        time.sleep(2)
 
         self._results["gtab_file"] = gtab_file
         self._results["B0_bet"] = B0_bet
         self._results["B0_mask"] = B0_mask
         self._results["dwi_file"] = self.inputs.dwi_file
 
+        tmp_files = [fbval_tmp_path, fbvec_tmp_path, dwi_file_tmp_path]
+        for j in tmp_files:
+            if j is not None:
+                os.remove(j)
+
         return runtime

pynets/core/interfaces.py

@@ -1569,9 +1569,8 @@ def save_ts_to_file(
                                               "%s" % ("%s%s%s" % ("hpass-",
                                                                   hpass,
                                                                   "Hz_") if hpass is not None else ""),
-                                              "%s" % ("%s%s%s" % ("extract-",
-                                                                  extract_strategy,
-                                                                  "_") if extract_strategy is not None else ""),
+                                              "%s" % ("%s%s" % ("extract-",
+                                                                  extract_strategy) if extract_strategy is not None else ""),
                                               ".npy",
                                               )
 

pynets/core/utils.py

@@ -1392,7 +1392,8 @@ def dmri_connectometry(
     base_dirname = f"dmri_connectometry_{ID}"
     dmri_connectometry_wf = pe.Workflow(name=base_dirname)
 
-    if template_name == "MNI152_T1" or template_name == "colin27" or template_name == "CN200":
+    if template_name == "MNI152_T1" or template_name == "colin27" or \
+        template_name == "CN200":
         template = pkg_resources.resource_filename(
             "pynets", f"templates/{template_name}_brain_{vox_size}.nii.gz"
         )

pynets/core/workflows.py

@@ -161,7 +161,7 @@ def evaluate_streamline_plausibility(dwi_data, gtab, mask_data, streamlines,
     import dipy.tracking.life as life
     import dipy.core.optimize as opt
     from dipy.tracking._utils import _mapping_to_voxel
-    from dipy.data import get_sphere
+    # from dipy.data import get_sphere
     from dipy.tracking import utils
     from dipy.tracking.streamline import Streamlines
 
@@ -202,11 +202,14 @@ def evaluate_streamline_plausibility(dwi_data, gtab, mask_data, streamlines,
                     dwi_data.shape).astype('bool') * dwi_data)
     # ! Remember this 4d masking function !
 
-    sphere = get_sphere(sphere)
     fiber_model = life.FiberModel(gtab)
     fiber_fit = fiber_model.fit(data_in_mask, streamlines_in_brain,
                                 affine=affine,
-                                sphere=sphere)
+                                sphere=False)
+    # sphere = get_sphere(sphere)
+    # fiber_fit = fiber_model.fit(data_in_mask, streamlines_in_brain,
+    #                             affine=affine,
+    #                             sphere=sphere)
     streamlines = list(np.array(streamlines_in_brain)[
                            np.where(fiber_fit.beta > 0)[0]])
     pruned_count = len(streamlines)
@@ -220,6 +223,8 @@ def evaluate_streamline_plausibility(dwi_data, gtab, mask_data, streamlines,
                           'LiFE model. Is your acquisition hemispheric? '
                           'Also check the gradient table for errors. \n'))
         return streamlines_in_brain
+    else:
+        del streamlines_in_brain
 
     model_predict = fiber_fit.predict()
     model_error = model_predict - fiber_fit.data