基于奇异值分解法的微震子台网监测能力分析

陈法兵 王 颖 任文涛 王传朋 王元杰 李 岩

1）中国北京 100013 煤炭科学研究总院开采设计研究分院

2）中国北京 100013 天地科技股份有限公司开采设计事业部

3）中国山东274718 郓城煤矿

4）中国山东272055 山东唐口煤业有限公司

0 引言

微震震源参数测定的精确性是合理分析微震数据和有效防治冲击地压灾害的基础。因此精确测定微震事件的震源位置和发震时刻是微震监测台网的首要任务。微震事件定位的精确性取决于多种因素，比如台站分布、速度模型、震相读取、走时区域异常、拾震器运行状态和环境噪音、台基情况等（焦远碧，1986；焦远碧等，1990；朱元清等，1997；姜长宁，2001；郭飙，2002），其中速度模型等因素是现阶段不可能完全解决的，震相读取误差具有一定的随机性。煤矿井下存在不规则的采空区和地质构造带，且采空区范围时刻处于动态变化当中，走时区域异常较为复杂，一般通过多点放炮试验得出不同源检距的平均波速，以此作为初始参数。进行拾震器选址时，要求候选点避开破碎岩体、大断层等。建造水泥台基，台基锚杆深度必须达到坚硬原岩。台站布置在硐室内，有效避开皮带机和转载机等大型用电设备的电磁干扰以及运输车辆和人类活动的噪音干扰。综上所述，只有测震台网布置是人为可以有效优化和改进的（许俊奇，1988 ；张玲等，2010）。本文使用奇异值分解法，从测震台网布置的角度对微震监测子台网的定位误差进行分析。

Kijko A（1977）、Nitzan 等（1990）、唐礼忠等（2006）、巩思园（2010）利用最佳实验设计理论中的D值优化法，得到震中和震源的定位误差。Yasuo Sato等（1965）通过蒙特卡洛法得到地震参数的误差值，从而计算出震源位置、起震时刻和波速的标准偏差。巩思园（2010）、唐礼忠等（2006）在分析测震台网监测能力时假设，对于任何一个微震事件，所有台站能接收到事件波形，但现实中小能量事件和强烈衰减的大能量事件可能无法触发所有台站甚至触发台站数达不到可以进行震源定位的最小台站数。对于某一震动事件，把能够接收到明显震动波形的台站组合称为子台网组合，即震源所能激发台站的集合，S（ S为全部台站的集合）。对于既定微震台网S，子台网组合是震源位置和震级的函数，即）。因此，针对不同位置、不同震级的微震事件，应该选择不同的子台网组合，对微震事件进行定位计算和误差分析（赵英萍等，2006；季爱东等，2011），因此，称为动态子台网组合。在合理确定子台网的基础上进行的定位误差分析，可精确反映现场实际情况，对现场微震监测具有指导意义。

1 定位子台网组合的判别方法

单台地震监测能力是指单台拾震器对于不同能量等级的微震事件所能监测的最大距离（单新建等，1996）。此处假设台站性能和读取震相的能力相同。单台地震监测能力决定了测震台网的密度大小和实际监测效能。

在假设煤矿地质和开采条件各个方向上无较大差别的情况下，可以近似把微震事件的能量和监测半径r联系起来（国家地震局地球物理研究所，1978），关系式为

式中，q接近2，μ与q都可以求出，不同煤矿取值不同。里氏震级公式为

式中，M为震级，E为能量，a、b为地震常数。联立式（1）、式（2）可得

按照波兰煤矿规程的规定，式中a=1.8，b=1.9。

对于任意给定震级和震源位置的微震事件，利用式（3）可以得到有效监测半径r，结合台站分布可唯一确定子台网组合。微震事件的震源位置和震级变化均会引起子台网组合的变化，这种动态的子台网组合合理反映了参与定位台站的实际情况。图1为唐口煤矿能量分别为1×103J，2×103J和5×103J的微震事件可激发的台站数目云图。

图1 不同能量等级微震事件的有效定位范围（a） E=103J；（b） E=2×103J；（c） E=5×103JFig.1 Effective location range of microseismic event with different energy

图2中菱形块折线代表子台网最大有效台站数，正方形块折线代表有效监测面积与矿井总面积的比值。微震事件的能量越大，有效监测半径越大，子台网包含的台站数目就越多。对同一能量的微震事件，发生在台网密度越大的区域，子台网包含的台站数目越多，对微震事件的监测能力越强。

