一种基于DEM和NDVI的气温降尺度方法

技术领域

本发明涉及气温降尺度方法，具体的说是涉及一种基于DEM(数字高程模型)和NDVI(归一化差分植被指数)的气温降尺度方法。

背景技术

近地面气温(一般为距地表2米高度处)在区域与全球尺度上的能量平衡和水循环平衡研究中是一个关键的物理参数，是地表能量平衡及陆表过程模拟的重要输入参数，已被广泛用于地表蒸散发估算、农业监测，气候变化、生态保护等各个方面。而目前绝大部分的气温产品，不管是用数学方法基于实测数据生成，还是基于数值预报模式生成，这些气温产品有一个共同点，空间分辨率基本都在0.1°及以上。低空间分辨率使得它们均在一定程度上丢失局部的地理特征，在很多以气温为输入参数的研究中，很难被直接应用。为了能将气温数据直接用到后续研究中就必须考虑提高现有气温数据的空间分辨率。

近地面气温是地表与大气进行水交换与能量交换的关键参数。得到地面上高空间分辨且较为准确的气温数据是保证相关研究精度的前提。目前得到高空间分辨率气温数据的方法主要分为三种，第一种是根据地面实测站点数据通过空间插值方法得到；第二种是基于卫星遥感数据反演得到；第三种是对气象模式输出的气温数据进行降尺度处理。

基于地面站点实测气温数据通过空间插值方法得到的气温数据不确定性大，对站点密度与站点分布要一定要求。而在实际情况中，站点布设需要耗费大量的人力物力，尤其是荒芜的山区与高原地区，通常情况是城区站点布设密集，山地站点稀少。因此通过插值方法得到的气温数据在山区不确定性更大，难以大规模的进行实际应用。

基于卫星遥感资料反演的地表温度，需要考虑的因素多，遥感资料本身存在误差，反演过程中存在误差累积，不确定性大。目前气温的反演方法也还在进一步完善中，且通过遥感数据反演得到的气温是卫星过境时刻气温，时间分辨率低。

而对气象模式输出的气温资料，运用统计降尺度方法得到的气温数据，虽然分辨率有所提高，在地市较为平坦的地方能达到精度要求，但是模式输出的数据同样具有较大的不确定性，在很多地方的适用性仍然很低。绝大部分都建立在统计降尺度基础上，更多的是考虑的从数学方法的角度去提升气温数据的分辨率，一些改进的统计降尺度方法，则需要依赖一些实测数据，很少考虑局部地形与植被影响的降尺度方法。蔡明勇、杨胜天等运用DEM与双线性内插的方法对全球陆地同化系统提供的3小时25公里空间分辨率地表大气温度数据进行空间降尺度得到1km近地面气温数据，但是只是单一的考虑DEM，而且对局地的气候与气温差异考虑不足，在海平面高度进行插值再反推到地面气温，丧失了气温产品本身的整体空间信息，没有考虑植被对温度调控作用。

目前有许多气温产品，时间分辨率达到三小时。例如中国科学院青藏高原研究所发布的中国区域地面气象要素数据集中提供的气温数据，具有较为可靠的适用性，但空间分辨率只有0.1°；欧洲中尺度天气预报中心提供的ERA-Interim再分析资料中的2米气温数据，经过生成模式的不断改进，在很多地方也具有较高的适用性，但空间分辨率最高也只有0.125°。而在地形复杂的山区以及高原地区，在网格大小为0.1°的范围内往往有地形的急剧起伏，以及植被覆盖的变化，这些地学因子都与气温有很大的相关性，尤其地表高程的变化，直接影响气温的升高或降低。

发明内容

本发明的目的是，针对上述问题，提出一种对现有气象模式生成的离地表2米高度处的粗分辨率气温数据进行降尺度的方法。

本发明的技术方案是：

针对现有气温数据的低空间分辨率，对局部地形与植被信息表达不足的特性，本方法引进DEM与NDVI对原始气温数据进行降尺度处理。根据一年中植被的生长情况，将一年划分为植被生长季和非植被生长季，分别对植被生长季和非植被生长季进行降尺度。

在植被生长季高程是气温变化的主导因素，根据相关研究，两者的平均相关系数为-0.6905，地表植被覆盖浓密，植被对近地面气温具有一定的调控作用。因此，在植被生长季，本文同时考虑DEM与NDVI对气温的影响。

对低分辨率的原始气温数据、DEM、NDVI，在整个研究区范围内进行筛选，找到NDVI存在的所有像元，将得到的这些低空间分辨率像元对应的地表气温与DEM、NDVI做多元回归，得到回归参数，再将高空间分辨率的DEM、NDVI数据带入回归方程，得到所有NDVI存在的高空间分辨率像元的降尺度气温值。对于NDVI不存在的像元，其气降尺度方法与非植被生长季一致。计算公式如公式(1)、(2)所示。

SAT_s＝A*DEM_s+B*NDVI_s+C>

SAT_h＝A*DEM_h+B*NDVI_h+C>

