地球物理学报  2016, Vol. 59 Issue (8): 3037-3045   PDF    
引用本文
王璐琛, 常旭, 王一博. 2016. 干涉走时微地震震源定位方法. 地球物理学报, 59(8): 3037-3045, doi: 10.6038/cjg20160826.  
WANG L C, CHANG X, WANG Y B. 2016. Locating micro-seismic events based on interferometric traveltime inversion. Chinese J. Geophys. (in Chinese), 59(8): 3037-3045, doi: 10.6038/cjg20160826.  
干涉走时微地震震源定位方法
王璐琛1,2 , 常旭1 , 王一博1     
1. 中国科学院页岩气与地质工程重点实验室, 中国科学院地质与地球物理研究所, 北京 100029;
2. 中国科学院大学, 北京 100049
收稿日期 2015-02-05, 2016-01-27 收修定稿
基金项目: 国家自然科学基金（41230317，41274112）资助.
作者简介: 王璐琛,男,1989年生,博士,主要从事地震波成像、微地震的反演和震源机制研究.E-mail:wangluchen@mail.iggcas.ac.cn
通讯作者: 常旭,研究员,主要从事地震波传播、地震成像和地球物理监测研究.E-mail:changxu@mail.iggcas.ac.cn
摘要: 本文基于地震波场干涉原理，建立了干涉走时微地震震源定位方法.该方法将两个接收点相对于一个微地震事件的走时差（称为干涉走时）的扰动作为残差函数，通过迭代求解最小残差函数，最终获得震源的空间位置.干涉走时震源定位方法利用两个接收点的到时差消除发震时刻未知和速度模型误差的影响，简化了震源定位算法.数值计算表明，本文提出的干涉走时定位方法在速度模型有误差的情况下仍然可以获得准确的微地震震源定位.
关键词: 微地震      干涉定位      震源反演      走时差     
Locating micro-seismic events based on interferometric traveltime inversion
WANG Lu-Chen1,2, CHANG Xu1, WANG Yi-Bo1     
1. Key Laboratory of Shale Gas and Geoengineering, Institute of Geology and Geophysics, Chinese Academy of Sciences, Beijing 100029, China;
2. University of Chinese Academy of Sciences, Beijing 100049, China
Abstract: To locate sources of micro-seismic events, a method based on seismic interferometric traveltime is proposed in this paper. This method builds on the concept of seismic wave-field interferometry. In this method, we define the interferometric traveltime as the difference of traveltimes from one seismic source to other two different receivers. The locations of micro-seismic events can be obtained by inverting seismic source's spatial position through iteration to obtain minimum interferometric traveltime residual. We set up a difference objective function to eliminate the seismic origin time in seismic arrival time. It can simplify the micro-seismic locating problem, because it reduces the velocity disturbance's influence on traveltime difference. The examples of the numerical calculation show that this method is stable under the velocity model disturbance, and has the advantage of interferometric locating in the multi-event problem with high accuracy..
Key words: Micro-seismic      Interferometric locating      Source inversion      Traveltime difference     
1 引言

本文针对非常规油气开发中水力压裂诱发的微地震事件的定位问题开展研究.当前对这一类诱发微地震事件的定位，主要沿用地震学研究中天然地震的定位方法.天然地震震源定位的传统方法是Geiger(1912)提出的单地震事件定位方法，该方法通过假定震源位置与发震时刻，并由此计算震源到各个检波点的计算到时.用所有检波点的计算到时与观测到时之差作为目标函数，求解目标函数的极小值即可获得震源的时空位置.地震学研究领域以Geiger的方法为基础派生了单地震事件定位的多种算法，例如：将Geiger法进行线性求解程序化的Hypo71系列(Lee and Lahr，1975)；提高Geiger法计算稳定性的阻尼最小二乘反演法(Lienert et al.，1986)；使用Backus-Gilbert反演理论的震源反演法(Klein，1978).近年来适用于复杂模型的多事件定位方法得到更广泛地研究，例如单纯形法定位(Prugger and Gendzwill，1988)、震源与台站校正反演法(Douglas，1976)、主事件定位法(Spence，1980)、双差定位法(Waldhauser and Ellsworth，2000)等.准确的震源定位是震源机制、震源几何构造和源区介质结构等地震学问题研究的基础(Nolen-Hoeksema and Ruff，2001; Aki and Richards，2002; Song and Toks z，2011)，准确的微地震震源定位是评价压裂效果和追踪压裂裂缝分布(Maxwell and Urbancic，2002; Daniels et al.，2007; Eisner et al.，2010)的重要依据.

