- Email: [email protected]
Using a time-frequency distribution to identify buried channels in reflection seismic data
Accepted Manuscript Using a time-frequency distribution to identify buried channels in reflection seismic data
Azita Nikoo, Amin Roshandel Kahoo, Hamid Hassanpour, Hamed Saadatnia
PII: DOI: Reference:
S1051-2004(16)30002-1 http://dx.doi.org/10.1016/j.dsp.2016.03.008 YDSPR 1928
To appear in:
Digital Signal Processing
Please cite this article in press as: A. Nikoo et al., Using a time-frequency distribution to identify buried channels in reflection seismic data, Digit. Signal Process. (2016), http://dx.doi.org/10.1016/j.dsp.2016.03.008
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Using a Time-Frequency Distribution to Identify Buried Channels in Reflection Seismic Data 1 Azita Nikoo a, േ , Amin Roshandel Kahoo b, Hamid Hassanpour c, Hamed Saadatnia d a
Department of Mining, Petroleum and Geophysics Engineering, University of Shahrood, Iran
b
Assistant Professor, Department of Mining, Petroleum and Geophysics Engineering, University of
Shahrood, Iran c
Professor, Department of Computer Engineering & IT, University of Shahrood, Iran
േ
d
Exploration Directorate, NIOC, Iran
Corresponding author. E-mail addresses: [email protected] (A. Nikoo), [email protected] (A. Roshandel Kahoo), [email protected] (H. Hassanpour), [email protected] (H. Saadatnia).
Abstract Geological events such as thin beds and channels cannot be easily revealed on seismic sections due to the interference of reflections from the top and bottom of the layer. Buried channel is one of the hydrocarbon traps, which is important in oil and gas exploration. Spectral decomposition can be used to indicate subtle changes in channel thickness. In using the Fourier transform only the frequency content of the changes is displayed; hence, the exact onset of the changes is missed. Time-frequency distributions are suitable approaches to display and interpret information embedded in non-stationary signals such as seismic signals. Spectral decomposition can be used for seismic attributes calculation, which are used for imaging of thin beds. The conventional spectral decompositions such as Short Time Fourier Transform (STFT) and Wigner-Ville Distribution (WVD) have limitations in terms of Heisenberg uncertainty principle and cross-terms artefacts, respectively. In this paper, we used the Reduced Interference Distribution (RID) for buried channel identification to overcome the mentioned limitations. We compared the obtained results of RID with that of Smoothed Pseudo WVD. The real seismic data is selected from one of the western oil fields of Iran. Our results indicate that a better resolution is achieved by RID in both vertical and lateral stratigraphic features. Keywords: Spectral decomposition, Seismic attributes, Wigner–Ville distribution, Reduced interference distribution.
1. Introduction Seismic waves are mechanical waves that travel through the Earth’s layers or along its surface. Seismic signals propagate into the Earth and reflect from boundary separating rocks with different properties [1]. Indeed, seismic reflection records carry important information, such as amplitude and frequency, associated with the Earth and its subsurface structure. Seismic waves are non-stationary because the Earth acts as a low-pass filter and changes the frequency content of seismic signals during the propagation [2]. Ordinary 3-D seismic data contains energy reflected from the subsurface at a wide range of frequencies, all of which are compounded in a typical seismic volume. Different frequencies of seismic reflection data can image geological features with frequency-dependent nature. Geological interpretations from seismic reflection data are mostly drawn from the time-amplitude domain. It is used to improve and increase the application of seismic data sets for geological interpretations [3]. Spectral decomposition can be helpful to see the amplitude or phase of reflections at particular frequencies. Besides, the spectral decomposition analysis allows the geophysicist to quantify amplitude variation with frequency and extract geological information from seismic data associated with frequency that normally cannot be viewed in conventional seismic data [4]. Fourier transform is not a suitable tool to analyze non-stationary signals. It cannot show the frequency variations of seismic signals over time. Time-frequency analysis is more appropriate for seismic data due to their non-stationary nature. Spectral decomposition has been widely used in seismic data processing and interpretation. We can explore stratigraphic settings, thin beds and buried channels using spectral decomposition of seismic sections [4-7]. The spectral decomposition can be employed to detect lowfrequency shadows associated with hydrocarbons [8-10]. Sinha et al. in [11] employed the continuous wavelet transform to computing a time-frequency map for seismic signal. Time-frequency map is used for seismic attributes computation. Askari and Siahkoohi in [12] attenuated the ground roll noise (a common type of coherent noise in seismic data) from seismic data using the S and x-f-k transforms. Roshandel Kahoo and Siahkoohi in [13] performed amplitude versus offset (AVO) analysis in the time-frequency
domain for detection of the gas reservoir. They used modified high resolution time–frequency transforms to suppress random noise from seismic data using peak filtering method [14]. The spectral decomposition can be used for time varying deconvolution [15, 16]. There are several methods in literature for spectral decomposition, such as short time Fourier transform (STFT) [17, 18], Wigner–Ville distribution (WVD) [19, 20], wavelet transform [21] and S–transform [22, 23]. Conventional spectral decomposition methods have some limitations, such as the Heisenberg uncertainty principle and cross–terms existence, which limit their applications in analyzing non-stationary signals. The conventional method of making a time-frequency map using the STFT limits time-frequency resolution by a predefined window length [11]. Therefore, the valuable information is lost or blurred in the time–frequency representations due to a trade-off between the time and frequency resolution. No window function is used to calculate the WVD of a signal. Therefore, WVD has high resolution both in time and frequency domain. However, the applications of this distribution are limited by cross–terms existence. Several methods such as Pseudo WVD (PWVD) and Smoothed Pseudo WVD (SPWVD) have been introduced to eliminate the cross–terms [24-27]. Smoothed versions of WVD use a time and frequency dependent function as a kernel to reduce cross terms. However, these methods reduce the energy distribution concentration in time-frequency plane [14, 28-30]. Cone shaped distribution (CSD) and adaptive spectrogram (ADS) are used for reducing cross-terms in biological signal [31]. In this paper, we used Reduced Interference Distribution (RID) [32] to calculate the time–frequency distribution which is found to be more appropriate for this application. Seismic signals are broadband (8-80 Hz) FM signals and RID is a suitable choice here. This method of distribution with slight impact on the time–frequency resolution increasingly reduces the cross terms. Analysis of thin–bed tuning on seismic reflectivity has been studied extensively [33-37]. During the last decade, spectral decomposition technique has proven to be an excellent tool to describe thin beds associated with channel sands, alluvial fans, and other similar structures [38]. In this paper, we focused on extracting seismic attributes based on the RID such as instantaneous center frequency, instantaneous
bandwidth frequency, root mean square (RMS) frequency, peak frequency and peak amplitude. These attributes are directly related to the thickness of the thin-bed. Since a buried channel is like a thin-bed, we used the extracted attributes for detection of buried channels.
2. Methodology In this section first we explained one of the effective time-frequency distributions and also a way in which to improve them. We then illustrated the nature of seismic data and the characteristics of thin layer on seismic sections, and proposed some useful seismic attributes to reveal buried channels.
2.1 Time-Frequency Distribution Resolution of time and frequency in STFT is limited by the Heisenberg uncertainty principle. One of the advantages of WVD over the STFT is to distribute energy density in time–frequency plane without losing time or frequency resolution. The WVD of signalݔሺݐሻis defined as [39]: ାஶ
ఛ
ఛ
ଶ
ଶ
ܹܸܦሺݐǡ ݂ሻ ൌ ିஶ ܺ ቀ ݐ ቁ ܺ כቀ ݐെ ቁ ݁ ିଶగఛ ݀߬ሺͳሻ where ܺሺݐሻ is the analytic signal of signal ݔሺݐሻ, ܺ כሺݐሻ denotes the complex conjugate of ܺሺݐሻ, ݂ is frequency, ݐis time and ߬ is the time delay. The analytic signal is a complex signal [40] and can be computed as: ܺሺݐሻ ൌ ݔሺݐሻ ݆ ܪሾݔሺݐሻሿሺʹሻ where ܪሾݔሺݐሻሿ is the Hilbert transform of ݔሺݐሻ and ݆ ൌ ξെͳ. ܭ ሺݐǡ ߬ሻ is defined as the instantaneous autocorrelation function of ܺሺݐሻ: ߬ ߬ ܭ ሺݐǡ ߬ሻ ൌ ܺ ቀ ݐ ቁ ܺ כቀ ݐെ ቁ ሺ͵ሻ ʹ ʹ By substitution of Eq. (3) into Eq. (1), we can rewrite Eq. (1) as:
ାஶ
ܹܸܦሺݐǡ ݂ሻ ൌ ିஶ ܭ ሺݐǡ ߬ሻ݁ ିଶగఛ ݀߬ ሺͶሻ Eq. (4) shows that WVD is equal to the Fourier transform of the instantaneous autocorrelation function of ܺሺݐሻ with respect to ߬. Because of bilinear characteristic, cross-terms interference is emerged for multi-component signals. For a signal composed of two componentsݔሺݐሻ ൌ ݔଵ ሺݐሻ ݔଶ ሺݐሻ, based on the quadratic superposition principle [20, 41], the WVD is defined as: ܹܸܦ௫ ሺݐǡ ݂ሻ ൌ ܹܸܦ ௫భ ሺݐǡ ݂ሻ ܹܸܦ௫మ ሺݐǡ ݂ሻ ܹܸܦ ௫భ ǡ௫మ ሺݐǡ ݂ሻ ܹܸܦ௫మ ǡ௫భ ሺݐǡ ݂ሻ ሺͷሻ ᇣᇧᇧᇧᇧᇧᇧᇧᇧᇤᇧᇧᇧᇧᇧᇧᇧᇧᇥ ᇣᇧᇧᇧᇧᇧᇧᇧᇧᇧᇧᇤᇧᇧᇧᇧᇧᇧᇧᇧᇧᇧᇥ Ǧ
Ǧ
The two last terms in right hand of Eq. (5) are cross–terms with the following characteristics [20, 41]: 1- They are found midway between the interacting components; 2- They oscillate proportionally to the inter-components' time-frequency distance; 3- They have a direction of oscillation orthogonal to the straight line connecting the components. PWVD and SPWVD are two well-known distributions with less cross-term compared with WVD. In these distributions, the cross-terms are reduced by smoothing the time-frequency representation due to their oscillation property. The smoothing reduces the cross-terms, but with the cost of reduction of resolution in time–frequency representation. Here we exploit the RID method for reducing cross–terms which has the least impact on the resolution of the time–frequency representation. This distribution belongs to the general Cohen class of distributions. Cohen’s distribution of a signal ݔሺݐሻ is defined as [42]: ାஶ ାஶ ାஶ
߬ ߬ ܴܶܨ௫ ሺݐǡ ݂ሻ ൌ න න න ݁ ଶగ௩ሺ௧ି௨ሻ ݃ሺݒǡ ߬ሻܺ ቀ ݐ ቁ ܺ כቀ ݐെ ቁ ݁ ିଶగఛ ݀߬݀ݑ݀ݒሺሻ ʹ ʹ ିஶ ିஶ ିஶ
where ݃ሺݒǡ ߬ሻis a 2D kernel and specifies the distribution. The kernel of RID emphasizes the auto-terms and deemphasizes the cross-terms. The RID keeps almost all the desirable properties of the WVD, but
with considerably reduced interference. In this paper, we used Binomial kernel, which is a particularly attractive discrete form which can be computed very efficiently [32]. To study the efficiency of the WVD, PWVD, SPWVD and RID, we tested them on a synthetic signal. The synthetic signal is composed of the sum of two sweep signals that have been repeated twice with total time of 4 s and sampling interval of 0.004 s. Two sweep signals have the bandwidth frequency of 10 to 40 Hz and 60 to 80 Hz respectively. The synthetic signal is shown in Figure 1. Signal 2 1.5 1 0.5 0 -0.5 -1 -1.5 -2
0
0.5
1
1.5
2 Time (s)
2.5
3
3.5
4
Figure 1. A synthetic signal comprised of the sum of two sweep signals that have been repeated twice.
The obtained time-frequency distributions using the aforementioned methods (WVD, PWVD, SPWVD and RID) are shown in Figure 2. As can be seen, WVD has serious cross–terms which are reduced in PWVD and SPWVD but the resolution is also reduced. The RID approach provides an advantage over the WVD, PWVD and SPWVD. The RID retains almost all of the desirable properties of the WVD, but with considerably reduced interference.
WVD
PWVD
0
18
20
0
18
20 16
40 60
14
80
16
40 60
14
80
100
12
120
100
12
120 0.5
1
1.5
2
2.5
3
3.5
4
0.5
10
1
1.5
SPWVD
8
0 20 Frequency (HZ)
2
2.5
3
3.5
4
10
(b)
(a)
6
40
RID
8
0 20
6
40
60
4
80
4
60 80
2
100 120 0.5
1
1.5
2 2.5 Time (s)
3
(c)
3.5
4
0
2 100 120 0.5
1
1.5
2
2.5
3
3.5
4
0
(d)
Figure 2. Time–Frequency distribution of synthetic signal (Figure 1) obtained by (a) WVD, (b) PWVD, (c) SPWVD and (d) RID.
2.2 The Nature of Seismic Reflection Data Reflection seismology is a method of exploration geophysics that uses the principles of seismology to estimate the properties of the Earth's subsurface from reflected seismic waves [43, 44]. The main objective of exploration seismology is identifying the subsurface structures associated with energy sources, especially hydrocarbons. To achieve this goal, mechanical waves are generated by seismic sources known as seismic source wavelet and travel in the Earth at a speed governed by the acoustic impedance of the Earth layers. The acoustic impedance of each layer is defined as the product of seismic wave velocity and the density of the Earth layer [45, 46]. Some of the emitted energy is reflected back when encountering an interface between two layers with different acoustic impedances, while some of the energy will be transmitted through the boundary. The reflected waves from different interfaces are recorded over a time series called the seismic trace or trace by land (geophone) or marine (hydrophone) receivers. The amplitude of the reflected wave is determined by multiplying the amplitude of the incident
wave by the seismic reflection coefficient, determined by the impedance contrast between the two materials. In fact, the Earth is a linear time-invariant system with an impulse response known as Earth reflection coefficient series. Therefore, a seismic trace is the convolution of seismic source wavelet with the Earth reflection coefficient series. At its most basic, the seismic reflection technique consists of generating seismic waves and measuring the time taken for the waves to travel from the source, reflect off an interface and be detected by an array of receivers at the surface [43]. Knowing the travel times from the source to various receivers, and the velocity of the seismic waves, a geophysicist then attempts to build up an image of the subsurface. Generally, the seismic data can be recorded in two and three dimensions. Due to various economic and scientific reasons, 3D recording of reflection seismic data is much better than 2D recording [47, 48]. A 3D recorded seismic data can be considered as pseudo 3D model of the Earth after special data processing of reflection seismology. Figure 3 shows conceptually what we are trying to achieve with such a recording.
Figure 3. Schematic showing a seismic data after processing relating to real-world geology (after [49]).
Seismic resolution is a key to extract stratigraphical details from seismic data. The resolution is defined as the ability to distinguish between two near interface (i.e. events) under the surface. The relationship between the layer thickness, the signal wavelength and the resolution of top and bottom horizon of the
ఒ
thin bed was first studied by Ricker [50]. At a layer thickness of less than (the tuning thickness), where ସ
ߣis the predominant wavelength of the seismic wavelet, which is calculated by the layer velocity, the reflective wavelets from top and bottom of the layer interfere constructively and produce a single event of high amplitude. It means two interfaces become indistinguishable in time. At spacings greater than that, two interfaces begin to be resolvable as two separate events, so the seismogram contains complete information for each interface. Therefore, seismic resolution is reduced by reduction in layer thickness and interpretation of these events has some problems. As mentioned before, we treat buried channels as thin layers. In this regard, a synthetic 2D seismic model of a thin layer with 30 Hz Ricker wavelet is designed. In Figure 4(a) the red line shows the top interface of the thin layer and the blue line shows the bottom interface; the thickness of the layer is between 0–135 m. The velocities and densities of different layers and the 2D seismic section of a thin bed are shown in Figure 4(a) and 4(b). To compare the results of RID to other time-frequency distributions and also to show the function of time-frequency distribution in seismology, we applied time-frequency distributions on every trace of the 2D synthetic seismic section respectively. A 2D time-frequency section is produced from each trace and a time-frequency cube is produced by placing the 2D time-frequency sections in succession. The different dimentions of timefrequency cube are time, frequency and trace number. Then, we extracted iso-frequency sections from the cube, which are 30 Hz in our case. By comparison with other sections in Figure 5, it is clear that RID can better separate top and bottom interfaces of a thin layer.
(a) Trace No. 0
20
40
60
80
100
120
140
160
180
200
Time Sample
20 40 60 80 100 120
(b) Figure 4. (a) Geological model of thin bed (b) synthetic seismic section.
-7
STFT 20
12
60 80
10
12
80
10
100 50
100
(a)
150
200
SPWVD
8
120
50
100
(b)
150
200
RID
6
20 Time Sample
40 60
100
8
6
20 4
40 60
4
40 60
2
80 100 120
x 10 14
20
40
120
-7
WVD
x 10 14
2
80 100
50
100 Trace No.
150
(c)
200
120
50
100
150
200
(d)
Figure 5. 30 Hz iso-frequency section of (a) STFT (b) WVD (c) SPWVD (d) RID.
2.3 Seismic Time-Frequency Attributes Seismic attribute is a quantity extracted or derived from seismic data that leads to a better geological or geophysical interpretation of the data [40]. The study of seismic attributes is an effective way to gather some important information about the geometry and the physical parameters of the subsurface such as
reservoir thickness, faults, porosity, permeability and hydrocarbon reservoir properties. There are a large number of seismic attributes, and each of them has a special application in practice [38]. Time-frequency attribute is an important category of seismic attributes, which is widely used in seismic data interpretation. In exploration seismology, a buried channel can be considered as a thin layer. Impulse response of a thin bed shown in Figure 6 can be defined as [51]: ܴሺݐሻ ൌ ܴଵ ݓሺݐሻ ܶଵ ܴଶ ݓሺ ݐ ߬ሻሺሻ where ܴሺݐሻ is the reflection from thin bed, ݓሺݐሻis the seismic source wavelet, ܴଵ and ܴଶ are the reflection coefficients of top and bottom boundary of the thin bed, respectively. ܶଵ is the transmission coefficient of top boundary of the thin bed and ߬ is the temporal thickness of the thin bed. The transfer function [52] can be computed by transforming the Eq. (7) into frequency domain using Fourier transform as: ܴሺ݆߱ሻ ൌ ൫ܴଵ ݁ ିఠ௧ ܶଵ ܴଶ ݁ ିఠሺ௧ାఛሻ ൯ܹሺ݆߱ሻ ൌ ܭሺ݆߱ሻܹሺ݆߱ሻሺͺሻ where ܹሺ݆߱ሻis the Fourier transform of߱ሺݐሻ and ܭሺ݆߱ሻ is the Fourier transform of the reflectivity series. It should be noted that the thin bed acts as a filter and changes the spectrum of incident waves. If ܴଵ ൌ െܴଶ ൌ ͳ and the transmission coefficient ܶଵ is close to 1, frequency characteristics of reflectivity series would be as Eq. (9) and its amplitude spectrum is defined as Eq. (10). ܭሺ݆߱ሻ ൌ ൫ͳ െ ݁ ିఠఛ ൯݁ ିఠ௧ ሺͻሻ ȁܭሺ݆߱ሻȁ ൌ ඥʹሾͳ െ ܿݏሺ߱߬ሻሿሺͳͲሻ This condition occurs when the velocity of the channel infills is larger than surrounding formation. In this case, destructive interference occurs at the edges of channels, so the minimum amplitude can be
expected. Due to reduction of the destructive interference while moving toward the center of the channel, the amplitude increases [51]. If the velocity of the channel infills is larger than the upper layer and smaller than the bottom layer, then ܴଵ ൌ ܴଶ ൌ ͳ and Eq. (9) and (10) can be rewritten as Eq. (11) and (12). ܭሺ݆߱ሻ ൌ ൫ͳ ݁ ିఠఛ ൯݁ ିఠ௧ ሺͳͳሻ ȁܭሺ݆߱ሻȁ ൌ ඥʹሾͳ ܿݏሺ߱߬ሻሿሺͳʹሻ In this case, unlike the previous case, constructive interference occurs at the edges of channels, so the maximum amplitude can be expected. Due to increasing of the destructive interference while moving toward the center of the channel, the amplitude decreases. In both cases, the reflected waves from the top and base of the channel interfer with each other in such a way that the dominant frequency will be maximum at the edge of channel. Also, dominant frequency decreases while moving toward the center of the channel.
Figure 6. Thin bed model.
Due to the high sensitivity of the seismic data frequency to channel thickness, we used seismic timefrequency attributes to reveal the buried channels which are one of the important oil traps. Amplitude of instantaneous centroid frequency, amplitude of instantaneous bandwidth frequency, amplitude of RMS
frequency, peak frequency and peak amplitude are effective time-frequency attributes to highlight stratigraphic features such as channels. These kinds of attributes reveal the variation of frequency content of seismic data. Therefor, the thickness of different parts of a channel can be simply detected. The instantaneous center frequency, instantaneous bandwidth frequency and RMS frequency are defined as Eq. (13) to Eq. (15), respectively [53]. ஶ
݂ ሺݐሻ ൌ
݂ ൈ ܴܶܨሺݐǡ ݂ሻ݂݀ ஶ
ܴܶܨሺݐǡ ݂ሻ݂݀
ሺͳ͵ሻ
ஶ
݂ଶ ሺݐሻ
ൌ
ሺ݂ െ ݂ ሺݐሻሻଶ ൈ ȁܴܶܨሺݐǡ ݂ሻȁଶ ݂݀ ஶ
ȁܴܶܨሺݐǡ ݂ሻȁଶ ݂݀
ሺͳͶሻ
ஶ
݂ଶ ሺݐሻ
ൌ
݂ ଶ ൈ ܴܶܨሺݐǡ ݂ሻ݂݀ ஶ
ܴܶܨሺݐǡ ݂ሻ݂݀
ሺͳͷሻ
where ܴܶܨሺݐǡ ݂ሻ is the time-frequency distribution. In this paper, the pointed out seismic time-frequency attributes are extracted from RID and SPWVD to identify the buried channels in synthetic and real seismic data. To identify the buried channel in seismic data, we must apply time-frequency decomposition on seismic data and then extract the aforementioned attributes.
3. Results To show the location of buried channel, we must apply time-frequency decomposition on seismic data and extract the aforementioned attributes. First, we applied the RID and SPWVD methods to each trace of 3D seismic data and computed considered attributes. Secondly, we located the resulting attributes in the same place of the considered trace; we repeated this action for the rest of traces. Therefore, we obtained 3D data of seismic attributes of the same size as the original 3D seismic data. Thirdly, we extracted time slices of a target horizon of 3D seismic attributes, which show the spread of channel. The related flow chart is illustrated in Figure 7.
Start
Take the seismic trace
Compute time-frequency decomposition
Compute seismic attribute
Place seismic attribute in same location as seismic trace
End
Figure 7. Flow chart of computing seismic attributes.
To show the efficiency of this method, we produced a synthetic sand-filled channel which is embedded in shale. The ratio of the reflection coefficients of top and bottom of the channel is considered as
ோభ ோమ
ൌ
െͳ. The 30 Hz Ricker wavelet is chosen for forward modeling as seismic source wavelet and the time sampling interval is 0.004 s. 3D seismic data of the channel and the time slice associated with the channel at time sample 41 are shown in Figure 8. The red oval shows the channel position on 3D seismic data. The channel spreads straightly through the Earth, which is clear in the displayed attributes. Peak amplitude and peak frequency attributes extracted from RID and SPWVD are shown in Figure 9. As can be seen in Figure 10 the amplitude of instantaneous center frequency, instantaneous bandwidth frequency and RMS frequency of RID and SPWVD clearly show the expansion of the channel, but there is an absolute difference between the resolution of RID and SPWVD.
X lin e No .
20 40 60 80 100 120 20
40
60 Inline No.
(a)
80
100
(b) Figure 8. (a) 3D synthetic seismic cube (b) Horizontal slice at time sample 41.
Magnitude
(b)
(c)
(d)
Frequency
(a)
Figure 9. Horizontal slice at the target zone (time sample 41) from (a) peak amplitude cube of RID and (b) peak amplitude cube of SPWVD and (c) peak frequency cube of RID and (d) peak frequency cube of SPWVD.
120
Amplitude of instantaneous center frequency of SPWVD 0.2 20 40 60 0.1 80 100 120 0 20 40 60 80 100 120
(a)
(b)
Amplitude of instantaneous bandwidth frequency of RID 20 40 60 80 100 120
0.2 0.1
20
40
60
(c)
80
100
120
0
Amplitude of instantaneous bandwidth frequency of SPWVD 20 40 60 80 100 120
Xline No.
0.1 0
(e)
40
60
80
100
120
0
Amplitude of RMS frequency of SPWVD
0.2
60 80 Inline No.
0.1
(d)
20 40 60 80 100 120 40
0.2
20
Amplitude of RMS frequency of RID
20
Magnitude
Amplitude of instantaneous center frequency of RID 0.2 20 40 60 0.1 80 100 120 0 20 40 60 80 100 120
100
120
20 40 60 80 100 120
0.2 0.1
20
40
60
80
100
120
0
(f)
Figure 10. Horizontal slice at the target zone (time sample 41) from instantaneous center frequency of (a) RID (b) SPWVD, Horizontal slice at the target zone (time sample 41) from instantaneous bandwidth frequency of (c) RID (d) SPWVD and Horizontal slice at the target zone (time sample 41) from RMS frequency of (e) RID (f) SPWVD.
Real seismic data is used to recognize the efficiency of RID to identify hidden channels by the above mentioned seismic time-frequency attributes. The real 3D seismic data is taken from one of the oil fields in the south-west of Iran with 400 inline and 600 crossline traces, with time sampling interval 0.004 s, which contains a channel at the 1.8 sec (time sample 10) time slice based on previous study [54]. In Figure 11, we show the 3D seismic cube and an extracted time slice at the 1.8 sec. As can be seen, it is difficult to recognize the buried channel in seismic time slice. The three seismic time-frequency attributes based on RID, the amplitude of instantaneous center frequency, amplitude of instantaneous bandwidth frequency and amplitude of RMS frequency attributes are shown in Figure 12. We extracted the time slices at 1.8 s from each of the three attributes and compared the result with that of SPWVD in Figures
13-15. Due to high resolution of RID, the different branches of buried channel (even thin branches) are more visible in the results of RID than SPWVD.
(a)
(b)
Figure 11. (a) 3D real seismic cube and (b) time slice at the 1.8 sec.
(b)
(a)
(c)
Figure 12. 3D cube of seismic time-frequency attributes extracted from RID (a) amplitude of instantaneous center frequency (b) amplitude of instantaneous bandwidth frequency (c) amplitude of RMS frequency.
Amplitude of instantaneous center frequency of RID 600
600
0.3 0.25
0.25
500
400
0.2
400
0.2
300
0.15
300
0.15
200
0.1
200
0.1
100
0.05
100
0.05
50
100
150
200 250 Xline No.
(a)
300
350
400
0
50
100
150
200
250
300
350
400
Magnitude
Inline No.
500
Amplitude of instantaneous center frequency of SPWVD 0.3
0
(b)
Figure 13. Time slice at the target zone (time sample 10) from amplitude of instantaneous center frequency extracted from (a) RID (b) SPWVD.
Amplitude of instantaneous bandwidth frequency of SPWVD
Amplitude of instantaneous bandwidth frequency of RID 600
0.4
600
500
0.35
500
0.4 0.35 0.3
0.3 Inline No.
400
0.25
0.25 300
0.2
300
0.2
200
0.15
200
0.15 0.1
0.1 100
Magnitude
400
100 0.05
0.05 50
100
150
200 250 Xline No.
300
350
400 0
50
100
150
200
250
300
350
0
400
(b)
(a)
Figure 14. Time slice at the target zone (time sample 10) from amplitude of instantaneous bandwidth frequency extracted from (a) RID (b) SPWVD.
Amplitude of RMS frequency of RID
Amplitude of RMS frequency of SPWVD
600
0.5
500
0.4
Inline No.
0.5
500
0.4
Magnitude
400
600
400
0.3
0.3 300
300
0.2
0.2 200
200 0.1
100
50
100
150
200 250 Xline No.
300
350
400
0
0.1
100
50
100
150
200
250
300
350
400
0
(b)
(a)
Figure 15. Time slice at the target zone (time sample 10) from amplitude of RMS frequency extracted form (a) RID (b) SPWVD.
We calculated the peak amplitude and peak frequency 3D seismic attributes using RID which are shown in Figure 16. Time slices at 1.8 s from two seismic attributes computed based on two time-frequency distribution (RID and SPWVD) are shown in Figure 17. As can be seen, low peak frequency corresponds to thicker part of the channel and vice versa. Conversely, high peak amplitude corresponds to thicker part of the channel. Therefore, we can identify the edge of channel and variation of its thickness.
(b)
(a)
Figure 16. (a) Peak amplitude (b) peak frequency attributes extracted from RID.
Peak Amplitude of RID 600 500
Peak Amplitude of SPWVD
1
600
0.8
500
0.8
400 0.6
300
0.6 300
0.4
200 100
0.2 50
100
150
200
250
300
350
200
0.4
100
0.2
400
50
100
150
(a)
250
300
350
400
(b) Peak Frequency of SPWVD
Peak Frequency of RID 600
60
600
60 50
500
50
500
400
40
400
40 30
30
300
200
20
200
20
100
10
100
10
300
50
100
150
200
250
Xline No.
(c)
300
350
400
0
50
100
150
200
250
300
350
Frequency
Inline No.
200
Magnitude
400
1
400
(d)
Fig. 17. Time slice at the target zone (time sample 10) from (a) peak amplitude cube extracted from and (b) peak amplitude cube extracted from SPWVD. (c) peak frequency cube extracted from and (d) peak frequency cube extracted from SPWVD.
4. Conclusions Frequency content of seismic signals carries valuable information of subtle changes, so analysis of this information is an effective way to survey subsurface structures. A spectral decomposition is a powerful tool in this field. Seismic signals are non-stationary multi-component signals and by employing conventional methods such as WVD, PWVD and SPWVD, some problems such as Heisenberg
uncertainty principle and cross–terms emerged. In this paper, we introduced the Reduced Interference Distribution to overcome the restrictions of the aforementioned methods. The RID keeps almost all the desirable properties of the WVD, but reduces considerably the interference cross–terms. We used this method of spectral decomposition to calculate the instantaneous central frequency, instantaneous bandwidth frequency, RMS frequency, peak frequency and peak amplitude attributes, which are used for imaging of the buried channel which is one of the oil traps. Actually, we wanted to improve the results of these attributes by RID. We tested this method of buried channel imaging on a synthetic and a real seismic data. We saw a reasonably good correlation between the channel's thickness in time and the mentioned seismic attributes. Therefore, these attributes can clearly delineate spreading channels vertically and horizontally, and correlate well with channel thickness.
References [1] [2] [3]
[4] [5]
[6] [7] [8] [9]
[10]
J. Pujol, Elastic Wave Propagation and Generation in Seismology, Cambridge University Press, 2003. Ö. Yilmaz, S.M. Doherty, Seismic Data Analysis: processing, inversion, and interpretation of seismic data, Society of Exploration Geophysicists, 2001. M.J. Duchesne, E.J. Halliday, J.V. Barrie, Analyzing seismic imagery in the time–amplitude and time–frequency domains to determine fluid nature and migration pathways: A case study from the Queen Charlotte Basin, offshore British Columbia, Journal of Applied Geophysics 73 (2011) 111-120. G. Partyka, J. Gridley, J. Lopez, Interpretational applications of spectral decomposition in reservoir characterization, The Leading Edge 18 (1999) 353-360. M.C. de Matos, Characterization of thin beds through joint time-frequency analysis applied to a turbidite reservoir in Campos Basin, Brazil, SEG Technical Program Expanded Abstracts 2005, pp. 1429-1432. J. Liu, K.J. Marfurt, Instantaneous spectral attributes to detect channels, Geophysics 72 (2007) P23-P31. R.H. Herrera, J. Han, M.v.d. Baan, Applications of the synchrosqueezing transform in seismic time-frequency analysis, GEOPHYSICS 79 (2014) V55-V64. J.P. Castagna, S. Sun, R.W. Siegfried, Instantaneous spectral analysis: Detection of low-frequency shadows associated with hydrocarbons, The Leading Edge 22 (2003) 120-127. M. Zarei, A. Roshandel Kahoo, H.R. Siahkoohi, Gas detection using deconvolutive short time Fourier transform, International Geophysical Conference and Oil & Gas Exhibition, Istanbul, Turkey, 2012, pp. 1-4. H. Sattari, A. Gholami, H.R. Siahkoohi, Seismic data analysis by adaptive sparse time-frequency decomposition, Geophysics 78 (2013) V207-V217.
[11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23]
[24]
[25]
[26] [27] [28] [29]
[30] [31]
[32] [33]
S. Sinha, P.S. Routh, P.D. Anno, J.P. Castagna, Spectral decomposition of seismic data with continuous-wavelet transform, Geophysics 70 (2005) P19-P25. R. Askari, H.R. Siahkoohi, Ground roll attenuation using the S and x-f-k transforms, Geophysical Prospecting 56 (2008) 105-114. A. Roshandel Kahoo, H.R. Siahkoohi, Gas detection from AVO analysis in time-frequency domain, 71st EAGE Conference & Exhibition, 2009. A. Roshandel Kahoo, H.R. Siahkoohi, Random noise suppression from seismic data using timefrequency peak filtering, 71st EAGE Conference & Exhibition, 2009. G.F. Margrave, M.P. Lamoureux, D.C. Henley, Gabor deconvolution: Estimating reflectivity by nonstationary deconvolution of seismic data, Geophysics 76 (2011) W15-W30. X.-K. Sun, Z. Sam, H.-W. Xie, Nonstationary sparsity-constrained seismic deconvolution, Applied Geophysics 11 (2014) 459-467. D. Gabor, Theory of communication. Part 1: The analysis of information, Electrical Engineers Part III: Radio and Communication Engineering, Journal of the Institution of 93 (1946) 429-441. A. Gholami, Sparse Time-Frequency Decomposition and Some Applications, Geoscience and Remote Sensing, IEEE Transactions on 51 (2013) 3598-3604. E. Wigner, On the quantum correction for thermodynamic equilibrium, Physical Review 40 (1932) 749. B. Boashash, Time Frequency Analysis, Elsevier Science, 2003. S. Mallat, A Wavelet Tour of Signal Processing: The Sparse Way, Elsevier Science, 2008. R.G. Stockwell, L. Mansinha, R. Lowe, Localization of the complex spectrum: the S transform, Signal Processing, IEEE Transactions on 44 (1996) 998-1001. M. Radad, A. Gholami, H.R. Siahkoohi, S-transform with maximum energy concentration: Application to non-stationary seismic deconvolution, Journal of Applied Geophysics 118 (2015) 155-166. F. Auger, P. Flandrin, Improving the readability of time-frequency and time-scale representations by the reassignment method, Signal Processing, IEEE Transactions on 43 (1995) 1068-1089. H.I. Choi, W.J. Williams, Improved time-frequency representation of multicomponent signals using exponential kernels, Acoustics, Speech and Signal Processing, IEEE Transactions on 37 (1989) 862-871. J. Jechang, W.J. Williams, Kernel design for reduced interference distributions, Signal Processing, IEEE Transactions on 40 (1992) 402-412. L. Cohen, Time-frequency Analysis, Prentice Hall PTR, 1995. M. Rauch-Davies, M. Ralston, Spectral decomposition–transform methods and fluid and reservoir prediction case study, 67th EAGE Conference & Exhibition, 2005. M. Ralston, M. Rauch-Davies, K. Li-Chun, X. Hui-Ping, Y. Di-Sheng, General method to reduce crossterm interference in the Wigner-Ville decomposition, SEG Technical Program Expanded Abstracts 2007, pp. 870-874. Y. Li, X. Zheng, Wigner-Ville distribution and its application in seismic attenuation estimation, Applied Geophysics 4 (2007) 245-254. Y. Hu, K.D.K. Luk, W.W. Lu, A. Holmes, J.C.Y. Leong, Comparison of time-frequency distribution techniques for analysis of spinal somatosensory evoked potential, Medical and Biological Engineering and Computing 39 (2001) 375-380. W.J. Williams, Reduced interference distributions: biological applications and interpretations, Proceedings of the IEEE 84 (1996) 1264-1280. M.B. Widess, How thin is a thin bed?, Geophysics 38 (1973) 1176-1180.
[34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48]
[49] [50]
[51] [52] [53] [54]
J.D. Robertson, H.H. Nogami, Complex seismic trace analysis of thin beds, Geophysics 49 (1984) 344-352. H. Zeng, A. John, G. Katherine, How thin is a thin bed? An alternative perspective, SEG Technical Program Expanded Abstracts, 2008, pp. 834-838. H. Zeng, Thin-bed prediction by geomorphology-constrained waveform analysis, SEG Technical Program Expanded Abstracts 2014, pp. 1395-1399. G. Li, M.D. Sacchi, S. Zhu, Y. Wang, Thin bed identification using zero-crossing-time stratal slices, SEG Technical Program Expanded Abstracts 2015, pp. 1678-1682. S. Chopra, K.J. Marfurt, Seismic Attributes for Prospect Identification and Reservoir Characterization, Society of Exploration Geophysicists, 2007. X. Wu, T. Liu, Spectral decomposition of seismic data with reassigned smoothed pseudo Wigner– Ville distribution, Journal of Applied Geophysics 68 (2009) 386-393. M.T. Taner, F. Koehler, R. Sheriff, Complex seismic trace analysis, Geophysics 44 (1979) 10411063. B. Boashash, Time-Frequency Signal Analysis and Processing: A Comprehensive Reference, Elsevier Science, 2015. H. Hassanpour, Time-frequency based detection of newborn EEG seizure, Queensland University of Technology, 2004. R.E. Sheriff, L.P. Geldart, Exploration Seismology, Cambridge University Press, 1995. M.R. Gadallah, R. Fisher, R.L. Fisher, Exploration Geophysics, Springer, 2008. W.M. Telford, L.P. Geldart, R.E. Sheriff, Applied Geophysics, Cambridge University Press, 1990. D.A. Herron, R.B. Latimer, First Steps in Seismic Interpretation, Society of Exploration Geophysicists, 2011. A. Cordsen, M. Galbraith, J. Peirce, B.A. Hardage, Planning Land 3-D Seismic Surveys, Society of Exploration Geophysicists, 2000. A.R. Brown, Interpretation of three-dimensional seismic data, 7 ed., Published jointly by American Association of Petroleum Geologists and the Society of Exploration Geophysicists, 2011. D.L. Baars, W.L. Watney, D.W. Steeples, E.A. Brostuen, Petroleum: a primer for Kansas, Educational Series, 1993. N. Ricker, Wavelet contraction, wavelet expansion, and the control of seismic resolution, Geophysics 18 (1953) 769–792. R.D. Han, Z.H. Wan, M.S. Chen, H.Y. Zhang, Application of time-frequency attributes based on generalized s-transform for thin bed indication, 73rd EAGE Conference & Exhibition, 2011. J.G. Proakis, D.G. Manolakis, Digital Signal Processing, Principles, Algorithms, and Applications, 4 ed., Pearson Prentice Hall, 2007. A.E. Barnes, Instantaneous spectral bandwidth and dominant frequency with applications to seismic reflection data, Geophysics 58 (1993) 419-428. R. Mohebian, M. Yari, M.A. Riahi, R. Ghanati, Channel detection using instantaneous spectral attributes in one of the SW Iran oil fields, Bollettino di Geofisica Teorica ed Applicata 54 (2013) 271-282.
Azita Nikoo received the B.S. degree in physics from Razi University of Kermanshah, Kermanshah, Iran, in 2009, the M.S. degree in Geophysics from Shahrood University of Technology, Shahrood, Iran, in 2013. She is the project implementers in Gity Kavosh Zaniar Consultants Company, Iran. Her research interests include application of time-frequency transform in exploration seismology, seismic attributes, and AVO. Hamed Saadatnia received the B.S. degree in Geology from University of Tehran, Tehran, Iran, in 2001 and the M.S. degrees in exploration seismology from the University of Tehran, Tehran, Iran, in 2004. Since 2007, he has been with the National Iranian Oil Company (NIOC), Tehran, Iran, where he is currently a geophysicist in seismic interpretation department. His research interests include reservoir characterization and seismic data interpretation.
Hamid Hassanpour received the B.S. degree in computer engineering from Iran University of Science and Technology, Tehran, Iran, in 1993, the M.S. degree in computer engineering from Amirkabir University of Technology, Tehran, Iran, in 1996, and the Ph.D. from the Queensland University of Technology, Brisbane, Australia, in 2004. He has a professor position in faculty of Computer Engineering & IT at the Shahrood University of Technology, Iran. His research interests include Image Processing, Signal Processing, time-frequency signal processing and analysis. Amin Roshandel Kahoo received the B.S. degree in mining engineering from University of Birjand, Birjand, Iran, and the M.S. and Ph.D. degrees in exploration seismology from the University of Tehran, Tehran, Iran, in 2005 and 2010 respectively. Since September 2010, he has been with the School of Mining, Petroleum & Geophysics Engineering, Shahrood University of Technology, Shahrood, Iran, where he is currently an Assistant Professor. His research interests include application of time-frequency transform in exploration seismology, seismic attributes, and seismic data processing and interpretation.
Azita Nikoo, Amin Roshandel Kahoo, Hamid Hassanpour, Hamed Saadatnia
PII: DOI: Reference:
S1051-2004(16)30002-1 http://dx.doi.org/10.1016/j.dsp.2016.03.008 YDSPR 1928
To appear in:
Digital Signal Processing
Please cite this article in press as: A. Nikoo et al., Using a time-frequency distribution to identify buried channels in reflection seismic data, Digit. Signal Process. (2016), http://dx.doi.org/10.1016/j.dsp.2016.03.008
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Using a Time-Frequency Distribution to Identify Buried Channels in Reflection Seismic Data 1 Azita Nikoo a, േ , Amin Roshandel Kahoo b, Hamid Hassanpour c, Hamed Saadatnia d a
Department of Mining, Petroleum and Geophysics Engineering, University of Shahrood, Iran
b
Assistant Professor, Department of Mining, Petroleum and Geophysics Engineering, University of
Shahrood, Iran c
Professor, Department of Computer Engineering & IT, University of Shahrood, Iran
േ
d
Exploration Directorate, NIOC, Iran
Corresponding author. E-mail addresses: [email protected] (A. Nikoo), [email protected] (A. Roshandel Kahoo), [email protected] (H. Hassanpour), [email protected] (H. Saadatnia).
Abstract Geological events such as thin beds and channels cannot be easily revealed on seismic sections due to the interference of reflections from the top and bottom of the layer. Buried channel is one of the hydrocarbon traps, which is important in oil and gas exploration. Spectral decomposition can be used to indicate subtle changes in channel thickness. In using the Fourier transform only the frequency content of the changes is displayed; hence, the exact onset of the changes is missed. Time-frequency distributions are suitable approaches to display and interpret information embedded in non-stationary signals such as seismic signals. Spectral decomposition can be used for seismic attributes calculation, which are used for imaging of thin beds. The conventional spectral decompositions such as Short Time Fourier Transform (STFT) and Wigner-Ville Distribution (WVD) have limitations in terms of Heisenberg uncertainty principle and cross-terms artefacts, respectively. In this paper, we used the Reduced Interference Distribution (RID) for buried channel identification to overcome the mentioned limitations. We compared the obtained results of RID with that of Smoothed Pseudo WVD. The real seismic data is selected from one of the western oil fields of Iran. Our results indicate that a better resolution is achieved by RID in both vertical and lateral stratigraphic features. Keywords: Spectral decomposition, Seismic attributes, Wigner–Ville distribution, Reduced interference distribution.
1. Introduction Seismic waves are mechanical waves that travel through the Earth’s layers or along its surface. Seismic signals propagate into the Earth and reflect from boundary separating rocks with different properties [1]. Indeed, seismic reflection records carry important information, such as amplitude and frequency, associated with the Earth and its subsurface structure. Seismic waves are non-stationary because the Earth acts as a low-pass filter and changes the frequency content of seismic signals during the propagation [2]. Ordinary 3-D seismic data contains energy reflected from the subsurface at a wide range of frequencies, all of which are compounded in a typical seismic volume. Different frequencies of seismic reflection data can image geological features with frequency-dependent nature. Geological interpretations from seismic reflection data are mostly drawn from the time-amplitude domain. It is used to improve and increase the application of seismic data sets for geological interpretations [3]. Spectral decomposition can be helpful to see the amplitude or phase of reflections at particular frequencies. Besides, the spectral decomposition analysis allows the geophysicist to quantify amplitude variation with frequency and extract geological information from seismic data associated with frequency that normally cannot be viewed in conventional seismic data [4]. Fourier transform is not a suitable tool to analyze non-stationary signals. It cannot show the frequency variations of seismic signals over time. Time-frequency analysis is more appropriate for seismic data due to their non-stationary nature. Spectral decomposition has been widely used in seismic data processing and interpretation. We can explore stratigraphic settings, thin beds and buried channels using spectral decomposition of seismic sections [4-7]. The spectral decomposition can be employed to detect lowfrequency shadows associated with hydrocarbons [8-10]. Sinha et al. in [11] employed the continuous wavelet transform to computing a time-frequency map for seismic signal. Time-frequency map is used for seismic attributes computation. Askari and Siahkoohi in [12] attenuated the ground roll noise (a common type of coherent noise in seismic data) from seismic data using the S and x-f-k transforms. Roshandel Kahoo and Siahkoohi in [13] performed amplitude versus offset (AVO) analysis in the time-frequency
domain for detection of the gas reservoir. They used modified high resolution time–frequency transforms to suppress random noise from seismic data using peak filtering method [14]. The spectral decomposition can be used for time varying deconvolution [15, 16]. There are several methods in literature for spectral decomposition, such as short time Fourier transform (STFT) [17, 18], Wigner–Ville distribution (WVD) [19, 20], wavelet transform [21] and S–transform [22, 23]. Conventional spectral decomposition methods have some limitations, such as the Heisenberg uncertainty principle and cross–terms existence, which limit their applications in analyzing non-stationary signals. The conventional method of making a time-frequency map using the STFT limits time-frequency resolution by a predefined window length [11]. Therefore, the valuable information is lost or blurred in the time–frequency representations due to a trade-off between the time and frequency resolution. No window function is used to calculate the WVD of a signal. Therefore, WVD has high resolution both in time and frequency domain. However, the applications of this distribution are limited by cross–terms existence. Several methods such as Pseudo WVD (PWVD) and Smoothed Pseudo WVD (SPWVD) have been introduced to eliminate the cross–terms [24-27]. Smoothed versions of WVD use a time and frequency dependent function as a kernel to reduce cross terms. However, these methods reduce the energy distribution concentration in time-frequency plane [14, 28-30]. Cone shaped distribution (CSD) and adaptive spectrogram (ADS) are used for reducing cross-terms in biological signal [31]. In this paper, we used Reduced Interference Distribution (RID) [32] to calculate the time–frequency distribution which is found to be more appropriate for this application. Seismic signals are broadband (8-80 Hz) FM signals and RID is a suitable choice here. This method of distribution with slight impact on the time–frequency resolution increasingly reduces the cross terms. Analysis of thin–bed tuning on seismic reflectivity has been studied extensively [33-37]. During the last decade, spectral decomposition technique has proven to be an excellent tool to describe thin beds associated with channel sands, alluvial fans, and other similar structures [38]. In this paper, we focused on extracting seismic attributes based on the RID such as instantaneous center frequency, instantaneous
bandwidth frequency, root mean square (RMS) frequency, peak frequency and peak amplitude. These attributes are directly related to the thickness of the thin-bed. Since a buried channel is like a thin-bed, we used the extracted attributes for detection of buried channels.
2. Methodology In this section first we explained one of the effective time-frequency distributions and also a way in which to improve them. We then illustrated the nature of seismic data and the characteristics of thin layer on seismic sections, and proposed some useful seismic attributes to reveal buried channels.
2.1 Time-Frequency Distribution Resolution of time and frequency in STFT is limited by the Heisenberg uncertainty principle. One of the advantages of WVD over the STFT is to distribute energy density in time–frequency plane without losing time or frequency resolution. The WVD of signalݔሺݐሻis defined as [39]: ାஶ
ఛ
ఛ
ଶ
ଶ
ܹܸܦሺݐǡ ݂ሻ ൌ ିஶ ܺ ቀ ݐ ቁ ܺ כቀ ݐെ ቁ ݁ ିଶగఛ ݀߬ሺͳሻ where ܺሺݐሻ is the analytic signal of signal ݔሺݐሻ, ܺ כሺݐሻ denotes the complex conjugate of ܺሺݐሻ, ݂ is frequency, ݐis time and ߬ is the time delay. The analytic signal is a complex signal [40] and can be computed as: ܺሺݐሻ ൌ ݔሺݐሻ ݆ ܪሾݔሺݐሻሿሺʹሻ where ܪሾݔሺݐሻሿ is the Hilbert transform of ݔሺݐሻ and ݆ ൌ ξെͳ. ܭ ሺݐǡ ߬ሻ is defined as the instantaneous autocorrelation function of ܺሺݐሻ: ߬ ߬ ܭ ሺݐǡ ߬ሻ ൌ ܺ ቀ ݐ ቁ ܺ כቀ ݐെ ቁ ሺ͵ሻ ʹ ʹ By substitution of Eq. (3) into Eq. (1), we can rewrite Eq. (1) as:
ାஶ
ܹܸܦሺݐǡ ݂ሻ ൌ ିஶ ܭ ሺݐǡ ߬ሻ݁ ିଶగఛ ݀߬ ሺͶሻ Eq. (4) shows that WVD is equal to the Fourier transform of the instantaneous autocorrelation function of ܺሺݐሻ with respect to ߬. Because of bilinear characteristic, cross-terms interference is emerged for multi-component signals. For a signal composed of two componentsݔሺݐሻ ൌ ݔଵ ሺݐሻ ݔଶ ሺݐሻ, based on the quadratic superposition principle [20, 41], the WVD is defined as: ܹܸܦ௫ ሺݐǡ ݂ሻ ൌ ܹܸܦ ௫భ ሺݐǡ ݂ሻ ܹܸܦ௫మ ሺݐǡ ݂ሻ ܹܸܦ ௫భ ǡ௫మ ሺݐǡ ݂ሻ ܹܸܦ௫మ ǡ௫భ ሺݐǡ ݂ሻ ሺͷሻ ᇣᇧᇧᇧᇧᇧᇧᇧᇧᇤᇧᇧᇧᇧᇧᇧᇧᇧᇥ ᇣᇧᇧᇧᇧᇧᇧᇧᇧᇧᇧᇤᇧᇧᇧᇧᇧᇧᇧᇧᇧᇧᇥ Ǧ
Ǧ
The two last terms in right hand of Eq. (5) are cross–terms with the following characteristics [20, 41]: 1- They are found midway between the interacting components; 2- They oscillate proportionally to the inter-components' time-frequency distance; 3- They have a direction of oscillation orthogonal to the straight line connecting the components. PWVD and SPWVD are two well-known distributions with less cross-term compared with WVD. In these distributions, the cross-terms are reduced by smoothing the time-frequency representation due to their oscillation property. The smoothing reduces the cross-terms, but with the cost of reduction of resolution in time–frequency representation. Here we exploit the RID method for reducing cross–terms which has the least impact on the resolution of the time–frequency representation. This distribution belongs to the general Cohen class of distributions. Cohen’s distribution of a signal ݔሺݐሻ is defined as [42]: ାஶ ାஶ ାஶ
߬ ߬ ܴܶܨ௫ ሺݐǡ ݂ሻ ൌ න න න ݁ ଶగ௩ሺ௧ି௨ሻ ݃ሺݒǡ ߬ሻܺ ቀ ݐ ቁ ܺ כቀ ݐെ ቁ ݁ ିଶగఛ ݀߬݀ݑ݀ݒሺሻ ʹ ʹ ିஶ ିஶ ିஶ
where ݃ሺݒǡ ߬ሻis a 2D kernel and specifies the distribution. The kernel of RID emphasizes the auto-terms and deemphasizes the cross-terms. The RID keeps almost all the desirable properties of the WVD, but
with considerably reduced interference. In this paper, we used Binomial kernel, which is a particularly attractive discrete form which can be computed very efficiently [32]. To study the efficiency of the WVD, PWVD, SPWVD and RID, we tested them on a synthetic signal. The synthetic signal is composed of the sum of two sweep signals that have been repeated twice with total time of 4 s and sampling interval of 0.004 s. Two sweep signals have the bandwidth frequency of 10 to 40 Hz and 60 to 80 Hz respectively. The synthetic signal is shown in Figure 1. Signal 2 1.5 1 0.5 0 -0.5 -1 -1.5 -2
0
0.5
1
1.5
2 Time (s)
2.5
3
3.5
4
Figure 1. A synthetic signal comprised of the sum of two sweep signals that have been repeated twice.
The obtained time-frequency distributions using the aforementioned methods (WVD, PWVD, SPWVD and RID) are shown in Figure 2. As can be seen, WVD has serious cross–terms which are reduced in PWVD and SPWVD but the resolution is also reduced. The RID approach provides an advantage over the WVD, PWVD and SPWVD. The RID retains almost all of the desirable properties of the WVD, but with considerably reduced interference.
WVD
PWVD
0
18
20
0
18
20 16
40 60
14
80
16
40 60
14
80
100
12
120
100
12
120 0.5
1
1.5
2
2.5
3
3.5
4
0.5
10
1
1.5
SPWVD
8
0 20 Frequency (HZ)
2
2.5
3
3.5
4
10
(b)
(a)
6
40
RID
8
0 20
6
40
60
4
80
4
60 80
2
100 120 0.5
1
1.5
2 2.5 Time (s)
3
(c)
3.5
4
0
2 100 120 0.5
1
1.5
2
2.5
3
3.5
4
0
(d)
Figure 2. Time–Frequency distribution of synthetic signal (Figure 1) obtained by (a) WVD, (b) PWVD, (c) SPWVD and (d) RID.
2.2 The Nature of Seismic Reflection Data Reflection seismology is a method of exploration geophysics that uses the principles of seismology to estimate the properties of the Earth's subsurface from reflected seismic waves [43, 44]. The main objective of exploration seismology is identifying the subsurface structures associated with energy sources, especially hydrocarbons. To achieve this goal, mechanical waves are generated by seismic sources known as seismic source wavelet and travel in the Earth at a speed governed by the acoustic impedance of the Earth layers. The acoustic impedance of each layer is defined as the product of seismic wave velocity and the density of the Earth layer [45, 46]. Some of the emitted energy is reflected back when encountering an interface between two layers with different acoustic impedances, while some of the energy will be transmitted through the boundary. The reflected waves from different interfaces are recorded over a time series called the seismic trace or trace by land (geophone) or marine (hydrophone) receivers. The amplitude of the reflected wave is determined by multiplying the amplitude of the incident
wave by the seismic reflection coefficient, determined by the impedance contrast between the two materials. In fact, the Earth is a linear time-invariant system with an impulse response known as Earth reflection coefficient series. Therefore, a seismic trace is the convolution of seismic source wavelet with the Earth reflection coefficient series. At its most basic, the seismic reflection technique consists of generating seismic waves and measuring the time taken for the waves to travel from the source, reflect off an interface and be detected by an array of receivers at the surface [43]. Knowing the travel times from the source to various receivers, and the velocity of the seismic waves, a geophysicist then attempts to build up an image of the subsurface. Generally, the seismic data can be recorded in two and three dimensions. Due to various economic and scientific reasons, 3D recording of reflection seismic data is much better than 2D recording [47, 48]. A 3D recorded seismic data can be considered as pseudo 3D model of the Earth after special data processing of reflection seismology. Figure 3 shows conceptually what we are trying to achieve with such a recording.
Figure 3. Schematic showing a seismic data after processing relating to real-world geology (after [49]).
Seismic resolution is a key to extract stratigraphical details from seismic data. The resolution is defined as the ability to distinguish between two near interface (i.e. events) under the surface. The relationship between the layer thickness, the signal wavelength and the resolution of top and bottom horizon of the
ఒ
thin bed was first studied by Ricker [50]. At a layer thickness of less than (the tuning thickness), where ସ
ߣis the predominant wavelength of the seismic wavelet, which is calculated by the layer velocity, the reflective wavelets from top and bottom of the layer interfere constructively and produce a single event of high amplitude. It means two interfaces become indistinguishable in time. At spacings greater than that, two interfaces begin to be resolvable as two separate events, so the seismogram contains complete information for each interface. Therefore, seismic resolution is reduced by reduction in layer thickness and interpretation of these events has some problems. As mentioned before, we treat buried channels as thin layers. In this regard, a synthetic 2D seismic model of a thin layer with 30 Hz Ricker wavelet is designed. In Figure 4(a) the red line shows the top interface of the thin layer and the blue line shows the bottom interface; the thickness of the layer is between 0–135 m. The velocities and densities of different layers and the 2D seismic section of a thin bed are shown in Figure 4(a) and 4(b). To compare the results of RID to other time-frequency distributions and also to show the function of time-frequency distribution in seismology, we applied time-frequency distributions on every trace of the 2D synthetic seismic section respectively. A 2D time-frequency section is produced from each trace and a time-frequency cube is produced by placing the 2D time-frequency sections in succession. The different dimentions of timefrequency cube are time, frequency and trace number. Then, we extracted iso-frequency sections from the cube, which are 30 Hz in our case. By comparison with other sections in Figure 5, it is clear that RID can better separate top and bottom interfaces of a thin layer.
(a) Trace No. 0
20
40
60
80
100
120
140
160
180
200
Time Sample
20 40 60 80 100 120
(b) Figure 4. (a) Geological model of thin bed (b) synthetic seismic section.
-7
STFT 20
12
60 80
10
12
80
10
100 50
100
(a)
150
200
SPWVD
8
120
50
100
(b)
150
200
RID
6
20 Time Sample
40 60
100
8
6
20 4
40 60
4
40 60
2
80 100 120
x 10 14
20
40
120
-7
WVD
x 10 14
2
80 100
50
100 Trace No.
150
(c)
200
120
50
100
150
200
(d)
Figure 5. 30 Hz iso-frequency section of (a) STFT (b) WVD (c) SPWVD (d) RID.
2.3 Seismic Time-Frequency Attributes Seismic attribute is a quantity extracted or derived from seismic data that leads to a better geological or geophysical interpretation of the data [40]. The study of seismic attributes is an effective way to gather some important information about the geometry and the physical parameters of the subsurface such as
reservoir thickness, faults, porosity, permeability and hydrocarbon reservoir properties. There are a large number of seismic attributes, and each of them has a special application in practice [38]. Time-frequency attribute is an important category of seismic attributes, which is widely used in seismic data interpretation. In exploration seismology, a buried channel can be considered as a thin layer. Impulse response of a thin bed shown in Figure 6 can be defined as [51]: ܴሺݐሻ ൌ ܴଵ ݓሺݐሻ ܶଵ ܴଶ ݓሺ ݐ ߬ሻሺሻ where ܴሺݐሻ is the reflection from thin bed, ݓሺݐሻis the seismic source wavelet, ܴଵ and ܴଶ are the reflection coefficients of top and bottom boundary of the thin bed, respectively. ܶଵ is the transmission coefficient of top boundary of the thin bed and ߬ is the temporal thickness of the thin bed. The transfer function [52] can be computed by transforming the Eq. (7) into frequency domain using Fourier transform as: ܴሺ݆߱ሻ ൌ ൫ܴଵ ݁ ିఠ௧ ܶଵ ܴଶ ݁ ିఠሺ௧ାఛሻ ൯ܹሺ݆߱ሻ ൌ ܭሺ݆߱ሻܹሺ݆߱ሻሺͺሻ where ܹሺ݆߱ሻis the Fourier transform of߱ሺݐሻ and ܭሺ݆߱ሻ is the Fourier transform of the reflectivity series. It should be noted that the thin bed acts as a filter and changes the spectrum of incident waves. If ܴଵ ൌ െܴଶ ൌ ͳ and the transmission coefficient ܶଵ is close to 1, frequency characteristics of reflectivity series would be as Eq. (9) and its amplitude spectrum is defined as Eq. (10). ܭሺ݆߱ሻ ൌ ൫ͳ െ ݁ ିఠఛ ൯݁ ିఠ௧ ሺͻሻ ȁܭሺ݆߱ሻȁ ൌ ඥʹሾͳ െ ܿݏሺ߱߬ሻሿሺͳͲሻ This condition occurs when the velocity of the channel infills is larger than surrounding formation. In this case, destructive interference occurs at the edges of channels, so the minimum amplitude can be
expected. Due to reduction of the destructive interference while moving toward the center of the channel, the amplitude increases [51]. If the velocity of the channel infills is larger than the upper layer and smaller than the bottom layer, then ܴଵ ൌ ܴଶ ൌ ͳ and Eq. (9) and (10) can be rewritten as Eq. (11) and (12). ܭሺ݆߱ሻ ൌ ൫ͳ ݁ ିఠఛ ൯݁ ିఠ௧ ሺͳͳሻ ȁܭሺ݆߱ሻȁ ൌ ඥʹሾͳ ܿݏሺ߱߬ሻሿሺͳʹሻ In this case, unlike the previous case, constructive interference occurs at the edges of channels, so the maximum amplitude can be expected. Due to increasing of the destructive interference while moving toward the center of the channel, the amplitude decreases. In both cases, the reflected waves from the top and base of the channel interfer with each other in such a way that the dominant frequency will be maximum at the edge of channel. Also, dominant frequency decreases while moving toward the center of the channel.
Figure 6. Thin bed model.
Due to the high sensitivity of the seismic data frequency to channel thickness, we used seismic timefrequency attributes to reveal the buried channels which are one of the important oil traps. Amplitude of instantaneous centroid frequency, amplitude of instantaneous bandwidth frequency, amplitude of RMS
frequency, peak frequency and peak amplitude are effective time-frequency attributes to highlight stratigraphic features such as channels. These kinds of attributes reveal the variation of frequency content of seismic data. Therefor, the thickness of different parts of a channel can be simply detected. The instantaneous center frequency, instantaneous bandwidth frequency and RMS frequency are defined as Eq. (13) to Eq. (15), respectively [53]. ஶ
݂ ሺݐሻ ൌ
݂ ൈ ܴܶܨሺݐǡ ݂ሻ݂݀ ஶ
ܴܶܨሺݐǡ ݂ሻ݂݀
ሺͳ͵ሻ
ஶ
݂ଶ ሺݐሻ
ൌ
ሺ݂ െ ݂ ሺݐሻሻଶ ൈ ȁܴܶܨሺݐǡ ݂ሻȁଶ ݂݀ ஶ
ȁܴܶܨሺݐǡ ݂ሻȁଶ ݂݀
ሺͳͶሻ
ஶ
݂ଶ ሺݐሻ
ൌ
݂ ଶ ൈ ܴܶܨሺݐǡ ݂ሻ݂݀ ஶ
ܴܶܨሺݐǡ ݂ሻ݂݀
ሺͳͷሻ
where ܴܶܨሺݐǡ ݂ሻ is the time-frequency distribution. In this paper, the pointed out seismic time-frequency attributes are extracted from RID and SPWVD to identify the buried channels in synthetic and real seismic data. To identify the buried channel in seismic data, we must apply time-frequency decomposition on seismic data and then extract the aforementioned attributes.
3. Results To show the location of buried channel, we must apply time-frequency decomposition on seismic data and extract the aforementioned attributes. First, we applied the RID and SPWVD methods to each trace of 3D seismic data and computed considered attributes. Secondly, we located the resulting attributes in the same place of the considered trace; we repeated this action for the rest of traces. Therefore, we obtained 3D data of seismic attributes of the same size as the original 3D seismic data. Thirdly, we extracted time slices of a target horizon of 3D seismic attributes, which show the spread of channel. The related flow chart is illustrated in Figure 7.
Start
Take the seismic trace
Compute time-frequency decomposition
Compute seismic attribute
Place seismic attribute in same location as seismic trace
End
Figure 7. Flow chart of computing seismic attributes.
To show the efficiency of this method, we produced a synthetic sand-filled channel which is embedded in shale. The ratio of the reflection coefficients of top and bottom of the channel is considered as
ோభ ோమ
ൌ
െͳ. The 30 Hz Ricker wavelet is chosen for forward modeling as seismic source wavelet and the time sampling interval is 0.004 s. 3D seismic data of the channel and the time slice associated with the channel at time sample 41 are shown in Figure 8. The red oval shows the channel position on 3D seismic data. The channel spreads straightly through the Earth, which is clear in the displayed attributes. Peak amplitude and peak frequency attributes extracted from RID and SPWVD are shown in Figure 9. As can be seen in Figure 10 the amplitude of instantaneous center frequency, instantaneous bandwidth frequency and RMS frequency of RID and SPWVD clearly show the expansion of the channel, but there is an absolute difference between the resolution of RID and SPWVD.
X lin e No .
20 40 60 80 100 120 20
40
60 Inline No.
(a)
80
100
(b) Figure 8. (a) 3D synthetic seismic cube (b) Horizontal slice at time sample 41.
Magnitude
(b)
(c)
(d)
Frequency
(a)
Figure 9. Horizontal slice at the target zone (time sample 41) from (a) peak amplitude cube of RID and (b) peak amplitude cube of SPWVD and (c) peak frequency cube of RID and (d) peak frequency cube of SPWVD.
120
Amplitude of instantaneous center frequency of SPWVD 0.2 20 40 60 0.1 80 100 120 0 20 40 60 80 100 120
(a)
(b)
Amplitude of instantaneous bandwidth frequency of RID 20 40 60 80 100 120
0.2 0.1
20
40
60
(c)
80
100
120
0
Amplitude of instantaneous bandwidth frequency of SPWVD 20 40 60 80 100 120
Xline No.
0.1 0
(e)
40
60
80
100
120
0
Amplitude of RMS frequency of SPWVD
0.2
60 80 Inline No.
0.1
(d)
20 40 60 80 100 120 40
0.2
20
Amplitude of RMS frequency of RID
20
Magnitude
Amplitude of instantaneous center frequency of RID 0.2 20 40 60 0.1 80 100 120 0 20 40 60 80 100 120
100
120
20 40 60 80 100 120
0.2 0.1
20
40
60
80
100
120
0
(f)
Figure 10. Horizontal slice at the target zone (time sample 41) from instantaneous center frequency of (a) RID (b) SPWVD, Horizontal slice at the target zone (time sample 41) from instantaneous bandwidth frequency of (c) RID (d) SPWVD and Horizontal slice at the target zone (time sample 41) from RMS frequency of (e) RID (f) SPWVD.
Real seismic data is used to recognize the efficiency of RID to identify hidden channels by the above mentioned seismic time-frequency attributes. The real 3D seismic data is taken from one of the oil fields in the south-west of Iran with 400 inline and 600 crossline traces, with time sampling interval 0.004 s, which contains a channel at the 1.8 sec (time sample 10) time slice based on previous study [54]. In Figure 11, we show the 3D seismic cube and an extracted time slice at the 1.8 sec. As can be seen, it is difficult to recognize the buried channel in seismic time slice. The three seismic time-frequency attributes based on RID, the amplitude of instantaneous center frequency, amplitude of instantaneous bandwidth frequency and amplitude of RMS frequency attributes are shown in Figure 12. We extracted the time slices at 1.8 s from each of the three attributes and compared the result with that of SPWVD in Figures
13-15. Due to high resolution of RID, the different branches of buried channel (even thin branches) are more visible in the results of RID than SPWVD.
(a)
(b)
Figure 11. (a) 3D real seismic cube and (b) time slice at the 1.8 sec.
(b)
(a)
(c)
Figure 12. 3D cube of seismic time-frequency attributes extracted from RID (a) amplitude of instantaneous center frequency (b) amplitude of instantaneous bandwidth frequency (c) amplitude of RMS frequency.
Amplitude of instantaneous center frequency of RID 600
600
0.3 0.25
0.25
500
400
0.2
400
0.2
300
0.15
300
0.15
200
0.1
200
0.1
100
0.05
100
0.05
50
100
150
200 250 Xline No.
(a)
300
350
400
0
50
100
150
200
250
300
350
400
Magnitude
Inline No.
500
Amplitude of instantaneous center frequency of SPWVD 0.3
0
(b)
Figure 13. Time slice at the target zone (time sample 10) from amplitude of instantaneous center frequency extracted from (a) RID (b) SPWVD.
Amplitude of instantaneous bandwidth frequency of SPWVD
Amplitude of instantaneous bandwidth frequency of RID 600
0.4
600
500
0.35
500
0.4 0.35 0.3
0.3 Inline No.
400
0.25
0.25 300
0.2
300
0.2
200
0.15
200
0.15 0.1
0.1 100
Magnitude
400
100 0.05
0.05 50
100
150
200 250 Xline No.
300
350
400 0
50
100
150
200
250
300
350
0
400
(b)
(a)
Figure 14. Time slice at the target zone (time sample 10) from amplitude of instantaneous bandwidth frequency extracted from (a) RID (b) SPWVD.
Amplitude of RMS frequency of RID
Amplitude of RMS frequency of SPWVD
600
0.5
500
0.4
Inline No.
0.5
500
0.4
Magnitude
400
600
400
0.3
0.3 300
300
0.2
0.2 200
200 0.1
100
50
100
150
200 250 Xline No.
300
350
400
0
0.1
100
50
100
150
200
250
300
350
400
0
(b)
(a)
Figure 15. Time slice at the target zone (time sample 10) from amplitude of RMS frequency extracted form (a) RID (b) SPWVD.
We calculated the peak amplitude and peak frequency 3D seismic attributes using RID which are shown in Figure 16. Time slices at 1.8 s from two seismic attributes computed based on two time-frequency distribution (RID and SPWVD) are shown in Figure 17. As can be seen, low peak frequency corresponds to thicker part of the channel and vice versa. Conversely, high peak amplitude corresponds to thicker part of the channel. Therefore, we can identify the edge of channel and variation of its thickness.
(b)
(a)
Figure 16. (a) Peak amplitude (b) peak frequency attributes extracted from RID.
Peak Amplitude of RID 600 500
Peak Amplitude of SPWVD
1
600
0.8
500
0.8
400 0.6
300
0.6 300
0.4
200 100
0.2 50
100
150
200
250
300
350
200
0.4
100
0.2
400
50
100
150
(a)
250
300
350
400
(b) Peak Frequency of SPWVD
Peak Frequency of RID 600
60
600
60 50
500
50
500
400
40
400
40 30
30
300
200
20
200
20
100
10
100
10
300
50
100
150
200
250
Xline No.
(c)
300
350
400
0
50
100
150
200
250
300
350
Frequency
Inline No.
200
Magnitude
400
1
400
(d)
Fig. 17. Time slice at the target zone (time sample 10) from (a) peak amplitude cube extracted from and (b) peak amplitude cube extracted from SPWVD. (c) peak frequency cube extracted from and (d) peak frequency cube extracted from SPWVD.
4. Conclusions Frequency content of seismic signals carries valuable information of subtle changes, so analysis of this information is an effective way to survey subsurface structures. A spectral decomposition is a powerful tool in this field. Seismic signals are non-stationary multi-component signals and by employing conventional methods such as WVD, PWVD and SPWVD, some problems such as Heisenberg
uncertainty principle and cross–terms emerged. In this paper, we introduced the Reduced Interference Distribution to overcome the restrictions of the aforementioned methods. The RID keeps almost all the desirable properties of the WVD, but reduces considerably the interference cross–terms. We used this method of spectral decomposition to calculate the instantaneous central frequency, instantaneous bandwidth frequency, RMS frequency, peak frequency and peak amplitude attributes, which are used for imaging of the buried channel which is one of the oil traps. Actually, we wanted to improve the results of these attributes by RID. We tested this method of buried channel imaging on a synthetic and a real seismic data. We saw a reasonably good correlation between the channel's thickness in time and the mentioned seismic attributes. Therefore, these attributes can clearly delineate spreading channels vertically and horizontally, and correlate well with channel thickness.
References [1] [2] [3]
[4] [5]
[6] [7] [8] [9]
[10]
J. Pujol, Elastic Wave Propagation and Generation in Seismology, Cambridge University Press, 2003. Ö. Yilmaz, S.M. Doherty, Seismic Data Analysis: processing, inversion, and interpretation of seismic data, Society of Exploration Geophysicists, 2001. M.J. Duchesne, E.J. Halliday, J.V. Barrie, Analyzing seismic imagery in the time–amplitude and time–frequency domains to determine fluid nature and migration pathways: A case study from the Queen Charlotte Basin, offshore British Columbia, Journal of Applied Geophysics 73 (2011) 111-120. G. Partyka, J. Gridley, J. Lopez, Interpretational applications of spectral decomposition in reservoir characterization, The Leading Edge 18 (1999) 353-360. M.C. de Matos, Characterization of thin beds through joint time-frequency analysis applied to a turbidite reservoir in Campos Basin, Brazil, SEG Technical Program Expanded Abstracts 2005, pp. 1429-1432. J. Liu, K.J. Marfurt, Instantaneous spectral attributes to detect channels, Geophysics 72 (2007) P23-P31. R.H. Herrera, J. Han, M.v.d. Baan, Applications of the synchrosqueezing transform in seismic time-frequency analysis, GEOPHYSICS 79 (2014) V55-V64. J.P. Castagna, S. Sun, R.W. Siegfried, Instantaneous spectral analysis: Detection of low-frequency shadows associated with hydrocarbons, The Leading Edge 22 (2003) 120-127. M. Zarei, A. Roshandel Kahoo, H.R. Siahkoohi, Gas detection using deconvolutive short time Fourier transform, International Geophysical Conference and Oil & Gas Exhibition, Istanbul, Turkey, 2012, pp. 1-4. H. Sattari, A. Gholami, H.R. Siahkoohi, Seismic data analysis by adaptive sparse time-frequency decomposition, Geophysics 78 (2013) V207-V217.
[11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23]
[24]
[25]
[26] [27] [28] [29]
[30] [31]
[32] [33]
S. Sinha, P.S. Routh, P.D. Anno, J.P. Castagna, Spectral decomposition of seismic data with continuous-wavelet transform, Geophysics 70 (2005) P19-P25. R. Askari, H.R. Siahkoohi, Ground roll attenuation using the S and x-f-k transforms, Geophysical Prospecting 56 (2008) 105-114. A. Roshandel Kahoo, H.R. Siahkoohi, Gas detection from AVO analysis in time-frequency domain, 71st EAGE Conference & Exhibition, 2009. A. Roshandel Kahoo, H.R. Siahkoohi, Random noise suppression from seismic data using timefrequency peak filtering, 71st EAGE Conference & Exhibition, 2009. G.F. Margrave, M.P. Lamoureux, D.C. Henley, Gabor deconvolution: Estimating reflectivity by nonstationary deconvolution of seismic data, Geophysics 76 (2011) W15-W30. X.-K. Sun, Z. Sam, H.-W. Xie, Nonstationary sparsity-constrained seismic deconvolution, Applied Geophysics 11 (2014) 459-467. D. Gabor, Theory of communication. Part 1: The analysis of information, Electrical Engineers Part III: Radio and Communication Engineering, Journal of the Institution of 93 (1946) 429-441. A. Gholami, Sparse Time-Frequency Decomposition and Some Applications, Geoscience and Remote Sensing, IEEE Transactions on 51 (2013) 3598-3604. E. Wigner, On the quantum correction for thermodynamic equilibrium, Physical Review 40 (1932) 749. B. Boashash, Time Frequency Analysis, Elsevier Science, 2003. S. Mallat, A Wavelet Tour of Signal Processing: The Sparse Way, Elsevier Science, 2008. R.G. Stockwell, L. Mansinha, R. Lowe, Localization of the complex spectrum: the S transform, Signal Processing, IEEE Transactions on 44 (1996) 998-1001. M. Radad, A. Gholami, H.R. Siahkoohi, S-transform with maximum energy concentration: Application to non-stationary seismic deconvolution, Journal of Applied Geophysics 118 (2015) 155-166. F. Auger, P. Flandrin, Improving the readability of time-frequency and time-scale representations by the reassignment method, Signal Processing, IEEE Transactions on 43 (1995) 1068-1089. H.I. Choi, W.J. Williams, Improved time-frequency representation of multicomponent signals using exponential kernels, Acoustics, Speech and Signal Processing, IEEE Transactions on 37 (1989) 862-871. J. Jechang, W.J. Williams, Kernel design for reduced interference distributions, Signal Processing, IEEE Transactions on 40 (1992) 402-412. L. Cohen, Time-frequency Analysis, Prentice Hall PTR, 1995. M. Rauch-Davies, M. Ralston, Spectral decomposition–transform methods and fluid and reservoir prediction case study, 67th EAGE Conference & Exhibition, 2005. M. Ralston, M. Rauch-Davies, K. Li-Chun, X. Hui-Ping, Y. Di-Sheng, General method to reduce crossterm interference in the Wigner-Ville decomposition, SEG Technical Program Expanded Abstracts 2007, pp. 870-874. Y. Li, X. Zheng, Wigner-Ville distribution and its application in seismic attenuation estimation, Applied Geophysics 4 (2007) 245-254. Y. Hu, K.D.K. Luk, W.W. Lu, A. Holmes, J.C.Y. Leong, Comparison of time-frequency distribution techniques for analysis of spinal somatosensory evoked potential, Medical and Biological Engineering and Computing 39 (2001) 375-380. W.J. Williams, Reduced interference distributions: biological applications and interpretations, Proceedings of the IEEE 84 (1996) 1264-1280. M.B. Widess, How thin is a thin bed?, Geophysics 38 (1973) 1176-1180.
[34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48]
[49] [50]
[51] [52] [53] [54]
J.D. Robertson, H.H. Nogami, Complex seismic trace analysis of thin beds, Geophysics 49 (1984) 344-352. H. Zeng, A. John, G. Katherine, How thin is a thin bed? An alternative perspective, SEG Technical Program Expanded Abstracts, 2008, pp. 834-838. H. Zeng, Thin-bed prediction by geomorphology-constrained waveform analysis, SEG Technical Program Expanded Abstracts 2014, pp. 1395-1399. G. Li, M.D. Sacchi, S. Zhu, Y. Wang, Thin bed identification using zero-crossing-time stratal slices, SEG Technical Program Expanded Abstracts 2015, pp. 1678-1682. S. Chopra, K.J. Marfurt, Seismic Attributes for Prospect Identification and Reservoir Characterization, Society of Exploration Geophysicists, 2007. X. Wu, T. Liu, Spectral decomposition of seismic data with reassigned smoothed pseudo Wigner– Ville distribution, Journal of Applied Geophysics 68 (2009) 386-393. M.T. Taner, F. Koehler, R. Sheriff, Complex seismic trace analysis, Geophysics 44 (1979) 10411063. B. Boashash, Time-Frequency Signal Analysis and Processing: A Comprehensive Reference, Elsevier Science, 2015. H. Hassanpour, Time-frequency based detection of newborn EEG seizure, Queensland University of Technology, 2004. R.E. Sheriff, L.P. Geldart, Exploration Seismology, Cambridge University Press, 1995. M.R. Gadallah, R. Fisher, R.L. Fisher, Exploration Geophysics, Springer, 2008. W.M. Telford, L.P. Geldart, R.E. Sheriff, Applied Geophysics, Cambridge University Press, 1990. D.A. Herron, R.B. Latimer, First Steps in Seismic Interpretation, Society of Exploration Geophysicists, 2011. A. Cordsen, M. Galbraith, J. Peirce, B.A. Hardage, Planning Land 3-D Seismic Surveys, Society of Exploration Geophysicists, 2000. A.R. Brown, Interpretation of three-dimensional seismic data, 7 ed., Published jointly by American Association of Petroleum Geologists and the Society of Exploration Geophysicists, 2011. D.L. Baars, W.L. Watney, D.W. Steeples, E.A. Brostuen, Petroleum: a primer for Kansas, Educational Series, 1993. N. Ricker, Wavelet contraction, wavelet expansion, and the control of seismic resolution, Geophysics 18 (1953) 769–792. R.D. Han, Z.H. Wan, M.S. Chen, H.Y. Zhang, Application of time-frequency attributes based on generalized s-transform for thin bed indication, 73rd EAGE Conference & Exhibition, 2011. J.G. Proakis, D.G. Manolakis, Digital Signal Processing, Principles, Algorithms, and Applications, 4 ed., Pearson Prentice Hall, 2007. A.E. Barnes, Instantaneous spectral bandwidth and dominant frequency with applications to seismic reflection data, Geophysics 58 (1993) 419-428. R. Mohebian, M. Yari, M.A. Riahi, R. Ghanati, Channel detection using instantaneous spectral attributes in one of the SW Iran oil fields, Bollettino di Geofisica Teorica ed Applicata 54 (2013) 271-282.
Azita Nikoo received the B.S. degree in physics from Razi University of Kermanshah, Kermanshah, Iran, in 2009, the M.S. degree in Geophysics from Shahrood University of Technology, Shahrood, Iran, in 2013. She is the project implementers in Gity Kavosh Zaniar Consultants Company, Iran. Her research interests include application of time-frequency transform in exploration seismology, seismic attributes, and AVO. Hamed Saadatnia received the B.S. degree in Geology from University of Tehran, Tehran, Iran, in 2001 and the M.S. degrees in exploration seismology from the University of Tehran, Tehran, Iran, in 2004. Since 2007, he has been with the National Iranian Oil Company (NIOC), Tehran, Iran, where he is currently a geophysicist in seismic interpretation department. His research interests include reservoir characterization and seismic data interpretation.
Hamid Hassanpour received the B.S. degree in computer engineering from Iran University of Science and Technology, Tehran, Iran, in 1993, the M.S. degree in computer engineering from Amirkabir University of Technology, Tehran, Iran, in 1996, and the Ph.D. from the Queensland University of Technology, Brisbane, Australia, in 2004. He has a professor position in faculty of Computer Engineering & IT at the Shahrood University of Technology, Iran. His research interests include Image Processing, Signal Processing, time-frequency signal processing and analysis. Amin Roshandel Kahoo received the B.S. degree in mining engineering from University of Birjand, Birjand, Iran, and the M.S. and Ph.D. degrees in exploration seismology from the University of Tehran, Tehran, Iran, in 2005 and 2010 respectively. Since September 2010, he has been with the School of Mining, Petroleum & Geophysics Engineering, Shahrood University of Technology, Shahrood, Iran, where he is currently an Assistant Professor. His research interests include application of time-frequency transform in exploration seismology, seismic attributes, and seismic data processing and interpretation.