当前位置:   article > 正文

MATLAB和scipy计算巴特沃斯通带滤波的比较_scipy滤波和matlab滤波结果不同

scipy滤波和matlab滤波结果不同

引言

最近从MATLAB转到python进行科学计算,来比较下MATLAB和SciPy计算巴特沃斯通带滤波的效果。

MATLAB

clear; close all; clc;

sampRat = 100.0; %采样频率

T = 10;
t = 0:1/sampRat:(T*sampRat-1)/sampRat;
x = 5*sin(2*pi*t*5)+10*sin(2*pi*t*25);


ff = fft(x);
ff = abs(ff);
ff = ff*2/T/sampRat;

f = 0.0:1/T:(sampRat-1/T);
figure()
plot(f, ff);
title('滤波前');

[b,a] = butter(4,[0.01/(sampRat/2),15/(sampRat/2)]);%设计butterworth滤波器(带通)

y =filter(b,a,x); 
fff = fft(y);
fff = abs(fff);
fff = fff*2/T/sampRat;
figure()
plot(f, fff)
title('滤波后');

SciPy

  1. # -*- coding: utf-8 -*-
  2. import scipy.signal as signal
  3. import matplotlib.pyplot as plt
  4. import numpy as np
  5. import math
  6. from scipy.signal import butter, lfilter
  7. def butter_bandpass_filter(data, lowcut, highcut, fs, order=5):
  8. b, a = butter_bandpass(lowcut, highcut, fs, order=order)
  9. y = lfilter(b, a, data)
  10. return y
  11. def butter_bandpass(lowcut, highcut, fs, order=5):
  12. nyq = 0.5 * fs
  13. low = lowcut / nyq
  14. high = highcut / nyq
  15. b, a = butter(order, [low, high], btype='bandpass')
  16. return b, a
  17. sampRat = 100
  18. dt = 1/sampRat
  19. T = 10
  20. t = np.linspace(0, T, T*sampRat, endpoint=False)
  21. y = 5*np.sin(2*math.pi*t*5)+10*np.sin(2*math.pi*t*25)
  22. filterY = butter_bandpass_filter(y, 0.01, 15, sampRat, 4)
  23. f = np.linspace(0, sampRat, T*sampRat, endpoint=False)
  24. ff = np.fft.fft(y)
  25. ff = np.abs(ff)*2/T/sampRat
  26. plt.figure()
  27. plt.plot(f, ff)
  28. plt.title('滤波前的频谱')
  29. plt.show()
  30. ff = np.fft.fft(filterY)
  31. ff = np.abs(ff)*2/T/sampRat
  32. plt.figure()
  33. plt.plot(f, ff)
  34. plt.title('滤波后的频谱')
  35. plt.show()




声明:本文内容由网友自发贡献,不代表【wpsshop博客】立场,版权归原作者所有,本站不承担相应法律责任。如您发现有侵权的内容,请联系我们。转载请注明出处:https://www.wpsshop.cn/w/IT小白/article/detail/106586
推荐阅读
相关标签
  

闽ICP备14008679号