ITK是一个功能很强大的医学图像处理公开库,搭配VTK用以显示图像,可以实现几乎所有医学图像处理的功能需要。ITK通常以C++包进行提供,当然也可以自己编译为Python包,不过编译过程比较繁琐耗时,而且很容易踩坑。但ITK官方进行的Python封装SimpleITK,则直接可以拿来使用;虽然有部分ITK的功能没有包含,但已基本够用了。我们在处理医学图像时,使用的基本都是SimpleITK。
本文就简单总结一下我们在处理这些图像时的经验,以便备忘,并为后来者参考。
1.读取文件
读取DICOM序列
医学图像中一个CT序列包含很多张图片,即一个case包含许多slice,使用SimpleITK可以直接读取一个序列,并方便地得到各种参数,将图像数据转换成numpy Array:
- import SimpleITK as sitk
- import numpy as np
- reader = sitk.ImageSeriesReader()
- dicom_names = reader.GetGDCMSeriesFileNames(case_path)
- reader.SetFileNames(dicom_names)
- image = reader.Execute()
- image_array = sitk.GetArrayFromImage(image) # z, y, x
- origin = image.GetOrigin() # x, y, z
- spacing = image.GetSpacing() # x, y, z
需要注意的是,SimpleITK读取的图像数据的坐标顺序为zyx,即从多少张切片到单张切片的宽和高;而据SimpleITK Image获取的origin和spacing的坐标顺序则是xyz。这些需要特别注意。
读取DICOM单张图片
可以将一个DICOM序列作为一个整体一次读入,也可以一张一张地读入每张切片:
- import SimpleITK as sitk
- import numpy as np
- image = sitk.ReadImage(slice_path)
- image_array = sitk.GetArrayFromImage(image) # z, y, x
这里需要注意的除了坐标顺序是zyx之外,还需注意,即使读取单张切片,所得到的结果也是3维的,只不过第一维是1。
读取mhd文件
涉及DICOM序列时,为了传输方便(从上百个dcm文件到一个mhd文件),很多情况下以mhd文件格式进行呈现,不过mhd文件只是一个很小的包含数据信息的文件,同时搭配的通常还有一个二进制的数据文件(格式为raw或zraw等)。使用SimpleITK读取这种文件也比较方便。
- import SimpleITK as sitk
- import numpy as np
- image = sitk.ReadImage(mhd_path)
- image_array = sitk.GetArrayFromImage(image) # z, y, x
- origin = image.GetOrigin() # x, y, z
- ...
有时候不想整个读取数据(因为比较大,读取处理比较慢),想要读取的只是一些基本信息,比如origin、spacing等等。这时可以只读取mhd文件,据此获取信息,读取方法比较简单,不再赘述。
2.处理文件
处理DICOM文件主要有插值等操作,可以直接使用SimpleITK(或者说是ITK)的相关函数,并通过pipeline结构进行处理。
插值
-
- import SimpleITK as sitk
- reader = sitk.ImageSeriesReader()
- dicom_names = reader.GetGDCMSeriesFileNames(case_path) reader.SetFileNames(dicom_names)
- image = reader.Execute()
- resample = sitk.ResampleImageFilter()
- resample.SetOutputDirection(image.GetDirection())
- resample.SetOutputOrigin(image.GetOrigin())
- newspacing = [1, 1, 1]
- resample.SetOutputSpacing(newspacing)
- newimage = resample.Execute(image)