CN107170019B - 一种快速低存储图像压缩感知方法 - Google Patents

一种快速低存储图像压缩感知方法 Download PDF

Info

Publication number
CN107170019B
CN107170019B CN201710380465.8A CN201710380465A CN107170019B CN 107170019 B CN107170019 B CN 107170019B CN 201710380465 A CN201710380465 A CN 201710380465A CN 107170019 B CN107170019 B CN 107170019B
Authority
CN
China
Prior art keywords
reconstruction
observation matrix
low
sampling
column
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201710380465.8A
Other languages
English (en)
Other versions
CN107170019A (zh
Inventor
王金铭
叶时平
徐振宇
陈超祥
蒋燕君
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang Shuren University
Original Assignee
Zhejiang Shuren University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang Shuren University filed Critical Zhejiang Shuren University
Priority to CN201710380465.8A priority Critical patent/CN107170019B/zh
Publication of CN107170019A publication Critical patent/CN107170019A/zh
Application granted granted Critical
Publication of CN107170019B publication Critical patent/CN107170019B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T9/00—Image coding
    • H—ELECTRICITY
    • H03—ELECTRONIC CIRCUITRY
    • H03M—CODING; DECODING; CODE CONVERSION IN GENERAL
    • H03M7/00—Conversion of a code where information is represented by a given sequence or number of digits to a code where the same, similar or subset of information is represented by a different sequence or number of digits
    • H03M7/30—Compression; Expansion; Suppression of unnecessary data, e.g. redundancy reduction
    • H03M7/3059—Digital compression and data reduction techniques where the original information is represented by a subset or similar information, e.g. lossy compression
    • H03M7/3062—Compressive sampling or sensing

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Multimedia (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Image Processing (AREA)
  • Medicines Containing Antibodies Or Antigens For Use As Internal Diagnostic Agents (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本发明公开了一种快速低存储图像压缩感知方法，属于图像处理领域。本发明采用低阶随机观测矩阵对原始信号进行全局采样和分块重构的方法，克服了传统压缩感知方法中观测矩阵需要占用大量存储空间和内存空间的缺点，同时克服了面对大尺寸图像重构实时性差的缺点，降低了随机观测矩阵所需的存储空间，降低了重构的计算复杂度，提升了重构的实时性。

Description

一种快速低存储图像压缩感知方法
技术领域
本发明属于图像处理领域，具体涉及一种快速低存储图像压缩感知方法。
背景技术
自压缩感知(Compressed Sensing，CS)理论提出以来，国内外研究人员的都在不断完善压缩感知的理论及拓展压缩感知的应用，但随着研究的深入，特别是针对压缩感知在大尺寸图像压缩采样和重构的研究中一直存在着亟待改善和解决的问题。
第一，在信号的非相关观测方面，虽然随机观测矩阵的投影方法具有理论上的完美特性，然而正由于其随机的特性，随机观测矩阵在硬件实现、存储分配和重构算法构造，都需要占用大量的存储空间和内存空间，在实际应用中受到很大的限制。为此，相关的科研人员针对观测矩阵占用存储空间大、内存利用率高等问题分别提出了各种解决方法，如分块压缩感知(Block Compressed Sensing，BCS)，该方法将原始图像划分成大小相等的若干块(如16×16)，根据采样率，设计相应大小的观测矩阵对原始信号进行分块采样。根据采样信息进行分块重构或全局重构。此外，Kronecker压缩感知方法、结构化观测矩阵、低秩观测矩阵、确定性观测矩阵等等都是能够有效降低观测矩阵存储空间的方法。
第二，在信号的优化重构过程中，随着图像尺寸规模的增加，重构算法的运算量将呈指数级上升，使整幅图像的重构过程变得非常耗时，并且重构过程中数据存储量和内存占用率将成倍增加，大大降低了压缩感知应用的实时性能。如设定采样率为0.5时，利用迭代重加权算法(Iterative Re-weighted Least Squares，IRLS)重构一幅256×256大小的图像，所需的时间约为60s左右，若重构一幅512×512大小的图像，则将需要约1000s的时间，但若重构一幅大小为1024×1024，其重构时间则是无法忍受的。此外，虽然有一些重构算法具备较快的重构速度，如正交匹配追踪算法(Orthogonal Matching Pursuit，OMP)，但其重构精度却稍劣于IRLS。所以找到一种即具备较高的重构精度，同时较快重构速度的重构算法仍是需要研究的问题。
发明内容
本发明的目的在于解决现有技术中存在的问题，并提供一种快速低存储图像压缩感知方法。
本发明所采用的具体技术方案如下：
快速低存储图像压缩感知方法，至少包括以下步骤：
步骤101：输入大小为N×N的测试图像，对测试图像进行小波变换，得到小波系数θN×N；
步骤102：根据设定的采样率M/N，生成低阶正交随机高斯观测矩阵Φ(M/t)×(N/t)，其大小为
Figure GDA0001332605030000021
其中，N表示原始测试图像小波系数任意列向量的长度，M表示针对任意列的测量数，t为分组数，
Figure GDA0001332605030000022
为正整数；
步骤103：结合半张量积运算方法，利用低阶高斯随机观测矩阵Φ(M/t)×(N/t)对小波系数进行压缩采样，得到相应的采样值yM×N，压缩采样方法为：
Figure GDA00013326050300000210
其中符号
Figure GDA00013326050300000211
表示向量的左半张量积运算；
步骤104：对大小为M×N的测量值，对其中的任一列yN×1进行如下方式的分组处理，将一列测量值划分成t组：
Figure GDA0001332605030000023
其中
Figure GDA0001332605030000024
表示进行分组后的第i组测量值，是长度为M/t的列向量，i＝1，2,3，…，t；y(M/t-1)t+1为列向量中的元素；
步骤105：对于每一组测量值
Figure GDA0001332605030000025
利用低阶高斯随机观测矩阵Φ(M/t)×(N/t)，设定原始信号长度N/t，进行分块重构，并得到每一组重构结果
Figure GDA0001332605030000026
其中
Figure GDA0001332605030000027
步骤106：根据步骤105所得的每一组重构结果
Figure GDA0001332605030000028
进行如下方式的排序，得到一列稀疏信号的最终重构结果
Figure GDA0001332605030000029
Figure GDA0001332605030000031
式中
Figure GDA0001332605030000032
表示
Figure GDA0001332605030000033
中第(jt+1)至((j+1)t)个元素，其中j＝1，2，…，N/t；
Figure GDA0001332605030000034
为第i列重构稀疏解中的第j个元素；
步骤107：对剩余的N-1列测量值进行步骤104、105所述方法重构，并按照步骤106所述对每一列重构结果进行排序，最终得原始图像稀疏系数的重构结果
Figure GDA0001332605030000035
步骤108：对重构稀疏系数
Figure GDA0001332605030000036
进行小波逆变换，得到大小为N×N的重构图像。
作为优选，所述步骤102中根据设定的采样率M/N，生成低阶正交随机高斯观测矩阵Φ(M/t)×(N/t)的具体步骤如下：
步骤1021：设定采样率M/N；
步骤1022：设定参数t，确定低阶观测矩阵大小(M/t)×(N/t)；
步骤1023：根据设定观测矩阵大小，系统生成正交高斯随机矩阵Φ(M/t)×(N/t)。
进一步的，所述步骤103中的利用低阶正交高斯随机观测矩阵Φ(M/t)×(N/t)对原始信号的小波系数进行全局压缩采样，包括如下步骤：
步骤1031：输入低阶正交高斯随机观测矩阵Φ(M/t)×(N/t)和原始信号小波系数θN×N；
步骤1032：根据设定参数t，将小波系数θN×N的每一列
Figure GDA0001332605030000037
分割成N/t个等长的块：
Figure GDA0001332605030000038
所有块均为t维列向量，得到分块后的小波系数：
Figure GDA0001332605030000039
步骤1033：取分块后的一列小波系数
Figure GDA00013326050300000310
利用已生成的低阶高斯随机观测矩阵Φ(M/t)×(N/t)与分块后的小波系数进行压缩采样，得到一列长度为M的采样值，压缩采样方法为：
Figure GDA0001332605030000041
其中φi,j∈Φ(M/t)×(N/t)，i＝1,2,…,M/t，j＝1,2,…,N/t；
步骤1034：利用低阶高斯随机观测矩阵完成每一列稀疏系数的压缩采样后，得到最后的压缩测量值yM×N。
进一步的，所述步骤105中对于每一组测量值
Figure GDA0001332605030000042
采用现有的基于lq-范数的迭代重加权方法进行重构，包括如下步骤：
步骤1051：输入一组测量值
Figure GDA0001332605030000043
低阶高斯随机观测矩阵Φ(M/t)×(N/t)，原始信号长度N/t；
步骤1052：初始化参数：令w(0)＝(1,…,1)N/t×1，ε0＝1，
Figure GDA0001332605030000044
j＝1，2，…，N/t；
步骤1053：根据设定的w(n)构建大小为(N/t)×(N/t)的对角矩阵D(n)；
步骤1054：根据如下方式计算并更新
Figure GDA0001332605030000045
Figure GDA0001332605030000046
步骤1055：判断
Figure GDA0001332605030000047
是否小于
Figure GDA0001332605030000048
若满足，则执行步骤1056，否则跳过步骤1056；
步骤1056：更新εn值，εn+1＝ρεn，其中0<ρ<1；
步骤1057：按照如下方式计算并更新w(n)：
Figure GDA0001332605030000049
其中
Figure GDA00013326050300000410
表示D(n)中第p个对角分量(p＝1,2,…,N/t)；q值的取值范围为(0，1)；
Figure GDA0001332605030000051
为第n次迭代重构的结果；
步骤1058：更新迭代次数，并判断εn是否满足条件，若不满足，返回步骤1053，继续迭代计算；若满足则退出迭代，并输出一组稀疏解
Figure GDA0001332605030000052
步骤1059：依次对剩余的每组测量值
Figure GDA0001332605030000053
进行如步骤1051～1058的重构，输出全部t组重构稀疏系数。
利用本发明的方法，构建低阶随机观测矩阵，利用该观测矩阵对原始信号进行全局采样和分块重构，经比较，利用本发明的采样和重构方法，在保证重构质量的同时，即可成倍降低观测矩阵的存储空间(最小可降为传统压缩感知方法中观测矩阵所需存储空间的1/1024)，又可大大提升重构的实时性，如一幅512×512大小的图像，利用本发明方法结合OMP重构算法，最短的重构时间仅需1.5s左右，比传统方法的提升了几十倍的速度。
附图说明
图1为本发明的流程图；
图2为不同大小高斯随机观测矩阵压缩采样和重构2维灰度图像Peppers(256×256)比较；
图3为不同大小高斯随机观测矩阵压缩采样和重构2维灰度图像Lena(512×512)比较；
图4为不同大小高斯随机观测矩阵压缩采样和重构2维灰度图像Mandrill(1024×1024)比较。
具体实施方式
下面结合附图和具体实施方式对本发明做进一步阐述和说明。本发明中各个实施方式的技术特征在没有相互冲突的前提下，均可进行相应组合。
如图1所示，一种快速低存储图像压缩感知方法，包括以下步骤：
步骤101：输入大小为N×N的测试图像，对测试图像进行小波变换，得到小波系数θN×N；
步骤102：根据设定的采样率M/N，生成低阶正交随机高斯观测矩阵Φ(M/t)×(N/t)，其大小为
Figure GDA0001332605030000054
其中，N表示原始测试图像小波系数任意列向量的长度，M表示针对任意列的测量数，t为分组数，
Figure GDA0001332605030000055
为正整数；
该步骤中根据设定的采样率M/N，生成低阶正交随机高斯观测矩阵Φ(M/t)×(N/t)的具体步骤如下：
步骤1021：设定采样率M/N；
步骤1022：设定参数t，确定低阶观测矩阵大小(M/t)×(N/t)；
步骤1023：根据设定观测矩阵大小，系统生成正交高斯随机矩阵Φ(M/t)×(N/t)。
相比传统压缩感知方法，本发明中对大小为N×N原始稀疏信号进行压缩采样(其采样率为M/N)，所需高斯随机观测矩阵的大小仅为(M/t)×(N/t)，对于相同的数据类型，Φ(M/t)×(N/t)所需的存储空间仅为ΦM×N的1/t2；如利用传统压缩感知方法在硬件上实现对一幅大小为1024×1024原始图像的采样和重构(假设采样率为25％)，若系统中存储的数据均采用单精度浮点数格式，则高斯随机观测矩阵需要1024KB的存储空间；若选取t＝16，采用本发明所述方法，其观测矩阵所需的存储空间仅为4KB；若选取t＝32，采用本发明所述方法，其观测矩阵所需的存储空间仅为1KB；由此，本发明所述方法可以大大降低观测矩阵所需的存储空间，进而降低重构过程的内存占用率，计算复杂度；此外，当有t＝1时，本发明所述的步骤103就退化为传统的压缩感知方法。
步骤103：结合半张量积运算方法，利用低阶高斯随机观测矩阵Φ(M/t)×(N/t)对小波系数进行压缩采样，得到相应的采样值yM×N，压缩采样方法为：
Figure GDA0001332605030000064
其中符号
Figure GDA0001332605030000065
表示向量的左半张量积运算；
该步骤中的利用低阶正交高斯随机观测矩阵Φ(M/t)×(N/t)对原始信号的小波系数进行全局压缩采样，包括如下步骤：
步骤1031：输入低阶正交高斯随机观测矩阵Φ(M/t)×(N/t)和原始信号小波系数θN×N；
步骤1032：根据设定参数t，将小波系数θN×N的每一列
Figure GDA0001332605030000066
分割成N/t个等长的块：
Figure GDA0001332605030000061
所有块均为t维列向量，得到分块后的小波系数：
Figure GDA0001332605030000062
步骤1033：取分块后的一列小波系数
Figure GDA0001332605030000063
利用已生成的低阶高斯随机观测矩阵Φ(M/t)×(N/t)与分块后的小波系数进行压缩采样，得到一列长度为M的采样值，压缩采样方法为：
Figure GDA0001332605030000071
其中φi,j∈Φ(M/t)×(N/t)，i＝1,2,…,M/t，j＝1,2,…,N/t；
步骤1034：利用低阶高斯随机观测矩阵完成每一列稀疏系数的压缩采样后，得到最后的压缩测量值yM×N。
步骤104：对大小为M×N的测量值，对其中的任一列yN×1进行如下方式的分组处理，将一列测量值划分成t组：
Figure GDA0001332605030000072
其中
Figure GDA0001332605030000073
表示进行分组后的第i组测量值，是长度为M/t的列向量，i＝1，2,3，…，t；y(M/t-1)t+1为列向量中的元素；
步骤105：对于每一组测量值
Figure GDA0001332605030000074
利用低阶高斯随机观测矩阵Φ(M/t)×(N/t)，设定原始信号长度N/t，可以采用现有的重构方法(如IRLS和OMP等)进行分块重构，并得到每一组重构结果
Figure GDA0001332605030000075
其中
Figure GDA0001332605030000076
该步骤中对于每一组测量值
Figure GDA0001332605030000077
采用现有的基于lq-范数的迭代重加权方法进行重构，包括如下步骤：
步骤1051：输入一组测量值
Figure GDA0001332605030000078
低阶高斯随机观测矩阵Φ(M/t)×(N/t)，原始信号长度N/t；
步骤1052：初始化参数：令w(0)＝(1,…,1)N/t×1，ε0＝1，
Figure GDA0001332605030000081
j＝1，2，…，N/t；
步骤1053：根据设定的w(n)构建大小为(N/t)×(N/t)的对角矩阵D(n)；
步骤1054：根据如下方式计算并更新
Figure GDA0001332605030000082
Figure GDA0001332605030000083
步骤1055：判断
Figure GDA0001332605030000084
是否小于
Figure GDA0001332605030000085
若满足，则执行步骤1056，否则跳过步骤1056；
步骤1056：更新εn值，εn+1＝ρεn，其中0<ρ<1；
步骤1057：按照如下方式计算并更新w(n)：
Figure GDA0001332605030000086
其中
Figure GDA0001332605030000087
表示D(n)中第p个对角分量(p＝1,2,…,N/t)；q值的取值范围为(0，1)；
Figure GDA0001332605030000088
为第n次迭代重构的结果；
步骤1058：更新迭代次数，并判断εn是否满足条件，若不满足，返回步骤1053，继续迭代计算；若满足则退出迭代，并输出一组稀疏解
Figure GDA0001332605030000089
步骤1059：依次对剩余的每组测量值
Figure GDA00013326050300000810
进行如步骤1051～1058的重构，输出全部t组重构稀疏系数。
步骤106：根据步骤105所得的每一组重构结果
Figure GDA00013326050300000811
进行如下方式的排序，得到一列稀疏信号的最终重构结果
Figure GDA00013326050300000812
Figure GDA00013326050300000813
式中
Figure GDA00013326050300000814
表示
Figure GDA00013326050300000815
中第(jt+1)至((j+1)t)个元素，其中j＝1，2，…，N/t；
Figure GDA00013326050300000816
为第i列重构稀疏解中的第j个元素；
步骤107：对剩余的N-1列测量值进行步骤104、105所述方法重构，并按照步骤106所述对每一列重构结果进行排序，最终得原始图像稀疏系数的重构结果
Figure GDA0001332605030000091
步骤108：对重构稀疏系数
Figure GDA0001332605030000092
进行小波逆变换，得到大小为N×N的重构图像。
对于一个长度为N的稀疏向量θN×1，采样数为M，若利用大小为M×N的观测矩阵对其进行采样，需要进行MN次乘法运算和M(N-1)次的加法运算，其相应的重构计算复杂度达到了
Figure GDA0001332605030000093
而利用大小为(M/t)×(N/t)的观测矩阵对θN×1进行采样，仅需(M/t)(N/t)次的乘法运算和(M/t)(N/t-1)次的加法运算，通过本发明所述的分块重构方法，其重构的计算复杂度仅约为
Figure GDA0001332605030000094
为此，本发明所述的一种快速低存储图像压缩感知方法，能够有效降低观测矩阵所需的存储空间，同时利用分块重构的方法也有益于降低重构的计算复杂度，提升重构的实时性。下面基于上述具体实施方式所述的方法，以实施例对本发明的技术效果进行进一步说明。方法的步骤如前，不再赘述，仅对具体的实现方法及效果进行展示。
实施例
为验证本发明所述方法的有效性，文中针对2维灰度图像进行了验证和比较。针对2为灰度图像信号，设计了三组验证实验，第1组设定不同的采样率，分别构建不同大小的高斯随机矩阵进行采样并利用lq-范数(0<q<1)的IRLS方法进行重构，并比较重构图像的峰值信噪比(PSNR)、结构相似度(SSIM)以及重构时间。第2组同样设定不同的采样率，利用OMP重构方法，验证比较利用本发明所述的分块重构算法对于图像重构实时性的提升性能。第3组设定不同的采样率，与BCS和Kronecker等的低存储压缩感知方法进行比较。
实验中采用的2维灰度图像分别为256×256的Peppers，512×512的Lena以及1024×1024的Mandrill。实验平台配备了Intel i7-4600CPU，2.1GHz主频，8GB内存，64位Windows 8操作系统，仿真软件采用Matlab R2010b。
重构算法采用了基于lq-范数(0<q<1)的IRLS算法，重构前，根据采样所得及观测矩阵参数t，将测量值(大小为M×N)进行分块，则每块大小为(M/t)×N，并利用lq-范数(0<q<1)的IRLS算法对每一块进行重构。重构过程中，令第i块第j列测量值为Yi,j，则每一列的重构解有如下表示：
Figure GDA0001332605030000095
式中
Figure GDA0001332605030000096
表示第i块第j列的第(n+1)次迭代结果，Dn为(N/t)×(N/t)大小的对角矩阵，Dn中第p个对角分量为
Figure GDA0001332605030000101
其中
Figure GDA0001332605030000102
为第n次迭代权重，定义如下：
Figure GDA0001332605030000103
从而，本文所涉及的lq-范数(0<q<1)的迭代重加权算法如下:
1)初始化:w(0)＝(1,…,1)，ε0＝1，
Figure GDA0001332605030000104
2)根据w(n)得到N×N对角矩阵Dn；
3)根据式(28)计算并更新
Figure GDA0001332605030000105
4)判断
Figure GDA0001332605030000106
是否小于
Figure GDA0001332605030000107
若满足，则执行第5)步，否则跳过5)；
5)更新εn值，εn+1＝ρεn(0<ρ<1)；
6)根据式(29)计算并更新w(n)值；
7)更新迭代次数，并判断εn是否满足条件，若不满足，返回第2)步，继续迭代计算；若满足则退出迭代，并输出稀疏解
Figure GDA0001332605030000108
实验时,设定算法的终止条件是εn<10-8，即在每次重构过程中，设ε0＝1，若经过若干次迭代后，有εn<10-8，则退出重构。完成每组重构任务后，执行步骤106所示操作，最终得到稀疏解
Figure GDA0001332605030000109
此外，对于lq-范数(0<q<1)的IRLS算法，现有报道中证明了在q值较大时，其恢复稀疏解的能力相对较强，因此，后续实验中均选择q＝0.8来验证半张量积压缩感知压缩采样方法对稀疏解的重构性能。
第1组实验中，首先对原始2维灰度图像进行小波变换，得到稀疏的小波系数；根据设定的采样率(M/N＝0.125，0.25，0.375，0.4375，0.5，0.5625，0.625，0.75)，构建6个不同大小的高斯随机观测矩阵Φ(t)(t＝1，2，4，8，16，32)，分别对原始图像的稀疏系数进行压缩采样，表1中列出了采样率为0.5时观测矩阵的大小。
表1:实验中所采用的不同大小的观测矩阵(M/N＝0.5)
Figure GDA00013326050300001010
如表1所示，当t＝1时，其压缩采样方式就退化为传统压缩感知的采样模式，故后续实验主要针对t＝1和t>1(如t＝2，4，8，16，32)时的重构性能比较。另如表1所示，当t＝32时，其观测矩阵仅为传统观测矩阵(t＝1时)大小的1/1024，大大降低了观测矩阵所需的存储空间。
得到对应的采样信号y后，利用本发明所述的lq-范数方法，对采样信号y进行重构，得到重构的小波系数；最后利用小波逆变换，得到重构图像。对于每幅灰度图像的稀疏系数，分别进行50次采样和重构，分别观测其重构图像和重构图像峰值信噪比、SSIM值及重构时间。图2-4列出了采样率为0.5时Peppers、Lena及Mandrill的重构结果，其重构图像为50次重构结果中随机选取的。
图2-4中的(a)-(b)分别表示了对2维原始图像的重构结果。从图2-4中可以较为直观地看到，对于三种不同大小的原始图像，当t＝2，4，8，16，32时，其相应的重构图像质量基本与传统压缩采样方法(t＝1)保持在同一个水平。为进一步对重构图像质量进行评估，表2-4分别列出了采用6种不同大小的高斯随机观测矩阵分别经50次采样和重构后，其重构图像峰值信噪比、结构相似度及重构时间的平均值。
根据表2和表3的PSNR和SSIM值所示，对于Peppers(256×256)而言，当观测矩阵减小至传统观测矩阵的1/64时(即t＝8时)，其重构图像质量仍与t＝1时保持同一水平，当持续减小观测矩阵后，其重构图像的PSNR和SSIM值稍有降低，但观测矩阵所需存储空间减小至传统的1/256或1/1024。然，对于Lena(512×512)和Mandrill(1024×1024)而言，随着观测矩阵的减小，其重构图像仍能与t＝1时保持一致水平，实验结果充分说明了可以采用低阶观测矩阵对原始图像信号进行压缩采样和精确重构，特别对于大尺寸图像，能达到更为高效的低存储化效果。需要指出的是，随着观测矩阵的持续减小，虽然Peppers重构质量的平均值略有下降，但这并不意味着观测矩阵的减小，其重构质量一定存在下降的趋势。在50次重构实验中，对于Peppers而言，仍存在若干次重构结果达到了t＝1时的质量水平，通过计算传感矩阵Acs(t)的相关性发现，若系统生成的低阶随机观测矩阵能够较好地满足RIP特性要求的，则依然能够得到较高的重构质量，反之则会造成重构质量的下降。
表2不同大小的高斯随机观测矩阵重构2维图像峰值信噪比PSNR(dB)
Figure GDA0001332605030000121
表3不同大小的高斯随机观测矩阵重构2维图像的SSIM值
Figure GDA0001332605030000122
表4列出了三种不同大小的原始图像经利用IRLS方法重构所需的时间。由表4可知，当t＝1时，对于全部采样率，Peppers重构所需的平均时间约为61.25s，当t＝2时，重构的平均时间约为17.875s，速度提升了3倍多；随着t继续增大，重构所需的时间持续降低，速度最快提升了11倍多(当t＝8时)。随着图像尺寸的增大，其重构时间呈指数型增加。当t＝1时，对于全部采样率，Lena重构所需的平均时间约为1579s，当t＝2时，重构的平均时间约为268s，速度提升了近6倍；随着t继续增大，重构所需的时间持续降低，速度最快提升了近69倍(当t＝16时，平均重构时间约为23s)。而对于Mandrill而言，当t＝1时，随着采样率的增加，其所需的重构时间是无法忍受的，而采用本发明所述分块重构方法，其重构所需的平均时间仅约为99s(当t＝32时)，其重构速度加快了近260倍，极大地提升了大尺寸图像重构的实时性。
由上所述，利用本发明所述的低阶压缩全局采样和分块重构的方法，在保证重构质量的前提下，即可成倍降低观测矩阵所需的存储空间，又可大大提升重构的实时性，这对于进一步拓展压缩感知的实际应用有非常积极的意义。
表4不同大小的高斯随机观测矩阵重构2维图像重构时间
Figure GDA0001332605030000131
为进一步检验本发明所述的分块重构方法对于重构实时性的提升作用，本发明进行了第2组实验。实验中再次利用低阶高斯随机观测矩阵对Lena小波系数进行压缩采样并进行分块处理后，利用OMP方法进行了重构，重构的结束条件为残差小于10-9，表5列出了采用6种不同大小的高斯随机观测矩阵分别经50次采样和重构后，其重构图像峰值信噪比及重构时间的平均值。
表5不同大小高斯观测矩阵采样后经OMP重构的性能(Lena 512×512)
Figure GDA0001332605030000141
对比表2和表4中Lena图像的重构时间，可以发现OMP算法的重构实时性要明显优于IRLS方法，随着观测矩阵的持续减小且当t＝32时，OMP算法重构Lena最短时间仅需1.46s(对应采样率为0.125时)；最长重构时间也仅为8.55s，比传统压缩采样的OMP重构速度提升了近145倍(对应采样率为0.75时)，大大提升了图像重构的实时性；此外，重构一幅大小为256×256的图像，其最短时间仅需0.5s左右(读者可以自行验证，选t＝16即可)，由此可知，采用本发明所述的低阶压缩采样和分块重构的方法，对于进一步拓展压缩感知理论的应用有非常积极的意义。
同时，对比表2和表5中Lena图像的PSNR值，虽然OMP方法的结果相比IRLS方法略低，但仍然能够得到较为稳定的重构质量(相比t＝1时)。
为进一步验证及比较本发明所述方法的重构性能，我们进行了第3组实验。对比实验主要与分块压缩感知(BCS)方法和Kronecker压缩感知方法进行了比较。BCS压缩采样时人为地将原始图像分成大小相等的若干块(如16×16)，随后根据设定的采样率，构建观测矩阵对原始图像进行分块采样和重构；Kronecker CS方法首先根据设定的采样率构建1个或若干个低阶观测矩阵，随后利用Kronecker运算构建一个大小与传统压缩感知方法一样的观测矩阵(M×N)进行压缩采样并进行全局重构；本发明所述的STP-CS方法根据设定的采样率首先构建一个较小的观测矩阵(t>1)对原始信号进行全局采样，随后进行分块重构。对于三种不同的低存储化压缩感知方法，从本质上来说就是采用的观测矩阵大小与压缩采样方式的区别。
实验中，对BCS方法设定图像块的大小为16×16，Kronecker方法中采用两个矩阵(一个大小为2×2，一个为M/2×N/2)进行张量积运算得到对应的观测矩阵，本发明方法中设定t＝2，4，8，16和32，采用的观测矩阵均为高斯随机观测矩阵，重构算法采用lq-范数(q＝0.8)的IRLS方法。对比实验中设定的采样率为M/N＝0.125，0.25，0.375，0.4375，0.5，0.5625，0.625及0.75。对Lena图像的小波系数分别进行50次采样和重构，取PSNR的平均值和重构时间的平均值。
表6与其他低存储压缩感知方法比较(Lena 512×512)
Figure GDA0001332605030000151
根据表6中PSNR值所示，对于相同的采样率，利用本发明所述方法，即便采用t＝32时的高斯随机观测矩阵，其所得重构图像的PSNR值均明显高于Kronecker方法所得结果。且相比BCS方法，本发明方法均能接近或超过BCS方法所得结果。此外，对于表6中所示的重构时间，本发明方法明显快于Kronecker方法，对于不同的采样率，本发明所述方法绝大部分的重构实时性要优于BCS方法，仅当采样率较低时(如0.125和0.25时)，相比BCS要略低。
此外，在对比实验中，对于BCS方法，其块大小的选择对重构质量有很大的影响，且对于Kronecker方法，其进行Kronecker运算的观测矩阵大小对重构质量也同样有着非常大的影响。而本发明所述方法，其选用的观测矩阵大小对重构质量的影响相比小了很多。
综上所述，本发明所述的STP-CS压缩感知模型可以利用低阶观测矩阵实现对原始信号的压缩采样，且在保证重构质量的前提下，极大地提升了重构的实时性，也成倍地降低了随机观测矩阵所需的存储空间。
以上所述的实施例只是本发明的一种较佳的方案，然其并非用以限制本发明。有关技术领域的普通技术人员，在不脱离本发明的精神和范围的情况下，还可以做出各种变化和变型。因此凡采取等同替换或等效变换的方式所获得的技术方案，均落在本发明的保护范围内。

