      subroutine infld(j)
c
c=======================================================================
c                                                                    ===
c  This subroutine reads in initial  and boundary conditions from    ===
c  initialization NetCDF file.  INFLD is called at initialization    ===
c  and re-start procedures by OCEAN and once per timestep by STEP    ===
c  to set boundary conditions. Boundary conditions are persistent    ===
c  or  interpolated linearly between boundary updates.  INFLD  is    ===
c  also called during an assimilation cycle, if appropriate.         ===
c                                                                    ===
#ifndef analytical
c  Calls:  NCDID, NCDINQ, NCOPN, NCPOPT, NCVGT, NCVGTC, NCVGT1,      ===
c          NCVID  (NetCDF library)                                   ===
c                                                                    ===
#else
c  Calls:  ANABNDRY, ANACOAST, ANADINFO, ANAMRHO, ANATOPO, ANAUVPT   ===
c          DPTHBNDRY  (analytical set-up)                            ===
#endif
c          EXITUS, LNBLK                                             ===
#ifdef coast
c          LANDSEA  (coastal set-up)                                 ===
c          MSK_BDY (coastal set-up)                                  ===
#endif
#ifdef oias
c          OIASSI   (OI assimilation set-up)                         ===
#endif
#ifndef analytical
c                                                                    ===
c  WARNING:   Character argument to NetCDF routines NCDID and NCVID  ===
c             is (upper/lower) case sensitive.                       ===
#endif
c                                                                    ===
c=======================================================================
c
c-----------------------------------------------------------------------
c  Define global data.
c-----------------------------------------------------------------------
c
#include <cdefs.h>
#include <param.h>
#include <pconst.h>
#include <fullwd.h>
#include <scalar.h>
#include <fields.h>
#ifndef analytical
# include <workspa.h>
#endif
#include <bndata.h>
#include <moddat.h>
#include <options.h>
#include <iounits.h>
#ifdef oias
# include <oiopts.h>
#endif
#include <hybrid.h>
#include <runid.h>
#include <rhomean.h>
#ifndef analytical
# include <netcdf.inc>
# include <pe_netcdf.h>
#endif
#ifndef gridold
# include <vertslabs.h>
#endif
c
c-----------------------------------------------------------------------
c  Define local data.
c-----------------------------------------------------------------------
c
#ifndef analytical
      integer ijmxkm
      parameter (ijmxkm=ijmx*km)
      integer i,ib,imax,j,jmax,jj,k,kmax,lenstr,m,ntrc,t0inpindx
      integer count(5),start(5)
      integer lnblk
      logical old
      FLOAT
     *      dayper,dt1,dt2,dt3
      FLOAT
     *      fld(ijmxkm)
      character*44 fldnam
      character*80 varnam
c
      parameter (dayper=c1e6,t0inpindx=1)
c
      data old /.false./
#else
      integer j
#endif
#ifndef analytical
c
c-----------------------------------------------------------------------
c  On first call, open initialization NetCDF file and read in model
c  configuration parameters and fields.
c-----------------------------------------------------------------------
c
      if(itt.lt.0) then
c
c  Open input NetCDF file.
c
        call ncpopt(ncverbos)
        lenstr=lnblk(inpname,len(inpname))
        ncinpid=ncopn(inpname(1:lenstr),ncnowrit,rcode)
        if(rcode.ne.0) then
          write(stdout,900) inpname(1:lenstr)
          call exitus('INFLD')
        endif
c
c  Read in titles.
c
        varid=ncdid(ncinpid,'lentit',rcode)
        if(rcode.eq.0) then
          call ncdinq(ncinpid,varid,varnam,lenstr,rcode)
          varid=ncvid(ncinpid,'titleic',rcode)
          if(rcode.eq.0) then
            start(1)=1
            count(1)=lenstr
            call ncvgtc(ncinpid,varid,start,count,titleic,lenstr,rcode)
          else
            write(stdout,910) 'titleic'
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'lentit'
          call exitus('INFLD')
        endif
        varid=ncdid(ncinpid,'lenrun',rcode)
        if(rcode.eq.0) then
          call ncdinq(ncinpid,varid,varnam,lenstr,rcode)
          varid=ncvid(ncinpid,'runic',rcode)
          if(rcode.eq.0) then
            start(1)=1
            count(1)=lenstr
            call ncvgtc(ncinpid,varid,start,count,runic,lenstr,rcode)
          else
            write(stdout,910) 'runic'
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'lenrun'
          call exitus('INFLD')
        endif
