Interoperable Physics Driver for NGGPS
grrad.f
Go to the documentation of this file.
1 
3 
122 ! ========================================================== !!!!!
123 ! 'module_radiation_driver' descriptions !!!!!
124 ! ========================================================== !!!!!
125 ! !
126 ! this is the radiation driver module. it prepares atmospheric !
127 ! profiles and invokes main radiation calculations. !
128 ! !
129 ! in module 'module_radiation_driver' there are twe externally !
130 ! callable subroutine: !
131 ! !
132 ! 'radinit' -- initialization routine !
133 ! input: !
134 ! ( si, NLAY, me ) !
135 ! output: !
136 ! ( none ) !
137 ! !
138 ! 'radupdate' -- update time sensitive data used by radiations !
139 ! input: !
140 ! ( idate,jdate,deltsw,deltim,lsswr, me ) !
141 ! output: !
142 ! ( slag,sdec,cdec,solcon ) !
143 ! !
144 ! 'grrad' -- setup and invoke main radiation calls !
145 ! input: !
146 ! ( prsi,prsl,prslk,tgrs,qgrs,tracer,vvl,slmsk, !
147 ! xlon,xlat,tsfc,snowd,sncovr,snoalb,zorl,hprim, !
148 ! alvsf,alnsf,alvwf,alnwf,facsf,facwf,fice,tisfc, !
149 ! sinlat,coslat,solhr,jdate,solcon, !
150 ! cv,cvt,cvb,fcice,frain,rrime,flgmin, !
151 ! icsdsw,icsdlw, ntcw,ncld,ntoz, NTRAC,NFXR, !
152 ! dtlw,dtsw, lsswr,lslwr,lssav, !
153 ! IX, IM, LM, me, lprnt, ipt, kdt,deltaq,sup,cnvw,cnvc, !
154 ! output: !
155 ! htrsw,topfsw,sfcfsw,dswcmp,uswcmp,sfalb,coszen,coszdg, !
156 ! htrlw,topflw,sfcflw,tsflw,semis,cldcov, !
157 ! input/output: !
158 ! fluxr !
159 ! optional output: !
160 ! htrlw0,htrsw0,htrswb,htrlwb !
161 ! !
162 ! !
163 ! external modules referenced: !
164 ! 'module physparam' in 'physparam.f' !
165 ! 'module funcphys' in 'funcphys.f' !
166 ! 'module physcons' in 'physcons.f' !
167 ! !
168 ! 'module module_radiation_gases' in 'radiation_gases.f' !
169 ! 'module module_radiation_aerosols' in 'radiation_aerosols.f' !
170 ! 'module module_radiation_surface' in 'radiation_surface.f' !
171 ! 'module module_radiation_clouds' in 'radiation_clouds.f' !
172 ! !
173 ! 'module module_radsw_cntr_para' in 'radsw_xxxx_param.f' !
174 ! 'module module_radsw_parameters' in 'radsw_xxxx_param.f' !
175 ! 'module module_radsw_main' in 'radsw_xxxx_main.f' !
176 ! !
177 ! 'module module_radlw_cntr_para' in 'radlw_xxxx_param.f' !
178 ! 'module module_radlw_parameters' in 'radlw_xxxx_param.f' !
179 ! 'module module_radlw_main' in 'radlw_xxxx_main.f' !
180 ! !
181 ! where xxxx may vary according to different scheme selection !
182 ! !
183 ! !
184 ! program history log: !
185 ! mm-dd-yy ncep - created program grrad !
186 ! 08-12-03 yu-tai hou - re-written for modulized radiations !
187 ! 11-06-03 yu-tai hou - modified !
188 ! 01-18-05 s. moorthi - NOAH/ICE model changes added !
189 ! 05-10-05 yu-tai hou - modified module structure !
190 ! 12-xx-05 s. moorthi - sfc lw flux adj by mean temperature !
191 ! 02-20-06 yu-tai hou - add time variation for co2 data, and !
192 ! solar const. add sfc emiss change !
193 ! 03-21-06 s. Moorthi - added surface temp over ice !
194 ! 07-28-06 yu-tai hou - add stratospheric vocanic aerosols !
195 ! 03-14-07 yu-tai hou - add generalized spectral band interp !
196 ! for aerosol optical prop. (sw and lw) !
197 ! 04-10-07 yu-tai hou - spectral band sw/lw heating rates !
198 ! 05-04-07 yu-tai hou - make options for clim based and modis !
199 ! based (h. wei and c. marshall) albedo !
200 ! 09-05-08 yu-tai hou - add the initial date and time 'idate' !
201 ! and control param 'ICTM' to the passing param list!
202 ! to handel different time/date requirements for !
203 ! external data (co2, aeros, solcon, ...) !
204 ! 10-10-08 yu-tai hou - add the ICTM=-2 option for combining !
205 ! initial condition data with seasonal cycle from !
206 ! climatology. !
207 ! 03-12-09 yu-tai hou - use two time stamps to keep tracking !
208 ! dates for init cond and fcst time. remove volcanic!
209 ! aerosols data in climate hindcast (ICTM=-2). !
210 ! 03-16-09 yu-tai hou - included sub-column clouds approx. !
211 ! control flags isubcsw/isubclw in initialization !
212 ! subroutine. passed auxiliary cloud control arrays !
213 ! icsdsw/icsdlw (if isubcsw/isubclw =2, it will be !
214 ! the user provided permutation seeds) to the sw/lw !
215 ! radiation calculation programs. also moved cloud !
216 ! overlapping control flags iovrsw/iovrlw from main !
217 ! radiation routines to the initialization routines.!
218 ! 04-02-09 yu-tai hou - modified surface control flag iems to !
219 ! have additional function of if the surface-air !
220 ! interface have the same or different temperature !
221 ! for radiation calculations. !
222 ! 04-03-09 yu-tai hou - modified to add lw surface emissivity !
223 ! as output variable. changed the sign of sfcnsw to !
224 ! be positive value denote solar flux goes into the !
225 ! ground (this is needed to reduce sign confusion !
226 ! in other part of model) !
227 ! 09-09-09 fanglin yang (thru s.moorthi) added QME5 QME6 to E-20!
228 ! 01-09-10 sarah lu - added gocart option, revised grrad for!
229 ! gocart coupling. calling argument modifed: ldiag3 !
230 ! removed; cldcov/fluxr sequence changed; cldcov is !
231 ! changed from accumulative to instant field and !
232 ! from input/output to output field !
233 ! 01-24-10 sarah lu - added aod to fluxr, added prslk and !
234 ! oz to setaer input argument (for gocart coupling),!
235 ! added tau_gocart to setaer output argument (for, !
236 ! aerosol diag by index of nv_aod) !
237 ! 07-08-10 s.moorthi - updated the NEMS version for new physics !
238 ! 07-28-10 yu-tai hou - changed grrad interface to allow all !
239 ! components of sw/lw toa/sfc instantaneous values !
240 ! being passed to the calling program. moved the !
241 ! computaion of sfc net sw flux (sfcnsw) to the !
242 ! calling program. merged carlos' nmmb modification.!
243 ! 07-30-10 s. moorthi - corrected some errors associated with !
244 ! unit changes !
245 ! 12-02-10 s. moorthi/y. hou - removed the use of aerosol flags !
246 ! 'iaersw' 'iaerlw' from radiations and replaced !
247 ! them by using the runtime variable laswflg and !
248 ! lalwflg defined in module radiation_aerosols. !
249 ! also replaced param nspc in grrad with the use of !
250 ! max_num_gridcomp in module radiation_aerosols. !
251 ! jun 2012 yu-tai hou - added sea/land madk 'slmsk' to the !
252 ! argument list of subrotine setaer call for the !
253 ! newly modified horizontal bi-linear interpolation !
254 ! in climatological aerosols schem. also moved the !
255 ! virtual temperature calculations in subroutines !
256 ! 'radiation_clouds' and 'radiation_aerosols' to !
257 ! 'grrad' to reduce repeat comps. renamed var oz as !
258 ! tracer to reflect that it carries various prog !
259 ! tracer quantities. !
260 ! - modified to add 4 compontents of sw !
261 ! surface downward fluxes to the output. (vis/nir; !
262 ! direct/diffused). re-arranged part of the fluxr !
263 ! variable fields and filled the unused slots for !
264 ! the new components. added check print of select !
265 ! data (co2 value for now). !
266 ! - changed the initialization subrution !
267 ! 'radinit' into two parts: 'radinit' is called at !
268 ! the start of model run to set up radiation related!
269 ! fixed parameters; and 'radupdate' is called in !
270 ! the time-loop to update time-varying data sets !
271 ! and module variables. !
272 ! sep 2012 h-m lin/y-t hou added option of extra top layer for !
273 ! models with low toa ceiling. the extra layer will !
274 ! help ozone absorption at higher altitude. !
275 ! nov 2012 yu-tai hou - modified control parameters through !
276 ! module 'physparam'. !
277 ! jan 2013 yu-tai hou - updated subr radupdate for including !
278 ! options of annual/monthly solar constant table. !
279 ! mar 2013 h-m lin/y-t hou corrected a bug in extra top layer !
280 ! when using ferrier microphysics. !
281 ! may 2013 s. mooorthi - removed fpkapx !
282 ! jul 2013 r. sun - added pdf cld and convective cloud water and!
283 ! cover for radiation !
284 ! aug 2013 s. moorthi - port from gfs to nems !
285 ! 13Feb2014 sarah lu - add aerodp to fluxr !
286 ! Apr 2014 Xingren Wu - add sfc SW downward fluxes nir/vis and !
287 ! sfcalb to export for A/O/I coupling !
288 ! jun 2014 y-t hou - revised code to include surface up and !
289 ! down spectral components sw fluxes as output. !
290 ! !
291 !!!!! ========================================================== !!!!!
292 !!!!! end descriptions !!!!!
293 !!!!! ========================================================== !!!!!
294 
295 
296 
297 !========================================!
299 !........................................!
300 !
301  use physparam
302  use physcons, only : eps => con_eps, &
303  & epsm1 => con_epsm1, &
304  & fvirt => con_fvirt &
305  &, rocp => con_rocp
306  use funcphys, only : fpvs
307 
308  use module_radiation_astronomy,only: sol_init, sol_update, coszmn
309  use module_radiation_gases, only : nf_vgas, getgases, getozn, &
310  & gas_init, gas_update
311  use module_radiation_aerosols,only : nf_aesw, nf_aelw, setaer, &
312  & aer_init, aer_update, &
313  & nspc1
314  use module_radiation_surface, only : nf_albd, sfc_init, setalb, &
315  & setemis
316  use module_radiation_clouds, only : nf_clds, cld_init, &
317  & progcld1, progcld2, progcld3,&
318  & diagcld1
319 
320  use module_radsw_parameters, only : topfsw_type, sfcfsw_type, &
321  & profsw_type,cmpfsw_type,nbdsw
322  use module_radsw_main, only : rswinit, swrad
323 
324  use module_radlw_parameters, only : topflw_type, sfcflw_type, &
325  & proflw_type, nbdlw
326  use module_radlw_main, only : rlwinit, lwrad
327 !
328  implicit none
329 !
330  private
331 
332 ! --- version tag and last revision date
333  character(40), parameter :: &
334  & VTAGRAD='NCEP-Radiation_driver v5.2 Jan 2013 '
335 ! & VTAGRAD='NCEP-Radiation_driver v5.1 Nov 2012 '
336 ! & VTAGRAD='NCEP-Radiation_driver v5.0 Aug 2012 '
337 
339 
341  real (kind=kind_phys) :: qmin
343  real (kind=kind_phys) :: qme5
345  real (kind=kind_phys) :: qme6
347  real (kind=kind_phys) :: epsq
348 ! parameter (QMIN=1.0e-10, QME5=1.0e-5, QME6=1.0e-6, EPSQ=1.0e-12)
349  parameter(qmin=1.0e-10, qme5=1.0e-7, qme6=1.0e-7, epsq=1.0e-12)
350 ! parameter (QMIN=1.0e-10, QME5=1.0e-20, QME6=1.0e-20, EPSQ=1.0e-12)
351 
353  real, parameter :: prsmin = 1.0e-6
354 
356  integer :: itsfc =0
357 
359  integer :: month0=0, iyear0=0, monthd=0
361  logical :: loz1st =.true.
362 
365  integer, parameter :: ltp = 0 ! no extra top layer
366 ! integer, parameter :: LTP = 1 ! add an extra top layer
367  logical, parameter :: lextop = (ltp > 0)
368 
369 ! --- publicly accessible module programs:
370 
371  public radinit, radupdate, grrad
372 
373 
374 ! =================
375  contains
376 ! =================
377 
384 !-----------------------------------
385  subroutine radinit( si, NLAY, me )
386 !...................................
387 
388 ! --- inputs:
389 ! & ( si, NLAY, me )
390 ! --- outputs:
391 ! ( none )
392 
393 ! ================= subprogram documentation block ================ !
394 ! !
395 ! subprogram: radinit initialization of radiation calculations !
396 ! !
397 ! usage: call radinit !
398 ! !
399 ! attributes: !
400 ! language: fortran 90 !
401 ! machine: wcoss !
402 ! !
403 ! ==================== definition of variables ==================== !
404 ! !
405 ! input parameters: !
406 ! si : model vertical sigma interface !
407 ! NLAY : number of model vertical layers !
408 ! me : print control flag !
409 ! !
410 ! outputs: (none) !
411 ! !
412 ! external module variables: (in module physparam) !
413 ! isolar : solar constant cntrol flag !
414 ! = 0: use the old fixed solar constant in "physcon" !
415 ! =10: use the new fixed solar constant in "physcon" !
416 ! = 1: use noaa ann-mean tsi tbl abs-scale with cycle apprx!
417 ! = 2: use noaa ann-mean tsi tbl tim-scale with cycle apprx!
418 ! = 3: use cmip5 ann-mean tsi tbl tim-scale with cycl apprx!
419 ! = 4: use cmip5 mon-mean tsi tbl tim-scale with cycl apprx!
420 ! iaerflg : 3-digit aerosol flag (abc for volc, lw, sw) !
421 ! a:=0 use background stratospheric aerosol !
422 ! =1 include stratospheric vocanic aeros !
423 ! b:=0 no topospheric aerosol in lw radiation !
424 ! =1 compute tropspheric aero in 1 broad band for lw !
425 ! =2 compute tropspheric aero in multi bands for lw !
426 ! c:=0 no topospheric aerosol in sw radiation !
427 ! =1 include tropspheric aerosols for sw !
428 ! ico2flg : co2 data source control flag !
429 ! =0: use prescribed global mean co2 (old oper) !
430 ! =1: use observed co2 annual mean value only !
431 ! =2: use obs co2 monthly data with 2-d variation !
432 ! ictmflg : =yyyy#, external data ic time/date control flag !
433 ! = -2: same as 0, but superimpose seasonal cycle !
434 ! from climatology data set. !
435 ! = -1: use user provided external data for the !
436 ! forecast time, no extrapolation. !
437 ! = 0: use data at initial cond time, if not !
438 ! available, use latest, no extrapolation. !
439 ! = 1: use data at the forecast time, if not !
440 ! available, use latest and extrapolation. !
441 ! =yyyy0: use yyyy data for the forecast time, !
442 ! no further data extrapolation. !
443 ! =yyyy1: use yyyy data for the fcst. if needed, do !
444 ! extrapolation to match the fcst time. !
445 ! ioznflg : ozone data source control flag !
446 ! =0: use climatological ozone profile !
447 ! =1: use interactive ozone profile !
448 ! ialbflg : albedo scheme control flag !
449 ! =0: climatology, based on surface veg types !
450 ! =1: modis retrieval based surface albedo scheme !
451 ! iemsflg : emissivity scheme cntrl flag (ab 2-digit integer) !
452 ! a:=0 set sfc air/ground t same for lw radiation !
453 ! =1 set sfc air/ground t diff for lw radiation !
454 ! b:=0 use fixed sfc emissivity=1.0 (black-body) !
455 ! =1 use varying climtology sfc emiss (veg based) !
456 ! =2 future development (not yet) !
457 ! icldflg : cloud optical property scheme control flag !
458 ! =0: use diagnostic cloud scheme !
459 ! =1: use prognostic cloud scheme (default) !
460 ! icmphys : cloud microphysics scheme control flag !
461 ! =1 zhao/carr/sundqvist microphysics scheme !
462 ! =2 brad ferrier microphysics scheme !
463 ! =3 zhao/carr/sundqvist microphysics+pdf cloud & cnvc,cnvw!
464 ! iovrsw : control flag for cloud overlap in sw radiation !
465 ! iovrlw : control flag for cloud overlap in lw radiation !
466 ! =0: random overlapping clouds !
467 ! =1: max/ran overlapping clouds !
468 ! isubcsw : sub-column cloud approx control flag in sw radiation !
469 ! isubclw : sub-column cloud approx control flag in lw radiation !
470 ! =0: with out sub-column cloud approximation !
471 ! =1: mcica sub-col approx. prescribed random seed !
472 ! =2: mcica sub-col approx. provided random seed !
473 ! lcrick : control flag for eliminating CRICK !
474 ! =t: apply layer smoothing to eliminate CRICK !
475 ! =f: do not apply layer smoothing !
476 ! lcnorm : control flag for in-cld condensate !
477 ! =t: normalize cloud condensate !
478 ! =f: not normalize cloud condensate !
479 ! lnoprec : precip effect in radiation flag (ferrier microphysics) !
480 ! =t: snow/rain has no impact on radiation !
481 ! =f: snow/rain has impact on radiation !
482 ! ivflip : vertical index direction control flag !
483 ! =0: index from toa to surface !
484 ! =1: index from surface to toa !
485 ! !
486 ! subroutines called: sol_init, aer_init, gas_init, cld_init, !
487 ! sfc_init, rlwinit, rswinit !
488 ! !
489 ! usage: call radinit !
490 ! !
491 ! =================================================================== !
492 !
493  implicit none
494 
495 ! --- inputs:
496  integer, intent(in) :: NLAY, me
497 
498  real (kind=kind_phys), intent(in) :: si(:)
499 
500 ! --- outputs: (none, to module variables)
501 
502 ! --- locals:
503 
504 !
505 !===> ... begin here
506 !
509  itsfc = iemsflg / 10 ! sfc air/ground temp control
510  loz1st = (ioznflg == 0) ! first-time clim ozone data read flag
511  month0 = 0
512  iyear0 = 0
513  monthd = 0
514 
515  if (me == 0) then
516 ! print *,' NEW RADIATION PROGRAM STRUCTURES -- SEP 01 2004'
517  print *,' NEW RADIATION PROGRAM STRUCTURES BECAME OPER. ', &
518  & ' May 01 2007'
519  print *, vtagrad !print out version tag
520  print *,' - Selected Control Flag settings: ICTMflg=',ictmflg, &
521  & ' ISOLar =',isolar, ' ICO2flg=',ico2flg,' IAERflg=',iaerflg, &
522  & ' IALBflg=',ialbflg,' IEMSflg=',iemsflg,' ICLDflg=',icldflg, &
523  & ' ICMPHYS=',icmphys,' IOZNflg=',ioznflg
524  print *,' IVFLIP=',ivflip,' IOVRSW=',iovrsw,' IOVRLW=',iovrlw, &
525  & ' ISUBCSW=',isubcsw,' ISUBCLW=',isubclw
526 ! write(0,*)' IVFLIP=',ivflip,' IOVRSW=',iovrsw,' IOVRLW=',iovrlw,&
527 ! & ' ISUBCSW=',isubcsw,' ISUBCLW=',isubclw
528  print *,' LCRICK=',lcrick,' LCNORM=',lcnorm,' LNOPREC=',lnoprec
529  print *,' LTP =',ltp,', add extra top layer =',lextop
530 
531  if ( ictmflg==0 .or. ictmflg==-2 ) then
532  print *,' Data usage is limited by initial condition!'
533  print *,' No volcanic aerosols'
534  endif
535 
536  if ( isubclw == 0 ) then
537  print *,' - ISUBCLW=',isubclw,' No McICA, use grid ', &
538  & 'averaged cloud in LW radiation'
539  elseif ( isubclw == 1 ) then
540  print *,' - ISUBCLW=',isubclw,' Use McICA with fixed ', &
541  & 'permutation seeds for LW random number generator'
542  elseif ( isubclw == 2 ) then
543  print *,' - ISUBCLW=',isubclw,' Use McICA with random ', &
544  & 'permutation seeds for LW random number generator'
545  else
546  print *,' - ERROR!!! ISUBCLW=',isubclw,' is not a ', &
547  & 'valid option '
548  stop
549  endif
550 
551  if ( isubcsw == 0 ) then
552  print *,' - ISUBCSW=',isubcsw,' No McICA, use grid ', &
553  & 'averaged cloud in SW radiation'
554  elseif ( isubcsw == 1 ) then
555  print *,' - ISUBCSW=',isubcsw,' Use McICA with fixed ', &
556  & 'permutation seeds for SW random number generator'
557  elseif ( isubcsw == 2 ) then
558  print *,' - ISUBCSW=',isubcsw,' Use McICA with random ', &
559  & 'permutation seeds for SW random number generator'
560  else
561  print *,' - ERROR!!! ISUBCSW=',isubcsw,' is not a ', &
562  & 'valid option '
563  stop
564  endif
565 
566  if ( isubcsw /= isubclw ) then
567  print *,' - *** Notice *** ISUBCSW /= ISUBCLW !!!', &
568  & isubcsw, isubclw
569  endif
570  endif
571 
587 ! Initialization
588 
589  call sol_init ( me ) ! --- ... astronomy initialization routine
590 
591  call aer_init ( nlay, me ) ! --- ... aerosols initialization routine
592 
593  call gas_init ( me ) ! --- ... co2 and other gases initialization routine
594 
595  call sfc_init ( me ) ! --- ... surface initialization routine
596 
597  call cld_init ( si, nlay, me) ! --- ... cloud initialization routine
598 
599  call rlwinit ( me ) ! --- ... lw radiation initialization routine
600 
601  call rswinit ( me ) ! --- ... sw radiation initialization routine
602 !
603  return
604 !...................................
605  end subroutine radinit
606 !-----------------------------------
608 
624 !-----------------------------------
625  subroutine radupdate( idate,jdate,deltsw,deltim,lsswr, me, &
626  & slag,sdec,cdec,solcon)
627 !...................................
628 
629 ! ================= subprogram documentation block ================ !
630 ! !
631 ! subprogram: radupdate calls many update subroutines to check and !
632 ! update radiation required but time varying data sets and module !
633 ! variables. !
634 ! !
635 ! usage: call radupdate !
636 ! !
637 ! attributes: !
638 ! language: fortran 90 !
639 ! machine: ibm sp !
640 ! !
641 ! ==================== definition of variables ==================== !
642 ! !
643 ! input parameters: !
644 ! idate(8) : ncep absolute date and time of initial condition !
645 ! (yr, mon, day, t-zone, hr, min, sec, mil-sec) !
646 ! jdate(8) : ncep absolute date and time at fcst time !
647 ! (yr, mon, day, t-zone, hr, min, sec, mil-sec) !
648 ! deltsw : sw radiation calling frequency in seconds !
649 ! deltim : model timestep in seconds !
650 ! lsswr : logical flags for sw radiation calculations !
651 ! me : print control flag !
652 ! !
653 ! outputs: !
654 ! slag : equation of time in radians !
655 ! sdec, cdec : sin and cos of the solar declination angle !
656 ! solcon : sun-earth distance adjusted solar constant (w/m2) !
657 ! !
658 ! external module variables: !
659 ! isolar : solar constant cntrl (in module physparam) !
660 ! = 0: use the old fixed solar constant in "physcon" !
661 ! =10: use the new fixed solar constant in "physcon" !
662 ! = 1: use noaa ann-mean tsi tbl abs-scale with cycle apprx!
663 ! = 2: use noaa ann-mean tsi tbl tim-scale with cycle apprx!
664 ! = 3: use cmip5 ann-mean tsi tbl tim-scale with cycl apprx!
665 ! = 4: use cmip5 mon-mean tsi tbl tim-scale with cycl apprx!
666 ! ictmflg : =yyyy#, external data ic time/date control flag !
667 ! = -2: same as 0, but superimpose seasonal cycle !
668 ! from climatology data set. !
669 ! = -1: use user provided external data for the !
670 ! forecast time, no extrapolation. !
671 ! = 0: use data at initial cond time, if not !
672 ! available, use latest, no extrapolation. !
673 ! = 1: use data at the forecast time, if not !
674 ! available, use latest and extrapolation. !
675 ! =yyyy0: use yyyy data for the forecast time, !
676 ! no further data extrapolation. !
677 ! =yyyy1: use yyyy data for the fcst. if needed, do !
678 ! extrapolation to match the fcst time. !
679 ! !
680 ! module variables: !
681 ! loz1st : first-time clim ozone data read flag !
682 ! !
683 ! subroutines called: sol_update, aer_update, gas_update !
684 ! !
685 ! =================================================================== !
686 !
687  implicit none
688 
689 ! --- inputs:
690  integer, intent(in) :: idate(:), jdate(:), me
691  logical, intent(in) :: lsswr
692 
693  real (kind=kind_phys), intent(in) :: deltsw, deltim
694 
695 ! --- outputs:
696  real (kind=kind_phys), intent(out) :: slag, sdec, cdec, solcon
697 
698 ! --- locals:
699  integer :: iyear, imon, iday, ihour
700  integer :: kyear, kmon, kday, khour
701 
702  logical :: lmon_chg ! month change flag
703  logical :: lco2_chg ! cntrl flag for updating co2 data
704  logical :: lsol_chg ! cntrl flag for updating solar constant
705 !
706 !===> ... begin here
707 !
710 ! --- ... time stamp at fcst time
711 
712  iyear = jdate(1)
713  imon = jdate(2)
714  iday = jdate(3)
715  ihour = jdate(5)
716 
717 ! --- ... set up time stamp used for green house gases (** currently co2 only)
718 
719  if ( ictmflg==0 .or. ictmflg==-2 ) then ! get external data at initial condition time
720  kyear = idate(1)
721  kmon = idate(2)
722  kday = idate(3)
723  khour = idate(5)
724  else ! get external data at fcst or specified time
725  kyear = iyear
726  kmon = imon
727  kday = iday
728  khour = ihour
729  endif ! end if_ictmflg_block
730 
731  if ( month0 /= imon ) then
732  lmon_chg = .true.
733  month0 = imon
734  else
735  lmon_chg = .false.
736  endif
737 
740  if (lsswr) then
741 
742  if ( isolar == 0 .or. isolar == 10 ) then
743  lsol_chg = .false.
744  elseif ( iyear0 /= iyear ) then
745  lsol_chg = .true.
746  else
747  lsol_chg = ( isolar==4 .and. lmon_chg )
748  endif
749  iyear0 = iyear
750 
751  call sol_update &
752 ! --- inputs:
753  & ( jdate,kyear,deltsw,deltim,lsol_chg, me, &
754 ! --- outputs:
755  & slag,sdec,cdec,solcon &
756  & )
757 
758  endif ! end_if_lsswr_block
759 
762  if ( lmon_chg ) then
763  call aer_update ( iyear, imon, me )
764  endif
765 
768  if ( monthd /= kmon ) then
769  monthd = kmon
770  lco2_chg = .true.
771  else
772  lco2_chg = .false.
773  endif
774 
775  call gas_update ( kyear,kmon,kday,khour,loz1st,lco2_chg, me )
776 
777  if ( loz1st ) loz1st = .false.
778 
780 ! call sfc_update ( iyear, imon, me )
781 
783 ! call cld_update ( iyear, imon, me )
784 !
785  return
786 !...................................
787  end subroutine radupdate
788 !-----------------------------------
790 
942 !-----------------------------------
943  subroutine grrad &
944  & ( prsi,prsl,prslk,tgrs,qgrs,tracer,vvl,slmsk, & ! --- inputs
945  & xlon,xlat,tsfc,snowd,sncovr,snoalb,zorl,hprim, &
946  & alvsf,alnsf,alvwf,alnwf,facsf,facwf,fice,tisfc, &
947  & sinlat,coslat,solhr,jdate,solcon, &
948  & cv,cvt,cvb,fcice,frain,rrime,flgmin, &
949  & icsdsw,icsdlw, ntcw,ncld,ntoz, ntrac,nfxr, &
950  & dtlw,dtsw, lsswr,lslwr,lssav, uni_cld,lmfshal,lmfdeep2, &
951  & ix, im, lm, me, lprnt, ipt, kdt, deltaq,sup,cnvw,cnvc, &
952  & htrsw,topfsw,sfcfsw,dswcmp,uswcmp,sfalb,coszen,coszdg, & ! --- outputs:
953  & htrlw,topflw,sfcflw,tsflw,semis,cldcov, &
954  & fluxr & ! --- input/output:
955  &, htrlw0,htrsw0,htrswb,htrlwb & !! --- optional outputs:
956  & )
958 ! ================= subprogram documentation block ================ !
959 ! !
960 ! this program is the driver of radiation calculation subroutines. * !
961 ! It sets up profile variables for radiation input, including * !
962 ! clouds, surface albedos, atmospheric aerosols, ozone, etc. * !
963 ! * !
964 ! usage: call grrad * !
965 ! * !
966 ! subprograms called: * !
967 ! setalb, setemis, setaer, getozn, getgases, * !
968 ! progcld1, progcld2, diagcds, * !
969 ! swrad, lwrad, fpvs * !
970 ! * !
971 ! attributes: * !
972 ! language: fortran 90 * !
973 ! machine: ibm-sp, sgi * !
974 ! * !
975 ! * !
976 ! ==================== definition of variables ==================== !
977 ! !
978 ! input variables: !
979 ! prsi (IX,LM+1) : model level pressure in Pa !
980 ! prsl (IX,LM) : model layer mean pressure Pa !
981 ! prslk (IX,LM) : exner function = (p/p0)**rocp !
982 ! tgrs (IX,LM) : model layer mean temperature in k !
983 ! qgrs (IX,LM) : layer specific humidity in gm/gm !
984 ! tracer(IX,LM,NTRAC):layer prognostic tracer amount/mixing-ratio !
985 ! incl: oz, cwc, aeros, etc. !
986 ! vvl (IX,LM) : layer mean vertical velocity in pa/sec !
987 ! slmsk (IM) : sea/land mask array (sea:0,land:1,sea-ice:2) !
988 ! xlon (IM) : grid longitude in radians, ok for both 0->2pi !
989 ! or -pi -> +pi ranges !
990 ! xlat (IM) : grid latitude in radians, default to pi/2 -> !
991 ! -pi/2 range, otherwise adj in subr called !
992 ! tsfc (IM) : surface temperature in k !
993 ! snowd (IM) : snow depth water equivalent in mm !
994 ! sncovr(IM) : snow cover in fraction !
995 ! snoalb(IM) : maximum snow albedo in fraction !
996 ! zorl (IM) : surface roughness in cm !
997 ! hprim (IM) : topographic standard deviation in m !
998 ! alvsf (IM) : mean vis albedo with strong cosz dependency !
999 ! alnsf (IM) : mean nir albedo with strong cosz dependency !
1000 ! alvwf (IM) : mean vis albedo with weak cosz dependency !
1001 ! alnwf (IM) : mean nir albedo with weak cosz dependency !
1002 ! facsf (IM) : fractional coverage with strong cosz dependen !
1003 ! facwf (IM) : fractional coverage with weak cosz dependency !
1004 ! fice (IM) : ice fraction over open water grid !
1005 ! tisfc (IM) : surface temperature over ice fraction !
1006 ! sinlat(IM) : sine of the grids' corresponding latitudes !
1007 ! coslat(IM) : cosine of the grids' corresponding latitudes !
1008 ! solhr : hour time after 00z at the t-stepe !
1009 ! jdate (8) : current forecast date and time !
1010 ! (yr, mon, day, t-zone, hr, min, sec, mil-sec) !
1011 ! solcon : solar constant (sun-earth distant adjusted) !
1012 ! cv (IM) : fraction of convective cloud !
1013 ! cvt, cvb (IM) : convective cloud top/bottom pressure in pa !
1014 ! fcice : fraction of cloud ice (in ferrier scheme) !
1015 ! frain : fraction of rain water (in ferrier scheme) !
1016 ! rrime : mass ratio of total to unrimed ice ( >= 1 ) !
1017 ! flgmin : minimim large ice fraction !
1018 ! icsdsw/icsdlw : auxiliary cloud control arrays passed to main !
1019 ! (IM) radiations. if isubcsw/isubclw (input to init) !
1020 ! are set to 2, the arrays contains provided !
1021 ! random seeds for sub-column clouds generators !
1022 ! lmfshal : mass-flux shallow conv scheme flag !
1023 ! lmfdeep2 : scale-aware mass-flux deep conv scheme flag !
1024 ! ntcw : =0 no cloud condensate calculated !
1025 ! >0 array index location for cloud condensate !
1026 ! ncld : only used when ntcw .gt. 0 !
1027 ! ntoz : =0 climatological ozone profile !
1028 ! >0 interactive ozone profile !
1029 ! NTRAC : dimension veriable for array oz !
1030 ! NFXR : second dimension of input/output array fluxr !
1031 ! dtlw, dtsw : time duration for lw/sw radiation call in sec !
1032 ! lsswr, lslwr : logical flags for sw/lw radiation calls !
1033 ! lssav : logical flag for store 3-d cloud field !
1034 ! IX,IM : horizontal dimention and num of used points !
1035 ! LM : vertical layer dimension !
1036 ! me : control flag for parallel process !
1037 ! lprnt : control flag for diagnostic print out !
1038 ! ipt : index for diagnostic printout point !
1039 ! kdt : time-step number !
1040 ! deltaq : half width of uniform total water distribution !
1041 ! sup : supersaturation in pdf cloud when t is very low!
1042 ! cnvw : layer convective cloud water !
1043 ! cnvc : layer convective cloud cover !
1044 ! !
1045 ! output variables: !
1046 ! htrsw (IX,LM) : total sky sw heating rate in k/sec !
1047 ! topfsw(IM) : sw radiation fluxes at toa, components: !
1048 ! (check module_radsw_parameters for definition) !
1049 ! %upfxc - total sky upward sw flux at toa (w/m**2) !
1050 ! %dnflx - total sky downward sw flux at toa (w/m**2) !
1051 ! %upfx0 - clear sky upward sw flux at toa (w/m**2) !
1052 ! sfcfsw(IM) : sw radiation fluxes at sfc, components: !
1053 ! (check module_radsw_parameters for definition) !
1054 ! %upfxc - total sky upward sw flux at sfc (w/m**2) !
1055 ! %dnfxc - total sky downward sw flux at sfc (w/m**2) !
1056 ! %upfx0 - clear sky upward sw flux at sfc (w/m**2) !
1057 ! %dnfx0 - clear sky downward sw flux at sfc (w/m**2) !
1058 ! dswcmp(IX,4) : dn sfc sw spectral components: !
1059 ! ( :, 1) - total sky sfc downward nir direct flux !
1060 ! ( :, 2) - total sky sfc downward nir diffused flux !
1061 ! ( :, 3) - total sky sfc downward uv+vis direct flux !
1062 ! ( :, 4) - total sky sfc downward uv+vis diff flux !
1063 ! uswcmp(IX,4) : up sfc sw spectral components: !
1064 ! ( :, 1) - total sky sfc upward nir direct flux !
1065 ! ( :, 2) - total sky sfc upward nir diffused flux !
1066 ! ( :, 3) - total sky sfc upward uv+vis direct flux !
1067 ! ( :, 4) - total sky sfc upward uv+vis diff flux !
1068 ! sfalb (IM) : mean surface diffused sw albedo !
1069 ! coszen(IM) : mean cos of zenith angle over rad call period !
1070 ! coszdg(IM) : daytime mean cosz over rad call period !
1071 ! htrlw (IX,LM) : total sky lw heating rate in k/sec !
1072 ! topflw(IM) : lw radiation fluxes at top, component: !
1073 ! (check module_radlw_paramters for definition) !
1074 ! %upfxc - total sky upward lw flux at toa (w/m**2) !
1075 ! %upfx0 - clear sky upward lw flux at toa (w/m**2) !
1076 ! sfcflw(IM) : lw radiation fluxes at sfc, component: !
1077 ! (check module_radlw_paramters for definition) !
1078 ! %upfxc - total sky upward lw flux at sfc (w/m**2) !
1079 ! %upfx0 - clear sky upward lw flux at sfc (w/m**2) !
1080 ! %dnfxc - total sky downward lw flux at sfc (w/m**2) !
1081 ! %dnfx0 - clear sky downward lw flux at sfc (w/m**2) !
1082 ! semis (IM) : surface lw emissivity in fraction !
1083 ! cldcov(IX,LM) : 3-d cloud fraction !
1084 ! tsflw (IM) : surface air temp during lw calculation in k !
1085 ! !
1086 ! input and output variables: !
1087 ! fluxr (IX,NFXR) : to save time accumulated 2-d fields defined as:!
1088 ! 1 - toa total sky upwd lw radiation flux !
1089 ! 2 - toa total sky upwd sw radiation flux !
1090 ! 3 - sfc total sky upwd sw radiation flux !
1091 ! 4 - sfc total sky dnwd sw radiation flux !
1092 ! 5 - high domain cloud fraction !
1093 ! 6 - mid domain cloud fraction !
1094 ! 7 - low domain cloud fraction !
1095 ! 8 - high domain mean cloud top pressure !
1096 ! 9 - mid domain mean cloud top pressure !
1097 ! 10 - low domain mean cloud top pressure !
1098 ! 11 - high domain mean cloud base pressure !
1099 ! 12 - mid domain mean cloud base pressure !
1100 ! 13 - low domain mean cloud base pressure !
1101 ! 14 - high domain mean cloud top temperature !
1102 ! 15 - mid domain mean cloud top temperature !
1103 ! 16 - low domain mean cloud top temperature !
1104 ! 17 - total cloud fraction !
1105 ! 18 - boundary layer domain cloud fraction !
1106 ! 19 - sfc total sky dnwd lw radiation flux !
1107 ! 20 - sfc total sky upwd lw radiation flux !
1108 ! 21 - sfc total sky dnwd sw uv-b radiation flux !
1109 ! 22 - sfc clear sky dnwd sw uv-b radiation flux !
1110 ! 23 - toa incoming solar radiation flux !
1111 ! 24 - sfc vis beam dnwd sw radiation flux !
1112 ! 25 - sfc vis diff dnwd sw radiation flux !
1113 ! 26 - sfc nir beam dnwd sw radiation flux !
1114 ! 27 - sfc nir diff dnwd sw radiation flux !
1115 ! 28 - toa clear sky upwd lw radiation flux !
1116 ! 29 - toa clear sky upwd sw radiation flux !
1117 ! 30 - sfc clear sky dnwd lw radiation flux !
1118 ! 31 - sfc clear sky upwd sw radiation flux !
1119 ! 32 - sfc clear sky dnwd sw radiation flux !
1120 ! 33 - sfc clear sky upwd lw radiation flux !
1121 !optional 34 - aeros opt depth at 550nm (all components) !
1122 ! 35 - aeros opt depth at 550nm for du component !
1123 ! 36 - aeros opt depth at 550nm for bc component !
1124 ! 37 - aeros opt depth at 550nm for oc component !
1125 ! 38 - aeros opt depth at 550nm for su component !
1126 ! 39 - aeros opt depth at 550nm for ss component !
1127 ! !
1128 ! optional output variables: !
1129 ! htrswb(IX,LM,NBDSW) : spectral band total sky sw heating rate !
1130 ! htrlwb(IX,LM,NBDLW) : spectral band total sky lw heating rate !
1131 ! !
1132 ! !
1133 ! definitions of internal variable arrays: !
1134 ! !
1135 ! 1. fixed gases: (defined in 'module_radiation_gases') !
1136 ! gasvmr(:,:,1) - co2 volume mixing ratio !
1137 ! gasvmr(:,:,2) - n2o volume mixing ratio !
1138 ! gasvmr(:,:,3) - ch4 volume mixing ratio !
1139 ! gasvmr(:,:,4) - o2 volume mixing ratio !
1140 ! gasvmr(:,:,5) - co volume mixing ratio !
1141 ! gasvmr(:,:,6) - cf11 volume mixing ratio !
1142 ! gasvmr(:,:,7) - cf12 volume mixing ratio !
1143 ! gasvmr(:,:,8) - cf22 volume mixing ratio !
1144 ! gasvmr(:,:,9) - ccl4 volume mixing ratio !
1145 ! !
1146 ! 2. cloud profiles: (defined in 'module_radiation_clouds') !
1147 ! --- for prognostic cloud --- !
1148 ! clouds(:,:,1) - layer total cloud fraction !
1149 ! clouds(:,:,2) - layer cloud liq water path !
1150 ! clouds(:,:,3) - mean effective radius for liquid cloud !
1151 ! clouds(:,:,4) - layer cloud ice water path !
1152 ! clouds(:,:,5) - mean effective radius for ice cloud !
1153 ! clouds(:,:,6) - layer rain drop water path !
1154 ! clouds(:,:,7) - mean effective radius for rain drop !
1155 ! clouds(:,:,8) - layer snow flake water path !
1156 ! clouds(:,:,9) - mean effective radius for snow flake !
1157 ! --- for diagnostic cloud --- !
1158 ! clouds(:,:,1) - layer total cloud fraction !
1159 ! clouds(:,:,2) - layer cloud optical depth !
1160 ! clouds(:,:,3) - layer cloud single scattering albedo !
1161 ! clouds(:,:,4) - layer cloud asymmetry factor !
1162 ! !
1163 ! 3. surface albedo: (defined in 'module_radiation_surface') !
1164 ! sfcalb( :,1 ) - near ir direct beam albedo !
1165 ! sfcalb( :,2 ) - near ir diffused albedo !
1166 ! sfcalb( :,3 ) - uv+vis direct beam albedo !
1167 ! sfcalb( :,4 ) - uv+vis diffused albedo !
1168 ! !
1169 ! 4. sw aerosol profiles: (defined in 'module_radiation_aerosols') !
1170 ! faersw(:,:,:,1)- sw aerosols optical depth !
1171 ! faersw(:,:,:,2)- sw aerosols single scattering albedo !
1172 ! faersw(:,:,:,3)- sw aerosols asymmetry parameter !
1173 ! !
1174 ! 5. lw aerosol profiles: (defined in 'module_radiation_aerosols') !
1175 ! faerlw(:,:,:,1)- lw aerosols optical depth !
1176 ! faerlw(:,:,:,2)- lw aerosols single scattering albedo !
1177 ! faerlw(:,:,:,3)- lw aerosols asymmetry parameter !
1178 ! !
1179 ! 6. sw fluxes at toa: (defined in 'module_radsw_main') !
1180 ! (topfsw_type -- derived data type for toa rad fluxes) !
1181 ! topfsw(:)%upfxc - total sky upward flux at toa !
1182 ! topfsw(:)%dnfxc - total sky downward flux at toa !
1183 ! topfsw(:)%upfx0 - clear sky upward flux at toa !
1184 ! !
1185 ! 7. lw fluxes at toa: (defined in 'module_radlw_main') !
1186 ! (topflw_type -- derived data type for toa rad fluxes) !
1187 ! topflw(:)%upfxc - total sky upward flux at toa !
1188 ! topflw(:)%upfx0 - clear sky upward flux at toa !
1189 ! !
1190 ! 8. sw fluxes at sfc: (defined in 'module_radsw_main') !
1191 ! (sfcfsw_type -- derived data type for sfc rad fluxes) !
1192 ! sfcfsw(:)%upfxc - total sky upward flux at sfc !
1193 ! sfcfsw(:)%dnfxc - total sky downward flux at sfc !
1194 ! sfcfsw(:)%upfx0 - clear sky upward flux at sfc !
1195 ! sfcfsw(:)%dnfx0 - clear sky downward flux at sfc !
1196 ! !
1197 ! 9. lw fluxes at sfc: (defined in 'module_radlw_main') !
1198 ! (sfcflw_type -- derived data type for sfc rad fluxes) !
1199 ! sfcflw(:)%upfxc - total sky upward flux at sfc !
1200 ! sfcflw(:)%dnfxc - total sky downward flux at sfc !
1201 ! sfcflw(:)%dnfx0 - clear sky downward flux at sfc !
1202 ! !
1203 !! optional radiation outputs: !
1204 !! 10. sw flux profiles: (defined in 'module_radsw_main') !
1205 !! (profsw_type -- derived data type for rad vertical profiles) !
1206 !! fswprf(:,:)%upfxc - total sky upward flux !
1207 !! fswprf(:,:)%dnfxc - total sky downward flux !
1208 !! fswprf(:,:)%upfx0 - clear sky upward flux !
1209 !! fswprf(:,:)%dnfx0 - clear sky downward flux !
1210 !! !
1211 !! 11. lw flux profiles: (defined in 'module_radlw_main') !
1212 !! (proflw_type -- derived data type for rad vertical profiles) !
1213 !! flwprf(:,:)%upfxc - total sky upward flux !
1214 !! flwprf(:,:)%dnfxc - total sky downward flux !
1215 !! flwprf(:,:)%upfx0 - clear sky upward flux !
1216 !! flwprf(:,:)%dnfx0 - clear sky downward flux !
1217 !! !
1218 !! 12. sw sfc components: (defined in 'module_radsw_main') !
1219 !! (cmpfsw_type -- derived data type for component sfc fluxes) !
1220 !! scmpsw(:)%uvbfc - total sky downward uv-b flux at sfc !
1221 !! scmpsw(:)%uvbf0 - clear sky downward uv-b flux at sfc !
1222 !! scmpsw(:)%nirbm - total sky sfc downward nir direct flux !
1223 !! scmpsw(:)%nirdf - total sky sfc downward nir diffused flux !
1224 !! scmpsw(:)%visbm - total sky sfc downward uv+vis direct flx !
1225 !! scmpsw(:)%visdf - total sky sfc downward uv+vis diff flux !
1226 ! !
1227 ! external module variables: !
1228 ! ivflip : control flag for in/out vertical indexing !
1229 ! =0 index from toa to surface !
1230 ! =1 index from surface to toa !
1231 ! icmphys : cloud microphysics scheme control flag !
1232 ! =1 zhao/carr/sundqvist microphysics scheme !
1233 ! =2 brad ferrier microphysics scheme !
1234 ! =3 zhao/carr/sundqvist microphysics +pdf cloud !
1235 ! !
1236 ! module variables: !
1237 ! itsfc : =0 use same sfc skin-air/ground temp !
1238 ! =1 use diff sfc skin-air/ground temp (not yet) !
1239 ! !
1240 ! ====================== end of definitions ======================= !
1241 !
1242  implicit none
1243 
1244 ! --- inputs: (for rank>1 arrays, horizontal dimensioned by IX)
1245  integer, intent(in) :: IX,IM, LM, NTRAC, NFXR, me,
1246  & ntoz, ntcw, ncld, ipt, kdt
1247  integer, intent(in) :: icsdsw(im), icsdlw(im), jdate(8)
1248 
1249  logical, intent(in) :: lsswr, lslwr, lssav, lprnt, uni_cld, &
1250  & lmfshal, lmfdeep2
1251 
1252  real (kind=kind_phys), dimension(IX,LM+1), intent(in) :: prsi
1253 
1254  real (kind=kind_phys), dimension(IX,LM), intent(in) :: prsl,
1255  & prslk, tgrs, qgrs, vvl, fcice, frain, rrime, deltaq, cnvw,
1256  & cnvc
1257  real (kind=kind_phys), dimension(IM), intent(in) :: flgmin
1258  real(kind=kind_phys), intent(in) :: sup
1259 
1260  real (kind=kind_phys), dimension(IM), intent(in) :: slmsk,
1261  & xlon, xlat, tsfc, snowd, zorl, hprim, alvsf, alnsf, alvwf,
1262  & alnwf, facsf, facwf, cv, cvt, cvb, fice, tisfc,
1263  & sncovr, snoalb, sinlat, coslat
1264 
1265  real (kind=kind_phys), intent(in) :: solcon, dtlw, dtsw, solhr,
1266  & tracer(ix,lm,ntrac)
1267 
1268  real (kind=kind_phys), dimension(IX,LM),intent(inout):: cldcov
1269 
1270 ! --- outputs: (horizontal dimensioned by IX)
1271  real (kind=kind_phys), dimension(IX,LM),intent(out):: htrsw,htrlw
1272 
1273  real (kind=kind_phys), dimension(IX,4), intent(out) :: dswcmp,
1274  & uswcmp
1275 
1276  real (kind=kind_phys), dimension(IM), intent(out):: tsflw,
1277  & sfalb, semis, coszen, coszdg
1278 
1279  type(topfsw_type), dimension(IM), intent(out) :: topfsw
1280  type(sfcfsw_type), dimension(IM), intent(out) :: sfcfsw
1281 
1282  type(topflw_type), dimension(IM), intent(out) :: topflw
1283  type(sfcflw_type), dimension(IM), intent(out) :: sfcflw
1284 
1285 ! --- variables are for both input and output:
1286  real (kind=kind_phys), intent(inout) :: fluxr(ix,nfxr)
1287 
1288 !! --- optional outputs:
1289  real (kind=kind_phys), dimension(IX,LM,NBDSW), optional,
1290  & intent(out) :: htrswb
1291  real (kind=kind_phys), dimension(IX,LM,NBDLW), optional,
1292  & intent(out) :: htrlwb
1293  real (kind=kind_phys), dimension(ix,lm), optional, &
1294  & intent(out) :: htrlw0
1295  real (kind=kind_phys), dimension(ix,lm), optional, &
1296  & intent(out) :: htrsw0
1297 
1298 ! --- local variables: (horizontal dimensioned by IM)
1299  real (kind=kind_phys), dimension(IM,LM+1+LTP):: plvl, tlvl
1300 
1301  real (kind=kind_phys), dimension(IM,LM+LTP) :: plyr, tlyr, qlyr, &
1302  & olyr, rhly, qstl, vvel, clw, prslk1, tem2da, tem2db, tvly
1303 
1304  real (kind=kind_phys), dimension(IM) :: tsfa, cvt1, cvb1, tem1d, &
1305  & sfcemis, tsfg, tskn
1306 
1307  real (kind=kind_phys), dimension(IM,LM+LTP,NF_CLDS) :: clouds
1308  real (kind=kind_phys), dimension(IM,LM+LTP,NF_VGAS) :: gasvmr
1309  real (kind=kind_phys), dimension(IM, NF_ALBD) :: sfcalb
1310  real (kind=kind_phys), dimension(IM, NSPC1) :: aerodp
1311  real (kind=kind_phys), dimension(IM,LM+LTP,NTRAC) :: tracer1
1312 
1313  real (kind=kind_phys), dimension(IM,LM+LTP,NBDSW,NF_AESW)::faersw
1314  real (kind=kind_phys), dimension(IM,LM+LTP,NBDLW,NF_AELW)::faerlw
1315 
1316  real (kind=kind_phys), dimension(IM,LM+LTP) :: htswc
1317  real (kind=kind_phys), dimension(IM,LM+LTP) :: htlwc
1318 
1319  real (kind=kind_phys), dimension(IM,LM+LTP) :: gcice, grain, grime
1320 
1321 !! --- may be used for optional sw/lw outputs:
1322 !! take out "!!" as needed
1323  real (kind=kind_phys), dimension(IM,LM+LTP) :: htsw0
1324 !! type (profsw_type), dimension(IM,LM+1+LTP) :: fswprf
1325  type(cmpfsw_type), dimension(IM) :: scmpsw
1326  real (kind=kind_phys), dimension(IM,LM+LTP,NBDSW) :: htswb
1327 
1328  real (kind=kind_phys), dimension(IM,LM+LTP) :: htlw0
1329 !! type (proflw_type), dimension(IM,LM+1+LTP) :: flwprf
1330  real (kind=kind_phys), dimension(IM,LM+LTP,NBDLW) :: htlwb
1331 
1332  real (kind=kind_phys) :: raddt, es, qs, delt, tem0d, cldsa(im,5)
1333 
1334  integer :: i, j, k, k1, lv, icec, itop, ibtc, nday, idxday(im), &
1335  & mbota(IM,3), mtopa(IM,3), LP1, nb, LMK, LMP, kd, lla, llb, &
1336  & lya, lyb, kt, kb
1337 
1338 ! --- for debug test use
1339 ! real (kind=kind_phys) :: temlon, temlat, alon, alat
1340 ! integer :: ipt
1341 ! logical :: lprnt1
1342 
1343 !
1344 !===> ... begin here
1345 !
1346  lp1 = lm + 1 ! num of in/out levels
1347 
1348 ! --- ... set local /level/layer indexes corresponding to in/out variables
1349 
1350  lmk = lm + ltp ! num of local layers
1351  lmp = lmk + 1 ! num of local levels
1352 
1353  if ( lextop ) then
1354  if ( ivflip == 1 ) then ! vertical from sfc upward
1355  kd = 0 ! index diff between in/out and local
1356  kt = 1 ! index diff between lyr and upper bound
1357  kb = 0 ! index diff between lyr and lower bound
1358  lla = lmk ! local index at the 2nd level from top
1359  llb = lmp ! local index at toa level
1360  lya = lm ! local index for the 2nd layer from top
1361  lyb = lp1 ! local index for the top layer
1362  else ! vertical from toa downward
1363  kd = 1 ! index diff between in/out and local
1364  kt = 0 ! index diff between lyr and upper bound
1365  kb = 1 ! index diff between lyr and lower bound
1366  lla = 2 ! local index at the 2nd level from top
1367  llb = 1 ! local index at toa level
1368  lya = 2 ! local index for the 2nd layer from top
1369  lyb = 1 ! local index for the top layer
1370  endif ! end if_ivflip_block
1371  else
1372  kd = 0
1373  if ( ivflip == 1 ) then ! vertical from sfc upward
1374  kt = 1 ! index diff between lyr and upper bound
1375  kb = 0 ! index diff between lyr and lower bound
1376  else ! vertical from toa downward
1377  kt = 0 ! index diff between lyr and upper bound
1378  kb = 1 ! index diff between lyr and lower bound
1379  endif ! end if_ivflip_block
1380  endif ! end if_lextop_block
1381 
1382  raddt = min(dtsw, dtlw)
1383 
1384 ! --- ... for debug test
1385 ! alon = 120.0
1386 ! alat = 29.5
1387 ! ipt = 0
1388 ! do i = 1, IM
1389 ! temlon = xlon(i) * 57.29578
1390 ! if (temlon < 0.0) temlon = temlon + 360.0
1391 ! temlat = xlat(i) * 57.29578
1392 ! lprnt1 = abs(temlon-alon) < 1.1 .and. abs(temlat-alat) < 1.1
1393 ! if ( lprnt1 ) then
1394 ! ipt = i
1395 ! exit
1396 ! endif
1397 ! enddo
1398 
1399 ! print *,' in grrad : raddt=',raddt
1400 
1403 
1404  if ( itsfc == 0 ) then ! use same sfc skin-air/ground temp
1405  do i = 1, im
1406  tskn(i) = tsfc(i)
1407  tsfg(i) = tsfc(i)
1408  enddo
1409  else ! use diff sfc skin-air/ground temp
1410  do i = 1, im
1411 !! tskn(i) = ta (i) ! not yet
1412 !! tsfg(i) = tg (i) ! not yet
1413  tskn(i) = tsfc(i)
1414  tsfg(i) = tsfc(i)
1415  enddo
1416  endif
1417 
1419 !
1420 ! convert pressure unit from pa to mb
1421  do k = 1, lm
1422  k1 = k + kd
1423  do i = 1, im
1424 ! plvl(i,k1) = 10.0 * prsi(i,k) ! cb (kpa) to mb (hpa)
1425 ! plyr(i,k1) = 10.0 * prsl(i,k) ! cb (kpa) to mb (hpa)
1426  plvl(i,k1) = 0.01 * prsi(i,k) ! pa to mb (hpa)
1427  plyr(i,k1) = 0.01 * prsl(i,k) ! pa to mb (hpa)
1428  tlyr(i,k1) = tgrs(i,k)
1429  prslk1(i,k1) = prslk(i,k)
1430 
1432 ! es = min( prsl(i,k), 0.001 * fpvs( tgrs(i,k) ) ) ! fpvs in pa
1433  es = min( prsl(i,k), fpvs( tgrs(i,k) ) ) ! fpvs and prsl in pa
1434  qs = max( qmin, eps * es / (prsl(i,k) + epsm1*es) )
1435  rhly(i,k1) = max( 0.0, min( 1.0, max(qmin, qgrs(i,k))/qs ) )
1436  qstl(i,k1) = qs
1437  enddo
1438  enddo
1439  do j = 1, ntrac
1440  do k = 1, lm
1441  k1 = k + kd
1442  do i = 1, im
1443  tracer1(i,k1,j) = tracer(i,k,j)
1444  enddo
1445  enddo
1446  enddo
1447 
1448  do i = 1, im
1449 ! plvl(i,LP1+kd) = 10.0 * prsi(i,LP1) ! cb (kpa) to mb (hpa
1450  plvl(i,lp1+kd) = 0.01 * prsi(i,lp1) ! pa to mb (hpa)
1451  enddo
1452 
1453  if ( lextop ) then ! values for extra top layer
1454  do i = 1, im
1455  plvl(i,llb) = prsmin
1456  if ( plvl(i,lla) <= prsmin ) plvl(i,lla) = 2.0*prsmin
1457  plyr(i,lyb) = 0.5 * plvl(i,lla)
1458  tlyr(i,lyb) = tlyr(i,lya)
1459 ! prslk1(i,lyb) = (plyr(i,lyb)*0.001) ** rocp ! plyr in hPa
1460  prslk1(i,lyb) = (plyr(i,lyb)*0.00001) ** rocp ! plyr in Pa
1461 
1462  rhly(i,lyb) = rhly(i,lya)
1463  qstl(i,lyb) = qstl(i,lya)
1464  enddo
1465 
1466  do j = 1, ntrac
1467  do i = 1, im
1468 ! --- note: may need to take care the top layer amount
1469  tracer1(i,lyb,j) = tracer1(i,lya,j)
1470  enddo
1471  enddo
1472  endif
1473 
1474 ! --- ... extra variables needed for ferrier's microphysics
1475 
1476  if (icmphys == 2) then
1477  do k = 1, lm
1478  k1 = k + kd
1479 
1480  do i = 1, im
1481  gcice(i,k1) = fcice(i,k)
1482  grain(i,k1) = frain(i,k)
1483  grime(i,k1) = rrime(i,k)
1484  enddo
1485  enddo
1486 
1487  if ( lextop ) then
1488  do i = 1, im
1489  gcice(i,lyb) = fcice(i,lya)
1490  grain(i,lyb) = frain(i,lya)
1491  grime(i,lyb) = rrime(i,lya)
1492  enddo
1493  endif
1494  endif ! if_icmphys
1495 
1497 
1498  if (ntoz > 0) then ! interactive ozone generation
1499 
1500  do k = 1, lmk
1501  do i = 1, im
1502  olyr(i,k) = max( qmin, tracer1(i,k,ntoz) )
1503  enddo
1504  enddo
1505 
1506  else ! climatological ozone
1507 
1508 ! print *,' in grrad : calling getozn'
1509  call getozn &
1510 ! --- inputs:
1511  & ( prslk1,xlat, &
1512  & im, lmk, &
1513 ! --- outputs:
1514  & olyr &
1515  & )
1516 
1517  endif ! end_if_ntoz
1518 
1520 
1521  call coszmn &
1522 ! --- inputs:
1523  & ( xlon,sinlat,coslat,solhr, im, me, &
1524 ! --- outputs:
1525  & coszen, coszdg &
1526  & )
1527 
1539 
1540 ! --- ... set up non-prognostic gas volume mixing ratioes
1541 
1542  call getgases &
1543 ! --- inputs:
1544  & ( plvl, xlon, xlat, &
1545  & im, lmk, &
1546 ! --- outputs:
1547  & gasvmr &
1548  & )
1549 
1551 
1552  do k = 2, lmk
1553  do i = 1, im
1554  tem2da(i,k) = log( plyr(i,k) )
1555  tem2db(i,k) = log( plvl(i,k) )
1556  enddo
1557  enddo
1558 
1559  if (ivflip == 0) then ! input data from toa to sfc
1560 
1561  do i = 1, im
1562  tem1d(i) = qme6
1563  tem2da(i,1) = log( plyr(i,1) )
1564  tem2db(i,1) = 1.0
1565  tsfa(i) = tlyr(i,lmk) ! sfc layer air temp
1566  tlvl(i,1) = tlyr(i,1)
1567  tlvl(i,lmp) = tskn(i)
1568  enddo
1569 
1570  do k = 1, lm
1571  k1 = k + kd
1572 
1573  do i = 1, im
1574  qlyr(i,k1) = max( tem1d(i), qgrs(i,k) )
1575  tem1d(i) = min( qme5, qlyr(i,k1) )
1576  tvly(i,k1) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k1)) ! virtual T (K)
1577  enddo
1578  enddo
1579 
1580  if ( lextop ) then
1581  do i = 1, im
1582  qlyr(i,lyb) = qlyr(i,lya)
1583  tvly(i,lyb) = tvly(i,lya)
1584  enddo
1585  endif
1586 
1587  do k = 2, lmk
1588  do i = 1, im
1589  tlvl(i,k) = tlyr(i,k) + (tlyr(i,k-1) - tlyr(i,k)) &
1590  & * (tem2db(i,k) - tem2da(i,k)) &
1591  & / (tem2da(i,k-1) - tem2da(i,k))
1592  enddo
1593  enddo
1594 
1595  else ! input data from sfc to toa
1596 
1597  do i = 1, im
1598  tem1d(i) = qme6
1599  tem2da(i,1) = log( plyr(i,1) )
1600  tem2db(i,1) = log( plvl(i,1) )
1601  tsfa(i) = tlyr(i,1) ! sfc layer air temp
1602  tlvl(i,1) = tskn(i)
1603  tlvl(i,lmp) = tlyr(i,lmk)
1604  enddo
1605 
1606  do k = lm, 1, -1
1607  do i = 1, im
1608  qlyr(i,k) = max( tem1d(i), qgrs(i,k) )
1609  tem1d(i) = min( qme5, qlyr(i,k) )
1610  tvly(i,k) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k)) ! virtual T (K)
1611  enddo
1612  enddo
1613 
1614  if ( lextop ) then
1615  do i = 1, im
1616  qlyr(i,lyb) = qlyr(i,lya)
1617  tvly(i,lyb) = tvly(i,lya)
1618  enddo
1619  endif
1620 
1621  do k = 1, lmk-1
1622  do i = 1, im
1623  tlvl(i,k+1) = tlyr(i,k) + (tlyr(i,k+1) - tlyr(i,k)) &
1624  & * (tem2db(i,k+1) - tem2da(i,k)) &
1625  & / (tem2da(i,k+1) - tem2da(i,k))
1626  enddo
1627  enddo
1628 
1629  endif ! end_if_ivflip
1630 
1632 
1633  nday = 0
1634  do i = 1, im
1635  if (coszen(i) >= 0.0001) then
1636  nday = nday + 1
1637  idxday(nday) = i
1638  endif
1639  enddo
1640 
1641 ! write(0,*)' plvl=',plvl(ipt,1:65)
1642 ! write(0,*)' plyr=',plyr(ipt,1:64)
1643 ! write(0,*)' tlyr=',tlyr(ipt,1:64)
1644 ! write(0,*)' tlvl=',tlvl(ipt,1:65)
1645 ! write(0,*)' qlyr=',qlyr(ipt,1:10)*1000
1646 
1649 
1650 !check print *,' in grrad : calling setaer '
1651 
1652  call setaer &
1653 ! --- inputs:
1654  & ( plvl,plyr,prslk1,tvly,rhly,slmsk,tracer1,xlon,xlat, &
1655  & im,lmk,lmp, lsswr,lslwr, &
1656 ! --- outputs:
1657  & faersw,faerlw,aerodp &
1658  & )
1659 
1671 
1672 ! --- ... obtain cloud information for radiation calculations
1673 
1674  if (ntcw > 0) then ! prognostic cloud scheme
1675 
1676  do k = 1, lmk
1677  do i = 1, im
1678  clw(i,k) = 0.0
1679  enddo
1680 
1681  do j = 1, ncld
1682  lv = ntcw + j - 1
1683  do i = 1, im
1684  clw(i,k) = clw(i,k) + tracer1(i,k,lv) ! cloud condensate amount
1685  enddo
1686  enddo
1687  enddo
1688 
1689  do k = 1, lmk
1690  do i = 1, im
1691  if ( clw(i,k) < epsq ) clw(i,k) = 0.0
1692  enddo
1693  enddo
1694 !
1695 ! --- add suspended convective cloud water to grid-scale cloud water
1696 ! only for cloud fraction & radiation computation
1697 ! it is to enhance cloudiness due to suspended convec cloud water
1698 ! for zhao/moorthi's (icmphys=1) &
1699 ! ferrier's (icmphys=2) microphysics schemes
1700 !
1701  if (icmphys == 1 .or. icmphys == 2) then
1702  do k = 1, lmk
1703  do i = 1, im
1704  clw(i,k) = clw(i,k) + cnvw(i,k)
1705  enddo
1706  enddo
1707  endif
1708 !
1709 
1710  if (icmphys == 1) then ! zhao/moorthi's prognostic cloud scheme
1711 
1712  call progcld1 &
1713 ! --- inputs:
1714  & ( plyr,plvl,tlyr,tvly,qlyr,qstl,rhly,clw, &
1715  & xlat,xlon,slmsk, im, lmk, lmp, &
1716  & uni_cld, lmfshal, lmfdeep2, cldcov(1:im,1:lm), &
1717 ! --- outputs:
1718  & clouds,cldsa,mtopa,mbota &
1719  & )
1720 
1721  elseif (icmphys == 2) then ! ferrier's microphysics
1722 
1723 ! print *,' in grrad : calling progcld2'
1724  call progcld2 &
1725 ! --- inputs:
1726  & ( plyr,plvl,tlyr,tvly,qlyr,qstl,rhly,clw, &
1727  & xlat,xlon,slmsk, gcice,grain,grime,flgmin, &
1728  & im, lmk, lmp, lmfshal, lmfdeep2, &
1729 ! --- outputs:
1730  & clouds,cldsa,mtopa,mbota &
1731  & )
1732 !
1733  elseif(icmphys == 3) then ! zhao/moorthi's prognostic cloud+pdfcld
1734 
1735  call progcld3 &
1736 ! --- inputs:
1737  & ( plyr,plvl,tlyr,tvly,qlyr,qstl,rhly,clw,cnvw,cnvc, &
1738  & xlat,xlon,slmsk, &
1739  & im, lmk, lmp, &
1740  & deltaq, sup,kdt,me, &
1741 ! --- outputs:
1742  & clouds,cldsa,mtopa,mbota &
1743  & )
1744 
1745  endif ! end if_icmphys
1746 
1747  else ! diagnostic cloud scheme
1748 
1749  do i = 1, im
1750 ! cvt1(i) = 10.0 * cvt(i)
1751 ! cvb1(i) = 10.0 * cvb(i)
1752  cvt1(i) = 0.01 * cvt(i)
1753  cvb1(i) = 0.01 * cvb(i)
1754 
1755  enddo
1756 
1757  do k = 1, lm
1758  k1 = k + kd
1759 
1760  do i = 1, im
1761 ! vvel(i,k1) = 10.0 * vvl(i,k)
1762  vvel(i,k1) = 0.01 * vvl(i,k)
1763  enddo
1764  enddo
1765 
1766  if ( lextop ) then
1767  do i = 1, im
1768  vvel(i,lyb) = vvel(i,lya)
1769  enddo
1770  endif
1771 
1772 ! --- compute diagnostic cloud related quantities
1773 
1774  call diagcld1 &
1775 ! --- inputs:
1776  & ( plyr,plvl,tlyr,rhly,vvel,cv,cvt1,cvb1, &
1777  & xlat,xlon,slmsk, &
1778  & im, lmk, lmp, &
1779 ! --- outputs:
1780  & clouds,cldsa,mtopa,mbota &
1781  & )
1782 
1783  endif ! end_if_ntcw
1784 
1785 ! --- ... start radiation calculations
1786 ! remember to set heating rate unit to k/sec!
1787 
1788  if (lsswr) then
1789 
1792 
1793 
1794  call setalb &
1795 ! --- inputs:
1796  & ( slmsk,snowd,sncovr,snoalb,zorl,coszen,tsfg,tsfa,hprim, &
1797  & alvsf,alnsf,alvwf,alnwf,facsf,facwf,fice,tisfc, &
1798  & im, &
1799 ! --- outputs:
1800  & sfcalb &
1801  & )
1802 
1804 
1805  do i = 1, im
1806  sfalb(i) = max(0.01, 0.5 * (sfcalb(i,2) + sfcalb(i,4)))
1807  enddo
1808 
1809  if (nday > 0) then
1810 
1812 ! print *,' in grrad : calling swrad'
1813 
1814  if ( present(htrswb) .and. present(htrsw0)) then
1815 
1816  call swrad &
1817 ! --- inputs:
1818  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
1819  & clouds,icsdsw,faersw,sfcalb, &
1820  & coszen,solcon, nday,idxday, &
1821  & im, lmk, lmp, lprnt, &
1822 ! --- outputs:
1823  & htswc,topfsw,sfcfsw &
1824 !! --- optional:
1825 !! &, HSW0=htsw0,FLXPRF=fswprf &
1826  &, hsw0=htsw0,hswb=htswb,fdncmp=scmpsw &
1827  & )
1828 
1829  do k = 1, lm
1830  k1 = k + kd
1831 
1832  do j = 1, nbdsw
1833  do i = 1, im
1834  htrswb(i,k,j) = htswb(i,k1,j)
1835  enddo
1836  enddo
1837  enddo
1838 
1839  else if ( present(htrswb) .and. .not. present(htrsw0) ) then
1840 
1841  call swrad &
1842 ! --- inputs:
1843  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
1844  & clouds,icsdsw,faersw,sfcalb, &
1845  & coszen,solcon, nday,idxday, &
1846  & im, lmk, lmp, lprnt, &
1847 ! --- outputs:
1848  & htswc,topfsw,sfcfsw &
1849 !! --- optional:
1850 !! &, hsw0=htsw0,flxprf=fswprf &
1851  &, hswb=htswb,fdncmp=scmpsw &
1852  & )
1853 
1854  do k = 1, lm
1855  k1 = k + kd
1856  do j = 1, nbdsw
1857  do i = 1, im
1858  htrswb(i,k,j) = htswb(i,k1,j)
1859  enddo
1860  enddo
1861  enddo
1862 
1863  else if ( present(htrsw0) .and. .not. present(htrswb) ) then
1864 
1865  call swrad &
1866 ! --- inputs:
1867  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
1868  & clouds,icsdsw,faersw,sfcalb, &
1869  & coszen,solcon, nday,idxday, &
1870  & im, lmk, lmp, lprnt, &
1871 ! --- outputs:
1872  & htswc,topfsw,sfcfsw &
1873 !! --- optional:
1874 !! &, hsw0=htsw0,flxprf=fswprf &
1875  &, hsw0=htsw0,fdncmp=scmpsw &
1876  & )
1877 
1878  else
1879 
1880  call swrad &
1881 ! --- inputs:
1882  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
1883  & clouds,icsdsw,faersw,sfcalb, &
1884  & coszen,solcon, nday,idxday, &
1885  & im, lmk, lmp, lprnt, &
1886 ! --- outputs:
1887  & htswc,topfsw,sfcfsw &
1888 !! --- optional:
1889 !! &, HSW0=htsw0,FLXPRF=fswprf,HSWB=htswb &
1890  &, fdncmp=scmpsw &
1891  & )
1892 
1893  endif
1894 
1895  do k = 1, lm
1896  k1 = k + kd
1897 
1898  do i = 1, im
1899  htrsw(i,k) = htswc(i,k1)
1900  enddo
1901  enddo
1902  if (present(htrsw0)) then
1903  do k = 1, lm
1904  k1 = k + kd
1905  do i = 1, im
1906  htrsw0(i,k) = htsw0(i,k1)
1907  enddo
1908  enddo
1909  endif
1910 
1911 ! --- surface down and up spectral component fluxes
1912 
1913  do i = 1, im
1914  dswcmp(i,1) = scmpsw(i)%nirbm
1915  dswcmp(i,2) = scmpsw(i)%nirdf
1916  dswcmp(i,3) = scmpsw(i)%visbm
1917  dswcmp(i,4) = scmpsw(i)%visdf
1918 
1919  uswcmp(i,1) = scmpsw(i)%nirbm * sfcalb(i,1)
1920  uswcmp(i,2) = scmpsw(i)%nirdf * sfcalb(i,2)
1921  uswcmp(i,3) = scmpsw(i)%visbm * sfcalb(i,3)
1922  uswcmp(i,4) = scmpsw(i)%visdf * sfcalb(i,4)
1923  enddo
1924 
1925  else ! if_nday_block
1926 
1927  do k = 1, lm
1928  do i = 1, im
1929  htrsw(i,k) = 0.0
1930  enddo
1931  enddo
1932 
1933  sfcfsw = sfcfsw_type( 0.0, 0.0, 0.0, 0.0 )
1934  topfsw = topfsw_type( 0.0, 0.0, 0.0 )
1935  scmpsw = cmpfsw_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 )
1936 
1937  do k = 1, 4
1938  do i = 1, im
1939  dswcmp(i,k) = 0.0
1940  uswcmp(i,k) = 0.0
1941  enddo
1942  enddo
1943 
1944 !! --- optional:
1945 !! fswprf= profsw_type( 0.0, 0.0, 0.0, 0.0 )
1946 
1947  if ( present(htrswb) ) then
1948  do j = 1, nbdsw
1949  do k = 1, lm
1950  do i = 1, im
1951  htrswb(i,k,j) = 0.0
1952  enddo
1953  enddo
1954  enddo
1955  endif
1956  if ( present(htrsw0) ) then
1957  do k = 1, lm
1958  do i = 1, im
1959  htrsw0(i,k) = 0.0
1960  enddo
1961  enddo
1962  endif
1963 
1964  endif ! end_if_nday
1965 
1966  endif ! end_if_lsswr
1967 
1968 ! write(0,*)' htrsw=',htrsw(ipt,1:64)*86400
1969  if (lslwr) then
1970 
1973 
1974  call setemis &
1975 ! --- inputs:
1976  & ( xlon,xlat,slmsk,snowd,sncovr,zorl,tsfg,tsfa,hprim, &
1977  & im, &
1978 ! --- outputs:
1979  & sfcemis &
1980  & )
1981 
1983 ! print *,' in grrad : calling lwrad'
1984 
1985  if ( present(htrlwb) .and. present(htrlw0) ) then
1986 
1987  call lwrad &
1988 ! --- inputs:
1989  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
1990  & clouds,icsdlw,faerlw,sfcemis,tsfg, &
1991  & im, lmk, lmp, lprnt, &
1992 ! --- outputs:
1993  & htlwc,topflw,sfcflw &
1994 !! --- optional:
1995 !! &, HLW0=htlw0,FLXPRF=flwprf &
1996  &, hlw0=htlw0 &
1997  &, hlwb=htlwb &
1998  & )
1999 
2000  do k = 1, lm
2001  k1 = k + kd
2002 
2003  do j = 1, nbdlw
2004  do i = 1, im
2005  htrlwb(i,k,j) = htlwb(i,k1,j)
2006  enddo
2007  enddo
2008  enddo
2009 
2010  else if ( present(htrlwb) .and. .not. present(htrlw0) ) then
2011 
2012  call lwrad &
2013 ! --- inputs:
2014  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
2015  & clouds,icsdlw,faerlw,sfcemis,tsfg, &
2016  & im, lmk, lmp, lprnt, &
2017 ! --- outputs:
2018  & htlwc,topflw,sfcflw &
2019 !! --- optional:
2020 !! &, hlw0=htlw0,flxprf=flwprf &
2021  &, hlwb=htlwb &
2022  & )
2023 
2024  do k = 1, lm
2025  k1 = k + kd
2026 
2027  do j = 1, nbdlw
2028  do i = 1, im
2029  htrlwb(i,k,j) = htlwb(i,k1,j)
2030  enddo
2031  enddo
2032  enddo
2033  else if ( present(htrlw0) .and. .not. present(htrlwb) ) then
2034 
2035  !print *,'call lwrad saving clear sky component'
2036  call lwrad &
2037 ! --- inputs:
2038  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
2039  & clouds,icsdlw,faerlw,sfcemis,tsfg, &
2040  & im, lmk, lmp, lprnt, &
2041 ! --- outputs:
2042  & htlwc,topflw,sfcflw &
2043 !! --- optional:
2044 !! &, hlw0=htlw0,flxprf=flwprf &
2045  &, hlw0=htlw0 &
2046  & )
2047 
2048  else
2049 
2050  call lwrad &
2051 ! --- inputs:
2052  & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
2053  & clouds,icsdlw,faerlw,sfcemis,tsfg, &
2054  & im, lmk, lmp, lprnt, &
2055 ! --- outputs:
2056  & htlwc,topflw,sfcflw &
2057 !! --- optional:
2058 !! &, HLW0=htlw0,FLXPRF=flwprf,HLWB=htlwb &
2059  & )
2060 
2061  endif
2062 
2063  do i = 1, im
2064  semis(i) = sfcemis(i)
2065 ! --- save surface air temp for diurnal adjustment at model t-steps
2066  tsflw(i) = tsfa(i)
2067  enddo
2068 
2069  do k = 1, lm
2070  k1 = k + kd
2071 
2072  do i = 1, im
2073  htrlw(i,k) = htlwc(i,k1)
2074  enddo
2075  enddo
2076 
2077  if (present(htrlw0)) then
2078  do k = 1, lm
2079  k1 = k + kd
2080  do i = 1, im
2081  htrlw0(i,k) = htlw0(i,k1)
2082  enddo
2083  enddo
2084  endif
2085 
2086  endif ! end_if_lslwr
2087 
2089 ! --- ... collect the fluxr data for wrtsfc
2090 
2091  if (lssav) then
2092 
2093  if ( lsswr ) then
2094  do i = 1, im
2095  fluxr(i,34) = fluxr(i,34) + dtsw*aerodp(i,1) ! total aod at 550nm
2096  fluxr(i,35) = fluxr(i,35) + dtsw*aerodp(i,2) ! DU aod at 550nm
2097  fluxr(i,36) = fluxr(i,36) + dtsw*aerodp(i,3) ! BC aod at 550nm
2098  fluxr(i,37) = fluxr(i,37) + dtsw*aerodp(i,4) ! OC aod at 550nm
2099  fluxr(i,38) = fluxr(i,38) + dtsw*aerodp(i,5) ! SU aod at 550nm
2100  fluxr(i,39) = fluxr(i,39) + dtsw*aerodp(i,6) ! SS aod at 550nm
2101  enddo
2102  endif
2103 ! --- save lw toa and sfc fluxes
2104 
2105  if (lslwr) then
2106  do i = 1, im
2107 ! --- lw total-sky fluxes
2108  fluxr(i,1 ) = fluxr(i,1 ) + dtlw * topflw(i)%upfxc ! total sky top lw up
2109  fluxr(i,19) = fluxr(i,19) + dtlw * sfcflw(i)%dnfxc ! total sky sfc lw dn
2110  fluxr(i,20) = fluxr(i,20) + dtlw * sfcflw(i)%upfxc ! total sky sfc lw up
2111 ! --- lw clear-sky fluxes
2112  fluxr(i,28) = fluxr(i,28) + dtlw * topflw(i)%upfx0 ! clear sky top lw up
2113  fluxr(i,30) = fluxr(i,30) + dtlw * sfcflw(i)%dnfx0 ! clear sky sfc lw dn
2114  fluxr(i,33) = fluxr(i,33) + dtlw * sfcflw(i)%upfx0 ! clear sky sfc lw up
2115  enddo
2116  endif
2117 
2118 ! --- save sw toa and sfc fluxes with proper diurnal sw wgt. coszen=mean cosz over daylight
2119 ! part of sw calling interval, while coszdg= mean cosz over entire interval
2120 
2121  if (lsswr) then
2122  do i = 1, im
2123  if (coszen(i) > 0.) then
2124 ! --- sw total-sky fluxes
2125 ! -------------------
2126  tem0d = dtsw * coszdg(i) / coszen(i)
2127  fluxr(i,2 ) = fluxr(i,2) + topfsw(i)%upfxc * tem0d ! total sky top sw up
2128  fluxr(i,3 ) = fluxr(i,3) + sfcfsw(i)%upfxc * tem0d ! total sky sfc sw up
2129  fluxr(i,4 ) = fluxr(i,4) + sfcfsw(i)%dnfxc * tem0d ! total sky sfc sw dn
2130 ! --- sw uv-b fluxes
2131 ! --------------
2132  fluxr(i,21) = fluxr(i,21) + scmpsw(i)%uvbfc * tem0d ! total sky uv-b sw dn
2133  fluxr(i,22) = fluxr(i,22) + scmpsw(i)%uvbf0 * tem0d ! clear sky uv-b sw dn
2134 ! --- sw toa incoming fluxes
2135 ! ----------------------
2136  fluxr(i,23) = fluxr(i,23) + topfsw(i)%dnfxc * tem0d ! top sw dn
2137 ! --- sw sfc flux components
2138 ! ----------------------
2139  fluxr(i,24) = fluxr(i,24) + scmpsw(i)%visbm * tem0d ! uv/vis beam sw dn
2140  fluxr(i,25) = fluxr(i,25) + scmpsw(i)%visdf * tem0d ! uv/vis diff sw dn
2141  fluxr(i,26) = fluxr(i,26) + scmpsw(i)%nirbm * tem0d ! nir beam sw dn
2142  fluxr(i,27) = fluxr(i,27) + scmpsw(i)%nirdf * tem0d ! nir diff sw dn
2143 ! --- sw clear-sky fluxes
2144 ! -------------------
2145  fluxr(i,29) = fluxr(i,29) + topfsw(i)%upfx0 * tem0d ! clear sky top sw up
2146  fluxr(i,31) = fluxr(i,31) + sfcfsw(i)%upfx0 * tem0d ! clear sky sfc sw up
2147  fluxr(i,32) = fluxr(i,32) + sfcfsw(i)%dnfx0 * tem0d ! clear sky sfc sw dn
2148  endif
2149  enddo
2150  endif
2151 
2152 ! --- save total and boundary layer clouds
2153 
2154  if (lsswr .or. lslwr) then
2155  do i = 1, im
2156  fluxr(i,17) = fluxr(i,17) + raddt * cldsa(i,4)
2157  fluxr(i,18) = fluxr(i,18) + raddt * cldsa(i,5)
2158  enddo
2159 
2160 ! --- save cld frac,toplyr,botlyr and top temp, note that the order
2161 ! of h,m,l cloud is reversed for the fluxr output.
2162 ! --- save interface pressure (pa) of top/bot
2163 
2164  do j = 1, 3
2165  do i = 1, im
2166  tem0d = raddt * cldsa(i,j)
2167  itop = mtopa(i,j) - kd
2168  ibtc = mbota(i,j) - kd
2169  fluxr(i, 8-j) = fluxr(i, 8-j) + tem0d
2170  fluxr(i,11-j) = fluxr(i,11-j) + tem0d * prsi(i,itop+kt)
2171  fluxr(i,14-j) = fluxr(i,14-j) + tem0d * prsi(i,ibtc+kb)
2172  fluxr(i,17-j) = fluxr(i,17-j) + tem0d * tgrs(i,itop)
2173  enddo
2174  enddo
2175  endif
2176 
2177  if (.not. uni_cld) then
2178  do k = 1, lm
2179  k1 = k + kd
2180  do i = 1, im
2181  cldcov(i,k) = clouds(i,k1,1)
2182  enddo
2183  enddo
2184  endif
2185 
2186 ! --- save optional vertically integrated aerosol optical depth at
2187 ! wavelenth of 550nm aerodp(:,1), and other optional aod for
2188 ! individual species aerodp(:,2:NSPC1)
2189 
2190 ! if ( laswflg ) then
2191 ! if ( NFXR > 33 ) then
2192 ! do i = 1, IM
2193 ! fluxr(i,34) = fluxr(i,34) + dtsw*aerodp(i,1) ! total aod at 550nm (all species)
2194 ! enddo
2195 
2196 ! if ( lspcodp ) then
2197 ! do j = 2, NSPC1
2198 ! k = 33 + j
2199 
2200 ! do i = 1, IM
2201 ! fluxr(i,k) = fluxr(i,k) + dtsw*aerodp(i,j) ! aod at 550nm for indiv species
2202 ! enddo
2203 ! enddo
2204 ! endif ! end_if_lspcodp
2205 ! else
2206 ! print *,' !Error! Need to increase array fluxr size NFXR ',&
2207 ! & ' to be able to output aerosol optical depth'
2208 ! stop
2209 ! endif ! end_if_nfxr
2210 ! endif ! end_if_laswflg
2211 
2212  endif ! end_if_lssav
2213 !
2214  return
2215 !...................................
2216  end subroutine grrad
2217 !-----------------------------------
2218 
2219 
2220 !
2222 !........................................!
2223  end module module_radiation_driver !
2224 !========================================!
integer, save iovrsw
cloud overlapping control flag for SW
Definition: physparam.f:192
integer, save isubcsw
sub-column cloud approx flag in SW radiation
Definition: physparam.f:224
real(kind=kind_phys) epsq
EPSQ=1.0e-12.
Definition: grrad.f:347
integer, save iaerflg
aerosol effect control flag
Definition: physparam.f:148
real(kind=kind_phys) qmin
QMIN=1.0e-10.
Definition: grrad.f:341
logical, save lcnorm
in-cld condensate control flag
Definition: physparam.f:198
logical, save lcrick
eliminating CRICK control flag
Definition: physparam.f:196
real(kind=kind_phys) qme5
QME5=1.0e-7.
Definition: grrad.f:343
real, parameter prsmin
toa pressure minimum value in mb (hPa)
Definition: grrad.f:353
integer, save ialbflg
surface albedo scheme control flag
Definition: physparam.f:209
logical loz1st
first-time climatological ozone data read flag
Definition: grrad.f:361
integer, save ico2flg
co2 data source control flag
Definition: physparam.f:165
integer, save icmphys
cloud micorphysics scheme control flag
Definition: physparam.f:190
integer, save icldflg
cloud optical property scheme control flag
Definition: physparam.f:188
integer, save iemsflg
surface emissivity scheme control flag
Definition: physparam.f:211
elemental real(krealfp) function, public fpvs(t)
Definition: funcphys.f:703
character(40), parameter vtagrad
Definition: grrad.f:333
logical, save lnoprec
precip effect on radiation flag (Ferrier microphysics)
Definition: physparam.f:200
integer, parameter ltp
optional extra top layer on top of low ceiling models LTP=0: no extra top layer ...
Definition: grrad.f:365
subroutine, public radinit(si, NLAY, me)
This subroutine is the initialization of radiation calculations.
Definition: grrad.f:386
integer, save ioznflg
ozone data source control flag
Definition: physparam.f:169
integer, save ictmflg
external data time/date control flag
Definition: physparam.f:167
real(kind=kind_phys) qme6
QME6=1.0e-7.
Definition: grrad.f:345
integer, save isolar
solar constant scheme control flag
Definition: physparam.f:134
integer itsfc
control flag for lw sfc air/ground interface temp setting
Definition: grrad.f:356
subroutine, public radupdate(idate, jdate, deltsw, deltim, lsswr, me, slag, sdec, cdec, solcon)
This subroutine calls many update subroutines to check and update radiation required but time varying...
Definition: grrad.f:627
integer, save isubclw
sub-column cloud approx flag in LW radiation
Definition: physparam.f:226
integer, save ivflip
vertical profile indexing flag
Definition: physparam.f:222
logical, parameter lextop
Definition: grrad.f:367
integer, save iovrlw
cloud overlapping control flag for LW
Definition: physparam.f:194
subroutine, public grrad(prsi, prsl, prslk, tgrs, qgrs, tracer, vvl, slmsk, xlon, xlat, tsfc, snowd, sncovr, snoalb, zorl, hprim, alvsf, alnsf, alvwf, alnwf, facsf, facwf, fice, tisfc, sinlat, coslat, solhr, jdate, solcon, cv, cvt, cvb, fcice, frain, rrime, flgmin, icsdsw, icsdlw, ntcw, ncld, ntoz, NTRAC, NFXR, dtlw, dtsw, lsswr, lslwr, lssav, uni_cld, lmfshal, lmfdeep2, IX, IM, LM, me, lprnt, ipt, kdt, deltaq, sup, cnvw, cnvc, htrsw, topfsw, sfcfsw, dswcmp, uswcmp, sfalb, coszen, coszdg, htrlw, topflw, sfcflw, tsflw, semis, cldcov, fluxr , htrlw0, htrsw0, htrswb, htrlwb )
This subroutine is the driver of radiation calculation subroutines. It sets up profile variables for ...
Definition: grrad.f:957