
Research Article


Increasing the Robustness of Homomorphic Deconvolution for Nonstationary Seismic Wavelet 

Muhammad Sajid
and
D.P. Ghosh



ABSTRACT

Study presents a novel algorithm to improve the frequency bandwidth of the
seismic data through removal of nonstationary wavelet in cepstrum domain. Algorithm
incorporates useful properties of both homomorphic deconvolution and spectral
decomposition by this way it estimates and filters the wavelet at each translation
of spectral decomposing Gaussian window which improve the resolution and tolerates
the nonstationary properties of seismic wavelet. This study describes the comprehensive
mathematical formulation of the algorithm and its testing on 1D, 2D synthetic
and real seismic section to confirm its applicability. Comparative testing shows
that the algorithm effectively removed the smearing effect of wavelet which
led to broader frequency spectrum and improved seismic resolution.





Received: September 03, 2013;
Accepted: November 08, 2013;
Published: April 10, 2014


INTRODUCTION Seismic data is considered to be a convolution of the nonstationary wavelet with the subsurface geology. Seismic wavelet changes its shape and frequency component continuously. Rather than the Fourier transform, Spectral decomposition contains both time and frequency information at the same instant and it is susceptible to local geology (Partyka et al., 1999). Its spectral decomposing window can be of different shapes and size. ShortTimeFourierTransform (STFT) with Gaussian analysis window is known as Gabor transform (Gabor, 1946). Cepstrum is defined as Inverse Fourier transform of the log of the magnitude of Fourier transform. It is mostly used in speech analysis and echo elimination (Oppenheim et al., 1968; Oppenheim, 1965). In seismic, it is being used for homomorphic deconvolution (Ulrych, 1971; Stoffa et al., 1974; Buhl et al., 1974; Tribolet et al., 1977) and for calculation of thin bed thickness with better accuracy (Hall, 2006; Liu et al., 2008).
We developed a new algorithm which we named as CepstrumWaveletFiltering (CWF)
which estimate and filter out the nonstationary wavelet at each translation
of the spectral decomposing window in the cepstrum domain by this way it removes
the stationary wavelet limitation of the homomorphic deconvolution algorithm
by incorporation of spectral decomposition algorithm in it. It is simple, fast
and perfectly reconstructs the seismic trace back from the cepstrum domain through
preservation the phase information in its workflow.
MATERIALS AND METHODS
As ShortTimeFourierTransform (STFT) implement the Fourier transform on each windowed version of the signal (Eq. 1) (Allen and Rabiner, 1977; Auger and Flandrin, 1995), so application of natural logarithm on each translation Fourier coefficient provides the logarithm of the magnitude of Fourier coefficients as a real part whereas the imaginary part contains the phase spectrum of the given translation as described in Eq. 2. Inverse Fourier transform of the real part provides the real cepstrum (Eq. 3) which is being used to recover the reflectivity by filtering out the wavelet whereas the phase spectrum of each translation is preserved for later perfect reconstruction of signal from the cepstrum domain as described in Eq. 4.
where, t is Time, f is frequency, x (t) is time domain signal, h (t) is decomposing window, τ is window translation along time axis, C is real cepstrongram, γ (t) is spectral reconstructiol window, n is total number of sample in trace, b is current location of window and w is half width of gaussian window. Gaussian window (Eq. 5) flanks are used to smoothly filter out the initial quefrencies from each translation (τ) of real cepstrogram. Figure 1 shows a graphical illustration of Gaussian filtering operation on 1D synthetic thin bed model (Fig. 1b), created by convolution of Ricker wavelet (f_{dom} = 35Hz) with reflectivity series (Fig. 1a), in a simplified form without using STFT in the algorithm. Synthetic synthetic trace comprises a set of 6 thin beds models for the following conditions:
•  Single event wavelet (f_{dom} = 35 Hz and Tuning thickness = 12 msec) 
• 
Two events are in between 0 to flat spot thickness (i.e., 3 samples or 6 msec thickness) 
• 
Two events are below flat spot thickness (i.e., 4 samples or 8 msec thickness) 
• 
Two events are at Flat spot thickness (i.e., 5 samples or 10 msec thickness) (Ricker’s Criterion) (Ricker, 1953) 
• 
Tuning thickness (i.e., 6 samples or 12 msec thickness) (Rayliegh’s Criterion) (Kallweit and Wood, 1982) 
• 
Greater than Tuning thickness (i.e., 7 samples or 14 msec thickness) 

