!
!##################################################################
!##################################################################
!######                                                      ######
!######                   PROGRAM ARPSTINTRP                 ######
!######                                                      ######
!######                     Developed by                     ######
!######     Center for Analysis and Prediction of Storms     ######
!######                University of Oklahoma                ######
!######                                                      ######
!##################################################################
!##################################################################
!


PROGRAM arpstintrp,100
!
!-----------------------------------------------------------------------
!
!  PURPOSE:
!
!  This program interpolates two ARPS history data on grid of the same
!  size to a time inbetween them. The output will be written into a new
!  history dump file.
!
!-----------------------------------------------------------------------
!
!  AUTHOR: Ming Xue
!
!  2/25/1999. Written based on ARPSTINTRP.
!
!  MODIFICATION HISTORY:
!
!  1999/10/13 (Gene Bassett)
!  Corrected a roundoff error problem for the history dump reference
!  times.  Made the history dump output characterists similar
!  to arpsintrp.
!
!  2001/06/18 (Gene Bassett)
!  Corrected error with absolute time (iabstinit variables).
!
!  2002/03/19 (Keith Brewster)
!  Corrected time calculations for the case when the user chooses
!  to have curtim be relative to initime of the input file.
!
!-----------------------------------------------------------------------
!
!-----------------------------------------------------------------------
!
!  Variable Declarations.
!
!-----------------------------------------------------------------------
!
  IMPLICIT NONE
!
!-----------------------------------------------------------------------
!
!  Dimension of the base grid (input data).
!
!-----------------------------------------------------------------------
!
  INTEGER :: nx,ny,nz
  INTEGER :: nstyps
!
!-----------------------------------------------------------------------
!
!  ARPS arrays. The last dimension is for the two sets of variables.
!
!-----------------------------------------------------------------------
!
  REAL, ALLOCATABLE :: u     (:,:,:,:)  ! Total u-velocity (m/s).
  REAL, ALLOCATABLE :: v     (:,:,:,:)  ! Total v-velocity (m/s).
  REAL, ALLOCATABLE :: w     (:,:,:,:)  ! Total w-velocity (m/s).
  REAL, ALLOCATABLE :: ptprt (:,:,:,:)  ! Perturbation potential temperature
                              ! from that of base state atmosphere (Kelvin).
  REAL, ALLOCATABLE :: pprt  (:,:,:,:)  ! Perturbation pressure from that
                              ! of base state atmosphere (Pascal).
  REAL, ALLOCATABLE :: qv    (:,:,:,:)  ! Water vapor mixing ratio (kg/kg).
  REAL, ALLOCATABLE :: qc    (:,:,:,:)  ! Cloud water mixing ratio (kg/kg).
  REAL, ALLOCATABLE :: qr    (:,:,:,:)  ! Rain water mixing ratio (kg/kg).
  REAL, ALLOCATABLE :: qi    (:,:,:,:)  ! Cloud ice mixing ratio (kg/kg).
  REAL, ALLOCATABLE :: qs    (:,:,:,:)  ! Snow mixing ratio (kg/kg).
  REAL, ALLOCATABLE :: qh    (:,:,:,:)  ! Hail mixing ratio (kg/kg).
  REAL, ALLOCATABLE :: tke   (:,:,:,:)  ! Turbulent Kinetic Energy ((m/s)**2)
  REAL, ALLOCATABLE :: kmh   (:,:,:,:)  ! Horizontal turb. mixing coef. for
                              ! momentum. ( m**2/s )
  REAL, ALLOCATABLE :: kmv   (:,:,:,:)  ! Vertical turb. mixing coef. for
                              ! momentum. ( m**2/s )

  INTEGER, ALLOCATABLE :: soiltyp (:,:,:)   ! Soil type
  REAL, ALLOCATABLE :: stypfrct(:,:,:)   ! Fraction of soil type
  INTEGER, ALLOCATABLE :: vegtyp  (:,:)          ! Vegetation type
  REAL, ALLOCATABLE :: roufns  (:,:)          ! Surface roughness
  REAL, ALLOCATABLE :: lai     (:,:)          ! Leaf Area Index
  REAL, ALLOCATABLE :: veg     (:,:)          ! Vegetation fraction

  REAL, ALLOCATABLE :: tsfc   (:,:,:,:)   ! Ground sfc. temperature (K)
  REAL, ALLOCATABLE :: tsoil  (:,:,:,:)   ! Deep soil temperature (K)
  REAL, ALLOCATABLE :: wetsfc (:,:,:,:)   ! Surface soil moisture
  REAL, ALLOCATABLE :: wetdp  (:,:,:,:)   ! Deep soil moisture
  REAL, ALLOCATABLE :: wetcanp(:,:,:,:)   ! Canopy water amount
  REAL, ALLOCATABLE :: snowdpth(:,:,:)           ! Snow cover

  REAL, ALLOCATABLE :: raing(:,:,:)         ! Grid supersaturation rain
  REAL, ALLOCATABLE :: rainc(:,:,:)         ! Cumulus convective rain
  REAL, ALLOCATABLE :: prcrate(:,:,:,:)     ! precipitation rate (kg/(m**2*s))
                                 ! prcrate(:,:,:,:) = total precip. rate
                                 ! prcrate(:,:,:,:) = grid scale precip. rate
                                 ! prcrate(:,:,:,:) = cumulus precip. rate
                                 ! prcrate(:,:,:,:) = microphysics precip. rate

  REAL, ALLOCATABLE :: radfrc(:,:,:,:) ! Radiation forcing (K/s)
  REAL, ALLOCATABLE :: radsw (:,:,:)    ! Solar radiation reaching the surface
  REAL, ALLOCATABLE :: rnflx (:,:,:)    ! Net radiation flux absorbed by surface

  REAL, ALLOCATABLE :: usflx (:,:,:)    ! Surface flux of u-momentum (kg/(m*s**2))
  REAL, ALLOCATABLE :: vsflx (:,:,:)    ! Surface flux of v-momentum (kg/(m*s**2))
  REAL, ALLOCATABLE :: ptsflx(:,:,:)    ! Surface heat flux (K*kg/(m*s**2))
  REAL, ALLOCATABLE :: qvsflx(:,:,:)    ! Surface moisture flux (kg/(m**2*s))

  REAL, ALLOCATABLE :: ubar  (:,:,:)   ! Base state u-velocity (m/s).
  REAL, ALLOCATABLE :: vbar  (:,:,:)   ! Base state v-velocity (m/s).
  REAL, ALLOCATABLE :: wbar  (:,:,:)   ! Base state w-velocity (m/s).
  REAL, ALLOCATABLE :: ptbar (:,:,:)   ! Base state potential temperature (K)
  REAL, ALLOCATABLE :: pbar  (:,:,:)   ! Base state pressure (Pascal).
  REAL, ALLOCATABLE :: rhobar(:,:,:)   ! Base state air density (kg/m**3)
  REAL, ALLOCATABLE :: qvbar (:,:,:)   ! Base state water vapor specific humidity
                                 ! (kg/kg,2).
  REAL, ALLOCATABLE :: x     (:)         ! The x-coord. of the physical and
                                 ! computational grid. Defined at u-point.
  REAL, ALLOCATABLE :: y     (:)         ! The y-coord. of the physical and
                                 ! computational grid. Defined at v-point.
  REAL, ALLOCATABLE :: z     (:)         ! The z-coord. of the computational grid.
                                 ! Defined at w-point on the staggered grid.
  REAL, ALLOCATABLE :: zp    (:,:,:)   ! The physical height coordinate defined at
                                 ! w-point on the staggered grid.
  REAL, ALLOCATABLE :: hterain(:,:)     ! The height of terrain.

  REAL, ALLOCATABLE :: uprt   (:,:,:,:)    ! Perturbation u-velocity (m/s)
  REAL, ALLOCATABLE :: vprt   (:,:,:,:)    ! Perturbation v-velocity (m/s)
  REAL, ALLOCATABLE :: qvprt  (:,:,:,:)    ! Perturbation water vapor specific humidity (kg/kg)

  REAL, ALLOCATABLE :: tem1  (:,:,:)       ! Temporary array
  REAL, ALLOCATABLE :: tem2  (:,:,:)       ! Temporary array
  REAL, ALLOCATABLE :: tem3  (:,:,:)       ! Temporary array
  REAL, ALLOCATABLE :: tem4  (:,:,:)       ! Temporary array
