CN104143043A - 一种多功能气候数据模型及其应用 - Google Patents

Abstract

本发明公开了一种多功能气候数据模型及其应用，该多功能气候数据模型以气温海拔调整方程、双线性距离加权插值、谐波拟合等方法为基础，实现任意尺度下气候变量的生成。通过这些数据，接合物种分布模型，可以实现对未来气候条件下树种的适宜性分布进行预测，为林业部门造林或再造林时的树种选择提供依据，从而提高林地生产力。本发明所产生的气候数据，还可以为森林生态系统生长模型提供长期的气候数据支持，从而提高生长模型的预测精度；通过本发明的降尺度方法，可以对政府间气候变化专门委员会的众多大气环流模型预测结果进行任意尺度转换，从而为多模型多情景下的相关研究提供数据支持。

Description

一种多功能气候数据模型及其应用

技术领域

本发明属于气候数据模型技术领域，具体涉及一种多功能气候数据模型的构建方法及其应用。

背景技术

1906～2005年间，全球平均气温上升了0.56～0.92℃；根据不同的气候情景模拟估计未来100年中，全球平均温度将继续上升1.4℃～5.8℃，极端气候事件发生的频率和强度也将持续增加。如此大幅度的气候变化将在大区域范围内影响植被的类型和分布，增加森林灾害(火灾、病虫害等)的发生频率与强度，对我国林业发展产生诸多影响。

预测全球气候变化对森林生态系统的影响，评价森林生态系统对全球气候变化的响应及适应，是全球变化研究的重要领域之一，也是针对森林生态系统制定有效的适应气候变化措施的基础。准确可靠的气候基础数据是以上研究的前提，其获取途径通常为：气象站直接观测数据、IPCC大气环流模型的模拟数据以及全球相关研究机构的长期监测气候数据库等。然而这些数据通常空间分辨率较低、数据格式复杂且变量类型不统一。利用统计方法插值气象站观测数据可以有效地提高数据空间分辨率，其中以ANUSPLIN和PRISM算法应用最广泛。ANUSPLIN是基于薄盘平滑样条法内插多变量数据的算法；钱永兰等于2010年将其应用于获取中国的气温和降水变量数据集，取得了较好的效果。PRISM由美国俄勒冈州立大学模型组设计的基于知识的气候模型系统，其数据成果主要为北美800m-4000m每月气候栅格表面，朱华忠等于2003年通过收集中国的气象站观测数据应用PRISM生成了中国的逐月气温和降水数据。范泽孟等于2012年使用曲面建模(HASM)的方法，对全球气温和中国降水量对大气环流模型HadCM3的三个情景数据进行了降尺度模拟。

综上所述，如何有效的获得高精度气候数据已经成为一个重要的研究问题，而上述研究所提供的方法仍然存在一定的缺点，首先，由于无法动态降尺度，静态生成的历年每月高空间分辨率气候变量栅格数据将占用大量存储空间，其数据往往只能应用于一个固定的尺度条件下。其次，上述方法均使用气象站观测值进行插值，因此需要大量的高精度观测值对其生成的结果进行精度保证，然而，很多时候由于部分气象站观测数据的缺失，将会直接影响局部插值的效果和数据准确性。再次，上述方法所产生的气候变量相对较少，往往一些对林木生长有较大影响的气候变量无法获得，因此很难满足林业相关研究的实际需要。最后，以上方法所获得的气候数据往往只针对某一个时间范围或几个特定的时期，无法在时间尺度上提供全面的数据支持。

发明内容

发明目的：针对现有技术中存在的不足，本发明的目的是提供一种多功能气候数据模型，以气温海拔调整方程、双线性距离加权插值、谐波拟合等方法为基础，实现任意尺度下气候变量的生成。通过这些数据，接合物种分布模型，可以实现对未来气候条件下树种的适宜性分布进行预测，为林业部门造林或再造林时的树种选择提供依据，从而提高林地生产力。本发明的另一目的是提供一种上述多功能气候数据模型的应用。

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

一种多功能气候数据模型，其特征在于，由以下各步骤构建而成：

(1)通过逐步回归的方法，构建基线气温栅格数据与基线数字高程模型按月的36个多元线性回归模型，并分别对回归模型求海拔的一阶偏导数，形成气温变量的海拔调整方程；其中，多元线性回归模型的因变量包括：Lat²，Lat³，Lon²，Lon³，Elev²，Lat×Lon，(Lat×Lon)²，(Lat×Lon)³，(Lon/Lat)²，(Lon/Lat)³，Lat²×Lon，Lon²×Lat，Lon²/Lat，Lon/Lat²，Lat×Lon×Elev，Lat×Lon×Elev²，(Lon/Lat)×Elev，(Lon/Lat)×Elev²，Lat×Elev，Lat×Elev²，Lon×Elev，Lon×Elev²，Lat为纬度、Lon为经度、Elev为海拔；

海拔调整方程如下：

△T＝△Elev×(a+b₁×2×Elev+b₂×Lat×Lon+b₃×2×Lat×Lon×Elev+b₄×Lon/Lat+b₅×2×Lon/Lat×Elev+b₆×Lat+b₇×2×Lat×Elev+b₈×Lon+b₉×2×Lon×Elev)式(1)；

式中，a,b_i(i＝1到9)为参数，值见表1和2，△T为气温差，△Elev为高程差；

表1 1-6海拔调整方程参数值

表2 7-12月海拔调整方程参数值

(2)利用双线性距离加重插值算法将原始基线降水栅格表面转换为任意空间尺度的基线降水数据P_d；结合步骤(1)的气温海拔调整方程及双线性距离加权插值算法，将原始基线气温栅格表面生成任意空间尺度的基线气温数据T_be；

(3)不同时期基础气候数据的生成；数据包括1901-2012年每年、以及大气环流模型预测的未来2010-2039、2040-2069、2070-2099三个时期的累年未来的每月距平基础气候数据；对于任意一点P，在距平基础气候数据上搜索其四个临近像元，并采用双线性距离加权插值获得A_n，即该点的距平值，与步骤(2)中获得的对应位置T_be或P_d相加，从而得到最终的基础气候变量模拟值；

(4)结合步骤(3)中生成的月平均基础气温变量，通过式(5)中的谐波拟合方程产生每日气温变量，生成三个新的派生气温变量，日均气温、日最高气温及日最低气温；

式中，t为时间；T为月分号；D为每月天数；m为当月的第几天；Y_T相关气候变量的月平均值；

$a_j=[\left(\frac{\mathrm{\πj}}{12}\right)/\sin \left(\frac{\mathrm{\πj}}{12}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\frac{2\mathrm{\πjT}}{12}\right)/6],$

$a_6=[\left(\frac{\mathrm{\π}}{2}\right)/\sin \left(\frac{\mathrm{\π}}{2}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\mathrm{\πT}\right)/12],$

b₆＝0，

t＝(T-0.5)+(m-0.5)/D。

(5)无霜期天数、霜冻天数及光合有效辐射的计算：

a)结合步骤(3)中获得的月均最低温，采用经验公式6计算每月的无霜期天数，再通过公式7间接计算每月的霜冻天数；

${\mathrm{NFFD}}_m=\min (N_m,0.2762+1.0114\×(\frac{1}{1+1.15\×e^{(-0.4\×T_{\min })}}\×N_m))$    式(6)

NFD_m＝N_m-NFFD_m   式(7)

式中，NFFD_m为无霜期天数，NFD_m为霜冻天数，T_min为月均最低气温，N_m为当月的天数；

b)结合步骤(4)所生成的日最低气温、日最高气温及样地纬度，采用式(8)来计算日光合有效辐射，计算结果为MJ；

${\mathrm{PAR}}_d=0.7\×(1-e^{-(0.036\×e^{-0.154\×{\mathrm{\Δ}}_{\mathrm{mT}}}){\mathrm{\Δ}}_{\mathrm{dT}}^{2.4}})\×Q_{\mathrm{ex}}$    式(8)

其中，

Q_ex＝(0.0864/π)×1367×(1+0.033×cos(2π×d_Y/365))×hd×sin(Dec)×sin(Lat)+cos(Dec)×cos(Lat)×sin(hd)，

Dec＝Arcsin(0.4×Sin(2π×(d_Y-82)/365))，

