CN109885808A

CN109885808A - 一种近地表气象要素计算方法

Info

Publication number: CN109885808A
Application number: CN201811597123.2A
Authority: CN
Inventors: 王智伟; 杨栋; 徐兰静; 刘艳峰; 岳泓辰
Original assignee: Xian University of Architecture and Technology
Current assignee: Xian University of Architecture and Technology
Priority date: 2018-12-26
Filing date: 2018-12-26
Publication date: 2019-06-14
Anticipated expiration: 2038-12-26
Also published as: CN109885808B

Abstract

本发明公开了一种近地表气象要素计算方法，主要包括：采集目标区域的气象相关数据；将目标区域划分为A区域、B区域、C区域和D区域四个区域；利用ANN方法确定最佳输入变量组合；依据该最佳输入变量组合确定最佳ANN映射模型；计算B区域的气象要素计算数据；计算D区域的气象要素计算数据；A区域、C区域的台站气象要素数据和B区域、D区域的计算气象要素估算数据共同组成目标区域的近地表气象要素数据。本发明结合遥感方法和空间插值方法，有效降低了台站稀疏导致的计算误差，并解决了遥感盲区近地表气象要素数据计算的问题。

Description

一种近地表气象要素计算方法

技术领域

本发明属于遥感应用技术领域，具体涉及一种近地表气象要素计算方法。

背景技术

空气温湿度等近地表气象要素是衡量地球气候变化、水文循环等的重要指标，也是开展各类气象服务的基础。在人类生产生活方面，这些气象数据可用于指导农牧业生产、空气调节***运行。在科学研究方面，可用于生态***建模、流行性疾病传播预警研究及气候变化等领域。在防灾减灾方面，也可用于林火探测等。

然而，受经济发展水平及人口分布特征等因素影响，部分地区气象台站分布稀疏且均匀性较差。现有台站气象数据不足以指导农牧业生产和支撑当地流行性疾病传播预警、气候变化、林火探测等研究。故探索适合复杂地形、台站稀疏区域准确获取气象数据的方法是非常必要的。

传统的高分辨率气象要素估算主要是使用当地及邻近区域内气象站数据进行空间插值，以实现数据点尺度到面尺度的转换。空间插值主要分为三类：非地统计学、地统计学和两者相结合的方法。如反距离加权法、回归模型法、薄板样条函数法属于非地统计学方法。地统计学方法主要包括克里格法及其改进算法等。结合方法则有梯度距离平方反比法、回归克里格等。研究人员对各种方法的适用性做了大量研究，结果表明空间插值方法对于台站密集且分布均匀区域的气象要素插值准确度较高。但由于受到台站数量及分布等因素的影响，在台站稀疏且分布均匀性较差的区域，插值误差较大，难以满足小尺度气象分布特征的研究。因此，随着遥感技术的迅速发展，使用遥感数据估算气象要素的方法开始形成并逐渐发展，以克服传统空间插值方法的上述缺陷。

气象要素遥感估算方法与传统的空间插值方法不同。它使用遥感获取的地-气***参数数据来估算近地表气象要素。卫星遥感提供具有空间连续覆盖、全球尺度的多种地-气***参数产品，如地表温度T_s、大气总可降水量TPW、土地覆盖类型LCT等。然而受到现有技术手段的限制，现在还不能直接通过遥感直接获取近地表气象数据。目前主流的方法主要是基于地表特征与近地表气象要素的耦合关系。因此，现有气象数据遥感估算方法在台站稀疏区域准确度比空间插值方法更高。但对于复杂地形区域，现有遥感估算方法准确度仍存在可提升的空间。此外，云覆盖、积雪、霾等造成的遥感观测盲区也是近地表气象要素计算的一大难点。

发明内容

为解决现有气象要素估算方法中存在的上述局限，本发明提供了一种基于遥感数据的近地表气象要素计算方法，解决现有的方法准确度低的问题。

为了实现上述目的，本发明采用如下技术方案予以实现：

一种近地表气象要素计算方法，包括以下步骤：

步骤1、采集目标区域的气象相关数据，其中气象相关数据包括从目标区域气象台站及其邻近区域气象台站获取的台站气象要素数据、从遥感平台获取的目标区域及其邻近区域的遥感数据；

