function [W,e,Xout] = NLMS(x,d,N,mu,alpha,Xin,Winit);
% [W,e,Xout] = NLMS(x,d,N,mu,alpha,Xin,Winit);
%
% Implementation of the normalized least mean square algorithm.
% For complex inputs, this algorithms return the complex conjugated 
% filter coefficients.
%
% Inputs:
% 	x	reference signal
%       d	desired signal
% 	N	# of taps in filter
%	mu	step size (0...1)
%	alpha   power bias term
% 	Xin	initialization of state vector
%	Winit	initialization of filter weights
%
% Outputs:
%	W	adapted filter weights (column vector)
%	e	error signal (row vector)
%	Xout	last state vector (column vector)
%
% St.Weiss, SIPI-USC, 16.2.1997

% ***** check input parameters for correctness *****
if (~exist('Winit')),
   W = zeros(N,1);
else
   if length(Winit)~=N,
      error('weight initialization must match filter length');
   end;
   W = Winit;
endif;

Lx = length(x);
[m,n] = size(x);
if n>m,
   x = x.';
endif;

if (~exist('Xin')),
   x = [zeros(N-1,1); x];
   p = alpha;
else
   if length(Xin)~=(N-1) 
      error('state initialization must match filter length minus one');
   endif;
   x = [Xin; x];
   p = Xin'*Xin+alpha;
endif; 
END

% ***** perform adaptive filtering *****
for k = 1:Lx,
   X = x(k+N-1:-1:k);
   y = W'*X;
   e(k) = d(k)-y;
   p = p+x(k+N-1)*conj(x(k+N-1));
   W = W+(2*mu*conj(e(k))/p)*X;
   p = p-x(k)*conj(x(k));
endfor;

% ***** prepare output parameter *****
Xout = x(Lx+1:Lx+N-1);