AUTOMATED STATIC CORRECTIONS*

The widespread use of common depth point techniques has emphasized the need for accurate static corrections. Manual interpretation methods can give excellent results, but a computer technique is desirable because of the great volumn of data recorded in common depth point shooting. The redundancy inherent in common depth point data may be used to compute a statistical estimate of the static corrections. The corrections are assumed to be time-invarient, surface-consistent, and independent of frequency. Surface consistency implies that all traces from a particular shot will receive the same shot static correction and all traces from a particular receiver position will receive the same receiver correction. Time shifts are computed for all input traces using crosscorrelation functions between common depth point traces. The time shift for each trace is composed of a shot static, a receiver static, residual normal moveout if present, and noise. Estimates of the shot and receiver static corrections are obtained by averaging different sets of the measured time shifts. Time shifts which are greatly in error are detected and removed from the computations. The method is useful for data which has a moderate to good signal to noise ratio. Residual normal moveout should be corrected before estimating the statics. The program estimates the statics for correctly stacking common depth point traces but it is not sensitive to constant or very slowly changing static errors.     

共接收点倾斜叠加波动方程偏移

共接收点倾斜叠加波动方程偏移,本质上是一种叠前偏移方法.每给定一个斜率P,对经过叠前（动校正前）常规处理的地震记录中的各共接收点道集,沿直线t=τ+px进行倾斜叠加,就形成一个共接收点倾斜叠加剖面.对之进行波动方程偏移,该偏移剖面将代表地下真实构造.对一系列的p,我们可以得到一系列这样的偏移剖面.对它们作共接收点叠加,偏移叠加剖面的信噪比将超过水平叠加剖面.本文导出了在均匀、水平层状及非均匀介质条件下的共接收点倾斜叠加波动方程偏移算法.     

AUTOMATIC EDITING OF NOISY SEISMIC DATA1

Seismic data often contain traces that are dominated by noise; these traces should be removed (edited) before multichannel filtering or stacking. Noise bursts and spikes should be edited before single channel filtering. Spikes can be edited using a running median filter with a threshold; noise bursts can be edited by comparing the amplitudes of each trace to those of traces that are nearby in offset-common midpoint space. Relative amplitude decay rates of traces are diagnostic of their signal-to-noise (S/N) ratios and can be used to define trace editing criteria. The relative amplitude decay rate is calculated by comparing the time-gated trace amplitudes to a control function that is the median trace amplitude as a function of time, offset, and common midpoint. The editing threshold is set using a data-adaptive procedure that analyses a histogram of the amplitude decay rates. A performance evaluation shows that the algorithm makes slightly fewer incorrect trace editing decisions than human editors. The procedure for threshold setting achieves a good balance between preserving the fold of the data and removing the noisiest traces. Tests using a synthetic seismic line show that the relative amplitude decay rates are diagnostic of the traces’S/N ratios. However, the S/N ratios cannot be accurately usefully estimated at the start of processing, where noisy-trace editing is most needed; this is the fundamental limit to the accuracy of noisy trace editing. When trace equalization is omitted from the processing flow (as in amplitude-versus-offset analysis), precise noisy-trace editing is critical. The S/N ratio of the stack is more sensitive to type 2 errors (failing to reject noisy traces) than it is to type 1 errors (rejecting good traces). However, as the fold of the data decreases, the S/N ratio of the stack becomes increasingly sensitive to type 1 errors.     

地震资料处理中的聚束滤波方法

聚束滤波方法的基本原理,包括信号与噪音模型、设计准则和所获得的聚束滤波器．本文所用的自适应聚束滤波器是根据具有参数化的动校正量、振幅及相位随偏移距变化（ＭＶＯ、ＡＶＯ及ＰＶＯ）的一次反射和多次反射信号模型而设计的．但在实际应用中ＰＶＯ通常不能得到证实．结果分析能为进一步的从ＡＶＯ及ＰＶＯ获取岩石类型信息的专题研究提供资料．人工合成资料的例子给出了参数化的自适应聚束滤波的实施细节和所设计滤波器的响应特性．实际资料的例子表明数据自适应聚束滤波在叠前共中点道集地震多次波消除上比Ｒａｄｏｎ变换方法更灵活、有效．类似于其它叠前处理过程,自适应聚束滤波的优越性在信噪比较高的资料上体现得最为明显．     