2 地震参数的误差估计

图2 最大有效台站数和监测面积百分比Fig.2 Max of the number of effective station and monitoring area ratio

假设微震台网包含n个台站，各台站坐标 β=（xi, yi, zi），i=1,2, …, n；震源坐标为 （x0, y0, z0），微震事件震源的发震时刻为 t0；震源参数θ=（x0, y0, z0, t0）。当发生一个微震事件时，该台网记录到m个观测到时tj，j=1,2, …,m。

第i个台站理论走时Ti由式（4）表示

第i个台站的到时残差ri由式（5）表示

在震源参数近似值θ*=（x*0, y*0, z*0, t*0）处对走时Ti应用一阶Taylor展开式，观测到时ti可由式（6）表示

联立式（5）、式（6）可得

使用矩阵表示如下

式中，A为m×n维微震波走时偏导数矩阵，δθ为n×1维震源位置改正量，r为m×1维矢量，表示观测走时与计算走时的残差。因为可以进行微震事件震源定位的前提条件是必须最少有4个台站可以判读震相初至，所以n≥4。

对系数矩阵进行奇异值分解，见式（9）

式中U为m×m维矩阵，列向量为AAT的m个正交归一本征矢量。V为n×n维矩阵，其列向量为ATA的正交归一本征矢量。T表示矩阵转置运算。S是m×n维对角奇异值矩阵。

奇异值按递减排列，即λ1≥λ2≥…≥λn≥0，其中非零奇异值的数目等于A的秩。

估计参量空间的协方差矩阵可由下式推出

其中 E 为数学期望，A+为 A 的 Moore-Penrose 广义逆。A+=VS+UTS+=diag （λ+1, λ+2,…, λ+n） ，其中λ+

i=1/λi。C（r）为观测值的协方差矩阵，主对角元素是各观测值的到时方差（Lee et al，1984；朱介寿等，1988）。可得到估计参量空间的协方差矩阵H

H的主对角元素是震源参数标准误差的估计值。震中和震源深度的标准误差为

协方差矩阵H的非对角元素可以用于估计震源参数间的线性相关性，而较高的震源参数线性相关性会导致震源定位迭代过程失稳，对计算结果和收敛速度影响很大（田玥等，2002；陈祥熊，2007）。因此，震源参数的线性相关性可以用来评价震源定位结果的优劣。相关系数公式为

β 的置信区间： β2≤ t2/ （t2+n2-2），其中 t=t（1－α/2）是在 1-α/2 概率水平上有 n-2 个自由度的student分布值（赵仲和，1983）。β在（0，1）内取值，0代表毫不相关，1代表相关系数的最大值。

经度—深度

纬度—深度

图3所示为震级M=1时，经度—深度方向线性相关系数βxz，纬度—深度方向线性相关系数βyz的等值线图。由图可知，在工作面开采区域范围内，相关系数βxz和βyz都较小，基本都分布在0.1范围之内，最小值达到0.06，属极低度线性相关范畴。因为震源定位精度与线性相关系数成反比，因此基于奇异值分解技术的微震台网定位精度分析是可靠的。

取P波平均波速为4 000 m/s，到时读取方差为0.05 s。利用奇异值分解法计算参量空间协方差矩阵，分别得到震级为0.8和1时的震中和震源深度定位误差等值线图4和图5，等值线范围则由子台网组合判别方法确定，显然，M=0.8时有效监测范围并没有覆盖整个分析区域。

由图4可知，微震台网包络范围内震中定位误差较小，一般在0—20 m范围内，但是台网边缘震中定位误差逐渐增大。震源定位误差比震中定位误差大，在采掘工作面范围内约100 m以内。当密集子台网处于台网边缘时，会产生误差剧烈变化带，对台网内定位精度造成影响。由图5可知，震级为1时，台网的监测范围已经覆盖分析区域，工作面内震中定位误差在10 m以内，震源定位误差在40 m以内，定位精度显然比M=0.8时高。

图3 线性相关系数等值线（M=1）（a） βxz；（b） βyzFig.3 Contour map of inearly dependent coeffi cient（M=1）

图4 震中和震源深度定位误差等值线图（M=0.8）Fig.4 Error contour map of epicenter and focal depth（M=0.8）

图5 震中和震源深度定位误差云图（M=1）

3 现场放炮数据验证

