da_get_innov_vector_sonde_sfc.inc

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