CN109557031A

CN109557031A - 一种岩心高光谱蚀变信息提取方法

Info

Publication number: CN109557031A
Application number: CN201811312872.6A
Authority: CN
Inventors: 田青林; 潘蔚; 陈雪娇; 余长发
Original assignee: Beijing Research Institute of Uranium Geology
Current assignee: Beijing Research Institute of Uranium Geology
Priority date: 2018-11-06
Filing date: 2018-11-06
Publication date: 2019-04-02

Abstract

本发明属于遥感地质调查技术领域，具体涉及一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法。本发明包括如下步骤：步骤1、岩心高光谱影像数据的读入；步骤2、数据预处理；步骤3、蚀变矿物端元光谱获取；步骤4、应用小波包变换分解每个像元的高光谱曲线；步骤5、计算蚀变矿物端元光谱和每个像元光谱的信息熵特征矢量；步骤6、利用权重光谱角制图算法计算信息熵特征矢量夹角；步骤7、输出岩心高光谱影像蚀变信息提取结果。本发明能够有效去除岩心高光谱数据的波段相关性，降低数据冗余，减小原始光谱波段数多、波段相关性强对信息提取精度的影响，适用范围广泛。

Description

一种岩心高光谱蚀变信息提取方法

技术领域

本发明属于遥感地质调查技术领域，具体涉及一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法。

背景技术

高光谱遥感具有光谱分辨率高，“图谱”合一的特点，是21世纪遥感领域重要的研究方向之一。近年来，高光谱技术在我国地质勘查领域应用越发广泛，利用高光谱数据进行矿物填图，提取蚀变矿物信息，为快速圈定成矿远景区和优选靶区发挥重要作用。与多光谱遥感数据相比，高光谱数据具有数据量大、波段数多且波段间存在强相关性特点，这为直接对高光谱数据进行处理带来一定的困难。

小波包分解是一种比小波分解更为精细的分解方法，能同时对低频和高频信息进行分解，而且这种分解既无冗余，也无疏漏。从物理意义上看，小波包变换通过伸缩和平移等运算对信号进行多尺度分析，继承发展了窗口傅里叶变换的局部化思想，同时克服了窗口大小不随频率变化等不足，提供了一个随频率改变的时频窗口，解决了傅里叶变换不能解决的许多问题。光谱角制图算法(spectral angle mapper,SAM)是一种基于光谱曲线整体相似性的算法，要求地物的整体光谱均参与计算，致使非特征光谱及噪声的比重会压制诊断性特征光谱对分类的影响。而权重光谱角制图算法(weight spectral angle mapper,WSAM)是在图像像元光谱与参考光谱差异较大的特征区间进行权重设置，再计算权重设置后的光谱角，能够突出局部特征差异。

发明内容

本发明解决的技术问题：本发明提供一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，能够有效去除岩心高光谱数据的波段相关性，降低数据冗余，减小原始光谱波段数多、波段相关性强对信息提取精度的影响，适用范围广泛。

本发明采用的技术方案：

一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，包括如下步骤：

步骤1、岩心高光谱影像数据的读入；步骤2、数据预处理；步骤3、蚀变矿物端元光谱获取；步骤4、应用小波包变换分解每个像元的高光谱曲线；步骤5、计算蚀变矿物端元光谱和每个像元光谱的信息熵特征矢量；步骤6、利用权重光谱角制图算法计算信息熵特征矢量夹角；步骤7、输出岩心高光谱影像蚀变信息提取结果。

所述步骤1中，对读入的任一高光谱像元X，X＝(x₁,x₂,…,x_k)，其中k是岩心高光谱数据的波段总数，x_i(i＝1，2，…，k)代表第i个波段的光谱响应值。

所述步骤2包括如下步骤：

步骤2.1、对Hyspex系统的定标参数对测量得到的岩心高光谱原始数据进行辐射定标，实现图像DN值到辐射值的转换；

步骤2.2、采用经验线性法进行反射率计算，实现辐射值到反射率的转换；

步骤2.3、运用小波变换方法对反射率数据进行光谱去噪，减弱或消除锯齿噪声的影响。

所述步骤3中，对步骤2中预处理后的数据进行MNF变换，通过降维去噪处理实现信号噪声分离；利用纯净像元指数方法提取图像中的“纯”像元；通过N维散度分析得到各“纯”像元的光谱曲线，对照标准矿物光谱库及理论分析确定最终的矿物端元光谱。

