CN101266216A - 标定双能ct系统的方法和图像重建方法 - Google Patents

Abstract

公开了一种标定双能CT系统的方法和图像重建方法，以计算被扫描物体的原子序数图像和特征密度图像，以及任意能量下的衰减系数图像。本发明还消除了因为射线硬化引起的杯状伪影的影响。该标定双能CT系统的方法包括步骤：选定至少两种不同的材料；探测双能射线对所述至少两种不同的材料在不同厚度组合下的透射射线，以获得投影值；按照所述不同厚度组合和相应的投影值相互对应的方式创建查找表。该图像重建方法包括步骤：利用双能射线扫描被检物体，以获取双能投影值；根据预先创建的查找表计算与所述双能投影值相对应的基材料系数投影值；以及基于所述基材料系数投影值重建基材料系数分布图像。从而，基于基材料系数分布图像可以计算被检物体的原子序数图像、特征密度图像和衰减系数图像。与现有技术相比，本发明提出的方法具有校正过程简单、计算精度高、不受X射线硬化影响的优点。

Description

标定双能CT系统的方法和图像重建方法

技术领域

本发明涉及辐射成像技术，具体涉及一种标定双能CT系统的方法和相应的图像重建方法，它能够消除因为射线硬化引起的杯状伪影。

背景技术

随着技术发展，计算机断层成像(CT)技术开始用到旅客行李检查系统中。常用的CT技术采用X光机作为射线源，产生能够连续分布的X射线。常规的图像重建方法所得到的图像表示物体的衰减系数分布，受到射线硬化的影响，出现了杯状伪影。

现有双能CT图像重建算法通常先采用常规的CT重建方法得到物体的高能和低能衰减系数图像，然后计算得到密度和原子序数图像。现有的这种处理方法由于不能消除射线硬化引起的杯状伪影，导致计算结果不准确并且降低了材料识别的准确度。

发明内容

鉴于上述问题，完成了本发明。本发明的一个目的是提出一种标定双能量CT系统的方法，以及一种图像重建方法。在本发明中，选定基材料，然后制作阶梯块和一系列厚度的矩形块，测量不同厚度组合下的投影值，就可以得到双能查找表，实现系统标定。另外，在利用两种基材料对双能量CT系统进行标定后，用双能CT重建算法可以获取物体的原子序数和密度图像，以及任意能量下的衰减系数图像。

根据本发明的一个方面，提出了一种标定双能CT系统的方法，包括步骤：选定至少两种不同的材料；探测双能射线对所述至少两种不同的材料在不同厚度组合下的透射射线，以获得投影值；按照所述不同厚度组合和相应的投影值相互对应的方式创建查找表。

优选地，所述双能射线包括低能射线和高能射线。

优选地，所述低能射线和所述高能射线是X射线。

优选地，所述至少两种不同的材料包括碳和铝。

根据本发明的另一方面，提出了一种图像重建方法，包括步骤：利用双能射线扫描被检物体，以获取双能投影值；根据预先创建的查找表计算与所述双能投影值相对应的基材料系数投影值；以及基于所述基材料系数投影值重建基材料系数分布图像。

优选地，所述的图像重建方法，还包括步骤：基于基材料系数分布图像计算被检物体的原子序数图像。

优选地，所述的图像重建方法，还包括步骤：基于基材料系数分布图像计算被检物体的特征密度图像。

优选地，所述的图像重建方法，还包括步骤：基于基材料系数分布图像计算被检物体的衰减系数图像。

优选地，如下创建所述查找表：选定至少两种不同的材料；探测双能射线对所述至少两种不同的材料在不同厚度组合下的透射射线，以获得投影值；按照所述不同厚度组合和相应的投影值相互对应的方式创建所述查找表。

与现有技术相比，本发明提出的方法具有校正过程简单、计算精度高、不受X射线硬化影响的优点。

根据本发明的图像重建方法所重建的些图像可以作为安全检查中判断物质属性的依据，以提高安全检查的准确度。

本发明的图像重建方法所取得的图像重建结果的精度高。根据模拟结果，原子序数和密度值的误差在1％以内。

根据本发明，可以获取被检物体在任意能量下的衰减系数图像，不受X射线能谱硬化的影响。

附图说明

