da_get_innov_vector_ssmt2.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_ssmt2 (it, grid, ob, iv)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 integer, intent(in) :: it ! External iteration.
10 type(domain), intent(in) :: grid ! first guess state.
11 type(y_type), intent(inout) :: ob ! Observation structure.
12 type(iv_type), intent(inout) :: iv ! O-B structure.
13
14 integer :: n ! Loop counter.
15 integer :: i, j, k ! Index dimension.
16 integer :: num_levs ! Number of obs levels.
17 real :: dx, dxm ! Interpolation weights.
18 real :: dy, dym ! Interpolation weights.
19 real, allocatable :: model_rh(:,:) ! Model value rh at ob location.
20
21 real :: v_h(kms:kme) ! Model value h at ob hor. location.
22 real :: v_p(kms:kme) ! Model value p at ob hor. location.
23
24 integer :: irh, irhf
25
26 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_ssmt2")
27
28 allocate (model_rh(1:iv%info(ssmt2)%max_lev,iv%info(ssmt2)%n1:iv%info(ssmt2)%n2))
29 model_rh(:,:) = 0.0
30
31 irh = 0 ; irhf = 0
32
33 do n=iv%info(ssmt2)%n1,iv%info(ssmt2)%n2
34 num_levs = iv%info(ssmt2)%levels(n)
35
36 if (num_levs < 1) cycle
37
38 ! [1.1] Get horizontal interpolation weights:
39
40 i = iv%info(ssmt2)%i(1,n)
41 j = iv%info(ssmt2)%j(1,n)
42 dx = iv%info(ssmt2)%dx(1,n)
43 dy = iv%info(ssmt2)%dy(1,n)
44 dxm = iv%info(ssmt2)%dxm(1,n)
45 dym = iv%info(ssmt2)%dym(1,n)
46
47 do k=kts,kte
48 v_h(k) = dym*(dxm*grid%xb%h(i,j ,k) + dx*grid%xb%h(i+1,j ,k)) &
49 + dy *(dxm*grid%xb%h(i,j+1,k) + dx*grid%xb%h(i+1,j+1,k))
50 v_p(k) = dym*(dxm*grid%xb%p(i,j ,k) + dx*grid%xb%p(i+1,j ,k)) &
51 + dy *(dxm*grid%xb%p(i,j+1,k) + dx*grid%xb%p(i+1,j+1,k))
52 end do
53
54 num_levs=0
55 do k=1, iv%info(ssmt2)%levels(n)
56 if (iv % ssmt2(n) % h(k) > 0.0) then
57 call da_to_zk(iv % ssmt2(n) % h(k), v_h, v_interp_h, iv%info(ssmt2)%zk(k,n))
58 else if (iv % ssmt2(n) % p(k) > 1.0) then
59 call da_to_zk(iv % ssmt2(n) % p(k), v_p, v_interp_p, iv%info(ssmt2)%zk(k,n))
60 end if
61
62 if (iv%info(ssmt2)%zk(k,n) > 0.0) then
63 num_levs=num_levs+1
64 iv%info(ssmt2)%zk(num_levs,n)=iv%info(ssmt2)%zk(k,n)
65
66 ob % ssmt2(n) % rh(num_levs) = ob % ssmt2(n) % rh(k)
67
68 iv % ssmt2(n) % rh(num_levs) % qc = iv % ssmt2(n) % rh(k) % qc
69
70 iv % ssmt2(n) % rh(num_levs) % error = iv % ssmt2(n) % rh(k) % error
71 end if
72 end do
73
74 iv%info(ssmt2)%levels(n) = num_levs
75 end do
76
77 call da_convert_zk (iv%info(ssmt2))
78
79 call da_interp_lin_3d (grid%xb%rh, iv%info(ssmt2), model_rh)
80
81 do n=iv%info(ssmt2)%n1, iv%info(ssmt2)%n2
82 do k = 1, iv%info(ssmt2)%levels(n)
83 iv % ssmt2(n) % rh(k) % inv = 0.0
84
85 !-----------------------------------------------------------------
86 ! [3.0] Interpolation:
87 !-----------------------------------------------------------------
88
89 if (ob % ssmt2(n) % rh(k) > missing_r .AND. &
90 iv % ssmt2(n) % rh(k) % qc >= obs_qc_pointer) then
91 iv % ssmt2(n) % rh(k) % inv = ob % ssmt2(n) % rh(k) - model_rh(k,n)
92 end if
93
94 ! write(122,'(2i4,i8,5f15.5)')n, k, iv%ssmt2(n)%height_qc(k), &
95 ! iv%ssmt2(n)%info%lat, iv%ssmt2(n)%info%lon, &
96 ! iv%ssmt2(n)%h(k), &
97 ! ob % ssmt2(n) % rh(k), model_rh(k,n)
98
99 end do
100 end do
101
102 !------------------------------------------------------------------------
103 ! [5.0] Perform optional maximum error check:
104 !------------------------------------------------------------------------
105
106 if (check_max_iv) call da_check_max_iv_ssmt2(iv, it, irh,irhf)
107
108 if (rootproc .and. check_max_iv_print) then
109 write(unit = check_max_iv_unit, fmt ='(A,i5,A)') &
110 'For outer iteration ',it,' Total Rejections for SSMI(T2) follows:'
111 write(unit = check_max_iv_unit, fmt = '(/,2(2(A,I6),/))') &
112 'Number of failed RH(T2) observations: ',irhf, ' on ',irh
113 end if
114
115 deallocate (model_rh)
116
117 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_ssmt2")
118
119 end subroutine da_get_innov_vector_ssmt2
120
121