da_find_layer_tl.inc
References to this file elsewhere.
1 subroutine da_find_layer_tl(layer,tv,pre,pre_ma,tv_ma,ks,ke,TGL_tv,TGL_pre_ma,TGL_tv_ma)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: tangent-linear routine for da_find_layer
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 integer, intent(in) :: ks, ke
10 integer, intent(out) :: layer
11 real, intent(in) :: tv_ma(ks-1:ke+1)
12 real, intent(inout) :: pre_ma(ks-1:ke+1)
13 real, intent(in) :: TGL_tv_ma(ks-1:ke+1)
14 real, intent(inout) :: TGL_pre_ma(ks-1:ke+1)
15 real, intent(in) :: pre
16 real, intent(out) :: tv
17 real, intent(out) :: TGL_tv
18
19 real :: TGL_alpha
20 integer :: k
21 real :: alpha, coef1, coef2
22
23 if (trace_use_frequent) call da_trace_entry("da_find_layer_tl")
24
25 ! coef1, coef2 are temporarily used in this routine
26
27 if (pre >= pre_ma(ks)) then
28 ! Below model bottom
29 layer = ks
30 coef1=log(pre/pre_ma(ks+1))/(pre_ma(ks)* &
31 (log(pre_ma(ks)/pre_ma(ks+1)))**2)
32 coef2=log(pre_ma(ks)/pre)/(pre_ma(ks+1)* &
33 (log(pre_ma(ks)/pre_ma(ks+1)))**2)
34 TGL_alpha = coef1 * TGL_pre_ma(ks) + coef2 * TGL_pre_ma(ks+1)
35 alpha = log(pre_ma(ks)/pre)/log(pre_ma(ks)/pre_ma(ks+1))
36
37 TGL_tv = (1.0-alpha)*TGL_tv_ma(ks) + &
38 (tv_ma(ks+1)-tv_ma(ks))*TGL_alpha + &
39 alpha*TGL_tv_ma(ks+1)
40 TGL_pre_ma(ks-1) = 0.0
41 tv = tv_ma(ks) * (1.0-alpha) + tv_ma(ks+1) * alpha
42 pre_ma(ks-1) = pre
43 else if (pre <= pre_ma(ke)) then
44 ! Above model top
45 layer = ke+1
46 coef1=log(pre/pre_ma(ke))/(pre_ma(ke-1)* &
47 (log(pre_ma(ke-1)/pre_ma(ke)))**2)
48 coef2=log(pre_ma(ke-1)/pre)/(pre_ma(ke)* &
49 (log(pre_ma(ke-1)/pre_ma(ke)))**2)
50 TGL_alpha = coef1 * TGL_pre_ma(ke-1) + coef2 * TGL_pre_ma(ke)
51 alpha = log(pre_ma(ke-1)/pre)/log(pre_ma(ke-1)/pre_ma(ke))
52
53 TGL_tv = (1.0-alpha)*TGL_tv_ma(ke-1) + &
54 (tv_ma(ke)-tv_ma(ke-1))*TGL_alpha + &
55 alpha*TGL_tv_ma(ke)
56 TGL_pre_ma(ke+1) = 0.0
57 tv = tv_ma(ke-1) * (1.0-alpha) + tv_ma(ke) * alpha
58 pre_ma(ke+1) = pre
59 else
60 ! Between model layers
61 do k=ks,ke-1
62 if (pre>=pre_ma(k+1) .and. pre<pre_ma(k)) then
63 layer = k+1
64 coef1=log(pre/pre_ma(k+1))/(pre_ma(k)* &
65 (log(pre_ma(k)/pre_ma(k+1)))**2)
66 coef2=log(pre_ma(k)/pre)/(pre_ma(k+1)* &
67 (log(pre_ma(k)/pre_ma(k+1)))**2)
68 TGL_alpha = coef1 * TGL_pre_ma(k) + coef2 * TGL_pre_ma(k+1)
69 alpha = log(pre_ma(k)/pre)/log(pre_ma(k)/pre_ma(k+1))
70 TGL_tv = (1.0-alpha)*TGL_tv_ma(k) + &
71 (tv_ma(k+1)-tv_ma(k))*TGL_alpha + &
72 alpha*TGL_tv_ma(k+1)
73 tv = tv_ma(k) * (1.0-alpha) + tv_ma(k+1) * alpha
74 exit
75 end if
76 end do
77 end if
78
79 if (trace_use_frequent) call da_trace_exit("da_find_layer_tl")
80
81 end subroutine da_find_layer_tl
82
83