要获得精确的定位结果，传统意义上的定位方法首先要求速度模型尽可能准确，其次要求准确地拾取走时间，从而获得纵横波不同相位的走时信息.但是这些条件在微地震定位中并不容易满足，因为通常微地震信号的频率高达数百至数千赫兹(Warpinski，1994)，数据的信噪比很低，并且介质的复杂程度很高.另外，微地震仅能传播很近的距离，因此我们仅能从有限空间中对其进行观测，这就增大了观测扰动对定位的精度的影响.近年来双差定位法的提出提高了震源定位方法的抗干扰性和精度，但是该方法中震源位置和发震时刻的反演是耦合在一起的，不利于以位置定位为主要目标的微地震震源反演.

除了以上介绍的走时类震源定位方法，另有研究使用地震波场干涉聚焦成像的方法进行震源成像定位.干涉成像有两种实现方式：一种是在数据域对地震道集进行干涉，生成虚拟道集(Schuster et al.，2004; Poliannikov et al.，2011)，进而用虚拟道集进行震源成像；另一种是在成像域使用干涉成像条件(Sava and Poliannikov，2008; 王晨龙等，2013; 李振春等，2014)，实现震源成像.采用成像的方式，避免了对记录的初至拾取，也能容易的处理发震时刻的影响(Artman，2010).Poliannikov等(2011)将干涉法应用于多震源相互定位上，李振春等(2014)的研究结果表明，干涉成像能很容易地推广至弹性波.但是成像方法的缺点在于，对于多震源同时成像问题，复杂波至会降低成像的分辨率，引入成像假象，降低震源定位的精度.另外，成像效果很大程度上取决于记录数据的丰富程度，对于观测记录较少的情况，成像效果不够理想.

本文利用同一震源事件对不同检波点的地震波走时干涉来消除震源发震时刻在反演中的影响，提高微地震事件的定位效果.本文构建了微地震干涉法走时差目标函数，以及相应的反演算法.通过数值计算并在模型中加入速度扰动，证明了干涉走时算法的稳定性和定位的精度，同时证明了该方法对速度模型存在误差时具有良好的适应性.

2 干涉走时的定义

在物理光学中，波场的干涉指的是对两个不同光波场重叠而产生新信号的这一过程.通过干涉，相同传播路径上的波场能够相互抵消，剩下的部分为光在不相同路径上传播的相位差.应用到地震学领域，就可以描述地震波在不同传播路径上的走时差.

地震波干涉法的基本思路如图 1所表示，由于路径上相近的部分传播的过程在运动学和动力学的角度上是几乎相同的，干涉后会相互抵消，所以剩下差的就是位置扰动所导致的结果.地震波波前的“震相到时”可以在频域中表示为

图 1 微地震事件的观测系统 Fig. 1 Observational system for a micro-seismic event

(1)

上式表示的是同一个震源到两个检波点各自的波场，其中表示的是频域中的震源子波，τsm和τsn分别表示震源到两个检波点的走时.将这两个波场进行干涉的话，就会得到相位iω(τsm－τsn)，这表示的恰好就是各自路径上到时的差值.这个时差可以通过道集互相关的方式或者两道初至拾取求差的方式获得，本文是采用后一种方式获得这个差值.

图 1所示的是一个简单的地震事件的观测系统，其中虚线表示从震源事件i发出的到检波点的地震射线.方程(2)所示的是针对震源事件i的微地震干涉法定位基本数学原理，公式为

(2)

其中TTmi表示震源事件i到检波点m的走时(Travel-Time)，因此TTmi－TTni表示的就是同一个震源事件到两个不同检波点的干涉走时(走时差).drmni是干涉走时的观测值与理论计算值之间的差，上标“obs”表示实际观测的干涉走时，上标“cal”表示理论计算的干涉走时.由于同一个震源事件的发震时刻是一定的，因此在方程(2)中，我们用同一个震源到不同检波点的到时Tmi来代替走时TTmi，计算求差函数Tmi－Tni，就可以消除事件i未知的发震时刻的影响，并且这个到时差也就等同于走时差TTmi－TTni.方程(2)就是本文干涉反演的目标函数，也可以说是干涉走时残差函数.

