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