da_tb_adj.inc
References to this file elsewhere.
1 subroutine da_tb_adj(ifreq,theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld, &
2 ! tbv,tbh, &
3 ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
4 ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh )
5
6 !-----------------------------------------------------------------------
7 ! Purpose: TBD
8 ! Output : ADJ_p0, ADJ_ta, ADJ_gamma, ADJ_sst, ADJ_wv, ADJ_hwv, ADJ_u
9 ! ADJ_alw, ADJ_zcld
10 ! Input : ADJ_tbv, ADJ_tbh, tbv, tbh
11 !-----------------------------------------------------------------------
12
13 implicit none
14
15 integer, intent(in ) :: ifreq
16 real , intent(in ) :: theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld
17 real , intent(inout) :: ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
18 ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld
19 real , intent(in ) :: ADJ_tbv,ADJ_tbh
20 ! real , intent(in ) :: tbv,tbh
21
22 real :: freq(4),ebiasv(4),ebiash(4),cf1(4),cf2(4),cg(4)
23
24 real :: f,costhet,gx2,tbup,tbdn,tauatm,sigma,remv,remh
25 real :: effangv,effangh,tbdnv,foam,foam_save,emissv,emissh
26 real :: dum
27 real :: refv,refh,semv,semh,scatv,scath,tbdnh
28 real :: ADJ_gx2,ADJ_tbup,ADJ_tbdn,ADJ_tauatm,ADJ_sigma
29 real :: ADJ_tremv,ADJ_remh,ADJ_effangv,ADJ_effangh
30 real :: ADJ_tbdnh,ADJ_dum,ADJ_foam,ADJ_emissv
31 real :: ADJ_emissh,ADJ_refv,ADJ_refh,ADJ_semv
32 real :: ADJ_semh,ADJ_scatv,ADJ_scath,ADJ_remv,ADJ_tbdnv
33 real :: ADJ_theta
34
35 real :: fem
36
37 data freq/19.35,22.235,37.0,85.5/
38
39 ! empirical bias corrections for surface emissivity
40
41 data ebiasv/0.0095, 0.005,-0.014, -0.0010/
42 data ebiash/0.004, 0.0,-0.023, 0.043/
43
44
45 data cf1/0.0015,0.004,0.0027,0.005/
46 data cg/4.50e-3, 5.200e-3, 5.5e-3, 7.2e-3 /
47
48 data cf2/0.0023,0.000,0.0002,-0.006/
49
50 ! 'foam' emissivity
51 data fem /1.0/
52
53 f=0.;costhet=0.;gx2=0.;tbup=0.;tbdn=0.;tauatm=0.
54 sigma=0.;remv=0.;remh=0.;effangv=0.;effangh=0.
55 tbdnv=0.;foam=0.;foam_save=0.;emissv=0.;emissh=0.
56 dum=0.;refv=0.;refh=0.;semv=0.;semh=0.;scatv=0.
57 scath=0.;tbdnh=0.;ADJ_gx2=0.;ADJ_tbup=0.;ADJ_tbdn=0.
58 ADJ_tauatm=0.;ADJ_sigma=0.;ADJ_tremv=0.;ADJ_remh=0.
59 ADJ_effangv=0.;ADJ_effangh=0.;ADJ_tbdnh=0.;ADJ_dum=0.
60 ADJ_foam=0.;ADJ_emissv=0.;ADJ_emissh=0.;ADJ_refv=0.
61 ADJ_refh=0.;ADJ_semv=0.;ADJ_semh=0.;ADJ_scatv=0.
62 ADJ_scath=0.;ADJ_remv=0.;ADJ_tbdnv=0.
63 ADJ_theta =0.
64
65
66 ! write (unit=stdout,fmt=*) 'ifreq',ifreq,theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld, &
67 ! tbv,tbh, &
68 ! ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
69 ! ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
70
71 f = freq(ifreq)
72 costhet = cos(theta*0.017453)
73
74 ! effective surface slope variance
75
76 gx2 = cg(ifreq)* u
77
78 ! get upwelling atmospheric brightness temperature
79
80 call tbatmos(ifreq,theta,p0,wv,hwv,ta,gamma,alw,zcld, &
81 tbup,tbdn,tauatm)
82
83 ! convert transmittance to optical depth
84
85 sigma = -alog(tauatm)*costhet
86
87 ! get rough surface emissivity
88
89 call roughem(ifreq,gx2,sst,theta,remv,remh)
90
91 ! get effective zenith angles for scattered radiation at surface
92
93 call effang(ifreq,theta,gx2,sigma,effangv,effangh)
94
95 ! get effective sky brightness temperatures for scattered radiation
96
97 call tbatmos(ifreq,effangv,p0,wv,hwv,ta,gamma,alw,zcld, &
98 dum,tbdnv,dum)
99
100 call tbatmos(ifreq,effangh,p0,wv,hwv,ta,gamma,alw,zcld, &
101 dum,tbdnh,dum)
102
103 ! compute 'foam' coverage
104
105 foam = cf1(ifreq)* u
106
107 if (u .gt. 5.0) then
108 foam_save = foam
109 foam = foam + cf2(ifreq)*( u-5.0)
110 end if
111
112 ! compute surface emissivities and reflectivity
113
114 emissv = foam*fem + (1.0 - foam)*(remv + ebiasv(ifreq))
115 emissh = foam*fem + (1.0 - foam)*(remh + ebiash(ifreq))
116 refv = 1.0 - emissv
117 refh = 1.0 - emissh
118
119 ! compute surface emission term
120
121 semv = sst*emissv
122 semh = sst*emissh
123
124 ! compute surface scattering term
125
126 scatv = refv*tbdnv
127 scath = refh*tbdnh
128
129 ! combine to get space-observed brightness temperature
130
131 ! tbv = tbup + tauatm*(semv + scatv)
132 ! tbh = tbup + tauatm*(semh + scath)
133
134
135 ! start
136 ! write (unit=stdout,fmt=*) 'ifreq 1',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
137 ! ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
138
139
140 ADJ_tbup = ADJ_tbh !!! first
141 ADJ_tauatm = ADJ_tbh*(semh + scath) !!! first
142 ADJ_semh = tauatm*ADJ_tbh !!! first
143 ADJ_scath = tauatm*ADJ_tbh !!! first
144
145 ADJ_tbup = ADJ_tbv + ADJ_tbup
146 ADJ_tauatm = ADJ_tbv*(semv + scatv) + ADJ_tauatm
147 ADJ_semv = tauatm*ADJ_tbv !!! first
148 ADJ_scatv = tauatm*ADJ_tbv !!! first
149
150 ADJ_refh = ADJ_scath*tbdnh !!! first
151 ADJ_tbdnh = refh*ADJ_scath !!! first
152 ADJ_refv = ADJ_scatv*tbdnv !!! first
153 ADJ_tbdnv = refv*ADJ_scatv !!! first
154 ADJ_sst = ADJ_semh*emissh + ADJ_sst
155 ADJ_emissh = sst*ADJ_semh !!! first
156 ADJ_sst = ADJ_semv*emissv + ADJ_sst
157 ADJ_emissv = sst*ADJ_semv !!! first
158
159 ADJ_emissh = - ADJ_refh + ADJ_emissh
160 ADJ_emissv = - ADJ_refv + ADJ_emissv
161
162 ADJ_foam = ADJ_emissh*fem !!! first
163 ADJ_foam = - ADJ_emissh*(remh + ebiash(ifreq)) + ADJ_foam
164 ADJ_remh = (1.0 - foam)*ADJ_emissh !!! first
165 ADJ_foam = ADJ_emissv*fem + ADJ_foam
166 ADJ_foam = - ADJ_emissv*(remv + ebiasv(ifreq)) + ADJ_foam
167 ADJ_remv = (1.0 - foam)*ADJ_emissv !!! first
168
169 if (u .gt. 5.0) then
170 ADJ_u = cf2(ifreq)*ADJ_foam + ADJ_u
171 foam=foam_save
172 end if
173 ADJ_u = cf1(ifreq)*ADJ_foam + ADJ_u
174
175 ADJ_dum = 0.0
176 dum = 0.0
177 ADJ_effangh = 0.0
178 call da_tbatmos_adj(ifreq,effangh,p0,wv,hwv,ta,gamma,alw, &
179 zcld,dum,tbdnh,dum, &
180 ADJ_effangh,ADJ_p0,ADJ_wv,ADJ_hwv, &
181 ADJ_ta,ADJ_gamma,ADJ_alw,ADJ_zcld, &
182 ADJ_dum,ADJ_tbdnh,ADJ_dum )
183 dum = 0.0
184 ADJ_dum = 0.0
185
186 ! write (unit=stdout,fmt=*) 'ifreq 2',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
187 ! ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
188
189 ADJ_effangv = 0.0
190 call da_tbatmos_adj(ifreq,effangv,p0,wv,hwv,ta,gamma,alw, &
191 zcld,dum,tbdnv,dum, &
192 ADJ_effangv,ADJ_p0,ADJ_wv,ADJ_hwv, &
193 ADJ_ta,ADJ_gamma,ADJ_alw,ADJ_zcld, &
194 ADJ_dum,ADJ_tbdnv,ADJ_dum )
195
196 ADJ_gx2=0.0
197 ADJ_sigma=0.0
198 ! write (unit=stdout,fmt=*) 'ifreq 3',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
199 ! ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
200
201 call da_effang_adj(ifreq,theta,gx2,sigma,effangv,effangh, &
202 ADJ_gx2,ADJ_sigma,ADJ_effangv,ADJ_effangh )
203
204 call da_roughem_adj(ifreq,gx2,sst,theta,remv,remh, &
205 ADJ_gx2,ADJ_sst,ADJ_remv,ADJ_remh )
206
207 ADJ_tauatm = - costhet*ADJ_sigma/tauatm + ADJ_tauatm
208
209
210 ! write (unit=stdout,fmt=*) 'ifreq 4',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
211 ! ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
212
213 call da_tbatmos_adj(ifreq,theta,p0,wv,hwv,ta,gamma,alw,zcld, &
214 tbup,tbdn,tauatm, &
215 ADJ_theta,ADJ_p0,ADJ_wv,ADJ_hwv,ADJ_ta,ADJ_gamma, &
216 ADJ_alw,ADJ_zcld,ADJ_tbup,ADJ_tbdn, &
217 ADJ_tauatm )
218
219 ADJ_theta=0.0 ! first
220
221 ADJ_u = cg(ifreq)*ADJ_gx2 + ADJ_u
222
223 ! write (unit=stdout,fmt=*) 'ifreq 5',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv, &
224 ! ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
225 ! end
226
227 end subroutine da_tb_adj
228
229