!
!-----------------------------------------------------------------------
!
!  Misc. local variables:
!
!-----------------------------------------------------------------------
!
  INTEGER :: i, j, k
  REAL :: amin, amax

  CHARACTER (LEN=80) :: basdmpfn
  INTEGER :: lbasdmpf
  CHARACTER (LEN=80) :: ternfn,sfcoutfl,soiloutfl,temchar
  INTEGER :: lternfn,lfn
  INTEGER :: iss,is

!wdt Copyright (c) 2001 Weather Decision Technologies, Inc.
  REAL :: zpmax
!
!-----------------------------------------------------------------------
!
!  Include files:
!
!-----------------------------------------------------------------------
!
  INCLUDE 'phycst.inc'
  INCLUDE 'globcst.inc'
  INCLUDE 'bndry.inc'
  INCLUDE 'indtflg.inc'
  INCLUDE 'grid.inc'
!wdt Copyright (c) 2001 Weather Decision Technologies, Inc.
  INCLUDE 'exbc.inc'

  INTEGER :: hinfmt,houtfmt, nchin, nchout
  CHARACTER (LEN=80) :: filename
  CHARACTER (LEN=80) :: grdbasfn
  INTEGER :: lenfil,lengbf

  INTEGER :: grdbas
  INTEGER :: ireturn

  REAL :: time
  INTEGER :: gboutcnt, vroutcnt

  DATA    gboutcnt, vroutcnt /0,0/

  INTEGER :: nfilemax
  PARAMETER (nfilemax=2)

  CHARACTER (LEN=80) :: hisfile(nfilemax)
  INTEGER :: nhisfile,nd, length, lenstr
  CHARACTER (LEN=80) :: timsnd
  CHARACTER (LEN=80) :: new_runname
  INTEGER :: tmstrln

  REAL :: times(nfilemax), outtime, alpha, beta
  INTEGER :: iabstinit,iabstinit1
  INTEGER :: ioffset
  INTEGER :: year1,month1,day1,hour1,minute1,second1,ioutabst,          &
          ioutabstinit
!
!-----------------------------------------------------------------------
!
!  namelist Declarations:
!
!-----------------------------------------------------------------------
!
  INTEGER :: use_data_t
  CHARACTER (LEN=19) :: initime  ! Real time in form of 'year-mo-dy:hr:mn:ss'

  NAMELIST /INPUT/hinfmt,nhisfile,grdbasfn,hisfile

!wdt Copyright (c) 2001 Weather Decision Technologies, Inc.: ngbrz,zbrdmp
  NAMELIST /output/ runname,use_data_t,initime,outtime,                 &
            dirname,exbcdmp,hdmpfmt,grbpkbit,hdfcompr,                  &
            grdout,basout,varout,mstout,rainout,prcout,iceout,          &
            tkeout, trbout,sfcout,landout,radout,flxout,                &
            qcexout,qrexout,qiexout,qsexout,qhexout,                    &
            totout,filcmprs,sfcdmp,soildmp,ngbrz,zbrdmp

  INTEGER :: nsize, nxy,nxyz
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
!  Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
  mgrid = 1
  nestgrd = 0
