% IMKLINVT Inverse of imaginary modified Kontorovich-Lebedev transform
%
function [Phi,err0,n]=IMKLinvT(fMKLinv,x,B,phiB0,eps0)
ab=r_jacobi01(200);
n=0; Phi=1; Phi0=0;
while abs(Phi-Phi0)>eps0*abs(Phi)
  n=n+1; Phi0=Phi;
  if n>200
    disp('n exceeds 200')
    break
  end
  xw=gauss(n,ab);
  Phi=0;
  for k=1:n
    xg=B*xw(k,1); wg=xw(k,2);
    f=feval(fMKLinv,xg);
    [Kr,Ki]=macdonald(x,1/2,xg);
    Phi=Phi+wg*Ki*f;
  end
  Phi=B*Phi;
end
err0=.482*x^(-3/4)*phiB0/abs(Phi);