地球物理学报  2013, Vol. 56 Issue (2): 522-530   PDF    
引用本文
王福昌, 万永革, 钱小仕 等. 由地震分布丛集性给出断层参数的一种新方法. 地球物理学报, 2013, 56(2): 522-530, doi: 10.6038/cjg20130216.  
Wang F C, Wan Y G, Qian X S, et al. A new method for determining fault planes parameters according to earthquake clustering. Chinese J. Geophys. (in Chinese), 2013, 56(2): 522-530, doi: 10.6038/cjg20130216.  
由地震分布丛集性给出断层参数的一种新方法
王福昌 , 万永革 , 钱小仕 , 张丽娟 , 张梅东     
中国地震局防灾科技学院, 河北三河 065201
收稿日期 2011-11-29, 2012-10-19 收修定稿
基金项目: 中央高校基本科研业务专项资金创新项目团队资助计划(ZY20110101)和河北省自然科学基金(A2011408006)资助
作者简介: 王福昌, 男, 1974年生, 副教授, 现主要从事多元统计分析和优化方法在地震学中应用研究.E-mail:fzmath@139.com
摘要: 由于大范围内地质构造的复杂性和介质的非均匀性, 发震断层面的几何形态一般十分复杂.如果大地震的破裂过程涉及多个断层的活动, 则发震断层并非是单一断层平面, 而是多个断层面的组合.利用地震空间位置分布丛集性, 即震源点成丛位于断层面附近的假设, 结合稳健扩充算法和主成分分析给出一种可以重构活断层网络三维空间结构的新方法.该方法每次从震源点集中处开始, 利用假设检验扩充子断层面, 并得到多个子断层面.接着按震源点属于最近断层面的准则把各子断层面内的震源点进行竞争, 并根据一定假设合并和删除一些子断层面, 最后用主成分分析确定每个子断层面参数.于是可根据地震事件目录给出一组矩形区域来描述断层面网络结构, 其中每个矩形断层面由其位置、走向和倾角确定.通过计算机模拟发现, 新方法可成功地重建模拟地震目录的断层面, 最后用于南加州1992年6月28日发生的Landers地震部分余震目录中, 得到各个子断层面参数与已知地质破裂或隐伏断层相当一致.
关键词: 断层面解      小震丛集      主成分分析      Landers地震     
A new method for determining fault planes parameters according to earthquake clustering
WANG Fu-Chang, WAN Yong-Ge, QIAN Xiao-Shi, ZHANG Li-Juan, ZHANG Mei-Dong     
Institute of Disaster Prevention of China Earthquake Administration, Hebei Sanhe 065201, China
Abstract: Because of the complex characteristics of the tectonic systems and medium heterogeneity in a wide region, the geometric shape of faults is usually complicated. If a major earthquake takes place, multiple faults are involved. And thus the fault planes are not a single fault plane but combined by several faults. We propose a new method that is able to reconstruct the three-dimensional structure of the active part of a fault network using the spatial location of earthquakes. Every time, the cluster consisting of neighboring several seismic events is selected randomly from the seismic events, and is gradually expanded by hypothesis testing, thus several clusters can be extracted. Subsequently, according to the criteria that any event point should be assigned to nearest cluster, the seismic events in extracted clusters are reassigned to different clusters. Some clusters are deleted or merged if they are too small or too close. Based on the remaining seismic event clusters, the parameters of fault planes are determined by the principal component analysis. So, given a catalog of seismic events, according to above procedure, we can find a series of rectangular-shaped fault planes which is characterized by its location, strike angle and dip angle to describe the spatial rupture structure of a large earthquake. The simulation results show that the algorithm successfully reconstructs the fault segments of synthetic earthquake catalogs. Applied to the real catalog constituted of a subset of the aftershock sequence of the 28 June 1992 Landers earthquake in southern California, the reconstructed plane segments fully agree with faults already known on geological maps or with blind faults..
Key words: Fault plane solution      Small earthquake clustering      Principal component analysis      Landers earthquake     
1 引言

