SUBROUTINE phys_rates
!!-----------------------------------------------------------------
! rates of :
! dilution and exchange between 2 boxes (box-specific)
! emission & deposition (box-specific) 
! NB: this subroutine can currently ONLY handle 2 boxes.
! Modification required for > 2 boxes.
!!-----------------------------------------------------------------

      USE flags_module,ONLY: emis_fg,depos_fg
      USE akparameter_module
      USE forcing_params_module
      USE module_data_gecko_main

      IMPLICIT NONE
!------------------------------------------------------
      !PRINT*,"starting physical rates"

      rem  = 0.
      rdep = 0.
      rdil = 0.
      rex  = 0.

!---------------
      SELECT CASE (ibox)
!---------------
        CASE(1)
! box 1 dilution & exchange
          IF (dilfix .eq. 1) THEN
            rdil = dilconst
          ELSE 
            IF (dhdt > 0.) THEN
              rdil =  dhdt/height
              IF(nbox.eq.2)THEN
              rex = conc(:,2)*dhdt/height
              ELSE
              rex = cbg(:)*dhdt/height
              ENDIF
            ENDIF
          ENDIF
! box 1 emission
        IF (emis_fg .GT. 0) rem = eflux/height

! box 1 deposition for gas, aerosol species
        IF (depos_fg .gt. 0) THEN
          DO i=1,ndepspe
            rdep(iddepspe(i)) = vd(i)/height
          ENDDO
          DO i=1, nsat
            rdep(idasat(i)) = inorg_aer(ibox)%vdepaer / height
          ENDDO
        ENDIF

! save conc (box 1, start of timestep) for later if needed by second box
        IF (nbox.GT.1) cbot_sav(:) = conc(:,1)

!---------------
        CASE(2)
! box 2 dilution
          IF (dilfix .eq. 1) THEN
            rdil = dilconst
          ELSE
            IF (dhdt >= 0.) THEN
              rdil =  0.
              rex  =  0.
            ELSE ! (dhdt < 0.) 
! box 2 needs conc(:,1) from the beginning of the timestep!
!             hence use cbot_sav
              rdil = -dhdt/(htop-height)
              rex  = -cbot_sav*dhdt/(htop-height)
            ENDIF
          ENDIF

! box 2 dilution, add the loss due to atmospheric subsidence
! i.e  gases are "pushed away" by subsiding air from the free trop.
          rdil = rdil + vs/(htop - height)

! box 2 emission is the source coming from atmospheric subsidence
! i.e  gases are slowly subsiding from the free trop.
          rem = cbg *vs / (htop - height)

!---------------
      END SELECT
!---------------

! -------------------------------------
      END SUBROUTINE phys_rates
! -------------------------------------