da_get_innov_vector_qscat.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_qscat( 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
25 real :: zk ! Interpolation vertical coordinator.
26
27 real, dimension(xp%kms:xp%kme) :: v_h ! Model value h at ob hor. location.
28
29 integer :: itu, itvv, ituf, itvvf
30
31 if (iv % num_qscat < 1) return
32
33 if (trace_use) call da_trace_entry("da_get_innov_vector_qscat")
34
35 do n=iv%ob_numb(iv%current_ob_time-1)%qscat + 1, iv%ob_numb(iv%current_ob_time)%qscat
36
37 itu = 0 ; itvv = 0 ; ituf = 0 ; itvvf = 0
38
39 ! [1.1] Get horizontal interpolation weights:
40
41 i = iv%qscat(n)%loc%i
42 j = iv%qscat(n)%loc%j
43 dx = iv%qscat(n)%loc%dx
44 dy = iv%qscat(n)%loc%dy
45 dxm = iv%qscat(n)%loc%dxm
46 dym = iv%qscat(n)%loc%dym
47
48 do k=xp%kts,xp%kte
49 v_h(k) = dym*(dxm*xb%h(i,j ,k) + dx*xb%h(i+1,j ,k)) &
50 + dy *(dxm*xb%h(i,j+1,k) + dx*xb%h(i+1,j+1,k))
51 end do
52
53 zk=missing_r
54
55 if (iv % qscat(n) % h > missing_r) then
56 call da_to_zk(iv % qscat(n) % h, v_h, xp, v_interp_h, zk)
57 if (zk < 1.0) then
58 zk = 1.0
59 end if
60 end if
61
62 iv%qscat(n)%zk=zk
63
64 !------------------------------------------------------------------------
65 ! [2.0] Initialise components of innovation vector:
66 !------------------------------------------------------------------------
67
68 iv % qscat(n) % u % inv = 0.0
69 iv % qscat(n) % v % inv = 0.0
70
71 if (zk < 0.0 .and. .not.anal_type_verify) then
72 iv % qscat(n) % u % qc = missing
73 iv % qscat(n) % v % qc = missing
74 else
75
76 ! [1.2] Interpolate horizontally:
77 call da_interp_obs_lin_2d( xb % u, xp, i, j, dx, dy, dxm, dym, &
78 model_u, iv%qscat(n)%zk)
79 call da_interp_obs_lin_2d( xb % v, xp, i, j, dx, dy, dxm, dym, &
80 model_v, iv%qscat(n)%zk)
81
82 !------------------------------------------------------------------------
83 ! [2.0] Initialise components of innovation vector:
84 !------------------------------------------------------------------------
85
86 !------------------------------------------------------------------------
87 ! [3.0] Fast interpolation:
88 !------------------------------------------------------------------------
89
90 if (ob % qscat(n) % u > missing_r .AND. &
91 iv % qscat(n) % u % qc >= obs_qc_pointer) then
92
93 iv % qscat(n) % u % inv = ob % qscat(n) % u - model_u
94 end if
95
96 if (ob % qscat(n) % v > missing_r .AND. &
97 iv % qscat(n) % v % qc >= obs_qc_pointer) then
98
99 iv % qscat(n) % v % inv = ob % qscat(n) % v - model_v
100 end if
101
102 !------------------------------------------------------------------------
103 ! [5.0] Perform optional maximum error check:
104 !------------------------------------------------------------------------
105
106 if (check_max_iv)&
107 call da_check_max_iv_qscat(it, iv % qscat(n), itu, itvv, ituf,itvvf)
108 end if
109 end do
110
111 if (rootproc .and. check_max_iv_print) then
112 write(unit = check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it,&
113 ', Total Rejections for qscat follows:'
114 write(unit = check_max_iv_unit, fmt = '(/,4(2(A,I6),/))') &
115 'Number of failed u-wind observations: ',ituf, ' on ',itu, &
116 'Number of failed v-wind observations: ',itvvf,' on ',itvv
117 end if
118
119 if (trace_use) call da_trace_exit("da_get_innov_vector_qscat")
120
121 end subroutine da_get_innov_vector_qscat
122
123