我正在尝试从采样率为 500Hz 的 10 分钟长 EEG 信号中过滤 theta 范围(3-8 Hz)。这是我的代码。请帮助我了解出了什么问题。现在过滤后的信号似乎被破坏了。非常感谢!
fs=500;
Wp = [3 8]/(fs/2); Ws = [2.5 8.5]/(fs/2);
Rp = 3; Rs = 40;
[n,Wn] = buttord(Wp,Ws,Rp,Rs);
[b,a] = butter(n,Wn,'bandpass');
fdata = filter(b,a,data);
x=0:ts:((length(data)/fs)-ts);
f=-fs/2:fs/(length(data)-1):fs/2;
subplot(2,2,1)
plot(x,data)
subplot(2,2,2)
z1=abs(fftshift(fft(data)));
plot(f,z1)
xlim([0 150]);
subplot(2,2,3)
plot(x,fdata)
subplot(2,2,4)
z=abs(fftshift(fft(fdata)));
plot(f,z);
xlim([0 150]);
最佳答案
您的代码(第 4 行)给出了一个滤波器阶数 n ，等于 37。我遇到了 数值精度 与如此大阶数的巴特沃斯滤波器的问题；即使订单低至 8。问题是 butter 为大订单提供了荒谬的 b 和 a 值。检查您的 b 和 a 向量，您会看到它们包含大约 1e21 (!)
解决方案是使用滤波器的 零极点表示 ，而不是系数 ( b , a ) 表示。您可以阅读有关此 here 的更多信息。特别是，
在您的情况下，您可以按照以下方式进行:
[z, p, k] = butter(n,Wn,'bandpass');
[sos,g] = zp2sos(z,p,k);
filt = dfilt.df2sos(sos,g);
fdata = filter(filt,data)
关于matlab - matlab中的脑电带通滤波器，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/23664631/