DOI: 10.11817/j.issn.1672-7207.2015.08.039

基于误差最小原理的微震震源参数反演

高永涛^{1, 2}，吴庆良^{1, 2}，吴顺川^{1, 2}，季毛伟^{1, 2}，杨凯^{1, 2}，刘超^{1, 2}

(1. 北京科技大学 金属矿山高效开采与安全教育部重点实验室，北京，100083；

2. 北京科技大学 土木与环境工程学院，北京，100083)

摘要：为减小矿山地质、岩性、施工等因素对震源参数测量精度的影响，提出一种基于误差最小原理的震源参数反演算法(均匀速度模型)，在目标函数取最小值条件下拟合出最优震源参数。研究结果表明：无论是对内场震源点还是对外场震源点进行定位分析，提出的反演方法比传统的定位算法具有更高的定位精度。克服传统方法由于现场实测波速与真实波速之间存在误差而影响到拟合精度的缺点，可满足实际工程监测需要，研究结果可为相关工程的微震准确定位提供参考和借鉴。

关键词：误差最小原理；定位算法；精度；微震；震源参数

中图分类号：TU46             文献标志码：A         文章编号：1672-7207(2015)08-3054-07

Source parameters inversion based on minimum error principle

GAO Yongtao^{1, 2}, WU Qingliang^{1, 2}, WU Shunchuan^{1, 2}, JI Maowei^{1, 2}, YANG Kai^{1, 2}, LIU Chao^{1, 2}

(1. State Key Laboratory of the Ministry of Education of China for High-Efficient Mining and Safety of Metal Mines,

University of Science and Technology Beijing, Beijing 100083, China;

2. School of Civil & Environment Engineering, University of Science and Technology Beijing, Beijing 100083, China)

Abstract: In order to reduce the impact of the mine geology, lithology, construction and other factors on the measurement accuracy of the source parameters, inversion algorithm was proposed based on minimum error principle using uniform velocity model, which fits the optimal source parameters through the search of the minimum value of the objective function. The results show that the proposed inversion method has a higher location accuracy than the traditional algorithm for both internal field source and external field source. The former overcomes the deficiency of the traditional methods for the fitting error resulting from the error between field measured velocity and true value, which can satisfy the field engineering demand. The research results may offer reference for location analysis in the similar projects using the micro-seismic observation system.

Key words: minimal error principle; location algorithm; accuracy; microseismic; source parameters

随着深部矿产资源的采掘和隧道工程的建设，地下工程不断走向深部，高应力条件下岩爆、地质缺陷活化等问题日渐突出，单纯靠传统的岩石力学手段已不能满足现实安全监测的需要，基于地球物理学发展起来的微震技术可以有效地监测岩石微破裂发生的位置，近年来，已广泛地应用于矿山、水电等地下工程监测。微震事件的空间分布规律包含着岩体不稳定区域的位置、范围和方向等大量信息，因此，提高微震事件的定位精度是岩体稳定性分析工作的重要内容之一，而有关定位算法的研究又是各项研究中的重点。国内外学者在震源定位算法研究方面做了大量工作并取得了丰硕成果^[1-6]。经典方法包括线性化地震定位方法，联合反演理论(JED)以及由JED算法发展而来的相对定位法(ATD)^[1-2]，在一定程度上都提高了定位精度。此外，还有一些著名的定位算法，如震源位置与速度结构的联合反演算法(SSH)、双重残差定位法(DDA)、震源扫描算法(SSA)等经典算法^[3-5]，至今仍被广泛采用。刘培洵等^[7]提出稳健估计方法，能够有效解决因个别随机误差偏大所引起精度下降的问题，得到较好的估计值。吕进国等^[8]基于单纯形法的优点，并结合模拟退火法的全局收敛性，提出采用稳健的模拟退火-单纯形法进行微震定位，提高了定位精度。朱权洁等^[9]采用 k-means 聚类方法对组合定位结果进行优化处理，减小定位结果的奇异性。董陇军等^[10-11]总结了前人的研究工作，提出3种无需预先测量速度的二维和三维震源定位的数学形式，按其因变量为到时、到时差、到时差商依次称为TT，TD和TDQ法。在实际工程中，由于地质构造的复杂性、岩体的各向异性、地下空区、采场、巷道的纵横交错、施工干扰等因素影响，造成实际参与定位计算的微震波波速、信号到时、检波器坐标等参数测量值与真实值存在误差(如检波器深孔安装时钻孔的角度和深度测量误差)，进而导致参与定位反演时的非线性方程组常数项并不是精确为零，严重影响了震源位置、发震时间等震源参数的拟合精度。为此，本文作者基于误差最小原理，提出一种微震震源参数反演算法，在目标函数取最小值条件下拟合出最优震源参数，以期为相关工程的微震准确定位提供参考。

