基于微震或声发射及震源机制约束的岩体裂隙表征方法

摘要

本发明提供一种基于微震或声发射及震源机制约束的岩体裂隙表征方法，涉及岩石裂隙表征技术领域。该方法首先读取微震或声发射事件发生的三维定位坐标以及包含微震或声发射事件能量、事件震级特征的点云数据信息；再对微震或声发射点云数据进行SOR降噪滤波，确定有效微震或声发射事件；然后对有效微震或声发射事件进行欧式聚类；并采用B样条曲面拟合算法将同一类有效微震或声发射事件拟合为裂隙面；最后采用矩张量反演优化裂隙面方位，完成岩体裂隙的表征。同时，还还对最终经反演修正方位的裂隙面进行可视化渲染展示。该方法实现了全自动化反演空间岩体内部裂隙面，提高了微震监测作为矿山安全监测的准确性、便捷性、可靠性。

说明书

技术领域

本发明涉及岩石裂隙表征技术领域，尤其涉及一种基于微震或声发射及震源机制约束的岩体裂隙表征方法。

背景技术

矿山发生大规模岩体失稳破裂前，岩体内部已经积累了大量的微破裂。随着积累的微破裂越来越多，内部形成大规模贯通式破裂，即形成贯通的宏观裂隙面。这种微破裂(以下简称“微震”)发生时，会产生具备一定特征的应力波形，经岩石传播至外界。由于传感器技术的发展，该种波形已经可以拾取，并可根据特定的特征进行分类。同时，使用相关原理已经可以利用波形对岩体内部发生破裂的位置进行定位，获取矿山围岩发生岩体破裂的时空分布，并可通过深层次的震源机制反演分析获取破裂机理。

在已经公开的中国专利(CN113050159A，一种煤岩水力压裂裂缝微震定位及扩展机理监测方法)中：首先根据实际情况布置微震传感器，采集压裂过程中的原始波形信号；在此基础上提出可更新长度的事件筛选窗口，利用长短时窗法扫描各通道数据，识别有效微震事件；通过定义波形敏感值筛选出煤岩体水力压裂诱发的有效微震事件，然后利用局部AIC法获得微震事件的高精度到时，并计算出破裂源的位置，同时根据P波幅值及初动信息反演震源机制；最后将微震定位和震源机制信息结合起来，可以揭示水力压裂裂缝的空间分布形态和破裂机制，确定煤层水力压裂影响范围，这对于压裂过程中参数实时动态调整、压裂半径确定、储层评价具有指导性意义。

该专利在原理上介绍并解释了基于提取微震波形信号反演岩体内部产生的裂缝的机理，同时对利用微震信号进岩体破裂位置进行定位的运算效率做了一定的改善。

中国专利(CN214544804U，煤矿井下微震拾震传感器定位系统)中公开了一种煤矿井下微震拾震传感器定位系统，包括：井下环网(10)，具有地面监控端(11)和与之相耦接的定位基站(13)；拾震传感器(20)，所述拾震传感器(20)适于与所述定位基站(13)之间通过无线通讯相耦接。该定位系统能够方便、及时、有效地测定拾震传感器的坐标，且节省了重新建立井下无线信号覆盖所需要的成本，实用性比较强；相较于传统的人员井下现场测量，节省了时间、提高了拾震传感器坐标测定的效率，同时保证测量的精度；基于煤矿现有的人员定位系统进行拾震传感器坐标测定，保证了井下无线信号的全覆盖。

该专利在实际工程应用上对使用微震进行矿山安全监测进行了诠释，同时改进了实际应用中发生传感器信号不稳定的问题。

以上两个专利分别从原理和应用的角度，对使用微震监测矿山围岩内部发生岩石破裂损伤的应用价值即可实施性进行了全面的介绍。矿山安全不仅关乎生产经营效益，更关系到矿工的切身利益，使用微震监测系统反演并可视化构造矿山围岩内部发生的破裂，对预知矿山可能发生的岩体失稳、边坡滑塌、巷道坍塌等危及矿山生产安全的事件，都具有重要意义。

