基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法

技术领域

本发明属于信号处理技术领域，尤其涉及一种基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法。

背景技术

压缩感知是由Donoho，Candès和Tao等学者于2004年提出的一种高效的信号处理理论。该理论充分挖掘了自然信号固有的稀疏或可压缩特性，在远低于Nyquist采样频率的条件下，利用非相关测量矩阵将原始高维信号投影到低维空间，得到低维的测量值，并通过求解非线性优化问题可以无失真地重建原始信号。不同于传统的基于Nyquist采样定理的信号处理方法，压缩感知在信号采样的同时就可对其进行适当的压缩，直接获取信号的压缩测量值，为实现高分辨率的信号采样提供了新的思路，在图像处理、生物传感、无线通信、模拟信息转换等诸多领域均有着广泛的应用前景。

信号重建是压缩感知理论的核心问题之一。信号重建的实质是在已知测量矩阵和测量值的基础上，反过来求解原始信号的问题。常见的求解模型有最小l0范数模型和最小l1范数模型，与之相对应的重建算法即为基于l0范数的匹配追踪类算法和基于l1范数的凸优化类算法。

匹配追踪类算法通过选择测量矩阵中与测量值最匹配的列向量来构建原始信号的稀疏逼近。凸优化类算法用l1范数替代l0范数并通过线性规划方法求解，但该类算法计算复杂度高，无法应用于数据量较大的信号重建。相比而言，匹配追踪类算法以其较低的计算复杂度成为了当前重建算法的研究热点。

现有压缩感知的信号重建算法研究大多是对一般匹配追踪(MP)类重建算法(如OMP算法、ROMP算法和CoSaMP算法等)进行直接改进，均未考虑压缩感知过程中测量矩阵的作用特点或影响，而在基于MP类算法的信号重建方法中测量矩阵直接参与了信号的重建，因此测量矩阵的作用与影响不应忽视。测量矩阵也是压缩感知研究的核心问题之一，现有的测量矩阵研究表明，基于分块哈达玛(Hardmard)矩阵构建的测量矩阵相对其他测量矩阵具有性能优势，如计算简单、便于硬件实现，信号重建质量较高等。

针对具有特定结构与性能优势的测量矩阵，设计与其结构特点相适应的重建算法，不仅能保持高的信号重建质量，还有利于简化算法运算、提高其实时性、降低算法的空间复杂度，这些对一般应用系统尤其是处理能力相对弱、存储空间小的资源有限系统而言尤为重要。

发明内容

本发明的目的是结合维数任意的分块哈达玛测量矩阵的结构特点，提供基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法。

本发明提供了一种基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法，包括：

步骤S1：根据待测稀疏信号X的维数N和采样率α的积向上取整得到X的测量个数M，并构造正交的分块哈达玛对角矩阵Φ_N×N，取Φ_N×N的前M行作为测量矩阵Φ；

其中，N＝Σk_i2^i-1,i＝1,2,…,k_i＝0或1，0.8≥α≥0.1，或1，且k_i不能全部为0，l_i不能全部为0；

步骤S2：利用所述测量矩阵Φ对待测稀疏信号X进行观测，得到待测稀疏信号X的M维测量值Y；

步骤S3：根据M和N的值和所述测量矩阵Φ中维数不超过M的哈达玛子块的具体维数选择由测量值Y求解待测稀疏信号X的重建信号的组合求解方法，并按选定的具体组合求解方法，求解待测稀疏信号X的重建信号

其中，该待测稀疏信号可以是可压缩信号经稀疏化处理得到的。

在本发明的方案中，k_i和l_i的值是由M和N确定的，只要M和N不为0，k_i和l_i也就不会都为0。

进一步的，步骤S1中分块哈达玛测量矩阵Φ_N×N的构造规则为：

当N＝2ⁿ,n＝2,3,...且M≥2^n-1时，将两个2^n-1阶哈达玛矩阵块沿主对角线排列并归一化得到块对角矩阵Φ_N×N；

