Skip to content

Diffusion Imaging Post Processing

Jump to bottom
Anthony Dick edited this page Sep 18, 2019 · 7 revisions

Run Diffusion Postprocessing

#sample command: #bash dwi_preproc_submitjobs_AHEAD.sh 1120012 S1 /home/data/dcnlab/AHEAD/dset>

sub=$1
#sub=4152024
echo $sub
ses=$2
#ses=S1
echo $ses
workdir=$3
#workdir=/home/data/dcnlab/AHEAD/dset
echo $workdir
MODULESHOME=/home/share/Modules/3.2.10
. $MODULESHOME/../global/profile.modules

module add miniconda/2.7
module load slicer/4.8.1
module load dtiprep/1.2.4
module load tortoise/3.1.1
module load afni/17.3.06
module load dsistudio/1.0
module load fsl/5.0.8
TORTDIR=/home/applications/tortoise/3.1.1/DIFFPREPV311/bin/bin
DRBUDDIDIR=/home/applications/tortoise/3.1.1/DRBUDDIV311/bin
SLICERDIR='/home/applications/slicer/4.8.1/lib/Slicer-4.8/cli-modules'
DIFFCALCDIR=/home/applications/tortoise/3.1.1/DIFFCALC/DIFFCALCV311

run=$(dir $workdir/sub-$sub/ses-$ses/dwi/sub-"$sub"_ses-"$ses"_run-*.nii.gz)
run=${run##*/}
run=${run##*run-}
run=${run%%_dwi*}

cd $workdir
mkdir -p $sub/ses-$ses/dwi_processed
cd $workdir/$sub/ses-$ses/dwi_processed
cp -rf $workdir/sub-$sub/ses-$ses/dwi $workdir/$sub/ses-$ses/dwi_processed
cp -rf $workdir/sub-$sub/ses-$ses/fmap $workdir/$sub/ses-$ses/dwi_processed
cp -rf $workdir/sub-$sub/ses-$ses/anat $workdir/$sub/ses-$ses/dwi_processed
newworkdir=$workdir/$sub/ses-$ses/dwi_processed
echo "New working directory is $newworkdir"
mkdir -p $newworkdir/dwi/dwi_preproc

Quality Control before preprocessing

#cd $workdir/sub-$sub/ses-$ses/dwi/
#cp sub-"$sub"_ses-"$ses"_run-"$run"_dwi.bval sub-"$sub"_ses-"$ses"_run-"$run"_dwi.bval
cd $workdir/scripts/QAscripts/QAscripts_version4/
bash qa_dti_v4.sh $newworkdir/dwi/sub-"$sub"_ses-"$ses"_run-01_dwi.nii.gz $newworkdir/dwi/sub-"$sub"_ses-"$ses"_run-01_dwi.bval $newworkdir/dwi/sub-"$sub"_ses-"$ses"_run-01_dwi.bvec 0 $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_preQC.csv

echo "End quality assurance script."

Convert the DWI NIFTI/bvec/bvals to NRRD format

$SLICERDIR/DWIConvert --inputVolume ../sub-"$sub"_ses-"$ses"_run-"$run"_dwi.nii.gz --inputBVectors ../sub-"$sub"_ses-"$ses"_run-"$run"_dwi.bvec --inputBValues ../sub-"$sub"_ses-"$ses"_run-"$run"_dwi.bval --conversionMode FSLToNrrd --outputVolume sub-"$sub"_ses-"$ses"_run-"$run"_dwi.nrrd

Run DTIPrep quality control

Convert the NRRD file to short format for the next conversion

Convert the NRRD file back to FSL format (NIFTI/bvec/bval)

Convert the T1 scan and prep for TORTOISE

    #These next two commands were supposed to make the T1 better for DRBUDDI but didn't
    #bet ../../anat/sub-"$sub"_ses-"$ses"_T1w.nii.gz ../../anat/sub-"$sub"_ses-"$ses"_T1w_bet.nii.gz -B -f 0.4 -g 0 #run FSL bet command
    #fslmaths ../../anat/sub-"$sub"_ses-"$ses"_T1w.nii.gz -sub ../../anat/sub-"$sub"_ses-"$ses"_T1w_bet.nii.gz -div 10 -add ../../anat/sub-"$sub"_ses-"$ses"_T1w_bet.nii.gz ../../anat/sub-"$sub"_ses-"$ses"_T1watt.nii.gz #run fsl maths command

    if [ -e ../../anat/sub-"$sub"_ses-"$ses"_T2w.nii.gz ]; then
        anatmod=T2w
    else
        anatmod=T1w
    fi

    $TORTDIR/ImportNIFTI -i ../../anat/sub-"$sub"_ses-"$ses"_"$anatmod".nii.gz --phase horizontal
    mv ../../anat/sub-"$sub"_ses-"$ses"_"$anatmod".nii_proc ./

Convert the DWI scan to TORTOISE format

    $TORTDIR/ImportNIFTI -i sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short.nii -b sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short.bval -v sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short.bvec -p horizontal

Convert the AP field map for DWI scan to TORTOISE format

    if [ -e $newworkdir/fmap/sub-"$sub"_ses-"$ses"_acq-dwi_dir-AP_run-*_epi.nii.gz ]; then
    	gunzip $newworkdir/fmap/sub-"$sub"_ses-"$ses"_acq-dwi_dir-AP_run-*_epi.nii.gz
	else
        fmapfile=$(dir $newworkdir/fmap/sub-"$sub"_ses-"$ses"_acq-dwi_dir-AP_run-*_epi.nii)
    fi
    fmapfile=${fmapfile%.*}
    tmpfmapfile=${fmapfile%.*}
    if [ -e $tmpfmapfile.bval ]; then
        rm $tmpfmapfile.bval
    fi
    if [ -e $tmpfmapfile.bvec ]; then
        rm $tmpfmapfile.bvec
    fi
    echo "0" > $tmpfmapfile.bval
    for ((i=1;i<=3;i++)); do
        echo "0" >> $tmpfmapfile.bvec
    done
    fmapfile=${fmapfile##*/}
    $TORTDIR/ImportNIFTI -i ../../fmap/$fmapfile --phase horizontal
    fmapfile=${fmapfile%.*}
    mv ../../fmap/"$fmapfile"_proc ./

Get voxel resolution

    vsize=($(3dinfo -ad3 sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short.nii))

echo "Begin DIFFPREP"
    #Do DIFFPREP’ing expecting DRBUDDI
    #$TORTDIR/DIFFPREP -i sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short.list --will_be_drbuddied 1 --keep_intermediate 0 --upsampling all --res ${vsize[0]} ${vsize[1]} ${vsize[2]}
	#Do DIFFPREP’ing without DRBUDDI
    $TORTDIR/DIFFPREP -i sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short.list --will_be_drbuddied 0 --keep_intermediate 0 --upsampling all --res ${vsize[0]} ${vsize[1]} ${vsize[2]}

cd $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc

Convert TORTOISE bmtxt table to separate bvec and bval files

echo "Convert TORTOISE bmtxt file"
    $TORTDIR/TORTOISEBmatrixToFSLBVecs sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bmtxt

#echo "Begin DRBUDDI"

#Do DRBUDDI’ing
#$DRBUDDIDIR/DR_BUDDI_withoutGUI --up_data sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc.list --down_data "$fmapfile"_proc/$fmapfile.list --structural sub-"$sub"_ses-"$ses"_"$anatmod".nii_proc/sub-"$sub"_ses-"$ses"_"$anatmod".nii --distortion_level medium --metric CCSK --metric MSJac --metric CCJac --restrict-deformation 1 --enforce_deformation_antisymmetry --n_DWIs 0 --step 0 --res ${vsize[0]} ${vsize[1]} ${vsize[2]}

#Convert TORTOISE bmtxt table to separate bvec and bval files
#$TORTDIR/TORTOISEBmatrixToFSLBVecs sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI.bmtxt

echo "Round bvals" perl -i -pe 's/(\d*.\d*)/int($1+0.5)/ge' $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bvals

echo "Begin quality assurance script after postprocessing" #Quality Control after preprocessing cd $workdir/scripts/QAscripts/QAscripts_version4/

#bash qa_dti_v4.sh $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI.nii $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI.bvals $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI.bvecs 0 $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_postQC.csv

#If not DRBUDDI bash qa_dti_v4.sh $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.nii $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bvals $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bvecs 0 $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_postQC.csv

#cd $workdir/sub-$sub/ses-$ses/dwi/dwi_preproc #echo "Run DIFFCALC for chi square map" #$DIFFCALCDIR/EstimateTensorWLLS -i /$workdir/sub-$sub/ses-$ses/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI.list #$DIFFCALCDIR/EstimateTensorNLLS -i /$workdir/sub-$sub/ses-$ses/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI.list --save_CS #$DIFFCALCDIR/ComputeFAMap /$workdir/sub-$sub/ses-$ses/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI_L0_DT.nii #$DIFFCALCDIR/ComputeTRMap /$workdir/sub-$sub/ses-$ses/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI_L0_DT.nii #$DIFFCALCDIR/ComputeDECMap /$workdir/sub-$sub/ses-$ses/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc_DRBUDDI_L0_DT.nii

#Generate DSIStudio Source files

echo "Generate DSIStudio Maps"

cd $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc

dsi_studio --action=src --source=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.nii --bval=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bvals --bvec=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bvecs --output=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.src.gz

#Run DSIStudio Reconstruction using GQI

dsi_studio --action=rec --source=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.src.gz --method=4 --param0=1.25

#Run DSIStudio Reconstruction using DTI

dsi_studio --action=rec --source=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.src.gz --method=1

echo "Perform topup"

cd $newworkdir/fmap/ cp $workdir/scripts/acquisition_parameters.txt ./ fslmerge -t $sub.merge.APPA.nii sub-"$sub"_ses-"$ses"_acq-dwi_dir-AP_run-02_epi.nii sub-"$sub"_ses-"$ses"_acq-dwi_dir-PA_run-01_epi.nii.gz topup --imain=$sub.merge.APPA.nii --datain=acquisition_parameters.txt --config=b02b0.cnf --out=my_topup_results --fout=my_field --iout=my_unwarped_images --verbose cp $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.nii ./ applytopup --imain=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.nii --inindex=1 --datain=acquisition_parameters.txt --topup=my_topup_results --out=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC_topup.nii --method=jac --verbose cp ./sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC_topup.nii.gz $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc/ cd $newworkdir/dwi/dwi_preproc/sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_proc gunzip sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC_topup.nii.gz

echo "Generate DSIStudio Maps After TopUp"

dsi_studio --action=src --source=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC_topup.nii --bval=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bvals --bvec=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC.bvecs --output=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC_topup.src.gz

#Run DSIStudio Reconstruction using GQI

dsi_studio --action=rec --source=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC_topup.src.gz --method=4 --param0=1.25

#Run DSIStudio Reconstruction using DTI

dsi_studio --action=rec --source=sub-"$sub"_ses-"$ses"_run-"$run"_dwi_QCed_short_DMC_topup.src.gz --method=1

Clone this wiki locally