logmapo.m
上传用户:look542
上传日期:2009-06-04
资源大小:784k
文件大小:3k
源码类别:

传真(Fax)编程

开发平台:

Matlab

  1. function L_all = logmapo(rec_s,g,L_a,ind_dec)
  2. % Copyright Nov 1998, Yufei Wu
  3. % MPRG lab, Virginia Tech.
  4. % for academic use only
  5. % Log_MAP algorithm using straightforward method to compute branch metrics
  6. % no approximation is used.
  7. % Can be simplified to Max-Log-MAP by using approximation ln(e^x+e^y) = max(x,y).
  8. % Input: rec_s: scaled received bits. 
  9. %               rec_s = 0.5 * L_c * yk = ( 2 * a * rate * Eb/N0 ) * yk
  10. %        g: code generator for the component RSC code, in binary matrix form.
  11. %        L_a: a priori info. for the current decoder, 
  12. %               scrambled version of extrinsic Inftyo. of the previous decoder.
  13. %        ind_dec: index of decoder. Either 1 or 2. 
  14. %               Encoder 1 is assumed to be terminated, while encoder 2 is open.
  15. %
  16. % Output: L_all: log-likelihood ratio of the symbols. Complete information.
  17. % Total number of bits: Inftyo. + tail
  18. L_total = length(rec_s)/2;
  19. [n,K] = size(g); 
  20. m = K - 1;
  21. nstates = 2^m;          % number of states in the trellis
  22. % Set up the trellis
  23. [next_out, next_state, last_out, last_state] = trellis(g);
  24. Infty = 1e10;
  25. % Initialization of Alpha
  26. Alpha(1,1) = 0; 
  27. Alpha(1,2:nstates) = -Infty*ones(1,nstates-1);
  28. % Initialization of Beta
  29. if ind_dec==1
  30.    Beta(L_total,1) = 0;
  31.    Beta(L_total,2:nstates) = -Infty*ones(1,nstates-1); 
  32. elseif ind_dec==2
  33.    Beta(L_total,1:nstates) = zeros(1,nstates);
  34. else
  35.    fprintf('ind_dec is limited to 1 and 2!n');
  36. end
  37. % Trace forward, compute Alpha
  38. for k = 2:L_total+1
  39.     for state2 = 1:nstates
  40.       gamma = -Infty*ones(1,nstates);
  41.       gamma(last_state(state2,1)) = (-rec_s(2*k-3)+rec_s(2*k-2)*last_out(state2,2))....
  42.            -log(1+exp(L_a(k-1)));
  43.       gamma(last_state(state2,2)) = (rec_s(2*k-3)+rec_s(2*k-2)*last_out(state2,4))....
  44.            +L_a(k-1)-log(1+exp(L_a(k-1)));
  45.       if(sum(exp(gamma+Alpha(k-1,:)))<1e-300)
  46.          Alpha(k,state2)=-Infty;
  47.       else
  48.          Alpha(k,state2) = log( sum( exp( gamma+Alpha(k-1,:) ) ) );  
  49.       end   
  50.     end
  51.     tempmax(k) = max(Alpha(k,:));
  52.     Alpha(k,:) = Alpha(k,:) - tempmax(k);
  53. end     
  54. % Trace backward, compute Beta
  55. for k = L_total-1:-1:1
  56.   for state1 = 1:nstates
  57.      gamma = -Infty*ones(1,nstates);
  58.      gamma(next_state(state1,1)) = (-rec_s(2*k+1)+rec_s(2*k+2)*next_out(state1,2))....
  59.            -log(1+exp(L_a(k+1)));
  60.      gamma(next_state(state1,2)) = (rec_s(2*k+1)+rec_s(2*k+2)*next_out(state1,4))....
  61.            +L_a(k+1)-log(1+exp(L_a(k+1)));
  62.      if(sum(exp(gamma+Beta(k+1,:)))<1e-300)
  63.         Beta(k,state1)=-Infty;
  64.      else
  65.         Beta(k,state1) = log(sum(exp(gamma+Beta(k+1,:))));
  66.      end   
  67.   end
  68.   tempmax_b(k) = max(Beta(k,:));
  69.   Beta(k,:) = Beta(k,:) - tempmax_b(k);
  70. end
  71. % Compute the soft output, log-likelihood ratio of symbols in the frame
  72. for k = 1:L_total
  73.   for state2 = 1:nstates
  74.      gamma0 = (-rec_s(2*k-1)+rec_s(2*k)*last_out(state2,2))....
  75.            -log(1+exp(L_a(k)));
  76.      gamma1 = (rec_s(2*k-1)+rec_s(2*k)*last_out(state2,4))...
  77.            +L_a(k)-log(1+exp(L_a(k)));
  78.      temp0(state2) = exp(gamma0 + Alpha(k,last_state(state2,1)) + Beta(k,state2));
  79.      temp1(state2) = exp(gamma1 + Alpha(k,last_state(state2,2)) + Beta(k,state2));
  80.   end
  81.   L_all(k) = log(sum(temp1)) - log(sum(temp0));
  82. end