张勤耘，金双根

(1.上海大学 通信与信息工程学院，上海 200444；2.中国科学院 上海天文台，上海 200030)

0 引言

火山喷发时，会喷射出大量的岩石块、火山灰、火山气体等物质。这些物质会严重影响建筑物的完整性，破坏陆地或空中交通安全，扼制植物的生长，影响人类的健康生活和安全。因此监测并预报火山灰十分重要。世界各地的火山活动可通过各种各样的地基或空基设备进行监测，但这2类设备各有优势与缺点：地基设备可以提供相对持续的高精度时空覆盖度，但是密度不足够大、覆盖度有限；相反地，卫星系统可以确保测量的全球空间密度以及所研究区域的概要图，但是空间覆盖率较低，精度不高。为了更好地降低火山烟流的危害，相关的监测预报系统和设备如意大利地震台网中心(Istituto Nazionale di Geofisica e Vulcanologia，INGV)、激光雷达(light detection and ranging，LiDAR)、多普勒雷达(doppler radar，DR)等已经逐渐得以完善[1-3]。除此之外，火山烟流也可能影响全球卫星导航系统(global navigation satellite system，GNSS)的卫星信号。

文献[4]利用无电离层延迟测量(ionosphere-free linear combination，LC)的异常，采用地震层析成像算法将折射率异常和推算的火山灰的温度映射在时空图中，分析了2000-08的日本三宅岛火山(Miyake Jima)相关情况，表明热火山灰造成了全球定位系统(global position system，GPS)信号路径的延迟[4]。文献[5]采用稳健型统计学方法，以2007-09-04和2007-09-05的意大利埃特纳火山(Etna)为例，探索大气中火山烟流是否会影响GPS非差拟合相位残差(post-fit phase residuals，PPR)，研究结果表明，火山烟流明显影响了GPS PPR值，且是日本三宅岛火山所得残差值的1/10，猜测残差值与火山爆发强度存在一定的相关性[5]。文献[6]以运动学定位方法为基础，以2009年阿拉斯加州的里道特火山(Redoubt)为例，获取全球卫星导航系统(global navigation satellite system，GNSS)数据的相位残差，将大规模、短周期的定位偏移与穿过稠密的火山烟流的“卫星-测站对”的相位延迟联系起来，并提供了一种测量火山烟流扩散方向的方法即沿着卫星在天空图中的运行轨迹表示出相位残差，由此猜想利用卫星与地表测站之间的GPS信号的相位残差可以监测火山，同时也指出由于火山烟流的浓度可能不足够高而不能对卫星信号造成显著影响，即这种方法并不适用于所有的火山烟流[6]。上述方法都是使用标准的地学定位处理技术，并将火山灰信号作为一种非典型的大气效应，需要分析GPS载波相位观测范围，还要评估站心坐标、轨道、时钟、大气延迟以及相位模糊度等影响。而且火山灰与GPS接收机随时间变化，“卫星-测站对”的路径并不一定就通过火山灰；需要利用最小二乘估计法解算，其他模型误差都会影响待估参数。文献[7]提出利用GPS信号信噪比数据(signal-to-noise ratio，SNR)估计2008年阿拉斯加州奥克莫克火山(Okmok)和2009年里道特火山(Redoubt)，研究表明在火山爆发后的一小段时间内，SNR会在有一定的突降之后恢复到之前的变化趋势上[7]。文献[8]使用了非差拟合相位残差方法PPR和SNR探测方法监测2012-08-06爆发的新西兰汤加里罗火山(Tongariro)，进一步表明使用GPS数据监测火山灰的可能性，有助于提高对火山的认识[8]。

但以上研究仍存在诸多问题，如不是所有GPS SNR观测都能监测火山灰[9-10]、信噪比受电离层闪烁的影响，以及火山灰传播速度不清楚。本文利用GPS SNR数据监测2016-10-08的日本阿苏山火山灰变化，给出GPS SNR异常值探测火山烟流的结果，并估计火山灰平均传播速度。

