welch2000.m

经典功率谱估计的matlab仿真。利用的是直接法

clear,t=[0:0.001:2-0.001];
sn=[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1];
xiang=[1.28,1.94,4.07,4.22,3.24,3.04,1.47,1.15,1.13,0.97,4.29,3.09,2.82,3.26,3.93,1.64,0.91,1.79,1.47,4.83,2.13,4.00,2.25,5.30,5.25,3.06,4.69,3.80];%相位值
fu=sn.*[1.46,1.46,1.46,1.46,1.46,1.46,1.46,1.46,1.46,1.66,1.85,1.95,2.15,2.34,2.44,2.63,2.73,2.93,3.12,3.22,3.42,3.61,3.71,3.91,4.10,4.20,4.39,13.96];%幅度值
fi=[0.88,1.52,2.16,2.80,3.44,4.08,4.72,5.36,6.00,6.64,7.28,7.92,8.56,9.20,9.84,10.48,11.12,11.76,12.40,13.04,13.68,14.32,14.96,15.60,16.24,16.88,17.52,25.68];%频率值
k=1:2000;
y(29,2000)=zeros;
m=0;
while m~=28,m=m+1;
y(m+1,k)=fu(m)*sin(2*pi*t(k)*fi(m)+xiang(m))+y(m,k);
end,y;
low=y(29,k);

Fs=1000;
nfft=1024*16;
window=boxcar(500);
%window=hamming(150)
noverlap=250;
pwelch(low,window,noverlap,nfft,Fs);
%axis([0,30,-20,25])