一种遥感影像变化识别方法及系统

技术领域

本发明涉及遥感数据分析技术领域，特别涉及一种遥感影像变化识别方法及系统。

背景技术

随着科学技术的发展，人类开发资源与改造自然的能力不断增强，自然界的变化和人类的各种活动，每天都在改变着地表景观及土地利用形式。由于遥感对地观测具有实时、快速、覆盖范围广、多光谱、周期性等特点，遥感技术已经成为变化检测最主要的技术手段。

遥感影像的变化识别就是从不同时间获取的同一地理区域的多时相遥感影像中，分析地表覆盖变化趋势及演化规律，在土地利用/覆盖、灾害检测、森林植被变化、城市扩展与布局等诸多领域发挥着积极重要的作用，具有迫切的科学应用需求和广泛的应用前景。

现有技术中使用PCA和K-Means对多时空卫星图像进行无监督的变化检测但是现有技术中识别变化的依据是k-means聚类中样本数量最少的一类，但在实际问题中，若两期影像全部发生变化，那么所识别的结果会不够准确。

发明内容

为了解决现有技术的问题，本发明提供了一种遥感影像变化识别方法及系统。

一方面，提供了一种遥感影像变化识别方法，所述方法包括：

裁剪前后两期对应边界区域的栅格数据影像，

对所述栅格数据影像进行主成分分析，得到灰度图像；

对所述灰度图像进行缩放，使得大小完全一致，并求差值，得到差值图像；

对所述差值图像进行主成分分析，增加区分度；

对分析后的所述差值图像进行聚类分析，选定阈值，挑选出高于所述阈值中心点所涵盖的所有样本点作为变化像素，其余为未变化像素，分别对所述变化像素和所述未变化像素赋值，得到二值图像；

将所述二值图像转换为矢量图，统计所述矢量图中赋值得到所述变化像素的面积和所述矢量图的总面积，得到变化率。

进一步地，分别对所述变化像素和所述未变化像素赋值时，所述变化像素赋值为255，所述未变化像素赋值为0。

进一步地，所述裁剪前后两期对应边界区域的栅格数据影像包括：

以要素图斑为单位，从矢量数据图层中导出每个图斑所对应的以要素ID命令的矢量边界数据，每一份矢量边界数据分别裁剪前后两期对应边界区域的栅格数据影像。

进一步地，对所述栅格数据影像进行主成分分析包括：

GDAL分别读取对应边界区域的前后两期栅格数据影像的波段数据以及坐标信息，将波段数据读入数组；

对GDAL读取到的波段数组进行变形和转置，处理成主成分分析所需要的正确数据格式，调用sklearn库提供的PCA接口对处理后的数组进行主成分分析，将原本的三维数据投影到一维特征上，去除图像冗余信息，将原本的多波段栅格数据影像转换为单波段的灰度图像。

进一步地，对所述差值图像进行主成分分析包括：

提取所述差值图像的特征向量，构建差值图像的特征向量数据集，作为主成分分析的数据集，调用sklearn的PCA接口对数据集执行均值归一化后确定其协方差矩阵，计算协方差矩阵的特征向量和特征值，特征向量按照特征值的降序排序；为差值图像的每个像素构造向量，构成差值图像的特征向量空间，将差值图像的特征向量空间投影到特征向量上。

进一步地，对分析后的所述差值图像进行聚类分析包括：

对投影后的差值图像的特征向量空间进行k-means算法聚类，在聚类的过程中，手动设置聚类中心数，随机选取相应个数的聚类中心，计算样本点到各聚类中心点的距离，根据所述距离，确定各个样本点所属类别，计算同类样本点的中心点，并将其设定为新的中心，重复此步骤，直到收敛。

另一方面，提供了一种遥感影像变化识别系统，所述系统包括：

裁剪模块：用于裁剪前后两期对应边界区域的栅格数据影像；

第一分析模块：用于对所述栅格数据影像进行主成分分析，得到灰度图像；

求差模块：用于对所述灰度图像进行缩放，使得大小完全一致，并求差值；

第二分析模块：用于对所述差值图像进行主成分分析，增加区分度；

聚类分析模块：用于对分析后的所述差值图像进行聚类分析；

