CN106353821A - 一种微地震事件定位方法 - Google Patents

Abstract

本发明提供了一种微地震事件定位方法，该方法利用任意两道间的波形互相关，计算拾取初至时差，同时利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数,对微地震事件进行优化迭代定位。该方法可以有效的克服低信噪比同向轴不连续的问题；同时利用“利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数”，对于包含N个建波器的观测系统可以提供N（N-1）/2种检波器的组合，从而提供了大量的信息冗余，提高微地震事件的定位精度。

Description

一种微地震事件定位方法

方法领域

本发明涉及水力压裂微地震监测技术领域，特别是微地震事件定位方法。

背景技术

油气田的微地震监测一种主要用于油气田开发的新的地震方法，它是利用水力压裂、油气开采、常规注水等石油工程作业而产生的地震波。常规的井中微震事件定位前提是，具有较高的信噪比，微地震事件的同向轴具有良好的连续性，在此前提下利用微地震事件初至的匹配程度作为目标函数进行定位。这存在两个问题：1)初至拾取的人为误差问题2)信噪比较低时微地震事件同向轴没有良好的连续性是无法拾取初至。

发明内容

本发明的目的是针对现有技术存在的问题，提供一种微地震事件定位方法。该方法利用任意两道间的波形互相关，计算拾取初至时差，同时利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数,对微地震事件进行优化迭代定位。

本发明通过如下方法措施来实现：

一种微地震事件定位方法,首先确定工区的层状速度模型和检波器的位置,确定微地震事件信号，利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数,对微地震事件进行优化迭代定位。

上述方案中优化迭代定位过程包含如下步骤：

1)确定初始的微地震位置p₀(x₀,y₀,z₀)，和迭代终止阈值d_th；

2)在可能的微地震事件位置范围内应用一定网格划分空间；

3)根据初始位置p₀进行优化定位，其定位结果为p₁(x₁,y₁,z₁)；

4)计算初始位置p₀与优化位置p₁之间距离d₀₁，其计算公式如下

d₀₁＝[(x₁-x₀)²+(y₁-y₀)²+(z₁-z₀)²]^1/2 公式(1)；

5)如果d₀₁>d_th,令p₀＝p₁，转到步骤3)；否则结束迭代微地震事件位置p_e＝p₁；

上述方案进一步包括：

步骤3)中的优化定位过程包含以下步骤：

以初始的微地震位置p₀和划分的网格点为基础，应用网格搜索算法针对每一个空间网格点(x_s,y_s,z_s)计算目标函数；

选取目标函数最小的网格点的位置作为定位结果p₁。

进一步，假设包含N级检波器的观测系统，任意两级检波器对为R_i,R_j，所述的空间网格点(x_s,y_s,z_s)的目标函数建立包含如下步骤：

计算空间网格点(x_s,y_s,z_s)对应的两级检波器间的拾取初至的差为Δt_p ^ij和计算初至的差为Δt_cs ^ij；

网格点(x_s,y_s,z_s)的目标函数应用下式计算

$F(x_s,y_s,z_s)=\underset{i=1}{\overset{N}{\Σ}}\underset{j>i}{\overset{N}{\Σ}}{\mathrm{\λ}}_{ij}|{{\mathrm{\Δt}}_p}^{ij}-{{\mathrm{\Δt}}_{cs}}^{ij}|$ 公式(2)

其中λ_ij为检波器对R_i,R_j对应的权值。

进一步，检波器对R_i,R_j间的拾取初至的差Δt_p ^ij和互相关系数r_ij应用两道间的微地震记录进行互相关运算的方法求取，其计算步骤如下：

分别选取检波器对R_i,R_j的微地震信号分别为Sⁱ(h₁)和S^j(h₂)，其中1≤h₁≤n₁，1≤h₂≤n₂采样点Sⁱ(1)和S^j(1)对应的时刻分别为t^p _i和t^p _j，设微地震信号的采样间隔为Δt；

分别将信号Sⁱ(h₁)和S^j(h₂)延拓为S^I(h₃)和S^J(h₄)，其中1≤h₃≤n₁+n₂-1，1≤h₄≤n₁+n₂-1，且

S^I(h₃)＝Sⁱ(h₃) 1≤h₃≤n₁

公式(3)

S^I(h₃)＝0 n₁+1≤h₃≤n₁+n₂-1

S^J(h₄)＝0 1≤h₄≤n₁-1

