随着涡度相关法被广泛用来测定地表与大气间物质和能量交换，通量贡献区分析的概念和方法也应运而生。在地气交换研究中，基本方法都是基于单个空间点的连续观测，并以此来获得感兴趣区域的湍流参数的定量描述［1］。对于面积足够大、下垫面均一的生态系统而言，来自各方向的通量是相同的，因此通过测定冠层上方的湍流通量可以很好地反映真实的地气交换情况，其涡度相关的观测值可以反映生态系统真实的净生态系统交换量。然而由于实际观测中并不能保证观测区域的均质性，不可避免会有下垫面不均一、地形相对复杂等情况，这就增加了观测难度和不确定性，给湍流数据分析、数据质量控制以及数据信息解读等带来困难。通量贡献区函数恰好反映了地表通量来源，即“源强”与观测点之间的函数关系，因此，可以利用通量贡献区分析定性或定量研究通量信息与地面状况之间的关系。 通量贡献区分析首要是求解通量贡献函数，目前适合不同大气条件和地表植被状况的通量贡献区模型主要包括解析模型、拉格朗日随机模型、大涡模拟模型［2-3］。其中FSAM(FluxSourceAreaModel)是由Schmid［4］提出的欧拉解析模型，该模型通过平流扩散的解析解建立起来，依据大气扩散学中的倒置烟羽扩散理论，假设风向沿X轴负方向、侧风向的扩散遵从高斯分布，不考虑垂直方向的通量扩散，同时忽略沿风向的扩散，各风向上的扩散相互独立，该模型假设通量观测位于常通量层并且通量均来自地面，且气流水平均匀［5］。由于该模型原理简单并考虑下垫面的不均匀性，因而得到了广泛的应用［6］。 本文利用FSAM模型［4，7］对位于北京市大兴区永定河沿河沙地的杨树人工林生态系统的涡度相关观测数据进行通量贡献区分析，其目的包括:(1)确定观测塔周围通量贡献区的时空分布特点;(2)分析和解释观测数据的空间代表性。 1研究区概况 研究区位于北京市大兴区榆垡镇林场(39°31'50″N，116°15'07″E)。该地区是海河水系永定河洪积冲积平原，海拔30m，地势平坦。研究区为杨树人工林，主要营造于1998、2001、2003年，总面积为0．8km2，平均胸径为14．5cm，平均树高为16．2m，平均栽植密度为3m×2m。林下植被种类较少，主要为紫苜蓿(MedicagosativaLinn)，尖头叶藜(ChenopodiumacuminatumWilld)，黄花蒿(ArtemisiaannuaL．)等。每年4—10月为生长季。 2研究方法 2．1通量与小气候测定系统 研究区观测仪器包括涡度相关观测系统和微气象梯度观测系统，观测塔高37m，主要观测仪器包括:开路式H2O/CO2红外气体分析仪(LI-7500，Li-Cor，NE)，三维超声风速仪(CSAT3，CS)，净辐射仪(Q7．1，REBS)，安装高度均为20m;气压计(CS105，CS，USA)，翻斗式自动雨量计(TE525-L，CS，USA)，安装高度分别为21m和22．5m;4层空气温湿度传感器(HMP45Cprobe，CS，USA)，安装高度为5、10、15m和20m;3个土壤温度传感器(TCAV107，CS，USA)和3个土壤热通量板(HFT3，Seattle，WA)，均置于地表以下5cm、10cm和20cm处;两个土壤水分观测仪TDR(CS616，CS，USA)置于地表以下5cm和20cm处。 2．2涡度相关通量观测数据处理 本文选取大兴观测站2009年1—12月连续的涡度相关观测数据和气象观测数据。塔上数据为10Hz实时数据和30min在线平均数据，利用数据采集器(产自美国，型号为CRR5000)采集，塔下常规气象数据为30min平均数据，用数据采集器(产自美国，型号为CR23X)进行采集。采集到的数据利用SAS软件进行后处理，处理方法包括:野点剔除、趋势去除、平面坐标拟合［8］、空气密度效应修正［9］、湍流脉动数据的稳态测试［10］和缺失数据插补，最后得到完整的通量序列。根据本文研究需要，在数据统计分析之前，还要进行如下的数据剔除:1)由于降雨时仪器响应问题会出现测量值的不准确或不合理的现象，所以剔除降雨量不为零的数据;2)由于风速较小时湍流混合发展不充分，因此剔除风速小于1m/s的数据;3)由于摩擦风速小于临界值时大气湍流混合程度弱，湍流发展不充分，因此剔除u*＜0．15m/s的数据;4)根据模型适用范围要求，剔除强稳定条件(ξ＞1)和强不稳定条件(ξ＜－1)下的数据。基于本文剔除数据的原则，约保留了60%的数据用于通量贡献区分析，一般情况下50%以上数据就很有代表性。 2．3通量贡献区模型 通量贡献区(footprint)是指对近地面层某一点所观测到的对湍流交换过程有贡献的有效源(汇)区域或者说对观测点的通量大小产生主要影响的表面区域［11］。通量贡献函数或源权重函数是描述表面源或汇空间分布与观测高度zm处观测到的通量信息之间关系的函数，其函数值可以解释为表面区域某一点对观测值的相对权重［7］。为了方便获得某一观测点测定的某一特定源权重的通量来源的最小可能范围，模型引入P水平概念，它表示在通量贡献达到P时的最小区域上的源权重函数的积分相对于整体的比值。本文根据Yang等［13］提出的最小化成本函数的计算公式，利用单层三维超声风速计测定数据和湍流相似性理论求解动力学粗糙度z0。随着涡度相关技术的发展，该方法已经为许多微气象学者所应用［14］。本文依据全年风向风速分布状况确定该年的主风方向，再计算不同风向和大气稳定度情况下的地表粗糙度。通常假设连续10d的z0值是常数，每10d算1个z0值，取连续3个z0的算术平均值作为该月的z0值，再根据需要将结果按生长季和非生长季计算出z0多月平均值。 3结果与分析 3．1风场分析及模型输入参数 2009年全年39．6%的风来自146．25°—213．75°方向，24．55%的风来自326.25°—33．75°风向(图2)，可视这两个方向为该年的主风方向，这主要是因为样地所在地区处于温带季风气候区，盛行风向随季节变化而变化，冬季主要受到来自高纬度内陆极地大陆气团影响，盛行偏北风，而夏季主要是受到变性热带海洋气团的影响，盛行偏南风。此外，19．06%和16．79%的风分别来自33．75°—146．25°和213．75°—326.25°两个方向，相对很小。全年1m/s附近的风速最多，其次是2m/s附近，个别时间风速较大，大约在11m/s，较大风速发生在偏北风方向上(图3)。年平均风速是1．86m/s。模型输入参数均在模型要求的变化范围之内(表1)。生长季时的zm/z0略小于非生长季，说明在生长季粗糙元以上的测量高度较小，即动力学粗糙度较大，虽然粗糙元高度和动力学粗糙度不完全相等，但是高大粗糙元和单位面积粗糙元数量多时，对应的动力学粗糙度就大［16］。在生长季，植物正处于生长阶段植物叶片茂盛，粗糙元平均垂直范围较高，而且单位面积上粗糙元的数量较多，导致动力学粗糙度较大。此外稳定条件下的zm/z0小于不稳定条件下，这可能是因为不稳定条件下风速较大，而一般而言，随着风速的增大，z0增大［16］。比较不同稳定度条件下的zm/L的绝对值可知，不稳定条件时其值较小，而L可视作湍流混合的垂直尺度，因此可见不稳定条件下湍流混合的垂直尺度较大。此外，146．25°—213．75°风向非生长季时和213．75°—326.25°风向时的σv/u*值较大，说明在这两种情况下横向风湍流对通量贡献区的影响比较大。#p#分页标题#e# 3．2通量贡献区空间分布 3．2．1主风方向上的分布 应用FSAM模型，按不同风向输入参数，再根据模型输出参数绘制P等值线分布图，其中等值线P由内向外依次为0．2、0．4、0．6和0．8。主风方向146．25°—213．75°上相同P水平下，非生长季贡献区空间范围整体上略小于生长季时的空间范围。稳定条件下迎风方向上各P水平贡献区的范围大于不稳定条件下的范围，生长季时尤其突出。例如80%通量贡献区，即P=0．8时，ea=333．18m，eb=429．28m，ec=343．14m，ed=488．74m，aa=58．85m，ab=68．83m，ac=55．05m，ad=64．93m，生长季稳定条件下的来流方向范围比不稳定条件下大135m左右，非生长季稳定条件下的来流方向范围比不稳定条件下大80m左右。此外，通量最大值所在位置依次为118．70、144．33、113．06m和141．98m。最大值位置随稳定条件增强而远离观测位置，非生长季比生长季离观测点更远(图4)。在326．25°—33．75°方向上，最大通量值依次位于来流方向115．80、133．1、116．76m和141．95m，距观测点的距离随稳定度增加而变大，且非生长季略小于生长季(图5)。在不稳定条件下，非生长季源区范围显著大于生长季的源区范围;稳定条件下，非生长季源区范围小于生长季源区范围。不稳定条件下的各P水平等值线包围的面积明显小于稳定条件下，而且生长季不稳定条件下的通量源区明显小于其他情况下通量源区，这可能是因为此时湍流扩散较快，涡度相关观测仪器捕捉不到较远处的湍流通量信息。可见，此时大气条件的变化对通量源区的影响较大。此外，两个方向上通量贡献区在稳定条件下均大于不稳定条件下。这主要是因为稳定条件下大气垂直运动弱，湍流扩散速度较慢，通量塔能观测到距离观测点较远位置的通量信息，所以源区面积较大。相比较而言，326．25°—33．75°方向上的源区面积略大于146．25°—213．75°方向上的源区面积。 3．2．2非主风方向上的分布 风向在33．75°—146．25°、213．75°—326．25°之间时，最大通量值依次位于来流方向110．89m、139．76m、121．82m和134．23m。稳定条件下的源区范围均大于不稳定条件下的源区范围，这与主风方向时的情况一致，主要是不同大气稳定度条件下湍流扩散强度不同造成的。当P=0．8时，风向在33．75°—146．25°时，横向范围分别为50—350m和50—450m，垂直方向分别为－75—75m和－100—100m;风向在213．75°—326.25°时，横向范围约为50—400m和50—550m，而垂直于风向的范围为－100—100m和－125—125m。此外，在各P水平上，风向213．75°—326．25°的通量源区范围大于33．75°—146．25°风向的通量源区范围(图6)。 3．3通量贡献区的时间分布 在146．25°—213．75°、213．75°—326．25°和326．25°—33．75°三个方向上，大气稳定条件出现的频率多于不稳定条件，在33．75°—146．25°方向上大气稳定条件出现的频率少于不稳定条件(图7)。结合风向风速在一年中出现的频率(图1、图2)，可以计算得到通量贡献区在时间尺度上的分布情况:在2009年1月—12月的1a里，主风方向146．25°—213．75°和326．25°—33．75°上，生长季和非生长季在稳定和不稳定条件下的贡献区分布比例为15．7%、10．64%、7．64%、8．58%和6．94%、4．37%、9．49%、5．94%，而在213．75°—326.25°，33．75°—146．25°方向上稳定和不稳定条件下的贡献区分布比例分别为7．59%、8．00%和8．32%、6.79%。综合以上数据分析可见，69．3%的通量信息来源于通量塔的偏北风与偏南风方向，其中42．56%的信息来自于偏南风方向上。因此，2009年通量贡献区主要变化发生在塔偏南方50—450m范围内，通量贡献区主要在偏南方80—170m的范围内。一般来说，在大气处于不稳定条件时，观测塔较少测到来自塔偏西方向的通量;而在稳定条件时，观测塔较少测到来自塔偏东方向的通量。另外，生长季时主风146．25°—213．75°方向稳定条件出现的概率最大，非生长季时326．25°—33．75°方向稳定条件出现的概率最大。因此，生长季时，通量数据主要来源于偏南方向50—500m，非生长季时，通量数据主要来源于偏北方向50—400m。 4讨论和结论 本文利用FSAM模型分析了2009年通量贡献区时空分布特点，从而确定观测塔测得的数据的空间代表性。从通量贡献区时空分布来看，本文中生长季的源区面积普遍大于非生长季，与赵晓松等［17］应用FSAM模型对长白山阔叶红松林进行通量贡献区分析的结果不同，而且，最大贡献区在塔偏南方80—170m范围，比赵晓松等的研究结果100—400m范围要小。这可能是因为本文林分高度和观测高度都比较小，而且z0的确定方法不同，赵晓松等［17］采用多层风速数据利用空气动力学方法计算出z0，本文采用的是单层风速数据和最小化成本函数法计算出z0。此外也可以用生态学方法［16］确定z0，但是这些方法都是非常简单的经验公式，多适用于静态研究，而且气象学意义也不明确。不同的确定方法可得到不同的z0值，因此得到的通量贡献区分布也有一定差异。此外，本文研究显示通量贡献区随大气稳定度的增强而增大，贡献区沿来流方向逐渐远离观测点，这与很多学者［17-18］的研究结果相同。观测的范围越远，观测值越有可能受到远处非林区的下垫面的影响，从而造成观测值与实际值的差异。 从通量贡献区模型来看，FSAM欧拉解析模型的算法简化了气流流动的物理过程，使得数学上很简单，减少了数值计算量。最终得到简化的二维通量贡献区，且有较好的准确定［19］。但是该模型没有考虑冠层内部的气流运动和沿气流方向的扩散，使得该模型在应用于高大的植被贡献区研究时会产生一定的误差。G?ckede等［20］研究表明向前拉格朗日随机模型更适用于高大植被观测区研究，但是并不能证明哪种模型的模拟结果更准确，这主要是因为没有参考值，而且研究区复杂性使得模型无法进行验证。此外两种模型都能很好的描述基本的湍流情况［20］。FokenandLeclearc［21］研究也指出尽管FSAM模型较简化而且有一定的局限性，但是依然能够提供有价值的通量信息。#p#分页标题#e# 在不考虑观测高度对通量贡献区影响的条件下得出如下结论:(1)随大气稳定度的增强，通量贡献区增大，通量最大值距观测点的距离增大，通量贡献最少的源区从塔偏西方转为塔偏东方;(2)全年通量贡献区分布较离散，基本在距离观测点500m以内，不稳定条件下较集中在400m以内。此外，生长季的通量贡献区面积普遍大于非生长季;(2)通量观测塔测得的数据有69．3%来自观测塔的偏北方与偏南方，其中42．56%来自于观测塔的偏南方。对观测值产生主要影响的下垫面位于观测塔偏南方80m到170m范围内。生长季时，通量观测数据主要来源于偏南方向50—500m，而偏北方向对观测数据的贡献最少;非生长季时，通量观测数据主要来源于偏北方向50—400m。致谢:感谢北京林业大学大兴课题组全体同学在数据采集和处理上的帮助，感谢北京林业大学同小娟老师、中国科学院地理科学与资源研究所郑宁博士在论文修改中给予的帮助。