留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

GPS时间序列中同震和震后形变的自动识别和估计

苏利娜 张勇

苏利娜, 张勇. GPS时间序列中同震和震后形变的自动识别和估计[J]. 武汉大学学报 ● 信息科学版, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016
引用本文: 苏利娜, 张勇. GPS时间序列中同震和震后形变的自动识别和估计[J]. 武汉大学学报 ● 信息科学版, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016
SU Lina, ZHANG Yong. Automatic Detection and Estimation of Coseismic and Postseismic Deformation in GPS Time Series[J]. Geomatics and Information Science of Wuhan University, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016
Citation: SU Lina, ZHANG Yong. Automatic Detection and Estimation of Coseismic and Postseismic Deformation in GPS Time Series[J]. Geomatics and Information Science of Wuhan University, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016

GPS时间序列中同震和震后形变的自动识别和估计

doi: 10.13203/j.whugis20170016
基金项目: 

中国地震局监测、预测、科研三结合课题 CEA-JC/3JH-172703

详细信息
    作者简介:

    苏利娜, 博士生, 主要从事GPS数据处理与地壳形变研究。sulinawhu@163.com

  • 中图分类号: P228

Automatic Detection and Estimation of Coseismic and Postseismic Deformation in GPS Time Series

Funds: 

the Combining of Monitoring/Forecast/Research Fund Program of China Earthquake Administration CEA-JC/3JH-172703

