Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/DOM/domqco.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/DOM/domqco.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/DOM/domqco.F90	(revision 12779)
@@ -114,8 +114,4 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) :: Kbb, Kmm, Kaa
-      !!----------------------------------------------------------------------
-      INTEGER ::   ji, jj, jk
-      INTEGER ::   ii0, ii1, ij0, ij1
-      REAL(wp)::   zcoef
       !!----------------------------------------------------------------------
       !
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/oce.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/oce.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/oce.F90	(revision 12779)
@@ -35,6 +35,4 @@
    !! ------------
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ssh, uu_b,  vv_b   !: SSH [m] and barotropic velocities [m/s]
-   !REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:,:) ::   ssh           !: SSH [m]
-   !REAL(wp), PUBLIC, ALLOCATABLE, SAVE        , DIMENSION(:,:,:) ::   uu_b,  vv_b   !: barotropic velocities [m/s]
    
    !! Arrays at barotropic time step:                   ! befbefore! before !  now   ! after  !
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/stpMLF.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/stpMLF.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/stpMLF.F90	(revision 12779)
@@ -245,14 +245,14 @@
       IF( lk_diamlr  )   CALL dia_mlr                           ! Update time used in multiple-linear-regression analysis
 
+      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
+      ! Now ssh filtering
+      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+                         CALL ssh_atf    ( kstp, Nbb, Nnn, Naa, ssh )            ! time filtering of "now" sea surface height
+                         CALL dom_qco_r3c( ssh(:,:,Nnn), r3t_f, r3u_f, r3v_f )   ! "now" ssh/h_0 ratio from filtrered ssh
 #if defined key_top
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
       ! Passive Tracer Model
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-                         ALLOCATE( zssh_f(jpi,jpj) )
-                         CALL ssh_atf    ( kstp, Nbb, Nnn, Naa , ssh, zssh_f )   ! time filtering of "now" sea surface height
-                         CALL dom_qco_r3c( zssh_f, r3t_f, r3u_f, r3v_f )         ! "now" ssh/h_0 ratio from filtrered ssh
-                         !
                          CALL trc_stp    ( kstp, Nbb, Nnn, Nrhs, Naa )           ! time-stepping
-                         DEALLOCATE( zssh_f )
 #endif
 
@@ -300,8 +300,6 @@
 !!    place.
 !!
-                         CALL zdyn_ts       ( Nnn, Naa, uu, vv )                    ! barotrope ajustment
+                         CALL mlf_baro_corr ( Nnn, Naa, uu, vv )                    ! barotrope ajustment
                          CALL finalize_sbc  ( kstp, Nbb, Naa, uu, vv, ts )          ! boundary condifions
-                         CALL ssh_atf       ( kstp, Nbb, Nnn, Naa, ssh )            ! time filtering of "now" sea surface height
-                         CALL dom_qco_r3c   ( ssh(:,:,Nnn), r3t_f, r3u_f, r3v_f )   ! "now" ssh/h_0 ratio from filtrered ssh
                          CALL tra_atf_qco   ( kstp, Nbb, Nnn, Naa, ts )             ! time filtering of "now" tracer arrays
                          CALL dyn_atf_qco   ( kstp, Nbb, Nnn, Naa, uu, vv  )        ! time filtering of "now" velocities and scale factors
@@ -377,7 +375,7 @@
 
 
-   SUBROUTINE zdyn_ts (Kmm, Kaa, puu, pvv)
-      !!----------------------------------------------------------------------
-      !!                  ***  ROUTINE zdyn_ts  ***
+   SUBROUTINE mlf_baro_corr (Kmm, Kaa, puu, pvv)
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE mlf_baro_corr  ***
       !!
       !! ** Purpose :   Finalize after horizontal velocity.
@@ -426,5 +424,5 @@
       ENDIF
       !
-   END SUBROUTINE zdyn_ts
+   END SUBROUTINE mlf_baro_corr
 
 
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OFF/dtadyn.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OFF/dtadyn.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OFF/dtadyn.F90	(revision 12779)
@@ -23,5 +23,9 @@
    USE c1d             ! 1D configuration: lk_c1d
    USE dom_oce         ! ocean domain: variables
+#if ! defined key_qco 
    USE domvvl          ! variable volume
+#else
+   USE domqco
+#endif
    USE zdf_oce         ! ocean vertical physics: variables
    USE sbc_oce         ! surface module: variables
