matlab程序中功率谱分析的经典常用方法

1、一、直接法clear;clc;close all;%清除变量；清屏；关闭当前图形窗口Fs=1000; t=0:1/Fs:1;nfft=2048;%改变nfft的值可对比不同采样值时的谱估计效果%*生成信号、噪声*%x1=cos(2*pi*40*t)+3*cos(2*pi*45*t);%信号x2=randn(size(t);%信号+噪声x3=x1+x2;Pxx,f=periodogram(x3,window,nfft,Fs); %plot(f,10*log10(Pxx); title(直接法 nfft=2048);set(gca,xlim,1 120);ylabel(Am/dB);xlabel(

2、Frequency/Hz);二、间接法Fs=1000;% 采样频率n=0:1/Fs:1;%产生含有噪声的序列x1=cos(2*pi*40*n)+3*cos(2*pi*45*n);%噪声直接法信号 x2=randn(size(n);%噪声-可编辑修改-x3=x1+x2;%信号+噪声 nfft=1024;CXk=fft(cxn); Pxx=abs(CXk);cxn=xcorr(x3);%计算序列的自相关函数index=0:round(nfft/2-1); f=index*Fs/nfft;plot_Pxx=10*log10(Pxx(index+1); figure (1)plot(f,plot_Px

3、x);title(间接法 nfft=1024);ylabel(Am/dB);set(gca,xlim,1 120); xlabel(Frequency/Hz);Bartlett 法clear;clc;close all;%清除变量；清屏；关闭当前图形窗口Fs=1000; t=0:1/Fs:1; nfft=1024;*生成信号、噪声*信号%噪声%信号+噪声%海明窗x仁 cos(2*pi*40*t)+3*cos(2*pi*45*t);% x2=randn(size(t);x3=x1+x2;window=hamming(512);noverlap=0; % 数据无重叠p=0.9; %置信概率 Pxx

4、,Pxxc=psd(x3,nfft,Fs,window,noverlap,p);index=0:round(nfft/2-1);k=index*Fs/nfft;plot_Pxx=10*log10(Pxx(index+1); plot_Pxxc=10*log10(Pxxc(index+1);figure(1)plot(k,plot_Pxx);title(Bartlett法海明窗);set(gca,xlim,1 120); ylabel(Am/dB); xlabel(Frequency/Hz);四、Welch法clear;clc;close all; % 清除变量；清屏；关闭当前图形窗口Fs=10

5、00;t=0:1/Fs:1; nfft=1024;0/ *x2=randn(size(t);%噪声x3=x1+x2;%信号+噪声window=hamming(512);%海明窗noverlap=128;range=onesided;%频率间隔为0 Fs/2，只计算一半的频生成信号、噪声x1=cos(2*pi*40*t)+3*cos(2*pi*45*t);%信号Pxx1,f=pwelch(x3,window,noverlap,nfft,Fs,range);plot_Pxx 1=10*log10(Pxx1); figure(1);plot(f,plot_Pxx1);title(Welch法海明窗)

6、;ylabel(Am/dB); set(gca,xlim,1 120); xlabel(Frequency/Hz);对所给的实测信号进行谱估计，本文采用了周期图法和Welch法。其程序如下所示:、周期图法clear;clc;close all;% 清除变量；清屏；关闭当前图形窗口 Fs=2000; % 采样频率load Chanel8Xia2Data.mat; x1=Chanel8Xia2Gray; x2=Chanel8Xia2Sti; Mlag=length(x1);(%*求功率谱密度 *%Pxx1,f=periodogram(x1,window,length(x1),Fs);%直接法Pxx

7、2,f=periodogram(x2,window,length(x1),Fs);plot_Pxx 1=10*log10(Pxx1);plot_Pxx2=10*log10(Pxx2);%*显示功率谱密度曲线 *%figure(1);plot(f,plot_Pxx1,b); axis(0,100,-50,40);grid on; title( set(gca,xlim,1 100);ylabel(Am/dB);hold on;plot(f,plot_Pxx2,r); axis(0,100,-50,40);grid on;title(直接法);xlabel(Frequency/Hz); ylabe

8、l(Am/dB);直接法);二、Welch 法clear;clc;close all;%清除变量；清屏；关闭当前图形窗口Fs=2000;load Chanel8Xia2Data.mat;x1=Chanel8Xia2Gray;x2=Chanel8Xia2Sti;%海明窗%数据无重叠window仁hamming(1024);noverlap=256;range=onesided;%频率间隔为0 Fs/2，只计算一半频、*求功率谱密度*Pxx1,f=pwelch(x1,window1,noverlap,length(x1),Fs,range);Pxx2,f=pwelch(x2,window1,noverlap,length(x1),Fs,range);plot_Pxx 1=10*log10(Pxx1);plot_Pxx2=10*log10(Pxx2);%*显示功率谱密度曲线*% figure(1);plot(f,plot_Pxx1,b);axis(0,100,-30,30);grid on; title(Welch法海明窗);ylabel(Am/dB);set(gca,xlim,1 100); hold