功率谱密度相关方法的MATLAB实现

xlabel('频率/Hz');ylabel('功率谱/dB');axis([0 4000 -20 0]) title('最大熵法 Order=20滤波后的信号功率谱');gridon %%

%%功率谱密度 %PSD_WELCH方法

Nfft=512;n=0:sigLength-1;t=n/Fs; %数据长度、时间序列 window=hanning(256); %选用的窗口

noverlap=128; %分段序列重叠的采样点数（长度） dflag='none'; %不做趋势处理

[Pxx,Pxxc,f]=psd(x,Nfft,Fs,window,noverlap,0.95); %功率谱估计,并以0.95的置信度给出置信区间，无返回值是绘制出置信区间 figure;subplot(211);

plot(f,10*log10(Pxx)); %绘制功率谱 xlabel('频率/Hz');ylabel('功率谱/dB');

grid on;title('PSD—Welch方法的原始信号功率谱') subplot(212)

[Pxx,Pxxc,f]=psd(y0,Nfft,Fs,window,noverlap,0.95); %功率谱估计,并以0.95的置信度给出置信区间，无返回值是绘制出置信区间 plot(f,10*log10(Pxx)); %绘制功率谱

xlabel('频率/Hz');ylabel('功率谱/dB');axis([0 4000 -30 0]) grid on;title('PSD—Welch方法的滤波后的信号功率谱') %%

%用多窗口法(MTM)

n=0:sigLength-1;t=n/Fs; %数据长度、分段数据长度，时间序列 [Pxx1,f]=pmtm(x,2,Nfft,Fs); %用多窗口法(NW=4)估计功率谱 figure;subplot(2,1,1),plot(f,10*log10(Pxx1)); %绘制功率谱 xlabel('频率/Hz');ylabel('功率谱/dB');

title('多窗口法(MTM) nw=2原始信号功率谱');gridon

[Pxx,f]=pmtm(y0,2,Nfft,Fs); %用多窗口法(NW=2)估计功率谱 subplot(2,1,2),plot(f,10*log10(Pxx)); %绘制功率谱

xlabel('频率/Hz');ylabel('功率谱/dB');axis([0 4000 -80 -20]) title('多窗口法(MTM) nw=2滤波后的信号功率谱');gridon %%

%采用Welch方法估计功率谱 noverlap=128; %重叠数据 dflag='none'; Nfft=1024; figure;subplot(2,1,1)

psd(x,Nfft,Fs,window,noverlap,dflag); %采用Welch方法估计功率谱 xlabel('频率/Hz');ylabel('功率谱/dB') title('Welch方法原始信号功率谱');gridon subplot(2,1,2)

psd(y0,Nfft,Fs,window,noverlap,dflag); %采用Welch方法估计功率谱 xlabel('频率/Hz');ylabel('功率谱/dB'),axis([0 4000 -30 0]) title('Welch方法滤波后的信号功率谱');gridon

pmusic(y0,[7,1.1],Nfft,Fs,32,16); %采用多信号分类法估计功率谱 xlabel('频率/Hz'); ylabel('功率谱/dB') title('通过MUSIC法估计的伪谱')