CN104821003A - 一种ct图像重建方法 - Google Patents

Abstract

本发明公开了一种CT图像重建方法，包括如下步骤：步骤一、获得原始CT图像；步骤二、获得当前加权惩罚因子；步骤三、构建加权全变分重建模型；步骤四、求解加权全变分的最小值；步骤五、获得更新CT图像；步骤六、判断获得的更新CT图像是否满足迭代终止条件。本发明提供的CT图像重建方法，通过应用自定义的加权惩罚因子函数构建加权全变分重建模型，利用求解加权全变分的最小值以更新CT图像，并通过不断迭代更新直至输出最终的重建CT图像；该CT图像重建方法解决了滤波反投影算法在数据采样不充分条件下存在的伪影问题，大大提高了重建CT图像的质量。

Description

一种CT图像重建方法

技术领域

本发明涉及CT成像领域，特别涉及一种CT图像重建方法。

背景技术

目前，CT扫描成像中主要采用基于滤波反投影的解析重建方法来生成被扫描物体的断层图像，首先通过对物体进行CT扫描获得原始投影数据，再采用一维斜坡滤波器对投影数据进行卷积处理，并将卷积处理后的投影数据乘上距离倒数加权的反投影因子用以更新待重建点的CT图像数据值。上述重建方法对于数据采样率比较充分的CT扫描数据能取得较好的重建结果，但当CT扫描的采样率比较低，不满足奈奎斯特采样定律的条件时，无法精确重建出原始信号，将导致重建图像存在比较明显的伪影，降低CT图像的质量。

发明内容

针对上述技术中存在的不足之处，本发明提供了一种可以解决滤波反投影算法在数据采样不充分条件下存在的伪影，并能够清晰、准确重建CT图像的CT图像重建方法。

本发明采用的技术方案是：一种CT图像重建方法，包括如下步骤：步骤一、获取CT扫描的原始投影数据，并根据原始投影数据获得原始CT图像；步骤二、根据当前CT图像获得当前加权惩罚因子；步骤三、根据当前加权惩罚因子构建加权全变分重建模型；步骤四、求解加权全变分的最小值；步骤五、根据求解获得的加权全变分的最小值，对当前CT图像进行重建，获得更新CT图像；步骤六、判断获得的更新CT图像是否满足迭代终止条件：若是，则直接输出获得的更新CT图像作为最终的重建CT图像；若否，则重复执行步骤二至步骤五，直至满足迭代终止条件，输出作为最终的重建CT图像；其中，将原始CT图像设置为初次迭代的当前CT图像，以第k-1次迭代重建出的更新CT图像作为当前CT图像进入第k次迭代，k≥2。

优选的，所述步骤一中的原始CT图像是对原始投影数据依次进行空气校正、卷积、反投影处理后获得的。

优选的，所述步骤二中的当前加权惩罚因子通过解析如下公式获得：

      

其中：a和δ是用于调整当前加权惩罚因子的参数，x_i，j表示当前CT图像的数据。

优选的，所述步骤三中的加权全变分重建模型通过解析如下公式获得：

       $\begin{array}{cc}{x}^{*}=\underset{x}{\arg \min }f\left(x\right)& \\ =\underset{x}{\arg \min }\underset{i,j}{\mathrm{\Σ}}\left\{w_{i,j}\right|\left|{\left(\mathrm{Dx}\right)}_{i,j}\right|\left|\right\}& s.t.\left|\right|\mathrm{Ac}-p\left|\right|\≤\mathrm{\ϵ},x_{i,j}\≥0\end{array},$

其中：A为当前CT图像的系统矩阵，p为投影数据，ε为投影误差。

优选的，所述步骤四中求解加权全变分的最小值通过解析如下公式获得：

       $\begin{array}{c}{x}^{*}=\underset{x}{\arg \min }L(x,\mathrm{\λ},\mathrm{\η})\\ =\underset{x}{\arg \min }\{f\left(x\right)+\mathrm{\λ}\left(\right||\mathrm{Ax}-p||-\mathrm{\ϵ})-\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_i\}\end{array},$