由于微地震事件是由于压裂所诱发的，其发生时间基本可以用压裂的时段来表示，不同时段压裂所诱发的微地震在时间记录上相隔较远易于分离，所以在微地震震源反演问题中，所需要反演的主要参数是地震事件的空间位置，因此去除未知的发震时刻能够减少反演变量，改善反演结果的分辨率和数据精度.

本文的干涉法同传统的双差定位方法(Waldhauser and Ellsworth，2000; Zhang and Thurber，2003; Zhou et al.，2010)不同，本文干涉法使用的是同一个震源对不同的检波点到时求差，而双差法是用不同的两个地震事件i和j对同一个检波点的到时求差，其目标函数定义为

(3)

它此时对求差事件存在的一个约束条件是，两个震源之间的距离要远小于震源到检波点的最大距离，也就是两个进行走时求差的震源在整个模型的尺度上要“足够的接近”.这种方法能去除大部分速度结构的影响，因为相近事件到同一个检波点的传播路径基本相同，所以公式(3)中的走时差异主要是源自于两个事件之间的介质结构的影响，从而转化为两个事件位置的相对距离.但是对于微地震震源定位而言，接近震源区域的介质速度结构由于缺少探查数据，是不确定性较大的参量.所以应用双差法的话，不能精确的将两个事件在震源区域介质中的关系(的相对距离)对应到二者的走时差上.

而由图 1所示的地震波传播路径可以看出，对相邻检波点应用本文的干涉法，则能去除大部分速度结构的影响.此时的走时差异主要是由靠近检波点部分的介质所引起的，而这部分的介质性质由于靠近检波点，所以能被测井等方式获得.因此本文干涉法受到介质不确定性的影响很小.

还要强调的是，本文的干涉法是针对独立的震源进行计算的，不需要震源事件存在空间位置的约束条件，得到的也是震源的实际位置而非相对位置.还有重要的一点是，从公式(3)易见，双差法是两个震源事件之间的走时求差，因此不能消除未知震源发震时刻的影响，这对微地震震源定位而言是一种缺陷.

3 干涉走时反演定位

根据公式(2)所表示的目标函数(亦称干涉走时残差函数)，通过提取不同事件的初至到时，然后求 差即可得到走时差(干涉走时).再通过反演的方式，求解给定范数的最小残差意义上的震源位置(Tarantola，1984a，1984b; Mora，1987).

从射线的角度来说，地震事件i到检波点m的到时Tmi可以表示为

(4)

其中τi表示的是震源的发震时刻，积分部分表示的是地震波在传播路径上的总走时，震源的空间扰动所导致地震波在传播时间上的扰动为

(5)

将这一扰动(称为到时残差)进行Taylor展开，可以得到

(6)

公式(6)说明，到时的变化是受到震源参数和速度模型的不确定影响的.

我们定义同一事件i到不同的检波点m，n的到时差为

(7)

这一定义同双差法(Zhang et al.，2009; Waldhauser and Ellsworth，2000)的差别在于公式(7)可消除发震时刻的不确定性对于到时差的影响，从而将未知变量缩减到空间坐标和介质模型上.同时考虑到微地震发震时刻的误差尺度远大于空间传播走时的误差尺度(微震的传播范围很小，因而空间扰动造成的走时误差尺度也会很小)，如果将发震时刻反演同位置反演结合在一起，则位置扰动所产生的到时偏差会由于太小而湮没在发震时刻误差所造成的到时偏差中，不利于对震源位置的反演，因此本文这样的处理应该是非常合适的.另外由图 1可以看到震源发出的射线在离开震源的时候路径上是比较接近的，因此在干涉的时候就能一定程度上消除这一部分介质的影响，所以模型介质的影响主要是体现在靠近检波点的区域内.而相对于震源区域，靠近检波点的介质物理性质是更容易得知的.

如果忽略速度模型不准确的影响，则到时残差可以表示为

(8)

另外对于观测数据和理论计算数据又可以将上面的残差函数表示为

(9)

其中(ΔTmni)obs表示观测数据的到时差，表示根据当前的事件空间位置所模拟出来计算数据到时差.这样根据公式(8)和公式(9)就可以建立出联合残差方程，简记为

(10)

