Skip to content

Commit

Permalink
Merge pull request #143 from khanlab/warps-rmfullgrid
Browse files Browse the repository at this point in the history
- adds rules for composite warps from unfold and remove map to full grid
- added a rule to create the ref space for warp files in native space
(similarly to what we already do for unfold space)
- moved ref vols to the seg_{modality} dir (were improperly sitting in the root
folder)
  • Loading branch information
akhanf authored Dec 26, 2021
2 parents 68fce70 + 42537ec commit 1d2924c
Show file tree
Hide file tree
Showing 7 changed files with 126 additions and 266 deletions.
1 change: 1 addition & 0 deletions hippunfold/workflow/Snakefile
Original file line number Diff line number Diff line change
Expand Up @@ -67,6 +67,7 @@ if 'hippb500' in config['modality']:
include: 'rules/preproc_hippb500.smk'

include: 'rules/autotop.smk'
include: 'rules/warps.smk'
include: 'rules/shape_inject.smk'
include: 'rules/gifti.smk'
include: 'rules/subfields.smk'
Expand Down
102 changes: 3 additions & 99 deletions hippunfold/workflow/rules/autotop.smk
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@ rule create_unfold_ref:
origin = 'x'.join(config['unfold_vol_ref']['origin']),
orient = config['unfold_vol_ref']['orient']
output:
nii = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards'])
nii = bids(root='work',datatype='seg_{modality}',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards'])
group: 'subj'
container: config['singularity']['autotop']
shell:
Expand All @@ -86,9 +86,9 @@ rule create_unfold_ref:
#this was unfold_phys_coords.nii in matlab implementation
rule create_unfold_coord_map:
input:
nii = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards'])
nii = bids(root='work',datatype='seg_{modality}',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards'])
output:
nii = bids(root='work',space='unfold',suffix='refcoords.nii.gz',**config['subj_wildcards'])
nii = bids(root='work',datatype='seg_{modality}',space='unfold',suffix='refcoords.nii.gz',**config['subj_wildcards'])
group: 'subj'
container: config['singularity']['autotop']
shell:
Expand All @@ -102,99 +102,3 @@ def get_laminar_coords(wildcards):
coords_io = bids(root='work',datatype='seg_{modality}',dir='IO',suffix='coords.nii.gz',desc='equivol',space='corobl',hemi='{hemi}', **config['subj_wildcards'])
return coords_io

rule create_warps:
input:
unfold_ref_nii = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards']),
unfold_phys_coords_nii = bids(root='work',space='unfold',suffix='refcoords.nii.gz',**config['subj_wildcards']),
coords_ap = bids(root='work',datatype='seg_{modality}',dir='AP',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='{hemi}', **config['subj_wildcards']),
coords_pd = bids(root='work',datatype='seg_{modality}',dir='PD',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='{hemi}', **config['subj_wildcards']),
coords_io = get_laminar_coords
params:
interp_method = 'linear'
resources:
mem_mb = 16000
output:
warp_unfold2native = bids(root='work',**config['subj_wildcards'],suffix='autotop/Warp_unfold2native.nii',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
warp_native2unfold= bids(root='work',**config['subj_wildcards'],suffix='autotop/Warp_native2unfold.nii',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
warpitk_unfold2native = bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_unfold2native.nii',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
warpitk_native2unfold= bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_native2unfold.nii',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
#warp_unfold2native_extrap = bids(root='work',**config['subj_wildcards'],suffix='autotop/Warp_unfold2native_extrapolateNearest.nii',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
group: 'subj'
log: bids(root='logs',**config['subj_wildcards'],hemi='{hemi,Lflip|R}',modality='seg-{modality}',suffix='create_warps.txt')
script: '../scripts/create_warps.py'

# TODO: add this extrapolation of the warp file, extrapolate_warp_unfold2native.m
# .. for now just using original warp
#rule extrapolate_warp_unfold2nii:


#full-grid correction of unfolded space
rule map_to_full_grid:
input:
coords_ap = bids(root='work',datatype='seg_{modality}',dir='AP',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='{hemi}', **config['subj_wildcards']),
coords_pd = bids(root='work',datatype='seg_{modality}',dir='PD',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='{hemi}', **config['subj_wildcards']),
coords_io = get_laminar_coords,
warpitk_native2unfold= bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_native2unfold.nii',desc='cropped',space='corobl',hemi='{hemi}',modality='{modality}'),
unfold_ref = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards'])
params:
dims = config['unfold_vol_ref']['dims'],
script = os.path.join(config['snakemake_dir'],'workflow','scripts','mapUnfoldToFullGrid.sh'),
warp_dir = bids(root='work',**config['subj_wildcards'],suffix='autotop',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}')
output:
warp_unfoldtemplate2unfold = bids(root='work',**config['subj_wildcards'],suffix='autotop/Warp_unfoldtemplate2unfold.nii',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
warp_unfold2unfoldtemplate = bids(root='work',**config['subj_wildcards'],suffix='autotop/Warp_unfold2unfoldtemplate.nii',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
warpitk_unfoldtemplate2unfold = bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_unfold2unfoldtemplate_0InverseWarp.nii.gz',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
warpitk_unfold2unfoldtemplate = bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_unfold2unfoldtemplate_0Warp.nii.gz',desc='cropped',space='corobl',hemi='{hemi,Lflip|R}',modality='{modality}'),
container: config['singularity']['autotop']
group: 'subj'
threads: 8
resources:
time = 15 #15min
log: bids(root='logs',**config['subj_wildcards'],space='corobl',hemi='{hemi,Lflip|R}',modality='seg{modality}',suffix='mapUnfoldToFullGrid.txt')
shell:
'SINGULARITYENV_ITK_GLOBAL_DEFAULT_NUMBER_OF_THREADS={threads} '
'{params.script} {input.coords_ap} {input.coords_pd} {input.coords_io} {input.unfold_ref} {params.warp_dir} {params.dims} &> {log}'


rule compose_warps_corobl2unfold_rhemi:
""" Compose corobl to unfold (unfold-template), for right hemi (ie no flip)"""
input:
native2unfold = bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_native2unfold.nii',desc='cropped',space='corobl',hemi='R',modality='{modality}'),
unfold2unfoldtemplate = bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_unfold2unfoldtemplate_0Warp.nii.gz',desc='cropped',space='corobl',hemi='R',modality='{modality}'),
ref = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards'])
output:
corobl2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='R',from_='corobl',to='unfold',mode='image')
container: config['singularity']['autotop']
group: 'subj'
shell: 'ComposeMultiTransform 3 {output.corobl2unfold} -R {input.ref} {input.unfold2unfoldtemplate} {input.native2unfold}'


rule compose_warps_corobl2unfold_lhemi:
""" Compose corobl to unfold (unfold-template), for left hemi (ie with flip)"""
input:
native2unfold = bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_native2unfold.nii',desc='cropped',space='corobl',hemi='Lflip',modality='{modality}'),
unfold2unfoldtemplate = bids(root='work',**config['subj_wildcards'],suffix='autotop/WarpITK_unfold2unfoldtemplate_0Warp.nii.gz',desc='cropped',space='corobl',hemi='Lflip',modality='{modality}'),
ref = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards']),
flipLR_xfm = os.path.join(config['snakemake_dir'],'resources','desc-flipLR_type-itk_xfm.txt')
output:
corobl2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='L',from_='corobl',to='unfold',mode='image')
container: config['singularity']['autotop']
group: 'subj'
shell: 'ComposeMultiTransform 3 {output.corobl2unfold} -R {input.ref} {input.unfold2unfoldtemplate} {input.native2unfold} {input.flipLR_xfm}'



rule compose_warps_t1_to_unfold:
""" Compose warps from T1w to unfold """
input:
corobl2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi}',from_='corobl',to='unfold',mode='image'),
ref = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards']),
t1w2corobl = bids(root='work',datatype='anat',**config['subj_wildcards'],suffix='xfm.txt',from_='T1w',to='corobl',desc='affine',type_='itk'),
output:
bids(root='results',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi}',from_='T1w',to='unfold',mode='image')
container: config['singularity']['autotop']
group: 'subj'
shell: 'ComposeMultiTransform 3 {output} -R {input.ref} {input.corobl2unfold} {input.t1w2corobl}'



22 changes: 20 additions & 2 deletions hippunfold/workflow/rules/common.smk
Original file line number Diff line number Diff line change
Expand Up @@ -93,7 +93,9 @@ def get_final_transforms():
else:
output_ref = 'corobl'

return expand(
xfms = []

xfms.extend(expand(
bids(
root='results',
datatype='seg_{modality}',
Expand All @@ -105,7 +107,23 @@ def get_final_transforms():
mode='image'),
space=output_ref,
hemi=config['hemi'],
allow_missing=True)
allow_missing=True))

xfms.extend(expand(
bids(
root='results',
datatype='seg_{modality}',
**config['subj_wildcards'],
suffix='xfm.nii.gz',
hemi='{hemi}',
from_='unfold',
to='{space}',
mode='image'),
space=output_ref,
hemi=config['hemi'],
allow_missing=True))

return xfms


def get_final_anat():
Expand Down
13 changes: 5 additions & 8 deletions hippunfold/workflow/rules/gifti.smk
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,9 @@ hemi_to_structure = {'L': 'CORTEX_LEFT', 'Lflip': 'CORTEX_LEFT', 'R': 'CORTEX_RI
surf_to_secondary_type = {'midthickness': 'MIDTHICKNESS', 'inner': 'PIAL', 'outer': 'GRAY_WHITE'}


rule warp_gii_unfoldtemplate2unfold:
"""warp from template space to subj unfolded"""
rule cp_template_to_unfold:
"""cp template unfold surf to subject"""
input:
warp = bids(root='work',**config['subj_wildcards'],suffix='autotop/Warp_unfoldtemplate2unfold.nii',desc='cropped',space='corobl',hemi='{hemi}',modality='{modality}'),
gii = os.path.join(config['snakemake_dir'],'resources','unfold_template','tpl-avg_space-unfold_den-{density}_{surfname}.surf.gii')
params:
structure_type = lambda wildcards: hemi_to_structure[wildcards.hemi],
Expand All @@ -18,7 +17,7 @@ rule warp_gii_unfoldtemplate2unfold:
container: config['singularity']['autotop']
group: 'subj'
shell:
'wb_command -surface-apply-warpfield {input.gii} {input.warp} {output.gii} && '
'cp {input.gii} {output.gii} && '
'wb_command -set-structure {output.gii} {params.structure_type} -surface-type {params.surface_type}'
' -surface-secondary-type {params.secondary_type}'

Expand Down Expand Up @@ -53,16 +52,14 @@ rule calc_unfold_template_coords:
"wb_command -set-structure {output.coords_gii} {params.structure_type}"







#subj unfolded surf might have a few vertices outside the bounding box.. this constrains all the vertices to the warp bounding box
rule constrain_surf_to_bbox:
input:
gii = bids(root='work',datatype='surf_{modality}',den='{density}',suffix='{surfname}.surf.gii', space='unfolded',hemi='{hemi}', **config['subj_wildcards']),
ref_nii = bids(root='work',space='unfold',suffix='refvol.nii.gz',**config['subj_wildcards']),
ref_nii = bids(root='work',space='unfold',datatype='seg_{modality}',suffix='refvol.nii.gz',**config['subj_wildcards']),
output:
gii = bids(root='work',datatype='surf_{modality}',den='{density}',suffix='{surfname}.surf.gii',desc='constrainbbox', space='unfolded',hemi='{hemi}', **config['subj_wildcards'])
group: 'subj'
Expand All @@ -71,7 +68,7 @@ rule constrain_surf_to_bbox:
#warp from subj unfolded to corobl
rule warp_gii_unfold2native:
input:
warp = bids(root='work',**config['subj_wildcards'],suffix='autotop/Warp_unfold2native.nii',desc='cropped',space='corobl',hemi='{hemi}',modality='{modality}'),
warp = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi}',from_='unfold',to='corobl',mode='surface'),
gii = bids(root='work',datatype='surf_{modality}',den='{density}',suffix='{surfname}.surf.gii', desc='constrainbbox',space='unfolded',hemi='{hemi}', **config['subj_wildcards'])
params:
structure_type = lambda wildcards: hemi_to_structure[wildcards.hemi],
Expand Down
94 changes: 94 additions & 0 deletions hippunfold/workflow/rules/warps.smk
Original file line number Diff line number Diff line change
@@ -0,0 +1,94 @@
rule create_native_coord_ref:
input:
coords_ap = bids(root='work',datatype='seg_{modality}',dir='AP',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='{hemi}', **config['subj_wildcards']),
output:
nii = bids(root='work',datatype='seg_{modality}',space='corobl',hemi='{hemi}',suffix='refcoords.nii.gz',**config['subj_wildcards'])
group: 'subj'
container: config['singularity']['autotop']
shell:
'c3d {input} -cmp -omc {output}'



rule create_warps:
input:
unfold_ref_nii = bids(root='work',space='unfold',datatype='seg_{modality}',suffix='refvol.nii.gz',**config['subj_wildcards']),
unfold_phys_coords_nii = bids(root='work',space='unfold',datatype='seg_{modality}',suffix='refcoords.nii.gz',**config['subj_wildcards']),
coords_ap = bids(root='work',datatype='seg_{modality}',dir='AP',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='{hemi}', **config['subj_wildcards']),
coords_pd = bids(root='work',datatype='seg_{modality}',dir='PD',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='{hemi}', **config['subj_wildcards']),
coords_io = get_laminar_coords,
native_ref_coords_nii = bids(root='work',datatype='seg_{modality}',space='corobl',hemi='{hemi}',suffix='refcoords.nii.gz',**config['subj_wildcards'])
params:
interp_method = 'linear'
resources:
mem_mb = 16000
output:
warp_unfold2native = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi,Lflip|R}',from_='unfold',to='corobl',mode='surface'),
warp_native2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi,Lflip|R}',from_='corobl',to='unfold',mode='surface'),
warpitk_unfold2native = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi,Lflip|R}',from_='unfold',to='corobl',mode='image'),

