Python Scipy Freqz: Analyze Signal Frequency Response

When working on signal processing projects, I often need to analyze how different systems respond to various frequencies. This is where the freqz function from SciPy’s signal module becomes invaluable. It helps me visualize and understand frequency responses of digital filters without diving into complex mathematical calculations.

In this article, I’ll share several practical ways to use scipy.signal.freqz for analyzing digital filters and signals. I’ve used these techniques in my projects, and they’ve saved me countless hours of manual calculations.

Let’s get into how you can leverage this powerful function in your own Python projects.

SciPy’s freqz Function

The freqz function computes the frequency response of a digital filter. Simply put, it shows how a filter affects different frequency components of a signal.

This function is handy when:

Designing audio filters for music applications
Working on communication systems
Analyzing speech or biomedical signals
Developing control systems

Let’s start with a basic example to understand how it works.

Read SciPy Signal

Method 1: Basic Usage of freqz with Filter Coefficients

The easiest way to use freqz is with filter coefficients. This method is perfect when you already have your filter defined by its numerator and denominator coefficients.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Define a simple low-pass filter
b = [0.25, 0.25, 0.25, 0.25]  # Numerator coefficients
a = [1]  # Denominator coefficients

# Calculate frequency response
w, h = signal.freqz(b, a)

# Convert to frequency in Hz (assuming 1000 Hz sample rate)
f = w * 1000 / (2 * np.pi)

# Plot the frequency response
plt.figure(figsize=(10, 6))
plt.plot(f, 20 * np.log10(abs(h)))
plt.title('Filter Frequency Response')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Amplitude (dB)')
plt.grid(True)
plt.show()

I executed the above example code and added the screenshot below

This example creates a simple moving average filter and visualizes its frequency response. I find this particularly useful when I need a quick assessment of how my filter will perform.

Check out SciPy Integrate

Method 2: Use freqz with Second-Order Sections (SOS)

When working with higher-order filters, I often use the SOS format to maintain numerical stability. The freqz_sos function is perfect for this scenario.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Create a 6th-order Butterworth bandpass filter in SOS format
sos = signal.butter(6, [0.1, 0.4], 'bandpass', output='sos')

# Calculate frequency response
w, h = signal.sosfreqz(sos, worN=2000)

# Plot magnitude response
plt.figure(figsize=(10, 6))
plt.plot(w / np.pi, 20 * np.log10(abs(h)))
plt.title('Bandpass Filter Frequency Response')
plt.xlabel('Normalized Frequency (×π rad/sample)')
plt.ylabel('Magnitude (dB)')
plt.grid(True)
plt.axvline(0.1, color='r', linestyle='--', alpha=0.3)
plt.axvline(0.4, color='r', linestyle='--', alpha=0.3)
plt.show()

I executed the above example code and added the screenshot below

This example creates a bandpass filter that could be used in applications like isolating speech frequencies in a noisy environment.

Read SciPy Misc

Method 3: Customize the Frequency Points

Sometimes I need to analyze specific frequency ranges in detail. The worN parameter lets me specify exactly which frequency points to evaluate.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Define a notch filter (to remove 60 Hz noise from a US power line)
fs = 1000  # Sampling frequency of 1 kHz
f0 = 60    # Notch frequency (60 Hz power line interference)
Q = 30     # Quality factor

# Design the filter
b, a = signal.iirnotch(f0, Q, fs)

# Create custom frequency points with more detail around 60 Hz
freq_points = np.concatenate([
    np.linspace(0, 50, 100),
    np.linspace(50, 70, 200),  # More points around 60 Hz
    np.linspace(70, 500, 100)
])

# Calculate the frequency response at custom points
w, h = signal.freqz(b, a, worN=2*np.pi*freq_points/fs)

# Convert to Hz for plotting
f = w * fs / (2 * np.pi)

# Plot magnitude response
plt.figure(figsize=(10, 6))
plt.plot(f, 20 * np.log10(abs(h)))
plt.title('Notch Filter for 60 Hz Power Line Interference')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude (dB)')
plt.grid(True)
plt.axvline(60, color='r', linestyle='--', label='Notch Frequency')
plt.legend()
plt.show()

I executed the above example code and added the screenshot below

This example creates a notch filter to remove 60 Hz power line interference, a common issue in audio recording and biomedical signal processing in the US.

Check out SciPy Stats

Method 4: Analyze Phase Response

When working on applications like audio processing or communication systems, I often need to analyze both magnitude and phase responses.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Design an 8th-order Chebyshev Type I filter
b, a = signal.cheby1(8, 1, 0.2)

# Calculate frequency response
w, h = signal.freqz(b, a)

# Plot magnitude and phase
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(10, 10))

# Magnitude response
ax1.plot(w / np.pi, 20 * np.log10(abs(h)))
ax1.set_title('Chebyshev Type I Filter Frequency Response')
ax1.set_xlabel('Normalized Frequency (×π rad/sample)')
ax1.set_ylabel('Magnitude (dB)')
ax1.grid(True)
ax1.set_ylim([-80, 5])

# Phase response
angles = np.unwrap(np.angle(h))
ax2.plot(w / np.pi, angles)
ax2.set_xlabel('Normalized Frequency (×π rad/sample)')
ax2.set_ylabel('Phase (radians)')
ax2.grid(True)
ax2.set_yticks(np.arange(-10, 11, 2) * np.pi/2)
ax2.set_yticklabels(['$-5\pi$', '$-4\pi$', '$-3\pi$', '$-2\pi$', '$-\pi$', 
                     '0', '$\pi$', '$2\pi$', '$3\pi$', '$4\pi$', '$5\pi$'])

plt.tight_layout()
plt.show()

This example shows both magnitude and phase responses of a Chebyshev filter, which is useful when phase distortion is a concern, such as in high-fidelity audio applications.

Read SciPy Convolve

Method 5: Compare Multiple Filters

When designing a system, I often need to compare different filter options to choose the best one for my needs.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Design three different filters
# Sample rate
fs = 44100  # Standard audio sample rate (Hz)
cutoff = 5000  # Cutoff frequency (Hz)
normalized_cutoff = cutoff / (fs / 2)  # Normalized cutoff frequency

# Design filters
b_butter, a_butter = signal.butter(4, normalized_cutoff, 'low')
b_cheby, a_cheby = signal.cheby1(4, 1, normalized_cutoff, 'low')
b_bessel, a_bessel = signal.bessel(4, normalized_cutoff, 'low')

# Calculate frequency responses
w, h_butter = signal.freqz(b_butter, a_butter)
_, h_cheby = signal.freqz(b_cheby, a_cheby)
_, h_bessel = signal.freqz(b_bessel, a_bessel)

# Convert to frequency in Hz for plotting
f = w * fs / (2 * np.pi)

# Plot all filters
plt.figure(figsize=(12, 7))
plt.plot(f, 20 * np.log10(abs(h_butter)), label='Butterworth')
plt.plot(f, 20 * np.log10(abs(h_cheby)), label='Chebyshev Type I')
plt.plot(f, 20 * np.log10(abs(h_bessel)), label='Bessel')
plt.axvline(cutoff, color='r', linestyle='--', alpha=0.3, label='Cutoff (5 kHz)')
plt.title('Comparison of 4th-Order Low-Pass Filters')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude (dB)')
plt.grid(True)
plt.legend()
plt.xlim([0, 15000])
plt.ylim([-80, 5])
plt.show()

This example compares three common filter types (Butterworth, Chebyshev, and Bessel) that could be used in an audio application to limit high-frequency content.

Read SciPy Ndimage Rotate

Method 6: Analyze Zeros, Poles, and Gain (ZPK) Format