公式(8)相当于公式(10)的左边项，公式(9)相当于公式(10)的右边项.如果我们考虑三维空间中存在N个待定位的地震事件(每个事件包含三方向的未知数)，并且检波点两两组合有M对，则公式(10)中的走时对于事件空间位置的梯度矩阵 G 为

(11)

相应的事件的位置参量和数据残差 d 分别为

(12)

(13)

所以只要求解公式(10)就可以得到微地震事件位置 的修正，通过迭代就可以反演出各个事件的空间位置.

对于微震事件理论走时的计算，本文中采用的是射线类方法中的弯曲方法(Thurber and Aki，1987; Koketsu and Sekine，1998).弯曲法相比于传统射线类方法中的试射法，具有更高的计算效率，尤其是在三维空间的传播射线计算上.弯曲法基本思想是循环调整连接震源点与检波点之间的任意初始路径(假设的射线路经)，直到成为接近真实的射线路径.本文所使用的是其中的伪弯曲射线法(Julian and Gubbins，1977; Um and Thurber，1987).该方法使用近似射线追踪确定一个估计的射线路径，然后用估计的射线路经作为射线追踪的起始点来确定射线在三维空间中的路径.

4 模型计算 4.1 已知真实速度模型的干涉走时定位

为了结合页岩气开发中的微地震研究区域的特点(监测范围较小、明显的层状地层和包含断层)，本文选择了国际上发布的Amoco理论模型中的一个小区域(图 2a中方框所示的部分)作为速度模型.所选择的是横向范围37255 m到39250 m、纵向范围2505 m到4500 m各400个网格的正方形区域，原始网格大小横向和纵向都是5 m.本文将此区域的网格大小设置为1 m，并采用模拟井中接收的方式，观测系统如图 2b所示，其中三角表示检波点，一共垂直分布21个，检波点间距18 m.

进行模拟的40个地震事件的真实位置如下面图 3a所示，实心圆点表示微地震震源位置，不同颜色的源点表明微震事件具有不同的发震时刻.图 3b是其中某一个地震事件的观测记录，对此记录进行初至拾取便可得到观测到时.另外，图中圆点表示反演迭代第20次时采用伪弯曲射线法计算的初至到时.由初至时刻便可计算该事件对不同检波器的干涉时差.实际反演中将所有的记录合并在一起计算 各事件的干涉时差，同时反演40个地震事件的位置.

图 2 (a)Amoco模型，所选择的区域为图中红框所示，(b)此区域放大后的示意图，黑色三角表示井中检波点 Fig. 2 (a) Amoco model. Selected area is shown by red box，(b)Enlarged computing area. Black triangle represents well geophone
图 3 (a)模拟的40个微地震事件的真实位置，(b)某一个微地震事件对应的记录，圆点表示计算的到时 Fig. 3 (a)True positions of simulated 40 micro-seismic events，(b)Seismic record of one micro-seismic event. Dots indicate calculated arrival times

使用本文的干涉反演方法，反演定位的结果如图 4a和图 4b所示.图 4中的空心圆圈表示微地震震源的真实位置，五角星表示事件位置的反演结果，其中相同颜色的空心圆与五角星对应于同一个地震事件，从而可以对比事件定位的效果.所有事件的迭代计算初始位置都是图 4a中间偏上部位的蓝色圆点处.

图 4 (a) 40个事件干涉走时反演的结果，(b)局部放大的效果 Fig. 4 (a) Interferometric inversion results of 40 events，(b)Local amplification of Fig. 4a

从图 4b中蓝色虚线圈出的部分可以看到：对于相对距离较大的离散事件(图 4b中圆圈圈出的地震事件)，反演的效果比较好，事件反演都较好的收敛到真实的位置.对于相对距离较近的事件(图 4b方框圈出的地震事件)，由于事件之间的相互影响，所以反演得到的震源位置会存在一定的偏差，但是依旧能维持整体的成群特征.只有图 4a中高速断层体的下部存在少部分事件点反演结果严重偏离自己的真实位置.

4.2 速度模型存在误差时的干涉走时定位

