function L_all = logmapo(rec_s,g,L_a,ind_dec)% Copyright Nov 1998, Yufei Wu% MPRG lab, Virginia Tech.% for academic use only% Log_MAP algorithm using straightforward method to compute branch metrics% no approximation is used.% Can be simplified to Max-Log-MAP by using approximation ln(e^x+e^y) = max(x,y).% Input: rec_s: scaled received bits. %               rec_s = 0.5 * L_c * yk = ( 2 * a * rate * Eb/N0 ) * yk%        g: code generator for the component RSC code, in binary matrix form.%        L_a: a priori info. for the current decoder, %               scrambled version of extrinsic Inftyo. of the previous decoder.%        ind_dec: index of decoder. Either 1 or 2. %               Encoder 1 is assumed to be terminated, while encoder 2 is open.%% Output: L_all: log-likelihood ratio of the symbols. Complete information.% Total number of bits: Inftyo. + tailL_total = length(rec_s)/2;[n,K] = size(g); m = K - 1;nstates = 2^m;          % number of states in the trellis% Set up the trellis[next_out, next_state, last_out, last_state] = trellis(g);for state=1:nstates    temp2=last_out(state,2:2:4);  for bit=0:1    b2_last_state(state,bit+1)=last_state(state,find(temp2==2*bit-1));  end   if (last_out(state,2)==-1)       b2_last_out(state,:)=last_out(state,:);   else       b2_last_out(state,1:2)=last_out(state,3:4);       b2_last_out(state,3:4)=last_out(state,1:2);   endendInfty = 1e100;% Initialization of AlphaAlpha(1,1) = 0; Alpha(1,2:nstates) = -Infty*ones(1,nstates-1);% Initialization of Betaif ind_dec==1   Beta(L_total,1) = 0;   Beta(L_total,2:nstates) = -Infty*ones(1,nstates-1); elseif ind_dec==2   Beta(L_total,1:nstates) = zeros(1,nstates);else   fprintf('ind_dec is limited to 1 and 2!\n');end% Trace forward, compute Alphafor k = 2:L_total+1    for state2 = 1:nstates      gamma = -Infty*ones(1,nstates);      gamma(last_state(state2,1)) = (-rec_s(2*k-3)+rec_s(2*k-2)*last_out(state2,2));%-log(1+exp(L_a(2*k-3)));      gamma(last_state(state2,2)) = (rec_s(2*k-3)+rec_s(2*k-2)*last_out(state2,4));% +L_a(2*k-3)-log(1+exp(L_a(2*k-3)));      if(sum(exp(gamma+Alpha(k-1,:)))<1e-300)         Alpha(k,state2)=-Infty;      else         Alpha(k,state2) = log( sum( exp( gamma+Alpha(k-1,:) ) ) );        end       end    tempmax(k) = log(sum(exp(Alpha(k,:))));    Alpha(k,:) = Alpha(k,:) - tempmax(k);end     % Trace backward, compute Betafor k = L_total-1:-1:1  for state1 = 1:nstates     gamma = -Infty*ones(1,nstates);     gamma(next_state(state1,1)) = (-rec_s(2*k+1)+rec_s(2*k+2)*next_out(state1,2));%-log(1+exp(L_a(2*k+1)));     gamma(next_state(state1,2)) = (rec_s(2*k+1)+rec_s(2*k+2)*next_out(state1,4));% +L_a(2*k+1)-log(1+exp(L_a(2*k+1)));     if(sum(exp(gamma+Beta(k+1,:)))<1e-300)        Beta(k,state1)=-Infty;     else        Beta(k,state1) = log(sum(exp(gamma+Beta(k+1,:))));     end     end  Beta(k,:) = Beta(k,:) - tempmax(k+1);end% Compute the soft output, log-likelihood ratio of symbols in the framefor k = 1:L_total  for state2 = 1:nstates     gamma0 = (-rec_s(2*k-1)+rec_s(2*k)*last_out(state2,2));%-log(1+exp(L_a(2*k-1)));     gamma1 = (rec_s(2*k-1)+rec_s(2*k)*last_out(state2,4)) ;%+L_a(2*k-1)-log(1+exp(L_a(2*k-1)));     temp0(state2) = exp(gamma0 + Alpha(k,last_state(state2,1)) + Beta(k,state2));     temp1(state2) = exp(gamma1 + Alpha(k,last_state(state2,2)) + Beta(k,state2));          gamma00 = (rec_s(2*k-1)*b2_last_out(state2,1)+rec_s(2*k)*b2_last_out(state2,2)) ;%-log(1+exp(L_a(2*k)));     gamma11 = (rec_s(2*k-1)*b2_last_out(state2,3)+rec_s(2*k)*b2_last_out(state2,4)) ;%+L_a(2*k)-log(1+exp(L_a(2*k)));          temp00(state2)= exp(gamma00 + Alpha(k,b2_last_state(state2,1)) + Beta(k,state2));     temp01(state2)= exp(gamma11 + Alpha(k,b2_last_state(state2,2)) + Beta(k,state2));  end  L_all(2*k-1) = log(sum(temp1)) - log(sum(temp0));  L_all(2*k) = log(sum(temp01)) - log(sum(temp00));end