公式(4)

S^J(h₄)＝Sⁱ(h₄) n₁≤h₄≤n₁+n₂-1

对Sⁱ(h₁)和S^j(h₂)进行互相关运算用下式表示

$\begin{array}{cc}{\mathrm{XR}}_{ij}\left(\mathrm{\τ}\right)=\underset{k=1}{\overset{\mathrm{\τ}}{\Σ}}{S}^I\left(k\right)S^I(n_1+n_2-\mathrm{\τ})& 1\≤\mathrm{\τ}\≤n_1+n_2-1\end{array}$ 公式(5)

考虑到微地震信号Sⁱ(h₁)和S^j(h₂)极性可能相反，对互相关函数取绝对值即：

XR_ij(τ)＝|XR_ij(τ)| 公式(6)

设τ＝τ_max时XR_ij(τ)取得最大值，则Δt_p ^ij按照下式计算

Δt_p ^ij＝t^p _j-t^p _i+Δt*(τ_max-n₁-n₂+1) j≥i

公式(7)

Δt_p ^ij＝t^p _i-t^p _j-Δt*(τ_max-n₁-n₂+1) j<i

两道间的最大互相关系数r_ij按照下式计算

$r_{ij}={\mathrm{XR}}_{ij}\left({\mathrm{\&tau;}}_{\max }\right)/\underset{k=1}{\overset{{\mathrm{\&tau;}}_{\max }}{\&Sigma;}}{S}^I{\left(k\right)}^2\underset{k=1}{\overset{{\mathrm{\&tau;}}_{\max }}{\&Sigma;}}{S}^J{(n_1+n_2-\mathrm{\&tau;})}^2$ 公式(8)

其中互相关系数r_ij代表检波器对R_i,R_j间的拾取初至的差Δt_p ^ij的可信程度。

空间网格点(x_s,y_s,z_s)对应检波器对R_i,R_j间的计算初至的差Δt_c ^ij应用以下步骤实现：

通过射线追踪方法空间网格点位置(x_s,y_s,z_s)到对应两级检波器对R_i,R_j的初至旅行时t_cs ⁱ和t_cs ^j；

应用以下公式计算检波器对间的计算初至的差

Δt_cs ^ij＝t_cs ⁱ-t_cs ^j j≥i

公式(9)

Δt_cs ^ij＝t_cs ⁱ-t_cs ^j j<i

检波器对R_i,R_j的权重系数λ_ij的确定过程包含如下步骤：

根据确定的层状速度模型，分别应用射线追踪方法计算微地震初始位置p₀到检波器对R_i,R_j位置的射线线路经Ray_i和Ray_j；

设射线路径Ray_i和Ray_j所经过的点分别为(pⁱ ₁,pⁱ ₂,pⁱ ₃..pⁱ _Ni)和(p^j ₁,p^j ₂,p^j ₃..p^j _Nj)，其中pⁱ ₁和p^j ₁微地震初始位置p₀，pⁱ _Ni和p^j _Nj分别为检波器对R_i,R_j所在的位置，(pⁱ ₂,pⁱ ₃..pⁱ _Ni-1)为射线路径Ray_i与层状介质的分界面的交点，(p^j ₂,p^j ₃..p^j _Nj-1)为射线路径Ray_j与层状介质的分界面的交点,则射线路径Ray_i和Ray_j在p₀位置处的射线方向用如下单位向量表示，

V_i＝[vⁱ _x vⁱ _y vⁱ _z]

公式(10)

V_j＝[v^j _x v^j _y v^j _z]

其中

