地球物理学报  2012, Vol. 55 Issue (1): 238-251   PDF    
引用本文
李文杰, 曲寿利, 魏修成, 等. 非零偏VSP弹性波叠前逆时深度偏移技术探讨. 地球物理学报, 2012, 55(1): 238-251, doi: 10.6038/j.issn.0001-5733.2012.01.023.  
LI W J, QU S L, WEI X C, et al. Discussion about reverse-time depth migration technology of prestack elastic wavefield for offset VSP data. Chinese J. Geophys. (in Chinese), 2012, 55(1): 238-251, doi: 10.6038/j.issn.0001-5733.2012.01.023.  
非零偏VSP弹性波叠前逆时深度偏移技术探讨
李文杰, 曲寿利, 魏修成, 宁俊瑞     
中国石油化工股份有限公司石油勘探开发研究院, 北京 100083
收稿日期 2010-07-01, 2011-11-01 收修定稿
基金项目: 国家高技术研究发展计划(863)(2006AA06Z203)资助
作者简介: 李文杰,男,1966年生,2004年毕业于中国石油大学, 博士,主要从事地震资料处理方法研究. E-mail:liwenjie2006beijing@yahoo.com.cn
摘要: 非零偏VSP地震资料是一种多分量资料,处理非零偏VSP资料,弹性波叠前逆时深度偏移技术无疑是最适合的处理技术.本文从二维各向同性介质的弹性波波动方程出发,研究了对非零偏VSP资料进行叠前逆时深度偏移的偏移算法,讨论了逆时传播过程中的边值问题和数值频散问题及其相应的解决方案;采用求解程函方程计算得到地下各点的地震波初至时间作为成像时间,实现了非零偏VSP资料的叠前逆时深度偏移.最后进行了模型试算和非零偏VSP地震资料的试处理,结果表明该方法不受地层倾角限制,较适用于高陡构造地区或介质横向速度变化较大地区的非零偏VSP地震资料处理.
关键词: 非零偏VSP      各向同性介质      弹性波场      纵波      转换波      叠前逆时深度偏移     
Discussion about reverse-time depth migration technology of prestack elastic wavefield for offset VSP data
LI Wen-Jie, QU Shou-Li, WEI Xiu-Cheng, NING Jun-Rui     
Exploration & Production Research Institute, SINOPEC, Beijing 100083, China
Abstract: Offset VSP seismic data is a kind of multicomponent data, the elastic prestack reverse-time depth migration technology is no doubt the best suitable for processing offset VSP data. This paper studied prestack reverse-time depth migration algorithms based on elastic wave equation in 2D isotropic media for offset VSP data, and discussed the boundary value problems and frequency dispersion problems and their relevant solution in wavefield reverse-time propagation, then used the solution of eikonal equation to get the one-way traveltime from the source to each point in the image space, which is regarded as the image time at that point, and complemented the prestack reverse-time depth migration for offset VSP data. Finally, we carried out model test and processed offset VSP data. The processing results proved that the scheme is not influenced by lay inclination and suited to process offset VSP data from the areas with steep subsurface structures and big variation of horizontal velocities.
Key words: Offset VSP      Isotropic media      Scalar elastic wavefield      P-wave      Converted wave      Prestack reverse-time depth migration     
1 引 言

