Research on microseismic event imaging based on waveform clustering analysis
-
摘要: 微地震数据的质量极大地影响震源定位结果的准确性，尤其在地面微地震数据中P波、S波初至不明显或难以识别的低信噪比信号，基于初至拾取的定位方法往往具有较大误差。由于诱发地震的震源机制不同，不同检波器接收到的信号特征也不尽相同，这些因素都会对微地震事件准确定位造成一定的影响。本文提出一种基于聚类分析的微地震数据优选方法，以提高微地震事件振幅叠加定位的成像质量。首先，通过聚类分析计算各数据道之间的欧式距离，通过聚类图去除距离较大和较小的数据道，即去除信号特征不明显或噪声干扰严重数据道; 其次，对各道信号进行互相关计算，提取信号特征相关性强数据道，去除P波特征不明显的数据道; 最后，对筛选数据道进行振幅叠加定位。与原始数据道的定位结果对比表明，微地震聚类分析定位效果更加聚焦，最大震源点位置更加清晰，定位结果得到明显改善。Abstract: The quality of microseismic data greatly affects the accuracy of the source location results, especially in surface microseismic data where the first arrivals of P waves and S waves are not obvious or the low signal-to-noise ratio signals are difficult to identify.Secondly, due to the different focal mechanisms that induce microseismic, the signal characteristics received by different geophones are not the same.All of these factors have a certain impact on the accurate location of microseismic events.This paper proposes a method for selecting microseismic data based on cluster analysis to improve the imaging quality of microseismic event amplitude superposition locating.First, the Euclidean distance between the trace channels is calculated by cluster analysis, and the trace channels with larger and smaller distances are removed through the clustering graph, that is, trace channels with insignificant signal characteristics or severe noise interference are removed.Second, the cross-correlation calculation is performed on the signals of each channel, the trace channels with strong signal characteristic correlation being extracted, and those with unobvious P wave characteristics removed.Finally, the screened trace channels are superimposed and located with amplitude superposition.The comparison of the locating results with the original trace channels shows that the locating effect of the microseismic cluster analysis is more focused, the location of the focal point is clearer, and the locating results are significantly improved.
-
Keywords:
- microseismic /
- cluster analysis /
- waveform classification /
- data denoising /
- events location
-
微地震事件定位面临两大难题，一是发震时刻难以确定，二是微地震事件的信号难以识别[1]。相对于井中采集的信号，地面检波器阵列通常信噪比较低，原因主要是震源到地面检波器距离较远，导致子波衰减严重，其次是地面观测阵列更容易受噪声(施工噪声、环境噪声等)影响[2]。
为了克服这两大难题，微地震震源定位发展了两类方法，一类是震源成像的方法，另一类是基于微地震事件振幅叠加方法[3-4]。振幅叠加或网格搜索方法，是将储层一定区域内划分为大小相同的网格，根据速度模型计算每个网格到所有检波器的走时，得到所有网格点到全部检波器的走时表[5]。根据走时表的时差关系或利用波形的互相关，提取走时差进行偏移叠加数据道波形，从而反演储层区域的成像结果，成像结果的最大值位置即为震源位置。该方法的优点是，无须识别、拾取初至即可进行对地下目标范围和地震记录进行扫描，完成震源的成像结果[6]。但通常叠加的数据不会进行筛选，导致噪声干扰数据参与运算，降低了成像效果。同时，不同的震源机制(ISO源、DC源、CLVD源)会辐射不同的波场特征，地面不同位置检波器观测到的P波和S波特征也不相同[7]。在叠加成像过程中，低信噪比和信号特征差异都会对成像结果造成影响。
本文应用了一种聚类分析方法消除这种影响。聚类分析又称为群分析，它是研究(样品和变量)分类问题的一种多元统计方法，并对相似属性特征的对象进行分组的数据挖掘方法[8]。正如统计学中样本之间有不同的分类方法，信号也可以看作是等间隔、连续采样的样本。通过聚类分析可以定义信号与信号之间的距离，距离越小，信号差异越小; 反之，差异越大。根据不同信号的特征对信号进行聚类计算，筛选出信号特征接近的数据集进行成像，改善成像的质量。相比于传统波形互相关计算，聚类分析对地震道噪声敏感性更低，不易受噪声影响。波形互相关计算需精确提取事件的完整波形才能计算出较为精确的相似性结果，而聚类分析只需包含事件本身就可对其进行分类。
1. 波形聚类方法
波形聚类方法主要包括数据道之间的距离计算方法和聚类算法。
学者提出多种距离计算方法，如欧氏(Euclidean)距离、切比雪夫(Chebychev)距离、兰氏(Canberra)距离、马氏(Mahalanobis)距离、夹角余弦(Cosine)距离以及汉明(Hamming)距离等[9]。每种距离计算适用不同的样本特点，尝试以上方法对微地震数据进行了测试，发现欧式距离能获得不错的聚类结果。欧氏距离计算方法如下：
$$ d_{i j}=\left(\sum\limits_{a=1}^{p}\left|x_{i a}-x_{j a}\right|^{2}\right)^{1 / 2} $$ (1) 式中，dij为数据道Xi与数据道Xj的距离; P为数据道中采样点总数; xia、xja分别为数据道Xi与Xj中第a个采样点位置的振幅值。
定义dij有如下性质：dij≥0，dij=dji且dij≤dik+dkj。欧式距离计算结果与波形数据的量纲有关，因此在计算前需要对各道数据进行归一化处理。通过式(1)可得距离矩阵D(0)：
$$ \boldsymbol{D}_{(0)}=\left[\begin{array}{cccc} d_{11} & d_{12} & \cdots & d_{1 n} \\ d_{21} & d_{22} & \cdots & d_{2 n} \\ \vdots & \vdots & & \vdots \\ d_{n 1} & d_{n 2} & \cdots & d_{n n} \end{array}\right] $$ (2) 式(2)中，d11=d22=…=dnn=0，D是一个对称矩阵，在聚类时只需计算上三角或下三角形部分即可。根据D(0)可对n道数据进行分类，距离近的数据归为一类，距离远的数据归为另一同类。
聚类也有多种方法，如最短距离法、最长距离法以及重心法等。本文采用最短距离法来说明聚类的过程：
(1) 确定距离的计算方法。计算波形数据之间的距离，得到距离矩阵D(0)，定义类与类之间的距离为Dij，开始时每道数据自成一类，显然此时Dij=dij。
(2) 找出D(0)的非对角线最小元素，设为dpq，则将波形数据道X(p)和X(q)合并成一个新的类，记为Gr={X(p)，X(q)}。
(3) 计算所有剩余类Gk与新类Gr的距离：
$$ {{D}_{kr}}={{\min }_{{{X}_{(i)}}\in {{G}_{k}},{{X}_{(j)}}\in {{G}_{r}}}}{{d}_{ij}} \\ =\min \left\{ {{\min }_{{{x}_{(i)}}\in {{G}_{k}},{{x}_{(j)}}\in {{G}_{p}}}}{{d}_{ij}},{{\min }_{{{x}_{(i)}}\in {{G}_{k}},{{X}_{(j)}}\in {{G}_{q}}}}{{d}_{ij}} \right\} \\ =\min \left\{ {{D}_{kp}},{{D}_{kq}} \right\} $$ (3) (4) 将D(0)中p、q行及p、q列删除，然后在D(0)的第一行第一列位置插入Gr，得到一个n-1维距离矩阵，记为D(1)。
(5) 对D(1)重复步骤(2)(3)(4)得到D(2)，直到所有的元素归为一类(距离矩阵中剩下一行一列)时，聚类结束。
2. 聚类分析的应用
2.1 正演数据聚类
现以合成波形数据为例，说明波形聚类的过程。图 1(a)(c)(e)(g)分别为Ricker子波、Morlet子波、偏移的Morlet子波以及以上3种子波叠加的波形。微地震事件在发震时频率很高，传播过程中衰减严重[10-12]，因此分别对4种子波叠加高斯白噪声模拟真实的信号，如图 1(b)(d)(f)(h)所示。可以看出，叠加噪声后的合成信号信噪比很低，波形特征较为接近。
![]()
图 1 4种类型子波与对应加噪后合成波形Figure 1. Four types of wavelets and the corresponding synthesized waveform with noise对4种不同的子波分别叠加3次高斯白噪声，如图 2所示。每次叠加噪声均为随机叠加，因此叠加噪声后波形各不相同，将12道合成数据随机打乱，其中1、5、9道为Ricker子波合成波形，2、3、8道为Morlet子波合成波形，4、6、10道为偏移Morlet子波合成波形，7、11、12道为合成子波合成波形。
聚类结果如图 3所示。根据式(1)计算各道之间欧式距离，得到距离矩阵，在初始状态下每一道各为一类。第1步，查找各类之间的最小距离，将1、9合并为一个类; 第2步，重复计算新类和剩余类的距离，将最小距离的5和1、9合并为新的类1、9、5。重复上述2个步骤后，2、8、3与4、10、6各分为同一类。分类结果中只有7没有与11、12直接归为同一类，这是由于7前面组合的新类较为相似。分类结果表明，聚类方法基本能将12道合成波形数据正确分类。
2.2 微地震信号聚类
将聚类分析方法在事件微地震监测数据中进行了应用。M井组位于山西省临汾市沁水县马壁乡(沁水盆地南部)，该井组为水平井，水平井段垂直深度在920~1 030 m之间。该地区主要煤层为3号煤，煤层深度由东南向西北方向逐渐变浅，煤层厚度2.4~7.5 m，水平井段附近无显著大断层发育[13-14]。根据声波测井曲线计算，储层P波速度为3 300~3 500 m/s，地表速度较高，P波速度为3 900~4 200 m/s。煤样测试表明，3号煤层渗透率在0.025~0.029 mD范围，整体渗透性较差。
M井组进行了水力压裂改造，以扩大煤储层的煤层气产量，同时地面进行了微地震监测。微地震监测观测系统地面投影如图 4所示，其中黄色图钉表示检波器位置，红色图钉为水平井段起始位置，绿色图钉为井底，红色实线为井轨迹。
![]()
图 4 M井组地面微地震观测阵列与井轨迹地面投影Figure 4. Geophone array of well M group and ground projection of well trajectoryM-1井第5段水力压裂地面微地震监测事件A波形如图 5所示，数据记录共43道，时窗内共800采样点，采样频率2 ms，波形显示时窗的长度为1.6 s。从波形图中可以看出，2、14、25为明显的噪声干扰道。此外，由于地面检波器分布位置不同，波形特征存在明显的差异。
对图 5的微地震进行聚类分析，结果如图 6所示。横轴坐标对应图 5中的数据序号，纵坐标表示分类距离。由于数据较多，为了能清晰地看出聚类结果，分类距离从左至右依次增大。根据分类结果中距离的变化趋势可以看出，矩形内数据道之间的距离差值较小，变化梯度低; 矩形框外的数据道相对距离变化较大，因此矩形框内的数据道可以看成一类事件。选取图 6中蓝色虚线矩形框内的数据进行叠加定位。保留的数据波形都具有相似的波形特征，如较强的P波或S波。去除数据道波形为噪声道或信号差异较大的数据波形。
从图 5的波形数据中提取3类典型特征的微地震信号，如图 7所示。
Ⅰ类信号特征具有明显P波和S波，且P波波动明显; Ⅱ类信号无明显P波，但S波振幅较强; Ⅲ类信号具有较强的P波，而无明显的S波。3类信号主要由震源机制与地面采集位置确定。由于基于振幅叠加的微地震事件成像方法与波形特征密切相关，噪声干扰道以及不同信号特征会影响最终的成像结果，因此在波形叠加过程中，一般采用绝对值振幅或平方振幅值[15]，Ⅰ类和Ⅲ类可避免极性相反造成的振幅相互抵消，但Ⅱ类在P波振幅叠加中不可避免地带来振幅削弱的影响。利用波形之间的相关性，去除了Ⅱ类数据，最终保留的数据道如图 8所示。
![]()
图 8 事件A聚类后保留数据与删除数据波形Figure 8. Retained traces and deleted traces after event A clustering2.3 信号聚类对定位成像的优化
分别对聚类筛选前后的数据道进行绕射叠加成像。绕射叠加的原理是，通过射线追踪算法计算同一个震源点到每一个检波器的时间偏移量，当不同检波器校准偏移量为同一发震时刻时，振幅叠加必然最大。其定位基本流程是，将目标区域按定位精度划分成多个小体元(假设震源点)，计算体元到每个检波器的偏移量，然后沿着偏移量进行绕射叠加，直至遍历所有体元。其成像公式[16]为
$$ \operatorname{Img}(x, \tau)=\sum\limits_{i=1}^{N} u_{i}\left[\tau+t_{i}(x)\right] $$ (4) 式中，N为检波器个数; ui为第i个检波器记录的地震记录，即原始波形信息，通常为波形振幅的绝对值或振幅平方; ti(x)为从体元到第i个检波器的理论到时; τ为发震时刻。
整个成像过程是将所有检波器波形记录上的振幅按照计算的走时进行叠加，在震源位置的网格点处叠加能量最强。
分别对聚类前后的数据进行绕射叠加成像。由于已知微地震事件所在位置，叠加时选取600 ms时窗长度进行绕射叠加。全部波形数据成像结果如图 9(a)所示，聚类后保留数据成像结果如图 9(b)所示。背景颜色为相似性叠加值，红色表示高相似性值，黄色、绿色依次降低，蓝色最低。根据叠加成像的原理，一个震源点释放的地震波信号由地面各站点接收后虽然有时差、振幅差或畸变[17]，但是具有内在相关性。绕射叠加成像结果中最大能量值位置即是震源发生的位置[18-19]。聚类分析的本质是增加事件的成像值，改善定位结果。
对比聚类前后可以看出，图 9(b)成像结果相比图 9(a)更加聚焦，红色区域过渡更集中，并且聚类保留数据的相似性值域明显提高(全部数据叠加相似性值域为0.33~0.47，筛选后数据叠加相似性值域为0.40~0.54)，从而提高了定位结果的可信度。在成像结果中，只需提取对应的最大相似性值位置或拓扑结构的中心，便可确定震源点[20]。因此，聚类筛选后的数据在成像效果和相似性值域均得到了明显改善，有利于微地震事件的进一步分析。
事件B与事件C的微地震事件成像结果如图 10、图 11所示。可以看出，聚类分析方法改善了成像的质量，提高了相似性值域范围。
3. 结论
本文将统计学中的聚类分析方法在微地震成像中进行了应用，将每道信号看作一个样本，根据信号的特征进行分类，再通过各道相关性去除P波振幅较弱的数据道。
聚类算法在实际应用中计算效率高，识别较为准确，能够快速识别和提取距离相近的数据道，无需人工初至拾取进行微地震事件成像，适合于微地震监测实时定位处理。
筛选后的地震数据成像结果明显改善，且叠加的相似性值域明显增大，提高了定位的可靠度。
致谢: 感谢深听(北京)科技有限责任公司DeepListen解释软件，为本文提供了数据处理的基础。 -
![]()
图 1 4种类型子波与对应加噪后合成波形
Figure 1. Four types of wavelets and the corresponding synthesized waveform with noise
![]()
图 4 M井组地面微地震观测阵列与井轨迹地面投影
Figure 4. Geophone array of well M group and ground projection of well trajectory
![]()
图 8 事件A聚类后保留数据与删除数据波形
Figure 8. Retained traces and deleted traces after event A clustering
-
[1] McClellan J H, Eisner L, Liu E T, et al. Array processing in microseismic monitoring: detection, enhancement, and localization of induced seismicity[J]. IEEE Signal Processing Magazine, 2018, 35(2): 99-111. DOI: 10.1109/MSP.2017.2776798
[2] 达姝瑾, 李学贵, 董宏丽, 等. 微地震震源定位方法综述[J]. 吉林大学学报: 地球科学版, 2020, 50(4): 1228-1239. Da Shujin, Li Xuegui, Dong Hongli, et al. Summary of microseismic location methods[J]. Journal of Jilin University: Earth Science Edition, 2020, 50(4): 1228-1239.
[3] 常旭, 李政, 王鹏, 等. 基于频率衰减补偿的微地震定位方法[J]. 地球物理学报, 2018, 61(1): 250-257. Chang Xu, Li Zheng, Wang Peng, et al. Micro-seismic location based on frequency attenuation compensation[J]. Chinese Journal of Geophysics, 2018, 61(1): 250-257.
[4] 刘劲松, 王赟, 姚振兴. 微地震信号到时自动拾取方法[J]. 地球物理学报, 2013, 56(5): 1660-1666. https://www.cnki.com.cn/Article/CJFDTOTAL-DQWX201305024.htm Liu Jinsong, Wang Yun, Yao Zhenxing. On micro-seismic first arrival identification: a case study[J]. Chinese Journal of Geophysics, 2013, 56(5): 1660-1666. https://www.cnki.com.cn/Article/CJFDTOTAL-DQWX201305024.htm
[5] 黄国娇, 巴晶, 钱卫. 层状TI介质中微地震定位和各向异性速度结构同时反演[J]. 地球物理学报, 2020, 63(7): 2846-2857. https://www.cnki.com.cn/Article/CJFDTOTAL-DQWX202007032.htm Huang Guojiao, Ba Jing, Qian Wei. Simultaneous inversion for anisotropic velocity structure and microseismic location in layered TI media[J]. Chinese Journal of Geophysics, 2020, 63(7): 2846-2857. https://www.cnki.com.cn/Article/CJFDTOTAL-DQWX202007032.htm
[6] 曾志毅, 张建中. 利用微地震记录互相关成像的震源定位方法[J]. 石油地球物理勘探, 2020, 55(2): 360-372, 388, 231. https://www.cnki.com.cn/Article/CJFDTOTAL-SYDQ202002015.htm Zeng Zhiyi, Zhang Jianzhong. Source location through microseismic cross-correlation imaging[J]. Oil Geophysical Prospecting, 2020, 55(2): 360-372, 388, 231. https://www.cnki.com.cn/Article/CJFDTOTAL-SYDQ202002015.htm
[7] Geller R J, Mueller C S. Four similar earthquakes in central California[J]. Geophysical Research Letters, 1980, 7(10): 821-824. http://www.onacademic.com/detail/journal_1000035778135710_07c0.html
[8] Hadi A S, Kaufman L, Rousseeuw P J. Finding groups in data: an introduction to cluster analysis[J]. Technometrics, 1992, 34(1): 111. http://esyuwapu.ru/nilyjas.pdf
[9] 任雪松, 于秀林. 多元统计分析[M]. 北京: 中国统计出版社, 2010. [10] 王璐琛, 常旭, 王一博. 干涉走时微地震震源定位方法[J]. 地球物理学报, 2016, 59(8): 3037-3045. https://www.cnki.com.cn/Article/CJFDTOTAL-DQWX201608026.htm Wang Luchen, Chang Xu, Wang Yibo. Locating micro-seismic events based on interferometric traveltime inversion[J]. Chinese Journal of Geophysics, 2016, 59(8): 3037-3045. https://www.cnki.com.cn/Article/CJFDTOTAL-DQWX201608026.htm
[11] 郑晶, 曹子原, 姜天琪, 等. 基于深度信念神经网络的微震波到时拾取方法[J]. 矿业科学学报, 2018, 3(6): 521-526. http://kykxxb.cumtb.edu.cn/article/id/180 Zheng Jing, Cao Ziyuan, Jiang Tianqi, et al. Deep belief neural network-based arrival picking for microseismic data[J]. Journal of Mining Science and Technology, 2018, 3(6): 521-526. http://kykxxb.cumtb.edu.cn/article/id/180
[12] 孙远, 杨峰, 郑晶, 等. 基于变分模态分解和小波能量熵的微震信号降噪[J]. 矿业科学学报, 2019, 4(6): 469-479. http://kykxxb.cumtb.edu.cn/article/id/248 Sun Yuan, Yang Feng, Zheng Jing, et al. Research on microseismic signal denoising based on variational mode decomposition and wavelet energy entropy[J]. Journal of Mining Science and Technology, 2019, 4(6): 469-479. http://kykxxb.cumtb.edu.cn/article/id/248
[13] 李德伟, 杨瑞召, 孟令斌, 等. 煤储层中微地震事件震级影响因素分析[J]. 煤炭学报, 2020, 45(11): 3895-3902. https://www.cnki.com.cn/Article/CJFDTOTAL-MTXB202011024.htm Li Dewei, Yang Ruizhao, Meng Lingbin, et al. Analysis of factors influencing the magnitude of microseismic events in coal reservoirs[J]. Journal of China Coal Society, 2020, 45(11): 3895-3902. https://www.cnki.com.cn/Article/CJFDTOTAL-MTXB202011024.htm
[14] 颜志丰, 琚宜文, 唐书恒, 等. 沁水盆地南部煤层气储层压裂过程数值模拟研究[J]. 地球物理学报, 2013, 56(5): 1734-1744. Yan Zhifeng, Jü Yiwen, Tang Shuheng, et al. Numerical simulation study of fracturing process in coalbed methane reservoirs in southern Qinshui basin[J]. Chinese Journal of Geophysics, 2013, 56(5): 1734-1744.
[15] 缪思钰, 张海江, 陈余宽, 等. 基于微地震定位和速度成像的页岩气水力压裂地面微地震监测[J]. 石油物探, 2019, 58(2): 262-271, 284. https://www.cnki.com.cn/Article/CJFDTOTAL-SYWT201902013.htm Miao Siyu, Zhang Haijiang, Chen Yükuan, et al. Surface microseismic monitoring of shale gas hydraulic fracturing based on microseismic location and tomography[J]. Geophysical Prospecting for Petroleum, 2019, 58(2): 262-271, 284. https://www.cnki.com.cn/Article/CJFDTOTAL-SYWT201902013.htm
[16] 田宵, 张雄, 张华, 等. 全干涉成像的微地震定位方法研究[J]. 地球物理学报, 2020, 63(8): 3105-3115. Tian Xiao, Zhang Xiong, Zhang Hua, et al. Full-interferometry imaging method for microseismic location[J]. Chinese Journal of Geophysics, 2020, 63(8): 3105-3115.
[17] 姜天琪, 裴烁瑾. 基于网格搜索-牛顿迭代法的微震震源定位算法[J]. 矿业科学学报, 2019, 4(6): 480-488. http://kykxxb.cumtb.edu.cn/article/id/249 Jiang Tianqi, Pei Shuojin. Micro-seismic event location based on Newton iteration method and grid-search method[J]. Journal of Mining Science and Technology, 2019, 4(6): 480-488. http://kykxxb.cumtb.edu.cn/article/id/249
[18] 王维波, 周瑶琪, 春兰. 地面微地震监测SET震源定位特性研究[J]. 中国石油大学学报: 自然科学版, 2012, 36(5): 45-50, 55. https://www.cnki.com.cn/Article/CJFDTOTAL-SYDX201205010.htm Wang Weibo, Zhou Yaoqi, Chun Lan. Characteristics of source localization by seismic emission tomgraphy for surface based on microseismic monitoring[J]. Journal of China University of Petroleum: Edition of Natural Science, 2012, 36(5): 45-50, 55. https://www.cnki.com.cn/Article/CJFDTOTAL-SYDX201205010.htm
[19] 张永成, 郝海金, 李兵, 等. 煤层气水平井微地震成像裂缝监测应用研究[J]. 煤田地质与勘探, 2018, 46(4): 67-71. https://www.cnki.com.cn/Article/CJFDTOTAL-MDKT201804011.htm Zhang Yongcheng, Hao Haijin, Li Bing, et al. Application of microseismic monitoring and imaging of fractures in horizontal CBM well[J]. Coal Geology & Exploration, 2018, 46(4): 67-71. https://www.cnki.com.cn/Article/CJFDTOTAL-MDKT201804011.htm
[20] 单仁亮, 宋永威, 白瑶, 等. 基于小波包变换的爆破信号能量衰减特征研究[J]. 矿业科学学报, 2018, 3(2): 119-128. http://kykxxb.cumtb.edu.cn/article/id/129 Shan Renliang, Song Yongwei, Bai Yao, et al. Research on the energy attenuation characteristics of blasting vibration signals based on wavelet packet transformation[J]. Journal of Mining Science and Technology, 2018, 3(2): 119-128. http://kykxxb.cumtb.edu.cn/article/id/129
下载:
计量
- 文章访问数: 892
- HTML全文浏览量: 240
- PDF下载量: 63
