viterbi_detector.m

小区初搜为GSM系统中的一个关键过程

function [ rx_burst ] = viterbi_detector(SYMBOLS,NEXT,PREVIOUS,START,STOPS,Y,Rhh)
%
% VITERBI_DETECTOR:
%           This matlab code does the actual detection of the
%           received sequence. As indicated by the name the algorithm
%           is the viterbi algorithm, which is a MLSE. At this time 
%           the approch is to use Ungerboecks modified algorithm, and
%           to return hard output only.
%
% SYNTAX:   [ rx_burst ] 
%           = 
%           viterbi_detector(SYMBOLS,NEXT,PREVIOUS,START,STOPS,Y,Rhh)
%
% INPUT:    SYMBOLS:  The table of symbols corresponding the the state-
%                     numbers. Format as made by make_symbols.m
%           NEXT:     A transition table containing the next legal
%                     states, as it is generated by the code make_next.m
%           PREVIOUS: The transition table describing the legal previous
%		      states as generated by make_previous.m
%           START:    The start state of the algorithm.
%           STOPS:    The legal stop states.
%           Y:        Complex baseband representation of the matched 
%                     filtered and down converted received signal, as it
%                     is returned by mf.m
%           Rhh:      The autocorrelation as estimated by mf.m
%
% OUTPUT:   rx_burst: The most likely sequence of symbols.
%
% SUB_FUNC: make_increment
%
% WARNINGS: None.
%
% TEST(S):  Tested with no noise, perfect syncronization, and channel 
%           estimation/filtering. (Refer to viterbi_ill.m)
%
% AUTOR:    Jan H. Mikkelsen / Arne Norre Ekstr鴐
% EMAIL:    [e-mail] / [e-mail]
%
% $Id: viterbi_detector.m,v 1.7 1997/11/18 12:41:26 aneks Exp $

% KNOWLEDGE OF Lh AND M IS NEEDED FOR THE ALGORITHM TO OPERATE
%
 
filename1='D:\MATLAB6p1\work\searchSCH\i_sch600.txt';

filename2='D:\MATLAB6p1\work\searchSCH\q_sch600.txt';

i = textread(filename1,'%s');

q = textread(filename2,'%s');

i = (hex2dec(i)-512)/512;

q = (hex2dec(q)-512)/512;

length_i=length(i);

length_q=length(q);

f=i+j*q; %注意 I Q的顺序反了的话对结 
r=f;
TRAINING = [1, 0, 1, 1, 1, 0, 0, 1, 0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 0, 1, 1];
T_SEQ = T_SEQ_gen(TRAINING);
Lh=4;
OSR=4;
% [Y, Rhh] = mafi_SCH_new(r,Lh,T_SEQ,OSR);
[ SYMBOLS ] = make_symbols(Lh);
[ M , Lh ] = size(SYMBOLS);

% THE NUMBER OF STEPS IN THE VITERBI
%
STEPS=length(Y);   %匹配滤波的输出，为148

% INITIALIZE TABLES (THIS YIELDS A SLIGHT SPEEDUP).
%
METRIC = zeros(M,STEPS);        %M行，148列，可以简化为M行2列
SURVIVOR = zeros(M,STEPS);      %幸存路径

% DETERMINE PRECALCULATABLE PART OF METRIC
%
[ NEXT ] = make_next(SYMBOLS);
[ PREVIOUS ] = make_previous(SYMBOLS);
[ INCREMENT ] = make_increment(SYMBOLS,NEXT,Rhh);
%INCREMENT=make_increment(SYMBOLS,NEXT,Rhh);%incemen，当所有的参数定了的时候，也是一个固定的表
                                           %照资料上的说法，这就是一个M*M的矩阵
% THE FIRST THING TO DO IS TO ROLL INTO THE ALGORITHM BY SPREADING OUT 
% FROM 	THE START STATE TO ALL THE LEGAL STATES. 
%
[ START ] = make_start(Lh,SYMBOLS);
PS=START;
[ PREVIOUS ] = make_previous(SYMBOLS);
[ STOPS ] = make_stops(Lh,SYMBOLS);
[ INCREMENT ] = make_increment(SYMBOLS,NEXT,Rhh);

