首页 > 解决方案 > Matlab:sign()函数的傅立叶系数正在振荡

问题描述

我正在编写一个接收时变信号并返回 FFT 的函数。它遵循 Matlab 文档:

它适用于非常简单的正弦曲线总和。

对于其他模拟信号,比如阶跃函数,我最终得到了一些 fubar,它以类似锯齿的方式振荡。这是一个使用重载阶跃函数的最小示例,它绘制了 FFT 权重的虚部:

%Create and plot sign function 
Fs = 1000;
dt = 1/Fs;
tvals = -1: dt: 1-dt;
yvals = sign(tvals); %don't use heaviside() it is very slow

subplot(2,1,1);
    plot(tvals, yvals, 'k', 'LineWidth', 2); hold on; 
    axis([-1 1 -1.1 1.1]);

%Calculate center-shifted FFT and plot imaginary part
fftInitial = fft(yvals);
n = length(fftInitial);
dF = Fs/n;
frequenciesShifted = -Fs/2: dF: Fs/2-dF;  %zero-centered frequency range
fftShifted = fftshift(fftInitial/length(yvals));

subplot(2,1,2);
plot(frequenciesShifted, imag(fftShifted), 'b', 'LineWidth', 2);hold on
xlim([-8 8])

这是结果图: 在此处输入图像描述

请注意,已知的想象解是2/jwj(-2/w)

请注意,除了锯齿,这是我最关心的问题,这些权重的包络似乎甚至没有遵循这一点。它实际上似乎围绕原点翻转。不知道我在这里做错了什么组合。

基于一些有用的反馈,人们指出了这个问题:
Analytical Fourier transform vs FFT of functions in Matlab

特别是,在时间数组中不包括零可能会导致问题,这是那里的主要问题。我的代码中的时间数组包含一个零,所以它看起来不像是重复的。我通过时移步骤将我的零推到了数组的前面(尽管坦率地说,我已经做了很多 fft() 没有这样做,而且它们看起来都很好,所以我不认为这是问题所在,但是我现在只想将其作为一个问题删除)。所以我们最终得到:

Fs = 1000;
dt = 1/Fs;
tvals = 0: dt: 2-dt;
yvals = sign(tvals-1); %don't use heaviside() it is very slow
zeroInd = find(yvals == 1, 1, 'first');
yvals(zeroInd-1) =0;
%Calculate center-shifted FFT and plot imaginary part
fftInitial = fft(yvals);
n = length(fftInitial);
dF = Fs/n;
frequenciesShifted = -Fs/2: dF: Fs/2-dF;  %zero-centered frequency range
fftShifted = fftshift(fftInitial/length(yvals));

%Plot stuff
subplot(2,1,1);
plot(tvals, yvals, 'k', 'LineWidth', 2); hold on; 
axis([0 2 -1.1 1.1]);  
subplot(2,1,2);
plot(frequenciesShifted, imag(fftShifted), 'b', 'LineWidth', 2);hold on
xlim([-8 8])
grid on;

我仍然可以使用相同的不正确信封获得相同的锯齿形功能。所以,虽然我确实看到我的问题和那个问题密切相关,但我不确定它们是否重复。而且我真的希望能够绘制这些组件值的中途合理的外观图(我确实处理了一些基本上是阶跃函数的生理信号(与我的测量仪器相比,它们移动得非常快,所以这不是对我来说只是一个学术练习)。

标签: matlabfft

解决方案


您的代码有两个问题:

  1. DFT(FFT 算法计算 DFT)将原点放在最左边的 bin 处。创建yvals这样的原点位于中间会导致输出频谱成为其长度偏移一半的信号的频谱。这会导致非常高频的振荡。修复方法是ifftshift在调用. 在另一个问题中fft查看更多信息

  2. DFT 假设(可以解释为假设)输入信号是周期性的。这导致了第二次大跃进。基本上,您的信号看起来像一个移位的框函数,因此您的变换看起来像一个修改相位的 sinc 函数。解决方案是调用fft. 参见例如this other question

修改您的代码如下:

yvals = sign(tvals);
yvals = yvals .* hanning(numel(yvals), 'periodic').'; % Apply windowing function

% ...

fftInitial = fft(ifftshift(yvals)); % Shift signal before calling FFT

这是您的代码现在给出的输出:

图形


推荐阅读