Estoy tratando de crear un filtro de paso bajo simple, pero obtuve lo que considero un resultado sorprendente al observar la respuesta de frecuencia de un filtro Butterworth simple.
He copiado gran parte del ejemplo a continuación de esta otra publicación . He agregado un código en la parte inferior del script para comparar los espectros de entrada y salida con la respuesta de frecuencia del filtro. Esperaría que el espectro de salida sea el producto del espectro de entrada y la respuesta de frecuencia :
Sin embargo, el gráfico a continuación muestra que el filtro en realidad aumenta algunos componentes de baja frecuencia; vea cómo la línea roja está por encima del verde debajo de .
¿Alguien podría explicar por qué es esto?
import numpy as np
from scipy.signal import butter, lfilter, freqz
import matplotlib.pyplot as plt
from scipy.fftpack import fft as fft
def butter_lowpass(cutoff, fs, order=5):
nyq = 0.5 * fs
normal_cutoff = cutoff / nyq
b, a = butter(order, normal_cutoff, btype='low', analog=False)
return b, a
def butter_lowpass_filter(data, cutoff, fs, order=5):
b, a = butter_lowpass(cutoff, fs, order=order)
y = lfilter(b, a, data)
return y
# Filter requirements.
order = 6
fs = 30.0 # sample rate, Hz
cutoff = 3.667 # desired cutoff frequency of the filter, Hz
# Get the filter coefficients so we can check its frequency response.
b, a = butter_lowpass(cutoff, fs, order)
# Plot the frequency response.
w, h = freqz(b, a, worN=8000)
plt.subplot(2, 1, 1)
plt.plot(0.5*fs*w/np.pi, np.abs(h), 'b')
plt.plot(cutoff, 0.5*np.sqrt(2), 'ko')
plt.axvline(cutoff, color='k')
plt.xlim(0, 0.5*fs)
plt.title("Lowpass Filter Frequency Response")
plt.xlabel('Frequency [Hz]')
plt.grid()
# Demonstrate the use of the filter.
# First make some data to be filtered.
T = 5.0 # seconds
n = int(T * fs) # total number of samples
t = np.linspace(0, T, n, endpoint=False)
# "Noisy" data. We want to recover the 1.2 Hz signal from this.
data = np.sin(1.2*2*np.pi*t) + 1.5*np.cos(9*2*np.pi*t) + 0.5*np.sin(12.0*2*np.pi*t)
# Filter the data, and plot both the original and filtered signals.
y = butter_lowpass_filter(data, cutoff, fs, order)
plt.subplot(2, 1, 2)
plt.plot(t, data, 'b-', label='data')
plt.plot(t, y, 'g-', linewidth=2, label='filtered data')
plt.xlabel('Time [sec]')
plt.grid()
plt.legend()
plt.subplots_adjust(hspace=0.35)
plt.show()
def calculateFFT(time,signal):
N=len(signal)
df=1/((time[-1]-time[0]))
frequencies=[i*df for i in range(int(N/2.0))]
fftValues = [2.0/N*abs(i) for i in fft(signal,N)[0:N/2.0] ]
return frequencies,fftValues
plt.subplot(2, 1, 1)
originalfreqs,originalFFT=calculateFFT(t,data)
plt.plot(originalfreqs,originalFFT,"g",label="original")
filteredfreqs,filteredFFT=calculateFFT(t,y)
plt.plot(filteredfreqs,filteredFFT,"r",label="filtered")
plt.legend()
fuente