2012年9月7日11时19分(北京时间)，云南省昭通市彝良县与贵州省毕节市交界发生M_S5.7地震，震中为27.51°N，103.97°E，震源深度14 km；12时16分，彝良县又发生M_S5.6地震，震中为27.56°N，104.03°E，震源深度14 km(数据来自中国地震台网正式地震目录，CSN).据云南地震台网统计，截至10月7日24时，震源区记录到震级大于1.0的余震共944次，最大余震为M_L4.3级.据民政部统计，至2012年9月10日9时，地震共造成81人死亡，74.4万人受灾，直接经济损失35.1亿元(王东坡等，2013).

彝良县位于小江深大断裂东侧，区内发育一系列活动断裂，地震频繁.震源区位于昭通—鲁甸断裂和其前缘主断裂——彝良—会泽断裂之间(闻学泽等，2013)(吕坚和王小娜等学者称之为石门断裂)，在小江地震带影响区之内(图 1a).彝良地震是典型的双震，震源机制为走滑型，兼有少量逆冲分量(吕坚等，2013)，其余震分布在半径约为30 km的区域内.小江地震带历史上曾发生多次强烈地震，昭通断裂也曾发生多次地震.自1948年以来，震区100 km范围发生5.0级以上地震达30余次，比较典型的有1948年10月9日发生在贵州威宁的5.8级地震，距离本次地震震中约11 km；1948年10月10日发生在云南大关的5.8级地震，距离本次地震震中约23 km；其中震级最大的是1974年5月11日发生在大关的7.1级地震，距离本次地震震中78 km，造成1423人死亡，1600多人受伤的惨痛灾难(王东坡等，2013).

图 1

Fig. 1

图 1 研究区域断裂(a)和台站(b)分布图 XSHF： 鲜水河断裂； ANHF： 安宁河断裂； ZMHF： 则木河断裂； PDHF： 普渡河断裂； XJF： 小江断裂； LFF： 莲峰断裂； ZTF： 昭通—鲁甸断裂； YLF： 彝良—会泽断裂； XWF： 寻甸—宣威断裂.Fig. 1 Faults (a) and stations (b) distribution map of the study area XSHF： Xianshuihe fault； ANHF： Anninghe fault； ZMHF： Zemuhe fault； PDHF： Puduhe fault； XJF： Xiaojiang fault； LFF： Lianfeng fault； ZTF： Zhaotong-Ludian fault； YLF： Yiliang-Huize fault； XWF： Xundian-Xuanwei fault.

地震定位的精度会受到台网布局、可用定位震相、地震波到时的读数精度、地壳速度模型等诸多因素的影响.2000年，Waldhauser等在前人的工作基础上提出了双差地震定位法(Waldhauser and Ellsworth，2000，2002)，此方法采用走时差数据作为输入数据，并可使用互相关技术获得精确的走时差，减弱了地壳速度模型和到时读数精度对定位结果的影响，使其特别适用于发生在断层附近的大地震和其余震序列的重定位，定位误差小，可以得到震源的精确的时空分布.双差定位法提出后在国内外得到了广泛应用(杨智娴等，2003；黄媛等，2008；王未来等，2012，2014；房立华等，2013，2014；苏金蓉等，2013；张广伟等，2014；Waldhauser and Ellsworth，2002; Waldhauser and Schaff，2008; Waldhauser and Tolstoy，2011；Schaff and Waldhauser，2005；Lin et al.，2007；Hauksson et al.，2012；Matoza et al.，2013).

云南地震台网给出的地震目录显示彝良地震的余震分布散乱，没有明显的趋势(图 2).本文提出时域多通道相关检测函数，结合云南地震台网记录的地震波形数据，计算得到了高精度的波形互相关走时差数据；使用双差定位法对彝良地震和之后1个月的余震序列进行了重定位.定位结果为研究震源区发震构造，分析余震活动的时空特性提供了重要参考.

本文采用云南地震台网提供的2012年9月7日11时至10月7日22时，即彝良地震发生后1个月内的观测报告数据和地震事件3分量波形数据(采样率为100 Hz).选取观测报告中震级在M_L1.0级以上的944个地震，其中2次M_L5.0~5.9，2次M_L4.0~4.9，25次M_L3.0~3.9，915次M_L1.0~2.9(图 2).观测报告给出的中M_S5.7主震的震中为27.55°N，103.99°E，震源深度10 km；M_S5.6主震的震中为27.44°N，103.99°E，震源深度2 km(图 2中的五角星).两个主震的震源位置与中国地震台网正式地震目录(以下简称CSN)的结果(图 2a中的实心圆)有较大偏差.尤其是M_S5.6主震明显偏离地震群，且震源深度过浅.但为了保持数据的一致性，本文仍然采用云南地震台网观测报告中的结果作为初始位置进行重定位.

图 2

Fig. 2

图 2 初始地震分布图(断层构造闻学泽等，2013) (a)震中分布图，黑色圆圈表示云南台网定位结果，黑色五角星表示云南台网定位结果，黑色实心圆表示中国地震台网正式地震目录结果，黑色圆点表示ClobalCMT结果.AA′剖面与两主震震源机制走向一致，BB′剖面与AA′剖面垂直；(b)所有地震在剖面AA′上的投影；(c)所有地震在剖面BB′上的投影.Fig. 2 Initial location of earthquakes(Fault structure quote from Wen et al.，2013) (a)The epicenter distribution，black circles represent the results of Yunnan Seismic Network，black stars represent the results of Yunnan seismic network，black solid rounds represent the results of CSN，black dots represent the location of two main earthquakes comes from GlobalCMT，profile AA′ consistent with the strike of the two main focal mechanism，profile BB′ vertical with profile AA′；(b)All the earthquakes along the profile AA′；(c)All the earthquakes along the profile BB′.

