% cost_update( col, v, vl) - update FIR and G for Markov time-varying ISI
%
function cost_update( col, v, vl)

global KNOWN;
global DET;
global GAUSS;
global GM;
global ML;
global KALEH;

global trans;
global states;
global FIR;
global FIR_0;
global FIR_save;
global FIR_g;
global FIR_length;
global Memory;
global N;

global G;
global GG;
global qq;
global alpha;
global alpha2;

global B;
global avg_h;
global Rh;
global inv_Rh;
global inv_Rh_d;
global var_noise;
global z;


% disp('cost_update'), col


% don't need cost_update for DET, KALEH or direct ML
%
if( DET | KALEH | ML | col >= N)

  return;

end


c1 = col + 1;

if( KNOWN)

  FIR = FIR_save(:,c1);

% no symbol estimates available yet
%
elseif( isempty(B) )

  if( GAUSS)

    % disp( 'cost_update: isempty(B) and GAUSS');

    % reasonable initialization
    %
     FIR = avg_h(:,c1);
     % G = Rh; -- G already set in cost_init()

    % junk initialization, just for testing
    %
    % original, unstable for high SNR:
    %
    % t = 2*(rand(FIR_length,1) > 0.5) - 1; % random [-1,+1]
    % G = inv( t*(t'/var_noise) + inv_Rh);
    % FIR = G*( t*(z(c1)/var_noise) + inv_Rh*avg_h(:,c1));
    %
    % new stable way, 26 June 2000
    %
    % t = 2*(rand(FIR_length,1) > 0.5) - 1; % random [-1,+1]
    % s2 = var_noise + t'*(Rh*t);
    % G = Rh - (Rh*t)*((t'*Rh)/s2);
    % FIR = avg_h(:,c1) + Rh*(t*((z(c1)-t'*avg_h(:,c1))/s2));

  elseif( GM)

     FIR = alpha * FIR;

  else

    error( 'cost_init type flag not set.');

  end

% have symbol estimates in B
%
else 

  if( GAUSS)

    % disp( 'cost_update: ~isempty(B) and GAUSS');

    % original, unstable for high SNR:
    %
    % t = B(c1,:)';
    % G = inv( t*(t'/var_noise) + inv_Rh);
    % FIR = G*( t*(z(c1)/var_noise) + inv_Rh*avg_h(:,c1));

    % new stable way, 24 June 2000
    %
    t = B(c1,:)';
    s2 = var_noise + t'*(Rh*t);
    G = Rh - (Rh*t)*((t'*Rh)/s2);
    FIR = avg_h(:,c1) + Rh*(t*((z(c1)-t'*avg_h(:,c1))/s2));

  elseif( GM)

    Fv = [1:FIR_length];

    t = B(col,:)';

    if( col == 1) 

      G = t*(t'/var_noise) + (1+alpha2)*inv_Rh;
      q = t*(z(col)/var_noise) + alpha*(inv_Rh*FIR_0);

    else

      G = t*(t'/var_noise) + inv_Rh;
      q = t*(z(col)/var_noise);

      CGk1 = inv_Rh * GG(:,Fv+(col-2)*FIR_length);
      G = G + alpha2*(inv_Rh - CGk1*inv_Rh);
      q = q + alpha*(CGk1*qq(:,col-1));

    end

    G = inv(G);

    GG(:,Fv+(col-1)*FIR_length) = G;
    qq(:,col) = q;

    % set next G and q
    %
    t = B(c1,:)';
    CGk1 = inv_Rh * G;
    q = t*(z(c1)/var_noise) + alpha*(CGk1*q);	% RHS G and q are from previous section
    G = t*(t'/var_noise) + inv_Rh - alpha2*(CGk1*inv_Rh);

    if( c1 < N)
      CGk1 = inv_Rh * GG(:,Fv+c1*FIR_length);
      G = G + alpha2*(inv_Rh - CGk1*inv_Rh);
      q = q + alpha*(CGk1*qq(:,c1+1));
    end

    G = inv(G);
    FIR = G*q;

    % FIR_g(:,c1) = FIR;

% [FIR, FIR_save(:,c1)]

  else

    error( 'cost_update bad type flag');

  end % if GAUSS or GM

end % isempty(B)

trans = FSM_trans( states, FIR);
