SciPy - FFTpack

SciPy - Fourier Transforms (scipy.fft)

Fourier Transformation is a powerful mathematical tool that decomposes a signal (often a time-domain signal) into its constituent frequencies. Analyzing a signal in the frequency domain can reveal characteristics hidden in the time domain. It’s widely used in signal processing, image analysis, audio engineering, and many other scientific and engineering fields. SciPy provides comprehensive Fourier Transform capabilities through the scipy.fft module (which is the modern replacement for the older scipy.fftpack).

Fast Fourier Transform (FFT)

The Fast Fourier Transform (FFT) is an efficient algorithm to compute the Discrete Fourier Transform (DFT) and its inverse.

One-Dimensional Discrete Fourier Transform

The DFT transforms a sequence of N complex numbers $x_0, x_1, \ldots, x_{N-1}$ into another sequence of N complex numbers $X_0, X_1, \ldots, X_{N-1}$, often called frequency bins. The fft() function in scipy.fft computes the DFT, and ifft() computes the inverse DFT.

import numpy as np
from scipy.fft import fft, ifft

# Create a simple input array (time-domain signal)
x = np.array([1.0, 2.0, 1.0, -1.0, 1.5])
print(f"Original signal x: {x}")

# Apply the FFT
y = fft(x)
print(f"FFT of x (frequency-domain signal y):\n{y}")

# Apply the Inverse FFT (IFFT)
y_inverted = ifft(y)
print(f"IFFT of y (should be close to original x):\n{y_inverted}")
# Note: IFFT result might have very small imaginary parts due to floating point precision.
# We can take the real part if the original signal was real.
print(f"Real part of IFFT of y: {np.real(y_inverted)}")

The output will be:

Original signal x: [ 1.   2.   1.  -1.   1.5]
FFT of x (frequency-domain signal y):
[ 4.50000000+0.j          2.08155948-1.65109876j -1.83155948+1.60822041j
 -1.83155948-1.60822041j  2.08155948+1.65109876j]
IFFT of y (should be close to original x):
[ 1.0+0.j   2.0+0.j   1.0+0.j  -1.0+0.j   1.5+0.j ]
Real part of IFFT of y: [ 1.   2.   1.  -1.   1.5]

Let’s consider a more practical example with a noisy sine wave.

import numpy as np
from scipy.fft import fft, fftfreq
import matplotlib.pyplot as plt # For plotting (optional)

# Signal parameters
sampling_rate = 100  # Hz (samples per second)
duration = 5.0       # seconds
frequency = 2.0      # Hz (frequency of the sine wave)

num_samples = int(sampling_rate * duration)
time_vector = np.linspace(0.0, duration, num_samples, endpoint=False)

# Generate a sine wave with some noise
signal = np.sin(2 * np.pi * frequency * time_vector) + 0.5 * np.random.randn(num_samples)

print(f"Number of samples: {signal.size}")

# Compute the FFT
signal_fft = fft(signal)

# Generate the frequency bins for the x-axis of the FFT plot
# fftfreq(number_of_samples, sample_spacing)
sample_spacing = 1.0 / sampling_rate
frequency_bins = fftfreq(num_samples, d=sample_spacing)

# Plotting (optional, requires matplotlib)
# plt.figure(figsize=(12, 6))
# plt.subplot(1, 2, 1)
# plt.plot(time_vector, signal)
# plt.title("Time-Domain Signal")
# plt.xlabel("Time (s)")
# plt.ylabel("Amplitude")

# plt.subplot(1, 2, 2)
# # We plot the magnitude of the FFT (abs value) and usually only the positive frequencies
# positive_freq_indices = frequency_bins >= 0
# plt.plot(frequency_bins[positive_freq_indices], np.abs(signal_fft[positive_freq_indices]))
# plt.title("Frequency-Domain Signal (FFT Magnitude)")
# plt.xlabel("Frequency (Hz)")
# plt.ylabel("Magnitude")
# plt.xlim(0, sampling_rate / 2) # Show up to Nyquist frequency
# plt.grid(True)
# plt.tight_layout()
# plt.show()

# For text output, let's show some FFT values and corresponding frequencies
print("\nSample FFT magnitudes and their frequencies:")
for i in range(5):
    print(f"Freq: {frequency_bins[i]:.2f} Hz, Magnitude: {np.abs(signal_fft[i]):.2f}")
# Find the dominant frequency component
dominant_freq_idx = np.argmax(np.abs(signal_fft[frequency_bins >= 0])) # consider positive freqs
print(f"Dominant positive frequency: {frequency_bins[dominant_freq_idx]:.2f} Hz")

This code generates a noisy sine wave. The fftfreq function calculates the frequencies corresponding to each point in the FFT output. The FFT output signal_fft contains complex numbers; np.abs(signal_fft) gives their magnitudes, which represent the strength of each frequency component. For real input signals, the FFT output is symmetric, so we often only analyze the positive frequencies.

A sample text output (exact values will vary due to noise):

Number of samples: 500

Sample FFT magnitudes and their frequencies:
Freq: 0.00 Hz, Magnitude: 1.89
Freq: 0.20 Hz, Magnitude: 6.09
Freq: 0.40 Hz, Magnitude: 6.37
Freq: 0.60 Hz, Magnitude: 3.90
Freq: 0.80 Hz, Magnitude: 10.47
Dominant positive frequency: 2.00 Hz

The dominant frequency found should be close to the original frequency of 2.0 Hz we set for the sine wave.

Discrete Cosine Transform (DCT)

The Discrete Cosine Transform (DCT) is related to the DFT but uses only cosine functions. It expresses a finite sequence of data points in terms of a sum of cosine functions oscillating at different frequencies. DCTs are widely used in signal and image compression (e.g., JPEG).

SciPy’s scipy.fft module provides dct for the forward DCT and idct for the inverse DCT (IDCT).

from scipy.fft import dct, idct
import numpy as np

# Input array
x_dct = np.array([4., 3., 5., 10., 5., 3.])
print(f"Original array for DCT: {x_dct}")

# Compute DCT (default is Type-II DCT)
dct_coeffs = dct(x_dct)
print(f"DCT coefficients:\n{dct_coeffs}")

# Reconstruct the original signal using IDCT (default is Type-III IDCT for Type-II DCT)
x_reconstructed = idct(dct_coeffs)
print(f"Reconstructed array using IDCT:\n{x_reconstructed}")

The output will be:

Original array for DCT: [ 4.  3.  5. 10.  5.  3.]
DCT coefficients:
[ 60.          -3.48476592 -13.85640646  11.3137085   6.
  -6.31319305]
Reconstructed array using IDCT:
[ 4.  3.  5. 10.  5.  3.]

There are different types of DCTs (Type-I, Type-II, Type-III, Type-IV). scipy.fft.dct and scipy.fft.idct allow specifying the type using the type parameter. The default dct is Type-II, and its corresponding inverse idct is Type-III. These are common choices for applications like compression.

The scipy.fft module also provides functions for 2-D (and N-D) FFTs and DCTs (e.g., fft2, dctn), which are essential for image processing and other multidimensional data analysis.