一种用于编码孔径快照成像光谱仪的系统响应标定方法

摘要

本发明公开的一种简便的系统响应标定方法，用于编码孔径快照成像光谱系统，属于计算摄像学领域。针对现有技术中编码孔径快照成像光谱系统响应函数标定方法费时、费力、成本高的技术问题，本发明公开的一种用于编码孔径快照成像光谱仪的系统响应标定方法，将响应函数的标定问题转化为线性方程组的求解问题，无需要引入额外的器件即能够较为精确地得到系统的响应函数，进而提高编码孔径快照成像光谱系统标定效率，节约标定成本。

说明书

技术领域

本发明涉及一种用于编码孔径快照成像光谱仪的、简便的系统响应标定方法，属于计算摄像学领域。

背景技术

光谱成像不同于传统的RGB成像或全色成像，其在更宽的光谱范围上对场景的光谱信号进行了更精细的采样，能够获取更为丰富的场景信息。光谱成像的这一特点使其在污水检测、植被研究、大气监测、医学诊断等领域相比传统成像技术具有得天独厚的优势，因此正在被越来越广泛地投入应用。

目前，光谱成像技术可以分为“完全采样”的传统光谱成像技术和“欠采样”的基于压缩感知的计算光谱成像技术。前者虽然能够以较高的精度获取空间-光谱信息，但是不可避免地要在空间、时间或者光谱维度上进行推扫。随着人们对获取动态场景的光谱信息的需求不断提升，传统的光谱成像技术显得越来越捉襟见肘。而基于压缩感知的计算光谱成像技术有望突破传统光谱成像技术中空间分辨率、光谱分辨率、时间分辨率不可兼得的困境，近年来已经受到越来越多的关注。

基于压缩感知的计算光谱成像技术通过对场景的空间信息和光谱信息进行调制，能够在单次曝光中获取全部的空间-光谱信息，后期利用特定的重构算法，就可以获得不同波段下场景的图像。与传统光谱成像技术相比，基于压缩感知的计算光谱成像技术能够胜任动态场景的观测。

Ashwin Wagadarikar等人设计的编码孔径快照成像光谱仪(Coded ApertureSnapshot Spectral Imager,CASSI)是最常用的计算光谱成像装置之一。在此基础上还衍生出了CASSI+RGB相机、CASSI+全色相机等变体。CASSI系统的核心部件包括一块编码板和一块色散棱镜。编码板对场景信息在空间维度上进行编码，色散棱镜把不同波长的光谱信号调制到不同的空间位置。精确地描述系统对空间-光谱信号的调制过程对后期的重构至关重要。

由于直接标定CASSI系统响应函数存在的困难，目前的工作大多假设系统响应为理想脉冲响应，仅仅考虑位置调制而忽略强度调制。这种不精确的描述在相当大的程度上损害了重构的精度。Henry Arguello等人试图用DMD(Digital Mirror Device)替代编码模板来对系统响应进行标定(详见Henry Arguello,Hoover Rueda,Yuehao Wu,DennisW.Prather,and Gonzalo R.Arce,"Higher-order computational model for codedaperture spectral imaging,"Appl.Opt.52,D12-D21(2013))，但是这种标定方法需要逐DMD单元进行操作，费时费力，并且DMD的引入极大地增加了系统的成本。

发明内容

针对现有技术中CASSI系统响应函数标定方法费时、费力、成本高的技术问题，本发明公开的一种用于编码孔径快照成像光谱仪的系统响应标定方法，将响应函数的标定问题转化为线性方程组的求解问题，无需要引入额外的器件即能够较为精确地得到系统的响应函数，进而提高CASSI系统标定效率，节约标定成本。

本发明公开的一种用于编码孔径快照成像光谱仪的系统响应标定方法，包括如下步骤：

步骤101：用不同波长的均匀单色光照射CASSI系统中的编码模板M，得到并记录不同波长的均匀单色光在探测器上成的像

CASSI系统由物镜、编码模板、中继镜、色散棱镜和灰度相机构成。目标场景为I(x,y,λ)，f(x,y)为探测器上成的像，公式(1)描述了CASSI系统的成像过程：

其中，h(x,y,λ)表示系统对波长为λ的单色光的响应，m(x,y)表示编码模板，取值为0、1。对(1)式进行离散化，得到：

