1 subroutine da_get_innov_vector_sonde_sfc( 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(:,:,:,:)
18 integer :: n ! Loop counter.
19 integer :: i, j, k ! Index dimension.
20 real :: dx, dxm ! Interpolation weights.
21 real :: dy, dym ! Interpolation weights.
22 real, allocatable :: model_u(:,:) ! Model value u at oblocation.
23 real, allocatable :: model_v(:,:) ! Model value v at oblocation.
24 real, allocatable :: model_t(:,:) ! Model value t at oblocation.
25 real, allocatable :: model_p(:,:) ! Model value p at oblocation.
26 real, allocatable :: model_q(:,:) ! Model value q at oblocation.
27 real, allocatable :: model_hsm(:,:)
29 real :: v_h(kms:kme) ! Model value h at ob hor. location.
30 real :: v_p(kms:kme) ! Model value p at ob hor. location.
35 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_sonde_sfc")
37 allocate (model_u(1,iv%info(sonde_sfc)%n1:iv%info(sonde_sfc)%n2))
38 allocate (model_v(1,iv%info(sonde_sfc)%n1:iv%info(sonde_sfc)%n2))
39 allocate (model_t(1,iv%info(sonde_sfc)%n1:iv%info(sonde_sfc)%n2))
40 allocate (model_p(1,iv%info(sonde_sfc)%n1:iv%info(sonde_sfc)%n2))
41 allocate (model_q(1,iv%info(sonde_sfc)%n1:iv%info(sonde_sfc)%n2))
42 allocate (model_hsm(1,iv%info(sonde_sfc)%n1:iv%info(sonde_sfc)%n2))
45 do n=iv%info(sonde_sfc)%n1,iv%info(sonde_sfc)%n2
46 if (iv%sonde_sfc(n)%u%qc == fails_error_max) iv%sonde_sfc(n)%u%qc = 0
47 if (iv%sonde_sfc(n)%v%qc == fails_error_max) iv%sonde_sfc(n)%v%qc = 0
48 if (iv%sonde_sfc(n)%t%qc == fails_error_max) iv%sonde_sfc(n)%t%qc = 0
49 if (iv%sonde_sfc(n)%p%qc == fails_error_max) iv%sonde_sfc(n)%p%qc = 0
50 if (iv%sonde_sfc(n)%q%qc == fails_error_max) iv%sonde_sfc(n)%q%qc = 0
54 if (sfc_assi_options == sfc_assi_options_1) then
55 do n=iv%info(sonde_sfc)%n1,iv%info(sonde_sfc)%n2
56 ! [1.1] Get horizontal interpolation weights:
58 i = iv%info(sonde_sfc)%i(1,n)
59 j = iv%info(sonde_sfc)%j(1,n)
60 dx = iv%info(sonde_sfc)%dx(1,n)
61 dy = iv%info(sonde_sfc)%dy(1,n)
62 dxm = iv%info(sonde_sfc)%dxm(1,n)
63 dym = iv%info(sonde_sfc)%dym(1,n)
67 iv%sonde_sfc(n)%p%inv = ob%sonde_sfc(n)%p
68 iv%sonde_sfc(n)%t%inv = ob%sonde_sfc(n)%t
69 iv%sonde_sfc(n)%q%inv = ob%sonde_sfc(n)%q
70 iv%sonde_sfc(n)%u%inv = ob%sonde_sfc(n)%u
71 iv%sonde_sfc(n)%v%inv = ob%sonde_sfc(n)%v
73 if (iv % sonde_sfc(n) % h > missing_r) then
75 v_h(k) = dym*(dxm*grid%xb%h(i,j ,k) + dx*grid%xb%h(i+1,j ,k)) &
76 + dy *(dxm*grid%xb%h(i,j+1,k) + dx*grid%xb%h(i+1,j+1,k))
79 hd = v_h(kts) - iv % sonde_sfc(n) % h
81 if (abs(hd) <= Max_StHeight_Diff .or. anal_type_verify ) then
82 if (iv % sonde_sfc(n) % h < v_h(kts)) then
83 iv%info(sonde_sfc)%zk(:,n) = 1.0+1.0e-6
84 call da_obs_sfc_correction(iv%info(sonde_sfc), iv%sonde_sfc(n), n, grid%xb)
87 call da_to_zk(iv % sonde_sfc(n) % h, v_h, v_interp_h, iv%info(sonde_sfc)%zk(1,n))
90 else if (ob % sonde_sfc(n) % p > 1.0) then
92 v_p(k) = dym*(dxm*grid%xb%p(i,j ,k) + dx*grid%xb%p(i+1,j ,k)) &
93 + dy *(dxm*grid%xb%p(i,j+1,k) + dx*grid%xb%p(i+1,j+1,k))
96 call da_to_zk(ob % sonde_sfc(n) % p, v_p, v_interp_p, iv%info(sonde_sfc)%zk(1,n))
98 if (iv%info(sonde_sfc)%zk(1,n) < 0.0 .and. .not. anal_type_verify) then
99 iv % sonde_sfc(n) % p % inv = missing_r
100 iv % sonde_sfc(n) % p % qc = missing_data
101 iv%info(sonde_sfc)%zk(:,n) = 1.0+1.0e-6
106 call da_convert_zk (iv%info(sonde_sfc))
108 if (.not.anal_type_verify ) then
109 do n=iv%info(sonde_sfc)%n1,iv%info(sonde_sfc)%n2
110 if (iv%info(sonde_sfc)%zk(1,n) < 0.0) then
111 iv % sonde_sfc(n) % u % qc = missing_data
112 iv % sonde_sfc(n) % v % qc = missing_data
113 iv % sonde_sfc(n) % t % qc = missing_data
114 iv % sonde_sfc(n) % q % qc = missing_data
115 iv % sonde_sfc(n) % p % qc = missing_data
120 ! [1.2] Interpolate horizontally:
122 call da_interp_lin_3d (grid%xb%u, iv%info(sonde_sfc), model_u,'u')
123 call da_interp_lin_3d (grid%xb%v, iv%info(sonde_sfc), model_v,'v')
125 call da_interp_lin_3d (grid%xb%u, iv%info(sonde_sfc), model_u)
126 call da_interp_lin_3d (grid%xb%v, iv%info(sonde_sfc), model_v)
128 call da_interp_lin_3d (grid%xb%t, iv%info(sonde_sfc), model_t)
129 call da_interp_lin_3d (grid%xb%q, iv%info(sonde_sfc), model_q)
130 call da_interp_lin_3d (grid%xb%p, iv%info(sonde_sfc), model_p)
131 else if (sfc_assi_options == sfc_assi_options_2) then
132 ! 1.2.1 Surface assimilation approach 2(10-m u, v, 2-m t, q, and sfc_p)
134 call da_interp_lin_2d (grid%xb%u10, iv%info(sonde_sfc), 1,model_u)
135 call da_interp_lin_2d (grid%xb%v10, iv%info(sonde_sfc), 1,model_v)
136 call da_interp_lin_2d (grid%xb%t2, iv%info(sonde_sfc), 1,model_t)
137 call da_interp_lin_2d (grid%xb%q2, iv%info(sonde_sfc), 1,model_q)
138 call da_interp_lin_2d (grid%xb%psfc, iv%info(sonde_sfc), 1,model_p)
140 do n=iv%info(sonde_sfc)%n1,iv%info(sonde_sfc)%n2
142 iv%sonde_sfc(n)%p%inv = ob%sonde_sfc(n)%p
143 iv%sonde_sfc(n)%t%inv = ob%sonde_sfc(n)%t
144 iv%sonde_sfc(n)%q%inv = ob%sonde_sfc(n)%q
145 iv%sonde_sfc(n)%u%inv = ob%sonde_sfc(n)%u
146 iv%sonde_sfc(n)%v%inv = ob%sonde_sfc(n)%v
148 if (iv%sonde_sfc(n)%p%qc >= 0) then
149 ! model surface p, t, q, h at observed site:
151 call da_interp_lin_2d_partial (grid%xb%terr, iv%info(sonde_sfc), 1, n, n, model_hsm(:,n))
153 ho = iv%sonde_sfc(n)%h
157 if (iv%sonde_sfc(n)%t%qc >= 0 .and. iv%sonde_sfc(n)%q%qc >= 0) then
158 to = ob%sonde_sfc(n)%t
159 qo = ob%sonde_sfc(n)%q
160 call da_sfc_pre(psfcm, model_p(1,n), model_t(1,n), model_q(1,n), model_hsm(1,n), ho, to, qo)
161 else if (iv%sonde_sfc(n)%t%qc >= 0 .and. iv%sonde_sfc(n)%q%qc < 0) then
162 to = ob%sonde_sfc(n)%t
163 call da_sfc_pre(psfcm, model_p(1,n), model_t(1,n), model_q(1,n), model_hsm(1,n), ho,to)
165 call da_sfc_pre(psfcm, model_p(1,n), model_t(1,n), model_q(1,n), model_hsm(1,n), ho)
168 ! Pressure at the observed height:
174 do n=iv%info(sonde_sfc)%n1,iv%info(sonde_sfc)%n2
175 !-----------------------------------------------------------------------
176 ! [3.0] Fast interpolation:
177 !-----------------------------------------------------------------------
179 if (ob % sonde_sfc(n) % u > missing_r .AND. iv % sonde_sfc(n) % u % qc >= obs_qc_pointer) then
180 iv % sonde_sfc(n) % u % inv = iv%sonde_sfc(n)%u%inv - model_u(1,n)
182 iv % sonde_sfc(n) % u % inv = 0.0
185 if (ob % sonde_sfc(n) % v > missing_r .AND. iv % sonde_sfc(n) % v % qc >= obs_qc_pointer) then
186 iv % sonde_sfc(n) % v % inv = iv%sonde_sfc(n)%v%inv - model_v(1,n)
188 iv % sonde_sfc(n) % v % inv = 0.0
191 if (ob % sonde_sfc(n) % p > 0.0 .AND. iv % sonde_sfc(n) % p % qc >= obs_qc_pointer) then
192 iv % sonde_sfc(n) % p % inv = iv%sonde_sfc(n)%p%inv - model_p(1,n)
194 iv % sonde_sfc(n) % p % inv = 0.0
197 if (ob % sonde_sfc(n) % t > 0.0 .AND. iv % sonde_sfc(n) % t % qc >= obs_qc_pointer) then
198 iv % sonde_sfc(n) % t % inv = iv%sonde_sfc(n)%t%inv - model_t(1,n)
200 iv % sonde_sfc(n) % t % inv = 0.0
203 if (ob % sonde_sfc(n) % q > 0.0 .AND. iv % sonde_sfc(n) % q % qc >= obs_qc_pointer) then
204 iv % sonde_sfc(n) % q % inv = iv%sonde_sfc(n)%q%inv - model_q(1,n)
206 iv % sonde_sfc(n) % q % inv = 0.0
210 !-----------------------------------------------------------------------
211 ! [5.0] Perform optional maximum error check:
212 !-----------------------------------------------------------------------
214 if ( check_max_iv ) &
215 call da_check_max_iv_sonde_sfc(iv,ob, it, num_qcstat_conv)
222 deallocate (model_hsm)
224 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_sonde_sfc")
226 end subroutine da_get_innov_vector_sonde_sfc