Python SciPy for Signals

Signal processing with Python SciPy — filter design, spectral analysis, convolution, and system simulation using scipy.signal.

Introduction

Python with NumPy and SciPy has become the dominant open-source platform for signal processing. The scipy.signal module provides comprehensive tools for filter design, spectral analysis, system simulation, and waveform generation. Combined with Matplotlib for visualization and NumPy for numerical computation, Python offers a complete signal processing environment that is free, widely supported, and increasingly used in both academia and industry.

Unlike MATLAB (proprietary, expensive), Python is freely available and integrates with the broader scientific computing ecosystem. Skills developed in Python for signal processing transfer directly to machine learning, data science, and embedded systems programming. The syntax closely mirrors MATLAB for most signal processing operations, making transitions between the two platforms straightforward.

This guide provides working code for every major signal processing task you will encounter in a signals and systems course, organized by workflow stage from signal creation through analysis and filtering.

Essential Imports

python example

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal
from scipy.fft import fft, fftfreq, ifft
from scipy.io import wavfile  # For audio file I/O

Generating Signals

Example

fs = 1000 # Sampling frequency t = np.arange(0, 1, 1/fs) # Time vector: 0 to 1 second N = len(t) # Number of samples # Sinusoid x_sin = np.sin(2*np.pi*50*t) # 50 Hz sine # Exponential decay x_exp = np.exp(-5*t) # Unit step x_step = np.heaviside(t - 0.2, 0.5) # Step at t=0.2 # Rectangular pulse x_rect = signal.square(2*np.pi*5*t, duty=0.3) # 5 Hz, 30% duty # Chirp (frequency sweep) x_chirp = signal.chirp(t, f0=10, t1=1, f1=200, method='linear') # Multi-tone signal x_multi = (np.sin(2*np.pi*100*t) + 0.5*np.sin(2*np.pi*250*t) + 0.3*np.sin(2*np.pi*400*t)) # Signal with additive noise x_noisy = x_sin + 0.5*np.random.randn(N) # Impulse signal x_impulse = np.zeros(N) x_impulse[0] = 1

Important conventions: NumPy arrays are zero-indexed. Use np.pi (not math.pi) for vectorized operations. The signal.chirp function handles linear, quadratic, and logarithmic frequency sweeps.

Convolution

Example

# Discrete convolution x = np.array([1, 2, 3, 2, 1]) h = np.array([1, -1, 0.5]) y = np.convolve(x, h) print(f"Output: {y}") # Length = 5 + 3 - 1 = 7 # Mode options control output length y_full = np.convolve(x, h, mode='full') # Length Lx + Lh - 1 y_same = np.convolve(x, h, mode='same') # Length max(Lx, Lh) y_valid = np.convolve(x, h, mode='valid') # Length |Lx - Lh| + 1 # Continuous-time approximation via discrete convolution dt = 0.001 t_ct = np.arange(0, 2, dt) x_ct = np.exp(-2*t_ct) h_ct = np.exp(-3*t_ct) y_ct = np.convolve(x_ct, h_ct) * dt # Scale by dt for integral approx # FFT-based convolution (faster for long signals) y_fft = signal.fftconvolve(x_ct, h_ct, mode='full') * dt

Performance note: np.convolve uses direct computation ($O(N \cdot M)$). For long signals, signal.fftconvolve uses FFT-based multiplication ($O(N \log N)$) and is significantly faster when $N > 500$.

Filter Design and Application

Example

# Butterworth low-pass filter order = 4 fc = 100 # Cutoff frequency in Hz b, a = signal.butter(order, fc, btype='low', fs=fs) # Apply filter (two methods) y_filtered = signal.filtfilt(b, a, x_noisy) # Zero-phase (non-causal) y_causal = signal.lfilter(b, a, x_noisy) # Causal filtering # FIR filter design (windowed sinc) numtaps = 51 fir_coeff = signal.firwin(numtaps, 100, fs=fs) y_fir = signal.lfilter(fir_coeff, 1, x_noisy) # Bandpass filter b_bp, a_bp = signal.butter(4, [80, 120], btype='band', fs=fs) y_bp = signal.filtfilt(b_bp, a_bp, x_noisy) # Notch filter (remove specific frequency) b_notch, a_notch = signal.iirnotch(60, Q=30, fs=fs) # Remove 60 Hz y_notch = signal.filtfilt(b_notch, a_notch, x_noisy) # Chebyshev Type I filter b_ch, a_ch = signal.cheby1(4, 0.5, 100, btype='low', fs=fs) # Elliptic filter (sharpest transition) b_el, a_el = signal.ellip(4, 0.5, 40, 100, btype='low', fs=fs)