其中，λ(||Ax-p||-ε)＝0，η_ix_i＝0是需满足的互松弛条件，λ≥0，η_i≥0是需满足的非负性约束条件，λ为预设经验参数，具体的，求解加权全变分的最小值可转化为求η的最小化问题。

优选的，所述η的最小化通过解析如下公式获得：

       ${\mathrm{\η}}^{*}=\arg \min \left\{\right||\▿f\left(x\right)+\mathrm{\λ}\frac{A^T(\mathrm{Ax}-p)}{\left|\right|\mathrm{Ax}-p\left|\right|}|\left|\right\}s.t.\mathrm{\λ}>0.$

优选的，所述迭代终止条件设置为：迭代次数k达到预先设定的最大次数k_max或两次迭代误差||x^(k+1)-x^(k)||小于给定阈值T。

本发明与现有技术相比，其有益效果是：本发明提供的CT图像重建方法，通过应用自定义的加权惩罚因子函数构建加权全变分重建模型，利用求解加权全变分的最小值以更新CT图像，并通过不断迭代更新直至输出最终的重建CT图像；该CT图像重建方法解决了滤波反投影算法在数据采样不充分条件下存在的伪影问题，大大提高了重建CT图像的质量。

附图说明

图1为本发明所述CT图像重建方法的流程示意图。

具体实施方式

下面结合附图对本发明做进一步的详细说明，以令本领域技术人员参照说明书文字能够据以实施。

如图1所示，本发明提供了一种CT图像重建方法，包括如下步骤：

步骤一、获取CT扫描的原始投影数据，并对原始投影数据依次进行空气校正、卷积、反投影处理，获得原始CT图像；

步骤二、根据当前CT图像获得当前加权惩罚因子，当前加权惩罚因子通过解析如下公式获得：

      

其中：a和δ是用于调整当前加权惩罚因子的参数，x_i，j表示当前CT图像的数据；

步骤三、根据当前加权惩罚因子构建加权全变分重建模型，加权全变分重建模型通过解析如下公式获得：

       $\begin{array}{cc}{x}^{*}=\underset{x}{\arg \min }f\left(x\right)& \\ =\underset{x}{\arg \min }\underset{i,j}{\mathrm{\Σ}}\left\{w_{i,j}\right|\left|{\left(\mathrm{Dx}\right)}_{i,j}\right|\left|\right\}& s.t.\left|\right|\mathrm{Ac}-p\left|\right|\≤\mathrm{\ϵ},x_{i,j}\≥0\end{array},$

其中：A为当前CT图像的系统矩阵，p为投影数据，ε为投影误差；

步骤四、求解加权全变分的最小值，求解加权全变分的最小值通过解析如下公式获得：

       $\begin{array}{c}{x}^{*}=\underset{x}{\arg \min }L(x,\mathrm{\λ},\mathrm{\η})\\ =\underset{x}{\arg \min }\{f\left(x\right)+\mathrm{\λ}\left(\right||\mathrm{Ax}-p||-\mathrm{\ϵ})-\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_i\}\end{array},$

其中，λ(||Ax-p||-ε)＝0，η_ix_i＝0是需满足的互松弛条件，λ≥0，η_i≥0是需满足的非负性约束条件，λ为预设经验参数，具体的，求解加权全变分的最小值可转化为求η的最小化问题，所述η的最小化通过解析如下公式获得：

       ${\mathrm{\η}}^{*}=\arg \min \left\{\right||\▿f\left(x\right)+\mathrm{\λ}\frac{A^T(\mathrm{Ax}-p)}{\left|\right|\mathrm{Ax}-p\left|\right|}|\left|\right\}s.t.\mathrm{\λ}>0;$

步骤五、根据求解获得的加权全变分的最小值，对当前CT图像进行重建，获得更新CT图像；

步骤六、判断获得的更新CT图像是否满足迭代终止条件，所述迭代终止条件设置为：迭代次数k达到预先设定的最大次数k_max或两次迭代误差||x^(k+1)-x^(k)||小于给定阈值T：若是，则直接输出获得的更新CT图像作为最终的重建CT图像；若否，则重复执行步骤二至步骤五，直至满足迭代终止条件，输出作为最终的重建CT图像；

