%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%題目：????????????????????????????????????????????????????????????????????%
%已知AR（4）過程：???????????????????????????????????????????????????????????%
%x（n）=2.760x（n-1）-3.809x（n-2）+2.654x（n-3）-0.924x（n-4）+w（n）????????????????%
%?試用Burg算法利用Matlab估計(jì)模型參數(shù)、階數(shù)和功率譜。?????????????????????????%
%實(shí)驗(yàn)報(bào)告要求：實(shí)驗(yàn)仿真100次，但不要求給出所有100次結(jié)果，?????????????????????%
%只給出估計(jì)結(jié)果的方差和偏差（用估計(jì)結(jié)果均值與真值的差來估計(jì)偏差）。????????????%
%實(shí)驗(yàn)報(bào)告要求寫清楚基本原理，提交程序并要有關(guān)鍵注釋。?????????????????????????%
%電子通信工程?研1205?于飛?2012020148????????????????????????????????????????%
%2012年10月11日%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear?allclose?allclc;

sum=zeros（5001）;%psd估計(jì)值寄存器
for?o=1:1:100
????%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
????%產(chǎn)生均值為零、方差為1、服從高斯（正態(tài)）分布的白噪聲信號(hào)u（n）%
????%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
????w=randn（150000）;
????w=w-mean（w）;
????power_u=var（w）;
????m_u=mean（w）;
????%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
????%使u（n）經(jīng)過系統(tǒng)?x（n）=2.760x（n-1）-3.809x（n-2）+2.654x（n-3）-0.924x（n-4）+w（n）%
????%得到所要研究的隨機(jī)過程x（n）%
????%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
????c=1;

????%d=[1-1.3521.338-0.6620.240];%寬帶過程
????d=[1-2.7603.809-2.6540.924];%窄帶過程
????y=filter（cdw）;
????%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
????%burg算法計(jì)算參數(shù)、階數(shù)、功率譜?????????????????????????????????????????????%
????%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
????N=256;
????x=y（25000:25000+N-1）;%截取長(zhǎng)度為256的樣本
????f=x;%初始條件，前向預(yù)測(cè)值
????g=x;%初始條件，后向預(yù)測(cè)值
????p=（1/N）*x*x‘;%初始誤差功率
????a=zeros（110）;%產(chǎn)生用來估計(jì)系數(shù)的向量a
????b=zeros（110）;%產(chǎn)生用來估計(jì)系數(shù)的向量b
????a（1）=1;%初始條件，初始系數(shù)a0=1;

????for?m=1:1:10%循環(huán)遞推，最多10階
????????num=-2*f（m+1:N）*g（m:N-1）‘;%計(jì)算反饋系數(shù)k的分子
????????den=f（m+1:N）*f（m+1:N）‘+g（m:N-1）*g（m:N-1）‘;%計(jì)算反射系數(shù)k的分母
????????k=num/den;%計(jì)算反射系數(shù)
????%%%%%%%%%%%%%%%%%%%%%%-----------計(jì)算系數(shù)a------------%%%%%%%%%%%%%%%%%%%%%%
????????for?n=1:1:m-1
???????????b（n）=a（m-n）;
????????end;
????????if?m==1
????????????a（1）=k;
????????else
????????????b（m）=1;
????????end;

????????a=a+k*b;
????%%%%%%%%%%%%%%%%%%%%%%------------------------------%%%%%%%%%%%%%%%%%%%
????%%%%%%%%%%%%%%%%%-----------計(jì)算前、后項(xiàng)預(yù)測(cè)值-----------%%%%%%%%%%%%%%%
????????ft=f;
????????gt=g;
????????f=ft+k*[0gt（1:N-1）];
????????g=k*ft+[0gt（1:N-1）];
????%%%%%%%%%%%%%%%%%%%%%%------------------------------%%%%%%%%%%%%%%%%%%%
????????pn=（1-k^2）*p;%計(jì)算誤差功率
????????if?p-pn?????????????break;
????????end;
????????p=pn;
????end;
????figure（1）
????a=[1a（1:m）];
????[HW]=freqz（1a500）;
????psd=power_u*（abs（H）.^2）;
????sum?=sum?+psd;
????%plot（（1/（2*pi））*W10*log10（psd））;%用估計(jì)的參數(shù)估計(jì)功率譜并繪出藍(lán)色
????hold?on;
????if?o==1
????????[HW]=freqz（cd500）;
????????psdn=power_u*（abs（H）.^2）;
????????plot（（

?屬性????????????大小?????日期????時(shí)間???名稱
-----------?---------??----------?-----??----
?????目錄???????????0??2012-11-10?10:51??burgyf\
?????文件????????3493??2012-10-14?14:24??burgyf\ex2.asv
?????文件????????3537??2012-10-14?14:30??burgyf\ex2.m
?????文件???????69120??2012-11-08?22:25??burgyf\實(shí)驗(yàn)報(bào)告.doc