CN113780307A

CN113780307A - 一种区域年度最大蓝绿空间信息提取方法

Info

Publication number: CN113780307A
Application number: CN202110974210.0A
Authority: CN
Inventors: 鲁军景; 孙雷刚; 左璐; 刘剑锋; 郝庆涛; 郭风华; 马晓倩; 柏会子; 黄亚云
Original assignee: Institute Of Geography Hebei Academy Of Sciences
Current assignee: Institute Of Geography Hebei Academy Of Sciences
Priority date: 2021-08-24
Filing date: 2021-08-24
Publication date: 2021-12-10

Abstract

本发明公开了一种区域年度最大蓝绿空间信息提取方法，包括：对区域同年度的所有遥感影像进行辐射定标和大气校正，得到区域多时相遥感影像结果；提取遥感影像有效范围；根据构建的去云指数、暗像元指数、NDVI、水体指数提取相应的专题信息；对单时相遥感影像进行去云和暗像元处理，并利用NDVI提取植被信息；在提取暗像元的基础上，利用海蓝、绿波段和水体指数逐步提取水体信息；对单时相植被和水体进行合成；对区域多时相蓝绿空间信息进行空间融合，最终得到区域年度最大蓝绿空间信息结果。本发明通过构建去云指数、暗像元指数和水体指数逐层剔除云层和阴影、季节差异等因素对植被和水体信息的影响，操作过程简单、灵活，易于在大尺度范围推广应用。

Description

一种区域年度最大蓝绿空间信息提取方法

技术领域

本发明涉及一种区域年度最大蓝绿空间信息提取方法，属于生态遥感、国土空间规划领域。

背景技术

坚持生态优先，推动绿色发展，是建设幸福宜居之城的必然趋势。蓝绿空间指区域范围内以提供生态系统服务为主的水体和绿地两种生态元素的总称，其建设质量直接影响着区域的生态宜居环境与生态建设水平。在生态、美丽、和谐、宜居的新时代城市发展建设的理念下，蓝绿空间是区域公共空间的主体和特色精华所在，为市民供给休闲游憩服务，改善提升城市环境风貌，并对优化城市空间结构、隔离防护各种危险和污染具有重要意义。近几年“国土开发格局”和“蓝绿空间”的要义在国家重要文件和重大规划中不断被深化。2018年《河北雄安新区规划纲要》明确要求雄安新区蓝绿空间占比稳定在70％；而且在国土规划新背景下，蓝绿空间规划已作为重要支撑纳入市县级国土空间规划体系，体现了国土空间规划体系对建设“生态文明、美丽中国”的高度重视。但随着城镇化、工业化进程的不断加快，2020年末中国常住人口城镇化率已经超过60％；在城镇扩展过程中的土地利用/土地覆盖变化(LUCC)客观上表现出生产、生活性用地与生态性用地间的双向流动特征，加剧了蓝绿生态格局的演变，使得蓝绿生态要素受损，引发环境污染，水土流失，自然生境侵蚀等系列生态问题，严重制约区域生态宜居环境建设。因此，急需快速、全面、准确摸清“蓝绿”家底，掌握区域蓝绿位置、数量等空间基础资料，为蓝绿空间系统规划与日常管理决策提供依据，确保区域范围内蓝绿空间的布局更加科学合理，进而推动区域的高质量可持续发展。