专利“一种煤岩水力压裂裂缝微震定位及扩展机理监测方法”(CN113050159A)，及专利“煤矿井下微震拾震传感器定位系统”(CN214544804U)对微震监测的原理和实现分别进行了分析与设计，但是其两项技术最终生成的仅为包含微震特征的三维空间点云，并没有从点云处理的角度进行进一步的、根据点云信息可视化的、岩体破裂面全自动的、专业生成与分析。目前矿山使用的多数微震监测系统，都对微震事件发生的位置进行了各种方式的描述，或可以导出微震事件位置，或可以对微震事件进行可视化三维展示，但是仅局限于对微震事件本身进行展示，并未有系统可对微震事件反应的岩体破裂面进行三维描述。由于各微震事件最终反映出来包含三维空间坐标信息，即最终生成了各种包含微震特征信息点的集合，但是现有技术并没有利用该集合进行全方位自动化的岩体裂隙面的分析反演与展示。

现有技术对空间点的集合(以下简称点云)进行处理有较多方式。在《点云库PCL从入门到精通》(机械工业出版社，郭浩，9787111615521)第378页，使用点云处理技术进行三维空间地形点云进行处理并进行可视化三维展示，该技术已经对经雷达测量产生的点云进行了全方位的处理，并对点云数据进行了专业的加工。

虽然在《点云库PCL从入门到精通》(机械工业出版社，郭浩，9787111615521)中对点云处理技术进行了详细阐释，但是并没有结合矿山围岩安全监测的微震监测方法进行空间三维的裂隙面可视化处理展示，更没有从点云生成的本质入手进行裂隙面的约束优化。这就导致缺乏矿业相关的点云处理技术，使微震或声发射定位坐标点的处理缺乏专业点云处理技术支撑。

发明内容

本发明要解决的技术问题是针对上述现有技术的不足，提供一种基于微震或声发射及震源机制约束的岩体裂隙表征方法，自动使用包含微震或声发射特征的点云数据结合震源机制分析反演矿山围岩内部岩体裂隙面，并全方位三维展示分析裂隙面的空间分布。

为解决上述技术问题，本发明所采取的技术方案是：基于微震或声发射及震源机制约束的岩体裂隙表征方法，读取微震或声发射点云数据；

对微震或声发射点云数据进行降噪滤波，确定有效微震或声发射事件；

对有效微震或声发射事件进行聚类；

将同一类有效微震或声发射事件拟合为裂隙面；

反演优化裂隙面方位，完成岩体裂隙的表征。

所述方法具体包括以下步骤：

步骤1：初始化微震或声发射事件分布处理空间，建立微震或声发射点云数据读取路径；

步骤2：读取微震或声发射点云数据；读取微震或声发射事件发生的三维定位坐标以及包含微震或声发射事件能量、事件震级特征的点云数据信息；

步骤3：对微震或声发射点云数据进行SOR降噪滤波，确定有效微震或声发射事件；

步骤4：对有效微震或声发射事件进行欧式聚类；

步骤5：采用B样条曲面拟合算法将同一类有效微震或声发射事件拟合为裂隙面；

步骤6：采用矩张量反演优化裂隙面方位，完成岩体裂隙的表征。

所述步骤3的具体方法为：

步骤3.1：初始化SOR滤波降噪算法，读取微震或声发射事件的三维空间信息，并生成微震或声发射事件数据集；

步骤3.2：读取第一个微震或声发射事件坐标，统计该微震或声发射事件距离其他所有微震或声发射事件的总距离；利用该微震或声发射事件距离其他全部微震或声发射事件的总距离，同微震或声发射事件数量相除，最终得出该微震或声发射事件距离其他全部微震或声发射事件的平均距离；遍历统计全部微震或声发射事件距离其他微震或声发射事件的平均距离；

步骤3.3：统计各微震或声发射事件的平均距离信息，生成平均距离的统计数据；

步骤3.4：根据平均距离的统计信息生成全部平均距离数据的平均值μ及方差σ；

步骤3.5：将平均距离d满足下面公式的微震或声发射事件确定为有效微震或声发射事件，其余为未产生成片微震或声发射事件并贯通为裂隙面的无效微震或声发射事件；

