ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

項目一:LDPC在碼刪除信道下譯碼

2019-07-10 20:03:34  阅读:209  来源: 互联网

标签:項目 end Hr tt else 下譯碼 刪除 排序 order


clc;clear;
for loop = 1:1000
    load PEGH_3_6_1024.mat
    [m,n] = size(H);
    %    H_simple = zhu_rank(H);                    %H化为行最简阵
    %***********************高斯消元开始****************************************
    a = H;
    for j = 1:m
        for i = j:m-1
            for k = i+1:m
                if a(i,j)==0 && a(k,j)==1
                    a([i,k],:)=a([k,i],:);
                end
            end
        end
        for i = j+1:m
            if a(i,j)==1 &&a(j,j)==1
                a(i,:)=xor(a(i,:),a(j,:));
            end
        end
    end
    %F_F=a
    %*********************调整顺序*****************
    %********统计******
    for i = 1:m
        for j = i:n     %最大值只能统计到n,不能统计到m+n,m为单位矩阵
            if a(i,j)==1
                order(i) = j;
                break;
            else
                order(i) = m+1;%全零视为越界
            end
        end
    end
    order_1 = order;%读取顺序      %未整理前每行第一个不为0数字所在位置组成的向量
    %*********第一次排序*****
    for i = 1:m-1     %i,j是冒泡前后排序,不是行列排序
        for j = i+1:m
            if order_1(i)>order_1(j)
                t = a(i,:);
                a(i,:) = a(j,:);
                a(j,:) = t;
                t = order_1(i);
                order_1(i) = order_1(j);
                order_1(j) = t;
            end
        end
    end
    order_2 = order_1;
    %调整顺序   % 整理后没行第一个不为0数字所在位置组成的向量
    a_begin = a;       %整理后的排序
    %***************第一次行最简化简*******************************
    for j = m:-1:1              %列元素排序,和上面行元素排序不同,不可替换
        for i = j:-1:1          %****确定位置***
            if a(i,j)~=0
                countt(j) = i;  %找到每一列最后一个非零元素所在的位置
                break;          %找到每一列从下至上第一个非零元素后,跳出本轮循环
            else
                countt(j) = m+1;
            end
        end
    end
    countt;
    for j = m:-1:1
        for i = countt(j)-1:-1:1  %全零列不影响
            if a(i,j)~=0
                a(i,:) = xor(a(i,:),a(countt(j),:));
            end
        end
    end
    a_simple = a;
    %********************************************************
    
    %********第二次统计&排序*******
    for i = 1:m
        for j = i:n     %最大值只能统计到n,不能统计到m+n,m为单位矩阵
            if a(i,j)==1
                count = j;
                break;
            else
                count = m+1;
            end
        end
        order(i)= count;
    end
    order_3 = order;  %读取顺序,便于行变换
    for i = 1:m-1    %i,j是冒泡前后排序,不是行列排序
        for j = i+1:m
            if order_3(i)>order_3(j)
                tt = a(i,:);
                a(i,:) = a(j,:);
                a(j,:) = tt;
                t = order_3(i);
                order_3(i) = order_3(j);
                order_3(j) = t;
            end
        end
    end
    order_4 = order_3;  %最终顺序
    a_third = a ;
    %*********第二次行最简化简:消除第一次化简的影响(非行阶梯)*******************
    for j = m:-1:1
        for i = j:-1:1
            if a(i,j)~=0
                record(j) = i;
                break;
            else
                record(j) = m+1;
            end
        end
    end
    record;
    for j = m:-1:1
        for i = record(j)-1:-1:1
            if a(i,j)~=0
                a(i,:) = xor(a(i,:),a(record(j),:));
            end
        end
    end
    a_last = a;
    %*********第二次化简完成************************************************
    %G_G = trace(a(1:m,1:m));
    BB= a(1:m,n-m+1:n);                           %变换后后半部分矩阵
    H_simple=a(1:m,1:n);                           %变换后前半部分矩阵
    %***********************高斯消元结束****************************************
    
    %***********************生成矩阵运行开始************************************
    Hs = H(1:m,1:n-m);
    Hr = H(1:m,n-m+1:n);
    H_H = [Hr Hs];       %只进行列变换
    G = zhu_rank(H_H);                        %%%%%%
    %*******先考虑列排序,再考虑行变换******
    for j = m:-1:1
        for i = j:-1:1
            if G(i,j)~=0
                G_order(j) = i;
                break;
            else
                G_order(j) = m+1;
            end
        end
    end
    G_order ;         %记录列向量中从下至上首个非零元素
    %****************预备工作完成*********************
    %****************开始第一次列变换***********************
    for i = 1:m-1
        for j = i+1:m
            if G(i,i)==0   %和冒泡排序的区别
                if G_order(i)>G_order(j) %G_order元素大小进行排序,对应G的列元素进行
                    t = G_order(i);         %顺序
                    G_order(i) = G_order(j);
                    G_order(j) = t;
                    tt = G(:,i);            %化简
                    G(:,i) = G(:,j);
                    G(:,j) = tt;
                    ttt = H_H(:,i);
                    H_H(:,i) = H_H(:,j);    %影子
                    H_H(:,j) = ttt;
                end
            end
        end
    end
    G_order_1 = G_order;     %从小至大依次排序
    G_1 = G;
    %**************第一次列变换完成*************************
    %******************消除1 2 2 3 4 6 6影响**********
    for i = 1:m
        if G_order_1(i)~=i
            for j = i:m-1
                for k = i+1:m
                    if G_order_1(k) == i
                        t = G_order_1 (i);
                        G_order_1 (i) = G_order_1 (k);
                        G_order_1 (k) = t;
                        tt = G_1(:,k);
                        G_1(:,k) = G_1(:,j);
                        G_1(:,j) = tt;
                        ttt = H_H(:,k);
                        H_H(:,k) = H_H(:,j);
                        H_H(:,j) = ttt;
                    end
                end
            end
        end
    end
    G_order_2 = G_order_1;      %重新排序
    G_2 = G_1;
    %******************前m列消除工作完成***********************
    %******************与m列后面的列向量进行置换开始************
    Hr_ji = trace(G_2(1:m,1:m));
    if Hr_ji < m
        for i = 1:m
            if G_2(i,i)==0                     %Hr 定位非满秩列
                for j = m+1:n                 %G矩阵进行变换
                    G_2(:,[i,j])=G_2(:,[j,i]);    %列变换
                    H_H(:,[i,j]) = H_H(:,[j,i]);
                    kkk = zhu_rank(G_2);
                    k = trace(kkk(1:m,1:m));
                    if k==Hr_ji+1
                        Hr_ji = Hr_ji+1 ;      %秩加1
                        break;
                    else
                        G_2(:,[i,j])=G_2(:,[j,i]); %换回来
                        H_H(:,[i,j]) = H_H(:,[j,i]);
                    end
                end
            end
        end
        H_H_finished = H_H;
        test = zhu_rank(H_H_finished);
    else
        %disp('原Hr矩阵满秩');
        test = G_2;
    end
    middle_G = test(1:m,m+1:n);
    G =[eye(n-m) middle_G'];%eye(n-m)
    %***********************生成矩阵运行结束************************************
    Cs = rand(1,m)<0.5;                        %随机生成信息位
    C = mod(Cs*G,2);                           %利用生成矩阵G生成完整码字
    W = zeros(8,n);
    %***********************randperm随机删除码字开始***********************************
% % %     for del_rate = 36:43
% % %         C1 = C;
% % %         rand_del = randperm(n,fix(del_rate*n*0.01));
% % %         for bb = 1:length(rand_del)
% % %             C1(rand_del(bb)) = 8;              %将删除码字赋值8
% % %         end
% % %         C_del = C1;  
    %***********************randperm随机删除码字结束************************************   
    %***********************概率随机删除码字開始************************************
    for del_rate = 40%36:43
        C1 = C;
        rand_del = rand(1,n)<del_rate*0.01;
        find_1 =find(rand_del==1);
        length(find_1);
        for bb = 1:length(find_1)
            C1(find_1) = 8;
        end
        C_del = C1;
    %***********************概率随机删除码字结束************************************
        %*********************删除信道译码开始**************************************
        HH = H;
        b = zeros(1,n);
        while(length(find (C_del(1,1:n)==8))~=0)                          %&&length(find(HH(1:m,1:n)==1))~=0)
            sect_first = zeros(1,m);
            for i = 1:m
                sect_first(i) = length(find(HH(i,1:n))==1);
            end
            for j = 1:n
                for i = 1:m
                    if HH(i,j)==1&&C_del(j)~=8
                        b(i) = xor(b(i),C_del(j));
                        HH(i,j) = 0;               % 刪除邊
                    end
                end
            end
            for i = 1:m
                hop = find(HH(i,1:n)==1);
                if length(hop)==1
                    if C_del(hop)==8
                        C_del(hop) = b(i);
                        HH(i,hop) = 0;
                    else
                        HH(i,hop) = 0;
                    end
                end
            end
            sect_last = zeros(1,m);
            for i = 1:m
                sect_last(i) = length(find(HH(i,1:n))==1);%%HH notH
            end
            %             %************判斷環開始*******************
            sect_check = (sect_first == sect_last);
            sum_jud = 0;
            for i = 1:m
                sum_jud = sum_jud + sect_check(i);
            end
            if sum_jud == m
                break;
            end
            %************判斷環結束*******************
        end
        %*********************删除信道译码完成**************************************
        W(del_rate-35,1:n)=C_del;                 %Y_1译码输出码字
    end
    %*************duo個碼字檢驗開始*****************
    for i = 1:8
        sum = 0;
        for j = 1:n
            if W(i,j)~=C(1,j)
                sum = sum +1;
            end
        end
        BER(i,loop) = sum/n;
    end
end
BER

for i=1:8
    sum = 0;
    for j = 1:loop
        sum = sum + BER(i,j);
    end
    BER_averg(i) = sum/n;
end
BER_averg

t = 0.36:0.01:0.43;
y = 0,0,0,0,0.0002,0.0098,0.0802,0.1964;
plot(t,y)
        
%*************duo個碼字檢驗結束*****************
    
   
%*********************删除信道译码结束**************************************
%***************************误码率计算开始**********************************
%         sum = 0;
%         for i = 1:n
%             if W(i) == C(i)
%                 count(i) = 0;
%             else
%                 count(i) = 1;
%             end
%             sum = sum + count(i);
%         end
%         Ber(x-19,loop) = sum/n;
%**************************误码率计算结束************************************
    %end
%end



% 
% %*************统计平均开始*****************
% for i = 1:8
%     BER = 0;
%     for j = 1:loop
%         BER = BER + Ber(i,j);
%     end
%     BER_1(i) = BER/loop;
% end
% BER_1
% %*************统计平均结束*****************




%**********************錯誤譯碼開始(無迭代)********************************
% % % for j = 1:n
% % %     for i = 1:m
% % %         if HH(i,j)==1&&(C_del(j)==0||C_del(j)==1)
% % %             b(i) = xor(b(i),C_del(j));             %更新
% % %             HH(i,j) = 8;                          %删除
% % %             C_yi(j) = C_del(j);                   %记录
% % %             k = find(HH(i,1:n)==1);               %统计
% % %             if length(k)==1                       %传送
% % %                 if C_del(k)==8
% % %                     C_del(k) = b(i);
% % %                     HH(i,k) = 8 ;                     %删除
% % %                     if k<j                            %子程序
% % %                         for i_son = 1:m
% % %                             if HH(i_son,k)==1
% % %                                 b(i_son) = xor(b(i_son),C_del(k)); %子更新
% % %                                 HH(i_son,k) = 8;
% % %                             end
% % %                         end
% % %                     end
% % %                 else
% % %                     HH(i,k)=8;
% % %                 end
% % %             end
% % %         end
% % %     end
% % % end
%****************************無迭代運算結束*********************************************















































