基于行波分离的起伏地表逆时偏移成像

摘要：传统的逆时偏移成像方法对起伏地形的适应性不强，无法对复杂地表条件下的地下构造准确成像。而全波场互相关成像条件会在偏移结果中产生低频噪声，以及在地下构造速度变化剧烈的区域，还可能产生偏移假象。为了解决这些问题，提出基于行波分离的起伏地表逆时偏移成像算法。在这种算法中，首先推导曲线坐标系下的声波方程，采用贴体网格和Lebedev网格实现了起伏地表下的高精度数值模拟；然后利用Poynting矢量对地震波场进行行波分离，并推导相应的互相关成像条件；最终选择合适的波场成分进行互相关成像。结果表明，该算法可以消除起伏地表的影响、压制低频噪声和偏移假象，从而获得高信噪比的成像结果。

关键词：行波分离； Poyting矢量； 贴体网格； 起伏地表； 逆时偏移

中图分类号：P 631.4""" 文献标志码：A

引用格式：卫二祥，李振春，黄建平，等. 基于行波分离的起伏地表逆时偏移成像［J］. 中国石油大学学报（自然科学版），2025，49（1）：72-80.

WEI Erxiang， LI Zhenchun， HUANG Jianping， et al. Reverse time migration based on wavefield decomposition with surface topography［J］. Journal of China University of Petroleum （Edition of Natural Science） ，2025，49（1）：72-80.

Reverse time migration based on wavefield decomposition

with surface topography

WEI Erxiang1，2， LI Zhenchun1，2， HUANG Jianping1，2， WANG Jiushuan3， YANG Jing3

（1.State Key Laboratory of Deep Oil and Gas， China University of Petroleum （East China）， Qingdao 266580， China;

2.School of Geosciences， China University of Petroleum （East China）， Qingdao 266580， China;

3.Bureau of Geophysical Prospecting INC， China National Petroleum Corporation， Zhuozhou 072751， China）

Abstract： Conventional reverse time migration struggles to adapt to undulating terrains， leading to inaccuracies in imaging subsurface structures under complex surface conditions. Additionally， full-wavefield imaging conditions often introduce low-frequency noise and migration artifacts， particularly in areas significant velocity variations in the underground structure. To address these challenges， a reverse time migration imaging algorithm for undulating surfaces based on wavefield decomposition is proposed. This algorithm begins by deriving acoustic wave equations in a curvilinear coordinate system， enabling high-precision numerical simulations for undulating surfaces using body-fitted and Lebedev grids. Subsequently， wavefield decomposition is performed utilizing the Poynting vector， and corresponding cross-correlation imaging conditions are established. Finally， appropriate wavefield components are selected for cross-correlation imaging. The results demonstrate that the proposed algorithm effectively eliminates the impact of undulating surfaces， suppresses low-frequency noise， and reduces migration artifacts. As a result， it achieves high signal-to-noise ratio imaging， providing more accurate and reliable subsurface structure imaging.

Keywords： wavefield decomposition; Poynting vector; body-fitted grid; surface topography; reverse time migration

收稿日期：2024-05-13

基金项目：“十四五”重大项目（2021QNLM020001）; 国家创新群体项目（41821002）; 国家自然科学基金项目（42074133）; 山东能源集团深层高温地热重大科技项目（SNKJ2022A06-R23）; 长庆油田公司复杂油田地震波成像技术的基础理论研究项目（2023-10502）; 中石油BGP项目（2023DJ8307）

第一作者：卫二祥（1992-），男，博士研究生，研究方向为地震波正演及偏移成像。E-mail： 976924258@qq.com。

通信作者：李振春（1963-），男，教授，博士，博士生导师，研究方向为地震波传播与成像。E-mail： leonli@upc.edu.cn。

文章编号：1673-5005（2025）01-0072-09""" doi：10.3969/j.issn.1673-5005.2025.01.007