!
!-----------------------------------------------------------------------
!
!  Set the default parameters
!
!-----------------------------------------------------------------------
!
  hinfmt    = 10
  grdbasfn  = 'X'
  nhisfile  = 1
  hisfile(1) = 'X'
  use_data_t = 1
  initime ='0000-00-00:00:00:00'

  runname   = 'runname_not_set'

  dirname   = './'
  exbcdmp   = 1
  hdmpfmt   = 1
  grbpkbit  = 16
  hdfcompr  = 0
  filcmprs  = 0
  basout    = 0
  grdout    = 0
  varout    = 1
  mstout    = 1
  iceout    = 1
  tkeout    = 1
  trbout    = 0
  rainout   = 0
  sfcout    = 0
  snowout   = 0
  landout   = 0
  qcexout   = 0
  qrexout   = 0
  qiexout   = 0
  qsexout   = 0
  qhexout   = 0
  sfcdmp    = 1
  soildmp   = 1
!wdt Copyright (c) 2001 Weather Decision Technologies, Inc.: ngbrz,zbrdmp
  ngbrz = 5
  zbrdmp = 10000.0

  WRITE(6,'(/9(/2x,a)/)')                                               &
     '###############################################################', &
     '###############################################################', &
     '###                                                         ###', &
     '###                Welcome to ARPSTINTRP                    ###', &
     '###                                                         ###', &
     '###############################################################', &
     '###############################################################'

  READ(5,INPUT, END=100)

  nhisfile = 2

  WRITE(6,'(/a/)') ' Input control parameters read in are:'

  IF( hinfmt == 5) THEN
    WRITE(6,'(2(2x,a))')                                                &
        'The Savi3D not supported as an INPUT file format',             &
        'Job stopped in ARPSTINTRP.'
    STOP 9102
  END IF

  WRITE(6,'(1x,a,i3)') ' hinfmt      =', hinfmt
  WRITE(6,'(1x,a,i3)') ' nhisfile    =', nhisfile

  length = LEN( grdbasfn )
  CALL strlnth(  grdbasfn, length )
  WRITE(6,'(1x,a,a)')  ' grdbasfn    =', grdbasfn(1:length)

  DO i=1,nhisfile
    length = LEN( hisfile(i) )
    CALL strlnth(  hisfile(i), length )
    WRITE(6,'(1x,a,i3,a,a)') ' hisfile(',i,')=',hisfile(i)(1:length)
  END DO

!
!-----------------------------------------------------------------------
!
!  Set the control parameters for output:
!
!-----------------------------------------------------------------------
!
  WRITE(6,'(/a/)')                                                      &
      ' Reading in control parameters for the output data files..'

  READ (5,output,END=100)

!wdt Copyright (c) 2001 Weather Decision Technologies, Inc.
  IF ( exbcdmp == 4 ) rayklow = -1

  WRITE(6,'(/2x,a,a)')                                                  &
      'The run name to be used for constructing output file name is ',  &
      runname

  new_runname = runname

  totout = 1

  CALL get_dims_from_data(hinfmt,trim(grdbasfn),                    &
       nx,ny,nz,nstyps, ireturn)

  IF (nstyps <= 0) nstyps = 1
  nstyp = nstyps
