#!/bin/sh
# sienax - Structural Image Evaluation, including Normalisation, of Atrophy (X-sectional)
#
# Stephen Smith, FMRIB Image Analysis Group
#
# Copyright (C) 1999-2006 University of Oxford
#
# SHCOPYRIGHT
Usage() {
echo ""
echo "Usage: sienax <input_fileroot> [-d] [-B \"bet options\"] [-2] [-t2] [-t <t>] [-b <b>] [-r] [-lm <lesion_mask>] [-S \"segmentation options\"]"
echo ""
echo "-d : debug (don't delete intermediate files)"
echo "-B \"betopts\" : options to pass to BET brain extraction (inside double-quotes), e.g. -B \"-f 0.3\""
echo "-2 : two-class segmentation (don't segment grey and white matter separately)"
echo "-t2 : T2-weighted input image (default T1-weighted)"
echo "-t <t> : ignore from t (mm) upwards in MNI152/Talairach space"
echo "-b <b> : ignore from b (mm) downwards in MNI152/Talairach space. b should probably be -ve"
echo "-r : regional - use standard-space masks to give peripheral cortex GM volume (3-class segmentation only) and ventricular CSF volume"
echo "-lm <mask> : use lesion (or lesion+CSF) mask to remove incorrectly labelled \"grey matter\" voxels"
echo "-S \"segopts\" : options to pass to FAST segmentation (inside double-quotes), e.g. -S \"-i 20\""
echo ""
exit
}
if [ _$FSLDIR = _ ] ; then
FSLDIR=/usr/local/fsl
export FSLDIR
fi
[ _$1 = _ ] && Usage
I=`${FSLDIR}/bin/remove_ext $1`;
[ `${FSLDIR}/bin/imtest $I` = 0 ] && Usage
echo "-----------------------------------------------------------------------" > ${I}.sienax
echo "" >> ${I}.sienax
echo " SIENA - Structural Image Evaluation, using Normalisation, of Atrophy" >> ${I}.sienax
echo " part of FSL www.fmrib.ox.ac.uk/fsl" >> ${I}.sienax
echo " running cross-sectional atrophy measurement: sienax version 2.4" >> ${I}.sienax
echo " sienax $@" >> ${I}.sienax
echo "" >> ${I}.sienax
shift
debug=0
regional=0
betopts=""
segopts=""
nseg=3
stdroi=""
origin3=37 # `avwval ${FSLDIR}/etc/standard/avg152T1 origin3`
pixdim3=2 # `avwval ${FSLDIR}/etc/standard/avg152T1 pixdim3`
imtype=-t1
while [ _$1 != _ ] ; do
if [ $1 = -d ] ; then
debug=1
shift
elif [ $1 = -r ] ; then
regional=1
shift
elif [ $1 = -B ] ; then
betopts=$2
shift 2
elif [ $1 = -S ] ; then
segopts=$2
shift 2
elif [ $1 = -2 ] ; then
nseg=2
shift
elif [ $1 = -t2 ] ; then
imtype=-t2
shift
elif [ $1 = -t ] ; then
stdt=`echo $2 | sed 's/-/_/g'`
stdt=`echo "10 k $stdt $pixdim3 / $origin3 + p" | dc -`
stdroi="$stdroi -roi 0 1000000 0 1000000 0 $stdt 0 1"
shift 2
elif [ $1 = -b ] ; then
stdb=`echo $2 | sed 's/-/_/g'`
stdb=`echo "10 k $stdb $pixdim3 / $origin3 + p" | dc -`
stdroi="$stdroi -roi 0 1000000 0 1000000 $stdb 1000000 0 1"
shift 2
elif [ $1 = -lm ] ; then
lm=$2
shift 2
else
Usage
fi
done
if [ $regional = 1 ] ; then
if [ $nseg = 2 ] ; then
echo "Can't do regional analysis with 2-class segmentation"
exit
fi
fi
echo "---------- extract brain --------------------------------------------" >> ${I}.sienax
${FSLDIR}/bin/bet $I ${I}_brain -s -m $betopts >> ${I}.sienax
${FSLDIR}/bin/avwmaths++ ${I}_brain -sub `$FSLDIR/bin/avwstats++ ${I}_brain -p 0` -mas ${I}_brain_mask ${I}_brain -odt float
echo "" >> ${I}.sienax
echo "---------- register to stdairach space using brain and skull --------" >> ${I}.sienax
echo "(do not worry about histogram warnings)" >> ${I}.sienax
${FSLDIR}/bin/pairreg ${FSLDIR}/etc/standard/avg152T1_brain ${I}_brain ${FSLDIR}/etc/standard/avg152T1_skull ${I}_brain_skull ${I}2std.mat >> ${I}.sienax 2>&1
${FSLDIR}/bin/avscale ${I}2std.mat ${FSLDIR}/etc/standard/avg152T1 > ${I}2std.avscale
xscale=`grep Scales ${I}2std.avscale | awk '{print $4}'`
yscale=`grep Scales ${I}2std.avscale | awk '{print $5}'`
zscale=`grep Scales ${I}2std.avscale | awk '{print $6}'`
vscale=`echo "10 k $xscale $yscale * $zscale * p"|dc -`
echo "VSCALING $vscale" >> ${I}.sienax
echo "" >> ${I}.sienax
echo "---------- mask with std mask ---------------------------------------" >> ${I}.sienax
${FSLDIR}/bin/convert_xfm -inverse -omat ${I}2std_inv.mat ${I}2std.mat
MASK=${FSLDIR}/etc/standard/avg152T1_brain_mask_dil
if [ "$stdroi" != "" ] ; then
${FSLDIR}/bin/avwmaths $MASK $stdroi ${I}_stdmaskroi
MASK=${I}_stdmaskroi
fi
${FSLDIR}/bin/flirt -in $MASK -ref ${I}_brain -out ${I}_stdmask -applyxfm -init ${I}2std_inv.mat
${FSLDIR}/bin/avwmaths ${I}_stdmask -thr 0.5 -bin ${I}_stdmask
${FSLDIR}/bin/avwmaths ${I}_brain -mask ${I}_stdmask ${I}_stdmaskbrain
if [ $regional = 1 ] ; then
${FSLDIR}/bin/flirt -in ${FSLDIR}/etc/standard/avg152T1_strucseg_periph -ref ${I}_brain -out ${I}_stdmask_segperiph -applyxfm -init ${I}2std_inv.mat
${FSLDIR}/bin/avwmaths ${I}_stdmask_segperiph -thr 0.5 -bin ${I}_stdmask_segperiph
${FSLDIR}/bin/avwmaths ${FSLDIR}/etc/standard/avg152T1_strucseg -thr 4.5 -bin ${I}_tmpmask
${FSLDIR}/bin/flirt -in ${I}_tmpmask -ref ${I}_brain -out ${I}_stdmask_segvent -applyxfm -init ${I}2std_inv.mat
${FSLDIR}/bin/avwmaths ${I}_stdmask_segvent -thr 0.5 -bin ${I}_stdmask_segvent
/bin/rm ${I}_tmpmask*
fi
echo "" >> ${I}.sienax
echo "---------- segment tissue into types --------------------------------" >> ${I}.sienax
if [ $nseg = 2 ] ; then
${FSLDIR}/bin/fast1 -c 2 $imtype -e -ov $segopts ${I}_stdmaskbrain >> ${I}.sienax 2>&1
echo "" >> ${I}.sienax
echo "---------- convert brain volume into normalised volume --------------" >> ${I}.sienax
S=`${FSLDIR}/bin/avwstats ${I}_stdmaskbrain_pve_1 -m -v`
xa=`echo $S | awk '{print $1}'`
xb=`echo $S | awk '{print $3}'`
ubrain=`echo "2 k $xa $xb * 1 / p" | dc -`
nbrain=`echo "2 k $xa $xb * $vscale * 1 / p" | dc -`
else
if [ _$lm != _ ] ; then
${FSLDIR}/bin/avwmaths_32R $lm -bin -mul -1 -add 1 -mul ${I}_stdmaskbrain ${I}_stdmaskbrain
fi
${FSLDIR}/bin/fast1 $imtype -e -ov $segopts ${I}_stdmaskbrain >> ${I}.sienax 2>&1
if [ _$lm != _ ] ; then
${FSLDIR}/bin/avwmaths_32R $lm -bin -max ${I}_stdmaskbrain_pve_2 ${I}_stdmaskbrain_pve_2
${FSLDIR}/bin/avwmaths_32SI $lm -bin -mul 3 -max ${I}_stdmaskbrain_seg ${I}_stdmaskbrain_seg
fi
echo "" >> ${I}.sienax
echo "---------- convert brain volume into normalised volume --------------" >> ${I}.sienax
echo "" >> ${I}.sienax
echo "tissue volume unnormalised-volume" >> ${I}.sienax
if [ $regional = 1 ] ; then
${FSLDIR}/bin/avwmaths_32R ${I}_stdmaskbrain_pve_1 -mas ${I}_stdmask_segperiph ${I}_stdmaskbrain_pve_1_segperiph
S=`${FSLDIR}/bin/avwstats ${I}_stdmaskbrain_pve_1_segperiph -m -v`
xa=`echo $S | awk '{print $1}'`
xb=`echo $S | awk '{print $3}'`
uxg=`echo "2 k $xa $xb * 1 / p" | dc -`
xg=`echo "2 k $xa $xb * $vscale * 1 / p" | dc -`
echo "pgrey $xg $uxg (peripheral grey)" >> ${I}.sienax
${FSLDIR}/bin/avwmaths_32R ${I}_stdmaskbrain_pve_0 -mas ${I}_stdmask_segvent ${I}_stdmaskbrain_pve_0_segvent
S=`${FSLDIR}/bin/avwstats ${I}_stdmaskbrain_pve_0_segvent -m -v`
xa=`echo $S | awk '{print $1}'`
xb=`echo $S | awk '{print $3}'`
uxg=`echo "2 k $xa $xb * 1 / p" | dc -`
xg=`echo "2 k $xa $xb * $vscale * 1 / p" | dc -`
echo "vcsf $xg $uxg (ventricular CSF)" >> ${I}.sienax
fi
S=`${FSLDIR}/bin/avwstats ${I}_stdmaskbrain_pve_1 -m -v`
xa=`echo $S | awk '{print $1}'`
xb=`echo $S | awk '{print $3}'`
ugrey=`echo "2 k $xa $xb * 1 / p" | dc -`
ngrey=`echo "2 k $xa $xb * $vscale * 1 / p" | dc -`
echo "GREY $ngrey $ugrey" >> ${I}.sienax
S=`${FSLDIR}/bin/avwstats ${I}_stdmaskbrain_pve_2 -m -v`
xa=`echo $S | awk '{print $1}'`
xb=`echo $S | awk '{print $3}'`
uwhite=`echo "2 k $xa $xb * 1 / p" | dc -`
nwhite=`echo "2 k $xa $xb * $vscale * 1 / p" | dc -`
echo "WHITE $nwhite $uwhite" >> ${I}.sienax
ubrain=`echo "2 k $uwhite $ugrey + 1 / p" | dc -`
nbrain=`echo "2 k $nwhite $ngrey + 1 / p" | dc -`
fi
echo "BRAIN $nbrain $ubrain" >> ${I}.sienax
${FSLDIR}/bin/overlay 1 1 -c ${I} -a ${I}_stdmaskbrain_seg 1.9 5 ${I}_render
if [ $regional = 1 ] ; then
${FSLDIR}/bin/overlay 1 1 -c ${I} -a ${I}_stdmaskbrain_pve_1_segperiph 0.3 0.7 ${I}_periph_render
${FSLDIR}/bin/overlay 1 1 -c ${I} -a ${I}_stdmaskbrain_pve_0_segvent 0.3 0.7 ${I}_vent_render
fi
if [ $debug = 0 ] ; then
/bin/rm -f `$FSLDIR/bin/imglob ${I}_brain* ${I}_stdmask*`
/bin/rm -f ${I}2std.avscale ${I}2std_inv.mat
fi
echo "$nbrain"