1 引言

全球定位系统（GPS）具有全天候、全时段获取地面测站高精度绝对位移的能力，因而，它已经广泛应用于地震学领域，成为地壳运动监测和地震监测的一种有力手段，并发展成一门新的学科-GPS地震学（GPS Seismology）^[1-2].最初，利用GPS技术获取地震同震位移（静态位移）^[3-4]，主要用于分析地壳形变特征与地震孕育、发生之间的关系^[5]；之后，随着GPS观测精度和数据处理方法的不断改进和提高，尤其是高频（1Hz）和超高频（20~50Hz）GPS技术的出现^[6]，GPS技术不仅能用于长期连续的地表位移观测，也能用于短期瞬时动态运动的监测^[7-12]，对以记录速度和加速度的地震仪/强震仪起到了重要的补充作用，为研究地震震源破裂过程、地震前后地壳形变短期变化过程等问题提供了多窗口的监测手段^[5]；近年来，全球或区域高频GPS连续观测网络的建立，结合高速通讯传输，使得GPS技术在地震学实时应用（Real-time GPS Seismology）方面成为了研究热点，其主要包括地震震级震源快速确定^[13-15]，断层模型快速反演^[16-17]，地震预警^[18-19]、海啸预警^[20-21]等多方面研究内容.但是受限于现有高采样率GPS单历元解算方法，软件和实时产品等诸因素的制约，在实现实时GPS地震应用这一目标前，仍然有诸多问题需要解决，而高采样GPS数据的实时处理又是其中最为基础和核心的问题.

目前，GPS高频数据处理主要采用单历元瞬时定位算法，获取测站瞬时位移，常用方法有两种：一种是利用相对定位技术，其动态定位模式采用双差模式，该方法需要一个或多个远离震中的台站作为参考站，利用震区监测台站与所选的参考站观测数据一起组成双差观测值，进而解算监测台站每个历元的瞬时坐标.该方法的不足处在于需要不受地震影响的参考站.在强震中，即使远离震中的台站也会在一定程度上受到地震的影响而产生位移，由此得到的震区监测台站的位移会受到参考站位移的影响，而不能真正反映监测台站的绝对真实位移.即使可以找到远离震中的参考站，长距离动态相对定位的精度也难以得到保证.近年来快速发展起来的另一种方法是精密单点定位技术（PPP技术）^[12]，该方法虽然避免了选择参考站的问题，但是该种方法需要从IGS网站下载精密星历和精密钟差，而这些精密产品目前还存在滞后性，难以满足实时性需求.尽管区域增强PPP技术（PPP-RTK技术）能够通过参考站网实时估算卫星轨道和钟差改正数实现实时精密定位^[22]，但是该方法仍面临数据通讯以及实时数据质量控制等诸多问题，此外，在PPP计算中，由于受到很多误差因素的影响，单站模糊度快速固定仍旧是一个还没有很好解决的复杂问题，且目前单站模糊度固定解的可靠性还有待于进一步提高^[23].

为了克服上述两种方法存在的固有局限性，Colosimo等学者从历元间位移变化入手，利用单站载波相位观测值确定测站历元间位移变化量，进而转化为瞬时位移，并用实测GPS数据验证了该方法的可行性^[24]，为GPS应用于实时地震监测提供了一种新的途径.基于这一思想，本文则从更具有物理意义的测站速度着手，利用单站GPS多普勒观测数据和GPS广播星历，基于GPS多普勒频移测速法实时确定台站运动速度，从而实时获取GPS监测站运动状态.同样，该方法也可以克服差分需要参考站和PPP需要IGS实时产品的局限性.

2 单站GPS测速方法

目前利用GPS技术获取运动物体高精度的速度主要是采用GPS接收机测量得到的多普勒观测值直接求解速度.多普勒测速同伪距单点定位类似，仅需广播星历即可，只是观测值不再是伪距观测值，而是多普勒频移观测值或导出多普勒观测值.

2.1 多普勒测速的观测模型

载波相位（L₁载波）的一般距离观测方程如公式（1）^[25-26]，卫星到地面GPS接收机的几何距离见公式（2）：

(1)

(2)