@@ -52,5 +56,7 @@
    PUBLIC   dta_dyn_sed        ! called by nemo_gcm
    PUBLIC   dta_dyn_atf        ! called by nemo_gcm
+#if ! defined key_qco
    PUBLIC   dta_dyn_sf_interp  ! called by nemo_gcm
+#endif
 
    CHARACTER(len=100) ::   cn_dir          !: Root directory for location of ssr files
@@ -128,26 +134,31 @@
       IF( l_ldfslp .AND. .NOT.lk_c1d )   CALL  dta_dyn_slp( kt, Kbb, Kmm )    ! Computation of slopes
       !
-      ts(:,:,:,jp_tem,Kmm) = sf_dyn(jf_tem)%fnow(:,:,:)  * tmask(:,:,:)    ! temperature
-      ts(:,:,:,jp_sal,Kmm) = sf_dyn(jf_sal)%fnow(:,:,:)  * tmask(:,:,:)    ! salinity
-      wndm(:,:)         = sf_dyn(jf_wnd)%fnow(:,:,1)  * tmask(:,:,1)    ! wind speed - needed for gas exchange
-      fmmflx(:,:)       = sf_dyn(jf_fmf)%fnow(:,:,1)  * tmask(:,:,1)    ! downward salt flux (v3.5+)
-      fr_i(:,:)         = sf_dyn(jf_ice)%fnow(:,:,1)  * tmask(:,:,1)    ! Sea-ice fraction
-      qsr (:,:)         = sf_dyn(jf_qsr)%fnow(:,:,1)  * tmask(:,:,1)    ! solar radiation
-      emp (:,:)         = sf_dyn(jf_emp)%fnow(:,:,1)  * tmask(:,:,1)    ! E-P
+      ts    (:,:,:,jp_tem,Kmm) = sf_dyn(jf_tem)%fnow(:,:,:) * tmask(:,:,:)    ! temperature
+      ts    (:,:,:,jp_sal,Kmm) = sf_dyn(jf_sal)%fnow(:,:,:) * tmask(:,:,:)    ! salinity
+      wndm  (:,:)              = sf_dyn(jf_wnd)%fnow(:,:,1) * tmask(:,:,1)    ! wind speed - needed for gas exchange
+      fmmflx(:,:)              = sf_dyn(jf_fmf)%fnow(:,:,1) * tmask(:,:,1)    ! downward salt flux (v3.5+)
+      fr_i  (:,:)              = sf_dyn(jf_ice)%fnow(:,:,1) * tmask(:,:,1)    ! Sea-ice fraction
+      qsr   (:,:)              = sf_dyn(jf_qsr)%fnow(:,:,1) * tmask(:,:,1)    ! solar radiation
+      emp   (:,:)              = sf_dyn(jf_emp)%fnow(:,:,1) * tmask(:,:,1)    ! E-P
       IF( ln_dynrnf ) THEN 
-         rnf (:,:)      = sf_dyn(jf_rnf)%fnow(:,:,1) * tmask(:,:,1)    ! E-P
-         IF( ln_dynrnf_depth .AND. .NOT. ln_linssh )    CALL  dta_dyn_hrnf(Kmm)
-      ENDIF
-      !
-      uu(:,:,:,Kmm)        = sf_dyn(jf_uwd)%fnow(:,:,:) * umask(:,:,:)    ! effective u-transport
-      vv(:,:,:,Kmm)        = sf_dyn(jf_vwd)%fnow(:,:,:) * vmask(:,:,:)    ! effective v-transport
-      ww(:,:,:)        = sf_dyn(jf_wwd)%fnow(:,:,:) * tmask(:,:,:)    ! effective v-transport
+         rnf(:,:)              = sf_dyn(jf_rnf)%fnow(:,:,1) * tmask(:,:,1)    ! E-P
+         IF( ln_dynrnf_depth .AND. .NOT.ln_linssh )   CALL dta_dyn_hrnf( Kmm )
+      ENDIF
+      !
+      uu(:,:,:,Kmm)            = sf_dyn(jf_uwd)%fnow(:,:,:) * umask(:,:,:)    ! effective u-transport
+      vv(:,:,:,Kmm)            = sf_dyn(jf_vwd)%fnow(:,:,:) * vmask(:,:,:)    ! effective v-transport
+      ww(:,:,:)                = sf_dyn(jf_wwd)%fnow(:,:,:) * tmask(:,:,:)    ! effective v-transport
       !
       IF( .NOT.ln_linssh ) THEN
          ALLOCATE( zemp(jpi,jpj) , zhdivtr(jpi,jpj,jpk) )