步骤2、将目标区域划分为四个区域，分别为：A区域：遥感数据可用的台站区域，B区域：遥感数据可用的非台站区域，C区域：遥感数据不可用的台站区域，D区域：遥感数据不可用的非台站区域；

步骤3、利用人工神经网络方法对A区域的遥感变量进行降维处理，获得降维后的遥感变量；

步骤4、将A区域与降维后的遥感变量对应的遥感数据作为输入，A区域的台站气象要素数据作为输出，利用人工神经网络方法构建遥感数据和台站气象要素数据之间的映射模型；

步骤5、将B区域中与降维后的遥感变量对应的遥感数据作为输入，利用步骤4确定的映射模型计算B区域的气象要素计算数据；

步骤6、将步骤2中A区域的台站气象要素数据、C区域的台站气象要素数据和步骤5获得的B区域的气象要素计算数据作为输入数据，使用空间插值法计算D区域的气象要素计算数据；

A区域、C区域的台站气象要素数据和B区域、D区域的气象要素计算数据共同组成目标区域的近地表气象要素数据。

进一步的，所述的步骤3中，将A区域的台站气象要素数据作为输出数据，A区域的遥感数据作为输入数据，假定人工神经网络映射模型的隐含层为某一确定节点数，然后不断增加可以最大程度降低该ANN映射模型误差的输入变量，直至所有输入变量筛选完成，选取ANN映射模型误差最小时的输入变量组合作为最佳输入变量组合，即降维后的遥感变量。

进一步的，所述的步骤4中，将步骤3确定的降维后的遥感变量对应的遥感数据作为输入数据，A区域的台站气象要素数据作为输出数据，分别对不同隐含层节点的人工神经网络映射模型进行训练，并计算训练结束后的每个人工神经网络映射模型的误差，选取误差最小时对应的人工神经网络映射模型作为映射模型。

进一步的，所述的步骤1中的气象相关数据还包括一些辅助变量数据，其中，辅助变量数据主要有经纬度、高程、地形起伏度或儒略日等，这些辅助变量数据和遥感数据组合作为ANN映射模型的输入变量。

进一步的，所述的步骤1中还包括对采集的气象相关数据进行预处理，具体为：对从遥感平台采集的遥感图像进行图像配准、几何纠正、正射纠正、裁剪、镶嵌、重采样至目标分辨率后，按照气象台站在遥感图像中的位置提取对应的遥感数据，并将台站气象要素数据、遥感数据进行匹配。

具体的，所述的步骤3中，在每次增加输入变量后，将增加变量后的遥感数据按照一定的数据量比例划分为训练集、测试集和验证集，使用训练集对ANN映射模型进行训练，使用验证集决定训练停止的时刻，使用测试集验证训练结束后的ANN映射模型的误差；最终得到多个误差，选取ANN映射模型误差最小时的输入变量组合作为最佳输入变量组合。

具体的，所述的步骤4中，将输入变量的数据比例划分为训练集、测试集和验证集，使用训练集对ANN映射模型进行训练，使用验证集决定训练停止的时刻，使用测试集验证训练结束后的ANN映射模型的误差；最终得到多个误差，选取最小时对应的ANN映射模型作为最佳ANN映射模型。

具体的，所述的步骤6中的空间插值法主要有反距离加权插值法、薄板样条函数法、克里格法或改进算法。

与现有技术相比，本发明的有益效果是：

本发明的方法将待测定气象要素的区域划分为不同的几个区域，再结合遥感方法和空间插值方法，即将人工神经网络(ANN)模型的估算结果作为空间插值方法的输入变量。有效降低了台站稀疏导致的插值误差，并且空间插值方法的应用有效避免了遥感盲区这一近地表气象要素数据计算问题。

附图说明

图1是本发明方法计算气象要素的流程(图中‘○’代表方法；‘’代表已知数据；‘□’代表估算数据；虚线箭头表示区域间模型或数据传递；实线箭头表示流程)。

图2是步骤3中最佳输入变量筛选流程图。

图3是川西藏区复杂地形及台站分布。

图4是实施例中ANN映射模型最佳输入变量筛选结果。

以下结合说明书附图和具体实施方式对本发明做具体说明。

