da_transform_xtoy_radar_adj.inc
References to this file elsewhere.
1 subroutine da_transform_xtoy_radar_adj(xb, iv, xp, jo_grad_y, jo_grad_x)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 !------------------------------------------------------------------------
8 ! This subroutine is the adjoint of Doppler radar observation operators.
9 !------------------------------------------------------------------------
10
11 implicit none
12
13 type (xb_type), intent(in) :: xb ! first guess state.
14 type (ob_type), intent(in) :: iv ! obs. inc vector (o-b).
15 type (xpose_type), intent(in) :: xp ! Domain decomposition vars.
16 type (y_type) , intent(inout) :: jo_grad_y ! grad_y(jo)
17 type (x_type) , intent(inout) :: jo_grad_x ! grad_x(jo)
18
19 integer :: n ! Loop counter.
20 integer :: i, j, k ! Index dimension.
21 integer :: num_levs
22 real :: dx, dxm !
23 real :: dy, dym !
24
25 real, dimension(max_ob_levels) :: model_p ! Model value p at ob location.
26 real, dimension(max_ob_levels) :: model_u ! Model value u at ob location.
27 real, dimension(max_ob_levels) :: model_v ! Model value v at ob location.
28 real, dimension(max_ob_levels) :: model_w ! Model value w at ob location.
29 real, dimension(max_ob_levels) :: model_rho ! Model rho at ob location.
30 real, dimension(max_ob_levels) :: model_qrn ! Model qrn at ob location.
31 real, dimension(max_ob_levels) :: model_qrnb! Model qrn at ob location.
32 real :: model_ps
33
34 real :: xr,yr,zr
35
36 ! W_HALF is vertical velocity at half-sigma levels.
37
38 if (iv%num_Radar > 0) then
39 do n=iv%ob_numb(iv%current_ob_time-1)%radar + 1, &
40 iv%ob_numb(iv%current_ob_time)%radar
41 do k = 1, max_ob_levels
42 model_u(k) = 0.
43 model_v(k) = 0.
44 model_w(k) = 0.
45 model_qrn(k) = 0.
46 end do
47
48 num_levs = iv % Radar(n) % info % levels
49
50 i = iv%Radar(n)%loc%i
51 j = iv%Radar(n)%loc%j
52 dx = iv%Radar(n)%loc%dx
53 dy = iv%Radar(n)%loc%dy
54 dxm = iv%Radar(n)%loc%dxm
55 dym = iv%Radar(n)%loc%dym
56
57 model_ps = iv%Radar(n)%model_ps
58
59 ! [1.7] Calculate rv and rf at OBS location
60
61 xr = xb%ds * (iv%Radar(n)%loc%x - iv%Radar(n)%stn_loc%x)
62 yr = xb%ds * (iv%Radar(n)%loc%y - iv%Radar(n)%stn_loc%y)
63
64 do k = 1, num_levs
65
66 model_qrnb(k) = iv%Radar(n)%model_qrn(k)
67 model_p(k) = iv%Radar(n)%model_p(k)
68
69 if (iv % Radar(n) % height_qc(k) /= below_model_surface .and. &
70 iv % Radar(n) % height_qc(k) /= above_model_lid) then
71
72 if (use_Radar_rf) then
73 if (iv % Radar(n) % rf(k) % qc >= obs_qc_pointer) then
74 call da_reflectivity_adj(jo_grad_y%Radar(n)%rf(k), &
75 model_qrn(k), model_qrnb(k))
76 end if
77 end if
78
79 if (use_Radar_rv) then
80 if (iv % Radar(n) % rv(k) % qc >= obs_qc_pointer) then
81 zr=iv%Radar(n)%height(k) - iv%Radar(n)%stn_loc%elv
82
83 call da_radial_velocity_adj(jo_grad_y%Radar(n)%rv(k), &
84 model_p(k), model_u(k), model_v(k), model_w(k), &
85 model_qrn(k), model_ps, xr, yr, zr, model_qrnb(k))
86 end if
87 end if
88 end if
89
90 jo_grad_y%Radar(n)%rv(k) = 0.0
91 jo_grad_y%Radar(n)%rf(k) = 0.0
92
93 end do
94
95 ! [1.6] Interpolate horizontally from crs points:
96
97 call da_interp_lin_3d_adj(jo_grad_x % wh, xp, &
98 i, j, dx, dy, dxm, dym, &
99 model_w, num_levs, iv%Radar(n)%zk, &
100 num_levs)
101 call da_interp_lin_3d_adj(jo_grad_x % qrn, xp, &
102 i, j, dx, dy, dxm, dym, &
103 model_qrn, num_levs, iv%Radar(n)%zk, &
104 num_levs)
105 call da_interp_lin_3d_adj(jo_grad_x % v, xp, i, j, dx, dy, &
106 dxm, dym, model_v, num_levs, iv%Radar(n)%zk, num_levs)
107 call da_interp_lin_3d_adj(jo_grad_x % u, xp, i, j, dx, dy, &
108 dxm, dym, model_u, num_levs, iv%Radar(n)%zk, num_levs)
109 end do
110 end if
111
112 end subroutine da_transform_xtoy_radar_adj
113
114