CN106815817A - 一种改进的高光谱图像去噪方法 - Google Patents

Abstract

本发明公开了一种改进的高光谱图像去噪方法，包括以下步骤：(1)将待去噪的高光谱图像转换成空间光谱联合的二维矩阵；(2)根据空间光谱联合的二维矩阵，采用基于欧式距离的近邻相似度计算策略计算高光谱图像的像素点与近邻的空间相似度；(3)结合像素点空间局域相似性和光谱间低秩性建立去噪模型，恢复出原始无噪数据；(4)采用原始无噪数据恢复出三维无噪高光谱图像。本发明较为显著地提高了去噪的效果,可以针对实时响应和高精度不同需求的应用场景给出可靠的高光谱图像分类。

Description

一种改进的高光谱图像去噪方法

技术领域

本发明涉及高光谱图像处理领域，特别涉及一种改进的高光谱图像去噪方法。

背景技术

近年来，高光谱图像以其拥有丰富的光谱信息和能够对地物进行识别的空间信息受到了越来越多的关注。并且其已成功地在考古发现、地质勘探、深林火灾检测、大气监控、军事作战等领域得到广泛的应用和产生深远的影响。然而高光谱图像在采集和传输的过程中，往往会受到多种不同类型噪声的污染，很大程度降低了数据的可靠性，而且还对后续的解混、分割和目标检测等造成严重的影响。因此，研究高光谱图像去噪技术在其应用当中是非常有必要的。

目前出现的基于光谱信号和二维图像去噪技术的高光谱图像去噪算法，都取得了不错的效果。但由于高光谱图像具有丰富的光谱信息和空间信息等特点，单一地利用光谱信息或者空间信息进行去噪，就其去噪效果而言是远远不够的。

发明内容

发明目的：本发明为了解决高光谱图像采集、传输过程等可能出现噪声的问题，提出了一种改进高光谱图像去噪方法,较为显著地提高了去噪的效果,可以针对实时响应和高精度不同需求的应用场景给出可靠的高光谱图像分类。

技术方案：本发明所述的改进的高光谱图像去噪方法，包括以下步骤：

(1)将待去噪的高光谱图像转换成空间光谱联合的二维矩阵；

(2)根据空间光谱联合的二维矩阵，采用基于欧式距离的近邻相似度计算策略计算高光谱图像的像素点与近邻的空间相似度；

(3)结合像素点空间局域相似性和光谱间低秩性建立去噪模型，恢复出原始无噪数据；

(4)采用原始无噪数据恢复出三维无噪高光谱图像。

进一步的，步骤(1)具体包括：

(1-1)获取待去噪的高光谱图像X∈R^m×n×b，其中，m、n分别是其空间结构的行数和列数，b是波段数；

(1-2)将高光谱图像X的第i个像素点在所有波段上的值记为一个向量d_i，i＝1,2,...,mn，mn是像素点数；

(1-3)将所有向量进行集合，即为空间光谱联合的二维矩阵D＝[d₁，d₂，...，d_mn]。

进一步的，步骤(2)具体包括：

(2-1)定义高光谱图像的模型为：D＝A+E+N，其中，A代表原始无噪数据，具有低秩特性，且A＝[α₁，α₂，...，α_mn]，向量α_i表示去噪后的第i个像素点在所有波段上的值，E代表稀疏噪声，N代表高斯噪声；

(2-2)计算像素点的空间相似度S：其中，σ表示标准差，tr(·)表示矩阵的迹，L是拉普拉斯矩阵。

进一步的，步骤(3)具体包括：

(3-1)结合空间像素点的相似性和光谱间低秩性，建立去噪模型为：

s.t. D＝A+E+N

式中，λ、γ、β分别为系数噪声项、高斯噪声项和空间邻域信息的折中因子；

(3-2)对去噪模型求解，得到原始无噪数据A。

进一步的，步骤(3-2)具体包括：

(3-2-1)将去噪模型转变成以下等价形式：

s.t. D＝A+E+N，A＝J

(3-2-2)其增广拉格朗日函数为：

式中，Y₁和Y₂为系数矩阵，μ为步长；

(3-2-3)固定其他项来更新J：

其中，针对核函数和F范数采用奇异值阈值的方法来求解，表示第k次迭代时的值；