随着非零偏VSP 勘探技术在地震勘探领域的推广，有关非零偏VSP地震资料处理技术的研究也越来越受到业界的重视，目前，对非零偏VSP 资料的处理主要还是采用常规的VSP-CDP 转换或Kirchhoff积分偏移，但这两种方法都有其局限性:一,VSP-CDP转换仅对水平层状介质或横向变速缓慢的介质才能取得比较好的处理结果；二，利用Kirchhoff积分偏移方法对VSP 数据成像时，要求速度场横向变化不能太大.由此可见，这两类方法只适用于横向变速较为缓慢的情况，如何寻找一种适用范围更广的非零偏VSP 地震资料处理技术一直成为研究人员探索的课题.近年来，弹性波叠前逆时深度偏移技术在多波多分量地震资料处理技术领域成为了人们研究的热门技术，这主要得益于逆时偏移算法具有不同于其他偏移算法以下的优点:首先，逆时算法不受地下界面倾角的限制，适宜于各种倾角地层的波场偏移；其次，这种算法能处理横向速度变化较大地区的地震波偏移问题；再者，逆时偏移算法精度较高，算法稳健；除此以外，逆时偏移方法在处理回转波方面也有着其他偏移算法无法比拟的优点.由此可知，对于处理属于多波多分量地震资料的非零偏VSP资料，弹性波叠前逆时深度偏移技术是最适合的处理技术.在弹性波叠前逆时深度偏移技术理论研究领域，国外起步较早，文献[1-6]展示了这方面的研究成果，在国内，也有不少研究者开始从事该领域以及相关领域的研究工作，得到了许多有价值的研究成果[7-12]，笔者在文献[12]中深入探讨了基于地面多波多分量地震资料的弹性波叠前逆时深度偏移方法技术.然而，虽然到目前为止，人们在弹性波叠前逆时深度偏移技术研究方面取得了一系列的进展，但总体来说，国内外在这一领域还处在起步阶段，还需要投入更多的精力来研究它.

本文研究了二维各向同性介质中两分量的非零偏VSP叠前逆时深度偏移方法，首先从弹性波波动方程出发，利用有限差分法推导了各向同性介质中弹性波逆时延拓的有限差分格式，并在反传播过程中进行纵、横波波场分离；并利用程函方程的有限差分法计算模型中各网格点的初至时间作为成像条件，实现了各向同性介质中非零偏VSP地震资料的叠前逆时深度偏移.

2 非零偏VSP 弹性波叠前逆时深度偏移方法原理

非零偏VSP叠前逆时深度偏移技术主要包括两部分内容:一是波场的逆时传播，二是偏移成像.逆时传播是地震波数值模拟的反问题，这时已知的是井中各接收点的波场值，假设t≥tL(tL为井中各接收点的最大记录时间)时成像区域内的各点波场值为零，由此可以将弹性波场的逆时传播转化为式(1)所示的边值问题和初值问题来解决:

(1)

其中，u(x，z，t)，w(x，z，t)分别为速度模型某个网格点(x，z)所对应的位移的水平分量和垂直分量；fu(z，t)，fw(z，t)分别为井中接收到的波场位移的水平分量值和垂直分量值；x，z分别为水平距离和垂直距离；x=0表示VSP接收井所在的位置.

2.1 弹性波场逆时传播算法

在各向同性介质条件下，假设X和Z是二维介质情况下沿X轴和Z轴的坐标，X轴的方向向右，Z轴的方向朝下，可以利用两个共轭的二阶偏微分方程来描述介质中P波的运动以及垂向偏振的SV波的运动，在这里不考虑沿水平方向偏振的SH 波，并假定U，W分别为水平方向与垂直方向的位移，ρ为介质的密度，t为时间，λ 和μ 为具体介质中的拉梅系数，那么可以得到关于各向同性非均匀介质的弹性波动方程[13]:

(2)

假定密度ρ 为常数，这样就可以把方程组(2)看成是随空间位置变化的P 波和SV 波速度的函数.其中λ 和μ 与纵、横波速度α(x，z)和β(x，z)的关系如下:

(3)

Kelly[13]给出了上述二维各向同性介质中的弹性波方程正演的差分格式，由此便可以推导出弹性波方程逆时传播的差分格式如下:

(4)

(5)

2.2 稳定性条件

通过平面谐波分析,(4)式和(5)式的稳定性条 件如下[4]:

(6)

其中，h= min(Δx，Δz)

2.3 边值问题与数值频散问题的处理

在逆时传播过程中，将井中接收到的地震波场作为波源，从最大时间向最小时间进行波场逆时传播，同样遇到边值问题，左边为接收井，波场值已知，其他三边的边值问题待解决，在此，利用通过旁轴近似推导的吸收边界方程(详细推导过程可参考文献[14])来解决逆时传播过程中所遇到的边界问题.由于逆时传播的传播方向为由外向内，因此四边的吸收边界方程与正演时有所不同，他们所对应的吸收边界方程如下:

速度模型顶部所对应的吸收边界方程为

(7a)

速度模型底部所对应的吸收边界方程为

(7b)

