📚 课程目录第3章图像预处理与三维重建前处理

3. 图像预处理与三维重建前处理

三维重建并不是把原始 DICOM 序列直接丢给算法就能得到可靠模型。CT/MR 设备输出的体数据往往带有噪声、层厚不均、层间距远大于层内像素间距、不同序列之间空间位置不一致等问题。如果跳过前处理，后续的面绘制、体绘制、阈值分割和几何测量都会继承这些误差。本节把前处理拆成五个相互衔接的环节：去噪、配准、插值、重采样与各向同性体素化，并说明每个环节的输入、关键参数、可观察结果与判断依据。

一、去噪：先分清噪声类型再选滤波器

CT 图像的噪声近似服从泊松分布，在低剂量扫描中表现为明显的颗粒状起伏；MR 图像则同时存在热噪声和 Rician 噪声，背景区域常出现偏亮的偏置。去噪的目标不是把图像“磨平”，而是在抑制噪声的同时保住边缘、细小血管和骨小梁结构。

常用方法有三类。第一类是线性平滑，如高斯滤波，实现简单但会模糊边缘，一般只用于对分辨率要求不高的预览。第二类是非线性边缘保持滤波，如中值滤波、双边滤波、各向异性扩散（Perona-Malik）。双边滤波同时考虑空间距离和灰度差，能在平滑同质区域的同时保留高对比边界，是 CT 骨组织前处理的常用选择。第三类是基于变换域或学习的方法，如小波阈值去噪、BM3D、DnCNN，适合噪声水平较高且对细节要求严格的场景。

以一幅 512×512、层厚 1 mm 的头部 CT 为例，若目视发现脑实质内颗粒噪声明显，可先尝试双边滤波：空间域 σ_s 取 1.5～2.5 像素，灰度域 σ_r 取 20～40 HU。σ_r 过大会把灰质与白质的边界一起抹掉，过小则去噪不足。判断依据是：同质区域的标准差下降，而脑沟、脑回和颅骨内板的边缘梯度基本不变。若使用各向异性扩散，迭代次数通常控制在 5～15 次，扩散系数 κ 取 20～50 HU，迭代过多会出现“块状”伪影。

去噪应在配准和插值之前完成，但不宜过度。过度去噪会改变后续阈值分割的灰度分布，使骨与软组织的分界发生漂移。

二、配准：把不同序列、不同时相拉到同一空间

配准解决的是“同一解剖结构在不同图像中位置不一致”的问题。三维重建中常见的配准场景包括：CT 与 MR 融合显示、PET 与 CT 衰减校正、同一患者不同时间点的随访对比、以及多模态图像引导手术。

配准由四个要素构成：变换模型、相似性测度、优化策略和插值方式。变换模型分刚体（6 自由度：3 平移 + 3 旋转）、仿射（12 自由度）、非线性（B 样条或 Demons）。相似性测度中，同模态常用均方误差或归一化互相关，跨模态常用互信息（MI）或归一化互信息（NMI）。优化常用梯度下降、Powell 或 L-BFGS。

实际操作中，建议采用“粗到精”的两级策略：先对降采样 4 倍的图像做刚体配准，得到初始变换；再在原分辨率上做仿射或 B 样条精配准。以 CT-MR 头部配准为例，刚体阶段互信息值通常从 0.2 左右提升到 0.6 以上，精配准后可达到 0.75～0.85。判断配准是否成功，不能只看测度值，还要用棋盘格叠加或红绿通道融合图目视检查脑室、颅骨和眼球等标志结构是否对齐。

三、插值：重采样时如何“造”出新体素

插值发生在重采样、配准变换和任意角度切面提取中。它的本质是：已知离散网格上的灰度值，估计非网格位置的值。常用方法有最近邻、线性、三次卷积（B 样条）和窗 sinc。

最近邻插值速度快、不产生新灰度值，适合标签图、分割掩膜和骨组织阈值图，但会产生阶梯状边缘。线性插值在 CT 中较常用，平滑度与速度折中。三次卷积在 MR 和需要高质量斜切面的场景中更合适，但会引入轻微振铃。窗 sinc 质量最高，计算量也最大，多用于科研级重采样。

