王同军 黄春梅 舒涛 戴应洪 甘利涛
1)重庆市地震局,重庆市渝北区红黄路339号 401147
2)四川省地震局,成都 610041
2013年4月之前,四川测震台网使用SeisDPS地震分析软件进行地震编目工作,地震日常定位方法为区域搜索法。2013年4月之后,编目工作使用的是Jopens5.2软件的地震分析模块Msdp,使用的速度模型是川滇三维走时表。该走时表针对川滇地区地震台站密度较低、壳幔速度结构横向变化剧烈的特点,在震相走时计算和地震定位时考虑了地球扁率、台站高程和地形起伏等的影响,其三维速度模型的震相理论到时显示功能,有助于震相分析人员在复杂介质结构区域正确识别震相。但因分层较多,该走时表在速报和编目时处理速度相对较慢。
2015年,四川测震台网人员参与了“全国区域一维速度模型”课题的研究工作,根据四川发生的地震和爆破震相等资料,通过数据拟合、波速稳定性分析、折合走时曲线拟合、Hyposat批处理搜索和结果检验等方法,建立了四川地区一维地壳速度模型(以下简称 2015模型),并在四川台网的地震速报和编目中选择应用。该模型从稳定性、定位精度以及PTD方法验证等方面都证明,在相同定位方法的前提下,由该模型定位的震源深度更接近于实际情况。
因四川测震台网接入了周边邻省共享台站数据,这使得九寨沟MS7.0地震成为网内地震。对于网内地震来说,在不同速度模型下,各种定位方法结果的经纬度、震级等要素基本变化不大,只是有的模型不适合当地,台站的震相到时残差可能会相对较大。但是不同模型下,很重要的地震参数——震源深度会变化较大。如果使用的速度模型浅于当地真实地壳结构,一些定位方法则定不出深度,显示为0。本文研究旨在利用四川2015模型,选择Msdp嵌入的几种常用定位方法,对九寨沟MS7.0地震进行重新定位,比较震源深度结果,并用PTD方法来验证合适的震源深度,以期得到较为可靠的深度结果。
四川测震台网工作人员选取四川测震台网2009年1月~2014年12月间,记录台站数大于30个(包括邻省台站)、震级ML≥1.0的地震事件4529次,震相共计211189个。利用波速拟合、走时折合等处理方法确定建立速度模型所需的各项参数的变化范围;设定各个参数的变化步长,利用 Hyposat地震定位程序,对选取的地震数据进行批处理,再根据处理结果,选取残差较小且上地壳和下地壳的深度符合前人研究结果的参数,初步建立新的速度模型;利用批量定位的残差、PTD方法定位的深度、人工地震数据、典型地震事件等对新建模型进行验证,分析新模型的合理性。由此确定的最优一维速度模型如表1所示。
表1 四川2015地壳速度模型
2013年 10月 21日 01:00和 01:10,芦山MS7.0地震科考项目进行了2次人工爆破,第 1次在芦山县双石镇(30.27°N,102.90°E),药量为 1.5t,第 2次 在 宝 兴 县 永 富 乡 (30.67°N,102.56°E),药量为 2.0t。根据四川测震台网记录到的地震波形,在Hyposat定位方法下,比较四川2015模型和四川分区速度模型(AA、BB模型)的定位结果,结果显示,AA模型定位宝兴永富人工爆破深度为35.4km,明显不符合实际深度;BB模型定位2次人工爆破深度均为0km,但震中误差相对较大;使用四川2015模型的震源深度均为0,震中误差最小,更接近真实情况。
利用PTD方法对各个模型进行了验证,结果表明,使用四川2015模型时的震相分布更服从于高斯正态分布,表明这一速度模型更为符合实际情况。使用四川2015模型与分区速度模型对2009年1月~2014年10月MS≥3.0的319次地震进行Hyposat批处理,分区速度模型(AA模型)的定位残差为0.52,分区速度模型(BB模型)的为0.71,四川2015模型的为0.48;水平误差和深度误差方面,相较于编目数据库结果而言,四川2015模型的水平误差和深度误差都很小,且趋势稳定。2015模型与分区速度模型(AA模型)的水平误差相差不大,但AA模型的深度误差出现了较多的异常值。综上所述,四川2015模型的定位结果相对较好。
2017年8月8日九寨沟发生MS7.0地震,四川测震台网可接收数据的全部108个台站(包括周边邻省共享台站)均清晰地记录了此次地震的波形,距震中最近的台站为四川九寨沟台(JZG),震中距约36km,最远的为云南大姚台(DAY),约864km。台网所有台站定位的空隙角为64°,该地震属于台站包围较好的网内地震。
九寨沟地震震中烈度为Ⅸ度,但在震中对地表的破坏性相对较小,远场区域震感相对较强烈,因此,该地震震源深度是非常关键的参数。在地震定位时,垂直方向没有台站,无法用垂直方向台站震相走时残差来约束震源深度(宋秀青等,2016),并且也没有震中距小于震源深度的台站(假定震源深度不超过33km),最近的九寨沟台(JZG)震中距为36km。所以,在各种速度模型下的各种定位方法中,震源深度均不是确定值。
以下为几个研究机构给出的初判震源深度:①中国地震台网中心利用国家台站波形数据,通过面波振幅谱波陷特征,初判此次地震震源深度为12km;②北京大学、中国地震局地球物理研究所给出的初判震源深度为10km,余震的震源深度为6~22km;③中国地震局地震预测研究所大震应急组初判震源深度为5km。
速报震源深度方面,四川测震台网和中国地震台网中心速报震源深度均为20km;甘肃测震台网13km;青海测震台网15km;陕西测震台网10km。四川测震台网中心编目的震源深度为23km。
统计了四川测震台网提交到中国地震台网中心的2010年以来九寨沟地区地震的编目深度结果发现,震源深度为1~5km的地震占14%;6~10km的占11%;11~15km的占41%;16~20km的占 27%;>20km的占 7%。从统计上看,该地区震源深度的优势范围为 11~15km。
以上各个深度结果差别较大,若要确定相对准确的震源深度,就需要适合当地的速度模型和定位程序,特别是合适的速度模型对于震源深度的准确测定有着非常重要的作用。
地震定位方法是产生更精确地震资料的基础,定位方法的正确选取,是产出可靠地震资料的保障,基于此才可为相关各领域的研究及应用提供更好的数据服务。
下面介绍MSDP地震分析软件下嵌入的几种定位方法。
(1)单纯型法是在n维空间中,用n+1个顶点构成1个多面体。依据单纯型运算规则,计算各个顶点的函数值,然后进行比较,确定顶点的优劣;然后计算新点,用好的顶点代替,不断地改变顶点,使单纯型朝着目标函数最小的方向移动,最终获得准确解(陈贵美等,2009)。
(2)HYP2000定位程序采用传统Geiger法的基本思路:将走时在初值附近做泰勒展开,在此基础上建立观测方程组,然后求解。在具体计算方法上,即先将观测方程组降维,不化成正规方程组,直接用奇异值分解最小二乘法方程组。HYP2000还可采用分区水平分层速度模型来为每个台站制定不同的速度模型。
(3)Hyposat计算原理同HYP2000,在定位初始化中可以设置初始深度、初始深度误差等,深度计算还可以采用固定和反演同时参与计算。本次计算设置的初始深度为15km,初始深度误差为15km。
(4)川滇3D(LOC3D)使用的速度模型为地震行业科技专项“川滇地区地震走时表编制”的研究成果“川滇地区三维速度模型”,使用的定位方法为单纯型法。
(5)利用Jopens5.2的Msdp定位程序,配置单纯型法、Hyposat方法、HYP2000方法的速度模型为四川2015模型,在将波形放大为1s的软件窗口下,重新精确标注了109个台站的12个初至Pg震相、97个初至Pn震相,几种定位方法的结果见表2,残差见图1。
表2 几种定位方法定位结果和各台网正式编目结果
图1 各种定位方法下台站初至波残差
由图1可见,九寨沟MS7.0地震因震级较大,发震时有一个破裂过程,震源对于近台来说,已经不能做为一个点源,这导致各个台站的初至波震相残差在近台(0~200km)相对较大,而在200~750km,残差相对较稳定,远台(大于750km)由于地震波已经传至其他地区,速度模型已经有了很大的变化,也导致残差相对较大。
台站震相走时残差反映了理论走时与工作人员标注的实际走时之间的差别,其系统偏差主要反映了理论走时与实际观测的不一致性。这种系统偏差包含了某一射线路径上实际地球结构与标准对称地球模型的差异(薛峰等,1998)。一般来说,走时残差是定位误差、震源区域介质以及浅层的不均匀性、台站附近介质速度异常及地震波的频散效应等多种因素引起的走时误差的积累。
走时系统偏差是将各个台站Pg或Pn震相残差取绝对值后,再计算所有残差的算术平均值。与平均误差相比,走时系统偏差由于离差被绝对值化,不会出现正负相抵消的情况,更好地反映了预测值误差的实际情况。
九寨沟MS7.0地震定位结果中的发震时刻、经纬度、震级等基本要素,在各个测震台网中心编目结果中差别不大,在此不作讨论。差别较大的是震源深度,为5~23km,没有统一的深度值。
由表2、图1还可见,几种定位方法计算出的系统偏差和深度差别较大,HYP2000定位方法深度最浅,系统偏差最大;Hyposat方法系统偏差最小,深度为7km;单纯法和川滇3D法系统偏差相近,深度相差了4km。
距离最近的九寨沟台(JZG)震中距约 36km,Sg-Pg时间差约 4.98s,按照虚波速度为8km/s的波速,则震源距为39.84km,根据震中距、震源距以及震源深度组成的勾股定理,可得震源深度约为15km。据此,四川2015模型下,单纯法的12km和川滇3D的16km深度结果较为可靠;各个测震台网的编目结果见表2,由表2可见,中国地震台网中心的10km、甘肃测震台网的13km,较为可靠。分析认为,这2个测震台网地震定位所使用的台站均能较好地包围震中,结果相对较合理。
张晁军等(2010)根据震源深度计算公式,在震中距、速度模型一定的情况下,得出震相走时残差对深度误差有较大影响的结论。在距震中100km处的台站,当震相走时残差为0.05s时,震源深度的误差达到±2.66km。震源深度误差也随着震中距的增大而增加,这也意味着确定震源深度必须要有近台参与定位。
杨成荣(2014)对2010年1月18日19时35分新疆呼图壁ML4.0地震,利用相同的速度模型(新疆3400走时表)、不同的定位方法进行定位,除了单纯法的深度结果为5km(经震中距、震源距以及震源深度组成的勾股定理验证为结果不可取)以外,其他方法得到的震源深度较一致,均为23km左右。
由以上的研究结果可以看出,只要选取适合当地的速度模型,排除明显错误的震源深度结果,则在理论范围内的深度结果均较一致。因此,速度模型较定位方法对深度的影响稍大。
综上所述认为,单纯型法的12km左右至川滇3D方法的16km左右的深度范围,为此次九寨沟MS7.0地震的初步定位深度范围。
PTD深度定位程序是根据朱元清等(1990)提出的利用Pn、Pg震相测定震源深度的方法而研制的程序。程序的实现过程中,在指定的深度范围内直接搜索震源深度,在该深度下,初至波理论计算值与实测值之差最小。该方法利用了初至震相Pg、Pn相对较清晰、信噪比较高、标注误差小的特点,以及莫霍界面的波速比较容易确定且比较稳定,并结合直达Pg波和折射Pn波的离源射线反向容易辨认的特征,能比较敏感地反映震源深度的变化。理论计算结果显示,当深度变化5km时,对震中距为100km的地震台,其到时差的变化可以达到0.70s。如果只用Pg波的走时,当震源深度变化10km时,Pg波的走时差只改变了约0.07s(宋秀青等,2014)。PTD方法的优点在于,采用不同台站的初至Pn、Pg到时差计算震源深度,可以避开发震时刻引起的误差,且参与计算的Pn-Pg数据对较多,如本研究采用了797个样本对,提高了区间的置信度。
针对上、下地壳的震源深度有多解的情况,对于程序搜索到的最小Pn-Pg走时差对应的深度值,若其走时比深度上下变化1km导致的局部走时梯度大3.3倍时,则认为该组数据无法搜到合适的震源深度,舍弃该深度;若在该深度处计算的直达波走时与对应的走时差大于此深度下变化1km导致的局部走时梯度的3.3倍时,则认为该深度值不符合直达波的走时要求,该组数据无法搜到合适的震源深度,舍弃该深度。这样的条件约束能有助于剔除多解情况下的虚值。
在PTD程序配置好四川2015模型后,其他的计算参数设置如下:当理论到时与实测到时之差大于0.5s时,舍弃该样本点;大于1km的深度定义为有效深度,参与最后的深度平均值计算;为了保证了Pg、Pn震相的准确性,震中距≤150km的Pg震相以及震中距≥150km的Pn震相才参与计算深度,不指定震源深度在哪一层界面内,在地壳厚度内搜索。
设置计算参数后,该PTD计算软件使用的数据为Hyposat定位方法下的结果文件数据,共97个Pn震相、12个 Pg震相,可用 Pn-Pg数据对 797个,Hyposat定位的震源深度为6.7km。图2为四川2015模型PTD定位方法震源深度结果。图2中标注的深度(D)后面的第2个值为PTD方法重新计算的深度结果,为12.3km。
图2 四川2015模型PTD定位方法震源深度结果
当理论深度为12km时,如果Pg、Pn的到时随机误差分别增加0.5s,则深度误差幅度为±2km(宋秀青等,2014)。计算结果显示,此次九寨沟MS7.0地震震源深度为12.3km,取整为12km,则误差范围为±2km,这些震相到时组合的计算结果均为独立的深度结果。由于理论模型与实际应用的地壳模型间存在差别,所以,在取得大量单个深度结果的情况下,深度结果能满足统计规律,结果也符合高斯分布,并且这797个独立的震源深度的集中度最大值与高斯分布极值耦合的最好,因此认为,本次计算结果(12±2)km为九寨沟MS7.0地震震源深度。
本文比较了各省级测震台网提交到中国地震台网中心的地震编目结果后发现,四川、甘肃测震台网是距震中最近的台网,其所用台站能较好地包围震中,而其他台网台站分布单边,地震波传播路径较远且传播路径上介质速度结构复杂,故PTD方法无法配置速度模型,所以只从中国地震台网中心调取了这2个测震台网的正式编目结果,取其观测报告中的Pn、Pg震相,分别配置当地速度模型(四川2015模型已在上文作了说明,甘肃2015模型为vPg=6.12km/s,vPb=6.60km/s,vPn=8.13km/s,vP/vS=1.72,康拉德界面深度 21km,莫霍界面深度52km),用PTD方法重新计算了震源深度,使用四川测震台网编目震相,计51个Pn震相、50个Pg震相,可用Pn-Pg数据对为412个,计算结果为 13.7km(图3(a));使用甘肃测震台网编目震相,计98个Pn震相、26个Pg震相,可用Pn-Pg数据对为1156个,计算结果为11.2km(图3(b))。这2个深度结果相差不大,深度的集中度最大值与高斯分布极值耦合的相当好,并与图2的12.3km非常接近,这3个深度取平均值为12.4km。故认为震源深度为(12±2)km较为合适。
图3 编目震相PTD定位方法结果——震源深度高斯分布
九寨沟地区地处青藏高原向四川盆地陡跌的两大地貌单元的过渡地带,该地区活动断裂发育,地震活动频繁,是我国的地震活跃区之一,历史上曾发生多次有感破坏性地震。该区域处于松潘-平武地震活动带内相对稳定区,到目前为止,此次九寨沟地震直接破坏作用不太大,但由它诱发的泥石流、滑坡、山崩等,则会对九寨沟自然景观造成一定影响(郭建强等,2000)。根据现场调查,九寨沟MS7.0地震的直接破坏作用不大,这也从侧面说明了该地震有一定的深度,震源深度相对小于10km的可能性较小。
综上所述,本文参照四川测震台网2010年以来九寨沟地区地震的深度统计优势范围,对比在四川2015模型下单纯法定位和川滇3D方法定位的残差和深度,以及按九寨沟台(JZG)计算震源深度的范围,初步判定震源深度范围,再根据四川、甘肃测震台网编目的震相数据,用PTD计算结果验证,最终判断此次九寨沟地震的震源深度为(12±2)km较为合适。
致谢:感谢审稿专家在论文修改中提出的宝贵意见,给出的修改建议丰富了论文的内容,提高了论文的水平;感谢“全国区域一维速度模型”课题组提供PTD深度定位程序。