反投影重建技术以数学理论为依托，加之计算机技术的支持，在重建领域取得了长足的进步，在临床扫描、工业零件、木乃伊等多领域应用广泛^[1]，常用解析法重建图像，然而在稀疏采样时，解析法重建的图像会出现明显的条状伪影^[2]。随着压缩感知理论的提出，使得高精度稀疏重建成为可能，总变差（total variation，TV）最小化模型是目前重建问题中使用最广泛的正则化方法^[3]之一。Sidky等^[4]提出TV最小化算法，在仿真模体上可以从高度欠采样的投影数据中高精度地重建图像，可以有效减少辐射剂量，提高扫描速度，在稀疏重建和含噪投影重建中发挥了重要作用，并且该算法可适用于多种层析成像模式，如计算机断层成像（computed tomography，CT）、磁共振成像（magnetic resonance imaging，MRI）、电子顺磁共振成像（electron paramagnetic resonance imaging，EPRI）中。随后，出现了许多TV变体以及将标量值图像重建扩展到矢量值图像重建的算法，如：HOTV^[5-6]、TpV^[7-8]、NLTV^[9]、VTV^[10]等，在一定程度上提高了重建精度。

但在实际应用中，各种成像模态均有多通道任务。在CT中有多能谱CT、对比剂增强CT、不同层厚CT等不同类型的多通道CT，并且多通道之间所形成的图像具有很强的结构相似性；在MRI和EPRI中，有T1或T2不同的加权项，图像之间也具有结构相似性。那么，多通道图像是否可以用TV类算法实现高精度重建成为备受关注的问题。最直接的方式是将各个通道的TV值求和，也称为逐通道TV（channel-by-channel TV，TVs），由于TVs重建图像是单独地惩罚每个通道，并未耦合其他通道的图像信息，忽略了各通道图像之间的相似性。于是，Rigie等^[11-12]提出一种TV的多通道泛化方法——TnV模型并将其应用于多通道光谱CT图像重建，通过鼓励不同通道图像具有相同的边缘结构，并且梯度向量指向同一个方向，更好的保留了边界信息，但其推导公式过于复杂。

此外，通常求解TV类最小化模型使用ASD-POCS（adaptive steepest descent- projection onto convex sets）算法^[13-14]，但此算法存在很大的缺点——算法参数需要人为凭经验设定；Chambolle-Pock（CP）算法^[15]是一种优秀的凸优化问题求解算法，它以表达形式简单、算法参数不需要人为设定和每个子最优化问题均有闭合解的优点受到一些学者的青睐；Pan团队系统地总结了CP算法，并推导了与CT相关的多个优化问题的CP算法实例并给出伪代码^[16]；Qiao等^[17]提出了一种新型数据约束的平衡TV——bTV模型，并使用CP算法详细推导了bTV-CP算法实例，推导步骤简单明了，通过引入算法参数ν提高了收敛速率，并对算法的正确性及算法参数对收敛速率的影响进行了详细研究。

鉴于此，本文面向多通道图像重建，以CT重建为研究范例，提出一种CP算法框架下的核TV多通道图像重建算法。基于上述bTV模型的CP算法求解步骤，重新推导核TV模型中较为复杂的部分，并引入算法参数ν提高收敛速率。通过在仿真模体和真实CT模体上的实验，验证算法和计算机实现的正确性；分析算法的收敛性；评估算法参数对收敛速率的影响；通过改变稀疏采样角度的个数，评估核TV模型的稀疏重建能力；在含有不同强度噪声的投影数据下重建图像，评估核TV模型的去噪能力；并用多个图像质量衡量标准对TVs模型和核TV模型的重建结果进行定性和定量分析。

1.   相关知识

1.1   成像模型

目前所有的CT成像模型均是离散的，即离散—离散（discrete to discrete，D2D）的成像模型。可表示为：

多通道图像可表示为：

将D2D模型扩展到多通道：

