利用高频GPS进行地震动态变形分析及地震定位*

吴继忠 吴文坛

(1)南京工业大学测量系，南京 210009

2)河北省测绘资料档案馆，石家庄050031)

利用高频GPS进行地震动态变形分析及地震定位*

吴继忠1)吴文坛2)

(1)南京工业大学测量系，南京 210009

2)河北省测绘资料档案馆，石家庄050031)

利用精密单点定位技术，得出了日本Mw 9.0地震中8个GPS站的动态变形序列。结果显示：GPS参考站震时最大形变量明显大于永久性变形，东向最大位移达5.9 m，北方向最大位移达1.8 m，水平位移总体朝向震中方向。根据各个参考站记录的地震波到达时间，得到了震中位置及S波传播速度。

高频GPS;日本Mw 9.0地震;精密单点定位;同震变形;地震定位

1 引言

随着GPS误差模型的不断精化、数据处理方法的改进以及高频GPS接收机的出现，利用GPS监测强震地面运动过程得到了越来越多的应用。GPS不仅可以观测到周期小于1秒的位移量，而且可以检测到超长周期的地壳运动。本文将以2011年3月11日日本Mw 9.0地震为例，利用高频GPS观测结果分析同震地表动态变形过程及该区域内8个GPS参考站震后短时运动结果，并利用GPS参考站的动态变形序列实现地震定位。

2 高频GPS数据处理方法

高频GPS数据处理方法主要有双差和非差两种模式。在现有的数据处理软件中，BERNESE软件同时具有用精密单点定位和双差定位获取单历元解的功能，且其中的双差解精度较高。GAMIT软件的动态定位模块TRACK具有用双差定位获取单历元解的功能，在地震GPS数据处理中已经得到广泛应用［1-3］。

双差解算的关键问题之一是模糊度的快速分解，例如在GAMIT软件中是先利用MW组合宽巷模糊度分解公式计算出宽巷模糊度，然后再利用宽巷模糊度对L1和L2模糊度进行分解。对地震观测的GPS数据进行双差解算的一个难点是如何寻找稳定的参考站，选择较远的参考站又不利于模糊度的固定，采用非差解算则可以避免这个问题。非差解算主要是精密单点定位，可以对单个测站的数据进行处理。精密单点定位一般均采用模糊度浮点解，其模糊度不能固定为整数的原因在于接收机与卫星存在一个初始的相位延迟(Uncalibrated Phase Delays，UPD)，致使非差或单差模糊度本身不具备整数性质，只有双差模糊度才具有整数性。在PANDA采用的模糊度固定新方法中，利用地面跟踪网数据估计单差的UPD(SD-UPD)，来去掉单差模糊度的小数部分，从而固定单差模糊度，利用此方法80%的独立模糊度可以被固定，地面点固定解相比浮点解精度提高30%左右［4］。Blewit［5］提出PPP算法中快速解算初始相位模糊整周数的方法——AMBIZAP，该算法逐点进行PPP定位、逐基线进行独立基线解算，定位精度比传统PPP法有较大提高，此算法与GIPSY软件相结合大大提高了GPS数据处理的速度，适合于GPS巨型网的解算，在新一版的GIPSY软件中将融合AMBIZAP方法。

3 日本Mw 9.0地震中GPS参考站的变形分析

日本列岛地处太平洋西北，是太平洋板块、北美板块、菲律宾板块和欧亚板块的交汇处。在该区域太平洋板块以约80 mm/a的速率相对欧亚板块向西运动，插入北美板块下部，与北美板块、欧亚大陆板块、菲律宾板块相互作用，形成了世界上最典型的大型俯冲带——西太平洋俯冲带。在日本北部，太平洋板块俯冲到北美板块与欧亚大陆板块之下，3月11日的Mw9.0地震即发生在太平洋板块与北美板块交界的俯冲带上。地震发生后多家研究机构给出的震源机制解表明地震发生在下插板片上边界，破裂面西倾9°～14°，释放地震矩3.6×1022Nm［6］。JPL和Caltech［6］的ARIA研究团队在震后迅速分析了覆盖日本全境的GEONET网络GPS观测数据，给出了日本境内同震位移场。结果显示，本州岛北部地区存在显著的东向位移，水平位移量在靠近震中区的东海岸达到极大，约为5.3 m，同时这一地区还存在1.1 m左右的沉降，在本州岛北部西海岸也发现约0.7 m的水平位移。本文选择8个参考站进行分析(图1)，这些站的震中距最大为576 km，最小为145 km。

