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

Abstract

本发明适用于医学和工业领域计算机断层扫描系统领域，提供了一种CT图像的重建方法，其特征在于，所述方法包括以下步骤；步骤S101，获取CT扫描所采集的投影数据；步骤S102，应用成像模型方程式G＝MX，采用

公式计算出重建图像；其中，G为CT扫描所采集的投影数据，M为系统矩阵，X为重建的CT图像，Xⁱ表示第i次迭代后得到的重建图像，λ表示收敛系数。本发明在CT扫描时，进行高欠采样数据采集可以有效缩短CT扫描时间、所采集的数据量较小，从而提高了图像重建速度、由于扫描时间短，从而减少病人运动所引起的伪影，能在保证CT图像质量的前提下大幅降低病人所承受的CT辐射剂量。

Description

一种CT图像的重建方法

技术领域

本发明属于医学和工业领域计算机断层扫描系统领域，尤其涉及一种CT图像的重建方法。

背景技术

计算机断层成像（CT）是通过无损方式获取物体内部结构信息的一种重要成像手段，它拥有高分辨率、高灵敏度以及多层次等众多优点，是我国装机量最大的医疗影像诊断设备之一，被广泛应用于各个医疗临床检查领域。然而，由于CT扫描过程中需要使用X射线，随着人们对辐射潜在危害的逐步了解，CT辐射剂量问题越来越受到人们的重视。合理使用低剂量（As Low AsReasonably Achievable，ALARA）原则要求在满足临床诊断的前提下，尽量降低对患者的辐射剂量。因此，研究和开发新的低剂量CT成像方法，既能保证CT成像质量又减少有害的辐射剂量，对于医疗诊断领域具有重要的科学意义和应用前景。

现有技术提出了一种稀疏角度CT图像的重建方法。该方法是：获取CT机的系统参数和一周期内的稀疏角度投影数据，并将所获得的稀疏角度投影数据按采集的时间顺序等分成若干组；应用重建模型采用辅助函数法求解得到要求的迭代运算方法进行循环迭代运算，当循环次数达到预设的次数时即停止迭代运算，并将当前所得到的迭代运算结果作为重建图像。

然而，上述现有技术方案由于要进行全采样，导致CT扫描时间长；由于全采样所采集的数据量较大，导致图像重建速度慢；由于扫描时间长，导致出现病人运动所引起的伪影；病人所承受的CT辐射剂量高。

发明内容

本发明的目的在于提供一种CT图像的重建方法，旨在解决现有技术的稀疏角度CT图像的重建方法由于要进行全采样，导致CT扫描时间长；由于全采样所采集的数据量较大，导致图像重建速度慢；由于扫描时间长，导致出现病人运动所引起的伪影；病人所承受的CT辐射剂量高的问题。

本发明是这样实现的，一种CT图像的重建方法，所述方法包括以下步骤；

步骤S101，获取CT扫描所采集的投影数据；

步骤S102，应用成像模型方程式G＝MX，采用 $X^{i+1}=X^i+(G_i-M_i\·X^i)\·M_i^T\·\mathrm{\λ}/(M_i\·M_i)$ 公式计算出重建图像；其中，G为CT扫描所采集的投影数据，M为系统矩阵，X为重建的CT图像，Xⁱ表示第i次迭代后得到的重建图像，λ表示收敛系数。

在本发明中，本发明在CT扫描时，进行高欠采样数据采集可以有效缩短CT扫描时间、所采集的数据量较小，从而提高了图像重建速度、由于扫描时间短，从而减少病人运动所引起的伪影，能在保证CT图像质量的前提下大幅降低病人所承受的CT辐射剂量。

附图说明

图1是本发明实施例提供的CT图像的重建方法的实现流程示意图。

具体实施方式

为了使本发明的目的、技术方案及有益效果更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

请参阅图1，为本发明实施例提供的CT图像的重建方法的实现流程，其包括以下步骤；

在步骤S101中，获取CT扫描所采集的投影数据；

在步骤S102中，应用成像模型方程式G＝MX，采用

公式计算出重建图像；其中，G为CT扫描所采集的投影数据，M为系统矩阵，X为重建的CT图像，Xⁱ表示第i次迭代后得到的重建图像，λ表示收敛系数。

具体的CT图像的重建方法如下：

1、定义成像模型方程如公式（1），其中G为CT扫描所采集的投影数据，M为系统矩阵，X为重建的CT图像。定义参数a＝0.2；b＝0.9，b一般取0-1之间的数值。特殊的，b也可直接取值为0，以用来简化程序；c＝0.96，λ＝1，d＝0.993。

G＝MX    （1）

2、将重建图像X的初始值赋为0，即X＝0；

3、进入主循环程序，将程序执行i=1,L,istepnum次，一般取istepnum为20到100之间的整数；

4、将X的值赋给X⁰，即X⁰＝X；

