da_final_write_y.inc

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