基于相对TV最小的CT图像重建算法

张家浩; 乔志伟

doi:10.15953/j.ctta.2022.190

基于相对TV最小的CT图像重建算法

doi: 10.15953/j.ctta.2022.190

张家浩,
乔志伟^,

山西大学计算机与信息技术学院, 太原030006

基金项目: 国家自然科学基金面上项目（模型与数据耦合驱动的快速四维EPRI肿瘤氧成像（62071281））；中央引导地方科技发展资金项目（新型TV和学习先验联合约束的快速四维EPRI成像方法（YDZJSX2021A003））；山西省回国留学人员科研资助项目（基于新型四维TV正则机理的快速EPRI肿瘤氧成像方法研究（2020-008））

详细信息

作者简介:
张家浩：男，山西大学计算机科学与技术专业硕士研究生，主要从事医学图像重建、图像处理等方面的研究，E-mail：821967394@qq.com

乔志伟：男，博士，山西大学计算机与信息技术学院教授、博士生导师，主要从事医学图像重建、信号处理、大规模最优化等方面的研究，E-mail：zqiao@sxu.edu.cn

通讯作者:
乔志伟*，

中图分类号: O 242；TP 391.41
计量
- 文章访问数: 90
- HTML全文浏览量: 30
- PDF下载量: 36
- 被引次数: 0
出版历程
- 收稿日期: 2022-09-25
- 修回日期: 2022-10-13
- 录用日期: 2022-10-17
- 网络出版日期: 2022-11-04
- 刊出日期: 2023-03-31

Computed Tomography Reconstruction Algorithm Based on Relative Total Variation Minimization

ZHANG Jiahao,
QIAO Zhiwei^,

School of Computer and Information Technology, Shanxi University, Taiyuan 030006, China

摘要

摘要: 总变差（TV）最小算法是一种有效的CT图像重建算法，可以对稀疏或含噪投影数据进行高精度重建。然而，在某些情况下，TV算法会产生阶梯状伪影。在图像去噪领域，相对TV算法展现了优于TV算法的性能。鉴于此，将相对TV模型引入图像重建，提出相对TV最小优化模型，并在自适应梯度下降-投影到凸集（ASD-POCS）框架下设计对应的求解算法，以进一步提升重建精度。以Shepp-Logan、FORBILD及真实CT图像仿真模体进行重建实验，验证了该算法的正确性并评估了算法的稀疏重建能力和抗噪能力。实验结果表明，相对TV算法可以实现逆犯罪，可以对稀疏或含噪投影数据进行高精度重建，与TV算法相比，该算法可以取得更高的重建精度。
- 计算机断层成像 /
- 相对总变差最小 /
- 稀疏重建 /
- ASD-POCS算法
Abstract: The total variation (TV) minimization algorithm is an effective CT image reconstruction algorithm that can reconstruct sparse or noisy projection data with high accuracy. However, in some cases, the TV algorithm produces stepped artifacts. The relative TV algorithm outperforms TV algorithm in the field of image denoising. In view of this, the relative TV model is introduced into image reconstruction, a relative TV minimum optimization model is proposed, and the corresponding solution algorithm is designed under the framework of adaptive gradient descent projection to the convex set (ASD-POCS) to further improve reconstruction accuracy. The reconstruction experiments were conducted with Shepp Logan, Forbild, and real CT image simulation models to verify the anti-crime ability of the algorithm and evaluate the sparse reconstruction and anti-noise abilities of the algorithm. The experimental findings reveal that the algorithm outperforms the TV method in terms of anti-crime capability and the ability to reconstruct sparse or noisy projection data with high precision. Compared with the TV algorithm, the algorithm can achieve higher reconstruction accuracy.
- computed tomography /
- minimum relative total variation /
- sparse reconstruction /
- ASD-POCS algorithm

HTML全文

图 1 噪声图像的偏导数特点

Figure 1. Partial derivative characteristics of noisy images

下载: 全尺寸图片幻灯片

图 2 RTV算法的正确性研究关于Shepp-Logan模体的重建结果

Figure 2. Correctness research of RTV algorithm on reconstruction results of Shepp-Logan phantom

下载: 全尺寸图片幻灯片