-         zhdivtr(:,:,:) = sf_dyn(jf_div)%fnow(:,:,:)  * tmask(:,:,:)    ! effective u-transport
+         zhdivtr(:,:,:) = sf_dyn(jf_div )%fnow(:,:,:) * tmask(:,:,:)    ! effective u-transport
          emp_b  (:,:)   = sf_dyn(jf_empb)%fnow(:,:,1) * tmask(:,:,1)    ! E-P
          zemp   (:,:)   = ( 0.5_wp * ( emp(:,:) + emp_b(:,:) ) + rnf(:,:) + fwbcorr ) * tmask(:,:,1)
-         CALL dta_dyn_ssh( kt, zhdivtr, ssh(:,:,Kbb), zemp, ssh(:,:,Kaa), e3t(:,:,:,Kaa) )  !=  ssh, vertical scale factor & vertical transport
+#if defined key_qco
+         CALL dta_dyn_ssh( kt, zhdivtr, ssh(:,:,Kbb), zemp, ssh(:,:,Kaa) )
+         CALL dom_qco_r3c( ssh(:,:,Kaa), r3t(:,:,Kaa), r3u(:,:,Kaa), r3v(:,:,Kaa) )
+#else
+         CALL dta_dyn_ssh( kt, zhdivtr, ssh(:,:,Kbb), zemp, ssh(:,:,Kaa), e3t(:,:,:,Kaa) )  !=  ssh, vertical scale factor
+#endif
          DEALLOCATE( zemp , zhdivtr )
          !                                           Write in the tracer restart file
@@ -329,7 +340,11 @@
         ENDIF
         !
+#if defined key_qco
+        CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb) )
+        CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm) )
+#else
         DO jk = 1, jpkm1
-           e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( 1._wp + ssh(:,:,Kmm) * tmask(:,:,1) / ( ht_0(:,:) + 1.0 - tmask(:,:,1) ) )
-        ENDDO
+           e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( 1._wp + ssh(:,:,Kmm) * r1_ht_0(:,:) * tmask(:,:,jk) )
+        END DO
         e3t(:,:,jpk,Kaa) = e3t_0(:,:,jpk)
 
@@ -342,5 +357,5 @@
         ! ------------------------------------
         CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w(:,:,:,Kmm), 'W' )
-  
+!!gm this should be computed from ssh(Kbb)  
         e3t(:,:,:,Kbb)  = e3t(:,:,:,Kmm)
         e3u(:,:,:,Kbb)  = e3u(:,:,:,Kmm)
@@ -366,5 +381,6 @@
         gdepw(:,:,:,Kbb) = gdepw(:,:,:,Kmm)
         !
-      ENDIF
+     ENDIF
+#endif
       !
       IF( ln_dynrnf .AND. ln_dynrnf_depth ) THEN       ! read depht over which runoffs are distributed
@@ -389,4 +405,5 @@
             ENDIF
          END_2D
+!!st pourquoi on n'utilise pas le gde3w ici plutôt que de faire une boucle ? 
          DO_2D_11_11
             h_rnf(ji,jj) = 0._wp
@@ -413,4 +430,5 @@
    END SUBROUTINE dta_dyn_init
 
+   
    SUBROUTINE dta_dyn_sed( kt, Kmm )
       !!----------------------------------------------------------------------
@@ -529,4 +547,5 @@
    END SUBROUTINE dta_dyn_sed_init
 
+   
    SUBROUTINE dta_dyn_atf( kt, Kbb, Kmm, Kaa )
      !!---------------------------------------------------------------------
@@ -551,5 +570,7 @@
       !
    END SUBROUTINE dta_dyn_atf
+
    
+#if ! defined key_qco    
    SUBROUTINE dta_dyn_sf_interp( kt, Kmm )
       !!---------------------------------------------------------------------
@@ -588,5 +609,7 @@
       !
    END SUBROUTINE dta_dyn_sf_interp
-
+#endif
+
+   
    SUBROUTINE dta_dyn_ssh( kt, phdivtr, psshb,  pemp, pssha, pe3ta )
       !!----------------------------------------------------------------------
@@ -606,5 +629,5 @@
       !!          The boundary conditions are w=0 at the bottom (no flux)
       !!
