turbo_qam_2a.m

Turbo编码16QAM在AWGN信道中的性能Turbo码份量码为[1,17/15],码率为1/2,采用随机交织器

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Turbo编码16QAM在AWGN信道中的性能
%% Turbo码份量码为(1,17/15),码率为1/2,采用随机交织器
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function    Turbo_QAM_2a
    
    INFO_LEN = 8000;                                                        %% 信息比特数目
    SYMBOL_LEN = INFO_LEN/2;                                                %% 16QAM符号数目
    CODE_RATE = 2.0;                                                        %% Turbo码的码率(bit/Symbol)
    AddrTable = randperm(INFO_LEN);                                         %% 生成随机交织器的地址
    Trellis = poly2trellis(4,[15,17],15);
    
    for BitSnr = 2:0.1:8,
        LinearSnr = 10^(BitSnr/10);        
        NoiseVar = 10/(2*CODE_RATE*LinearSnr);
        BitSum = 0;
        BitErrSum = 0;
        while BitErrSum < 200,
            %% Turbo 码的编码
            InfoData = randint(1,INFO_LEN);
            CodeData1 = convenc(InfoData,Trellis,0);                        %% 第一个RSC编码器  
            CheckData1 = CodeData1(2:2:end);                                %% 取校验比特序列
            IntrlvInfo = intrlv(InfoData,AddrTable);                        %% 信息序列交织
            CodeData2 = convenc(IntrlvInfo,Trellis,0);                      %% 第二个RSC编码器
            CheckData2 = CodeData2(2:2:end);                                %% 取校验比特序列
            for i =1:(INFO_LEN/2),                                          %% 对校验序列进行删余
                CheckData(i*2-1) = CheckData1(i*2-1);                       
                CheckData(i*2) = CheckData2(i*2);
            end
                
            %% 16QAM 调制
            StarTable = [-1 -3 1 3];
            MapData = InfoData*2 + CheckData;
            for i=1 : SYMBOL_LEN,
                IndexRe = MapData(i*2-1) + 1;
                IndexIm = MapData(i*2) + 1;
                Modulate(i) = StarTable(IndexRe) + j*StarTable(IndexIm);
            end

            %% 加上白噪声
            NoiseRe = sqrt(NoiseVar)*randn(1,SYMBOL_LEN);
            NoiseIm = sqrt(NoiseVar)*randn(1,SYMBOL_LEN);
            WhiteNoise = NoiseRe + j*NoiseIm;
            Receive = Modulate + WhiteNoise;
            
            %% 计算每个映射比特的度量
            ReceiveRe = real(Receive);
            ReceiveIm = imag(Receive);
            for i=1:SYMBOL_LEN,
                [Metric1,Metric2] = QamBitMetric(ReceiveRe(i),NoiseVar);
                SoftInfo(i*2-1) = Metric1;
                SoftCheck(i*2-1) = Metric2;
                [Metric1,Metric2] = QamBitMetric(ReceiveIm(i),NoiseVar);
                SoftInfo(i*2) = Metric1;
                SoftCheck(i*2) = Metric2;
            end
            for i = 1:(INFO_LEN/2),
                FirstCheck(i*2-1) =  SoftCheck(i*2-1);
                FirstCheck(i*2) = 0;
                SecondCheck(i*2-1) = 0;
                SecondCheck(i*2) = SoftCheck(i*2);
            end
            
            %% Turbo 码的迭代译码
            ExtInfo = zeros(1,INFO_LEN);                                        %% 初始化外信息
            for Iter = 1:8,
                FirstInput = SoftInfo + ExtInfo;                                %% 计算第一个译码器的输入
                FirstOutput = MapDecode(FirstInput,FirstCheck,INFO_LEN);        %% 第一个MAP译码
                FirstInfo = FirstOutput - ExtInfo;                              %% 减去外信息
                SecondInput = intrlv(FirstInfo,AddrTable);                      %% 输出结果进行交织
                SecondOutput = MapDecode(SecondInput,SecondCheck,INFO_LEN);     %% 第二个MAP译码器
                NewExtInfo = SecondOutput - SecondInput;                        %% 计算新的外信息
                ExtInfo = deintrlv(NewExtInfo,AddrTable);                       %% 解交织后作为下次迭代的外信息
            end        
            DetectData = int8(SecondOutput>=0);                                 %% 第二个MAP译码器输出结果判决
            DecodeData = deintrlv(DetectData,AddrTable);                        %% 解交织

            %% 统计误比特率
            BitErrCount = sum(InfoData ~= DecodeData);
            BitErrSum = BitErrSum + BitErrCount;
            BitSum = BitSum + INFO_LEN;
        end
        BitErrRate = BitErrSum/BitSum;
        sprintf('Eb/N0 = %f , BER = %e \n',BitSnr,BitErrRate)
    end
