for i=0:number_of_states-1
    for j=0:fanout-1
        [next_state,out_put]=G_func1(i,j,L,fanout);
        nextstate(i+1,j+1)=next_state;
        output(i+1,j+1)=out_put;
        input(i+1,next_state+1)=j;
    end
end 
% ********************************************************* %


% ********************維特比譯碼部分**********************%

depth_of_trellis=length(source_coded);%*******************
Eb=1;
% Es=Eb*BPH*(1/2);
Es=Eb*BPH*(kRS/nRS);% Es=Eb*BPH*(k/n)
% Ej0=(Eb*number_of_out*Q)/(Eb_to_Nj);% 每個(gè)多音干擾的能量Ej0
Ej0=(Es*number_of_out)/(BPH*Q*Eb_to_Nj); % 每跳時(shí)間內(nèi)每個(gè)多音干擾的能量Ej0
sgma=sqrt(Eb/(2*Eb_to_No));% AWGN的均方根
demod_input=zeros(number_of_out,depth_of_trellis);
f=zeros(1,depth_of_trellis);
rc=zeros(1,number_of_out);
rs=zeros(1,number_of_out);
D=0; % D 記錄網(wǎng)格圖的當(dāng)前狀態(tài),這里初始狀態(tài)是0狀態(tài)

% *******************信道和非相干解調(diào)部分:加多音干擾和噪聲,然后非相干解調(diào)****************** %
for i=1:depth_of_trellis % i表示網(wǎng)格圖的時(shí)間走勢(shì)
    f(i)=output(D+1,source_coded(i)+1); % f(i)是i時(shí)刻的分支轉(zhuǎn)移輸出,即i時(shí)刻的跳頻頻率號(hào),頻率號(hào)范圍是[0,number_of_out-1]而不是[1,number_of_out]
    theta_s=2*pi*rand;% 跳頻信號(hào)的相位
    theta_j=2*pi*rand;% 干擾音的相位
%   J=randint(1,Q,number_of_out); % J 矩陣中存放Q個(gè)干擾音所在的頻率號(hào),干擾音所在頻率號(hào)范圍也是[0,number_of_out-1]而不是[1,number_of_out]
    J=gen_multijammer(Q,number_of_out);% J 矩陣中存放Q個(gè)干擾音所在的頻率號(hào),干擾音所在頻率號(hào)范圍也是[0,number_of_out-1]而不是[1,number_of_out]
    for j=0:number_of_out-1
        if (j==f(i))
            rc(j+1)=sqrt(Es)*cos(theta_s)+sgma*randn;
            rs(j+1)=sqrt(Es)*sin(theta_s)+sgma*randn;
        else
            rc(j+1)=sgma*randn;
            rs(j+1)=sgma*randn;
        end
    end
    for k=1:Q
        for j=0:number_of_out-1
            if (j==J(k))
                rc(j+1)=rc(j+1)+sqrt(Ej0)*cos(theta_j);
                rs(j+1)=rs(j+1)+sqrt(Ej0)*sin(theta_j);
            end
        end
    end
    for j=0:number_of_out-1
        demod_input(j+1,i)=sqrt(rc(j+1)^2+rs(j+1)^2);
    end
    D=nextstate(D+1,source_coded(i)+1); 
end



state_metric=zeros(number_of_states,2);
survivor_state=zeros(number_of_states,depth_of_trellis+1);
for i=1:depth_of_trellis-L
    flag=zeros(1,number_of_states);
    if i<=L+1
        step=2^((L+1-i)*BPH);
    else
        step=1;
    end
    for j=0:step:number_of_states-1
        for m=0:fanout-1
            branch_metric=demod_input(output(j+1,m+1)+1,i);
            if((state_metric(nextstate(j+1,m+1)+1,2)<state_metric(j+1,1)...
                    +branch_metric)|flag(nextstate(j+1,m+1)+1)==0)
                state_metric(nextstate(j+1,m+1)+1,2)=state_metric(j+1,1)+branch_metric;
                survivor_state(nextstate(j+1,m+1)+1,i+1)=j;
                flag(nextstate(j+1,m+1)+1)=1;
            end
        end
    end
    state_metric=state_metric(:,2:-1:1);
end
for i=depth_of_trellis-L+1:depth_of_trellis
    flag=zeros(1,number_of_states);
    last_stop=number_of_states/(2^((i-depth_of_trellis+L-1)*BPH));
    for j=0:last_stop-1
        branch_metric=demod_input(output(j+1,m+1)+1,i);
        if((state_metric(nextstate(j+1,1)+1,2)<state_metric(j+1,1)...
            +branch_metric)|flag(nextstate(j+1,1)+1)==0)
            state_metric(nextstate(j+1,1)+1,2)=state_metric(j+1,1)+branch_metric;
            survivor_state(nextstate(j+1,1)+1,i+1)=j;
            flag(nextstate(j+1,1)+1)=1;
        end
    end
    state_metric=state_metric(:,2:-1:1);
