da_transform_vtox_adj.inc
References to this file elsewhere.
1 subroutine da_transform_vtox_adj(cv_size, xb, xbx, be, ep, xa, xp, vp, vv, cv)
2
3 !--------------------------------------------------------------------------
4 ! Purpose: Control variable transform Adjoint
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 perturbation structure.
14 type (x_type), intent(inout) :: xa ! grad_x(jo)
15 type (xpose_type), intent(inout) :: xp ! Dimensions and xpose buffers.
16 type (vp_type),intent(out) :: vp ! grdipt/level cv (local).
17 type (vp_type),intent(out) :: vv ! grdipt/eof cv (local).
18 real, intent(inout) :: cv(1:cv_size) ! control variables.
19
20 integer :: i, j, k
21 real :: sdmd, s1md, mu
22 real, dimension(xb%kms:xb%kme) :: p, mr_a, mr_b
23 real :: PU, PD, coeff
24
25 if (trace_use) call da_trace_entry("da_transform_vtox_adj")
26
27 !-------------------------------------------------------------------------
28 ! Polar treatment for Global
29 !-------------------------------------------------------------------------
30
31 if (global) then
32 ! Poles treatment for global WRFVAR
33 call da_get_avpoles(xa%u,xa%v,&
34 ids, ide, jds, jde, &
35 ims, ime, jms, jme, kms, kme, &
36 its, ite, jts, jte, kts, kte )
37 call da_get_aspoles(xa%t, &
38 ids, ide, jds, jde, &
39 ims, ime, jms, jme, kms, kme, &
40 its, ite, jts, jte, kts, kte )
41 call da_get_aspoles(xa%p, &
42 ids, ide, jds, jde, &
43 ims, ime, jms, jme, kms, kme, &
44 its, ite, jts, jte, kts, kte )
45 call da_get_aspoles(xa%q, &
46 ids, ide, jds, jde, &
47 ims, ime, jms, jme, kms, kme, &
48 its, ite, jts, jte, kts, kte )
49 call da_get_aspoles(xa%psfc, &
50 ids, ide, jds, jde, &
51 ims, ime, jms, jme, 1, 1, &
52 its, ite, jts, jte, 1, 1 )
53 end if
54
55 ! Compute w increments using Richardson's eqn.
56 if (Use_RadarObs) then
57
58 do k=xp%kts,xp%kte
59 do j=xp%jts,xp%jte
60 do i=xp%its,xp%ite
61 xa%w(i,j,k)=xa%w(i,j,k)+0.5*xa%wh(i,j,k)
62 xa%w(i,j,k+1)=xa%w(i,j,k+1)+0.5*xa%wh(i,j,k)
63 xa%wh(i,j,k)=0.
64 end do
65 end do
66 end do
67
68 call da_uvprho_to_w_adj(xb, xa, xp)
69
70 if (use_radarobs .and. Use_Radar_rf) then
71
72 ! Partition of hydrometeor increments via warm rain process
73
74 call da_moist_phys_adj( xb, xa, xp)
75
76 endif
77
78 end if
79
80 !-------------------------------------------------------------------------
81 ! If Testing_wrfvar = .true., not "XToY" transform needed to do here (YRG):
82
83 if (.not.testing_wrfvar) then
84
85 if (use_ssmt1obs .or. use_ssmt2obs .or. Use_GpspwObs .or. &
86 Use_SsmiTbObs .or. Use_SsmiRetrievalObs .or. use_GpsrefObs) then
87
88 if (use_ssmiTbobs) call da_transform_xtotb_adj(xa, xb)
89
90 if (use_ssmt1obs .or. use_ssmt2obs .or. &
91 Use_SsmiTbObs .or. Use_SsmiRetrievalObs) then
92 if (global) then
93 call da_error(__FILE__,__LINE__, &
94 (/"xb%speed is not available, see da_transfer_kmatoxb.inc"/))
95 end if
96 call da_transform_xtoseasfcwind_adj(xa, xb)
97 end if
98
99 ! GPS Refractivity:
100 if (use_GpsrefObs) &
101 call da_transform_xtogpsref_adj(xa, xb, xp)
102
103 ! Now for PW.
104 call da_transform_xtotpw_adj(xa, xb)
105
106 end if
107
108 if (sfc_assi_options == 2) call da_transform_xtowtq_adj(xp, xb, xa)
109
110 end if
111 !-------------------------------------------------------------------------
112
113 call da_pt_to_rho_adj(xb, xa)
114 do j=xb%jts,xb%jte
115 do i=xb%its,xb%ite
116 if (fg_format == fg_format_wrf) then
117 mu=0.0
118 s1md=0.0
119
120 p(:)=0.0
121
122 do k=xb%kts,xb%kte
123 mr_b(k) = xb%q(i,j,k)/(1.0 - xb%q(i,j,k))
124 s1md=s1md+(1.0+mr_b(k))*xb%dnw(k)
125
126 p(k) = p(k) + 0.5*xa%p(i,j,k)
127 p(k+1) = p(k+1) + 0.5*xa%p(i,j,k)
128
129 mu = mu - p(k)*(1.0+mr_b(k))*xb%dnw(k)
130
131 mr_a(k) = - p(k)*xb%psac(i,j)*xb%dnw(k)
132 p(k+1) = p(k+1) + p(k)
133 end do
134
135 xa%psfc(i,j) = xa%psfc(i,j) - mu/s1md
136 sdmd=-mu*xb%psac(i,j)/s1md
137
138 do k=xb%kts,xb%kte
139 mr_a(k) = mr_a(k) + sdmd*xb%dnw(k)
140 xa%q(i,j,k) = xa%q(i,j,k) + mr_a(k)/(1.0 - xb%q(i,j,k))**2
141 end do
142 else if (fg_format == fg_format_kma_global)then
143 do k=xb%kts,xb%kte
144 if (k == xb%kte) then
145 coeff = xb%KMA_B(K)/(xb%KMA_A(K)+xb%KMA_B(K)* &
146 xb%psfc(I,J)/100.0)
147 else
148 PU = xb%KMA_A(K+1) + xb%KMA_B(K+1)*xb%psfc(I,J)/100.0
149 PD = xb%KMA_A(K ) + xb%KMA_B(K )*xb%psfc(I,J)/100.0
150 coeff=xb%KMA_B(K)*1.0/(PD-PU)**2*(-PU*(LOG(PD)-LOG(PU))+PD-PU)&
151 + xb%KMA_B(K+1)*1.0/(PD-PU)**2*(PD*(LOG(PD)-LOG(PU))-PD+PU)
152 end if
153
154 xa%psfc(i,j) = xa % psfc(i,j) + &
155 xb%p(i,j,k) * xa % p(i,j,k)/100.0 * coeff
156 end do
157 end if
158 end do
159 end do
160
161 if (global) then
162 call da_set_boundary_xa(xa, xb)
163 end if
164
165 !-------------------------------------------------------------------------
166 ! [3.0]: Perform x = u_p vp transform::
167 !-------------------------------------------------------------------------
168
169 call da_zero_vp_type (vp)
170 call da_transform_vptox_adj(xb, xa, vp, be, ep, xp)
171
172 !-------------------------------------------------------------------------
173 ! [2.0]: Perform vp = u_v vv transform:
174 !-------------------------------------------------------------------------
175
176 call da_zero_vp_type (vv)
177
178 ! Uv for alpha is a null transform:
179 ! if (be % ne > 0) then
180 ! vv % alpha(its:ite,jts:jte,1:be%ne) = &
181 ! vp % alpha(its:ite,jts:jte,1:be%ne)
182 ! end if
183
184 if (vert_corr == 2) then
185 call da_vertical_transform('u_adj', be, &
186 xb % vertical_inner_product, &
187 vv, vp)
188 else
189 vv % v1(its:ite,jts:jte,kts:kte) = vp % v1(its:ite,jts:jte,kts:kte)
190 vv % v2(its:ite,jts:jte,kts:kte) = vp % v2(its:ite,jts:jte,kts:kte)
191 vv % v3(its:ite,jts:jte,kts:kte) = vp % v3(its:ite,jts:jte,kts:kte)
192 vv % v4(its:ite,jts:jte,kts:kte) = vp % v4(its:ite,jts:jte,kts:kte)
193 vv % v5(its:ite,jts:jte,kts:kte) = vp % v5(its:ite,jts:jte,kts:kte)
194 end if
195
196 !-------------------------------------------------------------------------
197 ! [1.0]: perform vv = u_h cv transform:
198 !-------------------------------------------------------------------------
199
200 if (global) then
201 call da_transform_vtovv_global_adj(cv_size, xbx, be, cv, vv)
202 else
203 call da_transform_vtovv_adj(cv_size, xb, be, cv, vv, xp)
204 end if
205
206 if (trace_use) call da_trace_exit("da_transform_vtox_adj")
207
208 end subroutine da_transform_vtox_adj
209
210