当N＝2ⁿ,n＝2,3,...且2^n-j≤M＜2^n-j+1,j＝1,2,3,...,n＞j时，将两个2^n-j阶哈达玛矩阵块和从2^n-j+1到2^n-1维的哈达玛矩阵块各一个，沿主对角线按从小到大顺序排列并归一化得到块对角矩阵Φ_N×N；

当2^n-1＜N＜2ⁿ,n＝2,3,...时，将对应k_i为1的Σk_i个不同维的哈达玛矩阵块沿主对角线排列并归一化得到块对角矩阵Φ_N×N，其中，2^k1+2^k2+…+2^kr＝N,r＝Σk_i。

其中，2^k维哈达玛矩阵块的表达式为：

进一步的，步骤S1中的分块哈达玛测量矩阵Φ_N×N为正交矩阵，其形式如下：

当N＝2ⁿ,n＝2,3,...且M≥2^n-1时，

当N＝2ⁿ,n＝2,3,...且2^n-j≤M＜2^n-j+1,j＝2,3,...,n＞j时，

当2^n-1＜N＜2ⁿ,n＝2,3,...时，

其中，分别代表维数为2^k1，2^k2，2^kr，2^n-1，2^n-j，2^n-j+1的哈达玛矩阵子块。

进一步的，步骤S2中，利用测量矩阵Φ对信号X进行观测的数学表达式为：Y＝ΦX。

进一步的，步骤S3中，由Y求解X的重建信号的求解方法的确定规则为：

当N＝2ⁿ,n＝2,3,...且M＝2^n-j,j＝1,2,3,...,n＞j时，根据可逆矩阵方程求解方法直接求得X的重建信号的前M维分量其中为Y的第1到第2^n-j维分量；其余维分量取0，即

当N＝2ⁿ,n＝2,3,...且2^n-j＜M＜2^n-j+1,j＝1,2,3,...,n＞j时，按所述M＝2^n-j时的求解方法求得的前2^n-j维分量，按一般匹配追踪类算法求解后N-2^n-j维分量，将前2^n-j维分量和后N-2^n-j维分量组合得到重建信号

当2^n-1＜N＜2ⁿ,n＝2,3,...且时，根据可逆矩阵方程求逆直接计算求解方法求得X的重建信号的前M维分量其中其余维的分量取0，即

其中，Y₁为Y的第1到第2ⁱ¹维分量，Y₂为Y的第2ⁱ¹+1到第2ⁱ²维分量，Y_m为Y的第2^m-1+1到第2^m维分量；

当2^n-1＜N＜2ⁿ,n＝2,3,...且时，按时的求解方法求的前维分量，按一般匹配追踪类算法求的后维分量，将前维分量和后维分量组合得到N维的重建信号

其中，步骤S3的确定规则中的一般匹配追踪类算法包括：正交匹配追踪算法(OMP)、正则化正交匹配追踪算法(ROMP)、压缩采样匹配追踪算法(CoSaMP)等。

有益效果

本发明提供的一种基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法，首先通过待测稀疏信号的维数和测量个数的取值范围来构造便于通过对哈达玛矩阵块求逆来直接求解的测量矩阵，从而可以根据稀疏信号的维数、对应采样率的测量个数和测量矩阵中的哈达玛子块的维数来选择对应的信号重建方法；然后当测量个数正好等于测量矩阵中哈达玛子块维数不超过测量个数的哈达玛子块的维数之和时，采用可逆矩阵方程来求解重建信号的全部非0分量；否则，当测量个数大于测量矩阵中某几个哈达玛子块的维数之和时，将压缩感知中的测量信号重建分解为两个部分，即：第一部分是求可逆矩阵方程的解、第二部分是求欠定方程组的稀疏解，前者通过可逆矩阵方程求逆直接求得重建信号中部分分量的精确解，后者根据一般匹配追踪类重建算法得到余下分量的最优解。与直接利用一般匹配追踪类算法求取最优解相比，由于本发明提出的方法中利用了可逆矩阵方程求得重建信号中部分分量的精确解，所以提高了信号的重建质量，尤其对于经过稀疏化处理得到的主要信息集中于低频段的稀疏信号；同时，因为减小了需要求欠定方程组的稀疏解的规模，且在由可逆矩阵方程求精确解的过程中省去了最小二乘求解过程，因而能显著降低计算复杂度，减少信号重建所需的运行时间和对存储空间的占用。

