1 引言

大气边界层(ABL)是人类活动和各项生态环境过程发生的主要大气层区域^[1]。与上层自由大气交界处存在的逆温层所对应高度被定义为大气边界层高度(ABLH)^[2]。ABL作为整层大气中最靠近下垫面的一层,除了受大气层结影响外,还受到下垫面的影响,在重霾天气的形成过程中发挥重要作用。因此,探测ABL的结构特征十分重要,目前无法直接测量ABLH,只能依赖于对其他资料的反演。常用的ABLH求解方法有位温法、梯度法、Richardson数法、曲线拟合法、小波变换法及干绝热法等^[3]。

激光雷达作为一种先进的探测手段,将其观测数据用于反演ABLH早已被认可。张婉春等^[4]采用梯度法利用激光雷达回波信号反演北京2013年1月重霾期间的ABLH,反演结果与探空数据具有很好的一致性,平均差异仅为86 m,同时发现严重污染时ABLH低于500 m,日平均高度约为424 m。吕立慧等^[5]利用2014年11月京津冀地区10台激光雷达数据,采用梯度法获取的ABLH结果和美国全球资料数据同化模式模拟的边界层高度相关性达到0.82,表明梯度法是稳定可靠的,同时发现京津冀地区ABLH在300~900 m之间,并且东南方向较高。王琳等^[6]评估了激光雷达数据反演边界层高度的4种方法,发现这4种方法各有优缺点和适用范围,其中梯度法、标准偏差法、小波协方差变换法比较相近,准确性高但不稳定,而曲线拟合法的稳定性好,但结构精度不高。

此外也有研究表明WRF(Weather Research and Forecasting)模式可以在一定程度上相对准确地模拟出ABLH, 陈炯等^[7]利用WRF模式对北京地区夏季ABLH日变化进行了模拟分析,认为使用高分辨率的模拟基本上能够反映出ABLH的日变化特征。赵世强等^[8]利用激光雷达结合WRF数值模式估算出兰州远郊榆中地区夏季ABLH,指出两种方法反演的ABLH具有较好的一致性,但选取的边界层方案不能同时兼顾热对流边界层和动力作用边界层的模拟。何心河等^[9]的模拟研究指出WRF模式在浓度低时模拟结果较好,而对浓度快速大幅上升和极值的模拟能力还不够。王丽霞等^[10]评估了WRF模式中YSU(Yonsei University Scheme)、MYJ(Mellor-Yamada-Janjic Scheme)和ACM2(Asymmetric Convective Model,Vision 2)三种不同的边界层方案模拟结果,指出三种不同的边界层方案均能够模拟出冬季边界层高度日变化的基本特征,其中,YSU、MYJ和ACM2模拟结果与位温廓线法所得结果对比的平均绝对误差分别为79.56%,124.64%和85.03%,均方根误差分别为110.17,200.5,122.25 m,消光系数分别为0.85、0.76和0.91。

2015年12月北京市先后启动两次重污染红色预警,这是《北京市空气重污染应急预案》自2013年通过以来北京及全国初次启动大气重污染红色预警^[11]。本文针对12月份重污染期间北京市的ABLH进行研究,利用激光雷达反演评估WRF模拟结果,同时分析ABLH与颗粒物质量浓度之间的关系。

2 研究方法

2.1 激光雷达反演ABLH

激光雷达是探测大气气溶胶物理及光学特性时空分布的有效手段,被广泛应用于大气气溶胶、大气成分和云的研究领域^[12-14]。激光雷达可长时间连续不断地自动进行探测工作,其探测数据具有高时空分辨率的特征,因此,可以更好地实时、连续测量ABL结构的时空分布^[15-18]。

2.1.1 激光雷达系统

激光雷达系统总体结构由三部分组成:激光发射单元、接收单元、数据采集和控制单元。激光发射单元由激光器和发射光学单元两部分组成^[19]。接收单元是激光雷达系统的核心部分,需要根据探测目标的不同来设计不同类型的接收望远镜和后继光学单元,以降低探测盲区,提高系统信噪比,减少对大气多次散射光子的接收。数据采集和控制单元是决定激光雷达数据质量的关键部分之一,采集电子学系统的响应速度、带宽、动态范围及信噪比直接决定雷达系统时空分辨能力和对弱信号的解析能力,表1为该激光雷达系统的主要参数。

