da_effang_adj.inc
References to this file elsewhere.
1 subroutine da_effang_adj(ifreq,theta,gx2,sigma,effangv,effangh, &
2 ADJ_gx2,ADJ_sigma,ADJ_effangv,ADJ_effangh )
3
4 implicit none
5
6 !-----------------------------------------------------------------
7 ! Purpose: Calculate the effective zenith angle of reflected microwave
8 ! radiation at SSM/I frequencies for vertical and horizontal polarization
9 !
10 ! Output: ADJ_gx2, ADJ_sigma
11 ! Input: ADJ_effangv,ADJ_effangh,effangv,effangh
12 !-----------------------------------------------------------------
13
14 integer, intent(in) :: ifreq
15 real, intent(in) :: theta,gx2,sigma
16 real, intent(inout) :: ADJ_gx2,ADJ_sigma
17 real, intent(inout) :: ADJ_effangv,ADJ_effangh
18 real, intent(out) :: effangv,effangh
19
20 real c19v,c19h,c22v,c22h,c37v,c37h,c85v,c85h
21 real s19v(6),s19h(6),s22v(6),s22h(6), &
22 s37v(6),s37h(6),s85v(6),s85h(6)
23
24 real :: alnsig,gg,ggg,xd,xx,xd_save,xx_save
25 real :: z1,z2,z3,z4,z5,z6
26 real :: y,dth,angh,angv,alnthv,alnthh,alnthv_save,alnthh_save
27 real :: ADJ_alnsig,ADJ_gg,ADJ_ggg,ADJ_xd
28 real :: ADJ_z1,ADJ_z2,ADJ_z3,ADJ_z4,ADJ_z5,ADJ_z6,ADJ_alnthv
29 real :: ADJ_y,ADJ_dth,ADJ_angh,ADJ_angv,ADJ_xx,ADJ_alnthh
30
31 data c19v,c19h,c22v,c22h,c37v,c37h,c85v,c85h &
32 /-.5108,.5306,-.5108,.5306,-.6931,.1823,-.9163,.3000/
33 data s19v /.225E+2,.698E+2,-.238E+2,-.648E+1,.402E+0,.262E+1/
34 data s19h /.743E+1,.507E+2,-.206E+2,-.191E+1,.648E-1,.291E+1/
35 data s22v /.249E+2,.701E+2,-.240E+2,-.714E+1,.405E+0,.256E+1/
36 data s22h /.498E+1,.442E+2,-.190E+2,-.129E+1,.803E-2,.277E+1/
37 data s37v /.215E+2,.573E+2,-.211E+2,-.670E+1,.443E+0,.253E+1/
38 data s37h /.869E+1,.571E+2,-.257E+2,-.302E+1,.237E+0,.386E+1/
39 data s85v /.116E+2,.263E+2,-.101E+2,-.358E+1,.270E+0,.175E+1/
40 data s85h /.736E+1,.568E+2,-.254E+2,-.248E+1,.196E+0,.387E+1/
41
42 if (trace_use) call da_trace_entry("da_effang_adj")
43
44 alnsig=0.0;gg=0.0;ggg=0.0;xd=0.0;xx=0.0;xd_save=0.0;xx_save=0.0
45 z1=0.0;z2=0.0;z3=0.0;z4=0.0;z5=0.0;z6=0.0;y=0.0;dth=0.0;angh=0.0
46 angv=0.0;alnthv=0.0;alnthh=0.0;alnthv_save=0.0;alnthh_save=0.0
47 ADJ_alnsig=0.0;ADJ_gg=0.0;ADJ_ggg=0.0;ADJ_xd=0.0
48 ADJ_z1=0.0;ADJ_z2=0.0;ADJ_z3=0.0;ADJ_z4=0.0;ADJ_z5=0.0
49 ADJ_z6=0.0;ADJ_alnthv=0.0;ADJ_y=0.0;ADJ_dth=0.0;ADJ_angh=0.0
50 ADJ_angv=0.0;ADJ_xx=0.0;ADJ_alnthh=0.0
51
52
53 if (gx2 .le. 0.0 .or. sigma .le. 0.0) then
54 effangv = theta
55 effangh = theta
56 return
57 end if
58 alnsig = alog(sigma)
59 gg = gx2*gx2
60 ggg = gg*gx2
61
62 if (ifreq .eq. 1) then
63 xd = alnsig - c19v
64 xx = xd*xd
65 z1 = xx*ggg
66 z2 = xd*ggg
67 z3 = xd*gg
68 z4 = xx*gg
69 z5 = xx*gx2
70 z6 = xd*gx2
71 alnthv = s19v(1)*z1 + s19v(2)*z2 + s19v(3)*z3 + &
72 s19v(4)*z4 + s19v(5)*z5 + s19v(6)*z6
73 alnthv_save = alnthv
74 alnthv = alnthv + 3.611
75
76 xd_save = xd
77 xx_save = xx
78
79 xd = alnsig - c19h
80 xx = xd*xd
81 z1 = xx*ggg
82 z2 = xd*ggg
83 z3 = xd*gg
84 z4 = xx*gg
85 z5 = xx*gx2
86 z6 = xd*gx2
87
88 alnthh = s19h(1)*z1 + s19h(2)*z2 + s19h(3)*z3 + &
89 s19h(4)*z4 + s19h(5)*z5 + s19h(6)*z6
90 alnthh_save = alnthh
91 alnthh = alnthh + 3.611
92
93 else if (ifreq .eq. 2) then
94 xd = alnsig - c22v
95 xx = xd*xd
96 z1 = xx*ggg
97 z2 = xd*ggg
98 z3 = xd*gg
99 z4 = xx*gg
100 z5 = xx*gx2
101 z6 = xd*gx2
102 alnthv = s22v(1)*z1 + s22v(2)*z2 + s22v(3)*z3 + &
103 s22v(4)*z4 + s22v(5)*z5 + s22v(6)*z6
104
105 alnthv_save = alnthv
106 alnthv = alnthv + 3.611
107
108 xd_save = xd
109 xx_save = xx
110
111 xd = alnsig - c22h
112 xx = xd*xd
113 z1 = xx*ggg
114 z2 = xd*ggg
115 z3 = xd*gg
116 z4 = xx*gg
117 z5 = xx*gx2
118 z6 = xd*gx2
119 alnthh = s22h(1)*z1 + s22h(2)*z2 + s22h(3)*z3 + &
120 s22h(4)*z4 + s22h(5)*z5 + s22h(6)*z6
121
122 alnthh_save = alnthh
123 alnthh = alnthh + 3.611
124
125 else if (ifreq .eq. 3) then
126 xd = alnsig - c37v
127 xx = xd*xd
128 z1 = xx*ggg
129 z2 = xd*ggg
130 z3 = xd*gg
131 z4 = xx*gg
132 z5 = xx*gx2
133 z6 = xd*gx2
134 alnthv = s37v(1)*z1 + s37v(2)*z2 + s37v(3)*z3 + &
135 s37v(4)*z4 + s37v(5)*z5 + s37v(6)*z6
136
137 alnthv_save = alnthv
138 alnthv = alnthv + 3.611
139
140 xd_save = xd
141 xx_save = xx
142
143 xd = alnsig - c37h
144 xx = xd*xd
145 z1 = xx*ggg
146 z2 = xd*ggg
147 z3 = xd*gg
148 z4 = xx*gg
149 z5 = xx*gx2
150 z6 = xd*gx2
151 alnthh = s37h(1)*z1 + s37h(2)*z2 + s37h(3)*z3 + &
152 s37h(4)*z4 + s37h(5)*z5 + s37h(6)*z6
153
154 alnthh_save = alnthh
155 alnthh = alnthh + 3.611
156
157 else if (ifreq .eq. 4) then
158 xd = alnsig - c85v
159 xx = xd*xd
160 z1 = xx*ggg
161 z2 = xd*ggg
162 z3 = xd*gg
163 z4 = xx*gg
164 z5 = xx*gx2
165 z6 = xd*gx2
166 alnthv = s85v(1)*z1 + s85v(2)*z2 + s85v(3)*z3 + &
167 s85v(4)*z4 + s85v(5)*z5 + s85v(6)*z6
168
169 alnthv_save = alnthv
170 alnthv = alnthv + 3.611
171
172 xd_save = xd
173 xx_save = xx
174
175 xd = alnsig - c85h
176 xx = xd*xd
177 z1 = xx*ggg
178 z2 = xd*ggg
179 z3 = xd*gg
180 z4 = xx*gg
181 z5 = xx*gx2
182 z6 = xd*gx2
183 alnthh = s85h(1)*z1 + s85h(2)*z2 + s85h(3)*z3 + &
184 s85h(4)*z4 + s85h(5)*z5 + s85h(6)*z6
185
186 alnthh_save = alnthh
187 alnthh = alnthh + 3.611
188 end if
189
190 angv = 90.0 - exp(alnthv)
191 angh = 90.0 - exp(alnthh)
192 y = 1.0 - 28.0*gx2
193 if (y .lt. 0.0) then
194 y = 0.0
195 end if
196
197 dth = (theta - 53.0)*y
198 effangv = angv + dth
199 effangh = angh + dth
200
201 ! start
202
203 if (gx2 .le. 0.0 .or. sigma .le. 0.0) then
204 ADJ_effangv = 0.0
205 ADJ_effangh = 0.0
206 return
207 end if
208
209 ADJ_angh = ADJ_effangh
210 ADJ_dth = ADJ_effangh
211
212 ADJ_angv = ADJ_effangv
213 ADJ_dth = ADJ_effangv + ADJ_dth
214
215 ADJ_y = (theta - 53.0)*ADJ_dth
216
217 if (y .lt. 0.0) then
218 ADJ_y = 0.0
219 end if
220
221 ADJ_gx2 = - 28.0*ADJ_y + ADJ_gx2
222
223 ADJ_alnthh = - ADJ_angh*exp(alnthh)
224 ADJ_alnthv = - ADJ_angv*exp(alnthv)
225
226 if (ifreq .eq. 1) then
227
228 alnthh = alnthh_save
229
230 ADJ_z1 = s19h(1)*ADJ_alnthh
231 ADJ_z2 = s19h(2)*ADJ_alnthh
232 ADJ_z3 = s19h(3)*ADJ_alnthh
233 ADJ_z4 = s19h(4)*ADJ_alnthh
234 ADJ_z5 = s19h(5)*ADJ_alnthh
235 ADJ_z6 = s19h(6)*ADJ_alnthh
236
237 ADJ_xd = ADJ_z6*gx2
238 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
239 ADJ_xx = ADJ_z5*gx2
240 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
241 ADJ_xx = ADJ_z4*gg + ADJ_xx
242 ADJ_gg = xx*ADJ_z4
243 ADJ_xd = ADJ_z3*gg + ADJ_xd
244 ADJ_gg = xd*ADJ_z3 + ADJ_gg
245 ADJ_xd = ADJ_z2*ggg + ADJ_xd
246 ADJ_ggg = xd*ADJ_z2
247 ADJ_xx = ADJ_z1*ggg + ADJ_xx
248 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
249
250 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
251 ADJ_alnsig = ADJ_xd
252
253 ! vertical
254 xd = xd_save
255 xx = xx_save
256
257 alnthv = alnthv_save
258
259 ADJ_z1 = s19v(1)*ADJ_alnthv
260 ADJ_z2 = s19v(2)*ADJ_alnthv
261 ADJ_z3 = s19v(3)*ADJ_alnthv
262 ADJ_z4 = s19v(4)*ADJ_alnthv
263 ADJ_z5 = s19v(5)*ADJ_alnthv
264 ADJ_z6 = s19v(6)*ADJ_alnthv
265
266 ADJ_xd = ADJ_z6*gx2 !+ ADJ_xd
267 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
268 ADJ_xx = ADJ_z5*gx2 !+ ADJ_xx
269 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
270 ADJ_xx = ADJ_z4*gg + ADJ_xx
271 ADJ_gg = xx*ADJ_z4 + ADJ_gg
272 ADJ_xd = ADJ_z3*gg + ADJ_xd
273 ADJ_gg = xd*ADJ_z3 + ADJ_gg
274 ADJ_xd = ADJ_z2*ggg + ADJ_xd
275 ADJ_ggg = xd*ADJ_z2 + ADJ_ggg
276 ADJ_xx = ADJ_z1*ggg + ADJ_xx
277 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
278
279 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
280
281 ADJ_alnsig = ADJ_xd + ADJ_alnsig
282
283 else if (ifreq .eq. 2) then
284
285 ADJ_z1 = s22h(1)*ADJ_alnthh
286 ADJ_z2 = s22h(2)*ADJ_alnthh
287 ADJ_z3 = s22h(3)*ADJ_alnthh
288 ADJ_z4 = s22h(4)*ADJ_alnthh
289 ADJ_z5 = s22h(5)*ADJ_alnthh
290 ADJ_z6 = s22h(6)*ADJ_alnthh
291
292 ADJ_xd = ADJ_z6*gx2
293 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
294 ADJ_xx = ADJ_z5*gx2
295 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
296 ADJ_xx = ADJ_z4*gg + ADJ_xx
297 ADJ_gg = xx*ADJ_z4 + ADJ_gg
298 ADJ_xd = ADJ_z3*gg + ADJ_xd
299 ADJ_gg = xd*ADJ_z3 + ADJ_gg
300 ADJ_xd = ADJ_z2*ggg + ADJ_xd
301 ADJ_ggg = xd*ADJ_z2 + ADJ_ggg
302 ADJ_xx = ADJ_z1*ggg + ADJ_xx
303 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
304
305 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
306 ADJ_alnsig = ADJ_xd
307
308 ! vertical
309 xd = xd_save
310 xx = xx_save
311
312 alnthv = alnthv_save
313
314 ADJ_z1 = s22v(1)*ADJ_alnthv
315 ADJ_z2 = s22v(2)*ADJ_alnthv
316 ADJ_z3 = s22v(3)*ADJ_alnthv
317 ADJ_z4 = s22v(4)*ADJ_alnthv
318 ADJ_z5 = s22v(5)*ADJ_alnthv
319 ADJ_z6 = s22v(6)*ADJ_alnthv
320
321 ADJ_xd = ADJ_z6*gx2
322 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
323 ADJ_xx = ADJ_z5*gx2
324 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
325 ADJ_xx = ADJ_z4*gg + ADJ_xx
326 ADJ_gg = xx*ADJ_z4 + ADJ_gg
327 ADJ_xd = ADJ_z3*gg + ADJ_xd
328 ADJ_gg = xd*ADJ_z3 + ADJ_gg
329 ADJ_xd = ADJ_z2*ggg + ADJ_xd
330 ADJ_ggg = xd*ADJ_z2 + ADJ_ggg
331 ADJ_xx = ADJ_z1*ggg + ADJ_xx
332 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
333 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
334 ADJ_alnsig = ADJ_xd + ADJ_alnsig
335
336 else if (ifreq .eq. 3) then
337
338 ADJ_z1 = s37h(1)*ADJ_alnthh
339 ADJ_z2 = s37h(2)*ADJ_alnthh
340 ADJ_z3 = s37h(3)*ADJ_alnthh
341 ADJ_z4 = s37h(4)*ADJ_alnthh
342 ADJ_z5 = s37h(5)*ADJ_alnthh
343 ADJ_z6 = s37h(6)*ADJ_alnthh
344
345 ADJ_xd = ADJ_z6*gx2
346 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
347 ADJ_xx = ADJ_z5*gx2
348 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
349 ADJ_xx = ADJ_z4*gg + ADJ_xx
350 ADJ_gg = xx*ADJ_z4 + ADJ_gg
351 ADJ_xd = ADJ_z3*gg + ADJ_xd
352 ADJ_gg = xd*ADJ_z3 + ADJ_gg
353 ADJ_xd = ADJ_z2*ggg + ADJ_xd
354 ADJ_ggg = xd*ADJ_z2 + ADJ_ggg
355 ADJ_xx = ADJ_z1*ggg + ADJ_xx
356 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
357 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
358 ADJ_alnsig = ADJ_xd
359
360 ! vertical
361 xd = xd_save
362 xx = xx_save
363
364 alnthv = alnthv_save
365
366 ADJ_z1 = s37v(1)*ADJ_alnthv
367 ADJ_z2 = s37v(2)*ADJ_alnthv
368 ADJ_z3 = s37v(3)*ADJ_alnthv
369 ADJ_z4 = s37v(4)*ADJ_alnthv
370 ADJ_z5 = s37v(5)*ADJ_alnthv
371 ADJ_z6 = s37v(6)*ADJ_alnthv
372
373 ADJ_xd = ADJ_z6*gx2
374 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
375 ADJ_xx = ADJ_z5*gx2
376 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
377 ADJ_xx = ADJ_z4*gg + ADJ_xx
378 ADJ_gg = xx*ADJ_z4 + ADJ_gg
379 ADJ_xd =ADJ_z3*gg + ADJ_xd
380 ADJ_gg = xd*ADJ_z3 + ADJ_gg
381 ADJ_xd = ADJ_z2*ggg + ADJ_xd
382 ADJ_ggg = xd*ADJ_z2 + ADJ_ggg
383 ADJ_xx = ADJ_z1*ggg + ADJ_xx
384 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
385 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
386 ADJ_alnsig = ADJ_xd + ADJ_alnsig
387
388 else if (ifreq .eq. 4) then
389
390 ADJ_z1 = s85h(1)*ADJ_alnthh
391 ADJ_z2 = s85h(2)*ADJ_alnthh
392 ADJ_z3 = s85h(3)*ADJ_alnthh
393 ADJ_z4 = s85h(4)*ADJ_alnthh
394 ADJ_z5 = s85h(5)*ADJ_alnthh
395 ADJ_z6 = s85h(6)*ADJ_alnthh
396
397 ADJ_xd = ADJ_z6*gx2
398 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
399 ADJ_xx = ADJ_z5*gx2
400 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
401 ADJ_xx = ADJ_z4*gg + ADJ_xx
402 ADJ_gg = xx*ADJ_z4 + ADJ_gg
403 ADJ_xd = ADJ_z3*gg + ADJ_xd
404 ADJ_gg = xd*ADJ_z3 + ADJ_gg
405 ADJ_xd = ADJ_z2*ggg + ADJ_xd
406 ADJ_ggg = xd*ADJ_z2 + ADJ_ggg
407 ADJ_xx = ADJ_z1*ggg + ADJ_xx
408 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
409 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
410 ADJ_alnsig = ADJ_xd
411
412 ! vertical
413 xd = xd_save
414 xx = xx_save
415
416 alnthv = alnthv_save
417
418 ADJ_z1 = s85v(1)*ADJ_alnthv
419 ADJ_z2 = s85v(2)*ADJ_alnthv
420 ADJ_z3 = s85v(3)*ADJ_alnthv
421 ADJ_z4 = s85v(4)*ADJ_alnthv
422 ADJ_z5 = s85v(5)*ADJ_alnthv
423 ADJ_z6 = s85v(6)*ADJ_alnthv
424
425 ADJ_xd = ADJ_z6*gx2
426 ADJ_gx2 = xd*ADJ_z6 + ADJ_gx2
427 ADJ_xx = ADJ_z5*gx2
428 ADJ_gx2 = xx*ADJ_z5 + ADJ_gx2
429 ADJ_xx = ADJ_z4*gg + ADJ_xx
430 ADJ_gg = xx*ADJ_z4 + ADJ_gg
431 ADJ_xd = ADJ_z3*gg + ADJ_xd
432 ADJ_gg = xd*ADJ_z3 + ADJ_gg
433 ADJ_xd = ADJ_z2*ggg + ADJ_xd
434 ADJ_ggg = xd*ADJ_z2 + ADJ_ggg
435 ADJ_xx = ADJ_z1*ggg + ADJ_xx
436 ADJ_ggg = xx*ADJ_z1 + ADJ_ggg
437 ADJ_xd = 2.0*xd*ADJ_xx + ADJ_xd
438 ADJ_alnsig = ADJ_xd + ADJ_alnsig
439 end if
440
441 ADJ_gg = ADJ_ggg*gx2 + ADJ_gg
442 ADJ_gx2 = gg*ADJ_ggg + ADJ_gx2
443 ADJ_gx2 = 2.0*gx2*ADJ_gg + ADJ_gx2
444
445 if (abs(sigma) .gt. 0.0) then
446 ADJ_sigma = ADJ_alnsig/sigma + ADJ_sigma
447 ! else
448 ! ADJ_sigma = 0.0
449 end if
450
451 if (trace_use) call da_trace_exit("da_effang_adj")
452
453 end subroutine da_effang_adj
454
455