从下面结合附图的详细描述中，本发明的上述特征和优点将更明显，其中：

图1是根据本发明实施例的双能CT系统的示意图；

图2是用来说明根据本发明实施例的标定双能CT系统的方法的示意图；

图3是用来描述根据本发明实施例的标定双能CT系统的方法和图像重建方法的流程图；

图4是有机玻璃瓶装水的横截面示意图；

图5示出了根据现有技术重建的结果图像和根据本发明的方法重建的结果图像，其中图5A是采用常规方法重建得到的低能衰减系数图像，图5B是采用本发明实施例的图像重建方法得到的60keV下的衰减系数图像，图5A和图5B的图像的显示灰度窗均为[0.12 0.21]；图5C是根据本发明实施例的方法重建的特征密度图像，显示灰度窗为[0.6 1.12]；图5D是根据本发明实施例的方法重建的原子序数图像，显示灰度窗为[68]；图5E和图5F分别表示沿着图5A与图5B所示的图像的中心线所取得的像素值的曲线；图5G和图5H分别表示沿着图5C和图5D所示的图像的中心线所取得的像素值的曲线以及真实的图像曲线，其中实线表示重建结果，虚线表示真实值。

具体实施方式

下面，参考附图详细说明本发明的优选实施方式。在附图中，虽然示于不同的附图中，但相同的附图标记用于表示相同的或相似的组件。为了清楚和简明，包含在这里的已知的功能和结构的详细描述将被省略，否则它们将使本发明的主题不清楚。

·CT数学原理

将二维分布u(x，y)沿着某个方向θ求线积分，便得到一维的函数p_θ(t)，该函数称为u(x，y)在θ角度的投影。如果能够得到各个方向的投影p_θ(t)，那么可以根据Radon变换精确计算得到二维分布u(x，y)。从投影得到二维分布的过程称为重建。这便是CT的数学原理。

实际应用中，X光机和探测器围绕物体旋转一圈，便测量得到物体的某个切片的衰减系数分布在各个方向的投影，从而可以根据CT原理重建得到物体切片的衰减系数二维分布。

·基材料分解模型

在小型X射线安全检查系统所涉及的能量范围内(＜200keV)，材料的线衰减系数可以用下面的解析表达式(1)来近似表示：

μ(E)＝a₁f_p(E)+a₂f_KN(E)                 (1)

$a_1=\frac{\mathrm{\ρZ}}{M}{Z}^n---\left(2\right)$

$a_2=\frac{\mathrm{\ρZ}}{M}---\left(3\right)$

其中，f_p(E)表示光电效应截面随着能量的变化，f_KN(E)表示康普顿散射截面随着能量的变化，f_p(E)与f_KN(E)均有已知的解析表达式。另外，常数a₁和a₂与材料的原子序数、质量数和密度有关，其表达式如(2)和(3)式所示，其中Z表示原子序数，M表示质量数，ρ表示密度(g/cm³)，n为常数。

由于每种材料的线衰减系数都可以被公式(1)中的两个系数a₁和a₂唯一确定，因此可以选取两种基材料，比如碳和铝，用基材料的线衰减系数的线性组合表示其它所有材料的线性衰减系数，如下式(4)所示：

μ(E)＝b₁μ₁(E)+b₂μ₂(E)                       (4)

其中，μ(E)为任意一种被检查物体的材料的线衰减系数，μ₁(E)和μ₂(E)为所选的两种基材料的线衰减系数，b₁和b₂称为基材料系数。

按照下面的公式(5)，定义原子序数的2倍与质量数的比值与密度的乘积为特征密度：

${\mathrm{\ρ}}^{*}=\mathrm{\ρ}\frac{2Z}{M}---\left(5\right)$

设两种基材料的原子序数和特征密度分别为(Z₁，ρ₁ ^*)和(Z₂，ρ₂ ^*)，那么根据上面的公式(1)～(4)可以推导出任意一种材料的原子序数和特征密度的表达式如公式(6)和(7)所示：

${\mathrm{\ρ}}^{*}=b_1{\mathrm{\ρ}}_1^{*}+b_2{\mathrm{\ρ}}_2^{*}---\left(6\right)$

$Z={\left(\frac{b_1{\mathrm{\ρ}}_1^{*}{Z}_1^n+b_2{\mathrm{\ρ}}_2^{*}{Z}_2^n}{b_1{\mathrm{\ρ}}_1^{*}+b_2{\mathrm{\ρ}}_2^{*}}\right)}^{1/n}---\left(7\right)$

·基材料投影模型

X光管产生的能谱通常为连续谱，探测器对X射线的能量响应函数也不是常数。设能谱N(E)与能量响应函数P_d(E)的乘积为S(E)，并且将S(E)归一化：

${\∫}_0^{E_m}S\left(E\right)\mathrm{dE}=1---\left(8\right)$

那么一条投影线上的投影值的表达式为如下的积分式：

$p=-\ln \frac{I}{I_0}=-\ln {\∫}_0^{E_m}S\left(E\right)\exp \left(\text{-}\underset{l}{\∫}\mathrm{\μ}(E,x,y)\mathrm{dl}\right)\mathrm{dE}---\left(9\right)$

其中，I₀和I分别表示射线被物体衰减前和衰减后的探测器读数值，E_m表示射线的最大能量，l表示射线穿过的路径。

公式(9)即为实际系统的测量投影值p与二维分布μ(x，y)的关系。可以看出，由于X射线多色性，公式(9)并不表示μ(x，y)沿着某条直线的线积分，因而并不满足CT数学原理的要求。常规重建算法忽略这个不一致性，重建得到的μ(x，y)图像含有杯状的伪影，称为硬化伪影。

现有的双能CT方法先采用常规重建算法得到两组μ(x，y)，然后计算原子序数和密度等信息，不能消除射线多色性的影响。本发明采用基材料分解的思路解决了这个问题。

将基材料分解模型代入到公式(9)中，可以得到基于基材料系数的投影值表达式：

$p=-\ln {\∫}_0^{E_m}S\left(E\right)\exp (-{\∫}_l[{\mathrm{\μ}}_1\left(E\right)b_1(x,y)+{\mathrm{\μ}}_2\left(E\right)b_2(x,y)\left]\mathrm{dl}\right)\mathrm{dE}---\left(10\right)$

将公式(10)中的沿着路径l的积分用下面的公式(11)和(12)表达：

∫₁b₁(x，y)dl＝B₁                (11)

∫₁b₂(x，y)dl＝B₂                (12)

这样，按照公式(11)和(12)的定义，称B₁和B₂为基材料系数投影值。设获取到了每个角度下的完整的基材料系数投影值，那么便可以根据CT重建理论重建得到基材料系数b₁和b₂的分布，从而根据基材料分解模型计算出物体的原子序数和特征密度分布，以及任意能量下的线衰减系数值。

·基材料系数投影值的求解

双能CT在两个能量下采集投影数据，便得到如下的双能投影数据：

$p_1(B_1,B_2)=-\ln {\∫}_0^{E_1}{S}_1\left(E\right)\exp [-B_1{\mathrm{\μ}}_1\left(E\right)-B_2{\mathrm{\μ}}_2\left(E\right)]\mathrm{dE}---\left(13\right)$

$p_2(B_1,B_2)=-\ln {\∫}_0^{E_2}{S}_2\left(E\right)\exp [-B_1{\mathrm{\μ}}_1\left(E\right)-B_2{\mathrm{\μ}}_2\left(E\right)]\mathrm{dE}---\left(14\right)$

理论上，在测量得到(p₁，p₂)后，可以根据方程(13)和(14)求解得到(B₁，B₂)。但是，由于上述两个方程均为对数积分方程，无法解析求解。另外，常用的非线性迭代求解方法计算量很大，并且不容易得到稳定解。

但是，发明人注意到，当射线穿过厚度为d₁和d₂的第一基材料和第二基材料后，测量到的双能投影数据如下面的表达式(15)和(16)：

$p_1=-\ln {\∫}_0^{E_1}{S}_1\left(E\right)\exp [-d_1{\mathrm{\μ}}_1\left(E\right)-d_2{\mathrm{\μ}}_2\left(E\right)]\mathrm{dE}---\left(15\right)$

$p_2=-\ln {\∫}_0^{E_2}{S}_2\left(E\right)\exp [-d_1{\mathrm{\μ}}_1\left(E\right)-d_2{\mathrm{\μ}}_2\left(E\right)]\mathrm{dE}---\left(16\right)$

对比公式(13)、(14)和公式(15)、(16)可以看出，如果测量到的投影数据对(p₁，p₂)与{p₁(B₁，B₂)，p₂(B₁，B₂)}相同，那么基材料系数投影值对(B₁，B₂)与基材料的厚度组合(d₁，d₂)是完全相同的。

因此，测量已知材料的不同厚度组合下的双能投影，就可以得到双能投影数据对(p₁，p₂)与基材料系数投影值对(B₁，B₂)的对应关系，从而形成查找表。

在图像重建过程中，根据对被检物体进行双能射线扫描来测量双能投影数据对(p₁，p₂)。然后，根据双能投影数据对(p₁，p₂)从查找表中查找对应的基材料系数投影值对(B₁，B₂)，或者在找到近似的基材料投影值对(B’₁，B’₂)的情况下，通过线性插值来实现对基材料系数投影值对(B₁，B₂)。上述的计算过程明显比通过求解对数方程的过程简单得多。

图1是根据本发明实施例的双能CT系统的示意图。如图1所示，射线源100在控制器500的控制下按照预定的时序产生能量连续分布的双能X射线。被检物体200放在承载机构300上。承载机构300在控制器500的控制下可以匀速旋转，也可以上升和下降。探测器阵列400设置在与射线源100相对的位置，在控制器500的控制下接收透射过被检物体200的透射射线，得到针对第一能量的探测信号和针对第二能量的探测信号。由探测器阵列400检测的信号被转换成数字信号，存储在计算机600中，用于后续的标定或者重建操作。

图2是用来说明根据本发明实施例的标定双能CT系统的方法的示意图。在图2中示出了探测器阵列所采用的双层结构，其中厚度较小的低能晶体410在前面，主要吸收X射线的低能量部分；厚度较大的高能晶体420在后面，主要吸收X射线的高能量部分。低能晶体410和高能晶体420探测的信号被诸如AD转换器这样的辅助电路转换成数字信号。因此，探测器阵列400能够分别输出高能信号和低能信号。

图3是用来描述根据本发明实施例的标定双能CT系统的方法和图像重建方法的流程图。图3的左侧流程图示出了标定双能CT系统的过程，而图3的右侧流程图示出了图像重建方法的详细过程。

如图2中的基材料X和Y所示，选定两种常见的材料作为基材料，比如碳X和铝Y(S110)。将其中一种基材料，例如碳X，制作成阶梯状，另一种基材料，例如铝Y，制作成各种厚度的长方体。根据每一对基材料的厚度值(d₁，d₂)都可以测量得到相应的低能和高能投影值对(p₁，p₂)。

按照图2所示的几何关系，当基材料从上往下经过辐射区域时，探测器阵列400便测量得到了基材料Y的某个厚度与基材料X的一系列厚度所对应的双能投影值。改变基材料Y的厚度，重复上面的测量，便可以得到X和Y在各个厚度组合下的双能投影值。所有的测量结果形成了双能投影值对与基材料厚度组合的对应关系(S120)。

在检查过程中，首先将被检查物体放在承载机构300上，然后控制器500控制射线源200发出双能X射线在各个角度下照射被检查物体，并由探测器阵列400得到双能投影值对(S210)。

接下来，利用上述创建的查找表计算双能投影值对所对应的厚度组合(d₁，d₂)。这样，就可以得到基材料系数投影值对(B₁，B₂)(S220)。然后，根据CT重建算法便得到了基材料系数分布b₁和b₂的分布图像(S230)。

另外，根据基材料系数图像计算得到物体的原子序数和特征密度图像，以及任意能量下的衰减系数图像(S240)。

为了验证上述重建方法，进行了如下的数值模拟实验。设X光机高压为140kV，探测器阵列采用CsI晶体，采用蒙特卡罗方法模拟得到能谱和探测器的能量响应函数，然后用公式(9)解析计算双能投影值。

