da_get_innov_vector_geoamv.inc

References to this file elsewhere.
1 subroutine da_get_innov_vector_geoamv( it, xb, xp, ob, iv)
2 
3    !-------------------------------------------------------------------------
4    ! Purpose: Calculates innovation vector does QC for  Geo. AMV's
5    !               
6    ! Update :
7    !     01/24/2007    Syed RH Rizvi
8    !     Updated for "VERIFY"       
9    !-------------------------------------------------------------------------
10 
11    implicit none
12 
13    integer, intent(in)           :: it      ! External iteration.
14    type(xb_type), intent(in)    :: xb      ! first guess state.
15    type(xpose_type), intent(in) :: xp      ! Domain decomposition vars.
16    type(y_type),  intent(in)    :: ob      ! Observation structure.
17    type(ob_type), intent(inout) :: iv      ! O-B structure.
18 
19    integer                        :: n        ! Loop counter.
20    integer                        :: i, j, k  ! Index dimension.
21    integer                        :: num_levs ! Number of obs levels.
22 
23    real                         :: dx, dxm  ! Interpolation weights.
24    real                         :: dy, dym  ! Interpolation weights.
25    real, dimension(1:max_ob_levels) :: model_u  ! Model value u at ob location.
26    real, dimension(1:max_ob_levels) :: model_v  ! Model value v at ob location.
27 
28    real, dimension(xp%kms:xp%kme) :: v_p      ! Model value p at ob hor. location.
29 
30    integer           :: itu,ituf,itvv,itvvf
31 
32    if (iv % num_geoamv > 0) then
33       itu   = 0
34       itvv  = 0
35       ituf  = 0
36       itvvf = 0
37 
38       do n=iv%ob_numb(iv%current_ob_time-1)%geoamv + 1, iv%ob_numb(iv%current_ob_time)%geoamv
39 
40          ! [1.3] Get horizontal interpolation weights:
41 
42          num_levs = iv % geoamv(n) % info % levels
43          if (num_levs < 1) cycle
44 
45          model_u(:) = 0.0
46          model_v(:) = 0.0
47 
48          i = iv%geoamv(n)%loc%i
49          j = iv%geoamv(n)%loc%j
50          dx = iv%geoamv(n)%loc%dx
51          dy = iv%geoamv(n)%loc%dy
52          dxm = iv%geoamv(n)%loc%dxm
53          dym = iv%geoamv(n)%loc%dym
54 
55          do k=xp%kts,xp%kte
56             v_p(k) = dym*(dxm*xb%p(i,j  ,k) + dx*xb%p(i+1,j  ,k)) &
57                    + dy *(dxm*xb%p(i,j+1,k) + dx*xb%p(i+1,j+1,k))
58          end do
59 
60          do k=1, iv % geoamv(n) % info % levels
61             iv%geoamv(n)%zk(k)=missing_r
62 
63             if (iv % geoamv(n) % p(k) > 1.0) then
64                call da_to_zk(iv % geoamv(n) % p(k), v_p, xp, v_interp_p, iv%geoamv(n)%zk(k))
65             end if
66 
67             if (iv%geoamv(n)%zk(k) < 0.0 .and.  .not.anal_type_verify) then
68                iv % geoamv(n) % u(k) % qc = missing
69                iv % geoamv(n) % v(k) % qc = missing
70             end if
71          end do
72 
73          call da_interp_lin_3d( xb % u, xp, i, j, dx, dy, dxm, dym, &
74                              model_u, max_ob_levels,iv%geoamv(n)%zk,num_levs)
75          call da_interp_lin_3d( xb % v, xp, i, j, dx, dy, dxm, dym, &
76                              model_v, max_ob_levels,iv%geoamv(n)%zk,num_levs)
77 
78          do k = 1, iv % geoamv(n) % info % levels
79             iv % geoamv(n) % u(k) % inv = 0.0
80             iv % geoamv(n) % v(k) % inv = 0.0
81             if (ob % geoamv(n) % u(k) > missing_r .AND. &
82                  iv % geoamv(n) % u(k) % qc >= obs_qc_pointer) then
83                iv % geoamv(n) % u(k) % inv = ob % geoamv(n) % u(k) - model_u(k)
84             end if
85 
86             if (ob % geoamv(n) % v(k) > missing_r .AND. &
87                  iv % geoamv(n) % v(k) % qc >= obs_qc_pointer) then
88                iv % geoamv(n) % v(k) % inv = ob % geoamv(n) % v(k) - model_v(k)
89             end if
90          end do
91 
92          !------------------------------------------------------------------------
93          ! Perform optional maximum error check:
94          !------------------------------------------------------------------------
95 
96          if (check_max_iv) then
97            call da_check_max_iv_geoamv(it, iv % geoamv(n),itu,ituf,itvv,itvvf)
98          end if
99       end do
100 
101       if (rootproc .and. check_max_iv_print) then
102          write(unit = check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it, &
103             ', Total Rejections for Geo. AMVs follows:'
104          write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
105             'Number of failed u-wind observations:     ',ituf, ' on ',itu,   &
106             'Number of failed v-wind observations:     ',itvvf,' on ',itvv,   &
107             'Finally Total Geo. AMVs rejections ',ituf+itvvf,' on ',itu+itvv
108       end if
109    end if
110 
111 end subroutine da_get_innov_vector_geoamv
112 
113