da_get_innov_vector_qscat.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_qscat (it, grid, ob, iv)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 integer, intent(in) :: it ! External iteration.
10 type(domain), intent(in) :: grid ! first guess state.
11 type(y_type), intent(inout) :: ob ! Observation structure.
12 type(iv_type), intent(inout) :: iv ! O-B structure.
13
14 integer :: n ! Loop counter.
15 integer :: i, j, k ! Index dimension.
16 real :: dx, dxm ! Interpolation weights.
17 real :: dy, dym ! Interpolation weights.
18
19 real, allocatable :: model_u(:,:) ! Model value u at ob location.
20 real, allocatable :: model_v(:,:) ! Model value v at ob location.
21
22 real :: v_h(kms:kme) ! Model value h at ob hor. location.
23
24 integer :: itu, itvv, ituf, itvvf
25
26 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_qscat")
27
28 itu = 0 ; itvv = 0 ; ituf = 0 ; itvvf = 0
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 do n=iv%info(qscat)%n1,iv%info(qscat)%n2
34
35 ! [1.1] Get horizontal interpolation weights:
36
37 i = iv%info(qscat)%i(1,n)
38 j = iv%info(qscat)%j(1,n)
39 dx = iv%info(qscat)%dx(1,n)
40 dy = iv%info(qscat)%dy(1,n)
41 dxm = iv%info(qscat)%dxm(1,n)
42 dym = iv%info(qscat)%dym(1,n)
43
44 do k=kts,kte
45 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))
46 end do
47
48 if (iv % qscat(n) % h > missing_r) then
49 call da_to_zk(iv % qscat(n) % h, v_h, v_interp_h, iv%info(qscat)%zk(1,n))
50 if (iv%info(qscat)%zk(1,n) < 1.0) then
51 iv%info(qscat)%zk(1,n) = 1.0
52 end if
53 end if
54 end do
55
56 call da_convert_zk (iv%info(qscat))
57
58 if (.not. anal_type_verify) then
59 do n=iv%info(qscat)%n1,iv%info(qscat)%n2
60 if (iv%info(qscat)%zk(1,n) < 0.0) then
61 iv%qscat(n)%u%qc = missing
62 iv%qscat(n)%v%qc = missing
63 end if
64 end do
65 end if
66
67 call da_interp_lin_3d (grid%xb%u, iv%info(qscat), model_u)
68 call da_interp_lin_3d (grid%xb%v, iv%info(qscat), model_v)
69
70 do n=iv%info(qscat)%n1,iv%info(qscat)%n2
71
72 !------------------------------------------------------------------------
73 ! [2.0] Initialise components of innovation vector:
74 !------------------------------------------------------------------------
75
76 !------------------------------------------------------------------------
77 ! [3.0] Fast interpolation:
78 !------------------------------------------------------------------------
79
80 if (ob % qscat(n) % u > missing_r .AND. &
81 iv % qscat(n) % u % qc >= obs_qc_pointer) then
82 iv % qscat(n) % u % inv = ob % qscat(n) % u - model_u(1,n)
83 end if
84
85 if (ob % qscat(n) % v > missing_r .AND. &
86 iv % qscat(n) % v % qc >= obs_qc_pointer) then
87 iv % qscat(n) % v % inv = ob % qscat(n) % v - model_v(1,n)
88 end if
89 end do
90
91 !------------------------------------------------------------------------
92 ! [5.0] Perform optional maximum error check:
93 !------------------------------------------------------------------------
94
95 if (check_max_iv) call da_check_max_iv_qscat(iv, it, itu, itvv, ituf,itvvf)
96
97 if (rootproc .and. check_max_iv_print) then
98 write(unit = check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it,&
99 ', Total Rejections for qscat follows:'
100 write(unit = check_max_iv_unit, fmt = '(/,4(2(A,I6),/))') &
101 'Number of failed u-wind observations: ',ituf, ' on ',itu, &
102 'Number of failed v-wind observations: ',itvvf,' on ',itvv
103 end if
104
105 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_qscat")
106
107 end subroutine da_get_innov_vector_qscat
108
109