附图说明

图1是本发明提出的基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法实现的流程图；

图2是不同测量值数目下，采用OMP算法、ROMP算法、CoSaMP算法与本发明提出的方法对大小为512×512的Lena图像压缩测量后重建的峰值信噪比的比较曲线图：其中，图2(a)为采用OMP算法与采用本发明方法对Lena图像重建的峰值信噪比比较曲线；图2(b)为采用ROMP算法与采用本发明方法对Lena图像重建的峰值信噪比比较曲线；图2(c)为采用CoSaMP算法与采用本发明方法对Lena图像重建的峰值信噪比比较曲线；

图3是测量值数目为200时，OMP算法、ROMP算法、CoSaMP算法与本发明方法对大小为512×512的Lena图像压缩测量后重建的效果图：其中，图3(a)为Lena原始图像；图3(b)、(c)分别为OMP算法与本发明方法对Lena重建的图像；图3(d)、(e)分别为ROMP算法与本发明方法对Lena重建的图像；图3(f)、(g)分别为ROMP算法与本发明方法对Lena重建的图像；

图4是测量值数目为300时，OMP算法、ROMP算法、CoSaMP算法与本发明方法对大小为512×512的Lena图像压缩测量后重建的效果图：其中，图4(a)为Lena原始图像；图4(b)、(c)分别为OMP算法与本发明方法对Lena重建的图像；图4(d)、(e)分别为ROMP算法与本发明方法对Lena重建的图像；图4(f)、(g)分别为ROMP算法与本发明方法对Lena重建的图像；

图5是不同测量值数目下，采用OMP算法、ROMP算法、CoSaMP算法与本发明提出的方法对大小为480×480的Lena图像压缩测量后重建的峰值信噪比的比较曲线图：其中，图5(a)为采用OMP算法与采用本发明方法对Lena图像重建的峰值信噪比比较曲线；图5(b)为采用ROMP算法与采用本发明方法对Lena图像重建的峰值信噪比比较曲线；图5(c)为采用CoSaMP算法与采用本发明方法对Lena图像重建的峰值信噪比比较曲线；

图6是测量值数目为200时，OMP算法、ROMP算法、CoSaMP算法与本发明提出的方法对大小为480×480的Lena图像压缩测量后重建的效果图：其中，图6(a)为Lena原始图像；图6(b)、(c)分别为OMP算法与本发明方法对Lena压缩测量后重建的图像；图6(d)、(e)分别为ROMP算法与本发明方法对Lena压缩测量后重建的图像；图6(f)、(g)分别为ROMP算法与本发明方法对Lena压缩测量后重建的图像；

图7是测量值数目为300时，OMP算法、ROMP算法、CoSaMP算法与本发明提出的方法对大小为480×480的Lena图像压缩测量后重建的效果图：其中，图7(a)为Lena原始图像；图7(b)、(c)分别为OMP算法与本发明方法对Lena压缩测量后重建的图像；图7(d)、(e)分别为ROMP算法与本发明方法对Lena压缩测量后重建的图像；图7(f)、(g)分别为ROMP算法与本发明方法对Lena压缩测量后重建的图像。

具体实施方式

下面结合附图对本发明的具体实施方式做进一步的描述。

图1示出了本发明提出的基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法的工作流程，具体为：

步骤S1：根据待测稀疏信号X的维数N(N＝Σk_i2^i-1,i＝1,2,…,k_i＝0或1)和采样率α(0.8≥α≥0.1)的积向上取整得到X的测量个数或1)，并构造正交的分块哈达玛对角矩阵Φ_N×N，取Φ_N×N的前M行作为测量矩阵Φ。

其中，分块哈达玛测量矩阵Φ_N×N的构造规则为：