式中，λ₁表示L₁载波波长（~19.04cm），φ_m^j表示载波相位观测值，t表示信号接收时刻，c表示光速，τ_m^j表示卫星j到接收机m之间的传播时间，δ_m^t和δt^j分别表示接收机钟误差和卫星钟误差，T^j _m表示对流层延迟量，I^j _m表示电离层延迟量，N^j _m为整周模糊度，φ ^j _m（0）为初始相位偏差，dM为多路径效应误差，dR为相对论效应误差，ε^j _m为接收机测量噪声，ρ _m^j表示信号发射时刻的卫星位置（x^s，y^s，z^s）到信号接收时刻的接收机位置（x，y，z）之间的几何距离，其中r代表三维位置向量，dr代表相应的位置误差向量.

对GPS距离观测方程进行微分可得多普勒观测方程，在时间间隔很短（＜1s）的条件下，卫星轨道、对流层、电离层、相对论效应等误差的变化率为微小量，可忽略其对测速结果的影响^[27]，因此可得简化的测速观测模型：

(3)

(4)

公式（3）（4）中，卫星的位置rj，卫星的速度ṙ ^j以及卫星的钟速（钟差一次变率）δṫ ^j都可以利用广播星历精确计算^[28]，接收机的位置rm利用伪距单点定位确定，是多普勒观测值.因此，公式（4）中只包含GPS测站速度ṙ _m和钟速δṫ _m四个未知参数，当GPS测站观测到4颗或以上有效卫星，便可以利用最小二乘解法求得GPS测站的速度ṙ_m.

2.2 多普勒观测值的选择

多普勒观测值可以采用接收机直接产生的原始多普勒观测值，也可以采用由载波相位观测值经过时间差分构造的导出多普勒观测值.现有研究表明原始多普勒观测值精度通常低于导出多普勒观测值的精度^[29].本文采用一阶中心差分法构造多普勒观测值（t），其形式为：

(5)

式中， _T为T时刻的导出多普勒观测值，φ_T+ΔT、φT-ΔT分别为T时刻后一历元和前一历元的载波相位观测值，Δt为GPS观测数据的采样间隔.

假设L₁（或者L₂）载波相位观测值精度为σ_φ，根据公式（5）和误差传播定律，可得导出多普勒观测值精度为.一般地，σ_φ取值为1~2mm，则取决于采样间隔，采样率越高，导出多普勒观测值的测量精度越低.若采用L₁载波相位观测值构造多普勒观测值，当采样率为1 Hz时，约为0.7~1.4mm/s，而当采样率为5Hz时，则为3.5~7mm/s.若采用双频消电离层组合观测值φ_IF，即

(6)

可得组合后的载波相位观测值精度为σ_φIF=2.99σ_φ，观测噪声放大.因此本文采用L₁载波相位观测值构造多普勒观测值进行速度求解.

此外，利用载波相位观测值基于一阶中心差分构造多普勒观测值不可避免会受到原始载波相位观测值中周跳的影响.因此，在数据处理时需要对周跳进行处理.本文中对周跳只进行探测，周跳探测方法选用目前最为成熟和有效的TurboEdit方法^[30].

除了受周跳影响外，对于采样率为1 Hz的数据，导出多普勒观测值还会受到电离层残余误差和多路径残余误差等周期性误差影响.恒星日滤波法是一种有效削弱周期性误差影响的方法，该方法利用GPS卫星轨道的恒星日重复性，来降低与卫星轨道和接收机环境相关的周期性误差影响^[31-34].可用于改善GPS测速的精度.需要指出的是，对于5 Hz采样率的数据，由于采样率时间极短，相邻历元间周期性误差几何强相关，一阶中心差分后基本完全削弱了其对构造多普勒观测值的影响，可以不进行该项改正.

3 应用实例

