rs_cnv_awgn.m

由RS码、卷积码构成的级联码

close all;
clear all;
input = randint(1,10^4,2); % Random data
%######################################################################
%将输入序列拆分成d部分，每部分含K*m个信息元；
%再将每部分的K*m个信息元分成K段，将每段的m个二进制信息元转换成一个十进制信息元
%这样输入信息就可看成有d部分，每部分含K个十进制码元的码组
m=8;
N=2^m-1;
t=16;
K=N-2*t;
pp = 285;    %primitive polynomial本原多项式
pt = poly2trellis(7,[155 117]); % Define trellis.
delay = 48;
%######################################################################
%查看是否需要补0，使输入input是K*m的整数倍
if rem(length(input),K*m)>0
    input=[input,zeros(size(1:K*m-rem(length(input),K*m)))];
end
%######################################################################
d= length(input)/(K*m);
SNRindB=0:10;
 %######################################################################
 %rs编译码
%  a=1;
%  while a==1
%      for w=1:length(SNRindB)
%          for i=1:d
%              for j=1:K
%                  x=input(1+(j-1)*m+(i-1)*(K*m):j*m+(i-1)*(K*m));
%                  input_rs_deci(1,j)=bin2deci(x);%二进制转换为十进制
%              end
%              input_rs = gf(input_rs_deci, m, pp);
%              code_rs= rsenc(input_rs, N, K);
%              code_in_rs= gf_double(code_rs);%gf域转换到double型
%              code_out_rs = awgn(code_in_rs, SNRindB(w),0);
%              for f=1:N
%              if code_out_rs(f)>N
%                 code_out_rs(f)=N;
%              end
%              if code_out_rs(f)<0
%                 code_out_rs(f)=0;
%              end
%              end
%             %######################################################################
%             decode_out_rs=gf(code_out_rs,m,pp);
%             [dec_rs0,cnumerr] = rsdec(decode_out_rs,N,K);
%             dec_rs=double(dec_rs0.x);
%             for j=1:K
%                 p=dec_rs(1,j);
%                 output_rs(1+(j-1)*m+(i-1)*(K*m):j*m+(i-1)*(K*m))=deci2bin(p,m);%十进制转换为二进制
%             end
%          end
%          pe_rs=length(find(output_rs~=input))/length(input);
%              err_rs(1,w)=pe_rs;
%      end
%   figure;
%      semilogy(SNRindB,err_rs,'b*-');
%     hold on
%     axis([0,10,1e-6,1]);
%     xlabel('SNR(dB)');
%     ylabel('Pe');
%     a=2;
%  end
 %######################################################################
 %卷积编码
 a=2;
 while a==2
     for w=1:length(SNRindB)
         code_cnv1 = convenc(input,pt); 
         code_cnv2 = convenc([input(length(input)-delay+1:end),zeros(1,delay)],pt); 
         %######################################################################
      
         code_out_cnv1 = awgn(code_cnv1, SNRindB(w),'measured');
         code_out_cnv2 = awgn(code_cnv2, SNRindB(w),'measured');
         %######################################################################
         %内码（卷积）译码
         %八电平量化
         qcode1 = quantiz(code_out_cnv1,[0.001,.1,.3,.5,.7,.9,.999]);
         qcode2 = quantiz(code_out_cnv2,[0.001,.1,.3,.5,.7,.9,.999]);
         decode_cnv_delay1 = vitdec(qcode1,pt,delay,'cont','soft',3); 
         decode_cnv_delay2 = vitdec(qcode2,pt,delay,'cont','soft',3);
         decode_cnv1=[decode_cnv_delay1(delay+1:end)];
         decode_cnv2=[decode_cnv_delay2(delay+1:end)];
         decode_cnv=[decode_cnv1,decode_cnv2];
         
         pe_cnv=length(find(decode_cnv~=input))/length(input);
         err_cnv(1,w)=pe_cnv;
     end

     semilogy(SNRindB,err_cnv,'g.-');
     hold on;
     axis([0,10,1e-6,1]);
     xlabel('SNR(dB)');
     ylabel('Pe');
     a=3;
 end 

 %######################################################################
  %rs+卷积
  a=3;
 while a==3
     for w=1:length(SNRindB)
         for i=1:d
             for j=1:K
                 x=input(1+(j-1)*m+(i-1)*(K*m):j*m+(i-1)*(K*m));
                 input_deci(1,j)=bin2deci(x);%二进制转换为十进制
             end
             %######################################################################
             %外码（rs）编码
             input_rs_cnv = gf(input_deci, m, pp);
             code= rsenc(input_rs_cnv, N, K);
             %######################################################################
             input_cnv_deci= gf_double(code);%gf域转换到double型
             for s=1:N
                 p=input_cnv_deci(1,s);
                 input_cnv(1+(s-1)*m:s*m)=deci2bin(p,m);%十进制转换为二进制
             end
             %######################################################################
             %内码（卷积）编码
             output_cnv1 = convenc(input_cnv,pt); 
             output_cnv2 = convenc([input_cnv(length(input_cnv)-delay+1:end),zeros(1,delay)],pt);  % Encode the data.
             %#############################################################
             % Add noise.
             code_out1 = awgn(output_cnv1, SNRindB(w),'measured');
             code_out2 = awgn(output_cnv2, SNRindB(w),'measured');
             %######################################################################
             %内码（卷积）译码
             %八电平量化
             qcode1 = quantiz(code_out1,[0.001,.1,.3,.5,.7,.9,.999]);
             qcode2 = quantiz(code_out2,[0.001,.1,.3,.5,.7,.9,.999]);
             decode_delay1 = vitdec(qcode1,pt,delay,'cont','soft',3); 
             decode_delay2 = vitdec(qcode2,pt,delay,'cont','soft',3);
             decode1=[decode_delay1(delay+1:end)];
             decode2=[decode_delay2(delay+1:end)];
             decode=[decode1,decode2];
             %######################################################################
             for s=1:N
                 x= decode (1+(s-1)*m:s*m);
                 decode_deci(1,s)=bin2deci(x);%二进制转换为十进制
             end
             %######################################################################
             decode_deci=gf(decode_deci,m,pp);
             [dec,cnumerr] = rsdec(decode_deci,N,K);
            dec=double(dec.x);%外码译码
             for j=1:K
                 p=dec(1,j);
                 output(1+(j-1)*m+(i-1)*(K*m):j*m+(i-1)*(K*m))=deci2bin(p,m);%十进制转换为二进制
             end
         end
         [num,pe] = biterr(output,input);
         err(1,w)=pe;
     end
     semilogy(SNRindB,err,'k*-');
     a=4;
 end
 %######################################################################
  %rs+卷积+交织
  a=4;
 while a==4
     for w=1:length(SNRindB)
         for i=1:d
             for j=1:K
                 x=input(1+(j-1)*m+(i-1)*(K*m):j*m+(i-1)*(K*m));
                 input_deci_jz(1,j)=bin2deci(x);%二进制转换为十进制
             end
             %######################################################################
             %外码（rs）编码
             input_rs_cnv_jz = gf(input_deci_jz, m, pp);
             code_jz= rsenc(input_rs_cnv_jz, N, K);
             %######################################################################
             input_cnv_deci_jz= gf_double(code_jz);%gf域转换到double型
             for s=1:N
                 p=input_cnv_deci_jz(1,s);
                 input_cnv_jz(1+(s-1)*m:s*m)=deci2bin(p,m);%十进制转换为二进制
             end
             %######################################################################
             %内码（卷积）编码
             output_cnv_jz01 = convenc(input_cnv_jz,pt); 
             output_cnv_jz02 = convenc([input_cnv_jz(length(input_cnv_jz)-delay+1:end),zeros(1,delay)],pt);  % Encode the data.
             %#############################################################
             %交织
             c=16;
             b1=length(output_cnv_jz01)/c;
             b2=length(output_cnv_jz02)/c;
             for q=1:c
                 for p=1:b1
                     output_cnv_jz1(q,p)=output_cnv_jz01(1,c*(p-1)+q);
                 end
             end
               for q=1:c
                 for p=1:b2
                     output_cnv_jz2(q,p)=output_cnv_jz02(1,c*(p-1)+q);
                 end
             end
             % Add noise.
             code_out_jz01 = awgn(output_cnv_jz1, SNRindB(w),'measured');
             code_out_jz02 = awgn(output_cnv_jz2, SNRindB(w),'measured');
              %解交织
             for p=1:b1
                 for q=1:c
                     code_out_jz1(1,c*(p-1)+q)=code_out_jz01(q,p);
                 end 
             end
              for p=1:b2
                 for q=1:c
                     code_out_jz2(1,c*(p-1)+q)=code_out_jz02(q,p);
                 end 
             end
             %######################################################################
             %内码（卷积）译码
             %八电平量化
             qcode_jz1 = quantiz(code_out_jz1,[0.001,.1,.3,.5,.7,.9,.999]);
             qcode_jz2 = quantiz(code_out_jz2,[0.001,.1,.3,.5,.7,.9,.999]);
             decode_delay_jz1 = vitdec(qcode_jz1,pt,delay,'cont','soft',3); 
             decode_delay_jz2 = vitdec(qcode_jz2,pt,delay,'cont','soft',3);
             decode_jz1=[decode_delay_jz1(delay+1:end)];
             decode_jz2=[decode_delay_jz2(delay+1:end)];
             decode_jz=[decode_jz1,decode_jz2];
             %#############################################################
             for s=1:N
                 x= decode_jz (1+(s-1)*m:s*m);
                 decode_deci_jz(1,s)=bin2deci(x);%二进制转换为十进制
             end
             %######################################################################
             decode_deci_jz=gf(decode_deci_jz,m,pp);
             [dec_jz,cnumerr] = rsdec(decode_deci_jz,N,K);
            dec_jz=double(dec_jz.x);%外码译码
             for j=1:K
                 p=dec_jz(1,j);
                 output_jz(1+(j-1)*m+(i-1)*(K*m):j*m+(i-1)*(K*m))=deci2bin(p,m);%十进制转换为二进制
             end
         end
         [num,pe] = biterr(output_jz,input);
         err_jz(1,w)=pe;
     end
     semilogy(SNRindB,err_jz,'ro-');
     a=5;
 end
  %######################################################################
  %未编码
    for w=1:length(SNRindB)
  input_ori=awgn(input,SNRindB(w),'measured');
   qcode_ori = quantiz(input_ori,[0.5]);
     pe_ori=length(find(qcode_ori~=input))/length(input);
     err_ori(1,w)=pe_ori;
    end
     semilogy(SNRindB,err_ori,'b.--');  
 legend('卷积','rs+卷积','rs+卷积+交织','未编码');