当N＝2ⁿ(n＝2,3,...)且M≥2^n-1时，将两个2^n-1阶哈达玛矩阵块沿主对角线排列并归一化得到块对角矩阵Φ_N×N；

当N＝2ⁿ,n＝2,3,...且2^n-j≤M＜2^n-j+1,j＝1,2,3,...,n＞j时，将两个2^n-j阶哈达玛矩阵块和从2^n-j+1到2^n-1维的哈达玛矩阵块各一个，沿主对角线按从小到大顺序排列并归一化得到块对角矩阵Φ_N×N；

当2^n-1＜N＜2ⁿ,n＝2,3,...时，将对应k_i为1的Σk_i个不同维的哈达玛矩阵块沿主对角线排列并归一化得到块对角矩阵Φ_N×N。

其中，2^k维哈达玛矩阵块的表达式为：

该分块哈达玛测量矩阵Φ_N×N具体为正交矩阵，其形式如下：

当N＝2ⁿ,n＝2,3,...且M≥2^n-1时，

当N＝2ⁿ,n＝2,3,...,且2^n-j≤M＜2^n-j+1,j＝1,2,3,...,n＞j时，

当2^n-1＜N＜2ⁿ,n＝2,3,...时，

其中，分别代表维数为2^k1，2^k2，2^kr，2^n-1，2^n-j，2^n-j+1的哈达玛矩阵子块。

步骤S2：利用分块哈达玛测量矩阵对待测稀疏信号X进行观测，得到待测稀疏信号X的M维测量值Y，数学表达式为：Y＝ΦX。

步骤S3：根据M和N的值和所述测量矩阵Φ中维数不超过M的哈达玛子块的具体维数选择由测量值Y求解待测稀疏信号X的重建信号的组合求解方法，并按选定的具体组合求解方法，求解待测稀疏信号X的重建信号

而求解方法的确定规则为：

当N＝2ⁿ,n＝2,3,...,且M＝2^n-j,j＝1,2,3,...,n＞j时，根据可逆矩阵方程求解方法直接求得X的重建信号的前M维分量其中，为Y的第1到第2^n-j维分量；其余维分量取0，即

当N＝2ⁿ,n＝2,3,...,且2^n-j＜M＜2^n-j+1,j＝1,2,3,...,n＞j时，按所述M＝2^n-j时的求解方法求得的前2^n-j维分量，按一般匹配追踪类算法(MP)求解后N-2^n-j维分量，将前2^n-j维分量和后N-2^n-j维分量组合得到N维的重建信号

当2^n-1＜N＜2ⁿ,n＝2,3,...且时，通过可逆矩阵方程求逆矩阵后直接计算求解的方法求得X的重建信号的前M维分量其中其余的分量取0，即其中，Y₁为Y的第1到第2ⁱ¹维分量，Y₂为Y的第2ⁱ¹+1到第2ⁱ²维分量，Y_m为Y的第2^m-1+1到第2^m维分量；

当2^n-1＜N＜2ⁿ,n＝2,3,...且时，按时的求解方法求的前维分量，按一般匹配追踪类算法求的后维分量，将前维分量和后维分量组合得到N维的重建信号

上述的一般匹配追踪类算法包括：正交匹配追踪(OMP)算法、正则化正交匹配追踪(ROMP)算法、压缩采样匹配追踪(CoSaMP)算法等。

下面通过具体实施例来验证本发明所提算法的优良性能。需要指出的是，该实施例只是示例性的，并不是要限制本发明的适用范围。

实施例一：

对大小为512×512的Lena、Baboon、Peppers和Plane灰度图像进行压缩重建，具体实施过程为：

一、采用离散余弦变换(DCT)基对上述图像信号进行稀疏化处理得到稀疏信号X。

二、构造维数任意的分块哈达玛测量矩阵。

因N＝512，对应大于0.5的采样率α，相应的M>256，例如测量值数目M取300时，将两个维数均为256的哈达玛矩阵子块H₂₅₆沿对角线排列并归一化得到矩阵Φ_512×512：

