siso_map_decoder_1.m

turbo 码matlab编译码程序

function [decision_c,extrin_information]=SISO_MAP_decoder_1(x,x1,CONS,prior_information)
% map algorithm,zero force RSC,in the first interation ,extrinc is zero sequence from decoder 1 to
% decoder 2
% dummy length is Npre
K=length(x);
m=3;
Npre=7*m;
for t=1:K
    for l=1:8                       % l is the trellis numbers
           p_exp(t)=exp(prior_information(t));
           p_1(t)=p_exp(t)/(1+p_exp(t));
           p_0(t)=1/(1+p_exp(t));
           
           [next_state_1,x_p]=nxt_stat1513(l-1,1);
           R_r_1(t,next_state_1+1)=p_1(t)*exp(CONS*((x(t)-1)^2+(x1(t)-(2*x_p-1))^2));
           [next_state_0,x_p]=nxt_stat1513(l-1,0);
           R_r_0(t,next_state_0+1)=p_0(t)*exp(CONS*((x(t)+1)^2+(x1(t)-(2*x_p-1))^2));
                
    end   
end
% SISO 1  
A_a_initial(1)=1;
A_a_initial(2:8)=0;
    
for t=1:K/4     %  虚拟长度为一个子帧的长度t=1:K/2 
    for ll=1:8
            A_a(t,ll)=0;
    end
    A_a_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);
           
            if t==1
                A_a(t,next_state_0+1)=A_a(t,next_state_0+1)+A_a_initial(ll)*R_r_0(t,next_state_0+1);
                A_a(t,next_state_1+1)=A_a(t,next_state_1+1)+A_a_initial(ll)*R_r_1(t,next_state_1+1);
            else
                A_a(t,next_state_0+1)=A_a(t,next_state_0+1)+A_a(t-1,ll)*R_r_0(t,next_state_0+1);
                A_a(t,next_state_1+1)=A_a(t,next_state_1+1)+A_a(t-1,ll)*R_r_1(t,next_state_1+1);

            end
             
    end 
    for ll=1:8
       A_a_totall(t)= A_a_totall(t)+A_a(t,ll);
    end
    for ll=1:8
        A_a(t,ll)=A_a(t,ll)/A_a_totall(t);
    end
end
% SISO 2 前向递推
dummy_A_a(K/4-Npre,:)=1/8;
for t=K/4-Npre+1:K/2
    for ll=1:8
           dummy_A_a(t,ll)=0;
    end
    dummy_A_a_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);
           
            
           
                dummy_A_a(t,next_state_0+1)=dummy_A_a(t,next_state_0+1)+dummy_A_a(t-1,ll)*R_r_0(t,next_state_0+1);
                dummy_A_a(t,next_state_1+1)=dummy_A_a(t,next_state_1+1)+dummy_A_a(t-1,ll)*R_r_1(t,next_state_1+1);

  
    end 
    for ll=1:8
       dummy_A_a_totall(t)= dummy_A_a_totall(t)+dummy_A_a(t,ll);
    end
    for ll=1:8
        dummy_A_a(t,ll)=dummy_A_a(t,ll)/dummy_A_a_totall(t);
    end
end
for t=K/4+1:K/2
    for ll=1:8
        A_a(t,ll)=dummy_A_a(t,ll);
    end
    A_a_totall(t)=dummy_A_a_totall(t);
end
% SISO 3 前向递推
dummy_A_a(K/2-Npre,:) =1/8;
for t=K/2-Npre+1:3*K/4
    for ll=1:8
           dummy_A_a(t,ll)=0;
    end
    dummy_A_a_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);
           

                dummy_A_a(t,next_state_0+1)=dummy_A_a(t,next_state_0+1)+dummy_A_a(t-1,ll)*R_r_0(t,next_state_0+1);
                dummy_A_a(t,next_state_1+1)=dummy_A_a(t,next_state_1+1)+dummy_A_a(t-1,ll)*R_r_1(t,next_state_1+1);

    end 
    for ll=1:8
       dummy_A_a_totall(t)= dummy_A_a_totall(t)+dummy_A_a(t,ll);
    end
    for ll=1:8
        dummy_A_a(t,ll)=dummy_A_a(t,ll)/dummy_A_a_totall(t);
    end
end
for t=K/2+1:3*K/4
    for ll=1:8
        A_a(t,ll)=dummy_A_a(t,ll);
    end
    A_a_totall(t)=dummy_A_a_totall(t);
end
% SISO 4 前向递推
dummy_A_a(3*K/4-Npre,:) =1/8;
for t=3*K/4-Npre+1:K
    for ll=1:8
           dummy_A_a(t,ll)=0;
    end
    dummy_A_a_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);
           
      
                dummy_A_a(t,next_state_0+1)=dummy_A_a(t,next_state_0+1)+dummy_A_a(t-1,ll)*R_r_0(t,next_state_0+1);
                dummy_A_a(t,next_state_1+1)=dummy_A_a(t,next_state_1+1)+dummy_A_a(t-1,ll)*R_r_1(t,next_state_1+1);

    end 
    for ll=1:8
       dummy_A_a_totall(t)= dummy_A_a_totall(t)+dummy_A_a(t,ll);
    end
    for ll=1:8
        dummy_A_a(t,ll)=dummy_A_a(t,ll)/dummy_A_a_totall(t);
    end
end
for t=3*K/4+1:K
    for ll=1:8
        A_a(t,ll)=dummy_A_a(t,ll);
    end
    A_a_totall(t)=dummy_A_a_totall(t);
end