长期以来，区域蓝绿空间分布主要依靠人工进行地面实测，这种解译方式直接利用了人脑的智慧和知识，解译结果一般具有很高的可靠性，但费时、费力，解译成本高，而且定量化的分析精度较差，不适合大尺度研究。遥感技术的快速发展，为更加精确地获取“蓝绿”空间数据提供了现实的技术支撑。2010年之前的植被和水体提取受限于遥感技术与影像质量的影响，只能进行简单基于像素波段间的计算；Rouse等利用红光和和近红外波段构建了归一化植被指数NDVI，Gao等根据植被指数的构建原理，基于TM影像利用绿光波段和近红外波段构建了归一化差异水体指数NDWI，这两个指数均取得了较好的效果，且广泛应用于生态研究的各个领域。期间随着计算机技术的发展，以阈值法为基础基于机器学习的决策树和支持向量机SVM方法陆续应用到植被和水体提取中。近十年来，随着三高遥感数据的出现以及植被和水体遥感定量与反演研究的深入，面向对象的技术逐步发展起来；崔齐等提出一种基于面向对象的矢量约束来实现高空间分辨率遥感影像水体提取的方法，该方法能够准确提取细小水体信息；蒲智基于QuickBird遥感数据，采用面向对象的技术对遥感影像进行阴影检测和补偿，然后对阴影校正后的影像进行绿地信息提取；面向对象的提取方法抗噪性强，但更适合地物纹理突出的高分辨率遥感影像。从2015年开始植被和水体提取的研究越来越广泛，其他技术相继出现，其中具有代表性的是基于机器学习的深度学习方法在植被和水体中的应用。但新型提取方法主要针对特殊地形或精度要求偏高的区域，算法的可泛化性较低，目前植被和水体的提取仍以阈值法和分类器法为主要识别方法。

目前大多数的“蓝绿”空间识别都是基于成像质量很高的单时相影像数据进行的单一要素的研究，而实际应用中，受天气、地形、传感器本身的影响，云层、阴影以及“同物异谱和异物同谱”现象都会影响地物的识别精度；地表植被和水体也会有季节性差异，因此需要充分利用多时相数据，准确获取区域的最大蓝绿空间分布；而且在实际城市建设过程中长久以来的蓝绿分治、水岸分建，造成了蓝绿断裂、水岸分割等诸多问题，使得城市生态的完整性方面存在“空缺”，相互融合不足，生态魅力特色难以彰显。

发明内容

针对现有研究和技术存在的不足，本发明提出了一种区域年度最大蓝绿空间信息提取方法，充分的考虑了地表植被和水体的季节性变化，以及云层、阴影和“同物异谱和异物同谱”现象的干扰，最大限度的挖掘影像数据中的植被和水体信息，该方法操作简单、易于大尺度推广应用。

为解决上述问题，本发明所采取的技术方案是：

一种区域年度最大蓝绿空间信息提取方法，包括以下步骤：

S1，对同一年度内植被生长季和水体涨水季的遥感影像进行预处理(辐射定标和大气校正)，得到不同时相代表地物真实反射率的遥感影像结果图；

S2，针对步骤S1所得遥遥感影像结果图，以任一时相遥感影像结果图为基准，对其他遥感影像结果图进行精确的几何配准，形成一个在地理空间坐标上相匹配的遥感影像结果图集合；

S3，针对步骤S2中得到的遥感影像结果图集合中的每一时相遥感影像结果图分别进行以下操作：

S3-1，针对该时相遥感影像结果图分别计算去云指数、暗像元指数、归一化植被指数和水体指数，得到相应的指数空间分布结果；

S3-2，对该时相遥感影像结果图进行波段组合和拉伸显示；

S3-3，对该时相遥感影像结果图进行植被覆盖区提取，详细操作步骤如下：

S3-3-1，对步骤S3-1中的归一化植被指数NDVI进行去云处理，操作步骤如下：

S3-3-1-1，在该时相遥感影像结果图云区随机选取样本点，统计样本点在去云指数值的分布情况，设定最小值为下限阈值，并生成0-1二值掩膜图，0代表其他区域，1代表初步提取的云区；

S3-3-1-2，在步骤S3-3-1-1初步提取的云区范围内随机选取样本点，统计样本点在热红外波段的分布情况，设定最大值为上限阈值，最后生成0-1二值掩膜图，0代表地表温度较高的区域，1代表云区地表温度较低的区域；

S3-3-1-3，将步骤S3-3-1-1和S3-3-1-2中的两个二值掩膜进行空间运算，得到最终的云覆盖区0-1二值图，0代表云区，1代表其他区域；

S3-3-2，对步骤S3-1中的归一化植被指数NDVI进行去暗像元处理，操作步骤如下：

