稀疏反演多次波去除策略与效果分析

多次波去除是海上勘探中的一项关键技术,在一定程度上决定着最终成像结果的质量.表面相关多次波去除方法(Surface-Related Multiple Elimination,SRME)是一种基于反馈迭代模型的数据驱动类方法,该方法利用地震数据本身预测多次波,进一步对其进行匹配相减获得一次波结果,这往往导致与多次波重叠的一次波被错误地去除.稀疏反演一次波估计方法(Estimation of Primaries by Sparsity Inversion,EPSI)在原理上同样基于反馈迭代模型,通过求解一次波和表面相关多次波总残差最小化问题直接获得一次波.该方法无多次波的匹配相减过程,可减少与多次波发生重叠的一次波产生损失.另外,EPSI可以很好地解决SRME无法解决的近道缺失问题.理论数据和实际资料处理对比表明SRME对数据完整性以及一次波和多次波的分布情况有较强的依赖性,而EPSI几乎不受这两个因素的影响.EPSI相比SRME方法计算成本高,尚未在工业界广泛应用,但随着高性能计算的普及,该方法的实用性有望得到突破.     

Using surface multiples to estimate primaries by sparse inversion from blended data

For data acquired with conventional acquisition techniques, surface multiples are usually considered as noise events that obscure the primaries. However, in this paper we demonstrate that for the situation of blended acquisition, meaning that different sources are shooting in a time‐overlapping fashion, multiples can be used to ‘deblend’ the seismic measurements. We utilize the recently introduced estimation of primaries by sparse inversion (EPSI) methodology, in which the primary impulse responses are considered to be the unknowns in a large‐scale inversion process. With some modifications the estimation of primaries by sparse inversion method can be used for blended seismic data. As output this process gives unblended primary impulse responses with point sources and receivers at the surface, which can be used directly in traditional imaging schemes. It turns out that extra information is needed to improve on the deblending of events that do not have much associated multiple energy in the data, such as steep events at large offsets. We demonstrate that this information can be brought in during acquisition and during processing. The methodology is illustrated on 2D synthetic data.     

Phase inversion deconvolution for long and short period multiples attenuationl1

A new approach to deconvolution has been developed to improve the attenuation of multiple energy. This approach to deconvolution is unique in that it not only eliminates the usual assumptions of a minimum phase lag wavelet and a random distribution of impulses, but also overcomes the noise limitation of the homomorphic deconvolution and its inherent instability to phase computation. We attempt to analyse the continuous alteration of the acoustic waveform during the propagation through a linear system. Based on the results of this analysis, the surface-related measurements are described as a convolution of the impulse response of the system with the non-stationary forward wavelet which includes all multiple terms generated within the system. The amplitude spectrum of the forward wavelet is recovered from the amplitude spectrum of the recorded signal, using the difference between the rate of decay of the source wavelet and the duration of the measurement. The phase spectrum of the forward wavelet is estimated using the Hilbert transform and the fact that the mixed phase lag wavelet can be presented as a convolution of the minimum and maximum phase lag wavelets. The multiples are discriminated from primaries by comparison of the phase spectrum of the seismic signal and the inverse of the forward wavelet. Therefore, the technique is called phase inversion deconvolution (PID). This approach requires no velocity information in order to recognize and attenuate multiple energy. Therefore, primary energy is recovered in the near-offset region where the velocity differential between primary and multiple energies is very small.     

逆子波域消除多次波方法研究

SRMA（与表面相关多次波的衰减）算法包含预测和相减两步.相减算法中,当多次波与反射波同相轴相交时,如何有效减去多次波、保留反射波,是面临的主要问题.通过分析非正交性对滤波器（逆子波）的影响,可以证明：逆子波因非正交性产生的误差呈近似的高斯分布.在此基础上,本文提出了在逆子波域（单道自适应相减滤波的滤波算子的集合）,利用其误差的概率分布特征,对逆子波进行估计,用逆子波的估计对逆子波进行校正来消除多次波的方法.其步骤为：首先用SRMA方法预测出表面多次波,并对每一单炮进行单道自适应相减,得到逆子波,形成逆子波域;其次,在逆子波域采用中值滤波,提取接近真实逆子波的逆子波估计;第三,在逆子波域用逆子波估计对畸变的逆子波进行校正;最后采用校正后的逆子波来衰减多次波.通过简单模型和SMARRT模型的测试,该方法不仅能够有效减去多次波,而且在相交的区域,能够保持反射波同相轴的连续性并恢复其正确的振幅.     

衰减干涉假象的海洋一次波与多次波同时成像方法

海水与空气间的强波阻抗差使得海洋地震资料普遍发育自由表面相关多次波,如何利用好多次波所携带的有效信息已成为提高海洋地震资料成像品质的新突破点.基于面炮偏移的一次波与多次波同时成像方法能够避免多次波预测精度的影响,但是,正向传播的震源子波与反向延拓的自由表面相关多次波所产生的干涉假象严重制约了该技术的应用,本文提出了一种基于单程波偏移算子,可在成像域压制干涉假象的一次波与多次波同时成像方法.其中包含了三个步骤:第一,传统单程波偏移成像方法中的震源子波替换为一次波、多次波与震源子波,初始上行延拓波场为一次波与多次波,基于单程波算子的波场延拓与互相关成像条件的应用得到包含干涉假象的一次波与多次波同时成像;第二,以子波为震源,自由表面相关多次波为记录,按照传统单程波偏移成像方法得到干涉假象;第三,基于最小二乘匹配滤波算法,将第一步的成像结果与第二步的干涉假象进行匹配相减,得到干涉假象衰减后的一次波与多次波同时成像,避开了由于实际资料子波无法准确提取而造成一次波与多次波对成像能量级的不一致性.Sigsbee2B模型测试验证了本方法的有效性,并在我国某探区深海实际资料处理中得到了成功应用,深层基底得到了清晰刻画,并且照明均衡度明显改善.     

基于保幅拉东变换的多次波衰减

为在去除多次波时有效保护地震一次反射波数据的AVO现象,给后续反演、解释提供准确的地震数据,本文提出了一种基于保幅拉东变换的多次波衰减方法,该方法是对常规抛物拉东变换的修改,把常规的稀疏拉东变换在拉东域分成两部分：一部分用于模拟零偏移距处的反射波能量,增加的另一部分用于模拟反射波振幅的AVO特性.该方法不仅考虑了反射波同相轴的形状,还考虑了反射波同相轴振幅幅度的变化,从而可把反射波信息进行有效转换,进而有利于多次波的消除,更好地恢复有效波的能量.在把地震数据由时间域转换到拉东域时,本文采用了IRLS算法实现保幅拉东算子的反演.模型数据和实际地震道集的试算分析表明,与常规拉东变换相比,保幅拉东变换在去除多次波的同时可有效保护一次反射波的AVO现象.     

VSP上下行反射波联合成像方法研究

VSP资料上下行波场发育丰富.本文在分析VSP直达波、上行反射波、下行反射波传播路径及其照明范围的基础上,指出了常规VSP波动方程偏移方法缺陷,进而通过修改波场延拓方式,提出了上下行反射波联合成像方法,并在高频近似下分析了该方法的成像原理.该方法不需要进行VSP上下行反射波场分离,能够同时对VSP资料中的一次反射波、自由表面多次波、层间多次波进行成像,比常规成像剖面具有更宽的成像范围和更好的成像效果.该方法能够对下行一次反射波进行成像,从而可以实现常规偏移方法难以处理的高陡倾角构造成像.模拟资料和实际资料处理证明了本文方法的正确性.     

Full‐wavefield migration: using surface and internal multiples in imaging

