# Design an IIR Highpass Butterworth Filter using Bilinear Transformation Method in Scipy – Python

• Last Updated : 07 Jan, 2022

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.

## Python3

 `# 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.

## Python3

 `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.

## Python3

 `# 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

## Python3

 `# 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

## Python3

 `# 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

## Python3

 `# 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.

## Python3

 `# 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

## Python3

 `# Call mfreqz function to plot the ` `# magnitude and phase response ` `mfreqz(z, p, Fs) `

Output:

Step 9: Plotting the impulse & step response

## Python3

 `# Call impz function to plot impulse and  ` `# step response of the filter ` `impz(z, p) `

Output:

Below is the implementation:

## Python3

 `# 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:

My Personal Notes arrow_drop_up
Recommended Articles
Page :