基于深度学习融合驱动的微震稀疏台网定位方法

1.本发明涉及一种微震定位方法，具体是一种基于深度学习融合驱动的微震稀疏台网定位方法，属于井下微震定位技术领域。


背景技术：

2.随着我国煤炭生产重心向深部转移，开采条件日益复杂，冲击地压灾害日趋严重，已成为煤矿开采的重大安全隐患之一。为最大限度减小冲击地压造成的危害，众多监测预警技术，如钻屑法、地音、电磁辐射以及应力在线等广泛应用于冲击地压的监测预警，其中基于微震系统的定位与监测分析技术是最基础与使用最广泛的冲击地压监测预警技术。微震事件的精准定位是开展后续监测预警分析的基础，但是受煤矿复杂的采掘环境与台网布设不足等影响，目前的微震定位方法存在一定的局限性。
3.微震事件的定位关键在于提高对震源中心的定位精度。如今在台站数量充足的情况下，可以较为可靠的实现微震事件的定位。其中主流的微震定位方法包括纵波首次进入时间法和相对定位法。其中，纵波首次进入时间法中的定位方程具有四个未知参数，至少需要四个台站的数据才可以定位，少于四个台站的情况无法定位，视为台站缺失；而传统相对定位法，例如双差定位法，对监测数据的完整程度要求较高，在台站缺乏情况下会极大限制其定位精度，并且在单台站情况下定位精度较差。因此在井下台网布设不理想的情况下，传统方法往往会出现信息缺失导致定位性能严重下降。


技术实现要素：

