MODULE time_mod
   use prec, only: rstd
   PRIVATE

   INTEGER, SAVE :: ncid
!$OMP THREADPRIVATE(ncid)
   INTEGER, SAVE :: time_counter_id
!$OMP THREADPRIVATE(time_counter_id)
   INTEGER, SAVE :: it
!$OMP THREADPRIVATE(it)

   INTEGER, SAVE :: itau0 = 0
!$OMP THREADPRIVATE(itau0)

   REAL(rstd), SAVE :: dt ! time step (s)
!$OMP THREADPRIVATE(dt)
   REAL(rstd), SAVE :: write_period
!$OMP THREADPRIVATE(write_period)
   INTEGER, SAVE    :: itau_out, itau_adv, itau_dissip, itau_physics, itaumax
!$OMP THREADPRIVATE(itau_out, itau_adv, itau_dissip, itau_physics, itaumax)

   INTEGER, SAVE    :: itau_check_conserv
!$OMP THREADPRIVATE(itau_check_conserv)

   INTEGER, SAVE :: day_step ! number of time steps per day
   REAL(rstd), SAVE :: ndays ! number of days to run
!$OMP THREADPRIVATE(day_step,ndays)

   CHARACTER(LEN=255) :: time_style
!$OMP THREADPRIVATE(time_style)

   PUBLIC create_time_counter_header, update_time_counter, close_time_counter, init_time, &
      dt, write_period, itau_out, itau_adv, itau_dissip, itau_physics, itaumax, &
      itau_check_conserv, &
      day_step, ndays, time_style, itau0

CONTAINS

   SUBROUTINE init_time
      USE earth_const, ONLY: daysec, scale_factor
      USE getin_mod, ONLY: getin
      USE mpipara, ONLY: is_mpi_root
      IMPLICIT NONE
      REAL(rstd) :: run_length

      ! @getin_doc daysec realm=time type=real default=86400
      ! Length of the day, in seconds.
      CALL getin('daysec', daysec)

      ! @getin_doc dt realm=time type=real default=0
      ! Time step for a full dynamics type step. Limited by horizontal resolution and speed of sound.
      ! If `dt==0` then `day_step` must be set to a non-zero value.
      dt = 0
      CALL getin("dt", dt)

      ! @getin_doc day_step realm=time type=integer default=0
      ! Number of full dynamics time steps per day (see `daysec`).
      ! Active if `dt==0`.
      day_step = 0
      CALL getin('day_step', day_step)
      IF (dt == 0) THEN
         ! use day_step instead
         IF (day_step == 0) THEN
            ! Error, neither dt nor day_step are set
            WRITE (*, *) "init_time: Error, neither dt nor day_step are set!"
            STOP
         ELSE
            dt = daysec/day_step
         END IF
      ELSE
         IF (day_step /= 0) THEN
            ! Error, both dt and day_step are set
            WRITE (*, *) "init_time: Error, both dt=", dt, " and ", &
               "day_step=", day_step, " are set!"
            STOP
         ELSE
            day_step = NINT(daysec/dt)
         END IF
      END IF

      ! @getin_doc write_period realm=time type=real default=0
      ! Time in seconds between outputs. For DYNAMICO, 'output' means sending data to XIOS.
      ! Thus XIOS receives from DYNAMICO data sampled at a time resolution
      ! of `write_period`. Any XIOS operation requested via XML configuration files,
      ! such as time averaging/min/max, operates on that data.
      ! Ultimately XIOS writes data at a frequency `output_freq` defined in XML files on a per-field basis.
      ! `write_period` should be a divisor of all XIOS output frequencies (see however `time_style`).
      write_period = 0
      CALL getin('write_period', write_period)

      ! @getin_doc time_style realm=time type=String default=none
      ! Valid values :
      ! * 'none' : normal time style
      ! * 'dcmip' : time parameters are rescaled (divided) by `scale factor` for DCMIP2012 small-planet experiments.
      !   This applies to `dt`, `write_period` and dissipation time scales.
      !   Time values visible in the DYNAMICO log and seen by XIOS are the rescaled values. Especially XML files
      !   must set 'output_freq' in rescaled time, i.e. multiples of `write_period/scale_factor`.
      !   Dividing `dt` by `scale_factor` may result in a non-integer time step which may cause issues with XIOS.
      ! FIXME : It seems `run_length` is not rescaled ?
      time_style = 'none'
      CALL getin('time_style', time_style)

      SELECT CASE (TRIM(time_style))
      CASE ('none') ! do nothing
      CASE ('dcmip') ! rescale time step for small-planet experiments
         dt = dt/scale_factor
         write_period = write_period/scale_factor
         IF (is_mpi_root) PRINT *, 'Output frequency (scaled) set to ', write_period
      CASE DEFAULT
         IF (is_mpi_root) PRINT *, "Bad selector for variable time_style >", TRIM(time_style), "> options are <none>, <dcmip>"
         STOP
      END SELECT

      itau_out = FLOOR(.5 + write_period/dt)
      IF (is_mpi_root) PRINT *, 'Output frequency itau_out = ', itau_out

      ! @getin_doc run_length realm=time type=real default=0
      ! Total length of the run, in seconds. Must be a multiple of `dt`.
      ! Active if `itaumax==0`, in which case either `run_length` or `ndays` must be set to a non-zero value, but not both.
      run_length = 0
      CALL getin('run_length', run_length)

      ! @getin_doc run_length realm=time type=integer default=0
      ! Total length of the run, in days (see `daysec`).
      ! Active if `itaumax==0`, in which case either `run_length` or `ndays` must be set to a non-zero value, but not both.
      ndays = 0
      CALL getin('ndays', ndays)

      ! @getin_doc itaumax realm=time type=integer default=0
      ! Total length of the run in units of dynamics time steps `dt`, i.e. the effect is equivalent
      ! to setting `run_length` to `dt*itaumax`. If set to zero (default), `run_length` or `ndays` apply.
      itaumax = 0
      CALL getin('itaumax', itaumax)

      IF (itaumax == 0) THEN
         ! use run_length or ndays
         IF (run_length == 0) THEN
            !look for ndays instead
            IF (ndays == 0) THEN
               ! Error, neither run_length nor ndays is set
               WRITE (*, *) "init_time: Error, neither itaumx nor run_length nor ndays are set!"
               STOP
            ELSE
               run_length = ndays*daysec
            END IF
         ELSE ! run_length/=0
            IF (ndays /= 0) THEN
               ! Error, both ndays and run_length are set
               WRITE (*, *) "init_time: Error, both ndays=", ndays, " and ", &
                  "run_length=", run_length, " are set!"
               STOP
            ELSE
               ndays = run_length/daysec
            END IF
         END IF
         itaumax = NINT(run_length/dt)
      ELSE ! itaumax/=0
         ! itaumax overrides run_length and ndays
         run_length = dt*itaumax
         ndays = run_length/daysec
      END IF

      itau_adv = 1
      CALL getin('itau_adv', itau_adv)

      itau_dissip = 0 ! set to zero which implies itau_dissip will be automatically computed (see init_dissip)
      CALL getin('itau_dissip', itau_dissip)

      itau_physics = 1
      CALL getin('itau_physics', itau_physics)
      if (itau_physics <= 0) itau_physics = HUGE(itau_physics)

      itau_check_conserv = HUGE(itau_check_conserv)
      CALL getin('itau_check_conserv', itau_check_conserv)

      IF (is_mpi_root) THEN
         WRITE (*, *) 'init_time: run for run_length=', run_length, ' (s)'
         WRITE (*, *) '           i.e. ndays=', ndays, ' (days)'
         WRITE (*, *) '              itaumax=', itaumax, ' (iterations)'
         WRITE (*, *) '           dt=', dt, ' (s)'
         WRITE (*, *) '           itau_adv=', itau_adv
         WRITE (*, *) '           itau_dissip=', itau_dissip
         WRITE (*, *) '           itau_physics=', itau_physics
      END IF

      CALL create_time_counter_header

   END SUBROUTINE init_time

   SUBROUTINE create_time_counter_header
      USE netcdf_mod
      USE prec
      USE getin_mod
      USE mpipara
      IMPLICIT NONE
      INTEGER :: status
      INTEGER :: timeid, dtid
      REAL(rstd) :: dt
      CHARACTER(LEN=255) :: time_frequency

      IF (no_io) RETURN

      dt = 0
      CALL getin("dt", dt)