More Information
    Author Bio:

    SU Lina, PhD candidate, specializes in GPS data processing and crustal deformation. E-mail: sulinawhu@163.com

  • 摘要: 地震发生时,附近GPS观测站会捕捉到同震和震后形变,同震和震后形变对地学研究和维持动态参考框架具有重要作用,识别和估计同震和震后形变也是GPS时间序列处理的重要内容。提出综合检校法识别和估计GPS时间序列中的同震和震后形变,首先利用该方法对GPS站建模,综合考虑多源的同震阶跃及兼容性和均方根RMS减小率来判断测站是否受到地震影响,然后用试错法搜索最优衰减常数τ后,利用最小二乘估计包含同震和震后形变的所有参数,最后以真实震例验证了该方法的有效性。
  • 图  1  程序流程图

    Figure  1.  Flow Chart of Program

    图  2  antc站GPS时间序列和模型拟合

    Figure  2.  GPS Time Series and Its Fitting of antc Station

    图  3  2010年Mw 8.7智利Maule地震GPS站点的同震和震后形变图

    Figure  3.  The Co/Post-seismic Deformation Maps of GPS Stations Affected by Mw 8.7 2010 Maule, Chile Earthquake

    表  1  2010年Mw 8.7智利Maule地震的GPS站点的综合检校法判定

    Table  1.   The Result of Comprehensive Inspection for GPS Stations in Mw 8.7 2010 Maule, Chile Earthquake

    站名 (R1-R2)/R2 (R1-R3)/R3 同震模型值/mm 7 d坐标差/mm 30 d坐标差/mm 判断
    N E N E N E N E N E
    unsa 0.07 0.33 0.04 0.20 -2.9 -9.9 -0.9 -6.0 -1.6 -4.5
    cont 207.20 226.84 0.22 0.29 -688.1 -3061.8 -685.9 -3053.8 -688.7 3050.8 同震+震后
    conz 187.28 236.16 0.20 0.24 -673.4 -3057.7 -673.7 -3042.8 -674.2 3056.3 同震+震后
    areq 0.04 0.10 0.03 0.09 -2.7 -1.5 2.8 1.6 1.0 1.4
    iqqe 0.02 0.01 0.01 0.00 0.7 -0.4 3.2 1.8 2.1 2.3
    copo 0.04 -0.02 0.01 0.00 -0.9 -1.2 7.3 0.8 5.7 1.6
    antc 41.47 74.60 0.38 0.39 90.3 -834.2 185.5 -815.5 192.7 -855.7 同震+震后
    ufpr 0.17 0.35 0.08 0.01 -0.9 -3.5 1.2 -1.6 -0.8 -1.8
    sant 43.93 40.65 0.29 0.40 -146.4 -250.6 -141.3 -256.1 -144.1 -263.6 同震+震后
    cfag 6.11 13.71 0.38 0.43 -18.7 -33.5 -15.0 -34.6 -17.9 -34.4 同震+震后
    lhcl 1.29 5.86 1.28 5.75 1.0 -0.4 2.7 -2.3 2.3 -6.2 同震+震后
    coyq 0.03 0.00 0.01 0.00 -4.6 0.0 -0.2 -0.3 -1.6 0.3
    tucu 0.23 0.16 0.04 0.03 -4.4 -3.5 -0.3 -3.4 -1.6 -3.1
    bue1 0.61 3.70 0.23 0.57 4.1 -18.4 0.0 -19.3 -3.4 -19.4 同震+震后
    bue2 0.90 6.85 0.49 0.54 -3.7 -21.6 0.1 -18.1 -3.0 -18.8 同震+震后
    lpgs 0.42 4.81 0.34 0.59 2.1 -15.9 0.5 -17.0 -0.9 -17.3 同震+震后
    parc 0.03 0.02 0.04 0.00 1.4 0.9 1.8 0.6 1.0 0.1
    rio2 0.01 0.02 0.01 0.00 -0.7 -0.5 2.3 0.1 0.7 0.9
    falk 0.04 0.14 0.09 0.00 1.4 -2.1 3.4 -1.4 1.0 -0.3
    autf 0.04 0.00 0.00 0.00 -1.1 -0.5 0.2 -2.4 0.8 -0.1
    vbca 1.65 26.32 1.44 0.58 -0.5 -34.0 0.6 -33.6 1.0 -35.5 同震+震后
    udec 6.68 135.68 6.84 54.49 0.0 0.0 0.0 0.0 0.0 0.0 震后
    tmco 1.93 71.60 1.92 48.40 0.0 0.0 0.0 0.0 0.0 0.0 震后
    sold 5.58 28.55 5.03 28.97 0.0 0.0 0.0 0.0 0.0 0.0 震后
    rgao 10.38 54.56 10.45 46.78 0.0 0.0 0.0 0.0 0.0 0.0 震后
    pclm 77.95 96.66 41.60 47.18 0.0 0.0 0.0 0.0 0.0 0.0 震后
    nihu 2.07 2.83 0.29 0.21 0.0 0.0 0.0 0.0 0.0 0.0 震后
    lmhs 30.24 98.10 25.28 49.02 0.0 0.0 0.0 0.0 0.0 0.0 震后
    crrl 5.74 56.48 5.62 43.93 0.0 0.0 0.0 0.0 0.0 0.0 震后
    chml 32.52 76.78 24.46 54.49 0.0 0.0 0.0 0.0 0.0 0.0 震后
    cbqc 24.51 58.25 19.01 41.55 0.0 0.0 0.0 0.0 0.0 0.0 震后
    arco 15.76 111.16 16.45 51.55 0.0 0.0 0.0 0.0 0.0 0.0 震后
    注:R1为没有同震和震后形变参数的RMS;R2为有同震和震后参数的RMS;R3为仅有同震形变参数的RMS
    下载: 导出CSV

    表  2  2010年Mw 8.7智利Maule地震GPS站点的同震和震后形变

    Table  2.   The Co/Post-seismic Deformations of GPS Stations Affected by Mw 8.7 2010 Maule, Chile Earthquake

    站名 同震形变/mm 震后形变振幅k/mm
    N E U N E U
    cont -688.226 -3 069.632 -36.141 -5.421 -103.423 -15.406
    conz -673.486 -3 063.147 -45.626 -3.103 -117.714 -7.920
    antc 192.831 -846.466 -43.891 33.176 -166.727 55.959
    sant -146.648 -252.708 -21.747 -6.724 -31.889 8.451
    cfag -19.137 -33.798 4.239 -1.537 -1.761 3.571
    lhcl 0.969 -0.442 0.874 5.877 -23.953 -2.131
    bue1 -4.419 -18.917 -4.079 -0.004 -3.959 -1.372
    bue2 -4.868 -22.130 -6.581 -0.573 -6.169 -4.017
    lpgs -2.191 -16.092 -3.818 -1.324 -5.980 -4.643
    vbca 0.319 -35.562 20.767 3.792 -16.190 -13.249
    udec - - - -5.250 -149.621 -32.032
    tmco - - - 2.528 -63.345 -10.037
    sold - - - 16.692 -104.595 -9.904
    rgao - - - -9.220 -53.935 -6.209
    pclm - - - 60.684 -105.512 -75.152
    nihu - - - -4.370 -131.004 -17.253
    lmhs - - - 24.841 -122.520 45.170
    crrl - - - 4.588 -46.147 -12.539
    chml - - - 25.821 -80.596 9.641
    cbqc - - - 29.967 -104.018 -12.760
    arco - - - -12.150 -118.309 -31.955
    下载: 导出CSV

    表  3  11个地震的自动探测受地震影响的GPS站点

    Table  3.   Automatic Detection Result of the Affected GPS Sites for Eleven Earthquakes

    地震 时间 震中地理位置 震级 影响站点
    San Simeon 2003.9767 35.75°, -121.15° Mw 6.6 p278 p576 p067 p526 p295 crbt lows rnch mee2 pkdb mee1 hogs land pomm masw mida uslo carh cand mnmc hunt tblp van1 p304
    Parkfield 2004.7445 35.92°, -120.54° Mw 6.0 rnch pkdb pomm land mida hogs cand carh mnmc lows masw hunt p281 tblp p282 crbt p283 mee2 mee1 p295 p532 p539 uslo
    Alaska 2002.8397 63.52°, -147.53° Mw 7.9 grnr eil1 clgo fair gnaa tlka atw2 cena pot3 tsea chi3 ken1 seld whit
    Izmit 1999.6260 40.70°, 29.98° Mw 7.4 ista tubi ankr
    Okhotsk Sea 2013.3959 54.87°, 153.28° Mw 8.3 pets
    Costa Rica 2012.6817 10.09°, -85.31° Mw 7.6 hua2 caba qsec epza ind1 lmnl saju lafe pne2 vera mana cn30
    Tokachi-oki 2003.7356 42.21°, 143.84° Mw 8.3 stk2 yssk khaj
    Tohoku-oki 2011.1904 38.32°, 142.37° Mw 9.0 khaj mizu tskb bjfs gmsd tsk2 mtka kgni usud yssk aira daej suwn chan mcil shao ksmv
    Padang 2007.6973 -3.78°, 100.99° Mw 8.4 bsat lais slbu mlkn pski msai ntus coco
    Indian Ocean 2012.2773 2.35°, 92.82° Mw 8.6 pbri pbr2 ntus cusv coco xmis ban2 dgar dgav
    Central California 2007.8315 37.37°, -121.78° Mw 5.6 Monb p226 lutz mhcb p218 p253 p227 p228
    下载: 导出CSV
  • [1] Bock Y, Agnew D C, Fang P, et al. Determination of Crustal Deformation from the Landers Earthquake Sequence Using Continuous Geodetic Mea-surement[J]. Nature, 1993, 361:338-340
    [2] Larson K M, Freymueller J T, Philipsen S. Global Plate Velocities from the Global Positioning System[J]. J Geophys Res, 1997, 102(B5):9961-9981 doi:  10.1029/97JB00514
    [3] Gan W J, Zhang P Z, Zheng K S, et al.Present-day Crustal Motion Within the Tibetan Plateau Inferred from GPS Measurements[J].J Geophys Res, 2007, 112(B8):1-14 doi:  10.1029-2005JB004120/
    [4] Yue H, Lay T.Inversion of High-rate (1 sps) GPS Data for Rupture Process of the 11 March 2011 Tohoku Earthquake (Mw 9.1)[J].Geophys Res Lett, 2011, 38:L00G09 https://www.deepdyve.com/lp/wiley/inversion-of-high-rate-1-sps-gps-data-for-rupture-process-of-the-11-EVB1QY7XkL
    [5] Paul S. Earthquake and Volcano Deformation[M]. USA:Princeton University Press, 2010
    [6] Bock Y, Melgar D.Physical Applications of GPS Geodesy:A Review[J]. Reportson Progress in Physics, 2016, 79(10):106801, doi: 10.1088/00344885
    [7] Fang P, Su L, Goldberg D, et al. Modeling Postseismic Deformation of GNSS Global Tracking Network Stations[C].International Symposium Asia Pacific Space Geodynamics, Moscow, Russia, 2015
    [8] Williams S D P. Offsets in Global Positioning System Time Series[J].J Geophys Res, 2003, 108:211-227 http://cn.bing.com/academic/profile?id=b88eb3304d55ea7cc882ee40649f2ac5&encoded=0&v=paper_preview&mkt=zh-cn
    [9] Gazeaux J, Williams S, King M, et al. Detecting Offsets in GPS Time Series:First Results from the Detection of Offsets in GPS Experiment[J]. J Geophys Res, 2013, 118:2397-2407 doi:  10.1002/jgrb.50152
    [10] Riel B, Simons M, Agram P, et al. Detecting Transient Signals in Geodetic Time Series Using Sparse Estimation Techniques[J]. J Geophys Res Solid Earth, 2014, 119(6):5140-5160 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=JJ0233732274
    [11] Zhang J. Continuous GPS Measurements of Crustal Deformation in Southern California[D]. San Diego: Univ of Calif, 1996
    [12] Nikolaidis R. Observation of Geodetic and Seismic Deformation with the Global Positioning System[D]. San Diego: Univ of Calif, 2002
    [13] Shen Z, Jackson D, Feng Y, et al. Postseismic Deformation Following the Landers Earthquake, California, 28 June 1992[J]. Bull Seismol Soc Am, 1994, 84(3):780-791 http://cn.bing.com/academic/profile?id=0a1d301f91d20f2d773ba96d596e4fc0&encoded=0&v=paper_preview&mkt=zh-cn
    [14] Marone C, Scholz C H, Bilham R. On the Mecha-nics of Earthquake Afterslip[J]. J Geophys Res, 1991, 96:8441-8452 doi:  10.1029/91JB00275
    [15] Savage J C, Svarc J. Postseismic Deformation Associated with the 1992 Mw=7.3 Landers Earthquake, Southern California[J].J Geophys Res, 1997, 102:7565-7577 http://cn.bing.com/academic/profile?id=e516a314c49c5888e11df930e9442a62&encoded=0&v=paper_preview&mkt=zh-cn
    [16] Savage J C, Svarc J. Postseismic Relaxation Follo-wing the 1992 M7.3 Landers and 1999 M7.1 Hector Mine Earthquakes, Southern California[J]. J Geophys Res, 2009, 114:B01401 http://cn.bing.com/academic/profile?id=13c3dd0571395b1f2f93a26c7d32b5ea&encoded=0&v=paper_preview&mkt=zh-cn
    [17] Savage J C, Langbein J. Postearthquake Relaxation After the 2004 M6 Parkfield, California, Earthquake and Rate-and-State Friction[J]. J Geophys Res, 2008, 113:B10407 doi:  10.1029/2008JB005723
    [18] Gonzalez-Ortega A, Fialko Y, Sandwell D, et al. El Mayor-Cucapah (Mw 7.2) Earthquake:Early Near-field Postseismic Deformation from InSAR and GPS Observations[J]. J Geophys Res, 2014, 119:1482-1497 doi:  10.1002/2013JB010193
    [19] Dong D, Fang P, Bock Y, et al. Spatiotemporal Filtering Using Principal Component Analysis and Karhunen-Loeve Expansion Approaches for Regional GPS Network Analysis[J]. J Geophys Res, 2006, 111:B03405 http://d.old.wanfangdata.com.cn/NSTLQK/NSTL_QKJJ02751244/
  • [1] 卞畏畏, 伍吉仓, 张磊, 高宇.  强震时空统计分析及InSAR同震形变场空间分布特征 . 武汉大学学报 ● 信息科学版, 2022, 47(6): 875-886. doi: 10.13203/j.whugis20220176
    [2] 李艳艳, 殷海涛, 韩林桥.  利用MSMPCA去噪的CEEMD方法监测高频GNSS同震形变 . 武汉大学学报 ● 信息科学版, 2022, 47(3): 352-360. doi: 10.13203/j.whugis20190356
    [3] 余鹏飞, 陈威, 乔学军, 赵斌, 李刚, 熊维.  基于多源SAR数据的2022年门源Ms 6.9地震同震破裂模型反演研究 . 武汉大学学报 ● 信息科学版, 2022, 47(6): 898-906. doi: 10.13203/j.whugis20220093
    [4] 柴海山, 陈克杰, 魏国光, 方荣新, 邹蓉, 祝会忠.  北斗三号与超高频GNSS同震形变监测:以2021年青海玛多Mw 7.4地震为例 . 武汉大学学报 ● 信息科学版, 2022, 47(6): 946-954. doi: 10.13203/j.whugis20220140
    [5] 李志才, 丁开华, 张鹏, 温扬茂, 赵利江, 陈建峰.  GNSS观测的2021年青海玛多地震(Mw 7.4)同震形变及其滑动分布 . 武汉大学学报 ● 信息科学版, 2021, 46(10): 1489-1497. doi: 10.13203/j.whugis20210301
    [6] 李宇磊, 张永志, 王帅, 刘泰.  2015年尼泊尔Mw7.8地震的震后形变机制—震后余滑效应和黏滞性松弛效应 . 武汉大学学报 ● 信息科学版, 2020, 45(10): 1578-1587. doi: 10.13203/j.whugis20180349
    [7] 蔡杰华, 张路, 董杰, 董秀军, 廖明生, 许强.  九寨沟震后滑坡隐患雷达遥感早期识别与形变监测 . 武汉大学学报 ● 信息科学版, 2020, 45(11): 1707-1716. doi: 10.13203/j.whugis20200263
    [8] 郭南男, 赵静旸.  一种改进的GPS区域叠加滤波算法 . 武汉大学学报 ● 信息科学版, 2019, 44(8): 1220-1225. doi: 10.13203/j.whugis20180049
    [9] 许才军, 邓长勇, 周力璇.  利用方差分量估计的地震同震滑动分布反演 . 武汉大学学报 ● 信息科学版, 2016, 41(1): 37-44. doi: 10.13203/j.whugis20150500
    [10] 杨国华, 朱爽, 梁洪宝, 杨博.  芦山Ms 7.0级地震震前及同震地表形变 . 武汉大学学报 ● 信息科学版, 2015, 40(1): 121-127.
    [11] 裴媛媛, 廖明生, 王寒梅.  时间序列SAR影像监测堤坝形变研究 . 武汉大学学报 ● 信息科学版, 2013, 38(3): 266-269.
    [12] 时间序列SAR影像监测堤坝形变研究 . 武汉大学学报 ● 信息科学版, 2013, 38(3): 266-269.
    [13] 李志才, 张鹏, 温扬茂, 廖瑛.  利用GPS和海底基准点观测形变反演日本大地震(Mw 9.0)同震断层滑动分布 . 武汉大学学报 ● 信息科学版, 2013, 38(1): 40-43.
    [14] 丁开华, 许才军, 温扬茂.  汶川地震震后形变的GPS反演 . 武汉大学学报 ● 信息科学版, 2013, 38(2): 131-135.
    [15] 许才军, 何平, 温扬茂, 张磊.  日本2011 Tohoku-Oki Mw 9.0级地震的同震形变及其滑动分布反演:GPS和InSAR约束 . 武汉大学学报 ● 信息科学版, 2012, 37(12): 1387-1391.
    [16] 董晓燕, 丁晓利, 李志伟, 王永哲.  一种新的SAR像素偏移量估计流程及其在同震形变监测中的应用 . 武汉大学学报 ● 信息科学版, 2011, 36(7): 789-792.
    [17] 刁法启, 熊熊, 李军, 郑勇.  利用GPS观测反演Tokachi-Oki地震(日本)及强余震的同震滑移分布 . 武汉大学学报 ● 信息科学版, 2010, 35(3): 270-273.
    [18] 方荣新, 施闯, 辜声峰.  基于PPP动态定位技术的同震地表形变分析 . 武汉大学学报 ● 信息科学版, 2009, 34(11): 1340-1343.
    [19] 邢乐林, 李建成, 李辉, 孙文科.  Sumatra-Adaman大地震同震和震后形变的GRACE卫星检测 . 武汉大学学报 ● 信息科学版, 2009, 34(9): 1080-1084.
    [20] 彭颖, 许才军, 刘洋.  联合地震位错模型和InSAR数据构建2017年九寨沟Mw6.5地震同震三维形变场 . 武汉大学学报 ● 信息科学版, 0, 0(0): 0-0. doi: 10.13203/j.whugis20200289
  • 加载中