-      !! ** action  :   ssh(:,:,Kaa) / e3t(:,:,:,Kaa) / ww
+      !! ** action  :   ssh(:,:,Kaa) / e3t(:,:,k,Kaa) / ww
       !!
       !! Reference  : Leclair, M., and G. Madec, 2009, Ocean Modelling.
@@ -630,9 +653,10 @@
       !                                                ! Sea surface  elevation time-stepping
       pssha(:,:) = ( psshb(:,:) - z2dt * ( r1_rho0 * pemp(:,:)  + zhdiv(:,:) ) ) * ssmask(:,:)
-      !                                                 ! 
-      !                                                 ! After acale factors at t-points ( z_star coordinate )
-      DO jk = 1, jpkm1
-        pe3ta(:,:,jk) = e3t_0(:,:,jk) * ( 1._wp + pssha(:,:) * tmask(:,:,1) / ( ht_0(:,:) + 1.0 - tmask(:,:,1) ) )
-      END DO
+      !
+      IF( PRESENT( pe3ta ) ) THEN                      ! After acale factors at t-points ( z_star coordinate )
+         DO jk = 1, jpkm1
+            pe3ta(:,:,jk) = e3t_0(:,:,jk) * ( 1._wp + pssha(:,:) * r1_ht_0(:,:) * tmask(:,:,jk) )
+         END DO
+      ENDIF
       !
    END SUBROUTINE dta_dyn_ssh
@@ -657,4 +681,5 @@
       !!----------------------------------------------------------------------
       !
+!!st code dupliqué même remarque que plus haut pourquoi ne pas utiliser gdepw ?
       DO_2D_11_11
          h_rnf(ji,jj) = 0._wp
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OFF/nemogcm.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OFF/nemogcm.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OFF/nemogcm.F90	(revision 12779)
@@ -28,4 +28,7 @@
    USE usrdef_nam     ! user defined configuration
    USE eosbn2         ! equation of state            (eos bn2 routine)
+#if defined key_qco
+   USE domqco         ! tools for scale factor         (dom_qco_r3c  routine)
+#endif
    !              ! ocean physics
    USE ldftra         ! lateral diffusivity setting    (ldf_tra_init routine)
@@ -117,7 +120,17 @@
                                 CALL dta_dyn    ( istp, Nbb, Nnn, Naa )       ! Interpolation of the dynamical fields
 #endif
+#if ! defined key_sed_off
+         IF( .NOT.ln_linssh ) THEN
+                                CALL dta_dyn_atf( istp, Nbb, Nnn, Naa )       ! time filter of sea  surface height and vertical scale factors
+# if defined key_qco
+                                CALL dom_qco_r3c( ssh(:,:,Kmm), r3t_f, r3u_f, r3v_f )
+# endif
+         ENDIF
                                 CALL trc_stp    ( istp, Nbb, Nnn, Nrhs, Naa ) ! time-stepping
-#if ! defined key_sed_off
-         IF( .NOT.ln_linssh )   CALL dta_dyn_atf( istp, Nbb, Nnn, Naa )       ! time filter of sea  surface height and vertical scale factors
+# if defined key_qco
+                                !r3t(:,:,Kmm) = r3t_f(:,:)                     ! update ssh to h0 ratio
+                                !r3u(:,:,Kmm) = r3u_f(:,:)
+                                !r3v(:,:,Kmm) = r3v_f(:,:)
+# endif
 #endif
          ! Swap time levels
@@ -127,7 +140,9 @@
          Naa = Nrhs
          !
-#if ! defined key_sed_off
+#if ! defined key_qco
+# if ! defined key_sed_off
          IF( .NOT.ln_linssh )   CALL dta_dyn_sf_interp( istp, Nnn )  ! calculate now grid parameters
-#endif
+# endif
+#endif         
                                 CALL stp_ctl    ( istp, indic )  ! Time loop: control and print
          istp = istp + 1
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zbio.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zbio.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zbio.F90	(revision 12779)
@@ -40,4 +40,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zche.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zche.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zche.F90	(revision 12779)
@@ -132,4 +132,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zfechem.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zfechem.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zfechem.F90	(revision 12779)
@@ -33,4 +33,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zprod.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zprod.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zprod.F90	(revision 12779)
@@ -48,4 +48,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zsink.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zsink.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p4zsink.F90	(revision 12779)
@@ -40,4 +40,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p5zprod.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p5zprod.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/P4Z/p5zprod.F90	(revision 12779)
@@ -52,4 +52,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/SED/oce_sed.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/SED/oce_sed.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/SED/oce_sed.F90	(revision 12779)
@@ -13,5 +13,8 @@
    USE dom_oce , ONLY :   glamt     =>   glamt          !: longitude of t-point (degre)
    USE dom_oce , ONLY :   gphit     =>   gphit          !: latitude  of t-point (degre)
-!!st USE dom_oce , ONLY :   e3t       =>   e3t            !: latitude  of t-point (degre)
+!!st 
+#if ! defined key_qco
+   USE dom_oce , ONLY :   e3t       =>   e3t            !: latitude  of t-point (degre)
+#endif
    USE dom_oce , ONLY :   e3t_1d    =>   e3t_1d         !: reference depth of t-points (m)
    USE dom_oce , ONLY :   gdepw_0   =>   gdepw_0        !: reference depth of t-points (m)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/trcwri_pisces.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/trcwri_pisces.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/PISCES/trcwri_pisces.F90	(revision 12779)
@@ -21,4 +21,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcadv.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcadv.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcadv.F90	(revision 12779)
@@ -58,5 +58,6 @@
    INTEGER, PARAMETER ::   np_UBS     = 4   ! 3rd order Upstream Biased Scheme
    INTEGER, PARAMETER ::   np_QCK     = 5   ! QUICK scheme
-   
+
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcatf.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcatf.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcatf.F90	(revision 12779)
@@ -31,5 +31,9 @@
    USE trd_oce
    USE trdtra
+# if defined key_qco
+   USE traatfqco
+# else
    USE traatf
+# endif
    USE bdy_oce   , ONLY: ln_bdy
    USE trcbdy          ! BDY open boundaries
@@ -50,4 +54,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
@@ -107,10 +112,10 @@
             DO jn = 1, jptra
                CALL trd_tra( kt, Kmm, Kaa, 'TRC', jn, jptra_zdfp, ztrdt(:,:,:,jn) )
-            ENDDO
+            END DO
          ENDIF
 
          ! total trend for the non-time-filtered variables. 
          zfact = 1.0 / rn_Dt
-         ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3t*T)/e3tn; e3tn cancel from ts(Kmm) terms
+         ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3ta*Ta)/e3tn; e3tn cancel from ts(Kmm) terms
          IF( ln_linssh ) THEN       ! linear sea surface height only
             DO jn = 1, jptra
@@ -129,5 +134,5 @@
          DO jn = 1, jptra
             CALL trd_tra( kt, Kmm, Kaa, 'TRC', jn, jptra_tot, ztrdt(:,:,:,jn) )
-         ENDDO
+         END DO
          !
          IF( ln_linssh ) THEN       ! linear sea surface height only
@@ -146,14 +151,19 @@
             DO jn = 1, jptra
                CALL trd_tra( kt, Kmm, Kaa, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) )
-            ENDDO
+            END DO
          END IF
          !
       ELSE     
          IF( .NOT. l_offline ) THEN ! Leap-Frog + Asselin filter time stepping