1 理论和方法

GPS信号是广大用户发送的用于导航定位的调制波，它包含处于L波段的载波、测距码、导航电文。常用载波有2种：L1和L2波段的载波频率分别为fL1=154f0=1 575.42 MHz、fL2=120f0=1 227.60 MHz；对应的波长为λ1=19.032 cm、λ2=24.420 cm。现代化后，又增加了L5载波，频率为fL5=115f0=1 176.45 MHz、λ5=25.480 cm。测距码包括C/A码和P码，前者只能调制在L1上，码率为1.023 MHz，后者的码率为10.23 MHz。

另外，GPS接收机通常记录SNR数据，这些数据是载波与噪声密度比(C/N0)，即信号功率与噪声功率的比值。SNR与C/N0之间的关系可以表示为SNR=(C/N0)/B，以分贝(dB)为单位；其中B为带宽。相对应载波L1和L2，信噪比类型主要有S1和S2，此外三频观测还有S5等。SNR数据不提供卫星发射机与接收机之间的距离信息，因此对定位解决方案没有直接贡献，也是由于这个原因，SNR数据通常被地理物理学家和大地测量工作者忽略。然而SNR数据可以直接测得信号阻塞。当GPS信号经过火山灰时，部分信号会发生衰减或者散射，这就意味着到达GPS接收机的信号功率会降低。相反地，大气中的水汽并不直接影响信号功率[7]。这种对水不敏感的特性正是建立L波段雷达观测灰云的原因。PPR方法既受水汽、温度的影响，同时也受火山灰的影响；而SNR仅受火山灰的影响：这也是利用GPS信噪比SNR探测火山灰的优点之一[10]。

信噪比SNR数据值非常依赖高度角。当卫星高度角从5°到65°变化时，SNR的值会缓慢增加，这主要是归功于天线增益的类型。L1 SNR数据在高度角大于25°时有高频噪声(系统噪声或随机噪声)。相反地，海拔更高的L2 SNR数据的高频噪声的级别更小。这种差异是由于C/A码比L2P或L2C短而更易遭受交流干扰。然而，L2P却具有比C/A或L2C更低的SNR值，这是因为接收机在检索过程中不能使用加密码。注意只有在2005年以后发射的卫星上才可以获得L2C码。此外除非用户要求记录L2C码，否则不会记录。

当高度角小于30°时，SNR数据的振荡是由地面反射导致的。反射信号比直接信号经过更长的距离并且与直接信号进行干涉，这将导致需要观测调制。高度角低于30°的SNR数据的频率与土壤湿度、雪深和海平面高度相关。这里，只考虑具有更高的高度角的SNR数据作为探测火山灰的兴趣因子。这主要是因为高度角高的L2P SNR数据具有更好的精确度，所以用它去分离火山灰的影响。

在利用GPS信噪比SNR数据探测火山烟流时，首先要对SNR数据进行预处理，舍弃高度角低于30°的数据，因为低于30°的SNR数据极容易受到地面反射的影响，与土壤湿度、雪深、海平面高度等因素相关[11-15]。通过多项式拟合得到拟合信噪比值，并与实际信噪比值相减，得到信噪比差异值。分析SNR异常值随时间变化的情况，进而利用SNR数据监测火山灰变化。此外，还要注意可以取火山爆发前30 min的SNR数据的标准差作为性能控制尺度，探测值需要是标准差的2.5倍；若探测值不满足这个条件，则该卫星记录的SNR数据不可用。

此外，在条件允许的情况下，需要验证信噪比是否受电离层闪烁的影响，即需要一些接收机数据作为对比，以证明SNR数据的改变不是由电离层闪烁产生的；要想更充分地分析利用SNR数据监测火山灰的优劣性，还需要对比雷达扫描结果，以确定灰尘粒子尺寸和密度对于使用SNR数据分析火山灰的限制。

