CN114241067A - 一种基于多项式再生采样核的电容层析图像重建方法 - Google Patents

Abstract

一种基于多项式再生采样核的电容层析图像重建方法，将ECT图像建模为一个有限新息率FRI信号，然后用多项式再生采样核滤波后均匀采样；FRI采样过程的目的是从输入信号中提取特征信息，将原灵敏度矩阵和FRI测量矩阵结合，合并后的灵敏度矩阵通过零向量扩展来重构，并随机排列新灵敏度矩阵的行向量，得到对管道横截面两相流体介电常数分布的稀疏测量矩阵；以同样方式展开的电容，得到对管道横截面两相流体介电常数分布的稀疏测量向量；构建一个稀疏化的测量方程，将电容层析图像重建问题转换为L0范数最小化问题，最后通过求解L0范数优化问题来重构原始图像。本发明效率较高，精度较高，为电容层析成像研究提供了一种新的方法和手段。

Description

一种基于多项式再生采样核的电容层析图像重建方法

技术领域

本发明涉及图像处理技术领域，具体涉及一种基于多项式再生采样核的电容层析图像重建方法。

背景技术

电容层析成像技术(Electrical Capacitance Tomography，ECT)是过程层析成像技术(Process Tomography)中的一种技术，可以用来检测流动参数。通常用于工业管道内的多相流检测，近年来ECT技术已经被广泛应用。ECT技术具有非侵入性、适用范围广、结构简单、安全、成本低等许多的优点，因此ECT技术为解决两相流的可视化参数检测提供了一条较好的途径。

电容层析成像的基本思想是：不同流体具有不同的物理化学特性，表现为介电常数不相同。管道内部混合流体的介电常数会随着管道内部的流体分布变化而变化，这会导致电容传感器所测量得电容测量值发生变化。数据采集系统将电容测量值上传到计算机，通过高实时性的图像重建算法将介电常数的分布情况转换为可视化的图像。典型的 ECT系统结构图如图1所示，100为计算机，200是检测线圈，300是屏蔽层，400是励磁线圈。

被测两相流管道截面的介电常数分布ε(x,y)与从电容传感器获取的不同测量电极间的电容向量C_i之间的关系可以用式(1)来描述：

C_i＝∫∫_DS_i(x,y)ε(x,y)dxdy (1)

其中，C_i(i＝1,2...m)为测量电容值，m为电容传感器的两两组合个数， m是一个整数。D是被测的区域，(x,y)∈D。ε(x,y)是被测区域内部的介电常数分布函数；S_i(x,y)是对应于C_i的灵敏度分布函数，

因为在实际的ECT系统中，所划分像素个数远大于电容传感器的数量，即n＞m，要在测量值不足的情况下求方程的解。在数学上，这是一个欠定的问题。欠定问题的解不具有唯一性。图像重构算法是电容层析成像的关键所在，目前已经提出了不少比较成熟的图像重构算法，传统算法有Landweber算法、线性反投影法(LBP)、Tikhonov 正则化算法等。LBP算法重建复杂流型图像的质量相对较低，但LBP 算法的优点是简单快速；Landweber算法在大多数情况下都是一种有效的算法，与LBP算法相比，具有更好的图像质量，近年来得到了广泛的应用，但是这种算法也有缺点，它是通过最速下降来确定搜索的方向，在步长参数选择不适当的情况下易陷入局部极值；Tikhonov 正则化算法可以输出较好的图像，但这种算法重建出的图像边缘不是特别理想。迄今为止，提高电容层析图像重建的质量和精度仍然是个关键技术问题。

发明内容

为了克服现有技术的不足，针对电容层析成像问题，本发明提出一种高效率和高精度的基于多项式再生采样核的电容层析图像重建方法，将ECT图像建模为一个有限新息率(Finite Rate of Innovation，FRI) 信号，然后用多项式再生采样核滤波后均匀采样；FRI采样过程的目的是从输入信号中提取特征信息，将原灵敏度矩阵和FRI测量矩阵结合，合并后的灵敏度矩阵通过零向量扩展来重构，并随机排列新灵敏度矩阵的行向量，得到对管道横截面两相流体介电常数分布的稀疏测量矩阵；以同样方式展开的电容，得到对管道横截面两相流体介电常数分布的稀疏测量向量；构建一个稀疏化的测量方程，将电容层析图像重建问题转换为L0范数最小化问题，最后通过求解L0范数优化问题来重构原始图像。