(3-2-4)固定其他项来更新A：

(3-2-5)固定其他项来更新E：

(3-2-6)固定其他项来更新N：

(3-2-7)固定其他项来更新Y₁&Y₂：

Y_1，k+1＝Y_1，k+μ_k(D-A_k+1-E_k+1-N_k+1)

Y_2，k+1＝Y_2，k+μ_k(A_k+1-J_k+1)

(3-2-8)固定其他项来更新μ_k：

式中，ρ表示迭代步长，ε₀表示迭代阈值；

(3-2-9)判断终止条件：

||D-A-E-N||_∞＜ε₁

||J-A||_∞＜ε₂

式中，ε₁表示设定的阈值，ε₂表示设定的阈值；

(3-2-10)迭代终止后，得到无噪数据A。

有益效果：本发明与现有技术相比，其显著优点是：首先，通过结合空间领域相似信息，较好地保留了空间细节信息；再者，就是利用改进的低秩表示理论，针对不同噪声给出了统一的解决方法，在模型中，针对椒盐噪声、条纹噪声等稀疏噪声，专门给出了稀疏项E,而对高斯噪声也给出了高斯项N，能较为有效地去除了稀疏与高斯的混合噪声。本发明将两种方法建立在一个统一的模型中，去噪效果显著，具备较高的使用价值。

附图说明

图1是本发明的一个实施例的流程示意图；

图2是步骤(3-2)的具体流程示意图。

具体实施方式

如图1所示，本实施例的改进的高光谱图像去噪方法，包括以下步骤：

步骤1：将待去噪的高光谱图像转换成空间光谱联合的二维矩阵。

该步骤具体包括：

(1-1)获取待去噪的高光谱图像X∈R^m×n×b，其中，m、n分别是其空间结构的行数和列数，b是波段数；

(1-2)将高光谱图像X的第i个像素点在所有波段上的值记为一个向量d_i，i＝1,2,...,mn，mn是像素点数；

(1-3)将所有向量进行集合，即为空间光谱联合的二维矩阵D＝[d₁，d₂，...，d_mn]。

步骤2：根据空间光谱联合的二维矩阵，采用基于欧式距离的近邻相似度计算策略计算高光谱图像的像素点与近邻的空间相似度。

该步骤具体包括：

(2-1)定义高光谱图像的模型为：D＝A+E+N，其中，A代表原始无噪数据，具有低秩特性，且A＝[a₁，a₂，...，a_mn]，向量α_i表示去噪后的第i个像素点在所有波段上的值，E代表稀疏噪声，N代表高斯噪声；

(2-2)计算像素点的空间相似度S：其中，σ表示标准差，tr(·)表示矩阵的迹，L是拉普拉斯矩阵。

具体的，为了尽量保留原始图像空间的局域结构信息，认为如果d_i与其在空间结构的邻域内的d_j是非常相似的(例如它们是来自同一片草地的像元)，则与他们相对应的α_i与α_j也是非常相似的。W代表权重矩阵，则d_i与它空间邻域内的k-近邻d_j的权重为：

为了能够保留空间局部结构信息，可以通过最小化下面的函数来达到：

式中，L＝P-W,P是对角矩阵。

步骤3：结合像素点空间局域相似性和光谱间低秩性建立去噪模型，恢复出原始无噪数据。该步骤的方法定义为改进RPCA方法(Spatial Neighboring Similarity andImprove RPCA,S_IRPCA)。

步骤(3)具体包括：

(3-1)结合空间像素点的相似性和光谱间低秩性，建立去噪模型为：

s.t.D＝A+E+N

式中，λ、γ、β分别为系数噪声项、高斯噪声项和空间邻域信息的折中因子；

(3-2)对去噪模型求解，得到原始无噪数据A。

如图2所示，去噪模型的求解步骤具体为：

(3-2-1)将去噪模型转变成以下等价形式：

s.t.D＝A+E+N，A＝J

(3-2-2)其增广拉格朗日函数为：

式中，Y₁和Y₂为系数矩阵，μ为步长；

(3-2-3)固定其他项来更新J：

其中，针对核函数和F范数采用奇异值阈值的方法来求解，表示第k次迭代时的值；

(3-2-4)固定其他项来更新A：

