www.pudn.com > matlab.rar > viterbi.m, change:2015-04-19,size:6573b


function [code,channel_output,survivor_state,cumulated_metric]=viterbi(G,k,right_channel_output) 
%right_channel_output = [  1     1     1     0     0     0     1     0     0     1     0     1     0     0     0     1     1 0 1 1]; 
p = length(right_channel_output); 
for i = 1:p 
    if rem(i,8)==0 
        channel_output(i) = bitxor(right_channel_output(i),1); 
    else 
        channel_output(i) = right_channel_output(i); 
    end 
end 
G=[1 1 1 1 0 0 1;1 0 1 1 0 1 1];  
k = 1; 
 
%VITERBI        卷积码的维特比解码器 
 %[decoder_ouput,survivor_state,cumulated_metric]=viterbi(G,k,channel_output) 
 %               G是一个n*Lk矩阵,该矩阵的每一行确 
 %               定了从移位记错器到第n个输出间的连接, 
 %               是码速率。 
 %               survivor_state是表示通过网络的最佳路径的矩阵。 
 %               量度在另一个函数metric(x,y)中给出,而且可根据 
 %               硬判决和软判决来指定。 
 %               该算法最小化了量度而不是最大化似然 
 n=size(G,1);    %取出矩阵G的一维大小,即得出输出端口 
     %   检查大小 
 if rem(size(G,2),k)~=0      %当G列数不是k的整数倍时 
    error('Size of G and k do not agree')    %发出出错信息 
 end 
 if rem(size(channel_output,2),n)~=0         %当输出量元素个数不是输出端口的整数倍时 
    error('Channel output not of the right size')    %发出出错信息 
 end 
 L=size(G,2)/k;          %得出移位数,即寄存器的个数 
     %       由于L-1个寄存器的状态即可表示出输出状态, 
     %       所以总的状态数number_of_states可由前L-1个 
     %       寄存器的状态组合来确定 
 number_of_states=2^((L-1)*k); 
     %   产生状态转移矩阵、输出矩阵和输入矩阵 
 for j=0:number_of_states-1      %j表示当前寄存器组的状态因为状态是从零 
                                     %开始的,所以循环从0到number_of_states-1 
                                
       for l=0:2^k-1              %l为从k个输入端的信号组成的状态,总的状 
                                      %态数为2^k,所以循环从0到2^k-11 
                                      
       %  nxt_stat完成从当前的状态和输入的矢量得出下寄存器组的一个状态 
       [next_state,memory_contents]=nxt_stat(j,l,L,k); 
             % input数组值是用于记录当前状态到下一个状态所要的输入信号矢量 
       %  input数组的维数:    一维坐标x=j+1指当前状态的值 
       %                      二维坐标y=next_state+1指下一个状态的值 
       %  由于Matlab中数组的下标是从1开始的,而状态值 
       %  是从0开始的,所以以上坐标值为:状态值+1 
       input(j+1,next_state+1)=l;  
       %  branch_output用于记录在状态j下输入l时的输出 
       branch_output=rem(memory_contents*G',2); 
       %  nextstate数组记录了当前状态j下输入l时的下一个状态 
       nextstate(j+1,l+1)=next_state; 
       %  output数组记录了当前状态j下输入l时的输出(十进制) 
       output(j+1,l+1)=bin2deci(branch_output); 
    end 
 end 
 %       state_metric数组用于记录译码过程在每状态时的汉明距离 
 %       state_metric大小为number_of_states2,(:,1)当前 
 %       状态位置的汉明距离,为确定值,而(:,2)为当前状态加输入 
 %       得到的下一个状态汉明距离,为临时值 
 state_metric=zeros(number_of_states,2); 
 %       depth_of_trellis用于记录网格图的深度 
 depth_of_trellis=length(channel_output)/n; 
 %       输出矩阵,每一列为一个输出状态 
 channel_output_matrix=reshape(channel_output,n,depth_of_trellis); 
 %       survivor_state描述译码过程中在网格图中的路径 
 survivor_state=zeros(number_of_states,depth_of_trellis+1); 
 %开始无尾信道输出的解码 
 for i=1:depth_of_trellis    %i指示网格图的深度 
    %    flag矩阵用于记录网格图中的某一列是否被访问过 
    flag=zeros(1,number_of_states); 
    if i<=L 
       step=2^((L-i)*k);      %在网格图的开始处,并不是所有的状态都取 
    else                         %用step来说明这个变化 
       step=1;                    %状态数从1→2→4→...→number_of_states 
    end 
    for j=0:step:number_of_states-1      %j表示寄存器的当前状态 
       for l=0:2^k-1                  %l为当前的输入 
          branch_metric=0;        %用于记录码间距离 
          %       将当前状态下输入状态l时的输出output转为n位二进制,以便 
          %       计算码间距离(说明:数组坐标大小变化同上)。 
          binary_output=deci2bin(output(j+1,l+1),n); 
          %       计算实际的输出码同网格图中此格某种输出的码间距离 
          for ll=1:n 
             branch_metric=branch_metric+metric(channel_output_matrix(ll,i),binary_output(ll)); 
          end 
          %       选择码间距离较小的那条路径 
          %       选择方法: 
          %               当下一个状态没有被访问时就直接赋值,否则,用比它小的将其覆盖 
          if((state_metric(nextstate(j+1,l+1)+1,2)>=state_metric(j+1,1)+branch_metric)|flag(nextstate(j+1,l+1)+1)==0) 
             % 下一个状态的汉明距离(临时值)=当前状态的汉明距离(确定值)+ 码间距离 
             state_metric(nextstate(j+1,l+1)+1,2)=state_metric(j+1,1)+branch_metric; 
             % survivor_state数组的一维坐标为下一个状态值,二维坐标为此状态 
             % 在网格图中的列位置,记录的数值为当前状态,这样就可以从网格图中 
             % 某位置的某个状态得出其对应上一个列位置的状态,从而能很方便的完 
             % 成译码过程。 
             survivor_state(nextstate(j+1,l+1)+1,i+1)=j; 
             flag(nextstate(j+1,l+1)+1)=1;     %指示该状态已被访问过 
          end 
       end 
    end 
    state_metric=state_metric(:,2:-1:1);     %移动state_metric,将临时值移为确定值 
    i = i+1; 
 end 
 %开始尾部信道输出解码 
  
 %   从最佳路径中产生解码 
 %   译码过程可从数组survivor_state的最后一个位置向前逐级译码 
 state_sequence=zeros(1,depth_of_trellis+1); 
 %   survivor_state数组的最后的输出状态肯定是“0” 
  
 %   逐级译码过程 
 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); 
 %  state_sequence(1,i+1) = survivor_state(state_sequence(1,i+1)+1,i+1); 
 end 
 decorder_output_matrix=zeros(k,depth_of_trellis); 
 for i=1:depth_of_trellis 
    %    根据数组input的定义来得出从当前状态到下一个状态的输入信号矢量 
    dec_output_deci=input(state_sequence(1,i)+1,state_sequence(1,i+1)+1); 
    %    转成二进制信号 
    dec_output_bin=deci2bin(dec_output_deci,k); 
    %    将一次译码存入译码输出矩阵decoder_output_matrix相应的位置 
    decoder_output_matrix(:,i)=dec_output_bin(k:-1:1)'; 
 end 
 %       按照一维序列形式重新组织输出 
decoder_output=reshape(decoder_output_matrix,1,k*(depth_of_trellis)); 
 
 %       state_metric为网格图最后一个列位置中“0”状态位置的汉明距 
 %       离,这个值就是整个译码过程中的汉明距离。 
 cumulated_metric=state_metric(1,1); 
 %卷积码的维特比译码函数 
  
i = (p-12)/2; 
for j = 1:i 
  code(j) = decoder_output(j); 
end