逆时偏移是一种基于双程波方程进行波场延拓的高精度成像方法，相比于Krichhoff偏移、高斯束偏移以及单程波偏移等方法，具有不受地层倾角限制，对高陡构造准确成像的优势［1-2］。传统的偏移成像方法大多基于水平地表假设，对起伏地形的适应性不强，容易导致复杂地表条件下的构造成像质量不高［3-4］。在起伏地表条件下，地震波场的准确模拟是逆时偏移成像的基础。为此，国内外研究人员对起伏地表数值模拟算法进行了大量的研究和改进。一些学者在笛卡尔坐标系中采用矩形网格离散计算区域［5-7］，即使在起伏地表附近采用细网格进行离散，仍会出现阶梯状网格，无法从根本上消除虚假散射。 Tessmer等［8］、 Hestholm等［9］、董良国等［10］、Tarrass等［11］、Qu等［12］采用垂向网格坐标变换法将起伏边界转换为水平边界，然后在变换后的矩阵网格中进行数值模拟，基本消除了阶梯近似产生的虚假散射，但是这种方法在地表起伏剧烈的情况可能会不稳定。贴体网格能够贴合起伏边界，对起伏地表有着很好的适应性。因此一些研究学者将其引入到起伏地表数值模拟中，并取得了很好的模拟效果。Komatitsch等［13］采用贴体网格实现了起伏地表伪谱法数值模拟。Zhang等［14］首先将贴体网格应用在有限差分中，采用同位网格实现了一阶速度-应力方程的正演模拟。此后，Appel等［15］、兰海强等［16］将贴体网格有限差分算法推广到二阶弹性波方程中，实现了稳定、显式的二阶精度数值模拟。然而，采用同位网格需要更小的网格步长才能达到交错网格的差分精度，导致了计算成本的增加。丘磊等［17］推导了曲线坐标系下的一阶速度-应力方程，并使用标准交错网格实现了正演模拟。但是在曲线坐标系中，方程的速度、应力分量不再满足交错分布，使用标准交错网格求解时需要复杂的插值计算，不仅增加了计算量，还降低了模拟精度。为此，李庆洋等［3］、杨宇等［18］将Lebedev网格引入到贴体网格有限差分中，避免了标准交错网格的插值误差，降低了算法的复杂度，实现了起伏地表的高精度数值模拟。逆时偏移由于波动方程双向传播的特点，地震波在传播过程中遇到波阻抗界面就会发生反射，在成像过程中该反射波与正常传播的波互相关而产生低频噪声，严重影响成像质量和进一步的处理和解释精度［19-20］。早期的低频噪声的压制方法主要是在使用波动方程进行波场延拓时尽量避免发生背向反射，但是这类方法无法从根本上压制偏移噪声。为此，一些学者基于低频噪声的产生原理，发展了基于波场分离的成像方法。Liu等［21］基于Hilbert变换实现了上下行波波场分离，并且只保留方向相反的波场分量进行互相关成像。Fei等［22］提出只保留震源波场下行波和检波点波场上行波的互相关成像条件，避免了在速度剧烈变化的区域产生偏移假象。王一博等［23］进一步将波场分离为左上、左下、右上、右下4个组分，得到了16种成像结果，有效地压制了低频噪声和串扰噪声。以上方法都是基于Hilbert变换，在波场分离过程中需要进行大量Fourier变换，因而需要较高的计算成本和储存成本。而基于Poynting矢量的波场分离方法凭借其算法简单、计算量小的优势被广泛应用于逆时偏移成像。Yoon等［24］首次利用Poynting矢量确定波场的传播方向，实现了波场分离互相关成像。此后，很多学者开始基于Poynting矢量进行波场分解，针对不同的地质体选择合适的波场组分进行互相关成像［25-26］，取得了很好的成像效果。郭旭等［27］进一步将该方法推广到各向异性VTI介质逆时偏移成像中，有效地压制了低频噪声。王晓毅等［20］在最小二乘逆时偏移中使用Poynting矢量进行波场分解，压制了梯度计算时的低频噪声，提高了算法的收敛速度。笔者首先基于贴体网格和Lebedev网格实现起伏地表有限差分算法；然后利用Poynting矢量确定地震波场的传播方向，对震源波场和检波点波场进行方向行波分离，并推导相应的互相关成像条件；最终实现基于行波分离的起伏地表逆时偏移成像算法。