选择原则可以概括为：标签和掩膜用最近邻，灰度图像用线性或三次卷积，需要保持定量 HU 值时优先线性并检查插值后骨皮质的峰值是否被削低。

四、重采样与各向同性体素化

重采样是改变体素网格的空间分辨率或几何形状。各向同性体素化是重采样的一个特例：把原本层间距大于层内间距的体数据，重采样为三个方向间距相等的体素，例如把 0.7 mm×0.7 mm×5 mm 的 CT 变成 0.7 mm×0.7 mm×0.7 mm 或 1 mm×1 mm×1 mm。

为什么必须做？因为面绘制和体绘制中的梯度计算、等值面提取、距离变换都假设体素在三个方向上的物理尺寸一致。若 z 方向间距是 x 方向的 7 倍，梯度算子会严重各向异性，重建出的骨表面在上下方向被“压扁”或出现阶梯。

实现上，先读取 DICOM 中的 PixelSpacing（行、列间距）和 SliceThickness/SpacingBetweenSlices，计算目标网格。例如原始尺寸 512×512×120，像素间距 0.7 mm×0.7 mm，层间距 5 mm，目标各向同性 1 mm。则新网格为：x 方向 512×0.7/1≈358，y 方向 358，z 方向 120×5/1=600。然后用线性插值在 z 方向补出 600 层。若目标设为 0.7 mm，则 z 方向需补到约 857 层，数据量增加约 7 倍，需权衡内存与精度。

# 以 SimpleITK 为例的各向同性重采样
import SimpleITK as sitk

img = sitk.ReadImage("ct_series.mha")
orig_spacing = img.GetSpacing()          # 例如 (0.7, 0.7, 5.0)
orig_size    = img.GetSize()             # 例如 (512, 512, 120)

target_spacing = (1.0, 1.0, 1.0)
target_size = [
    int(round(orig_size[i] * orig_spacing[i] / target_spacing[i]))
    for i in range(3)
]

resampler = sitk.ResampleImageFilter()
resampler.SetOutputSpacing(target_spacing)
resampler.SetSize(target_size)
resampler.SetOutputDirection(img.GetDirection())
resampler.SetOutputOrigin(img.GetOrigin())
resampler.SetInterpolator(sitk.sitkLinear)
resampler.SetDefaultPixelValue(-1000)    # CT 背景填充空气

iso = resampler.Execute(img)
sitk.WriteImage(iso, "ct_iso_1mm.mha")

重采样后要检查三件事：一是体素间距是否真的相等，可用 iso.GetSpacing() 确认；二是骨皮质的 CT 值峰值是否明显下降，若下降超过 100 HU，说明插值过度平滑，可改用最近邻或降低目标分辨率；三是矢状面和冠状面重建图是否还保留层间阶梯，若仍有，说明 z 方向插值不足或原始层厚过大。

五、把五个环节串成可复现的流程

实际工程中，这五个环节不是孤立执行的，而是有明确的先后与依赖关系。推荐顺序是：去噪 → 配准 → 插值 → 重采样 → 各向同性体素化。其中配准和重采样都依赖插值，因此插值方式的选择会同时影响两个环节；各向同性体素化本质上是重采样的目标设定，可以合并执行。

一个可复现的前处理记录至少应包含：原始图像的行列间距与层间距、去噪方法与参数、配准的变换矩阵或位移场、插值方式、目标 spacing 与 size、输出图像的 spacing 校验值，以及骨皮质峰值 HU 的前后对比。这些字段直接决定后续三维重建结果能否被他人复现和比较。

若原始层厚远大于层内间距（如 5 mm 层厚、0.7 mm 像素间距），单纯插值到各向同性并不能恢复真实的 z 方向分辨率，只能改善几何一致性。此时应在报告中注明“z 方向为插值生成”，避免把插值结果当作真实采集精度用于诊断测量。

完成前处理后，体数据才具备进入分割、表面重建和体积测量的条件。判断前处理是否合格，可以回到三个可观察指标：同质区域噪声是否下降、多模态标志结构是否对齐、骨皮质峰值是否保持在可接受范围。三者同时满足，才说明这一轮前处理没有引入新的系统性偏差。