Diffraction imaging by uniform asymptotic theory and double exponential fitting

Seismic diffracted waves carry valuable information for identifying geological discontinuities. Unfortunately, the diffraction energy is generally too weak, and standard seismic processing is biased to imaging reflection. In this paper, we present a dynamic diffraction imaging method with the aim of enhancing diffraction and increasing the signal‐to‐noise ratio. The correlation between diffraction amplitudes and their traveltimes generally exists in two forms, with one form based on the Kirchhoff integral formulation, and the other on the uniform asymptotic theory. However, the former will encounter singularities at geometrical shadow boundaries, and the latter requires the computation of a Fresnel integral. Therefore, neither of these methods is appropriate for practical applications. Noting the special form of the Fresnel integral, we propose a least‐squares fitting method based on double exponential functions to study the amplitude function of diffracted waves. The simple form of the fitting function has no singularities and can accelerate the calculation of diffraction amplitude weakening coefficients. By considering both the fitting weakening function and the polarity reversal property of the diffracted waves, we modify the conventional Kirchhoff imaging conditions and formulate a diffraction imaging formula. The mechanism of the proposed diffraction imaging procedure is based on the edge diffractor, instead of the idealized point diffractor. The polarity reversal property can eliminate the background of strong reflection and enhance the diffraction by same‐phase summation. Moreover,the fitting weakening function of diffraction amplitudes behaves like an inherent window to optimize the diffraction imaging aperture by its decaying trend. Synthetic and field data examples reveal that the proposed diffraction imaging method can meet the requirement of high‐resolution imaging, with the edge diffraction fully reinforced and the strong reflection mostly eliminated.     

Exact time-domain solutions for acoustic diffraction by a half plane

We derive exact time-domain solutions for scattering of acoustic waves by a half plane by inverse Fourier transforming the frequency-domain integral solutions. The solutions consist of a direct term, a reflected term and two diffraction terms. The diffracting edge induces step function discontinuities in the direct and reflected, terms at two shadow boundries. At each boundary, the associated diffraction term reaches a maximum amplitude of half the geometrical optics term and has a signum function discontinuity so that the total field remains continuous. We evaluate solutions for practical point source configurations by numerically convolving the impulse diffraction responses with a wavelet. We solve the associated problems of convolution with a singular, truncated diffraction operator by analytically derived correction techniques. We produce a zero offset section and compare it to a Kirchhoff integral solution. Our exact diffraction hyperbola exhibits noticeable asymmetry, with higher amplitudes on the reflector side of the edge. Near the apex of the hyperbola the Kirchhoff solution approximates the exact diffraction term symmetric in amplitude about the reflection shadow boundary, but omits the other low amplitude term necessary to ensure continuity at the direct shadow boundary.     

THE PERFECTLY REFLECTING WEDGE USED AS A CONTROL MODEL IN SEISMIC DIFFRACTION MODELLING*

The numerical modelling of seismic diffraction, e.g., at faults and other discontinuities, generally requires the use of fast approximate methods. The geophysicist responsible for the development of such numerical methods has a real need of exact solutions to certain ideal geometries to check the accuracy of his calculations. One such exact solution, which is available, is the acoustic wave solution to the perfectly reflecting wedge. The solution is three-dimensional and the source is an explosive point source. This model is ideal for seismic diffraction; the solution has the advantage of being exact, truly three-dimensional and of being in the convenient form of the temporal and spatial impulse response. More complicated sources which are extended in either space or time can, therefore, be modelled exactly by numerical integration. This paper presents some examples of the use of the perfectly reflecting wedge as a control model for an asymptotic high frequency diffraction modelling method. This control model has revealed that certain survey and wedge configurations can yield significant disagreement with, e.g., the Kirchhoff approximation. Such configurations could occur during VSP modelling when the survey lies in the near field or in the shadow zone of a high contrast fault. This control model has also been instructive in demonstrating why the high frequency, asymptotic, approximation is generally very good and has indicated a possible improvement to the Kirchhoff approximation for edge diffraction.     

Generalized non‐hyperbolic approximation for qP‐wave relative geometrical spreading in a layered transversely isotropic medium with a vertical symmetry axis