When to use which filter:

filtfilt: stored data analysis (zero phase distortion, doubles the order effectively)
lfilter: real-time or causal simulation (introduces group delay)
sosfilt/sosfiltfilt: numerically stable for high-order filters (use output='sos' in design)

Second-Order Sections (Recommended for High Orders)

python example

# For orders > 6, use SOS format to avoid numerical instability
sos = signal.butter(10, 100, btype='low', fs=fs, output='sos')
y_sos = signal.sosfiltfilt(sos, x_noisy)

Spectral Analysis

Example

# FFT-based spectrum X = fft(x_sin) freqs = fftfreq(N, 1/fs) plt.figure() plt.plot(freqs[:N//2], 2*np.abs(X[:N//2])/N) plt.xlabel('Frequency (Hz)') plt.ylabel('Magnitude') plt.title('Single-Sided Amplitude Spectrum') # Power Spectral Density (Welch method) f_psd, Pxx = signal.welch(x_noisy, fs=fs, nperseg=256) plt.figure() plt.semilogy(f_psd, Pxx) plt.xlabel('Frequency (Hz)') plt.ylabel('PSD (V²/Hz)') plt.title('Power Spectral Density') # Spectrogram (time-frequency analysis) f_spec, t_spec, Sxx = signal.spectrogram(x_chirp, fs=fs, nperseg=128) plt.figure() plt.pcolormesh(t_spec, f_spec, 10*np.log10(Sxx), shading='gouraud') plt.ylabel('Frequency (Hz)') plt.xlabel('Time (s)') plt.colorbar(label='Power (dB)') plt.title('Spectrogram') # Windowed FFT with explicit window window = signal.windows.hann(N) X_win = fft(x_sin * window) plt.plot(freqs[:N//2], 2*np.abs(X_win[:N//2])/np.sum(window))

System Analysis

Example

# Define transfer function: H(s) = 1/(s^2 + 2s + 5) num = [1] den = [1, 2, 5] sys = signal.TransferFunction(num, den) # Impulse response t_imp, y_imp = signal.impulse(sys) plt.figure() plt.plot(t_imp, y_imp) plt.title('Impulse Response') plt.xlabel('Time (s)') plt.grid(True) # Step response t_step, y_step = signal.step(sys) # Frequency response (Bode plot) w, mag, phase = signal.bode(sys) fig, (ax1, ax2) = plt.subplots(2, 1, sharex=True) ax1.semilogx(w, mag); ax1.set_ylabel('Magnitude (dB)') ax2.semilogx(w, phase); ax2.set_ylabel('Phase (degrees)') ax2.set_xlabel('Frequency (rad/s)') # Pole-zero analysis zeros, poles, gain = signal.tf2zpk(num, den) print(f"Poles: {poles}") print(f"Zeros: {zeros}") print(f"System is {'stable' if all(np.real(poles) < 0) else 'unstable'}")

Z-Transform and Digital Systems

Example

# Discrete-time system: H(z) = (1 + 0.5z^-1)/(1 - 0.8z^-1) b_d = [1, 0.5] # Numerator coefficients [b0, b1, ...] a_d = [1, -0.8] # Denominator coefficients [a0, a1, ...] # Frequency response w, H = signal.freqz(b_d, a_d, worN=512, fs=fs) plt.figure() plt.subplot(2, 1, 1) plt.plot(w, 20*np.log10(np.abs(H))) plt.ylabel('Magnitude (dB)') plt.subplot(2, 1, 2) plt.plot(w, np.degrees(np.angle(H))) plt.ylabel('Phase (degrees)') plt.xlabel('Frequency (Hz)') # Pole-zero analysis z, p, k = signal.tf2zpk(b_d, a_d) print(f"Zeros: {z}") print(f"Poles: {p}") print(f"Stable: {all(np.abs(p) < 1)}") # Impulse response of discrete system imp = np.zeros(50) imp[0] = 1 h_discrete = signal.lfilter(b_d, a_d, imp) plt.stem(h_discrete)

