9.2 Fourier analysis and filtering

Fourier Analysis in Geophysics

Principles and Applications

Fourier analysis works on a simple but powerful idea: any complex signal can be broken down into a sum of sinusoidal (sine and cosine) waves, each with its own frequency, amplitude, and phase. In geophysics, the signals you're working with are rarely simple. A seismic trace, for instance, contains reflections from multiple layers, noise from the surface, and instrument artifacts all stacked on top of each other. Fourier analysis lets you pull those apart by moving from the time domain (amplitude vs. time) into the frequency domain (amplitude vs. frequency).

The Fourier transform is the mathematical operation that performs this conversion. It comes as a pair:

• The forward Fourier transform converts from time to frequency domain
• The inverse Fourier transform converts back from frequency to time domain

This back-and-forth capability is what makes Fourier analysis so useful. You can transform your data into the frequency domain, manipulate specific frequency components (remove noise, isolate a signal band), and then transform back to get a cleaner time-domain signal.

Fourier analysis applies across nearly every branch of geophysical data processing:

• Signal processing of seismic, gravity, magnetic, and electromagnetic data
• Noise reduction by identifying and removing unwanted frequency components
• Spectral analysis to characterize the frequency content of a dataset
• Data compression by retaining only the most significant frequency components

Mathematical Foundations

The continuous forward Fourier transform converts a time-domain function $f(t)$ into its frequency-domain representation $F(\omega)$:

$F(\omega) = \int_{-\infty}^{\infty} f(t) \, e^{-i\omega t} \, dt$

The inverse Fourier transform recovers the original signal:

$f(t) = \frac{1}{2\pi} \int_{-\infty}^{\infty} F(\omega) \, e^{i\omega t} \, d\omega$

Here, $\omega$ is angular frequency (in radians per second), $i$ is the imaginary unit, and the complex exponential $e^{-i\omega t}$ encodes both sine and cosine components via Euler's formula.

Several properties make the Fourier transform especially practical:

• Linearity: The transform of a sum of signals equals the sum of their individual transforms.
• Scaling: Compressing a signal in time spreads it out in frequency, and vice versa.
• Convolution theorem: Convolution in the time domain becomes simple multiplication in the frequency domain. This is why filtering is so much more efficient in the frequency domain. Instead of performing a computationally expensive convolution, you just multiply two spectra together.

Time vs. Frequency Domain Conversion

Discrete Fourier Transform (DFT)

Geophysical data is recorded as discrete samples, not continuous functions. The discrete Fourier transform (DFT) handles this by operating on a finite set of $N$ samples.

The forward DFT is:

$X[k] = \sum_{n=0}^{N-1} x[n] \, e^{-i\frac{2\pi}{N}kn}$

The inverse DFT is:

$x[n] = \frac{1}{N} \sum_{k=0}^{N-1} X[k] \, e^{i\frac{2\pi}{N}kn}$

where $x[n]$ is the $n$-th time-domain sample and $X[k]$ is the $k$-th frequency-domain coefficient.

Computing the DFT directly requires $O(N^2)$ operations, which becomes prohibitively slow for large datasets. The fast Fourier transform (FFT) is an algorithm that exploits symmetries in the DFT to reduce this to $O(N \log N)$. For a seismic survey with millions of samples, this difference is the reason frequency-domain processing is feasible at all. Most FFT implementations work most efficiently when $N$ is a power of 2, so data is often zero-padded to the next power of 2 before transformation.

Sampling and Aliasing

Before any Fourier analysis can happen, the continuous geophysical signal must be properly sampled and digitized. The Nyquist-Shannon sampling theorem sets the fundamental rule: the sampling frequency $f_s$ must be at least twice the highest frequency component $f_{max}$ present in the signal.

$f_s \geq 2 f_{max}$

The frequency $f_s / 2$ is called the Nyquist frequency. If the signal contains energy above this frequency, those components get "folded back" into lower frequencies during digitization. This is aliasing, and it's irreversible once the data is recorded. A 120 Hz signal sampled at 200 Hz, for example, would appear as an 80 Hz signal in the digitized data.

To prevent aliasing, anti-aliasing filters (analog low-pass filters) are applied before digitization to remove frequency content above the Nyquist frequency.

Power and Phase Spectra

The Fourier transform $F(\omega)$ is complex-valued, meaning it contains both magnitude and phase information. These are typically analyzed separately:

