Month: November 2011

RNA Interpolation

November 21, 2011 Documentation No comments

A new paper is added to the collection of reproducible documents:
Seismic data interpolation beyond aliasing using regularized nonstationary autoregression

Fourier finite differences

November 16, 2011 Documentation No comments

A new paper is added to the collection of reproducible documents:
Fourier finite-difference wave propagation

Program of the month: sfenvelope

November 5, 2011 Programs No comments

Complex trace attributes were introduced into geophysics by the paper

Taner, M. T., F. Koehler, and R. E. Sheriff, 1979, Complex seismic trace analysis: Geophysics, 44, 1041-1063.

If $s(t)$ is the input seismic trace, then the analytical trace is defined as the complex-valued signal $a(t) = s(t)+i h(t)$, where $h(t)$ is the Hilbert transform of $s(t)$

$$h(t) =\displaystyle \frac{1}{\pi} \int \frac{s(\tau)}{t-\tau} d\tau\;.$$

The signal envelope is the positive signal $e(t)=\sqrt{s^2(t)+h^2(t)}$. A phase-rotated seismic signal is $p(t)=s(t)\,\cos{\phi} +h(t)\,\sin{\phi}$ where $\phi$ is the phase of rotation. By default, sfenvelope computes the signal envelope. It can also produce a phase-rotated signal if given hilb=y and phase=. If phase=90 (the default value), the phase-rotated signal will be simply the Hilbert transform of the input.

The following figure from book/rsf/rsf/sfenvelope illustrates an application of sfenvelope:

Computing the discrete Hilbert transform is not a trivial task. In the Fourier domain, the continuous Hilbert transform is given by

$$\displaystyle H(\omega) = i\,\operatorname{sgn}(\omega)\,S(\omega)$$

where $\operatorname{sgn}$ is the sign function. The discontinuity of the sign function in the frequency domain at $\omega=0$ is related to the slow $1/t$ decay of the filter impulse response in the time domain. The discontinuity at the Nyquist frequency creates additional oscillations. Different practical implementations shorten the filter impulse response by effectively smoothing the Fourier-domain discontinuities. The Madagascar implementation of the discrete Hilbert transform follows the algorithm described in

Pei, S.-C., and P.-H. Wang, 2001, Closed-form design of maximally flat FIR Hilbert transformers, differentiators, and fractional delayers by power series expansion: IEEE Trans. on Circuits and Systems, v. 48, No. 4, 389-398.

The accuracy/cost trade-off is controlled by two parameters: order= and ref=. The following figures frombook/rsf/rsf/sfenvelope illustrate the effect of the order= parameter:

The Seismic Unix implementation (suhilb program) applies a Hamming window in the time domain. For some reason, it has the filter polarity reversed:

A multidimensional analog of the Hilbert transform is the Riesz transform. It is implemented in the sfriesz program.


November 5, 2011 Documentation No comments

New users often ask for a brief comprehensible introduction to Madagascar. Jeff Gowdin contributes A brief introduction to Madagascar, a tutorial book designed to provide such an introduction.

The purpose of this document is to teach new users how to use the powerful tools in Madagascar to: process data, produce documents and build your own Madagascar programs.

You can access the tutorial in the following forms: