1 引言(Introduction)

全球范围内，虽然城市面积只占到陆地面积的2.4%(Potere and Schneider， 2007)，但是全球约75%~80%的二氧化碳排放都与城市有关(O′Meara，1999)．城市化已成为环境变化和改变碳循环的重要组成部分．因此了解城市发展与碳循环之间关系对于评估未来大气中温室气体浓度变化和帮助减轻气候变化带来的影响十分重要．

涡度相关法(EC)作为一种有效的可直接测定湍流通量的方法已经广泛应用在生态系统CO₂和能量通量观测中(Baldocchi，2003)．然而涡度相关法在城市通量观测中的应用还很稀少，而且已有的观测只覆盖了很少的城市类型(Velasco and Roth，2010)．在城市区域展开EC观测，需要注意的一个问题是城市异质的下垫面特征对通量观测的影响(Grimmond et al., 2004)．例如北京325 m气象铁塔不同方向的地表粗糙度z₀在2.1~6.3 m之间(Al-Jiboori and Hu， 2005)；伦敦市区z₀在0.27~0.87 m，平均建筑高度z_h则在5.6~8.8 m 之间(Wood et al., 2010)，城市下垫面的异质性对EC观测高度和传感器方向的选择有着重要的影响(Velasco and Roth， 2010).过去10年中大多数EC观测集中在发达国家城市，在发展中国家的观测非常少见，只有Velasco et al(2005)在墨西哥城，Burri et al(2009)在开罗的短期观测，以及Song and Wang(2012)、 Liu et al(2012)、Song et al(2013)在北京城市不同高度内进行的年尺度的观测．目前城市环境下的 EC 通量观测仍显不足但开展迅速，准确评估城市CO₂通量还需适用于城市条件下的观测和数据处理方案．Velasco and Roth(2010)指出了城市CO₂通量观测需要注意的问题：一定的观测高度以避免由单一建筑、树木等导致的微尺度变异，建议涡度相关传感器应高于惯性副层；合适观测平台以避免其结构对流场的改变从而导致附加的通量贡献．Kotthaus and Grimmond(2012)针对大部分城市通量观测高度偏低，通量易受微尺度(10¹~10² m)短时高强度的浓度信号的影响，发展了微尺度人为源辨别技术从而将微尺度排放源从局地尺度(10²~10⁴ m)中分离出来.Nordbo等(2012)发展了城市环境下热通量数据处理方案，首次给出了城市环境下显热和潜热通量的涡度相关方法检测限．上述对城市EC观测方法的研究加深了我们对此方法在城市适用性方面的理解，但仍需对城市复杂环境下通量的计算和质量控制进行系统分析以应对不断开展的城市通量观测．本研究选取北京325 m气象铁塔展开涡度相关梯度观测．分析探讨城市复杂下垫面和大气湍流特性对通量计算、校正以及对质量控制的影响，对比不同的计算、校正方案在城市条件下的适用性，为城市通量观测提供参考．

2 观测实验(Observation experiment)

本研究所述观测于2008—2009年在中国科学院大气物理研究所325 m气象铁塔上进行(图 1)．铁塔位于北京市北三环路与北四环路之间(39°9′N，116°4′E)，海拔高度为48.63 m.在观测铁塔2 km半径内，地势相对平坦，地表覆盖以建筑为主，其地表覆盖率λ_b为0.65，其次为植被和道路，地表覆盖率λ_v和λ_r分别为0.23和0.12(Song and Wang， 2012)．区域内建筑高度分布不均，南边平均建筑高度为50 m而其他方向为20 m(李倩等，2003).在观测铁塔20 km半径内，平均建筑高度为18.3 m，地表覆盖率λ_b、λ_v和λ_r分别为0.65、0.21和0.14(Miao et al., 2012)．

