%自適應濾波器

clear;

v=1;               %信號方差
l=100;             %樣本個數
u=0.1;            %步長
hi=-0.8;           %hi
h=zeros(l,1);
h(1)=0;            %h(0)
mp=0;              %平均超調量
ct=1;              %統計次數

for count=1:ct
    x=randn(1,l);
    meanx=sum(x)/l;
    covx=sum((x-meanx)*(x-meanx)')/l;
    x=0.5*(x-meanx)/sqrt(covx);
    w=randn(1,l);
    meanw=sum(w)/l;
    covw=sum((w-meanw)*(w-meanw)')/l;
    w=sqrt(v)*(w-meanw)/sqrt(covw);
    
    s=hi*x;
    y=s+w;
    
    t=0:1:l-1;
    %plot(t,s,t,y);
    
    cstep=l;   %converge step
    r=x*x'/l;
    for i=2:l
        h(i)=hi+((1-2*u*r)^(i-1))*(h(1)-hi);     %理論
        if abs((h(i)-hi)/hi)<=0.05 cstep=i;break,end %收斂速度
    end
    
    for i=cstep+1:l
        h(i)=hi+((1-2*u*r)^(i-1))*(h(1)-hi);    %理論
    end
    
    hm=zeros(1,l);
    sm=zeros(1,l);
    hm(1)=h(1);    
    em=0;
    
    for i=1:l-1
        sm(i)=hm(i)*x(i);
        hm(i+1)=hm(i)+2*u*(y(i)-hm(i)*x(i))*x(i);
        if i>=cstep em=em+(y(i)-sm(i))^2;end
    end
    
    sm(l)=hm(l)*x(l);
    em=em/(l-1-cstep);
    m=(em-v)/v;
    mp=mp+m;
end
mp=mp/ct


plot(t+1,h,':b',t+1,hm,'r');
legend('理論值','估計值',0);
title('h(n)的理論值和估計值');
xlabel('n');

figure;
plot(t+1,s,':b',t+1,sm,'r');
legend('s(n)','估計值',0);
title('s(n)及其估計值');
xlabel('n');