% Demonstration of LMS algorithm for noise cancellation.
% Rich Kozick, Spring 1997

% Desired signal
s = auread('/home/pleiades/ACM/snd/clint_eastwood.au');
Ls = length(s);

% Interference + random noise
n = auread('/home/pleiades/ACM/snd/spacemusic.au');
Sign = 0.1;
n = n(1:Ls) + Sign*randn(Ls,1);

y = s + n;

Sigx = 0.01;
bx = [1];	% bx and ax are filtering on n to produce x
ax = [1 0.9];
Dx = 10;	% Delay of n that appears in x
x = filter(bx, ax, n);
x = [x(Dx:(Ls)); zeros(Dx-1,1)] + Sigx*randn(Ls,1);

N = 20;		% Length of adaptive filter

% LMS algorithm for adaptive noise cancellation

w = zeros(N,1);
mu = 1/(10*N*var(x));
for k=N:Ls
  xk = x(k:-1:(k-N+1));
  nhat(k) = w'*xk;
  e(k) = y(k) - nhat(k);
  w = w + 2*mu*e(k)*xk;
end

% The signal estimate is in the vector e