CN103106676B - 一种基于低剂量投影数据滤波的x射线ct图像重建方法 - Google Patents

Abstract

一种基于低剂量投影数据滤波的X射线CT图像重建方法,包括：(1)获取CT设备的***参数和低剂量扫描协议下的投影数据；(2)对步骤（1）中所获取的投影数据进行方差估计；(3)进行滤波器设计，对步骤（1）中获取的投影数据进行滤波，获得滤波后的投影数据；(4)对步骤（3）中获得的滤波后的投影数据进行CT图像重建，得到最终的CT图像。步骤（3）具体包括：进行滤波器形式设计（I）；将步骤（1）中的投影数据代入公式（I）循环执行公式（I），达到预设的次数时停止，所得到的迭代运算结果为最终滤波后的投影数据。该方法可在大幅减少图像噪声和伪影的同时较好地保持原有图像的分辨率，最终实现CT图像的优质重建。

Description

一种基于低剂量投影数据滤波的X射线CT图像重建方法

技术领域

本发明涉及一种医学影像的图像重建方法，具体涉及一种基于投影数据滤波的X射线CT图像重建方法。

背景技术

X射线CT扫描已经广泛应用于临床医学影像诊断，但是CT扫描过程中过高的X射线辐射剂量会存在致癌风险。为了降低对使用者的损害,如何最大限度地降低X射线使用剂量已经成为医学CT成像领域研究的关键技术之一。

为了降低X射线辐射剂量，现有技术中使用的最简便的途径就是降低CT扫描过程中的管电流和扫描时间。然而，由于降低了管电流和扫描时间，就会使得投影数据中含有大量的噪声，基于传统的滤波反投影方法重建的图像质量存在严重的退化现象，难以满足临床诊断需要。

为了在保证图像质量的前提下大幅降低X射线辐射剂量，诸多基于降低管电流和扫描时间的低剂量CT图像重建方法相继提出，例如基于统计模型的迭代重建方法和基于投影数据滤波的解析重建方法。其中，基于统计模型的迭代重建方法，通过对采集的投影数据的噪声以及成像***进行图像重建模型构建，可以实现低剂量CT图像优质重建；基于投影数据滤波的解析重建方法，是通过对采集的投影数据的噪声以及成像***进行数据滤波建模，再通过解析重建方法实现快速且优质的低剂量CT图像重建。

在重建同样大小的CT图像时,基于统计模型的迭代重建方法需要对目标函数进行几十甚至上百次的反复迭代求解，导致图像重建时间大幅增加。所花费的时间远远大于传统的解析重建方法，不能满足临床的所需要的CT实时显像要求。

而现有技术中的基于投影数据滤波的解析重建方法，在投影数据降噪过程中不可避免地会导致图像原有细节信息的丢失，从而导致相应CT图像分辨率的下降。

因此，针对现有技术不足，提供一种基于低剂量投影数据滤波的X射线CT图像重建方法,能够针对降低管电流和扫描时间的低剂量CT投影数据进行图像重建,所重建的图像质量好,能够保持原有图像的分辨率。

发明内容

本发明提供一种基于低剂量投影数据滤波的X射线CT图像重建方法,能够针对降低管电流和扫描时间的低剂量CT投影数据进行图像重建,所重建的图像质量好,能够保持原有图像的分辨率。

本发明的上述目的通过如下技术手段实现。

一种基于低剂量投影数据滤波的X射线CT图像重建方法，依次包括如下步骤：

(1)获取CT设备的***参数和低剂量扫描协议下的投影数据y_raw；

(2)对步骤(1)中所获取的投影数据y_raw进行逐个数据点上的方差估计；

(3)进行滤波器设计，对步骤(1)中获取的投影数据y_raw进行滤波，获得滤波后的投影数据y_filtered；

(4)对步骤(3)中获得的滤波后的投影数据y_filtered进行CT图像重建，得到最终的CT图像。

上述步骤(3)具体包括：

(3.1)进行滤波器形式设计，为：

$y_{filtered}^{n+1}=w\·y_{raw}+(1-w)\·F\left(y_{filtered}^n\right)---\left(I\right);$