式中，SATs是低分辨率原始像元气温值；DEMs是低分辨率DEM像元值；NDVIs是低分辨率NDVI像元值；A为低分辨率时DEM的回归系数；B为低分辨率时NDVI的回归系数；C是回归方程常数项；SAT_h降尺度后高分辨率气像元气温值；DEM_h为与降尺度后气温分辨率一致的高分辨率DEM像元值；NDVI_h为与降尺度后气温分辨率一致的高分辨率NDVI像元值。

在非植被生长季，绝大部分植被将处于凋敝期，地表植被覆盖明显减少，裸露地表将增多，此时植被对气温的调节作用将明显较少,若此时再与植被生长季一样考虑植被对气温的调节作用，则在降尺度的过程中必将带来较大误差。此时在本文中将只考虑在同一时刻影响气温变化的主导因素即地表高程变化的影响。

在非植被生长季，降尺度分为两步，我们首先求得每一个降尺度后高分辨率像元与所在的与原始气温数据分辨率一致的低分辨率的大像元的高程差，再结合气温垂减率，得到每一个降尺度后高分辨率小像元与对应的原始低分辨率气温数据之间的气温补偿值，以每一个小像元对应大像元的原始气温值作为小像元的背景值，加上补偿值，得到降尺度后的每一个小像元的值。如公式3所示。

SAT_Dh＝SAT_Bh+λ*ΔH_h>

式中SAT_Dh是降尺度后高分辨率像元加上气温补偿值后；SAT_Bh是高分辨率像元的原始背景气温值；λ为当地气温直减率；△H_h是高分辨率像元与对应低分辨率大像元之间的高程差。

基于DEM降尺度时我们以每一个大像元值作为大像元范围内所有小像元的背景值，这样就会导致降尺度后在一些原始像元的分界处保持了与原始数据一致的差异性，气温数值出现跳跃，出现明显的分界。因此我们对基于DEM降尺度后的数据，考虑气温具有空间相关性和空间连续性，做平滑处理。我们用一个3*3的窗口的平均值去代替这个窗口最中间的像元值，一定程度上去减弱降尺度后保持的与原始像元一致的边界效应，得到最终的连续性更强的降尺度结果。如公式(4)所示：

式中，SAT_h(i,j)是平滑后第i列j行个小像元值；SAT_Dh(m,n)是加上气温补偿值后第m列n行个小像元值。

本发明的有益效果是：针对现有方法的不足，对已有低空间分辨率的气温数据集，考虑数字高程模型(DEM)变化，加入局部地理信息对气温影响，同时在植被生长季考虑植被覆盖对近地面气温的调控作用，将原始数据有效降尺度到高空间分辨率。

附图说明

图1本发明的技术流程图

图2青藏高原地区海拔高度

图3青藏高原地区2010年2月1日0时中国区域地面气象要素数据集原始气温数据

图4青藏高原地区2010年2月1日0时中国区域地面气象要素数据集气温数据降尺度后

图5局部地区2010年2月1日0时中国区域地面气象要素数据集气温数据降尺度前后对比

图6局部地区2010年2月1日0时ERA-Interim2米气温数据降尺度前后对比

图7青藏高原地区2010年6月1日0时中国区域地面气象要素数据集原始气温数据

图8青藏高原地区2010年6月1日0时中国区域地面气象要素数据集气温数据降尺度后

图9局部地区2010年6月1日0时中国区域地面气象要素数据集气温数据降尺度前后对比

图10局部地区2010年2月1日0时ERA-Interim2米气温数据降尺度前后对比

图11玛曲站6、7月11时、23时降尺度后站点温度，原始温度，降尺度温度对比折线图

具体实施方式

下面结合附图和实例详细描述本发明的技术方案。

实验区选择为中国境内广大的青藏高原第三极地区(E73°_E106°，N40°_N23°)。在这一地区包括了西南河流源区，高程变化急剧，如图(2)所示在整个青藏高原地区内海拔变化从60米到8000多米，地形的复杂多变对气温变化有一定程度影响，尤其是山顶，山谷，河谷等地方。在植被生长季植被覆盖各异，不同的植被覆盖对气温也有大小各异的影响。