For complex filter analysis, I sometimes prefer working with zeros, poles, and gain representation rather than transfer function coefficients.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Define a filter using zeros, poles, and gain
zeros = [np.exp(1j * np.pi / 4), np.exp(-1j * np.pi / 4)]
poles = [0.8 * np.exp(1j * np.pi / 8), 0.8 * np.exp(-1j * np.pi / 8)]
gain = 1.0

# Calculate frequency response
w, h = signal.freqz_zpk(zeros, poles, gain)

# Plot magnitude response
plt.figure(figsize=(10, 6))
plt.plot(w / np.pi, 20 * np.log10(abs(h)))
plt.title('ZPK Filter Frequency Response')
plt.xlabel('Normalized Frequency (×π rad/sample)')
plt.ylabel('Magnitude (dB)')
plt.grid(True)
plt.show()

# Plot zeros and poles
plt.figure(figsize=(8, 8))
unit_circle = plt.Circle((0, 0), 1, fill=False, color='black', linestyle='--', alpha=0.3)
plt.gca().add_patch(unit_circle)
plt.scatter(np.real(zeros), np.imag(zeros), marker='o', color='blue', label='Zeros')
plt.scatter(np.real(poles), np.imag(poles), marker='x', color='red', label='Poles')
plt.axis('equal')
plt.grid(True)
plt.title('Pole-Zero Plot')
plt.xlabel('Real Part')
plt.ylabel('Imaginary Part')
plt.legend()
plt.show()

This example is particularly useful for understanding filter stability and behavior from a theoretical perspective.

Read SciPy Stats Z-score

Method 7: Practical Application – Create an Equalizer

Let’s build something practical – a simple three-band equalizer that could be used in a music application.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Sample rate
fs = 44100  # Hz

# Define bands
low_cutoff = 500  # Hz
high_cutoff = 4000  # Hz

# Normalize frequencies
low_norm = low_cutoff / (fs / 2)
high_norm = high_cutoff / (fs / 2)

# Create filters for each band (4th order)
b_low, a_low = signal.butter(4, low_norm, 'lowpass')
b_mid, a_mid = signal.butter(4, [low_norm, high_norm], 'bandpass')
b_high, a_high = signal.butter(4, high_norm, 'highpass')

# Calculate frequency responses
worN = 4096
w, h_low = signal.freqz(b_low, a_low, worN=worN)
_, h_mid = signal.freqz(b_mid, a_mid, worN=worN)
_, h_high = signal.freqz(b_high, a_high, worN=worN)

# Convert to Hz for plotting
f = w * fs / (2 * np.pi)

# Apply gain to each band (in dB)
low_gain = 3  # Boost bass by 3 dB
mid_gain = 0  # Leave midrange flat
high_gain = 2  # Boost treble by 2 dB

# Convert dB to linear scale for multiplication
h_low_boosted = h_low * (10 ** (low_gain / 20))
h_mid_boosted = h_mid * (10 ** (mid_gain / 20))
h_high_boosted = h_high * (10 ** (high_gain / 20))

# Combined response (sum of all bands)
h_combined = h_low_boosted + h_mid_boosted + h_high_boosted

# Plot
plt.figure(figsize=(12, 8))
plt.semilogx(f, 20 * np.log10(abs(h_low_boosted)), label='Low Band (+3 dB)')
plt.semilogx(f, 20 * np.log10(abs(h_mid_boosted)), label='Mid Band (0 dB)')
plt.semilogx(f, 20 * np.log10(abs(h_high_boosted)), label='High Band (+2 dB)')
plt.semilogx(f, 20 * np.log10(abs(h_combined)), 'k--', label='Combined Response')

plt.title('Three-Band Equalizer Frequency Response')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude (dB)')
plt.grid(True, which="both", ls="--")
plt.legend()
plt.xlim([20, 20000])  # Standard audio frequency range
plt.axvline(low_cutoff, color='r', linestyle='--', alpha=0.3)
plt.axvline(high_cutoff, color='r', linestyle='--', alpha=0.3)
plt.show()

