一种高光谱地物分类识别方法与流程

本发明涉及地物分类识别技术领域，具体地讲涉及一种高光谱地物分类识别方法。

背景技术：

高光谱遥感技术利用先进的光谱成像技术，带来上百个波段高达分米级别分辨率的遥感数据，提供了异常丰富的光谱波段信息，为诸多民用领域，如森林土地资源调查,自然灾害观测、气象预报等提供了全新的解决方法。

然而，高光谱遥感数据相对普通数据而言，数据量大大增加，使得处理过程中耗费大量时间和服务器资源；同时，多达上百个波段的光谱数据造成高光谱遥感数据普遍存在“同物异谱，同谱异物”现象，使得地物分类很容易出现错分和漏分现象。

技术实现要素：

根据现有技术中存在的问题，本发明提供了一种高光谱地物分类识别方法，能够快速的对高光谱遥感图像的数据进行预处理，并高精确度的对高光谱遥感图像进行地物分类。

为实现上述发明目的，本发明提供了一种高光谱地物分类识别方法，包括如下步骤：

S1，数据校正，对已知分类情况的遥感图像依次进行大气校正和几何校正，得到校正后的遥感图像；

S2，去噪处理，对校正后的遥感图像进行去噪处理，得到去噪后的遥感图像；

S3，重组分批，对去噪后的遥感图像进行重组分批，得到重组后的遥感图像；

S4，特征提取，对重组后的遥感图像使用构建的多层深度网络进行训练，提取出重组后的遥感图像的深层特征；

S5，分类识别，将待分类识别图像进行步骤S1、步骤S2和步骤S3，再通过步骤S4中所述的网络进行训练，得到待分类识别图像的深层特征，将待分类识别图像的深层特征与步骤S4中得到的深层特征进行比较，得到待分类识别图像的最终分类结果。

优选的，步骤S1中，所述遥感图像通过高光谱传感器获取，依次对遥感图像进行大气校正和几何校正；所述大气校正是指在ENVIClassic软件中对遥感图像的ConvertFileParameter参数进行设置使其辐射亮度值符合FLAASH要求，再通过MultispectralSettings进行多光谱设置，得到大气校正后的遥感图像；所述几何校正是指在ENVIClassic软件中通过调整UTM中的Zone的数值和Datum的参数，再在RegistrationParameters的对话框中对WrapMethod、Resampling和Background进行设置，完成对大气校正后的遥感图像的几何校正，得到校正后的遥感图像。

进一步优选的，步骤S2中，对校正后的遥感数图像进行去噪处理的公式如下：

其中c(x,y,t)为扩散系数；

当时：

当时：

当时：

其中，I(x,y,t)是偏微分方程的解，也就是去噪后的遥感图像，x、y表示去噪后的遥感图像的像素点的坐标，t表示噪声的方差；I表示校正后的遥感图像，Gσ是高斯滤波，σ是高斯核的标准偏差，div和分别表示散度算子和梯度算子，||表示模值，*表示卷积；k是校正后的遥感图像数据的梯度阈值,k＞0，Δk是k的偏置,0＜Δk≤5；Ixx和Iyy是校正后的遥感图像I的二阶偏导数，扩散系数c(x,y,t)是关于梯度的函数，α和β是连贯系数，α+β＝1。

更进一步优选的，步骤S3的具体步骤如下：

S31，利用matalab软件内置的imread函数读取去噪后的遥感图像的A个波段数据，将它们都转化为A个大小为M×M的图片，每个图片相当于一个像素点，再把A个大小为M×M的图片拉伸为A列M²×1大小的矩阵，此矩阵即是由遥感图像转化为的A×M²个像素点；

S32，再将每列矩阵平均分为C份，每份含有B个元素,B＝M²/C，即将A列M²×1的矩阵数据分成C份大小为B×A的数据，最终得到A个含有C份大小为B×A的图像，即重组后的遥感图像。

更进一步优选的，步骤S4中，利用限制玻尔兹曼机构建多层网络对重组后的遥感图像进行训练，提取出重组后的遥感图像的深层特征，具体过程如下：