325 m铁塔涡度相关通量观测系统为架设在47 m(2008年1月开始观测)、140 m和280 m(2008年6月开始观测)高度的开路涡度相关系统(OPEC).OPEC系统主要有开路非色散红外气体CO₂-水汽分析仪(IRGA，Li7500，Li-CorInc，Lincoln Nebraska，USA)、三维超声风速仪(CSAT-3，Campbell，Utah，USA)和数据采集器(Model CR5000，Campbell，Utah，USA)组成.超声和红外气体分析仪的采样频率为10 Hz，实时观测的数据记录在数据采集器.为了尽量减小通量塔身对来流方向空气运动的扭曲，将红外水汽-CO₂分析仪和三维超声风速仪安装在2.5 m 长的支臂上，并朝向盛行风向.除了通量观测系统外，在铁塔8、15、32、47、65、100、120、140、160、180、200、240和280 m 同步开展了温度、湿度、风速，风向等气象因素观测.2009年6月开始在47、140和280 m高度辅以净辐射观测(CNR1，Kipp&Zonen，Delft，Netherlands)．

水汽及CO₂通量通过垂直风速脉动(w^′(m · s^-1))和标量物质浓度(c^′(μmol · m^-3))脉动的乘积(协方差)在一定时间间隔内的平均值来计算．

式中，上横线表示时间平均，N表示时间平均内的采样数目，脉动值为瞬时观测值与平均值的差(例如:c′_i=c_i－).

3 结果与讨论(Results and discussion) 3.1 通量计算 3.1.1 平均周期

基于涡度相关技术观测地气间CO₂/H₂O 通量，需要考虑通量变化特征和微气象的相关原则来确定适宜的数据平均周期，其原则是尽可能的包含各种频率的湍流成分，在这个意义上，数据平均周期越长越好，但是CO₂/H₂O 通量可能有着明显的日变化，因此平均周期还必须可以分辨CO₂/H₂O通量的日变化特征；其次平均时间周期还必须能分辨短周期的零星事件的影响.与自然生态系统通常取半小时平均周期一样，目前绝大多数城市CO₂/H₂O 通量也都取此平均周期(Velasco and Roth， 2010)．然而，目前我们并不确定对通量有影响的所有低频成分的来源，因此还需要评价低频成分对城市通量观测结果的影响．假定大气处于稳态，利用O_give函数可以确定合适的平均周期(Foken and Wichura， 1996).O_give函数协谱从高频到低频的累积积分(Stull，1988)，可表示如下

式中，f₀为可分辨的最低频率，f_nyquist为nyquist频率．所有频率对应的O_give函数值等于相应时间系列的协方差值Co．O_give函数在某个频率处集中为常数，此频率就可以传换为湍流通量计算的平均周期．图 2显示了47 m高度1 d内不同时刻CO₂，水汽和显热通量的O_give函数，最长周期为90 min，图中每一条曲线代表一天内每90 min的采样周期．可以看出渐近线的频率界限值在1 d内变化不大，大部分O_give函数在接近15 min平均周期处O_give函数就趋向一个稳定值，且30 min至45 min的平均周期内基本没有波动，但部分CO₂、水汽O_give函数在更长的平均周期内有比较明显的发散，这与在自然生态系统观测到的有明显不同，自然生态系统中可观察到120 min平均周期内O_give函数仍然保持不变(Sun et al., 2006)．与CO₂、水汽稍有不同的是，温度的O_give函数在较长的时间周期内相对平稳，这可能是由于温度的日变化主要受太阳辐射的影响，城市人为热释放的时间变异对温度的影响并不显著.CO₂和水汽O_give函数的发散发生在城市交通的早晚高峰，说明在城市这种源排放时间变异大的系统，过长的通量平均周期并不适合，因此在通量计算中取30 min平均周期．从上图还可以看出，在观测高度10 Hz的采样频率足以分辨出湍流输送的高频部分．

3.1.2 坐标旋转

