Add --fs-no-resume option to reuse existing freesurfer outputs without resuming (eg. longitudinal pipeline base)

.circleci/config.yml

@@ -241,7 +241,7 @@ jobs:
           command: |
             mkdir -p /tmp/fslicense
             cd /tmp/fslicense
-            echo "cHJpbnRmICJrcnp5c3p0b2YuZ29yZ29sZXdza2lAZ21haWwuY29tXG41MTcyXG4gKkN2dW12RVYzelRmZ1xuRlM1Si8yYzFhZ2c0RVxuIiA+IGxpY2Vuc2UudHh0Cg==" | base64 -d | sh
+              echo "cHJpbnRmICJrcnp5c3p0b2YuZ29yZ29sZXdza2lAZ21haWwuY29tXG41MTcyXG4gKkN2dW12RVYzelRmZ1xuRlM1Si8yYzFhZ2c0RVxuIiA+IGxpY2Vuc2UudHh0Cg==" | base64 -d | sh
       - persist_to_workspace:
           root: /tmp
           paths:
@@ -544,6 +544,24 @@ jobs:
       - store_artifacts:
           path: /tmp/ds005/derivatives
           destination: fasttrack
+      - run:
+          name: Check fs-reuse-base using existing freesufer output
+          no_output_timeout: 5m
+          command: |
+            bash /tmp/src/smriprep/.circleci/ds005_run.sh --fs-reuse-base
+      - run:
+          name: Clean working directory
+          when: on_fail
+          command: |
+            rm -rf /tmp/ds005/work/smriprep_wf/fsdir_run_*/
+            find /tmp/ds005/work \( -name "*.nii.gz" -or -name "*.nii" -or -name "*.gii" -or -name "*.h5" \) \
+                -exec sh -c 'rm -f {}; touch {}' \;
+      - store_artifacts:
+          path: /tmp/ds005/work
+          destination: fs_reuse_base
+      - store_artifacts:
+          path: /tmp/ds005/derivatives
+          destination: fs_reuse_base
   ds054:
     <<: *machine_defaults
     environment:

smriprep/cli/run.py

@@ -214,6 +214,12 @@ def get_parser():
         help='Path to existing FreeSurfer subjects directory to reuse. '
         '(default: OUTPUT_DIR/freesurfer)',
     )
