Skip to content

WIP: Automatically reshape FreeSurfer ico7 niftis #315

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 4 commits into from
Jun 16, 2015
Merged

WIP: Automatically reshape FreeSurfer ico7 niftis #315

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
d26215e
NF: Automatically reshape FreeSurfer ico7 niftis
effigies Jun 10, 2015
c6526f7
TEST: FreeSurfer ico7 surface NIfTI hack
effigies Jun 10, 2015
3477a13
TEST: Update nitest-freesurfer submodule
effigies Jun 15, 2015
51cc7de
STY: Remove blank lines from test_freesurfer_ico7_hack
effigies Jun 16, 2015
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion nibabel-data/nitest-freesurfer
Open in desktop
Submodule nitest-freesurfer updated from 011086 to df9ccc
38 changes: 29 additions & 9 deletions nibabel/nifti1.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -698,16 +698,26 @@ def get_data_shape(self):
Allows for freesurfer hack for large vectors described in
https://github.com/nipy/nibabel/issues/100 and
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/save_nifti.m?spec=svn5022&r=5022#77

Allows for freesurfer hack for 7th order icosahedron surface described
in
https://github.com/nipy/nibabel/issues/309
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/load_nifti.m?r=8776#86
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/save_nifti.m?r=8776#50
'''
shape = super(Nifti1Header, self).get_data_shape()
# Apply freesurfer hack for vector
if shape != (-1, 1, 1): # Normal case
if shape == (-1, 1, 1):
vec_len = int(self._structarr['glmin'])
if vec_len == 0:
raise HeaderDataError('-1 in dim[1] but 0 in glmin; '
'inconsistent freesurfer type header?')
return (vec_len, 1, 1)
# Apply freesurfer hack for ico7 surface
elif shape == (27307, 1, 6):
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do I understand right - that this is one particular image ('ico7')? So the hack only applies to that image?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

To surfaces with that number of vertices. I can't find a good reference to the arguments as to why they do this[0], but a reasonably common FreeSurfer surface is a warped icosahedron, and in particular they recommend using the 7th-order icosahedron. This is the number of vertices used in their fsaverage group template for comparisons across subjects.

Their MATLAB code indicates that this specific surface shape is hard-coded:
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/save_nifti.m?r=8776#50
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/load_nifti.m?r=8776#86

Edit: Nipype's documentation indicates that icosahedra up to the 7th order are supported. 6th order and under have few enough vertices to not require a hack, which explains why it's so specific to ico7.

[0] I believe the use of the icosahedron is to allow for near constant neighbor-neighbor distances.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I noticed this:

"""
By default, mri_surf2surf will save the output as multiple 'slices'; has no effect for paint/w output format. For ico, the output will appear to be a 'volume' with Nv/R colums, 1 row, R slices and Nf frames, where Nv is the number of vertices on the surface. For icosahedrons, R=6. For others, R will be the prime factor of Nv closest to 6. Reshaping is for logistical purposes (eg, in the analyze format the size of a dimension cannot exceed 2^15). Use this flag to prevent this behavior. This has no effect when the output type is paint.
"""

at https://surfer.nmr.mgh.harvard.edu/fswiki/mri_surf2surf - comment for --noreshape.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The easiest way to clarify the situation is to simply try all of the icosahedra:

$ for ORDER in {1..7}; do
    mri_surf2surf --hemi rh --srcsubject fsaverage \
        --srcsurfval $SUBJECTS_DIR/fsaverage/surf/rh.orig.avg.area.mgh \
        --trgsubject ico --trgsurfval rh.ico$ORDER.nii \
        --trgicoorder $ORDER &> /dev/null \
    && python -c "from nibabel import load;\
          print ($ORDER, load('rh.ico$ORDER.nii').header._structarr['dim'][1:4])" \
    && rm rh.ico$ORDER.nii
  done 
(1, array([42,  1,  1], dtype=int16))
(2, array([162,   1,   1], dtype=int16))
(3, array([642,   1,   1], dtype=int16))
(4, array([2562,    1,    1], dtype=int16))
(5, array([10242,     1,     1], dtype=int16))
(6, array([-1,  1,  1], dtype=int16))
(7, array([27307,     1,     6], dtype=int16))

So whatever their ideal of what ought to be done for icosahedra with >2^15 vertices, in practice, they use the (27307, 1, 6) shape for ico7, the large vector hack for ico6 and just (x, 1, 1) for lesser icosahedra. (ico8 does not work.)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there any possibility that someone will try to use this with a volume with some other number of vertices? Would they expect it to work?

return (163842, 1, 1)
else: # Normal case
return shape
vec_len = int(self._structarr['glmin'])
if vec_len == 0:
raise HeaderDataError('-1 in dim[1] but 0 in glmin; inconsistent '
'freesurfer type header?')
return (vec_len, 1, 1)

def set_data_shape(self, shape):
''' Set shape of data
Expand All @@ -725,12 +735,22 @@ def set_data_shape(self, shape):
Applies freesurfer hack for large vectors described in
https://github.com/nipy/nibabel/issues/100 and
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/save_nifti.m?spec=svn5022&r=5022#77

