CN110261900A - 一种基于速度信息的地下浅层微震定位系统 - Google Patents

Abstract

本发明公开了一种基于速度信息的地下浅层微震定位系统，包括震动传感器阵列、信号预处理模块、定位参数提取模块、定位模型建模模块和定位模型解算模块；本发明采用地下震动场波场互易原理，通过传感器阵列数据进行波场逆时数值模拟，将地层速度场信息进行网格划分，利用交错网格有限差分法实现了时‑空场全波形逆时外推，生成了传感器作为“震源”下，群波阵面三维空间叠加图，建立了含有震源信息的定位模型。本发明建立了更加精细化的定位模型，采用能量聚焦的方法，在建立目标函数的基础上，利用QPSO在震源模型中进行全局快速搜索，实现了大范围内快速震源解算。实现了大区域、多尺度、多线程的快速精准扫描定位。

Description

一种基于速度信息的地下浅层微震定位系统

技术领域

本发明属于爆破震动测试领域，具体涉及一种基于速度信息的地下浅层微震定位系统。

背景技术

地下浅层震源是指地下震源深度不超过100m的空间内发生的事件，其定位技术是解决军事领域中高价值弹药地下炸点定位以及侵彻轨迹测量等主要途径；是实现民用领域中地质监测，工程爆破、文物防盗监测、煤矿勘察、地表结构组成成分分析、地质结构探索、地下稀有矿物质勘探、石油探测挖掘等重要手段。

与深层地震、煤矿深层开采、石油勘探等大区域、大当量、大深度、长时间震源定位相比较，该类震源定位具有以下特点：(1)观测数据有限，在未知地层模型前提下，仅具有单次爆炸产生的一次震动数据；(2)地下震源深度较浅，一般不超过100m，浅层地质结构复杂且未知，无法借鉴深层地壳结构模型建立浅层速度场模型；(3)震动波波群混叠复杂，爆炸近场土壤本构特性呈弹塑性，弹性波受地面反射、折射影响较大，震相特征不明显；(4)定位精度要求较高，在100m的散布区域内立体定位误差小于1m，不同于大深度地震定位中百米定位误差要求，属于小区域高精度定位问题。

目前，在地下空间定位中，主要针对深层地震定位系统，没有直接用于所述浅层小区域定位系统，而深层地震定位系统无法直接应用到浅层定位中，主要存在如下问题：

(1)深层地震定位主要方法之一——基于Geiger的走时定位方法，该方法重点解决在地质均匀且速度已知条件下的起爆点定位问题，但对于浅层复杂地质条件，速度变化较大，在采用该方法时，速度参量仅能设定为平均值，定位误差较大；

(2)深层地震定位主要方法之二——偏振角度定位方法，该方法采用少量观测节点即可实现起爆点快速定位，但在复杂地质条件下，特别存在强反射界面时，表层入射角度无法表征震源和传感器接收点之间的真实射线路径，因此存在定位假象。

发明内容

本发明的目的是提供一种基于速度信息的地下浅层微震定位系统，解决了地下浅层震源无法定位的问题，还解决了地下浅层震源无法精细化定位的问题。

本发明的技术方案为：一种基于速度信息的地下浅层微震定位系统，其特征在于：包括震动传感器阵列、信号预处理模块、定位参数提取模块、定位模型建模模块和定位模型解算模块；震动传感器阵列设置在地表，用于采集震源产生的震动信号并发送给信号预处理模块；信号预处理模块对震动信号中的地噪声进行预处理后发送给定位参数提取模块；定位参数提取模块用于获取震动信号中的初至波到时信息、高精度速度信息、地下浅层速度信息和震源激励信息；震源激励信息包括震源激励频率以及震源子波类型；定位模型建模模块用于形成地下震源能量场模型，包括虚拟震源激励模块、波阵面提取模块和波阵面融合模块；虚拟震源激励模块在高精度速度信息的基础上，将每个传感器节点作为虚拟震源，以其对应的震源子波类型作为震源子波，共用主频作为震源激励频率，利用交错网格有限差分法进行地下波场三维逆时数值模拟，使每个传感器节点形成一组三维波阵面序列；波阵面提取模块根据每个传感器节点对应的初至波到时信息，得到初至波到时时刻所对应的波阵面；波阵面融合模块对所有的传感器节点对应的波阵面图像进行融合，形成地下震源能量场模型；定位模型解算模块对融合后的震源定位模型进行解算，得到多个波阵面的相交点，即为震源位置震源位置。