%Hr_ji = trace(G);                 %利用最简阵的迹的和判断Hr是否满秩
% ber_1 = 0;
% ber_2 = zeros(1,m);
% for x = 1:7
%     for xx = 1:1000        
%         %*******************************************************************
%         %                        Hr通过列变换变为满秩矩阵
%         %*******************************************************************       
%         if Hr_ji < m
%             for i = 1:m
%                 if G(i,i)~=1                    %Hr 定位非满秩列
%                     for j = 1:m                  %Hs 扫列
%                         tt  = Hs(:,j);           %列交换
%                         Hs(:,j) = Hr(:,i);
%                         Hr(:,i) =  tt;
%                         k = trace(zhu_rank(Hr));
%                         if k==Hr_ji+1
%                             Hr_ji = Hr_ji+1;      %退步
%                             break                 %返回上一个循环
%                         end
%                     end
%                 end
%             end
%             CC = zhu_rank(Hr);                     %Hr最简矩阵
%             Hr_ji_1 = trace(CC );        %Hr经过列变换后的迹
%         end
%         %**********************************************************************
%         Hr_1 = zhu_nijuzhen(Hr);           %Hr_1为Hr的逆矩阵
%         yz = mod(Hr_1*Hr,2);               %验证Hr确实是Hr的逆矩阵
%         Cs = rand(1,m)<0.5;               %系统码字,rand生成
%         %*********Cr = Hr_1*Hs*Cs*********************************************
%         HH = mod(Hr_1*Hs,2);
%         Cr = mod(HH*Cs',2);                %Cr校验码字计算得到
%         C = [Cs Cr'];                      %生成码字C
%         %**********************************************************************
%         %****************此时H矩阵已经发生变化
%         H = [Hs,Hr];                       %校验矩阵H  
%         %*********************************************************************
%         a = zeros(1,m);                       % 生成验证码字a
%         k = fix(n*((36+x)/100));
% %         k = fix(0.1*x*n);
%         suiji_shan = randperm(n,k);           %生成随机删除码字的地址,是一个向量
%         for bb = 1:length(suiji_shan)
%             C(suiji_shan(bb)) = 8;             %将删除码字赋值8
%         end
%         Y = C;                                %经过删除信道删除后的码字
%         k = delete_decode(H,a,Y);             %经过删除信道译码后的码字
%         ber(x,xx) = sum(xor(C,k));               %错误译码数量
%       end
% %     Ber = ber(x,:);
% end

标签:項目,end,Hr,tt,else,下譯碼,刪除,排序,order
来源: https://blog.csdn.net/zujiusheng/article/details/95371163

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有