warpitk_native2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi,Lflip|R}',from_='corobl',to='unfold',mode='image')
group: 'subj'
log: bids(root='logs',**config['subj_wildcards'],hemi='{hemi,Lflip|R}',modality='seg-{modality}',suffix='create_warps.txt')
script: '../scripts/create_warps.py'



rule compose_warps_corobl2unfold_lhemi:
""" Compose corobl to unfold (unfold-template), for left hemi (ie with flip)"""
input:
native2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='Lflip',from_='corobl',to='unfold',mode='image'),
ref = bids(root='work',space='unfold',datatype='seg_{modality}',suffix='refvol.nii.gz',**config['subj_wildcards']),
flipLR_xfm = os.path.join(config['snakemake_dir'],'resources','desc-flipLR_type-itk_xfm.txt')
output:
corobl2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='L',from_='corobl',to='unfold',mode='image')
container: config['singularity']['autotop']
group: 'subj'
shell: 'ComposeMultiTransform 3 {output.corobl2unfold} -R {input.ref} {input.native2unfold} {input.flipLR_xfm}'


#consider renaming to state this composition is to "un-flip"
rule compose_warps_unfold2corobl_lhemi:
""" Compose unfold (unfold-template) to corobl, for left hemi (ie with flip)"""
input:
unfold2native = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='Lflip',from_='unfold',to='corobl',mode='image'),
flipLR_xfm = os.path.join(config['snakemake_dir'],'resources','desc-flipLR_type-itk_xfm.txt'),
unfold_ref = bids(root='work',space='unfold',datatype='seg_{modality}',suffix='refvol.nii.gz',**config['subj_wildcards']),
ref = bids(root='work',datatype='seg_{modality}',dir='AP',suffix='coords.nii.gz',desc='laplace',space='corobl',hemi='L', **config['subj_wildcards']),
output:
unfold2corobl = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='L',from_='unfold',to='corobl',mode='image')
container: config['singularity']['ants']
group: 'subj'
shell: 'antsApplyTransforms -o [{output.unfold2corobl},1] -r {input.ref} -t {input.flipLR_xfm} -t {input.unfold2native} -i {input.unfold_ref} -v'