c
c  Read in and check main dimension parameters.
c
        varid=ncvid(ncinpid,'imt',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,imax,rcode)
          if(imax.ne.imt) then
            write(stdout,930) 'IMT, IMAX: ',imt,imax
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'imt'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'jmt',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,jmax,rcode)
          if(jmax.ne.jmt) then
            write(stdout,930) 'JMT, JMAX: ',jmt,jmax
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'jmt'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'km',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,kmax,rcode)
          if(kmax.ne.km) then
            write(stdout,930) 'KM, KMAX: ',km,kmax
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'km'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'nt',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,ntrc,rcode)
          if(ntrc.ne.nt) then
            write(stdout,930) 'NT, NTRC: ',nt,ntrc
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'nt'
          call exitus('INFLD')
        endif
        varid=ncdid(ncinpid,'tslev',rcode)
        if(rcode.eq.0) then
          call ncdinq(ncinpid,varid,varnam,nprof,rcode)
          if(nprof.gt.mprof) then
            write(stdout,930) 'MPROF, NPROF: ',mprof,nprof
            call exitus('INFLD')
          endif
        else
          write(stdout,940) 'tslev'
          call exitus('INFLD')
        endif
c
c  Read in flags and switches, IFLAG.
c
        start(1)=1
        count(1)=10
        varid=ncvid(ncinpid,'iflag',rcode)
        if(rcode.eq.0) then
          call ncvgt(ncinpid,varid,start,count,iflag,rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'iflag'
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'iflag'
          call exitus('INFLD')
        endif
c
c  Read in model domain parameters.
c
        varid=ncvid(ncinpid,'coord',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,coord,rcode)
          if((rcode.eq.0).and.((coord.lt.0).or.(coord.gt.2))) then
            write(stdout,1000) 'coord',coord
            call exitus('INFLD')
          endif
         else
          coord=0
        endif
        varid=ncvid(ncinpid,'rlngd',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,rlngd,rcode)
        else
          write(stdout,920) 'rlngd'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'rlatd',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,rlatd,rcode)
        else
          write(stdout,920) 'rlatd'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'delx',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,delx,rcode)
        else
          old = .true.
        endif
        varid=ncvid(ncinpid,'dely',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,dely,rcode)
        else
          old = .true.
        endif
        varid=ncvid(ncinpid,'thetad',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,thetad,rcode)
        else
          write(stdout,920) 'thetad'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'gridx',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,gridx,rcode)
        else
          write(stdout,920) 'gridx'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'gridy',rcode)
        if(rcode.eq.0) then
          call ncvgt1(ncinpid,varid,1,gridy,rcode)
        else
          write(stdout,920) 'gridy'
          call exitus('INFLD')
        endif
c
c  Ensure backward compatibility
c
        if (old) then
          if (coord.ne.1) then
            delx = c0
            dely = c0
            write (stdout,1020)
           else
            delx  = rlngd
            dely  = rlatd
            rlngd = c0
            rlatd = c0
            write (stdout,1030) rlngd,rlatd
          endif
        endif
