Radiation Scheme in CCPP
radiation_gases.f
Go to the documentation of this file.
1 ! ========================================================== !!!!!
2 ! 'module_radiation_gases' description !!!!!
3 ! ========================================================== !!!!!
4 ! !
5 ! set up ozone climatological profiles and other constant gas !
6 ! profiles, such as co2, ch4, n2o, o2, and those of cfc gases. All !
7 ! data are entered as mixing ratio by volume, except ozone which is !
8 ! mass mixing ratio (g/g). !
9 ! !
10 ! in the module, the externally callabe subroutines are : !
11 ! !
12 ! 'gas_init' -- initialization !
13 ! input: !
14 ! ( me ) !
15 ! output: !
16 ! ( none ) !
17 ! !
18 ! 'gas_update' -- read in data and update with time !
19 ! input: !
20 ! ( iyear, imon, iday, ihour, loz1st, ldoco2, me ) !
21 ! output: !
22 ! ( none ) !
23 ! !
24 ! 'getozn' -- setup climatological ozone profile !
25 ! input: !
26 ! ( prslk,xlat, !
27 ! IMAX, LM ) !
28 ! output: !
29 ! ( o3mmr ) !
30 ! !
31 ! 'getgases' -- setup constant gas profiles for LW and SW !
32 ! input: !
33 ! ( plvl, xlon, xlat, !
34 ! IMAX, LMAX ) !
35 ! output: !
36 ! ( gasdat ) !
37 ! !
38 ! external modules referenced: !
39 ! 'module machine' in 'machine.f' !
40 ! 'module funcphys' in 'funcphys.f' !
41 ! 'module physcons' in 'physcons.f !
42 ! 'module module_iounitdef' in 'iounitdef.f' !
43 ! !
44 ! unit used for radiative active gases: !
45 ! ozone : mass mixing ratio (g/g) !
46 ! co2 : volume mixing ratio (p/p) !
47 ! n2o : volume mixing ratio (p/p) !
48 ! ch4 : volume mixing ratio (p/p) !
49 ! o2 : volume mixing ratio (p/p) !
50 ! co : volume mixing ratio (p/p) !
51 ! cfc11 : volume mixing ratio (p/p) !
52 ! cfc12 : volume mixing ratio (p/p) !
53 ! cfc22 : volume mixing ratio (p/p) !
54 ! ccl4 : volume mixing ratio (p/p) !
55 ! cfc113: volume mixing ratio (p/p) !
56 ! !
57 ! !
58 ! program history: !
59 ! may 2003 - y-t hou create rad_module.f that collectively !
60 ! combines several radiation computation supporting !
61 ! programs into fortran 90 module structure (gases !
62 ! and aerosols, etc.) !
63 ! apr 2004 - y-t hou modified to add astronomy and surface !
64 ! module components. !
65 ! feb 2005 - y-t hou rewrite the component modules into !
66 ! separate individule modules for thier corresponding !
67 ! tasks. here as radiation_gases.f !
68 ! mar 2006 - y-t hou add initialization subroutine to co2 and !
69 ! other gases. historical 2-d co2 data are added. !
70 ! sep 2008 - y-t hou add parameter ictm to control the input !
71 ! data time at the model initial condition. !
72 ! oct 2008 - y-t hou modify the initialization code to add the !
73 ! option of superimposing climatology seasonal cycle !
74 ! to the initial condition data (currently co2 only) !
75 ! nov 2008 - y-t hou fix bugs in superimposing climatology !
76 ! seasonal cycle calculations !
77 ! aug 2011 - y-t hou fix a bug in subr getgases doing vertical !
78 ! co2 mapping. (for iflip=0 case, not affact opr). !
79 ! aug 2012 - y-t hou modified subr getozn. moved the if-first !
80 ! block to subr gas_init to ensure threading safe in !
81 ! climatology ozone applications. (not affect gfs) !
82 ! also changed the initialization subr into two parts:!
83 ! 'gas_init' is called at the start of run to set up !
84 ! module parameters; and 'gas_update' is called within!
85 ! the time loop to check and update data sets. defined!
86 ! the climatology ozone parameters k1oz,k2oz,facoz as !
87 ! module variables and are set in subr 'gas_update' !
88 ! nov 2012 - y-t hou modified control parameters thru module !
89 ! 'physparam'. !
90 ! jan 2013 - z. janjic/y. hou modified ilon (longitude index) !
91 ! computing formula in subroutine getgases to work !
92 ! properly for models with either of 0->360 or !
93 ! -180->180 zonal grid directions. !
94 ! !
95 ! !
96 !!!!! ========================================================== !!!!!
97 !!!!! end descriptions !!!!!
98 !!!!! ========================================================== !!!!!
99 
100 
105 !========================================!
107 !........................................!
108 !
109  use physparam, only : ico2flg, ictmflg, ioznflg, ivflip, &
112  & kind_phys, kind_io4
113  use funcphys, only : fpkapx
114  use physcons, only : con_pi
115  use ozne_def, only : jmr => latsozc, loz => levozc, &
116  & blte => blatc, dlte=> dphiozc, &
117  & timeozc => timeozc
118  use module_iounitdef, only : nio3clm, nico2cn
119 !
120  implicit none
121 !
122  private
123 
124 ! --- version tag and last revision date
125  character(40), parameter :: &
126  & VTAGGAS='NCEP-Radiation_gases v5.1 Nov 2012 '
127 ! & VTAGGAS='NCEP-Radiation_gases v5.0 Aug 2012 '
128 
129 ! --- parameter constants
131  integer, parameter, public :: nf_vgas = 10 ! number of gas species
133  integer, parameter :: imxco2 = 24 ! input co2 data lon points
135  integer, parameter :: jmxco2 = 12 ! input co2 data lat points
137  integer, parameter :: minyear = 1957 ! earlist year 2-d co2 data
138  ! available
139 
140  real (kind=kind_phys), parameter :: resco2=15.0 ! horiz res in degree
141  real (kind=kind_phys), parameter :: raddeg=180.0/con_pi ! rad->deg conversion
142  real (kind=kind_phys), parameter :: prsco2=788.0 ! pres lim for 2-d co2 (mb)
143  real (kind=kind_phys), parameter :: hfpi =0.5*con_pi ! half of pi
144 
145 ! --- parameter constants for gas volume mixing ratioes
146  real (kind=kind_phys), parameter :: co2vmr_def = 350.0e-6
147  real (kind=kind_phys), parameter :: n2ovmr_def = 0.31e-6
148  real (kind=kind_phys), parameter :: ch4vmr_def = 1.50e-6
149  real (kind=kind_phys), parameter :: o2vmr_def = 0.209
150  real (kind=kind_phys), parameter :: covmr_def = 1.50e-8
151  real (kind=kind_phys), parameter :: f11vmr_def = 3.520e-10 ! aer 2003 value
152  real (kind=kind_phys), parameter :: f12vmr_def = 6.358e-10 ! aer 2003 value
153  real (kind=kind_phys), parameter :: f22vmr_def = 1.500e-10 ! aer 2003 value
154  real (kind=kind_phys), parameter :: cl4vmr_def = 1.397e-10 ! aer 2003 value
155  real (kind=kind_phys), parameter :: f113vmr_def= 8.2000e-11 ! gfdl 1999 value
156 
157 ! --- ozone seasonal climatology parameters defined in module ozne_def
158 ! - 4x5 ozone data parameter
159 ! integer, parameter :: JMR=45, LOZ=17
160 ! real (kind=kind_phys), parameter :: blte=-86.0, dlte=4.0
161 ! - geos ozone data
162 ! integer, parameter :: JMR=18, LOZ=17
163 ! real (kind=kind_phys), parameter :: blte=-85.0, dlte=10.0
164 
165 ! --- module variables to be set in subroutin gas_init and/or gas_update
166 
167 ! - variables for climatology ozone (ioznflg = 0)
168  real (kind=kind_phys), allocatable :: pkstr(:), o3r(:,:,:)
169  integer :: k1oz = 0, k2oz = 0
170  real (kind=kind_phys) :: facoz = 0.0
171 
172 ! - arrays for co2 2-d monthly data and global mean values from observed data
173  real (kind=kind_phys), allocatable :: co2vmr_sav(:,:,:)
174  real (kind=kind_phys), allocatable :: co2cyc_sav(:,:,:)
175 
176  real (kind=kind_phys) :: co2_glb = co2vmr_def
177  real (kind=kind_phys) :: gco2cyc(12)
178  data gco2cyc(:) / 12*0.0 /
179 
180  integer :: kyrsav = 0
181  integer :: kmonsav = 1
182 
183 ! --- public interfaces
184 
186 
187 
188 ! =================
189  contains
190 ! =================
191 
214 !-----------------------------------
215  subroutine gas_init &
216 !...................................
218 ! --- inputs:
219  & ( me )
220 ! --- outputs: ( none )
221 
222 ! =================================================================== !
223 ! !
224 ! gas_init sets up ozone, co2, etc. parameters. if climatology ozone !
225 ! then read in monthly ozone data. !
226 ! !
227 ! inputs: dimemsion !
228 ! me - print message control flag 1 !
229 ! !
230 ! outputs: (to the module variables) !
231 ! ( none ) !
232 ! !
233 ! external module variables: (in physparam) !
234 ! ico2flg - co2 data source control flag !
235 ! =0: use prescribed co2 global mean value !
236 ! =1: use input global mean co2 value (co2_glb) !
237 ! =2: use input 2-d monthly co2 value (co2vmr_sav) !
238 ! ictmflg - =yyyy#, data ic time/date control flag !
239 ! = -2: same as 0, but superimpose seasonal cycle !
240 ! from climatology data set. !
241 ! = -1: use user provided external data for the fcst !
242 ! time, no extrapolation. !
243 ! = 0: use data at initial cond time, if not existed!
244 ! then use latest, without extrapolation. !
245 ! = 1: use data at the forecast time, if not existed!
246 ! then use latest and extrapolate to fcst time.!
247 ! =yyyy0: use yyyy data for the forecast time, no !
248 ! further data extrapolation. !
249 ! =yyyy1: use yyyy data for the fcst. if needed, do !
250 ! extrapolation to match the fcst time. !
251 ! ioznflg - ozone data control flag !
252 ! =0: use climatological ozone profile !
253 ! >0: use interactive ozone profile !
254 ! ivflip - vertical profile indexing flag !
255 ! co2usr_file- external co2 user defined data table !
256 ! co2cyc_file- external co2 climotology monthly cycle data table !
257 ! !
258 ! internal module variables: !
259 ! pkstr, o3r - arrays for climatology ozone data !
260 ! !
261 ! usage: call gas_init !
262 ! !
263 ! subprograms called: none !
264 ! !
265 ! =================================================================== !
266 !
267  implicit none
268 
269 ! --- inputs:
270  integer, intent(in) :: me
271 
272 ! --- output: ( none )
273 
274 ! --- locals:
275  real (kind=kind_phys), dimension(IMXCO2,JMXCO2) :: co2dat
276  real (kind=kind_phys) :: co2g1, co2g2
277  real (kind=kind_phys) :: pstr(loz)
278  real (kind=kind_io4) :: o3clim4(jmr,loz,12), pstr4(loz)
279 
280  integer :: imond(12), ilat(jmr,12)
281  integer :: i, j, k, iyr, imo
282  logical :: file_exist, lextpl
283  character :: cline*100, cform*8
284  data cform / '(24f7.2)' / !! data format in IMXCO2*f7.2
285 !
286 !===> ... begin here
287 !
288  if ( me == 0 ) print *, vtaggas ! print out version tag
289 
290  kyrsav = 0
291  kmonsav = 1
292 
293 ! --- ... climatology ozone data section
294 
295  if ( ioznflg > 0 ) then
296  if ( me == 0 ) then
297  print *,' - Using interactive ozone distribution'
298  endif
299  else
300  if ( timeozc /= 12 ) then
301  print *,' - Using climatology ozone distribution'
302  print *,' timeozc=',timeozc, ' is not monthly mean', &
303  & ' - job aborting in subroutin gas_init!!!'
304  stop
305  endif
306 
307  allocate (pkstr(loz), o3r(jmr,loz,12))
308  rewind nio3clm
309 
310  if ( loz == 17 ) then ! For the operational ozone climatology
311  do k = 1, loz
312  read (nio3clm,15) pstr4(k)
313  15 format(f10.3)
314  enddo
315 
316  do imo = 1, 12
317  do j = 1, jmr
318  read (nio3clm,16) imond(imo), ilat(j,imo), &
319  & (o3clim4(j,k,imo),k=1,10)
320  16 format(i2,i4,10f6.2)
321  read (nio3clm,20) (o3clim4(j,k,imo),k=11,loz)
322  20 format(6x,10f6.2)
323  enddo
324  enddo
325  else ! For newer ozone climatology
326  read (nio3clm)
327  do k = 1, loz
328  read (nio3clm) pstr4(k)
329  enddo
330 
331  do imo = 1, 12
332  do k = 1, loz
333  read (nio3clm) (o3clim4(j,k,imo),j=1,jmr)
334  enddo
335  enddo
336  endif ! end if_LOZ_block
337 !
338  do imo = 1, 12
339  do k = 1, loz
340  do j = 1, jmr
341  o3r(j,k,imo) = o3clim4(j,k,imo) * 1.655e-6
342  enddo
343  enddo
344  enddo
345 
346  do k = 1, loz
347  pstr(k) = pstr4(k)
348  enddo
349 
350  if ( me == 0 ) then
351  print *,' - Using climatology ozone distribution'
352  print *,' Found ozone data for levels pstr=', &
353  & (pstr(k),k=1,loz)
354 ! print *,' O3=',(o3r(15,k,1),k=1,LOZ)
355  endif
356 
357  do k = 1, loz
358  pkstr(k) = fpkapx(pstr(k)*100.0)
359  enddo
360  endif ! end if_ioznflg_block
361 
362 ! --- ... co2 data section
363 
365 
366  lab_ico2 : if ( ico2flg == 0 ) then
367 
368  if ( me == 0 ) then
369  print *,' - Using prescribed co2 global mean value=', &
370  & co2vmr_def
371  endif
372 
373  else lab_ico2
374 
375  lab_ictm : if ( ictmflg == -1 ) then ! input user provided data
376 
377  inquire (file=co2usr_file, exist=file_exist)
378  if ( .not. file_exist ) then
379  print *,' Can not find user CO2 data file: ',co2usr_file, &
380  & ' - Stopped in subroutine gas_init !!'
381  stop
382  else
383  close (nico2cn)
384  open(nico2cn,file=co2usr_file,form='formatted',status='old')
385  rewind nico2cn
386  read (nico2cn, 25) iyr, cline, co2g1, co2g2
387  25 format(i4,a94,f7.2,16x,f5.2)
388  co2_glb = co2g1 * 1.0e-6
389 
390  if ( ico2flg == 1 ) then
391  if ( me == 0 ) then
392  print *,' - Using co2 global annual mean value from', &
393  & ' user provided data set:',co2usr_file
394  print *, iyr,cline(1:94),co2g1,' GROWTH RATE =', co2g2
395  endif
396  elseif ( ico2flg == 2 ) then
397  allocate ( co2vmr_sav(imxco2,jmxco2,12) )
398 
399  do imo = 1, 12
400  read (nico2cn,cform) co2dat
401 !check print cform, co2dat
402 
403  do j = 1, jmxco2
404  do i = 1, imxco2
405  co2vmr_sav(i,j,imo) = co2dat(i,j) * 1.0e-6
406  enddo
407  enddo
408  enddo
409 
410  if ( me == 0 ) then
411  print *,' - Using co2 monthly 2-d data from user', &
412  & ' provided data set:',co2usr_file
413  print *, iyr,cline(1:94),co2g1,' GROWTH RATE =', co2g2
414 
415  print *,' CHECK: Sample of selected months of CO2 data'
416  do imo = 1, 12, 3
417  print *,' Month =',imo
418  print *, co2vmr_sav(1,:,imo)
419  enddo
420  endif
421  else
422  print *,' ICO2=',ico2flg,' is not a valid selection', &
423  & ' - Stoped in subroutine gas_init!!!'
424  stop
425  endif ! endif_ico2flg_block
426 
427  close (nico2cn)
428  endif ! endif_file_exist_block
429 
430  else lab_ictm ! input from observed data
431 
432  if ( ico2flg == 1 ) then
433  if ( me == 0 ) then
434  print *,' - Using observed co2 global annual mean value'
435  endiF
436  elseif ( ico2flg == 2 ) then
437  allocate ( co2vmr_sav(imxco2,jmxco2,12) )
438 
439  if ( me == 0 ) then
440  print *,' - Using observed co2 monthly 2-d data'
441  endif
442  else
443  print *,' ICO2=',ico2flg,' is not a valid selection', &
444  & ' - Stoped in subroutine gas_init!!!'
445  stop
446  endif
447 
448  if ( ictmflg == -2 ) then
449  inquire (file=co2cyc_file, exist=file_exist)
450  if ( .not. file_exist ) then
451  if ( me == 0 ) then
452  print *,' Can not find seasonal cycle CO2 data: ', &
453  & co2cyc_file,' - Stopped in subroutine gas_init !!'
454  endif
455  stop
456  else
457  allocate( co2cyc_sav(imxco2,jmxco2,12) )
458 
459 ! --- ... read in co2 2-d seasonal cycle data
460  close (nico2cn)
461  open (nico2cn,file=co2cyc_file,form='formatted', &
462  & status='old')
463  rewind nico2cn
464  read (nico2cn, 35) cline, co2g1, co2g2
465  35 format(a98,f7.2,16x,f5.2)
466  read (nico2cn,cform) co2dat ! skip annual mean part
467 
468  if ( me == 0 ) then
469  print *,' - Superimpose seasonal cycle to mean CO2 data'
470  print *,' Opened CO2 climatology seasonal cycle data',&
471  & ' file: ',co2cyc_file
472 !check print *, cline(1:98), co2g1, co2g2
473  endif
474 
475  do imo = 1, 12
476  read (nico2cn,45) cline, gco2cyc(imo)
477  45 format(a58,f7.2)
478 !check print *, cline(1:58),gco2cyc(imo)
479  gco2cyc(imo) = gco2cyc(imo) * 1.0e-6
480 
481  read (nico2cn,cform) co2dat
482 !check print cform, co2dat
483  do j = 1, jmxco2
484  do i = 1, imxco2
485  co2cyc_sav(i,j,imo) = co2dat(i,j) * 1.0e-6
486  enddo
487  enddo
488  enddo
489 
490  close (nico2cn)
491  endif ! endif_file_exist_block
492  endif
493 
494  endif lab_ictm
495  endif lab_ico2
496 
497  return
498 !
499 !...................................
500  end subroutine gas_init
501 !-----------------------------------
502 
505 
506 !!\param[in] iyear integer, 1, year of the requested data for fcst
507 !!\param[in] imon integer, 1, month of the year
508 !!\param[in] iday integer, 1, day of the month
509 !!\param[in] ihour integer, 1, hour of the day
510 !!\param[in] loz1st logical, 1, clim ozone 1st time update control flag
511 !!\param[in] ldoco2 logical, 1, co2 update control flag
512 !!\param[in] me integer, 1, print message control flag
513 !!\param[out] NONE
514 !!\section external External Module Variables:
515 !!\n physparam::ico2flg - co2 data source control flag
516 !!\n =0: use prescribed co2 global mean value
517 !!\n =1: use input global mean co2 value (co2_glb)
518 !!\n =2: use input 2-d monthly co2 value (co2vmr_sav)
519 !!\n physparam::ictmflg - =yyyy#, data ic time/date control flag
520 !!\n = -2: same as 0, but superimpose seasonal cycle from climatology data set.
521 !!\n = -1: use user provided external data for the fcst time, no extrapolation.
522 !!\n = 0: use data at initial condition time, if not existed then use latest, without extrapolation.
523 !!\n = 1: use data at the forecast time, if not existed then use latest and extrapolate to fcst time.
524 !!\n =yyyy0: use yyyy data for the forecast time, no further data extrapolation.
525 !!\n =yyyy1: use yyyy data for the fcst. if needed, do extrapolation to match the fcst time.
526 !!\n physparam::ioznflg - ozone data control flag
527 !!\n =0: use climatological ozone profile
528 !!\n >0: use interactive ozone profile
529 !!\n physparam::ivflip - vertical profile indexing flag
530 !!\n physparam::co2usr_file - external co2 user defined data table
531 !!\n physparam::co2cyc_file - external co2 climotology monthly cycle data table
532 !-----------------------------------
533  subroutine gas_update &
534 !...................................
536 ! --- inputs:
537  & ( iyear, imon, iday, ihour, loz1st, ldoco2, me )
538 ! --- outputs: ( none )
539 
540 ! =================================================================== !
541 ! !
542 ! gas_update reads in 2-d monthly co2 data set for a specified year. !
543 ! data are in a 15 degree lat/lon horizontal resolution. !
544 ! !
545 ! inputs: dimemsion !
546 ! iyear - year of the requested data for fcst 1 !
547 ! imon - month of the year 1 !
548 ! iday - day of the month 1 !
549 ! ihour - hour of the day 1 !
550 ! loz1st - clim ozone 1st time update control flag 1 !
551 ! ldoco2 - co2 update control flag 1 !
552 ! me - print message control flag 1 !
553 ! !
554 ! outputs: (to the module variables) !
555 ! ( none ) !
556 ! !
557 ! external module variables: (in physparam) !
558 ! ico2flg - co2 data source control flag !
559 ! =0: use prescribed co2 global mean value !
560 ! =1: use input global mean co2 value (co2_glb) !
561 ! =2: use input 2-d monthly co2 value (co2vmr_sav) !
562 ! ictmflg - =yyyy#, data ic time/date control flag !
563 ! = -2: same as 0, but superimpose seasonal cycle !
564 ! from climatology data set. !
565 ! = -1: use user provided external data for the fcst !
566 ! time, no extrapolation. !
567 ! = 0: use data at initial cond time, if not existed!
568 ! then use latest, without extrapolation. !
569 ! = 1: use data at the forecast time, if not existed!
570 ! then use latest and extrapolate to fcst time.!
571 ! =yyyy0: use yyyy data for the forecast time, no !
572 ! further data extrapolation. !
573 ! =yyyy1: use yyyy data for the fcst. if needed, do !
574 ! extrapolation to match the fcst time. !
575 ! ioznflg - ozone data control flag !
576 ! =0: use climatological ozone profile !
577 ! >0: use interactive ozone profile !
578 ! ivflip - vertical profile indexing flag !
579 ! co2dat_file- external co2 2d monthly obsv data table !
580 ! co2gbl_file- external co2 global annual mean data table !
581 ! !
582 ! internal module variables: !
583 ! co2vmr_sav - monthly co2 volume mixing ratio IMXCO2*JMXCO2*12 !
584 ! co2cyc_sav - monthly cycle co2 vol mixing ratio IMXCO2*JMXCO2*12 !
585 ! co2_glb - global annual mean co2 mixing ratio !
586 ! gco2cyc - global monthly mean co2 variation 12 !
587 ! k1oz,k2oz,facoz !
588 ! - climatology ozone parameters 1 !
589 ! !
590 ! usage: call gas_update !
591 ! !
592 ! subprograms called: none !
593 ! !
594 ! =================================================================== !
595 !
596  implicit none
597 
598 ! --- inputs:
599  integer, intent(in) :: iyear, imon, iday, ihour, me
600 
601  logical, intent(in) :: loz1st, ldoco2
602 
603 ! --- output: ( none )
604 
605 ! --- locals:
606  real (kind=kind_phys), dimension(IMXCO2,JMXCO2) :: co2dat, co2ann
607  real (kind=kind_phys) :: co2g1, co2g2, rate
608 
609  integer :: i, id, j, l, iyr, imo, iyr1, iyr2, jyr, idyr
610  integer, save :: mdays(13), midmon=15, midm=15, midp=45
611 ! --- number of days in a month
612  data mdays / 31,28,31,30,31,30,31,31,30,31,30,31,30 /
613 
614  logical :: file_exist, lextpl, change
615  character :: cline*100, cform*8, cfile1*26
616  data cform / '(24f7.2)' / !! data format in IMXCO2*f7.2
617 !
618 !===> ... begin here
619 !
620 ! --- ... ozone data section
621 
622  if ( ioznflg == 0 ) then
623  midmon = mdays(imon)/2 + 1
624  change = loz1st .or. ( (iday==midmon) .and. (ihour==0) )
625 !
626  if ( change ) then
627  if ( iday < midmon ) then
628  k1oz = mod(imon+10, 12) + 1
629  midm = mdays(k1oz)/2 + 1
630  k2oz = imon
631  midp = mdays(k1oz) + midmon
632  else
633  k1oz = imon
634  midm = midmon
635  k2oz = mod(imon, 12) + 1
636  midp = mdays(k2oz)/2 + 1 + mdays(k1oz)
637  endif
638  endif
639 !
640  if (iday < midmon) then
641  id = iday + mdays(k1oz)
642  else
643  id = iday
644  endif
645 
646  facoz = float(id - midm) / float(midp - midm)
647  endif
648 
649 ! --- ... co2 data section
650 
651  if ( ico2flg == 0 ) return ! use prescribed global mean co2 data
652  if ( ictmflg ==-1 ) return ! use user provided co2 data
653  if ( .not. ldoco2 ) return ! no need to update co2 data
654 
655  if ( ictmflg < 0 ) then ! use user provided external data
656  lextpl = .false. ! no time extrapolation
657  idyr = iyear ! use the model year
658  else ! use historically observed data
659  lextpl = ( mod(ictmflg,10) == 1 ) ! flag for data extrapolation
660  idyr = ictmflg / 10 ! year of data source used
661  if ( idyr == 0 ) idyr = iyear ! not specified, use model year
662  endif
663 
664 ! --- ... auto select co2 2-d data table for required year
665 
666  kmonsav = imon
667  if ( kyrsav == iyear ) return
668  kyrsav = iyear
669  iyr = iyear
670 
671 ! --- ... for data earlier than MINYEAR (1957), the data are in
672 ! the form of semi-yearly global mean values. otherwise,
673 ! data are monthly mean in horizontal 2-d map.
674 
675  lab_if_idyr : if ( idyr < minyear .and. ictmflg > 0 ) then
676 
677  if ( me == 0 ) then
678  print *,' Requested CO2 data year',iyear,' earlier than', &
679  & minyear
680  print *,' Which is the earliest monthly observation', &
681  & ' data available.'
682  print *,' Thus, historical global mean data is used'
683  endif
684 
685 ! --- ... check to see if requested co2 data file existed
686 
687  inquire (file=co2gbl_file, exist=file_exist)
688  if ( .not. file_exist ) then
689  print *,' Requested co2 data file "',co2gbl_file, &
690  & '" not found - Stopped in subroutine gas_update!!'
691  stop
692  else
693  close(nico2cn)
694  open (nico2cn,file=co2gbl_file,form='formatted',status='old')
695  rewind nico2cn
696 
697  read (nico2cn, 24) iyr1, iyr2, cline
698  24 format(i4,4x,i4,a48)
699 
700  if ( me == 0 ) then
701  print *,' Opened co2 data file: ',co2gbl_file
702 !check print *, iyr1, iyr2, cline(1:48)
703  endif
704 
705  if ( idyr < iyr1 ) then
706  iyr = iyr1
707 !check if ( me == 0 ) then
708 ! print *,' Using earlist available co2 data, year=',iyr1
709 !check endif
710  endif
711 
712  i = iyr2
713  lab_dowhile1 : do while ( i >= iyr1 )
714 ! read (NICO2CN,26) jyr, co2g1, co2g2
715 ! 26 format(i4,4x,2f7.2)
716  read (nico2cn, *) jyr, co2g1, co2g2
717 
718  if ( i == iyr .and. iyr == jyr ) then
719  co2_glb = (co2g1+co2g2) * 0.5e-6
720  if ( ico2flg == 2 ) then
721  do j = 1, jmxco2
722  do i = 1, imxco2
723  co2vmr_sav(i,j,1:6) = co2g1 * 1.0e-6
724  co2vmr_sav(i,j,7:12) = co2g2 * 1.0e-6
725  enddo
726  enddo
727  endif
728 
729  if ( me == 0 ) print *,' Co2 data for year',iyear, &
730  & co2_glb
731  exit lab_dowhile1
732  else
733 !check if ( me == 0 ) print *,' Skip co2 data for year',i
734  i = i - 1
735  endif
736  enddo lab_dowhile1
737 
738  close ( nico2cn )
739  endif ! end if_file_exist_block
740 
741  else lab_if_idyr
742 
743 ! --- ... set up input data file name
744 
745  cfile1 = co2dat_file
746  write(cfile1(19:22),34) idyr
747  34 format(i4.4)
748 
749 ! --- ... check to see if requested co2 data file existed
750 
751  inquire (file=cfile1, exist=file_exist)
752  if ( .not. file_exist ) then
753 
754  lab_if_ictm : if ( ictmflg > 10 ) then ! specified year of data not found
755  if ( me == 0 ) then
756  print *,' Specified co2 data for year',idyr, &
757  & ' not found !! Need to change namelist ICTM !!'
758  print *,' *** Stopped in subroutine gas_update !!'
759  endif
760  stop
761  else lab_if_ictm ! looking for latest available data
762  if ( me == 0 ) then
763  print *,' Requested co2 data for year',idyr, &
764  & ' not found, check for other available data set'
765  endif
766 
767  lab_dowhile2 : do while ( iyr >= minyear )
768  iyr = iyr - 1
769  write(cfile1(19:22),34) iyr
770 
771  inquire (file=cfile1, exist=file_exist)
772  if ( me == 0 ) then
773  print *,' Looking for CO2 file ',cfile1
774  endif
775 
776  if ( file_exist ) then
777  exit lab_dowhile2
778  endif
779  enddo lab_dowhile2
780 
781  if ( .not. file_exist ) then
782  if ( me == 0 ) then
783  print *,' Can not find co2 data source file'
784  print *,' *** Stopped in subroutine gas_update !!'
785  endif
786  stop
787  endif
788  endif lab_if_ictm
789  endif ! end if_file_exist_block
790 
791 ! --- ... read in co2 2-d data for the requested month
792 
793  close(nico2cn)
794  open (nico2cn,file=cfile1,form='formatted',status='old')
795  rewind nico2cn
796  read (nico2cn, 36) iyr, cline, co2g1, co2g2
797  36 format(i4,a94,f7.2,16x,f5.2)
798 
799  if ( me == 0 ) then
800  print *,' Opened co2 data file: ',cfile1
801  print *, iyr, cline(1:94), co2g1,' GROWTH RATE =', co2g2
802  endif
803 
804 ! --- ... add growth rate if needed
805  if ( lextpl ) then
806 ! rate = co2g2 * (iyear - iyr) ! rate from early year
807 ! rate = 1.60 * (iyear - iyr) ! avg rate over long period
808  rate = 2.00 * (iyear - iyr) ! avg rate for recent period
809  else
810  rate = 0.0
811  endif
812 
813  co2_glb = (co2g1 + rate) * 1.0e-6
814  if ( me == 0 ) then
815  print *,' Global annual mean CO2 data for year', &
816  & iyear, co2_glb
817  endif
818 
819  if ( ictmflg == -2 ) then ! need to calc ic time annual mean first
820 
821  if ( ico2flg == 1 ) then
822  if ( me==0 ) then
823  print *,' CHECK: Monthly deviations of climatology ', &
824  & 'to be superimposed on global annual mean'
825  print *, gco2cyc
826  endif
827  elseif ( ico2flg == 2 ) then
828  co2ann(:,:) = 0.0
829 
830  do imo = 1, 12
831  read (nico2cn,cform) co2dat
832 !check print cform, co2dat
833 
834  do j = 1, jmxco2
835  do i = 1, imxco2
836  co2ann(i,j) = co2ann(i,j) + co2dat(i,j)
837  enddo
838  enddo
839  enddo
840 
841  do j = 1, jmxco2
842  do i = 1, imxco2
843  co2ann(i,j) = co2ann(i,j) * 1.0e-6 / float(12)
844  enddo
845  enddo
846 
847  do imo = 1, 12
848  do j = 1, jmxco2
849  do i = 1, imxco2
850  co2vmr_sav(i,j,imo) = co2ann(i,j)+co2cyc_sav(i,j,imo)
851  enddo
852  enddo
853  enddo
854 
855  if ( me==0 ) then
856  print *,' CHECK: Sample of 2-d annual mean of CO2 ', &
857  & 'data used for year:',iyear
858  print *, co2ann(1,:)
859  print *,' CHECK: AFTER adding seasonal cycle, Sample ', &
860  & 'of selected months of CO2 data for year:',iyear
861  do imo = 1, 12, 3
862  print *,' Month =',imo
863  print *, co2vmr_sav(1,:,imo)
864  enddo
865  endif
866  endif ! endif_icl2flg_block
867 
868  else ! no need to calc ic time annual mean first
869 
870  if ( ico2flg == 2 ) then ! directly save monthly data
871  do imo = 1, 12
872  read (nico2cn,cform) co2dat
873 !check print cform, co2dat
874 
875  do j = 1, jmxco2
876  do i = 1, imxco2
877  co2vmr_sav(i,j,imo) = (co2dat(i,j) + rate) * 1.0e-6
878  enddo
879  enddo
880  enddo
881 
882  if ( me == 0 ) then
883  print *,' CHECK: Sample of selected months of CO2 ', &
884  & 'data used for year:',iyear
885  do imo = 1, 12, 3
886  print *,' Month =',imo
887  print *, co2vmr_sav(1,:,imo)
888  enddo
889  endif
890  endif ! endif_ico2flg_block
891 
892  do imo = 1, 12
893  gco2cyc(imo) = 0.0
894  enddo
895  endif ! endif_ictmflg_block
896  close ( nico2cn )
897 
898  endif lab_if_idyr
899 
900  return
901 !
902 !...................................
903  end subroutine gas_update
904 !-----------------------------------
905 
930 !-----------------------------------
931  subroutine getgases &
932 !...................................
934 ! --- inputs:
935  & ( plvl, xlon, xlat, &
936  & imax, lmax, &
937 ! --- outputs:
938  & gasdat &
939  & )
940 
941 ! =================================================================== !
942 ! !
943 ! getgases set up global distribution of radiation absorbing gases !
944 ! in volume mixing ratio. currently only co2 has the options from !
945 ! observed values, all other gases are asigned to the climatological !
946 ! values. !
947 ! !
948 ! inputs: !
949 ! plvl(IMAX,LMAX+1)- pressure at model layer interfaces (mb) !
950 ! xlon(IMAX) - grid longitude in radians, ok both 0->2pi or !
951 ! -pi -> +pi arrangements !
952 ! xlat(IMAX) - grid latitude in radians, default range to !
953 ! pi/2 -> -pi/2, otherwise see in-line comment !
954 ! IMAX, LMAX - horiz, vert dimensions for output data !
955 ! !
956 ! outputs: !
957 ! gasdat(IMAX,LMAX,NF_VGAS) - gases volume mixing ratioes !
958 ! (:,:,1) - co2 !
959 ! (:,:,2) - n2o !
960 ! (:,:,3) - ch4 !
961 ! (:,:,4) - o2 !
962 ! (:,:,5) - co !
963 ! (:,:,6) - cfc11 !
964 ! (:,:,7) - cfc12 !
965 ! (:,:,8) - cfc22 !
966 ! (:,:,9) - ccl4 !
967 ! (:,:,10) - cfc113 !
968 ! !
969 ! external module variables: (in physparam) !
970 ! ico2flg - co2 data source control flag !
971 ! =0: use prescribed co2 global mean value !
972 ! =1: use input global mean co2 value (co2_glb) !
973 ! =2: use input 2-d monthly co2 value (co2vmr_sav) !
974 ! ivflip - vertical profile indexing flag !
975 ! !
976 ! internal module variables used: !
977 ! co2vmr_sav - saved monthly co2 concentration from sub gas_update !
978 ! co2_glb - saved global annual mean co2 value from gas_update !
979 ! gco2cyc - saved global seasonal variation of co2 climatology !
980 ! in 12-month form !
981 ! ** note: for lower atmos co2vmr_sav may have clim monthly deviations !
982 ! superimposed on init-cond co2 value, while co2_glb only !
983 ! contains the global mean value, thus needs to add the !
984 ! monthly dglobal mean deviation gco2cyc at upper atmos. for !
985 ! ictmflg/=-2, this value will be zero. !
986 ! !
987 ! usage: call getgases !
988 ! !
989 ! subprograms called: none !
990 ! !
991 ! =================================================================== !
992 !
993  implicit none
994 
995 ! --- input:
996  integer, intent(in) :: IMAX, LMAX
997  real (kind=kind_phys), intent(in) :: plvl(:,:), xlon(:), xlat(:)
998 
999 ! --- output:
1000  real (kind=kind_phys), intent(out) :: gasdat(:,:,:)
1001 
1002 ! --- local:
1003  integer :: i, k, ilat, ilon
1004 
1005  real (kind=kind_phys) :: xlon1, xlat1, tmp
1006 
1007 !===> ... begin here
1008 
1009 ! --- ... assign default values
1010 
1011  do k = 1, lmax
1012  do i = 1, imax
1013  gasdat(i,k,1) = co2vmr_def
1014  gasdat(i,k,2) = n2ovmr_def
1015  gasdat(i,k,3) = ch4vmr_def
1016  gasdat(i,k,4) = o2vmr_def
1017  gasdat(i,k,5) = covmr_def
1018  gasdat(i,k,6) = f11vmr_def
1019  gasdat(i,k,7) = f12vmr_def
1020  gasdat(i,k,8) = f22vmr_def
1021  gasdat(i,k,9) = cl4vmr_def
1022  gasdat(i,k,10)= f113vmr_def
1023  enddo
1024  enddo
1025 
1026 ! --- ... co2 section
1027 
1028  if ( ico2flg == 1 ) then
1029 ! --- use obs co2 global annual mean value only
1030 
1031  do k = 1, lmax
1032  do i = 1, imax
1033  gasdat(i,k,1) = co2_glb + gco2cyc(kmonsav)
1034  enddo
1035  enddo
1036 
1037  elseif ( ico2flg == 2 ) then
1038 ! --- use obs co2 monthly data with 2-d variation at lower atmos
1039 ! otherwise use global mean value
1040 
1041  tmp = raddeg / resco2
1042  do i = 1, imax
1043  xlon1 = xlon(i)
1044  if ( xlon1 < 0.0 ) xlon1 = xlon1 + con_pi ! if xlon in -pi->pi, convert to 0->2pi
1045  xlat1 = hfpi - xlat(i) ! if xlat in pi/2 -> -pi/2 range
1046 !note xlat1 = xlat(i) ! if xlat in 0 -> pi range
1047 
1048  ilon = min( imxco2, int( xlon1*tmp + 1 ))
1049  ilat = min( jmxco2, int( xlat1*tmp + 1 ))
1050 
1051  if ( ivflip == 0 ) then ! index from toa to sfc
1052  do k = 1, lmax
1053  if ( plvl(i,k) >= prsco2 ) then
1054  gasdat(i,k,1) = co2vmr_sav(ilon,ilat,kmonsav)
1055  else
1056  gasdat(i,k,1) = co2_glb + gco2cyc(kmonsav)
1057  endif
1058  enddo
1059  else ! index from sfc to toa
1060  do k = 1, lmax
1061  if ( plvl(i,k+1) >= prsco2 ) then
1062  gasdat(i,k,1) = co2vmr_sav(ilon,ilat,kmonsav)
1063  else
1064  gasdat(i,k,1) = co2_glb + gco2cyc(kmonsav)
1065  endif
1066  enddo
1067  endif
1068  enddo
1069  endif
1070 
1071 !
1072  return
1073 !...................................
1074  end subroutine getgases
1075 !-----------------------------------
1076 
1083 !-----------------------------------
1084  subroutine getozn &
1085 !...................................
1087 ! --- inputs:
1088  & ( prslk,xlat, &
1089  & imax, lm, &
1090 ! --- outputs:
1091  & o3mmr &
1092  & )
1093 
1094 ! =================================================================== !
1095 ! !
1096 ! getozn sets up climatological ozone profile for radiation calculation!
1097 ! !
1098 ! this code is originally written By Shrinivas Moorthi !
1099 ! !
1100 ! inputs: !
1101 ! prslk (IMAX,LM) - exner function = (p/p0)**rocp !
1102 ! xlat (IMAX) - latitude in radians, default to pi/2 -> -pi/2 !
1103 ! range, otherwise see in-line comment !
1104 ! IMAX, LM - horizontal and vertical dimensions !
1105 ! !
1106 ! outputs: !
1107 ! o3mmr (IMAX,LM) - output ozone profile in mass mixing ratio (g/g)!
1108 ! !
1109 ! module variables: !
1110 ! k1oz, k2oz - ozone data interpolation indices !
1111 ! facoz - ozone data interpolation factor !
1112 ! ivflip - control flag for direction of vertical index !
1113 ! !
1114 ! usage: call getozn !
1115 ! !
1116 ! =================================================================== !
1117 !
1118  implicit none
1119 
1120 ! --- inputs:
1121  integer, intent(in) :: IMAX, LM
1122 
1123  real (kind=kind_phys), intent(in) :: prslk(:,:), xlat(:)
1124 
1125 ! --- outputs:
1126  real (kind=kind_phys), intent(out) :: o3mmr(:,:)
1127 
1128 ! --- locals:
1129  real (kind=kind_phys) :: o3i(imax,loz), wk1(imax), deglat, elte, &
1130  & tem, tem1, tem2, tem3, tem4, temp
1131 
1132  integer :: i, j, k, l, j1, j2, ll
1133 !
1134 !===> ... begin here
1135 !
1136  elte = blte + (jmr-1)*dlte
1137 
1138  do i = 1, imax
1139  deglat = xlat(i) * raddeg ! if xlat in pi/2 -> -pi/2 range
1140 ! deglat = 90.0 - xlat(i)*raddeg ! if xlat in 0 -> pi range
1141 
1142  if (deglat > blte .and. deglat < elte) then
1143  tem1 = (deglat - blte) / dlte + 1
1144  j1 = tem1
1145  j2 = j1 + 1
1146  tem1 = tem1 - j1
1147  elseif (deglat <= blte) then
1148  j1 = 1
1149  j2 = 1
1150  tem1 = 1.0
1151  elseif (deglat >= elte) then
1152  j1 = jmr
1153  j2 = jmr
1154  tem1 = 1.0
1155  endif
1156 
1157  tem2 = 1.0 - tem1
1158  do j = 1, loz
1159  tem3 = tem2*o3r(j1,j,k1oz) + tem1*o3r(j2,j,k1oz)
1160  tem4 = tem2*o3r(j1,j,k2oz) + tem1*o3r(j2,j,k2oz)
1161  o3i(i,j) = tem4*facoz + tem3*(1.0 - facoz)
1162  enddo
1163  enddo
1164 
1165  do l = 1, lm
1166  ll = l
1167  if (ivflip == 1) ll = lm -l + 1
1168 
1169  do i = 1, imax
1170  wk1(i) = prslk(i,ll)
1171  enddo
1172 
1173  do k = 1, loz-1
1174  temp = 1.0 / (pkstr(k+1) - pkstr(k))
1175 
1176  do i = 1, imax
1177  if (wk1(i) > pkstr(k) .and. wk1(i) <= pkstr(k+1)) then
1178  tem = (pkstr(k+1) - wk1(i)) * temp
1179  o3mmr(i,ll) = tem * o3i(i,k) + (1.0 - tem) * o3i(i,k+1)
1180  endif
1181  enddo
1182  enddo
1183 
1184  do i = 1, imax
1185  if (wk1(i) > pkstr(loz)) o3mmr(i,ll) = o3i(i,loz)
1186  if (wk1(i) < pkstr(1)) o3mmr(i,ll) = o3i(i,1)
1187  enddo
1188  enddo
1189 !
1190  return
1191 !...................................
1192  end subroutine getozn
1193 !-----------------------------------
1194 
1195 !
1196 !........................................!
1197  end module module_radiation_gases !
1198 !========================================!
real(kind=kind_phys), dimension(:), allocatable pkstr
real(kind=kind_phys), parameter ch4vmr_def
real(kind=kind_phys), parameter n2ovmr_def
subroutine, public gas_init
This subroutine sets up ozone, co2, etc. parameters. if climatology ozone then read in monthly ozone ...
integer, save ioznflg
ozone dta source control flag
Definition: physparam.f:175
real(kind=kind_phys), parameter co2vmr_def
real(kind=kind_phys), parameter o2vmr_def
integer, parameter imxco2
input co2 dat lon points
character, save co2usr_file
external co2 user defined data table
Definition: physparam.f:181
integer, parameter, public nf_vgas
number of gas species
character, save co2gbl_file
external co2 global annual mean data table
Definition: physparam.f:179
real(kind=kind_phys), parameter hfpi
real(kind=kind_phys), parameter raddeg
subroutine, public gas_update
This subroutine reads in 2-d monthly co2 data set for a specified year. data are in a 15 degree lat/l...
This module contains some the most frequently used math and physics constants for gcm models...
Definition: physcons.f:29
real(kind=kind_phys), parameter f12vmr_def
character(40), parameter vtaggas
subroutine, public getgases
This subroutine sets up global distribution of radiation absorbing gases in volume mixing ratio...
This module defines commonly used control variables/parameters in physics related programs...
Definition: physparam.f:27
integer, save ico2flg
co2 data source control flag
Definition: physparam.f:171
real(kind=kind_phys), dimension(12) gco2cyc
real(kind=kind_phys), dimension(:,:,:), allocatable o3r
character, save co2dat_file
external co2 2d monthly obsv data table
Definition: physparam.f:177
real(kind=kind_phys), parameter f113vmr_def
integer, save ivflip
vertical profile indexing flag
Definition: physparam.f:224
real(kind=kind_phys) co2_glb
real(kind=kind_phys), parameter covmr_def
character, save co2cyc_file
external co2 clim monthly cycle data table
Definition: physparam.f:183
real(kind=kind_phys), dimension(:,:,:), allocatable co2vmr_sav
integer, parameter jmxco2
input co2 data lat points
real(kind=kind_phys), parameter resco2
real(kind=kind_phys), parameter prsco2
real(kind=kind_phys), dimension(:,:,:), allocatable co2cyc_sav
integer, parameter minyear
earlist year 2-d co2 data available
real(kind=kind_phys), parameter con_pi
Definition: physcons.f:41
real(kind=kind_phys), parameter f11vmr_def
subroutine, public getozn
This subroutine sets up climatological ozone profile for radiation calculation this code is originall...
This module sets up ozone climatological profiles and other constant gas profiles, such as co2, ch4, n2o, o2, and those of cfc gases. All data are entered as mixing ratio by volume, except ozone which is mass mixing ratio (g/g).
real(kind=kind_phys), parameter f22vmr_def
real(kind=kind_phys) facoz
integer, save ictmflg
external data time/date control flag
Definition: physparam.f:173
real(kind=kind_phys), parameter cl4vmr_def