da_get_innov_vector_airsr.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_airsr( it, grid, 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 !----------------------------------------------------------------------
9
10 implicit none
11
12 integer, intent(in) :: it ! External iteration.
13 type(domain), intent(in) :: grid ! first guess state.
14 type(y_type), intent(inout) :: ob ! Observation structure.
15 type(iv_type), intent(inout) :: iv ! O-B structure.
16
17 integer :: n ! Loop counter.
18 integer :: i, j, k ! Index dimension.
19 integer :: num_levs ! Number of obs levels.
20 real :: dx, dxm ! Interpolation weights.
21 real :: dy, dym ! Interpolation weights.
22
23 ! real :: model_h(1:max_ob_levels) ! Model value h at ob location.
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, allocatable :: model_q(:,:) ! Model value q at ob location.
29 real, allocatable :: model_t(:,:) ! Model value t at ob location.
30
31 integer :: itt,ittf,itqv,itqvf
32
33 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_airsr")
34
35 itt = 0; itqv = 0;
36 ittf = 0; itqvf = 0;
37
38 allocate (model_t(iv%info(airsr)%max_lev,iv%info(airsr)%n1:iv%info(airsr)%n2))
39 allocate (model_q(iv%info(airsr)%max_lev,iv%info(airsr)%n1:iv%info(airsr)%n2))
40
41 model_t(:,:) = 0.0
42 model_q(:,:) = 0.0
43
44 do n=iv%info(airsr)%n1, iv%info(airsr)%n2
45 num_levs = iv%info(airsr)%levels(n)
46
47 if (num_levs < 1) cycle
48
49 ! [1.1] Get horizontal interpolation weights:
50
51 i = iv%info(airsr)%i(1,n)
52 j = iv%info(airsr)%j(1,n)
53 dx = iv%info(airsr)%dx(1,n)
54 dy = iv%info(airsr)%dy(1,n)
55 dxm = iv%info(airsr)%dxm(1,n)
56 dym = iv%info(airsr)%dym(1,n)
57
58 do k=kts,kte
59 v_h(k) = dym*(dxm*grid%xb%h(i,j,k) + dx*grid%xb%h(i+1,j,k)) + dy *(dxm*grid%xb%h(i,j+1,k) + dx*grid%xb%h(i+1,j+1,k))
60 v_p(k) = dym*(dxm*grid%xb%p(i,j,k) + dx*grid%xb%p(i+1,j,k)) + dy *(dxm*grid%xb%p(i,j+1,k) + dx*grid%xb%p(i+1,j+1,k))
61 end do
62
63 do k=1, num_levs
64 if (iv%airsr(n)%p(k) > 1.0) then
65 call da_to_zk(iv%airsr(n)%p(k), v_p, v_interp_p, iv%info(airsr)%zk(k,n))
66 else if (iv%airsr(n)%h(k) > 0.0) then
67 call da_to_zk(iv%airsr(n)%h(k), v_h, v_interp_h, iv%info(airsr)%zk(k,n))
68 end if
69
70 if (iv%info(airsr)%zk(k,n) < 0.0 .and. .not. anal_type_verify) then
71 iv%airsr(n)%t(k)%qc = missing
72 iv%airsr(n)%q(k)%qc = missing
73 end if
74 end do
75 end do
76
77 call da_convert_zk (iv%info(airsr))
78
79 ! [1.2] Interpolate horizontally to ob:
80 call da_interp_lin_3d (grid%xb%t, iv%info(airsr), model_t)
81 call da_interp_lin_3d (grid%xb%q, iv%info(airsr), model_q)
82
83 do n=iv%info(airsr)%n1,iv%info(airsr)%n2
84 !------------------------------------------------------------------------
85 ! [2.0] Initialise components of innovation vector:
86 !------------------------------------------------------------------------
87
88 do k = 1, iv%info(airsr)%levels(n)
89 iv%airsr(n)%t(k)%inv = 0.0
90 iv%airsr(n)%q(k)%inv = 0.0
91
92 !------------------------------------------------------------------------
93 ! [3.0] Interpolation:
94 !------------------------------------------------------------------------
95
96 if (ob%airsr(n)%t(k) > missing_r .AND. iv%airsr(n)%t(k)%qc >= obs_qc_pointer) then
97 iv%airsr(n)%t(k)%inv = ob%airsr(n)%t(k) - model_t(k,n)
98 end if
99
100 if (ob%airsr(n)%q(k) > missing_r .AND. iv%airsr(n)%q(k)%qc >= obs_qc_pointer) then
101 iv%airsr(n)%q(k)%inv = ob%airsr(n)%q(k) - model_q(k,n)
102 end if
103 end do
104 end do
105
106 !------------------------------------------------------------------------
107 ! [5.0] Perform optional maximum error check:
108 !------------------------------------------------------------------------
109
110 if (check_max_iv) call da_check_max_iv_airsr(iv, it, itt,ittf,itqv,itqvf)
111
112 if (rootproc .and. check_max_iv_print) then
113 write(unit = check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it, &
114 ', Total rejections for airsr follows:'
115
116 write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
117 'Number of failed temperature observations:',ittf, ' on ',itt, &
118 'Number of failed mixing ratio observations:',itqvf,' on ',itqv, &
119 'Final total airsr rejections: ',ittf+itqvf,' on ',itt +itqv
120 end if
121
122 deallocate (model_t)
123 deallocate (model_q)
124
125 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_airsr")
126
127 end subroutine da_get_innov_vector_airsr
128
129