         program stats_loc
!
!******** compute statistics from an I/O API file
!         compute max of 8hr moving average
!         output goes to I/O API file
!         modified from stats.F by George Pouliot
!         by Dongming Hwang

      implicit none
      
      include 'parms3.ext'      ! I/O API constants
      include 'fdesc3.ext'      ! I/O API file description data structure
      include 'iodecl3.ext'     ! I/O API function declarations 

      INTEGER, EXTERNAL  :: ENVINT ! ENVINT
      INTEGER      STATUS          ! ENVINT status

       integer               ::  maxwords
       integer               ::  maxwordlength, maxlinelength
       integer               ::  nwords

       parameter (maxwords = 100,maxlinelength=100,maxwordlength=16)

       character(len=maxwordlength) words(maxwords)
       character(len=maxlinelength) line
       character(len=maxwordlength) keyword

       character(len=maxwordlength) , allocatable :: spec_names(:)
       character(len=maxwordlength) , allocatable :: species(:)

       character(len=maxwordlength) , allocatable :: statschars(:)  

       character(len=maxlinelength)     :: species2op, default
       character(len=maxlinelength)     :: stats_list
       integer, allocatable  :: spec_lookup(:)
       integer, allocatable  :: stats_lookup(:)

       integer logdev

       character(len=8)     :: buf1, buf2, buf3
       character(len=8)     :: buf4, buf5
       integer              :: osteps, o_day
!
!******* file header information for input & output file
!     
       integer                   :: sdate, stime,tstep,mxrec,nvars
       integer                   :: ncols, nrows,nlays, nsteps
       integer                   :: nthik, gdtype 
       character(len=16)         :: gdname 
!
!******* file header information for input file
! 
       integer            vtypei( MXVARS3 ) ! variable type:  M3(INT|REAL|DBLE)
       character(len=16)  vnamei( MXVARS3 ) ! variable names (length MXDLEN3=80)
       character(len=16)  unitsi( MXVARS3 ) !   "   units or 'none' (MXDLEN3=80)
       character(len=80)  vdesci( MXVARS3 ) !   "      descriptions (MXDLEN3=80)

!
!******* file header information for output file
!     
       integer            vtypeo( MXVARS3 ) ! variable type:  M3(INT|REAL|DBLE)
       character(len=16)  vnameo( MXVARS3 ) ! variable names (length MXDLEN3=80)
       character(len=16)  unitso( MXVARS3 ) !   "   units or 'none' (MXDLEN3=80)
       character(len=80)  vdesco( MXVARS3 ) !   "      descriptions (MXDLEN3=80)
!
!******* environment variables
!
        integer                   :: infil_id
        character(2)  infil_ch
        character(len=16) input, output, this_program,sainput,saoutput
        character(len=16) species_list, stats_env

        character(len=60) zone1, zone2, zone3, zone4
!
!******* data storage
!        
        real, allocatable ::        data( : , : , : )
        real, allocatable ::    data_buf( : , : , : , : )

        real, allocatable ::      sa_data( : , : , : )
        real, allocatable ::  sa_data_buf( : , : , : , : )

        real, allocatable :: eighthr_max( : , : , : )
        real, allocatable :: eighthr_sum( : , : , : )
        real, allocatable :: eighthr_dat( : , : , : , : )

        real, allocatable :: sa_eighthr_max( : , : , : )
        real, allocatable :: sa_eighthr_sum( : , : , : )
        real, allocatable :: sa_eighthr_dat( : , : , : , : )

        integer, allocatable ::  zone_col( : , : )
        integer, allocatable ::  zone_row( : , : )

        integer                   :: zone_ptr, hr_id
        integer                   :: num_zn(4), i
        integer                   :: z1, z2, z3, z4
        integer                   :: jdate, jtime
        integer                   ::  tdate, ttime
        integer                   :: num_spe, num_stats

        integer                   :: rdate, rtime, wdate, wtime
!
!******* indexing
!
        integer                   :: oindex, iindex
!
!******* loop variables
!
        integer                   :: v_loop, t_loop, i_loop, j_loop
        integer                   :: l_loop
        integer                   :: k_loop, s_loop
!
!******* error handling
!
        integer istatus
        character(len=80) mesg
!
!******** external functions
!
        integer trimlen
        external trimlen