hd＝Arccos(-tan(Lat)×tan(Dec)，

Lat＝π×siteLat/180；

式中，PAR_d为日光合有效辐射，Q_ex为忽略大气层影响情况下的太阳辐射，hd为日照时数，Dec为太阳倾斜角，d_Y为一年中的第几天，△_mT为30天平均最高气温与最低气温的温差，△_dT为当日的最高温与最低温的温差，siteLat为样地纬度；

(6)其它派生气候变量计算

根据步骤(3)生成基础气候数据，按年或季节进行统计，快速形成每年或不同季节的平均气温、平均最高气温、平均最低气温，总降水量；

通过步骤(4)所生成的日气温变量值，形成按月或年的各种积温变量；

通过统计月均气温大于5℃的月数生成温暖指数；

统计大于5℃、小于30℃的月均气温之和的年均值得到生物温度。

步骤2)中，对于任意的数据点P，首先在原始的基线气温和基线数字高程模型栅格表面上分别搜索其临近的四个像元；并采用双线性距离加权插值获得插值温度T_bil及插值海拔Elev_bil；然后将高程差△Elev＝Elev-Elev_bil代入气温海拔调整方程获得温差△T；最后将△T与T_bil相加获得T_be，即P点基线温度；其次，在原始基线降水栅格表面上搜索临近的四个像元，通过双线性距离加重插值得到P_d，即P点的基线降水量；双线性距离加权插值算法见公式2、3、4及图1。

P＝T₁×d₄+T₂×d₃   式(2)

T₁＝p₁×d₂+p₄×d₁   式(3)

T₂＝p₂×d₂+p₃×d₁   式(4)。

有益效果：与现有技术相比，本发明以气温海拔调整方程、双线性距离加权插值、谐波拟合等方法为基础，实现任意尺度下气候变量的生成。通过这些数据，接合物种分布模型，可以实现对未来气候条件下树种的适宜性分布进行预测，为林业部门造林或再造林时的树种选择提供依据，从而提高林地生产力；本发明所产生的气候数据，还可以为森林生态系统生长模型提供长期的气候数据支持，从而提高生长模型的预测精度；通过本发明的降尺度方法，可以对政府间气候变化专门委员会的众多大气环流模型预测结果进行任意尺度转换，从而为多模型多情景下的相关研究提供数据支持。

附图说明

图1是双线性距离加权原理图；

图2是日平均气温实测值与拟合值比较结果图；

图3是日最高气温实测值与拟合值比较结果图；

图4是日最低气温实测值与拟合值比较结果图；

图5是月太阳有效辐射实测值与估计值比较图；

图6是1980-2009时期气候条件下林分干枝生物量、根生物量、叶生物量及平均胸径随林龄的变化结果图；

图7是林分实测值与预测值比较结果图，图A为林分干枝生物量(t/ha)，图B为林分根生物量(t/ha)，图C为林分叶生物量(t/ha)，图D为林分平均胸径(cm)；

图8是在1980-2009及2040-2069时期气候条件下林分生长预测变量比较结果图；

图9是1980-2009时期和2040－2069时期温度适宜度变化图，图A为高精度数字高程模型，图B为2040－2069时期温度适宜度减去1980-2009时期温度适宜度的差值，图C为1980-2009时期温度适宜度，图D为2040-2069时期温度适宜度。

图10是1980-2009时期和2040－2069时期的年均气温分布图。

具体实施方式

下面结合具体实施例对本发明做进一步的说明。

实施例1