其中，$ {\boldsymbol{g}}_{i} $是大小为M的列向量，表示第i个通道的投影数据；$ {\boldsymbol{u}}_{i} $是大小为N的列向量，表示第i个通道的图像；$ {\boldsymbol{A}}_{i} $是第i个通道的系统矩阵，大小为M×N，由于各个通道的成像条件一致，因此各个通道的系统矩阵是相同的；$ {{A}}_{i(m,n)} $为矩阵$ {\boldsymbol{A}}_{i} $中的元素，表示第i个通道中第$ n $个像素对第$ m $个投影数据测量值的贡献；$ L $为通道数。

该模型求解实质上是求线性方程组解的过程，但在实际情况中，该方程组会因规模大、扫描范围覆盖不足和噪声影响导致求解变得异常困难。为了解决上述问题，需进一步将其建模为一个最优化问题。

1.2   TVs最优化模型

为求解上述大规模病态方程组，可以引入稀疏先验等先验知识设计最优化模型。数据约束的TVs最小化模型可表示为：

其中，

$ {\boldsymbol{u}}^{*} $表示满足约束条件的图像，即最优化模型的解；$ {\parallel \boldsymbol{g}-\boldsymbol{A}\boldsymbol{u}\parallel }_{2} $表示生成投影数据与真实投影数据的$\ell_{2}$范数；$ \varepsilon $为数据容差，大小与噪声相关。经过分析可知，求$ \mathrm{T}\mathrm{V}\mathrm{s}\left(\boldsymbol{u}\right) $的最优解即求每个通道$ \mathrm{T}\mathrm{V}\left(\boldsymbol{u}\right) $的最优解，$ \mathrm{T}\mathrm{V}\left(\boldsymbol{u}\right) $表示图像$ \boldsymbol{u} $的TV范数，即图像梯度大小变化$ |\boldsymbol{D}\boldsymbol{u}{|}_{\mathrm{m}\mathrm{a}\mathrm{g}} $的$\ell_{1}$范数：

$ |\cdot {|}_{\mathrm{m}\mathrm{a}\mathrm{g}} $表示求模运算，本文中为$\ell_{2}$范数模，即：

$ \boldsymbol{D} $表示图像梯度变化，以大小为N×N的二维图像为例，$ \boldsymbol{D} $是一个大小为2N×N的矩阵：

$ {\boldsymbol{D}}_{1} $、$ {\boldsymbol{D}}_{2} $的大小均为N×N，分别表示图像沿$ x $轴和$ y $轴的梯度变化，设$ {u}_{s,t} $表示图像$ \boldsymbol{u} $中一个特定位置的像素点，$ s\in [1,N] $,$ t\in [1,N] $，则$ {\boldsymbol{D}}_{1} $、$ {\boldsymbol{D}}_{2} $可表示为：

那么，式（6）可表示为：

1.3   CP算法求解TVs模型

CP算法是一种原始—对偶算法，可以求解如（12）式的最优化模型：

其中，$ \boldsymbol{x} $，$ \boldsymbol{y} $分别是空间$ \boldsymbol{X} $和$ \boldsymbol{Y} $中的有限维向量；$ \boldsymbol{K} $为一个矩阵，表示$ \boldsymbol{X} $到$ \boldsymbol{Y} $的线性变化；$ F $和$ G $是凸函数，但不要求平滑。

则式（4）可表示为：

其中，$\delta_{\mathrm{b}\mathrm{a}\mathrm{l}\mathrm{l}\left(\varepsilon \right)}\left(\boldsymbol{x}\right)$为指示函数：

根据CP算法作以下变化，

那么，式（18）的最邻近映射为：

其中：$ {F}^{\mathrm{*}} $表示函数$ F $的凸共轭函数，$ \mathrm{p}\mathrm{r}\mathrm{o}{\mathrm{x}}_{\sigma }\left[F\right] $表示函数$ F $的最邻近映射；$ \mathrm{\sigma } $为非负常数；$ l $为通道序数。