Compensation for geometrical spreading along the ray‐path is important in amplitude variation with offset analysis especially for not strongly attenuative media since it contributes to the seismic amplitude preservation. The P‐wave geometrical spreading factor is described by a non‐hyperbolic moveout approximation using the traveltime parameters that can be estimated from the velocity analysis. We extend the P‐wave relative geometrical spreading approximation from the rational form to the generalized non‐hyperbolic form in a transversely isotropic medium with a vertical symmetry axis. The acoustic approximation is used to reduce the number of parameters. The proposed generalized non‐hyperbolic approximation is developed with parameters defined by two rays: vertical and a reference rays. For numerical examples, we consider two choices for parameter selection by using two specific orientations for reference ray. We observe from the numerical tests that the proposed generalized non‐hyperbolic approximation gives more accurate results in both homogeneous and multi‐layered models than the rational counterpart.     

Common‐shot Fresnel beam migration based on wave‐field approximation in effective vicinity under complex topographic conditions

Gaussian beam migration is a versatile imaging method for geologically complex land areas, which overcomes the limitation of Kirchhoff migration in imaging multiple arrivals and has no steep‐dip limits of one‐way wave‐equation migration. However, its imaging accuracy depends on the geometry of Gaussian beam that is determined by the initial parameter of dynamic ray tracing. As a result, its applications in exploration areas with strong variations in topography and near‐surface velocity are limited. Combined with the concept of Fresnel zone and the theory of wave‐field approximation in effective vicinity, we present a more robust common‐shot Fresnel beam imaging method for complex topographic land areas in this paper. Compared with the conventional Gaussian beam migration for irregular topography, our method improves the beam geometry by limiting its effective half‐width with Fresnel zone radius. Moreover, through a quadratic travel‐time correction and an amplitude correction that is based on the wave‐field approximation in effective vicinity, it gives an accurate method for plane‐wave decomposition at complex topography, which produces good imaging results in both shallow and deep zones. Trials of two typical models and its application in field data demonstrated the validity and robustness of our method.     

3-D acoustic modelling of edge diffractions — revisited

An accurate, fast, and simple algorithm for 3-D acoustic modelling of seismic edge diffractions, originally developed in the 1980s, is revisited in this paper. The main objective is to reintroduce this simple approach to edge-diffraction modelling and for the first time give the details of the theory in the open literature. The method is based on a combination of Kirchhoff theory and uniform asymptotic techniques developed within a high-frequency assumption. The diffraction contributions are then computed at stationary edge points only, by analogy with the geometrical ray contributions associated with internal stationary points or specular points. To be able to handle sampling inaccuracies of the critical edge points, a modified algorithm is proposed. Its robustness is verified in case of scattering from a circular edge. Also the extension from rigid or free boundary conditions to the case of edges defined by two penetrable surfaces is discussed in this paper. Both experimental and synthetic 3-D data are presented to demonstrate the potential of this edge-diffraction modelling technique. Since all parameters needed in the computations are obtained from dynamic ray tracing, the algorithm can readily be incorporated in existing software packages for 3-D seismic ray modelling.     

Geometrical spreading factor of head waves in cased and open boreholes

Factors (coefficients) of geometrical spreading of compressional and shear head waves are calculated for an impulse multipole source of elastic oscillations in boreholes. It is shown that the length of the logging tool (i.e., the distance between the source and the nearest receiver) used for sonic measurements and the velocities of elastic waves in the medium both contribute to the factor of geometrical spreading. For a high-velocity formation (the shear wave velocity in the rock is higher than the compressional wave velocity in the fluid that fills the borehole) and a sufficiently long sonic tool with a monopole source, the coefficient of geometrical spreading is approximated by asymptotic formula 1/Z [Roever et al., 1974; Krauklis and Krauklis, 1976], where Z is the length of the tool; i.e., the amplitude of the compressional head wave decreases proportionally to the distance between the source and the receiver. In acoustically soft formations, this approximation is inapplicable even for long tools with length Z > 4 m. Waveforms in cased boreholes have a significant frequency dispersion even in case of good-quality cementing, and the factor of geometrical spreading there depends considerably on the length of the tool and the elastic properties of the rocks.     

复杂地表条件下叠前菲涅尔束偏移方法