图1 GPS参考站的分布Fig.1 Distribution of GPS reference stations

表1 GPS参考站在UTC 5：40—5：55时的同震位移Tab.1 Coseismic displacements of GPS reference stations during UTC 5：40—5：55

图1中8个GPS参考站在地震发生15分钟内的同震位移量见表1［6］。对比可以发现，参考站0550的震中距最小，其位移量最大，东方向位移达5 m，北方向和高程方向的位移量也超过1m，明显高于其他参考站的位移量。所有参考站位移方向总体向东，高程方向为明显的下沉。

为了进一步分析此次地震的震时地表动态形变过程，本文利用精密单点定位技术(无模糊度固定)获取了8个GPS参考站在地震发生15分钟(UTC5：40—5：55)内的单历元动态解，得到900个历元的动态坐标序列，软件处理过程见文献［7］。图2(a)～(h)显示的是8个GPS参考站在震时水平方向的运动轨迹，从图中可以看出，图2(a)～(h)的共同特点是地震过程中变形的最大幅度远超过其震后变形的幅度，东方向位移量明显大于北方向位移量，各个站南北方向和东西方向在震时均有明显的弹性形变，例如参考站0017在UTC 05：50：09东方向最大变形达到0.39 m，在UTC 5：49：16北方向最大变形达0.25 m，而最终东方向变形量在0.06 m左右，北方向变形量在0.07 m左右，其他参考站也出现类似的弹性变形过程，都是先朝震中方向运动，然后发生回弹。东方向的振幅明显高于北方向，各站点时序的振幅随震中距减弱，显示出地震波能量随距离衰减。就站点的初始位置和最后位置而言，站点0017、0151、0166、0550位于震源的西北方向，其水平运动方向为东偏南，且较为一致;站点0201、0587、0223、3052位于震源的西南方向，其水平运动方向大体朝西，南北向略有差异，动态解显示的站点最后位置位移量与高精度静态解(图2中圆圈的位置)吻合较好，这说明动态解的结果是可靠的。

4 利用GPS参考站进行地震定位

设n个台站的观测到时为t1，t2，…，tn，求震源(x0，y0，z0)及发震时刻t0，并使目标函数最小，即

图2 GPS站点在震时水平方向的运动轨迹Fig.2 Horizontal tracks of GPS stations during earthquake

其中ri为到时残差，ri=ti-t0-Ti(x0，y0，z0)，Ti为震源到第i个台站的计算走时。由于GPS参考站的位置已知，首先需要确定地震波达到各参考站的时刻，由于已用PPP方法获得了参考站在地震发生15分钟(UTC 5：40—5：55)内的单历元动态解，可以通过分析其坐标序列的变化来确定地震波的到时，本文采用文献［3］的方法，即根据其坐标序列计算中误差，取2倍中误差作为阈值，阈值的计算公式为：

式中N、E、U分别为北、东、高坐标分量，n为历元总数(文中n=900)，k为历元序号，σN、σE、σU分别为北、东、高方向相邻历元坐标差的阈值。获得了每个GPS参考站的地震波到时后，对式(1)进行解算，得到相应参数。

式(1)是一个非线性模型，由于待定参数的初始值不易确定，适合采用迭代算法来解算。常见的迭代算法有牛顿法、高斯-牛顿法、信赖域法、Levenberg-Marquardt算法(简称L-M算法)等，其中L-M算法是目前应用较为广泛的一种无条件约束的优化方法，对参数初始值不敏感，能有效地处理奇异和非正定矩阵，本文采用L-M算法进行非线性模型参数估计，收敛准则为式(1)前后两次迭代结果之差小于10-4，经过51次迭代后收敛。参数估计结果得到此次地震的震中位置为38.224°N，142.343°E，S波平均传播速度为2.91 km/s，表2显示了各个参考站的地震波到时及其计算残差。由计算结果可知，震中位置的估计结果与美国地质调查局提供的结果比较接近，S波波速与以往研究结果也基本一致。

表2 地震波到时及其残差(单位：s)Tab.2 Arrival time and calculating residuals of seismic wave(unit：s)

5 结论