A COMPARISON OF STATISTICAL AND DETERMINISTIC WIENER DECONVOLUTION OF DEEP-TOW SEISMIC DATA*

Wiener ‘spiking’ deconvolution of seismic traces in the absence of a known source wavelet relies upon the use of digital filters, which are optimum in a least-squares error sense only if the wavelet to be deconvolved is minimum phase. In the marine environment in particular this condition is frequently violated, since bubble pulse oscillations result in source signatures which deviate significantly from minimum phase. The degree to which the deconvolution is impaired by such violation is generally difficult to assess, since without a measured source signature there is no optimally deconvolved trace with which the spiked trace may be compared. A recently developed near-bottom seismic profiler used in conjunction with a surface air gun source produces traces which contain the far-field source signature as the first arrival. Knowledge of this characteristic wavelet permits the design of two-sided Wiener spiking and shaping filters which can be used to accurately deconvolve the remainder of the trace. In this paper the performance of such optimum-lag filters is compared with that of the zero-lag (one-sided) operators which can be evaluated from the reflected arrival sequence alone by assuming a minimum phase source wavelet. Results indicate that the use of zero-lag operators on traces containing non-minimum phase wavelets introduces significant quantities of noise energy into the seismic record. Signal to noise ratios may however be preserved or even increased during deconvolution by the use of optimum-lag spiking or shaping filters. A debubbling technique involving matched filtering of the trace with the source wavelet followed by optimum-lag Wiener deconvolution did not give a higher quality result than can be obtained simply by the application of a suitably chosen Wiener shaping filter. However, cross correlation of an optimum-lag spike filtered trace with the known ‘actual output’ of the filter when presented with the source signature is found to enhance signal-to-noise ratio whilst maintaining improved resolution.     

3D simultaneous joint PP‐PS prestack seismic inversion at Schiehallion field,United Kingdom Continental Shelf

A workflow for simultaneous joint PP‐PS prestack inversion of data from the Schiehallion field on the United Kingdom Continental Shelf is presented and discussed. The main challenge, describing reasonable PS to PP data registration before any prestack or joint PP‐PS inversion, was overcome thanks to a two‐stage process addressing the signal envelope, then working directly on the seismic data to estimate appropriate time‐variant time‐shift volumes. We evaluated the benefits of including PS along with PP prestack seismic data in a joint inversion process to improve the estimated elastic property quality and also to enable estimation of density compared with other prestack and post‐stack inversion approaches. While the estimated acoustic impedance exhibited a similar quality independent of the inversion used (PP post‐stack, PP prestack or joint PP‐PS prestack inversion) the shear impedance estimation was noticeably improved by the joint PP‐PS prestack inversion when compared to the PP prestack inversion. Finally, the density estimated from joint PP and PS prestack data demonstrated an overall good quality, even where not well‐controlled. The main outcome of this study was that despite several data‐related limitations, inverting jointly correctly processed PP and PS data sets brought extra value for reservoir delineation as opposed to PP‐only or post‐stack inversion.     

IMPROVED LOW-FREQUENCY DECAY ESTIMATION USING THE MULTITAPER SPECTRAL ANALYSIS METHOD1

Spectral analysis is one of the most ubiquitous signal processing tools used in exploration geophysics. Among many applications, it is used simply to look at the frequency content of seismic traces, to find notches, to estimate wavelets under the minimum-phase assumption, and to match broadband synthetic seismograms to seismic data. Seismic spectra exhibit very large dynamic ranges, particularly at low frequencies. Estimation of low-frequency decay is very important for accurate modelling. However, when using traditional spectral estimates incorporating smoothing windows, too much sidelobe energy leaks from high power into low power areas, spoiling our ability to estimate low-frequency spectral decay. The multitaper method of spectral analysis due to D. Thomson does not employ just a single window, but rather a set of orthogonal data tapers. It is possible to have much less sidelobe contamination, while maintaining a stable estimate. The trace is tapered by each of a subset of the orthogonal tapers, and a raw spectral estimate produced in each case. These are combined to produce a final spectral estimate. The technique can be made adaptive by applying different weights to the different raw spectra at different frequencies. A comparison of seismic spectral estimation using this multitaper technique with a traditional approach having the same analysis bandwidth and stability demonstrates the very different estimates of spectral decay in the areas of high dynamic range. The multitaper approach provides estimates with much reduced sidelobe leakage, and hence is a very appealing method for reflection seismology.     

