da_transform_vtox.inc
References to this file elsewhere.
1 subroutine da_transform_vtox(grid, cv_size, xbx, be, ep, cv, vv, vp)
2
3 !--------------------------------------------------------------------------
4 ! Purpose: Control variable transform x' = Uv.
5 !--------------------------------------------------------------------------
6
7 implicit none
8
9 type(domain), intent(inout) :: grid
10 integer, intent(in) :: cv_size ! Size of cv array.
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
18 integer :: i, j, k
19 real :: sdmd, s1md, mu
20 real :: p(kms:kme)
21 real :: mr_a(kms:kme)
22 real :: mr_b(kms:kme)
23 real :: PU, PD, coeff
24
25 if (trace_use) call da_trace_entry("da_transform_vtox")
26
27 call da_zero_x (grid%xa)
28
29 !----------------------------------------------------------------------
30 ! [1.0]: Perform vv = u_h cv transform:
31 !----------------------------------------------------------------------
32
33 if (global) then
34 call da_transform_vtovv_global(cv_size, xbx, be, cv, vv)
35 else
36 call da_transform_vtovv(grid, cv_size, be, cv, vv)
37 end if
38
39 !----------------------------------------------------------------------
40 ! [2.0]: Perform vp = u_v vv transform:
41 !----------------------------------------------------------------------
42
43 if (vert_corr == vert_corr_2) then
44 call da_vertical_transform('u', be, grid%xb % vertical_inner_product, vv, vp)
45 else
46 vp % v1(its:ite,jts:jte,kts:kte) = vv % v1(its:ite,jts:jte,kts:kte)
47 vp % v2(its:ite,jts:jte,kts:kte) = vv % v2(its:ite,jts:jte,kts:kte)
48 vp % v3(its:ite,jts:jte,kts:kte) = vv % v3(its:ite,jts:jte,kts:kte)
49 vp % v4(its:ite,jts:jte,kts:kte) = vv % v4(its:ite,jts:jte,kts:kte)
50 vp % v5(its:ite,jts:jte,kts:kte) = vv % v5(its:ite,jts:jte,kts:kte)
51 end if
52
53 ! WHY?
54 ! if (be % ne > 0) then
55 ! vp % alpha(its:ite,jts:jte,1:be%ne) = vv%alpha(its:ite,jts:jte,1:be%ne)
56 ! end if
57
58 !----------------------------------------------------------------------
59 ! [3.0]: Perform x = u_p vp transform::
60 !----------------------------------------------------------------------
61
62 call da_transform_vptox(grid, vp, be, ep)
63
64 !----------------------------------------------------------------------
65 ! [4.0]: Move the following:
66 !----------------------------------------------------------------------
67
68 do j=jts,jte
69 do i=its,ite
70 if (fg_format == fg_format_wrf) then
71 sdmd=0.0
72 s1md=0.0
73 do k=kts,kte
74 mr_a(k) = grid%xa%q(i,j,k)/(1.0 - grid%xb%q(i,j,k))**2
75 mr_b(k) = grid%xb%q(i,j,k)/(1.0 - grid%xb%q(i,j,k))
76
77 sdmd=sdmd+mr_a(k)*grid%xb%dnw(k)
78 s1md=s1md+(1.0+mr_b(k))*grid%xb%dnw(k)
79 end do
80
81 mu=-(grid%xa%psfc(i,j)+grid%xb%psac(i,j)*sdmd)/s1md
82
83 p(kte+1)=0.0
84
85 do k=kte,kts,-1
86 p(k)=p(k+1)-(mu*(1.0+mr_b(k)) &
87 + grid%xb%psac(i,j)*mr_a(k))*grid%xb%dnw(k)
88
89 grid%xa%p(i,j,k)=0.5*(p(k)+p(k+1))
90 end do
91 else if (fg_format == fg_format_kma_global) then
92 do k=kts,kte
93 if (k == kte) then
94 coeff=grid%xb%KMA_B(K)/(grid%xb%KMA_A(K)+grid%xb%KMA_B(K)*grid%xb%psfc(I,J)/100.0)
95 else
96 PU = grid%xb%KMA_A(K+1) + grid%xb%KMA_B(K+1)*grid%xb%psfc(I,J)/100.0
97 PD = grid%xb%KMA_A(K ) + grid%xb%KMA_B(K )*grid%xb%psfc(I,J)/100.0
98 coeff=grid%xb%KMA_B(K) *1.0/(PD-PU)**2*(-PU*(LOG(PD)-LOG(PU)) &
99 + PD-PU)&
100 + grid%xb%KMA_B(K+1)*1.0/(PD-PU)**2*(PD*(LOG(PD)-LOG(PU))-PD+PU)
101 end if
102 ! Here since grid%xa%psfc holds value in Pa. dlnp -> dp
103 grid%xa%p(i,j,k) = grid%xb%p(i,j,k) * grid%xa%psfc(I,J)/100.0 * coeff
104
105 end do
106 end if
107 end do
108 end do
109
110 call da_pt_to_rho_lin(grid)
111
112 ! If test_wrfvar = .true., not "XToY" transform needed to do here:
113
114 if (.not.test_wrfvar) then
115 ! Exchange grid%xa halo region.
116 #ifdef DM_PARALLEL
117 #include "HALO_XA.inc"
118 #endif
119
120 if (sfc_assi_options == 2) then
121 call da_transform_xtowtq (grid)
122 ! Exchange grid%xa (surface variable) halo region.
123 #ifdef DM_PARALLEL
124 #include "HALO_SFC_XA.inc"
125 #endif
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(grid)
133
134 ! GPS Refractivity:
135 if (use_gpsrefobs) &
136 call da_transform_xtogpsref_lin(grid)
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 (/"grid%xb%speed is not available, see da_transfer_kmatoxb.inc"/))
143 end if
144 call da_transform_xtoseasfcwind_lin(grid)
145 end if
146 if (use_ssmitbobs) call da_transform_xtotb_lin (grid)
147
148 ! Exchange grid%xa halo region.
149 #ifdef DM_PARALLEL
150 #include "HALO_SSMI_XA.inc"
151 #endif
152 end if
153 end if
154
155 ! Compute w increments using Richardson's eqn.
156
157 if (use_radarobs) then
158 call da_uvprho_to_w_lin(grid)
159
160 do k=kts,kte
161 do j=jts,jte
162 do i=its,ite
163 grid%xa%wh(i,j,k)=0.5*(grid%xa%w(i,j,k)+grid%xa%w(i,j,k+1))
164 end do
165 end do
166 end do
167
168 #ifdef DM_PARALLEL
169 #include "HALO_RADAR_XA_W.inc"
170 #endif
171
172 if (use_radarobs .and. use_radar_rf) then
173 ! Partition of hydrometeor increments via warm rain process
174 call da_moist_phys_lin(grid)
175 end if
176 end if
177
178 !---------------------------------------------------------------
179 ! Polar treatment for Global
180 !---------------------------------------------------------------
181
182 if (global) then
183 call da_get_vpoles(grid%xa%u,grid%xa%v, &
184 ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
185 call da_get_spoles(grid%xa%t, &
186 ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
187 call da_get_spoles(grid%xa%p, &
188 ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
189 call da_get_spoles(grid%xa%q, &
190 ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
191 call da_get_spoles(grid%xa%psfc, &
192 ids, ide, jds, jde, ims, ime, jms, jme, 1, 1, its, ite, jts, jte, 1, 1)
193 call da_set_boundary_xa(grid)
194 end if
195
196 if (trace_use) call da_trace_exit("da_transform_vtox")
197
198 end subroutine da_transform_vtox
199
200