2020, Vol. 63
2. 中国科学技术大学地球与空间科学学院, 合肥 230026;
3. 东方地球物理公司, 涿州 072751
2. School of Earth and Space Sciences, University of Science and Technology of China, Hefei 230026, China;
3. BGP Inc., China national petroleum corporation, Zhuozhou 072751, China
微地震监测是近年来发展的一种新的油气田开发技术.早在20世纪20年代, 人们就已经注意到油气开采过程会诱发或者触发地震.一直到近70年代, 微震监测技术才发展起来用于油气田开发.水力压裂产生的微震事件能量很弱, 震级一般在-3至0级之间, 大多数微震事件的频率范围在200~1500 Hz之间, 且持续时间比较短(小于1 s)(Warpinski, 2009; Maxwell, 2010；赵博雄等, 2014; 梁北援等, 2015).微地震监测使用不同布局的检波器持续监测水力压裂过程产生的微地震事件, 检波器可以放置在地面也可以放置在观测井中(Eisner et al., 2009, 2010; 朱亚东洋等, 2017).其中，井中监测是最常用的微地震监测方式, 压裂井和监测井之间的距离一般不大于800 m, 可以采用单井或多井监测系统.微地震监测方法一般通过微震事件的分布来监测压裂区域裂缝的发育情况和几何形状, 评价压裂效果, 因此获得准确的震源位置是微地震监测的基础(Grechka, 2010; Maxwell, 2010; Zimmer, 2011; Liang et al., 2016).
目前常用的定位方法有很多种(毛庆辉等，2019), 比如基于概率统计的贝叶斯方法(Zhang et al., 2017)、基于走时的反演方法(Geiger, 1912; Jones et al., 2014; Zhou et al., 2015; Tian et al., 2016; Chen et al., 2017)及基于波形的定位方法(王晨龙等, 2013; 王璐琛等, 2016; Zheng et al., 2016; Li et al., 2016, 2017；王向腾等, 2019)等.基于走时的定位方法首先需要拾取P波和S波的初至, 通过找到与初至最吻合的空间位置来进行定位, 通常一个破裂过程可以产生上千个微震事件, 因此手动拾取会耗费大量的时间.其次, 由于井下监测的信号较弱, 低信噪比导致的拾取误差会影响定位精度.后来, 人们借鉴地震勘探中的偏移成像思想发展了不需要拾取走时的偏移定位方法, 该方法将微震事件的能量聚焦到空间网格点上, 能量最强的网格点就是震源位置, 低信噪比的数据通过叠加可以压制噪声的影响, 因此近年来被广泛应用在微地震定位中(Kao and Shan, 2004; Baker et al., 2005).
偏移定位方法主要基于绕射叠加或者克希霍夫叠加, 依据成像点到接收点的理论走时叠加波形中相关的信息(Kao and Shan, 2004; Gharti et al., 2010).绕射叠加使用P波和S波的波形信息进行成像时需要搜索震源的发震时刻, 所有时间切片中能量最大的切片对应的时间被认为是该地震的发震时刻, 该切片中能量最强的网格点为震源位置.然而搜索发震时刻会耗费大量的时间，且发震时刻的误差也会影响最后的定位精度.Li等(2015)借鉴地震干涉成像的思想对震源进行定位, 利用波形互相关提取的走时差信息可避免搜索发震时刻, 同时根据微地震井下监测数据的特点对不同的走时差信息进行加权来提高定位精度.王璐琛等(2016)构建干涉走时差目标函数求取震源位置, 利用相同震源到不同检波器的走时差信息可以避免反演发震时刻, 同时也降低了震源区域速度模型误差对定位结果的影响.
本文在Li等(2015)提出的加权干涉成像方法的基础上提出一种全干涉成像方法, 即同时使用地震记录的互相关和自相关道集.对于某个地震事件, 将所有检波器接收到的波形信息做互相关可以提取震源在不同检波器的4类走时差信息, 而地震记录的自相关可以得到震源到相同检波器的S波和P波到时差, 简称为S-P走时差.全干涉成像方法添加的S-P走时差可以提高震源成像精度.本文首先介绍了绕射叠加和干涉成像两种方法的基本原理；然后采用一个井下监测的理论模型分析了三种方法的成像特点、发震时刻误差对绕射叠加定位结果的影响、采用蒙特卡洛随机方法添加噪音进而研究每种方法定位结果的不确定性；最后, 将三种方法应用到单井监测的水力压裂实际数据中.
1 基本原理 1.1 绕射叠加绕射叠加思想是根据成像点到接收点的理论走时叠加波形中相关的信息(Kao and Shan, 2004; Gharti et al., 2010).该方法首先把研究区域细分成很多网格点, 每个网格点都可能是震源位置, 然后用射线追踪计算每一个网格点到检波器的理论走时, 最后根据理论走时计算每个网格点的叠加能量作为成像值，公式为
|
(1) |
式中, N为检波器个数, ui为第i个检波器记录到的地震记录, 可以用原始波形信息, 也可以使用不同的波形函数(Trojanowski and Eisner, 2017), 如波形振幅的绝对值, 波形振幅的平方, 波形的包络, 相似系数, 长短时窗的比值等.ti(x)是从网格点x到第i个检波器的理论到时, τ是发震时刻.整个成像过程是将所有检波器波形记录上的振幅按照计算的走时进行叠加, 在震源位置的网格点处叠加能量最强, 该方法需要搜索发震时刻.
1.2 干涉成像为了避免搜索发震时刻, 基于走时的网格搜索方法使用S-P走时差和检波器之间的走时差来进行定位(Zhou et al., 2015, 2016).Li等(2015)提出的加权弹性波干涉成像法, 利用互相关提取检波器之间的走时差可以避免搜索发震时刻.本文在Li等(2015)的基础上, 增加了检波器道集自相关提取的S-P走时差来提高震源成像精度，表达式为
|
(2) |
|
(3) |
|
(4) |
式(2)中Ci, j代表第i个检波器和第j个检波器的波形互相关道集, 包含了4种走时差信息, 即P波和S波走时差, P波走时差, S波走时差以及S波和P波走时差, 如图 1所示.由于井下监测时S波信号强于P波, 图 1显示的互相关记录中P波走时差难以识别, 因此本文将P波走时差权重赋予零, 只使用其他三种走时差信息进行定位.公式(3)中的AiSP是第i个检波器的波形自相关记录, 包含震源在该检波器上的S-P走时差.成像时仅用自相关道集的右半边数据, 并且删掉零时移处的大振幅值避免干扰成像结果.本文把公式(2)和公式(3)相乘作为最后的成像结果(相乘可以突出不同成像结果中共同的强叠加振幅值), 如公式(4)所示.由于该方法同时使用了互相关和自相关道集, 包含了所有组合的走时差信息, 因此本文称为全干涉成像方法(Full-Interferometry Imaging Method), 简称为FIIM.为了方便对比, Li等(2015)提出的加权干涉成像方法(Weighted-Interferometry Imaging Method)在本文中简称为WIIM.
|
图 1 全干涉成像方法示意图 Fig. 1 Sketch of the full-interferometry imaging method |
如图 2a所示, 首先设计一个二维的井下监测理论模型.包括一个震源和12个检波器, 震源子波为雷克子波, 中心频率为60 Hz, 震源位置为(400, 325)m, 检波器间距为15 m分布在垂直井中.速度模型为水平层状、各向同性介质, 大小为500×500 m, 两个分界面的深度分别为100 m和200 m, 本文采用二维有限差分算法计算震源到检波器的波场(Graves, 1996; Zhang and Chen, 2006), 图 2b是正演的波形(包括X和Z分量), 横坐标代表记录到的走时, 纵坐标代表检波器序号, 从图中可以清晰地看到P波和S波的初至, 且S波的能量比P波强很多.为了提高计算效率, 先用射线追踪法计算目标区域内的所有网格点到检波器的走时表, 网格间距为2 m.
|
图 2 (a) 井中监测的速度模型及检波器-震源位置；(b)正演波形.上图为X分量, 下图为Z分量 Fig. 2 (a) Velocity model and source-receiver geometry of the downhole monitoring well; (b) Forward modeling waveforms, including the X (upper) and Z (lower) componentse |
假设震源的发震时刻已知, 用图 2b中的P波和S波信息分别进行成像, 图 3a、b、c分别为P波、S波以及绕射叠加(即两者相乘)的成像结果, 为了方便对比, 图中展示的为振幅值均一化后的成像结果.绕射叠加成像结果表明：能量聚焦的区域和真实的震源位置(加号)是吻合的.
|
图 3 (a) P波、(b)S波和(c)绕射叠加的成像结果(振幅均一化处理) Fig. 3 Normalized imaging results of the P wave (a), S wave (b) and diffraction (c) stacking |
以上的测试假设震源的发震时刻已知, 然而在处理实际数据时, 震源的发震时刻是未知的, 因此成像的同时需按一定的步长搜索发震时刻, 搜索范围内所有时间切片中能量最大的切片对应的时间为发震时刻, 该切片中能量最强的网格点为震源位置.为了研究发震时刻精度对绕射叠加方法定位的影响, 将真实的发震时刻设置为0.1 s, 不同的发震时刻对成像结果的影响如图 4所示.由图 4分析可知, 当发震时刻为0.05 s时, 定位误差高达200 m；当发震时刻为0.11 s时, 定位误差约为30 m；只有当发震时刻误差在1 ms左右时, 定位误差可以控制在10 m以内.由此可见, 发震时刻的精度对定位结果的误差影响非常大.然而, 发震时刻的精度依赖于搜索步长，步长越小，获得的发震时刻越准确，但是运算成本也越高.因此实际应用中选择合适的搜索步长非常重要, 以均衡运算成本和成像结果的准确性.
|
图 4 不同的发震时刻对应的定位误差 Fig. 4 Location errors with different origin times |
接下来使用干涉成像方法进行定位.首先, 对原始波形数据进行相关计算.检波器记录的波形自相关结果如图 5a所示, 波形自相关系数在零时移处有最大的相关系数值, 但是有用的信息是P波和S波的时间延迟(如图 5a中黑色箭头所指), 因此在成像时需要截断零时移对应的数据避免干扰定位结果.图 5b是所有检波器记录的波形互相关, 12个检波器对应有12×(12-1)/2=66条互相关道集, 从图 5b中可以清晰地看到检波器之间的S波和P波的走时差, S波走时差以及P波和S波走时差, 而P波走时差由于能量太弱难以识别.利用上述互相关和自相关道集分别进行S-P走时差、加权干涉和全干涉成像, 结果如图 6所示.
|
图 5 (a) 波形自相关的结果(X分量)；(b)检波器之间波形互相关的结果(X分量) Fig. 5 (a) Auto-correlograms of the X component; (b) Cross-correlograms of X component |
|
图 6 (a) S-P走时差、(b)加权弹性波干涉成像方法(WIIM)以及(c)全干涉成像方法(FIIM)振幅均一化的成像结果；(d)沿着图(c)中AB和CD方向WIIM和FIIM的振幅值均一化后的对比结果 Fig. 6 Normalized imaging results of the S-P differential times(a), weighted-interferometry(b), and full-interferometry imaging methods(c); (d) Comparison of the normalized amplitudes along the lines AB and CD |
图 6a为S-P走时差的成像结果, 与Zhou等(2015)基于拾取的网格搜索方法中S-P走时差的走时残差分布相似.图 6b为加权干涉成像结果, 可以看出该方法很好的将能量聚焦在真实的震源附近.图 6c为全干涉成像方法的成像结果, 即将S-P走时差(图 6a)和加权干涉(图 6b)的成像结果相乘.图 6d对比了加权干涉和全干涉成像结果沿着AB和CD方向(图 6c)的能量值, 两种方法的能量值均做了均一化处理.由于S-P走时差的成像结果主要分布在垂直于震源-检波器方向, 相比于加权干涉成像方法, 全干涉成像方法可以降低震源-检波器方向(AB)上伪震源的能量, 而垂直于震源-检波器方向上(CD)两种方法的结果相似.因此, 全干涉成像方法得到准确的震源位置的同时可以降低伪震源的干扰.值得注意的是该方法虽然省去了搜索发震时刻, 但是成像公式中添加了自相关项, 因此全干涉成像法的计算效率低于加权干涉成像.
2.3 噪声测试接下来对图 2b的理论波形添加噪声进行测试, 在原始波形上分别添加噪声水平为最高振幅20%、30%、35%及40%的高斯随机噪声(Li et al., 2020).本文采用蒙特卡洛方法研究定位的不确定性(Maxwell, 2009; Tian et al., 2016, 2017).图 7给出了100次蒙特卡洛模拟随机噪声的绕射叠加成像(绿色点), 加权干涉成像(蓝色点)和全干涉成像方法(红色点)的定位结果, 其中黑色加号为真实的震源位置.从图 7可以看出绕射叠加法的误差分布成条带状, 而干涉成像方法的误差分布成椭圆形.相比于绕射叠加方法, 加权干涉成像方法的定位误差随着噪声水平的加大而明显增加, 尤其在沿着震源-检波器方向, 因此加权干涉成像方法的抗噪能力较弱.而全干涉成像方法中添加的S-P走时差信息可以降低震源-检波器方向的定位误差, 尤其在35%和40%的高噪声水平下, 相比于加权干涉成像方法提高了定位精度.另外, 从20%和30%噪声水平下的定位结果分析, 两种干涉成像方法的精确度均高于绕射叠加方法, 本文猜测这是由于绕射叠加方法使用N个数据定位时, 干涉成像方法利用互相关可以得到N×(N-1)/2个数据, 而全干涉法可以使用N+N×(N-1)/2个数据, 数据量的增加可以一定程度地提高成像精度.
|
图 7 100次蒙特卡洛模拟随机噪声 对原始波形分别添加最高振幅的20%、30%、35%以及40%的高斯随机噪声.蓝色、红色和绿色分别代表加权干涉成像, 全干涉成像和绕射叠加方法的定位结果.其中黑色加号为真实的震源位置. Fig. 7 Monte Carlo simulation by 100 times with noises Random Gauss noise of 20%, 30%, 35%, and 40% of the maximum amplitude are added into original waveforms, respectively. The blue, red, and green dots indicate the weighted interferometry imaging, the full-interferometry imaging and diffraction imaging results, respectively. The black plus sign is the true source. |
本节把上述讨论的三种方法应用到单井监测的水力压裂实际数据中, 该套数据由东方地球物理公司(BGP)提供.在进行定位之前需建立一个速度模型, 图 8中灰色实线为声波测井曲线, 从中提取一个包括P波和S波的五层速度模型(黑色实线).检波器分布深度范围为2800~3020 m, 间距为20 m.图 9a为监测到的一个微震事件的三分量波形图, Z分量上可以同时看到P波和S波, 而X和Y分量上的S波信号较强, P波信号较弱, 因此我们选择Z分量做波形自相关提取S-P走时差.Z分量波形自相关和三分量数据波形互相关的结果如图 9b所示, 其中互相关道集包含66条记录而自相关道集包含12条记录.
|
图 8 测井曲线以及提取的速度模型(包括P波和S波) 三角形代表检波器的深度坐标. Fig. 8 Sonic-logging data and the extracted velocity model(including P and S waves) Triangles represent the depth of receivers. |
|
图 9 (a) 单井监测到的三分量数据；(b) X、Y和Z分量互相关道集以及Z分量自相关道集 Fig. 9 (a) Three-component data obtained by single well monitoring; (b) The cross-correlograms of X, Y and Z components and the auto-correlograms of the Z components |
在使用绕射叠加方法时, P波成像使用Z分量的数据而S波成像使用X和Y分量的数据.正如2.1小节中讨论的, 由于发震时刻未知, 绕射叠加方法需要搜索发震时刻.野外数据采集时检波器记录到的波形图是连续的, 数据处理时需要进行事件检测并对原始记录图进行截断, 因此发震时刻的搜索需要考虑正负值.发震时刻为负值说明在事件记录(截断后)的零时刻之前事件就已经发生了.实际数据的发震时刻搜索步长一般为1~2个采样率.为了节省搜索时间, 本文采用“先粗略, 再精细”的搜索方式.第一步的搜索范围为-0.1~0.1 s, 步长为10 ms, 不同发震时刻对应的成像结果如图 10a所示.结果表明当发震时刻为-0.06 s时, 对应切面的叠加能量最强.考虑到搜索效率以及1 ms的采样率, 第二步的搜索范围将以-0.06 s为中心左右展开, 范围为-0.069 s到-0.051 s, 步长减小为2 ms以提高发震时刻搜索精度.通过对比成像结果的幅值, 可计算出当发震时刻为-0.065 s时对应的叠加能量最大.图 10b为两次搜索过程中不同的发震时刻对应的成像切面的最大振幅值, 黑色实心点为每次搜索时最优的发震时刻.
|
图 10 绕射叠加方法中发震时刻的搜索过程 (a)成像结果的3D时间切片；(b)不同的发震时刻对应的成像切面的最大振幅值. Fig. 10 Process of searching the origin time during the diffraction imaging (a) 3D slice view of imaging profiles; (b) The curve of different origin times and their corresponding maximum amplitudes. |
从压裂段选取30个微震事件进行定位.首先使用图 8中提取的速度模型计算每个网格点到检波器的走时表, 然后采用三种方法分别成像, 最后从成像剖面中选取叠加能量最大的点确定为震源位置.三种方法在单井监测时只能得到二维的定位结果, 可使用BGP提供的方位角信息将二维定位结果转化到三维空间, 如图 11所示, 图中绿色点, 蓝色点和红色点分别代表绕射叠加, 加权干涉和全干涉成像法的定位结果, 黑色点为BGP提供的基于走时的网格搜索方法的定位结果.结果表明三种方法得到的事件分布很相似, 基本都聚集在压裂井(黑色实线)的附近, 因此从整体上来说, 三种方法的定位结果是合理的.但注意到加权干涉成像在第28个地震事件上的定位结果(图 11中箭头所指的)和其他两种方法的结果相差很大.图 12给出了第28个地震的加权干涉和全干涉的二维成像结果, 可以看到加权干涉成像结果在(278, 3130)m处出现了一个比合理位置能量更强的伪震源, 从而导致定位结果的偏差.由于全干涉成像方法添加了S-P走时差信息, 在震源-检波器方向上增加了约束, 使得叠加能量最终落在了较为合理的位置上(402, 3198)m.从计算效率(Intel Core i7-7500 CPU)来看, 绕射叠加方法耗时约3.2 h, 加权干涉成像方法耗时约1.0 h, 全干涉成像法耗时约1.5 h.实际数据的应用结果表明全干涉成像方法相比于加权干涉法可以降低伪震源的干扰, 提供更合理的定位结果, 同时计算效率高于绕射叠加方法.
|
图 11 不同方法的定位结果 Fig. 11 Location results with different methods |
|
图 12 第28个微震事件振幅均一化的成像结果 (a) 全干涉成像法的结果.加号代表成像结果中能量最强的点；(b) 加权干涉成像法的结果. Fig. 12 The normalized imaging results of the 28th micro-seismic events (a) Full-interferometry imaging. The black plus sign represents the strongest stacking energy,which is regarded as the final location；(b) Weighted-interferometry imaging. |
本文深入探讨了绕射叠加方法, 加权干涉成像方法和本文提出的全干涉成像方法在微地震单井监测中的应用.绕射叠加方法直接使用P波和S波进行成像的同时需要搜索发震时刻, 计算效率低且发震时刻误差会增加定位结果的不确定性.加权干涉成像方法利用波形互相关提取了事件在不同检波器的走时差信息, 提高了计算效率, 但是在高噪声水平下沿着震源-检波器方向易出现高能量的伪震源.添加了S-P走时差的全干涉成像方法能够有效地降低震源-检波器方向的定位误差.理论模型测试和实际数据结果表明全干涉成像方法可以提供准确的震源位置, 同时避免了搜索发震时刻, 可以广泛地应用于微地震单井监测.
致谢 感谢GeoTomo公司提供实际数据处理软件MiVuTM，感谢张伟教授为本研究提供的有限差分正演程序.
Baker T, Granat R, Clayton R W. 2005. Real-time earthquake location using Kirchhoff reconstruction. Bulletin of the Seismological Society of America, 95(2): 699-707. DOI:10.1785/0120040123 |
Chen Y K, Zhang H J, Miao Y Y, et al. 2017. Back azimuth constrained double-difference seismic location and tomography for downhole microseismic monitoring. Physics of the Earth and Planetary Interiors, 264: 35-46. DOI:10.1016/j.pepi.2016.10.003 |
Eisner L, Duncan P M, Heigl W M, et al. 2009. Uncertainties in passive seismic monitoring. The Leading Edge, 28(6): 648-655. DOI:10.1190/1.3148403 |
Eisner L, Hulsey B J, Duncan P M, et al. 2010. Comparison ofsurface and borehole locations of induced seismicity. Geophysical Prospecting, 58(5): 809-820. DOI:10.1111/j.1365-2478.2010.00867.x |
Geiger L. 1912. Probability method for the determination of earthquake epicenters from the arrival time only. Bull St Louis Univ, 8(1): 56-71. |
Gharti H N, Oye V, Roth M, et al. 2010. Automated microearthquake location using envelope stacking and robust global optimization. Geophysics, 75(4): MA27-MA46. DOI:10.1190/1.3432784 |
Graves R W. 1996. Simulating seismic wave propagation in 3D elastic media using staggered-grid finite differences. Bulletin of the Seismological Society of America, 86(4): 1091-1106. |
Grechka V. 2010. Data-acquisition design for microseismic monitoring. The Leading Edge, 29(3): 278-282. DOI:10.1190/1.3353723 |
Jones G A, Kendall J M, Bastow I D, et al. 2014. Locating microseismic events using borehole data. Geophysical Prospecting, 62(1): 34-49. DOI:10.1111/1365-2478.12076 |
Kao H, Shan S J. 2004. The source-scanning algorithm:Mapping the distribution of seismic sources in time and space. Geophysical Journal International, 157(2): 589-594. DOI:10.1111/j.1365-246X.2004.02276.x |
Li G, Liu X Q, Tang J T, et al. 2020. Improved shift-invariant sparse coding for noise attenuation of magnetotelluric data. Earth, Planets and Space, 72(1): 4311-4322. DOI:10.1186/s40623-020-01173-7 |
Li L, Chen H, Wang X M. 2015. Weighted-elastic-wave interferometric imaging of microseismic source location. Applied Geophysics, 12(2): 221-234. DOI:10.1007/s11770-015-0479-z |
Li L, Chen H, Wang X M. 2016. Relative elastic interferometric imaging for microseismic source location. Journal of Geophysics and Engineering, 13(5): 733-744. DOI:10.1088/1742-2132/13/5/733 |
Li L, Becker D, Chen H, et al. 2017. A systematic analysis of correlation-based seismic location methods. Geophysical Journal International, 212(1): 659-678. |
Liang B Y, Shen C, Leng C B, et al. 2015. Development of microseismic monitoring for hydro-fracturing. Progress in Geophysics (in Chinese), 30(1): 401-410. DOI:10.6038/pg20150158 |
Liang C T, Yu Y Y, Yang Y H, et al. 2016. Joint inversion of source location and focal mechanism of microseismicity. Geophysics, 81(2): KS41-KS49. |
Mao Q H, Wang P, Zeng J. 2019. Review of hydro-fracturing microseismic event location methods. Progress in Geophysics (in Chinese), 34(5): 1878-1886. DOI:10.6038/pg2019CC0240 |
Maxwell S. 2010. Microseismic:Growth born from success. The Leading Edge, 29(3): 338-343. DOI:10.1190/1.3353732 |
Maxwell S C. 2009. Microseismic location uncertainty. CSEG Recorder, 34(4): 41-46. |
Tian X, Zhang W, Zhang J. 2016. Cross double-difference inversion for microseismic event location using data from a single monitoring well. Geophysics, 81(5): KS183-KS194. DOI:10.1190/geo2016-0198.1 |
Tian X, Zhang W, Zhang J. 2017. Cross double-difference inversion for simultaneous velocity model update and microseismic event location. Geophysical Prospecting, 65(5): 259-273. |
Trojanowski J, Eisner L. 2017. Comparison of migration-based location and detection methods for microseismic events. Geophysical Prospecting, 65(1): 47-63. DOI:10.1111/1365-2478.12366 |
Wang C L, Cheng J B, Yin C, et al. 2013. Microseismic events location of surface and borehole observation with reverse-time focusing using interferometry technique. Chinese Journal of Geophysics (in Chinese), 56(9): 3184-3196. DOI:10.6038/cjg20130931 |
Wang L C, Chang X, Wang Y B. 2016. Locating micro-seismic events based on interferometric traveltime inversion. Chinese Journal of Geophysics (in Chinese), 59(8): 3037-3045. DOI:10.6038/cjg20160826 |
Wang X T, Ni S D, Zhou Y, et al. 2019. Topography effects for focal depth inversion of shallow earthquakes based on waveforms:A case study for the M6. 3 event on 3 September 2017 in Democratic People's Republic of Korea. Chinese Journal of Geophysics (in Chinese), 62(12): 4684-4695. DOI:10.6038/cjg2019M0474 |
Warpinski N. 2009. Microseismic monitoring inside and out. Journal of Petroleum Technology, 61(11): 80-85. DOI:10.2118/118537-JPT |
Zhang W, Chen X F. 2006. Traction image method for irregular free surface boundaries in finite difference seismic wave simulation. Geophysical Journal International, 167(1): 337-353. DOI:10.1111/j.1365-246X.2006.03113.x |
Zhang Z S, Rector J W, Nava M J. 2017. Simultaneous inversion of multiple microseismic data for event locations and velocity model with Bayesian inference. Geophysics, 82(3): KS27-KS39. |
Zhao B X, Wang Z R, Liu R, et al. 2014. Review of Microseismic monitoring technology research. Progress in Geophysics (in Chinese), 29(4): 1882-1888. DOI:10.6038/pg20140454 |
Zheng Y K, Wang Y B, Chang X. 2016. Wave equation based microseismic source location and velocity inversion. Physics of the Earth and Planetary Interiors, 261: 46-53. DOI:10.1016/j.pepi.2016.07.003 |
Zhou H, Zhang W, Zhang J. 2016. Downhole microseismic monitoring for low signal-to-noise ratio events. Journal of Geophysics and Engineering, 13(5): 805-816. DOI:10.1088/1742-2132/13/5/805 |
Zhou W, Wang L S, Guan L P, et al. 2015. Microseismic event location using an inverse method of joint P-S phase arrival difference and P-wave arrival difference in a borehole system. Journal of Geophysics and Engineering, 12(2): 220-226. |
Zhu Y D Y, Wang J L, Sun F, et al. 2017. Micro-seismic monitoring and instrument for hydraulic fracturing in the low-permeability oilfield. Chinese Journal of Geophysics (in Chinese), 60(11): 4282-4293. DOI:10.6038/cjg20171116 |
Zimmer U. 2011. Microseismic design studies. Geophysics, 76(6): WC17-WC25. DOI:10.1190/geo2011-0004.1 |
梁北援, 沈琛, 冷传波, 等. 2015. 微地震压裂监测技术研发进展. 地球物理学进展, 30(1): 401-410. DOI:10.6038/pg20150158 |
毛庆辉, 王鹏, 曾隽. 2019. 水力压裂微地震事件定位方法综述. 地球物理学进展, 34(5): 1878-1886. DOI:10.6038/pg2019CC0240 |
王晨龙, 程玖兵, 尹陈, 等. 2013. 地面与井中观测条件下的微地震干涉逆时定位算法. 地球物理学报, 56(9): 3184-3196. DOI:10.6038/cjg20130931 |
王璐琛, 常旭, 王一博. 2016. 干涉走时微地震震源定位方法. 地球物理学报, 59(8): 3037-3045. DOI:10.6038/cjg20160826 |
王向腾, 倪四道, 周勇, 等. 2019. 地形起伏对基于地震波形的浅源地震深度反演影响——以2017年9月3日朝鲜M6. 3事件为例.地球物理学报, 62(12): 4684-4695. DOI:10.6038/cjg2019M0474 |
赵博雄, 王忠仁, 刘瑞, 等. 2014. 国内外微地震监测技术综述. 地球物理学进展, 29(4): 1882-1888. DOI:10.6038/pg20140454 |
朱亚东洋, 王金磊, 孙峰, 等. 2017. 水力压裂微地震井地联合监测系统及仪器. 地球物理学报, 60(11): 4282-4293. DOI:10.6038/cjg20171116 |