基于波场传播方向的波动方程偏移有效去噪方法(英文)

本文研究分析了双程波波动方程偏移成像中广泛存在的三种主要噪声,特别是针对过去研究中没能很好解决的存在于高速盐丘悬垂边界附近的射线状噪声,提出了基于优化成像条件的有效去噪方法。射线状噪声主要来自于震源一侧波场的下行透射波分量和接收阵列一侧波场的上行散射波互相关成像。这一部分能量具有较强的互相关性,但并不携带真实的反射面信息。它广泛存在叠前偏移成像中,与信号在强度上同量级。多数情况下偏移成像中的相关噪声由方向性传播的波场能量产生。利用波场梯度得到的波场传播角度,可以分离出噪声对应的波场能量,并在成像条件中减去。采用这一方法可以有效地去除多种噪声,包括直达波噪声、散射波噪声和射线状噪声。该去噪方法不依赖波场外推算子,在需要时可以方便地运用到几乎所有的波动方程偏移中去。并且该去噪方法针对噪声的物理根源,对信号的损害很小。对去噪后的偏移成像结果额外地进行波数域滤波处理,可以进一步提高叠加图像的质量。这一去噪方法在超广角单程波偏移成像中取得良好效果,我们同时期待其在其他双程波波动方程偏移特别是逆时偏移(RTM)中的成功运用。     

Gabor‐frame‐based Gaussian packet migration

We present a Gaussian packet migration method based on Gabor frame decomposition and asymptotic propagation of Gaussian packets. A Gaussian packet has both Gaussian‐shaped time–frequency localization and space–direction localization. Its evolution can be obtained by ray tracing and dynamic ray tracing. In this paper, we first briefly review the concept of Gaussian packets. After discussing how initial parameters affect the shape of a Gaussian packet, we then propose two Gabor‐frame‐based Gaussian packet decomposition methods that can sparsely and accurately represent seismic data. One method is the dreamlet–Gaussian packet method. Dreamlets are physical wavelets defined on an observation plane and can represent seismic data efficiently in the local time–frequency space–wavenumber domain. After decomposition, dreamlet coefficients can be easily converted to the corresponding Gaussian packet coefficients. The other method is the Gabor‐frame Gaussian beam method. In this method, a local slant stack, which is widely used in Gaussian beam migration, is combined with the Gabor frame decomposition to obtain uniform sampled horizontal slowness for each local frequency. Based on these decomposition methods, we derive a poststack depth migration method through the summation of the backpropagated Gaussian packets and the application of the imaging condition. To demonstrate the Gaussian packet evolution and migration/imaging in complex models, we show several numerical examples. We first use the evolution of a single Gaussian packet in media with different complexities to show the accuracy of Gaussian packet propagation. Then we test the point source responses in smoothed varying velocity models to show the accuracy of Gaussian packet summation. Finally, using poststack synthetic data sets of a four‐layer model and the two‐dimensional SEG/EAGE model, we demonstrate the validity and accuracy of the migration method. Compared with the more accurate but more time‐consuming one‐way wave‐equation‐based migration, such as beamlet migration, the Gaussian packet method proposed in this paper can correctly image the major structures of the complex model, especially in subsalt areas, with much higher efficiency. This shows the application potential of Gaussian packet migration in complicated areas.     

Imaging artefacts of artificial diving waves in reverse time migration: cause analysis in the angle domain and an effective removal strategy

In areas with strong velocity gradients, traditional reverse time migration based on cross-correlation imaging condition not only produces low-frequency noise but also generates diving wave artefacts. The artefacts caused by diving waves have no typical low-frequency characteristics and cannot be eliminated by simple high-pass filtering approaches. We apply the wave-field decomposition imaging condition to analyse the causes of false images in reverse time migration by decomposing the full wave-field into up-going and down-going components in the angle domain. We find that artificial diving wave imaging artefacts, which are generated by the cross-correlation between the up-going source and down-going receiver wave-fields in areas with strong velocity gradients, arise at large angles. We propose an efficient strategy by means of the wavelength-dependent smoothing operator to eliminate artefacts from artificial diving waves in reverse time migration. Specifically, the proposed method provides more reasonable down-going wave-fields in areas with sharp velocity constructs by considering the factor of varying seismic wavelengths during wave propagation, and the artificial components of diving waves are eliminated in a straightforward manner. Meanwhile, the other wave-field components that contribute to true subsurface images are minimally affected. Benefiting from a smoothed velocity, the proposed method can be adapted to the traditional reverse time migration imaging frame, which reveals significant implementation potential for the seismic exploration industry. A salt model is designed and included to demonstrate the effectiveness of our approach.     

