da_transform_xtoy_radar_adj.inc
References to this file elsewhere.
1 subroutine da_transform_xtoy_radar_adj(grid, iv, 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 (domain), intent(in) :: grid
14 type (iv_type), intent(in) :: iv ! obs. inc vector (o-b).
15 type (y_type) , intent(inout) :: jo_grad_y ! grad_y(jo)
16 type (x_type) , intent(inout) :: jo_grad_x ! grad_x(jo)
17
18 integer :: k ! Index dimension.
19
20 integer :: n
21
22 real, allocatable :: model_p(:,:)
23 real, allocatable :: model_u(:,:)
24 real, allocatable :: model_v(:,:)
25 real, allocatable :: model_w(:,:)
26 real, allocatable :: model_qrn(:,:)
27 real, allocatable :: model_qrnb(:,:)
28 real, allocatable :: model_ps(:)
29
30 real :: xr,yr,zr
31
32 real :: alog10
33
34 if (trace_use) call da_trace_entry("da_transform_xtoy_radar_adj")
35
36 alog10= alog(10.0)
37
38 allocate (model_p(iv%info(radar)%max_lev,iv%info(radar)%n1:iv%info(radar)%n2))
39 allocate (model_u(iv%info(radar)%max_lev,iv%info(radar)%n1:iv%info(radar)%n2))
40 allocate (model_v(iv%info(radar)%max_lev,iv%info(radar)%n1:iv%info(radar)%n2))
41 allocate (model_w(iv%info(radar)%max_lev,iv%info(radar)%n1:iv%info(radar)%n2))
42 allocate (model_qrn(iv%info(radar)%max_lev,iv%info(radar)%n1:iv%info(radar)%n2))
43 allocate (model_qrnb(iv%info(radar)%max_lev,iv%info(radar)%n1:iv%info(radar)%n2))
44 allocate (model_ps(iv%info(radar)%n1:iv%info(radar)%n2))
45
46 ! Needed
47 model_u = 0.0
48 model_v = 0.0
49 model_w = 0.0
50 model_qrn = 0.0
51
52 ! W_HALF is vertical velocity at half-sigma levels.
53
54 model_ps(iv%info(radar)%n1:iv%info(radar)%n2) = iv%radar(iv%info(radar)%n1:iv%info(radar)%n2)%model_ps
55
56 do n=iv%info(radar)%n1,iv%info(radar)%n2
57
58 ! [1.7] Calculate rv and rf at OBS location
59
60 xr = grid%xb%ds * (iv%info(radar)%x(1,n) - iv%radar(n)%stn_loc%x)
61 yr = grid%xb%ds * (iv%info(radar)%y(1,n) - iv%radar(n)%stn_loc%y)
62
63 model_qrnb(1:iv%info(radar)%levels(n),n) = iv%radar(n)%model_qrn(:)
64 model_p (1:iv%info(radar)%levels(n),n) = iv%radar(n)%model_p(:)
65
66 do k = 1,iv%info(radar)%levels(n)
67 if (iv % radar(n) % height_qc(k) /= below_model_surface .and. &
68 iv % radar(n) % height_qc(k) /= above_model_lid) then
69
70 if (use_radar_rf) then
71 if (iv % radar(n) % rf(k) % qc >= obs_qc_pointer) then
72
73 model_qrn(k,n) = model_qrn(k,n) + leh2/(model_qrnb(k,n)*alog10) * jo_grad_y%radar(n)%rf(k)
74
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,n), model_u(k,n), model_v(k,n), model_w(k,n), &
84 model_qrn(k,n), model_ps(n), xr, yr, zr, model_qrnb(k,n))
85
86 end if
87 end if
88 end if
89 end do
90 jo_grad_y%radar(n)%rv(:) = 0.0
91 jo_grad_y%radar(n)%rf(:) = 0.0
92 end do ! n
93
94 ! [1.6] Interpolate horizontally from crs points:
95
96 call da_interp_lin_3d_adj (jo_grad_x % wh, iv%info(radar), model_w)
97 call da_interp_lin_3d_adj (jo_grad_x % qrn, iv%info(radar), model_qrn)
98 call da_interp_lin_3d_adj (jo_grad_x % v, iv%info(radar), model_v)
99 call da_interp_lin_3d_adj (jo_grad_x % u, iv%info(radar), model_u)
100
101 deallocate (model_p)
102 deallocate (model_u)
103 deallocate (model_v)
104 deallocate (model_w)
105 deallocate (model_qrn)
106 deallocate (model_qrnb)
107 deallocate (model_ps)
108
109 if (trace_use) call da_trace_exit("da_transform_xtoy_radar_adj")
110
111 end subroutine da_transform_xtoy_radar_adj
112
113