图 3 RTV算法的正确性研究关于FORBILD模体的重建结果

Figure 3. Correctness research of RTV algorithm on reconstruction results of FORBILD phantom

下载: 全尺寸图片幻灯片

图 4 RTV算法重建Shepp-Logan模体收敛行为的展示

Figure 4. Display of convergence behavior when reconstructing Shepp-Logan phantom by RTV algorithm

下载: 全尺寸图片幻灯片

图 5 RTV算法重建FORBILD模体收敛行为的展示

Figure 5. Display of convergence behavior when reconstructing FORBILD phantom by RTV algorithm

下载: 全尺寸图片幻灯片

图 6 不同模体在20个投影角度下使用RTV和TV最小化重建算法的RMSE趋势曲线的比较

Figure 6. Comparison of RMSE trend curves of the same phantom using RTV and TV minimization reconstruction algorithms at 20 projection angles

下载: 全尺寸图片幻灯片

图 7 RTV算法和TV算法对于FORBILD模体重建结果的比较：图像上方的数字表示投影个数；左边文字表示使用的算法

Figure 7. Comparison of FORBILD phantom reconstruction results between the RTV algorithm and the TV algorithm: the number above the image represents the number of projections; the text on the left indicates the algorithm used

下载: 全尺寸图片幻灯片

图 8 RTV算法和TV算法对于Shepp-Logan模体重建结果的比较：图像上方的数字表示投影个数；左边文字表示使用的算法

Figure 8. Comparison of Shepp-Logan phantom reconstruction results between the RTV algorithm and the TV algorithm: the number above the image represents the number of projections; the text on the left indicates the algorithm used

下载: 全尺寸图片幻灯片

图 9 RTV算法和TV算法对于真实CT图像重建结果的比较：图像上方的数字表示投影个数；左边文字表示使用的算法

Figure 9. Comparison of CT-phantom reconstruction results between RTV algorithm and TV algorithm: the number above the image represents the number of projections; The text on the left indicates the algorithm used

下载: 全尺寸图片幻灯片

图 10 在50个投影角度并于投影数据中加入方差为0.05的高斯白噪声条件下分别使用RTV和TV重建算法进行重建的实验结果

Figure 10. Under the condition of 50 projection angles and adding Gaussian white noise with a variance of 0.05 to the projection data, RTV and TV minimization reconstruction algorithms are used to reconstruct the results, respectively

下载: 全尺寸图片幻灯片

图 11 不同模体在50个投影角度并于投影数据中加入方差为0.05的高斯白噪声条件下使用RTV和TV最小化重建算法重建结果的RMSE趋势曲线的比较

Figure 11. Comparison of RMSE trend curves of different phantom reconstruction results the using RTV and TV minimization reconstruction algorithms with 50 projection angles and adding Gaussian white noise with variance of 0.05 added to the projection data

下载: 全尺寸图片幻灯片

图 12 λ不同取值的实验结果和其各自对应的中心线波形图

Figure 12. Experimental results with different λ values and their corresponding centerline waveforms

下载: 全尺寸图片幻灯片

图 13 ε不同取值的实验结果和其各自对应的中心线波形图

Figure 13. Experimental results with different ε values and their corresponding centerline waveforms

下载: 全尺寸图片幻灯片

