c**** RDCHMRT c subroutine rdchmrt(version) use filunit use grid use ptemiss use chmstry use pigsty use rtracchm use tracer c c----CAMx v7Beta6 190902 c c----------------------------------------------------------------------- c Description: c----------------------------------------------------------------------- c c This routine reads the chemical parameters file for RTRAC c technology of the source apportionment algorithm. c c Copyright 1996 - 2018 c Ramboll c c Argument description: c Input: c version C model version character string c c----------------------------------------------------------------------- c LOG: c----------------------------------------------------------------------- c c 01/16/02 --gwilson-- Original development c 11/12/03 --cemery-- Added secondary RTRAC species to primary c RTRAC species c 10/29/09 --cemery-- Added RTRAC surface model parameters c c----------------------------------------------------------------------- c Include files: c----------------------------------------------------------------------- c include 'camx.prm' include 'flags.inc' include 'rtracsrf.inc' c c----------------------------------------------------------------------- c External functions: c----------------------------------------------------------------------- c integer istrln c c----------------------------------------------------------------------- c Argument declaration: c----------------------------------------------------------------------- c character*20 version c c----------------------------------------------------------------------- c Local variables: c----------------------------------------------------------------------- c integer mvecedge real rkoh, rko3, rkno3, adj character*200 line character*80 action character*10 camxv, camxvin, sptmp c character*6 flag(MXTRSP) character*6 sectmp(MXTRSP) c c----------------------------------------------------------------------- c Data statements: c----------------------------------------------------------------------- c data camxv /'VERSION7BE'/ c c----------------------------------------------------------------------- c Entry point: c----------------------------------------------------------------------- c c --- initalize some variables --- c lrtgas = .TRUE. lrtaero = .TRUE. write(idiag,'(//,A,//)',ERR=7009) & ' ******* Reading RTRAC chemical definitions file.' c c --- open the file ---- c open(unit=iorchm,file=chmfil,ERR=7010,status='OLD') c c --- read and verify the version number --- c action = 'Reading version number in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line camxvin = line(21:30) write(idiag,'(A)',ERR=7009) line(:istrln(line)) call jstlft( camxvin ) call toupper( camxvin ) if (camxvin .NE. camxv ) then write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(/,A)') & ' CAMx version in RTRAC chemistry file INVALID' write(iout,'(2A)') ' Expecting: ',camxv write(iout,'(2A)') ' Found: ',camxvin write(iout,'(1X,3A)') 'Make sure you are using an RTRAC ', & 'file designed for ',version call camxerr() endif c c --- ignore the description record ---- c action = 'Reading description record in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) c c --- read the number of gas species ---- c action = 'Reading number of gas species in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line read(line,'(20X,I10)',ERR=7005) nrtgas write(idiag,'(A,I5)',ERR=7009) line(1:20), nrtgas if( nrtgas .LE. 0 ) lrtgas = .FALSE. c c --- read the number of aero species ---- c action = & 'Reading number of aerosol species in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line read(line,'(20X,I10)',ERR=7005) nrtaero write(idiag,'(A,I5)',ERR=7009) line(1:20), nrtaero if (lparttn .AND. nrtaero.lt.nrtgas) then write(iout,'(A)') 'RTRAC gas-aerosol partitioning is selected.' write(iout,'(A,A)') 'You must list at least the same number', & ' of aerosol tracers as gas tracers.' write(iout,'(A,A)') 'NOTE: gas tracer 1 will partition to', & ' aerosol tracer 1, etc.' write(iout,'(A)') 'CAMx is stopping.' call camxerr() endif if( nrtaero .LE. 0 ) lrtaero = .FALSE. c c --- read the number of photolysis reactions ---- c action = & 'Reading number of photolysis reactions in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line read(line,'(20X,I10)',ERR=7005) nrtphot write(idiag,'(A,I5)',ERR=7009) line(1:20), nrtphot c c --- read the number of thermal reactions ---- c action = & 'Reading number of thermal reactions in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line read(line,'(20X,I10)',ERR=7005) nrtherm write(idiag,'(A,I5)',ERR=7009) line(1:20), nrtherm c c --- make sure that there is enough array space --- c ntotsp = nrtgas + nrtaero nrtrac = nrtgas + nrtaero nsaspc = ntotsp notimespc = ntotsp c c --- make sure the arrays are allocated large enough --- c if( ntotsp .GT. MXTRSP ) then write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,*) 'A parameter in the camx.prm is not ', & 'sufficiently large.' write(iout,*) 'Please change the value for parameter: MXTRSP' write(iout,*) 'It should be set to a value of at least: ',ntotsp call flush(iout) call camxerr() endif c c --- call routine to allocate the arrays --- c mvecedge = MAX(ncol(1),nrow(1)) call alloc_tracer_specs(ngrid,ncol,nrow,nlay,nlayers_ems, & tectyp,mvecedge,iout) call alloc_rtracchm(ngrid,ncol,nrow,nspec,ntotsp) call alloc_tracer_vdep(ngrid,ncol,nrow,ntotsp) c c --- if doing sampling grids, allocate the arrays --- c if( lsmptrc ) & call alloc_tracer_sample(nsample,ncolsmp,nrowsmp,nrtsmpcels) c --- initialize the variables --- c do i = 1,ntotsp jnum(i) = 0 rtjfact(i) = 0. aoh(i) = 0. eaoh(i) = 0. boh(i) = 0. troh(i) = 0. ano3(i) = 0. eano3(i) = 0. bno3(i) = 0. trno3(i) = 0. ao3(i) = 0. eao3(i) = 0. bo3(i) = 0. tro3(i) = 0. flag(i) = ' ' sectmp(i) = ' ' eqkoa(i) = 0. khydro(i) = 0. kleach(i) = 0. kpen(i) = 0. kphot(i) = 0. enddo c c --- read the gas species identification records --- c if( nrtgas .GT. 0 ) then action = 'Reading gas species identification'// & ' record in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) do i = 1,nrtgas read(iorchm,'(A)',ERR=7000,END=7001) line c c --- parse the line ---- c sptmp = line(6:15) call jstlft( sptmp ) flag(i) = line(16:21) call jstlft( flag(i) ) call toupper( flag(i) ) if( flag(i) .EQ. 'PRIM ' ) then lsecnd(i) = .FALSE. else if( flag(i) .EQ. 'SEC ' ) then lsecnd(i) = .TRUE. sectmp(i) = line(22:27) call jstlft( sectmp(i) ) if( ipigflg .EQ. IRONPIG ) then idx = 0 lreg(i) = .FALSE. do ispc=1,i-1 if( ptname(ispc)(1:6) .EQ. sectmp(i) ) idx = ispc enddo if( idx .EQ. 0 ) goto 7011 ksec(i) = idx else idx = 0 lreg(i) = .FALSE. do ispc=1,nspec if( spname(ispc)(1:6) .EQ. sectmp(i) ) then idx = ispc lreg(i) = .TRUE. ksec(i) = idx endif enddo if( .NOT. lreg(i) ) then do ispc=1,i-1 if( ptname(ispc)(1:6) .EQ. sectmp(i) ) then idx = ispc ksec(i) = idx endif enddo endif if( idx .EQ. 0 ) goto 7004 endif else goto 7003 endif c c ---- passes all checks, load into arrays --- c ptname(i) = sptmp lsamap(i) = i lsagas(i) = .TRUE. loutsa(i) = .TRUE. read(line,'(27X,6F10.0)',ERR=7005) rtlbnd(i), & rthlaw(i), rttfact(i), rtdrate(i), rtreact(i), rtscale(i) enddo c c --- echo the data to the diag file ---- c do i = 1,nrtgas write(idiag,'(I3,2X,A10,2A6,2E10.3,4F10.3)',ERR=7009) i, & ptname(i),flag(i),sectmp(i),rtlbnd(i), rthlaw(i), & rttfact(i), rtdrate(i),rtreact(i),rtscale(i) enddo endif c c --- read the soil model data records c if (lsrfmodrt) then action = 'Reading surface model species records'// & ' in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) do i = 1,nrtgas read(iorchm,'(A)',ERR=7000,END=7001) line sptmp = line(6:15) call jstlft( sptmp ) if (sptmp .NE. ptname(i)) goto 7020 read(line,'(15X,5F10.0)',ERR=7005) & eqkoa(i),khydro(i),kleach(i),kpen(i),kphot(i) write(idiag,'(I3,2X,A10,5E10.3)',ERR=7009) i, & ptname(i),eqkoa(i),khydro(i),kleach(i),kpen(i),kphot(i) if (eqkoa(i).lt.0. .or. khydro(i).lt.0. .or. kleach(i).lt.0. & .or. kpen(i).lt.0. .or. kphot(i).lt.0.) then write(iout,'(A)') 'Rates for this species are <0!' write(iout,'(A)') 'Only positive rates are allowed --- Stopping' call camxerr() endif enddo endif c c --- read the aerosol species identification records --- c if( nrtaero .GT. 0 ) then action = 'Reading aerosol species identification'// & ' record in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) do j = 1,nrtaero i = nrtgas+j read(iorchm,'(A)',ERR=7000,END=7001) line c c --- parse the line ---- c sptmp = line(6:15) call jstlft( sptmp ) ptname(i) = sptmp lsamap(i) = i lsagas(i) = .FALSE. loutsa(i) = .TRUE. read(line,'(16X,4F10.0)',ERR=7005) rtlbnd(i), & rtdens(i), rtlcut(i), rtucut(i) enddo c c --- echo the data to the diag file ---- c do j = 1,nrtaero i = nrtgas+j write(idiag,'(I3,2X,A10,E10.3,3F10.4)',ERR=7009) i, & ptname(i), rtlbnd(i), rtdens(i), & rtlcut(i), rtucut(i) rtdens(i) = rtdens(i)*1.e6 enddo endif c c --- read the phololysis reactions --- c if( nrtphot .GT. 0 ) then action = 'Reading photolysis reaction'// & ' record in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) do irec=1,nrtphot read(iorchm,'(A)',ERR=7000,END=7001) line sptmp = line(1:10) call jstlft( sptmp ) idx = 0 do i = 1,nrtgas if( ptname(i) .EQ. sptmp ) idx = i enddo if( idx .EQ. 0 ) goto 7006 read(line,'(10X,I10,F10.0)',ERR=7005) jnum(idx), & rtjfact(idx) iphot = 0 do i = 1,nphot1 if( jnum(idx) .EQ. idphot1(i) ) iphot = i enddo do i = 1,nphot2 if( jnum(idx) .EQ. idphot2(i) .OR. & jnum(idx) .EQ. idphot3(i) ) iphot = i enddo if( iphot .EQ. 0 ) goto 7007 enddo c c --- echo data to diag file ---- c do i = 1,nrtgas write(idiag,'(A10,I5,F10.4)',ERR=7009) ptname(i), & jnum(i), rtjfact(i) enddo endif c c --- read the thermal reactions --- c if( nrtherm .GT. 0 ) then action = 'Reading thermal reaction'// & ' record in RTRAC Chemistry file.' read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) read(iorchm,'(A)',ERR=7000,END=7001) line write(idiag,'(A)',ERR=7009) line(:istrln(line)) do irec=1,nrtherm read(iorchm,'(A)',ERR=7000,END=7001) line sptmp = line(1:10) call jstlft( sptmp ) idx = 0 do i = 1,nrtgas if( ptname(i) .EQ. sptmp ) idx = i enddo if( idx .EQ. 0 ) goto 7006 c sptmp = line(11:16) call jstlft( sptmp ) if( sptmp .EQ. 'OH ' ) then read(line,'(16X,2F12.0,2F10.0)',ERR=7005) aoh(idx), & eaoh(idx), boh(idx), troh(idx) else if( sptmp .EQ. 'NO3 ' ) then read(line,'(16x,2F12.0,2F10.0)',ERR=7005) ano3(idx), & eano3(idx), bno3(idx), trno3(idx) else if( sptmp .EQ. 'O3 ' ) then read(line,'(16x,2F12.0,2F10.0)',ERR=7005) ao3(idx), & eao3(idx), bo3(idx), tro3(idx) else goto 7008 endif enddo c c --- echo data to diag file ---- c do i = 1,nrtgas write(idiag,'(A10,A6,2E12.4,2F10.4)',ERR=7009) ptname(i), & ' OH ',aoh(i), eaoh(i), boh(i), troh(i) write(idiag,'(A10,A6,2E12.4,2F10.4)',ERR=7009) ptname(i), & ' NO3 ',ano3(i), eano3(i), bno3(i), trno3(i) write(idiag,'(A10,A6,2E12.4,2F10.4)',ERR=7009) ptname(i), & ' O3 ',ao3(i), eao3(i), bo3(i), tro3(i) enddo c c-----Provide diagnostic info for checking rate expressions c See CHEMRT for explanation of function ARRHEN c write(idiag,'(/,a)') & 'Diagnostic info for checking RTRAC rate expressions' write(idiag,'(/,a,/,a)') & 'Rates at 298 K and 1013 mbar in ppm-n min-1', & ' Species OH rate O3 rate NO3 rate' adj = 1. do i = 1,nrtgas rkoh = arrhen(aoh(i),eaoh(i),boh(i),troh(i),298.)*adj rko3 = arrhen(ao3(i),eao3(i),bo3(i),tro3(i),298.)*adj rkno3 = arrhen(ano3(i),eano3(i),bno3(i),trno3(i),298.)*adj write(idiag,'(I5,1X,A10,1P3E12.4)') i, ptname(i), & rkoh, rko3, rkno3 enddo write(idiag,'(/,a,/,a)') & 'Rates at 273 K and 1013 mbar in ppm-n min-1', & ' Species OH rate O3 rate NO3 rate' adj = 298./273. do i = 1,nrtgas rkoh = arrhen(aoh(i),eaoh(i),boh(i),troh(i),273.)*adj rko3 = arrhen(ao3(i),eao3(i),bo3(i),tro3(i),273.)*adj rkno3 = arrhen(ano3(i),eano3(i),bno3(i),trno3(i),273.)*adj write(idiag,'(I5,1X,A10,1P3E12.4)') i, ptname(i), & rkoh, rko3, rkno3 enddo endif c c --- return to calling routine --- c goto 9999 c c----------------------------------------------------------------------- c Error messages: c----------------------------------------------------------------------- c 7000 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) call camxerr() c 7001 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') 'Premature end-of-file reached.' write(iout,'(1X,A)') action(:istrln(action)) call camxerr() c 7003 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,2A)') 'Invalid keyword for species type: ', & line(16:21) write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) call camxerr() c 7004 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,3A)') 'Specified primary species: ',line(22:27), & ' is not in regular model species list' write(iout,'(1X,A)') 'nor in the RTRAC species list read so far.' write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) write(iout,'(1X,2A)') 'Any child RTRAC tracer must be listed after ', & 'its parent RTRAC tracer.' call camxerr() c 7005 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) call camxerr() c 7006 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,3A)') 'Species listed not included in ', & 'list of gas species definitions.' write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) call camxerr() c 7007 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,2A,I3)') 'Reaction # listed is not a ', & 'photolysis reaction. Rxn Id: ',jnum(idx) write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) call camxerr() c 7008 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,3A)') 'Invalid species for reaction: ', & line(11:16) write(iout,'(1X,A)') 'Must be one of the following: ' write(iout,'(10X,A)') 'OH' write(iout,'(10X,A)') 'NO3' write(iout,'(10X,A)') 'O3' write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) call camxerr() c 7009 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') 'Cannot write to the diag file.' call camxerr() c 7010 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') 'Cannot open file: ',chmfil(:istrln(chmfil)) call camxerr() c 7011 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,3A)') 'Specified primary species: ',line(22:27), & ' is not in the RTRAC species list read so far.' write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) write(iout,'(1X,2A)') 'Any child RTRAC tracer must be listed after', & ' its parent RTRAC tracer.' write(iout,'(1X,2A,/,1X,A)') 'Secondary RTRAC tracers can only map to', & ' other RTRAC tracers if IRON PiG',' is invoked.' call camxerr() c 7020 continue write(iout,'(//,A)') 'ERROR in RDCHMRT:' write(iout,'(1X,A)') action(:istrln(action)) write(iout,'(1X,A,A,A)') 'Surface model gas species: ',sptmp, & ' does not match the RTRAC species list.' write(iout,'(1X,A,/,A)') 'Line read: ',line(:istrln(line)) write(iout,'(1X,A,A)') 'Expecting species: ',ptname(i) write(iout,'(1X,2A,/,1X,A)') 'The surface model parameters must be', & ' given in the same order as the RTRAC',' gas species list' call camxerr() c c----------------------------------------------------------------------- c Return point: c----------------------------------------------------------------------- c 9999 continue c write(idiag,'(//,A,//)',ERR=7009) & ' ******* Finished RTRAC chemical definitions file.' close(iorchm) return end