* Reads in a 3D dataset and a DICOM-based affine matrix
* Outputs the average and maximum displacement that the matrix produces
when applied to the edge voxels of the 3D dataset's automask.
* The motivation for this program was to check if two
affine transformation matrices are 'close' -- but of course,
you can use this program for anything else you like.
* Suppose you have two affine transformation matrices that
transform a dataset XX.nii to MNI space, stored in files
AA.aff12.1D and BB.aff12.1D
and they aren't identical but they are close. How close?
* If these matrices are from 3dAllineate (-1Dmatrix_save),
then each matrix transforms DICOM-order coordinates
in XX.nii to MNI-space.
* So Inverse(AA) transforms MNI-space to XX-space
* So Inverse(AA)*BB transforms MNI-space to MNI-space,
going back to XX-space via matrix Inverse(AA) and then forward
to MNI-space via BB.
* This program (3dmaxdisp) can compute the average and maximum
displacement of Inverse(AA)*BB over the edges of the MNI template,
which will give you the answer to 'How close?' are the matrices.
If these displacements are on the order of a voxel size
(e.g., 1 mm), then the two matrices are for practical purposes
'close enough' (in Zhark's opinion).
* How to do this?
cat_matvec AA.aff12.1D -I BB.aff12.1D > AinvB.aff12.1D
3dmaxdisp -dset ~/abin/MNI152_2009_template_SSW.nii.gz'' -matrix AinvB.aff12.1D
* Results are sent to stdout, two numbers per row (average and maximum),
one row of output for each matrix row given. Usually you will want to
capture stdout to a file with '>' or '| tee', depending on your further plans.
-inset ddd }= The input dataset is 'ddd', which is used only to form
*OR* }= the mask over which the displacements will be computed.
-dset ddd }=
-matrix mmm = File 'mmm' has 12 numbers per row, which are assembled
into the 3x4 affine transformation matrix to be applied
to the coordinates of the voxels in the dataset mask.
* As a special case, you can use the word 'IDENTITY'
for the matrix filename, which should result in
a max displacement of 0 mm.
* If there is more than 1 row in 'mmm', then each row
is treated as a separate matrix, and the max displacement
will be computed separately for each matrix.
-verb = Print a few progress reports (to stderr).
Author: Zhark the Displacer (AKA Bob the Inverted) -- June 2021