close to final working version

- only one atlas from CLI instead of multiples (was using atlas[0] some
places)
- remove final transforms (not used with nativesurf workflow)
- giftis only for atlas (needs updated multihist7 files - workflow to
generate these in resources/atlas-v2 for now)
- use unfoldedreg when the atlas has the metrics
   - option for dentate-based unfoldreg too (unused right now)

TODO:
 - update atlas files (use another location for testing)
 - check unfoldreg warps
 - make nativesurf the default

hippunfold/config/snakebids.yml

@@ -187,11 +187,10 @@ parse_args:
       - 'mouse'
     default: 
       - 'multihist7'
-    nargs: '+'
     help: 'Select the atlas (unfolded space) to use for subfield labels. (default: %(default)s)'
 
   --no_unfolded_reg:
-    help: 'Do not perform unfolded space (2D) registration based on thickness, curvature, and gyrification for closer alignment to the reference atlas. NOTE: only multihist7 has these features currently, so this unfolded_reg is automatically skipped if a different atlas is chosen. (default: %(default)s)'
+    help: 'Do not perform unfolded space (2D) registration based on surface metrics (e.g. thickness, curvature, and gyrification) for closer alignment to the reference atlas. NOTE: only multihist7 has these features currently, so this unfolded_reg is automatically skipped if a different atlas is chosen. (default: %(default)s)'
     default: False
     action: store_true
 
@@ -450,12 +449,21 @@ template_files:
 
 atlas_files:
   multihist7:
-    label_nii: sub-maxprob_label-hipp_desc-manualsubfieldsunfoldaligned_dseg.nii.gz
-    lut: desc-subfields_atlas-multihist7_dseg.tsv
-    label_list: labellist.txt
-    thickness: thickness.nii.gz
-    curvature: curvature.nii.gz
-    gyrification: gyrification.nii.gz
+    lut: desc-subfields_atlas-multihist7_dseg.tsv #<-- perhaps generate this from the cifti instead of having user provide it
+    surf_gii: hemi-{hemi}_space-unfold_label-{label}_{surf_type}.surf.gii
+    metric_gii: hemi-{hemi}_label-{label}_{metric}.shape.gii
+    label_gii: hemi-{hemi}_label-{label}_dseg.label.gii
+    label_wildcards:
+     - hipp
+    hemi_wildcards:
+     - L
+     - R
+    metric_wildcards:
+     - thickness
+     - gyrification
+     - curvature
+
+
   bigbrain:
     label_nii: sub-bigbrain_hemi-{hemi}_label-hipp_desc-manualsubfields_dseg.nii.gz
     lut: desc-subfields_atlas-bigbrain_dseg.tsv
@@ -536,6 +544,7 @@ resource_urls:
     synthseg_v0.1: 'zenodo.org/record/8184230/files/trained_model.3d_fullres.Task102_synsegGenDetailed.nnUNetTrainerV2.model_best.tar'
     synthseg_v0.2: 'zenodo.org/record/8184230/files/trained_model.3d_fullres.Task203_synthseg.nnUNetTrainerV2.model_best.tar'
   atlas:
+    fakeatlas: 'this.is.a/fake/url'
     multihist7: 'files.ca-1.osf.io/v1/resources/v8acf/providers/osfstorage/65395b782827451220b86dd8/?zip='
     bigbrain: 'files.ca-1.osf.io/v1/resources/v8acf/providers/osfstorage/65395b8b13d27b123094c96f/?zip='
     magdeburg: 'files.ca-1.osf.io/v1/resources/v8acf/providers/osfstorage/65395b8013d27b122f94c938/?zip='
@@ -691,7 +700,6 @@ skip_inject_template_labels: False
 force_nnunet_model: False
 t1_reg_template: False
 generate_myelin_map: False
-no_unfolded_reg: False
 root: results
 use_template_seg: False
 template_seg_smoothing_factor: 2

