da_get_innov_vector_airsr.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_airsr( it, xb, xp, ob, iv)
2
3 !----------------------------------------------------------------------
4 ! Purpose: a) Rcomputes innovation vecrot at
5 ! AIRS retrieval locations
6 ! b) Does quality control check on innovation vector
7 ! if required.
8 ! Update :
9 ! 01/24/2007 Syed RH Rizvi
10 ! Updated for "VERIFY"
11 !----------------------------------------------------------------------
12
13 implicit none
14
15 integer, intent(in) :: it ! External iteration.
16 type(xb_type), intent(in) :: xb ! first guess state.
17 type(xpose_type), intent(in) :: xp ! Domain decomposition vars.
18 type(y_type), intent(inout) :: ob ! Observation structure.
19 type(ob_type), intent(inout) :: iv ! O-B structure.
20
21 integer :: n ! Loop counter.
22 integer :: i, j, k ! Index dimension.
23 integer :: num_levs ! Number of obs levels.
24 real :: dx, dxm ! Interpolation weights.
25 real :: dy, dym ! Interpolation weights.
26
27 real, dimension(1:max_ob_levels) :: model_t ! Model value t at ob location.
28 real, dimension(1:max_ob_levels) :: model_q ! Model value q at ob location.
29 ! real, dimension(1:max_ob_levels) :: model_h ! Model value h at ob location.
30
31 real, dimension(xp%kms:xp%kme) :: v_h ! Model value h at ob hor. location.
32 real, dimension(xp%kms:xp%kme) :: v_p ! Model value p at ob hor. location.
33
34 integer :: itt,ittf,itqv,itqvf
35
36 if (iv%num_airsr < 1) return
37
38 if (trace_use) call da_trace_entry("da_get_innov_vector_airsr")
39
40 itt = 0; itqv = 0;
41 ittf = 0; itqvf = 0;
42
43 do n=iv%ob_numb(iv%current_ob_time-1)%airsr + 1, iv%ob_numb(iv%current_ob_time)%airsr
44 num_levs = iv%airsr(n)%info%levels
45
46 if (num_levs < 1) cycle
47
48 model_t(:) = 0.0
49 model_q(:) = 0.0
50 ! model_h(:) = 0.0
51
52 ! [1.1] Get horizontal interpolation weights:
53
54 i = iv%airsr(n)%loc%i
55 j = iv%airsr(n)%loc%j
56 dx = iv%airsr(n)%loc%dx
57 dy = iv%airsr(n)%loc%dy
58 dxm = iv%airsr(n)%loc%dxm
59 dym = iv%airsr(n)%loc%dym
60
61 do k=xp%kts,xp%kte
62 v_h(k) = dym*(dxm*xb%h(i,j ,k) + dx*xb%h(i+1,j ,k)) &
63 + dy *(dxm*xb%h(i,j+1,k) + dx*xb%h(i+1,j+1,k))
64 v_p(k) = dym*(dxm*xb%p(i,j ,k) + dx*xb%p(i+1,j ,k)) &
65 + dy *(dxm*xb%p(i,j+1,k) + dx*xb%p(i+1,j+1,k))
66 end do
67
68 do k=1, num_levs
69 iv%airsr(n)%zk(k)=missing_r
70
71 if (iv%airsr(n)%p(k) > 1.0) then
72 call da_to_zk(iv%airsr(n)%p(k), v_p, xp, v_interp_p, iv%airsr(n)%zk(k))
73 else if (iv%airsr(n)%h(k) > 0.0) then
74 call da_to_zk(iv%airsr(n)%h(k), v_h, xp, v_interp_h, iv%airsr(n)%zk(k))
75 end if
76
77 if (iv%airsr(n)%zk(k) < 0.0 .and. .not.anal_type_verify) then
78 iv%airsr(n)%t(k)%qc = missing
79 iv%airsr(n)%q(k)%qc = missing
80 end if
81 end do
82
83 ! [1.2] Interpolate horizontally to ob:
84
85 call da_interp_lin_3d( xb%t, xp, i, j, dx, dy, dxm, dym, &
86 model_t, max_ob_levels, iv%airsr(n)%zk, num_levs)
87 call da_interp_lin_3d( xb%q, xp, i, j, dx, dy, dxm, dym, &
88 model_q, max_ob_levels, iv%airsr(n)%zk, num_levs)
89
90 ! call da_interp_lin_3d( xb%h, xp, i, j, dx, dy, dxm, dym, &
91 ! model_h, max_ob_levels, iv%airsr(n)%zk, num_levs)
92
93 !------------------------------------------------------------------------
94 ! [2.0] Initialise components of innovation vector:
95 !------------------------------------------------------------------------
96
97 do k = 1, iv%airsr(n)%info%levels
98 iv%airsr(n)%t(k)%inv = 0.0
99 iv%airsr(n)%q(k)%inv = 0.0
100
101 !------------------------------------------------------------------------
102 ! [3.0] Interpolation:
103 !------------------------------------------------------------------------
104
105 if (ob%airsr(n)%t(k) > missing_r .AND. &
106 iv%airsr(n)%t(k)%qc >= obs_qc_pointer) then
107 iv%airsr(n)%t(k)%inv = ob%airsr(n)%t(k) - model_t(k)
108 end if
109
110 if (ob%airsr(n)%q(k) > missing_r .AND. &
111 iv%airsr(n)%q(k)%qc >= obs_qc_pointer) then
112 iv%airsr(n)%q(k)%inv = ob%airsr(n)%q(k) - model_q(k)
113 end if
114 end do
115
116 !------------------------------------------------------------------------
117 ! [5.0] Perform optional maximum error check:
118 !------------------------------------------------------------------------
119
120 if (check_max_iv) &
121 call da_check_max_iv_airsr(it, iv % airsr(n),itt,ittf,itqv,itqvf)
122 end do
123
124 if (rootproc .and. check_max_iv_print) then
125 write(unit = check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it, &
126 ', Total Rejections for airsr follows:'
127
128 write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
129 'Number of failed temperature observations:',ittf, ' on ',itt, &
130 'Number of failed mixing ratio observations:',itqvf,' on ',itqv, &
131 'Finally Total airsr rejections: ',ittf+itqvf,' on ',itt +itqv
132 end if
133
134 if (trace_use) call da_trace_exit("da_get_innov_vector_airsr")
135
136 end subroutine da_get_innov_vector_airsr
137
138