ACTA Scientiarum Naturalium Universitatis Pekinensis
亚洲地区1918—1970年7级以上强震重新定位与震源机制测量及结果的不确定性分析
韩甲源 张献兵† 周仕勇†
北京大学理论与应用地球物理研究所, 北京大学地球与空间科学学院, 北京 100871; † 通信作者, E-mail: xbzhang@pku.edu.cn (张献兵), zsy@pku.edu.cn (周仕勇)
摘要 通过光学字符识别(OCR)软件, 整理世界地震台网1918—1970年的世界台站震相报告, 对该时期发生在亚洲及周边地区7级以上强震的震源位置及部分强震的震源机制(断层面解)进行测量。使用1918—1970年《国际地震资料汇编》(ISS)和EHB Bulletin的震相到时和到时差, 对该时段338个7级以上地震事件进行重新定位, 得到其中316个地震事件的可靠定位。使用1933—1970年ISS和EHB Bulletin的P波初动数据, 对该时段244个地震事件求解震源机制解, 得到其中209个地震事件的震源机制解。研究结果丰富了亚洲地区基于全球地震台网观测的1918—1970年仪器记录早期近现代强震目录震源机制资料, 并提供一种研究近现代地震的可行方法, 对开展其他地区或其他震级范围的仪器记录早期近现代地震研究具有启发作用。关键词 光学字符识别(OCR); 近现代地震; 重定位; 文字识别; 震源机制解
868
区。仅在过去一百年间, 中国发生的死亡人数超过5万的地震就有3次, 日本也发生过几次危害巨大的地震。对于亚洲地区的强震研究, 仅仅依靠近40多年的现代地震台网记录是远远不够的。因此, 基于 1918—1970年现代仪器测量的地震资料开展亚洲地区强震发震位置和震源机制的研究, 对于认识强震发震构造和开展地震危险性分析都具有十分重要的意义。
为了保障地震定位的精度, 需要有足够多的地震台站记录。考虑到20世纪初亚洲地区地震台网监测力度很弱, 本研究只对发生在那个时期的7级以上强震开展定位和震源机制测量。本文做如下设定: 1) 研究区范围为 26°3′E—169°40′W, 11°7′S— 77°43′N; 2) 研究时段设定为1918年1月1日—1970年 12 月 31 日; 3) 定位地震的震级范围为Ms≥7.0。本研究构建的含震源机制解的地震目录包括地震发生的时间、经度、纬度、震级和震源深度, 以及发震断层的走向、倾角和滑动角。
2 地震定位与震源机制测量2.1 地震资料概述
自1882年意大利地震学家Milne, Ewing和Gray研制的摆式地震仪在日本东京大学成功地记录并绘制一个地方震的地面振动图以来, 20世纪开始, 现代地震仪在世界各地陆续安装并投入观测, 全球范围的地震仪器记录越来越多地被收集和发布。1897
年, Milne开始收集世界各地60多个地震台站的记录, 并发布 Shide Circulars (BAASSC, 1900‒1912)。1913 年 Milne 逝世后, 牛津大学 Turner教授继续进行全球地震台站记录资料的收集和发布(BAASSC, 1913‒1917)。1918 年, 这份全球地震报告被更名为《国际地震资料汇编》(International Seismological Summery, 1918‒1963), 本研究的基础资料主要取自《国际地震资料汇编》报告, 该报告的资料来源于图 1所示全球基准地震台网记录(从图1可见, 欧洲、日本和美国的地震台网较为密集, 南半球的地震台网较为稀疏)。
EHB Bulletin 包含 1960—2008年的地震资料,是ISC (International Seismological Centre) Bulletin 的
[6]精简版。该目录使用Engdahl 等 的算法, 显著地改善之前ISS等地震目录的震中定位精度。同时,通过对地震波走时数据的残差校正, 反演得到全球速度结构。EHB Bulletin目录持续更新至2008 年,直到新的ISC定位算法[9]在2009年投入使用。
2.2研究方法
使用最新的地球平均速度结构模型Crust5.1 计算地震震相的理论走时, 能够有效地提高地震的定位精度。我们使用光学字符识别(optical character recognition, OCR)软件, 识别并整理1918—1970年期间ISS和EHB Bulletin中全球地震台记录的P波、S波及其他震相到时与到时差; 使用 Schweitzer[10]的 HYPOSAT程序, 对该时段的338个地震事件重
新定位; 使用俞春泉等[11]改进的格点尝试法, 基于1933—1970年期间ISS 的P波初动数据, 求解该时段244个地震事件的震源机制解。由于文档地震数据的扫描效果较差, 因此在OCR识别过程中, 选择合适的识别软件极为重要。在对多个软件的识别结果进行测试后, 本研究采用Adobe公司的ABBYY软件进行OCR识别。图2通过 1937年8月 11日发生在 6.8233°S, 116.1884°E 的地震数据, 展示 ABBYY的识别效果。ABBYY 对数字的识别准确率可达95%以上, 再经过简单的人工修正, 即能达到数据处理要求。即使有些到时数据的错误没有被纠正, 这些数据也会因为与正演结果的偏差太大而从数据集中剔除。虽然ABBYY对数字的识别准确率较高, 但对台站名称的识别效果较差, 原因可能是数字的结构相对简单, 可选择的范围较小。对于台站名称识别错误的情况, 本文尽量采用批量处理的形式给予纠正。这是由于台站名称的识别错误情况往往是相同的, 比如“b”被识别成“6”。这种相同的错误可以通过编程来批量地解决, 具体步骤如下: 1) 浏览各个台站的识别结果,对于单个台站, 列出常见的几种识别错误; 2) 通过建立“字典”数据结构, 将常见的错误结果与正确结果相对应。
在 OCR识别之后, 要将台站名称与定位程序
中的代号对应起来, 这样才能与反演程序中的台站代号相对应, 完成反演步骤。本文采用字典的数据结构进行操作, 自动将台站名称转化为台站代号。在转化的过程中, 可能遇到多个台站名称相同的情况, 此时就要舍弃数据或者精确度。我们对这个问题的解决方法如下: 如果一个台站名称对应的几个台站相距较远(大于 0.5°), 则把这一台站名称对应的数据全部舍弃; 如果一个台站名称对应的几个台站相距较近(小于 0.5°), 就将这一台站名称对应的数据全部归属其中一个台站。这样做的好处是, 既可以保证一定的精确度, 也不会因舍弃的数据过多对结果造成影响。接下来, 是数据的预处理步骤。这部分比较琐碎, 因为除需要对识别错误的数字进行更正外, 还需要对格式进行调整, 以期满足反演程序的输入要求。格式调整后, 对于数据中没有的参数, 除将时间的标准差设为1 s外, 其他参数(如反方位角和离源角等)则没有使用。经过以上处理后, 每个地震事件大约能得到100~500个有效的P, S 震相。
本文将IASP91与 CRUST5.1模型相结合, 作为震源参数的反演模型。相对于其他标准模型(如PREM, AK135), IASP91模型的优势在于对主要震相到时的计算较为准确。对于那些对局部速度结构敏感的震相(Pn, ss和 pp等地表反射震相), 我们采
用 CRUST5.1代替IASP91的地壳速度结构, 对于其他震相, 则采用 IASP91进行计算。对地震进行初次定位时, 先使用固定深度定位。本文使用0和 30 km作为初始深度, 进行两次定位。之所以选这两个深度, 是因为7级以上地震的震源深度中位数大约为30 km, 而以0 km作为初始震源深度可以更好地反演浅源地震的震中位置。HYPOSAT程序定位的时候进行多次迭代, 迭代后进行正演计算残差,去掉残差较大的几组数据, 使结果更具可靠性。确定震中位置后, 反演程序进行进一步的迭代, 更精确地确定震源深度, 并对震中的位置做微小的调整。我们在几次反演结果中选择水平方向不确定度最小的作为最终结果, 记录发震位置和发震时间等。由于本文研究的时段(1900—1970年)全球地震台网稀疏, 台间距为数百公里, 因此缺乏对震源深度的约束。将震源深度设为0, 10 或 30 km进行测试, 结果表明对定位结果(震中位置)影响不大。
2.3 定位结果
本文使用ISS和EHB Bulletin的震相到时数据和到时差, 借助 IASPEI91 模型和 CRUST5.1 模型,对 1918—1970年期间的338个地震事件进行定位,得到其中316个可靠性好、精度高的定位结果。我们从定位结果中可以发现, 震中位置的反演精度与台站记录数有关。对于台站记录数较多的地震事件
(震级大或发震年份晚), 震中位置的误差较小, 在0.1°之内; 对于台站记录数较少的地震事件(震级小或发震年份早), 震中位置的误差较大, 在 1°之内。
由于震源深度的反演结果存在较大误差, 且近现代强震的震相到时数据中鲜有pp或其他对震源深度敏感的震相, 因此对近现代强震震源深度的测定较为困难。由于本文将初始震源深度设为0和30 km, 进行两次深度计算, 取水平方向不确定度最小的作为最终结果, 因此反演结果中仅有部分地震事件得到非0或非 30 km且稳定的震源深度。其中,震源深度为0~15 km的地震有108 个, 15~50 km的有 148 个, 50~200 km的有43 个, 200 km以上的有17个。
在地震定位的精确度方面, 95%的置信区间内,经度的不确定度平均为0.199°, 纬度的不确定度平均为 0.168°; 除确定深度的定位结果外, 深度的不确定度平均为15.59 km; 定位的走时均方根残差为6.67 s, 走时差的均方根值为9.47 s。
从图 3 可以清楚地看出, 尽管由于台站稀疏,台网方位角覆盖度低, 导致板块边缘地震的定位误差大, 但本文确定的50 km以上的中‒深源地震绝大多数发生在太平洋或菲律宾板块向欧亚板块俯冲的岛弧一带, 清晰地揭示出这3个构造板块间的相互接触关系。我们注意到, 印度尼西亚北部苏维拉西
海附近地震的定位精度较差, 原因可能是该地区在本文研究时段(1900—1970年)内几乎没有地震台站(图1)。在中国大陆松潘‒甘孜块体内部也有一个地震的定位误差较大, 经查证, 此地震是发生于1923年3月24日的四川道孚地震, 由于当时全球地震台网建立不久, 数据较少, 数据质量较差, 因此定位结果不够精确。图3显示, 由于陆地上地震台站的方位角分布较好, 因此定位结果通常比发生在大洋或大洋边缘地区的地震好。
将本文定位结果与文献[3‒4]中相同地震的定位结果(二者都是基于现代地震仪记录的定位)进行比较, 可以发现二者基本上一致, 证明了本文结果的稳定性与可靠性。图4展示本文定位结果对《中国历史强震目录》[1]中 7级以上地震震中位置的修正, 由于《中国历史强震目录》主要参考现场灾害考察的文献档案, 因此对震中的确定有较大的人为因素导致的不确定性。
2.4 基于 P波初动的震源机制解
利用ISS 中 1933—1970年期间244个地震事件的P波初动数据, 使用俞春泉等[11]改进的格点尝试法, 求得209个地震事件的震源机制解(图5)。根据P波初动对断层面解的约束程度, 我们将解分为A, B, C三类。A类解为准确可靠的震源机制解, P波初动在震源球投影图上分布比较均匀, 对得到的两
个节面有较好的限制; B类解为较为准确可靠的震源机制解, 可靠性比A类解略差, P波初动在震源球投影图上分布比较合理, 基本上能够约束住两个节面; Ｃ类解的可靠性相对更差,P波初动在震源球投影图上分布不均匀, 不能很好地约束两个节面, 但基本上能判断出震源机制解的性质。对于C类解,其矛盾比相对较高, 且存在多组可能的断层面解。我们根据2.3节的定位结果, 推断是哪个构造带或断裂上的地震, 参照地震所在构造带或断裂的几何参数, 选取与其最匹配的解。209个震源机制的多解表达图及详细解在附录1中展示, 相应的震源机制节面参数信息见附录2。附录2中131个C类解是根据亚洲及邻区的地质资料[12‒14]确定的距图3所示地震最近的构造断层的走向和倾角, 并在格点尝试法程序自动给出的所有断层面解中选取与该断层走向最接近的一组作为最终的断层面解(黑粗线), 同时兼顾解的矛盾比, 即在断层面参数接近的情况下,以矛盾比更低作为确定最终解的原则。
为了研究100多年来亚洲地区构造运动的稳定性, 并论证本文测定的强震震源机制解的合理性,我们下载了哈佛大学有关本文研究区1976 年1月1日至 2018年 10月1日的7级以上强震CMT (central moment tensor)解(https://www.globalcmt.org)(图 6)。对比图5与图 6可以看到, 近半个世纪以来, 哈佛
大学测定的亚洲地区7级以上强震的震源机制解数目与我们挖掘历史文档记录测定的近现代7级以上强震的震源机制解的数目大体上相当, 我们测定的强震震源机制解对中国大陆部分地区(如南北地震带)以及天山构造带等重要活动构造带的强震资料有明显的补充。还可以看到, 太平洋或菲律宾板块向欧亚板块俯冲的岛弧一带地震的震源机制几乎都是逆冲型, 反映太平洋和菲律宾板块向欧亚板块俯冲的动力学过程; 天山构造带发生的强震破裂机制也是逆冲型, 反映天山南缘塔里木盆地等刚性地块向天山造山带下插的动力学过程[15]。
3 结束语
本研究在前人工作的基础上, 通过地质资料、
国内外地震台记录和近现代地震目录, 利用震相走时定位和P波初动求解震源机制解, 得到亚洲地区1918—1970年间316个7级以上历史地震的重定位地震目录以及1933—1970年间209个历史地震包含震源机制解的近现代地震目录, 对亚洲地区近现代地震资料进行了较好的补充和完善。
由于震相到时存在误差, 定位结果中某些近现代地震的发震位置及震源深度误差较大, 尤其是发震时间早及台站记录少的地震, 发震位置及震源深度存在更大的不确定性。此外, 对于基于P波初动求解的断层面解, 如果P波初动的台站记录数据不够多, 数据点在震源球上的投影覆盖范围就不够大,导致即使使用更精确的反演方法, 也难以约束断层面解。由于历史原因造成的观测限制, 我们尚未找到更好的对缺少观测地震事件的震源参数进行精确测量的方法, 本文结果仅能对缺乏强震台站记录地区开展强震发震构造等相关研究提供可用的参考资料, 对于结果的不确定性以及对相关研究的可能影响, 建议参考图3(定位结果的不确定性)进行详细的分析。
现代地震台网记录的历史很短, 地震地质、工程地震及地震危险性分析的相关研究常陷入地震目录太少或存在不确定性的困境。因此, 充分挖掘历史强震的信息, 并对其中部分错误进行修正, 是一项非常有意义的工作。本研究的另一个意义在于提供一种研究历史地震的可行方法, 可以对其他地区或其他震级范围的历史地震研究起到启发作用。
参考文献
[1] 国家地震局震害防御司. 中国历史强震目录. 北京:地震出版社, 1995 [2] 中国地震局震害防御司. 中国近代地震目录. 北京:中国科学技术出版社, 1999 [3] 贾科, 周仕勇. 西南地区 1900—1970年历史地震震源参数推断及结果的不确定性分析. 地球物理学报, 2012, 55(9): 2948‒2962 [4] 高彬, 贾科, 周仕勇. 华北地区 1900—1970 年 5级以上地震发震位置与震源机制研究. 地球物理学报, 2016, 59(11): 4089‒4099 [5] 周仕勇, 许忠淮. 现代地震学教程. 北京: 北京大学出版社, 2010 [6] Engdahl E R, Hilst R V D, Buland R. Global teleseismic earthquake relocation with improved travel times and procedures for depth determination. Bull Seismol Soc Am, 1998, 88(3): 722‒743 [7] Villaseñor A, Engdahl E R. A digital hypocenter catalog for the international seismological summary. Seismological Research Letters, 2005, 76(5): 554‒559 [8] 许力生, 陈运泰. 用数字化宽频带波形资料反演共和地震的震源参数. 地震学报, 1997, 19(2): 113‒128 [9] Giacomo D D, Engdahl B, Bondar I, et al. A new reference global instrumental earthquake catalogue (1900‒2009). AGU Fall Meeting Abstracts, 2011, 158 (1): 19‒26 [10] Schweitzer J. How can the ISC location procedures be improved?. Physics of the Earth & Planetary Interiors, 2006, 158(1): 19‒26 [11] 俞春泉, 陶开, 崔效锋, 等. 用格点尝试法求解P波初动震源机制解及解的质量评价. 地球物理学报, 2009, 52(5): 1402‒1411 [12] Yin A. Cenozoic tectonic evolution of Asia: a preliminary synthesis. Tectonophysics, 2010, 488(1): 293‒325 [13] Holm R J, Rosenbaum G, Richards S W. Post 8 Ma reconstruction of Papua New Guinea and Solomon Islands: microplate tectonics in a convergent plate boundary setting. Earth-science Reviews, 2016, 156: 66‒81 [14] Wiens D A, Demets C, Gordon R G, et al. A diffuse plate boundary model for Indian Ocean tectonics. Geophysical Research Letters, 2013, 12(7): 429‒432 [15] 周仕勇, 许忠淮, 陈晓非. 伽师强震群震源特征及震源机制力学成因分析. 地球物理学报, 2001, 44 (5): 654‒662