!
!-----------------------------------------------------------------------
!
!  Allocate the arrays.
!
!-----------------------------------------------------------------------
!
  ALLOCATE(  ALLOCATE(  ALLOCATE(  ALLOCATE(ptprt(nx,ny,nz,2))
  ALLOCATE(pprt(nx,ny,nz,2))
  ALLOCATE(qv(nx,ny,nz,2))
  ALLOCATE(qc(nx,ny,nz,2))
  ALLOCATE(qr(nx,ny,nz,2))
  ALLOCATE(qi(nx,ny,nz,2))
  ALLOCATE(qs(nx,ny,nz,2))
  ALLOCATE(qh(nx,ny,nz,2))
  ALLOCATE(tke(nx,ny,nz,2))
  ALLOCATE(kmh(nx,ny,nz,2))
  ALLOCATE(kmv(nx,ny,nz,2))

  ALLOCATE(soiltyp(nx,ny,nstyps))
  ALLOCATE(stypfrct(nx,ny,nstyps))
  ALLOCATE(vegtyp(nx,ny))
  ALLOCATE(roufns(nx,ny))
  ALLOCATE(lai(nx,ny))
  ALLOCATE(veg(nx,ny))

  ALLOCATE(tsfc(nx,ny,0:nstyps,2))
  ALLOCATE(tsoil(nx,ny,0:nstyps,2))
  ALLOCATE(wetsfc(nx,ny,0:nstyps,2))
  ALLOCATE(wetdp(nx,ny,0:nstyps,2))
  ALLOCATE(wetcanp(nx,ny,0:nstyps,2))
  ALLOCATE(snowdpth(nx,ny,2))

  ALLOCATE(raing(nx,ny,2))
  ALLOCATE(rainc(nx,ny,2))
  ALLOCATE(prcrate(nx,ny,4,2))

  ALLOCATE(radfrc(nx,ny,nz,2))
  ALLOCATE(radsw(nx,ny,2))
  ALLOCATE(rnflx(nx,ny,2))

  ALLOCATE(usflx(nx,ny,2))
  ALLOCATE(vsflx(nx,ny,2))
  ALLOCATE(ptsflx(nx,ny,2))
  ALLOCATE(qvsflx(nx,ny,2))

  ALLOCATE(ubar(nx,ny,nz))
  ALLOCATE(vbar(nx,ny,nz))
  ALLOCATE(wbar(nx,ny,nz))
  ALLOCATE(ptbar(nx,ny,nz))
  ALLOCATE(pbar(nx,ny,nz))
  ALLOCATE(rhobar(nx,ny,nz))
  ALLOCATE(qvbar(nx,ny,nz))
  ALLOCATE(  ALLOCATE(  ALLOCATE(  ALLOCATE(zp(nx,ny,nz))
  ALLOCATE(hterain(nx,ny))

  ALLOCATE(uprt(nx,ny,nz,2))
  ALLOCATE(vprt(nx,ny,nz,2))
  ALLOCATE(qvprt(nx,ny,nz,2))

  ALLOCATE(tem1(nx,ny,nz))
  ALLOCATE(tem2(nx,ny,nz))
  ALLOCATE(tem3(nx,ny,nz))
  ALLOCATE(tem4(nx,ny,nz))

  nxyz = nx*ny*nz
  nxy  = nx*ny

  CALL flzero(u     , nxyz*2)
  CALL flzero(v     , nxyz*2)
  CALL flzero(uprt  , nxyz*2)
  CALL flzero(vprt  , nxyz*2)
  CALL flzero(w     , nxyz*2)
  CALL flzero(ptprt , nxyz*2)
  CALL flzero(pprt  , nxyz*2)
  CALL flzero(qvprt , nxyz*2)
  CALL flzero(qc    , nxyz*2)
  CALL flzero(qr    , nxyz*2)
  CALL flzero(qi    , nxyz*2)
  CALL flzero(qs    , nxyz*2)
  CALL flzero(qh    , nxyz*2)
  CALL flzero(tke   , nxyz*2)
  CALL flzero(kmh   , nxyz*2)
  CALL flzero(kmv   , nxyz*2)
  CALL flzero(radfrc, nxyz*2)

  CALL flzero(ubar  , nxyz)
  CALL flzero(vbar  , nxyz)
  CALL flzero(wbar  , nxyz)
  CALL flzero(ptbar , nxyz)
  CALL flzero(pbar  , nxyz)
  CALL flzero(rhobar, nxyz)
  CALL flzero(qvbar , nxyz)

  CALL flzero(x, nx)
  CALL flzero(y, ny)
  CALL flzero(z, nz)
  CALL flzero(zp    , nxyz)
  CALL flzero(hterain,nx*ny)

  nsize = nx*ny*(1+nstyps)*2
  CALL flzero(tsfc   ,nsize)
  CALL flzero(tsoil  ,nsize)
  CALL flzero(wetsfc ,nsize)
  CALL flzero(tsfc   ,nsize)
  CALL flzero(wetdp  ,nsize)
  CALL flzero(wetcanp,nsize)

  DO iss=1,nstyps
    DO j=1,ny
      DO i=1,nx
        soiltyp (i,j,iss) = 0
        stypfrct(i,j,iss) = 0.0
      END DO
    END DO
  END DO

  DO is=1,2
    DO j=1,ny
      DO i=1,nx
        snowdpth(i,j,is) = 0

        vegtyp (i,j) = 0
        lai    (i,j) = 0.0
        roufns (i,j) = 0.0
        veg    (i,j) = 0.0

        raing  (i,j,is) = 0.0
        rainc  (i,j,is) = 0.0
        prcrate(i,j,1,is) = 0.0
        prcrate(i,j,2,is) = 0.0
        prcrate(i,j,3,is) = 0.0
        prcrate(i,j,4,is) = 0.0
        radsw  (i,j,is) = 0.0
        rnflx  (i,j,is) = 0.0
        usflx  (i,j,is) = 0.0
        vsflx  (i,j,is) = 0.0
        ptsflx (i,j,is) = 0.0
        qvsflx (i,j,is) = 0.0
      END DO
    END DO
  END DO


  ldirnam=LEN(dirname)
  CALL strlnth( dirname , ldirnam)

  lengbf=LEN(grdbasfn)
  CALL strlnth( grdbasfn, lengbf)
  WRITE(6,'(/a,a)')' The grid/base name is ', grdbasfn(1:lengbf)
!
!-----------------------------------------------------------------------
!
!  Loop over data files
!
!-----------------------------------------------------------------------
!
  ireturn = 0

  DO nd=1,2

    filename = hisfile(nd)

    lenfil=LEN(filename)
    CALL strlnth( filename, lenfil)
    WRITE(6,'(/a,a,a)')                                                 &
        ' Data set ', filename(1:lenfil) ,' to be processed.'
!
!-----------------------------------------------------------------------
!
!  Read all input data arrays
!
!-----------------------------------------------------------------------
!
    CALL dtaread(nx,ny,nz,nstyps,                                       &
        hinfmt,nchin,grdbasfn(1:lengbf),lengbf,                         &
        filename(1:lenfil),lenfil,time, x,y,z,zp,                       &
        uprt(1,1,1,nd),vprt (1,1,1,nd),w  (1,1,1,nd),ptprt(1,1,1,nd),   &
        pprt(1,1,1,nd),qvprt(1,1,1,nd),qc (1,1,1,nd),qr   (1,1,1,nd),   &
        qi  (1,1,1,nd),qs   (1,1,1,nd),qh (1,1,1,nd),tke  (1,1,1,nd),   &
        kmh (1,1,1,nd),kmv  (1,1,1,nd),                                 &
        ubar,vbar,wbar,ptbar,pbar,rhobar,qvbar,                         &
        soiltyp,stypfrct,vegtyp,lai,roufns,veg,                         &
        tsfc(1,1,0,nd),tsoil(1,1,0,nd),wetsfc(1,1,0,nd),                &
        wetdp(1,1,0,nd),wetcanp(1,1,0,nd),snowdpth(1,1,nd),             &
        raing(1,1,nd),rainc(1,1,nd),prcrate(1,1,1,nd),                  &
        radfrc(1,1,1,nd),radsw(1,1,nd),rnflx(1,1,nd),                   &
        usflx(1,1,nd),vsflx(1,1,nd),ptsflx(1,1,nd),qvsflx(1,1,nd),      &
        ireturn, tem1, tem2, tem3)

    WRITE(6,'(/a,/2(a,i2),a,i4,a,/3(a,i2),a,f13.3,a/)')                 &
        'History data read in for time: ',                              &
        'month=',month,',   day=',   day,',   year=',year,',',          &
        'hour =',hour ,',minute=',minute,', second=',second,            &
        ', time=',time,'(s)'

    CALL ctim2abss(year,month,day,hour,minute,second,iabstinit)

    IF(nd == 1) THEN  ! Save the values for data set 1
      year1  = year
      month1 = month
      day1   = day
      hour1  = hour
      minute1= minute
      second1= second
      iabstinit1 = iabstinit
    END IF

    times(nd) = time + int(iabstinit - iabstinit1)

  END DO

  IF( hinfmt == 9 .AND. ireturn == 2 ) THEN
    WRITE(6,'(/1x,a/)') 'The end of GrADS file was reached.'
    CLOSE ( nchin )
    CALL retunit( nchin )
    GO TO 9001
  END IF

  IF( ireturn /= 0 ) GO TO 9002            ! Read was unsuccessful
!

  IF( use_data_t == 1) THEN ! Use init time in input file
    ioutabstinit=iabstinit1
    year  = year1
    month = month1
    day   = day1
    hour  = hour1
    minute= minute1
    second= second1
  ELSE
    READ(initime,'(i4.4,1x,i2.2,1x,i2.2,1x,i2.2,1x,i2.2,1x,i2.2)')      &
                   year,month,day,hour,minute,second
    CALL ctim2abss(year,month,day,hour,minute,second,ioutabstinit)
  END IF

  PRINT*,'ioutabstinit=',ioutabstinit
  PRINT*,'iabstinit1=',iabstinit1

  ioffset = ioutabstinit - iabstinit1
  times(1) = times(1) - ioffset
  times(2) = times(2) - ioffset
  curtim = outtime 

  WRITE(6,'(/a,/2(a,i2),a,i4,a,/3(a,i2),a,f13.3,a/)')                   &
      'In output file, the reference time is',                          &
      'month=',month,',   day=',   day,',   year=',year,',',            &
      'hour =',hour ,',minute=',minute,', second=',second,              &
      ', & the time relative to this reference =',curtim

  IF ( curtim > MAX(times(1),times(2)) .OR.                             &
         curtim < MIN(times(1),times(2))) THEN
    WRITE (*,*) "WARNING: Performing extrapolation.  Desired time ",    &
                "is outside the range of the reference files."
  END IF

  IF( times(2) == times(1)) THEN
    WRITE (*,*) "ERROR: times in reference files are the same, ",       &
                "can't perform interpolation."
    STOP 1
  END IF

  alpha = (times(2)-curtim)/(times(2)-times(1))
  beta = 1.0-alpha

  WRITE (*,*)
  WRITE (*,*) "Relative weights: file 1",alpha,"   file 2",beta

!-----------------------------------------------------------------------
!
!  Calculate total fields from that for base state and perturbations
!
!-----------------------------------------------------------------------
!
  DO k=1,nz
    DO j=1,ny
      DO i=1,nx
        uprt  (i,j,k,1)=alpha*uprt  (i,j,k,1)+beta*uprt  (i,j,k,2)
        vprt  (i,j,k,1)=alpha*vprt  (i,j,k,1)+beta*vprt  (i,j,k,2)
        w     (i,j,k,1)=alpha*w     (i,j,k,1)+beta*w     (i,j,k,2)
        ptprt (i,j,k,1)=alpha*ptprt (i,j,k,1)+beta*ptprt (i,j,k,2)
        pprt  (i,j,k,1)=alpha*pprt  (i,j,k,1)+beta*pprt  (i,j,k,2)
        qvprt (i,j,k,1)=alpha*qvprt (i,j,k,1)+beta*qvprt (i,j,k,2)
        qc    (i,j,k,1)=alpha*qc    (i,j,k,1)+beta*qc    (i,j,k,2)
        qr    (i,j,k,1)=alpha*qr    (i,j,k,1)+beta*qr    (i,j,k,2)
        qi    (i,j,k,1)=alpha*qi    (i,j,k,1)+beta*qi    (i,j,k,2)
        qs    (i,j,k,1)=alpha*qs    (i,j,k,1)+beta*qs    (i,j,k,2)
        qh    (i,j,k,1)=alpha*qh    (i,j,k,1)+beta*qh    (i,j,k,2)
        tke   (i,j,k,1)=alpha*tke   (i,j,k,1)+beta*tke   (i,j,k,2)
        kmh   (i,j,k,1)=alpha*kmh   (i,j,k,1)+beta*kmh   (i,j,k,2)
        kmv   (i,j,k,1)=alpha*kmv   (i,j,k,1)+beta*kmv   (i,j,k,2)
        u     (i,j,k,1)=uprt (i,j,k,1)+ubar (i,j,k)
        v     (i,j,k,1)=vprt (i,j,k,1)+vbar (i,j,k)
        qv    (i,j,k,1)=qvprt(i,j,k,1)+qvbar(i,j,k)
        radfrc(i,j,k,1)=alpha*radfrc(i,j,k,1)+beta*radfrc(i,j,k,2)
      END DO
    END DO
  END DO

  DO is=0,nstyp
    DO j=1,ny-1
      DO i=1,nx-1
        tsfc   (i,j,is,1)=alpha*tsfc   (i,j,is,1)                       &
                          +beta*tsfc   (i,j,is,2)
        tsoil  (i,j,is,1)=alpha*tsoil  (i,j,is,1)                       &
                          +beta*tsoil  (i,j,is,2)
        wetsfc (i,j,is,1)=alpha*wetsfc (i,j,is,1)                       &
                          +beta*wetsfc (i,j,is,2)
        wetdp  (i,j,is,1)=alpha*wetdp  (i,j,is,1)                       &
                          +beta*wetdp  (i,j,is,2)
        wetcanp(i,j,is,1)=alpha*wetcanp(i,j,is,1)                       &
                          +beta*wetcanp(i,j,is,2)
      END DO
    END DO
  END DO

  DO j=1,ny-1
    DO i=1,nx-1
      snowdpth(i,j,1)=alpha*snowdpth(i,j,1)+beta*snowdpth(i,j,2)
      raing  (i,j,1)=alpha*raing  (i,j,1)+beta*raing  (i,j,2)
      rainc  (i,j,1)=alpha*rainc  (i,j,1)+beta*rainc  (i,j,2)

      prcrate(i,j,1,1)=alpha*prcrate(i,j,1,1)+beta*prcrate(i,j,1,2)
      prcrate(i,j,2,1)=alpha*prcrate(i,j,2,1)+beta*prcrate(i,j,2,2)
      prcrate(i,j,3,1)=alpha*prcrate(i,j,3,1)+beta*prcrate(i,j,3,2)
      prcrate(i,j,4,1)=alpha*prcrate(i,j,4,1)+beta*prcrate(i,j,4,2)

      radsw  (i,j,1)=alpha*radsw  (i,j,1)+beta*radsw  (i,j,2)
      rnflx  (i,j,1)=alpha*rnflx  (i,j,1)+beta*rnflx  (i,j,2)
      usflx  (i,j,1)=alpha*usflx  (i,j,1)+beta*usflx  (i,j,2)
      vsflx  (i,j,1)=alpha*vsflx  (i,j,1)+beta*vsflx  (i,j,2)
      ptsflx (i,j,1)=alpha*ptsflx (i,j,1)+beta*ptsflx (i,j,2)
      qvsflx (i,j,1)=alpha*qvsflx (i,j,1)+beta*qvsflx (i,j,2)
    END DO
  END DO

!
!-----------------------------------------------------------------------
!
!  Print out the max/min of output varaibles.
!
!-----------------------------------------------------------------------
!
  WRITE(6,'(/1x,a/)')                                                   &
      'Min. and max. of data interpolated to the new time:'

  CALL a3dmax0(x,1,nx,1,nx,1,1,1,1,1,1,1,1, amax,amin)
  WRITE(6,'(/1x,2(a,e13.6))') 'xmin    = ', amin,',  xmax    =',amax

  CALL a3dmax0(y,1,ny,1,ny,1,1,1,1,1,1,1,1, amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'ymin    = ', amin,',  ymax    =',amax

  CALL a3dmax0(z,1,nz,1,nz,1,1,1,1,1,1,1,1, amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'zmin    = ', amin,',  zmax    =',amax

  CALL a3dmax0(zp,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz,                    &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'zpmin   = ', amin,', zpmax    =',amax

  CALL a3dmax0(ubar,1,nx,1,nx,1,ny,1,ny-1,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'ubarmin = ', amin,',  ubarmax =',amax

  CALL a3dmax0(vbar,1,nx,1,nx-1,1,ny,1,ny,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'vbarmin = ', amin,',  vbarmax =',amax

  CALL a3dmax0(ptbar,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,               &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'ptbarmin= ', amin,',  ptbarmax=',amax

  CALL a3dmax0(pbar,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'pbarmin = ', amin,',  pbarmax =',amax

  CALL a3dmax0(rhobar,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,              &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'rhobarmin=', amin,', rhobarmax=',amax

  CALL a3dmax0(qvbar,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,               &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qvbarmin= ', amin,',  qvbarmax=',amax

  CALL a3dmax0(uprt,1,nx,1,nx,1,ny,1,ny-1,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'uprtmin = ', amin,',  uprtmax =',amax

  CALL a3dmax0(vprt,1,nx,1,nx-1,1,ny,1,ny,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'vprtmin = ', amin,',  vprtmax =',amax

  CALL a3dmax0(w,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz,                     &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'wmin    = ', amin,',  wmax    =',amax

  CALL a3dmax0(ptprt,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,               &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'ptprtmin= ', amin,',  ptprtmax=',amax

  CALL a3dmax0(pprt,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'pprtmin = ', amin,',  pprtmax =',amax

  CALL a3dmax0(qvprt,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,               &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qvprtmin= ', amin,',  qvprtmax=',amax

  CALL a3dmax0(qc,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qcmin   = ', amin,',  qcmax   =',amax

  CALL a3dmax0(qr,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qrmin   = ', amin,',  qrmax   =',amax

  CALL a3dmax0(qi,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qimin   = ', amin,',  qimax   =',amax

  CALL a3dmax0(qs,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qsmin   = ', amin,',  qsmax   =',amax

  CALL a3dmax0(qh,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                  &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qhmin   = ', amin,',  qhmax   =',amax

  CALL a3dmax0(tke,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                 &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'tkemin  = ', amin,',  tkemax  =',amax

  CALL a3dmax0(kmh,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                 &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'kmhmin  = ', amin,',  kmhmax  =',amax

  CALL a3dmax0(kmv,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,                 &
              amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'kmvmin  = ', amin,',  kmvmax  =',amax

  CALL a3dmax0(raing,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'raingmin= ', amin,',  raingmax=',amax

  CALL a3dmax0(rainc,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'raincmin= ', amin,',  raincmax=',amax

  CALL a3dmax0(prcrate(1,1,1,1),1,nx,1,nx-1,1,ny,1,ny-1,                &
               1,1,1,1,amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'prcr1min= ', amin,',  prcr1max=',amax

  CALL a3dmax0(prcrate(1,1,2,1),1,nx,1,nx-1,1,ny,1,ny-1,                &
               1,1,1,1,amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'prcr2min= ', amin,',  prcr2max=',amax

  CALL a3dmax0(prcrate(1,1,3,1),1,nx,1,nx-1,1,ny,1,ny-1,                &
               1,1,1,1,amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'prcr3min= ', amin,',  prcr3max=',amax

  CALL a3dmax0(prcrate(1,1,4,1),1,nx,1,nx-1,1,ny,1,ny-1,                &
               1,1,1,1,amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'prcr4min= ', amin,',  prcr4max=',amax

  DO iss = 0, nstyp

    CALL a3dmax0(tsfc(1,1,iss,1),1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,       &
                 amax,amin)
    WRITE(6,'(1x,2(a,e13.6),a,i3)')                                     &
        'tsfcmin = ', amin,',  tsfcmax =',amax,' for soil type=',iss

    CALL a3dmax0(tsoil(1,1,iss,1),1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,      &
                 amax,amin)
    WRITE(6,'(1x,2(a,e13.6),a,i3)')                                     &
        'tsoilmin= ', amin,',  tsoilmax=',amax,' for soil type=',iss

    CALL a3dmax0(wetsfc(1,1,iss,1),1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,     &
                 amax,amin)
    WRITE(6,'(1x,2(a,e13.6),a,i3)')                                     &
        'wetsmin = ', amin,',  wetsmax =',amax,' for soil type=',iss

    CALL a3dmax0(wetdp(1,1,iss,1),1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,      &
                 amax,amin)
    WRITE(6,'(1x,2(a,e13.6),a,i3)')                                     &
        'wetdmin = ', amin,',  wetdmax =',amax,' for soil type=',iss

    CALL a3dmax0(wetcanp(1,1,iss,1),                                    &
                 1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,amax,amin)
    WRITE(6,'(1x,2(a,e13.6),a,i3)')                                     &
        'wetcmin = ', amin,',  wetcmax =',amax,' for soil type=',iss

  END DO

  CALL a3dmax0(roufns,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                  &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'roufnmin =', amin,', roufnmax =',amax

  CALL a3dmax0(veg,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                     &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'vegmin  = ', amin,',   vegmax =',amax

  CALL a3dmax0(radfrc,1,nx,1,nx-1,1,ny,1,ny-1,1,nz,1,nz-1,              &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'radfnmin =', amin,', radfnmax =',amax

  CALL a3dmax0(radsw,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                   &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'radswmin =', amin,', radswmax =',amax

  CALL a3dmax0(rnflx,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                   &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'rnflxmin =', amin,', rnflxmax =',amax

  CALL a3dmax0(usflx,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                   &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'usflxnmin =', amin,', usflxmax =',amax

  CALL a3dmax0(vsflx,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                   &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'vsflxmin =', amin,', vsflxmax =',amax

  CALL a3dmax0(ptsflx,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                  &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'ptflxmin =', amin,', ptflxmax =',amax

  CALL a3dmax0(qvsflx,1,nx,1,nx-1,1,ny,1,ny-1,1,1,1,1,                  &
               amax,amin)
  WRITE(6,'(1x,2(a,e13.6))') 'qvflxmin =', amin,', qvflxmax =',amax

  DO k=1,nz
    DO j=1,ny
      DO i=1,nx
        tem1 (i,j,k)= 0.0        ! To be put in place of wbar
      END DO
    END DO
  END DO
!
!-----------------------------------------------------------------------
!
!  Data dump of the model grid and base state arrays:
!
!  First find a unique name basdmpfn(1:lbasdmpf) for the grid and
!  base state array dump file
!
!  If grid/base state data has been written out once, skip
!  the following writing block. Also no need to write out
!  separate data for Savi3D dump. The same for GrADS dump.
!
!-----------------------------------------------------------------------
!
  WRITE (cmnt(nocmnt),'(a,i4,a,i4,a,i4)')                               &
      ' nx =',nx,', ny =',ny,', nz =',nz

  runname = new_runname
  houtfmt = hdmpfmt
  grbpkbit = 16

  CALL gtlfnkey(runname, lfnkey)

  IF(houtfmt /= 9 ) THEN

    IF( gboutcnt == 1 ) GO TO 500 ! If done already, skip this part.

    CALL gtbasfn(runname(1:lfnkey),dirname,ldirnam,hdmpfmt,             &
                 1,0,basdmpfn,lbasdmpf)

    PRINT*
    PRINT*,'Output grid/base state file is ', basdmpfn(1:lbasdmpf)

    grdbas = 1      ! Dump out grd and base state arrays only

    CALL dtadump(nx,ny,nz,nstyps,hdmpfmt,nchout,                        &
                 basdmpfn(1:lbasdmpf),grdbas,filcmprs,                  &
                 u,v,w,ptprt,pprt,qv,qc,qr,qi,qs,qh,                    &
                 tke,kmh,kmv,                                           &
                 ubar,vbar,tem1,ptbar,pbar,rhobar,qvbar,                &
                 x,y,z,zp,hterain, tem1,tem1,tem1,                      &
                 soiltyp,stypfrct,vegtyp,lai,roufns,veg,                &
                 tsfc,tsoil,wetsfc,wetdp,wetcanp,snowdpth,              &
                 raing,rainc,prcrate,                                   &
                 radfrc,radsw,rnflx,                                    &
                 usflx,vsflx,ptsflx,qvsflx,                             &
                 tem2,tem3,tem4)

    gboutcnt = 1

    500     CONTINUE

  END IF
!
!-----------------------------------------------------------------------
!
!  Then the time dependent fields:
!
!-----------------------------------------------------------------------
!
  IF( .NOT. (houtfmt == 9 .AND. vroutcnt == 1) ) THEN
!
!-----------------------------------------------------------------------
!
!  Reconstruct the file name using the specified directory name
!
!-----------------------------------------------------------------------
!
    CALL gtdmpfn(runname(1:lfnkey),dirname,                             &
                 ldirnam,curtim,hdmpfmt,1,0, hdmpfn, ldmpf)

  END IF

  WRITE(6,'(a,a)') 'Writing t-dependent variable history dump ',        &
                    hdmpfn(1:ldmpf)
  grdbas = 0

  CALL dtadump(nx,ny,nz,nstyps,hdmpfmt,nchout,                          &
               hdmpfn(1:ldmpf),grdbas,filcmprs,                         &
               u,v,w,ptprt,pprt,qv,qc,qr,qi,qs,qh,                      &
               tke,kmh,kmv,                                             &
               ubar,vbar,tem1,ptbar,pbar,rhobar,qvbar,                  &
               x,y,z,zp,hterain, tem1,tem1,tem1,                        &
               soiltyp,stypfrct,vegtyp,lai,roufns,veg,                  &
               tsfc,tsoil,wetsfc,wetdp,wetcanp,snowdpth,                &
               raing,rainc,prcrate,                                     &
               radfrc,radsw,rnflx,                                      &
               usflx,vsflx,ptsflx,qvsflx,                               &
               tem2,tem3,tem4)
!
!-----------------------------------------------------------------------
!
!  Write out soil model variable file
!
!-----------------------------------------------------------------------
!
  IF ( sfcin == 1 ) THEN

    CALL cvttsnd( curtim, timsnd, tmstrln )

    soiloutfl = runname(1:lfnkey)//".soilvar."//timsnd(1:tmstrln)
    lfn = lfnkey + 9 + tmstrln

    IF( dirname /= ' ' ) THEN

      temchar = soiloutfl
      soiloutfl = dirname(1:ldirnam)//'/'//temchar
      lfn  = lfn + ldirnam + 1

    END IF

    CALL fnversn(soiloutfl, lfn)

    IF (soildmp > 0) THEN
      PRINT *, 'Writing soil data to ',soiloutfl(1:lfn)

      CALL wrtsoil(nx,ny,nstyps, soiloutfl(1:lfn),dx,dy,                &
               mapproj,trulat1,trulat2,trulon,sclfct,ctrlat,ctrlon,     &
               1,1,1,1,1,1,                                             &
               tsfc,tsoil,wetsfc,wetdp,wetcanp,snowdpth,soiltyp)

      IF (soildmp == 1) CALL soilcntl(nx,ny, soiloutfl(1:lfn),          &
                1,1,1,1,1,1, x,y)
    END IF

  END IF       ! sfcin.eq.1

!-----------------------------------------------------------------------
!
!  Write out surface property data file: sfcoutfl .
!
!-----------------------------------------------------------------------
!
  IF ( landin == 1 ) THEN

    sfcoutfl = runname(1:lfnkey)//".sfcdata"
    lfn = lfnkey + 8

    IF( dirname /= ' ' ) THEN

      temchar = sfcoutfl
      sfcoutfl = dirname(1:ldirnam)//'/'//temchar
      lfn  = lfn + ldirnam + 1

    END IF

    CALL fnversn(sfcoutfl, lfn)

    IF (sfcdmp > 0) THEN
      PRINT *, 'Write surface property data in ',sfcoutfl(1:lfn)

      CALL wrtsfcdt(nx,ny,nstyps,sfcoutfl(1:lfn), dx,dy,                &
               mapproj,trulat1,trulat2,trulon,sclfct,ctrlat,ctrlon,     &
                  1,1,1,1,1,0,                                          &
                  soiltyp,stypfrct,vegtyp,lai,roufns,veg,veg)

      IF (sfcdmp == 1) CALL sfccntl(nx,ny, sfcoutfl(1:lfn),             &
                 1,1,1,1,1,0, x,y, tem1,tem2)

    END IF

  END IF       ! landin.eq.1

  STOP 0

  100   WRITE(6,'(a)')                                                  &
          'Error reading NAMELIST file. Program ARPSTINTRP stopped.'
  STOP 9104

  9001  CONTINUE

  WRITE(6,'(/2x,a)')'For the output grid:'
  WRITE(6,'(2x,a,f12.4)')                                               &
      'The latitude  of the output grid center, ctrlat=',ctrlat
  WRITE(6,'(2x,a,f12.4/)')                                              &
      'The longitude of the output grid center, ctrlon=',ctrlon
  WRITE(6,'(2x,a/2x,a,2f15.4,a)')                                       &
      'The SW corner (i,j)=(2,2) of the grid is located at ',           &
      '(',xgrdorg,ygrdorg,') of the input grid.'

  STOP 9001

  9002  CONTINUE
  WRITE(6,'(1x,a,i2,/1x,a)')                                            &
      'Data read was unsuccessful. ireturn =', ireturn,                 &
      'Job stopped in ARPSTINTRP.'

  STOP 9002
END PROGRAM arpstintrp