限制玻尔兹曼机的网络结构由可见层和隐藏层组成；n表示可见层中神经元的个数，m表示隐藏层中神经元的个数，v表示可见层的状态向量，vj∈v，vj＝(v1,v2,...,vn)，vj表示第j个可见层中神经元的状态向量；h表示隐藏层的状态向量，hi∈h，hi＝(h1,h2,...,hm)，hi表示第i个隐藏层中神经元的状态向量；a表示可见层的偏置向量，aj∈a，aj＝(a1,a2,...,an)，aj表示第j个可见层中神经元的偏置向量；b表示隐藏层的偏置向量，bi∈b，bi＝(b1,b2,...,bm)，bi表示第i个隐藏层中神经元的偏置向量；W＝(wi,j)∈R^n×m表示可见层和隐藏层之间的权值矩阵，wij表示可见层中第i个神经元和隐藏层中第j个神经元之间的权值，R^n×m表示所有n行m列的矩阵；

将重组后的遥感图像输入到限制玻尔兹曼机的网络结构的可见层，可见层有a个节点，将A列中每行的A个像素点分别输入到可见层的一个节点中，a＝A，求得隐藏层中第i个神经元的激活概率为：

P(hi＝1|v)＝σ(bi+∑vjWij)

可见层中第j个神经元的激活概率为：

其中，P()表示激活概率函数，为sigmoid激活函数；

将可见层的状态设置为一个训练样本v⁽⁰⁾:＝v，随机初始化W，执行k步Gibbs采样，计算所有隐藏单元的激活概率，确定可见层中第j个神经元vj取值为1的概率，随着a列中每行的像素点的不断输入，进而产生可见层的一个重构，按下式迭代更新：

Δbi＝Δbi+[P(hi＝1|v⁽⁰⁾)-P(hi＝1|v^(k))]

其中，表示可见层中第j个神经元的初始值，表示可见层中第j个神经元迭代k次的值，Δwij示Wij的偏置值，Δaj表示aj的偏置值，Δbj表示bj的偏置值；等号左边的Δwij表示等号右边的Δwij的下一个迭代值，等号左边的Δaj表示等号右边的Δaj的下一个迭代值，等号左边的Δbj表示等号右边的Δbj的下一个迭代值；

得到限制玻尔兹曼机训练算法的参数刷新公式：

其中，η是训练网络的学习率，η＞0，取值越大，参数更新的速度越快；ns表示训练样本的个数，即M²×A个；W,a和b代表一个类别的深层特征，等号左边的W表示等号右边的W的下一个刷新值，等号左边的a表示等号右边的a的下一个刷新值，等号左边的b表示等号右边的b的下一个刷新值；

将W,a和b与标签进行对比并计算出相应的差值，若差值小于或等于设定的loss率大小，则停止迭代和刷新W、a、b，则此时的W,a和b代表着该标签的类别；

如果差值大于设定的loss率大小，就继续迭代及刷新W、a、b，直到W,a和b与标签的差值小于或等于设定的loss率大小，最终的W,a和b代表着该标签对应类别的特征。

更进一步优选的，步骤S5的具体过程如下：先将待分类识别图像依次进行步骤S1所述的大气校正、几何校正及步骤S2所述的去噪处理，再将处理后的待分类识别图像输入到步骤S5所述的限制玻尔兹曼机的网络结构中，得到处理后的待分类识别图像的W,a和b的值，将处理后的待分类识别图像的W,a和b的值与步骤S5中得到的W,a和b值进行比较，两者最接近的那一组步骤S5中得到的W,a和b值对应的类别就是处理后的待分类识别图像的最终分类结果。

更进一步优选的，所述标签是指对实验场地进行实地测量，将每个不同的物种区域都进行了不同的标注，每种标注对应着一组W,a、b。

本发明的有益效果在于：

1)本发明的高光谱地物分类识别方法中的特征提取利用了限制玻尔兹曼机构建多层训练网络的算法，通过深度学习对遥感图像的数据进行深入挖掘，进而提取到遥感图像数据的深层特征，最终通过将深层特征的参数与标记的样本进行对比，从而确定遥感图像的各像素点的分类，提高了遥感图像的高光谱地物分类的准确度。

2)本发明的高光谱地物分类识别方法中的去噪过程利用了自适应扩散系数对高光谱遥感图像的不同噪声区域进行不同程度的扩散平滑，不仅有效的去除噪声干扰，而且又能很好地保留了边缘信息，便于后续的分类。

3)本发明的高光谱地物分类识别方法中，由于高光谱遥感图像一般都包含上百个波段，每个波段相当于一张图片，如果将这上百幅图片直接训练，将会耗费大量时间和服务器资源。因此，本发明方法的重组分批过程中将多波段的高光谱遥感图像数据进行单波段剥离，将上百个波段的数据分别拉伸为一个列向量，存放在一个大矩阵中，再将其分割成若干行，每次按照网络的节点数分批输入数据进行后续的训练，提高了光谱遥感图像的数据预处理速度。