有益效果：本发明采用地下震动场波场互易原理，通过传感器阵列数据进行波场逆时数值模拟，将地层速度场信息进行网格划分，利用交错网格有限差分法实现了时-空场全波形逆时外推，生成了传感器作为“震源”下，群波阵面三维空间叠加图，建立了含有震源信息的定位模型。本发明充分利用了浅层地质结构中震动波波场传输特性，以地层速度场网格划分方式，建立了更加精细化的定位模型，为地下浅层震源定位提供一种新的震源模型的建模方法。采用能量聚焦的方法，在建立目标函数的基础上，利用QPSO在震源模型中进行全局快速搜索，实现了大范围内快速震源解算。实现了大区域、多尺度、多线程的快速精准扫描定位。

具体实施方式

为使本发明的目的、内容和优点更加清楚，对本发明的具体实施方式作进一步详细描述。

本发明提供一种基于速度信息的地下浅层微震定位系统，包括震动传感器阵列、信号预处理模块、定位参数提取模块、定位模型建模模块、定位模型解算模块；震动传感器阵列等间距设置在地表，用于采集震源产生的震动信号并发送给信号预处理模块；信号预处理模块一方面用于对震动信号中的地噪声进行预处理，另一方面用于剔除传感器阵列无效数据，可以采用马氏距离等数据一致性评价方法实现；然后将处理后的震动信号发送给定位参数提取模块；定位参数提取模块用于获取震动信号中的初至波到时信息、地下浅层速度信息和震源激励信息；定位参数提取模块首先利用长短时窗法(STA/LTA)提取每个传感器中的震源传播至传感器节点的初至波到时信息，然后利用初至波到时信息，采用初至波走时层析法提取地下浅层速度信息，通过抛物插值最短路径射线追踪法(PTISPR)对速度信息进行修正，得到高精度速度信息；震源激励信息包括震源激励频率以及震源子波类型，定位参数提取模块对传感器阵列数据进行FFT变换，采用相关法对传感器阵列的谱特征进行一致性分析，获取共用主频f_max作为震源激励频率；对传感器阵列数据进行同态滤波处理，得出子波类型；

定位模型建模模块用于形成地下震源能量场模型，包括虚拟震源激励模块、波阵面提取模块、波阵面融合模块；虚拟震源激励模块在高精度速度信息的基础上，将每个传感器节点作为虚拟震源，以其对应的震源子波类型作为震源子波，共用主频f_max作为震源激励频率，利用交错网格有限差分法进行地下波场三维逆时数值模拟，使每个传感器节点形成一组三维波阵面序列；波阵面提取模块根据每个传感器节点对应的初至波到时信息，得到初至波到时时刻所对应的波阵面；波阵面融合模块对所有的传感器节点对应的波阵面图像，采用加权平均空域融合方法进行融合，形成震源定位模型，多个波阵面的相交点即为震源位置；

定位模型解算模块对融合后的震源定位模型进行解算，得到震源位置：该模块可以通过网格搜索法对震源定位模型中的聚焦点进行扫描，实现震源的快速定位；也可以采用卷积核等图像特征识别方法，得到交叉点坐标，即震源位置。

进一步地，为了实现精细化定位，还可以利用能量聚焦原理，构建能流目标函数，采用量子粒子群(QPSO)算法对能流目标函数进行解算，不断迭代寻优，输出震源位置，实现精细化定位，具体包括以下步骤：

S1、利用能量聚焦原理，构建能流目标函数。

能流目标函数定位的主要原理是通过计算给定的震源位置处能量的聚焦程度来实现定位。

首先定义能流矢量：

