%LSeg2
clear%清理工作間變量
L=15;% M序列的周期
y1=1;y2=1;y3=1;y4=0;%四個移位積存器的輸出初始值
for i=1:L;%開始循環，長度為L
    x1=xor(y3,y4);%第一個移位積存器的輸入是第3個與第4個移位積存器的輸出的“或”
    x2=y1;%第二個移位積存器的輸入是第3個移位積存器的輸出
    x3=y2;%第三個移位積存器的輸入是第2個移位積存器的輸出
    x4=y3;%第四個移位積存器的輸入是第3個移位積存器的輸出
    y(i)=y4;%取出第四個移位積存器幅值為"0"和"1"的輸出信號，
    if y(i)>0.5,u(i)=-0.03;%如果M序列的值為"1"時,辨識的輸入信號取“-0.03”
    else u(i)=0.03;%當M序列的值為"0"時,辨識的輸入信號取“0.03”
    end%小循環結束
    y1=x1;y2=x2;y3=x3;y4=x4;%為下一次的輸入信號做準備
end%大循環結束，產生輸入信號u 
z(2)=0;z(1)=0;%取z的前兩個初始值為零
for k=3:15;%循環變量從3到15   
    z(k)=1.5*z(k-1)-0.7*z(k-2)+u(k-1)+0.5*u(k-2);%給出理想的辨識輸出采樣信號    
end
c0=[0.001 0.001 0.001 0.001]';%直接給出被辨識參數的初始值,即一個充分小的實向量
p0=10^6*eye(4,4);%直接給出初始狀態P0，即一個充分大的實數單位矩陣
E=0.000000005;%相對誤差E=0.000000005
c=[c0,zeros(4,14)];%被辨識參數矩陣的初始值及大小
e=zeros(4,15);%相對誤差的初始值及大小
for k=3:15; %開始求K 
    h1=[-z(k-1),-z(k-2),u(k-1),u(k-2)]'; x=h1'*p0*h1+1; x1=inv(x); %開始求K(k)
    k1=p0*h1*x1;%求出K的值
    d1=z(k)-h1'*c0; c1=c0+k1*d1;%求被辨識參數c
    e1=c1-c0;%求參數當前值與上一次的值的差值
    e2=e1./c0;%求參數的相對變化
    e(:,k)=e2; %把當前相對變化的列向量加入誤差矩陣的最后一列       
    c0=c1;%新獲得的參數作為下一次遞推的舊參數
    c(:,k)=c1;%把辨識參數c 列向量加入辨識參數矩陣的最后一列 
    p1=p0-k1*k1'*[h1'*p0*h1+1];%求出 p(k)的值
    p0=p1;%給下次用
    if e2<=E break;%若參數收斂滿足要求，終止計算
    end%小循環結束
end%大循環結束
c%顯示被辨識參數
e%顯示辨識結果的收斂情況
%分離參數
a1=c(1,:); a2=c(2,:); b1=c(3,:); b2=c(4,:); ea1=e(1,:); ea2=e(2,:); eb1=e(3,:); eb2=e(4,:);