CN105787905A - 一种基于动态电流的锥束ct环状伪影校正方法 - Google Patents

Abstract

本发明公开了一种基于动态电流的锥束CT环状伪影校正方法，克服了现有技术中，重建图像中有环状伪影残留问题。该发明含有以下步骤：步骤1、利用动态电流下探元响应与管电流是否满足线性关系，将探元分为坏点和响应不一致探元两类；步骤2、判断不同管电流下单个探元响应的增量是否为零，对第一类坏点进行检测；步骤3、计算每一个探元的输出响应与管电流的相关系数，利用相关性分析的方法对第二类坏点进行检测；步骤4、计算探元响应和管电流间的线性回归方程，以单个管电流下所有探元响应的均值为基准，计算探元的一致性校正参数矩阵。本发明解决了现有坏点检测方法阈值确定困难的问题，本发明方法对环状伪影校正效果较好，通用性较强。

Description

一种基于动态电流的锥束CT环状伪影校正方法

技术领域

该发明涉及一种锥束CT环状伪影校正方法，特别是涉及一种基于动态电流的锥束CT环状伪影校正方法。

背景技术

锥束CT(ComputedTomography)***主要由X射线源、平板探测器、精密机械和计算机组成，能够在非接触、不破坏的条件下得到物体内部的三维结构信息，已广泛应用于无损检测和医学诊断等领域。实际应用中，受平板探测器制造工艺的影响，在相同射线能量照射下，各探元响应之间存在不同程度的差异，容易导致CT重建图像中出现一系列以重建中心为圆心的同心圆环。这些表现为圆环形状的环状伪影严重干扰了后续对图像的分析、处理等工作。因此研究环状伪影校正方法具有重要的实际意义。

由于平板探测器的探元在一定时期内具有稳定性，已有环状伪影校正方法中探测器校正法的实用性最强。体现在：将探测器的校正参数生成模板，然后每次成像时加载模板进行校正，使用一段时间后更新模板即可。检测探测器坏点进而生成坏点分布模板是探测器校正法的一个关键步骤。

坏点检测方法通常先采集一定曝光水平下未放置任何物体时的探测器输出图像(以下简称曝光图像)，然后计算曝光图像的某些参数，再基于这些参数设置一定的阈值来识别坏点。有方法是直接计算曝光图像中探元响应的均值和方差，而后认为响应的方差超过一定范围的探元即为坏点；还有方法是计算探元响应的增益系数，之后将增益系数超过一定范围的探元标记为坏点；上述两种方法容易受曝光不均匀的影响，且较难确定合适的阈值。Tang等人提出的一种基于小波分析的坏点检测方法效果较好，不受曝光不均匀的影响，但需根据坏点分布情况，将呈点状、簇状或线状分布的坏点分开检测，检测结果显然受坏点分布形式影响且不能准确定位坏点的位置。本发明根据探元响应与管电流线性相关的特点，提出了一种基于动态电流的锥束CT环状伪影校正方法。先对不同管电流下响应恒定的探元进行检测；然后计算每一个探元响应与管电流的相关系数，对响应随机变化和迅速达到最大的探元进行检测。将检测出的探元均作为坏点记入坏点模板，并进行插值校正，将每个管电流下探元响应的均值作为基准，计算同一组管电流下各探元响应曲线方程与基准曲线方程的转化关系，得到全部探元的一致性校正参数矩阵。利用上述方法得到的坏点模板和一致性校正参数，在投影采集时对探元的输出响应进行校正。

发明内容

本发明克服了现有技术中，重建图像中有环状伪影残留问题，提供一种效果显著的基于动态电流的锥束CT环状伪影校正方法。

本发明的技术解决方案是，提供一种具有以下步骤的基于动态电流的锥束CT环状伪影校正方法：含有以下步骤：步骤1、利用动态电流下探元响应与管电流是否满足线性关系，将探元分为坏点和响应不一致探元两类；步骤2、判断不同管电流下单个探元响应的增量是否为零，对第一类坏点进行检测；步骤3、计算每一个探元的输出响应与管电流的相关系数，利用相关性分析的方法对第二类坏点进行检测；步骤4、计算探元响应和管电流间的线性回归方程，以单个管电流下所有探元响应的均值为基准，计算探元的一致性校正参数矩阵。