Claims (1)

1.一种快速低存储图像压缩感知方法，其特征在于：至少包括以下步骤：
步骤101：输入大小为N×N的测试图像，对测试图像进行小波变换，得到小波系数θN×N；
步骤102：根据设定的采样率M/N，生成低阶正交随机高斯观测矩阵Φ(M/t)×(N/t)，其大小为
Figure FDA0002378988120000011
其中，N表示原始测试图像小波系数任意列向量的长度，M表示针对任意列的测量数，t为分组数，t>1，
Figure FDA0002378988120000012
为正整数；
步骤103：结合半张量积运算方法，利用低阶高斯随机观测矩阵Φ(M/t)×(N/t)对小波系数进行压缩采样，得到相应的采样值yM×N，压缩采样方法为：
yM×N＝Φ(M/t)×(N/t)⋉θN×N
其中符号“⋉”表示向量的左半张量积运算；
步骤104：对大小为M×N的测量值，对其中的任一列yN×1进行如下方式的分组处理，将一列测量值划分成t组：
Figure FDA0002378988120000013
其中
Figure FDA0002378988120000014
表示进行分组后的第i组测量值，是长度为M/t的列向量，i＝1，2,3，…，t；y(M/t-1)t+1为列向量中的元素；
步骤105：对于每一组测量值
Figure FDA0002378988120000015
利用低阶高斯随机观测矩阵Φ(M/t)×(N/t)，设定原始信号长度N/t，进行分组重构，并得到每一组重构结果
Figure FDA0002378988120000016
其中
Figure FDA0002378988120000017
步骤106：根据步骤105所得的每一组重构结果
Figure FDA0002378988120000018
进行如下方式的排序，得到一列稀疏信号的最终重构结果
Figure FDA0002378988120000019
Figure FDA0002378988120000021
式中
Figure FDA0002378988120000022
表示
Figure FDA0002378988120000023
中第(jt+1)至((j+1)t)个元素，其中j＝1，2，…，N/t；
Figure FDA0002378988120000024
为第i列重构稀疏解中的第j个元素；
步骤107：对剩余的N-1列测量值进行步骤104、105所述方法重构，并按照步骤106所述对每一列重构结果进行排序，最终得原始图像稀疏系数的重构结果
Figure FDA0002378988120000025
步骤108：对重构稀疏系数
Figure FDA0002378988120000026
进行小波逆变换，得到大小为N×N的重构图像；
所述步骤102中根据设定的采样率M/N，生成低阶正交随机高斯观测矩阵Φ(M/t)×(N/t)的具体步骤如下：
步骤1021：设定采样率M/N；
步骤1022：设定参数t，确定低阶观测矩阵大小(M/t)×(N/t)；
步骤1023：根据设定观测矩阵大小，系统生成正交高斯随机矩阵Φ(M/t)×(N/t)；
所述步骤103中的利用低阶正交高斯随机观测矩阵Φ(M/t)×(N/t)对原始信号的小波系数进行全局压缩采样，包括如下步骤：
步骤1031：输入低阶正交高斯随机观测矩阵Φ(M/t)×(N/t)和原始信号小波系数θN×N；
步骤1032：根据设定参数t，将小波系数θN×N的每一列θi N×1分割成N/t个等长的块：θi (1)，θi (2)，…，θi (N/t),所有块均为t维列向量，得到分块后的小波系数：
Figure FDA0002378988120000027
步骤1033：取分块后的一列小波系数θi (1)，θi (2)，…，θi (t)，利用已生成的低阶高斯随机观测矩阵Φ(M/t)×(N/t)与分块后的小波系数进行压缩采样，得到一列长度为M的采样值，压缩采样方法为：
Figure FDA0002378988120000031
其中φi,j∈Φ(M/t)×(N/t)，i＝1,2,…,M/t，j＝1,2,…,N/t；
步骤1034：利用低阶高斯随机观测矩阵完成每一列稀疏系数的压缩采样后，得到最后的压缩测量值yM×N；
所述步骤105中对于每一组测量值
Figure FDA0002378988120000032
采用现有的基于lq-范数的迭代重加权方法进行重构，包括如下步骤：
步骤1051：输入一组测量值
Figure FDA0002378988120000033
低阶高斯随机观测矩阵Φ(M/t)×(N/t)，原始信号长度N/t；
步骤1052：初始化参数：令w(0)＝(1,…,1)N/t×1，ε0＝1，
Figure FDA0002378988120000034
Figure FDA0002378988120000035
步骤1053：根据设定的w(n)构建大小为(N/t)×(N/t)的对角矩阵D(n)；
步骤1054：根据如下方式计算并更新
Figure FDA0002378988120000036
Figure FDA0002378988120000037
步骤1055：判断
Figure FDA0002378988120000038
是否小于
Figure FDA0002378988120000039
若满足，则执行步骤1056，否则跳过步骤1056；
步骤1056：更新εn值，εn+1＝ρεn，其中0<ρ<1；
步骤1057：按照如下方式计算并更新w(n)：
Figure FDA00023789881200000310
其中
Figure FDA00023789881200000311
表示D(n)中第p个对角分量(p＝1,2,…,N/t)；q值的取值范围为(0，1)；
Figure FDA0002378988120000041
为第n次迭代重构的结果；
步骤1058：更新迭代次数，并判断εn是否满足条件，若不满足，返回步骤1053，继续迭代计算；若满足则退出迭代，并输出一组稀疏解
Figure FDA0002378988120000042
步骤1059：依次对剩余的每组测量值
Figure FDA0002378988120000043
进行如步骤1051～1058的重构，输出全部t组重构稀疏系数。
CN201710380465.8A 2017-05-25 2017-05-25 一种快速低存储图像压缩感知方法 Active CN107170019B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710380465.8A CN107170019B (zh) 2017-05-25 2017-05-25 一种快速低存储图像压缩感知方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710380465.8A CN107170019B (zh) 2017-05-25 2017-05-25 一种快速低存储图像压缩感知方法

