source: NEMO/branches/2019/dev_r11879_ENHANCE-05_SimonM-Harmonic_Analysis/src/OCE/DIA/diadetide.F90 @ 12097

MODULE diadetide
   !!======================================================================
   !!                      ***  MODULE diadetide  ***
   !! Computation of weights for daily detided model diagnostics
   !!======================================================================
   !! History :       !  2019  (S. Mueller)
   !!----------------------------------------------------------------------
   USE par_oce        , ONLY :   wp, jpi, jpj
   USE in_out_manager , ONLY :   lwp, numout
   USE iom            , ONLY :   iom_put
   USE dom_oce        , ONLY :   rdt, nsec_day
   USE phycst         , ONLY :   rpi
   USE tide_mod       , ONLY :   tide_harmo, jpmax_harmo, Wave
   USE xios

   IMPLICIT NONE
   PRIVATE

   LOGICAL, PUBLIC                               ::   lk_diadetide
   INTEGER                                       ::   ndiadetide
   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:)     ::   tdiadetide

   PUBLIC ::   dia_detide_init, dia_detide

   !!----------------------------------------------------------------------
   !! NEMO/OCE 4.0 , NEMO Consortium (2019)
   !! $Id$
   !! Software governed by the CeCILL license (see ./LICENSE)
   !!----------------------------------------------------------------------
CONTAINS

   SUBROUTINE dia_detide_init
      !!----------------------------------------------------------------------
      !!               ***  ROUTINE dia_detide_init  ***
      !!
      !! ** Purpose : initialisation of the weight computation for daily
      !!              detided diagnostics (currently M2-detiding only)
      !!
      !!----------------------------------------------------------------------

      LOGICAL                             ::   llxatt_enabled
      INTEGER, ALLOCATABLE, DIMENSION(:)  ::   itide                                    ! Tidal-constituent index
      REAL(wp), ALLOCATABLE, DIMENSION(:) ::   ztide_omega, ztide_u, ztide_v, ztide_f   ! Tidal-constituent parameters
      REAL(wp)                            ::   zdt
      INTEGER                             ::   jn

      ! Enquire detiding activation state (test for presence of detiding-related
      ! weights field and output file group)
      IF ( xios_is_valid_field( "diadetide_weight" ).AND.xios_is_valid_filegroup( "diadetide_files" ) ) THEN
         llxatt_enabled = .TRUE.
      ELSE
         llxatt_enabled = .FALSE.
      END IF
      lk_diadetide = llxatt_enabled

      IF (lwp) THEN
         WRITE (numout, *)
         WRITE (numout, *) 'dia_detide_init : weight computation for daily detided model diagnostics'
         WRITE (numout, *) '~~~~~~~~~~~~~~~'
         WRITE (numout, *) '                  lk_diadetide = ', lk_diadetide
      END IF

      IF (lk_diadetide) THEN
         ! Retrieve information about M2 tidal constituent
         ALLOCATE( ztide_omega(1), ztide_v(1), ztide_u(1), ztide_f(1), itide(1) )
         DO jn = 1, jpmax_harmo
            IF (TRIM( Wave(jn)%cname_tide ) == 'M2') itide(1) = jn
         END DO
         CALL tide_harmo( ztide_omega, ztide_v, ztide_u, ztide_f, itide, 1 )
         ! For M2, twice the tidal period spans slightly more than one full
         ! day. Compute the maximum number of equal intervals that span exactly
         ! twice the tidal period *and* whose mid-points fall within a 24-hour
         ! period from midnight to midnight.
         zdt = Wave(itide(1))%nt * 2.0_wp * rpi / ztide_omega(1)
         ndiadetide = FLOOR( zdt / ( zdt - 86400.0_wp ) )
         DEALLOCATE( ztide_omega, ztide_v, ztide_u, ztide_f, itide )
         ! Compute mid-points of the intervals to be included in the detided
         ! average
         ALLOCATE ( tdiadetide(ndiadetide) )
         DO jn = 1, ndiadetide
            tdiadetide(jn) = ( REAL( jn, KIND=wp) - 0.5_wp ) * zdt / REAL( ndiadetide, KIND=wp ) - ( zdt - 86400.0_wp ) * 0.5_wp
         END DO
      END IF

   END SUBROUTINE dia_detide_init

   SUBROUTINE dia_detide( kt )
      !!----------------------------------------------------------------------
      !!                  ***  ROUTINE dia_detide  ***
      !!
      !! ** Purpose : weight computation for daily detided model diagnostics
      !!----------------------------------------------------------------------

      INTEGER, INTENT(in)          ::   kt
      REAL(wp), DIMENSION(jpi,jpj) ::   zwght_2D
      REAL(wp)                     ::   zwght, ztmp
      INTEGER                      ::   jn

      ! Compute detiding weight at the current time-step; the daily total weight
      ! is one, and the daily summation of a diagnosed field multiplied by this
      ! weight should provide daily detided averages
      zwght = 0.0_wp
      DO jn = 1, ndiadetide
         ztmp = ( tdiadetide(jn) - REAL( nsec_day, KIND=wp ) ) / rdt
         IF ( ( ztmp < 0.5_wp ).AND.( ztmp >= -0.5_wp ) ) THEN
            zwght = zwght + 1.0_wp / REAL( ndiadetide, KIND=wp )
         END IF
      END DO
      zwght_2D(:,:) = zwght
      CALL iom_put( "diadetide_weight", zwght_2D)

   END SUBROUTINE dia_detide

END MODULE diadetide