S3-3-2-1，在该时相遥感影像结果图暗像元区(水域、阴影)随机选取样本点，统计样本点在暗像元指数值的分布情况，设定最大值为上限阈值，最后生成0-1二值掩膜图，0代表暗像元区，1代表其他区域；

S3-3-3，将步骤S3-3-1-3和S3-3-2-1中二值图与S3-1中的归一化植被指数NDVI进行单时相空间信息运算，最后得到消除云和暗像元噪声的归一化植被指数NDVI；

S3-3-4，在该时相遥感影像结果图植被覆盖区随机选取样本点，统计样本点在步骤S3-3-3中去除噪声后的NDVI值的分布情况，设定最小值为下限阈值，最终得到消除噪声后的植被分布结果二值图，0代表其他区域，1代表植被覆盖区；

S3-4，利用决策树方法对该时相遥感影像结果图进行提取水体处理，操作步骤如下：

S3-4-1，对步骤S3-3-2-1中的二值掩膜图进行空间运算，生成新的0-1二值掩膜图，0代表其他区域，1代表暗像元区；

S3-4-2，在步骤S3-4-1的基础上，对海蓝波段和绿波段进行阈值设定，去除云阴影覆盖下城镇用地和山体阴影，最后生成0-1二值掩膜图，0代表剔除云阴影覆盖下城镇用地和山体阴影，1代表其他暗像元区；

S3-4-3，对步骤S3-4-2中提取的其他暗像元区范围内的水体覆盖区域随机选取样本点，统计样本点在水体指数值的分布情况，设定最小值为下限阈值；最后生成0-1二值掩膜图，0代表其他区域，1代表初步提取的水体；

S3-4-4，将步骤S3-4-1、S3-4-2和S3-4-3中的三个二值掩膜进行空间运算，得到最终的水体区域0-1二值图，0代表其他区域，2代表水体；

S3-5，对步骤S3-3中植被二值分布图和步骤3-4中的水体二值分布图进行单时相信息空间合成，生成该时相最终的蓝绿空间信息结果图；

S4，基于步骤S3所得各时相的蓝绿空间信息结果图，将区域内年度所有时相的最终蓝绿空间信息结果图进行融合，形成一个包含多个波段的影像数据；

S5，在步骤S4基础之上，对覆盖相同区域的最终蓝绿空间信息结果进行最大值合成，得到该年度区域最大的蓝绿空间信息分布图。

进一步的，步骤S2还包括针对遥感影像结果图集合所有波段提取重合范围，得到该遥感影像结果图集合有效的遥感影像数据。

进一步的，步骤S3-1中构建的去云指数QYI，其计算公式为：QYI＝R_B+R_G+R_R，R_B为蓝波段的光谱反射率值，R_G为绿波段的光谱反射率值，R_R为红波段的光谱反射率值。

进一步的，步骤S3-1中构建的暗像元指数DPI，其计算公式为：DPI＝R_NIR+R_SWIR1+R_SWIR2，R_NIR为近红外波段的光谱反射率值，R_SWIR1为1.57-1.65μm波长范围的短波红外光谱反射率值，R_SWIR2为2.11-2.29μm波长范围的短波红外光谱反射率值。

进一步的，步骤S3-1中构建的归一化植被指数NDVI，其计算公式为：

RNIR为近红外波段的光谱反射率值，R_R为红波段的光谱反射率值。

进一步的，步骤S3-1中构建的水体指数WI，其计算公式为：

其中，

R_G为绿波段的光谱反射率值，R_MIR为中红外波段的光谱反射率值，R_NIR为近红外波段的光谱反射率值。

进一步的，步骤S3-3-1-3中的空间运算的计算公式为：YC_掩膜＝((YC1*TIRS1)eq0)*1，其中YC_掩膜为云区二值掩膜图，YC1是步骤S3-3-1-1中初步识别的云区二值掩膜图，TIRS1是步骤S3-3-1-2中的地表温度分区的二值掩膜图。