-            IF( ln_linssh ) THEN   ;   CALL tra_atf_fix( kt, Kbb, Kmm, Kaa, nittrc000,         'TRC', ptr, jptra )                     !     linear ssh
-            ELSE                   ;   CALL tra_atf_vvl( kt, Kbb, Kmm, Kaa, nittrc000, rn_Dt, 'TRC', ptr, sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh
+# if defined key_qco
+            IF( ln_linssh ) THEN   ;   CALL tra_atf_fix_lf( kt, Kbb, Kmm, Kaa, nittrc000,        'TRC', ptr, jptra )                     !     linear ssh
+            ELSE                   ;   CALL tra_atf_qco_lf( kt, Kbb, Kmm, Kaa, nittrc000, rn_Dt, 'TRC', ptr, sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh
+# else
+            IF( ln_linssh ) THEN   ;   CALL tra_atf_fix   ( kt, Kbb, Kmm, Kaa, nittrc000,        'TRC', ptr, jptra )                     !     linear ssh
+            ELSE                   ;   CALL tra_atf_vvl   ( kt, Kbb, Kmm, Kaa, nittrc000, rn_Dt, 'TRC', ptr, sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh
+# endif
             ENDIF
          ELSE
-                                       CALL trc_atf_off( kt, Kbb, Kmm, Kaa, ptr )       ! offline 
+                                       CALL trc_atf_off   ( kt, Kbb, Kmm, Kaa, ptr )       ! offline
          ENDIF
          !
@@ -182,5 +192,5 @@
    END SUBROUTINE trc_atf
 
-
+# if ! defined key_qco
    SUBROUTINE trc_atf_off( kt, Kbb, Kmm, Kaa, ptr )
       !!----------------------------------------------------------------------
@@ -198,8 +208,8 @@
       !!                This can be summurized for tempearture as:
       !!             ztm = ( e3t_n*tn + rbcp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )   ln_dynhpg_imp = T
-      !!                  /( e3t(:,:,:,Kmm)    + rbcp*[ e3t(:,:,:,Kbb)    - 2 e3t(:,:,:,Kmm)    + e3t(:,:,:,Kaa)    ] )   
+      !!                  /( e3t(:,:,jk,Kmm)    + rbcp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] )   
       !!             ztm = 0                                                       otherwise
       !!             tb  = ( e3t_n*tn + rn_atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )
-      !!                  /( e3t(:,:,:,Kmm)    + rn_atfp*[ e3t(:,:,:,Kbb)    - 2 e3t(:,:,:,Kmm)    + e3t(:,:,:,Kaa)    ] )
+      !!                  /( e3t(:,:,jk,Kmm)    + rn_atfp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] )
       !!             tn  = ta 
       !!             ta  = zt        (NB: reset to 0 after eos_bn2 call)
@@ -257,5 +267,79 @@
       !
    END SUBROUTINE trc_atf_off
-
+# else
+   SUBROUTINE trc_atf_off( kt, Kbb, Kmm, Kaa, ptr )
+      !!----------------------------------------------------------------------
+      !!                   ***  ROUTINE tra_atf_off  ***
+      !!
+      !!          !!!!!!!!!!!!!!!!! REWRITE HEADER COMMENTS !!!!!!!!!!!!!!
+      !!
+      !! ** Purpose :   Time varying volume: apply the Asselin time filter  
+      !! 
+      !! ** Method  : - Apply a thickness weighted Asselin time filter on now fields.
+      !!              - save in (ta,sa) a thickness weighted average over the three 
+      !!             time levels which will be used to compute rdn and thus the semi-
+      !!             implicit hydrostatic pressure gradient (ln_dynhpg_imp = T)
+      !!              - swap tracer fields to prepare the next time_step.
+      !!                This can be summurized for tempearture as:
+      !!             ztm = ( e3t_n*tn + rbcp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )   ln_dynhpg_imp = T
+      !!                  /( e3t(:,:,jk,Kmm)    + rbcp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] )   
+      !!             ztm = 0                                                       otherwise
+      !!             tb  = ( e3t_n*tn + rn_atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )
+      !!                  /( e3t(:,:,jk,Kmm)    + rn_atfp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] )
+      !!             tn  = ta 
+      !!             ta  = zt        (NB: reset to 0 after eos_bn2 call)
+      !!
+      !! ** Action  : - (tb,sb) and (tn,sn) ready for the next time step
+      !!              - (ta,sa) time averaged (t,s)   (ln_dynhpg_imp = T)
+      !!----------------------------------------------------------------------
+      INTEGER                                   , INTENT(in   ) ::  kt            ! ocean time-step index
+      INTEGER                                   , INTENT(in   ) ::  Kbb, Kmm, Kaa ! time level indices
+      REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) ::  ptr           ! passive tracers
+      !!     
+      INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices
+      REAL(wp) ::   ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar
+      REAL(wp) ::   ze3t_b, ze3t_n, ze3t_a, ze3t_f           !   -      -
+      !!----------------------------------------------------------------------
+      !
+      IF( kt == nittrc000 )  THEN
+         IF(lwp) WRITE(numout,*)
+         IF(lwp) WRITE(numout,*) 'trc_atf_off : Asselin time filtering'
+         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
+         IF( .NOT. ln_linssh ) THEN
+            rfact1 = rn_atfp * rn_Dt
+            rfact2 = rfact1 / rho0
+         ENDIF
+        !  
+      ENDIF
+      !
+      DO jn = 1, jptra      
+         DO_3D_11_11( 1, jpkm1 )
+            ze3t_b = 1._wp + r3t(ji,jj,Kbb) * tmask(ji,jj,jk)
+            ze3t_n = 1._wp + r3t(ji,jj,Kmm) * tmask(ji,jj,jk)
+            ze3t_a = 1._wp + r3t(ji,jj,Kaa) * tmask(ji,jj,jk)
+            !                                         ! tracer content at Before, now and after
+            ztc_b  = ptr(ji,jj,jk,jn,Kbb) * ze3t_b
+            ztc_n  = ptr(ji,jj,jk,jn,Kmm) * ze3t_n
+            ztc_a  = ptr(ji,jj,jk,jn,Kaa) * ze3t_a
+            !
+            ztc_d  = ztc_a  - 2. * ztc_n  + ztc_b
+            !
+            ze3t_f = 1._wp + r3t_f(ji,jj)*tmask(ji,jj,jk)
+            ztc_f  = ztc_n  + rn_atfp * ztc_d
+            !
+            IF( .NOT. ln_linssh .AND. jk == mikt(ji,jj) ) THEN           ! first level 
+               ztc_f  = ztc_f  - rfact1 * ( sbc_trc(ji,jj,jn) - sbc_trc_b(ji,jj,jn) )
+            ENDIF
+
+            ze3t_f = 1.e0 / ze3t_f
+            ptr(ji,jj,jk,jn,Kmm) = ztc_f * ze3t_f     ! time filtered "now" field
+            !
+         END_3D
+         ! 
+      END DO
+      !
+   END SUBROUTINE trc_atf_off
+# endif
+   
 #else
    !!----------------------------------------------------------------------
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcdmp.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcdmp.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcdmp.F90	(revision 12779)
@@ -45,4 +45,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcldf.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcldf.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcldf.F90	(revision 12779)
@@ -44,4 +44,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcsbc.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcsbc.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trcsbc.F90	(revision 12779)
@@ -30,4 +30,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
@@ -49,5 +50,5 @@
       !!            The surface freshwater flux modify the ocean volume
       !!         and thus the concentration of a tracer as :
