Dear Sarah and Dianne,
(and FSL and SPM mailing lists, since it's a bit of a cross-over
question; sorry to those of you who read both!)

Starting with the easy bit... The reference volume is usually the
target of the registration, for FLIRT this is -ref, for SPM's Coreg
GUI it's "reference", but for the spm_coreg function (at the MATLAB
command line) it's VG. With non-rigid (though see below), for SPM's
normalise it's "template" , for the unified segmentation it's the
tissue probability maps (TPMs), and for DARTEL it's (any one of) the
group-wise average templates.

For converting rigid/affine matrices, I'm afraid it depends how you're
doing the registration in SPM. If you're calling spm_coreg directly,
then as described in the help, the output relates the voxel
coordinates in the images via:
  VF.mat\spm_matrix(x)*VG.mat
or in other words, M=spm_matrix(x) maps from target world/mm
coordinates to source world/mm coordinates. Note that VF.mat\blah
means inv(VF.mat)*blah, which maps from (source) world/mm coordinates
to source voxel coordinates.

You could then convert M to a FLIRT matrix (from source to target
"scaled" coordinates) using worldmat2flirtmat in the package of
snippets here:
  http://www.nitrc.org/snippet/detail.php?type=package&id=1

If you run Coreg from the SPM GUI, things are a bit more complicated,
since it modifies the voxel-world mapping in the source image, VF.mat,
in-place. Noting that VF.mat\M*VG.mat can be written
inv(inv(M)*VF.mat)*VG.mat, we see that VF.mat gets premultiplied by
the inverse of the output from spm_coreg (this happens in either
spm8/config/spm_run_coreg_estimate.m or spm5/spm_config_coreg.m --
search for spm_get_space, which is the function used to get/set the
transformation in the NIfTI headers).

Two things to note about this. First, Sarah, if you've coregistered
one image twice (to two different targets, or using different
parameters) then unless you've done the registrations on two different
copies of the original source image, the second coreg will have
started from where the first finished, which probably isn't what you
want, if you're trying to compare them fairly. Second, you need to
have kept a copy of the original source image (or its VF.mat) before
the registration to recover the matrix M from the result after the
registration, if you want to convert this for use with rmsdiff. You
could then do:
  after = spm_get_space('after.nii')
  before = spm_get_space('before.nii')
  M = before / after % or before*inv(after)
You could then proceed as above.

If you're doing non-rigid (non-affine) spatial normalisation, then you
can't convert to a matrix, and probably can't use FSL's rmsdiff
(though FSL might have an equivalent for displacement fields... I'm
not sure...). Probably the easiest thing in that case would be to
compute RMS differences yourself in MATLAB/SPM. You should be able to
convert blah_sn.mat files (or blah_seg_sn.mat from unified
segmentation, and also u_blah flow fields from DARTEL) to deformation
fields (y_blah) using the Deformations utility. The y_ files are
five-dimensional NIfTI images with dimensions [X Y Z 1 3], where X, Y
and Z are the numbers of voxels in these dimensions of the target
image, and the three components in the fifth dimension (the fourth was
reserved for time) are the x, y, and z coordinates of where each
target voxel maps to in source world/mm space. You can therefore
compare two different mappings (to the same target; things will be
more complicated if you wanted to compare mappings to two different
targets...) by directly computing the RMS difference between these
fields. E.g. the following might work, though I haven't tested it for
typos etc...
  A = nifti('y_a.nii');
  B = nifti('y_b.nii');
  D = A.dat(:, :, :, 1, :) - B.dat(:, :, :, 1, :);
  % If the above throws an error, check y_a and y_b relate to the same target
  err2 = sum(D.^2, 5); % sum over 5th dimension
  % err2 is a three-dimensional image of squared errors (in mm^2)
  rms = sqrt(mean(err2(:)))

This would also let you compute maximum error in mm, etc. etc.
  max(sqrt(err2(:)))
  hist(sqrt(err2(:))) % plot histogram

Note that the above mean and max are over all voxels in the target; it
might make more sense to use only voxels within a mask (defined over
the same space as the target image). E.g. if your y_ deformation comes
from unified segmentation, you should be able to do something like:
  G = spm_vol(fullfile(spm('dir'), 'TPM', 'grey.nii'));
  G = spm_read_vols(G) > 0.5;
then e.g.
  rms = sqrt(mean(err2(G)))
to compute RMS error over just the voxels with more than 50% (prior)
probability of being GM.

I hope (some of?!) that helps!

Ged


2009/10/5 Dianne Patterson <[log in to unmask]>:
> Dear Group,
> Has anyone answered this question? What is the relationship between spm
> registration matrices and fsl registration matrices?
> I'm looking forward to hearing more about this.
> -Dianne
>
> On Fri, Oct 2, 2009 at 3:35 PM, Sarah Andersen <[log in to unmask]>
> wrote:
>>
>> Dear Group,
>> I just learned about the rmsdiff utility and would really like to use it
>> to compare two
>> different normalizations of the same image. From the usage it seems that I
>> need to input
>> the two matrix files, one from each of the normalizations. However, I did
>> these
>> normalizations in SPM and FSL doesn't seem to recognize the SPM output
>> matrix. Is there a
>> way to convert the SPM matrix into something that FSL can read?
>> Also is the reference volume the initial brain before either of the
>> normalizations?
>>
>> Thanks!
>> Sarah
>>
>
>
>
> --
> Dianne Patterson, Ph.D.
> [log in to unmask]
> University of Arizona
> SLHS 328
> 621-5105
>