利用超声风速仪观测风速时，仪器的设置不当(而这几乎是无法避免的)会产生一个明显的垂直分量，每一度的倾斜误差会造成大约3%~4%的标量通量误差和大约14%的动量通量误差(Dyer et al., 1982)．因此在地形起伏较大的森林或者地表粗糙程度变化较大的城市进行涡度相关观测时，必须考虑坐标旋转来进行倾斜校正．利用3层高度(47、140和280 m)的通量观测系统，我们分析了旋转角度和平面拟合对CO₂通量观测的影响．考虑到观测塔周围不同的地表粗糙度，按照每30°划分为12个流场坐标，设每个流场坐标内的平均风场在一个与地面平行的平面内，并确定此平面，然后将每一时次的水平风速u、v，垂直风速w旋转到此平面.垂直旋转角Φ为旋转后平均w=0 时所旋转的角度.图 3显示了不同高度的垂直旋转角度Φ随不同风向的变化，可以看出在47 m，垂直旋转角度Φ(8.3±2.1)°较140和280 m高，其最大值(23.5°)所对应的风向(180°)的平均建筑高度为50 m(图 1)，这种大的旋转角度可能是由于高的建筑群改变了流场的方向并且产生了一个人为的旋转角度所致．随着观测高度的增加，城市建筑对流场的影响减弱，角度Φ也随之降低，140和280 m平均Φ值分别为(4.5±0.9)°和(7.3±1.7)°．值得注意的是在280 m高度风向为300°时存在一个较大Φ值(15°)，其原因还需进一步研究．一种可能的原因是此方向风速较大导致塔体顶端摆动，使得超声风速仪记录了明显的垂直速度分量．

城市复杂地形可能导致不同坐标旋转方案计算的通量有明显差异．图 4显示了二次坐标旋转和平面拟合方法计算的2009年1月平均CO₂通量随风向的变化．整体上采用平面拟合计算的CO₂通量比二次旋转稍低．对应于垂直旋转角度Φ的变化，两种方法所计算通量的差异随观测高度变化，其中最大的差异出现在47 m处建筑物分布发生明显变化的方向(180°~200°)．分析表明这是由于受建筑物的影响，在此高度平面拟合方法计算的垂直风速脉动较二次旋转计算的低，从而导致对湍流交换通量的系统性低估．对于本研究所述的城市下垫面，选择不同的两种坐标旋转方案导致的差异在15%以内，这种差异在城市通量计算中不能被忽视，且应在评估城市CO₂收支的不确定性时给予考虑．鉴于平面拟合方法受地面建筑物的影响尤其在观测高度较低时，因此建议在城市通量计算中使用二次或者三次坐标旋转．如果采用平面拟合方案，则须按照风向分别拟合(至少划分为两个平面，其中一个包含主要的建筑群)，而不能简单的假设为均一下垫面而采用一个平面来拟合，尤其在使用诸如Eddypro(Li-CorInc，Lincoln Nebraska，USA)和TK2(Mauder et al., 2008)等软件直接计算通量时需特别注意．

3.2 通量数据评价 3.2.1 稳态检验

通量的涡度相关观测需要稳态条件，即是指所有的统计参数不随时间发生变化(Panofsky and Dutton， 1984；Kaimal and Finnigan， 1994)，非稳态对湍流通量观测影响十分严重．Foken和Wichura(1996)提出的一种稳态检验方法在涡度相关技术测定的湍流通量数据质量评价中得到广泛应用(Foken and Wichura， 1996；Aubinet et al., 2003)．这种方法通过比较一个时间序列内统计参数值和这个时间序列里的几个时间间隔内统计参数的平均值来进行稳态检验．例如，对于一个垂直风速w和变量x(水平风速或者标量浓度)协方差时间序列(30 min)，其数据长度为N(10 Hz的采样频率，N=18000)，将其分割为M=N/L=4~8个短时间间隔内序列的数据(通常选取N/L=6，即5 min的10 Hz数据则L=3000)．L个短时间间隔内信号w和x的协方差可用方程(3)表示：

为进行稳态检验，需要计算N/L个短时间间隔信号w和x的协方差的算术平均值，可以用方程(4)表示：

对于全周期时间序列信号w和x的协方差可以利用方程(5)表示：

因此，可以用方程(6)定义稳态测试比值为：

如果方程(6)的稳态测试比值小于20%~30%，可以认为湍流通量处于稳态条件(Foken and Wichura， 1996；Aubinet et al.，2000)，如果比值小于60%，认为湍流通量数据质量可以接受(Aubinet et al., 2000)．

