PRO calc_phi,lam
;
; Pre-calculation of the Lorentzian and Gaussian contributions phi
; to scattered light for pre-set parameters b (") and
; distances r (") from the disk center.
; INPUT : lam - working wavelength
; OUTPUT: IDL save file 'phi_val_lam' containing arrays of
; 	  phi(r) for each parameter b.
;
; To be used for fitting the scattering parameters in scatter.pro
; Calls functions PHIL.PRO, PHIG.PRO, CLV.PRO, QROMB.PRO, BESELI.PRO

; Distances r = 0", 930", 931", ... , 990" from disk center
r=findgen(62)+930.		; step 1", limb at 960"
r=shift(r,1)
r(0)=0.
n=n_elements(r)

; Lorentzian parameter b = 10" (fixed)
print,'calculating Lorentz contribution for B=10" and',lam
l10=fltarr(n)
for i=0,n-1 do l10(i)=PHIL(r(i),10,lam)

; Gaussian parameters
b=[3.,5.,10.,15.,20.,40.,50.,60.,70.,80.,90.]

g=fltarr(n,11)
for j=0,10 do begin
print,'calculating Gaussian contribution for b=',b(j)
  for i=0,n-1 do g(i,j)=PHIG(r(i),b(j),lam)
endfor

g3=reform(g(*,0))
g5=reform(g(*,1))
g5(0)=g(0,2)		;correction of phi(0,5")
g3(0:4)=g5(0:4)		;correction of phi(0:4,3")
g10=reform(g(*,2))
g15=reform(g(*,3))
g20=reform(g(*,4))
g40=reform(g(*,5))
g50=reform(g(*,6))
g60=reform(g(*,7))
g70=reform(g(*,8))
g80=reform(g(*,9))
g90=reform(g(*,10))

; Saving
print,'Saving in phi_val_'+strtrim(fix(lam),2)
save,l10,g3,g5,g10,g15,g20,g40,g50,g60,g70,g80,g90, $
	file='phi_val_'+strtrim(fix(lam),2)

END