CN110501072B - 一种基于张量低秩约束的快照式光谱成像系统的重构方法 - Google Patents

一种基于张量低秩约束的快照式光谱成像系统的重构方法 Download PDF

Info

Publication number
CN110501072B
CN110501072B CN201910793484.2A CN201910793484A CN110501072B CN 110501072 B CN110501072 B CN 110501072B CN 201910793484 A CN201910793484 A CN 201910793484A CN 110501072 B CN110501072 B CN 110501072B
Authority
CN
China
Prior art keywords
tensor
imaging system
hyperspectral image
dimensional
snapshot
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
CN201910793484.2A
Other languages
English (en)
Other versions
CN110501072A (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.)
Beijing Institute of Technology BIT
Original Assignee
Beijing Institute of Technology BIT
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 Beijing Institute of Technology BIT filed Critical Beijing Institute of Technology BIT
Priority to CN201910793484.2A priority Critical patent/CN110501072B/zh
Publication of CN110501072A publication Critical patent/CN110501072A/zh
Application granted granted Critical
Publication of CN110501072B publication Critical patent/CN110501072B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
    • G01J3/00—Spectrometry; Spectrophotometry; Monochromators; Measuring colours
    • G01J3/12—Generating the spectrum; Monochromators
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
    • G01J3/00—Spectrometry; Spectrophotometry; Monochromators; Measuring colours
    • G01J3/28—Investigating the spectrum
    • G01J3/2823—Imaging spectrometer
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
    • G01J3/00—Spectrometry; Spectrophotometry; Monochromators; Measuring colours
    • G01J3/28—Investigating the spectrum
    • G01J3/2823—Imaging spectrometer
    • G01J2003/2826—Multispectral imaging, e.g. filter imaging
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
    • G01J3/00—Spectrometry; Spectrophotometry; Monochromators; Measuring colours
    • G01J3/28—Investigating the spectrum
    • G01J2003/283—Investigating the spectrum computer-interfaced
    • G01J2003/284—Spectral construction

Abstract

本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，属于计算摄像学领域。本发明应用于编码孔径快照式光谱成像系统和基于全色相机的双相机光谱系统两种快照式光谱成像系统，首先利用高光谱图像的非局部相似性构建三维张量，然后使用维度可区分的张量低秩约束模型(DLTR,dimension‑discriminative low‑rank tensor regularization)挖掘三维张量的结构特性，包括空间自相似性、光谱相关性和空‑谱联合相关性；最后交替更新、迭代求解，从而完成高精度的高光谱图像重构。本发明能够更好地传递高光谱图像的高维物理特性，更好地挖掘其内在结构特性，大幅度提高快照式光谱成像系统的重构质量，具有重构精度高的优点。

Description