其中y_raw为步骤(1)中采集的投影数据；与分别为第n次和第n+1次滤波后的投影数据，n为自然数，表示迭代运算的次数；w为权重因子；

(3.2)将步骤(1)中的投影数据y_raw作为初始值代入公式(I)中进行迭代运算，得到滤波后的投影数据然后以同样的方法依次进行下一次迭代运算，依次不断地循环执行公式(I)，当循环次数达到预设的次数时即停止迭代运算，并将所得到的迭代运算结果作为最终的滤波后的投影数据y_filtered。

其中，迭代运算的次数n范围设置为大于等于5而小于等于30。

优选的，迭代运算的次数n设置为9或者10。

其中，权重因子w为方差的负相关函数。

权重因子w与方差的关系为或者或者其中β和ε表示正实数。

上述步骤(1)中CT设备的***参数包括X射线入射光子强度、***电子噪声的方差

上述步骤(2)采用基于小邻域图像的局部方差估计方法或者基于CT投影数据噪声特性的方差估计方法对步骤(1)获取的投影数据y_raw进行逐个数据点上的方差估计。

具体的，方差的估计公式为：其中，I₀为待计算方差的第i个数据点的X射线入射光子强度；为***电子噪声的方差。

上述步骤(4)具体是采用滤波反投影法或者迭代重建法进行CT图像重建。

本发明提供一种基于低剂量投影数据滤波的X射线CT图像重建方法,依次包括如下步骤：(1)获取CT设备的***参数和低剂量扫描协议下的投影数据y_raw；(2)对步骤(1)中所获取的投影数据y_raw进行逐个数据点上的方差估计；(3)进行滤波器设计，对步骤(1)中获取的投影数据y_raw进行滤波，获得滤波后的投影数据y_filtered；(4)对步骤(3)中获得的滤波后的投影数据y_filtered进行CT图像重建，得到最终的CT图像。该方法通过对CT投影数据的噪声设计滤波(I)式中的权重因子，将滤波后投影数据与原始采集的投影数据进行加权融合，有效解决了一般性降噪技术进行投影数据滤波所导致的噪声去除与分辨率保持的矛盾性问题，可以在大幅减少图像噪声和伪影的同时较好地保持原有图像的分辨率，最终实现CT图像的优质重建。

附图说明

利用附图对本发明作进一步的说明，但附图中的内容不构成对本发明的任何限制。

图1为本发明方法的流程示意图。

图2为本发明实施例2中仿真所使用的Shepp-Logan体模数据示意图。

图3为本发明实施例2中模拟生成的低剂量CT投影数据；其中，(a)为模拟生成的低剂量CT投影数据y_raw图像；(b)为本发明方法滤波后的投影数据图像；(c)为图像(a)与(b)的差图。

图4为本发明实施例2中由滤波处理后的低剂量投影数据FBP重建的图像：其中，(a)为本发明方法滤波后投影数据的重建结果；(b)为一般的均值滤波方法滤波后投影数据的重建结果，均值滤波窗为3×3。

图5为由模拟生成的低剂量投影数据直接通过FBP法重建的图像。

具体实施方式

结合以下实施例对本发明作进一步描述。

实施例1。

一种基于低剂量投影数据滤波的X射线CT图像重建方法，依次包括如下步骤。

(1)获取CT设备的***参数和低剂量扫描协议下的投影数据y_raw。

(2)对步骤(1)中所获取的投影数据y_raw进行逐个数据点上的方差估计。

(3)进行滤波器设计，对步骤(1)中获取的投影数据y_raw进行滤波，获得滤波后的投影数据y_filtered。

(4)对步骤(3)中获得的滤波后的投影数据y_filtered进行CT图像重建，得到最终的CT图像。

具体的，步骤(1)中CT设备的***参数包括X射线入射光子强度I₀、***电子噪声的方差等，通常这些数值可以直接从测试设备读取或者得到。

步骤(2)采用基于小邻域图像的局部方差估计方法或者基于CT投影数据噪声特性的方差估计方法对步骤(1)获取的投影数据y_raw进行逐个数据点上的方差估计。

方差的估计公式可以为：其中，I₀为待计算方差的第i个数据点的X射线入射光子强度；为***电子噪声的方差。

