1 引言(Introduction)

近年来, 随着航运事业的高速发展, 中国已经成为世界上首屈一指的海事国家(Wee et al., 2018).然而, 在港口和航运业繁荣发展的同时, 船舶航运也带来了大量的大气污染物排放.研究表明, 2014年环渤海地区、长江三角洲地区(长三角)和珠江三角洲地区(珠三角)的船舶氮氧化物(NO_x)排放分别占NO_x总排放量的约13%、28%和17%(Chen et al., 2017), 船舶已成为我国沿海地区重要的排放源.为了控制船舶排放, 2018年底中国交通运输部发布了更为严格的《船舶大气污染物排放控制区实施方案》, 旨在设立船舶大气污染物排放控制区, 降低船舶硫氧化物、氮氧化物、颗粒物和挥发性有机物等大气污染物的排放, 这也表明了我国政府管理部门对加强船舶排放控制与精细化管理的决心.而船舶大气污染物排放控制区的规划和精细化船舶管理方案的发布需要依赖船舶排放清单及时空分布特征研究.

对于船舶排放空间表征, 传统基于港口或航道通航能力的船舶排放空间表征, 难以准确表征船舶排放空间分布；使用美国船只自动化互助救助系统(Automated Mutual-Assistance Vessel Rescue System, AMVER)或国际综合海洋大气数据集(International Comprehensive Ocean Atmosphere Data Set, ICOADS)等数据库作为船舶排放空间表征(Endresen et al., 2003；Wang et al., 2007), 由于船舶区域分布差异大, 对于亚洲地区并不适用.船舶自动识别系统(Automatic Identification System, AIS)具有覆盖区域广、船舶数据量大等优势, 能发射并接收船舶轨迹点经纬度、航行速度与航行时间等动态信息, 可更好地表征船舶排放空间分布.目前, 采用动力法结合AIS动态信息及船舶登记系统存储的静态参数计算船舶排放, 能反映船舶的瞬时排放量化和时空动态变化特征, 已经成为船舶排放清单研究的主流方法(刘欢等, 2017).Liu等(2016)、朱倩茹等(2017)和尹佩玲等(2017)基于AIS数据与动力法分别建立高分辨率的东亚地区远洋船舶排放清单、珠三角船舶排放清单和宁波-舟山港船舶排放清单.相较于传统统计方法, 基于AIS数据的动力法在时空特征表征方面更为精细、时效性更高.然而, 受AIS工作模式的制约, 由AIS网络通信问题导致的轨迹点漂移与缺失(吴建华等, 2017), 在一定程度上会造成船舶排放量低估和船舶空间分布特征表征错误, 进而影响船舶排放的空气质量成效评估和船舶排放控制区的划分.此外, 船舶静态参数的缺失也会大大影响船舶排放的计算.因此, 如何修复AIS轨迹、补充船舶静态参数对于精细化船舶管理尤为重要.

目前, 对于AIS轨迹修复的研究主要采用基于航行信息插值的轨迹修复方法和结合数值插值的轨迹修复方法.基于航行信息插值的轨迹修复方法主要是考虑航速航向的AIS航迹插值方法(王超等, 2015), 能较好地把控船舶航向角度变化, 但对航行信息的依赖程度相对较大.结合数值插值的轨迹修复方法主要包括传统线性插值法(王晨熙等, 2009)、三次样条插值法(刘立群等, 2015)和多维阵列插值法(Wang et al., 2017).其中, 传统线性插值法利用轨迹点经纬度与航行时间进行线性插值, 但对于非线性轨迹修复误差大；多维阵列插值法利用多种船舶信息构建船舶时空共生模式分析补充船舶轨迹, 对航行信息依赖程度大且较为复杂；相反, 三次样条插值对航行信息的依赖程度较小, 并且对非线性轨迹也能较好地修复.然而, 目前少有研究将AIS轨迹修复技术和动力法相结合, 应用于船舶排放清单表征改进.随着大数据的发展和船舶排放管理的精细化, 有必要探索AIS轨迹修复技术对船舶排放表征的改进效果.