Multiple scattering is usually ignored in migration algorithms, although it is a genuine part of the physical reflection response. When properly included, multiples can add to the illumination of the subsurface, although their crosstalk effects are removed. Therefore, we introduce full‐wavefield migration. It includes all multiples and transmission effects in deriving an image via an inversion approach. Since it tries to minimize the misfit between modeled and observed data, it may be considered a full waveform inversion process. However, full‐wavefield migration involves a forward modelling process that uses the estimated seismic image (i.e., the reflectivities) to generate the modelled full wavefield response, whereas a smooth migration velocity model can be used to describe the propagation effects. This separation of modelling in terms of scattering and propagation is not easily achievable when finite‐difference or finite‐element modelling is used. By this separation, a more linear inversion problem is obtained. Moreover, during the forward modelling, the wavefields are computed separately in the incident and scattered directions, which allows the implementation of various imaging conditions, such as imaging reflectors from below, and avoids low‐frequency image artefacts, such as typically observed during reverse‐time migration. The full wavefield modelling process also has the flexibility to image directly the total data (i.e., primaries and multiples together) or the primaries and the multiples separately. Based on various numerical data examples for the 2D and 3D cases, the advantages of this methodology are demonstrated.     

表面多次波最小二乘逆时偏移成像

使用相同的炮记录,多次波偏移能提供比反射波偏移更广的地下照明和更多的地下覆盖但是同时产生很多的串声噪声.相比传统逆时偏移,最小二乘逆时偏移反演的反射波成像结果具有更高的分辨率和更均衡的振幅.我们主要利用最小二乘逆时偏移压制多次波偏移产生的串声噪声.多次波最小二乘逆时偏移通常需要一定的迭代次数以较好地消除串声噪声.若提前将一阶多次波从所有阶数的多次波中过滤出来,使用相同的迭代次数,一阶多次波的最小二乘逆时偏移能够得到具有更高信噪比的成像剖面,而且能够提供与多次波最小二乘逆时偏移相似的有效地下结构成像.     

Passive seismic data primary estimation and noise removal via focal‐denoising closed‐loop surface‐related multiple elimination based on 3D L1‐norm sparse inversion

Passive seismic has recently attracted a great deal of attention because non‐artificial source is used in subsurface imaging. The utilization of passive source is low cost compared with artificial‐source exploration. In general, constructing virtual shot gathers by using cross‐correlation is a preliminary step in passive seismic data processing, which provides the basis for applying conventional seismic processing methods. However, the subsurface structure is not uniformly illuminated by passive sources, which leads to that the ray path of passive seismic does not fit the hyperbolic hypothesis. Thereby, travel time is incorrect in the virtual shot gathers. Besides, the cross‐correlation results are contaminated by incoherent noise since the passive sources are always natural. Such noise is kinematically similar to seismic events and challenging to be attenuated, which will inevitably reduce the accuracy in the subsequent process. Although primary estimation for transient‐source seismic data has already been proposed, it is not feasible to noise‐source seismic data due to the incoherent noise. To overcome the above problems, we proposed to combine focal transform and local similarity into a highly integrated operator and then added it into the closed‐loop surface‐related multiple elimination based on the 3D L1‐norm sparse inversion framework. Results proved that the method was capable of reliably estimating noise‐free primaries and correcting travel time at far offsets for a foresaid virtual shot gathers in a simultaneous closed‐loop inversion manner.     

Marine seismic wavefield measurement to remove sea-surface multiples