式中，|A|表示传感器节点接收到的某震源事件的能量，是矢量的长度，表示传感器处波传播的方向余弦。对位于(x,y,z)处的震源，在第k个传感器处的能流矢量可表示为T表示传感器节点所记录接收到的震动信号的时间，表示该传感器记录的T时刻的信号振幅，表示传播方向。因此分别表示三分量记录中震动信号的振幅值。构建如下能流目标函数：

式中，表示上述有限差分逆时数值模拟所计算出的第k个传感器位置处的地震波传播的方向余弦，m代表地层速度信息，(x,y,z)代表震源坐标，|·|表示矢量点积，p表示2范数，N代表传感器个数。

S2、采用量子粒子群(QPSO)算法对能流目标函数进行解算，不断迭代寻优，输出震源位置，具体包括以下步骤：

a、设置粒子种群规模为30，维度为3，在S10利用网格搜索法得到震源初始解(X₀,Y₀,Z₀)的基础上，以网格尺寸scale的10倍为范围，随机生成维度为3的30个粒子，形成初始粒子群(X_i,Y_i,Z_i),其中i代表第i代粒子，初始i＝1；

b、将S11.1的能流目标函数作为粒子的适应度值，即使下式取最大值。

c、根据目标函数计算粒子适应度值，计算每一代粒子适应度得到局部最优位置P_ibest，其中i代表第i代粒子，种群所有代中最大适应度对应粒子作为全局最优位置G_best，粒子群每代平均最优位置

d、更新粒子位置(X_i+1,n,Y_i+1,n,Z_i+1,n)

X_i+1,n＝p_i,n,x±β*|M_best-X_i,n|*ln(1/u) u＝rand(0,1) (4)

Y_i+1,n＝p_i,n,y±β*|M_best-Y_i,n|*ln(1/u) u＝rand(0,1) (5)

Z_i+1,n＝p_i,n,z±β*|M_best-Z_i,n|*ln(1/u) u＝rand(0,1) (6)

其中，β是系数创造力，调节它的值能控制算法的收敛速度，通常情况下，β为0.5。

其中，和是介于0和1之间的随机函数。

e、更新每个粒子每一维的局部最优点P_ibest和群体的全局最优点G_best，即当某个粒子的适应度值优于其经历的最优点P_ibest时，则将其作为新的局部最优点；类似的，若某个粒子的适应度值优于所有粒子所经历的全局最优点G_best，则将其作为新的全局最优点；

f、判断是否达到结束条件，即适应度值是否在误差允许范围内或更新代数达到设定的上限，若满足结束条件则停止更新，若不满足到则跳转至步骤c进行下一次迭代计算。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和变形，这些改进和变形也应视为本发明的保护范围。

Claims (10)

1.一种基于速度信息的地下浅层微震定位系统，其特征在于：包括震动传感器阵列、信号预处理模块、定位参数提取模块、定位模型建模模块和定位模型解算模块；

震动传感器阵列设置在地表，用于采集震源产生的震动信号并发送给信号预处理模块；

信号预处理模块对震动信号中的地噪声进行预处理后发送给定位参数提取模块；

定位参数提取模块用于获取震动信号中的初至波到时信息、高精度速度信息、地下浅层速度信息和震源激励信息；震源激励信息包括震源激励频率以及震源子波类型；

定位模型建模模块用于形成地下震源能量场模型，包括虚拟震源激励模块、波阵面提取模块和波阵面融合模块；虚拟震源激励模块在高精度速度信息的基础上，将每个传感器节点作为虚拟震源，以其对应的震源子波类型作为震源子波，共用主频作为震源激励频率，利用交错网格有限差分法进行地下波场三维逆时数值模拟，使每个传感器节点形成一组三维波阵面序列；波阵面提取模块根据每个传感器节点对应的初至波到时信息，得到初至波到时时刻所对应的波阵面；波阵面融合模块对所有的传感器节点对应的波阵面图像进行融合，形成地下震源能量场模型；

定位模型解算模块对融合后的震源定位模型进行解算，得到多个波阵面的相交点，即为震源位置震源位置。

2.根据权利要求1所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：震动传感器阵列等间距设置在地表。

