# The Leading Edge

- © 2014 by The Society of Exploration Geophysicists

## Abstract

Statistical white reflectivity, minimum-phase wavelets, and stationarity are regarded to be three major shortcomings in dealing with the Robinson convolution model. Modern reflectivity inversion methods (e.g., constrained reflectivity inversion) mainly attempt to reduce side effects of the first two assumptions, but most of them overlook the important fact that seismic signals are typically nonstationary. Numerical tests and practical applications have verified that nonstationarity does have great influence on reflectivity estimation and probably would result in misplaced positions or wrongly restored amplitudes for estimated reflectivity. Nonstationarity can be tackled from a perspective of hyperbolic smoothing in Gabor deconvolution. Specially, the energy relationship of hyperbolic strips in a log spectrum can be taken as a quantitative indicator in balancing nonstationarity and conditioning seismic traces to the assumption of unchanging wavelets. Applications to marine seismic data show that balancing nonstationarity finally helps to recover subtle reflectivity information and promotes better displays of geologic features of the subsurface.

Statistical white reflectivity, minimum-phase wavelets, and stationarity are regarded to be three major shortcomings in dealing with the Robinson convolution model. Modern reflectivity inversion methods (e.g., constrained reflectivity inversion) mainly attempt to reduce side effects of the first two assumptions, but most of them overlook the important fact that seismic signals are typically nonstationary. Numerical tests and practical applications have verified that nonstationarity does have great influence on reflectivity estimation and probably would result in misplaced positions or wrongly restored amplitudes for estimated reflectivity. Nonstationarity can be tackled from a perspective of hyperbolic smoothing in Gabor deconvolution. Specially, the energy relationship of hyperbolic strips in a log spectrum can be taken as a quantitative indicator in balancing nonstationarity and conditioning seismic traces to the assumption of unchanging wavelets. Applications to marine seismic data show that balancing nonstationarity finally helps to recover subtle reflectivity information and promotes better displays of geologic features of the subsurface.

## Introduction

It is a common practice that the assumptions of minimum-phase wavelet and statistical white reflectivity are used to enable feasible solutions for the Robinson convolution model (Robinson, 1954). However, it has become widely accepted that wavelets are often of mixed phase, whereas reflectivity seldom behaves according to the Gaussian distribution (Walden and Hosken, 1986). For this reason, blind deconvolution was proposed, which is implemented as a constrained inversion problem by using more practical assumptions on reflectivity.

In spite of this, most methods are still established on the basis of stationarity, in which a wavelet is always assumed to be invariant. However, it seldom makes sense in practice; the wavelet always changes during its propagation as a result of many types of physical effects of earth filtering. In other words, nonstationarity is a common phenomenon that shakes the basis of most reflectivity inversion methods.

As a significant improvement on the Robinson convolution model, the nonstationary convolution model was proposed in the late 1990s and forms the basis of Gabor deconvolution. This model allows the wavelet to evolve with time (Margrave, 1998; Margrave and Lamoureux, 2001) by incorporating attenuation effects, which is always regarded as a major contributor to nonstationarity. Therefore, Gabor deconvolution not only can achieve better resolution than stationary deconvolution types but also is more capable of addressing issues of nonstationarity than either time-varying deconvolution or inverse-*Q* filtering. That is mainly because time-varying deconvolution is sensitive to empirical parameters (e.g., window size) (Griffiths et al., 1977; Koehler and Taner, 1985), and inverse-*Q* filtering might suffer from instability (Bickel and Natarajan, 1985; Hargreaves and Calvert, 1991; Wang, 2002).

To begin with, influences of nonstationarity are first examined with actual examples in this paper. More important, we attempt to seek a feasible way to quantify and balance nonstationarity from the perspective of hyperbolic smoothing that is the kernel of Gabor deconvolution.

## Constrained reflectivity inversion

It is widely accepted that actual reflectivity sequences seldom are different from white-noise sequences. Their stochastic behaviors are not white or Gaussian noise (Canadas, 2002). Observations in various geographic zones have verified that their distributions were symmetrical but had a narrow central peak and tail that decayed much more slowly than Gaussian distribution (Walden and Hosken, 1986).

For this reason, the Cauchy constraint is used frequently in reflectivity inversion because it generalizes the stochastic properties of real reflectivity sequences well. In this case, the objective function for reflectivity inversion can be established as (Sacchi, 1997; Canadas, 2002)
(1)
where ** r** and