确定某一地区地震发震断层面的位置、走向角、倾角等参数在地质构造和地球动力学研究中具有重要意义[1].对于浅源地震, 一般认为主震发生后, 大量余震在断层面及其附近发生, 因而震源点位置的空间分布可以约束断层面的形状和位置.假设用一个矩形平面来模拟发震断层, 且大多数地震发生在断层面的附近, 则可通过地震震源位置来求解发震断层的走向角、倾角及位置参数.王鸣和王培德(1992)[2]根据这种原理, 用设在极震区的一小孔径临时地震台网的资料分析研究了1989年10月18日大同-阳高地震部分余震的空间分布和综合断层面解.然而他们利用余震分布求断层面最小二乘解时使用了Gauss-Newton法, 该方法在迭代过程中对初始参数选择十分敏感, 而且未给出断层走向角和倾角的误差估计及断层面位置参数估计.王福昌等(2008)[3]基于L1范数给出能够减小离群值影响的稳健数学模型和相应的粒子群优化算法.王福昌等(2010)[4]用主成分分析思想和SVD分解方法求解断层面走向角和倾角, 避免了算法迭代过程, 加快了计算速度, 且可估计出矩形断层面四个顶点坐标位置, 但含有噪声震源点时算法不稳健, 即离群震源点存在可能严重影响参数估计值.上述方法的一种不足是都不能处理同时拥有多个子断层面的复杂破裂面.万永革等(2008)[1]在数学模型中引入定位误差信息并用模拟退火算法, 结合反演理论和应力场等地震学知识给出唐山地震几个主要断层面的走向、倾角和滑动角参数的估计值和误差, 每一个断层面分支的位置, 得到了较好的结果, 但每个子断层面划分仍不能自动给出.AnsariA等(2009)[5]使用模糊聚类对伊朗的历史地震目录进行自动模糊聚类分析, 与地质学上得到结果较为一致, 但没有利用深度信息, 没给出走向角、倾角和位置等信息, Ouillon G等(2008)[6]用Optimal Anisotropic Dynamic Clustering (OADC)方法对美国加州1992年Landers地震的断层面结构进行了较好的估计, 但是当数据包含异常值时不稳健, 断层面分支的个数确定和算法计算速度与稳定性等问题仍需进一步完善.张广伟等(2011)[7]也曾利用小震空间分布估算过唐山地震断层面解, 证实利用地震丛集性确定断层面参数方法的有效性.

一般来讲, 大多数地震均沿着活断层发生[8].近年来随着双差定位法[9]和波形互相关技术[10-12]在地震学中的应用, 对小地震的定位精度和数量越来越高, 使得采用小地震的丛集性确定深部活断层几何形态成为可能[13-14].然而, 上述研究者采用小震精确定位的丛集性描述断层几何形态并非定量描述, 因此迫切需要定量化方法来对深部断层进行更为精确的描述.

由于大范围内介质性质存在非均匀性和地质构造的多变性, 发震断层面的几何形态往往十分复杂, 并非是一个简单的平面.特别是, 如果一次特大地震的破裂过程涉及多个断层的活动, 则发震断层并非是单一的断层平面, 而是由多个断层面组合而成.本文利用稳健扩充算法, 使用震源点离断层面中心的马氏(Mahalanois)距离最小的准则把震源点进行分类, 以使得震源点隶属于最近的断层面, 在震源点集中的区域提取一系列可能断层面, 并将满足一定假设的断层面合并, 约简断层面个数, 最后利用主成分分析给出每个断层面走向角、倾角和矩形断层面四个顶点估计值.

先用计算机模拟验证本文方法的有效性, 再把该方法用于Landers地震目录中, 重建结果与已知文献结果比较吻合.

2 求断层面结构的数学模型 2.1 提取一个断层面的稳健扩充算法

地震震源点地理空间位置分布有时并不十分集中, 在提取断层面时, 这些散乱分布的震源点会严重影响断层面结构参数的确定.假设在断层面附近地震震源点是密集的, 而将散乱分布的震源点作为离群值(outliers).根据该假设给出一种稳健扩充算法, 只要震源点分布数据满足假设, 算法就能正确估计出模型参数.该算法主要思路是, 从一个密集子集出发, 逐渐扩充, 并在扩充过程中保持震源点分布集合的协方差矩阵的行列式值尽可能小, 即保持集合尽可能"紧密".因为从几何上看, 按照主成分分析原理, 正定协方差矩阵的行列式值正比于数据分布空间椭球体的体积, 行列式值越小, 含有同样多震源点的置信椭球体的体积越小, 从而数据分布越密集.

