da_get_innov_vector_synop.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_synop( 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 ! model 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 :: v_h(kms:kme) ! Model value h at ob hor. loc
20 real :: v_p(kms:kme) ! Model value p at ob hor. loc
21
22 real :: hd, psfcm
23
24 real :: ho, to, qo
25
26 integer :: itu,ituf,itvv,itvvf,itp,itpf,itt,ittf,itqv,itqvf
27
28 real, allocatable :: model_u(:,:)
29 real, allocatable :: model_v(:,:)
30 real, allocatable :: model_t(:,:)
31 real, allocatable :: model_q(:,:)
32 real, allocatable :: model_p(:,:)
33 real, allocatable :: model_hsm(:,:)
34
35 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_synop")
36
37 itu = 0; itvv = 0; itp = 0; itt = 0; itqv = 0;
38 ituf = 0; itvvf = 0; itpf = 0; ittf = 0; itqvf = 0;
39
40 allocate (model_u (1,iv%info(synop)%n1:iv%info(synop)%n2))
41 allocate (model_v (1,iv%info(synop)%n1:iv%info(synop)%n2))
42 allocate (model_t (1,iv%info(synop)%n1:iv%info(synop)%n2))
43 allocate (model_q (1,iv%info(synop)%n1:iv%info(synop)%n2))
44 allocate (model_p (1,iv%info(synop)%n1:iv%info(synop)%n2))
45 allocate (model_hsm(1,iv%info(synop)%n1:iv%info(synop)%n2))
46
47 if (sfc_assi_options == sfc_assi_options_1) then
48
49 do n=iv%info(synop)%n1,iv%info(synop)%n2
50 ! [1.1] Get horizontal interpolation weights:
51 i = iv%info(synop)%i(1,n)
52 j = iv%info(synop)%j(1,n)
53 dx = iv%info(synop)%dx(1,n)
54 dy = iv%info(synop)%dy(1,n)
55 dxm = iv%info(synop)%dxm(1,n)
56 dym = iv%info(synop)%dym(1,n)
57
58 ! Surface correction
59
60 iv%synop(n)%p%inv = ob%synop(n)%p
61 iv%synop(n)%t%inv = ob%synop(n)%t
62 iv%synop(n)%q%inv = ob%synop(n)%q
63 iv%synop(n)%u%inv = ob%synop(n)%u
64 iv%synop(n)%v%inv = ob%synop(n)%v
65
66 if (iv % synop(n) % h > missing_r) then
67 do k=kts,kte
68 v_h(k) = dym*(dxm*grid%xb%h(i,j ,k) + dx*grid%xb%h(i+1,j ,k)) &
69 + dy *(dxm*grid%xb%h(i,j+1,k) + dx*grid%xb%h(i+1,j+1,k))
70 end do
71
72 hd = v_h(kts) - iv % synop(n) % h
73 if (abs(hd) <= Max_StHeight_Diff .or. anal_type_verify ) then
74 if (iv % synop(n) % h < v_h(kts)) then
75 iv%info(synop)%zk(:,n) = 1.0+1.0e-6
76 call da_obs_sfc_correction(iv%info(synop), iv%synop(n), n, grid%xb)
77
78 ! To keep the original "ob" with no change for multiple
79 ! outer-loops use:
80 ! ob%synop(n)%p = iv%synop(n)%p%inv
81 ! ob%synop(n)%t = iv%synop(n)%t%inv
82 ! ob%synop(n)%q = iv%synop(n)%q%inv
83 ! ob%synop(n)%u = iv%synop(n)%u%inv
84 ! ob%synop(n)%v = iv%synop(n)%v%inv
85 else
86 call da_to_zk(iv % synop(n) % h, v_h, v_interp_h, iv%info(synop)%zk(1,n))
87 end if
88 end if
89 else if (ob % synop(n) % p > 1.0) then
90 do k=kts,kte
91 v_p(k) = dym*(dxm*grid%xb%p(i,j ,k) + dx*grid%xb%p(i+1,j ,k)) &
92 + dy *(dxm*grid%xb%p(i,j+1,k) + dx*grid%xb%p(i+1,j+1,k))
93 end do
94
95 call da_to_zk(ob % synop(n) % p, v_p, v_interp_p, iv%info(synop)%zk(1,n))
96 if (iv%info(synop)%zk(1,n) < 0.0 .and. .not.anal_type_verify ) then
97 iv % synop(n) % p % inv = missing_r
98 iv % synop(n) % p % qc = missing
99 iv%info(synop)%zk(:,n) = 1.0+1.0e-6
100 end if
101 end if
102 end do
103
104 call da_convert_zk (iv%info(synop))
105
106 if (.not.anal_type_verify ) then
107 do n=iv%info(synop)%n1,iv%info(synop)%n2
108 if (iv%info(synop)%zk(1,n) < 0.0) then
109 iv % synop(n) % u % qc = missing
110 iv % synop(n) % v % qc = missing
111 iv % synop(n) % t % qc = missing
112 iv % synop(n) % q % qc = missing
113 iv % synop(n) % p % qc = missing
114 end if
115 end do
116 end if
117
118 ! [1.2] Interpolate horizontally:
119 call da_interp_lin_3d (grid%xb%u, iv%info(synop), model_u)
120 call da_interp_lin_3d (grid%xb%v, iv%info(synop), model_v)
121 call da_interp_lin_3d (grid%xb%t, iv%info(synop), model_t)
122 call da_interp_lin_3d (grid%xb%q, iv%info(synop), model_q)
123 call da_interp_lin_3d (grid%xb%p, iv%info(synop), model_p)
124 else if (sfc_assi_options == sfc_assi_options_2) then
125 ! Surface data assmiilation approach 2
126 !------------------------------------
127
128 ! 1.2.1 Surface assmiilation approach 2(10-m u, v, 2-m t, q, and sfc_p)
129
130 call da_interp_lin_2d (grid%xb%u10, iv%info(synop), 1, model_u(1,:))
131 call da_interp_lin_2d (grid%xb%v10, iv%info(synop), 1, model_v(1,:))
132 call da_interp_lin_2d (grid%xb%t2, iv%info(synop), 1, model_t(1,:))
133 call da_interp_lin_2d (grid%xb%q2, iv%info(synop), 1, model_q(1,:))
134 call da_interp_lin_2d (grid%xb%psfc, iv%info(synop), 1, model_p(1,:))
135 call da_interp_lin_2d (grid%xb%terr, iv%info(synop), 1, model_hsm(1,:))
136
137 do n=iv%info(synop)%n1,iv%info(synop)%n2
138 if (iv%synop(n)%p%qc >= 0) then
139 ho = iv%synop(n)%h
140 to = -888888.0
141 qo = -888888.0
142
143 if (iv%synop(n)%t%qc >= 0 .and. iv%synop(n)%q%qc >= 0) then
144 to = ob%synop(n)%t
145 qo = ob%synop(n)%q
146 call da_sfc_pre (psfcm, model_p(1,n), model_t(1,n), model_q(1,n), model_hsm(1,n), ho, to, qo)
147 else if (iv%synop(n)%t%qc >= 0 .and. iv%synop(n)%q%qc < 0) then
148 to = ob%synop(n)%t
149 call da_sfc_pre (psfcm, model_p(1,n), model_t(1,n), model_q(1,n), model_hsm(1,n), ho, to)
150 else
151 call da_sfc_pre (psfcm, model_p(1,n), model_t(1,n), model_q(1,n), model_hsm(1,n), ho)
152 end if
153
154 ! Pressure at the observed height:
155 model_p(1,n) = psfcm
156 end if
157 end do
158 end if
159
160 do n=iv%info(synop)%n1,iv%info(synop)%n2
161
162 !--------------------------------------------------------------------
163 ! [3.0] Fast interpolation:
164 !--------------------------------------------------------------------
165
166 if (ob % synop(n) % u > missing_r .AND. iv % synop(n) % u % qc >= obs_qc_pointer) then
167 iv % synop(n) % u % inv = iv%synop(n)%u%inv - model_u(1,n)
168 else
169 iv % synop(n) % u % inv = 0.0
170 end if
171
172 if (ob % synop(n) % v > missing_r .AND. iv % synop(n) % v % qc >= obs_qc_pointer) then
173 iv % synop(n) % v % inv = iv%synop(n)%v%inv - model_v(1,n)
174 else
175 iv % synop(n) % v % inv = 0.0
176 end if
177
178 if (ob % synop(n) % p > 0.0 .AND. iv % synop(n) % p % qc >= obs_qc_pointer) then
179 iv % synop(n) % p % inv = iv%synop(n)%p%inv - model_p(1,n)
180 else
181 iv % synop(n) % p % inv = 0.0
182 end if
183
184 if (ob % synop(n) % t > 0.0 .AND. iv % synop(n) % t % qc >= obs_qc_pointer) then
185 iv % synop(n) % t % inv = iv%synop(n)%t%inv - model_t(1,n)
186 else
187 iv % synop(n) % t % inv = 0.0
188 end if
189
190 if (ob % synop(n) % q > 0.0 .AND. iv % synop(n) % q % qc >= obs_qc_pointer) then
191 iv % synop(n) % q % inv = iv%synop(n)%q%inv - model_q(1,n)
192 else
193 iv % synop(n) % q % inv = 0.0
194 end if
195 end do
196
197 !--------------------------------------------------------------------
198 ! [5.0] Perform optional maximum error check:
199 !--------------------------------------------------------------------
200
201 if (check_max_iv) then
202 call da_check_max_iv_synop(iv, it, itu,ituf,itvv,itvvf,itp,itpf,itt,ittf,itqv,itqvf)
203 end if
204
205 if (rootproc .and. check_max_iv_print) then
206 write(unit= check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it, &
207 ', Total Rejections for Synop follows:'
208
209 write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
210 'Number of failed u-wind observations: ',ituf, ' on ',itu, &
211 'Number of failed v-wind observations: ',itvvf,' on ',itvv, &
212 'Number of failed pressure observations: ',itpf, ' on ',itp, &
213 'Number of failed temperature observations:',ittf, ' on ',itt, &
214 'Number of failed mixing ratio observations:',itqvf,' on ',itqv, &
215 'Finally Total Synop rejections ',ituf+itvvf+itpf+ittf+itqvf,' on ',&
216 itu +itvv +itp +itt +itqv
217 end if
218
219 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_synop")
220
221 end subroutine da_get_innov_vector_synop
222
223