[发明专利]一种基于GlobeLand30的历史地表覆盖快速重建方法

权利要求书

1.一种基于GlobeLand30的历史地表覆盖快速重建方法，包括以下步骤：

第一步、准备地表覆盖产品GlobeLand30、遥感影像Landsat数据以及DEM数据，选择感兴趣区域并生成感兴趣区域的矢量文件，对GlobeLand30、DEM与遥感影像Landsat数据进行裁剪；

第二步、基于形态学开运算构建对地表覆盖产品土地斑块的优化从而获得斑块单元：将地表覆盖产品GlobeLand30按照类别划分为若干二值图层，同时对每一层二值图像进行窗口大小递增的形态学开运算，针对不同窗口大小下的计算结果，当前背景像元占比大于用户设定阈值时，对二值图层进行合并，获得优化后的斑块单元；

第三步、计算影像归一化差值光谱向量NDSV：通过遍历所有波段反射率的两两组合，对每一对波段组合计算归一化差值指数，获得对应的归一化差值光谱向量NDSV，所有的归一化差值光谱向量NDSV构成了影像光谱特征集；

第四步、将每一个斑块单元作为聚类计算单元，采用K-means算法进行无监督聚类，采用K-means++算法逐个初始化每个斑块单元的聚类中心：对于任意斑块单元，从中随机选择一个像元的归一化差值光谱向量NDSV作为第一个簇的聚类中心，当目前有k个聚类中心时，计算当前斑块单元内其他像元的归一化差值光谱向量NDSV到第k个聚类中心μ_k的欧式距离，并根据该欧式距离设置对应像元的抽样权重，根据抽样权重随机抽取第k+1个聚类中心，重复该过程直到簇的数量满足用户预设要求；

第五步、采用方差比准则VRC进行每个斑块单元聚类簇数的寻优：定义一个聚类簇数集合，对于任意斑块单元，将所有像元的归一化差值光谱向量NDSV输入K-means算法进行聚类，根据聚类结果计算方差比准则VRC结果，通过遍历簇数集合内的所有元素，获得最大方差比准则VRC对应的簇数K_j作为当前斑块单元内聚类的最优簇数；

第六步、根据第四、第五步方法，对所有斑块单元执行K-means方法，根据最优簇数聚类结果获得全局影像伪样本集：将第三步获得的光谱影像特征集作为K-means输入特征，遍历第二步获得的所有斑块单元，对各个斑块单元内的所有像元使用第四步方法获得聚类中心，采用第五步方法获得最优簇数，在此基础上获得聚类结果，令斑块单元内像元数量占比最大的簇继承斑块单元的类别信息，获得伪样本集；

第七步、根据地表覆盖产品确定地表覆盖类别体系，对每一个地类类别选取对应的关键特征，构建多个二类高斯混合模型：计算全局影像的改进归一化水体指数MNDWI作为水体类型的关键特征，计算非水体区域的归一化植被指数NDVI分别作为林地、草地、人工地表类型的关键特征，计算植被区域的坡度作为耕地类型的关键特征，将每一个关键特征的特征值分布视为两个高斯分布的混合，构建如下高斯混合模型：

其中为概率密度函数，y表示关键特征的特征值，α₁、α₂分别为第一和第二个高斯分布的混合系数、μ₁、μ₂分别为两个分布的均值向量，∑₁、∑₂分布为两个分布的协方差向量，p(y|μ₁,∑₁)与p(y|μ₂,∑₂)分别为第一与第二个高斯分布密度函数。以上{(α_j,μ_j,∑_j)|j＝1,2}统称为混合成分参数；

第八步、采用迭代期望最大化EM算法估算每个高斯混合模型中的混合成分参数，完成对高斯混合模型的分解，从而完成对伪样本集的优化，获得训练样本集；

第九步、利用变化检测与分类方法完成目标历史时相的地表覆盖重建：判断目标时相与地表覆盖产品是否属于同一年份，如果属于同一年份，则将第八步中优化后的训练样本输入随机森林RF分类器完成训练，获得目标时相地表覆盖分类结果；如果属于不同年份，则采用双时相影像变化检测方法获得不变区域，传递不变区域样本标签，结合目标时相影像形成训练样本，训练样本输入随机森林RF分类器完成训练，获得目标历史时相地表覆盖分类结果。

2.根据权利要求1所述的基于GlobeLand30的历史地表覆盖快速重建方法，其特征在于：第二步所述的土地斑块优化方法，其中形态学结构元形状为正方形，领域大小范围为[3,21],步长为2，背景像素比例ξ设置阈值为5％。