to_nibabel() produces wrong transforms about antspy

[cookpa]

ITK uses LPS coordinates and NIFTI RAS. This is handled in the I/O, when ITK converts qform into a direction cosine matrix or writes a cosine matrix back out as a qform.

https://github.com/InsightSoftwareConsortium/ITK/blob/master/Modules/IO/NIFTI/src/itkNiftiImageIO.cxx

See methods SetImageIOOrientationFromNIfTI and SetNIfTIOrientationFromImageIO.

According to this page

https://nipy.org/nibabel/coordinate_systems.html

NiBabel always wants RAS coordinates so I think it should be safe to apply the same conversion steps in the to_nibabel and from_nibabel functions.

Some more on this:

https://itk.org/Doxygen/html/classitk_1_1ImageBase.html#a23ebc76de3b60bb1eeabf66cd4dabc48

This is the base ITK class that handles voxel to physical space conversion. It references this page

https://itk.org/Wiki/Proposals:Orientation#Some_notes_on_the_DICOM_convention_and_current_ITK_usage

[lmoesch]

This issue is not nibabel related, since the qform is saved within the nifti header and should be parsed the same way in all implementations. I checked the values within SPM and MRIcroN and they are the same as the nibabel output, so there has to be an issue by extracting the values via the ReadImageInformation().

As reading and writing through ITK methods results in correct qforms, my guess is that ITK stores the quatern_b, quatern_c and quatern_d values from the header as well to correct the directions.

[stnava]

thx - i am leaving this open until a test is added to help prevent this from coming up in the future.

[stnava]

consider 1e73b82

[stnava]

very nice find. a pull request would be helpful. it's probably just changing some indexing somewhere in the function get_qform.

[lmoesch]

As for the pull request, I'm not really sure where to fix this issue (prolly somewhere in the itkNiftiImageIO.cxx). As far as I understand, the direction and origin are parsed from the q_form without paying respect to quatern_b, quatern_c and quatern_d.

For my current code I wrote a wrapper to replace ants.from_nibabel()

def nii_to_ants(nii_img : nib.Nifti1Image):
    q_form = nii_img.get_qform()
    spacing = nii_img.header['pixdim'][1:4]
    
    origin = q_form[:3,3]
    direction = q_form[:3,:3] / spacing

    return ants.from_numpy(data=nii_img.get_data().astype(np.float), origin=origin.tolist(), spacing=spacing.tolist(), direction=direction)

This however assumes 3D data and Nifti1 headers.

[stnava]

indeed, i am not sure where you would fix this either. when i search:

https://github.com/ANTsX/ANTsPy/search?q=get_qform&unscoped_q=get_qform

i dont find any code in ANTsPy. thoughts?

[lmoesch]

When I understand the code correctly, the direction, origin and spacing is derived from the ITK function itk::ImageIOBase::ReadImageInformation() in ANTsPy/ants/lib/LOCAL_antsImageHeaderInfo.cxx So this has to be an issue in the used ITK version.

[stnava]

a reproducible / self-contained example of this problem:

import ants
import nibabel as nib
fn = ants.get_data( "ch2" )
ants_mni = ants.image_read( fn )
nii_mni = nib.load( fn )
ants_mni = ants_mni.to_nibabel()
print(ants_mni.get_qform())
print(nii_mni.get_qform())

[stnava]

the problem may be resolvable by alterations to this bit of code:

ANTsPy/ants/utils/convert_nibabel.py

Lines 18 to 21 in ad7fcd2

	array_data = image.numpy()
	affine = np.hstack([image.direction*np.diag(image.spacing),np.array(image.origin).reshape(3,1)])
	affine = np.vstack([affine, np.array([0,0,0,1.])])
	new_img = nib.Nifti1Image(array_data, affine)

[lmoesch]

The problem is, that the metadata already mismatches upon loading the image via ants.image_read().

    nii_img = nib.load('MNI152.nii.gz')
    ants_img = ants.from_nibabel(nii_img) # This is the same as ants_img = ants.image_read('MNI152.nii.gz') 

    print(nii_img.header.get_qform()[:3,3].tolist())
    print(ants_img.origin) 

(90.0, -126.0, -72.0)
(-90.0, 126.0, -72.0)

[stnava]

this particular difference may be the difference between LPS vs whatever nibabel uses.

[stnava]

can you share the image you are using ( upload to issue ) so we can be on the same page.

[lmoesch]

https://github.com/nilearn/nilearn/raw/master/nilearn/datasets/data/avg152T1_brain.nii.gz

[stnava]

upon further reflection, i think that this is just an image orientation difference - ie a difference in coordinate systems (itk uses LPS). @cookpa - can you briefly comment and / or suggest a relevant link ?

[stnava]

more discussion here: https://discourse.itk.org/t/nifti-orientation-issues/431/2 it doesnt "convert images" ... rather, the reference frame is interpreted as LPS. the issue here is, i think, that the ANTsPy "to_nibabel" function is probably wrong.

[lmoesch]

To make a long story short: ITK assumes NIFTI images to be in RAS space and simply flips the first two axes upon loading (code). Therefore it is sufficient to do the same upon to_nibabel().

