da_transform_xtopsfc.inc
References to this file elsewhere.
1 subroutine da_transform_xtopsfc(xb, xa, xp, synop, y_synop)
2
3 !---------------------------------------------------------------------
4 ! Purpose: TBD
5 !---------------------------------------------------------------------
6
7 implicit none
8
9 type (x_type), intent(in) :: xa ! gridded analysis increment.
10 type (xb_type), intent(in) :: xb ! first guess state.
11 type (xpose_type), intent(in) :: xp ! domain decomposition vars.
12
13 type (synop_type), intent(in) :: synop
14 type (residual_synop_type), intent(out) :: y_synop
15
16 integer :: i, j ! index dimension.
17 real :: dx, dxm !
18 real :: dy, dym !
19
20 real :: hsm, to, qo
21 real :: tsm, qsm, psm
22 real :: psm_prime
23
24 ! 1.0 Get cross pt. horizontal interpolation weights:
25
26 i = synop%loc%i
27 dy = synop%loc%dy
28 dym = synop%loc%dym
29 j = synop%loc%j
30 dx = synop%loc%dx
31 dxm = synop%loc%dxm
32
33 ! 2.0 Surface assimilation approach 2 (2-m t and q, 10-m u and v)
34
35 call da_interp_lin_2d(xa%u10, xp%ims, xp%ime, xp%jms, xp%jme, &
36 i, j, dx, dy, dxm, dym, y_synop%u)
37 call da_interp_lin_2d(xa%v10, xp%ims, xp%ime, xp%jms, xp%jme, &
38 i, j, dx, dy, dxm, dym, y_synop%v)
39 call da_interp_lin_2d(xa%psfc, xp%ims, xp%ime, xp%jms, xp%jme, &
40 i, j, dx, dy, dxm, dym, psm_prime)
41 call da_interp_lin_2d(xa%t2, xp%ims, xp%ime, xp%jms, xp%jme, &
42 i, j, dx, dy, dxm, dym, y_synop%t)
43 call da_interp_lin_2d(xa%q2, xp%ims, xp%ime, xp%jms, xp%jme, &
44 i, j, dx, dy, dxm, dym, y_synop%q)
45 if (synop%p%qc < 0) then
46 y_synop%p = 0.0
47 return
48 end if
49
50 ! 3.0 The pressure at the observed height:
51
52 ! 3.1 Constants:
53
54 to = -888888.0
55 qo = -888888.0
56
57 ! Terrain height at the observed site:
58
59 ! 3.2 model background surface p, t, q, h at observed site:
60
61 call da_interp_lin_2d(xb%terr, xp%ims, xp%ime, xp%jms, xp%jme, &
62 i, j, dx, dy, dxm, dym, hsm)
63 call da_interp_lin_2d(xb%t2, xp%ims, xp%ime, xp%jms, xp%jme, &
64 i, j, dx, dy, dxm, dym, tsm)
65 call da_interp_lin_2d(xb%q2, xp%ims, xp%ime, xp%jms, xp%jme, &
66 i, j, dx, dy, dxm, dym, qsm)
67 call da_interp_lin_2d(xb%psfc, xp%ims, xp%ime, xp%jms, xp%jme, &
68 i, j, dx, dy, dxm, dym, psm)
69
70 ! 3.3 perturbations t, q, p at the model surface:
71
72 ! 4.0 Compute the perturbation of the surface pressure perturbation
73 ! at the observed height
74
75 if (synop%t%qc >= 0 .and. synop%q%qc >= 0) then
76 ! 4.1 Observed value = background + innovation: both t, q available:
77 ! ----------------------------------------
78
79 to = tsm + synop%t%inv
80 qo = qsm + synop%q%inv
81 call da_sfc_pre_lin(y_synop%p, psm_prime, y_synop%t, y_synop%q, &
82 psm, tsm, qsm, hsm, synop%h, to, qo)
83 else if (synop%t%qc >= 0 .and. synop%q%qc < 0) then
84
85 ! 4.2 Observed value = background + innovation: only t available
86 ! ----------------------------------------
87
88 to = tsm + synop%t%inv
89 call da_sfc_pre_lin(y_synop%p, psm_prime, y_synop%t, y_synop%q, &
90 psm, tsm, qsm, hsm, synop%h, to)
91 else
92 ! 4.3 No observed t and q available:
93 ! -----------------------------
94
95 call da_sfc_pre_lin(y_synop%p, psm_prime, y_synop%t, y_synop%q, &
96 psm, tsm, qsm, hsm, synop%h)
97 end if
98
99 end subroutine da_transform_xtopsfc
100
101