da_transform_xtotb_lin.inc
References to this file elsewhere.
1 subroutine da_transform_xtotb_lin(xa, xb)
2
3 !-------------------------------------------------------------------------
4 ! Purpose: TBD
5 !-------------------------------------------------------------------------
6
7 implicit none
8
9 type (xb_type), intent(in) :: xb ! first guess state.
10 type (x_type), intent(inout):: xa ! gridded analysis increment.
11
12 integer :: ks,ke ! Range of 3rd dimension of arrays.
13 integer :: i,j,k
14 real :: dum1, dum2, zrhom, TGL_zrhom
15
16 real :: psfc,ta,gamma,sst,htpw,speed,alw,zcld,tpw
17 real :: TGL_psfc,TGL_ta,TGL_gamma,TGL_sst,TGL_tpw
18 real :: TGL_htpw,TGL_speed,TGL_alw,TGL_zcld
19
20 ks = xb%kts
21 ke = xb%kte
22
23 ! WHY
24 ! call da_transform_xtotpw(xa, xb)
25
26 do j=xb%jts,xb%jte
27 do i=xb%its,xb%ite
28 ! surface pressure (mb) (940 -1030)
29
30 psfc = 0.01*xb%psfc(i,j)
31 TGL_psfc = 0.01*xa%psfc(i,j)
32
33 ! sea surface temperature (k) (273 - 303) (doesnot change)
34
35 sst = xb%tgrn(i,j)
36 ! TGL_sst = xa%tgrn(1,1)
37 TGL_sst = 0.0
38
39 ! effective surface air temperature (263 - 303)
40
41 ta = xb%tgrn(i,j) + &
42 (xb%t(i,j,ks)-xb%tgrn(i,j))*log(2.0/0.0001)/ &
43 log((xb%h(i,j,ks)-xb%terr(i,j))/0.0001)
44
45 TGL_ta = (xa%t(i,j,ks)-0.)*log(2.0/0.0001)/ &
46 log((xb%h(i,j,ks)-xb%terr(i,j))/0.0001)
47
48 ! gamma is an emperical formula and zcld is given for now
49
50 gamma = (ta-270)*0.023 + 5.03 ! effective lapse rate (km) (4.0 - 6.5)
51 zcld = 1 ! effective cloud height (km)
52 ! = 1 if no cloud infomation
53 TGL_gamma = TGL_ta*0.023
54 TGL_zcld = 0.0
55
56
57 ! total precipitable water in (kg/m**2) (0 - 70)
58
59 tpw = xb%tpw(i,j)*10.0
60 TGL_tpw = xa%tpw(i,j)*10.0
61
62 ! speed, surface wind speed (m/sec) (0 - 30) , take 10 m wind
63
64 speed = xb%speed(i,j)
65 TGL_speed = xa%speed(i,j)
66
67 ! Column liquid water (kg/m**2) (0 - 0.5) (no data now. So, do it later.)
68
69 alw = 0.0
70 TGL_alw = 0.0
71
72 ! Column height weighted moisture density on the grid locally
73
74 zrhom = 0.0
75 do k=ks,ke
76 zrhom=zrhom+(xb%hf(i,j,k+1)-xb%hf(i,j,k))*xb%h(i,j,k)* &
77 xb%q(i,j,k)*xb%rho(i,j,k)
78 end do
79
80 TGL_zrhom = 0.0
81 do k = ks,ke
82 TGL_zrhom = (xb%hf(i,j,k+1)-xb%hf(i,j,k))*xb%h(i,j,k)* &
83 (xb%q(i,j,k)*xa%rho(i,j,k) + &
84 xa%q(i,j,k)*xb%rho(i,j,k)) + TGL_zrhom
85 end do
86
87 ! WHY
88 ! call da_transform_xtozrhoq(xb, i, j, zh, zf, zrhom)
89 ! call da_transform_xtozrhoq_lin(xb, xa, i, j, zh, zf, TGL_zrhom)
90
91 ! Column moisture density on the grid locally
92
93 htpw = zrhom/tpw/1000.
94 TGL_htpw = TGL_zrhom/tpw/1000. &
95 - TGL_tpw/tpw*htpw
96
97 dum1 = 0.0
98
99 call da_tb_tl(1,53.,psfc,ta,gamma,sst,tpw,htpw,speed,alw,zcld, &
100 ! xb%tb19v(i,j),xb%tb19h(i,j), &
101 TGL_psfc,TGL_ta,TGL_gamma,TGL_sst, &
102 TGL_tpw,TGL_htpw,TGL_speed,TGL_alw, &
103 TGL_zcld,xa%tb19v(i,j),xa%tb19h(i,j) )
104
105 call da_tb_tl(2,53.,psfc,ta,gamma,sst,tpw,htpw,speed,alw,zcld, &
106 ! xb%tb22v(i,j),dum1, &
107 TGL_psfc,TGL_ta,TGL_gamma,TGL_sst, &
108 TGL_tpw,TGL_htpw,TGL_speed,TGL_alw, &
109 TGL_zcld,xa%tb22v(i,j),dum2 )
110
111 call da_tb_tl(3,53.,psfc,ta,gamma,sst,tpw,htpw,speed,alw,zcld, &
112 ! xb%tb37v(i,j),xb%tb37h(i,j), &
113 TGL_psfc,TGL_ta,TGL_gamma,TGL_sst, &
114 TGL_tpw,TGL_htpw,TGL_speed,TGL_alw, &
115 TGL_zcld,xa%tb37v(i,j),xa%tb37h(i,j) )
116
117 call da_tb_tl(4,53.,psfc,ta,gamma,sst,tpw,htpw,speed,alw,zcld, &
118 ! xb%tb85v(i,j),xb%tb85h(i,j), &
119 TGL_psfc,TGL_ta,TGL_gamma,TGL_sst, &
120 TGL_tpw,TGL_htpw,TGL_speed,TGL_alw, &
121 TGL_zcld,xa%tb85v(i,j),xa%tb85h(i,j) )
122 end do
123 end do
124
125 end subroutine da_transform_xtotb_lin
126
127