Diseñe un filtro elíptico de paso de banda IIR usando Scipy-Python

IIR significa Infinite Impulse Response, es una de las características sorprendentes de muchos sistemas invariantes de tiempo lineal que se caracterizan por tener una respuesta de impulso h(t)/h(n) que no llega a 0 en ninguna etapa sino que persiste indefinidamente. .

¿Qué es el filtro elíptico de paso de banda IIR?

El filtro elíptico es un tipo especial de filtro que se utiliza en el procesamiento de señales digitales cuando se necesita una transición rápida de la banda de paso a la de parada. 

Las especificaciones son las siguientes:  

  • Frecuencia de banda de paso: 1400-2100 Hz
  • Frecuencia de banda de parada: 1050-24500 Hz
  • Ondulación de banda de paso: 0.4dB
  • Atenuación de la banda de parada: 50 dB
  • Frecuencia de muestreo: 7 kHz
  • Trazaremos la respuesta de magnitud y fase del filtro.

Enfoque paso a paso:

Paso 1: Importación de todas las bibliotecas necesarias.

Python3

# import required library
import numpy as np
import scipy.signal as signal
import matplotlib.pyplot as plt

Paso 2: Definición de funciones definidas por el usuario mfreqz() e impz() . El mfreqz es una función para la gráfica de magnitud y fase y el impz es una función para impulso y respuesta de paso].

Python3

# Function to depict magnitude
# and phase plot
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)
 
    # 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()

Paso 3: Definir variables con las especificaciones dadas del filtro.

Python3

# Given specification
 
# Sampling frequency in Hz
Fs = 7000
 
# Pass band frequency in Hz
fp = np.array([1400, 2100])
 
# Stop band frequency in Hz
fs = np.array([1050, 2450])
 
# Pass band ripple in dB
Ap = 0.4
 
# Stop band attenuation in dB
As = 50

Paso 4: Calcule la frecuencia de corte

Python3

# Compute pass band and stop band edge frequencies
 
# Normalized passband edge
# frequencies w.r.t. Nyquist rate
wp = fp/(Fs/2)
 
# Normalized stopband
# edge frequencies
ws = fs/(Fs/2)

Paso 5: Calcular el orden del filtro digital Elliptic Bandpass.

Python3

# Compute order of the elliptic filter
# using signal.ellipord
N, wc = signal.ellipord(wp, ws, Ap, As)
 
# Print the order of the filter and
# cutoff frequencies
print('Order of the filter=', N)
print('Cut-off frequency=', wc)

Paso 6: Diseñe un filtro de paso de banda elíptico digital.

Python3

# Design digital elliptic bandpass filter
# using signal.ellip function
z, p = signal.ellip(N, Ap, As, wc, 'bandpass')
 
 
# Print numerator and denomerator
# coefficients of the filter
print('Numerator Coefficients:', z)
print('Denominator Coefficients:', p)

Paso 7: Grafique la magnitud y la respuesta de fase.

Python3

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

Paso 8: Trace el impulso y la respuesta al escalón del filtro.

Python3

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

A continuación se muestra la implementación completa del enfoque paso a paso anterior:

Python3

# Import required library
import numpy as np
import scipy.signal as signal
import matplotlib.pyplot as plt
 
 
# Function to depict magnitude
# and phase plot
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)
 
    # 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
 
# Sampling frequency in Hz
Fs = 7000
 
# Pass band frequency in Hz
fp = np.array([1400, 2100])
 
# Stop band frequency in Hz
fs = np.array([1050, 2450])
 
# Pass band ripple in dB
Ap = 0.4
 
# Stop band attenuation in dB
As = 50
 
# Compute pass band and
# stop band edge frequencies
# Normalized passband edge frequencies
# w.r.t. Nyquist rate
wp = fp/(Fs/2)
 
# Normalized stopband edge frequencies
ws = fs/(Fs/2)
 
# Compute order of the elliptic filter
# using signal.ellipord
N, wc = signal.ellipord(wp, ws, Ap, As)
 
# Print the order of the filter and cutoff frequencies
print('Order of the filter=', N)
print('Cut-off frequency=', wc)
 
# Design digital elliptic bandpass filter
# using signal.ellip() function
z, p = signal.ellip(N, Ap, As, wc, 'bandpass')
 
 
# Print numerator and denomerator coefficients of the filter
print('Numerator Coefficients:', z)
print('Denominator Coefficients:', p)
 
 
# Depicting visualizations
 
# Call mfreqz 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)

Producción:

Publicación traducida automáticamente

Artículo escrito por sagnikmukherjee2 y traducido por Barcelona Geeks. The original can be accessed here. Licence: CCBY-SA

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *