program format93 *_______________________________________________________________________ * date: 2/17/93 * author : Bill Ridgway (ridgway@climate.gsfc.nasa.gov, 301-286-9138) * purpose: to covert gfdl, gsfc, or kiae data formats to fm93 * routines : * rgfdl : read data in gfdl format from unit 5 * rgsfc : read data in gsfc format from unit 5 * rkiae : read data in kiae format from unit 5 * rfm93 : read data in fm93 format from unit 5 * wfm93 : write data in fm93 format to unit 6 *_______________________________________________________________________ * proceedure: 1. ingest data at 5-10 cm-1 resolution into program memory * 2. write data in newly specified format to output file *_______________________________________________________________________ * parameters: * mlayer: maximum number of model layers * mlevel: maximum number of model levels * mbands: maximum number of spectral bands *_______________________________________________________________________ * internal variables containing all model data * header(10): 10-line file header containing 80 characters per line * atmlbl : 3-character atmospheric descriptor * 'trp' tropical * 'mls' mid-latitude summer * 'mlw' mid-latitude winter * 'sas' sub-arctic summer * 'saw' sub-arctic winter * '200' 200 (K) isothermal * '250' 250 (K) isothermal * '300' 300 (K) isothermal * iphase : icrccm phase 1 or 2 (see report) * icase : icrccm case number (see report) * nlayer : actual number of model layers * nlevel : actual number of model levels = nlayer+1 * nbands : actual number of spectral bands * psurf : surface pressure (mb) * tsurf : surface temperature (K) * ptrop : tropopause pressure (mb) * flgctm : flag to indicate h2o continuum is/is not (1/0) included * ppmco2 : co2 mixing ratio (ppm) * ppmch4 : ch4 mixing ratio (ppm) * ppmn2o : n2o mixing ratio (ppm) * pthair : total gas column for dry air (#air molecules/cm**2) * pthh2o : total h2o column (#h2o molecules/cm**2) * pthco2 : co2 path (#co2 molecules/cm**2) * ptho3 : ozone path (#ozone molecules/cm**2) * pthch4 : ch4 path (#ch4 molecules/cm**2) * pthn2o : n2o path (#n2o molecules/cm**2) * bandv1 : wavenumber lower limit of each band * bandv2 : wavenumber upper limit of each band * plevel : pressure levels starting at top (mb) * tlevel : temperature at pressure levels (K) * player : pressure at mid-point of each layer (mb) * tlayer : temperature characteristic of each layer (K) * alayer : dry air mass of layer (#air molecules/cm**2) * wlayer : water vapor volume mixing ratio of each layer (ppm) * olayer : ozone mixing ratio of each layer (ppm) * fluxup : radiative flux upward by level and band (w/m**2) * fluxdn : radiative flux downward by level and band (w/m**2) * fluxnt : net radiative flux upward by level and band (w/m**2) * heatrt : layer heating rates indexed by layer and band (K/day) * tropup : tropopause radiative flux upward by band (w/m**2) * tropdn : tropopause radiative flux downward by band (w/m**2) * tropnt : tropopause net radiative flux upward by band (w/m**2) *_______________________________________________________________________ * units and conversions: *_______________________ * dry air mass per layer = [1/G] * [1/Mair] * dP * (molecules/cm**2) [g/cm**2/mb] (#air/g) (mb) * H2O absorber amount = RH2O * [1/G] * [1/Mh2o] * dP * (molecules/cm**2) (g/g) [g/cm**2/mb] (#h2o/g) (mb) * CO2 amount = 10.e-6 * FCO2 * [1/G] * [1/Mair] * dP * (molecules/cm**2) (ppmv) [g/cm**2/mb] (#air/g) (mb) * O3 absorber amount = RO3 * [1/G] * [1/Mo3] * dP * (molecules/cm**2) (g/g) [g/cm**2/mb] (#o3/g) (mb) * water: RH2O is mass mixing ratio in (g/g) * co2: FCO2 is molecular mixing ratio in (ppmv) * ozone: RO3 is mass mixing ratio in (g/g) * [1/G] = 1.02 g/cm**2/mb * [1/Mh2o] = 6.0256e+23/18. #h2o/g * [1/Mo3] = 6.0256e+23/48. #o3/g * [1/Mair] = 6.0256e+23/28.97 #air/g (dry air) * dP = layer thickness mb * tropopause pressures for standard atmospheres: * trp=tropical 17 km or 93.7 mb * mls=mid-latitude summer 13 km or 179.0 mb * mlw=mid-latitude winter 10 km or 256.8 mb * sas=subarctic summer 10 km or 267.7 mb * saw=subarctic winter 9 km or 282.9 mb *_______________________________________________________________________ * parameters defined in driver and passed to all subroutines: integer iounit,mlayer,mlevel,mbands parameter (mlayer=150,mlevel=151,mbands=600) * variables passed to or created by subroutines: character*80 header(10) character*3 atmlbl integer iphase,icase,nlayer,nlevel,nbands real psurf,tsurf,ptrop,flgctm real ppmco2,ppmch4,ppmn2o real pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o real bandv1(0:mbands),bandv2(0:mbands) real plevel(mlevel),tlevel(mlevel) real player(mlayer),tlayer(mlayer) real alayer(mlayer),wlayer(mlayer),olayer(mlayer) real fluxup(mlevel,0:mbands) real fluxdn(mlevel,0:mbands) real fluxnt(mlevel,0:mbands) real heatrt(mlayer,0:mbands) real tropup( 0:mbands) real tropdn( 0:mbands) real tropnt( 0:mbands) *_______________________________________________________________________ * identify file format as gsfc, gfdl, etc from first 4 characters * ingest data by calling appropriate subroutine * re-write data in 93 format using subroutine wfm93 (write format-93) character*4 first4 iounit=5 read(iounit,'(a)') first4 rewind iounit if(first4.eq.'GFDL') then call rgfdl + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) else if(first4.eq.' ...') then call rgsfc + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) else if(first4.eq.'JUNE') then call rkiae + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) else if(first4.eq.'NASA') then call rfm93 + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) else write(6,'(a)') 'data file format not recognized.' stop endif iounit=6 call wfm93 + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) stop end subroutine rgfdl + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) *_______________________________________________________________________ * purpose: rgfdl - read gfdl data in earlier 1988 data format *_______________________________________________________________________ * parameters defined in driver and passed to all subroutines: integer iounit,mlayer,mlevel,mbands * variables passed to or created by subroutines: character*80 header(10) character*3 atmlbl integer iphase,icase,nlayer,nlevel,nbands real psurf,tsurf,ptrop,flgctm real ppmco2,ppmch4,ppmn2o real pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o real bandv1(0:mbands),bandv2(0:mbands) real plevel(mlevel),tlevel(mlevel) real player(mlayer),tlayer(mlayer) real alayer(mlayer),wlayer(mlayer),olayer(mlayer) real fluxup(mlevel,0:mbands) real fluxdn(mlevel,0:mbands) real fluxnt(mlevel,0:mbands) real heatrt(mlayer,0:mbands) real tropup( 0:mbands) real tropdn( 0:mbands) real tropnt( 0:mbands) *_______________________________________________________________________ * icrccm phase, case #, atmosphere, etc. entered from separate file open(1,file='temp_header') read(1,*) iphase read(1,*) icase read(1,'(a3)') atmlbl read(1,*) flgctm read(1,*) ppmco2 read(1,*) ppmch4 read(1,*) ppmn2o close(1) * blank out header initially do k=1,10 write(header(k),'(80a1)') (' ',j=1,80) end do * read 3 lines of header info from gfdl file, number of levels and bands read( iounit,'(a80)') header(1) read( iounit,'(a80)') header(3) read( iounit,'(a80)') header(5) read( iounit,*) nlevel,nbands nlayer=nlevel-1 * read spectral sums (k=0th elements) do k=0,0 read( iounit,*) bandv1(k),bandv2(k) read( iounit,*) (fluxnt(l,k), l = 1,nlevel ) read( iounit,*) (fluxup(l,k), l = 1,nlevel ) read( iounit,*) (fluxdn(l,k), l = 1,nlevel ) read( iounit,*) (heatrt(l,k), l = 1,nlayer ) read( iounit,*) ptrop,tropnt(k),tropup(k),tropdn(k) end do * read atmospheric profile information read( iounit,*) ( plevel(l), l = 1,nlevel ) read( iounit,*) ( tlevel(l), l = 1,nlevel ) read( iounit,*) ( player(l), l = 1,nlayer ) read( iounit,*) ( tlayer(l), l = 1,nlayer ) read( iounit,*) ( wlayer(l), l = 1,nlayer ) read( iounit,*) ( olayer(l), l = 1,nlayer ) * ingest flux and heating rate data band-by-band do k=1,nbands read( iounit,*) bandv1(k),bandv2(k) read( iounit,*) (fluxnt(l,k), l = 1,nlevel ) read( iounit,*) (fluxup(l,k), l = 1,nlevel ) read( iounit,*) (fluxdn(l,k), l = 1,nlevel ) read( iounit,*) (heatrt(l,k), l = 1,nlayer ) read( iounit,*) ptrop,tropnt(k),tropup(k),tropdn(k) end do * set up gfdl scaling factors for h2o and o3 if(iphase.eq.1.and.icase.eq.17) then wfactor=0.75 else if(iphase.eq.1.and.icase.eq.18) then wfactor=0.75 else if(iphase.eq.1.and.icase.eq.21) then wfactor=1.25 else if(iphase.eq.1.and.icase.eq.22) then wfactor=1.25 else if(iphase.eq.2.and.icase.eq.15) then wfactor=0.75 else if(iphase.eq.2.and.icase.eq.16) then wfactor=0.75 else if(iphase.eq.2.and.icase.eq.17) then wfactor=1.25 else if(iphase.eq.2.and.icase.eq.18) then wfactor=1.25 else wfactor=1. end if if(icase.eq.24.and.header(3)(30:40).eq.'0.75 X CONC') then ofactor=0.75 else if(icase.eq.24.and.header(3)(30:40).eq.'1.25 X CONC') then ofactor=1.25 else ofactor=1. end if * re-define units as ppm, molecules/cm**2 do j=1,nlayer alayer(j) = 1.02 * (6.0256e+23/28.97) * (plevel(j+1)-plevel(j)) wlayer(j) = wfactor * 1.e+6 * wlayer(j) * 28.97/18.00 olayer(j) = ofactor * 1.e+6 * olayer(j) * 28.97/48.00 end do * sum column amounts for dry air and 5 gases pthair = 0. pthh2o = 0. ptho3 = 0. do j=1,nlayer pthair = pthair + alayer(j) pthh2o = pthh2o + 1.e-6 * wlayer(j) * alayer(j) ptho3 = ptho3 + 1.e-6 * olayer(j) * alayer(j) end do pthco2 = 1.e-6 * ppmco2 * pthair pthch4 = 1.e-6 * ppmch4 * pthair pthn2o = 1.e-6 * ppmn2o * pthair * miscellaneous psurf = plevel(nlevel) tsurf = tlevel(nlevel) return end subroutine rkiae + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) *_______________________________________________________________________ * purpose: rkiae - read kiae icrccm data in earlier 1988 data format *_______________________________________________________________________ * parameters defined in driver and passed to all subroutines: integer iounit,mlayer,mlevel,mbands * variables passed to or created by subroutines: character*80 header(10) character*3 atmlbl integer iphase,icase,nlayer,nlevel,nbands real psurf,tsurf,ptrop,flgctm real ppmco2,ppmch4,ppmn2o real pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o real bandv1(0:mbands),bandv2(0:mbands) real plevel(mlevel),tlevel(mlevel) real player(mlayer),tlayer(mlayer) real alayer(mlayer),wlayer(mlayer),olayer(mlayer) real fluxup(mlevel,0:mbands) real fluxdn(mlevel,0:mbands) real fluxnt(mlevel,0:mbands) real heatrt(mlayer,0:mbands) real tropup( 0:mbands) real tropdn( 0:mbands) real tropnt( 0:mbands) *_______________________________________________________________________ * icrccm phase, case #, atmosphere, etc. entered from separate file open(1,file='temp_header') read(1,*) iphase read(1,*) icase read(1,'(a3)') atmlbl read(1,*) flgctm read(1,*) ppmco2 read(1,*) ppmch4 read(1,*) ppmn2o close(1) * blank out header initially do k=1,10 write(header(k),'(80a1)') (' ',j=1,80) end do * read 3 lines of header info from gfdl file, number of levels and bands read( iounit,'(a80)') header(1) read( iounit,'(a80)') header(3) write(header(6),'(a)') + ' *************************************************************' write(header(7),'(a)') + ' PLEASE NOTE: WATER VAPOR AND OZONE PROFILES ARE NOT PROVIDED.' write(header(8),'(a)') + ' HEATING RATES ARE REPORTED AT FLUX LEVELS. ' write(header(9),'(a)') + ' *************************************************************' read( iounit,*) nlevel,nbands nlayer=nlevel-1 * read spectral sums (k=0th elements) do k=0,0 read( iounit,*) bandv1(k),bandv2(k) read( iounit,*) (fluxnt(l,k), l = 1,nlevel ) read( iounit,*) (fluxup(l,k), l = 1,nlevel ) read( iounit,*) (fluxdn(l,k), l = 1,nlevel ) read( iounit,*) (heatrt(l,k), l = 1,nlayer ) read( iounit,*) tropnt(k),tropup(k),tropdn(k) end do * read atmospheric profile information read( iounit,'(1x,7e11.5)') ( plevel(l), l = 1,nlevel ) read( iounit,'(1x,7e11.5)') ( tlevel(l), l = 1,nlevel ) read( iounit,'(1x,7e11.5)') ( player(l), l = 1,nlayer ) read( iounit,'(1x,7e11.5)') ( tlayer(l), l = 1,nlayer ) read( iounit,'(1x,7e11.5)') ( tlayer(l), l = 1,nlayer ) read( iounit,'(1x,7e11.5)') ( tlayer(l), l = 1,nlayer ) * find two levels adjoining tropopause if(atmlbl.eq.'trp') ptrop= 93.7 if(atmlbl.eq.'mls') ptrop=179.0 if(atmlbl.eq.'mlw') ptrop=256.8 if(atmlbl.eq.'sas') ptrop=267.7 if(atmlbl.eq.'saw') ptrop=282.9 do j=1,nlevel if(plevel(j).lt.ptrop.and.ptrop.le.plevel(j+1)) then jtrop1=j jtrop2=j+1 endif end do * interpolate fluxes to tropopause do k=0,0 tropnt(k)=( (plevel(jtrop2)-ptrop) * fluxnt(jtrop1,k) + + (ptrop-plevel(jtrop1)) * fluxnt(jtrop2,k)) + / (plevel(jtrop2)-plevel(jtrop1)) tropup(k)=( (plevel(jtrop2)-ptrop) * fluxup(jtrop1,k) + + (ptrop-plevel(jtrop1)) * fluxup(jtrop2,k)) + / (plevel(jtrop2)-plevel(jtrop1)) tropdn(k)=( (plevel(jtrop2)-ptrop) * fluxdn(jtrop1,k) + + (ptrop-plevel(jtrop1)) * fluxdn(jtrop2,k)) + / (plevel(jtrop2)-plevel(jtrop1)) end do * ingest flux and heating rate data band-by-band do k=1,nbands read( iounit,*) bandv1(k),bandv2(k) read( iounit,*) (fluxnt(l,k), l = 1,nlevel ) read( iounit,*) (fluxup(l,k), l = 1,nlevel ) read( iounit,*) (fluxdn(l,k), l = 1,nlevel ) read( iounit,*) (heatrt(l,k), l = 1,nlayer ) read( iounit,*) tropnt(k),tropup(k),tropdn(k) tropnt(k)=( (plevel(jtrop2)-ptrop) * fluxnt(jtrop1,k) + + (ptrop-plevel(jtrop1)) * fluxnt(jtrop2,k)) + / (plevel(jtrop2)-plevel(jtrop1)) tropup(k)=( (plevel(jtrop2)-ptrop) * fluxup(jtrop1,k) + + (ptrop-plevel(jtrop1)) * fluxup(jtrop2,k)) + / (plevel(jtrop2)-plevel(jtrop1)) tropdn(k)=( (plevel(jtrop2)-ptrop) * fluxdn(jtrop1,k) + + (ptrop-plevel(jtrop1)) * fluxdn(jtrop2,k)) + / (plevel(jtrop2)-plevel(jtrop1)) end do * re-define units as ppm, molecules/cm**2 do j=1,nlayer alayer(j) = 1.02 * (6.0256e+23/28.97) * (plevel(j+1)-plevel(j)) wlayer(j) = 1.e+6 * wlayer(j) * 28.97/18.00 olayer(j) = 1.e+6 * olayer(j) * 28.97/48.00 end do * sum column amounts for dry air and 5 gases pthair = 0. pthh2o = 0. ptho3 = 0. do j=1,nlayer pthair = pthair + alayer(j) pthh2o = pthh2o + 1.e-6 * wlayer(j) * alayer(j) ptho3 = ptho3 + 1.e-6 * olayer(j) * alayer(j) end do pthco2 = 1.e-6 * ppmco2 * pthair pthch4 = 1.e-6 * ppmch4 * pthair pthn2o = 1.e-6 * ppmn2o * pthair * miscellaneous psurf = plevel(nlevel) tsurf = tlevel(nlevel) return end subroutine rgsfc + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) *_______________________________________________________________________ * purpose: rgsfc - read gsfc icrccm data in earlier 1988 data format *_______________________________________________________________________ * parameters defined in driver and passed to all subroutines: integer iounit,mlayer,mlevel,mbands * variables passed to or created by subroutines: character*80 header(10) character*3 atmlbl integer iphase,icase,nlayer,nlevel,nbands real psurf,tsurf,ptrop,flgctm real ppmco2,ppmch4,ppmn2o real pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o real bandv1(0:mbands),bandv2(0:mbands) real plevel(mlevel),tlevel(mlevel) real player(mlayer),tlayer(mlayer) real alayer(mlayer),wlayer(mlayer),olayer(mlayer) real fluxup(mlevel,0:mbands) real fluxdn(mlevel,0:mbands) real fluxnt(mlevel,0:mbands) real heatrt(mlayer,0:mbands) real tropup( 0:mbands) real tropdn( 0:mbands) real tropnt( 0:mbands) *_______________________________________________________________________ character*80 csave * blank out header initially do k=1,10 write(header(k),'(80a1)') (' ',j=1,80) end do * read 10 lines of header info from gsfc file read (iounit,'(16x,a80)') header do k=1,10 csave=header(k) csave(1:77)=csave(4:80) header(k)=csave end do * read key parameters which identify this case iphase=1 read (iounit,910) icase,atmlbl,nlayer,nlevel,nbands, + psurf,tsurf,ptrop, + pthair,pthh2o,ptho3,flgctm,ppmco2,ppmch4,ppmn2o 910 format(16x,i16, + /,16x,13x,a3, + 3(/,16x,i16), + 10(/,16x,e16.6)) * ingest spectrally summed results (k=0) and all spectral band data do k=0,nbands read (iounit,'(16x,f16.2)') bandv1(k),bandv2(k) read (iounit,930) fluxup(nlevel,k),fluxdn(nlevel,k), + tropup( k),tropdn( k), + fluxup( 1,k) 930 format(16x,e16.6,4(/,16x,e16.6)) tropnt(k)=tropup(k)-tropdn(k) read (iounit,'(1x)') read (iounit,'(1x)') do j=1,nlayer read (iounit,'(1p,e16.6,3p,e16.4,1p,3e16.6)') + player(j),tlayer(j),wlayer(j),olayer(j),heatrt(j,k) end do read (iounit,'(1x)') read (iounit,'(1x)') do j=1,nlevel read (iounit,'(1p,e16.6,3p,e16.4,1p,3e16.6)') + plevel(j),tlevel(j),fluxup(j,k),fluxdn(j,k),fluxnt(j,k) end do end do * re-define units as ppm, molecules/cm**2 do j=1,nlayer alayer(j) = 1.02 * (6.0256e+23/28.97) * (plevel(j+1)-plevel(j)) wlayer(j) = 1.e+6 * wlayer(j) * 28.97/18.00 olayer(j) = 1.e+6 * olayer(j) * 28.97/48.00 end do * sum column amounts for dry air and 5 gases pthair = 0. pthh2o = 0. ptho3 = 0. do j=1,nlayer pthair = pthair + alayer(j) pthh2o = pthh2o + 1.e-6 * wlayer(j) * alayer(j) ptho3 = ptho3 + 1.e-6 * olayer(j) * alayer(j) end do pthco2 = 1.e-6 * ppmco2 * pthair pthch4 = 1.e-6 * ppmch4 * pthair pthn2o = 1.e-6 * ppmn2o * pthair * icrccm phase, case # over-ridden using separate file open(1,file='temp_header') read(1,*) iphase read(1,*) icase read(1,'(a3)') atmlbl close(1) return end subroutine rfm93 + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) *_______________________________________________________________________ * purpose: rfm93 - read 1993-format icrccm data file, return data arrays *_______________________________________________________________________ * parameters defined in driver and passed to all subroutines: integer iounit,mlayer,mlevel,mbands * variables passed to or created by subroutines: character*80 header(10) character*3 atmlbl integer iphase,icase,nlayer,nlevel,nbands real psurf,tsurf,ptrop,flgctm real ppmco2,ppmch4,ppmn2o real pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o real bandv1(0:mbands),bandv2(0:mbands) real plevel(mlevel),tlevel(mlevel) real player(mlayer),tlayer(mlayer) real alayer(mlayer),wlayer(mlayer),olayer(mlayer) real fluxup(mlevel,0:mbands) real fluxdn(mlevel,0:mbands) real fluxnt(mlevel,0:mbands) real heatrt(mlayer,0:mbands) real tropup( 0:mbands) real tropdn( 0:mbands) real tropnt( 0:mbands) *_______________________________________________________________________ * blank out header initially do k=1,10 write(header(k),'(80a1)') (' ',j=1,80) end do * read 10 lines of header info from fm93 file read (iounit,'(a80)') header * read key parameters which identify this case read (iounit,910) iphase,icase,atmlbl,nlayer,nlevel,nbands, + psurf,tsurf,ptrop,flgctm, + ppmco2,ppmch4,ppmn2o, + pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o 910 format(24x,i16, + /,24x,i16, + /,24x,13x,a3, + 3(/,24x,i16), + 3(/,24x,f16.2), + 1(/,24x,f16.0), + 3(/,24x,f16.2), + 6(/,24x,f16.4)) * ingest spectrally summed results (k=0 only) k=0 read (iounit,'(24x,f16.2)') bandv1(k),bandv2(k) read (iounit,'(24x,e16.4)') fluxup(nlevel,k),fluxdn(nlevel,k), + tropup( k),tropdn( k), + fluxup( 1,k) tropnt(k)=tropup(k)-tropdn(k) read (iounit,'(1x)') read (iounit,'(1x)') do j=1,nlayer read (iounit,'(6e12.4)') + player(j),tlayer(j),alayer(j),wlayer(j),olayer(j),heatrt(j,k) end do read (iounit,'(1x)') read (iounit,'(1x)') do j=1,nlevel read (iounit,'(5e12.4)') + plevel(j),tlevel(j),fluxup(j,k),fluxdn(j,k),fluxnt(j,k) end do * read band-by-band flux summaries, fluxes, heating rates do k=1,nbands read(iounit,'(2f15.2)') bandv1(k),bandv2(k) read(iounit,'(5e15.6)') fluxup(nlevel,k),fluxdn(nlevel,k), + tropup( k),tropdn( k), + fluxup( 1,k) read(iounit,'(5e15.6)') (fluxup(j,k),j=1,nlevel) read(iounit,'(5e15.6)') (fluxdn(j,k),j=1,nlevel) read(iounit,'(5e15.6)') (fluxnt(j,k),j=1,nlevel) read(iounit,'(5e15.6)') (heatrt(j,k),j=1,nlayer) end do return end subroutine wfm93 + (iounit,mlayer,mlevel,mbands + ,header,atmlbl,iphase,icase,nlayer,nlevel,nbands + ,psurf,tsurf,ptrop,flgctm,ppmco2,ppmch4,ppmn2o + ,pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o + ,bandv1,bandv2,plevel,tlevel + ,player,tlayer,alayer,wlayer,olayer + ,fluxup,fluxdn,fluxnt,heatrt + ,tropup,tropdn,tropnt) *_______________________________________________________________________ * purpose: wfm93 - write icrccm data in new standard (1993) format *_______________________________________________________________________ * parameters defined in driver and passed to all subroutines: integer iounit,mlayer,mlevel,mbands * variables passed to or created by subroutines: character*80 header(10) character*3 atmlbl integer iphase,icase,nlayer,nlevel,nbands real psurf,tsurf,ptrop,flgctm real ppmco2,ppmch4,ppmn2o real pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o real bandv1(0:mbands),bandv2(0:mbands) real plevel(mlevel),tlevel(mlevel) real player(mlayer),tlayer(mlayer) real alayer(mlayer),wlayer(mlayer),olayer(mlayer) real fluxup(mlevel,0:mbands) real fluxdn(mlevel,0:mbands) real fluxnt(mlevel,0:mbands) real heatrt(mlayer,0:mbands) real tropup( 0:mbands) real tropdn( 0:mbands) real tropnt( 0:mbands) *_______________________________________________________________________ * print file header with 10 lines of info in any format write(iounit,'(a80)') header * print descriptive parameters write(iounit,910) iphase,icase,atmlbl,nlayer,nlevel,nbands, + psurf,tsurf,ptrop,flgctm, + ppmco2,ppmch4,ppmn2o, + pthair,pthh2o,pthco2,ptho3,pthch4,pthn2o 910 format('.....icrccm.phase.number',i16 + ,/,'......icrccm.case.number',i16 + ,/,'.......atmospheric.label',a16 + ,/,'.....number.model.layers',i16 + ,/,'.....number.model.levels',i16 + ,/,'...number.spectral.bands',i16 + ,/,'........surface.pressure',0p,f16.2,' mb' + ,/,'.....surface.temperature',0p,f16.2,' k' + ,/,'.....tropopause.pressure',0p,f16.2,' mb' + ,/,'......h2o.continuum.flag',0p,f16.0,' (0=off,1=included)' + ,/,'........co2.mixing.ratio',0p,f16.2,' ppm' + ,/,'........ch4.mixing.ratio',0p,f16.2,' ppm' + ,/,'........n2o.mixing.ratio',0p,f16.2,' ppm' + ,/,'....total.dry.air.column',1p,e16.4,' molecules/cm**2' + ,/,'........total.h2o.column',1p,e16.4,' molecules/cm**2' + ,/,'........total.co2.column',1p,e16.4,' molecules/cm**2' + ,/,'.........total.o3.column',1p,e16.4,' molecules/cm**2' + ,/,'........total.ch4.column',1p,e16.4,' molecules/cm**2' + ,/,'........total.n2o.column',1p,e16.4,' molecules/cm**2') * print spectrally itegrated summary of icrccm fluxes k=0 write(iounit,920) bandv1(k),bandv2(k) 920 format('....spectral.region.from',0p,f16.2,' cm-1 ' + ,/,'......................to',0p,f16.2,' cm-1 ') write(iounit,930) fluxup(nlevel,k),fluxdn(nlevel,k), + tropup( k),tropdn( k), + fluxup( 1,k) 930 format('.......flux.up.@.surface',1p,e16.4,' watts/m**2' + ,/,'.....flux.down.@.surface',1p,e16.4,' watts/m**2' + ,/,'....flux.up.@.tropopause',1p,e16.4,' watts/m**2' + ,/,'..flux.down.@.tropopause',1p,e16.4,' watts/m**2' + ,/,'..flux.up.@.top.of.atmos',1p,e16.4,' watts/m**2') * print table of atmospheric profile data and layer heating rates write(iounit,'(2a)') + ' layer >> p t air mass' +,' h2o o3 heating' +,' (mb) (K) (#/cm**2)' +,' (ppm) (ppm) (deg/day)' do j=1,nlayer write(iounit,'(1p,e12.4,3p,e12.4,1p,4e12.4)') + player(j),tlayer(j),alayer(j),wlayer(j),olayer(j),heatrt(j,k) end do * print table of computed fluxes write(iounit,'(a)') + ' level >> p t flux up flux down net flux' +,' (mb) (k) (w/m**2) (w/m**2) (w/m**2)' do j=1,nlevel write(iounit,'(1p,e12.4,3p,e12.4,1p,3e12.4)') + plevel(j),tlevel(j),fluxup(j,k),fluxdn(j,k),fluxnt(j,k) end do * print band-by-band flux summaries, fluxes, heating rates do k=1,nbands write(iounit,'(0p,2f15.2)') bandv1(k),bandv2(k) write(iounit,'(1p,5e15.6)') fluxup(nlevel,k),fluxdn(nlevel,k), + tropup( k),tropdn( k), + fluxup( 1,k) write(iounit,'(1p,5e15.6)') (fluxup(j,k),j=1,nlevel) write(iounit,'(1p,5e15.6)') (fluxdn(j,k),j=1,nlevel) write(iounit,'(1p,5e15.6)') (fluxnt(j,k),j=1,nlevel) write(iounit,'(1p,5e15.6)') (heatrt(j,k),j=1,nlayer) end do return end