设D={xi, i=1, 2,…, n}为震源点的观测样本集, 其中xi为p维行向量, 将xi看作p维空间中一个点.在实际问题中, p取2或3.当p=2时对应只利用震源点的经纬度信息, 当p=3则再加入深度信息.

首先, 从某个初始位置x0开始, 寻找离x0最近的p+1个震源点, 构成初始子集S0.然后, 逐步扩充该子集.在每一步中, 寻找到子集Si的马氏距离最小的震源点, 并利用卡方检验验证它是否与当前集合中的点属于同一模型, 如果是, 则加入子集中, 继续下一步, 否则就停止扩充.

这里使用马氏距离而不用欧氏距离是因为虽然欧氏距离简单易求, 但它只能有效地描述球形分布, 而马氏距离可以比较好地刻画椭球形分布的高维正态分布.当椭球充分扁平时, 可以用它很好地近似表示平面.震源点x到集合Si的马氏距离定义为

其中, 分别为震源点集Si的均值和未修正协方差估计量, ni为Si中震源点数.

对于每个使用稳健扩充算法提取出来的震源点子集Si, 当p=2时利用主成分分析可以确定断层面走向和估计断层区域在地表的矩形投影, 当p=3时利用主成分分析可以确定断层面走向、倾角和矩形断层面四个顶点位置[4].

设z为任意震源点, 由于d2(z, Si)服从X2(p)分布, 设Xα2(p)为自由度为p的卡方分布X2(p)的α分位数, 若d2(z, Si) < Xα2(p), 则将z作为属于断层面的震源点, 添加到当前震源点子集Si中, 并进行下一步迭代, 否则就停止扩充震源点子集.

在震源点集合扩充的过程中, 特别是扩充的前几步, 由于估计的误差较大, 添加到震源点集合的数据点未必和该集合中的震源点属于同一断层面, 同时初始震源点子集中的数据点也未必属于同一断层面.这些错误划分的震源点会导致高估协方差矩阵行列式值.为提高方法的稳健性, 在每次迭代中, 使用Rousseeuw提出的C-Step方法[15]对当前数据集进行调整.其方法的思想为:设Si={x1, x2,…, xm}(m < n)为D的m元子集.将D中所有点按照到Si的马氏距离进行排序, 得到d(x'1, Si) < d(x'2, Si) <… < d(x'n, Si), 用排序后的前m个点构成新集合Sinew={x'1, x'2,…, x'm}来取代原来的子集Si.

为防止算法早熟, 对于前Nmin步中的震源点, 不做假设检验, 直接将最近(马氏距离最小)的震源点加入Si中.这里的Nmin是算法所能识别的断层面模型集Si应包含的最小震源点数, 可根据实际数据分布情形确定, 本文中Nmin取为30~50.

综上, 单个子断层面结构的稳健扩充算法详细步骤如下:

(1)设震源点x0为起始位置, 以距离x0最近的p+1个震源点作为初始子断层面数据集合S0, 令i=0;

(2)在第i步, 对子断层面数据集合Si={x1, x2,…, xni }, 计算其均值向量和协方差矩阵; 在p=3时, 利用最小特征值对应特征向量为子断层面法向量, 可以算出子断层面的走向角和倾角, 最小特征值对应震源点分布在子断层面法向量上的投影方差, 从而估计出断层面"厚度"; 若最小特征值相对于其余两个特征值太大, 则是由于定位精度太低或者不能被认为能构成一个子断层面;

