da_transform_vtox.inc
References to this file elsewhere.
1 subroutine da_transform_vtox(cv_size, xb, xbx, be, ep, cv, vv, vp, xp, xa)
2
3 !--------------------------------------------------------------------------
4 ! Purpose: Control variable transform x' = Uv.
5 !--------------------------------------------------------------------------
6
7 implicit none
8
9 integer, intent(in) :: cv_size ! Size of cv array.
10 type (xb_type), intent(in) :: xb ! first guess.
11 type (xbx_type),intent(in) :: xbx ! For header & non-grid arrays.
12 type (be_type), intent(in) :: be ! background errors.
13 type (ep_type), intent(in) :: ep ! Ensemble perturbations.
14 real, intent(in) :: cv(1:cv_size) ! control variables.
15 type (vp_type),intent(out) :: vv ! grdipt/eof cv (local).
16 type (vp_type),intent(inout) :: vp ! grdipt/level cv (local).
17 type (xpose_type), intent(inout) :: xp ! Dimensions and xpose buffers.
18 type (x_type), intent(out) :: xa ! model space analysis (local).
19
20 integer :: i, j, k
21
22 real :: sdmd, s1md, mu
23 real, dimension(xb%kms:xb%kme) :: p, mr_a, mr_b
24 real :: PU, PD, coeff
25
26 if (trace_use) call da_trace_entry("da_transform_vtox")
27
28 call da_zero_x (xa)
29
30 !----------------------------------------------------------------------
31 ! [1.0]: Perform vv = u_h cv transform:
32 !----------------------------------------------------------------------
33
34 if (global) then
35 call da_transform_vtovv_global(cv_size, xbx, be, cv, vv)
36 else
37 call da_transform_vtovv(cv_size, xb, be, cv, vv, xp)
38 end if
39
40 !----------------------------------------------------------------------
41 ! [2.0]: Perform vp = u_v vv transform:
42 !----------------------------------------------------------------------
43
44 if (vert_corr == 2) then
45 call da_vertical_transform('u', be, &
46 xb % vertical_inner_product, &
47 vv, vp)
48
49 else
50 vp % v1(its:ite,jts:jte,kts:kte) = vv % v1(its:ite,jts:jte,kts:kte)
51 vp % v2(its:ite,jts:jte,kts:kte) = vv % v2(its:ite,jts:jte,kts:kte)
52 vp % v3(its:ite,jts:jte,kts:kte) = vv % v3(its:ite,jts:jte,kts:kte)
53 vp % v4(its:ite,jts:jte,kts:kte) = vv % v4(its:ite,jts:jte,kts:kte)
54 vp % v5(its:ite,jts:jte,kts:kte) = vv % v5(its:ite,jts:jte,kts:kte)
55 end if
56
57 ! WHY
58 ! if (be % ne > 0) then
59 ! vp % alpha(its:ite,jts:jte,1:be%ne) = vv%alpha(its:ite,jts:jte,1:be%ne)
60 ! end if
61
62 !----------------------------------------------------------------------
63 ! [3.0]: Perform x = u_p vp transform::
64 !----------------------------------------------------------------------
65
66 call da_transform_vptox(xb, vp, xp, xa, be, ep)
67
68 !----------------------------------------------------------------------
69 ! [4.0]: Move the following:
70 !----------------------------------------------------------------------
71
72 do j=xb%jts,xb%jte
73 do i=xb%its,xb%ite
74 if (fg_format == fg_format_wrf) then
75 sdmd=0.0
76 s1md=0.0
77 do k=xb%kts,xb%kte
78 mr_a(k) = xa%q(i,j,k)/(1.0 - xb%q(i,j,k))**2
79 mr_b(k) = xb%q(i,j,k)/(1.0 - xb%q(i,j,k))
80
81 sdmd=sdmd+mr_a(k)*xb%dnw(k)
82 s1md=s1md+(1.0+mr_b(k))*xb%dnw(k)
83 end do
84
85 mu=-(xa%psfc(i,j)+xb%psac(i,j)*sdmd)/s1md
86
87 p(xb%kte+1)=0.0
88
89 do k=xb%kte,xb%kts,-1
90 p(k)=p(k+1)-(mu*(1.0+mr_b(k)) &
91 + xb%psac(i,j)*mr_a(k))*xb%dnw(k)
92
93 xa%p(i,j,k)=0.5*(p(k)+p(k+1))
94 end do
95 else if (fg_format == fg_format_kma_global) then
96 do k=xb%kts,xb%kte
97 if (k == xb%kte) then
98 coeff=xb%KMA_B(K)/(xb%KMA_A(K)+xb%KMA_B(K)*xb%psfc(I,J)/100.0)
99 else
100 PU = xb%KMA_A(K+1) + xb%KMA_B(K+1)*xb%psfc(I,J)/100.0
101 PD = xb%KMA_A(K ) + xb%KMA_B(K )*xb%psfc(I,J)/100.0
102 coeff=xb%KMA_B(K) *1.0/(PD-PU)**2*(-PU*(LOG(PD)-LOG(PU)) &
103 + PD-PU)&
104 + xb%KMA_B(K+1)*1.0/(PD-PU)**2*(PD*(LOG(PD)-LOG(PU))-PD+PU)
105 end if
106 ! Here since xa%psfc holds value in Pa. dlnp -> dp
107 xa%p(i,j,k) = xb%p(i,j,k) * xa%psfc(I,J)/100. * coeff
108
109 end do
110 end if
111 end do
112 end do
113
114 call da_pt_to_rho_lin(xb, xa, xp)
115
116 ! If Testing_wrfvar = .true., not "XToY" transform needed to do here:
117
118 if (.not.testing_wrfvar) then
119 ! Exchange XA halo region.
120 call wrf_dm_halo(xp%domdesc,xp%comms,xp%halo_id4)
121
122 if (sfc_assi_options == 2) then
123 call da_transform_xtowtq (xp, xb, xa)
124 ! Exchange XA (surface variable) halo region.
125 call wrf_dm_halo(xp%domdesc,xp%comms,xp%halo_id6)
126 end if
127
128 if (use_ssmt1obs .or. use_ssmt2obs .or. Use_GpspwObs .or. &
129 Use_SsmiTbObs .or. Use_SsmiRetrievalObs .or. Use_GpsrefObs) then
130
131 ! Now do something for PW
132 call da_transform_xtotpw(xa, xb)
133
134 ! GPS Refractivity:
135 if (use_GpsrefObs) &
136 call da_transform_xtogpsref_lin(xa, xb, xp)
137
138 if (use_ssmt1obs .or. use_ssmt2obs .or. &
139 Use_SsmiTbObs .or. Use_SsmiRetrievalObs) then
140 if (global) then
141 call da_error(__FILE__,__LINE__, &
142 (/"xb%speed is not available, see da_transfer_kmatoxb.inc"/))
143 end if
144 call da_transform_xtoseasfcwind_lin(xa, xb)
145 end if
146 if (Use_SsmiTbObs) call da_transform_xtotb_lin (xa, xb)
147
148 ! Exchange XA halo region.
149 call wrf_dm_halo(xp%domdesc,xp%comms,xp%halo_id8)
150 end if
151 end if
152
153 ! Compute w increments using Richardson's eqn.
154
155 if (Use_RadarObs) then
156 call da_uvprho_to_w_lin(xb, xa, xp)
157
158 do k=xp%kts,xp%kte
159 do j=xp%jts,xp%jte
160 do i=xp%its,xp%ite
161 xa%wh(i,j,k)=0.5*(xa%w(i,j,k)+xa%w(i,j,k+1))
162 end do
163 end do
164 end do
165
166 call wrf_dm_halo(xp%domdesc,xp%comms,xp%halo_id13)
167
168 if (use_radarobs .and. Use_Radar_rf) then
169 ! Partition of hydrometeor increments via warm rain process
170
171 call da_moist_phys_lin( xb, xa, xp)
172
173 endif
174 end if
175
176 !---------------------------------------------------------------
177 ! Polar treatment for Global
178 !---------------------------------------------------------------
179
180 if (global) then
181 call da_get_vpoles(xa%u,xa%v, &
182 ids, ide, jds, jde, &
183 ims, ime, jms, jme, kms, kme, &
184 its, ite, jts, jte, kts, kte )
185 call da_get_spoles(xa%t, &
186 ids, ide, jds, jde, &
187 ims, ime, jms, jme, kms, kme, &
188 its, ite, jts, jte, kts, kte )
189 call da_get_spoles(xa%p, &
190 ids, ide, jds, jde, &
191 ims, ime, jms, jme, kms, kme, &
192 its, ite, jts, jte, kts, kte )
193 call da_get_spoles(xa%q, &
194 ids, ide, jds, jde, &
195 ims, ime, jms, jme, kms, kme, &
196 its, ite, jts, jte, kts, kte )
197 call da_get_spoles(xa%psfc, &
198 ids, ide, jds, jde, &
199 ims, ime, jms, jme, 1, 1, &
200 its, ite, jts, jte, 1, 1 )
201 call da_set_boundary_xa(xa, xb)
202 end if
203
204 if (trace_use) call da_trace_exit("da_transform_vtox")
205
206 end subroutine da_transform_vtox
207
208