基于SART算法的CL硬化伪影校正方法研究

曹大泉，王雅霄，阙介民，孙翠丽，魏存峰，魏 龙

(1.中国科学院 高能物理研究所 核分析技术重点实验室，北京 100049；

2.北京市射线成像技术与装备工程技术研究中心，北京 100049)

X射线薄板层析成像(CL)系统是一种非同轴扫描的计算机断层成像系统，专门针对板状构件进行结构成像和缺陷检测。CL与常规锥束CT的差异在于转台的转轴不再垂直于X射线源焦点与探测器中心的连线(称为光轴)，而是成一定的倾斜角度。CL的扫描方式解决了CT扫描板状构件时的两个难题：1) 某些角度穿透困难，图像信噪比低；2) 高空间分辨率下无法完成扫描。另一方面，这种特殊的扫描结构使得CL的投影数据不满足数据完备性条件，因此，传统的解析重建方法容易出现不同断层间图像信息混叠。

实际数据获取系统与图像重建数学模型之间的差异导致产生各种图像伪影。由于CL的射线源同样采用X光管，而非重建模型中的单色光，因此导致产生射束硬化伪影。硬化伪影的校正方法很多，如滤片预硬化法、多项式拟合法、蒙特卡罗法、双能校正法等。这些校正方法均对硬化伪影有一定的校正效果，不过都有一定的限制条件。CL检测的目标样品形状多不规则(如电路板的焊点)，不易提取准确边界，不适合采用重投影多项式拟合法[1]；材质多样且多为复合材料，缺乏先验知识，不宜使用蒙特卡罗法模拟校正[2]；滤片预硬化法简单实用，但校正效果不及其他方法理想[3]；双能校正法成本高，在实验室不易实现。考虑到射束硬化伪影广泛存在于CL重建图像中，严重降低了图像质量，因此需要一种不同于上述射束硬化伪影校正的方法来校正CL图像的伪影。文献[4]提出一种校正常规CT图像射束硬化伪影的方法，该方法在SART重建算法中加入射束硬化物理过程，通过准确表达射线投影获取模型来校正CT图像的硬化伪影。本工作进一步研究该校正方法，将该方法用于CL图像的硬化伪影校正。

1 校正方法及数据获取

1.1 CL硬化伪影校正方法

SART算法是一种对ART算法进行改进的代数迭代重建算法[5]。与ART对每条射线依次更新不同，SART对每个投影角度下所有的射线同时更新，因此重建速度更快。SART迭代重建首先将被重建区域赋初值0，对该区域做前投影，将前投影与真实投影的差值进行权重的反投影。然后将反投影的结果迭加到初始被重建区域，其体素更新过程可表示为：

λ

(1)

考虑到SART在更新每个重建点像素时，将每条射线上的衰减值按权重进行累加，权重值即为射线在每个重建网格中的交线长度，因此可将射束硬化校正的模型引入SART重建的过程，称该种图像重建算法为SART-BHC算法。SART-BHC使用射线在重建点的有效衰减系数累加代替衰减系数的直接累加，可更准确地反映多色光谱投影数据获取的过程，进而在重建图像时抑制射束硬化伪影。有效衰减系数与多色光谱、射线在物体中穿过的长度及物体本身材质相关，可表示为[6]：

(2)

其中：α为校正强度系数；β为能量相关系数；Eeff为有效衰减系数。

将式(2)代入SART迭代公式可得SART-BHC算法更新公式[4]：

(3)

其中，N为重建体素个数。

(4)

(5)

其中：SO为光源到物体旋转中心的距离；OD为物体旋转中心到探测器中心的距离；φ为中心射线与物体旋转轴的夹角；θ为物体绕z″的旋转角；di为第i条射线与探测器的交点。

在物体坐标系下，得到每个投影角度下射线源及每个探测器像素的位置坐标，就可利用射线驱动的Siddon算法[7]计算系统矩阵ω。