进一步的，步骤S3-3-3中的空间运算的计算公式为：NDVI_去噪＝NDVI*YC_掩膜*DPI_掩膜，其中NDVI_去噪为消除云层和暗像元噪声后的植被指数信息结果，NDVI是归一化植被指数，YC_掩膜是步骤S3-3-1-3中的云区二值掩膜图，DPI_掩膜是步骤S3-3-2-1中的暗像元二值掩膜图。

进一步的，步骤S3-4-4中的空间运算的计算公式为：Water_掩膜＝DPI_掩膜1*YY1*W11，其中Water_掩膜为水体二值掩膜图，DPI_掩膜1是步骤S3-4-1中的暗像元二值掩膜图，YY1是步骤S3-4-2中剔除云阴影覆盖下城镇用地和山体阴影的暗像元区的二值图，WI1是步骤S3-4-3中初步提取的水体的二值掩膜图。

进一步的，步骤S3-5中的单时相信息空间合成方法的计算公式为：Blue_Green_单时相＝Veg_掩膜+Water_掩膜，其中Blue_Green_单时相为单时相蓝绿空间信息结果图，Veg_掩膜是步骤S3-3-4中的植被分布二值掩膜图，Water_掩膜是步骤S3-4-4中的水体二值掩膜图。

本发明的有益效果在于：

(1)本发明考虑了云层、阴影等因素，构建了去云指数、暗像元指数和水体指数，消除了云、阴影等噪声对区域蓝绿空间分布提取结果的影响，实现了区域年度最大蓝绿空间信息提取，操作过程简单、灵活，易于在大尺度范围椎广应用。

(2)本发明通过QYI和热红外波段组合有效剔除城镇高亮地物与云层信息存在的“异物同谱”现象，提高了云层的提取结果精度。

(3)本发明通过DPI、海蓝、绿波段和WI的组合利用决策树的方法有效剔除了山体阴影、云阴影及其与城镇等地物共同作用后同水体存在的“异物同谱”现象，提前高了水体的提取结果精度。

(4)相比现有的技术方法，本发明从“蓝绿”结合的角度出发，充分考虑了植被和水体的季节性差异，利用遥感影像数据的多时相特征，提出了一种区域年度最大蓝绿空间信息提取方法，解决了长久以来的蓝绿断裂现象，弥补了区域年度最大蓝绿空间提取的空白，可为快速获取国土空间规划中的蓝绿空间占比提供技术支撑。

附图说明

图1为本发明的流程图；

图2为四个时相的遥感影像图；

图3为三个时相的云层提取结果图；

图4为四个时相的暗像元提取结果图；

图5为四个时相经过去云和去阴影处理后的最终的蓝绿空间信息结果图；

图6为区域年度最大蓝绿空间信息分布图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。

本实施例以某区域2020年度5月6日、5月22日、8月10日和9月11日四个不同时相的Landsat8遥感影像为例进行说明。由于受恶劣天气条件影响，部分遥感影像数据出现缺失，但考虑到短时间植被和水体不会突然改变，可以用邻近时相的数据进行补充，因此本实例使用上述四个不同时相的数据提取本区域最大年度蓝绿空间信息。

参考图1，图1为本发明一种区域年度最大蓝绿空间信息提取方法的流程图，包括以下步骤：

S1，收集该区域同一年度植被生长季和水体涨水季内所有遥感影像，分别对其进行精确的辐射定标和大气辐射校正，得到年度内所有时相的代表地物真实反射率的遥感影像结果图。如图2为该区域四个时相的遥感影像图，图2中第一幅图为5月6日的遥感影像图，第二幅图为5月22日的遥感影像图，第三幅图为8月10日的遥感影像图，第四幅图为9月11日的遥感影像图；除了季节性的差异外，同时都有不同程度的云层和阴影等背景噪声的影响，任一单时相遥感数据提取结果均不能完全代表该区域本年度的蓝绿空间分布状况。