一种基于张量低秩约束的快照式光谱成像系统的重构方法
技术领域
本发明专利涉及一种用于快照式光谱成像系统的高光谱图像重构方法，尤其是能够高精度地重构高光谱图像的方法，属于计算摄像学领域。
背景技术
高光谱图像能够获取目标场景的连续光谱采样。与RGB彩色图像或者全色图像相比，高光谱图像能够提供丰富的光谱信息和细节结构。目前，高光谱图像已经被应用到环境监控、地质勘探、航空航天和医学检测等多个领域。
高光谱图像被称为三维数据立方体，包含二维的空间信息和一维的光谱信息。为了获取高光谱图像，传统成像光谱系统均采用扫描的方式，一次获取一个点、一条线或者一个谱带的光谱信息。这种方式成像速度慢，无法进行动态场景的光谱成像。近年来，随着计算成像技术和压缩感知理论的发展，快照式光谱成像系统打破了传统成像方式的壁垒，能够在远低于奈奎斯特采样数下完成目标场景的压缩光谱采样，并利用重构算法完成高光谱图像的获取。最具代表性的快照式光谱成像系统是杜克大学David Brady等人提出的编码孔径快照光谱成像系统(CASSI，coded aperture snapshot spectral imaging)。该系统通过对光谱信息编码、色散和积分等步骤，实现三维信息的二维压缩采样，能够通过单次曝光来获取目标场景的光谱信息，成像速度快。在此基础上，Wang等人提出的基于全色相机的双路系统(DCCHI，dual-camera compressive hyperspectral imaging)将CASSI和一个普通全色相机结合，实现光谱信息和空间信息的对称互补，提高重构质量。而如何从二维压缩采样中重构出三维高光谱图像是此类快照式光谱成像系统需要解决的重要问题。
从二维压缩采样中重构出三维高光谱图像是一个欠定难题。传统的重构算法从压缩感知原理出发，利用高光谱图像的先验信息构建带约束的目标方程，然后进行优化求解。目前，最常用的先验约束是稀疏约束，包括基于分段平滑的全变差约束法(Totalvariation，TV)，基于正交稀疏变换的梯度投影法(Gradient Projection for SparseReconstruction，GPSR)等都被用来重构高光谱图像。但此类算法容易引起图像过平滑现象，丢失了原有的细节纹理。此外，基于高光谱图像的空间-光谱相关性的矩阵低秩约束(Low-Rank Matrix Approximation，LRMA)也被用来重构高光谱图像，重构质量有一定提高。但这种算法与稀疏约束类重构算法相同，均把三维高光谱图像向量化成一维向量，破坏了高光谱图像本身的高维结构特性，重构质量仍不够理想，极大地限制了快照式光谱成像系统在实际中的应用。
发明内容
针对现有重构算法忽略高光谱图像高维结构特性、重构精度不理想的问题。本发明要解决的技术问题是提供一种基于张量低秩约束的快照式光谱成像系统的重构方法，能够更好地传递高光谱图像的高维物理特性，更好地挖掘其内在结构特性，大幅度提高快照式光谱成像系统的重构质量，具有重构精度高的优点。
为达到以上目的，本发明采用以下技术方案。
本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，应用于编码孔径快照式光谱成像系统和基于全色相机的双相机光谱系统两种快照式光谱成像系统，首先利用高光谱图像的非局部相似性构建三维张量，然后使用维度可区分的张量低秩约束模型(DLTR，dimension-discriminative low-rank tensor regularization)挖掘三维张量的结构特性，包括空间自相似性、光谱相关性和空-谱联合相关性；最后交替更新、迭代求解，从而完成高精度的高光谱图像重构。
本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，包括以下步骤：
步骤101：建立快照式光谱成像系统模型，并输入快照式光谱成像系统的采样图像Y、标定后的成像系统前向响应矩阵H、正则化系数τ、权重w、重构迭代次数Imax。
步骤101中所述快照式光谱成像系统是编码孔径快照式光谱成像系统(CASSI，coded aperture snapshot spectral imaging)和基于全色相机的双相机光谱系统(DCCHI，dual-camera compressive hyperspectral imaging)两种快照式光谱成像系统。CASSI系统主要由物镜、编码模版、中继镜、色散棱镜和全色相机等部件构成。目标场景的高光谱图像F大小为M×N×Ω，高光谱图像F上任意一点的像素值为f(i，j，λ)，且1≤i≤M，1≤j≤N，1≤λ≤Ω。其中，M×N表示高光谱图像的空间分辨率，Ω表示高光谱图像的谱段数。入射光经物镜到达编码模版进行随机的0-1编码。经编码后的图像到达色散棱镜后，不同频段的图像会沿着竖直方向进行偏移。最后所有频段的图像到达全色相机后进行叠加，得到压缩的二维混叠图像。CASSI系统的数学模型为：
Figure BDA0002178814280000021
其中ω(λ)表示CCD相机的频谱响应函数，Cu(i，j)表示编码模版函数，φ(λ)表示色散棱镜的波段偏移函数，gc(i，j)为二维混叠观测图像。将上述CASSI系统的成像响应模型写成矩阵形式为：
Gc＝HcF (2)
其中F为高光谱图像，大小为M×N×Ω；Hc表示CASSI系统的前向响应矩阵，大小为M(N+Ω-1)×MNΩ，包含ω(λ)、Cu(i，j)、φ(λ)以及积分的共同作用；Gc表示CASSI二维压缩观测图像，大小为M×(N+Ω-1)。
DCCHI系统由一个分光镜、一个CASSI系统和一个全色相机这两个分支构成。进入DCCHI系统的入射光先经过分光镜被一分为二，一半进入CASSI系统，一半直接到达全色相机。进入CASSI系统的成像过程如前所述。进入全色相机分支的入射光会直接到达灰度相机得到目标场景的二维灰度投影，DCCHI系统模型为：
Figure BDA0002178814280000031
将上式写成矩阵形式，得到全色相机分支的成像响应模型：
Gp＝HpF (4)
其中Hp表示全色相机分支的前向响应矩阵，大小为MN×MNΩ，包含ω(λ)以及积分的共同作用；Gp表示全色相机的观测图像，大小为M×N。将公式(2)和公式(4)联立，得到DCCHI系统的成像响应模型：
Figure BDA0002178814280000032
快照式光谱成像系统的成像响应模型如下：
G＝HF (6)
则对于CASSI系统，G＝Gc，H＝Hc；对于DCCHI系统，G＝[Gc；Gp]，H＝[Hc；Hp]。
并输入快照式光谱成像系统的采样图像Y、标定后的成像系统前向响应矩阵H、正则化系数τ、权重w、重构迭代次数Imax，用于后续步骤。
步骤102：初始化重构高光谱图像F0，初始化迭代次数t＝0。
步骤102所述重构高光谱图像F0的初始化方式为：
F0＝HTG (7)
其中HT表示快照式光谱成像系统的前向响应H的转置，即从二维压缩观测反投影到三维数据立方体的过程。
步骤103：对高光谱图像进行空间维的重叠块采样，空间块大小为s×s，步长为Δ，采样得到大小为s×s×Ω的三维立方体。三维立方体块的总数为L＝((M-s)/(s-Δ)+1)×((N-s)/(s-Δ)+1)。
步骤104：对步骤103中得到的每一个三维立方体，先将其空间维向量化，得到大小为s2×Ω的二维空-谱块。然后分别以每一个二维空-谱块为中心，在大小为W×W的窗口里，使用最近邻算法搜索k个与所述中心空-谱块最接近的空-谱块。
步骤104所述最近邻算法的具体内容为，先计算窗口内每一个空-谱块与中心空-谱块的欧几里得距离，然后将得到的所有距离进行从小到大排序，则前k个最小距离对应的空-谱块就是要搜索的k个最接近的空-谱块。
步骤105：将每一个空-谱块和其在步骤104中搜索得到的k个最接近的空-谱块进行数据整合，得到大小为s2×Ω×k的三维张量P，则该三维张量P能够体现高光谱图像的非局部相似性。用符号R表示步骤103至步骤105的采样、搜索和整合过程，则对于第l个(1≤l≤L)三维张量Pl表示为：
Pl＝RlF (8)
步骤105所述张量是矩阵的推广，维度大于等于三的数据均被称为张量。
步骤105所述高光谱图像的非局部相似性包含三个方面，分别是空间相似性、光谱相似性和联合相关性。其中，空间相似性是指图像在空间纹理上相似，光谱相关性是指光谱谱线相似，联合相关性则是指中心空-谱块与k个最接近的空-谱块的距离上非常接近，故其在空间维和光谱维均有较高的相关性。
步骤106：从l＝1到l＝L，即对每一个三维张量使用基于维度可区分的张量低秩约束模型进行张量低秩复原，获得去噪后的三维张量。
步骤106中所述维度可区分的张量低秩约束模型为：
Figure BDA0002178814280000041
其中，运算符
Figure BDA0002178814280000042
表示Frobenius范数的平方，τ为正则化系数，
Figure BDA0002178814280000043
表示对张量Pl沿第n个维度进行模展开以后的矩阵，
Figure BDA0002178814280000044
表示
Figure BDA0002178814280000045
的第r个奇异值，Dn表示张量Pl的第n个维度的长度，ε为一个正的小数。在公式(9)中，wn为平衡不同维度的低秩程度的系数。由于张量Pl在三个维度分别具有空间相似性、光谱相似性和联合相关性，因此Pl在三个维度上进行模展开后的矩阵均具有低秩特性，但其低秩程度有所不同。系数wn是用来衡量所述三维张量Pl的三个维度的低秩程度，能够更好地表达高光谱图像在不同维度上的物理属性的差异性，进而提高张量低秩复原的精度。
步骤106所述张量低秩复原方法的目的是利用张量低秩约束从带噪声的RlF中恢复出去噪后的
Figure BDA0002178814280000046
根据公式(9)的维度可区分的张量低秩约束模型，复原
Figure BDA0002178814280000047
的优化目标方程为：
Figure BDA0002178814280000048
采用逐维度更新的方式来求解上述优化问题。令
Figure BDA0002178814280000049
逐维度更新的优化目标方程为：
Figure BDA0002178814280000051
其中
Figure BDA0002178814280000052
上式中
Figure BDA0002178814280000053
的解为：
Figure BDA0002178814280000054
其中
Figure BDA0002178814280000055
Figure BDA0002178814280000056
为对矩阵
Figure BDA0002178814280000057
进行奇异值分解的结果，diag(·)表示以括号内元素构成的对角矩阵，Sα，ε(σ)表示奇异值收缩算子，定义为：
Figure BDA0002178814280000058
其中c0＝|σ|-ε，c1＝(c0)2-4(α-ε|σ|)。
最后，令
Figure BDA0002178814280000059
即完成张量
Figure BDA00021788142800000510
的低秩复原，获得去噪后的三维张量，其中运算符fold3(·)表示将括号内的矩阵沿第3个维度逆向变换至张量形式。
步骤107：使用步骤106计算得到的张量
Figure BDA00021788142800000511
更新高光谱图像F。
步骤107所述更新高光谱图像F的优化目标方程为：
Figure BDA00021788142800000512
公式(14)中高光谱图像F的解为：
Figure BDA00021788142800000513
由于矩阵H的规模很大，无法直接求其解析解，因此需利用共轭梯度下降法求高光谱图像F的近似解，从而完成高光谱图像F的更新。
步骤108：更新参数，当前迭代次数t＝t+1，并转至步骤103进行迭代，直至t＝Imax，即使用基于维度可区分的张量低秩约束模型进行张量低秩复原，然后更新高光谱图像进行迭代求解，实现更好地传递高光谱图像的高维物理特性，更好地挖掘其内在结构特性，大幅度提高快照式光谱成像系统的重构质量，从而完成快照式光谱成像系统的高光谱图像重构。
有益效果：
1、本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，使用张量来挖掘高光谱图像的高维结构特性，能够克服现有高光谱重构方法中使用向量化对其结构特性的破坏问题，进而提高高光谱图像的重构精度。
2、本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，使用维度可区分的张量低秩约束模型来同时表达高光谱图像在多个维度的相关性以及物理差异性，有利于提高高光谱图像的重构精度。
3、本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，使用基于维度可区分的张量低秩约束模型来更新并完成高光谱图像的重构，能够大幅度提高重构精度。
4、本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，使用交替更新、迭代求解的思路完成高光谱图像的重构，能够保证收敛性和鲁棒性。
5、本发明公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，能够应用到多个快照式光谱成像系统，包括编码孔径快照式光谱成像系统以及基于全色相机的双相机光谱成像系统，具有良好的扩展性。
附图说明
图1是本发明中用于编码孔径快照光谱成像系统的结构图；
图2是本发明中用于基于全色相机的双相机光谱成像系统的结构图；
图3是本发明公开的基于张量低秩约束的快照式光谱成像系统的重构方法的流程图；
图4是本发明和对比方法在两种快照式光谱成像系统下，对测试图像1进行仿真重构后波长为600nm时的结果图，其中图4(a)为参考图像，图4(b)、图4(d)、图4(f)和图4(h)分别为TV、GPSR、LRMA和本发明在CASSI系统下的重建结果，图4(c)、图4(e)、图4(g)和图4(i)分别为TV、GPSR、LRMA和本发明在DCCHI系统下的重建结果，
图5是本发明和对比方法在两种快照式光谱成像系统下，对测试图像9进行仿真重构后波长为600nm时的结果图，其中图5(a)为参考图像，图5(b)、图5(d)、图5(f)和图5(h)分别为TV、GPSR、LRMA和本发明在CASSI系统下的重建结果，图5(c)、图5(e)、图5(g)和图5(i)分别为TV、GPSR、LRMA和本发明在DCCHI系统下的重建结果，
具体实施方式
为了更好地说明本发明的目的和优点，下面结合附图和实例对发明内容做进一步说明。
实施例1：
本实施例公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，应用于编码孔径快照光谱成像系统(Coded Aperture Snapshot Spectral Imager，CASSI)(详见Wagadarikar A，John R，Willett R，Brady D.Single disperser design for codedaperture snapshot spectral imaging[J].Applied optics.2008,47(10)：B44-B51.)和基于全色相机的双相机光谱成像系统(DCCHI，dual-camera compressive hyperspectralimaging)(详见Wang L，Xiong Z，Gao D，et al.Dual-camera design for coded aperturesnapshot spectral imaging[J].Applied Optics，2015，54(4)：848-58.)。其中，CASSI系统使用编码孔径和色散介质来对三维高光谱图像分别进行空间和光谱调制，通过探测器获取二维混叠投影。DCCHI系统是在CASSI系统的基础上增加了一个全色相机分支得到。DCCHI在获取CASSI的二维混叠投影的同时，使用全色相机来获取目标场景的二维未混叠的全色投影，能够提高重构质量。这两种系统能够通过单次曝光来获取目标场景的光谱信息，成像速度快，近年来得到该领域广泛的关注和研究。
如何从二维压缩采样中重构出三维高光谱图像是快照式光谱成像系统需要解决的重要问题。传统的重构算法从压缩感知原理出发，利用高光谱图像的先验信息构建带约束的目标方程，然后进行优化求解。目前，最常用的先验约束是稀疏约束，包括基于分段平滑的全变差约束法(Total variation，TV)(详见Chambolle A.An Algorithm for TotalVariation Minimization and Applications[M].KluwerAcademic Publishers，2004.)，基于正交稀疏变换的梯度投影法(Gradient Projection for Sparse Reconstruction，GPSR)(详见M.A.T.Figueiredo，R.D.Nowak，and S.J.Wright.Gradient projection forsparse reconstruction：Application to compressed sensing and other inverseproblems.In IEEE Journal of Selected Topics in Signal Processing，pages 586-597，2007.)等都被用来重构高光谱图像。但此类算法容易引起图像过平滑现象，丢失了原有的细节纹理。此外，基于高光谱图像的空间-光谱相关性的矩阵低秩约束(Low-RankMatrix Approximation，LRMA)(详见Y.Fu，Y.Zheng，I.Sato，and Y.Sato.Exploitingspectral-spatial correlation for coded hyperspectral image restoration.InIEEE Conference on Computer Vision and Pattern Recognition，pages3727-3736，June 2016.)也被用来重构高光谱图像，重构质量有一定提高。但上述所有算法在处理的过程中，均把三维高光谱图像向量化成一维向量，破坏了高光谱图像本身的高维结构特性，重构质量仍不够理想，极大地限制了快照式光谱成像系统在实际中的应用。
针对现有重构算法忽略高光谱图像高维结构特性、重构精度不理想的问题。本实施例提供一种基于张量低秩约束的快照式光谱成像系统的重构方法，能够更好地传递高光谱图像的高维物理特性，更好地挖掘其内在结构特性，大幅度提高快照式光谱成像系统的重构质量，具有重构精度高的优点。
如图3所示，本实施例公开的一种基于张量低秩约束的快照式光谱成像系统的重构方法，具体实现步骤如下：
步骤101：建立快照式光谱成像系统模型，并输入快照式光谱成像系统的采样图像Y、标定后的成像系统前向响应矩阵H、正则化系数τ、权重w、重构迭代次数Imax。
步骤101中所述快照式光谱成像系统是编码孔径快照式光谱成像系统(CASSI，coded aperture snapshot spectral imaging)和基于全色相机的双相机光谱系统(DCCHI，dual-camera compressive hyperspectral imaging)两种快照式光谱成像系统。CASSI系统主要由物镜、编码模版、中继镜、色散棱镜和全色相机等部件构成。目标场景的高光谱图像F大小为M×N×Ω，高光谱图像F上任意一点的像素值为f(i，j，λ)，且1≤i≤M，1≤j≤N，1≤λ≤Ω。其中，M×N表示高光谱图像的空间分辨率，Ω表示高光谱图像的谱段数。入射光经物镜到达编码模版进行随机的0-1编码。经编码后的图像到达色散棱镜后，不同频段的图像会沿着竖直方向进行偏移。最后所有频段的图像到达全色相机后进行叠加，得到压缩的二维混叠图像。CASSI系统的数学模型为：
Figure BDA0002178814280000081
其中ω(λ)表示CCD相机的频谱响应函数，Cu(i，j)表示编码模版函数，φ(λ)表示色散棱镜的波段偏移函数，gc(i，j)为二维混叠观测图像。将上述CASSI系统的成像响应模型写成矩阵形式为：
Gc＝HcF (2)
其中F为高光谱图像，大小为M×N×Ω；Hc表示CASSI系统的前向响应矩阵，大小为M(N+Ω-1)×MNΩ，包含ω(λ)、Cu(i，j)、φ(λ)以及积分的共同作用；Gc表示CASSI二维压缩观测图像，大小为M×(N+Ω-1)。
DCCHI系统由一个分光镜、一个CASSI系统和一个全色相机这两个分支构成。进入DCCHI系统的入射光先经过分光镜被一分为二，一半进入CASSI系统，一半直接到达全色相机。进入CASSI系统的成像过程如前所述。进入全色相机分支的入射光会直接到达灰度相机得到目标场景的二维灰度投影，DCCHI系统模型为：
Figure BDA0002178814280000082
将上式写成矩阵形式，得到全色相机分支的成像响应模型：
Gp＝HpF (4)
其中Hp表示全色相机分支的前向响应矩阵，大小为MN×MNΩ，包含ω(λ)以及积分的共同作用；Gp表示全色相机的观测图像，大小为M×N。将公式(2)和公式(4)联立，得到DCCHI系统的成像响应模型：
Figure BDA0002178814280000083
快照式光谱成像系统的成像响应模型如下：
G＝HF (6)
则对于CASSI系统，G＝Gc，H＝Hc；对于DCCHI系统，G＝[Gc；Gp]，H＝[Hc；Hp]。
并输入快照式光谱成像系统的采样图像Y、标定后的成像系统前向响应矩阵H、正则化系数τ、权重w、重构迭代次数Imax，用于后续步骤。
步骤102：初始化重构高光谱图像F0，初始化迭代次数t＝0。
步骤102所述重构高光谱图像F0的初始化方式为：
F0＝HTG (7)
其中HT表示快照式光谱成像系统的前向响应H的转置，即从二维压缩观测反投影到三维数据立方体的过程。
步骤103：对高光谱图像进行空间维的重叠块采样，空间块大小为s×s，步长为Δ，采样得到大小为s×s×Ω的三维立方体。三维立方体块的总数为L＝((M-s)/(s-Δ)+1)×((N-s)/(s-Δ)+1)。
步骤104：对步骤103中得到的每一个三维立方体，先将其空间维向量化，得到大小为s2×Ω的二维空-谱块。然后分别以每一个二维空-谱块为中心，在大小为W×W的窗口里，使用最近邻算法搜索k个与所述中心空-谱块最接近的空-谱块。
步骤104所述最近邻算法的具体内容为，先计算窗口内每一个空-谱块与中心空-谱块的欧几里得距离，然后将得到的所有距离进行从小到大排序，则前k个最小距离对应的空-谱块就是要搜索的k个最接近的空-谱块。
步骤105：将每一个空-谱块和其在步骤104中搜索得到的k个最接近的空-谱块进行数据整合，得到大小为s2×Ω×k的三维张量P，则该三维张量P能够体现高光谱图像的非局部相似性。用符号R表示步骤103至步骤105的采样、搜索和整合过程，则对于第l个(1≤l≤L)三维张量Pl表示为：
Pl＝RlF (8)
步骤105所述张量是矩阵的推广，维度大于等于三的数据均被称为张量。
步骤105所述高光谱图像的非局部相似性包含三个方面，分别是空间相似性、光谱相似性和联合相关性。其中，空间相似性是指图像在空间纹理上相似，光谱相关性是指光谱谱线相似，联合相关性则是指中心空-谱块与k个最接近的空-谱块的距离上非常接近，故其在空间维和光谱维均有较高的相关性。
步骤106：从l＝1到l＝L，即对每一个三维张量使用基于维度可区分的张量低秩约束模型进行张量低秩复原，获得去噪后的三维张量。
步骤106中所述维度可区分的张量低秩约束模型为：
Figure BDA0002178814280000101
其中，运算符
Figure BDA0002178814280000102
表示Frobenius范数的平方，τ为正则化系数，
Figure BDA0002178814280000103
表示对张量Pl沿第n个维度进行模展开以后的矩阵，
Figure BDA0002178814280000104
表示
Figure BDA0002178814280000105
的第r个奇异值，Dn表示张量Pl的第n个维度的长度，ε为一个正的小数。在公式(9)中，wn为平衡不同维度的低秩程度的系数。由于张量Pl在三个维度分别具有空间相似性、光谱相似性和联合相关性，因此Pl在三个维度上进行模展开后的矩阵均具有低秩特性，但其低秩程度有所不同。系数wn是用来衡量所述三维张量Pl的三个维度的低秩程度，能够更好地表达高光谱图像在不同维度上的物理属性的差异性，进而提高张量低秩复原的精度。
步骤106所述张量低秩复原方法的目的是利用张量低秩约束从带噪声的RlF中恢复出去噪后的
Figure BDA0002178814280000106
根据公式(9)的维度可区分的张量低秩约束模型，复原
Figure BDA0002178814280000107
的优化目标方程为：
Figure BDA0002178814280000108
采用逐维度更新的方式来求解上述优化问题。令
Figure BDA0002178814280000109
逐维度更新的优化目标方程为：
Figure BDA00021788142800001010
其中
Figure BDA00021788142800001011
上式中
Figure BDA00021788142800001012
的解为：
Figure BDA00021788142800001013
其中
Figure BDA00021788142800001014
Figure BDA00021788142800001015
为对矩阵
Figure BDA00021788142800001016
进行奇异值分解的结果，diag(·)表示以括号内元素构成的对角矩阵，Sα，ε(σ)表示奇异值收缩算子，定义为：
Figure BDA00021788142800001017
其中c0＝|σ|-ε，c1＝(c0)2-4(α-ε|σ|)。
最后，令
Figure BDA00021788142800001018
即完成张量
Figure BDA00021788142800001019
的低秩复原，获得去噪后的三维张量，其中运算符fold3(·)表示将括号内的矩阵沿第3个维度逆向变换至张量形式。
步骤107：使用步骤106计算得到的张量
Figure BDA0002178814280000111
更新高光谱图像F。
步骤107所述更新高光谱图像F的优化目标方程为：
Figure BDA0002178814280000112
公式(14)中高光谱图像F的解为：
Figure BDA0002178814280000113
由于矩阵H的规模很大，无法直接求其解析解，因此需利用共轭梯度下降法求高光谱图像F的近似解，从而完成高光谱图像F的更新。
步骤108：更新参数，当前迭代次数t＝t+1，并转至步骤103进行迭代，直至t＝Imax，即使用基于维度可区分的张量低秩约束模型进行张量低秩复原，然后更新高光谱图像进行迭代求解，实现更好地传递高光谱图像的高维物理特性，更好地挖掘其内在结构特性，大幅度提高快照式光谱成像系统的重构质量，从而完成快照式光谱成像系统的高光谱图像重构。
为说明本发明的效果，本实施例在两种快照式光谱成像系统上仿真实验并进行对比分析。
1.实验条件
本实验的硬件测试条件为：Inter i7 6700，内存16G，Matlab 2015b。测试所用高光谱图片来自于CAVE数据集(详见F.Yasuma，T.Mitsunaga，D.Iso，andS.K.Nayar.Generalized assorted pixel camera：postcapture control ofresolution，dynamic range，and spectrum.IEEE Transactions on Image Processing，19(9)：2241-53，2010.)。CASSI中编码孔径模版为p＝0.5的贝努利随机矩阵，色散棱镜的色散为线性等距色散。对比方法为基于全变差约束的TV算法、基于稀疏梯度投影的GPSR算法和基于矩阵低秩约束的LAMA算法。实验参数设置为τ＝1，w＝[0.1，0.1，1]。为了定量的衡量重构结果的质量，使用峰值信噪比(Peak signal to noise ratio，PSNR)衡量重构结果的空间质量，使用相对无量纲全局误差(erreurrelative globale adimensionnelledesynth`ese，ERGAS)(详见L.Wald.Data fusion：definitions and architectures：fusionof images of different spatial resolutions.Presses des MINES，2002.)衡量重构结果的光谱保真度。其中，PSNR值越大，ERGAS值越小，重构精度越高。
2.实验结果
为了验证本发明在高光谱重构精度上的提升，在CASSI系统和DCCHI系统两个快照式光谱成像系统上测试本发明公开的方法和对比方法的重构精度。表1和表2分别是CASSI系统下和DCCHI系统下，对CAVE数据集上10幅图像进行重构的结果。
表1 CASSI系统下的重构精度对比
Figure BDA0002178814280000121
表2 DCCHI系统下的重构精度对比
Figure BDA0002178814280000122
从表1和表2的结果可以看出，TV算法和GPSR算法均是基于稀疏先验，获得的重构结果也非常相近。LRMA算法是基于矩阵低秩约束，能够获得较好的重构质量，说明低秩约束比稀疏约束更加利于高光谱图像的重构。而本发明公开的基于张量低秩约束的重构算法能够获得最好的重构结果，空间质量和光谱精度均明显优于对比算法，从而说明张量低秩约束重构方法的有效性。
图4为图像1在两种快照式光谱成像系统下的600nm处的重构结果，图5则为图像9在两种快照式光谱成像系统下的600nm处的重构结果。从图4和图5的结果均可以看出，TV算法和GPSR算法的重构图像表现的过于平滑，丢失细节，且引入很多噪声。而本发明公开方法的重构结果最好，空间细节和纹理信息都得到非常准确的复原，再次说明本发明公开的重构方法的有效性。
以上所述的具体描述，对发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (8)

1.一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：包括以下步骤，
步骤101：建立快照式光谱成像系统模型，并输入快照式光谱成像系统的采样图像Y、标定后的成像系统前向响应矩阵H、正则化系数τ、权重w、重构迭代次数Imax；
步骤102：初始化重构高光谱图像F0，初始化迭代次数t＝0；
步骤103：对高光谱图像进行空间维的重叠块采样，空间块大小为s×s，步长为Δ，采样得到大小为s×s×Ω的三维立方体；三维立方体块的总数为L＝((M-s)/(s-Δ)+1)×((N-s)/(s-Δ)+1)；其中F为高光谱图像，大小为M×N×Ω；
步骤104：对步骤103中得到的每一个三维立方体，先将其空间维向量化，得到大小为s2×Ω的二维空-谱块；然后分别以每一个二维空-谱块为中心，在大小为W×W的窗口里，使用最近邻算法搜索k个与中心空-谱块最接近的空-谱块；
步骤105：将每一个空-谱块和其在步骤104中搜索得到的k个最接近的空-谱块进行数据整合，得到大小为s2×Ω×k的三维张量P，则该三维张量P能够体现高光谱图像的非局部相似性；用符号R表示步骤103至步骤105的采样、搜索和整合过程，则对于第l个(1≤l≤L)三维张量Pl表示为：
Pl＝RlF (8)
步骤106：从l＝1到l＝L，即对每一个三维张量使用基于维度可区分的张量低秩约束模型进行张量低秩复原，获得去噪后的三维张量；
步骤107：使用步骤106计算得到的张量
Figure FDA0002445012220000011
更新高光谱图像F；
步骤108：更新参数，当前迭代次数t＝t+1，并转至步骤103进行迭代，直至t＝Imax，即使用基于维度可区分的张量低秩约束模型进行张量低秩复原，然后更新高光谱图像进行迭代求解，实现更好地传递高光谱图像的高维物理特性，更好地挖掘其内在结构特性，大幅度提高快照式光谱成像系统的重构质量，从而完成快照式光谱成像系统的高光谱图像重构。
2.如权利要求1所述的一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：步骤101中所述快照式光谱成像系统是编码孔径快照式光谱成像系统(CASSI，coded aperture snapshot spectral imaging)和基于全色相机的双相机光谱系统(DCCHI,dual-camera compressive hyperspectral imaging)两种快照式光谱成像系统；CASSI系统主要由物镜、编码模版、中继镜、色散棱镜和全色相机等部件构成；目标场景的高光谱图像F大小为M×N×Ω，高光谱图像F上任意一点的像素值为f(i,j,λ)，且1≤i≤M，1≤j≤N，1≤λ≤Ω；其中，M×N表示高光谱图像的空间分辨率，Ω表示高光谱图像的谱段数；入射光经物镜到达编码模版进行随机的0-1编码；经编码后的图像到达色散棱镜后，不同频段的图像会沿着竖直方向进行偏移；最后所有频段的图像到达全色相机后进行叠加，得到压缩的二维混叠图像；CASSI系统的数学模型为：
Figure FDA0002445012220000021
其中ω(λ)表示CCD相机的频谱响应函数，Cu(i,j)表示编码模版函数，φ(λ)表示色散棱镜的波段偏移函数，gc(i,j)为二维混叠观测图像；将上述CASSI系统的成像响应模型写成矩阵形式为：
Gc＝HcF (2)
其中F为高光谱图像，大小为M×N×Ω；Hc表示CASSI系统的前向响应矩阵，大小为M(N+Ω-1)×MNΩ，包含ω(λ)、Cu(i,j)、φ(λ)以及积分的共同作用；Gc表示CASSI二维压缩观测图像，大小为M×(N+Ω-1)；
DCCHI系统由一个分光镜、一个CASSI系统和一个全色相机这两个分支构成；进入DCCHI系统的入射光先经过分光镜被一分为二，一半进入CASSI系统，一半直接到达全色相机；进入CASSI系统的成像过程如前所述；进入全色相机分支的入射光会直接到达灰度相机得到目标场景的二维灰度投影，DCCHI系统模型为：
Figure FDA0002445012220000022
将上式写成矩阵形式，得到全色相机分支的成像响应模型：
Gp＝HpF (4)
其中Hp表示全色相机分支的前向响应矩阵，大小为MN×MNΩ，包含ω(λ)以及积分的共同作用；Gp表示全色相机的观测图像，大小为M×N；将公式(2)和公式(4)联立，得到DCCHI系统的成像响应模型：
Figure FDA0002445012220000023
快照式光谱成像系统的成像响应模型如下：
G＝HF (6)
则对于CASSI系统，G＝Gc，H＝Hc；对于DCCHI系统，G＝[Gc；Gp]，H＝[Hc；Hp]；
并输入快照式光谱成像系统的采样图像Y、标定后的成像系统前向响应矩阵H、正则化系数τ、权重w、重构迭代次数Imax，用于后续步骤。
3.如权利要求2所述的一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：步骤102所述重构高光谱图像F0的初始化方式为：
F0＝HTG (7)
其中HT表示快照式光谱成像系统的前向响应H的转置，即从二维压缩观测反投影到三维数据立方体的过程。
4.如权利要求3所述的一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：步骤104所述最近邻算法的具体实现方法为，先计算窗口内每一个空-谱块与中心空-谱块的欧几里得距离，然后将得到的所有距离进行从小到大排序，则前k个最小距离对应的空-谱块就是要搜索的k个最接近的空-谱块。
5.如权利要求4所述的一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：步骤105所述张量是矩阵的推广，维度大于等于三的数据均被称为张量；
步骤105所述高光谱图像的非局部相似性包含三个方面，分别是空间相似性、光谱相似性和联合相关性；其中，空间相似性是指图像在空间纹理上相似，光谱相关性是指光谱谱线相似，联合相关性则是指中心空-谱块与k个最接近的空-谱块的距离上非常接近，故其在空间维和光谱维均有较高的相关性。
6.如权利要求5所述的一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：步骤106中所述维度可区分的张量低秩约束模型为：
Figure FDA0002445012220000031
其中，运算符
Figure FDA0002445012220000032
表示Frobenius范数的平方,τ为正则化系数，
Figure FDA0002445012220000033
表示对张量Pl沿第n个维度进行模展开以后的矩阵，
Figure FDA0002445012220000034
表示
Figure FDA0002445012220000035
的第r个奇异值，Dn表示张量Pl的第n个维度的长度，ε为一个正的小数；在公式(9)中，wn为平衡不同维度的低秩程度的系数；由于张量Pl在三个维度分别具有空间相似性、光谱相似性和联合相关性，因此Pl在三个维度上进行模展开后的矩阵均具有低秩特性，但其低秩程度有所不同；系数wn是用来衡量所述三维张量Pl的三个维度的低秩程度，能够更好地表达高光谱图像在不同维度上的物理属性的差异性，进而提高张量低秩复原的精度。
7.如权利要求6所述的一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：步骤106所述张量低秩复原方法的目的是利用张量低秩约束从带噪声的RlF中恢复出去噪后的
Figure FDA0002445012220000036
根据公式(9)的维度可区分的张量低秩约束模型，复原
Figure FDA0002445012220000037
的优化目标方程为：
Figure FDA0002445012220000038
采用逐维度更新的方式来求解优化问题；令
Figure FDA0002445012220000039
逐维度更新的优化目标方程为：
Figure FDA00024450122200000310
其中
Figure FDA0002445012220000041
上式中
Figure FDA0002445012220000042
的解为：
Figure FDA0002445012220000043
其中
Figure FDA0002445012220000044
Figure FDA0002445012220000045
为对矩阵
Figure FDA0002445012220000046
进行奇异值分解的结果，diag(·)表示以括号内元素构成的对角矩阵，Sα,ε(σ)表示奇异值收缩算子，定义为：
Figure FDA0002445012220000047
其中c0＝|σ|-ε，c1＝(c0)2-4(α-ε|σ|)；
最后，令
Figure FDA0002445012220000048
即完成张量
Figure FDA0002445012220000049
的低秩复原，获得去噪后的三维张量，其中运算符fold3(·)表示将括号内的矩阵沿第3个维度逆向变换至张量形式。
8.如权利要求7所述的一种基于张量低秩约束的快照式光谱成像系统的重构方法，其特征在于：步骤107所述更新高光谱图像F的优化目标方程为：
Figure FDA00024450122200000410
公式(14)中高光谱图像F的解为：
Figure FDA00024450122200000411
由于矩阵H的规模很大，无法直接求其解析解，因此需利用共轭梯度下降法求高光谱图像F的近似解，从而完成高光谱图像F的更新。
CN201910793484.2A 2019-08-26 2019-08-26 一种基于张量低秩约束的快照式光谱成像系统的重构方法 Active CN110501072B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910793484.2A CN110501072B (zh) 2019-08-26 2019-08-26 一种基于张量低秩约束的快照式光谱成像系统的重构方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910793484.2A CN110501072B (zh) 2019-08-26 2019-08-26 一种基于张量低秩约束的快照式光谱成像系统的重构方法

Publications (2)

Publication Number Publication Date
CN110501072A CN110501072A (zh) 2019-11-26
CN110501072B true CN110501072B (zh) 2020-07-24

Family

ID=68589721

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910793484.2A Active CN110501072B (zh) 2019-08-26 2019-08-26 一种基于张量低秩约束的快照式光谱成像系统的重构方法

Country Status (1)

Country Link
CN (1) CN110501072B (zh)

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110926611A (zh) * 2019-12-12 2020-03-27 北京理工大学 一种应用于压缩感知光谱成像系统的噪声抑制方法
CN111174912B (zh) * 2020-01-03 2021-02-23 南京大学 一种快照型解色散模糊的高光谱成像方法
CN111397733B (zh) * 2020-04-23 2021-03-02 湖南大学 一种单/多帧快照式光谱成像方法、系统及介质
CN112241937B (zh) * 2020-07-22 2023-10-13 西安电子科技大学 一种基于神经网络的高光谱图像重构方法
CN112556848B (zh) * 2020-11-28 2022-09-06 南京理工大学 双相机压缩测量协同张量表示的融合计算成像方法
CN112989593B (zh) * 2021-03-09 2022-08-23 南京理工大学 基于双相机的高光谱低秩张量融合计算成像方法
CN115791628A (zh) * 2021-09-10 2023-03-14 北京与光科技有限公司 光谱恢复方法、光谱恢复装置和电子设备
CN114119426B (zh) * 2022-01-26 2022-07-01 之江实验室 非局部低秩转换域与全连接张量分解图像重构方法及装置
CN115950837B (zh) * 2023-03-10 2023-05-23 湖南师范大学 基于即插即用先验的快照式光谱成像方法、系统及介质
CN116563174B (zh) * 2023-07-11 2023-09-29 江西师范大学 一种图像的重建方法、装置和计算机存储介质

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107764797A (zh) * 2017-09-21 2018-03-06 天津大学 一种基于低秩张量算法的拉曼光谱图像数据预处理方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10436871B2 (en) * 2017-04-24 2019-10-08 Cedars-Sinai Medical Center Low-rank tensor imaging for multidimensional cardiovascular MRI

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107764797A (zh) * 2017-09-21 2018-03-06 天津大学 一种基于低秩张量算法的拉曼光谱图像数据预处理方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于低秩张量分析的高光谱图像降维与分类;陈昭等;《红外与毫米波学报》;20131231;全文 *

Also Published As

Publication number Publication date
CN110501072A (zh) 2019-11-26

Similar Documents

Publication Publication Date Title
CN110501072B (zh) 一种基于张量低秩约束的快照式光谱成像系统的重构方法
Cai et al. Mst++: Multi-stage spectral-wise transformer for efficient spectral reconstruction
CN111260576B (zh) 一种基于去噪三维卷积自编码网络的高光谱解混算法
CN107525588B (zh) 一种基于gpu的双相机光谱成像系统的快速重构方法
Rao et al. A residual convolutional neural network for pan-shaprening
CN109697697B (zh) 基于优化启发的神经网络的光谱成像系统的重构方法
CN109146787B (zh) 一种基于插值的双相机光谱成像系统的实时重建方法
CN113160294B (zh) 图像场景深度的估计方法、装置、终端设备和存储介质
CN109883548B (zh) 基于优化启发的神经网络的光谱成像系统的编码优化方法
CN111080567A (zh) 基于多尺度动态卷积神经网络的遥感图像融合方法及系统
He et al. Fast hyperspectral image recovery of dual-camera compressive hyperspectral imaging via non-iterative subspace-based fusion
CN114119444B (zh) 一种基于深度神经网络的多源遥感图像融合方法
CN105761234A (zh) 一种基于结构稀疏表示的遥感影像融合方法
CN109886898B (zh) 基于优化启发的神经网络的光谱成像系统的成像方法
CN111174912B (zh) 一种快照型解色散模糊的高光谱成像方法
CN109887050B (zh) 一种基于自适应字典学习的编码孔径光谱成像方法
CN108288256A (zh) 一种多光谱马赛克图像复原方法
CN110717947B (zh) 一种基于外部和内部训练的高质量光谱重构方法
CN109087262B (zh) 一种多视图光谱图像的重建方法、存储介质
CN109978802A (zh) 基于nsct和pcnn的压缩感知域内的高动态范围图像融合方法
CN110926611A (zh) 一种应用于压缩感知光谱成像系统的噪声抑制方法
CN112989593B (zh) 基于双相机的高光谱低秩张量融合计算成像方法
Zeng et al. U-net-based multispectral image generation from an rgb image
CN111397733B (zh) 一种单/多帧快照式光谱成像方法、系统及介质
CN112556848B (zh) 双相机压缩测量协同张量表示的融合计算成像方法

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
GR01 Patent grant
GR01 Patent grant