所述步骤1中动态电流的利用方法如下：固定X射线发生器的管电压，在探测器的最佳动态响应范围内调整管电流，采集当前管电流下的多张探测器输出图像，并取其平均值作为该电流下的探测器输出图像；最终连续设置多个不同的管电流，可以获得一个与该电流序列对应的探测器输出图像序列。

所述步骤1中坏点和响应不一致探元的分类方法如下：由辐射探测器的工作原理可知，在一定范围内探测器产生的信号与入射射线强度成正比；其中，X射线的强度E₀可以用如下公式表示

E₀＝kZIV^m,公式(1)

式中，k是比例系数，Z是靶材料的原子序数，I、V分别表示X射线管电流和管电压，参数m约等于2；由此可以得出固定管电压下，一定曝光时间内探元响应与管电流I间具有线性相关关系

g(x,y)＝a(x,y)I+b(x,y),公式(2)

其中，g(x,y)代表探元(x,y)处的响应，a(x,y)为关系曲线的斜率，b(x,y)为关系曲线的截距，且其中包括了无射线照射时探元的暗电流输出；探元响应与管电流关系偏离线性相关关系的探元即为坏点；响应与管电流关系满足线性相关关系的探元为响应不一致探元。

所述步骤2中第一类坏点的特征如下：探元的输出值恒定，包括当射线管电流增加时输出值恒定为0、最大值或其他常数的探元，这类探元的响应不随管电流的变化而变化。

所述步骤3中第二类坏点的特征如下：探元的输出值变化异常，包括射线管电流连续增加时输出值增加速度过快响应迅速达到最大值和输出值随机变化的探元，这类探元的响应-管电流关系不满足公式(2)。

所述步骤3中第二类坏点的检测方法具体如下：设探测器上探元数量为M×N，X射线管电流的连续递增序列为I＝(I₁,I₂,…,I_n)，生成一组长度为m的伪随机数r对I重新排列得到随机序列I_r＝(I_r1,I_r2,…,I_rm)，则探元(x,y)对应的响应序列接着用皮尔逊相关系数ρ(x,y)来度量I_r与之间的相关程度，具体如下：

式中，是序列I_r的平均值，是序列的平均值，利用式(3)可以得到所有探元对应的相关系数矩ρ，再利用公式(4)进行坏点检测，

式中，q是坏点检测的阈值，矩阵D中取值为0的点即为坏点，根据D中坏点的位置，随后在投影图像中利用插值法对坏点进行校正，步骤2中检测出的第一类坏点也记入模板D。

所述步骤4中探元的一致性校正参数矩阵求解过程如下：对全部坏点进行插值校正后，全部探元响应-管电流关系均满足公式(2)中的方程形式，由式(2)推导可知，任意两探元对应的响应-管电流关系曲线g(x₁,y₁)＝a(x₁,y₁)I+b(x₁,y₁)和g(x₂,y₂)＝a(x₂,y₂)I＋b(x₂,y₂)也具有如下转化关系

g(x₂,y₂)＝g(x₁,y₁)P₁＋Q₁,公式(5)

其中，P₁＝a(x₂,y₂)/a(x₁,y₁)，Q₁＝b(x₂,y₂)-b(x₁,y₁)·a(x₂,y₂)/a(x₁,y₁)，由于全部探元响应的均值代表了探元响应的平均水平，故在一致性校正中以不同管电流下所有探元响应的平均值作为基准，令某电流下所有探元响应的平均值

则序列I对应的探元响应的平均值序列

可以计算得到与I之间的线性回归方程以为基准序列，将其它探元的响应输出g_I(x,y)利用公式(5)向做一致性校正，校正后的值

g′I(x,y)＝g_I(x,y)P(x,y)+Q(x,y),公式(8)

其中，P(x,y)＝a₀/a(x,y)，Q(x,y)＝b₀-b(x,y)·a₀/a(x,y)，这样全部探元响应均被校正到均值水平，由此可以得到所有探元的一致性校正系数矩阵P和矩阵Q。

与现有技术相比，本发明基于动态电流的锥束CT环状伪影校正方法具有以下优点：本发明利用探元响应与X射线管电流之间线性相关这一特点，将探元进行分类检测和校正，提出了一种基于动态电流的环状伪影去除方法，解决了现有坏点检测方法阈值确定困难的问题。在不同探测器平台上开展实际实验表明，本发明方法对环状伪影校正效果较好，通用性较强。