2.   本文方法

2.1   核TV最优化模型

由于TVs模型重建多通道图像时没有耦合其他通道信息，并未有效利用通道间图像相似性的特点。因此，本文提出一种数据约束的核TV最小化模型，可以表示为：

其中：$ {\boldsymbol{u}}^{\mathrm{*}} $表示模型的最优解，即重建图像；$ \mathrm{T}\mathrm{n}\mathrm{V}\left(\boldsymbol{u}\right) $表示$ \boldsymbol{u} $的核TV范数，即$\mathrm{T}\mathrm{n}\mathrm{V}\left(\boldsymbol{u}\right)= \big\|\;\boldsymbol{u}\;\big\|_{\mathrm{T}\mathrm{n}\mathrm{V}}= \big\|\;\boldsymbol{J}\boldsymbol{u}\;\big\|_{1,\mathrm{*}}={\sum }_{i,j}\big\|\;\left(\boldsymbol{J}\boldsymbol{u}\right)(i,j)\;\big\|_{\mathrm{*}}$，对特定位置$ (i,j) $像素的核范数^[18-20]等于其雅克比矩阵的奇异值之和，即$\big\|\;\boldsymbol{J}\boldsymbol{u}(i,j)\;\big\|_{\mathrm{*}}=\big\|\;\overrightarrow{\boldsymbol{\sigma }}(i,j)\;\big\|_{1}$；${\big\| \boldsymbol{g}-\boldsymbol{A}\boldsymbol{u}\big\| }_{2}$为生成投影数据与实际投影数据的$\ell_{2}$范数；$ \varepsilon $为数据容差，与噪声大小有关。

特定位置$ (i,j) $的雅克比矩阵可表示为：

其中，

与此同时，引入平衡参数$ \mathrm{\nu } $以调整数据保真项与正则项的相对大小，在不影响模型解的情况下可以有效提高收敛速率，相关参数选取参考文献^[17,21-23]，得到新的TnV模型：

2.2   TnV-CP算法实例推导

式（25）可以表示为CP算法通用的求解形式：

其中：$ {\boldsymbol{u}}^{\mathrm{*}} $表示满足约束条件的图像，即被重建图像；$ \mathrm{\nu } $为非零常数；$\delta_{\mathrm{b}\mathrm{a}\mathrm{l}\mathrm{l}\left(\varepsilon \right)}(\boldsymbol{g}-\boldsymbol{A}\boldsymbol{u})$为指示函数，表示当向量$ (\boldsymbol{g}-\boldsymbol{A}\boldsymbol{u}) $在半径为$ \varepsilon $球的内部时，函数值为0，否则为$ \mathrm{\infty } $，可以表示为：

根据CP算法做以下变化：

其中，$ \boldsymbol{J} $为图像像素点的雅可比矩阵。由凸共轭函数和最邻近映射的定义可得$ {F}_{1}\left(\boldsymbol{p}\right) $和$ {F}_{2}\left(\boldsymbol{z}\right) $的凸共轭函数及$ {F}_{1}^{*} $、$ {F}_{2}^{*} $和$ G $的最邻近映射为：

其中，式（33）和式（34）是求取数据保真项和正则项的凸共轭函数。$ \varepsilon $为数据容差，与噪声大小有关；$ \overrightarrow{\boldsymbol{\sigma }}(i,j) $是$ (i,j) $位置像素雅克比矩阵$ \boldsymbol{z}(i,j) $的奇异值。将式（33）和式（34）代入最邻近映射函数中，可得式（30）～式（32）的最邻近映射为式（35）～式（37）。$ \sigma $为非零常数，$ \boldsymbol{z}(i,j) $的奇异值分解为$ \boldsymbol{U}\sum {\boldsymbol{V}}^{\mathrm{T}} ，$$ {\sum }_{\mathrm{P}}\mathrm{表} $示奇异值的软阈值缩放，即${\sum }_{\mathrm{P}}=\mathrm{m}\mathrm{i}\mathrm{n}\Big(1,|\sum |\Big)$。