为了验证本文所提方法在实际地震监测应用中的可行性和效果，本文以2010年4月4日发生在墨西哥BajaCalifornia北部的M_w7.2级EI-Mayor-Cucapah地震为例开展研究，该地震发生于22:40:42 UTC（GPSTime-UTC=15s），震中位置为32.259°N，115.287°W，震源深度为10km.选取震中邻近区域200km范围内由UNAVCO-PBO提供的5 Hz高频GPS数据（如图 1所示）进行处理，处理时间范围为22:00:00至23:00:00（GPS Time）.考虑到GPS数据处理结果对比分析方便以及篇幅限制等因素，在后续正文中根据震中距由近及远选择了7个具有邻近强震仪数据的监测站（P500，P496，P744，P499，P487，P485，P484）进行重点分析.有关GPS各站信息具体列于表 1，同时表 1中还列出了与GPS站并置的强震仪测站信息，强震仪数据由USGS-NSMP（美国地质调查局下设的国家强震运动工程）提供.图 1为震中及各站分布示意图.

3.1 地震期间GPS测站速度估计

采用本文提出的方法对GPS测站数据进行模拟实时处理，研究测站在地震期间的瞬时运动状态.卫星轨道采用GPS广播星历，卫星截止高度角设置为7°，最终的速度结果转换至N/E/U方向上进行分析.图 2所示为测站P496速度序列，横轴代表由GPS周秒表示的时间，时间范围从22:00:00到23:00:00（地震发震时刻22:40:57，81657s），纵轴为相应的速度，三幅子图从上至下分别对应北向、东向和垂向速度.从图中可知，震前测站速度分量为零均值的白噪声序列，以22:00:00到22:40:00这段时间测站速度计算，可得N、E和U各方向速度精度分别为6.2mm/s，3.4mm/s和16.8mm/s.在地震波到达测站期间，受地震波影响，测站剧烈晃动，各向速度发生变化，其中N、E方向最大速度接近0.6m/s，U方向最大速度约0.15m/s；随着地震波的消退，测站逐渐趋于静止，速度再次表现为零均值随机噪声序列.由结果可得，采用GPS测速法获得的测站瞬时速度，能够清晰、完整地记录整个地震过程，而其所捕获的地震波信号，可以为地震学后续研究提供数据依据.

图 2给出了测站P496的时域运动状态变化序列，为了进一步分析地震发生前后测站运动信号的频域特征，图 3给出了测站P496N/E/U各方向的速度功率谱密度图，实线表示由地震发震时刻前1200个历元（22:36:57至22:40:56）所得结果，而虚线表示由地震发震时刻后1200个历元（22:40:57至22:44:56）所得结果.图 3显示，地震发生后，测站水平向速度各频率功率谱密度均变大，尤其在低频部分（N/E＜1Hz；U＜0.5 Hz），速度结果在N/E/U三个方向上的功率谱密度都显著增加，由此说明，采用GPS测速所得结果能够明显地记录地震波信号低频部分的信息.从图中还可知，在地震发生之前，U方向的功率谱较之N/E方向的大，说明U方向的速度序列振幅较大，这反映出该方向的测量精度较之其它两方向的低，这与前面所得结论一致.此外，对比U方向震前和震时速度结果功率谱密度，可以发现当频率值超过0.5 Hz时，两者十分接近，在该频率范围内的地震波信号已经无法与测量噪声区分，从另一方面也反映出U方向的测量精度较差，一定程度上限制了U方向的应用.

图 4给出了7个GPS测站在地震期间的瞬时速度，横轴代表GPS时间，纵轴代表相应的速度，图中黑色竖向加粗虚线表示地震发生时刻，由于在未发生地震期间，各站速度趋近于零，为直观展示地震期间各站运动状态，这里仅列出了地震发生后200s内的速度结果.需要说明的是，在图中，为清晰描述所有测站三个方向的速度变化，对各测站各方向的速度序列在图上标示时进行了相应的平移.由图 4结果可知，在地震发生后，地震波从震中传出，各站各向分量在地震到达时均有不同程度的振动，相比垂直方向而言，水平方向的振动更为明显，不仅振幅大，而且持续时间长.产生该现象的原因，一方面是由于地震波携带的能量大部分集中于水平向传播的横波（S波）中，另一方面原因则是由于利用GPS进行测速，其在垂直方向的误差较水平方向更大，垂向信号淹没于测量噪声中.图 4中测站P496在N/E/U三个方向上振动均最为明显，这是因为测站P496处于地震带破裂的方向，由此说明，测站所受地震波的影响程度不仅取决于距离震中的距离，而且还与测站是否处于地震破裂带方向有关.此外，从图中所示结果还可知，各站发生振动的先后顺序与震中距有关，距离震中近的测站先感知（P500，P496，P744，P499），距离震中较远的测站后感知（P487，P485，P484），这是由于地震波在同介质区域传播波速近似相同所致.依此物理特性，当距离震中较近的GPS测站监测到地震波来临，即可向远距离区域发出预警警报，实施异地预警.