基于波场分离的弹性波逆时偏移

尽管叠前逆时偏移成像精度高,但仅针对单一纵波的成像也可能形成地下介质成像盲区,由于基于弹性波方程的逆时偏移成像可形成多波模式的成像数据,因此弹性波逆时偏移成像可提供更为丰富的地下构造信息.本文依据各向同性介质的一阶速度-应力方程组构建震源和检波点矢量波场,再利用Helmholtz分解提取纯纵波和纯横波波场,使用震源归一化的互相关成像条件获得纯波成像,避免了直接使用坐标分量成像而引起的纵横波串扰问题.针对转换波成像的极性反转问题,文中提出一种共炮域极性校正方法.为有效节约存储成本,也提出一种适用于弹性波逆时偏移的震源波场逆时重建方法,在震源波场正传过程中,仅保存PML边界内若干层的速度分量波场,进而逆时重建出所有分量的震源波场.本文分别对地堑模型和Marmousi2模型进行了弹性波逆时偏移成像测试,结果表明:所提出的共炮域极性校正方法正确有效,基于波场分离的弹性波逆时偏移成像的纯波数据能够对复杂地下构造准确成像.     

多波多分量高斯束叠前深度偏移

本文对基于弹性波动理论的多波多分量高斯束偏移进行了完整且详细的分析和公式推导,实现了3D空间多分量(矢量)波场的直接成像.由于当前多数基于弹性波动方程的偏移方法只是假设应力边界条件为自由地表边界条件,这种假设不符合垂直地震剖面(VSP)和海底电缆(OBC)等地震数据.为此本文详细分析了实际应用中常见的三种弹性各向同性介质模型的应力边界条件:自由空间、海底和自由地表模型.在上行传播假设情况下,获得了应力边界条件与位移边界条件的关系式.在此基础上,准确推导了3D多波多分量高斯束波场延拓和偏移成像公式,并在偏移过程中实现了完整的多波型自动分离.由于常规的互相关成像条件不适用于矢量波场成像,本文引用了散度/旋度互相关成像条件.通过约定PS转换波的正向旋转方向解决了3D空间PS成像极性翻转问题.利用2D和3D模型数据偏移成像验证了我们所提出的多波多分量高斯束偏移方法的可行性.     

特征高斯波包叠前深度偏移方法

高斯波包（Gaussian packet）传播算子可在局部时空域高效地计算局部波包的传播.高斯波包叠前深度偏移的基础是在Gabor变换域描述观测数据,再利用高斯波包传播算子计算炮点波场和检波点波场,两者相关即可得到偏移结果.利用炮道集的局部τ-p特征在Gabor变换域表达观测数据,可以仅关注部分高斯波包框架函数上的数据投影,这样既实现了波场的压缩存储,同时可利用高斯波包传播算子反传框架函数以实现整个炮道集的快速反传.这些综合了观测数据局部τ-p特征的高斯波包函数称为特征高斯波包（characteristic Gaussian packet,CGP）,相应的波场反传称为特征高斯波包反传.理论及数值分析证明了上述特征高斯波包反传方法是有效且快速的.炮点正传波场也利用高斯波包传播算子模拟.利用互相关成像条件可实现特征高斯波包叠前深度偏移（characteristic Gaussian packet pre-stack depth migration,CGPM）.由于高斯波包传播算子描述了局部方向及局部空间的波场,所以CGPM可以自然地提取角度域成像道集（ADCIG）,并易于实现面向目标叠前深度偏移,从而作为偏移引擎为偏移速度分析（MVA）服务.数值实验证明了CGPM和面向目标CGPM的有效性和实用性.     

基于解析时间波场外推与波场分解的逆时偏移方法研究

