信息处理 - 一个行为如此奇怪的简单低通滤波器 - 吾爱随笔录

一个行为如此奇怪的简单低通滤波器

信息处理 matlab 频谱过滤器设计低通滤波器有限脉冲响应

2022-01-30 07:23:54

我对过滤器设计了解不多。在阅读了几篇关于 wikipedia 的文章后，我正在玩过滤器设计，我看到了一些奇怪的东西，并希望得到一些帮助来理解这里到底发生了什么。所以我有一个信号在时间上均匀测量（每 30 秒），总共有 1006 次测量，我想最终估计它的功率谱密度。但是我只想查看某些高频，所以我决定通过一个低通滤波器运行它以获得零通道和相邻小频率的功率。然后我将减去滤波后的信号并查看较高频率的功率残差。

首先，我刚刚编写了自己的天真理想低通滤波器，我在其中进行 FFT，将所有高于 0.8 mHz 的频率归零，然后进行逆 FFT，这就是我得到的图片。理想的过滤器

这张照片对我来说非常有意义，因为滤波后的信号只是围绕原始信号振荡，而这些就是我看到的低频。

接下来，我使用 MATLAB 的 DSP System Toolbox 中包含的基于 GUI 的滤波器构建器来设计我自己的滤波器。因此，我使用等波纹方法制作了一个低通、有限脉冲响应、单速率滤波器，通过的频率幅度为 1dB，停止频率幅度衰减 60dB。由于范围必须以严格的正数开始，我使用了 10e-9 Hz 到 0.0008 Hz 的允许范围。滤波器实现为直接形式的 FIR。这是频率响应。设计滤波器的频率响应

我的问题是，当我通过这个设计的过滤器运行相同的数据时，过滤后的信号看起来很奇怪。设计的过滤数据

有人可以告诉我这里发生了什么吗？由于我使用的是 MATLAB 自己的实现，我认为数学或计算没有任何问题，所以这一定是我自己的理解。这是预期的吗？是不是任何地方都没有问题，而图片就是它本来的样子？我应该查看其他过滤器吗？我期待这张照片看起来或多或少像第一张照片，其中一些带有背景功率的低频只是被覆盖了。

我在想要查看的特定频率范围内使用高通滤波器甚至带通滤波器，但它们的行为都很奇怪。跟信号本身有关系吗？背景功率是巨大的，所以它是否会在整个频谱范围内泄漏并弄乱一切？谢谢！

编辑 1

感谢大家的意见。我的原始信号是那个蓝色的 U 形信号。如果你们想要的话，我可以在我弄清楚如何制作可折叠部分后立即在此处发布测量结果。最初的问题是估计 0.8mHz 和 8mHz 之间的功率谱密度。我的第一直觉是只有一个低通滤波器，然后减去它。然后对残差使用multitaper方法来估计PSD。由于去趋势，我想到了减法。由于通常会去除数据的趋势，因此我决定使用过滤后的信号去除趋势。

首先，我编写了自己的简单数字滤波器，在频率空间中进行滤波。这是 MATLAB 代码。

Fs = 1/30;         % Sampling Frequency - per second
L = length(data);  % Length of the signal
fftb = fft(data);
f = Fs/2*linspace(0,1,L/2+1);
f(end+1:2*end-1)=fliplr(f(2:end));
fftb(f(:)>0.0008)=0;
naivefiltered = real(ifft(fftb));

只需采用 FFT，将我不想要的频率归零，然后采用逆 FFT，有时由于某种原因很复杂，所以我采用实部，我得到了上面的第一张图片。

然后我使用 MATLAB 的滤波器设计工具来设计这个滤波器，这里是代码。

Fpass = 3.3135e-09;  % Passband Frequency
Fstop = 0.0008;      % Stopband Frequency
Apass = 1;           % Passband Ripple (dB)
Astop = 60;          % Stopband Attenuation (dB)
Fs    = 1/30;        % Sampling Frequency
h = fdesign.lowpass('fp,fst,ap,ast', Fpass, Fstop, Apass, Astop, Fs);
Hd = design(h, 'equiripple','MinOrder', 'any','StopbandShape', 'flat');

