Python scipy利用快速傅里叶变换实现滤波

import numpy as np
from scipy import fft

PI = np.pi*2
fs = 60     #采样频率
T = 100     #采样周期数
N = fs*T    #采样点
t = np.linspace(0, T, N)
noise = 2 * np.random.randn(*t.shape)
s = 2 * np.sin(PI * t) + 3 * np.sin(22 * PI * t) + noise
F = fft.fft(s)
f = fft.fftfreq(N, 1.0/fs)

import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(2,2,1)
ax.plot(t, s)
ax.set_title("t vs s")
f_abs = np.abs(F)
ax = fig.add_subplot(2,2,2)
ax.plot(f, f_abs)
ax.set_title("fs vs |F|")

xlims = [[0,2], [21,23]]
for i, xlim in enumerate(xlims):
    ax = fig.add_subplot(2,2,3+i)
    ax.stem(f, f_abs)
    ax.set_title("fs vs |F|")
    ax.set_xlim(xlim)

plt.show()

fig = plt.figure(1)
f_filt = F * (np.abs(f) < 2)
s_filt = fft.ifft(f_filt)
ax = fig.add_subplot()
ax.plot(t, s, lw=0.2)
ax.plot(t, s_filt.real, lw=2)
ax.set_title("threshold=2")
ax.set_xlim([0,10])
plt.show()