实际上火山爆发时，喷发出的物质除了呈固态小颗粒的火山灰，还有含有水蒸气和硫化物的火山气体。在海拔较高的地方，火山灰和水蒸气会分离。这是因为在海拔较高的位置，水蒸气依然可以保持充足的密度和体积，但是火山灰却会沉降并向水平方向的扩散，这也导致在水平方向和垂直方向上离火山越远火山灰密度越低；所以距离火山较近的测站受火山灰影响大，而其他距离较远的测站受火山灰影响微弱，即火山灰密度在信噪比SNR的衰减现象中扮演着十分重要的角色。

2 结果与分析

阿苏山(Mount Aso)是日本最大的活火山，也是世界大型的活火山之一，位于九州岛熊本县东北部，最高点海拔达到1 592 m。2016-10-08，当地时间凌晨1：46，日本气象厅(Japan Meterological Agency，JMA)报道阿苏火山中岳火山口发生爆发式喷发。由于当日环境比较恶劣，影响了视觉观察，但是卫星影像显示火山爆发当日火山灰到达了11 km的高度。火山沉积物向西北侧面扩散了1.6 km，向东南侧扩散了1 km，在东北侧也比较多。

我们利用的是日本国土地理院(Geographical Information Authority of Japan，GSI)所建立的GNSS地球观测网(GNSS Earth observation network system，GEONET)数据，这些GNSS数据的采样间隔为30或1 s，采样间隔依赖于接收机所处位置的环境、GNSS接收机的类型和遥测系统。根据文献[6]的研究结果，15到30 s的标准采样率就足以分析经过烟流的信号；本次研究所使用的数据间隔是30 s。

本文选取阿苏山火山30 km范围以内的测站，其测站位置与阿苏山的关系如图1所示，其中较粗曲线部分表示第281天16：00到18：00的卫星轨迹，较细曲线部分表示观测窗口全时段内的卫星轨迹。

针对各个测站提取火山爆发时刻前后16：00—18：00(UTC)时间段内的GNSS SNR数据，分别分析SNR与时间以及卫星高度角之间的变化关系，探索SNR与时间以及卫星高度角之间的相关性，并通过二项式拟合得到对应于各个测站的SNR异常值，然后统计信噪比异常发生的概率。经过对上述9个测站的试验与分析，最终呈现PRN32卫星与阿苏山10 km范围以内的3个测站即0701、0703、0704之间的“测站-卫星对”的结果，其中0703测站的情况如图2所示。

从图2可知：对于0703测站，第279天、第280天、第282天、第283天所对应的信噪比SNR随高度角的增大而逐渐增大，而变化趋势则逐渐减小，没有任何变化异常；而对应于第281天在卫星高度角40°左右时的信噪比SNR明显下降3 dB左右(而对于0701测站和0704测站，第279天、第280天、第281天、第282天、第283天，SNR没有任何异常变化现象，且各天SNR的变化趋势同一，变化曲线基本同一)。故得到初步结果：火山灰仅对0703-PRN32的“测站-卫星对”有影响，这与火山灰的扩散方向和卫星测站对之间的路径是否穿过火山灰等因素有关。

由于接收机观测数据中并未记录到对应于L1的信噪比数据S1，但记录了与L2和L5对应的信噪比数据S2和S5；所以对于有明显下降现象的测站0703，除了考虑上述S2外，又考虑了S5(见图3)。综合图2与图3可以看出，同0703测站的S2数据相比，整体来看S5的值要比S2的值高15 dB左右，而且在高度角40°左右时，第281天的S5也有相对明显的异常下降现象。

为进一步研究火山灰与信噪比异常变化之间的关系，可以利用4阶多项式拟合得到期望SNR值，并与实际SNR值相减，得到SNR差异值。如图4为第281天针对0703测站的PRN32卫星的观测值信噪比差异值随高度角变化的情况。

由图4可知：S2的RMSE为0.683 5，S5的RMSE为0.393 7；在火山爆发前，SNR差异值稳定在[—1，1]这个区间；在火山爆发后，SNR异常差异值振荡现象大，其中对应于S2的最大差异绝对值为3 dB，对应于S5的最大差异绝对值为1.9 dB。

为进一步研究信噪比与火山灰之间的关系，可以从“信噪比的下降现象是否在火山爆发后发生”这一方面考虑(如图5所示)。

由图5可知，信噪比异常下降现象的确发生在火山爆发时刻之后2 min左右，而且异常变化间隔大概在3 min左右，变化间隔较短。这可能一方面是因为在火山爆发瞬间喷发出大量高浓度的火山烟流，而烟流快速上升的同时火山灰也逐渐向侧方向扩散，火山烟流的浓度降低；另一方面因为卫星与测站的相对位置也在发生变化，卫星的直射信号路径随着时间可能穿过浓度较低的火山烟流或者不穿过火山烟流。

经过上述针对各测站的卫星高度角和时间与信噪比之间的变化关系的“控制变量”对比研究分析，即针对同一测站同一卫星研究不同天的SNR变化分析，以及针对不同测站同一卫星同一天的SNR变化分析，排除了由于测站自身因素、测站周边地况、卫星自身因素以及其他天气情况等因素，得到0703测站所对应的PRN32卫星在第281天(火山爆发当天)的SNR异常的确是由于火山灰引起。

通过大量数据综合考虑0703测站所有卫星的信噪比差异值(实际信噪比值与拟合信噪比值相减得到)，并统计各信噪比差异值的概率情况如图6所示。

由图6可知：信噪比差异值出现的概率小于0.12，而差异值在[—3，—1.75]区间内的概率小于0.02；故可以判断信噪比异常现象的出现概率是极其小的，若没有外在因素的影响，信噪比异常现象不会出现，即火山灰的确对GPS信号造成了影响。阿苏山周边9个测站对应的PRN32信噪比SNR2的时间序列变化分析结果如图7所示。

由图7可知：对于PRN32卫星，整体来看，除了0703测站对应的信噪比值在16：50左右有一个突然的降低之外，各测站的信噪比SNR没有骤降而是随时间的增长而逐渐降低而且降低的幅度越来越大；在17：20之前，除了0703测站对应的信噪比差异值在16：50左右最大骤降幅度达到了3 dB Hz，其他测站对应的信噪比差异值绝大多数在[—1，1]之间；而在17：20之后，所有测站的信噪比差异值振荡幅度增大，而考虑到在高度角小于30°时测站接收机受多路径效应的影响，这种大幅度振荡现象也可以得到解释。图8再次验证在“卫星-测站对”的信号路径未经过火山灰时，信噪比差异值绝对值小于1 dB，也从另一个角度反映了测站到火山的距离并不是唯一影响信噪比差异值的因素。

综合火山爆发后火山灰到达卫星测站之间的连线(下文将此相交点称为交点)所经过的距离，以及火山爆发到信噪比发生变化时所经过的时间，大致可以估计出火山灰传播的平均速度。当信噪比发生变化时，已知变化时间dt，卫星高度角ele，测站与火山口之间的距离s，估计的交点到火山口之间连线与铅垂线之间的夹角β(β可以是一个区间角)，可以求出卫星-测站连线与交点-火山口连线之间的夹角α，最后可以根据正弦定理求得交点与火山口之间的距离ρ，进而求出火山灰传播的平均速度v(30 m/s左右)。

3 结束语