附图说明

图1为本发明的地物分类识别方法的流程框图。

图2为本发明的地物分类识别方法中对去噪后的遥感图像进行重组分批示意图。

图3本发明的地物分类识别方法中的限制玻尔兹曼机构建多层训练网络的示意图。

图4本发明的地物分类识别方法的流程示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

如图1、图4所示，一种高光谱地物分类识别方法，包括如下步骤：

S1，数据校正，对已知分类情况的遥感图像依次进行大气校正和几何校正，得到校正后的遥感图像；

具体的，大气校正的步骤：启动ENVI Classic软件，将RGB Color设置为353，接着选择Basic Tools→Preprocessing→Calibration Utilities→Landset Calibration，选择已知分类情况的遥感图像，Landset定标工具会从元数据文件中自动获取相关的参数信息，包括成像日期，定标参数等，选择Calibration Type：Radiance，将文件重命名为“11”；在主菜单中，选择Basic Tools→Band Math，在Enter an expression输入b1/10.0，点击Add to list，OK，接着在Variables to Bands Pairings对话框中，点击“map variable to input file”，选择“11”，然后点击Choose，文件命名为“22”；在主菜单中选择Basic Tools—Convert Data(BSQ、BIL、BIP)，打开“22”，OK，按钮打开Convert File Parameter，选择Output Interleave：BIL，Convert In Place：Yes，确定，就获得了符合FLAASH要求的辐射亮度值；在主菜单中，选择Spectral→FLAASH，打开FLAASH功能。单击Input Radiance Image按钮，打开“22”，接着选择Use single factor for all bands(single scale factor:10.00000)，选择输出文件名和路径，命名为“33，然后点击Multispectral Settings进行多光谱设置，高级设置默认即可；最后单击Apply按钮即可，得到大气校正后的遥感图像“33”。

几何校正的步骤如下：在ENVI Classic中打开得到大气校正后的遥感图像“33”，多光谱影像选择5,4,3波段打开；选择file→open vector file，打开矢量数据.shp文件，选择Map→registration→select gcps:image to map,选择UTM，Zone数值改为16，Datum改为NAD83，点击OK得到Ground Control Points Selection；在主影像窗口中，选择Tools→Pixel Locator，将影像光标移到地物的某一个位置，点击Apply；在矢量窗口中，用鼠标左键点击并拖曳矢量光标，当十字形光标位于所需要的道路交叉时，松开鼠标左键，把矢量光标放置在道路的交叉口上，点击鼠标右键，在弹出的窗口中选择Export Map Location，新的地图坐标就会出现在Ground Control Points Selection对话框中，点击AddPoint，添加控制点，并观察RMS误差的变化；选择完控制点后，在Control Points Selection对话框中选择Options→Wrap Files，在Input Wrap Image对话框中选择要校正的影像，点击OK，得到Registration Parameters 对话框，将Wrap Method改为RST，Resampling改为Cubic Convolution，Background设为0；输出校正后的影像名，点击OK,开始进行影像到地形图的校正，结束后最终得到校正后的遥感图像。

S2，去噪处理，对校正后的遥感图像进行去噪处理，得到去噪后的遥感图像；

具体的，对校正后的遥感数图像进行去噪处理的公式如下：

其中c(x,y,t)为扩散系数；

当时：

当时：

当时：

其中，I(x,y,t)是偏微分方程的解，也就是去噪后的遥感图像，x、y表示去噪后的遥感图像的像素点的坐标，t表示噪声的方差；I表示校正后的遥感图像，Gσ是高斯滤波，σ是高斯核的标准偏差，div和分别表示散度算子和梯度算子，||表示模值，*表示卷积；k是校正后的遥感图像数据的梯度阈值,k＞0，Δk是k的偏置,0＜Δk≤5；Ixx和Iyy是校正后的遥感图像I的二阶偏导数，扩散系数c(x,y,t)是关于梯度的函数，α和β是连贯系数，α+β＝1。

S3，重组分批，对去噪后的遥感图像进行重组分批，得到重组后的遥感图像；

具体的，如图2所示，首先利用matalab软件内置的imread函数读取去噪后的遥感图像的A个波段数据，将它们都转化为A个大小为M×M的图片，每个图片相当于一个像素点，再把A个大小为M×M的图片拉伸为A列M²×1大小的矩阵，此矩阵即是由遥感图像转化为的A×M²个像素点；