需要说明的是，方差的估计公式并不局限于本实施例中的这种情况，也可以根据需要灵活设置或者选择。

步骤(3)具体包括：

(3.1)进行滤波器形式设计，为：

$y_{filtered}^{n+1}=w\·y_{raw}+(1-w)\·F\left(y_{filtered}^n\right)---\left(I\right);$

其中y_raw为步骤(1)中采集的投影数据；与分别为第n次和第n+1次滤波后的投影数据，n为自然数，表示迭代运算的次数；w为权重因子。F(·)为传统的滤波器，例如均值滤波、中值滤波，非局部平均滤波等。

(3.2)将步骤(1)中的投影数据y_raw作为初始值代入公式(I)中进行迭代运算，得到滤波后的投影数据然后以同样的方法依次进行下一次迭代运算，依次不断地循环执行公式(I)，当循环次数达到预设的次数时即停止迭代运算，并将所得到的迭代运算结果作为最终的滤波后的投影数据y_filtered。

通常，本发明的方法迭代运算的次数n范围设置为大于等于5而小于等于30，优选设置为9次或者10次。不仅能够有效满足图像处理的要求，而且迭代次数少，图像重建时间短。

权重因子w为方差的负相关函数，二者之间的关系为或者或者其中β和ε表示正实数。

步骤(4)具体是采用滤波反投影法或者迭代重建法进行CT图像重建。

本发明所公开方法，相比一般的投影数据滤波方法，充分考虑了CT投影数据的噪声特性，能够准确进行建模，并根据噪声方差设计权重将滤波后投影数据与原始采集的投影数据进行加权融合，解决了一般性单纯滤波降噪强度过大引起的图像细节丢失。实践证明，通过本发明的方法所实现的CT重建图像，在大幅减少图像噪声和伪影的同时能够较好地保持原有图像的分辨率。

实施例2。

以计算机仿真的数字体模数据为例来描述本发明的基于低剂量投影数据滤波的X射线CT图像重建方法的具体实施过程。

如图1所示，本发明一种基于低剂量投影数据滤波的X射线CT图像重建方法，具体包括下列步骤。

采用图2所示的Shepp-Logan数字体模图像作为本发明的计算机仿真实验对象。体模图像像素矩阵大小设为512×512，模拟CT机的X射线源到旋转中心和探测器的距离分别为570mm和1040mm，旋转角在[0,2π]间，采样值为1160，每个采样角对应672个探测器单元，探测器单元的大小为1.407mm。

首先进行步骤(1)，通过CT***仿真生成大小为1160×672的投影数据y_raw，如图3(a)所示，其中X射线的入射光子强度I₀为2.5×10⁵，***电子噪声的方差为11.0。在实际的CT数据采集中，投影数据和***参数即入射光子强度和***电子噪声的方差均可以直接获取。

(2)接着，对步骤1中模拟生成的CT投影数据y_raw和***参数I₀和进行逐个数据点上的方差估计。

方差的估计公式为：其中，I₀为第i个数据点的X射线入射光子强度，即I₀＝2.5×10⁵；为***电子噪声的方差，即 ${\mathrm{\σ}}_e^2=\mathrm{11.0.}$

(3)接着利用本发明所示的滤波方法进行滤波。

根据迭代公式(I)： $y_{filtered}^{n+1}=w\·y_{raw}+(1-w)\·F\left(y_{filtered}^n\right);$

权重因子参数β＝0.001；滤波器F(·)为3×3邻域的均值滤波，其形式为N_i为第i个数据点周围3×3的邻域。

然后将步骤1中的投影数据y_raw作为初始值代入公式(I)中进行迭代运算，得到滤波后的投影数据然后以同样的方法依次进行下一次迭代运算，依次不断地循环迭代运算，当循环次数达到10次时，得到最终的滤波后的投影数据如图3(b)所示。

(4)最后对对滤波的投影数据使用传统扇形束滤波反投影(FilteredBack-Projection,FBP)算法进行重建，得到重建图像，如图4(a)所示。

为了对比本发明所示方法的效果，对模拟生成的同组CT投影数据，采用一般的均值滤波技术进行滤波处理，然后对滤波的投影数据采用传统扇形束FBP算法进行重建，得到重建图像，如图4(b)所示。

