da_get_innov_vector_polaramv.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_polaramv( it, grid, ob, iv)
2
3 !----------------------------------------------------------------------
4 ! Purpose: Calculates innovation vector, does QC for polaramv
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(in) :: ob ! Observation structure.
12 type(iv_type), intent(inout) :: iv ! O-B structure.
13
14 integer :: n, k ! Loop counter.
15 integer :: i (kms:kme)
16 integer :: j (kms:kme)
17 real :: dx (kms:kme)
18 real :: dxm(kms:kme)
19 real :: dy (kms:kme)
20 real :: dym(kms:kme)
21 real, allocatable :: model_u(:,:) ! Model value u at ob location.
22 real, allocatable :: model_v(:,:) ! Model value v at ob location.
23
24 real :: v_p(kms:kme) ! Model value p at ob hor. location.
25
26 integer :: itu,ituf,itvv,itvvf
27
28 if (trace_use_dull) call da_trace_entry("da_get_innov_vector_polaramv")
29
30 itu = 0; itvv = 0;
31 ituf = 0; itvvf = 0;
32
33 allocate (model_u(iv%info(polaramv)%max_lev,iv%info(polaramv)%n1:iv%info(polaramv)%n2))
34 allocate (model_v(iv%info(polaramv)%max_lev,iv%info(polaramv)%n1:iv%info(polaramv)%n2))
35
36 model_u(:,:) = 0.0
37 model_v(:,:) = 0.0
38
39 do n=iv%info(polaramv)%n1, iv%info(polaramv)%n2
40 if (iv%info(polaramv)%levels(n) < 1) cycle
41
42 ! [1.1] Get horizontal interpolation weights:
43
44 if (position_lev_dependant) then
45 i(:) = iv%info(polaramv)%i(:,n)
46 j(:) = iv%info(polaramv)%j(:,n)
47 dx(:) = iv%info(polaramv)%dx(:,n)
48 dy(:) = iv%info(polaramv)%dy(:,n)
49 dxm(:) = iv%info(polaramv)%dxm(:,n)
50 dym(:) = iv%info(polaramv)%dym(:,n)
51 do k=kts,kte
52 v_p(k) = dym(k)*(dxm(k)*grid%xb%p(i(k),j(k),k) + dx(k)*grid%xb%p(i(k)+1,j(k),k)) &
53 + dy(k) *(dxm(k)*grid%xb%p(i(k),j(k)+1,k) + dx(k)*grid%xb%p(i(k)+1,j(k)+1,k))
54 end do
55 else
56 i(1) = iv%info(polaramv)%i(1,n)
57 j(1) = iv%info(polaramv)%j(1,n)
58 dx(1) = iv%info(polaramv)%dx(1,n)
59 dy(1) = iv%info(polaramv)%dy(1,n)
60 dxm(1) = iv%info(polaramv)%dxm(1,n)
61 dym(1) = iv%info(polaramv)%dym(1,n)
62
63 v_p(kts:kte) = dym(1) * (dxm(1)*grid%xb%p(i(1),j(1),kts:kte) + dx(1)*grid%xb%p(i(1)+1,j(1),kts:kte)) &
64 + dy(1) * (dxm(1)*grid%xb%p(i(1),j(1)+1,kts:kte) + dx(1)*grid%xb%p(i(1)+1,j(1)+1,kts:kte))
65 end if
66
67 do k=1, iv%info(polaramv)%levels(n)
68 if (iv%polaramv(n)%p(k) > 1.0) then
69 call da_to_zk (iv%polaramv(n)%p(k), v_p,v_interp_p, iv%info(polaramv)%zk(k,n))
70 end if
71 end do
72 end do
73
74 call da_convert_zk(iv%info(polaramv))
75
76 if (.not. anal_type_verify) then
77 do n=iv%info(polaramv)%n1,iv%info(polaramv)%n2
78 do k=1, iv%info(polaramv)%levels(n)
79 if (iv%info(polaramv)%zk(k,n) < 0.0) then
80 iv%polaramv(n)%u(k)%qc = missing
81 iv%polaramv(n)%v(k)%qc = missing
82 end if
83 end do
84 end do
85 end if
86
87 ! [1.2] Interpolate horizontally to ob:
88
89 call da_interp_lin_3d (grid%xb%u, iv%info(polaramv), model_u)
90 call da_interp_lin_3d (grid%xb%v, iv%info(polaramv), model_v)
91
92 do n=iv%info(polaramv)%n1,iv%info(polaramv)%n2
93 do k = 1, iv%info(polaramv)%levels(n)
94 iv % polaramv(n) % u(k) % inv = 0.0
95 iv % polaramv(n) % v(k) % inv = 0.0
96
97 if (ob % polaramv(n) % u(k) > missing_r .AND. &
98 iv % polaramv(n) % u(k) % qc >= obs_qc_pointer) then
99 iv % polaramv(n) % u(k) % inv = ob % polaramv(n) % u(k) - model_u(k,n)
100 end if
101
102 if (ob % polaramv(n) % v(k) > missing_r .AND. &
103 iv % polaramv(n) % v(k) % qc >= obs_qc_pointer) then
104 iv % polaramv(n) % v(k) % inv = ob % polaramv(n) % v(k) - model_v(k,n)
105 end if
106 end do
107 end do
108
109 !------------------------------------------------------------------------
110 ! Perform optional maximum error check:
111 !------------------------------------------------------------------------
112
113 if (check_max_iv) call da_check_max_iv_polaramv(iv, it, itu,ituf,itvv,itvvf)
114
115 if (rootproc .and. check_max_iv_print) then
116 write(unit = check_max_iv_unit, fmt ='(A,i5,A)')'For outer iteration ',it, &
117 ', Total Rejections for Polar AMVs follows:'
118 write(unit = check_max_iv_unit, fmt = '(/,10(2(A,I6),/))') &
119 'Number of failed u-wind observations: ',ituf, ' on ',itu, &
120 'Number of failed v-wind observations: ',itvvf,' on ',itvv, &
121 'Finally Total Polar AMVs rejections ',ituf+itvvf,' on ',itu+itvv
122 end if
123
124 deallocate (model_u)
125 deallocate (model_v)
126
127 if (trace_use_dull) call da_trace_exit("da_get_innov_vector_polaramv")
128
129 end subroutine da_get_innov_vector_polaramv
130
131