This example creates a three-band equalizer typically found in audio applications. It’s perfect for adjusting bass, midrange, and treble frequencies in music production or playback systems.

Read SciPy Find Peaks

Method 8: Analyze FIR Filter Design with Window Methods

When designing FIR filters, I often need to analyze how different window functions affect the frequency response.

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Sample rate and filter parameters
fs = 44100  # Hz
cutoff = 1000  # Hz
normalized_cutoff = cutoff / (fs / 2)

# Design FIR filters with different windows (same filter length)
numtaps = 101
windows = ['hamming', 'blackman', 'hann', 'boxcar']
responses = []

plt.figure(figsize=(12, 8))

for window in windows:
    # Design filter
    b = signal.firwin(numtaps, normalized_cutoff, window=window)
    
    # Calculate frequency response
    w, h = signal.freqz(b, [1], worN=2000)
    f = w * fs / (2 * np.pi)
    
    # Plot
    plt.semilogx(f, 20 * np.log10(abs(h)), label=f'{window.capitalize()} Window')
    responses.append(h)

plt.title('FIR Low-Pass Filter with Different Windows')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude (dB)')
plt.grid(True, which="both", ls="--")
plt.legend()
plt.axvline(cutoff, color='r', linestyle='--', label='Cutoff Frequency')
plt.xlim([20, 20000])
plt.ylim([-120, 10])
plt.show()

This example is particularly useful when designing audio filters where the transition bandwidth and stopband attenuation requirements vary.

Check out Python SciPy Chi-Square Test

Method 9: Implement a Real-world Example – Telephone Channel

Let’s simulate a typical telephone channel filter that limits audio to the standard telephone bandwidth (300-3400 Hz):

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# Sample rate
fs = 44100  # Hz

# Define telephone bandwidth
low_freq = 300   # Hz
high_freq = 3400  # Hz

# Normalize frequencies
low_norm = low_freq / (fs / 2)
high_norm = high_freq / (fs / 2)

# Design bandpass filter
b, a = signal.butter(4, [low_norm, high_norm], 'bandpass')

# Calculate frequency response
worN = 4096
w, h = signal.freqz(b, a, worN=worN)
f = w * fs / (2 * np.pi)

# Plot magnitude response
plt.figure(figsize=(10, 6))
plt.semilogx(f, 20 * np.log10(abs(h)))
plt.title('Telephone Channel Filter Response')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude (dB)')
plt.grid(True, which="both", ls="--")
plt.axvline(low_freq, color='r', linestyle='--', label='Lower Cutoff (300 Hz)')
plt.axvline(high_freq, color='r', linestyle='--', label='Upper Cutoff (3400 Hz)')
plt.legend()
plt.xlim([20, 20000])
plt.ylim([-80, 5])

# Add annotations
plt.annotate('Voice Frequencies', xy=(1000, -5), xytext=(1000, -20),
             arrowprops=dict(facecolor='black', shrink=0.05))
plt.show()

This example simulates the limited frequency response of a traditional telephone channel, which is useful when developing communication systems or audio processing applications.

Check out Python SciPy Exponential

When Should You Use freqz?

From my experience, scipy.signal.freqz is most useful in these situations:

When designing digital filters and you need to visualize their response
When comparing different filter designs to select the best one for your application
When you need to understand how a filter will affect specific frequency components
When debugging filter implementations that aren’t behaving as expected
When teaching or learning about digital signal processing concepts

Understanding frequency responses has helped me immensely in developing better audio processing, communication systems, and sensor data analysis applications.

The scipy.signal library offers a comprehensive set of tools for signal processing, and freqz is one of its most valuable functions for filter analysis. By mastering this function, you can develop more effective signal processing solutions for various applications.

Whether you’re building audio processors, communication systems, or analyzing biomedical signals, understanding frequency responses through tools like freqz will make your signal processing journey much smoother.