da_effang_tl.inc
References to this file elsewhere.
1 subroutine da_effang_tl(ifreq,theta,gx2,sigma,effangv,effangh, &
2 TGL_gx2,TGL_sigma,TGL_effangv,TGL_effangh )
3
4 !--------------------------------------------------------------------------
5 ! Purpose : Calculate the effective zenith angle of reflected microwave
6 ! radiation at SSM/I frequencies for vertical and horizontal polarization
7 !
8 ! Input :: TGL_gx2, TGL_sigma
9 ! Output :: TGL_effangv,TGL_effangh,effangv,effangh
10 !--------------------------------------------------------------------------
11
12 implicit none
13
14 integer, intent(in) :: ifreq
15 real, intent(in) :: theta,gx2,sigma,TGL_gx2, TGL_sigma
16 real, intent(out) :: TGL_effangv,TGL_effangh,effangv,effangh
17
18 real c19v,c19h,c22v,c22h,c37v,c37h,c85v,c85h
19 real s19v(6),s19h(6),s22v(6),s22h(6), &
20 s37v(6),s37h(6),s85v(6),s85h(6)
21
22 real :: alnsig,gg,ggg,xd,xx
23 real :: z1,z2,z3,z4,z5,z6,alnthv
24 real :: y,dth,angh,angv,alnthh
25 real :: TGL_alnsig,TGL_gg,TGL_ggg,TGL_xd
26 real :: TGL_z1,TGL_z2,TGL_z3,TGL_z4,TGL_z5,TGL_z6,TGL_alnthv
27 real :: TGL_y,TGL_dth,TGL_angh,TGL_angv,TGL_xx,TGL_alnthh
28
29 data c19v,c19h,c22v,c22h,c37v,c37h,c85v,c85h &
30 /-.5108,.5306,-.5108,.5306,-.6931,.1823,-.9163,.3000/
31 data s19v /.225E+2,.698E+2,-.238E+2,-.648E+1,.402E+0,.262E+1/
32 data s19h /.743E+1,.507E+2,-.206E+2,-.191E+1,.648E-1,.291E+1/
33 data s22v /.249E+2,.701E+2,-.240E+2,-.714E+1,.405E+0,.256E+1/
34 data s22h /.498E+1,.442E+2,-.190E+2,-.129E+1,.803E-2,.277E+1/
35 data s37v /.215E+2,.573E+2,-.211E+2,-.670E+1,.443E+0,.253E+1/
36 data s37h /.869E+1,.571E+2,-.257E+2,-.302E+1,.237E+0,.386E+1/
37 data s85v /.116E+2,.263E+2,-.101E+2,-.358E+1,.270E+0,.175E+1/
38 data s85h /.736E+1,.568E+2,-.254E+2,-.248E+1,.196E+0,.387E+1/
39
40 if (trace_use) call da_trace_entry("da_effang_tl")
41
42 if (gx2 .le. 0.0 .or. sigma .le. 0.0) then
43 effangv = theta
44 TGL_effangv = 0.0
45 effangh = theta
46 TGL_effangh = 0.0
47 return
48 end if
49 alnsig = alog(sigma)
50 if (abs(sigma) .gt. 0.0) then
51 TGL_alnsig = TGL_sigma/sigma
52 else
53 TGL_alnsig = 0.0
54 end if
55 gg = gx2*gx2
56 TGL_gg = 2.0*gx2*TGL_gx2
57 ggg = gg*gx2
58 TGL_ggg = TGL_gg*gx2 + gg*TGL_gx2
59
60 if (ifreq .eq. 1) then
61
62 xd = alnsig - c19v
63 TGL_xd = TGL_alnsig
64 xx = xd*xd
65 TGL_xx = 2.0*xd*TGL_xd
66 z1 = xx*ggg
67 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
68 z2 = xd*ggg
69 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
70 z3 = xd*gg
71 TGL_z3 = TGL_xd*gg + xd*TGL_gg
72 z4 = xx*gg
73 TGL_z4 = TGL_xx*gg + xx*TGL_gg
74 z5 = xx*gx2
75 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
76 z6 = xd*gx2
77 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
78 alnthv = s19v(1)*z1 + s19v(2)*z2 + s19v(3)*z3 + &
79 s19v(4)*z4 + s19v(5)*z5 + s19v(6)*z6
80 TGL_alnthv = s19v(1)*TGL_z1 + s19v(2)*TGL_z2 + s19v(3)*TGL_z3 + &
81 s19v(4)*TGL_z4 + s19v(5)*TGL_z5 + s19v(6)*TGL_z6
82 alnthv = alnthv + 3.611
83
84 xd = alnsig - c19h
85 TGL_xd = TGL_alnsig
86 xx = xd*xd
87 TGL_xx = 2.0*xd*TGL_xd
88 z1 = xx*ggg
89 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
90 z2 = xd*ggg
91 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
92 z3 = xd*gg
93 TGL_z3 = TGL_xd*gg + xd*TGL_gg
94 z4 = xx*gg
95 TGL_z4 = TGL_xx*gg + xx*TGL_gg
96 z5 = xx*gx2
97 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
98 z6 = xd*gx2
99 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
100
101 alnthh = s19h(1)*z1 + s19h(2)*z2 + s19h(3)*z3 + &
102 s19h(4)*z4 + s19h(5)*z5 + s19h(6)*z6
103 TGL_alnthh = s19h(1)*TGL_z1 + s19h(2)*TGL_z2 + s19h(3)*TGL_z3 + &
104 s19h(4)*TGL_z4 + s19h(5)*TGL_z5 + s19h(6)*TGL_z6
105
106 alnthh = alnthh + 3.611
107
108 else if (ifreq .eq. 2) then
109 xd = alnsig - c22v
110 TGL_xd = TGL_alnsig
111 xx = xd*xd
112 TGL_xx = 2.0*xd*TGL_xd
113 z1 = xx*ggg
114 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
115 z2 = xd*ggg
116 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
117 z3 = xd*gg
118 TGL_z3 = TGL_xd*gg + xd*TGL_gg
119 z4 = xx*gg
120 TGL_z4 = TGL_xx*gg + xx*TGL_gg
121 z5 = xx*gx2
122 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
123 z6 = xd*gx2
124 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
125 alnthv = s22v(1)*z1 + s22v(2)*z2 + s22v(3)*z3 + &
126 s22v(4)*z4 + s22v(5)*z5 + s22v(6)*z6
127 TGL_alnthv = s22v(1)*TGL_z1 + s22v(2)*TGL_z2 + s22v(3)*TGL_z3 + &
128 s22v(4)*TGL_z4 + s22v(5)*TGL_z5 + s22v(6)*TGL_z6
129 alnthv = alnthv + 3.611
130 ! TGL_alnthv = TGL_alnthv
131
132 xd = alnsig - c22h
133 TGL_xd = TGL_alnsig
134 xx = xd*xd
135 TGL_xx = 2.0*xd*TGL_xd
136 z1 = xx*ggg
137 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
138 z2 = xd*ggg
139 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
140 z3 = xd*gg
141 TGL_z3 = TGL_xd*gg + xd*TGL_gg
142 z4 = xx*gg
143 TGL_z4 = TGL_xx*gg + xx*TGL_gg
144 z5 = xx*gx2
145 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
146 z6 = xd*gx2
147 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
148 alnthh = s22h(1)*z1 + s22h(2)*z2 + s22h(3)*z3 + &
149 s22h(4)*z4 + s22h(5)*z5 + s22h(6)*z6
150 TGL_alnthh = s22h(1)*TGL_z1 + s22h(2)*TGL_z2 + s22h(3)*TGL_z3 + &
151 s22h(4)*TGL_z4 + s22h(5)*TGL_z5 + s22h(6)*TGL_z6
152 alnthh = alnthh + 3.611
153 ! TGL_alnthh = TGL_alnthh
154 else if (ifreq .eq. 3) then
155 xd = alnsig - c37v
156 TGL_xd = TGL_alnsig
157 xx = xd*xd
158 TGL_xx = 2.0*xd*TGL_xd
159 z1 = xx*ggg
160 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
161 z2 = xd*ggg
162 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
163 z3 = xd*gg
164 TGL_z3 = TGL_xd*gg + xd*TGL_gg
165 z4 = xx*gg
166 TGL_z4 = TGL_xx*gg + xx*TGL_gg
167 z5 = xx*gx2
168 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
169 z6 = xd*gx2
170 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
171 alnthv = s37v(1)*z1 + s37v(2)*z2 + s37v(3)*z3 + &
172 s37v(4)*z4 + s37v(5)*z5 + s37v(6)*z6
173 TGL_alnthv = s37v(1)*TGL_z1 + s37v(2)*TGL_z2 + s37v(3)*TGL_z3 + &
174 s37v(4)*TGL_z4 + s37v(5)*TGL_z5 + s37v(6)*TGL_z6
175 alnthv = alnthv + 3.611
176 ! TGL_alnthv = TGL_alnthv
177
178 xd = alnsig - c37h
179 TGL_xd = TGL_alnsig
180 xx = xd*xd
181 TGL_xx = 2.0*xd*TGL_xd
182 z1 = xx*ggg
183 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
184 z2 = xd*ggg
185 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
186 z3 = xd*gg
187 TGL_z3 = TGL_xd*gg + xd*TGL_gg
188 z4 = xx*gg
189 TGL_z4 = TGL_xx*gg + xx*TGL_gg
190 z5 = xx*gx2
191 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
192 z6 = xd*gx2
193 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
194 alnthh = s37h(1)*z1 + s37h(2)*z2 + s37h(3)*z3 + &
195 s37h(4)*z4 + s37h(5)*z5 + s37h(6)*z6
196 TGL_alnthh = s37h(1)*TGL_z1 + s37h(2)*TGL_z2 + s37h(3)*TGL_z3 + &
197 s37h(4)*TGL_z4 + s37h(5)*TGL_z5 + s37h(6)*TGL_z6
198 alnthh = alnthh + 3.611
199 ! TGL_alnthh = TGL_alnthh
200 else if (ifreq .eq. 4) then
201 xd = alnsig - c85v
202 TGL_xd = TGL_alnsig
203 xx = xd*xd
204 TGL_xx = 2.0*xd*TGL_xd
205 z1 = xx*ggg
206 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
207 z2 = xd*ggg
208 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
209 z3 = xd*gg
210 TGL_z3 = TGL_xd*gg + xd*TGL_gg
211 z4 = xx*gg
212 TGL_z4 = TGL_xx*gg + xx*TGL_gg
213 z5 = xx*gx2
214 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
215 z6 = xd*gx2
216 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
217 alnthv = s85v(1)*z1 + s85v(2)*z2 + s85v(3)*z3 + &
218 s85v(4)*z4 + s85v(5)*z5 + s85v(6)*z6
219 TGL_alnthv = s85v(1)*TGL_z1 + s85v(2)*TGL_z2 + s85v(3)*TGL_z3 + &
220 s85v(4)*TGL_z4 + s85v(5)*TGL_z5 + s85v(6)*TGL_z6
221 alnthv = alnthv + 3.611
222 ! TGL_alnthv = TGL_alnthv
223
224 xd = alnsig - c85h
225 TGL_xd = TGL_alnsig
226 xx = xd*xd
227 TGL_xx = 2.0*xd*TGL_xd
228 z1 = xx*ggg
229 TGL_z1 = TGL_xx*ggg + xx*TGL_ggg
230 z2 = xd*ggg
231 TGL_z2 = TGL_xd*ggg + xd*TGL_ggg
232 z3 = xd*gg
233 TGL_z3 = TGL_xd*gg + xd*TGL_gg
234 z4 = xx*gg
235 TGL_z4 = TGL_xx*gg + xx*TGL_gg
236 z5 = xx*gx2
237 TGL_z5 = TGL_xx*gx2 + xx*TGL_gx2
238 z6 = xd*gx2
239 TGL_z6 = TGL_xd*gx2 + xd*TGL_gx2
240 alnthh = s85h(1)*z1 + s85h(2)*z2 + s85h(3)*z3 + &
241 s85h(4)*z4 + s85h(5)*z5 + s85h(6)*z6
242 TGL_alnthh = s85h(1)*TGL_z1 + s85h(2)*TGL_z2 + s85h(3)*TGL_z3 + &
243 s85h(4)*TGL_z4 + s85h(5)*TGL_z5 + s85h(6)*TGL_z6
244 alnthh = alnthh + 3.611
245 ! TGL_alnthh = TGL_alnthh
246 end if
247 angv = 90.0 - exp(alnthv)
248 TGL_angv = - TGL_alnthv*exp(alnthv)
249 angh = 90.0 - exp(alnthh)
250 TGL_angh = - TGL_alnthh*exp(alnthh)
251 y = 1.0 - 28.0*gx2
252 TGL_y = - 28.0*TGL_gx2
253 if (y .lt. 0.0) then
254 y = 0.0
255 TGL_y = 0.0
256 end if
257 dth = (theta - 53.0)*y
258 TGL_dth = (theta - 53.0)*TGL_y
259 effangv = angv + dth
260 TGL_effangv = TGL_angv + TGL_dth
261 effangh = angh + dth
262 TGL_effangh = TGL_angh + TGL_dth
263
264 if (trace_use) call da_trace_exit("da_effang_tl")
265
266 end subroutine da_effang_tl
267
268