1. Repeat main loop		13. $for{\text{ } }i = 1{\text{:} }{N_{{\rm{grad}}} }{\text{ do RTV - ASD loop} }$：
2. $ {{\boldsymbol{f}}}_{p}={{\boldsymbol{f}}}^{(n)}$		${\boldsymbol{d} }{ {\boldsymbol{f} }^{(k)} } = {\nabla _f}\left\\| { { {\boldsymbol{f} }^{(n)} } } \right\\|_{{\rm{RTV}}}^{(k)}$
3. ${\text{for } }i = 1{\text{ : } }{N_d}{\text{ do } }{\boldsymbol{f} } = {\boldsymbol{f} } + \beta { {\boldsymbol{A} }_{\boldsymbol{i} } }\displaystyle \frac{ { {g_i} - { {\boldsymbol{A} }_{\boldsymbol{i} } }{\boldsymbol{ \times f} } } }{ { { {\boldsymbol{A} }_{\boldsymbol{i} } }{\boldsymbol{ \times } }{ {\boldsymbol{A} }_{\boldsymbol{j} } } } }$		$ {\boldsymbol{d}}{{\boldsymbol{f}}^{(k)}} = {{{\boldsymbol{d}}{{\boldsymbol{f}}^{(k)}}} \mathord{\left/ {\vphantom {{{\boldsymbol{d}}{{\boldsymbol{f}}^{(k)}}} {{{\left\\| {{\boldsymbol{d}}{{\boldsymbol{f}}^{(k)}}} \right\\|}_2}}}} \right. } {{{\left\\| {{\boldsymbol{d}}{{\boldsymbol{f}}^{(k)}}} \right\\|}_2}}} $
4. ${\text{for } }i = 1{\text{ : } }{N_i}{\text{ do if } }{{\boldsymbol{f}}_i} < 0{\text{ then } }{{\boldsymbol{f}}_i} = 0$		$ {{\boldsymbol{f}}^{(n)}} = {{\boldsymbol{f}}^{(n)}} - {d_\alpha }*{\boldsymbol{d}}{{\boldsymbol{f}}^{(k)}} $
5. ${ {\boldsymbol{f} }_{{\rm{res}}} } = { {\boldsymbol{f} }^{(n)} }$		end for
6. ${{\boldsymbol{g}}^{(n)}} = {\boldsymbol{A}} \times {{\boldsymbol{f}}^{(n)}}$		14. $ \nabla {{\boldsymbol{f}}_{{\text{RTV}}}} = {\left\\| {{{\boldsymbol{f}}^{(n)}} - {{\boldsymbol{f}}_p}} \right\\|_2} $
7. ${d_p} = {\left\\| {{{\boldsymbol{g}}^{(n)}}{\boldsymbol{ - }}{{\boldsymbol{g}}_{\boldsymbol{0}}}} \right\\|_2}$		15. ${\rm{if}}{\text{ } }\nabla { {\boldsymbol{f} }_{ {\text{RTV} } } } > {r_{\max} }*\nabla { {\boldsymbol{f} }_{ {\text{POCX} } } }{\text{ } }{\rm{and}}{\text{ } }{d_p} > \in {\text{ then} }$
10. $\nabla {{\boldsymbol{f}}_{{\text{POCX}}}} = {\left\\| {{{\boldsymbol{f}}^{(n)}} - {{\boldsymbol{f}}_p}} \right\\|_2}$		${d_\alpha } = {d_\alpha }*{\alpha _{{\rm{red}}} }$
11. ${\text{if first iterarion , then }}{d_\alpha } = \alpha \cdot \nabla {{\boldsymbol{f}}_{{\text{POCX}}}}$		16. $\beta = \beta \cdot {\beta _{{\rm{red}}} }$
12. ${{\boldsymbol{f}}_p} = {{\boldsymbol{f}}^{(n)}}$		17. ${\text{until}}\left\{ {{\text{stopping criteria}}} \right\}$
		18. ${\text{return } }\;{{\boldsymbol{f}}_{ {\rm{res} } } }$

下载: 导出CSV

表 1 RTV算法和TV算法重建FORBILD模体的RMSE和SSIM比较

Table 1. RMSE and SSIM comparison of RTV and TV algorithms for FORBILD phantom reconstruction

项目	算法	投影个数
项目	算法	20	30	40	50
RMSE	RTV	30.0×10^-5	9.96×10^-5	7.11×10^-5	5.60×10^-5
RMSE	TV	5130×10^-5	2530×10^-5	620×10^-5	590×10^-5
SSIM	RTV	0.9587	0.9792	0.9972	0.9992
SSIM	TV	0.9293	0.9583	0.9943	0.9986

下载: 导出CSV

表 2 RTV算法和TV算法重建Shepp-Logan模体的RMSE和SSIM比较

Table 2. RMSE and SSIM comparison of the RTV and TV algorithms for Shepp-Logan phantom reconstruction

项目	算法	投影个数
项目	算法	20	30	40	50
RMSE	RTV	8.010×10^-5	4.791×10^-5	3.349×10^-5	2.586×10^-5
RMSE	TV	1610×10^-5	530×10^-5	160×10^-5	80×10^-5
SSIM	RTV	0.9866	0.9948	0.9978	0.9995
SSIM	TV	0.9812	0.9923	0.9956	0.9991

下载: 导出CSV

表 3 RTV算法和TV算法重建真实CT模体的RMSE和SSIM比较

Table 3. RMSE and SSIM comparison of the RTV and TV algorithms for CT phantom reconstruction

项目	算法	投影个数
项目	算法	20	30	40	50
RMSE	RTV	0.0403	0.0294	0.0210	0.0120
RMSE	TV	0.0467	0.0302	0.0218	0.0168
SSIM	RTV	0.9213	0.9546	0.9743	0.9956
SSIM	TV	0.8977	0.9387	0.9542	0.9842

下载: 导出CSV

表 5 RTV算法和TV算法在不同等级噪声下重建Shepp-Logan模体的RMSE和SSIM比较

Table 5. Comparison of RMSE and SSM of the RTV and TV algorithms for reconstructing Shepp-Logan phantom under different levels of noise

项目	算法	噪声等级
项目	算法	0.01	0.02	0.03	0.04	0.05
RMSE	RTV	0.0015	0.0028	0.0035	0.0048	0.0055
RMSE	TV	0.0065	0.0081	0.0102	0.0126	0.0148
SSIM	RTV	0.9943	0.9823	0.9781	0.9653	0.9526
SSIM	TV	0.9914	0.9814	0.9727	0.9594	0.9512

下载: 导出CSV

表 4 RTV算法和TV算法在不同等级噪声下重建FORBILD模体的RMSE和SSIM比较

Table 4. Comparison of RMSE and SSM of the RTV and TV algorithms for reconstructing FORBILD phantom under different levels of noise

项目	算法	噪声等级
项目	算法	0.01	0.02	0.03	0.04	0.05
RMSE	RTV	0.0021	0.0029	0.0032	0.0043	0.0050
RMSE	TV	0.0085	0.0107	0.0134	0.0162	0.0184
SSIM	RTV	0.9923	0.9843	0.9812	0.9756	0.9712
SSIM	TV	0.9836	0.9793	0.9754	0.9698	0.9652

下载: 导出CSV

表 6 RTV算法和TV算法在不同等级噪声下重建真实CT模体的RMSE和SSIM比较

Table 6. Comparison of RMSE and SSM of the RTV and TV algorithms for reconstructing CT phantom under different levels of noise

项目	算法	噪声等级
项目	算法	0.01	0.02	0.03	0.04	0.05
RMSE	RTV	0.0188	0.0198	0.0205	0.0225	0.0233
RMSE	TV	0.0194	0.0215	0.0233	0.0247	0.0286
SSIM	RTV	0.9258	0.9042	0.8938	0.8814	0.8715
SSIM	TV	0.9123	0.8926	0.8626	0.8612	0.8523

下载: 导出CSV

参考文献(22)

