da_read_omb_tmp.inc
References to this file elsewhere.
1 subroutine da_read_omb_tmp(filename,unit_in,num,obs_type_in, nc)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 !-------------------------------------------------------------------------
8 ! read diagnostics written to temporary file by WRFVAR
9 !-------------------------------------------------------------------------
10
11 implicit none
12
13 integer ,intent (in) :: unit_in
14 integer ,intent (inout) :: num
15 character*(*),intent (in) :: obs_type_in, filename
16 integer ,intent (in) :: nc
17
18 integer :: num_obs, ios
19 character*20 :: ob_type
20 logical :: if_write
21
22 character*5 :: stn_id
23 integer :: n, k, kk, l, levels, dummy_i
24 real :: lat, lon, press, height, dummy
25 real :: tpw_obs, tpw_inv, tpw_err, tpw_inc
26 real :: u_obs, u_inv, u_error, u_inc, &
27 v_obs, v_inv, v_error, v_inc, &
28 t_obs, t_inv, t_error, t_inc, &
29 p_obs, p_inv, p_error, p_inc, &
30 q_obs, q_inv, q_error, q_inc, &
31 ref_obs, ref_inv, ref_error, ref_inc, &
32 spd_obs, spd_inv, spd_err, spd_inc
33 integer :: u_qc, v_qc, t_qc, p_qc, q_qc, tpw_qc, spd_qc, ref_qc
34
35 if (trace_use) call da_trace_entry("da_read_omb_tmp")
36
37 open(unit=unit_in,file=trim(filename),form='formatted',status='old',iostat=ios)
38 if (ios /= 0) Then
39 call da_error(__FILE__,__LINE__, &
40 (/"Cannot open file"//trim(filename)/))
41 end if
42
43 reports: do
44
45 read(unit_in,'(a20,i8)', end = 999, err = 1000) ob_type,num_obs
46 if_write = .false.
47 if (index(ob_type,OBS_type_in(1:nc)) > 0) if_write = .true.
48
49 SELECT CASE (trim(adjustl(ob_type)))
50
51 CASE ('synop' , 'metar' , 'ships' , 'buoy' , 'sonde_sfc' )
52 if (num_obs > 0) then
53 do n = 1, num_obs
54 read(unit_in,'(i8)')levels
55 if (if_write) then
56 write(omb_unit,'(i8)')levels
57 num = num + 1
58 end if
59 do k = 1, levels
60 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
61 kk,l, stn_id, & ! Station
62 lat, lon, press, & ! Lat/lon, pressure
63 u_obs, u_inv, u_qc, u_error, u_inc, &
64 v_obs, v_inv, v_qc, v_error, v_inc, &
65 t_obs, t_inv, t_qc, t_error, t_inc, &
66 p_obs, p_inv, p_qc, p_error, p_inc, &
67 q_obs, q_inv, q_qc, q_error, q_inc
68 if (if_write) &
69 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
70 num, k, stn_id, & ! Station
71 lat, lon, press, & ! Lat/lon, pressure
72 u_obs, u_inv, u_qc, u_error, u_inc, &
73 v_obs, v_inv, v_qc, v_error, v_inc, &
74 t_obs, t_inv, t_qc, t_error, t_inc, &
75 p_obs, p_inv, p_qc, p_error, p_inc, &
76 q_obs, q_inv, q_qc, q_error, q_inc
77 end do
78 end do
79 end if
80 if (if_write) exit reports
81 cycle reports
82
83 CASE ('geoamv' , 'polaramv' )
84
85 if (num_obs > 0) then
86 do n = 1, num_obs
87 read(unit_in,'(i8)')levels
88 if (if_write) then
89 write(omb_unit,'(i8)')levels
90 num = num + 1
91 end if
92 do k = 1, levels
93 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
94 kk, l, stn_id, & ! Station
95 lat, lon, press, & ! Lat/lon, pressure
96 u_obs, u_inv, u_qc, u_error, u_inc, &
97 v_obs, v_inv, v_qc, v_error, v_inc
98 if (if_write) &
99 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
100 num, k, stn_id, & ! Station
101 lat, lon, press, & ! Lat/lon, pressure
102 u_obs, u_inv, u_qc, u_error, u_inc, &
103 v_obs, v_inv, v_qc, v_error, v_inc
104
105 end do
106 end do
107 end if
108
109 if (if_write) exit reports
110 cycle reports
111
112 CASE ('gpspw' )
113 if (num_obs > 0) then
114 do n = 1, num_obs
115 read(unit_in,'(i8)')levels
116 if (if_write) then
117 write(omb_unit,'(i8)')levels
118 num = num + 1
119 end if
120 do k = 1, levels
121 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
122 kk,l, stn_id, & ! Station
123 lat, lon, dummy, & ! Lat/lon, dummy
124 tpw_obs, tpw_inv, tpw_qc, tpw_err, tpw_inc
125 if (if_write) &
126 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
127 num, k, stn_id, & ! Station
128 lat, lon, dummy, & ! Lat/lon, dummy
129 tpw_obs, tpw_inv, tpw_qc, tpw_err, tpw_inc
130 end do
131 end do
132 end if
133
134 if (if_write) exit reports
135 cycle reports
136
137 CASE ('sound' )
138
139 if (num_obs > 0) then
140 do n = 1, num_obs
141 read(unit_in,'(i8)')levels
142 if (if_write) then
143 write(omb_unit,'(i8)')levels
144 num = num + 1
145 end if
146 do k = 1, levels
147 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
148 kk,l, stn_id, & ! Station
149 lat, lon, press, & ! Lat/lon, dummy
150 u_obs, u_inv, u_qc, u_error, u_inc, &
151 v_obs, v_inv, v_qc, v_error, v_inc, &
152 t_obs, t_inv, t_qc, t_error, t_inc, &
153 q_obs, q_inv, q_qc, q_error, q_inc
154 if (if_write) &
155 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
156 num, k, stn_id, & ! Station
157 lat, lon, press, & ! Lat/lon, dummy
158 u_obs, u_inv, u_qc, u_error, u_inc, &
159 v_obs, v_inv, v_qc, v_error, v_inc, &
160 t_obs, t_inv, t_qc, t_error, t_inc, &
161 q_obs, q_inv, q_qc, q_error, q_inc
162 end do
163 end do
164 end if
165 if (if_write) exit reports
166 cycle reports
167
168 CASE ('airep' )
169
170 if (num_obs > 0) then
171 do n = 1, num_obs
172 read(unit_in,'(i8)') levels
173 if (if_write) then
174 write(omb_unit,'(i8)')levels
175 num = num + 1
176 end if
177 do k = 1, levels
178 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
179 kk,l, stn_id, & ! Station
180 lat, lon, press, & ! Lat/lon, dummy
181 u_obs, u_inv, u_qc, u_error, u_inc, &
182 v_obs, v_inv, v_qc, v_error, v_inc, &
183 t_obs, t_inv, t_qc, t_error, t_inc
184 if (if_write) &
185 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
186 num, k, stn_id, & ! Station
187 lat, lon, press, & ! Lat/lon, dummy
188 u_obs, u_inv, u_qc, u_error, u_inc, &
189 v_obs, v_inv, v_qc, v_error, v_inc, &
190 t_obs, t_inv, t_qc, t_error, t_inc
191 end do
192 end do
193 end if
194
195 if (if_write) exit reports
196 cycle reports
197
198 CASE ('pilot' , 'profiler' )
199 if (num_obs > 0) then
200 do n = 1, num_obs
201 read(unit_in,'(i8)')levels
202 if (if_write) then
203 write(omb_unit,'(i8)')levels
204 num = num + 1
205 end if
206 do k = 1, levels
207 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
208 kk,l, stn_id, & ! Station
209 lat, lon, press, & ! Lat/lon, dummy
210 u_obs, u_inv, u_qc, u_error, u_inc, &
211 v_obs, v_inv, v_qc, v_error, v_inc
212 if (if_write) &
213 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
214 num, k, stn_id, & ! Station
215 lat, lon, press, & ! Lat/lon, dummy
216 u_obs, u_inv, u_qc, u_error, u_inc, &
217 v_obs, v_inv, v_qc, v_error, v_inc
218 end do
219 end do
220 end if
221 if (if_write) exit reports
222 cycle reports
223
224 CASE ('ssmir' )
225 if (num_obs > 0) then
226 do n = 1, num_obs
227 read(unit_in,'(i8)')levels
228 if (if_write) then
229 write(omb_unit,'(i8)')levels
230 num = num + 1
231 end if
232 do k = 1, levels
233 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
234 kk,l, stn_id, & ! Station
235 lat, lon, dummy, & ! Lat/lon, dummy
236 spd_obs, spd_inv, spd_qc, spd_err, spd_inc
237 if (if_write) &
238 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
239 num, k, stn_id, & ! Station
240 lat, lon, dummy, & ! Lat/lon, dummy
241 spd_obs, spd_inv, spd_qc, spd_err, spd_inc
242 end do
243 end do
244 end if
245
246 if (if_write) exit reports
247 cycle reports
248
249
250 CASE ('ssmiT' )
251 if (num_obs > 0) then
252 do n = 1, num_obs
253 read(unit_in,'(i8)')levels
254 if (if_write) then
255 write(omb_unit,'(i8)')levels
256 num = num + 1
257 end if
258 do k = 1, levels
259 read(unit_in,'(2i8,a5,2f9.2,f17.7,7(2f17.7,i8,2f17.7))', err= 1000)&
260 kk,l, stn_id, & ! Station
261 lat, lon, dummy, & ! Lat/lon, dummy
262 dummy, dummy, dummy_i, dummy, dummy, &
263 dummy, dummy, dummy_i, dummy, dummy, &
264 dummy, dummy, dummy_i, dummy, dummy, &
265 dummy, dummy, dummy_i, dummy, dummy, &
266 dummy, dummy, dummy_i, dummy, dummy, &
267 dummy, dummy, dummy_i, dummy, dummy, &
268 dummy, dummy, dummy_i, dummy, dummy
269 if (if_write) &
270 write(omb_unit,'(2i8,a5,2f9.2,f17.7,7(2f17.7,i8,2f17.7))', err= 1000)&
271 num,k,stn_id, & ! Station
272 lat, lon, dummy, & ! Lat/lon, dummy
273 dummy, dummy, dummy_i, dummy, dummy, &
274 dummy, dummy, dummy_i, dummy, dummy, &
275 dummy, dummy, dummy_i, dummy, dummy, &
276 dummy, dummy, dummy_i, dummy, dummy, &
277 dummy, dummy, dummy_i, dummy, dummy, &
278 dummy, dummy, dummy_i, dummy, dummy, &
279 dummy, dummy, dummy_i, dummy, dummy
280 end do
281 end do
282 end if
283 if (if_write) exit reports
284 cycle reports
285
286 CASE ('satem' )
287 if (num_obs > 0) then
288 do n = 1, num_obs
289 read(unit_in,'(i8)') levels
290 if (if_write) then
291 write(omb_unit,'(i8)')levels
292 num = num + 1
293 end if
294 do k = 1, levels
295 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
296 kk,l, stn_id, & ! Station
297 lat, lon, dummy, & ! Lat/lon, dummy
298 dummy,dummy, dummy_i, dummy, dummy
299 if (if_write) &
300 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
301 num,k,stn_id, & ! Station
302 lat, lon, dummy, & ! Lat/lon, dummy
303 dummy,dummy, dummy_i, dummy, dummy
304 end do
305 end do
306 end if
307
308 if (if_write) exit reports
309 cycle reports
310
311 CASE ('ssmt1' , 'ssmt2' )
312 if (num_obs > 0) then
313 do n = 1, num_obs
314 read(unit_in,'(i8)') levels
315 if (if_write) then
316 write(omb_unit,'(i8)')levels
317 num = num + 1
318 end if
319 do k = 1, levels
320 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
321 kk,l, stn_id, & ! Station
322 lat, lon, dummy, & ! Lat/lon, dummy
323 dummy,dummy, dummy_i, dummy, dummy
324 if (if_write) &
325 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
326 num,k,stn_id, & ! Station
327 lat, lon, dummy, & ! Lat/lon, dummy
328 dummy,dummy, dummy_i, dummy, dummy
329 end do
330 end do
331 end if
332
333 if (if_write) exit reports
334 cycle reports
335
336 CASE ('qscat' )
337 if (num_obs > 0) then
338 do n = 1, num_obs
339 read(unit_in,'(i8)') levels
340 if (if_write) then
341 write(omb_unit,'(i8)')levels
342 num = num + 1
343 end if
344 do k = 1, levels
345 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
346 kk,l, stn_id, & ! Station
347 lat, lon, press, & ! Lat/lon, dummy
348 u_obs, u_inv, u_qc, u_error, u_inc, &
349 v_obs, v_inv, v_qc, v_error, v_inc
350 if (if_write) &
351 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
352 num,k,stn_id, & ! Station
353 lat, lon, press, & ! Lat/lon, dummy
354 u_obs, u_inv, u_qc, u_error, u_inc, &
355 v_obs, v_inv, v_qc, v_error, v_inc
356 end do
357 end do
358 end if
359
360 if (if_write) exit reports
361 cycle reports
362
363 CASE ('bogus' )
364 ! TC Bogus data is written in two records
365 ! 1st record holds info about surface level
366 ! 2nd is for upper air
367
368 if (num_obs > 0) then
369 do n = 1, num_obs
370 read(unit_in,'(i8)') levels
371 if (if_write) then
372 write(omb_unit,'(i8)')levels
373 num = num + 1
374 end if
375 do k = 1, levels
376 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
377 kk,l, stn_id, & ! Station
378 lat, lon, press, & ! Lat/lon, dummy
379 u_obs, u_inv, u_qc, u_error, u_inc, &
380 v_obs, v_inv, v_qc, v_error, v_inc
381 if (if_write) &
382 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
383 num,l,stn_id, & ! Station
384 lat, lon, press, & ! Lat/lon, dummy
385 u_obs, u_inv, u_qc, u_error, u_inc, &
386 v_obs, v_inv, v_qc, v_error, v_inc
387 end do
388 read(unit_in,'(i8)') levels
389 if (if_write) then
390 write(omb_unit,'(i8)')levels
391 end if
392 do k = 1, levels
393 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
394 kk,l, stn_id, & ! Station
395 lat, lon, press, & ! Lat/lon, dummy
396 u_obs, u_inv, u_qc, u_error, u_inc, &
397 v_obs, v_inv, v_qc, v_error, v_inc
398 if (if_write) &
399 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
400 num,l,stn_id, & ! Station
401 lat, lon, press, & ! Lat/lon, dummy
402 u_obs, u_inv, u_qc, u_error, u_inc, &
403 v_obs, v_inv, v_qc, v_error, v_inc
404 end do
405 end do
406 end if
407 if (if_write) exit reports
408 cycle reports
409
410 CASE ('airsr' )
411 if (num_obs > 0) then
412 do n = 1, num_obs
413 read(unit_in,'(i8)') levels
414 if (if_write) write(omb_unit,'(i8)')levels
415 num = num + 1
416 do k = 1, levels
417 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
418 kk,l, stn_id, & ! Station
419 lat, lon, press, & ! Lat/lon, dummy
420 t_obs, t_inv, t_qc, t_error, t_inc, &
421 q_obs, q_inv, q_qc, q_error, q_inc
422 if (if_write) &
423 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
424 num,k,stn_id, & ! Station
425 lat, lon, press, & ! Lat/lon, dummy
426 t_obs, t_inv, t_qc, t_error, t_inc, &
427 q_obs, q_inv, q_qc, q_error, q_inc
428 end do
429 end do
430 end if
431
432 if (if_write) exit reports
433 cycle reports
434
435 CASE ('gpsref' )
436 if (num_obs > 0) then
437 do n = 1, num_obs
438 read(unit_in,'(i8)') levels
439 if (if_write) write(omb_unit,'(i8)')levels
440 num = num + 1
441 do k = 1, levels
442 read(unit_in,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
443 kk,l, stn_id, & ! Station
444 lat, lon, height, & ! Lat/lon, height
445 ref_obs, ref_inv, ref_qc, ref_error, ref_inc
446 if (if_write) &
447 write(omb_unit,'(2i8,a5,2f9.2,f17.7,5(2f17.7,i8,2f17.7))', err= 1000)&
448 num,k,stn_id, & ! Station
449 lat, lon, height, & ! Lat/lon, height
450 ref_obs, ref_inv, ref_qc, ref_error, ref_inc
451 end do
452 end do
453 end if
454
455 if (if_write) exit reports
456 cycle reports
457 case default;
458 !
459 write(unit=message(1), fmt='(a,a20,a,i3)') &
460 'Got unknown obs_type string:', trim(ob_type),' on unit ',unit_in
461 call da_error(__FILE__,__LINE__,message(1:1))
462 END SELECT
463 END DO reports
464
465 999 continue
466 close (unit_in)
467
468 if (trace_use) call da_trace_exit("da_read_omb_tmp")
469 return
470
471 1000 continue
472 write(unit=message(1), fmt='(a,i3)') &
473 'read error on unit: ',unit_in
474 call da_warning(__FILE__,__LINE__,message(1:1))
475
476 end subroutine da_read_omb_tmp