The most significant is subclassingMGHHeader fromLabeledWrapStruct. It appeared to have been a reimplementation ofLabeledWrapStruct in many places without theendianness option. I have simply raised an exception on attempts to instantiate as little-endian.

The next most significant change is moving fromMdc (matrix of direction cosines) andPxyz_c (center point in world coordinates) to the more explicitx/y/z/c_ras column vectors. These are consistent with most FreeSurfer literature and program outputs, and correctly encode the shape of the direction cosines as column-major rather than row-major (needing to be transposed when reading/writing).

Somewhat confusingly, there are two ways of referring to the voxel and world coordinates, andXYZ means world in one, and voxel in the other, and in both cases the alternative includes "R" and "S".

Voxel	World
ColumnRowSlice	XYZ
XYZ	RightAnteriorSuperior

As CRS/XYZ appears to be FreeSurfer-unique, I'd prefer to settle on XYZ/RAS, and perhaps consider non-XYZ variable names when they can reduce confusion.

TODO:

Update tests with new header fields
Annotate header dtype with descriptions of fields (possibly further updating field names, shapes)
Improve footer handling (not targeting full features, but access to TR would be good)
Fixget_data_shape handling (seeENH: Add image slicing #550)

effigies force-pushed therf/mgh branch 3 times, most recently from5241b1b to832f584Compare

October 28, 2017 14:43

effigies changed the title~~[WIP] RF: Simplify MGHImage~~RF: Simplify MGHImage and add footer fields

Oct 28, 2017

Copy link

coveralls commentedOct 28, 2017

Coverage increased (+0.02%) to 96.27% when pulling8870068 on effigies:rf/mgh intob1bf5e7 on nipy:master.

Copy link

codecov-io commentedOct 28, 2017•
edited
Loading

Codecov Report

Merging#569 intomaster willincrease coverage by0.09%.
The diff coverage is97.83%.

@@            Coverage Diff             @@##           master     #569      +/-   ##==========================================+ Coverage   94.32%   94.41%   +0.09%==========================================  Files         177      177                Lines       24697    24884     +187       Branches     2640     2656      +16     ==========================================+ Hits        23295    23495     +200+ Misses        925      913      -12+ Partials      477      476       -1

Impacted Files	Coverage Δ
nibabel/spatialimages.py	`96.05% <100%> (ø)`	⬆️
nibabel/tests/test_spatialimages.py	`99.7% <100%> (ø)`	⬆️
nibabel/freesurfer/tests/test_mghformat.py	`98.44% <97.31%> (-1.56%)`	⬇️
nibabel/freesurfer/mghformat.py	`95.53% <98.42%> (+7.03%)`	⬆️

Continue to review full report at Codecov.

Legend -Click here to learn more
Δ = absolute <relative> (impact),ø = not affected,? = missing data
Powered byCodecov. Last update2817d1b...8d720fe. Read thecomment docs.

Copy link

coveralls commentedOct 30, 2017

Coverage increased (+0.06%) to 96.311% when pulling21cbbcf on effigies:rf/mgh intob1bf5e7 on nipy:master.

Copy link

MemberAuthor

effigies commentedOct 30, 2017

This is ready for review, if anybody's got a few minutes.

matthew-brett reviewed

Oct 30, 2017

View reviewed changes

Copy link

Member

matthew-brett left a comment

Choose a reason for hiding this comment

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

Main question - is it really worth renaming Mdc and the other fields?

nibabel/freesurfer/mghformat.py Outdated


		def__init__(self,
		binaryblock=None,
		endianness='>',

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

Does it make sense to allow passing this flag? It also changes the API.

Copy link

MemberAuthor

effigiesOct 30, 2017

Choose a reason for hiding this comment

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

I did this because of places whereLabeledWrapStruct calls__init__ assuming the type signature of the base class. I can instead make sure to override any methods that make such calls, just figured this would be a little shorter.

Copy link

Member

matthew-brettOct 31, 2017

Choose a reason for hiding this comment

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

I played with the idea of making WrapStruct have a parent or similar that had forced endianness, but I seem to remember it got a bit messy.

nibabel/freesurfer/mghformat.py

		returnklass(hdr_str+ftr_str,check=check)

		defget_affine(self):
		''' Get the affine transform from the header information.

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

Blank line after first line (sorry, I know that wasn't present before either).

nibabel/freesurfer/mghformat.py Outdated

		MdcD=np.hstack((hdr['x_ras'],hdr['y_ras'],hdr['z_ras']))*hdr['voxelsize']
		vol_center=MdcD.dot(hdr['dims'][:3].reshape(-1,1))/2
		affine[:3, :3]=MdcD
		affine[:3, [3]]=hdr['c_ras']-vol_center

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

3 rather than[3] ?

Copy link

MemberAuthor

effigiesOct 30, 2017

Choose a reason for hiding this comment

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

The RHS is a column, so the[3] keeps the LHS expecting a column. Is there a trick I don't know here?

nibabel/freesurfer/mghformat.py

		hdr=self._header_data
		zooms=hdr['delta']
		returntuple(zooms[:])
		# Do not return time zoom (TR) if 3D image

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

Maybe worth a note in the docstring about the time zoom and 3D?

nibabel/freesurfer/mghformat.py Outdated

		('dims','>i4', (4,)),# 4; width, height, depth, nframes
		('type','>i4'),# 20; data type
		('dof','>i4'),# 24; degrees of freedom
		('ras_good','>i2'),# 28; *_ras fields valid

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

I guess these renamings are API changes - our our Freesurfer colleagues happy with these name changes? Any way to offer the old names as alternatives for a while?

Copy link

MemberAuthor

effigiesOct 30, 2017

Choose a reason for hiding this comment

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

Yeah, I can re-add a_header_data structure (which I've replaced withLabeledWrapStruct._structarr) that will expose the old fields and raise aDeprecationWarning.

nibabel/freesurfer/mghformat.py Outdated

		('dof','>i4'),# 24; degrees of freedom
		('ras_good','>i2'),# 28; *_ras fields valid
		('voxelsize','>f4', (3,)),# 30; zooms (X, Y, Z)
		('x_ras','>f4', (3,1)),# 42; X direction cosine column

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

(3, 1) rather than(3,)?

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

Is this split into three columns really useful? Compared to the native 3, 3 size?

Copy link

MemberAuthor

effigiesOct 30, 2017

Choose a reason for hiding this comment

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

By including all three columns as a single, row-orderedMdc matrix, we're exposing a transposed matrix. It was being transposed correctly in and out, but if someone does decide to use it directly, they need to know this. I prefer explicit to silently misleading...

I could use(3,) if that's important. But it will mean usersknowing that they're columns and not row vectors (can add to the docs).

For context, this is the comparable portion of the C data structure, which is just 16 individualfloat fields:https://github.com/freesurfer/freesurfer/blob/5367eec84621b331271e97f4e1f1502c3f6d0d87/include/mri.h#L184-L191

Copy link

Member

matthew-brettOct 31, 2017

Choose a reason for hiding this comment

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

OK - fine to leave as is (3,1) etc.

nibabel/freesurfer/mghformat.py Outdated

		hdr=self._structarr
		MdcD=np.hstack((hdr['x_ras'],hdr['y_ras'],hdr['z_ras']))*hdr['voxelsize']
		vol_center=MdcD.dot(hdr['dims'][:3].reshape(-1,1))/2
		affine[:3, :3]=MdcD

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

return nib.affines.from_matvec(MdcD, hdr['c_ras'] - vol_center)

nibabel/freesurfer/mghformat.py Outdated


		Ignores byte order; always big endian
		'''
		structarr=super(MGHHeader,klass).default_structarr()

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

Here might be worth checking and error if endian not big.

nibabel/freesurfer/mghformat.py Outdated


		# Assign after we've had a chance to raise exceptions
		hdr['voxelsize']=voxelsize
		hdr['x_ras']=Mdc[:, [0]]

Copy link

Member

matthew-brettOct 30, 2017

Choose a reason for hiding this comment

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

hdr['x_ras'], hdr['y_ras'], hdr['z_ras'] = Mdc.T

Copy link

MemberAuthor

effigies commentedOct 30, 2017

Re: Renaming, the biggest motivation is that a lot of the field names seem only to align with one slideshow that documents the spaces/transforms, and not at all with either the FreeSurfer source code or the outputs ofmri_info or other FreeSurfer tools, let alone how anybody outside FreeSurfer refers to these spaces.

So my strong preference is to move to a more consistent naming scheme (I'm open to alternative suggestions), but I'm also perfectly happy to provide a backwards-compatible interface to access the data using the old field names to get back the old shapes.

matthew-brett reviewed

Oct 31, 2017

View reviewed changes

Copy link

Member

matthew-brett left a comment

Choose a reason for hiding this comment

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

I guess it would be good to get more feedback from Freesurferers out there.

nibabel/freesurfer/mghformat.py Outdated

		('dof','>i4'),# 24; degrees of freedom
		('ras_good','>i2'),# 28; *_ras fields valid
		('voxelsize','>f4', (3,)),# 30; zooms (X, Y, Z)
		('x_ras','>f4', (3,1)),# 42; X direction cosine column

Copy link

Member

matthew-brettOct 31, 2017

Choose a reason for hiding this comment

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

OK - fine to leave as is (3,1) etc.

nibabel/freesurfer/mghformat.py Outdated


		def__init__(self,
		binaryblock=None,
		endianness='>',

Copy link

Member

matthew-brettOct 31, 2017

Choose a reason for hiding this comment

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

I played with the idea of making WrapStruct have a parent or similar that had forced endianness, but I seem to remember it got a bit messy.

Copy link

MemberAuthor

effigies commentedNov 1, 2017

Okay. I've pinged the Python neuroimaging list and the FreeSurfer support list. Hopefully we can find somebody willing to share their opinions on the Internet. :-)

Copy link

Member

mwaskom commentedNov 1, 2017•
edited
Loading

I don't have a strong opinion here, but I would guess it's more likely that a user has encountered the fscoordinates slide deck than read through the freesurfer source code.

Copy link

Member

satra commentedNov 1, 2017

@effigies - perhaps chatting with doug would be a good idea as well. also folks at corticometrics may have some suggestions as well.

however, i do not have a strong opinion either. since i rarely do much with those files than reading into an array, doing some processing, and writing. or just visualization.

Copy link

MemberAuthor

effigies commentedNov 1, 2017•
edited
Loading

@satra Can I expect them to respond to the post on the FreeSurfer list, or is there a preferred way of contacting them more directly?

Copy link

Member

satra commentedNov 1, 2017

@effigies - i hadn't seen that the mail was on the freesurfer list. hopefully they respond - otherwise it's your call :)

Copy link

MemberAuthor

effigies commentedNov 7, 2017

Okay, it sounds like there's a slight preference from@mwaskom to follow the FS documentary slideshow and a preference (of unknown strength) from@matthew-brett to avoid changing names without a Very Good Reason. I'll defer to this for now and revertMdc andPxyz_c, because I don't think it's critical to my more pressing goals here of getting zooms for 4D images to include TR (to match the behavior of Analyze-like images) and enforcing the 3/4D constraint.

However, I do think I would like to move to more "standard" FreeSurfer naming, because this construct ("volume geometry", or dimensions + x/y/z/c_ras) shows up in several places, includingm3z warp files, LTA source and destination geometries (see#565), and in a text appendage to surface files:

nibabel/nibabel/freesurfer/io.py

Lines 50 to 73 in2139ce0

	def_read_volume_info(fobj):
	"""Helper for reading the footer from a surface file."""
	volume_info=OrderedDict()
	head=np.fromfile(fobj,'>i4',1)
	ifnotnp.array_equal(head, [20]):# Read two bytes more
	head=np.concatenate([head,np.fromfile(fobj,'>i4',2)])
	ifnotnp.array_equal(head, [2,0,20]):
	warnings.warn("Unknown extension code.")
	returnvolume_info

	volume_info['head']=head
	forkeyin ['valid','filename','volume','voxelsize','xras','yras',
	'zras','cras']:
	pair=fobj.readline().decode('utf-8').split('=')
	ifpair[0].strip()!=keyorlen(pair)!=2:
	raiseIOError('Error parsing volume info.')
	ifkeyin ('valid','filename'):
	volume_info[key]=pair[1].strip()
	elifkey=='volume':
	volume_info[key]=np.array(pair[1].split()).astype(int)
	else:
	volume_info[key]=np.array(pair[1].split()).astype(float)
	# Ignore the rest
	returnvolume_info

nibabel/nibabel/freesurfer/io.py

Lines 462 to 488 in2139ce0

	def_serialize_volume_info(volume_info):
	"""Helper for serializing the volume info."""
	keys= ['head','valid','filename','volume','voxelsize','xras','yras',
	'zras','cras']
	diff=set(volume_info.keys()).difference(keys)
	iflen(diff)>0:
	raiseValueError('Invalid volume info: %s.'%diff.pop())

	strings=list()
	forkeyinkeys:
	ifkey=='head':
	ifnot (np.array_equal(volume_info[key], [20])ornp.array_equal(
	volume_info[key], [2,0,20])):
	warnings.warn("Unknown extension code.")
	strings.append(np.array(volume_info[key],dtype='>i4').tostring())
	elifkeyin ('valid','filename'):
	val=volume_info[key]
	strings.append('{0} = {1}\n'.format(key,val).encode('utf-8'))
	elifkey=='volume':
	val=volume_info[key]
	strings.append('{0} = {1} {2} {3}\n'.format(
	key,val[0],val[1],val[2]).encode('utf-8'))
	else:
	val=volume_info[key]
	strings.append('{0} = {1:0.10g} {2:0.10g} {3:0.10g}\n'.format(
	key.ljust(6),val[0],val[1],val[2]).encode('utf-8'))
	returnb''.join(strings)

But I think it does make sense to postpone the discussion of field naming to an actual attempt to provide a common interface to volume geometries (which is actually re-implemented in several FreeSurfer data structures).

Copy link

coveralls commentedNov 9, 2017

Coverage increased (+0.004%) to 96.257% when pulling1bad727 on effigies:rf/mgh intob1bf5e7 on nipy:master.

Copy link

Member

matthew-brett commentedNov 10, 2017

Chris - sorry - I should have said - I don't care too much about the name changes, especially if there's a deprecation route. I should have said before, but if you think these names are more canonical, that's OK - we can always explain in the docstring and comments what the relationship is to the slide deck.

Copy link

MemberAuthor

effigies commentedNov 10, 2017

Okay. I think I'll still push it off for now. If/when we add support for LTA/M3Z etc, it will be time enough to consider a consistent naming scheme, rather than deprecating the old one now and finding later that we would have preferred something slightly different.

I will add docs though noting thatMdc is a transposed matrix. (It's a little annoying that I can't specify row/column order in the dtype.)

effigies force-pushed therf/mgh branch from1bad727 to79a5876Compare

November 23, 2017 01:17

Copy link

coveralls commentedNov 23, 2017

Coverage increased (+0.003%) to 96.258% when pulling79a5876 on effigies:rf/mgh intod420dcf on nipy:master.

effigies added8 commits

November 23, 2017 09:09

FIX/ENH: get/set_data_shape more robustly

bb56000

RF: Subclass MGHHeader from LabeledWrapStruct

d7b472a

RF: Remove one-line private method

415ee10

FIX: goodRASFlag check, update

4e520ac

RF: Remove unused _load_cache

a83f80b

FIX: Label x/y/z_ras correctly

ae5c085

RF: Make direction cosines explicit, not transposed

6b6e512

RF: Pxyz_c to c_ras

971c98a

effigies added4 commits

November 23, 2017 09:09

ENH: Check for expected byte order

16c8bb6

ENH: Re-add _header_data with deprecation warning

743db3b

FIX: Return on get, not on set

495bbce

REVERT renamed/reshaped header fields

821b6b6

effigies force-pushed therf/mgh branch from79a5876 to821b6b6Compare

November 23, 2017 14:11

Copy link

coveralls commentedNov 23, 2017

Coverage increased (+0.003%) to 96.258% when pulling821b6b6 on effigies:rf/mgh into2817d1b on nipy:master.

RF/TEST: Check endianness rather than ignoring

fede660

Copy link

coveralls commentedDec 9, 2017

Coverage increased (+0.004%) to 96.258% when pullingfede660 on effigies:rf/mgh into2817d1b on nipy:master.

RF: Remove endianness, override copy, as_byteswapped, diagnose_binary…

0d4cdd0

…block

Copy link

coveralls commentedDec 9, 2017

Coverage decreased (-0.02%) to 96.236% when pulling0d4cdd0 on effigies:rf/mgh into2817d1b on nipy:master.

effigies added2 commits

December 9, 2017 11:21

FIX: Do not double-set attributes

93f3c70

TEST: TestWrapStruct features, set_zooms, deprecated fields

664f11c

Copy link

coveralls commentedDec 9, 2017

Coverage increased (+0.04%) to 96.293% when pulling664f11c on effigies:rf/mgh into2817d1b on nipy:master.

effigies added2 commits

December 9, 2017 13:00

TEST: Add and test version check

3f62445

TEST: Test MGHHeader._HeaderData fully

fa73ad4

effigies force-pushed therf/mgh branch from30fc9ab tofa73ad4Compare

December 9, 2017 18:28

Copy link

MemberAuthor

effigies commentedDec 9, 2017

@matthew-brett Got back to this. I think I've addressed all of your comments, and I've fleshed out some more tests. Ready for another review.

Copy link

coveralls commentedDec 9, 2017•
edited
Loading

Coverage increased (+0.08%) to 96.331% when pullingfa73ad4 on effigies:rf/mgh into2817d1b on nipy:master.

Copy link

MemberAuthor

effigies commentedDec 18, 2017

Hi@matthew-brett, just bugging you about this. No worries if you're busy.

matthew-brett reviewed

Dec 18, 2017

View reviewed changes

Copy link

Member

matthew-brett left a comment

Choose a reason for hiding this comment

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

Just a few small comments.

nibabel/freesurfer/mghformat.py Outdated

		''' Set zooms into header fields

		See docstring for ``get_zooms`` for examples
		Sets the spaing of voxels in the x, y, and z dimensions.

Copy link

Member

matthew-brettDec 18, 2017

Choose a reason for hiding this comment

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

Typo 'spacing'

nibabel/freesurfer/mghformat.py


		See docstring for ``get_zooms`` for examples
		Sets the spaing of voxels in the x, y, and z dimensions.
		For four-dimensional files, a temporal zoom (repetition time, or TR, in