1  基于误差最小原理的目标函数  构建

假设震源点的空间坐标为(x₀，y₀，z₀)，编号为i(i=1，2，…，n)的检波器空间坐标记为(x_i，y_i，z_i)，d_i为震源点与编号i(i=1，2，…，n)的检波器之间的距离，则

     (1)

记v₀为震源到各检波器之间岩层中的地震波传播速度(均匀速度模型)，t₀为震源产生的时刻，t_i为检波器接收到信号的时间，为第i个检波器到收到信号与震源产生的时间差，则

                (2)

即

  (3)

由于矿山现场地质地形复杂、资源开采时人类活动干扰等原因，严重影响了震源参数的测量精度，最终导致参与定位计算时的非线性方程组常数项不为0，本文将误差项记为λ_i(为不失一般性，此处λ_i可以取任意值)。

  (4)

若震源参数未知，则式(4)可简写为，显然，当|λ_i|越小时，求解的误差越小，反演出的震源参数越接近实际情况。为消除数值正负对求解结果的影响，构造目标函数如下：

         (5)

2  目标函数的求解与讨论

引入向量值函数F(X)，将式(5)转化成二次函数形式，定义函数

 (6)

则函数的极小点X^*就是目标函数(5)的最优解，

为求解X^*，令(为梯度散子)，由于

因此，

     (7)

采用Gauss-Newton或者Levenbery-Marguardt迭代法将函数F(X)线性化即可求解。

根据误差项λ_i的取值不同做如下讨论：

1) 当λ_i全部为0时，任选m个(m≤n，m为待求变量个数)方程参与定位计算，均可反演出震源参数的精确解。

2) 当λ_i不全部为0时(假设有k项不为0)，则取该k项构造目标函数并进行最优化计算，而将其他n-k项(n-k＜m，否则可取该n-k项参照情况1)进行反演)补充为约束条件，即

    (8)

约束条件为

此类问题也可用λ_i为0的n-k项对目标函数进行降维处理，进而提高求解精度。

3) 当λ_i全部不为0时，为提高求解精度，可先根据数据采集情况剔除误差较大项，然后构造目标函数，反演出震源参数的最优解。

此外，为提高计算速度，可根据现场实际情况补充约束条件，尤其是在情况2)和3)时，可以缩小搜索范围，加快寻优速度。约束条件形式如下：

s.t. ；

式中：，，，和分别为x₀，y₀，z₀，v₀和t₀的上下界。

3  算例分析

假定在空间区域O内，6个传感器位于球体的6个方位点，坐标(单位均为m)分别是A(0，500，500)，B(500，0，500)，C(1 000，500，500)，D(500，1 000，500)，E(500，500，1 000)，F(500，500，0)，假定微震源G(625，156，578)，H(232，604，746)位于内场区域，微震源I(673，700，1 200)，J(912，675，934)位于外场区域，波在介质中传播的等效波速v为5.2 m/ms，发震时间统一为某日10:00:20，传感器与震源空间分布位置如图1所示，正常情况下各传感器收到各震源P波的时间见表1。

图1  传感器与震源空间分布位置示意图

Fig. 1  Distribution of sensors and microseismic sources

表1  传感器收到各震源P波信号的时间

