function [L_all,start_state]=logmap(rec,state_map,L_a_1row,iter,start_state,L1)%rec中包含A,B,parity1三部分,等長,因此網格圖中經歷的狀態轉移次數L_total=rec長度的1/3L_total = length(rec)/3;nstates = 8;          % number of states in the trellisL_all=zeros(1,2*L_total);% Set up the trellisstate_map = trellis;%無窮大Infty = 1e10;%前項概率alpha,三維分別對應4種輸入,網格圖中L_total+1個時刻,狀態總數Alpha=zeros(4,L_total+1,nstates);% Initialization of Alpha%前L1次等概初始化%L1+1次迭代開始,用L1次迭代估計的編碼器初始狀態start_state來初始化if iter<=L1    Alpha(:,1,:) = -log(nstates); else    Alpha(:,1,:)=-Infty;    Alpha(:,1,start_state+1)=0;end%后項概率beta,三維分別對應4種輸入,網格圖中L_total+1個時刻,狀態總數Beta=zeros(4,L_total,nstates);% Initialization of Betaif iter<=L1    Beta(:,L_total,:) = -log(nstates); else    Beta(:,L_total,:)=-Infty;    Beta(:,L_total,start_state+1)=0;end% L_a is seperated into two rows:各代表A,B的邊信息L_a=zeros(2,L_total);L_a(1,:)=L_a_1row(1:2:end);L_a(2,:)=L_a_1row(2:2:end);tempmax=zeros(1,L_total);% Trace forward, compute Alpha%計算狀態轉移概率gamma,對4種輸入分別計算,將平方計算簡化為乘積運算%rec(3*k-5),rec(3*k-4),rec(3*k-3)對應3比特組,前兩位是信息,末位是校驗.%gamma(1,i)=rec(3*k-5)*(-1)+rec(3*k-4)*(-1)+rec(3*k-3)*output(00,i)+log(1/(1+exp(L_a(1,k-1))))+log(1/(1+exp(L_a(2,k-1))))%output(00,i)是第i個狀態下輸入為00時的輸出,通過查表state_map得到for k = 2:L_total+1    gamma = -Infty*ones(4,nstates);    for state2 = 1:nstates        gamma(1,state2) = (-rec(3*k-5)-rec(3*k-4)+rec(3*k-3)*state_map(state2,11))...            -log(1+exp(L_a(1,k-1)))-log(1+exp(L_a(2,k-1)));        gamma(2,state2) = (-rec(3*k-5)+rec(3*k-4)+rec(3*k-3)*state_map(state2,14))...            -log(1+exp(L_a(1,k-1)))+L_a(2,k-1)-log(1+exp(L_a(2,k-1)));        gamma(3,state2) = (rec(3*k-5)-rec(3*k-4)+rec(3*k-3)*state_map(state2,17))...            +L_a(1,k-1)-log(1+exp(L_a(1,k-1)))-log(1+exp(L_a(2,k-1)));        gamma(4,state2) = (rec(3*k-5)+rec(3*k-4)+rec(3*k-3)*state_map(state2,20))...            +L_a(1,k-1)-log(1+exp(L_a(1,k-1)))+L_a(2,k-1)-log(1+exp(L_a(2,k-1)));    end%計算Alpha:Alpha(i,k,j)=sum(exp(Alpha(m,k-1,nextstate(i,j))+gamma(i,j)))--對m求和,m=0,1,2,3%nextstate(i,j)是j狀態下輸入為i時轉移到的狀態    for i=1:4        for j=1:nstates            Alpha(i,k,j)=exp(Alpha(1,k-1,state_map(j,i))+gamma(i,j))+exp(Alpha(2,k-1,state_map(j,i))+gamma(i,j))...                +exp(Alpha(3,k-1,state_map(j,i))+gamma(i,j))+exp(Alpha(4,k-1,state_map(j,i))+gamma(i,j));            if(Alpha(i,k,j)<1e-300)                Alpha(i,k,j)=-Infty;            else                Alpha(i,k,j) = log(Alpha(i,k,j));              end          end    end%歸一化    tempmax(k)=Alpha(1,k,1);    temp(1) = max(Alpha(1,k,:));    temp(2) = max(Alpha(2,k,:));    temp(3) = max(Alpha(3,k,:));    temp(4) = max(Alpha(4,k,:));    temp=sort(temp);    tempmax(k) = temp(4);    Alpha(1,k,:) = Alpha(1,k,:) - tempmax(k);    Alpha(2,k,:) = Alpha(2,k,:) - tempmax(k);    Alpha(3,k,:) = Alpha(3,k,:) - tempmax(k);    Alpha(4,k,:) = Alpha(4,k,:) - tempmax(k);end     % Trace backward, compute Beta%計算狀態轉移概率gamma,對4種輸入分別計算,將平方計算簡化為乘積運算%rec(3*k+1),rec(3*k+2),rec(3*k+3)對應3比特組,前兩位是信息,末位是校驗.%gamma(1,i)=rec(3*k+1)*(-1)+rec(3*k+2)*(-1)+rec(3*k+3)*output(00,i)+log(1/(1+exp(L_a(1,k-1))))+log(1/(1+exp(L_a(2,k-1))))%output(00,i)是第i個狀態下輸入為00時的輸出,通過查表state_map得到for k = L_total-1:-1:1    gamma = -Infty*ones(4,nstates);    for state1 = 1:nstates        gamma(1,state1) = (-rec(3*k+1)-rec(3*k+2)+rec(3*k+3)*state_map(state1,23))....            -log(1+exp(L_a(1,k+1)))-log(1+exp(L_a(2,k+1)));        gamma(2,state1) = (-rec(3*k+1)+rec(3*k+2)+rec(3*k+3)*state_map(state1,26))....            -log(1+exp(L_a(1,k+1)))+L_a(2,k+1)-log(1+exp(L_a(2,k+1)));        gamma(3,state1) = (rec(3*k+1)-rec(3*k+2)+rec(3*k+3)*state_map(state1,29))....            +L_a(1,k+1)-log(1+exp(L_a(1,k+1)))-log(1+exp(L_a(2,k+1)));        gamma(4,state1) = (rec(3*k+1)+rec(3*k+2)+rec(3*k+3)*state_map(state1,32))....            +L_a(1,k+1)-log(1+exp(L_a(1,k+1)))+L_a(2,k+1)-log(1+exp(L_a(2,k+1)));    end%計算Beta:Beta(i,k,j)=sum(exp(Beta(m,k+1,laststate(i,j))+gamma(i,j)))--對m求和,m=0,1,2,3%laststate(i,j)是輸入為i時轉移到j狀態的那個狀態    for i=1:4        for j=1:nstates            Beta(i,k,j)=exp(Beta(1,k+1,state_map(j,i+4))+gamma(i,j))+exp(Beta(2,k+1,state_map(j,i+4))+gamma(i,j))...                +exp(Beta(3,k+1,state_map(j,i+4))+gamma(i,j))+exp(Beta(4,k+1,state_map(j,i+4))+gamma(i,j));            if(Beta(i,k,j)<1e-300)                Beta(i,k,j)=-Infty;            else                Beta(i,k,j) = log(Beta(i,k,j));              end          end    end%歸一化    Beta(1,k,:) = Beta(1,k,:) - tempmax(k+1);    Beta(2,k,:) = Beta(2,k,:) - tempmax(k+1);    Beta(3,k,:) = Beta(3,k,:) - tempmax(k+1);    Beta(4,k,:) = Beta(4,k,:) - tempmax(k+1);end% Compute the soft output, log-likelihood ratio of symbols in the frame%//delta(1,i)~delta(4,i):分別代表第i個時刻下信息為00,01,10,11的概率for k = 1:L_total    for state2 = 1:nstates        gamma0 = (-rec(3*k-2)-rec(3*k-1)+rec(3*k)*state_map(state2,11))....            -log(1+exp(L_a(1,k)))-log(1+exp(L_a(2,k)));        gamma1 = (-rec(3*k-2)+rec(3*k-1)+rec(3*k)*state_map(state2,14))....            -log(1+exp(L_a(1,k)))+L_a(2,k)-log(1+exp(L_a(2,k)));        gamma2 = (rec(3*k-2)-rec(3*k-1)+rec(3*k)*state_map(state2,17))....            +L_a(1,k)-log(1+exp(L_a(1,k)))-log(1+exp(L_a(2,k)));        gamma3 = (rec(3*k-2)+rec(3*k-1)+rec(3*k)*state_map(state2,20))....            +L_a(1,k)-log(1+exp(L_a(1,k)))+L_a(2,k)-log(1+exp(L_a(2,k)));  %temp=delta*Alpha(k-1)*Beta(k)           temp0(state2) = sum(exp(gamma0 + Alpha(:,k,state_map(state2,1)) + Beta(1,k,state2)))...            +sum(exp(gamma0 + Alpha(:,k,state_map(state2,1)) + Beta(2,k,state2)))...            +sum(exp(gamma0 + Alpha(:,k,state_map(state2,1)) + Beta(3,k,state2)))...            +sum(exp(gamma0 + Alpha(:,k,state_map(state2,1)) + Beta(4,k,state2)));        temp1(state2) = sum(exp(gamma1 + Alpha(:,k,state_map(state2,2)) + Beta(2,k,state2)))...            +sum(exp(gamma1 + Alpha(:,k,state_map(state2,2)) + Beta(1,k,state2)))...            +sum(exp(gamma1 + Alpha(:,k,state_map(state2,2)) + Beta(3,k,state2)))...            +sum(exp(gamma1 + Alpha(:,k,state_map(state2,2)) + Beta(4,k,state2)));        temp2(state2) = sum(exp(gamma2 + Alpha(:,k,state_map(state2,3)) + Beta(3,k,state2)))...            +sum(exp(gamma2 + Alpha(:,k,state_map(state2,3)) + Beta(1,k,state2)))...            +sum(exp(gamma2 + Alpha(:,k,state_map(state2,3)) + Beta(2,k,state2)))...            +sum(exp(gamma2 + Alpha(:,k,state_map(state2,3)) + Beta(4,k,state2)));        temp3(state2) = sum(exp(gamma3 + Alpha(:,k,state_map(state2,4)) + Beta(4,k,state2)))...            +sum(exp(gamma3 + Alpha(:,k,state_map(state2,4)) + Beta(1,k,state2)))...            +sum(exp(gamma3 + Alpha(:,k,state_map(state2,4)) + Beta(2,k,state2)))...            +sum(exp(gamma3 + Alpha(:,k,state_map(state2,4)) + Beta(3,k,state2)));    end    L_all(2*k-1)=log(sum(temp2)+sum(temp3))-log(sum(temp0)+sum(temp1));    L_all(2*k)  =log(sum(temp1)+sum(temp3))-log(sum(temp0)+sum(temp2));end%狀態估計,用Alpha的N時刻狀態概率和Beta的0時刻狀態概率聯合估計if iter==L1        % determine the start and end state,and return as start_state    p_endstate=zeros(1,nstates);    p_startstate=zeros(1,nstates);    % p_endstate%用Alpha的最后一個時刻的狀態估計編碼器的終止狀態,計算每個狀態概率    p2=zeros(1,nstates);    for i=1:nstates        p2(i)=sum(exp(Alpha(:,L_total+1,i)));    end    p1=sum(p2);    p_endstate=log(p2)-log(p1);    % p_startstate%用Beta的第0個時刻的狀態估計編碼器的起始狀態,計算每個狀態概率    p2=zeros(1,nstates);    beta=zeros(4,nstates);    delta = -Infty*ones(4,nstates);    for state1 = 1:nstates        delta(1,state1) = (-rec(1)-rec(2)+rec(3)*state_map(state1,23))-log(1+exp(L_a(1,1)))-log(1+exp(L_a(2,1)));        delta(2,state1) = (-rec(1)+rec(2)+rec(3)*state_map(state1,26))-log(1+exp(L_a(1,1)))+L_a(2,1)-log(1+exp(L_a(2,1)));        delta(3,state1) = (rec(1)-rec(2)+rec(3)*state_map(state1,29))+L_a(1,1)-log(1+exp(L_a(1,1)))-log(1+exp(L_a(2,1)));        delta(4,state1) = (rec(1)+rec(2)+rec(3)*state_map(state1,32))+L_a(1,1)-log(1+exp(L_a(1,1)))+L_a(2,1)-log(1+exp(L_a(2,1)));    end    for i=1:4        for j=1:nstates            beta(i,j)=exp(Beta(1,1,state_map(j,i+4))+delta(i,j))+exp(Beta(2,1,state_map(j,i+4))+delta(i,j))...                +exp(Beta(3,1,state_map(j,i+4))+delta(i,j))+exp(Beta(4,1,state_map(j,i+4))+delta(i,j));            if(beta(i,j)<1e-300)                beta(i,j)=-Infty;            else                beta(i,j) = log(beta(i,j));              end          end    end    p2=zeros(1,nstates);    for i=1:nstates        p2(i)=sum(exp(beta(:,i)));    end    p1=sum(p2);    p_startstate=log(p2)-log(p1);    % Sc%用Alpha和Beta聯合估計編碼器的起始和終止狀態,選出概率最大的一個    p=zeros(1,nstates);    start_state=0;    for i=1:nstates        p(i)=p_endstate(i)+p_startstate(i);        if p(i)>p(start_state+1)            start_state=i-1;        end    endend