1 subroutine da_get_innov_vector_pilot( it,num_qcstat_conv, grid, ob, iv)
3 !-----------------------------------------------------------------------
5 ! Updated for Analysis on Arakawa-C grid
6 ! Author: Syed RH Rizvi, MMM/ESSL/NCAR, Date: 10/22/2008
7 !-----------------------------------------------------------------------
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(:,:,:,:)
17 integer :: n ! Loop counter.
18 integer :: i, j, k ! Index dimension.
20 real :: dx, dxm ! Interpolation weights.
21 real :: dy, dym ! Interpolation weights.
23 real, allocatable :: model_u(:,:) ! Model value u at ob location.
24 real, allocatable :: model_v(:,:) ! Model value v at ob location.
26 real :: v_h(kms:kme) ! Model value h at ob hor. location.
27 real :: v_p(kms:kme) ! Model value p at ob hor. location.
29 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_pilot")
31 allocate (model_u(iv%info(pilot)%max_lev,iv%info(pilot)%n1:iv%info(pilot)%n2))
32 allocate (model_v(iv%info(pilot)%max_lev,iv%info(pilot)%n1:iv%info(pilot)%n2))
38 do n=iv%info(pilot)%n1,iv%info(pilot)%n2
39 do k=1, iv%info(pilot)%levels(n)
40 if (iv%pilot(n)%u(k)%qc == fails_error_max) iv%pilot(n)%u(k)%qc = 0
41 if (iv%pilot(n)%v(k)%qc == fails_error_max) iv%pilot(n)%v(k)%qc = 0
46 do n=iv%info(pilot)%n1,iv%info(pilot)%n2
47 ! [1.3] Get horizontal interpolation weights:
49 i = iv%info(pilot)%i(1,n)
50 j = iv%info(pilot)%j(1,n)
51 dx = iv%info(pilot)%dx(1,n)
52 dy = iv%info(pilot)%dy(1,n)
53 dxm = iv%info(pilot)%dxm(1,n)
54 dym = iv%info(pilot)%dym(1,n)
57 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))
58 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))
61 do k=1, iv%info(pilot)%levels(n)
62 if (iv % pilot(n) % p(k) > 1.0) then
63 call da_to_zk(iv % pilot(n) % p(k), v_p, v_interp_p, iv%info(pilot)%zk(k,n))
64 else if (iv%pilot(n) % h(k) > missing_r) then
65 call da_to_zk(iv % pilot(n) % h(k), v_h, v_interp_h, iv%info(pilot)%zk(k,n))
68 if (iv%info(pilot)%zk(k,n) < 0.0 .and. .not.anal_type_verify) then
69 iv % pilot(n) % u(k) % qc = missing_data
70 iv % pilot(n) % v(k) % qc = missing_data
75 call da_convert_zk (iv%info(pilot))
77 ! [1.4] Interpolate horizontally:
79 call da_interp_lin_3d (grid%xb%u, iv%info(pilot), model_u,'u')
80 call da_interp_lin_3d (grid%xb%v, iv%info(pilot), model_v,'v')
82 call da_interp_lin_3d (grid%xb%u, iv%info(pilot), model_u)
83 call da_interp_lin_3d (grid%xb%v, iv%info(pilot), model_v)
86 do n=iv%info(pilot)%n1,iv%info(pilot)%n2
87 !------------------------------------------------------------------------
88 ! [2.0] Initialise components of innovation vector:
89 !------------------------------------------------------------------------
91 do k = 1, iv%info(pilot)%levels(n)
92 iv % pilot(n) % u(k) % inv = 0.0
93 iv % pilot(n) % v(k) % inv = 0.0
95 !------------------------------------------------------------------------
96 ! [4.0] Fast interpolation:
97 !------------------------------------------------------------------------
99 if (ob % pilot(n) % u(k) > missing_r .AND. iv % pilot(n) % u(k) % qc >= obs_qc_pointer) then
100 iv % pilot(n) % u(k) % inv = ob % pilot(n) % u(k) - model_u(k,n)
103 if (ob % pilot(n) % v(k) > missing_r .AND. iv % pilot(n) % v(k) % qc >= obs_qc_pointer) then
104 iv % pilot(n) % v(k) % inv = ob % pilot(n) % v(k) - model_v(k,n)
109 !------------------------------------------------------------------------
110 ! [5.0] Perform optional maximum error check:
111 !------------------------------------------------------------------------
113 if ( check_max_iv ) &
114 call da_check_max_iv_pilot(iv, it, num_qcstat_conv)
119 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_pilot")
121 end subroutine da_get_innov_vector_pilot