几种读nii图像方法的轴序比较

读 .nii / .nii.gz 图像并转成 numpy 可用 medpy.io、nibabel、itk、SimpleITK 几种方法，然而几种方法读出来的轴序有出入，本篇比较此几种方法。

Datum

所用数据来自 verse，经 iTomxy/data/verse/preprocess.py 预处理，朝向和轴序是 RAI。三视图（可视化：show_nii）：

其中最短的是 LR = 35，次短 AP = 118，最长 SI = 214，所以如果维持 RAI 轴序，读出来的形状应是 [35, 118, 214]，给下面代码结果作参考。

Test

• 模拟 preprocess_mmwhs.py，itk 兼测与 numpy、simpleitk 交叉使用（保存、重读）

import os
import numpy as np
import itk
import SimpleITK as sitk
import nibabel as nib
import medpy.io as medio

# 所选的一个 verse 数据
f = os.path.expanduser("~/data/verse/processed-verse19/test/sub-verse012_image.nii.gz")

# medpy
lab_medio, _ = medio.load(f)
print("medpy:", lab_medio.shape)

# nibabel
lab_nib = nib.load(f) # <class 'nibabel.nifti1.Nifti1Image'>
print("nibabel:", lab_nib.shape)
lab_nib_np = lab_nib.get_fdata()
print("nibabel -> numpy:", lab_nib_np.shape)

# itk（兼测与 numpy、sitk 交叉使用）
lab_itk = itk.imread(f)
print("itk:", lab_itk.shape)
lab_itk_np = itk.GetArrayViewFromImage(lab_itk)
print("itk -> numpy:", lab_itk_np.shape)
# numpy、sitk 保存
np.save("lab_itk_np.npy", lab_itk_np)
lab_itk_np_sitk = sitk.GetImageFromArray(lab_itk_np)
print("lab_itk_np_sitk:", lab_itk_np_sitk.GetSize())
sitk.WriteImage(lab_itk_np_sitk, "lab_itk_np_sitk.nii.gz")
# numpy 重读
lab_itk_np_re = np.load("lab_itk_np.npy")
print("lab_itk_np_re:", lab_itk_np_re.shape)
# itk 重读
lab_itk_np_sitk_itk = itk.imread("lab_itk_np_sitk.nii.gz")
print("lab_itk_np_sitk_itk:", lab_itk_np_sitk_itk.shape)
lab_itk_np_sitk_itk_np = itk.GetArrayViewFromImage(lab_itk_np_sitk_itk)
print("lab_itk_np_sitk_itk -> numpy:", lab_itk_np_sitk_itk_np.shape)
# sitk 重读
lab_itk_np_sitk_sitk = sitk.ReadImage("lab_itk_np_sitk.nii.gz")
print("lab_itk_np_sitk_sitk:", lab_itk_np_sitk_sitk.GetSize())
lab_itk_np_sitk_sitk_np = sitk.GetArrayFromImage(lab_itk_np_sitk_sitk)
print("lab_itk_np_sitk_itk -> numpy:", lab_itk_np_sitk_sitk_np.shape)

# sitk
lab_sitk = sitk.ReadImage(f)
print("sitk:", lab_sitk.GetSize())
lab_sitk_np = sitk.GetArrayFromImage(lab_sitk)
print("sitk -> numpy:", lab_sitk_np.shape)

输出：

medpy: (35, 118, 214)
nibabel: (35, 118, 214)
nibabel -> numpy: (35, 118, 214)
itk: (214, 118, 35)
itk -> numpy: (214, 118, 35)
lab_itk_np_sitk: (35, 118, 214)					# numpy 转 sitk（保存前）
lab_itk_np_re: (214, 118, 35)					# numpy 重读
lab_itk_np_sitk_itk: (214, 118, 35)				# itk 重读
lab_itk_np_sitk_itk -> numpy: (214, 118, 35)	# itk 重读后转 numpy
lab_itk_np_sitk_sitk: (35, 118, 214)			# sitk 重读
lab_itk_np_sitk_itk -> numpy: (214, 118, 35)	# sitk 重读后转 numpy
sitk: (35, 118, 214)
sitk -> numpy: (214, 118, 35)

结论：

• medpy.io、nibabel 读不改轴序；
• itk 读会逆转轴序（本例即 RAI 变 IAR）；
• sitk 读不会，但在转 numpy 时逆了，参考 [1,2]。

Copying Meta Information

用 SimpleITK 读 nii、转成 numpy，做一些预处理之后，再转回 SimpleITK 格式，meta information 会变，因为 numpy 本身不会存。为避免 meta information 丢失，可用 CopyInformation 函数从原 SimpleITK 格式图像中复制 meta information（origin、spacing、direction）：

Help on method CopyInformation in module SimpleITK.SimpleITK:

CopyInformation(srcImage) method of SimpleITK.SimpleITK.Image instance
    CopyInformation(Image self, Image srcImage)
    
    
    Copy common meta-data from an image to this one.
    
    
    Copies the Origin, Spacing, and Direction from the source image to
    this image. The meta-data dictionary is not copied.
    
    It is required for the source Image's dimension and size to match, this image's attributes, otherwise an
    exception will be generated.

例程：

import numpy as np
import SimpleITK as sitk

img = sitk.ReadImage("sub-verse012_ct.nii.gz")
print(img.GetSize(), img.GetDirection()) # 原尺寸、朝向

# 转成 numpy
img_np = sitk.GetArrayFromImage(img)
print(img_np.shape)
print(img_np.min(), img_np.max())

# 截 CT 窗口
window_level = 300 # in [250, 300]
window_width = 250 # in [1000, 1500]
img_clip = np.clip(img_np, window_level - window_width, window_level + window_width)

# 转回 SimpleITK 格式
img_clip_sitk = sitk.GetImageFromArray(img_clip)
print(img_clip_sitk.GetSize(), img_clip_sitk.GetDirection()) # 恢复 meta info 前
# 恢复 meta information
img_clip_sitk.CopyInformation(img)
print(img_clip_sitk.GetSize(), img_clip_sitk.GetDirection()) # 恢复 meta info 后
# 保存
sitk.WriteImage(img_clip_sitk, f"clip-{window_level}-{window_width}.nii")

输出：

# 原
(194, 252, 35) (-0.08921860326264591, 0.12383144150247688, -0.988284277732299, 0.9959593006554869, 0.0008781031215354981, -0.08980145398859996, -0.01025242750260076, -0.9923028786768314, -0.12340942123145082)
# 恢复  meta info 前
(194, 252, 35) (1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0)
# 恢复  meta info 后
(194, 252, 35) (-0.08921860326264591, 0.12383144150247688, -0.988284277732299, 0.9959593006554869, 0.0008781031215354981, -0.08980145398859996, -0.01025242750260076, -0.9923028786768314, -0.12340942123145082)

References

1. Order of images returned by GetArrayFromImage - where is this this documented? #1088
2. Conversion between numpy and SimpleITK