Publications (2)

Publication Number Publication Date
CN107170019A CN107170019A (zh) 2017-09-15
CN107170019B true CN107170019B (zh) 2020-10-23

Family

ID=59821860

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710380465.8A Active CN107170019B (zh) 2017-05-25 2017-05-25 一种快速低存储图像压缩感知方法

Country Status (1)

Country Link
CN (1) CN107170019B (zh)

Families Citing this family (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109194968B (zh) * 2018-09-13 2020-12-25 天津大学 一种融合信源信道译码的图像压缩感知方法
CN109559357B (zh) * 2018-09-17 2020-01-24 衡阳师范学院 一种基于小波包阈值的图像块压缩感知重构方法
CN109756740B (zh) * 2019-01-18 2020-02-18 北京邮电大学 基于最优测量矩阵的半张量图像压缩方法和图像恢复方法
CN110390622A (zh) * 2019-07-18 2019-10-29 江西财经大学 一种基于半张量积压缩感知的加密方法
CN111193925B (zh) * 2019-12-25 2021-12-07 杭州中威电子股份有限公司 一种基于分块向量内积的图像压缩感知编码及归一化方法
CN113378143B (zh) * 2021-07-06 2023-06-23 江西财经大学 基于半张量压缩感知的加密域可逆信息隐藏与认证方法
CN114125447B (zh) * 2021-11-02 2023-11-14 浙江树人学院（浙江树人大学） 一种基于分块及换位算法的压缩感知快速重构方法
CN114418866A (zh) * 2021-11-23 2022-04-29 北京邮电大学 一种低采样率下的深度学习图像恢复方法
CN117880434B (zh) * 2024-03-07 2024-05-24 江西财经大学 基于双2d压缩感知融合的可逆缩略图保留加密方法与系统

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2017039651A1 (en) * 2015-09-02 2017-03-09 Siemens Healthcare Gmbh Fast sparse computed tomography image reconstruction from few views

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2017039651A1 (en) * 2015-09-02 2017-03-09 Siemens Healthcare Gmbh Fast sparse computed tomography image reconstruction from few views

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
低存储化压缩感知;王金铭，叶时平，徐振宇，蒋燕君;《中国图像图形学报》;20160905;第21卷(第7期);第3页第1部分至第6页第5部分 *

Also Published As

Publication number Publication date
CN107170019A (zh) 2017-09-15

Similar Documents

Publication Publication Date Title
CN107170019B (zh) 一种快速低存储图像压缩感知方法
CN109671029B (zh) 基于伽马范数最小化的图像去噪方法
Li et al. A fast algorithm for learning overcomplete dictionary for sparse representation based on proximal operators
CN107146251B (zh) 一种三维人脸模型的对称性分析方法
CN105488767B (zh) 一种基于最小二乘优化的压缩感知图像快速重建方法
CN110569961A (zh) 一种神经网络训练方法、装置及终端设备
CN109840585A (zh) 一种面向稀疏二维卷积的运算方法和系统
CN109446478B (zh) 一种基于迭代和可重构方式的复协方差矩阵计算系统
CN105630741A (zh) 一种改进的按位替换法求矩阵逆矩阵模块
WO2020223850A1 (en) System and method for quantum circuit simulation
CN111401413A (zh) 一种基于优化理论的带规模约束的并行聚类方法
CN112288046B (zh) 一种用于神经网络的基于混合粒度的联合稀疏方法
Mohanty et al. Design and performance analysis of fixed-point jacobi svd algorithm on reconfigurable system
CN116170601B (zh) 基于四列列向量分块奇异值分解的图像压缩方法
CN115443657A (zh) 一种最近邻搜索方法、编码器、解码器及存储介质
CN106707243A (zh) 一种重构雷达信号的广义正则化正交匹配追踪方法
CN112087235A (zh) 基于伪逆感知字典的稀疏度自适应doa估计方法及系统
CN107154061A (zh) 一种分块压缩感知的正则化解码方法
CN109088796B (zh) 基于网络层析成像技术的网络流量矩阵预测方法
Wang et al. Hardware efficient architectures of improved Jacobi method to solve the eigen problem
CN108875547B (zh) 一种面向资源受限环境的目标跟踪方法及装置
CN114125447B (zh) 一种基于分块及换位算法的压缩感知快速重构方法
Li et al. High-performance tensor decoder on GPUs for wireless camera networks in IoT
CN115983343A (zh) 基于FPGA的YOLOv4卷积神经网络轻量化方法
GB2571818A (en) Selecting encoding options

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
CB02 Change of applicant information
CB02 Change of applicant information

Address after: 312030 Jiangxia Road, Yang Xun Qiao Town, Keqiao District, Shaoxing City, Zhejiang Province, No. 2016

Applicant after: Zhejiang Shuren University

Address before: Hangzhou City, Zhejiang province 310015 Gongshu District shurenuniversity Road No. 8

Applicant before: Zhejiang Shuren University

GR01 Patent grant
GR01 Patent grant