附图说明

图1是本发明基于动态电流的锥束CT环状伪影校正方法中探测器校正流程图；

图2是本发明基于动态电流的锥束CT环状伪影校正方法中不同探测器管电流下探元响应灰度变化曲线图，如PerkinElmer公司型号为XRD1621xNES的探测器；

图3是本发明基于动态电流的锥束CT环状伪影校正方法中不同探测器管电流下探元响应灰度变化曲线图，如Varian公司型号为PaxScan4030E的探测器；

图4是本发明基于动态电流的锥束CT环状伪影校正方法中PerkinElmer公司型号为XRD1621xNES的探测器上的坏点分布；

图5是本发明基于动态电流的锥束CT环状伪影校正方法中Varian公司型号为PaxScan4030E的探测器上的坏点分布；

图6是本发明基于动态电流的锥束CT环状伪影校正方法中PerkinElmer公司型号为XRD1621xNES探测器实验铝柱的重建结果，(a)为传统方法校正；(b)为本发明方法校正；

图7是本发明基于动态电流的锥束CT环状伪影校正方法中Varian公司型号为PaxScan4030E的探测器实验塑料薄膜的重建结果，(a)为传统方法校正后；(b)为本发明方法校正后。

具体实施方式

下面结合附图和具体实施方式对本发明基于动态电流的锥束CT环状伪影校正方法作进一步说明：

探元响应与管电流间的关系：由辐射探测器的工作原理可知，在一定范围内探测器产生的信号与入射射线强度成正比。其中，X射线的强度E₀可以用如下公式表示

E₀＝kZIV^m,公式(1)

式中，k是比例系数，Z是靶材料的原子序数，I、V分别表示X射线管电流和管电压，参数m约等于2。由此可以得出固定管电压下，一定曝光时间内探元响应与管电流I间具有线性相关关系

g(x,y)＝a(x,y)I+b(x,y),公式(2)

其中，g(x,y)代表探元(x,y)处的响应，a(x,y)为关系曲线的斜率，b(x,y)为关系曲线的截距，且其中包括了无射线照射时探元的暗电流输出。在实际应用中，暗电流输出作为噪声，需从投影图像中扣除。

探元的分类：根据探元响应与管电流间的关系分析表明，正常情况下探元响应与管电流线性相关。然而，受平板探测器制造工艺的影响，个别探元对X射线响应存在异常，本发明认为探元响应-管电流关系偏离线性关系的都是坏点，并根据响应的不同将坏点分为如下两类：

(1)输出值恒定的探元，包括当射线管电流增加时输出值恒定为0、最大值或其他常数的探元，这类探元的响应不随管电流的变化而变化；

(2)输出值变化异常的探元，包括射线管电流增加时输出曲线斜率过大响应迅速达到最大值和输出值随机变化的探元，这类探元的响应-管电流关系曲线不满足公式(2)。

此外，还有几种情况探元响应与管电流满足公式(2)的线性关系，但是关系方程不同，这类探元称为响应不一致探元。不一致探元的出现是由于探元存在个体差异，对强度相同的X射线响应不一致，所以表示各探元响应与管电流线性关系曲线的方程也具有随机性。

坏点和响应不一致探元都会导致重建图像中出现环状伪影，但是二者的输出响应与管电流关系不同，需要分别加以校正。

探测器校正：前面对探元进行了响应特点分析，并将探元分成坏点和不一致探元两类，接下来介绍基于坏点的分类检测和一致性校正的探测器校正方法，本发明提出的方法流程见附图1。方法需要对定义的两类坏点分开检测，得到坏点模板D，并根据D中的位置先对坏点进行校正，再计算探元的一致性校正参数矩阵P、Q对响应不一致探元进行校正。

坏点检测及校正：

对于第一类坏点，可以通过判断不同管电流下单个探元响应的增量是否为零进行检测，并在坏点记录模板D中将坏点所在位置记为0。对于第二类坏点，其响应与管电流不满足公式(2)的线性相关关系，针对这类坏点本发明提出一种基于相关性分析的检测方法。

设探测器上探元数量为M×N，X射线管电流的连续递增序列为I＝(I₁,I₂,…,I_n)，生成一组长度为m的伪随机数r对I重新排列得到随机序列I_r＝(I_r1,I_r2,…,I_rm)，则探元(x,y)对应的响应序列