双程波方程逆时深度偏移是复杂介质高精度成像的有效技术,但其结果中通常包含成像方法引起的噪音和假象,一般的滤波方法会破坏成像剖面上的振幅,其中的假象也会给后续地质解释带来困扰.将波场进行方向分解然后实现入射波与反射波的相关成像能够有效地消除这类成像噪音,并提高逆时偏移成像质量.波传播方向的分解通常在频率波数域实现,它会占用大量的存储和计算资源,不便于在沿时间外推的逆时深度偏移中应用.本文提出解析时间波场外推方法,可以在时间外推的每个时间片上实现波传播方向的显式分解,逆时深度偏移中利用分解后的炮检波场进行对应的相关运算,实现成像噪音和成像信号的分离.在模型和实际数据上的测试表明,相比于常规互相关逆时偏移成像结果,本文方法能够有效地消除低频成像噪音和特殊地质构造导致的成像假象.     

基于Poynting矢量行波分离的最小二乘逆时偏移

因为在逆时偏移中基于双程波动方程构建震源波场和检波器波场,所以在波场延拓过程中地震波遇到波阻抗界面时,背向发育的反射波会与正常传播的波场互相关产生较强振幅的低频噪声.这一特点使得以逆时偏移为基础的最小二乘偏移方法在梯度计算时同样存在着低频噪声的干扰,从而导致反演收敛的速度减慢.考虑到计算量和存储成本的因素,本文借助Po...     

叠前地震数据特征波场分解、偏移成像与层析反演

本文提出了一套叠前地震数据稀疏表达(特征波场合成)、深度偏移成像和层析成像的处理流程.不同于传统的变换域中的数据稀疏表达理论,本文利用局部平面波的传播方向(慢度矢量),在中心炮检点处同时进行波束合成,从而将地震数据投影到局部平面波域(高维空间)中.由于波束合成后的地震数据描述了局部平面波的方向特征,因此称之为特征波场.然而波束合成算法需要估计局部平面波的慢度矢量.当地震数据受噪声干扰时,难以在常规τ-p谱中自动估计局部平面波的射线参数(慢度矢量).本文提出了基于反演理论的特征波场合成方法,可以同时反演局部平面波及其传播方向,从而提高特征波合成的自动化程度并保持方法的稳健性.通过特征波场合成,可以将地震数据分解为单独的震相(波形).这样的数据可以直接用来成像及反演.在局部平面波域中,由于局部平面波的入射与出射射线参数已知,传统的Kirchhoff叠前深度偏移(PSDM)和高斯束/控制束PSDM可以实现从"沿等时面的画弧"到"向反射点(段)的直接投影"的转变,叠前偏移的效率以及成像质量可以同时提高.此外,特征波场与地下反射点(段)的一对一映射关系使得叠前深度偏移与层析成像融为一体,可以极大地提高速度反演的效率.数值试验证明了特征波场合成、叠前深度成像以及层析反演的有效性.     

隧道工程勘察中的超前预报成像方法研究

针对地下工程领域隧道超前预报地震波波场传播与成像中存在的问题,通过数值模拟,构建二维含低速异常的隧道介质模型,研究隧道弹性波场传播规律和异常体边界成像准确性.首先,利用一阶速度-应力波动方程和高阶交错网格有限差分计算方法,导出隧道超前预报数值模拟的稳定性条件和边界条件,对上述隧道模型进行数值模拟,识别波场特征;其次,利用叠前逆时偏移成像方法,对压制噪音干扰后的波场在互相关成像条件下,对隧道模型中的异常体边界进行逆时偏移成像.研究结果表明：采用高阶交错网格有限差分正演获得异常体边界清晰的反射波和角点产生的散射波;逆时偏移算法获得隧道内异常体准确成像结果,从而大大提高隧道超前预报的分辨率与准确性;靠近掌子面单一震源、多道接收观测系统对异常体成像效果最佳,为隧道内高效数据采集提供理论依据.     

Converted‐wave beam migration with sparse sources or receivers

