LT_BP_performance

clc;clear;
m = 500;
K = 500;
overhead = [550 600 700 800];
i_ii =10000;
for iii = 1:4
    n = overhead(iii);
    for i_i = 1:i_ii
        if mod(i_i,100)==0
            i_i
        end
        c = 0.04;
        o = 0.5;
        x = ones(1,m);
        x_1 = rand(1,m)<0.5 ;           %信源碼字
        Pgenerate = rand(1,n);          %隨機生成n個數
        degree = [1:K];
        %********************************************
        %         for d = 1:K
        %             index_d = 2:K;
        %             if find(d==index_d)
        %                 PP = 1/d/(d-1);
        %             elseif d==1
        %                 PP = 1/K;
        %             end
        %             P(d) = PP;
        %         end
        %         P;
        %*********************************************
        S = c*sqrt(K)*log(K/o);
        Z = 0;
        for d = 1:K
            index_d = 2:K;
            if find(d==index_d)
                PP = 1/d/(d-1);
            elseif d==1
                PP = 1/K;
            end
            P_P(d) = PP;
            d_index = 1:round(K/S)-1;
            if (find(d==d_index))
                tt = S/K/d;
            elseif d==round(K/S)
                tt = S/K*log(S/o);
            else
                tt = 0;
            end
            t(d) = tt;
            Z = Z + P_P(d)+t(d);
        end
        P_P;
        t;
        Z;
        for d = 1:K
            P(d) = (P_P(d)+t(d))/Z;
        end
        P;
        sum = 0;
        for i = 1:K
            sum = sum + P(i);
            sum_1(i) = sum;
        end
        sum_1;
        for j = 1:length(Pgenerate)
            for i = 1:length(sum_1)
                if Pgenerate(j)<sum_1(i)
                    label(j) = i;
                    break;
                end
            end
        end
        label;
        for j = 1:n
            a = zeros(1,m);
            vector = randperm(m,label(j));
            for i = 1:length(vector)
                a(vector(i)) = 1;
            end
            G(1:m,j) = a;
        end
        G;
        y = mod(x_1*G,2);
        while(1)
            %***********判斷準備開始1*****************
            for i = 1:n
                front_judge(i) = length(find(G(:,i)==1));
            end
            %***********判斷準備結束1*****************
            count = zeros(1,n);
            for j = 1:n
                index_1 = find(G(:,j)==1);
                if length(index_1)==1     %條件滿足
                    count(j) = index_1;
                    x(index_1) = y(j);
                    G(index_1,j) = 0;
                end
            end
            t = 1;
            countt = zeros(1,length(find(count~=0)));
            for j = 1:n
                if count(j)~=0
                    countt(t) = count(j);
                    t = t + 1;
                end
            end
            countt;
            for i = 1:length(countt)         %橫坐標count
                for j = 1:n
                    if G(countt(i),j) == 1
                        y(j) = xor(x(countt(i)),y(j));
                        G(countt(i),j) = 0;
                    end
                end
            end
            %***********判斷準備開始2*****************
            for i = 1:n
                last_judge(i) = length(find(G(:,i)==1));
            end
            %***********判斷準備結束2*****************
            if last_judge==front_judge
                break;
            end
        end
        x_2 = x;
        check = 0;
        for i = 1:length(x_1)
            if x_1(i)~=x_2(i)
                check = check + 1;
            end
        end
        check_1(i_i) = check/length(x_1);
    end
    check_2(iii,1:i_ii) = check_1(1:i_ii);
end
for i = 1:4
    berr = 0;
    for j = 1:i_ii
        berr = berr + check_2(i,j);
    end
    ber(i) = berr/i_ii;    
end
ber
for i = 1:4
    ferr = 0;
    for j = 1:i_ii
        if check_2(j)~=0
            ferr = ferr+1;
        end
    end
    fer(i) = ferr/i_ii;
end
fer







%         while(1)
%             %***********判斷準備開始1*****************
%             for i = 1:n
%                 front_judge(i) = length(find(G(:,i)==1));
%             end
%             %***********判斷準備結束1*****************
%             %             for j = 1:n
%             %                 index_y = find(G(:,j)==1);
%             %                 if length(index_y)==1 %度1
%             %                     x(index_y) = y(j);%賦值
%             %                     G(index_y,j) = 0; %刪邊
%             %                     for jj = 1:n
%             %                         if G(index_y,jj)==1
%             %                             y(jj) = xor(x(index_y),y(jj));
%             %                             G(index_y,jj)=0;
%             %                         end
%             %                     end
%             %                 end
%             %             end
%             jj = 1;
%             for j = 1:n
%                 index_1 = find(G(:,j)==1);
%                 if length(index_1)==1     %條件滿足
%                     count(jj) = index_1;
%                     jj = jj + 1;
%                     x(index_1) = y(j);
%                     G(index_1,j) = 0;
%                 end
%             end
%             G;
%             count; %對於心願節點來説是1，3
%             for i = 1:length(count)         %橫坐標count
%                 index_2 = find(G(i,:)==1);  %縱坐標index_2
%                 for j = 1:length(index_2)
%                     y(index_2(j)) = xor(count(i),index_2(j));
%                     G(count(i),index_2(j)) = 0;
%                 end
%             end
%             %***********判斷準備開始2*****************
%             for i = 1:n
%                 last_judge(i) = length(find(G(:,i)==1));
%             end
%             %***********判斷準備結束2*****************
%             if last_judge==front_judge
%                 break;
%             end
%         end


%*********************************************
% degree = [1:K];
% for i = 1:K-1
%     for j = i:K
%         if P(i)>P(j)
%             t = P(i);
%             P(i) = P(j);
%             P(j) = t;
%             tt = degree(i);
%             degree(i) = degree(j);
%             degree(j) = tt;
%         end
%     end
% end
% P1 =[  P]              %從小到達排序
% degree01 = degree ;  %對應度值
 
% for j = 1:length(Pgenerate)
%     for i = 1:length(P1)
%         if Pgenerate(j)<=P1(i)
%             label(j) = degree01(i);
%             break;
%         else
%             label(j) = degree01(length(P1));
%             break;
%         end
%     end
% end
% label;
%