Allows for freesurfer hack for 7th order icosahedron surface described
in
https://github.com/nipy/nibabel/issues/309
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/load_nifti.m?r=8776#86
https://code.google.com/p/fieldtrip/source/browse/trunk/external/freesurfer/save_nifti.m?r=8776#50
'''
# Apply freesurfer hack for vector
hdr = self._structarr
shape = tuple(shape)
if (len(shape) == 3 and shape[1:] == (1, 1) and
shape[0] > np.iinfo(hdr['dim'].dtype.base).max): # Freesurfer case

# Apply freesurfer hack for ico7 surface
if shape == (163842, 1, 1):
shape = (27307, 1, 6)
# Apply freesurfer hack for vector
elif (len(shape) == 3 and shape[1:] == (1, 1) and
shape[0] > np.iinfo(hdr['dim'].dtype.base).max):
try:
hdr['glmin'] = shape[0]
except OverflowError:
Expand Down
30 changes: 29 additions & 1 deletion nibabel/tests/test_nifti1.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -24,8 +24,11 @@
Nifti1Pair, Nifti1Extension, Nifti1Extensions,
data_type_codes, extension_codes, slice_order_codes)

from ..freesurfer import load as mghload

from .test_arraywriters import rt_err_estimate, IUINT_TYPES
from .test_helpers import bytesio_filemap, bytesio_round_trip
from .nibabel_data import get_nibabel_data, needs_nibabel_data

from numpy.testing import (assert_array_equal, assert_array_almost_equal,
assert_almost_equal)
Expand Down Expand Up @@ -242,7 +245,7 @@ def test_magic_offset_checks(self):
'file nifti1; setting to minimum value '
'of ' + str(hdr.single_vox_offset))

def test_freesurfer_hack(self):
def test_freesurfer_large_vector_hack(self):
# For large vector images, Freesurfer appears to set dim[1] to -1 and
# then use glmin for the vector length (an i4)
HC = self.header_class
Expand Down Expand Up @@ -288,6 +291,31 @@ def test_freesurfer_hack(self):
hdr.set_data_shape(constructor(shape))
assert_equal(hdr.get_data_shape(), shape)

@needs_nibabel_data('nitest-freesurfer')
def test_freesurfer_ico7_hack(self):
HC = self.header_class
hdr = HC()
# Test that using ico7 shape automatically uses factored dimensions
hdr.set_data_shape((163842, 1, 1))
assert_array_equal(hdr._structarr['dim'][1:4], np.array([27307, 1, 6]))
# Test consistency of data in .mgh and mri_convert produced .nii
nitest_path = os.path.join(get_nibabel_data(), 'nitest-freesurfer')
mgh = mghload(os.path.join(nitest_path, 'fsaverage', 'surf',
'lh.orig.avg.area.mgh'))
nii = load(os.path.join(nitest_path, 'derivative', 'fsaverage', 'surf',
'lh.orig.avg.area.nii'))
assert_equal(mgh.shape, nii.shape)
assert_array_equal(mgh.get_data(), nii.get_data())
assert_array_equal(nii.header._structarr['dim'][1:4],
np.array([27307, 1, 6]))
# Test writing produces consistent nii files
with InTemporaryDirectory():
nii.to_filename('test.nii')
nii2 = load('test.nii')
assert_equal(nii.shape, nii2.shape)
assert_array_equal(nii.get_data(), nii2.get_data())
assert_array_equal(nii.get_affine(), nii2.get_affine())

def test_qform_sform(self):
HC = self.header_class
hdr = HC()
Expand Down
6 changes: 5 additions & 1 deletion nibabel/tests/test_nifti2.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,11 @@ class _Nifti2Mixin(object):
sizeof_hdr = Nifti2Header.sizeof_hdr
quat_dtype = np.float64

def test_freesurfer_hack(self):
def test_freesurfer_large_vector_hack(self):
# Disable this check
pass

def test_freesurfer_ico7_hack(self):
# Disable this check
pass

Expand Down