将式（35）～式（37）代入CP算法框架中，可得如表1所示的TnV-CP算法实例的伪代码。

表  1  TnV-CP算法伪代码

Table  1.  The TnV-CP algorithm pseudocode

输入：	$ \boldsymbol{g},\boldsymbol{A},\boldsymbol{D},\epsilon ;\mathrm{\nu } $
1）	$L={\parallel \frac{\boldsymbol{A} }{\mathrm{\nu }\boldsymbol{D} }\parallel }_{\mathrm{S}\mathrm{V} },\;\sigma =\frac{1}{L},\;\tau =\frac{1}{L},\;\theta =1,\;n=0$
2）	${\boldsymbol{u} }_{0}=0,\;{\stackrel{-}{\boldsymbol{u} } }_{0}=0,\;{\boldsymbol{p} }_{0}=0$
3）	$ \mathrm{r}\mathrm{e}\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{t} $
4）	$ {\boldsymbol{p}}_{n+1}=\mathrm{max}(\left|\right|{\boldsymbol{p}}_{n}+\sigma (\boldsymbol{A}\overline{\boldsymbol{u}}-\boldsymbol{g})|{|}_{2}-\sigma \epsilon ,0)\cdot ({\boldsymbol{p}}_{n}+\sigma (\boldsymbol{A}\overline{\boldsymbol{u}}-\boldsymbol{g}\left)\right)/\left|\right|{\boldsymbol{p}}_{n}+\sigma (\boldsymbol{A}\overline{\boldsymbol{u}}-\boldsymbol{g})|{|}_{2} $
5）	$ {\mathbf{z}(\mathrm{i},\mathrm{j})}_{\mathrm{n}+1}=\boldsymbol{U}{\sum }_{\mathrm{P}}{\boldsymbol{V}}^{\mathrm{T}} $
6）	$\theta =\sqrt{\mathrm{\nu }/(\mathrm{\nu }+2\tau )},\;\tau =\theta *\tau ,\;\sigma =\sigma /\theta$
7）	$ {\boldsymbol{u}}_{n+1}={\boldsymbol{u}}_{n}-\tau {\boldsymbol{A}}^{\mathrm{T}}{\boldsymbol{p}}_{n+1}-\tau \mathrm{\upsilon }{\boldsymbol{D}}^{\mathrm{T}}{\boldsymbol{z}}_{n+1} $
8）	$ {\stackrel{-}{\boldsymbol{u}}}_{n+1}={\boldsymbol{u}}_{n+1}+\theta ({\boldsymbol{u}}_{n+1}-{\boldsymbol{u}}_{n}) $
9）	$ n=n+1 $
10）	$\mathrm{u}\mathrm{n}\mathrm{t}\mathrm{i}\mathrm{l}\;n\ge N$

下载: 导出CSV 

| 显示表格

其中，N为迭代次数；$ \parallel {\cdot \parallel }_{\mathrm{S}\mathrm{V}} $为矩阵的最大奇异值，具体求解步骤参考文献^[16]；$ \boldsymbol{A} $为系统矩阵，本文采用像素驱动法求取系统矩阵$ \boldsymbol{A} $^[24]；$ \boldsymbol{g} $为获取的投影数据；具体参数设置参考文献^{[17,21-23,25-26]}。

2.3   图像质量评价标准

图像质量最直接的评价方式是主观评价，即观测者肉眼观察重建图像与原图的区别，但此过程比较耗时，也需要观测者的耐心和丰富的经验。为了全面评估重建图像的质量，本文采用主客观结合的评价方式对重建图像的质量进行定性和定量分析，其中客观评价指标采用以下3种评价方法：

（1）图像误差即均方根误差（root mean square error，RMSE）

其中，$ {\boldsymbol{f}} $为重建图像；$ {\boldsymbol{u}} $为真实图像；$ {L}^{2} $为图像总像素个数。

