CN103810733B - 一种稀疏角度x射线ct图像的统计迭代重建方法 - Google Patents
一种稀疏角度x射线ct图像的统计迭代重建方法 Download PDFInfo
- Publication number
- CN103810733B CN103810733B CN201410071978.7A CN201410071978A CN103810733B CN 103810733 B CN103810733 B CN 103810733B CN 201410071978 A CN201410071978 A CN 201410071978A CN 103810733 B CN103810733 B CN 103810733B
- Authority
- CN
- China
- Prior art keywords
- data
- projection
- image
- ray
- raw
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Landscapes
- Apparatus For Radiation Diagnosis (AREA)
Abstract
本发明公开一种稀疏角度X射线CT图像的统计迭代重建方法，其包括步骤：获取CT设备的系统参数及稀疏角度低剂量X射线CT扫描的投影数据yraw，对投影数据yraw进行基于中值先验约束的数据恢复处理，得到恢复后的投影数据yrestored；对投影数据yraw及恢复后的投影数据yrestored进行自适应加权处理，得到加权后的投影数据yweight；对加权后的投影数据yweight采用统计迭代重建方法进行图像重建，得到重建后的低剂量CT图像。本发明通过基于中值先验约束的投影数据恢复以及自适应加权处理，对采集到的稀疏低剂量CT投影数据有效的降噪处理，最终实现稀疏低剂量CT图像的优质重建，重建后的CT图像能够有效地去除重建图像中噪声引起的斑块效应，CT图像质量有着明显的改善。
Description
技术领域
本发明涉及医学影像处理技术，尤其是涉及一种稀疏角度X射线CT图像的统计迭代重建方法。
背景技术
X射线CT扫描已经广泛应用于临床医学影像诊断，但是CT扫描过程中过高的X射线辐射剂量会存在致癌风险。为了降低对使用者的损害，如何最大限度地降低X射线使用剂量已经成为医学CT成像领域研究的关键技术之一。
为了降低X射线辐射剂量，现有技术可以通过以下两种方式实现低剂量CT扫描：一是降低CT扫描过程中的管电流和扫描时间，即低mAs（milliampere-seconds）扫描方式，二是减少CT扫描过程中X射线球管旋转一周内的曝光次数，即稀疏角度扫描方式。目前，将上述两种扫描方式结合在一起，即在稀疏角度扫描的同时适量地降低mAs，成为临床应用研究的热点。同时降低CT扫描的mAs以及球管曝光角度，将会导致采集到的投影数据中不仅含有大量的噪声，还会由于稀疏角度扫描使得采集的数据具有不完备性，从而使得常规的图像重建方法，如传统的滤波反投影方法（filtered back-projection,FBP），重建的图像质量出现严重的退化现象，难以满足临床诊断需要。
为了在保证图像质量的前提下大幅降低X射线辐射剂量，诸多基于统计模型的低剂量CT图像迭代重建方法相继提出。通过对采集的投影数据的噪声以及成像系统进行图像重建模型构建，基于统计模型的迭代重建方法可以实现低剂量CT图像优质重建。同时考虑到稀疏采集数据的不完备性，为了保证统计迭代重建方法的鲁棒性，通常会引入一定的先验信息进行约束。目前研究中，基于重建图像分段光滑假设的全变分（total variation,TV）约束被广泛引入到统计迭代重建方法。其中代表性方法为Sidky等人提出的基于TV的凸集投影（projection onto convex sets，POCS）重建算法，简称TV-POCS方法，对于稀疏角度采集的数据具有较好地图像重建效果，可以有效地去除图像中的噪声和条形伪影。
而现有技术中的针对稀疏角度扫描的统计迭代重建方法，均是默认高mAs的情况，没有考虑稀疏角度扫描与低mAs的结合情况。在同时存在稀疏角度扫描和低mAs的情况，TV-POCS方法将会受到低mAs情况下投影数据大量噪声的影响，使得重建图像产生噪声引起的斑块效应。
发明内容
针对现有技术不足，本发明提供一种稀疏角度X射线低剂量CT图像的统计迭代重建方法，能够针对同时降低管电流和扫描时间以及减少球管曝光角度情况下的稀疏角度低剂量射线CT图像进行图像重建，重建后的图像能够有效地去除重建图像中噪声引起的斑块效应。
本发明采用如下技术方案实现：一种稀疏角度X射线CT图像的统计迭代重建方法，其包括步骤：
获取CT设备的系统参数及稀疏角度低剂量X射线CT扫描的投影数据yraw，对投影数据yraw进行基于中值先验约束的数据恢复处理，得到恢复后的投影数据yrestored；
对投影数据yraw及恢复后的投影数据yrestored进行自适应加权处理，得到加权后的投影数据yweight；
对加权后的投影数据yweight采用统计迭代重建方法进行图像重建，得到重建后的低剂量CT图像。
其中，对投影数据yraw进行基于中值先验约束的数据恢复处理的步骤包括：
设置投影数据yraw进行恢复处理的数据恢复模型为y表示待恢复的投影数据，Σ为一对角矩阵，其对角线元素为对应投影数据yraw的方差R(y)为先验约束项，λ>0为正则化调节参数；
通过梯度下降法或者高斯-塞德尔优化算法对数据恢复模型进行迭代求解，得到恢复后的投影数据yrestored。
其中，先验约束项R(y)的设计采用了基于中值滤波推衍出的中值先验形式，即R(y)=(y-median(y))2/2，其中，median(y)表示投影数据y的中值滤波结果值。
其中，其中I0为投影数据yraw中第i个数据点的X射线入射光子强度，为CT设备的系统电子噪声的方差。
其中，自适应加权处理的计算公式为yweight=w·yraw+(1-w)·yrestored，其中权重因子w是根据投影数据yraw的方差计算得到。
其中，或者其中δ、β和ε均为预设的非负实数。
其中，阈值参数δ取值为0.01。
其中，对加权后的投影数据yweight采用统计迭代重建方法进行图像重建的步骤具体包括：
设置统计迭代重建方法的图像重建模型为
其中，μ为待重建的CT图像，H为CT设备中成像系统矩阵，W为一对角矩阵，T表示矩阵的转置运算，R(μ)表示统计迭代重建中的先验约束项，γ为正则化调节参数；
对上述图像重建模型求解，得到重建后的低剂量CT图像。
其中，对上述图像重建模型求解方法为TV-POCS求解方法、梯度下降算法或者共轭梯度优化算法。
与现有技术相比，本发明具有如下有益效果：
本发明技术方案能够针对同时降低管电流和扫描时间以及减少球管曝光角度情况下的稀疏角度低剂量CT投影数据，通过基于中值先验约束的投影数据恢复以及自适应加权处理，对采集到的稀疏低剂量CT投影数据有效的降噪处理，同时较好地保持了原有数据的细节信息，有效解决了现有技术（如TV-POCS技术）直接从采集到的未经处理的稀疏低剂量CT投影数据进行图像重建所产生的图像块斑效应，最终实现稀疏低剂量CT图像的优质重建，重建后的CT图像能够有效地去除重建图像中噪声引起的斑块效应，CT图像质量有着明显的改善。
附图说明
图1为本发明方法的流程示意图。
图2为本发明实施例中仿真所使用的XCAT体模数据示意图。
图3a、图3b和图3c为本发明实施例中由模拟生成的稀疏低剂量投影数据（58个投影角度）分别通过FBP重建方法、TV-POCS重建方法本发明方法的重建结果示意图。
图4a、图4b和图4c为本发明实施例中由模拟生成的稀疏低剂量投影数据（116个投影角度）分别通过FBP重建方法、TV-POCS重建方法本发明方法的重建结果示意图。
具体实施方式
本发明涉及一种医学影像的图像重建方法，具体来说涉及一种基于稀疏角度低剂量X射线CT图像的统计迭代重建方法。
如图1所示，本发明一个优选实施例包括如下实现步骤：
步骤S1、获取CT设备的系统参数和低剂量X射线扫描协议下（低mAs和稀疏角度扫描）的投影数据yraw。
CT设备的系统参数主要包括X射线入射光子强度I0，系统电子噪声的方差σe2等。
步骤S2、对步骤S1中获取的投影数据yraw进行数据恢复处理，得到恢复后的投影数据yrestored。
对投影数据yraw进行恢复处理的数据恢复模型为
其中，y表示待恢复的投影数据，Σ为一对角矩阵，其对角线元素为对应投影数据yraw的方差R(y)为先验约束项，λ>0为正则化调节参数，用于刻画先验约束的强度。
考虑到中值滤波在对数据噪声进行平滑处理的同时，对图像数据中的某些细节也能起到很好的保护作用，先验约束项R(y)的设计采用了基于中值滤波推衍出的中值先验形式，即R(y)=(y-median(y))2/2，其中，median(y)表示投影数据y的中值滤波结果值。
具体的，投影数据yraw的方差通过对投影数据yraw的逐个数据点估计得到，即通过公式估计得到，其中I0为第i个数据点的X射线入射光子强度，为系统电子噪声的方差；也可以通过局部邻域方差估计等其他方式得到。
具体的，式（I）中数据恢复模型的求解，可以通过梯度下降法或者高斯-塞德尔优化算法进行迭代求解，得到恢复后的投影数据yrestored。
步骤S3、对恢复后的投影数据yrestored和步骤S1中获取的投影数据yraw进行自适应加权处理，得到加权后的投影数据yweight。
自适应加权处理的计算公式为yweight=w·yraw+(1-w)·yrestored，其中权重因子w是根据投影数据yraw的方差计算得到。
具体的，权重因子w的计算公式可以为二值形式，如也可以为非二值形式，如或者其中δ、β和ε均为非负实数。比如，在一个优选实施例中，阈值参数δ取值为0.01。
步骤S4、对步骤S3中得到的加权后的投影数据yweight采用统计迭代重建方法进行图像重建，得到重建后的低剂量CT图像。
统计迭代重建方法的图像重建模型为：
其中，μ为待重建的CT图像，H为CT设备中成像系统矩阵，可由步骤S1获取的系统参数确定，W为一对角矩阵，其对角线元素为权重因子，不考虑权重的情况下，为单位对角矩阵，T表示矩阵的转置运算。R(μ)表示统计迭代重建中的先验约束项，可以为全变分约束，即R(μ)=TV(μ)，也可以为其他形式的先验约束，如Huber先验、非局部均值先验等。γ为正则化调节参数，为一正实数，用于控制先验约束的强度。
具体的，式（II）中图像重建模型的求解，可以采用现有文献（E.Sidky,C.Kao,andX.Pan,“Accurate image reconstruction from few-views and limited-angle data indivergent beam CT,”J.X-Ray Sci.Technol.,vol.14,pp.119-39,2006）所给的TV-POCS求解方法，也可以采用梯度下降算法或者共轭梯度优化算法进行求解，得到重建后最终的CT图像。
采用图2所示的XCAT数字体模图像作为本发明的计算机仿真实验对象。体模图像像素矩阵大小设为512×512，模拟CT机的X射线源到旋转中心和探测器的距离分别为570mm和1040mm，旋转角在[0,2π]间，采样值分别为58和116，每个采样角对应672个探测器单元，探测器单元的大小为1.407mm。通过CT设备仿真分别生成58个投影角度和116个投影角度的两组稀疏角度低剂量CT投影数据yraw，其中X射线的入射光子强度I0为3.0×104，系统电子噪声的方差为11.0。在实际的CT数据采集中，投影数据和系统参数即入射光子强度I0和系统电子噪声的方差均可以直接获取。
对步骤S1中模拟生成的低剂量CT投影数据yraw进行数据恢复处理。具体的，投影数据yraw的方差通过公式估计得到，数据恢复模型的求解采用高斯-塞德尔优化算法求解，得到恢复后的投影数据yrestored。
对步骤S2中得到的恢复后投影数据yrestored和步骤S1中获取的投影数据yraw进行自适应加权处理，得到加权后的投影数据yweight。具体的，权重因子计算采用二值形式，即其中，阈值参数δ取值为0.01。
对步骤S3中得到的加权后的投影数据yweight采用统计迭代重建方法（比如TV-POCS方法）进行图像重建，得到重建后最终的低剂量CT图像。
为了对比本发明所示方法的效果，对模拟生成的两组CT投影数据：图3a、图3b和图3c为本发明实施例中由模拟生成的稀疏低剂量投影数据（58个投影角度）分别通过FBP重建方法、TV-POCS重建方法本发明方法的重建结果示意图；图4a、图4b和图4c为本发明实施例中由模拟生成的稀疏低剂量投影数据（116个投影角度）分别通过FBP重建方法、TV-POCS重建方法本发明方法的重建结果示意图。
将三种重建方法对应的重建结果进行比较可见，相比FBP重建方法，统计迭代重建方法可以有效地减少重建图像中的噪声和伪影。而本发明方法相比TV-POCS重建方法对于噪声引起的斑块效应具有更好的去除效果，尤其在较少的投影角度情况下，效果更加明显，如图3c所示。
同时，为了说明本发明方法相比TC-POCS重建方法对重建图像质量的改善，我们计算了图2中部所示的方框区域内的两种方法重建结果的相对均方根误差。在58个投影角度的情况下，本发明方法和TC-POCS重建方法对应的相对均方根误差值分别为0.015和0.029；在116个投影角度的情况下，本发明方法和TC-POCS重建方法对应的相对均方根误差值分别为0.012和0.014。可以看出本发明方法相比TC-POCS重建方法对重建图像质量有着明显的改善，尤其在较少的投影角度情况下，相对均方根误差可以减少一半。
综上，本发明能够针对同时降低管电流和扫描时间以及减少球管曝光角度情况下的稀疏角度低剂量X射线CT图像的投影数据，通过基于中值先验约束的投影数据恢复以及自适应加权处理，对采集到的稀疏角度低剂量CT投影数据有效的降噪处理，同时较好地保持了原有数据的细节信息，有效解决了现有技术（如TV-POCS技术）直接从采集到的未经处理的稀疏角度低剂量X射线CT投影数据进行图像重建所产生的图像块斑效应，最终实现稀疏角度低剂量X射线CT图像的优质重建。
以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。
Claims (7)
1.一种稀疏角度X射线CT图像的统计迭代重建方法，其特征在于包括步骤：
获取CT设备的系统参数及稀疏角度低剂量X射线CT扫描的投影数据yraw，对投影数据yraw进行基于中值先验约束的数据恢复处理，得到恢复后的投影数据yrestored；
对投影数据yraw及恢复后的投影数据yrestored进行自适应加权处理，得到加权后的投影数据yweight；
对加权后的投影数据yweight采用统计迭代重建方法进行图像重建，得到重建后的低剂量CT图像；
其中，对投影数据yraw进行基于中值先验约束的数据恢复处理的步骤包括：设置投影数据yraw进行恢复处理的数据恢复模型为y表示待求解的恢复后的投影数据，Σ为一对角矩阵，其对角线元素为对应投影数据yraw的方差R(y)为先验约束项，λ＞0为正则化调节参数；通过梯度下降法或者高斯‐塞德尔优化算法对数据恢复模型进行迭代求解，得到恢复后的投影数据yrestored；
其中，自适应加权处理的计算公式为yweight＝w·yraw+(1-w)·yrestored，其中权重因子w是根据投影数据yraw的方差计算得到。
2.根据权利要求1所述一种稀疏角度X射线CT图像的统计迭代重建方法，其特征在于，先验约束项R(y)的设计采用了基于中值滤波推衍出的中值先验形式，即R(y)＝(y-median(y))2/2，其中，median(y)表示投影数据y的中值滤波结果值。
3.根据权利要求1所述一种稀疏角度X射线CT图像的统计迭代重建方法，其特征在于，其中I0为投影数据yraw中第i个数据点的X射线入射光子强度，为CT设备的系统电子噪声的方差。
4.根据权利要求1所述一种稀疏角度X射线CT图像的统计迭代重建方法，其特征在于，或者其中δ、β和ε均为预设的非负实数。
5.根据权利要求4所述一种稀疏角度X射线CT图像的统计迭代重建方法，其特征在于，阈值参数δ取值为0.01。
6.根据权利要求1所述一种稀疏角度X射线CT图像的统计迭代重建方法，其特征在于，对加权后的投影数据yweight采用统计迭代重建方法进行图像重建的步骤具体包括：
设置统计迭代重建方法的图像重建模型为
其中，μ为待重建的CT图像，H为CT设备中成像系统矩阵，W为一对角矩阵，T表示矩阵的转置运算，R(μ)表示统计迭代重建中的先验约束项，γ为正则化调节参数；
对上述图像重建模型求解，得到重建后的低剂量CT图像。
7.根据权利要求6所述一种稀疏角度X射线CT图像的统计迭代重建方法，其特征在于，对上述图像重建模型求解方法为TV‐POCS求解方法、梯度下降算法或者共轭梯度优化算法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410071978.7A CN103810733B (zh) | 2014-02-28 | 2014-02-28 | 一种稀疏角度x射线ct图像的统计迭代重建方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410071978.7A CN103810733B (zh) | 2014-02-28 | 2014-02-28 | 一种稀疏角度x射线ct图像的统计迭代重建方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN103810733A CN103810733A (zh) | 2014-05-21 |
| CN103810733B true CN103810733B (zh) | 2017-04-05 |
Family
ID=50707451
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201410071978.7A Active CN103810733B (zh) | 2014-02-28 | 2014-02-28 | 一种稀疏角度x射线ct图像的统计迭代重建方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN103810733B (zh) |
Families Citing this family (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104103086B (zh) * | 2014-06-06 | 2017-02-15 | 华南理工大学 | 一种稀疏采样角度下基于变分不等式的ct图像重建方法 |
| CN104240210B (zh) * | 2014-07-21 | 2018-08-10 | 南京邮电大学 | 基于压缩感知的ct图像迭代重建方法 |
| WO2016091594A1 (en) * | 2014-12-10 | 2016-06-16 | Koninklijke Philips N.V. | Statistically weighted regularization in multi-contrast imaging |
| CN110136218B (zh) * | 2019-03-28 | 2022-10-28 | 中国人民解放军战略支援部队信息工程大学 | 基于噪声生成机制与数据驱动紧框架的ct投影去噪重建方法及装置 |
| CN110751702B (zh) | 2019-10-29 | 2023-06-16 | 上海联影医疗科技股份有限公司 | 图像重建方法、系统、装置及存储介质 |
| CN111080736B (zh) * | 2019-12-11 | 2023-09-08 | 电子科技大学 | 一种基于稀疏变换的低剂量ct图像重建方法 |
| CN111652951B (zh) * | 2020-05-07 | 2023-06-06 | 中国工程物理研究院材料研究所 | 一种稀疏角度快中子ct成像方法 |
| CN116664707B (zh) * | 2022-02-17 | 2024-02-13 | 北京朗视仪器股份有限公司 | 一种ct图像迭代重建方法、装置及计算机设备 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103136773A (zh) * | 2013-02-05 | 2013-06-05 | 南方医科大学 | 一种稀疏角度x射线ct成像方法 |
| CN103136772A (zh) * | 2012-10-25 | 2013-06-05 | 南方医科大学 | 基于加权阿尔法散度约束的x射线低剂量ct图像重建方法 |
-
2014
- 2014-02-28 CN CN201410071978.7A patent/CN103810733B/zh active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103136772A (zh) * | 2012-10-25 | 2013-06-05 | 南方医科大学 | 基于加权阿尔法散度约束的x射线低剂量ct图像重建方法 |
| CN103136773A (zh) * | 2013-02-05 | 2013-06-05 | 南方医科大学 | 一种稀疏角度x射线ct成像方法 |
Non-Patent Citations (1)
| Title |
|---|
| Iterative Image Reconstruction for Ultra-low-dose CT with a Combined Low-mAs and Sparse-view Protocol;Yunwan Zhang;《35th Annual International Conference of the IEEE EMBS》;20130707;5107-5110 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN103810733A (zh) | 2014-05-21 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN103810733B (zh) | 一种稀疏角度x射线ct图像的统计迭代重建方法 | |
| CN111325686B (zh) | 一种基于深度学习的低剂量pet三维重建方法 | |
| Jia et al. | GPU-based iterative cone-beam CT reconstruction using tight frame regularization | |
| CN103810734B (zh) | 一种低剂量x射线ct投影数据恢复方法 | |
| CN102737392B (zh) | 一种低剂量x线ct图像的非局部正则化先验重建方法 | |
| CN102314698B (zh) | 基于阿尔法散度约束的全变分最小化剂量ct重建方法 | |
| CN101980302A (zh) | 投影数据恢复导引的非局部平均低剂量ct重建方法 | |
| CN103810735A (zh) | 一种低剂量x射线ct图像统计迭代重建方法 | |
| CN112102213B (zh) | 低剂量ct图像处理方法、扫描系统及计算机存储介质 | |
| Park et al. | Characterization of metal artifacts in X‐ray computed tomography | |
| CN103106676B (zh) | 一种基于低剂量投影数据滤波的x射线ct图像重建方法 | |
| Peng et al. | An irregular metal trace inpainting network for x‐ray CT metal artifact reduction | |
| CN103679706B (zh) | 一种基于图像各向异性边缘检测的ct稀疏角度重建方法 | |
| Zhang et al. | REDAEP: Robust and enhanced denoising autoencoding prior for sparse-view CT reconstruction | |
| Fang et al. | Removing ring artefacts for photon-counting detectors using neural networks in different domains | |
| Zhang et al. | Iterative metal artifact reduction for x‐ray computed tomography using unmatched projector/backprojector pairs | |
| WO2020019684A1 (zh) | 低剂量pet图像重建方法、装置、设备及存储介质 | |
| Liang et al. | Guest editorial low-dose CT: what has been done, and what challenges remain? | |
| CN106056644A (zh) | Ct扫描的数据处理方法及装置 | |
| Yang et al. | Shading correction assisted iterative cone-beam CT reconstruction | |
| CN103793890A (zh) | 一种能谱ct图像的恢复处理方法 | |
| Zhang et al. | Computed tomography sinogram inpainting with compound prior modelling both sinogram and image sparsity | |
| CN110458908A (zh) | 基于有限角度迭代重建超视野ct图像的方法 | |
| CN105844678A (zh) | 基于全广义变分正则化的低剂量x射线ct图像重建方法 | |
| CN109472836A (zh) | 一种ct迭代重建中伪影校正方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |