www.pudn.com > ln.rar > mainp1test4.asv

 
clear;clc;close all; 
 
display('This program is for the paper titled:A class of novel blind source extraction algorithms based on a linear predictor') 
display('for the whitened case and updated by the eqn.(21)') 
load ABio7 
sig=benchcicbar; 
s0=sig(2,:);s1=sig(3,:); 
s2=sig(4,:);s3=sig(5,:); 
%%%%%%%%%%------beginning of fiuring the source signal-----%%%%%%%%% 
figure 
title('The source signals') 
subplot(4,1,1);plot(s0);ylabel('s_1'); 
subplot(4,1,2);plot(s1);ylabel('s_2'); 
subplot(4,1,3);plot(s2);ylabel('s_3'); 
subplot(4,1,4);plot(s3);ylabel('s_4'); 
 
%%%%%%%%%%------end of fiuring the source signal-----%%%%%%%%% 
A=[-0.9370 0.9448 0.3651 -0.2617; 
    -0.8263 -0.3478 -0.2872 0.6279; 
    0.3730 0.1975 0.9953 0.8325; 
    0.8955 0.3578 0.5649 0.6257]; 
s=[s0;s1;s2;s3]; 
x=A*s; 
figure 
title('The mixed signals') 
subplot(4,1,1);plot(x(1,:));ylabel('x_1'); 
subplot(4,1,2);plot(x(2,:));ylabel('x_2'); 
subplot(4,1,3);plot(x(3,:));ylabel('x_3'); 
subplot(4,1,4);plot(x(4,:));ylabel('x_4'); 
 
 
%%%%%%%%  whiten the mixed signal 
covarianceMatrix=cov(x',1);     % 产生协方差阵 
[E,D]=eig(covarianceMatrix);           % 求出协方差的特征向量和特征值 
 
[x1,WhitenMatrix]=whitenSig(x,E,D);     % 对信号进行白化 
figure 
title('The whitened signals') 
subplot(4,1,1);plot(x1(1,:));ylabel('x1_1'); 
subplot(4,1,2);plot(x1(2,:));ylabel('x1_2'); 
subplot(4,1,3);plot(x1(3,:));ylabel('x1_3'); 
subplot(4,1,4);plot(x1(4,:));ylabel('x1_4'); 
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 
[m,n]=size(x1); 
P=30;%%%%%-----length of the linear predicor 
%b=rand(1,P)'-0.5;%%%------------linear predictor coefficient vector---------- 
 
%%%%%%%%%%%%%%--------------update the weight vector--------------------- 
xx1=x1;%%%%%xx1 is the copy of x1 
SS=[]; 
mu=0.0025;%%%%%%%%%%%%-------learn stepsize--- 
w=rand(m,1); 
wold=w; 
 
%for sn=1:m 
for sn=1:1 
   w=rand(m,1); 
   wold=w; 
   %P=P-(sn-1)*5; 
   b=rand(1,P)'; 
   E=[]; 
    
   for n=P+1:n 
    y(n)=w'*x1(:,n); 
    yv=0; 
    for p=1:P 
        yv=yv+b(p)*y(n-p); 
     end 
    e=y(n)-yv; 
    w=w-mu*e*x1(:,n); 
    w=w/norm(w,'fro'); 
    wold=w; 
   end 
   y=w'*x1; 
   SS=[SS;y]; 
    
   P=5; 
    
    Q=10*eye(P); 
     for i=P+1:n 
       y=w'*x1; 
       Y=[]; 
       X=[]; 
       for j=1:P 
          Y=[Y;y(i-j)]; 
          X=[X,x1(:,i-j)]; 
       end 
       K=Q*Y*inv(Y'*Q*Y+1); 
       v=y(i)-Y'*b; 
       b=b+K*v; 
       Q=Q-K*Y'*Q; 
       hatx=zeros(m,1); 
       for j=1:P 
          hatx=hatx+b(j)*X(:,j); 
       end 
       hatx=x1(:,i)-hatx;      
       e=y(i)-b'*Y; 
       E=[E,e]; 
       w=wold-mu*e*hatx;   
       w=w/norm(w); 
       wold=w; 
   end 
   
   y=w'*x1; 
   SS=[SS;y]; 
   %hatw=rand(m,1); 
   %error=1; 
   %%%update hatw 
 
  %Eyx=x1*diag(y); 
  %Eyx=sum(Eyx,2); 
  %x1=x1-Eyx/norm(y)^2*y; 
    
   % while(error>0.01) 
       %x2=x1-hatw*y; 
       %hatw2=hatw; 
       %hatw=hatw+0.0001*x2*y'; 
       %error=norm(hatw2-hatw,'fro'); 
   %end 
   %x1=x2; 
end 
figure 
title('The extracted signals') 
subplot(2,1,1);plot(SS(1,:));ylabel('(a)'); 
subplot(2,1,2);plot(SS(2,:));ylabel('(b)'); 
%subplot(4,1,1);plot(SS(1,:));ylabel('s_1'); 
%subplot(4,1,2);plot(SS(2,:));ylabel('s_2'); 
%subplot(4,1,3);plot(SS(3,:));ylabel('s_3'); 
%subplot(4,1,4);plot(SS(4,:));ylabel('s_4');