其中，

速度模型右边所对应的吸收边界方程为

(7c)

其中，式(7c)所对应的参数项为

方程(7a-7c)为非零偏VSP 地震波场逆时传播过程中其他三边所对应的吸收边界方程.

在波场反传播数值模拟中，除了边值问题外，还需要解决数值频散的问题，数值频散现象也称为网格发散，它是速度随频率而发生的一种正常变化，相对于低频地震波来说，高频地震波被衰减了，地震波波形出现明显的拖尾现象.一般来说，数值频散问题随着网格间距的增大而变得更加严重.到目前为止，解决数值频散问题的方法有多种[15-18]，本文采用文献[18]中的通量传输校正技术消除逆时传播过程中的数值频散现象取得了较好的效果.

2.4 纵、横波分离

一般来说，井中接收到的波场记录是弹性波场的水平分量和垂直分量(即水平位移和垂直位移)，其中水平分量和垂直分量中都包含了纵波成分和横波成分，为了能得到真正的纵、横波叠前深度偏移剖面，需要进行波场分解，本文采用文献[6]所描述的方法来进行波场分解.根据纵、横波所具有的动力学特征，纵波可以利用弹性波场位移的散度来表示:

(8)

其中u，w分别表示弹性波场的水平位移和垂直位移.

而转换波则可以用弹性波场位移的旋度来表示:

(9)

2.5 成像条件

偏移成像条件的计算是弹性波逆时偏移的一个重要组成部分，所谓成像条件是指在地震波场逆时传播过程中地下各点的成像时间.成像条件可以分为以下三种:(1)零时间成像条件[19]；(2)激发时间成像条件[6]；(3)波阻抗成像条件[20].通常，弹性波叠前逆时深度偏移采用激发时间成像条件.

激发时间成像条件最初是由Claerbout在文献[21]中提出的，它的核心思想是:下行波的到达时间等于上行波的出发时间.这个概念是建立在反射波成像的基础上的，它对多次波、转换波成像都是适用的，求取激发时间成像条件就是求取地下各个网格点的初至时间，求取地下各个网格点的初至波到达时间有多种途径[22-26]，人们通常采用的方法主要有解析法和射线追踪法，但利用这两种方法计算初至时间经常存在以下几方面的局限性:(1)解析法只能适用于少数特殊的速度分布，对于一些较复杂的速度分布则难以应用；(2)利用射线追踪法求取初至时间有时会出现盲区，特别是针对一些较复杂的速度分布区域尤为如此；(3)利用以上方法计算的初至时间可能只是局部极小时间，并非真正的初至时间；(4)利用以上方法比较耗时，特别是在计算多炮检距时往往需要大量的时间.为了克服以上方法求取初至时间的局限性，有研究者提出了一种利用有限差分近似程函方程的方法来计算初至时间，利用该方法计算初至时间不会出现盲区，同时可以大大降低计算量，并具有较高的计算精度，后来人们对此方法进行了不同程度的改进，本文所用的计算初至时间方法则为其中的一种，该方法利用有限差分近似求解程函方程并与局部算法相结合来计算初至时间，文献[26]中对该方法进行了详细的描述，数值模拟结果表明:该方法具有较高的精度，能适合于任意复杂的速度模型.

3 数值模拟

为了了解前面所述的叠前逆时深度偏移方法在非零偏VSP资料处理中的应用效果，设计了三个各向同性弹性介质速度模型，利用弹性波场正演的方法得到模型的非零偏VSP的合成记录.假设震源只激发纵波，震源子波为主频30Hz的雷克子波.激发源位于地表以下10 m 的深度处，接收井位于速度模型的左边，时间采样间隔为1 ms，空间采样间隔为5m，三个模型的偏移距分别为450m，1000m，1000m，并假定反射横波仅为P-SV 转换波.从合成记录中分离出上行波场(即反射波场)和下行波场，然后利用前面所述的叠前逆时深度偏移的方法，对模型的上行波场进行叠前逆时深度偏移，分别得到偏移后的纵、横波深度偏移剖面.