!
!  read input file name ...
!
        input = 'CMAQ_FILE'
        infil_id = 1
        write(infil_ch, '(I2.2)') infil_id
        write(*,*)'infil_ch ', infil_ch, 'infil_id ', infil_id
        input(10:11) = infil_ch
        write(*,*)"input ", input

        sainput = 'SA_FILE'
!        infil_id = 1
!        write(infil_ch, '(I2.2)') infil_id
!        write(*,*)'infil_ch ', infil_ch, 'infil_id ', infil_id
        sainput(8:9) = infil_ch
        write(*,*)"sainput ", sainput

        output = 'STATS_FILE'
        saoutput = 'SA_STATS_FILE'
        this_program = 'IOAPI_STATS'
        species_list = 'SPECIES'
        stats_env = 'STATS'

        logdev = init3()	!  initialization returns unit # for log
!
!******** open input file and get header information
!
        if ( .not. open3( input, FSREAD3, this_program ) ) 
     &     then
              MESG = 'Could not open file "' //
     &           input( 1: TRIMLEN( input) ) //
     &            '" for input'
              CALL M3EXIT( this_program, 0, 0, MESG, 2 )
        end if

        if ( .not. desc3(input)) then
           mesg = 'Could not get description info for file "' //
     &           input( 1: TRIMLEN( input) ) //'"'
           call m3exit( this_program, 0, 0, mesg, 2 )     
        end if

        ncols = NCOLS3D
        nrows = NROWS3D
        nvars = NVARS3D
        nsteps = MXREC3D
        sdate = SDATE3D
        stime = STIME3D
        tstep = TSTEP3D
	mxrec  = MXREC3D
	nlays  = NLAYS3D
	nthik  = NTHIK3D
	gdtype = GDTYP3D      
        gdname = GDNAM3D

!*!*!*! Open SA file

        if ( .not. open3( sainput, FSREAD3, this_program ) ) 
     &     then
              MESG = 'Could not open file "' //
     &           input( 1: TRIMLEN( sainput) ) //
     &            '" for input'
              CALL M3EXIT( this_program, 0, 0, MESG, 2 )
        end if

        if ( .not. desc3(sainput)) then
           mesg = 'Could not get description info for file "' //
     &           input( 1: TRIMLEN( sainput) ) //'"'
           call m3exit( this_program, 0, 0, mesg, 2 )     
        end if

!
!******* store variable information from input; construct default variable list
!
        default(1:1) = ' '
	do v_loop = 1,nvars
	vtypei(v_loop) = VTYPE3D(v_loop)
	vnamei(v_loop) = VNAME3D(v_loop)
        unitsi(v_loop) = UNITS3D(v_loop)
	vdesci(v_loop) = VDESC3D(v_loop)
        if (v_loop .eq. 1) then
            default = 
     &         vnamei(v_loop)(1:trimlen(vnamei(v_loop)))
        else
            default = 
     &         default(1:trimlen(default))//
     &         vnamei(v_loop)(1:trimlen(vnamei(v_loop)))
        endif
	      enddo

!
!***** get variable list from environment (use default list if not set by user)
!
        call envstr(species_list,'species to output',                
     &      default,species2op,istatus)

!
!***** parse the variable list
!
        call parseline(maxlinelength,maxwords, nwords,maxwordlength,
     &      species2op,words,keyword)     	 
         
	num_spe = nwords
        write (logdev,*)                                              
     &      'Number of species to put in output:',num_spe
        write (logdev,*)                                               
     &      'Number of species to put in output:',num_spe

!
!********* allocate storage for variable list for output
!
	if (num_spe .gt. 0) then	    	 
	   allocate (species(num_spe))
           allocate (spec_lookup(num_spe))
	   do i_loop = 1,num_spe
              species(i_loop) =                                        
     &           words(i_loop)(1:trimlen(words(i_loop)))
           enddo
	endif

!
!********** create cross reference indexing for output variable list 
!********** from input variable list
!
        do i_loop = 1,num_spe
           do v_loop = 1,nvars
              if (vnamei(v_loop) .eq. species(i_loop)) then
                 spec_lookup(i_loop) = v_loop                    
              endif
           enddo
        enddo

!
!******** get list of statistics from environtment 
!******** (use all stats if not supplied by user)
!
        call envstr(stats_env,'stats to output',                     
     &      '8HRMX',stats_list,istatus)
!
!******** parse the stats character string
!           
        call parseline(maxlinelength,maxwords, nwords,maxwordlength,
     &     stats_list,words,keyword)     	 
         
	num_stats = nwords
        write (logdev,*)                                            
     &        'Number of stats to put in output:',num_stats
        write (logdev,*)                                          
     &        'Number of stats to put in output:',num_stats