Table 1  P-arrival time of seismic source recorded by sensors                  ms

将本文提出的定位方法(MEP)与传统的“复合四四定位方法”(Four-Four Combination)^[12]的定位精度进行比较分析。一些矿山在微震监测系统安装完毕后采用校验炮测量波速，但由于地质构造复杂，岩体各项异性等，造成某一次校验炮所测波速并不能代表微震波在矿山各个方位的实际传播速度。因此，给微震P波速度一个较小的相对误差浮动(±1%，±3%，±5%)，即在波速分别取5.252，5.148，5.356，5.044，5.460，4.940 m/ms时，分析2种定位方法的定位结果，见表2。

从表2可以看出：采用MEP反演算法可准确求解出任意位置的震源参数(文中解的精度设定为e^-7)。采用传统的FFC算法，在已知区域的真实速度时也可求出任意位置的震源参数准确解，然而当给定速度与真实速度有误差浮动时，其定位误差随着给定波速误差的增大而增大，本算例中在x，y，z 3个方向的坐标误差最大达到53.386 6 m，而绝对距离误差最大达到53.722 7 m。由于地质构造的复杂多样、岩体的各向异性、地下空区、采场、巷道的纵横交错等因素的影响，在实际工程中很难准确测出某个监测区域的真实波速，甚至在某些复杂条件下，难以保证测量波速与真实波速的相对误差浮动控制在3%以内。本算例中当相对误差浮动为±3%时，传统方法定位的绝对距离误差最大达到31.932 m，以此说明采用传统方法进行定位时，其精度很难满足实际工程需要。

表2  采用新方法和传统方法分别对内外场震源点定位求解的误差比较

Table 2  Error comparisons of location results for both internal and external field source using new method and traditional method

图2所示分别为采用2种算法定位求解时，G，H，I和J 4个震源点的绝对距离误差随速度误差(量纲一)变化的对比曲线。从图2可以看出：采用传统的FFC算法，绝对距离误差与给定波速误差成正比，并随着其上下浮动而成近似对称分布。而采用MEP反演算法进行定位分析时，定位精度不随波速误差浮动而变化，只需保证算法中解的精度足够小即可。

图2  新方法和传统方法定位求解的绝对距离误差随速度误差变化曲线

Fig. 2  Comparison curves of positioning absolute distance error with velocity error by new method and traditional method

图3给出了分别采用2种算法定位求解时，绝对距离误差随检波器与震源点的平均距离变化的对比曲线。从图3可看出：采用MEP反演算法进行定位分析时，定位精度不受检波器到震源点的平均距离长短影响，同样只需保证算法中解的精度足够小即可，故而理论上，采用此方法定位求解时精度与震源点的内外场位置无关。而采用传统的FFC算法，无论是内场震源点还是外场震源点，绝对距离误差均与检波器到震源点的平均距离成正比。由于检波器到外场震源点的平均距离一定大于其到内场震源点的平均距离，故而为提高定位精度，要尽量确保微震点在检波器阵列之内。

图3  2种方法定位求解的绝对距离误差随检波器与震源点的平均距离变化曲线

Fig. 3  Comparison curves of positioning absolute distance error with average distance between detectors and microseismic sources using two methods

4  某磷矿现场爆破试验及结果分析

某磷矿矿体顶板直接分布在第四系中更新统含黏性土砾石含水层中，厚度较大，透水性较强，矿区资源开采过程中，顶板存在破裂和突水淹井的危险。为防止Q2含水层地下水涌入矿坑，降低开发风险，保护地质环境，合理开发利用矿产资源，并检验当前采矿方法对顶板含水层与隔水层的扰动影响情况，安装北京科技大学微震监测系统对其进行监测研究。系统配有12个传感器，主要目的是监测西区第2-12号勘探线425水平以上区域顶板和矿柱的稳定性。检波器及爆破点的坐标和时间参数见表3。

