%?預報誤差松弛法算法辨識模型參數
%?辨識模型：z（k）-2.725*z（k-1）+2.469*z（k-2）-0.7441*z（k-3）=0.000861*u（k-1）+
%?0.003201*u（k-2）+0.0007421*u（k-3）+e（k）
%?
%?單輸入單輸出
%============================================
close?all
clear
clc
%==========產生M序列作為輸入===============
x1=[0?1?0?1?1?0?1?1?1];??%initial?value
n=504;?%n為脈沖數目
M1=[];??%存放M序列
for?i=1:n
temp=xor（x1（4）x1（9））;
M1（i）=x1（9）;
for?j=9:-1:2
x1（j）=x1（j-1）;
end
x1（1）=temp;
end
u=M1;

%產生均值為0，方差為0.1的高斯白噪聲
E=randn（1504）;

%帶辨識系數矩陣---兩輸入兩輸出
A1=-2.725;
A2=2.469;
A3=-0.7441;
B1=0.000861;
B2=0.003201;
B3=0.0007421;
D1=0;
%產生觀測矩陣
z=zeros（1504）;
z（:1）=-1;
z（:2）=-3;
z（:3）=2;
for?i=4:504
????z（i）=-A1*z（i-1）-A2*z（i-2）-A3*z（i-3）+B1*u（i-1）+B2*u（i-2）+B3*u（i-3）+E（i）;?
end
%=============開始辨識=======================
%-----首先估計A?B系數--------
COV=0;
for?k=1:500
COV=COV+E（k）*E（k）‘;
end
disp（‘噪聲協方差矩陣：‘）