!
!******** allocate storage to store the list of stats to compute
!
	if (num_stats .gt. 0) then	    	 
	   allocate (statschars(num_stats))
           allocate (stats_lookup(num_stats))
	   do i_loop = 1,num_stats
              statschars(i_loop) =                              
     &           words(i_loop)(1:trimlen(words(i_loop)))
	   enddo
	endif

!
!********* create cross-reference index for stats list
!
        do i_loop = 1,num_stats
           if (    statschars(i_loop) .eq. 'HRMAX') then
              stats_lookup(i_loop) = 1
           elseif (statschars(i_loop) .eq. '8HRMX') then
              stats_lookup(i_loop) = 2
           elseif (statschars(i_loop) .eq. 'HRAVG') then
              stats_lookup(i_loop) = 5
           elseif (statschars(i_loop) .eq. 'DYSUM') then
              stats_lookup(i_loop) = 6
           endif
        enddo

!
!********** setup of header information for output file
!
        NCOLS3D = ncols
        NROWS3D = nrows
        NVARS3D = num_stats*num_spe
c...        MXREC3D = (nsteps-offset)/daylength
        TSTEP3D = tstep
        NLAYS3D = nlays
        GDTYP3D = gdtype
        GDNAM3D = gdname
        NTHIK3D = nthik

!
!******** set up output variable name descriptions and units
!         
        do v_loop = 1,num_spe   ! loop over variable names
           oindex = num_stats*(v_loop-1)
           iindex = spec_lookup(v_loop)

           do j_loop = 1,num_stats   ! loop for each type of statistic
              VNAME3D(oindex+j_loop) = vnamei(iindex)(1:
     &          trimlen(vnamei(iindex)))//statschars(j_loop)
              VTYPE3D(num_stats*(v_loop-1)+j_loop) = M3REAL
              UNITS3D(num_stats*(v_loop-1)+j_loop) = unitsi(iindex)
              VDESC3D(num_stats*(v_loop-1)+j_loop) = vdesci(iindex)
           enddo
        enddo

        do i_loop = 1, num_spe*num_stats
           vnameo(i_loop) = VNAME3D(i_loop)
           vtypeo(i_loop) = VTYPE3D(i_loop)
           unitso(i_loop) = UNITS3D(i_loop)
           vdesco(i_loop) = VDESC3D(i_loop)
        enddo
!
!******** open output file
!
        TSTEP3D = 240000
        if ( .not. open3( output, FSCREA3, this_program ) ) 
     &     then
              MESG = 'Could not open file "' //
     &            output( 1: TRIMLEN( output) ) //
     &            '" for output'
              CALL M3EXIT( this_program, 0, 0, MESG, 2 )
        end if

!*!*!*!*! Open SA output file

        if ( .not. open3( saoutput, FSCREA3, this_program ) ) 
     &     then
              MESG = 'Could not open file "' //
     &            output( 1: TRIMLEN( saoutput) ) //
     &            '" for output'
              CALL M3EXIT( this_program, 0, 0, MESG, 2 )
        end if
           
!
!******** allocate storage for data
!           
       allocate          (data(ncols,nrows,nlays))
       allocate       (sa_data(ncols,nrows,nlays))

       allocate   (eighthr_dat(ncols,nrows,nlays,8))
       allocate   (eighthr_max(ncols,nrows,nlays))
       allocate   (eighthr_sum(ncols,nrows,nlays))

       allocate   (sa_eighthr_dat(ncols,nrows,nlays,8))
       allocate   (sa_eighthr_max(ncols,nrows,nlays))
       allocate   (sa_eighthr_sum(ncols,nrows,nlays))

       allocate   (data_buf(ncols,nrows,nlays,4))
       allocate   (sa_data_buf(ncols,nrows,nlays,4))
       allocate   (zone_col(4,ncols*nrows))
       allocate   (zone_row(4,ncols*nrows))

c input zone grid

