IIR stands for Infinite Impulse Response, It is one of the striking features of many linear-time invariant systems that are distinguished by having an impulse response h(t)/h(n) which does not become zero after some point but instead continues infinitely.
What is IIR Highpass Butterworth ?
It basically behaves just like an ordinary digital Highpass Butterworth Filter with an infinite impulse response.
The specifications are as follows:
- Pass band frequency: 2-4 kHz
- Stop band frequency: 0-500 Hz
- Pass band ripple: 3dB
- Stop band attenuation: 20 dB
- Sampling frequency: 8 kHz
- We will plot the magnitude, phase, impulse, step response of the filter.
Step-by-step Approach:
Step 1: Importing all the necessary libraries.
# import required library import numpy as np
import scipy.signal as signal
import matplotlib.pyplot as plt
|
Step 2: Defining user-defined functions mfreqz() and impz(). mfreqz is a function for magnitude and phase plot & impz is a function for impulse and step response.
def mfreqz(b, a, Fs):
# Compute frequency response of the filter
# using signal.freqz function
wz, hz = signal.freqz(b, a)
# Calculate Magnitude from hz in dB
Mag = 20 * np.log10( abs (hz))
# Calculate phase angle in degree from hz
Phase = np.unwrap(np.arctan2(np.imag(hz), np.real(hz))) * ( 180 / np.pi)
# Calculate frequency in Hz from wz
Freq = wz * Fs / ( 2 * np.pi) # START CODE HERE ### (≈ 1 line of code)
# Plot filter magnitude and phase responses using subplot.
fig = plt.figure(figsize = ( 10 , 6 ))
# Plot Magnitude response
sub1 = plt.subplot( 2 , 1 , 1 )
sub1.plot(Freq, Mag, 'r' , linewidth = 2 )
sub1.axis([ 1 , Fs / 2 , - 100 , 5 ])
sub1.set_title( 'Magnitude Response' , fontsize = 20 )
sub1.set_xlabel( 'Frequency [Hz]' , fontsize = 20 )
sub1.set_ylabel( 'Magnitude [dB]' , fontsize = 20 )
sub1.grid()
# Plot phase angle
sub2 = plt.subplot( 2 , 1 , 2 )
sub2.plot(Freq, Phase, 'g' , linewidth = 2 )
sub2.set_ylabel( 'Phase (degree)' , fontsize = 20 )
sub2.set_xlabel(r 'Frequency (Hz)' , fontsize = 20 )
sub2.set_title(r 'Phase response' , fontsize = 20 )
sub2.grid()
plt.subplots_adjust(hspace = 0.5 )
fig.tight_layout()
plt.show()
# Define impz(b,a) to calculate impulse response # and step response of a system input: b= an array # containing numerator coefficients,a= an array containing #denominator coefficients def impz(b, a):
# Define the impulse sequence of length 60
impulse = np.repeat( 0. , 60 )
impulse[ 0 ] = 1.
x = np.arange( 0 , 60 )
# Compute the impulse response
response = signal.lfilter(b, a, impulse)
# Plot filter impulse and step response:
fig = plt.figure(figsize = ( 10 , 6 ))
plt.subplot( 211 )
plt.stem(x, response, 'm' , use_line_collection = True )
plt.ylabel( 'Amplitude' , fontsize = 15 )
plt.xlabel(r 'n (samples)' , fontsize = 15 )
plt.title(r 'Impulse response' , fontsize = 15 )
plt.subplot( 212 )
step = np.cumsum(response) # Compute step response of the system
plt.stem(x, step, 'g' , use_line_collection = True )
plt.ylabel( 'Amplitude' , fontsize = 15 )
plt.xlabel(r 'n (samples)' , fontsize = 15 )
plt.title(r 'Step response' , fontsize = 15 )
plt.subplots_adjust(hspace = 0.5 )
fig.tight_layout()
plt.show()
|
Step 3:Define variables with the given specifications of the filter.
# Given specification Fs = 8000 # Sampling frequency in Hz
fp = 2000 # Pass band frequency in Hz
fs = 500 # Stop Band frequency in Hz
Ap = 3 # Pass band ripple in dB
As = 20 # Stop band attenuation in dB
# Compute Sampling parameter Td = 1 / Fs
|
Step 4:Computing the cut-off frequency
# Compute cut-off frequency in radian/sec wp = 2 * np.pi * fp # pass band frequency in radian/sec
ws = 2 * np.pi * fs # stop band frequency in radian/sec
|
Step 5: Pre-wrapping the cut-off frequency
# Prewarp the analog frequency Omega_p = ( 2 / Td) * np.tan(wp * Td / 2 ) # Prewarped analog passband frequency
Omega_s = ( 2 / Td) * np.tan(ws * Td / 2 ) # Prewarped analog stopband frequency
|
Step 6: Computing the Butterworth Filter
# Compute Butterworth filter order and cutoff frequency N, wc = signal.buttord(Omega_p, Omega_s, Ap, As, analog = True )
# Print the values of order and cut-off frequency print ( 'Order of the filter=' , N)
print ( 'Cut-off frequency=' , wc)
|
Output:
Step 7: Design analog Butterworth filter using N and wc by signal.butter() function.
# Design analog Butterworth filter using N and # wc by signal.butter function b, a = signal.butter(N, wc, 'high' , analog = True )
# Perform bilinear Transformation z, p = signal.bilinear(b, a, fs = Fs)
# Print numerator and denomerator coefficients # of the filter print ( 'Numerator Coefficients:' , z)
print ( 'Denominator Coefficients:' , p)
|
Output:
Step 8: Plotting the Magnitude & Phase Response
# Call mfreqz function to plot the # magnitude and phase response mfreqz(z, p, Fs) |
Output:
Step 9: Plotting the impulse & step response
# Call impz function to plot impulse and # step response of the filter impz(z, p) |
Output:
Below is the implementation:
# import required library import numpy as np
import scipy.signal as signal
import matplotlib.pyplot as plt
# User defined functions mfreqz for # Magnitude & Phase Response def mfreqz(b, a, Fs):
# Compute frequency response of the filter
# using signal.freqz function
wz, hz = signal.freqz(b, a)
# Calculate Magnitude from hz in dB
Mag = 20 * np.log10( abs (hz))
# Calculate phase angle in degree from hz
Phase = np.unwrap(np.arctan2(np.imag(hz), np.real(hz))) * ( 180 / np.pi)
# Calculate frequency in Hz from wz
Freq = wz * Fs / ( 2 * np.pi) # START CODE HERE ### (≈ 1 line of code)
# Plot filter magnitude and phase responses using subplot.
fig = plt.figure(figsize = ( 10 , 6 ))
# Plot Magnitude response
sub1 = plt.subplot( 2 , 1 , 1 )
sub1.plot(Freq, Mag, 'r' , linewidth = 2 )
sub1.axis([ 1 , Fs / 2 , - 100 , 5 ])
sub1.set_title( 'Magnitude Response' , fontsize = 20 )
sub1.set_xlabel( 'Frequency [Hz]' , fontsize = 20 )
sub1.set_ylabel( 'Magnitude [dB]' , fontsize = 20 )
sub1.grid()
# Plot phase angle
sub2 = plt.subplot( 2 , 1 , 2 )
sub2.plot(Freq, Phase, 'g' , linewidth = 2 )
sub2.set_ylabel( 'Phase (degree)' , fontsize = 20 )
sub2.set_xlabel(r 'Frequency (Hz)' , fontsize = 20 )
sub2.set_title(r 'Phase response' , fontsize = 20 )
sub2.grid()
plt.subplots_adjust(hspace = 0.5 )
fig.tight_layout()
plt.show()
# Define impz(b,a) to calculate impulse # response and step response of a system # input: b= an array containing numerator # coefficients,a= an array containing #denominator coefficients def impz(b, a):
# Define the impulse sequence of length 60
impulse = np.repeat( 0. , 60 )
impulse[ 0 ] = 1.
x = np.arange( 0 , 60 )
# Compute the impulse response
response = signal.lfilter(b, a, impulse)
# Plot filter impulse and step response:
fig = plt.figure(figsize = ( 10 , 6 ))
plt.subplot( 211 )
plt.stem(x, response, 'm' , use_line_collection = True )
plt.ylabel( 'Amplitude' , fontsize = 15 )
plt.xlabel(r 'n (samples)' , fontsize = 15 )
plt.title(r 'Impulse response' , fontsize = 15 )
plt.subplot( 212 )
step = np.cumsum(response) # Compute step response of the system
plt.stem(x, step, 'g' , use_line_collection = True )
plt.ylabel( 'Amplitude' , fontsize = 15 )
plt.xlabel(r 'n (samples)' , fontsize = 15 )
plt.title(r 'Step response' , fontsize = 15 )
plt.subplots_adjust(hspace = 0.5 )
fig.tight_layout()
plt.show()
# Given specification Fs = 8000 # Sampling frequency in Hz
fp = 2000 # Pass band frequency in Hz
fs = 500 # Stop Band frequency in Hz
Ap = 3 # Pass band ripple in dB
As = 20 # Stop band attenuation in dB
# Compute Sampling parameter Td = 1 / Fs
# Compute cut-off frequency in radian/sec wp = 2 * np.pi * fp # pass band frequency in radian/sec
ws = 2 * np.pi * fs # stop band frequency in radian/sec
# Prewarp the analog frequency Omega_p = ( 2 / Td) * np.tan(wp * Td / 2 ) # Prewarped analog passband frequency
Omega_s = ( 2 / Td) * np.tan(ws * Td / 2 ) # Prewarped analog stopband frequency
# Compute Butterworth filter order and cutoff frequency N, wc = signal.buttord(Omega_p, Omega_s, Ap, As, analog = True )
# Print the values of order and cut-off frequency print ( 'Order of the filter=' , N)
print ( 'Cut-off frequency=' , wc)
# Design analog Butterworth filter using N and # wc by signal.butter function b, a = signal.butter(N, wc, 'high' , analog = True )
# Perform bilinear Transformation z, p = signal.bilinear(b, a, fs = Fs)
# Print numerator and denomerator coefficients of the filter print ( 'Numerator Coefficients:' , z)
print ( 'Denominator Coefficients:' , p)
# Call mfreqz function to plot the magnitude # and phase response mfreqz(z, p, Fs) # Call impz function to plot impulse and step # response of the filter impz(z, p) |
Output: