哈尔滨市冬半年污染传输特征及潜在源区分析
2021-04-20李树岭隋亚男徐盛荣景学义
李树岭, 赵 玲, 花 丛, 隋亚男， 张 鑫， 徐盛荣, 景学义
(1.黑龙江省龙云气象科技有限责任公司气象院士工作站，哈尔滨 150028； 2.哈尔滨市气象台，哈尔滨 150028；3.黑龙江省气象台，哈尔滨 150030； 4.国家气象中心，北京 100081； 5.哈尔滨市环境监测中心站，哈尔滨 150076)
引 言
近年来，随着社会经济的快速发展、城市化进程的加快，大气污染已经成为城市面临的突出的环境问题之一[1-4]，尤其是我国的北方地区，秋冬季节环境空气质量受大气污染影响更为严重，频发的重污染天气已对社会公众的健康构成威胁[5-7]。
一个地区空气污染程度主要取决于气象条件、当地污染源排放量及外来污染物的传输影响。已有研究[5,8-9]表明，风速小(或静风)、高湿、逆温偏强等气象条件对污染天气的形成有利，而外来污染物的输送影响，需结合气象条件、局地地形条件及区域污染物分布等进行分析。HYSPLIT模式和GIS技术结合的TrajStat软件提供的轨迹聚类方法，是分析污染物传输轨迹的重要工具。许多学者[10-13]以此为基础,对国内大气污染较重的城市和地区的污染物传输问题进行了研究。王喜全等[14]通过对京津冀城市群和辽宁中部城市群先后出现的灰霾污染天气个例分析指出，辽宁中部城市群出现的灰霾过程一定程度上是京津冀地区灰霾污染跨区域输送影响的结果，并指出了跨区域影响的污染物传输通道。王芳等[15]研究发现，广州市大气污染较重时段主要受特殊气象条件和珠江三角洲地区周边城市排放源的影响，本地源排放与周边城市污染物输送叠加使大气污染加重。严晓瑜等[16]分析发现，位于新疆、甘肃、内蒙古、青海的西北源区及四川、陕西的东南源区，是影响银川市PM10和PM2.5浓度的两个主要潜在源区，西北气流轨迹和东南气流轨迹也是影响银川市污染物浓度的最重要的输送路径。在张磊等[17]研究的基础上，花丛[18-19]等通过对京津冀地区冬半年污染物传输特征分析，构建了传输指数，并对传输指数进行改进，改进后的传输指数能够较好地反映外来传输作用对本地污染物浓度的影响。
哈尔滨是我国纬度最高的省会城市。在每年的秋冬季节，由于受城市供暖燃煤排放的增加、周边秸秆焚烧及外来污染源传输的影响，哈尔滨已经成为我国污染最严重的省会城市之一。目前对影响哈尔滨市空气质量的外来源地及输送路径规律的研究还较少。
本文利用HYSPLIT后向轨迹模式及潜在源区贡献法，结合哈尔滨市2014-2017年冬半年PM2.5浓度监测资料，对主要污染源传输特征进行统计分析，总结出对哈尔滨市PM2.5浓度影响较大的传输路径和污染物潜在源区，同时结合污染物监测数据构建传输指数，为哈尔滨市空气污染预报预警、空气污染成因评估及防治调控的决策提供参考。
1 资料与方法
1.1 资料来源
本文所用PM2.5浓度监测数据来源于哈尔滨市环境监测中心站。选取哈尔滨市辖区11个国控点[20](图1)监测数据，取其平均值代表其空气质量状况。气象数据来源于中国气象局，为地面及探空观测资料。后向轨迹模拟采用NCEP(National Center for Environmental Prediction, 美国国家环境预报中心) 全球资料同化系统GDAS(Global Data Assimilation System) 提供的1°×1°分析场资料作为初始场。
图1 哈尔滨市环境监测站分布示意图
1.2 后向轨迹模式介绍
后向轨迹模拟采用美国国家海洋和大气管理局空气资源实验室研发的HYSPLIT轨迹模式。该模式是欧拉和拉格朗日混合型扩散模式，具有处理多种气象要素输入场、多种物理过程和不同类型污染物排放源功能的较为完整的输送、扩散和沉降模式 , 已经被广泛地应用于多种污染物在各个地区的传输和扩散的研究中。
本文利用HYSPLIT后向轨迹模式对哈尔滨(45.75°N、126.77°E)2014-2017年冬半年(2014年10月-2015年3月、2015年10月-2016年3月、2016年10月-2017年3月、2017年10月-2018年3月)共计24个月的逐6 h气团进行后向轨迹模拟，轨迹计算起始高度为500 m，每条轨迹模拟时长为24 h,时间分辨率为1 h。取35°-60°N、100°-135°E的区域作为研究范围，将其划分为0.25°×0.25°的网格。
1.3 后向轨迹聚类分析
聚类分析是按照样本的相似性、亲疏度将样本进行分型划类的一种统计分析方法。后向轨迹聚类是根据气团轨迹的空间相似度，即传输速度和方向，对所有到达模式受点的气团轨迹进行分组聚类，以判断受点不同时段主导气流方向和污染物潜在来源。本文利用TrajStat软件中的角距离聚类算法，对冬半年到达哈尔滨的气流轨迹进行聚类分析，得到冬半年影响哈尔滨的不同输送气流类型。
1.4 潜在源区贡献法
本文在分析污染物潜在源区时，采用了PSCF(Potential Source Contribution Function)即潜在源区贡献法[21]。这是一种应用广泛的污染物源地统计方法，通过结合气团轨迹和某污染物要素值，给出可能的排放源位置。PSCF函数表征经过某一区域的气团到达观测点时对应的要素值超过设定阈值的条件概率。由于哈尔滨冬半年重污染天气中首要污染物以PM2.5为主，结合《HJ 633-2012环境空气质量指数(AQI)技术规定(试行)》[22]，设定当轨迹对应的PM2.5浓度>150 μg/m3时，定义该轨迹为重污染轨迹，其值越大表示经过该网格污染轨迹的比例越高。取35°-60°N、100°-135°E的区域作为研究范围，将其划分为0.25°×0.25°的网格。经过网格(i,j)的重污染轨迹数为mij，经过网格(i,j)的总轨迹数为nij，则PSCF函数的表达式如式(1)所示。
由于PSCF是一种条件概率，为避免nij较小的网格出现非常高的PSCFij，通过引入经验权重函数W(nij)对其进行降误差处理[21]，表达式如式(2)所示。
式中，nave为研究区域内每个网格的平均轨迹端点数。经统计，本文中nave取5。
(1)
(2)
1.5 传输指数计算
气团的传输轨迹中包含了大量对重污染天气有指示意义的信息。张磊等[17]结合轨迹模式的输出结果和污染物排放源强度，设计了一种可评估污染物平流输送强度的客观参数，并证明其可较为合理地表征传输强度。花丛[18-19]等在此基础上引入潜在源贡献因子，结合传输轨迹和PM2.5排放源清单构建传输强度，计算公式为
Tl(i,j)=Rl(i,j)×E(i,j)×Wd(i,j)×Wtl(i,j)×PSCF
(3)
式中，T为传输强度，R为输送概率，E为PM2.5排放源强度，Wd为距离权重函数，Wtl为时间权重函数，PSCF为潜在源贡献因子。下标l和(i，j)是对应的轨迹和网格。其中R、Wd、Wtl的具体计算方法可参考花丛等[18]的研究。检验结果表明，传输强度与PM2.5浓度有较好的一致性，可以作为分析重污染天气成因的有效工具，在预报评估中发挥支撑作用。
该计算方法中，采用的是清华大学PM2.5排放源，E不随时间变化。考虑到哈尔滨地区冬半年多受秸秆燃烧等突发性大气污染事件影响，为更真实地反映由近地面PM2.5浓度实时变化及排放强度突变引起的输送强度变化，计算轨迹模拟时段内各监测站点的PM2.5平均浓度，作为网格点下垫面污染物强度。改进后的E(i,j)计算公式为
(4)
其中,EPM2.5(s)为监测站点s在轨迹模拟时段内的PM2.5平均浓度，sn为网格点(i,j)周边0.5°×0.5°范围内监测站点数。
由式(3)计算得出各网格点传输强度T1(i,j)后，求和即得到轨迹l的传输指数Tall：
(5)
其中，m、n分别代表经向和纬向上的网格数。
2 结果与分析
2.1 哈尔滨市大气污染特征
参照《HJ 633-2012环境空气质量指数(AQI)技术规定(试行)》[22]，同时结合哈尔滨市冬半年大气污染首要污染物为PM2.5的实际情况，本文定义0
表1 2014-2017年冬半年哈尔滨PM2.5浓度分级
统计2014-2017年冬半年每日间隔6 h(02时、08时、14时、20时)的PM2.5浓度逐月分布情况(图略)。在统计时段，各月PM2.5浓度中位数分别为40、74、86、82、59、48，第三个四分位数分别为87、141、135、141、98、76；其中11月、12月、1月空气质量状况较其他月份偏差，为重污染天气相对高发时段，近25%的时次PM2.5浓度≥150 μg/m3。
上述这些特征均表明哈尔滨冬半年清洁天气和轻度污染天气占主导地位，程度极重的污染天气较少。但是每个月均出现过PM2.5浓度>150 μg/m3的重污染天气。统计时段内共有160个时次PM2.5浓度超过500 μg/m3，占总数的0.9%，均出现在10-12月，其中10-11月占95.7%；共有9个时次PM2.5浓度超过1000 μg/m3，均出现在10月下旬至11月上旬。10-11月是哈尔滨周边农区秸秆集中焚烧期，受秸秆焚烧的影响，PM2.5浓度超过500 μg/m3时次多出现在这一时期。2017年11月1日16时哈尔滨PM2.5浓度曾达到1543 μg/m3的极高值。
另外，一天中污染程度最重的时次多出现在夜间至早晨时段，即20时、02时、08时，其中02时和08时最多，空气污染程度有日变化特点。
2.2 污染传输特征分析
2.2.1 后向轨迹聚类分析
利用Trajstat软件[10-13]模拟哈尔滨2014-2017年冬半年每日逐6 h气团后向轨迹，每条轨迹模拟时长为24 h，时间分辨率为1 h，共获得2916条气团后向轨迹。为体现近地层传输特征，同时减少地面摩擦影响，模拟起始高度选为500 m。由于本研究仅关注研究点(哈尔滨)的后向气流轨迹输送的水平来向，因此采用Trajstat软件中角距离聚类算法对冬半年到达哈尔滨的气流轨迹进行聚类分析，得到冬半年影响哈尔滨不同输送气流类型，在此基础上对每组气流所对应的污染物浓度特征进行统计分析。
对后向轨迹模式模拟的2916条轨迹进行聚类分析(图2)，得到4条主要传输路径。将模拟起始时刻哈尔滨市PM2.5浓度与轨迹进行对应，经过数据质量控制，共得到含有PM2.5浓度信息的轨迹2894条。其中，当PM2.5浓度>150 μg/m3时，定义该轨迹为重污染轨迹，共有460条，占总轨迹的15.9%。
图2 2014-2017年冬半年哈尔滨24 h后向轨迹分析及聚类结果
第1类轨迹(西南路径)共617条，占轨迹总数的21.3%，主要为来自内蒙古东部或辽宁中西部经吉林中西部地区进入到哈尔滨市的西南路径气流。其中，重污染轨迹为147条，占此类轨迹总数量的23.8%、占轨迹总数量的5.1%。
中堂画的边上还有大相框，和对联紧挨着，相框里有不同时期的全家福，也有老年人的单人照，新人的结婚照，婴儿的百日照。黑白泛黄的，彩色的，都在里面，相似的五官，相似的神情，只有服装的样式是不一样的，发型也不一样。
第2类轨迹(西北偏西路径)共1298条，占轨迹总数的44.8%，主要为来自内蒙古东北部经齐齐哈尔、绥化南部等地进入到哈尔滨市的西北偏西路径气流，高空为西北气流，近地面以偏西气流为主。其中，重污染轨迹为200条，占此类轨迹总数量的15.4%、占轨迹总数量的6.9%。
第3类轨迹(偏北路径)共821条，占轨迹总数的28.4%，主要来自黑龙江省大兴安岭、黑河经绥化等地进入到哈尔滨市的偏北路径气流，从高空到地面均为偏北气流。其中，重污染轨迹为71条，占此类轨迹总数量的8.6%、占轨迹总数量的2.5%。
第4类轨迹(偏东路径)共158条，占轨迹总数的5.5%。这类传输路径多出现在地面高压南侧或南来地面低压北侧的弱辐合区内，一般为气团沿小兴安岭西南侧经过哈尔滨北部遇大青山脉阻挡及在引导气流的影响下向西折移到哈尔滨市城区，所经区域在秋末冬初季节属于秸秆燃烧高排放量区域。其中，重污染轨迹为42条，占此类轨迹总数量的26.6%、占轨迹总数量的1.5%。此类轨迹为4类轨迹中数目最少、轨迹最短的。
第2类轨迹和第3类轨迹长度明显长于第1类和第4类的。在模拟时间相同的前提下，说明西南路径和偏东路径气团的移动速度明显小于西北偏西路径和偏北路径气团的移速，其中偏东路径气团(第4类轨迹)移动速度尤其缓慢，说明其主导时段的天气较为静稳，为污染物的积累创造了有利的气象条件。
西南路径气流和偏东路径气流的平均PM2.5浓度较大，分别为128.3 μg/m3和115.1 μg/m3。重污染轨迹在各自类别轨迹数中所占比例分别为23.8%和26.6%，重污染轨迹的平均PM2.5浓度分别为292 μg/m3和275 μg/m3。相比之下，西北偏西路径和偏北路径的平均PM2.5浓度分别为68.7 μg/m3和87.4 μg/m3，重污染轨迹在各自类别轨迹中出现概率分别为8.6%和15.4%，说明来自这两个方向的气团较为清洁。
2.2.2 不同传输轨迹与主要天气系统之间的关系
西南路径(第1类轨迹)多出现在弱低压系统前和高压西侧或西南侧。这种传输路径，有利于暖湿空气和吉林省中西部地区的污染物向哈尔滨地区输送，有利于污染物吸湿增长及近地面逆温的形成和加强；并且在这些区域容易出现弱的气旋性弱辐合，有利于污染物积聚；西南路径气流还可以使污染物在小兴安岭地形的阻挡下积聚在松嫩平原。所以西南路径传输轨迹容易导致哈尔滨冬半年出现重污染天气。
以2017年10月18-20日重污染过程为例。2017年10月18日16时至20日15时，哈尔滨连续48 h出现重度及以上级别污染，严重污染达到37 h，首要颗粒物为PM2.5，平均浓度为438 μg/m3。 19日00-01时污染程度最重，PM2.5平均浓度达到837 μg/m3，其中平房东轻厂监测站PM2.5峰值浓度高达1487 μg/m3。
在500 hPa平均高度场和距平场上(图3)，东亚中高纬地区以平直的偏西气流为主，范围较小的冷涡位于黑龙江以北的俄罗斯大陆上，哈尔滨位于槽底部正负距平交界处。在近地面1000 hPa平均高度场和距平场上(图4)，黑龙江省位于弱的低压区内，等压线稀疏，哈尔滨位于低压前西南气流里，风力小，逆温强，不利于污染物扩散。18日夜间到19日上午及19日夜间到20日上午，哈尔滨分别在两个低压系统前的西南气流里，空气污染十分严重。19日下午在前一个弱低压后弱偏北气流的控制下，空气质量短时间好转。单站地面三线图和风的变化及地面水平风场特征显示(图略)，在污染较重阶段对应着随时间或空间有西南风向偏东风、再向偏北风气旋式辐合变化的过程，有利于污染物积聚。
图3 2017年10月18-20日500 hPa平均高度场和距平场
图4 2017年10月18-20日1000 hPa平均高度场和距平场
以2018年2月2-6日为例分析西北偏西路径(第2类轨迹)与天气系统的关系。在500 hPa平均高度场和距平场上(图5)，乌拉尔山高压脊强度偏强、范围偏大，黑龙江省位于高压脊前的西北气流里，为正高度距平。在近地面1000 hPa平均高度场和距平场上(图6) ，哈尔滨位于强大的蒙古高压前干冷的西北气流里，西北气流向哈尔滨输送清洁的干冷空气，空气扩散条件和空气质量较好。
图5 2018年2月2-6日500 hPa平均高度场和距平场
图6 2018年2月2-6日1000 hPa平均高度场和距平场
以2016年2月13-14日为例分析偏北路径(第3类轨迹)与天气系统的关系。在500 hPa平均高度场和距平场上(图略)，2016年2月中下旬东亚大槽位置偏西，日本及以东的太平洋海域原东亚大槽平均位置被高压脊取代，乌拉尔山高压脊强度偏强、范围偏大。中国北方在乌拉尔山高压脊前东亚大槽的影响下，冷空气势力偏强，黑龙江省为较强的西北气流。其中2月13-14日，在近地面1000 hPa平均高度场和距平场上(图略)，哈尔滨位于较强的蒙古高压前和较强日本海低压之间干冷的偏北气流里，等压线密集，风力较大，近地面无逆温，空气扩散条件和空气质量较好。
偏东路径(第4类轨迹)多出现在高压南侧或低压北侧低压倒槽区，风速常较小。这种传输路径，有利于暖湿空气和哈尔滨阿城、宾县等地污染物向哈尔滨市区输送，有利于污染物吸湿增长及近地面逆温的形成和加强。在这些区域，还容易出现弱的气旋性弱辐合，有利于污染物积聚。偏东传输轨迹路径短，表明大气中、低层无明显影响系统活动或影响系统弱，天气形势静稳度高，气团移动速度缓慢。该路径传输的空气经过较高的排放强度地区，与哈尔滨城区周边近地面PM2.5浓度较高的污染物充分混合，使污染物浓度增大。因此偏东路径传输轨迹容易导致哈尔滨市区冬半年出现重污染天气。
以2016年11月3-5日为例分析偏东路径(第4类轨迹)与天气系统的关系。在500 hPa平均高度场和距平场上(图略)，东亚大槽位置偏西偏强，槽的经向度大，冷涡中心位于雅库茨克附近，黑龙江省位于冷涡西南侧西北气流里，冷空气势力偏强。在近地面1000 hPa高度场和距平场上(图略)，黑龙江省位于高压底部的偏东气流里，逐6 h天气形势变化显示，哈尔滨先后位于弱高压南部、鞍形场及较强高压南部，均为偏东暖湿气流控制，等压线稀疏，风力小，再加之周边农民集中焚烧秸秆，致使哈尔滨城区出现严重的空气污染。当哈尔滨位于鞍形场靠近蒙古低压东北侧的偏东气流里时，污染浓度最高，4日16时至23时，连续8 h PM2.5浓度≥1000 μg/m3。
2.2.3 污染气团潜在源区分析
PSCF函数表征经过某一区域的气团到达观测点时对应的要素值超过设定阈值的条件概率。高PSCF值所对应的区域即污染物的潜在源区，经过该区域的轨迹为污染物传输的主要路径。利用TrajStat软件对哈尔滨PM2.5进行污染气流潜在源区分析。计算结果(图7)显示，影响哈尔滨城区PM2.5浓度的主要潜在贡献源区为松嫩平原和辽河平原一带。其中，PSCF函数值在0.3以上的高影响区域出现在黑龙江哈尔滨地区、吉林中部、辽宁北部及内蒙古通辽等地。这些较强的潜在贡献源区也是污染物西南传输路径或偏东传输路径的必经之地。上述区域为我国东北主要农区，秋末冬初秸秆焚烧多，且冬季供暖燃煤多，污染物排放强度大，气团携带的污染物浓度相对较高。
2.3 污染物传输指数与PM2.5浓度之间的关系
参考花丛[18-19]等传输指数计算方法，计算2014-2017年冬半年哈尔滨逐6 h传输指数，结果最大值为6.9，最小值为0.03。对6 h间隔传输指数与相应时次PM2.5浓度做相关性分析，其相关系数为0.5，通过α=0.01的显著性检验。传输指数越大，对应PM2.5浓度>150 μg/m3的概率越大，传输指数越小，对应空气污染的概率越小。分析发现，当传输指数≥2.5时， 70% 时次的PM2.5浓度≥ 200 μg/m3，80% 时次的PM2.5浓度≥ 150 μg/m3；当传输指数≤0.2时，仅有5%时次的PM2.5浓度≥ 100 μg/m3，1%时次的PM2.5≥ 150 μg/m3，无PM2.5浓度≥ 200 μg/m3的情况发生。这些都说明传输指数可以在一定程度上描述哈尔滨城区污染物浓度变化，解释污染成因。但是传输作用仅是污染天气的形成条件之一，所以传输指数还不能完全描述污染物的浓度变化，污染物浓度还受污染源排放量及扩散条件等影响。
图7 哈尔滨潜在源区贡献因子(PSCF)分布
花丛等[18]对京津冀地区传输指数的研究分析发现，平流输送对污染物的传输导致了PM2.5浓度的变化，PM2.5浓度相对传输条件的变化有一定滞后响应时间。通过分析发现，哈尔滨PM2.5浓度变化和传输指数的变化之间也存在一定滞后性。为计算滞后响应时间的长度，对2014-2017年冬半年每日逐6 h传输指数和PM2.5浓度进行时滞相关分析，结果显示，在PM2.5滞后传输指数6 h的情况下，两者相关性也较为明显，相关系数达到0.44，通过α=0.01的显著性检验。说明传输指数对污染物浓度变化有大约6 h的预报提前量，可以此作为参照，估算污染物浓度的变化趋势。
2.4 传输强度特征
将所计算时段内各网格点上的传输强度进行叠加，即可得到传输强度的空间分布状况，即不同区域对目标城市的综合传输影响。对于周边区域而言，黑龙江省西南部、吉林中部、辽宁北部传输强度大于0.1(见图8)，说明上述地区的污染物可通过区域传输对哈尔滨产生较明显影响。对于哈尔滨本地来说，哈尔滨市城郊及其双城、呼兰、阿城、五常西部等县(市)距离哈尔滨市区近，冬半年供暖燃煤多，秸秆焚烧多，排放强度较高，传输强度超过10甚至超过100，传输作用最为明显。距离哈尔滨城区越近， 传输强度越大。哈尔滨市城区污染源主要来自本地和附近县市排放，在对污染源进行调控和减排时，可参考传输强度的空间分布进行决策，尽量以小的经济成本实现最好的区域联防联控效果。
图8 2014-2017年冬半年哈尔滨传输强度分布
3 结 论
(1) 哈尔滨冬半年清洁天气和轻度污染天气占主导地位，PM2.5浓度>150 μg/m3的重污染时次和重污染轨迹仅占15.1%和15.9%,但是每月都有重污染时次出现，其中11月、12月、1月空气质量比其他月份偏差，有近25%的时次PM2.5浓度接近150 μg/m3。一天中污染程度最重的时次多出现在夜间至凌晨时段。
(2) 后向轨迹聚类分析结果表明，影响哈尔滨冬半年空气质量的传输路径有4条。偏东路径和西南路径为主要污染路径，偏北路径和西北偏西路径为相对清洁路径。重污染轨迹数在其所属类别传输路径轨迹数中占比由高到低依次为：偏东路径(26.6%)、西南路径(23.8%)、西北偏西路径(14.5%)、偏北路径(8.6%)。
(3) 西南路径气流多出现在弱低压系统前和高压西侧或西南侧，偏东路径气流多出现在高压南侧或低压北侧的低压倒槽区。这两条路径的气流有利于污染物和暖湿空气向哈尔滨输送及污染物的积聚，容易导致空气重污染。偏北路径和西北偏西路径，多出现在较强的蒙古高压前和低压后西北或偏北气流里，有利于向哈尔滨输送干冷空气和清洁空气，有利于污染物的沉降和扩散。
(4) 松嫩平原和辽河平原是影响哈尔滨的主要污染物潜在源区。哈尔滨市城郊及双城、呼兰、阿城、五常西部等地传输强度超过10,甚至超过100，传输作用最为明显。哈尔滨冬半年逐6 h间隔的传输指数与PM2.5浓度存在一定的正相关关系。传输指数对PM2.5浓度变化有大约6 h的预报提前量，可以此作为参照，估算污染物浓度的变化趋势。