return


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 4AM映射符号的比特度量的计算
%%  映射方式 bit1 bit2  00 -> -1  01 -> -3
%%                      10 -> 1   11 -> 3 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function  [Metric1,Metric2] = QamBitMetric(Recv,NoiseVar)
    Prob1 = exp(-(Recv-3)^2/(2*NoiseVar)) + exp(-(Recv-1)^2/(2*NoiseVar));
    Prob0 = exp(-(Recv+3)^2/(2*NoiseVar)) + exp(-(Recv+1)^2/(2*NoiseVar));
    Metric1 = log(Prob1/Prob0);
    
    Prob1 = exp(-(Recv-3)^2/(2*NoiseVar)) + exp(-(Recv+3)^2/(2*NoiseVar));
    Prob0 = exp(-(Recv-1)^2/(2*NoiseVar)) + exp(-(Recv+1)^2/(2*NoiseVar));
    Metric2 = log(Prob1/Prob0);
return


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 概率域RSC(1,17/15)卷积码的MAP译码,分量码最后状态不归0
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function    SoftOutput = MapDecode(Info,Check,Length)

    ST0 = 1; ST1 = 2;  ST2 = 3; ST3 = 4;                                      %% 定义状态号
    ST4 = 5; ST5 = 6;  ST6 = 7; ST7 = 8; 
    Alpha0 = zeros(8,Length);
    Alpha1 = zeros(8,Length);

    %% 接收第1个数据
    Alpha0(ST0,1) = 1;
    Alpha1(ST0,1) = exp(Info(1)+Check(1));

    LSum = Alpha0(ST0,1) + Alpha1(ST0,1);
    Alpha0(ST0,1) = Alpha0(ST0,1)/LSum;
    Alpha1(ST0,1) = Alpha1(ST0,1)/LSum;
    
    %% 接收第2个数据
    Alpha0(ST0,2) = Alpha0(ST0,1);
    Alpha1(ST0,2) = Alpha0(ST0,1)*exp(Info(2)+Check(2));
    Alpha0(ST4,2) = Alpha1(ST0,1);
    Alpha1(ST4,2) = Alpha1(ST0,1)*exp(Info(2)+Check(2));
    
    LSum = sum(Alpha0(:,2)) + sum(Alpha1(:,2));
    Alpha0(:,2) = Alpha0(:,2)/LSum;
    Alpha1(:,2) = Alpha1(:,2)/LSum;
    
    %% 接收第3个数据
    Alpha0(ST0,3) = Alpha0(ST0,2);
    Alpha1(ST0,3) = Alpha0(ST0,2)*exp(Info(3)+Check(3));
    Alpha0(ST2,3) = Alpha1(ST4,2)*exp(Check(3));
    Alpha1(ST2,3) = Alpha1(ST4,2)*exp(Info(3));
    Alpha0(ST4,3) = Alpha1(ST0,2);
    Alpha1(ST4,3) = Alpha1(ST0,2)*exp(Info(3)+Check(3));
    Alpha0(ST6,3) = Alpha0(ST4,2)*exp(Check(3));
    Alpha1(ST6,3) = Alpha0(ST4,2)*exp(Info(3));
    
    LSum = sum(Alpha0(:,3)) + sum(Alpha1(:,3));
    Alpha0(:,3) = Alpha0(:,3)/LSum;
    Alpha1(:,3) = Alpha1(:,3)/LSum;
    
    %% 接收第4个数据
    Alpha0(ST0,4) = Alpha0(ST0,3);
    Alpha1(ST0,4) = Alpha0(ST0,3)*exp(Info(4)+Check(4));
    Alpha0(ST1,4) = Alpha0(ST2,3);
    Alpha1(ST1,4) = Alpha0(ST2,3)*exp(Info(4)+Check(4));
    Alpha0(ST2,4) = Alpha1(ST4,3)*exp(Check(4));
    Alpha1(ST2,4) = Alpha1(ST4,3)*exp(Info(4));
    Alpha0(ST3,4) = Alpha1(ST6,3)*exp(Check(4));
    Alpha1(ST3,4) = Alpha1(ST6,3)*exp(Info(4));
    
    Alpha0(ST4,4) = Alpha1(ST0,3);
    Alpha1(ST4,4) = Alpha1(ST0,3)*exp(Info(4)+Check(4));
    Alpha0(ST5,4) = Alpha1(ST2,3);
    Alpha1(ST5,4) = Alpha1(ST2,3)*exp(Info(4)+Check(4));
    Alpha0(ST6,4) = Alpha0(ST4,3)*exp(Check(4));
    Alpha1(ST6,4) = Alpha0(ST4,3)*exp(Info(4));
    Alpha0(ST7,4) = Alpha0(ST6,3)*exp(Check(4));
    Alpha1(ST7,4) = Alpha0(ST6,3)*exp(Info(4));
    
    LSum = sum(Alpha0(:,4)) + sum(Alpha1(:,4));
    Alpha0(:,4) = Alpha0(:,4)/LSum;
    Alpha1(:,4) = Alpha1(:,4)/LSum;

    %% Alpha的迭代计算
    for i = 5:Length,
        Alpha0(ST0,i) = (Alpha0(ST0,i-1)+Alpha1(ST1,i-1));
        Alpha1(ST0,i) = (Alpha0(ST0,i-1)+Alpha1(ST1,i-1)) * exp(Info(i)+Check(i));
        Alpha0(ST1,i) = (Alpha0(ST2,i-1)+Alpha1(ST3,i-1));
        Alpha1(ST1,i) = (Alpha0(ST2,i-1)+Alpha1(ST3,i-1)) * exp(Info(i)+Check(i));
        Alpha0(ST2,i) = (Alpha1(ST4,i-1)+Alpha0(ST5,i-1)) * exp(Check(i));
        Alpha1(ST2,i) = (Alpha1(ST4,i-1)+Alpha0(ST5,i-1)) * exp(Info(i));
        Alpha0(ST3,i) = (Alpha1(ST6,i-1)+Alpha0(ST7,i-1)) * exp(Check(i));
        Alpha1(ST3,i) = (Alpha1(ST6,i-1)+Alpha0(ST7,i-1)) * exp(Info(i));
        
        Alpha0(ST4,i) = (Alpha1(ST0,i-1)+Alpha0(ST1,i-1));
        Alpha1(ST4,i) = (Alpha1(ST0,i-1)+Alpha0(ST1,i-1)) * exp(Info(i)+Check(i));
        Alpha0(ST5,i) = (Alpha1(ST2,i-1)+Alpha0(ST3,i-1));
        Alpha1(ST5,i) = (Alpha1(ST2,i-1)+Alpha0(ST3,i-1)) * exp(Info(i)+Check(i));
        Alpha0(ST6,i) = (Alpha0(ST4,i-1)+Alpha1(ST5,i-1)) * exp(Check(i));
        Alpha1(ST6,i) = (Alpha0(ST4,i-1)+Alpha1(ST5,i-1)) * exp(Info(i));
        Alpha0(ST7,i) = (Alpha0(ST6,i-1)+Alpha1(ST7,i-1)) * exp(Check(i));
        Alpha1(ST7,i) = (Alpha0(ST6,i-1)+Alpha1(ST7,i-1)) * exp(Info(i));
        
        LSum = sum(Alpha0(:,i)) + sum(Alpha1(:,i));
        Alpha0(:,i) = Alpha0(:,i)/LSum;
        Alpha1(:,i) = Alpha1(:,i)/LSum;
    end
       
    %% Beta的初始化
    Beta0 = zeros(8,Length);
    Beta1 = zeros(8,Length);
    Beta0(:,Length) = 1/16;
    Beta1(:,Length) = 1/16;
    
    %% Beta 的迭代处理
    for i = Length:-1:5,
        Beta0(ST0,i-1) = Beta0(ST0,i) + Beta1(ST0,i)*exp(Info(i)+Check(i));
        Beta1(ST0,i-1) = Beta0(ST4,i) + Beta1(ST4,i)*exp(Info(i)+Check(i));
        Beta0(ST1,i-1) = Beta0(ST4,i) + Beta1(ST4,i)*exp(Info(i)+Check(i));
        Beta1(ST1,i-1) = Beta0(ST0,i) + Beta1(ST0,i)*exp(Info(i)+Check(i));
        Beta0(ST2,i-1) = Beta0(ST1,i) + Beta1(ST1,i)*exp(Info(i)+Check(i));
        Beta1(ST2,i-1) = Beta0(ST5,i) + Beta1(ST5,i)*exp(Info(i)+Check(i));
        Beta0(ST3,i-1) = Beta0(ST5,i) + Beta1(ST5,i)*exp(Info(i)+Check(i));
        Beta1(ST3,i-1) = Beta0(ST1,i) + Beta1(ST1,i)*exp(Info(i)+Check(i));

        Beta0(ST4,i-1) = Beta0(ST6,i)*exp(Check(i)) + Beta1(ST6,i)*exp(Info(i));
        Beta1(ST4,i-1) = Beta0(ST2,i)*exp(Check(i)) + Beta1(ST2,i)*exp(Info(i));
        Beta0(ST5,i-1) = Beta0(ST2,i)*exp(Check(i)) + Beta1(ST2,i)*exp(Info(i));
        Beta1(ST5,i-1) = Beta0(ST6,i)*exp(Check(i)) + Beta1(ST6,i)*exp(Info(i));
        Beta0(ST6,i-1) = Beta0(ST7,i)*exp(Check(i)) + Beta1(ST7,i)*exp(Info(i));
        Beta1(ST6,i-1) = Beta0(ST3,i)*exp(Check(i)) + Beta1(ST3,i)*exp(Info(i));
        Beta0(ST7,i-1) = Beta0(ST3,i)*exp(Check(i)) + Beta1(ST3,i)*exp(Info(i));
        Beta1(ST7,i-1) = Beta0(ST7,i)*exp(Check(i)) + Beta1(ST7,i)*exp(Info(i));
                
        LSum = sum(Beta0(:,i-1)) + sum(Beta1(:,i-1));
        Beta0(:,i-1) = Beta0(:,i-1)/LSum;
        Beta1(:,i-1) = Beta1(:,i-1)/LSum;
    end
    
    %% 接收第4个数据
    Beta0(ST0,3) = Beta0(ST0,4) + Beta1(ST0,4)*exp(Info(4)+Check(4));
    Beta1(ST0,3) = Beta0(ST4,4) + Beta1(ST4,4)*exp(Info(4)+Check(4));
    Beta0(ST2,3) = Beta0(ST1,4) + Beta1(ST1,4)*exp(Info(4)+Check(4));
    Beta1(ST2,3) = Beta0(ST5,4) + Beta1(ST5,4)*exp(Info(4)+Check(4));
    Beta0(ST4,3) = Beta0(ST6,4)*exp(Check(4)) + Beta1(ST6,4)*exp(Info(4));
    Beta1(ST4,3) = Beta0(ST2,4)*exp(Check(4)) + Beta1(ST2,4)*exp(Info(4));
    Beta0(ST6,3) = Beta0(ST7,4)*exp(Check(4)) + Beta1(ST7,4)*exp(Info(4));
    Beta1(ST6,3) = Beta0(ST3,4)*exp(Check(4)) + Beta1(ST3,4)*exp(Info(4)) ;
                
    LSum = sum(Beta0(:,3)) + sum(Beta1(:,3));
    Beta0(:,3) = Beta0(:,3)/LSum;
    Beta1(:,3) = Beta1(:,3)/LSum;

    %% 接收第3个数据
    Beta0(ST0,2) = Beta0(ST0,3) + Beta1(ST0,3)*exp(Info(3)+Check(3));
    Beta1(ST0,2) = Beta0(ST4,3) + Beta1(ST4,3)*exp(Info(3)+Check(3));
    Beta0(ST4,2) = Beta0(ST6,3)*exp(Check(3)) + Beta1(ST6,3)*exp(Info(3));
    Beta1(ST4,2) = Beta0(ST2,3)*exp(Check(3)) + Beta1(ST2,3)*exp(Info(3));

    LSum = sum(Beta0(:,2)) + sum(Beta1(:,2));
    Beta0(:,2) = Beta0(:,2)/LSum;
    Beta1(:,2) = Beta1(:,2)/LSum;

    %% 接收第2个数据
    Beta0(ST0,1) = Beta0(ST0,2) + Beta1(ST0,2)*exp(Info(2)+Check(2));
    Beta1(ST0,1) = Beta0(ST4,2) + Beta0(ST4,2)*exp(Info(2)+Check(2));
    
    LSum = sum(Beta0(:,1)) + sum(Beta1(:,1));
    Beta0(:,1) = Beta0(:,1)/LSum;
    Beta1(:,1) = Beta1(:,1)/LSum;
        
    %% 输出软信息
    for i = 1:Length,
        LSum0 = Alpha0(:,i)' * Beta0(:,i);
        LSum1 = Alpha1(:,i)' * Beta1(:,i);
        SoftOutput(i) = log(LSum1/LSum0);
    end
return