为验证上述理论分析的合理性，以该矿某2次爆破信号为例，根据检波器的坐标和相应的到时(见表3，已剔除信噪比不理想或者后波先达的数据)，采用本文提出的MEP反演算法和传统的FFC方法进行定位计算，并与实测位置比较分析。表4为2种方法的误差对比。

从表4可以看出，无论是内场震源点还是外场震源点，采用基于误差最小原理的反演算法的定位精度都比传统的“复合四四定位方法”的高，对于内场震源点，MEP算法的x，y，z 3个坐标方向的平均误差为6.74 m，最大为10.05 m，绝对距离误差为12.51 m；而FFC算法3个坐标轴方向的平均误差为13.20 m，最大达到17.58 m，绝对距离误差为25.13 m；对于外场震源点，MEP算法的x，y，z 3个坐标方向的平均误差为14.60 m，最大达到18.61 m，绝对距离误差为26.61 m；而FFC算法3个坐标轴方向的平均误差为17.59 m，最大的达到24.12 m，绝对距离误差为31.72 m，见图4。

表3  检波器及爆破点的坐标和时间参数

Table 3  Coordinate and time parameters of sensors and shot points

表4  新方法和传统方法定位误差比较

Table 4  Comparisons of location error using new method and traditional method

此外，从图4可以看出，无论是对于内场震源点还是外场震源点，采用MEP算法定位时地非线性方程组平均拟合误差都比传统的“复合四四定位方法”的小，对于外场震源点，MEP算法的平均拟合误差为2.41 m，最大为4.58 m；而FFC算法的平均拟合误差为4.76 m，最大为10.2 m；对于外场震源点，MEP算法的平均拟合误差为1.44 m，最大为2.91 m；而FFC算法的平均拟合误差为1.82 m，最大为2.92 m。

以上充分说明MEP算法比FFC法优越，定位精度较高，尽管MEP算法也是基于均质速度模型，但它是通过全局寻优，不断迭代搜索能最准确反应各传感器坐标、到时和波速之间关系的函数值，即使得震源参数之间的非线性拟合误差最小。克服了传统方法由于现场实测波速与真实波速之间存在误差而影响到拟合精度的缺点。理论上，采用MEP算法的震源拟合精度不受内外场空间位置影响，而传统方法的内场震源拟合精度应高于外场震源拟合精度。但本算例的实际结果与理论相悖，其原因主要是由于地质、岩性、施工等因素，微震波在传播过程中发生多次反射、透射，造成数据异常(如信噪比不理想或者远波先达等)，而本文在外场震源定位时已剔除大量异常数据，大大减小了误差。

图4  新方法和传统方法的定位误差和拟合误差比较

Fig. 4  Comparisons of location error and fitting error using new method and traditional method

以上算例分析及爆破试验均证实了新方法MEP的科学性、合理性及正确性，可以在实际工程中推广使用。

5  结论

1) 为减小因矿山现场地质地形复杂、资源开采时人类活动干扰等原因对震源参数测量精度的影响，提出一种基于误差最小原理的震源参数反演算法(均匀速度模型)，在目标函数取最小值条件下拟合出最优震源参数。

2) 本文提出的MEP反演方法比传统的FFC算法具有更高的的定位精度。理论上在保证求解精度条件下(文中解的精度设定为e^-7)，MEP算法的定位精度不随波速误差浮动而变化，且不受检波器到震源点的平均距离长短影响；而传统的FFC定位算法在给定速度与真实速度的误差浮动为±3%时，其精度难以满足实际工程需要，且定位结果的绝对距离误差与检波器到震源点的平均距离大小成正比。

3) 无论是内场震源点还是外场震源点进行定位分析，采用MEP反演算法的定位精度都高于传统FFC算法，但由于现场地质、岩性、施工等因素复杂，可能会影响微震波在地层中的传播路径，故而作者建议在微震监测台网布置时应尽量保证监测区域在检波器阵列之内，且在定位计算时应剔除异常数据。

参考文献：

[1] Douglas A. Joint epicenter determination[J]. Nature, 1976, 215: 45-48.

