da_get_innov_vector_sound.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_sound( 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
24 real, dimension(1:max_ob_levels) :: model_u ! Model value u at ob location.
25 real, dimension(1:max_ob_levels) :: model_v ! Model value v at ob location.
26 real, dimension(1:max_ob_levels) :: model_t ! Model value t at ob location.
27 real, dimension(1:max_ob_levels) :: model_q ! Model value q at ob location.
28 ! real, dimension(1:max_ob_levels) :: model_h ! Model value h at ob location.
29
30 real, dimension(xp%kms:xp%kme) :: v_h ! Model value h at ob hor. location.
31 real, dimension(xp%kms:xp%kme) :: v_p ! Model value p at ob hor. location.
32
33 integer :: itu,ituf,itvv,itvvf,itt,ittf,itqv,itqvf
34
35
36
37 if (iv%num_sound < 1) return
38
39 if (trace_use) call da_trace_entry("da_get_innov_vector_sound")
40
41 itu = 0; itvv = 0; itt = 0; itqv = 0;
42 ituf = 0; itvvf = 0; ittf = 0; itqvf = 0;
43
44 do n=iv%ob_numb(iv%current_ob_time-1)%sound + 1, &
45 iv%ob_numb(iv%current_ob_time)%sound
46 num_levs = iv%sound(n)%info%levels
47
48 if (num_levs < 1) cycle
49
50 model_u(:) = 0.0
51 model_v(:) = 0.0
52 model_t(:) = 0.0
53 model_q(:) = 0.0
54 ! model_h(:) = 0.0
55
56 ! [1.1] Get horizontal interpolation weights:
57
58 i = iv%sound(n)%loc%i
59 j = iv%sound(n)%loc%j
60 dx = iv%sound(n)%loc%dx
61 dy = iv%sound(n)%loc%dy
62 dxm = iv%sound(n)%loc%dxm
63 dym = iv%sound(n)%loc%dym
64
65 do k=xp%kts,xp%kte
66 v_h(k) = dym*(dxm*xb%h(i,j ,k) + dx*xb%h(i+1,j ,k)) &
67 + dy *(dxm*xb%h(i,j+1,k) + dx*xb%h(i+1,j+1,k))
68 v_p(k) = dym*(dxm*xb%p(i,j ,k) + dx*xb%p(i+1,j ,k)) &
69 + dy *(dxm*xb%p(i,j+1,k) + dx*xb%p(i+1,j+1,k))
70 end do
71
72 do k=1, num_levs
73 iv%sound(n)%zk(k)=missing_r
74
75 if (iv%sound(n)%p(k) > 1.0) then
76 call da_to_zk(iv%sound(n)%p(k), v_p, xp, v_interp_p, &
77 iv%sound(n)%zk(k))
78 else if (iv%sound(n)%h(k) > 0.0) then
79 call da_to_zk(iv%sound(n)%h(k), v_h, xp, v_interp_h, &
80 iv%sound(n)%zk(k))
81 end if
82
83 if (iv%sound(n)%zk(k) < 0.0 .and. .not.anal_type_verify) then
84 iv%sound(n)%u(k)%qc = missing
85 iv%sound(n)%v(k)%qc = missing
86 iv%sound(n)%t(k)%qc = missing
87 iv%sound(n)%q(k)%qc = missing
88 end if
89 end do
90
91 ! [1.2] Interpolate horizontally to ob:
92
93 call da_interp_lin_3d( xb%u, xp, i, j, dx, dy, dxm, dym, &
94 model_u, max_ob_levels, iv%sound(n)%zk, num_levs)
95 call da_interp_lin_3d( xb%v, xp, i, j, dx, dy, dxm, dym, &
96 model_v, max_ob_levels, iv%sound(n)%zk, num_levs)
97 call da_interp_lin_3d( xb%t, xp, i, j, dx, dy, dxm, dym, &
98 model_t, max_ob_levels, iv%sound(n)%zk, num_levs)
99 call da_interp_lin_3d( xb%q, xp, i, j, dx, dy, dxm, dym, &
100 model_q, max_ob_levels, iv%sound(n)%zk, num_levs)
101
102 ! call da_interp_lin_3d( xb%h, xp, i, j, dx, dy, dxm, dym, &
103 ! model_h, max_ob_levels, iv%sound(n)%zk, num_levs)
104
105 !----------------------------------------------------------------------
106 ! [2.0] Initialise components of innovation vector:
107 !----------------------------------------------------------------------
108
109 do k = 1, iv%sound(n)%info%levels
110 iv%sound(n)%u(k)%inv = 0.0
111 iv%sound(n)%v(k)%inv = 0.0
112 iv%sound(n)%t(k)%inv = 0.0
113 iv%sound(n)%q(k)%inv = 0.0
114
115 !-------------------------------------------------------------------
116 ! [3.0] Interpolation:
117 !-------------------------------------------------------------------
118
119 if (ob%sound(n)%u(k) > missing_r .AND. &
120 iv%sound(n)%u(k)%qc >= obs_qc_pointer) then
121 iv%sound(n)%u(k)%inv = ob%sound(n)%u(k) - model_u(k)
122 end if
123
124 if (ob%sound(n)%v(k) > missing_r .AND. &
125 iv%sound(n)%v(k)%qc >= obs_qc_pointer) then
126 iv%sound(n)%v(k)%inv = ob%sound(n)%v(k) - model_v(k)
127 end if
128
129
130 if (ob%sound(n)%t(k) > missing_r .AND. &
131 iv%sound(n)%t(k)%qc >= obs_qc_pointer) then
132 iv%sound(n)%t(k)%inv = ob%sound(n)%t(k) - model_t(k)
133 end if
134
135 if (ob%sound(n)%q(k) > missing_r .AND. &
136 iv%sound(n)%q(k)%qc >= obs_qc_pointer) then
137 iv%sound(n)%q(k)%inv = ob%sound(n)%q(k) - model_q(k)
138 end if
139 end do
140
141 !----------------------------------------------------------------------
142 ! [5.0] Perform optional maximum error check:
143 !----------------------------------------------------------------------
144
145 if (check_max_iv) then
146 call da_check_max_iv_sound(it, iv % sound(n), &
147 itu,ituf,itvv,itvvf,itt,ittf,itqv,itqvf)
148 end if
149 end do
150
151 if (rootproc .and. check_max_iv_print) then
152 write(unit = check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it, &
153 ', Total Rejections for Sound follows:'
154
155 write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
156 'Number of failed u-wind observations: ',ituf, ' on ',itu, &
157 'Number of failed v-wind observations: ',itvvf,' on ',itvv, &
158 'Number of failed temperature observations:',ittf, ' on ',itt, &
159 'Number of failed mixing ratio observations:',itqvf,' on ',itqv, &
160 'Finally Total Sound rejections: ',ituf+itvvf+ittf+itqvf,' on ',&
161 itu +itvv+itt +itqv
162 end if
163
164 if (trace_use) call da_trace_exit("da_get_innov_vector_sound")
165
166 end subroutine da_get_innov_vector_sound
167
168