# Go to T2 contrast
cd data/t2 
# Spinal cord segmentation.
sct_deepseg_sc -i t2.nii.gz -c t2 -qc ~/qc_singleSubj

##——————————————————————————————————
#Register to template (PAMP50):
# Vertebral labeling
sct_label_vertebrae -i t2.nii.gz -s t2_seg.nii.gz -c t2 -qc ~/qc_singleSubj

# Create labels at C3 and T2 mid-vertebral levels. These labels are needed for template registration.
sct_label_utils -i t2_seg_labeled.nii.gz -vert-body 3,9 -o t2_labels_vert.nii.gz
# Generate a QC report to visualize the two selected labels on the anatomical image
sct_qc -i t2.nii.gz -s t2_labels_vert.nii.gz -p sct_label_utils -qc ~/qc_singleSubj

# Register t2->template.
sct_register_to_template -i t2.nii.gz -s t2_seg.nii.gz -l t2_labels_vert.nii.gz -c t2 -qc ~/qc_singleSubj

#Warp template objects 
sct_warp_template -d t2.nii.gz -w warp_template2anat.nii.gz -a 0 -qc ~/qc_singleSubj

#——————————————————————————————————
#Segment GM and WM 
cd ../t2s
# Segment gray matter (check QC report afterwards)
sct_deepseg_gm -i t2s.nii.gz -qc ~/qc_singleSubj

# Full spinal cord segmentation
sct_deepseg_sc -i t2s.nii.gz -c t2s -qc ~/qc_singleSubj

# Subtract GM segmentation from cord segmentation to obtain WM segmentation
sct_maths -i t2s_seg.nii.gz -sub t2s_gmseg.nii.gz -o t2s_wmseg.nii.gz

GM-informed template registration 
# Register template->t2 (using warping field generated from template<->t2 registration)
# Tips: Here we use the WM seg for the iseg/dseg fields in order to account for both the cord and the GM shape.

sct_register_multimodal -i "${SCT_DIR}/data/PAM50/template/PAM50_t2s.nii.gz" -iseg "${SCT_DIR}/data/PAM50/template/PAM50_wm.nii.gz" -d t2s.nii.gz -dseg t2s_wmseg.nii.gz -initwarp ../t2/warp_template2anat.nii.gz -initwarpinv ../t2/warp_anat2template.nii.gz -owarp warp_template2t2s.nii.gz -owarpinv warp_t2s2template.nii.gz -param step=1,type=seg,algo=rigid:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -qc ~/ qc_singleSubj

# Warp template
sct_warp_template -d t2s.nii.gz -w warp_template2t2s.nii.gz -qc ~/qc_singleSubj

——————————————————————————————————
#Calculate DTI Metrics:

cd ../dmri 
# Compute mean dMRI from dMRI data
sct_dmri_separate_b0_and_dwi -i dmri.nii.gz -bvec bvecs.txt

# Segment SC on mean dMRI data
# Note: This segmentation does not need to be accurate-- it is only used to create a mask around the cord
sct_deepseg_sc -i dmri_dwi_mean.nii.gz -c dwi -qc ~/qc_singleSubj

# Create mask (to be used in the subsequent motion correction step)
sct_create_mask -i dmri_dwi_mean.nii.gz -p centerline,dmri_dwi_mean_seg.nii.gz -size 35mm

#Motion correction for dMRI 

# Motion correction (moco)
sct_dmri_moco -i dmri.nii.gz -m mask_dmri_dwi_mean.nii.gz -bvec bvecs.txt -qc ~/qc_singleSubj -qc-seg dmri_dwi_mean_seg.nii.gz

#Register template to dMRI 

# Segment SC on motion-corrected mean dwi data (check results in the QC report)
sct_deepseg_sc -i dmri_moco_dwi_mean.nii.gz -c dwi -qc ~/qc_singleSubj

# Register template->dwi via t2s to account for GM segmentation
# Tips: In general for DWI we use the PAM50_t1 contrast, which is close to the dwi contrast (although here we are not using type=im in -param, so it will not make a difference).
# Note: the flag "-initwarpinv" provides a transformation dmri->template, in case you would like to bring all your DTI metrics in the PAM50 space (e.g. group averaging of FA maps)
sct_register_multimodal -i "${SCT_DIR}/data/PAM50/template/PAM50_t1.nii.gz" -iseg "${SCT_DIR}/data/PAM50/template/PAM50_cord.nii.gz" -d dmri_moco_dwi_mean.nii.gz -dseg dmri_moco_dwi_mean_seg.nii.gz -initwarp ../t2s/warp_template2t2s.nii.gz -initwarpinv ../t2s/warp_t2s2template.nii.gz -owarp warp_template2dmri.nii.gz -owarpinv warp_dmri2template.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -qc ~/qc_singleSubj

# Warp template
sct_warp_template -d dmri_moco_dwi_mean.nii.gz -w warp_template2dmri.nii.gz -qc ~/qc_singleSubj

#Compute DTI 

# Compute DTI metrics using dipy [1]
sct_dmri_compute_dti -i dmri_moco.nii.gz -bval bvals.txt -bvec bvecs.txt

# Tips: the flag "-method restore" estimates the tensor with robust fit (RESTORE method [2])
# Compute FA within the white matter from individual level 2 to 5
sct_extract_metric -i dti_FA.nii.gz -f label/atlas -l 51 -method map -vert 2:9 -vertfile label/template/PAM50_levels.nii.gz -perlevel 1 -o m_FA_in_wm.csv

sct_extract_metric -i dti_MD.nii.gz -f label/atlas -l 51 -method map -vert 2:9 -vertfile label/template/PAM50_levels.nii.gz -perlevel 1 -o m_MD_in_wm.csv

sct_extract_metric -i dti_AD.nii.gz -f label/atlas -l 51 -method map -vert 2:9 -vertfile label/template/PAM50_levels.nii.gz -perlevel 1 -o m_AD_in_wm.csv

sct_extract_metric -i dti_RD.nii.gz -f label/atlas -l 51 -method map -vert 2:9 -vertfile label/template/PAM50_levels.nii.gz -perlevel 1 -o m_RD_in_wm.csv