% SISO 4  后向递推
% for l=1:8
%     B_b_K(1:8)=A_a(K,l);
% end
B_b_K=[1,0,0,0,0,0,0,0];    
for t=K-1:-1:3*K/4+1
    for ll=1:8
            B_b(t,ll)=0;
    end
    B_b_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);
           
            if t==K-1
                B_b(t,ll)=B_b_K(next_state_0+1)*R_r_0(t+1,next_state_0+1)+B_b_K(next_state_1+1)*R_r_1(t+1,next_state_1+1);
            else
                B_b(t,ll)=B_b(t+1,next_state_0+1)*R_r_0(t+1,next_state_0+1)+B_b(t+1,next_state_1+1)*R_r_1(t+1,next_state_1+1);
            end
             
    end 
    for ll=1:8
        B_b_totall(t)=A_a_totall(t+1);
    end
    for ll=1:8
        B_b(t,ll)=B_b(t,ll)/B_b_totall(t);
    end
end

% SISO 3 后向递推
dummy_B_b(3*K/4+Npre,:)=A_a(3*K/4+Npre,:);
for t=3*K/4+Npre-1:-1:K/2+1
    for ll=1:8
            dummy_B_b(t,ll)=0;
    end
    dummy_B_b_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);

            dummy_B_b(t,ll)=dummy_B_b(t+1,next_state_0+1)*R_r_0(t+1,next_state_0+1)+dummy_B_b(t+1,next_state_1+1)*R_r_1(t+1,next_state_1+1);
             
    end 
    for ll=1:8
        dummy_B_b_totall(t)=A_a_totall(t+1);
    end
    for ll=1:8
        dummy_B_b(t,ll)=dummy_B_b(t,ll)/dummy_B_b_totall(t);
    end
end
for t=3*K/4:-1:K/2+1
    for ll=1:8
        B_b(t,ll)=dummy_B_b(t,ll);
    end
end
% SISO 2后向递推
dummy_B_b(K/2+Npre,:)=A_a(K/2+Npre,:);
    
for t=K/2+Npre-1:-1:K/4+1
    for ll=1:8
            dummy_B_b(t,ll)=0;
    end
    dummy_B_b_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);
           
            dummy_B_b(t,ll)=dummy_B_b(t+1,next_state_0+1)*R_r_0(t+1,next_state_0+1)+dummy_B_b(t+1,next_state_1+1)*R_r_1(t+1,next_state_1+1);
      
    end 
    for ll=1:8
        dummy_B_b_totall(t)=A_a_totall(t+1);
    end
    for ll=1:8
        dummy_B_b(t,ll)=dummy_B_b(t,ll)/dummy_B_b_totall(t);
    end
end
for t=K/2:-1:K/4+1
    for ll=1:8
        B_b(t,ll)=dummy_B_b(t,ll);
    end
end

% SISO 1后向递推
dummy_B_b(K/4+Npre,:)=A_a(K/4+Npre,:);

    
for t=K/4+Npre-1:-1:1
    for ll=1:8
            dummy_B_b(t,ll)=0;
    end
    dummy_B_b_totall(t)=0;
    for ll=1:8
        
            [next_state_0,x_p]=nxt_stat1513(ll-1,0); 
            [next_state_1,x_p]=nxt_stat1513(ll-1,1);
           
    
                dummy_B_b(t,ll)=dummy_B_b(t+1,next_state_0+1)*R_r_0(t+1,next_state_0+1)+dummy_B_b(t+1,next_state_1+1)*R_r_1(t+1,next_state_1+1);
 
    end 
    for ll=1:8
        dummy_B_b_totall(t)=A_a_totall(t+1);
    end
    for ll=1:8
        dummy_B_b(t,ll)=dummy_B_b(t,ll)/dummy_B_b_totall(t);
    end
end
for t=K/4:-1:1
    for ll=1:8
        B_b(t,ll)=dummy_B_b(t,ll);
    end
end




 for t=1:K
     ln_1(t)=0;
     ln_0(t)=0;
     for ll=1:8
         if t==1
             
             [next_state_1,x_p]=nxt_stat1513(ll-1,1);
             ln_1(t)=ln_1(t)+A_a_initial(ll)*R_r_1(t,next_state_1+1)*B_b(t,next_state_1+1);
             [next_state_0,x_p]=nxt_stat1513(ll-1,0);
             ln_0(t)=ln_0(t)+A_a_initial(ll)*R_r_0(t,next_state_0+1)*B_b(t,next_state_0+1);
             
         else if t>1&t<K
             
             [next_state_1,x_p]=nxt_stat1513(ll-1,1);
             ln_1(t)=ln_1(t)+A_a(t-1,ll)*R_r_1(t,next_state_1+1)*B_b(t,next_state_1+1);
             [next_state_0,x_p]=nxt_stat1513(ll-1,0);
             ln_0(t)=ln_0(t)+A_a(t-1,ll)*R_r_0(t,next_state_0+1)*B_b(t,next_state_0+1);
         else
             [next_state_1,x_p]=nxt_stat1513(ll-1,1);
             ln_1(t)=ln_1(t)+A_a(t-1,ll)*R_r_1(t,next_state_1+1)*B_b_K(next_state_1+1);
             [next_state_0,x_p]=nxt_stat1513(ll-1,0);
             ln_0(t)=ln_0(t)+A_a(t-1,ll)*R_r_0(t,next_state_0+1)*B_b_K(next_state_0+1);
           
         end
     end
     
     ln_c(t)=log(ln_1(t)/ln_0(t));
     if ln_c(t)>0
         decision_c(t)=1;
     else
         decision_c(t)=0;
     end
     
 end
end 

extrin_information=ln_c+4*CONS*x-prior_information;