% [dmin,Ne,N1,N2,N3,Nb,Nb1,Nb2,Nb3,LD,LD1,LD2,LD3,L] = ...
%   dmin_main(nstates,ninput,trellis_transition,trellis_output)
% Function finding the minimum distance (dmin), the number of
% neighbors at minimum distance (Ne), the number of neighbors
% corresponding to next to minimum distance paths (N1,N2,N3), the
% number of bit errors for the minimum distance paths (Nb), the
% number of bit errors corresponding to next to minimum distance
% paths (Nb1,Nb2,Nb3), the maximum path length for a minimum
% distance error event (LD), and the maximum path lengths
% corresponding to next to minimum distance error events
% (LD1,LD2,LD3). 
%
% The arguments are nstates (the number of states), ninput (the
% number of input bits), trellis_transition (a matrix describing
% the trellis state transitions) and trellis_output (a matrix
% describing the trellis outputs).
% The number of columns of the matrices equal the number of inputs,
% and the number of rows equal the number of trellis states. 
%
% e.g. the 4-state convolutional code of fig. 8.6 has:
% nstates = 4
% ninput = 1
% trellis_transition=[1 2; 3 4; 1 2; 3 4]
% trellis_output=[0 3; 2 1; 3 0; 1 2]
%
% George Ginis, May 2001
%
% Many thanks to Hichan Moon for the several helpful discussions, and
% the sharing of ideas contributing to the improvement of the algorithm. 
% 

function [dmin,Ne,N1,N2,N3,Nb,Nb1,Nb2,Nb3,LD,LD1,LD2,LD3,L] = dmin_main(nstates,ninput,trellis_transition,trellis_output)

% The list L contains information about the variables 
% dmin, Ne, N1, N2, N3, Nb, Nb1, Nb2, Nb3, LD, LD1, LD2, LD3 
% corresponding to each state.   
  
% upper bound on minimum distance
% (must be initialized to a finite large value)
dmin_ub=1e6;

% neighbor depth indicates the distance depth above minimum distance
% e.g. ndepth = 3 says that only paths up to dmin+3 are of interest
ndepth = 3;

% initialize L
L=dmin_init(nstates,ninput,trellis_transition,trellis_output,ndepth,dmin_ub);

% initialize search length
slen = 1;

% initialize vector keeping information about state 1 
L1=zeros(1,1+3*(ndepth+1));
L1(1)=dmin_ub;  % initialize distance to the upper bound value 

% define maximum search length
max_slen = 200; 

% search continues (trellis "propagates") until either the maximum
% search length has been reached, or the distances of all states
% have exceeded the distance of the (ndepth)-th path, which implies
% that all paths of interest have already been discovered.  
cont = 1;    % acts like a flag 
while (cont == 1 & slen <= max_slen) 
  
  % increment search length
  slen = slen + 1;
  
  % "propagate" paths one step forward 
  % returns updated information about states and flag indicating
  % whether all paths of interest have already been discovered
  [L,cont]=...
    dmin_prop(L,nstates,ninput,...
	      trellis_transition,trellis_output,ndepth,dmin_ub,L1);
   
  % save information about state 1 from last search
  L1old = L1;
  
  % update number of neighbors information
  L1(1:ndepth+2) = polyadd( L1old(1:ndepth+2), L(1,1:ndepth+2) );

  % update number of bits in error information
  poly_bits_add = polyadd(...
    [L1old(1) L1old(ndepth+3:2*ndepth+3)], [L(1,1) L(1,ndepth+3:2*ndepth+3)] );
  L1(ndepth+3:2*ndepth+3) = poly_bits_add(2:ndepth+2);      
  
  % update path length quantities
  if (L1old(1) > L1(1)) 
    % if path of lower minimum distance is found, path length
    % information must be updated
    L1(2*ndepth+4:3*ndepth+4) = slen;
    
  elseif (L1old(1) == L1(1)) 
    
    for k=0:ndepth
  
    if (L1old(2+k) < L1(2+k)) 
      % if new paths were discovered in the last search, path
      % length information is updated 
      L1(2*ndepth+4+k) = slen;
    end
  
  end
 
  end
  
end

% define quantities of interest
dmin=L1(1);
Ne=L1(2);
N1=L1(3);
N2=L1(4);
N3=L1(5);
Nb=L1(6);
Nb1=L1(7);
Nb2=L1(8);
Nb3=L1(9);
LD=L1(10);
LD1=L1(11);
LD2=L1(12);
LD3=L1(13);