本发明解决其技术问题所采用的技术方案是：

一种基于多项式再生采样核的电容层析图像重建方法，包括以下步骤：

步骤一，初始化，获取并归一化数据；

当介电常数相近时，电容向量C_i与介电常数分布函数ε(x,y)之间的关系可近似为线性，再对其进行离散处理，最后归一化，得到ECT 矩阵形式的数学模型：

C＝Sg (2)

其中，C是电容传感器得到的大小为m×1的电容测量值，是一个归一化后的向量；g是大小为n×1的灰度向量，对应于管道内部两相流介电常数分布；S是大小为m×n的灵敏度矩阵，是归一化后的矩阵；m 为测量极板间两两组合的个数，m为整数，n为所剖分的像素个数，n 为整数；由此公式(2)可得，ECT图像重建的问题转化为通过测量电容向量C求解灰度向量g的过程；

步骤二，有限新息率FRI信号建模

由于灰度向量是离散的，是一个典型的离散有限信息率FRI信号，因此把Tikhonov正则化算法得到的灰度向量g_v建模为狄拉克脉冲序列信号g(x)：

g_v＝(S^T·S+α·I)^-1·S^TC (3)

其中，g_v是一个大小为n×1的灰度向量，α是正则化参数，I是一个大小为n×n的单位阵，因为S^T·S仍然有不适定的情况，所以加入正则化参数α，使得S^T·S+α·I适定；通过求解广义逆矩阵，有效处理病态逆问题；

其中，x＝1,2...,N是像素的位置，L是g(x)非零灰度值的个数，a_l∈[0,1] 是像素x_l的灰度值。x_l∈{1,2,....,N}是非零灰度值的像素位置；

步骤三，有限新息率FRI采样

用多项式再生采样核

对狄拉克脉冲序列信号g(x)滤波后均匀采样，得到样本y_k(k＝1,2,...,K)；

均匀采样后得到K个样本y_k：

其中，y(x)＝g(x)*h(x)是FRI信号g(x)滤波后的结果，＜·,·＞表示内积的过程，y_k是滤波后均匀采样的样本值。T是采样间隔，K是FRI采样值的个数；

步骤四，从FRI采样样本中获取测量值，构建观测矩阵，用样条系数C_m,k对FRI采样值y_k进行加权求和，可以得到M个测量值,M个测量值可以由K个样本y_k计算得到：

假设u_m＝τ_m·T^m，公式(6)写成矩阵形式：

因为像素位置x_l属于集合{1,2,...,N}，a_l∈[0,1]是像素x_l的灰度值,测量值u_m表示为集合{1,2,...,N}中所有元素的完备线性组合，所以公式(7) 写成：

U＝Ag (9)

其中，U＝[u₀,u₁,...,u_M-1]^T∈R^M×1是测量向量，A是由集合

构成的M×N测量矩阵，g＝[g₁,g₂,...,g_N]^T是ECT的 N×1灰度向量；

步骤五，数据融合，构建具有稀疏性的测量方程，过程如下：

步骤5.1：结合FRI观测矩阵，由于公式(2)和(9)都是对g的一个多维观测，所以可以进行数据融合，得到一个新的测量方程，结合公式(2)和(9)得到新的测量方程：

步骤5.2：通过补零的方式增加虚拟电极，通过补零的方式增加虚拟电极，获得更多的测量电容，由于任意两个虚拟电极之间形成的灵敏度矩阵为0，所有虚拟电极之间的测量电容也为0，故ECT系统的数学模型重写为：

λ₀＝S₀g (11)

其中，λ₀是在数据融合后加入一些虚拟电极的归一化电容向量；S₀是 ECT系统在数据融合后加入一些虚拟电极的归一化灵敏度矩阵；

步骤5.3：随机重组λ₀和S₀的行向量，过程如下：加上虚拟电极后，对灵敏度矩阵S₀的每一行都进行随机排列，从而S₀成为一个新的灵敏度矩阵S_new，将电容向量和灵敏度矩阵以相同的顺序随机排列，相应的测量电容λ₀转化为λ_new，使其具有一定的随机性；这样，随机观测矩阵和稀疏基在一定程度上不相关，因此ECT系统的数学模型重写为：