GlobalCMT给出了彝良双震的最佳双力偶解(表 1)和矩张量解(表 2).本文将彝良双震的矩张量解绘制于图 2a中(采用下半球投影)，黑色五角星表示云南台网定位结果，黑色实心圆表示中国地震台网正式地震目录结果，由虚线和震源球连接.GlobalCMT给出的定位结果由实线和震源球连接.

表 1(Table 1)

表 1 GlobalCMT最佳双力偶解 Table 1 Best double-couple solution of GlobalCMT 震级 节面1 节面2 走向/(°) 倾角/(°) 滑动角/(°) 走向/(°) 倾角/(°) 滑动角/(°) MS5.7 236 63 136 350 51 35 MS5.6 234 58 155 338 69 34

表 1 GlobalCMT最佳双力偶解 Table 1 Best double-couple solution of GlobalCMT

表 2(Table 2)

表 2 GlobalCMT矩张量解 Table 2 Moment tensor solution of GlobalCMT 震级 Mrr Mtt Mpp Mrt Mrp Mtp 矩张量的指数部分/(N·m) MS5.7 1.700 0.694 -2.390 -0.972 -0.231 -1.780 1017 MS5.6 0.464 0.553 -1.020 -0.381 0.298 -0.464 1017

表 2 GlobalCMT矩张量解 Table 2 Moment tensor solution of GlobalCMT

我们绘制了地震事件震相走时曲线图(图 3)，实心点表示震相的观测走时，实线表示由速度模型M1(图 4)计算而得的理论走时.可见，当震中距大于250 km时，观测走时与理论走时相差较大，S波尤其明显，并且，约97%的观测走时数据的震中距都在250 km以内.因此，本文选择震中距小于250 km的走时资料，共11572条，其中P波5238条，S波6334条，每个地震均保证至少有4条走时记录.用于重定位的地震台站共32个，均来自云南台网.在震源区附近有6个近台(图 1b)，震中距均在100 km以内.近台数据能极大地提高定位精度，尤其能够提高震源深度的精度(张天中等，2007；黄媛等，2008).但由于地震发生在云贵交界，台站的整体分布在震源区东部有所欠缺.

本文采用双差定位法(Waldhauser and Ellsworth，2000)对彝良地震序列进行重定位.双差定位法是一种相对地震定位方法，反演的是一个地震群中的每个地震相对于该地震群的质心的相对位置，因而不需要主事件，适用于大地震群的空间跨度(杨智娴等，2003).双差定位法的优点主要体现在以下几个方面：(1)观测走时差(双差)可以用地震目录的震相到时数据计算，也可以用波形互相关技术得到，其中互相关技术获得的走时差可以达到时间序列采样的精度甚至更小，这也是双差定位法最大的优点.波形互相关技术的应用也为应用S波走时差精确定位注入了新的活力.对于十分相似的信号，采样率为100 Hz的数字化数据的P波和S波到时差可精确到1 ms(Poupinet et al.，1984).(2)双差定位法要求地震对到台站的距离须远远大于地震对中两地震之间的距离(一般要求距离比大于10).由于地震对到同一台站的路径相近，两者的观测走时差主要由地震对中两地震之间的地壳结构决定，不依赖于地震对与台站之间的地壳结构.通过引入地震对观测走时差，能够有效地减小由于对地壳结构了解不够精细而引起的误差.(3)将一定半径范围内的每个地震与其相邻的地震组成地震对建立连接.在定位过程中，每个地震对走时差方程共同组成双差方程组，联立求解，使得这些地震互相约束，并在反演迭代加权过程中舍弃与其他地震连接不紧密的地震，最终能够得到十分精确的相对位置.和传统地震定位法相比，双差定位法的定位结果比地震目录数据提高了一个量级(Waldhauser and Ellsworth，2000).

尽管双差定位法能够减弱地壳速度结构对定位的影响，能够得到震源间十分精确的相对位置，但速度结构仍然会影响地震群的绝对位置和尺度(黄媛等，2008).本次彝良地震以及2014年8月3日鲁甸地震之后，许多学者对这一区域进行了研究，我们在相关文献中查找到此区域较广泛使用的4种不同一维速度模型(图 4)：王未来等(2014)和房立华等(2014)用于鲁甸地震研究的速度模型，记为M1；赵小艳和孙楠(2014)用于鲁甸地震研究的速度模型，记为M2；吕坚等(2013)用于彝良地震研究的速度模型，基于朱介寿等(2005)的研究成果记为M3；赵旭等(2014)用于鲁甸地震研究的速度模型，基于王椿镛等(2002)的研究成果，记为M4.M1和M2都基于熊绍柏等(1993)的人工地震测深的研究成果.

图 4

Fig. 4

图 4 速度模型Fig. 4 Velocity model