(3)计算所有震源点到子断层面数据集合Si的马氏距离函数d(xk, Si), k=1, 2,…, n, 并由小到大进行排序, 得d(x'1, Si) < d(x'2, Si) <… < d(x'n, Si), 把Si更新为具有相同震源点数的更致密的子断层面数据集合S'i={x'1, x'2,…, x'ni};

(4)计算新子集S'i的均值向量μ'i和协方差矩阵Σ'i;

(5)设S'i的最邻近点为z, 满足, 计算z关于当前模型的马氏距离的平方值;

(6)或, 令, 转(2);

(7)否则停止迭代, 并输出最后的子断层面参数估计.

2.2 利用Mean Shift算法选取子断层面初始位置

尽管上述稳健扩充算法具有很强的抵抗离群值的能力, 但它仍为一种局部搜索算法, 对初始位置的选取有较强的依赖性.当初始点位置距离子断层面较远时, 算法仍有可能过高地估计断层面厚度而得出错误结果.为避免这种情况发生, 本文使用基于核密度估计的Mean Shift算法[16-17]过程, 可以保证初始点总是位于密度较大的区域, 即初始子集总是子断层面震源点数据, 于是算法可以更大概率得到正确结果.

给定p维空间的n个样本数据{xj, j=1, 2,…, n}, Mean Shift算法迭代公式为:

其中称为核函数, H为p×p维带宽矩阵.本文取核函数为高斯核函数, 带宽矩阵取对角矩阵, 其参数的估计利用Scott(1987)[18]公式:, σj是在第j个维度上数据分布的标准偏差, 使用如下估计得到.对应的Mean Shift算法可以描述为:

(1)设置初始值y0, 结束条件ε;

(2)计算yt+1=MH(yt);

(3)判断是否满足, 如果满足, 则退出; 否则, 用yt+1替代yt, 转至(2).

由于Mean Shift向量MH(x)正比于概率密度函数的梯度, 所以它总是指向概率密度增加最大的方向.因此Mean Shift算法思想是通过迭代使每个待处理的点"漂移"到概率密度函数的局部极大值点处.我们从样本集中随机地选择N0个数据点, 并对每个点运行Mean Shift算法程序直至收敛, 这样得到N0个模式点.然后从这些模式点中选择使核密度估计值达最大的点作为稳健扩充算法的起始点.

使用Mean Shift过程和稳健扩充算法可从震源点数据中提取子断层面并估计其参数.然而, 为了把上述算法推广到提取多个子断层面参数估计, 还需要做如下工作.

2.3 提取多个断层面及其参数估计

针对含有多断层面污染数据的多断层面估计过程有时会产生二义性.在没有其它先验知识的情况下, 多断层面参数估计的策略只能是基于震源点数据本身, 可使用统计假设检验来消除多个子断层面估计过程中可能产生的二义性.

多个子断层面参数估计的过程如图 1所示.

图 1 估计多个子断层面的过程 Fig. 1 Process of estimating fault segments

利用断层面厚度的比较, 可以确定某个断层面的合法性.设新提取的断层面Mnew含有vnew个震源点其断层面"厚度"为.设参考断层面的"厚度"为, 震源点数v0.计算如下统计量.

设Fα(vnew, v0)为自由度为(vnew, v0)的F分布的α分位数.如果T2≤Fα(vnew, v0), 则接受新断层面为合法断层面, 并在剩余的震源点集上继续提取.否则, 认为震源点集中不再包含有意义的合法断层面, 终止断层面提取.在实际地震数据处理中, 可以根据地震定位误差给出断层面的置信厚度, 从而简化检验过程.

每次提取一个新合法断层面后, 都在新断层面和已提取断层面之间进行竞争.设新断层面及其样本集为(Pnew, Snew), 已提取断层面及其震源点集为(Pk, Sk), k=1, 2,…, L.对任意数据点, 比较马氏距离和.将x分配给马氏距离最小的断层面.在模型竞争后, 检查所有断层面的震源点集, 若某个断层面所含震源点数为零或少于某一定值, 则认为该断层面不再有效, 并将其从已经提取的断层面中删除.

当所有合法断层面被抽取后, 还需进行断层面的合并, 基本思想是合并要满足两个条件:(1)两个属于子断层面的集合S 1, S 2有交集, 即存在满足或满足; (2)合并后方差不发生较大改变, 即对于任意两个断层面, 先将两震源点集合并, 估计合并后的震源点集的参数和方差, 将合并后断层面方差与原断层面方差作比较, 计算下面的统计量:

其中, σmerge2为合并后的方差, 即为合并后协方差矩阵的最小特征值, σ12, σ22分别为原来两个模型的方差.若一个断层面中心位于另一个断层面内, 且分别是原来两个断层面的震源点数目, 则将两个断层面合并, 否则不进行合并操作.完成所有合并操作后, 所得到的断层面集合就是最后的断层面结构估计的结果.

2.4 计算机模拟

采用计算机模拟来验证方法的有效性.

在三维空间中模拟产生三个子断层面, 假设x, y, z分别代表经度、纬度和深度并统一单位为km, 其方程分别为P1:y=-5, -10≤x≤10, 0≤z≤10, P2:y=5, -10≤x≤10, 0≤z≤10, P3:x=0, -10≤y≤10, 0≤z≤10.先在每个子断层面上随机取100个点模拟真实震源点, 然后围绕这些点进行随机扰动来模拟定位误差, 每个扰动分量服从均匀分布U(0, 0.1), 构成一个含有300震源点的测试样本集, 在图 2a中用黑点表示.同时为了检验算法抗干扰能力, 在空间区域内随机均匀产生300个噪声点, 在图 2b中用黑×"表示.使用稳健扩充算法、模型竞争和合并过程对这600个震源点数据进行断层面的提取.图 2c为按本文算法在含噪声数据中提取出来的子断层面及其震源点分布, 可见与图 2a有较好的对应.图 2d中三个矩形平面表示提取出来的三个子断层面在含噪声数据中的分布, 可以看出提取出来的子断层与原始的发震子断层非常接近, 算法很好地给出了子断层位置.

图 2 (a)模拟的断层面附近原始震源点; (b)在震源点混入噪声点; (c)算法提取出来的断层面及其震源点; (d)断层面在原始含噪声点数据中分布 Fig. 2 (a) Simulated hypocenters distribution around fault segments; (b) Mixing hypocenters distribution with outliers (c) Extracted fault segments and their hypocenters distribution; (d) Extracted fault segments in hypocenters distribution with outliers

还可以验证每一个提取出来的断层面与原始平面夹角都小于0.5°, 每个断层面的厚度都略小于0.02.算法从数据集中正确地拟合出三个断层面, 偏差可被定位误差所解释.模拟数据和得到的结果验证了本文方法的有效性, 考虑断层面的走向、倾角和长度等几何特征, 每种情形都准确地得到断层结构的空间范围.这里的测试数据是假定地震震源点均匀地分布在断层面附近, 为模拟更真实的地震目录, 将对非均匀分布的复杂断层结构进一步测试.

3 在Landers地震余震序列分析中的应用

由于新方法用于模拟数据中得到了正确结果, 促使我们把它应用到真实地震目录中, 例如1992年发生于美国南加州Landers地震余震序列数据.我们使用Lin G等(2007)[19]的地震定位数据, 出于计算速度考虑, 取5000个余震震源点.由于断层面的三维空间结构比较复杂, 这里只给出一系列拟合子断层面在地表上的投影.

3.1 余震数据处理

设(xi, yi, zi)为第i(i=1, 2,…, n)个余震震源点在北东下地理坐标系中的空间位置.每个精确定位余震震源点数据格式为东经(φ)、北纬(θ)和深度(h)组成, 记为(φi, θi, hi), 其中东经、北纬的单位是角度, 深度单位是km, 使用时需要转化为单位都是km的坐标(xi, yi, zi), 近似计算公式如下:

设地球上1弧度长近似为cl=111.199km, 纬度、经度和深度平均值分别为和, 则可得

3.2 子断层面参数的确定

对转换后的数据进行断层面提取与合并, 得到一系列子断层面.下面给出确定子断层面参数的方法.

令xj=(xj, yj, zj), j=1, 2,…, ni, 设Si={x1, x2,…, xni}为属于第i个断层面的震源点数据, 对Si进行主成分分析:计算, 把对角化后得3个特征值λ1≥λ2≥λ3≥0和相应单位特征向量u1, u2, u3.最大特征值λ1提供了断层面的最大长度, 对应的u1提供长度的方向,λ2和u2给出断层面宽度信息, λ3和u3给出断层面厚度信息.如果假设震源点在断层面附近,特征值和特征向量提供了断层面的长度和方位, 其中第三特征向量u3为断层面法向量.

特征值和断层面的长度及相应的方向之间的关系依赖于数据点在平面附近具体的分布.易证, 若ξ为均匀分布在区间[0, L]上的一维随机变量, 则ξ的方差为, 所以.因此, 如果地震均匀分布在一个长为L宽为W的平面附近, 则有,.第三个特征值λ3的平方根是地震定位标准差.如果断层面代表了断层和相关的地震, λ3应为定位的不确定性.程序以最自然的方式给出了类似于断层面的对象.对每一个断层面, 向量u3可确定断层面的走向角和倾角.

根据万永革等(2000)[20], 如果断层面的走向角为ψ(0≤ψ≤2π)、倾角为δ(0≤δ≤π/2), 则断层面的一个单位法向量为n0=(sinψsinδ, -cosψsinδ, cosδ), 于是u3//n0, 当u33 > 0时, u3=n0, u33 < 0时-u3=n0.若求出u3的第三个分量u33 < 0, -u3仍为λ3对应的特征向量.不妨假定u33 > 0, 则由u3=n0可确定走向角ψ, 倾角δ, 公式如下:

在求出断层面法向量u3后即可确定走向角和倾角, 均值μi即为断层面中心, λ1, λ2分别确定断层面长度与宽度, u1, u2确定断层面长度和宽度方向, 从而可确定矩形断层面四个顶点位置(空间坐标).

求出断层面顶点坐标后, 再把相应坐标分量x、y和z的值转化为φ、θ和h值.φ0、θ0和h0定义同上, 则可得θ=θ0+y/cl, φ=φ0+x/(clcosθr), h=z+h0.

3.3 数据分析结果

按照本文方法, 根据Lin G等(2007)[19]精确定位数据进行了分析, 拟合出Landers地震一系列子断层结果, 并通过绘制地表投影图识别出已知的断层面(图 3a).提取出的震源点用"·"表示, 离群震源点用"+"表示.提出离群点后, 对剩余震源点进行主成分分析来确定子断层面参数.得到的各矩形子断层面在地表投影为四边形, 用大写英文字母标识.可以看出这时断层面范围与余震分布吻合很好.图 3b给出剩余震源点及相应子断层面在地表投影, 为显示得更清楚, 图 3c给出本文方法确定断层面的地表投影.为便于比较, 在图 3d中给出文献[6]的结果, 通过对比可以看出, 两种方法得到的结果十分类似, 都较好地拟合出了符合实际的断层面.文献中原图即为图 3d, 单位与本文不一致, 但形状是类似的.本文方法还可识别一些未知的隐伏断层, 这对于地质构造研究是有意义的.

图 3 (a)Landers地震余震和11个拟合子断层面的地表投影; (b)剔除离群震源点后的结果; (c)11个拟合子断层面的地表投影; (d)文献[6]中得到的16个子断层面的结果 Fig. 3 (a) Surface projection of aftershocks and 11 fault segments of Landers earthquake; (b) Results after deleting outlier aftershocks; (c) Surface projection of 11 faults segments of Landers earthquake; (d) Results of 16 faults segments in Ref.[6]

本方法可能产生缺少任何构造解释意义的伪断层面, 产生的原因是有些地区地震分布不符合平面假设.在Landers余震序列中, 可以通过拟合断层面与已知的Landers地区断层面倾角大都接近于垂直的事实相比, 倾角太小的则为伪断层面.

表 1列出了本文算法给出的11个断层面的参数, 这些参数确定了断层面的大小、位置、走向角和倾角等.

表 1 由本文方法确定的Landers余震序列拟合平面参数及其与在同一区域已知断层的对应 Table 1 Correspondence between fitting planes parameters of the Landers aftershock sequences obtained by this method and kno、vn faults in the same area

通过图 3a和图 3d的对比可以看出图 3a的B和图 3b的C, 图 3a的E和图 3d的D对应较好.相对应断层面的参数如表 2、表 3.

表 2 由本文方法得到的部分结果 Table 2 The partial results by the new method
表 3 Ouillon等(2008)[6]得到的部分结果 Table 3 The partial results by Ouillon et al.(2008)[6]

本文得到的B和E断层面与Ouillon等(2008)的C和D相互对应, 通过对比表 2和表 3可以发现, 本文方法与文献结果比较接近, 说明了本文方法的有效性.

4 结论与讨论

本文提出了一种新方法, 可以把三维震源点数据划分为可以被解释为断层面的明显的类.在模拟含噪声地震数据中, 该方法以很高精度恢复了存在断层面, 于是将该方法用于Landers地震余震序列中, 并能够通过绘制地表映射图识别出已知的断层面.另外, 该方法还允许识别一些未知的隐伏断层.虽然地震学家可以通过肉眼对付同样的情况, 然而, 本方法的优点是可以自动地刻画断层面的特征和属于断层面的震源点而不需要人工指定它们.这种方法提供了断层的方位、地震序列的活动范围, 为约束断层面提供了一种新的约束[21].

相对于经度和纬度值而言, 深度值对断层面倾角估计有重要影响.由于在地震定位时, 深度值的精度相对较低, 从而导致确定出的走向精度较高, 倾角精度偏低.参数精度可以使用万永革等(2008)[1]的反演方法, 也可以使用Bootstrap方法来确定.

震源点空间分布决定了断层面的各个参数, 它的分布规律是否符合丛集性假设是应用本方法能否成功的关键.当数据假设不满足时, 参数精度会受到较大影响.例如当深度定位不够精确, 且接近均匀分布时, 得到的倾角往往较小, 甚至断层面几乎与地面平行.

考虑到地震预测和震灾评估, 大地震附近区域关于断层网络结构的知识也会有助于检验关于应力触发作用的不同假设[22-24].至今有两种方法用来检验主破裂断层几何形态和余震空间位置之间的联系.一是考虑余震所有平行于主断层的破裂断层.二是假定破裂结构是优化定位的.没有一种方法被证明是真实的 [25].因此准确地给出断层网络三维结构成为解释静态和动态地震触发和理解以几十年为尺度的断裂机制的重要参考.

参考文献
[1] 万永革, 沈正康, 刁桂苓, 等. 利用小震分布和区域应力场确定大震断层面参数方法及其在唐山地震序列中的应用. 地球物理学报 , 2008, 51(3): 793–804. Wan Y G, Shen Z K, Diao G L, et al. An algorithm of fault parameter determination using distribution of small earthquakes and parameters of regional stress field and its application to Tangshan earthquake sequence. Chinese J. Geophys (in Chinese) , 2008, 51(3): 793-804.
[2] 王鸣, 王培德. 1989年10月18日大同-阳高地震的震源机制和发震构造. 地震学报 , 1992, 14(4): 407–415. Wang M, Wang P D. The mechanism and seismogenic structure of Datong-Yanggao earthquake which occurred on Oct. 18, 1989.. Acta Seismologica Sinica (in Chinese) , 1992, 14(4): 407-415.
[3] 王福昌, 万永革, 胡顺田. 粒子群算法在主震断层面参数估计中的应用. 地震研究 , 2008, 31(2): 149–154. Wang F C, Wan Y G, Hu S T. Application of particle swarm optimization to the estimation of mainshock fault plane parameters. Journal of Seismological Research (in Chinese) , 2008, 31(2): 149-154.
[4] 王福昌, 曹慧荣, 万永革. 线性Errors-in-Variables模型在确定汶川地震主震断层面中的应用. 数理统计与管理 , 2010, 29(3): 381–390. Wang F C, Cao H R, Wan Y G. Application of linear errors-in-variables model for determination main earthquake's fault parameters of Wenchuan earthquake. Journal of Applied Statistics and Management (in Chinese) , 2010, 29(3): 381-390.
[5] Ansari A, Noorzad A, Zafarani H. Clustering analysis of the seismic catalog of Iran. Computer & Geosciences , 2009, 35(1): 475-486.
[6] Ouillon G, Ducorbier C, Sornette D. Automatic reconstruction of fault networks from seismicity catalogs:Three-dimensional optimal anisotropic dynamic clustering. J. Geophys. Res. , 2008, 113(B1). DOI:10.1029/2007JB005032
[7] 张广伟, 雷建设, 谢富仁, 等. 华北地区小震精定位及构造意义. 地震学报 , 2011, 33(6): 699–714. Zhang G W, Lei J S, Xie F R, et al. Precise relocation of small earthquakes occurred in North China and its tectonic implication. Acta Seismologica Sinica (in Chinese) , 2011, 33(6): 699-714.
[8] 杨智娴, 陈运泰, 张宏志. 张北-尚义地震序列的重新定位和发震构造. 地震学报 , 2002, 24(4): 366–377. Yang Z X, Chen Y T, Zhang H Z. Relocation and seismogenic structure of the 1998 Zhangbei-Shangyi earthquake sequence. Acta Seismologica Sinica (in Chinese) , 2002, 24(4): 366-377.
[9] Waldhause F, Ellsworth W L. A double difference earthquake location algorithm:method and application to the Northern Hayward fault, California. Bull. Seism. Soc. Amer. , 2002, 24(4): 366-377.
[10] Schaff P D, Waldhauser F. Waveform cross-correlation-based differential travel-time measurements at the Northern California Seismic Network. Bull. Seism. Soc. Amer. , 2005, 95(6): 2446-2461. DOI:10.1785/0120040221
[11] Lin G Q, Shearer P D, Hauksson E. Applying a three-dimensional velocity model, waveform cross correlation, and cluster analysis to locate southern California seismicity from 1981 to 2005. J. Geophys. Res. , 2007, 112: B12309. DOI:10.1029/2007JB004986
[12] 吕鹏, 丁志峰, 朱露培. 结合波形互相关的双差定位方法在2008年汶川地震余震序列中的应用. 地震学报 , 2011, 33(4): 407–419. Lü P, Ding Z F, Zhu L P. Application of double difference relocation technique to aftershocks of 2008 Wenchuan earthquake using waveform cross-correlation. Acta Seismologica Sinica (in Chinese) , 2011, 33(4): 407-419.
[13] 朱艾斓, 徐锡伟, 胡平, 等. 首都圈地区小震重新定位及其在地震构造研究中的应用. 地质评论 , 2005, 51(3): 268–274. Zhu A L, Xu X W, Hu P, et al. Relocation of small earthquakes in Beijing area and its implication to seismotectonics. Geological Review (in Chinese) , 2005, 51(3): 268-274.
[14] 朱艾斓, 徐锡伟, 周永胜, 等. 川西地区小震重新定位及其活动构造意义. 地球物理学报 , 2005, 48(3): 629–636. Zhu A L, Xu X W, Zhou Y S, et al. Relocation of small earthquakes in western Sichuan, China and its implications for active tectonics. Chinese J. Geophys (in Chinese) , 2005, 48(3): 629-636.
[15] Rousseeuw P J, Van Driessen K. A fast algorithm for the minimum covariance determinant estimator. Technometrics , 1999, 41(3): 212-223. DOI:10.1080/00401706.1999.10485670
[16] 李乡儒, 吴福朝, 胡占义. 均值漂移算法的收敛性. 软件学报 , 2005, 16(3): 365–374. Li X R, Wu F C, Hu Z Y. Convergence of a mean shift algorithm. Journal of Software (in Chinese) , 2005, 16(3): 365-374. DOI:10.1360/jos160365
[17] 周芳芳, 樊晓平, 叶榛. 均值漂移算法的研究与应用. 控制与决策 , 2007, 22(8): 841–847. Zhou F F, Fan X P, Ye Z. Mean shift research and applications. Control and Decision (in Chinese) , 2007, 22(8): 841-847.
[18] Scott D W, Terrell G R. Biased and unbiased cross validation in density estimation. Journal of the American Statistical Association , 1987, 82(400): 1131-1146. DOI:10.1080/01621459.1987.10478550
[19] Lin G. LSH:an earthquake relocation catalog using southern California pick and waveform data, (2007-5-30). http://www.rsmas.miami.edu/personal/glin/LSH.html.
[20] 万永革, 吴忠良, 周公威, 等. 根据震源的两个节面的走向角和倾角求滑动角. 地震地磁观测与研究 , 2000, 21(5): 26–30. Wan Y G, Wu Z L, Zhou G W, et al. How to get rake angle of the earthquake fault from known strike and dip of the two nodal planes. Seismological and Geomagnetic Observation and Research (in Chinese) , 2000, 21(5): 26-30.
[21] 方颖, 江在森. 断裂带变形分析方法研究. 大地测量与地球动力学 , 2011, 31(2): 53–55. Fang Y, Jiang Z S. On fault zone deformation analysis methods. Journal of Geodesy and Geodynamics (in Chinese) , 2011, 31(2): 53-55.
[22] Wan Y G, Wu Z L, Zhou G W, et al. Global test of seismic static stress triggering model. Acta Seismologica Sinica , 2002, 15(3): 318-332. DOI:10.1007/s11589-002-0065-3
[23] Wan Y G, Wu Z L, Zhou G W. Focal mechanism dependence of static stress triggering of earthquakes. Tectonophysics , 2004, 390(1-4): 235-243. DOI:10.1016/j.tecto.2004.03.028
[24] Wan Y G, Shen Z K. Coulomb failure stress changes on faults caused by the 2008Mw7.9 Wenchuan earthquake, China. Tectonophysics , 2010, 491(1-4): 105-118. DOI:10.1016/j.tecto.2010.03.017
[25] Steacy S, Nalbant S S, McCloskey J, et al. Onto what planes should Coulomb stress perturbations be resolved?. J. Geophys. Res , 2005, 110: B05S15. DOI:10.1029/2004JB003356