本文利用GEONET的SNR数据探测了2016-10-08日本阿苏山火山灰变化，验证利用SNR数据的变化监测火山灰的可行性。在阿苏山火山爆发期间，火山北侧的0703测站所接收到的PRN32卫星的SNR数据受到较大影响，信噪比的值在较短时间内会有明显的降低，在16：50左右最大骤降幅度达到了3 dB；其他卫星几乎没有影响，而其他附近测站SNR数据几乎也没有明显异常现象，主要是这些测站与卫星之间信号没有穿过一定浓度的火山灰。根据火山灰到达卫星测站之间的距离，以及火山爆发到信噪比发生变化时所经过的时间，可以估计出火山灰传播平均速度约为30 m/s。将来可以利用格洛纳斯卫星导航系统(global navigation satellite system,GLONASS)、伽利略卫星导航系统(Galileo satellite navigation system，Galileo)和中国北斗卫星导航系统(BeiDou navigation satellite system，BDS)的SNR数据对火山灰进行探测，并与GPS数据的结果进行对比，进一步研究利用GNSS SNR数据探测火山灰的可行性。此外，还可以进一步考虑火山烟流的成分对GNSS SNR数据的具体影响。

[1] SCOLLO S，PRESTIFILIPPO M，SPATA G，et al.Monitoring and forecasting Etna volcanic plumes[J].Natural Hazards & Earth System Sciences，2009，9(5)：11993.

[2] SCOLLO S，PRESTIFILIPPO M，PECORA E，et al.Height estimation of the 2011—2013 Etna lava fountains[J].Annals of Geophysics，2014，57(2)：1-25.

[3] GOUHIER M，HARRIS A，CALVARI S，et al.Erratum to：Lava discharge during Etna’s January 2011 fire fountain tracked using MSG-SEVIRI[J].Bulletin of Volcanology，2012，74(5)：1261-1261.

[4] HOULIÉ N，BRIOLE P，NERCESSIAN A，et al.Sounding the plume of the 18 August 2000 eruption of Miyakejima volcano(Japan)using GPS[J].Nature，2005，414(5)：583-583.

[5] ARANZULLA M，CANNAVF，SCOLLO S，et al.Volcanic ash detection by GPS signal[J].GPS Solutions，2013，17(4)：485-497.

[6] GRAPENTHIN R，FREYMUELLER J T，KAUFMAN A M.Geodetic observations during the 2009 eruption of Redoubt Volcano，Alaska[J].Journal of Volcanology & Geothermal Research，2013，259(6)：115-132.

[7] LARSON K M.A new way to detect volcanic plumes[J].Geophysical Research Letters，2013，40(11)：2657-2660.

[8] FOURNIER N，JOLLY A D.Detecting complex eruption sequence and directionality from high-rate geodetic observations：The August 6，2012 Te Maari eruption，Tongariro，New Zealand[J].Journal of Volcanology & Geothermal Research，2014(286)：387-396.

[9] ARANZULLA M，CANNAVO F，Scollo S.Detection of Volcanic Plumes by GPS：the 23 November 2013 Episode on Mt.Etna[J].Annals of geophysics(Annali di geofisica)，2014，57(Fast Track 2).

[10] OHTA Y，IGUCHI M.Advective diffusion of volcanic plume captured by dense GNSS network around Sakurajima volcano：a case study of the vulcanian eruption on July 24，2012[J].Earth，Planets and Space，2015，67(1)：157.

[11] LARSON K M，BRAUN J J，SMALL E E，et al.GPS multipath and its relation to near-surface soil moisture content[J].IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing，2010，3(1)：91-99.

[12] NIEVINSKI F G，LARSON K M.Inverse modeling of GPS multipath for snow depth estimation—part I：formulation and simulations[J].IEEE Transactions on Geoscience & Remote Sensing，2014，52(10)：6555-6563.

[13] NIEVINSKI F G，LARSON K M.Inverse modeling of GPS multipath for snow depth estimation—part II：application and validation[J].IEEE Transactions on Geoscience & Remote Sensing，2014，52(10)：6564-6573.

[14] JIN S G，NAJIBI N.Sensing snow height and surface temperature variations in Greenland from GPS reflected signals[J].Advances in Space Research，2014，53(11)：1623-1633.

[15] JIN S G，QIAN X D，KUTOGLU H.Snow depth variations estimated from GPS-reflectometry：a case study in Alaska from L2P SNR data[J].Remote Sensing，2016，8(1)：63.