所述步骤4中，选定小波基函数ψ及分解层数j，应用ψ小波基对每个高光谱像元X进行j层小波包变换，得到p个子分量，p＝0，1，…，2^j-1，它们对应的能量分解系数分别为d_i ^p，i为信号在对应子分量中的分解系数编号。

所述步骤5中，首先，根据步骤4方法得到每个矿物端元光谱曲线和像元光谱曲线的小波包分解系数，计算每个子分量的能量E_p，可由小波包分解系数的平方和计算得到，即

然后，计算信号总能量E，为所有子分量的能量之和，即

其次，计算相对小波能量l_p，为每个子分量能量与信号总能量之比，即

l_p＝E_p/E (3)

每个子分量上小波包信息熵定义为

W_p＝-l_plogl_p (4)

通过计算得到每个像元光谱曲线X所对应的p个子分量小波包信息熵，构成p维信息熵特征矢量T，即

T＝(W₁,W₂,…,W_p) (5)

所述步骤6中，光谱曲线包含的波段总数记为N，假设差异较大的特征区间包含有N₁个波段，对应的像元光谱和参考光谱反射率值为x_aj和y_aj，j∈[1,N₁]，差异较小的普通区间包含有N₂个波段，对应的像元光谱和参考光谱的反射率值为x_bi和y_bi，i∈[1,N₂]，N₁+N₂＝N，将差异较大的特征区间添加权重系数γ，γ≥1，γ＝1时即为传统光谱角制图算法，权重光谱角制图算法的计算公式如下：

找到步骤5计算得到的蚀变矿物端元光谱和每个像元光谱的信息熵特征矢量差异较大的区间，对其设置权重，采用公式(6)计算权重设置后的信息熵特征矢量夹角。

所述步骤7中，通过对步骤6计算得到的矿物端元光谱曲线与像元光谱曲线所对应的信息熵特征矢量夹角设定阈值，提取矿化蚀变，最终输出岩心高光谱影像蚀变信息提取结果。

本发明的有益效果：

本发明提供的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，能够有效去除岩心高光谱数据的波段相关性，降低数据冗余，减小原始光谱波段数多、波段相关性强对信息提取精度的影响，同时增大了类间地物的可区分性，更利于分类和信息提取，适用范围广泛。

附图说明

图1为本发明所提供的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法的流程图。

具体实施方式

下面结合附图和具体实施例对本发明作进一步详细说明。

本发明所提供的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，包括如下步骤：

步骤1，读取试验区Hyspex岩心高光谱数据，总计256个波段。对任一高光谱像元X，记为X＝(x₁,x₂,…,x₂₅₆)，

其中，x_i(i＝1，2，…，256)代表第i个波段的光谱响应值。

步骤2，数据预处理：

首先，利用ENVI软件对Hyspex系统的定标参数对测量得到的岩心高光谱原始数据进行辐射定标，实现图像DN值到辐射值的转换；然后，在ENVI软件中采用经验线性法(empirical line calibration，ELC)进行反射率计算，实现辐射值到反射率的转换；最后，在matlab软件中运用小波变换方法对反射率数据进行光谱去噪，减弱或消除锯齿噪声的影响。

步骤3，蚀变矿物端元光谱获取：

在ENVI软件中对预处理后的数据进行MNF变换，通过降维去噪处理实现信号噪声分离；利用纯净像元指数(pixel purity index，PPI)方法提取图像中的“纯”像元；通过N维散度分析得到各“纯”像元的光谱曲线，对照标准矿物光谱库及理论分析确定最终的矿物端元光谱。具体包括高岭石、地开石、绿泥石和伊利石4种矿物端元光谱曲线。

步骤4，应用小波包变换分解每个像元的高光谱曲线：

选定小波基函数ψ为db3小波，分解层数j为5层，则将X＝(x₁,x₂,…,x₂₅₆)进行5层小波包分解，得到32个子分量，它们对应的能量分解系数分别为i为信号在对应子分量中的分解系数编号。

步骤5，计算蚀变矿物端元光谱和每个像元光谱的信息熵特征矢量：

首先，根据步骤4方法得到每个矿物端元光谱曲线和像元光谱曲线的小波包分解系数，计算每个子分量的能量E_p，可由小波包分解系数的平方和计算得到，即

然后，计算信号总能量E，为所有子分量的能量之和，即

l_p＝E_p/E (3)

每个子分量上小波包信息熵定义为

W_p＝-l_plogl_p (4)

通过计算得到每个像元光谱曲线X所对应的32个子分量小波包信息熵，构成p维信息熵特征矢量T，即

T＝(W₁,W₂,…,W₃₂) (5)

步骤6，利用权重光谱角制图算法计算信息熵特征矢量夹角：