vⁱ _x＝(xⁱ ₂-x₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2

vⁱ _y＝(yⁱ ₂-y₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2 公式(11)

vⁱ _z＝(zⁱ ₂-z₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2

v^j _x＝(x^j ₂-x₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2

v^j _y＝(y^j ₂-y₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2 公式(12)

v^j _z＝(z^j ₂-z₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2

其中(x₀,y₀,z₀)为p₀所对应的坐标位置,(xⁱ ₂，yⁱ ₂，zⁱ ₂)为p^j ₂点所对应的坐标位置,(x^j ₂，y^j ₂，z^j ₂)为p^j ₂点所对应的坐标位置；

检波器对R_i,R_j对应的射线线路径Ray_i和Ray_j在p₀处的射线差异可以用L_ΔVij表示，其计算公式如下

$L_{{\mathrm{\&Delta;V}}_{ij}}=\left|\right|{\mathrm{\&Delta;V}}_{ij}|{|}_2$ 公式(13)

其中

ΔV_ij＝V_i-V_j 公式(14)

||ΔV_ij||₂代表的是向量ΔV_ij的第二范数，

射线差异代表定位结果对R_i,R_j间的拾取初至差与计算初至差的差的误差的敏感程度，越大敏感程度越大，反之亦然；

检波器对R_i,R_j对应的权值λ_ij应用互相关系数r_ij与射线路径差异的模表示，用下式计算

${\mathrm{\&lambda;}}_{ij}=r_{ij}*L_{{\mathrm{\&Delta;V}}_{ij}}$ 公式(15)。

本发明利用任意两道间的波形互相关，计算拾取初至时差，同时利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数,对微地震事件进行优化迭代定位。该方法可以有效的克服低信噪比同向轴不连续的问题；同时利用“利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数”，对于包含N个建波器的观测系统可以提供N(N-1)/2种检波器的组合，从而提供了大量的信息冗余，提高微地震事件的定位精度。合理的权重系数选择对于定位结果至关重要，本发明中主要从以下两方面考虑一是：计算拾取初至时差时的最大互相关时的相关系数，用来体现该时差的可信程度，相关系数越高，可信度越高；二是：利用两级检波器的在“微震事件点位置处”射线路径的差异，射线路径的差异可以体现道间时差对定位结果的敏感程度，射线路径差异越大，敏感程度越低，权重越大。

附图说明

图1本发明处理流程图。

图2正演观测系统及定位过程和结果。

其中图2-a正演观测系统及定位过程和结果；标注1为检波器，标注2为网格扫描过过程中的第一个点(-200,200,2800)，标注3为可能的微地震区域，标注4为初始微地震位置p₀(100,100,2900)。

图2-b是图2-a震源点附近放大图；标注5为首次迭代p₀为(100,100,2900)时的定位结果(0，0，2990)，标注6为最终定位结果(0,0,2998)，标注7为事件点的真实位置，标注8为应用常规定位方法的定位结果(-12,10,3014)。

图3正演速度模型。

图4三种不同子波。

图5正演微地震记录及拾取的初至；标注9为拾取的初至。

图6利用互相关法得到的任意检波器对拾取初至的差。

图7利用互相关法计算任意检波器对拾取初至的差时的互相关系数。

图8网格点(-200,200,2800)对应的任意两级检波器间的计算初至的差。

图9相对于初始点位置p₀计算得到的任意两级检波器射线路径的差异。

图10相对于首次迭代初始点位置p₀为(100,100,2900)时计算得到的任意两级检波器权重系数。

图11初次迭代初始点位置p₀为(100,100,2900)时计算所有网格点的目标函数和最小目标函数值点；标注10为目标函数最小值点。

具体实施方式

下文特举出较佳实施例，对本发明作详细说明如下。

结合附图1，本实施例的具体流程为：

一种微地震事件定位方法,首先确定工区的层状速度模型和检波器的位置,确定微地震事件信号，利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数,对微地震事件进行优化迭代定位。其中，

优化迭代定位过程包含如下步骤：

1)确定初始的微地震位置p₀(x₀,y₀,z₀)，和迭代终止阈值d_th；

2)在可能的微地震事件位置范围内应用一定网格划分空间；

3)根据初始位置p₀进行优化定位，其定位结果为p₁(x₁,y₁,z₁)；

4)计算初始位置p₀与优化位置p₁之间距离d₀₁，其计算公式如下

d₀₁＝[(x₁-x₀)²+(y₁-y₀)²+(z₁-z₀)²]^1/2 公式(1)；

5)如果d₀₁>d_th,令p₀＝p₁，转到步骤3)；否则结束迭代微地震事件位置p_e＝p₁；

步骤3)中的优化定位过程包含以下步骤：

以初始的微地震位置p₀和划分的网格点为基础，应用网格搜索算法针对每一个空间网格点(x_s,y_s,z_s)计算目标函数；

选取目标函数最小的网格点的位置作为定位结果p₁。

进一步，假设包含N级检波器的观测系统，任意两级检波器对为R_i,R_j，所述的空间网格点(x_s,y_s,z_s)的目标函数建立包含如下步骤：