高斯束偏移虽然克服了Kirchhoff偏移不能处理多波至和单程波动方程偏移不能对陡倾构造准确成像的问题,但在复杂地表条件下其偏移精度取决于所选择的初始束宽度,即当初始宽度较小时,近地表成像精度较高,但此时中深层成像质量较差;反之当初始宽度较大时,中深层成像质量提高,但近地表成像精度降低.针对高斯束偏移中深层和浅层成像精度的矛盾,本文发展了一种适用于陆地复杂地表条件的叠前菲涅尔束偏移方法.基于惠更斯-菲涅尔原理,本文首先给出了菲涅尔束的概念及其表征的格林函数,并采用有效邻域波场近似理论和反褶积成像条件,导出了复杂地表条件下叠前保幅深度偏移公式.最后,针对常规旁轴射线追踪中的数值噪音,给出了一种压制策略.同高斯束偏移相比,本文方法不仅解决了中深层和浅层成像精度的矛盾,而且提高了复杂地表条件下平面波的分解精度,使得偏移结果更加准确可靠.典型的模型算例验证了本文方法的有效性和稳健性.     

Second-order interpolation of later-arrival traveltimes

The performance of a 3D prestack migration of the Kirchhoff type can be significantly enhanced if the computation of the required stacking surface is replaced by an efficient and accurate method for the interpolation of diffraction traveltimes. Thus, input traveltimes need only be computed and stored on coarse grids, leading to considerable savings in CPU time and computer storage. However, interpolation methods based on a local approximation of the traveltime functions fail in the presence of triplications of the wavefront or later arrivals. This paper suggests a strategy to overcome this problem by employing the coefficients of a hyperbolic traveltime expansion to locate triplications and correct for the resulting errors in the interpolated traveltime tables of first and later arrivals.     

Geometrical constraints for the configuration of the Vrancea (Romania) intermediate-depth seismicity nest

For a seismogenic area like Vrancea (Romania) with well-defined geometrical features of the seismicity production in space and time, the numerical simulation of the earthquake generation process (e.g. cellular automaton) looks highly attractive. The delimitation, as accurately as possible, of the geometrical features of the seismically active system in the Vrancea subcrustal zone is essential to constrain the simulation modeling. As a first approximation, the seismicity pattern is close to a fault plane NE–SW oriented, extended roughly vertically between 60 and 170 km depth. A characteristic median plane is defined by minimizing the distance of hypocenters. The average distance of the hypocenters to the median plane is around 5 km. However, a more detailed investigation of the geometrical configuration of seismicity indicates a fragmentation of the active body located in the upper mantle in two segments. The seismicity pattern is well approximated by a planar distribution in each segment. In the transition zone, between the upper and lower segment, the hypocenter distribution is more dispersed and shows a disruption among the two planar segments, measured by about 9 km apart laterally one relative to the other. The two segments hosted the major Vrancea events recorded in the last two centuries (for which we have available location of acceptable accuracy). The narrow transition zone at about 100 km depth is interpreted as a weaker segment, possibly caused by a dehydration process or by an infiltration of asthenosphere material from the back side of the South-Eastern Carpathian arc system. It is still debatable if fragmentation in two seismically active segments reflects the existence of two neighbouring separate blocks (upper, continental and lower, oceanic block) or a consequence of a breaking process separating a continental block into two parts. The segmentation of the descending lithosphere and the edge effects are apparently stationary, at least for the time interval since 1985 to the present, for which the earthquake catalogue is reliable (homogeneous).     

Numerical tests of 3D true-amplitude zero-offset migration1

An amplitude-preserving migration aims at imaging compressional primary (zero-or) non-zero-offset reflections into 3D time or depth-migrated reflections so that the migrated wavefield amplitudes are a measure of angle-dependent reflection coeffcients. The principal objective is the removal of the geometrical-spreading factor of the primary reflections. Various migration/inversion algorithms involving weighted diffraction stacks proposed recently are based on Born or Kirchhoff approximations. Here, a 3D Kirchhoff-type zero-offset migration approach, also known as a diffraction-stack migration, is implemented in the form of a time migration. The primary reflections of the wavefield to be imaged are described a priori by the zero-order ray approximation. The aim of removing the geometrical- spreading loss can, in the zero-offset case, be achieved by not applying weights to the data before stacking them. This case alone has been implemented in this work. Application of the method to 3D synthetic zero-offset data proves that an amplitude-preserving migration can be performed in this way. Various numerical aspects of the true-amplitude zero-offset migration are discussed.     

ATTENUATION OF RANDOM NOISE BY 2-D AND 3-D CDP STACKING AND KIRCHHOFF MIGRATION*