Practical Tips

Use filtfilt for zero-phase filtering when processing stored data (not real-time)
Use lfilter for causal filtering in real-time or simulation scenarios
Always check filter stability with np.all(np.abs(np.roots(a)) < 1)
Use SOS format for high-order IIR filters to avoid numerical issues: pass output='sos'
Use signal.welch instead of raw FFT for PSD estimation of noisy signals
Normalize frequency correctly: pass fs= parameter to avoid manual Nyquist division
Use signal.fftconvolve for long signals — it is $O(N\log N)$ vs $O(N^2)$
Window before FFT to reduce spectral leakage (Hanning window is a good default)

Key Takeaways

scipy.signal provides comprehensive signal processing: filters, spectra, systems
NumPy handles array operations; Matplotlib handles visualization; SciPy provides DSP algorithms
np.convolve for discrete convolution; signal.fftconvolve for efficient long convolutions
Filter design: butter, cheby1, ellip, firwin — all available in scipy.signal
System analysis: impulse, step, bode, freqz — matching MATLAB functionality
Use SOS format (output='sos') for numerically stable high-order IIR filtering
Python is free, widely adopted, and skills transfer directly to ML and data science

Introduction

Essential Imports

python example

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal
from scipy.fft import fft, fftfreq, ifft
from scipy.io import wavfile  # For audio file I/O

Generating Signals

Example

Convolution

Example

Filter Design and Application

Example

When to use which filter:

filtfilt: stored data analysis (zero phase distortion, doubles the order effectively)
lfilter: real-time or causal simulation (introduces group delay)
sosfilt/sosfiltfilt: numerically stable for high-order filters (use output='sos' in design)

Second-Order Sections (Recommended for High Orders)

python example

# For orders > 6, use SOS format to avoid numerical instability
sos = signal.butter(10, 100, btype='low', fs=fs, output='sos')
y_sos = signal.sosfiltfilt(sos, x_noisy)

Spectral Analysis

Example

System Analysis

Example

Z-Transform and Digital Systems

Example

Practical Tips

Use filtfilt for zero-phase filtering when processing stored data (not real-time)
Use lfilter for causal filtering in real-time or simulation scenarios
Always check filter stability with np.all(np.abs(np.roots(a)) < 1)
Use SOS format for high-order IIR filters to avoid numerical issues: pass output='sos'
Use signal.welch instead of raw FFT for PSD estimation of noisy signals
Normalize frequency correctly: pass fs= parameter to avoid manual Nyquist division
Use signal.fftconvolve for long signals — it is $O(N\log N)$ vs $O(N^2)$
Window before FFT to reduce spectral leakage (Hanning window is a good default)

Key Takeaways

scipy.signal provides comprehensive signal processing: filters, spectra, systems
NumPy handles array operations; Matplotlib handles visualization; SciPy provides DSP algorithms
np.convolve for discrete convolution; signal.fftconvolve for efficient long convolutions
Filter design: butter, cheby1, ellip, firwin — all available in scipy.signal
System analysis: impulse, step, bode, freqz — matching MATLAB functionality
Use SOS format (output='sos') for numerically stable high-order IIR filtering
Python is free, widely adopted, and skills transfer directly to ML and data science

Python SciPy for Signals

Introduction

Essential Imports

Generating Signals

Convolution

Filter Design and Application

Second-Order Sections (Recommended for High Orders)

Spectral Analysis

System Analysis

Z-Transform and Digital Systems

Practical Tips

Key Takeaways

Continue learning this concept

Python SciPy for Signals

Introduction

Essential Imports

Generating Signals

Convolution

Filter Design and Application

Second-Order Sections (Recommended for High Orders)

Spectral Analysis

System Analysis

Z-Transform and Digital Systems

Practical Tips

Key Takeaways

Continue learning this concept