Matlab中fft函数的用法及关键问题详解
FFT是Fast Fourier Transform(快速傅里叶变换)的简称,这种算法可以减少计算DFT(离散傅里叶变换,关于此更详细的说明见后文)的时间,大大提高了运算效率,并曾经一度被认为是信号分析技术划时代的进步,其重要性由此可见一斑。闲话少叙,言归正传。
基于FFT在信号分析中的重要性,其必然会成为MATLAB的座上宾。FFT算法在MATLAB中实现的函数是Y=fft(x,n)。刚接触频谱分析用到FFT时,几乎都会对MATLAB的fft函数产生一些疑惑,本文本着从问题出发的原则,主要着手对一下几个问题进行解释:
(一)fft函数计算得到的Y是输入信号x的频谱吗?如果不是还要经过怎样的变换?为什么要除以N。
(二)如何计算Y对应的频率f,并绘制(f,Y)频谱图?
(三)如何根据离散信号的长度确定n的数值?
下面以MATLAB帮助文档中的例子来一一看这几个问题。
- Fs = 1000; % Sampling frequency
- T = 1/Fs; % Sample time
- L = 1000; % Length of signal
- t = (0:L-1)*T; % Time vector
- % Sum of a 50 Hz sinusoid and a 120 Hz sinusoid
- x = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t);
- y = x + 2*randn(size(t)); % Sinusoids plus noise
- subplot(2,1,1)
- plot(Fs*t(1:50),x(1:50))
- title('Sinusoids Signal')
- xlabel('time (milliseconds)')
- subplot(2,1,2)
- plot(Fs*t(1:50),y(1:50))
- title('Signal Corrupted with Zero-Mean Random Noise')
- xlabel('time (milliseconds)')
- NFFT = 2^nextpow2(L); % Next power of 2 from length of y
- Y = fft(y,NFFT)/L; % (1)
- f = Fs/2*linspace(0,1,NFFT/2+1);
- % Plot single-sided amplitude spectrum.
- figure
- plot(f,2*abs(Y(1:NFFT/2+1))) % (2)
- title('Single-Sided Amplitude Spectrum of y(t)')
- xlabel('Frequency (Hz)')
- ylabel('|Y(f)|')
1、关于问题(一)
程序(1)处为何要除以信号的采样长度L?
由Fourier变换对
可知,fft函数直接计算得到的X(k)并不是频谱幅值。
在x(j)的Fourier级数(3)中,谐波分量
对应的幅值为X(k)/N,因此必须对fft得到的结果除以离散信号的长度N才能得到频谱幅值。
2、关于问题(三)
程序(1)处,fft函数的第二个参数NFFT为何取值
- NFFT = 2^nextpow2(L);
由(2)可知幅值谱只取了前半部分,并且还要乘以倍数2。也就是说频谱幅值是
- 2*abs(Y(1:NFFT/2+1))
这是因为信号的频谱是前后对称的,而且一般负频率没有物理意义。也就是说绝对值相等的频率,其对应的频率幅值是相等的,所以把正频率对应的频率幅值的两倍作为频谱幅值。
3、关于问题(二)
由离散傅里叶变换(DFT)的推导过程可知,如果信号的时间长度是T0,则用DFT进行频谱分析的频率分辨率是1/T0。由于fft函数取信号的长度是NFFT,采样频率是Fs,由此可知fft的频率分辨率是Fs/NFFT。因此与频谱幅值
- 2*abs(Y(1:NFFT/2+1))
- f=1/(NFFT*T)*(0:NFFT/2);
- f = Fs/2*linspace(0,1,NFFT/2+1);