3D post‐stack interval velocity analysis with effective use of datuming

Interval velocity analysis using post‐stack data has always been a desire, mainly for 3D data sets. In this study we present a method that uses the unique characteristics of migrated diffractions to enable interval velocity analysis from three‐dimensional zero‐offset time data. The idea is to perform a standard three‐dimensional prestack depth migration on stack cubes and generate three‐dimensional common image gathers that show great sensitivity to velocity errors. An efficient ‘top‐down’ scheme for updating the velocity is used to build the model. The effectiveness of the method is related to the incorporation of wave equation based post‐stack datuming in the model building process. The proposed method relies on the ability to identify diffractions along redatumed zero‐offset data and to analyse their flatness in the migrated local angle domain. The method can be considered as an additional tool for a complete, prestack depth migration based interval velocity analysis.     

Robust multiple suppression using adaptive beamforming

Minimum variance unbiased (MVU) beamforming is a type of multichannel filtering which extracts coherent signals without distortion, whilst minimizing residual noise power. Adaptive beamforming estimates signal and noise characteristics as part of the extraction process. The adaptive beamformer used here is designed from models of primary and multiple reflection signals having parametrically specified moveout and amplitude variation with offset (MVO and AVO). Phase variation with offset (PVO) can also be included but it is not usually justified in practice. The resulting analysis provides data for input into AVO and PVO schemes for obtaining lithological information. Synthetic data examples illustrate details of implementation of parametric adaptive MVU beamforming and the response characteristics of the resultant design. Real data examples show that data-adaptive beamforming is more flexible and more effective in attenuating multiples in prestack common-midpoint seismic data than Radon transform methods. In common with other prestack multichannel processes, the advantages of beamforming are shown to best effect in data with a good signal-to-noise ratio.     

Handling dip discrimination phenomenon in common‐reflection‐surface stack via combination of output‐imaging‐scheme and migration/demigration

In the application of a conventional common‐reflection‐surface (CRS) stack, it is well‐known that only one optimum stacking operator is determined for each zero‐offset sample to be simulated. As a result, the conflicting dip situations are not taken into account and only the most prominent event contributes to any a particular stack sample. In this paper, we name this phenomenon caused by conflicting dip problems as ‘dip discrimination phenomenon’. This phenomenon is not welcome because it not only leads to the loss of weak reflections and tips of diffractions in the final zero‐offset‐CRS stacked section but also to a deteriorated quality in subsequent migration. The common‐reflection‐surface stack with the output imaging scheme (CRS‐OIS) is a novel technique to implement a CRS stack based on a unified Kirchhoff imaging approach. As far as dealing with conflicting dip problems is concerned, the CRS‐OIS is a better option than a conventional CRS stack. However, we think the CRS‐OIS can do more in this aspect. In this paper, we propose a workflow to handle the dip discrimination phenomenon based on a cascaded implementation of prestack time migration, CRS‐OIS and prestack time demigration. Firstly, a common offset prestack time migration is implemented. Then, a CRS‐OIS is applied to the time‐migrated common offset gather. Afterwards, a prestack time demigration is performed to reconstruct each unmigrated common offset gather with its reflections being greatly enhanced and diffractions being well preserved. Compared with existing techniques dealing with conflicting dip problems, the technique presented in this paper preserves most of the diffractions and accounts for reflections from all possible dips properly. More importantly, both the post‐stacked data set and prestacked data set can be of much better quality after the implementation of the presented scheme. It serves as a promising alternative to other techniques except that it cannot provide the typical CRS wavefield attributes. The numerical tests on a synthetic Marmousi data set and a real 2D marine data set demonstrated its effectiveness and robustness.     