hippunfold/resources/atlas-v2/README.md

@@ -0,0 +1,7 @@
+# WIP - workflow for updating atlases to new v2 specs
+
+This example is run from the existing multihist atlas folder.
+
+Still needs updates for cifti structure labels etc (will add this in once we have an updated wb_command/wb_view in our deps)
+
+We currently have templates and atlases in hippunfold (selected by --template and --atlas), where the former selects the volume-based template used for volumetric registration, and the latter selects the unfolded-space subfields (and metrics for unfolded registration). I am proposing the latter become purely surface-based, so instead of providing labels and metrics in the volumetric unfolded space, we use giftis for everything, and this now becomes surface-templates instead of atlases. 

hippunfold/resources/atlas-v2/Snakefile

@@ -0,0 +1,90 @@
+
+hippunfold_dir='/localscratch/hippunfold_khanlab'
+hipp_dseg='sub-maxprob_label-hipp_desc-manualsubfieldsunfoldaligned_dseg.nii.gz'
+labellist='labellist_withdg.txt'
+
+hemis=['L','R']
+metrics=['thickness','curvature','gyrification']
+surftypes=['midthickness','inner','outer']
+
+hemi_to_structure = {"L": "CORTEX_LEFT", "R": "CORTEX_RIGHT"}
+surf_to_secondary_type = {
+    "midthickness": "MIDTHICKNESS",
+    "inner": "PIAL",
+    "outer": "GRAY_WHITE",
+}
+
+
+
+rule all:
+    input:
+        expand('hemi-{hemi}_label-{label}_dseg.label.gii',label=['hipp','dentate'],hemi=hemis),
+        expand('hemi-{hemi}_label-{label}_{metric}.shape.gii',label=['hipp'],hemi=hemis,metric=metrics),
+        expand('hemi-{hemi}_space-unfold_label-{label}_{surftype}.surf.gii',label=['hipp','dentate'],hemi=hemis,surftype=surftypes),
+
+
+
+rule hipp_dseg:
+    input:
+        hipp_dseg
+    output:
+        temp('hipp_dseg.nii.gz')
+    shell:
+        'cp {input} {output}'
+
+rule dentate_dseg: 
+    input:
+        hipp_dseg
+    params:
+        dg_label=6
+    output:
+        temp('dentate_dseg.nii.gz')
+    shell:
+        'c3d {input} -scale 0 -shift {params.dg_label} -o {output}'
+
+rule surf_gifti:
+    input:
+        surf_gii=Path(hippunfold_dir) / 'hippunfold' / 'resources' / 'unfold_template_{label}' / 'tpl-avg_space-unfold_den-unfoldiso_{surftype}.surf.gii'
+    params:
+        structure_type=lambda wildcards: hemi_to_structure[wildcards.hemi],
+        secondary_type=lambda wildcards: surf_to_secondary_type[wildcards.surftype],
+        surface_type="FLAT",
+    output:
+        surf_gii='hemi-{hemi}_space-unfold_label-{label}_{surftype}.surf.gii'
+    shell:
+        'cp {input} {output} && '
+        "wb_command -set-structure {output} {params.structure_type} -surface-type {params.surface_type}"
+        " -surface-secondary-type {params.secondary_type}"
+
+
+rule label_gifti:
+    input:
+        dseg_nii='{label}_dseg.nii.gz',
+        lut=labellist,
+        surf_gii='hemi-{hemi}_space-unfold_label-{label}_midthickness.surf.gii'
+    params:
+        structure_type=lambda wildcards: hemi_to_structure[wildcards.hemi],
+    output:
+        label_gii='hemi-{hemi}_label-{label}_dseg.label.gii'
+    shadow: 'minimal'
+    shell:
+        'wb_command -volume-label-import {input.dseg_nii} {input.lut} label_vol.nii && '
+        'wb_command -volume-label-to-surface-mapping label_vol.nii {input.surf_gii} {output.label_gii} && '
+        "wb_command -set-structure {output} {params.structure_type}"
+
+
+
+rule metric_gifti:
+    input:
+        metric_nii='{metric}.nii.gz',
+        surf_gii='hemi-{hemi}_space-unfold_label-{label}_inner.surf.gii' #existing metric niftis are on inner space
+    params:
+        structure_type=lambda wildcards: hemi_to_structure[wildcards.hemi],
+    output:
+        metric_gii='hemi-{hemi}_label-{label}_{metric}.shape.gii'
+    shadow: 'minimal'
+    shell:
+        'wb_command -volume-to-surface-mapping {input.metric_nii} {input.surf_gii} {output.metric_gii} -trilinear && '
+        "wb_command -set-structure {output} {params.structure_type}"
+
+

