da_tprh_to_q_adj.inc
References to this file elsewhere.
1 subroutine da_tprh_to_q_adj (grid)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Adjoint of da_tprh_to_q_adj.
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 type (domain), intent(inout) :: grid
10
11 integer :: is, ie ! 1st dim. end points.
12 integer :: js, je ! 2nd dim. end points.
13 integer :: ks, ke ! 3rd dim. end points.
14 integer :: i, j, k ! Loop counter.
15 real :: qs_prime_over_qs(its:ite,jts:jte,kts:kte) ! Temp.
16
17 if (trace_use) call da_trace_entry("da_tprh_to_q_adj")
18
19 !---------------------------------------------------------------------------
20 ! [1.0] initialise:
21 !---------------------------------------------------------------------------
22
23 is = its; ie = ite
24 js = jts; je = jte
25 ks = kts; ke = kte
26
27 if (test_wrfvar) then
28 is = its-1
29 js = jts-1
30
31 ie = ite+1
32 je = jte+1
33
34 if ( is < ids ) is = ids
35 if ( js < jds ) js = jds
36
37 if ( ie > ide ) ie = ide
38 if ( je > jde ) je = jde
39 end if
40
41 !---------------------------------------------------------------------------
42 ! [2.0] Calculate relative humidity increment:
43 !---------------------------------------------------------------------------
44
45 do k = ks, ke
46 do j = js, je
47 do i = is, ie
48 qs_prime_over_qs(i,j,k) = grid%xb % q(i,j,k) * grid%xa % q(i,j,k)
49
50 grid%xa % rh(i,j,k) = grid%xa % rh(i,j,k) + qs_prime_over_qs(i,j,k) / &
51 grid%xb % rh(i,j,k)
52 end do
53 end do
54 end do
55
56 !---------------------------------------------------------------------------
57 ! [2.0] Calculate saturation specific humidity ratio qs~/qs:
58 !---------------------------------------------------------------------------
59
60 call da_tp_to_qs_adj (grid, qs_prime_over_qs )
61
62 if (trace_use) call da_trace_exit("da_tprh_to_q_adj")
63
64 end subroutine da_tprh_to_q_adj
65
66 !subroutine da_tprh_to_q_adj( t, p, es, q, rh, &
67 ! t_prime, p_prime, rh_prime, q_prime, n )
68
69 !---------------------------------------------------------------------------
70 ! Purpose: Adjoint of da_tprh_to_q_adj.
71 !---------------------------------------------------------------------------
72
73 ! implicit none
74
75 ! integer i, n
76 ! real t, es, p, q, rh,t_prime, p_prime, rh_prime, q_prime
77 ! dimension t (n) ! Temperature.
78 ! dimension es (n) ! Saturation vapour pressure.
79 ! dimension p (n) ! Pressure.
80 ! dimension q (n) ! Specific humidity.
81 ! dimension rh (n) ! Relative Humidity.
82 ! dimension t_prime (n) ! Temperature increment.
83 ! dimension p_prime (n) ! Pressure increment.
84 ! dimension rh_prime(n) ! Pressure increment.
85 ! dimension q_prime (n) ! Pressure increment.
86
87 ! real temp, qs_prime_over_qs ! Temporary storage.
88 ! dimension qs_prime_over_qs(n) ! qs~/qs.
89
90 ! do i = 1,n
91 ! temp = q(i) * q_prime(i)
92
93 !---------------------------------------------------------------------------
94 ! [2.0] Calculate relative humidity increment:
95 !---------------------------------------------------------------------------
96
97 ! rh_prime(i) = rh_prime(i) + temp / rh(i)
98 ! qs_prime_over_qs(i) = temp
99 ! end do
100
101 !---------------------------------------------------------------------------
102 ! [1.0] Calculate saturation specific humidity ratio qs~/qs:
103 !---------------------------------------------------------------------------
104
105 ! call da_tp_to_qs_adj( t, p, es, t_prime, p_prime, qs_prime_over_qs, n )
106
107 !end subroutine da_tprh_to_q_adj
108
109