wrfv2_fire/var/da/da_pseudo/da_ao_stats_pseudo.inc

   1 subroutine da_ao_stats_pseudo  (stats_unit, iv, re)
   2
   3    !-----------------------------------------------------------------------
   4    ! Purpose: TBD
   5    !-----------------------------------------------------------------------
   6
   7    implicit none
   8
   9    integer,        intent(in)    :: stats_unit    ! Output unit for stats.
  10    type (iv_type), intent(inout) :: iv            ! O-B
  11    type  (y_type), intent(in)    :: re            ! O-A
  12
  13    type (stats_pseudo_type) :: stats
  14    integer                  :: nu, nv, nt, np, nq
  15    integer                  :: n
  16    real                     :: o_minus_b, o_minus_a, sigma_o, sigma_b
  17
  18    if (trace_use_dull) call da_trace_entry("da_ao_stats_pseudo")
  19
  20    nu = 0
  21    nv = 0
  22    nt = 0
  23    np = 0
  24    nq = 0
  25
  26    stats%maximum%u = maxmin_type (missing_r, 0, 0)
  27    stats%maximum%v = maxmin_type (missing_r, 0, 0)
  28    stats%maximum%t = maxmin_type (missing_r, 0, 0)
  29    stats%maximum%p = maxmin_type (missing_r, 0, 0)
  30    stats%maximum%q = maxmin_type (missing_r, 0, 0)
  31    stats%minimum%u = maxmin_type(-missing_r, 0, 0)
  32    stats%minimum%v = maxmin_type(-missing_r, 0, 0)
  33    stats%minimum%t = maxmin_type(-missing_r, 0, 0)
  34    stats%minimum%p = maxmin_type(-missing_r, 0, 0)
  35    stats%minimum%q = maxmin_type(-missing_r, 0, 0)
  36
  37    stats%average = residual_pseudo1_type(0.0, 0.0, 0.0, 0.0, 0.0)
  38    stats%rms_err = stats%average
  39
  40    do n=1, iv%info(pseudo)%nlocal
  41       if (iv%info(pseudo)%proc_domain(1,n)) then
  42          call da_stats_calculate (n, 0, iv%pseudo(n)%u%qc, &
  43             re%pseudo(n)%u, nu, &
  44             stats%minimum%u, stats%maximum%u, &
  45             stats%average%u, stats%rms_err%u)
  46          call da_stats_calculate (n, 0, iv%pseudo(n)%v%qc, &
  47             re%pseudo(n)%v, nv, &
  48             stats%minimum%v, stats%maximum%v, &
  49             stats%average%v, stats%rms_err%v)
  50          call da_stats_calculate (n, 0, iv%pseudo(n)%t%qc, &
  51             re%pseudo(n)%t, nt, &
  52             stats%minimum%t, stats%maximum%t, &
  53             stats%average%t, stats%rms_err%t)
  54          call da_stats_calculate (n, 0, iv%pseudo(n)%p%qc, &
  55             re%pseudo(n)%p, np, &
  56             stats%minimum%p, stats%maximum%p, &
  57             stats%average%p, stats%rms_err%p)
  58          call da_stats_calculate (n, 0, iv%pseudo(n)%q%qc, &
  59             re%pseudo(n)%q, nq, &
  60             stats%minimum%q, stats%maximum%q, &
  61             stats%average%q, stats%rms_err%q)
  62
  63          if (nu > 0) then
  64             o_minus_b = iv%pseudo(n)%u%inv
  65             o_minus_a = re%pseudo(n)%u
  66             sigma_o   = iv%pseudo(n)%u%error
  67          else if (nv > 0) then
  68             o_minus_b = iv%pseudo(n)%v%inv
  69             o_minus_a = re%pseudo(n)%v
  70             sigma_o   = iv%pseudo(n)%v%error
  71          else if (nt > 0) then
  72             o_minus_b = iv%pseudo(n)%t%inv
  73             o_minus_a = re%pseudo(n)%t
  74             sigma_o   = iv%pseudo(n)%t%error
  75          else if (np > 0) then
  76             o_minus_b = iv%pseudo(n)%p%inv
  77             o_minus_a = re%pseudo(n)%p
  78             sigma_o   = iv%pseudo(n)%p%error
  79          else if (nq > 0) then
  80             o_minus_b = iv%pseudo(n)%q%inv
  81             o_minus_a = re%pseudo(n)%q
  82             sigma_o   = iv%pseudo(n)%q%error
  83          end if
  84
  85          write(stats_unit,'(/a,a1,a,e15.5)')' pseudo ', pseudo_var, ' o-b: ', o_minus_b
  86          write(stats_unit,' (a,a1,a,e15.5)')' pseudo ', pseudo_var, ' o-a: ', o_minus_a
  87          write(stats_unit,' (a,a1,a,e15.5)')' pseudo ', pseudo_var, ' sigma_o: ', sigma_o
  88
  89          ! calculate equivalent sigma_b using o-a=(o-b)*sigma_o/(sigma_o+sigma_b)
  90          sigma_b = sqrt ((o_minus_b - o_minus_a) / o_minus_a) * sigma_o
  91          write(stats_unit,'(a,a1,a,e15.5)')' pseudo ', pseudo_var, ' sigma_b: ', sigma_b
  92       end if
  93    end do
  94
  95    ! do inter-processor communication to gather statistics.
  96    call da_proc_sum_int (nu)
  97    call da_proc_sum_int (nv)
  98    call da_proc_sum_int (nt)
  99    call da_proc_sum_int (np)
 100    call da_proc_sum_int (nq)
 101    iv%nstats(pseudo) = nu + nv + nt + np + nq
 102
 103    call da_proc_stats_combine(stats%average%u, stats%rms_err%u, &
 104       stats%minimum%u%value, stats%maximum%u%value, &
 105       stats%minimum%u%n, stats%maximum%u%n, &
 106       stats%minimum%u%l, stats%maximum%u%l)
 107    call da_proc_stats_combine(stats%average%v, stats%rms_err%v, &
 108       stats%minimum%v%value, stats%maximum%v%value, &
 109       stats%minimum%v%n, stats%maximum%v%n, &
 110       stats%minimum%v%l, stats%maximum%v%l)
 111    call da_proc_stats_combine(stats%average%t, stats%rms_err%t, &
 112       stats%minimum%t%value, stats%maximum%t%value, &
 113       stats%minimum%t%n, stats%maximum%t%n, &
 114       stats%minimum%t%l, stats%maximum%t%l)
 115    call da_proc_stats_combine(stats%average%p, stats%rms_err%p, &
 116       stats%minimum%p%value, stats%maximum%p%value, &
 117       stats%minimum%p%n, stats%maximum%p%n, &
 118       stats%minimum%p%l, stats%maximum%p%l)
 119    call da_proc_stats_combine(stats%average%q, stats%rms_err%q, &
 120       stats%minimum%q%value, stats%maximum%q%value, &
 121       stats%minimum%q%n, stats%maximum%q%n, &
 122       stats%minimum%q%l, stats%maximum%q%l)
 123
 124    if (rootproc) then
 125       if (nu /= 0 .or. nv /= 0 .OR. nt /= 0 .or. np /= 0 .OR. nq /= 0) then
 126          write(unit=stats_unit, fmt='(/a/)') ' O-A Diagnostics for pseudo'
 127          call da_print_stats_pseudo (stats_unit, nu, nv, nt, np, nq, stats)
 128       end if
 129    end if
 130
 131    if (trace_use_dull) call da_trace_exit("da_ao_stats_pseudo")
 132
 133 end subroutine da_ao_stats_pseudo
 134
 135