【Python】import SimpleITK as sitk

1、操作系统接口

(1) os模块提供了不少与操作系统相关联的函数。

import os
print(os.getcwd())   #输出当前的工作目录
os.chdir('path')     #修改当前工作目录到指定路径
print(os.getcwd())  # 从输出可以看出当前路径已经改变

os.system('mkdir today')       # 执行系统命令 mkdir 

每一条system函数执行时，会创建一个子进程在系统上执行命令行，子进程的执行结果无法影响主进程。
mkdir:在当前路径下创建名为today的新目录

建议使用 “import os” 风格而非 “from os import *”。这样可以保证随操作系统不同而有所变化的 os.open() 不会覆盖内置函数 open()。
在使用 os 这样的大型模块时内置的 dir() 和 help() 函数非常有用:

dir(os)
help(os)

(2) 针对日常的文件和目录管理任务，:mod:shutil 模块提供了一个易于使用的高级接口:

import shutil
shutil.copyfile('file1','file2')  #将文件1拷贝给文件2
shutil.move('file1','file2')  #将文件1移动到文件2目录下

2、文件通配符

glob模块提供了一个函数用于从目录通配符搜索中生成文件列表:

import glob
print(glob.glob('path/*.py'))     #输出指定路径下的所有.py文件

3、命令行参数

通用工具脚本经常调用命令行参数。这些命令行参数以链表形式存储于 sys 模块的 argv 变量。

当程序执行时，Python 从命令行获取所有值并将它们存储在 sys.argv 列表中。列表的第一个元素 sys.argv[0] 是脚本的完整路径(或脚本名称——取决于具体操作系统)。列表的第二个元素是脚本的第一个命令行参数，即 sys.argv[1]，依此类推。

（1） sys.argv
（2） getopt、getopt.getopt 方法、Exception getopt.GetoptError
（3） argparse

4、使用SimpleITK读取、保存、处理nii文件

(1) nii格式

nii格式就是后缀名为.nii或.nii.gz的文件，该格式又叫NIfTI-1。MRI图像或者CT图像通常会以这种格式保存。

这种格式的作用，简单来理解就是将索引坐标映射到体素坐标。在计算机中的数据都是离散的，比如一个数组可以按照索引取对应的值，但是在现实生活中，距离都是连续的，那如何给计算机中离散的点（索引坐标）赋予一个连续的坐标（即体素坐标）呢？

在nii格式中，为了将索引坐标(数组下标)映射到体素坐标（空间坐标），除了保存图像的数据外，即那一个个离散的像素，还保存了一些额外信息，比如每个像素间的距离，原点坐标，方向等等，这样根据像素间的距离就可以计算某一像素在空间中真正的坐标了。

在做深度学习时，训练过程中我们并不关心空间坐标，只需要把nii文件中保存的像素值拿出来转化为numpy或tensor输入网络就可以了，空间距离对我们也没什么用。但当我们做对比实验时，会发现将不同网络输出的预测结果保存为nii文件，然后用3DSlicer打开后，二者的预测结果可能完全不一样，这是因为起始坐标和原始label是对不上的，所以需要设置成一样的。

(2) 读取nii成numpy格式

只需要使用SimpleITK先读取图像，然后从图像中获取数组就可以了。最后image就是一个numpy数组，但是需要注意的是nii文件默认保存数据的顺序是[x, y, z]，但是numpy数组保存数据的顺序是[z, y, x]，刚好是反过来的.。因为在训练时一般是[x, y, z]，所以我们在dataset中需要将图像的坐标轴转换一下，即做一个transpose(2, 1, 0)操作。

import SimpleITK as sitk
image_path = ''
image = sitk.ReadImage(image_path)
image = sitk.GetArrayFromImage(image)

(3) 将numpy格式保存成nii