+    g_fs.add_argument(
+        '--fs-reuse-base',
+        action='store_true',
+        dest='fs_reuse_base',
+        help='Reuse freesurfer base template (from longitudinal preprocessing)',
+    )
     g_fs.add_argument(
         '--cifti-output',
         nargs='?',
@@ -607,6 +613,7 @@ def build_workflow(opts, retval):
         freesurfer=opts.run_reconall,
         fs_subjects_dir=opts.fs_subjects_dir,
         hires=opts.hires,
+        fs_reuse_base=opts.fs_reuse_base,
         layout=layout,
         longitudinal=opts.longitudinal,
         low_mem=opts.low_mem,

smriprep/workflows/anatomical.py

@@ -111,6 +111,7 @@ def init_anat_preproc_wf(
     cifti_output: ty.Literal['91k', '170k', False] = False,
     name: str = 'anat_preproc_wf',
     skull_strip_fixed_seed: bool = False,
+    fs_reuse_base: bool = False,
 ):
     """
     Stage the anatomical preprocessing steps of *sMRIPrep*.
@@ -183,6 +184,10 @@ def init_anat_preproc_wf(
         Do not use a random seed for skull-stripping - will ensure
         run-to-run replicability when used with --omp-nthreads 1
         (default: ``False``).
+    fs_reuse_base : bool
+        Adjust pipeline to reuse base template
+        of an existing longitudinal freesurfer output
+        (default: ``False``).
 
     Inputs
     ------
@@ -275,6 +280,7 @@ def init_anat_preproc_wf(
         sloppy=sloppy,
         omp_nthreads=omp_nthreads,
         skull_strip_fixed_seed=skull_strip_fixed_seed,
+        fs_reuse_base=fs_reuse_base,
     )
     template_iterator_wf = init_template_iterator_wf(spaces=spaces)
     ds_std_volumes_wf = init_ds_anat_volumes_wf(
@@ -458,6 +464,7 @@ def init_anat_fit_wf(
     sloppy: bool = False,
     name='anat_fit_wf',
     skull_strip_fixed_seed: bool = False,
+    fs_reuse_base: bool = False,
 ):
     """
     Stage the anatomical preprocessing steps of *sMRIPrep*.
@@ -1040,6 +1047,7 @@ def init_anat_fit_wf(
         name='surface_recon_wf',
         omp_nthreads=omp_nthreads,
         hires=hires,
+        fs_reuse_base=fs_reuse_base,
         precomputed=precomputed,
     )
     if t2w or flair:

smriprep/workflows/base.py

@@ -46,6 +46,7 @@ def init_smriprep_wf(
     freesurfer,
     fs_subjects_dir,
     hires,
+    fs_reuse_base,
     layout,
     longitudinal,
     low_mem,
@@ -88,6 +89,7 @@ def init_smriprep_wf(
                 freesurfer=True,
                 fs_subjects_dir=None,
                 hires=True,
+                fs_reuse_base=False,
                 layout=BIDSLayout('.'),
                 longitudinal=False,
                 low_mem=False,
@@ -178,6 +180,7 @@ def init_smriprep_wf(
             freesurfer=freesurfer,
             derivatives=derivatives,
             hires=hires,
+            fs_reuse_base=fs_reuse_base,
             layout=layout,
             longitudinal=longitudinal,
             low_mem=low_mem,
@@ -214,6 +217,7 @@ def init_single_subject_wf(
     derivatives,
     freesurfer,
     hires,
+    fs_reuse_base,
     layout,
     longitudinal,
     low_mem,
@@ -258,6 +262,7 @@ def init_single_subject_wf(
                 freesurfer=True,
                 derivatives=[],
                 hires=True,
+                fs_reuse_base=False,
                 layout=BIDSLayout('.'),
                 longitudinal=False,
                 low_mem=False,
@@ -286,6 +291,9 @@ def init_single_subject_wf(
         Enable FreeSurfer surface reconstruction (may increase runtime)
     hires : :obj:`bool`
         Enable sub-millimeter preprocessing in FreeSurfer
+    fs_reuse_base : bool
+        Adjust pipeline to reuse base template
+        of an existing longitudinal freesurfer output
     layout : BIDSLayout object
         BIDS dataset layout
     longitudinal : :obj:`bool`
@@ -418,6 +426,7 @@ def init_single_subject_wf(
         precomputed=deriv_cache,
         freesurfer=freesurfer,
         hires=hires,
+        fs_reuse_base=fs_reuse_base,
         longitudinal=longitudinal,
         msm_sulc=msm_sulc,
         name='anat_preproc_wf',

smriprep/workflows/surfaces.py

@@ -64,6 +64,7 @@ def init_surface_recon_wf(
     *,
     omp_nthreads: int,
     hires: bool,
+    fs_reuse_base: bool,
     precomputed: dict,
     name='surface_recon_wf',
 ):
@@ -129,14 +130,21 @@ def init_surface_recon_wf(
             :simple_form: yes
 
             from smriprep.workflows.surfaces import init_surface_recon_wf
-            wf = init_surface_recon_wf(omp_nthreads=1, hires=True, precomputed={})
+            wf = init_surface_recon_wf(
+                omp_nthreads=1,
+                hires=True,
+                fs_reuse_base=False,
+                precomputed={})
 
     Parameters
     ----------
     omp_nthreads : int
         Maximum number of threads an individual process may use
     hires : bool
         Enable sub-millimeter preprocessing in FreeSurfer
+    fs_reuse_base : bool
+        Adjust pipeline to reuse base template
+        of an existing longitudinal freesurfer output
 
     Inputs
     ------
@@ -238,39 +246,82 @@ def init_surface_recon_wf(
         name='sync',
     )
 
+    if not fs_reuse_base:
+        recon_config = pe.Node(FSDetectInputs(hires_enabled=hires), name='recon_config')
+
+        fov_check = pe.Node(niu.Function(function=_check_cw256), name='fov_check')
+        fov_check.inputs.default_flags = ['-noskullstrip', '-noT2pial', '-noFLAIRpial']
+
+        autorecon1 = pe.Node(
+            ReconAll(directive='autorecon1', openmp=omp_nthreads),
+            name='autorecon1',
+            n_procs=omp_nthreads,
+            mem_gb=5,
+        )
+        autorecon1.interface._can_resume = False
+        autorecon1.interface._always_run = True
+
+        skull_strip_extern = pe.Node(FSInjectBrainExtracted(), name='skull_strip_extern')
+
+        autorecon_resume_wf = init_autorecon_resume_wf(omp_nthreads=omp_nthreads)
+
+        # fmt:off
+        workflow.connect([
+            # Configuration
+            (inputnode, recon_config, [('t1w', 't1w_list'),
+                                       ('t2w', 't2w_list'),
+                                       ('flair', 'flair_list')]),
+            # Passing subjects_dir / subject_id enforces serial order
+            (inputnode, autorecon1, [('subjects_dir', 'subjects_dir'),
+                                     ('subject_id', 'subject_id')]),
+            (autorecon1, skull_strip_extern, [('subjects_dir', 'subjects_dir'),
+                                              ('subject_id', 'subject_id')]),
+            (skull_strip_extern, autorecon_resume_wf, [('subjects_dir', 'inputnode.subjects_dir'),
+                                                       ('subject_id', 'inputnode.subject_id')]),
+            # Reconstruction phases
+            (inputnode, autorecon1, [('t1w', 'T1_files')]),
+            (inputnode, fov_check, [('t1w', 'in_files')]),
+            (fov_check, autorecon1, [('out', 'flags')]),
+            (recon_config, autorecon1, [('t2w', 'T2_file'),
+                                        ('flair', 'FLAIR_file'),
+                                        ('hires', 'hires'),
+                                        # First run only (recon-all saves expert options)
+                                        ('mris_inflate', 'mris_inflate')]),
+            (inputnode, skull_strip_extern, [('skullstripped_t1', 'in_brain')]),
+            (recon_config, autorecon_resume_wf, [('use_t2w', 'inputnode.use_T2'),
+                                                 ('use_flair', 'inputnode.use_FLAIR')]),
+            # Generate mid-thickness surfaces
+            (autorecon_resume_wf, get_surfaces, [
+                ('outputnode.subjects_dir', 'subjects_dir'),
+                ('outputnode.subject_id', 'subject_id'),
+            ]),
+            (autorecon_resume_wf, save_midthickness, [
+                ('outputnode.subjects_dir', 'base_directory'),
+                ('outputnode.subject_id', 'container'),
+            ]),
+        ])
+        # fmt:on
+    else:
+        fs_base_inputs = autorecon1 = pe.Node(nio.FreeSurferSource(), name='fs_base_inputs')
+
+        # fmt:off
+        workflow.connect([
+            (inputnode, fs_base_inputs, [('subjects_dir', 'subjects_dir'),
+                                         ('subject_id', 'subject_id')]),
+            # Generate mid-thickness surfaces
+            (inputnode, get_surfaces, [
+                ('subjects_dir', 'subjects_dir'),
+                ('subject_id', 'subject_id'),
+            ]),
+            (inputnode, save_midthickness, [
+                ('subjects_dir', 'base_directory'),
+                ('subject_id', 'container'),
+            ]),
+        ])
+        # fmt:on
+
+    # fmt:off
     workflow.connect([
-        # Configuration
-        (inputnode, recon_config, [('t1w', 't1w_list'),
-                                   ('t2w', 't2w_list'),
-                                   ('flair', 'flair_list')]),
-        # Passing subjects_dir / subject_id enforces serial order
-        (inputnode, autorecon1, [('subjects_dir', 'subjects_dir'),
-                                 ('subject_id', 'subject_id')]),
-        (autorecon1, skull_strip_extern, [('subjects_dir', 'subjects_dir'),
-                                          ('subject_id', 'subject_id')]),
-        (skull_strip_extern, autorecon_resume_wf, [('subjects_dir', 'inputnode.subjects_dir'),
-                                                   ('subject_id', 'inputnode.subject_id')]),
-        # Reconstruction phases
-        (inputnode, autorecon1, [('t1w', 'T1_files')]),
-        (inputnode, fov_check, [('t1w', 'in_files')]),
-        (fov_check, autorecon1, [('out', 'flags')]),
-        (recon_config, autorecon1, [('t2w', 'T2_file'),
-                                    ('flair', 'FLAIR_file'),
-                                    ('hires', 'hires'),
-                                    # First run only (recon-all saves expert options)
-                                    ('mris_inflate', 'mris_inflate')]),
-        (inputnode, skull_strip_extern, [('skullstripped_t1', 'in_brain')]),
-        (recon_config, autorecon_resume_wf, [('use_t2w', 'inputnode.use_T2'),
-                                             ('use_flair', 'inputnode.use_FLAIR')]),
-        # Generate mid-thickness surfaces
-        (autorecon_resume_wf, get_surfaces, [
-            ('outputnode.subjects_dir', 'subjects_dir'),
-            ('outputnode.subject_id', 'subject_id'),
-        ]),
-        (autorecon_resume_wf, save_midthickness, [
-            ('outputnode.subjects_dir', 'base_directory'),
-            ('outputnode.subject_id', 'container'),
-        ]),
         (get_surfaces, midthickness, [
             ('white', 'in_file'),
             ('graymid', 'graymid'),