[1]	PAN X, SIDKY E Y, VANNIER M. Why do commercial CT scanners still employ traditional, filtered back-projection for image reconstruction?[J]. Inverse Problems, 2008, 25(12): 1230009.
[2]	SIDKY E Y, PAN X. Image reconstruction in circular cone-beam computed tomography by constrained, total-variation minimization[J]. Physics in Medicine & Biology, 2008, 53(17): 4777−4807.
[3]	YU L, YU Z, SIDKY E Y, et al. Region of interest reconstruction from truncated data in circular cone-beam CT[J]. IEEE Transactions on Medical Imaging, 2006, 25(7): 869−881. doi: 10.1109/TMI.2006.872329
[4]	CANDES E, ROMBERG J, TAO T. Robust uncertainty principles: Exact signal reconstruction from highly incomplete frequency information[J]. IEEE Transactions on Information Theory, 2006, 52(2): 489-509.
[5]	VOGEL C R. Iterative method for total variation denoising[J]. SIAM Journal on Scientific Computing, 1996, 17(1): 227−238. doi: 10.1137/0917016
[6]	SIDKY E Y, EMIL Y, KAO C M, et al. Accurate image reconstruction from few-views and limited-angle data in divergent-beam CT[J]. Journal of X-ray Science & Technology, 2006, 14: 119−139.
[7]	DONOHO D L. Compressed sensing[J]. IEEE Transactions on Information Theory, 2006, 52: 1289-1306.
[8]	LIU B, KATSEVICH A, YU H. Interior tomography with curvelet-based regularization[J]. Journal of X-ray Science and Technology, 2016, 25(1): 1−13.
[9]	SIDKY E Y, CHARTRAND R, BOONE J M, et al. Constrained TpV minimization for enhanced exploitation of gradient sparsity: Application to CT image reconstruction[J]. IEEE Journal of Translational Engineering in Health & Medicine, 2014, 2(6): 1−18.
[10]	CHAMBOLLE A, POCK T. An introduction to continuous optimization for imaging[J]. Acta Numerica, 2016, 25: 161-319.
[11]	MBOLLE A, POCK T. A first-order primal-dual algorithm for convex problems with applications to imaging[J]. Journal of Mathematical Imaging and Vision, 2011, 40(1): 120−145. doi: 10.1007/s10851-010-0251-1
[12]	SIDKY E Y, JRGENSEN J H, PAN X. Convex optimization problem prototyping for image reconstruction in computed tomography with the Chambolle-Pock algorithm[J]. Physics in Medicine & Biology, 2011, 57(10): 3065−3091.
[13]	XU L, YAN Q, XIA Y, et al. Structure extraction from texture via relative total variation[J]. ACM Transactions on Graphics, 2012, 31(6): 139.
[14]	RIGIE D S, RIVIèRE P L. Joint reconstruction of multi-channel, spectral CT data via constrained total nuclear variation minimization[J]. Physics in Medicine & Biology, 2015, 60(5): 1741.
[15]	QIAO Z, REDLER G, EPEL B, et al. A balanced total-variation-Chambolle-Pock algorithm for EPR imaging[J]. Journal of Magnetic Resonance, 2021, 328: 107009. doi: 10.1016/j.jmr.2021.107009
[16]	闫慧文, 乔志伟. 基于ASD-POCS框架的高阶TpV图像重建算法[J]. CT理论与应用研究, 2021,30(3): 279−289. DOI: 10.15953/j.1004-4140.2021.30.03.01. YAN H W, QIAO Z W. High order TPV image reconstruction algorithm based on ASD-POCS framework[J]. CT Theory and Applications, 2021, 30(3): 279−289. DOI: 10.15953/j.1004-4140.2021.30.03.01. (in Chinese).
[17]	XIN J, LIANG L, CHEN Z, et al. Anisotropic total variation minimization method for limited-angle CT reconstruction[C]//Proc. SPIE 8506, Developments in X-ray Tomography VIII, 85061C, 2012. https://doi.org/10.1117/12.930339.
[18]	YU H, WANG G. A soft-threshold filtering approach for reconstruction from a limited number of projections[J]. Physics in Medicine & Biology, 2010, 55(13): 3905−3916.
[19]	乔志伟. 总变差约束的数据分离最小图像重建模型及其Chambolle-Pock求解算法[J]. 物理学报, 2018,67(19): 198701. DOI: 10.7498/aps.67.20180839. QIAO Z W. Total variation constrained data separation minimum image reconstruction model and its Chambolle-Pock algorithm[J]. Acta Physica Sinica, 2018, 67(19): 198701. DOI: 10.7498/aps.67.20180839. (in Chinese).
[20]	YU Z, NOO F, DENNERLEIN F, et al. Simulation tools for two-dimensional experiments in X-ray computed tomography using the FORBILD head phantom[J]. Physics in Medicine & Biology, 2012, 57(13): 237−252.
[21]	SIDDON R L. Fast calculation of the exact radiological path for a three-dimensional CT array[J]. Medical Physics, 1985, 12(2): 252−255. doi: 10.1118/1.595715
[22]	GORDON R, BENDER R, HERMAN G T. Algebraic reconstruction techniques (ART) for three-dimensional electron microscopy and X-ray photography[J]. Journal of Theoretical Biology, 1970, 29(3): 471−481. doi: 10.1016/0022-5193(70)90109-8