[2] Spence W. Relative epicenter determination using P-wavearrival- time differences[J]. Bulletin of the Seismological Society of American, 1980, 70(1): 171-183.

[3] Romney C. Seismic waves from the Dixie Valley Fairview Peak earthquakes[J]. Bulletin of the Seismological Society of American, 1957, 47(4): 301-319.

[4] Waldhauser F, Ellsworth W L. A double-difference earthquake location algorithm: Method and application to the Northern Hayward Fault, California[J]. Bulletin of the Seismological Society of American, 2000, 90(6): 1353-1368.

[5] Crosson R S. Crustal structure modeling of earthquake data, 1, Simultameous least squares estimation of hypocenter and velocity parameters[J]. Journal of Geophysical Research, 1976, 81(17): 3036-3046.

[6] 林峰, 李庶林, 薛云亮, 等. 基于不同初值的微震源定位方法[J]. 岩石力学与工程学报, 2010, 29(5): 996-1002.

LIN Feng, LI Shulin, XUE Yunliang, et al. Microseismic sources location methods based on different initial values[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(5): 996-1002.

[7] 刘培洵, 刘力强, 黄元敏, 等. 声发射定位的稳健算法[J]. 岩石力学与工程学报, 2009, 28(增刊1): 2760-2764.

LIU Peixun, LIU Liqiang, HUANG Yuanmin, et al. Robust arithmetic for acoustic emission location[J]. Chinese Journal of Rock Mechanics and Engineering, 2009, 28(Suppl 1): 2760-2764.

[8] 吕进国, 姜耀东, 赵毅鑫, 等. 基于稳健模拟退火-单纯形混合算法的微震定位研究[J]. 岩土力学, 2013, 34(8): 2195-2203.

L Jinguo, JIANG Yaodong, ZHAO Yixin, et al. Study of microseismic positioning based on steady simulated annealing-simplex hybrid algorithm[J]. Rock and Soil Mechanics, 2013, 34(8): 2195-2203.

[9] 朱权洁, 姜福兴, 缪华祥, 等. 基于聚类分析的微震定位二次优化研究[J]. 采矿与安全工程学报, 2014, 31(2): 196-202.

ZHU Quanjie, JIANG Fuxing, MIAO Huaxiang, et al. Quadratic optimization of mining microseismic source positioning based on cluster analysis[J]. Journal of Mining and Safety Engineering, 2014, 31(2): 196-202.

[10] 董陇军, 李夕兵, 唐礼忠, 等. 无需预先测速的微震震源定位的数学形式及震源参数确定[J]. 岩石力学与工程学报, 2011, 30(10): 2057-2067.

DONG Longjun, LI Xibing, TANG Lizhong, et al. Mathematical functions and parameters for microseismic source location without pre-measuring speed[J]. Chinese Journal of Rock Mechanics and Engineering, 2011, 30(10): 2057-2067.

[11] Li X B, Dong L J. Comparison of two methods in acoustic emission source location using four sensors without measuring sonic speed[J]. Sensor Letters, 2011, 9(5): 1501-1505.

[12] 李会义, 姜福兴, 杨淑华. 基于Matlab的岩层微地震破裂定位求解及其应用[J]. 煤炭学报, 2006, 31(2): 154-158.

LI Huiyi, JIANG Fuxing, YANG Shuhua. Research and application of microseismic monitoring location of strata fracturing based on Matlab[J]. Journal of China Coal Society, 2006, 31(2): 154-158.

(编辑  赵俊)

收稿日期：2014-08-08；修回日期：2014-11-18

基金项目(Foundation item)：长江学者和创新团队发展计划项目(IRT0950)；中央直属高校基本科研业务费资助项目(FRF-SD-12-002A)(Project (IRT0950) supported by Program for Changjiang Scholars and Innovative Research Team in University；Project (FRF-SD-12-002A) supported by the Fundamental Research Funds for the Central Universities)

通信作者：吴庆良，博士研究生，从事岩土和采矿工程方面的研究；E-mail：wuqingliangmy@163.com