da_get_innov_vector_airep.inc

References to this file elsewhere.
1 subroutine da_get_innov_vector_airep( 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    integer                        :: num_levs ! Number of obs levels.
21    real                           :: dx, dxm  ! Interpolation weights.
22    real                           :: dy, dym  ! Interpolation weights.
23    real, dimension(1:max_ob_levels) :: model_u  ! Model value u at ob location.
24    real, dimension(1:max_ob_levels) :: model_v  ! Model value v at ob location.
25    real, dimension(1:max_ob_levels) :: model_t  ! Model value t at ob location.
26 
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    integer           :: itu,ituf,itvv,itvvf,itt,ittf
32 
33    if (iv % num_airep < 1) return
34    
35    if (trace_use) call da_trace_entry("da_get_innov_vector_airep")
36 
37    itu   = 0; itvv    = 0; itt  = 0;
38    ituf  = 0; itvvf   = 0; ittf = 0;
39 
40    do n=iv%ob_numb(iv%current_ob_time-1)%airep + 1, iv%ob_numb(iv%current_ob_time)%airep
41       num_levs = iv % airep(n) % info % levels
42       if (num_levs < 1) cycle
43 
44       model_u(:) = 0.0
45       model_v(:) = 0.0
46       model_t(:) = 0.0
47 
48       ! [1.1] Get horizontal interpolation weights:
49 
50       i = iv%airep(n)%loc%i
51       j = iv%airep(n)%loc%j
52       dx = iv%airep(n)%loc%dx
53       dy = iv%airep(n)%loc%dy
54       dxm = iv%airep(n)%loc%dxm
55       dym = iv%airep(n)%loc%dym
56 
57       do k=xp%kts,xp%kte
58          v_h(k) = dym*(dxm*xb%h(i,j  ,k) + dx*xb%h(i+1,j  ,k)) &
59                  + dy *(dxm*xb%h(i,j+1,k) + dx*xb%h(i+1,j+1,k))
60          v_p(k) = dym*(dxm*xb%p(i,j  ,k) + dx*xb%p(i+1,j  ,k)) &
61                  + dy *(dxm*xb%p(i,j+1,k) + dx*xb%p(i+1,j+1,k))
62       end do
63 
64       num_levs=0
65       do k=1, iv % airep(n) % info % levels
66          iv%airep(n)%zk(k)=missing_r
67          if (iv % airep(n) % p(k) > 1.0) then
68             call da_to_zk(iv % airep(n) % p(k), v_p, xp, v_interp_p, &
69                iv%airep(n)%zk(k))
70          else if (iv % airep(n) % h(k) > 0.0) then
71             call da_to_zk(iv % airep(n) % h(k), v_h, xp, v_interp_h, &
72                iv%airep(n)%zk(k))
73          end if
74 
75          if (iv%airep(n)%zk(k) < 0.0 .and. .not.anal_type_verify ) then
76             iv % airep(n) % u(k) % qc = missing
77             iv % airep(n) % v(k) % qc = missing
78             iv % airep(n) % t(k) % qc = missing
79          end if
80       end do
81 
82       num_levs = iv % airep(n) % info % levels
83 
84       ! [1.2] Interpolate horizontally:
85       call da_interp_lin_3d( xb % u, xp, i, j, dx, dy, dxm, dym, &
86          model_u, max_ob_levels, iv%airep(n)%zk, num_levs)
87 
88       call da_interp_lin_3d( xb % v, xp, i, j, dx, dy, dxm, dym, &
89          model_v, max_ob_levels, iv%airep(n)%zk, num_levs)
90 
91       call da_interp_lin_3d( xb % t, xp, i, j, dx, dy, dxm, dym, &
92           model_t, max_ob_levels, iv%airep(n)%zk, num_levs)
93 
94       !-------------------------------------------------------------------
95       ! [2.0] Initialise components of innovation vector:
96       !-------------------------------------------------------------------
97 
98       do k = 1, iv % airep(n) % info % levels
99          iv % airep(n) % u(k) % inv = 0.0
100          iv % airep(n) % v(k) % inv = 0.0
101          iv % airep(n) % t(k) % inv = 0.0
102 
103          !----------------------------------------------------------------
104          ! [3.0] Fast interpolation:
105          !----------------------------------------------------------------
106 
107          if (ob % airep(n) % u(k) > missing_r .AND. &
108               iv % airep(n) % u(k) % qc >= obs_qc_pointer) then
109             iv % airep(n) % u(k) % inv = ob % airep(n) % u(k) - model_u(k)
110          end if
111 
112          if (ob % airep(n) % v(k) > missing_r .AND. &
113               iv % airep(n) % v(k) % qc >= obs_qc_pointer) then
114             iv % airep(n) % v(k) % inv = ob % airep(n) % v(k) - model_v(k)
115          end if
116 
117          if (ob % airep(n) % t(k) > missing_r .AND. &
118               iv % airep(n) % t(k) % qc >= obs_qc_pointer) then
119             iv % airep(n) % t(k) % inv = ob % airep(n) % t(k) - model_t(k)
120          end if
121       end do
122 
123       !-------------------------------------------------------------------
124       ! [5.0] Perform optional maximum error check:
125       !-------------------------------------------------------------------
126 
127       if (check_max_iv)then
128          call da_check_max_iv_airep(it, iv % airep(n),itu,ituf,itvv, &
129             itvvf,itt,ittf)
130       end if
131    end do
132 
133    if (rootproc .and. check_max_iv_print) then
134       write(unit = check_max_iv_unit, fmt ='(A,i5,A)') &
135          'For outer iteration ',it, &
136          ', Total Rejections for Airep follows:'
137       write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
138          'Number of failed u-wind observations:     ',ituf, ' on ',itu,   &
139          'Number of failed v-wind observations:     ',itvvf,' on ',itvv,   &
140          'Number of failed temperature observations:',ittf, ' on ',itt,   &
141          'Finally Total Airep rejections ',ituf+itvvf+ittf,' on ', &
142          itu +itvv +itt
143    end if
144    
145    if (trace_use) call da_trace_exit("da_get_innov_vector_airep")
146 
147 end subroutine da_get_innov_vector_airep
148 
149