-3

我有一个包含 2048 个 x 和 y 值对的 ascii 文件。我只想知道如何在 MATLAB 中绘制 y 的 fft。我正在编写以下 MATLAB 代码,但无法找到合适的结果。

我怎样才能做到这一点?这是我尝试过的:

I = load('data1.asc');

for i = 1:2048
    y = I(:,2);
end

plot(x)

Fs = 40000;                    
T = 1/Fs;                   
L = 2000;     
NFFT = 2^nextpow2(L);
Y = abs(fft(y,NFFT))/L;
f = Fs/2*linspace(0,1,NFFT/2+1);

figure, plot(f,2*abs(Y(1:NFFT/2+1))) 
axis([0 40000 0 40])
xlabel('Frequency (Hz)')
ylabel('|Y(f)|')
4

2 回答 2

4

与其直接进入 MATLAB 的 FFT 例程,不如考虑使用periodogram函数。当人们说“FFT”时,他们通常指的是 PSD 或周期图,即使用适当的窗口化和 FFT 化样本绘制的功率谱密度图。MATLAB 中的周期图函数会为您处理所有这些细节,即应用窗口函数、计算 FFT、从 FFT 输出导出幅度、适当缩放轴,甚至在需要时绘图。

注意:periodogram在 MATLAB Signal Processing Toolbox 中 - 如果您无权访问它,那么您可以考虑使用具有克隆的Octave(免费 MATLAB 克隆)periodogram,否则您需要自己将各种构建块放在一起:

  • 窗口函数
  • 快速傅里叶变换
  • 计算幅度
  • 注意频率和幅度值的缩放
  • 绘图 PSD
于 2013-07-26T10:16:01.470 回答
1

代码的fft一部分对我来说看起来不错。但是,这一点没有多大意义:

for i = 1:2048
    y = I(:,2);
end

你想在这里做什么?您根本没有i在 for 循环中使用循环索引。

另外,我假设y长度为 2000,你能确认一下吗?否则L = 2000应改为L = length(y). 同样,我假设数据的采样频率是40kHz,否则Fs = 40000是不正确的。

编辑以下评论中的讨论:

使用您提供的数据,我得到相同的结果。我唯一做的就是在分析中排除最后一个数据点,当它下降到零时。你读取数据的方式对我来说仍然没有意义。注意我使用的是 Octave,而不是 MATLAB,但代码应该在 MATLAB 中给出相同的结果。

load('ascii_value.txt')
y = ascii_value(1:end-1,2);
plot(y)
L=length(y);
Fs = 40000;
T = 1/Fs;
NFFT = 2^nextpow2(L);
Y = abs(fft(y,NFFT))/L;
f = Fs/2*linspace(0,1,NFFT/2+1);
figure, plot(f,2*abs(Y(1:NFFT/2+1)))
axis([0 40000 0 40])
xlabel('Frequency (Hz)')
ylabel('|Y(f)|')

信号如下所示:

在此处输入图像描述

和这样的 FFT:

在此处输入图像描述

注意:如果您以 40 kHz 采样,您的 FFT 只能达到 20kHz(奈奎斯特频率)。

于 2013-07-24T08:38:16.203 回答