|d-μ|≤σ。

所述步骤4的具体方法为：

步骤4.1：读取经过步骤3确定的有效微震或声发射事件的三维坐标信息及微震特征；

步骤4.2：将全部有效微震或声发射事件按照一定的空间矢量方向生成索引代码，同一个方向梯度上对应的全部有效微震或声发射事件按照距离矢量方向起点由远到近依次生成事件索引；

步骤4.3：将全部生成事件索引的有效微震或声发射事件的微震特征生成微震或声发射事件特征权重，并将特征权重与索引相对应；

步骤4.4：统计并生成有效微震或声发射事件的特征权重由大至小的索引值序列，并提取特征权重最大的微震或声发射事件索引，置为初始微震或声发射事件；

步骤4.5：读取初始微震或声发射事件索引对应的坐标；

步骤4.6：判断该初始微震或声发射事件是否经过读取操作，若经过读取操作，则读取剩余有效微震或声发射事件中特征权重最大的微震或声发射事件索引，并转入步骤4.5；若未经过搜索操作则进行下一步；

步骤4.7：判断当前微震或声发射事件设定距离内是否存在其他微震或声发射事件，若存在则读取该微震或声发射事件索引，并进入下一步；若不存在则执行步骤4.10；

步骤4.8：根据微震或声发射事件索引读取该微震或声发射事件的坐标，并将该微震或声发射事件与初始微震或声发射事件置为同一类，并为这类微震或声发射事件设置类标签；

步骤4.9：判断是否有剩余未设置标签的微震或声发射事件，若有，则提取其中特征权重最大的微震或声发射事件的索引，将该微震或声发射事件置为初始微震或声发射事件，并进入步骤4.5；若没有，则将已经归类的全部微震或声发射事件按照各自的类标签进行存储。

所述步骤5的具体方法为：

步骤5.1：将已经分类好的有效微震或声发射事件，读取出同一个类中的全部有效微震或声发射事件，微震或声发射事件包括微震定位及微震特征的全部信息；

步骤5.2：根据各同类有效微震或声发射事件的特征权重及三维空间位置信息，使用主成分分析法，将该类中特征权重最大的两个微震或声发射事件作为主成分分析的两个初始点，利用两个初始点结合同类的全部有效微震或声发射事件张成最终裂隙面两基本向量即最终裂隙面的基本方向；

步骤5.3：利用已经张成的裂隙面基本方向，采用B样条曲线拟合方法对同一类有效微震或声发射事件进行面边界线条拟合，形成闭环曲线作为最终裂隙面的外轮廓线；

步骤5.4：对已经形成的裂隙面的外轮廓线进行曲面的三角化，对最终的裂隙面内部进行三角网格划分；

步骤5.5：读取全部同一类的有效微震或声发射事件，结合已经建立的三角网格对最终要生成的矿山岩体裂隙曲面进行最终拟合，形成该同一类有效微震或声发射事件的最终裂隙面；

步骤5.6：判断是否还有剩余的经过步骤4归纳的有效微震或声发射事件类没有进行最终拟合裂隙面，若有则读取该类全部有效微震或声发射事件，进入步骤5.2；若没有，则整理各类有效微震或声发射事件生成的裂隙面进行存储。

所述步骤6的具体方法为：

步骤6.1：读取已经生成裂隙面的微震或声发射事件类中全部有效微震或声发射事件，读取有效微震或声发射事件对应的波形信号；

步骤6.2：对有效微震或声发射事件的波形信号进行预处理，调和波形振幅及频率，匹配微震或声发射事件的微震或声发射事件特征；

步骤6.3：对有效微震或声发射事件的震源进行矩张量反演；

步骤6.4：使用主应力约束判断得出矩张量反演的裂隙三维方位信息；

步骤6.5：根据裂隙三维方位信息，调用微震或声发射事件特征权重对已经反演的裂隙方位进行二次加权修正；

步骤6.6：根据经过修正的裂隙方位，拟合并调整步骤5得到的有效微震或声发射事件生成的裂隙面；