1" 方法原理

1.1" 曲线坐标系下的波动方程

当起伏地表存在时，数值模拟所用的离散网格需要很好地贴合起伏边界来避免阶梯离散产生的虚假散射，这样的网格就是贴体网格，网格内部的光滑性和边界的正交性保证了数值模拟的稳定和精度［14］。贴体网格生成的本质是通过坐标变换将规则的计算域转换到不规则的物理域（图1）。

贴体网格生成以后，便可以获得曲线坐标系中各个网格点（ξ，η） 与笛卡尔坐标系中的网格点（x，y） 的一一 对应关系：

x=x（ξ，η），z=z（ξ，η）.（1）

在笛卡尔坐标系中，各向同性介质中一阶速度-应力声波方程可以表示为

vxt=1ρpx ，

vzt=1ρpz ，

pt=ρV2pvxx+vzz.（2）

式中，vx、vz分别为质点的水平和垂直方向的振动速度；p为声压；ρ 为介质密度；Vp为声波速度。

domain and grid in computational domain

对方程（2）采用链式法则，可推导出曲线坐标系下的声波波动方程：

vxt=1ρ（ξxpξ+ηxpη），

vzt=1ρ（ξzpξ+ηzpη），

pt=ρV2pξxvxξ+ξzvzξ+ηxvxη+ηzvzη. （3）

1.2" Lebedev网格有限差分方法

由方程（3）可以看出，方程中的每个变量都在同一个网格的水平和垂直方向上存在空间偏导数，此时标准交错网格不再满足变量交错分布的特点，求解不在交错网格上的变量偏导数需要复杂的插值运算，增加了计算量，降低了模拟精度。为此，选用Lebedev网格来求解方程（3）。Lebedev网格的核心思想是速度和声压的不同分量交错定义在同一网格上（图2）。

从图2中的变量分布可看出，每个变量都定义在同一网格的不同位置，因此需要分别进行更新计算。以vx 为例，其时间2阶、2M阶空间差分精度的离散形式为

vxn+12i+12，j=vxn-12i+12，j+Δtρi+12，jξxi+12，jFξpni+12，j+ηxi+12，jFηpni+12，j，

Fξpni+12，j=∑Lm=1cmΔξpni+12+2m-12，j-pni+12-2m-12，j，

Fηpni+12，j=∑Lm=1cmΔηpni+12，j+2m-12-pni+12，j-2m-12. （4）

vxn+12i，j+12=vxn-12i，j+12+Δtρi，j+12ξxi，j+12Fξpni，j+12+ηxi，j+12Fηpni，j+12，

Fξpni，j+12=∑Lm=1cmΔξpni+2m-12，j+12-pni-2m-12，j+12，

Fηpni，j+12=∑Lm=1cmΔηpni，j+12+2m-12-pni，j+12-2m-12. （5）

式中，Δt 为时间采样间隔；Δξ、Δη 分别为水平和垂直方向的网格间距；cm 为有限差分系数。方程（4）和（5）用来分别更新不同位置的vx 变量，同理可以得到其他变量的离散形式。

1.3" Poynting矢量波场分离

Poynting矢量的概念最早出现在电磁学领域，由Yoon等［24］引入到地震勘探领域，并给出Poynting矢量的计算表达式为

P=-ptp. （6）