具体实施方式

本发明通过人工神经网络(ANN)方法获取遥感数据可用(即晴朗天气)区域的气象要素数据，然后再以遥感数据可用区域估算的气象要素数据作为空间插值方法的输入，获取遥感数据不可用(即非晴朗天气)区域的气象要素数据，流程见图1。具体的，本发明的近地表气象要素计算方法主要包括以下步骤：

步骤1、对于待测定的区域，通过现场实测或者中国气象数据共享网络平台搜集的数据作为该区域的气象相关数据，其中气象相关数据至少包括从目标区域气象台站及其邻近区域(一般为目标区域1～3个经纬度之内的区域)气象台站获取的台站气象要素数据、从遥感平台获取的目标区域及其邻近区域(目标区域1～3个经纬度之内区域)的遥感数据，其中，台站气象要素数据一般指气温、相对湿度、风速等。遥感数据一般指地表温度、植被指数、地表反射率、大气可降水量、地表覆盖类型等，主要来源于MODIS、PSS、METOP、FY、NOAA、Landsat等遥感平台。具体根据实际需要选择所需的变量及其对应的数据。

进一步的，气象相关数据还可以包括一些辅助变量数据，其中，辅助变量数据主要有经纬度、高程、地形起伏度或儒略日等，以增加后期气象要素计算数据的准确性，在后期使用过程中，这些辅助变量数据和遥感数据组合作为ANN映射模型的输入变量。

为了进一步增加上述数据的准确性，需对上述数据进行预处理，具体为：对从遥感平台采集的遥感图像进行图像配准、几何纠正、正射纠正、裁剪、镶嵌、重采样至目标分辨率后，按照气象台站在遥感图像中的位置提取对应的遥感数据，并将台站数据、遥感数据以及辅助变量数据按照经纬度和观测日期匹配。

步骤2、将上述目标区域划分为四个区域，分别为：A区域：遥感数据可用的台站区域，B区域：遥感数据可用的非台站区域，C区域：遥感数据不可用的台站区域，D区域：遥感数据不可用的非台站区域；根据实际统计结果，一般不同区域面积之间存在B>D>>A>C的关系。

步骤3、利用ANN方法确定最佳输入变量组合，即对A区域的遥感变量进行降维处理；

将A区域的台站气象要素数据作为ANN映射模型的输出变量，A区域的遥感数据(或遥感数据和辅助数据的组合数据)作为ANN映射模型的输入变量，此时ANN映射模型的节点数不确定，假定ANN映射模型的隐含层为某一确定节点数(比如节点数为100)，通过ANN方法筛选出适合于输出变量计算的最佳输入变量组合，具体筛选过程为：(1)、首先选择可以估算台站气象数据的最佳输入变量；(2)、在已筛选的最佳输入变量的基础上增加一个其它变量，将该变量组合作为新的输入变量，计算ANN映射模型误差(如均方根误差RMSE和平均绝对误差MAE)；(3)、重复上一步骤(2)，直至所有变量筛选完成，然后选取ANN映射模型误差最小时的输入变量组合作为最佳输入变量组合。

在每次增加输入变量后，将增加变量后的遥感数据(或遥感数据和辅助数据的组合数据)按照一定的数据量比例划分为训练集、测试集和验证集，该比例关系使用者可自行选择，一般为70％、15％和15％。使用训练集对ANN映射模型进行训练，使用验证集决定训练停止的时刻，使用测试集验证训练结束后的ANN映射模型的性能；最终得到多个误差，选取ANN映射模型误差最小时的输入变量组合作为最佳输入变量组合；具体流程见图2所示，图2中，“X”为ANN输入变量数据集，“VN”为X中变量数目，“T”为存储上依次循环中集合X的中间变量，“K”为预测因子(即遥感变量和辅助变量的组合变量)变量数目，“Predictors’”为估算的气象要素的最佳变量组合。

步骤4、确定最佳ANN映射模型：

将步骤3确定的最佳输入变量对应的A区域遥感数据(或遥感数据和辅助数据的组合数据)作为ANN映射模输入数据，A区域的台站气象要素数据作为输出数据，分别对不同隐含层节点的ANN映射模型进行训练，如隐含层节点为可选择10～100(也可根据变量数目需求选择其他节点数范围)，总共91个节点，然后计算每个节点对应的ANN映射模型的误差，最后选取误差最小的ANN映射模型作为最佳ANN映射模型。

