Merge pull request #344 from CoBrALab/no_BSpline2

No b spline2

rabies/analysis_pkg/utils.py

@@ -17,7 +17,7 @@ def resample_prior_maps(in_file, ref_file, transforms = [], inverses = []):
         warped_vol_fname = os.path.abspath(
             "deformed_volume" + str(x) + ".nii.gz")
         warped_volumes.append(warped_vol_fname)
-        exec_applyTransforms(transforms=transforms, inverses=inverses, input_image=volumes_list[x], ref_image=ref_file, output_image=warped_vol_fname, mask=False)
+        exec_applyTransforms(transforms=transforms, inverses=inverses, input_image=volumes_list[x], ref_image=ref_file, output_image=warped_vol_fname, interpolation='Linear')
 
     sample_volume = sitk.ReadImage(
         warped_volumes[0])

rabies/parser.py

@@ -324,7 +324,7 @@ def get_parser():
         )
     g_registration.add_argument(
         '--anat_robust_inho_cor', type=str,
-        default='apply=false,masking=false,brain_extraction=false,template_registration=SyN',
+        default='apply=false,masking=false,brain_extraction=false,keep_mask_after_extract=false,template_registration=SyN',
         help=
             "When selecting this option, inhomogeneity correction is executed twice to optimize \n"
             "outcomes. After completing an initial inhomogeneity correction step, the corrected outputs \n"
@@ -342,6 +342,10 @@ def get_parser():
             " combined masks from inhomogeneity correction. This will enhance brain edge-matching, but \n"
             " requires good quality masks. This must be selected along the 'masking' option.\n"
             " *** Specify 'true' or 'false'. \n"
+            "* keep_mask_after_extract: If using brain_extraction, use the mask to compute the registration metric \n"
+            " within the mask only. Choose to prevent stretching of the images beyond the limit of the brain mask \n"
+            " (e.g. if the moving and target images don't have the same brain coverage).\n"
+            "*** Specify 'true' or 'false'. \n"
             "* template_registration: Specify a registration script for the alignment of the \n"
             " dataset-generated unbiased template to a reference template for masking.\n"
             "*** Rigid: conducts only rigid registration.\n"
@@ -361,15 +365,15 @@ def get_parser():
         )
     g_registration.add_argument(
         '--bold_robust_inho_cor', type=str,
-        default='apply=false,masking=false,brain_extraction=false,template_registration=SyN',
+        default='apply=false,masking=false,brain_extraction=false,keep_mask_after_extract=false,template_registration=SyN',
         help=
             "Same as --anat_robust_inho_cor, but for the EPI images.\n"
             "(default: %(default)s)\n"
             "\n"
         )
     g_registration.add_argument(
         '--commonspace_reg', type=str,
-        default='masking=false,brain_extraction=false,template_registration=SyN,fast_commonspace=false',
+        default='masking=false,brain_extraction=false,keep_mask_after_extract=false,template_registration=SyN,fast_commonspace=false',
         help=
             "Specify registration options for the commonspace registration.\n"
             "* masking: Combine masks derived from the inhomogeneity correction step to support \n"
@@ -380,6 +384,10 @@ def get_parser():
             " combined masks from inhomogeneity correction. This will enhance brain edge-matching, but \n"
             " requires good quality masks. This must be selected along the 'masking' option.\n"
             "*** Specify 'true' or 'false'. \n"
+            "* keep_mask_after_extract: If using brain_extraction, use the mask to compute the registration metric \n"
+            " within the mask only. Choose to prevent stretching of the images beyond the limit of the brain mask \n"
+            " (e.g. if the moving and target images don't have the same brain coverage).\n"
+            "*** Specify 'true' or 'false'. \n"
             "* template_registration: Specify a registration script for the alignment of the \n"
             " dataset-generated unbiased template to the commonspace atlas.\n"
             "*** Rigid: conducts only rigid registration.\n"
@@ -411,7 +419,7 @@ def get_parser():
         )
     g_registration.add_argument(
         "--bold2anat_coreg", type=str, 
-        default='masking=false,brain_extraction=false,registration=SyN',
+        default='masking=false,brain_extraction=false,keep_mask_after_extract=false,registration=SyN',
         help=
             "Specify the registration script for cross-modal alignment between the EPI and structural\n"
             "images. This operation is responsible for correcting EPI susceptibility distortions.\n"
@@ -422,6 +430,10 @@ def get_parser():
             " inhomogeneity correction. This will enhance brain edge-matching, but requires good quality \n"
             " masks. This must be selected along the 'masking' option.\n"
             "*** Specify 'true' or 'false'. \n"
+            "* keep_mask_after_extract: If using brain_extraction, use the mask to compute the registration metric \n"
+            " within the mask only. Choose to prevent stretching of the images beyond the limit of the brain mask \n"
+            " (e.g. if the moving and target images don't have the same brain coverage).\n"
+            "*** Specify 'true' or 'false'. \n"
             "* registration: Specify a registration script.\n"
             "*** Rigid: conducts only rigid registration.\n"
             "*** Affine: conducts Rigid then Affine registration.\n"
@@ -467,6 +479,14 @@ def get_parser():
             "(default: %(default)s)\n"
             "\n"
         )
