da_get_innov_vector_sonde_sfc.inc

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