We propose a new method for removing sea-surface multiples from marine seismic reflection data in which, in essence, the reflection response of the earth, referred to a plane just above the sea-floor, is computed as the ratio of the plane-wave components of the upgoing wave and the downgoing wave. Using source measurements of the wavefield made during data acquisition, three problems associated with earlier work are solved: (i) the method accommodates source arrays, rather than point sources; (ii) the incident field is removed without simultaneously removing part of the scattered field; and (iii) the minimum-energy criterion to find a wavelet is eliminated. Pressure measurements are made in a horizontal plane in the water. The source can be a conventional array of airguns, but must have both in-line and cross-line symmetry, and its wavefield must be measured and be repeatable from shot to shot. The problem is formulated for multiple shots in a two-dimensional configuration for each receiver, and for multiple receivers in a two-dimensional configuration for each shot. The scattered field is obtained from the measurements by subtracting the incident field, known from measurements at the source. The scattered field response to a single incident plane wave at a single receiver is obtained by transforming the common-receiver gather to the frequency–wavenumber domain, and a single component of this response is obtained by Fourier transforming over all receiver coordinates. Each scattered field component is separated into an upgoing wave and a downgoing wave using the zero-pressure condition at the water-surface. The upgoing wave may then be expressed as a reflection coefficient multiplied by the incident downgoing wave plus a sum of scattered downgoing plane waves, each multiplied by the corresponding reflection coefficient. Keeping the upgoing scattered wave fixed, and using all possible incident plane waves for a given frequency, yields a set of linear simultaneous equations for the reflection coefficients which are solved for each plane wave and for each frequency. To create the shot records that would have been measured if the sea-surface had been absent, each reflection coefficient is multiplied by complex amplitude and phase factors, for source and receiver terms, before the five-dimensional Fourier transformation back to the space–time domain.     

Estimating primaries by sparse inversion of the 3D curvelet transform and the L1-norm constraint

In this paper, we built upon the estimating primaries by sparse inversion (EPSI) method. We use the 3D curvelet transform and modify the EPSI method to the sparse inversion of the biconvex optimization and L1-norm regularization, and use alternating optimization to directly estimate the primary reflection coefficients and source wavelet. The 3D curvelet transform is used as a sparseness constraint when inverting the primary reflection coefficients, which results in avoiding the prediction subtraction process in the surface-related multiples elimination (SRME) method. The proposed method not only reduces the damage to the effective waves but also improves the elimination of multiples. It is also a wave equationbased method for elimination of surface multiple reflections, which effectively removes surface multiples under complex submarine conditions.     

Closed-loop SRME based on 3D L1-norm sparse inversion

In many situations, the quality of seismic imaging is largely determined by a proper multiple attenuation as preprocessing step. Despite the widespread application of surface-related multiple elimination (SRME) and estimation of primaries by sparse inversion (EPSI) for the removal of multiples, there still exist some limitations in the process of prediction and subtraction (SRME) or inversion (EPSI), which make the efficiency of multiple attenuation less satisfactory. To solve these problems, a new fully data-driven method called closed-loop SRME was proposed, which combines the robustness of SRME and the multi-dimensional inversion strategy of EPSI. Due to the selection of inversion approach and constraint, primary estimation by closed-loop SRME may fall into a local optimum during the solving process, which lowers the accuracy of deep information and weakens the continuity of seismic events. To avoid these shortcomings, we first modified the solving method for closed-loop SRME to an L1 norm-based bi-convex optimization method, which stabilizes the solution. Meanwhile, in the L1 norm constraint-based optimization process, the 3D sparsifying transform, being a 2D Curvelet-1D wavelet transform, is brought in as a 3D sparse constraint. In the 3D sparsifying domain, the data become sparser, thus making the result of optimization more accurate, the information of seismic events more continuous and the resolution higher. Examples on both synthetic and field data demonstrate that the method proposed in this paper, compared with the traditional SRME and closed-loop SRME, have an excellent effect on primary estimation and suppress multiples effectively.     

Source wavelet correction for practical Marchenko imaging: A sub-salt field-data example from the Gulf of Mexico