THE ESTIMATION OF SIGNAL SPECTRA AND RELATED QUANTITIES BY MEANS OF THE MULTIPLE COHERENCE FUNCTION *

A seismic trace recorded with suitable gain control can be treated as a stationary time series. Each trace, χ_j(t), from a set of traces, can be broken down into two stationary components: a signal sequence, α_j(t) *s(t—τ_j), which correlates from trace to trace, and an incoherent noise sequence, n_j(t), which does not correlate from trace to trace. The model for a seismic trace used in this paper is thus χ_j(t) =α_j(t) * s(t—τ_j) +n_j(t) where the signal wavelet α_j(t), the lag (moveout) of the signal τ_j, and the noise sequence n_j(t) can vary in any manner from trace to trace. Given this model, a method for estimating the power spectra of the signal and incoherent noise components on each trace is presented. The method requires the calculation of the multiple coherence function γ_j(f) of each trace. γ_j(f) is the fraction of the power on traced at frequency f that can be predicted in a least-square error sense from all other traces. It is related to the signal-to-noise power ratio ρ_j(f) by where K_j(f) can be computed and is in general close to 1.0. The theory leading to this relation is given in an Appendix. Particular attention is paid to the statistical distributions of all estimated quantities. The statistical behaviour of cross-spectral and coherence estimates is complicated by the presence of bias as well as random deviations. Straightforward methods for removing this bias and setting up confidence limits, based on the principle of maximum likelihood and the Goodman distribution for the sample multiple coherence, are described. Actual field records differ from the assumed model mainly in having more than one correctable component, components other than the required sequence of reflections being lumped together as correlated noise. When more than one correlatable component is present, the estimate for the signal power spectrum obtained by the multiple coherence method is approximately the sum of the power spectra of the correlatable components. A further practical drawback to estimating spectra from seismic data is the limited number of degrees of freedom available. Usually at least one second of stationary data on each trace is needed to estimate the signal spectrum with an accuracy of about 10%. Examples using synthetic data are presented to illustrate the method.     

Common-reflection-surface (CRS) stack for common offset

We provide a data-driven macro-model-independent stacking technique that migrates 2D prestack multicoverage data into a common-offset (CO) section. We call this new process the CO common-reflection-surface (CRS) stack. It can be viewed as the generalization of the zero-offset (ZO) CRS stack, by which 2D multicoverage data are stacked into a well-simulated ZO section. The CO CRS stack formula can be tailored to stack P-P, S-S reflections as well as P-S or S-P converted reflections. We point out some potential applications of the five kinematic data-derived attributes obtained by the CO CRS stack for each stack value. These include (i) the determination of the geometrical spreading factor for reflections, which plays an important role in the construction of the true-amplitude CO section, and (ii) the separation of the diffractions from reflection events. As a by-product of formulating the CO CRS stack formula, we have also derived a formula to perform a data-driven prestack time migration.     

Prestack nonstationary deconvolution based on variable-step sampling in the radial trace domain

The conventional nonstationary convolutional model assumes that the seismic signal is recorded at normal incidence. Raw shot gathers are far from this assumption because of the effects of offsets. Because of such problems, we propose a novel prestack nonstationary deconvolution approach. We introduce the radial trace （RT） transform to the nonstationary deconvolution, we estimate the nonstationary deconvolution factor with hyperbolic smoothing based on variable-step sampling （VSS） in the RT domain, and we obtain the high-resolution prestack nonstationary deconvolution data. The RT transform maps the shot record from the offset and traveltime coordinates to those of apparent velocity and traveltime. The ray paths of the traces in the RT better satisfy the assumptions of the convolutional model. The proposed method combines the advantages of stationary deconvolution and inverse Q filtering, without prior information for Q. The nonstationary deconvolution in the RT domain is more suitable than that in the space-time （XT） domain for prestack data because it is the generalized extension of normal incidence. Tests with synthetic and real data demonstrate that the proposed method is more effective in compensating for large-offset and deep data.     

Geomechanical property estimation of unconventional reservoirs using seismic data and rock physics