与解析重建算法(如滤波反投影(FBP)重建算法)需要完备的投影数据相比，SART-BHC在投影数据不完备的情况下也能得到较好的重建结果，因此更适用于CL这种无法获取完备投影数据的扫描系统。SART-BHC的射束硬化伪影校正方法通过将有效衰减系数引入迭代过程，不需被测物体和光谱的信息，不需对被测物体进行刻度实验以建立硬化校正模型。与标准SART重建相比，不用增加额外的迭代次数。此外，由于以往锥形束CT射束硬化伪影校正方法无法适用于CL扫描结构[8]，该校正方法可很好地应用于CL这种特殊的扫描结构。由于很多工业部件无法获取准确的成分信息，或无法进行刻度实验以建立硬化伪影校正模型，本文提出的方法可很好满足工业部件硬化伪影校正的需要。

图1 CL系统几何示意图

1.2 模拟数据获取

模拟数据不受实验硬件条件不足的限制，又可排除实验中各种因素的干扰，可方便准确地验证校正算法。本文采用蒙特卡罗方法在Geant4环境下搭建CL扫描系统模型，获取投影数据以验证算法的有效性。

首先，在Geant4环境中利用电子打靶获取X射线光谱，模拟参数列于表1。模拟空间中的几何关系示于图2。

表1 模拟电子打靶主要参数

图2 模拟电子打靶的几何关系

高速电子(160 keV)轰击在钨靶上，大部分电子与钨靶作用后散射出去，极少部分电子与钨靶作用后产生X射线，将探测器收集到的X射线绘制成能谱，结果示于图3。

图3 160 keV电子打钨靶产生的X光谱

图4 CL扫描模型

按图3光谱的分布发射一定数量的X光子，就相当于在模拟环境中建立一连续谱分布的X射线源。将该光谱作为X射线源，在Gate环境中建立CL扫描模型(图4)，模型参数列于表2。

表2 CL投影模型主要参数

利用上述模型，样品绕转轴旋转360°，获取间隔均匀的60幅投影图，图5示出其中1幅。

图5 模拟CL扫描的投影图

1.3 实验数据获取

a——FXE-225.99型微焦点X光机；b——PaxScan® 2520V型探测器；c——Huber转台；d——偏转转台

本实验在1台225 kV高分辨的CL系统上进行，实验平台示于图6。CL系统搭建在光学平台上，主要包括1台FXE-225.99型微焦点X光机，1台PaxScan® 2520V型探测器，1台加装1个很薄的有机玻璃圆筒的Huber转台，以及位于Huber转台下方用于控制偏转角的偏转转台。为了与模拟数据更好符合，加工一铝质圆盘，圆盘直径30 mm，高10 mm。调整CL系统的参数(表3)，对铝盘进行CL扫描，得到图7所示的投影图。

表3 实验参数

图7 实验投影图

2 结果分析

2.1 模拟数据校正结果

为证明SART-BHC校正方法的可行性，作为对比，将蒙特卡罗模拟获得的投影利用FBP算法重建成断层图。FBP是一种最普通的重建算法，由于其算法简单、重建速度快，广泛应用于CT图像重建领域。FBP的重建参数列于表4，重建的CL结果示于图8。由图8可看出，CL扫描重建出的断层图像有明显的杯状伪影。灰度曲线的底部较平，分析可知，这个杯状伪影不仅来自于硬化，还有一部分来自于层间混叠效应，即CL的倾斜扫描方式使得图像重建时不同层之间的信息相互干扰。层间混叠效应在本实验中体现在把本属于铝盘内部体素对投影的贡献，重建图像时一部分丢失在铝盘外部，这样产生层间混叠伪影，同时降低了内部的灰度，加重了杯状伪影。

表4 模拟数据的FBP算法重建参数

a——铝盘的断层图像；b——断层图像灰度曲线

利用SART-BHC算法重建蒙特卡罗模拟的投影数据，重建参数列于表5，SART-BHC算法重建结果示于图9。由图9可见，模体的断层图像不再存在杯状伪影，灰度曲线的趋势变为平直，与实际情况相符。因此，应用SART-BHC算法对CL扫描方式所产生的杯状伪影有很好的校正效果，证明对该算法改进的正确性及该算法对CL扫描方式的杯状伪影校正的有效性。

