利用近震及远震波形反演的北京mw=51级地震震源深度

利用近震及远震波形反演的北京mw=51级地震震源深度

1地震震源机制的确定—引言以中亚地区为主体的新构造运动非常活跃。作为一个特殊的区域，它在大陆动力学研究中发挥着非常重要的作用。首都圈地区是中国政治、经济和文化的中心。历史上有多次大地震，造成了巨大的经济损失和人员伤亡。华北克拉通的研究一直是我们研究的重点和热点，也得到了地震和地质构造的各种结果。吴京通过切割波的分离，获得了首都圈西北部地壳环境的各向异性。当湖南伟等人通过16048条p波第一次到达时，他们获得了北京天津唐区3维p波结构（经度112125；纬度3545）。从1973年1月1日到2007年7月1日，usgs提供了地震、主要断裂带、文安地震和近地震站的分布信息。圆圈代表地震的级别小于4，波形代表地震的级别大于4.0。不同颜色代表不同的地震深度。黑色带显示了中国北方的几个断裂带。历史上著名的大地震主要集中在这些大断裂带上。黑色五角星列出了2006-07-04年文安地震的位置。这表明地震没有位于断裂带，因此引起了我们的极大兴趣。黑色三角是本文用来反演和计算的近地震站的分布。该站覆盖的方位角良好，我们可以获得该地震的源机制，并提供强有力的数据支持和机会。地震震源机制的确定,对于地震本身的研究、地震的孕震机理的解释及震后应力的分布,具有十分重要的意义.为此,对于震源机制理论和方法的研究,一直是地震学研究的热点,并取得了大量重要成果.1980年,Helmbergeretal.利用地震P波(P&&Pnl)模拟,通过理论和实际波形的对比,采用格点搜索的方法,得到了研究地震的震源机制.随后Wallaceet.al通过对震源机制进行方位角、倾角、滑移角的参数化,并利用最小二乘反演结合互相关的方法,得到了地震的震源机制.DziewonskiandWoodhouse给出了同时考虑近震波形和全球波形P波初到确定地震震源机制的方法,并提出利用正态模态叠加的方法研究了地震的震源机制.L.E.Jones同时考虑直达P和地表多次反射的P波得到了地震的震源机制和破裂过程.随着矩张量方法的出现,利用长周期地震矩张量反演的方法被广泛应用[17,18,19,20,21,22,23,24,25,26,27].而其他一些反演方法,主要针对近震波形的直达波和Moho面的多次反射波[28,29,30,31,32,33,34,35,36].近年来,随着计算性能提高及大量高质量宽频地震数据的使用,利用地震波形反演震源机制成为可能.Zhao等通过分割波形记录为Pnl和Snl部分,分别赋予不同的权重,通过格点搜索的方法进行地震震源机制的反演.最近,利用远震波形和近震波形相结合的方法被广泛应用于地震震源机制、地震定位等的研究中.杨智娴等利用双差的方法,利用经验格林函数的方法得到了地震的震源参数.本文在处理震中距小于600km的近震波形时,采用了Zhao(1994)提出的CAP方法,在处理30°～90°远震波形的时候,我们主要考虑直达P波、源区地表反射的sP和pP波,来同时确定地震的震级、深度和震源机制的方位角、倾角、滑移角等物理量.并通过了与哈佛大学基于中心矩张量方法(CMT)得到的地震震源机制解,利用P波初动得到的由CCDSN提供的震源机制解(P-initialmethod)进行了比较,三者较为相似;进一步通过余震分布、理论波形和实际观测对比发现,通过综合考虑近震和远震波形能快速准确地反演得到地震的震源机制和震源深度.2近震与远震的反演方法2.1格点搜索方法基于研究区域的地壳速度模型,利用近震波形进行震源机制反演,其反演原理为:设μ(t)是台站记录到的去除仪器响应后的地震波形,s(t)是相应的理论计算出的波形,如等式(1)所示:sj(t)=Μ03∑i=1Aij(ϕ-θ,δ,λ)Gij(h,Δ,t),(1)sj(t)=M0∑i=13Aij(ϕ−θ,δ,λ)Gij(h,Δ,t),(1)其中,j=1,2,3对应于垂向、径向、切向分量,Gij是对应各个方向的格林函数,Aij是辐射的衰减系数,M0是地震的矩张量,ϕ和Δ地震的方位角和震中距.需要反演得到的地震深度(h)、方位角(θ)、倾角(δ)、滑移角(λ)可以通过解等式(2)得到:u(t)=s(t)‚(2)其中μ(t)和s(t)分别为实际地震记录和理论地震图.为此,我们可以根据格点搜索的方法,同时搜索可能的地震深度(h)、方位角(θ)、倾角(δ)、滑移角(λ),同时给定误差函数作为测量标准,从而得到最佳的搜索出的震源机制.对于近震和远震波形反演,我们采用不同方法计算格林函数和误差测量函数.对于近震情况,近震波形受地壳横向各向异性(散射体的存在)影响较为明显,导致不同的频率范围内的波形差别较大,在这种情况下,我们采用频率-波数方法(F-K)来计算各台站处的格林函数.另外,在波形反演过程中,反演结果容易受到波形中较强部分的影响,如果采用全波形反演,长周期的面波将对反演结果起主导地位.为了克服这一反演中存在的困难,我们在进行近震反演的时候,利用的是Zhao提出的Cut-And-Paste(CAP)方法.该方法把整个波形分为P波部分(Pnl)和面波部分(Sur),对两部分的三分量共5部分(Pnl不存在切向分量)给定不同的权重进行反演,这样可以充分考虑各部分波形对于震源机制的贡献.我们选择格点搜索的方法,并选取等式(3)所示的误差测量函数.e=∥u(t)-s(t)∥,(3)对于远震波形,地震波的图像主要受震源机制、源时间函数、台站附近的速度结构影响,为此在计算各台站下方的格林函数的时候,我们采用射线理论.考虑到不同台站震中距差别较大,记录到的波形存在数量级的差别,为此我们采用归一化的误差测量函数,如等式(4)所示.通过格点搜索方法,在适当范围内循环地震深度、方位角、倾角、滑移角,得到相对误差最小时的地震的震源机制和地震深度.e=∥u(t)-s(t)∥∥u(t)∥⋅∥s(t)∥,(4)3数据分析与反演结果3.1独立强度生长记录本文用到的数据包括两部分:近震数据主要来自于首都圈台网、国家台网提供的震中距小于600km的波形记录,台站信息如表1中左栏所示;远震数据主要来自IRIS提供的、震中距在30°～90°之间的、信噪比(SNR)较高的波形数据,台站信息如表1中右栏所示.3.2近震波形反演的震源机制一般说来,人工试验反射波的方法能够得到准确的炮点位置和台站分布资料,通过偏移的方法,可以得到较为清晰的速度剖面,本文在选择源区模型的时候,采用了反射波得到的结果.由于文安地区地壳速度结构变化较为明显,呈现出由东向西递增趋势,为了得到更为准确的格林函数,对于不同方位角上的台站,我们采用了不同的地壳速度模型,如表2所示.其中东部、南部、北部的台站采用表2给出的模型,沉积层厚度在2km,中地壳存在一个14.5km的低速层;西部区域台站采用表2给出的速度模型,沉积层厚度为1km,中地壳低速体厚度为7km.通过频率-波数域(F-K)的方法,计算了基于上述模型、在不同震中距台站的理论格林函数.并将实际数据截断为Pnl部分和Sur部分,为了去除噪声的影响,采用了对信号进行滤波后比较的方法.滤波范围分别为:对于Pnl,保留5～15s的波形信号,对于Sur部分则保留10s以上的信号.对于理论计算得到的波形,采用相同的滤波范围,利用格点搜索及互相关的方法,计算二者的互相关系数、振幅和到时差.利用近震CAP反演得到的地震震源机制解、理论波形和实际观测对比、误差值、台站投影到震源球上的位置如图2所示.图中灰色波形为计算出的理论波形,黑色波形为实际观测的波形.基于表2中的地壳速度模型,利用F-K方法计算得到的理论波形,在不同方位角台站上都能较好地与实际记录到的波形匹配,且其绝对振幅大小也较为相似.互相关系数多数都在60%以上,且波形的三分量及Pnl和Sur部分都有较好的对应关系.同时也应该看到,二者在到时上存在一些差距,造成这些偏差的原因,可能是由于我们在选择模型时,P、S波的速度值与真实值之间存在差别,同时一维模型不能考虑到北京地区地壳和上地幔的横向不均匀性和存在的各向异性.由近震波形反演得到的震源机制为:方位角210°、倾角80°、滑移角-150°,地震为典型的走滑型.我们进一步比较了误差随不同地震试验深度变化的关系,图3给出了在计算不同震源深度,利用网格搜索得到的震源机制最优解的结果,其中横轴为不同的地震深度,纵轴为该深度计算理论波形和实际观测的最小二乘误差值.可以看出:(1)反演得到的震源机制较为稳定,随震源深度的变化不是很明显.(2)最佳的震源深度出现在14.0km左右,但是深度范围在11.5km至15.5km的偏差值较小.为此,本文引入了远震波形尤其是深度震相(pP)进行震源机制的反演,该方法首先可以进一步验证本文近震波形反演的震源机制的结果,更为重要的是提高地震震源深度的分辨率.3.3远震台定位特征分析本文基于近震波形,利用CAP方法反演得到了地震的震源机制如图3.为了验证所得震源机制的可靠性,利用该机制计算在远震台站上的理论波形,与对应远震台站上的波形进行了比较.当震中距小于30°时,会出现P波的三重值现象,不利于波形的比较和识别;当震中距大于90°时,来自CMB的绕射波及其前驱波(precussor)也会对我们需要的直达P波、sP、pP波带来影响.为此,本文用到的远震波形的台站震中距在30°～90°之间,分布如图4所示.由于在反演过程中,需要P、sP、pP波的波形的信噪比高,震相清晰,我们从78个远震台站中,仅仅找到8个信噪比高、震相清晰的台站,台站分布如图4所示.黑色三角代表台站的分布,右边为对应的台站的名称,黑色五角星为文安地震的震中分布.通过两步格点搜索的方法,即搜索过程分两次完成.第一步,对参数化的方位角、倾角、滑移角进行全值搜索,给定较大的搜索步长,可以在不影响结果的情况下,节约时间.全值搜索的方位角的范围为0°到360°,倾角范围为0°到90°,滑移角范围为-180°到180°,步长为10°;第二步,在第一步最小值附近,变化参数范围,选择较小的搜索步长,从而得到更加准确的震源机制解.实际观测波形中存在较强的背景噪声干扰,本文通过频谱分析发现,噪声的主要频率范围6～10s之间.为此我们在比较远震波形的时候,采用了滤波比较的方法,本文采用的滤波范围为:带通滤波0.8～3Hz.如图5所示,本文将搜索得到的最佳的震源机制解,所对应的理论波形与实际的波形数据滤波,并按照绝对时间对齐,可见对于我们找到的8个台站,直达P波、sP、pP在到时、极性、振幅大小各个方面都较为一致.其中灰色的波形为计算得到的理论波形,黑色的波形为实际观测的数据波形.由于加入了远震波形,可以使我们反演得到的地震深度的分辨率大大提高.为了显示不同深度的地震波形在远震记录上的差别,我们进行了深度分辨率测试,如图6所示.本文计算了以文安地震为震中,以搜索出的最佳震源机制解为震源机制的地震,在台站COLA上记录到的理论波形.测试的地震深度从10～18km.通过理论和实际波形的比较,可以看出地震深度在14.0km时,直达P、sP、pP三个波形最为一致,而地震深度在13.0km或者15.0km时,上述三波形就存在较为明显的提前和延后性.可见远震波形中的深度震相对于地震的深度具有较高分辨率.4文安地震国内及周边区域地震采集对比我们比较了本文综合近震和远震波形反演的震源机制、Harvard的CMT、CCDSN的结果,如图7所示,可以看出我们的结果与Harvard提供的结果较为接近,都是接近纯走滑型震源机制,与CCDSN的震源机制存在一定的差别,但是地震的大体类型走滑型地震还是一致的,河北省地震局通过实地调查及震后形变、地震灾害分布等情况也得到了地震为沿北东向走滑的地震类型.与此同时,本文还分析了国家地震局台网中心提供的文安地震发生后文安及周边区域(经度114°～118°,纬度38°～41°)2个月内的余震目录,如图7所示.同时,可以看到本文通过波形反演得到的地震的震级为5.1,与当地实地地震调查的结果5.1级一致.需要强调的是,本文采用了直达P波和源区地表反射的P波(pP)进行震源机制解的反演,能够对地震深度进行较好的限制.对于远震地震覆盖(方位角)分布良好时,本文采用的利用P、sP、pP波形反演的方法是一种较为独立、可靠的确定地震震源深度的方法.本文反演得到的地震深度,与中国国家台网中心得到的较为一致,与哈佛大学给出的地震深度存在较大的差距.究其原因可能是,哈佛大学在求解地震震源机制解时,主要利用的是直达P波波形数据信息,没有利用近源地表反射的深度震相(pP、sP)信息.同时从本文给出的不同远震深度给出的理论波形与实际波形的比较图(图6)可以看出,地震的深度在14～15km之间,哈佛大学给出的22.7km,计算得到理论的pP、sP波的到时,远远大于实际观测数据的pP、sP波的到时.根据本文由波形方法反演出的震源深度可以看出,此次文安地震发生在中地壳的上部,该区域岩石层为花岗岩,在该温度和压力情况下会出现部分熔融和蠕变活动,造成中地壳上部应力的积累和大地震的发生.同时还可以看到,本文利用波形反演方法得到的地震震源机制与利用P波初到等方法得到的震源机制较为吻合,且具有对于地震深度更高的分辨率.5基于速度模型的地震采集本文通过国家台网、首都圈台网提供的近震波形数据,及IRIS提供的远震波形数据,反演得到了2006/07/04文安地震的震级为Mw=5.1,震源机制:走向角为210°、倾角为80°、滑移角为-150°,地震深度为14.0～15.0km.华北地区为非均匀多层地壳结构,上地壳主要为中、新生代沉积层,电阻率较低,而中地壳主要是花岗岩层,受温度和压力的影响会出现部分熔融和蠕变活动,造成应力的积累和大地震的发生,这可能是此次文安地震发生的机制之一.同时考虑近震及远震波形进行波形反演,为研究中等震级的地震研究提供了更为准确的方法,尤其是地震深度的确定.从图3也可以看出,当近震台站分布情况较为理想时,利用CAP方法得到的震源机制较为稳定,随深度的变化影响较小;而采用远震波形的sP、pP波,与直达P波的出射方向在深度方向相反,通过它们之间的到时差,可以很好地确定地震的深度.二者结合,可以得到较为理想的地震的震源机制和地震深度.与此同时,我们在计算过程中,采用的波形的绝对大小、振幅与地震释放的能量相关,所以能够同时得到较为准确的地震震级.本文在计算格林函数的时候,考虑到格林函数与物质本身的性质有关,对于不同方位角上的台站,根据不同的地质构造背景,给出了两种不同的速度模型,这样有助于得到更为真实的台站和震中之间的格林函数.速度模型对于震源机制反演结果的影响视情况而定,其中关键因素为反演时关注的波形的频率范围:(1)利用较低频的地震信号反演震源机制的时候,地震记录初动极性(polarity)的正负,决定着反演的震源机制的结果,对于小的速度结构的变化不是很敏感,此时震源机制对地球模型的依赖较小.(2)利用较高频的地