4.本发明的目的是提供一种基于深度学习融合驱动的微震稀疏台网定位方法，该方法能够对少于四台站的井下微震位置进行准确定位。
5.为了实现上述目的，本发明提供一种基于深度学习融合驱动的微震稀疏台网定位方法，包括以下步骤：
6.步骤1：将台站安装在井下开采巷道中，并收集长期微震波形数据，将目标检测区域发生过的微震的波形数据进行记录；
7.步骤2：在目标检测区域中，按微震发生的位置给微震波形数据标记标签，并将发生过的微震分为不同的地理集群；使用k-means算法将微震聚类，获得u个聚类中心，根据u个聚类中心定义聚类区域，聚类中心则为聚类区域的质心，划分u个多边形，划分多边形区域的方法学名叫voronoi多边形，即划分的多边形内的每个点到该多边形质心的距离都小于到其他多边形质心的距离，地图上的所有微震发生位置都以点表示，分配给最接近的聚类区域，由此将微震的波形数据分类为包含n个类：0类对应于没有任何微震的井下噪声，1类到n类对应于来自相应地理区域的微震簇；
8.步骤3：利用生成对抗网络的方法生成额外新波形窗口的方法，将步骤2中在1类到n类中的数据明显少于其他类数据的微小数据集生成额外数据窗口以扩大数据集；
9.步骤4：将步骤2的数据和步骤3生成的数据分为两个独立测试集和训练集，其中测
试集和训练集的比例为2：8；
10.步骤5：采用微震波形数据作为输入，将频率为100hz的微震波形数据切割成10s一个的窗口，则10s的数据窗口可生成1000个波形采样样本，将每1hz波形数据窗口的波形图进行输入，设计微震区域检测模型，所述微震区域检测模型包括：波形数据输入、数据规则化、卷积运算、激活、池化、全连接与softmax分类操作，微震区域检测模型在训练后还进行初始化权重以及偏置，通过反向传播算法得到优化模型，运行后微震区域检测模型输出一个n维向量，该n维向量的第n位为对应改微震属于位置u的概率，微震区域检测模型计算概率分布数据的方法为：
[0011][0012]
式中：w是所有权重的集合；
[0013]
b是所有偏差的集合；
[0014]
zc为每一类的计算分数；
[0015]
pc为每类最后的计算概率；
[0016]
c为类的标签；
[0017]
z0为波形数据；
[0018]
zk为求和函数中每一类的计算分数；
[0019]
n为类；
[0020]
k为求和函数的内部参数；
[0021]
步骤6：在单台站情况下，可对同一水平面微震事件进行区域定位，双台站和三台站情况下可进行精准空间定位，首先由台站d来计算微震a和微震b之间的间隔距离|r
ab
|，利用台站d采集井下开采巷道附近产生的微震信号，获取p波和s波到达台站d的微震速度v
p
和vs，将v
p
和vs进行计算得到速度系数kv；
[0022]kv
＝v
pvs
/(v
p-vs)
ꢀꢀꢀ
(1)
[0023]
步骤7：获取p波和s波到达台站d的时间，根据p波和s波传播到台站d的时间不同来确定微震位置的方法，当p波和s波传播到台站d时，就会触发台站d内的检波器记录下接收到p波和s波到达检波器的时刻，根据到时不同理论，可得出微震分别到达台站d的距离|r
ed
|：
[0024][0025]
步骤8：如果|r
ab
|远小于这些微震和台站d之间的距离，那么微震间隔距离用以下等式近似，其中是微震a、b之间距离的近似值，|r
ad
|和|r
bd
|分别是微震a和微震b到单台站d的距离：
[0026][0027][0028]
式中：是s波到达微震a的时间；
[0029]
是p波到达微震a的时间；
[0030]
是s波到达微震b的时间；
[0031]
是p波到达微震b的时间；
[0032]
步骤9：若拥有两个台站d1和d2且同样满足远离两个震源这一特点的话，分别对两个台站d1和d2的数据进行计算，然后使用勾股定理以获得更好的近似：
[0033][0034]
步骤10：若拥有三个台站d1、d2和d3，则按三个台站d1、d2和d3距离震源位置比较，选取位置更远的两个台站执行步骤9；
[0035]
步骤11：获取步骤5中得出的区域中l个已知的非共面微震的震源位置，通过台站计算这些已知震源的位置和待测震源之间的间距并建立距离方程组来得出目标位置。首先根据区域内已有微震簇位置构建微震坐标系，在一个由l 1次微震组成的微震簇，其中x1,x2,
…
x
l 1
为每次微震发生时的震源中心坐标，在三维坐标系r3中定义，获取每对微震i和j之间的欧氏距离，|r
ij
|是微震i和j之间的距离，|xi|和|xj|分别是微震i和j与参考系的原点的距离(i,j＝1,2,
…
,l 1)，可得：
[0036][0037]
步骤12：重构坐标系，利用微震簇的内部结构不随平移和旋转而改变的特性(即可以平移和旋转它，但其中的内部结构不改变)，因此将第l 1次微震发生位置设计为原点，即x
l 1
＝(0,0,0),此时|r
il 1
|＝|xi|，|r
jl 1
|＝|xj|，因此将步骤11转化为：
[0038][0039]
步骤13:通过公式(7)得出x，其中g＝v(1:k,1:3)，∑＝s(1:3,1:3)，l:m为沿着一个轴从第l个位置到第m个定义一个l
×
3的矩阵x用来放置步骤11中设计的原始坐标系的微震震源坐标，再定义一个l
×
3的矩阵xa用来放置步骤12中新坐标系的微震震源坐标，其中x
l 1
＝(0,0,0)，矩阵xa、x的每一行为震源的xyz三维坐标，首先使用平移向量将xa转换为x，将xa中每个坐标都进行偏移转换，其中分别是原始参考系和新参考系中微震簇的几何中心的位置，
[0040]
x＝x
i i＝1,2,
…
,l
ꢀꢀꢀ
(8)
[0041]
步骤14：定义一个l
×
l大小的距离矩阵r，存放步骤11中每对微震之间的距离，第i行第j列即为微震i和j之间距离，下面为每一位的计算公式；
[0042][0043]
步骤15：矩阵r可由公式(7)和公式(9)得
[0044]
r＝xx
t r
l
×
l
ꢀꢀ
(10)
[0045]
则现在可得r的奇异值分解，其中v是一个正交矩阵，s(l
×
l)是r的奇异值对角矩阵
[0046][0047]
步骤16：通过此公式得出x中的参数，其中g＝v(1:l,1:3)，∑＝s(1:3,1:3)，h:m表示沿着一个轴从第h个位置到第m个位置提取矩阵的元素，
[0048]
x＝g∑
1/2
ꢀꢀ
(12)
[0049]
与现有技术相比，本发明由传统的相对定位方法配合机器学习算法改进而来，提出了一种基于深度学习融合驱动的微震稀疏台网定位方法，可以在台站数量不足的情况下也能精准定位微震位置；首先，研究基于聚类算法的微震定位数据标注方法，构建定位波形数据集，在此阶段利用收集到的井下台站记录的长期微震波形数据，通过聚类算法计算波形数据，将井下划分成多个可能的微震区域，从而构建微震定位数据集，其次，基于已构建的微震定位数据集，研究台网稀疏情况下基于深度学习的区域定位模型，在此阶段设计机器学习模型进行预测，该模型通过对波形数据的处理，将输出待预测微震发生位置概率最大的区域。最后，在区域定位的基础之上，本发明进一步拟设计矿山微震事件精准定位方法，利用已经监测出的区域内所有曾发生过的微震位置和其震波数据，通过震波的速度和到时计算获得与台站的距离，并由此来计算每一对微震之间的位置，最后经由建立坐标系构建和矩阵计算的方法获得目标微震的位置；本发明能够在井下台站数量少于四个的情况下进行准确定位。
附图说明
[0050]
图1为本发明微震定位系统中微震区域确定的流程图；
[0051]
图2为本发明微震定位系统中微震区域微震精确定位的流程图；
[0052]
图3为本发明微震区域划分和井下台站设计示意图；
[0053]
图4为本发明单台站测间距示意图；
[0054]
图5为本发明双台站测间距示意图；
[0055]
图6为本发明三台站测间距示意图。
具体实施方式
[0056]
下面结合附图对本发明做进一步说明。
[0057]
如图1-图3所示，一种基于深度学习融合驱动的微震稀疏台网定位方法，包括以下步骤：
[0058]
步骤1：将台站安装在井下开采巷道中，并收集长期微震的波形数据，将目标检测区域发生过的微震的波形数据进行记录；
[0059]
步骤2：在目标检测区域中，按微震发生的位置给微震波形数据标记标签，并将发生过的微震分为不同的地理集群；使用k-means算法将微震聚类，获得u个聚类中心，根据u个聚类中心定义聚类区域，聚类中心则为聚类区域的质心，划分u个多边形，划分多边形区域的方法学名叫voronoi多边形，即划分的多边形内的每个点到该多边形质心的距离都小于到其他多边形质心的距离，地图上的所有微震发生位置都以点表示，分配给最接近的聚类区域，由此将微震的波形数据分类为包含n个类：0类对应于没有任何微震的井下噪声，1类到n类对应于来自相应地理区域的微震簇；
[0060]
步骤3：利用生成对抗网络的方法生成额外新波形窗口的方法，将步骤2中在1类到n类中的数据明显少于其他类数据的微小数据集生成额外数据窗口以扩大数据集；
[0061]
步骤4：将步骤2的数据和步骤3生成的数据分为两个独立测试集和训练集，其中测试集和训练集的比例为2：8；
[0062]
步骤5：采用微震波形数据作为输入，将频率为100hz的微震波形数据切割成10s一个的窗口，则10s的数据窗口可生成1000个波形采样样本，将每1hz波形数据窗口的波形图进行输入，设计微震区域检测模型，所述微震区域检测模型包括：波形数据输入、数据规则化、卷积运算、激活、池化、全连接与softmax分类操作，微震区域检测模型在训练后还进行初始化权重以及偏置，通过反向传播算法得到优化模型，运行后微震区域检测模型输出一个n维向量，该n维向量的第n位为对应改微震属于位置u的概率，微震区域检测模型计算概率分布数据的方法为：
[0063][0064]
式中：w是所有权重的集合；
[0065]
b是所有偏差的集合；
[0066]
zc为每一类的计算分数；
[0067]
pc为每类最后的计算概率；
[0068]
c为类的标签；
[0069]
z0为波形数据；
[0070]
zk为求和函数中每一类的计算分数；
[0071]
n为类；
[0072]
k为求和函数的内部参数；
[0073]
步骤6：在单台站情况下，可对同一水平面微震事件进行区域定位，双台站和三台站情况下可进行精准空间定位，首先由台站d来计算微震a和微震b之间的间隔距离|r
ab
|，利用台站d采集井下开采巷道附近产生的微震信号，获取p波和s波到达台站d的微震速度v
p
和vs，将v
p
和vs进行计算得到速度系数kv：
[0074]kv
＝v
pvs
/(v
p-vs)
ꢀꢀ
(1)
[0075]
步骤7：获取p波和s波到达台站d的时间，根据p波和s波传播到台站d的时间不同来确定微震位置的方法，当p波和s波传播到台站d时，就会触发台站d内的检波器记录下接收到p波和s波到达检波器的时刻，根据到时不同理论，可得出微震分别到达台站d的距离|r
ed
|：
[0076][0077]
步骤8：如果|r
ab
|远小于这些微震和台站d之间的距离，那么微震间隔距离用以下等式近似，其中是微震a、b之间距离的近似值，|r
ad
|和|r
bd
|分别是微震a和微震b到单台站d的距离，如图4所示：
[0078]
[0079][0080]
式中：是s波到达微震a的时间；
[0081]
是p波到达微震a的时间；
[0082]
是s波到达微震b的时间；
[0083]
是p波到达微震b的时间；
[0084]
步骤9：如图5所示，若拥有两个台站d1和d2且同样满足远离两个震源这一特点的话，分别对两个台站d1和d2的数据进行计算，然后使用勾股定理以获得更好的近似：
[0085][0086]
步骤10：如图6所示，若拥有三个台站d1、d2和d3，则按台站距离震源位置比较，选取位置更远的两个台站执行步骤9；
[0087]
步骤11：获取步骤5中得出的区域中l个已知的非共面微震的震源位置，通过台站计算这些已知震源的位置和待测震源之间的间距并建立距离方程组来得出目标位置。首先根据区域内已有微震簇位置构建微震坐标系，在一个由l 1次微震组成的微震簇，其中x1,x2,
…
x
l 1
为每次微震发生时的震源中心坐标，在三维坐标系r3中定义，获取每对微震i和j之间的欧氏距离，|r
ij
|是微震i和j之间的距离，|xi|和|xj|分别是微震i和j与参考系的原点的距离(i,j＝1,2,
…
,l 1)，可得：
[0088][0089]
步骤12：重构坐标系，利用微震簇的内部结构不随平移和旋转而改变的特性(即可以平移和旋转它，但其中的内部结构不改变)，因此将第l 1次微震发生位置设计为原点，即x
l 1
＝(0,0,0),此时|r
il 1
|＝|xi|，|r
jl 1
|＝|xj|，因此将步骤11转化为：
[0090][0091]
步骤13:通过公式(7)得出x，其中g＝v(1:k,1:3)，∑＝s(1:3,1:3)，l:m为沿着一个轴从第l个位置到第m个定义一个l
×
3的矩阵x用来放置步骤11中设计的原始坐标系的微震震源坐标，再定义一个l
×
3的矩阵xa用来放置步骤12中新坐标系的微震震源坐标，其中x
l 1
＝(0,0,0)，矩阵xa、x的每一行为震源的xyz三维坐标，首先使用平移向量将xa转换为x，将xa中每个坐标都进行偏移转换，其中分别是原始参考系和新参考系中微震簇的几何中心的位置，
[0092]
x＝x
i i＝1,2,
…
,l
ꢀꢀꢀ
(8)
[0093]
步骤14：定义一个l
×
l大小的距离矩阵r，存放步骤11中每对微震之间的距离，第i行第j列即为微震i和j之间距离，下面为每一位的计算公式：
[0094][0095]
步骤15：矩阵r可由公式(7)和公式(9)得
[0096]
r＝xx
t r
l
×
l
ꢀꢀ
(10)
[0097]
则现在可得r的奇异值分解，其中v是一个正交矩阵，s(l
×
l)是r的奇异值对角矩阵
[0098][0099]
步骤16：通过此公式得出x中的参数，其中g＝v(1:l,1:3)，∑＝s(1:3,1:3)，h:m表示沿着一个轴从第h个位置到第m个位置提取矩阵的元素，
[0100]
x＝g∑
1/2
ꢀꢀ
(12)。
[0101]
本发明主要使用的是台站，可以输出微震的波形，安装在有微震风险的井下开采巷道,如图1所示(步骤1-5)，首先，研究基于聚类算法的微震定位数据标注方法，构建定位波形数据集。在此阶段利用收集到的井下台站记录的长期微震波形数据，通过聚类算法计算波形数据，将井下划分成多个可能的微震区域，从而构建微震定位数据集。其次，基于已构建的微震定位数据集，研究台网稀疏情况下基于深度学习的区域定位模型，在此阶段设计机器学习模型进行预测，该模型通过对波形数据的处理，将输出待预测微震发生位置概率最大的区域,最后，如图2所示(步骤6-16)，在区域定位的基础之上，本发明进一步拟设计矿山微震事件精准定位方法，利用已经监测出的区域内所有曾发生过的微震位置和其震波数据，通过震波的速度和到时计算获得与台站的距离，并由此来计算每一对微震之间的位置，最后经由建立坐标系构建和矩阵计算的方法获得目标微震的位置。本发明能够在井下台站数量少于四个的情况下进行准确定位。