二值图像生成模块：用于选定阈值，挑选出高于所述阈值中心点所涵盖的所有样本点作为变化像素，其余为未变化像素，分别对所述变化像素和所述未变化像素赋值，得到二值图像；

转换模块：用于将所述二值图像转换为矢量图；

第一计算模块：用于统计所述矢量图中赋值得到所述变化像素的面积和所述矢量图的总面积，得到变化率。

进一步地，所述第一分析模块包括：

读取模块：用于GDAL分别读取对应边界区域的前后两期栅格数据影像的波段数据以及坐标信息，将波段数据读入数组；

变形和转置模块：用于对GDAL读取到的波段数组进行变形和转置，处理成主成分分析所需要的正确数据格式；

主成分分析模块：用于调用sklearn库提供的PCA接口对处理后的数组进行主成分分析，将原本的三维数据投影到一维特征上，去除图像冗余信息，将原本的多波段栅格数据影像转换为单波段的灰度图像。

进一步地，所述第二分析模块包括：

提取模块：用于提取所述差值图像的特征向量，构建差值图像的特征向量数据集；

均值归一化模块：用于调用sklearn的PCA接口对数据集执行均值归一化后确定其协方差矩阵；

第二计算模块：用于计算协方差矩阵的特征向量和特征值；

构造模块：用于为差值图像的每个像素构造向量，以此构成差值图像的特征向量空间；

投影模块：用于将差值图像的特征向量空间投影到特征向量上。

进一步地，所述聚类分析模块包括：

第三计算模块：用于随机选取相应个数的聚类中心，计算样本点到各聚类中心的距离；

类别确定模块：用于根据所述距离，确定各个样本点所属类别；

第四计算模块：用于计算同类样本点的中心点，并将其设定为新的中心。

本发明实施例提供的技术方案带来的有益效果是：本发明采用GDAL库读取遥感影像数据，使用主成分分析法，将多波段的遥感影像数据转换为灰度图像，利用k-means无监督聚类的原理，对同一地理空间下不同时间采集到的遥感影像进行分析、对比以及变化识别，最终计算出同一区域遥感影像在前后两个时间点内发生变化的比率，提高了影像变化识别的准确度，优化图像处理的速度。

附图说明

为了更清楚地说明本发明实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明实施例一提供的一种遥感影像变化识别方法的流程图；

图2是本发明实施例二提供的一种遥感影像变化识别系统的示意图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方式作进一步地详细描述。

实施例一

一种遥感影像变化识别方法，参见图1，该方法包括以下步骤：

步骤(1)：裁剪前后两期对应边界区域的栅格数据影像，具体包括：新建一个文件夹，以要素图斑为单位，从矢量数据图层中导出每个图斑所对应的矢量边界数据，每份矢量边界数据都以图斑ID为名称进行保存；新建两个文件夹，用每一份矢量边界数据分别裁剪前后两期对应边界区域的栅格数据影像，同一期的栅格数据影像存放在一个文件夹中。

步骤(2)：对所述栅格数据影像进行主成分分析，具体包括：GDAL分别读取对应边界区域的前后两期栅格数据影像的波段数据以及坐标信息，将波段数据读入数组；对GDAL读取到的波段数组进行变形和转置，处理成主成分分析所需要的正确数据格式，调用sklearn库提供的PCA接口对处理后的数组进行主成分分析，将原本的三维数据投影到一维特征上，去除图像冗余信息，将原本的多波段栅格数据影像转换为单波段的灰度图像。

需要说明的是，GDAL是读写大量的栅格空间数据格式的开源库，GDAL可以读取栅格影像的波段数据以及相关地理坐标信息参数。

步骤(3)：对灰度图像进行数组缩放，使得大小完全一致，之后逐像素求差值，得到前后两期影像的差值图像。

步骤(4)：对差值图像进行主成分分析，增加区分度，具体包括：提取差值图像的特征向量，构建差值图像的特征向量数据集，作为主成分分析的数据集，调用sklearn的PCA接口对数据集执行均值归一化后确定其协方差矩阵，计算协方差矩阵的特征向量和特征值，特征向量按照特征值的降序排序；为差值图像的每个像素构造向量，构成差值图像的特征向量空间，将差值图像的特征向量空间投影到特征向量上，有助于后续聚类分析将像素分为变化像素和未变化像素。