通过在人地系统主题数据库(http://www.data.ac.cn/index.asp)获得可下载公用的气象站观测数据，包括累年每月平均气温、累年每月平均最低气温、累年每月平均最高气温及累年每月降水量；形成气象站观测值。将各气象站的经纬度值代入双线性距离加权插值算法(如图1)，在1961-1990累年基线气候栅格表面和基线数字高程栅格表面上(基线气温栅格表面包括月均气温、月均最低气温、月均最高气温，空间分辨率为0.5°，来自于http://www.cru.uea.ac.uk/cru/data/hrg/的CRU TS3.21，基线降水栅格表面空间分辨率为1km，来自于人地系统主题数据库http://www.data.ac.cn/index.asp，基线数字高程栅格表面空间分辨率为0.5°，来自于http://srtm.csi.cgiar.org/)获得插值气候变量值(插值气温T_bil和插值降水P_d)和插值海拔(Elev_bil)，再利用各气象站的海拔(Elev)和插值海拔计算高程差(△Elev)，将各气象站的经纬度和高程差代入下式计算得到温差(△T)；

△T＝△Elev×(a+b₁×2×Elev+b₂×Lat×Lon+b₃×2×Lat×Lon×Elev+b₄×Lon/Lat+b₅×2×Lon/Lat×Elev+b₆×Lat+b₇×2×Lat×Elev+b₈×Lon+b₉×2×Lon×Elev)式(1)；

式中，a,b_i(i＝1到9)为参数，值见表1和2，△T为气温差，△Elev为高程差。

表1 1-6海拔调整方程参数值

表2 7-12月海拔调整方程参数值

最终各气象站的基线气温预测值(T_be)为T_bil与△T相加，基线降水量预测值为P_d。将气象站观测值与预测值进行1：1线拟合，其精度R²值见表3，数据分析表明，所生成的任意空间尺度基线数据能够很好的反应气象观测值，特别是经过海拔调整以后的基线气温数据，其拟合精度均达到0.8以上。

表3 任意空间尺度基线数据与气象站观测值拟合1:1线的R²值

实施例2

通过收集全国32个城市的累年日均气温观测数据，通过按月统计得到月均气温、月均最高气温及月均最低气温，将这三个气温变量的月平均值分别代入下式计算日气温数据拟合值，

式中，t为时间；T为月分号(1到12)；D为每月天数；m为当月的第几天；Y_T相关气候变量的月平均值；

$a_j=[\left(\frac{\mathrm{\πj}}{12}\right)/\sin \left(\frac{\mathrm{\πj}}{12}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\frac{2\mathrm{\πjT}}{12}\right)/6],$

$a_6=[\left(\frac{\mathrm{\π}}{2}\right)/\sin \left(\frac{\mathrm{\π}}{2}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\mathrm{\πT}\right)/12],$

b₆＝0，

t＝(T-0.5)+(m-0.5)/D。

定义公式中的月分号T＝1到12，1到12月每月天数(D)分别为31、28、31、30、31、30、31、31、30、31、30、31天，气候变量Y_T为月均气温、月均最高气温和月均最低气温；对日均气温数据观测值与日均气温拟合值通过拟合方程y＝ax+b进行了精度评价，结果见表4、图2-4。数据分析表明，三个日气温变量的拟合精度(R²)值均达到0.98以上，32个拟合方程的斜率(a)都在1左右，其均值达到0.994，拟合方程的截距(b)也都在1以下，表明方程拟合值与实际观测值的误差没有显著性差异。

表4 日观测气温与日拟合气温值验证结果

实施例3

使用一组来自于2006年澳洲东南部观测的日最高气温、日最低气温及日光合有效辐射数据，样地纬度-42.9°。首先按月统计平均得到月均最低气温和月均最高气温值；将这两个气温变量的每月平均值分别代入下式计算每日的日最低气温和日最高气温；

式中，t为时间；T为月分号；D为每月天数；m为当月的第几天；Y_T相关气候变量的月平均值；

$a_j=[\left(\frac{\mathrm{\πj}}{12}\right)/\sin \left(\frac{\mathrm{\πj}}{12}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\frac{2\mathrm{\πjT}}{12}\right)/6],$

$a_6=[\left(\frac{\mathrm{\π}}{2}\right)/\sin \left(\frac{\mathrm{\π}}{2}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\mathrm{\πT}\right)/12],$

b₆＝0，

t＝(T-0.5)+(m-0.5)/D。

定义公式中的月分号T＝1到12，1到12月每月天数(D)分别为31、28、31、30、31、30、31、31、30、31、30、31天，气候变量Y_T为月均最高气温和月均最低气温；再将计算获得的日最高气温、日最低气温结合样地纬度代入下式计算获得每日光合有效辐射估计值；

${\mathrm{PAR}}_d=0.7\×(1-e^{-(0.036\×e^{-0.154\×{\mathrm{\Δ}}_{\mathrm{mT}}}){\mathrm{\Δ}}_{\mathrm{dT}}^{2.4}})\×Q_{\mathrm{ex}}$    式(8)

其中，

Q_ex＝(0.0864/π)×1367×(1+0.033×cos(2π×d_Y/365))×hd×sin(Dec)×sin(Lat)+cos(Dec)×cos(Lat)×sin(hd)，

Dec＝Arcsin(0.4×Sin(2π×(d_Y-82)/365))，

hd＝Arccos(-tan(Lat)×tan(Dec)，

Lat＝π×siteLat/180；

式中，PAR_d为日光合有效辐射，Q_ex为忽略大气层影响情况下的太阳辐射，hd为日照时数，Dec为太阳倾斜角，d_Y为一年中的第几天，△_mT为30天平均最高气温与最低气温的温差，△_dT为当日的最高温与最低温的温差，siteLat为样地纬度；

最后对观测和估计的日光合有效辐射值按月统计；最终形成每月日均光合有效辐射观测值和估计值，并进行拟合分析(图5)，拟合精度(R²)达到0.992，斜率接近于1，估计值与观测值间并没有显著性差异。

实施例4样地尺度模拟气候变量在林分生长模型中的应用

使用赵梅芳在湖南会同获得的杉木林生长数据(Zhao M,Xiang W,Peng C,TianD.Simulating age-related changes in carbon storage and allocation in a Chinese fir plantation growing in southern China using the3-PG model.Forest Ecology and Management2009；257(6):1520-1531)，地理位置为北纬26°50'，东经109°45'，平均海拔858m。

首先，利用样地经纬度信息通过双线性距离加权插值算法(图1)，在1961-1990累年基线气候栅格表面和基线数字高程栅格表面上获得样地的插值气候变量值(插值气温T_bil和插值降水P_d)和插值海拔(Elev_bil)，其中，基线气温栅格表面包括月均气温、月均最低气温、月均最高气温，其空间分辨率为0.5°，来自于http://www.cru.uea.ac.uk/cru/data/hrg/的CRU TS3.21，基线降水栅格表面空间分辨率为1km，来自于人地系统主题数据库http://www.data.ac.cn/index.asp，基线数字高程栅格表面空间分辨率为0.5°，来自于http://srtm.csi.cgiar.org/；再利用样地海拔(Elev)和插值海拔计算高程差(△Elev)，将样地的经纬度和高程差代入下式计算得到温差(△T)，样地的基线气温变量值(T_be)为T_bil与△T相加，基线降水量值为P_d。

△T＝△Elev×(a+b₁×2×Elev+b₂×Lat×Lon+b₃×2×Lat×Lon×Elev+b₄×Lon/Lat+b₅×2×Lon/Lat×Elev+b₆×Lat+b₇×2×Lat×Elev+b₈×Lon+b₉×2×Lon×Elev)式(1)；

式中，a,b_i(i＝1到9)为参数，值见表1和2，△T为气温差，△Elev为高程差；

其次，利用样地经纬度在1980-2009每年和2040-2069累年的距平气候栅格表面(1980-2009每年距平数据通过来自于http://www.cru.uea.ac.uk/cru/data/hrg/的绝对气候栅格表面CRU TS3.21与基线气候栅格表面相减获得，2040-2069累年距平数据来自于http://www.ipcc-data.org/，模型NIES:MIROC3_2-HI的SRB1情景)进行双线性距离加权插值，获得样地的距平值(A_n)，最终1980-2009每年和2040-2069累年每月的月均气温、月均最高气温和月均最低气温由T_be和A_n相加形成，月降温水量由P_d和A_n相加形成。

然后，将月均最高气温和月均最低气温分别代入下式计算每日的日最低气温和日最高气温，定义公式(5)中的月分号T＝1到12，1到12月每月天数(D)分别为31、28、31、30、31、30、31、31、30、31、30、31天，气候变量Y_T为月均最高气温和月均最低气温；

式中，t为时间；T为月分号；D为每月天数；m为当月的第几天；Y_T相关气候变量的月平均值；

$a_j=[\left(\frac{\mathrm{\πj}}{12}\right)/\sin \left(\frac{\mathrm{\πj}}{12}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\frac{2\mathrm{\πjT}}{12}\right)/6],$

$a_6=[\left(\frac{\mathrm{\π}}{2}\right)/\sin \left(\frac{\mathrm{\π}}{2}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\mathrm{\πT}\right)/12],$

b₆＝0，

t＝(T-0.5)+(m-0.5)/D。

并将计算获得的日最高气温、日最低气温结合样地纬度代入下面公式计算获得每日光合有效辐射估计值，对日光合有效辐射值按月平均，形成1980-2009每年和2040-2069累年每月的日均光合有效辐射值。

${\mathrm{PAR}}_d=0.7\×(1-e^{-(0.036\×e^{-0.154\×{\mathrm{\Δ}}_{\mathrm{mT}}}){\mathrm{\Δ}}_{\mathrm{dT}}^{2.4}})\×Q_{\mathrm{ex}}$    式(8)

其中，

Q_ex＝(0.0864/π)×1367×(1+0.033×cos(2π×d_Y/365))×hd×sin(Dec)×sin(Lat)+cos(Dec)×cos(Lat)×sin(hd)，

Dec＝Arcsin(0.4×Sin(2π×(d_Y-82)/365))，

hd＝Arccos(-tan(Lat)×tan(Dec)，

Lat＝π×siteLat/180；

式中，PAR_d为日光合有效辐射，Q_ex为忽略大气层影响情况下的太阳辐射，hd为日照时数，Dec为太阳倾斜角，d_Y为一年中的第几天，△_mT为30天平均最高气温与最低气温的温差，△_dT为当日的最高温与最低温的温差，siteLat为样地纬度；

最后，将月均最低气温代入公式(6)计算1980-2009每年和2040-2069累年每月无霜期天数(NFFD_m)，再将其代入公式(7)中计算得到1980-2009每年和2040-2069累年每月霜冻天数(NFD_m)，其中1到12月每月天数(N_m)分别为31、28、31、30、31、30、31、31、30、31、30、31天。

${\mathrm{NFFD}}_m=\min (N_m,0.2762+1.0114\×(\frac{1}{1+1.15\×e^{(-0.4\×T_{\min })}}\×N_m))$    式(6)

NFD_m＝N_m-NFFD_m   式(7)

式中，NFFD_m为无霜期天数，NFD_m为霜冻天数，T_min为月均最低气温，N_m为当月的天数；

将1980-2009每年的值按年平均，最终生成的1980-2009累年和2040-2069累年每月气候数据结果见表5和表6。

表5 1980-2009累年平均气候数据

月分	月均最高气温	月均最低气温	月均气温	月降水量	光合有效辐射	霜冻天数
							1月	8.69	2.37	5.42	32.15	10.55	9.03
2月	10.46	4.38	7.36	50.77	12.50	4.12
							3月	14.71	7.73	11.16	73.57	16.55	0.93
4月	21.26	13.27	17.21	100.00	20.39	0
							5月	25.92	17.55	21.71	157.20	22.45	0
6月	29.00	21.19	25.04	207.84	22.04	0
							7月	31.34	23.03	27.14	232.52	22.54	0
8月	31.33	22.47	26.89	166.01	21.90	0
							9月	28.29	18.87	23.56	121.04	20.28	0

10月	22.58	14.14	18.23	78.07	16.38	0
							11月	17.06	8.68	12.83	38.79	13.63	0.43
12月	11.69	3.82	7.71	20.93	11.91	5.63

根据梅芳在湖南会同的杉木林生长参数，将表5和表6的气候数据分别带入林分生长模型3-PG(花利忠.雷州半岛桉树人工林3PG模型的研究[硕士论文].福建:福建农林大学；2004；赵梅芳.基于3-PG机理模型的杉木林碳固定及蒸散量模拟研究[硕士论文].长沙:中南林业科技大学；2008.)中进行模拟。模型的输出林分生长变量包括：林分干枝生物量、林分根生物量、林分叶生物量和林分平均胸径。

表6 2040-2069累年平均气候数据

月分	月均最高气温	月均最低气温	月均气温	月降水量	光合有效辐射	霜冻天数
							1月	11.49	3.74	7.27	28.61	12.43	5.79
2月	11.78	4.92	7.96	41.62	13.61	3.33
							3月	17.07	9.23	12.60	66.87	17.71	0.24
4月	23.07	14.60	18.62	124.29	20.96	0
							5月	28.05	19.50	23.04	177.95	22.55	0
6月	31.43	22.61	26.35	201.94	23.39	0
							7月	34.18	24.91	29.05	186.34	23.53	0
8月	34.15	24.46	28.62	176.08	22.66	0
							9月	30.35	20.79	25.28	136.35	20.28	0
10月	24.94	15.88	19.73	94.99	16.92	0
							11月	18.77	10.52	14.26	53.15	13.50	0
12月	13.76	5.60	9.26	25.28	12.09	2.79

杉木林1980-2009时期30年林分生长模拟结果如图6所示，将模拟的4个林分生长变量与观测的杉木林生长数据进行拟合精度评价，如图7所示。结果数据分析表明，4个林分生长变量的预测值与实测值的拟合精度(R²)均在0.97以上，完全能够达到3-PG生长模型应用的要求。杉木林2040-2069时期的林分生长模拟结果与1980-2009时期的对比如图8所示，可见未来气候条件下，将会对杉木林在这一样地的生长产生一定的不利影响；2040-2069时期的林分生长变量在林分生长初期与1980-2009时期没有明显差异，但随着林龄的增加，这种差异将变得逐渐明显。

实施例5区域尺度模拟气温变量在树种温度适宜性中的应用

通过在SRTM网站(http://srtm.csi.cgiar.org/)上下载得到的中国东南部地区的高精度数字高程模型(图9A)，在其每一个像元点上获得纬度、经度及海拔信息。

首先，利用每个像元点的经纬度信息采用双线性距离加权插值算法(图1)，在1961-1990累年基线月均气候栅格表面和基线数字高程栅格表面上获得每个像元的插值月均气温变量值(插值气温T_bil)和插值海拔(Elev_bil)，其中，基线月均气候栅格表面的空间分辨率为0.5°，http://www.cru.uea.ac.uk/cru/data/hrg/的CRU TS3.21，基线数字高程栅格表面空间分辨率为0.5°，来自于http://srtm.csi.cgiar.org/；再利用每个像元的海拔(Elev)和插值海拔计算高程差(△Elev)，将每个像元的经纬度和高程差代入下式计算得到温差(△T)；每个像元的基线月均气温变量值(T_be)为T_bil与△T相加。

△T＝△Elev×(a+b₁×2×Elev+b₂×Lat×Lon+b₃×2×Lat×Lon×Elev+b₄×Lon/Lat+b₅×2×Lon/Lat×Elev+b₆×Lat+b₇×2×Lat×Elev+b₈×Lon+b₉×2×Lon×Elev)式(1)；

式中，a,b_i(i＝1到9)为参数，值见表1和2，△T为气温差，△Elev为高程差；

其次，利用每个像元的经纬度在1980-2009每年和2040-2069累年的距平月均气温栅格表面(1980-2009每年距平数据通过来自于http://www.cru.uea.ac.uk/cru/data/hrg/的绝对气候栅格表面CRU TS3.21与基线气候栅格表面相减获得，2040-2069累年距平数据来自于http://www.ipcc-data.org/，模型NIES:MIROC3_2-HI的SRB1情景)进行双线性距离加权插值，获得每个像元的距平值(A_n)，最终1980-2009每年和2040-2069累年每月的月均气温由T_be和A_n相加形成。

最后，对已经生成的1980-2009每年月均气温，对每个像元12个月的值进行平均，获得1980-2009每年的年平均气温，再对每个像元30年的值进行平均，获得1980-2009累年的年平均气温；对已经生成的2040-2069累年每月的月均气温，对每个像元12个月的值进行平均，形成2040-2069累年的年平均气温；两个时期的年均气温分布如图10。

根据杉木的生长三基点温度：生长最高温(T_max)、生长最低温(T_min)及生长最适温(T_opt)，使用温度适宜性模型(公式9)进行杉木温度适宜性研究。

同实施例4数据来源一样，根据赵梅芳在湖南会同获得的杉木林生长三基点温度，T_max＝40℃，T_min＝0℃，T_opt＝17.5℃；将上述所生成的两个时期的区域年均气温栅格数据的每一个像元点气温(T)带入公式(9)，从而形成杉木林1980-2009时期(图9C)和2040-2069时期(图9D)的温度适宜度分布图；两期温度适宜度图的差值如图9B所示，可见，未来气温在这一区域的上升，会使得这一区域温度适宜度在北部有增加的均势，而在南部有降低的趋势，即杉木的生长温度适宜区有一定的向北偏移现象。

从以上结果可以看了，本发明不仅可以为样地尺度气候变化研究提供基础数据，还可以为区域尺度林业相关研究提供数据支持，是其他方法难以比拟的。

Claims (2)

1.一种多功能气候数据模型，其特征在于，由以下各步骤构建而成：

(1)通过逐步回归的方法，构建基线气温栅格数据与基线数字高程模型按月的36个多元线性回归模型，并分别对回归模型求海拔的一阶偏导数，形成气温变量的海拔调整方程；其中，多元线性回归模型的因变量包括：Lat²，Lat³，Lon²，Lon³，Elev²，Lat×Lon，(Lat×Lon)²，(Lat×Lon)³，(Lon/Lat)²，(Lon/Lat)³，Lat²×Lon，Lon²×Lat，Lon²/Lat，Lon/Lat²，Lat×Lon×Elev，Lat×Lon×Elev²，(Lon/Lat)×Elev，(Lon/Lat)×Elev²，Lat×Elev，Lat×Elev²，Lon×Elev，Lon×Elev²，Lat为纬度、Lon为经度、Elev为海拔；

海拔调整方程如下：

△T＝△Elev×(a+b₁×2×Elev+b₂×Lat×Lon+b₃×2×Lat×Lon×Elev+b₄×Lon/Lat+b₅×2×Lon/Lat×Elev+b₆×Lat+b₇×2×Lat×Elev+b₈×Lon+b₉×2×Lon×Elev)式(1)；

式中，a,b_i(i＝1到9)为参数，值见表1和2，△T为气温差，△Elev为高程差；

表1 1-6海拔调整方程参数值

表2 7-12月海拔调整方程参数值

(2)利用双线性距离加重插值算法将原始基线降水栅格表面转换为任意空间尺度的基线降水数据P_d；结合步骤(1)的气温海拔调整方程及双线性距离加权插值算法，将原始基线气温栅格表面生成任意空间尺度的基线气温数据T_be；

(3)不同时期基础气候数据的生成；数据包括1901-2012年每年、以及大气环流模型预测的未来2010-2039、2040-2069、2070-2099三个时期的累年未来的每月距平基础气候数据；对于任意一点P，在距平基础气候数据上搜索其四个临近像元，并采用双线性距离加权插值获得A_n，即该点的距平值，与步骤(2)中获得的对应位置T_be或P_d相加，从而得到最终的基础气候变量模拟值；

(4)结合步骤(3)中生成的月平均基础气温变量，通过式(5)中的谐波拟合方程产生每日气温变量，生成三个新的派生气温变量，日均气温、日最高气温及日最低气温；

式中，t为时间；T为月分号；D为每月天数；m为当月的第几天；Y_T相关气候变量的月平均值；

$a_j=[\left(\frac{\mathrm{\πj}}{12}\right)/\sin \left(\frac{\mathrm{\πj}}{12}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\frac{2\mathrm{\πjT}}{12}\right)/6],$

$a_6=[\left(\frac{\mathrm{\π}}{2}\right)/\sin \left(\frac{\mathrm{\π}}{2}\right)]\underset{T}{\mathrm{\Σ}}[Y_T\cos \left(\mathrm{\πT}\right)/12],$

b₆＝0，

t＝(T-0.5)+(m-0.5)/D。

(5)无霜期天数、霜冻天数及光合有效辐射的计算：

a)结合步骤(3)中获得的月均最低温，采用经验公式6计算每月的无霜期天数，再通过公式7间接计算每月的霜冻天数；

${\mathrm{NFFD}}_m=\min (N_m,0.2762+1.0114\×(\frac{1}{1+1.15\×e^{(-0.4\×T_{\min })}}\×N_m))$    式(6)

NFD_m＝N_m-NFFD_m   式(7)

式中，T_min为月均最低气温，N_m为当月的天数；

b)结合步骤(4)所生成的日最低气温、日最高气温及样地纬度，采用式(8)来计算日光合有效辐射，计算结果为MJ；

${\mathrm{PAR}}_d=0.7\×(1-e^{-(0.036\×e^{-0.154\×{\mathrm{\Δ}}_{\mathrm{mT}}}){\mathrm{\Δ}}_{\mathrm{dT}}^{2.4}})\×Q_{\mathrm{ex}}$    式(8)

其中，

Q_ex＝(0.0864/π)×1367×(1+0.033×cos(2π×d_Y/365))×hd×sin(Dec)×sin(Lat)+cos(Dec)×cos(Lat)×sin(hd)，

Dec＝Arcsin(0.4×Sin(2π×(d_Y-82)/365))，

hd＝Arccos(-tan(Lat)×tan(Dec)，

Lat＝π×siteLat/180；

式中，PAR_d为日光合有效辐射，Q_ex为忽略大气层影响情况下的太阳辐射，hd为日照时数，Dec为太阳倾斜角，d_Y为一年中的第几天，△_mT为30天平均最高气温与最低气温的温差，△_dT为当日的最高温与最低温的温差，siteLat为样地纬度；

(6)其它派生气候变量计算

根据步骤(3)生成基础气候数据，按年或季节进行统计，快速形成每年或不同季节的平均气温、平均最高气温、平均最低气温，总降水量；

通过步骤(4)所生成的日气温变量值，形成按月或年的各种积温变量；

通过统计月均气温大于5℃的月数生成温暖指数；

统计大于5℃、小于30℃的月均气温之和的年均值得到生物温度。

2.根据权利要求1所述的多功能气候数据模型，其特征在于，步骤2)中，对于任意的数据点P，首先在原始的基线气温和基线数字高程模型栅格表面上分别搜索其临近的四个像元；并采用双线性距离加权插值获得插值温度T_bil及插值海拔Elev_bil；然后将高程差△Elev＝Elev-Elev_bil代入气温海拔调整方程获得温差△T；最后将△T与T_bil相加获得T_be，即P点基线温度；其次，在原始基线降水栅格表面上搜索临近的四个像元，通过双线性距离加重插值得到P_d，即P点的基线降水量；双线性距离加权插值算法见公式2、3、4及图1；

P＝T₁×d₄+T₂×d₃   式(2)

T₁＝p₁×d₂+p₄×d₁   式(3)

T₂＝p₂×d₂+p₃×d₁   式(4)。