式中，P为Poynting矢量；p为声压场， 为梯度算子。由于用方程（3）进行逆时偏移成像过程中时，p/t，vx/t，vz/t都已经计算过，所以求取Poynting矢量时，基本不会额外增加计算量。现在通过式（6）来计算地震波的传播方向，用一个均匀介质模型来验证波场分离效果。模型速度为3000 m/s， 模型大小为301×301，网格间距均为10 m， 时间采样间隔为1 ms， 采用主频为20 Hz的雷克子波作为震源，震源位置为模型网格坐标（151，151） 处。图3展示了400 ms 时刻的波场快照和Poynting矢量，图3（b）和（c）分别是Poynting矢量的水平分量和垂直分量，图中黑色为正值，白色为负值，代表了地震波的传播方向。波场分离得到的上、下、左、右行波快照如图4所示，可以看到，根据Poynting矢量的方向实现了行波分离。

1.4" 互相关成像条件

常规逆时偏移采用的是震源波场和检波点波场，即不进行波场分离的全波场互相关成像条件［28］，其表达式为

I（x，z）=∫Tmax0S（x，z，t）R（x，z，t）dt.（7）

式中，I（x，z） 为成像结果；t 为时间；Tmax 为波场最大传播时间；S（x，z，t）、R（x，z，t） 分别为正向传播的震源波场和反向传播的检波点波场。

在逆时偏移过程中，震源波场正向传播和检波点波场反向延拓都包含了上行波和下行波，分别表示为

S（x，z，t）=Sd（x，z，t）+Su（x，z，t），R（x，z，t）=Rd（x，z，t）+Ru（x，z，t）. （8）

式中，下标d，u分别代表下行波和上行波。

将式（8）代入式（7）中，可以得到波场分离的互相关成像条件：

I（x，z）=∫Tmax0Sd（x，z，t）Ru（x，z，t）dt+∫Tmax0Su（x，z，t）Rd（x，z，t）dt+∫Tmax0Sd（x，z，t）Rd（x，z，t）dt+∫Tmax0Su（x，z，t）Ru（x，z，t）dt=

Idu（x，z）+Iud（x，z）+Idd（x，z）+Iuu（x，z）.（9）

其中，Liu等指出，Idd、Iuu 表示的是同向波场互相关成像，是低频噪声产生的主要原因，不应该包含在成像条件中。Fei等［22］进一步指出，Iud 虽然包含有效信息，但在偏移速度变化剧烈的区域容易产生偏移假象干扰成像质量，也应该从互相关成像条件中去除，最终可得到波场分离的互相关成像条件：

I（x，z）=Idu（x，z）=∫Tmax0Sd（x，z，t）Ru（x，z，t）dt.（10）

式（10）表示震源波场的下行波和检波点波场的上行波互相关成像条件。

与互相关成像条件相比，震源归一化成像条件可以获得更好的成像振幅［29］，因此采用震源归一化互相关成像条件：

I（x，z）=∫Tmax0Sd（x，z，t）Ru（x，z，t）dt∫Tmax0S（x，z，t）S（x，z，t）dt . （11）

2" 模型测试

2.1" 三层模型

为了直观展示波场分离的互相关成像结果，验证该算法的正确性，选用一个三层模型来进行测试，如图5所示，图5（a）和（b）分别为物理域和计算域速度模型。模型的起伏边界由高斯函数y（x）=100exp（-（x-81）2/400）-100exp（-（x-221）2/400）+120

生成，模型大小为301×201，水平和垂直方向的网格间距均为10 m。观测系统为全速度场接收，共301道接收，道间距为10 m，共激发1炮，炮点位置在网格（151，1）处，震源采用主频为20 Hz的雷克子波，时间采样间隔为1 ms， 总的记录时长为1.5 s。基于该模型生成用于逆时偏移的炮记录，该过程是使用方程（3）进行有限差分正演模拟得到。

正演模拟得到的波场快照与炮记录如图6和7所示。从图6可以看到，波场中没有出现虚假散射和数值频散，说明本文正演算法在处理起伏地表问题的优势；而从图7可以发现，受起伏地表的影响，反射波同相轴发生了扭曲，在后续的逆时偏移过程中需要消除这种影响才可以得到准确的成像结果。

