function [g_draw] = draw_g(e,g_prior,ind_e,r_m,r_s,i_init)

% Draw g from Stochastic Volatility Model using chi-squared mixtureapprox;
%   e = data;
%   g = standard deviation of innovation for ln(volatilty);
%   ind_e = mixture indicators;
%   r_m = mixture means
%   r_s = mixture sd
%   sigma_e_min = minimum value of sigma_e
%   i_init = 1 if vague prior used; = 0 is ln(sig2) = 0 is used;
%;

big = 1.0e+6;  % large number     
T=size(e,1);
c=0.001;
lnres2=log(e.^2 + c);      % c = 0.001 factor from ksv, restud(1998), page 370) 
mean_t = ind_e*r_m;
sigma_t = ind_e*r_s;
var_t = sigma_t.^2;

 % Kalman Filtering 
 ye=lnres2-mean_t;
 p1t=zeros(T+1,1);
 p2t=zeros(T+1,1);
 x1t=zeros(T+1,1);
 x2t=zeros(T+1,1);
 x3_draw=zeros(T+1,1);
 
 % Compute log-likelihood values for each g value
 n_g = size(g_prior,1);
 g_values = g_prior(:,1);
 p_values = g_prior(:,2);
 llf_vec = NaN*zeros(n_g,1);
 
 for ig = 1:n_g;
 	gam = g_values(ig)^2; 
  % -- Compute Covariance Matrix  -- 
  x1=0;
  if i_init == 1;
    p1=big;
  elseif i_init == 0;
    p1 = 0;
  end;
  llf = 0;
  for t=1:T;
 	 x2=x1;
 	 p2=p1+gam;
 	 if isnan(ye(t)) == 1;
 	  x1 = x2;
 	  p1 = p2;
   else;
    e = ye(t)-x2;
 	  h=p2+var_t(t);
 	  k=p2/h;
 	  p1=p2-k*p2;
 	  x1=x2+k*e;
 	  llf = llf - 0.5*(log(h) + (e*e)/h); 
 	 end;
  end;
  llf_vec(ig) = llf;
 end;
 llf_vec = llf_vec - max(llf_vec);   % normalize so largest value is 0 ;
 lf_vec = exp(llf_vec);
 lf_marg = lf_vec'*p_values;
 g_post = (lf_vec.*p_values)./lf_marg;
 cum_prob = cumsum(g_post);
 u = rand(1,1);
 aa = u < cum_prob;
 bb = n_g+1-sum(aa);
 g_draw = g_values(bb);
 
 %{
 g_post
 cum_prob
 u
 bb
 g_draw
 %}
 
end