本发明选择用皮尔逊相关系数ρ(x,y)来度量I_r与之间的相关程度，具体如下：

式中，是序列I_r的平均值，是序列的平均值。利用公式(4)可以得到所有探元对应的相关系数矩ρ，再利用公式(5)进行坏点检测。

公式中，q是坏点检测的阈值，实验表明q取0.9能够满足多种不同类型探测器的坏点检测需要。矩阵D中取值为0的点即为坏点，根据D中坏点的位置，在投影图像中利用插值法对坏点进行校正。对全部坏点进行插值校正后，全部探元响应-管电流关系均满足公式(2)中的方程形式，由公式(2)推导可知，任意两探元对应的响应-管电流关系曲线g(x₁,y₁)＝a(x₁,y₁)I+b(x₁,y₁)和g(x₂,y₂)＝a(x₂,y₂)I+b(x₂,y₂)也具有如下转化关系

g(x₂,y₂)＝g(x₁,y₁)P₁+Q₁,公式(5)

其中，P₁＝a(x₂,y₂)/a(x₁,y₁)，Q₁＝b(x₂,y₂)-b(x₁,y₁)·a(x₂,y₂)/a(x₁,y₁)。

一致性校正：坏点校正后全部探元响应-管电流关系均满足公式(2)中的方程形式，由于全部探元响应的均值代表了探元响应的平均水平，故在一致性校正中以不同管电流下所有探元响应的平均值作为基准。令某电流下所有探元响应的平均值

则序列I对应的探元响应的平均值序列

代入公式(2)可知以为基准序列，将其它探元的响应输出g_I(x,y)利用公式(5)向做一致性校正，校正后的值

g′_I(x,y)＝g_I(x,y)P(x,y)+Q(x,y),公式(8)

其中，P(x,y)＝a₀/a(x,y)，Q(x,y)＝b₀-b(x,y)·a₀/a(x,y)，这样全部探元响应均被校正到均值水平，由此可以得到所有探元的一致性校正系数矩阵P和矩阵Q。

为验证本发明方法的有效性和适用性，在几款不同类型探测器上开展了实验，并与基于偏置、增益和探测器坏点模板的传统校正方法分别进行了比较。具体实验参数见表1。

表1实验参数

每个管电流下采集50张曝光图像取平均，得到几款探测器的探元响应随管电流变化关系曲线如附图2所示。可以看出，在设置的各组实验参数下，探元响应的均值与X射线管电流间都满足良好的线性关系。

在坏点检测效果方面，附图3和附图4分别展示了在不同***中利用本发明方法检测得到的几款探测器的坏点分布，图中的白色‘*’表示利用本发明方法检测到的坏点在探测器上的位置。可以看出几款探测器上都存在着大量本发明所述的坏点，且其分布形式多种多样，然而本发明方法对坏点的检测能力并不受坏点分布形式的影响。

为了更直观地观察本发明方法的校正效果，在表1中的实验参数下，实验1、2分别在2.0mA电流下对直径4cm铝柱进行了CT扫描，在200μA电流下对塑料薄膜进行了CT扫描，实验结果见附图5和附图6。

对比附图5和附图6中利用传统方法校正得到的重建结果(a)和利用本发明方法校正得到的重建结果(b)，可以看出本发明方法对由响应不一致探元造成的轻微环状伪影(见附图5)和由坏点造成的明显环状伪影(见附图6(a))均能达到很好的去除效果，相比于传统校正方法对环状伪影校正得更加彻底。

Claims (7)

1.一种基于动态电流的锥束CT环状伪影校正方法，其特征在于：含有以下步骤：

步骤1、利用动态电流下探元响应与管电流是否满足线性关系，将探元分为坏点和响应不一致探元两类；

步骤2、判断不同管电流下单个探元响应的增量是否为零，对第一类坏点进行检测；

步骤3、计算每一个探元的输出响应与管电流的相关系数，利用相关性分析的方法对第二类坏点进行检测；

步骤4、计算探元响应和管电流间的线性回归方程，以单个管电流下所有探元响应的均值为基准，计算探元的一致性校正参数矩阵。

2.根据权利要求1所述的基于动态电流的锥束CT环状伪影校正方法，其特征在于：所述步骤1中动态电流的利用方法如下：固定X射线发生器的管电压，在探测器的最佳动态响应范围内调整管电流，采集当前管电流下的多张探测器输出图像，并取其平均值作为该电流下的探测器输出图像；最终连续设置多个不同的管电流，可以获得一个与该电流序列对应的探测器输出图像序列。

