WIP: Automatically reshape FreeSurfer ico7 niftis

[effigies]

Closes #309.

Tried to follow in the footsteps of f53c168. Presents programmer with data object equivalent to an MGH, leaving header intact for writing to file.

Tests planned:

1. mri_convert a .mgh file in nibabel-data/nitest-freesurfer/fsaverage/surf to '.nii' and verify img.get_data() is equal
2. Create new image with ico7 shape, verify img.header._structarr['dim'][1:4] == (27307, 1, 6)

[effigies]

@matthew-brett How does one actually run the test_nifti1 tests?

$ python -c 'import nibabel.tests.test_nifti1 as t; t.TestNifti1Image()' 
Traceback (most recent call last):
  File "<string>", line 1, in <module>
  File "/usr/lib/python2.7/unittest/case.py", line 191, in __init__
    (self.__class__, methodName))
ValueError: no such test method in <class 'nibabel.tests.test_nifti1.TestNifti1Image'>: runTest

[matthew-brett]

Run the tests locally you mean?

You could start with git submodule update --init to get the test data and then nosetests nibabel\tests\test_nifti1.py from the nibabel root directory (containing setup.py).

[effigies]

Thanks. Wrote tests that pass locally, with change in https://bitbucket.org/nipy/nitest-freesurfer/pull-request/1/. If that PR goes through, I'll push the tests.

[matthew-brett]

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

[effigies]

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.

[matthew-brett]

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.

[effigies]

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.)

[matthew-brett]

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?

[matthew-brett]

Ben - I noticed your email to the freesurfer list, and I think it got no reply. What do you think we should do here?

[effigies]

I figured I'd at least wait until Monday, since I sent this email out on a Friday in the summer.

If we don't hear anything, my feeling is we should follow the lead of load/save_nifti.m and only handle these two special cases until there's evidence of other special cases in the wild.

Though the ico6 case is indicated in the documentation, being unable to produce it with the tool the documentation was written about makes me skeptical that they've bothered implementing it. (They recommend everybody use ico7, anyway.) We could still add the case, as it's an unlikely shape to appear by accident, but then that introduces the problem of saving files: Do we use the (-1, 1, 1) or the (6827, 1, 6)?

Assuming the decision is to restrict ourselves to the existing cases, is there anything you need from me regarding the nitest-freesurfer PR?

[matthew-brett]

Just needs tests - I guess these are waiting on the nitest-freesurfer PR

[effigies]

Yup. In case BitBucket didn't alert you, I updated the PR over there.

[matthew-brett]

Merged the data commit, thanks for that.

[effigies]

Seems not to be fetching the new nitest-freesurfer. Is that just a thing that requires time to sync? https://travis-ci.org/nipy/nibabel/jobs/66921373

[matthew-brett]

Ah no - you need to update the submodule. The submodule points to one exact commit, not a branch - https://matthew-brett.github.io/pydagogue/git_submodules.html

[matthew-brett]

Great - thanks very much for this - in it goes.