-      !!            tr(Krhs) = tr(Krhs) + emp * tr(Kmm) / e3t   for k=1
+      !!            tr(Krhs) = tr(Krhs) + emp * tr(Kmm) / e3t_   for k=1
       !!         where emp, the surface freshwater budget (evaporation minus
       !!         precipitation ) given in kg/m2/s is divided
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trdmxl_trc.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trdmxl_trc.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/TRP/trdmxl_trc.F90	(revision 12779)
@@ -51,4 +51,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcbc.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcbc.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcbc.F90	(revision 12779)
@@ -48,4 +48,5 @@
    !! * Substitutions
 #  include "do_loop_substitute.h90"
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcdta.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcdta.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcdta.F90	(revision 12779)
@@ -41,4 +41,5 @@
    !! Substitutions
 #include "do_loop_substitute.h90"
+#include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
@@ -206,5 +207,5 @@
                      ztp(jk) = ptrcdta(ji,jj,jpkm1)
                   ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1
-                     DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1)
+                     DO jkk = 1, jpkm1                                  ! when  gdept_1d(jkk) < zl < gdept_1d(jkk+1)
                         IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN
                            zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk))
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcini.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcini.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcini.F90	(revision 12779)
@@ -30,5 +30,6 @@
    
    PUBLIC   trc_init   ! called by opa
-
+   
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcrst.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcrst.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcrst.F90	(revision 12779)
@@ -33,4 +33,5 @@
    PUBLIC   trc_rst_cal
 
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcstp.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcstp.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcstp.F90	(revision 12779)
@@ -36,5 +36,6 @@
    REAL(wp) ::   rsecfst, rseclast       ! ???
    REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE ::   qsr_arr   ! save qsr during TOP time-step
-
+   
+#  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/TOP 4.0 , NEMO Consortium (2018)
Index: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcwri.F90
===================================================================
--- NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcwri.F90	(revision 12761)
+++ NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/TOP/trcwri.F90	(revision 12779)
@@ -60,7 +60,9 @@
        CALL iom_put( "e3v_0", e3v_0(:,:,:) )
        !
+#if ! defined key_qco
        CALL iom_put( "e3t" , e3t(:,:,:,Kmm) )
        CALL iom_put( "e3u" , e3u(:,:,:,Kmm) )
        CALL iom_put( "e3v" , e3v(:,:,:,Kmm) )
+#endif 
        !
       ENDIF
