Hi Joe,

Sorry about the slow reply.

> what is vol1.mat and vol2.mat?

I intended vol1 and vol2 to mean the volume structures returned by 
spm_vol, e.g.
   vol1 = spm_vol('image1.img');
in which case vol1.mat is the 4-by-4 affine transformation matrix 
mapping from voxel to world/mm coordinates for this volume.

If two images are affinely (or rigidly) coregistered (e.g. from 
realign, or coreg) then the world coordinates correspond, so voxel 
correspondence from image 1 to image 2 is given by
   xyz2 = inv(vol2.mat) * vol1.mat * xyz1;

> after coreg using SPM, all I have is sn.mat, 
> which upon being opened, shows Affine matrix and a couple of log kind of
> structs. Am I missing something here? (using SPM5).

It's crucial to distinguish the different types of registration. My 
original code snippet (and the stuff about vol1 and vol2 above) is for 
rigid/affine registered images. If you have an sn.mat then this is 
from non-rigid spatial-normalisation using the Normalise button (or 
unified segmentation and spatial normalisation in SPM5) not from the 
Coreg button.

You can use the sn.mat to map from target to source coordinates (which 
usually means from the template to the original individual image). I 
originally suggested converting it to a deformation field first, but 
that isn't necessary.

If you want to go in the other direction, and you don't already have 
the inverse from seg_inv_sn.mat, then you will still need to use the 
Deformations utility to write out the inverse as a deformation field.

For the actual mechanics of doing this...

>> Maybe I should just code up something to do this

...I have now written some code which can handle the three cases: 
rigid/affine, spatial normalisation, deformation field. I think the 
code is pretty well documented (I hope you agree -- but please ask if 
anything is unclear!).

   http://www.cs.ucl.ac.uk/staff/gridgway/vbm/map_coords.m

For just the case of spatial normalisation, there is also a function 
in the SPM5 directory -- but it has a slight bug if the sn.mat comes 
from unified segmentation. A bug-fixed (I hope) version is available here:

   http://www.cs.ucl.ac.uk/staff/gridgway/vbm/get_orig_coord5.m

Note that get_orig_coord5 requires its input coordinates to be in 
world/mm space, and can return either world or voxel space results. 
map_coords usually expects its input to be in voxel space (except the 
special case mentioned in the help) and can return world and voxel 
space results. They give the same results after accounting for this; 
use whichever is more convenient for you.

Best,
Ged.