在工作面开采过程中，曾多次放卸压炮和断顶炮，微震监测系统记录到这些爆破事件的波形，见图6。统计18个能量大于1×104J的典型放炮事件，见表1。对比实际的震中和震源误差与本文计算方法得出的震中和震源误差，已知微震监测系统震中误差允许范围为±20m，震源误差允许范围为±50m，则震中吻合率为77.8%，震源吻合率为94.4%。可知本文对定位精度的评价方法可信度较高，可以为现场提供可靠参考。

利用拾震器单台监测能力分析和奇异值分解技术相结合的方法对微震台网的定位动态子台网进行定位误差计算，可以预先估计台网的震中和震源定位精度，对台网布置方案的制定和优化具有重要的指导意义。

图6 爆破震动波形

表1 计算误差与实际误差对比Table 1 Contrast of calculation error and actual error

4 结论

（1）微震事件的定位子台网是随着微震事件的能量以及微震台网布置的变化而变化的，因此定位子台网是动态的。在台站分布已知情况下，通过对单台监测能力的研究可以确定定位子台网的组合。

（2）平面震源参数与垂直震源参数的线性相关系数较小，故震源参数相关性小，表明唐口煤矿微震监测台网布置较好，定位精度较高。由此可见奇异值分解法分析结果的可靠性。

（3）奇异值具有较好的稳定性，通过奇异值分解方法对定位方程组进行处理，可以排除干扰因素，得到稳定精确解。实际计算误差与奇异值分解法所得误差在误差允许范围内吻合率较高，证明利用奇异值分解法进行误差评估是有效的。

（4）本文所述基于微震定位动态子台网组合的判别和奇异值分解技术的地震参数定位精度估计，可以对微震台网的地震定位精度进行预评价。奇异值分解法对微震台网的建立及后续台站的增减和调整具有指导作用。

陈祥熊. 地震精确定位问题综述[J].福建地震，2007，23（1-2）.

巩思园，窦林名，曹安业，等. 煤矿微震监测台网优化布设研究[J]. 地球物理学报，2010，53（2）：457-465.

郭飙. 首都圈数字地震台网的微震定位实验[J]. 地震地质，2002，24（3）：453-460.

国家地震局地球物理研究所.近震分析[ M] .北京：地震出版社，1978.

季爱东，钟普浴，于澄，等，山东数字地震台网定位子台网选取对定位结果的影响[J]. 地震地磁观测与研究，2011，32 （4）：62-66.

姜长宁.江苏数字地震台网地震定位能力分析[J].地震学刊，2001，21（2）：1-5.

焦远碧，等.华北地区测震台网的监测能力[ J] .地震，1986，6（1）：26-31 .

焦远碧，吴开统，等.我国地震台网监测能力及台网观测条件质量评定[J]. 中国地震，1990，6（4）：1-7.

单新建，韩京，许静.新疆测震台网历史监测能力及现状[J].内陆地震，1996，10（1）.

唐礼忠，杨承祥，潘长良.大规模深井开采微震监测系统站网布置优化[J]．岩石力学与工程学报，2006，25（10）：2 036-2 042.

田玥，陈晓非.地震定位研究综述[J]. 地球物理学进展，2002，17（1）.

许俊奇，测震台网布局对震源参数精度的影响[J].地震地磁观测与研究，1988，3.

张玲，梁向军，董春丽等.山西数字地震台网监测能力分析[J].山西地震，2010，4（144）.

赵英萍，高景春，等. 定位子台网的选取对定位结果影响研究[J]. 地震，2006，26（4）：37-44.

赵仲和.区域地震台网地震定位能力分析[J].地震学报，1983，5（4）：467-476.

朱介寿，等.地震学中的计算方法[M].北京：地震出版社，1988：35-47.

朱元清，赵仲和.提高地震定位精度新方法的研究[J]. 地震地磁观测与研究，1997，18（5）.

Kijko A. An Algorithm for the Optimum Distribution of a Regional Seismic Network-I.Pageoph.，1977，（115）：999-1 009.

Kijko A. An Algorithm for the Optimum Distribution of a Regional SeismicⅡ[J]. An analysis of accuracy of location of local earthquake depending on the number of seismic station[J]. Pageoph.，1977，（115）：1011-1021.

Lee，S W 斯图尔特.微震台网的原理及应用[M]. 北京：地震出版社，1984：113-115.

Nitzan Rabinowitz，David M. Steinberg. Optimal confi guration of a seismographic network：a statistical approach[J].

Sato Y，Skoko D. Opimum distribution of seismic observation point Ⅱ[J].Bull of earthquake Res Inst Tokyo Univ，1965，（43）：451-457.