其中，Ω代表探测器上(i,j)像素点对应的场景中的区域，其大小由响应函数决定；h(i,j,λ)为离散化的响应函数，描述场景中不同的点对探测器上一点响应强度的贡献。

用不同波长的均匀单色光照射CASSI系统中的编码模板M，则式(1)、式(2)中所述场景为均匀单色场景，场景I(x,y,λ)与探测器上的像f(x,y)由式(1)、式(2)联系。记录不同波长的均匀单色光照射下编码模板在探测器上成的像

步骤102：选取上一点(i,j)，找到其在编码模板M上对应的区域Ω，构造用于估计CASSI系统响应的线性方程如公式(3)所示，

m₁₁·h₁₁+m₁₂·h₁₂+…+m_kl·h_kl＝M_ij(3)

所述线性方程左边为区域Ω内各点的亮度的加权求和，加权求和权重即为待求的离散化后的系统响应；方程右边为(i,j)点的亮度。

步骤103：选取像上点(i,j)相邻的点(i′,j′)，找到相邻的点(i′,j′)在编码模板上对应的区域Ω′，所述区域Ω′通过对Ω进行平移得到，平移量即为(i′-i,j′-j)。构造用于估计CASSI系统响应的线性方程如公式(4)所示，

m′₁₁·h₁₁+m′₁₂·h₁₂+…+m′_kl·h_kl＝M_i′j′(4)

所述线性方程左边为区域Ω′内各点的亮度的加权求和，加权求和权重即为待求的离散化后的系统响应；方程右边为(i′,j′)点的亮度。

步骤104：重复步骤103，得到N个估计CASSI系统响应的线性方程。这N个线性方程的矩阵表示形式为：

AH＝Y(5)

其中，A为0-1系数矩阵，H为待求的权重，Y为上选取的各点的亮度。

为了保证解的鲁棒性，步骤104所述N应大于待求的离散化后的权重的个数。

步骤105：用最小二乘法求解步骤104得到的N个估计CASSI系统响应的线性方程得到离散形式的CASSI系统响应H，即实现点(i,j)的CASSI系统响应函数标定。

利用最小二乘法求解公式，得到H的最小二乘解为：

H＝(A^TA)^-1A^TY(4)

H即为离散形式的CASSI系统响应，即实现单点CASSI系统响应函数标定。

还包括步骤106：考虑系统响应的空间不一致性，则在像上选取不同的点(i,j)，重复上述步骤101-105，进行多点CASSI系统响应标定，得到不同点的CASSI系统响应，即完成不同点的CASSI系统响应函数标定。多点标定能够捕捉到CASSI系统响应的空间不一致性，进而能够更精确描述CASSI系统响应。

有益效果：

1、本发明公开的一种用于编码孔径快照成像光谱仪的系统响应标定方法，选取标定图像上的点及其在编码模板M上对应的区域，构造用于估计CASSI系统响应的线性方程，即将响应函数的标定问题转化为求解线性方程组的问题，无需要引入额外的器件即能够较为精确地得到系统的响应函数，进而提高CASSI系统响应函数标定效率，节约标定成本。

2、本发明公开的一种编码孔径光谱成像系统的重构方法，效率高、成本低，可以推广使用，进而能够实现精确地描述CASSI系统的调制过程，有利于光谱图像重构精度的提升。

附图说明

图1是本发明中提到的编码孔径快照成像光谱仪的系统结构图；

图2是本发明公开的一种用于编码孔径快照成像光谱仪的系统响应标定方法的流程图；

图3是本发明实施例中用到的编码模板；

图4是编码模板在450nm均匀单色光照射下在探测器上成的像；

图5是本发明实施例中利用最小二乘拟合得到的数值与真实值的对比；

图6是本发明实施例中考虑系统响应和不考虑系统响应对balloons进行仿真重构得到的结果(500nm)的对比，其中：图6-a为考虑系统响应的结果，图6-b为不考虑系统响应的结果；

图7是本发明实施例中考虑系统响应和不考虑系统响应对feathers进行仿真重构得到的结果(500nm)的对比，其中：图7-a为考虑系统响应的结果，图7-b为不考虑系统响应的结果；

图8是本发明实施例中考虑系统响应和不考虑系统响应对cd进行仿真重构得到的结果(500nm)的对比，其中：图8-a为考虑系统响应的结果，图8-b为不考虑系统响应的结果。

具体实施方式

