function [w,h,WM,E,mu,M,J] = lms_fir(muv);

% Inlämningsuppgift - LMS 
%
% Kurs: Projektkurs i adaptiv signalbehandling för F4
%
% 20020315 Jonas Rutström
%  
%
%  Mål  : Ökad förståelse för LMS-algoritmen. 
%  
%  Del 1: Teori
%  
%         Läsa kap.6 i Widrow.
%
%  Del 2: Matlab
%  
%  
%  -------- Exempel på implementering av LMS-algoritmen -----------------
% 

  close all
    
  nfir = 2;
  npts = 200;
  WM   = zeros(npts,nfir);
  E    = zeros(length(muv),npts);
  
% Färgande filter för insignalen.
% Graden av färgning betäms genom brytfrekvensen fc.
% fc kan ligga i intervallet [0-1]. Detta filter kommer att bestämma utseendet på
% autokorrelationsmatrisen R. 

  fc      = 1;
  [Bf,Af] = butter(2,fc);

% Insignal
  x = randn(npts,1);      % Vitt brus 
  x = filter(Bf,Af,x);    % Färgat brus ( Vitt -> Bf/Af -> Färgat )
  
  R      = xcorrm(x,2)    % Autokorrelationsmatrisen för insignalen
  R_inv  = inv(R)
  lam = mean(eig(R));     
 
% Överföringsfunktion - FIR-filter
  
  H = [0.4 0.5];          % Nollställen

% Önskad signal
    
  d = filter(H,1,x);   

% LMS algoritm

% Yttre loop för att beräkna mse för olika värden på mu

  for j = 1:length(muv);
     
     mu = muv(j);           % plocka ut aktuellt mu.
     w = zeros(nfir,1);     % startvärden för vikterna
     y = zeros(npts,1);     % startvärden för utsignalen
     
    % ---- LMS loopen ----------
    for i = nfir:length(x)-1;
       
      WM(i,:) = w';        % Spara undan vikerna för varje iteration i en matris.
       
      y(i)    = w'*x(i:-1:i-nfir+1);  
      e       = d(i) - y(i);
      
      if e^2 > 1
        E(j,i) = 1;
      else
        E(j,i) = (e^2);
      end
      
      mu = mu;
      
      w    = w + 2*mu*e*x(i:-1:i-nfir+1);

    end

    % ---- slut LMS loop ----------
       
 end
 % Plottar
 subplot(211), plot([y d]); xlabel('iterationer');
 subplot(212), bar(w);      xlabel('vikter');
 
 figure;
 
 subplot(211), plot(WM);
 subplot(212), 
 
  % Här beräknas felytan
    w1 = [-0.5:.01:1.2]';
    w2 = [-0.5:.01:1.2]';

    [W1 W2]  = meshgrid(w1,w2);
    [W10 W20]= meshgrid(w1-H(1),w2-H(2));
    
    J = R(1,1)*W1.*W1 + (2*R(1,2))*W2.*W1 + R(2,2)*W2.*W2;
    J = -(J - 2 * (R(1,1)*W1.*W10 + R(1,2)*W1.*W20 + R(2,1)*W2.*W10 + R(2,2)*W2.*W20));
    
 contour(w1,w2,J,20);
 hold
 plot(WM(:,1),WM(:,2),'-');
 plot(H(1),H(2),'rx');

 %mesh([1:npts],muv,E),xlabel('iterationer');ylabel('mu');



  
  
  
  
  