步骤6.7：判断是否还存在未经矩张量反演修正方位的裂隙面，若有则读取该裂隙面对应的有效微震或声发射事件类信息，并进入步骤6.2；若没有，则将全部经过矩张量反演修正方位的裂隙面数据导出。

所述方法还创建中转缓存空间，用于中转留存反演生成裂隙面过程产生的点云降噪、点云聚类、点云分割、点云拟合结果以及矩张量反演修正的中间结果的运行数据。

所述方法还对最终经反演修正方位的裂隙面进行可视化渲染，并将渲染结果输出到可视化图形展示界面。

采用上述技术方案所产生的有益效果在于：本发明提供的基于微震或声发射及震源机制约束的岩体裂隙表征方法，使用了SOR降噪滤波算法对离群的三维空间微震或声发射事件进行了过滤，实现了对微震或声发射事件三维空间分布点的离群点去除；对三维空间微震或声发射事件进行了带有能量权重的欧式聚类，使经过降噪的有效微震或声发射事件划分出各自生成裂隙面的范围，对各有效微震或声发射事件进行了有效的空间归类；对三维空间微震或声发射事件进行了带有能量权重的B样条曲面拟合，使经过聚类的有效微震或声发射事件生成了空间裂隙面；针对已经拟合得出的空间裂隙面进行了矩张量反演，根据裂隙形成机理修正裂隙方向。

本发明方法融合了三维空间点云处理技术及地球物理中矩张量反演理论，实现了由包含矿山微震三维空间坐标、波形、能量权重等特征的微震或声发射事件信息反演出三维空间矿山围岩内部裂隙面，并最终实现了三维可视化综合综合展示，为矿山围岩稳定性分析，加强矿山安全，提高矿企生产安全性等方面提供可靠依据。

与现有技术相比，本发明实现了全自动化反演空间岩体内部裂隙面，而不仅仅是提供了微震或声发射事件的点状空间展示，从而进一步提高了微震和声发射监测作为矿山安全监测的准确性、便捷性、可靠性。

附图说明

图1为本发明实施例提供的基于微震或声发射及震源机制约束的岩体裂隙表征方法的流程图；

图2为本发明实施例提供的SOR降噪的流程图；

图3为本发明实施例提供的欧式聚类的流程图；

图4为本发明实施例提供的B样条曲面拟合裂隙的流程图；

图5为本发明实施例提供的矩张量反演修正裂隙面方位的流程图；

图6为本发明实施例提供的未经处理的微震或声发射事件点云图；

图7为本发明实施例提供的微震或声发射事件反演裂隙面局部图。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。

本实施例中，基于微震或声发射及震源机制约束的岩体裂隙表征方法，如图1所示，包括以下步骤：

步骤1：初始化微震或声发射事件分布处理空间，建立微震或声发射点云数据读取路径，调整点云数据反演裂隙面的可视化展示界面为相应合适大小；

步骤2：读取微震或声发射点云数据；读取微震或声发射事件发生的三维定位坐标以及包含微震或声发射事件能量、事件震级特征的点云数据信息；同时创建中转缓存，用于中转留存后续反演生成裂隙面过程产生的点云降噪、点云聚类、点云分割、点云拟合结果、矩张量反演修正的中间结果的运行数据空间；

步骤3：对微震或声发射点云数据进行SOR降噪滤波，确定有效微震或声发射事件；在处理微震或声发射点云时，由于实际矿山围岩在矿山开采过程中产生了各种由于开采扰动的内部裂隙，某些微震或声发射事件实际并未引起裂隙破裂贯通最终生成裂隙面。故此，需要根据各微震或声发射事件定位坐标的离散性及统计学特征，分析并自动将这种扰动事件排除在后续分析的微震或声发射事件点云的范围，同时可以提高微震反演裂隙面的精确度。主要流程如图2所示，具体包括以下步骤：

步骤3.1：初始化SOR滤波降噪算法，读取已经存入中转缓存的微震或声发射事件三维空间信息，并生成微震或声发射事件数据集；

