CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_radiation_surface.F90
1
5
7
8 use machine, only: kind_phys
9
10 contains
11
17 subroutine gfs_radiation_surface_init (me, ialb, iems, semis_file, con_pi, errmsg, errflg)
18
20
21 implicit none
22
23 integer, intent(in) :: me, ialb, iems
24 character(len=26), intent(in) :: semis_file
25 real(kind_phys), intent(in) :: con_pi
26 character(len=*), intent(out) :: errmsg
27 integer, intent(out) :: errflg
28
29 ! Initialize CCPP error handling variables
30 errmsg = ''
31 errflg = 0
32
33 if ( me == 0 ) then
34 print *,'In GFS_radiation_surface_init, before calling sfc_init'
35 print *,'ialb=',ialb,' iems=',iems
36 end if
37
38 ! Call surface initialization routine
39 call sfc_init ( me, ialb, iems, semis_file, con_pi, errmsg, errflg )
40
41 end subroutine gfs_radiation_surface_init
42
43
47 subroutine gfs_radiation_surface_run ( &
48 ialb, im, nf_albd, frac_grid, lslwr, lsswr, lsm, lsm_noahmp, &
49 lsm_ruc, xlat, xlon, slmsk, lndp_type, n_var_lndp, sfc_alb_pert,&
50 lndp_var_list, lndp_prt_list, landfrac, snodl, snodi, sncovr, &
51 sncovr_ice, fice, zorl, hprime, tsfg, tsfa, tisfc, coszen, &
52 cplice, min_seaice, min_lakeice, lakefrac, use_lake_model, &
53 alvsf, alnsf, alvwf, alnwf, facsf, facwf, &
54 semis_lnd, semis_ice, semis_wat, snoalb, use_cice_alb, con_ttp, &
55 albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, &
56 albdvis_ice, albdnir_ice, albivis_ice, albinir_ice, &
57 semisbase, semis, sfcalb, sfc_alb_dif, errmsg, errflg)
58
60 epsln, &
62
63 implicit none
64
65 integer, intent(in) :: im, nf_albd, ialb
66 logical, intent(in) :: frac_grid, lslwr, lsswr, use_cice_alb, cplice
67 integer, intent(in) :: lsm, lsm_noahmp, lsm_ruc, lndp_type, n_var_lndp
68 real(kind=kind_phys), intent(in) :: min_seaice, min_lakeice, con_ttp
69 integer, dimension(:), intent(in) :: use_lake_model
70
71 real(kind=kind_phys), dimension(:), intent(in) :: xlat, xlon, slmsk, &
72 sfc_alb_pert, &
73 landfrac, lakefrac, &
74 snodl, snodi, sncovr, &
75 sncovr_ice, fice, zorl, &
76 hprime, tsfg, tsfa, tisfc, &
77 coszen, alvsf, alnsf, alvwf, &
78 alnwf, facsf, facwf, snoalb
79 real(kind=kind_phys), dimension(:), intent(in), optional :: lndp_prt_list
80 character(len=3) , dimension(:), intent(in), optional :: lndp_var_list
81 real(kind=kind_phys), dimension(:), intent(in), optional :: albdvis_ice, albdnir_ice, &
82 albivis_ice, albinir_ice
83
84 real(kind=kind_phys), dimension(:), intent(inout) :: albdvis_lnd, albdnir_lnd, &
85 albivis_lnd, albinir_lnd, &
86 semis_lnd, semis_ice, semis_wat
87 real(kind=kind_phys), dimension(:), intent(inout) :: semisbase, semis
88 real(kind=kind_phys), dimension(:,:), intent(inout) :: sfcalb
89 real(kind=kind_phys), dimension(:), intent(inout) :: sfc_alb_dif
90
91 character(len=*), intent(out) :: errmsg
92 integer, intent(out) :: errflg
93
94 ! Local variables
95 integer :: i
96 real(kind=kind_phys) :: lndp_alb
97 real(kind=kind_phys), dimension(im) :: cimin, fracl, fraci, fraco
98 logical, dimension(im) :: icy
99
100 ! Initialize CCPP error handling variables
101 errmsg = ''
102 errflg = 0
103
104 ! Return immediately if neither shortwave nor longwave radiation are called
105 if (.not. lsswr .and. .not. lslwr) return
106
107 do i=1,im
108 if (lakefrac(i) > f_zero) then
109 cimin(i) = min_lakeice
110 else
111 cimin(i) = min_seaice
112 endif
113 enddo
114
115 ! Set up land/ice/ocean fractions for emissivity and albedo calculations
116 if (.not. frac_grid) then
117 do i=1,im
118 if (slmsk(i) == 1) then
119 fracl(i) = f_one
120 fraci(i) = f_zero
121 fraco(i) = f_zero
122 icy(i) = .false.
123 else
124 fracl(i) = f_zero
125 fraco(i) = f_one
126 if(fice(i) < cimin(i)) then
127 fraci(i) = f_zero
128 icy(i) = .false.
129 else
130 fraci(i) = fraco(i) * fice(i)
131 icy(i) = .true.
132 endif
133 fraco(i) = max(f_zero, fraco(i)-fraci(i))
134 endif
135 enddo
136 else
137 do i=1,im
138 fracl(i) = landfrac(i)
139 fraco(i) = max(f_zero, f_one - fracl(i))
140 if(fice(i) < cimin(i)) then
141 fraci(i) = f_zero
142 icy(i) = .false.
143 else
144 fraci(i) = fraco(i) * fice(i)
145 icy(i) = .true.
146 endif
147 fraco(i) = max(f_zero, fraco(i)-fraci(i))
148 enddo
149 endif
150
151 if (lslwr) then
154 call setemis (lsm, lsm_noahmp, lsm_ruc, frac_grid, cplice, &
155 use_lake_model, lakefrac, xlon, xlat, slmsk, &
156! frac_grid, min_seaice, xlon, xlat, slmsk, &
157 snodl, snodi, sncovr, sncovr_ice, zorl, tsfg, &
158 tsfa, hprime, semis_lnd, semis_ice, semis_wat,&
159 im, fracl, fraco, fraci, icy, & ! --- inputs
160 semisbase, semis) ! --- outputs
161 endif
162
163 if (lsswr) then
165 lndp_alb = -999.
166 if (lndp_type==1) then
167 do i =1,n_var_lndp
168 if (lndp_var_list(i) == 'alb') then
169 lndp_alb = lndp_prt_list(i)
170 endif
171 enddo
172 endif
173
176
177 call setalb (slmsk, lsm, lsm_noahmp, lsm_ruc, use_cice_alb, snodi, sncovr, sncovr_ice, &
178 snoalb, zorl, coszen, tsfg, tsfa, hprime, frac_grid, lakefrac, &
179 alvsf, alnsf, alvwf, alnwf, facsf, facwf, fice, tisfc, &
180 albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, &
181 albdvis_ice, albdnir_ice, albivis_ice, albinir_ice, &
182 im, nf_albd, sfc_alb_pert, lndp_alb, fracl, fraco, fraci, icy, ialb, &
183 con_ttp, & ! --- inputs
184 sfcalb ) ! --- outputs
185
187 sfc_alb_dif(:) = max(0.01, 0.5 * (sfcalb(:,2) + sfcalb(:,4)))
188 endif
189
190 end subroutine gfs_radiation_surface_run
191
192 end module gfs_radiation_surface
real(kind=kind_phys), parameter, public f_zero
real(kind=kind_phys), parameter, public f_one
subroutine, public sfc_init(me, ialbflg, iemsflg, semis_file, con_pi, errmsg, errflg)
This subroutine is the initialization program for surface radiation related quantities (albedo,...
subroutine, public setemis(lsm, lsm_noahmp, lsm_ruc, frac_grid, cplice, use_lake_model, lakefrac, xlon, xlat, slmsk, snodl, snodi, sncovr, sncovr_ice, zorlf, tsknf, tairf, hprif, semis_lnd, semis_ice, semis_wat, imax, fracl, fraco, fraci, icy, semisbase, sfcemis)
This subroutine computes surface emissivity for LW radiation.
subroutine, public setalb(slmsk, lsm, lsm_noahmp, lsm_ruc, use_cice_alb, snodi, sncovr, sncovr_ice, snoalb, zorlf, coszf, tsknf, tairf, hprif, frac_grid, lakefrac, alvsf, alnsf, alvwf, alnwf, facsf, facwf, fice, tisfc, lsmalbdvis, lsmalbdnir, lsmalbivis, lsmalbinir, icealbdvis, icealbdnir, icealbivis, icealbinir, imax, nf_albd, albppert, pertalb, fracl, fraco, fraci, icy, ialbflg, con_ttp, sfcalb)
This subroutine computes four components of surface albedos (i.e., vis-nir, direct-diffused) accordin...
real(kind=kind_phys), parameter, public epsln
This module sets up surface albedo for SW radiation and surface emissivity for LW radiation.