// 5tap LMS algorithm simulation

clear;                                     // clear all variable

SIMLEN=300;                                // simulation length
TAPLEN=5;                                  // FIR filter length
MU=0.01;                                   // step size parameter
h=[1, 0, -3, -2, 3];                       // unkown system FIR filter coefficient

xvec=loadwave('male.wav');                 // male voice (fs=8kHz)
//xvec=loadwave('female.wav');             // female voice (fs=8kHz)
SIMLEN=length(xvec);                       // get simulation data length

evec=zeros(1:SIMLEN);
wvec1=zeros(1:SIMLEN); wvec2=zeros(1:SIMLEN); wvec3=zeros(1:SIMLEN);
wvec4=zeros(1:SIMLEN); wvec5=zeros(1:SIMLEN);
w=zeros(1:TAPLEN);
xbuf=zeros(1:TAPLEN);

for i=1:SIMLEN                             // simulation loop 
  x=xvec(i);                               // read input signal
  xbuf=[x, xbuf(1:TAPLEN-1)];              // shift register
  d=h*xbuf';                               // calculate desired signal
  y=w*xbuf';                               // adaptive FIR filter
  e=d-y;                                   // calculate error signal
  
  w=w+2.0*MU*e*xbuf;                       // LMS algorithm
  //w=w+2.0*MU*e*xbuf/(xbuf*xbuf'+1e-9);   // NLMS (Normalized LMS) algorithm
    
  wvec1(i)=w(1); wvec2(i)=w(2);  wvec3(i)=w(3);
  wvec4(i)=w(4); wvec5(i)=w(5);
  evec(i)=e;
end

scf(0); clf; plot(wvec1); plot(wvec2); plot(wvec3); plot(wvec4); plot(wvec5);
scf(1); clf; plot(evec);
scf(2); clf; plot(10*log10(evec+1e-9));