计算空间网格点(x_s,y_s,z_s)对应的两级检波器间的拾取初至的差为Δt_p ^ij和计算初至的差为Δt_cs ^ij；

网格点(x_s,y_s,z_s)的目标函数应用下式计算

$F(x_s,y_s,z_s)=\underset{i=1}{\overset{N}{\&Sigma;}}\underset{j<i}{\overset{N}{\&Sigma;}}{\mathrm{\&lambda;}}_{ij}|{{\mathrm{\&Delta;t}}_p}^{ij}-{{\mathrm{\&Delta;t}}_{cs}}^{ij}|$ 公式(2)

其中λ_ij为检波器对R_i,R_j对应的权值。

进一步，检波器对R_i,R_j间的拾取初至的差Δt_p ^ij和互相关系数r_ij应用两道间的微地震记录进行互相关运算的方法求取，其计算步骤如下：

分别选取检波器对R_i,R_j的微地震信号分别为Sⁱ(h₁)和S^j(h₂)，其中1≤h₁≤n₁，1≤h₂≤n₂采样点Sⁱ(1)和S^j(1)对应的时刻分别为t^p _i和t^p _j，设微地震信号的采样间隔为Δt；

分别将信号Sⁱ(h₁)和S^j(h₂)延拓为S^I(h₃)和S^J(h₄)，其中1≤h₃≤n₁+n₂-1，1≤h₄≤n₁+n₂-1，且

S^I(h₃)＝Sⁱ(h₃) 1≤h₃≤n₁

公式(3)

S^I(h₃)＝0 n₁+1≤h₃≤n₁+n₂-1

S^J(h₄)＝0 1≤h₄≤n₁-1

公式(4)

S^J(h₄)＝Sⁱ(h₄) n₁≤h₄≤n₁+n₂-1

对Sⁱ(h₁)和S^j(h₂)进行互相关运算用下式表示

$\begin{array}{cc}{\mathrm{XR}}_{ij}\left(\mathrm{\&tau;}\right)=\underset{k=1}{\overset{\mathrm{\&tau;}}{\&Sigma;}}{S}^I\left(k\right)S^I(n_1+n_2-\mathrm{\&tau;})& 1\&le;\mathrm{\&tau;}\&le;n_1+n_2-1\end{array}$ 公式(5)

考虑到微地震信号Sⁱ(h₁)和S^j(h₂)极性可能相反，对互相关函数取绝对值即：

XR_ij(τ)＝|XR_ij(τ)| 公式(6)

设τ＝τ_max时XR_ij(τ)取得最大值，则Δt_p ^ij按照下式计算

Δt_p ^ij＝t^p _j-t^p _i+Δt*(τ_max-n₁-n₂+1) j≥i

公式(7)

Δt_p ^ij＝t^p _i-t^p _j-Δt*(τ_max-n₁-n₂+1) j<i

两道间的最大互相关系数r_ij按照下式计算

$r_{ij}={\mathrm{XR}}_{ij}\left({\mathrm{\&tau;}}_{\max }\right)/\underset{k=1}{\overset{{\mathrm{\&tau;}}_{\max }}{\&Sigma;}}{S}^I{\left(k\right)}^2\underset{k=1}{\overset{{\mathrm{\&tau;}}_{\max }}{\&Sigma;}}{S}^J{(n_1+n_2-\mathrm{\&tau;})}^2$ 公式(8)

其中互相关系数r_ij代表检波器对R_i,R_j间的拾取初至的差Δt_p ^ij的可信程度。

空间网格点(x_s,y_s,z_s)对应检波器对R_i,R_j间的计算初至的差Δt_c ^ij应用以下步骤实现：

通过射线追踪方法空间网格点位置(x_s,y_s,z_s)到对应两级检波器对R_i,R_j的初至旅行时t_cs ⁱ和t_cs ^j；

应用以下公式计算检波器对间的计算初至的差

Δt_cs ^ij＝t_cs ⁱ-t_cs ^j j≥i

公式(9)

Δt_cs ^ij＝t_cs ⁱ-t_cs ^j j<i

检波器对R_i,R_j的权重系数λ_ij的确定过程包含如下步骤：

根据确定的层状速度模型，分别应用射线追踪方法计算微地震初始位置p₀到检波器对R_i,R_j位置的射线线路经Ray_i和Ray_j；

