pro ciret2_expo_ret_sub_analyt, dbz, vdq2

common pwr_law_coeff, az, av, am, bz, bv, bm

common micro_params_expo22, iwc,dist,conc,mass_median_length, count_median_length, $
	count_mean_length, mass_mean_length

common micro_params_expo, lambda,N0,rho_air

common micro_params_expo23, v_mass_mean

; dbz, vdq2 assumed to come in mm6/m3 and m/s

vdq=vdq2*100.

; calculate lambda

term1=((1.d)+(bv/((6.d)+bz)))*av
gamma1=gamma((6.d)+bz+bv)
gamma2=gamma((6.d)+bz)

lambda=(vdq/(term1*(gamma1/gamma2)))^((-1.d)/bv)	; vdq should be cm/s and lambda will be returned in 1/cm

; calculate N0

ze=1.e-6*(10.d^(dbz/10.d))	; cm^6/m^3

denom=az*((6.d)+bz)*(gamma((6.d)+bz))
numer=ze*(lambda^((6.d)+bz+(1.d)))

N0=numer/denom		; units should be 1./m^3/cm

;calculate the IWC

iwc=am*N0*bm*(gamma(bm))/(lambda^(bm+1))	; will return g/m^3 if am is cgs, N0 is (1/m^3/cm), lambda is (1/cm^(bm+1))


; calcualte the concentration

conc=N0/lambda

count_mean_length=1./lambda

mass_mean_length=(bm+1.d)/lambda

count_median_length=abs((alog(0.5))/(lambda))

v_mass_mean=av*(1+(bv/bm))*(lambda^(-bv))*(gamma(bm+bv)/gamma(bm))


; now the mass median length


;dm=count_median_length & ratio=0.0 & count=0
;iwc_total=-am*N0*expon_integral_series_7(0.2,bm,(-1.*lambda))
;while ratio lt 0.5 and count lt 1000 do begin

;	dm=dm+0.5e-4

;	iwc_int=-am*N0*expon_integral_series_7(dm,bm,(-1.*lambda))

;	count=count+1
;	ratio=iwc_int/iwc_total

;endwhile

mass_median_length=0.0		;dm

dist=-9999.

end