图(3) / 表(3)
计量
  • 文章访问数:  1022
  • HTML全文浏览量:  78
  • PDF下载量:  492
  • 被引次数: 0
出版历程
  • 收稿日期:  2017-01-14
  • 刊出日期:  2018-10-05

GPS时间序列中同震和震后形变的自动识别和估计

doi: 10.13203/j.whugis20170016
    基金项目:

    中国地震局监测、预测、科研三结合课题 CEA-JC/3JH-172703

    作者简介:

    苏利娜, 博士生, 主要从事GPS数据处理与地壳形变研究。sulinawhu@163.com

  • 中图分类号: P228

摘要: 地震发生时,附近GPS观测站会捕捉到同震和震后形变,同震和震后形变对地学研究和维持动态参考框架具有重要作用,识别和估计同震和震后形变也是GPS时间序列处理的重要内容。提出综合检校法识别和估计GPS时间序列中的同震和震后形变,首先利用该方法对GPS站建模,综合考虑多源的同震阶跃及兼容性和均方根RMS减小率来判断测站是否受到地震影响,然后用试错法搜索最优衰减常数τ后,利用最小二乘估计包含同震和震后形变的所有参数,最后以真实震例验证了该方法的有效性。

English Abstract

苏利娜, 张勇. GPS时间序列中同震和震后形变的自动识别和估计[J]. 武汉大学学报 ● 信息科学版, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016
引用本文: 苏利娜, 张勇. GPS时间序列中同震和震后形变的自动识别和估计[J]. 武汉大学学报 ● 信息科学版, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016
SU Lina, ZHANG Yong. Automatic Detection and Estimation of Coseismic and Postseismic Deformation in GPS Time Series[J]. Geomatics and Information Science of Wuhan University, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016
Citation: SU Lina, ZHANG Yong. Automatic Detection and Estimation of Coseismic and Postseismic Deformation in GPS Time Series[J]. Geomatics and Information Science of Wuhan University, 2018, 43(10): 1504-1510. doi: 10.13203/j.whugis20170016
  • 长期以来,GPS连续观测手段被广泛应用于板块和断裂运动、地震破裂等地球物理领域的研究,对于揭示地震孕育发生起到了重要作用[1-4]。地震孕育、发生的过程是地壳运动产生应力应变积累、集中,最终导致地壳局部破裂,释放弹性应变能的过程[5]。一个完整的地震周期包括了地壳的应力积累引起的震间形变,断层突然破裂产生的同震形变,以及地壳和上地幔逐步恢复到稳态的震后形变过程[6],GPS连续观测在监测地震周期方面发挥了重要作用。同震和震后形变信息的提取,一方面对于研究地震破裂的物理机制和地壳流变学起着至关重要的作用;另一方面,改正同震和震后形变的影响,得到更准确的震间速度,是研究板块运动、应力应变积累和地震孕育,乃至地震预报的基础。此外,提供同震和震后形变参数也是定义和维持动态国际地球参考框架(international terrestrial reference frame,ITRF)的重要部分[7]。因此,GPS时间序列中的同震和震后形变的研究显得十分重要。

    国内外发生了多次显著地震,震中周围上百千米甚至上千千米的GPS站点都会受到影响[6],震后形变持续几年甚至几十年[5-6]。筛选出受影响的GPS站点,建立合适的模型是GPS时间序列分析处理的重要工作之一。为避免冗繁的人工挑选,有学者提出了自动探测阶跃[8-9]和瞬态形变[10]的方法,但未知时刻阶跃和瞬态形变的探测难度很高,容易漏判或者误判,准确度和探测幅度都不如人工方法,因此仍以人工识别为主。针对同震和震后形变的识别,一般是在地震发生后,人工检查震中一定范围内所有GPS站的时间序列,判断是否受到同震和震后形变影响,或者是根据GPS时间序列的形态查找阶跃和非线性形变对应的时刻和区域是否曾发生地震。人工筛选法主观性强,工作量大,有时也会有疏漏和错误。本文提出自动探测和估计同震和震后形变方法(automatic detection and estimation of coseismic and postseismic deformation, ADECP),利用综合检校法自动判断GPS站点是否有同震和震后形变发生,整体估计包括同震和震后形变的所有参数。最后,基于SOPAC(Scripps Orbit and Permanent Array Center)提供的全球跟踪站和加州GPS连续站时间序列数据,以智利Maule地震和加州Hector Mine地震等10余个地震为例,对该方法的有效性进行了验证。

    • 原始GPS时间序列包含了站点速度、天线移动或设备更换引起的阶跃、同震和震后形变、年周期和半年周期的季节项,还有一些未模型化的瞬态形变和噪声。考虑到历元间分量的相关性非常小[11]ti历元的坐标分量(ΔN, ΔE, ΔU)时间序列可以模型化为[12]

      $$ \begin{array}{l} y\left( {{t_i}} \right) = a + b{t_i} + c\;{\rm{sin}}\left( {2{\rm{ \mathsf{ π} }}{t_i}} \right) + d\;{\rm{cos}}\left( {2{\rm{ \mathsf{ π} }}{t_i}} \right) + e\;{\rm{sin}}\left( {4{\rm{ \mathsf{ π} }}{t_i}} \right) + f\;{\rm{cos}}\left( {4{\rm{ \mathsf{ π} }}{t_i}} \right) + \sum\limits_{j = 1}^{ng} {{g_j}} H\left( {{t_i}- {T_g}_{_j}} \right){t_i} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\sum\limits_{j = 1}^{nh} {{h_j}} H\left( {{t_i}- {T_h}_{_j}} \right){t_i} + \sum\limits_{j = 1}^{nk} {{k_j}} \left( {1- {{\rm{e}}^{\left[{-\left( {\frac{{{t_i}-T{k_j}}}{{\tau j}}} \right)} \right]}}} \right)H\left( {{t_i} -{T_{{k_j}}}} \right) + \varepsilon \\ \end{array} $$ (1)

      式中,ti表示以t0起的历元,单位为a; a是初始位移值; 线性速度b表示震间速度; cdef表示年周期和半年周期项系数;gTg历元由于天线或设备改变等原因引起的阶跃;hTh历元发生的同震形变;指数模型${k_j}\left( {1- {{\rm{e}}^{\left[{-\left( {\frac{{{t_i}-T{k_j}}}{{\tau j}}} \right)} \right]}}} \right)$表示Tk历元发生振幅为k、衰减常数为τ的震后形变;ε为观测噪声;H代表Heaviside阶梯函数:

      $$ H = \left\{ \begin{array}{l} 0, {t_i}-T < 0\\ 1, {t_i}-T \ge 0 \end{array} \right. $$

      常用的震后形变模型除基于地壳线性粘弹性流变机制的指数模型[13]外,还有基于余滑的对数模型[14]:${k_j}{\rm{log}}\left( {1 + \frac{{{t_i}-{T_{{k_j}}}}}{{{\tau _j}}}} \right)$,均已广泛应用于1992年Mw 7.2 Landers、1999年Mw 7.1 Hector Mine、2004年Mw 6.0 Parkfield、2010年Mw 7.2 El Mayor-Cucapah等震后的形变研究[13, 15-18]

    • 本文的程序ADECP分为识别和估计两个模块。识别模块首先计算站点到震中的距离,然后根据震级选择较大范围内的所有GPS站作为可疑站点,利用综合检校法判断可疑站点是否受地震影响。如果是,则增加同震和震后形变参数,进入估计模块。估计模块中首先利用trail-and-error试错法搜索最优衰减常数τ[15, 19],然后将τ作为已知值,用最小二乘方法估算各参数大小和不确定度,具体程序流程见图 1

      图  1  程序流程图

      Figure 1.  Flow Chart of Program

    • 当地震发生时,受影响的GPS时间序列会有明显的同震和震后形变,模型拟合时增加同震和震后形变参数可以有效减小RMS值,RMS值减小越多,说明发生同震和震后形变越大,基于此,本文利用RMS值减小率来判断GPS时间序列是否存在同震和震后形变。有些站可能存在一些未模型化的非线性形变,为了防止将这些非线性形变误判为震后效应,还需计算同震形变和不同来源的同震形变的一致性来进一步检验, 增加可靠性,即为考虑同震阶跃和RMS减小率的综合检校法。

      具体而言,综合以下3条准则进行判断:

      1) 分别计算有、无同震和震后形变参数的模型值RMS,对比RMS减小率,以此判断是否可能存在同震和震后形变。

      2) 分别计算有、无震后形变参数的模型值RMS,对比RMS减小率,以此判断是否可能只存在同震形变而没有震后形变。

      3) 多种方法计算同震阶跃,第1种是模型计算同震形变大小,第2种是通过坐标差计算发震时刻的同震阶跃,分别计算发震时刻较长时间和较短时间坐标差,这里采用一个月前后坐标差和一周前后的坐标差,时间段可修改。

      分析3种同震阶跃的大小和彼此的一致性,以此排除其他未模型化的非线性形变的干扰,减少误判。如果是其他类似非线性变化,增加参数,RMS势必会减小,会产生虚假的同震模型值,但这些非线性变化在地震时刻没有阶跃,利用坐标差可以排除这种情况。值得注意的是,此处计算的3种来源的同震形变仅作为自动探测的依据,而最终同震形变大小取决于模型整体估算结果。

      若1)、2)、3)条准则均成立,则认为存在同震和震后形变,在模型中增加同震和震后形变参数;若1)、3)成立,2)不成立,则认为地震影响较小,仅存在同震形变,在模型中只增加同震形变参数;若1)、2)准测成立,表明站点震后开始观测,仅受到震后形变影响;否则认为该点不受到地震影响。

      在进行RMS对比时:

      1) 首先要对数据进行预处理,移除速度、天线或设备改变引起的阶跃,年周期和半年周期项,以及其他地震的同震和震后形变,防止没有同震和震后形变参数时其他参数补偿拟合值而导致RMS减小率不准确。

      2) 统计有无同震和震后形变参数RMS减小率时,只统计震前、后各一年的RMS,防止数据太长,RMS减小效果被降低。

      3) 由于高程精度低,本方法只检校N、E方向,一旦任一水平方向检测出有同震和震后形变发生,则N、E、U这3个方向需同时增加同震和震后的形变参数。

    • 从式(1)可以看出GPS时间序列模型关于x =[a b c d e f g h k τ]T是非线性的,如果假定衰减常数τ已知,模型可以转化为关于x =[a b c d e f g h k]T的线性模型,就可以利用线性最小二乘求得各参数x。本文采用试错法搜索最优τ,即试算一定范围内所有值逐一作为τ,利用最小二乘进行拟合,产生最小卡方值(观测值和模型值残差平方和)对应的τ值为最优τ,该方法可以避免非线性拟合陷入局部最优。一旦最优τ确定,再采用最小二乘重新估计各参数和不确定度。值得注意的是,在同震和震后形变识别时,最优衰减常数τ还未确定,可以根据经验设定一个先验τ,该值大小不影响自动识别结果。

    • 本文选取2010-02-27 Mw 8.7智利Maule地震和1999-09-12 Mw 7.1美国加州Hector Mine地震为例来验证程序的可靠性。GPS数据为SOPAC提供的原始的GPS站时间序列,超过3倍四分位距(IQR)[12],并与模型值差水平50 mm、垂直100 mm视为粗差并剔除。

    • 程序ADECP将距震中3 000 km内共32个GPS站作为可疑站点,识别出其中21个GPS站受到该次地震的影响,cont等10个站存在同震和震后形变,udec等11个站建于震后,仅发生震后形变。表 1列出了综合检校法的各项值以及判断结果。该识别结果与SOPAC人工筛查的结果一致,表明综合检校法可以有效探测出受该次地震影响的站点。

      表 1  2010年Mw 8.7智利Maule地震的GPS站点的综合检校法判定

      Table 1.  The Result of Comprehensive Inspection for GPS Stations in Mw 8.7 2010 Maule, Chile Earthquake

      站名 (R1-R2)/R2 (R1-R3)/R3 同震模型值/mm 7 d坐标差/mm 30 d坐标差/mm 判断
      N E N E N E N E N E
      unsa 0.07 0.33 0.04 0.20 -2.9 -9.9 -0.9 -6.0 -1.6 -4.5
      cont 207.20 226.84 0.22 0.29 -688.1 -3061.8 -685.9 -3053.8 -688.7 3050.8 同震+震后
      conz 187.28 236.16 0.20 0.24 -673.4 -3057.7 -673.7 -3042.8 -674.2 3056.3 同震+震后
      areq 0.04 0.10 0.03 0.09 -2.7 -1.5 2.8 1.6 1.0 1.4
      iqqe 0.02 0.01 0.01 0.00 0.7 -0.4 3.2 1.8 2.1 2.3
      copo 0.04 -0.02 0.01 0.00 -0.9 -1.2 7.3 0.8 5.7 1.6
      antc 41.47 74.60 0.38 0.39 90.3 -834.2 185.5 -815.5 192.7 -855.7 同震+震后
      ufpr 0.17 0.35 0.08 0.01 -0.9 -3.5 1.2 -1.6 -0.8 -1.8
      sant 43.93 40.65 0.29 0.40 -146.4 -250.6 -141.3 -256.1 -144.1 -263.6 同震+震后
      cfag 6.11 13.71 0.38 0.43 -18.7 -33.5 -15.0 -34.6 -17.9 -34.4 同震+震后
      lhcl 1.29 5.86 1.28 5.75 1.0 -0.4 2.7 -2.3 2.3 -6.2 同震+震后
      coyq 0.03 0.00 0.01 0.00 -4.6 0.0 -0.2 -0.3 -1.6 0.3
      tucu 0.23 0.16 0.04 0.03 -4.4 -3.5 -0.3 -3.4 -1.6 -3.1
      bue1 0.61 3.70 0.23 0.57 4.1 -18.4 0.0 -19.3 -3.4 -19.4 同震+震后
      bue2 0.90 6.85 0.49 0.54 -3.7 -21.6 0.1 -18.1 -3.0 -18.8 同震+震后
      lpgs 0.42 4.81 0.34 0.59 2.1 -15.9 0.5 -17.0 -0.9 -17.3 同震+震后
      parc 0.03 0.02 0.04 0.00 1.4 0.9 1.8 0.6 1.0 0.1
      rio2 0.01 0.02 0.01 0.00 -0.7 -0.5 2.3 0.1 0.7 0.9
      falk 0.04 0.14 0.09 0.00 1.4 -2.1 3.4 -1.4 1.0 -0.3
      autf 0.04 0.00 0.00 0.00 -1.1 -0.5 0.2 -2.4 0.8 -0.1
      vbca 1.65 26.32 1.44 0.58 -0.5 -34.0 0.6 -33.6 1.0 -35.5 同震+震后
      udec 6.68 135.68 6.84 54.49 0.0 0.0 0.0 0.0 0.0 0.0 震后
      tmco 1.93 71.60 1.92 48.40 0.0 0.0 0.0 0.0 0.0 0.0 震后
      sold 5.58 28.55 5.03 28.97 0.0 0.0 0.0 0.0 0.0 0.0 震后
      rgao 10.38 54.56 10.45 46.78 0.0 0.0 0.0 0.0 0.0 0.0 震后
      pclm 77.95 96.66 41.60 47.18 0.0 0.0 0.0 0.0 0.0 0.0 震后
      nihu 2.07 2.83 0.29 0.21 0.0 0.0 0.0 0.0 0.0 0.0 震后
      lmhs 30.24 98.10 25.28 49.02 0.0 0.0 0.0 0.0 0.0 0.0 震后
      crrl 5.74 56.48 5.62 43.93 0.0 0.0 0.0 0.0 0.0 0.0 震后
      chml 32.52 76.78 24.46 54.49 0.0 0.0 0.0 0.0 0.0 0.0 震后
      cbqc 24.51 58.25 19.01 41.55 0.0 0.0 0.0 0.0 0.0 0.0 震后
      arco 15.76 111.16 16.45 51.55 0.0 0.0 0.0 0.0 0.0 0.0 震后
      注:R1为没有同震和震后形变参数的RMS;R2为有同震和震后参数的RMS;R3为仅有同震形变参数的RMS

      一旦确定站点受到地震影响,相应地增加同震或震后形变参数,进入参数估计模块。以sant站为例,利用trial-and-error方法历数0~500之间每一个值为τ,计算相应卡方值。τ=113时卡方值最小,则最优τ为113,然后用最小二乘方法重新估计各参数和不确定度。图 2为antc站去除速度趋势的时间序列图,黑色圆点为观测值,红色线为模型拟合值,蓝色线为同震发生时刻。增加同震和震后形变参数,模型拟合的RMS分别为2.8 mm、4.57 mm和6.28 mm,拟合结果很好。用同样的方法计算其他20个站,利用时间序列模型整体估计得到所有站的同震和震后形变大小(见表 2图 3)。图 3中黑色实心点为所有GPS站点,红色圆圈表示受到地震影响的站点,蓝色实心表示时间序列中存在同震和震后形变,蓝色三角表示该站建于震后,仅存在震后形变,箭头显示同震(图 3(a))和震后形变振幅k(图 3(b))。从表 2图 3可以看出同震和震后形变主要集中在30°~40°S。本次地震为逆冲性地震,Nazca板块俯冲于南美洲板块之下,因此东西向比南北向受影响范围大。cont和conz位于震中西侧且距离最近,发生最大同震形变,而位于震中东侧的antc发生最大震后形变。

      图  2  antc站GPS时间序列和模型拟合

      Figure 2.  GPS Time Series and Its Fitting of antc Station

      表 2  2010年Mw 8.7智利Maule地震GPS站点的同震和震后形变

      Table 2.  The Co/Post-seismic Deformations of GPS Stations Affected by Mw 8.7 2010 Maule, Chile Earthquake

      站名 同震形变/mm 震后形变振幅k/mm
      N E U N E U
      cont -688.226 -3 069.632 -36.141 -5.421 -103.423 -15.406
      conz -673.486 -3 063.147 -45.626 -3.103 -117.714 -7.920
      antc 192.831 -846.466 -43.891 33.176 -166.727 55.959
      sant -146.648 -252.708 -21.747 -6.724 -31.889 8.451
      cfag -19.137 -33.798 4.239 -1.537 -1.761 3.571
      lhcl 0.969 -0.442 0.874 5.877 -23.953 -2.131
      bue1 -4.419 -18.917 -4.079 -0.004 -3.959 -1.372
      bue2 -4.868 -22.130 -6.581 -0.573 -6.169 -4.017
      lpgs -2.191 -16.092 -3.818 -1.324 -5.980 -4.643
      vbca 0.319 -35.562 20.767 3.792 -16.190 -13.249
      udec - - - -5.250 -149.621 -32.032
      tmco - - - 2.528 -63.345 -10.037
      sold - - - 16.692 -104.595 -9.904
      rgao - - - -9.220 -53.935 -6.209
      pclm - - - 60.684 -105.512 -75.152
      nihu - - - -4.370 -131.004 -17.253
      lmhs - - - 24.841 -122.520 45.170
      crrl - - - 4.588 -46.147 -12.539
      chml - - - 25.821 -80.596 9.641
      cbqc - - - 29.967 -104.018 -12.760
      arco - - - -12.150 -118.309 -31.955

      图  3  2010年Mw 8.7智利Maule地震GPS站点的同震和震后形变图

      Figure 3.  The Co/Post-seismic Deformation Maps of GPS Stations Affected by Mw 8.7 2010 Maule, Chile Earthquake

    • 程序ADECP将距震中约600 km内共172个GPS站作为可疑站点,逐个判断其中131个GPS站受到该次地震的影响。综合检校法的各项值以及判断结果因篇幅有限未列出,结果表明综合检校法可以有效探测出受该次地震影响的站点。与SOPAC人工筛选结果相比,离震中较近的测站,发生形变较大,判断结果完全一致,而whc1等共计16个站点判别上存在差异,这些站点离震中较远,形变非常微小,是否受到同震和震后形变本身具有模糊性,人工判断通过观察GPS时间序列形态来判断,而ADECP通过分析RMS减小率以及同震形变来判断, 人工判断和ADECP方法存在主客观差异。但这些站是否增加同震或震后形变参数影响非常小,两者判别差异可接受。

      基于时间序列模型整体估计可得到所有站的同震和震后形变大小, 因篇幅有限未列出。在所有观测站点中,最大同震形变发生在ldes站,N、E、U方向同震位移分别为191.4 mm、62.449 mm和28.621 mm。而最大震后形变发生在oaes站,震后形变振幅N、E、U方向分别为-49.193 mm、41.808 mm和29.039 mm, 而ldes振幅为43.022 mm、35.36 mm和30.514 mm。

    • 本文提出的综合检校法成功识别1999-10-16 Mw 7.1美国加州Hector Mine地震和2010-02-27 Mw 8.7智利Maule地震影响的GPS测站,该方法在2011年日本Mw 9.0地震等其他11个地震中也得到成功应用,表 3列出了自动探测的存在同震或震后形变的GPS站。需要说明的是,本文所分析的GPS时间序列为SOPAC提供的全球和加州连续站,全球站不包含各国自有的站点。总体而言,自动化识别和估计同震和震后形变程序准确而高效,免除了繁冗的手动筛选过程,尤其是密集的GPS网络,比如1999年Mw 7.1 Hector Mine地震影响了加州百余个GPS连续站。

      表 3  11个地震的自动探测受地震影响的GPS站点

      Table 3.  Automatic Detection Result of the Affected GPS Sites for Eleven Earthquakes

      地震 时间 震中地理位置 震级 影响站点
      San Simeon 2003.9767 35.75°, -121.15° Mw 6.6 p278 p576 p067 p526 p295 crbt lows rnch mee2 pkdb mee1 hogs land pomm masw mida uslo carh cand mnmc hunt tblp van1 p304
      Parkfield 2004.7445 35.92°, -120.54° Mw 6.0 rnch pkdb pomm land mida hogs cand carh mnmc lows masw hunt p281 tblp p282 crbt p283 mee2 mee1 p295 p532 p539 uslo
      Alaska 2002.8397 63.52°, -147.53° Mw 7.9 grnr eil1 clgo fair gnaa tlka atw2 cena pot3 tsea chi3 ken1 seld whit
      Izmit 1999.6260 40.70°, 29.98° Mw 7.4 ista tubi ankr
      Okhotsk Sea 2013.3959 54.87°, 153.28° Mw 8.3 pets
      Costa Rica 2012.6817 10.09°, -85.31° Mw 7.6 hua2 caba qsec epza ind1 lmnl saju lafe pne2 vera mana cn30
      Tokachi-oki 2003.7356 42.21°, 143.84° Mw 8.3 stk2 yssk khaj
      Tohoku-oki 2011.1904 38.32°, 142.37° Mw 9.0 khaj mizu tskb bjfs gmsd tsk2 mtka kgni usud yssk aira daej suwn chan mcil shao ksmv
      Padang 2007.6973 -3.78°, 100.99° Mw 8.4 bsat lais slbu mlkn pski msai ntus coco
      Indian Ocean 2012.2773 2.35°, 92.82° Mw 8.6 pbri pbr2 ntus cusv coco xmis ban2 dgar dgav
      Central California 2007.8315 37.37°, -121.78° Mw 5.6 Monb p226 lutz mhcb p218 p253 p227 p228

      筛选阈值的设定也是综合检校法的重点,既要把受影响的站点尽可能挑选出来,又要避免误选,可以根据实际情况和站点质量调整。10余个地震的结果表明ADECP程序可以基本筛选出受影响的测站,对一些观测质量较差的站点也可以准确判断。在发生微小形变的情况,同震和震后形变具有一定模糊性,由于主客观差异,与人工探测结果会不一致,但这些形变非常微小,对于时间序列的拟合和速度的影响非常小。

      其中,ADECP的识别模块分析震前、震后各一年的数据判断GPS连续站是否存在同震和震后形变,而参数估计模块则利用全部时段数据整体最小二乘估计各参数。参数估计的准确性一方面取决于数据质量,另一方面取决于同震和震后形变是否被正确识别,其他阶跃是否被准确标注时间。

      本文提出的自动探测方法基于已知的发震时刻和数学模型,判断GPS站是否存在同震和震后形变,因此准确度较高,但不适用于探测未知时刻的阶跃和其他瞬态形变。设备故障或人为干扰等原因引起的阶跃需预先标注时间,ADECP增加相应参数并估计大小,阶跃时间是否正确标注间接影响同震和震后形变的识别和参数估计。如果GPS观测期间有多次地震发生,基于准确的发震时刻,ADECP可以逐次判断GPS时间序列是否受到地震的影响。针对每一次地震,先假定非目标地震的同震和震后形变存在,估计并改正(如果不存在,相应参数会很小);再利用综合检校法判断是否受目标地震的影响。以此类推,待所有地震判断完毕,再重新整体估计所有参数。

      综合检校法的有效性和GPS站时间序列的同震和震后形变信号强度直接相关,同震和震后形变越大,自动识别有效性越高,可信度也越高。而GPS站的同震和震后形变大小和震级、震中距、震源深度、地质构造等多种因素有关。本文实验的地震最小震级为Mw 5.6,地表GPS站有明显的同震和震后形变。同时,还检查了2000-2016年加州地区Mw 5.0至Mw 5.5共16个地震附近GPS站数据,发现GPS站较难反映Mw 5.5以下地震的同震和震后形变,主要与地表发生形变较小有关,还可能与测站和震中距离等因素有关。

      衰减常数τ的确定是选取卡方值最小即拟合度最好的一个值作为最优τ,是数学方法上的最优,却忽略了地球物理学上的意义。此外,参数计算时采取单一站各自估计,忽略了各站间同震和震后形变的空间相关性。实际上,断层破裂引起的同震形变和余滑或粘弹性松弛等原因引起的震后形变在空间上必然相关,如果考虑这种空间相关性,可以避免其他干扰,参数估计会更合理稳健,这也是本文将要改进的地方。

参考文献 (19)

目录

    /

    返回文章
    返回