为了更好地说明本发明的目的和优点，下面结合附图和实例对发明内容做进一步说明。

实施例1：

本实施例公开的一种用于编码孔径快照成像光谱仪的系统响应标定方法，应用于CASSI系统，将响应函数的标定问题转化为一个线性方程组的求解问题，最终通过计算最小二乘解的方式得到系统响应函数。本实施例的流程图如图2所示。

本实施例公开的一种用于编码孔径快照成像光谱仪的系统响应标定方法，包含以下步骤：

步骤101：分别用400nm,410nm,…,700nm波长的均匀单色光照射CASSI系统的编码模板，记录编码模板在探测器上成的像

步骤101中所述400nm-700nm均匀单色光由单色仪发出。成像过程中调整成像系统放大率，使编码模板上每个单元所成像的大小等于4×4个像素，随后将这16个像素合并为1个像素。编码模板为0、1编码模板。

步骤102：选取上一点(i,j)，找到其在编码模板M上对应的区域Ω，构造一个线性方程，如式(6)所示，

m₁₁·h₁₁+m₁₂·h₁₂+…+m₅₅·h₅₅＝M_ij(6)

所选区域Ω的大小为5×5个编码模板单元。所述线性方程左边为区域Ω内各点的亮度的加权求和，加权求和权重即为待求的离散化后的系统响应；方程右边为(i,j)点的亮度。

步骤103：选取上点(i,j)相邻的点(i′,j′)，找到其在编码模板上对应的区域Ω′，所述区域Ω′通过对Ω进行平移得到，平移量即为(i′-i,j′-j)。构造用于估计CASSI系统响应的线性方程如公式(7)所示，

m′₁₁·h₁₁+m′₁₂·h₁₂+…+m′₅₅·h₅₅＝M_i′j′(7)

所述线性方程左边为区域Ω′内各点的亮度的加权求和，加权求和权重即为待求的离散化后的系统响应；方程右边为(i′,j′)点的亮度。

步骤104：重复步骤103，得到49个方程。所选(i′,j′)为(i,j)最近邻的48个点。这49个线性方程的矩阵表示形式为：

AH＝Y(5)

其中，A为0-1系数矩阵，H为待求的权重，Y为上选取的各点的亮度。

步骤105：用最小二乘法求解步骤104得到的49个估计CASSI系统响应的线性方程得到离散形式的CASSI系统响应H，即实现点(i,j)的CASSI系统响应函数标定。

利用最小二乘法求解公式，得到H的最小二乘解为：

H＝(A^TA)^-1A^TY(4)

H即为离散形式的CASSI系统响应，即实现单点CASSI系统响应函数标定。

步骤106：在像上依次选取不同的点(i,j)，重复上述步骤101-105，完成所有点的响应函数标定。

本实施例将从两方面说明本发明的效果，一是系统响应估计的准确性，二是在描述CASSI系统调制过程时考虑系统响应后对重建精度带来的提升。

1.实验条件

本实验的硬件测试条件为：Inter i7 6700，8G RAM，Matlab 2016a。标定所用编码模板及其在探测器上成的像如图3、图4所示。仿真所用高光谱图片来自于哥伦比亚大学的Cave数据集。输入的压缩采样图像根据前面标定得到的系统响应生成，大小为512×542。用TwIST算法(详见Bioucas-Dias JM,Figueiredo MAT.A New TwIST:Two-Step IterativeShrinkage/Thresholding Algorithms for Image Restoration[J].IEEE Transactionson Image Processing.2007,16(12):2992-3004.)进行重构，其中采样矩阵分别为传统的假定系统响应为理想脉冲响应的采样矩阵和考虑系统响应的采样矩阵。

2.实验结果

为了验证系统响应估计的准确性，计算最小二乘的均方根误差(RMSE)。

归一化后计算得到的RMSE为0.0082。图5直观地展示了真实值与拟合值的差异。从RMSE以及图5来看，本发明能够较为精确地估计系统响应。

为了说明本发明对重建精度带来的有益提升，比较不同采样矩阵下TwIST的重建结果。几种指标如下表所示。

表1重建精度对比

图6、图7、图8分别为balloons、feathers、cd重建结果(500nm)对比。

从表1的结果以及图6、图7、图8可以看出，在考虑真实系统响应之后，重建精度有了较大的提升，说明了本发明能够带来的有益效果。

以上所述的具体描述，对发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。