（2）结构相似性（structural similarity，SSIM）

其中，${{\mu }_{f}，\mu }_{u}，{\sigma }_{fu}{，\mu }_{f}^{2}，{\mu }_{u}^{2}$分别为图像$ {\boldsymbol{f}} $和$ {\boldsymbol{u}} $的平均值、协方差和方差；为了避免分式分母为0，$ {C}_{1} $和$ {C}_{2} $为常数。

（3）峰值信噪比（peak signal to noise ratio，PSNR）

其中，$ n $为每个采样点的比特数；$ ({2}^{n}-1) $为信号的最大值；$ \mathrm{M}\mathrm{S}\mathrm{E} $为均方误差。

3.   实验结果与分析

本节包括5个研究点：①算法正确性验证——验证 TnV模型及求解过程的正确性；②使用仿真模体分析 TnV-CP算法的收敛性；③在仿真模体及真实CT模体上分析算法参数对收敛速率的影响；④在不同稀疏角度下重建图像，评估TnV模型的稀疏重建能力；⑤在固定角度不同强度噪声投影下重建图像，评估TnV模型的去噪能力。

本文所有实验均使用3个通道，第1通道的图像为模体的真值图像，第2和第3通道的图像为真值图像经过非线性变化后的图像（图1）。为了方便起见，后续所有实验结果均只展示第1通道的相关图像和实验数据。

图  1  三种模体在3个通道中的图像

Figure  1.  Images of three motifs in three channels

下载: 全尺寸图片 幻灯片

3.1   算法正确性验证——验证TnV模型及求解过程的正确性

算法正确性验证是在投影数据完备的前提下，如果重建图像与真值图像的误差充分小，说明设计的重建模型、求解算法及计算机实现均是正确的。在本节中，使用Shepp-Logan和FORBILD仿真模体对TnV模型和求解过程进行正确性验证。模体大小均为256×256，旋转中心位于图像中心[128,128]的位置，探元长度为1，探测器中探元的个数与模体的长度保持一致，在[0,π]范围内均匀采集360个角度的投影数据进行重建。$ \mathrm{\nu } $的取值为0.1，其余参数均通过计算得到，并将收敛条件设置为图像误差RMSE$ \le $10^-4。

图2是TnV模型在Shepp-Logan模体的重建结果图。图2（a）和图2（b）分别是真值图像和重建图像，两幅图像几乎一致，用肉眼无法区分。图2（c）是真值图像与重建图像的水平中心线波形比较，可见两条水平中心线波形基本重合，在主观评价上，TnV模型取得了高精度重建。

图  2  TnV模型在Shepp-Logan模体的重建结果

Figure  2.  Reconstruction results of TnV model in Shepp-Logan phantom

下载: 全尺寸图片 幻灯片

图3是TnV模型在FORBILD模体的重建结果图。图3（a）和图3（b）分别是真值图像和重建图像，两幅图像基本一致，用肉眼难以区分。图3（c）是真值图像与重建图像的水平中心线波形比较，可见两条水平中心线波形基本重合，在主观评价上，TnV模型取得了很高的重建精度。

图  3  TnV模型在FORBILD模体的重建结果

Figure  3.  Reconstruction results of TnV model in FORBILD phantom

下载: 全尺寸图片 幻灯片

客观评价采用RMSE评价标准对重建图像质量进行评估。理论上，当RMSE的数值小于10^-4时，显示器便已无法区分重建图像与真值图像，可认为算法正确并收敛。图4（a）和图4（b）分别为TnV模型在Shepp-Logan和FORBILD模体重建过程中RMSE走势图，可以看出两种模体重建的误差均小于10^-4，如果继续迭代至收敛状态，RMSE可分别达到10^-7和10^-6。两种模体的RMSE走势图均说明TnV模型及求解过程是正确的，并且取得了高精度重建。

图  4  TnV模型重建过程中图像误差趋势图

