;
;jeff snider, univeristy of wyoming, march 2006
;
;The measurement of c_twomey is at a state which has a specific volume = specific_volume_meas
;
pro twomey, w_local, tk_local, p_local, r_n_twomey
;
common  constants, gascon, rho_h2o, mw_nh42so4, rho_nh42so4, sigma_0,      $
        sigma_t, dsigmadmolality_nh42so4, dsigmadmolality_nacl, mw_h2o,    $
        mw_nacl, rho_nacl, tkmelt, c_nh42so4, a_nh42so4, c_nacl, a_nacl,   $
        cp_air_o, cp_h2o_o, mw_air, gravity, p0, alpha, beta_local,        $
        gascon_h2o, cw_h2o_o, tk_o, lv_o, ew_o, tk_aerosol, salt_type,     $
        epsilon, c_twomey, k_twomey, aerosol_type, right_tail,             $
        gascon_air, tc_base, tk_base, p_base, h_start, hmax,               $
        r0, dt, rmin, rmax, coef, epsilon_aerosol, specific_volume_meas,   $
        nchan, mw_salt, rho_salt, a_salt, c_salt, specific_volume_base,    $
        dsigmadmolality_salt, rho_insoluble, volume_insoluble_to_soluble,  $
        mixrat_tot_1,r_1,geo_sigma_1,mixrat_tot_2,r_2,geo_sigma_2,   $
        tk_start,p_start,sratio_start,mixrat_tot_3,r_3,geo_sigma_3
;
ew            = ew_jeff(tk_local)
r_v           = epsilon * ew / (p_local - ew)
cp_tot        = cp_air_o + r_v*cp_h2o_o
lv_t          = lv_o + (cw_h2o_o - cp_h2o_o)*(tk_o - tk_local)
diffusivity   = dif_h2o(p_local,tk_local)
conductivity  = cond_air(tk_local-tkmelt)
q1_johnson    = lv_t*gravity/(cp_air_o*gascon_h2o*tk_local^2.) - gravity/(gascon_air*tk_local)
q2_johnson    = p_local/(epsilon*ew) + lv_t^2./(cp_air_o*gascon_h2o*tk_local^2.)
e_johnson     = ew*conductivity*diffusivity*gascon_h2o*tk_local^2. / (rho_h2o*(conductivity*gascon_h2o^2.*tk_local^3.+diffusivity*lv_t^2.*ew))
alpha_johnson = 100.^k_twomey*c_twomey*specific_volume_meas
zeta_johnson  = (0.159/(rho_h2o*q2_johnson))*(q1_johnson/e_johnson)^(3./2.)
r_n_twomey    = (alpha_johnson)^(2./(k_twomey+2.))*(zeta_johnson*w_local^(3./2.) / $
                (k_twomey*beta(k_twomey/2.,1.5)))^(k_twomey/(k_twomey+2.))
;
end