因此, 本研究以2013年广东省AIS船舶数据为例, 采用基于时间和经纬度的三次样条轨迹修复方法对AIS船舶轨迹进行轨迹修复, 并结合动力法计算船舶排放, 分析对比AIS轨迹修复前后船舶排放表征的差异, 并利用空气质量模型CMAQ和卫星观测评估AIS轨迹修复对船舶排放表征和广东沿海空气质量模拟的改进效果.

2 研究方法与资料(Research methods and data) 2.1 基于时间和经纬度的三次样条轨迹修复方法

三次样条插值方法是通过求解三弯矩方程组得出曲线函数组的方法, 最终得到通过一系列形值点的光滑曲线.由于AIS船舶轨迹灵活多变且数据量大, 因而导致部分AIS轨迹航向角等数据常有缺失.基于时间和经纬度的三次样条插值方法, 在运用数值插值方法进行插值的过程中考虑了船舶航行情况, 既能保持数值插值的灵活性与低计算成本, 又能反映船舶航行的实际状况, 是较好的船舶轨迹修复方法.

本研究会在轨迹修复前后对船舶轨迹点进行清洗, 将分布异常、速度异常等轨迹点剔除后再进行船舶排放计算.由于在船舶排放计算过程中不同航行工况需采用不同的排放因子, 不同航行工况的排放差异较大, 因此, 需要区分航行工况进行研究.基于2013年广东省AIS数据分析发现, AIS上报时间间隔也会随着航行工况改变而改变.为更好地捕捉不同航行工况下的船舶轨迹变化情况来进行轨迹修复, 本研究根据船舶航行速度区分不同航行工况, 包括正常航行、慢速航行、机动航行及靠泊, 分别取各航行工况下AIS平均上报时间间隔作为进行轨迹修复的判断依据和轨迹修复插值时间间隔, 工况划分限值与平均上报时间间隔如表 1所示.

基于时间和经纬度的三次样条轨迹修复方法, 取[t_i-2, t_i-1, t_i, t_i+1]4个船舶轨迹点位置及速度信息, 对于[t_i-1, t_i]时间段分别作时间-经度、时间-纬度和时间-速度的三次样条插值.

首先建立三次样条插值模型, 设i点的时间t_i、经度x_i、纬度y_i, 三次样条函数S(t)在节点t_i处的二阶导数为M_i, 可得到三次样条表达式为：

(1)

式中, M_i为未知参数；h_i=t_i+1-t_i.通过对S(t)进行求导得出S′(t), 利用S′t_i+0=S′t_i-0可得：

(2)

(3)

(4)

(5)

由于方程式比未知数少, 因此, 根据第一类边界条件补充i=1和i=n的端点方程, 最后可得到关于参数M₀, M₁, …, M_n的n+1阶线性方程组, 其三弯矩阵方程为：

(6)

根据上式求出M₀, M₁, …, M_n的值, 即可求出时间-经度、时间-纬度的三次样条插值, 通过合并上述两条三次样条曲线得到船舶的航行轨迹.对于三次样条轨迹修复方法的误差验证参照文献(刘立群等, 2015).

2.2 船舶排放计算及静态参数回归补充

本文主要针对2013年广东省海域已识别的AIS船舶进行排放计算与研究, 包括远洋船舶273艘、沿海船舶215艘, 共计488艘, 船舶类型包括集装箱船、普通货船、油轮、客船、拖船及其他船等14种, 结合船舶静态信息数据与AIS轨迹数据计算船舶排放, 每个船舶轨迹点污染物排放量计算模型如下：

(7)

(8)

(9)

式中,m =1、2、3、4，代表靠泊、机动、慢速、正常航行4种航行工况；p=1、2、3, 代表主机、辅机和锅炉3种动力设备；EM为污染物排放量(t)；MCR为发动机额定功率(kW)；LF为负载率；T_i和T_i-1分别为第i条和第i-1条AIS动态信息上报时间(h)；EF为排放因子(g·kW^-1·h^-1)；EF₀为基础排放因子(g·kW^-1·h^-1)；FCF为燃料校正因子；LLA为低负载校正因子(当船舶负载率低于20%时燃油燃烧效率会下降, 排放量会相应上升, 此时引入低负载校正因子对排放因子进行细化修订).

利用公式(7)计算船舶排放需要完备的船舶静态信息数据库作参数支撑.结合劳氏数据库与船舶登记系统数据建立船舶静态信息数据库, 但船舶静态信息数据库中的部分参数仍不完整, 部分船舶档案常缺失船舶主辅机功率、最大航速等重要参数.因此, 本研究利用不同类型船舶载重吨位拟合船舶主机功率(表 2), 通过分析不同类型船舶辅机与主机比例来计算补充辅机功率(表 3)；此外, 对于缺失的最大航速则采用该类型船舶最大航速推荐值(表 4), 逐步补充完善船舶静态信息.

2.3 空气质量模拟平台搭建与排放清单

为验证AIS轨迹修复对船舶排放表征的改进作用, 本研究搭建以WRF-SMOKE-CMAQ为核心的空气质量模拟平台, 评估AIS轨迹修复对船舶排放空间表征和海域空气质量的影响.其中, CMAQ是在我国应用广泛且成熟的空气质量模式, WRF为CMAQ提供气象模拟数据, SMOKE对排放清单做时间分配、空间分配、化学物种分配和烟羽抬升等, 最终输出供CMAQ模式使用的网格化、小时化、物种化排放源清单.本研究采用的模拟时段为2014年7月, 由于7月广东省主要以南风为主, 因此, 船舶排放对广东沿海的影响更为突出.模型模拟区域以Lambert投影坐标系, 中心经度为114°E, 中心纬度为28.5°N, 两条平行标准纬度分别为15°N和40°N.设置两层嵌套模拟网格, 分辨率分别设置为27 km×27 km(D1)和9 km×9 km(D2), 其中, 第一层网格(D1)覆盖中国大部分地区和东亚及东南亚部分国家, 第二层网格(D2)覆盖华南大部分地区, 包括广东省及周边相邻部分省份与南海部分海域, 具体如图 1所示.模拟区域垂直方向共设置18个气压层, 层间距自下而上逐渐增大.化学机制选用改进的CB05机理(CB05CL-AE5-AQ)(Byun et al., 2006), 该机理同时考虑了液相和气溶胶化学, 采用EBI(Euler Backward Iterative)算法求解.WRF模型所使用的物理过程参数化方案见文献(王延龙等, 2018).

本文模型模拟采用2014年广东省人为排放清单(范小莉等, 2017；卞雅慧等, 2018；王肖丽等, 2019)作为输入.中国区域采用清华大学开发的2012MEIC清单(Center for Earth System Science, Tsinghua University. MEIC Model, http://www.meicmodel.org/), 其他区域采用亚洲REAS清单(Zhang et al., 2009), 天然源排放采用天然源排放估算模型MEGAN v2.10(Guenther et al., 2012)的计算结果.由于AIS数据包含的船舶类型有限, 只包括沿海与远洋部分类型船舶, 未能涵盖广东省海域所有船舶.因此, 未涵盖的船舶类型采用Li等(2018)的船舶排放清单计算方法进行补充, 得到2014年广东省海域船舶排放清单总量(表 5).本研究基于动力法计算已识别的AIS船舶排放量, 将船舶轨迹点排放量叠加到网格, 每个网格船舶排放量占比作为船舶排放空间分配因子, 应用于SMOKE排放清单模型当中, 得到空间网格化清单用于下一步模型模拟.

本文利用广东省沿海区域NO₂小时监测浓度验证模拟平台的可靠性.根据广东沿海区域分布, 将沿海区域对出海域及珠江水系划分为粤东海域、珠江口海域及粤西海域(图 1).图 1同时给出了各站点NO₂模拟值与监测值的标准化平均偏差(Normalized Mean Bias, NMB)分布, 总体标准差为-11%, 能较好地反映沿海站点NO₂浓度变化.通常将模型模拟浓度与监测浓度的偏差在-30%~30%内的情况视为模拟效果良好(王延龙等, 2018).由站点NMB分布及研究区域划分可以看出(图 1), 除珠江口部分站点高估外, 其他区域尤其是沿海地区的整体模拟偏差在可接受范围内, 表明模拟平台具备较好的稳定性和可靠性.

本文继续采用OMI卫星遥感反演的分辨率为0.25°的对流层NO₂垂直柱浓度(Vertical Column Density, VCD)(https://giovanni.gsfc.nasa.gov/giovanni/)与模拟浓度进行对比, 来验证AIS轨迹修复对广东沿海地区空气质量模拟的改进作用.由图 2可知, 观测反演区域与模拟D1网格范围一致, 垂直层高度及反演时段与模拟高度及时段一致.

3 分析与讨论(Results analysis and discussion) 3.1 AIS轨迹修复的排放验证

为验证轨迹修复对船舶排放计算的影响, 本研究抽取3段具有不同航行状态的完整轨迹, 通过删除部分轨迹点以模拟4种常见的AIS轨迹缺失情景, 然后结合动力法计算对比原始完整轨迹、模拟缺失轨迹和修复后轨迹的NO_x排放.由于船舶航行速度和航行工况都会影响船舶排放计算, 因此, 本研究抽取的3段完整船舶轨迹具有不同的航行工况.其中, 轨迹1为船舶维持在正常航行工况(正常轨迹), 轨迹2为船舶由慢速航行加速转换为正常航行工况(加速轨迹), 轨迹3为船舶由正常航行减速转换为慢速航行工况(减速轨迹).此外, 为对比不同轨迹缺失情景对船舶排放计算的差异, 本研究模拟4种轨迹缺失情景, 包括中间缺失100个轨迹点的中段缺失情景、在轨迹两端各缺失50个轨迹点的首尾缺失情景、每10个轨迹点缺10个的连续多段缺情景及每3个轨迹点缺3个的零散缺失情景.这3段轨迹的12种模拟缺失情景见图 3, 其NO_x排放变化见表 6.

在4种不同轨迹缺失情景中, 中段缺失情景的船舶排放误差最为明显, 特别是对于船舶工况变化的轨迹.对于加速轨迹缺失段, 由于动力法默认采用缺失段时间间隔和时间间隔后一轨迹点的速度值进行排放量计算, 忽略了缺失轨迹的加速过程, 容易导致船舶排放高估.相反, 对于减速轨迹缺失段, 会导致船舶排放低估.在本案例中, 加速轨迹在中段缺失的情景下, 排放计算偏差可达14%左右, 其中, 缺失段排放量比修复后排放量高44%；在减速轨迹中段缺失的情景下, 造成了18%的船舶排放低估, 缺失段排放量比修复后排放量低42%.因此, 在这种轨迹缺失情景下进行轨迹修复, 能降低10%以上的船舶排放计算误差, 同时弥补空间上的排放缺失.在首尾缺失情景中, 由于每段缺失的点数相对较少且不涉及船舶航行工况的转变, 偏差相对较小.但一旦涉及船舶航行工况的转变, 如减速轨迹首尾缺失情景图 3i, 排放误差增大, 可造成7%以上的船舶排放高估.其余轨迹缺失情景中, 随着每段缺失轨迹的点数减少, 排放误差也在减小.因此, 轨迹修复通过补充船舶速度变化过程, 对于缺失轨迹点数越多的船舶轨迹排放计算改善程度越大, 尤其是涉及船舶航行工况改变的情况.

AIS轨迹修复能减少船舶排放量错误分配, 提高船舶排放的空间表征.基于动力法计算船舶轨排放量时, 会默认将两个轨迹点之间的排放量分配到后一轨迹点上.这就导致在轨迹缺失情景下, 会将缺失段整段的排放量都分配到缺失段时间间隔后一轨迹点, 即缺失段后端点, 由于缺失时间长与速度变化大等原因, 将使得该轨迹点排放量较修复后轨迹点排放量高估, 甚至高出数百倍.如在加速轨迹中段缺失情景, 缺失段后端点排放量比修复后轨迹点排放量高197倍, 而减速轨迹中段缺失情景, 缺失段后端点排放量也比修复后轨迹点排放量高101倍.

3.2 基于轨迹修复的广东省船舶排放表征改进

本研究将AIS轨迹修复方法应用于广东省船舶排放表征改进.根据2.3节中的广东省各研究海域划分, 表 7统计了广东省各研究海域船舶轨迹点数, 轨迹修复点数差异的空间分布见图 4a.由表 7可见, 基于时间和经纬度的三次样条轨迹修复方法能够大幅度提高广东省AIS数据的轨迹点数量, 轨迹点总数由4685773个增至5746664个.特别是在粤东海域, 修复后轨迹点增加高达88%, 粤西海域也有明显的增加, 约36%.在粤东海域, 轨迹点缺失集中在揭阳和汕尾对出海域(图 4a), 可能与该地区陆地基站覆盖范围有关.图 4c中框出的粤东重点修复海域轨迹修复后, 船舶轨迹点数急增20倍左右.粤西海域在距离大陆越远的海域轨迹缺失越严重, 其原因可能更多的与船舶AIS信号随距离减弱或陆地基站覆盖范围有关.80%以上船舶轨迹点集中于珠江口海域, 船舶离岸距离相对较近, AIS接收信号较好, 轨迹点缺失的比例最低.因此, 基于原始AIS数据进行船舶排放空间分布, 会造成粤东与粤西海域的排放缺失, 轨迹修复补充了缺失轨迹的船舶空间分布, 使船舶分布更连贯.

本研究继续利用动力法计算对比轨迹修复前后的广东省海域船舶NO_x排放, 结果见表 8.对比图 4d和图 4e可知, 总体上, 修复后的船舶排放在空间分布上更为密集、连贯, 主要表现在粤东揭阳与汕尾对出海域及粤西阳江与茂名对出海域.2014年广东省海域有72%的船舶排放集中于珠江口, 尽管修复后轨迹点数增加近60万个, 但珠三角区域NO_x排放量仅降低0.9%.这可能是因为在珠江口海域船舶航速较低, 同时轨迹缺失主要以零散缺失为主, 轨迹修复细化了珠江口的船舶排放空间分布.类似的, 粤西海域由于轨迹修复细化补充了船舶点航行速度变化过程, 船舶排放量稍有减小.但在轨迹点增加最为明显的粤东海域, 修复后船舶NO_x排放增加了6%, 粤东重点修复海域船舶NO_x排放增加了2.7倍, 轨迹修复大大补充了该区域的船舶排放(图 4c).总体而言, 轨迹修复对区域船舶NO_x总排放量的影响较小, 但对局部区域船舶排放与空间表征有明显的改进.

3.3 轨迹修复对模拟船舶排放污染影响

为将上述研究结果应用于空气质量模型模拟, 本文研究了轨迹修复前后船舶排放空间分配因子的差异, 结果如图 5a所示.由3.1节船舶排放分析可知, 在轨迹缺失情况下, 基于动力法计算船舶排放会把整段缺失轨迹的排放量分配到缺失段后端点, 导致错误分配.结合船舶在广东海域双向航行, 若按原始轨迹进行空间分配将造成图 5a中框出的粤东重点修复海域网格排放量低估73%, 而该海域两端排放量高估26%.

将轨迹修复前后船舶排放空间分配因子应用到空气质量模拟平台, 并利用模型模拟与卫星观测对比的手段, 验证轨迹修复对沿海区域空气质量的影响, 图 5c~5e展示了轨迹修复前后船舶排放NO₂模拟浓度及变化分布.其中, 卫星观测使用OMI反演的NO₂垂直柱浓度数据, 模拟结果月均柱浓度与卫星观测月均柱浓度的相关系数r约为0.79, 相关性较强.为分析NO₂模拟柱浓度与卫星反演NO₂柱浓度的差值情况, 绘制了原始轨迹NO₂模拟柱浓度与卫星反演NO₂柱浓度差值的空间分布(图 5b), 并分别选定各海域模拟浓度变化集中的网格做箱线图(图 6)进行分析对比.由图 5b与图 6可知, 在广东海域, 总体上NO₂模拟浓度低估约29%, 其可能与船舶排放低估有关, 在珠江口部分区域高估与珠江口附近城市排放集中也有一定关系.对比轨迹修复前后NO₂模拟浓度空间分布(图 5c和图 5d)可见, 原始轨迹模型模拟的NO₂浓度分布在粤东揭阳和汕尾对出海域出现明显的中断, 造成该处船舶排放对沿海区域的影响被低估；轨迹修复后补充了该处船舶排放及对沿海区域的影响, 使得NO₂模拟浓度的空间分布更加连贯.同时, 轨迹修复也在一定程度上提高了区域NO₂模拟浓度, 粤东揭阳和汕尾对出海域NO₂模拟浓度增加20%以上, 在粤东重点修复海域NO₂模拟浓度的增幅达40%.由此可见, 轨迹修复能有效弥补轨迹缺失海域的船舶排放, 改进船舶排放模拟效果.从图 5b中框出的各海域模拟浓度变化集中的网格与图 6可知, 在轨迹修复后, 粤东海域与粤西海域NO₂模拟柱浓度均有一定程度增大, 与卫星观测浓度的差距由47%降至40%, 在珠江口海域NO₂模拟柱浓度稍有降低, 与卫星观测的差距由22%降至21%, 尤其在粤东重点修复区域NO₂模拟柱浓度与观测浓度的差距由51%减至6%, 差距减少45%, 总体上模拟结果更接近卫星观测结果, 说明轨迹修复能有效提高模拟准确度.

3.4 讨论

本研究采用基于时间与经纬度的三次样条轨迹修复方法, 总体上能满足轨迹修复的精度要求, 但对于非线性轨迹的修复具有一定误差.此外, 三次样条插值函数的建立需要一定量的轨迹历史数据, 并且受到插值轨迹参照点的影响较大.由此看出, 本研究采用的轨迹修复方法仍有一定的改善空间.在AIS数据质量较好、航向角等信息较为齐全的情况下, 增加航向角在轨迹修复中的考虑是较为有效的改进方向.

结合轨迹修复研究船舶排放清单是目前船舶排放清单空间表征研究的重点.张帆等(2019)采用考虑船舶航速航向的AIS航迹插值方法对船舶轨迹进行插值修复, 建立了动态船舶排放清单.而本研究基于时间与经纬度的三次样条轨迹修复方法, 虽然需要一定的轨迹历史数据来提供插值轨迹参照点, 但对于船舶航行信息的依赖程度相对较小, 在AIS数据质量较差的情况下也能适用.此外, 本研究结合模型模拟与卫星观测体现轨迹修复对船舶排放清单模拟的改善.在大数据发展趋势下, 下一步研究工作应该探讨如何将轨迹修复应用于动态船舶排放清单, 寻找更好的动态船舶排放清单建立方法, 提高动态船舶排放清单的准确度.

4 结论(Conclusions)

1) 本研究通过分析轨迹修复的排放验证得出, 轨迹修复对于轨迹缺失点数越多的轨迹排放计算的改善程度越大, 尤其是涉及船舶航行工况改变的情况, 能降低约10%以上的船舶排放计算误差, 在空间上能补充船舶排放空间分布和减少船舶排放量错误分配.

2) 轨迹修复后广东省船舶轨迹点数增加约100万, 船舶NO_x排放增加了0.6%.粤东海域是轨迹点与排放缺失最集中的海域, 轨迹修复后轨迹点增幅为88%, NO_x排放在广东省船舶排放的占比提升至22%, 其中, 粤东重点修复海域的NO_x排放量增加2.7倍.原始轨迹在广东省海域较为稀疏, 在粤东海域出现明显轨迹缺失；轨迹修复后广东省海域船舶轨迹更为密集, 粤东海域船舶轨迹得以补充, 船舶排放空间分布更连贯.

3) 与卫星观测对比, 广东省海域船舶排放模拟结果以低估为主, 在珠江口部分区域出现高估.轨迹修复后NO₂模拟空间分布更连贯, 粤东、粤西海域NO₂模拟浓度均有一定程度的提高, 珠江口海域NO₂模拟浓度出现一定程度的降低, 轨迹修复后粤东重点修复海域船舶排放模拟浓度与观测结果的差距减少45%, 模拟结果总体上更接近卫星反演, 从而提高了船舶排放的模拟准确度.