需要说明的是，主成分分析是一种基于变量协方差矩阵对数据进行压缩降维、去噪的有效方法。

步骤(5)：对分析后的差值图像进行聚类分析，具体包括：对投影后的差值图像的特征向量空间进行k-means算法聚类，在聚类的过程中，手动设置聚类中心数，随机选取相应个数的聚类中心，计算样本点到各聚类中心的距离，根据距离，确定各个样本点所属类别，计算同类样本点的中心点，并将其设定为新的中心，重复此步骤，直到收敛，其原则是：类内的距离应尽可能小，类外的距离应尽可能大，由于发生变化的区域在差值图像中的像素高于未变化的像素，故发生变化区域像素的聚类平均值更高，在二值图像中，像素的变化范围在0～255之间，因此，设定一个像素平均值的阈值，高于该阈值的聚类平均值所涵盖的像素则被视为发生变化的像素，其余为未变化像素，变化像素的像素值赋为255，显示白色，未变化像素的像素值赋为0，显示黑色，将原本的灰度图像转换为二值图像。

步骤(6)：将二值图像转换为矢量图，统计矢量图中所有白色像素所占面积和矢量图的总面积，计算出变化区域所占的百分比，得到变化率。

实施例二

参见图2，一种遥感影像变化识别系统，包括：裁剪模块、第一分析模块、求差模块、第二分析模块、聚类分析模块、二值图像生成模块、转换模块和第一计算模块；裁剪模块用于裁剪前后两期对应边界区域的栅格数据影像；第一分析模块：用于对栅格数据影像进行主成分分析，得到灰度图像，第一分析模块包括：读取模块：用于GDAL分别读取对应边界区域的前后两期栅格数据影像的波段数据以及坐标信息，将波段数据读入数组；变形和转置模块：用于对GDAL读取到的波段数组进行变形和转置，处理成主成分分析所需要的正确数据格式；主成分分析模块：用于调用sklearn库提供的PCA接口对处理后的数组进行主成分分析，将原本的三维数据投影到一维特征上，去除图像冗余信息，将原本的多波段栅格数据影像转换为单波段的灰度图像。

求差模块用于对灰度图像进行缩放，使得大小完全一致，并求差值；第二分析模块：用于对差值图像进行主成分分析，增加区分度，第二分析模块具体包括：提取模块：用于提取差值图像的特征向量，构建差值图像的特征向量数据集，均值归一化模块：用于调用sklearn的PCA接口对数据集执行均值归一化后确定其协方差矩阵，第二计算模块：用于计算协方差矩阵的特征向量和特征值，构造模块：用于为差值图像的每个像素构造向量，以此构成差值图像的特征向量空间，投影模块：用于将差值图像的特征向量空间投影到特征向量上。

聚类分析模块用于对分析后的差值图像进行聚类分析，聚类分析模块具体包括：第三计算模块：用于随机选取相应个数的聚类中心，计算样本点到各聚类中心的距离；类别确定模块：用于根据距离，确定各个样本点所属类别；第四计算模块：用于计算同类样本点的中心点，并将其设定为新的中心。

二值图像生成模块：用于选定阈值，挑选出高于阈值中心点所涵盖的所有样本点作为变化像素，其余为未变化像素，变化像素的像素值赋为255，未变化像素的像素值赋为0，得到二值图像；转换模块：用于将二值图像转换为矢量图；第一计算模块：用于统计矢量图中所有像素值为255的面积和总面积，得到变化率。

值得说明的是，本发明采用GDAL库读取遥感影像数据，使用主成分分析法，将多波段的遥感影像数据转换为灰度图像，利用k-means无监督聚类的原理，对同一地理空间下不同时间采集到的遥感影像进行分析、对比以及变化识别，最终计算出同一区域遥感影像在前后两个时间点内发生变化的比率，提高了影像变化识别的准确度，优化图像处理的速度。

上述本发明实施例序号仅仅为了描述，不代表实施例的优劣。

以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。