1 subroutine da_get_innov_vector_qscat (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.
19 real :: dx, dxm ! Interpolation weights.
20 real :: dy, dym ! Interpolation weights.
22 real, allocatable :: model_u(:,:) ! Model value u at ob location.
23 real, allocatable :: model_v(:,:) ! Model value v at ob location.
25 real :: v_h(kms:kme) ! Model value h at ob hor. location.
27 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_qscat")
29 allocate (model_u(iv%info(qscat)%max_lev,iv%info(qscat)%n1:iv%info(qscat)%n2))
30 allocate (model_v(iv%info(qscat)%max_lev,iv%info(qscat)%n1:iv%info(qscat)%n2))
33 do n=iv%info(qscat)%n1,iv%info(qscat)%n2
34 if (iv%qscat(n)%u%qc == fails_error_max) iv%qscat(n)%u%qc = 0
35 if (iv%qscat(n)%v%qc == fails_error_max) iv%qscat(n)%v%qc = 0
39 do n=iv%info(qscat)%n1,iv%info(qscat)%n2
41 ! [1.1] Get horizontal interpolation weights:
43 i = iv%info(qscat)%i(1,n)
44 j = iv%info(qscat)%j(1,n)
45 dx = iv%info(qscat)%dx(1,n)
46 dy = iv%info(qscat)%dy(1,n)
47 dxm = iv%info(qscat)%dxm(1,n)
48 dym = iv%info(qscat)%dym(1,n)
51 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))
54 if (iv % qscat(n) % h > missing_r) then
55 call da_to_zk(iv % qscat(n) % h, v_h, v_interp_h, iv%info(qscat)%zk(1,n))
56 if (iv%info(qscat)%zk(1,n) < 1.0) then
57 iv%info(qscat)%zk(1,n) = 1.0
62 call da_convert_zk (iv%info(qscat))
64 if (.not. anal_type_verify) then
65 do n=iv%info(qscat)%n1,iv%info(qscat)%n2
66 if (iv%info(qscat)%zk(1,n) < 0.0) then
67 iv%qscat(n)%u%qc = missing_data
68 iv%qscat(n)%v%qc = missing_data
74 call da_interp_lin_3d (grid%xb%u, iv%info(qscat), model_u,'u')
75 call da_interp_lin_3d (grid%xb%v, iv%info(qscat), model_v,'v')
77 call da_interp_lin_3d (grid%xb%u, iv%info(qscat), model_u)
78 call da_interp_lin_3d (grid%xb%v, iv%info(qscat), model_v)
81 do n=iv%info(qscat)%n1,iv%info(qscat)%n2
83 !------------------------------------------------------------------------
84 ! [2.0] Initialise components of innovation vector:
85 !------------------------------------------------------------------------
87 !------------------------------------------------------------------------
88 ! [3.0] Fast interpolation:
89 !------------------------------------------------------------------------
91 if (ob % qscat(n) % u > missing_r .AND. &
92 iv % qscat(n) % u % qc >= obs_qc_pointer) then
93 iv % qscat(n) % u % inv = ob % qscat(n) % u - model_u(1,n)
96 if (ob % qscat(n) % v > missing_r .AND. &
97 iv % qscat(n) % v % qc >= obs_qc_pointer) then
98 iv % qscat(n) % v % inv = ob % qscat(n) % v - model_v(1,n)
102 !------------------------------------------------------------------------
103 ! [5.0] Perform optional maximum error check:
104 !------------------------------------------------------------------------
106 if ( check_max_iv ) &
107 call da_check_max_iv_qscat(iv, it, num_qcstat_conv)
112 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_qscat")
114 end subroutine da_get_innov_vector_qscat