c   inpu zone definition
       read(10,*) buf1, buf2, buf3, buf4, buf5
       write(*,*)' line ', buf1, buf2, buf3, buf4, buf5
       do i = 1, 4
          num_zn(i) = 0
       enddo
       z1 = 1
       z2 = 1
       z3 = 1
       z4 = 1
       do i = 1, ncols*nrows
           read(10,*,end=112) i_loop, j_loop, k_loop, buf1, t_loop
           if ( t_loop .eq. -8 ) then
              zone_col(4, z4) = i_loop
              zone_row(4, z4) = j_loop
              z4 = z4+1
              num_zn(4) = num_zn(4)+1
           else if ( t_loop .eq. -7 ) then
              zone_col(3, z3) = i_loop
              zone_row(3, z3) = j_loop
              z3 = z3+1
              num_zn(3) = num_zn(3)+1
           else if ( t_loop .eq. -6 ) then
              zone_col(2, z2) = i_loop
              zone_row(2, z2) = j_loop
              z2 = z2+1
              num_zn(2) = num_zn(2)+1
           else if ( t_loop .eq. -5 ) then
              zone_col(1, z1) = i_loop
              zone_row(1, z1) = j_loop
              z1 = z1+1
              num_zn(1) = num_zn(1)+1
           end if
        enddo

112     continue

        write(*,*)'zone summary: Eastern:Centeral:Mountain:Pacific',
     *           (num_zn(i),i=1,4)
        write(*,*)' z1 ', num_zn(1), 'col row', 
     *          zone_col(1,num_zn(1)), zone_row(1,num_zn(1))
        write(*,*)' z2 ', num_zn(2), 'col row', 
     *          zone_col(2,num_zn(2)), zone_row(2,num_zn(2))
        write(*,*)' z3 ', num_zn(3), 'col row', 
     *          zone_col(3,num_zn(3)), zone_row(3,num_zn(3))
        write(*,*)' z4 ', num_zn(4), 'col row', 
     *          zone_col(4,num_zn(4)), zone_row(4,num_zn(4))


        o_day = ENVINT( 'OUT_DAY', ' ', 0, STATUS )
        osteps = o_day*24
        write(*,*)'total output # of day:hour ', o_day, osteps

c loop on species

        do v_loop = 1,num_spe

c  get jdate/jtime : output date/time
           jdate = sdate
           jtime = stime
c  get rdate/rtime : reading date/time
           rdate = jdate
           rtime = stime + 50000     ! 00 eastern time = GMT + 5
           wdate = jdate
           wtime = 0                 ! local time

c-----------------------------------------------------------
c      fetch 7 hours data
c      first hour data: pre-fetch 4 hours data
        hr_id = 1
        do t_loop = 1, 4
           if ( .not. 
     &        read3(input,vnamei(spec_lookup(v_loop)),
     &        ALLAYS3,rdate,rtime,data(1,1,1))
     &        ) then

              infil_id = infil_id+1
              write(infil_ch, '(I2.2)') infil_id
              input(10:11) = infil_ch
              sainput(8:9) = infil_ch

              if ( .not. open3( input, FSREAD3, this_program ) ) 
     &           then
                 MESG = 'Could not open file "' //
     &              input( 1: TRIMLEN( input) ) //
     &              '" for input'
                    CALL M3EXIT( this_program, 0, 0, MESG, 2 )
                 end if

                 if ( .not. 
     &                read3(input,vnamei(spec_lookup(v_loop)),
     &                ALLAYS3,rdate,rtime,data(1,1,1))
     &                ) then
                    mesg = 'Error reading '//
     &                  vnamei(spec_lookup(v_loop))//
     &                 'from file '//
     &                 input( 1: TRIMLEN( input ) )
                    call m3exit( this_program, 0, 0, mesg, 2 )
                 end if                                    

               end if
               do i_loop = 1, ncols
               do j_loop = 1, nrows
               do k_loop = 1, nlays
                    data_buf(i_loop, j_loop, k_loop, t_loop) = 
     &                    data(i_loop, j_loop, k_loop)
               enddo
               enddo
               enddo

!*!*!*!*!*! SA 

           if ( .not. 
     &        read3(sainput,vnamei(spec_lookup(v_loop)),
     &        ALLAYS3,rdate,rtime,sa_data(1,1,1))
     &        ) then

