1#!/bin/tcsh 2 3## help the pitiful luser? 4 5set dohelp = 0 6if( $#argv == 0 ) set dohelp = 1 7if( dohelp == 0 )then 8 if( $argv[1] == "-help" ) set dohelp = 1 9endif 10 11if( $dohelp )then 12 echo " " 13 echo "Usage: @banded_clustsim list_of_directories" 14 echo " " 15 echo "The input to this script is a list of afni_proc.py results directories." 16 echo "The steps followed are:" 17 echo " a) determine the collective passband for the stimuli from the" 18 echo " X.nocensor.xmat.1D matrix file in each directory (program stimband)." 19 echo " b) bandpass each directory's pre-whitened residual file, and determine" 20 echo " its spatial ACF = AutoCorrelation Function (program 3dFWHMx)" 21 echo " c) average the spatial ACF parameters across directories" 22 echo " d) use the resulting smoothness to estimate the cluster-size" 23 echo " thresholds (program 3dClustSim)" 24 echo "The goal is to provide cluster significance levels for group tests." 25 echo " " 26 echo "** This script is currently in development and is experimental -- RWCox **" 27 exit 0 28endif 29 30## scan inputs and assemble the lists of files to process 31 32set nbad = 0 33set matlist = ( ) 34set errlist = ( ) ; set nerrWH = 0 ; set nerrOT = 0 35set msklist = ( ) ; set nmskGM = 0 ; set nmskGR = 0 36 37## We prefer whitened residuals (from 3dREMLfit -Rwherr), but 38## will process non-whitened residuals if that's all there is 39 40set errts_template_WH1 = 'errts.\*.WH+tlrc.HEAD' 41set errts_template_WH2 = 'errts_whitened.\*+tlrc.HEAD' 42set errts_template_OTHER = 'errts.\*+tlrc.HEAD' 43 44echo "========== Scanning for input files ==========" 45 46foreach aaa ( $argv ) 47 48 if( ! -d $aaa )then 49 echo "** ERROR: directory $aaa does not exist" ; @ nbad++ ; continue 50 endif 51 52 if( ! -f $aaa/X.nocensor.xmat.1D )then 53 echo "** ERROR: matrix $aaa/X.nocensor.xmat.1D does not exist" ; @ nbad++ 54 else 55 set matlist = ( $matlist $aaa/X.nocensor.xmat.1D ) 56 endif 57 58 unset mmm 59 if( -f $aaa/mask_GM_resam+tlrc.HEAD )then 60 set mmm = $aaa/mask_GM_resam+tlrc.HEAD ; @ nmskGM++ 61 else if( -f $aaa/mask_group+tlrc.HEAD )then 62 set mmm = $aaa/mask_group+tlrc.HEAD ; @ nmskGR++ 63 endif 64 if( $?mmm == 0 )then 65 echo "** ERROR: neither mask $aaa/mask_GM_resam+tlrc.HEAD nor $aaa/mask_group+tlrc.HEAD is found" 66 @ nbad++ 67 else 68 set msklist = ( $msklist $mmm ) 69 endif 70 71 unset fff 72 if( $?fff == 0 )then 73 set eee = ( `find $aaa -name errts\*WH\*+tlrc.HEAD` ) 74 if( $#eee > 0 )then 75 set fff = $eee[1] 76 @ nerrWH++ 77 endif 78 endif 79 80 if( $?fff == 0 )then 81 set eee = ( `find $aaa -name errts\*whitened\*+tlrc.HEAD` ) 82 if( $#eee > 0 )then 83 set fff = $eee[1] 84 @ nerrWH++ 85 endif 86 endif 87 88 if( $?fff == 0 )then 89 set eee = ( `find $aaa -name errts\*+tlrc.HEAD` ) 90 if( $#eee > 0 ) 91 set fff = $eee[1] 92 @ nerrOT++ 93 endif 94 endif 95 96 if( $?fff == 0 )then 97 echo "** ERROR: no file of form $errts'*'+tlrc.HEAD in directory $aaa" 98 @ nbad++ 99 else 100 set errlist = ( $errlist $fff ) 101 endif 102 103end 104 105if( ! -f $argv[1]/mask_group+tlrc.HEAD )then 106 echo "** ERROR: mask file $argv[1]/mask_group+tlrc.HEAD does not exist" ; nbad++ 107endif 108 109if( $nbad > 0 )then 110 echo "** FATAL ERROR: Can't continue after such problems :-(" 111 exit 1 112endif 113 114echo "-- Using $nerrWH prewhitened errts datasets and $nerrOT other errts datasets" 115echo "-- Using $nmskGM gray-matter mask datasets and $nmskGR group mask datasets" 116 117## echo "-- list of errts datasets:" 118## echo " $errlist" 119## echo "-- list of mask datasets:" 120## echo " $msklist" 121 122## process all the matrices at once to get the collective passband 123 124echo "========== Computing stimulus passband ==========" 125 126set bbb = ( `stimband $matlist` ) 127 128if( $#bbb < 2 || "$bbb[1]" == "$bbb[2]" )then 129 echo "program stimband failed on the matrices for some reason :-(" 130 exit 1 131endif 132 133echo "=" 134echo "Stimulus passband = $bbb[1] to $bbb[2] Hz" 135echo "=" 136 137## bandpass and 3dFWHMx each noise file separately 138 139set rrr = `3dnewid -fun` 140 141echo "========== Bandpassing and computing ACF smoothness ==========" 142 143echo "# 3dFWHMx output" > BPtemp.$rrr.1D 144 145foreach nnn ( `count -dig 1 1 $#errlist` ) 146 147 set eee = $errlist[$nnn] 148 set mmm = $msklist[$nnn] 149 150 3dTproject -input $eee -mask $mmm -prefix BPtemp.$rrr.nii -passband $bbb[1] $bbb[2] 151 152# the 'tail -1' is to keep only the 2nd (ACF) line from 3dFWHMx 153 154 3dFWHMx -acf BPtemp.$nnn.$rrr.1D -mask $mmm BPtemp.$rrr.nii | tail -1 >> BPtemp.$rrr.1D 155 156 \rm BPtemp.$rrr.nii 157 158end 159 160## compute the average of the ACF outputs 161 162set bbb = ( `1dsum -mean BPtemp.$rrr.1D` ) 163 164cp BPtemp.$rrr.1D BPacf.1D 165 166echo "Average ACF FWHM = $bbb[4] mm" 167 168## run 3dClustSim 169 170echo "========== Running 3dClustSim ==========" 171 1723dClustSim -acf $bbb[1] $bbb[2] $bbb[3] \ 173 -mask $argv[1]/mask_group+tlrc.HEAD \ 174 -both -MEGA -sumup -sumup -prefix BPclust 175 176\rm BPtemp*${rrr}* 177 178echo "========== e Finito ==========" 179exit 0 180