hippunfold/resources/atlas-v2/labellist_withdg.txt

@@ -0,0 +1,12 @@
+Subiculum
+1 0 0 255 255
+CA1
+2 133 222 255 255
+CA2
+3 0 255 170 255
+CA3
+4 255 162 0 255
+CA4
+5 255 0 0 255
+DG
+6 255 255 0 255

hippunfold/workflow/Snakefile

@@ -61,10 +61,6 @@ if config["skip_inject_template_labels"]:  # TODO do we need this?
 
 limit_to_list.add(config["modality"])
 
-# if atlas doesn't contain alignment features
-if not config["atlas"] == ["multihist7"]: 
-    config["no_unfolded_reg"] = True
-
 
 inputs = snakebids.generate_inputs(
     bids_dir=config["bids_dir"],

hippunfold/workflow/rules/common.smk

@@ -155,67 +155,6 @@ def get_final_coords():
     return coords
 
 
-def get_final_transforms():
-    xfms = []
-
-    xfms.extend(
-        inputs[config["modality"]].expand(
-            bids(
-                root=root,
-                datatype="warps",
-                **inputs.subj_wildcards,
-                label="{autotop}",
-                suffix="xfm.nii.gz",
-                hemi="{hemi}",
-                from_="{space}",
-                to="unfold",
-                mode="image",
-            ),
-            space=ref_spaces,
-            autotop=config["autotop_labels"],
-            hemi=config["hemi"],
-            allow_missing=True,
-        )
-    )
-
-    xfms.extend(
-        inputs[config["modality"]].expand(
-            bids(
-                root=root,
-                datatype="warps",
-                **inputs.subj_wildcards,
-                label="{autotop}",
-                suffix="xfm.nii.gz",
-                hemi="{hemi}",
-                from_="unfold",
-                to="{space}",
-                mode="image",
-            ),
-            space=ref_spaces,
-            autotop=config["autotop_labels"],
-            hemi=config["hemi"],
-            allow_missing=True,
-        )
-    )
-
-    xfms.extend(
-        inputs[config["modality"]].expand(
-            bids(
-                root=root,
-                datatype="warps",
-                **inputs.subj_wildcards,
-                label="{autotop}",
-                suffix="refvol.nii.gz",
-                space="unfold",
-            ),
-            autotop=config["autotop_labels"],
-            allow_missing=True,
-        )
-    )
-
-    return xfms
-
-
 def get_final_anat():
     anat = []
     if "T1w" in ref_spaces or "T2w" in ref_spaces:
@@ -335,11 +274,10 @@ def get_final_qc():
 def get_final_subj_output():
     subj_output = []
     subj_output.extend(get_final_spec())
-   # subj_output.extend(get_final_subfields())
-  #  subj_output.extend(get_final_coords())
-   # subj_output.extend(get_final_transforms())
-   # subj_output.extend(get_final_anat())
-   # subj_output.extend(get_final_qc())
+    subj_output.extend(get_final_subfields())
+    subj_output.extend(get_final_coords())
+    subj_output.extend(get_final_anat())
+    subj_output.extend(get_final_qc())
     return subj_output
 
 

hippunfold/workflow/rules/gifti.smk

@@ -214,98 +214,6 @@ rule affine_gii_to_native:
     shell:
         "wb_command -surface-apply-affine {input.gii} {input.xfm} {output.gii}"
 
-rule resample_atlas_to_refvol:
-    """this is just done in case the atlas has a different unfolded config than the current run"""
-    input:
-        refvol=bids(
-            root=root,
-            space="unfold",
-            label="hipp",
-            datatype="warps",
-            suffix="refvol.nii.gz",
-            **inputs.subj_wildcards
-        ),
-        atlas_dir=lambda wildcards: Path(download_dir) / "atlas" / wildcards.atlas,
-    params:
-        atlas=lambda wildcards, input: Path(input.atlas_dir)
-        / config["atlas_files"][wildcards.atlas]["label_nii"].format(**wildcards),
-    output:
-        label_nii=bids(
-            root=work,
-            datatype="anat",
-            suffix="subfields.nii.gz",
-            space="unfold",
-            hemi="{hemi}",
-            label="hipp",
-            atlas="{atlas}",
-            **inputs.subj_wildcards
-        ),
-    log:
-        bids(
-            root="logs",
-            suffix="resamplesubfieldrefvol",
-            space="unfold",
-            hemi="{hemi}",
-            label="hipp",
-            atlas="{atlas}",
-            **inputs.subj_wildcards
-        ),
-    container:
-        config["singularity"]["autotop"]
-    group:
-        "subj"
-    shell:
-        "antsApplyTransforms -d 3 -n MultiLabel -i {params.atlas} -r {input.refvol} -o {output.label_nii} -v &> {log}"
-
-
-rule nii_to_label_gii:
-    input:
-        label_nii=bids(
-            root=work,
-            datatype="anat",
-            suffix="subfields.nii.gz",
-            space="unfold",
-            hemi="{hemi}",
-            label="hipp",
-            atlas="{atlas}",
-            **inputs.subj_wildcards
-        ),
-        surf=os.path.join(
-            workflow.basedir,
-            "..",
-            "resources",
-            "unfold_template_hipp",
-            "tpl-avg_space-unfold_den-{density}_midthickness.surf.gii",
-        ),
-        atlas_dir=lambda wildcards: Path(download_dir) / "atlas" / wildcards.atlas,
-    params:
-        label_list=lambda wildcards, input: Path(input.atlas_dir)
-        / config["atlas_files"][wildcards.atlas]["label_list"].format(**wildcards),
-        structure_type=lambda wildcards: hemi_to_structure[wildcards.hemi],
-    output:
-        label_gii=bids(
-            root=root,
-            datatype="surf",
-            den="{density}",
-            suffix="subfields.label.gii",
-            space="{space}",
-            hemi="{hemi}",
-            label="hipp",
-            atlas="{atlas}",
-            **inputs.subj_wildcards
-        ),
-    group:
-        "subj"
-    container:
-        config["singularity"]["autotop"]
-    shadow:
-        "minimal"
-    shell:
-        "wb_command -volume-to-surface-mapping {input.label_nii} {input.surf} temp.shape.gii -enclosing  && "
-        "wb_command -metric-label-import temp.shape.gii {params.label_list} {output.label_gii} && "
-        "wb_command -set-structure {output.label_gii} {params.structure_type}"
-
-
 def get_cmd_cifti_metric(wildcards, input, output):
     cmd = f"wb_command  -cifti-create-dense-scalar {output}"
     if "L" in config["hemi"]:
@@ -505,8 +413,7 @@ rule create_spec_file_hipp:
                 **inputs.subj_wildcards
             ),
             surfname=["midthickness"],
-            space="{space}",
-            #space=["{space}", "unfold"],
+            space=["{space}", "unfold","unfoldreg"],
             allow_missing=True,
         ),
         cifti_metrics=lambda wildcards: inputs[config["modality"]].expand(