da_rttov_tl.inc
References to this file elsewhere.
1 #ifdef RTTOV
2 subroutine da_rttov_tl(inst, nchanl, nprofiles, con_vars, aux_vars, &
3 con_vars_tl, aux_vars_tl, tb)
4
5 !---------------------------------------------------------------------------
6 ! PURPOSE: interface to the tangent linear subroutine of RTTOV
7 !---------------------------------------------------------------------------
8
9 implicit none
10
11 #include "rttov_tl.interface"
12
13 integer , intent (in) :: inst, nchanl, nprofiles
14 type (con_vars_type), intent (in) :: con_vars (nprofiles)
15 type (con_vars_type), intent (in) :: con_vars_tl (nprofiles)
16 type (aux_vars_type), intent (in) :: aux_vars (nprofiles)
17 type (aux_vars_type), intent (in) :: aux_vars_tl (nprofiles)
18 real , intent (out) :: tb(nchanl,nprofiles)
19
20 ! local variables
21 integer :: n, nc
22 Integer :: alloc_status(140)
23
24 ! RTTOV input parameters
25 integer :: nfrequencies, nchannels, nbtout
26 integer :: nchan(nprofiles)
27 integer , pointer :: lprofiles(:)
28 type(rttov_coef) :: coef
29 type(profile_type) :: profiles(nprofiles), profiles_tl(nprofiles)
30 logical :: addcloud
31 real , allocatable :: surfem(:)
32 integer , pointer :: channels (:), polarisations(:,:)
33 logical, allocatable :: calcemis (:)
34
35 ! RTTOV out parameters
36 integer :: errorstatus(nprofiles)
37
38 ! RTTOV inout parameters
39 real , pointer :: emissivity (:), emissivity_tl (:)
40 type (radiance_type) :: radiance, radiance_tl
41 type (transmission_type) :: transmission, transmission_tl
42
43 call da_trace_entry("da_rttov_tl")
44
45 nchan (:) = nchanl
46 coef = coefs(inst)
47 addcloud = .false.
48 alloc_status(:) = 0
49
50 do n = 1, nprofiles
51 profiles(n) % nlevels = con_vars(n) % nlevels
52 allocate (profiles(n)%p(profiles(n) % nlevels), stat=alloc_status(1))
53 allocate (profiles(n)%t(profiles(n) % nlevels), stat=alloc_status(2))
54 allocate (profiles(n)%q(profiles(n) % nlevels), stat=alloc_status(3))
55 allocate (profiles(n)%o3(profiles(n) % nlevels), stat=alloc_status(4))
56 allocate (profiles(n)%co2(profiles(n) % nlevels), stat=alloc_status(5))
57 allocate (profiles(n)%clw(profiles(n) % nlevels), stat=alloc_status(6))
58 if (any(alloc_status /= 0)) then
59 WRITE(UNIT=message(1),FMT='(A,I5)') &
60 "mem allocation error to for profiles",n
61 call da_error(__FILE__,__LINE__,message(1:1))
62 end if
63
64 profiles(n) % ozone_data = .false.
65 profiles(n) % co2_data = .false.
66 profiles(n) % clw_data = .false.
67
68 profiles(n) % p(:) = coef%ref_prfl_p(:)
69 profiles(n) % t(:) = con_vars(n)%t(:)
70 profiles(n) % q(:) = con_vars(n)%q(:)
71 profiles(n) % o3(:) = 0.0 !con_vars(n)%o3(:)
72 profiles(n) % co2(:) = 0.0 !con_vars(n)%co2(:)
73 profiles(n) % clw(:) = 0.0 !con_vars(n)%clw(:)
74
75 profiles(n) % skin % surftype = aux_vars (n) % surftype
76 profiles(n) % skin % t = aux_vars (n) % surft
77 profiles(n) % skin % fastem (:) = 0.0 ! aux_vars (n) % fastem (:)
78
79 profiles(n) % s2m % t = aux_vars (n) % t2m
80 profiles(n) % s2m % q = aux_vars (n) % q2m
81 profiles(n) % s2m % o = 0.0 !aux_vars (n) % o3
82 profiles(n) % s2m % p = con_vars (n) % ps
83 profiles(n) % s2m % u = aux_vars (n) % u10
84 profiles(n) % s2m % v = aux_vars (n) % v10
85
86 profiles(n) % zenangle = aux_vars (n) % satzen
87 profiles(n) % azangle = aux_vars (n) % satazi
88
89 profiles(n) % ctp = 500.0
90 profiles(n) % cfraction = 0.0
91
92 profiles_tl(n) % nlevels = con_vars_tl(n) % nlevels
93 allocate (profiles_tl(n)%p(profiles_tl(n) % nlevels), stat=alloc_status(1))
94 allocate (profiles_tl(n)%t(profiles_tl(n) % nlevels), stat=alloc_status(2))
95 allocate (profiles_tl(n)%q(profiles_tl(n) % nlevels), stat=alloc_status(3))
96 allocate (profiles_tl(n)%o3(profiles_tl(n) % nlevels), stat=alloc_status(4))
97 allocate (profiles_tl(n)%co2(profiles_tl(n) % nlevels), stat=alloc_status(5))
98 allocate (profiles_tl(n)%clw(profiles_tl(n) % nlevels), stat=alloc_status(6))
99 if (any(alloc_status /= 0)) then
100 WRITE(UNIT=message(1),FMT='(A,I5)') &
101 "mem allocation error to for profiles_tl",n
102 call da_error(__FILE__,__LINE__,message(1:1))
103 end if
104
105 profiles_tl(n) % ozone_data = .false.
106 profiles_tl(n) % co2_data = .false.
107 profiles_tl(n) % clw_data = .false.
108
109 profiles_tl(n) % p(:) = 0.0
110 profiles_tl(n) % t(:) = con_vars_tl(n)%t(:)
111 profiles_tl(n) % q(:) = con_vars_tl(n)%q(:)
112 profiles_tl(n) % o3(:) = 0.0 !con_vars(n)%o3(:)
113 profiles_tl(n) % co2(:) = 0.0 !con_vars(n)%co2(:)
114 profiles_tl(n) % clw(:) = 0.0 !con_vars(n)%clw(:)
115
116 profiles_tl(n) % skin % surftype = -1
117 profiles_tl(n) % skin % t = 0.0 !aux_vars_tl (n) % surft
118 profiles_tl(n) % skin % fastem (:) = 0.0 ! aux_vars (n) % fastem (:)
119
120 profiles_tl(n) % s2m % t = 0.0 !aux_vars_tl (n) % t2m
121 profiles_tl(n) % s2m % q = 0.0 !aux_vars_tl (n) % q2m
122 profiles_tl(n) % s2m % o = 0.0 !aux_vars_tl (n) % o3
123 profiles_tl(n) % s2m % p = con_vars_tl (n) % ps
124 profiles_tl(n) % s2m % u = 0.0 !aux_vars_tl (n) % u10
125 profiles_tl(n) % s2m % v = 0.0 !aux_vars_tl (n) % v10
126
127 profiles_tl(n) % zenangle = -1
128 profiles_tl(n) % azangle = -1
129
130 profiles_tl(n) % ctp = 0.0 !500.
131 profiles_tl(n) % cfraction = 0.0
132
133 end do
134
135 #ifdef RTTOV
136 call rttov_setupchan(nprofiles, nchan, coef, & ! in
137 nfrequencies, nchannels, nbtout) ! out
138 #endif
139
140
141 allocate (lprofiles(nfrequencies), stat = alloc_status(31))
142 allocate (channels (nfrequencies), stat = alloc_status(32))
143 allocate (polarisations(nchannels, 3), stat = alloc_status(33))
144 allocate (emissivity(nchannels), stat = alloc_status(34))
145 allocate (emissivity_tl(nchannels), stat = alloc_status(134))
146 allocate (calcemis(nchannels), stat = alloc_status(35))
147 allocate (surfem (nchannels), stat = alloc_status(36))
148
149 ! allocate transmittance structure
150 allocate (transmission % tau_surf (nchannels) ,stat= alloc_status(8))
151 allocate (transmission % tau_layer (coef % nlevels, nchannels) ,stat= alloc_status(9))
152 allocate (transmission % od_singlelayer(coef % nlevels, nchannels ),stat= alloc_status(10))
153
154 allocate (transmission_tl % tau_surf (nchannels) ,stat= alloc_status(108))
155 allocate (transmission_tl % tau_layer (coef % nlevels, nchannels) ,stat= alloc_status(109))
156 allocate (transmission_tl % od_singlelayer(coef % nlevels, nchannels ),stat= alloc_status(110))
157
158
159 ! allocate radiance results arrays with number of channels
160 allocate (radiance % clear (nchannels) ,stat= alloc_status(11))
161 allocate (radiance % cloudy (nchannels) ,stat= alloc_status(12))
162 allocate (radiance % total (nchannels) ,stat= alloc_status(13))
163 allocate (radiance % bt (nchannels) ,stat= alloc_status(14))
164 allocate (radiance % bt_clear (nchannels) ,stat= alloc_status(15))
165 allocate (radiance % upclear (nchannels) ,stat= alloc_status(16))
166 allocate (radiance % dnclear (nchannels) ,stat= alloc_status(17))
167 allocate (radiance % reflclear(nchannels) ,stat= alloc_status(18))
168 allocate (radiance % overcast (coef % nlevels, nchannels) ,stat= alloc_status(19))
169 ! allocate the cloudy radiances with full size even
170 ! if not used
171 allocate (radiance % downcld (coef % nlevels, nchannels) ,stat= alloc_status(20))
172
173 allocate (radiance % out (nbtout) ,stat= alloc_status(121))
174 allocate (radiance % out_clear(nbtout) ,stat= alloc_status(122))
175 allocate (radiance % total_out(nbtout) ,stat= alloc_status(123))
176 allocate (radiance % clear_out(nbtout) ,stat= alloc_status(124))
177
178 ! allocate radiance results arrays with number of channels
179 allocate (radiance_tl % clear (nchannels) ,stat= alloc_status(111))
180 allocate (radiance_tl % cloudy (nchannels) ,stat= alloc_status(112))
181 allocate (radiance_tl % total (nchannels) ,stat= alloc_status(113))
182 allocate (radiance_tl % bt (nchannels) ,stat= alloc_status(114))
183 allocate (radiance_tl % bt_clear (nchannels) ,stat= alloc_status(115))
184 allocate (radiance_tl % upclear (nchannels) ,stat= alloc_status(116))
185 allocate (radiance_tl % dnclear (nchannels) ,stat= alloc_status(117))
186 allocate (radiance_tl % reflclear(nchannels) ,stat= alloc_status(118))
187 allocate (radiance_tl % overcast (coef % nlevels, nchannels) ,stat= alloc_status(119))
188 ! allocate the cloudy radiances with full size even
189 ! if not used
190 allocate (radiance_tl % downcld (coef % nlevels, nchannels) ,stat= alloc_status(120))
191
192 allocate (radiance_tl % out (nbtout) ,stat= alloc_status(121))
193 allocate (radiance_tl % out_clear(nbtout) ,stat= alloc_status(122))
194 allocate (radiance_tl % total_out(nbtout) ,stat= alloc_status(123))
195 allocate (radiance_tl % clear_out(nbtout) ,stat= alloc_status(124))
196
197 if (any(alloc_status /= 0)) then
198 call da_error(__FILE__,__LINE__, &
199 (/"mem allocation error prior to rttov_tl"/))
200 end if
201
202 surfem (:) = 0.0
203 #ifdef RTTOV
204 call rttov_setupindex(nchan, nprofiles, nfrequencies, & ! in
205 nchannels, nbtout, coef, surfem, & ! in
206 lprofiles, channels, polarisations, & ! out
207 emissivity ) ! out
208 #endif
209
210 nc = nchannels/nprofiles
211
212 if (coef%id_sensor == 1) then ! infrared sensor
213 calcemis (1:nchannels) = .true.
214 emissivity (1:nchannels) = 0.0
215 emissivity_tl (1:nchannels) = 0.0
216 else if (coef%id_sensor == 2) then ! microwave sensor
217 do n = 1, nprofiles
218 if (profiles(n) % skin % surftype == 1) then ! sea
219 calcemis ((n-1)*nc+1:n*nc) = .true.
220 emissivity ((n-1)*nc+1:n*nc) = 0.0
221 emissivity_tl ((n-1)*nc+1:n*nc) = 0.0
222 else ! 0:land ; 2:sea-ice
223 calcemis ((n-1)*nc+1:n*nc) = .false.
224 emissivity ((n-1)*nc+1:n*nc) = 0.9
225 emissivity_tl ((n-1)*nc+1:n*nc) = 0.0
226 end if
227 end do
228 end if
229
230 #ifdef RTTOV
231 call rttov_tl(&
232 & errorstatus, &! out
233 & nfrequencies, &! in
234 & nchannels, &! in
235 & nbtout, &! in
236 & nprofiles, &! in
237 & channels, &! in
238 & polarisations, &! in
239 & lprofiles, &! in
240 & profiles, &! in
241 & coef, &! in
242 & addcloud, &! in
243 & calcemis, &! in
244 & emissivity, &! inout
245 & profiles_tl, &! in
246 & emissivity_tl, &! inout
247 & transmission, &! inout
248 & transmission_tl, &! inout
249 & radiance, &! inout
250 & radiance_tl) ! inout
251 #endif
252
253 ! rttov87 generates warnings we want to ignore
254 if (any(errorstatus(:) == errorstatus_fatal)) then
255 write (message(1),*) 'rttov_direct error code = ', errorstatus(:)
256 write (message(2),*) 'nfrequencies = ', nfrequencies
257 write (message(3),*) 'nchannels = ', nchannels
258 write (message(4),*) 'nbtout = ', nbtout
259 write (message(5),*) 'nprofiles = ', nprofiles
260 write (message(6),*) 'channels = ', channels
261 write (message(7),*) 'polarisations = ', polarisations
262 write (message(8),*) 'lprofiles = ', lprofiles
263 write (message(9),*) 'addcloud = ', addcloud
264 write (message(10),*) 'calcemis = ', calcemis
265 write (message(11),*) 'profiles%s2m = ', profiles(1)%s2m
266 write (message(12),*) 'profiles%skin = ', profiles(1)%skin
267 write (message(13),*) 'profiles%zenangle = ', profiles(1)%zenangle
268 write (message(14),*) 'profiles%azangle = ', profiles(1)%azangle
269 write (message(15),*) 'profiles%p = ', profiles(1)%p
270 write (message(16),*) 'profiles%t = ', profiles(1)%t
271 write (message(17),*) 'profiles%q = ', profiles(1)%q
272 write (message(18),*) 'emissivity = ', emissivity
273 write (message(19),*) 'radiance = ', radiance%out_clear
274 write (message(20),*) 'profiles_tl%s2m = ', profiles_tl(1)%s2m
275 write (message(21),*) 'profiles_tl%skin = ', profiles_tl(1)%skin
276 write (message(22),*) 'profiles_tl%zenangle = ', profiles_tl(1)%zenangle
277 write (message(23),*) 'profiles_tl%azangle = ', profiles_tl(1)%azangle
278 write (message(24),*) 'profiles_tl%p = ', profiles_tl(1)%p
279 write (message(25),*) 'profiles_tl%t = ', profiles_tl(1)%t
280 write (message(26),*) 'profiles_tl%q = ', profiles_tl(1)%q
281 write (message(27),*) 'emissivity_tl = ', emissivity_tl
282 write (message(28),*) 'radiance_tl = ', radiance_tl%out_clear
283 call da_warning(__FILE__,__LINE__,message(1:28))
284 end if
285
286 nc = nbtout / nprofiles
287 do n = 1, nprofiles
288 tb(1:nc,n) = radiance_tl % out_clear((n-1)*nc+1:n*nc)
289 end do
290
291 deallocate (lprofiles)
292 deallocate (channels)
293 deallocate (polarisations)
294 deallocate (emissivity)
295 deallocate (emissivity_tl)
296 deallocate (calcemis)
297 deallocate (surfem)
298 do n = 1, nprofiles
299 deallocate (profiles(n)%p)
300 deallocate (profiles(n)%t)
301 deallocate (profiles(n)%q)
302 deallocate (profiles(n)%o3)
303 deallocate (profiles(n)%co2)
304 deallocate (profiles(n)%clw)
305
306 deallocate (profiles_tl(n)%p)
307 deallocate (profiles_tl(n)%t)
308 deallocate (profiles_tl(n)%q)
309 deallocate (profiles_tl(n)%o3)
310 deallocate (profiles_tl(n)%co2)
311 deallocate (profiles_tl(n)%clw)
312 end do
313
314 ! deallocate transmittance structure
315 deallocate (transmission % tau_surf ,stat= alloc_status(6))
316 deallocate (transmission % tau_layer ,stat= alloc_status(7))
317 deallocate (transmission % od_singlelayer,stat= alloc_status(8))
318
319 ! deallocate transmittance structure
320 deallocate (transmission_tl % tau_surf ,stat= alloc_status(106))
321 deallocate (transmission_tl % tau_layer ,stat= alloc_status(107))
322 deallocate (transmission_tl % od_singlelayer,stat= alloc_status(108))
323
324 ! deallocate radiance results arrays with number of channels
325 deallocate (radiance % clear ,stat=alloc_status(9))
326 deallocate (radiance % cloudy ,stat=alloc_status(10))
327 deallocate (radiance % total ,stat=alloc_status(11))
328 deallocate (radiance % bt ,stat=alloc_status(12))
329 deallocate (radiance % bt_clear ,stat=alloc_status(13))
330 deallocate (radiance % upclear ,stat=alloc_status(14))
331 deallocate (radiance % dnclear ,stat=alloc_status(15))
332 deallocate (radiance % reflclear,stat=alloc_status(16))
333 deallocate (radiance % overcast ,stat=alloc_status(17))
334 deallocate (radiance % downcld ,stat=alloc_status(18))
335 deallocate (radiance % out ,stat= alloc_status(19))
336 deallocate (radiance % out_clear ,stat= alloc_status(20))
337 deallocate (radiance % total_out ,stat= alloc_status(21))
338 deallocate (radiance % clear_out ,stat= alloc_status(22))
339
340 deallocate (radiance_tl % clear ,stat=alloc_status(109))
341 deallocate (radiance_tl % cloudy ,stat=alloc_status(110))
342 deallocate (radiance_tl % total ,stat=alloc_status(111))
343 deallocate (radiance_tl % bt ,stat=alloc_status(112))
344 deallocate (radiance_tl % bt_clear ,stat=alloc_status(113))
345 deallocate (radiance_tl % upclear ,stat=alloc_status(114))
346 deallocate (radiance_tl % dnclear ,stat=alloc_status(115))
347 deallocate (radiance_tl % reflclear,stat=alloc_status(116))
348 deallocate (radiance_tl % overcast ,stat=alloc_status(117))
349 deallocate (radiance_tl % downcld ,stat=alloc_status(118))
350 deallocate (radiance_tl % out ,stat= alloc_status(119))
351 deallocate (radiance_tl % out_clear ,stat= alloc_status(120))
352 deallocate (radiance_tl % total_out ,stat= alloc_status(121))
353 deallocate (radiance_tl % clear_out ,stat= alloc_status(122))
354
355
356 if (any(alloc_status /= 0)) then
357 call da_error(__FILE__,__LINE__, &
358 (/"mem deallocation error"/))
359 end if
360
361 call da_trace_exit("da_rttov_tl")
362
363 end subroutine da_rttov_tl
364 #endif