|Series||Geophysical References Series|
|Title||Digital Imaging and Deconvolution: The ABCs of Seismic Exploration and Processing|
|Author||Enders A. Robinson and Sven Treitel|
|Store||SEG Online Store|
What is wavelet processing? To attenuate multiples and to increase the signal-to-noise ratio, the field trace x generally undergoes such operations as signature deconvolution, predictive deconvolution (either spike or gapped or both), frequency and wavenumber filtering, stacking, and other ancillary operations. The essential data-processing operations along with any extra signal-enhancing techniques transform the field trace x into the gross trace .
The gross wavelet g does not have most of the undesirable effects, such as distortions of the source signal, instruments, absorption, and multiple reflections, all of which are present in the field wavelet. Thus, the gross wavelet is an improvement. As a rule, however, the gross wavelet will be neither minimum delay nor causal. Often it will be nonsymmetric and frequently quite spread out. Thus, the gross wavelet is not entirely desirable for a seismic interpreter.
Accordingly, various filtering operations, which come under the heading of wavelet processing, are used to improve the quality and interpretability of final displays. Wavelet processing creates a visual picture of the subsurface, which then can be interpreted meaningfully by geologists and geophysicists. At this final stage, the physical properties of transmitted wave motion are no longer important or even needed. Instead, a nonphysical wavelet, called the interpreter wavelet h, is introduced solely for its visual convenience.
The interpreter wavelet h has no basis in physics. The purpose of wavelet processing is to convert the gross wavelet g into a sharp zero-phase interpreter wavelet h. A two-sided shaping filter f can be designed to transform g at the input to h at the output. We write this operation as . Application of this shaping filter to the gross trace yields as output a trace that we might call the interpreter trace z, which is given by
The interpreter trace provides us with a clearer picture than does the gross trace. The interpreter trace is the final seismic display trace — that is, the ultimate processed seismic trace. It consists of an interpreter wavelet attached to each of the estimated reflection coefficients in the reflectivity function. At this point, all wave-propagation features that were present in the original received trace — such as source and receiver effects, inelastic absorption, and multiple reflections — should have been removed as well as possible.
We might say that the wavelet-processing step in the overall processing sequence changes a physical seismic wavelet into an artificial interpreter wavelet h. The criterion for choosing a particular form of interpreter wavelet has to do with resolution. Interpreter wavelet h should be as sharp as is possible while remaining consistent with the data. In other words, the interpreter wavelet should have as large a bandwidth as possible.
A wavelet’s sharpness depends on its phase spectrum as well as on its amplitude spectrum. For a given amplitude spectrum, the greatest sharpness occurs when the phase spectrum is zero — that is, when the wavelet is zero phase. In summary, an interpreter wavelet should have a broad amplitude spectrum and a zero-phase spectrum. Such a wavelet is symmetric and damps with the greatest degree of rapidity that is consistent with its amplitude spectrum.
In summary, wavelet processing involves filtering the seismic trace to change the shape of the underlying seismic wavelet in advantageous ways. This filtering operation is implemented with a shaping filter. The physical field wavelet is a one-sided (or causal) signal whose initial nonzero value occurs at the arrival time of the reflected event. Because the earth necessarily acts as a causal filter, clearly it is impossible for the physical wavelet to have a precursor (i.e., nonzero amplitudes) arriving before the wavelet’s actual arrival time. However, the seismic wavelet takes some time (a few milliseconds or so) to build up to its peak amplitude, and thus this peak amplitude occurs after the wavelet’s arrival time. Almost always, some secondary peaks follow the main wavelet peak. The objective of wavelet processing is to convert the physical field seismic wavelet into an artificial interpreter wavelet, which is more easily pickable and hence more pleasing to the eye of an interpreter. Like great art, this artifact does not have to be a wavelet that satisfies physical laws, but it must be esthetically pleasing.
Over the years, seismic interpreters have decided that the most pleasing wavelet is one having a zero-phase-lag characteristic — namely, the familiar zero-phase wavelet. A wavelet is zero phase if and only if the wavelet is an autocorrelation function. In other words, the concept of a zero-phase wavelet and the concept of an autocorrelation function are identical. An autocorrelation is symmetric about time zero, it has its maximum value at time zero, and it rapidly damps off in both time directions. Thus, in wavelet processing, the actual seismic wavelet, which has its beginning at the arrival time of the reflection, is replaced (as well as possible) by a zero-phase wavelet with its center (i.e., its maximum value) coinciding with the reflection’s arrival time.
Perhaps at this point we can insert a few further comments about nature and nurture. Seismology is concerned with the physics of wave propagation through the earth. The actual recorded seismic waves, before any processing, contain the greatest amount of information about this physical phenomenon. A petroleum geologist, however, is interested in a picture of a cross section of the earth. The painting of a good picture requires a significant element of art. The whole purpose of data processing is to discard the strictly physical picture of wave propagation and replace it by a more artistic rendition of the geology. If the artwork is done properly, interpreters can look at this picture and see the geology with their own eyes. The ultimate goal of digital processing is to relieve the geologist from having to worry about wave propagation effects; these effects should be corrected as much as possible to produce a clearer picture of the subsurface. Nature gives us the recorded seismic traces, but nurture, in the form of deconvolution, stacking, wavelet processing, and migration, provides the final geologic picture. The discipline known as seismic stratigraphy attests to the efficacy of seismic-processing methods that are designed to satisfy these criteria.
There are many approaches to wavelet processing, and they embrace a wide range of linear and nonlinear methods. However, as time goes on, we expect that the various piecemeal elements employed in current seismic-processing systems ultimately will meld with the 3D inversion schemes now under development. The result could be an overall inversion technique that uses geophysical as well as geologic effects to give a 3D geometric and lithologic picture of the subsurface.
Wavelet-processing methods in geophysics usually entail use of shaping and deconvolution filters. Wavelet processing is a two-part process. The first and most important part is the specification of a mathematical model that adequately represents the geophysical phenomenon under consideration. The second part of wavelet processing is implementation of the numerical procedures to carry out the needed filtering operation. The two parts go hand in hand. If one specifies key model components in a statistical sense, then a statistical wavelet-processing method results. Because specification precedes implementation, a large specification error can overwhelm any achievable reduction in the implementation error. Thus, research that provides a better understanding of the geophysical principles can pay off handsomely when it comes to computer implementation. In short, better science leads to better computer usage.
Let us give some examples of shaping filters. Figures 5 through 10 contain four diagrams each. The top two diagrams show the input signal and the desired output signal, whereas the bottom two diagrams show the filter and the actual output signal. The filter is always a two-sided least-squares operator designed to act on the input signal to produce (in the least-squares sense) the desired output signal.
A two-sided filter with coefficients consists of two components. The causal component contains the coefficients for time indices . The anticausal component contains the coefficients for time indices . In other words, the left side of a two-sided filter represents the anticausal (or anticipation) part, whereas the right side (including the coefficient for time index 0) represents the causal (or memory) part. The actual output signal is the convolution of the filter with the input signal. How well the filter succeeds can be appreciated by comparing the actual output signal with the desired output signal.
The following entries must be provided: (1) the sampling rate (in milliseconds) of the required signal; (2) the length (in milliseconds) of the required signal; and (3) the time origin of the required signal. For the time-origin parameter, one of two choices usually is made: front (which positions time zero at the first sample point) or center (which positions time zero at the center sample point). In Figures 8, 9, and 10, we see that wrong choices of the time-origin parameter lead to poor results.
Figure 5 shows an example in which the intent is to shape a mixed-phase input wavelet into a more compressed mixed-phase desired output wavelet. Note that the anticipation part of the filter does most of the work to accomplish this task.
In Figure 6, the input wavelet is the minimum-phase equivalent of the input wavelet of Figure 5, and the desired output wavelet is the minimum-phase equivalent of the desired output wavelet in Figure 5. In theory, the memory part of the filter should do the work, whereas the anticipation part ideally should be zero. The reason the anticipation part is not strictly zero is because the data windows are of finite rather than infinite length.
Figure 7 illustrates the shaping of a wide symmetric wavelet into a narrower symmetric wavelet. We observe that both the filter and the actual output are also symmetric. Note the prominent unwanted side lobes on the actual output. These side lobes are the stiff price that is paid for the desire to reduce the width of the input signal.
Figure 8 illustrates a case for which neither the input nor the desired output was centered. In other words, for each of these signals, the center sample point was not placed at time zero. Shaping a delayed zero-phase wavelet into a delayed spike never should be attempted. Figure 8 shows the poor result. Instead, the time-origin parameter for both input and desired output should have been the center sample point of each signal.
It is not advisable either to shape a delayed zero-phase wavelet into a delayed narrow pulse. Such an attempt could make matters even worse. Figure 9 shows the poor result. Instead, the time-origin parameter for both input and desired output should have been the center sample point.
The approach shown in Figure 10 does not work well either. The figure shows what happens in an attempt to shape the delayed zero-phase wavelet into a delayed shorter-duration minimum-phase wavelet. Once more, the result is poor. Instead, the time-origin parameter for the input should have been the center sample point, and the time-origin parameter for the desired output should have been its front sampling point.
Let us describe the characteristics of the required input signal. The input signal required for the design of a shaping filter can be obtained in various ways. The input signal could be the recorded impulse response of the recording system or it could be the far-field signature or it could be a correlated vibroseis sweep signal. Alternatively, the input signal can be created in the computer by using its mathematical specifications.
Any signal can be altered in different ways. For example, the signal can be convolved with a theoretical geophone response or an inelastic absorption operator or the calculated impulse response of the recording filter or any combination of such components. Further processing is also possible. For example, the signal thus created can be band-pass-filtered, integrated, differentiated, deconvolved, or spectrally balanced, among other possibilities. Any combination of these possibilities could produce the required input signal.
An instrument response can be generated as a sequence of convolutions of different components, such as a low-cut filter, a notch filter, a high-cut filter, a ghost-removal filter, a geophone response, and an absorption response. The user must give individual specifications for these components. For example, for a low-cut filter component, one must enter the low-frequency roll-off point and the low-frequency slope. A notch-filter component could have a 50- or 60-Hz notch, for example. For a high-cut filter, one must enter the high-frequency roll-off point and the high-frequency slope. For either a marine-source ghost response or a marine-receiver ghost response, one must specify the two-way travel-time as well as the reflection coefficient of the water surface. For a calculated geophone response, one must specify the natural resonant frequency of the geophone string as well as its damping factor.
Next let us describe characteristics of the required shaping filter. After the required input signal has been either recorded or obtained, the characteristics of the shaping filter must be specified. The filter design can be used to determine a variety of shaping filters for the given input signal. Design options can include computation of inverse, band-pass, dephasing, minimum-phase, or spiking filters. An ideal inverse filter is the full amplitude and phase inverse of the given input signal. A band-pass filter can convert the input signal into a zero-phase band-pass signal or a minimum-phase band-pass signal.
A spiking filter converts the input wavelet into an approximation of a perfect spike signal. Wavelet spiking generally is implemented in the time domain; other options can be obtained conveniently in either domain. In the frequency domain, the minimum-phase wavelet spectrum can be found with the Hilbert transform. In particular, this minimum-phase spectrum can be computed as the Hilbert transform of the wavelet’s log-amplitude spectrum. In other words, the minimum-phase spectrum and the log-amplitude spectrum are a Hilbert transform pair, as shown in Silvia and Robinson (1979).
To design spiking filters, one usually must apply prewhitening, which involves a small percentage increase in the magnitude of the zero-lag wavelet autocorrelation coefficient. A prewhitening level of about 1% is usually sufficient. For higher signal-to-noise ratios, a lower prewhitening value, say 0.1%, often is satisfactory.
If we want to filter in the time domain, the minimum-phase counterpart of the given wavelet can be found by taking the least-squares inverse of the least-squares inverse of a given wavelet. In other words, this method finds the minimum-phase wavelet in two steps. First, the minimum-delay inverse filter is computed from the wavelet autocorrelation. It can be shown that this least-squares inverse filter is a minimum-delay filter (see for, example, Robinson, 1967a). Second, the minimum-phase inverse to the previously determined minimum-phase inverse filter is computed. This second minimum-phase inverse to the first minimum-phase inverse is the then required filter. One must specify the length (in milliseconds) of the least-squares inverse filter used for this minimum-delay calculation. One also must specify the prewhitening level — that is, the percentage boost of the zero-lag autocorrelation value used in the least-squares filter design.
On occasion, it might be desirable to apply a band-pass operator to the wavelet filter. In that case, the following information would have to be supplied to an appropriate program: (1) Should the operator be zero phase or minimum phase? and (2) Which are the four corner frequencies that specify the passband? Note that the four corner frequencies represent, respectively, the 0% and 100% points of the low-cut ramp and the 100% and 0% points of the high-cut ramp, in hertz. Such ramps could be formed by Hanning (cosine) tapers in the frequency domain (Robinson and Silvia, 1978). For example, the entry: 4-12-50-75 would generate a filter with a 12- to 50-Hz passband, along with an 8-Hz-wide low-cut ramp and a 25-Hz-wide high-cut ramp.
Often, we must apply a symmetric taper to the output filter to reduce spurious edge effects that can occur with an optimal finite-length filter. Taper choices include (1) a Hanning (cosine) ramp or (2) a Bartlett (linear) ramp (Robinson and Silvia, 1978). The length (in milliseconds) of the taper window as well as the percentage of the calculated filter to be held flat must be given. For example, if this percentage is set to 80%, then the filter will have to be tapered over the first 10% and the last 10% of its length. A value of 100% is equivalent to no taper — that is, it is a boxcar truncation.
Once the filter is designed, it is convolved with the input to yield the actual output. Throughout the seismic section, continual examination of the actual outputs with the desired outputs provides a means for monitoring the design process in action and for evaluating changes in the parameters. Such a diagnostic option would be helpful for determining filter effectiveness as well as for selecting the best possible filter parameters.
- ↑ Silvia, M. T., and E. A. Robinson, 1979, Deconvolution of geophysical time series in the exploration for oil and natural gas: Elsevier.
- ↑ Robinson, E. A., 1967a, Statistical communication and detection: Hafner Publishing Co.
- ↑ Robinson, E. A., and M. T. Silvia, 1978, Digital signal processing and time series analysis: Holden-Day.
|Previous section||Next section|
|White convolutional model||All-pass filter|
|Previous chapter||Next chapter|
Also in this chapter
- The shaping filter
- Spiking filter
- White convolutional model
- All-pass filter
- Convolutional model
- Nonminimum-delay wavelet
- Signature deconvolution
- Dual-sensor wavelet estimation
- Deconvolution: Einstein or predictive?
- Appendix I: Exercises