! module restart use params,only: nlev use flds_hist,only: nf_hist,f_hist,nf_restart,fnames_restart, | nf_diaghist,fnames_diaghist implicit none #include ! ! Uars solstice flux data (formerly in uarsflx.h): integer,parameter :: nbin_uars=170 real,dimension(nbin_uars) :: | wv_low, wv_high, sminflx, smedflx, smaxflx ! ! Formerly in sftrf.h: integer,parameter :: nsftrf=158 real,dimension(nsftrf) :: | vl1 ,vl2 ,sfx ,sgry ,sgo2hz ,sgo31 ,sgo32 integer,parameter :: nsfdim1=17, nsfdima=9, nsfdimb=5 real,dimension(nsfdim1,nsfdima) :: aij real,dimension(nsfdim1,nsfdimb) :: bij real,dimension(nsfdima) :: ad00, ad10 real,dimension(nsfdimb) :: bd00, bd10 ! ! Formerly in atmosp.h: integer, parameter :: natm=71,nfatm=8 real :: atm(natm,nfatm) ! real,dimension(nlev) :: qcool contains !----------------------------------------------------------------------- subroutine read_source ! ! Read source history netcdf file: ! use input,only: source,start_year,start_day use params,only: nlev use flds_ionatm,only: xiop,xiop4s,xiop2p,xiop2d use flds_atmos,only: xncl,xnclo,xno2,xnn2,xno,tn,xno21d,xno21s use flds_sodium,only: xnapco2,xnaph2o,xnaop use flds_modelz,only: zp implicit none ! ! Local: integer :: i,k,kk,ncid,istat,ndims,nvars,ngatts,idunlim,len, | ihist integer :: ilev=0,ibin_uars=0,isftrf=0,isfdim1=0,isfdima=0, | isfdimb=0,iatm=0,ifatm=0 character(len=80) :: varname character(len=120) :: char120 logical :: exists real :: xmm real,dimension(nlev) :: alamm,presmb,un,vn ! real :: rmcl(nlev),rmclo(nlev) real :: smcl(45) = | (/1.7E-13, 1.4E-12, 4.0E-12, 6.8E-12, 2.1E-11, 3.0E-11, 3.8E-11, | 5.2E-11, 6.6E-11, 7.5E-11, 8.3E-11, 8.9E-11, 8.2E-11, 9.0E-11, | 1.E-10, 2.1E-10, 3.3E-10, 1.4E-9, 1.8E-9, 2.2E-9, 2.3E-9, | 2.4E-9, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, | 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, | 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, | 1.0e-20, 1.0e-20, 1.0e-20/) real :: smclo(45) = | (/1.9E-10, 6.0E-10, 7.1E-10, 8.6E-10, 5.2E-10, 3.5E-10, 2.0E-10, | 7.4E-11, 3.0E-11, 2.0E-11, 8.0E-12, 3.6E-12,8.3E-13, 6.0E-13, | 1.9E-13, 1.0E-13, 8.0E-14, 2.3E-14, 2.0E-14,2.9E-15, 2.0E-15, | 1.2E-15, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, | 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, | 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, 1.0e-20, | 1.0e-20, 1.0e-20, 1.0e-20/) ! ! File must be on the disk (not going to mss at this time): write(6,"(/,72('-'))") write(6,"('Enter read_source: source=',a)") trim(source) call expand_path(source) ! expand any env vars imbedded in source path inquire(file=source,exist=exists) if (.not.exists) then write(6,"('>>> read_source: cannot find SOURCE history file ', | a)") trim(source) stop 'source' endif ! ! Open netcdf file: istat = nf_open(source,NF_NOWRITE,ncid) if (istat /= NF_NOERR) then write(char120,"('Error return from nf_open for netcdf ', | 'file ',a,' istat=',a)") trim(source),istat call handle_ncerr(istat,char120) stop 'readfile' endif write(6,"('read_source: Opened netcdf file ',a,' ncid=',i6)") | trim(source),ncid ! ! Get number of dims, vars, atts, and id of unlimited dimension. istat = nf_inq(ncid,ndims,nvars,ngatts,idunlim) ! write(6,"('readfile: ndims=',i3,' nvars=',i3,' ngatts=',i3, ! | ' idunlim=',i3)") ndims,nvars,ngatts,idunlim ! ! Loop over number of dimension on the file (may be only one): len = -1 do i=1,ndims istat = nf_inq_dim(ncid,i,varname,len) ! ! Check number of levels dimension against nlev: if (trim(varname)=='lev') then if (len /= nlev) then write(6,"('>>> readfile: nlev=',i3,' but lev dimension ', | 'read from file = ',i3)") nlev,len stop 'nlev' endif ilev = 1 endif ! ! Check additional dimensions: if (trim(varname)=='nbin_uars') then if (len /= nbin_uars) then write(6,"('>>> readfile: len on file=',i4, | ' but nbin_uars=',i4)") len,nbin_uars stop 'nbin_uars' endif ibin_uars = 1 endif if (trim(varname)=='nsftrf') then if (len /= nsftrf) then write(6,"('>>> readfile: len on file=',i4, | ' but nsftrf=',i4)") len,nsftrf stop 'nsftrf' endif isftrf = 1 endif if (trim(varname)=='nsfdim1') then if (len /= nsfdim1) then write(6,"('>>> readfile: len on file=',i4, | ' but nsfdim1=',i4)") len,nsfdim1 stop 'nsfdim1' endif isfdim1 = 1 endif if (trim(varname)=='nsfdima') then if (len /= nsfdima) then write(6,"('>>> readfile: len on file=',i4, | ' but nsfdima=',i4)") len,nsfdima stop 'nsfdima' endif isfdima = 1 endif if (trim(varname)=='nsfdimb') then if (len /= nsfdimb) then write(6,"('>>> readfile: len on file=',i4, | ' but nsfdimb=',i4)") len,nsfdimb stop 'nsfdimb' endif isfdimb = 1 endif enddo ! i=1,ndims if (ilev==0) then write(6,"('>>> readfile: could not find dimension ''lev''', | ' on source file ',a)") trim(source) endif if (ibin_uars==0) then write(6,"('>>> readfile: could not find dimension ', | '''nbin_uars'' on source file ',a)") trim(source) endif if (isftrf==0) then write(6,"('>>> readfile: could not find dimension ', | '''nsftrf'' on source file ',a)") trim(source) endif if (isfdim1==0) then write(6,"('>>> readfile: could not find dimension ', | '''nsfdim1'' on source file ',a)") trim(source) endif if (isfdima==0) then write(6,"('>>> readfile: could not find dimension ', | '''nsfdima'' on source file ',a)") trim(source) endif if (isfdimb==0) then write(6,"('>>> readfile: could not find dimension ', | '''nsfdimb'' on source file ',a)") trim(source) endif if (ilev==0.or.ibin_uars==0.or.isftrf==0.or.isfdim1==0.or. | isfdima==0.or.isfdimb==0) stop 'dimensions' ! call findhist(ncid,start_year,start_day,ihist) if (ihist < 0) then write(6,"('Note read_source: error from findhist -- ', | 'assuming ihist = 1.')") ihist = 1 elseif (ihist == 0) then write(6,"('>>> read_source: error searching for year=', | i4,' day=',i4,' on history file ',a)") start_year, | start_day,trim(source) stop 'findhist' else ! write(6,"('read_source: Found start_year=',i4,' start_day=', ! | i4,' ihist=',i4)") start_year,start_day,ihist endif ! ! Read restart fields: do i=1,nf_restart call rdvar(fnames_restart(i),ncid,i,ihist) enddo ! i=1,nf_hist ! ! Read history diags: ! (fields on the history that are not calculated in the model) do i=1,nf_diaghist call rdvar(fnames_diaghist(i),ncid,i,ihist) enddo ! i=1,nf_hist ! rmcl = 0. ! whole-array init to avoid indef (not in old code) rmclo = 0. ! whole-array init to avoid indef (not in old code) do k=1,45 ! was do 46 in old restart kk = 1+(k-1)*2 rmcl(kk) = smcl(k) rmclo(kk) = smclo(k) enddo ! k=1,45 do k=1,44 ! was do 47 in old restart kk = 2+(k-1)*2 rmcl(kk) = (smcl (k+1)+smcl (k))*0.5 rmclo(kk) = (smclo(k+1)+smclo(k))*0.5 enddo ! k=1,44 ! ! Present day temperatue at 30 mb ! tn(1)=225. ! CO2 doubled yr=2000 temperature at 30 mb tn(1)=215. ! CO2 halved ice age temperature at 30 mb ! tn(1)=235. ! tn(1)=tnms(1)+10. ! tn(1)=tnms(1)-10. ! tn(1)=tnms(1) ! do k=1,nlev ! was do 944 in old restart xiop4s(k)=xiop(k) xiop2p(k)=1.e-20 xiop2d(k)=1.e-20 xmm=xno2(k)+xnn2(k)+xno(k) ! Base case xncl(k) =rmcl(k) *xmm*1.0 xnclo(k)=rmclo(k)*xmm*1.0 ! 2100 GLOBAL CHANGE CASE ! xncl(k) =rmcl(k) *xmm*1.5 ! xnclo(k)=rmclo(k)*xmm*1.5 ! ! Below inits were in old glbmean.F (after read source): xno21d(k) = 1.e-20 ! flds_atmos xno21s(k) = 1.e-20 ! flds_atmos qcool(k) = 1.e-20 ! local alamm(k) = 0. ! local presmb(k) = 5.e-7*exp(-zp(k)) ! local un(k) = 0. vn(k) = 0. xnapco2(k)=1.e-20 ! sodium xnaph2o(k)=1.e-20 ! sodium xnaop(k)=1.e-20 ! sodium enddo ! ! Read additional fields. In previous versions of glbmean, these were ! read from stdin. In this version, they were added to the netcdf ! history file (see ~foster/tgcm/glbmean/appendhist). ! ! Get additional fields: write(6,"(' ')") call get1dvar('wv_low' ,wv_low ,nbin_uars,ncid) call get1dvar('wv_high',wv_high,nbin_uars,ncid) call get1dvar('sminflx',sminflx,nbin_uars,ncid) call get1dvar('smedflx',smedflx,nbin_uars,ncid) call get1dvar('smaxflx',smaxflx,nbin_uars,ncid) call get1dvar('vl1',vl1,nsftrf,ncid) call get1dvar('vl2',vl2,nsftrf,ncid) call get1dvar('sfx',sfx,nsftrf,ncid) call get1dvar('sgry',sgry,nsftrf,ncid) call get1dvar('sgo2hz',sgo2hz,nsftrf,ncid) call get1dvar('sgo31',sgo31,nsftrf,ncid) call get1dvar('sgo32',sgo32,nsftrf,ncid) call get2dvar('aij',aij,nsfdim1,nsfdima,ncid) call get2dvar('bij',bij,nsfdim1,nsfdimb,ncid) call get1dvar('ad00',ad00,nsfdima,ncid) call get1dvar('ad10',ad10,nsfdima,ncid) call get1dvar('bd00',bd00,nsfdimb,ncid) call get1dvar('bd10',bd10,nsfdimb,ncid) call get2dvar('atm',atm,natm,nfatm,ncid) ! istat = nf_close(ncid) write(6,"('Finished reading source history ',a)") trim(source) write(6,"(72('-'),/)") end subroutine read_source !------------------------------------------------------------------- subroutine get1dvar(fname,array,idim1,ncid) ! ! Args: character(len=*),intent(in) :: fname integer,intent(in) :: idim1,ncid real,intent(out) :: array(idim1) ! ! Local: integer :: istat,start(1),count(1),idvar character(len=120) :: char120 ! ! Get var id: istat = nf_inq_varid(ncid,fname,idvar) if (istat /= NF_NOERR) then write(6,"('Error return from nf_inq_varid for var ',a )") | trim(fname) stop 'get1dvar' endif ! ! Read var: start(1) = 1 count(1) = idim1 istat = nf_get_vara_double(ncid,idvar,start,count,array) if (istat /= NF_NOERR) then write(char120,"('>>> get1dvar: Error return from ', | 'nf_get_vara_double for field ',a)") trim(fname) call handle_ncerr(istat,char120) return endif write(6,"('Read 1-d field ',a,' min,max=',2e12.4)") | fname,minval(array),maxval(array) end subroutine get1dvar !------------------------------------------------------------------- subroutine get2dvar(fname,array,idim1,idim2,ncid) ! ! Args: character(len=*),intent(in) :: fname integer,intent(in) :: idim1,idim2,ncid real,intent(out) :: array(idim1,idim2) ! ! Local: integer :: istat,start(2),count(2),idvar character(len=120) :: char120 ! ! Get var id: istat = nf_inq_varid(ncid,fname,idvar) if (istat /= NF_NOERR) then write(6,"('Error return from nf_inq_varid for var ',a )") | trim(fname) stop 'get2dvar' endif ! ! Read var: start(1) = 1 start(2) = 1 count(1) = idim1 count(2) = idim2 istat = nf_get_vara_double(ncid,idvar,start,count,array) if (istat /= NF_NOERR) then write(char120,"('>>> get2dvar: Error return from ', | 'nf_get_vara_double for field ',a)") trim(fname) call handle_ncerr(istat,char120) return endif write(6,"('Read array ',a,' min,max=',2e12.4)") | fname,minval(array),maxval(array) end subroutine get2dvar !------------------------------------------------------------------- subroutine rdvar(fname,ncid,ifld,ihist) ! ! Read field fname(nlev) from netcdf file attached to ncid. ! use cons,only: zpcons use flds_modelz use flds_msisatm use flds_atmos use flds_ionatm use flds_ionzrt use flds_dissoc use flds_heat use flds_sodium implicit none ! ! Args: character(len=*),intent(in) :: fname integer,intent(in) :: ncid,ifld,ihist ! ! Local: integer :: istat,start(1),count(1),idvar,k,ndims integer :: start2d(2),count2d(2) character(len=120) :: char120 real :: rdfield(nlev),rdfield_2d(nlev,1) real :: edyrd(nlev) ! edy is set by init_eddy (cons.F) ! ! Get variable id: istat = nf_inq_varid(ncid,fname,idvar) if (istat /= NF_NOERR) then write(char120,"('Error return from nf_inq_varid', + ' for field ',a)") trim(fname) call handle_ncerr(istat,char120) return endif istat = nf_inq_varndims(ncid,idvar,ndims) if (ndims < 1.or.ndims > 2) then write(6,"('>>> rdvar: not reading field ',a,' because ndims=', | i3)") trim(fname),ndims return endif ! ! Get variable (nlev) from the file: if (ndims==1) then call get1dvar(fname,rdfield,nlev,ncid) else ! ! Get 2-d variable (nlev,itime) for current time: start2d(1) = 1 count2d(1) = nlev start2d(2) = ihist ! index to desired starting time count2d(2) = 1 istat =nf_get_vara_double(ncid,idvar,start2d,count2d,rdfield_2d) if (istat /= NF_NOERR) then write(char120,"('>>> rdvar: Error return from ', | 'nf_get_vara_double for 2-d field ',a,' nlev=',i3,' ihist=', | i4)") trim(fname),nlev call handle_ncerr(istat,char120) return endif write(6,"('Read 2-d field ',i4,': ',a,' (ihist=',i4,')', | ' min,max=',2e12.4)") ifld,fname,ihist, | minval(rdfield_2d),maxval(rdfield_2d) rdfield(:) = rdfield_2d(:,1) endif ! if (trim(fname)=='TE') then ! write(6,"('rdvar: ihist=',i4,' TE=',/,(4e16.8))") ! | ihist,rdfield ! endif ! ! Transfer field read to appropriate array. ! Arrays are use-associated from the appropriate fields module ! (see fmodules.F). ! select case(fname) ! Field name Array Module case ('ZP ') ; zp = rdfield ! flds_modelz case ('ZPMS ') ; zpms = rdfield ! flds_msisatm case ('ZPHT ') ; zpht = rdfield ! flds_modelz case ('ZPMHT ') ; zpmht = rdfield ! flds_msisatm case ('TN ') ; tn = rdfield ! flds_atmos case ('TNMS ') ; tnms = rdfield ! flds_msisatm case ('XNO ') ; xno = rdfield ! flds_atmos case ('XOM ') ; xom = rdfield ! flds_msisatm case ('XNO1D ') ; xno1d = rdfield ! flds_atmos case ('XNN2 ') ; xnn2 = rdfield ! flds_atmos case ('XN2M ') ; xn2m = rdfield ! flds_msisatm case ('XNO2 ') ; xno2 = rdfield ! flds_atmos case ('XO2M ') ; xo2m = rdfield ! flds_msisatm case ('SHT ') ; sht = rdfield ! flds_atmos case ('SHTMS ') ; shtms = rdfield ! flds_msisatm case ('AMAS ') ; amas = rdfield ! flds_atmos case ('AMASS ') ; amass = rdfield ! flds_msisatm case ('RHO ') ; rho = rdfield ! flds_atmos case ('RHOMS ') ; rhoms = rdfield ! flds_msisatm case ('TE ') ; te = rdfield ! flds_ionatm case ('TI ') ; ti = rdfield ! flds_ionatm case ('XNE ') ; xne = rdfield ! flds_ionatm case ('XNEEE ') ; xneee = rdfield ! flds_ionatm case ('XNPI ') ; xnpi = rdfield ! flds_ionatm case ('XNNI ') ; xnni = rdfield ! flds_ionatm case ('XIOP ') ; xiop = rdfield ! flds_ionatm case ('XINOP ') ; xinop = rdfield ! flds_ionatm case ('XINP ') ; xinp = rdfield ! flds_ionatm case ('XIN2P ') ; xin2p = rdfield ! flds_ionatm case ('XIO2P ') ; xio2p = rdfield ! flds_ionatm case ('XN2D ') ; xn2d = rdfield ! flds_atmos case ('XN4S ') ; xn4s = rdfield ! flds_atmos case ('XN4SM ') ; xn4sm = rdfield ! flds_msisatm case ('XNNO ') ; xnno = rdfield ! flds_atmos case ('XNNO2 ') ; xnno2 = rdfield ! flds_atmos case ('XNNOZ ') ; xnnoz = rdfield ! flds_atmos case ('XNHE ') ; xnhe = rdfield ! flds_atmos case ('XNHEM ') ; xnhem = rdfield ! flds_msisatm case ('XNARG ') ; xnarg = rdfield ! flds_atmos case ('XNARGM ') ; xnargm = rdfield ! flds_msisatm case ('XNH2O ') ; xnh2o = rdfield ! flds_atmos case ('XNH2 ') ; xnh2 = rdfield ! flds_atmos case ('XNCH4 ') ; xnch4 = rdfield ! flds_atmos case ('XNOH ') ; xnoh = rdfield ! flds_atmos case ('XNHO2 ') ; xnho2 = rdfield ! flds_atmos case ('XNH2O2 ') ; xnh2o2 = rdfield ! flds_atmos case ('XNH ') ; xnh = rdfield ! flds_atmos case ('XNHM ') ; xnhm = rdfield ! flds_msisatm case ('XNHTOT ') ; xnhtot = rdfield ! flds_atmos case ('XNCO ') ; xnco = rdfield ! flds_atmos case ('XNCO2 ') ; xnco2 = rdfield ! flds_atmos case ('QNTOT ') ; qntot = rdfield ! flds_heat case ('XIRCOOL ') ; xircool= rdfield ! flds_heat case ('DO2T ') ; do2t = rdfield ! flds_dissoc case ('DO3T ') ; do3t = rdfield ! flds_dissoc case ('SPED ') ; sped = rdfield ! flds_ionatm case ('SHAL ') ; shal = rdfield ! flds_ionatm case ('PARCD ') ; parcd = rdfield ! flds_ionatm case ('QTIN ') ; qtin = rdfield ! flds_ionzrt case ('XNO3 ') ; xno3 = rdfield ! flds_atmos case ('XNAS ') ; xnas = rdfield ! flds_sodium case ('XNAT ') ; xnat = rdfield ! flds_sodium case ('EDY ') ; edyrd = rdfield ! local case ('QNN ') ; qnn = rdfield ! flds_heat case ('QNC ') ; qnc = rdfield ! flds_heat case ('QSRC ') ; qsrc = rdfield ! flds_heat case ('QSRB ') ; qsrb = rdfield ! flds_heat case ('QN1D ') ; qn1d = rdfield ! flds_heat case ('XLEN ') ; xlen = rdfield ! flds_heat (totheat) case ('QO3PR ') ; qo3pr = rdfield ! flds_heat case ('QIAN ') ; qian = rdfield ! flds_ionzrt case ('QJOUL ') ; qjoul = rdfield ! flds_heat case ('QNO3 ') ; qno3 = rdfield ! flds_heat (totheat) case ('QNLYA ') ; qnlya = rdfield ! flds_heat case ('XNOC ') ; xnoc = rdfield ! flds_heat (totheat) case ('XCO2C ') ; xco2c = rdfield ! flds_heat (totheat) case ('XO3PC ') ; xo3pc = rdfield ! flds_heat (totheat) case ('HCM ') ; hcm = rdfield ! flds_heat (totheat) case ('HCE ') ; hce = rdfield ! flds_heat (totheat) case ('DICHM ') ; dichm = rdfield ! flds_dissoc case ('DNCHM ') ; dnchm = rdfield ! flds_dissoc case ('DO2LYA ') ; do2lya = rdfield ! flds_dissoc case ('DO2SRC ') ; do2src = rdfield ! flds_dissoc case ('DO2SRB ') ; do2srb = rdfield ! flds_dissoc case ('QNIGHT ') ; qnight = rdfield ! flds_heat case ('QOP ') ; qop = rdfield ! flds_ionzrt case ('QNOP ') ; qnop = rdfield ! flds_ionzrt case ('QO2P ') ; qo2p = rdfield ! flds_ionzrt case ('QN2P ') ; qn2p = rdfield ! flds_ionzrt case ('QNP ') ; qnp = rdfield ! flds_ionzrt case ('QNSPE ') ; qnspe = rdfield ! flds_ionzrt case ('GWHEAT ') ; gwheat = rdfield ! flds_heat (totheat) case ('QPRO ') ; qpro = rdfield ! flds_ionzrt case ('QPROH ') ; qproh = rdfield ! flds_heat case ('QCR ') ; qcr = rdfield ! flds_ionzrt case ('QIXRAY ') ; qixray = rdfield ! flds_ionzrt case ('XNN2O ') ; xnn2o = rdfield ! flds_atmos case ('XNCH4L ') ; xnch4l = rdfield ! flds_atmos case ('XNNOY ') ; xnnoy = rdfield ! flds_atmos case ('XNHNO3 ') ; xnhno3 = rdfield ! flds_atmos case ('XNN2O5 ') ; xnn2o5 = rdfield ! flds_atmos case ('XNH2O2L ') ; xnh2o2l= rdfield ! flds_atmos case ('XNNO3 ') ; xnno3 = rdfield ! flds_atmos case ('XN4SL ') ; xn4sl = rdfield ! flds_atmos case ('XNO21D ') ; xno21d = rdfield ! flds_atmos case ('XNO21S ') ; xno21s = rdfield ! flds_atmos case default write(6,"('rdvar: Field ',a,' not transferred to array.')") | fname end select ! if (trim(fname)=='XNHTOT') ! | write(6,"('readhist: xnhtot=',/,(6e12.4))") xnhtot ! if (trim(fname)=='XNH2O') ! | write(6,"('readhist: xnh2o=',/,(6e12.4))") xnh2o ! if (trim(fname)=='XNH2O2') ! | write(6,"('readhist: xnh2o2=',/,(6e12.4))") xnh2o2 ! if (trim(fname)=='XNH2') ! | write(6,"('readhist: xnh2=',/,(6e12.4))") xnh2 ! if (trim(fname)=='XNOH') ! | write(6,"('readhist: xnoh=',/,(6e12.4))") xnoh ! if (trim(fname)=='XNHO2') ! | write(6,"('readhist: xnho2=',/,(6e12.4))") xnho2 ! if (trim(fname)=='XNH') ! | write(6,"('readhist: xnh=',/,(6e12.4))") xnh ! if (trim(fname)=='XNCH4') ! | write(6,"('readhist: xnch4=',/,(6e12.4))") xnch4 ! ! zp read from history must be same as zp set by cons (zpcons): ! if (trim(fname)=='ZP') then istat = 0 do k=1,nlev if (zp(k) /= zpcons(k)) then istat = 1 endif enddo ! k=1,nlev if (istat==1) then write(6,"(/,'>>> restart: zp not equal to zpcons:')") write(6,"('zp read from history = ',/,(6e12.4))") zp write(6,"('zpcons set by cons = ',/,(6e12.4))") zpcons write(6,"('This means you are trying to run the model ',/, | 'at a different vertical zp resolution than the source', | ' history.',/)") call shutdown('zp') endif endif end subroutine rdvar !------------------------------------------------------------------- subroutine findhist(ncid,iyear,iday,ihist) implicit none ! ! Args: integer,intent(in) :: ncid integer,intent(inout) :: iyear,iday integer,intent(out) :: ihist ! ! Local: integer :: i,istat,id_year,id_day,id_ut,id_time,ntimes real :: year,day,ut,time ! ihist = 0 ! ! Get number of histories (times) on the file: ! If this nf_inq_unlimdim fails, this is probably an old glbmean ! history, which did not have a time dimension (only one history ! per file). In this case, go ahead and read the single history, ! and assume it is for iyear,iday (start_year, start_day). ! istat = nf_inq_unlimdim(ncid,id_time) ! id of unlimited record var if (istat /= NF_NOERR .or. id_time < 0) then write(6,"(/,'>>> findhist: Error getting id of unlimited ', | 'dimension.')") write(6,"('This is probably an old history, without the time ', | 'dimension.',/)") ihist = -1 return else write(6,"('findhist: id_time=',i3)") id_time endif istat = nf_inq_dimlen(ncid,id_time,ntimes) ! length of time var if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error getting length of time record variable') ! ! Get variable ids for year, day, ut: ! (if any of these fail, then this is probably an old history file ! that did not have time histories. In this case, simply return ! ihist == 1) ! istat = nf_inq_varid(ncid,"year",id_year) if (istat /= NF_NOERR) write(6,"('>>> findhist: error getting', | ' id for var year.')") istat = nf_inq_varid(ncid,"day",id_day) if (istat /= NF_NOERR) write(6,"('>>> findhist: error getting', | ' id for var day.')") istat = nf_inq_varid(ncid,"ut",id_ut) if (istat /= NF_NOERR) write(6,"('>>> findhist: error getting', | ' id for var ut.')") ! write(6,"(/,'Searching for year,day=',2i5,' ntimes=',i4)") | iyear,iday,ntimes do i=1,ntimes istat = nf_get_var1_double(ncid,id_year,i,year) istat = nf_get_var1_double(ncid,id_day ,i,day ) istat = nf_get_var1_double(ncid,id_ut ,i,ut ) istat = nf_get_var1_double(ncid,id_time,i,time) ! write(6,"('findhist: i=',i4,' year,day,ut,time=',4f10.2)") ! | i,year,day,ut,time if (year==float(iyear).and.day==float(iday).and.ut==0.) then ihist = i write(6,"('Found history for iyear,iday =', | 2i5,' (ut==0): ihist=',i4)") iyear,iday,ihist return endif enddo ! i=1,ntimes ! ! If search fails, default to last history on the file: ! write(6,"('>>> findhist: could not find history for year=', | i5,' day=',i5,': searched ',i4,' histories.')") | iyear,iday,ntimes ihist = ntimes istat = nf_get_var1_double(ncid,id_year,ntimes,year) istat = nf_get_var1_double(ncid,id_day ,ntimes,day ) write(6,"('Will read last history: ihist=',i3,' year=',i4, | ' day=',i3)") ntimes,int(year),int(day) ! ! Setting iyear and iday here will override namelist start_year, ! start_day with year,day of history being read. For now, honor ! user's start day and year, even tho that history was not found ! on the source file. ! ! iyear = int(year) ! iday = int(day) write(6,"(' ')") end subroutine findhist !------------------------------------------------------------------- end module restart