3.根据权利要求1所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：信号预处理模块一方面用于对震动信号中的地噪声进行预处理，另一方面用于剔除传感器阵列无效数据。

4.根据权利要求1所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：定位参数提取模块利用长短时窗法提取每个传感器中的震源传播至传感器节点的初至波到时信息，然后利用初至波到时信息，采用初至波走时层析法提取地下浅层速度信息，通过抛物插值最短路径射线追踪法对速度信息进行修正，得到高精度速度信息。

5.根据权利要求4所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：定位参数提取模块对传感器阵列数据进行FFT变换，采用相关法对传感器阵列的谱特征进行一致性分析，获取共用主频f_max作为震源激励频率；对传感器阵列数据进行同态滤波处理，得出子波类型。

6.根据权利要求5所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：波阵面融合模块采用加权平均空域融合方法对所有的传感器节点对应的波阵面图像进行融合。

7.根据权利要求5所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：定位模型解算模块可以通过网格搜索法对震源定位模型中的聚焦点进行扫描，实现震源的快速定位；也可以采用图像特征识别方法，得到交叉点坐标，即震源位置。

8.根据权利要求1-7任一项所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：得到震源位置后，构建能流目标函数，采用量子粒子群算法，对能流目标函数进行解算，根据给定的震源位置处能量的聚焦程度输出精确震源位置。

9.根据权利要求8所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：能流目标函数构建如下：

首先定义能流矢量：

式中，|A|表示传感器节点接收到的某震源事件的能量，是矢量的长度，表示传感器处波传播的方向余弦；对位于(x,y,z)处的震源，在第k个传感器处的能流矢量可表示为T表示传感器节点所记录接收到的震动信号的时间，表示该传感器记录的T时刻的信号振幅，表示传播方向；因此 分别表示三分量记录中震动信号的振幅值；

构建如下能流目标函数：

式中，表示上述有限差分逆时数值模拟所计算出的第k个传感器位置处的地震波传播的方向余弦，m代表地层速度信息，(x,y,z)代表震源坐标，|·|表示矢量点积，p表示2范数，N代表传感器个数。

10.根据权利要求9所述的一种基于速度信息的地下浅层微震定位系统，其特征在于：采用量子粒子群算法对能流目标函数进行解算，不断迭代寻优，输出震源位置，具体包括以下步骤：

a、设置粒子种群规模和维度，在得到的震源位置(X₀,Y₀,Z₀)的基础上，以网格尺寸的倍数为范围，随机生成设置维度下的粒子，形成初始粒子群(X_i,Y_i,Z_i),其中i代表第i代粒子，初始i＝1；

b、将能流目标函数作为粒子的适应度值，即使下式取最大值：

c、根据目标函数计算粒子适应度值，计算每一代粒子适应度得到局部最优位置P_ibest，其中i代表第i代粒子，种群所有代中最大适应度对应粒子作为全局最优位置G_best，粒子群每代平均最优位置为：

d、更新粒子位置(X_i+1,n,Y_i+1,n,Z_i+1,n)

X_i+1,n＝p_i,n,x±β*|M_best-X_i,n|*ln(1/u) u＝rand(0,1) (4)

Y_i+1,n＝p_i,n,y±β*|M_best-Y_i,n|*ln(1/u) u＝rand(0,1) (5)

Z_i+1,n＝p_i,n,z±β*|M_best-Z_i,n|*ln(1/u) u＝rand(0,1) (6)

其中，β是系数创造力，调节它的值能控制算法的收敛速度：

其中，和是介于0和1之间的随机数；

e、更新每个粒子每一维的局部最优点P_ibest和群体的全局最优点G_best，即当某个粒子的适应度值优于其经历的最优点P_ibest时，则将其作为新的局部最优点；类似的，若某个粒子的适应度值优于所有粒子所经历的全局最优点G_best，则将其作为新的全局最优点；

f、判断是否达到结束条件，即适应度值是否在误差允许范围内或更新代数达到设定的上限，若满足结束条件则停止更新，若不满足到则跳转至步骤c进行下一次迭代计算。