title "non-central chi-square inverse distribution function(cinv)".
data list free/ beta df lambda.
begin data
0.05  1 6.6468
0.95  1 6.6468 
0.20  1 6.6468 
end data.

compute cmean=df+2*lambda.
compute start=cmean.
compute step=start/20.
compute prec=.00000001.

loop i=0 to 1000.
   compute cval=start-step*i.
end loop if ((beta)>ncdf.chisq(cval,df,lambda)).
do if (i=0).
loop i=0 to 1000.
   compute cval=start+step*i.
end loop if ((beta)<ncdf.chisq(cval,df,lambda)).
compute i=i*(-1).
end if.
compute betaa=ncdf.chisq(cval,df,lambda).

*begin loop for interval halving.
compute dir=i/abs(i).
compute beta=beta.
compute newn=cval.
compute lambda=lambda.
compute intv=step*2.
compute betadiff=abs(betaa-beta).
loop if (betadiff>=prec).
  compute  intv=intv*.5.
  compute  newn=newn +dir*intv*.5.
*compute new betaa.
*compute betaa .
  compute betaa=ncdf.chisq(newn,df,lambda).
  compute betadiff=abs(betaa-beta).
  do if (betaa<beta).
    compute  dir=1. 
    else .
    compute dir= -1.
  end if.
end loop .
compute  betaa=ncdf.chisq(newn,df,lambda).
compute cval=newn.
format betaa,cval,lambda (f10.5).
format i,df (f8.0).
*print / beta df  betaa cval lambda i.
*execute.
report format=list/var= 
        beta 'beta' 'desired'(8)
       df 'df'(5)
        betaa 'beta' 'aquired'(10)
       lambda 'lambda' 'noncentral' 'parameter'
        cval 'cval' 'cinv' (10)
       i 'initial' 'loop' (8).