wrfv2_fire/var/da/da_minimisation/da_write_diagnostics.inc

   1 subroutine da_write_diagnostics(it, grid,num_qcstat_conv, ob, iv, re, y, j)
   2
   3    !---------------------------------------------------------------------------
   4    ! Purpose: Output data assmiilation diagnostics
   5    !---------------------------------------------------------------------------
   6
   7    implicit none
   8
   9    integer,        intent(in)    :: it
  10    type (domain),  intent(in)    :: grid
  11    integer,        intent(inout) :: num_qcstat_conv(:,:,:,:)
  12    type (y_type),  intent(in)    :: ob      ! Observation structure.
  13    type (iv_type), intent(inout) :: iv      ! innovation vector.
  14    type (y_type),  intent(inout) :: re      ! residual vector.
  15    type (y_type),  intent(in)    :: y       ! y = H(x_inc) structure.
  16    type (j_type),  intent(inout) :: j       ! Cost function.
  17
  18    character(len=13) filename
  19
  20    if (trace_use) call da_trace_entry("da_write_diagnostics")
  21
  22    if (rootproc .and. print_detail_outerloop) then
  23       write(filename,'(a,i2.2)') "statistics_",it
  24       call da_get_unit (stats_unit)
  25       open(unit=stats_unit,file=filename,status="replace")
  26    endif
  27
  28    iv%nstats(:) = 0
  29
  30    !---------------------------------------------------------------------------
  31    ! [1.0] Calculate innovation vector (O-B) statistics:
  32    !---------------------------------------------------------------------------
  33
  34    if (iv%info(sound)%ntotal > 0) then
  35       call da_oi_stats_sound     (stats_unit, iv)
  36       call da_oi_stats_sonde_sfc (stats_unit, iv)
  37    end if
  38    if (iv%info(mtgirs)%ntotal   > 0) call da_oi_stats_mtgirs   (stats_unit, iv)
  39    if (iv%info(tamdar)%ntotal   > 0) call da_oi_stats_tamdar   (stats_unit, iv)
  40    if (iv%info(tamdar)%ntotal   > 0) call da_oi_stats_tamdar_sfc(stats_unit, iv)
  41    if (iv%info(synop)%ntotal    > 0) call da_oi_stats_synop    (stats_unit, iv)
  42    if (iv%info(geoamv)%ntotal   > 0) call da_oi_stats_geoamv   (stats_unit, iv)
  43    if (iv%info(polaramv)%ntotal > 0) call da_oi_stats_polaramv (stats_unit, iv)
  44    if (iv%info(airep)%ntotal    > 0) call da_oi_stats_airep    (stats_unit, iv)
  45    if (iv%info(pilot)%ntotal    > 0) call da_oi_stats_pilot    (stats_unit, iv)
  46    if (iv%info(metar)%ntotal    > 0) call da_oi_stats_metar    (stats_unit, iv)
  47    if (iv%info(ships)%ntotal    > 0) call da_oi_stats_ships    (stats_unit, iv)
  48    if (iv%info(gpspw)%ntotal    > 0) call da_oi_stats_gpspw    (stats_unit, iv)
  49    if (iv%info(gpsref)%ntotal   > 0) call da_oi_stats_gpsref   (stats_unit, iv)
  50    if (iv%info(ssmi_tb)%ntotal  > 0) call da_oi_stats_ssmi_tb  (stats_unit, iv)
  51    if (iv%info(ssmi_rv)%ntotal  > 0) call da_oi_stats_ssmi_rv  (stats_unit, iv)
  52    if (iv%info(satem)%ntotal    > 0) call da_oi_stats_satem    (stats_unit, iv)
  53    if (iv%info(ssmt1)%ntotal    > 0) call da_oi_stats_ssmt1    (stats_unit, iv)
  54    if (iv%info(ssmt2)%ntotal    > 0) call da_oi_stats_ssmt2    (stats_unit, iv)
  55    if (iv%info(qscat)%ntotal    > 0) call da_oi_stats_qscat    (stats_unit, iv)
  56    if (iv%info(profiler)%ntotal > 0) call da_oi_stats_profiler (stats_unit, iv)
  57    if (iv%info(buoy)%ntotal     > 0) call da_oi_stats_buoy     (stats_unit, iv)
  58    if (iv%info(radar)%ntotal    > 0) call da_oi_stats_radar    (stats_unit, iv)
  59    if (iv%info(bogus)%ntotal    > 0) call da_oi_stats_bogus    (stats_unit, iv)
  60    if (iv%info(airsr)%ntotal    > 0) call da_oi_stats_airsr    (stats_unit, iv)
  61
  62 #if defined(RTTOV) || defined(CRTM)
  63    if (iv%num_inst > 0) call da_oi_stats_rad (stats_unit, iv)
  64 #endif
  65
  66    if (num_pseudo > 0) call da_oi_stats_pseudo (stats_unit, iv)
  67
  68    !---- ----------------------------------------------------------------------
  69    ! [2.0] Calculate residual vector (O-A) statistics:
  70    !---------------------------------------------------------------------------
  71 if (.not. anal_type_verify) then
  72
  73    if (iv%info(sound)%ntotal > 0) then
  74       call da_ao_stats_sound (stats_unit, iv, re)
  75       call da_ao_stats_sonde_sfc (stats_unit, iv, re)
  76    end if
  77    if (iv%info(mtgirs)%ntotal   > 0) call da_ao_stats_mtgirs   (stats_unit, iv, re)
  78    if (iv%info(tamdar)%ntotal   > 0) call da_ao_stats_tamdar   (stats_unit, iv, re)
  79    if (iv%info(tamdar)%ntotal   > 0) call da_ao_stats_tamdar_sfc(stats_unit, iv, re)
  80    if (iv%info(synop)%ntotal    > 0) call da_ao_stats_synop    (stats_unit, iv, re)
  81    if (iv%info(geoamv)%ntotal   > 0) call da_ao_stats_geoamv   (stats_unit, iv, re)
  82    if (iv%info(polaramv)%ntotal > 0) call da_ao_stats_polaramv (stats_unit, iv, re)
  83    if (iv%info(airep)%ntotal    > 0) call da_ao_stats_airep    (stats_unit, iv, re)
  84    if (iv%info(pilot)%ntotal    > 0) call da_ao_stats_pilot    (stats_unit, iv, re)
  85    if (iv%info(metar)%ntotal    > 0) call da_ao_stats_metar    (stats_unit, iv, re)
  86    if (iv%info(ships)%ntotal    > 0) call da_ao_stats_ships    (stats_unit, iv, re)
  87    if (iv%info(gpspw)%ntotal    > 0) call da_ao_stats_gpspw    (stats_unit, iv, re)
  88    if (iv%info(gpsref)%ntotal   > 0) call da_ao_stats_gpsref   (stats_unit, iv, re)
  89    if (iv%info(ssmi_tb)%ntotal  > 0) call da_ao_stats_ssmi_tb  (stats_unit, iv, re)
  90    if (iv%info(ssmi_rv)%ntotal  > 0) call da_ao_stats_ssmi_rv  (stats_unit, iv, re)
  91    if (iv%info(satem)%ntotal    > 0) call da_ao_stats_satem    (stats_unit, iv, re)
  92    if (iv%info(ssmt1)%ntotal    > 0) call da_ao_stats_ssmt1    (stats_unit, iv, re)
  93    if (iv%info(ssmt2)%ntotal    > 0) call da_ao_stats_ssmt2    (stats_unit, iv, re)
  94    if (iv%info(qscat)%ntotal    > 0) call da_ao_stats_qscat    (stats_unit, iv, re)
  95    if (iv%info(profiler)%ntotal > 0) call da_ao_stats_profiler (stats_unit, iv, re)
  96    if (iv%info(buoy)%ntotal     > 0) call da_ao_stats_buoy     (stats_unit, iv, re)
  97    if (iv%info(radar)%ntotal    > 0) call da_ao_stats_radar    (stats_unit, iv, re)
  98    if (iv%info(bogus)%ntotal    > 0) call da_ao_stats_bogus    (stats_unit, iv, re)
  99    if (iv%info(airsr)%ntotal    > 0) call da_ao_stats_airsr    (stats_unit, iv, re)
 100
 101 #if defined(CRTM) || defined(RTTOV)
 102    if (iv%num_inst > 0) call da_ao_stats_rad (stats_unit, iv, re)
 103 #endif
 104
 105    if (num_pseudo > 0) call da_ao_stats_pseudo (stats_unit, iv, re)
 106
 107    !---------------------------------------------------------------------------
 108    ! [3.0] Calculate analysis increment (A-B) statistics:
 109    !---------------------------------------------------------------------------
 110
 111    call da_analysis_stats (grid, stats_unit)
 112
 113    !---------------------------------------------------------------------------
 114    ! [4.0] Write VAR diagnostic :
 115    !------------------------------------------------------------------------------
 116
 117    call da_get_var_diagnostics ( iv, j)
 118
 119 end if
 120
 121    !---- -------------------------------------------------------------------------
 122    !  [5.0] Write observation data (O, O-B, O-A, y=hx'):
 123    !------------------------------------------------------------------------------
 124
 125    call da_write_obs(it, ob, iv, re)
 126
 127    ! Write ETKF observation files if required (note - 1PE only at present):
 128    if (anal_type_verify) then
 129       call da_write_obs_etkf(ob, iv, re)
 130    end if
 131
 132    call da_final_write_obs(it, iv)
 133
 134 if (.not. anal_type_verify) then
 135    call da_write_y(iv, y)
 136
 137    call da_final_write_y(iv)
 138    call da_print_qcstat(it, iv, num_qcstat_conv)
 139 end if
 140
 141    if (rootproc .and. print_detail_outerloop) then
 142        close(stats_unit)
 143        call da_free_unit (stats_unit)
 144    end if
 145
 146    if (trace_use) call da_trace_exit("da_write_diagnostics")
 147
 148 end subroutine da_write_diagnostics
 149
 150