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