Back to the main page.
Bug 2540 - ft_volumesegment: segmented mri incorrectly inherits units from input mri
| Status | CLOSED FIXED |
| Reported | 2014-04-16 11:34:00 +0200 |
| Modified | 2019-08-10 12:29:46 +0200 |
| Product: | FieldTrip |
| Component: | core |
| Version: | unspecified |
| Hardware: | PC |
| Operating System: | Windows |
| Importance: | P5 normal |
| Assigned to: | |
| URL: | |
| Tags: | |
| Depends on: | |
| Blocks: | |
| See also: |
Roemer van der Meij - 2014-04-16 11:34:32 +0200
Created attachment 613 downsampled mri The attached (downsampled) mri has units = 'mm', verifiable by with ft_sourceplot. If I convert the units to 'cm', then segment the volume, the resulting segmentedmri has units = 'mm', even though in reality it is in 'cm', verifiable with ft_sourceplot. % segment mri load('mri.mat') mri = ft_convert_units(mri,'cm'); cfg = []; cfg.write = 'no'; cfg.coordsys = 'ctf'; segmentedmri = ft_volumesegment(cfg, mri); cfg = []; cfg.funparameter = 'gray'; ft_sourceplot(cfg,segmentedmri);
Jan-Mathijs Schoffelen - 2014-04-16 11:43:16 +0200
Well, it seems that the units in the input MRI are not pushed through correctly, right? Can you look into this Roemer?
Roemer van der Meij - 2014-04-16 11:49:21 +0200
No, it's a bit more complicated I think. Before the MRI is passed to spm functions, it is transformed again to mm at line 279. This is done correctly (verified with ft_sourceplot). The units of this converted mri are then later (line 426) correctly passed through to the output segmentedmri, being mm's. However, for some reason, after all the spm segmentation, the units are actually in cm's again. (verified with ft_sourceplot). Perhaps it lies somewhere in the transform matrix not being correctly converted in any of these ft_convert_units steps.
Jan-Mathijs Schoffelen - 2014-04-16 11:52:01 +0200
sqrt(sum(mri.transform(1:3,1:3).^2)) will give a clue as to the units. If you have a 1 mm voxel-spaced volume, if it's in mm you should get 1, if it's in cm you should get 0.1
Roemer van der Meij - 2014-04-16 11:59:05 +0200
Ah, gotcha. At line 278 the input mri.transform is saved as original.transform (being cm in this case). This is later put back into the output segmentedmri. However, the units are not saved. The units that are put in the output segmentedmri are those of the converted mri (for spm processing), which is always in mm. I can save the original unit, and put it back later with the original transform. That should solve it. However, there's an isfield on mri.unit at line 426... There are cases when the input mri did not have units? Or is this a deprecated case? If not, perhaps we can do something like this. Check if input mri had units, if so, save transform, coordsys and unit at line 277 in original.transform etc. Then, if orginal was not empty, put back all 3 fields in the output segmentedmri. If original is empty, use these fields from the converted mri (generated at line 280).
Roemer van der Meij - 2014-04-16 12:08:29 +0200
I've done it slightly different then I suggested. I now save the original unit when it was present, and put it back later if it was. The units from the converted mri are now never used. This is valid if the following is always true: * out = ft_convert_unit(mri) will always have a field .unit, if the input had a field .unit. Would you agree JM? Then I'll commit.
Jan-Mathijs Schoffelen - 2014-04-16 12:15:29 +0200
if you call out=ft_convert_units(mri); the function will estimate the units from mri-structure, if it was missing before. does that answer your question?
Jan-Mathijs Schoffelen - 2014-04-16 12:26:54 +0200
I guess the cleanes fix would be: at the top of the function: if ~isfield(mri, 'unit') mri = ft_convert_units(mri); end inputunit = mri.unit; at the bottom of the function: segmentedmri = ft_convert_units(segmentedmri, inputunit);
Roemer van der Meij - 2014-04-16 12:29:19 +0200
Ah, I see now that my assertion is not true (sorry for the quick updates ;)). I see that ft_convert_units tries to automatically determine the units if they are not present. So it's possible they were not present before the call to ft_volumesegment. I changed it a bit. If no units are present, it will try to determine them automatically by mri = ft_convert_units(mri), and then use it later on. It seems pretty safe now, I'll commit the changes. (mid-air collision there, but we ended up with the same idea :)) I didn't add the conversion at the bottom, as the transform, coordsys and units are now saved and added back, so there is no need for it.
Roemer van der Meij - 2014-04-16 12:43:43 +0200
setting as fixed
Roemer van der Meij - 2014-04-16 12:44:16 +0200
setting as fixed (mid-air collision with myself??)