function ext_coeff_expo_ZV_bvbzbm, Z, Vdq, bm, bv, bz, am, av, az, 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.

;parameter inputs are bv, bz, and bm.

; 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.7d-5)*(10.d)     ; the 1000/100 converts from kg/m*s to g/cm*s
kin_visc=dyn_visc/rho_air
g=981.d


;ab=(0.04394+0.06049)/2.d & 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

; c_sigma first

term1=(4.d*g*((ab*kin_visc)^(1./bb)))/(rho_air*(kin_visc^2))
term2=bm+3.-((1.+bv)/bb)
c_sigma=term1*am*(av^(-1.d/bb))*gamma(term2)

; cv

term1=(gamma(6.+bv+bz))/(gamma(6.+bz))
term2=1.+(bv/(6.+bz))
cv=av*term1*term2

; cv2

term1=(4.-bm+bz+((1.+bv)/bb))/bv
cv2=cv^term1

;cz

cz=az*(6.+bz)*(gamma(6.+bz))

; the exponent to vdq is term1

ext_coeff=(c_sigma*cv2/cz)*Z/(Vdq^term1)


return, ext_coeff
end