然后选取Φ_512×512的前M行并归一化作为测量矩阵Φ；

对应小于0.5且大于0.25的采样率α(对应256>M>128)，例如测量值数目M为200时，将两个维数为128的哈达玛矩阵子块H₁₂₈、一个维数为256的哈达玛矩阵子块H₂₅₆按维数从小到大沿对角线排列并归一化得到矩阵Φ_512×512：

然后选取Φ_512×512的前M行并归一化作为测量矩阵Φ。

三、利用维数任意的分块哈达玛测量矩阵Φ对图像稀疏化处理后的信号X进行压缩测量。

四、分别利用本发明提出的方法和一般匹配追踪类算法(OMP算法、ROMP算法和CoSaMP算法)对压缩后图像进行重建。

附图2给出了不同测量值数目下，OMP算法、ROMP算法、CoSaMP算法与本发明提出的方法对512×512的Lena图像压缩测量后重建信号的峰值信噪比对比曲线。依据本发明方法，先要根据M或采样率α值确定分块哈达玛对角阵中的子块。当N＝2ⁿ即2的整数次幂时，若M≥2^n-1，则块对角阵由两个相同的2^n-1维块对角排列组成；若2^n-2≤M<2^n-1,则将上一情况中左上角的2^n-1维块分解为两个2^n-2维块，即对角阵由两个2^n-2维和一个2^n-1维块组成；若2^n-3≤M<2^n-2,则要将上一情况中对角阵左上角的2^n-2维块分解为两个2^n-3维块；若M更小，则继续将上一情况中左上角的最小子块分解为两个低一级的子块。实际中，对应N＝512的最小子块的维数不宜小于32，因为此时的采样率已低于0.1，采样率过低将无法保证信号重建成功。在本例的峰值信噪比的对比曲线中，当M取值在150～250这一段时，本发明方法的峰值信噪比曲线上升趋势不明显，其原因是：对应N＝512的情况，在上述M值范围内的测量矩阵中子块维数最大为128，也即能实现精确求解重建信号的分量为128个，超出128的部分只能采用一般MP类方法近似求解。M值在300～400范围时曲线缓变的原因也相同。但无论M取何值，都能看出本发明方法相对基于一般MP类方法的重建算法性能优势明显。

图3和图4分别给出了测量值数目分别为200和300时，OMP算法、ROMP算法、CoSaMP算法与本发明方法对512×512的Lena图像压缩测量后重建的效果图。容易看出，本发明方法的重建图像比一般MP类算法的重建图像更清晰。

为不失一般性，表1给出了测量值数目分别为200、300和400时，512×512的Baboon、Peppers和Plane图像在本发明方法和一般匹配追踪类算法下，重建图像的峰值信噪比。

表1本发明方法与一般MP类算法对512×512的三幅典型图像进行压缩测量后重建信号的峰值信噪比

从表1中可看出，本发明方法的重建图像的峰值信噪比明显高于一般匹配追踪类算法的重建图像峰值信噪比。这跟前面对Lena图像重建所得结论一致。

表2给出了不同测量值数目下，三个典型MP类算法与本发明方法对512×512的Lena图像重建的运行时间比较。

表2OMP算法、ROMP算法、CoSaMP算法与本发明方法对512×512的Lena图像的压缩测量信号重建的运行时间

(单位：s)，CPU：Inter(R)Core(TM)i5-2400，运行环境：Matlab R2014a

从表2可看出，本发明方法的运行时间明显少于一般MP类算法。表2的数据对比表明，利用基于维数任意的分块哈达玛测量矩阵的结构与性能特点，将信号重建分两部分，即：采用可逆矩阵方程解析求解重建信号中的大部分(不低于50％)和采用MP方法近似求解重建信号中的小部分(低于50％)，因为可逆矩阵方程直接求解的方式具有简单、高效和精确的特点以及通过减小MP方式求解的规模能降低运算量，从而使本发明方法可显著缩短信号重建算法运行时间。

实施例2：