仍需要注意的是，可能从网络中输出预测结果的坐标顺序为[x, y, z]，但是保存时numpy数组轴的顺序一定要是[z, y, x ]，做一下transpose(2, 1, 0)或permute(2, 1, 0)操作，这样才是正确的。当然，这一步操作取决于dataset类是否交换了维度，如果交换了维度，网络输出后自然需要交换回来。

# image是一个三维numpy数组，image_path是要保存的路径
sitk.WriteImage(image, image_path)

(4) 什么是origin、Direction、Spacing，以及如何设置它们

origin就是原点的坐标，direction大概是轴的方向，spacing就是每个像素间所代表的真实世界的距离。对这些有一个大概的认识后，实际使用中，我们只需要以原图为基准，将不同网络输出预测结果的元数据属性设置成和原图一样就可以了，这样用3D Slicer打开时位置就能对应起来，需要注意的是保存前的numpy轴的顺序必须得是[z, y, x]，否则需要先交换再设置元数据属性，这样才是正确的！代码如下：

# 读取原图像
origin_path = ''
origin = sitk.ReadImage(origin_path)

# 读取预测图像
pred_path = ''
pred = sitk.ReadImage(pred_path)

# 将预测结果的元数据属性设置成和原图像一样
pred.SetDirection(origin.GetDirection())
pred.SetOrigin(origin.GetOrigin())
pred.SetSpacing(origin.GetSpacing())

# 保存处理后的
sitk.WriteImage(pred, pred_path)

(5) 示例

做了四个对比试验，nnunet，nnformer，unetr，our，根据nii的ID来依次处理，将他们的元数据都设成和原始image的元数据一样。

import SimpleITK as sitk

def setMetaMessage(target, origin):

    target.SetDirection(origin.GetDirection())
    target.SetOrigin(origin.GetOrigin())
    target.SetSpacing(origin.GetSpacing())
    return target


def process(id):
    image_path = 'C:/Users/result'

    image = sitk.ReadImage(os.path.join(image_path, str(id) + '_image.nii.gz'))
    print('image size:', image.GetSize())
    sitk.WriteImage(image, os.path.join('C:/Users/process result', str(id) + '_image.nii.gz'))


    label = sitk.ReadImage(os.path.join(image_path, str(id) + '_label.nii.gz'))
    print('label size:', label.GetSize())
    label = setMetaMessage(label, image)
    sitk.WriteImage(label, os.path.join('C:/Users/process result', str(id) + '_label.nii.gz'))


    nnunet = sitk.ReadImage(os.path.join(image_path, str(id) + '_nnunet.nii.gz'))
    print('nnunet size:', nnunet.GetSize())
    nnunet = setMetaMessage(nnunet, image)
    sitk.WriteImage(nnunet, os.path.join('C:/Users/process result', str(id) + '_nnunet.nii.gz'))


    nnformer = sitk.ReadImage(os.path.join(image_path, str(id) + '_nnformer.nii.gz'))
    print('nnformer size:', nnformer.GetSize())
    nnformer = setMetaMessage(nnformer, image)
    sitk.WriteImage(nnformer, os.path.join('C:/Users/process result', str(id) + '_nnformer.nii.gz'))


    unetr = sitk.ReadImage(os.path.join(image_path, str(id) + '_unetr.nii.gz'))
    unetr = sitk.GetArrayFromImage(unetr)
    unetr = unetr.transpose(2, 1, 0)  # 换轴
    unetr = sitk.GetImageFromArray(unetr)
    print('unetr size:', unetr.GetSize())
    unetr = setMetaMessage(unetr, image)
    sitk.WriteImage(unetr, os.path.join('C:/Users/process result', str(id) + '_unetr.nii.gz'))


    our = sitk.ReadImage(os.path.join(image_path, str(id) + '_our.nii.gz'))
    our = sitk.GetArrayFromImage(our)
    our = our.transpose(2, 1, 0)  # 换轴
    our = sitk.GetImageFromArray(our)
    print('our size:', our.GetSize())
    our = setMetaMessage(our, image)
    sitk.WriteImage(our, os.path.join('C:/Users/process result', str(id) + '_our.nii.gz'))



