DSP - 使用 FFT 獲得功率頻譜密度估計
阿新 • • 發佈:2020-08-09
具有指定取樣率的偶數長度輸入
對於一個取樣率為 1 kHz 的偶數長度訊號,分別使用fft
和periodogram
獲得其週期圖。比較二者的結果。
建立一個含N(0,1) 加性噪聲的 100 Hz 正弦波訊號。取樣頻率為 1 kHz。訊號長度為 1000 個取樣點。使用隨機數生成器的預設設定以獲得可重現的結果。
rng default Fs = 1000; t = 0:1/Fs:1-1/Fs; x = cos(2*pi*100*t) + randn(size(t));
使用fft
獲取週期圖。訊號是偶數長度的實數值訊號。由於訊號是實數值訊號,您只需要對正負頻率之一進行功率估計。
為了保持總功率不變,將同時在兩組(正頻率和負頻率)中出現的所有頻率乘以因子 2。零頻率 (DC) 和 Nyquist 頻率不會出現兩次。繪製結果。
N = length(x); xdft = fft(x); xdft = xdft(1:N/2+1); psdx = (1/(Fs*N)) * abs(xdft).^2; psdx(2:end-1) = 2*psdx(2:end-1); freq = 0:Fs/length(x):Fs/2;
計算並使用periodogram
繪製週期圖。二者的結果相同。
periodogram(x,rectwin(length(x)),length(x),Fs)
具有歸一化頻率的輸入
通過fft
為使用歸一化頻率的輸入生成周期圖。建立一個帶N(0,1) 加性噪聲的正弦波訊號。該正弦波的角頻率為π/4弧度/取樣點。使用隨機數生成器的預設設定以獲得可重現的結果。
rng default n = 0:999; x = cos(pi/4*n) + randn(size(n));
使用fft
獲取週期圖。訊號是偶數長度的實數值訊號。由於訊號是實數值訊號,您只需要對正負頻率之一進行功率估計。為了保持總功率不變,將同時在兩組(正頻率和負頻率)中出現的所有頻率乘以因子 2。零頻率 (DC) 和 Nyquist 頻率不會出現兩次。繪製結果。
N = length(x); xdft = fft(x); xdft = xdft(1:N/2+1); psdx = (1/(2*pi*N)) * abs(xdft).^2; psdx(2:end-1) = 2*psdx(2:end-1); freq = 0:(2*pi)/N:pi;
計算並使用periodogram
繪製週期圖。二者的結果相同。
periodogram(x,rectwin(length(x)),length(x))
具有歸一化頻率的複數值輸入
使用fft
為具有歸一化頻率的複數值輸入生成周期圖。採用一個帶N(0,1) 復噪聲的復指數訊號,角頻率為π/4弧度/取樣點。採用隨機數生成器的預設設定,以獲得可重現的結果。
rng default n = 0:999; x = exp(1j*pi/4*n) + [1 1j]*randn(2,length(n))/sqrt(2);
使用fft
獲得週期圖。由於輸入是複數值,此處求[0,2π)弧度/取樣點區間內的週期圖。繪製結果。
N = length(x); xdft = fft(x); psdx = (1/(2*pi*N)) * abs(xdft).^2; freq = 0:(2*pi)/N:2*pi-(2*pi)/N;
使用periodogram
獲取並繪製週期圖。比較 PSD 估計。
periodogram(x,rectwin(length(x)),length(x),'twosided')
Reference
1. MathWorks