纪小芳,龚 元,郑 翔,鲁建兵,冯 明,庄家尧,叶立新,刘胜龙,方万力,王 丹,何雪凯,姜 姜,*
1 南方现代林业协同创新中心,南京林业大学林学院,南京 210037 2 南京林业大学生物与环境学院,南京 210037 3 凤阳山-百山祖国家级自然保护区凤阳山管理处,龙泉 323700
在全球气候变化的背景下,2002—2011年间CO2浓度的增加速率达到了(2.0±0.1 )ppm/a,其对温室效应的贡献大约78%[1]。森林生态系统作为主要的碳汇之一,在陆地生态系统中扮演着重要的角色。过去几十年大量的研究围绕森林生态系统的碳循环展开,试图了解和掌握其碳排放动态,为分析和解决温室效应等全球环境问题寻找途径。
涡度协方差(Eddy covariance,EC)方法常用来测量生态系统尺度上地表与大气的CO2气体交换。涡度协方差也称涡度相关法是一种多用于测量植被下垫面与大气间热量、物质和动量交换的工具,目前为国际上主流的基于微气象理论的通量观测技术,20世纪90年代以来,随着涡度相关法的发展和应用为直接观测不同生态系统CO2排放和吸收提供了技术方法[2]。早期的涡度相关技术多应用于农田、湿地、草地和森林等几个关键陆地生态系统与大气之间的CO2交换研究[3],但受观测高度,大气边界层高度,大气稳定度等环境因素的影响,通量塔传感器所测得的通量值仅代表下垫面的一定区域[4],因此在应用涡度协方差法时,需要确定其空间代表性和通量源区[5]。目前通量源区的模型主要有解析模型、大涡模拟模型、拉格朗日模型和闭合模型等几类。解析模型的假设基础较多,理论上仅适用于下垫面平缓的区域,主要采用梯度扩散理论、二维平流扩散方程以及相似理论得出通量贡献区,较著名的有Kormann and Meixner(KM)模型、FSAM模型和Horst-Weil模型[6- 9];大涡模型最初用于大气和环境科学的研究,但其物理机理、计算比较复杂繁琐,耗时,消耗大量的存储空间,不适宜于长期通量观测数据下的计算[10- 11];拉格朗日模型是基于拉格朗日粒子扩散的数值模拟,理论上严格考虑了扩散的均匀分布约束,可以正确的反映出非均匀湍流中的扩散,较为著名的有Hsieh模型以及Kljun模型[12- 13]。2015年由Kljun[14]等提出的三维通量足迹模型更加适用于长期连续且多时间序列的通量塔观测源区计算,运算速度更快。目前应用FSAM模型[5,9,15- 23]和KM模型[24- 27]的研究较多,但是对Kljun模型[28- 29]实际进行运用的研究较少。对通量源区的主要研究结果发现生长季的源区分布在任何状态下均小于非生长季[5,24- 25,30- 31];但是也有研究表明在大气不稳定状态下时,非生长季的源区面积大于生长季[32];在稳定大气条件下的通量源区长度要显著大于不稳定条件下[33];不同时间其通量源区的大小有所不同[16],随着研究的深入与实际需要针对下垫面不均一的森林生态系统碳通量的监测方法也逐渐成为研究的热点[34]。Schmid H P 等[35]对位于尼日尔萨赫勒地区的虎灌木(粗砂和灌木交错分布)进行研究发现从单一位置进行测量不能代表生态系统的空间湍流通量。Sogachev A等[36]利用SCADIS模型对山脊不同位置的通量值进行研究发现,相比较于同等高度,背山脊处的通量所受到的干扰最大,而通量塔的最佳建设位置在山脊顶部。而目前使用的最多的为ART (Agroscope Reckenholz Tanikon) Footprint Tool模型(基于KM模型)来实现对空间异质性不同土地类型对所测通量所占比例的确定[12,27]。
本研究利用Kljun通量足迹模型(后简称Kljun模型)和ART Footprint Tool(KM模型),对浙江凤阳山通量塔2017年全年的观测数据进行分析,探讨该森林生态系统在不同时间跨度、不同大气条件下的通量源区变化情况以及不同土地类型的对通量观测值的贡献,对研究区通量观测值的空间代表性作出解释,为其他森林生态系统CO2通量研究和类似下垫面的通量源区评估提供服务和参考。
研究区位于中国东南部的浙江丽水龙泉市凤阳山自然保护区内,该区建于1975年,现有管理面积15171hm2[37- 38]。位于东经119°06′—119°15′、北纬27°46′—27°58′。保护区森林覆盖率达到了90.8%,属于中亚热带温暖湿润气候区,同时受海洋性气候和季风影响较大,研究区内雨量充沛,多年平均降水量达到2400mm,降水集中于4—6月,占全年的80%,湿度大,雾多;年均气温12.8℃,极端高低温分别为30.2℃、-12.5℃,月均气温最高28℃左右,最低6℃—13℃,年蒸发量达到1170 mm以上,无霜期275 d,有效积温约6500℃[39- 40]。四季分明,一般3月底4月初入春;7月初入夏;8月中旬入秋;11月中旬入冬,光照资源丰富[41]。土壤类型主要为黄壤土,分布在海拔800 m之上的高山坡地。地形以山地为主,地势较高,沟壑交错,保护区的黄茅尖,海拔高达1928 m,为江浙第一高峰。保护区内的通量塔于2016年10月搭建完成,位于研究区中心,其所在下垫面为阴坡24°,森林遍山野岭,风浪区无限大,树冠高度大部分分布均匀,符合Kljun模型的使用条件[42],周围林地以木荷、杉木和木楠混交为主,平均冠层高度约15 m,平均林龄约40 a。植物资源丰富,形成多优势种结构特征。人工杉木林和人工柳杉林分布在海拔1400—1500 m左右,木荷在海拔300 m,900 m和1500 m处均有分布,黄山松在海拔600 m和1300 m处均有分布,多脉青冈和黄山木兰在海拔1500 m左右有分布[43]。
图1 凤阳山研究区土地利用类型地形图Fig.1 Topographic map and land use types in Fengyangshan research area
研究区中心建有高40m的观测铁塔,塔上装有通量和梯度系统。通量观测系统的传感器有数据采集器(CSI CR300,Campbell Scientific Inc.,USA),三维超声及CO2/H2O分析仪(IRGASON,Campbell Scientific Inc.,USA),四分量辐射表(CNR4,Campbell Scientific Inc.,USA),空气温湿度(HMP155A,Campbell Scientific Inc.,USA),土壤热通量板(HFP01,Campbell Scientific Inc.,USA)。配套的气象观测系统为6层梯度观测(空间高度2、8、16、24、32、40m,土壤深度10、20、30、40、60、90cm),梯度系统测量大气中不同高度的温度、湿度、风速、风向以及不同深度土壤含水量和土壤温度,适用于不同的下垫面和大气条件,是边界层气象、农林气象、大气环境监测最普遍运用和最基本的观测手段。主要传感器有土壤温度传感器(TCAV)、土壤水分传感器(Campbell CS616,Campbell Scientific Inc.,USA)、风向传感器(MetOne 020C,Campbell Scientific Inc.,USA)、数据采集器(CSI CR300,Campbell Scientific Inc.,USA)、6层风速传感器(MetOne 010C,Campbell Scientific Inc.,USA)、HOBO雨量传感器(Onset RG3-M,Campbell Scientific Inc.,USA)和空气温湿度传感器(HMP155A,Campbell Scientific Inc.,USA)。整套观测系统于2016年10月修建完成。研究人员每隔三个月去现场进行数据下载和维护仪器等工作。本研究选取凤阳山通量塔2017全年的连续通量观测数据。采集到的数据经过在线程序(EasyFluxTM-PC, Campbell Scientific Instruments, USA, https://www.campbellsci.com/)处理,EasyFlux程序中所用代码均与EddyPro (Li-COR,NE, USA, https://www.licor.com/)软件所用源代码一致。主要处理过程包括:野点去除、坐标旋转、WPL校正、建立数据质量等级指标(1—9)等最后得到30 min时间间隔的通量和微气象数据序列[44- 45]。之后对数据进行筛选,主要是去除降雨及降雨前后1h数据;剔除红外分析仪信号强度低于0.8的数据;去除质量控制等级标记为“9”的数据[46]。基于本文剔除数据的原则,2017年的通量数据有效率为51%,来进行通量贡献区分析,一般通量数据有效率超过50%则具有代表性[26,28]。
依据凤阳山当地的实际气候状况,季节划分以3、4、5月为春季,6、7、8为夏季,9、10、11为秋季,1、2、12月为冬季[47]。
一般通量观测塔上红外分析仪所测得的通量值为某一时刻迎风方向上对观测值产生影响的下垫面空间代表区域的碳源或碳汇强度,对观测点通量值有贡献的下垫面区域即为通量贡献区[48- 49]。以观测点通量塔为原点(0,0)建立坐标,x轴的正值方向代表上风距离:
(1)
式中Zm为有效观测点高度,f为碳源或碳汇的转换函数即footprint函数,代表了表面上某点(x,y)对Zm处观测值的贡献率密度,Qc代表的是表面碳源或碳汇的点源强度,Fc为在Zm处所测得的通量值,R代表对通量值有贡献的下垫面区域[7]。
Zm=Zreceptor-Zd
(2)
式中Zreceptor为通量塔相对于地面的观测高度,Zd为零平面位移高度,采用森林冠层高度的2/3[50- 51]。
而侧风积分足迹函数fiy则主要受到上风距离(x),有效动力学高度(Zm),大气边界层高度(h),摩擦风速(u*)和垂直风速脉动的标准差(σw)等几个参数决定。以上参数除了h需要单独计算,主要与L(Obukhov,奥布霍夫长度)有关,详细计算方法Kljun文章中给出[14],其他的都可以通过EasyFluxTM在线程序处理获取[46]。由量纲分析(Π定理),将以上变量有可能构成的无量纲参数组如下所示:
Π1=fiyZm
(3)
(4)
(5)
(6)
侧风积分足迹的无量纲函数F*可以由为上风距离X*表示,即F*=φ(X*),由X*=Π2Π3-1Π4,F*=Π1Π3-1Π4可得到:
(7)
式中的a,b,c,d为通过后向拉格朗日随机粒子扩散模型(LPDM-B)不断的实验计算出的各拟合参数,均与粗糙度(z0)有关。用Zm/L来判别大气稳定度,其中L为Obukhov长度(可由EasyFlux在线程序计算获取),当Zm/L>0时,大气为稳定状态,而当Zm/L<0时,大气为不稳定状态[21]。一般在分析中常用贡献率密度的P水平等值线所包围的区域(一般取P=0.8或0.9)来表示EC系统的观测范围[21],这里参考Kim等[42]基于Kljun模型的森林生态系统通量源区的研究方法,本研究选取P=0.8[42],Kljun等[42]还提供了该模型的在线数据处理网站(http://www.footprint.kljun.net/index.php),该研究下载了该模型的开源Matlab (MathWorks,USA,https://www.mathworks.com/)函数代码应用于浙江凤阳山针阔混交林森林生态系统的通量足迹研究。
从图2和表1可以清晰的看到凤阳山森林全年和四季的风速和风向分布频率,除秋季的主风向为东北风之外,其余季节与全年的主风向均为西南风。结合图1中的地形可以看出,风向受到地形的很大影响,由于通量塔位于山坡上,所以其西北方向的风受到山坡高度的影响,很难传输到通量塔所在位置,而西南方向海拔从通量塔位置开始逐渐增高,因此不会完全阻挡风,给风的传输提供了可能,东北和西南方向也是如此。同时从表1中的白天和黑夜不同风频数据可以看出,全年昼夜风向相反,白天和黑夜的主风向分别为东北风和西南风,可能受到山谷风的影响[52]。
图2 研究区风向风速玫瑰图Fig.2 The wind rose map in the study area
表1 不同时间不同风向的风频
利用Kljun模型计算出30 min时间间隔的90%的通量贡献率的通量足迹,并绘制出其所对应的上风距离分布图和累计频率概率图。从图中可以明显的看出,通量贡献区最远点的分布比较集中,距离通量塔观测点最远的距离为7000m,主要的通量测量值来自于东北和西南方向,结合全年的风向风速图2可以发现,风向和风速是影响通量贡献区的主要因素。
对全年的贡献率在10—80%的源区分布进行分析发现,在不稳定条件下的源区面积明显小于稳定条件。当大气处于不稳定条件时,其通量贡献率在80%的源区面积为0.89 km2,源区长度在323.62—839.62 m间,0—90°方向上源区长度达到最大,此时大气做垂直剧烈运动,物质不稳定输送速度很快,涡流更易形成和传播,EC系统所测得的通量值来自上风向更短距离,导致其源区面积较小;而当大气处于稳定条件时,其通量贡献率在80%的源区面积为4.88 km2,源区长度在672.69—2390.49 m之间,在180—270°方向上源区长度达到最大,空气的湍流运动较弱,物质上下垂直扩散速度比较缓慢,EC系统可以测得来自较远距离的涡流,因此所代表的源区面积较大[16,24,53]。
图6表示在40m观测高度Kljun模型评价浙江凤阳山针阔混交林分别在春、夏、秋、冬季节在大气条件稳定和不稳定条件通量贡献率为80%时的源区分布图。可以看出,通量贡献区存在季节差异,不同季节通量贡献区的面积和主要分布方向都有差异,稳定条件下的源区面积要明显大于不稳定大气条件下,与全年的源区分布规律相似,且源区的方向和形状均与风向风速图相一致。春季在不稳定和稳定条件的源区面积分别为0.64 km2和4.89 km2,源区长度分别在277.08—707.41 m和622.73—1887.58 m之间,分别在0—90°和180—270°方向源区长度达到最大。夏季在不稳定和稳定条件下的源区面积分别为1.07 km2和2.58 km2,源区长度分别在442.04—1822.53 m和329.40—1880.01 m之间,源区长度最大值均出现在180—270°方向。秋季在不稳定和稳定条件下的源区面积分别为1.42 km2和5.37 km2,源区长度分别在359.69—1212.98 m和441.73—2717.48 m之间,源区长度最大值出现在0—90°和180—270°风向。冬季不稳定和稳定条件下的源区面积分别为0.70 km2和5.54 km2,源区长度分别在122.43—756.86 m和650.90—2468.66 m之间,源区长度最大值风向为0—90°和180—270°风向。明显可以看出,源区长度的最大值出现方向主要在东北(0—90°)和西南(180—270°)方向,这完全符合全年源区的分布规律。在大气条件不稳定的条件下,源区长度分布不超过2000 m,源区范围从大到小排列为:秋季、夏季、冬季和春季,可以看出生长季的贡献区在大气不稳定条件下较大;当大气条件稳定时,源区长度分布不超过3000 m,源区范围从大到小排列为:冬季、秋季、春季和夏季,此时生长季的源区分布最小。
图3 通量足迹最远点分布图Fig.3 The distribution of the flux footprint furthest point
图4 通量足迹最远点累计频率图Fig.4 The probability chart of the flux footprint furthest point
图5 全年源区分布Fig.5 Flux footprint through 2017(0,0)点为观测点,东南西北方向分别为x正轴、y负轴、x负轴和y正轴
图6 不同季节的源区分布Fig.6 Flux footprint during different seasons
图7 侧风积分函数图Fig.7 Crosswind integral function diagram
从图7中可以看出通量贡献率峰值均位于传感器附近,不超过150m,且随着大气条件从不稳定到稳定,其在上风向的倾斜更少。在不稳定条件下,其湍流强度较高导致化合物向上传输和较短的传输距离与传输时间。在稳定条件下侧风向距离是明显高于不稳定条件下,该结论与Kim[42]等研究结果一致,与前文通量贡献区在不稳定条件下的源区分布范围小于大气稳定条件的结果一致。Zm/L<0,表示大气不稳定状态,当Zm/L>0,表示大气稳定状态[22]。从表2可以看出在春、夏、秋、冬季,当大气条件不稳定时,Z0、h、sigmav和u*值均大于稳定条件下的,在稳定条件下,夏季的z0、h、sigma和u*最高,冬季的L值最高;在不稳定条件下,夏季的z0、sigmav、u*和L的绝对值最高,冬季的h最高。
下垫面的土地利用类型对通量值有重要影响,所以能将两者进行结合,可以对源区不同植被对通量观测值的贡献有更加直观的认识。本研究选取1120×1180 m的研究区范围利用ART Footprint Tool[6,54- 56]进行定量化分析,将该范围分割成3304个20×20m的方格,确定其格内的土地利用类型并求出其对通量值的贡献率,最后叠加得到不同植被类型在全年对通量贡献的百分占比(图8)所选范围全年总贡献率达90%。从图中可以看出源区贡献从大到小依次为针阔混交林、阔叶林、建筑用地及道路、杉木林、毛竹林、柳杉林、黄山松林。结合图1可以看出通常土地类型的面积越大,其占比就较大。毛竹林的面积是柳杉林的两倍,但是其贡献百分比却是柳杉林的3.7倍,可以看出毛竹林在对该研究区的碳预算功能有重要作用。该区域内的建筑用地主要为生态定位站的实验楼及标准气象场,其本身的CO2通量理论上可以忽略不计。
表2 侧风积分函数输入参数
图8 不同土地利用类型的源区贡献率占比Fig.8 Proportion of source area contribution by land use type
通量贡献区主要受到风向和风速的影响,但是自然条件中的风向风速几乎时刻在发生变化,为了更好的对贡献区进行分析,本文分析了不同大气层结条件和不同季节下的贡献区分布。如结果(图6)中所示,通量源区的位置在不同季节有相应的变化,主要受到风向风速的影响,因为夏季的风向风速变化幅度较大,所以其对应的源区分布范围无论是在大气稳定还是不稳定的条件下,都分布较广,冬季时的大气较为稳定从图6中的冬b可以看出,该季节在西南方向的源区分布范围明显较其他季节大,而在大气稳定的条件下夏季的源区范围较小,主要是夏季植物茂盛生长,导致湍流形成比较迅速,从而导致源区的范围较小,该结论与朱明佳[53]等的研究结果一致,在大气不稳定的条件下源区范围从大到小排列为:夏季、秋季、春季和冬季,可以看出生长季的贡献区在大气不稳定条件下较大,该结论与龚笑飞[16]等研究结论一致;无论从全年还是季节的尺度上来看,在大气不稳定的情况下,其源区分布范围均小于大气稳定条件,主要是由于当大气条件不稳定时,空气上下交换较快,垂直运动强烈,扩散也较快,最终使通量值所代表的源区范围较小,然而在大气稳定条件下,空气垂直产生湍流的速度减缓,扩散减弱,源区范围可以追溯到上风方向更远处,该结论与前人[5,16,24- 25,31- 32,57- 58]研究结论均一致。本研究区的源区距离最远可达到7000m(图3),这可能也与下垫面、冠层高度和仪器的观测高度有关,本研究中观测高度为40m,冠层高度大约15m,超过了冠层高度的2.5倍。与魏远[5]等人所得出源区最远点距离为3500m结果不同,主要其仪器高度为25m,冠层高度为16m,由于仪器测量高度不同从而导致源区最远点的距离有差异;而顾永剑[17]等对湿地生态系统进行源区分析得出其源区长度不超过300m,主要由于湿地的植被均为较低矮的草类植被,仪器观测高度仅为4.8m,因此造成其源区长度小于本研究。而下垫面会对摩擦速率产生影响,从而对通量贡献区产生影响,在未来的研究中,特别是在地形及其复杂的情况下,探索研究出能够考虑到下垫面,地形变化等方面的源区模型并加以运用。本研究发现建筑用地及道路的源区分布贡献率也会达到5.48%,仅次于常绿阔叶林。虽然这部分对CO2通量理论上可以忽略不计,但人为活动产生的CO2排放会低估对该区域的森林碳汇功能。
本文利用2017年全年浙江凤阳山针阔混交林生态系统通量塔观测资料,在保证数据质量的前提下,对其所测通量值所代表的源区分布进行分析,从而对研究区通量观测值的空间代表性作出解释,为其他森林生态系统CO2通量研究和类似下垫面的通量源区评估提供服务和参考,研究发现:
研究区内全年盛行东北风和西南风,风向在东北方向的占总风频的34.16%,风向在西南方向的占总风频的46.96%;而风向是影响源区分布的主要因素,研究区的源区分布主要也在东北和西南方向,从源区的季节分布也可以看出,各季节的源区分布方向和形状基本与该季节的风向一致;当通量贡献率达到90%和80%时,通量观测值分别来自距通量塔7000m和3000m范围,且各季节的侧风向积分函数峰值均位于传感器附近,不超过150m。在大气层结不稳定条件下,湍流强度较高;稳定时,源区范围大于大气不稳定条件,主要由于大气不稳定条件直接影响了大气运动的强度,湍流形成速度与传播距离;季节尺度上,当大气条件不稳定时,源区范围从大到小排列为:秋季、夏季、冬季和春季;当大气条件稳定时源区范围从大到小排列为:冬季、秋季、春季和夏季源区贡献从大到小依次为针阔混交林、阔叶林、建筑用地及道路、杉木林、毛竹林、柳杉林、黄山松林。