将图4中(a)和(b)两种滤波方法对应的重建结果进行比较可见，本发明方法能够在保持图像分辨率的同时，更好地减少图像中的噪声和伪影，如图4箭头所示。

同时，为了说明投影数据滤波对传统扇形束FBP算法重建图像质量的影像，图5给出了由CT投影数据直接采用传统扇形束FBP算法重建得到的图像。比较图4和图5中的重建结果，可以看出本发明的投影数据滤波技术可以有效地减少重建图像中的噪声和伪影。

本发明的基于低剂量投影数据滤波的X射线CT图像重建方法,通过对CT投影数据的噪声设计滤波(I)式中的权重因子，将滤波后投影数据与原始采集的投影数据进行加权融合，有效解决了一般性降噪技术进行投影数据滤波所导致的噪声去除与分辨率保持的矛盾性问题，可以在大幅减少图像噪声和伪影的同时较好地保持原有图像的分辨率，最终实现CT图像的优质重建。

需要说明的是，本发明中一般的高斯型滤波器均可选为滤波器F(g)，如均值滤波、中值滤波，非局部平均滤波等。

需要说明的是，权重因子w构造为方差的负相关函数，除了本发明实例中给出的反比例函数外，也可以为或者等其他函数形式，其中β,ε表示正实数。

最后应当说明的是，以上实施例仅用以说明本发明的技术方案而非对本发明保护范围的限制，尽管参照较佳实施例对本发明作了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明技术方案的实质和范围。

Claims (3)

1.一种基于低剂量投影数据滤波的X射线CT图像重建方法，其特征在于：依次包括如下步骤：

(1)获取CT设备的***参数和低剂量扫描协议下的投影数据y_raw；

(2)对步骤(1)中所获取的投影数据y_raw进行逐个数据点上的方差估计；

(3)进行滤波器设计，对步骤(1)中获取的投影数据y_raw进行滤波，获得滤波后的投影数据y_filtered；

(4)对步骤(3)中获得的滤波后的投影数据y_filtered进行CT图像重建，得到最终的CT图像；

所述步骤(2)采用基于小邻域图像的局部方差估计方法或者基于CT投影数据噪声特性的方差估计方法对步骤(1)获取的投影数据y_raw进行逐个数据点上的方差估计；

方差的估计公式为： ${\mathrm{\σ}}_{y_{raw}}^2=\frac{1}{I_0}{e}^{y_{raw}}(1+\frac{{\mathrm{\σ}}_e^2-1.25}{I_0}{e}^{y_{raw}}),$ 其中，I₀为待计算方差的第i个数据点的X射线入射光子强度；为***电子噪声的方差；

所述步骤(3)具体包括：

(3.1)进行滤波器形式设计，为：

$y_{filtered}^{n+1}=w\·y_{raw}+(1-w)\·F\left(y_{filtered}^n\right)---\left(I\right);$

其中y_raw为步骤(1)中采集的投影数据；与分别为第n次和第n+1次滤波后的投影数据，n为自然数，表示迭代运算的次数；w为权重因子；

(3.2)将步骤(1)中的投影数据y_raw作为初始值代入公式(I)中进行迭代运算，得到滤波后的投影数据然后以同样的方法依次进行下一次迭代运算，依次不断地循环执行公式(I)，当循环次数达到预设的次数时即停止迭代运算，并将所得到的迭代运算结果作为最终的滤波后的投影数据y_filtered；

权重因子w为方差的负相关函数；

权重因子w与方差的关系为 $w=1/\left({\mathrm{\β\σ}}_{y_{raw}}^2\right)$ 或者 $w=e^{-{\mathrm{\β\σ}}_{y_{raw}}^2}$ 或者其中β和ε表示正实数；

所述步骤(4)具体是采用滤波反投影法或者迭代重建法进行CT图像重建。

2.根据权利要求1所述的基于低剂量投影数据滤波的X射线CT图像重建方法，其特征在于：迭代运算的次数n范围设置为大于等于5而小于等于30。

3.根据权利要求2所述的基于低剂量投影数据滤波的X射线CT图像重建方法，其特征在于：迭代运算的次数n设置为9或10。