CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_rrtmg_post.F90
1
7
9 contains
10
21 subroutine gfs_rrtmg_post_run (im, km, kmp1, lm, ltp, kt, kb, kd, nspc1, &
22 nfxr, nday, lsswr, lslwr, lssav, fhlwr, fhswr, raddt, coszen, &
23 coszdg, prsi, tgrs, aerodp, cldsa, mtopa, mbota, clouds1, &
24 cldtaulw, cldtausw, sfcflw, sfcfsw, topflw, topfsw, scmpsw, &
25 fluxr, total_albedo, errmsg, errflg)
26
27 use machine, only: kind_phys
31
32 implicit none
33
34 ! Interface variables
35 integer, intent(in) :: im, km, kmp1, lm, ltp, kt, kb, kd, &
36 nspc1, nfxr, nday
37 logical, intent(in) :: lsswr, lslwr, lssav
38 real(kind=kind_phys), intent(in) :: raddt, fhlwr, fhswr
39
40 real(kind=kind_phys), dimension(im), intent(in) :: coszen, coszdg
41
42 real(kind=kind_phys), dimension(im,kmp1), intent(in) :: prsi
43 real(kind=kind_phys), dimension(im,km), intent(in) :: tgrs
44
45 real(kind=kind_phys), dimension(im,NSPC1), intent(in) :: aerodp
46 real(kind=kind_phys), dimension(im,5), intent(in) :: cldsa
47 integer, dimension(im,3), intent(in) :: mbota, mtopa
48 real(kind=kind_phys), dimension(im,lm+LTP), intent(in) :: clouds1
49 real(kind=kind_phys), dimension(im,lm+LTP), intent(in) :: cldtausw
50 real(kind=kind_phys), dimension(im,lm+LTP), intent(in) :: cldtaulw
51 real(kind=kind_phys), dimension(im), intent(inout) :: total_albedo
52
53 type(sfcflw_type), dimension(im), intent(in) :: sfcflw
54 type(sfcfsw_type), dimension(im), intent(in) :: sfcfsw
55 type(cmpfsw_type), dimension(im), intent(in) :: scmpsw
56 type(topflw_type), dimension(im), intent(in) :: topflw
57 type(topfsw_type), dimension(im), intent(in) :: topfsw
58
59 real(kind=kind_phys), dimension(im,nfxr), intent(inout) :: fluxr
60
61 character(len=*), intent(out) :: errmsg
62 integer, intent(out) :: errflg
63
64 ! Local variables
65 integer :: i, j, k, k1, itop, ibtc
66 real(kind=kind_phys) :: tem0d, tem1, tem2
67
68 ! Initialize CCPP error handling variables
69 errmsg = ''
70 errflg = 0
71
72 if (.not. (lsswr .or. lslwr)) return
73
74! - For time averaged output quantities (including total-sky and
75! clear-sky SW and LW fluxes at TOA and surface; conventional
76! 3-domain cloud amount, cloud top and base pressure, and cloud top
77! temperature; aerosols AOD, etc.), store computed results in
78! corresponding slots of array fluxr with appropriate time weights.
79
80! --- ... collect the fluxr data for wrtsfc
81
82 if (lssav) then
83 if (lsswr) then
84 do i=1,im
85! fluxr(i,34) = fluxr(i,34) + fhswr*aerodp(i,1) ! total aod at 550nm
86! fluxr(i,35) = fluxr(i,35) + fhswr*aerodp(i,2) ! DU aod at 550nm
87! fluxr(i,36) = fluxr(i,36) + fhswr*aerodp(i,3) ! BC aod at 550nm
88! fluxr(i,37) = fluxr(i,37) + fhswr*aerodp(i,4) ! OC aod at 550nm
89! fluxr(i,38) = fluxr(i,38) + fhswr*aerodp(i,5) ! SU aod at 550nm
90! fluxr(i,39) = fluxr(i,39) + fhswr*aerodp(i,6) ! SS aod at 550nm
91 fluxr(i,34) = aerodp(i,1) ! total aod at 550nm
92 fluxr(i,35) = aerodp(i,2) ! DU aod at 550nm
93 fluxr(i,36) = aerodp(i,3) ! BC aod at 550nm
94 fluxr(i,37) = aerodp(i,4) ! OC aod at 550nm
95 fluxr(i,38) = aerodp(i,5) ! SU aod at 550nm
96 fluxr(i,39) = aerodp(i,6) ! SS aod at 550nm
97 enddo
98 endif
99
100! --- save lw toa and sfc fluxes
101 if (lslwr) then
102 do i=1,im
103! --- lw total-sky fluxes
104 fluxr(i,1 ) = fluxr(i,1 ) + fhlwr * topflw(i)%upfxc ! total sky top lw up
105 fluxr(i,19) = fluxr(i,19) + fhlwr * sfcflw(i)%dnfxc ! total sky sfc lw dn
106 fluxr(i,20) = fluxr(i,20) + fhlwr * sfcflw(i)%upfxc ! total sky sfc lw up
107! --- lw clear-sky fluxes
108 fluxr(i,28) = fluxr(i,28) + fhlwr * topflw(i)%upfx0 ! clear sky top lw up
109 fluxr(i,30) = fluxr(i,30) + fhlwr * sfcflw(i)%dnfx0 ! clear sky sfc lw dn
110 fluxr(i,33) = fluxr(i,33) + fhlwr * sfcflw(i)%upfx0 ! clear sky sfc lw up
111 enddo
112 endif
113
114! --- save sw toa and sfc fluxes with proper diurnal sw wgt. coszen=mean cosz over daylight
115! part of sw calling interval, while coszdg= mean cosz over entire interval
116 if (lsswr) then
117 do i = 1, im
118 if (coszen(i) > 0.) then
119! --- sw total-sky fluxes
120! -------------------
121 tem0d = fhswr * coszdg(i) / coszen(i)
122 fluxr(i,2 ) = fluxr(i,2) + topfsw(i)%upfxc * tem0d ! total sky top sw up
123 fluxr(i,3 ) = fluxr(i,3) + sfcfsw(i)%upfxc * tem0d ! total sky sfc sw up
124 fluxr(i,4 ) = fluxr(i,4) + sfcfsw(i)%dnfxc * tem0d ! total sky sfc sw dn
125! --- sw uv-b fluxes
126! --------------
127 fluxr(i,21) = fluxr(i,21) + scmpsw(i)%uvbfc * tem0d ! total sky uv-b sw dn
128 fluxr(i,22) = fluxr(i,22) + scmpsw(i)%uvbf0 * tem0d ! clear sky uv-b sw dn
129! --- sw toa incoming fluxes
130! ----------------------
131 fluxr(i,23) = fluxr(i,23) + topfsw(i)%dnfxc * tem0d ! top sw dn
132! --- sw sfc flux components
133! ----------------------
134 fluxr(i,24) = fluxr(i,24) + scmpsw(i)%visbm * tem0d ! uv/vis beam sw dn
135 fluxr(i,25) = fluxr(i,25) + scmpsw(i)%visdf * tem0d ! uv/vis diff sw dn
136 fluxr(i,26) = fluxr(i,26) + scmpsw(i)%nirbm * tem0d ! nir beam sw dn
137 fluxr(i,27) = fluxr(i,27) + scmpsw(i)%nirdf * tem0d ! nir diff sw dn
138! --- sw clear-sky fluxes
139! -------------------
140 fluxr(i,29) = fluxr(i,29) + topfsw(i)%upfx0 * tem0d ! clear sky top sw up
141 fluxr(i,31) = fluxr(i,31) + sfcfsw(i)%upfx0 * tem0d ! clear sky sfc sw up
142 fluxr(i,32) = fluxr(i,32) + sfcfsw(i)%dnfx0 * tem0d ! clear sky sfc sw dn
143 endif
144 enddo
145 endif
146
147! --- save total and boundary layer clouds
148
149 if (lsswr .or. lslwr) then
150 do i=1,im
151 fluxr(i,17) = fluxr(i,17) + raddt * cldsa(i,4)
152 fluxr(i,18) = fluxr(i,18) + raddt * cldsa(i,5)
153 enddo
154
155! --- save cld frac,toplyr,botlyr and top temp, note that the order
156! of h,m,l cloud is reversed for the fluxr output.
157! --- save interface pressure (pa) of top/bot
158
159 do j = 1, 3
160 do i = 1, im
161 tem0d = raddt * cldsa(i,j)
162 itop = mtopa(i,j) - kd
163 ibtc = mbota(i,j) - kd
164 fluxr(i, 8-j) = fluxr(i, 8-j) + tem0d
165 fluxr(i,11-j) = fluxr(i,11-j) + tem0d * prsi(i,itop+kt)
166 fluxr(i,14-j) = fluxr(i,14-j) + tem0d * prsi(i,ibtc+kb)
167 fluxr(i,17-j) = fluxr(i,17-j) + tem0d * tgrs(i,itop)
168 enddo
169 enddo
170
171! Anning adds optical depth and emissivity output
172 if (lsswr .and. (nday > 0)) then
173 do j = 1, 3
174 do i = 1, im
175 tem0d = raddt * cldsa(i,j)
176 itop = mtopa(i,j) - kd
177 ibtc = mbota(i,j) - kd
178 tem1 = 0.
179 do k=ibtc,itop
180 tem1 = tem1 + cldtausw(i,k) ! approx .55 um channel
181 enddo
182 fluxr(i,43-j) = fluxr(i,43-j) + tem0d * tem1
183 enddo
184 enddo
185 endif
186
187 if (lslwr) then
188 do j = 1, 3
189 do i = 1, im
190 tem0d = raddt * cldsa(i,j)
191 itop = mtopa(i,j) - kd
192 ibtc = mbota(i,j) - kd
193 tem2 = 0.
194 do k=ibtc,itop
195 tem2 = tem2 + cldtaulw(i,k) ! approx 10. um channel
196 enddo
197 fluxr(i,46-j) = fluxr(i,46-j) + tem0d * (1.0-exp(-tem2))
198 enddo
199 enddo
200 endif
201
202 endif
203
204 endif ! end_if_lssav
205
206! --- The total sky (with clouds) shortwave albedo
207 total_albedo = 0.0
208 if (lsswr) then
209 where(topfsw(:)%dnfxc>0) total_albedo(:) = topfsw(:)%upfxc/topfsw(:)%dnfxc
210 endif
211!
212 end subroutine gfs_rrtmg_post_run
214 end module gfs_rrtmg_post
subroutine gfs_rrtmg_post_run(im, km, kmp1, lm, ltp, kt, kb, kd, nspc1, nfxr, nday, lsswr, lslwr, lssav, fhlwr, fhswr, raddt, coszen, coszdg, prsi, tgrs, aerodp, cldsa, mtopa, mbota, clouds1, cldtaulw, cldtausw, sfcflw, sfcfsw, topflw, topfsw, scmpsw, fluxr, total_albedo, errmsg, errflg)
This module contains LW band parameters set up.
Definition radlw_param.f:61
This module is for specifying the band structures and program parameters used by the RRTMG-SW scheme.
Definition radsw_param.f:62
derived type for LW fluxes at surface
Definition radlw_param.f:87
derived type for LW fluxes at top of atmosphere
Definition radlw_param.f:78
derived type for special components of surface SW fluxes
derived type for SW fluxes at surface
Definition radsw_param.f:89
derived type for SW fluxes at TOA
Definition radsw_param.f:79