2

我很难用 'buttord' 和 'butter' 函数计算巴特沃斯系数。我的目标是从我构建的时间序列中过滤掉噪音。时间序列有一个红噪声分量和一个频率为 0.3 Hz 的正弦信号。时间序列的采样频率为 10 Hz。

按照关于“buttord”的文档http://www.mathworks.com/help/signal/ref/buttord.html我计算了规范的 [n, Wn]（遵循文档的示例 1）：

Wp = 0.33/(sfreq/2); Ws = 0.37/(sfreq/2);
passripp = 0.1; stopatten = 40;
[n,Wn] = buttord(Wp,Ws,passripp,stopatten);
[b,a] = butter(n,Wn);
y_butter = filter(b,a,timeseries(:,2));

用时间绘制 y_butter 只会让我到处都是零！

我尝试使用“freqz”来检查滤波器的频率响应（使用 512 个样本）：

freqz(b,a,512,sfreq)

其绘图显示过渡带在 1 到 4 Hz 之间！

我对过滤器的理解是：

  • 0.3 Hz 信号
  • >> 0.3 Hz 时的噪声
  • 通过从 0 到 0.33 Hz 的所有内容
  • 衰减 0.36 Hz 以上的一切

您的帮助将不胜感激！

数据可以在这里下载：http ://dl.dropbox.com/u/1918592/detrendedTS.mat 'ts' 的第 1 列是时间变量，第 2 列是数据变量

我去除了线性拟合（Matlab'detrend'）以消除一些走开的红噪声行为。

4

1 回答 1

3

您发布的代码会生成第 57 个 (!!!!) 订单过滤器来执行您想要执行的操作。该freqz命令确实表明该过滤器是有效的，但是尝试直接实现这样一个高阶过滤器而不首先使用将其转换为二阶部分tf2sos无疑会导致严重的数值错误。这可能就是您只看到零的原因。如果将滤波器系数转换为二阶部分，然后级联滤波器，则实际上应该使用该freqz命令获得观察到的滤波器输出。

之所以首先获得如此高的滤波器阶数，是因为您将开始带和停止带放置得非常接近。buttord对于简单的低通巴特沃斯滤波器，完全不需要使用该函数。只需使用...

[b,a] = butter(n,cut/(sfreq/2));

...只需选择n给您想要的滚降（6,12,18 dB 等...）。当使用双精度数据表示时，您只会遇到n大约 10 个问题，但您可能会发现一个相对较小的过滤器顺序可以完成您想要的工作。

于 2012-11-08T10:42:12.787 回答