图 1为一大倾角倾斜层各向同性弹性介质速度模型及该模型的各个网格点初至走时等值线图，图 2 (a,b)分别为该模型正演波场的水平分量和垂直分量，图 2 (c,d)分别为经过波场分离后的反射波场的水平分量和垂直分量，图 3(a,b)分别为对上行波场进行偏移成像所得到的纵、横波逆时深度偏移剖面.

图 1 倾斜层介质速度模型（a)及其地震波 走时初至等时线图（b)(等值线单位：ms) Fig. 1 Slant layer velocity model and its seismic first arrival traveltime contoursCtime^s unit： ms)
图 2 由倾斜层速度模型得到的非零偏VSP弹性波波场分量及波场分离后的上行波场 a水平分量；（b)垂直分量；（)上行波水平分量；（d)上行波垂直分量. Fig. 2 The offset VSP elastic wavefield components of velocity model in Fig. 1 and corresponding upgoing wavefield components (a)X-component； (b) Z-component； () X-component of upgoing wavefield；d) Z-component of upgoing wavefield.
图 3 对图 2中的反射波场进行逆时 深度偏移所得到的深度剖面 a反射p波偏移剖面；（b)反射psv波偏移剖面. Fig. 3 Depth sections obtained by applying reverse time depth migration to the reflected wavefield in Fig. 2 (a)P wave depth section ； (b) PSV wave depth section.

图 4为三层水平层状各向同性弹性介质速度模型及该模型的各个网格点初至走时等值线图，图 5(a,b)分别为该模型正演波场的水平分量和垂直分量，图 5 (c,d)分别为经过波场分离后的反射波场的水平分量和垂直分量，图 6(a,b)分别为对上行波场进行偏移成像所得到的纵、横波逆时深度偏移剖面.

图 4 三层水平层状介质速度模型（a)及其地震波走时初至等时线图（b)(等值线单位：ms) Fig. 4 Three layers velocity model and its seismic first arrival traveltime contours (time's unit: ms)
图 5 由水平层速度模型得到的非零偏VSP弹性波波场分量及波场分离后的上行波场 (a)水平分量；（b)垂直分量；（c)上行波水平分量；（d)上行波垂直分量. Fig. 5 The offset VSP elastic wavefield components of velocity model in Fig. 4 and corresponding upgoing wavefield components (a)X-component； (b) Z-component； (c) X-component of upgoing wavefield；d) Z-component of upgoing wavefield.
图 6 对图 5中的反射波场进行逆时深度偏移所得到的深度剖面 (a)反射p波偏移剖面；(b)反射psv波偏移剖面. Fig. 6 Depth sections obtained by applying reverse time depth migration to the reflected waveiield in Fig. 5 (a)P wave depth section ; (b) PSV wave depth section.

图 7为各向同性弹性介质断层速度模型(模型第一层与第二层之间为一个正断层，第二层与第三层的分界面为一水平界面)及模型的各个网格点初至走时等值线图，图 8(a,b)分别为图 7所示的断层模型正演波场的水平分量和垂直分量，图 8(c,d)为波场分离后的上行波波场分量，图 9(a,b)为上行波场进行偏移成像所得到的纵、横波深度剖面.

图 7 断层模型（a)及其地震波走时初至等时线图（b)(等值线单位：ms) Fig. 7 Fault model (a) and its seismic first arrival traveltime contours (b) (time’s unit: ms)
图 8 由断层速度模型得到的非零偏VSP弹性波波场分量及波场分离后的上行波场 ()水平分量；（b)垂直分量；（)上行波水平分量；（d)上行波垂直分量. Fig. 8 The offset VSP elastic wavefield components of velocity model in Fig. 7 and corresponding upgoing wavefield components (a)X-component; (b) Z-component; (c) X-component of upgoing wavefield;(d) Z-component of upgoing wavefield.
图 9 对图 8中的反射波场进行逆时深度偏移所得到的深度剖面 (a)PP波偏移剖面；(b)PSV波偏移剖面. Fig. 9 Depth sections obtained by applying reverse time depth migration to the reflected wavefield (a)PP wave depth section；(b)PSV wave depth section.

从偏移结果上可以看到，反射纵波与转换波的逆时深度偏移剖面能正确地反映模型的构造形态和精确位置，但根据地震波的传播规律，在VSP 测井过程中，只能接收到来自地下地层某部分的地震信息，因此只能对某部分地层进行成像，这也是VSP测井技术本身的缺陷.