rule compose_warps_t1_to_unfold:
""" Compose warps from T1w to unfold """
input:
corobl2unfold = bids(root='work',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi}',from_='corobl',to='unfold',mode='image'),
ref = bids(root='work',space='unfold',datatype='seg_{modality}',suffix='refvol.nii.gz',**config['subj_wildcards']),
t1w2corobl = bids(root='work',datatype='anat',**config['subj_wildcards'],suffix='xfm.txt',from_='T1w',to='corobl',desc='affine',type_='itk'),
output:
bids(root='results',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi}',from_='T1w',to='unfold',mode='image')
container: config['singularity']['ants']
group: 'subj'
shell: 'ComposeMultiTransform 3 {output} -R {input.ref} {input.corobl2unfold} {input.t1w2corobl}'


rule compose_warps_unfold_to_cropt1:
""" Compose warps from unfold to cropT1w """
input:
unfold2corobl = bids(root='results',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi}',from_='unfold',to='corobl',mode='image'),
ref = bids(root='work',datatype='seg_{modality}',suffix='cropref.nii.gz', space='T1w',hemi='{hemi}', **config['subj_wildcards']),
t1w2corobl = bids(root='work',datatype='anat',**config['subj_wildcards'],suffix='xfm.txt',from_='T1w',to='corobl',desc='affine',type_='itk'),
unfold_ref = bids(root='work',space='unfold',datatype='seg_{modality}',suffix='refvol.nii.gz',**config['subj_wildcards']),
output:
unfold2cropt1 = bids(root='results',datatype='seg_{modality}',**config['subj_wildcards'],suffix='xfm.nii.gz',hemi='{hemi}',from_='unfold',to='T1w',mode='image')
container: config['singularity']['ants']
group: 'subj'
shell: 'antsApplyTransforms -o [{output.unfold2cropt1},1] -r {input.ref} -t [{input.t1w2corobl},1] -t {input.unfold2corobl} -i {input.unfold_ref} -v'




5 changes: 3 additions & 2 deletions hippunfold/workflow/scripts/create_warps.py
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,7 @@ def summary(name, array):
#load unfolded coordinate map
#unfold_ref_nib = nib.load(snakemake.input.unfold_ref_nii)
unfold_phys_coords_nib = nib.load(snakemake.input.unfold_phys_coords_nii)
native_ref_coords_nib = nib.load(snakemake.input.native_ref_coords_nii)

#load laplace coords
coord_ap_nib = nib.load(snakemake.input.coords_ap)
Expand Down Expand Up @@ -250,8 +251,8 @@ def summary(name, array):
#and save ITK warp too
f = convert_warp_to_itk(native_to_unfold)
warpitk_unfold2native_nib = nib.Nifti1Image(f.astype(dt),
coord_ap_nib.affine,
coord_ap_nib.header)
native_ref_coords_nib.affine,
native_ref_coords_nib.header)
warpitk_unfold2native_nib.to_filename(snakemake.output.warpitk_unfold2native)


Loading

0 comments on commit 1d2924c

Please sign in to comment.