+    g_resampling.add_argument(
+        "--interpolation", type=str, default='Linear',
+        help=
+            "Select the interpolator which will be used by antsApplyTransforms (e.g. 'Linear' or 'BSpline[5]') \n"
+            "to resample preprocessed timeseries. \n"
+            "(default: %(default)s)\n"
+            "\n"
+        )    
 
     g_stc = preprocess.add_argument_group(
         title='STC Options', 
@@ -1052,24 +1072,34 @@ def read_parser(parser, args):
             name='bold_inho_cor')
 
         opts.commonspace_reg = parse_argument(opt=opts.commonspace_reg, 
-            key_value_pairs = {'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 
+            key_value_pairs = {'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 'keep_mask_after_extract':['true', 'false'], 
                 'template_registration':['Rigid', 'Affine', 'SyN', 'no_reg'], 'fast_commonspace':['true', 'false']},
             name='commonspace_reg')
 
         opts.bold2anat_coreg = parse_argument(opt=opts.bold2anat_coreg, 
-            key_value_pairs = {'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 
+            key_value_pairs = {'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 'keep_mask_after_extract':['true', 'false'], 
                 'registration':['Rigid', 'Affine', 'SyN', 'no_reg']},
             name='bold2anat_coreg')
 
         opts.anat_robust_inho_cor = parse_argument(opt=opts.anat_robust_inho_cor, 
-            key_value_pairs = {'apply':['true', 'false'], 'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 
+            key_value_pairs = {'apply':['true', 'false'], 'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 'keep_mask_after_extract':['true', 'false'], 
                 'template_registration':['Rigid', 'Affine', 'SyN', 'no_reg']},
             name='anat_robust_inho_cor')
 
         opts.bold_robust_inho_cor = parse_argument(opt=opts.bold_robust_inho_cor, 
-            key_value_pairs = {'apply':['true', 'false'], 'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 
+            key_value_pairs = {'apply':['true', 'false'], 'masking':['true', 'false'], 'brain_extraction':['true', 'false'], 'keep_mask_after_extract':['true', 'false'], 
                 'template_registration':['Rigid', 'Affine', 'SyN', 'no_reg']},
             name='bold_robust_inho_cor')
+
+        # check that masking/extraction options are well set
+        for name,opt in zip(['--commonspace_reg','--bold2anat_coreg','--anat_robust_inho_cor','--bold_robust_inho_cor'],
+                            [opts.commonspace_reg, opts.bold2anat_coreg, opts.anat_robust_inho_cor, opts.bold_robust_inho_cor]):
+            if opt['brain_extraction']:
+                if not opt['masking']:
+                    raise ValueError(f"For {name}, 'masking' must be set to 'true' if 'brain_extraction' is set to 'true'.")
+            if opt['keep_mask_after_extract']:
+                if not opt['brain_extraction']:
+                    raise ValueError(f"For {name}, 'brain_extraction' must be set to 'true' if 'keep_mask_after_extract' is set to 'true'.")
 
     elif opts.rabies_stage == 'confound_correction':
         opts.frame_censoring = parse_argument(opt=opts.frame_censoring, 

rabies/preprocess_pkg/commonspace_reg.py

@@ -11,7 +11,7 @@
 from .preprocess_visual_QC import PlotOverlap,template_masking
 
 
-def init_commonspace_reg_wf(opts, commonspace_masking, brain_extraction, template_reg, fast_commonspace, inherit_unbiased, output_folder, transforms_datasink, num_procs, output_datasinks, joinsource_list, name='commonspace_reg_wf'):
+def init_commonspace_reg_wf(opts, commonspace_masking, brain_extraction, keep_mask_after_extract, template_reg, fast_commonspace, inherit_unbiased, output_folder, transforms_datasink, num_procs, output_datasinks, joinsource_list, name='commonspace_reg_wf'):
     # commonspace_wf_head_start
     """
     This workflow handles the alignment of all MRI sessions to a common space. This is conducted first by generating
@@ -42,6 +42,10 @@ def init_commonspace_reg_wf(opts, commonspace_masking, brain_extraction, templat
                                 combined masks from inhomogeneity correction. This will enhance brain edge-matching, but 
                                 requires good quality masks. This should be selected along the 'masking' option.
                                 *** Specify 'true' or 'false'. 
+                                * keep_mask_after_extract: If using brain_extraction, use the mask to compute the registration metric
+                                within the mask only. Choose to prevent stretching of the images beyond the limit of the brain mask
+                                (e.g. if the moving and target images don't have the same brain coverage).
+                                *** Specify 'true' or 'false'.
                                 * template_registration: Specify a registration script for the alignment of the 
                                 dataset-generated unbiased template to the commonspace atlas.
                                 *** Rigid: conducts only rigid registration.
@@ -53,13 +57,14 @@ def init_commonspace_reg_wf(opts, commonspace_masking, brain_extraction, templat
                                 template_registration. This option can be faster, but may decrease the quality of 
                                 alignment between subjects. 
                                 *** Specify 'true' or 'false'. 
-                                (default: masking=false,brain_extraction=false,template_registration=SyN,fast_commonspace=false)
+                                (default: masking=false,brain_extraction=false,keep_mask_after_extract=false,template_registration=SyN,fast_commonspace=false)
                         
     Workflow:
         parameters
             opts: command line interface parameters
             commonspace_masking: whether masking is applied during template generation and registration
             brain_extraction: whether brain extraction is applied for template registration
+            keep_mask_after_extract: whether to keep using mask to delineate metric computation after brain extraction.
             template_reg: registration method
             fast_commonspace: whether the template generation step is skipped and instead each scan is registered directly in commonspace
             output_folder: specify a folder to execute the workflow and store important outputs
@@ -124,7 +129,7 @@ def init_commonspace_reg_wf(opts, commonspace_masking, brain_extraction, templat
                                                     'inverse_warp', 'warped_image']),
                                             name="atlas_reg_inherited")
     else:
-        atlas_reg = pe.Node(Function(input_names=['reg_method', 'brain_extraction', 'moving_image', 'moving_mask', 'fixed_image', 'fixed_mask', 'rabies_data_type'],
+        atlas_reg = pe.Node(Function(input_names=['reg_method', 'brain_extraction', 'keep_mask_after_extract', 'moving_image', 'moving_mask', 'fixed_image', 'fixed_mask', 'rabies_data_type'],
                                         output_names=['affine', 'warp',
                                                     'inverse_warp', 'warped_image'],
                                         function=run_antsRegistration),
@@ -136,6 +141,7 @@ def init_commonspace_reg_wf(opts, commonspace_masking, brain_extraction, templat
                 brain_extraction=False
 
         atlas_reg.inputs.brain_extraction = brain_extraction
+        atlas_reg.inputs.keep_mask_after_extract = keep_mask_after_extract
         atlas_reg.inputs.rabies_data_type = opts.data_type
         atlas_reg.plugin_args = {
             'qsub_args': f'-pe smp {str(3*opts.min_proc)}', 'overwrite': True}
@@ -261,7 +267,7 @@ def resample_unbiased_mask(unbiased_template, template_mask,to_atlas_affine,to_a
             transform_list=[to_atlas_affine,to_atlas_inverse_warp]
             inverse_list=[1,0]
             exec_applyTransforms(transforms = transform_list, inverses = inverse_list, 
-                input_image = template_mask, ref_image = unbiased_template, output_image = unbiased_mask, mask=True)
+                input_image = template_mask, ref_image = unbiased_template, output_image = unbiased_mask, interpolation='GenericLabel')
             return unbiased_mask
 
         resample_unbiased_mask_node = pe.Node(Function(input_names=['unbiased_template', 'template_mask','to_atlas_affine','to_atlas_inverse_warp'],
@@ -361,13 +367,14 @@ def resample_unbiased_mask(unbiased_template, template_mask,to_atlas_affine,to_a
                                                                             'unbiased_to_atlas_inverse_warp', 'warped_unbiased']),
                                                 name="inherit_unbiased_inputnode")
 
-            inherit_unbiased_reg_node = pe.Node(Function(input_names=['reg_method', 'brain_extraction', 'moving_image', 'moving_mask', 'fixed_image', 'fixed_mask', 'rabies_data_type'],
+            inherit_unbiased_reg_node = pe.Node(Function(input_names=['reg_method', 'brain_extraction', 'keep_mask_after_extract', 'moving_image', 'moving_mask', 'fixed_image', 'fixed_mask', 'rabies_data_type'],
                                             output_names=['affine', 'warp',
                                                         'inverse_warp', 'warped_image'],
                                             function=run_antsRegistration),
                                 name='inherit_unbiased_reg', mem_gb=2*opts.scale_min_memory)
             inherit_unbiased_reg_node.inputs.reg_method = 'Rigid' # this is the registration modelbuild conducts
             inherit_unbiased_reg_node.inputs.brain_extraction = False # brain extraction is not applied during template building
+            inherit_unbiased_reg_node.inputs.keep_mask_after_extract = False # brain extraction is not applied during template building
             inherit_unbiased_reg_node.inputs.rabies_data_type = opts.data_type
 
             if commonspace_masking: # this parameter should be inherited from previous run
@@ -544,7 +551,7 @@ def prep_commonspace_transform(native_ref, atlas_mask, native_to_unbiased_affine
     from rabies.utils import exec_applyTransforms
     # resample the atlas brain mask to native space
     exec_applyTransforms(transforms = commonspace_to_native_transform_list, inverses = commonspace_to_native_inverse_list, 
-        input_image = atlas_mask, ref_image = native_ref, output_image = native_mask, mask=True)
+        input_image = atlas_mask, ref_image = native_ref, output_image = native_mask, interpolation='GenericLabel')
 
     return native_mask, native_to_commonspace_transform_list,native_to_commonspace_inverse_list,commonspace_to_native_transform_list,commonspace_to_native_inverse_list
 

rabies/preprocess_pkg/hmc.py

@@ -239,7 +239,7 @@ def register_slice(fixed_image, moving_image):
     registration_method.SetMetricSamplingStrategy(registration_method.NONE)
     # registration_method.SetMetricSamplingPercentage(0.01)
 
-    registration_method.SetInterpolator(sitk.sitkBSpline)
+    registration_method.SetInterpolator(sitk.sitkLinear)
 
     # Optimizer settings.
     registration_method.SetOptimizerAsGradientDescent(
@@ -278,7 +278,7 @@ def slice_specific_registration(i, ref_file, timeseries_file):
 
         final_transform = register_slice(fixed_image, moving_image)
         moving_resampled = sitk.Resample(moving_image, fixed_image, final_transform,
-                                         sitk.sitkBSplineResamplerOrder4, 0.0, moving_image.GetPixelID())
+                                         sitk.sitkLinear, 0.0, moving_image.GetPixelID())
 
         resampled_slice = sitk.GetArrayFromImage(moving_resampled)
         volume_array[:, j, :] = resampled_slice