end
state_sequence=zeros(1,depth_of_trellis+1);
for i=1:depth_of_trellis
    state_sequence(1,depth_of_trellis-i+1)=survivor_state((state_sequence(1,depth_of_trellis+2-i)...
        +1),depth_of_trellis-i+2);
end
% decoder_output=zeros(1,BPH*(depth_of_trellis-L));% length(decoder_output)=BPH*2N
decoder_output=zeros(1,depth_of_trellis-L);% depth_of_trellis=length(source_coded)=3002=3000+L,存放3000個(gè)fanout=4進(jìn)制符號(hào)流
for i=1:depth_of_trellis-L
    dec_output_deci=input(state_sequence(1,i)+1,state_sequence(1,i+1)+1);% 輸出的十進(jìn)制符號(hào)寄存器dec_output_deci(其實(shí)是fanout=4進(jìn)制的符號(hào))
                                                                         % 因?yàn)橄旅嬉痪鋎ec_output_bin=deci2change(dec_output_deci,BPH,2)顯然是
                                                                         % 把dec_output_deci轉(zhuǎn)化為 BPH 位二進(jìn)制的,所以dec_output_deci
                                                                         % 應(yīng)該是2^BPH=fanout進(jìn)制而不是十進(jìn)制的
    decoder_output(i)=dec_output_deci; % decoder_output矩陣存放的是G函數(shù)維特比譯碼后得到的3000個(gè)普通域的fanout進(jìn)制的符號(hào)
%     if(BPH~=1)
%         dec_output_bin=deci2change(dec_output_deci,BPH,2);
%         decoder_output((i-1)*BPH+1:i*BPH)=dec_output_bin; %decoder_output是維特比譯碼輸出的二進(jìn)制序列,length(decoder_output)=BPH*2N
%     else
%         decoder_output(1,i)=dec_output_deci;
%     end
end

% 解交織和譯碼之前要先將普通域的fanout=4進(jìn)制的符號(hào)轉(zhuǎn)化為GF(16)域的16進(jìn)制符號(hào),準(zhǔn)備進(jìn)行基于GF(16)域的16進(jìn)制符號(hào)的解交織和RS譯碼
% decoder_output_hex=zeros(1,15*N/7);% 存放普通域的16進(jìn)制符號(hào)
decoder_output_hex=zeros(1,nRS*N/kRS);
% decoder_output_GF=zeros(1,15*N/7);% 存放GF域的16進(jìn)制符號(hào)
decoder_output_GF=zeros(1,nRS*N/kRS);
for i=1:(depth_of_trellis-L)/2   % 即 i=1:15*N/7
    w1=decoder_output(2*(i-1)+1); % decoder_output中有3000個(gè)fanout=4進(jìn)制符號(hào),把兩位4進(jìn)制符號(hào)轉(zhuǎn)化為1個(gè)普通域的16進(jìn)制符號(hào)
    w2=decoder_output(2*i);       % 這段"把兩位4進(jìn)制符號(hào)轉(zhuǎn)化為1個(gè)普通域的16進(jìn)制符號(hào)"的程序通過測(cè)試,正確!
    w=4*w1+1*w2;
    for j=0:2^M-1
        if (w==j)
            decoder_output_hex(i)=j;
        end
    end
end
decoder_output_GF=gf(decoder_output_hex,M);% 把普通域的16進(jìn)制符號(hào)轉(zhuǎn)化為GF(16)域的16進(jìn)制符號(hào),總共1500個(gè)符號(hào),準(zhǔn)備解交織后進(jìn)行RS譯碼

% % ***********************隨機(jī)解交織*********************%
% deint_output_bin=deinterleave(decoder_output,alpha); % BPH*2N個(gè)比特進(jìn)行解交織,2N個(gè)符號(hào)即BPH*2N個(gè)比特
% 
% % ****************************************************%
% % *********************解交織部分(塊交織)(基于比特的解交織)****************%
% % 僅適用于BPH=2且N=1000時(shí)
% deint_output_bin=zeros(1,2*BPH*N);% 解交織輸出deint_output_bin
% B1=zeros(200,200);
% for i=1:2*BPH*N
%     B1(i)=decoder_output(i);
% end
% B=B1.';
% for i=1:2*BPH*N
%     deint_output_bin(i)=B(i);% deint_output_bin為解交織后輸出的二進(jìn)制比特序列,行矢量
% end
% % ************************************************************************%