需要用到的数据有中国高时空分辨率地面气象要素驱动数据集中的近地面瞬时气温数据，欧洲中期天气预报中心的ERA-interim的2米气温数据。SRTM(Shuttle RadarTopography Mission)90米分辨率DEM(Digital Elevation Model)数据，由美国宇航局(NASA)和国防部国家测绘局(NIMA)联合测量，数据来源于中国科学院计算机网络信息中心国际科学数据镜像网站(http://www.gscloud.cn)。MODIS的MYD11A2植被指数产品，该产品包含多种植被指数数据，本文主要用到其中的NDVI(Normalized Difference VegetationIndex)数据，空间分辨率为1KM，时间分辨率为16天合成，时间范围为每一年的植被生长季，大致时间为每年的第129天到257天之间的数据。

(1)数据进预处理

第一步：将90米空间分辨率DEM数据采用双线性内插法在ACGIS中重采样为0.01°空间分辨率，将得到0.01°空间分辨率的DEM重采样至0.1°空间分辨率，重采样方法为求得0.1°网格内100个0.01°像元的平均值，作为重采样0.1°像元的DEM值。

第二步：NDVI数据采用MRT批量从1KM分辨率重采样至0.01°分辨率并裁剪出与研究区相匹配的区域，因为1KM和0.01°在空间的大小很接近，所以重采样方法采用最邻近法。将第一步处理得到的0.01°分辨率NDVI进一步处理得到0.1°分辨率的NDVI，处理方法为，找到每一个0.1°大像元对应的100个0.01°像元，并判断其中存在缺失值的像元个数，如果在这个0.1°的大像元中缺失的小像元的个数小于30个，则用余下的数据存在的小像元的平均值作为最终的0.1°分辨率像元的NDVI值，否则对应的0.1°像元值为缺失值。

第三步：两种气温数据从CDF格式读取为tiff文件，并按照实验区域范围进行裁剪。最后对ERA-interim气温数据重采样至0.1°分辨率，与中国高时空分辨率地面气象要素驱动数据集中的近地面瞬时气温数据分辨率一致。

(2)植被生长季气温降尺度

对0.1°分辨率的原始气温数据、DEM、NDVI，在整个研究区范围内进行筛选，找到NDVI存在的所有像元，将得到的这些像元对应的地表气温与DEM、NDVI做多元回归，得到回归参数，再将0.01°分辨率的DEM、NDVI数据带入回归方程，得到所有NDVI存在的0.01°分辨率像元的降尺度气温值。如公式(5)、(6)所示。

SAT_0.1＝A*DEM_0.1+B*NDVI_0.1+C>

SAT_0.01＝A*DEM_0.01+B*NDVI_0.01+C>

式中SAT_0.1是0.1°分辨率原始像元气温值；DEM0.1是0.1°分辨率DEM像元值；NDVI_0.1是0.1°分辨率NDVI像元值；A是0.1°分辨率时DEM的回归系数；B是0.1°分辨率时NDVI的回归系数；C是回归方程常数项；SAT_0.01是0.01°分辨率降尺度后像元气温值；DEM_0.01是0.01°分辨率DEM像元值；NDVI_0.01是0.01°分辨率NDVI像元值。

(3)非植被生长季降尺度

在非植被生长季，首先基于DEM进行降尺度处理。我们求得每一个0.01°分辨率像元与所在的0.1°分辨率的大像元的高程差，再结合青藏高原地区的气温垂减率，得到每一个0.01°小像元与0.1°大像元的气温补偿值，以每一个小像元对应大像元的原始气温值作为小像元的背景值，加上补偿值，得到降尺度后的每一个小像元的值。如公式(7)所示。

SAT_D0.01＝SAT_B0.01+λ*DH_0.01>

式中，SATD_0.01是0.01°像元加上气温补偿值后；SATB_0.01是0.01°像元的原始背景气温值；λ是青藏高原地区气温直减率，根据江净超、翁笃鸣等对我国气温直减率的研究，得到我国青藏高寒区气温直减率，文中取值为0.5K/100m；△H_0.01是0.01°像元与对应0.1度大像元之间的高程差。

基于DEM降尺度时我们以每一个大像元值作为大像元范围内所有小像元的背景值，这样就会导致降尺度后在一些原始像元的分界处保持了与原始数据一致的差异性，气温数值出现跳跃，出现明显的分界。因此我们对基于DEM降尺度后的数据，考虑气温具有空间相关性和空间连续性，做平滑处理。我们用一个3*3的窗口的平均值去代替这个窗口最中间的像元值，一定程度上去减弱降尺度后保持的与原始像元一致的边界效应，得到最终的连续性更强的降尺度结果。平滑过程如公式(9)所示。

式中，SAT_0.01(i,j)为平滑后第i列j行个小像元值；SATD_0.01(m,n)为加上气温补偿值后第m列n行个小像元值。

结果

1.在非植被生长季，地表起伏很大的地方，一个0.1°像元的气温值不能很好的反应地表的真实气温，高程的急剧变化会导致气温的升高或降低。只有更的高空间分辨率才能较为准确的表达地表的地形起伏，在考虑地形对气温的影响之后，本方法能比较准确地描述这些区域地表的气温变化，在保留原始数据整体空间属性的基础上能得到高空间分辨率，连续性更强的气温数据。；

2.在植被生长季，除了高程是气温变化的主导因子以外，植被覆盖对气温有不可忽视的影响，引入DEM与NDVI对气温产品进行降尺度，不仅考虑地形因子的影响，也考虑了植被指数的影响，能更为准确的得到降尺度结果。本方法在植被生长季引入植被指数以后进行降尺度，结果较为理想。

3.本方法在整个一年时间内，不管是植被生长季还是非植被生长季降尺度度结果都能较为有效的提高原始数据的空间分辨率，相比原始气温数据能更好的表达局部气温的变，更加符合客观规律，能很好的达到提高原始气温数据空间分辨率的目的。