图8是不同成像条件得到的偏移结果，图8（a）和（b）分别为全波场成像结果、Laplace滤波后的全波场成像结果，可以看出滤波后能压制大部分低频噪声，但是在浅层噪声并未被完全压制干净，而且在第一个层界面附近出现了偏移假象。图8（c）～（e）分别是震源上行波与检波点上行波成像结果、震源下行波与检波点下行波成像结果、震源上行波与检波点下行波成像结果，可以看出这些结果几乎不包含有效信息，是低频噪声产生的原因，在成像时应该舍去。图8（f）是本文所用的震源下行波与检波点上行波互相关成像条件得到的成像结果，对比图8（a）、（b）和（f）可得，提出的成像算法可以消除起伏地表的影响，有效压制低频噪声，去除偏移假象，得到准确的偏移成像结果。

2.2" 修改的SEG/EAGE盐丘模型

为了进一步验证提出的基于行波分离的起伏地表逆时偏移算法对复杂模型的适用性，采用修改的SEG/EAGE盐丘模型进行测试，如图9所示，图9（a）和（b）分别为物理域和计算域速度模型。该盐丘模型存在速度变化剧烈的区域和高陡的断层构造，其起伏边界由正弦函数y（x）=50sin（2π/150x）+60

生成，模型网格大小为676×225，水平和垂直方向的网格间距均为

10 m。采用中间激发，两边接收的观测系统，90炮均匀地布设在地表，炮检距为50 m，第1炮位于网格坐标（113 ， 1）处。每炮有225个检波点接收，道间距为10 m。震源为20 Hz的雷克子波，时间采样间隔为0.8 ms，总的采样时长为2.0 s， 总的时间采样点为2500个。分别采用基于矩形网格、贴体网格的全波场互相关成像算法以及本文行波分离成像算法对合成地震数据进行偏移成像测试，得到的成像结果如图10所示，图10（a）和（b）分别为基于矩形网格和贴

体网格的全波场成像Laplace滤波后的结果，图10（c）为本文行波分离成像结果。由图10（a）和（b）可以看出，受起伏地表的影响，在基于矩形网格的成像结果中，浅层存在散射噪音，地层同相轴发生了一定程度的扭曲，尤其是2 km深处的水平地层最为明显；而基于贴体网格的成像结果消除了这些影响。对比图10（b）和（c）可以发现，全波场成像算法得到的剖面经过Laplace滤波后，在浅层仍然存在低频噪声，而本文行波分离的成像结果中，偏移噪声得到了有效压制，成像剖面整体上更加清晰，成像质量更高。为了更清晰地展示本文算法的优势，将图10（b）和（c）红框中的区域进行放大显示，放大结果如图10（d）和（e）。通过对比可以看出，全波场成像结果在盐丘顶部位置存在偏移假象，如图10（d）中黑色箭头所指位置，会给后续的解释工作带来困难。而本文行波分离成像方法所得结果（图10（e））在对应位置有效地消除了偏移假象，得到了准确的成像结果。

3" 结束语

基于贴体网格和Lebedev网格的有限差分算法消除了阶梯离散造成的虚假散射，避免了标准交错网格的插值误差，可用于起伏地表下的正演模拟及逆时偏移成像。同时在逆时偏移的震源波场正向传播和检波点波场逆时延拓过程中，利用Poynting矢量对波场实现了行波分离，并选择有效波场成分进行互相关成像，从而实现了基于行波分离的起伏地表逆时偏移成像。数值模拟测试表明，与全波场互相关成像方法相比，提出的成像方法在基本不会额外增加计算量的前提下，不仅能够压制逆时偏移成像中的低频噪声，而且可以消除偏移假象，得到更高分辨率、更准确的成像结果。

参考文献：

［1］" BAYSAL E， KOSLOFF D D， SHERWOOD J W C." Reverse time migration［J］. Geophysics， 1983，48（11）：1514-1524.

