[发明专利]基于空谱三维稀疏先验的高光谱解混压缩感知方法

权利要求书

1.一种基于空谱三维稀疏先验的高光谱解混压缩感知方法，其特征在于包括以下步骤：

步骤一、对于包含n_p个像素，n_b个波段的高光谱图像其中每一个像素的光谱x_i表示成所有端元的线性组合，其中，为对应的丰度值向量；X表示成丰度值矩阵和端元矩阵W的乘积，如下：

X＝WH               (1)

H中，行方向是空间维，代表不同像素的光谱在同一个端元上的投影；列方向是光谱维，代表某个像素的光谱在不同端元上的投影；

步骤二、使用满足高斯随机分布，列归一化的随机采样矩阵对场景的光谱进行随机采样，得到压缩数据其中，m表示n_b个波段压缩后的长度，m＜n_b；

F＝AX＝AWH                   (2)

步骤三、面积有限的成像场景，通常仅包含有限个端元；根据场景信息在USGS光谱库中抽取n_e个光谱组成端元矩阵W；

步骤四、为了通过模型(2)精确求解H，综合考虑H空间维和光谱维上的稀疏性，设计了空谱三维稀疏先验，约束H的求解；具体过程如下：

(a)构建H的空谱三维稀疏先验：

minHλ1||H||1,1+λ2||H||TV,1---(3)]]>

其中，λ₁和λ₂表示两种约束的比重，λ₁＝λ₂＝1；‖·‖₁表示l₁范数，||H||TV,1=Σj=1ne||D(h‾j)||1=Σj=1ne||D(ϵjTH)||1,]]>e_j和ε_j分别表示和空间中的第j个单位向量，表示H中的第j行，D_i表示对信号第i个分量进行梯度运算的线性操作符；

(b)构建丰度值的其他先验；根据线性混合模型，引入丰度值矩阵非负且光谱维全和为1的限制，其中，和是全部元素为1，长度分别为n_e和n_p的向量；仅引入H＞0的限制；

(c)构建H的重建模型；结合式(2)、式(3)及H＞0得到H的重建模型：

H=argminHλ1||H||1,1+λ2||H||TV,1,s.t.F=AWH,H>0---(4)]]>

为了消除约束条件和数据保真项(2)之间的耦合作用，简化后续的优化过程，向式(4)中引入分离变量v_j＝He_j,得到：

minH,υij,uijΣj=1npλ1||vj||1+Σj=1neΣi=1npλ2|uij|,---(5)]]>

s.t.vj=Hej,∀j;uij=Di(ϵjTH),∀i,j;F=AWH,H>0]]>

(d)求解模型(5)得到丰度值矩阵H的估计具体求解过程如下：

①根据模型(5)构建关于H,v_j,u_ij的增广拉格朗日方程

其中，α＝2⁵,κ＝2⁵,β＝2¹³,γ＝2⁵为二次项惩罚系数，ψ_j,π_ij,Π,υ为对应的拉格朗日乘子，乘子的所有元素初始为0，‖·‖₂表示l₂范数，初始化外循环计数器C₀＝1；

②固定拉格朗日乘子和H，更新分离变量v_j,u_ij；形式如下：

vj=max{|Hej-ψiα|-λ1α,0}sgn(Hej-ψjα)]]>

                                  (7)

uij=max{|Di(ϵjTH)-πijκ|-λ2κ,0}sgn[Di(ϵjTH)-πijκ]]]>

③固定拉格朗日乘子和分离变量v_j,u_ij，采用梯度下降法更新H；将当前拉格朗日乘子下第一次H设为H⁰；假定第k+1次更新，由H^k得到H^k+1，形式如下：

其中，为关于H一阶导数，v_j,u_ij已知，形式如下：

τ为梯度下降步长；在初始化过程中，当第一次更新H⁰时，τ采用最速下降法进行初始化；之后更新H^k,k＝1,2,...时，对τ采用两点步长梯度法进行初始化；两点步长梯度法具体形式如下：

其中，tr(·)表示矩阵的迹；τ的细化过程具体如下：

(a)代入初始化的τ，根据式(8)得到H^k+1，设置参数δ＝3.2×10^-4，η＝0.6和计数器c＝0；

(b)判断H^k+1是否满足如下的条件：

如果不满足，更新计数器c＝c+1；

如果c＜5，缩小步长τ＝τ·η，继续循环判断是否满足式(11)；

否则τ由最速下降法确定，然后由式(8)得到更新的H^k+1；

否则，得到更新的H^k+1；

④循环②、③步直至‖H^k+1-H^k‖₂/‖H^k‖₂≤10^-3，得到当前固定的拉格朗日乘子下最优的丰度值矩阵H′；

如果‖H′-H⁰‖₂/‖H⁰‖₂≤10^-3或者，外循环超过300次，即C₀＞300，则最终估计的丰度值矩阵结束循环；

否则，更新拉格朗日乘子：

ψjk+1+ψjk-α(Hej-vj),πijk+1=πijk-κ[Di(ϵjTH)-uij],---(12)]]>

Π^k+1＝Π^k-β(AWH-F)

外循环计数器C₀＝C₀+1，重新设置H⁰＝H′,k＝0并循环②、③和④步；

步骤五、依据线性混合模型(1)，引入选择的端元矩阵W，得到重建的光谱数据

X^=WH^---(13).]]>