5、进入子循环程序I，按公式（2）计算重建图像X。其中Xⁱ表示第i次迭代后得到的重建图像；λ表示收敛系数，一般取小于1的正数，按实际经验一般在0.9-1.0范围内取值。

$X^{i+1}=X^i+(G_i-M_i\·X^i)\·M_i^T\·\mathrm{\λ}/(M_i\·M_i)---\left(2\right)$

子循环程序I结束；

6、进入子循环程序II，对重建图像X内的所有点依次进行如下操作

子循环程序II结束；

7、由公式（4），计算此时X与步骤4中X⁰的差异参数η

$\mathrm{\η}=\sqrt{\underset{s,t,v}{\mathrm{\Σ}}{(X_{s,t,v}-X_{s,t,v}^0)}^2}---\left(4\right)$

8、进行如下判断操作，即如果主循环程序执行的是第一次循环，则给ω赋初值ω＝a×η；

9、将X的值赋给X⁰，即X⁰＝X；

10、进入子循环III，

$\stackrel{i}{X}=(1-b)\·{\▿}_{X_{s,t,v}}\left(\right|\left|X_{s,t,v}\right|\left|\right)+b\·{\▿}_{X_{s,t,v}}\left(\right||X_{s,t,v}-X_{s,t,v}^p|\left|\right)---\left(6\right)$

$\hat{X}=\frac{\stackrel{i}{X}}{\sqrt{\underset{s,t,v}{\mathrm{\Σ}}{\left({\stackrel{i}{X}}_{s,t,v}\right)}^2}}---\left(7\right)$

$X=X-\mathrm{\ω}\·\hat{X}---\left(8\right)$

子循环程序III结束；

子循环III一般执行5-40次。其中：

为先验图像，一般为全扫描重建得到的CT图像。而||Xs，t，v，||可由公式（9）得到，

可由公式（10）得到，也可参考公式（10）用类似方法得到，

$\left|\right|X_{s,t,v}\left|\right|=\underset{s,t,v}{\mathrm{\Σ}}\sqrt{{(X_{s,t,v}-X_{s-1,t,v})}^2+{(X_{s,t,v}-X_{s,t-1,v})}^2+{(X_{s,t,v}-X_{s,t,v-1})}^2}---\left(9\right)$

${\▿}_{X_{s,t,v}}\left(\right|\left|X_{s,t,v}\right|\left|\right)=\frac{(X_{s,t,v}-X_{s-1,t,v})+(X_{s,t,v}-X_{s,t-1,v})+(X_{s,t,v}-X_{s,t,v-1})}{\sqrt{{(X_{s,t,v}-X_{s-1,t,v})}^2+{(X_{s,t,v}-X_{s,t-1,v})}^2+{(X_{s,t,v}-X_{s,t,v-1})}^2}}$

$-\frac{(X_{s+1,t,v}-X_{s,t,v})}{\sqrt{{(X_{s+1,t,v}-X_{s,t,v})}^2+{(X_{s+1,t,v}-X_{s+1,t-1,v})}^2+{(X_{s+1,t,v}-X_{s+1,t,v-1})}^2}}$

$-\frac{(X_{s,t+1,v}-X_{s,t,v})}{\sqrt{{(X_{s,t+1,v}-X_{s-1,t+1,v})}^2+{(X_{s,t+1,v}-X_{s,t,v})}^2+{(X_{s,t+1,v}-X_{s,t+1,v-1})}^2}}$

$-\frac{(X_{s,t,v+1}-X_{s,t,v})}{\sqrt{{(X_{s,t,v+1}-X_{s-1,t,v+1})}^2+{(X_{s,t,v+1}-X_{s,t-1,v+1})}^2+{(X_{s,t,v+1}-X_{s,t,v})}^2}}---\left(10\right)$

11、再次给ω和λ赋值

ω＝ω×c    （11）

λ＝λ×d    （12）

12、返回主循环程序。

综上所述，本发明实施例在CT扫描时，进行高欠采样数据采集可以有效缩短CT扫描时间、所采集的数据量较小，从而提高了图像重建速度、由于扫描时间短，从而减少病人运动所引起的伪影，能在保证CT图像质量的前提下大幅降低病人所承受的CT辐射剂量。

本领域普通技术人员可以理解实现上述实施例方法中的全部或部分步骤是可以通过程序来指令相关的硬件来完成，所述的程序可以存储于一计算机可读取存储介质中，所述的存储介质，如ROM/RAM、磁盘、光盘等。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种CT图像的重建方法，其特征在于，所述方法包括以下步骤；

步骤S101，获取CT扫描所采集的投影数据；

步骤S102，应用成像模型方程式G＝MX，采用

公式计算出重建图像；其中，G为CT扫描所采集的投影数据，M为系统矩阵，X为重建的CT图像，Xⁱ表示第i次迭代后得到的重建图像，λ表示收敛系数。