(3-2-5)固定其他项来更新E：

(3-2-6)固定其他项来更新N：

(3-2-7)固定其他项来更新Y₁&Y₂：

Y_1，k+1＝Y_1，k+μ_k(D-A_k+1-E_k+1-N_k+1)

Y_2，k+1＝Y_2，k+μ_k(A_k+1-J_k+1)

(3-2-8)固定其他项来更新μ_k：

式中，ρ表示迭代步长，ε₀表示迭代阈值；

(3-2-9)判断终止条件：

||D-A-E-N||_∞式ε₁

||J-A||_∞＜ε₂

式中，ε₁表示设定的阈值，ε₂表示设定的阈值；

(3-2-10)迭代终止后，得到无噪数据A。

步骤4：采用原始无噪数据A恢复出三维无噪高光谱图像。

以上所揭露的仅为本发明一种较佳实施例而已，不能以此来限定本发明之权利范围，因此依本发明权利要求所作的等同变化，仍属本发明所涵盖的范围。

Claims (5)

1.一种改进的高光谱图像去噪方法，其特征在于包括以下步骤：

(1)将待去噪的高光谱图像转换成空间光谱联合的二维矩阵；

(2)根据空间光谱联合的二维矩阵，采用基于欧式距离的近邻相似度计算策略计算高光谱图像的像素点与近邻的空间相似度；

(3)结合像素点空间局域相似性和光谱间低秩性建立去噪模型，恢复出原始无噪数据；

(4)采用原始无噪数据恢复出三维无噪高光谱图像。

2.根据权利要求1所述的改进的高光谱图像去噪方法，其特征在于：步骤(1)具体包括：

(1-1)获取待去噪的高光谱图像X∈R^m×n×b，其中，m、n分别是其空间结构的行数和列数，b是波段数；

(1-2)将高光谱图像X的第i个像素点在所有波段上的值记为一个向量d_i，i＝1,2,...,mn，mn是像素点数；

(1-3)将所有向量进行集合，即为空间光谱联合的二维矩阵D＝[d₁，d₂，...，d_mn]。

3.根据权利要求2所述的改进的高光谱图像去噪方法，其特征在于：步骤(2)具体包括：

(2-1)定义高光谱图像的模型为：D＝A+E+N，其中，A代表原始无噪数据，具有低秩特性，且A＝[α₁，α₂，...，α_mn]，向量α_i表示去噪后的第i个像素点在所有波段上的值，E代表稀疏噪声，N代表高斯噪声；

(2-2)计算像素点的空间相似度S：其中，σ表示标准差，tr(·)表示矩阵的迹，L是拉普拉斯矩阵。

4.根据权利要求3所述的改进的高光谱图像去噪方法，其特征在于：步骤(3)具体包括：

(3-1)结合空间像素点的相似性和光谱间低秩性，建立去噪模型为：

$\underset{A,E,N}{\min }\left|\right|A|{|}_{*}+\mathrm{\λ}|\left|E\right|{|}_{2,1}+\mathrm{\γ}\left|\right|N|{|}_F^2+\mathrm{\β}tr\left({\mathrm{ALA}}^T\right)$

s.t.D＝A+E+N

式中，λ、γ、β分别为系数噪声项、高斯噪声项和空间邻域信息的折中因子；

(3-2)对去噪模型求解，得到原始无噪数据A。

5.根据权利要求4所述的改进的高光谱图像去噪方法，其特征在于：步骤(3-2)具体包括：

(3-2-1)将去噪模型转变成以下等价形式：

$\underset{A,E,N}{\min }\left|\right|A|{|}_{*}+\mathrm{\λ}|\left|E\right|{|}_{2,1}+\mathrm{\γ}\left|\right|N|{|}_F^2+\mathrm{\β}tr\left({\mathrm{JLJ}}^T\right)$

s.t.D＝A+E+N，A＝J

(3-2-2)其增广拉格朗日函数为：

$\begin{array}{c}L=\left|\right|A|{|}_{*}+\mathrm{\λ}|\left|E\right|{|}_{2,1}+\mathrm{\γ}\left|\right|N|{|}_F^2+\mathrm{\β}tr\left({\mathrm{JLJ}}^T\right)+tr\left(Y_1^T\right(D-A-E-N\left)\right)\\ +tr\left(Y_2^T\right(A-J\left)\right)+\frac{\mathrm{\μ}}{2}\left(\right||D-A-E-N|{|}_F^2+\left|\right|A-J\left|{|}_F^2\right)\end{array}$