光谱曲线包含的波段总数记为N，假设差异较大的特征区间包含有N₁个波段，对应的像元光谱和参考光谱反射率值为x_aj和y_aj，j∈[1,N₁]，差异较小的普通区间包含有N₂个波段，对应的像元光谱和参考光谱的反射率值为x_bi和y_bi，i∈[1,N₂]，N₁+N₂＝N，将差异较大的特征区间添加权重系数γ，γ≥1，γ＝1时即为传统光谱角制图算法，权重光谱角制图算法的计算公式如下：

步骤7，输出岩心高光谱影像蚀变提取结果：

通过对步骤6计算得到的矿物端元光谱曲线与像元光谱曲线所对应的信息熵特征矢量夹角设定阈值，提取矿化蚀变，最终输出高岭石、地开石、绿泥石和伊利石4种蚀变信息提取结果。

Claims

1.一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：包括如下步骤：

步骤(1)、岩心高光谱影像数据的读入；步骤(2)、数据预处理；步骤(3)、蚀变矿物端元光谱获取；步骤(4)、应用小波包变换分解每个像元的高光谱曲线；步骤(5)、计算蚀变矿物端元光谱和每个像元光谱的信息熵特征矢量；步骤(6)、利用权重光谱角制图算法计算信息熵特征矢量夹角；步骤(7)、输出岩心高光谱影像蚀变信息提取结果。

2.根据权利要求1所述的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：

所述步骤(1)中，对读入的任一高光谱像元X，X＝(x₁,x₂,…,x_k)，其中k是岩心高光谱数据的波段总数，x_i(i＝1，2，…，k)代表第i个波段的光谱响应值。

3.根据权利要求2所述的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：

所述步骤(2)包括如下步骤：

步骤(2.1)、对Hyspex系统的定标参数对测量得到的岩心高光谱原始数据进行辐射定标，实现图像DN值到辐射值的转换；

步骤(2.2)、采用经验线性法进行反射率计算，实现辐射值到反射率的转换；

步骤(2.3)、运用小波变换方法对反射率数据进行光谱去噪，减弱或消除锯齿噪声的影响。

4.根据权利要求3所述的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：

所述步骤(3)中，对步骤(2)中预处理后的数据进行MNF变换，通过降维去噪处理实现信号噪声分离；利用纯净像元指数方法提取图像中的“纯”像元；通过N维散度分析得到各“纯”像元的光谱曲线，对照标准矿物光谱库及理论分析确定最终的矿物端元光谱。

5.根据权利要求4所述的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：

所述步骤(4)中，选定小波基函数ψ及分解层数j，应用ψ小波基对每个高光谱像元X进行j层小波包变换，得到p个子分量，p＝0，1，…，2^j-1，它们对应的能量分解系数分别为i为信号在对应子分量中的分解系数编号。

6.根据权利要求5所述的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：

所述步骤(5)中，首先，根据步骤(4)方法得到每个矿物端元光谱曲线和像元光谱曲线的小波包分解系数，计算每个子分量的能量E_p，可由小波包分解系数的平方和计算得到，即

然后，计算信号总能量E，为所有子分量的能量之和，即

l_p＝E_p/E (3)

每个子分量上小波包信息熵定义为

W_p＝-l_plogl_p (4)

T＝(W₁,W₂,…,W_p) (5)。

7.根据权利要求6所述的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：

所述步骤(6)中，光谱曲线包含的波段总数记为N，假设差异较大的特征区间包含有N₁个波段，对应的像元光谱和参考光谱反射率值为x_aj和y_aj，j∈[1,N₁]，差异较小的普通区间包含有N₂个波段，对应的像元光谱和参考光谱的反射率值为x_bi和y_bi，i∈[1,N₂]，N₁+N₂＝N，将差异较大的特征区间添加权重系数γ，γ≥1，γ＝1时即为传统光谱角制图算法，权重光谱角制图算法的计算公式如下：

找到步骤(5)计算得到的蚀变矿物端元光谱和每个像元光谱的信息熵特征矢量差异较大的区间，对其设置权重，采用公式(6)计算权重设置后的信息熵特征矢量夹角。

8.根据权利要求7所述的一种基于小波包变换和权重光谱角制图的岩心高光谱蚀变信息提取方法，其特征在于：

所述步骤(7)中，通过对步骤(6)计算得到的矿物端元光谱曲线与像元光谱曲线所对应的信息熵特征矢量夹角设定阈值，提取矿化蚀变，最终输出岩心高光谱影像蚀变信息提取结果。