An extension of a previously developed rock physics model is made that quantifies the relationship between the ductile fraction of a brittle/ductile binary mixture and the isotropic seismic reflection response. By making a weak scattering (Born) approximation and plane wave (eikonal) approximation, with a subsequent ordering according to the angles of incidence, singular value decomposition analyses are performed to understand the stack weightings, number of stacks, and the type of stacks that will optimally estimate two fundamental rock physics parameters – the ductile fraction and the compaction and/or diagenesis. It is concluded that the full PP stack, i.e., sum of all PP offset traces, and the “full” PS stack, i.e., linear weighted sum of PS offset traces, are the two optimal stacks needed to estimate the two rock physics parameters. They dominate over both the second‐order amplitude variation offset “gradient” stack, which is a quadratically weighted sum of PP offset traces that is effectively the far offset traces minus the near offset traces, and the higher order fourth order PP stack (even at large angles of incidence). Using this result and model‐based Bayesian inversion, the seismic detectability of the ductile fraction (shown by others to be the important rock property for the geomechanical response of unconventional reservoir fracking) is demonstrated on a model characteristic of the Marcellus shale play.     

WEIGHTED STACKING OF SEISMIC DATA USING AMPLITUDE-DECAY RATES AND NOISE AMPLITUDES1

The signal-to-noise (S/N) ratio of seismic reflection data can be significantly enhanced by stacking. However, stacking using the arithmetic mean (straight stacking) does not maximize the S/N ratio of the stack if there are trace-to-trace variations in the S/N ratio. In this case, the S/N ratio of the stack is maximized by weighting each trace by its signal amplitude divided by its noise power, provided the noise is stationary. We estimate these optimum weights using two criteria: the amplitude-decay rate and the measured noise amplitude for each trace. The amplitude-decay rates are measured relative to the median amplitude-decay rate as a function of midpoint and offset. The noise amplitudes are measured using the data before the first seismic arrivals or at late record times. The optimum stacking weights are estimated from these two quantities using an empirical equation. Tests with synthetic data show that, even after noisy-trace editing, the S/N ratio of the weighted stack can be more than 10 dB greater than the S/N ratio of the straight stack, but only a few decibels more than the S/N ratio of the trace equalized stack. When the S/N ratio is close to 0 dB, a difference of 4 dB is clearly visible to the eye, but a difference of 1 dB or less is not visible. In many cases the S/N ratio of the trace-equalized stack is only a few decibels less than that of the optimum stack, so there is little to be gained from weighted stacking. However, when noisy-trace editing is omitted, the S/N ratio of the weighted stack can be more than 10 dB greater than that of the trace-equalized stack. Tests using field data show that the results from straight stacking, trace-equalized stacking, and weighted stacking are often indistinguishable, but weighted stacking can yield slight improvements on isolated portions of the data.     

RESIDUAL STATICS ESTIMATION BY STACK-POWER MAXIMIZATION IN THE FREQUENCY DOMAIN1

Traditionally, residual static corrections are based on timeshifts estimated for individual CMP sorted traces, which are later resolved into surface-consistent statics. This is a stable and attractive procedure because the data flow is simple and the memory storage required is limited. An alternative station-oriented method maximizing the stack-power estimates surface-consistent static corrections directly. The statics evaluation in this method involves several CMP gathers, which should improve the prediction of statics on noise-contaminated data. In this paper the performance of the above methods will be compared using synthetic as well as real seismic data. Neither method is capable of estimating large statics compared to the dominating period, because local optimization might fail. Global Monte Carlo search by, for instance, simulated annealing has been used to overcome the cycle-skipping problems when proper field statics are missing. Although this procedure is computationally very heavy, it may be the only way to deal with large residual statics. In order to enlarge the operational field for local optimization, it is suggested that the stack-power in the frequency domain is maximized. This makes it easy to change the frequency band during the optimization. Making use of the frequency domain will also normally be faster than the traditional time-domain optimization even for a limited number of iterations. Moreover, the main memory storage required can be significantly reduced, since it is only necessary to keep the frequency band in the memory, where the signal-to-noise ratio is good.     

CRS-stack-based residual static correction

