generate_QA_stats

#!/bin/bash

usage()
{
	echo "Usage: `basename $0` <subject id> <session no>"
	echo
	echo "<subject id> format should be numbers only: e.g. 100044 or 100054"
	echo
	echo "<session no> should be a single number: e.g. 1 or 2"
	echo
}

if [ $# -lt 2 ]; then
	echo
	echo "./generate_QA_data <subject id> <session no>"
	echo "e.g. ./generate_QA_data 100044 1"
	echo
	exit
fi

## define paths and constants
SUBJECTDIR=/mnt/panuc/udallp2/subjects

SUBID=$1
SESSION=$2
CURRDIR=$SUBJECTDIR/$SUBID/session$SESSION
GROUP=`cat $CURRDIR/0_group`
OUTPUTFILE=$CURRDIR/QA/${SUBID}_QA_stats.csv

# note: missing files are automatically flagged
# the cutoff points for when measures are flagged
skullstrip_cutoff=0.9			# skullstripDICE < cutoff gets flagged
reg_warn_cutoff=0.75			# regDICE < cutoff gets flagged
absdisp_cutoff=2.0			# mean absdisp > cutoff gets flagged
reldisp_cutoff=2.0			# mean reldisp > cutoff gets flagged
percent_outliers_cutoff=10.0		# percent_outliers > cutoff gets flagged
fdvals_cutoff=0.5			# fdvals > cutoff gets flagged
num_vols_cutoff=240			# num_vols != 240 gets flagged

# expected scan parameters, actual != expected gets flagged 
slices=37
echoes=3
dyn_scans=240
FOV_RL=224
FOV_AP=224
FOV_FH=129.5
TR=2500.001
TE=9.5
voxel_size=3.5
slice_thickness=3.5
slice_orientation=1			# note: slice_orientation of 1 in par file is "transverse"

if [[ $GROUP == "PATIENT" ]]; then
	tasks=( rest_on rest_off axcpt_on axcpt_off )
else
	tasks=( rest_on axcpt_on )
fi

> $OUTPUTFILE

# write header to output file
echo "measure,data,flag" >> $OUTPUTFILE

########################### get and flag task-INDEPENDENT measures  #############################
# flag scan parameters
for file in $CURRDIR/parrec/*WIP_ME*.PAR
do (
	shopt -s nocasematch
	if [[ $file == *"ME_RS"* ]]; then
		if [[ $file == *"off"* ]]; then
			scan_type=RS_off
		else
			scan_type=RS_on
		fi
	elif [[ $file == *"ME_Task"* ]]; then
		if [[ $file == *"off"* ]]; then	
			scan_type=Task_off
		else
			scan_type=Task_on
		fi
	fi
	slices_actual=`echo $(grep slices $file | sed 's/.$//' | awk '{print $NF}')`
	echoes_actual=`echo $(grep echoes $file | sed 's/.$//' | awk '{print $NF}')`
	dyn_scans_actual=`echo $(grep dynamics $file | sed 's/.$//' | awk '{print $NF}')`
	FOV=`grep FOV $file`
	FOV_RL_actual=`echo $FOV | sed 's/.$//' | awk '{print $NF}'`
	FOV_AP_actual=`echo $FOV | awk '{print $(NF-2)}'`
	FOV_FH_actual=`echo $FOV | awk '{print $(NF-1)}'`
	TR_actual=`echo $(grep Repetition $file | awk '{print $(NF-1)}')`
	parameter_line=`sed "101q;d" $file`
	parameter_line=`echo $parameter_line | sed 's/.$//'`
	params=($parameter_line)
	TE_actual=${params[30]}
	voxel_size_actual=${params[28]}
	slice_thickness_actual=${params[22]}
	slice_orientation_actual=${params[25]}
	echo "${scan_type}_slices,${slices_actual},"`awk -v n1=$slices -v n2=$slices_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_echoes,${echoes_actual},"`awk -v n1=$echoes -v n2=$echoes_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_dyn_scans,${dyn_scans_actual},"`awk -v n1=$dyn_scans -v n2=$dyn_scans_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_FOV_RL,${FOV_RL_actual},"`awk -v n1=$FOV_RL -v n2=$FOV_RL_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_FOV_AP,${FOV_AP_actual},"`awk -v n1=$FOV_AP -v n2=$FOV_AP_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_FOV_FH,${FOV_FH_actual},"`awk -v n1=$FOV_FH -v n2=$FOV_FH_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_TR,${TR_actual},"`awk -v n1=$TR -v n2=$TR_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_TE,${TE_actual},"`awk -v n1=$TE -v n2=$TE_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_voxel_size,${voxel_size_actual},"`awk -v n1=$voxel_size -v n2=$voxel_size_actual \
		'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_slice_thickness,${slice_thickness_actual},"`awk -v n1=$slice_thickness \
		-v n2=$slice_thickness_actual 'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
	echo "${scan_type}_slice_orientation,${slice_orientation_actual},"`awk -v n1=$slice_orientation \
		-v n2=$slice_orientation_actual 'BEGIN {if (n1==n2) print "0"; else print "1";}'` >> $OUTPUTFILE
) done
flagged_scan_params=`grep ",1$" $OUTPUTFILE`
flagged_scan_params=`echo $flagged_scan_params | tr " " "\n"`
RS_off_flagged=`grep "RS_off" <<< "$flagged_scan_params"`
echo ME_RS_off,NA,`[[ -z $RS_off_flagged ]] && echo 0 || echo 1` >> $OUTPUTFILE
RS_on_flagged=`grep "RS_on" <<< "$flagged_scan_params"`
echo ME_RS_on,NA,`[[ -z $RS_on_flagged ]] && echo 0 || echo 1` >> $OUTPUTFILE
Task_off_flagged=`grep "Task_off" <<< "$flagged_scan_params"`
echo ME_Task_off,NA,`[[ -z $Task_off_flagged ]] && echo 0 || echo 1` >> $OUTPUTFILE
Task_on_flagged=`grep "RS_off" <<< "$flagged_scan_params"`
echo ME_Task_on,NA,`[[ -z $Task_on_flagged ]] && echo 0 || echo 1` >> $OUTPUTFILE

# get skullstripDICE 
if [[ -f $CURRDIR/mprage/skullstripDICE ]]; then
	ssDICE=`grep 'value="' $CURRDIR/mprage/skullstripDICE | awk '{print $7}' | grep -o '[0-9].[0-9]*'`
	echo -ne "skullstripDICE,$ssDICE," >> $OUTPUTFILE
	echo $ssDICE'<'$skullstrip_cutoff | bc -l >> $OUTPUTFILE
else
	echo "skullstripDICE,NA,1" >> $OUTPUTFILE
fi

########################## get and flag task-DEPENDENT measures  ###############################
for TASK in "${tasks[@]}"
do (
	# flag registration DICE coefficients
	if [[ -f $CURRDIR/$TASK/${TASK}_regDICE ]]; then
		regDICE=`grep 'value="' $CURRDIR/$TASK/${TASK}_regDICE | awk '{print $7}' | grep -o '[0-9].[0-9]*'`
		echo -n "${TASK}_regDICE,$regDICE," >> $OUTPUTFILE
		echo $regDICE'<'$reg_warn_cutoff | bc -l >> $OUTPUTFILE
	else
		echo "${TASK}_regDICE,NA,1" >> $OUTPUTFILE
	fi
	
	# flag motion mean metrics
	if [[ -f $CURRDIR/mcflirt_data/${TASK}_data_abs_mean.rms ]]; then
		absdisp_mean=`cat $CURRDIR/mcflirt_data/${TASK}_data_abs_mean.rms`
		echo -n "${TASK}_abs_mean_displacement,$absdisp_mean," >> $OUTPUTFILE
		echo $absdisp_mean'>'$absdisp_cutoff | bc -l >> $OUTPUTFILE
	else
		echo "${TASK}_abs_mean_displacement,NA,1" >> $OUTPUTFILE
	fi
	if [[ -f $CURRDIR/mcflirt_data/${TASK}_data_rel_mean.rms ]]; then
		reldisp_mean=`cat $CURRDIR/mcflirt_data/${TASK}_data_rel_mean.rms`
		echo -n "${TASK}_rel_mean_displacement,$reldisp_mean," >> $OUTPUTFILE
		echo $reldisp_mean'>'$reldisp_cutoff | bc -l >> $OUTPUTFILE
	else
		echo "${TASK}_rel_mean_displacement,NA,1" >> $OUTPUTFILE
	fi

	# flag outliers
	if [[ -f $CURRDIR/$TASK/${TASK}_e002_outliers_volumes.txt ]]; then
		total_outliers=`cat $CURRDIR/$TASK/${TASK}_e002_outliers_volumes.txt | wc -l`
		echo "${TASK}_total_outliers,$total_outliers,0" >> $OUTPUTFILE
	else 
		echo "${TASK}_total_outliers,NA,1" >> $OUTPUTFILE
	fi
	if [[ -f $CURRDIR/mcflirt_data/${TASK}_data.nii.gz ]]; then
		num_vols=`fslval $CURRDIR/mcflirt_data/${TASK}_data.nii.gz dim4` 
		num_vols=`echo $num_vols | xargs`
		echo -n "${TASK}_num_vols,$num_vols," >> $OUTPUTFILE
		if [ "$num_vols" == "$num_vols_cutoff" ]; then
			echo 0 >> $OUTPUTFILE
		else
			echo 1 >> $OUTPUTFILE
		fi
	else
		echo "${TASK}_num_vols,NA,1" >> $OUTPUTFILE
	fi
	if [[ -f $CURRDIR/$TASK/${TASK}_e002_outliers_volumes.txt ]] \
			&& [[ -f $CURRDIR/mcflirt_data/${TASK}_data.nii.gz ]]; then
		percent_outliers=`bc <<< "scale = 10; (( 100*$total_outliers/$num_vols ))"`
		echo -ne "${TASK}_percent_outliers,$percent_outliers," >> $OUTPUTFILE
		echo $percent_outliers'>'$percent_outliers_cutoff | bc -l >> $OUTPUTFILE
	else 
		echo "${TASK}_percent_outliers,NA,1" >> $OUTPUTFILE
	fi
	

	# flag fdvals means
	if [[ -f $CURRDIR/$TASK/${TASK}_e002_fd_vals ]]; then
		fdvals_e002_mean=`awk '{ s+=$1 } END {print s/NR}' RS="\n" $CURRDIR/$TASK/${TASK}_e002_fd_vals`	
		echo -n "${TASK}_fdvals_e002_mean,$fdvals_e002_mean," >> $OUTPUTFILE
		echo $fdvals_e002_mean'>'$fdvals_cutoff | bc -l >> $OUTPUTFILE
	else
		echo "${TASK}_fdvals_e002_mean,NA,1" >> $OUTPUTFILE
	fi
	if [[ -f $CURRDIR/$TASK/${TASK}_e00213_tsoc_fd_vals ]]; then
		fdvals_tsoc_mean=`awk '{ s+=$1 } END {print s/NR}' RS="\n" $CURRDIR/$TASK/${TASK}_e00213_tsoc_fd_vals`	
		echo -n "${TASK}_fdvals_tsoc_mean,$fdvals_tsoc_mean," >> $OUTPUTFILE
		echo $fdvals_tsoc_mean'>'$fdvals_cutoff | bc -l >> $OUTPUTFILE
	else
		echo "${TASK}_fdvals_tsoc_mean,NA,1" >> $OUTPUTFILE
	fi
	if [[ -f $CURRDIR/$TASK/${TASK}_e00213_medn_fd_vals ]]; then
		fdvals_medn_mean=`awk '{ s+=$1 } END {print s/NR}' RS="\n" $CURRDIR/$TASK/${TASK}_e00213_medn_fd_vals`	
		echo -n "${TASK}_fdvals_medn_mean,$fdvals_medn_mean," >> $OUTPUTFILE
		echo $fdvals_medn_mean'>'$fdvals_cutoff | bc -l >> $OUTPUTFILE
	else
		echo "${TASK}_fdvals_medn_mean,NA,1" >> $OUTPUTFILE
	fi

	# flag dvars means
	thresh_exists=1
	if [[ -f $CURRDIR/$TASK/${TASK}_e002_dvars_thresh ]]; then
		dvarsvals_e002_thresh=`cat $CURRDIR/$TASK/${TASK}_e002_dvars_thresh`
		echo "${TASK}_dvarsvals_e002_thresh,$dvarsvals_e002_thresh,0" >> $OUTPUTFILE
	else
		echo "${TASK}_dvarsvals_e002_thresh,NA,1" >> $OUTPUTFILE
		thresh_exists=0
	fi
	if [[ -f $CURRDIR/$TASK/${TASK}_e002_dvars_vals ]]; then 
		dvarsvals_e002_mean=`awk '{ total += $1; count++ } END {print total/count}' $CURRDIR/$TASK/${TASK}_e002_dvars_vals`	
		echo -n "${TASK}_dvarsvals_e002_mean,$dvarsvals_e002_mean," >> $OUTPUTFILE
		if [[ $thresh_exists -eq 1 ]]; then
			echo $dvarsvals_e002_mean'>'$dvarsvals_e002_thresh | bc -l >> $OUTPUTFILE
		else
			echo 1 >> $OUTPUTFILE
		fi
	else
		echo "${TASK}_dvarsvals_e002_mean,NA,1" >> $OUTPUTFILE
	fi
	thresh_exists=1
	if [[ -f $CURRDIR/$TASK/${TASK}_e00213_tsoc_dvars_thresh ]]; then
		dvarsvals_tsoc_thresh=`cat $CURRDIR/$TASK/${TASK}_e00213_tsoc_dvars_thresh`
		echo "${TASK}_dvarsvals_tsoc_thresh,$dvarsvals_tsoc_thresh,0" >> $OUTPUTFILE
	else
		echo "${TASK}_dvarsvals_tsoc_thresh,NA,1" >> $OUTPUTFILE
		thresh_exists=0
	fi
	if [[ -f $CURRDIR/$TASK/${TASK}_e00213_tsoc_dvars_vals ]]; then 
		dvarsvals_tsoc_mean=`awk '{ total += $1; count++ } END {print total/count}' $CURRDIR/$TASK/${TASK}_e00213_tsoc_dvars_vals`	
		echo -n "${TASK}_dvarsvals_tsoc_mean,$dvarsvals_tsoc_mean," >> $OUTPUTFILE
		if [[ $thresh_exists -eq 1 ]]; then
			echo $dvarsvals_tsoc_mean'>'$dvarsvals_tsoc_thresh | bc -l >> $OUTPUTFILE
		else
			echo 1 >> $OUTPUTFILE
		fi
	else
		echo "${TASK}_dvarsvals_tsoc_mean,NA,1" >> $OUTPUTFILE
	fi
	thresh_exists=1
	if [[ -f $CURRDIR/$TASK/${TASK}_e00213_medn_dvars_thresh ]]; then
		dvarsvals_medn_thresh=`cat $CURRDIR/$TASK/${TASK}_e00213_medn_dvars_thresh`
		echo "${TASK}_dvarsvals_medn_thresh,$dvarsvals_medn_thresh,0" >> $OUTPUTFILE
	else
		echo "${TASK}_dvarsvals_medn_thresh,NA,1" >> $OUTPUTFILE
		thresh_exists=0
	fi
	if [[ -f $CURRDIR/$TASK/${TASK}_e00213_medn_dvars_vals ]]; then 
		dvarsvals_medn_mean=`awk '{ total += $1; count++ } END {print total/count}' $CURRDIR/$TASK/${TASK}_e00213_medn_dvars_vals`
		echo -n "${TASK}_dvarsvals_medn_mean,$dvarsvals_medn_mean," >> $OUTPUTFILE
		if [[ $thresh_exists -eq 1 ]]; then
			echo $dvarsvals_medn_mean'>'$dvarsvals_medn_thresh | bc -l >> $OUTPUTFILE
		else
			echo 1 >> $OUTPUTFILE
		fi
	else
		echo "${TASK}_dvarsvals_medn_mean,NA,1" >> $OUTPUTFILE
	fi
) done