2.1.2 梯度法

ABL与上部自由大气层之间的覆盖逆温作用使边界层内富含大量的气溶胶粒子,而自由大气层中只含有少量的气溶胶粒子,使边界层和自由大气层中的气溶胶粒子质量浓度发生明显变化,根据气溶胶质量浓度在空间垂直分布的梯度变化规律来确定ABLH。而激光雷达探测到的后向散射回波信号与气溶胶粒子质量浓度呈正相关,因此,可以用激光雷达后向散射回波信号的梯度变化来反演ABLH^[20]。

对激光雷达接收到的米散射后向回波信号进行处理得到的激光雷达距离平方校正回波信号(PRR)为

P(R)R2=CY(R)β(R,λ)T(R,λ),(1)

式中, R为探测高度,C为激光雷达的系统因子,Y(R)为与探测高度R相关的几何因子,λ为激光波长,β(R,λ)为探测高度R处的大气后向散射系数,T(R,λ)为大气透过率。PRR廓线强度对应于相应高度大气气溶胶质量浓度的大小,PRR廓线梯度变化代表大气气溶胶质量浓度的变化,这个梯度变化最大的位置对应ABLH。

梯度法采用回波信号随高度衰减速率大小作为判断ABLH的依据,定义PRR梯度为

D(R)=d[P(R)R2]/dR,(2)

反映PRR在不同高度R处的衰减变化情况。ABLH对应于PRR衰减最快的位置,也就是D(R)取最大值时所对应的高度值^[21]。

在利用梯度法计算边界层高度之前,需对激光雷达回波信号进行校正,以保证准确反演边界层高度。首先,选取16~20 km处回波信号进行背景噪声扣除;其次,将各层高度上的回波信号强度乘上距离的平方再进行距离平方校正;然后,采用基于最小二乘的卷积拟合算法的Savitzky-Golay滤波器对回波信号进行滑动窗口平滑降噪;最后,对底层大气回波信号进行几何重叠因子校正,并将结果与高层大气回波信号进行拼接,得到几何因子校正后的回波信号强度。

2.2 WRF数值模式模拟ABLH

WRF数值模式是由NCAR(National Center for Atmospheric Research)及NOAA (National Oceanic and Atmospheric Administration)联合相关大学和研究机构开发的新一代中尺度预报模式同化系统,给理想化的动力学研究、全物理过程的天气预报、空气质量预报以及区域气候模拟提供了一个公用的模式框架^[7,22]。

2.2.1 模型设置

本研究所用的数值模式为WRF3.7.1,研究区域覆盖京津冀地区。模拟采用双向反馈的三重嵌套网格,模拟中心点为(36.228°N,111.278°E),水平方向采用Lambert投影,水平分辨率分别为36,12,4 km,第一层网格数为100×100,第二层网格数为103×103,第三层网格数为151×127。垂直层为101层,为了更细致地描述边界层结果,边界层内设置较高的分辨率,2 km以下设置了34层,垂直层顶气压为50 hPa。模拟时间段为2015年11月20日至31日(北京时间),其中11月份的数据为预积分数据。模型所需的气象输入数据采用美国国家环境预测中心(NCEP)发布的FNL (final operational global analysis)再分析数据作为模式初始场,网格分辨率为1°×1°,每6 h 提供一次边界条件。表2为模拟所采用的物理参数化方案,模式微物理过程选择Lin方案^[23],长波辐射选择RRTM (rapid radiative transfer model)方案^[24],短波辐射选择RRTMG (rapid radiative transfer model for general circulation models)方案^[24],陆面过程选择Noah方案^[25],边界层采用YSU方案^[26],积云对流采用Grell-Freitas方案^[27],地形数据采用美国USGS 30s的LANDUSE资料,该资料有24种土地利用分类,可为高分辨率数值模拟提供必要的数据基础。

2.2.2 理查森数法

理查森数法(Richardson,简称Ri)^[28]是从WRF模式中计算ABLH的标准方法,总体理查森数的计算公式为

