!*==sbemit.spg  processed by SPAG 6.70Rc at 15:33 on 10 Dec 2012
!*------------------   SPAG Configuration Options   --------------------
!*--0323,76 000101,-1 000000102011332010100002000000210211210,136 10 --
!*--1100000011112111000000000000,10,10,10,10,10,10,900,100 200000000 --
!*--000000010000000000000,72,72 73,42,38,33 00011112110000100000000  --
!*----------------------------------------------------------------------
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!$Author: joelevin $
!$Date: 2011-03-24 11:33:26 -0500 (Thu, 24 Mar 2011) $
!$Revision: 11724 $
!$HeadURL: https://eweru-dev1.eweru.ksu.edu/svn/code/weps1/branches/WEPS_F90_update/weps.src/src/lib_erosion/sbemit.for $
!**********************************************************************
!     subroutine sbemit
!**********************************************************************
      subroutine sbemit(ounit,ws,hhr)
      use i_p1werm
      use i_m1sim
      use i_m1flag
      use i_file
      use i_m2geo
      use i_e2erod
      use c_datetime
      use s_caldatw
      implicit none
!*--SBEMIT25
!
!*** Start of declarations rewritten by SPAG
!
! Dummy arguments
!
      real :: hhr,ws
      integer :: ounit
      intent (in) hhr,ounit,ws
!
! Local variables
!
      real :: aegt,aegt10,aegtss,emit10,emitss,emittot,tt
      real,save :: aegt10p,aegtp,aegtssp,tims
      integer,save :: da,initflg,mo,prev_erosion_jday,yr
      integer :: i,j
!
!*** End of declarations rewritten by SPAG
!
 
!     To calc the emissions for each time step of the input wind speed
!     The emissions for EPA are the suspension component
!      with units kg m-2 s-1.
!     To write out a file in the format:
!      12 blank col, yr, mo, day, hr, soucename, emissionrate
!
!     Instructions & logic:
!     To get ntstep period emissions output on erosion days:
!       user sets am0efl = 3 in WEPS configuration screen
!          subroutine openfils creates output file emit.out
!          EROSION calls sbemit to write heading in emit.out file,
!          & sets am0efl to 98, then calls sbemit
!          to print (hourly) Weps emissions on erosion days.
!       or
!       user sets ae0efl (print flg)=4 in stand_alone input file
!           EROSION opens emit.out file, calls sbemit to write headings
!           & sets  ae0efl to 99, then calls sbemit
!            to print period emissions for an erosion day.
!
!
!     +++ ARGUMENT DECLARATIONS +++
!     +++ ARGUMENT DEFINITIONS +++
!
!
!     + + + GLOBAL COMMON BLOCKS + + +
!
!
!     + + + LOCAL COMMON BLOCKS + + +
!
!
!     +++ PARAMETERS +++
!
!     +++ LOCAL VARIABLES +++
 
 !     real    hr
 !     save    hr
 
 
!
!     +++ LOCAL VARIABLE DEFINITIONS +++
!
!
!     +++ OUTPUT FORMATS +++
!
 1000 format (1x,'  yr  mo  day     hr  ws  emission (kg m-2 s-1)')
                             !Unit number for detail grid erosion
 1100 format (22x,'        total    salt/creep    susp      PM10')
 1200 format (1x,3(i4),f7.3,f6.2,1x,4(f11.8))
!
!     +++ END SPECIFICATIONS +++
!
!     set initial conditions
 
      if (initflg==0) then
         initflg = initflg + 1
 
         prev_erosion_jday = am0jd - 1   !init to previous day
 
       !  aegtp   = 0.0
       !  aegtssp = 0.0
       !  aegt10p = 0.0
         tims = 3600*24/ntstep  !seconds in each emission period
         call caldatw(da,mo,yr)     !Set day, month and year
         write (0,*) 'First ntstep is: ',ntstep,tims,tims/3600
!          hr = 0
 
         write (ounit,*) 'SBEMIT output'
!          write (ounit,*) 'Suspended emissions < 0.10 mm dia.'
         write (ounit,*)
 
          ! Print date of Run
 1300    format (1x,'Date of run: ',a3,' ',i2.2,', ',i4,' ',&i2.2,':',  &
               & i2.2,':',i2.2)
         write (ounit,1300) mstring,dt(3),dt(1),dt(5),dt(6),dt(7)
         write (ounit,*)
 
         write (ounit,1000)
         write (ounit,1100)
         write (ounit,*)
      end if
    ! else
 
         !init prev erosion hr values to zero if this is new erosion day
      if (prev_erosion_jday/=am0jd) then
         prev_erosion_jday = am0jd
         aegtp = 0.0
         aegtssp = 0.0
         aegt10p = 0.0
      end if
 
      write (0,*) 'Subsequent ntstep is: ',ntstep,tims,tims/3600
!          if (hr .ge. 24) then
!             hr = tims/3600
      call caldatw(da,mo,yr)
!          else
!             hr = hr + tims/3600
!          endif
!     calculate averages of inner grid points
      aegt = 0.0
      aegtss = 0.0
      aegt10 = 0.0
!
      do j = 1,jmax - 1
         do i = 1,imax - 1
            aegt = aegt + egt(i,j)
            aegtss = aegtss + egtss(i,j)
            aegt10 = aegt10 + egt10(i,j)
         end do
      end do
      tt = (imax-1)*(jmax-1)
      aegt = -aegt/tt               ! change signs to positive=emission
      aegtss = -aegtss/tt
      aegt10 = -aegt10/tt
 
!  Commented out so that emission results don't get summed from hr to hr
!  when erosion stops.  - LEW 5/26/06
!     Hourly (or peroid) emission rate (kg m-2 s-1)
   !      if (aegtssp .gt. aegtss) then  !main set egt arrays to zero
   !          aegtp   = 0.0
   !          aegtssp = 0.0
   !          aegt10p = 0.0
   !      endif
!
      emittot = (aegt-aegtp)/tims
      emitss = (aegtss-aegtssp)/tims
      emit10 = (aegt10-aegt10p)/tims
!
!     Save prior hour average emission
      aegtp = aegt
      aegtssp = aegtss
      aegt10p = aegt10
!
!     Write to emit.out file
      write (ounit,1200) yr,mo,da,hhr,ws,&emittot,emittot - emitss,     &
                       & emitss,emit10
!
   !   endif
      end subroutine sbemit