在每次ANN映射模型训练过程中，同上述步骤3记载的过程一样，将输入数据按照一定比例划分为训练集、测试集和验证集，使用训练集对ANN映射模型进行训练，使用验证集决定训练停止的时刻，使用测试集验证训练结束后的ANN映射模型的误差；最终得到多个误差，根据上述举例的节点数，得到91个误差。选取误差最小时对应的ANN映射模型作为最佳ANN映射模型，即确定了输入变量数目、模型的隐含节点数和输出变量数目，具体流程见图2所示。

另外，当最佳输入变量组合中部分变量不可用(如检测不到)时，应用其他可用变量组合重复上述步骤4，生成相应的最佳ANN映射模型。

步骤5、确定B区域的气象要素计算数据：

将B区域中与步骤3获得的最佳输入变量对应的遥感数据作为输入变量，使用步骤4确定的最佳ANN映射模型计算B区域的气象要素计算数据。

步骤6、将步骤2中A区域的台站气象要素数据、C区域的台站气象要素数据和步骤5获得的B区域的气象要素计算数据作为输入数据，使用空间插值法计算D区域的气象要素计算数据。本发明的空间插值法主要有反距离加权插值法、薄板样条函数法、克里格法或改进算法，可根据操作复杂度或准确度等选择合适的方法，优选反距离加权插值法。

最后，A区域、C区域的台站气象要素数据和B区域、D区域的计算气象要素估算数据共同组成目标区域的近地表气象要素数据。

以下给出本发明的具体实施例，需要说明的是本发明并不局限于以下具体实施例中，凡在本申请技术方案基础上做的等同变换均落入本发明的保护范围。

实施例

以川西藏区近地表日最高和日最低气温估算为例。如附图3所示，川西藏区位于中国四川省西部，由甘孜州、阿坝州、雅安市宝兴县、石棉县和凉山州冕宁县等地区组成。总面积约238000km²，位于东经97.2°～104.5°，北纬27.7°～34.4°之间。地形复杂，海拔771～6983m，岭谷高差在500米以上。人员主要集中在县城等少数区域，仅有16个具有连续气象观测数据的国家基准气象站。地域广阔、地形复杂、人员居住分散及台站分布稀疏是川西藏区的主要特征。北部为川西北高原，地势由西向东倾斜，分为西北方向的丘状高原和东北方向的高平原。山高谷深、河流众多是南部川西山地的特征，根据切割深浅可分为西南方向的高山原和南部中域的高山峡谷区。川西藏区属于高原亚温带气候区。

步骤1：数据搜集及预处理

1)台站气象要素数据

川西藏区及其邻近区域41个国家级基准站2012年1月至2012年12月的日值数据被用于本发明方法的验证。数据包括日最高气温T_a,max、日最低气温T_a,min。气象要素数据均来源于中国气象数据网的中国地面气候资料日值数据集(V3.0)。

2)遥感数据

本实施例使用四种MODIS数据产品：(1)地表温度(MOD11A1和MYD11A1)，(2)植被指数(MOD13A2)，(3)土地覆盖(MCD12Q1)，(4)陆地反射(MOD09GA)，所有数据均来自于美国国家航空航天局网站。遥感数据产品均进行图像配准、几何纠正、正射纠正、裁剪、镶嵌、重采样等处理，至1km×1km空间分辨率，以匹配地表温度产品的空间分辨率。通过台站坐标提取遥感图像中对应位置的数据。

遥感数据包括四种T_s(Terra日间地表温度T_{s,Terra,daytime}、Terra夜间地表温度T_{s,Terra,nighttime}、Aqua日间地表温度T_{s,Aqua,daytime}、Aqua夜间地表温度T_{s,Aqua,nighttime})、两种植被指数(近红外波段归一化差分植被指数NDVI、近红外波段增强型植被指数EVI)、土地覆盖类型LCT、七个波段的陆地反射数据(Band₁、Band₂、Band₃、Band₄、Band₅、Band₆、Band₇)和太阳天顶角SZA。由于遥感产品时间分辨率的差异，需对数据进行处理，以匹配地表温度和台站观测数据的时间分辨率(一天一次)。例如植被指数产品每16天获取一次数据。假设植被在16天之内未有变化，每天的植被指数等于观测周期(16天)内一次观测的植被指数。使用MODIS Reprojection Tool对遥感图像进行提取和重采样；使用ArcGIS对遥感图像进行分析、数据提取、空间插值等；使用MATLAB对数据进行分析和构建模型。

