Thanks to help from Wolfgang and Yingli (and Bob too), I've put together two scripts to do the rotation of the gradient vectors. The first script, I'll include in this message and a second script in the following message. The first script rotates all the gradient vectors based on the rotation of a single-brick dataset (or the first brick of a dataset). The script used 3dvolreg rotation history in the header if it is available. If not, it will try 3dWarpDrive instead.
If you save it as it is named here, you can use it like this:
tcsh @rotatevectorsallx1 volreg_dset gradients.1D
where volreg_dset is the name of your volume registered or 3dWarpDrive aligned dataset and gradients.1D is the name of the gradient vector file to be transformed. The output gradient vector file is named rotatedgrads.1D.
#!/bin/tcsh
# @rotatevectorsallx1
#
# rotate vectors from tensor.1D based on first brick of a volume
# registered dataset
#
# use the matrix vector rotation fields from the header attributes
# if 3dvolreg was used, they will be named VOLREG_MATVEC_nnnnnn
# if 3dWarpDrive was used, they will be named
WARPDRIVE_MATVEC_FOR_nnnnnn
if ("$1" == '' || "$1" == '-help' || "$1" == '-h') then
goto HELP
endif
if ("$#" < 2) then
echo "Not enough parameters"
goto END
endif
set dset = $1
set tensor1D = $2
1dcat $2 > rotatedgrads.1D
set ngrads = `cat rotatedgrads.1D|wc -l`
echo "Number of gradients is $ngrads"
@ ngrads--
echo "Dataset: $dset"
echo "Gradient vectors: $tensor1D"
echo "Rotating vectors..."
rm -f rotatedgrads.1D
# put the matrix vector on three rows to be compatible with Vecwarp
( cat_matvec "$dset"'::VOLREG_MATVEC_000000' > tempmatvec.1D ) >
& /dev/null
# if there are no attributes from 3dvolreg in header, try 3dWarpDrive
rotations
set check = `wc tempmatvec.1D`
if ($check[1] == 0) then
echo "Looking for WARPDRIVE attribute instead"
( cat_matvec "$dset"'::WARPDRIVE_ROTMAT_FOR_000000' >
tempmatvec0.1D ) > & /dev/null
set check = `wc tempmatvec0.1D`
if ($check[1] == 0) then
echo "No rotation matrix in header"
goto END
endif
else
# replace the translation vector with 0s
1dcat tempmatvec.1D'[0..2]' '1D:3@0' > tempmatvec0.1D
endif
foreach grad ( `count 0 $ngrads`)
# get the original vector for the current gradient
1dcat "$tensor1D{$grad}" > tempvector.1D
# rotate the vector for the current gradient
Vecwarp -matvec tempmatvec0.1D -backward -input tempvector.1D -
force -output temprotgrad.1D > & /dev/null
# append the rotated gradients to a single file
cat temprotgrad.1D >> rotatedgrads.1D
rm temprotgrad.1D
rm tempvector.1D
end
rm tempmatvec*.1D
echo "Output gradient vector file is rotatedgrads.1D"
goto END
HELP:
echo "Usage: `basename $0` volreg_dset gradients.1D"
echo "Rotate DTI gradient vectors based on rotations from volume
registration"
echo "Current script assumes 3dvolreg or 3dWarpDrive corrected
dataset with"
echo "preference to 3dvolreg rotation parameters"
echo "Output file is rotatedgrads.1D"
echo ""
echo "Daniel Glen SSCC/NIMH/NIH"
echo ""
goto END
END:
