LM Algoritması ile FIR Filtre Parametre Kestirimi

Posted in Optimizasyon, Yazılar by On Temmuz 5, 2015. No comments

Aşağıda FIR filtre parametre kestirim işlemlerini gerçekleştiren Matlab kodları verilmiştir.

  1. clc, clear all, close all;
  2.  
  3. g_SAMPLE = 1000;
  4. g_ITERATION = 2;
  5. g_PARAMETER = 20;
  6.  
  7. f_1 = 10;
  8. f_2 = 500;
  9. f_3 = 1000;
  10. f_s = 5000;
  11. n_s = 1 : g_SAMPLE;
  12. INPUT = cos(2 * pi * n_s * f_1 / f_s) + cos(2 * pi * n_s * f_2 / f_s) + ...
  13.         cos(2 * pi * n_s * f_3 / f_s);
  14. figure, plot(INPUT);
  15. [n_a n_b] = fir1(g_PARAMETER - 1, 200 / f_s, 'low');
  16. DESIRED_OUTPUT = filter_conv(INPUT, n_a)
  17. figure, plot(DESIRED_OUTPUT);
  18.  
  19. lamda = 0.01;
  20.  
  21. p_a = rand(1, g_PARAMETER);
  22.  
  23. for n_i = 1 : g_ITERATION
  24.     E(n_i) = 0;
  25.    
  26.     for n_j = 1 : g_SAMPLE - g_PARAMETER + 1;
  27.         OUTPUT_TRAIN(n_j) = filter_conv(INPUT(n_j : n_j + g_PARAMETER - 1), p_a);
  28.         ERROR(n_j) = DESIRED_OUTPUT(n_j) - OUTPUT_TRAIN(n_j);
  29.         E(n_i) = E(n_i) + 0.5 * ERROR(n_j) ^ 2;
  30.         JACOBIAN(n_j, : ) = INPUT(n_j : n_j + g_PARAMETER - 1);
  31.     end
  32.  
  33.     kontrol = 1;
  34.  
  35.     while(kontrol == 1)
  36.        
  37.          J_new = inv(JACOBIAN' * JACOBIAN + lamda * eye(g_PARAMETER,g_PARAMETER)) * JACOBIAN' * ERROR';
  38.          p_a_new = p_a + J_new';
  39.  
  40.          E_NEW(n_i) = 0;
  41.  
  42.          for n_j = 1 : g_SAMPLE - g_PARAMETER + 1;
  43.             OUTPUT_TRAIN(n_j) = filter_conv(INPUT(n_j : n_j + g_PARAMETER - 1), p_a_new);
  44.             ERROR_NEW(n_j) = DESIRED_OUTPUT(n_j) - OUTPUT_TRAIN(n_j);
  45.             E_NEW(n_i) = E_NEW(n_i) + 0.5 * ERROR_NEW(n_j) ^ 2;
  46.            
  47.         end
  48.        
  49.        if E_NEW(n_i) < E(n_i)
  50.           p_a = p_a_new;
  51.           lamda = lamda * .5;
  52.            kontrol = 0;
  53.        else
  54.             lamda = lamda / .5;
  55.        end;
  56.      end  
  57. end;
  58.  
  59. TEST_OUTPUT = filter_conv(INPUT, p_a)
  60. figure, plot(TEST_OUTPUT);
  1. function FILTER_OUT = filter_conv(INPUT, n_a)
  2.  
  3. for n_i = 1 : length(INPUT) - length(n_a) + 1
  4.     FILTER_OUT(n_i) = sum(INPUT(n_i : n_i + length(n_a) - 1).* n_a);
  5.    
  6. end