为什么通过清零 FFT 箱进行过滤是个坏主意?

信息处理 fft 过滤器
2022-01-06 20:42:22

通过对信号执行 FFT,将一些 bin 归零,然后执行 IFFT 来过滤信号非常容易。例如:

t = linspace(0, 1, 256, endpoint=False)
x = sin(2 * pi * 3 * t) + cos(2 * pi * 100 * t)
X = fft(x)
X[64:192] = 0
y = ifft(X)

高频分量被这个“砖墙”FFT 滤波器完全去除。

但我听说这不是一个好方法。

  • 为什么这通常是一个坏主意?
  • 在某些情况下这是一个好的或好的选择吗?

[根据 pichenettes 的建议]

3个回答

在频域中将 bin 归零与在频域中乘以矩形窗口相同。在频域中乘以一个窗口与在时域中对该窗口进行变换的循环卷积相同。矩形窗口的变换是 Sinc 函数 (sin(ωt)/ωt)。请注意,Sinc 函数有许多大波纹和延伸时域孔径全宽的波纹。如果可以输出所有这些波纹(振铃)的时域滤波器是一个“坏主意”,那么归零箱也是如此。

对于 FFT 孔径宽度中“区间”或非整数周期的任何光谱内容,这些波纹将是最大的。因此,如果您的原始 FFT 输入数据是该窗口中有些非周期性的任何数据的窗口(例如,大多数非同步采样的“真实世界”信号),那么这些特定的伪影将由归零箱产生。

另一种看待它的方式是,每个 FFT 结果 bin 代表时域中正弦波的某个频率。因此,将 bin 归零将产生与减去该正弦波相同的结果,或者等效地,添加具有精确 FFT bin 中心频率但相位相反的正弦波。请注意,如果时域中某些内容的频率在 FFT 宽度上不是纯整数周期的,那么尝试通过添加一个完全整数周期正弦波的倒数来消除非整数周期信号将产生,而不是静音,而是看起来更像“节拍”音符(不同频率的 AM 调制正弦波)。同样,可能不是想要的。

相反,如果您的原始时域信号只是几个纯未调制的正弦波,它们在 FFT 孔径宽度中都是完全整数周期的,那么归零 FFT 箱将去除指定的无伪影。

这个问题也让我困惑了很久。@hotpaw2 的解释很好。您可能对使用 matlab 进行的简单实验感兴趣。

https://poweidsplearningpath.blogspot.com/2019/04/dftidft.html


更新信息。

为了验证这个事实很简单,我们只需要仔细观察理想(?)带通滤波器的脉冲响应频谱,该滤波器只是将 FFT 箱清零。为什么我需要添加副词“谨慎”?如果我们只是用同样大小的 FFT 来观察脉冲的响应,就会被欺骗,如图1所示。尽管如此,如果我们在观察滤波器的输出时添加 DFT 的阶数,即对脉冲响应进行零填充,我们可以发现所谓的吉布斯现象,即频域波纹,如图2所示。

结果实际上来自窗口效应。如果你想完全理解问题,请参考DSP圣经(1)的第7.6章和第10.1-10.2章。总而言之,这里指出了三个关键点。

  1. 窗口大小和 DFT(FFT) 的阶数完全独立。不要将它们混合在一起。
  2. 窗口的属性(类型/尺寸)决定了 DTFT 的形状。(例如,更宽的主瓣导致频率响应中的瞬态带更宽。)
  3. DFT只是频域中DTFT的采样。此外,DFT 的阶数越高,DFT 的频谱越密集。

因此,借助图2中更密集的光谱,我们可以看到理想(假)带通滤波器的掩模。

在此处输入图像描述欺骗频率。回复。

在此处输入图像描述频率中的吉布斯现象。回复。

(1) Alan V. Oppenheim 和 Ronald W. Schafer。2009. 离散时间信号处理(第 3 版)。Prentice Hall Press,上马鞍河,新泽西州,美国。

fps = 15;

LPF = 1;
HPF = 2;

n = -511:512;
n0 = 0;
imp = (n==n0);

NyquistF = 1/2*fps;

%% Ideal BPF
tmp_N = 512;
tmp_n = 0:1:tmp_N-1;
freq = ( n .* fps) ./ tmp_N;
F = fft(imp, tmp_N);  
F_bpf = IdealBandpassFilter(F, fps, LPF, HPF);
imp_rep =[real(ifft(F_bpf))'];

% Zero padding.
imp_rep2 =[zeros(1,2048) real(ifft(F_bpf))' zeros(1,2048)];

N = 2^nextpow2(length(imp_rep));
F = fft(imp_rep,N);
freq_step = fps/N;
freq = -fps/2:freq_step:fps/2-freq_step;
freq = freq(N/2+1:end)';

figure;
plot(freq,abs(F(1:N/2)));
xlabel('freq(Hz)');
ylabel('mag');
title('Mis leading Freq Response');


N = 2^nextpow2(length(imp_rep2));
F = fft(imp_rep2,N);
freq_step = fps/N;
freq = -fps/2:freq_step:fps/2-freq_step;
freq = freq(N/2+1:end)';

figure;
plot(freq,abs(F(1:N/2)));
xlabel('freq(Hz)');
ylabel('mag');
title('Zero Padding (DFT) with more points');

%% Function
function filered_signal = IdealBandpassFilter(input_signal, fs, w1, w2)

    N = length(input_signal);
    n = 0:1:N-1;
    freq = ( n .* fs) ./ N;

    filered_signal = zeros(N, 1);

    for i = 1:N
        if freq(i) > w1 & freq(i) < w2
            filered_signal(i) = input_signal(i);
        end

    end
end

FFT 给出的时间分辨率很差,即它没有给出特定频率在什么时间存在的信息。它提供了给定信号持续时间的现有频率分量的信息。

通过在 FFT 中将 bin 归零,在时域中的 IFFT 之后分辨率会很差。