Merge remote-tracking branch 'origin/release-v4.5.2'
[WRF.git] / var / da / da_ssmi / tb.inc
blobef523f1180bfcbb7b17e14a8306801f09b1c49fb
1 !--------------------------------------------------------------------------
2 !Input :   ifREQ = (1,2,3, or 4 ) for (19.35, 22.235, 37, or 85.5 GHz)
3 !          THETA = incidence angle (deg.)
4 !                                                 (approximate  range)
5 !          P0    = surface pressure (mb)                    (940 -1030)
6 !          WV    = precipitable water (kg/m**2)               (0 - 70)
7 !          HWV   = water vapor density scale height (km)    (0.5 - 3.0)
8 !          TA, GAMMA = "effective" surface air temperature 
9 !                        and lapse rate (K ; km)
10 !               T(z) = Ta - gamma*z              (263 - 303 ; 4.0 - 6.5)
12 !          SST = sea surface temperature (K)     (273 - 303)
13 !          U  = wind speed (m/sec)               (0 - 30)
14 !          ALW = column liquid water (kg/m**2)   (0 - 0.5)
15 !              note: ALW > 0.5 usually indicates precipitation
17 !          ZCLD = effective cloud height (km)    
18 !              (Ta - GAMMA*ZCLD)  should be warmer than 243 K in order
19 !               to be reasonably realistic.
21 !Output:   TBV, TBH = vertically and horizontally polarized brightness
22 !               temperature (K) at frequency given by ifREQ  
24 !----------------------- tbmod_ssmi.f ---  cut here ----------------------
25 ! VERSION:  July 1997
26 ! This has been carefully calibrated against cloud-free pixels over the 
27 ! globe for the month of Nov 87.   It's not perfect, but don't change it 
28 ! unless you come up with 
29 ! a better calibration data set.  One remaining task it to add splines for
30 ! foam between 5 and 12 m/sec
32 !  Reference:
34 !       Petty, 1990: "On the Response of the Special Sensor Microwave/Imager
35 !         to the Marine Environment - Implications for Atmospheric Parameter
36 !         Retrievals"  Ph.D. Dissertation, University of Washington, 291 pp.
37 !         (See below for J.Atmos.Ocean.Tech. articles in press or pending)
39 !   coded in a quick-and-dirty fashion, and without any warranty, by:
40 !   Grant W. Petty
41 !   Earth and Atmospheric Sciences Dept.
42 !   Purdue University
43 !   West Lafayette, in  47907
44 !   USA
46 !   Tel. No. (317) 494-2544
47 !   Internet address : gpetty@rain.atms.purdue.edu
49 !----------------------------------------------------------------------
50       subroutine tb(ifreq,theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld, &
51          tbv,tbh)
52       implicit none
53       integer ifreq
55       real   , intent(in   ) :: theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld
56       real   , intent(  out) :: tbv,tbh
57       real freq(4),ebiasv(4),ebiash(4),cf1(4),cf2(4),cg(4)
59       real :: f,costhet,gx2,sigma,remv,remh,tbup,tbdn,tauatm
60       real :: effangv,effangh,tbdnv,dum,foam,emissv,emissh
61       real :: refv,refh,semv,semh,scatv,scath,tbdnh
63       data freq/19.35,22.235,37.0,85.5/
65 ! empirical bias corrections for surface emissivity
68       data ebiasv/0.0095, 0.005,-0.014, -0.0010/
69       data ebiash/0.004,   0.0,-0.023, 0.043/
72       data cf1/0.0015,0.004,0.0027,0.005/
73       data cg/4.50e-3, 5.200e-3, 5.5e-3, 7.2e-3 /
75       data cf2/0.0023,0.000,0.0002,-0.006/
77 ! 'foam' emissivity
78       real :: fem
79       data fem /1.0/
81       f = freq(ifreq)
82       costhet = cos(theta*0.017453)
84 ! effective surface slope variance
86       gx2 = cg(ifreq)*u
88 ! get upwelling atmospheric brightness temperature
90       call tbatmos(ifreq,theta,p0,wv,hwv,ta,gamma,alw,zcld, &
91                                    tbup,tbdn,tauatm)
93 ! convert transmittance to optical depth
95       sigma = -alog(tauatm)*costhet
97 ! get rough surface emissivity
99       call roughem(ifreq,gx2,sst,theta,remv,remh)
101 ! get effective zenith angles for scattered radiation at surface
103       call effang(ifreq,theta,gx2,sigma,effangv,effangh)
105 ! get effective sky brightness temperatures for scattered radiation
107       call tbatmos(ifreq,effangv,p0,wv,hwv,ta,gamma,alw,zcld, &
108                    dum,tbdnv,dum)
110       call tbatmos(ifreq,effangh,p0,wv,hwv,ta,gamma,alw,zcld, &
111                    dum,tbdnh,dum)
112 ! compute 'foam' coverage
113       foam = cf1(ifreq)*u
115       if (u .gt. 5.0) then
116         foam = foam + cf2(ifreq)*(u-5.0)
117       end if
119 ! compute surface emissivities and reflectivity
120       emissv = foam*fem + (1.0 - foam)*(remv + ebiasv(ifreq))
121       emissh = foam*fem + (1.0 - foam)*(remh + ebiash(ifreq))
122       refv = 1.0 - emissv
123       refh = 1.0 - emissh
124 ! compute surface emission term
125       semv = sst*emissv
126       semh = sst*emissh
127 ! compute surface scattering term
128       scatv = refv*tbdnv
129       scath = refh*tbdnh
130 ! combine to get space-observed brightness temperature
131       tbv = tbup + tauatm*(semv + scatv)
132       tbh = tbup + tauatm*(semh + scath)
135       end subroutine tb