《场地响应及路径衰减联合反演方法软件实现及应用.doc》由会员分享，可在线阅读，更多相关《场地响应及路径衰减联合反演方法软件实现及应用.doc（15页珍藏版）》请在三一文库上搜索。
1、场地响应及路径衰减联合反演方法软件实现及应用朱新运 摘要：地震波衰减及场地响应研究是地震学研究领域重要的前沿课题。基于谱比法，联合反演求解场地响应及地震波衰减参数,这种算法理论上是可行的，据此，基于MATLAB开发平台开发研制了地震波衰减参数与台站场地响应联合反演软件。本文介绍软件设计的思路、理论原理，方法实现及资料处理过程，简述软件模块功能及使用方法。该软件计算过程稳定、可靠，界面设计友好，可用于日常地震资料的分析，也可作为地震科学研究的工具。关键词 场地响应 衰减参数 联合反演 软件 1 引言在对耦合了路径、场地、仪器、地震辐射的方向性等多种响应成分的地震记录（陈运泰等，2000）进行必要
2、限定的情况下计算波谱比，可以扣除耦合的几何扩散、仪器响应、地震辐射方向性，但是场地响应与非均匀性、非弹性路径衰减则难以简单地扣除的。为了获得较为准确的衰减系数，Chun等(1987)提出的“逆向双台谱比法”。该方法是计算Lg波衰减参数比较可靠的方法（Chung,2002），它通过两台站连线两端延长线位置上发生的成对地震的记录谱比计算衰减参数。由于这种方法在资料选择上要求比较苛刻，影响了该方法的推广应用。Chung等（2005）对“逆向双台谱比法”进行了改进，改进后的“逆向双台谱比法”不再严格要求地震对以及台站对处在一条直线上，只要求地震对及台站对之间的距离满足一定的条件即可。这一改进扩大了地震
3、资料的选择范围，但对资料的应用仍然受到一定的限制。场地响应描述的是近地表地层介质对地震地面运动的放大作用，研究台站场地响应（Bonilla et al, 1997; Riepl et al, 1998;Langston, 1979; Nath et al, 2000,2001）的方法很多，但基本上是通过路径衰减项简化求得场地响应的。根据Shearer（1999）的研究结果，场地响应值与近地表地层介质的介质密度及地震波速度有关，与地震波传播距离无关，而地震衰减参数是地震波路径衰减的描述,仅与距离有关，由于地震Lg波的衰减对地震辐射模式不敏感(Mitra et al，2006)，因而使用Lg波作为
4、基本资料进行相关研究时，台站方位可以忽略。据此，朱新运等（朱新运等，2007，ZHU Xin_Yun,et al,2007;2012;2013;2014）以地震Lg波台站记录谱为观测数据，在资料选择上，使用Ms5级以上地震，为保证通过谱比消去震源因子时不致引入较大误差，限定台站对地震方位角变化不大，同时为保证衰减参数稳定，限定地震波穿过构造相对单一的区域，利用路径衰减与台站震中距有关、而场地响应与震中距无关这一特点，通过联合反演将场地响应与路径衰减互相分开，得到了较好的结果。本文中作者对场地响应与路径衰减参数联合反演方法论证，同时介绍了在MATLAB平台下研制的路径衰减与场地响应联合反演软件，
5、及使用该软件研究了华北、华东等区域的地震波衰减及场地响应参数成果。软件设计构件合理、界面美观、计算功能完善、使用方便，可以作为地震波相关研究的实用工具。2 理论与方法台站记录的地震波的频谱可以表示为（Hasegawa,1985）： （1）式中，表示频率；表示震源因子项；表示震源辐射图型因子；表示震中距；表示仪器效应；表示场地响应；表示几何扩散因子， （2）是表示几何扩散特性的参数，对于面波波，体波；表示介质衰减项，它可以表示为： （3）式中，是衰减系数， （4）在式（4）中，表示波速度；。假定参与运算的地震数为，需要评估的台站数为，两个台站（第台与第台，）对同一个地震（第个地震）的记录谱的比可
6、表示为： （5）在上式中，对于同一地震，震源因子项之比可以消去；震源辐射图型因子与震源机制解有关，对较大地震选取资料时可考虑震中距离与台站间距离相比很大，以至于震源对台站之间方位角变化很小，则：，对较小地震则可看作点源，震源辐射因子影响不大，从而，在此情况下震源辐射图型因子项之比可以消去；仪器效应与是已知的；是台站记录谱；由式（2）给出。将式（3）代入式（5）整理后得： （6） 对式（6）两边取自然对数可以得： （7）式中， （8） (9) （10）令： （11） （12）则： （13）选取资料时考虑研究时间段内，在研究区即地震波射线穿过的区域内，没有发生大的构造活动，衰减参数保持恒定；台站场
7、地响应仅与频率有关，两台站对地震效应的幅度的比值仅与频率有关；在参与评估的台站中，对应任何频率总有一个台站的效应是理想的，即：对应该频率该台站效应为理想基岩场地响应，效应幅度为一个单位。首先以第一个台为参考台站，则： （14）式中，各角标意义与前述各式意义相同，写成矩阵形式，则为：(14)(15)在式（15）中， 参与评估台站数，分析频率个，频率代码（相关角标），则对应的有个未知数，不同频率下各台站场地响应对参考台站场地响应比值有，共个未知数，方程组的个数为。解方程（15）可得相应未知数，其中对应衰减部分为，利用（4）式，由可得到值；其余部分为场地响应比值数据，参考台站对自身在不同频率下的比值
8、均为1，对不同分析频率在各台站场地响应比值数据中逐一搜索比值最小值，以该最小值为相应频率在理想基岩场地响应的幅度单位，把该频率其余台站对应结果与相应最小值对比，得到该频率在相应台站的最终效应幅度，进而得到所有台站对所有分析频率的场地响应。3 软件实现软件设计包括文件操作、参数设置及获取、各个计算过程实现等模块,程序设计依据MATLAB平台特征，数据运算尽可能矩阵化，共有代码尽可能函数化。主窗口及窗口控件初始化直接使用GUI生成，程序控件一般属性，包括控件类别、控件在主窗口中的位置等通过GUI模板设计后进行相应赋值，主界面见图1。 图1 程序主界面（群速度高端太小，没有截取到合理数据段、终止计算
9、）3.1 参数设置软件对各关键控制参数均提供了数据接口，其中高、低频率及频率间隔数限定软件计算从低频到高频以等间隔分隔的各频率结果，考虑到在付氏变换过程中记录谱的突跳等影响数据稳定，软件设计平滑功能，在平滑窗键入相关参数，在数据计算过程中进行数据平滑。群速度参数用于控制地震波有效成分截取，该参数包括低端速度和高端速度，软件根据该参数截取对应时间窗的地震波形，为了有效截取波形，在准备地震文件时应提供较为准确的地震震源位置。几何扩散参数用于对波形进行几何扩散校正，一般面波校正参数为0.5，体波校正参数为1。截取有效波形后加窗进行付氏变换，软件提供汉宁窗、汉明窗等六种窗函数及窗宽大小编辑窗，可设计计
10、算过程中使用的具有相关属性的窗函数。3.2 结果显示及处理结果处理区包括了图像查看、结果存储等控件，图像查看包括过程图像查看和结果图像查看，过程图像查看相当于计算过程的回放；通过查看关键步骤产生数据图像及结果图像，可以对整个计算过程及计算结果是否合理有一个全面了解。在“波、谱图”区域，用户可通过键入地震记录条数并选择“波形图”或“波谱图”查看本条地震记录波形及在速度窗控制下波形被截取或波谱计算情况，在“衰减显示”下方选择衰减分向查看该分向衰减图形（图2）或选择“三向综合”查看所有三个分向计算的衰减结果，在“综合场地响应”下方选择分向名，绘图区绘制在同一坐标系下各台场地响应，在“单台场地响应”下
11、选择台站代码，绘图区绘制该台站场地响应结果。结果存储包括数据结果和图像结果存储两部分，“当前图像存储”按钮完成当前绘图区图像存储，为便于把多台场地响应结果绘制在同一页面，该软件专门设计了“多台场地响应”存储区，键入单页面纵、横参数，点击“存储”按钮完成多台站场地响应图件结果存储。图2、软件显示图形实例（1）.N_S向衰减结果4 算例：首都圈区域Lg波衰减及场地响应研究4.1 资料条件首都圈数字地震台网由107个地震台站组成，均为3分量数字地震仪，采样率为50Hz。其中44个宽频带台站的仪器为在20s20 Hz速度平坦的FBS一3A地震仪；58个短周期台站的仪器为在120 Hz速度平坦的FSS-
12、3地震仪，其余5个为在120s20 Hz之间平坦的CTS一1甚宽频带地震仪。甚宽带和宽频带地震仪大多分布在山区的地表台站，短周期地震仪大多分布在平原地区，5台安装在地表，其余均为安装在第四纪黄土层中的井下。上述仪器的性能如图3所示。图3 首都圈数字地震台网地震仪器性能首都圈地区地质构造复杂，尤其是山区和平原区地壳厚度、地质构造差异较大，东南为凹陷盆地，而西北部为山区。考虑到首都圈地区地质构造的复杂性，我们根据地震射线传播路径将评估衰减参数的资料分成相互独立的四组进行评估；第一组：西北部山区台站对同区域地震的记录；第二组：东南部凹陷盆地台站对同区域地震记录；第三组：西北部山区台站对东南部凹陷盆地
13、地震的记录；第四组：东南部凹陷盆地台站对西北部山区地震的记录。其中，第一、第二两组的地震射线分别在各自相同的区域传播，第三、四组均沿相反路径穿过地壳陡变区域。在2006年01月至2008年底记录的404次震级ML1.75.4的地震中，作者挑选了信噪比高的1000多条波形，涉及的地震台达81个，地震射线复盖了首都圈西北山区和东南盆地，地震、地震台及地震的射线分布情况如图2所示。图4 地震台、地震射线及地震分布图。 (a)隆起区台站对同区域地震的记录射线（粗虚线）东南部凹陷盆地台站对同区域地震记录（细实线）;(b)西北部山区台站对东南部凹陷盆地区域地震的记录射线（粗虚线），东南部凹陷盆地台站对西北
14、部山区地震记录（细实线）4.2 资料处理在本研究中，Lg波段的截取采用了固定速度窗法。利用清晰的地震波震相到时求得g波的速度上界3.65km/s，取g波速度下界为2.6km/s (Murphy et al,1982; Cappillo et al, 1985)，以该速度窗截取波形（图5），对截取的g波形段补零扩展至最近的2的整数次方倍，为防止谱泄漏，两端加2的余弦窗进行富利叶变换、然后做仪器校正、几何衰减校正，并以3个单位频率步长进行平滑（图5），使用数据段信号谱大于噪声谱（图6）。根据式（14）进行拟合（图7）。图5 Lg波段的截取（a）及波谱图（b）实例图6 信号与噪声谱比较图7 数据拟合
15、实例4.3 数据结果4.3.1 衰减系数从同一组地震射线中三个分向计算结果可以发现在6Hz以上总趋势发生拐弯现象，且在三个分向出现明显差别（图8）。这个结果与Shin等(1987)及Raúl等（2008）的研究结果一致，Shin等(1987)在研究Lg波衰减时发现在7Hz以上，Lg波的衰减发生了较大变化，并采取了一定的办法进行校正，可是并没有获得理想结果，Shin等(1987)对地震波衰减不一致的情况做了解释，认为Lg波的高频部分受到了Sn及Pn波的尾波干扰，因此建议研究Lg波衰减可在7Hz以下或者挑选震中距小于600km的地震的记录。Raúl等（2008）在对墨西哥东北索
16、诺拉州地震波衰减研究中也出现6或7Hz以上衰减不一致的情况，但没有进一步解释。为此，在本研究中我们仅计算16Hz之间衰减结果并拟合相关关系。图8 地震波1Hz至10Hz的衰减关系第一组射线包括西北部山区台站对同区域地震的86条记录图（No.1），第二组射线包括东南部凹陷盆地台站对同区域地震的1036条记录图（No.2），第三组射线包括西北部山区台站对东南部凹陷盆地地震的记录182条（No.3），第四组射线包括东南部凹陷盆地台站对西北部山区地震记录196条（No.4）。四组衰减关系中第一组最大，第四组最小，其余两组居中（图9仅画出垂直向）。图9 地震波衰减关系(No.1、No.2、No.3、No
17、.4分别代表第一、第二、第三、第四组数据)4.3.2 场地响应根据资料情况，首都圈区域东北地震资料较少，无法满足数据运算条件。在西北山区及东南凹陷盆地区域部分台站前后仪器响应不一致甚至差别较大，有些在资料选择时间段没有足够地震记录或有较多地震记录但在选择Lg波段时资料被剔除等原因，实际上共有81个台站得到稳定的场地响应。在参与评估的台站中，山区台站基底均为硬基底类地表台；平原区选择的台站中有8个台站为地表基岩台，l0个台站为井下泥岩台，其余台站均位于井下黄土沉积层中。基岩的场地响应多数在16 Hz频带内比较平坦，部分在频率高端有抬升（56Hz）。在低频段（12Hz），黄土沉积的场地响应普遍高于
18、硬基底的基岩场地（图10）。场地响应描述的是近地表地层介质对地震地面运动的放大作用，是衡量台站附近地下介质对地震波的阻尼影响程度的物理量。根据Shearer（Shearer，1999）有关的研究结果，松软的沉积场地对地震动的放大作用较大；坚硬的基岩场地则反之，本研究结果与此结论一致。图10 场地响应图。实线表示垂直向场地响应数据结果，点线代表南北向结果，虚线表示东西向结果5 讨论与结论使用该软件除了华北区域衰减及台站参数进行研究之外，也处理了华东区域地震波衰减及场地响应参数，并对这些结果与其他研究者在全球不同区域使用不同方法所获得的区域Lg波衰减参数，笔者对这些搜集并与之进行了比较（图11）。
19、本研究数据结果较好的反映了研究该区域地震活动水平。图11 不同区域Lg波衰减参数比较本研究以谱比法为基础，利用路径衰减与台站震中距有关、而场地响应与震中距无关这一特点，通过数据计算实现衰减参数项与场地响应项分离。软件设计以MATLAB为平台开发，设计上充分考虑开发平台特征，发挥其优势，以节省内存、运行速度快、流程可控、软件可读性强及操作简便为基本要求，功能上集数据预处理、参数计算、结果处理于一体，应用可视化操作，界面友好、功能完善，数据输入、输出及参数配置操作方便实用。参考文献陈运泰，吴忠良，王培德等.数字地震学北京：地震出版社：114朱新运，于俊谊，2008，基于MATLAB的小震震源参数计
20、算软件研制，2008，西北地震学报，30（6）：380384.朱新运，陈运泰，用Lg波资料反演场地响应与地震波衰减参数，2007，地震学报，29（6）：569580.朱新运，2013，数字地震波分析系统，浙江大学出版社。朱新运,2016,华北盆地Lg 波衰减及场地响应特征, 地球科学， 41(12):2109-2117Borcherdt R D. 1970. Effects of local geology on ground motion near San Francisco Bay. Bull Seism Soc Amer, 60:2961.Brune J N.1970.Tectonic
21、stress and the spectra of seismic shear waves from earthquakeJ.J Geophys Res,75(26)：4 9975 009Brune J N.1971.CorrectionFJ.J Geophys Res.76：5002Chun K-K, Gogdon F W, Richard J K, et al.1987. A novel technique for measuring Lg attenuation-results from eastern Canada between 1 to10. Bull Seism Soc Amer
22、,77:398419.Garccíf D, Singh S K, Herráiz M, et al.,2004.Inslab Earthquakes of Central Mexico: Source spectra, and stress drop.Bull Seism Soc Amer.,94:782802.Hartzell S H. 1992. Site response estimation from earthquake data. Bull Seism Soc Amer, 82:23082327.Hasegawa H S. 1985. Attenuation o
23、f Lg waves in the Canada Shield. Bull Seism Soc Amer, 75:15691582.Herrmann R B, Kijko A.1983.Modeling some empirical vertical component Lg relation, DARPA/AFGL Seismic Research Symposium,6-8 May,1985,U.S. Air Force Academy, Colorado Springs, ColoradoKing J L, Tucker B E. 1984, Observed variations of
24、 earthquakes motion across a sediment-filled valley. Bull Seism Soc Amer., 74:137152.Langston C A. 1979. Structure under Mount Rainier, Washington, inferred from teleseismic body waves. J Geop-hys Res. 84:47494762.Moya C A，Aguirre J，Irikura K2000Inversion of source parameters and site effects from s
25、trong ground motion records using genetic algorithms-JBull Seism Soc Amer,90(4)：977992Nath S K, Biswas N N, Dravinski M, et al. 2001. Determination of S-wave site response in Anchorage, Alaska in the 19 Hz frequency band, Pure Appl Geophys. 159:26732698.Nath S K, Sengupta P, Sengupta S, et al. 2000.
26、 Site response estimation using strong motion network: a step towards microzonation of Sikkim Himalayas. Seismology 2000, Curr Sci, 79:13161326.Tucker B E, King J L, Hatzfeld D. 1984. Observations of hard-rock effects. Bull Seism Soc Amer, 74:121136. Xie J, Mitchell B J. 1990. Attenuation of multiph
27、ase surface waves in the basin and range province. I. Lg and Lg coda. Geophys J Int, 102: 121137.ZHU Xin-yun, CHEN Yun-tai, An inversion of site response and Lg attenuation using Lg waveform,2007, Acta Seismologica Sinica，20（6）：605616.Zhu, X. Y. and Y. T. Chen, 2012, An Inversion of Lg-wave Attenuation and Site Response in the North China Region, Bull. Seismol. Soc. Am. 102(6), 2600-2610, doi:10.1785/0120110294.2012.Zhu XinYun, 2014, An Inversion of Lg-Wave Attenuation and Site Response from Seismic Spectral Ratios in the Eastern China Region, Bull. Seismol. Soc. Am 104(3):1389-1399.