Gaussian beam depth migration overcomes the single‐wavefront limitation of most implementations of Kirchhoff migration and provides a cost‐effective alternative to full‐wavefield imaging methods such as reverse‐time migration. Common‐offset beam migration was originally derived to exploit symmetries available in marine towed‐streamer acquisition. However, sparse acquisition geometries, such as cross‐spread and ocean bottom, do not easily accommodate requirements for common‐offset, common‐azimuth (or common‐offset‐vector) migration. Seismic data interpolation or regularization can be used to mitigate this problem by forming well‐populated common‐offset‐vector volumes. This procedure is computationally intensive and can, in the case of converted‐wave imaging with sparse receivers, compromise the final image resolution. As an alternative, we introduce a common‐shot (or common‐receiver) beam migration implementation, which allows migration of datasets rich in azimuth, without any regularization pre‐processing required. Using analytic, synthetic, and field data examples, we demonstrate that converted‐wave imaging of ocean‐bottom‐node data benefits from this formulation, particularly in the shallow subsurface where regularization for common‐offset‐vector migration is both necessary and difficult.     

基于时空局域化dreamlet单程波算子的观测系统沉降法偏移

Dreamlet偏移成像目的是探索一类能够对地震波场和单程波传播算子同时分解和压缩的理论和方法,也即实现在压缩域的传播与成像、地震数据在时间和空间的非平稳性质,决定了要实现地震数据的有效稀疏表示,分解方法必须在时间和空间上同时具有局域化性质.Dreamlet由时间和空间局部分解原子的张量积构成,可以看作一种脉冲-小波束形式的波场分解原子.时空局域化的dreamlet单程波传播算子在对波场沿深度方向延拓时,地震数据在时间轴上总是向同一方向流动.随着深度的增加,部分用于成像浅层结构的数据归位至其空间位置后被dreamlet算子丢弃,波场的有效记录时间变短,每一步用于波场延拓的计算量也相应下降.为了充分发挥这一优势,本文介绍dreamlet传播算子与观测系统沉降法偏移相结合的理论与方法.观测系统沉降法偏移每一步都将记录到的所有数据向下延拓,沉降后的波场等效于把源和检波器都放置在目标深度所能接收的反射数据.Dreamlet观测系统沉降过程只保留用于成像观测系统下部地质结构的有效数据,自动丢弃已经用于成像观测系统上部而对下部成像没有贡献的信号.本文通过二维SEG/EAGE叠后和Marmousi叠前数据算例展示了dreamlet传播算子应用于观测系统沉降法偏移的这一特点.数值算例结果显示,在不影响成像质量的前提下,该偏移方法能够有效减少传播数据量,为发展一种快速高效的偏移方法提供了新的思路.     

Directional‐oriented wavefield imaging: a new wave‐based subsurface illumination imaging condition for reverse time migration

The key objective of an imaging algorithm is to produce accurate and high‐resolution images of the subsurface geology. However, significant wavefield distortions occur due to wave propagation through complex structures and irregular acquisition geometries causing uneven wavefield illumination at the target. Therefore, conventional imaging conditions are unable to correctly compensate for variable illumination effects. We propose a generalised wave‐based imaging condition, which incorporates a weighting function based on energy illumination at each subsurface reflection and azimuth angles. Our proposed imaging kernel, named as the directional‐oriented wavefield imaging, compensates for illumination effects produced by possible surface obstructions during acquisition, sparse geometries employed in the field, and complex velocity models. An integral part of the directional‐oriented wavefield imaging condition is a methodology for applying down‐going/up‐going wavefield decomposition to both source and receiver extrapolated wavefields. This type of wavefield decomposition eliminates low‐frequency artefacts and scattering noise caused by the two‐way wave equation and can facilitate the robust estimation for energy fluxes of wavefields required for the seismic illumination analysis. Then, based on the estimation of the respective wavefield propagation vectors and associated directions, we evaluate the illumination energy for each subsurface location as a function of image depth point and subsurface azimuth and reflection angles. Thus, the final directional‐oriented wavefield imaging kernel is a cross‐correlation of the decomposed source and receiver wavefields weighted by the illuminated energy estimated at each depth location. The application of the directional‐oriented wavefield imaging condition can be employed during the generation of both depth‐stacked images and azimuth–reflection angle‐domain common image gathers. Numerical examples using synthetic and real data demonstrate that the new imaging condition can properly image complex wave paths and produce high‐fidelity depth sections.     

矢量黏弹性衰减补偿高斯束偏移

