最小二乘法作谱估计的MATLAB程序

close all;
clear all;

N=512;
m=1:N;
xn=sqrt（20）*sin（2*pi*0.2*m）+sqrt（2）*sin（2*pi*0.213*m）+randn（1N）; %产生含有噪声的序列xn
Rx=xcorr（xn）;
q=200;
%以下p1、p2分别是AR的阶数
p1=4;
p2=100;
%p1=4时：
p=p1;
R=zeros（pp）;
r=zeros（1p）;
for i=1:p
    for j=1:p
        R（ij）=R（ij）+Rx（1N+q-p+i+j-1）;
    end
end
for m=1:p
    r（m）=Rx（1N+q+m）;
end
A=-（inv（R‘*R））*R‘*r‘;
A（p+1）=1;
figure（1）
[H1w1]=freqz（1A）;
f1=w1/（2*pi）;%归一化
figure（1）
plot（f1