olsblms.m

《自适应滤波算法与实现》（第二版）源码

%OLSBLMS Open-Loop Subband LMS algorithm
%
%   'ifile.mat' - input file containing:
%      Nr - members of ensemble
%      N - iterations
%      Sx - standard deviation of input
%      B - coefficient vector of plant (numerator)
%      A - coefficient vector of plant (denominator)
%      Sn - standard deviation of measurement noise
%      u - convergence factor
%      gamma - small constant to prevent the updating factor from getting too large
% 
%      a - small factor for the updating equation for the signal energy in the subbands
%      M - number of subbands
%      hk - analysis filterbank, with the filter coefficients in the rows
%      fk - synthesis filterbank, with the filter coefficients in the rows
%      Nw - order of the adaptive filter in the subbands
%      
%   'ofile.mat' - output file containing:
%      MSE - mean-square error
%      MSEsub - mean-square error in the subbands

clear all		% clear memory
load ifile;		% read input variables
L=M;                    % interpolation/decimation factor

for l=1:Nr
    disp(['Run: ' num2str(l)])
    nc=sqrt(Sn)*randn(1,N*M);     % noise at system output 
    xin=(rand(1,N*M)-.5);
    xin=sqrt(Sx)*sqrt(12)*xin;    % input signal
    d=filter(B,A,xin);            % desired signal
    d=d+nc;                       % unknown system output

    % Analysis of desired and input signals
    for k=1:M
	xaux=filter(hk(k,:),1,d);
	dsb(k,:)=xaux(find(mod((1:length(xaux))-1,L)==0)); % desired signal split in subbands
	xaux=filter(hk(k,:),1,xin);
	xsb(k,:)=xaux(find(mod((1:length(xaux))-1,L)==0)); % input signal split in subbands
    end;
  
    for m=1:M
	w_ol(m,:)=zeros(1,Nw+1);      % initial coefficient vector for
                                      % subband m
	x_ol(m,:)=[zeros(1,Nw+1)];    % initial input vector for subband m
	sig_ol(m)=0;
    end

    for k=1:N
      if mod(k,200)==0 disp(['Iteration: ' num2str(k)]), end;
      for m=1:M 
	  % Open-loop
	  % --------------------------------------------------
	  x_ol(m,:)=[xsb(m,k) x_ol(m,1:Nw)]; % new input vector for subband m
	  xaux=shiftdim(x_ol(m,:),1);
	  waux=shiftdim(w_ol(m,:),1);
	  e_ol(m,l,k)=dsb(m,k)-waux'*xaux;   % error at subband m
	  sig_ol(m)=(1-a)*sig_ol(m)+a*(xsb(m,k))^2;
	  unlms_ol=2*u/(gamma+(Nw+1)*sig_ol(m));
	  w_ol(m,:)=w_ol(m,:)+unlms_ol*e_ol(m,l,k)*shiftdim(xaux,-1); % new coefficient vector for subband m
      end;
    end;
    wk_ol(l,:,:)=w_ol;
end;

for m=1:M
    mse_ol(m,:)=mean(e_ol(m,:,:).^2,2);    % MSE at subband m
end;
MSE=mean(mse_ol,1);                 % overall MSE
MSEsub = mse_ol;

save ofile MSE MSEsub;