In the case of onshore data sets, the acquired reflection events can be strongly impaired due to rough top‐surface topography and inhomogeneities in the uppermost low‐velocity layer, the so‐called weathering layer. Without accounting for these influences, the poor data quality will make data processing very difficult. Usually, the correction for the top‐surface topography is not perfect. The residuals from this correction and the influence of the weathering layers lead to small distortions along the reflection events. We integrated a residual static correction method into our data‐driven common‐reflection‐surface‐stack‐based imaging workflow to further eliminate such distortions. The moveout‐corrected traces and the stacked pilot trace are cross‐correlated to determine a final estimate of the surface‐consistent residual statics in an iterative manner. As the handling of top‐surface topography within the common‐reflection‐surface stack is discussed in a separate paper in this special issue, the corresponding residual static correction will be explained in more detail. For this purpose, the results obtained with a data set from the Arabian Peninsula will be presented.     

Random noise attenuation using a structure-oriented weighted singular value decomposition

Singular value decomposition (SVD) is a useful method for random noise suppression in seismic data processing. A structure-oriented SVD (SOSVD) approach which incorporates structure prediction to the SVD filter is effcient in attenuating noise except distorting seismic events at faults and crossing points. A modified SOSVD approach using a weighted stack, called structure-oriented weighted SVD (SOWSVD), is proposed. In this approach, the SVD filter is used to attenuate noise for prediction traces of a primitive trace which are produced via the plane-wave prediction. A weighting function related to local similarity and distance between each prediction trace and the primitive trace is applied to the denoised prediction traces stacking. Both synthetic and field data examples suggest the SOWSVD performs better than the SOSVD in both suppressing random noise and preserving the information of the discontinuities for seismic data with crossing events and faults.

     

Velocity analysis after migration

The double‐square‐root (DSR) equation used in pre‐stack migration is formulated in terms of velocity‐dependent and velocity‐independent terms. The velocity‐dependent term is shown to be the hyperbolic normal moveout (NMO) correction, whereas the velocity‐independent term is related to the recording geometry only. This separation of the velocity‐dependent term offers a means of applying vertical corrections to an initial migration velocity field. Using this concept, procedures are described both for velocity determination and for achieving improved structural imaging.
This decoupling is accurate both for constant‐velocity media and for media whose velocity varies as a function of depth. In media whose velocity varies as a function of both space and depth, a procedure is described for building velocity models through common‐image gather (CIG) stacking following prestack depth migration (PSDM) and time conversion (TC). This so‐called PSDM‐TC stack procedure provides a means of (a) incorporating both vertical and lateral velocity updates into an initial velocity model, (b) obtaining improved structural imaging by using a non‐optimal velocity model for the prestack depth migration, and (c) updating velocity by flattening CIGs and maximizing stack energy. The procedure can be applied to both P‐P wave and P‐SV wave migration.     

L'EGALISATION SPECTRALE: UN MOYEN D'AMELIORER LA QUALITE DES DONNEES SISMIQUES*

In land seismic surveys spectrum equalization can increase the quality of seismic data in a selected frequency band. The power of lower frequencies in the spectrum of input traces is generally greater than that of higher frequencies, particularly in land seismic surveys because of ground roll. In order to improve the quality of seismic data it is necessary to raise the energy of higher frequencies to the same level as that of lower frequencies, without alteration of the phases. The first step of the method is to compute the amplitude spectrum of each input trace to determine a weighting function which is then applied to the amplitude spectrum in order to balance it. The function is the inverse of the short wavelength variation of the amplitude spectrum. The short wavelength variation can be obtained by interpolation between average values of the modulus of the amplitude spectrum computed in narrow bands within a selected band of frequencies. Another way of obtaining the short wavelength variation is to apply a low-pass filter to the amplitude spectrum. The calculations are readily performed in the frequency domain by the Fourier transform. Spectrum equalization is automatically adjusted to each trace and does not modify the average amplitude in the time domain. However, as the frequency band and energy of the ground roll both vary according to the distance from the shot, spectrum equalization tends to make the spectrum of output traces independent of the offset distance. The use of spectrum equalization before any two-dimensional filtering improves ground roll elimination. Continuity and resolution of horizons are also increased by spectrum equalization before CDP stack. Several examples of applications of spectrum equalization to seismic land and marine surveys are shown.