S2，针对步骤S1所得遥遥感影像结果图，以任一时相遥感影像结果图为基准，对其他遥感影像结果图进行精确的几何配准，形成一个在地理空间坐标上相匹配的遥感影像结果图集合。由于图中有部分是不重合，而为了保证本发明的准确性，还需要针对遥感影像结果图集合所有波段提取重合范围，得到该遥感影像结果图集合有效的遥感影像数据。在本实施例中，利用ENVI软件小熊工具箱提取所有波段重合范围，得到每一时相有效的遥感影像数据。

S3-1，针对该时相遥感影像结果图分别计算去云指数、暗像元指数、归一化植被指数和水体指数，得到相应的指数空间分布结果。

卫星遥感在成像的过程经常会受到天气的影响，云层遮挡等问题会造成遥感影像上地物信息的损失，而且在提取植被时云层还会导致植被指数出现异常值等问题。本发明在分析云层与植被、水体、城镇、裸土等典型地物的光谱差异的基础上，利用云层的高亮度特性构建去云指数QYI，计算公式为：QYI＝R_B+R_G+R_R，R_B为蓝波段的光谱反射率值，R_G为绿波段的光谱反射率值，R_R为红波段的光谱反射率值。

云层除了本身对地物信息的遮挡外，其阴影也会造成地物信息的损失，而且阴影常与水体在影像上显示低亮度的暗像元存在“异物同谱”现象；考虑到阴影和水体整体暗像元提取对植被的提取精度不会造成影响，还可以为水体的提取提供初步的范围，因此本发明从暗像元的角度将水体和阴影作为整体构建了暗像元指数DPI，其计算公式为：DPI＝R_NIR+R_SWIR1+R_SWIR2，R_NIR为近红外波段的光谱反射率值，R_SWIR1为1.57-1.65μm波长范围的短波红外光谱反射率值，R_SWIR2为2.11-2.29μm波长范围的短波红外光谱反射率值。

归一化植被指数NDVI已广泛应用于植被遥感监测的各个领域，其最基础的特性就是可以根据NDVI值的大小提取植被覆盖区域，因此本发明采用NDVI阈值法对植被覆盖区域进行提取，其计算公式为：

水体是遥感影像上的典型地物，但容易与云阴影和山体阴影出现“异物同谱”现象，从而无法将阴影和水体较好的分开。目前水体的识别常使用NDWI和MNDWI两个指数，MNDWI利用短波红外代替近红外波段对NDWI进行了改进，但在本实验中两者均未取得很好的效果，因此，本发明在MNDWI的基础上添加了近红外波段构建了水体指数WI，其计算公式为：

其中，

R_G为绿波段的光谱反射率值，R_MIR为中红外波段的光谱反射率值，R_NIR为近红外波段的光谱反射率值。在本实施例中MNDWI公式中的R_MIR使用的是Landsat8数据中的SWIR1(1.57-1.65μm)波段的反射率值。

S3-2，对该时相遥感影像结果图进行NIR-Red-Green波段组合和Logarithmic拉伸显示。

S3-3-1-1，在该时相遥感影像结果图云区随机选取样本点，统计样本点在去云指数值的分布情况，设定最小值为下限阈值，并生成0-1二值掩膜图，0代表其他区域，1代表初步提取的云区。

S3-3-1-2，去云指数QYI主要考虑了云的高亮度特性，但实际提取的云区常混有城镇高亮区，从云区和城镇地表温度差异的角度出发，在步骤S3-3-1-1初步提取的云区范围内随机选取样本点，统计样本点在热红外波段的分布情况，设定最大值为上限阈值，最后生成0-1二值掩膜图，0代表地表温度较高的区域，1代表云区地表温度较低的区域；对云层覆盖区进行提纯，去除城镇高亮区异物同谱现象的影响。

S3-3-1-3，将步骤S3-3-1-1和S3-3-1-2中的两个二值掩膜进行空间运算，得到最终的云覆盖区0-1二值图，0代表云区，1代表其他区域。其空间运算的计算公式为：

YC_掩膜＝((YC1*TIRS1)eq 0)*1

其中YC_掩膜为云区二值掩膜图，YC1是步骤S3-3-1-1中初步识别的云区二值掩膜图，TIRS1是步骤S3-3-1-2中的地表温度分区的二值掩膜图。