!              infil_id = infil_id+1
              write(infil_ch, '(I2.2)') infil_id
              input(10:11) = infil_ch
              sainput(8:9) = infil_ch

              if ( .not. open3( sainput, FSREAD3, this_program ) ) 
     &           then
                 MESG = 'Could not open file "' //
     &              input( 1: TRIMLEN( sainput) ) //
     &              '" for input'
                    CALL M3EXIT( this_program, 0, 0, MESG, 2 )
                 end if

                 if ( .not. 
     &                read3(sainput,vnamei(spec_lookup(v_loop)),
     &                ALLAYS3,rdate,rtime,sa_data(1,1,1))
     &                ) then
                    mesg = 'Error reading '//
     &                  vnamei(spec_lookup(v_loop))//
     &                 'from file '//
     &                 input( 1: TRIMLEN( sainput ) )
                    call m3exit( this_program, 0, 0, mesg, 2 )
                 end if                                    

               end if
               do i_loop = 1, ncols
               do j_loop = 1, nrows
               do k_loop = 1, nlays
                    sa_data_buf(i_loop, j_loop, k_loop, t_loop) = 
     &                    sa_data(i_loop, j_loop, k_loop)
               enddo
               enddo
               enddo
!*!*!*!

               call nextime(rdate,rtime,tstep)
           end do   ! end t_loop for first hour data
           write(*,*)' complete first four hour data fetch ...'

           zone_ptr = 1
           call zone_merge(eighthr_dat, data_buf, 
     &                     ncols, nrows, nlays,
     *                     zone_ptr, 
     *                     zone_col, zone_row, num_zn, hr_id)

           call sa_zone_merge(sa_eighthr_dat, sa_data_buf, 
     &                     ncols, nrows, nlays,
     *                     zone_ptr, 
     *                     zone_col, zone_row, num_zn, hr_id)

             write(*,*)' rdate:rtime ------------- ',rdate, rtime


c        get next 6 hour data
        do l_loop = 2, 7

           if ( .not. 
     &          read3(input,vnamei(spec_lookup(v_loop)),
     &          ALLAYS3,rdate,rtime,data(1,1,1))
     &          ) then

c    ...  end of this day, continue next day ...
              infil_id = infil_id+1
              write(infil_ch, '(I2.2)') infil_id
              input(10:11) = infil_ch
              sainput(8:9) = infil_ch
              if ( .not. open3( input, FSREAD3, this_program ) ) 
     &           then
                 MESG = 'Could not open file "' //
     &              input( 1: TRIMLEN( input) ) //
     &              '" for input'
                 CALL M3EXIT( this_program, 0, 0, MESG, 2 )
              end if

              if ( .not. 
     &             read3(input,vnamei(spec_lookup(v_loop)),
     &             ALLAYS3,rdate,rtime,data(1,1,1))
     &             ) then
                 mesg = 'Error reading '//
     &               vnamei(spec_lookup(v_loop))//
     &              'from file '//
     &              input( 1: TRIMLEN( input ) )
                 call m3exit( this_program, 0, 0, mesg, 2 )
              end if                                    
           end if               

           do i_loop = 1, ncols
           do j_loop = 1, nrows
           do k_loop = 1, nlays
                data_buf(i_loop, j_loop, k_loop, zone_ptr) = 
     &                data(i_loop, j_loop, k_loop)
           enddo
           enddo
           enddo

!*!*!*!*! SA part

           if ( .not. 
     &          read3(sainput,vnamei(spec_lookup(v_loop)),
     &          ALLAYS3,rdate,rtime,sa_data(1,1,1))
     &          ) then

c    ...  end of this day, continue next day ...
!              infil_id = infil_id+1
              write(infil_ch, '(I2.2)') infil_id
              input(10:11) = infil_ch
              sainput(8:9) = infil_ch
              if ( .not. open3( sainput, FSREAD3, this_program ) ) 
     &           then
                 MESG = 'Could not open file "' //
     &              input( 1: TRIMLEN( sainput) ) //
     &              '" for input'
                 CALL M3EXIT( this_program, 0, 0, MESG, 2 )
              end if

              if ( .not. 
     &             read3(sainput,vnamei(spec_lookup(v_loop)),
     &             ALLAYS3,rdate,rtime,sa_data(1,1,1))
     &             ) then
                 mesg = 'Error reading '//
     &               vnamei(spec_lookup(v_loop))//
     &              'from file '//
     &              input( 1: TRIMLEN( sainput ) )
                 call m3exit( this_program, 0, 0, mesg, 2 )
              end if                                    
           end if               

           do i_loop = 1, ncols
           do j_loop = 1, nrows
           do k_loop = 1, nlays
                sa_data_buf(i_loop, j_loop, k_loop, zone_ptr) = 
     &                sa_data(i_loop, j_loop, k_loop)
           enddo
           enddo
           enddo

