CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
radiation_tools.F90
1
3
5 use machine, only: &
6 kind_phys ! Working type
7 implicit none
8
9 real(kind_phys) :: &
10 rrtmgp_minp, & ! Minimum pressure allowed in RRTMGP
11 rrtmgp_mint ! Minimum temperature allowed in RRTMGP
12contains
13
15 subroutine cmp_tlev(nCol,nLev,minP,p_lay,t_lay,p_lev,tsfc,t_lev)
16 ! Inputs
17 integer, intent(in) :: &
18 nCol,nLev
19 real(kind_phys),intent(in) :: &
20 minP
21 real(kind_phys),dimension(nCol),intent(in) :: &
22 tsfc
23 real(kind_phys),dimension(nCol,nLev),intent(in) :: &
24 p_lay,t_lay
25 real(kind_phys),dimension(nCol,nLev+1),intent(in) :: &
26 p_lev
27
28 ! Outputs
29 real(kind_phys),dimension(nCol,nLev+1),intent(out) :: &
30 t_lev
31
32 ! Local
33 integer :: iCol,iLay, iSFC, iTOA
34 logical :: top_at_1
35 real(kind_phys), dimension(nCol,nLev) :: tem2da, tem2db
36
37 top_at_1 = (p_lev(1,1) .lt. p_lev(1, nlev))
38 if (top_at_1) then
39 isfc = nlev
40 itoa = 1
41 else
42 isfc = 1
43 itoa = nlev
44 endif
45
46 if (itoa .eq. 1) then
47 tem2da(1:ncol,2:isfc) = log(p_lay(1:ncol,2:isfc))
48 tem2db(1:ncol,2:isfc) = log(p_lev(1:ncol,2:isfc))
49 do icol = 1, ncol
50 tem2da(icol,1) = log(p_lay(icol,1) )
51 tem2db(icol,1) = log(max(minp, p_lev(icol,1)) )
52 tem2db(icol,isfc) = log(p_lev(icol,isfc) )
53 enddo
54 t_lev(1:ncol,1) = t_lay(1:ncol,itoa)
55 do ilay = 2, isfc
56 do icol = 1, ncol
57 t_lev(icol,ilay) = t_lay(icol,ilay) + (t_lay(icol,ilay-1) - t_lay(icol,ilay))&
58 * (tem2db(icol,ilay) - tem2da(icol,ilay)) &
59 / (tem2da(icol,ilay-1) - tem2da(icol,ilay))
60 enddo
61 enddo
62 t_lev(1:ncol,isfc+1) = tsfc(1:ncol)
63 else
64 tem2da(1:ncol,2:itoa) = log(p_lay(1:ncol,2:itoa))
65 tem2db(1:ncol,2:itoa) = log(p_lev(1:ncol,2:itoa))
66 do icol = 1, ncol
67 tem2da(icol,1) = log(p_lay(icol,1))
68 tem2db(icol,1) = log(p_lev(icol,1))
69 tem2db(icol,itoa) = log(max(minp, p_lev(icol,itoa)) )
70 enddo
71
72 t_lev(1:ncol,1) = tsfc(1:ncol)
73 do ilay = 1, itoa-1
74 do icol = 1, ncol
75 t_lev(icol,ilay+1) = t_lay(icol,ilay) + (t_lay(icol,ilay+1) - t_lay(icol,ilay))&
76 * (tem2db(icol,ilay+1) - tem2da(icol,ilay)) &
77 / (tem2da(icol,ilay+1) - tem2da(icol,ilay))
78 enddo
79 enddo
80 t_lev(1:ncol,itoa+1) = t_lay(1:ncol,itoa)
81 endif
82
83 end subroutine cmp_tlev
84
86 subroutine check_error_msg(routine_name, error_msg)
87 character(len=*), intent(in) :: &
88 error_msg, routine_name
89
90 if(error_msg /= "") then
91 print*,"ERROR("//trim(routine_name)//"): "
92 print*,trim(error_msg)
93 return
94 end if
95 end subroutine check_error_msg
96
97end module radiation_tools