inv_ebwES_jd

0001 function [xlambda,xYval,xYflag] = inv_ebwES_jd(xq,N,xabr,xguess);
0002 %
0003 global q_ebw_ES N_ebw_ES abr_ebw_ES
0004 %
0005 h = waitbar(0,'Calculates ES Roots');
0006 %
0007 for ix = 1:length(xq);
0008 %
0009     q_ebw_ES = xq(ix);
0010     N_ebw_ES = N;
0011     abr_ebw_ES = xabr(ix);
0012     guess_ES = xguess(ix);
0013 %
0014     try
0015         [xlambda(ix),xYval(ix),xYflag(ix)] = fzero(@F_jd,guess_ES);
0016     catch
0017         xlambda(ix) = NaN;
0018         xYval(ix) = NaN;
0019         xYflag(ix) = 0;
0020     end
0021 %
0022     waitbar(ix/length(xq),h);
0023 end
0024 close(h)
0025 %
0026 % -----------------------------------------------
0027 %
0028 function Y = F_jd(lambda);
0029 %
0030 global q_ebw_ES N_ebw_ES abr_ebw_ES
0031 %
0032 alpha = lambda/abr_ebw_ES^2;
0033 Xlambda = repmat(lambda,[N_ebw_ES,1]);
0034 I = besseli(1:N_ebw_ES,lambda)'; % takes the Bessel function on lambda with orders 1:N
0035 Y = alpha - 2*sum(exp(-Xlambda).*I.*(1:N_ebw_ES)'.^2./(q_ebw_ES^2 - (1:N_ebw_ES)'.^2));
0036 
0037 
0038 
0039 
0040 
0041 
0042 
0043