SimpleITK 获取CT spacing 底层原理
一、层厚、层间距概念
层厚: CT扫描机扫描出来的断层的层的厚度, 通常用Slice thickness表示, 比如 5mm
层间隔:一般用 Slice interval 或 Slice increment来表示,比如 5mm。在 CT 中,层间隔表示相邻两层中心之间的距离
如上图所示, 在CT中, 相邻两层CT中心之间的距离代表层间距, 当层厚==层间距时,无重叠, 代表CT断层之间是连续的, 那么同理, 当层厚!=层间距时, 代表断层之间是有距离或者重叠的
二、pydicom中的层厚以及spacing
spacing:
一张二维的医学图像,图像尺寸为[6x7],每个像素之间的物理空间距离(spacing),x轴、y轴分别是(20,30)。扩展到三维医学图像,就出现【体素】,就是三维的“像素”。相比像素有三个轴(x、y、z),同样每个轴都有spacing。 意思相同大小尺寸的医学图片,spacing也是不同。
如下是pydicom以及SimpleITK获取spacing:
import pydicom as dicom
import numpy as np
from pathlib import Path
PathDicom = 'dcm序列所在文件夹路径'
data_dir = Path(PathDicom)
dicom_files = list(data_dir.glob('*.dcm')) # 获取目录下所有 .dcm 文件
dicom_files.sort(key=lambda x: int(x.stem.rsplit('.', 1)[1])) #将 .dcm文件 按序号排列 注意: sorted() 方法会按字符排列, 此处需要按末尾数字排列
slices = [dicom.read_file(dicom_file, force=True) for dicom_file in dicom_files]
# 采用以下方法返回的SliceThickness是错误的
print(slices[0].SliceThickness)
# 正确的层厚获取方法
thickness = np.abs(slices[0].SliceLocation - slices[1].SliceLocation) # 获取 dcm 文件的层厚
# 获取整个 CT 序列的 spacing, 等同于 SimpleITK 的 GetSpacing() 方法
spacing_z_list = []
for i in range(1, len(slices)):
spacing_z_list.append(slices[i].SliceLocation - slices[i - 1].SliceLocation)
spacing_z = abs(sum(spacing_z_list) / len(spacing_z_list)) # 获取Z轴的spacing
print(f"spacing_z : {spacing_z}")
如上代码所示, 因为dicom文件层次之间可能会存在重叠或者间隔, 故z轴的spacing应该取所有相邻序列层厚的均值, 即abs(sum(spacing_z_list) / len(spacing_z_list))
。
注意:在临床CT数据中, 有时候会遇到一些数据在ITK-SNAP软件中显示如下, 有多个component:
这种情况下, 使用SimpleITK获取的Z轴spacing就是利用我上述代码的逻辑计算的, 并且单独的component之间数据的层厚是相同的, 只是因为两个component连接处的层厚大大超出单独component内的层厚, 故不注意的话会造成一些问题。这个component信息可以通过pydicom包中AcquisitionNumber获取
singleDicomFile # 单独的一个dcm文件 print(singleDicomFile.AcquisitionNumber) # 查看某个dcm文件所属component
参考链接:
【CTMR技术】CT与MR中的层厚&层间隔一样吗?
DICOM数据字段说明
关于DICOM文件的SliceThickness参数错误问题
【碎碎念】医学图像处理工具SimpleITK与DataSet