!$OMP BARRIER
!$OMP MASTER
      IF (is_mpi_root) THEN
         status = NF90_CREATE('time_counter.nc', NF90_CLOBBER, ncid)
         status = NF90_DEF_DIM(ncid, 'time_counter', NF90_UNLIMITED, timeid)
         status = NF90_DEF_VAR(ncid, 'time_counter', NF90_DOUBLE, (/timeid/), time_counter_id)
         status = NF90_PUT_ATT(ncid, time_counter_id, "long_name", "time")
         status = NF90_PUT_ATT(ncid, time_counter_id, "units", "seconds since 2000-01-01 00:00:00")
         status = NF90_PUT_ATT(ncid, time_counter_id, "calendar", "noleap")
         status = NF90_DEF_VAR(ncid, 'mdt', NF90_DOUBLE, varid=dtid)

         WRITE (time_frequency, *) write_period
         PRINT *, TRIM(time_frequency)
         status = NF90_PUT_ATT(ncid, NF90_GLOBAL, "time_frequency", TRIM(time_frequency)//"s")

         status = NF90_ENDDEF(ncid)

         status = NF90_PUT_VAR(ncid, dtid, dt)
      END IF
      it = 0
!$OMP END MASTER
!$OMP BARRIER

   END SUBROUTINE create_time_counter_header

   SUBROUTINE update_time_counter(time)
      USE netcdf_mod
      USE mpipara
      USE prec
      IMPLICIT NONE
      REAL(rstd), INTENT(IN) ::time
      INTEGER :: status
      REAL(rstd) ::time_array(1)

      IF (no_io) RETURN

!$OMP BARRIER
!$OMP MASTER
      time_array(1) = time

      it = it + 1
      IF (is_mpi_root) THEN
         status = NF90_PUT_VAR(ncid, time_counter_id, time_array, start=(/it/), count=(/1/))
         status = NF90_SYNC(ncid)
      END IF
!$OMP END MASTER
!$OMP BARRIER

   END SUBROUTINE update_time_counter

   SUBROUTINE close_time_counter
      USE netcdf_mod
      USE mpipara
      IMPLICIT NONE
      INTEGER :: status

      IF (no_io) RETURN

!$OMP BARRIER
!$OMP MASTER
      IF (is_mpi_root) status = NF90_CLOSE(ncid)
!$OMP END MASTER
!$OMP BARRIER

   END SUBROUTINE close_time_counter

END MODULE time_mod