为研究不同稳态大气湍流过程对城市涡度相关湍流通量数据的影响，根据式(6)对2008年1月的47 m通量数据进行了稳态检验．结果表明47 m处有22%的通量数据满足稳态条件(稳态测试比值<30%)，48%的通量数据质量可以接受(稳态测试比值<60%)．图 5比较了不同稳定度条件下30 min平均通量和5个6 min平均通量的差异，图中ξ为稳定度参数，ξ<0表示不稳定的大气状态，ξ<0表示稳定的大气状态．可以看出非稳态发生并没有明显的日变化趋势，在各个时间段出现的频率基本相同．这与自然生态系统非稳态倾向于夜间出现的结果(宋涛，2007)并不相同，这可能是在城市环境下，白天易受局地源短时水平平流影响所致．这种水平平流由短时间内点源释放的烟羽或者风向改变传输的烟羽产生.值得注意的是，来源于同一风向上大于30 min的平流不能被稳态检验所识别．

3.2.2 湍流积分统计特性检验

近地层湍流方差相似性检验可以判断大气湍流是否能够很好的发展与形成，是否符合近地层大气湍流运动的莫宁-奥布霍夫相似性理论．通过湍流方差相似性关系检验可以发现非均质地形条件下的一些典型效应(Aubinet et al., 2000)，如果由于障碍物或者仪器自身导致的机械湍流，则湍流方差相似性关系观测值会大于模型预测值，这是因为湍流的机械产生项永远促进湍流的发展．与稳态检验类似，如果湍流方差相似性关系的观测值与模拟值相差不大，可以认为湍流通量数据质量是可以接受的，因此湍流方差相似性检验也可以作为湍流通量数据分析和控制的标准．

由于目前缺少CO₂和H₂O通量的方差相似性函数，因而方差相似性检验主要集中在风速分量w和u．垂直风速湍流方差相似性关系预测值和基于莫宁-奥布霍夫相似理论的预测值的比较可以用式(7)表示：

式中，(σ_w/u_*)_measured和(σ_w/u_*)_model分别为实际观测和预测值.理论模型采用Panofsky等(1977)给出的相似性函数：

式中，c₁、c₂为常数，L为莫宁-奥布霍夫长度(m)，Z为观测高度(m).

为研究不同的大气湍流过程对城市涡度相关湍流通量数据的影响，使用2008年1月47 m处通量数据进行了方差相似性检验．图 6为方差相似性检验数据质量分级分布频率，由图可知与稳态检验结果相似，较差质量等级的数据所占比例较少．高标志位数据出现的时间大多在夜间和下午(图略)，大多数异常数据都能够用质量分级标志位来解释．但是值得注意的是并非所有质量分级标志位高的数据都为异常数据，这说明根据质量分级标志位剔除数据而后进行数据插补(gap filling)时需谨慎，因为将标志位高但质量好数据剔除，然后利用日平均法(Liu et al., 2012)或者神经网络插补(Song and Wang， 2012)会引进新的不确定性．在夜间大气稳定层结的条件下，几乎所有的涡度相关技术应用上的限制都会发生，包括仪器本身和气象因素的限制(Massman and Lee， 2002)，这些限制是导致夜间通量数据质量不高的主要原因．

研究者总是希望有尽可能多的高质量数据用于计算交换通量的季节和年总量.利用摩擦风速u_*阈值剔除的弱湍流交换数据是广泛应用一种数据质量控制和评价的方法(Saito et al., 2005)，为了评价质量控制(稳态和湍流方差相似性检验)后的数据表现，将两种质量控制后的结果进行了比较．之所以选用摩擦风速进行数据质量控制，是因为它代表了湍流的机械产生项，它的大小代表了湍流混合的强度.考虑到城市高的冠层高度，这里选取0.3 m · s^-1作为摩擦风速阈值.图 7显示通量数据随u_*的变化，可以看出，夜晚在低于u_*阈值情况下，有约64%的高质量数据，如若按照u_*阈值而剔除，将使通量数据序列有很大的缺失，而使年总量的估计带来不确定性．此外在高于u_*阈值情况下，由于近地层处于非稳态，有相当数量的通量数据没有通过稳态和湍流方差相似性检验而被剔除．这说明仅仅采用u_*检验来控制观测的通量数据似乎是不够的．因此建议在城市通量观测中，通量数据的质量评价应以稳态和湍流方差相似性检验为主．

3.2.3 平流对通量观测的影响