其中，将原始CT图像设置为初次迭代的当前CT图像，以第k-1次迭代重建出的更新CT图像作为当前CT图像进入第k次迭代，k≥2。

实施例1：

步骤一、读取待处理的CT扫描原始投影数据，并对原始投影数据依次进行空气校正、卷积、反投影处理，获得原始CT图像；初始化迭代参数：k_max＝5，λ＝2；

步骤二、根据当前CT图像获得当前加权惩罚因子，当前加权惩罚因子通过解析如下公式获得：

      

其中：1≤i，j≤N，其中N是当前CT图像的大小，a和δ是用于调整当前加权惩罚因子的参数，x_i，j表示当前CT图像的数据，a越大，δ越小，加权惩罚因子的效果越理想，选择过大的a或过小的δ易导致重建算法收敛于局部极值，使得重建CT图像中出现局部亮点或者暗点，为了更好的保持重建CT图像中的边缘信息，可以选择小于边缘梯度的δ和适中的a，在本实施例中，δ＝10^-3，a＝5；

步骤三、根据当前加权惩罚因子构建加权全变分重建模型，加权全变分重建模型通过解析如下公式获得：

       $\begin{array}{cc}{x}^{*}=\underset{x}{\arg \min }f\left(x\right)& \\ =\underset{x}{\arg \min }\underset{i,j}{\mathrm{\Σ}}\left\{w_{i,j}\right|\left|{\left(\mathrm{Dx}\right)}_{i,j}\right|\left|\right\}& s.t.\left|\right|\mathrm{Ac}-p\left|\right|\≤\mathrm{\ϵ},x_{i,j}\≥0\end{array},$

其中：A为当前CT图像的系统矩阵，其实现可通过距离驱动的三维投影与反投影方法，p为投影数据，ε为投影误差，与光子散射、电子噪声、受检物体等多种因素相关；

步骤四、求解加权全变分的最小值，加权全变分重建模型根据最优化理论中的(KKT)条件可以转为其拉格朗日问题求解：

       $\begin{array}{c}{x}^{*}=\underset{x}{\arg \min }L(x,\mathrm{\λ},\mathrm{\η})\\ =\underset{x}{\arg \min }\{f\left(x\right)+\mathrm{\λ}\left(\right||\mathrm{Ax}-p||-\mathrm{\ϵ})-\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_i\}\end{array},$

其中，λ(||Ax-p||-ε)＝0，η_ix_i＝0是需满足的互松弛条件，λ≥0，η_i≥0是需满足的非负性约束条件，λ为预设经验参数，可在迭代过程中更新；

在本实施例中，λ^(k+1)＝0.95*λ^(k)，如果要满足投影项约束条件，则λ＞0，x_i为待重建的数据，若要求得x_i的非零解，则其对应的η_i→0，因此上述求解过程可进一步转化为求η的最小化问题：

       ${\mathrm{\η}}^{*}=\arg \min \left\{\right||\▿f\left(x\right)+\mathrm{\λ}\frac{A^T(\mathrm{Ax}-p)}{\left|\right|\mathrm{Ax}-p\left|\right|}|\left|\right\}s.t.\mathrm{\λ}>0,$

对上述η的最小化问题求解，可以通过交替使用凸集投影法和自适应梯度下降过程求解；

步骤五、根据求解获得的加权全变分的最小值，对当前CT图像进行重建，获得更新CT图像；

步骤六、判断获得的更新CT图像是否满足迭代终止条件，所述迭代终止条件设置为：迭代次数k达到预先设定的最大次数k_max或两次迭代误差||x^(k+1)-x^(k)||小于给定阈值T：若是，则直接输出获得的更新CT图像作为最终的重建CT图像；若否，则重复执行步骤二至步骤五，直至满足迭代终止条件，输出作为最终的重建CT图像；

其中，将原始CT图像设置为初次迭代的当前CT图像，以第k-1次迭代重建出的更新CT图像作为当前CT图像进入第k次迭代，k≥2。

