Diffusivity angle script

This is a function used to compute diffusivity angle of infrared spectra for given atmospheric conditions.

Reference: Feng J, Huang Y. Diffusivity factor approximation for spectral outgoing longwave radiation, Journal of Atmospheric Sciences (submitted)

function [diff,alpha]=dfa(opts,Ht,H,gamma,wavenum,ts,mode_geo);

%%%% opts: total optical depth, nadir view

%%%% Ht: scale height (km) of optical depth

%%%% gamma: temperature lapse rate (K/km); poistively defined

%%%% ts: surface temperature (K)

%%%% mode_geo: 0 if plane-parral atmosphere; 1 if curved earth surface

R       =6371.23; %%% Earth radius (km)

hkB     =1.438833;%%% h/kB

%H      =100;     %%% Height of TOA (km)

r       =0.5772156;%% Euler’s constant

alpha   =(hkB*wavenum/ts^2)*gamma*Ht*(exp(hkB*wavenum/ts))/(exp(hkB*wavenum/ts)1);

%alpha  =(hkB*wavenum/ts^2/2+1/ts)*gamma*Ht;

if opts<1.45

if mode_geo==0

        scale=1;

        tmpmin=acos(exp(0.5))/pi*180;

elseif mode_geo==1

        scale=(1H/R)^2;

else

        display(‘mode_geo must be 0 or 1!’)

end

        a=(12*alpha*opts0.5*alpha*opts^2*(r+log(opts)2)+alpha/12*opts^3)*scale(1alpha/3*opts^3);

        b=0.5*alpha*opts^3+alpha*opts;

        c=alpha/4*opts^2+alpha/4*opts^3;

        diff=acos((b-sqrt(b^24*a*c))/2/a)/pi*180;

        if mode_geo==1

        tmp=(b-sqrt(b^24*a*c))/2/a;

        diff=acos(sqrt(1(1tmp^2)*scale))/pi*180;

        end

end

if opts>=1.45

if mode_geo==0

        diff=acos(exp(0.5))/pi*180;

elseif mode_geo==1

       scale=(1H/R)^2;

       diff=acos(sqrt(1(1-exp(0.5)^2)*scale))/pi*180;

end

end