var/da/da_statistics/da_analysis_stats.inc

   1 subroutine da_analysis_stats  (grid, stats_unit, stats_unit2)
   2
   3    !------------------------------------------------------------------------
   4    ! Purpose: Calculate min, max, mean and RMS of input 1d field.
   5    !------------------------------------------------------------------------
   6 #if (WRF_CHEM == 1)
   7 use module_state_description, only : num_chem, PARAM_FIRST_SCALAR                   &
   8                      ,p_chem_ic_p25,  p_chem_ic_sulf                                &
   9                      ,p_chem_ic_bc1, p_chem_ic_bc2, p_chem_ic_oc1, p_chem_ic_oc2    &
  10                      ,p_chem_ic_dust_1, p_chem_ic_dust_2, p_chem_ic_seas_1, p_chem_ic_seas_2  &
  11                      ,p_chem_ic_bc_a01, p_chem_ic_bc_a02, p_chem_ic_bc_a03    &
  12                      ,p_chem_ic_oc_a01, p_chem_ic_oc_a02, p_chem_ic_oc_a03    &
  13                      ,p_chem_ic_so4_a01, p_chem_ic_so4_a02, p_chem_ic_so4_a03    &
  14                      ,p_chem_ic_no3_a01, p_chem_ic_no3_a02, p_chem_ic_no3_a03    &
  15                      ,p_chem_ic_nh4_a01, p_chem_ic_nh4_a02, p_chem_ic_nh4_a03    &
  16                      ,p_chem_ic_cl_a01, p_chem_ic_cl_a02, p_chem_ic_cl_a03    &
  17                      ,p_chem_ic_na_a01, p_chem_ic_na_a02, p_chem_ic_na_a03    &
  18                      ,p_chem_ic_oin_a01, p_chem_ic_oin_a02, p_chem_ic_oin_a03
  19 #endif
  20
  21    implicit none
  22
  23    type (domain), intent (in) :: grid
  24    integer,       intent (in) :: stats_unit ! Output unit for stats.
  25    integer, optional,  intent (in) :: stats_unit2 ! Output unit for chem stats.
  26
  27    integer :: i, j, k
  28    integer :: ij_g, ijk_g   ! ij, ijk for global domain.
  29    integer :: kdim  ! k range
  30
  31    real    :: um, vm, tm, pm, qm , qcwm, qrnm ! On local domain.
  32    real    :: qcim, qsnm, qgrm
  33    real    :: rij_g, rijk_g      ! On global domain.
  34
  35    type (maxmin_field_type) :: max_u(kts:kte), max_v(kts:kte),     &
  36                                max_t(kts:kte), max_p(kts:kte),     &
  37                                max_q(kts:kte),                     &
  38                                min_u(kts:kte), min_v(kts:kte),     &
  39                                min_t(kts:kte), min_p(kts:kte),     &
  40                                min_q(kts:kte)
  41
  42    type (maxmin_field_type) :: max_qcw(kts:kte), max_qrn(kts:kte), &
  43                                max_qci(kts:kte), max_qsn(kts:kte), &
  44                                max_qgr(kts:kte),                   &
  45                                min_qcw(kts:kte), min_qrn(kts:kte), &
  46                                min_qci(kts:kte), min_qsn(kts:kte), &
  47                                min_qgr(kts:kte)
  48 #if (WRF_CHEM == 1)
  49    type (maxmin_field_type) :: max_chem(kts:kte,num_chem), min_chem(kts:kte,num_chem)
  50    real    :: chemm(num_chem), chemv(kts:kte,num_chem)
  51    integer :: ic
  52 #endif
  53
  54    real                     :: uv(kts:kte), vv(kts:kte),     &
  55                                tv(kts:kte), pv(kts:kte),     &
  56                                qv(kts:kte)
  57    real                     :: qcwv(kts:kte), qrnv(kts:kte), &
  58                                qciv(kts:kte), qsnv(kts:kte), &
  59                                qgrv(kts:kte)
  60
  61    call da_trace_entry("da_analysis_stats")
  62
  63    kdim = kte-kts+1
  64
  65    ij_g = mix * mjy
  66    ijk_g = ij_g * mkz
  67
  68    rij_g  = 1.0/real(ij_g)
  69    rijk_g = 1.0/real(ijk_g)
  70
  71    if (rootproc) then
  72       write(unit=stats_unit, fmt='(/a/)') ' Minimum of gridded analysis increments'
  73       select case ( cloud_cv_options )
  74       case ( 0 )
  75          write(unit=stats_unit, fmt='(6a/)') &
  76               ' Lvl         ', &
  77               'u     i    j          ', &
  78               'v     i    j          ', &
  79               't     i    j          ', &
  80               'p     i    j          ', &
  81               'q     i    j'
  82       case ( 1 )
  83          write(unit=stats_unit, fmt='(8a/)') &
  84               ' Lvl         ', &
  85               'u     i    j          ', &
  86               'v     i    j          ', &
  87               't     i    j          ', &
  88               'p     i    j          ', &
  89               'q     i    j          ', &
  90               'qcw   i    j          ', &
  91               'qrn   i    j'
  92       case ( 2, 3 )
  93          write(unit=stats_unit, fmt='(11a/)') &
  94               ' Lvl         ', &
  95               'u     i    j          ', &
  96               'v     i    j          ', &
  97               't     i    j          ', &
  98               'p     i    j          ', &
  99               'q     i    j          ', &
 100               'qcw   i    j          ', &
 101               'qrn   i    j          ', &
 102               'qci   i    j          ', &
 103               'qsn   i    j          ', &
 104               'qgr   i    j'
 105       end select
 106
 107 #if (WRF_CHEM == 1)
 108       write(unit=stats_unit2, fmt='(/a/)') ' Minimum of gridded analysis increments'
 109       select case ( chem_cv_options )
 110        case ( 10 )
 111           write(unit=stats_unit2, fmt='(11a/)') &
 112                ' Lvl         ', &
 113                'p25      i    j          ', &
 114                'sulf     i    j          ', &
 115                'bc1      i    j          ', &
 116                'bc2      i    j          ', &
 117                'oc1      i    j          ', &
 118                'oc2      i    j          ', &
 119                'dust_1   i    j          ', &
 120                'dust_2   i    j          ', &
 121                'seas_1   i    j          ', &
 122                'seas_2   i    j'
 123        case ( 108 )
 124           write(unit=stats_unit2, fmt='(9a/)') &
 125                ' Lvl, so4aj i j, so4ai i j, nh4aj i j, nh4ai i j, no3aj i j,', &
 126                'no3ai  i j, naaj i j, naai i j, claj i j, clai i j,', &
 127                'asoa1j i j, asoa1i i j, asoa2j i j, asoa2i i j, ', &
 128                'asoa3j i j, asoa3i i j, asoa4j i j, asoa4i i j,', &
 129                'bsoa1j i j, bsoa1i i j, bsoa2j i j, bsoa2i i j, ', &
 130                'bsoa3j i j, bsoa3i i j, bsoa4j i j, bsoa4i i j,', &
 131                'orgpaj i j, orgpai i j, ecj i j,    eci i j,', &
 132                'p25j   i j,   p25i i j, antha  i j, seas   i j, soila i j,', &
 133                'so2 i j,    no2 i j,    o3 i j,     co i j'
 134       end select
 135 #endif
 136    end if
 137
 138    call da_maxmin_in_field(grid%xa%u(its:ite,jts:jte,kts:kte), max_u, min_u)
 139    call da_proc_maxmin_combine(kdim, max_u, min_u)
 140    call da_maxmin_in_field(grid%xa%v(its:ite,jts:jte,kts:kte), max_v, min_v)
 141    call da_proc_maxmin_combine(kdim, max_v, min_v)
 142    call da_maxmin_in_field(grid%xa%t(its:ite,jts:jte,kts:kte), max_t, min_t)
 143    call da_proc_maxmin_combine(kdim, max_t, min_t)
 144    call da_maxmin_in_field(grid%xa%p(its:ite,jts:jte,kts:kte), max_p, min_p)
 145    call da_proc_maxmin_combine(kdim, max_p, min_p)
 146    call da_maxmin_in_field(grid%xa%q(its:ite,jts:jte,kts:kte), max_q, min_q)
 147    call da_proc_maxmin_combine(kdim, max_q, min_q)
 148
 149    if ( cloud_cv_options >= 1 ) then
 150       call da_maxmin_in_field(grid%xa%qcw(its:ite,jts:jte,kts:kte), max_qcw, min_qcw)
 151       call da_proc_maxmin_combine(kdim, max_qcw, min_qcw)
 152       call da_maxmin_in_field(grid%xa%qrn(its:ite,jts:jte,kts:kte), max_qrn, min_qrn)
 153       call da_proc_maxmin_combine(kdim, max_qrn, min_qrn)
 154    end if
 155
 156    if ( cloud_cv_options >= 2 ) then
 157       call da_maxmin_in_field(grid%xa%qci(its:ite,jts:jte,kts:kte), max_qci, min_qci)
 158       call da_proc_maxmin_combine(kdim, max_qci, min_qci)
 159       call da_maxmin_in_field(grid%xa%qsn(its:ite,jts:jte,kts:kte), max_qsn, min_qsn)
 160       call da_proc_maxmin_combine(kdim, max_qsn, min_qsn)
 161       call da_maxmin_in_field(grid%xa%qgr(its:ite,jts:jte,kts:kte), max_qgr, min_qgr)
 162       call da_proc_maxmin_combine(kdim, max_qgr, min_qgr)
 163    end if
 164
 165 #if (WRF_CHEM == 1)
 166    if ( chem_cv_options >= 10 ) then
 167    do ic=PARAM_FIRST_SCALAR, num_chem
 168      call da_maxmin_in_field(grid%xachem%chem_ic(its:ite,jts:jte,kts:kte, ic), max_chem(:,ic), min_chem(:,ic))
 169      call da_proc_maxmin_combine(kdim, max_chem(:,ic), min_chem(:,ic))
 170    end do
 171    end if
 172 #endif
 173
 174    um = 999999.0
 175    vm = 999999.0
 176    tm = 999999.0
 177    pm = 999999.0
 178    qm = 999999.0
 179
 180    qcwm = 999999.0
 181    qrnm = 999999.0
 182    qcim = 999999.0
 183    qsnm = 999999.0
 184    qgrm = 999999.0
 185
 186 #if (WRF_CHEM == 1)
 187    chemm = 999999.0
 188 #endif
 189
 190    do k = kts, kte
 191       if (rootproc) then
 192          if ( abs(min_q(k)%value) < 1.e-30 ) min_q(k)%value = 0.0
 193          select case ( cloud_cv_options )
 194          case ( 0 )
 195             write(unit=stats_unit, fmt='(i4,4(f12.4,2i5),e12.4,2i5)')    k, &
 196                min_u(k), min_v(k), min_t(k), min_p(k), min_q(k)
 197          case ( 1 )
 198             write(unit=stats_unit, fmt='(i4,4(f12.4,2i5),3(e12.4,2i5))') k, &
 199                min_u(k), min_v(k), min_t(k), min_p(k), min_q(k),            &
 200                min_qcw(k), min_qrn(k)
 201          case ( 2, 3 )
 202             write(unit=stats_unit, fmt='(i4,4(f12.4,2i5),6(e12.4,2i5))') k, &
 203                min_u(k), min_v(k), min_t(k), min_p(k), min_q(k),            &
 204                min_qcw(k), min_qrn(k), min_qci(k), min_qsn(k), min_qgr(k)
 205          end select
 206
 207 #if (WRF_CHEM == 1)
 208          select case ( chem_cv_options )
 209          case ( 10 )
 210             write(unit=stats_unit2, fmt='(i4,4(f12.4,2i5),6(e12.4,2i5))') k, &
 211                min_chem(k, p_chem_ic_p25), min_chem(k, p_chem_ic_sulf), &
 212                min_chem(k, p_chem_ic_bc1), min_chem(k, p_chem_ic_bc2), min_chem(k, p_chem_ic_oc1), min_chem(k, p_chem_ic_oc2), &
 213                min_chem(k, p_chem_ic_dust_1), min_chem(k, p_chem_ic_dust_2), min_chem(k, p_chem_ic_seas_1), min_chem(k, p_chem_ic_seas_1)
 214          case ( 20 )
 215             write(unit=stats_unit2, fmt='(i4,24(f12.4,2i5))') k, &
 216                min_chem(k, p_chem_ic_bc_a01), min_chem(k, p_chem_ic_bc_a02), min_chem(k, p_chem_ic_bc_a03), &
 217                min_chem(k, p_chem_ic_oc_a01), min_chem(k, p_chem_ic_oc_a02), min_chem(k, p_chem_ic_oc_a03), &
 218                min_chem(k, p_chem_ic_so4_a01), min_chem(k, p_chem_ic_so4_a02), min_chem(k, p_chem_ic_so4_a03), &
 219                min_chem(k, p_chem_ic_no3_a01), min_chem(k, p_chem_ic_no3_a02), min_chem(k, p_chem_ic_no3_a03), &
 220                min_chem(k, p_chem_ic_nh4_a01), min_chem(k, p_chem_ic_nh4_a02), min_chem(k, p_chem_ic_nh4_a03), &
 221                min_chem(k, p_chem_ic_cl_a01), min_chem(k, p_chem_ic_cl_a02), min_chem(k, p_chem_ic_cl_a03), &
 222                min_chem(k, p_chem_ic_na_a01), min_chem(k, p_chem_ic_na_a02), min_chem(k, p_chem_ic_na_a03), &
 223                min_chem(k, p_chem_ic_oin_a01), min_chem(k, p_chem_ic_oin_a02), min_chem(k, p_chem_ic_oin_a03)
 224          case ( 108 )     ! racm_soa_vbs_da
 225           write(unit=stats_unit2,fmt='(A,i4,39(f8.3,2i4))') 'k=',k, (min_chem(k,ic),ic=PARAM_FIRST_SCALAR,num_chem)
 226          end select
 227 #endif
 228       end if ! rootproc
 229
 230       um=minval(min_u(:)%value)
 231       vm=minval(min_v(:)%value)
 232       tm=minval(min_t(:)%value)
 233       pm=minval(min_p(:)%value)
 234       qm=minval(min_q(:)%value)
 235       if ( cloud_cv_options >= 1 ) then
 236          qcwm=minval(min_qcw(:)%value)
 237          qrnm=minval(min_qrn(:)%value)
 238       end if
 239       if ( cloud_cv_options >= 2 ) then
 240          qcim=minval(min_qci(:)%value)
 241          qsnm=minval(min_qsn(:)%value)
 242          qgrm=minval(min_qgr(:)%value)
 243       end if
 244    end do
 245
 246 #if (WRF_CHEM == 1)
 247    do ic=PARAM_FIRST_SCALAR, num_chem
 248       chemm(ic)=minval(min_chem(:,ic)%value)
 249    end do
 250 #endif
 251
 252    if (rootproc) then
 253       select case ( cloud_cv_options )
 254       case ( 0 )
 255          write(unit=stats_unit, fmt='(a,4(f12.4,10x),e12.4)') ' ALL', &
 256             um, vm, tm, pm, qm
 257       case ( 1 )
 258          write(unit=stats_unit, fmt='(a,4(f12.4,10x),3(e12.4,10x))') ' ALL', &
 259             um, vm, tm, pm, qm, qcwm, qrnm
 260       case ( 2, 3 )
 261          write(unit=stats_unit, fmt='(a,4(f12.4,10x),6(e12.4,10x))') ' ALL', &
 262             um, vm, tm, pm, qm, qcwm, qrnm, qcim, qsnm, qgrm
 263       end select
 264
 265 #if (WRF_CHEM == 1)
 266       select case ( chem_cv_options )
 267       case ( 10 )
 268           write(unit=stats_unit2, fmt='(a,4(f12.4,10x),6(e12.4,10x))') ' ALL', &
 269              chemm(p_chem_ic_p25), chemm(p_chem_ic_sulf), chemm(p_chem_ic_bc1), chemm(p_chem_ic_bc2), chemm(p_chem_ic_oc1), &
 270              chemm(p_chem_ic_oc2), chemm(p_chem_ic_dust_1), chemm(p_chem_ic_dust_2), chemm(p_chem_ic_seas_1), chemm(p_chem_ic_seas_2)
 271       case ( 20 )
 272           write(unit=stats_unit2, fmt='(a,24(f12.4,10x))') ' ALL', &
 273              chemm(p_chem_ic_bc_a01), chemm(p_chem_ic_bc_a02), chemm(p_chem_ic_bc_a03), &
 274              chemm(p_chem_ic_oc_a01), chemm(p_chem_ic_oc_a02), chemm(p_chem_ic_oc_a03), &
 275              chemm(p_chem_ic_so4_a01), chemm(p_chem_ic_so4_a02), chemm(p_chem_ic_so4_a03), &
 276              chemm(p_chem_ic_no3_a01), chemm(p_chem_ic_no3_a02), chemm(p_chem_ic_no3_a03), &
 277              chemm(p_chem_ic_nh4_a01), chemm(p_chem_ic_nh4_a02), chemm(p_chem_ic_nh4_a03), &
 278              chemm(p_chem_ic_cl_a01), chemm(p_chem_ic_cl_a02), chemm(p_chem_ic_cl_a03), &
 279              chemm(p_chem_ic_na_a01), chemm(p_chem_ic_na_a02), chemm(p_chem_ic_na_a03), &
 280              chemm(p_chem_ic_oin_a01), chemm(p_chem_ic_oin_a02), chemm(p_chem_ic_oin_a03)
 281      case ( 108 )     ! racm_soa_vbs_da
 282           write(unit=stats_unit2,fmt='(A,39f8.3)') 'Column_Min', (chemm(ic),ic=PARAM_FIRST_SCALAR,num_chem)
 283       end select
 284 #endif
 285
 286       write(unit=stats_unit, fmt='(/a/)') &
 287          ' Maximum of gridded analysis increments'
 288
 289 #if (WRF_CHEM == 1)
 290       write(unit=stats_unit2, fmt='(/a/)') ' Maximum of gridded analysis increments'
 291       select case ( chem_cv_options )
 292        case ( 10 )
 293           write(unit=stats_unit2, fmt='(11a/)') &
 294                ' Lvl         ', &
 295                'p25      i    j          ', &
 296                'sulf     i    j          ', &
 297                'bc1      i    j          ', &
 298                'bc2      i    j          ', &
 299                'oc1      i    j          ', &
 300                'oc2      i    j          ', &
 301                'dust_1   i    j          ', &
 302                'dust_2   i    j          ', &
 303                'seas_1   i    j          ', &
 304                'seas_2   i    j'
 305       end select
 306 #endif
 307
 308       select case ( cloud_cv_options )
 309       case ( 0 )
 310          write(unit=stats_unit, fmt='(6a/)') &
 311             ' Lvl         ', &
 312             'u     i    j          ', &
 313             'v     i    j          ', &
 314             't     i    j          ', &
 315             'p     i    j          ', &
 316             'q     i    j'
 317       case ( 1 )
 318          write(unit=stats_unit, fmt='(8a/)') &
 319             ' Lvl         ', &
 320             'u     i    j          ', &
 321             'v     i    j          ', &
 322             't     i    j          ', &
 323             'p     i    j          ', &
 324             'q     i    j          ', &
 325             'qcw   i    j          ', &
 326             'qrn   i    j'
 327       case ( 2, 3 )
 328          write(unit=stats_unit, fmt='(11a/)') &
 329             ' Lvl         ', &
 330             'u     i    j          ', &
 331             'v     i    j          ', &
 332             't     i    j          ', &
 333             'p     i    j          ', &
 334             'q     i    j          ', &
 335             'qcw   i    j          ', &
 336             'qrn   i    j          ', &
 337             'qci   i    j          ', &
 338             'qsn   i    j          ', &
 339             'qgr   i    j'
 340       end select
 341    end if !rootproc
 342
 343    do k = kts, kte
 344       if (rootproc) then
 345          if ( abs(max_q(k)%value) < 1.e-30 ) max_q(k)%value = 0.0
 346          select case ( cloud_cv_options )
 347          case ( 0 )
 348             write(unit=stats_unit, fmt='(i4,4(f12.4,2i5),e12.4,2i5)') k,    &
 349                max_u(k), max_v(k), max_t(k), max_p(k), max_q(k)
 350          case ( 1 )
 351             write(unit=stats_unit, fmt='(i4,4(f12.4,2i5),3(e12.4,2i5))') k, &
 352                max_u(k), max_v(k), max_t(k), max_p(k), max_q(k),            &
 353                max_qcw(k), max_qrn(k)
 354          case ( 2, 3 )
 355             write(unit=stats_unit, fmt='(i4,4(f12.4,2i5),6(e12.4,2i5))') k, &
 356                max_u(k), max_v(k), max_t(k), max_p(k), max_q(k),            &
 357                max_qcw(k), max_qrn(k), max_qci(k), max_qsn(k), max_qgr(k)
 358          end select
 359
 360 #if (WRF_CHEM == 1)
 361          select case ( chem_cv_options )
 362          case ( 10 )
 363             write(unit=stats_unit2, fmt='(i4,4(f12.4,2i5),6(e12.4,2i5))') k, &
 364                max_chem(k, p_chem_ic_p25), max_chem(k, p_chem_ic_sulf), &
 365                max_chem(k, p_chem_ic_bc1), max_chem(k, p_chem_ic_bc2), max_chem(k, p_chem_ic_oc1), max_chem(k, p_chem_ic_oc2), &
 366                max_chem(k, p_chem_ic_dust_1), max_chem(k, p_chem_ic_dust_2), max_chem(k, p_chem_ic_seas_1), max_chem(k, p_chem_ic_seas_1)
 367          case ( 20 )
 368             write(unit=stats_unit2, fmt='(i4,24(f12.4,2i5))') k, &
 369                max_chem(k, p_chem_ic_bc_a01), max_chem(k, p_chem_ic_bc_a02), max_chem(k, p_chem_ic_bc_a03), &
 370                max_chem(k, p_chem_ic_oc_a01), max_chem(k, p_chem_ic_oc_a02), max_chem(k, p_chem_ic_oc_a03), &
 371                max_chem(k, p_chem_ic_so4_a01), max_chem(k, p_chem_ic_so4_a02), max_chem(k, p_chem_ic_so4_a03), &
 372                max_chem(k, p_chem_ic_no3_a01), max_chem(k, p_chem_ic_no3_a02), max_chem(k, p_chem_ic_no3_a03), &
 373                max_chem(k, p_chem_ic_nh4_a01), max_chem(k, p_chem_ic_nh4_a02), max_chem(k, p_chem_ic_nh4_a03), &
 374                max_chem(k, p_chem_ic_cl_a01), max_chem(k, p_chem_ic_cl_a02), max_chem(k, p_chem_ic_cl_a03), &
 375                max_chem(k, p_chem_ic_na_a01), max_chem(k, p_chem_ic_na_a02), max_chem(k, p_chem_ic_na_a03), &
 376                max_chem(k, p_chem_ic_oin_a01), max_chem(k, p_chem_ic_oin_a02), max_chem(k, p_chem_ic_oin_a03)
 377          case ( 108 )     ! racm_soa_vbs_da
 378             write(unit=stats_unit2,fmt='(A,i4,39(f8.3,2i4))') 'k=',k, (max_chem(k,ic),ic=PARAM_FIRST_SCALAR,num_chem)
 379          end select
 380 #endif
 381       end if
 382
 383       um=maxval(max_u(:)%value)
 384       vm=maxval(max_v(:)%value)
 385       tm=maxval(max_t(:)%value)
 386       pm=maxval(max_p(:)%value)
 387       qm=maxval(max_q(:)%value)
 388       if ( cloud_cv_options >= 1 ) then
 389          qcwm=maxval(max_qcw(:)%value)
 390          qrnm=maxval(max_qrn(:)%value)
 391       end if
 392       if ( cloud_cv_options >= 2 ) then
 393          qcim=maxval(max_qci(:)%value)
 394          qsnm=maxval(max_qsn(:)%value)
 395          qgrm=maxval(max_qgr(:)%value)
 396       end if
 397    end do
 398
 399 #if (WRF_CHEM == 1)
 400    do ic=PARAM_FIRST_SCALAR, num_chem
 401       chemm(ic)=maxval(max_chem(:,ic)%value)
 402    end do
 403 #endif
 404
 405    if (rootproc) then
 406       select case ( cloud_cv_options )
 407       case ( 0 )
 408         write(unit=stats_unit, fmt='(a,4(f12.4,10x),e12.4)') ' ALL', &
 409            um, vm, tm, pm, qm
 410       case ( 1 )
 411         write(unit=stats_unit, fmt='(a,4(f12.4,10x),3(e12.4,10x))') ' ALL', &
 412            um, vm, tm, pm, qm, qcwm, qrnm
 413       case ( 2, 3 )
 414         write(unit=stats_unit, fmt='(a,4(f12.4,10x),6(e12.4,10x))') ' ALL', &
 415            um, vm, tm, pm, qm, qcwm, qrnm, qcim, qsnm, qgrm
 416       end select
 417
 418 #if (WRF_CHEM == 1)
 419       select case ( chem_cv_options )
 420       case ( 10 )
 421           write(unit=stats_unit2, fmt='(a,4(f12.4,10x),6(e12.4,10x))') ' ALL', &
 422              chemm(p_chem_ic_p25), chemm(p_chem_ic_sulf), chemm(p_chem_ic_bc1), chemm(p_chem_ic_bc2), chemm(p_chem_ic_oc1), &
 423              chemm(p_chem_ic_oc2), chemm(p_chem_ic_dust_1), chemm(p_chem_ic_dust_2), chemm(p_chem_ic_seas_1), chemm(p_chem_ic_seas_2)
 424       case ( 20 )
 425           write(unit=stats_unit2, fmt='(a,24(f12.4,10x))') ' ALL', &
 426              chemm(p_chem_ic_bc_a01), chemm(p_chem_ic_bc_a02), chemm(p_chem_ic_bc_a03),  &
 427              chemm(p_chem_ic_oc_a01), chemm(p_chem_ic_oc_a02), chemm(p_chem_ic_oc_a03),  &
 428              chemm(p_chem_ic_so4_a01), chemm(p_chem_ic_so4_a02), chemm(p_chem_ic_so4_a03),  &
 429              chemm(p_chem_ic_no3_a01), chemm(p_chem_ic_no3_a02), chemm(p_chem_ic_no3_a03),  &
 430              chemm(p_chem_ic_nh4_a01), chemm(p_chem_ic_nh4_a02), chemm(p_chem_ic_nh4_a03),  &
 431              chemm(p_chem_ic_cl_a01), chemm(p_chem_ic_cl_a02), chemm(p_chem_ic_cl_a03),  &
 432              chemm(p_chem_ic_na_a01), chemm(p_chem_ic_na_a02), chemm(p_chem_ic_na_a03),  &
 433              chemm(p_chem_ic_oin_a01), chemm(p_chem_ic_oin_a02), chemm(p_chem_ic_oin_a03)
 434       case ( 108 )     ! racm_soa_vbs_da
 435           write(unit=stats_unit2,fmt='(A,39f8.3)') 'Column_Max', (chemm(ic),ic=PARAM_FIRST_SCALAR,num_chem)
 436       end select
 437 #endif
 438
 439       write(unit=stats_unit, fmt='(/a/)') ' Mean of gridded analysis increments'
 440
 441 #if (WRF_CHEM == 1)
 442       write(unit=stats_unit2, fmt='(/a/)') ' Mean of gridded analysis increments'
 443       select case ( chem_cv_options )
 444        case ( 10 )
 445           write(unit=stats_unit2, fmt='(a/)') &
 446             ' Lvl        p25           sulf           bc1           bc2           oc1           oc2         dust_1         dust_2         seas_1         seas_2'
 447       end select
 448 #endif
 449
 450       select case ( cloud_cv_options )
 451       case ( 0 )
 452          write(unit=stats_unit, fmt='(a/)') &
 453             ' Lvl        u           v           t           p           q'
 454       case ( 1 )
 455          write(unit=stats_unit, fmt='(a/)') &
 456             ' Lvl        u           v           t           p           q           qcw         qrn'
 457       case ( 2, 3 )
 458          write(unit=stats_unit, fmt='(a/)') &
 459             ' Lvl        u           v           t           p           q           qcw         qrn         qci         qsn         qgr'
 460       end select
 461    end if !rootproc
 462
 463    um = 0.0
 464    vm = 0.0
 465    tm = 0.0
 466    pm = 0.0
 467    qm = 0.0
 468
 469 #if (WRF_CHEM == 1)
 470    chemm = 0.0
 471    if ( chem_cv_options >= 10 ) then
 472    do ic=PARAM_FIRST_SCALAR, num_chem
 473      do k = kts, kte
 474         chemv(k,ic) = sum(grid%xachem%chem_ic(its:ite,jts:jte,k,ic))
 475      end do
 476    call da_proc_sum_real (chemv(:,ic))
 477    end do
 478    end if
 479 #endif
 480
 481    do k = kts, kte
 482       uv(k) = sum(grid%xa%u(its:ite,jts:jte,k))
 483       vv(k) = sum(grid%xa%v(its:ite,jts:jte,k))
 484       tv(k) = sum(grid%xa%t(its:ite,jts:jte,k))
 485       pv(k) = sum(grid%xa%p(its:ite,jts:jte,k))
 486       qv(k) = sum(grid%xa%q(its:ite,jts:jte,k))
 487    end do
 488    call da_proc_sum_real (uv)
 489    call da_proc_sum_real (vv)
 490    call da_proc_sum_real (tv)
 491    call da_proc_sum_real (pv)
 492    call da_proc_sum_real (qv)
 493
 494    if ( cloud_cv_options >= 1 ) then
 495       qcwm = 0.0
 496       qrnm = 0.0
 497        do k = kts, kte
 498           qcwv(k) = sum(grid%xa%qcw(its:ite,jts:jte,k))
 499           qrnv(k) = sum(grid%xa%qrn(its:ite,jts:jte,k))
 500        end do
 501       call da_proc_sum_real (qcwv)
 502       call da_proc_sum_real (qrnv)
 503    end if
 504    if ( cloud_cv_options >= 2 ) then
 505       qcim = 0.0
 506       qsnm = 0.0
 507       qgrm = 0.0
 508        do k = kts, kte
 509           qciv(k) = sum(grid%xa%qci(its:ite,jts:jte,k))
 510           qsnv(k) = sum(grid%xa%qsn(its:ite,jts:jte,k))
 511           qgrv(k) = sum(grid%xa%qgr(its:ite,jts:jte,k))
 512        end do
 513       call da_proc_sum_real (qciv)
 514       call da_proc_sum_real (qsnv)
 515       call da_proc_sum_real (qgrv)
 516    end if
 517
 518    if (rootproc) then
 519       do k = kts, kte
 520          select case ( cloud_cv_options )
 521          case ( 0 )
 522             write(unit=stats_unit, fmt='(i4,4f12.4,e12.4)')  k, &
 523                uv(k)*rij_g, vv(k)*rij_g, tv(k)*rij_g, &
 524                pv(k)*rij_g, qv(k)*rij_g
 525          case ( 1 )
 526             write(unit=stats_unit, fmt='(i4,4f12.4,3e12.4)') k, &
 527                uv(k)*rij_g, vv(k)*rij_g, tv(k)*rij_g, &
 528                pv(k)*rij_g, qv(k)*rij_g,                        &
 529                qcwv(k)*rij_g, qrnv(k)*rij_g
 530          case ( 2, 3 )
 531             write(unit=stats_unit, fmt='(i4,4f12.4,6e12.4)') k, &
 532                uv(k)*rij_g, vv(k)*rij_g, tv(k)*rij_g, &
 533                pv(k)*rij_g, qv(k)*rij_g,                        &
 534                qcwv(k)*rij_g, qrnv(k)*rij_g, qciv(k)*rij_g,     &
 535                qsnv(k)*rij_g, qgrv(k)*rij_g
 536          end select
 537
 538 #if (WRF_CHEM == 1)
 539          select case ( chem_cv_options )
 540          case ( 10 )
 541             write(unit=stats_unit2, fmt='(i4,10e14.4)') k, &
 542                chemv(k, p_chem_ic_p25)*rij_g, chemv(k, p_chem_ic_sulf)*rij_g, &
 543                chemv(k, p_chem_ic_bc1)*rij_g, chemv(k, p_chem_ic_bc2)*rij_g, chemv(k, p_chem_ic_oc1)*rij_g, chemv(k, p_chem_ic_oc2)*rij_g, &
 544                chemv(k, p_chem_ic_dust_1)*rij_g, chemv(k, p_chem_ic_dust_2)*rij_g, chemv(k, p_chem_ic_seas_2)*rij_g, chemv(k, p_chem_ic_seas_2)*rij_g
 545          case ( 20 )
 546             write(unit=stats_unit2, fmt='(i4,24e14.4)') k, &
 547                chemv(k, p_chem_ic_bc_a01)*rij_g, chemv(k, p_chem_ic_bc_a02)*rij_g, chemv(k, p_chem_ic_bc_a03)*rij_g,  &
 548                chemv(k, p_chem_ic_oc_a01)*rij_g, chemv(k, p_chem_ic_oc_a02)*rij_g, chemv(k, p_chem_ic_oc_a03)*rij_g,  &
 549                chemv(k, p_chem_ic_so4_a01)*rij_g, chemv(k, p_chem_ic_so4_a02)*rij_g, chemv(k, p_chem_ic_so4_a03)*rij_g,  &
 550                chemv(k, p_chem_ic_no3_a01)*rij_g, chemv(k, p_chem_ic_no3_a02)*rij_g, chemv(k, p_chem_ic_no3_a03)*rij_g,  &
 551                chemv(k, p_chem_ic_nh4_a01)*rij_g, chemv(k, p_chem_ic_nh4_a02)*rij_g, chemv(k, p_chem_ic_nh4_a03)*rij_g,  &
 552                chemv(k, p_chem_ic_cl_a01)*rij_g, chemv(k, p_chem_ic_cl_a02)*rij_g, chemv(k, p_chem_ic_cl_a03)*rij_g,  &
 553                chemv(k, p_chem_ic_na_a01)*rij_g, chemv(k, p_chem_ic_na_a02)*rij_g, chemv(k, p_chem_ic_na_a03)*rij_g,  &
 554                chemv(k, p_chem_ic_oin_a01)*rij_g, chemv(k, p_chem_ic_oin_a02)*rij_g, chemv(k, p_chem_ic_oin_a03)*rij_g
 555          case ( 108 )     ! racm_soa_vbs_da
 556           write(unit=stats_unit2,fmt='(A,i4,39f8.3)') 'k=',k, &
 557               (chemv(k,ic)*rij_g, ic=PARAM_FIRST_SCALAR,num_chem)
 558          end select
 559
 560          do ic=PARAM_FIRST_SCALAR, num_chem
 561            chemm(ic)=chemm(ic)+chemv(k,ic)
 562          end do
 563 #endif
 564
 565          um=um+uv(k)
 566          vm=vm+vv(k)
 567          tm=tm+tv(k)
 568          pm=pm+pv(k)
 569          qm=qm+qv(k)
 570          if ( cloud_cv_options >= 1 ) then
 571             qcwm = qcwm + qcwv(k)
 572             qrnm = qrnm + qrnv(k)
 573          end if
 574          if ( cloud_cv_options >= 2 ) then
 575             qcim = qcim + qciv(k)
 576             qsnm = qsnm + qsnv(k)
 577             qgrm = qgrm + qgrv(k)
 578          end if
 579       end do !k loop
 580    end if !rootproc
 581
 582    if (rootproc) then
 583       select case ( cloud_cv_options )
 584       case ( 0 )
 585          write(unit=stats_unit, fmt='(a,4f12.4,e12.4)')  ' ALL',   &
 586             um*rijk_g, vm*rijk_g, tm*rijk_g, pm*rijk_g, qm*rijk_g
 587       case ( 1 )
 588          write(unit=stats_unit, fmt='(a,4f12.4,3e12.4)') ' ALL',   &
 589             um*rijk_g, vm*rijk_g, tm*rijk_g, pm*rijk_g, qm*rijk_g, &
 590             qcwm*rijk_g, qrnm*rijk_g
 591       case ( 2, 3 )
 592          write(unit=stats_unit, fmt='(a,4f12.4,6e12.4)') ' ALL',   &
 593             um*rijk_g, vm*rijk_g, tm*rijk_g, pm*rijk_g, qm*rijk_g, &
 594             qcwm*rijk_g, qrnm*rijk_g, qcim*rijk_g, qsnm*rijk_g, qgrm*rijk_g
 595       end select
 596
 597 #if (WRF_CHEM == 1)
 598       select case ( chem_cv_options )
 599          case ( 10 )
 600          write(unit=stats_unit2, fmt='(a,10e14.4)') ' ALL',   &
 601             chemm(p_chem_ic_p25)*rijk_g, chemm(p_chem_ic_sulf)*rijk_g, &
 602             chemm(p_chem_ic_bc1)*rijk_g, chemm(p_chem_ic_bc2)*rijk_g, chemm(p_chem_ic_oc1)*rijk_g, chemm(p_chem_ic_oc2)*rijk_g, &
 603             chemm(p_chem_ic_dust_1)*rijk_g, chemm(p_chem_ic_dust_2)*rijk_g, chemm(p_chem_ic_seas_2)*rijk_g, chemm(p_chem_ic_seas_2)*rijk_g
 604          case ( 20 )
 605          write(unit=stats_unit2, fmt='(a,24e14.4)') ' ALL',   &
 606             chemm(p_chem_ic_bc_a01)*rijk_g, chemm(p_chem_ic_bc_a02)*rijk_g, chemm(p_chem_ic_bc_a03)*rijk_g, &
 607             chemm(p_chem_ic_oc_a01)*rijk_g, chemm(p_chem_ic_oc_a02)*rijk_g, chemm(p_chem_ic_oc_a03)*rijk_g, &
 608             chemm(p_chem_ic_so4_a01)*rijk_g, chemm(p_chem_ic_so4_a02)*rijk_g, chemm(p_chem_ic_so4_a03)*rijk_g, &
 609             chemm(p_chem_ic_no3_a01)*rijk_g, chemm(p_chem_ic_no3_a02)*rijk_g, chemm(p_chem_ic_no3_a03)*rijk_g, &
 610             chemm(p_chem_ic_nh4_a01)*rijk_g, chemm(p_chem_ic_nh4_a02)*rijk_g, chemm(p_chem_ic_nh4_a03)*rijk_g, &
 611             chemm(p_chem_ic_cl_a01)*rijk_g, chemm(p_chem_ic_cl_a02)*rijk_g, chemm(p_chem_ic_cl_a03)*rijk_g, &
 612             chemm(p_chem_ic_na_a01)*rijk_g, chemm(p_chem_ic_na_a02)*rijk_g, chemm(p_chem_ic_na_a03)*rijk_g, &
 613             chemm(p_chem_ic_oin_a01)*rijk_g, chemm(p_chem_ic_oin_a02)*rijk_g, chemm(p_chem_ic_oin_a03)*rijk_g
 614          case ( 108 )     ! racm_soa_vbs_da
 615           write(unit=stats_unit2,fmt='(A,39f8.3)') 'ALL', (chemm(ic)*rijk_g,ic=PARAM_FIRST_SCALAR,num_chem)
 616       end select
 617 #endif
 618
 619
 620       write(unit=stats_unit, fmt='(/a/)') ' RMSE of gridded analysis increments'
 621
 622 #if (WRF_CHEM == 1)
 623       write(unit=stats_unit2, fmt='(/a/)') ' RMSE of gridded analysis increments'
 624       select case ( chem_cv_options )
 625        case ( 10 )
 626           write(unit=stats_unit2, fmt='(a/)') &
 627             ' Lvl        p25           sulf           bc1           bc2           oc1           oc2         dust_1         dust_2         seas_1         seas_2'
 628       end select
 629 #endif
 630
 631       select case ( cloud_cv_options )
 632       case ( 0 )
 633          write(unit=stats_unit, fmt='(a/)') &
 634             ' Lvl        u           v           t           p           q'
 635       case ( 1 )
 636          write(unit=stats_unit, fmt='(a/)') &
 637             ' Lvl        u           v           t           p           q           qcw         qrn'
 638       case ( 2, 3 )
 639          write(unit=stats_unit, fmt='(a/)') &
 640             ' Lvl        u           v           t           p           q           qcw         qrn         qci          qsn         qgr'
 641       end select
 642    end if !rootproc
 643
 644    um = 0.0
 645    vm = 0.0
 646    tm = 0.0
 647    pm = 0.0
 648    qm = 0.0
 649    uv = 0.0
 650    vv = 0.0
 651    tv = 0.0
 652    pv = 0.0
 653    qv = 0.0
 654
 655 #if (WRF_CHEM == 1)
 656    chemm = 0.0
 657    chemv = 0.0
 658 #endif
 659
 660    qcwv = 0.0
 661    qrnv = 0.0
 662    qciv = 0.0
 663    qsnv = 0.0
 664    qgrv = 0.0
 665    qcwm = 0.0
 666    qrnm = 0.0
 667    qcim = 0.0
 668    qsnm = 0.0
 669    qgrm = 0.0
 670
 671    do k = kts, kte
 672       do j=jts,jte
 673          do i=its,ite
 674             uv(k) = uv(k) + grid%xa%u(i,j,k) * grid%xa%u(i,j,k)
 675             vv(k) = vv(k) + grid%xa%v(i,j,k) * grid%xa%v(i,j,k)
 676             tv(k) = tv(k) + grid%xa%t(i,j,k) * grid%xa%t(i,j,k)
 677             pv(k) = pv(k) + grid%xa%p(i,j,k) * grid%xa%p(i,j,k)
 678             qv(k) = qv(k) + grid%xa%q(i,j,k) * grid%xa%q(i,j,k)
 679          end do
 680       end do
 681    end do
 682    call da_proc_sum_real (uv)
 683    call da_proc_sum_real (vv)
 684    call da_proc_sum_real (tv)
 685    call da_proc_sum_real (pv)
 686    call da_proc_sum_real (qv)
 687
 688 #if (WRF_CHEM == 1)
 689    do ic=PARAM_FIRST_SCALAR, num_chem
 690    do k = kts, kte
 691       do j=jts,jte
 692          do i=its,ite
 693             chemv(k,ic) = chemv(k,ic) + grid%xachem%chem_ic(i,j,k,ic) * grid%xachem%chem_ic(i,j,k,ic)
 694          end do
 695       end do
 696    call da_proc_sum_real (chemv(:,ic))
 697    end do
 698    end do
 699 #endif
 700
 701    if ( cloud_cv_options >= 1 ) then
 702       do k = kts, kte
 703          do j=jts,jte
 704             do i=its,ite
 705                qcwv(k) = qcwv(k) + grid%xa%qcw(i,j,k) * grid%xa%qcw(i,j,k)
 706                qrnv(k) = qrnv(k) + grid%xa%qrn(i,j,k) * grid%xa%qrn(i,j,k)
 707             end do
 708          end do
 709       end do
 710       call da_proc_sum_real (qcwv)
 711       call da_proc_sum_real (qrnv)
 712    end if
 713
 714    if ( cloud_cv_options >= 2 ) then
 715       do k = kts, kte
 716          do j=jts,jte
 717             do i=its,ite
 718                qciv(k) = qciv(k) + grid%xa%qci(i,j,k) * grid%xa%qci(i,j,k)
 719                qsnv(k) = qsnv(k) + grid%xa%qsn(i,j,k) * grid%xa%qsn(i,j,k)
 720                qgrv(k) = qgrv(k) + grid%xa%qgr(i,j,k) * grid%xa%qgr(i,j,k)
 721             end do
 722          end do
 723       end do
 724       call da_proc_sum_real (qciv)
 725       call da_proc_sum_real (qsnv)
 726       call da_proc_sum_real (qgrv)
 727    end if
 728
 729    if (rootproc) then
 730       do k = kts, kte
 731          select case ( cloud_cv_options )
 732          case ( 0 )
 733             write(unit=stats_unit, fmt='(i4,4f12.4,e12.4)')  k, &
 734                sqrt(uv(k)*rij_g),   &
 735                sqrt(vv(k)*rij_g),   &
 736                sqrt(tv(k)*rij_g),   &
 737                sqrt(pv(k)*rij_g),   &
 738                sqrt(qv(k)*rij_g)
 739          case ( 1 )
 740             write(unit=stats_unit, fmt='(i4,4f12.4,3e12.4)') k, &
 741                sqrt(uv(k)*rij_g),   &
 742                sqrt(vv(k)*rij_g),   &
 743                sqrt(tv(k)*rij_g),   &
 744                sqrt(pv(k)*rij_g),   &
 745                sqrt(qv(k)*rij_g),   &
 746                sqrt(qcwv(k)*rij_g), &
 747                sqrt(qrnv(k)*rij_g)
 748          case ( 2, 3 )
 749             write(unit=stats_unit, fmt='(i4,4f12.4,6e12.4)') k, &
 750                sqrt(uv(k)*rij_g),   &
 751                sqrt(vv(k)*rij_g),   &
 752                sqrt(tv(k)*rij_g),   &
 753                sqrt(pv(k)*rij_g),   &
 754                sqrt(qv(k)*rij_g),   &
 755                sqrt(qcwv(k)*rij_g), &
 756                sqrt(qrnv(k)*rij_g), &
 757                sqrt(qciv(k)*rij_g), &
 758                sqrt(qsnv(k)*rij_g), &
 759                sqrt(qgrv(k)*rij_g)
 760          end select
 761
 762 #if (WRF_CHEM == 1)
 763          select case ( chem_cv_options )
 764          case ( 10 )
 765             write(unit=stats_unit2, fmt='(i4,10e14.4)') k, &
 766                sqrt(chemv(k, p_chem_ic_p25)*rij_g),   &
 767                sqrt(chemv(k, p_chem_ic_sulf)*rij_g),   &
 768                sqrt(chemv(k, p_chem_ic_bc1)*rij_g),   &
 769                sqrt(chemv(k, p_chem_ic_bc2)*rij_g),   &
 770                sqrt(chemv(k, p_chem_ic_oc1)*rij_g),   &
 771                sqrt(chemv(k, p_chem_ic_oc2)*rij_g), &
 772                sqrt(chemv(k, p_chem_ic_dust_1)*rij_g), &
 773                sqrt(chemv(k, p_chem_ic_dust_2)*rij_g), &
 774                sqrt(chemv(k, p_chem_ic_seas_1)*rij_g), &
 775                sqrt(chemv(k, p_chem_ic_seas_2)*rij_g)
 776          case ( 20 )
 777             write(unit=stats_unit2, fmt='(i4,24e14.4)') k, &
 778                sqrt(chemv(k, p_chem_ic_bc_a01)*rij_g), sqrt(chemv(k, p_chem_ic_bc_a02)*rij_g), sqrt(chemv(k, p_chem_ic_bc_a03)*rij_g),  &
 779                sqrt(chemv(k, p_chem_ic_oc_a01)*rij_g), sqrt(chemv(k, p_chem_ic_oc_a02)*rij_g), sqrt(chemv(k, p_chem_ic_oc_a03)*rij_g),  &
 780                sqrt(chemv(k, p_chem_ic_so4_a01)*rij_g), sqrt(chemv(k, p_chem_ic_so4_a02)*rij_g), sqrt(chemv(k, p_chem_ic_so4_a03)*rij_g),  &
 781                sqrt(chemv(k, p_chem_ic_no3_a01)*rij_g), sqrt(chemv(k, p_chem_ic_no3_a02)*rij_g), sqrt(chemv(k, p_chem_ic_no3_a03)*rij_g),  &
 782                sqrt(chemv(k, p_chem_ic_nh4_a01)*rij_g), sqrt(chemv(k, p_chem_ic_nh4_a02)*rij_g), sqrt(chemv(k, p_chem_ic_nh4_a03)*rij_g),  &
 783                sqrt(chemv(k, p_chem_ic_cl_a01)*rij_g), sqrt(chemv(k, p_chem_ic_cl_a02)*rij_g), sqrt(chemv(k, p_chem_ic_cl_a03)*rij_g),  &
 784                sqrt(chemv(k, p_chem_ic_na_a01)*rij_g), sqrt(chemv(k, p_chem_ic_na_a02)*rij_g), sqrt(chemv(k, p_chem_ic_na_a03)*rij_g),  &
 785                sqrt(chemv(k, p_chem_ic_oin_a01)*rij_g), sqrt(chemv(k, p_chem_ic_oin_a02)*rij_g), sqrt(chemv(k, p_chem_ic_oin_a03)*rij_g)
 786          case ( 108 )     ! racm_soa_vbs_da
 787             write(unit=stats_unit2,fmt='(A,i4,39f8.3)') 'k=',k, &
 788               (sqrt(chemv(k,ic)*rij_g), ic=PARAM_FIRST_SCALAR,num_chem)
 789          end select
 790          do ic=PARAM_FIRST_SCALAR, num_chem
 791            chemm(ic)=chemm(ic)+chemv(k,ic)
 792          end do
 793 #endif
 794
 795          um=um+uv(k)
 796          vm=vm+vv(k)
 797          tm=tm+tv(k)
 798          pm=pm+pv(k)
 799          qm=qm+qv(k)
 800          if ( cloud_cv_options >= 1 ) then
 801             qcwm=qcwm+qcwv(k)
 802             qrnm=qrnm+qrnv(k)
 803          end if
 804          if ( cloud_cv_options >= 2 ) then
 805             qcim=qcim+qciv(k)
 806             qsnm=qsnm+qsnv(k)
 807          end if
 808       end do !k loop
 809    end if !rootproc
 810
 811    if (rootproc) then
 812       select case ( cloud_cv_options )
 813       case ( 0 )
 814          write(unit=stats_unit, fmt='(a,4f12.4,e12.4)')  ' ALL', &
 815             sqrt(um*rijk_g), sqrt(vm*rijk_g), sqrt(tm*rijk_g),   &
 816             sqrt(pm*rijk_g), sqrt(qm*rijk_g)
 817       case ( 1 )
 818          write(unit=stats_unit, fmt='(a,4f12.4,3e12.4)') ' ALL', &
 819             sqrt(um*rijk_g), sqrt(vm*rijk_g), sqrt(tm*rijk_g),   &
 820             sqrt(pm*rijk_g), sqrt(qm*rijk_g),                    &
 821             sqrt(qcwm*rijk_g), sqrt(qrnm*rijk_g)
 822       case ( 2, 3 )
 823          write(unit=stats_unit, fmt='(a,4f12.4,6e12.4)') ' ALL', &
 824             sqrt(um*rijk_g), sqrt(vm*rijk_g), sqrt(tm*rijk_g),   &
 825             sqrt(pm*rijk_g), sqrt(qm*rijk_g),                    &
 826             sqrt(qcwm*rijk_g), sqrt(qrnm*rijk_g), sqrt(qcim*rijk_g), &
 827             sqrt(qsnm*rijk_g), sqrt(qgrm*rijk_g)
 828       end select
 829
 830 #if (WRF_CHEM == 1)
 831          select case ( chem_cv_options )
 832          case ( 10 )
 833          write(unit=stats_unit2, fmt='(a,10e14.4)') ' ALL', &
 834          sqrt(chemm(p_chem_ic_p25)*rijk_g), sqrt(chemm(p_chem_ic_sulf)*rijk_g),   &
 835          sqrt(chemm(p_chem_ic_bc1)*rijk_g), sqrt(chemm(p_chem_ic_bc2)*rijk_g), sqrt(chemm(p_chem_ic_oc1)*rijk_g), sqrt(chemm(p_chem_ic_oc2)*rijk_g), &
 836          sqrt(chemm(p_chem_ic_dust_1)*rijk_g), sqrt(chemm(p_chem_ic_dust_2)*rijk_g), &
 837          sqrt(chemm(p_chem_ic_seas_1)*rijk_g), sqrt(chemm(p_chem_ic_seas_2)*rijk_g)
 838          case ( 20 )
 839          write(unit=stats_unit2, fmt='(a,24e14.4)') ' ALL', &
 840          sqrt(chemm(p_chem_ic_bc_a01)*rijk_g), sqrt(chemm(p_chem_ic_bc_a02)*rijk_g), sqrt(chemm(p_chem_ic_bc_a03)*rijk_g), &
 841          sqrt(chemm(p_chem_ic_oc_a01)*rijk_g), sqrt(chemm(p_chem_ic_oc_a02)*rijk_g), sqrt(chemm(p_chem_ic_oc_a03)*rijk_g), &
 842          sqrt(chemm(p_chem_ic_so4_a01)*rijk_g), sqrt(chemm(p_chem_ic_so4_a02)*rijk_g), sqrt(chemm(p_chem_ic_so4_a03)*rijk_g), &
 843          sqrt(chemm(p_chem_ic_no3_a01)*rijk_g), sqrt(chemm(p_chem_ic_no3_a02)*rijk_g), sqrt(chemm(p_chem_ic_no3_a03)*rijk_g), &
 844          sqrt(chemm(p_chem_ic_nh4_a01)*rijk_g), sqrt(chemm(p_chem_ic_nh4_a02)*rijk_g), sqrt(chemm(p_chem_ic_nh4_a03)*rijk_g), &
 845          sqrt(chemm(p_chem_ic_cl_a01)*rijk_g), sqrt(chemm(p_chem_ic_cl_a02)*rijk_g), sqrt(chemm(p_chem_ic_cl_a03)*rijk_g), &
 846          sqrt(chemm(p_chem_ic_na_a01)*rijk_g), sqrt(chemm(p_chem_ic_na_a02)*rijk_g), sqrt(chemm(p_chem_ic_na_a03)*rijk_g), &
 847          sqrt(chemm(p_chem_ic_oin_a01)*rijk_g), sqrt(chemm(p_chem_ic_oin_a02)*rijk_g), sqrt(chemm(p_chem_ic_oin_a03)*rijk_g)
 848          case ( 108 )     ! racm_soa_vbs_da
 849           write(unit=stats_unit2,fmt='(A,39f8.3)') 'ALL', (sqrt(chemm(ic)*rijk_g),ic=PARAM_FIRST_SCALAR,num_chem)
 850          end select
 851 #endif
 852    end if
 853
 854    call da_trace_exit("da_analysis_stats")
 855
 856 contains
 857
 858 #include "da_maxmin_in_field.inc"
 859
 860 end subroutine da_analysis_stats