对大小为480×480(N＝480＝2⁸+2⁷+2⁶+2⁵≠2ⁿ)的Lena、Baboon、Peppers和Plane灰度图像进行压缩重建。具体实施过程为：

一、采用离散余弦变换基对上述图像进行稀疏化处理后得到稀疏信号X。

二、构造维数任意的分块哈达玛测量矩阵。

对应N＝480的情况，将维数分别为32、64、128和256的哈达玛矩阵子块H₃₂、H₆₄、H₁₂₈和H₂₅₆按维数从小到大沿对角线排列并归一化得到矩阵Φ_480×480：

然后选取Φ_480×480的前M行并归一化作为测量矩阵Φ。

三、利用维数任意的分块哈达玛测量矩阵Φ对稀疏表示后图像进行压缩测量。

四、分别利用本发明提出的方法和一般匹配追踪类算法(OMP算法、ROMP算法和CoSaMP算法)对压缩后图像进行重建。

附图5给出了不同测量值数目下，OMP算法、ROMP算法、CoSaMP算法与本发明提出的方法对480×480的Lena图像压缩测量后重建的峰值信噪比比较曲线。这里，每组曲线中，本发明方法用到的一般匹配追踪类算法仍然与对比算法相同。从这三组曲线可以看出，本发明方法的重建图像的峰值信噪比明显高于一般匹配追踪类算法的重建图像的峰值信噪比。

附图6和附图7给出了测量值数目分别为200和300时，OMP算法、ROMP算法、CoSaMP算法与本发明方法对480×480的Lena图像压缩测量后重建的效果图。可以看出，本发明方法的重建图像比一般MP类算法的重建图像清晰度高。

为了不失一般性，表3给出了测量值数目为200、300和400时，480×480的Baboon、Peppers和Plane图像在本发明方法和一般匹配追踪类算法下，重建图像的峰值信噪比。

表3本发明方法和一般匹配追踪类算法对480×480的其他图像压缩测量后重建的峰值信噪比

从表3可以看出，本发明方法的重建图像的峰值信噪比比一般匹配追踪类算法的平均高出1～4dB，这跟Lena图像得到的结论一致。

表4给出了不同测量值数目下，本发明方法和三种典型一般MP类算法对480×480的Lena图像重建的运行时间。

表4本发明方法和三个典型MP类算法对480×480的Lena图像重建的运行时间(单位：s)

从表4可看出，与一般匹配追踪类算法的运行时间比较，本发明方法的运行时间显著减少。

综上所述，本发明提供的一种基于维数任意的分块哈达玛测量矩阵的压缩感知信号重建方法，首先通过待测稀疏信号的维数和测量个数的取值范围来构造便于通过对哈达玛快矩阵块求逆来直接求解的测量矩阵，从而可以根据稀疏信号的维数、对应采样率的测量个数和测量矩阵中的哈达玛子块的维数来选择对应的信号重建方法；然后当测量个数正好等于测量矩阵中哈达玛子块维数不超过测量个数的哈达玛子块的维数之和时，采用可逆矩阵方程来求解重建信号的全部非0分量；否则，当测量个数大于测量矩阵中某几个哈达玛子块的维数之和时，将压缩感知中的测量信号重建分解为两个部分，即：第一部分是求可逆矩阵方程的解、第二部分是求欠定方程组的稀疏解，前者通过可逆矩阵方程求逆直接求得重建信号中部分分量的精确解，后者根据一般匹配追踪类重建算法得到余下分量的最优解。与直接利用一般匹配追踪类算法求取最优解相比，由于本发明提出的方法中利用了可逆矩阵方程求得重建信号中部分分量的精确解，从而提高了信号的重建质量；同时，因为减小了需要求欠定方程组的稀疏解的规模，在由可逆矩阵方程求精确解的过程中省去了最小二乘求解过程，因而能显著降低计算复杂度，减少信号重建所需的运行时间和对存储空间的占用。

以上所述仅为本发明的实施例而已，并不用以限制本发明，凡在本发明精神和原则之内，所作任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。