基于弹性波动理论的多波多分量高斯束偏移具有计算效率高和成像准确等优点.但是目前此方法没有考虑实际地下介质的黏弹性对地震波传播的影响,从而无法补偿能量衰减和校正相位畸变,这使得该方法对一些含高黏弹性地层的成像效果不佳.针对衰减区域的成像问题,本文提出一种黏弹性衰减补偿高斯束偏移方法,该方法以多波多分量矢量波场弹性高斯束偏移方法为基础,在偏移过程中沿射线路径通过引入品质因子Q来考虑黏弹性影响并进行衰减补偿.该方法能够在偏移过程中实现PP波和PS波的自动分离及分别成像.同时,本文给出了在矢量波场偏移过程中提取角度域共成像点道集的方法,以便用于成像质量控制,并为后续速度和黏弹性参数反演提供所需的数据.本文利用2D层状模型和洼陷模型进行了方法测试,其成像结果验证了本文所提出的黏弹性衰减补偿高斯束偏移方法的可行性和有效性.     

True amplitude elastic Gaussian beam imaging of multicomponent walkaway vertical seismic profiling data

We develop the true‐amplitude prestack migration of multicomponent data based on the use of elastic Gaussian beams for walkaway vertical seismic profile (VSP) acquisition systems. It consists in a weighted summation of multishot data with specific weights, computed by tracing elastic Gaussian beams from each imaging point of the target area towards the sources and receivers. Each pair of beams may be connected with either a pair of P‐rays (PP‐image) or the P‐ray towards sources and the S‐ray to receivers (PS‐image) and is uniquely determined by dip (the angle of the bisector between the rays and the vertical direction) and opening (the angle between the rays) angles. Shooting from the bottom towards the acquisition system helps to avoid well‐known troubles, in particular multipathing for the imaging conditions in complex velocity models. The ability to fix the dip angle and implement summation over opening angles leads to the so‐called selective images that contain mostly interfaces with desired slopes. On the other hand, a set of images computed for a range of opening angles by summation over all available dip angles is used as input of an AVO‐like inversion procedure for the recovery of elastic parameters. The feasibility of this imaging procedure is verified by synthetic data for 2D realistic elastic models.     

基于激发振幅成像条件的探地雷达逆时偏移成像

针对基于互相关成像条件的探地雷达（GPR）逆时偏移计算效率低、存储量大及易产生低频假象的不足,本文将激发振幅成像条件应用于GPR逆时偏移成像中.通过在源点电磁波场正向传播过程计算每个网格点的能量密度,并保存最大能量密度的时刻和相应的电磁波场值;在接收点电磁波场逆向传播过程提取每个网格点最大能量密度时刻及对应的电磁波场值,并利用保存的最大能量源点电磁波场及走时做归一化,从而获得了依赖反射系数成像剖面,避免了源点正向传播电磁波场的存储和重建.此外,为了提高电磁波场的模拟精度,采用了基于三角形剖分的时间域有限单元法（FETD）计算电磁波正向和逆向传播过程.最后通过模型试算表明：激发振幅成像条件相比于归一化互相关成像条件,成像结果低频噪声更弱,空间分辨率更高,计算效率提高了近2倍.     

Global-scale seismic interferometry: theory and numerical examples

