da_get_var_diagnostics.inc

References to this file elsewhere.
1 subroutine da_get_var_diagnostics( iv, j)
2 
3    !-----------------------------------------------------------------------
4    ! Purpose: TBD
5    !-----------------------------------------------------------------------
6 
7    implicit none
8 
9    type(ob_type), intent(inout):: iv      ! innovation vector.
10    type(j_type), intent(inout) :: j       ! Cost function.
11 
12    integer                      :: num_obs_tot
13    integer                      :: i,k
14    real                         :: jo_radiance
15    real                         :: temp(63)
16 
17    if (trace_use) call da_trace_entry("da_get_var_diagnostics")
18 
19    !--------------------------------------------------------------------------
20    ! [1.0] Sum up Jo across processors:
21    !--------------------------------------------------------------------------
22 
23    num_obs_tot = num_gpspw_tot + num_synop_tot + num_metar_tot + &
24                  num_pilot_tot + num_ssmi_tot  + num_ssmt1_tot + &
25                  num_ssmt2_tot + num_satem_tot + num_airep_tot + &
26                  num_geoamv_tot+ num_polaramv_tot+ num_ships_tot + &
27                  num_sound_tot + num_qscat_tot + num_profiler_tot+ &
28                  num_buoy_tot + num_gpsref_tot + num_bogus_tot   + &
29                  num_radiance_tot + num_airsr_tot
30 
31    temp(1)  = j % jo % synop_u
32    temp(2)  = j % jo % synop_v
33    temp(3)  = j % jo % synop_t
34    temp(4)  = j % jo % synop_p
35    temp(5)  = j % jo % synop_q
36    temp(6)  = j % jo % metar_u
37    temp(7)  = j % jo % metar_v
38    temp(8)  = j % jo % metar_t
39    temp(9)  = j % jo % metar_p
40    temp(10) = j % jo % metar_q
41    temp(11) = j % jo % ships_u
42    temp(12) = j % jo % ships_v
43    temp(13) = j % jo % ships_t
44    temp(14) = j % jo % ships_p
45    temp(15) = j % jo % ships_q
46    temp(16) = j % jo % geoamv_u
47    temp(17) = j % jo % geoamv_v
48    temp(18) = j % jo % polaramv_u
49    temp(19) = j % jo % polaramv_v      
50    temp(20) = j % jo % gpspw_tpw       
51    temp(21) = j % jo % gpsref_ref      
52    temp(22) = j % jo % sound_u         
53    temp(23) = j % jo % sound_v         
54    temp(24) = j % jo % sound_t         
55    temp(25) = j % jo % sound_q         
56    temp(26) = j % jo % sonde_sfc_u     
57    temp(27) = j % jo % sonde_sfc_v     
58    temp(28) = j % jo % sonde_sfc_t     
59    temp(29) = j % jo % sonde_sfc_p     
60    temp(30) = j % jo % sonde_sfc_q     
61    temp(31) = j % jo % airep_u         
62    temp(32) = j % jo % airep_v         
63    temp(33) = j % jo % airep_t         
64    temp(34) = j % jo % pilot_u         
65    temp(35) = j % jo % pilot_v         
66    temp(36) = j % jo % bogus_u         
67    temp(37) = j % jo % bogus_v         
68    temp(38) = j % jo % bogus_t         
69    temp(39) = j % jo % bogus_q         
70    temp(40) = j % jo % bogus_slp       
71    temp(41) = j % jo % ssmir_speed     
72    temp(42) = j % jo % ssmir_tpw       
73    temp(43) = j % jo % ssmi_tb19v      
74    temp(44) = j % jo % ssmi_tb19h      
75    temp(45) = j % jo % ssmi_tb22v      
76    temp(46) = j % jo % ssmi_tb37v      
77    temp(47) = j % jo % ssmi_tb37h      
78    temp(48) = j % jo % ssmi_tb85v      
79    temp(49) = j % jo % ssmi_tb85h      
80    temp(50) = j % jo % satem_thickness 
81    temp(51) = j % jo % ssmt1_t         
82    temp(52) = j % jo % ssmt2_rh        
83    temp(53) = j % jo % qscat_u         
84    temp(54) = j % jo % qscat_v         
85    temp(55) = j % jo % profiler_u      
86    temp(56) = j % jo % profiler_v      
87    temp(57) = j % jo % buoy_u          
88    temp(58) = j % jo % buoy_v          
89    temp(59) = j % jo % buoy_t          
90    temp(60) = j % jo % buoy_p          
91    temp(61) = j % jo % buoy_q          
92    temp(62) = j % jo % airsr_t         
93    temp(63) = j % jo % airsr_q         
94 
95    call da_proc_sum_real(temp(:))
96 
97    j % jo % synop_u         = temp(1)  
98    j % jo % synop_v         = temp(2)  
99    j % jo % synop_t         = temp(3)  
100    j % jo % synop_p         = temp(4)  
101    j % jo % synop_q         = temp(5)  
102    j % jo % metar_u         = temp(6)  
103    j % jo % metar_v         = temp(7)  
104    j % jo % metar_t         = temp(8)  
105    j % jo % metar_p         = temp(9)  
106    j % jo % metar_q         = temp(10) 
107    j % jo % ships_u         = temp(11) 
108    j % jo % ships_v         = temp(12) 
109    j % jo % ships_t         = temp(13) 
110    j % jo % ships_p         = temp(14) 
111    j % jo % ships_q         = temp(15) 
112    j % jo % geoamv_u        = temp(16) 
113    j % jo % geoamv_v        = temp(17) 
114    j % jo % polaramv_u      = temp(18) 
115    j % jo % polaramv_v      = temp(19) 
116    j % jo % gpspw_tpw       = temp(20) 
117    j % jo % gpsref_ref      = temp(21) 
118    j % jo % sound_u         = temp(22) 
119    j % jo % sound_v         = temp(23) 
120    j % jo % sound_t         = temp(24) 
121    j % jo % sound_q         = temp(25) 
122    j % jo % sonde_sfc_u     = temp(26) 
123    j % jo % sonde_sfc_v     = temp(27) 
124    j % jo % sonde_sfc_t     = temp(28) 
125    j % jo % sonde_sfc_p     = temp(29) 
126    j % jo % sonde_sfc_q     = temp(30) 
127    j % jo % airep_u         = temp(31) 
128    j % jo % airep_v         = temp(32) 
129    j % jo % airep_t         = temp(33) 
130    j % jo % pilot_u         = temp(34) 
131    j % jo % pilot_v         = temp(35) 
132    j % jo % bogus_u         = temp(36) 
133    j % jo % bogus_v         = temp(37) 
134    j % jo % bogus_t         = temp(38) 
135    j % jo % bogus_q         = temp(39) 
136    j % jo % bogus_slp       = temp(40) 
137    j % jo % ssmir_speed     = temp(41) 
138    j % jo % ssmir_tpw       = temp(42) 
139    j % jo % ssmi_tb19v      = temp(43) 
140    j % jo % ssmi_tb19h      = temp(44) 
141    j % jo % ssmi_tb22v      = temp(45) 
142    j % jo % ssmi_tb37v      = temp(46) 
143    j % jo % ssmi_tb37h      = temp(47) 
144    j % jo % ssmi_tb85v      = temp(48) 
145    j % jo % ssmi_tb85h      = temp(49) 
146    j % jo % satem_thickness = temp(50) 
147    j % jo % ssmt1_t         = temp(51) 
148    j % jo % ssmt2_rh        = temp(52) 
149    j % jo % qscat_u         = temp(53) 
150    j % jo % qscat_v         = temp(54) 
151    j % jo % profiler_u      = temp(55) 
152    j % jo % profiler_v      = temp(56) 
153    j % jo % buoy_u          = temp(57) 
154    j % jo % buoy_v          = temp(58) 
155    j % jo % buoy_t          = temp(59) 
156    j % jo % buoy_p          = temp(60) 
157    j % jo % buoy_q          = temp(61) 
158    j % jo % airsr_t         = temp(62) 
159    j % jo % airsr_q         = temp(63) 
160 
161    if (use_rad) then
162       jo_radiance = 0.0
163       do i = 1, iv%num_inst                 ! loop for sensor
164          call da_proc_sum_ints(j % jo % rad(i)% num_ichan(:))
165          call da_proc_sum_real(j % jo % rad(i) % jo_ichan(:))
166          jo_radiance = jo_radiance + sum(j % jo % rad(i) % jo_ichan(:))
167       end do
168    end if
169 
170    !-----------------------------------------------------------------------------
171    ! [2.0] Print out VAR diagnostics:
172    !-----------------------------------------------------------------------------
173 
174    if (rootproc) then
175 
176       write(unit=stdout,fmt=*) ' '
177       write(unit=stdout,fmt='(A)') 'Diagnostics'
178       write(unit=stdout,fmt='(A,F12.2)')   '   Final cost function J       = ', j % total
179       write(unit=stdout,fmt=*) ' '
180 
181       write(unit=stdout,fmt='(a,i8)')    '   Total number of obs.        = ', num_obs_tot
182       write(unit=stdout,fmt='(a,f15.5)') '   Final value of J            = ', j % total
183       write(unit=stdout,fmt='(a,f15.5)') '   Final value of Jo           = ', j % jo % total
184       write(unit=stdout,fmt='(a,f15.5)') '   Final value of Jb           = ', j % jb
185       write(unit=stdout,fmt='(a,f15.5)') '   Final value of Je           = ', j % je
186       if (num_obs_tot > 0) &
187          write(unit=stdout,fmt='(a,f15.5)') '   Final J / total num_obs     = ', j % total / &
188                                                           real(num_obs_tot)
189       write(unit=stdout,fmt='(a,f15.5)') '   Jb factor used(1)           = ', var_scaling1
190       write(unit=stdout,fmt='(a,f15.5)') '   Jb factor used(2)           = ', var_scaling2
191       write(unit=stdout,fmt='(a,f15.5)') '   Jb factor used(3)           = ', var_scaling3
192       write(unit=stdout,fmt='(a,f15.5)') '   Jb factor used(4)           = ', var_scaling4
193       write(unit=stdout,fmt='(a,f15.5)') '   Jb factor used(5)           = ', var_scaling5
194       write(unit=stdout,fmt='(a,f15.5)') '   Jb factor used              = ', jb_factor
195       write(unit=stdout,fmt='(a,f15.5)') '   Je factor used              = ', je_factor
196       write(unit=stdout,fmt=*) ' '
197 
198       if (use_rad) then
199          write(unit=stdout,fmt='(a,i8)')    '   Total number of radiances    = ', num_radiance_tot
200          write(unit=stdout,fmt='(a,f15.5)') '   Cost function for radiances  = ', jo_radiance
201          write(unit=stdout,fmt=*) ' '
202       end if
203 
204       ! [4.2] Output components of Jo:
205 
206       if (iv % num_synop_glo > 0) then
207          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    synop obs, Jo(actual)  = ', &
208                                   iv % num_synop_glo, num_synop_tot, &
209                                   j % jo % synop_u, iv % synop_ef_u, &
210                                   j % jo % synop_v, iv % synop_ef_v, &
211                                   j % jo % synop_t, iv % synop_ef_t, &
212                                   j % jo % synop_p, iv % synop_ef_p, &
213                                   j % jo % synop_q, iv % synop_ef_q
214 
215       end if
216 
217       if (trace_use) call da_trace("da_get_var_diagnostics", &
218          message="Memory increase from internal write")
219 
220       if (iv % num_metar_glo > 0) then
221          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    metar obs, Jo(actual)  = ', &
222                                iv % num_metar_glo, num_metar_tot, &
223                                j % jo % metar_u, iv % metar_ef_u, &
224                                j % jo % metar_v, iv % metar_ef_v, &
225                                j % jo % metar_t, iv % metar_ef_t, &
226                                j % jo % metar_p, iv % metar_ef_p, &
227                                j % jo % metar_q, iv % metar_ef_q    
228       end if
229 
230       if (iv % num_ships_glo > 0) then
231          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    ships obs, Jo(actual)  = ', &
232                                iv % num_ships_glo, num_ships_tot, &
233                                j % jo % ships_u, iv % ships_ef_u, &
234                                j % jo % ships_v, iv % ships_ef_v, &
235                                j % jo % ships_t, iv % ships_ef_t, &
236                                j % jo % ships_p, iv % ships_ef_p, &
237                                j % jo % ships_q, iv % ships_ef_q                                
238       end if
239 
240 
241       if (iv % num_geoamv_glo > 0) then
242          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    geoamv ob, Jo(actual)  = ', &
243                                iv % num_geoamv_glo, num_geoamv_tot, &
244                                j % jo % geoamv_u, iv % geoamv_ef_u, &
245                                j % jo % geoamv_v, iv % geoamv_ef_v, &
246                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0
247       end if
248 
249       if (iv % num_polaramv_glo > 0) then
250          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    polaramv,  Jo(actual)  = ', &
251                                iv % num_polaramv_glo, num_polaramv_tot, &
252                                j % jo % polaramv_u, iv % polaramv_ef_u, &
253                                j % jo % polaramv_v, iv % polaramv_ef_v, &
254                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0
255       end if
256 
257 
258       if (iv % num_gpspw_glo > 0) then
259          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    gpspw obs, Jo(actual)  = ', &
260                                iv % num_gpspw_glo, num_gpspw_tot, &
261                                j % jo % gpspw_tpw, iv % gpspw_ef_tpw, &
262                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0
263       end if
264 
265       if (iv % num_gpsref_glo > 0) then
266          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    gpsref obs, Jo(actual)  = ', &
267                                iv % num_gpsref_glo, num_gpsref_tot, &
268                                j % jo % gpsref_ref, iv % gpsref_ef_ref, &
269                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0
270       end if
271 
272       if (iv % num_sound_glo > 0) then
273          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    sound obs, Jo(actual)  = ', &
274                                iv % num_sound_glo, num_sound_tot, &
275                                j % jo % sound_u, iv % sound_ef_u, &
276                                j % jo % sound_v, iv % sound_ef_v, &
277                                j % jo % sound_t, iv % sound_ef_t, &
278                                j % jo % sound_q, iv % sound_ef_q, 0.0, 1.0 
279          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    sonde obs, Jo(actual)  = ', &
280                                iv % num_sound_glo, num_sound_tot, &
281                                j % jo % sonde_sfc_u, iv % synop_ef_u, &
282                                j % jo % sonde_sfc_v, iv % synop_ef_v, &
283                                j % jo % sonde_sfc_t, iv % synop_ef_t, &
284                                j % jo % sonde_sfc_p, iv % synop_ef_p, &
285                                j % jo % sonde_sfc_q, iv % synop_ef_q
286       end if
287 
288       if (iv % num_airep_glo > 0) then
289          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    airep obs, Jo(actual)  = ', &
290                                iv % num_airep_glo, num_airep_tot, &
291                                j % jo % airep_u, iv % airep_ef_u, &
292                                j % jo % airep_v, iv % airep_ef_v, &
293                                j % jo % airep_t, iv % airep_ef_t, &
294                                0.0, 1.0, 0.0, 1.0
295       end if
296 
297       if (iv % num_bogus_glo > 0) then
298          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    bogus obs, Jo(actual)  = ', &
299                                iv % num_bogus_glo, num_bogus_tot, &
300                                j % jo % bogus_u, iv % bogus_ef_u, &
301                                j % jo % bogus_v, iv % bogus_ef_v, &
302                                j % jo % bogus_t, iv % bogus_ef_t, &
303                                j % jo % bogus_q, iv % bogus_ef_q, &
304                                j % jo % bogus_slp, iv % bogus_ef_slp
305       end if
306 
307       if (iv % num_pilot_glo > 0) then
308          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    pilot obs, Jo(actual)  = ', &
309                                iv % num_pilot_glo, num_pilot_tot, &
310                                j % jo % pilot_u, iv % pilot_ef_u, &
311                                j % jo % pilot_v, iv % pilot_ef_v, &
312                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0
313       end if
314 
315       if (iv % num_ssmi_retrieval_glo > 0) then
316          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    ssmir obs, Jo(actual) = ', &
317                                   iv % num_ssmi_retrieval_glo, num_ssmi_tot, &
318                                   j % jo % ssmir_speed, iv % ssmir_ef_speed, &
319                                   j % jo % ssmir_tpw, iv % ssmir_ef_tpw, &
320                                   0.0, 1.0, 0.0, 1.0, 0.0, 1.0
321       end if
322 
323       if (iv % num_satem_glo > 0) then
324          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    satem obs, Jo(actual)  = ', &
325                                iv % num_satem_glo, num_satem_tot, &
326                                j % jo % satem_thickness, iv % satem_ef_thickness, &
327                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0
328       end if
329 
330       if (iv % num_ssmt1_glo > 0) then
331          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    ssmt1 obs, Jo(actual)  = ', &
332                                iv % num_ssmt1_glo, num_ssmt1_tot, &
333                                j % jo % ssmt1_t, iv % ssmt1_ef_t, &
334                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0
335       end if
336 
337       if (iv % num_ssmt2_glo > 0) then
338          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    ssmt2 obs, Jo(actual)  = ', &
339                                iv % num_ssmt2_glo, num_ssmt2_tot, &
340                                j % jo % ssmt2_rh, iv % ssmt2_ef_rh, &
341                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0
342       end if
343 
344       if (iv % num_qscat_glo > 0) then
345          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    qscat obs, Jo(actual)  = ', &
346                                iv % num_qscat_glo, num_qscat_tot, &
347                                j % jo % qscat_u, iv % qscat_ef_u, &
348                                j % jo % qscat_v, iv % qscat_ef_v, &
349                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0
350       end if
351 
352       if (iv % num_buoy_glo > 0) then
353          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    buoy  obs, Jo(actual)  = ', &
354                                iv % num_buoy_glo, num_buoy_tot, &
355                                j % jo % buoy_u, iv % buoy_ef_u, &
356                                j % jo % buoy_v, iv % buoy_ef_v, &
357                                j % jo % buoy_t, iv % buoy_ef_t, &
358                                j % jo % buoy_p, iv % buoy_ef_p, &
359                                j % jo % buoy_q, iv % buoy_ef_q                                
360       end if
361 
362       if (iv % num_profiler_glo > 0) then
363          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    profiler,  Jo(actual)  = ', &
364                                iv % num_profiler_glo, num_profiler_tot, &
365                                j % jo % profiler_u, iv % profiler_ef_u, &
366                                j % jo % profiler_v, iv % profiler_ef_v, &
367                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0
368       end if
369       if (iv % num_airsr_glo > 0) then
370          write(unit=jo_unit,fmt='(a30,2i8,10f15.5)')'    airsr obs, Jo(actual)  = ', &
371                                iv % num_airsr_glo, num_airsr_tot, &
372                                j % jo % airsr_t, iv % airsr_ef_t, &
373                                j % jo % airsr_q, iv % airsr_ef_q, &
374                                0.0, 1.0, 0.0, 1.0, 0.0, 1.0
375       end if
376       do i = 1, iv%num_inst                 ! loop for sensor
377        do k = 1, iv%instid(i)%nchan
378           if (j % jo % rad(i) % num_ichan(k) > 0) then
379             write(unit=jo_unit,fmt='(a30,a16,i5,i10,8f15.5)')'    radiance,  Jo(actual)  = ', &
380                                iv%instid(i)%rttovid_string, k , &
381                                j % jo % rad(i) % num_ichan(k), &
382                                j % jo % rad(i) % jo_ichan(k)
383           end if
384        end do
385       end do
386    end if
387 
388    if (trace_use) call da_trace_exit("da_get_var_diagnostics")
389       
390 end subroutine da_get_var_diagnostics
391 
392