1)通过对GPS参考站的静态解和动态解的对比，表明动态解显示的站点最后位置位移量与静态解吻合较好，但动态解能够获取地震所造成的同震位移场时空演化过程，因而高频GPS观测结果更具有应用价值。由于无法获得地震仪的观测数据，因而未能将GPS与地震仪二者结果进行对比。

2)日本Mw 9.0地震中，静态解和动态解均表明GPS参考站东向位移非常显著，总体朝向震中方向，同时从动态解可以看出震时最大地表形变量明显大于永久性形变。

3)根据GPS参考站的位置及相应的地震波达到时间，反演得到震中位置与美国地质调查局提供的结果比较接近，反演得到的S波波速与以往研究结果也基本一致。

致谢 JPL和Caltech的ARIA团队提供了GPS同震位移资料的下载，NGDS提供了高频GPS数据资料，在此一并表示感谢!

1 Wang G Q，et al.Comparisons of ground motions from colocated and closely spaced one-sample-persecond Global Positioning System and accelerograph recordings of the 2003 M 6.5 San Simeon，California，earthquake in the Parkfield region［J］.Bull Seismol Soc Am.，2007，97(1B)：76-90.

2 吴继忠.基于GPS观测的Baja California地震地壳变形分析［J］.武汉大学学报(信息科学版)，2011，36(4)：437 -440.(Wu Jizhong.Crustal deformation analysis derived from GPS observations during Baja California earthquake［J］.Geomatics and Information Science of Wuhan University，2011，36(4)：437-440)

3 徐韶光，等.汶川地震同震形变的静态和动态分析［J］.大地测量与地球动力学，2010，(3)：27-30.(Xu Shaoguang，et al.Static and kinematic analysis od coseismic deformation of wenchuan Ms8.0 earthquake［J］.Journal of Geodesy and Geodynamics，2010，(3)：27-30)

4 施闯，等.卫星导航系统综合分析处理软件PANDA及研究进展［J］.航天器工程，2009，18(4)：64-69.(Shi Chuang，et al.PANDA：comprehensive processing software for satellite navigation systems and its research progress［J］.Spacecraft Engineering，2009，18(4)：64-69)

5 Blewitt，G.Fixed point theorems of GPS carrier phase ambiguity resolution and their application to massive network processing：Ambizap［J］J.Geophys.Res.，2008，113( B12410)：1-12.

6 王敏，等.全球定位系统测定的2011年日本宫城Mw9.0级地震远场同震位移［J］.科学通报，2011，56(20)：1 593 -1 596.(Wang Ming，et al.Far-field coseismic displacements associated with the 2011 Tohoku-oki earthquake in Japan observed by Global Positioning System［J］.Chinese Sci Bull，2011，56(20)：1 593-1 596)

7 方荣新，等.GPS地震仪：PANDA软件测试结果与验证［J］.武汉大学学报(信息科学版)，2011，36(4)：453-456.(Fang Rongxin，et al.GPS Seismometer：PANDA software testing results and validation［J］.Geomatics and Information Science of Wuhan University，2011，36(4)：453 -456)

KINEMATIC ANALYSIS OF COSEISMIC DEFORMATION AND SEISMIC LOCATION USING HIGH-RATE GPS

Wu Jizhong1)and Wu Wentan2)

(1)Department of Surveying，Nanjing Industrial University，Nanjing 210009 2)Hebei Provincial Archives of Surveying＆Mapping，Shijiazhuang050031)

The high-rate GPS data from 8 stations in the 2011 Mw9.0 Tohoku earthquake were processed by using precise point positioning，and the displacement time series at each GPS site were analyzed to derive displacement waveforms for this event.The results indicate that the maximum displacement of GPS stations during the earthquake is larger than the permanent deformation.The kinematic solutions show that the maximum amplitude of east and north displacement is up to 5.9 m，1.8 m respectively，and horizontal displacements were towards the epicenter on the whole.At last，the epicenter and the S wave velocity is derived from the arrival time recorded by the GPS stations in the region.

high-rate GPS;Mw9.0 Tohoku earthquake;precise point positioning;coseismic displacement;seismic location

1671-5942(2012)02-0020-04

2011-11-18

精密工程与工业测量信息局重点实验室开放基金(PF2011-9);地球空间环境与大地测量教育部重点实验室开放研究基金(10-01-12);江西省数字国土重点实验室开放基金(DLLJ201105)

吴继忠，男，1981年生，博士，主要研究方向为GPS卫星定位技术应用.E-mail：jzwumail@163.com

P228.1

A