da_search_obs.inc
References to this file elsewhere.
1 subroutine da_search_obs (ob_type_string, unit_in, num_obs, nth, iv, found_flag)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: Search obs. in gts_omb.000
5 !-----------------------------------------------------------------------
6
7 !-------------------------------------------------------------------------
8 ! read iv=O-B structure written by WRFVAR
9 !-------------------------------------------------------------------------
10
11 implicit none
12
13 type (iv_type), intent(inout) :: iv ! O-B structure.
14 integer, intent(in) :: unit_in, nth, num_obs
15 character(len=20), intent(in) :: ob_type_string
16 logical, intent(out) :: found_flag
17
18 character*5 :: stn_id
19 real :: lat, lon
20 integer :: n, n_dummy, k, levels
21
22 if (trace_use) call da_trace_entry("da_search_obs")
23
24 SELECT CASE (trim(adjustl(ob_type_string)))
25
26 CASE ('synop')
27
28 do n = 1, num_obs
29 read(unit_in,'(i8,a5,2E22.13)') n_dummy, stn_id, lat, lon
30 if ( iv%info(synop)%id(nth) == stn_id .and. &
31 iv%info(synop)%lat(1,nth) == lat .and. &
32 iv%info(synop)%lon(1,nth) == lon ) then
33
34 read(unit_in,'(E22.13,5(E22.13,i8,E22.13))')&
35 iv%synop(nth)%h, &
36 iv%synop(nth)%u, &! O-B u
37 iv%synop(nth)%v, &! O-B v
38 iv%synop(nth)%t, &! O-B t
39 iv%synop(nth)%p, &! O-B p
40 iv%synop(nth)%q ! O-B q
41 found_flag = .true.
42 rewind (unit_in)
43 read(unit_in,*)
44 return
45 else
46 read(unit_in,*)
47 endif
48 enddo
49 found_flag = .false.
50
51 CASE ('metar')
52
53 do n = 1, num_obs
54 read(unit_in,'(i8,a5,2E22.13)') n_dummy, stn_id, lat, lon
55 if ( iv%info(metar)%id(nth) == stn_id .and. &
56 iv%info(metar)%lat(1,nth) == lat .and. &
57 iv%info(metar)%lon(1,nth) == lon ) then
58
59 read(unit_in,'(E22.13,5(E22.13,i8,E22.13))')&
60 iv%metar(nth)%h, &
61 iv%metar(nth)%u, &! O-B u
62 iv%metar(nth)%v, &! O-B v
63 iv%metar(nth)%t, &! O-B t
64 iv%metar(nth)%p, &! O-B p
65 iv%metar(nth)%q ! O-B q
66 found_flag = .true.
67 rewind (unit_in)
68 read(unit_in,*)
69 return
70 else
71 read(unit_in,*)
72 endif
73 enddo
74 found_flag = .false.
75
76 CASE ('ships')
77
78 do n = 1, num_obs
79 read(unit_in,'(i8,a5,2E22.13)') n_dummy, stn_id, lat, lon
80 if ( iv%info(ships)%id(nth) == stn_id .and. &
81 iv%info(ships)%lat(1,nth) == lat .and. &
82 iv%info(ships)%lon(1,nth) == lon ) then
83
84 read(unit_in,'(E22.13,5(E22.13,i8,E22.13))')&
85 iv%ships(nth)%h, &
86 iv%ships(nth)%u, &! O-B u
87 iv%ships(nth)%v, &! O-B v
88 iv%ships(nth)%t, &! O-B t
89 iv%ships(nth)%p, &! O-B p
90 iv%ships(nth)%q ! O-B q
91 found_flag = .true.
92 rewind (unit_in)
93 read(unit_in,*)
94 return
95 else
96 read(unit_in,*)
97 endif
98 enddo
99 found_flag = .false.
100
101 CASE ('sonde_sfc')
102
103 do n = 1, num_obs
104 read(unit_in,'(i8,a5,2E22.13)') n_dummy, stn_id, lat, lon
105 if ( iv%info(sonde_sfc)%id(nth) == stn_id .and. &
106 iv%info(sonde_sfc)%lat(1,nth) == lat .and. &
107 iv%info(sonde_sfc)%lon(1,nth) == lon ) then
108
109 read(unit_in,'(E22.13,5(E22.13,i8,E22.13))')&
110 iv%sonde_sfc(nth)%h, &
111 iv%sonde_sfc(nth)%u, &! O-B u
112 iv%sonde_sfc(nth)%v, &! O-B v
113 iv%sonde_sfc(nth)%t, &! O-B t
114 iv%sonde_sfc(nth)%p, &! O-B p
115 iv%sonde_sfc(nth)%q ! O-B q
116 found_flag = .true.
117 rewind (unit_in)
118 read(unit_in,*)
119 return
120 else
121 read(unit_in,*)
122 endif
123 enddo
124 found_flag = .false.
125
126 CASE ('sound')
127
128 do n = 1, num_obs
129 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
130 if ( iv%info(sound)%id(nth) == stn_id .and. &
131 iv%info(sound)%lat(1,nth) == lat .and. &
132 iv%info(sound)%lon(1,nth) == lon ) then
133
134 do k = 1, levels
135 read(unit_in,'(2E22.13,4(E22.13,i8,E22.13))')&
136 iv % sound(nth) % h(k), &
137 iv % sound(nth) % p(k), & ! Obs Pressure
138 iv%sound(nth)%u(k), &! O-B u
139 iv%sound(nth)%v(k), &! O-B v
140 iv%sound(nth)%t(k), &! O-B t
141 iv%sound(nth)%q(k) ! O-B q
142 enddo
143 found_flag = .true.
144 rewind (unit_in)
145 read(unit_in,*)
146 return
147 else
148 do k = 1, levels
149 read(unit_in,*)
150 enddo
151 endif
152 enddo
153 found_flag = .false.
154
155 CASE ('buoy')
156
157 do n = 1, num_obs
158 read(unit_in,'(i8,a5,2E22.13)') n_dummy, stn_id, lat, lon
159 if ( iv%info(buoy)%id(nth) == stn_id .and. &
160 iv%info(buoy)%lat(1,nth) == lat .and. &
161 iv%info(buoy)%lon(1,nth) == lon ) then
162
163 read(unit_in,'(E22.13,5(E22.13,i8,E22.13))')&
164 iv%buoy(nth)%h, &
165 iv%buoy(nth)%u, &! O-B u
166 iv%buoy(nth)%v, &! O-B v
167 iv%buoy(nth)%t, &! O-B t
168 iv%buoy(nth)%p, &! O-B p
169 iv%buoy(nth)%q ! O-B q
170 found_flag = .true.
171 rewind (unit_in)
172 read(unit_in,*)
173 return
174 else
175 read(unit_in,*)
176 endif
177 enddo
178 found_flag = .false.
179
180 CASE ('geoamv')
181
182 do n = 1, num_obs
183 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
184 if ( iv%info(geoamv)%id(nth) == stn_id .and. &
185 iv%info(geoamv)%lat(1,nth) == lat .and. &
186 iv%info(geoamv)%lon(1,nth) == lon ) then
187
188 do k = 1, levels
189 read(unit_in,'(E22.13,2(E22.13,i8,E22.13))')&
190 iv % geoamv(nth) % p(k), & ! Obs Pressure
191 iv%geoamv(nth)%u(k), &! O-B u
192 iv%geoamv(nth)%v(k)
193 enddo
194 found_flag = .true.
195 rewind (unit_in)
196 read(unit_in,*)
197 return
198 else
199 do k = 1, levels
200 read(unit_in,*)
201 enddo
202 endif
203 enddo
204 found_flag = .false.
205
206 CASE ('gpspw')
207
208 do n = 1, num_obs
209 read(unit_in,'(i8,a5,2E22.13)') n_dummy, stn_id, lat, lon
210 if ( iv%info(gpspw)%id(nth) == stn_id .and. &
211 iv%info(gpspw)%lat(1,nth) == lat .and. &
212 iv%info(gpspw)%lon(1,nth) == lon ) then
213
214 read(unit_in,'(E22.13,i8,E22.13)')&
215 iv%gpspw(nth)%tpw
216 found_flag = .true.
217 rewind (unit_in)
218 read(unit_in,*)
219 return
220 else
221 read(unit_in,*)
222 endif
223 enddo
224 found_flag = .false.
225
226 CASE ('ssmir')
227
228 do n = 1, num_obs
229 read(unit_in,'(i8,2E22.13)') n_dummy, lat, lon
230 if ( &
231 iv%info(ssmi_rv)%lat(1,nth) == lat .and. &
232 iv%info(ssmi_rv)%lon(1,nth) == lon ) then
233
234 read(unit_in,'(2(E22.13,i8,E22.13))')&
235 iv%ssmi_rv(nth)%speed, & ! O-B speed
236 iv%ssmi_rv(nth)%tpw ! O-BA tpw
237 found_flag = .true.
238 rewind (unit_in)
239 read(unit_in,*)
240 return
241 else
242 read(unit_in,*)
243 endif
244 enddo
245 found_flag = .false.
246
247 CASE ('airep')
248
249 do n = 1, num_obs
250 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
251 if ( iv%info(airep)%id(nth) == stn_id .and. &
252 iv%info(airep)%lat(1,nth) == lat .and. &
253 iv%info(airep)%lon(1,nth) == lon ) then
254
255 do k = 1, levels
256 read(unit_in,'(2E22.13,3(E22.13,i8,E22.13))')&
257 iv % airep(nth) % h(k), &
258 iv % airep(nth) % p(k), & ! Obs pressure
259 iv%airep(nth)%u(k), &! O-B u
260 iv%airep(nth)%v(k), &! O-B v
261 iv%airep(nth)%t(k)
262 enddo
263 found_flag = .true.
264 rewind (unit_in)
265 read(unit_in,*)
266 return
267 else
268 do k = 1, levels
269 read(unit_in,*)
270 enddo
271 endif
272 enddo
273 found_flag = .false.
274
275 CASE ('polaramv')
276
277 do n = 1, num_obs
278 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
279 if ( iv%info(polaramv)%id(nth) == stn_id .and. &
280 iv%info(polaramv)%lat(1,nth) == lat .and. &
281 iv%info(polaramv)%lon(1,nth) == lon ) then
282
283 do k = 1, levels
284 read(unit_in,'(E22.13,2(E22.13,i8,E22.13))')&
285 iv % polaramv(nth) % p(k), & ! Obs Pressure
286 iv%polaramv(nth)%u(k), &! O-B u
287 iv%polaramv(nth)%v(k)
288 enddo
289 found_flag = .true.
290 rewind (unit_in)
291 read(unit_in,*)
292 return
293 else
294 do k = 1, levels
295 read(unit_in,*)
296 enddo
297 endif
298 enddo
299 found_flag = .false.
300
301 CASE ('pilot')
302
303 do n = 1, num_obs
304 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
305 if ( iv%info(pilot)%id(nth) == stn_id .and. &
306 iv%info(pilot)%lat(1,nth) == lat .and. &
307 iv%info(pilot)%lon(1,nth) == lon ) then
308
309 do k = 1, levels
310 read(unit_in,'(E22.13,2(E22.13,i8,E22.13))')&
311 iv % pilot(nth) % p(k), & ! Obs Pressure
312 iv%pilot(nth)%u(k), &! O-B u
313 iv%pilot(nth)%v(k)
314 enddo
315 found_flag = .true.
316 rewind (unit_in)
317 read(unit_in,*)
318 return
319 else
320 do k = 1, levels
321 read(unit_in,*)
322 enddo
323 endif
324 enddo
325 found_flag = .false.
326
327 CASE ('ssmi_tb')
328
329 do n = 1, num_obs
330 read(unit_in,'(i8,2E22.13)') n_dummy, lat, lon
331 if ( &
332 iv%info(ssmi_tb)%lat(1,nth) == lat .and. &
333 iv%info(ssmi_tb)%lon(1,nth) == lon ) then
334
335 read(unit_in,'(7(E22.13,i8,E22.13))')&
336 iv%ssmi_tb(nth)%tb19h, & ! O-B Tb19h
337 iv%ssmi_tb(nth)%tb19v, & ! O-B Tb19v
338 iv%ssmi_tb(nth)%tb22v, & ! O-B Tb22v
339 iv%ssmi_tb(nth)%tb37h, & ! O-B Tb37h
340 iv%ssmi_tb(nth)%tb37v, & ! O-B Tb37v
341 iv%ssmi_tb(nth)%tb85h, & ! O-B Tb85h
342 iv%ssmi_tb(nth)%tb85v ! O-B Tb85v
343 found_flag = .true.
344 rewind (unit_in)
345 read(unit_in,*)
346 return
347 else
348 read(unit_in,*)
349 endif
350 enddo
351 found_flag = .false.
352
353 CASE ('satem')
354
355 do n = 1, num_obs
356 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
357 if ( iv%info(satem)%id(nth) == stn_id .and. &
358 iv%info(satem)%lat(1,nth) == lat .and. &
359 iv%info(satem)%lon(1,nth) == lon ) then
360
361 do k = 1, levels
362 read(unit_in,'(E22.13,(E22.13,i8,E22.13))')&
363 iv % satem(nth) % p(k), & ! Obs Pressure
364 iv%satem(nth)%thickness(k)
365 enddo
366 found_flag = .true.
367 rewind (unit_in)
368 read(unit_in,*)
369 return
370 else
371 do k = 1, levels
372 read(unit_in,*)
373 enddo
374 endif
375 enddo
376 found_flag = .false.
377
378 CASE ('ssmt1')
379
380 do n = 1, num_obs
381 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
382 if ( iv%info(ssmt1)%id(nth) == stn_id .and. &
383 iv%info(ssmt1)%lat(1,nth) == lat .and. &
384 iv%info(ssmt1)%lon(1,nth) == lon ) then
385
386 do k = 1, levels
387 read(unit_in,'(E22.13,(E22.13,i8,E22.13))')&
388 iv % ssmt1(nth) % h(k), & ! Obs Pressure
389 iv%ssmt1(nth)%t(k)
390 enddo
391 found_flag = .true.
392 rewind (unit_in)
393 read(unit_in,*)
394 return
395 else
396 do k = 1, levels
397 read(unit_in,*)
398 enddo
399 endif
400 enddo
401 found_flag = .false.
402
403 CASE ('ssmt2')
404
405 do n = 1, num_obs
406 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
407 if ( iv%info(ssmt2)%id(nth) == stn_id .and. &
408 iv%info(ssmt2)%lat(1,nth) == lat .and. &
409 iv%info(ssmt2)%lon(1,nth) == lon ) then
410
411 do k = 1, levels
412 read(unit_in,'(E22.13,(E22.13,i8,E22.13))')&
413 iv % ssmt2(nth) % h(k), & ! Obs Pressure
414 iv%ssmt2(nth)%rh(k)
415 enddo
416 found_flag = .true.
417 rewind (unit_in)
418 read(unit_in,*)
419 return
420 else
421 do k = 1, levels
422 read(unit_in,*)
423 enddo
424 endif
425 enddo
426 found_flag = .false.
427
428 CASE ('qscat')
429
430 do n = 1, num_obs
431 read(unit_in,'(i8,a5,2E22.13)') n_dummy, stn_id, lat, lon
432 if ( iv%info(qscat)%id(nth) == stn_id .and. &
433 iv%info(qscat)%lat(1,nth) == lat .and. &
434 iv%info(qscat)%lon(1,nth) == lon ) then
435
436 read(unit_in,'(E22.13,2(E22.13,i8,E22.13))')&
437 iv % qscat(nth) % h, & ! Obs height
438 iv%qscat(nth)%u, &! O-B u
439 iv%qscat(nth)%v ! O-B v
440 found_flag = .true.
441 rewind (unit_in)
442 read(unit_in,*)
443 return
444 else
445 read(unit_in,*)
446 endif
447 enddo
448 found_flag = .false.
449
450 CASE ('profiler')
451
452 do n = 1, num_obs
453 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
454 if ( iv%info(profiler)%id(nth) == stn_id .and. &
455 iv%info(profiler)%lat(1,nth) == lat .and. &
456 iv%info(profiler)%lon(1,nth) == lon ) then
457
458 do k = 1, levels
459 read(unit_in,'(E22.13,2(E22.13,i8,E22.13))')&
460 iv % profiler(nth) % p(k), & ! Obs Pressure
461 iv%profiler(nth)%u(k), &! O-B u
462 iv%profiler(nth)%v(k) ! O-B v
463 enddo
464 found_flag = .true.
465 rewind (unit_in)
466 read(unit_in,*)
467 return
468 else
469 do k = 1, levels
470 read(unit_in,*)
471 enddo
472 endif
473 enddo
474 found_flag = .false.
475
476 CASE ('bogus')
477
478 do n = 1, num_obs
479 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
480 if ( iv%info(bogus)%id(nth) == stn_id .and. &
481 iv%info(bogus)%lat(1,nth) == lat .and. &
482 iv%info(bogus)%lon(1,nth) == lon ) then
483
484 read(unit_in,'(E22.13,i8,E22.13)')&
485 iv%bogus(nth)%slp
486 do k = 1, levels
487 read(unit_in,'(2E22.13,4(E22.13,i8,E22.13))')&
488 iv % bogus(nth) % h(k), &
489 iv % bogus(nth) % p(k), & ! Obs Pressure
490 iv%bogus(nth)%u(k), &! O-B u
491 iv%bogus(nth)%v(k), &! O-B v
492 iv%bogus(nth)%t(k), &! O-B t
493 iv%bogus(nth)%q(k) ! O-B q
494 enddo
495 found_flag = .true.
496 rewind (unit_in)
497 read(unit_in,*)
498 return
499 else
500 read(unit_in,*)
501 do k = 1, levels
502 read(unit_in,*)
503 enddo
504 endif
505 enddo
506 found_flag = .false.
507
508 CASE ('airsr')
509
510 do n = 1, num_obs
511 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
512 if ( iv%info(airsr)%id(nth) == stn_id .and. &
513 iv%info(airsr)%lat(1,nth) == lat .and. &
514 iv%info(airsr)%lon(1,nth) == lon ) then
515
516 do k = 1, levels
517 read(unit_in,'(E22.13,2(E22.13,i8,E22.13))')&
518 iv % airsr(nth) % p(k), & ! Obs Pressure
519 iv%airsr(nth)%t(k), &! O-B t
520 iv%airsr(nth)%q(k) ! O-B q
521 enddo
522 found_flag = .true.
523 rewind (unit_in)
524 read(unit_in,*)
525 return
526 else
527 do k = 1, levels
528 read(unit_in,*)
529 enddo
530 endif
531 enddo
532 found_flag = .false.
533
534 CASE ('gpsref')
535
536 do n = 1, num_obs
537 read(unit_in,'(2i8,a5,2E22.13)') n_dummy, levels, stn_id, lat, lon
538 if ( iv%info(gpsref)%id(nth) == stn_id .and. &
539 iv%info(gpsref)%lat(1,nth) == lat .and. &
540 iv%info(gpsref)%lon(1,nth) == lon ) then
541
542 do k = 1, levels
543 read(unit_in,'(E22.13,(E22.13,i8,E22.13))')&
544 iv % gpsref(nth) % h(k), & ! Obs Height
545 iv%gpsref(nth)%ref(k) ! O-B ref
546 enddo
547 found_flag = .true.
548 rewind (unit_in)
549 read(unit_in,*)
550 return
551 else
552 do k = 1, levels
553 read(unit_in,*)
554 enddo
555 endif
556 enddo
557 found_flag = .false.
558
559 case default;
560
561 write(unit=message(1), fmt='(a,a20,a,i3)') &
562 'Got unknown obs_type string:', trim(ob_type_string),' on unit ',unit_in
563 call da_error(__FILE__,__LINE__,message(1:1))
564
565 END SELECT
566
567 if (trace_use) call da_trace_exit("da_search_obs")
568 return
569 end subroutine da_search_obs