!*!*!*!

           call nextime(rdate,rtime,tstep)

c  update ptr
           zone_ptr = mod(zone_ptr+1, 4)
           if (zone_ptr .eq. 0) zone_ptr = 4

           hr_id = hr_id + 1
           if (hr_id .gt. 8) hr_id = mod(hr_id,8)
              call zone_merge(eighthr_dat, data_buf, 
     &                        ncols, nrows, nlays,
     *                        zone_ptr, 
     *                        zone_col, zone_row, num_zn, hr_id)

              call sa_zone_merge(sa_eighthr_dat, sa_data_buf, 
     &                        ncols, nrows, nlays,
     *                        zone_ptr, 
     *                        zone_col, zone_row, num_zn, hr_id)

        end do    ! end 7 hour fetch

        write(*,*)'complete 7 hour fetch, begin 8hr processing ...'

c-----------------------------------------------------------
c begin 8 hr process
        do t_loop = 1,osteps
	   write (logdev,*) t_loop, jdate,jtime
           if ( .not. 
     &          read3(input,vnamei(spec_lookup(v_loop)),
     &          ALLAYS3,rdate,rtime,data(1,1,1))
     &          ) then
c    ...  end of this day, continue next day ...
              infil_id = infil_id+1
              write(infil_ch, '(I2.2)') infil_id
              input(10:11) = infil_ch
              if ( .not. open3( input, FSREAD3, this_program ) ) 
     &            then
                  MESG = 'Could not open file "' //
     &               input( 1: TRIMLEN( input) ) //
     &               '" for input'
                  CALL M3EXIT( this_program, 0, 0, MESG, 2 )
              end if

              if ( .not. 
     &             read3(input,vnamei(spec_lookup(v_loop)),
     &             ALLAYS3,rdate,rtime,data(1,1,1))
     &             ) then
                 mesg = 'Error reading '//
     &                vnamei(spec_lookup(v_loop))//
     &               'from file '//
     &               input( 1: TRIMLEN( input ) )
                 call m3exit( this_program, 0, 0, mesg, 2 )
              end if                                    
            end if                                    

            do i_loop = 1, ncols
            do j_loop = 1, nrows
            do k_loop = 1, nlays
                 data_buf(i_loop, j_loop, k_loop, zone_ptr) = 
     &                 data(i_loop, j_loop, k_loop)
            enddo
            enddo
            enddo

!*!*!*!*!  SA part

           if ( .not. 
     &          read3(sainput,vnamei(spec_lookup(v_loop)),
     &          ALLAYS3,rdate,rtime,sa_data(1,1,1))
     &          ) then
c    ...  end of this day, continue next day ...
!              infil_id = infil_id+1
              write(infil_ch, '(I2.2)') infil_id
              input(10:11) = infil_ch
              sainput(8:9) = infil_ch
              if ( .not. open3( sainput, FSREAD3, this_program ) ) 
     &            then
                  MESG = 'Could not open file "' //
     &               input( 1: TRIMLEN( sainput) ) //
     &               '" for input'
                  CALL M3EXIT( this_program, 0, 0, MESG, 2 )
              end if

              if ( .not. 
     &             read3(sainput,vnamei(spec_lookup(v_loop)),
     &             ALLAYS3,rdate,rtime,sa_data(1,1,1))
     &             ) then
                 mesg = 'Error reading '//
     &                vnamei(spec_lookup(v_loop))//
     &               'from file '//
     &               input( 1: TRIMLEN( sainput ) )
                 call m3exit( this_program, 0, 0, mesg, 2 )
              end if                                    
            end if                                    

            do i_loop = 1, ncols
            do j_loop = 1, nrows
            do k_loop = 1, nlays
                 sa_data_buf(i_loop, j_loop, k_loop, zone_ptr) = 
     &                 sa_data(i_loop, j_loop, k_loop)
            enddo
            enddo
            enddo
!*!*!*!*!

            call nextime(rdate,rtime,tstep)

            zone_ptr = mod(zone_ptr+1, 4)
            if (zone_ptr .eq. 0) zone_ptr = 4
            hr_id = hr_id + 1
            if (hr_id .gt. 8) hr_id = mod(hr_id,8)
            call zone_merge(eighthr_dat, data_buf, 
     &                      ncols, nrows, nlays,
     *                      zone_ptr, 
     *                      zone_col, zone_row, num_zn, hr_id)

            call sa_zone_merge(sa_eighthr_dat, sa_data_buf, 
     &                      ncols, nrows, nlays,
     *                      zone_ptr, 
     *                      zone_col, zone_row, num_zn, hr_id)