S3-3-2-1，在该时相遥感影像结果图暗像元区(水域、阴影)随机选取样本点，统计样本点在暗像元指数值的分布情况，设定最大值为上限阈值，最后生成0-1二值掩膜图，0代表暗像元区，1代表其他区域。

S3-3-3，将步骤S3-3-1-3和S3-3-2-1中二值图与S3-1中的归一化植被指数NDVI进行单时相空间信息运算，最后得到消除云和暗像元噪声的归一化植被指数NDVI。其空间运算的计算公式为：

NDVI_去噪＝NDVI*YC_掩膜*DPI_掩膜

其中NDVI_去噪为消除云层和暗像元噪声后的植被指数信息结果，NDVI是归一化植被指数，YC_掩膜是步骤S3-3-1-3中的云区二值掩膜图，DPI_掩膜是步骤S3-3-2-1中的暗像元二值掩膜图。

S3-3-4，在该时相遥感影像结果图植被覆盖区随机选取样本点，统计样本点在步骤S3-3-3中去除噪声后的NDVI值的分布情况，设定最小值为下限阈值，最终得到消除噪声后的植被分布结果二值图，0代表其他区域，1代表植被覆盖区。

云层和暗像元提取结果是植被和水体精确提取的基础，由于5月22日遥感影像没有云层遮挡，因此最终给出了三个时相遥感影像的云层(图3)和四个时相遥感影像的暗像元提取结果图(图4)。如图3所示，依次为5月06日、8月10日和9月11日云层的提取结果。如图4所示，依次为5月06日、5月22日、8月10日和9月11日暗像元的提取结果。

S3-4-1，对步骤S3-3-2-1中的二值掩膜图进行空间运算，生成新的0-1二值掩膜图，0代表其他区域，1代表暗像元区。

S3-4-2，在步骤S3-4-1的基础上，对海蓝波段和绿波段设置上限阈值0，利用公式(b1 ge 0or b2 ge 0)*1去除云阴影覆盖下城镇用地和山体阴影，最后生成0-1二值掩膜图，0代表剔除云阴影覆盖下城镇用地和山体阴影，1代表其他暗像元区。

S3-4-3，对步骤S3-4-2中提取的其他暗像元区范围内的水体覆盖区域随机选取样本点，统计样本点在水体指数值的分布情况，设定最小值为下限阈值；最后生成0-1二值掩膜图，0代表其他区域，1代表初步提取的水体。

S3-4-4，将步骤S3-4-1、S3-4-2和S3-4-3中的三个二值掩膜进行空间运算，得到最终的水体区域0-1二值图，0代表其他区域，2代表水体。其空间运算的计算公式为：

Water_掩膜＝DPI_掩膜1*YY1*WI1

其中Water_掩膜为水体二值掩膜图，DPI_掩膜1是步骤S3-4-1中的暗像元二值掩膜图，YY1是步骤S3-4-2中剔除云阴影覆盖下城镇用地和山体阴影的暗像元区的二值图，WI1是步骤S3-4-3中初步提取的水体的二值掩膜图。

S3-5，对步骤S3-3中植被二值分布图和步骤3-4中的水体二值分布图进行单时相信息空间合成，生成该时相最终的蓝绿空间信息结果图。其单时相信息空间合成方法的计算公式为：

Blue-Green_单时相＝Veg_掩膜+Water_掩膜

其中Blue_Green_单时相为单时相蓝绿空间信息结果图，Veg_掩膜是步骤S3-3-4中的植被分布二值掩膜图，Water_掩膜是步骤S3-4-4中的水体二值掩膜图。

图5为四个单时相经过去云和去阴影处理后的最终的蓝绿空间信息结果图，其中图5中第一幅图为5月06日的区域蓝绿空间信息结果图，第二幅图为5月22日的区域蓝绿空间信息结果图，第三幅图为8月10日的区域蓝绿空间信息结果图，第四幅图为9月11日的区域蓝绿空间信息结果图。