表5 模拟数据的SART-BHC算法重建参数

a——铝盘的断层图像；b——断层图像灰度曲线

2.2 实验数据校正结果

为进一步证明SART-BHC校正方法在CL中的实用性，对CL实验平台上获取的实验数据进行校正研究。首先将在CL实验平台上获取的实验数据利用FBP算法重建断层图像，重建参数与表4所列参数相同，重建结果示于图10。为与模拟数据直观对照，从整幅断层图像的正中选取一部分，并压缩至与模拟图像同等大小。从图10可看出，实验数据的FBP重建结果与图8所示的模拟数据结果一致，同样有严重的杯状伪影，严重降低了图像质量，应予以校正。

a——铝盘的断层图像；b——断层图像灰度曲线

利用SART-BHC算法对CL的实验投影数据进行重建，主要重建参数列于表6，校正结果示于图11。从图11的重建结果可看出，铝质圆盘的断层图像不再存在杯状伪影，灰度曲线较图10所示的结果更加均匀，与实际相符。可见，SART-BHC应用在CL中可有效地校正杯状伪影，得到更好的图像质量。

表6 实验数据的SART-BHC算法重建参数

a——铝盘的断层图像；b——断层图像灰度曲线

3 结论

本文利用蒙特卡罗模拟及CL实验平台分别获取了CL扫描方式的模拟和实验的投影数据，并利用基于SART迭代重建的硬化伪影校正算法(SART-BHC算法)对两组投影数据分别进行重建，均得到了硬化伪影消除的断层图像。具体结论如下。

1) CL成像系统同样存在射束硬化现象，重建图像存在杯状伪影。杯状伪影部分来自射束硬化，也有一部分来自层间混叠效应。

2) 通过对SART-BHC的研究及对CL系统矩阵的推导，实现了应用于CL的改进SART-BHC算法。该算法通过将有效衰减系数引入迭代过程，在不增加重建迭代次数的情况下，可以很好消除CL图像的射束硬化伪影。

3) 与以往的硬化伪影校正算法不同，SART-BHC在校正过程中不需被测物体和光谱信息，不需对被测物体进行刻度实验以建立硬化校正模型。由于很多工业部件无法获取准确的成分信息，或无法进行刻度实验，SART-BHC可很好地满足工业部件硬化伪影校正的需要。

4) SART-BHC在校正射束硬化伪影的同时，能很好地适用于CL这种投影数据采样不完备的扫描几何，因此可较好地抑制CL层间混叠。

参考文献：

[1] ZHANG Q H, LU H N, YANG M. Beam hardening correction using polynomial fitting based on reprojection[J]. Optical Technique, 2005, 31(4): 633-638.

[2] 曾刚，郁忠强，阎永廉. 基于蒙特卡罗模拟的射束硬化校正方法[J]. 高能物理与核物理，2006，30(2)：178-182.

ZENG Gang, YU Zhongqiang, YAN Yonglian. Beam hardening correction based on Monte Carlo simulation[J]. High Energy Physics and Nuclear Physics, 2006, 30(2): 178-182(in Chinese).

[3] 赵维. 锥束CT射束硬化伪影校正[D]. 北京：中国科学院高能物理研究所，2012.

[4] BRABANT L, PAUWELS E, DIERICK M, et al. A novel beam hardening correction method requiring no prior knowledge, incorporated in an iterative reconstruction algorithm[J]. NDT & E International, 2012, 51: 68-73.

[5] KAK A C, SLANEY M. Principles of computerized tomographic imaging[M]. New York: IEEE Press, 1988.

[6] QUE J M, CAO D Q, ZHAO W, et al. Computed laminography and the reconstruction algorithm[J]. Chinese Phys C, 2012, 36(8): 777-783.

[7] SIDDON R. Fast calculation of the exact radiological path length for a three-dimensional CT array[J]. Med Phys, 1985, 12: 252-255.

[8] ZHAO W, FU G T, SUN C L, et al. Beam hardening correction for a cone-beam CT system and its effect on spatial resolution[J]. Chinese Phys C, 2011, 35(10): 978-985.