**are vectors of reflectivity sequence and seismic trace, respectively, and**

*d***represents the wavelet matrix for the convolution process. Furthermore, impedance constraint also is incorporated to ensure a more reliable estimation, i.e., (2) Here,**

*W***stands for the causal matrix,**

*C***is the vector of the acoustic-impedance trend, and**

*p**μ*and

*β*are trade-off parameters for the sparsity constraint and impedance constraint, respectively. By equating the gradient of equation 2 to zero, the normal equation for the constrained reflectivity inversion finally formulates as (Velis, 2008) (3) for which the iterative scheme usually is used as an efficient solution (Sacchi, 1997).

To begin with, we would like to examine the performance of the constrained reflectivity inversion with an actual reflectivity sequence (Figure 1b) from a carbonate zone. A synthetic trace is generated by convolving the reflectivity sequence with a minimum-phase wavelet, the dominant frequency of which falls at 30 Hz. Then, as shown in Figure 1a, this trace is contaminated by 5% white noise.

We attempt to simulate a stationary situation here because the wavelet is invariant from top to bottom. The initial guess on reflectivity and wavelet is provided by spiking deconvolution. Processing parameters finally are determined as *σ** _{r}* = 0.054,

*μ*= 0.005,

*β*= 0.03. The iteration process is terminated when the relative changing ratio of reflectivity is smaller than 10

^{–6}.

Figure 1c shows estimated reflectivity. For a better view, we show only the most important contents of reflectivity, amplitudes of which should be no less than 30% of the average. Recovered reflectivity (Figure 1c) is almost identical to actual reflectivity (Figure 1b), implying robustness of this type of constrained inversion algorithm under stationary situations.

## Nonstationarity and hyperbolic smoothing in Gabor deconvolution

Stationary situations seldom exist in the real world. Because of earth filtering, a wavelet changes all the time during its propagation in the subsurface. In fact, it is a common phenomenon that challenges the basic assumption of the unchanging wavelet for many inversion methods built on the traditional convolution model (Robinson, 1954), which fails to accommodate this type of physical effect. Nonstationarity might deflect us from the true answer if its influences are overlooked. In this paper, the solution to this problem begins with the hyperbolic smoothing in Gabor deconvolution.

In the Gabor-transformed domain, an actual trace should be modeled more reasonably as (Margrave, 1998; Margrave and Lamoureux, 2001)
(4)
where *ŵ*(*f*) is the time-independent frequency spectrum of the source wavelet and *α*(*t*, *f*) and *Ȓ*_{g}(*t*, *f*) stand for the time-frequency spectrum of attenuation and reflectivity, respectively. As a common effect in subsurface media, attenuation is always regarded as a major contributor to nonstationarity in seismic traces. Furthermore, the nonstationary convolution model in equation 4 would be changed into a linear one by logarithmic transformation, i.e.,
(5)
where *L* stands for logarithm operation on its subscript. This equation casts new light on the factorization of the nonstationary convolution model because its gradient is related closely to *Q*, whereas the intercept consists of source wavelet and reflectivity.

Hyperbolic smoothing (Margrave et al., 2011) usually is recommended to ensure better factorization. In this process, the *t-f* log spectrum needs to be divided into hyperbolic strips, and an iterative process proceeds to seek better estimationsof attenuation by suppressing the influence of the source wavelet until the energy of each hyperbolic stripe is almost balanced.

The main advantages of hyperbolic smoothing in the log spectrum are its accuracy and efficiency because separation in the log spectrum is done largely by subtractions instead of divisions and thus no longer needs prewhitening factors to ensure numerical stability (Sun et al., 2012; Sun et al., 2013). In most cases, only a few iterations are needed to complete the processing.

Here, more tests are conducted on the actual reflectivity sequence in Figure 1b to show a feasible way of balancing nonstationarity, which is extracted from the hyperbolic smoothing in the log spectrum. The magnitude spectrum of the actual reflectivity sequence in Figure 2a shows that its energy tends to concentrate in the high-frequency zone.

This special type of spectrum is called blue in reference to the visible spectrum, which is a common property for practical reflectivity sequences. If we investigate rms values of the 20 hyperbolic strips in the log spectrum, we might find that the resultant energy curve (i.e., the black curve in Figure 2d) is nearly flat under stationary situations, despite the fact that it is a blue spectrum instead of a white one.