为了寻找最适于数据和研究区域地质结构的速度模型，我们使用相同的数据和上述4个速度模型分别进行了地震重定位.我们从定位地震数量，地震群展布形态，主震位置，数据残差(RMS)和被定位到空气中的地震数量(AQ)等方面综合评价模型的优劣(表 3).表 3中，M4的结果最优，M1次之.但重定位地震分布显示，M4的主震位置与CSN的结果偏差较大，M1的主震位置与CSN的结果很接近，因此，最终我们选择M1用于本次彝良地震序列的双差定位研究.

表 3(Table 3)

表 3 4个模型的定位结果 Table 3 The location results of 4 models 模型名称 定位地震数 RMS/ms AQ M1 652 16 4 M2 616 23 10 M3 642 13 9 M4 661 13 0

表 3 4个模型的定位结果 Table 3 The location results of 4 models

本文在进行地震对匹配时，限定地震对最大震源距离为10 km(Waldhauser and Ellsworth，2000; Waldhauser and Schaff，2008；张广伟等，2014)，震中距最大为250 km.

选取研究区域内符合要求的地震，绘制并查看震相走时曲线(图 3)，去除明显偏离走时曲线的震相数据，制成特定格式的地震目录.对某一地震事件，在半径10 km的区域内寻找能够和它组成地震对的地震事件(称为邻居)，地震对中的相同震相的走时之差即为这个地震对的目录走时差，记为dt_cat.经过地震配对，我们共得到目录走时差数据1462360个，其中P波走时差数据572830个，S波走时差数据889530个.

图 3

Fig. 3

图 3 震相走时曲线(实心点表示观测走时，实线表示理论走时)Fig. 3 Travel time curve(the solid points represent observed travel times，the lines represent calculated travel times)

云南地震台网提供的地震波形数据是EVT格式的，我们首先将其转换为了SAC格式，并编写程序根据地震目录提取了相应的地震事件，将地震位置、发震时刻、观测走时等信息写入SAC头文件.对没有观测走时的地震事件SAC文件，我们用表 1中的速度模型计算了理论走时，再写入SAC头文件(Waldhauser and Schaff，2008).之后，对地震事件SAC文件依次进行了去均值，去线性趋势，去仪器响应和带通滤波，频带为1~6 Hz.检视地震波形，挑选信噪比高的数据进行下一步计算.

经预处理得到的地震事件SAC文件包括三个分量，分别为北向(N)，东向(E)和垂直向(Z).P波和S波在N，E，Z三个方向都存在分量，分量大小取决于地震和台站的方位角以及地震波的视出射角.在以往的研究中，一般用垂直分量(Z)计算P波互相关走时差，用水平分量(N、E)计算S波互相关走时差(Waldhauser and Ellsworth，2000；Matoza et al.，2013).同一震相不同分量的走时差和相关系数一般是不相同的，研究中选取大相关系数对应的走时差作为互相关走时差，相关系数作为权重.但是，相关系数只是对变量之间线性相关程度的定性描述，相关系数的小差别不能说明相关性一定存在明显差异.研究中，我们发现，用不同分量计算得到的相关系数大多差别不大，但走时差往往可以相差1 s以上，并且可能出现正负不一致的情况.同时，对于小震中距的浅源地震，P波在水平向的分量更大，不应只用垂直分量来计算互相关走时差.

为了解决上面的问题，本文提出时域多通道相关检测函数并用于计算波形互相关走时差数据(式(1)—(5)).

在时间域计算每个分量的相关检测函数C_m(式(1))，与传统的互相关函数不同的是，在两列数据错动的过程中，相关检测函数补充的是真实数据而不是0(Schaff and Waldhauser，2005).

计算归一化系数c_m:

将多个分量的相关检测函数和归一化系数叠加，得到多通道相关检测函数C(Yang et al.，2009):

函数C中的最大值(记为相关系数，用C_f表示)对应的值τ为相关改正量:

进一步地，通过地震和台站的经纬度可以计算地震到台站的方位角，进行坐标变换，将北向(N)和东向(E)分量转换为径向(R)和横向(T)分量.在层状均匀介质中，地震波在径向和垂直向所在平面内传播，P波能量主要集中在Z和R分量，SV波能量主要集中在Z和R分量，SH波能量主要集中在T分量(图 5).但是SV波会受到P波尾波的影响，导致初动不清晰.因此，本文用地震波的Z和R分量计算P波互相关走时差，用T分量计算S波的互相关走时差.

图 5

Fig. 5

图 5 BJT台站记录到的地震波形Fig. 5 The earthquake waveform recorded by station BJT

表 4和表 5给出了某个地震对的互相关走时差和相关系数计算结果，算法来自Schaff和Waldhauser(2005).这个地震对的间距为5.45 km，震级分别为M_S5.6和M_L4.3.ZAT台站到两个地震的距离分别为37.34 km和39.5 km，平均方位角约为223°；BJT台站到两个地震的距离分别为138.5 km和141.88 km，平均方位角约为105°.ZAT台站P波Z、N、R三个分量的走时差基本相同，平均相关系数大于0.9，R分量的相关系数大于Z分量的相关系数，说明小震中距的地震P波能量主要集中在水平分量；E，T两个分量的走时差与其他分量差别较大，同时相关系数也略小.由于震中距加大，数据的信噪比降低，BJT台站各分量的相关系数均小于ZAT台站各分量的相关系数.表 5中T分量的相关系数都在0.9以上，并且大于Z分量和R分量的相关系数，说明T分量的S波信噪比较高.

