How to Use FFT in NumPy
The Fast Fourier Transform is one of the most important algorithms in signal processing. It takes a signal that varies over time and decomposes it into its constituent frequencies. Think of it as...
Key Insights
- NumPy’s FFT transforms time-domain signals into frequency-domain representations, revealing the underlying frequency components that compose any signal
- The
fftfreq()function is essential for interpreting FFT output—without it, you’re just looking at meaningless array indices - For real-valued signals (which covers most practical applications), use
rfft()instead offft()for better performance and cleaner output
Introduction to FFT
The Fast Fourier Transform is one of the most important algorithms in signal processing. It takes a signal that varies over time and decomposes it into its constituent frequencies. Think of it as reverse-engineering a chord to find the individual notes being played.
NumPy provides a complete FFT implementation in its numpy.fft module. It’s not the fastest option available (libraries like FFTW or Intel MKL are faster), but it’s good enough for most applications and requires zero additional dependencies. If you’re doing exploratory data analysis, prototyping, or working with moderately-sized datasets, NumPy’s FFT is the pragmatic choice.
Understanding the Basics: Time vs. Frequency Domain
Every signal can be represented in two ways. The time domain shows amplitude over time—what you see on an oscilloscope. The frequency domain shows amplitude at each frequency—what you see on a spectrum analyzer.
The FFT converts between these representations. The forward transform goes from time to frequency; the inverse transform goes back. Both are lossless—you can round-trip without losing information.
Here’s how to generate and visualize a simple sine wave:
import numpy as np
import matplotlib.pyplot as plt
# Sampling parameters
sample_rate = 1000 # Hz
duration = 1.0 # seconds
frequency = 5 # Hz
# Generate time array and sine wave
t = np.linspace(0, duration, int(sample_rate * duration), endpoint=False)
signal = np.sin(2 * np.pi * frequency * t)
# Plot the time-domain signal
plt.figure(figsize=(10, 4))
plt.plot(t[:200], signal[:200]) # First 200ms
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')
plt.title('5 Hz Sine Wave')
plt.grid(True)
plt.show()
This creates a clean 5 Hz sine wave sampled at 1000 Hz. The signal completes 5 full cycles per second, and we’re capturing 1000 samples per second—more than enough to represent it accurately.
Core NumPy FFT Functions
NumPy’s FFT module has four functions you’ll use constantly:
numpy.fft.fft(a) - Computes the discrete Fourier transform. Returns a complex array of the same length as the input. The magnitude of each element represents amplitude at that frequency; the angle represents phase.
numpy.fft.ifft(a) - The inverse transform. Takes frequency-domain data and reconstructs the time-domain signal.
numpy.fft.fftfreq(n, d) - Returns the frequencies corresponding to each FFT bin. The parameter n is the signal length, and d is the sample spacing (1/sample_rate).
numpy.fft.fftshift(a) - Shifts the zero-frequency component to the center of the array. Makes visualization more intuitive.
Let’s apply these to our sine wave:
# Compute FFT
fft_result = np.fft.fft(signal)
frequencies = np.fft.fftfreq(len(signal), 1/sample_rate)
# Get magnitude (absolute value of complex numbers)
magnitude = np.abs(fft_result)
# Plot frequency spectrum (positive frequencies only)
positive_mask = frequencies >= 0
plt.figure(figsize=(10, 4))
plt.plot(frequencies[positive_mask], magnitude[positive_mask])
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude')
plt.title('Frequency Spectrum')
plt.xlim(0, 50) # Focus on 0-50 Hz range
plt.grid(True)
plt.show()
# Find the peak frequency
peak_idx = np.argmax(magnitude[:len(magnitude)//2])
print(f"Dominant frequency: {frequencies[peak_idx]:.1f} Hz")
The output shows a sharp spike at 5 Hz—exactly where our sine wave’s frequency is. The magnitude at that frequency is proportional to the amplitude of the original signal.
Working with Real-World Signals
Real signals rarely contain just one frequency. Let’s create a composite signal and use FFT to decompose it:
# Create a signal with three frequency components
freq1, freq2, freq3 = 10, 25, 60 # Hz
amp1, amp2, amp3 = 1.0, 0.5, 0.3 # Amplitudes
signal = (amp1 * np.sin(2 * np.pi * freq1 * t) +
amp2 * np.sin(2 * np.pi * freq2 * t) +
amp3 * np.sin(2 * np.pi * freq3 * t))
# Compute FFT
fft_result = np.fft.fft(signal)
frequencies = np.fft.fftfreq(len(signal), 1/sample_rate)
magnitude = np.abs(fft_result) / len(signal) # Normalize
# Plot
fig, axes = plt.subplots(2, 1, figsize=(10, 6))
# Time domain
axes[0].plot(t[:500], signal[:500])
axes[0].set_xlabel('Time (s)')
axes[0].set_ylabel('Amplitude')
axes[0].set_title('Composite Signal (Time Domain)')
# Frequency domain
positive_mask = frequencies >= 0
axes[1].plot(frequencies[positive_mask], magnitude[positive_mask] * 2)
axes[1].set_xlabel('Frequency (Hz)')
axes[1].set_ylabel('Amplitude')
axes[1].set_title('Frequency Spectrum')
axes[1].set_xlim(0, 100)
plt.tight_layout()
plt.show()
# Identify peaks
threshold = 0.1
peaks = frequencies[positive_mask][magnitude[positive_mask] * 2 > threshold]
print(f"Detected frequencies: {peaks}")
The FFT cleanly separates the three frequencies. The relative heights of the peaks correspond to the original amplitudes. This is the power of frequency-domain analysis—patterns that are invisible in time-domain data become obvious.
Practical Considerations
Sampling Rate and Nyquist Theorem
You can only detect frequencies up to half your sampling rate (the Nyquist frequency). A 1000 Hz sample rate means you can detect frequencies up to 500 Hz. Anything higher will alias—appearing as a lower frequency in your results.
Windowing
FFT assumes your signal is periodic. If it’s not (and real signals never are), you get spectral leakage—energy smearing across frequencies. Windowing functions taper the signal edges to reduce this:
# Generate a signal that doesn't complete full cycles
freq = 10.5 # Non-integer frequency causes leakage
signal = np.sin(2 * np.pi * freq * t)
# Apply Hanning window
window = np.hanning(len(signal))
signal_windowed = signal * window
# Compare FFT results
fft_no_window = np.abs(np.fft.fft(signal))
fft_with_window = np.abs(np.fft.fft(signal_windowed))
frequencies = np.fft.fftfreq(len(signal), 1/sample_rate)
fig, axes = plt.subplots(1, 2, figsize=(12, 4))
positive_mask = frequencies >= 0
axes[0].plot(frequencies[positive_mask], fft_no_window[positive_mask])
axes[0].set_xlim(0, 30)
axes[0].set_title('Without Windowing (Spectral Leakage)')
axes[0].set_xlabel('Frequency (Hz)')
axes[1].plot(frequencies[positive_mask], fft_with_window[positive_mask])
axes[1].set_xlim(0, 30)
axes[1].set_title('With Hanning Window (Clean Peak)')
axes[1].set_xlabel('Frequency (Hz)')
plt.tight_layout()
plt.show()
The windowed version shows a much cleaner peak with less energy bleeding into adjacent frequencies.
Use rfft() for Real Signals
If your input is real-valued (no imaginary component), use numpy.fft.rfft() instead. It’s faster and returns only the positive frequencies, since the negative frequencies are just mirror images:
# rfft is more efficient for real-valued signals
rfft_result = np.fft.rfft(signal)
rfft_frequencies = np.fft.rfftfreq(len(signal), 1/sample_rate)
print(f"fft output length: {len(np.fft.fft(signal))}")
print(f"rfft output length: {len(rfft_result)}")
Common Applications and Use Cases
FFT appears everywhere in signal processing. Audio analysis uses it for spectrum visualization and pitch detection. Image processing uses the 2D variant (numpy.fft.fft2) for filtering and compression. Communications systems use it for modulation and demodulation.
Here’s a practical example—a simple low-pass filter that removes high-frequency noise:
# Create a noisy signal
clean_signal = np.sin(2 * np.pi * 5 * t)
noise = 0.5 * np.random.randn(len(t))
noisy_signal = clean_signal + noise
# FFT
fft_result = np.fft.fft(noisy_signal)
frequencies = np.fft.fftfreq(len(noisy_signal), 1/sample_rate)
# Zero out frequencies above cutoff
cutoff_freq = 20 # Hz
fft_filtered = fft_result.copy()
fft_filtered[np.abs(frequencies) > cutoff_freq] = 0
# Inverse FFT to get filtered signal
filtered_signal = np.fft.ifft(fft_filtered).real
# Plot comparison
fig, axes = plt.subplots(3, 1, figsize=(10, 8))
axes[0].plot(t[:500], noisy_signal[:500])
axes[0].set_title('Noisy Signal')
axes[1].plot(t[:500], filtered_signal[:500])
axes[1].set_title('Filtered Signal (Low-Pass at 20 Hz)')
axes[2].plot(t[:500], clean_signal[:500])
axes[2].set_title('Original Clean Signal')
for ax in axes:
ax.set_xlabel('Time (s)')
ax.set_ylabel('Amplitude')
plt.tight_layout()
plt.show()
This filter works by transforming to frequency domain, zeroing out unwanted frequencies, and transforming back. It’s not the most sophisticated approach (a proper filter would use a smoother transition), but it demonstrates the principle.
Summary and Further Resources
NumPy’s FFT module gives you everything needed for basic frequency analysis. The core workflow is straightforward: generate or load your signal, apply fft() or rfft(), use fftfreq() to interpret the results, and optionally use ifft() to transform back.
Key functions to remember:
fft()/ifft()- Forward and inverse transformsrfft()/irfft()- Optimized versions for real-valued signalsfftfreq()/rfftfreq()- Get frequency values for each binfftshift()- Center the zero-frequency component
For more advanced signal processing—filter design, spectrograms, wavelets—look at SciPy’s scipy.signal and scipy.fft modules. They build on NumPy’s foundation with additional algorithms and convenience functions.
The NumPy FFT documentation covers edge cases and additional functions not discussed here. For theoretical background, any digital signal processing textbook will explain the mathematics behind these transforms.