%**************解交織部分(基于16進(jìn)制GF(16)域符號(hào)(即RS碼元)的解交織)***************%
% deint_output=zeros(N/7,15);
deint_output=zeros(N/kRS,nRS);
deint_output_GF=gf(deint_output,M); % 存放解交織后的GF(16)域的16進(jìn)制符號(hào)矩陣(尺寸為100*15),是 RS 譯碼器的輸入矩陣
% for i=1:15*N/7
for i=1:15*N/7
    deint_output_GF(i)=decoder_output_GF(i);% 對(duì)!必須將解調(diào)后的序列按列存放為100行15列的矩陣
end
%*****************************************************************************%
% % ***********(8,4)譯碼部分***********%
% deco_input=zeros(BPH*2*N,1); % 列矢量
% 
% deco_input(1:BPH*2*N)=deint_output_bin(1:BPH*2*N);
% deco_output=decode(deco_input,8,4,'linear',G)';% (8,4)譯碼器的譯碼輸出(長為BPH*N)個(gè)比特,行矢量)
% 
% % ***********************************%
%**************** RS 碼的譯碼部分(N,K,Dmin)=(15,7,9)的 RS 碼*****************%
deco_output_GF=rsdec(deint_output_GF,15,7);% RS 譯碼輸出矩陣GF(16)域(尺寸100*7)
                                           % 注意這個(gè)deco_output_GF矩陣對(duì)應(yīng)于RS編碼前的source矩陣,矩陣source的第一行存放的是原始信息符號(hào)流source_sequence
                                           % 序列的第1到7個(gè)符號(hào),第二行存放的是source_sequence序列的第8到14個(gè)符號(hào),依此類推.而這個(gè)deco_output_GF矩陣同樣如此,
                                           % 它的每個(gè)位置上的元素與前面source矩陣相同位置上的元素應(yīng)該是一一對(duì)應(yīng)的!
                                           % 下一步應(yīng)該把這個(gè)deco_output_GF矩陣的元素按行取出構(gòu)成一個(gè)1*N的序列,與原始信息符號(hào)流source_sequence進(jìn)行
                                           % 比較,計(jì)算符號(hào)錯(cuò)誤率
%****************************************************************************%

deco_output_GF1=deco_output_GF';% 轉(zhuǎn)置一下變成7行100列矩陣,每列為一個(gè)信息符號(hào)組    
deco_output_GF2=reshape(deco_output_GF1,1,N); % 把譯碼輸出reshape(reshape按列取出)成1行N=700列的GF(16)域符號(hào)矢量,準(zhǔn)備計(jì)算符號(hào)錯(cuò)誤個(gè)數(shù)
 
for i=1:(M/BPH)*15*N/7   % 計(jì)算維特比譯碼后,RS譯碼之前的4進(jìn)制符號(hào)錯(cuò)誤數(shù) 
    if (decoder_output(i)~=codedsequence_quart(i))
        vita_symbol_err(pp,rep)=vita_symbol_err(pp,rep)+1;% 這時(shí)vita_symbol_err是維特比譯碼后,RS譯碼之前的2^BPH=4進(jìn)制符號(hào)錯(cuò)誤數(shù)
    end
end

for i=1:N    % N=700,700個(gè)16進(jìn)制符號(hào),計(jì)算最終RS譯碼后的16進(jìn)制符號(hào)錯(cuò)誤數(shù)
    if(deco_output_GF2(i)~=source_sequence(i))
       num_of_err(pp,rep)=num_of_err(pp,rep)+1; % 這時(shí)num_of_err是2^M=16進(jìn)制符號(hào)錯(cuò)誤數(shù)
   end
end

waitbar(rep/times,WTbar)
end % 與最外層"多少遍rep" for循環(huán)對(duì)應(yīng)的end
close(WTbar)
waitbar(pp/length(Eb_to_Nj_in_dB),WTbarpp)
end % 對(duì)應(yīng)于最最外面pp信干比循環(huán)的end
close(WTbarpp)

PS_all=sum(num_of_err,2)/(N*times) % 計(jì)算出16進(jìn)制信息符號(hào)錯(cuò)誤率,其中N 是一遍的符號(hào)流長度,總共進(jìn)行times遍,共N*times個(gè)符號(hào)
PB_all=2^(M-1)/(2^M-1)*PS_all

symerr_num_biterr=zeros(2*length(Eb_to_Nj_in_dB),times);% symerr_num_biterr矩陣中奇數(shù)行是某一信干比下維特比譯碼后,RS譯碼之前的2^BPH=4進(jìn)制符號(hào)錯(cuò)誤數(shù)
                                                        % 偶數(shù)行是對(duì)應(yīng)的最終RS譯碼后的2^M=16進(jìn)制符號(hào)錯(cuò)誤數(shù)
for i=1:length(Eb_to_Nj_in_dB)
    symerr_num_biterr(2*(i-1)+1,:)=vita_symbol_err(i,:);
    symerr_num_biterr(2*i,:)=num_of_err(i,:);
end