c
c  Read in thicknesses for the flat reference levels.
c
        start(1)=1
        count(1)=km
        varid=ncvid(ncinpid,'hz',rcode)
        if(rcode.eq.0) then
          call ncvgt(ncinpid,varid,start,count,hz,rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'hz'
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'hz'
          call exitus('INFLD')
        endif
#ifdef dblsigma
c
c  Read in vertical transfomation (double) parameters.
c
        if(iflag(9).eq.1) then
          varid=ncvid(ncinpid,'zc1',rcode)
          if(rcode.eq.0) then
            call ncvgt1(ncinpid,varid,1,zc1,rcode)
          else
            write(stdout,920) 'zc1'
            call exitus('INFLD')
          endif
          varid=ncvid(ncinpid,'zc2',rcode)
          if(rcode.eq.0) then
            call ncvgt1(ncinpid,varid,1,zc2,rcode)
          else
            write(stdout,920) 'zc2'
            call exitus('INFLD')
          endif
          varid=ncvid(ncinpid,'zref',rcode)
          if(rcode.eq.0) then
            call ncvgt1(ncinpid,varid,1,zref,rcode)
          else
            write(stdout,920) 'zref'
            call exitus('INFLD')
          endif
          varid=ncvid(ncinpid,'zslope',rcode)
          if(rcode.eq.0) then
            call ncvgt1(ncinpid,varid,1,zslope,rcode)
          else
            write(stdout,920) 'zslope'
            call exitus('INFLD')
          endif
        endif
        write(stdout,1010) zc1,zc2,zslope,zref
#endif
c
c  Read in mean T-S profile.
c
        start(1)=1
        count(1)=nprof
        varid=ncvid(ncinpid,'tclima',rcode)
        if(rcode.eq.0) then
          call ncvgt(ncinpid,varid,start,count,tinit(1,1),rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'tclima'
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'tclima'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'sclima',rcode)
        if(rcode.eq.0) then
          call ncvgt(ncinpid,varid,start,count,tinit(1,2),rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'sclima'
            call exitus('INFLD')
          endif
        else
          write(stdout,920) 'sclima'
          call exitus('INFLD')
        endif
        varid=ncvid(ncinpid,'zclima',rcode)
        if(rcode.eq.0) then
          call ncvgt(ncinpid,varid,start,count,tinit(1,3),rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'zclima'
            call exitus('INFLD')
          endif
          do 10 k=1,nprof
            tinit(k,3)=abs(tinit(k,3))*m2cm
  10      continue
        else
          write(stdout,920) 'zclima'
          call exitus('INFLD')
        endif
c
c  Inquire time dimension and get input time ID and initialize
c  input time counter.
c
        varid=ncdid(ncinpid,'time0',rcode)
        if(rcode.eq.0) then
          call ncdinq(ncinpid,varid,varnam,tinpmax,rcode)
        else
          write(stdout,940) 'time0'
          call exitus('INFLD')
        endif
        varid=ncdid(ncinpid,'time',rcode)
        if(rcode.eq.0) then
          call ncdinq(ncinpid,varid,varnam,tinpmax,rcode)
        else
          write(stdout,940) 'time'
          call exitus('INFLD')
        endif
c
        t0inpid=ncvid(ncinpid,'time0',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'time0'
          call exitus('INFLD')
        endif
        tinpid=ncvid(ncinpid,'time',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'time'
          call exitus('INFLD')
        endif
        tinpindx=1
c
c  Get initialization time (usually, modified Julian day) for extraneous uses.
c
          start(1)=t0inpindx
          count(1)=1
          call ncvgt(ncinpid,t0inpid,start,count,dstart,rcode)
          if(rcode.ne.0) then
            write(stdout,990) 'time0',t0inpindx
            call exitus('INFLD')
          endif
c
c  Get variables ID from NetCDF initialization file.
c
        piniid=ncvid(ncinpid,'pbar',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'pbar'
          call exitus('INFLD')
        endif
        pbryid=ncvid(ncinpid,'pbry',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'pbry'
          call exitus('INFLD')
        endif
        qbryid=ncvid(ncinpid,'qbry',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'qbry'
          call exitus('INFLD')
        endif
        do 20 m=1,nt
          fldnam=tname(1,m)
          lenstr=lnblk(fldnam,len(fldnam))
          tiniid(m)=ncvid(ncinpid,fldnam(1:lenstr),rcode)
          if(rcode.ne.0) then
            write(stdout,920) fldnam(1:lenstr)
            call exitus('INFLD')
          endif
  20    continue
        tbryid=ncvid(ncinpid,'tbry',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'tbry'
          call exitus('INFLD')
        endif
        viniid=ncvid(ncinpid,'vclin',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'vclin'
          call exitus('INFLD')
        endif
        vbryid=ncvid(ncinpid,'vbry',rcode)
        if(rcode.ne.0) then
          write(stdout,920) 'vbry'
          call exitus('INFLD')
        endif
# ifdef coast
c
c-----------------------------------------------------------------------
c  Read in coastal segments and land/sea mask at the tracer points.
c-----------------------------------------------------------------------
c
        call landsea(ncinpid)
# endif
#else
c
c-----------------------------------------------------------------------
c  Set in model domain configuration.
c-----------------------------------------------------------------------
c
      if(itt.lt.0) then
        call anadinfo
# ifdef coast
c
c-----------------------------------------------------------------------
c  Setup coastal segments and land/sea mask at the tracer points.
c-----------------------------------------------------------------------
c
        call anacoast
# endif
#endif
#ifndef analytical
c
c-----------------------------------------------------------------------
# ifdef gridold
c  Read in bottom topography (meters) at the tracer points, HR.  Scale
# else
c  Read in bottom topography (meters) at the tracer points, HD.  Scale
# endif
c  to centimeters.
c-----------------------------------------------------------------------
c
        varid=ncvid(ncinpid,'tbath',rcode)
        if(rcode.eq.0) then
          start(1)=1
          count(1)=imt
          start(2)=1
          count(2)=jmt
# ifdef gridold
          call ncvgt(ncinpid,varid,start,count,hr,rcode)
# else
          call ncvgt(ncinpid,varid,start,count,hd,rcode)
# endif
          if(rcode.ne.0) then
            write(stdout,910) 'tbath'
            call exitus('INFLD')
          endif
          do 30 jj=1,jmt
          do 30 i=1,imt
# ifdef gridold
            hr(i,jj)=m2cm*abs(hr(i,jj))
# else
            hd(i,jj)=m2cm*abs(hd(i,jj))
# endif
  30      continue
        else
          write(stdout,920) 'tbath'
          call exitus('INFLD')
        endif
# ifndef gridold
c
c-----------------------------------------------------------------------
c  Read in bottom topography (meters) at the velocity points, HV.  Scale
c  to centimeters.
c-----------------------------------------------------------------------
c
        varid=ncvid(ncinpid,'vbath',rcode)
        if(rcode.eq.0) then
          start(1)=1
          count(1)=imt
          start(2)=1
          count(2)=jmt
          call ncvgt(ncinpid,varid,start,count,hv,rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'vbath'
            call exitus('INFLD')
          endif
          do 3001 jj=1,jmt
          do 3001 i=1,imt
            hv(i,jj)=m2cm*abs(hv(i,jj))
 3001    continue
        else
          write(stdout,920) 'vbath'
          call exitus('INFLD')
        endif
#  ifdef dblsigma
c
c-----------------------------------------------------------------------
c  Read in interface (meters) at the tracer points, tinter.  Scale
c  to centimeters.
c-----------------------------------------------------------------------
c
        varid=ncvid(ncinpid,'interface',rcode)
        if(rcode.eq.0) then
          start(1)=1
          count(1)=imt
          start(2)=1
          count(2)=jmt
          call ncvgt(ncinpid,varid,start,count,tinter,rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'interface'
            call exitus('INFLD')
          endif
          do 3002 jj=1,jmt
          do 3002 i=1,imt
            tinter(i,jj)=m2cm*abs(tinter(i,jj))
 3002    continue
        else
          write(stdout,920) 'interface'
          call exitus('INFLD')
        endif
c
c-----------------------------------------------------------------------
c  Read in interface (meters) at the velocity points, vinter.  Scale
c  to centimeters.
c-----------------------------------------------------------------------
c
        varid=ncvid(ncinpid,'viface',rcode)
        if(rcode.eq.0) then
          start(1)=1
          count(1)=imt
          start(2)=1
          count(2)=jmt
          call ncvgt(ncinpid,varid,start,count,vinter,rcode)
          if(rcode.ne.0) then
            write(stdout,910) 'vinter'
            call exitus('INFLD')
          endif
          do 3003 jj=1,jmt
          do 3003 i=1,imt
            vinter(i,jj)=m2cm*abs(vinter(i,jj))
 3003    continue
        else
          write(stdout,920) 'vinter'
          call exitus('INFLD')
        endif
#  endif
# endif
#else
c
c-----------------------------------------------------------------------
c  Set Analytical bottom topography at the tracer points.
c-----------------------------------------------------------------------
c
        call anatopo
#endif
      elseif(itt.eq.0) then
#ifndef analytical
c
c-----------------------------------------------------------------------
c  Read in initial conditions for row J.
c-----------------------------------------------------------------------
c
c  Read in internal mode velocity components.
c
        start(1)=xindx
        count(1)=1
        start(2)=1
        count(2)=km
        start(3)=1
        count(3)=imt
        start(4)=j
        count(4)=1
        start(5)=t0inpindx
        count(5)=1
        call ncvgt(ncinpid,viniid,start,count,fld,rcode)
        if(rcode.ne.0) then
          write(stdout,950) 'vclin x-component',j
          call exitus('INFLD')
        endif
        call flip(u,imt,km,fld,km,imt)
        start(1)=yindx
        count(1)=1
        call ncvgt(ncinpid,viniid,start,count,fld,rcode)
        if(rcode.ne.0) then
          write(stdout,950) 'vclin y-component',j
          call exitus('INFLD')
        endif
        call flip(v,imt,km,fld,km,imt)
c
c  Read in tracers (temperature, salinity, and other tracers if any).
c
        start(1)=1
        count(1)=km
        start(2)=1
        count(2)=imt
        start(3)=j
        count(3)=1
        start(4)=t0inpindx
        count(4)=1
        do 40 m=1,nt
          call ncvgt(ncinpid,tiniid(m),start,count,fld,rcode)
          if(rcode.ne.0) then
            write(stdout,960) 'trc',j,m
            call exitus('INFLD')
          endif
          call flip(t(1,1,m),imt,km,fld,km,imt)
  40    continue
c
c  Scale salinity to gain accuracy.
c
        do 50 k=1,km
        do 50 i=1,imt
          t(i,k,2)=t(i,k,2)-smean
  50    continue
c
c  Read in transport streamfunction.
c
        start(1)=1
        count(1)=imt
        start(2)=j
        count(2)=1
        start(3)=t0inpindx
        count(3)=1
        call ncvgt(ncinpid,piniid,start,count,p(1,j),rcode)
        if(rcode.ne.0) then
          write(stdout,950) 'pbar',j
          call exitus('INFLD')
        endif
#else
c
c-----------------------------------------------------------------------
c  Set Analytical initial conditions for row J.
c-----------------------------------------------------------------------
c
        call anauvpt(j)
#endif
c
        if(j.eq.jmt) then
c
c-----------------------------------------------------------------------
c  If last row of initialization calls, calculate mean background
c  density field.
c-----------------------------------------------------------------------
c
#ifndef analytical
          call meanrho
#else
          call anamrho
c
c-----------------------------------------------------------------------
c  If last row of initialization calls, set depths of boundary points.
c-----------------------------------------------------------------------
c
          call dpthbndry
#endif
c
c-----------------------------------------------------------------------
c  If last row of initialization calls, set length of boundary
c  condition vectors.
c-----------------------------------------------------------------------
c
          imx(1)=jmt
          imx(2)=imt
          imx(3)=jmt
          imx(4)=imt
#ifndef analytical
c
c-----------------------------------------------------------------------
c  If last row of initialization calls, read in first set of boundary
c  conditions.
c-----------------------------------------------------------------------
c
          do 80 ib=1,4
            start(2)=1
            count(2)=km
            start(3)=1
            count(3)=imx(ib)
            start(4)=ib
            count(4)=1
            start(5)=tinpindx
            count(5)=1
c
c  Internal mode velocity components.
c
            start(1)=xindx
            count(1)=1
            call ncvgt(ncinpid,vbryid,start,count,fld,rcode)
            if(rcode.ne.0) then
              write(stdout,970) 'vbry x-component',ib
              call exitus('INFLD')
            endif
            call flip(u1(1,1,ib),ijmx,km,fld,km,imx(ib))
            start(1)=yindx
            count(1)=1
            call ncvgt(ncinpid,vbryid,start,count,fld,rcode)
            if(rcode.ne.0) then
              write(stdout,970) 'vbry y-component',ib
              call exitus('INFLD')
            endif
            call flip(v1(1,1,ib),ijmx,km,fld,km,imx(ib))
c
c  Tracers (temperature, salinity, and other tracers if any).
c
            do 60 m=1,nt
              start(1)=m
              count(1)=1
              call ncvgt(ncinpid,tbryid,start,count,fld,rcode)
              if(rcode.ne.0) then
                write(stdout,980) 'tbry',ib,m
                call exitus('INFLD')
              endif
              call flip(t1(1,1,ib,m),ijmx,km,fld,km,imx(ib))
  60        continue
c
c  Scale salinity to gain accuracy.
c
            do 70 k=1,km
            do 70 i=1,imx(ib)
              t1(i,k,ib,2)=t1(i,k,ib,2)-smean
  70        continue
c
c  Transport streamfunction and time rate of change of vorticity.
c
            start(1)=1
            count(1)=imx(ib)
            start(2)=ib
            count(2)=1
            start(3)=tinpindx
            count(3)=1
            call ncvgt(ncinpid,pbryid,start,count,p1(1,1,ib),rcode)
            if(rcode.ne.0) then
              write(stdout,970) 'pbry',ib
              call exitus('INFLD')
            endif
            call ncvgt(ncinpid,qbryid,start,count,p1(1,2,ib),rcode)
            if(rcode.ne.0) then
              write(stdout,970) 'qbry',ib
              call exitus('INFLD')
            endif
  80      continue
c
c  Initialization time (usually, modified Julian day).
c
          start(1)=t0inpindx
          count(1)=1
          call ncvgt(ncinpid,t0inpid,start,count,ttbdyo,rcode)
          if(rcode.ne.0) then
            write(stdout,990) 'time0',t0inpindx
            call exitus('INFLD')
          endif
          ttbdy1=ttbdyo
#else
c
c-----------------------------------------------------------------------
c  If last row, set initial Analytical boundary conditions.
c-----------------------------------------------------------------------
c
          call anabndry
#endif
#ifdef coast
          call msk_bdy
#endif          
        endif
      else
#ifndef analytical
c
c=======================================================================
c  Set boundary conditions for this timestep.
c=======================================================================
c
c  1st, save last read set of boundary conditions
c
        ttbdy=(ttbdy1-ttbdyo)*day2sec
        if(ttbdy.lt.ttsec) then
          do 120 ib=1,4
            do 110 i=1,imx(ib)
              po(i,1,ib)=p1(i,1,ib)
              po(i,2,ib)=p1(i,2,ib)
              do 100 k=1,km
                uo(i,k,ib)=u1(i,k,ib)
                vo(i,k,ib)=v1(i,k,ib)
                do 90 m=1,nt
                  to(i,k,ib,m)=t1(i,k,ib,m)
  90            continue
 100          continue
 110        continue
 120      continue
c
c-----------------------------------------------------------------------
c  2nd, if not persistent boundary conditions, read in next set of
c  boundary conditions
c-----------------------------------------------------------------------
c
          if(iflag(3).ne.1) then
c
c  Increase counter for input time.
c
            tinpindx=tinpindx+1
            if(tinpindx.gt.tinpmax) then
              ttbdy1=ttbdyo+dayper
            else
c
c  Read in four boundary edges (1: west, 2: south, 3: east, 4: north).
c
              do 150 ib=1,4
                start(2)=1
                count(2)=km
                start(3)=1
                count(3)=imx(ib)
                start(4)=ib
                count(4)=1
                start(5)=tinpindx
                count(5)=1
c
c  Internal mode velocity components.
c
                start(1)=xindx
                count(1)=1
                call ncvgt(ncinpid,vbryid,start,count,fld,rcode)
                if(rcode.ne.0) then
                  write(stdout,970) 'vbry',ib
                  call exitus('INFLD')
                endif
                call flip(u1(1,1,ib),ijmx,km,fld,km,imx(ib))
                start(1)=yindx
                count(1)=1
                call ncvgt(ncinpid,vbryid,start,count,fld,rcode)
                if(rcode.ne.0) then
                  write(stdout,970) 'vbry',ib
                  call exitus('INFLD')
                endif
                call flip(v1(1,1,ib),ijmx,km,fld,km,imx(ib))
c
c  Tracers (temperature, salinity, and other tracers if any).
c
                do 130 m=1,nt
                  start(1)=m
                  count(1)=1
                  call ncvgt(ncinpid,tbryid,start,count,fld,rcode)
                  if(rcode.ne.0) then
                    write(stdout,980) 'tbry',ib,m
                    call exitus('INFLD')
                  endif
                  call flip(t1(1,1,ib,m),ijmx,km,fld,km,imx(ib))
 130            continue
c
c  Scale salinity to gain accuracy.
c
                do 140 k=1,km
                do 140 i=1,imx(ib)
                  t1(i,k,ib,2)=t1(i,k,ib,2)-smean
 140            continue
c
c  Transport streamfunction and time rate of change of vorticity.
c
                start(1)=1
                count(1)=imx(ib)
                start(2)=ib
                count(2)=1
                start(3)=tinpindx
                count(3)=1
                call ncvgt(ncinpid,pbryid,start,count,p1(1,1,ib),rcode)
                if(rcode.ne.0) then
                  write(stdout,970) 'pbry',ib
                  call exitus('INFLD')
                endif
                call ncvgt(ncinpid,qbryid,start,count,p1(1,2,ib),rcode)
                if(rcode.ne.0) then
                  write(stdout,970) 'qbry',ib
                  call exitus('INFLD')
                endif
 150          continue
c
c  Time since initialization.
c
              start(1)=tinpindx
              count(1)=1
              call ncvgt(ncinpid,tinpid,start,count,ttbdy1,rcode)
              if(rcode.ne.0) then
                write(stdout,990) 'time',tinpindx
                call exitus('INFLD')
              endif
            endif
          else
            ttbdy1=ttbdyo+dayper
          endif
        endif
c
c  3rd, if not persistent boundary conditions, interpolate boundary
c  conditions between updates.
c
        if(iflag(3).ne.1) then
          ttbdy=(ttbdy1-ttbdyo)*day2sec
          dt1=ttbdy-ttsec
          dt2=dtts
          dt3=dt1+dt2
          do 190 ib=1,4
            do 180 i=1,imx(ib)
              po(i,1,ib)=(po(i,1,ib)*dt1+p1(i,1,ib)*dt2)/dt3
              po(i,2,ib)=(po(i,2,ib)*dt1+p1(i,2,ib)*dt2)/dt3
              do 170 k=1,kmax
                uo(i,k,ib)=(uo(i,k,ib)*dt1+u1(i,k,ib)*dt2)/dt3
                vo(i,k,ib)=(vo(i,k,ib)*dt1+v1(i,k,ib)*dt2)/dt3
                do 160 m=1,nt
                  to(i,k,ib,m)=(to(i,k,ib,m)*dt1+t1(i,k,ib,m)*dt2)/dt3
 160            continue
 170          continue
 180        continue
 190      continue
        endif
#else
c
c=======================================================================
c  Set Analytical boundary conditions for this timestep.
c=======================================================================
c
        call anabndry
#endif
#ifdef coast
        call msk_bdy
#endif        
      endif
#ifdef oias
c
c=======================================================================
c  Begin intermittent OI assimilation.
c=======================================================================
c
      call oiassi(j)
#endif
c
#ifndef analytical
 900  format(/,' INFLD - unable to open initialization NetCDF file: ',a)
 910  format(/,' INFLD - error while reading variable: ',a)
 920  format(/,' INFLD - cannot find variable: ',a,2x,
     *         ' in initialization NetCDF file.')
 930  format(/,' INFLD - inconsistent dimension parameters, ',a,2i4)
 940  format(/,' INFLD - cannot find dimension: ',a,2x,
     *         ' in initialization NetCDF file.')
 950  format(/,' INFLD - error while reading variable: ',a,2x,
     *         ' at row J = ',i5)
 960  format(/,' INFLD - error while reading variable: ',a,2x,
     *         ' at row J = ',i5,2x,' and TRACER = ',i2)
 970  format(/,' INFLD - error while reading variable: ',a,2x,
     *         ' at boundary = ',i1)
 980  format(/,' INFLD - error while reading variable: ',a,2x,
     *         ' at boundary = ',i1,' and TRACER = ',i2)
 990  format(/,' INFLD - error while reading variable: ',a,2x,
     *         ' at index = ',i4)
1000  format(/,' INFLD - invalid value ',a,'=',i10)
#  ifdef dblsigma
1010  format(/'Double-Sigma Parameters:'//
     *  f10.2,'  ZC1     minimum depth of the coordinate interface ',
     *                  '(cm)',/,
     *  f10.2,'  ZC2     maximum depth of the coordinate interface ',
     *                  '(cm)',/,
     *  f10.4,'  ZSLOPE  slope parameter of the coordinate interface',/,
     *  f10.2,'  ZREF    reference depth for the coordinate interface ',
     *                  '(cm)'/)
#  endif
1020  format (/'+++Warning:  INFLD - old input file.'/13x,
     &        'Replacing domain values (DELX,DELY) by (0,0)')
1030  format (/'+++Warning:  INFLD - old input file.'/13x,
     &        'Replacing domain values (DELX,DELY) by (',f8.3,', ',f8.3,
     &        ')'/13x,'and domain values (RLNGD,RLATD) by (0,0)')
#endif
      return
      end