4 处理实例

为了了解本文提出的非零偏VSP 弹性波叠前逆时深度偏移方法在实际生产应用中的效果，利用该方法对某探区一口非零偏VSP 测井的地震资料进行了偏移试处理.

该井接收点的分布范围为1000~5640 m 之间，偏移距为2913m，共设置278个接收点，中间出现变观，观测段5640~4800 m 之间的道间距为10m，观测段4800~1000m 之间的道间距为20m，记录的时间采样间隔为1 ms，仅选用X分量和Z分量进行处理(如图 10 所示)，在实际资料处理中，假设震源只激发纵波，横波仅为转换波.

图 10 X分量(a)和Z分量(b)非零偏VSP原始地震资料 Fig. 10 Orignal offset VSP seismic data

为了实现利用弹性波叠前逆时深度偏移方法对以上的非零偏VSP 资料进行处理，设计了图 11 所示的简易非零偏VSP资料偏移处理流程.

图 11 利用弹性波叠前逆时深度偏移方法 处理非零偏VSP资料简易流程 Fig. 11 Simple flow of processing offset VSP data using elastic prestack reverse Time depth migration method

图 10为某井的非零偏VSP资料的X分量和Z分量的原始记录，从图上可以看到:记录初至起跳干脆，根据图 11所示的处理流程，整个非零偏VSP 资料的逆时深度偏移主要包括以下几个步骤:

(1)资料预处理:该处理过程主要包括给零偏和非零偏VSP原始记录加观测系统、剔除零偏和非零偏

VSP资料中的不正常道和空道，初至切除，静校正，去噪.(2)初至拾取:拾取预处理后的零偏VSP 地震记录中各道的初至波到达时间.

(3) 层速度计算:利用零偏VSP 地震记录中各道的初至时间计算VSP井位置所在地层的纵、横波层速度(如图 12所示).

图 12 利用零偏VSP资料的初至时间计算出的地层的纵、横波层速度 Fig. 12 P-wave and S-wave interval velocities computed using first arrival time picked with zero-offset VSP data

(4) 建立偏移速度模型:利用VSP 井所在位置地层的纵、横波层速度并结合过井测线剖面建立叠前逆时深度偏移所用的纵、横波速度模型.

(5) 波场分离:利用中值滤波法或F-K滤波法对非零偏VSP资料的X分量和Z分量记录进行波场分离得到逆时深度偏移所用的上行弹性波场.图 13为利用F-K滤波法对非零偏VSP原始波场进行波场分离后所得到的上行波场的水平分量和垂直分量.

图 13 X分量(a)和Z分量(b)非零偏VSP地震资料的上行波场 Fig. 13 Upward wavefield of offset VSP data

(6)计算成像条件:利用偏移速度模型，计算从激发点到速度模型中各个网格点的初至波到达时间，并以此作为偏移成像的成像条件.

(7) 波场反传播:利用前面所述的波场反传播方法和偏移速度模型对上行波场进行波场反传播.

(8) 根据计算的成像条件，进行偏移成像，得到纵、横波偏移深度剖面.图 14 分别为利用弹性波叠前逆时深度偏移方法得到的反射纵波、转换波的偏移深度剖面.

图 14 PP波（a)和PSV波（b)深度偏移剖面（道间距：10 m) Fig. 14 PP and VSP migrated depth sections(the interval between traces s 10 m)

(9) 对偏移深度剖面进行切除得到切除后的偏移深度剖面(如图 15所示)，然后将切除后的深度剖面进行插值和深-时转换，得到插值后的时间剖面(如图 16所示，图 17为某公司处理的反射纵波和转换波的VSPCDP剖面).

图 15 PP波（a)和PSV波（b)经过切除后的深度偏移剖面（道间距：0 m) Fig. 15 The cutted VSP migrated depth sections(the interval between traces is 10 m)
图 16 PP波（a)和PSV波（b)经过插值和时深转换后的时间剖面（道间距：m) Fig. 16 The intepolated VSP migrated sections after depth to time conversion(the interval between traces is 10 m)
图 17 某公司处理的PP波（a)和PSV波（b)反射纵波与转换波VSPCDP剖面 Fig. 17 PP-wave and PSV-wave VSPCDP section processed by a company

(10)将偏移后的时间剖面嵌入二维过井剖面，图 18、图 19、图 20分别为二维过井剖面，嵌入PP波偏移剖面的过井剖面，嵌入PSV波偏移剖面的过井剖面.通过上面的处理实例，可以看出，利用弹性波叠

图 18 过井测线时间偏移剖面 Fig. 18 2D surface seismic cross-well section
图 19 嵌人PP偏移剖面的过井测线时间偏移剖面 Fig. 19 PP-wave migrated section pasted in 2D surface seismic cross-well section
图 20 嵌人PSV偏移剖面的过井测线时间偏移剖面 Fig. 20 PSV-wave migrated section pasted in 2D surface seismic cross-well section

前逆时深度偏移方法处理非零偏VSP资料，其偏移结果较好地反映了地下介质的构造特征，可以在一定程度上改善非零偏VSP资料的处理质量，提高地震勘探精度.

5 结 论

本文研究了非零偏VSP 叠前弹性波逆时深度偏移方法，该方法不受地层倾角限制，较适用于高陡构造地区或介质横向速度变化较大地区的非零偏VSP地震资料处理；从数值模拟结果可以看出该方法能使断层速度模型精确成像，实际资料的试处理也表明该方法与其他处理方法相比具有较高的分辨率和较高的信噪比.