设射线路径Ray_i和Ray_j所经过的点分别为(pⁱ ₁,pⁱ ₂,pⁱ ₃..pⁱ _Ni)和(p^j ₁,p^j ₂,p^j ₃..p^j _Nj)，其中pⁱ ₁和p^j ₁微地震初始位置p₀，pⁱ _Ni和p^j _Nj分别为检波器对R_i,R_j所在的位置，(pⁱ ₂,pⁱ ₃..pⁱ _Ni-1)为射线路径Ray_i与层状介质的分界面的交点，(p^j ₂,p^j ₃..p^j _Nj-1)为射线路径Ray_j与层状介质的分界面的交点,则射线路径Ray_i和Ray_j在p₀位置处的射线方向用如下单位向量表示，

V_i＝[vⁱ _x vⁱ _y vⁱ _z]

公式(10)

V_j＝[v^j _x v^j _y v^j _z]

其中

vⁱ _x＝(xⁱ ₂-x₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2

vⁱ _y＝(yⁱ ₂-y₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2 公式(11)

vⁱ _z＝(zⁱ ₂-z₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2

v^j _x＝(x^j ₂-x₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2

v^j _y＝(y^j ₂-y₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2 公式(12)

v^j _z＝(z^j ₂-z₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2

其中(x₀,y₀,z₀)为p₀所对应的坐标位置,(xⁱ ₂，yⁱ ₂，zⁱ ₂)为p^j ₂点所对应的坐标位置,(x^j ₂，y^j ₂，z^j ₂)为p^j ₂点所对应的坐标位置；

检波器对R_i,R_j对应的射线线路径Ray_i和Ray_j在p₀处的射线差异可以用L_ΔVij表示，其计算公式如下

$L_{{\mathrm{\&Delta;V}}_{ij}}=\left|\right|{\mathrm{\&Delta;V}}_{ij}|{|}_2$ 公式(13)

其中

ΔV_ij＝V_i-V_j 公式(14)

||ΔV_ij||₂代表的是向量ΔV_ij的第二范数，

射线差异代表定位结果对R_i,R_j间的拾取初至差与计算初至差的差的误差的敏感程度，越大敏感程度越大，反之亦然；

检波器对R_i,R_j对应的权值λ_ij应用互相关系数r_ij与射线路径差异的模表示，用下式计算

${\mathrm{\&lambda;}}_{ij}=r_{ij}*L_{{\mathrm{\&Delta;V}}_{ij}}$ 公式(15)。

下面结合附图2-11，对上述实施例的实施方案做进一步说明。

建立如图2所示观测系统，应用如图3所示的速度模型，正演模拟过程中每级检波器随机选用不同的震源子波如图4所示，同时添加按照信噪比为1添加随机噪音最终获得如图5所示的微地震记录；

利用已知的速度模型(图3所示)和检波器位置(图2-a所示)对该微震记录进行重新定位

1)以偏离微震事件点173米的位置(100,100,2900)为初始位置p₀(如图2-a标注4)，设置迭代终止阈值d_th＝3米；

2)以震源位置为中心将其周围200*200*200三维空间(如图2-a标注3)作为可能的微地震事件范围，以2*2*2米的网格进行划分；

3)利用化相关法计算出任意检波器对R_i,R_j间的拾取初至的差Δt_p ^ij(如图6所示)，同时计算出任意检波器对互相关系数r_ij(如图7所示)；

4)根据射线追踪方法，利用公式9计出步骤1)划分网格中任意一点(x_s,y_s,z_s)相对于检波器对R_i,R_j间的计算初至的差Δt_cs ^ij(如图8所示为网格扫描过程中第一个网格点(-200,-200,2800)(如图2-a标注2)对应的所有检波器对计算初至的差)；

5)利用公式13计算任意检波器对R_i,R_j在p₀位置处的射线路径差的差异(如图9所示为首次迭代p₀为(100,100,2900)时，所有检波器对在p₀位置处的射线路径差的差异)；

6)根据公式15，利用步骤3)中计算的任意检波器对互相关系数和步骤5)中计算的计算任意检波器对R_i,R_j在p₀位置处的射线路径差的差异计算任意检波器对的权重系数λ_ij(如图10所示为首次迭代p₀为(100,100,2900)时，任意检波器对R_i,R_j的权重系数)。