步骤3.2：读取第一个微震或声发射事件坐标，统计该微震或声发射事件距离其他所有微震或声发射事件的总距离；利用该微震或声发射事件距离其他全部微震或声发射事件的总距离，同微震或声发射事件数量相除，最终得出该微震或声发射事件距离其他全部微震或声发射事件的平均距离；遍历统计输入的全部微震或声发射事件距离其他微震或声发射事件的平均距离，将微震或声发射事件定位坐标点的标号及其距离其他微震或声发射事件的平均距离的信息存入中转缓存；

步骤3.3：统计中转缓存内部各微震或声发射事件的平均距离信息，生成平均距离的统计数据；

步骤3.4：根据平均距离的统计信息生成全部平均距离数据的平均值μ及方差σ；

步骤3.5：将平均距离d满足下面公式的微震或声发射事件确定为有效微震或声发射事件，其余为未产生成片微震或声发射事件并贯通为裂隙面的无效微震或声发射事件；并将有效微震或声发射事件及无效微震或声发射事件的定位信息及微震特征信息分别存入中转缓存：

|d-μ|≤σ

步骤4：对有效微震或声发射事件进行欧式聚类；在实际的矿山围岩中，并不是全部发生的微震或声发射事件仅贯通为一个裂隙面，因而需要对距离相近的微震或声发射事件归为同一类。在定义为一类的微震或声发射事件时，需要加入微震特征，例如微震或声发射事件能量、事件能级、事件定位误差等特征。因为微震或声发射事件积累并贯通为一个裂隙面的过程中存在一个能量或者震级最大的微震或声发射事件，其余微震或声发射事件皆围绕该微震或声发射事件环绕发展，最终形成贯通。因而，在实际根据微震或声发射事件生成裂隙面的反演过程可以取能量权重最大者为该裂隙面的初始信息点。本发明在实际发明阶段亦考虑了其他微震特征因素，最终可以根据权重公式将微震或声发射事件能量、事件能级、事件定位误差这些特征一并加入，并不仅考虑能量影响。欧式聚类过程的主要流程如图3所示，具体为：

步骤4.1：初始化中转缓存的缓存空间，读取经过步骤3确定的有效微震或声发射事件的三维坐标信息及微震特征，并放入中转缓存；

步骤4.2：将全部有效微震或声发射事件按照一定的空间矢量方向生成索引代码，同一个方向梯度上对应的全部有效微震或声发射事件按照距离矢量方向起点由远到近依次生成事件索引；

步骤4.3：将全部生成事件索引的有效微震或声发射事件的微震特征按照前述的方式生成微震或声发射事件特征权重，并将特征权重与索引相对应；至此，各有效微震或声发射事件包括索引值、定位坐标、特征权重三方面信息；

步骤4.4：统计并生成有效微震或声发射事件的特征权重由大至小的索引值序列，将经过排序后的有效微震或声发射事件序列存入中转缓存，清空前述生成的无序有效微震或声发射事件，并提取特征权重最大的微震或声发射事件索引，置为初始微震或声发射事件；

步骤4.5：读取初始微震或声发射事件索引对应的坐标；

步骤4.6：判断该初始微震或声发射事件是否经过读取操作，若经过读取操作，则读取剩余有效微震或声发射事件中特征权重最大的微震或声发射事件索引，并转入步骤4.5；若未经过读取操作则进行下一步；

步骤4.7：判断当前微震或声发射事件设定距离内是否存在其他微震或声发射事件，若存在则读取该微震或声发射事件索引，并进入下一步；若不存在则执行步骤4.10；

步骤4.8：根据微震或声发射事件索引读取该微震或声发射事件的坐标，并将该微震或声发射事件与初始微震或声发射事件置为同一类，并为这类微震或声发射事件设置类标签；

步骤4.9：判断是否有剩余未设置标签的微震或声发射事件，若有，则提取其中特征权重最大的微震或声发射事件的索引，将该微震或声发射事件置为初始微震或声发射事件，并进入步骤4.5；若没有，则将已经归类的全部微震或声发射事件按照各自的类标签分别存入中转缓存；由于本步骤产生的初始微震或声发射事件可能不止一个，因而本步骤产生初始微震或声发射事件后欧式聚类过程进入多初始微震或声发射事件并行状态；

