! module wacxg5_module ! ! Read t,u,v,z from wacxg5_ncfiles (WACCM-X history files), for use in ! nudging timegcm lower boundary. Data are interpolated to current model ! time. Nudging takes place with 100% at lb, then reducing the percentage ! of waccmx data upward for wacxg5_nlev levels (see dt.F and duv.F). If ! wacxg5_nlev==1, then use only 100% at the lb (no nudging in this case). ! ! If wacxg5_daily==1, then daily-average the waccm data before interpolation. ! If wacxg5_zm==1, then zonally average before nudging. wacxg5_daily and ! wacxg5_zm are mutually exclusive. ! use params_module,only: nlon,nlonp4,glon1,nlat,dlat,dlon,nlevp1, | dlev,zibot use nchist_module,only: handle_ncerr use addfld_module,only: addfld use input_module ,only: | wacxg5_nlev, ! number of levels to nudge from lb | wacxg5_daily, ! 0/1 for whether or not to daily-average | wacxg5_zm, ! 0/1 for whether or not to take zonal means | nfiles =>wacxg5_nfiles, ! number of wacxg5 files | ncfiles=>wacxg5_ncfiles, ! wacxg5 file names | mxlen_filename use init_module,only: iyear,istep use cons_module,only: pi implicit none integer :: ncid=0 ! netcdf file id of currently open wacxg5 data file integer :: ntime=0 ! number of histories on open wacxg5 file ! integer,parameter :: mxtime=48 ! max number of times per file #include ! ! Year, day of year, and hour from wacxg5 file: integer,allocatable,dimension(:),save :: | wacxg5_date, ! YYYYMMDD | wacxg5_datesec, ! secs of current date | wacxg5_date0, ! YYYYMMDD | wacxg5_datesec1 ! secs of current date ! ! T,Z,U,V from wacxg5_ncfile, interpolated to current time. real,dimension(nlonp4,nlat,nlevp1) :: ! global 3d fields | t_wacxg5, u_wacxg5, v_wacxg5, z_wacxg5 real,dimension(nlat,nlevp1) :: ! zonal means | tzm_wacxg5, uzm_wacxg5, vzm_wacxg5, zzm_wacxg5 ! ! Percentage of wacxg5 data to nudge with from lbc to wacxg5_nlev real,allocatable,save :: wacxg5_nudge(:) ! (wacxg5_nlev) ! ! Current wacxg5_ncfile being used for nudging (for histories): character(len=mxlen_filename) :: wacxg5_ncfile_curr ! contains !----------------------------------------------------------------------- subroutine wacxg5_bndry(modeltime) implicit none ! ! Bracket wacxg5 data (on multiple files) and interpolate to modeltime. ! Return interpolated data to t_wacxg5, et.al. arrays module data above. ! ! Arg: integer,intent(in) :: modeltime(4) ! ! Local: integer :: istat,k real :: zp(nlevp1) ! ! Calculate nudging factor once per run (wacxg5_nlev is from namelist): if (istep==1) then allocate(wacxg5_nudge(wacxg5_nlev),stat=istat) do k=1,nlevp1 zp(k) = zibot+(k-1)*dlev enddo do k=1,wacxg5_nlev if (wacxg5_nlev==1) then ! 100% lb only, no nudging wacxg5_nudge(k) = 0. else wacxg5_nudge(k) = | cos(pi/2.*((zp(k)-zp(1))/(zp(wacxg5_nlev)-zp(1))))**2 endif enddo write(6,"('wacxg5_bndry: wacxg5_nlev=',i4,' wacxg5_nudge=', | /,(5f15.5))") wacxg5_nlev,wacxg5_nudge endif ! ! Interpolate to current model time, taking daily mean, zonal mean, ! or simple hourly according to namelist flags. ! (zonal means are calculated in get_wacxg5_hourly if wacxg5_zm > 0) ! if (wacxg5_daily > 0) then call get_wacxg5_daily(modeltime) else call get_wacxg5_hourly(modeltime) endif end subroutine wacxg5_bndry !----------------------------------------------------------------------- subroutine get_wacxg5_hourly(modeltime) ! ! Arg: integer,intent(in) :: modeltime(4) ! ! Local: integer,save :: ifile=1,itime0=1,itime1=2 real,save :: yfrac, ! year fraction at modeltime | yfrac0,yfrac1 ! year fraction at data itime0,1 logical,save :: isread(mxtime)=.false. integer :: istat real,dimension(nlon,nlat,nlevp1) :: ! full column | t0, u0, v0, z0, ! data at time index itime0 | t1, u1, v1, z1 ! data at time index itime1 100 continue if (ncid==0) call readfile(ncfiles(ifile),ifile,ncid) yfrac = yfrac_tgcm(iyear,modeltime) if (itime0 < itime1) yfrac0 = yfrac_wacx(itime0) yfrac1 = yfrac_wacx(itime1) if (yfrac >= yfrac0 .and. yfrac < yfrac1) then if (.not.isread(itime0)) then call readdata(ncid,ifile,itime0,t0,u0,v0,z0,0) isread(itime0) = .true. endif if (.not.isread(itime1)) then call readdata(ncid,ifile,itime1,t1,u1,v1,z1,1) isread(itime1) = .true. endif ! ! Interpdata saves interpolated data in t,u,v,z_wacxg5 (module data above) call interpdata(modeltime,yfrac,yfrac0,yfrac1,nlevp1, | t0,u0,v0,z0,t1,u1,v1,z1) wacxg5_ncfile_curr = ncfiles(ifile) ! for history file ! ! Calculate zonal means of interpolated data if necessary: if (wacxg5_zm > 0) call calc_zm else t0 = t1 ; u0 = u1 ; v0 = v1 ; z0 = z1 if (itime0 > itime1) then itime0 = 1 itime1 = 2 isread(:) = .false. else itime0 = itime1 itime1 = itime1+1 endif ! ! End of current file. Retain t0 data at yfrac0 from this last time, ! reset itime1=1, and go back up to open new file. Note that now ! itime0 > itime1 when interpolating between last time on previous ! file and first time on new file, see conditionals above. ! if (itime1 > ntime) then ifile = ifile+1 if (ifile > nfiles) stop 'out of wacxg5 files' istat = nf_close(ncid) ncid = 0 itime1 = 1 ! now itime0==ntime, so itime0 > itime1 isread(itime1) = .false. endif goto 100 endif end subroutine get_wacxg5_hourly !----------------------------------------------------------------------- subroutine calc_zm use mpi_module,only: mytid ! ! Calculate zonal means of t,u,v,z_wacxg5, returning tzm,uzm,vzm,zzm_wacxg5 ! These arrays are in module data above. ! ! Local: integer :: i,j,k real :: rlon real :: tzm_diag(nlevp1,nlonp4,nlat) real :: uzm_diag(nlevp1,nlonp4,nlat) real :: vzm_diag(nlevp1,nlonp4,nlat) rlon = 1./real(nlon) ! note 3->nlon+2 = nlon tzm_wacxg5=0. ; uzm_wacxg5=0. ; vzm_wacxg5=0. ; zzm_wacxg5=0. do k=1,nlevp1 ! ! Output zonal means of interpolated data: do j=1,nlat do i=3,nlon+2 ! data only (not periodic points) tzm_wacxg5(j,k) = tzm_wacxg5(j,k)+t_wacxg5(i,j,k) uzm_wacxg5(j,k) = uzm_wacxg5(j,k)+u_wacxg5(i,j,k) vzm_wacxg5(j,k) = vzm_wacxg5(j,k)+v_wacxg5(i,j,k) zzm_wacxg5(j,k) = zzm_wacxg5(j,k)+z_wacxg5(i,j,k) enddo ! i=3,nlon+2 tzm_wacxg5(j,k) = tzm_wacxg5(j,k)*rlon uzm_wacxg5(j,k) = uzm_wacxg5(j,k)*rlon vzm_wacxg5(j,k) = vzm_wacxg5(j,k)*rlon zzm_wacxg5(j,k) = zzm_wacxg5(j,k)*rlon enddo ! j=1,nlat enddo ! k=1,nlevp1 ! ! Save zm fields to secondary history: ! do j=1,nlat ! do i=3,nlon+2 ! tzm_diag(:,i,j) = tzm_wacxg5(j,:) ! redundant in lon ! uzm_diag(:,i,j) = uzm_wacxg5(j,:) ! redundant in lon ! vzm_diag(:,i,j) = vzm_wacxg5(j,:) ! redundant in lon ! enddo ! tzm_diag(:,1:2,j) = tzm_diag(:,nlon+1:nlon+2,j) ! tzm_diag(:,nlon+3:nlon+4,j) = tzm_diag(:,3:4,j) ! uzm_diag(:,1:2,j) = uzm_diag(:,nlon+1:nlon+2,j) ! uzm_diag(:,nlon+3:nlon+4,j) = uzm_diag(:,3:4,j) ! vzm_diag(:,1:2,j) = vzm_diag(:,nlon+1:nlon+2,j) ! vzm_diag(:,nlon+3:nlon+4,j) = vzm_diag(:,3:4,j) ! call addfld('tzm_diag',' ',' ',tzm_diag(:,:,j), ! | 'lev',1,nlevp1,'lon',1,nlonp4,j) ! call addfld('uzm_diag',' ',' ',uzm_diag(:,:,j), ! | 'lev',1,nlevp1,'lon',1,nlonp4,j) ! call addfld('vzm_diag',' ',' ',vzm_diag(:,:,j), ! | 'lev',1,nlevp1,'lon',1,nlonp4,j) ! enddo end subroutine calc_zm !----------------------------------------------------------------------- subroutine get_wacxg5_daily(modeltime) ! ! Read lower boundary wacxg5 data (t,u,v,z) from a series of hourly 1-day ! WACCM-X history files. Daily-average data across the 24 hours, and ! interpolate to current model time (day,hour,min). The averaged and ! interpolated data (t_wacxg5, et.al. in this module) is added to the model ! lower boundary in module lbc (lbc.F). ! ! Arg: integer,intent(in) :: modeltime(4) ! ! Local: integer,save :: | ifile,ifile0,ifile1, ! indices to wacxg5 files ncfiles(nfiles) | ncid0,ncid1 ! left,right netcdf file ids real,save :: | hour_prev ! decimal model hour at previous step real :: hour_curr, ! decimal hours of current model time | hour_int ! hour to interpolate logical,save :: first=.true. integer :: i,istat real,dimension(nlon,nlat,1) :: ! note lb only | t0, u0, v0, z0, ! data at day0 ! data at left interp boundary | t1, u1, v1, z1 ! data at day1 ! data at right interp boundary ! ! At the first model step, find and open wacxg5 file (subs findfile and ! readfile) for model day of first modeltime. If first model decimal ! hour is <= 12, then open the previous day file, otherwise open the ! next day file. Read and average data in both days (sub readdata_daily), ! and linearly interpolate to the current model hour (sub interpdata). ! Averaged data are assumed to be at ut 12, for purpose of interpolation. ! if (first) then call findfile(iyear,modeltime,ifile,ncid) ! returns ifile,ncid hour_curr = real(modeltime(2))+real(modeltime(3))/60. if (hour_curr <= 12.) then ! will need previous day ifile1 = ifile ncid1 = ncid ifile0 = ifile-1 if (ifile0 < 1) call shutdown('need previous waccx file') call readfile(ncfiles(ifile0),ifile0,ncid0) ! read previous file elseif (hour_curr > 12.) then ! will need next file ifile0 = ifile ncid0 = ncid ifile1 = ifile+1 if (ifile1 > nfiles) call shutdown('out of waccmx files') call readfile(ncfiles(ifile1),ifile1,ncid1) ! read next file endif ! ! Read and average data for both days: call readdata_daily(ncid0,t0,u0,v0,z0,0) call readdata_daily(ncid1,t1,u1,v1,z1,1) hour_prev = hour_curr endif ! first call ! ! Current model decimal hour: ! hour_curr = real(modeltime(2))+real(modeltime(3))/60. ! ! Advance and read the next wacxg5 file if model time has crossed the ! 12 ut boundary. Transfer indices and data from t1 to t0, and read ! new t1 data. ! if (hour_prev <= 12..and.hour_curr > 12.) then istat = nf_close(ncid0) ifile0 = ifile1 ncid0 = ncid1 t0 = t1 ; u0 = u1 ; v0 = v1 ; z0 = z1 ifile1 = ifile0+1 if (ifile1 > nfiles) call shutdown('out of waccmx files') call readfile(ncfiles(ifile1),ifile1,ncid1) ! open next file with ncid1 ! ! Read and average data of next day for new right-hand interpolation data: ! call readdata_daily(ncid1,t1,u1,v1,z1,1) ! write(6,"('wacxg5_daily: Advance to ifile1 ',i3, ! | ' day=',i3,' ut=',f9.4,' wacxg5 file=',/,a)") ! | ifile1,modeltime(1),hour_curr,trim(ncfiles(ifile1)) endif ! crossed ut12 boundary ! ! Adjust hour to interpolate to (hour_int), so interpolation ! boundaries can be 0->24. ! if (hour_curr <= 12) then hour_int = hour_curr+12. if (hour_int > 24.or.hour_int < 12.) | call shutdown('wacxg5_daily: bad hour_int > 24 or < 12') else hour_int = hour_curr-12. if (hour_int < 0.or.hour_int > 12.) | call shutdown('wacxg5_daily: bad hour_int < 0 or > 12') endif ! ! Interpolate to hour_int between 0 and 24. Interpdata defines ! wacxg5 module data t_wacxg5, etc. This is added to the lower ! boundary by the lbc module (lbc.F). ! call interpdata(modeltime,hour_int,0.,24.,1, | t0,u0,v0,z0,t1,u1,v1,z1) ! write(6,"('wacxg5_daily: iyear=',i4,' day=',i4,' hour_curr=', ! | f10.4,' hour_int=',f10.4,' ifile0,1=',2i4)") ! | iyear,modeltime(1),hour_curr,hour_int,ifile0,ifile1 ! write(6,"('wacxg5_daily: t0=',2e14.6,' u0=',2e14.6,' v0=',2e14.6, ! | ' z0=',2e14.6)") minval(t0),maxval(t0),minval(u0),maxval(u0), ! | minval(v0),maxval(v0),minval(z0),maxval(z0) ! write(6,"('wacxg5_daily: t1=',2e14.6,' u1=',2e14.6,' v1=',2e14.6, ! | ' z1=',2e14.6)") minval(t1),maxval(t1),minval(u1),maxval(u1), ! | minval(v1),maxval(v1),minval(z1),maxval(z1) first = .false. hour_prev = hour_curr end subroutine get_wacxg5_daily !----------------------------------------------------------------------- subroutine findfile(iyear,modeltime,ifile,ncid) ! ! Find current model year and day at hour 0 in wacxg5 data files. ! Return index ifile and logical unit ncid. ! ! Args: integer,intent(in) :: iyear,modeltime(4) integer,intent(out) :: ifile,ncid ! ! Local: logical :: found integer :: i,istat,wyear,wday,whour real :: yfrac ifile = 1 ncid = 0 100 continue if (ncid==0) call readfile(ncfiles(ifile),ifile,ncid) found = .false. searchloop: do i=1,ntime yfrac = yfrac_wacx(i,wyear,wday,whour) if (wyear==iyear.and.wday==modeltime(1).and.whour==0) then found = .true. exit searchloop endif enddo searchloop if (found) then return else ! advance to next file ifile = ifile+1 if (ifile > nfiles) then write(6,"('>>> findfile: could not find file with iyear=', | i4,' day=',i4)") iyear,modeltime(1) call shutdown('out of wacxg5 files') endif istat = nf_close(ncid) ncid = 0 goto 100 endif end subroutine findfile !----------------------------------------------------------------------- subroutine readfile(ncfile,ifile,ncid) ! ! Open wacxg5 file (WACCM-X history file written by Valery Yudin). ! Verify dimensions and coordinates, and read date and times (lb data ! is not read at this time). ! implicit none ! ! Args: character(len=*),intent(in) :: ncfile integer,intent(in) :: ifile integer,intent(out) :: ncid ! ! Local: integer :: istat,ndims,idunlim,ngatts,nvars integer :: nlat_rd,nlon_rd,nlev_rd integer :: id_lat,id_lon,id_time,id_lev integer :: idv_lev,idv_date,idv_datesec integer :: istart1(1),icount1(1) character(len=1024) :: dskfile real,allocatable :: wacxg5_levs(:) ! dskfile = ' ' call getfile(ncfile,dskfile) istat = nf_open(ncfile,NF_NOWRITE,ncid) istat = nf_inq(ncid,ndims,nvars,ngatts,idunlim) write(6,"('Opened wacxg5 file ',a)") trim(ncfile) ! ! Get and verify latitude dimension. ! (read nlat_rd, compare with nlat, which is in params.h): istat = nf_inq_dimid(ncid,'lat',id_lat) istat = nf_inq_dimlen(ncid,id_lat,nlat_rd) if (nlat_rd /= nlat) then write(6,"(/,'>>> readfile: bad nlat_rd=',i3, | ' -- should be nlat=',i3)") nlat_rd,nlat call shutdown('readfile') endif ! ! Get and verify longitude dimension. ! (read nlon_rd, compare with nlon, which is in params.h): istat = nf_inq_dimid(ncid,'lon',id_lon) istat = nf_inq_dimlen(ncid,id_lon,nlon_rd) if (nlon_rd /= nlon) then write(6,"(/,'>>> readfile: bad nlon_rd=',i3, | ' -- should be nlon=',i3)") nlon_rd,nlon call shutdown('readfile') endif ! ! Get and verify levels dimension. ! (read nlev_rd, compare with nlevp1, which is in params.h): istat = nf_inq_dimid(ncid,'lev',id_lev) istat = nf_inq_dimlen(ncid,id_lev,nlev_rd) if (nlev_rd /= nlevp1) then write(6,"(/,'>>> readfile: bad nlev_rd=',i3, | ' -- should be nlevp1=',i3)") nlev_rd,nlevp1 call shutdown('readfile') endif ! ! ilev dimension should be same as lev: istat = nf_inq_dimid(ncid,'ilev',id_lev) istat = nf_inq_dimlen(ncid,id_lev,nlev_rd) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'readfile: error getting levels dimension') if (nlev_rd /= nlevp1) then write(6,"(/,'>>> readfile: bad nlev_rd (ilev)=',i3, | ' -- should be nlevp1=',i3)") nlev_rd,nlevp1 call shutdown('readfile') endif ! ! Get number of times (value of unlimited dimension) on the file: istat = nf_inq_unlimdim(ncid,id_time) ! id of unlimited record var istat = nf_inq_dimlen(ncid,id_time,ntime) ! length of time var if (ntime > mxtime) then write(6,"('>>> Reading wacxg5 data file: ntime=',i4,' > ', | 'mxtime=',i4)") ntime,mxtime call shutdown('wacxg5 ntime') endif ! ! Allocate date and datesec vars: if (allocated(wacxg5_date)) deallocate(wacxg5_date) if (allocated(wacxg5_datesec)) deallocate(wacxg5_datesec) allocate(wacxg5_date(ntime),stat=istat) allocate(wacxg5_datesec(ntime),stat=istat) istart1(1) = 1 icount1(1) = ntime ! ! date (YYYYMMDD): istat = nf_inq_varid(ncid,"date",idv_date) istat = nf_get_vara_int(ncid,idv_date,istart1,icount1, | wacxg5_date) ! ! datesec (seconds in current day): istat = nf_inq_varid(ncid,"datesec",idv_datesec) istat = nf_get_vara_int(ncid,idv_datesec,istart1,icount1, | wacxg5_datesec) ! write(6,"('readfile: file=',a,' date=',/,(6i10))") ! | trim(ncfile),wacxg5_date ! write(6,"('readfile: datesec=',/,(6i10))") wacxg5_datesec end subroutine readfile !----------------------------------------------------------------------- real function yfrac_tgcm(iyear,modeltime) ! ! Convert timegcm modeltime (day,hour,min) to fraction of a year. ! ! Args: integer,intent(in) :: iyear,modeltime(4) ! ! Local: integer :: iday,isecs iday = modeltime(1) isecs = modeltime(2)*3600+modeltime(3)*60+modeltime(4) yfrac_tgcm=float(iyear)+ | (float(iday)-1.+float(isecs)/(3600.*24.))/365. end function yfrac_tgcm !----------------------------------------------------------------------- real function yfrac_wacx(it,wacx_year,wacx_day,wacx_hour) ! ! Convert waccm date and time to fraction of a year. ! ! Args: integer,intent(in) :: it integer,optional :: wacx_year,wacx_day,wacx_hour ! ! Local: integer :: wyear,wmonth,wday,wsec,m,wjday integer :: ndmon(12) = ! J F M A M J J A S O N D | (/31,28,31,30,31,30,31,31,30,31,30,31/) wyear = wacxg5_date(it)/10000 ! waccm year wmonth= (wacxg5_date(it)-wyear*10000)/100 ! waccm month of year wday = wacxg5_date(it)-(wyear*10000+wmonth*100) ! waccm day of month wsec = wacxg5_datesec(it) ! waccm secs into day wjday = 0 ! waccm julian day do m=1,wmonth if (m > 1) wjday = wjday + ndmon(m-1) enddo wjday = wjday+wday yfrac_wacx=float(wyear)+ ! waccm year+fraction | (float(wjday)-1.+float(wsec)/(3600.*24.))/365. if (present(wacx_year)) wacx_year = wyear if (present(wacx_day)) wacx_day = wjday if (present(wacx_hour)) wacx_hour = wsec/3600 end function yfrac_wacx !----------------------------------------------------------------------- subroutine readdata(ncid,ifile,itime,t,u,v,z,index) implicit none ! ! Read 3d data t,u,v,z from wacxg5 file. ! The file is already open under ncid. ! ! Args: integer,intent(in) :: ncid,ifile,itime,index real,dimension(nlon,nlat,nlevp1),intent(out) :: t,u,v,z ! ! Local: integer :: istat,id,start(4),count(4) t=0. ; u=0. ; v=0. ; z=0. start(1:3) = 1 start(4) = itime count(1) = nlon count(2) = nlat count(3) = nlevp1 ! full column for diag count(4) = 1 istat = nf_inq_varid(ncid,"TA",id) istat = nf_get_vara_double(ncid,id,start,count,t) istat = nf_inq_varid(ncid,"UA",id) istat = nf_get_vara_double(ncid,id,start,count,u) istat = nf_inq_varid(ncid,"VA",id) istat = nf_get_vara_double(ncid,id,start,count,v) istat = nf_inq_varid(ncid,"Z3A",id) istat = nf_get_vara_double(ncid,id,start,count,z) ! write(6,"('Read data: ifile ',i3,': itime',i1,'=',i3, ! | ' min,max T=',2e12.4,' U=',2e12.4,' V=',2e12.4,' Z=',2e12.4)") ! | ifile,index,itime,minval(t),maxval(t),minval(u),maxval(u), ! | minval(v),maxval(v),minval(z),maxval(z) end subroutine readdata !----------------------------------------------------------------------- subroutine readdata_daily(ncid,t,u,v,z,index) ! ! Read 2d (lower boundary only) data t,u,v,z from wacxg5 file. ! Read 24 hours (ntime) from file attached to ncid, and take daily ! average. The file is already open under ncid. ! ! Args: integer,intent(in) :: ncid,index real,dimension(nlon,nlat,1),intent(out) :: t,u,v,z ! lb only ! ! Local: integer :: i,j,istat,id,start(4),count(4) integer :: idt,idu,idv,idz real,dimension(nlon,nlat) :: wt,wu,wv,wz ! start(1:3) = 1 count(1) = nlon count(2) = nlat count(3) = 1 ! lb only count(4) = 1 istat = nf_inq_varid(ncid,"TA",idt) istat = nf_inq_varid(ncid,"UA",idu) istat = nf_inq_varid(ncid,"VA",idv) istat = nf_inq_varid(ncid,"Z3A",idz) t=0. ; u=0. ; v=0. ; z=0. do i=1,ntime start(4) = i istat = nf_get_vara_double(ncid,idt,start,count,wt) istat = nf_get_vara_double(ncid,idu,start,count,wu) istat = nf_get_vara_double(ncid,idv,start,count,wv) istat = nf_get_vara_double(ncid,idz,start,count,wz) t(:,:,1) = t(:,:,1)+wt(:,:) u(:,:,1) = u(:,:,1)+wu(:,:) v(:,:,1) = v(:,:,1)+wv(:,:) z(:,:,1) = z(:,:,1)+wz(:,:) enddo t(:,:,1) = t(:,:,1)/real(ntime) u(:,:,1) = u(:,:,1)/real(ntime) v(:,:,1) = v(:,:,1)/real(ntime) z(:,:,1) = z(:,:,1)/real(ntime) end subroutine readdata_daily !----------------------------------------------------------------------- subroutine interpdata(modeltime,yfrac,yfrac0,yfrac1,nlev, | t0,u0,v0,z0,t1,u1,v1,z1) ! ! Interpolate to yfrac between t,u,v,z0 at yfrac0 and t,u,v,z1 at yfrac1. ! Return module data t,u,v,z_wacxg5 w/ unit conversion where necessary. ! Functioni time_interp is in util.F. Number of levels is passed in. ! use mpi_module,only: lat0,lat1,lon0,lon1 implicit none ! ! Args: integer,intent(in) :: modeltime(4),nlev real,intent(in) :: yfrac,yfrac0,yfrac1 real,dimension(nlon,nlat,nlev),intent(in) :: | t0,u0,v0,z0,t1,u1,v1,z1 ! ! Local: integer :: i,j,k real,dimension(nlev,lon0:lon1) :: twacx,uwacx,vwacx,zwacx ! ! External: real,external :: time_interp ! util.F do k=1,nlev do j=1,nlat do i=3,nlon+2 t_wacxg5(i,j,k)= time_interp(t0(i-2,j,k),t1(i-2,j,k),yfrac0, | yfrac1,yfrac) u_wacxg5(i,j,k)= time_interp(u0(i-2,j,k),u1(i-2,j,k),yfrac0, | yfrac1,yfrac) v_wacxg5(i,j,k)= time_interp(v0(i-2,j,k),v1(i-2,j,k),yfrac0, | yfrac1,yfrac) z_wacxg5(i,j,k)= time_interp(z0(i-2,j,k),z1(i-2,j,k),yfrac0, | yfrac1,yfrac) enddo ! i=1,nlon ! ! Periodic points: t_wacxg5(1:2,j,k) = t_wacxg5(nlon+1:nlon+2,j,k) t_wacxg5(nlon+3:nlon+4,j,k) = t_wacxg5(3:4,j,k) u_wacxg5(1:2,j,k) = u_wacxg5(nlon+1:nlon+2,j,k) u_wacxg5(nlon+3:nlon+4,j,k) = u_wacxg5(3:4,j,k) v_wacxg5(1:2,j,k) = v_wacxg5(nlon+1:nlon+2,j,k) v_wacxg5(nlon+3:nlon+4,j,k) = v_wacxg5(3:4,j,k) z_wacxg5(1:2,j,k) = z_wacxg5(nlon+1:nlon+2,j,k) z_wacxg5(nlon+3:nlon+4,j,k) = z_wacxg5(3:4,j,k) enddo ! j=1,nlat enddo ! 1,nlev ! ! Convert from m/s to cm/s and m to cm: u_wacxg5 = u_wacxg5*100. v_wacxg5 = v_wacxg5*100. z_wacxg5 = z_wacxg5*100. ! do k=1,nlev ! write(6,"('Interpolated global wacxg5 at modeltime=',4i4, ! | ': k=',i3,' yfrac=',f16.10,' T min,max=',2e12.4,' U min,max=', ! | 2e12.4,' V min,max=',2e12.4,' Z min,max=',2e12.4)") modeltime, ! | k,yfrac, ! | minval(t_wacxg5(:,:,k)),maxval(t_wacxg5(:,:,k)), ! | minval(u_wacxg5(:,:,k)),maxval(u_wacxg5(:,:,k)), ! | minval(v_wacxg5(:,:,k)),maxval(v_wacxg5(:,:,k)), ! | minval(z_wacxg5(:,:,k)),maxval(z_wacxg5(:,:,k)) ! enddo ! ! Save 3d interpolated fields to secondary history: ! twacx=0. ; uwacx=0. ; vwacx=0. ; zwacx=0. ! do j=lat0,lat1 ! do i=lon0,lon1 ! do k=1,nlev ! twacx(k,i) = t_wacxg5(i,j,k) ! uwacx(k,i) = u_wacxg5(i,j,k) ! vwacx(k,i) = v_wacxg5(i,j,k) ! zwacx(k,i) = z_wacxg5(i,j,k) ! enddo ! k=1,nlev ! enddo ! i=lon0,lon1 ! call addfld('t_wacxg5',' ','deg K',twacx, ! | 'lev',1,nlev,'lon',lon0,lon1,j) ! call addfld('u_wacxg5',' ','cm/s',uwacx, ! | 'lev',1,nlev,'lon',lon0,lon1,j) ! call addfld('v_wacxg5',' ','cm/s',vwacx, ! | 'lev',1,nlev,'lon',lon0,lon1,j) ! call addfld('z_wacxg5',' ','cm',zwacx, ! | 'lev',1,nlev,'lon',lon0,lon1,j) ! enddo ! j=lat0,lat1 end subroutine interpdata !----------------------------------------------------------------------- end module wacxg5_module