积雪是气候系统重要的组成部分,是冰冻圈中地理分布最广泛、季节和年际变化最显著的成员。雪盖通过表面不同的反照率、地表通量,调制了陆面-大气间的能量交换,极大地影响了大气非绝热加热过程( 孙菽芬和李敬阳,2002 )。 肖林和车涛(2015) 在进行青藏高原积雪对气候反馈的初步研究时发现积雪对气候存在正的辐射强迫,冬季积雪的异常分布会导致东亚冬季风环流的异常( 陈海山等,2003 )。同时,春季积雪异常改变了夏季的土壤湿度和温度分布( 李震坤等,2009 )。 朱玉祥和丁一汇(2007) 则认为高原积雪是一种重要的陆面强迫因子,青藏高原冬季多雪,是引起中国东部夏季降水出现“南涝北旱”的一个重要原因。近年来全球气候异常,可能与前期的积雪异常有很大的关系,应关注积雪影响天气气候的问题( Li and Wang, 2011 ; 左志燕和张人禾,2012 ; 许立言和武炳义,2012a ; 2012b ; 张若楠等,2015 ; 张人禾等,2016 ; 杨凯等,2017 )。
对积雪物理过程的认识随着观测资料的不断丰富而逐渐深入,不同精度的反映积雪内部水、热物理过程的模型不断发展( Jordan,1991 ; Yang et al,1997 ; Sun et al,1999 ; Dai et al,2003 )。用于气候研究的绝大多数GCMs所采用的积雪参数化方案( 李伟平等,2009 )或模型则过于简单( Wu and Wu, 2004 ),这些简单的模型对积雪的变化难以满足当代气候研究的需要( 陈海山和孙照渤,2004 ),因此许多学者开展了积雪模拟研究, Wei and Dong(2015) 进行了CMIP5青藏高原积雪深度模拟, Yuan et al(2016) 在“SnowWind”模型基础上提高了对青藏高原积雪的估计, Xie and Hu(2016) , Xie et al(2018) 对比了两种积雪覆盖率的参数化方案,以及对积雪分布以及收支平衡的影响,认为积雪变率的表达对青藏高原积雪模拟效果的提高十分重要。 李弘毅和王建(2013) 则针对我国青藏高原积雪特征,提炼出现阶段空间分布式积雪水文模拟中的关键问题;模式结构对模拟结果至关重要,同时, Zhang(2015) 也发现影响陆面过程模型积雪模拟精度的不确定性中,驱动的不确定性最大,尤其是驱动中的降水和辐射,其次是模式结构以及参数化方案。因此,对于陆面过程模拟,驱动场和模式结构是决定模式模拟效果至关重要的因素,不同的驱动场会对模拟结果产生不同的影响。
国家气象信息中心提供的CLDAS(CMA Land Data Assimilation System)驱动数据,尽管在夏季的陆面过程模拟有着良好的表现( 潘旸等,2012 ),但是由于原降水驱动融合的背景场和站点观测在冬季存在一些问题,在冬季降水的刻画能力上有所欠缺,因此CLDAS研究团队对此进行了改进。本文分别选用了改进前后的两种CLDAS降水驱动陆面模式进行模拟,对积雪覆盖率(snow cover fraction,SCF)、雪深(snow depth,SD)、雪水当量(snow water equivalent,SWE)进行评估。一方面对比不同驱动对模式模拟结果所带来的影响,一方面分析改进的降水驱动的积雪模拟结果在中国区域的改进效果和具体表现。
1 研究区域和数据 1.1 研究区域中国积雪的分布自南向北逐渐增厚,由西向东,明显减少。中国积雪分布范围广泛,其中新疆、东北和青藏高原是三个大值区,主要分布在东北和内蒙古自治州东部地区、新疆北部和西部地区以及青藏高原地区。依据 李培基和米德生(1983) 应用中国气象台站积雪观测资料,以年累计积雪日数大于60 d的地区为稳定积雪区的定义,中国的东北区域、新疆地区及青藏高原地区均为稳定积雪区。本文选取东北区域(35°~61°N、115°~141°E)、新疆区域(34°~48°N、73°~96°E),以及青藏高原区域(26°~41°N、73°~105°E)这三个稳定积雪区和中国区域(15°~55°N、70°~140°E)对积雪模拟的效果进行评估。
1.2 CLDAS驱动数据 1.2.1 CLDAS驱动数据“CLDAS驱动数据”产品包含气温、气压、湿度、风速、降水、短波辐射六个要素,各要素为覆盖亚洲区域(0°~65°N、60°~160°E)的逐小时、等经纬度(0.0625°×0.0625°)网格产品,时间长度从2008年至今。气温、气压、湿度、风速产品以ECMWF数值分析/预报产品为背景场,通过多重网格三维变分技术融合地面自动站观测数据而形成;短波辐射产品基于DISORT辐射传输模型, 利用FY-2卫星全圆盘标称图VIS通道数据反演而形成;降水由FY-2/EMSIP降水产品与中国区域降水融合产品叠加拼接而成。其中中国区域的融合降水是由CMORPH(Climate Prediction Center Morphing)卫星降水以及自动站观测降水融合而成。但是尽管该降水( 潘旸等,2012 )在夏季取得了较好的表现, 但是它仍存在些许不足,一方面作为背景场的CMORPH卫星监测冬季降水的能力差,导致冬季部分地区降水数据缺失,另一方面,自动站在冬季降水观测不准确,导致该降水驱动在冬季与观测的一致性较差,进一步可能导致陆面模拟中积雪变量的低估。本文选取驱动的时间范围为2008—2014年。
1.2.2 CLDAS-Prcp驱动数据基于驱动数据降水存在的问题,CLDAS研究团队对此进行了改进,制作了CLDAS-Prcp驱动数据。使用MERRA2(Modern Era Retrospective-analysis for Res)再分析降水资料对CMORPH卫星降水资料进行了补充,将两者作为联合背景场,然后结合MERRA2再分析资料、日降水人工观测资料、气温资料以及自动站逐小时降水观测,得到新的逐小时降水资料,弥补冬季结冰期区域逐小时站点降水观测资料,制作了1998—2015年时间分辨率为逐小时、空间分辨率为0.0625°的中国区域融合降水产品。本文选用了采用CLDAS-Prcp驱动陆面模式进行积雪变量的模拟,时间范围为2008—2014年。
1.3 评估数据 1.3.1 MODIS积雪覆盖率数据MODIS是美国EOS(Earth Observing System)系统中主要的传感器之一,它具有较高的空间分辨率、光谱分辨率和时间分辨率,在积雪、植被研究等方面都得到了广泛应用( Yuan et al,2016 ; 李洛晞等,2016 )。MODIS/AQUA(TERRA)包括MYD10C1(MOD10C1)、MYD10C2(MOD10C2)和MYD-10CM(MOD10CM)积雪覆盖率产品,MYD10C1(MOD10C1)是逐日雪盖产品,是将分辨率为500 m的MYD10A1观测依据地表类型处理成0.05℃MG的格网数据产品。MODIS积雪产品的精度较高,在不受积雪条件和地表覆盖类型的影响,准确率大约可达90%( Hall et al,2002 ; 2010 ; Hall and Riggs, 2007 )。MOD10A2的绝对精度可以达到93%( Hall and Riggs, 2007 ),MOD10A1在加拿大地区雪季160 d的平均一致率可达到93%( Simic et al,2003 )。同时,MODIS积雪产品在中国地区的质量也得到广泛验证,新疆北部晴天时MOD10A1产品精度很高,总精度可达到98.5%,积雪分类精度为98.2%( 张学通等,2008 )。东北地区MOD10A2积雪产品则倾向于低估积雪覆盖面积,MOD10C2的积雪产品监测的积雪面积与较为相近( 仲桂新等,2010 )。因此本文采用MODIS/AQUA的数据集MYD10C1作为检验数据,数据采用等经纬度投影,空间分辨率为0.05°×0.05°,本评估选取的数据时间为2008—2014年。
1.3.2 站点雪深数据雪深的检验数据是由中国气象局国家气象信息中心提供的气象数据统一服务接口(Meteorological Unified Service Interface Community,MUSIC)的全国综合气象信息共享系统(CIMISS)数据。MUSIC是面向气象业务和科研,提供全国统一、标准、丰富的数据访问服务和应用编程接口。本文从CIMISS获取中国地面逐日积雪深度资料,该接口提供1951年1月1日00:00:00时至今的积雪深度资料(单位: cm),评估选用的时间是2008—2014年。
1.3.3 MERRA2雪水当量数据MERRA是新一代全球高分辨率卫星时代再分析资料,产品相比其他全球再分析资料在青藏高原地面气温表达方面具有一定的优势,在能量平衡与水循环等诸多领域得到了广泛应用( 除多等,2016 );在中亚地区的适用性评估结果表明,MERRA的降水资料比ERA-Interim具有更高的精度( 胡增运等,2013 ; 秦艳慧等,2015 )。MERRA的资料在中纬度密集植被区显示出比AMSRE更好的优势( Yi et al,2011 );而MERRA2数据是NASA发布的大气再分析数据,融合了最新的GEOS-5(Goddard Earth Observing System Data Assimilation System Version 5)数据。本文采用MERRA2中Land数据作为雪水当量的检验数据,Land数据包括一系列陆表变量,例如积雪覆盖率、积雪深度、雪水当量、叶面积指数以及地面总降水,空间分辨率为0.5°×0.625°,评估采用雪水当量(单位:kg·m -2 )作为模式模拟的检验数据,时间分辨率为1 h。
2 试验设计和方法 2.1 模式介绍Noah陆面过程模式基于在俄勒冈州立大学OSU模式,由美国环境预报中心多个机构共同发展、优化以及改进,可以同WRF模式大气以及其他模式耦合。Noah刻画了积雪的积累和消融过程,基于物理参数化方案模拟陆气之间能量、动量、水的交换,模拟积雪的积累、升华、融化和热交换处理,在积雪处理方面考虑了积雪老化、冰盖液态水再结冰、冻土以及积雪压实过程的动态雪深和密度变化过程,模式采用积雪-植被-土壤混合层。模式可以输出SWE、SD以及积雪密度,并且通过区域的经验消融曲线计算出积雪覆盖率(SCF)。Noah3.3-3.6采用的关于SCF参数化方案是( Koren et al,1999 )针对NECP气候模型雪盖和冻土参数化方案提出的,这个参数化方案强调了积雪覆盖率和雪水当量的关系。很多学者也对该模型进行了进一步的研究与评估( Ben,2010 ; 田静等,2011 ; 张果等,2016 )。本文陆面模拟采用Noah3.6。Noah中积雪融化受向下长波辐射、向下短波辐射以及感热通量和潜热通量影响,在这几个变量的影响下,如果地表温度低于0℃,则积雪不融化,否则积雪蒸发、升华以及融化可能会同时发生。当处于降雪或冰冻状态下,模式利用降水对雪水当量进行初始化,雪深受雪水当量、积雪密度、积雪表面温度以及土壤表面温度影响,由雪水当量和积雪密度计算得出。当雪水当量超过积雪为全覆盖时雪水当量阈值时,积雪覆盖率为100%,否则依据当时雪水当量、雪水当量阈值,调整系数重新计算积雪覆盖率,计算关系如下:
\begin{array}{l} RSNOW = SNEQV/SNUP{\rm{ }}\\ SNCOVR = {\rm{ }}1 - [EXP\left({ - SALP \times RSNOW} \right) - {\rm{ }}\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;RSNOW \times EXP\left({ - SALP} \right)] \end{array} 2.2 试验设计和评估指标试验采用Noah3.6陆面模式,利用CLDAS和CLDAS-Prcp分别进行模拟,对模式中积雪变量,包括积雪覆盖率、雪深、雪水当量,在东北区域、新疆区域、青藏高原区域三大稳定积雪区,以及中国区域2008—2014年积雪模拟的表现进行评估。模式输出模拟结果的频率为1 d,模式验证数据分别为MODIS/AQUA积雪覆盖率日数据MYD10C1,国家级自动站CIMISS雪深日数据以及MERRA2雪水当量小时数据。在利用MERRA2雪水当量小时数据进行比对时,首先要进行时间尺度变换,将小时数据处理成日平均雪水当量。为了进行对比评估,检验过程中将格点的数据统一处理成模式输出的分辨率0.0625°。本文采取的评估指标包括季节平均值、偏差(BIAS)、均方根误差(RMSE)以及相关系数(CORR)具体公式如下:
\begin{array}{l} BIAS = {N^{ - 1}}\sum\limits_{i = 1}^N {[{X_{\rm{m}}}\left(i \right) - {X_{\rm{o}}}\left(i \right)]} \\ RMSE = \sqrt {{N^{ - 1}}\sum\limits_{i = 1}^N {{{[{X_{\rm{m}}}\left(i \right) - {X_{\rm{o}}}\left(i \right)]}^2}{\rm{ }}} } \\ CORR = \frac{{\sum\limits_{i = 1}^N {\left[ {{X_{\rm{m}}}\left(i \right) - {N^{ - 1}}\sum\limits_{i = 1}^N {{X_{\rm{m}}}\left(i \right)} } \right]\left[ {{X_{\rm{o}}}\left(i \right) - {N^{ - 1}}\sum\limits_{i = 1}^N {{X_{\rm{o}}}\left(i \right)} } \right]} }}{{\sqrt {\sum\limits_{i = 1}^N {{{\left[ {{X_{\rm{m}}}\left(i \right) - {N^{ - 1}}\sum\limits_{i = 1}^N {{X_{\rm{m}}}\left(i \right)} } \right]}^2}} } \sqrt {\sum\limits_{i = 1}^N {{{\left[ {{X_{\rm{o}}}\left(i \right) - {N^{ - 1}}\sum\limits_{i = 1}^N {{X_{\rm{o}}}\left(i \right)} } \right]}^2}} } }} \end{array} 图 1 2008—2014年(a,c,e)积雪期(NDJ),(b,d,f)融雪期(FMA)中国区域季节平均积雪覆盖率(单位:%) (a, b)MODIS的观测数据,(c, d)CLDAS模拟结果,(e, f)CLDAS-Prcp模拟结果 Fig. 1 Seasonal MEAN of SCF during (a, c, e) snow-accumulating period (NDJ), and (b, d, f) snow-melting period (FMA) in China from 2008 to 2014 (unit: %) (a, b) MODIS observation, (c, d) CLDAS model result, (e, f) CLDAS-Prcp model result 图 2 2008—2014年(a, c)积雪期(NDJ),(b, d)融雪期(FMA)中国区域积雪覆盖率季节BIAS(单位:%) (a, b)CLDAS模拟结果,(c, d)CLDAS-Prcp模拟结果 Fig. 2 Seasonal BIAS of SCF during (a, c) snow-accumulating period (NDJ), and (b, d) snow-melting period (FMA) in China from 2008 to 2014 (unit: %) (a, b) CLDAS model result, (c, d) CLDAS-Prcp model result由2008—2014年中国区域积雪覆盖率逐日CORR( 图 3a , 3b )对比发现CLDAS的模拟结果( 图 3a )CORR较低,普遍在0.4以下,青藏高原和中国南部甚至出现负相关。在东北地区和新疆北部CORR普遍在0.2~0.5,新疆南部及青藏高原地区的CORR基本在0.1以下,说明CLDAS的模拟结果在时间变率上与卫星观测一致性较差;而CLDAS-Prcp的积雪覆盖率模拟结果( 图 3b )在整个中国区域CORR都得到提高,青藏高原区域的负相关基本消失,上升至0.3左右,东北地区的CORR最高可达0.6~0.7,新疆南部地区的CORR提升到0.5左右。虽然CLDAS-Prcp模拟结果的BIAS及RMSE上改进不是很明显,但是与观测相关上改进效果明显,说明改进后的模拟结果和观测一致性较好。
图 3 2008—2014年中国区域积雪覆盖率逐日CORR对比图 (a) CLDAS模拟结果,(b) CLDAS-Prcp模拟结果 Fig. 3 CORR of SCF in China from 2008 to 2014 (a) CLDAS model result, (b) CLDAS-Prcp model result在上述积雪覆盖率季节平均值的对比中,CLDAS-Prcp的表现在东北和新疆区域略高于观测,在BIAS和RMSE上的改进不是十分明显,这有可能是由于MODIS在东北森林地区反演精度有关。东北地区具有广大的山区,森林覆盖率高,蕴育着丰富的森林。 Su et al(2010) 发现,MODIS在森林地区以及山区对雪监测的能力不够,可能会对积雪覆盖率形成低估;MOD10A1在加拿大地区雪季160 d的平均的一致率可达到93%;其中,常绿阔叶林的一致性最低,只有80%,尤其是在融雪期( Simic et al,2003 );中国东北地区MOD10A2积雪产品倾向于低估积雪覆盖面积,林地植被的积雪产品积雪分类精度值变化最大,对于光学遥感,由于林冠对积雪信息的部分屏蔽效应,林下积雪覆盖特征很难被检测到,使得通常的线性光谱混合模型在林区很难得到较高的结果精度,其次云量和积雪深度也是影响影像精度的重要原因( 仲桂新等,2010 )。因此,CLDAS-Prcp模拟值在东北和新疆区域略高于MODIS,有可能是受MODIS本身低估积雪覆盖率的影响。
3.2 对雪深(SD)模拟影响的评估图 4 为2008—2014年中国区域季节平均积雪深度,其中雪深观测来自于CIMISS逐日积雪数据,发现三个积雪区中CIMISS站点观测的积雪主要集中在东北区域和新疆区域,雪深大概在0.05 m以上,其他区域雪深基本处在0.02 m以下。从两种驱动的积雪覆盖率模拟结果中可以看出,CLDAS-Prcp模拟结果的基本形态与CIMISS站点观测数据相吻合,CLDAS模拟结果偏低。
图 4 2008—2014年中国区域季节平均积雪深度(单位:m) (a)CIMISS雪深观测数据,(b)CLDAS模拟结果,(c)CLDAS-Prcp模拟结果 Fig. 4 Seasonal Mean of SD in China from 2008 to 2014 (unit: m) (a) CIMISS station observation, (b) CLDAS model result, (c) CLDAS-Prcp model result从2008—2014年中国区域雪水当量逐日CORR对比( 图 10 )可以看出CLDAS的模拟结果( 图 10a )CORR较低,在东北地区和新疆北部地区CORR普遍在0.2~0.5,新疆的南部以及青藏高原地区的CORR基本在0.1以下,青藏高原部分地区甚至为负相关;CLDAS-Prcp的雪水当量的模拟结果( 图 10b )在整个中国区域内的相关系数都得到了提高,在青藏高原大部分地区的CORR提高到0.7以上,东北地区可达0.9以上,而在新疆南部地区可提升到0.6~0.9左右,在整个中国地区CORR得到很大的提升,负相关基本消失。也可以看出CLDAS-Prcp的模拟结果在青藏高原和新疆的部分地区已由负相关转为正相关。
3.4 对突发降雪模拟能力的评估2008年1月我国南部地区发生了百年不遇的雪灾,河南、湖北、江苏和江西西北部、浙江北部出现大到暴雪,湖南、贵州、安徽南部和江西等地出现冻雨或冰冻天气,4次过程主要的发生时间为:1月10—16日,18—22日,25—29日,31日至2月2日( 沈永平等,2009 )。从2008年1月10日至2月3日CLDAS-Prcp模式平均积雪覆盖率(SCF)与国家气候中心发布的全国冰冻日数( 丁一汇,2010 )对比结果可以看出,模式模拟的积雪覆盖率与积雪日数空间分布一致。在这段时间,长沙及周边地区的积雪覆盖率达到了70%以上,这与全国冰冻日数15~20 d的地区相对应,而西宁、兰州、西安合肥等地处在冰冻日数5~10 d的地区,在CLDAS-Prcp的积雪覆盖率图上可以看到明显对应的空间分布,积雪覆盖率基本在10%~30%。除了新疆北部以及东北地区高估,CLDAS-Prcp的模式模拟结果与全国冰冻日数空间分布基本一致。说明CLDAS-Prcp不仅可以较好地模拟积雪变量,同时能够刻画出强度较大的突发降雪事件。
上述试验采用CLDAS和改进后的降水驱动(CLDAS-Prcp)分别驱动Noah3.6陆面模式对积雪变量进行模拟,对比了降水驱动对模拟结果所带来的影响,评估结果说明CLDAS-Prcp改善了原有驱动在冬季由于低估降水所造成的积雪量偏少的情况,并在几个主要积雪区和中国区域都有明显的改进,同时也印证了驱动数据的质量直接决定着模式对积雪的模拟效果。CLDAS-Prcp驱动改进了模式积雪变量的模拟效果,在不同的区域以及不同变量,改进效果不尽相同。从评估效果可以看出,CLDAS-Prcp驱动模式模拟结果能够较好地反映青藏高原的积雪变量,其积雪覆盖率、雪深、雪水当量的RMSE和BIAS分别为13.1%、2.2%,0.69×10 -2 m、-0.11×10 -2 m,2.49 kg·m -2 、0.447 kg·m -2 ,东北区域模拟与观测时间变率最为一致,积雪覆盖率、雪深、雪水当量的CORR分别为0.416,0.775,0.931,而雪水当量的改进效果最明显,RMSE和BIAS分别减小了54.8%和83.1%,CORR提高了0.47。而从三个积雪变量比较,三者与观测的一致性则表现为:雪水当量>雪深>积雪覆盖率。
通过2008年1月的中国南方暴雪事例,也可以看出,CLDAS-Prcp不仅可以较好地反映积雪变量的时间变率,同时能够较准确地反映突发的强度较大的降雪事件。
积雪覆盖率在东北和新疆区域模拟略高于MODIS,可能是作为检验数据的MODIS在森林地区以及山区对雪监测能力不够,对SCF形成低估造成的( Su et al, 2010 )。因此,后续工作可以通过地表类型进行区域划分,进行针对不同地表类型区域来评估,进一步深刻探讨新驱动在不同地表类型的表现。同时,未来工作可以考虑通过数据融合以及数据同化( 朱国富,2015a ; 2015b )方法对积雪覆盖率SCF的数据质量进行改进,减少驱动和模式结构的不确定性对积雪模拟造成的影响。
Hall D K, Riggs G A, Salomonson V V, et al, 2002. MODIS snow-cover products[J]. Remote Sens Environ, 83(1/2): 181-194. Jordan R, 1991. A one-dimensional temperature model for a snow cover[R]. Technical Documentation for SNTHERM. 89. Special Report 91-16. Hanover: USA Cold Regions Research and Engineering Laboratory: 49. Wu Tongwen, Wu Guoxiong, 2004. An empirical formula to compute snow cover fraction in GCMs[J]. Adv Atmos Sci, 21(4): 529-535. DOI:10.1007/BF02915720 Xie Zhipeng, Hu Zeyong, 2016. Comparison of the snow simulations in community land model using two snow cover fraction parameterizations[C]. Vienna: EGU General Assembly. Yuan Wenping, Xu Wenfang, Ma Minna, et al, 2016. Improved snow cover model in terrestrial ecosystem models over the Qinghai-Tibetan plateau[J]. Agric For Meteor, 218-219: 161-170. DOI:10.1016/j.agrformet.2015.12.004 Zhang Yongfei, 2015. Multivariate land snow data assimilation in the Northern Hemisphere: development, evaluation and uncertainty quantification of the extensible data assimilation system[D]. Austin: The University of Texas at Austin: 1-5.