Analysis of Accuracy Related Factors in Photon Counting X-ray CT Energy Imaging
-
摘要:
以光子计数探测器为核心的多能谱成像,成为X射线成像领域人们研究的焦点。本文利用Geant 4仿真软件生成射线源能谱,通过设置不同的能量阈值和探测器能量分辨率,获取不同范围的高低能能谱,扫描物体得到高低能投影数据进行双能重建,比较各组别下的等效原子序数和电子密度标准差。实验结果表明,探测器能量分辨率越高、高低能能谱区分度越大,成像效果越好,而能量阈值对成像精度的影响最终是通过改变能谱重合度,进而对成像精度造成影响。同时,随着能谱重合度变小、能谱宽度变窄,噪声也会加大,因此需要平衡能谱重合度和噪声的关系,才能达到更好的成像效果。
Abstract:Multi-energy spectrum imaging with photon counting detectors has become the focus of X-ray imaging. In this study, the Geant4 simulation software was used to generate the energy spectrum of the ray source, obtaining high-low energy spectra of different ranges by setting different detector energy resolutions and changing the energy threshold. High-low energy projection data of scan objects were obtained for dual-energy reconstruction. The standard deviation of the effective atomic number and electron density under each group were compared. The results show that the higher the energy resolution and the greater the difference between high and low energy spectra, the better is the imaging quality. The energy threshold influenced the imaging accuracy when the spectrum overlap ratio was changed. When the spectrum overlap became smaller and the spectral width became narrower, noise increased, making it necessary to balance the relationship between spectrum overlap and noise to achieve better imaging results.
-
伦琴于1895年发现X射线,而后英国工程师Hounsfield基于X射线发明出计算机断层成像(computed tomography,CT)机。自CT机问世几十年来,无论是在硬件性能(X光机、探测器)、扫描方式,还是重建算法上都获得了长足发展,在医疗、工业和安全检查等领域发挥着重要作用。更短的扫描和重建时间,以及更好的成像质量成为该领域研究人员追求的目标[1]。
双能CT是由Alvarez[2]在1976年首次提出,可以获取物体的等效原子序数和电子密度信息,从而更加准确地识别物质。多能谱CT是在双能CT的基础上,增加了参与扫描和重建的能区数量,引入了更多能谱信息,重建结果更加准确。
随着高计数率电子元器件的迅速发展,以光子计数探测器为核心的多能谱X射线成像技术成为射线检测领域的重要发展方向[3-4]。与传统能量积分探测器相比,光子计数探测器可以通过设置能量阈值来选择不同能区的信号,同时它将每个光子当作独立事件进行统计,一次扫描便能得到多段能谱,减小了辐射剂量[5]。但在实际的应用上也存在多能段光子数减少、信噪比下降,探测器性能不足造成的噪声等问题[6],对影响光子计数X射线CT成像精度的因素及其规律进行基础性研究十分必要。
本文对系统能谱的获取、模体的扫描和重建过程进行仿真实验,以重建得到的等效原子序数和电子密度的标准差作为成像质量的评价标准,得出探测器能量分辨率、能谱阈值和能谱范围等因素对成像精度影响的有关结论,为后续相关技术的工程应用化提供参考。
1. 双能重建基本原理
本文主要研究探测器能量分辨率、能谱阈值和能谱范围等对重建图像精度的影响规律。以双能重建为例,根据光子计数探测原理[7],进行双能成像的模拟。双能CT成像的重点研究内容之一是双能分解方法,常见的双能分解方法包含投影域前处理、图像域后处理和迭代方法[8]。
考虑到未来的工程化需求,本文选择的双能重建方法为基于基材料模型的双能投影预处理分解方法,基材料选取的是碳和铝。
材料模型公式:
$$ \mu \left(E\right)={{b}}_{1}{\mu }_{1}\left(E\right)+{{b}}_{2}{\mu }_{2}\left(E\right) \text{,} $$ (1) $ {\mu }_{1} $ (E)、$ {\mu }_{2} $ (E)分别为两种基材料的线性衰减系数;$ {b}_{1}$ 、$ {b}_{2}$ 对应两种基材料的分解系数,对于某一固定的物质,$ {b}_{1}$ 、$ {b}_{2}$ 是两个常数。双能分解即求解如下的双能投影积分方程组:
$$ \left\{\begin{aligned} {P}_{L}=\;& -{\ln}\bigg(\int \mathit{{S}_{L}} (E)\exp \Big({-{B}}_{1}{\mu }_{1}({E})-\\& {{B}}_{2}{\mu }_{2}({E})\Big){\rm d}{E}\bigg)+{\ln}\int \mathit{{S}_{L}}(E){\rm d}E\\ {P}_{H}=\;&-{\ln}\bigg(\int \mathit{{S}_{H}}(E) \exp \Big(-{{B}}_{1}{\mu }_{1}({E})-\\ & {{B}}_{2}{\mu }_{2}({E})\Big){\rm d}{E}\bigg)+{\ln}\int \mathit{{S}_{H}}(E){\rm d}E\end{aligned}\right. \text{,} $$ (2) 其中,SH、SL分别为高低能能谱,纵坐标记录的是光子数,PH、PL分别为高低能投影,
$ {{B}}_{1}=\displaystyle \int {{b}}_{1}\mathrm{d}{l}, {{B}}_{2}=\displaystyle \int {{b}}_{2}\mathrm{d}{l} $ 。求解
$ {B}_{1}$ 、$ {B}_{2}$ 的这个过程被称为投影分解过程。由于$ {B}_{1}$ 、$ {B}_{2}$ 为$ {b}_{1}$ 、$ {b}_{2}$ 的线积分投影值,求解出$ {B}_{1}$ 、$ {B}_{2}$ 后,根据CT重建的原理,利用滤波反投影重建算法[9],便可计算出$ {b}_{1}$ 和$ {b}_{2}$ ,由此可以计算材质的等效原子序数和电子密度信息,以完成材料的探测识别。计算公式:$$ {Z}_{\rm{eff}}={\left(\frac{{b}_{1}{\rho }_{e1}{Z}_{1}^{n}+{b}_{2}{\rho }_{e2}{Z}_{2}^{n}}{{b}_{1}{\rho }_{e1}+{b}_{2}{\rho }_{e2}}\right)}^{\tfrac{1}{n}} \text{,} $$ (3) $$ {\rho }_{e}={{b}_{1}\rho }_{e1}+{{b}_{2}\rho }_{e2} , $$ (4) 式(3)和式(4)中Z1、Z2分别为两种基材料的原子序数,ρe1、ρe2分别为两种基材料的电子密度,n一般取在3~4之间。
2. 光子计数CT成像系统能谱仿真
在Ubuntu系统下利用Geant 4软件,模拟高速电子撞击靶材射出X射线的过程。真空管内靶材为1 cm厚的钨,靶面倾角为25°,射线源射出能量为160 keV的电子,滤波片为2.5 mm厚的铝,射线穿过区域皆为真空,光子穿过虚拟探测器后被记录能量,得到射线源入射能谱。图1为能谱发射的几何示意,图2为得到的入射能谱。
通过对入射能谱与基于能量分辨率高斯函数的卷积来模拟探测器能量分辨率的影响,得到系统能谱[10]。后续实验所需的高低能能谱均基于系统能谱得到。探测器能量分辨率R通常用半高全宽(full width at half maximum,FWHM)W表示,定义:
$$ R=\frac{W}{E}\text{,} $$ (5) 其中,E是相对应的能量。W可以通过标准差σ与高斯分布的关系计算得到:
$$ W=2.355\;\sigma 。 $$ (6) 因此,标准差可以表示为:
$$ \sigma =\frac{R\times E}{2.355} 。 $$ (7) 根据探测器指定能量下的标准能量分辨率R0和能量E0,可由式(8)获得其他能量下的能量分辨率。由式(7)获得每个能量点下的标准差σ,并使用该标准差对原始入射能谱进行高斯模糊[11]。
$$ R={R}_{0}\sqrt{\frac{{E}_{0}}{E}} 。 $$ (8) 3. 投影数据加噪仿真
由于光子计数的特性遵循泊松分布[12-13],所以本文为投影数据添加泊松噪声,用于模拟真实环境。泊松分布的概率函数见式(9),该函数表示在给定时间间隔或空间区域内发生k次事件的概率,其中,
$ \lambda $ 是泊松分布的均值和方差。$$ P\Big(X=k\Big)=\frac{{\lambda }^{k}{\exp}{(-\lambda )}}{k!} 。 $$ (9) 具体加噪过程。
(1)光子穿过模体后的衰减规律:
$$ {n}=\int S\left(E\right){\exp}{\left(\int -l\mu \left(E\right)\mathrm{d}l\right)}\mathrm{d}E \text{,} $$ (10) 其中,n为入射后探测器接收到的光子数,S(E)为系统能谱,exp为自然常数,l表示射线穿过模体的厚度,μ(E)为水模的线性衰减系数。
将各个单能(E1、E2、E3、···)下的光子数ni=
$ S\left({E}_{\mathrm{i}}\right){\exp}{\left(\int -l\mu \left({E}_{\mathrm{i}}\right)\mathrm{d}l\right)} $ 分别代入式(9)中的$ \lambda $ ,将对应泊松分布下生成的随机数作为该能量下加噪后的光子数ni_noise。(2)将能谱各单能段加噪后的光子数ni_noise累加,并与入射前探测器接收到的光子数n0按下式进行变换,便可得到加噪后的投影Pnoise:
$$ {P}_{\mathrm{n}\mathrm{o}\mathrm{i}\mathrm{s}\mathrm{e}}=\mathrm{l}\mathrm{n}\left(\frac{{n}_{0}}{{n}_{\mathrm{n}\mathrm{o}\mathrm{i}\mathrm{s}\mathrm{e}}}\right) \text{,} $$ (11) 其中,入射前探测器接收到的光子数
$ {n}_{0}= \displaystyle \int S\left(E\right)\mathrm{d}E $ 。4. 实验仿真和结果分析
4.1 实验模体和成像几何
设置的仿真模体为一圆柱状物体,材料为钙(原子序数20,密度1.55 g/cm3),半径为3 cm,模体截面如图3所示。
进行扫描仿真实验的探测系统几何结构如图4所示。射线源为扇形束,扇束张角为22.6°;探测器线性等距排列,共256个探测单元,总长512 mm;射线源、模体、探测器中心位于同一条直线;射线源到模体中心和探测器中心的距离分别为640 mm和
1280 mm;每隔1° 扫描一次模体,每次扫描发射的光子数为105。4.2 能量分辨率的影响
市场上性能相对较好的光子计数探测器材料为碲锌镉(CZT),其在60~100 keV下的能量分辨率大约是10%。因此参考碲锌镉材料的分辨率作为最优值,从10% @ 100 keV到40% @ 100 keV依次递减能量分辨率,间隔10%取4组能量分辨率进行实验,并配以理想能量分辨率(0%)作为对比,高低能能量阈值选为80 keV,各实验组参数见表1。各实验组高斯模糊化之后的高低能能谱见图5。
表 1 能量分辨率实验对照表Table 1. Experimental comparison table of energy resolution编号 能量分辨率 低能区光子数 高能区光子数 1 0% @ 100 keV 82978 17286 2 10% @ 100 keV 82964 17271 3 20% @ 100 keV 82964 17218 4 30% @ 100 keV 82963 17132 5 40% @ 100 keV 82932 17005 穿过模体的光子由所设能量阈值被划入不同的能量通道,从而得到高低能投影数据。为投影数据添加泊松噪声后,进行双能重建[14-15]。
上述实验组的成像结果见表2。不同能量分辨率下的重建结果标准差曲线如图6所示。
表 2 能量分辨率实验成像结果表Table 2. Imaging results of energy resolution experiment编
号能量
分辨率等效原子
序数等效原子
序数标准差电子密度 电子密
度标准差1 0% @
100 keV19.538647 0.268186 1.518764 0.039321 2 10% @
100 keV19.549584 0.291602 1.517544 0.042370 3 20% @
100 keV19.544331 0.295689 1.518037 0.042150 4 30% @
100 keV19.544904 0.313434 1.518503 0.044175 5 40% @
100 keV19.570808 0.350212 1.515818 0.048496 上述结果表明,随着探测器能量分辨率的提高,成像结果受噪声影响程度下降,等效原子序数相较电子密度而言,受噪声影响程度更大。不难分析,光子计数探测器能量分辨率越低,表示其对目标能区光子的统计区域越大,这会使得能谱重叠范围增加,不利于能量分解。但研究发现,能量分辨率在30%甚至40%的情况下,密度和原子序数的标准差变化并不十分明显。
图7是不同探测器能量分辨率的重建结果对比图,电子密度图的灰度显示窗口为[1.44 1.60],等效原子序数图的灰度显示窗口为[18.3 20.7]。图8为同一排像素灰度曲线比较。
4.3 能量阈值及能量范围的影响
本部分考察的是高低能能谱重合程度对于成像的影响,在探测器能量分辨率相同的条件下,通过选取不同能量阈值,改变能量范围,可以得到重合度不一的高低能能谱,重合度的定义为:
$$ {\text{重合度}} = \frac{\text{能谱重叠部分面积}}{\text{全能谱范围面积}}。 $$ (12) 探测器能量分辨率选用30% @ 100 keV,具体对照实验见表3。各实验组能谱对比见图9。
表 3 能量范围实验对照表Table 3. Experimental comparison table of energy range编号 低能范围EL 高能范围EH 重合度Overlap_ratio 低能段光子数 高能段光子数 1 0~70 keV 110~160 keV 0.0025 74978 4419 2 0~80 keV 100~160 keV 0.0134 82323 7438 3 0~80 keV 80~160 keV 0.0648 82963 17132 4 0~100 keV 60~160 keV 0.1526 62561 37533 实验结果见表4,图10为重建结果的均值变化。图11为重建切片图,电子密度图的灰度显示窗口为[1.41 1.61],等效原子序数图的灰度显示窗口为[18.6 20.5]。图12为同一排像素的灰度比较。
表 4 能量范围实验结果Table 4. Results of energy range experiment编号 重合度Overlap_ratio 等效原子序数 等效原子序数标准差 电子密度 电子密度标准差 1 0.0025 19.550129 0.316059 1.517544 0.049246 2 0.0134 19.551358 0.298083 1.517176 0.043675 3 0.0648 19.544904 0.313434 1.518503 0.044175 4 0.1526 19.554059 0.409670 1.517345 0.055390 从结果可以看出,设置能谱阈值改变能谱范围,使得高低能能谱重合度发生改变,会对成像质量造成影响。盲目减小高低能能谱重合度,并不一定使重建结果更好,因为当光子数减少时,噪声会增大。如何平衡能谱重合度和噪声是下一步的研究工作。
5. 小结
光子计数探测器作为多能谱成像的关键器件,成为研究能谱成像的一个重要内容。本文通过实验,模拟光子计数探测器在不同能量分辨率、不同能谱阈值和能量范围下的成像实验。通过比较重建的等效原子序数图和电子密度图的标准差,评价成像结果的好坏,得到上述因素对光子计数X射线CT成像精度的影响规律。
(1)针对安检CT中常用的160 kV射线源能量,光子计数探测器能量分辨率越高,成像效果越好,但是能量分辨率对于材料的分辨差异并不显著。
(2)在探测器能量分辨率确定的情况下,恰当选取能谱阈值对于成像很有帮助。通过设置合适的阈值,使高低能能谱重合度减小,往往可以得到更好的成像结果。但是随着能谱变窄、重合度下降,噪声也会加大,需要平衡这两种因素,否则更低的能谱重合度也有可能带来更差的成像结果。
本文的研究成果将有助于能谱CT系统的总体设计,如从双能成像材料分辨的角度,并不一定过高追求能量分辨率,这对于能量成像CT系统的硬件选型具有积极的指导意义。
-
表 1 能量分辨率实验对照表
Table 1 Experimental comparison table of energy resolution
编号 能量分辨率 低能区光子数 高能区光子数 1 0% @ 100 keV 82978 17286 2 10% @ 100 keV 82964 17271 3 20% @ 100 keV 82964 17218 4 30% @ 100 keV 82963 17132 5 40% @ 100 keV 82932 17005 表 2 能量分辨率实验成像结果表
Table 2 Imaging results of energy resolution experiment
编
号能量
分辨率等效原子
序数等效原子
序数标准差电子密度 电子密
度标准差1 0% @
100 keV19.538647 0.268186 1.518764 0.039321 2 10% @
100 keV19.549584 0.291602 1.517544 0.042370 3 20% @
100 keV19.544331 0.295689 1.518037 0.042150 4 30% @
100 keV19.544904 0.313434 1.518503 0.044175 5 40% @
100 keV19.570808 0.350212 1.515818 0.048496 表 3 能量范围实验对照表
Table 3 Experimental comparison table of energy range
编号 低能范围EL 高能范围EH 重合度Overlap_ratio 低能段光子数 高能段光子数 1 0~70 keV 110~160 keV 0.0025 74978 4419 2 0~80 keV 100~160 keV 0.0134 82323 7438 3 0~80 keV 80~160 keV 0.0648 82963 17132 4 0~100 keV 60~160 keV 0.1526 62561 37533 表 4 能量范围实验结果
Table 4 Results of energy range experiment
编号 重合度Overlap_ratio 等效原子序数 等效原子序数标准差 电子密度 电子密度标准差 1 0.0025 19.550129 0.316059 1.517544 0.049246 2 0.0134 19.551358 0.298083 1.517176 0.043675 3 0.0648 19.544904 0.313434 1.518503 0.044175 4 0.1526 19.554059 0.409670 1.517345 0.055390 -
[1] 陈志强, 张丽, 金鑫. X射线安全检查技术研究新进展[J]. 科学通报, 2017, 62(13): 1350-1364. DOI: 10.1360/N972016-00698. CHEN Z Q, ZHANG I, JIN X. New progress in X-ray safety inspection technology[J]. Chinese Science Bulletin, 2017, 62(13): 1350-1364. DOI: 10.1360/N972016-00698. (in Chinese).
[2] ALVAREZ R E. Extraction of energy-dependent information in radiography[D]. California: Stanford University, 1976.
[3] 郝佳, 张丽, 陈志强, 等. 多能谱X射线成像技术及其在CT中的应用[J]. CT理论与应用研究, 2011, 20(1): 141-150. HAO J, ZHANG L, CHEN Z Q, et al. Multi-energy X-ray imaging technique and its application in computed tomography[J]. CT Theory and Applications, 2011, 20(1): 141-150. (in Chinese).
[4] 邸云霞, 孔慧华, 牛晓伟. 基于主成分分析的多能谱CT图像分析方法研究[J]. CT理论与应用研究, 2022, 31(6): 749-760. DOI: 10.15953/j.ctta.2022.077. DI Y X, KONG H H, NIU X W. Research on image analysis method of spectral CT based on principal component analysis[J]. CT Theory and Applications, 2022, 31(6): 749-760. DOI: 10.15953/j.ctta.2022.077. (in Chinese).
[5] 余子丽. 多能光子技术X-CT的图像重建方法研究[D]. 南京: 南京航空航天大学, 2016. YU Z L. A research on image reconstruction method for multi-energy photon counting X-CT[D]. Nangjing: Nanjing University of Aeronautics and Astronautics, 2016. (in Chinese).
[6] 张煜林. 基于光子计数探测器多能谱CT重建算法及应用研究[D]. 太原: 中北大学, 2017. ZHANG Y L. Multi-spectral CT reconstruction algorithm based on photon count detector and its application[D]. Taiyuan: North University of China, 2017. (in Chinese).
[7] HATEM ALKADHI, ANDRÉ EULER, DAVID MAINTZ, et al. Spectral Imaging[M]. Switzerland: Springer, 2022.
[8] 邸云霞. 基于能谱CT的材料分解方法研究[D]. 太原: 中北大学, 2023. DI Y X. Study on material decomposition method based on energy spectrum CT[D]. Taiyuan: North University of China, 2023. (in Chinese).
[9] AVINASH C K, MALCOLM S. Principles of computerized tomographic imaging[M]. IEEE Press, 1988.
[10] 张斯龙. 提高半导体探测器能量分辨率的信号获取方法研究[D]. 湘潭: 湘潭大学, 2021. ZHANG S L. Research on signal acquisition methodto improve high energy resolution of semiconductor detector[D]. Xiangtan: Xiangtan University, 2021. (in Chinese).
[11] OpenGATE Collaboration. Introduction of GATE [EB/OL]. https://opengate.readthedocs.io/en/latest/index.html, 2024.
[12] 寇松峰, 陈钱, 顾国华, 等. 光子计数成像研究[J]. 应用光学, 2008, 29(5): 675-678. DOI: 10.3969/j.issn.1002-2082.2008.05.003. KOU S F, CHEN Q, GU G H, et al. Research on photon counting imaging technology[J]. Journal of Applied Optics, 2008, 29(5): 675-678. DOI: 10.3969/j.issn.1002-2082.2008.05.003. (in Chinese).
[13] 陈松懋, 郝伟, 苏秀琴, 等. 光子计数成像算法的研究进展[J]. 激光与光电子学进展, 2021, 58(18): 1-16. DOI: 10.3788/LOP202158.1811010. CHEN S M, HAO W, SU X Q, et al. Research progress on photon counting imaging algorithms[J]. Laser & Optoelectronics Progress, 2021, 58(18): 1-16. DOI: 10.3788/LOP202158.1811010. (in Chinese).
[14] 李保磊, 张萍宇, 李斌, 等. X射线双能计算机层析成像投影分解的优化迭代方法[J]. 光学学报, 2017, 37(10): 1-10. DOI: 10.3788/AOS201737.1034001. LI B L, ZHANG P Y, LI B, et al. Optimized iterative method for projection decompositi on of X-ray dual-energy computed tomography[J]. Acta Optica Sinica, 2017, 37(10): 1-10. DOI:10.3788/AOS201737.1034001. (in Chinese).
[15] 庄天戈. CT原理与算法[M]. 上海: 上海交通大学出版社, 1992. ZHUANG T G. CT principle and algorithm[M]. Shanghai: Shanghai Jiaotong University Press, 1992. (in Chinese).