3)辅助变量数据

直接或者间接影响气温估算的变量也应被用于模型的构建，如经度Longitude、纬度Latitude、高程Elevation、儒略日JD、地形起伏度R等。气象站经、纬度数据来源于中国气象数据网；高程数据来源于地理空间数据云平台SRTMDEMUTM 90m分辨率数字高程数据产品，并经图像预处理过程至1km×1km空间分辨率。儒略日为1～366，分别对应于2012年366天。应用高程栅格数据计算地形起伏度。台站数据，遥感数据和辅助数据按照经纬度、儒略日进行匹配。

步骤2：区域划分：将川西藏区域分为四个区域，分别为：A区域：遥感T_s数据可用的台站区域，B区域：遥感T_s数据可用的非台站区域，C区域：遥感T_s数据不可用的台站区域，D区域：遥感T_s数据不可用的非台站区域。统计结果证实，区域B和D的面积之比为5.17：1。

步骤3：变量筛选：变量筛选的目的为选择适合于气温估算的ANN映射模型输入变量，详细流程见附图2。本实施例的筛选结果见附图4。当输入变量数VN小于8时，模型误差随输入变量数的增加而迅速下降；当VN大于8，小于12时，模型误差持续下降但速度较为缓慢；当VN大于12时，模型误差增加。故当VN为12时，模型性能达到最佳，即RMSE和MAE最小。故川西藏区气温估算时，ANN映射模型最佳输入变量组合包括：T_{s,Terra,nighttime}、T_{s,Aqua,nighttime}、JD、Latitude、Elevation、T_{s,Terra,daytime}、NDVI、T_{s,Aqua,daytime}、R、Band₆、Band₁、Band₅。可以发现Band₂、Band₃、Band₄、Band₇、EVI和LCT并未出现于最佳输入变量组合中。原因可能为：对于川西藏区而言，Band₂、Band₃、Band₄、Band₇、EVI和LCT反映的地表信息与NDVI或其它波段重复。SZA和Latitude被排除，可能是因为SZA和经度对气温分布的影响有限。

步骤4：模型构建：在步骤3筛选的最佳输入变量基础上，构建区域A遥感数据和台站气象要素数据之间的ANN映射模型。步骤3筛选的最佳输入变量对应的区域A的数据作为ANN的输入，台站气象要素数据作为ANN的输出。训练ANN以获取输入与输出变量的最佳映射关系，即为气象要素估算的ANN映射模型。应用筛选后的变量构建晴朗天气(遥感T_s数据可用)条件下估算T_a,max和T_a,min的ANN模型。由于四种地表温度数据可用性随时间变化而不同，因此采用可用的地表温度数据，辅以步骤3中筛选后的变量中除地表温度外的变量对应的数据构建不同条件下的近地表气温估算模型，模型汇总见表1。

表1.晴朗天气条件下估算近地表气温的ANN模型汇总

上述模型分别是可用地表温度组合对应的最佳ANN映射模型，在使用过程中，根据该地段区域能够获取的地表温度类型选择对应的模型。

步骤5：将步骤4表1中的ANN映射模型应用于区域B，以步骤3筛选的最佳输入变量对应的区域B的数据作为步骤4表1中ANN映射模型的输入来估算区域B近地表气象要素数据，即日最高温度和日最低温度。

步骤6：以A区域台站气象要素数据、C区域的台站气象要素数据和步骤5中B区域的气象要素计算数据作为输入，使用反距离加权插值法得到D区域的气象要素计算数据。结果见附图5。

图5中(a)、(b)、(c)、(d)分别为川西藏区在2012年2月4日(冬季)、5月14日(春季)、7月13日(夏季)、10月11日(秋季)的最高温度分布情况；图6中(a)、(b)、(c)、(d)分别为川西藏区在2012年2月4日(冬季)、5月14日(春季)、7月13日(夏季)、10月11日(秋季)的最低温度分布情况。

