1 引言

青藏高原(简称高原)被誉为地球“第三极”(Jane, 2008), 平均海拔超过4000 m, 占我国领土的四分之一, 位于我国西南部, 南邻副热带, 北至中纬度, 为世界最高的高原(马耀明等, 2006)。由于其特殊的地理位置及其地貌特征, 使其在全球气候系统尤其在亚洲季风的形成、爆发、持续时间及强度等方面扮演着重要的角色(杨梅学等, 1998; 张永生等, 1998; 吴国雄等, 1999, 2004, 2005; Ye et al, 1998; Wu et al, 1998; Tao et al, 1998; 樊威伟等, 2018)。大量研究表明(吴晓鸣等, 2013; 王介民, 1999; Ma et al, 2002; 严晓强等, 2018; 谢琰等, 2018), 高原陆气相互作用对区域及全球气候都有着重要的影响。王介民(1999)指出高纬度地区和山地地区是气候变化的敏感区, 特别是位于中纬度高原地区(含大面积冻土)的地气相互作用, 对区域气候和全球气候有着极其重要的作用。高原季风与南亚高压的强度和位置变化密切相关。动态高原季风指数与中上层经向风的相关性分析得出其夏季存在遥相关性, 表明遥相关性是高原季风与东亚季风和南亚季风之间的联系(Xun et al, 2012)。近47年来高原在总体上呈显著增暖趋势, 降水量呈现出增多趋势(李林等, 2010)。高原热力作用不仅对亚洲季风和降水变率有着重要影响, 而且还通过激发类似于亚洲-太平洋涛动的大尺度遥相关, 影响着北美和欧洲以及南印度洋的大气环流和气候(周秀骥等, 2009)。

高原对大气的各种热力和动力效应是通过大气边界层影响自由大气的, 所以对于高原大气边界层尤其是近地层微气象特征的分析研究显得尤为重要(李家伦等, 2000)。受高原独特下垫面状况影响, 当地加热表现出年际和年代际变率(周秀骥等, 2009)。陆气相互作用是指发生在陆地表面的热力、动力、水文以及生物物理、生物化学等一系列复杂过程, 以及这些过程与大气的相互作用过程(苏炳凯等, 1999)。而不同下垫面也具有不同的多样性, 构成复杂, 分布不均, 同时, 也让它们所构成的陆地表面多样且具有一定的复杂性。而地表作为大气能量的主要输入方, 极大地影响着地面与大气之间进行的各种相互作用, 在局部地区或特定时段甚至起着关键性决定作用(Guan et al, 2009)。正因如此, 高原地气相互作用的研究成为近30多年来的研究热点之一, 受到了国内外越来越多的关注。

气候系统对于陆面状态的异常变化十分敏感。而地表特征参数(动力学粗糙度、热力学粗糙度等)在陆面过程中起着十分重要的作用, 是引起气候变化的重要因子(贾立等, 2000)。高原的下垫面存在不同程度的起伏, 这为深入理解高原的地气相互作用带来了一定障碍。起伏的地表可能改变地表的粗糙元及其排列方式, 并导致地表粗糙度的改变。粗糙度的变化也会影响其他地表湍流输送特征参数的变化, 这也可能导致地表通量的变化。并且地表粗糙度是陆面过程模式和气候模式中的重要参数, 其大小在一定程度上反映了近地表气流与下垫面之间的物质和能量交换、传输强度及它们之间的相互作用大小(刘静等, 2007; Irannejad et al, 1998; Shao, 2000; 张强等, 2003; 谢志鹏等, 2017; 杨丽薇等, 2017)。由此可见, 计算地表粗糙度并总结其变化特征, 能够更加深刻地认识陆气相互作用。