我们依旧使用图 2b的速度模型作为正演模拟，并使用图 5a的扰动速度模型进行干涉走时震源位置反演测试.图 5a的速度模型左侧1/4靠近检波点处的速度模型同图 2b的真实速度模型相近，右侧3/4处速度模型是图 2b右侧对应位置速度进行平滑的结果，以此表示右侧部分(微地震事件发生区域)的介质速度存在扰动.速度扰动相对真实值的百分比见图 5b.40个模拟微地震事件的真实位置保持不变(见图 3a).

图 5 (a) 反演模型，(b)速度扰动量的百分比 Fig. 5 Inversion model with perturbation，(b)Percentages of velocity perturbation

从图 5b可以看出，本文所用扰动模型相对于真实模型最大有25%的速度扰动，平均扰动量为4.1%.用此扰动速度模型反演的结果见图 6.

对比图 6与图 4a可以看出，当背景速度不准的时候，震源的定位结果(图 6)存在一定的误差，但震源位置总体的聚集程度依旧同真实情况一致，说明干涉法能够抵抗一定程度的模型不准确所造成的干扰.定位误差的统计结果见图 7b.图 7的统计结果说明干涉法在速度不准的时候对多震源事件反演的效果(图 7b)同真实速度反演的效果(图 7a)差别不大，证明本文方法对模型存在速度扰动时仍然有良好的适应性.

图 6 40个事件的定位结果 Fig. 6 Locating results of 40 events
图 7 干涉法反演统计结果 (a)真实速度；(b)有扰动速度. Fig. 7 Statistical result of interferometric inversion method (a)True velocity model；(b)Velocity model with perturbation.

表 1是使用准确速度模型与有误差速度模型进行反演的结果统计信息，从表里可以看出：本文方法在模型精确的条件下最小的定位偏差仅为厘米量级，表明本方法的最高精度可以完美定位事件的位 置；两种模型下反演的最大偏差比较接近，都在25 m的量级上，说明本文的干涉法反演方法不会因为模型存在误差而导致反演结果出现较大的恶化，也就是鲁棒性良好；两种模型的平均偏差都在5 m附近(有误差模型的总平均误差稍大)，说明多事件所造成的定位扰动还是比较明显，需要进一步的提高迭代求解的算法分辨率.

表 1 两种模型反演的结果统计信息 Table 1 Statistical information of inversion results by using two models
5 结论

本文提出了基于地震波场干涉原理的走时震源定位方法.运用震源事件的到时差构建目标函数，一方面避免了发震时刻的误差对于微地震定位的影响，一方面减小了地下震源区域介质模型的不确定性对定位结果的影响，提高了定位结果的准确度.针对于层内分布的微地震事件，本文干涉法反演的震源位置比较集中，与预设震源位置的误差在允许范围内.相对位置距离较近的地震事件对于反演的精度有较大的影响(个别点甚至会偏差很远)，这就要求在反演求解器上进行更深入的研究，从而提高反演算法本身的算法分辨率.速度不准确的定位算例表明，本文的干涉法抗模型扰动的能力较好.