然后对于第三张图片，我正在使用

smbtotal = filter(Hd,data);

但是在 endolith 发表评论之后（我什至不知道 filtfilt）我读到了 filtfilt 并看到“filter()”引入了巨大的时间延迟，所以当我使用

smbtotal = filtfilt(Hd.Numerator,1,data);

图片确实看起来好多了。这里是。

在此处输入图像描述

所以看起来这个时间/相移是问题所在。所以现在我有三个问题。

1.对于您在此处看到的数据类型（U 形范围超过几个数量级），如果我想估计特定频段（0.8mHZ 到 8mHz）的功率，哪种技术更好？我应该减去低通滤波信号吗？我应该对信号进行高通滤波吗？或者我应该对它进行带通滤波，然后估计 PSD？

2.设计过滤器并非易事。但是为什么应用过滤器并不重要？为什么使用“filter()”会引入这个时间/相移？

3.我是在物理空间中卷积还是在频域中相乘有关系吗？哪一个更好？在我看来，如果我在频域中相乘（就像我用我的朴素滤波器做的那样），时间/相移没有问题。但是，如果我在物理空间中进行卷积（例如在 MATLAB 中使用“filter()”），那么您会遇到时间/相移问题吗？那么在频率空间中相乘总是更好吗？

感谢大家的时间。欣赏它！

2个回答

实时应用的真实滤波器通过截断和加窗无限脉冲响应来逼近理想滤波器，以产生有限脉冲响应；应用该滤波器需要将信号延迟一段适中的时间，让计算“看到”一点未来。这种延迟表现为相移。

您可以手动计算相移（即按照以下步骤手动计算）。

计算滤波器的傅里叶变换。
计算傅里叶变换的角度（或θ）。这意味着 F(w) = |F(w)|(theta)
theta 是信号的相移。

然后，您可以通过将相移角乘以采样（？）频率来计算时间延迟。

如果您想在没有时间延迟（或相移）的情况下工作，您应该查看可以通过双向滤波形成的零相位滤波器。

当您在频域中进行乘法运算时，您会忽略时域信号的某些方面，因此无法获得滤波器的真实表示。

例如：在频域中，我可以轻松定义 F(w)，使得 0:Fcutoff 为 0，Fcutoff:infinity 为 1。

频域中非理想高通滤波器的示例算法为：

 // Return RC high-pass filter output samples, given input samples,
 // time interval dt, and time constant RC
 function highpass(real[0..n] x, real dt, real RC)
   var real[0..n] y
   var real α := RC / (RC + dt)
   y[0] := x[0]
   for i from 1 to n
     y[i] := α * y[i-1] + α * (x[i] - x[i-1])
   return y

因此避免单独在傅里叶域中设计高通滤波器。这里 alpha 表示您必须人为引入的延迟。

这将是一个理想的 1 级高通滤波器，由于采样带来的复杂性，这在现实生活中是无法实现的。

您提到您想询问信号的频率响应，特别是您想知道信号高频分量的电平。（来自：“但是我只想查看某些高频”）。

在这种情况下，如果您只是使用 MATLAB 的内置函数计算信号的 DFT，这将告诉您信号的频率响应。例如，您可以这样做：

freqz(signal, 1, 2^nextpow2(length(signal)), your_sampling_frequency_here);

据我了解，这应该可以解决您的问题。从那里，您可以简单地研究信号的高频分量。

作为旁注，您通常不能过滤信号，然后减去过滤后的版本并查看残差。这是因为过滤信号的行为会引入时间/相位延迟，这意味着过滤后的信号和原始信号之间的减法会相应地失真。您可以使用一些技术来考虑这些延迟，但作为一般情况，请注意不要做这些事情。

其它你可能感兴趣的问题

上一篇使用 FFT 进行时域插值，最后填充零下一篇在 MATLAB 中估计 16-QAM 的最佳采样点