将湍流通量与存储通量之和作为涡度相关法计算的地气间交换通量是目前常用的方法.由于忽略了平流或者其他因素，这种建立在水平均质假设基础之上的方法可能会导致通量观测的系统误差(Finnigan et al., 2003)．为了将平流项应用到涡度相关观测中，Yi等(2000)在Wisconsin森林447 m高塔上采用多层通量观测，发现30～122 m间森林CO₂平流(水平加上垂直平流)损失高达27%，采用线性外推方法粗略的估计30 m处CO₂平流损失约为10%．采用与Yi等(2000)相似的方法，我们观察了不同高度涡度相关可观测项(湍流通量+存储通量)之差，并且假设各层通量之差是由总平流F_ad(水平平流F_adh +垂直平流F_adv)引起的，这里所述的垂直平流F_adv是指垂直的大气运动对标量物质的垂直输送，有别于湍流运动对标量物质的垂直输送．一般，在复杂地形或非均匀下垫面时，热力或地形强迫往往会引起局地环流或中小尺度运动，产生量级较小但作用十分关键的垂直运动. 而这些垂直运动往往会输送一定热量或者标量物质，形成垂直平流通量.尽管在城市区域无法完全排除由各个高度对应的通量贡献区内不同的CO₂排放源所导致的差异，但是通过选取合适的数据仍然可以就平流输送对通量观测影响得出一些有用的参考．下面的分析只选用2008年8月份数据，这是因为8月处于夏季没有供热用高架源，另外，此期间严格的减排措施尽可能的减少了工厂高架源．

图 8显示了140 m与280 m间-ΔF_ad的差异，可以看出最大-ΔF_ad出现在夜间约为15 μmol · m^-2 · s^-1，这种明显的差异表示夜间280 m平流对于通量贡献比140 m大，这种平流贡献主要是垂直平流，因为夜间大的浓度梯度的存在．例如，假设夜间存在0.05 m · s^-1垂直平均速度，那么只需夜间 5×10^-6的浓度差就能产生约10 μmol · m^-2 · s^-1垂直平流贡献，事实上夜间经常存在大于5×10^-6浓度差异(图略)．白天由于浓度梯度很小，垂直平流可以忽略，因此白天140 m与280 m间-ΔF_ad的差异可认为主要是水平平流引起的，类似的变化模态也存在47 m与140 m之间(图略)．2008年8月140 m与280 m累计的ΔF_ad占280 m处累计通量的约10%，而同期47 m与140 m累计的ΔF_ad占140 m处累计通量的约33%，这说明在上层，标量物质交换观测更为可靠．

为了评价各层之间总平流差异ΔF_ad有多少与水平平流相联系，我们首先定义水平风速在高度Z₁到Z₂上的积分(Yi et al., 2000)：

如果我们假设水平浓度梯度是常数，那么两层之间水平平流通量差为(Yi et al., 2000)：

可以看出，如果两层平流通量差ΔF_ad与平均风速有关，那么可以认为水平平流对平流通量有贡献．从图 8可以看出在280~140 m之间，IU与ΔF_ad有着较好的变化趋势．回归分析表明两者可决系数R²为0.48(p<0.05)．白天，垂直方向上CO₂混合均一，垂直平流通量很小(图 9)，但是水平梯度依然存在，IU与ΔF_ad较好的相关关系说明白天水平平流而不是垂直平流占据平流通量的主导．

垂直平流可写为式(11)：

式中，是整个控制体积内CO₂的平均浓度，w_r和c_r 为高度r处的平均垂直风速(m · s^-1)和CO₂浓度(×10^-6).为了定量地描述垂直平流与水平平流相互关系，我们取平均垂直风速为常量-0.01 m · s^-1，为天气尺度下沉气流特征速度(Yi et al., 2000)．这样就可以通过垂直浓度梯度和式(11)估算垂直平流F_adv，图 9显示了140 m和280 m间的垂直平流F_adv对于总平流的贡献，可以看出夜间垂直平流贡献最大，而白天在零值附近．真实的垂直平流的日变化模态可能与图 9所示有所不同，因为280 m高度处平均垂直速度可能比底层大．此外白天对流的存在，垂直速度可能会比较大而且更加波动，会使垂直平流日变化模态更加复杂，但也因为对流的存在，垂直速度的增加并不会明显增加白天的垂直平流，因为对流混合使浓度梯度变小．

