var/da/da_bogus/da_get_innov_vector_bogus.inc

   1 subroutine da_get_innov_vector_bogus(it,num_qcstat_conv, grid, ob, iv)
   2
   3    !------------------------------------------------------------------------------
   4    ! Purpose: calculate the innovations for the bogus data.
   5    !    Updated for Analysis on Arakawa-C grid
   6    !    Author: Syed RH Rizvi,  MMM/ESSL/NCAR,  Date: 10/22/2008
   7    !------------------------------------------------------------------------------
   8
   9    implicit none
  10
  11    integer,          intent(in)    :: it       ! External iteration.
  12    type(domain),     intent(in)    :: grid     ! first guess state.
  13    type(y_type),     intent(in)    :: ob       ! Observation structure.
  14    type(iv_type),    intent(inout) :: iv       ! O-B structure.
  15    integer,          intent(inout) :: num_qcstat_conv(:,:,:,:)
  16
  17
  18    integer :: n        ! Loop counter.
  19    integer :: i, j, k  ! Index dimension.
  20    integer :: num_levs ! Number of obs levels.
  21
  22    real    :: dx, dxm  ! Interpolation weights.
  23    real    :: dy, dym  ! Interpolation weights.
  24    real    :: v_h(kms:kme)      ! Model value h at ob hor. location.
  25    real    :: v_p(kms:kme)      ! Model value p at ob hor. location.
  26
  27    real, allocatable :: model_u(:,:)  ! Model value u at ob location.
  28    real, allocatable :: model_v(:,:)  ! Model value v at ob location.
  29    real, allocatable :: model_t(:,:)  ! Model value t at ob location.
  30    real, allocatable :: model_q(:,:)  ! Model value t at ob location.
  31    real    :: model_slp                  ! Model value slp at ob location.
  32
  33
  34    if (trace_use_dull) call da_trace_entry("da_get_innov_vector_bogus")
  35
  36    allocate (model_u(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
  37    allocate (model_v(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
  38    allocate (model_t(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
  39    allocate (model_q(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
  40    model_u(:,:) = 0.0
  41    model_v(:,:) = 0.0
  42    model_t(:,:) = 0.0
  43    model_q(:,:) = 0.0
  44
  45    if ( it > 1 ) then
  46       do n=iv%info(bogus)%n1,iv%info(bogus)%n2
  47          do k=1, iv%info(bogus)%levels(n)
  48             if (iv%bogus(n)%u(k)%qc == fails_error_max) iv%bogus(n)%u(k)%qc = 0
  49             if (iv%bogus(n)%v(k)%qc == fails_error_max) iv%bogus(n)%v(k)%qc = 0
  50             if (iv%bogus(n)%t(k)%qc == fails_error_max) iv%bogus(n)%t(k)%qc = 0
  51             if (iv%bogus(n)%q(k)%qc == fails_error_max) iv%bogus(n)%q(k)%qc = 0
  52          end do
  53       end do
  54    end if
  55
  56    do n=iv%info(bogus)%n1,iv%info(bogus)%n2
  57       num_levs = iv%info(bogus)%levels(n)
  58
  59       if (num_levs < 1) cycle
  60
  61       i   = iv%info(bogus)%i(1,n)
  62       j   = iv%info(bogus)%j(1,n)
  63       dx  = iv%info(bogus)%dx(1,n)
  64       dy  = iv%info(bogus)%dy(1,n)
  65       dxm = iv%info(bogus)%dxm(1,n)
  66       dym = iv%info(bogus)%dym(1,n)
  67
  68       do k=kts,kte
  69          v_h(k) = dym*(dxm*grid%xb%h(i,j,k) + dx*grid%xb%h(i+1,j,k)) + dy *(dxm*grid%xb%h(i,j+1,k) + dx*grid%xb%h(i+1,j+1,k))
  70          v_p(k) = dym*(dxm*grid%xb%p(i,j,k) + dx*grid%xb%p(i+1,j,k)) + dy *(dxm*grid%xb%p(i,j+1,k) + dx*grid%xb%p(i+1,j+1,k))
  71       end do
  72
  73       do k=1, iv%info(bogus)%levels(n)
  74          if (iv % bogus(n) % p(k) > 1.0) then
  75             call da_to_zk(iv % bogus(n) % p(k), v_p, v_interp_p, iv%info(bogus)%zk(k,n))
  76          else if (iv % bogus(n) % h(k) > 0.0) then
  77             call da_to_zk(iv % bogus(n) % h(k), v_h, v_interp_h, iv%info(bogus)%zk(k,n))
  78          end if
  79
  80          if (iv%info(bogus)%zk(k,n) < 0.0 .and.  .not.anal_type_verify) then
  81             iv % bogus(n) % u(k) % qc = missing_data
  82             iv % bogus(n) % v(k) % qc = missing_data
  83             iv % bogus(n) % t(k) % qc = missing_data
  84             iv % bogus(n) % q(k) % qc = missing_data
  85          end if
  86       end do
  87    end do
  88
  89    call da_convert_zk (iv%info(bogus))
  90
  91    ! [1.4] Interpolate horizontally:
  92
  93 #ifdef A2C
  94    call da_interp_lin_3d (grid%xb%u, iv%info(bogus), model_u,'u')
  95    call da_interp_lin_3d (grid%xb%v, iv%info(bogus), model_v,'v')
  96 #else
  97    call da_interp_lin_3d (grid%xb%u, iv%info(bogus), model_u)
  98    call da_interp_lin_3d (grid%xb%v, iv%info(bogus), model_v)
  99 #endif
 100    call da_interp_lin_3d (grid%xb%t, iv%info(bogus), model_t)
 101    call da_interp_lin_3d (grid%xb%q, iv%info(bogus), model_q)
 102
 103    do n=iv%info(bogus)%n1,iv%info(bogus)%n2
 104       num_levs = iv%info(bogus)%levels(n)
 105
 106       if (num_levs < 1) cycle
 107
 108       i   = iv%info(bogus)%i(1,n)
 109       j   = iv%info(bogus)%j(1,n)
 110       dx  = iv%info(bogus)%dx(1,n)
 111       dy  = iv%info(bogus)%dy(1,n)
 112       dxm = iv%info(bogus)%dxm(1,n)
 113       dym = iv%info(bogus)%dym(1,n)
 114
 115       model_slp = dym*(dxm*grid%xb%slp(i,j)   + dx*grid%xb%slp(i+1,j)) &
 116          + dy *(dxm*grid%xb%slp(i,j+1) + dx*grid%xb%slp(i+1,j+1))
 117
 118       !------------------------------------------------------------------------
 119       ! [2.0] Initialise components of innovation vector:
 120       !------------------------------------------------------------------------
 121
 122       iv % bogus(n) % slp % inv = 0.0
 123
 124       if (ABS(ob % bogus(n) % slp - missing_r) > 1.0 .AND. &
 125            iv % bogus(n) % slp % qc >= obs_qc_pointer) then
 126         iv % bogus(n) % slp % inv = ob % bogus(n) % slp - model_slp
 127       end if
 128
 129       do k = 1, iv%info(bogus)%levels(n)
 130          iv % bogus(n) % u(k) % inv = 0.0
 131          iv % bogus(n) % v(k) % inv = 0.0
 132          iv % bogus(n) % t(k) % inv = 0.0
 133          iv % bogus(n) % q(k) % inv = 0.0
 134
 135          !------------------------------------------------------------------------
 136          ! [4.0] Fast interpolation:
 137          !------------------------------------------------------------------------
 138
 139          if (ob % bogus(n) % u(k) > missing_r .AND. iv % bogus(n) % u(k) % qc >= obs_qc_pointer) then
 140            iv % bogus(n) % u(k) % inv = ob % bogus(n) % u(k) - model_u(k,n)
 141          end if
 142
 143          if (ob % bogus(n) % v(k) > missing_r .AND. iv % bogus(n) % v(k) % qc >= obs_qc_pointer) then
 144            iv % bogus(n) % v(k) % inv = ob % bogus(n) % v(k) - model_v(k,n)
 145          end if
 146
 147          if (ob % bogus(n) % t(k) > missing_r .AND. iv % bogus(n) % t(k) % qc >= obs_qc_pointer) then
 148             ! only for global Bogus(YRG 07/15/2005):
 149             if (iv%info(bogus)%platform(n)(8:12) /= 'TCBOG') then
 150                iv % bogus(n) % t(k) % inv = ob % bogus(n) % t(k) - model_t(k,n)
 151             else
 152                iv % bogus(n) % t(k) % inv = missing_r
 153                iv % bogus(n) % t(k) % qc  = missing_data
 154             end if
 155          end if
 156
 157          if (ob % bogus(n) % q(k) > missing_r .AND. iv % bogus(n) % q(k) % qc >= obs_qc_pointer) then
 158             ! only for global Bogus(YRG 07/15/2005):
 159             if (iv%info(bogus)%platform(n)(8:12) /= 'TCBOG') then
 160                iv % bogus(n) % q(k) % inv = ob % bogus(n) % q(k) - model_q(k,n)
 161             else
 162               iv % bogus(n) % q(k) % inv = missing_r
 163               iv % bogus(n) % q(k) % qc  = missing_data
 164             end if
 165          end if
 166       end do
 167    end do
 168
 169    !------------------------------------------------------------------------
 170    ! [5.0] Perform optional maximum error check:
 171    !------------------------------------------------------------------------
 172
 173    if ( check_max_iv ) &
 174       call da_check_max_iv_bogus(iv,ob, it, num_qcstat_conv)
 175
 176     deallocate (model_u)
 177     deallocate (model_v)
 178     deallocate (model_t)
 179     deallocate (model_q)
 180
 181    if (trace_use_dull) call da_trace_exit("da_get_innov_vector_bogus")
 182
 183 end subroutine da_get_innov_vector_bogus
 184
 185