参考文献
Aki K, Richards P G. 2002. Quantitative Seismology. 2nd ed. Sausalito, CA:University Science Books.
Artman B, Podladtchikov I, Witten B. 2010. Source location using time-reverse imaging. Geophysical Prospecting , 58(5): 861–873.
Daniels J L, Waters G A, Le Calvez J H, et al. 2007. Contacting more of the barnett shale through an integration of real-time microseismic monitoring, petrophysics, and hydraulic fracture design.//SPE Annual Technical Conference and Exhibition. Anaheim, California, U.S.A.:SPE.
Douglas A. 1976. Joint epicentre determination. Nature , 215(5096): 47–48.
Eisner L, Williams-Stroud S, Hill A, et al. 2010. Beyond the dots in the box:Microseismicity-constrained fracture models for reservoir simulation. The Leading Edge , 29(3): 326–333.
Geiger L. 1912. Probability method for the determination of earthquake epicenters from the arrival time only. Bull. St. Louis. Univ. , 8: 60–71.
Julian B R, Gubbins D. 1977. Three-dimensional seismic ray tracing. J. Geophys. , 43(1): 95–113.
Klein F W. 1978. Hypocenter location program HYPOINVERSE:Part I. Users guide to versions 1, 2, 3 and 4. Part Ⅱ. Source listings and notes. Open-File Report, 78-694. U. S. Geological Survey.
Koketsu K, Sekine S. 1998. Pseudo-bending method for three-dimensional seismic ray tracing in a spherical earth with discontinuities. Geophysical Journal International , 132(2): 339–346.
Lee W H K, Lahr J C. 1975. HYPO71:A computer program for determining hypocenter, magnitude, and first motion pattern of local earthquakes. Open-File Report, 75-311. U. S. Geological Survey.
Li Z C, Sheng G Q, Wang W B, et al. 2014. Time-reverse microseismic hypocenter location with interferometric imaging condition based on surface and downhole multi-components. Oil Geophysical Prospecting (in Chinese) , 49(4): 661–666.
Lienert B R, Berg E, Frazer L N. 1986. Hypocenter:An earthquake location method using centered, scaled, and adaptively damped least squares. Bull. Seismol. Soc. Am. , 76(3): 771–783.
Nolen-Hoeksema R C, Ruff L J. 2001. Moment tensor inversion of microseisms from the B-sand propped hydrofracture, M-site, Colorado. Tectonophysics , 336(1-4): 163–181.
Maxwell S C, Urbancic T I. 2002. Real-time 4D reservoir characterization using passive seismic data.//SPE Annual Technical Conference and Exhibition. San Antonio, Texas:SPE.
Mora P. 1987. Nonlinear two-dimensional elastic inversion of multioffset seismic data. Geophysics , 52(9): 1211–1228.
Poliannikov O V, Malcolm A E, Djikpesse H, et al. 2011. Microseismicity:Interferometric hydrofracture microseism localization using neighboring fracture. Geophysics , 76(6): WC27–WC36.
Prugger A F, Gendzwill D J. 1988. Microearthquake location:A nonlinear approach that makes use of a simplex stepping procedure. Bull. Seismol. Soc. Am. , 78(2): 799–815.
Sava P, Poliannikov O. 2008. Interferometric imaging condition for wave-equation migration. Geophysics , 73(2): S47–S61.
Schuster G T, Yu J, Sheng J, et al. 2004. Interferometric/daylight seismic imaging. Geophysical Journal International , 157(2): 838–852.
Song F, Toks z M N. 2011. Full-waveform based complete moment tensor inversion and source parameter estimation from downhole microseismic data for hydrofracture monitoring. Geophysics, 76(6):WC103-WC116.
Spence W. 1980. Relative epicenter determination using P-wave arrival-time differences. Bull. Seismol. Soc. Am. , 70(1): 171–183.
Tarantola A. 1984a. Inversion of seismic reflection data in the acoustic approximation. Geophysics , 49(8): 1259–1266.
Tarantola A. 1984b. Linearized inversion of seismic reflection data. Geophysical Prospecting , 32(6): 998–1015.
Thurber C H, Aki K. 1987. Three-dimensional seismic imaging. Annual Review of Earth and Planetary Sciences , 15: 115–139.
Um J, Thurber C. 1987. A fast algorithm for two-point seismic ray tracing. Bull. Seismol. Soc. Am. , 77(3): 972–986.
Waldhauser F, Ellsworth W L. 2000. A double-difference earthquake location algorithm:Method and application to the northern Hayward fault, California. Bull. Seismol. Soc. Am. , 90(6): 1353–1368.
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 J. Geophys. , 56(9): 3184–3196. doi: 10.6038/cjg20130931.
Warpinski N R. 1994. Interpretation of hydraulic fracture mapping experiments.//University of Tulsa Centennial Petroleum Engineering Symposium. Tulsa, Oklahoma, 291-300.
Zhang H J, Sarkar S, Toks z M N, et al. 2009. Passive seismic tomography using induced seismicity at a petroleum field in Oman. Geophysics , 74(6): WCB57–WCB69.
Zhang H J, Thurber C H. 2003. Double-difference tomography:the method and its application to the Hayward Fault, California. Bull. Seismol. Soc. Am. , 93(5): 1875–1889.
Zhou R M, Huang L J, Rutledge J. 2010. Microseismic eventlocation for monitoring injection using double-difference tomography. The Leading Edge , 29(2): 208–214.
李振春, 盛冠群, 王维波, 等. 2014. 井地联合观测多分量微地震逆时干涉定位. 石油地球物理勘探 , 49(4): 661–666.
王晨龙, 程玖兵, 尹陈, 等. 2013. 地面与井中观测条件下的微地震干涉逆时定位算法. 地球物理学报 , 56(9): 3184–3196.