Three-dimensional seismic surveys, in general, do not need the same high degree of CDP coverage as 2-D surveys to achieve a certain signal-to-noise ratio after migration. This can be shown theoretically for Kirchhoff migration and laterally uncorrelated noise. More precisely, there exists a formal relationship between the multiplicity of CDP coverage of a 3-D survey and that of a 2-D survey with the same signal-to-uncorrelated-noise ratio. Frequency and aperture are parameters in the corresponding expression. Heuristically the relationship can be obtained by applying the concept of the Fresnel zone. Though the mathematics in this paper refer to laterally uncorrelated noise, the underlying concepts can probably also be used for weakly correlated noise, e.g., for multiple reflections and for the low-frequency remnants of surface waves.     

A new approximation of the velocity-depth distribution and its application to the computation of seismic wave fields

Summary A new approximation of the velocity-depth distribution in a vertically inhomogeneous medium is suggested. This approximation guarantees the continuity of velocity and of its first and second derivatives and does not generate false low-velocity zones. It is very suitable for the computations of seismic wave fields in vertically inhomogeneous media by ray methods and its modifications, as it removes many false anomalies from the travel-time and amplitude-distance curves of seismic body waves. The ray integrals can be evaluated in a closed form; the resulting formulae for rays, travel times and geometrical spreading are very simple. They do not contain any transcendental functions (such asln (x) orsin ^–1, (x)) like other approximations; only the evaluation of one square root and of certain simple arithmetic expressions for each layer is required. From a computational point of view, the evaluation of ray integrals and of geometrical spreading is only slightly slower than for a system of homogeneous parallel layers and even faster than for a piece-wise linear approximation.     

Velocity model updating in prestack Kirchhoff time migration for PS converted waves: Part I – Theory

A velocity model updating approach is developed based on moveout analysis of the diffraction curve of PS converted waves in prestack Kirchhoff time migration. The diffraction curve can be expressed as a product of two factors: one factor depending on the PS converted‐wave velocity only, and the other factor depending on all parameters. The velocity‐dependent factor represents the hyperbolic behaviour of the moveout and the other is a scale factor that represents the non‐hyperbolic behaviour of the moveout. This non‐hyperbolic behaviour of the moveout can be corrected in prestack Kirchhoff time migration to form an inverse normal‐moveout common‐image‐point gather in which only the hyperbolic moveout is retained. This hyperbolic moveout is the moveout that would be obtained in an isotropic equivalent medium. A hyperbolic velocity is then estimated from this gather by applying hyperbolic moveout analysis. Theoretical analysis shows that for any given initial velocity, the estimated hyperbolic velocity converges by an iterative procedure to the optimal velocity if the velocity ratio is optimal or to a value closer to the optimal velocity if the velocity ratio is not optimal. The velocity ratio (V_P/V_S) has little effect on the estimation of the velocity. Applying this technique to a synthetic seismic data set confirms the theoretical findings. This work provides a practical method to obtain the velocity model for prestack Kirchhoff time migration.     

基于倾角-偏移距域道集的绕射波成像

地震绕射波是地下非连续性地质体的地震响应,绕射波成像对地下断层、尖灭和小尺度绕射体的识别具有重要的意义.在倾角域共成像点道集中,反射波同相轴表现为一条下凸曲线,能量主要集中在菲涅耳带内,绕射波能量则比较发散.由于倾角域菲涅耳带随偏移距变化而存在差异,因此本文提出一种在倾角-偏移距域道集中精确估计菲涅耳带的方法,在各偏移距的倾角域共成像点道集中实现菲涅耳带的精确切除,从而压制反射波.在倾角-偏移距域道集中还可以分别实现绕射波增强,绕射波同相轴相位校正,因此能量弱的绕射波可以清晰地成像.在倾角域共成像点道集中,反射波同相轴的最低点对应于菲涅耳带估计所用的倾角,因此本文提出一种在倾角域共成像点道集中直接自动拾取倾角场的方法.理论与实际资料试算验证了本文绕射波成像方法的有效性.     

地震勘探空间分辨力分析

未偏移地震剖面的空间分辨力研究通常采用第一菲涅尔带的解释方法,它是以几何地震学的射线理论为基础,仅考虑了地震波的运动学特征,不能详细描述地震波的动力学特征.本文从物理地震学的广义绕射叠加观点出发,综合分析地震波的运动学和动力学特征,由绕射叠加效应推导出了零炮检距地震道空间分辨力的计算公式.该公式的计算结果与第一菲涅尔带半径在数值上相当,并通过正演模拟验证了计算公式的合理性.因此利用物理地震学的原理研究地震勘探空间分辨力更符合地震波客观实际的传播特征.     