再将每列矩阵平均分为C份，每份含有B个元素,B＝M²/C，即将A列M²×1的矩阵数据分成C份大小为B×A的数据，最终得到A个含有C份大小为B×A的图像，即重组后的遥感图像。

S4，特征提取，对重组后的遥感图像使用构建的多层深度网络进行训练，提取出重组后的遥感图像的深层特征；

具体的，如图3所示，图3为限制玻尔兹曼机构建多层网络的示意图，利用限制玻尔兹曼机构建多层网络对重组后的遥感图像进行训练，提取出重组后的遥感图像的深层特征，具体过程如下：

限制玻尔兹曼机的网络结构由可见层和隐藏层组成；n表示可见层中神经元的个数，m表示隐藏层中神经元的个数，v表示可见层的状态向量，vj∈v，vj＝(v1,v2,...,vn)，vj表示第j个可见层中神经元的状态向量；h表示隐藏层的状态向量，hi∈h，hi＝(h1,h2,...,hm)，hi表示第i个隐藏层中神经元的状态向量；a表示可见层的偏置向量，aj∈a，aj＝(a1,a2,...,an)，aj表示第j个可见层中神经元的偏置向量；b表示隐藏层的偏置向量，bi∈b，bi＝(b1,b2,...,bm)，bi表示第i个隐藏层中神经元的偏置向量；W＝(wi,j)∈R^n×m表示可见层和隐藏层之间的权值矩阵，wij表示可见层中第i个神经元和隐藏层中第j个神经元之间的权值，R^n×m表示所有n行m列的矩阵；

将重组后的遥感图像输入到限制玻尔兹曼机的网络结构的可见层，可见层有a个节点，将A列中每行的A个像素点分别输入到可见层的一个节点中，a＝A，求得隐藏层中第i个神经元的激活概率为：

P(hi＝1|v)＝σ(bi+∑vjWij)

第j个可见单元的激活概率为：

其中，P()表示激活概率函数，为sigmoid激活函数；

将可见层的状态设置为一个训练样本v⁽⁰⁾:＝v，随机初始化W，执行k步Gibbs采样，计算所有隐藏单元的激活概率，确定可见层中第j个神经元vj取值为1的概率，随着a列中每行的像素点的不断输入，进而产生可见层的一个重构，按下式迭代更新：

Δbi＝Δbi+[P(hi＝1|v⁽⁰⁾)-P(hi＝1|v^(k))]

其中，表示可见层中第j个神经元的初始值，表示可见层中第j个神经元迭代k次的值，Δwij示Wij的偏置值，Δaj表示aj的偏置值，Δbj表示bj的偏置值；等号左边的Δwij表示等号右边的Δwij的下一个迭代值，等号左边的Δaj表示等号右边的Δaj的下一个迭代值，等号左边的Δbj表示等号右边的Δbj的下一个迭代值；

得到限制玻尔兹曼机训练算法的参数刷新公式：

其中，η是训练网络的学习率，η＞0，取值越大，参数更新的速度越快；ns表示训练样本的个数，即M²×A个；W,a和b代表一个类别的深层特征，等号左边的W表示等号右边的W的下一个刷新值，等号左边的a表示等号右边的a的下一个刷新值，等号左边的b表示等号右边的b的下一个刷新值；

将W,a和b与标签进行对比并计算出相应的差值，若差值小于或等于设定的loss率大小，则停止迭代和刷新W、a、b，则此时的W,a和b代表着该标签的类别；

如果差值大于设定的loss率大小，就继续迭代及刷新W、a、b，直到W,a和b与标签的差值小于或等于设定的loss率大小，最终的W,a和b代表着该标签对应类别的特征。

所述标签是指对实验场地进行实地测量，将每个不同的物种区域(松树区域、草地区域、杨树区域、沙滩区域、水源区域等)都进行了不同的标注，每种标注对应着一组W,a、b。

S5，先将待分类识别图像依次进行步骤S1所述的大气校正、几何校正及步骤S2所述的去噪处理，再将处理后的待分类识别图像输入到步骤S5所述的限制玻尔兹曼机的网络结构中，得到处理后的待分类识别图像的W,a和b的值，将处理后的待分类识别图像的W,a和b的值与步骤S5中得到的W,a和b值进行比较，两者最接近的那一组步骤S5中得到的W,a和b值对应的类别就是处理后的待分类识别图像的最终分类结果。

综上所述，本发明的地物分类识别方法能够快速的对高光谱遥感图像的数据进行预处理，并高精确度的对高光谱遥感图像进行地物分类。