S4，基于步骤S3所得四个单时相的蓝绿空间信息结果图，将区域内年度四个单时相的最终蓝绿空间信息结果图进行融合，形成一个包含多个波段的影像数据。在本实例中使用ENVI软件中的Layer Stacking工具将步骤S3所得四个单时相的蓝绿空间信息结果图合成得到一个在空间上相匹配的影像。

S5，在步骤S4基础之上，对覆盖相同区域的最终蓝绿空间信息结果进行最大值合成，得到该区域年度最大的蓝绿空间信息分布图。在本实例中使用ENVI软件中的band math工具利用公式“b1 ge b2 ge b3 ge b4”对步骤S4中的影像数据的每个像元进行空间运算，最终得到该区域年度最大的蓝绿空间信息分布图，如图6所示。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围。

Claims

1.一种区域年度最大蓝绿空间信息提取方法，其特征在于，包括以下步骤：

S1，对同一年度内植被生长季和水体涨水季的遥感影像进行预处理，得到不同时相代表地物真实反射率的遥感影像结果图；

S3-2，对该时相遥感影像结果图进行波段组合和拉伸显示；

S3-3-2-1，在该时相遥感影像结果图暗像元区随机选取样本点，统计样本点在暗像元指数值的分布情况，设定最大值为上限阈值，最后生成0-1二值掩膜图，0代表暗像元区，1代表其他区域；

2.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S2还包括针对遥感影像结果图集合所有波段提取重合范围，得到该遥感影像结果图集合有效的遥感影像数据。

3.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-1中构建的去云指数QYI，其计算公式为：QYI＝R_B+R_G+R_R，R_B为蓝波段的光谱反射率值，R_G为绿波段的光谱反射率值，R_R为红波段的光谱反射率值。

4.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-1中构建的暗像元指数DPI，其计算公式为：DPI＝R_NIR+R_SWIR1+R_SWIR2，R_NIR为近红外波段的光谱反射率值，R_SWIR1为1.57-1.65μm波长范围的短波红外光谱反射率值，R_SWIR2为2.11-2.29μm波长范围的短波红外光谱反射率值。

5.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-1中构建的归一化植被指数NDVI，其计算公式为：

R_NIR为近红外波段的光谱反射率值，R_R为红波段的光谱反射率值。

6.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-1中构建的水体指数WI，其计算公式为：

其中，

_G为绿波段的光谱反射率值，R_MIR为中红外波段的光谱反射率值，R_NIR为近红外波段的光谱反射率值。

7.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-3-1-3中的空间运算的计算公式为：YC_掩膜＝((YC1*TIRS1)0*1，其中YC_掩膜为云区二值掩膜图，YC1是步骤S3-3-1-1中初步识别的云区二值掩膜图，TIRS1是步骤S3-3-1-2中的地表温度分区的二值掩膜图。

8.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-3-3中的空间运算的计算公式为：NDVI_去噪＝NDVI*YC_掩膜*PI_掩膜，其中NDVI_去噪为消除云层和暗像元噪声后的植被指数信息结果，NDVI是归一化植被指数，YC_掩膜是步骤S3-3-1-3中的云区二值掩膜图，DPI_掩膜是步骤S3-3-2-1中的暗像元二值掩膜图。

9.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-4-4中的空间运算的计算公式为：Water_掩膜＝DPI_掩膜1*YY1*WI1，其中Water_掩膜为水体二值掩膜图，DPI_掩膜1是步骤S3-4-1中的暗像元二值掩膜图，YY1是步骤S3-4-2中剔除云阴影覆盖下城镇用地和山体阴影的暗像元区的二值图，WI1是步骤S3-4-3中初步提取的水体的二值掩膜图。

10.如权利要求1所述的一种区域年度最大蓝绿空间信息提取方法，其特征在于，步骤S3-5中的单时相信息空间合成方法的计算公式为：Blue__单时相＝Veg_掩膜+Water_掩膜，其中Blue__单时相为单时相蓝绿空间信息结果图，Veg_掩膜是步骤S3-3-4中的植被分布二值掩膜图，Water_掩膜是步骤S3-4-4中的水体二值掩膜图。