var/da/da_qscat/da_get_innov_vector_qscat.inc

   1 subroutine da_get_innov_vector_qscat (it,num_qcstat_conv, grid, ob, iv)
   2
   3    !-----------------------------------------------------------------------
   4    ! Purpose: TBD
   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(inout) :: ob      ! Observation structure.
  14    type(iv_type),    intent(inout) :: iv      ! O-B structure.
  15    integer,          intent(inout) :: num_qcstat_conv(:,:,:,:)
  16
  17    integer :: n        ! Loop counter.
  18    integer :: i, j, k  ! Index dimension.
  19    real    :: dx, dxm  ! Interpolation weights.
  20    real    :: dy, dym  ! Interpolation weights.
  21    real    :: speed, direction
  22
  23    real, allocatable :: model_u(:,:)  ! Model value u at ob location.
  24    real, allocatable :: model_v(:,:)  ! Model value v at ob location.
  25
  26    real    :: v_h(kms:kme)      ! Model value h at ob hor. location.
  27
  28    if (trace_use_dull) call da_trace_entry("da_get_innov_vector_qscat")
  29
  30    allocate (model_u(iv%info(qscat)%max_lev,iv%info(qscat)%n1:iv%info(qscat)%n2))
  31    allocate (model_v(iv%info(qscat)%max_lev,iv%info(qscat)%n1:iv%info(qscat)%n2))
  32
  33    if ( it > 1 ) then
  34       do n=iv%info(qscat)%n1,iv%info(qscat)%n2
  35          if (iv%qscat(n)%u%qc == fails_error_max) iv%qscat(n)%u%qc = 0
  36          if (iv%qscat(n)%v%qc == fails_error_max) iv%qscat(n)%v%qc = 0
  37       end do
  38    end if
  39
  40    do n=iv%info(qscat)%n1,iv%info(qscat)%n2
  41
  42       ! [1.1] Get horizontal interpolation weights:
  43
  44       i   = iv%info(qscat)%i(1,n)
  45       j   = iv%info(qscat)%j(1,n)
  46       dx  = iv%info(qscat)%dx(1,n)
  47       dy  = iv%info(qscat)%dy(1,n)
  48       dxm = iv%info(qscat)%dxm(1,n)
  49       dym = iv%info(qscat)%dym(1,n)
  50
  51       do k=kts,kte
  52          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))
  53       end do
  54
  55       if (iv % qscat(n) % h > missing_r) then
  56          call da_to_zk(iv % qscat(n) % h, v_h, v_interp_h, iv%info(qscat)%zk(1,n))
  57          if (iv%info(qscat)%zk(1,n) < 1.0) then
  58             iv%info(qscat)%zk(1,n) = 1.0
  59          end if
  60       end if
  61    end do
  62
  63    call da_convert_zk (iv%info(qscat))
  64
  65    if (.not. anal_type_verify) then
  66       do n=iv%info(qscat)%n1,iv%info(qscat)%n2
  67          if (iv%info(qscat)%zk(1,n) < 0.0) then
  68             iv%qscat(n)%u%qc = missing_data
  69             iv%qscat(n)%v%qc = missing_data
  70          end if
  71       end do
  72    end if
  73
  74 #ifdef A2C
  75    call da_interp_lin_3d (grid%xb%u, iv%info(qscat), model_u,'u')
  76    call da_interp_lin_3d (grid%xb%v, iv%info(qscat), model_v,'v')
  77 #else
  78    call da_interp_lin_3d (grid%xb%u, iv%info(qscat), model_u)
  79    call da_interp_lin_3d (grid%xb%v, iv%info(qscat), model_v)
  80 #endif
  81
  82    do n=iv%info(qscat)%n1,iv%info(qscat)%n2
  83
  84       !------------------------------------------------------------------------
  85       ! [2.0] Initialise components of innovation vector:
  86       !------------------------------------------------------------------------
  87
  88       !------------------------------------------------------------------------
  89       ! [3.0] Fast interpolation:
  90       !------------------------------------------------------------------------
  91           if (wind_sd_qscat) then
  92               call da_ffdduv_model (speed,direction,model_u(1,n), model_v(1,n), convert_uv2fd)
  93
  94               if (ob%qscat(n)%u > missing_r .AND. iv%qscat(n)%u%qc >= obs_qc_pointer) then
  95                   iv%qscat(n)%u%inv = ob%qscat(n)%u - speed
  96               end if
  97
  98               if (ob%qscat(n)%v > missing_r .AND. iv%qscat(n)%v%qc >= obs_qc_pointer) then
  99                   iv%qscat(n)%v%inv = ob%qscat(n)%v - direction
 100                   if (iv%qscat(n)%v%inv > 180.0 ) iv%qscat(n)%v%inv = iv%qscat(n)%v%inv - 360.0
 101                   if (iv%qscat(n)%v%inv < -180.0 ) iv%qscat(n)%v%inv = iv%qscat(n)%v%inv + 360.0
 102               end if
 103           else
 104                if (ob % qscat(n) % u > missing_r .AND. &
 105                    iv % qscat(n) % u % qc >= obs_qc_pointer) then
 106                    iv % qscat(n) % u % inv = ob % qscat(n) % u - model_u(1,n)
 107                end if
 108
 109                if (ob % qscat(n) % v > missing_r .AND. &
 110                    iv % qscat(n) % v % qc >= obs_qc_pointer) then
 111                    iv % qscat(n) % v % inv = ob % qscat(n) % v - model_v(1,n)
 112                end if
 113           end if
 114    end do
 115
 116    !------------------------------------------------------------------------
 117    ! [5.0] Perform optional maximum error check:
 118    !------------------------------------------------------------------------
 119
 120    if ( check_max_iv ) &
 121       call da_check_max_iv_qscat(iv, it, num_qcstat_conv)
 122
 123    deallocate (model_u)
 124    deallocate (model_v)
 125
 126    if (trace_use_dull) call da_trace_exit("da_get_innov_vector_qscat")
 127
 128 end subroutine da_get_innov_vector_qscat
 129
 130