da_transform_xtopsfc_adj.inc
References to this file elsewhere.
1 subroutine da_transform_xtopsfc_adj(grid, iv, obs_index, synop, j_synop_y, jo_grad_x)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 type (domain), intent(in) :: grid
10 type (iv_type), intent(in) :: iv
11 integer, intent(in) :: obs_index
12 type (synop_type), intent(in) :: synop(:)
13 type (residual_synop_type), intent(inout) :: j_synop_y(:) ! grad_y(jo)
14 type (x_type), intent(inout) :: jo_grad_x ! grad_x(jo)
15
16
17 integer :: n
18
19 real :: to, qo
20 real, allocatable :: hsm(:,:)
21 real, allocatable :: tsm(:,:)
22 real, allocatable :: qsm(:,:)
23 real, allocatable :: psm(:,:)
24 real, allocatable :: psm_prime(:,:)
25 real, allocatable :: u(:,:)
26 real, allocatable :: v(:,:)
27 real, allocatable :: t(:,:)
28 real, allocatable :: q(:,:)
29
30 if (trace_use) call da_trace_entry("da_transform_xtopsfc_adj")
31
32 allocate (hsm(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
33 allocate (tsm(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
34 allocate (qsm(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
35 allocate (psm(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
36 allocate (psm_prime(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
37 allocate (u(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
38 allocate (v(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
39 allocate (t(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
40 allocate (q(1,iv%info(obs_index)%n1:iv%info(obs_index)%n2))
41
42 psm_prime = 0.0
43
44 ! 2.1 Terrain height at the observed site (xj, yi):
45
46 call da_interp_lin_2d(grid%xb%terr, iv%info(obs_index), 1, hsm)
47 call da_interp_lin_2d(grid%xb%t2, iv%info(obs_index), 1, tsm)
48 call da_interp_lin_2d(grid%xb%q2, iv%info(obs_index), 1, qsm)
49 call da_interp_lin_2d(grid%xb%psfc, iv%info(obs_index), 1, psm)
50
51 do n=iv%info(obs_index)%n1,iv%info(obs_index)%n2
52 if (synop(n)%p%qc >= 0) then
53 !------------------------------------------------------------------------
54 ! 2.0 Calculate gradient with respect the pressure at the observed height:
55 !------------------------------------------------------------------------
56
57 to = -888888.0
58 qo = -888888.0
59
60 ! 2.3 Zero out the adjoint variables:
61
62 !----------------------------------------------------------------
63 ! 3.0 Surface pressure gradient at the observed height
64 !----------------------------------------------------------------
65
66 ! 4.0 Get the surface pressure gradient at the model surface height (hsm)
67 ! 4.1 Both observed to and qo available:
68 if (synop(n)%t%qc >= 0 .and. synop(n)%q%qc >= 0) then
69 to = tsm(1,n) + synop(n)%t%inv
70 qo = qsm(1,n) + synop(n)%q%inv
71 call da_sfc_pre_adj (j_synop_y(n)%p, psm_prime(1,n), j_synop_y(n)%t, j_synop_y(n)%q, &
72 psm(1,n), tsm(1,n), qsm(1,n), hsm(1,n), synop(n)%h, to, qo)
73
74 ! 4.2 only observed to available:
75 else if (synop(n)%t%qc >= 0 .and. synop(n)%q%qc < 0) then
76 to = tsm(1,n) + synop(n)%t%inv
77 call da_sfc_pre_adj (j_synop_y(n)%p, psm_prime(1,n), j_synop_y(n)%t, j_synop_y(n)%q, &
78 psm(1,n), tsm(1,n), qsm(1,n), hsm(1,n), synop(n)%h, to)
79
80 ! 4.3 Both observed to and qo missing:
81 else
82 call da_sfc_pre_adj (j_synop_y(n)%p, psm_prime(1,n), j_synop_y(n)%t, j_synop_y(n)%q, &
83 psm(1,n), tsm(1,n), qsm(1,n), hsm(1,n), synop(n)%h)
84 end if
85 end if
86 t(1,n)=j_synop_y(n)%t
87 q(1,n)=j_synop_y(n)%q
88 u(1,n)=j_synop_y(n)%u
89 v(1,n)=j_synop_y(n)%v
90 end do
91
92 ! 2.2 convert the jo_grad_y to the model grids (t2, q2, u10, v10, and psfc)
93 call da_interp_lin_2d_adj(jo_grad_x%t2, iv%info(obs_index), 1, t)
94 call da_interp_lin_2d_adj(jo_grad_x%q2, iv%info(obs_index), 1, q)
95 call da_interp_lin_2d_adj(jo_grad_x%u10, iv%info(obs_index), 1, u)
96 call da_interp_lin_2d_adj(jo_grad_x%v10, iv%info(obs_index), 1, v)
97 call da_interp_lin_2d_adj(jo_grad_x%psfc, iv%info(obs_index), 1, psm_prime)
98
99 deallocate (hsm)
100 deallocate (tsm)
101 deallocate (qsm)
102 deallocate (psm)
103 deallocate (psm_prime)
104 deallocate (u)
105 deallocate (v)
106 deallocate (t)
107 deallocate (q)
108
109 if (trace_use) call da_trace_exit("da_transform_xtopsfc_adj")
110
111 end subroutine da_transform_xtopsfc_adj
112
113