实验功率谱估计
实验目的:
1、掌握最大熵谱估计的基本原理。

2、了解最终预测误差（FPE）准则。

3、掌握周期图谱估计的基本原理。

4、掌握传统谱估计中直接法与间接法之间的关系。

5、复习快速傅里叶变换与离散傅里叶变换之间关系。

实验内容:
1、设两正弦信号的归一化频率分别为0.175和0.20，用最大熵法编程计算信噪比S/N=30dB、N=32点时该信号的最大熵谱估计结果。

2、用周期图法编程计算上述信号的谱估计结果。

程序示例:
1、最大熵谱估计
clc;
N=32;
SNR=30;
fs=1;
t=1:N;
t=t/fs;
y=sin(2*pi*0.175*t)+sin(2*pi*0.20*t);
x = awgn(y,SNR);
M=1;
P(M)=0;
Rx(M)=0;
for n=1:N
P(M)=P(M)+(abs(x(n)))^2;
ef(1,n)=x(n);
eb(1,n)=x(n);
end
P(M)=P(M)/N;
Rx(M)=P(M);
M=2;
A=0;
D=0;
for n=M:N
A=A+ef(M-1,n)*eb(M-1,n-1);
D=D+(abs(ef(M-1,n)))^2+(abs(eb(M-1,n-1)))^2; end
xishu=-2*A/D;
a(M-1,M-1)=-2*A/D;
P(M)=P(M-1)*(1-(abs(xishu))^2);
FPE(M-1)=P(M)*(N+M)/(N-M);
TH=FPE(M-1);
for n=M:N
ef(M,n)=ef(M-1,n)+xishu*eb(M-1,n-1);
eb(M,n)=eb(M-1,n-1)+xishu*ef(M-1,n);
end
M=M+1;
A=0;
D=0;
for n=M:N
A=A+ef(M-1,n)*eb(M-1,n-1);
D=D+(abs(ef(M-1,n)))^2+(abs(eb(M-1,n-1)))^2;
end
xishu=-2*A/D;
a(M-1,M-1)=-2*A/D;
P(M)=P(M-1)*(1-(abs(xishu))^2);
FPE(M-1)=P(M)*(N+M)/(N-M);
for m=1:M-2
a(M-1,m)=a(M-2,m)+xishu*a(M-2,M-1-m);
end
while FPE(M-1)<TH
TH=FPE(M-1);
for n=M:N
ef(M,n)=ef(M-1,n)+xishu*eb(M-1,n-1);
eb(M,n)=eb(M-1,n-1)+xishu*ef(M-1,n);
end
M=M+1;
A=0;
D=0;
for n=M:N
A=A+ef(M-1,n)*eb(M-1,n-1);
D=D+(abs(ef(M-1,n)))^2+(abs(eb(M-1,n-1)))^2;
end
xishu=-2*A/D;
a(M-1,M-1)=-2*A/D;
P(M)=P(M-1)*(1-(abs(xishu))^2);
FPE(M-1)=P(M)*(N+M)/(N-M);
for m=1:M-2
a(M-1,m)=a(M-2,m)+xishu*a(M-2,M-1-m);
end
end
T=1/fs;sum1=0;
f=0.01:0.01:0.5;
for m=1:M-1;
sum1=sum1+a(M-1,m)*exp(-j*2*pi*m*f*T);
end
s1=(abs(1+sum1)).^2;
s=P(M)*T./s1;
plot(f,10*log10(s),'k');
xlabel('f/fs');
ylabel('功率谱/dB');
2、周期图谱估计
clc;clear;
N=32;
SNR=30;
fs=1;
t=1:N;
t=t/fs;
y=sin(2*pi*0.175*t)+sin(2*pi*0.20*t);
x = awgn(y,SNR);
sum1=0;
f=0.05:0.01:0.5;
for m=1:N
sum1=sum1+x(m)*exp(-j*2*pi*m*f);
end
s=(abs(sum1)).^2/N;
plot(f,10*log10(s),'k');
xlabel('f/fs');
ylabel('功率谱/dB');
实验结果：
1、最大熵法估计结果：
2、周期图法估计结果：。