Progress in the imaging of the mantle and core is partially limited by the sparse distribution of natural sources; the earthquake hypocenters are mainly along the active lithospheric plate boundaries. This problem can be approached with seismic interferometry. In recent years, there has been considerable progress in the development of seismic interferometric techniques. The term seismic interferometry refers to the principle of generating new seismic responses by cross‐correlating seismic observations at different receiver locations. The application of interferometric techniques on a global scale could create sources at locations where no earthquakes occur. In this way, yet unknown responses would become available for the application of travel‐time tomography and surface‐wave dispersion studies. The retrieval of a dense‐enough sampling of source gathers would largely benefit the application of reflection imaging. We derive new elastodynamic representation integrals for global‐scale seismic interferometry. The relations are different from other seismic interferometry relations for transient sources, in the sense that they are suited for a rotating closed system like the Earth. We use a correlation of an observed response with a response to which free‐surface multiple elimination has been applied to account for the closed system. Despite the fact that the rotation of the Earth breaks source‐receiver reciprocity, the seismic interferometry relations are shown to be valid. The Coriolis force is included without the need to evaluate an extra term. We synthesize global‐scale earthquake responses and use them to illustrate the acoustic versions of the new interferometric relations. When the sampling of real source locations is dense enough, then both the responses with and without free‐surface multiples are retrieved. When we do not take into account the responses from the sources in the direct neighborhood of the seismic interferometry‐constructed source location, the response with free‐surface multiples can still be retrieved. Even when only responses from sources at a certain range of epicentral distances are available, some events in the Green's function between two receiver locations can still be retrieved. The retrieved responses are not perfect, but the artefacts can largely be ascribed to numerical errors. The reconstruction of internal events – the response as if there was a source and a receiver on (major) contrasts within the model – could possibly be of use for imaging. With modelling it is possible to discover in which region of the correlation panel stationary phases occur that contribute to the retrieval of events. This knowledge opens up a new way of filtering out undesired events and of discovering whether specific events could be retrieved with a given source‐receiver configuration.     

叠前逆时偏移影响因素分析

反射地震勘探中的偏移成像技术是获取地下介质构造形态最有效的手段之一.在叠前深度域偏移方法中,目前工业界采用的方法包括基于射线理论的波动方程积分解法和基于波动理论的微分波动方程单程波解法,这两类方法难以处理地震波横向速度变化剧烈的高陡倾角构造成像问题.近年来勘探地震学研究领域发展起来的叠前逆时偏移采用了双程波求解微分波动方程的算法,这种方法具有相位准确、不受介质横向速度变化和高陡倾角构造的影响、成像精度高、可以利用回转波正确成像等优点,从理论上弥补了当前工业界常规地震偏移所面临的成像缺陷.然而,叠前逆时偏移成像方法从理论走向实用尚需解决如下问题:计算速度和数据存储空间的节省、初始速度模型的建立、震源子波的选择、数值模型边界条件的定义和假像的消除等等.对于计算速度和存储量大的问题,随着计算机硬件的快速发展,将会不断得到改善,同时可以采取一些计算技术和存储策略来加以缓解.本文主要针对初始速度模型的建立、震源子波的选择、数值模型边界条件的定义和假像的消除这些因素,利用简单模型进行了分析.对于反射波造成的传播路径上的假像,给出了一种振幅补偿滤波方法.对勘探地球物理学界给出的SEG/EAGE二维盐丘模型、Marmousi模型和本研究设计的崎岖海底模型进行了叠前逆时偏移成像,均取得了较好的成像效果.     

利用三维高斯射线束成像进行地震定位

常规的地震定位方法通常需要拾取地震记录的初至,当初至不明显或被较高水平的噪声淹没时精度较低.本文采用基于三维高斯射线束的偏移成像方法对震源进行定位,较好地解决了该问题.通过三维高斯射线束对台站记录进行偏移归位,并将各台站成像结果的交点作为地震能量释放的中心位置;当各台站成像结果不能交于一点时,采用三维空间高斯滤波方法可实现震源位置的自动获取.提出的变网格计算方案极大地减少了计算量,显著地提高了成像精度和计算效率.利用首都圈地震台网数据,对涿鹿、滦县以及房山三个地震事件进行试算,结果表明:基于变网格三维高斯束偏移成像的地震定位方法自动化程度很高,而且具有较好的抗噪能力,特别适合处理低信噪比资料的地震定位问题.     

基于一阶弹性波方程的能量互相关成像条件

基于地震波场能量构建的能量互相关成像条件,具有易实现、物理意义明确及背向散射压制效果明显等优势.但是,目前构建的能量互相关成像条件仅适用于二阶弹性波方程,难以直接应用于一阶弹性波方程.为此,本文针对一阶弹性波方程,基于能量守恒定理及能量密度,构建以速度-应力为参数的能量范数以表征弹性波场能量,将速度-应力能量范数拓展为能量内积以提取弹性波场反射能量.震源端与检波端的基矢量正方向保持一致的基础上,构建得到可有效压制背向散射的弹性波能量成像条件.数值模拟结果表明:该成像条件可以得到背向散射压制、振幅有效保持的能量成像结果.