Fig. 1(ag): 
(a) Thin bed reflectivity model, (b) Synthetic seismic trace created by convolution with Ricker wavelet f_{dom }= 35 Hz, (c) Real cepstrum of the synthetic trace, (d) Phase spectrum of the trace, (e) Initial quefrencies are zeroed out by point wise multiplication of the smooth negative of Gaussian window flank, (f) Real cepstrum after filtering and (g) Improved resolution after the application of the CWF algorithm (reconstructed black trace) as compared to original (red curve) 
Figure 1c and d shows the real cepstrum and phase spectrum respectively whereas Fig. 1e shows the application of Gaussian filtering operation on real cepstrum by point wise multiplication of negative of the Gaussian window flanks with the initial quefrencies of the real cepstrum which smoothly reduces its values whereas Fig. 1f shows real cepstrum after filtering. Figure 1g shows the comparison of original synthetic trace (red curve) with CWF algorithm reconstructed trace which shows improved resolution.
APPLICATION Application on thin bed model: Figure 2 shows the CWF algorithm implementation on synthetic thin bed model which is created by using the above mention parameters. Figure 2b shows its real cepstrogram whereas Fig. 2c shows the Gaussian filtering window for each translation of cepstrogram which is used to reduce initial quefrencies. Figure 2d and e shows the filtered real cepstrogram and its reconstructed seismic trace respectively. By comparison of the original synthetic trace (Fig. 2a) with CWF reconstructed trace, it can be observe that the resolution of events is effectively improved.

Fig. 2(ae): 
(a) Original synthetic thin bed model, (b) Real cepstrogram, obtained by calculation of real cepstrum at each time step of spectrogram, (c) Gaussian filtering window for each time step of cepstrogram, (d) Filtered cepstrogram after implementation of the Gaussian window filtering and (e) Reconstructed synthetic trace from filtered cepstrogram which shows improved resolution 
Application on synthetic wedge model: The 2D synthetic model is created from the geological model of the stratigraphic layer of known thickness and elastic properties. In the synthetic seismic (Fig. 3a), trace interval is 2 m, sample rate is 2 msec, the wedge layer velocity is 2700 m sec^{1} and Ricker wavelet with the predominant frequency of 35 Hz:
is used.
As CWF algorithm filter out the wavelet at each translation (τ) without being the influence of user predefined spectral broadening window (Yilmaz, 2008) so the low frequencies of the synthetic seismic section remains whereas the higher frequencies become apparent because of the removal of the wavelet interference effect as can be observed in Fig. 3a. In comparison of improvement in seismic resolution, it can be observed that resolution of the wedge model is effectively improved after the application of CWF as can be observe through Fig. 3b and c. An event which started to merge at trace No. 22 is separable till trace No. 32 which shows a 10 m increase in resolution as highlighted by colour lines.
Application on real seismic data: Seismic data from Malay Basin are usually of high quality with good signal to noise ratio.

Fig. 3(ac): 
(a) Amplitude spectrum of synthetic wedge model with almost same trend but with higher values, (b) Original synthetic wedge model with S/N ratio of 3 and trace interval of 2 m and (c) Improved resolution, events which started to merge at trace No. 22 still separable till trace No. 30 (8 m improve in resolution) 

Fig. 4(ab): 
(a) Real seismic section contains a number of hidden features and (b) Improved seismic resolution and number of hidden features are recovered after the application of CWF algorithm, some are marked with arrow signs 

