matlab求单边功率谱

求功率谱有好几种方法，本文例两种，一种是先求信号频谱，然后取平方后平均；第二种是用welch法（可指定各种窗）。

rng default;
fs = 1000;   % 采样率
t = 0:1/fs:1-1/fs;
% x = cos(2*pi*50*t)+ randn(size(t));
x = cos(2*pi*50*t);
nfft = length(x);  %指定傅里叶变换的点数（点数越大频率分辨率越高）

figure(1);  %画单边频谱图
sig_dft = abs(fft(x, nfft));
sig_dft(2:end-1) = 2*sig_dft(2:end-1);  %因为只显示单边谱，所以需要把数值*2，折起来
sig_amp = mapminmax(sig_dft/nfft);  %归一化
x1 = (0:nfft-1)*fs/nfft; 
plot(x1(1:250),sig_amp(1:250));  
xlabel('Frequency(Hz)');
ylabel('Amplitude');
title('Single sided Amplitude Spectrum');

figure(2);  %常规方法画单边功率谱
xdft = fft(x);
xdft = xdft(1:nfft/2+1);
psdx = (1/(fs*nfft)) * abs(xdft).^2;
psdx(2:end-1) = 2*psdx(2:end-1);
freq = 0:fs/nfft:fs/2;
semilogy(freq,psdx);  %semilogy以对数方式绘制纵轴

figure(3);  %修正周期图功率谱密度估计方法
[p,f]=pwelch(x, hanning(nfft), nfft/2, nfft, fs);  %指定不同的窗会影响功率谱结果
semilogy(f,