例如，选定碳和铝作为基材料，原子序数分别为6和13，质量数为12.011和26.9815，密度为1g/cm³和2.7g/cm³，特征密度为0.999084g/cm³和2.601783g/cm³。碳材料的厚度序列为0～10cm，间隔为1cm；铝材料的厚度序列为0～1cm，间隔0.1cm。不同厚度组合下的双能投影值采用公式(15)和(16)计算得到。这样碳和铝材料均有11个厚度，形成11×11的查找表。

图4是有机玻璃瓶装水的横截面示意图。被检物体为装水的有机玻璃瓶，瓶壁厚度为5mm，瓶外直径为160mm，瓶内直径为150mm。有机玻璃的原子序数为6.56，密度为0.8g/cm³，特征密度为0.863g/cm³。水的原子序数为7.51，密度为1.0g/cm³，特征密度为1.11g/cm³。

在CT扫描过程中，采用平行束扫描，投影角度数为720，探测器个数为512，重建图像大小为521×512。

图5示出了根据现有技术重建的结果图像和根据本发明的方法重建的结果图像，其中图5A是采用常规方法重建得到的低能衰减系数图像，图5B是采用本发明实施例的图像重建方法得到的60keV下的衰减系数图像，图5A和图5B的图像的显示灰度窗均为[0.12 0.21]；图5C是根据本发明实施例的方法重建的特征密度图像，显示灰度窗为[0.6 1.12]；图5D是根据本发明实施例的方法重建的原子序数图像，显示灰度窗为[68]；图5E和图5F分别表示沿着图5A与图5B所示的图像的中心线所取得的像素值的曲线；图5G和图5H分别表示沿着图5C和图5D所示的图像的中心线所取得的像素值的曲线以及真实的图像曲线，其中实线表示重建结果，虚线表示真实值。

比较图5A和5B、图5C和5D，可以看到本发明实施例的方法得到的衰减系数图像可以消除射线硬化引起的杯状伪影。另外，从图5G和图5H中可以看到重建结果与真实值的差异非常小，这表明本方法可以获得很高的精度。

如上所述，根据本发明实施例的标定方法降低了复杂度。具体来说，选定多个基材料，然后制作阶梯块和一系列厚度的矩形块。通过测量不同厚度组合下的投影值，就可以得到双能查找表，实现系统标定。

另外，根据本发明实施例的图像重建方法所取得的图像重建结果精度高。根据模拟结果，原子序数和密度值的误差在1％以内。

另外，根据本发明的方法，可以获取被检物体在任意能量下的衰减系数图像，不受X射线能谱硬化的影响。

上面的描述仅用于实现本发明的实施方式，本领域的技术人员应该理解，在不脱离本发明的范围的任何修改或局部替换，均应该属于本发明的权利要求来限定的范围，因此，本发明的保护范围应该以权利要求书的保护范围为准。

Claims (9)

1. 一种标定双能CT系统的方法，包括步骤：

选定至少两种不同的材料；

探测双能射线对所述至少两种不同的材料在不同厚度组合下的透射射线，以获得投影值；

按照所述不同厚度组合和相应的投影值相互对应的方式创建查找表。

2. 如权利要求1所述的方法，其特征在于，所述双能射线包括低能射线和高能射线。

3. 如权利要求2所述的方法，其特征在于，所述低能射线和所述高能射线是X射线。

4. 如权利要求1所述的方法，其特征在于，所述至少两种不同的材料包括碳和铝。

5. 一种图像重建方法，包括步骤：

利用双能射线扫描被检物体，以获取双能投影值；

根据预先创建的查找表计算与所述双能投影值相对应的基材料系数投影值；以及

基于所述基材料系数投影值重建基材料系数分布图像。

6. 如权利要求5所述的图像重建方法，还包括步骤：

基于基材料系数分布图像计算被检物体的原子序数图像。

7. 如权利要求5所述的图像重建方法，还包括步骤：

基于基材料系数分布图像计算被检物体的特征密度图像。

8. 如权利要求5所述的图像重建方法，还包括步骤：

基于基材料系数分布图像计算被检物体的衰减系数图像。

9. 如权利要求5～8之一所述的图像重建方法，如下创建所述查找表：

选定至少两种不同的材料；

探测双能射线对所述至少两种不同的材料在不同厚度组合下的透射射线，以获得投影值；

按照所述不同厚度组合和相应的投影值相互对应的方式创建所述查找表。

Images