#!/bin/tcsh -f
PARSE:
echo "Parsing ..."
set Narg = $#
set tmp = ___tmp_WO_
set pd_default = 0
set cnt = 1
set ref_in = ''
set anat_in = ''
set KeepTmp = 0
set verb = 0
set pd = $pd_default
set suf = ""
set ss_opt = ''
set dxyz_list = ( 1 1 1 )
if ("$1" == '') goto HELP
while ($cnt <= $Narg)
set donext = 1;
if ($donext && "$argv[$cnt]" == "-help" || "$argv[$cnt]" == "-h") then
goto HELP
endif
if ($donext && "$argv[$cnt]" == "-dxyz") then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need voxel size after -dxyz"
goto END
else
@ cnt ++
set dxyz_list[1] = "$argv[$cnt]"
set dxyz_list[2] = "$argv[$cnt]"
set dxyz_list[3] = "$argv[$cnt]"
set donext = 0
endif
endif
if ($donext && "$argv[$cnt]" == "-dx") then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need voxel size after -dx"
goto END
else
@ cnt ++
set dxyz_list[1] = "$argv[$cnt]"
set donext = 0
endif
endif
if ($donext && "$argv[$cnt]" == "-dy") then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need voxel size after -dy"
goto END
else
@ cnt ++
set dxyz_list[2] = "$argv[$cnt]"
set donext = 0
endif
endif
if ($donext && "$argv[$cnt]" == "-dz") then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need voxel size after -dz"
goto END
else
@ cnt ++
set dxyz_list[3] = "$argv[$cnt]"
set donext = 0
endif
endif
if ($donext && "$argv[$cnt]" == "-suffix") then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need a string after -suffix"
goto END
else
@ cnt ++
set suf = "$argv[$cnt]"
set donext = 0
endif
endif
if ($donext && "$argv[$cnt]" == "-base") then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need template volume after -base"
goto END
else
@ cnt ++
set ref_in = "$argv[$cnt]"
set donext = 0
endif
endif
if ($donext && ( "$argv[$cnt]" == "-pad_base" || "$argv[$cnt]" == "-pad_input" ) ) then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need value after -pad_base (or -pad_input) "
goto END
else
@ cnt ++
set pd = "$argv[$cnt]"
if ($pd > 100 || $pd < 0) then
echo "Error: -base_pad should be between 0 and 100"
echo "I have $pd"
goto END
endif
set donext = 0
endif
endif
if ($donext && "$argv[$cnt]" == "-input") then
set pLoc = $cnt
if ($pLoc == $Narg) then
echo "Need volume after -input"
goto END
else
@ cnt ++
set anat_in = "$argv[$cnt]"
set donext = 0
endif
endif
if ($donext && "$argv[$cnt]" == "-keep_tmp") then
set KeepTmp = 1;
set donext = 0
endif
if ($donext && "$argv[$cnt]" == "-clean") then
rm -f ${tmp}_*.* >& /dev/null
goto END
set donext = 0
endif
if ($donext && "$argv[$cnt]" == "-verb") then
set verb = 1;
set donext = 0
endif
if ($donext == 1) then
echo "Error: Option or parameter '$argv[$cnt]' not understood"
goto END
endif
@ cnt ++
end
if ("$anat_in" == "" || "$ref_in" == "") then
echo "Error: Need both -base and -input parameters"
goto END
endif
set base_parsed = `@parse_afni_name $ref_in`
set anat_parsed = `@parse_afni_name $anat_in`
set base_o = $base_parsed[2]
set anat_o = $anat_parsed[2]
set base = $base_o
set anat = $anat_o
goto MOVE_BASE
ALIGN_CENTERS:
#Not needed, assuming there are no clipping problems.
#Better do this for all oblique sets outside of this script if necessary
rm -f ${anat}_shft*.???? >& /dev/null
@Align_Centers -base ${base}+orig -dset ${anat}+orig
set anat = ${anat}_shft
MOVE_BASE:
#Try different large angles, moving base because it is much larger
rm -f ${tmp}_${base}_rs*.???? >& /dev/null
3dresample -dxyz $dxyz_list[1] $dxyz_list[2] $dxyz_list[3] -inset ${base_parsed[1]}/${base}+orig -prefix ./${tmp}_${base}_rs
set base = ${tmp}_${base}_rs
rm -f ${tmp}_${anat}_rs*.???? >& /dev/null
3dresample -inset ${anat_parsed[1]}/${anat}+orig -prefix ./${tmp}_${anat}_rs -master ${base}+orig
set anat = ${tmp}_${anat}_rs
set anat_i = $anat
set base_i = $base
set olst = ''
set ilst = ''
set RMSmin = 9000;
set ANGmin = -1;
set CNTmin = 0;
set cnt = 0
foreach ang ( 0 20 40 )
echo "Trying ${ang}R 0A 0I...(`date`)"
@ cnt ++
set anat = ${anat_i}
rm -f ${base}_R${ang}*.???? >& /dev/null
3drotate -rotate ${ang}R 0A 0I -prefix ${base}_R${ang} ${base}+orig
set ilst = ($ilst ${base}_R${ang}+orig)
rm -f ${anat}_b${ang}_aa*.???? >& /dev/null
3dvolreg -twopass -cubic -prefix ${anat}_b${ang}_aa -dfile ${anat}_b${ang}_aa.log -base ${base}_R${ang}+orig ${anat}+orig.
#get the rms result
set dfile = `cat ${anat}_b${ang}_aa.log`
if ( `ccalc -form int -expr step"($RMSmin - $dfile[9])"` ) then
set ANGmin = $ang;
set RMSmin = $dfile[9];
set CNTmin = $cnt;
endif
set anat = ${anat}_b${ang}_aa
set olst = ($olst ${anat}+orig)
end
echo $olst
echo $ilst
echo "`date`:Optimal registration found at test number $CNTmin ANGmin = $ANGmin, RMSmin = $RMSmin" >> ${tmp}_LeLog.txt
tail ${tmp}_LeLog.txt
goto THE_FINAL_GIG
STORE_RESULTS:
#Not necessary unless you are debugging
rm -f ${anat_i}_aa_all*.???? ${base_i}_all*.???? >& /dev/null
3dTcat -prefix ${anat_i}_aa_all $olst
3dTcat -prefix ${base_i}_all $ilst
#afni -noplugouts ${base}+orig.HEAD ${anat_i}_aa_all+orig.HEAD ${base_i}_all+orig.HEAD &
APPLY_BEST_MOVE_BASE:
#get inverse transform of best alignment
cat_matvec "$olst[$CNTmin]"::VOLREG_MATVEC_000000 -I > ${tmp}_iVR_MV.1D
#apply it to the pre-moved base
3drotate -matvec_dicom ${tmp}_iVR_MV.1D -prefix `@GetAfniPrefix $ilst[$CNTmin]`_alnd_${anat_i} $ilst[$CNTmin]
echo "Run: afni -noplugouts `@GetAfniPrefix $ilst[$CNTmin]`_alnd_${anat_i}+orig.HEAD ${anat_i}+orig.HEAD"
THE_FINAL_GIG:
if ("$suf" == "") then
set suf = "_alnd_${anat_o}"
endif
#apply full xform to initial base
cat_matvec ${base}_R${ANGmin}+orig::ROTATE_MATVEC_000000 "$olst[$CNTmin]"::VOLREG_MATVEC_000000 -I > ${tmp}_BaseToAnat.1D
#apply to initial base
rm -f ${base_o}${suf}*.???? >& /dev/null
3drotate -matvec_dicom ${tmp}_BaseToAnat.1D -prefix ./${base_o}${suf} ${base_parsed[1]}/${base_o}+orig
echo "Run: afni -noplugouts ${base_o}_alnd_${anat_o}+orig.HEAD ${anat_o}+orig.HEAD"
goto CLEANUP
FINE_REG:
#Cannot do finer registration because scaling is not allowed, partial coverage...
if (0) then
3dWarpDrive -coarserot -base `@GetAfniPrefix $ilst[$CNTmin]`_alnd_${anat_i}+orig\
-prefix Anat_alnd_MovedBase -affine_general \
-input ${anat_i}+orig -quintic -verb
endif
CLEANUP:
if ($KeepTmp == 0) then
rm -f ${tmp}_*.* >& /dev/null
endif
goto END
HELP:
if (1) then
echo "Usage 1: A script to align a full coverage T1 weighted non-oblique dataset"
echo " to match a partial coverage T1 weighted non-oblique dataset "
echo " Alignment is done with a rotation and shift (6 parameters) transform"
echo " only."
echo ""
echo " Script is still in testing phase"
echo ""
echo " `basename $0` [options] <-base FullCoverageT1> <-input PartialCoverageObliqueT1>"
echo " Mandatory parameters:"
echo " -base FullCoverageT1: Reference anatomical full coverage volume."
echo ""
echo " -input PartialCoverageObliqueT1: The name says it all."
echo ""
echo " Optional parameters:"
#echo " -pad_base MM : Pad the base dset by MM mm in each directions."
#echo " That is needed to make sure that datasets"
#echo " requiring wild rotations do not get cropped."
#echo " Default is MM = $pd_default."
#echo " If your output dataset is clipped, try increasing"
#echo " MM to `ccalc -eval $pd_default+10` or `ccalc -eval $pd_default+20`."
echo " -suffix SUF : Output dataset name is formed by adding SUF to"
echo " the prefix of the base dataset."
echo " The default suffix is _alnd_PartialCoverageObliqueT1"
echo " -keep_tmp : Keep temporary files."
echo " -clean : Clean all temp files, likely left from -keep_tmp"
echo " option then exit."
echo " -dxyz MM : Cubic voxel size of output DSET in TLRC"
echo " space Default MM is 1. If you do not"
echo " want your output voxels to be cubic"
echo " Then use the -dx, -dy, -dz options below."
echo " -dx MX : Size of voxel in the x direction"
echo " (Right-Left). Default is 1mm."
echo " -dy MY : Size of voxel in the y direction"
echo " (Anterior-Posterior). Default is 1mm."
echo " -dz MZ : Size of voxel in the z direction."
echo " (Inferior-Superior).Default is 1mm."
echo " Example:"
echo " `basename $0` -base ah_SurfVol+orig. -input ah_T1W_anat+orig."
echo ""
echo ""
echo "Written by Ziad S. Saad, for Ikuko (ziad@nih.gov)"
echo " SSCC/NIMH/NIH/DHHS"
echo ""
goto END
endif
END:
syntax highlighted by Code2HTML, v. 0.9.1