da_final_write_obs.inc
References to this file elsewhere.
1 subroutine da_final_write_obs(iv)
2
3 !-------------------------------------------------------------------------
4 ! Purpose: Writes full diagnostics for O, (O-B) & OMA togather
5 !-------------------------------------------------------------------------
6
7 implicit none
8
9 type (iv_type), intent(in) :: iv ! O-B structure.
10
11 integer :: n, k, iunit
12 integer :: ios ! Error code from MPI routines.
13 integer :: num_obs, num_sound_obs
14 character(len=filename_len), allocatable :: filename(:)
15
16 if (trace_use) call da_trace_entry("da_final_write_obs")
17
18 #ifdef DM_PARALLEL
19 ! Wait to ensure all temporary files have been written
20 call mpi_barrier(comm, ierr)
21 #endif
22
23 if (rootproc) then
24 call da_get_unit(iunit)
25 allocate (filename(0:num_procs-1))
26 do k = 0,num_procs-1
27 write(unit=filename(k),fmt ='(a,i3.3)')'gts_omb_oma.',k
28 end do
29 call da_get_unit(omb_unit)
30 open(unit=omb_unit,file='gts_omb_oma',form='formatted', &
31 status='replace', iostat=ios)
32 if (ios /= 0) call da_error(__FILE__,__LINE__, &
33 (/"Cannot open file gts_omb_oma"/))
34 end if
35
36 num_obs = 0
37 if (iv%info(synop)%nlocal > 0) then
38 do n = 1, iv%info(synop)%nlocal
39 if(iv%info(synop)%proc_domain(1,n)) num_obs = num_obs + 1
40 end do
41 end if
42 call da_proc_sum_int(num_obs)
43 if (num_obs > 0 .and. rootproc) then
44 write(omb_unit,'(a20,i8)')'synop', num_obs
45 num_obs = 0
46 do k = 0,num_procs-1
47 call da_read_omb_tmp(filename(k),iunit,num_obs,'synop',5)
48 end do
49 end if
50
51 !------------------------------------------------------------------
52 ! [2] writing Metar
53 !------------------------------------------------------------------
54
55 num_obs = 0
56 if (iv%info(metar)%nlocal > 0) then
57 do n = 1, iv%info(metar)%nlocal
58 if (iv%info(metar)%proc_domain(1,n)) num_obs = num_obs + 1
59 end do
60 end if
61 call da_proc_sum_int(num_obs)
62 if (num_obs > 0 .and. rootproc) then
63 write(omb_unit,'(a20,20i8)')'metar', num_obs
64 num_obs = 0
65 do k = 0,num_procs-1
66 call da_read_omb_tmp(filename(k),iunit,num_obs,'metar',5)
67 end do
68 end if
69
70 !------------------------------------------------------------------
71 ! [3] writing Ships
72 !------------------------------------------------------------------
73
74 num_obs = 0
75 if (iv%info(ships)%nlocal > 0) then
76 do n = 1, iv%info(ships)%nlocal
77 if(iv%info(ships)%proc_domain(1,n)) num_obs = num_obs + 1
78 end do
79 end if
80 call da_proc_sum_int(num_obs)
81 if (num_obs > 0 .and. rootproc) then
82 write(omb_unit,'(a20,i8)')'ships', num_obs
83 num_obs = 0
84 do k = 0,num_procs-1
85 call da_read_omb_tmp(filename(k),iunit,num_obs,'ships',5)
86 end do
87 end if
88
89 !------------------------------------------------------------------
90 ! [4] writing GeoAMV
91 !------------------------------------------------------------------
92
93 num_obs = 0
94 if (iv%info(geoamv)%nlocal > 0) then
95 do n = 1, iv%info(geoamv)%nlocal
96 if (iv%info(geoamv)%proc_domain(1,n)) num_obs = num_obs + 1
97 end do
98 end if
99 call da_proc_sum_int(num_obs)
100 if (num_obs > 0 .and. rootproc) then
101 write(omb_unit,'(a20,i8)')'geoamv', num_obs
102 num_obs = 0
103 do k = 0,num_procs-1
104 call da_read_omb_tmp(filename(k),iunit,num_obs,'geoamv',6)
105 end do
106 end if
107
108 !------------------------------------------------------------------
109 ! [5] writing PolarAMV
110 !------------------------------------------------------------------
111
112 num_obs = 0
113 if (iv%info(polaramv)%nlocal > 0) then
114 do n = 1, iv%info(polaramv)%nlocal
115 if (iv%info(polaramv)%proc_domain(1,n)) num_obs = num_obs + 1
116 end do
117 end if
118 call da_proc_sum_int(num_obs)
119 if (num_obs > 0 .and. rootproc) then
120 write(omb_unit,'(a20,i8)')'polaramv', num_obs
121 num_obs = 0
122 do k = 0,num_procs-1
123 call da_read_omb_tmp(filename(k),iunit,num_obs,'polaramv',8)
124 end do
125 end if
126
127 !------------------------------------------------------------------
128 ! [5] writing GPSPW
129 !------------------------------------------------------------------
130
131 num_obs = 0
132 if (iv%info(gpspw)%nlocal > 0) then
133 do n = 1, iv%info(gpspw)%nlocal
134 if(iv%info(gpspw)%proc_domain(1,n)) num_obs = num_obs + 1
135 end do
136 end if
137 call da_proc_sum_int(num_obs)
138 if (num_obs > 0 .and. rootproc) then
139 write(omb_unit,'(a20,i8)')'gpspw', num_obs
140 num_obs = 0
141 do k = 0,num_procs-1
142 call da_read_omb_tmp(filename(k),iunit,num_obs,'gpspw',5)
143 end do
144 end if
145
146 !------------------------------------------------------------------
147 ! [6] writing Sonde
148 !------------------------------------------------------------------
149
150 num_obs = 0
151 if (iv%info(sound)%nlocal > 0) then
152 do n = 1, iv%info(sound)%nlocal
153 if (iv%info(sound)%proc_domain(1,n)) num_obs = num_obs + 1
154 end do
155 end if
156 call da_proc_sum_int(num_obs)
157 if (num_obs > 0 .and. rootproc) then
158 write(omb_unit,'(a20,i8)')'sound', num_obs
159 num_sound_obs = num_obs
160 num_obs = 0
161 do k = 0,num_procs-1
162 call da_read_omb_tmp(filename(k),iunit,num_obs,'sound',5)
163 end do
164
165 ! writing Sonde_sfc
166 write(omb_unit,'(a20,i8)')'sonde_sfc', num_sound_obs
167 num_obs = 0
168 do k = 0,num_procs-1
169 call da_read_omb_tmp(filename(k),iunit,num_obs,'sonde_sfc',9)
170 end do
171 end if
172
173 !------------------------------------------------------------------
174 ! [7] writing Airep
175 !------------------------------------------------------------------
176
177 num_obs = 0
178 if (iv%info(airep)%nlocal > 0) then
179 do n = 1, iv%info(airep)%nlocal
180 if(iv%info(airep)%proc_domain(1,n)) num_obs = num_obs + 1
181 end do
182 end if
183 call da_proc_sum_int(num_obs)
184 if (num_obs > 0 .and. rootproc) then
185 write(omb_unit,'(a20,i8)')'airep', num_obs
186 num_obs = 0
187 do k = 0,num_procs-1
188 call da_read_omb_tmp(filename(k),iunit,num_obs,'airep',5)
189 end do
190 end if
191
192 !------------------------------------------------------------------
193 ! [8] writing Pilot
194 !------------------------------------------------------------------
195
196 num_obs = 0
197 if (iv%info(pilot)%nlocal > 0) then
198 do n = 1, iv%info(pilot)%nlocal
199 if(iv%info(pilot)%proc_domain(1,n)) num_obs = num_obs + 1
200 end do
201 end if
202 call da_proc_sum_int(num_obs)
203 if (num_obs > 0 .and. rootproc) then
204 write(omb_unit,'(a20,i8)')'pilot', num_obs
205 num_obs = 0
206 do k = 0,num_procs-1
207 call da_read_omb_tmp(filename(k),iunit,num_obs,'pilot',5)
208 end do
209 end if
210
211 !------------------------------------------------------------------
212 ! [9] writing ssmi_rv
213 !------------------------------------------------------------------
214
215 num_obs = 0
216 if (iv%info(ssmi_rv)%nlocal > 0) then
217 do n = 1, iv%info(ssmi_rv)%nlocal
218 if(iv%info(ssmi_rv)%proc_domain(1,n)) num_obs = num_obs + 1
219 end do
220 end if
221 call da_proc_sum_int(num_obs)
222 if (num_obs > 0 .and. rootproc) then
223 write(omb_unit,'(a20,i8)')'ssmir', num_obs
224 num_obs = 0
225 do k = 0,num_procs-1
226 call da_read_omb_tmp(filename(k),iunit,num_obs,'ssmir',5)
227 end do
228 end if
229
230 !------------------------------------------------------------------
231 ! [10] writing SSMITB
232 !------------------------------------------------------------------
233
234 num_obs = 0
235 if (iv%info(ssmi_tb)%nlocal > 0) then
236 do n = 1, iv%info(ssmi_tb)%nlocal
237 if (iv%info(ssmi_tb)%proc_domain(1,n)) num_obs = num_obs + 1
238 end do
239 end if
240 call da_proc_sum_int(num_obs)
241 if (num_obs > 0 .and. rootproc) then
242 write(omb_unit,'(a20,i8)')'ssmiT', num_obs
243 num_obs = 0
244 do k = 0,num_procs-1
245 call da_read_omb_tmp(filename(k),iunit,num_obs,'ssmiT',5)
246 end do
247 end if
248
249 !------------------------------------------------------------------
250 ! [11] writing SATEM
251 !------------------------------------------------------------------
252
253 num_obs = 0
254 if (iv%info(satem)%nlocal > 0) then
255 do n = 1, iv%info(satem)%nlocal
256 if(iv%info(satem)%proc_domain(1,n)) num_obs = num_obs + 1
257 end do
258 end if
259 call da_proc_sum_int(num_obs)
260 if (num_obs > 0 .and. rootproc) then
261 write(omb_unit,'(a20,i8)')'satem', num_obs
262 num_obs = 0
263 do k = 0,num_procs-1
264 call da_read_omb_tmp(filename(k),iunit,num_obs,'satem',5)
265 end do
266 end if
267
268 !------------------------------------------------------------------
269 ! [12] writing SSMT1
270 !------------------------------------------------------------------
271
272 num_obs = 0
273 if (iv%info(ssmt1)%nlocal > 0) then
274 do n = 1, iv%info(ssmt1)%nlocal
275 if(iv%info(ssmt1)%proc_domain(1,n)) num_obs = num_obs + 1
276 end do
277 end if
278 call da_proc_sum_int(num_obs)
279 if (num_obs > 0 .and. rootproc) then
280 write(omb_unit,'(a20,i8)')'ssmt1', num_obs
281 num_obs = 0
282 do k = 0,num_procs-1
283 call da_read_omb_tmp(filename(k),iunit,num_obs,'ssmt1',5)
284 end do
285 end if
286
287 !------------------------------------------------------------------
288 ! [13] writing SSMT2
289 !------------------------------------------------------------------
290
291 num_obs = 0
292 if (iv%info(ssmt2)%nlocal > 0) then
293 do n = 1, iv%info(ssmt2)%nlocal
294 if(iv%info(ssmt2)%proc_domain(1,n)) num_obs = num_obs + 1
295 end do
296 end if
297 call da_proc_sum_int(num_obs)
298 if (num_obs > 0 .and. rootproc) then
299 write(omb_unit,'(a20,i8)')'ssmt2', num_obs
300 num_obs = 0
301 do k = 0,num_procs-1
302 call da_read_omb_tmp(filename(k),iunit,num_obs,'ssmt2',5)
303 end do
304 end if
305
306 !------------------------------------------------------------------
307 ! [14] writing QSCAT
308 !------------------------------------------------------------------
309
310 num_obs = 0
311 if (iv%info(qscat)%nlocal > 0) then
312 do n = 1, iv%info(qscat)%nlocal
313 if(iv%info(qscat)%proc_domain(1,n)) num_obs = num_obs + 1
314 end do
315 end if
316 call da_proc_sum_int(num_obs)
317 if (num_obs > 0 .and. rootproc) then
318 write(omb_unit,'(a20,i8)')'qscat', num_obs
319 num_obs = 0
320 do k = 0,num_procs-1
321 call da_read_omb_tmp(filename(k),iunit,num_obs,'qscat',5)
322 end do
323 end if
324
325 !------------------------------------------------------------------
326 ! [15] writing Profiler
327 !------------------------------------------------------------------
328
329 num_obs = 0
330 if (iv%info(profiler)%nlocal > 0) then
331 do n = 1, iv%info(profiler)%nlocal
332 if(iv%info(profiler)%proc_domain(1,n)) num_obs = num_obs + 1
333 end do
334 end if
335 call da_proc_sum_int(num_obs)
336 if (num_obs > 0 .and. rootproc) then
337 write(omb_unit,'(a20,i8)')'profiler', num_obs
338 num_obs = 0
339 do k = 0,num_procs-1
340 call da_read_omb_tmp(filename(k),iunit,num_obs,'profiler',8)
341 end do
342 end if
343
344 !------------------------------------------------------------------
345 ! [16] writing Buoy
346 !------------------------------------------------------------------
347
348 num_obs = 0
349 if (iv%info(buoy)%nlocal > 0) then
350 do n = 1, iv%info(buoy)%nlocal
351 if(iv%info(buoy)%proc_domain(1,n)) num_obs = num_obs + 1
352 end do
353 end if
354 call da_proc_sum_int(num_obs)
355 if (num_obs > 0 .and. rootproc) then
356 write(omb_unit,'(a20,i8)')'buoy', num_obs
357 num_obs = 0
358 do k = 0,num_procs-1
359 call da_read_omb_tmp(filename(k),iunit,num_obs,'buoy',4)
360 end do
361 end if
362
363 !------------------------------------------------------------------
364 ! [17] writing Bogus
365 !------------------------------------------------------------------
366
367 num_obs = 0
368 if (iv%info(bogus)%nlocal > 0) then
369 do n = 1, iv%info(bogus)%nlocal
370 if(iv%info(bogus)%proc_domain(1,n)) num_obs = num_obs + 1
371 end do
372 end if
373 call da_proc_sum_int(num_obs)
374 if (num_obs > 0 .and. rootproc) then
375 write(omb_unit,'(a20,i8)')'bogus', num_obs
376 num_obs = 0
377 do k = 0,num_procs-1
378 call da_read_omb_tmp(filename(k),iunit,num_obs,'bogus',5)
379 end do
380 end if
381
382 !------------------------------------------------------------------
383 ! writing AIRS retrievals:
384 !------------------------------------------------------------------
385
386 num_obs = 0
387 if (iv%info(airsr)%nlocal > 0) then
388 do n = 1, iv%info(airsr)%nlocal
389 if(iv%info(airsr)%proc_domain(1,n)) num_obs = num_obs + 1
390 end do
391 end if
392 call da_proc_sum_int(num_obs)
393 if (num_obs > 0 .and. rootproc) then
394 write(omb_unit,'(a20,i8)')'airsr', num_obs
395 num_obs = 0
396 do k = 0,num_procs-1
397 call da_read_omb_tmp(filename(k),iunit,num_obs,'airsr',5)
398 end do
399 end if
400
401 !------------------------------------------------------------------
402 ! writing GPS refractivity
403 !------------------------------------------------------------------
404
405 num_obs = 0
406 if (iv%info(gpsref)%nlocal > 0) then
407 do n = 1, iv%info(gpsref)%nlocal
408 if(iv%info(gpsref)%proc_domain(1,n)) num_obs = num_obs + 1
409 end do
410 end if
411 call da_proc_sum_int(num_obs)
412 if (num_obs > 0 .and. rootproc) then
413 write(omb_unit,'(a20,i8)')'gpsref', num_obs
414 num_obs = 0
415 do k = 0,num_procs-1
416 call da_read_omb_tmp(filename(k),iunit,num_obs,'gpsref',6)
417 end do
418 end if
419
420 if (rootproc) then
421 close(iunit)
422 close(omb_unit)
423 call da_free_unit(iunit)
424 call da_free_unit(omb_unit)
425 deallocate (filename)
426 end if
427
428 if (trace_use) call da_trace_exit("da_final_write_obs")
429
430 end subroutine da_final_write_obs
431
432