2.如权利要求1所述的方法，其特征在于，步骤S102具体包括以下步骤：步骤1、定义成像模型方程如公式（1），其中G为CT扫描所采集的投影数据，M为系统矩阵，X为重建的CT图像；

G＝MX    （1）

步骤2，将重建图像X的初始值赋为0，即X＝0；

步骤3、进入主循环程序，将程序执行i=1,L,istepnum次；

步骤4、将X的值赋给X⁰，即X⁰＝X；

步骤5、进入子循环程序I，按公式（2）计算重建图像X；其中Xi表示第i次迭代后得到的重建图像；λ表示收敛系数；

$X^{i+1}=X^i+(G_i-M_i\·X^i)\·M_i^T\·\mathrm{\λ}/(M_i\·M_i)---\left(2\right)$

子循环程序I结束；

步骤6、进入子循环程序II，对重建图像X内的所有点依次进行如下操作

子循环程序II结束；

步骤7、由公式（4），计算此时X与步骤4中X⁰的差异参数η；

$\mathrm{\η}=\sqrt{\underset{s,t,v}{\mathrm{\Σ}}{(X_{s,t,v}-X_{s,t,v}^0)}^2}---\left(4\right)$

步骤8、进行如下判断操作，即如果主循环程序执行的是第一次循环，则给ω赋初值ω＝a×η；

步骤9、将X的值赋给X⁰，即X⁰＝X；

步骤10、进入子循环III，

$\stackrel{i}{X}=(1-b)\·{\▿}_{X_{s,t,v}}\left(\right|\left|X_{s,t,v}\right|\left|\right)+b\·{\▿}_{X_{s,t,v}}\left(\right||X_{s,t,v}-X_{s,t,v}^p|\left|\right)---\left(6\right)$

$\hat{X}=\frac{\stackrel{i}{X}}{\sqrt{\underset{s,t,v}{\mathrm{\Σ}}{\left({\stackrel{i}{X}}_{s,t,v}\right)}^2}}---\left(7\right)$

$X=X-\mathrm{\ω}\·\hat{X}---\left(8\right)$

子循环程序III结束；

其中：

为先验图像，一般为全扫描重建得到的CT图像；而||Xs，t，v||可由公式（9）得到，

可由公式（10）得到，可由公式（10）得到，

$\left|\right|X_{s,t,v}\left|\right|=\underset{s,t,v}{\mathrm{\Σ}}\sqrt{{(X_{s,t,v}-X_{s-1,t,v})}^2+{(X_{s,t,v}-X_{s,t-1,v})}^2+{(X_{s,t,v}-X_{s,t,v-1})}^2}---\left(9\right)$

${\▿}_{X_{s,t,v}}\left(\right|\left|X_{s,t,v}\right|\left|\right)=\frac{(X_{s,t,v}-X_{s-1,t,v})+(X_{s,t,v}-X_{s,t-1,v})+(X_{s,t,v}-X_{s,t,v-1})}{\sqrt{{(X_{s,t,v}-X_{s-1,t,v})}^2+{(X_{s,t,v}-X_{s,t-1,v})}^2+{(X_{s,t,v}-X_{s,t,v-1})}^2}}$

$-\frac{(X_{s+1,t,v}-X_{s,t,v})}{\sqrt{{(X_{s+1,t,v}-X_{s,t,v})}^2+{(X_{s+1,t,v}-X_{s+1,t-1,v})}^2+{(X_{s+1,t,v}-X_{s+1,t,v-1})}^2}}$

$-\frac{(X_{s,t+1,v}-X_{s,t,v})}{\sqrt{{(X_{s,t+1,v}-X_{s-1,t+1,v})}^2+{(X_{s,t+1,v}-X_{s,t,v})}^2+{(X_{s,t+1,v}-X_{s,t+1,v-1})}^2}}$

$-\frac{(X_{s,t,v+1}-X_{s,t,v})}{\sqrt{{(X_{s,t,v+1}-X_{s-1,t,v+1})}^2+{(X_{s,t,v+1}-X_{s,t-1,v+1})}^2+{(X_{s,t,v+1}-X_{s,t,v})}^2}}---\left(10\right)$

步骤11、再次给ω和λ赋值

ω＝ω×c    （11）

λ＝λ×d    （12）

步骤12、返回主循环程序。

3.如权利要求2所述的方法，其特征在于，定义参数a＝0.2；b取0-1之间的数值；c＝0.96，λ＝1，d＝0.993。

4.如权利要求3所述的方法，其特征在于，所述b＝0.9。

5.如权利要求2所述的方法，其特征在于，所述取istepnum为20到100之间的整数。

6.如权利要求2所述的方法，其特征在于，所述λ取小于1的正数。

7.如权利要求6所述的方法，其特征在于，所述λ在0.9-1.0范围内取值。

Images