if __name__ == '__main__':

    process(447)

(6) 重采样

重采样可以这样理解，现在我们有一个2m2m2m的正方体，它的像素分辨率为100100100，即每个方向都存了100个离散的像素点，现在我们保持正方体的尺寸不变，还是2m2m2m，但是像素分辨率插值变为150150150，这样就缩小了每个像素间的space。所以这个函数也是需要我们传入一个目标space。可以简单理解为，体积不变，密度增大了，类比二维的resize操作~

重采样改变的是每个像素间所代表的物理距离，在重采样代码中，我们是重采样前像素个数为x，每个像素间的距离(即Spacing)为y；重采样后的像素个数为z，每个像素间的距离为m；我们为了保持物理体积不变(即xy = zm)，所以当改变m(即Spacing)后，相应的像素个数z也会通过插值改变。

在神经网络中，我们除了通过重采样改变像素个数外，还可以通过Resize操作来改变像素点。只不过为了最终可视化和原label大小一样，最好先把两者的Spacing调整为一样的。否则会出现长宽比发生变形的情况。重采样代码如下：

def resampleVolume(outspacing, vol, type):
    """
    将体数据重采样的指定的spacing大小\n
    paras：
    outpacing：指定的spacing，例如[1,1,1]
    vol：sitk读取的image信息，这里是体数据
    type：指定插值方法，一般对image采取线性插值，对label采用最近邻插值
    return：重采样后的数据
    """
    outsize = [0, 0, 0]
    # 读取文件的size和spacing信息
    inputsize = vol.GetSize()
    inputspacing = vol.GetSpacing()
 
    transform = sitk.Transform()
    transform.SetIdentity()
    # 计算改变spacing后的size，用物理尺寸/体素的大小
    outsize[0] = round(inputsize[0] * inputspacing[0] / outspacing[0])
    outsize[1] = round(inputsize[1] * inputspacing[1] / outspacing[1])
    outsize[2] = round(inputsize[2] * inputspacing[2] / outspacing[2])
 
    # 设定重采样的一些参数
    resampler = sitk.ResampleImageFilter()
    resampler.SetTransform(transform)
    # 图像使用线性插值，标签使用最近邻插值
    if type == 'linear':
        resampler.SetInterpolator(sitk.sitkLinear)
        resampler.SetOutputPixelType(sitk.sitkFloat32)  # image用float32存
    else:
        resampler.SetInterpolator(sitk.sitkNearestNeighbor)
        resampler.SetOutputPixelType(sitk.sitkUInt8)  # 标签用int8存储
    resampler.SetOutputOrigin(vol.GetOrigin())
    resampler.SetOutputSpacing(outspacing)
    resampler.SetOutputDirection(vol.GetDirection())
    resampler.SetSize(outsize)
    newvol = resampler.Execute(vol)
    return newvol

(7) simpleITK读取nii文件并显示

import SimpleITK as sitk
from matplotlib import pyplot as plt
 
def showNii(img):
    for i in range(img.shape[0]):
        plt.imshow(img[i,:,:],cmap='gray')
        plt.show()
 
itk_img = sitk.ReadImage('C:\\Users\\86472\\Desktop\\1552282517.831928.nii')
img = sitk.GetArrayFromImage(itk_img)
showNii(img)

(8) SimpleITK读写nii.gz文件

## using simpleITK to load and save data.
import SimpleITK as sitk
itk_img = sitk.ReadImage('./nifti.nii.gz')
img = sitk.GetArrayFromImage(itk_img)
print("img shape:",img.shape)
 
## save 
out = sitk.GetImageFromArray(img)
# # out.SetSpacing(itk_img.GetSpacing())
# # out.SetOrigin(itk_img.GetOrigin())
sitk.WriteImage(out,'simpleitk_save.nii.gz')

参考：
https://blog.csdn.net/qq_43705697/article/details/124928070