da_get_innov_vector_airep.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_airep( it, xb, xp, ob, iv)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 ! Update :
6 ! 01/24/2007 Syed RH Rizvi
7 ! Updated for "VERIFY"
8 !-----------------------------------------------------------------------
9
10 implicit none
11
12 integer, intent(in) :: it ! External iteration.
13 type(xb_type), intent(in) :: xb ! first guess state.
14 type(xpose_type), intent(in) :: xp ! Domain decomposition vars.
15 type(y_type), intent(inout) :: ob ! Observation structure.
16 type(ob_type), intent(inout) :: iv ! O-B structure.
17
18 integer :: n ! Loop counter.
19 integer :: i, j, k ! Index dimension.
20 integer :: num_levs ! Number of obs levels.
21 real :: dx, dxm ! Interpolation weights.
22 real :: dy, dym ! Interpolation weights.
23 real, dimension(1:max_ob_levels) :: model_u ! Model value u at ob location.
24 real, dimension(1:max_ob_levels) :: model_v ! Model value v at ob location.
25 real, dimension(1:max_ob_levels) :: model_t ! Model value t at ob location.
26
27
28 real, dimension(xp%kms:xp%kme) :: v_h ! Model value h at ob hor. location.
29 real, dimension(xp%kms:xp%kme) :: v_p ! Model value p at ob hor. location.
30
31 integer :: itu,ituf,itvv,itvvf,itt,ittf
32
33 if (iv % num_airep > 0) then
34
35 itu = 0; itvv = 0; itt = 0;
36 ituf = 0; itvvf = 0; ittf = 0;
37
38 do n=iv%ob_numb(iv%current_ob_time-1)%airep + 1, iv%ob_numb(iv%current_ob_time)%airep
39 num_levs = iv % airep(n) % info % levels
40 if (num_levs < 1) cycle
41
42 model_u(:) = 0.0
43 model_v(:) = 0.0
44 model_t(:) = 0.0
45
46 ! [1.1] Get horizontal interpolation weights:
47
48 i = iv%airep(n)%loc%i
49 j = iv%airep(n)%loc%j
50 dx = iv%airep(n)%loc%dx
51 dy = iv%airep(n)%loc%dy
52 dxm = iv%airep(n)%loc%dxm
53 dym = iv%airep(n)%loc%dym
54
55 do k=xp%kts,xp%kte
56 v_h(k) = dym*(dxm*xb%h(i,j ,k) + dx*xb%h(i+1,j ,k)) &
57 + dy *(dxm*xb%h(i,j+1,k) + dx*xb%h(i+1,j+1,k))
58 v_p(k) = dym*(dxm*xb%p(i,j ,k) + dx*xb%p(i+1,j ,k)) &
59 + dy *(dxm*xb%p(i,j+1,k) + dx*xb%p(i+1,j+1,k))
60 end do
61
62 num_levs=0
63 do k=1, iv % airep(n) % info % levels
64 iv%airep(n)%zk(k)=missing_r
65 if (iv % airep(n) % p(k) > 1.0) then
66 call da_to_zk(iv % airep(n) % p(k), v_p, xp, v_interp_p, &
67 iv%airep(n)%zk(k))
68 else if (iv % airep(n) % h(k) > 0.0) then
69 call da_to_zk(iv % airep(n) % h(k), v_h, xp, v_interp_h, &
70 iv%airep(n)%zk(k))
71 end if
72
73 if (iv%airep(n)%zk(k) < 0.0 .and. .not.anal_type_verify ) then
74 iv % airep(n) % u(k) % qc = missing
75 iv % airep(n) % v(k) % qc = missing
76 iv % airep(n) % t(k) % qc = missing
77 end if
78 end do
79
80 num_levs = iv % airep(n) % info % levels
81
82 ! [1.2] Interpolate horizontally:
83 call da_interp_lin_3d( xb % u, xp, i, j, dx, dy, dxm, dym, &
84 model_u, max_ob_levels, iv%airep(n)%zk, num_levs)
85
86 call da_interp_lin_3d( xb % v, xp, i, j, dx, dy, dxm, dym, &
87 model_v, max_ob_levels, iv%airep(n)%zk, num_levs)
88
89 call da_interp_lin_3d( xb % t, xp, i, j, dx, dy, dxm, dym, &
90 model_t, max_ob_levels, iv%airep(n)%zk, num_levs)
91
92 !-------------------------------------------------------------------
93 ! [2.0] Initialise components of innovation vector:
94 !-------------------------------------------------------------------
95
96 do k = 1, iv % airep(n) % info % levels
97 iv % airep(n) % u(k) % inv = 0.0
98 iv % airep(n) % v(k) % inv = 0.0
99 iv % airep(n) % t(k) % inv = 0.0
100
101 !----------------------------------------------------------------
102 ! [3.0] Fast interpolation:
103 !----------------------------------------------------------------
104
105 if (ob % airep(n) % u(k) > missing_r .AND. &
106 iv % airep(n) % u(k) % qc >= obs_qc_pointer) then
107 iv % airep(n) % u(k) % inv = ob % airep(n) % u(k) - model_u(k)
108 end if
109
110 if (ob % airep(n) % v(k) > missing_r .AND. &
111 iv % airep(n) % v(k) % qc >= obs_qc_pointer) then
112 iv % airep(n) % v(k) % inv = ob % airep(n) % v(k) - model_v(k)
113 end if
114
115 if (ob % airep(n) % t(k) > missing_r .AND. &
116 iv % airep(n) % t(k) % qc >= obs_qc_pointer) then
117 iv % airep(n) % t(k) % inv = ob % airep(n) % t(k) - model_t(k)
118 end if
119 end do
120
121 !-------------------------------------------------------------------
122 ! [5.0] Perform optional maximum error check:
123 !-------------------------------------------------------------------
124
125 if (check_max_iv)then
126 call da_check_max_iv_airep(it, iv % airep(n),itu,ituf,itvv, &
127 itvvf,itt,ittf)
128 end if
129 end do
130
131 if (rootproc .and. check_max_iv_print) then
132 write(unit = check_max_iv_unit, fmt ='(A,i5,A)') &
133 'For outer iteration ',it, &
134 ', Total Rejections for Airep follows:'
135 write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
136 'Number of failed u-wind observations: ',ituf, ' on ',itu, &
137 'Number of failed v-wind observations: ',itvvf,' on ',itvv, &
138 'Number of failed temperature observations:',ittf, ' on ',itt, &
139 'Finally Total Airep rejections ',ituf+itvvf+ittf,' on ', &
140 itu +itvv +itt
141 end if
142 end if
143
144 end subroutine da_get_innov_vector_airep
145
146