Figure  4.  Trend chart of image error during TnV model reconstruction

下载: 全尺寸图片 幻灯片

3.2   TnV-CP算法收敛性分析

为了更直观的分析TnV-CP算法的收敛性，引入下列两个衡量指标：

其中，${{M}}_{1}$为重建图像生成投影与原始投影数据的$\ell_{1}$范数距离，也称为数据残差；${{M}}_{2}$为重建图像TV值的相对误差。

图5（a）和图5（c）分别为TnV算法在两种模体的数据残差走势图，可见两幅图的数据残差呈现明显的下降趋势，直到收敛状态。图5（b）和图5（d）分别为TnV算法在两种模体TV相对误差的迭代走势图，两幅图的TV相对误差均小于10^-4。由于TV值是图像所有像素点的梯度大小变化之和，将TV相对误差分配到每个像素点上，每个像素点的相对误差小于10^-7，显示器已无法区分图像的变化，但此时，TV相对误差还在持续下降至收敛状态。

图  5  TnV-CP算法在两种仿真模体的收敛性

Figure  5.  Convergence of TnV-CP algorithm in two simulation models

下载: 全尺寸图片 幻灯片

从图4和图5的3种迭代走势图分析可知，3种图像质量衡量指标均已达到所允许的收敛状态。

3.3   算法参数对收敛速率的影响

为了探究算法参数$ \mathrm{\nu } $对收敛速率的影响，本节在两种仿真模体和真实CT图像模体上使用5种不同的参数组合，在10个稀疏角度下观察不同参数组合重建图像的RMSE下降趋势，比较算法参数对收敛速率的影响（图6）。图6（a）～图6（c）分别为不同参数组合在仿真模体和真实CT模体上对收敛速率的影响曲线图，可见3种模体的收敛速率不完全相同，参数选择依赖具体的重建对象，但对参数$ \mathrm{\nu }=0.1 $，3种模体的收敛速率与重建精度都是最佳的。因此，在后续实验中，算法参数均取$ \mathrm{\nu }=0.1 $。

图  6  不同算法参数对收敛速率的影响

Figure  6.  Influence of different algorithm parameters on convergence rate

下载: 全尺寸图片 幻灯片

3.4   稀疏重建能力评估

为了评估TnV模型的稀疏重建能力，本节采用3.3节中的3种模体在10、20、30、40和50个稀疏角度下分别使用TVs和TnV两种算法重建各通道图像，从主观和客观两个方面评估TnV算法的稀疏重建能力。在实验中，两种算法参数均保持一致。图7、图9和图11分别展示了在仿真模体和真实CT模体上两种算法在不同稀疏角度下均达到收敛状态时的重建图像。

图  7  两种算法在Shepp-Logan模体上稀疏重建结果比较

Figure  7.  Comparison of sparse reconstruction results between two algorithms on Shepp-Logan phantoms

下载: 全尺寸图片 幻灯片

图  9  两种算法在FORBILD模体上稀疏重建结果比较

Figure  9.  Comparison of sparse reconstruction results between two algorithms on FORBILD phantoms

下载: 全尺寸图片 幻灯片

图  11  两种算法在真实CT模体上稀疏重建结果比较

Figure  11.  Comparison of sparse reconstruction results between two algorithms on real CT phantoms

下载: 全尺寸图片 幻灯片

图7为两种算法在Shepp-Logan模体的重建结果图，可以看出两种算法在20个稀疏角度下重建图像与真值图像已无明显区别，稀疏重建能力强。

为了更直观的观察两种算法的重建结果，图8为在10个稀疏角度重建结果的放大图和水平中心线波形图。图8（a）和图8（b）为两种算法重建图像感兴趣区域的放大图，可以看出TVs算法重建局部图在边缘交汇处十分模糊，难以区分不同的边缘，而TnV算法重建局部图可以明显看出两个结构的交汇部分及单个结构的主要边缘，表明TnV算法的边缘保持能力强于TVs算法。图8（c）和图8（d）为两种算法重建图像水平中心线波形图，图8（d）中箭头所指的重建波形更接近真值波形，说明在采样角度极其稀疏的情况下，TnV算法的重建能力强于TVs算法。