本发明提供的CT图像重建方法，通过应用自定义的加权惩罚因子函数构建加权全变分重建模型，利用求解加权全变分的最小值以更新CT图像，并通过不断迭代更新直至输出最终的重建CT图像；该CT图像重建方法解决了滤波反投影算法在数据采样不充分条件下存在的伪影问题，大大提高了重建CT图像的质量。

尽管本发明的实施方案已公开如上，但其并不仅限于说明书和实施方式中所列运用，它完全可以被适用于各种适合本发明的领域，对于熟悉本领域的人员而言，可容易地实现另外的修改，因此在不背离权利要求及等同范围所限定的一般概念下，本发明并不限于特定的细节和这里示出与描述的图例。

Claims (7)

1.一种CT图像重建方法，其特征在于，包括如下步骤：

步骤一、获取CT扫描的原始投影数据，并根据原始投影数据获得原始CT图像；

步骤二、根据当前CT图像获得当前加权惩罚因子；

步骤三、根据当前加权惩罚因子构建加权全变分重建模型；

步骤四、求解加权全变分的最小值；

步骤五、根据求解获得的加权全变分的最小值，对当前CT图像进行重建，获得更新CT图像；

步骤六、判断获得的更新CT图像是否满足迭代终止条件：若是，则直接输出获得的更新CT图像作为最终的重建CT图像；若否，则重复执行步骤二至步骤五，直至满足迭代终止条件，输出作为最终的重建CT图像；

其中，将原始CT图像设置为初次迭代的当前CT图像，以第k-1次迭代重建出的更新CT图像作为当前CT图像进入第k次迭代，k≥2。

2.如权利要求1所述的CT图像重建方法，其特征在于，所述步骤一中的原始CT图像是对原始投影数据依次进行空气校正、卷积、反投影处理后获得的。

3.如权利要求1所述的CT图像重建方法，其特征在于，所述步骤二中的当前加权惩罚因子通过解析如下公式获得：

其中：a和δ是用于调整当前加权惩罚因子的参数，x_i，j表示当前CT图像的数据。

4.如权利要求3所述的CT图像重建方法，其特征在于，所述步骤三中的加权全变分重建模型通过解析如下公式获得：

$\begin{array}{ccc}{x}^{*}=\underset{x}{\arg \min }f\left(x\right)& & \\ =\underset{x}{\arg \min }\underset{i,j}{\mathrm{\Σ}}\left\{w_{i,j}\right|\left|{\left(\mathrm{Dx}\right)}_{i,j}\right|\left|\right\}& s.t.\left|\right|\mathrm{Ax}-p\left|\right|\≤\mathrm{\ϵ},& x_{i,j}\≥0\end{array},$

其中：A为当前CT图像的系统矩阵，p为投影数据，ε为投影误差。

5.如权利要求4所述的CT图像重建方法，其特征在于，所述步骤四中求解加权全变分的最小值通过解析如下公式获得：

$\begin{array}{c}{x}^{*}=\underset{x}{\arg \min }L(x,\mathrm{\λ},\mathrm{\η})\\ =\underset{x}{\arg \min }\{f\left(x\right)+\mathrm{\λ}\left(\right||\mathrm{Ax}-p||-\mathrm{\ϵ})-\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_i\},\end{array}$

其中，λ(||Ax-p||-ε)＝0，η_ix_i＝0是需满足的互松弛条件，λ≥0，η_i≥0是需满足的非负性约束条件，λ为预设经验参数，具体的，求解加权全变分的最小值可转化为求η的最小化问题。

6.如权利要求5所述的CT图像重建方法，其特征在于，所述η的最小化通过解析如下公式获得：

$\begin{array}{cc}{\mathrm{\η}}^{*}=\arg \min \left\{\right||\▿f\left(x\right)+\mathrm{\λ}\frac{A^T(\mathrm{Ax}-p)}{\left|\right|\mathrm{Ax}-p\left|\right|}|\left|\right\}& s.t.\mathrm{\λ}>0\end{array}.$

7.如权利要求6所述的CT图像重建方法，其特征在于，所述迭代终止条件设置为：迭代次数k达到预先设定的最大次数k_max或两次迭代误差||x^(k+1)-x^(k)||小于给定阈值T。