da_transform_xtoy_bogus.inc
References to this file elsewhere.
1 subroutine da_transform_xtoy_bogus (xa, xb, iv, xp, y)
2
3 !---------------------------------------------------------------------
4 ! Purpose: the linearized bogus observation operators.
5 !---------------------------------------------------------------------
6
7 implicit none
8
9 type (xb_type), intent(in) :: xb ! first guess state.
10 type (x_type), intent(in) :: xa ! gridded analysis increment.
11 type (ob_type), intent(in) :: iv ! Innovation vector (O-B).
12 type (xpose_type), intent(in) :: xp ! Domain decomposition vars.
13 type (y_type), intent(inout) :: y ! y = h (xa) (linear)
14
15 integer :: n ! Loop counter.
16 integer :: i, j, k ! Index dimension.
17 integer :: num_levs ! Number of obs levels.
18 real :: dx, dxm !
19 real :: dy, dym !
20 real, dimension(xp%kts:xp%kte) :: model_t,model_t9
21 real, dimension(xp%kts:xp%kte) :: model_q,model_q9
22 real, dimension(xp%kts:xp%kte) :: model_p_c,model_p_c9
23 real :: psfc_loc,terr_loc,psfc_loc9
24
25 if (iv%num_bogus > 0) then
26 do n = iv%ob_numb(iv%current_ob_time-1)%bogus+1, iv%ob_numb(iv%current_ob_time)%bogus
27
28 !xyh y%bogus(n)%slp = 0.0
29 !xyh y%bogus(n)%u(:) = 0.0
30 !xyh y%bogus(n)%v(:) = 0.0
31 !xyh y%bogus(n)%t(:) = 0.0
32 !xyh y%bogus(n)%q(:) = 0.0
33
34 num_levs = iv % bogus(n) % info % levels
35
36 if (num_levs < 1) cycle
37
38 ! [1.1] Get cross pt. horizontal interpolation weights:
39
40 i = iv%bogus(n)%loc%i
41 dy = iv%bogus(n)%loc%dy
42 dym = iv%bogus(n)%loc%dym
43 j = iv%bogus(n)%loc%j
44 dx = iv%bogus(n)%loc%dx
45 dxm = iv%bogus(n)%loc%dxm
46
47 ! [1.2] Interpolate horizontally from cross points:
48
49 do k = xp%kts, xp%kte
50 model_t9(k) = dym*(dxm*xa%t(i,j,k) + dx*xa%t(i+1,j,k)) &
51 + dy *(dxm*xa%t(i,j+1,k) + dx*xa%t(i+1,j+1,k))
52 model_t(k) = dym*(dxm*xb%t(i,j,k) + dx*xb%t(i+1,j,k)) &
53 + dy *(dxm*xb%t(i,j+1,k) + dx*xb%t(i+1,j+1,k))
54 model_q9(k) = dym*(dxm*xa%q(i,j,k) + dx*xa%q(i+1,j,k)) &
55 + dy *(dxm*xa%q(i,j+1,k) + dx*xa%q(i+1,j+1,k))
56 model_q(k) = dym*(dxm*xb%q(i,j,k) + dx*xb%q(i+1,j,k)) &
57 + dy *(dxm*xb%q(i,j+1,k) + dx*xb%q(i+1,j+1,k))
58 model_p_c9(k) = dym*(dxm*xa%p(i,j,k) + dx*xa%p(i+1,j,k)) &
59 + dy *(dxm*xa%p(i,j+1,k) + dx*xa%p(i+1,j+1,k))
60 model_p_c(k) = dym*(dxm*xb%p(i,j,k) + dx*xb%p(i+1,j,k)) &
61 + dy *(dxm*xb%p(i,j+1,k) + dx*xb%p(i+1,j+1,k))
62 end do
63
64 terr_loc = dym*(dxm*xb%terr(i,j) + dx*xb%terr(i+1,j)) &
65 + dy *(dxm*xb%terr(i,j+1) + dx*xb%terr(i+1,j+1))
66 psfc_loc = dym*(dxm*xb%psfc(i,j) + dx*xb%psfc(i+1,j)) &
67 + dy *(dxm*xb%psfc(i,j+1) + dx*xb%psfc(i+1,j+1))
68 psfc_loc9 = dym*(dxm*xa%psfc(i,j) + dx*xa%psfc(i+1,j)) &
69 + dy *(dxm*xa%psfc(i,j+1) + dx*xa%psfc(i+1,j+1))
70
71 call da_tpq_to_slp_lin (model_t, model_q, model_p_c, &
72 terr_loc, psfc_loc, &
73 model_t9, model_q9, &
74 model_p_c9, psfc_loc9, y%bogus(n)%slp, xp)
75
76 ! [1.4] Interpolate horizontally from dot points:
77 call da_interp_lin_3d(xa % u, xp, i, j, dx, dy, dxm, dym, &
78 y%bogus(n)%u, num_levs, iv%bogus(n)%zk, num_levs)
79 call da_interp_lin_3d(xa % v, xp, i, j, dx, dy, dxm, dym, &
80 y%bogus(n)%v, num_levs, iv%bogus(n)%zk, num_levs)
81 call da_interp_lin_3d(xa % t, xp, i, j, dx, dy, dxm, dym, &
82 y%bogus(n)%t, num_levs, iv%bogus(n)%zk, num_levs)
83 call da_interp_lin_3d(xa % q, xp, i, j, dx, dy, dxm, dym, &
84 y%bogus(n)%q, num_levs, iv%bogus(n)%zk, num_levs)
85 end do
86 end if
87
88 end subroutine da_transform_xtoy_bogus
89
90