Fig. 5(ab): 
(a) Real seismic from Malay Basin show gas seepage and ( b) CWF algorithm improved resolution of seismic section without being much adversely affected by low frequency zones 
Figure 4a shows the real seismic section which contains number of hidden features due to the interference effect of the seismic wavelet. Figure 4b shows the same seismic section after the application of propose algorithm, it can be observe that the CWF algorithm extracted a number of hidden features as some of them are marked with arrow sign and effectively improved the seismic resolution. In another example shown in Fig. 5, shows the application of CWF algorithm on seismic section containing gas clouds. As gas clouds produce large high frequency attenuation which leads to shadowing effect to underlying layers. Figure 5a and b show the original seismic and seismic after the proposed CWF algorithm. It can be easily observed that the proposed algorithm has considerably improved the resolution of seismic without being much adversely affected by gas cloud in original data.
CONCLUSION New proposed CepstrumWaveletFiltering (CWF) algorithm increase the robustness of homomorphic deconvolution by introduction of ShortTimeFourierTransform (STFT) in its algorithm. This incorporation leads to estimate and filter the seismic wavelet at each translation of spectral decomposing Gaussian window which improve the resolution of seismic data. This improvement in seismic resolution which leads to better thickness estimation, correlation to well logs and thereby improves the stratigraphic interpretation. These results are achieved without being constrained by well log information and starting model, which leads to unbiased data driven better resolution. As evident from its testing on 1D, 2D synthetic and on real seismic section, the proposed algorithm has potential use in seismic interpretive data processing. ACKNOWLEDGMENTS We thank Universiti Teknologi PETRONAS (UTP) providing the resources for research and Petronas for support and permission to show seismic data of Malay Basin. M. Sajid thanks University COMSATS for providing the opportunity to carry out research at UTP. We also thank Arthur E. Barnes from Petronas for valuable discussion and suggestions. We used the CREWES free software package for the synthetic seismic generation.

REFERENCES 
1: Allen, J.B. and L. Rabiner, 1977. A unified approach to shorttime Fourier analysis and synthesis. Proc. IEEE, 65: 15581564. CrossRef 
2: Auger, F. and P. Flandrin, 1995. Improving the readability of timefrequency and timescale representations by the reassignment method. IEEE Trans. Signal Process., 43: 10681089. CrossRef 
3: Buhl, P., P.L. Stoffa and G.M. Bryan, 1974. The application of homomorphic deconvolution to shallowwater marine seismologyPart II: Real data. Geophysics, 39: 417426. CrossRef 
4: Gabor, D., 1946. Theory of communication. Part 1: The analysis of information. Electr. Eng. Part III: Radio Communi. Eng. J. Inst., 93: 429441. Direct Link 
5: Hall, M., 2006. Predicting bed thickness with cepstral decomposition. Leading Edge, 25: 199204. CrossRef 
6: Kallweit, R.S. and L.C. Wood, 1982. The limits of resolution of zerophase wavelets. Geophysics, 47: 10351046. CrossRef 
7: Liu, C.Y., S.F. Xu, J.H. Yue and X.C. Wei, 2008. Application of cepstral analysis to thin beds. J. China Univ. Mining Technol., 18: 232236. CrossRef 
8: Oppenheim, A.V., R.W. Schafer and T.G. Stockham Jr., 1968. Nonlinear filtering of multiplied and convolved signals. IEEE Trans. Audio Electroacoustics, 16: 437466. CrossRef 
9: Oppenheim, A.V., 1965. Superposition in a class of nonlinear systems. Ph.D. Thesis, Research Laboratory of Electronics, MIT, Cambridge.
10: Partyka, G., J. Gridley and J. Lopez, 1999. Interpretational applications of spectral decomposition in reservoir characterization. Leading Edge, 18: 353360. CrossRef 
11: Ricker, N., 1953. Wavelet contraction, wavelet expansion and the control of seismic resolution. Geophysics, 18: 769792. CrossRef  Direct Link 
12: Stoffa, P., P. Buhl and G. Bryan, 1974. The application of homomorphic deconvolution to shallowwater marine seismologyPart I: Models. Geophysics, 39: 401416. CrossRef 
13: Tribolet, J., T.F. Quatieri and A.V. Oppenheim, 1977. Shorttime homomorphic analysis. Proceedings of the IEEE International Conference on Acoustics, Speech and Signal Processing, Volume 2, May 911, 1977, Connecticut, USA., pp: 716722 CrossRef 
14: Ulrych, T., 1971. Application of homomorphic deconvolution to seismology. Geophysics, 36: 650660. CrossRef 
15: Yilmaz, O., 2008. Seismic Data Analysis. Vol. 1, SEG Publishing Co., Tuls, pp: 231233