图  8  两种算法在10个稀疏角度下重建局部图

Figure  8.  Two algorithms for local graph reconstruction under 10 sparse angles

下载: 全尺寸图片 幻灯片

图9为两种算法在FORBILD模体的重建结果图，两种算法在30个角度下重建图像与真值图像在肉眼上难以区分，已实现高精度重建。在10和20个稀疏角度下，用肉眼可以观察到两幅图像有明显区别，TnV算法重建的边缘模糊部分更少，边缘保持能力更强。为了更直观的观察两种算法的重建结果，图10给出了两种算法在20个角度重建图像的水平中心线波形图，可以看出TnV算法重建图像的中心线波形更接近真值波形，重建精度更高。

图  10  20个稀疏角度下重建图像的水平中心线波形图

Figure  10.  Horizontal centerline waveform of reconstructed images unader 20 sparse angles

下载: 全尺寸图片 幻灯片

图11为两种算法在真实CT模体的重建结果图。在30～50个稀疏角度下，两种算法的重建结果在肉眼观察上几乎一致。在10和20个角度重建结果的红色框局部区域，TVs算法并未重建出真实模体的细小结构，而TnV算法可以简单重建出微小的结构，重建精度更高。

客观评价使用RMSE和SSIM两个评价指标对重建图像进行定量分析。表2～表4分别给出两种算法在3种模体上均达到收敛状态时的RMSE和SSIM值。结果表明，不论在哪种模体上，TnV算法的重建精度均高于TVs算法。

表  2  两种算法重建Shepp-Logan模体在不同稀疏角度下评估参数的比较

Table  2.  Comparison of evaluation parameters of Shepp-Logan phantom reconstructed by two algorithms under different sparsity angles

评估参数	算法	投影个数
		10	20	30	40	50
RMSE	TVs	0.079	0.002	6.626×10-7	4.881×10-7	4.764×10-7
	TnV	0.067	1.214×10-6	6.584×10-7	4.651×10-7	4.585×10-7
SSIM	TVs	0.798	0.999	1.000	1.000	1.000
	TnV	0.831	1.000	1.000	1.000	1.000

下载: 导出CSV 

| 显示表格

表  3  两种算法重建FORBILD模体在不同稀疏角度下评估参数的比较

Table  3.  Comparison of evaluation parameters between two algorithms in reconstructing FORBILD motifs under different sparse angles

评估参数	算法	投影个数
		10	20	30	40	50
RMSE	TVs	0.101	0.022	2.022×10-6	1.311×10-6	1.129×10-6
	TnV	0.090	0.014	1.484×10-6	1.256×10-6	1.034×10-6
SSIM	TVs	0.700	0.981	1.000	1.000	1.000
	TnV	0.736	0.992	1.000	1.000	1.000

下载: 导出CSV 

| 显示表格

表  4  两种算法重建真实CT模体在不同稀疏角度下评估参数的比较

Table  4.  Comparison of evaluation parameters between two algorithms in reconstructing real CT phantom under different sparsity angles

评估参数	算法	投影个数
		10	20	30	40	50
RMSE	TVs	0.104	0.072	0.054	0.043	0.035
	TnV	0.099	0.069	0.052	0.041	0.034
SSIM	TVs	0.593	0.748	0.840	0.894	0.923
	TnV	0.617	0.764	0.850	0.901	0.927

下载: 导出CSV 

| 显示表格

3.5   去噪能力评估

为了评估TnV模型的去噪能力，本节在60个固定投影角度下，在3个通道的投影数据中依次加入不同强度的噪声，通过比较两种算法的重建结果，定性定量的分析TnV模型的去噪能力。其余算法参数设置与3.4节保持一致。由于在其他强度噪声条件下，两种算法重建图像的规律是一致的，图12仅展示第1通道3种模体在30 dB噪声下的重建结果。