表 4(Table 4)

表 4 P波互相关走时差和相关系数表 Table 4 Cross-correlation travel-time difference and correlation coefficient of P wave 台站分量 ZAT BJT dtcc/s Cf dtcc/s Cf Z -1.42 0.8859 0.86 0.8728 N -1.44 0.9881 0.28 0.6950 E -2.47 0.8660 0.08 0.8044 R -1.44 0.9582 0.08 0.8253 T -2.18 0.7739 0.8 0.6802 本文 -1.44 0.9366 0.85 0.8607

表 4 P波互相关走时差和相关系数表 Table 4 Cross-correlation travel-time difference and correlation coefficient of P wave

表 5(Table 5)

表 5 S波互相关走时差和相关系数表 Table 5 Cross-correlation travel-time difference and correlation coefficient of S wave 台站分量 ZAT BJT dtcc/s Cf dtcc/s Cf Z -2.57 0.8109 -1.0 0.8003 N -2.61 0.8389 -1.51 0.7978 E -2.60 0.9495 -1.84 0.8874 R -3.01 0.7197 -1.84 0.7605 T -2.61 0.9208 -1.52 0.9126 本文 -2.61 0.9208 -1.52 0.9126

表 5 S波互相关走时差和相关系数表 Table 5 Cross-correlation travel-time difference and correlation coefficient of S wave

采用本文算法得到的BJT台站的P波走时差和相关系数与Z分量的结果接近，同时Z分量的相关系数在所有分量中也是最大的.两个台站不同分量的走时差都具有较好的一致性，但是我们也发现，BJT台站的P波走时差与S波走时差为一正一负，P波相关系数相比S波相关系数只有轻微的减小.ZAT台站距第一个地震较近，因此一般情况下，ZAT台站的走时差应为负值，这与表 4和表 5的结果相符合；BJT台站同样距第一个地震较近，走时差也应为负值.因此，BJT台站的P波走时差可能存在问题，需要进一步的检验(在后续的数据筛选中，BJT台站的P波走时差被去除了).上述结果说明了波形互相关走时差的计算结果严重依赖于数据的信噪比(由相关系数体现出来)，信噪比轻微的减小，就会使得走时差结果不稳定，因此必须对计算得到的互相关走时差数据进行筛选.

图 6是ZAT台站Z、R、T三分量的波形互相关改正图，分别对应P波和S波.两个地震波形先按震相对齐(实线代表震相，两虚线之间代表相关计算窗)，然后对虚线进行平移，移动量为相关改正量τ.平移后，两个地震震相的主波形基本重合了，说明本文的算法是正确的.

图 6

Fig. 6

图 6 ZAT台站波形互相关改正图 (实线表示MS5.6地震，虚线表示ML4.3地震)Fig. 6 Waveform cross-correlation correction of station STA (solid line represent MS5.6 earthquake，dotted line represent ML4.3 earthquake)

本文计算波形互相关走时差数据的具体步骤为：

(1)数据格式转化，提取事件波形，计算理论走时.

(2)数据预处理：去仪器响应，去均值，去线性趋势，带通滤波.

(3)人工挑选信噪比高的波形，分量变换.

(4)将相同台站记录的两个地震波形数据按震相(或理论震相)对齐，均放置滑动窗.

(5)设定合理的震相窗和滑动量.震相窗宽度应包含震相优势周期的2~3个波长，本文将P波窗设为P波到时前0.5 s至P波到时后0.5 s，将S波窗设为S波到时前1 s至S波到时后1 s.滑动量的大小是对观测走时的误差大小的估计(注意P和S震相不要混在一个时窗中)，滑动窗的初始位置为震相窗位置，本文设定了两个滑动量，分别为±1 s和±1.5 s(Waldhauser and Schaff，2008；Matoza et al.，2013).

(6)计算时域多通道相关检测函数.当相关系数C_f大于0.5时，保存相关改正量τ和波形互相关走时差dt_cc.由于设定了两个滑动量，会得到两组数据.

(7)若两组数据的相关改正量τ相等(同时两组数据的相关系数和互相关走时差也会相等)，则接受这组数据.如果互相关走时差是稳定的，那么两组数据的相关改正量应该相等或相差很小，这个方法可以判断互相关走时差是否正确(Waldhauser and Schaff，2008; Waldhauser and Tolstoy，2011).

目录走时差与相关改正量的和为波形互相关走时差，可以认为互相关走时差是在目录走时差的基础上进行了大小为相关改正量的误差修正而得到的.因此，在统计意义上，相关改正量τ应接近高斯分布(图 7).

图 7

Fig. 7

图 7 相关改正量柱状图(灰色表示原始数据，黑色表示经筛选的数据) (a)P波相关改正量； (b)S波相关改正量.Fig. 7 Histograms of correlative correction(gray represent original data，black represent selected data) (a)P wave correlative correction； (b)S wave correlative correction.

经过上述计算过程，我们共得到原始互相关走时差数据2436028个(图 7中以灰色表示)，其中P波走时差数据1217900个，S波走时差数据1218128个；经步骤(7)后，走时差数据剩余1804650个，占原始数据的74.08%，其中P波走时差数据824492个，S波走时差数据980158个.约26%的错误数据被去除了，保证了数据质量.