Imaging a target zone below a salt body can be challenging because large velocity contrasts in the overburden between the salt and surrounding sediments generate internal multiples, which interfere with primary reflections from the target level in the imaging process. This can lead to an erroneous interpretation of reflections in the sub-salt area if multiples are misinterpreted as primaries. The Marchenko redatuming method may enable imaging of the sub-salt target area where the effect of the multiply-scattering overburden is removed. This is achieved by creating a redatumed reflection response where virtual sources and receivers are located below the overburden using a macromodel of the velocity field and the surface reflection data. The accuracy of the redatumed data and the associated internal multiple removal, however, depends on the accurate knowledge of the source wavelet of the acquired reflection data. For the first time, we propose a method which can accurately and reliably correct the amplitudes of the reflection response in field data as required by the Marchenko method. Our method operates by iteratively and automatically updating the source function so as to cancel the most artefact energy in the focusing functions, which are also generated by the Marchenko method. We demonstrate the method on a synthetic dataset and successfully apply it to a field dataset acquired in a deep-water salt environment in the Gulf of Mexico. After the successful source wavelet estimation for the field dataset, we create sub-salt target-oriented images with Marchenko redatumed data. Marchenko images using the proposed source wavelet estimation show clear improvements, such as increased continuity of reflectors, compared to surface-based images and to conventional Marchenko images computed without the inverted source wavelet. Our improvements are corroborated by evidence in the literature and our own synthetic results.     

Analytical wavefront curvature correction to plane‐wave reflection coefficients for a weak‐contrast interface

Most amplitude versus offset (AVO) analysis and inversion techniques are based on the Zoeppritz equations for plane‐wave reflection coefficients or their approximations. Real seismic surveys use localized sources that produce spherical waves, rather than plane waves. In the far‐field, the AVO response for a spherical wave reflected from a plane interface can be well approximated by a plane‐wave response. However this approximation breaks down in the vicinity of the critical angle. Conventional AVO analysis ignores this problem and always utilizes the plane‐wave response. This approach is sufficiently accurate as long as the angles of incidence are much smaller than the critical angle. Such moderate angles are more than sufficient for the standard estimation of the AVO intercept and gradient. However, when independent estimation of the formation density is required, it may be important to use large incidence angles close to the critical angle, where spherical wave effects become important. For the amplitude of a spherical wave reflected from a plane fluid‐fluid interface, an analytical approximation is known, which provides a correction to the plane‐wave reflection coefficients for all angles. For the amplitude of a spherical wave reflected from a solid/solid interface, we propose a formula that combines this analytical approximation with the linearized plane‐wave AVO equation. The proposed approximation shows reasonable agreement with numerical simulations for a range of frequencies. Using this solution, we constructed a two‐layer three‐parameter least‐squares inversion algorithm. Application of this algorithm to synthetic data for a single plane interface shows an improvement compared to the use of plane‐wave reflection coefficients.     

L1范数约束被动源数据稀疏反演一次波估计

对于被动源地震数据,运用常规的互相关算法得到的虚拟炮记录中,不仅含有一次波反射信息,还包括了表面相关多次波.然而,通过传统的被动源数据稀疏反演一次波估计(EPSI)方法,可以求得只含有一次波,不含表面相关多次波的虚拟炮记录.本文改进了传统的被动源数据稀疏反演一次波估计问题的求解方法,将被动源稀疏反演一次波估计求解问题转化为双凸L1范数约束的最优化求解问题,避免了在传统的稀疏反演一次波估计过程中用时窗防止反演陷入局部最优化的情况.在L1范数约束最优化的求解过程中,又结合了2DCurvelet变换和小波变换,在2DCurvelet-wavelet域中,数据变得更加稀疏,从而使求得的结果更加准确,成像质量得到了改善.通过简单模型和复杂模型,验证了本文提出方法的有效性.     

子波相位不准对反演结果的影响(英文)