上述结果对高塔观测中平流项的重要性进行了评估，然而平流产生的原因复杂，非平坦地形，地表覆盖类型的改变以及中尺度环流均能导致平流，使得平流的估计存在大的不确定性，因此上述的结果只是一个粗略的估计并且在一定程度上是325 m观测塔的特定结果.

3.4 通量收支的不确定性评价

由于随机和系统误差，自然生态系统涡度相关长期观测的不确定性在±20%，这种不确定性有可能导致长期观测的通量结果有源和汇的差异(Goulden et al.，1996).目前鲜有针对城市区域涡度相关长期观测不确定性的研究，只有Velasco等(2005)对CO₂通量随机和系统误差，以及Nordbo等(2012)对显热和潜热通量的随机误差的分析.尽管对通量观测的误差传播很难精确量化，但误差对日平均通量影响的敏感性分析可以帮助量化通量结果的不确定性(Velasco et al., 2005)，如果误差是随机的，其大小随着用来计算平均通量数据量的增加而以N^-1/2的形式减少，其中N为通量观测的次数．这样，随机误差大小可以通过观察所计算的净通量的收敛程度来估计．与随机误差不同，涡度相关观测的系统误差，并不会随着观测数据的增加而改变，这是因为作为对观测数据持续性的补偿或者放大，系统误差以简单的线性形式增加(Moncrieff et al., 1996)，而且系统误差通常难以检测．

图 10a显示了2009年冬季每0.5 h数据加上不同比例随机误差p_r后的结果，其中水平实线代表观测的日平均CO₂通量31.1 μmol · m^-2 · s^-1，在水平实线上下的对称曲线束表示平均通量的整体随机误差，可以看出这些曲线的大小依赖于p_r的大小和给定的天数，其中随机误差随着数据集的增加而减小．使用这种技术可以获得通量最佳估算所需的观测天数．例如，假设每0.5 h通量有+20%的随机误差，可以看出1 d数据不足以分辨出31.1 μmol · m^-2 · s^-1的净排放，另一方面，10 d的数据可以分辨出2.56 μmol · m^-2 · s^-1的净排放，100 d可以分辨出0.77 μmol · m^-2 · s^-1的净排放．图 10b显示了每半个小时数据加上不同比例系统误差误差p_s后的结果，与图 10a一样，水平实线代表观测的日平均通量.虚线代表了对每半小时通量给定不同系统误差(横坐标)后的“实际”通量(纵坐标)．对于观测的日平均通量31.1μmol · m^-2 · s^-1，假设每半个小时通量有+20%的系统误差，那么实际通量则为37.4 μmol · m^-2 · s^-1．如前所述，涡度相关观测通量包含随机和系统两种误差，如果分别考虑每半个小时通量数据存在着± 20%的随机和系统误差，则冬季CO₂通量的随机和系统误差分别是±0.8 μmol · m^-2 · s^-1和±6.22 μmol · m^-2 · s^-1．将各误差分量的平方求和再开方后(Moncrieff et al., 1996)总误差约为±6.3 μmol · m^-2 · s^-1，约为日平均排放通量的20%．对于长期的通量观测，使用上述的方法进行误差分析时应避免所选取的数据集有明显的季节变化，因此建议在讨论城市CO₂收支时给出每个季节CO₂排放的不确定性．

4 结论(Conclusions)

1)采用涡度相关方法观测城市CO₂和能量通量时，不同的坐标旋转方法所计算通量差异在15%以内，这种差异随着观测高度的增加而减小．平面拟合方法受地面建筑物的影响明显，建议在城市通量计算中使用二次或者三次坐标旋转．

2)采用摩擦风速阈值的方法判定城市夜间通量数据质量会导致大量数据缺失，给年总量的估计带来不确定性．城市通量观测中，通量数据的质量评价应以稳态和湍流方差相似性检验为主，辅以功率谱和协谱分析．

3)城市通量观测中，平流项的贡献不能忽略，平流输送的贡献可占至观测通量的33%．