参考文献
[1] Kuo J T, Dai T F. Kirchholf elastic wave migration for the case of noncoincident source and receiver. Geophysics , 1984, 49: 1223-1238. DOI:10.1190/1.1441751
[2] Dai T F, and Kuo J T. Real data results of Kirchholf elastic wave migration. Geophysics , 1986, 51: 1006-1011. DOI:10.1190/1.1442139
[3] Sun R, McMechan G A. Prestack reverse-time migration for elastic waves with application to synthetic offset vertical seismic profiles. Proc.IEEE , 1986, 74: 457-465. DOI:10.1109/PROC.1986.13486
[4] Chang W F, McMechan G A. Elastic reverse-time migration. Geophysics , 1987, 52: 1365-1375. DOI:10.1190/1.1442249
[5] Chang W F, McMechan G A. 3D elastic prestack reverse-time depth migration. Geophysics , 1994, 59: 597-609. DOI:10.1190/1.1443620
[6] Sun R, McMechan G A. Scalar reverse-time depth migration of prestack elastic seismic data. Geophysics , 2001, 66: 1519-1527. DOI:10.1190/1.1487098
[7] 秦福浩. 弹性波克希霍夫积分偏移法. 地球物理学报 , 1988, 31(5): 56–78. Qin F H. Prestack Kirchhoff integral migration of elastic wave. Chinese J. Geophys. (in Chinese) , 1988, 31(5): 56-78.
[8] 张关泉, 周洪波. 地面记录的纵横波分离. 石油地球物理勘探 , 1995, 30(2): 181–191. Zhang G Q, Zhou H B. Separating compressional wave from shear wave in Surface seismic records. Oil Geophysical Prospecting (in Chinese) , 1995, 30(2): 181-191.
[9] 魏修成, 刘洋. 高精度转换波速度分析. 石油勘探与开发 , 2000, 27(2): 57–59. Liu Y, Wei X C. High precision velocity analysis of converted wave. Petroleum Exploraton and Development (in Chinese) , 2000, 27(2): 57-59.
[10] 刘洋, 魏修成. 转换波三参数速度分析及动校正方法. 石油地球物理勘探 , 2005, 40(5): 504–509. Liu Y, Wei X C. Three-parameter velocity analysis of converted wave and NMO method. Oil Geophysical Prospecting (in Chinese) , 2005, 40(5): 504-509.
[11] 李文杰, 魏修成, 宁俊瑞, 等. 叠前弹性波逆时深度偏移及波场分离技术探讨. 物探化探计算技术 , 2008, 30(6): 447–456. Li W J, Wei X C, Ning J R, et al. The discussion of reverse time depth migration and wavefield separ-ation for prestack elastic data. Computing Techniques for Geophysical and Geochemical Exploration (in Chinese) , 2008, 30(6): 447-456.
[12] 李文杰, 魏修成, 宁俊瑞, 等. 弹性波叠前逆时深度偏移技术及其应用. 中国石油大学学报(自然科学版) , 2009, 33(4): 52–58. Li W J, Wei X C, Ning J R, et al. Reverse-time depth migration technology of prestack elastic wave-field and its application. Journal of China University of Petroleum(Edition of Natural Science) (in Chinese) , 2009, 33(4): 52-58.
[13] Kelly K R, Ward R W, Sven Treitel, et al. Synthetic seismograms: a finite difference approach. Geophysics , 1976, 41(1): 2-27. DOI:10.1190/1.1440605
[14] 李文杰, 魏修成, 宁俊瑞. 一种新的弹性波数值模拟吸收边界条件. 石油地球物理勘探 , 2009, 44(4): 501–507. Li W J, Wei X C, Ning J R. A new numeric simulating absorbing boundary condition of elastic wave. Oil Geophysical Prospecting (in Chinese) , 2009, 44(4): 501-507.
[15] Boris J P, Book D L. Flux-corrected transport.I.SHASTA:A fluid transport algorithm that works. Comput.Phys , 1973, 11: 38-69. DOI:10.1016/0021-9991(73)90147-2
[16] Book D l, Boris J p, Hain K. Flux-corrected transport.II: generalization of the method. Comput.Phys , 1975, 18: 248-283. DOI:10.1016/0021-9991(75)90002-9
[17] 杨顶辉, 滕吉文. 各向异性介质中三分量地震记录的FCT有限差分模拟. 石油地球物理勘探 , 1997, 32(2): 181–190. Yang D H, Teng J W. FCT finite difference modeling of three-component seismic records in anisotropic medium. Oil Geophysical Prospecting (in Chinese) , 1997, 32(2): 181-190.
[18] 李文杰, 张改兰, 姜大建. 利用FCT方法压制弹性数值模拟中的数值频散. 物探化探计算技术 , 2011, 33(3): 248–251. Li W J, Zhang G L, Jiang D J. Suppressing numerical dispersion in finite difference modeling based on acoustic and elastic equation. Computing Techniques for Geophysical and Geochemical Exploration (in Chinese) , 2011, 33(3): 248-251.
[19] McMechan G A. Migration by extrapolation of time-dependent boundary values. Geophysical Prospecting , 1983, 31(2): 413-420.
[20] Steve T, Hildebrand. Reverse-time depth migration:Impedance imaging condition. Geophysics , 1987, 52(8): 1060-1064. DOI:10.1190/1.1442371
[21] Claerbout J F. Toward a unified theory of reflector mapping. Geophysics , 1971, 36(2): 467-481.
[22] 张霖斌, 刘迎曦. 有限差分法射线追踪. 石油地球物理勘探 , 1994, 29(5): 673–677. Zhang L B, Liu Y X. Ray tracing of finite difference method. Oil Geophysical Prospecting (in Chinese) , 1994, 29(5): 673-677.
[23] John Vidale. Finite difference calculation of travel times. Bull.Seis.Soc.Am , 1988, 78(6): 2062-2076.
[24] Fuhao Qin, Yi Luo, Kim B. Finite-difference solution of the eikonal equation along expanding wavefronts. Geophysics , 1992, 57(3): 478-487. DOI:10.1190/1.1443263
[25] 刘洪, 孟繁林, 李幼铭. 计算最少走时和射线路径的界面网全局方法. 地球物理学报 , 1995, 38(6): 823–831. Liu H, Meng F L, Li Y M. The interface grid method for seeking global minimum travel-time and the correspondent raypath. Chinese J. Geophys. (in Chinese) , 1995, 38(6): 823-831.
[26] 李文杰, 魏修成, 宁俊瑞, 等. 一种改进的利用程函方程计算旅行时的方法. 石油地球物理勘探 , 2008, 43(5): 589–594. Li W J, Wei X C, Ning J R, et al. A method using improved eikonal equation to computer travel time. Oil Geophysical Prospecting (in Chinese) , 2008, 43(5): 589-594.