经过地震配对和走时差计算，我们得到了大量数据.如果将这些数据直接用于定位，一方面会使得双差方程组过大，增加计算量和计算时间；另一方面，数据中可能存在错误和自相矛盾的情况，会加大双差方程组的病态，导致得不到稳定解.因此，必须对原始的走时差数据进行筛选.

本研究中从以下几个方面对数据进行筛选：

(1)互相关数据的相关系数必须大于0.7.

研究发现，相关系数小于0.7的两个地震波形数据的相关性不高，它们的互相关走时差不可靠(Waldhauser and Schaff，2008).相关系数的大小与震相窗的宽度有关，因而并不是严格的，一般在0.6~0.8之间选取.如，Waldhauser和Schaff(2008)，Waldhauser和Tolstoy(2011)选用C_f≥0.7的数据，Lin等(2007)选用C_f≥0.6的数据.

(2)走时差数据(包括目录走时差和互相关走时差)必须小于最大预测走时差(式6)，

图 8

Fig. 8

图 8 走时差数据筛选过程图 (保留实线(P波)和虚线(S波)下方的数据)Fig. 8 The selection process of travel-time difference data (keep the data under the solid line for P wave and the dash line for S wave)

经过两次筛选，目录走时差数据剩余1414739个，互相关走时差剩余1518216个(图 7中以黑色表示).

(3)互相关数据为主，目录数据为辅，地震对走时差数据至少为8个.

在以往的双差定位法的应用研究中，同时使用目录走时差和互相关走时差数据进行定位(Waldhauser and Ellsworth，2000).若地震对在同一台站既有目录走时差也有互相关走时差，且两者不相等，这样的数据参与定位会使得双差方程组中存在矛盾的约束，从而造成解的不稳定.解决这个问题的一个方法是，在反演中对两种数据设置不同的距离参数(WDCC，WDCT)和权重(WTCCP，WTCCS，WTCTP，WTCTS)，将相对的大权重(一般相差100倍)分配给大间距地震对的目录数据和小间距地震对的互相关数据，即让目录数据决定地震群的大尺度分布，让互相关数据决定临近地震之间的相对位置.这种定位过程相当于先使用目录数据进行定位，然后以目录数据定位结果为震源的初始位置，再用互相关数据进行定位.因此，每个地震的最终位置同时受到目录数据和互相关数据影响.

假设观测走时数据的读数误差均为σ₀，那么目录数据的误差为 ；由于进行了互相关改正，互相关数据的误差约为一个采样间隔，即σ_cc≈0.01s << σ₀.因此，与互相关数据相比，目录数据误差大，精度低，但是易于获取，数据量大，对研究区域的覆盖好；互相关数据误差小，精度高，但计算费时，对数据质量要求高，因而数据量小，一般情况下对研究区域的覆盖没有目录数据好.

许多研究中只使用互相关数据进行定位(Lin et al.，2007；Waldhauser and Tolstoy，2011；Hauksson et al.，2012；Matoza et al.，2013).为了同时兼顾目录数据和互相关数据的优点，本文在数据选取时，以互相关数据为主，目录数据为辅：若地震对在同一台站既有目录数据也有互相关数据，则选用互相关数据；地震对须至少有5个互相关数据，这是为了保证两个地震波形足够相似，互相关走时差稳定可靠；若地震对的互相关数据少于8个，则用目录数据补充至8个.

⑷限定邻居数目，优先选取近距离的邻居.

将包含某一地震的所有地震对按地震对间距从小到大排列，优先选取小间距的地震对，并限定邻居数不超过100个.

走时差数据经过筛选，即保证数据的精度和一致性，也控制了数据量，在一定程度上保证了双差方程组的稳定.

地震目录中含有944个地震事件，经筛选后用于定位的数据集中共含有704个地震事件、55203个地震对、502530个走时差数据，其中P波走时差数据235229个，S波走时差数据267301个.数据集中共有475370个互相关数据，占所有互相关数据的41.4%，27160个目录数据，占所有目录数据的8.4%.互相关数据与目录数据的比例约为17.5∶1.地震对的平均距离为3.75 km，平均每个地震对拥有9个走时差数据，平均每个地震拥有78个邻居.

以相关系数作为互相关数据的权重，目录数据的权重P波设为0.4，S波设为0.2，统一作为互相关数据输入双差定位程序hypoDD中.反演算法采用LSQR，迭代分为5组，每组5次.如果在迭代中某个地震被定位到空气中(AQ)，则将这个地震的震源深度重新设为上一次迭代的深度，继续进行下一步迭代.

在3.3节中，要求地震对走时差数据至少为8个，与OBSCC的取值相同，因此，定位数据集中的704个地震事件仅形成一个地震群，这个地震群包含所有的地震.调整迭代参数，使得迭代结束时双差数据的标准差(RMS)基本不再变化以及平均震源改变量(DX，DY，DZ，DT)在预想的误差范围之内.RMS从迭代开始时的418 ms减小为迭代结束时的16 ms，并且在最后两次迭代时，RMSCC已经不再变化；DX，DY，DZ，DT分别从迭代开始时的1014 m，797 m，1812 m，207 ms减小为迭代结束时的6 m，4 m，160 m，1 ms.由于地震事件波形数据的采样间隔是0.01 s，用于反演的速度模型M1中的最小波速是3.18 km·s^-1，0.01 s的走时误差对应的距离误差应在几十米至几百米的量级.迭代结束时DX，DY，DZ，DT的值已在预想的误差范围之内.AQ数量从迭代开始时的20个减小为迭代结束时的4个，说明迭代越来越稳定.

通过上述几点，我们判断迭代结束时的定位结果是合理可信的，应为最优解.

采用双差法对704个地震事件进行重定位，最终重定位地震652个，重定位率达到92.6%.经重定位后，本次M_S5.7主震的震中为27.516°N，103.951°E，震源深度6.9 km；M_S5.6主震的震中为27.543°N，104.023°E，震源深度7.27 km.重定位后，主震的震中位置与CSN的结果几乎相同，震中位置相差小于2 km，震源深度有一定差别(表 6).值得一提的是M_S5.6主震的定位结果.重定位前，M_S5.6主震明显偏离地震群，且震源深度仅为2 km.重定位后，M_S5.6主震位于地震群内，震源深度加深，与M_S5.7主震的连线与震源机制走向(AA′)基本一致.两个主震的重定位结果说明本文得到的结果是可靠的.

表 6(Table 6)

表 6 彝良双震震源位置 Table 6 Yiliang earthquakes source location 震级 云南地震台网 中国地震台网 本文结果 MS5.7 27.553°N,103.993°E,10 km 27.51°N,103.97°E,14 km 27.516°N,103.951°E,6.9 km MS5.6 27.442°N,103.988°E,2 km 27.56°N,104.03°E,14 km 27.543°N,104.023°E,7.27 km

表 6 彝良双震震源位置 Table 6 Yiliang earthquakes source location

图 9为彝良地震余震序列分布图，与图 2比较可以看出，重定位前，地震群分布比较松散，趋势性不明显.M_S5.7主震位于地震群内部，M_S5.6主震位于地震群外.重定位后，地震群紧缩为条带状，并沿主震震源机制走向(约为235°)分布，勾勒出与附近断层(闻学泽等，2013)相近的展布形态.M_S5.7主震位于地震群左端，M_S5.6主震位于地震群右端，两者连线与震源机制走向(AA′)一致.水平分布上，余震向北东向扩展，整个扩展仅不到20 km，并沿BB′剖面有凸起，与东南侧的彝良—会泽断裂的形态对应.深度上，震源分布较重定位前变浅并且更加集中，震源主要分布在5~8 km之间，两个主震仅相差不到400 m，深度都在7 km左右，与吕坚等(2013)的研究结果接近.在垂直于主震震源机制走向的剖面上，地震群出现轻微倾斜，宽度约为8 km，与主震震源机制倾向(约60°)基本一致，余震集中分布在主震附近，且大致勾勒出断层的产状.地震群整体位于彝良—会泽断裂的上盘，即北西盘.

图 9

Fig. 9

图 9 重定位地震分布图.剖面和符号意义同图2Fig. 9 Relocation earthquakes distribution. The profiles and symbols are the same with that in Fig.2

双差定位程序给出的东西、南北及深度三方向平均相对误差分别为15.1 m、12.6 m和28.1 m.由于定位中采用LSQR算法，协方差矩阵的对角线元素只做了近似计算，并且严格依赖内部迭代时恰当的收敛，上述误差不具有实际意义，仅供参考(Paige and Saunders，1982；Waldhauser and Ellsworth，2000; Waldhauser，2001).

本文定位中使用的数据主要以互相关数据为主，许多研究中只使用地震目录数据，为了比较两者的区别，我们同样使用了目录数据进行定位，结果如图 10所示.

图 10

Fig. 10

图 10 目录数据定位地震分布图.剖面和符号意义同图2Fig. 10 Relocation earthquakes distribution of catalog data. The profiles and symbols are the same with that in Fig.2

只使用目录数据时，重定位地震数682个，迭代结束时RMS和AQ分别为107 ms和41.重定位的地震数量虽然大于以互相关数据为主时的652个，但RMS和AQ也分别大于本文的以互相关数据为主(以下简称互相关数据)时的16ms和4，说明目录数据虽然在数据量上占有优势，但数据精度和稳定性都要比互相关数据差.

M_S5.7主震的震中为27.507°N，103.952°E，震源深度21.54 km；M_S5.6主震的震中为27.551°N，104.035°E，震源深度14.27 km.两个主震的位置与CSN和本文的结果基本相同，相差不到2 km.M_S5.7主震的深度加深，越过了彝良—会泽断裂的断层面，M_S5.7主震的深度与CSN的结果一致，位于断层面上.地震群的震中分布较定位前有所紧缩，虽不如互相关数据的结果紧密，但形态基本一致，都沿BB′剖面有凸起；地震群的深度分布较定位前有所加深，部分震源越过了断层面.地震群在垂直于断层的剖面上有轻微的倾斜，这一特点与互相关数据的结果相同.

总体上说，目录数据的结果较原始地震分布紧密，而互相关数据的结果又较目录数据的结果紧密.两者的震中分布比较接近，深度分布相差较大，可能是因为在地震定位中深度误差一般要大于震中误差的缘故.

LSQR算法给出的东西、南北及深度三方向平均相对误差分别为160.2 m、146.5 m和321 m(同样仅供参考)，比以互相关数据为主的结果高出一个量级.

本文采用自助(bootstrap)法(Efron，1982；Billings，1994)来估计重定位结果的误差.以定位过程中迭代结束时的双差数据的加权标准差(RMS)16 ms为误差上限，对输入走时差数据加在[-0.016 s，0.016 s]间均匀分布的随机误差(对目录数据为[-0.107 s，0.107 s])，加误差数据的权重都设为1，用同样的反演参数进行定位(Waldhauser and Ellsworth，2000; Waldhauser and Schaff，2008; Waldhauser and Tolstoy，2011).上述过程重复200次，分析重定位次数大于160次的地震事件.计算同一地震多次定位结果的标准差，作为这个地震事件的定位误差，同时分析震源参数之间的相关性.

两种数据的定位误差结果如表 7和图 11所示，其中Catalog表示目录数据的结果，CC表示本文的结果.在分析目录数据定位结果的误差时，M_S5.7主震被重定位了118次，M_S5.6主震被重定位了49次，而对互相关数据，两个主震都被重定位了200次.这可能是导致两个主震的目录数据的定位误差要远远大于互相关数据的原因，同样也说明了目录数据的数据质量要差于互相关数据.对整个地震群，目录数据的平均定位误差是互相关数据的5~6倍，并且东西向的误差略大于南北向的，原因可能是在震源区东部缺少台站；深度误差是震中误差的3~4倍，说明地震定位中，震源深度的确定仍然是比较困难的.由于重定位次数过少，目录数据的两个主震误差不具有代表性；对互相关数据，M_S5.6主震的误差大于M_S5.7主震是因为其初始位置离地震群较远，数据上缺少约束.

表 7(Table 7)

表 7 两种数据和主震的定位误差 Table 7 The location error of two data and main earthquakes 震源参数数据 东西向/m 南北向/m 深度/m 发震时刻/s Catalog CC Catalog CC Catalog CC Catalog CC 平均误差 280.3 55.2 223.4 43.0 690.1 186.7 0.05 0.01 MS5.7 2640.5 44.6 2011.7 29.3 1739.1 111.3 0.36 0.01 MS5.6 883.0 164.7 1182.1 119.3 3514.0 475.7 0.23 0.03

表 7 两种数据和主震的定位误差 Table 7 The location error of two data and main earthquakes

图 11

Fig. 11

图 11 目录数据(a)和互相关数据(b)相对定位误差直方图Fig. 11 Histograms of relative location error of catalog data (a) and cross-correlation data(b) Fig.10 Histograms of relative location error

两种数据的震源参数相关性如图 12和图 13所示.目录数据的结果显示，东西向与南北向分量，深度分量和发震时刻之间都存在弱相关性；南北向与深度和发震时刻之间存在中等的相关性；深度与发震时刻之间存在强相关性.互相关数据的结果显示，经过互相关改正后，南北向与东西向分量和深度的相关性有所减弱；发震时刻与南北向分量和深度的相关性得到极大改善.东西向分量与深度仍然存在弱相关性，原因可能是在震源区东部缺少台站.震源参数相关性的结果表明，好的数据质量能够改善震中、深度和发震时刻三者之间的相关性，使得定位结果更加可靠.

图 12

Fig. 12

图 12 目录数据的震源参数相关性柱状图Fig. 12 Histograms of the source parameter correlation of catalog data

图 13

Fig. 13

图 13 互相关数据的震源参数相关性柱状图Fig. 13 Histograms of the source parameter correlation of cross-correlation data

上述误差和相关性的分析表明，本文的以互相关数据为主的定位结果要优于目录数据的定位结果，是可信的.

在地震定位过程中，如何精确确定震源深度一直是个难题.彝良地震发生后，不同研究机构给出了各自的震源深度.中国地震台网中心(CENC)确定震源深度均为14 km；美国国家地震信息中心(NEIC)确定震源深度均为10 km；GlobalCMT反演的深度均为12 km；吕坚等(2013)采用区域地震台网数据和CAP方法反演了两个主震的震源机制解和震源深度，并利用远震记录深度震相(P、pP、sP)进一步确定震源深度均为6 km；本文采用区域地震台网数据和双差定位法得到的深度分别为6.9 km和7.27 km，与吕坚等的结果基本一致，原因可能是采用的数据均为区域地震台网数据.

我们在选择合适的速度模型用于定位时注意到，不同的速度模型的定位结果，地震群的震中分布相差不大，但深度分布有明显的差别(图 14和图 15).

图 14

Fig. 14

图 14 速度模型M2的重定位地震分布图.剖面和符号意义同图2Fig. 14 Relocation earthquakes distribution of velocity model M2. The profiles and symbols are the same with that in Fig.2

图 15

Fig. 15

图 15 速度模型M3的重定位地震分布图.剖面和符号意义同图2Fig. 15 Relocation earthquakes distribution of velocity model M3. The profiles and symbols are the same with that in Fig.2