c 8 hr max
            eighthr_sum(1:ncols,1:nrows,1:nlays) = 0.0
            sa_eighthr_sum(1:ncols,1:nrows,1:nlays) = 0.0
            do l_loop = 1,8
            do k_loop = 1,nlays
            do j_loop = 1,nrows
            do i_loop = 1,ncols
                eighthr_sum(i_loop, j_loop, k_loop) = 
     &             eighthr_sum(i_loop, j_loop, k_loop) +
     &             eighthr_dat(i_loop, j_loop, k_loop,l_loop) 
                sa_eighthr_sum(i_loop, j_loop, k_loop) = 
     &             sa_eighthr_sum(i_loop, j_loop, k_loop) +
     &             sa_eighthr_dat(i_loop, j_loop, k_loop,l_loop) 
            enddo
            enddo
            enddo
            enddo

            do k_loop = 1, nlays
            do j_loop = 1,nrows
            do i_loop = 1,ncols
                eighthr_sum(i_loop, j_loop, k_loop) = 
     &             eighthr_sum(i_loop, j_loop, k_loop)/8.0
                sa_eighthr_sum(i_loop, j_loop, k_loop) = 
     &             sa_eighthr_sum(i_loop, j_loop, k_loop)/8.0
            enddo
            enddo
            enddo       

c ...              starting of the day
            if (jtime .eq. 0) then
                eighthr_max(1:ncols,1:nrows,1:nlays) 
     &              = eighthr_sum(1:ncols, 1:nrows, 1:nlays)
                sa_eighthr_max(1:ncols,1:nrows,1:nlays) 
     &              = sa_eighthr_sum(1:ncols, 1:nrows, 1:nlays)
            endif

            do k_loop = 1, nlays
            do j_loop = 1,nrows
            do i_loop = 1,ncols
                if (eighthr_sum(i_loop, j_loop, k_loop) .gt.
     &              eighthr_max(i_loop, j_loop, k_loop)) then
                       eighthr_max(i_loop, j_loop, k_loop) =
     &                    eighthr_sum(i_loop, j_loop, k_loop)

!! keep SA 8hr max at same time as regular model output time 8hr max
                     sa_eighthr_max(i_loop, j_loop, k_loop) =
     &                  sa_eighthr_sum(i_loop, j_loop, k_loop)

                endif
            enddo
            enddo
            enddo       

            write(*,*)' rdate:rtime output ',rdate, rtime 

            if (mod(t_loop,24) .eq. 0) then
                if ( .not. 
     &             write3(output,
     &                vnameo(num_stats*(v_loop-1)+1),
     &                wdate,wtime,
     &                eighthr_max(1,1,1))) then
                      mesg = 'Error writing 8 hourly_max to file '//
     &                   output( 1: TRIMLEN( output ) ) 
                         call m3exit( this_program, 0, 0, mesg, 2 )
                end if

                if ( .not. 
     &             write3(saoutput,
     &                vnameo(num_stats*(v_loop-1)+1),
     &                wdate,wtime,
     &                sa_eighthr_max(1,1,1))) then
                      mesg = 'Error writing 8 hourly_max to file '//
     &                   output( 1: TRIMLEN( saoutput ) ) 
                         call m3exit( this_program, 0, 0, mesg, 2 )
                end if

                call nextime(wdate, wtime, 24*tstep)   ! next time step
            end if
            call nextime(jdate, jtime, tstep)   ! next time step

        enddo
        enddo  
!
!******** must shut down IO/API before exiting (necessary when writing to IO/API files)
!

        if ( .not. shut3() ) then
	        write (*,*) 'Error shutting down IOAPI'
	        stop 'ERROR: Abnormal Termination'
	endif

        write (logdev,*) this_program,': Normal Completion'

        stop
990     end

        subroutine updatemax(ncols,nrows, nlays,
     &             maxvals,data)

        implicit none
        integer ncols, nrows, nlays
        real maxvals(ncols,nrows,nlays)
        real data(ncols, nrows, nlays)

        integer i_loop, j_loop, k_loop

                do k_loop = 1, nlays
                   do j_loop = 1, nrows
                      do i_loop = 1, ncols
                         if (data(i_loop, j_loop, k_loop) .gt.
     &                       maxvals(i_loop,j_loop, k_loop)) then
                                  maxvals(i_loop, j_loop, k_loop)=
     &                            data(i_loop, j_loop, k_loop)
                         endif
                      enddo
                   enddo
               enddo


        return
        end
      
       subroutine parseline(maxlinelength,maxwords, nwords,
     &      maxwordlength,
     &      line,words,keyword)
