I’m using SimpleITK 1.2.4 to read and write NIfTI images. NIfTI docs say that it “should not occur” that qfac is 0. qfac is the first value of the pixdim array in the header and is expected to be -1 or 1.
When I create an image using GetImageFromArray, the written header is fine. If I read and write and image with SimpleITK, the header is incorrect. Here’s some code to show what I mean, and the output:
import gzip
import struct
from pathlib import Path
from tempfile import NamedTemporaryFile
import numpy as np
import nibabel as nib
import SimpleITK as sitk
def get_qfac(image_path):
image_path = Path(image_path)
with NamedTemporaryFile(suffix='.nii') as f:
if image_path.suffix == '.gz':
with gzip.open(image_path, 'rb') as f_in:
with open(f.name, 'wb') as f_out:
shutil.copyfileobj(f_in, f_out)
image_path = f.name
with open(image_path, 'rb') as f:
fmt = 8 * 'f'
size = struct.calcsize(fmt)
f.seek(76)
chunk = f.read(size)
pixdim = struct.unpack(fmt, chunk)
qfac = pixdim[0]
return qfac
def check_qfac(image_path):
qfac = get_qfac(image_path)
if qfac in (-1, 1):
print('qfac is ok:', qfac)
else:
print(f'qfac is {qfac} in {image_path}')
def main():
filepath_nib = '/tmp/test_nib.nii'
filepath_itk = '/tmp/test_itk.nii'
array = np.random.rand(10,10,10)
affine = np.eye(4)
print('Written with NiBabel')
nib.Nifti1Image(array, affine).to_filename(filepath_nib)
check_qfac(filepath_nib)
print('\nGetImageFromArray, written with ITK')
im = sitk.GetImageFromArray(array)
sitk.WriteImage(im, filepath_itk)
check_qfac(filepath_itk)
print('\nRead and write with ITK')
im = sitk.ReadImage(filepath_nib)
sitk.WriteImage(im, filepath_itk)
check_qfac(filepath_itk)
if __name__ == "__main__":
main()
Output:
Written with NiBabel
qfac is ok: 1.0
GetImageFromArray, written with ITK
qfac is ok: 1.0
Read and write with ITK
qfac is 0.0 in /tmp/test_itk.nii
When you first read an image in ITK, the Image may contain additional meta data tags from the reader. I did a test with NifTi and there were quite a lot of NifTi tags:
Try deleting all MetaData elements. The default behavior of an ITK filter is not to propagate the tags so “I’m += 0” would be the easiest was to accomplish this.
Instead of nifti, can you use a general-purpose 3D image file format, such as nrrd or metaimage?
Nifti makes sense for neuroimaging, but for all other fields, it is so inconvenient - there are many complications due to that directions can be defined using several different approaches, header is not human-readable and not easy to extend, there are lots of irrelevant (neuroimaging specific) fields but nothing for any other specialty, etc.
If there is a bug it would be in the itk:NifTiImageIO class, so that is an ITK bug.
I was able to reproduce the same behavior. However, if all NifTi metadata tags are removed before rewriting the file then the error does not occur. As @lassoan says, the NifTi spec is a bit complicated and it’s not clear what the valid states are for the redundant data.
There is something funny going on here that hopefully can be addressed. However, this is not as much of a “hack” as one might think. It is best practice to inspect the meta-data dictionary before writing an itk::Image. Since the meta-data dictionary is not propagate between filter generally this step is omitted. For this case the reader provides every detail of the NifTi header in the meta-data dictionary. These meta-data fields should not be written out as meta-data, but the fields are represented by then itk::Image’s geometric information which gets written in the NifTi header.
If you would like code that directly does the operations consider:
def strip_metadata(img):
for k in img.GetMetaDataKeys():
img.EraseMetaData(k)
I’ve looked into way to customize how spacing is defined q-form versus s-form. Anyway, I can instruct the NiftiIO to prefer the one or the other upon saving? For instance on load the trick of overwriting s-form/q-form meta data works. That is I can force ITK to prefer s-form instead of the default q-form but now way to achieve the same upon saving
ITK’s NIFTIO always writes out “Method 3” (sform > 0) NIFTI files, from what I can tell. This I believe was chosen because sform can store coordinates with more precision.
qfac is only a valid code for “Method 2”, it is not used for “Method 3” reading, which means that while invalid, it shouldn’t be used by other software reading the file.
The spec is “should not” vs. must not. But we know how well in general the NIFTI spec has been followed…