图  12  两种算法在不同模体上含噪重建结果比较

Figure  12.  Comparison of noisy reconstruction results between two algorithms on different phantoms

下载: 全尺寸图片 幻灯片

从图12可直观的看出TnV算法重建的模糊部分更少，边缘保持能力和去噪能力更强。表5、表6和表7详细给出了3种模体在不同强度噪声下重建图像的质量评价指标值，均表明TnV模型较TVs模型有着更高的重建精度，去噪能力更强。

表  5  两种算法重建Shepp-Logan模体在不同噪声强度下评估参数的比较

Table  5.  Comparison of evaluation parameters of Shepp Logan phantom reconstructed by two algorithms under different noise intensity

评估参数	算法	噪声等级
		30 dB	35 dB	40 dB	45 dB	50 dB
RMSE	TVs	0.035	0.025	0.018	0.012	0.008
	TnV	0.032	0.023	0.017	0.011	0.007
SSIM	TVs	0.958	0.980	0.989	0.994	0.997
	TnV	0.967	0.982	0.990	0.995	0.998
PSNR	TVs	29.090	31.910	34.840	38.460	42.270
	TnV	29.930	32.590	35.470	39.110	42.920

下载: 导出CSV 

| 显示表格

表  6  两种算法重建FORBILD模体在不同噪声强度下评估参数的比较

Table  6.  Comparison of evaluation parameters between two algorithms in reconstructing FORBILD phantom under different noise intensity

评估参数	算法	噪声等级
		30 dB	35 dB	40 dB	45 dB	50 dB
RMSE	TVs	0.073	0.060	0.043	0.030	0.019
	TnV	0.069	0.056	0.038	0.026	0.017
SSIM	TVs	0.875	0.913	0.947	0.974	0.988
	TnV	0.878	0.916	0.954	0.978	0.990
PSNR	TVs	22.700	24.390	27.330	30.570	34.230
	TnV	23.260	25.080	28.400	31.610	35.160

下载: 导出CSV 

| 显示表格

表  7  两种算法重建真实CT模体在不同噪声强度下评估参数的比较

Table  7.  Comparison of evaluation parameters between two algorithms in reconstructing real CT phantom under different noise intensity

评估参数	算法	噪声等级
		30 dB	35 dB	40 dB	45 dB	50 dB
RMSE	TVs	0.056	0.046	0.039	0.034	0.031
	TnV	0.052	0.044	0.037	0.033	0.030
SSIM	TVs	0.842	0.884	0.915	0.935	0.943
	TnV	0.854	0.896	0.923	0.938	0.947
PSNR	TVs	25.010	26.760	28.280	29.510	30.260
	TnV	25.610	27.100	28.630	29.690	30.520

下载: 导出CSV 

| 显示表格

4.   结语

本文面向多通道图像重建，以CT重建为研究范例，提出一种适用于多通道图像重建的核TV模型，并使用CP算法求解。在仿真模体上验证了算法及计算机实现的正确性，分析算法的收敛性，并在仿真模体与真实CT模体上探究算法参数对收敛速率的影响，评估了算法的稀疏重建能力与含噪重建能力。实验结果表明该算法能够考虑通道间图像的相似性，可以实现高精度图像重建，是一种高精度的多通道图像重建算法，并且可应用于各种成像模态的多通道重建场合。

本文提出的模型虽然能够利用通道间图像相似性的特征有效提高了重建精度，但在执行核范数最小化的过程中对不同奇异值有着相同的软阈值收缩过程，而不同奇异值具有不同的价值，较大的奇异值代表越多的结构信息，在重建过程中，对不同奇异值区别处理，调整相应的收缩幅度尽可能的保留更多的结构信息，精度可能会得到进一步提高，这将是接下来要重点研究的内容。