7)根据公式2，利用步骤3)中计算出的任意检波器对R_i,R_j间的拾取初至的差Δt_p ^ij和步骤4)中计算出的划分网格中任意一点相对于检波器对R_i,R_j间的计算初至的差Δt_cs ^ij以及步骤6)中计算得到的任意检波器对R_i,R_j的权重系数λ_ij，计算出所有扫描网格点的目标函数值(如图11所示为首次迭代初始点位置p₀为(100,100,2900)时所有网格点的目标函数值)

8)选取步骤7中目标函数值中最小的点对应的位置为当前迭代定位结果p₁。(如图2-b所示标注5为首次迭代p₀为(100,100,2900)时的定位结果(0，0，2990))；

9)根据公式1计算初始位置p₀与当前迭定位结果p₁的距离d₀₁(首次迭代p₀为(100,100,2900)时d₀₁＝10米)；

10)如果d₀₁>d_th(首次迭代d₀₁(10米)>d_th(3米))，则令p₀＝p₁转到步骤5)继续计算；否则结束迭代过程，微震事件的定位结果p_ms＝p₁(本次定位共迭代3次，最终定位结果为p_ms为(0,0,2998)，如图2-b中标注6所示)；

同时根据微地震事件拾取的初至(如图5标注9所示)根据常规的震源定位方法(以计算初至和拾取初至的匹配程度作为目标函数利用网格搜索算法进行定位)计算得到的定位结果为(-12,10,3014)(如图2-b标注8所示)，本发明方法能够很好的克服初至拾取困难造成的定位不准确问题，提高了定位精度。

Claims (7)

1.一种微地震事件定位方法,其特征在于首先确定工区的层状速度模型和检波器的位置,确定微地震事件信号，然后利用所有检波器对间的拾取初至差与计算初至差间的差的加权叠加作为目标函数,对微地震事件进行优化迭代定位。

2.根据权利要求1所述的微地震事件定位方法,其特征在于，微地震事件的优化迭代定位过程包含如下步骤：

1)确定初始的微地震位置p₀(x₀,y₀,z₀)，和迭代终止阈值d_th；

2)在可能的微地震事件位置范围内应用一定网格划分空间；

3)根据初始位置p₀进行优化定位，其定位结果为p₁(x₁,y₁,z₁)；

4)计算初始位置p₀与优化位置p₁之间距离d₀₁，其计算公式如下

d₀₁＝[(x₁-x₀)²+(y₁-y₀)²+(z₁-z₀)²]^1/2 公式(1)；

5)如果d₀₁>d_th,令p₀＝p₁，转到步骤3)；否则结束迭代微地震事件位置p_e＝p₁。

3.根据权利要求2所述的微地震事件定位方法,其特征在于，步骤3)中的优化定位过程包含以下步骤：

以初始的微地震位置p₀和划分的网格点为基础，应用网格搜索算法针对每一个空间网格点(x_s,y_s,z_s)计算目标函数；

选取目标函数最小的网格点的位置作为定位结果p₁。

4.根据权利要求3所述的微地震事件定位方法,其特征在于，假设包含N级检波器的观测系统，任意两级检波器对为R_i,R_j，所述的空间网格点(x_s,y_s,z_s)的目标函数建立包含如下步骤：

计算空间网格点(x_s,y_s,z_s)对应的两级检波器间的拾取初至的差为Δt_p ^ij和计算初至的差为Δt_cs ^ij；

网格点(x_s,y_s,z_s)的目标函数应用下式计算

公式(2)

其中λ_ij为检波器对R_i,R_j对应的权值。

5.根据权利要求4所述的微地震事件定位方法,其特征在于，检波器对R_i,R_j间的拾取初至的差Δt_p ^ij和互相关系数r_ij应用两道间的微地震记录进行互相关运算的方法求取，其计算步骤如下：

分别选取检波器对R_i,R_j的微地震信号分别为Sⁱ(h₁)和S^j(h₂)，其中1≤h₁≤n₁，1≤h₂≤n₂采样点Sⁱ(1)和S^j(1)对应的时刻分别为t^p _i和t^p _j，设微地震信号的采样间隔为Δt；

分别将信号Sⁱ(h₁)和S^j(h₂)延拓为S^I(h₃)和S^J(h₄)，其中1≤h₃≤n₁+n₂-1，1≤h₄≤n₁+n₂-1，且

S^I(h₃)＝Sⁱ(h₃) 1≤h₃≤n₁

公式(3)

S^I(h₃)＝0 n₁+1≤h₃≤n₁+n₂-1

S^J(h₄)＝0 1≤h₄≤n₁-1

