da_transform_xtoy_bogus_adj.inc
References to this file elsewhere.
1 subroutine da_transform_xtoy_bogus_adj(xb, iv, xp, jo_grad_y, jo_grad_x)
2
3 !---------------------------------------------------------------------
4 ! Purpose: the adjoint of bogus observation operators.
5 !---------------------------------------------------------------------
6
7 implicit none
8
9 type (xb_type), intent(in) :: xb ! first guess state.
10 type (ob_type), intent(in) :: iv ! obs. inc vector (o-b).
11 type (xpose_type), intent(in) :: xp ! Domain decomposition vars.
12 type (y_type) , intent(inout) :: jo_grad_y ! grad_y(jo)
13 type (x_type) , intent(inout) :: jo_grad_x ! grad_x(jo)
14
15 integer :: n ! Loop counter.
16 integer :: i, j, k ! Index dimension.
17 integer :: num_levs
18 real :: dx, dxm !
19 real :: dy, dym !
20
21 real, dimension(xp%kms:xp%kme) :: model_t,tt
22 real, dimension(xp%kms:xp%kme) :: model_q,qq
23 real, dimension(xp%kms:xp%kme) :: model_p_c,pp
24 real :: psfc_loc,terr_loc,ppsfc
25
26 if (iv%num_bogus > 0) then
27 do n= iv%ob_numb(iv%current_ob_time-1)%bogus+1, iv%ob_numb(iv%current_ob_time)%bogus
28
29 num_levs = iv % bogus(n) % info % levels
30
31 if (num_levs < 1) cycle
32
33 ! [1.1] Get cross pt. horizontal interpolation weights:
34
35 i = iv%bogus(n)%loc%i
36 dy = iv%bogus(n)%loc%dy
37 dym = iv%bogus(n)%loc%dym
38 j = iv%bogus(n)%loc%j
39 dx = iv%bogus(n)%loc%dx
40 dxm = iv%bogus(n)%loc%dxm
41
42 ! [1.1] Adjoint feedback from Y to X for u and v:
43
44 call da_interp_lin_3d_adj(jo_grad_x % u, xp, i, j, dx, dy, dxm, dym, &
45 jo_grad_y%bogus(n)%u, num_levs, &
46 iv%bogus(n)%zk, num_levs)
47 call da_interp_lin_3d_adj(jo_grad_x % v, xp, i, j, dx, dy, dxm, dym, &
48 jo_grad_y%bogus(n)%v, num_levs, &
49 iv%bogus(n)%zk, num_levs)
50 call da_interp_lin_3d_adj(jo_grad_x % t, xp, i, j, dx, dy, dxm, dym, &
51 jo_grad_y%bogus(n)%t, num_levs, &
52 iv%bogus(n)%zk, num_levs)
53 call da_interp_lin_3d_adj(jo_grad_x % q, xp, i, j, dx, dy, dxm, dym, &
54 jo_grad_y%bogus(n)%q, num_levs, &
55 iv%bogus(n)%zk, num_levs)
56
57 ! [1.2] Compute the feedback from SLP to t and q:
58
59 ! 1.2.1 Background at the observation location:
60
61 do k = xp%kts, xp%kte
62 model_t(k) = dym*(dxm*xb%t(i,j,k) + dx*xb%t(i+1,j,k)) &
63 + dy *(dxm*xb%t(i,j+1,k) + dx*xb%t(i+1,j+1,k))
64 model_q(k) = dym*(dxm*xb%q(i,j,k) + dx*xb%q(i+1,j,k)) &
65 + dy *(dxm*xb%q(i,j+1,k) + dx*xb%q(i+1,j+1,k))
66 model_p_c(k) = dym*(dxm*xb%p(i,j,k) + dx*xb%p(i+1,j,k)) &
67 + dy *(dxm*xb%p(i,j+1,k) + dx*xb%p(i+1,j+1,k))
68 end do
69
70 terr_loc = dym*(dxm*xb%terr(i,j) + dx*xb%terr(i+1,j)) &
71 + dy *(dxm*xb%terr(i,j+1) + dx*xb%terr(i+1,j+1))
72 psfc_loc = dym*(dxm*xb%psfc(i,j) + dx*xb%psfc(i+1,j)) &
73 + dy *(dxm*xb%psfc(i,j+1) + dx*xb%psfc(i+1,j+1))
74
75 ! 1.2.2 Compute the feedback from SLP to p, t, q, and psfc
76 ! at the observed location:
77
78 call da_tpq_to_slp_adj(model_t, model_q ,model_p_c, &
79 terr_loc, psfc_loc, &
80 tt, qq, pp, ppsfc, &
81 jo_grad_y%bogus(n)%slp, xp)
82
83 ! 1.2.3 feedback from the observed location to grid space:
84
85 ! 1.2.3.1 for Psfc
86
87 jo_grad_x % psfc(i,j) = jo_grad_x % psfc(i,j) + dym*dxm*ppsfc
88 jo_grad_x % psfc(i+1,j) = jo_grad_x % psfc(i+1,j)+ dym*dx*ppsfc
89 jo_grad_x % psfc(i,j+1) = jo_grad_x % psfc(i,j+1)+ dy*dxm*ppsfc
90 jo_grad_x % psfc(i+1,j+1)=jo_grad_x % psfc(i+1,j+1)+dy*dx*ppsfc
91
92 ! 1.2.3.2 for t, q, p
93
94 do k = xp%kts, xp%kte
95 jo_grad_x % t(i,j,k) = jo_grad_x % t(i,j,k)+dym*dxm*tt(k)
96 jo_grad_x % t(i+1,j,k) = jo_grad_x % t(i+1,j,k)+ dym*dx*tt(k)
97 jo_grad_x % t(i,j+1,k) = jo_grad_x % t(i,j+1,k)+ dy*dxm*tt(k)
98 jo_grad_x % t(i+1,j+1,k)=jo_grad_x % t(i+1,j+1,k)+dy*dx*tt(k)
99 jo_grad_x % q(i,j,k) = jo_grad_x % q(i,j,k)+dym*dxm*qq(k)
100 jo_grad_x % q(i+1,j,k) = jo_grad_x % q(i+1,j,k)+ dym*dx*qq(k)
101 jo_grad_x % q(i,j+1,k) = jo_grad_x % q(i,j+1,k)+ dy*dxm*qq(k)
102 jo_grad_x % q(i+1,j+1,k)=jo_grad_x % q(i+1,j+1,k)+dy*dx*qq(k)
103 jo_grad_x % p(i,j,k) = jo_grad_x % p(i,j,k)+dym*dxm*pp(k)
104 jo_grad_x % p(i+1,j,k) = jo_grad_x % p(i+1,j,k)+ dym*dx*pp(k)
105 jo_grad_x % p(i,j+1,k) = jo_grad_x % p(i,j+1,k)+ dy*dxm*pp(k)
106 jo_grad_x % p(i+1,j+1,k)=jo_grad_x % p(i+1,j+1,k)+dy*dx*pp(k)
107 end do
108 end do
109 end if
110
111 end subroutine da_transform_xtoy_bogus_adj
112
113