• The power spectrum is the squared magnitude $|F(\omega)|^2$. It shows how energy is distributed across frequencies. A seismic trace dominated by low-frequency reflections will show most of its power at low frequencies, while high-frequency noise will appear as elevated power at the high end.
• The phase spectrum is the argument (angle) of $F(\omega)$. It tells you the relative timing of each frequency component. Phase information is critical for preserving the shape and timing of waveforms during processing.

Together, the power and phase spectra fully characterize a signal's frequency content. You can reconstruct the original signal perfectly from these two spectra using the inverse Fourier transform.

Frequency-Domain Filtering

Types of Filters

Filtering in the frequency domain means selectively modifying specific frequency components of a signal. The four main filter types are:

• Low-pass filters preserve frequencies below a cutoff and attenuate those above it. These are used to remove high-frequency noise or smooth data. In seismic processing, a low-pass filter might remove random noise above 80 Hz while keeping the reflection signal below that.
• High-pass filters preserve frequencies above a cutoff and attenuate those below it. These remove low-frequency trends such as DC offsets or long-period drift. A high-pass filter at 5 Hz, for example, could remove slow instrument drift from a seismic record.
• Band-pass filters preserve frequencies within a specified range and attenuate everything outside it. These are the most commonly used filters in seismic processing. A typical band-pass might keep frequencies between 10 and 60 Hz, isolating the useful reflection bandwidth.
• Notch filters attenuate a narrow frequency band while preserving everything else. The classic application is removing 50 or 60 Hz power line interference from electromagnetic or seismic data.

Filter Design and Implementation

Designing a frequency-domain filter involves several choices:

1. Define the filter response by specifying cutoff frequencies, the width of the transition band (how sharply the filter rolls off), and the required stopband attenuation.
2. Choose a filter type. Common options include Butterworth (maximally flat passband), Chebyshev (sharper rolloff but with passband ripple), and elliptic filters (sharpest rolloff but with ripple in both passband and stopband).
3. Apply the filter. In practice, you multiply the signal's Fourier transform $F(\omega)$ by the filter's frequency response $H(\omega)$ to get the filtered spectrum $G(\omega) = F(\omega) \cdot H(\omega)$.
4. Transform back. Apply the inverse Fourier transform to $G(\omega)$ to obtain the filtered time-domain signal.

This multiplication-in-frequency approach is a direct application of the convolution theorem: multiplying spectra in the frequency domain is equivalent to convolving the signal with the filter's impulse response in the time domain, but computationally much faster.

Effects of Filtering on Geophysical Data

Signal-to-Noise Ratio Improvement

The primary goal of filtering is to improve the signal-to-noise ratio (SNR) by suppressing unwanted components:

• Low-pass filtering removes high-frequency random noise and measurement errors
• High-pass filtering removes low-frequency trends, DC offsets, or long-period drift
• Band-pass filtering isolates the frequency range where the desired signal is strongest, attenuating noise on both sides
• Notch filtering targets specific interference sources like power line noise or ground roll in seismic data

Potential Artifacts and Distortions

Every filter involves trade-offs. Applying a filter always changes the data, and not always in desirable ways:

• Low-pass filtering can blur sharp features and edges, reducing spatial or temporal resolution. Thin-bed reflections in seismic data, for instance, may be smeared together.
• High-pass filtering can amplify high-frequency noise, making the data appear grainy or noisy.
• Band-pass filtering can introduce Gibbs phenomenon (ringing artifacts) near sharp transitions or discontinuities, especially with steep filter rolloffs.
• Notch filtering can remove desired signal energy if the notch is too wide or poorly centered on the interference frequency.

Interpretation Considerations

Filtered data should always be compared against the original unfiltered data so you can assess what the filter actually did. A few practical guidelines:

• Choose filter parameters based on the known frequency content of both the signal and the noise. Examining the power spectrum before filtering helps you make informed choices.
• Consider the effects on amplitude, phase, and waveform shape. Zero-phase filters preserve event timing but aren't always achievable in practice; minimum-phase filters are causal but introduce phase distortion.
• Avoid over-filtering. Applying too many filters or using overly aggressive parameters can remove real geological signal along with the noise.
• Document and communicate your filter choices. Anyone interpreting the processed data needs to know what filters were applied, with what parameters, and what artifacts might have been introduced.