λ_new＝S_newg (12)

其中，λ_new是管道横截面两相流体介电常数分布的测量向量；S_new是管道横截面两相流体介电常数分布的测量矩阵；

步骤六，求观测方程的稀疏解，过程如下：

步骤6.1：稀疏变换，只有当输入信号具有稀疏性时，即非零元素个数较少时才能进行压缩感知，ECT系统中的大多数流型的信号稀疏程度不能满足压缩感知的要求，为了使ECT系统的输入信号能够适应压缩感知，需对输入信号进行正交变换，将其转化为稀疏信号，正交变换的过程描述为：

g＝ψs (13)

其中，大小为n×n的ψ为g的稀疏基。原始信号g在稀疏基上的投影为大小为n×1的稀疏向量s；

步骤6.2：求稀疏解，由于s是稀疏的，所以通过求解L0范数最小化问题，获得稀疏解

其中，L0范数||s||₀表示向量s中非零系数的个数，通过正交匹配追踪算法(OMP)来求解，ψ为信号g的稀疏基，是n×n的矩阵；

步骤七，图像重建。原始图像信号被估计为：

进一步，所述步骤二中，把Tikhonov正则化算法得到的灰度向量g_v建模为狄拉克脉冲序列信号g(x)，输入信号可以替换成LBP算法、 Landweber算法得到的灰度向量，不局限于Tikhonov正则化算法。

再进一步，所述步骤三中，用多项式再生采样核

对狄拉克脉冲序列信号g(x)滤波后均匀采样，得到样本y_k＝(1,2,...,K)；

M-1阶的多项式再生采样核的公式如下：

其中，M-1＝5是再生采样核的阶数，M-1阶再生采样核可以再生出M 个多项式，样条系数C_m,k由下式得出：

其中，

是

的对偶函数，属于准正交函数。

更进一步，所述步骤六中，为了获得稀疏解

使用OMP算法求解L0范数最小化问题。OMP算法可以替换成匹配追踪算法、分段匹配追踪算法、迭代阈值算法等算法，不局限于OMP算法。

本发明的有益效果主要表现在：效率较高，精度较高，为电容层析成像研究提供了一种新的方法和手段。

附图说明

图1是典型的ECT系统结构图。

图2是基于多项式再生采样核的电容层析图像重建方法的流程图。

图3是经典的FRI采样结构图。

图4是二种流型的原始图像，其中，(a)是双球流，(b)是环状流。

图5是二种流型的LBP算法图像重构图，其中，(a)是双球流， (b)是环状流。

图6是二种流型的Landweber算法图像重构图，其中，(a)是双球流，(b)是环状流。

图7是二种流型的Tikhonov正则化算法图像重构图，其中，(a) 是双球流，(b)是环状流。

图8是二种流型的LBP算法图像重构图，其中，(a)是双球流， (b)是环状流。

具体实施方式

下面结合附图对本发明作进一步描述。

参照图1～图8，一种基于多项式再生采样核的电容层析图像重建方法，包括以下步骤：

步骤一，初始化，获取并归一化数据；

当介电常数相近时，电容向量C_i与介电常数分布函数ε(x,y)之间的关系可近似为线性，再对其进行离散处理，最后归一化，得到ECT 矩阵形式的数学模型：

C＝Sg (2)

其中，C是电容传感器得到的大小为m×1的电容测量值，是一个归一化后的向量；g是大小为n×1的灰度向量，对应于管道内部两相流介电常数分布；S是大小为m×n的灵敏度矩阵，是归一化后的矩阵；m 为测量极板间两两组合的个数，m为整数，n为所剖分的像素个数，n 为整数；由此公式(2)可得，ECT图像重建的问题转化为通过测量电容向量C求解灰度向量g的过程；

步骤二，有限新息率FRI信号建模

由于灰度向量是离散的，是一个典型的离散有限信息率FRI信号，因此把Tikhonov正则化算法得到的灰度向量g_v建模为狄拉克脉冲序列信号g(x)：

g_v＝(S^T·S+α·I)^-1·S^TC (3)

其中，g_v是一个大小为n×1的灰度向量，α是正则化参数，I是一个大小为n×n的单位阵，因为S^T·S仍然有不适定的情况，所以加入正则化参数α，使得S^T·S+α·I适定；通过求解广义逆矩阵，有效处理病态逆问题；

其中，x＝1,2...,N是像素的位置，L是g(x)非零灰度值的个数，a_l∈[0,1] 是像素x_l的灰度值。x_l∈{1,2,....,N}是非零灰度值的像素位置；

步骤三，有限新息率FRI采样

用多项式再生采样核

对狄拉克脉冲序列信号g(x)滤波后均匀采样，得到样本y_k(k＝1,2,...,K)；

经典的FRI采样结构如图3所示，均匀采样后得到K个样本y_k：

其中，y(x)＝g(x)*h(x)是FRI信号g(x)滤波后的结果，＜·,·＞表示内积的过程，y_k是滤波后均匀采样的样本值。T是采样间隔，K是FRI采样值的个数；

步骤四，从FRI采样样本中获取测量值，构建观测矩阵，用样条系数C_m,k对FRI采样值y_k进行加权求和，可以得到M个测量值,M个测量值可以由K个样本y_k计算得到：

假设u_m＝τ_m·T^m，公式(6)写成矩阵形式：

因为像素位置x_l属于集合{1,2,...,N}，a_l∈[0,1]是像素x_l的灰度值,测量值u_m表示为集合{1,2,...,N}中所有元素的完备线性组合，所以公式(7) 写成：

U＝Ag (9)

其中，U＝[u₀,u₁,...,u_M-1]^T∈R^M×1是测量向量，A是由集合

构成的M×N测量矩阵，g＝[g₁,g₂,...,g_N]^T是ECT的 N×1灰度向量；

步骤五，数据融合，构建具有稀疏性的测量方程，过程如下：

步骤5.1：结合FRI观测矩阵，由于公式(2)和(9)都是对g的一个多维观测，所以可以进行数据融合，得到一个新的测量方程，结合公式(2)和(9)得到新的测量方程：

步骤5.2：通过补零的方式增加虚拟电极，通过补零的方式增加虚拟电极，获得更多的测量电容，由于任意两个虚拟电极之间形成的灵敏度矩阵为0，所有虚拟电极之间的测量电容也为0，故ECT系统的数学模型重写为：

λ₀＝S₀g (11)

其中，λ₀是在数据融合后加入一些虚拟电极的归一化电容向量；S₀是 ECT系统在数据融合后加入一些虚拟电极的归一化灵敏度矩阵；

步骤5.3：随机重组λ₀和S₀的行向量，过程如下：加上虚拟电极后，对灵敏度矩阵S₀的每一行都进行随机排列，从而S₀成为一个新的灵敏度矩阵S_new，将电容向量和灵敏度矩阵以相同的顺序随机排列，相应的测量电容λ₀转化为λ_new，使其具有一定的随机性；这样，随机观测矩阵和稀疏基在一定程度上不相关，因此ECT系统的数学模型重写为：

λ_new＝S_newg (12)

其中，λ_new是管道横截面两相流体介电常数分布的测量向量；S_new是管道横截面两相流体介电常数分布的测量矩阵；

步骤六，求观测方程的稀疏解，过程如下：

步骤6.1：稀疏变换，只有当输入信号具有稀疏性时，即非零元素个数较少时才能进行压缩感知，ECT系统中的大多数流型的信号稀疏程度不能满足压缩感知的要求，为了使ECT系统的输入信号能够适应压缩感知，需对输入信号进行正交变换，将其转化为稀疏信号，正交变换的过程描述为：

g＝ψs (13)

其中，大小为n×n的ψ为g的稀疏基。原始信号g在稀疏基上的投影为大小为n×1的稀疏向量s；

步骤6.2：求稀疏解，由于s是稀疏的，所以通过求解L0范数最小化问题，获得稀疏解

其中，L0范数||s||₀表示向量s中非零系数的个数，通过正交匹配追踪算法(OMP)来求解，ψ为信号g的稀疏基，是n×n的矩阵；

步骤七，图像重建。原始图像信号被估计为：

进一步，所述步骤二中，把Tikhonov正则化算法得到的灰度向量g_v建模为狄拉克脉冲序列信号g(x)，输入信号可以替换成LBP算法、 Landweber算法得到的灰度向量，不局限于Tikhonov正则化算法。

再进一步，所述步骤三中，用多项式再生采样核

