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