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_qrn ! Model qrn at ob location.
30 real, dimension(max_ob_levels) :: model_qrnb! Model qrn at ob location.
31 real :: model_ps
32
33 real :: xr,yr,zr
34
35 ! W_HALF is vertical velocity at half-sigma levels.
36
37 if (iv%num_Radar > 0) then
38 do n=iv%ob_numb(iv%current_ob_time-1)%radar + 1, &
39 iv%ob_numb(iv%current_ob_time)%radar
40 do k = 1, max_ob_levels
41 model_u(k) = 0.
42 model_v(k) = 0.
43 model_w(k) = 0.
44 model_qrn(k) = 0.
45 end do
46
47 num_levs = iv % Radar(n) % info % levels
48
49 i = iv%Radar(n)%loc%i
50 j = iv%Radar(n)%loc%j
51 dx = iv%Radar(n)%loc%dx
52 dy = iv%Radar(n)%loc%dy
53 dxm = iv%Radar(n)%loc%dxm
54 dym = iv%Radar(n)%loc%dym
55
56 model_ps = iv%Radar(n)%model_ps
57
58 ! [1.7] Calculate rv and rf at OBS location
59
60 xr = xb%ds * (iv%Radar(n)%loc%x - iv%Radar(n)%stn_loc%x)
61 yr = xb%ds * (iv%Radar(n)%loc%y - iv%Radar(n)%stn_loc%y)
62
63 do k = 1, num_levs
64
65 model_qrnb(k) = iv%Radar(n)%model_qrn(k)
66 model_p(k) = iv%Radar(n)%model_p(k)
67
68 if (iv % Radar(n) % height_qc(k) /= below_model_surface .and. &
69 iv % Radar(n) % height_qc(k) /= above_model_lid) then
70
71 if (use_Radar_rf) then
72 if (iv % Radar(n) % rf(k) % qc >= obs_qc_pointer) then
73 call da_reflectivity_adj(jo_grad_y%Radar(n)%rf(k), &
74 model_qrn(k), model_qrnb(k))
75 end if
76 end if
77
78 if (use_Radar_rv) then
79 if (iv % Radar(n) % rv(k) % qc >= obs_qc_pointer) then
80 zr=iv%Radar(n)%height(k) - iv%Radar(n)%stn_loc%elv
81
82 call da_radial_velocity_adj(jo_grad_y%Radar(n)%rv(k), &
83 model_p(k), model_u(k), model_v(k), model_w(k), &
84 model_qrn(k), model_ps, xr, yr, zr, model_qrnb(k))
85 end if
86 end if
87 end if
88
89 jo_grad_y%Radar(n)%rv(k) = 0.0
90 jo_grad_y%Radar(n)%rf(k) = 0.0
91
92 end do
93
94 ! [1.6] Interpolate horizontally from crs points:
95
96 call da_interp_lin_3d_adj(jo_grad_x % wh, xp, &
97 i, j, dx, dy, dxm, dym, &
98 model_w, num_levs, iv%Radar(n)%zk, &
99 num_levs)
100 call da_interp_lin_3d_adj(jo_grad_x % qrn, xp, &
101 i, j, dx, dy, dxm, dym, &
102 model_qrn, num_levs, iv%Radar(n)%zk, &
103 num_levs)
104 call da_interp_lin_3d_adj(jo_grad_x % v, xp, i, j, dx, dy, &
105 dxm, dym, model_v, num_levs, iv%Radar(n)%zk, num_levs)
106 call da_interp_lin_3d_adj(jo_grad_x % u, xp, i, j, dx, dy, &
107 dxm, dym, model_u, num_levs, iv%Radar(n)%zk, num_levs)
108 end do
109 end if
110
111 end subroutine da_transform_xtoy_radar_adj
112
113