da_check_rh.inc
References to this file elsewhere.
1 subroutine da_check_rh(grid)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Remove RH over 100% or under 10%
5 ! Made the modification to those levels, which RH are less than 95%
6 !---------------------------------------------------------------------------
7
8 implicit none
9
10 type (domain), intent(inout) :: grid
11
12 integer :: imod(kts:kte)
13 real :: rhtol(kts:kte)
14 real :: x_qs(kts:kte)
15 real :: dz(kts:kte)
16
17 integer :: i, j, k
18 real :: tol_adjust_moist, tol_moist, oldrha, each_moist, es, weight
19 real :: upper_modify_rh, lower_modify_rh
20
21 if (trace_use) call da_trace_entry("da_check_rh")
22
23
24 upper_modify_rh = 95.0
25 lower_modify_rh = 11.0
26
27 ! To get the grid%xa%rh for the moisture adjustments
28
29 if (cv_options_hum == cv_options_hum_specific_humidity) then
30 call da_tpq_to_rh_lin (grid)
31 end if
32
33 do i=its,ite
34 do j=jts,jte
35
36 tol_adjust_moist = 0.0
37 tol_moist = 0.0
38
39 do k=kts,kte
40 dz(k)=grid%xb%hf(i,j,k+1)-grid%xb%hf(i,j,k)
41
42 imod(k) = 0
43 x_qs(k) = 0.0
44 rhtol(k) = grid%xb%rh(i,j,k)+grid%xa%rh(i,j,k)
45 end do
46
47 do k=kts,kte
48 if (rhtol(k) .gt. maximum_rh) then
49 oldrha = grid%xa%rh(i,j,k)
50 ! modify grid%xa%rh
51 grid%xa%rh(i,j,k) = maximum_rh - grid%xb%rh(i,j,k)
52
53 call da_tp_to_qs(grid%xb%t(i,j,k)+grid%xa%t(i,j,k), &
54 grid%xb%p(i,j,k)+grid%xa%p(i,j,k), es, x_qs(k))
55
56 ! calculate grid%xa%q
57 call da_tprh_to_q_lin1(grid%xb%t(i,j,k), grid%xb%p(i,j,k), &
58 grid%xb%es(i,j,k), grid%xb%q(i,j,k), grid%xb%rh(i,j,k), grid%xa%t(i,j,k), &
59 grid%xa%p(i,j,k), grid%xa%rh(i,j,k), grid%xa%q(i,j,k))
60
61 tol_adjust_moist = tol_adjust_moist + x_qs(k)*(oldrha - &
62 grid%xa%rh(i,j,k))* dz(k)*(grid%xb%rho(i,j,k)+grid%xa%rho(i,j,k))
63 imod(k)=-1
64 end if
65
66 if (rhtol(k) .lt. minimum_rh) then
67 oldrha = grid%xa%rh(i,j,k)
68 grid%xa%rh(i,j,k) = minimum_rh - grid%xb%rh(i,j,k)
69 call da_tp_to_qs(grid%xb%t(i,j,k)+grid%xa%t(i,j,k), &
70 grid%xb%p(i,j,k)+grid%xa%p(i,j,k), es, x_qs(k))
71
72 call da_tprh_to_q_lin1(grid%xb%t(i,j,k), grid%xb%p(i,j,k), &
73 grid%xb%es(i,j,k), grid%xb%q(i,j,k), grid%xb%rh(i,j,k), grid%xa%t(i,j,k), &
74 grid%xa%p(i,j,k), grid%xa%rh(i,j,k), grid%xa%q(i,j,k))
75
76
77 tol_adjust_moist = tol_adjust_moist + x_qs(k)*(oldrha - &
78 grid%xa%rh(i,j,k))* dz(k)*(grid%xb%rho(i,j,k)+grid%xa%rho(i,j,k))
79 imod(k)=-1
80 end if
81 end do
82
83 if (tol_adjust_moist .gt. 0.0) then
84 do k=kts,kte
85 if (rhtol(k) .lt. upper_modify_rh .and. imod(k) .eq. 0) then
86 call da_tp_to_qs(grid%xb%t(i,j,k)+grid%xa%t(i,j,k), &
87 grid%xb%p(i,j,k)+grid%xa%p(i,j,k),es,x_qs(k))
88
89 each_moist = rhtol(k)*x_qs(k)* &
90 dz(k)*(grid%xb%rho(i,j,k)+grid%xa%rho(i,j,k))
91 tol_moist = tol_moist + each_moist
92 imod(k) = 1
93 end if
94 end do
95 end if
96
97 if (tol_adjust_moist .lt. 0.0) then
98 do k=kts,kte
99 if (rhtol(k) .gt. lower_modify_rh .and. imod(k) .eq. 0) then
100 call da_tp_to_qs(grid%xb%t(i,j,k)+grid%xa%t(i,j,k), &
101 grid%xb%p(i,j,k)+grid%xa%p(i,j,k), es, x_qs(k))
102
103 each_moist = rhtol(k)*x_qs(k)* &
104 dz(k)*(grid%xb%rho(i,j,k)+grid%xa%rho(i,j,k))
105 tol_moist = tol_moist + each_moist
106 imod(k) = 1
107 end if
108 end do
109 end if
110
111 if (tol_moist > 0) then
112 weight = tol_adjust_moist/tol_moist
113 do k=kts,kte
114 if (imod(k) .eq. 1) then
115 grid%xa%rh(i,j,k) = grid%xa%rh(i,j,k)+(grid%xb%rh(i,j,k)+grid%xa%rh(i,j,k))*weight
116 call da_tprh_to_q_lin1(grid%xb%t(i,j,k), grid%xb%p(i,j,k), grid%xb%es(i,j,k), &
117 grid%xb%q(i,j,k), grid%xb%rh(i,j,k), grid%xa%t(i,j,k), &
118 grid%xa%p(i,j,k), grid%xa%rh(i,j,k), grid%xa%q(i,j,k))
119
120 end if
121 end do
122 end if
123 end do
124 end do
125
126 if (trace_use) call da_trace_exit("da_check_rh")
127
128 end subroutine da_check_rh
129
130