% function gamma=SEgamma_noisy(M,src_var,SNR)
%
% The function returns the value of gamma to be used 
% for SE-CMA by the formula in the paper by Ghosh.

% A. Klein  5/01/98

function gamma=SEgamma_noisy(M,src_var,SNR)

% Declare accuracy threshold
THRESH=1e-10*SNR^2;

sigw=sqrt(src_var*10^(-SNR/10));
a=(2*[1:M]-M-1);
a=a/std(a,1)*sqrt(src_var);

% Get a rough upper/lower bound for the correct value of gamma within +/- 1
gamma2=0;
y=1;
while y>0
  gamma2=gamma2+1;
  y=sum(a.*(Q(a/sigw)-2*Q((sqrt(gamma2)+a)/sigw)));
end;
gamma1=gamma2-1;

%use binary search to find correct gamma
gamma=(gamma1+gamma2)/2;
y=sum(a.*(Q(a/sigw)-2*Q((sqrt(gamma)+a)/sigw)));
while abs(y)>THRESH
  if y>0
    gamma1=gamma;
  else
    gamma2=gamma;
  end;
  gamma=(gamma1+gamma2)/2;
  y=sum(a.*(Q(a/sigw)-2*Q((sqrt(gamma)+a)/sigw)));
end