3.2 GPS和强震仪并置站速度比较

为了验证GPS估计的测站瞬时速度的可靠性，采用NSMP（NationalStrong-MotionProject）提供的与测站P496和P744相距百米左右的NP5058和NP5028强震仪结果进行对比分析（其余5个并置站由于两站之间距离较远，所受到的地震影响并不相同，作为检验GPS测速的参考值不及距离较近的并置站合适，因此文中并未列出其余并置站对比结果）.本文所采用的强震仪结果为经由CGS/CSMIP（加州地质勘探局/强震观测项目）校正后的速度，因而可以直接作为参考值使用.强震仪数据采样率为200Hz，数据长度为270s.需要注意的是，强震仪数据中的时间为UTC，在与GPS测速结果比较前需要转为GPST，同时单位也需要由cm/s转为m/s.图 5为GPS测站P496直接获取的速度结果同其并置的强震仪测站积分得到的速度结果对比图.为了清晰反映两者的符合程度，仅截取从22:41:05至22:42:20GPST共计75s的结果.在图中，横轴为由GPS周秒表示的时间，纵轴为各方向对应的速度结果.从图 5可知，GPS结果与强震仪结果整体吻合得很好，特别是在地震波引起的地面晃动振幅达到最大及之后时间段内，两者所测地面速度运动状态几乎重合，此时地震波主要表现为长周期信号，由此说明GPS能够较好地捕获地震波低频信号.但是，从地震波初至到测站速度振幅达到最大这段时间内，强震仪记录较之GPS记录抖动更加频繁，说明强震仪捕捉到了地震波高频信号，而GPS却没有记录到.由此说明，GPS受限于采样率低的限制，其结果还无法较好地记录到同强震仪一致的诸多短周期高频地震波信号.图 6给出了测站P744与其并置的强震仪测站NP5028速度对比结果，可得出相同的结论，这里不再详述.

4 结论和展望

本文以GPS测速为基础，寻找到利用GPS进行实时地震监测的新途径.该方法与相对定位或精密单点定位方法求解测站瞬时位移相比，具有显著的优点：一是它避免了单站求解整周模糊度的难题，整体计算模型变得简单；二是它采用的导出多普勒观测值进行了历元间差分，因而影响GPS观测值的低阶误差项得到了削弱和消除；三是它所需的观测信息（GPS观测数据和GPS广播星历），可由GPS接收机自身提供，无需借助参考站网或者其它外部数据产品，适合实时处理.

实验及分析也进一步验证了该方法的可行性.通过GPS测站与其并置的强震仪测站速度结果对比分析，可得两者具有较好的一致性，尤其是长周期低频地震波信号能够很好的被GPS所捕获.此外，对比两者速度最大振幅，十分接近，由此可将GPS测速结果直接用于地震震度确定，并结合多测站速度结果快速生成震度图，用于推估震中位置和预估地震灾情等.但是从两者结果中也发现，GPS受限于采样率低的限制，还无法记录到同强震仪一致的诸多高频信号，今后有必要提高GPS采样率，扩充GPS捕获地震波信号的频带范围.

今后，随着GPS技术的实时性、连续性、精确性、完好性的不断完备，采用高频GPS观测数据进行地震等地质灾害监测具有广泛的应用前景.本文所提方法也还存在一定问题需要解决，比如如何能使速度积分成的位移更加准确，如何能够提高垂向速度精度等.此外，还需充分发掘利用GPS实时速度或者位移结果进行地震预警、地震学反演等相关问题的研究，使GPS技术能更好地应用于地震学领域.

致谢

衷心感谢PBO（Plate Boundary Observatory），IGS（International GNSS Service），USGS（United States Geological Survey），CGS（California Geological Survey）为本文提供的高频GPS数据以及强震仪校正结果文件.本项目研究在武汉大学地球空间环境与大地测量教育部重点实验室完成.