% NOTE THAT THE START STATE IS REFERRED TO AS STATE TO TIME 0
% AND THAT IT HAS NO METRIC.
%
S=NEXT(START,1);
METRIC(S,1)=real(conj(SYMBOLS(S,1))*Y(1))-INCREMENT(PS,S);
SURVIVOR(S,1)=START;        %survivor幸存路径，要保留每个状态每个时刻的PM，为最后的回溯准备
                            %保存的是，此时的这个状态是由前一时刻的哪个状态转移而来
S=NEXT(START,2);
METRIC(S,1)=real(conj(SYMBOLS(S,1))*Y(1))-INCREMENT(PS,S);%这个乘法也是复数运算
SURVIVOR(S,1)=START;

PREVIOUS_STATES=NEXT(START,:); %next为M×2

% MARK THE NEXT STATES AS REAL. N.B: COMPLEX INDICATES THE POLARITY
% OF THE NEXT STATE, E.G. STATE 2 IS REAL.
%
COMPLEX=0;%这个是一个奇偶标志量，表明是从奇状态还是偶状态开始，不同的状态开始的话，有点不同
       
for N = 2:Lh,
  if COMPLEX,   %COMPLEX这个量怎么指示的还是不清楚
    COMPLEX=0;
  else
    COMPLEX=1;
  end
  STATE_CNTR=0;
  for PS = PREVIOUS_STATES,     %刚开始，previous_states这个状态只有两个,但是这个是一个变量，外循环
    STATE_CNTR=STATE_CNTR+1;    %完了后，每次都重新变化
    S=NEXT(PS,1);
    METRIC(S,N)=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N))-INCREMENT(PS,S);
    SURVIVOR(S,N)=PS;
    USED(STATE_CNTR)=S;
    STATE_CNTR=STATE_CNTR+1;
    S=NEXT(PS,2);
    METRIC(S,N)=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N))-INCREMENT(PS,S);
    SURVIVOR(S,N)=PS;    
    USED(STATE_CNTR)=S;     %总感觉没有将所有的状态都遍历
  end
  PREVIOUS_STATES=USED;      %对开始时的一个初始化，注意它是一个变的矩阵
end
% AT ANY RATE WE WILL HAVE PROCESSED Lh STATES AT THIS TIME
%
PROCESSED=Lh;

% WE WANT AN EQUAL NUMBER OF STATES TO BE REMAINING. THE NEXT LINES ENSURE
% THIS.
%

if ~COMPLEX,               %这个的作用还是不理解，前面做了一次，后面又恢复，将所有的都变成了复数
  COMPLEX=1;               %可能是I（n）的关系，一个是实数，一个为纯虚数
  PROCESSED=PROCESSED+1;   %这个单独来一次起什么作用，为什么不和下面的一块进行
  N=PROCESSED;             %如果Lh取2，现在N＝3
  for S = 2:2:M,           %1到8的话，为偶
    PS=PREVIOUS(S,1);
    M1=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N)-INCREMENT(PS,S));
    PS=PREVIOUS(S,2);
    M2=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N)-INCREMENT(PS,S)); 
    if M1 > M2,
      METRIC(S,N)=M1;
      SURVIVOR(S,N)=PREVIOUS(S,1);
    else
      METRIC(S,N)=M2;
      SURVIVOR(S,N)=PREVIOUS(S,2);    %这个也是遍历状态的一部分，但是为什么不放到一起做
    end                               %具体的理论还不清楚
  end
end

