How to write a color 3D NIfTI with NiBabel?

Question

I've had no problem writing out 3D grayscale .nii files with nibabel and opening them in NIfTI viewers (Mango, MCIcron). However I haven't been able to write out 3D color, as each RGB plane is being interpreted as a different volume. E.g. the output from this:

import nibabel as nib
import numpy as np
nifti_path = "/my/local/path"
test_stack = (255.0 * np.random.rand(20, 201, 202, 3)).astype(np.uint8)
ni_img = nib.Nifti1Image(test_stack, np.eye(4))
nib.save(ni_img, nifti_path)

is seen as 3 separate 20x201x202 volumes. I also tried putting the color planes in the first axis (i.e. np.random.rand(3, 20, 201, 202)), but get the same issue. Looking around a bit it seems there is a "dataset" field that needs to be set to 128 for 24-bit RGB planar images. One nice thing about nibabel is how it automatically sets up the header based on the numpy array fed to it. However this is an ambiguous case and if I print the header information I can see it setting the datatype to 2 (uint8), which presumably is why viewers are interpreting it as separate volumes, not RGB24. I don't see any official support in the API for setting the datatype, but the documentation does mention access to the raw fields for those with "great courage". Doing this, i.e.

hdr = ni_img.header
raw = hdr.structarr
raw['datatype'] = 128

works in changing the header value

print(hdr)

gives "datatype : RGB" but when writing

nib.save(ni_img, nifti_path)

I get an error:

File "<python path>\lib\site-packages\nibabel\arraywriters.py", line 126, in scaling_needed
raise WriterError('Cannot cast to or from non-numeric types')
nibabel.arraywriters.WriterError: Cannot cast to or from non-numeric types

The exception is raised if some arr_dtype != out_dtype, so presumably my hacking of the raw header is causing some inconsistency down the line.

So, is there a proper way to do this?

Answer 1

Thanks to matthew.brett on the Neuroimaging analysis mailing list, I'm able to write out 3-d color NIfTI like so:

# ras_pos is a 4-d numpy array, with the last dim holding RGB
shape_3d = ras_pos.shape[0:3]
rgb_dtype = np.dtype([('R', 'u1'), ('G', 'u1'), ('B', 'u1')])
ras_pos = ras_pos.copy().view(dtype=rgb_dtype).reshape(shape_3d)  # copy used to force fresh internal structure
ni_img = nib.Nifti1Image(ras_pos, np.eye(4))
nib.save(ni_img, output_path)

Answer 2

Using the suggested method works, viewing the RGB volume, with e.g. ITK-SNAP is no problem,

# ras_pos is a 4-d numpy array, with the last dim holding RGB
shape_3d = ras_pos.shape[0:3]
rgb_dtype = np.dtype([('R', 'u1'), ('G', 'u1'), ('B', 'u1')])
ras_pos = ras_pos.copy().view(dtype=rgb_dtype).reshape(shape_3d)  # copy used 
#to force fresh internal structure
ni_img = nib.Nifti1Image(ras_pos, np.eye(4))
nib.save(ni_img, output_path)

but, reloading the image only works with get_data()

ni_img = nib.load(output_path)

# this will result in error
data = img.get_fdata()

# this will work fine, but get_data() is going to be removed.
data = img.get_data()

This is critical, because get_data() will be removed in a future release. Recommended is only using get_fdata(). There's currently an error in the get_fdata() method, which can not cast the RGB data in a valid numpy array.