M1和M2都是基于熊绍柏等(1993)的人工地震测深结果.M3来源于天然地震面波层析成像的研究结果，是四川、云南和贵州交界地区的平均结果(朱介寿等，2005).比较M1、M2、M3三个模型，M1浅层速度最高，M2浅层速度最低，M3介于两者之间.在保持走时不变的前提下，波速减小会使得地震波路径变短，最终表现为震源抬升和震中距减小.M2、M3的结果相比M1震源都有明显的上升.M2在5 km处P波速度由5.1 km·s^-1变为6.1 km·s^-1，导致震源深度分布在约5 km处出现分界；M3在3 km处P波速度由5.169 km·s^-1变为5.996 km·s^-1，导致震源深度分布在约3 km处出现分界，并且M_S5.6主震因过于偏离地震群，在迭代中被舍弃了.M1在20 km以上速度变化比较和缓，震源深度分布很集中.在8~16 km的深度，M1、M2、M3速度接近.三者震源深度分布的共同特点是震源都分布在15 km以上，并表现出轻微的倾斜.三个模型定位结果表明，速度模型与震源深度的相关性要大于与震中的相关性.因此，在地震定位中，速度模型与震源深度有较大的相关性，应尽量采用与研究区域地质结构相近的速度模型(Michelini and Lomax，2004).

王小娜等(2014)研究发现，彝良地震序列震源深度多集中在15 km以上，彝良地区的速度结构在0~10 km范围内出现明显的高速异常，而高速体下方10~30 km则为低速异常，彝良M_S5.7和M_S5.6主震位于高速异常区两侧的速度过渡区.M1在24~32 km处存在明显的低速层，与王小娜等的研究结果一致.在迭代过程中，M1表现的更加稳定，最终的RMS和AQ数较小，重定位地震数量较多，主震位置正确，因此，本文认为速度模型M1的定位结果更加可信.

根据重定位结果(图 9)可以确定M_S5.7、M_S5.6地震震源机制解中的NE走向节面为发震面，本次地震序列的发震构造应为一条NE走向的活动断裂，因此位于震源区北部的昭通断裂和位于震源区南部的彝良—会泽断裂均可作为发震构造.

根据GlobalCMT的结果，M_S5.7、M_S5.6地震的震源机制解的平均走向为235°，倾角为60°，滑动角为146°，因此，其发震断层运动应以右旋走滑为主，兼有逆冲分量；断层倾向应为NW，上盘为北西盘，下盘为南东盘.这样的地震序列应主要位于地表可见断裂的北西侧.彝良—会泽断裂是一条走向NE，倾向NW的右旋走滑兼逆冲断裂(闻学泽等，2011，2013)，且彝良地震序列位于它的北西侧.因此，彝良—会泽断裂应为本次彝良地震序列的发震构造，这一点与吕坚等(2013)结论相同.

根据地震地质构造和地震分布特点，川滇地区可以分成16个主要地震活动带，本次彝良地震位于马边—大关—昭通地震带上(苏有锦等，2003).彝良地区西侧是以鲜水河断裂、安宁河断裂、则木河断裂和小江断裂为东界，以金沙江断裂和红河断裂为西界的“川滇菱形块体”，东侧是相对稳定的南华活动地块区(张培震等，2003；邓启东等，2002).自新生代以来，川滇块体向SSE向运动，并受到南华地块的阻挡(闻学泽等，2011，2013).附近区域现代构造应力场以水平作用为主，主压应力优势方位为SSE-SE，主张应力优势方位为NE-NEE(林向东，2009；谢富仁等，2004).彝良—会泽断裂位于块体边界带上，不断受到川滇块体SSE向的挤出运动，表现出以右旋走滑为主、兼挤压逆冲的运动特征，应力不断积累最终破裂导致了本次彝良地震序列的发生.

本文采用双差定位法和波形互相关数据对2012年9月7日彝良地震和余震序列进行重定位，得出以下结论.

(1)本文提出并使用时域多通道相关检测函数来计算波形互相关走时差数据，之后对数据进行了除错和筛选，以互相关数据为主进行定位.同时分析了目录数据和互相关数据定位结果的误差和震源参数相关性，表明互相关数据的结果要优于目录数据的结果，好的数据质量能够降低参数之间的相关性.

(2)M_S5.7主震的震中为27.516°N，103.951°E，震源深度6.9 km；M_S5.6主震的震中为27.543°N，104.023°E，震源深度7.27 km.地震群紧缩为条带状，并沿主震震源机制走向分布，趋势性明显；深度总体分布较重定位前变浅并且更加集中，集中分布在5~8 km之间，两个主震仅相差不到400 m.地震群出现轻微倾斜，与主震震源机制倾向基本一致，震源整体位于彝良—会泽断裂的北西侧.东西向、南北向、深度分量和发震时刻的平均相对定位误差分别为55.2 m，43.0 m，186.7 m和0.01 s.

(3)地震定位中，震源深度与速度模型存在较大的相关性，即震源深度的变化可以用速度模型的变换来补偿，这一点在研究中需要特别注意.

基于Pg、Sg震相的地震定位只有在近台数量足够和台站方位覆盖较好的情况下才能获得较高精度的震源深度(Mori，1991).对于较大地震(M≥5)可以利用远震记录中可能存在的pP或sP来确定震源深度.对于中小地震(M < 5)可以利用近震深度震相(sPg、sPmP、sPn)以及它们的参考震相(Pg、PmP、Pn)来获得高精度的震源深度(孙茁等，2014).同时应尽可能地采用适合研究区域的速度模型(Michelini and Lomax，2004).

(4)分析重定位后震源分布情况和主震震源机制解的结果以及区域现代构造应力场，确定彝良—会泽断裂为本次彝良地震序列的发震构造.