% NOW THAT WE HAVE MADE THE RUN-IN THE REST OF THE METRICS ARE
% CALCULATED IN THE STRAIGHT FORWARD MANNER. OBSERVE THAT ONLY
% THE RELEVANT STATES ARE CALCULATED, THAT IS REAL FOLLOWS COMPLEX
% AND VICE VERSA.
%
N=PROCESSED+1;                %已经将所有的遍历进行完了
while N <= STEPS,             %奇偶分开，理由？是不是因为虚实交替原因。现在N＝4
  for S = 1:2:M-1,
    PS=PREVIOUS(S,1);
    M1=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N)-INCREMENT(PS,S));
    PS=PREVIOUS(S,2);
    M2=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N)-INCREMENT(PS,S)); 
    if M1 > M2,
      METRIC(S,N)=M1;
      SURVIVOR(S,N)=PREVIOUS(S,1);
    else
      METRIC(S,N)=M2;
      SURVIVOR(S,N)=PREVIOUS(S,2);
    end   
  end
  N=N+1;
  for S = 2:2:M,
    PS=PREVIOUS(S,1);
    M1=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N)-INCREMENT(PS,S));       %只返回实部
    PS=PREVIOUS(S,2);
    M2=METRIC(PS,N-1)+real(conj(SYMBOLS(S,1))*Y(N)-INCREMENT(PS,S)); 
    if M1 > M2,
      METRIC(S,N)=M1;
      SURVIVOR(S,N)=PREVIOUS(S,1);
    else
      METRIC(S,N)=M2;
      SURVIVOR(S,N)=PREVIOUS(S,2);
    end   
  end
  N=N+1;
end  


% HAVING CALCULATED THE METRICS, THE MOST PROBABLE STATESEQUENCE IS
% INITIALIZED BY CHOOSING THE HIGHEST METRIC AMONG THE LEGAL STOP 
% STATES.
%
BEST_LEGAL=0;            %判断什么情况下回溯的结果才是合法的
for FINAL = STOPS,       %stops不是唯一的，对Lh＝4的情况来说是，有两个截止状态
  if METRIC(FINAL,STEPS) > BEST_LEGAL,
    S=FINAL;
    BEST_LEGAL=METRIC(FINAL,STEPS);
  end
end

% UNCOMMENT FOR TEST OF METRIC
%
% METRIC
% BEST_LEGAL
% S
% pause

% HAVING FOUND THE FINAL STATE, THE MSK SYMBOL SEQUENCE IS ESTABLISHED
%
IEST(STEPS)=SYMBOLS(S,1);
N=STEPS-1;
while N > 0,
  S=SURVIVOR(S,N+1);
  IEST(N)=SYMBOLS(S,1);
  N=N-1;
end

% THE ESTIMATE IS NOW FOUND FROM THE FORMULA:
% IEST(n)=j*rx_burst(n)*rx_burst(n-1)*IEST(n-1)
% THE FORMULA IS REWRITTEN AS:
% rx_burst(n)=IEST(n)/(j*rx_burst(n-1)*IEST(n-1))
% FOR INITIALIZATION THE FOLLOWING IS USED:
% IEST(0)=1 og rx_burst(0)=1
%
rx_burst(1)=IEST(1)/(j*1*1);
for n = 2:STEPS,
  rx_burst(n)=IEST(n)/(j*rx_burst(n-1)*IEST(n-1)); %这个在DSP中应该怎么体现？结果很有可能是1，也就是说不变，除法没什么意义
end

% rx_burst IS POLAR (-1 AND 1), THIS TRANSFORMS IT TO
% BINARY FORM (0 AND 1).
%
rx_burst=(rx_burst+1)./2;    %完成差分解码

rx_data=[ rx_burst(4:42) , rx_burst(107:145) ];
rx_enc=rx_data;

[rx_block,FLAG_SS,PARITY_CHK] = channel_dec_SCH(rx_enc);
g = [1 0 1 0 1 1 1 0 1 0 1];
d = [rx_block 0 0 0 0 0 0 0 0 0 0];
[q,r] = deconv(d,g);
% ADJUST RESULT TO BINARY REPRESENTATION
L = length(r);
out = abs(r(L-9:L));
for n = 1:length(out),
  if ceil(out(n)/2) ~= floor(out(n)/2)
    out(n) = 1;
  else
    out(n) = 0;
  end
end
out2=out