Series approximation of the Kirchhoff method for wavefield computation on a rough subsurface

ABSTRACT Computation of the wavefield due to reflection from an irregular surface is carried out for subsurfaces with large radii of curvature. The Kirchhoff approximation is proved to be sufficiently accurate provided that the acoustic wavelength is sufficiently small with respect to the asperities of the rough surface. For cases where the irregular surface does not fulfil this condition, a series solution is proposed. The first term of this series appears to be the result obtained by conventional Kirchhoff approximation. The series, initially developed in the space–wavenumber domain by Meecham, is transformed into the space–time domain, and the general expression for the series is obtained by calculation of the normal derivative of the field function. The series solution, restricted to the first two terms, is illustrated by application to three synthetic examples. Applications show that the series approximation obtained by the Kirchhoff method contributes significantly to the modelling of narrow, steep and deep structures and consequently it appears that the second term in the series cannot be ignored in the computation of the wavefields arising from a rough surface.     

Wavefield interpolation in 3D large‐step Fourier wavefield extrapolation

Extrapolating wavefields and imaging at each depth during three‐dimensional recursive wave‐equation migration is a time‐consuming endeavor. For efficiency, most commercial techniques extrapolate wavefields through thick slabs followed by wavefield interpolation within each thick slab. In this article, we develop this strategy by associating more efficient interpolators with a Fourier‐transform‐related wavefield extrapolation method. First, we formulate a three‐dimensional first‐order separation‐of‐variables screen propagator for large‐step wavefield extrapolation, which allows for wide‐angle propagations in highly contrasting media. This propagator significantly improves the performance of the split‐step Fourier method in dealing with significant lateral heterogeneities at the cost of only one more fast Fourier transform in each thick slab. We then extend the two‐dimensional Kirchhoff and Born–Kirchhoff local wavefield interpolators to three‐dimensional cases for each slab. The three‐dimensional Kirchhoff interpolator is based on the traditional Kirchhoff formula and applies to moderate lateral velocity variations, whereas the three‐dimensional Born–Kirchhoff interpolator is derived from the Lippmann–Schwinger integral equation under the Born approximation and is adapted to highly laterally varying media. Numerical examples on the three‐dimensional salt model of the Society of Exploration Geophysicists/European Association of Geoscientists demonstrate that three‐dimensional first‐order separation‐of‐variables screen propagator Born–Kirchhoff depth migration using thick‐slab wavefield extrapolation plus thin‐slab interpolation tolerates a considerable depth‐step size of up to 72 ms, eventually resulting in an efficiency improvement of nearly 80% without obvious loss of imaging accuracy. Although the proposed three‐dimensional interpolators are presented with one‐way Fourier extrapolation methods, they can be extended for applications to general migration methods.     

Estimation of the rocks statistical parameters from traveltime measurements

In this paper the method for estimating the statistical parameters of the medium from traveltime measurements of refracted waves is applied to study the statistical characteristics of crystalline rocks at the Multifunctional Station Faido (Gotthard Base Tunnel, Switzerland). The method is based on the geometrical optics (GO) approximation. A covariance function for traveltime fluctuations has been obtained by considering quasihomogeneous fluctuations of sound velocity in a plain-stratified medium. Strongly anisometric (having unequal dimensions in different directions) random inhomogeneities were embedded in this medium. To estimate the statistical parameters around the tunnel, the traveltime fluctuations are calculated. It is assumed that each observation of traveltime-distance relation for a given shot-receiver group corresponds to a particular realization of a medium statistical ensemble. By calculating the variance and the zero cross intervals of the first derivative of traveltime fluctuations, the standard deviation of the velocity fluctuations and the characteristic horizontal scale of the inhomogeneities are estimated. Although the method allows to obtain the characteristic lengths of the inhomogeneities in vertical as well as in horizontal direction, the limited offset of the field data made it only possible to measure the latter. The estimated horizontal characteristic scale is about 13 m, which is reasonably close to the direct geological measurements in the studied region, where quartz lenses are dominant among the inhomogeneities. The standard deviation of the velocity is estimated as 4.5%, which might be caused by the fractured structure around the tunnel and also by the fault zone near the study area.