3.根据权利要求1所述的基于动态电流的锥束CT环状伪影校正方法，其特征在于：所述步骤1中坏点和响应不一致探元的分类方法如下：由辐射探测器的工作原理可知，在一定范围内探测器产生的信号与入射射线强度成正比；其中，X射线的强度E₀可以用如下公式表示

E₀＝kZIV^m,公式(1)

式中，k是比例系数，Z是靶材料的原子序数，I、V分别表示X射线管电流和管电压，参数m约等于2；由此可以得出固定管电压下，一定曝光时间内探元响应与管电流I间具有线性相关关系

g(x,y)＝a(x,y)I+b(x,y),公式(2)

其中，g(x,y)代表探元(x,y)处的响应，a(x,y)为关系曲线的斜率，b(x,y)为关系曲线的截距，且其中包括了无射线照射时探元的暗电流输出；探元响应与管电流关系偏离线性相关关系的探元即为坏点；响应与管电流关系满足线性相关关系的探元为响应不一致探元。

4.根据权利要求1所述的基于动态电流的锥束CT环状伪影校正方法，其特征在于：所述步骤2中第一类坏点的特征如下：探元的输出值恒定，包括当射线管电流增加时输出值恒定为0、最大值或其他常数的探元，这类探元的响应不随管电流的变化而变化。

5.根据权利要求1所述的基于动态电流的锥束CT环状伪影校正方法，其特征在于：所述步骤3中第二类坏点的特征如下：探元的输出值变化异常，包括射线管电流连续增加时输出值增加速度过快响应迅速达到最大值和输出值随机变化的探元，这类探元的响应-管电流关系不满足公式(2)。

6.根据权利要求1所述的基于动态电流的锥束CT环状伪影校正方法，其特征在于：所述步骤3中第二类坏点的检测方法具体如下：

设探测器上探元数量为M×N，X射线管电流的连续递增序列为I＝(I₁,I₂,…,I_n)，生成一组长度为m的伪随机数r对I重新排列得到随机序列I_r＝(I_r1,I_r2,…,I_rm)，则探元(x,y)对应的响应序列接着用皮尔逊相关系数ρ(x,y)来度量I_r与之间的相关程度，具体如下：

公式中，是序列I_r的平均值，是序列的平均值，利用公式(3)可以得到所有探元对应的相关系数矩ρ，再利用公式(4)进行坏点检测，

公式中，q是坏点检测的阈值，矩阵D中取值为0的点即为坏点，根据D中坏点的位置，随后在投影图像中利用插值法对坏点进行校正，步骤2中检测出的第一类坏点也记入模板D。

7.根据权利要求1所述的基于动态电流的锥束CT环状伪影校正方法，其特征在于：所述步骤4中探元的一致性校正参数矩阵求解过程如下：对全部坏点进行插值校正后，全部探元响应-管电流关系均满足公式(2)中的方程形式，由公式(2)推导可知，任意两探元对应的响应-管电流关系曲线g(x₁,y₁)＝a(x₁,y₁)I+b(x₁,y₁)和g(x₂,y₂)＝a(x₂,y₂)I+b(x₂,y₂)也具有如下转化关系

g(x₂,y₂)＝g(x₁,y₁)P₁+Q₁,公式(5)

其中，P₁＝a(x₂,y₂)/a(x₁,y₁)，Q₁＝b(x₂,y₂)-b(x₁,y₁)·a(x₂,y₂)/a(x₁,y₁)，由于全部探元响应的均值代表了探元响应的平均水平，故在一致性校正中以不同管电流下所有探元响应的平均值作为基准，令某电流下所有探元响应的平均值

则序列I对应的探元响应的平均值序列

可以计算得到与I之间的线性回归方程以为基准序列，将其它探元的响应输出g_I(x,y)利用公式(5)向做一致性校正，校正后的值

g′_I(x,y)＝g_I(x,y)P(x,y)+Q(x,y),公式(8)

其中，P(x,y)＝a₀/a(x,y)，Q(x,y)＝b₀-b(x,y)·a₀/a(x,y)，这样全部探元响应均被校正到均值水平，由此可以得到所有探元的一致性校正系数矩阵P和矩阵Q。