［2］" MCMECHAN G A. Migration by extrapolation of time-dependent boundary values［J］. Geophysical Prospecting， 1983，31（3）：413-420.

［3］" 李庆洋，黄建平，李振春，等.起伏地表贴体全交错网格仿真型有限差分正演模拟［J］.石油地球物理勘探，2015，50（4）：633-642.

LI Qingyang， HUANG Jianping， LI Zhenchun， et al. Undulating surface body-fitted grid seismic modeling based on fully staggered-grid mimetic finite difference［J］. Oil Geophysical Prospecting， 2015，50（4）：633-642.

［4］" 孙建国.复杂地表条件下地球物理场数值模拟方法评述［J］.世界地质，2007，26（3）：345-362.

SUN Jianguo. Methods for numerical modeling of geophysical fields under complex topographical conditions： a critical review［J］. Global Geology， 2007，26（3）：345-362.

［5］" ROBERTSSON J O A. A numerical free-surface condition for elastic/viscoelastic finite-difference modeling in the presence of topography［J］. Geophysics， 1996，61（6）：1921-1934.

［6］" OHMINATO T， CHOUET B A. A free-surface boundary condition for including 3D topography in the finite-difference method［J］. Bulletin of the Seismological Society of America， 1997，87（2）：494-515.

［7］" HAYASHI K， BURNS D R， TOKSZ M N. Discontinuous-grid finite-difference seismic modeling including surface topography［J］. Bulletin of the Seismological Society of America， 2001，91（6）：1750-1764.

［8］" TESSMER E， KOSLOF D， BEHLE A. Elastic wave propagation simulation in the presence of surface topography［J］. Geophysical Journal International， 1992，108：621-632.

［9］" HESTHOLM S， RUUD B. 2D finite-difference elastic wave modelling including surface topography1［J］. Geophysical Prospecting， 1994，42（5）：371-390.

［10］" 董良国，郭晓玲，吴晓丰，等.起伏地表弹性波传播有限差分法数值模拟［J］.天然气工业，2007，27（10）：38-41.

DONG Liangguo， GUO Xiaoling， WU Xiaofeng， et al. Finite difference numerical simulation for elastic wave propagation in rugged topography［J］. Natural Gas Industry， 2007，27（10）：38-41.

［11］" TARRASS I， GIRAUD L， THORE P. New curvilinear scheme for elastic wave propagation in presence of curved topography［J］. Geophysical Prospecting， 2011，59（5）：889-906.

［12］" QU Y， HUANG J， LI Z， et al. A hybrid grid method in an auxiliary coordinate system for irregular fluid-solid interface modelling［J］. Geophysical Journal International， 2017，208（3）：1540-1556.

［13］" KOMATITSCH D， COUTE F， MORA P. Tensorial formulation of the wave equation for modelling curved interfaces［J］. Geophysical Journal International， 1996，127（1）：156-168.

［14］" ZHANG W， CHEN X F. Traction image method for irregular free surface boundaries in finite difference seismic wave simulation［J］. Geophysical Journal International， 2006，167（1）：337-353.

［15］" APPEL D， PETERSSON N A. A stable finite difference method for the elastic wave equation on complex geometries with free surfaces［J］. Communications in Computational Physics，2009，4（1）：84-107.

［16］" 兰海强，刘佳，白志明.VTI介质起伏地表地震波场模拟［J］.地球物理学报，2011，54（8）：2072-2084.

LAN Haiqiang， LIU Jia， BAI Zhiming. Wave-field simulation in VTI media with irregular free surface［J］. Chinese Journal of Geophysics， 2011，54（8）：2072-2084.

［17］" 丘磊，田钢，石战结，等.起伏地表条件下有限差分地震波数值模拟：基于广义正交曲线坐标系［J］.浙江大学学报，2012，46（10）：1923-1931.

QIU Lei， TIAN Gang， SHI Zhanjie， et al. Finite-difference method for seismic wave numerical simulation in presence of topography—in generally orthogonal curvilinear coordinate system［J］. Journal of Zhejiang University， 2012，46（10）：1923-1931.