本文重点讨论在振幅谱估计准确的情况下,采用不同相位谱子波作为实际估计子波进行线性最小二乘反演,并对结果进行分析。除子波相位外,所有其它影响反演结果的因素均忽略。稀疏反射系数模型(块状波阻抗模型)反演结果表明:(1)使用不同相位谱子波进行反演,其反演结果合成的记录与原始记录都非常匹配,但反演的反射系数和声波阻抗结果与真实模型有差异;(2)反演结果的可靠程度主要与不同相位子波z变换的根的分布有关,当估计子波与真实子波Z变换的根的分布仅在单位圆附近有差异时,反演的反射系数和声波阻抗与真实模型很接近;(3)尽管反演前后地震记录都匹配了,并且评价反演结果好坏的柯西准则或改进柯西准则(反演参数没有进行自适应处理)已经达到了最优(最小),但反演结果与真实模型仍存在较大差异。最后,针对子波相位估计不准可能导致反演效果较差这个问题,我们提出采用求L1范数、丰度、变分、柯西准则(反演参数进行了自适应处理)或/和改进柯西准则(反演参数进行了自适应处理)的最优值或次优值作为评价准则的一种解决办法,理论上得到了好的效果。     

Viscoacoustic least-squares reverse-time migration of different-order free-surface multiples

Multiples have longer propagation paths and smaller reflection angles than primaries for the same source–receiver combination, so they cover a larger illumination area. Therefore, multiples can be used to image shadow zones of primaries. Least-squares reverse-time migration of multiples can produce high-quality images with fewer artefacts, high resolution and balanced amplitudes. However, viscoelasticity exists widely in the earth, especially in the deep-sea environment, and the influence of Q attenuation on multiples is much more serious than primaries due to multiples have longer paths. To compensate for Q attenuation of multiples, Q-compensated least-squares reverse-time migration of different-order multiples is proposed by deriving viscoacoustic Born modelling operators, adjoint operators and demigration operators for different-order multiples. Based on inversion theory, this method compensates for Q attenuation along all the propagation paths of multiples. Examples of a simple four-layer model, a modified attenuating Sigsbee2B model and a field data set suggest that the proposed method can produce better imaging results than Q-compensated least-squares reverse-time migration of primaries and regular least-squares reverse-time migration of multiples.     

Permeability inversion from low‐frequency seismoelectric logs in fluid‐saturated porous formations

In this paper, a method is proposed to invert permeability from seismoelectric logs in fluid‐saturated porous formations. From the analysis of both the amplitude and the phase of simulated seismoelectric logs, we find that the Stoneley wave amplitude of the ratio of the converted electric field to the pressure (REP) is sensitive to porosity rather than permeability while the tangent of the REP's phase is sensitive to permeability. The REP's phase reflects the phase discrepancy between the electric field and the pressure at the same location in the borehole. We theoretically derive the frequency‐dependent expression of the REP of the low‐frequency Stoneley wave and find that the tangent of the REP's argument is approximately in inverse proportion to permeability. We obtain an inversion formula and present the permeability inversion method by using the tangent of the REP's phase. To test this method, the permeabilities of different sandstones are inverted from the synthetic full‐waveform data of seismoelectric logs. A modified inversion process is proposed based on the analysis of inversion errors, by which the relative errors are controlled below 25% and they are smaller than those of the permeability inversion from the Stoneley wave of acoustic logs.     

地震多次波成像

地震资料含有各种类型多次波,而传统成像方法仅利用地震一次反射波成像,在地震成像前需将多次波去除.然而,多次波携带了丰富的地下结构信息,多次波偏移能够提供除反射波外的额外地下照明.修改传统逆时偏移方法,用包含一次反射波和多次波的原始记录代替震源子波,将SRME方法预测的表面多次波代替一次反射波作为输入数据,可将表面多次波成像.多次波成像的挑战和困难在于大量串扰噪声的产生,针对表面多次波成像中的成像噪声问题,将最小二乘逆时偏移方法与多次波分阶思想结合起来,发展可控阶数的表面多次波反演成像方法,有望初步实现高精度的表面多次波成像.在消除原始记录中的表面多次波后,通过逆散射级数方法预测得到层间多次波,将层间多次波作为逆时偏移方法的输入数据可将其准确归位到地下反射位置.数值实验表明,多次波成像能够有效地为地下提供额外照明,而可控阶表面多次波最小二乘逆时偏移成像方法几乎完全避免成像噪声.