However, I tried to utilize ants.reorient_image2() in order to achieve this and watched again a strange behavior. When loading an NIFTI image in RAS orientation

[[1. 0. 0.]
 [0. 1. 0.]
 [0. 0. 1.]]

this results in a direction

[[-1. 0. 0.]
 [0. -1. 0.]
 [0. 0. 1.]]

as the first axes are flipped. But when applying ants.reorient_image2(orientation='RAS') to it, it results in

[[1. 0. 0.]
 [0. 1. 0.]
 [0. 0. -1.]]

where the Z-Axes is flipped, which looks like a bug within itkOrientImageFilter.hxx to me.

[stnava]

i do not think that is likely. behavior looks correct here:

>>> mni = ants.image_read(ants.get_data('mni'))

>>> mni
ANTsImage (RPI)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (182, 218, 182)
	 Spacing    : (1.0, 1.0, 1.0)
	 Origin     : (-90.0, 126.0, -72.0)
	 Direction  : [ 1.  0.  0.  0. -1.  0.  0.  0.  1.]

>>> mniRAS
ANTsImage (RAS)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (182, 218, 182)
	 Spacing    : (1.0, 1.0, 1.0)
	 Origin     : (-90.0, -91.0, 109.0)
	 Direction  : [ 1.  0.  0.  0.  1.  0.  0.  0. -1.]

>>> mniLPS
ANTsImage (LPS)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (182, 218, 182)
	 Spacing    : (1.0, 1.0, 1.0)
	 Origin     : (91.0, 126.0, 109.0)
	 Direction  : [-1.  0.  0.  0. -1.  0.  0.  0. -1.]

[stnava]

however, your comment:

"To make a long story short: ITK assumes NIFTI images to be in RAS space and simply flips the first two axes upon loading (code). Therefore it is sufficient to do the same upon to_nibabel()."

is probably right. i guess the question is what happens in 2D and 4D.

[lmoesch]

i do not think that is likely. behavior looks correct here:

>>> mni = ants.image_read(ants.get_data('mni'))

>>> mni
ANTsImage (RPI)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (182, 218, 182)
	 Spacing    : (1.0, 1.0, 1.0)
	 Origin     : (-90.0, 126.0, -72.0)
	 Direction  : [ 1.  0.  0.  0. -1.  0.  0.  0.  1.]

>>> mniRAS
ANTsImage (RAS)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (182, 218, 182)
	 Spacing    : (1.0, 1.0, 1.0)
	 Origin     : (-90.0, -91.0, 109.0)
	 Direction  : [ 1.  0.  0.  0.  1.  0.  0.  0. -1.]

>>> mniLPS
ANTsImage (LPS)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (182, 218, 182)
	 Spacing    : (1.0, 1.0, 1.0)
	 Origin     : (91.0, 126.0, 109.0)
	 Direction  : [-1.  0.  0.  0. -1.  0.  0.  0. -1.]

According to nibabel documentation the z-axis value should be positive in RAS, or am I missing a point here?

[lmoesch]

[...] the question is what happens in 2D and 4D.

let me answer with the code from ITK:

const int max_defined_orientation_dims = ( dims > 3 ) ? 3 : dims;

It flips the first 2 dimensions in every case. And every dimension above the third is simply discarded.

[stnava]

ok - here is a fix:

>>> import ants
>>> import nibabel as nib
>>> fn = "avg152T1_brain.nii.gz"
>>> ants_img = ants.image_read( fn )
>>> nii_mni = nib.load( fn )
>>> ants_mni = ants_img.to_nibabel()
>>> print(ants_mni.get_qform())
[[  -2.    0.    0.   90.]
 [   0.    2.    0. -126.]
 [   0.    0.    2.  -72.]
 [   0.    0.    0.    1.]]
>>> print(nii_mni.get_qform())
[[  -2.    0.    0.   90.]
 [   0.    2.    0. -126.]
 [   0.    0.    2.  -72.]
 [   0.    0.    0.    1.]]
>>> ants.from_nibabel( nii_mni )
ANTsImage (RPI)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (91, 109, 91)
	 Spacing    : (2.0, 2.0, 2.0)
	 Origin     : (-90.0, 126.0, -72.0)
	 Direction  : [ 1.  0.  0.  0. -1.  0.  0.  0.  1.]

>>> ants_img 
ANTsImage (RPI)
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (91, 109, 91)
	 Spacing    : (2.0, 2.0, 2.0)
	 Origin     : (-90.0, 126.0, -72.0)
	 Direction  : [ 1.  0.  0.  0. -1.  0.  0.  0.  1.]

[stnava]

b15e286

[stnava]

apologies for pushing this to master but it seems to be correct .... agree?

[lmoesch]

I dont think the last-ind has to be flipped. But I'll check into this.

[lmoesch]

I'll checked into the fix and as far as I see it is not quite correct. ITK changes the signs of the first two elements of all vectors. Fortunately, this simplifies the whole fix to:

affine[:2,:] *= -1

UPDATE: The way the affine is created, you fix works nevertheless. When considering using the real transformation in the future, tho above approach is the way to go.