那曲地区位于高原腹地, 地处西藏北部, 是长江、怒江、拉萨河、易贡河等大江大河的源头。整个地形呈西高东低倾斜, 平均海拔4500 m以上, 这里海拔较高, 热量不足, 气候严寒干旱, 含氧量仅为海平面的一半, 藏北高原是西藏气候条件最恶劣的地区之一, 是典型的亚寒带气候区。高寒缺氧, 气候干燥, 昼夜温差大, 多大风天气(李凡等, 2004)。但这里有着广阔的天然草原, 所以建立了以那曲高寒气候环境观测研究站为中心的一套完备的中尺度观测网络(吴晓鸣等, 2013), 获得了大量宝贵的观测资料, 为更加准确地描述高原陆气相互作用和大气边界层结构提供了充足的依据。那曲地区下垫面为广阔的高地平原, 附着植被为高寒草地, 但下垫面存在不同程度的起伏, 有一定的复杂性, 为深刻研究高原的地气相互作用带来了一定的困难。

利用中分辨率成像光谱仪(MODIS)获得卫星资料, 通过那曲地区归一化植被指数(NDVI)资料进行研究, 得到动力学地表粗糙度。并且利用该区域3个观测站, 2008, 2010和2012年的大气湍流观测资料和自动气象站观测资料, 应用一种独立的确定地表动力学粗糙度(Z_0m)的方法(陈家宜等, 1993), 分析Z_0m的时间尺度动态变化规律及不同下垫面的不同结果。比对研究地表粗糙度的值是否可靠, 以期为高原地区的陆气相互作用研究提供准确的地表特征参数, 提高高原地区近地层理论研究水平。

2 资料选取和方法介绍 2.1 资料选取

所选区域是以位于高原中部的藏北高原中国科学院西北生态环境资源研究所那曲高寒气候与环境观测研究站(简称BJ)站为中心的200×200 km²的区域。在该区域上布设有北PAM(Portable Automated Meso-net)自动气象观测站(简称NPAM)、纳木错多圈层综合观测研究站(简称NAMC)、BJ三个气象观测站(图 1)。观测点周围的下垫面在较小空间尺度较为平坦, 在较大空间尺度存在一定起伏。使用的数据包括大气湍流观测资料和自动气象站观测资料。

BJ站位于31.37°N、91.90°E, 海拔4509 m。BJ观测点位于季节性冻土区, 植被为高寒草地。站点观测设备有超声风温仪(CAST3)、CO₂/H₂O红外开路分析仪(LI 7500)以及自动气象观测系统等(马耀明等, 2006), 使用BJ站2008和2012年全年的资料。NPAM站位于31°56′N, 91°43′E, 海拔约4700 m。实验场地面平坦、四周开阔。地面被草高15 cm的高原草甸覆盖。该实验站观测设备有超声风温仪和湿度脉动仪, 还包括2层气温和气湿、气压、平均风速、平均风向、地表辐射温度、土壤热通量、6层土壤湿度和温度及辐射等数据(马耀明等, 2001), 使用NPAM站2012年全年的资料。NAMC站位于30°46.44′N, 90°59.31′E, 海拔4730 m, 处于西藏自治区当雄县纳木错乡纳木错湖东南岸, 背靠念青唐古拉山脉, 下垫面为高寒草甸, 属典型的半干旱高原季风气候区。使用NAMC站2010年全年的资料。

选用的土地覆盖数据为GLC2009(Arino et al, 2010)数据。该数据为Envisat卫星2009年反演数据, 空间分辨率为300 m, 分类标准为Land Cover Classification System(LCCS), 该分类标准将全球地表分为23种不同的类型, 所选择区域中包括了其中的14种。由于所选地区的实际情况与GLC2009数据部分不符, 如所选地区并没有农田等下垫面。所以根据除多等(2010)所得的实际土地覆盖类型进行比较, 将灌溉农田、旱作农田、多耕地植被混合、多植被耕地混合、封闭草地和开放草地6个与实际不符的名称, 分别替换为山地草原、灌丛草甸、山地草甸、高寒草原、高寒草甸和高山稀疏植被6个与实际相符的名称。由于树木组成的下垫面占整体区域的比例少仅为0.36%, 所以将不再对常绿针叶林、混交林、多林地草地混合以及多草地林地混合四类下垫面进行研究。

卫星资料来源于MODIS, 该传感器是美国地球观测系统计划所发射的卫星TERRA和AQUA中重要的传感器。MODIS传感仪的波段覆盖可见光到热红外的全光谱, 可以实现对地表温度、地表植被覆盖、大气可降水量、云顶温度等地表和大气状况的探测, 最大空间分辨率为250 m。利用MODIS探测得到的归一化植被指数是来自MODIS的MYD13Q1产品, 每16天提供空间分辨率为250 m全球影像。选用那曲地区2008年、2010年和2012年全年73个资料。

2.2 研究方法 2.2.1 站点资料计算地表粗糙度的方法

利用单一高度超声风温仪的平均风速和湍流通量等测定值, 选择陈家宜等(1993)提出的地表粗糙度计算方案, 并得出其动态变化规律。该方案具体如下(杨耀先等, 2014):

根据Monin-Obukhov相似理论(Monin et al, 1954; Stull, 1991), 含层结稳定度订正函数(Panosky et al, 1984)的风廓线公式为

$ U(z)=\frac{u_{*}}{k}\left[\ln \frac{z-d}{Z_{0 \mathit{\boldsymbol{m}}}}-\psi_{m}(\zeta)\right], $

(1)

$ \left\{ {\begin{array}{*{20}{c}} {{\psi _{\rm{m}}}(\zeta) = 2\ln \left({\frac{{1 + x}}{2}} \right) + \ln \left({\frac{{1 + {x^2}}}{2}} \right) - {{\tan }^{ - 1}}(x) + \frac{\pi }{2}, }&{\zeta < 0}\\ {{\psi _{\rm{m}}}(\zeta) = - 5\zeta, }&{\zeta > 0} \end{array}} \right. $

(2)

其中: Z_0m为动力学地表粗糙度; z为观测高度; d为零平面位移, 取为0.03 m(杨耀先等, 2014), 由平均植被高度0.045 m计算得出; U为平均风速; u^*为摩擦速度; k=0.35为卡曼常数; L=$-\frac{u^{3}}{\left(k \frac{g}{\overline{\theta}}\right) \overline{\theta^{\prime} \omega}}$(Monin et al, 1954)为Monin-Obukhov长度; x=(1-16ζ)^1/4; ζ=(z-d)/L为大气稳定度参数。由式(1)可得:

$ \ln \frac{z-d}{Z_{0 \mathit{\boldsymbol{m}}}}=\frac{k U}{u_{*}}+\psi_{m}(\zeta), $

(3)

利用式(2)~(3)即可通过单一高度的$\frac{k U}{u_{*}}$的观测值和ζ来拟合确定Z_0m。

2.2.2 卫星资料计算地表粗糙度方法

考虑到研究区域主要下垫面为草地, 本文选用Massman模型(Massman, 1997, 1999)反演藏北那曲地区Z_0m。Massman模型计算方法如下:

$ \gamma=C_{1}+C_{2} \cdot \exp \left(-C_{3} \cdot C_{d} \cdot L A I\right), $

(4)

$ n_{e c}=\frac{C_{d} \cdot L A I}{2 \cdot \gamma^{2}}, $

(5)

$ d_{h}=1-\frac{\left[1-\exp \left(-2 \cdot n_{e c}\right)\right]}{2 \cdot n_{e c}}, $

(6)

$ \frac{Z_{0 \mathit{\boldsymbol{m}}}}{h}=\left[1-d_{h}\right] \cdot \exp \left(-\frac{k}{\gamma}\right), $

(7)

式中: C₁=0.32, C₂=0.26, C₃=15.1为模型中的常数, 与体积表面拖曳系数有关; LAI是叶面积指数; C_d=0.2为液面拖曳系数; n_ec为植被冠层内风速廓线消减系数; h为植被高度。在较早的很多研究中得出, 由于高原的高海拔环境, 造成了研究范围所属区域的低温, 以及对植被高度及稀疏程度产生一定影响。根据以往的研究理论基础, 本文认为藏北地区植被高度与归一化植被指数(NDVI)及海拔有关(Chen et al, 2013), 并在原基础上引入了海拔校正因子, 以下为计算公式:

$ \begin{aligned} H=& h_{\min }+\left(\frac{h_{\max }-h_{\min }}{N D V I_{\max }-N D V I_{\min }}\right)(N D V I\\ &-N D V I_{\min }) , \end{aligned} $

(8)

$ h=a c f \cdot H, $

(9)

式中: h_min和h_max为观测点观测到植被高度的最小和最大值; NDVI_max和NDVI_min为观测点的NDVI的最大值和最小值; H是基于植被高度与NDVI成正比这一假设得到植被高度; x为海拔, 由ASTER的DEM产品获得; acf为海拔校正因子(孙根厚, 2016), 用来表征藏北地区海拔对植被高度的影响, 它的形式是:

$ a c f=\left\{\begin{array}{cc}{0.149, } & {x>4800} \\ {11.809-0.0024 x, } & {4300<x<4800} \\ {1.49, } & {x<4300}\end{array}\right. $

(10)

所采用的LAI由MODIS的NDVI计算得到。计算公式为:

$ L A I=\left(\frac{N D V I *(1+N D V I)}{1-N D V I}\right)^{0.5}, $

(11)

3 结果分析 3.1 实测资料结果分析

图 2给出了那曲地区不同年份各站点的地表粗糙度的时间变化特征。从图 2中可以看出, Z_0m自2持续增加, 到7月、8月达到最大值, 其中2012年的BJ站与NPAM站显示出7月较8月稍大, 而2010年NAMC站和2008年BJ站8月比7月大。BJ站2008年的7月可能产生了降水过程, 导致7月Z_0m值小于6月。8月之后Z_0m开始减小, 12月的值大约与1月的值持平。总体来说, 各站点的Z_0m度变化呈春季到夏季逐月增大, 夏季到冬季逐月减小的变化特征。(其中NAMC站无2010年2月和10月实测数据)

3.2 MODIS反演地表粗糙度的时空变化特征

图 3为2008年BJ站点周围的200×200 km²的地表粗糙度分布。从图 3中可以看出, 2月较1月的Z_0m有所减小, 可能是由于降雪、温度等原因, 导致Z_0m很小, 并且持续减小。由于积雪融化温度回升等原因, 从2—5月Z_0m呈现逐月缓慢增大的变化特征, 6—8月Z_0m呈现逐月迅速增大的变化特征。图 3中从6月开始呈现出大量地表纹理, 说明下垫面的复杂性, 不论是地表或是植被都会对Z_0m产生重要影响。3—8月由于高原夏季风所带来的湿度、温度、气压等明显的变化, 促使地表植被的高度及覆盖度都在增大, Z_0m于8月达到峰值。尤其是5—6月的变化十分显著, 可能是由于6月对应高原夏季风的开始, 导致降水的增强, 从而植被的生长加快, Z_0m的迅速上升。6—8月由于持续不断的降水以及气温的不断上升等原因, 导致植被生长很旺盛, 但生长至成熟后不再发生变化, 对应为图中Z_0m最大值不变, 但由于充足的水分, 在这3个月Z_0m的大值区逐渐扩大, 8月达到最大范围。9—12月, 随着高原夏季风减退, 温度、湿度等都逐渐减小, 与高原夏季风盛行之时相比, 不再适宜植被的生长, 植被对Z_0m的贡献削弱, 地表植被高度及覆盖度逐渐减小, 因而, Z_0m不断减小, 并且大值区范围也逐渐减小。

3.3 MODIS反演地表粗糙度的时空

图 4和图 5分别给出了2010年和2012年BJ地区逐月地表粗糙度反演结果, 从图 4和图 5中可以看出, 2010年与2012年那曲地区的Z_0m同样呈现出了1—2月减小, 自2月始, Z_0m开始增大, 6月起, Z_0m增长迅速, 并于8月达到全年度峰值, 之后, Z_0m开始减小。

从图 3~5中也可看出, Z_0m不仅随时间尺度变化, 同时也随空间尺度变化。Z_0m在不同下垫面上呈现出的不同变化趋势。值得注意的是2008年11月, 那曲地区Z_0m整体较小, 普遍低至1 cm。据历史资料得知此时那曲地区有大范围降雪过程, 大范围降雪造成草甸下垫面均质平坦, 而降雪落地后, 容易形成块状、散乱、不连续分布的下垫面。并且在后文可以得到, 以冰雪为下垫面的区域地表粗糙度不超过1 cm, 这与历史天气过程相符合。所以认为降雪导致了11月的Z_0m非常小。并且11—12月Z_0m呈现出增长的趋势, 这是可能是由于温度、土壤未冻结或是其他原因, 导致积雪的融化, 进而地表粗糙度呈现出增长趋势(周扬等, 2017)。

3.4 卫星资料反演结果的评估

所选3个站点的下垫面都是高寒草甸, NPAM站点资料计算结果全年都大于其卫星资料反演结果, 只有9月与10月非常接近, 变化趋势相似, 站点资料计算结果最大值达到5 cm而卫星资料反演结果最大4.5 cm, 最大差异出现在5月, 达1.7 cm(图 6)。NAMC站点资料计算结果在4-11月与卫星资料反演结果非常接近, 但1月、3月以及12月的卫星资料反演结果却明显大于站点资料计算结果, 最大差异出现在在1月, 达1.5 cm。2008年BJ站点资料计算结果全年都大于其卫星资料反演结果, 其中站点资料计算结果在1—4月、7—11月与卫星资料反演结果非常接近, 但是在5月、6月和12月出现了较大的差异, 最大差异出现在5月, 达1.8 cm。2012年BJ的站点资料计算结果全年与其卫星资料反演结果相一致, 但在3—6月站点资料计算结果大于其卫星资料反演结果, 其他时间则是站点资料计算结果小于卫星资料反演结果, 最大差异出现在6月, 达1.1 cm。从图 6中可以看出, 整体情况下, 1月、2月、3月、11月和12月, 站点资料计算结果季节变化趋势与卫星资料反演结果一致。但4—10月的站点资料计算结果要大于卫星资料反演结果。从图 6还看出, 站点观测得到数据计算出的Z_0m普遍大于卫星资料计算出的结果, 这可能是由于平均的平滑作用导致的。值得注意的是, 2—7月, 根据确定地表粗糙度的独立方法计算而得单点Z_0m值显著增加, 而利用卫星数据所得结果却无显著增加。卫星结果显示出1—5月、11月、12月基本稳定在2 cm以下, 6—10月才有变化; 这与卫星资料反演的平均滑动作用有关。总的来说, 站点计算得到的结果普遍大于卫星反演得到的结果。

从Z_0m散点分布(图 7)中可以看出, 卫星数据与站点数据计算得到的地表粗糙度有明显的正相关关系。观测结果与反演结果的相关系数除了2010年NAMC站偏小, 其他的都要大, 说明下垫面的平均结果平滑了不同地区相同下垫面结果, 进一步说明该地卫星反演结果与站点计算的结果较为相似。但NAMC站点的结果却与其他站点的结果不同, 其与下垫面的结果比较相关系数为0.83288, 与该地卫星反演结果比较相关系数为0.61554。或是因为纳木错观测站离湖泊较近(1 km)受湖陆风等地方性小气候影响较大的原因。图 7结果均通过P=0.05的F检验。表明站点数据计算结果与卫星资料反演结果没有显著差异。

4 不同下垫面地表粗糙度的变化特征

根据地面覆盖资料GLC2009结合当地实际情况将那曲地区200×200 km²区域分成了10种不同的下垫面(Arino et al, 2010):山地草原、灌丛草甸、山地草甸、高寒草原、高寒草甸、高山稀疏植被、城镇用地、裸地、水体、冰盖与雪盖。

从那曲地区不同下垫面的Z_0m的月际变化(图 8)中可以看出, 可将那曲地区14种不同下垫面可以分为4类。第一类为城镇用地, 占该区域的0.07%。全年所有时期都比其他的下垫面的Z_0m度大, 而其Z_0m的变化也非常大, 可能是由于所选区域的城镇用地下垫面非常少以及人为原因造成的各种不规律变化。第二类为茂盛草地类, 包括灌丛草甸、山地草原、高寒草原和山地草甸4种下垫面, 占该区域的62.49%。这4种下垫面Z_0m的变化曲线较为相似, 与其他的下垫面相比, 只小于城镇用地。第三类为稀疏草地类, 包括高山稀疏植被、高寒草甸和裸地3种下垫面, 占该区域的33.74%。这3种下垫面的Z_0m又较为相似, 处于中等高度, 裸土的Z_0m位于这几个下垫面Z_0m的最低点, 高寒草甸的Z_0m变化较平稳, 其Z_0m受外界影响小。第四类为冰雪类, 包括冰盖与雪盖、水体2种下垫面, 占该区域的3.7%。这三种下垫面的Z_0m呈现出另一种现象, 全年的变化幅度都比较小, 而且这几种下垫面的Z_0m本来也小, 全年只有6月中旬时超过1 cm, 其他时间Z_0m都小于1 cm。从Z_0m的多年平均的季节变化可明显看出, 由于地表粗糙度大小的不同可将下垫面分为四类。并且15月Z_0m的变化是非常小的, 从5—8月达到峰值, 之后11—12月又降低到之前1—5月的水平。由于2008年11月的降雪, 可能导致了11月偏低的结果, 事实上往年11月地表粗糙度应高于12月。

从图 8中也可以看出, 那曲地区除了第四类下垫面的地区外, 其他地区的Z_0m变化为12月减小, 2月后开始增长, 至8月达到峰值, 8月后开始减小。但可从图 6中明显看到有几处Z_0m产生明显变化的阶段, 4月初、5月中旬、7月初、8月底、9月下旬几个时期。其中8月底的变化最为明显。各下垫面Z_0m平均变化超过2 cm。9月下旬的变化幅度也比较大, 平均变化超过1.5 cm。另需指出, 7月初变化比较特殊, 因为此次变化导致了水体和冰雪的Z_0m的明显增大。

值得注意的是, 因为云遮蔽等因素的影响, 2008年的5月、8月、11月; 2010年的8月、9月; 2012年的4月、7月, 整体Z_0m都有明显变化产生。由此得到了三年平均图上的8月、11月两次非常明显变化。其中11月根据其他的气象资料得出是由于降雪造成的, 而由2008年与2010年两年共同造成的8月变化的原因, 通过分析水体和冰雪类下垫面Z_0m突然增大, 所以认为是由于降水导致的。结合几次Z_0m的变化可分析得出降水、降雪、积雪会使下垫面变得较为均一平坦, 并使Z_0m也相对减小。

5 结论与讨论

通过对藏北高原那曲地区地表粗糙度的计算并分析其变化特征, 同时对比了卫星资料的反演和实测资料的计算结果, 得出以下主要结论:

(1) 卫星资料的反演结果和实测数据的计算结果基本一致, 均表明地表粗糙度整体呈现自28月持续增大, 8月达到峰值后开始减小, 至次年2月达到最低值。地表粗糙度的这种月际变化规律与高原夏季风带来的气象要素的变化存在诸多的联系。其中卫星反演结果2-5月的地表粗糙度增长缓慢。

(2) 通过卫星资料反演的地表粗糙度变化特征, 可根据地表粗糙度的大小将下垫面分为四类, 由大到小分别为城镇类、茂盛草地类、稀疏草地类和冰雪类。其中茂盛草地类占区域62.49%, Z_0m最大可达6 cm, 稀疏草地类占区域33.74%, Z_0m最大可达4 cm左右, 冰雪类占区域3.7%, Z_0m最大不超过1 cm。

(3) 比较卫星反演计算结果、站点卫星反演结果和实测资料计算结果, 得到几种结果均有明显正相关关系, 且站点卫星反演结果与实测计算结果有较好的拟合性。由于反演的平均滑动作用, 卫星反演资料较实测计算结果偏小。可以使用该方法计算区域的地表粗糙度结果并为模式提供真值进行模拟。

通过对地表粗糙度的研究, 有利于获得区域上的地表特征参数值, 为模式输入提供地面真值, 提高模式在高原地区的模拟水平, 加深对陆气相互作用过程的认识。