［18］" 黄建平，杨宇，李振春，等.基于M-PML边界的Lebedev网格起伏地表正演模拟方法及稳定性分析［J］.中国石油大学学报（自然科学版），2016，40（4）：47-56.

HUANG Jianping， YANG Yu， LI Zhenchun， et al. Lebedev grid finite-difference modeling for irregular free-surface and stability analysis based on M-PML boundary condition［J］. Journal of China University of Petroleum （Edition of Natural Science）， 2016，40（4）：47-56.

［19］" 杜启振，朱钇同，张明强，等.叠前逆时深度偏移低频噪声压制策略研究［J］.地球物理学报，2013，56（7）：2391-2401.

DU Qizhen， ZHU Yitong， ZHANG Mingqiang， et al. A study on the strategy of low wavenumber noise for prestack reverse-time depth migration［J］. Chinese Journal of Geophysics， 2013，56（7）：2391-2401.

［20］" 王晓毅，张江杰，许宏桥，等.基于Poynting矢量行波分离的最小二乘逆时偏移［J］.地球物理学报，2021，64（2）：645-655.

WANG Xiaoyi， ZHANG Jiangjie， XU Hongqiao， et al. Least-squares reverse time migration with wavefield decomposition based on the Poynting vector［J］. Chinese Journal of Geophysics， 2021，64（2）：645-655.

［21］" LIU F， ZHANG G， MORTON S A， et al. An effective imaging condition for reverse-time migration using wavefield decomposition［J］. Geophysics， 2011，76（1）：S29-S39.

［22］" FEI T W， LUO Y， YANG J， et al. Removing 1 images in reverse time migration： the concept of de-primary［J］. Geophysics， 2015，80（6）：S237-S244.

［23］" 王一博，郑忆康，薛清峰，等.基于Hilbert变换的全波场分离逆时偏移成像［J］.地球物理学报，2016，59（11）：4200-4211.

WANG Yibo， ZHENG Yikang， XUE Qingfeng， et al. Reverse time migration with Hilbert transform based full wavefield decomposition［J］. Chinese Journal of Geophysics， 2016，59（11）：4200-4211.

［24］" YOON K， MARFURT K J. Reverse-time migration using the Poynting vector［J］. Exploration Geophysics， 2006，37（1）：102-107.

［25］" CHEN T， HE B. A normalized wavefield separation cross-correlation imaging condition for reverse time migration based on Poynting vector［J］. Applied Geophysics， 2014，11（2）：158-166.

［26］" 韩文功，王江川，李振春，等.基于Poynting矢量和方向行波波场分离的逆时偏移成像条件［J］.地球物理学进展，2015，30（6）：2660-2665.

HAN Wengong， WANG Jiangchuan， LI Zhenchun， et al. Reverse-time migration imaging conditions based on Poynting vector and directional one-way wave field separating［J］. Progress in Geophysics， 2015，30（6）：2660-2665.

［27］" 郭旭，黄建平，李振春，等.基于行波分离的VTI介质逆时偏移［J］.物探与化探，2019，43（1）：100-109.

GUO Xu， HUANG Jianping， LI Zhenchun， et al. Reverse time migration in VTI media based on wavefield decomposition［J］. Geophysical and Geochemical Exploration， 2019，43（1）：100-109.

［28］" CLAERBOUT JON F. Toward a unified theory of reflector mapping［J］. Geophysics， 1971，36（3）：467.

［29］" 慕鑫茹，黄建平，黎国龙，等.黏声TTI介质中纯qP波逆时偏移成像方法［J］.中国石油大学学报（自然科学版），2023，47（2）：44-52.

MU Xinru， HUANG Jianping， LI Guolong， et al. Reverse time migration in viscoacoustic TTI media based on pure qP wave equation［J］. Journal of China University of Petroleum （Edition of Natural Science）， 2023，47（2）：44-52.

（编辑" 修荣荣）