对狄拉克脉冲序列信号g(x)滤波后均匀采样，得到样本y_k＝(1,2,...,K)；

M-1阶的多项式再生采样核的公式如下：

其中，M-1＝5是再生采样核的阶数，M-1阶再生采样核可以再生出M 个多项式，样条系数C_m,k由下式得出：

其中，

是

的对偶函数，属于准正交函数。

更进一步，所述步骤六中，为了获得稀疏解

使用OMP算法求解L0范数最小化问题。OMP算法可以替换成匹配追踪算法、分段匹配追踪算法、迭代阈值算法等算法，不局限于OMP算法。

为了验证本发明方法的性能，进行了仿真实验验证。利用 COMSOL Multiphysics软件建立三维仿真模型并生成仿真数据m×1的电容向量C和m×n的灵敏度矩阵S并导出数据。利用Matlab R2019a 仿真工具对导出的数据进行处理，用图像重构算法重建ECT的图像。三维仿真模型建立如下:管道为正方形，8个电极传感器组成阵列，两相流管道的截面被分割为20×20个像素，较密介质的介电常数为3，较稀介质的介电常数为1。流型设置为双球流和环状流。多项式再生采样核阶数为5。即M-1＝5。

二种流型的原始图像如图4所示，分别为双球流和环状流，三种传统算法和本发明方法重建的效果如图5、图6、图7和图8所示。从图形的比较可以看出，用本发明方法重建的图像非常接近原始图像，图像的边缘较平滑。本发明方法重建的图像质量明显优于其他三种传统算法。

为了对图像重建方法的性能进行数值评价，分析重建的图像与原图像的误差和相关系数。误差可以定义为重建的图像灰度向量与原始图像灰度向量的差异程度，误差的公式为:

其中，g是原始图像的灰度向量。

是g的估计，表示为重建后图像的灰度向量。

相关系数可以定义为重建后的图像灰度向量与原始图像灰度向量之间的线性相关度，相关系数的公式为:

其中,

是向量g的平均值，而

是向量

的平均值。

表1给出了三个传统算法和本发明方法在双球流和环状流两种流型下的重建图像误差和相关系数。

表1

由表1可得，本发明方法重建的图像误差比三种传统算法小得多。同样，本发明方法重建图像的相关系数优于三种传统算法。这意味着用本发明方法重建的图像与原始图像最相似。显然，本发明方法效果良好，获得了较好的图像。

Claims (4)

1.一种基于多项式再生采样核的电容层析图像重建方法，其特征在于，所述方法包括以下步骤：

步骤一，初始化，获取并归一化数据；

当介电常数相近时，电容向量C_i与介电常数分布函数ε(x,y)之间的关系可近似为线性，再对其进行离散处理，最后归一化，得到ECT矩阵形式的数学模型：

C＝Sg (2)

其中，C是电容传感器得到的大小为m×1的电容测量值，是一个归一化后的向量；g是大小为n×1的灰度向量，对应于管道内部两相流介电常数分布；S是大小为m×n的灵敏度矩阵，是归一化后的矩阵；m为测量极板间两两组合的个数，m为整数，n为所剖分的像素个数，n为整数；由此公式(2)可得，ECT图像重建的问题转化为通过测量电容向量C求解灰度向量g的过程；

步骤二，有限新息率FRI信号建模

由于灰度向量是离散的，是一个典型的离散有限信息率FRI信号，因此把Tikhonov正则化算法得到的灰度向量g_v建模为狄拉克脉冲序列信号g(x)：

g_v＝(S^T·S+α·I)^-1·S^TC (3)其中，g_v是一个大小为n×1的灰度向量，α是正则化参数，I是一个大小为n×n的单位阵，因为S^T·S仍然有不适定的情况，所以加入正则化参数α，使得S^T·S+α·I适定；通过求解广义逆矩阵，有效处理病态逆问题；

其中，x＝1,2...,N是像素的位置，L是g(x)非零灰度值的个数，a_l∈[0,1]是像素x_l的灰度值。x_l∈{1,2,....,N}是非零灰度值的像素位置；

步骤三，有限新息率FRI采样

用多项式再生采样核

对狄拉克脉冲序列信号g(x)滤波后均匀采样，得到样本y_k(k＝1,2,...,K)；

均匀采样后得到K个样本y_k：

