2. 吉林大学仪器科学与电气工程学院, 长春 130026
2. College of Instrument Science and Electrical Engineering, Jilin University, Changchun 130026, China
0 引言
页岩油藏及煤层气是两种近几年才开始引人注目的新能源。两者的开采都需要对岩层进行压裂以释放油、气。另外,常规油田为了提高产量、榨出最后一滴油,也要进行压裂施工。因此,水力压裂成为常规手段,对压裂效果的监测尤为必要[1-4]。
井地ERT (electrical resistance tomography)、井下测斜仪、微地震监测、放射性示踪剂等方法被广泛应用于水力压裂效果评价[5-6]。其中,微地震监测方法精度更高,监测获得的裂缝评价参数更多,监测结果更具实际指导意义,因此,微地震监测在目前水力压裂效果监测中是最有效、应用最广泛的方法[7-9]。
微震监测方法可探测压裂过程中产生破裂的位置,并由此推断裂缝的走向、形状和大小。微震监测通常有两种方式:井中监测和地面监测。井中监测是大型地球物理服务公司的通用手段,技术成熟,在监测井周围几百米的范围内能达到较高的水平和垂直定位精度。地面监测是近年来开始逐渐推广的手段,在地表布设检波器阵列,对阵列范围内深度不超过2 000 m的微震有比较好的水平分辨率,而垂直定位精度则随深度增加而降低。
近些年来,微震监测的理论研究成果较多,但实际数据例子较少。本文通过对一个地表微震监测实验采集到的实际数据的处理过程来介绍地表微震监测的整个过程。
1 方法原理 1.1 微地震信号预处理在微地震事件发生时,可以按采集数据将其分为优秀的微地震事件、较差的微地震事件和无微地震事件几种情况;在大多数情况下,采集到的微地震事件以较差的微地震事件为主。从地震数据记录可以发现:优秀的微地震事件呈现为清晰的脉冲信号,而能量较弱的微地震事件则显示频率高、持续时间短的特点,并且,能量越小,地下岩层破裂的裂缝长度越短;较差的微地震事件由于能量较弱,有可能会被环境噪声所淹没,在进行地面微地震监测时,环境噪声有可能对能量产生较大的干扰。因此,在识别微地震事件时,如何区分地震记录是微地震信号还是随机干扰至关重要。同时,针对较差的微地震事件,还需要对微地震波的到时进行精确的提取与记录[10-11]。
滤波处理在地震数据预处理中是最为常用的一种处理手段,但是此种方法主要针对的对象是勘探地震中的反射地震波,由于微地震波能量十分微弱,因此很多情况下这种滤波处理方法并不能适应微地震波预处理要求[12]。多道互相关检测技术在微地震波处理中能够获得较好的滤波效果,有效排除随机干扰,从而在信噪比较低的情况下提取出有效的微地震信号[13],可有效解决上述问题。
检波器采集的地震信号由地下微地震事件所引起,源信号在地下传播过程中会受到许多干扰噪声的影响,但是由于地震波传播路径、地下介质不同等因素,通常噪声与源信号的频率并不相同[14],即这些噪声干扰与微地震事件源信号的相关性很小。那么,将源信号与检波器所采集到的地震数据进行互相关计算之后,由于微地震信号与源信号频率相同,可以将微地震信号从信噪比较低的地震数据中识别出来,从而可以计算出微地震信号由震源点到检波器之间的传播时间。由于不同检波器采集到的微地震信号具有相关性,那么将两个检波器记录的信号进行互相关计算,就能够准确得到微地震事件在不同检波器之间的到时差。
对于随机过程,两个随机信号X=(x1, x2, …,xn),Y=(y1, y2, …, yn)的相关函数RXY(τ)与相关系数RXY(0)定义如下。
1) 相关函数:
式中:n为信号采样点数;τ为信号X的时移量。式(1)代表将X移动τ值后与Y进行相关计算。若存在τm,使τ=τm时RXY(τm)为最大值,那么信号X、Y之间的时移量为τm。
2) 相关系数:
归一化相关系数:
当
微地震信号子波受选取窗口大小、中心频率、频带宽度、信噪比以及信号在传播过程中发生畸变与衰减等因素的影响,相关检测方法获得的到时差并不精确,因此,在将数据进行相关处理后,通过长短时窗能量计算方式对微震波的到时波相进行识别[15]。
1.2 微震波相对到时波相自动判别在对微地震数据进行相关计算之后,还需要对相关的效果进行评价,从而识别微地震事件的准确到时。一般来说,地震信号都存在或多或少的噪声,但是由于真实的微地震信号与干扰噪声在频率与震幅上并不相同,因此,通过相关计算之后,加强了相关性较强的真实信号,并且有效地压制了噪声信号。通过计算在特定频率范围内包含时间变量的信噪比(SNR),可以识别微地震事件的准确到时[16-17]。针对检波器采集到地震数据能量的大小,分别对长时窗能量平均值与短时窗能量平均值进行计算。设长时窗能量平均值为LTA,代表地震信号能量缓慢变化;短时窗能量平均值为STA,代表地震信号突变能量,那么SNR=STA/LTA。
长短视窗比的基本原理是利用STA与LTA的比值来表示地震信号的能量变化。当一个地震信号到达时,由于STA的变化比LTA剧烈,因此SNR会明显变大,当SNR大于预先设定的阈值时,就可以认为该次事件是有效的微地震事件,从而将此刻的时间记录为地震波在各检波器上的相对到时。
在各种STA/LTA算法中,基于递归原理的STA/LTA算法计算速度较快,并且其结果平滑,在实际中应用非常普遍[18-19]。本文就采用递归的方法来计算STA/LTA值。设检波器采集到的地震信号为X=(x1, x2, ..., xn),则
式中:m为短时窗内的采样点数;l为长时窗内的采样点数。设采样率为dt,短时窗长度为tm,长时窗长度为tl,则tm=dt·m,tl=dt·l。对实际数据进行处理时,在地震信号序列X上滑动长短时窗,如果在某一个采样点上STA/LTA值发生突变,并大于预先设定的阈值时,就可以看作在这一时刻发生了一次微地震事件,并且将这一时刻记录为该微震事件的到时时刻。
从前面的分析可以看出,基于STA/LTA方法的微地震事件拾取精度不但与地震数据信噪比有关,并且受到阈值设定与时窗长度选取等因素的影响。
下面通过一个简单的例子(图 1)来计算经过相关处理后地震数据的STA/LTA值序列,从而说明长短时窗法如何自动检测短时间内检波器接收到的地震波能量。设定阈值为6,当SNR>6时,认为有微地震发生。通过图 1可以快速得出微地震波的到达时间为0.3 s。
在震相识别的计算过程中,阈值的设定对识别精度的影响很大:阈值设定过小,就不能很好地区分噪声干扰与微地震信号;而阈值设定过大会造成对能量较弱的微地震事件识别敏感度不够,从而发生漏检或者微地震信号到时误判的情况[20-21]。因此,在此种方法的基础上引入了一种以判别频率估计算法为基础的模式识别方法。
1.3 基于判别频率估计算法的模式识别方法频率估计(frequency estimate, FE)算法是朴素贝叶斯网络在分类识别中的一种具体应用, 茫現E算法只需要对训练数据集中的实例进行一次扫描,计算实例属性不同取值在不同分类类别中的条件概率,根据训练数据集中的数据构建条件概率表格;当有一个新实例需要进行分类识别的时候,基于已经生成的条件概率表格,使用朴素贝叶斯公式就可以计算出新实例在不同分类中的条件概率。
判别频率估计(discriminant frequency estima-tion, DFE)算法是一种基于贝叶斯网络分类的判别参数学习算法。FE算法在训练学习的过程中只是简单地将条件概率表格中参数的出现频率加1,没有考虑计算实例分类效果对计数的影响。实际上,在不断修改条件概率表格的过程中每一个步骤都有一个分类基础,即应用当前条件概率表格中参数的出现频率计算得出的局部概率分布。因此,DFE算法应用当前条件概率表格对计算实例进行分类,然后根据分类结果对条件概率表格进行修改[22]。
在通过模式识别方法进行分类识别之前,需要对识别对象进行特征拾取。特征属性的提取对识别对象描述越精确,分类识别的精度越高。微地震识别可以看作是在时间域上对地震波形的模式识别,在数据采集的过程中,实时判断微地震事件是否到来。有效提高微地震事件的识别精度,就能够在定位微地震事件时具有更高的效率。
模式识别方法以STA/LTA值作为特征参数,通过分类器识别微地震事件是否到来。与传统的STA/LTA方法相比,模式识别不需要预先设定阈值,降低了人工操作的干扰。在本文中,选用STA/LTA时间窗以1 s为步长、在地震数据上连续滑动100步所生成的能量比值序列作为微地震事件的特征属性。设描述微地震事件的特征向量为A,那么A={a1, a2, ..., a100},其中ai表示STA/LTA时间窗滑动到第i步时的能量比值。在进行微地震识别时,可能的分类结果只有微地震事件和噪声干扰两种结果。在进行分类识别之前,需要利用已有的地震数据按照设定的特征参数属性与分类结果作为训练数据集来生成分类器,然后采用分类器对实时采集到的地震信号进行识别。因此,在数据采集过程中,只对地震数据的特征属性参数进行计算,分类过程中的计算量并不大,分类效率很高,能够满足实时监测的要求。
使用一组实际采集的微地震数据对本方法的可行性进行验证。实验中的训练数据集由103条地震数据记录组成,其中真实的微地震事件记录为50条,剩余的地震数据记录为噪声数据或环境干扰数据。训练数据集中, 每一条记录的特征属性由一个长度为100的STA/LTA序列组成,包含初次拾取到发震时刻的前10个STA/LTA值和后90个STA/LTA值。图 2为从训练数据集中选取的4组地震数据波形及其对应的STA/LTA结果,短时窗长度为1 s,长时窗长度为20 s。从图 2可以看出,高信噪比地震数据的STA/LTA峰值在20以上,低信噪比地震数据的STA/LTA峰值只稍大于2,而干扰信号和背景噪音的STA/LTA峰值都大于2;因此,传统的STA/LTA方法已经不能将这4个信号准确地识别出来了。然而,使用DFE算法对这4个信号进行识别则能取得较好的结果。表 1为使用基于DFE算法的模式识别方法与传统STA/LTA方法的分类结果比较。从表 1可以看出,基于DFE算法的模式识别方法分类精度高于传统的STA/LTA方法10%以上,并且不需要人工设定阈值。
为了更加精确地对微地震事件的震源位置进行定位,需要在监测压裂施工开始前对地面采集站布置进行设计。本文采取振幅叠加定位方法对震源位置进行对位,需要对所有采集站的数据根据地震波的到时进行偏移叠加;为了取得更好的叠加效果,采用大张角大尺度排列方式,最大可能地接收到地下微地震信号信息。
传统的微地震震源定位方法当中,不论是纵横波时差法还是同型波时差法,都需要对地震信号的初至事件进行精确地拾取。在理论上,这些定位方法对震源点定位的精度较高;但是由于这些方法对于每一个微地震事件都需要在每一道地震数据上进行精确的初至拾取,因此计算量很大,处理效率较低。并且,当进行地面监测时,由于地震波在传播过程中衰减严重,而且微地震信号能量很弱,因此地面检波器所采集到的微地震信号常常无法对微地震事件的初至进行有效拾取。这就导致了传统的定位方法并不适用于地面微地震监测。
本项目采用的震幅叠加方法将地震波看作常规地震勘探中的绕射波,借用逆时偏移的方法对微地震事件的震源位置进行定位。这种方法不需要对微地震事件的初至进行精确拾取,因此提高了定位效率;并且多道叠加以后,能够使原本十分微弱的微地震事件得到加强,提高了定位的可靠性。在微地震观测实验中,数据采集台网能够将微地震信号记录在一段时间窗内,微地震事件震源位置与发震时间是未知的。将有可能发生微地震事件的监视目标区域按定位精度要求划分成多个小体元,每一个体元都被认为是一个潜在的微地震事件震源点。将微地震震源点看作勘探地震中的一个绕射点,根据逆时偏移原理,将所有数据道上收到的信号按照时间方向逆向传播,必然会在微地震震源点位置发生汇聚。因此,需要遍历监视区域内的所有体元,通过震幅叠加的方法,寻找到叠加后能量极大值的位置作为微地震事件的震源点位置。具体操作如图 3所示。
遍历每一个体元,把每一道数据按照相对应的时间偏移量进行偏移,然后将所有道数据进行叠加。假设某体元为实际地震震源点,那么偏移后各道数据到时相同,震幅得到加强;而非震源点各道数据到时不同,叠加后震幅抵消。
本项目通过扫描目标区域全部体元的办法,从震幅叠加后的能量聚焦情况来判断微地震事件的震源位置与发震时刻。假设第i个体元中心坐标为(xi, yi, zi),那么将所有数据道按照第i个体元到每一数据道的传播时间进行偏移后,通过叠加得到的能量总和E(xi, yi, zi)为
式中:S(xi, yi, zi, k, j)为第i个体元中心到第k个检波器在j时刻的震幅大小;M为检波器数量;T为时窗长度。
1.5 能量聚焦方法验证通过基于长短时窗分析的判别估计方法可以对微地震事件进行识别,但是其中有可能会包含一些“假”微地震事件;因此,采用能量聚集的方法对提取出的微地震事件再次进行验证分析,以区别出微地震事件和干扰信号。
微地震震源定位方法借用常规地震勘探的逆时偏移原理,通过对采集到的地震数据波形进行震幅叠加定位。而压裂地面监测过程中的“假”微地震事件基本都来源于地面上车辆、风等的干扰,在聚焦过程中将压裂目标区域划定成小的体元,干扰信号在聚焦结束后不会出现在目标区域内,因此可以将其判定为假事件。
首先对真实微地震事件和“假”微地震事件进行能量聚焦分析,结果如图 4所示。在俯视图中,A和B两个事件均达到归一化的最大能量值(图 4a),仅通过xy方向无法分辨出哪一个是真实的事件;在侧视图中(xz、yz方向),事件A依旧存在,并且在压裂的预期范围以内,而事件B相应位置没有出现能量聚焦点(图 4b、c)。因此,判定事件B为“假”微地震事件,在后期计算中将其去除。
2 方法验证笔者通过对河北北部一水平井的压裂监测对本文方法进行验证。据地质资料,待压裂目的煤层深约1 000 m,煤层厚度2~3 m。水平井施工沿煤层顶板上部数米进行。压裂监测目的在于了解煤层是否被压开而产生裂缝,从而形成煤层气传输通道。本次压裂井区构造位置位于宣后向斜直西北翼,在低山丘陵与平原交界处,西北向冲沟发育。
本次压裂检测共进行了6阶段(6 d)的压裂。压裂井分段压裂过程中,采用微地震监测技术监测压裂裂缝几何参数,确定裂缝是水平缝或垂直缝。对于水平裂缝,确定裂缝的空间形状、裂缝在空间各方向的伸展范围;对于垂直裂缝,确定裂缝方位和裂缝长度。当天给出监测结果,以便指导下一步压裂工作的实施方案。
2.1 监测方法采用井中和地面三维空间采集方式进行信号采集。本次监测共布置了6条测线,并标注了采集站位置。采集站间距260 m,每个采集站配备4个垂直分量检波器,检波器间距20 m,沿测线向远离阵列中心的方向排列。在每条放射线测线上,第一个检波器距离阵列中心点100 m。此监测阵列的设定旨在适应1 000 m目标区域深度,需要有足够宽的视角(从不同角度接收地震信号)才能比较准确地对震源进行定位。具体采集方案见图 5,井中检波器放置于图 5水文井内,检波器放置深度为800 m。
2.2 监测成果利用每一天的监测数据对压裂进行定位计算。根据测井数据,取地震波在地层中的传播速度v=2 000 m/s。图 6为一个典型微震事件的原始数据记录,图 7为其定位计算结果。从图 7可以看出,本次压裂监测获得能量较大的微地震事件都集中于井筒附近,并且裂缝的主要延展方向为基本垂直于井筒的东北方向。这一结果与当地地质构造基本一致,符合当地地应力发育情况。由于井中检波器数量较少,因此本次监测结果在竖直方向上精度较差。
图 8为所有6 d的定位计算结果。可以看出,每天压裂效果都很理想,监测到的微地震事件数量比较丰富。每天监测获得的裂缝发育情况与压裂设计情况相符合,具有较高的可信度。
2.3 计算结果分析根据以上监测成果可以得到每天压裂的施工情况,包括裂缝的缝高和缝长(表 2)。图 9为本次压裂施工设计单位提供的压裂裂缝模拟结果示意图。在图 9中,压裂模拟缝长为74.2 m,缝高为33.2 m。可以看出,表 2中的监测结果基本与模拟结果一致,其中:第1天、第4天、第6天压裂监测到的裂缝长度略小于模拟结果,第2天、第3天、第5天压裂监测到的裂缝长度略大于模拟结果;除第1天以外,其他几次压裂监测在裂缝高度上也与模拟结果基本一致。但是,总体上来看,与模拟结果相比,监测结果的误差在垂直方向上明显大于水平方向;这也是井中检波器数量较少、垂直方向上视角不足造成的结果。
1) 本文针对微地震数据进行预处理,并基于长时窗和短时窗判别频率估计算法对微地震事件进行提取,同时针对提取事件中包含的“假”微地震事件进行能量聚焦分析,将所有干扰事件进行剔除,最终实现了真实的微地震事件拾取。
2) 采用振幅叠加对微地震事件进行定位的方法将监测目标区域划分为多个体元,利用能量聚焦对微地震事件进行定位。该方法易于并行计算,并且对地震信号质量要求低、噪声容忍能力强。通过此种定位方法可以快速有效地进行野外压裂数据的现场处理,并通过后期数据成像方法将压裂裂缝的走向、方位、长度等信息快速通过计算机展现给操作人员,实现对压裂施工的现场指导。
致谢: 感谢所有参加此次野外实验的吉林大学仪器科学与电气工程学院的师生们。[1] | Song F, Toks z M N. Full-Waveform Based Complete Moment Tensor Inversionand Source Parameter Estimation from Downhole Microseismic Data for Hydrofracture Monitoring[J]. Geophysics, 2011, 76 : 103-116. |
[2] | De Pater C J, Groenenboom J. Active Seismic Moni-toring of Hydraulic Fractures in Laboratory Experiments[J]. International Journal of Rock Mechanics and Mining Sciences, 2001, 38 (6) : 777-785. DOI:10.1016/S1365-1609(01)00042-9 |
[3] | Ge M. Efficient Mine Microseismic Monitoring[J]. International Journal of Coal Geology, 2005, 6 (1/2) : 44-56. |
[4] | Warpinski N R. Microseismic Monitoring:Inside and out[J]. Journal of Petroleum Technology, 2009, 61 : 80-85. |
[5] | 刘财, 杨宝俊, 冯晅, 等. 论油气资源的多元勘探[J]. 吉林大学学报(地球科学版), 2016, 46 (4) : 1208-1220. Liu Cai, Yang Baojun, Feng Xuan, et al. Multivariate Exploration Technology of Hydrocarbon Resources[J]. Journal of Jilin University (Earth Science Edition), 2016, 46 (4) : 1208-1220. |
[6] | Bulant P, Eisner L, Pšencík I, et al. Importance of Borehole Deviation Surveys for Monitoring of Hydraulic Fracturing Treatments[J]. Geophysical Prospecting, 2007, 55 : 891-899. DOI:10.1111/gpr.2007.55.issue-6 |
[7] | Chambers K, Kendall J M, Brandsberg-Dahl S, et al. Testing the Ability of Surface Arrays to Monitor Microseismic Activity[J]. Geophysical Prospecting, 2010, 58 : 821-830. DOI:10.1111/j.1365-2478.2010.00893.x |
[8] | Li Junlun, Sadi Kuleli H, Zhang Haijiang, et al. Focal Mechanism Determination of Induced Microearth-quakes in an Oil Field Using Full Waveforms from Shallow and Deep Seismic Networks[J]. Geophysics, 2011, 76 : 87-101. |
[9] | Maxwell S C, Rutledge J, Jones R, et al. Petroleum Reservoir Characterization Using Downhole Microseismic Monitoring[J]. Geophysics, 2010, 75 : 129-137. DOI:10.1190/1.3507237 |
[10] | Zhang Shan, Liu Qinglin, Zhao Qun, et al. Aplica-tion of Microseismic Monitoring Technology in Development of Oil Field[J]. Geophysical Prospecting for Petroleum, 2002, 41 (2) : 226-231. |
[11] | Rutledge J T, Phillips W S. Hydraulic Stimulations of Natural Fracture as Revealed by Induced Microearthquakes Carthage Cotton Valley Gas Field East Texas[J]. Geophysics, 2003, 68 : 441-452. DOI:10.1190/1.1567214 |
[12] | Moriya H. Phase-only Correlation of Time-Varying Spectral Representations of Microseismic Data for Identification of Similar Seismic Events[J]. Geophysics, 2011, 76 : 37-45. |
[13] | Grechka V, Pech A, Tsvankin I. P-Wave Stacking Velocity Tomography for VTI Media[J]. Geophysical Prospecting, 2002, 50 : 151-168. DOI:10.1046/j.1365-2478.2002.00307.x |
[14] | Gei D, Eisner L, Suhadolc P. Feasibility of Estimating Vertical Transverse Isotropy from Microseismic Data Recorded by Surface Monitoring Arrays[J]. Geophysics, 2011, 76 : 117-126. |
[15] | 张山, 刘清林, 赵群, 等. 微地震监测技术在油田开发中的应用[J]. 石油物探, 2002, 41 (2) : 226-231. Zhang Shan, Liu Qinglin, Zhao Qun, et al. Microseismic Monitoring Technology in Oil Field Development[J]. Geophysical Prospecting for Petroleum, 2002, 41 (2) : 226-231. |
[16] | Ren T X, Reddish D J. Styles P Numerical Modelling and Microseismic Monitoring to Improve Strata Behaviour[J]. Computer Applications in the Minerals Industries, 2001, 5 : 1-10. |
[17] | Maxwell S C, Urbancic T I. The Role of Passive Microseismic Monitoring in the Instrumented Oil Field[C]//SPE Annual Technical Conference and Exhibition. San Antonio:[s. n.], 2002:85-87. |
[18] | 叶根喜, 姜福兴, 杨淑华. 时窗能量特征法拾取微地震波初始到时的可行性研究[J]. 地球物理学报, 2008, 51 (5) : 1574-1581. Ye Genxi, Jiang Fuxing, Yang Shuhua. Window Energy Feature Pickup Microseismic Wave First Arrival of Feasibility Study[J]. Chinese Journal of Geophysics, 2008, 51 (5) : 1574-1581. |
[19] | 高淑芳, 李山有, 武东坡, 等. 一种改进的STA/LTA震相自动识别方法[J]. 世界地震工程, 2008, 24 (2) : 37-41. Gao Shufang, Li Shanyou, Wu Dongpo, et al. An Improved STA/LTA Automatic Recognition Method of Seismic Phase[J]. World Earthquake Engineering, 2008, 24 (2) : 37-41. |
[20] | Li Zhengguang, Yang Runhai, Zhao Jinming, et al. Rock-Braking Test Research on Earthquake Sequence Type[J]. Journal of Seismological Research, 2005, 28 (4) : 388-393. |
[21] | 胡光书. 数字信号处理:理论、算法与实现[M]. 北京: 清华大学出版社, 1997 . Hu Guangshu. Digital Signal Processing:Theory, Algorithm and Implementation[M]. Beijing: Tsinghua University Press, 1997. |
[22] | 吕昊.基于油田压裂微地震监测的震相识别与震源定位方法研究[D].长春:吉林大学, 2012. Lü Hao. Research on the Seismic Phase Identification and the Source Location Based on Oilfield Fracture Microseismic Monitoring[D]. Changchun:Jilin University, 2012. |