This time, we simulate a nonstationarity situation by attenuating the synthetic trace with a three-layered *Q* model (right side of Figure 2b). Figure 2b shows the corresponding magnitude spectrum, in which the influence of the source wavelet is removed to make an analog with Figure 2a. At this stage, the energy curve in the log spectrum (i.e., the blue curve in Figure 2d) is no longer flat but exhibits an obvious gradient, which we suggest indicates the existence of nonstationarity.

Attenuation, source wavelet, and reflectivity can be estimated successively during hyperbolic smoothing. As a result, the magnitude spectrum of reflectivity in Figure 2c is basically restored. Meanwhile, the energy curve (i.e., the red curve in Figure 2d) is almost corrected to be as horizontal as the black curve.

If the constrained reflectivity inversion is applied directly to the attenuated trace without prior investigation of its nonstationarity, it might lead to inaccurate estimations on reflectivity, e.g., misplaced positions or wrongly restored amplitudes. As shown in Figure 3b, deep features become blurred or lost because of ignorance of nonstationarity. In contrast, reflectivity finally is well recovered after balancing nonstationary first.

During the whole process, Gabor deconvolution plays a role in preconditioning seismic traces to be modeled easily by unchanging wavelets, which is a basic assumption for constrained inversion. Moreover, because this method depends largely on analyzing the energy relationship in the log spectrum, it addresses practical issues without requiring knowledge of *Q* values beforehand.

## A field example

In this example, we will mainly discuss issues of nonstationarity in migrated marine data. Taking one seismic trace as an example, dominant frequency of the shallow half is higher and the frequency band is broader than that of the deep half (Figure 4a).

Furthermore, as shown in Figure 4b, its energy curve consisting of rms values of 100 hyperbolic strips in the log spectrum is a decreasing curve, which indicates that the seismic trace is under the influence of nonstationarity. The iteration process of hyperbolic smoothing in Gabor deconvolution can work to correct and balance the energy differences on this curve. Therefore, the deep spectrum (Figure 4c) is boosted, and its attenuation curve (Figure 4d) also is well corrected after application of Gabor deconvolution.

In practical processing, nonstationarity is a most overlooked issue, but it can impact the ultimate estimations of reflectivity greatly. Figure 5a and Figure 5b show the results obtained by implementing conventional constrained reflectivity inversion. Results here are displayed according to positive and negative reflectivity. However, quality of obtained reflectivity is improved further if nonstationarity is addressed and balanced before implementation of reflectivity inversion, with results shown in Figure 5c and Figure 5d.

In a routine manner, we directly applied constrained reflectivity inversion in our study area and then convolved the obtained reflectivity with a Ricker wavelet of a 60-Hz dominant frequency. As a result, Figure 6a presents an extracted amplitude slice at 2400 ms. In contrast, if we conduct Gabor deconvolution prior to the reflectivity inversion and do the same convolution operation, the corresponding amplitude slice is as Figure 6b shows.

It is found that many ambiguities shown in Figure 6a finally obtain clear characterizations and promote better definitions after nonstationarity is corrected. This finding reminds us that a certain amount of useful information might be veiled because of the influences of nonstationarity if we apply reflectivity inversion blindly without careful investigation.

## Summary

Nonstationarity is such a common phenomenon in the real world that it somewhat questions the validity of unchanging wavelets and thus poses great challenges for many techniques built on this model. Without prior investigation on nonstationarity, direct applications of constrained reflectivity inversion could lead to misplaced positions or wrongly restored amplitudes for estimated reflectivity.

In this paper, a feasible method that addresses nonstationarity is proposed from implementation of Gabor deconvolution. In this method, the energy curve consisting of the rms value of hyperbolic strips in the log spectrum can serve as a quantitative measurement of nonstationarity. Under this guidance, nonstationarity can be regarded as balanced until the attenuation curve is corrected to be almost flat. Balancing nonstationarity conditions seismic traces to be modeled readily by unchanging wavelets.

Joint applications with constrained reflectivity inversion can promote better characterizations for geologic structures. Otherwise, a certain amount of useful information probably would be veiled if the influence of nonstationarity is overlooked, as it is by many routine techniques.

## Acknowledgments

This research is funded partly by the National Basic Research Program of China (Grant No. 2011CB201103). The authors would like to thank the Laboratory for Integration of Geology and Geophysics at China University of Petroleum (Beijing) for permission to publish this work.