步骤5：采用B样条曲面拟合算法将同一类有效微震或声发射事件拟合为裂隙面；经过前述步骤，全部读取的微震或声发射事件已经经过无效事件排除，生成同一裂隙面的有效微震或声发射事件归为同一类。下面进行将同一类有效微震或声发射事件拟合为裂隙面。B样条曲面拟合算法过程的主要步骤流程图如图4所示，其具体步骤为：

步骤5.1：将已经分类好的有效微震或声发射事件，从中转缓存中读取出同一个类中的全部有效微震或声发射事件，微震或声发射事件包括微震定位及微震特征的全部信息；

步骤5.2：根据各同类有效微震或声发射事件的特征权重及三维空间位置信息，使用主成分分析法(PCA)，将该类中特征权重最大的两个微震或声发射事件作为主成分分析的两个初始点，利用两个初始点结合同类的全部有效微震或声发射事件张成最终裂隙面的两基本向量即最终裂隙面的基本方向；

步骤5.3：利用已经张成的裂隙面基本方向，采用B样条曲线拟合方法对同一类有效微震或声发射事件进行面边界线条拟合，形成闭环曲线作为最终裂隙面的外轮廓线；

步骤5.4：对已经形成的裂隙面的外轮廓线进行曲面的三角化，对最终的裂隙面内部进行三角网格划分；

步骤5.5：读取全部同一类的有效微震或声发射事件，结合已经建立的三角网格对最终要生成的矿山岩体裂隙曲面进行最终拟合，形成该同一类有效微震或声发射事件的最终裂隙面，并将生成的裂隙面存入中转缓存；

步骤5.6：判断是否还有剩余的经过步骤4归纳的有效微震或声发射事件类没有进行最终拟合裂隙面，若有则读取该类全部有效微震或声发射事件，进入步骤5.2；若没有，则整理各类有效微震或声发射事件生成的裂隙面，存入中间缓存；

步骤6：采用矩张量反演优化裂隙面方位；由于经过上述步骤生成的微震裂隙面仅使用了微震或声发射事件的特征权重及点云处理算法，并未对微震或声发射事件的波形进行处理及分析，即未考虑微震或声发射事件的形成机制。通过对微震或声发射事件的波形进行矩张量反演，可得出该事件的破裂类型以及该处震源的运动特性，继而可实现对上述得出裂隙面的方位优化。使用矩张量反演对已经生成的裂隙面进行进一步优化为本发明的创新点之一。矩张量反演优化裂隙面方位算法的主要流程如图5所示，其具体步骤为：

步骤6.1：读取已经生成裂隙面的微震或声发射事件类中全部有效微震或声发射事件，读取有效微震或声发射事件对应的波形信号；

步骤6.2：对有效微震或声发射事件波形信号进行预处理，调和波形振幅及频率，匹配微震或声发射事件的微震或声发射事件特征；

步骤6.3：对有效微震或声发射事件的震源进行矩张量反演；

步骤6.4：使用主应力约束判断得出矩张量反演的裂隙三维方位信息，将并存入中转缓存；

步骤6.5：根据裂隙三维方位信息，调用微震或声发射事件特征权重对已经反演的裂隙方位进行二次加权修正，提高反演裂隙面的位置准确性；

步骤6.6：根据经过修正的裂隙方位，拟合并调整中转缓存中步骤5得到的有效微震或声发射事件生成的裂隙面；

步骤6.7：判断中转缓存中是否还存在未经矩张量反演修正方位的裂隙面，若有则读取该裂隙面对应的有效微震或声发射事件类信息，并进入步骤6.2；若没有，则将全部经过矩张量反演修正方位的裂隙面数据导出，并存入内存；

步骤7：读取内存中存入的最终经反演修正方位的裂隙面进行可视化渲染，并将渲染结果输出到可视化图形展示界面。

本实施例，通过如图6所示的微震或声发射事件点云图及其对应的事件信息成功反演并精确修正了其矿山围岩内部如图7所示的裂隙面，成功地达到了实际为矿山提供生产安全决策依据的水平，为矿山安全生产真正地提供了精确的重点关注区域。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明权利要求所限定的范围。