对本实施例对估算川西藏区的气温误差进行计算，如表1所示，B区域气温估算的MAE为1.20～1.63℃，应用区域A的台站气象数据、区域C的台站气象数据以及区域B的气象要素估算数据以交叉验证方式评价本发明中反距离加权插值方法的准确度，反距离加权插值方法的MAE为1.59～3.23℃。根据误差传递原则，D区域气温估算的MAE为2.79～4.83℃。由于D区域占总区域面积的16.2％，因此应用本发明方法估算目标区域气温的MAE为16.2％·(2.79～4.83)+83.8％·(1.20～1.63)＝(1.45～2.15)℃

本发明还通过直接应用台站气象数据及反距离加权插值方法对川西藏区的最高温度和最低温度进行估算，结果证实该方法的平均绝对误差MAE为1.47～8.51℃。可以看出，本发明方法相对于传统的反距离加权插值方法具有更高的准确度。

需要说明的是，本发明并不局限于上述实施例。例如，除上述案例外，本发明所提出的方法还可应用于空气湿度等多种近地表气象要素的估算。在本发明公开的技术方案的基础上，本领域的技术人员根据所公开的技术内容，不需要创造性的劳动就可以对其中的一些技术特征作出一些替换和变形，这些替换和变形均在本发明的保护范围内。

Claims

1.一种近地表气象要素计算方法，其特征在于，包括以下步骤：

2.如权利要求1所述的近地表气象要素计算方法，其特征在于，所述的步骤3中，将A区域的台站气象要素数据作为输出数据，A区域的遥感数据作为输入数据，假定人工神经网络映射模型的隐含层为某一确定节点数，然后不断增加可以最大程度降低该ANN映射模型误差的输入变量，直至所有输入变量筛选完成，选取ANN映射模型误差最小时的输入变量组合作为最佳输入变量组合，即降维后的遥感变量。

3.如权利要求1所述的近地表气象要素计算方法，其特征在于，所述的步骤4中，将步骤3确定的降维后的遥感变量对应的遥感数据作为输入数据，A区域的台站气象要素数据作为输出数据，分别对不同隐含层节点的人工神经网络映射模型进行训练，并计算训练结束后的每个人工神经网络映射模型的误差，选取误差最小时对应的人工神经网络映射模型作为映射模型。

4.如权利要求1所述的近地表气象要素计算方法，其特征在于，所述的步骤1中的气象相关数据还包括一些辅助变量数据，其中，辅助变量数据主要有经纬度、高程、地形起伏度或儒略日等，这些辅助变量数据和遥感数据组合作为ANN映射模型的输入变量。

5.如权利要求1所述的近地表气象要素计算方法，其特征在于，所述的步骤1中还包括对采集的气象相关数据进行预处理，具体为：对从遥感平台采集的遥感图像进行图像配准、几何纠正、正射纠正、裁剪、镶嵌、重采样至目标分辨率后，按照气象台站在遥感图像中的位置提取对应的遥感数据，并将台站气象要素数据、遥感数据进行匹配。

6.如权利要求1所述的近地表气象要素计算方法，其特征在于，所述的步骤3中，在每次增加输入变量后，将增加变量后的遥感数据按照一定的数据量比例划分为训练集、测试集和验证集，使用训练集对ANN映射模型进行训练，使用验证集决定训练停止的时刻，使用测试集验证训练结束后的ANN映射模型的误差；最终得到多个误差，选取ANN映射模型误差最小时的输入变量组合作为最佳输入变量组合。

7.如权利要求1所述的近地表气象要素计算方法，其特征在于，所述的步骤4中，将输入变量的数据比例划分为训练集、测试集和验证集，使用训练集对ANN映射模型进行训练，使用验证集决定训练停止的时刻，使用测试集验证训练结束后的ANN映射模型的误差；最终得到多个误差，选取最小时对应的ANN映射模型作为最佳ANN映射模型。

8.如权利要求1所述的近地表气象要素计算方法，其特征在于，所述的步骤6中的空间插值法主要有反距离加权插值法、薄板样条函数法、克里格法或改进算法。