αRig(hi)=g(hi-h1)θv(h1)×θv(hi)-θv(h1)uhi)-u(h1)]2+[v(hi)-v(h1)]2+βu*2,(3)

式中,α_Rig为理查森数值,h_i为模式输出的各层高度,θ_v(h_i)为模式各层的虚位温,u(h_i)、v(h_i)分别为模式各层维向风、经向风,g为重力加速度,β为系数,u_*为摩擦系数。

理查森数法定义ABLH为从第一层开始,理查森数值α_Rig首次大于给定临界理查森数α_Ric的高度。YSU方案将ABLH定义为临界理查森数为零的高度,也即定义为在逆温层通量最小的层所在高度^[29]。

3 结果与讨论

3.1 激光雷达回波信号时空分布与空气质量指数变化趋势

地面监测数据来自全国城市空气质量实时发布平台中北京市奥体中心环境监测站点(40.0031° N,116.407 °E)数据(含下文的地面PM2.5数据),利用位于距离监测站点约4 km的北京市中国科学院大气物理研究所铁塔分部(39.9748°N,116.370°E) 的激光雷达进行连续立体监测。

新实施的环境空气质量指数将空气质量划分成优、良、轻度污染、中度污染、重度污染和严重污染共6个等级。对研究期间每小时的空气质量指数进行统计分析,所用共728条数据(丢失16条数据)。优、良、轻度污染、中度污染、重度污染和严重污染天气所占数据条数分别为123、121、87、60、165和172。其中,33.5%的时间为清洁天气(优、良天气),而66.5%的时间内出现不同程度的天气污染,其中重度污染和严重污染所占比例最大,分别为22.7%、23.6%。

激光雷达回波信号的时空分布与空气质量指数变化趋势如图1所示,激光雷达接收到的返回信号是每微秒光子数,得到最终的回波信号需经过光电转化、背景信号扣除、几何因子校正、距离平方校正及信号平滑滤波等一系列的预处理工作,其中LST表示当地标准时间。修正后的回波信号以不同颜色显示,反映大气中颗粒物在垂直方向上的分布状况,数值越大表明颗粒物含量越多。从观测结果来看,污染物主要集中累积在300 m以下,在高空并未观测到漂浮的颗粒物或输送过来的颗粒物污染源。第一次红色预警(12月8日7时至12月10日12时),污染物迅速清除;而第二次红色预警(12月19日7时至12月22日24时),污染物清除缓慢,持续时间长。

图 1. 激光雷达回波信号的时空分布和空气质量指数变化趋势

Fig. 1. Spatial and temporal distribution of the echo signal monitored by lidar and variation trend of air pollution index

下载图片 查看所有图片

3.2 ABLH特征分析

3.2.1 ABLH演变特征

图2为激光雷达和WRF模式给出的ABLH演变特征。从对比结果可以看出,两者结果趋于一致。WRF模拟与激光雷达反演结果的相关性为0.70,均方根误差为433 m,平均偏差为-264 m。清洁天气情况下ABLH较高,最高可到1000 m左右,而污染期间ABLH低,最低可能不足100 m。两者的极值并没有很好地吻合。

12月份共有两个较长的无污染时间段,第一个阶段为2日19时至6日11时,共持续89 h,激光雷达反演的平均ABLH为637 m(220~968 m),而WRF模拟的平均ABLH为583 m(27~1675 m),相对于激光雷达反演,WRF模拟的平均高度误差为8.5%。第二个阶段为15日10时至17日16时,共持续55 h,激光雷达反演的平均ABLH为660 m(448~880 m),而WRF模拟的平均ABLH为548 m(55~1182 m),平均高度误差为17.0%。

此外,12月份更多的是污染时间段,第一次红色预警前后污染过程(6日12时至11日02时)共持续111 h,激光雷达反演的平均ABLH为287 m(140~488 m),而WRF模拟的平均ABLH为146 m(27~489 m),平均高度误差为49.1%。第二次红色预警前后污染过程(18日04时至24日17时)共持续158 h,激光雷达反演的平均ABLH为370 m(136~850 m),而WRF模拟的平均ABLH

图 2. 激光雷达反演与WRF模拟的ABLH结果对比

Fig. 2. Comparison of ABLH results between lidar inversion and WRF simulation

下载图片 查看所有图片

为277 m(27~902 m),平均高度误差为25.1%。

3.2.2 ABLH日变化

图3为激光雷达和WRF模式给出的ABLH日变化特征曲线图。整体上,与雷达反演结果相比,WRF模拟的相关性为0.76,均方根误差为163 m,平均偏差为-61 m。因为WRF模拟的夜间ABLH相对较低,大多情况下高度都不到100 m,而相同时刻激光雷达反演的结果为300 m左右,从激光雷达观测结果来看,污染物主要累积在近地面300 m左右的大气层,因此根据质量浓度梯度求解ABLH的方法更加符合实际的边界层变化规律。

3.2.3 ABLH平均日变化

图4为ABLH平均日变化结果,整体上,与雷达

反演结果相比,WRF模拟的相关性为0.77,均方根误差为94 m,平均偏差为-69 m。WRF模拟的ABLH具有明显的日变化,凌晨至上午10时ABLH维持在(283±38) m,随后ABLH逐渐上升,于下午14时达到最高值569 m,随后ABLH逐渐降低,直到下午18时,边界层又趋于稳定。激光雷达观测到的ABLH不存在明显的日变化,ABLH从凌晨的382 m缓慢上升,同样在下午14时达到最高值486 m,随后又缓慢降低,在晚上21时达到最小值374 m。

3.3 ABLH与颗粒物质量浓度之间的关系

图5为激光雷达反演ABLH与PM2.5质量浓度的随时间变化的曲线图。可以看出,PM2.5质量浓度与ABLH相关性达到-0.85,表明PM2.5质量

图 3. ABLH日变化

Fig. 3. Daily variation of ABLH

下载图片 查看所有图片

图 4. ABLH平均日变化

Fig. 4. Daily average variation of ABLH

下载图片 查看所有图片

图 5. 激光雷达反演ABLH与PM2.5质量浓度变化

Fig. 5. Variation of ABLH by lidar inversion and PM2.5 mass concentration

下载图片 查看所有图片

浓度与ABLH存在很强的负相关性。如5日ABLH从881 m迅速下降到350 m左右,导致颗粒物质量浓度迅速升高,出现重污染天气。12月份连续出现多次重污染的原因除了北京独特的地理环境位置之外,根据国家气候中心的分析,12月份平均风速不到2 m/s,较常年偏小5%,小风日数(风速小于2 m/s)较常年同期偏多,风速低、小风日数多使污染物快速在近地面累积。其次,北京地区平均相对湿度约为71%,较常年同期偏高36%^[30],颗粒物的吸湿性增长使相对湿度大的天气情况下也会出现污染物质量浓度增加的情况。

4 结论

利用激光雷达对2015年12月份重霾天气下的北京进行连续性探测,采用梯度法对数据进行反演获取ABLH。同时,对WRF数值模式模拟的ABLH结果进行数据评估,发现虽然两种方法的结果具有较好的一致性,但是极值并没有很好地吻合。对于WRF模拟日变化结果来说,相关性达到0.76,均方根误差为163 m,平均偏差为-61 m。同时发现在清洁天气下,WRF模拟的准确性要高于污染天气下的模拟结果。

WRF模式的误差来源主要包括:模拟采用的USGS地形资料与实际地形高度存在一定的差异,使得下垫面与真实情况存在差异;初始场数据虽然采用是经过同化的FNL再分析资料,但气象场具有一定的复杂度,不可能完全用模拟来替代实际的大气状况;同时模式中存在较多的假设条件,这些都会给WRF模拟带来一定的误差。此外,需要注意的是,WRF模式中参数化方案不同对模拟结果具有较大的影响,不同参数化方案的性能分析也是下一步需要探讨的。

正是WRF模拟误差的存在,使得必须通过修正途径来提高WRF的模拟能力,而激光雷达具备长期持续观测能力,通过严格数据质量控制,可以有效保证数据精度和准确度。相对于其他方法,激光雷达反演可以更加方便准确地得到ABLH的时间序列变化,因此,可以利用激光雷达反演的ABLH对WRF等数值模式模拟的ABLH进行诊断,通过调整模型参数化方案,可以在一定程度上提高数值模式模拟的准确性。