其中，y(x)＝g(x)*h(x)是FRI信号g(x)滤波后的结果，<·,·>表示内积的过程，y_k是滤波后均匀采样的样本值。T是采样间隔，K是FRI采样值的个数；

步骤四，从FRI采样样本中获取测量值，构建观测矩阵，用样条系数C_m,k对FRI采样值y_k进行加权求和，可以得到M个测量值,M个测量值可以由K个样本y_k计算得到：

假设u_m＝τ_m·T^m，公式(6)写成矩阵形式：

因为像素位置x_l属于集合{1,2,...,N}，a_l∈[0,1]是像素x_l的灰度值,测量值u_m表示为集合{1,2,...,N}中所有元素的完备线性组合，所以公式(7)写成：

U＝Ag (9)

其中，U＝[u₀,u₁,...,u_M-1]^T∈R^M×1是测量向量，A是由集合

构成的M×N测量矩阵，g＝[g₁,g₂,...,g_N]^T是ECT的N×1灰度向量；

步骤五，数据融合，构建具有稀疏性的测量方程，过程如下：

步骤5.1：结合FRI观测矩阵，由于公式(2)和(9)都是对g的一个多维观测，所以可以进行数据融合，得到一个新的测量方程，结合公式(2)和(9)得到新的测量方程：

步骤5.2：通过补零的方式增加虚拟电极，通过补零的方式增加虚拟电极，获得更多的测量电容，由于任意两个虚拟电极之间形成的灵敏度矩阵为0，所有虚拟电极之间的测量电容也为0，故ECT系统的数学模型重写为：

λ₀＝S₀g (11)

其中，λ₀是在数据融合后加入一些虚拟电极的归一化电容向量；S₀是ECT系统在数据融合后加入一些虚拟电极的归一化灵敏度矩阵；

步骤5.3：随机重组λ₀和S₀的行向量，过程如下：加上虚拟电极后，对灵敏度矩阵S₀的每一行都进行随机排列，从而S₀成为一个新的灵敏度矩阵S_new，将电容向量和灵敏度矩阵以相同的顺序随机排列，相应的测量电容λ₀转化为λ_new，使其具有一定的随机性；这样，随机观测矩阵和稀疏基在一定程度上不相关，因此ECT系统的数学模型重写为：

λ_new＝S_newg (12)

其中，λ_new是管道横截面两相流体介电常数分布的测量向量；S_new是管道横截面两相流体介电常数分布的测量矩阵；

步骤六，求观测方程的稀疏解，过程如下：

步骤6.1：稀疏变换，只有当输入信号具有稀疏性时，即非零元素个数较少时才能进行压缩感知，ECT系统中的大多数流型的信号稀疏程度不能满足压缩感知的要求，为了使ECT系统的输入信号能够适应压缩感知，需对输入信号进行正交变换，将其转化为稀疏信号，正交变换的过程描述为：

g＝ψs (13)

其中，大小为n×n的ψ为g的稀疏基。原始信号g在稀疏基上的投影为大小为n×1的稀疏向量s；

步骤6.2：求稀疏解，由于s是稀疏的，所以通过求解L0范数最小化问题，获得稀疏解

其中，L0范数||s||₀表示向量s中非零系数的个数，通过正交匹配追踪算法(OMP)来求解，ψ为信号g的稀疏基，是n×n的矩阵；

步骤七，图像重建。原始图像信号被估计为：

2.如权利要求1所述的一种基于多项式再生采样核的电容层析图像重建方法，其特征在于，所述步骤二中，把Tikhonov正则化算法得到的灰度向量g_v建模为狄拉克脉冲序列信号g(x)。

3.如权利要求1或2所述的一种基于多项式再生采样核的电容层析图像重建方法，其特征在于，所述步骤三中，用多项式再生采样核

对狄拉克脉冲序列信号g(x)滤波后均匀采样，得到样本y_k＝(1,2,...,K)；

M-1阶的多项式再生采样核的公式如下：

其中，M-1＝5是再生采样核的阶数，M-1阶再生采样核可以再生出M个多项式，样条系数C_m,k由下式得出：

其中，

是

的对偶函数，属于准正交函数。

4.如权利要求1或2所述的一种基于多项式再生采样核的电容层析图像重建方法，其特征在于，所述步骤六中，为了获得稀疏解

使用OMP算法求解L0范数最小化问题。

Images