吕作勇 马晓静 房立华
1)广东省地震局,广州市先烈中路81号大院1号 510070
2)中国地震局地震监测与减灾技术重点实验室,广州市先烈中路81号大院1号 510070
3)广东省地震预警与重大工程安全诊断重点实验室,广州市先烈中路81号大院1号 510070
4)中国地震局地球物理研究所,北京 100081
震级作为地震参数的基础数据,除了用来划分地震大小,还在地球物理学、地震工程学以及地震社会学等领域中被广泛使用。自从Richter(1935)针对南加利福尼亚地区、伍德-安德森地震仪给出ML震级的定义及其量规函数,科学家们在此基础上还发展了一系列震级标度(如面波震级MS、体波震级MB、Mb等)。ML震级标度简单、实用的特性使其至今仍是测定中小地震最为常见的震级标度。
在我国,1959年李善邦根据Richter震级标度建立了中国近震ML震级的测定公式(李善邦,1981)。1971年全国地震会议中,针对我国常用的短周期和中长周期地震仪确定了近震ML震级量规函数R1和R2,1977年将其列入我国《地震台站观测规范》(试行)(国家地震局,1978)。我国区域台网实现数字化后,为了保证ML震级标度的继承性和统一性,根据《测震台网运行管理细则》①中国地震局监测预报司,2011,测震台网运行管理细则,中震测函(98号),近震ML震级需仿真成伍德-安德森地震仪记录后,在位移记录上量取水平向最大振幅来计算震级。目前,我国各区域台网都使用由广东省地震局开发的地震分析软件JOPENS-MSDP进行日常地震分析工作,该软件使用的量规函数为R1,因此各区域台网测定的地震ML震级参数都是据R1计算获得的。大量的观测事实和研究结果表明(严尊国等,1983、1992),ML震级量规函数存在区域性特点,故在我国使用统一的量规函数计算ML震级,会使结果具有离散性和不稳定性。对此,国内已有许多研究结果发表(高国英等,1987;阎志德等,1989;樊星等,1990;薛志照,1992;赵明淳等,2005;孟晓琴等,2008;项月文等,2010;袁卫红等,2012),这在一定程度上改善了中国ML震级计算的精度。为了解决震级测定的问题,IASPEI成立了震级工作组,那些与南加利福尼亚地区地壳结构相似的区域可以使用工作组推荐的ML震级量规函数,否则需要重新厘定(IASPEI,2011)。
目前国内对ML震级量规函数的研究是在假定新确立的量规函数计算的ML震级与原测定值之间不存在系统偏差的基础上,采用均匀震级系统的思想,利用震级残差统计来修正原有量规函数,进而获取新的量规函数(陈培善等,1983;严尊国等,1983;陈继峰等,2013)。国外对ML震级的量规函数的研究则是在ML震级定义的基础上,仿真成伍德-安德森地震仪的波形来测量振幅,利用大量的地震观测数据(振幅值与震源距)代入定义公式建立方程组,通过求解方程组获得待定系数,从而获取研究区的量规函数(Langston et al,1998;Wu et al,2005;Keir et al,2006;Miao et al,2007;Askari et al,2009;Bobbio et al,2009;Nguyen et al,2011;Ottemöller et al,2013;Saunders et al,2013;Bagh et al,2014)。
随着近年来宽频带数字地震观测的飞速发展,刘瑞丰等(2011)强调中国应该采用IASPEI推荐的新震级标度,以发挥宽频带数字地震资料的优势,使中国震级测定与国际接轨。广东地区一直使用R1量规函数,其适用性和准确性需要重新研究。本文在ML震级原始定义的基础上,采用国际通用的方法,通过收集大量的地震观测资料来研究广东地区ML震级量规函数。在仿真伍德-安德森地震仪波形时,本文通过对利用频率域和时间域2种仿真技术所获得的数据分别进行研究,得到了相应的量规函数。
中国地震观测台网 “十一·五”计划项目实施后,广东数字地震台网的台站数从16个增加到44个,其中5个为国家数字测震台,39个为区域数字测震台,各台站都配备宽频带地震计,台站分布如图1所示。本文收集了2007~2013年的地震波形资料和相应的震相资料,首先将地震波形资料仿真为伍德-安德森地震仪波形,然后自动量取水平向最大振幅。
地震波形仿真方法有频率域和时间域2种,其中时间域仿真能实现实时ML震级测定,在地震预警和地震自动速报中应用较多。频率域仿真伍德-安德森地震仪波形的步骤为:对地震波形先去均值和倾斜,经傅氏变换后除以各台站仪器响应,然后乘以伍德-安德森地震仪响应,最后再变换为时间域。本文采用Kanamori等(1999)提出的方法在时间域仿真为伍德-安德森地震仪波形。采用的伍德-安德森地震仪的固有周期为0.8s,阻尼系数为0.7,放大倍数为2080(Uhrhammer et al,1990)。图2给出了新丰江地震台记录的一个地震的波形以及仿真为伍德-安德森地震仪波形图。振幅的量取是在参考震相资料的基础上,在S波震相之后的时间窗内进行量取,该时间窗的大小与震中距大小成正比。参考的震相资料为广东地震台网的观测报告,其震相拾取都是由有经验的地震编目人员来处理的,因此能够确保自动量取振幅的精确性。通过上述处理过程,为保证数据的质量,本文按以下原则对数据进行挑选:①地震震级ML≥1.5,②每个地震至少8个台站的数据,③仿真后测量的最大振幅值≥30μm。按上述原则挑选了1590个地震的波形进行分析(地震分布如图1)。频率域仿真处理后,最终挑选了14703个振幅数据;时间域仿真处理后,最终挑选了13809个振幅数据。振幅数据为两水平向最大峰峰值的均值的一半(单位:mm)。振幅数据与震级和震源距的分布情况如图3所示,大多数数据分布在震源距500km以内及ML4.0以下。
图1 广东地区测震台站与地震分布图
图2 LTK地震台记录的某一地震速度波形(a)、时域(b)和频域(c)仿真伍德-安德森地震仪波形
Richter(1935)定义ML震级的计算公式为
其中,A为伍德-安德森地震仪记录的最大振幅(mm);S为台站校正,-lgA0为量规函数,其取值使震中距100km处记录A为1μm的地震ML等于零。根据式(1)以及规定的零级地震基准,Bakun等(1984)将量规函数写为
图3 振幅数据的震级和震源距的分布图
其中,R为震源距(km);n、k为与地震波衰减有关的系数,与地壳结构有关,为待定系数。Hutton等(1987)在对比了几个不同地区的ML震级标度后,认为用震中距100km的伍德-安德森地震仪记录定义的零级地震在不同地区存在差异,建议以震源距17km的伍德-安德森地震仪记录10mm的地震为3级地震,并以此作为ML震级标度的基准。因此本文利用下式来计算广东地区的量规函数
其中,Mi为第i个地震的震级,Aij为第i个地震的第j个台站记录的振幅,Rij为第i个地震的第j个台站的震源距。
将通过上述挑选的地震资料代入式(3),建立线性方程组(Miao et al,2007)
利用Press等(1986)编写的SVD算法来求解方程(4),即可获得地震的ML值以及系数n、k,进而获得广东地区量规函数。
通过对单个台站计算的ML震级与所有台站ML震级的均值差值的统计分析,从而获得各台站的台站校正S,即利用下式来计算
其中,Sj为第j个台站的校正值,Mij为第i个地震第j个台站计算的震级,Mi为第i个地震震级,N为台站总数。
对选出的1590个地震的波形样本,分别采用频率域和时间域的仿真技术得到的数据,利用上述方法进行处理,获得了广东地区的量规函数为
其中,R为震源距(km)。式(6)来源于由频率域仿真获得的数据拟合,而式(7)由时间域仿真获得的数据拟合得到。由于在时间域仿真伍德-安德森地震仪波形是利用 Kanamori等(1999)设计的递归滤波器公式,导致最终拟合得到的量规函数式(6)与式(7)存在很小的差异(系数差别在10%以内)。量规函数式(6)常用于广东地震台网的地震分析和编目工作中,而式(7)主要是针对地震实时处理,能够应用于地震自动速报、地震超快速报等方面。
由(6)、式(7)带入ML定义公式得到计算ML震级的公式为
其中,A为仿真伍德-安德森地震仪波形的振幅最大峰峰值的一半(mm),R为震源距(km)。
图4给出了本文获得的量规函数与我国大陆区域量规函数R1和IASPEI推荐的南加州地区的量规函数(Hutton et al,1987)的对比。在将振幅均以 mm计且仪器放大倍数都为2080的基础上,从图4中可以看出,针对2种仿真技术获得的量规函数差别很小,而与R1和IASPEI推荐的量规函数相比,在震源距>250km时随距离的衰减明显偏小,这应是由区域构造差异决定的,也进一步证明不同区域不能使用同一个量规函数,各区域台网需要更符合本台网构造特征的量规函数。中国大陆区域量规函数R1与IASPEI推荐的量规函数也存在差别,在震源距约500km处存在拐点,<500km时R1值大于IASPEI推荐的量规函数值,反之则R1值小于ASPEI推荐的量规函数值。对0~150km震源距内量规函数图形进行详细分析可知,本文结果能够较好地符合原始ML震级标度的基准以及Hutton等(1987)建议的ML震级标度的基准,而R1量规函数在震源距50km内则明显偏小。近年来各区域台网在地震编目工作中也发现存在近台计算的ML震级明显偏小的现象,这可能是当震中距较小时R1不够精确造成的。
图4 不同区域量规函数对比
在计算量规函数的同时,还得到了参与计算的地震的ML值,图5给出了频率域和时间域仿真处理得到的ML震级与由R1计算的ML震级的关系图。2种资料处理方法获得的新ML震级可以等价(图5(c))。与原始资料的ML震级比较,新计算的ML震级在约小于2.0时明显偏大,大于4.0时则明显偏小,这可能与大于4级的地震资料偏少有关。
图5 频率域仿真获得的M L震级(a)、时间域仿真获得的M L震级(b)与由 R1计算的M L震级关系图(黑色直线为最佳拟合直线,其公式见图上方,虚线的斜率为1且过原点),(c)为频率域和时间域仿真获得的M L震级对比图。新M L为本研究量规函数计算结果,带上标R的为时域仿真结果,其余为频域仿真结果,旧M L为R1量规函数计算结果
在区域台网地震分析工作中,经常遇到由于各种原因导致水平向地震波形数据丢失,地震分析人员往往会直接在垂直向量取振幅来计算ML震级,这与 Richter(1935)的ML震级定义不一致。由于台基对水平向地震波形的影响往往比垂直向更大,对不同台基类型利用各台站垂直向波形来计算的ML震级一致性更好(Havskov et al,2010)。为验证利用垂直向最大振幅来计算ML震级的可行性,本文统计了地震垂直向与水平向最大振幅的关系(如图6),其振幅比约为0.95,这与 Alsaker等(1991)和 Saunders等(2013)的研究结果一致。因此,利用垂直向波形来计算ML震级是可行的,且对结果影响较小。
某次地震的ML震级是取每个台站计算的ML值的算术平均值,由于震源能量辐射、传播路径以及台基状况的差异,各台站计算的ML值往往与均值存在差别(Richter,1958)。本文通过统计分析各台站ML震级的残差来获取台站校正(公式(5))。图7分别给出了频率域和时间域仿真时,各台站的ML震级的残差分布图,各台站的残差分布都满足正态分布。表1汇总了各台站的台站校正,其空间分布见图8。由式(2)可知,负的台站校正表示台站记录的波形被放大了,反之则表示被缩小了。
分别对上述2种仿真技术处理的计算结果进行统计,两者获得的台站校正基本一致,差别较小。从图8可以发现,广东地区中部区域多数台站校正为明显正值,其它台站校正大多为负值,西部区域多数台站校正为明显正值,东部区域台站校正有正有负,校正量值偏小。台站校正的分布特征与各台站台基状况以及仪器安装情况可能有关。有明显台站校正的台站其ML残差直方图的顶点也明显偏离0(如图7)。
图6 地震垂直向与水平向最大振幅关系图
图7 广东地震台网44个台站的M L残差统计直方图
由图7及表1可知,每个台站的ML残差数都大于30且满足正态分布,因此本文利用Z检验(Z-test)来验证各台站的台站校正与0的差异性是否显著。假设各台站的校正值都为0,通过下式来计算各台站的Z检验分数
式中,S为台站校正,σ为台站校正的标准差,n为样本数。当|Z|≥1.96时,在5%的显著性水平上,台站校正与0差异明显,即台站校正不为0,反之则表示假设成立,即台站校正为0。表1给出了广东台网44个台站的Z检验分数。|Z|值较小的台站,其台站校正也较小,|Z|<1.96的台站则不需要台站校正。
表1 广东地震台网的台站校正
图8 台站校正分布图
针对不同的仿真处理方法,本文获得了2个计算ML震级的公式((8)、式(9))。为了检验它们的精度以及是否比用R1量规函数计算ML震级有改进,定义各台站计算的ML值与均值的差为ML震级残差,并对其进行深入分析。
图9给出了观测数据中ML震级残差与利用式(8)、(9)计算的ML残差分布图。从统计直方图来看,所有的残差分布都在0附近,但是它们的离散程度还是存在差别的。通过上文获得的量规函数计算的震级残差分布基本一致,且离散度都比原有量规函数计算的震级离散度要小。通过上文获得的新量规函数以及台站校正,利用式(8)、(9)来计算ML震级分别能够降低离散程度约30%和29%。在广东地震台网日常编目工作中,利用式(8)计算ML震级,而式(9)则可以应用于地震自动速报、地震超快速报工作中。
图9 M L残差统计分布图
本文挑选了广东地震台网2007~2013年记录的1590个地震样本,采用频率域和时间域2种仿真技术来实现宽频带速度波形仿真为伍德-安德森地震仪波形,并在仿真后的波形上自动量取地震的最大振幅,在ML震级定义公式以及Hutton等(1987)建议的ML震级标度基准的基础上,分别利用上述2种仿真技术获得的振幅数据以及地震的震源距来研究获取广东地区近震ML震级量规函数,其中公式(8)为频率域仿真伍德-安德森地震仪波形来测定ML震级,常用于广东地震台网中小地震的分析和编目工作中,而公式(9)则为时间域仿真伍德-安德森地震仪波形来测定ML震级,能应用于地震自动速报和地震超快速报等工作中。结果表明,本研究获得的广东地区ML震级量规函数以及台站校正的应用,有效提高了广东地震台网测定ML震级的精度,可使ML震级测定的离散程度降低约30%。
针对2种仿真技术获得的广东地区量规函数差别很小,但在震源距>250km时本文所得量规函数则明显小于我国区域量规函数R1以及美国南加州地区的量规函数,这表明广东地区用于计算ML震级的地震波(S波或Lg波)随距离的衰减要弱于我国大陆区域以及美国南加州地区,这也进一步验证了在中国大陆区域不宜使用统一的量规函数,而需要考虑区域性特点。通过对比新量规函数计算的ML震级与R1计算的ML震级,发现小于2级时新ML震级偏大,这与观测早期因地震台站较为稀疏而统计R1量规函数时使用的地震样本震级较大有关,而大于4级时则新ML震级偏小,这可能是由本研究ML4.0以上地震数据偏少造成的。由于在区域台网地震分析工作中经常需要利用垂直向的波形来计算ML震级,本文统计得到水平向与垂直向最大振幅比约为0.95,这对最终ML震级的测定影响较小。通过统计分析单个台站计算的ML震级与该地震所有台站计算的ML均值的残差,本文得到了广东地震台网44个台站的台站校正。同一台站利用2种仿真技术得到的台站校正基本一致。Z检验结果表明,44个台站中有38个台站需要加入台站校正,这与仪器安装和台基的状况有一定关系。
致谢:感谢匿名审稿专家提出的宝贵意见和建议。