Phase is a useful underlying property of the analytic trace model of seismic data that can be used as both an interpretation aid and a means to calibrate and check interpretations on a given seismic dataset. We introduce the analytical trace model and demonstrate some of its usages. We provide working code in python for computation of the Hilbert Transform using a robust FFT-based method and explore 2 use cases for such computed quantities. Jupyter notebooks used for computation and generation of the figures are included in this project.

Introduction¶

The concept of phase permeates seismic data processing and signal processing in general, but it can be awkward to understand, and manipulating it directly can lead to surprising results. It doesn't help that the word phase is used to mean a variety of things, depending on whether we refer to the propagating wavelet, the observed wavelet, post-stack seismic attributes, or an entire seismic data set. Several publications have discussed the concepts and ambiguities Chant & Hastie, 1992Roden & Sepúlveda, 1999Liner, 2002.

The data we are using are the Penobscot 3D Purves, 2014 survey from the Open Seismic Repository, which is openly licensed CC-BY-SA by the Canada Nova Scotia Offshore Petroleum Board, and dGB Earth Sciences.

In this tutorial, we will focus on aspects of phase relevant to the interpreter. We will look at how to manipulate the phase of a seismic trace by manipulating the phase of its Fourier transform and will use that idea to generate the well-known instantaneous phase post-stack attribute. We will also check how close to zero phase our test data set is.

The plots included in this tutorial were created using standard python libraries and the code to reproduce them is available here on Curvenote and on the SEG tutorials GitHub repo.

The Hilbert Transform¶

The Fourier transform is complex. Taking the transform of any real signal will result in a set of complex coefficients. Complex numbers are essentially 2D vectors, meaning they have two components: magnitude and phase angle. Most of the time when dealing with Fourier transforms, we concentrate on the magnitude, which tells us about the distribution of signal energy through frequency. However, every signal also has a phase spectrum, and the phase encodes the signal's structure --- the distribution of the signal energy through time. We do not often examine phase spectrum because it is difficult to interpret, but we can manipulate the Fourier phase to change the structure of our signal without affecting its amplitude spectrum.

The Hilbert transform is a linear operator that produces a 90˚ phase shift in a signal, and it is a good first step in our exploration of phase. It is also commonly used in post-stack seismic analysis to generate the analytic signal from which we can compute the standard complex trace attributes such as envelope, instantaneous phase, and instantaneous frequency

The definition of the Hilbert transform is rather cryptic; it is much easier to consider in terms of its Fourier transform definition. The Hilbert transform $H$ of a signal $u$ is related to the Fourier transform $F$ like this:

F(H(u))(\omega) = \sigma_H(\omega) \cdot F(u)(\omega)

(1)

where:

\sigma_{\mathrm{H}}(\omega)=\left\{\begin{array}{l}i \text { for } \omega<0 \\ 0 \text { for } \omega=0 \\ -i \text { for } \omega>0\end{array}\right.

(2)

So we apply a Hilbert transform by multiplying all negative frequencies by $i$ and all positive frequencies by $-i$ , leaving any DC component untouched. This is perhaps not intuitive, but we can gain some additional insight by thinking about it geometrically.

The multiplication of two complex numbers $a$ (such as a Fourier coefficient) and $b$ (such as $\sigma_H$ above) can be thought of as a rotation in the complex plane. Euler's formula helps us to see that to rotate any complex number by a specific angle a, we must multiply it by the complex number $\mathrm{e}^{i \alpha}=\cos \alpha+i \sin \alpha$ .

We see that multiplication by $i$ alone, as in the equation above, is equivalent to a rotation by 90°. When we take the inverse Fourier transform, the result is a phase-shifted version of our signal.

So, we can generalize the definition of the Hilbert above to produce a phase shift to any angle, $\alpha$ :

\sigma_{n}(\omega)=\left\{\begin{array}{l}\mathrm{e}^{i \alpha} \text { for } \omega<0 \\ 0 \text { for } \omega=0 \\ \mathrm{e}^{-i \alpha} \text { for } \omega>0\end{array}\right.

(3)

Phase shifting in Python¶

Let us use this principle to phase-shift a seismic trace using open source python.

NOTE: If you are working in GNU Octave check back to fftshifter.m on github.

We will take the fast Fourier transform (FFT) of our signal, manipulate the coefficients, and apply the inverse FFT.

from math import ceil,exp
from scipy.fftpack import fft, ifft

def fftshifter(x, phase_shift_in_radians):
    # Create an array to apply the coefficient to 
    # positive and negative frequencies appropriately
    N = len(x);
    R0 = np.exp(-1j*phase_shift_in_radians);
    R = np.ones(x.shape, dtype=np.complex);
    R[0]=0.; 
    R[1:N//2] = R0;
    R[N//2:] = np.conj(R0);

    # Apply the phase shift in the frequency domain
    Xshifted = R*fft(x);
    
    # Recover the shifted time domain signal
    y = np.real(ifft(Xshifted));
    
    return y

Complex Trace Attributes¶

Now that we have this ability to phase-shift, what can we do with it?

We can now compute some seismic attributes. A whole set of commonly used attributes depends on the Hilbert transform and the code below generates the Hilbert transform, envelope, and phase traces for a trace from the Penobscot data set.

from scipy.signal import hilbert as hilbert_transform
hilbert = np.zeros_like(seismic)

for x in range(line.shape[0]):
    hilbert[x,:] = np.imag(hilbert_transform(line[x,:])[500:750])

z = seismic + 1j*hilbert
env = np.abs(z)
phase = np.angle(z)

See the accompanying notebook for full plotting code. A section of this trace is shown in Figure 2 below.

(2000.0, 3000.0)

Figure %s shows the same attributes now computed over the whole section. These are useful interpretation aids. Once you have the complex trace representation, z, there are many additional attributes that you can compute, such as instantaneous frequency, phase acceleration, and envelope-weighted frequency, the list goes on\dots merely by playing around with the complex trace.

Checking the phase of your seismic¶

The complex attributes generated from a seismic trace and its Hilbert transform are useful in seismic interpretation. However, the quality of our seismic interpretation depends on another aspect of phase, that is, the phase of our seismic data set and whether the observed wavelet has been corrected to a zero-phase response during processing. Roden & Sepúlveda (1999) give a practical way to assess the phase of seismic data as a quality check for seismic interpretation. They point out that where we have a strong reflection from a sharp geologic interface that we can confidently expect to produce a zero-phase response, we can measure the actual response in the seismic data set on that horizon to determine the amount of phase error present.

The principle that Roden & Sepúlveda (1999) apply is that at strong reflections in zero phase data, we would expect peaks (or troughs) in the trace to align with peaks in the envelope. Therefore, any horizon that we pick on zero-phase data could also be picked on the envelope of that data.

Figure %s shows a segment of our trace from Figure %s with envelope peaks marked as spikes (see find\_peaks.m). In Figure 2b, I have isolated those events and plotted the actual signal phase against the phase we would expect for a zero-phase data set.

The difference between the two could be systematic residual phase. If we think a particular reflector has minimal interference and no AVO or other effects, then the rest is residual phase. The strong trough at the Abenaki reflector in Figure suggests that there could be 0.64 radians, or about 37°, of residual phase in the data. That is something that we could now correct for by phase-shifting all our traces.

We can also look at this phase error at envelope peaks across our dataset as shown in Figure 5 below.

The Hilbert transform opens up a world of seismic attributes, some of which have everyday application for the interpreter. To run or build on the code clone the repo at http://github.com/seg/tutorials-2014.

The official version of this tutorial was published in the SEG's Leading Edge Magazine in October 2014.

References¶

Chant, I. J., & Hastie, L. M. (1992). Time-frequency analysis of magnetotelluric data. Geophysical Journal International, 111(2), 399–413. 10.1111/j.1365-246x.1992.tb00586.x
Roden, R., & Sepúlveda, H. (1999). The significance of phase to the interpreter: Practical guidelines for phase analysis. The Leading Edge, 18(7), 774–777. 10.1190/1.1438375
Liner, C. L. (2002). Phase, phase, phase. The Leading Edge, 21(5), 456–457. 10.1190/1.1885500
Purves, S. (2014). Phase and the Hilbert transform. The Leading Edge, 33(10), 1164–1166. 10.1190/tle33101164.1
Purves, S. (2014). Phase and the Hilbert transform. The Leading Edge, 33(10), 1164–1166. 10.1190/tle33101164.1

geoscience

phase_and_hilbert_transform.ipynb

Phase and the Hilbert Transform

Introduction¶

The Hilbert Transform¶

Phase shifting in Python¶

Complex Trace Attributes¶

Checking the phase of your seismic¶