Sparse View CT Reconstruction Algorithm Based on Non-Local Generalized Total Variation Regularization
-
摘要:
基于广义总变分（TGV）正则化的CT图像重建算法可以有效克服总变分（TV）正则化的阶梯效应，从而能保护重建图像过渡区域的结构特征。尽管TGV重建方法优于TV重建方法，但它仍然忽略了非局部自相似先验信息在恢复CT图像细节方面的显著作用。为了克服TGV重建方法的上述局限性，本文引入一种非局部广义总变分（NLTGV）正则项，并提出基于NLTGV正则化的稀疏角度CT重建算法。该方法不仅可以利用不同阶的非局部变分信息来保护图像结构特征，而且还可以利用非局部自相似性来恢复重建图像的细节。由于重建模型包含双非光滑项，难以直接求解，因此提出基于凸集投影的优化算法，将其分解为几个简单子问题实现有效求解。仿真和实验结果表明，与其他变分正则化重建方法相比，本文重建方法可以有效提高CT图像重建质量。
-
关键词:
- X射线CT /
- 稀疏角度采样 /
- 非局部广义全变分 /
- 凸集投影算法 /
- 分裂Bregman算法
Abstract:CT image reconstruction algorithm based on generalized total variation (TGV) can overcome the staircase effect of total variation (TV) regularization, thereby protecting the structural features of the reconstructed image transition region. Although the TGV reconstruction method is superior to the TV reconstruction method, it still ignores the role of non-local self-similar prior information in restoring CT image details. To overcome the aforementioned limitations of TGV reconstruction method, we introduce a non-local TGV (NLTGV) regularization term and propose a sparse view CT reconstruction algorithm based on NLTGV regularization. The proposed method can not only utilize non-local variational information of different orders to protect image structural features but can also utilize non-local self-similarity to restore the details of the reconstructed image. Owing to the inclusion of dual non-smooth terms in the reconstruction model, solving it directly is difficult. Therefore, we proposed an optimization algorithm based on convex set projection, which decomposes the problem into several sub-problems to be solved. The simulation and experimental results show that the proposed NLTGV regularization reconstruction method can effectively improve the quality of reconstructed images compared with other variational reconstruction methods.
-
X射线计算机断层扫描（computed tomography，CT）作为一种先进的无损成像技术，通过重建断层图像来可视化物体的内部结构，在医学影像诊断、工业无损检测、安测和逆向工程等领域有着广泛的应用。重建断层图像质量的好坏对后续的分析处理起到关键性作用。传统的解析重构方法如FBP算法需要足够多的投影数据来获得高质量的断层图像。众所周知，X射线对人体有辐射危害，在医学成像领域降低辐射剂量具有重要意义。一种常用的降低辐射剂量方法是将传统的密集投影角度采集模式转变为稀疏角度采样模式[1]。然而，利用稀疏角度采样投影数据重建CT图像是一个病态的逆问题，解析算法重建图像会出现条状伪影，显著降低CT图像质量，影响图像的可用性。在稀疏角度采样条件下，如何提高CT图像重建质量对于实际应用具有重要意义。
为了提高稀疏角度CT图像的重建质量，通常需引入先验信息，采用基于正则化的迭代算法进行重建。压缩感知理论[2]表明，利用稀疏性先验信息可有效减少重构原始信号所需的采样量。虽然CT图像本身不一定满足稀疏性，但通过适当的稀疏变换，可将其转换为稀疏信号。由于变分正则化在提高梯度域稀疏性方面的显著效果，在稀疏角度CT重建中得到了非常广泛的应用。
Sidky等[3-4]提出了自适应梯度下降−凸集投影算法（adaptive steepest descent-projection onto convex sets，ADS-POCS），求解全变分（total variation，TV）[5]最小化问题，用于解决投影数据不完全的CT重建问题。与传统的ART[6]和EM[7]相比，TV方法能有效抑制条状伪影，从而显著提高CT图像重建质量。
尽管TV方法取得了较好的效果，但其基于分段常量假设，重建图像易产生块状伪影，不利于保护重建图像中一些平滑过渡的结构特征。为了克服TV方法的上述弱点，有必要在正则化项中引入重建图像的高阶变分信息[8]。2014年，Niu等[9]基于全广义变分（total generalized variation，TGV）提出了惩罚加权最小二乘与TGV相结合（penalized weighted least-squares，PWLS-TGV）的重建算法，相比于TV方法能更好地保留图像细节信息。2015年，Sun等[10]提出利用Hessian矩阵的Frobenius范数作为正则项进行CT重建，可以有效抑制TV惩罚项的阶梯效应。2016年，Zhang等[11]提出了一种基于Lp范数保留边缘信息的p阶广义全变分正则化模型，能避免少角度投影情况下TV方法重建CT图像中的阶梯效应。2020年，Xi等[12]提出了一种高阶全变分（high order total variation，HOTV）算法，其研究表明高阶变分项能有效压制TV模型引入的阶梯伪影，同时能够用较少的投影数据重建出精度更高的图像。2021年，闫慧文等[13]引入高阶TpV（high-order total p-variation，HOTpV）正则化重建算法，利用Lp范数进一步提升了HOTV的重建效果。最近，Xi等[14]研究了自适应高阶全变分（adaptive-weighted high order TV，awHOTV），其研究表明与TV算法相比能取得更准确的重建结果。
除了利用不同阶数变分正则项的稀疏性之外，利用图像的非局部相似性被证明是进一步提升稀疏性的有效方法[15]。Kim等[16-17]提出一种基于非局部全变分（nonlocal total variation，NLTV）[18]的CT重建算法，取得了比TV方法更好的重建效果，然而仍然没有考虑到高阶变分信息。Ranftl等[19]提出了一种非局部广义总变分（nolocal total generalized total variation，NLTGV）正则项，能同时建模高阶变分的稀疏性与非局部自相似性，在计算机视觉深度任务估计中取得优良的效果。然而，目前在CT重建领域仍然缺少基于NLTGV的高阶非局部变分正则化图像重建方法研究。
为了克服现有高阶变分正则化CT重建方法未考虑非局部自相似性的问题，本文基于上述的NLTGV正则项，提出一种非局部高阶变分正则化稀疏角度CT图像重建算法，并且设计一种基于凸集投影与分裂布雷格曼方法（split Bregman）算法相结合的高效求解方案。由于NLTGV能够利用非局部高阶变分，因此所提方法可有效克服阶梯效应，并且可利用非局部相似性进一步提升重建图像的细节质量。仿真和实际数据研究表明，与现有的TV、TGV、NLTV变分正则化方法相比，该方法在噪声抑制和细节保持方面具有优势，可有效提高CT图像重建质量。
1. CT成像模型
CT成像的基本原理是通过测量X射线在不同方向上穿过物体的衰减量来重建物体内部的截面图像。假设X射线是单色的，X射线沿着传播路径L穿过被检测物体的衰减规律在数学上由比尔定律描述，
$$ I = {I_0}\exp \left({{{ - \displaystyle\int\limits_L {u(l){\rm d}l} }} }\right)\text{，} $$ (1) 其中，
$ {I_0} $ 是X射线源发射的初始强度，$ I $ 是X射线穿过物体后的强度，$ u(l) $ 是沿着传播路径的X射线衰减系数。为了获得重建断层图像，需要将待重建区域离散化为
$ {N_1} \times N{}_2 $ 个矩形网格，每个矩形网格中的X射线衰减系数可视为常数，定义为待重建图像的像素，第j个像素的数值表示为$ {u_j} $ 。因此，待重建图像可矢量化为$ {\boldsymbol{u}} = ({u_1},{u_2},\cdots,{u_j},\cdots,{u_N}) $ ，其中$ N = {N_1} \times N{}_2 $ 。式（1）可以离散化为代数方程，第i条 X射线传播路径所对应的投影
$ {p_i} $ 表示如下：$$ \sum\limits_{j = 1}^N {{a_{i,j}}} {u_j} = {p_i} = - \ln \left(\frac{{{I_i}}}{I}\right) \text{，} $$ (2) 其中，
$ {a_{i,j}} $ 表示像素$ {u_j} $ 的加权系数。将X射线穿过物体的所有路径均考虑在内，可得CT成像模型为：
$$ {\boldsymbol{Au}} = {\boldsymbol{p}} \text{，} $$ (3) 其中，
$ {\boldsymbol{A}} \in {R^{M \times N}} $ 是投影矩阵，$ {a_{i,j}} \geq 0 $ 为矩阵的第i行第j列位置处元素。矢量$ {\boldsymbol{p}} = ({p_1},{p_2}, \cdot \cdot \cdot , {p_i}, \cdot \cdot \cdot , {p_M})^{\mathrm{T}} $ 代表所有的投影数据，$ M $ 为投影测量数。CT重建问题可归结为基于式（3）估计CT断层图像$ {\boldsymbol{u}} $ 。2. 基于NLTGV正则化CT重建算法
2.1 NLTGV正则化
TV基于分段常数假设，容易导致CT重建图像存在块状伪影，广义总变分（TGV）进一步利用高阶变分信息来缓解上述问题，但仍然忽视了CT图像的非局部自相似信息。作为TGV的扩展，Ranftl等[19]提出非局部广义总变分正则项，定义为：
$$ \big|\big|{\boldsymbol{u}}\big|{\big|_{{\text{NLTGV}}}} = \mathop {\min }\limits_{{\boldsymbol{u}},{\boldsymbol{g}}}\Biggr( \int\limits_\Omega \int\limits_\Omega {w_1}(x,y)\big|u(x) - u(y) -\langle g(x),x - y \rangle {\mathrm{d}}x{\mathrm{d}}y \big|+ \sum\limits_{i = 1}^2 {\int\limits_\Omega {\int\limits_\Omega {{w_2}(x,y)\big|{g^i}(x) - {g^i}(y)\big|{\mathrm{d}}x{\mathrm{d}}y} } } \Biggr) ， $$ (4) 其中，
$ {\boldsymbol u}：\Omega \to R，g = {({g^1},{g^2})^{\mathrm{T}}}:\Omega \to {R^2} $ 分别表示图像和非局部变分，符号$ \left\langle , \right\rangle $ 表示内积运算。$ {w_1},{w_2}： \Omega \times \Omega \to {R^ + } $ 是支撑权重，度量任意两个像素的相似性，具体表示如下：$$ \left\{\begin{aligned} & {w_1}(x,y) = \frac{1}{{{C^i}}}\exp \left( - \frac{{{{\big\| {P(x) - P(y)} \big\|}_2}}}{{{h_u}}}\right)\times\\ &\qquad\qquad\exp \left( - \frac{{{{\big\| {x - y} \big\|}_2}}}{{{h_l}}}\right) \\ & {w_2}(x,y) = \gamma {w_1}(x,y) \end{aligned}\right. \text{，} $$ (5) 其中
$ {C^i} $ 是归一化因子，可以确保$ 0 \lt {w_1}(x,y) \leq 1 $ 。$ P(x) $ 和$ P(y) $ 分别代表以像素$ u(x) $ 和像素$ u(y) $ 为中心的图像块。对于任一像素，在搜索窗口Q（以该像素为中心的正方形窗口）中搜索其相似的像素，示意图如图1。权重$ {w_1}、{w_2} $ 的衰减快慢由正数的参数$ {h_u}、{h_l} $ 控制，权重$ {w_1} $ 构成类似于著名的双边滤波器[20]的滤波权重组成，$ {h_u} $ 控制权重$ {w_1} $ 随两个像素数值接近程度的变化情况，参数$ {h_l} $ 控制权重$ {w_1} $ 受两个像素空间距离的影响程度，可设置为半个搜索窗口的大小[20]。参数$ \gamma \gt 0 $ 是一个比例因子，对NLTGV的两个加权项起到平衡作用。![]()
图 1 搜索非局部相似图像块原理Figure 1. Schematic illustrating the searching of non-local similarity patchesNLTGV正则项有利于促进分段仿射数值解，因此NLTGV正则项有利于保护CT图像中的过渡区域特征，并且NLTGV使用权重
$ {w_1}、{w_2} $ 来引入非局部自相似先验信息。表1给出了NLTGV与其他广泛使用的变分正则项的特性对比。表 1 不同变分正则项特性的比较Table 1. Comparison of characteristics with different variation regularization算法 TV TGV NLTV NLTGV 邻域 局部 局部 非局部 非局部 变分阶次 一阶 高阶 一阶 高阶 综上，NLTGV是一种非局部的混合变分正则化项，能借助不同阶数的变分项保护CT图像渐变区域的特征。此外，NLTGV还充分考虑了CT图像的非局部相似性先验信息，能有效利用非局部结构相似性还原图像细节。
2.2 NLTGV正则化重建模型及求解
本文基于NLTGV正则项构建重建模型，并通过求解如下最小化问题来重建CT图像：
$$\left\{ \begin{aligned} &\mathop {\min }\limits_{\boldsymbol{u}} \big|\big|{\boldsymbol{u}}\big|{\big|_{{\text{NLTGV}}}} \\ &\begin{aligned} & {{\text{s}}{\text{.t.}}}\quad\begin{aligned} & {\boldsymbol{Au}} = {\boldsymbol{p}} \\ & {\boldsymbol{u}} \geq {\boldsymbol{0}} \end{aligned} \end{aligned} \end{aligned} \right. 。 $$ (6) 通过使用NLTGV正则项，所提重建模型可有效地同时利用高阶TV的稀疏性与非局部自相似冗余先验信息来提高噪声抑制和结构细节恢复的能力。重建模型（6）中优化目标函数是不可导的，直接求解较为困难，可利用基于变量分裂技术的优化方法，将其分裂为几个较为简单的子问题进行求解。受到ADS-POCS算法的启发，本文将问题（6）分裂为
$ {{\text{P}}_1} $ 和$ {{\text{P}}_2} $ 两个子问题进行交替求解，如式（7）。$ {{\text{P}}_1} $ 步骤可以看做凸集投影步骤，使重建图像满足成像模型与像素值非负约束,而$ {{\text{P}}_2} $ 步骤将经（simultaneous algebraic reconstruction technique，SART）联合代数重建技术迭代[21]与非负约束后的估计结果作为初始值，通过NLTGV正则化来抑制伪影与噪声，$$\left\{ \begin{aligned} &{{\text{P}}_1}:\left\{ \begin{aligned} &{\boldsymbol{f}} = {\text{SART}}({{\boldsymbol{u}}^{k - 1}}) \\ &{{\boldsymbol{u}}^{{\text{pos}}}} = {\text{pos}}({\boldsymbol{f}}) \end{aligned} \right. \\ & {{\text{P}}_2}:{{\boldsymbol{u}}^k} = \mathop {\arg \min }\limits_{\boldsymbol{u}} \left(\big|\big|{\boldsymbol{u}}\big|{\big|_{{\text{NLTGV}}}} + \frac{\mu }{2}{\big\| {{\boldsymbol{u}} - {{\boldsymbol{u}}^{{\text{pos}}}}} \big\|^2} \right) \end{aligned}\right. , $$ (7) 其中，SART代表联合代数算法迭代步骤，符号
$ {\text{pos}}( \cdot ) $ 代表对像素值取非负运算，$ {{\boldsymbol{u}}^{{\text{pos}}}} $ 为经SART算法迭代并进行非负约束后的图像，参数$ \mu $ 用来控制NLTGV抑制伪影和噪声的强度。具体求解过程如下：
（1）SART迭代
为了加快凸集投影算法的收敛速度，本文采用SART迭代步骤替代ADS-POCS算法中的代数重建技术（algebraic reconstruction technique，ART）迭代步骤，求解一致性约束，得到中间图像
$ {\boldsymbol{f}} $ ，令本轮SART算法迭代的初值$ {{\boldsymbol{f}}^0} = {{\boldsymbol{u}}^{k - 1}} $ ，即$ {{\boldsymbol{f}}^0} $ 为上一轮迭代的NLTGV稀疏正则化问题的解，则有：$$ \left\{\begin{aligned} &f_j^t = f_j^{t - 1} + \frac{\lambda }{{{A_{ + ,j}}}}\sum\limits_{i = 1}^M {\frac{{{a_{i,j}}}}{{{A_{i, + }}}}} \big({p_i} - {{\bar p}_i}\big) \\ &{A_{i, + }} = \sum\limits_{j = 1}^N {{a_{i,j}}} ,\,i = 1,2,\cdots,M, \\ &{A_{ + ,j}} = \sum\limits_{i = 1}^M {{a_{i,j}}} ,\,j = 1,2,\cdots,N, \\ &{{\bar p}_i} = {A_i}{f^{t - 1}} \end{aligned}\right. , $$ (8) 其中，迭代次数
$ t = 1,2,3, \cdots ,{t_{\max }} $ ，$ {t_{\max }} $ 为SART算法的最大迭代次数，$ {f_j} $ 为$ {\boldsymbol{f}} $ 的第j个分量，$ 0 \lt \lambda \leq 1 $ 为松弛因子，$ {A_{i, + }} $ 为矩阵$ {\boldsymbol{A}} $ 的第i行元素之和，$ {A_{ + ,j}} $ 为矩阵$ {\boldsymbol{A}} $ 的第j列元素之和。（2）非负约束
将重建的中间图像向非负空间投影，即：
$$ u_j^{{\text{pos}}} = \left\{ \begin{aligned} & f_j^{{t_{\max }}}, &&{\mathrm{if}}\;\;{f_j^{{t_{\max }}} \geq 0} \\ &0,&&{{\mathrm{if}}}\;\;{f_j^{{t_{\max }}} \leq 0} \end{aligned} \right. 。 $$ (9) （3）求解NLTGV稀疏正则化问题
通过对上述（7）中
$ {{\text{P}}_2} $ 问题的离散化，得下式：$$\mathop {\min }\limits_{{\boldsymbol{u}},{\boldsymbol{g}}}\Bigg( \sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {{u_i} - {u_j} - l_{ij}^T{g_i}} \right\|_1} +\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_2}(i,j)} } {\left\| {{g_i} - {g_j}} \right\|_1} + \frac{\mu }{2}{\big\| {{\boldsymbol{u}} - {{\boldsymbol{u}}^{{\text{pos}}}}} \big\|^2}\Bigg) \text{，} $$ (10) 其中，
$ {{\boldsymbol{l}}_{ij}} = {\left(l_{ij}^1,l_{ij}^2\right)^{\mathrm{T}}} $ ，$ l_{i,j}^1 = {x_i} - {x_j} $ 与$ l_{i,j}^2 = {y_i} - {y_j} $ 分别表示像素$ {u_i} $ 和$ {u_j} $ 间的水平距离和垂直距离，$ {Q_i} $ 代表以像素点$ {u_i} $ 为中心的搜索窗内像素索引。式（10）可进一步写成：
$$ \mathop {\min }\limits_{{\boldsymbol{u}},{\boldsymbol{g}}} \Bigg({\alpha _1}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {{u_i} - {u_j} - l_{ij}^T{g_i}} \right\|_1}+ {\alpha _2}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {{g_i} - {g_j}} \right\|_1} + \frac{1}{2}{\big\| {{\boldsymbol{u}} - {{\boldsymbol{u}}^{{\text{pos}}}}} \big\|^2}\Bigg) \text{，}$$ (11) 其中，
$ {\alpha _1} = 1/\mu $ ,$ {\alpha _2} = \gamma /\mu $ 为正则化参数。式（11）属于L1范数优化问题，由于L1范数非光滑不可导，难以采用经典的最速下降算法，共轭梯度算法等方法直接求解，一般可以采用交替乘子法，原始对偶法，临近梯度算法等基于分裂的优化方法求解。原始对偶法需要写出NLTGV的对偶形式，较为复杂，而临近梯度算法收敛速度较慢。因此为了有效求解式（11），本文采用Bregman算法[22]求解，属于交替乘子法的一种特殊迭代形式，专门为求解大型L1范数等非光滑最优化问题设计，其核心思想是将原问题分裂为一系列子问题依次求解，从而降低整个问题的难度。为了将Bregman算法应用于问题（11），首先引入变量s和r，将问题（11）转为如下等价的形式：
$$\left\{ \begin{aligned} &\mathop {\min }\limits_{{\boldsymbol{u}},{\boldsymbol{g}}} \Bigg({\alpha _1}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {{s_{ij}}} \right\|_1} +\\ &\qquad {\alpha _2}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {{r_{ij}}} \right\|_1} + \\ &\qquad\frac{1}{2}{\big|\big| {{\boldsymbol{u}} - {{\boldsymbol{u}}^{_{{\text{pos}}}}}} \big|\big|^2} \Bigg) \\ &{\text{s}}{\text{.t}}{\text{.}}\quad\begin{aligned} &{{s_{i,j}} = {u_i} - {u_j} - {\boldsymbol l}_{ij}^{\rm T}{{\boldsymbol g}_i}} \\&{{{\boldsymbol r}_{ij}} = {{\boldsymbol g}_i} - {{\boldsymbol g}_j}} \end{aligned} \end{aligned}\right. 。 $$ (12) 根据文献[22]中的Bregman算法迭代步骤，求解问题（12）的迭代式为：
$$ \begin{aligned} \;& {{\boldsymbol{u}}^q},{{\boldsymbol{g}}^q} = \mathop {\arg \min }\limits_{{\boldsymbol{u}},{\boldsymbol{g}}} \Biggr({\alpha _1}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {s_{_{i,j}}^{q - 1}} \right\|_1} + {\alpha _2}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {r_{_{i,j}}^{q - 1}} \right\|_1} + \frac{1}{2}{\big\| {{\boldsymbol{u}} - {{\boldsymbol{u}}^{_{{\text{pos}}}}}} \big\|^2}+ \\ & \;\;\;\;\frac{\rho }{2}\big\|s_{i,j}^{q - 1} - (u_i^{} - u_j^{} - l_{ij}^Tg_i^{q - 1}) + b_{i,j}^{q - 1}\big\|^2 +\frac{\rho }{2}\big\|r_{_{i,j}}^{q - 1} - (g_i^{} - g_j^{}) + d_{i,j}^{q - 1}\big\|^2\Biggr) \end{aligned} \text{，} $$ (13) $$ s_{i,j}^q = \mathop {\arg \min }\limits_{{{{\boldsymbol s}_{{\boldsymbol{i}},{\boldsymbol{j}}}}}} \Bigg({\alpha _1}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {{s_{i,j}}} \right\|_1} +\frac{\rho }{2}\big\|{s_{i,j}} - (u_i^q - u_j^q - l_{ij}^Tg_i^q) + b_{i,j}^{q - 1}\big\|^2\Bigg) \text{，} $$ (14) $$r_{i,j}^q = \mathop {\arg \min }\limits_{{{{\boldsymbol s}_{{\boldsymbol{i}},{\boldsymbol{j}}}}}}\Bigg( {\alpha _2}\sum\limits_{i = 1}^N {\sum\limits_{j \in {Q_i}}^{} {{w_1}(i,j)} } {\left\| {{r_{i,j}}} \right\|_1} +\frac{\rho }{2}\big\|{r_{i,j}} - (g_i^q - g_j^q) + d_{i,j}^{q - 1}\big\|^2\Bigg) \text{，} $$ (15) $$ b_{i,j}^q = b_{i,j}^{q - 1} + s_{i,j}^q - \left(u_i^q - u_j^q - l_{ij}^{\rm T}g_i^q\right) \text{，} $$ (16) $$ d_{i,j}^q = d_{i,j}^{q - 1} + r_{i,j}^q - \left(g_i^q - g_j^q\right) \text{，} $$ (17) 其中，
$ {b_{i,j}} $ ，$ {d_{i,j}} $ 为辅助变量，$ \rho \gt 0 $ 为罚系数。式（13）可采用共轭梯度算法求解，式（14）与式（15）具有解析解，即收缩阈值操作，表达式如下：
$$ s_{i,j}^q = {\text{shrinakge}}\left(\frac{\left(u_i^q - u_j^q - {\boldsymbol{l}}_{ij}^{\rm T}{\boldsymbol{g}}_i^q\right) - b_{i,j}^{q - 1}} {{\alpha _1}{w_1}(i,j)/\rho} \right)\text{，} $$ (18) $$ {\boldsymbol{r}}_{i,j}^q = {\text{shrinakge}}\left(\frac{\bigg(\left({\boldsymbol{g}}_i^q - {\boldsymbol{g}}_j^q\right) - {\boldsymbol{d}}_{i,j}^{q - 1}} {{\alpha _2}{w_1}(i,j)/\rho} \right) \text{，} $$ (19) 其中，shrinage(x, t)=sign(x)·max(|x|−t, 0)，t为阈值，sign(·)为符号函数，max(·)表述取两个变量中的最大值。
综上所述，基于NLTGV的迭代CT重建算法流程总结如表2：
表 2 基于NLTGV正则化的稀疏角度CT图像重建算法步骤Table 2. Procedures of a sparse-view CT image reconstruction algorithm based on NLTGV regularization算法1 基于NLTGV的稀疏角度CT重建算法 （1） 初始化: 给定初值$ {{{\boldsymbol u}}^0}， $$ {t_{\max }} ，$令初始化: 给定初值
$ {{{\boldsymbol u}}^0} ，$设置$ {t_{\max }} $,$ {q_{\max }} $,$ {k_{\max }} $,令$ k = 1 $（2） while （不满足停止准则） do （3） 令$ {{{\boldsymbol f}}^0} = {{{\boldsymbol u}}^{k - 1}} $，t=1 （4） while（t$ \leq $tmax） do （5） for j=1,2,···,N （6） $ f_j^t = f_j^{t - 1} + \displaystyle\frac{\lambda }{{{A_{ + ,j}}}}\sum\limits_{i = 1}^M {\frac{{{a_{i,j}}}}{{{A_{i, + }}}}} （{p_i} - {\bar p_i}） $ （7） $ {A_{i, + }} = \displaystyle\sum\limits_{j = 1}^N {{a_{i,j}}} ,\ i = 1,2,\cdots,M $ （8） $ {A_{ + ,j}} = \displaystyle\sum\limits_{i = 1}^M {{a_{i,j}}} ,\ j = 1,2,\cdots,N $ （9）$ {\bar p_i} = {A_i}{f^{t - 1}} $ （10）$ t = t + 1 $ （11） end for loop （12） end while loop （13） 非负约束，得到$ {{{\boldsymbol u}}^{{\text{pos}}}} $ （14） 令$ q = 1, $ b=0, d=0 （15） while（$ q \le {q^{\max }} $） do （16） 共轭梯度法求解u, g （17） 根据公式（18）与（19）求解s与r子问题 （18） 根据公式（16）与（17）更新变量b与d （19） $ q = q + 1 $ （20） end while loop （21） $ {{{\boldsymbol u}}^k} = {{{\boldsymbol u}}^{{q_{\max }}}} $ （22） $ k = k + 1 $ （23） end while loop （24） 输出: 重建图像$ {{{\boldsymbol u}}^*} $ 算法1的停止准则设置为相邻两次迭代的重建图像
$ {\boldsymbol{u}} $ 变化足够小或者迭代次数k达到最大迭代次数kmax。3. 数值仿真和实际实验
为了评估所提NLTGV重建算法的性能，本节采用模拟数据和实际投影数据进行图像重建，并与基于TV[3]、TGV[13]以及 NLTV[16]的重建算法进行比较。
3.1 定量评价指标
本文采用峰值信噪比（peak signal to noise ratio，PSNR）评估重建图像的准确性，结构相似性指数（structure similarity index measure，SSIM）[23]评估重构图像与参考图像之间的结构相似性，均方误差（root mean square error，RMSE）验证算法的收敛性。PSNR、SSIM及RMSE的定义：
$$ {\text{PSNR}} =10{\text{lo}}{{\text{g}}_{10}}\left(\frac{{{{(\max ({{\boldsymbol{u}}_{{\rm{ref}}}}))}^2}}}{{\displaystyle \sum\limits_{{i_1} = 1}^{{N_1}} {\sum\limits_{{i_2} = 1}^{{N_2}} {\frac{1}{{{N_1}{N_2}}}{{\Big({{\boldsymbol{u}}_{\rm{ref}}}({i_1},{i_2}) - {{\boldsymbol{u}}_{\rm{rec}}}({i_1},{i_2})\Big)}^2}} } }}\right) \text{，}$$ (20) $$ {\text{SSIM}} = \frac{{2{{\bar u}_1}{{\bar u}_2}\Big(2{\sigma _{12}} + {c_2}\Big)}}{{\Big(\bar u_1^2 + \bar u_2^2 + {c_1}\Big)\Big(\sigma _1^2 + \sigma _2^2 + {c_2}\Big)}} \text{，} $$ (21) $$ {\text{RMSE}} = \left( {\frac{1}{{{N_1}{N_2}}}{{\left\| {{{\boldsymbol{u}}_{ref}} - {{\boldsymbol{u}}_{rec}}} \right\|}^2}}\right)^{\tfrac{1}{2}} \text{，} $$ (22) 其中，
$ {{\boldsymbol{u}}_{\rm{ref}}} $ 和$ {{\boldsymbol{u}}_{\rm{rec}}} $ 分别是参考图像和重建图像。$ {\bar u_1}、 {\bar u_2}、{\sigma _1}、{\sigma _2}、{\sigma _{12}} $ 分别是$ {{\boldsymbol{u}}_{\rm{ref}}} $ 和$ {{\boldsymbol{u}}_{\rm{rec}}} $ 的均值、方差和协方差。正常数$ {c_1}、{c_2} $ 用于确保SSIM的计算数值稳定性。3.2 实验设置
在扇形束扫描模式下CT的相关成像参数：①X射线源到探测器的距离为
1400 mm；②X射线源到旋转轴的距离为1000 mm；③探测器数量为1024 ；④重建图像大小为 256像素×256像素，重建像素尺寸为0.5 mm$ \times $ 0.5 mm。本文实验均在Matlab 2014a环境下完成，计算机操作系统为Window 10，CPU为Intel i7-7700，内存为32 GB。度量两个像素的相似性图像块大小设置为3×3，搜索相似图像块的窗口设置为7×7，以平衡计算量和重建性能。算法的最大迭代次数kmax设置为20次，SART算法最大迭代次数
$ {t_{\max }} $ 设置为200次，Bregman算法最大迭代次数$ {q_{\max }} $ 设置为50次。权重$ {w_1}、{w_2} $ 中参数$ {h_l} $ 控制权重$ {w_1} $ 受两个像素空间距离的影响程度，可设置为半个搜索窗口的大小[19]，即设置为4，$ {h_u} $ 通过实验确定为0.5效果较好。为公平比较，在仿真和实际实验中各对比方法的正则化参数都进行精调，使各算法重建图像的RMSE都达到最小值。3.3 仿真数据重建
在仿真实验中，重建使用的参考图像如图2（a）所示，以等角度间隔在
$ [{0^\circ },{360^\circ }) $ 内采样90个角度投影数据。泊松高斯混合噪声模型[24]用于生成含噪声的投影数据，则沿着第i条X射线传输路径的投影测量值$ {I_i} $ 可以通过下式计算：$$ {I_i} = {\text{Possion}}\Big({I_0} e^{ - {A_i}u}\Big) + N\Big(0,{\sigma ^2}\Big),\;\;i = 1,2,\cdots,M \text{，} $$ (23) 其中，
$ {I_0} $ 为每条X射线的初始强度，$ {\sigma ^2} $ 为高斯噪声方差。$ \boldsymbol{A}_i $ 是投影矩阵A的第i行。在仿真实验中，$ {I_0} =50\ 000 $ 和$ {\sigma ^2} = 10 $ 。图2给出了各对比算法的重建结果。为了清楚地观察不同重建结果的细节，图3给出了图2中选取的感兴趣区域（region of interest，ROI）的放大图。从图2可看出，TV方法重建的图像中有明显的阶梯效应，图3中对应的ROI出现明显的块状伪影，如箭头所示。从图2中可观察到，TGV方法可以在一定程度上缓解块效状伪影，但仍存在明显的过平滑问题，相应的图3中ROI表明其重建图像的细节模糊较为严重。
![]()
图 3 仿真数据重建图像ROI放大图Figure 3. Zoom-in display of ROI of the images reconstructed from the simulated data如图3中NLTV重建图像的ROI放大图所示，NLTV方法重建图像的结构特征得到了一定恢复，但边界呈现锯齿状。这是因为NLTV方法虽然能利用非局部自相似先验信息，但是类似TV方法，只利用了一阶变分信息。与其他方法相比，所提出的NLTGV方法可以更好地保护图像中的结构特征，图像的结构特征更加清晰，更加接近参考图像，如图3中NLTGV重建图像的ROI所示。这可以归因于NLTGV方法不仅可以利用高阶变分对图像的过渡部分进行建模，而且还可以利用非局部先验信息来保留重建图像的细节特征。
为了客观评估不同方法的重建结果，表3列出各方法重建图像在ROI区域的PSNR和SSIM数值，其表明NLTGV方法在定量评估方面优于其他方法。
表 3 仿真数据重建图像的PSNR和SSIM结果Table 3. PSNR and SSIM for reconstructing images from simulated data采样角度数目 指标 算法 TV TGV NLTV NLTGV ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 50 PSNR 26.16 25.98 26.19 26.12 26.15 26.18 26.59 26.83 SSIM×10-2 80.02 79.97 80.37 80.52 78.73 79.71 80.48 81.78 70 PSNR 26.80 26.10 26.90 25.84 27.38 26.98 28.28 27.82 SSIM×10-2 81.96 80.90 79.78 80.16 83.05 83.53 85.10 84.93 90 PSNR 27.07 26.20 27.19 26.10 27.63 26.51 28.51 28.56 SSIM×10-2 81.64 82.95 82.37 81.80 82.62 83.34 85.11 87.29 110 PSNR 27.07 27.11 27.77 27.00 27.86 27.59 28.71 28.78 SSIM×10-2 84.24 84.74 82.71 84.03 83.20 84.72 85.47 87.36 不同采样角度数目下重建图像的RMSE随着迭代次数的变化情况如图4所示，其中横轴表示迭代次数，纵轴表示相应的RMSE值。从图中曲线可看出，在各采样角度下RSME数值均逐渐减小并趋于稳定，从而验证了算法的收敛性。
![]()
图 4 仿真数据重建图像的RMSE随迭代次数变化曲线Figure 4. RMSE of the image reconstructed from simulated data versus iterations3.4 实际数据重建
在真实数据研究中，本文选择一核桃样本来重建其断层图像。X射线管电压和电流分别设置为90 kVp和0.5 mA，探测器曝光时间设置为500 ms，在
$ [{0^\circ },{360^\circ }) $ 内等角度间隔采样90个投影数据。参考图像来自高剂量下全角度范围内360个等角度间隔投影数据的重建结果（X射线管电压90 kVp，电流1.2 mA，探测器曝光时间500 ms）。重建图像的像素大小为0.5 mm$ \times $ 0.5 mm，其他成像参数同模拟实验。实际投影数据重建图像如图5所示，表明不同算法抑制噪声和保留细节结构的能力不同。图5（a）显示了FBP算法利用上述提到的高剂量投影数据重建的参考图像。如图6中ROI1和ROI2的（b1）、（c1）、（d1）、（e1）所示，TV方法会导致重建图像产生块状伪影，并模糊细小的结构特征。这可解释为TV正则化项基于分段常量假设，因此在消除噪声的同时易产生阶梯效应。TGV方法可以有效地避免块状伪影，但它仍然对细节过于平滑，图6（b）中的（b2）、（c2）、（d2）、（e2）尤其明显。这是由于TGV正则化方法只使用局限于像素邻域的先验信息，因此在消除噪声的同时不可避免地对细节进行过度平滑。
![]()
图 6 实际数据重建图像ROI放大图Figure 6. Zoom-in display of ROI of the images reconstructed from the real data相比上述两种方法，NLTV正则化方法使用了空间域中的非局部自相似信息，可以有效去除噪声并保留重建图像的细节，如图5中NLTV重建图像所示。然而，图6中相应的ROI表明，NLTV方法仍然导致过渡区域产生了阶梯效应。图5中NLTGV重建图像与图6中对应的ROI表明NLTGV方法不仅可以有效抑制噪声与伪影，并且可更好地保护重建图像中平滑变化的特征，保留了更多的细节，进一步提高重建图像的对比度，在视觉上优于其他方法。
表4给出了不同采样角度下各算法重建图像被选取的ROI定量评估结果，其表明本文方法在定量评估方面优于其他对比算法，即重建精度更高。
表 4 真实数据重建图像的PSNR和SSIM结果Table 4. PSNR and SSIM for reconstructing images from real data采样角度数目 指标 算法 TV TGV NLTV NLTGV ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 50 PSNR 25.38 26.17 25.85 26.45 26.18 26.77 26.83 26.98 SSIM×10-2 79.97 79.98 79.52 81.12 79.71 82.74 81.78 82.70 70 PSNR 26.10 27.02 26.34 27.02 26.98 27.39 27.82 28.94 SSIM×10-2 80.90 80.00 80.16 83.67 83.53 86.91 84.93 88.18 90 PSNR 26.20 27.16 26.30 27.32 26.51 27.82 28.36 29.21 SSIM×10-2 81.95 84.68 82.80 84.27 83.34 84.21 87.29 89.04 110 PSNR 27.11 27.51 27.00 27.54 27.19 28.64 28.78 29.67 SSIM×10-2 83.74 84.42 84.03 84.67 84.72 86.22 87.36 89.56 为了进一步验证NLTGV算法的收敛性，图7显示了各采样角度下重建图像的RMSE随迭代次数变化情况，可观察到所提算法最终达到收敛状态。
![]()
图 7 实际数据重建图像的RMSE随迭代次数变化曲线Figure 7. RMSE of the image reconstructed from real data versus iterations4. 结论
本文提出一种NLTGV正则化方法来重建稀疏角度CT图像，其优势在于能同时利用高阶变分的稀疏性先验与非局部自相似性先验，有利于保护重建图像的结构特征。由于重建模型含有多个非光滑正则项，本文提出了一种基于凸集投影优化框架的交替求解策略。仿真和实际数据实验验证了算法的正确性与有效性。此外，本文NLTGV正则化方法很容易扩展到其他CT重建任务中，如内部重建、多能谱CT重建等。
与TGV等局部变分正则化方法相比，本文方法涉及到各像素点的非邻域变分运算，因而计算量更大，下一步工作重点将利用图形处理单元（GPU）并行计算技术加速重建过程。为了进一步提升重建效果，将Lp（
$ 0 \lt p \lt 1 $ ）范数与所提模型相结合实现非凸压缩感知重建也是未来的研究重点。 -
![]()
图 1 搜索非局部相似图像块原理
Figure 1. Schematic illustrating the searching of non-local similarity patches
![]()
图 3 仿真数据重建图像ROI放大图
Figure 3. Zoom-in display of ROI of the images reconstructed from the simulated data
![]()
图 4 仿真数据重建图像的RMSE随迭代次数变化曲线
Figure 4. RMSE of the image reconstructed from simulated data versus iterations
![]()
图 6 实际数据重建图像ROI放大图
Figure 6. Zoom-in display of ROI of the images reconstructed from the real data
![]()
图 7 实际数据重建图像的RMSE随迭代次数变化曲线
Figure 7. RMSE of the image reconstructed from real data versus iterations
表 1 不同变分正则项特性的比较
Table 1 Comparison of characteristics with different variation regularization
算法 TV TGV NLTV NLTGV 邻域 局部 局部 非局部 非局部 变分阶次 一阶 高阶 一阶 高阶 
下载: 导出CSV
表 2 基于NLTGV正则化的稀疏角度CT图像重建算法步骤
Table 2 Procedures of a sparse-view CT image reconstruction algorithm based on NLTGV regularization
算法1 基于NLTGV的稀疏角度CT重建算法 （1） 初始化: 给定初值$ {{{\boldsymbol u}}^0}， $$ {t_{\max }} ，$令初始化: 给定初值
$ {{{\boldsymbol u}}^0} ，$设置$ {t_{\max }} $,$ {q_{\max }} $,$ {k_{\max }} $,令$ k = 1 $（2） while （不满足停止准则） do （3） 令$ {{{\boldsymbol f}}^0} = {{{\boldsymbol u}}^{k - 1}} $，t=1 （4） while（t$ \leq $tmax） do （5） for j=1,2,···,N （6） $ f_j^t = f_j^{t - 1} + \displaystyle\frac{\lambda }{{{A_{ + ,j}}}}\sum\limits_{i = 1}^M {\frac{{{a_{i,j}}}}{{{A_{i, + }}}}} （{p_i} - {\bar p_i}） $ （7） $ {A_{i, + }} = \displaystyle\sum\limits_{j = 1}^N {{a_{i,j}}} ,\ i = 1,2,\cdots,M $ （8） $ {A_{ + ,j}} = \displaystyle\sum\limits_{i = 1}^M {{a_{i,j}}} ,\ j = 1,2,\cdots,N $ （9）$ {\bar p_i} = {A_i}{f^{t - 1}} $ （10）$ t = t + 1 $ （11） end for loop （12） end while loop （13） 非负约束，得到$ {{{\boldsymbol u}}^{{\text{pos}}}} $ （14） 令$ q = 1, $ b=0, d=0 （15） while（$ q \le {q^{\max }} $） do （16） 共轭梯度法求解u, g （17） 根据公式（18）与（19）求解s与r子问题 （18） 根据公式（16）与（17）更新变量b与d （19） $ q = q + 1 $ （20） end while loop （21） $ {{{\boldsymbol u}}^k} = {{{\boldsymbol u}}^{{q_{\max }}}} $ （22） $ k = k + 1 $ （23） end while loop （24） 输出: 重建图像$ {{{\boldsymbol u}}^*} $
下载: 导出CSV
表 3 仿真数据重建图像的PSNR和SSIM结果
Table 3 PSNR and SSIM for reconstructing images from simulated data
采样角度数目 指标 算法 TV TGV NLTV NLTGV ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 50 PSNR 26.16 25.98 26.19 26.12 26.15 26.18 26.59 26.83 SSIM×10-2 80.02 79.97 80.37 80.52 78.73 79.71 80.48 81.78 70 PSNR 26.80 26.10 26.90 25.84 27.38 26.98 28.28 27.82 SSIM×10-2 81.96 80.90 79.78 80.16 83.05 83.53 85.10 84.93 90 PSNR 27.07 26.20 27.19 26.10 27.63 26.51 28.51 28.56 SSIM×10-2 81.64 82.95 82.37 81.80 82.62 83.34 85.11 87.29 110 PSNR 27.07 27.11 27.77 27.00 27.86 27.59 28.71 28.78 SSIM×10-2 84.24 84.74 82.71 84.03 83.20 84.72 85.47 87.36
下载: 导出CSV
表 4 真实数据重建图像的PSNR和SSIM结果
Table 4 PSNR and SSIM for reconstructing images from real data
采样角度数目 指标 算法 TV TGV NLTV NLTGV ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 ROI1 ROI2 50 PSNR 25.38 26.17 25.85 26.45 26.18 26.77 26.83 26.98 SSIM×10-2 79.97 79.98 79.52 81.12 79.71 82.74 81.78 82.70 70 PSNR 26.10 27.02 26.34 27.02 26.98 27.39 27.82 28.94 SSIM×10-2 80.90 80.00 80.16 83.67 83.53 86.91 84.93 88.18 90 PSNR 26.20 27.16 26.30 27.32 26.51 27.82 28.36 29.21 SSIM×10-2 81.95 84.68 82.80 84.27 83.34 84.21 87.29 89.04 110 PSNR 27.11 27.51 27.00 27.54 27.19 28.64 28.78 29.67 SSIM×10-2 83.74 84.42 84.03 84.67 84.72 86.22 87.36 89.56
下载: 导出CSV
-
[1] PARCERO E, FLORES L, SÁNCHEZ M G, et al. Impact of view reduction in CT on radiation dose for patients[J]. Radiation Physics and Chemistry, 2017, 137(8): 173-175.
[2] CANDES E J, 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. DOI: 10.1109/TIT.2005.862083.
[3] SIDKY E Y, KAO C M, PAN X. Accurate image reconstruction from few-views and limited-angle data in divergent-beam CT[J]. Journal of X Ray Science & Technology, 2009, 14(2): 119-139.
[4] 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.
[5] RUDIN L I, OSHER S, FATEMI E. Nonlinear total variation based noise removal algorithms[J]. Physica D: Nonlinear Phenomena, 1992, 60(1/4): 259-268.
[6] 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.
[7] LANGE K, CARSON R E. EM reconstruction algorithms for emission and transmission tomography[J]. Journal of Computer Assisted Tomography, 1984, 8(2): 306-316.
[8] BREDIES, K, KUNISCH, K, POCK, T. Total genera-lized variation[J]. SIAM Journal on Imaging Sciences, 2010, 3(3): 492-526. DOI: 10.1137/090769521.
[9] NIU S, GAO Y, BIAN Z, et al. Sparse-view X-ray CT reconstruction via total generalized variation regularization[J]. Physics in Medicine & Biology, 2014, 59(12): 2997-3017.
[10] SUN T, SUN N, WANG J, et al. Iterative CBCT reconstruction using Hessian penalty[J]. Physics in Medicine and Biology, 2015, 60(5): 1965-1987. DOI: 10.1088/0031-9155/60/5/1965.
[11] ZHANG H M, WANG L Y, YAN B, et al. Constrained total generalized p-variation minimization for few-view X-ray computed tomography image reconstruction[J]. Plos One, 2016, 11(2): e0149899.
[12] XI Y R, QIAO Z R, WANG W J, et al. Study of CT image reconstruction algorithm based on high order total variation[J]. Optik, 2020, 204: 163814. DOI: 10.1016/j.ijleo.2019.163814.
[13] 闫慧文, 乔志伟. 基于 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 algorithms 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).
[14] XI Y, ZHOU P, YU H, et al. Adaptive-weighted high order TV algorithm for sparse-view CT reconstruction[J]. Medical Physics, 2023, 50: 5568-5584. DOI: 10.1002/mp.16371.
[15] HUANG J, MA J, LIU N, et al. Sparse angular CT reconstruction using non-local means based iterative-correction POCS[J]. Computers in Biology and Medicine, 2011, 41(4): 195-205. DOI: 10.1016/j.compbiomed.2011.01.009.
[16] KIM H, CHEN J, WANG A, et al. Non-local total-variation (NLTV) minimization combined with reweighted L1-norm for compressed sensing CT reconstruction[J]. Physics in Medicine & Biology, 2016, 61(18): 6878.
[17] KIM K, EL F G, LI Q. Low-dose CT reconstruction using spatially encoded nonlocal penalty[J]. Physics in Medicine & Biology, 2018, 63(3): 035045.
[18] GILBOA G, OSHER S. Nonlocal operators with applications to image processing[J]. Multiscale Modeling & Simulation, 2009, 7(3): 1005-1028.
[19] RANFTL R, BREDIES K, POCK T. Non-local total generalized variation for optical flow estimation[C]//ECCV, 2014, (6): 439-454.
[20] TOMASI C, MANDUCHI R. Bilateral filtering for gray and color images[C]//ICCV, 1998: 839-846.
[21] JIANG M, WANG G. Convergence of the simultaneous algebraic reconstruction technique (SART)[J]. IEEE Transactions on Image Processing, 2003, 12(8): 957-961. DOI: 10.1109/TIP.2003.815295.
[22] GOLDSTEIN T, OSHER S. The split Bregman method for L1-regularized problems[J]. SIAM Journal on Imaging Sciences, 2009, 2(2): 323-343. DOI: 10.1137/080725891.
[23] WANG Z, BOVIK A C, SHEIKH H R, et al. Image quality assessment: From error visibility to structural similarity[J]. IEEE Transactions on Image Processing, 2004, 13(4): 600-612. DOI: 10.1109/TIP.2003.819861.
[24] WANG J, LI T, XING L. Iterative image reconstruction for CBCT using edge-preserving prior[J]. Medical Physics, 2009, 36(1): 252-260. DOI: 10.1118/1.3036112.
计量
- 文章访问数: 141
- HTML全文浏览量: 56
- PDF下载量: 49