!       
!assumed delimiter between words is one or more blanks ' '
!       
       integer maxwords, maxwordlength, nwords, maxlinelength
       
       character(len=maxwordlength) words(maxwords)
       character(len=maxlinelength) line
       character(len=maxwordlength) keyword
 
!
!******** locals
!
         integer tlength,lbks, istart, iend, i_loop
         logical inword            !true if in a word, false otherwise
         character*1 char,prevchar
!
!******** external function calls
!         
         external lblank, trimlen
         integer trimlen, lblank
!         
!
!******* take a string of text (line) and split it into words
!        subject to the following limitations:
!        the maximum number of words in the string is maxwords
!        the maximum lenght of any word is maxwordlength
!        the maximum length of line is maxlinelength

         tlength = trimlen(line(1:maxlinelength))
         
         if (tlength .eq. maxlinelength) then
            write (*,*) 'ERROR: MAXLINELENGTH REACHED'
            write (*,*) 'increase the maximum length of the record'
            stop
         endif
         lbks = lblank(line(1:maxlinelength))
 
         nwords = 0
         inword = .true.

         istart = lbks+1
         iend = lbks+1
         prevchar = ' '
         
         do i_loop = lbks+1,tlength+1
            char = line(i_loop:i_loop)
            if ((char .eq. ' ') .and. (prevchar .ne. ' ')) then
               inword = .false.
               nwords = nwords + 1                              
               words(nwords) = line(istart:iend)
            else
                if((char .eq. ' ') .and. (prevchar .eq. ' ')) then
!
!************ do nothing
!
                else                
                   if (inword) then
                      iend = iend + 1
                   else
                      inword = .true.
                      istart = i_loop
                      iend = i_loop
                   endif                      
                endif
            endif
            prevchar = char
         enddo                  
         
         keyword = words(1)
         keyword = keyword(2:trimlen(keyword))
         
         return 
         end	      

c       merge data
        subroutine zone_merge(eighthr_dat, data_buf, 
     &                          ncols, nrows, nlays,
     *                          zone_ptr, 
     *                          zone_col, zone_row, num_zn, hr_id)
        real eighthr_dat(ncols, nrows, nlays, 8)
        real data_buf(ncols, nrows, nlays,4)
        integer   zone_col(4, ncols*nrows)
        integer   zone_row(4, ncols*nrows)
        integer   zone_ptr, hr_id
        integer   num_zn(4)

        integer i_loop, j_loop, k_loop
        integer i, zn
        integer c, r

c       merge data
        do i = 0, 3
           zn = mod(zone_ptr+i-1, 4) + 1

           do k_loop = 1, nlays
              do j_loop = 1, num_zn(i+1)
                    c = zone_col(i+1, j_loop)
                    r = zone_row(i+1, j_loop)
                    eighthr_dat(c, r, k_loop, hr_id) =
     &                  data_buf(c, r, k_loop, zn)
              enddo
           enddo

        enddo
         
        return 
        end	      



        subroutine sa_zone_merge(sa_eighthr_dat, sa_data_buf, 
     &                          ncols, nrows, nlays,
     *                          zone_ptr, 
     *                          zone_col, zone_row, num_zn, hr_id)
        real sa_eighthr_dat(ncols, nrows, nlays, 8)
        real sa_data_buf(ncols, nrows, nlays,4)
        integer   zone_col(4, ncols*nrows)
        integer   zone_row(4, ncols*nrows)
        integer   zone_ptr, hr_id
        integer   num_zn(4)

        integer i_loop, j_loop, k_loop
        integer i, zn
        integer c, r

c       merge data
        do i = 0, 3
           zn = mod(zone_ptr+i-1, 4) + 1

           do k_loop = 1, nlays
              do j_loop = 1, num_zn(i+1)
                    c = zone_col(i+1, j_loop)
                    r = zone_row(i+1, j_loop)
                    sa_eighthr_dat(c, r, k_loop, hr_id) =
     &                  sa_data_buf(c, r, k_loop, zn)
              enddo
           enddo

        enddo
         
        return 
        end	      