式中，Y₁和Y₂为系数矩阵，μ为步长；

(3-2-3)固定其他项来更新J：

$\begin{array}{c}{J}_{k+1}=\arg \underset{J}{\min }\mathrm{\β}tr\left({\mathrm{JLJ}}^T\right)+tr\left(Y_2^T\right(A-J\left)\right)+\frac{\mathrm{\μ}}{2}\left|\right|A-J|{|}_F^2\\ =\arg \underset{J}{\min }\frac{1}{\mathrm{\μ}}\left|\right|J|{|}_{*}+\frac{\mathrm{\η}}{2}||J-(J_k-\frac{1}{\mathrm{\η}}\[J_k(\frac{2\mathrm{\β}}{\mathrm{\μ}}L+I)-(A+\frac{Y_2}{\mathrm{\μ}})\])|{|}_F^2\end{array}$

其中，针对核函数和F范数采用奇异值阈值的方法来求解，⊙_k表示第k次迭代时⊙的值；

(3-2-4)固定其他项来更新A：

$\begin{array}{c}{A}_{k+1}=\arg \underset{A}{\min }tr\left(Y_1^T\right(D-A-E-N\left)\right)+tr\left(Y_2^T\right(A-J\left)\right)+\frac{\mathrm{\μ}}{2}\left(\right||D-A-E-N|{|}_F^2+\left|\right|A-J\left|{|}_F^2\right)\\ =\frac{1}{2}(D-E-N+J+\frac{Y_1}{\mathrm{\μ}}-\frac{Y_2}{\mathrm{\μ}})\end{array}$

(3-2-5)固定其他项来更新E：

$E_{k+1}=\arg \underset{E}{\min }\mathrm{\λ}\left|\right|E|{|}_{1,1}+tr\left(Y_1^T\right(D-A-E-N\left)\right)+\frac{\mathrm{\μ}}{2}||D-A-E-N|{|}_F^2$

$=\arg \underset{E}{\min }\mathrm{\λ}\left|\right|E|{|}_{2,1}+\frac{\mathrm{\μ}}{2}||E-(D-A-N+\frac{Y_1}{\mathrm{\μ}})|{|}_F^2$

(3-2-6)固定其他项来更新N：

$\begin{array}{c}{N}_{k+1}=\arg \underset{N}{\min }\mathrm{\γ}\left|\right|N|{|}_F^2+tr\left(Y_1^T\right(D-A-E-N\left)\right)+\frac{\mathrm{\μ}}{2}||D-A-E-N|{|}_F^2\\ =\frac{\mathrm{\μ}}{2\mathrm{\γ}+\mathrm{\μ}}(D-A-E+\frac{Y_1}{\mathrm{\μ}})\end{array}$

(3-2-7)固定其他项来更新Y₁&Y₂：

Y_1，k+1＝y_1，k+μ_k(D-A_k+1-E_k+1-N_k+1)

Y_2，k+1＝Y_2，k+μ_k(A_k+1-J_k+1)

(3-2-8)固定其他项来更新μ_k：

${\mathrm{\&mu;}}_{k+1}=\left\{\begin{array}{cc}{\mathrm{\&rho;\&mu;}}_k,& if\frac{{\mathrm{\&mu;}}_k\left|\right|E_{k+1}-E_k|{|}_F}{\left|\right|D|{|}_F}<{\mathrm{\&epsiv;}}_{0}or\frac{{\mathrm{\&mu;}}_k\left|\right|N_{k+1}-N_k|{|}_F}{\left|\right|D|{|}_F}<{\mathrm{\&epsiv;}}_0\\ {\mathrm{\&mu;}}_k,& otherwise\end{array}\right.$

式中，ρ表示迭代步长，ε₀表示迭代阈值；

(3-2-9)判断终止条件：

||D-A-E-N||_∞＜ε₁

||J-A||∞＜ε₂

式中，ε₁表示设定的阈值，ε₂表示设定的阈值；

(3-2-10)迭代终止后，得到无噪数据A。