公式(4)

S^J(h₄)＝Sⁱ(h₄) n₁≤h₄≤n₁+n₂-1

对Sⁱ(h₁)和S^j(h₂)进行互相关运算用下式表示

公式(5)

考虑到微地震信号Sⁱ(h₁)和S^j(h₂)极性可能相反，对互相关函数取绝对值即：

XR_ij(τ)＝|XR_ij(τ)| 公式(6)

设τ＝τ_max时XR_ij(τ)取得最大值，则Δt_p ^ij按照下式计算

Δt_p ^ij＝t^p _j-t^p _i+Δt*(τ_max-n₁-n₂+1) j≥i

公式(7)

Δt_p ^ij＝t^p _i-t^p _j-Δt*(τ_max-n₁-n₂+1) j<i

两道间的最大互相关系数r_ij按照下式计算

公式(8)

其中互相关系数r_ij代表检波器对R_i,R_j间的拾取初至的差Δt_p ^ij的可信程度。

6.根据权利要求5所述的微地震事件定位方法,其特征在于，空间网格点(x_s,y_s,z_s)对应检波器对R_i,R_j间的计算初至的差Δt_c ^ij应用以下步骤实现：

通过射线追踪方法空间网格点位置(x_s,y_s,z_s)到对应两级检波器对R_i,R_j的初至旅行时t_cs ⁱ和t_cs ^j；

应用以下公式计算检波器对间的计算初至的差

Δt_cs ^ij＝t_cs ⁱ-t_cs ^j j≥i

公式(9)

Δt_cs ^ij＝t_cs ⁱ-t_cs ^j j<i。

7.根据权利要求6所述的微地震事件定位方法,其特征在于，检波器对R_i,R_j的权重系数λ_ij的确定过程包含如下步骤：

根据确定的层状速度模型，分别应用射线追踪方法计算微地震初始位置p₀到检波器对R_i,R_j位置的射线线路经Ray_i和Ray_j；

设射线路径Ray_i和Ray_j所经过的点分别为(pⁱ ₁,pⁱ ₂,pⁱ ₃..pⁱ _Ni)和(p^j ₁,p^j ₂,p^j ₃..p^j _Nj)，其中pⁱ ₁和p^j ₁微地震初始位置p₀，pⁱ _Ni和p^j _Nj分别为检波器对R_i,R_j所在的位置，(pⁱ ₂,pⁱ ₃..pⁱ _Ni-1)为射线路径Ray_i与层状介质的分界面的交点，(p^j ₂,p^j ₃..p^j _Nj-1)为射线路径Ray_j与层状介质的分界面的交点,则射线路径Ray_i和Ray_j在p₀位置处的射线方向用如下单位向量表示，

V_i＝[vⁱ _x vⁱ _y vⁱ _z]

公式(10)

V_j＝[v^j _x v^j _y v^j _z]

其中

vⁱ _x＝(xⁱ ₂-x₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2

vⁱ _y＝(yⁱ ₂-y₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2 公式(11)

vⁱ _z＝(zⁱ ₂-z₀)/[(xⁱ ₂-x₀)²+(yⁱ ₂-y₀)²+(zⁱ ₂-z₀)²]^1/2

v^j _x＝(x^j ₂-x₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2

v^j _y＝(y^j ₂-y₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2 公式(12)

v^j _z＝(z^j ₂-z₀)/[(x^j ₂-x₀)²+(y^j ₂-y₀)²+(z^j ₂-z₀)²]^1/2

其中(x₀,y₀,z₀)为p₀所对应的坐标位置,(xⁱ ₂，yⁱ ₂，zⁱ ₂)为p^j ₂点所对应的坐标位置,(x^j ₂，y^j ₂，z^j ₂)为p^j ₂点所对应的坐标位置；

检波器对R_i,R_j对应的射线线路径Ray_i和Ray_j在p₀处的射线差异可以用表示，其计算公式如下

公式(13)

其中

ΔV_ij＝V_i-V_j 公式(14)

||ΔV_ij||₂代表的是向量ΔV_ij的第二范数，

射线差异代表定位结果对R_i,R_j间的拾取初至差与计算初至差的差的误差的敏感程度，越大敏感程度越大，反之亦然；

检波器对R_i,R_j对应的权值λ_ij应用互相关系数r_ij与射线路径差异的模表示，用下式计算

公式(15)。