da_vertical_transform.inc
References to this file elsewhere.
1 subroutine da_vertical_transform(string, be, vertical_wgt, vv, vp)
2
3 !---------------------------------------------------------------------
4 ! Purpose: TBD
5 !---------------------------------------------------------------------
6
7 implicit none
8
9 character (len=*), intent(in) :: string ! Character operation
10 type (be_type), intent(in) :: be ! Background error structure.
11 real, intent(in) :: vertical_wgt(ims:ime,jms:jme,kms:kme) ! Weighting.
12 type (vp_type), intent(inout) :: vv ! CV in gridpt/EOF space.
13 type (vp_type), intent(inout) :: vp ! CV in gridpt/level space.
14
15 integer :: j, m ! Loop counters.
16 real :: alpha_stddev_inv ! 1/ sigma_alpha
17
18 if (trace_use) call da_trace_entry("da_vertical_transform")
19
20 select case(string)
21
22 case ('u');
23
24 !-------------------------------------------------------------------
25 ! [1.0] Perform vp(i,j,k) = E L^{1/2} vv(i,j,m) transform:
26 !-------------------------------------------------------------------
27
28 if (be % v1 % mz > 0) then
29 call da_transform_vvtovp(&
30 be % v1 % evec, &
31 be % v1 % val, &
32 vertical_wgt, &
33 vv % v1, vp % v1, &
34 be % v1 % mz, kte)
35 else
36 vp % v1(its:ite,jts:jte,kts:kte) = 0.0
37 end if
38
39 if (be % v2 % mz > 0) then
40 call da_transform_vvtovp(&
41 be % v2 % evec, &
42 be % v2 % val, &
43 vertical_wgt, &
44 vv % v2, vp % v2, &
45 be % v2 % mz, kte)
46 else
47 vp % v2(its:ite,jts:jte,kts:kte) = 0.0
48 end if
49
50 if (be % v3 % mz > 0) then
51 call da_transform_vvtovp(&
52 be % v3 % evec, &
53 be % v3 % val, &
54 vertical_wgt, &
55 vv % v3, vp % v3, &
56 be % v3 % mz, kte)
57 else
58 vp % v3(its:ite,jts:jte,kts:kte) = 0.0
59 end if
60
61 if (be % v4 % mz > 0) then
62 call da_transform_vvtovp(&
63 be % v4 % evec, &
64 be % v4 % val, &
65 vertical_wgt, &
66 vv % v4, vp % v4, &
67 be % v4 % mz, kte)
68 else
69 vp % v4(its:ite,jts:jte,kts:kte) = 0.0
70 end if
71
72 if (be % v5 % mz > 0) then
73 if (global) then
74 vp % v5(its:ite,jts:jte,1) = vv % v5(its:ite,jts:jte,1)
75 else
76 call da_transform_vvtovp(&
77 be % v5 % evec, &
78 be % v5 % val, &
79 vertical_wgt, &
80 vv % v5, vp % v5, &
81 be % v5 % mz, kts)
82 end if
83 else
84 vp % v5(its:ite,jts:jte,kts:kts) = 0.0
85 end if
86
87 if (be % ne > 0) then
88 do m = 1, be % ne
89 do j = jts, jte
90 vp % alpha(its:ite,j,m) = vv % alpha(its:ite,j,m) * be % alpha % val(j,m)
91 end do
92 end do
93 end if
94
95 case ('u_inv');
96
97 !-------------------------------------------------------------------
98 ! [2.0] Perform vv(i,j,m) = L^{-1/2} E^T vp(i,j,k) transform:
99 !-------------------------------------------------------------------
100
101 if (be % v1 % mz > 0) then
102 call da_transform_vptovv(&
103 be % v1 % evec, &
104 be % v1 % val, &
105 vertical_wgt, &
106 vp % v1, vv % v1, &
107 be % v1 % mz, &
108 kds,kde, &
109 ims,ime, jms,jme, kms,kme, &
110 its,ite, jts,jte, kts,kte)
111 end if
112
113 if (be % v2 % mz > 0) then
114 call da_transform_vptovv(&
115 be % v2 % evec, &
116 be % v2 % val, &
117 vertical_wgt, &
118 vp % v2, vv % v2, &
119 be % v2 % mz, &
120 kds,kde, &
121 ims,ime, jms,jme, kms,kme, &
122 its,ite, jts,jte, kts,kte)
123 end if
124
125 if (be % v3 % mz > 0) then
126 call da_transform_vptovv(&
127 be % v3 % evec, &
128 be % v3 % val, &
129 vertical_wgt, &
130 vp % v3, vv % v3, &
131 be % v3 % mz, &
132 kds,kde, &
133 ims,ime, jms,jme, kms,kme, &
134 its,ite, jts,jte, kts,kte)
135 end if
136
137 if (be % v4 % mz > 0) then
138 call da_transform_vptovv(&
139 be % v4 % evec, &
140 be % v4 % val, &
141 vertical_wgt, &
142 vp % v4, vv % v4, &
143 be % v4 % mz, &
144 kds,kde, &
145 ims,ime, jms,jme, kms,kme, &
146 its,ite, jts,jte, kts,kte)
147 end if
148
149 if (be % v5 % mz > 0) then
150 if(global) then
151 vv % v5(its:ite,jts:jte,1) = vp % v5(its:ite,jts:jte,1)
152 else
153 call da_transform_vptovv(&
154 be % v5 % evec, &
155 be % v5 % val, &
156 vertical_wgt, &
157 vp % v5, vv % v5, &
158 be % v5 % mz, &
159 kds,kde, &
160 ims,ime, jms,jme, kms,kme, &
161 its,ite, jts,jte, kts,kte)
162 end if
163 end if
164
165 if (be % ne > 0) then
166 do m = 1, be % ne
167 do j = jts, jte
168 alpha_stddev_inv = 1.0 / be % alpha % val(j,m)
169 vv % alpha(its:ite,j,m) = vp % alpha(its:ite,j,m) * alpha_stddev_inv
170 end do
171 end do
172 end if
173
174 case ('u_adj');
175
176 !-------------------------------------------------------------------
177 ! [3.0] Perform vv_adj = U_{v}^{T} vp_adj transform:
178 !-------------------------------------------------------------------
179
180 if (be % v1 % mz > 0) then
181 call da_transform_vvtovp_adj(&
182 be % v1 % evec, &
183 be % v1 % val, &
184 vertical_wgt, &
185 vp % v1, &
186 vv % v1, be % v1 % mz, kte)
187 end if
188
189 if (be % v2 % mz > 0) then
190 call da_transform_vvtovp_adj(&
191 be % v2 % evec, &
192 be % v2 % val, &
193 vertical_wgt, &
194 vp % v2, &
195 vv % v2, be % v2 % mz, kte)
196 end if
197
198 if (be % v3 % mz > 0) then
199 call da_transform_vvtovp_adj(&
200 be % v3 % evec, &
201 be % v3 % val, &
202 vertical_wgt, &
203 vp % v3, &
204 vv % v3, be % v3 % mz, kte)
205 end if
206
207 if (be % v4 % mz > 0) then
208 call da_transform_vvtovp_adj(&
209 be % v4 % evec, &
210 be % v4 % val, &
211 vertical_wgt, &
212 vp % v4, &
213 vv % v4, be % v4 % mz, kte)
214 end if
215
216 if (be % v5 % mz > 0) then
217 if (global) then
218 vv % v5(its:ite,jts:jte,1) = vp % v5(its:ite,jts:jte,1)
219 else
220 call da_transform_vvtovp_adj(&
221 be % v5 % evec, &
222 be % v5 % val, &
223 vertical_wgt, &
224 vp % v5, &
225 vv % v5, be % v5 % mz, kts)
226 end if
227 end if
228
229 if (be % ne > 0) then
230 do m = 1, be % ne
231 do j = jts, jte
232 vv % alpha(its:ite,j,m) = vp % alpha(its:ite,j,m) * be % alpha % val(j,m)
233 end do
234 end do
235 end if
236
237 case default;
238
239 call da_error(__FILE__,__LINE__, &
240 (/"Invalid da_vertical_transform option "//trim(string)/))
241
242 end select
243
244 if (trace_use) call da_trace_exit("da_vertical_transform")
245
246 end subroutine da_vertical_transform
247
248