function ext_coeff_expo_single_rex, N0, mlambda, am, bm, av, bv, rho_air

; assume inputs are in cgs units

; This function calculates the extinction coeffecient assuming that Qext is 2 for large particles.  The formulation is based
; on the parameterization by Mitchell (1996) in his equation 22.  The cross sectional area is solved for and inserted
; into an expression which integrates over the size distribution.  This extinction coeffecient is consistent with the assumed
; characteristics of the ice crystal as represented by the mass and fall speed power laws.

; begin by choosing which value of the best number is appropriate.  Use a generic area relation.  This should only cause some
; error around the boundaries between coeffecients

dyn_visc=(1.7e-5)*(10.)		; the 1000/100 converts from kg/m*s to g/cm*s
kin_visc=dyn_visc/rho_air
g=981.
;ext_coeff_inc=fltarr(4)

;ab=[0.04394, 0.06049, 0.2072, 1.0865]
;bb=[0.970, 0.831, 0.638, 0.499]
;xx=[10.0,585.,1.e5,1.e8]

ab=(0.04394+0.06049)/2. & bb=(0.970+0.831)/2.	; this make the most sense from a physical perspective
;ab=(0.04394) & bb=(0.970)
;ab=0.06049 &
;bb=0.8


;ab=0.195 & bb=0.6		; this is the best set of coeffecients for 9/26 compared to Fu using Fu's published aspect ratios.

aa=0.22 & ba=1.75	; generic area coeffecients

;find the boundaries of the integrals

dd=0.001	; 10 microns in cm

;for j=0,3 do begin

;X=(2.*am*981.*rho_air*(dd^(bm+2.-ba)))/(aa*(dyn_visc^2))
;	if X le xx[j] and X lt xx[n_elements(xx)-1] and dd lt 0.25 then begin
;	lh_d=dd
;		while X lt xx[j] do begin
;			dd=dd+0.0001
;			X=(2.*am*981.*rho_air*(dd^(bm+2.-ba)))/(aa*(dyn_visc^2))
;		endwhile
;			dd=dd-0.0001 & rh_d=dd

; now integrate from lh_d to rh_d

;term1=exp_integral_soln(bm+2.-((1.+bv)/bb),-mlambda, 0.001, 0.25)
term1=(gamma(bm+3.-((1.+bv)/bb)))/(mlambda^(bm+3.-((1.+bv)/bb)))
term2=4.*g*N0*am*(av^(-1./bb))*((ab*kin_visc)^(1./bb))/(rho_air*(kin_visc^2))

ext_coeff=term1*term2

;	endif

;endfor

;ratio=ext_coeff_inc[0]/ext_coeff_inc[1]
;print, ext_coeff_inc
return, ext_coeff
end