2021, Vol. 64
2. 地质灾害防治与地质环境保护国家重点实验室(成都理工大学), 成都 610059
2. State Key Laboratory of Geohazard Prevention and Geoenviroment Protection, (Chengdu University of Technology), Chengdu 610059, China
早在20世纪初，学术界就已经发现小地震的发生频率远远高于大地震(Utsu, 1999).在当代，该规律已总结为著名的古登堡-理查德关系(GR关系)：log10N=a-bM，其中N表示震级大于M的地震的数量(Shearer, 2009).b通常在0.6到1.1之间(Utsu, 2002).b如果取其典型值1(Crampin and Gao, 2015)，则意味着地震震级减小1级，数量增至10倍.因此，识别小地震可以获得更加完备的地震目录.并且，对小地震的研究可以用于确定断层面的几何结构(Shearer, 2002; Yang et al., 2009)和理解大地震的诱发与成核机制(Dieterich, 1994; Miller et al., 2004).然而，小地震的能量弱，用常规方法识别有难度.因此，地震学迫切需要搜索小地震的可靠技术.
匹配滤波(matched filter)，或者说模板匹配(template matching)是检测小地震的一类方法(Van Trees, 1968; Whalen, 1971).在Gibbons和Ringdal(2006)将其用于密集台阵后，受到中外地震学界关注(李璐等, 2017).匹配滤波已经成功应用于多种研究领域，比如前、余震序列(Peng and Zhao, 2009; 谭毅培等, 2014; Kato and Nakagawa, 2014; 侯金欣和王宝善, 2017; Meng et al., 2018; Sianipar, 2020)、深源地震(姜金钟等, 2019)、诱发地震(induced earthquakes)(Yang et al., 2017; Cochran et al., 2018; Skoumal et al., 2020)、触发地震(triggered earthquakes)(Shelly et al., 2007; Aiken and Peng, 2014; Yao et al., 2015)、重复地震(Schmittbuhl et al., 2016; Yao et al., 2017; Uchida and Bürgmann, 2019)、塌陷地震(杨慧等, 2018)、滑坡事件(盛敏汉等, 2018)、核试验监测(Gibbons and Ringdal, 2012; Yao et al., 2018)、颤震(tremors)(Shelly et al., 2007)、火山相关的地震活动(Zhang and Wen, 2015a)、主动源事件(唐杰等, 2008)和声发射事件(侯金欣等, 2020).
匹配滤波方法是通过在连续地震观测记录中寻找与模板波形相似的信号来检测地震的.多数匹配滤波方法用滑动互相关评价信号间的相似程度(Gibbons and Ringdal, 2006).匹配滤波方法的最一般流程如下：首先，从地震目录中挑选有良好波形记录的已知地震事件用作模板事件.模板波形与各台站各分量的连续记录做滑动互相关(Yang et al., 2009)，得到各分量的互相关波形.然后，这些互相关波形叠加起来得到叠加互相关波形(Shelly et al., 2007)，或者进一步除以分量总数得到平均互相关波形(Peng and Zhao, 2009).当叠加互相关波形或者平均互相关波形的某振幅大于预设的阈值时，即获得一个候选检测.因为一个地震很可能被多个模板事件多次识别，匹配滤波方法会设定一个固定时窗.默认在这个固定时窗内只会发生一次地震，仅选择具备最大互相关值的候选检测为最终结果.
近年，有学者提出一些不以互相关为基础的匹配滤波方法.Skoumal等(2016)开发了一种专门用于探测重复地震的方法—Repeating Signal Detector(RSD).该方法使用欧几里得距离(Euclidian distance)来评判信号的相似性.如果两个信号完全相同，欧几里得距离是0.RSD只能处理单台数据.FAST方法把信号的特征信息提取为“指纹”，以匹配“指纹”的方式在单道波形记录中寻找相似信号(Yoon et al., 2015).匹配“指纹”比滑动互相关计算快很多.
随着近年机器学习的发展，机器学习也被用于地震检测(贾佳等, 2019).机器学习是实现人工智能的一种方式，即用普通计算机程序来模拟人类的智能.机器学习从数据中自动分析获得规律，并利用规律对未知数据进行预测.因为其中涉及大量的统计学理论，与推断统计学联系密切，机器学习也被称为统计学习理论.ConvNetQuake是一种基于机器学习的地震检测方法(Perol et al., 2018).该方法用非线性滤波器对波形进行分析和特征提取.在训练阶段，滤波器不断优化，进而能分辨地震信号和噪声.ConvNetQuake使用的具体技术是卷积神经网络(人工神经网络的一种，也称为深度学习).在ConvNetQuake之后，一些学者灵活应用卷积神经网络研发了很多新的机器学习地震检测方法.例如基于两个神经网络提取P波的到时和极性的方法(Ross et al., 2018a)、能区分P、S波的GPD(Ross et al., 2018b)和基于深度神经网络的PhaseNet(Zhu and Beroza, 2018).机器学习检测地震的最大优势是不会因为模板的数量增加而增加计算时间.但是目前的机器学习检测地震的方法都没有对来自不同观测台站的信号进行叠加，这不利于在噪声中检测地震.
多数匹配滤波方法默认检测到的地震和模板地震的震源位置相同，且不能检测到距离模板地震较远的小地震.这是因为当目标地震距离模板地震太远的时候，各个台站的互相关波形中的相关信号在叠加时会错开.Zhang和Wen(2015b)的Match & Locate方法可以同时对地震进行检测和定位.该方法对震源位置进行网格搜索，并在叠加互相关波形之前进行相对到时的校正.Match & Locate相对于其他多数匹配滤波方法，允许目标地震和模板地震有更远的距离，并提供相对定位.传统的匹配滤波(Peng and Zhao, 2009)是Match & Locate的特例.Liu等(2020)研发了基于GPU(Graphics Processing Unit，图形处理器)计算的Match & Locate.
上述的匹配滤波方法使用真实地震的真实波形作为模板波形，可以称为真实波形匹配滤波方法.这类方法依赖真实地震，存在一定的局限.例如，一次受关注的大地震发生后，地震学家和地震监测机构往往会在附近区域布设加密的观测台网.地震台站不可能记录到它布设好之前的地震.加密台网架设好之前的地震，包括主震，均不能作为这个加密台网的模板地震.真实波形匹配滤波方法无法立即应用在该加密台网上.Chamberlain和Townend (2018)提出了以虚拟地震的理论波形作为模板的理论波形匹配滤波方法.因为虚拟地震的波形是根据波动方程合成的，模板地震不足的问题便不再存在.然而，他们的研究仅仅使用了地震序列中主震的震源位置和震源机制解来合成模板波形.这会造成远离主震的小地震因各台站的互相关波形中的相关信号错开而漏检.而且，单一的震源机制解并不能代表所有地震的震源机制，可能会导致部分事件丢失.总之，不搜索震源位置和震源机制，虚拟地震的优势不能得到充分发挥.
本研究提出了一种以虚拟地震的理论格林函数为模板的匹配滤波方法(GFMF).该方法的大致流程如下：将研究区域划分为经度、纬度和深度三个维度的网格.每一个网格点为一个虚拟地震的震源点.各网格点上包括不同震源机制解的虚拟地震.所用的震源机制解离散地覆盖整个震源机制解空间.首先，依据这些虚拟地震的位置和台站位置合成格林函数.然后，假定虚拟地震在某一个网格点上，将观测数据和相应的理论格林函数滑动互相关得到互相关波形(CC).接着，假定虚拟地震具备某一个震源机制解，加权叠加各个互相关波形(加权系数在方法部分给出)，得到平均互相关波形(mean CC).用与传统匹配滤波方法一样的方式从平均互相关波形中得到候选检测.对所有震源位置和震源机制解重复上述步骤，得到所有的候选检测.最后，依据设定的固定时窗，从所有候选检测中选择出最终检测.观测数据与理论格林函数滑动互相关的这种方式让GFMF在保证计算结果不变的前提下，减少了为遍历震源机制解进行的滑动互相关计算的次数，从而减少计算时间.该优化借鉴自CAP方法(Zhao and Helmberger, 1994; Zhu and Helmberger, 1996; Zhu and Ben-Zion, 2013)中对于滑动互相关的处理方式，其数学基础是地震波形的互相关计算在本质上是时间序列的点积，遵守分配律.GFMF的细节和详细流程在方法部分介绍.
1 方法本节将先介绍GFMF的理论推导，再介绍阈值的设定方法，最后说明本研究使用的合成理论格林函数的具体技术.合成理论格林函数的具体技术并不属于GFMF方法本身，因为任何相关的具体技术都可以用于实现GFMF方法.在公式推导中，数学符号继承Gibbons和Ringdal (2006)，将地震波形和格林函数表示为离散的时间序列.
1.1 GFMF的理论推导在这一节，首先从最基本的互相关出发，逐步解释为什么GFMF可以在保证结果不变的情况下减少滑动互相关计算的次数，然后说明GFMF的整个流程.其中，严格区分互相关系数和互相关波形.
假定时间序列uN, Δt(0) 是某台站的某分量记录到的一段观测波形，其第一个采样点的时刻为0，采样间隔为Δt，具有N个采样点(下式中的T表示这是一个时间序列)：
|
(1) |
uM, Δt(t0)是uN, Δt(0) 的一个片段，其第一个采样点的时刻为t0，具有M个采样点(M≤N):
|
(2) |
vM, Δt(t1) 是一个模板波形，和uM, Δt(t0) 有相同的采样间隔和采样个数，但第一个采样点的时刻为t1:
|
(3) |
时间序列uM, Δt(t0)和vM, Δt(t1) 之间的点积或者说互相关系数定义为:
|
(4) |
进一步，uM, Δt(t0)和vM, Δt(t1) 的归一化点积或者说归一化互相关系数为：
|
(5) |
多数匹配滤波方法评价信号间的相似性依据的是公式(5).在程序实现上，可转为频率域的互相关计算.如果目标地震的观测波形和模板波形完全一样，它们的归一化互相关系数等于1.如果模板地震的发震时刻被修改为0，目标地震的发震时刻应为t0-t1.
如果模板地震是虚拟的，其波形可以表示为格林函数的加权叠加(安艺敬一和理查兹, 1986):
|
(6) |
其中，ai(i∈(1, n))是格林函数
以虚拟地震的理论波形为模板的匹配滤波方法(Chamberlain and Townend, 2018)首先依据公式(6)合成模板波形，然后依据公式(5)做模板波形与观测记录之间的互相关来评价信号间的相似性.在理论波形匹配滤波方法中，如果使用的震源机制有Nfc种，一个震源位置-台站对需要的归一化滑动互相关计算次数为3×Nfc(3为一个台站的分量数).
GFMF则不同.为了达到减少互相关计算次数的目的，GFMF将公式(6)代入公式(5)，得到公式(7)：
|
(7) |
〈vM, Δt(t1), vM, Δt(t1)〉可以表示为对
在搜索之前，研究区域需要划分为经度、纬度和深度三个维度的网格.每一个网格为一个震源位置.根据这些震源位置和台站位置，基于一维速度模型计算格林函数(图 1中的步骤一).为了在整个连续波形中搜索地震，我们需要从连续数据的开头到结尾，计算N-M+1次两M个采样点的时间序列(即连续观测波形的一个片段和一个格林函数)的归一化互相关系数，以得到这些归一化互相关系数组成的采样点个数为N-M+1，采样间隔为Δt的滑动归一化互相关波形CN-M+1, Δt(0)：
|
(8) |
|
图 1 GFMF工作流程图 Fig. 1 Workflow of the GFMF method |
公式(8)中的C′N-M+1, Δt(0) 只和震源参数(震源位置和震源机制)中的震源位置有关，而和震源机制无关.所以，GFMF先依据假定的虚拟地震的位置计算C′N-M+1, Δt(0) (本文称为CC)，即图 1中的步骤二.在假定了震源机制解后，对各CC加权叠加(权重系数为
|
(9) |
当mean CC的某个振幅值达到阈值时，即获得一个候选检测(图 1中的步骤四).GFMF用mean CC的MAD (median absolute deviation，中值绝对偏差)倍数做阈值.下一节会详细介绍.对所有震源位置和震源机制解重复上述步骤，得到所有的候选检测.经过上述搜索过程，可以得到大量的候选检测.因为同一个地震事件可能被多次检测到，我们假定在一个固定时间窗内只可能有一个地震发生，仅选择具备最大mean CC振幅的候选检测为最终检测，即新探测到的地震事件(图 1中的步骤五).
1.2 阈值的设定如果降低阈值，任何一种地震检测方法得到的检测数量都可以迅速增长.但是，这种增长伴随着更多的错误检测(下简称误检).地震检测都需要一个合理的阈值.GFMF用MAD的倍数划定阈值.Shelly等(2007)最先在匹配滤波法中用MAD的倍数划定阈值，而后为大量的匹配滤波方法相关研究所采用(Peng and Zhao, 2009; Aiken and Peng, 2014; Yao et al., 2015, 2017, 2018; Zhang and Wen, 2015a, b; Chamberlain and Townend, 2018).但是，很多学者都没有详细说明采用该方式的理由.本文对这个问题进行详细解释，以慰读者.
相对而言，噪声源是普遍的，而地震震源是罕见的.因此，地震观测波形中多数是噪声，而不是地震的信号.mean CC与观测波形密切相关，所以它的振幅也是更多和噪声有关，而不是地震.噪声是随机的.随机事件服从正态分布.所以，mean CC中的振幅的概率分布是正态分布.例如，图 2是后文利用实际数据进行的测试中的一天12个台站的8个虚拟震源的mean CC的振幅分布图.不难看出这就是正态分布典型的钟形曲线.
|
图 2 实际测试中的一天12个台站的8个虚拟震源的mean CC的振幅分布图 Fig. 2 Histograms for amplitudes values in mean CCs from 1-day real observed data of 12 seismic stations and 8 virtual-template events |
因为正态分布反映了噪声的随机性，所以它可以被用来估计噪声引起的误检数量(没有噪声就不会有误检).在正态分布中，大于某一个阈值的概率可以用累计分布函数计算.得到概率后乘以所有mean CC的采样总数就是估计的误检数量.如果阈值恰好是标准差的倍数，上述计算就可以大为化简.原因是每一个正态分布都可以经标准正态分布变换得到，而标准正态分布的相关计算则非常容易.Mean CC中振幅大于阈值(p)的采样点是误检的概率为P：
|
(10) |
其中，μ和σ分别是平均数和标准差.Φ是标准正态分布的累积分布函数，可用数值计算软件计算.标准差等于1.4826倍MAD(Rousseeuw and Croux, 1993).假设μ=0，并将p改写为MAD的倍数: 
|
(11) |
可得误检的数量NFD：
|
(12) |
Ncc是所有mean CC的采样总数.t和Δt分别是观测数据的时间长度和采样间隔. Nepicenter、Ndepth和NFM分别是震中、深度和震源机制的网格点的数目.
1.3 虚拟地震的理论格林函数的合成理论格林函数依据波动方程，可以用解析方法合成，也可以用数值方法合成.在本文中，我们采用一种频率-波数法(Zhu and Rivera, 2002).该研究提供了称为FK的程序包，可以合成限定频带的格林函数.FK在柱坐标系(ez, er, eθ)下，基于一维层状速度模型、震源深度和震中距合成三种基本断层的9个格林函数.这三种基本断层分别是45°倾滑断层、垂直倾滑断层和垂直走滑断层.FK计算的格林函数是阶跃函数在特定方位角产生的速度场.在FK中，公式(6)可在三个分量上分别具体化为：
|
(13) |
|
(14) |
|
(15) |
其中，Gz0、Gr0和Gθ0分别是45°倾滑断层在三分量上的格林函数(在一维速度层状模型中，Gθ0始终为0)；Gz1、Gr1和Gθ1分别是垂直倾滑断层在三分量上的格林函数；Gz2、Gr2和Gθ2分别是垂直走滑断层在三分量上的格林函数.因为Gθ0始终为0.这意味着，对于一个震源位置-台站对，GFMF需要做互相关的格林函数可进一步减少为8个.公式(6)中的ai，即在公式(11)-(13)中的格林函数的加权系数，是断层倾角δ、顺时针旋转的断层走向λ和顺时针旋转的台站方位角ϕ的函数.
2 应用测试 2.1 区域概况本文的应用测试选择在南加州的一个区域(图 3).选择该地区的原因有二.一是该地区的波形数据丰富.在小震中距的条件下，该地区有来自AZ(Frank Vernon and University California San Diego, 1982)和CI(California Institute of and United States Geological Survey, 1926)台网的12个宽频带台站可供使用(图 3).二是该地区有高品质的SCSN地震目录(Hutton et al., 2010)可以用来对GFMF的结果进行检验.据其介绍(https://scedc.caltech.edu/about/)，SCSN地震目录是美国完备性最好的地震目录，并且对本文涉及到的地震的定位结果均达到其说明文件(https://scedc.caltech.edu/eq-catalogs/docs/scec_dc.html)中定义的A等质量：震中误差在1 km内，深度误差在2 km内.
|
图 3 研究区域地图 红色网格点表示搜索过程中潜在的震中位置，经纬度方向间距0.02°.蓝色的沙滩球表示MW4.4主震.蓝色点是SCSN目录中的地震事件.三角形是地震台站.左下方的小图中的红框标明了研究区域的位置和大小. Fig. 3 Map of the study area The red dots are the potential epicenters in the grid search, whose grid intervals are 0.02° in longitude and latitude. The blue beachball represents the MW4.4 mainshock. The blue dots are the earthquakes in the SCSN catalog. The triangles are seismic stations. The red open square in the lower-inset marks the location of the study area. |
应用测试将搜索一次MW4.4地震前后14天的前、余震序列，并利用SCSN地震目录对搜索结果的时空分布特征、GR关系(b值和完备震级)和定位误差进行对比分析.然后与Match & Locate的搜索结果进行对比.最后，用GFMF给出中等地震的震源参数.
2.2 数据处理和主要的正反演参数本研究使用南加州地震数据中心(SCEDC)提供的web service服务获取波形数据和仪器响应文件.数据处理主要为使用SAC(Seismic Analysis Code，地震波分析软件包)(Goldstein et al., 2003; Goldstein and Snoke, 2005)对波形去毛刺、去均值、去线性趋势、波形歼灭、去仪器响应、重采样和滤波.搜索MW4.4地震的前震、余震序列的采样间隔为0.05 s，滤波频带为1~9 Hz.反演大于3级地震的采样间隔为0.2 s，滤波频带为0.01~0.2 Hz.
研究区域被划分为经度、纬度和深度三个维度的网格.经纬线方向上的网格间距为0.02°(图 3).震源深度为2.5、7.5、12.5和17.5 km，以平均覆盖0到20 km的深度范围.之所以特别地选择非整数的震源深度是因为FK在合成格林函数时要求震源深度和速度结构的界面深度不同.而本研究使用的当地的速度结构模型CVM-S4.26模型(Lee et al., 2014)是以整数深度为界面的.依据上述参数和可能的震中距计算格林函数，并进行滤波和时窗截取.滤波参数和对观测数据的一致.时窗截取范围为P波前1 s到S波后10 s.固定时窗为10 s，阈值为11倍MAD.
2.3 理论测试依次加入不同程度的干扰因素进行了四项理论测试，参见表 1.
|
表 1 各理论测试中不同干扰因素得到的不同的结果 Table 1 Different result obtained by different interference factors in theoretical tests |
假定要探测的地震位于116.82°W，33.5°N，深度12.5 km.发震时刻为60.0 s.震源机制的三个断层参数(strike，dip，rake)依次为150°、90°和0°.使用FK程序包的SYN命令，依据上述参数，为每一个台站合成120 s长度的波形数据.在这种情况下，GFMF没有任何干扰因素，可以正确得到输入的震源参数，并且mean CC最大值可达100.0%.能正确反演出FK的合成波形的震源参数说明了GFMF程序的正确性.
2.3.2 震源时间函数测试相对于基本测试，使用FK程序包的SYN命令在合成波形数据的时候，加入0.2 s的震源时间函数.这时，GFMF得到的发震时刻延迟了0.1 s.这一点并不意外，因为SYN命令加入的震源时间函数是一个等腰三角形.另外，平均互相关值降低为84.7%.这说明，震源时间函数对于反演各震源参数的影响完全是可以接受的.
2.3.3 噪声测试在震源时间函数测试的基础上，我们对FK合成的观测数据加入噪声进行了噪声测试.噪声是从真实的观测数据中选择的原始波形.叠加时，使地震信号与噪声间的信噪比为0.5.地震信号与噪声间的信噪比用P波到时前1 s到S波到时后10 s的方均根振幅量化.最后再对波形进行1~9 Hz的滤波，作为“观测数据”使用.经过GFMF搜索，除发震时刻延迟了0.1 s为60.1 s外，得到震源参数和假定的要探测的地震完全一致.发震时刻延迟的原因和前面的震源时间函数测试相同.此外，平均互相关值降低为34.0%.噪声测试说明GFMF凭借对互相关波形的叠加具有较好的对抗噪声的能力.
2.3.4 复杂模型测试最后，我们进行了复杂模型测试，以模拟实际应用中GFMF所用的速度模型和真实地球介质有较大差异的情形.复杂模型测试相对于上面的噪声测试的不同点是：FK依旧使用CVM-S4.26模型合成波形数据，而GFMF用的格林函数则是用PREM模型(preliminary reference Earth model)(Durek and Ekström, 1996)计算的.也就是说，复杂模型测试是在加入了震源时间函数和噪声的基础上，用简单的PREM模型搜复杂的CVM-S4.26模型.这样，我们就模拟了实际情况.这种情况下GFMF得到的发震时刻提前了1.5 s，得到的震中位置为输入位置的东南方的相邻网格点.深度则比输入值深5 km(即一个网格间距).另外，震源机制得到三个断层参数(strike，dip，rake)是150°、60°和0°，与输入值也有所不同.出现的偏差主要是速度模型差异造成的.最大平均互相关值是8.06%(图 4).尽管两个模型有非常大的差异，GFMF得到的发震时刻的误差依旧很小.这说明GFMF在复杂模型的情况下，虽然得到的震源参数不完全正确，检测事件的目的依旧可以达成.
|
图 4 复杂模型测试的波形拟合情况 (a) 三列波形分别为R、T和Z分量.参考时刻为发震时刻.台站名称标示在左侧.黑色波形为观测波形.红色波形是依据GFMF得到的震源参数合成的.每个波形的右上角显示互相关系数.(b) 平均互相关波形. Fig. 4 Waveform fitting of the complex model test (a) Traces in the three columns are the R, T, and Z components, respectively. The reference time is the origin time. Seismic station names are shown on the left. The black traces are observed waveforms. The red synthetic waveforms are calculated based on the source parameters obtained by the GFMF method. The upper right corner of each waveform displays the cross-correlation coefficient. (b) Mean CC. |
在对MW4.4地震的前震、余震序列的搜索中，GFMF识别到3273个地震事件.为完成这一数据长度为14天，采样间隔为0.05 s的搜索，一台戴尔®7070mf商务微型电脑(Intel® i7-9700处理器)需要11小时10分钟.GFMF在计算效率上可以实现实时处理.
图 5展示了一个被GFMF搜索到但未被SCSN地震目录收录的地震的波形.尽管该地震的波形(黑色)和相应虚拟地震的波形(红色)未能良好拟合，但该地震的波形清晰，证实是一次真实存在的地震.其中一些分量的互相关系数是非常低的，甚至有负数的情况.但是，平均互相关波形中却出现非常尖锐的一个突起.这个例子说明在高频带下，一维速度结构对于真实地球介质过度简化，使得GFMF无法良好地拟合观测波形，但得益于对互相关波形的叠加，不足以妨碍对地震事件的识别.Gibbons和Ringdal(2006)也得出过识别地震并不需要模板波形与目标波形有很高的相似度的结论.
|
图 5 一个被GFMF识别到，但未被SCSN地震目录收录的地震的波形拟合情况 (a) 三列波形分别为R、T和Z分量.参考时刻为发震时刻.台站名称标示在左侧.黑色波形为观测波形.红色波形是依据GFMF得到的震源参数合成的.每个波形的右上角显示互相关系数；(b) 平均互相关波形.该事件的发震时刻是2018-08-23T08∶42∶02.690. Fig. 5 Waveform fitting of one example event newly detected by the GFMF method but not existing in the SCSN catalog (a) Traces in the three columns are the R, T, and Z components, respectively. The reference time is the origin time. Seismic station names are shown on the left. The black traces are observed waveforms. The red synthetic waveforms are calculated based on the source parameters obtained by the GFMF method. The upper right corner of each waveform displays the cross-correlation coefficient. (b) Mean CC. The origin time of this detection is 2018-08-23T08:42:02.690. |
图 6展示搜索到的事件的时空特征和SCSN地震目录的对比.整体来看，二者的地震活动性有相似的时间变化和空间分布特征.二者的地震数量都是在8月12日到14日相对较少，而在15日主震发生后，地震数量骤然增加(图 6a).随后，二者的地震数量均按余震的一般规律，逐渐减少.GFMF的搜索结果和SCSN地震目录的空间分布也是相似的.在平面图中，SCSN地震目录中的地震事件主要分布在主震所在的及邻近的东、北和东北方向的四个网格中(图 6b).GFMF搜索结果也在这四个网格中有很多事件(图 6b).另外，有一些SCSN地震事件以主震为起点，近直线地向东北方向分布(图 6b).GFMF搜索结果也有这一特征.在剖面图中，二者均显示地震序列以主震为起点，向东逐渐向深部发展(图 6c).而且，GFMF搜索结果集中分布的倾角与SCSN目录事件的延伸倾角非常接近.可见，GFMF的搜索结果正确揭示了该地震序列的时间变化和空间分布特征.这说明，利用GFMF可以正确地认识该地震序列的时空特征.
|
图 6 GFMF搜索结果和SCSN地震目录中地震事件的时间变化和空间分布，及二者震源参数的差异分布 (a) 地震数量时间分布图.时间跨度为2018年8月12日到25日.每一个直方跨度6 h，其高度表示这6 h内的地震数量.红色和蓝色直方分别为GFMF搜索结果和SCSN地震目录.虚线表示2018年8月15日的MW4.4主震的发震时刻.(b) 地震事件空间分布平面图.蓝色点表示SCSN地震目录中的事件.网格区域的颜色表示该网格的GFMF搜索结果的地震数量.五角星为主震.(c)和(b)类似，但为地震事件空间分布的剖面图.(d, e, f) GFMF搜索结果和SCSN地震目录在发震时刻、震中位置和震源深度的差异分布图.图中虚线标明各自中值. Fig. 6 Time evolution and spatial distribution of GFMF detections and the SCSN catalog, and their difference of source parameters (a) Histograms for the number of events in 6-hour bins from August 12 to 25, 2018. The red and blue histograms are for GFMF detections and the SCSN catalog, respectively. The dashed line indicates the origin time of the MW4.4 mainshock on August 15, 2018. (b) Map view for the spatial distribution. The blue dots represent the earthquakes in the SCSN catalog. The color represents the number of GFMF detections in each square. The star represents the MW4.4 mainshock. (c) Similar to the (b), but in cross-section view. (d, e, f) Histograms of the difference respectively in the origin time, epicenter, and source depth between the GFMF detections and SCSN catalog. The dashed lines indicate the median differences. |
GFMF对SCSN地震目录中的地震事件的识别率是97.1%.我们在统计识别率时，采用了很高的标准.确认SCSN目录中的某事件A被GFMF成功搜索到，GFMF的搜索结果中就必须有事件B满足三个条件：(1) A与B的发震时刻的时差不大于固定时窗，(2) A是SCSN地震目录中，发震时刻上和B最接近的.(3) B是GFMF搜索结果中，发震时刻上A最接近的.在1056个SCSN目录事件中，1025个被GFMF搜索到，比例为97.1%.
我们借助SCSN地震目录对GFMF搜索结果的定位误差进行了简单估计.GFMF本身并不能给出定位误差，因为在本质上GFMF是从众多的网格点中依据互相关系数选择一个最优结果.也就是说，GFMF只能提供不同网格点的差异(可以用互相关系数的方差或者标准差量化).这和大家关心的定位误差并不是一回事.实际上，目前没有匹配滤波方法可以给出定位误差，包括Match & Locate(Zhang and Wen, 2015b).此处，我们借鉴Zhang和Wen(2015b)的作法，用其他的认为更准确的定位结果来估计定位误差.上文已经说明，SCSN地震目录中事件的震中误差在1 km内，深度误差在2 km内.所以我们用SCSN地震目录中已经有的1056个事件的定位结果对GFMF的结果进行了误差估计.我们直接默认SCSN的定位是准确的.这样.GFMF在发震时刻、震中和震源深度上的误差的中值分别为0.50 s、4.2 km和4.2 km(图 6d, e, f).
2.4.3 GR关系分析我们利用SCSN地震目录中的震级信息简要分析了GFMF搜索结果的GR关系(图 7).首先，利用SCSN地震目录的震级信息估计GFMF搜索结果的震级：依据SCSN地震目录中的事件在一个台站的竖直分量的最大振幅的对数(以10为底)，用直线拟合的方式，建立波形振幅和震级的经验关系.显然，不同台站的经验关系式是不同的.选用残差最小的CI_PLM台站的经验关系式计算搜索结果的震级(图 7a).最终，GR分析的结果显示，GFMF搜索结果的b值和完备性震级分别为1.12和-0.2(图 7b).相比，SCSN地震目录的b值和完备性震级分别为1.09和0.3.
|
图 7 震级、波形振幅和地震数量的关系 (a) 震级和波形振幅的关系.蓝色点表示SCSN地震目录中的事件，其横坐标为相应事件在CI_PLM台站的竖直分量的最大振幅值的对数(10为底)，纵坐标是SCSN地震目录中的震级.紫色直线为拟合的经验关系.(b) 震级和地震数量的关系(GR关系).横坐标为震级，纵坐标为不小于该震级的事件的数量的对数(10为底).红色圆点和蓝色三角分别表示GFMF搜索结果和SCSN地震目录.红色和蓝色虚线为分别对GFMF搜索结果和SCSN地震目录的拟合.竖直的两条黑色虚线指出了-0.2和0.3的震级，以便读图. Fig. 7 Relationship between magnitude, waveform amplitudes and the number of earthquakes (a) The relationship between magnitude and waveform amplitude. The blue dots represent SCSN cataloged earthquakes. The abscissa is the logarithm of the maximum amplitudes in the vertical component at CI_PLM station, and the ordinate is the magnitude from the SCSN catalog. The purple line is the empirical relationship of the fitting. (b) The relationship between magnitude and number of earthquakes (GR relationship). The abscissa is the magnitude, and the ordinate is the logarithm of the number of events not less than the magnitude. The red dots and blue triangles represent the GFMF detections and the SCSN cataloged events, respectively. The red and blue dashed lines are the fittings of the GFMF detections and the SCSN cataloged events, respectively. The two vertical black dashed lines indicate the magnitudes of -0.2 and 0.3 for easy reading. |
为了对比GFMF和真实波形匹配滤波方法的优缺点，我们选择真实波形匹配滤波方法中最好的Match & Locate搜索了这条地震序列.
2.5.1 主要的数据处理和反演参数说明Match & Locate测试中的多数数据处理参数与GFMF一致.阈值与GFMF的测试一样，也用11倍.模板波形和连续数据的滤波参数等与GFMF的测试一样.Match & Locate做到时校正用的一维速度模型与GFMF测试中计算格林函数的模型一样.所用的模板是该地区自2000年以来的大于3级的地震.
Match & Locate和GFMF的测试的参数差异主要是搜索网格.Match & Locate的测试中使用的搜索网格是Match & Locate 2.0版本程序包提供的默认值：搜索中心为相应模板的震中.搜索范围在经纬度方向上最大0.05°，网格间隔0.005°.深度上不进行搜索.搜索中心、范围和间隔与GFMF测试中不同的原因是GFMF进行网格搜索的目的是增加模板数量，进而增加检测数量，并筛选出最好的结果.Match & Locate进行网格搜索的目的主要是为了在模板事件与目标事件的距离远小于震中距的前提下进行相对定位.Match & Locate在检测小地震时通常是不搜索深度的，因为小地震只有近台才有信号(Zhang and Wen, 2015b).Match & Locate中不同网格点的模板波形没有波形形状的差异，只有到时的差异.如果只使用到时，除非既有近台又有远台，反演震源深度和发震时刻是耦合的(Billings et al., 1994; Zhang et al., 2014).而GFMF中，不同深度的模板的波形形状是不同的.GFMF对深度的搜索不仅基于到时，也基于波形形状.
2.5.2 测试结果Match & Locate一共检测到3213个地震事件(GFMF为3273个)，对1056个SCSN地震目录事件识别到1018个，识别率为96.4%(GFMF为97.1%).Match & Locate和GFMF的检测总数和识别率相当.
我们将Match & Locate分别与SCSN地震目录和GFMF进行更细致的对比.Match & Locate的搜索结果(图 8a中灰色直方)与SCSN地震目录(图 8a中蓝色直方)在单位时间的检测数量变化是同步的.先前在对比GFMF搜索结果和SCSN地震目录时，我们看到过同样的现象.Match & Locate得到的发震时刻与SCSN的差异的中值为0.26 s(图 8b)，震中差异的中值为1.8 km(图 8c).Match & Locate的定位精度明显好于GFMF(GFMF的发震时刻误差和震中误差分别为0.50 s和4.2 km).图 9中显示了Match & Locate(黑色实心曲线)和GFMF(红色曲线)的累积检测数量，可以看到Match & Locate与GFMF得到的地震的总数非常接近，且均明显多于SCSN地震目录(蓝色曲线).二者的数量除了在主震发生的8月15日明显激增外，在8月21日也有一个共同的小幅突变.经对比检查，GFMF和Match & Locate有2114个检测是重合的.可见GFMF与Match & Locate可以相互印证.这2114个事件的发震时刻的差异的中值为0.72 s，定位差异的中值为5.6 km(图 10).上述测试说明Match & Locate也能识别到大量事件.
|
图 8 Match & Locate搜索结果和SCSN地震目录中地震事件的时间变化，及二者震源参数的差异分布 (a) 地震数量时间分布图.时间跨度为2018年8月12日到25日.每一个直方跨度6 h，其高度表示这6 h内的地震数量.灰色和蓝色直方分别为Match & Locate搜索结果和SCSN地震目录.虚线表示2018年8月15日的MW4.4主震的发震时刻.红色的阶梯直方图表示仅使用部分模板的结果.(b, c) Match & Locate搜索结果和SCSN地震目录在发震时刻和震中位置的差异分布图.虚线标明各自中值.红色阶梯直方图表示了仅使用部分模板的结果.红色虚线表示仅使用部分模板的中值. Fig. 8 Time evolution of Match & Locate and the SCSN catalog, and their difference of source parameters (a) Histograms for the number of events in 6-hour bins from August 12 to 25, 2018. The gray and blue histograms are for Match & Locate detections and the SCSN catalog, respectively. The dashed line indicates the origin time of the MW4.4 mainshock on August 15, 2018. The red ladder histogram shows the results of using part templates. (b, c) Histograms of the difference respectively in the origin time and epicenter, between the Match & Locate detections and SCSN catalog. The dashed lines indicate the median differences. The red ladder histogram shows the results of using part templates. The red dashed line indicates the median value of part templates used. |
|
图 9 地震数量的时间累积图 红色、蓝色曲线分别表示GFMF的搜索结果和SCSN地震目录的结果.黑色的实线和虚线分别表示Match & Locate用全部和部分模板的结果. Fig. 9 Time cumulative number of earthquakes The red and blue curves represent the GFMF detections and the SCSN catalog, respectively. The black solid and the dashed line represent the Match & Locate detections with all and part templates, respectively. |
|
图 10 GFMF和Match & Locate搜索结果在发震时刻(a)和震中位置(b)的差异分布图 虚线表示中值. Fig. 10 Histograms of the difference respectively in the origin time (a) and epicenter (b) between the Match & Locate detections and SCSN catalog The dashed lines indicate the median differences. |
未来的地震是不能做模板的.我们假设学者进行科研工作是在2018年8月26日(该地震序列的最后一天是25日).这时候，该地区只有6个3级以上地震.我们只用这部分模板做测试，检查Match & Locate对模板的依赖程度.经试验，只使用这部分模板，检测总数从3213个减少为2422个，对SCSN地震目录只能识别963个.图 8a的红色阶梯直方表示了只是用部分模板的情形，可看到减少模板会明显减少Match & Locate的检测数量.图 9的黑色虚线也明显看到减少模板后的地震数量明显比使用全部模板时更少.另外，在模板减少后，定位误差也略微变大.发震时刻的误差增大至0.27 s(图 8b)，震中误差为1.9 km(图 8c).
2.6 中等地震的震源参数GFMF可以准确地联合反演该地区的中等地震(大于3级)的发震时刻、震源位置和震源机制.该地区自2000年以来震级大于3级的地震一共15个.表 2展示了这15个地震的GFMF反演结果和SCSN地震目录的震源参数差异的统计值，从中可以看到各种震源参数和SCSN地震目录都非常接近.发震时刻的差异的中值为1.25 s.震中和震源深度的误差中值均在2 km以内.震源机制解的P、T和N轴在三维空间内的夹角的中值和平均值小于15°.这表明，GFMF搜索结果和SCSN地震目录的震源机制的差异非常小.图 11展示了每一个地震的对比，从中可以得到更加直观的认识.SCSN地震目录反演震源机制使用的是初动和振幅比的方法.经上述对比可见，GFMF与这些传统方法获得的震源机制具有很好的一致性.另外，GFMF自身的波形拟合情况也较好.图 12展示了其中一个地震的各台站的波形拟合情况.可以看到绝大多数分量的拟合情况很好.部分分量如AZ_RDM的T分量和CI_LKH的Z分量的观测波形和理论波形在到时上存在错位.考虑到绝大多数分量已经拟合得很好，这些错位的可能原因是局部的速度异常.传统的震源机制反演程序，如CAP类的方法(Zhao and Helmberger, 1994; Zhu and Helmberger, 1996; 陈伟文等, 2012; Zhu and Ben-Zion, 2013; Chen et al., 2015; Bai et al., 2020)，可以对单分量进行平移，得到更好的拟合.GFMF需要对发震时刻、震源位置和震源机制进行联合反演，不能对波形进行平移.尽管如此，GFMF在低频滤波(0.01~0.2 Hz)下的拟合情况已经足够支持震源机制反演的结果.综上所述，与SCSN地震目录的对比和波形拟合情况可以证实GFMF在低频滤波的条件下可以正确得到中等地震包括震源机制在内的震源参数.
|
|
表 2 大于3级地震的GFMF搜索结果和SCSN地震目录的震源参数差异的统计值 Table 2 Statistical values of the difference in source parameters between GFMF detections and SCSN cataloged earthquakes with magnitude greater than 3 |
|
图 11 全部大于3级地震的GFMF反演结果 每一副子图对应于一个地震.蓝色和红色的沙滩球分别显示SCSN地震目录和GFMF反演的震源机制解.沙滩球上方是该地震的发震时刻和震级.沙滩球下第一排数字分别是GFMF反演结果和SCSN地震目录在发震时刻、震中和震源深度上的差异，第二排数字分别是两个震源机制解的P轴、T轴和N轴的夹角. Fig. 11 Inversion results of GFMF for earthquakes (Mag>3) Each sub-picture shows an earthquake′s information. The blue and red beachballs belong to the SCSN catalog and GFMF detections, respectively. The origin time and the magnitude are above the beachballs. The first row of numbers below the beachballs are the differences between the GFMF detections and the SCSN cataloged earthquakes in the origin time, epicenter, and source depth. The numbers in the second row are the angles of the P-axis, T-axis and N-axis between the two focal mechanism solutions. |
|
图 12 地震事件20181109145809的波形拟合情况 (a) 三列波形分别为R、T和Z分量.参考时刻为发震时刻.台站名称标示在左侧.黑色波形为观测波形.红色波形是依据GFMF得到的震源参数合成的.蓝色和红色竖线是理论P、S到时.(b) 平均互相关波形.每个波形的右上角显示互相关系数.该地震震级3.78级.采样间隔为0.2 s.滤波频带为0.01~0.2 Hz. Fig. 12 Waveform fitting of seismic event 20181109145809 (a) Traces in the three columns are the R, T, and Z components, respectively. The reference time is the origin time. Seismic station names are shown on the left. The black traces are observed waveforms. The red synthetic waveforms are calculated based on the source parameters obtained by the GFMF method. The blue and red vertical lines are the theoretical P and S time. (b) Mean CC. The magnitude of the earthquake is 3.78. The sampling interval is 0.2 s. The filtering frequency band is 0.01~0.2 Hz. |
合理划定阈值对所有地震检测方法都至关重要.不能为了获得更多的检测数量而无限制地降低阈值.可以依据公式(12)估计候选检测中的误检数.在本文的测试中，当阈值达到11倍MAD时，候选检测中的误检数约为0.08，远远小于1(表 3).即按照统计学规律可以认为结果中没有误检，这就是本文选择11倍MAD作为阈值的原因.如果进一步假设误检被筛选出来而成为最终检测的概率和正确检测一样，那么误检率就等于候选检测中的误检数除以候选检测的总数.在众多的采用了MAD划定阈值的匹配滤波方法的研究中，阈值一般是8倍MAD(Shelly et al., 2007)到13倍MAD(Yao et al., 2018).这个阈值范围内的GFMF检测总数和对SCSN地震目录的识别率见表 4.
|
|
表 3 不同阈值的误检估计 Table 3 False detection estimation of different thresholds |
|
|
表 4 不同阈值的检测总数和识别率 Table 4 Total number of detections and recognition rate of different thresholds |
格林函数做模板让GFMF在遍历震源机制解的部分节省了大量时间.无论是哪一种匹配滤波方法，最耗时的计算是滑动归一化互相关.一个虚拟地震在一个台站的格林函数和震源机制无关，恒定为8个(见1.3节).因此，在GFMF中，无论震源机制解的个数为多少，一个震源位置-台站对的互相关次数恒定为8次.因为理论波形和震源机制有关，理论波形匹配滤波法(Chamberlain and Townend, 2018)的滑动归一化互相关次数是3×Nfc(Nfc为震源机制解个数).也就是说，GFMF让滑动归一化互相关次数从3×Nfc减少为了8次.如果虚拟地震的震源位置个数和台站数量等其他因素均相同，只考虑震源机制解为自变量，GFMF将互相关部分的时间复杂度从线性级别降低到了常数级别，即计算时间和震源机制解个数原本为线性关系，降低为了常数关系.所以，GFMF可以节省大量时间.
3.3 GFMF与真实波形匹配滤波方法的比较经与Match & Locate的对比测试可以看出，GFMF与真实波形匹配滤波方法互有优势.首先，二者的结果中很大一部分是重合的，可以相互印证.另外，Match & Locate的定位精度明显比GFMF更为准确，其原因是真实模板波形包含了真实地球的介质信息.然而，GFMF的优势也很明显.使用GFMF不需要寻找真实地震做模板，不需要人工检视模板事件的波形数据质量.经上文的测试，若减少模板，Match & Locate检测到的地震数量会明显下降.真实波形匹配滤波方法要求研究区域有足够的模板事件，且需要人工挑选具有良好地震记录的波形数据.因为不需要真实波形做模板，GFMF不受这些限制.这使得一些场合更加适合GFMF而不是真实波形匹配滤波法，例如：(1)地震空区：地震空区的地震活动性不高，可以做模板事件的地震少.(2)短期运行的台网：一些临时台网往往运行时间较短.在这较短的时间内，可以做模板事件的地震可能很少.(3)实时监控：运用GFMF需要的人工工作量比真实匹配滤波方法少很多.例如，真实波形匹配滤波方法需要人工检视模板事件的波形数据，挑选出其中记录良好的台站和分量.另外，在台网运行的初期，真实模板地震的数量是不足的.需要等待真实模板地震的出现.而GFMF在新台网架设完毕后，可以立即开展工作.
3.4 GFMF对震源位置和震源机制的遍历GFMF对震源位置和震源机制进行了网格搜索，增加了模板数量，能搜索到更多小地震.正如在应用部分已经展示的，凭借对震源位置的网格搜索，GFMF能正确揭示地震序列的空间分布特征.对震源机制解的网格搜索使得GFMF能有更多的模板地震，从而识别到更多的地震事件.经测试，如果只使用单一的震源机制解，无论使用哪一个震源机制，得到的地震事件总数和识别率都远远低于完整震源机制解空间(图 13).单一震源机制解得到的事件总数最多只有965个，而全震源机制解是3273个；单一震源机制解得到的识别率最高只有60.4%，而全震源机制解是97.1%.各个单独的震源机制解得到的事件总数和识别率的中值分别为649和42.6%.主震的震源机制解得到的事件总数和识别率分别为531和35.7%，低于众震源机制解的平均水平(图 13).如果进行低频滤波，尽管不能像CAP类的方法平移波形以得到更好地拟合，GFMF依然可以得到中等地震(大于3级)的震源机制解.
|
图 13 完整震源机制解空间和单个震源机制解识别到的地震总数(a)和识别率(b)的对比 (a) 纵坐标是检测到的地震的数量，横坐标是震源机制的序号.顶部的红色虚线表示完整震源机制解空间的检测数量.每一个黑色点表示相应序号的震源机制的检测数量.震源机制解的序号按检测数量的降序排列.蓝色沙滩球表示了MW4.4主震的震源机制的位置.(b) 与(a)类似，但纵坐标是对SCSN地震目录事件的识别率. Fig. 13 Total detection number (a) and detection rate (b) of virtual-template events with varying and fixed focal mechanism(s) (a) The ordinate is the detection number, and the abscissa is the index of the focal mechanisms. The red dashed line at the top represents the detection number of varying focal mechanisms. Each black point represents the detection number of the corresponding focal mechanism. The index numbers of focal mechanisms are in descending order of the detection number. The blue beachball indicates the index of the focal mechanism of the MW4.4 mainshock. (b) Similar to (a), but for the recognition rate of the SCSN catalog. |
GFMF是一种以虚拟地震的理论格林函数为模板的匹配滤波方法.该方法是把观测波形和虚拟地震的理论格林函数做互相关，再将各互相关波形进行加权叠加得到平均互相关波形，最后依据平均互相关波形的振幅值识别地震.以格林函数为模板让GFMF在计算结果不变的前提下，减少互相关计算次数，节省计算时间.因为对震源位置和震源机制的网格搜索，GFMF可以探测到整个研究区域内所有可能位置的地震事件.和其他任何地震学方法一样，GFMF的应用也有适用的频带.在小地震的高频频带(1~9 Hz)下，GFMF可以对小地震进行识别，并凭借对位置和震源机制的遍历增加模板的数量而检测到更多地震事件和获得更高的检测率.在低频频带(0.01~0.2 Hz)下，GFMF可以可靠地反演中等地震(大于3级)的震源机制解等震源参数.
未来对GFMF的改进可以有两个方面的考虑，一是进一步提高计算速度(如支持GPU计算等).二是增加对完整矩张量的搜索.技术上需要再增加三个格林函数(Zhu and Rivera, 2002).自然界存在一些非双力偶的地震(Foulger and Julian, 2014)，而GFMF完全不对非双力偶成分进行搜索，会增加丢失这部分地震的风险.GFMF会在未来增加相关功能.目前，GFMF已经在网络公开，请访问www.seiswave.cn.
致谢 两位匿名审稿专家和《地球物理学报》编辑部提供了大量的修改建议.田冬冬博士给予了大量的指导和帮助.孟晓波博士、张淼教授、彭志刚教授、邓凯教授和王芳博士提供了意见建议.在此对他们深表感谢.除流程图外，本文其他图件使用GMT6 (Wessel et al., 2019)绘制.我们感谢GMT开发人员的无私奉献.
Aiken C, Peng Z G. 2014. Dynamic triggering of microearthquakes in three geothermal/volcanic regions of California. Journal of Geophysical Research: Solid Earth, 119(9): 6992-7009. DOI:10.1002/2014JB011218 |
Aki K, Richards P G. 1986. Quantitative Seismology Theory and Methods (in Chinese). Beijing: Seismological Press: 59.
|
Bai Q P, Ni S D, Chu R S, et al. 2020. gCAPjoint, a software package for full moment tensor inversion of moderately strong earthquakes with local and teleseismic waveforms. Seismological Research Letters, 91(6): 3550-3562. DOI:10.1785/0220200031 |
Billings S D, Kennett B L N, Sambridge M S. 1994. Hypocentre location: genetic algorithms incorporating problem-specific information. Geophysical Journal International, 118(3): 693-706. DOI:10.1111/j.1365-246X.1994.tb03994.x |
California Institute of Technology, United States Geological Survey Pasadena. 1926. Southern California Seismic Network. International Federation of Digital Seismograph Networks. DOI:10.7914/SN/CI |
Chamberlain C J, Townend J. 2018. Detecting real earthquakes using artificial earthquakes: on the use of synthetic waveforms in matched-filter earthquake detection. Geophysical Research Letters, 45(21): 11641-11649. DOI:10.1029/2018GL079872 |
Chen W W, Ni S D, Wang Z J, et al. 2012. Joint inversion with both local and teleseismic waveforms for source parameters of the 2010 Kaohsiung earthquake. Chinese Journal of Geophysics (in Chinese), 55(7): 2319-2328. DOI:10.6038/j.issn.0001-5733.2012.07.017 |
Chen W W, Ni S D, Kanamori H, et al. 2015. CAPjoint, a computer software package for joint inversion of moderate earthquake source parameters with local and teleseismic waveforms. Seismological Research Letters, 86(2A): 432-441. DOI:10.1785/0220140167 |
Cochran E S, Ross Z E, Harrington R M, et al. 2018. Induced earthquake families reveal distinctive evolutionary patterns near disposal wells. Journal of Geophysical Research: Solid Earth, 123(9): 8045-8055. DOI:10.1029/2018JB016270 |
Crampin S, Gao Y. 2015. The physics underlying Gutenberg-Richter in the earth and in the moon. Journal of Earth Science, 26(1): 134-139. DOI:10.1007/s12583-015-0513-3 |
Dieterich J. 1994. A constitutive law for rate of earthquake production and its application to earthquake clustering. Journal of Geophysical Research: Solid Earth, 99(B2): 2601-2618. DOI:10.1029/93JB02581 |
Durek J J, Ekström G. 1996. A radial model of anelasticity consistent with long-period surface-wave attenuation. Bulletin of the Seismological Society of America, 86(1A): 144-158. |
Foulger G R, Julian B R. 2014. Non-double-couple earthquakes. //Beer M, Kougioumtzoglou I A, Patelli E eds. Encyclopedia of Earthquake Engineering. Berlin Heidelberg: Springer, 1-31.
|
Gibbons S J, Ringdal F. 2006. The detection of low magnitude seismic events using array-based waveform correlation. Geophysical Journal International, 165(1): 149-166. DOI:10.1111/j.1365-246X.2006.02865.x |
Gibbons S J, Ringdal F. 2012. Seismic monitoring of the North Korea nuclear test site using a multichannel correlation detector. IEEE Transactions on Geoscience and Remote Sensing, 50(5): 1897-1909. DOI:10.1109/TGRS.2011.2170429 |
Goldstein P, Dodge D, Firpo M, et al. 2003. SAC2000: Signal processing and analysis tools for seismologists and engineers. International Geophysics, 81: 1613-1614. |
Goldstein P, Snoke A. 2005. SAC availability for the IRIS community. Incorporated Research Institutions for Seismology Newsletter, 7(1): 11. |
Harrington R M, Brodsky E E. 2009. Source duration scales with magnitude differently for earthquakes on the San Andreas fault and on secondary faults in parkfield, California. Bulletin of the Seismological Society of America, 99(4): 2323-2334. DOI:10.1785/0120080216 |
Hou J X, Wang B S. 2017. Temporal evolution of seismicity before and after the 2014 Ludian MS6.5 earthquake. Chinese Journal of Geophysics, 60(4): 1446-1456. DOI:10.6038/cjg20170418 |
Hou J X, Xie F, Ren Y Q, et al. 2020. Detection of acoustic emissions associated with the stick-slips of a meter-scale fault in laboratory by using the matched filter technique. Chinese Journal of Geophysics (in Chinses), 63(4): 1630-1641. DOI:10.6038/cjg2020N0102 |
Hutton K, Woessner J, Hauksson E. 2010. Earthquake monitoring in Southern California for Seventy-Seven Years (1932—2008). Bulletin of the Seismological Society of America, 100(2): 423-446. DOI:10.1785/0120090130 |
Jia J, Wang F Y, Wu Q J. 2019. Review of the application of machine learning in seismic detection and phase identification. China Earthquake Engineering Journal (in Chinese), 41(6): 1419-1425. DOI:10.3969/j.issn.1000-0844.2019.06.1419 |
Jiang J Z, Chen Q F, Li J. 2019. Waveform detection and location of deep earthquakes in the subduction zone beneath Northeast China. Chinese Journal of Geophysics (in Chinese), 62(8): 2930-2945. DOI:10.6038/cjg2019M0210 |
Kato A, Nakagawa S. 2014. Multiple slow-slip events during a foreshock sequence of the 2014 Iquique, Chile MW8.1 earthquake. Geophysical Research Letters, 41(15): 5420-5427. DOI:10.1002/2014GL061138 |
Lee E J, Chen P, Jordan T H, et al. 2014. Full-3-D tomography for crustal structure in Southern California based on the scattering-integral and the adjoint-wavefield methods. Journal of Geophysical Research: Solid Earth, 119(8): 6421-6451. DOI:10.1002/2014JB011346 |
Li L, Wang B S, Hou J X. 2017. Applications of matched filter technique in seismic data processing. Earthquake Research in China (in Chinese), 33(1): 14-22. DOI:10.3969/j.issn.1001-4683.2017.01.002 |
Liu M, Li H Y, Zhang M, et al. 2020. Graphics Processing Unit-Based Match and Locate (GPU-M & L): an improved match and locate method and its application. Seismological Research Letters, 91(2A): 1019-1029. DOI:10.1785/0220190241 |
Meng X F, Yang H F, Peng Z G. 2018. Foreshocks, b value map, and aftershock triggering for the 2011 MW5.7 Virginia earthquake. Journal of Geophysical Research: Solid Earth, 123(6): 5082-5098. DOI:10.1029/2017JB015136 |
Miller S A, Collettini C, Chiaraluce L, et al. 2004. Aftershocks driven by a high-pressure CO2 source at depth. Nature, 427(6976): 724-727. DOI:10.1038/nature02251 |
Peng Z G, Zhao P. 2009. Migration of early aftershocks following the 2004 Parkfield earthquake. Nature Geoscience, 2(12): 877-881. DOI:10.1038/ngeo697 |
Perol T, Gharbi M, Denolle M. 2018. Convolutional neural network for earthquake detection and location. Science Advances, 4(2): e1700578. DOI:10.1126/sciadv.1700578 |
Ross Z E, Meier M A, Hauksson E. 2018a. P wave arrival picking and first-motion polarity determination with deep learning. Journal of Geophysical Research: Solid Earth, 123(6): 5120-5129. DOI:10.1029/2017JB015251 |
Ross Z E, Meier M A, Hauksson E, et al. 2018b. Generalized seismic phase detection with deep learning. Bulletin of the Seismological Society of America, 108(5A): 2894-2901. DOI:10.1785/0120180080 |
Rousseeuw P J, Croux C. 1993. Alternatives to the median absolute deviation. Journal of the American Statistical Association, 88(424): 1273-1283. DOI:10.1080/01621459.1993.10476408 |
Schmittbuhl J, Karabulut H, Lengliné O, et al. 2016. Long-lasting seismic repeaters in the Central Basin of the Main Marmara Fault. Geophysical Research Letters, 43(18): 9527-9534. DOI:10.1002/2016GL070505 |
Shearer P M. 2002. Parallel fault strands at 9-km depth resolved on the Imperial Fault, Southern California. Geophysical Research Letters, 29(14): 19-1. DOI:10.1029/2002GL015302 |
Shearer P M. 2009. Introduction to Seismology. 2nd ed. Cambridge: Cambridge University Press.
|
Shelly D R, Beroza G C, Ide S. 2007. Non-volcanic tremor and low-frequency earthquake swarms. Nature, 446(7133): 305-307. DOI:10.1038/nature05666 |
Sheng M H, Chu R S, Wei Z G, et al. 2018. Study of microseismicity caused by Xishancun Landslide deformation in Li county, Sichuan Province. Chinese Journal of Geophysics (in Chinese), 61(1): 171-182. DOI:10.6038/cjg2018L0367 |
Sianipar D. 2020. Immediate foreshocks activity preceding the 2018 MW7.5 Palu earthquake in Sulawesi, Indonesia. Pure and Applied Geophysics, 177(6): 2421-2436. DOI:10.1007/s00024-020-02520-1 |
Skoumal R J, Brudzinski M R, Currie B S. 2016. An efficient repeating signal detector to investigate earthquake swarms. Journal of Geophysical Research: Solid Earth, 121(8): 5880-5897. DOI:10.1002/2016JB012981 |
Skoumal R J, Brudzinski M R, Currie B S, et al. 2020. Temporal patterns of induced seismicity in Oklahoma revealed from multi-station template matching. Journal of Seismology, 24(5): 921-935. DOI:10.1007/s10950-019-09864-9 |
Tan Y P, Cao J Q, Bian Z F, et al. 2014. Missing earthquakes detection and seismogenic structure of the Yuxian earthquake swarm in August of 2013. Acta Seismologica Sinica (in Chinese), 36(6): 1022-1031. DOI:10.3969/j.issn.0253-3782.2014.06.004 |
Tang J, Wang B S, Ge H K, et al. 2008. Study on weak signal detection of small shot in regional scale deep exploration. Chinese Journal of Geophysics (in Chinese), 51(6): 1810-1818. |
Uchida N, Bürgmann R. 2019. Repeating earthquakes. Annual Review of Earth and Planetary Sciences, 47: 305-332. DOI:10.1146/annurev-earth-053018-060119 |
Utsu T. 1999. Representation and analysis of the earthquake size distribution: a historical review and some new approaches. Pure and Applied Geophysics, 155(2): 509-535. DOI:10.1007/s000240050276 |
Utsu T. 2002. Statistical features of seismicity. International Geophysics, 81: 719-732. |
Van Trees H L. 1968. Detection Estimation and Modulation Theory. New York: John Wiley and Sons, Inc.
|
Vernon F, University California San Diego. 1982. ANZA Regional Network. International Federation of Digital Seismograph Networks. DOI:10.7914/SN/AZ |
Wessel P, Luis J F, Uieda L, et al. 2019. The generic mapping tools version 6. Geochemistry, Geophysics, Geosystems, 20(11): 5556-5564. DOI:10.1029/2019GC008515 |
Whalen A D. 1971. Detection of Signals in Noise. New York: Academic Press.
|
Yang H, Zhu L, Chu R. 2009. Fault-plane determination of the 18 April 2008 mount Carmel, Illinois, earthquake by detecting and relocating aftershocks. Bulletin of the Seismological Society of America, 99(6): 3413-3420. DOI:10.1785/0120090038 |
Yang H, Chu R S, Sheng M H. 2018. Source characterization of the 2015 gypsum mining earthquake in Pingyi, Shandong. Progress in Geophysics (in Chinese), 33(1): 125-132. DOI:10.6038/pg2018AA0509 |
Yang H F, Liu Y J, Wei M, et al. 2017. Induced earthquakes in the development of unconventional energy resources. Science China Earth Sciences, 60(9): 1632-1644. DOI:10.1007/s11430-017-9063-0 |
Yao D D, Peng Z G, Meng X F. 2015. Remotely triggered earthquakes in South-Central Tibet following the 2004 MW9.1 Sumatra and 2005 MW8.6 Nias earthquakes. Geophysical Journal International, 201(2): 543-551. DOI:10.1093/gji/ggv037 |
Yao D D, Walter J I, Meng X F, et al. 2017. Detailed spatiotemporal evolution of microseismicity and repeating earthquakes following the 2012 MW7.6 Nicoya earthquake. Journal of Geophysical Research: Solid Earth, 122(1): 524-542. DOI:10.1002/2016JB013632 |
Yao J Y, Tian D D, Lu Z, et al. 2018. Triggered seismicity after North Korea′s 3 September 2017 nuclear test. Seismological Research Letters, 89(6): 2085-2093. DOI:10.1785/0220180135 |
Yoon C E, O′Reilly O, Bergen K J, et al. 2015. Earthquake detection through computationally efficient similarity search. Science Advances, 1(11): e1501057. DOI:10.1126/sciadv.1501057 |
Zhang M, Tian D D, Wen L X. 2014. A new method for earthquake depth determination: stacking multiple-station autocorrelograms. Geophysical Journal International, 197(2): 1107-1116. DOI:10.1093/gji/ggu044 |
Zhang M, Wen L X. 2015a. Earthquake characteristics before eruptions of Japan′s Ontake volcano in 2007 and 2014. Geophysical Research Letters, 42(17): 6982-6988. DOI:10.1002/2015GL065165 |
Zhang M, Wen L X. 2015b. An effective method for small event detection: match and locate (M & L). Geophysical Journal International, 200(3): 1523-1537. DOI:10.1093/gji/ggu466 |
Zhao L S, Helmberger D V. 1994. Source estimation from broadband regional seismograms. Bulletin of the Seismological Society of America, 84(1): 91-104. |
Zhu L P, Helmberger D V. 1996. Advancement in source estimation techniques using broadband regional seismograms. Bulletin of the Seismological Society of America, 86(5): 1634-1641. |
Zhu L P, Rivera L A. 2002. A note on the dynamic and static displacements from a point source in multilayered media. Geophysical Journal International, 148(3): 619-627. DOI:10.1046/j.1365-246X.2002.01610.x |
Zhu L P, Ben-Zion Y. 2013. Parametrization of general seismic potency and moment tensors for source inversion of seismic waveform data. Geophysical Journal International, 194(2): 839-843. DOI:10.1093/gji/ggt137 |
Zhu W Q, Beroza G C. 2018. PhaseNet: a deep-neural-network-based seismic arrival-time picking method. Geophysical Journal International, 216(1): 261-273. DOI:10.1093/gji/ggy423 |
安艺敬一, 理查兹P G. 1986. 定量地震学理论和方法. 北京: 地震出版社: 59.
|
陈伟文, 倪四道, 汪贞杰, 等. 2012. 2010年高雄地震震源参数的近远震波形联合反演. 地球物理学报, 55(7): 2319-2328. DOI:10.6038/j.issn.0001-5733.2012.07.017 |
侯金欣, 王宝善. 2017. 2014年鲁甸MS6.5地震前后地震活动性. 地球物理学报, 60(4): 1446-1456. DOI:10.6038/cjg20170418 |
侯金欣, 谢凡, 任雅琼, 等. 2020. 利用模板匹配技术检测米尺度岩石断层黏滑实验中的声发射事件. 地球物理学报, 63(4): 1630-1641. DOI:10.6038/cjg2020N0102 |
贾佳, 王夫运, 吴庆举. 2019. 机器学习在地震检测与震相识别的应用综述. 地震工程学报, 41(6): 1419-1425. DOI:10.3969/j.issn.1000-0844.2019.06.1419 |
姜金钟, 陈棋福, 李姣. 2019. 中国东北的深源地震波形匹配检测及定位. 地球物理学报, 62(8): 2930-2945. DOI:10.6038/cjg2019M0210 |
李璐, 王宝善, 侯金欣. 2017. 模板匹配滤波技术在地震数据处理中的应用. 中国地震, 33(1): 14-22. DOI:10.3969/j.issn.1001-4683.2017.01.002 |
盛敏汉, 储日升, 危自根, 等. 2018. 四川省理县西山村滑坡运动变形过程中的微震研究. 地球物理学报, 61(1): 171-182. DOI:10.6038/cjg2018L0367 |
谭毅培, 曹井泉, 卞真付, 等. 2014. 2013年8月河北蔚县小震群遗漏地震检测与发震构造分析. 地震学报, 36(6): 1022-1031. DOI:10.3969/j.issn.0253-3782.2014.06.004 |
唐杰, 王宝善, 葛洪魁, 等. 2008. 小当量激发的远距离信号检测研究. 地球物理学报, 51(6): 1810-1818. DOI:10.3321/j.issn:0001-5733.2008.06.022 |
杨慧, 储日升, 盛敏汉. 2018. 2015山东平邑石膏矿塌陷地震震源参数测定. 地球物理学进展, 33(1): 125-132. DOI:10.6038/pg2018AA0509 |