使用 Python 中的带阻滤波器滤除频率范围并使用傅立叶变换 FFT 确认

Question

假设我有以下信号： y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(100.0 * 2.0*np.pi*x) + 0.2*np.sin(200 * 2.0*np.pi*x) 如何使用 Python 中的带阻滤波器滤除示例 100Hz？在此信号中，峰值位于 50Hz、100Hz 和 200Hz。如果可以使用 FFT 将其可视化以确认该频率已被正确过滤，这将很有帮助。

Answer 1

基于以下答案： Plotting a Fast Fourier Transform in Python 和： Bandstop filter 我写了以下代码：

import pandas as pd
import time
from scipy.signal import lfilter
import matplotlib.pyplot as plt
import scipy
import numpy as np

# In the below lines data are being filtered using Bandstop filter
print("Filtering using Bandstop filter...")
start_filtering_bandstop = time.time()

# Define filtering parameters:
order = 2
fs = 800.0 # sample rate, Hz
lowcut = 90 # desired cutoff frequency of the filter, Hz
highcut = 110

# Define plots
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 7))

# Number of samplepoints
N = 600
# sample spacing
T = 1.0 / fs
x = np.linspace(0.0, N*T, N)
y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(100.0 * 2.0*np.pi*x) + 0.2*np.sin(200 * 2.0*np.pi*x) # You can put there pandas series too...
ax1.plot(x, y, label='Signal before filtering')
print("Calculating FFT, please wait...")
yf = scipy.fftpack.fft(y)
xf = np.linspace(0.0, 1.0/(2.0*T), N/2)


ax1.set_title('Signal')
ax2.set_title('FFT')
ax2.plot(xf, 2.0/N * np.abs(yf[:N//2]), label='Before filtering')

def butter_bandstop_filter(data, lowcut, highcut, fs, order):
    nyq = 0.5 * fs
    low = lowcut / nyq
    high = highcut / nyq
    b, a = scipy.signal.butter(order, [low, high], btype='bandstop')#, fs, )
    y = lfilter(b, a, data)
    return y

print("Filtering signal, please wait...")
signal_filtered = butter_bandstop_filter(y, lowcut, highcut, fs, order)
ax1.plot(x, signal_filtered, label='Signal after filtering')
ax1.set(xlabel='X', ylabel='Signal values')
ax1.legend() # Don't forget to show the legend
ax1.set_xlim([0,0.8])
ax1.set_ylim([-1.5,2])

# Number of samplepoints
N = len(signal_filtered)
# sample spacing
T = 1.0 / fs
x = np.linspace(0.0, N*T, N)
y = signal_filtered
print("Calculating FFT after filtering, please wait...")
yf = scipy.fftpack.fft(y)
xf = np.linspace(0.0, 1.0/(2.0*T), N/2)

# Plot axes...
ax2.plot(xf, 2.0/N * np.abs(yf[:N//2]), label='After filtering')
ax2.set(xlabel='Frequency [Hz]', ylabel='Magnitude')
ax2.legend() # Don't forget to show the legend
plt.savefig('FFT_after_bandstop_filtering.png', bbox_inches='tight', dpi=300) # If dpi isn't set, the script execution will be faster
# Alternatively for immediate showing of plot:
# plt.show()
plt.close()

end_filtering_bandstop = time.time()
print("Data filtered using Bandstop filter in",round(end_filtering_bandstop - start_filtering_bandstop,2),"seconds!")

并获得了以下地块：

正如我们所见，100Hz 已使用带阻滤波器滤除。

Answer 2

为什么在快速傅立叶变换后频率 50 Hz 的幅度从 1 降低到 0.7？

使用 Python 中的带阻滤波器滤除频率范围并使用傅立叶变换 FFT 确认

Filter out range of frequencies using band stop filter in Python and confirm it using Fourier Transform FFT

python

numpy

fft

scipy