da_transform_xtoy_bogus.inc
References to this file elsewhere.
1 subroutine da_transform_xtoy_bogus (grid, iv, y)
2
3 !---------------------------------------------------------------------
4 ! Purpose: the linearized bogus observation operators.
5 !---------------------------------------------------------------------
6
7 implicit none
8
9 type (domain), intent(in) :: grid
10 type (iv_type), intent(in) :: iv ! Innovation vector (O-B).
11 type (y_type), intent(inout) :: y ! y = h (grid%xa) (linear)
12
13 integer :: n ! Loop counter.
14 integer :: i, j, k ! Index dimension.
15 integer :: num_levs ! Number of obs levels.
16 real :: dx, dxm
17 real :: dy, dym
18 real :: model_t(kts:kte)
19 real :: model_t9(kts:kte)
20 real :: model_q(kts:kte)
21 real :: model_q9(kts:kte)
22 real :: model_p_c(kts:kte)
23 real :: model_p_c9(kts:kte)
24 real :: psfc_loc,terr_loc,psfc_loc9
25
26 real, allocatable :: u(:,:)
27 real, allocatable :: v(:,:)
28 real, allocatable :: t(:,:)
29 real, allocatable :: q(:,:)
30
31 if (trace_use_dull) call da_trace_entry("da_transform_xtoy_bogus")
32
33 do n = iv%info(bogus)%n1,iv%info(bogus)%n2
34 num_levs = iv%info(bogus)%levels(n)
35
36 if (num_levs < 1) cycle
37
38 ! [1.1] Get cross pt. horizontal interpolation weights:
39
40 i = iv%info(bogus)%i(1,n)
41 dy = iv%info(bogus)%dy(1,n)
42 dym = iv%info(bogus)%dym(1,n)
43 j = iv%info(bogus)%j(1,n)
44 dx = iv%info(bogus)%dx(1,n)
45 dxm = iv%info(bogus)%dxm(1,n)
46
47 ! [1.2] Interpolate horizontally from cross points:
48
49 do k = kts, kte
50 model_t9(k) = dym*(dxm*grid%xa%t(i,j,k) + dx*grid%xa%t(i+1,j,k)) &
51 + dy *(dxm*grid%xa%t(i,j+1,k) + dx*grid%xa%t(i+1,j+1,k))
52 model_t(k) = dym*(dxm*grid%xb%t(i,j,k) + dx*grid%xb%t(i+1,j,k)) &
53 + dy *(dxm*grid%xb%t(i,j+1,k) + dx*grid%xb%t(i+1,j+1,k))
54 model_q9(k) = dym*(dxm*grid%xa%q(i,j,k) + dx*grid%xa%q(i+1,j,k)) &
55 + dy *(dxm*grid%xa%q(i,j+1,k) + dx*grid%xa%q(i+1,j+1,k))
56 model_q(k) = dym*(dxm*grid%xb%q(i,j,k) + dx*grid%xb%q(i+1,j,k)) &
57 + dy *(dxm*grid%xb%q(i,j+1,k) + dx*grid%xb%q(i+1,j+1,k))
58 model_p_c9(k) = dym*(dxm*grid%xa%p(i,j,k) + dx*grid%xa%p(i+1,j,k)) &
59 + dy *(dxm*grid%xa%p(i,j+1,k) + dx*grid%xa%p(i+1,j+1,k))
60 model_p_c(k) = dym*(dxm*grid%xb%p(i,j,k) + dx*grid%xb%p(i+1,j,k)) &
61 + dy *(dxm*grid%xb%p(i,j+1,k) + dx*grid%xb%p(i+1,j+1,k))
62 end do
63
64 terr_loc = dym*(dxm*grid%xb%terr(i,j) + dx*grid%xb%terr(i+1,j)) &
65 + dy *(dxm*grid%xb%terr(i,j+1) + dx*grid%xb%terr(i+1,j+1))
66 psfc_loc = dym*(dxm*grid%xb%psfc(i,j) + dx*grid%xb%psfc(i+1,j)) &
67 + dy *(dxm*grid%xb%psfc(i,j+1) + dx*grid%xb%psfc(i+1,j+1))
68 psfc_loc9 = dym*(dxm*grid%xa%psfc(i,j) + dx*grid%xa%psfc(i+1,j)) &
69 + dy *(dxm*grid%xa%psfc(i,j+1) + dx*grid%xa%psfc(i+1,j+1))
70
71 call da_tpq_to_slp_lin (model_t, model_q, model_p_c, terr_loc, psfc_loc, &
72 model_t9, model_q9, model_p_c9, psfc_loc9, y%bogus(n)%slp)
73 end do
74
75 allocate (u(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
76 allocate (v(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
77 allocate (t(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
78 allocate (q(iv%info(bogus)%max_lev,iv%info(bogus)%n1:iv%info(bogus)%n2))
79 u(:,:)=0.0
80 v(:,:)=0.0
81 t(:,:)=0.0
82 q(:,:)=0.0
83
84 call da_interp_lin_3d (grid%xa%u, iv%info(bogus), u)
85 call da_interp_lin_3d (grid%xa%v, iv%info(bogus), v)
86 call da_interp_lin_3d (grid%xa%t, iv%info(bogus), t)
87 call da_interp_lin_3d (grid%xa%q, iv%info(bogus), q)
88
89 do n=iv%info(bogus)%n1,iv%info(bogus)%n2
90 y%bogus(n)%u(:) = u(1:size(y%bogus(n)%u),n)
91 y%bogus(n)%v(:) = v(1:size(y%bogus(n)%v),n)
92 y%bogus(n)%t(:) = t(1:size(y%bogus(n)%t),n)
93 y%bogus(n)%q(:) = q(1:size(y%bogus(n)%q),n)
94 end do
95
96 deallocate (u)
97 deallocate (v)
98 deallocate (t)
99 deallocate (q)
100
101 if (trace_use_dull) call da_trace_exit("da_transform_xtoy_bogus")
102
103 end subroutine da_transform_xtoy_bogus
104
105