Index: NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/DIA/diaptr.F90
===================================================================
--- NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/DIA/diaptr.F90	(revision 13360)
+++ NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/DIA/diaptr.F90	(revision 14037)
@@ -22,4 +22,5 @@
    USE oce              ! ocean dynamics and active tracers
    USE dom_oce          ! ocean space and time domain
+   USE domain, ONLY : dom_tile
    USE phycst           ! physical constants
    !
@@ -32,20 +33,22 @@
    PRIVATE
 
+   INTERFACE ptr_sum
+      MODULE PROCEDURE ptr_sum_3d, ptr_sum_2d
+   END INTERFACE
+
    INTERFACE ptr_sj
       MODULE PROCEDURE ptr_sj_3d, ptr_sj_2d
    END INTERFACE
 
-   PUBLIC   ptr_sj         ! call by tra_ldf & tra_adv routines
-   PUBLIC   ptr_sjk        ! 
-   PUBLIC   dia_ptr_init   ! call in memogcm
    PUBLIC   dia_ptr        ! call in step module
    PUBLIC   dia_ptr_hst    ! called from tra_ldf/tra_adv routines
 
-   !                                  !!** namelist  namptr  **
-   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hstr_adv, hstr_ldf, hstr_eiv   !: Heat/Salt TRansports(adv, diff, Bolus.)
-   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hstr_ove, hstr_btr, hstr_vtr   !: heat Salt TRansports(overturn, baro, merional)
-
-   LOGICAL , PUBLIC ::   l_diaptr        !: tracers  trend flag (set from namelist in trdini)
-   INTEGER, PARAMETER, PUBLIC ::   nptr = 5  ! (glo, atl, pac, ind, ipc)
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   hstr_adv, hstr_ldf, hstr_eiv   !: Heat/Salt TRansports(adv, diff, Bolus.)
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   hstr_ove, hstr_btr, hstr_vtr   !: heat Salt TRansports(overturn, baro, merional)
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   pvtr_int, pzon_int             !: Other zonal integrals
+
+   LOGICAL, PUBLIC    ::   l_diaptr       !: tracers  trend flag
+   INTEGER, PARAMETER ::   jp_msk = 3
+   INTEGER, PARAMETER ::   jp_vtr = 4
 
    REAL(wp) ::   rc_sv    = 1.e-6_wp   ! conversion from m3/s to Sverdrup
@@ -56,9 +59,6 @@
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk34 ! mask out Southern Ocean (=0 south of 34°S)
 
-   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:)   :: p_fval1d
-   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:,:) :: p_fval2d
-
-   LOGICAL ::   ll_init = .TRUE.        !: tracers  trend flag (set from namelist in trdini)
-   
+   LOGICAL ::   ll_init = .TRUE.        !: tracers  trend flag
+
    !! * Substitutions
 #  include "do_loop_substitute.h90"
@@ -77,117 +77,338 @@
       INTEGER                         , INTENT(in)           ::   kt     ! ocean time-step index     
       INTEGER                         , INTENT(in)           ::   Kmm    ! time level index
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   pvtr   ! j-effective transport
+      REAL(wp), DIMENSION(A2D(nn_hls),jpk)    , INTENT(in), OPTIONAL ::   pvtr   ! j-effective transport
+      !!----------------------------------------------------------------------
+      !
+      IF( ln_timing )   CALL timing_start('dia_ptr')
+
+      IF( kt == nit000 .AND. ll_init )   CALL dia_ptr_init    ! -> will define l_diaptr and nbasin
+      !
+      IF( l_diaptr ) THEN
+         ! Calculate zonal integrals
+         IF( PRESENT( pvtr ) ) THEN
+            CALL dia_ptr_zint( Kmm, pvtr )
+         ELSE
+            CALL dia_ptr_zint( Kmm )
+         ENDIF
+
+         ! Calculate diagnostics only when zonal integrals have finished
+         IF( ntile == 0 .OR. ntile == nijtile ) CALL dia_ptr_iom(kt, Kmm, pvtr)
+      ENDIF
+
+      IF( ln_timing )   CALL timing_stop('dia_ptr')
+      !
+   END SUBROUTINE dia_ptr
+
+
+   SUBROUTINE dia_ptr_iom( kt, Kmm, pvtr )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE dia_ptr_iom  ***
+      !!----------------------------------------------------------------------
+      !! ** Purpose : Calculate diagnostics and send to XIOS
+      !!----------------------------------------------------------------------
+      INTEGER                         , INTENT(in)           ::   kt     ! ocean time-step index
+      INTEGER                         , INTENT(in)           ::   Kmm    ! time level index
+      REAL(wp), DIMENSION(A2D(nn_hls),jpk)    , INTENT(in), OPTIONAL ::   pvtr   ! j-effective transport
       !
       INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
-      REAL(wp) ::   zsfc,zvfc               ! local scalar
       REAL(wp), DIMENSION(jpi,jpj)     ::  z2d   ! 2D workspace
-      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zmask   ! 3D workspace
-      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  z3d    ! 3D workspace
-      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) ::  zts   ! 3D workspace
       REAL(wp), DIMENSION(jpj)      ::  zvsum, ztsum, zssum   ! 1D workspace
       !
       !overturning calculation
-      REAL(wp), DIMENSION(jpj,jpk,nptr) :: sjk, r1_sjk, v_msf  ! i-mean i-k-surface and its inverse
-      REAL(wp), DIMENSION(jpj,jpk,nptr) :: zt_jk, zs_jk ! i-mean T and S, j-Stream-Function
-
-      REAL(wp), DIMENSION(jpi,jpj,jpk,nptr)  :: z4d1, z4d2
-      REAL(wp), DIMENSION(jpi,jpj,nptr)      :: z3dtr ! i-mean T and S, j-Stream-Function
-      !!----------------------------------------------------------------------
-      !
-      IF( ln_timing )   CALL timing_start('dia_ptr')
-
-      IF( kt == nit000 .AND. ll_init )   CALL dia_ptr_init
-      !
-      IF( .NOT. l_diaptr )   RETURN
+      REAL(wp), DIMENSION(:,:,:  ), ALLOCATABLE ::   sjk, r1_sjk, v_msf  ! i-mean i-k-surface and its inverse
+      REAL(wp), DIMENSION(:,:,:  ), ALLOCATABLE ::   zt_jk, zs_jk        ! i-mean T and S, j-Stream-Function
+
+      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   z4d1, z4d2
+      REAL(wp), DIMENSION(:,:,:  ), ALLOCATABLE ::   z3dtr
+      !!----------------------------------------------------------------------
+      !
+      ALLOCATE( z3dtr(jpi,jpj,nbasin) )
 
       IF( PRESENT( pvtr ) ) THEN
          IF( iom_use( 'zomsf' ) ) THEN    ! effective MSF
-            DO jn = 1, nptr                                    ! by sub-basins
-               z4d1(1,:,:,jn) =  ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) )  ! zonal cumulative effective transport excluding closed seas
-               DO jk = jpkm1, 1, -1 
+            ALLOCATE( z4d1(jpi,jpj,jpk,nbasin) )
+            !
+            DO jn = 1, nbasin                                    ! by sub-basins
+               z4d1(1,:,:,jn) =  pvtr_int(:,:,jp_vtr,jn)                  ! zonal cumulative effective transport excluding closed seas
+               DO jk = jpkm1, 1, -1
                   z4d1(1,:,jk,jn) = z4d1(1,:,jk+1,jn) - z4d1(1,:,jk,jn)    ! effective j-Stream-Function (MSF)
                END DO
-               DO ji = 1, jpi
+               DO ji = 2, jpi
                   z4d1(ji,:,:,jn) = z4d1(1,:,:,jn)
                ENDDO
             END DO
             CALL iom_put( 'zomsf', z4d1 * rc_sv )
-         ENDIF
+            !
+            DEALLOCATE( z4d1 )
+         ENDIF
+         IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN
+            ALLOCATE( sjk(jpj,jpk,nbasin), r1_sjk(jpj,jpk,nbasin), v_msf(jpj,jpk,nbasin),   &
+               &      zt_jk(jpj,jpk,nbasin), zs_jk(jpj,jpk,nbasin) )
+            !
+            DO jn = 1, nbasin
+               sjk(:,:,jn) = pvtr_int(:,:,jp_msk,jn)
+               r1_sjk(:,:,jn) = 0._wp
+               WHERE( sjk(:,:,jn) /= 0._wp )   r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn)
+               ! i-mean T and S, j-Stream-Function, basin
+               zt_jk(:,:,jn) = pvtr_int(:,:,jp_tem,jn) * r1_sjk(:,:,jn)
+               zs_jk(:,:,jn) = pvtr_int(:,:,jp_sal,jn) * r1_sjk(:,:,jn)
+               v_msf(:,:,jn) = pvtr_int(:,:,jp_vtr,jn)
+               hstr_ove(:,jp_tem,jn) = SUM( v_msf(:,:,jn)*zt_jk(:,:,jn), 2 )
+               hstr_ove(:,jp_sal,jn) = SUM( v_msf(:,:,jn)*zs_jk(:,:,jn), 2 )
+               !
+            ENDDO
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sophtove', z3dtr )
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sopstove', z3dtr )
+            !
+            DEALLOCATE( sjk, r1_sjk, v_msf, zt_jk, zs_jk )
+         ENDIF
+
+         IF( iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN
+            ! Calculate barotropic heat and salt transport here 
+            ALLOCATE( sjk(jpj,1,nbasin), r1_sjk(jpj,1,nbasin) )
+            !
+            DO jn = 1, nbasin
+               sjk(:,1,jn) = SUM( pvtr_int(:,:,jp_msk,jn), 2 )
+               r1_sjk(:,1,jn) = 0._wp
+               WHERE( sjk(:,1,jn) /= 0._wp )   r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn)
+               !
+               zvsum(:) =    SUM( pvtr_int(:,:,jp_vtr,jn), 2 )
+               ztsum(:) =    SUM( pvtr_int(:,:,jp_tem,jn), 2 )
+               zssum(:) =    SUM( pvtr_int(:,:,jp_sal,jn), 2 )
+               hstr_btr(:,jp_tem,jn) = zvsum(:) * ztsum(:) * r1_sjk(:,1,jn)
+               hstr_btr(:,jp_sal,jn) = zvsum(:) * zssum(:) * r1_sjk(:,1,jn)
+               !
+            ENDDO
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sophtbtr', z3dtr )
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sopstbtr', z3dtr )
+            !
+            DEALLOCATE( sjk, r1_sjk )
+         ENDIF
+         !
+         hstr_ove(:,:,:) = 0._wp       ! Zero before next timestep
+         hstr_btr(:,:,:) = 0._wp
+         pvtr_int(:,:,:,:) = 0._wp
+      ELSE
+         IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' )  ) THEN    ! i-mean i-k-surface
+            ALLOCATE( z4d1(jpi,jpj,jpk,nbasin), z4d2(jpi,jpj,jpk,nbasin) )
+            !
+            DO jn = 1, nbasin
+               z4d1(1,:,:,jn) = pzon_int(:,:,jp_msk,jn)
+               DO ji = 2, jpi
+                  z4d1(ji,:,:,jn) = z4d1(1,:,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'zosrf', z4d1 )
+            !
+            DO jn = 1, nbasin
+               z4d2(1,:,:,jn) = pzon_int(:,:,jp_tem,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 )
+               DO ji = 2, jpi
+                  z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'zotem', z4d2 )
+            !
+            DO jn = 1, nbasin
+               z4d2(1,:,:,jn) = pzon_int(:,:,jp_sal,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 )
+               DO ji = 2, jpi
+                  z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'zosal', z4d2 )
+            !
+            DEALLOCATE( z4d1, z4d2 )
+         ENDIF
+         !
+         !                                ! Advective and diffusive heat and salt transport
+         IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN  
+            ! 
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sophtadv', z3dtr )
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sopstadv', z3dtr )
+         ENDIF
+         !
+         IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN  
+            ! 
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sophtldf', z3dtr )
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sopstldf', z3dtr )
+         ENDIF
+         !
+         IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN  
+            ! 
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sophteiv', z3dtr )
+            DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sopsteiv', z3dtr )
+         ENDIF
+         !
+         IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN
+             DO jn = 1, nbasin
+                z3dtr(1,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
+                DO ji = 2, jpi
+                   z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+                ENDDO
+             ENDDO
+             CALL iom_put( 'sophtvtr', z3dtr )
+             DO jn = 1, nbasin
+               z3dtr(1,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
+               DO ji = 2, jpi
+                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
+               ENDDO
+            ENDDO
+            CALL iom_put( 'sopstvtr', z3dtr )
+         ENDIF
+         !
+         IF( iom_use( 'uocetr_vsum_cumul' ) ) THEN
+            IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 )         ! Use full domain
+            CALL iom_get_var(  'uocetr_vsum_op', z2d ) ! get uocetr_vsum_op from xml
+            z2d(:,:) = ptr_ci_2d( z2d(:,:) )  
+            CALL iom_put( 'uocetr_vsum_cumul', z2d )
+            IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile )   ! Revert to tile domain
+         ENDIF
+         !
+         hstr_adv(:,:,:) = 0._wp       ! Zero before next timestep
+         hstr_ldf(:,:,:) = 0._wp
+         hstr_eiv(:,:,:) = 0._wp
+         hstr_vtr(:,:,:) = 0._wp
+         pzon_int(:,:,:,:) = 0._wp
+      ENDIF
+      !
+      DEALLOCATE( z3dtr )
+      !
+   END SUBROUTINE dia_ptr_iom
+
+
+   SUBROUTINE dia_ptr_zint( Kmm, pvtr )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE dia_ptr_zint ***
+      !!----------------------------------------------------------------------
+      !! ** Purpose : i and i-k sum operations on arrays
+      !!
+      !! ** Method  : - Call ptr_sjk (i sum) or ptr_sj (i-k sum) to perform the sum operation
+      !!              - Call ptr_sum to add this result to the sum over tiles
+      !!
+      !! ** Action  : pvtr_int - terms for volume streamfunction, heat/salt transport barotropic/overturning terms
+      !!              pzon_int - terms for i mean temperature/salinity
+      !!----------------------------------------------------------------------
+      INTEGER                     , INTENT(in)           :: Kmm          ! time level index
+      REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(in), OPTIONAL :: pvtr         ! j-effective transport
+      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE            :: zmask        ! 3D workspace
+      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE          :: zts          ! 4D workspace
+      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE            :: sjk, v_msf   ! Zonal sum: i-k surface area, j-effective transport
+      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE            :: zt_jk, zs_jk ! Zonal sum: i-k surface area * (T, S)
+      REAL(wp)                                           :: zsfc, zvfc   ! i-k surface area
+      INTEGER  ::   ji, jj, jk, jn                                       ! dummy loop indices
+      !!----------------------------------------------------------------------
+
+      IF( PRESENT( pvtr ) ) THEN
+         ! i sum of effective j transport excluding closed seas
+         IF( iom_use( 'zomsf' ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN
+            ALLOCATE( v_msf(A1Dj(nn_hls),jpk,nbasin) )
+
+            DO jn = 1, nbasin
+               v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) )
+            ENDDO
+
+            CALL ptr_sum( pvtr_int(:,:,jp_vtr,:), v_msf(:,:,:) )
+
+            DEALLOCATE( v_msf )
+         ENDIF
+
+         ! i sum of j surface area, j surface area - temperature/salinity product on V grid
          IF(  iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR.   &
             & iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN
-            ! define fields multiplied by scalar
+            ALLOCATE( zmask(A2D(nn_hls),jpk), zts(A2D(nn_hls),jpk,jpts), &
+               &      sjk(A1Dj(nn_hls),jpk,nbasin), &
+               &      zt_jk(A1Dj(nn_hls),jpk,nbasin), zs_jk(A1Dj(nn_hls),jpk,nbasin) )
+
             zmask(:,:,:) = 0._wp
             zts(:,:,:,:) = 0._wp
+
             DO_3D( 1, 0, 1, 1, 1, jpkm1 )
                zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
                zmask(ji,jj,jk)      = vmask(ji,jj,jk)      * zvfc
-               zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc  !Tracers averaged onto V grid
+               zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc !Tracers averaged onto V grid
                zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc
             END_3D
-         ENDIF
-         IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN
-            DO jn = 1, nptr
-               sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) )
-               r1_sjk(:,:,jn) = 0._wp
-               WHERE( sjk(:,:,jn) /= 0._wp )   r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn)
-               ! i-mean T and S, j-Stream-Function, basin
-               zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn)
-               zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn)
-               v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) ) 
-               hstr_ove(:,jp_tem,jn) = SUM( v_msf(:,:,jn)*zt_jk(:,:,jn), 2 )
-               hstr_ove(:,jp_sal,jn) = SUM( v_msf(:,:,jn)*zs_jk(:,:,jn), 2 )
-               !
-            ENDDO
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sophtove', z3dtr )
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sopstove', z3dtr )
-         ENDIF
-
-         IF( iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN
-            ! Calculate barotropic heat and salt transport here 
-            DO jn = 1, nptr
-               sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) )
-               r1_sjk(:,1,jn) = 0._wp
-               WHERE( sjk(:,1,jn) /= 0._wp )   r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn)
-               !
-               zvsum(:) = ptr_sj( pvtr(:,:,:), btmsk34(:,:,jn) )
-               ztsum(:) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) )
-               zssum(:) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) )
-               hstr_btr(:,jp_tem,jn) = zvsum(:) * ztsum(:) * r1_sjk(:,1,jn)
-               hstr_btr(:,jp_sal,jn) = zvsum(:) * zssum(:) * r1_sjk(:,1,jn)
-               !
-            ENDDO
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sophtbtr', z3dtr )
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sopstbtr', z3dtr )
-         ENDIF 
-         !
+
+            DO jn = 1, nbasin
+               sjk(:,:,jn)   = ptr_sjk( zmask(:,:,:)     , btmsk(:,:,jn) )
+               zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) )
+               zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) )
+            ENDDO
+
+            CALL ptr_sum( pvtr_int(:,:,jp_msk,:), sjk(:,:,:)   )
+            CALL ptr_sum( pvtr_int(:,:,jp_tem,:), zt_jk(:,:,:) )
+            CALL ptr_sum( pvtr_int(:,:,jp_sal,:), zs_jk(:,:,:) )
+
+            DEALLOCATE( zmask, zts, sjk, zt_jk, zs_jk )
+         ENDIF
       ELSE
-         !
-         zmask(:,:,:) = 0._wp
-         zts(:,:,:,:) = 0._wp
-         IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' )  ) THEN    ! i-mean i-k-surface 
+         ! i sum of j surface area - temperature/salinity product on T grid
+         IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' )  ) THEN
+            ALLOCATE( zmask(A2D(nn_hls),jpk), zts(A2D(nn_hls),jpk,jpts), &
+               &      sjk(A1Dj(nn_hls),jpk,nbasin), &
+               &      zt_jk(A1Dj(nn_hls),jpk,nbasin), zs_jk(A1Dj(nn_hls),jpk,nbasin) )
+
+            zmask(:,:,:) = 0._wp
+            zts(:,:,:,:) = 0._wp
+
             DO_3D( 1, 1, 1, 1, 1, jpkm1 )
                zsfc = e1t(ji,jj) * e3t(ji,jj,jk,Kmm)
@@ -196,88 +417,24 @@
                zts(ji,jj,jk,jp_sal) = ts(ji,jj,jk,jp_sal,Kmm) * zsfc
             END_3D
-            !
-            DO jn = 1, nptr
-               zmask(1,:,:) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) )
-               z4d1(:,:,:,jn) = zmask(:,:,:)
-            ENDDO
-            CALL iom_put( 'zosrf', z4d1 )
-            !
-            DO jn = 1, nptr
-               z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) &
-                  &            / MAX( z4d1(1,:,:,jn), 10.e-15 )
-               DO ji = 1, jpi
-                  z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'zotem', z4d2 )
-            !
-            DO jn = 1, nptr
-               z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) &
-                  &            / MAX( z4d1(1,:,:,jn), 10.e-15 )
-               DO ji = 1, jpi
-                  z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'zosal', z4d2 )
-            !
-         ENDIF
-         !
-         !                                ! Advective and diffusive heat and salt transport
-         IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN  
-            ! 
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sophtadv', z3dtr )
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sopstadv', z3dtr )
-         ENDIF
-         !
-         IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN  
-            ! 
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sophtldf', z3dtr )
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sopstldf', z3dtr )
-         ENDIF
-         !
-         IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN  
-            ! 
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sophteiv', z3dtr )
-            DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sopsteiv', z3dtr )
-         ENDIF
-         !
+
+            DO jn = 1, nbasin
+               sjk(:,:,jn)   = ptr_sjk( zmask(:,:,:)     , btmsk(:,:,jn) )
+               zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) )
+               zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) )
+            ENDDO
+
+            CALL ptr_sum( pzon_int(:,:,jp_msk,:), sjk(:,:,:)   )
+            CALL ptr_sum( pzon_int(:,:,jp_tem,:), zt_jk(:,:,:) )
+            CALL ptr_sum( pzon_int(:,:,jp_sal,:), zs_jk(:,:,:) )
+
+            DEALLOCATE( zmask, zts, sjk, zt_jk, zs_jk )
+         ENDIF
+
+         ! i-k sum of j surface area - temperature/salinity product on V grid
          IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN
+            ALLOCATE( zts(A2D(nn_hls),jpk,jpts) )
+
             zts(:,:,:,:) = 0._wp
+
             DO_3D( 1, 0, 1, 1, 1, jpkm1 )
                zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
@@ -285,33 +442,12 @@
                zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc
             END_3D
-             CALL dia_ptr_hst( jp_tem, 'vtr', zts(:,:,:,jp_tem) )
-             CALL dia_ptr_hst( jp_sal, 'vtr', zts(:,:,:,jp_sal) )
-             DO jn = 1, nptr
-                z3dtr(1,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
-                DO ji = 1, jpi
-                   z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-                ENDDO
-             ENDDO
-             CALL iom_put( 'sophtvtr', z3dtr )
-             DO jn = 1, nptr
-               z3dtr(1,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram !  (conversion in Gg)
-               DO ji = 1, jpi
-                  z3dtr(ji,:,jn) = z3dtr(1,:,jn)
-               ENDDO
-            ENDDO
-            CALL iom_put( 'sopstvtr', z3dtr )
-         ENDIF
-         !
-         IF( iom_use( 'uocetr_vsum_cumul' ) ) THEN
-            CALL iom_get_var(  'uocetr_vsum_op', z2d ) ! get uocetr_vsum_op from xml
-            z2d(:,:) = ptr_ci_2d( z2d(:,:) )  
-            CALL iom_put( 'uocetr_vsum_cumul', z2d )
-         ENDIF
-         !
+
+            CALL dia_ptr_hst( jp_tem, 'vtr', zts(:,:,:,jp_tem) )
+            CALL dia_ptr_hst( jp_sal, 'vtr', zts(:,:,:,jp_sal) )
+
+            DEALLOCATE( zts )
+         ENDIF
       ENDIF
-      !
-      IF( ln_timing )   CALL timing_stop('dia_ptr')
-      !
-   END SUBROUTINE dia_ptr
+   END SUBROUTINE dia_ptr_zint
 
 
@@ -320,5 +456,5 @@
       !!                  ***  ROUTINE dia_ptr_init  ***
       !!                   
-      !! ** Purpose :   Initialization, namelist read
+      !! ** Purpose :   Initialization
       !!----------------------------------------------------------------------
       INTEGER ::  inum, jn           ! local integers
@@ -327,18 +463,18 @@
       !!----------------------------------------------------------------------
 
-      l_diaptr = .FALSE.
-      IF(   iom_use( 'zomsf'    ) .OR. iom_use( 'zotem'    ) .OR. iom_use( 'zosal'    ) .OR.  &
-         &  iom_use( 'zosrf'    ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR.  &
-         &  iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) .OR. iom_use( 'sophtadv' ) .OR.  &
-         &  iom_use( 'sopstadv' ) .OR. iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) .OR.  & 
-         &  iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) .OR. iom_use( 'sopstvtr' ) .OR.  &
-         &  iom_use( 'sophtvtr' ) .OR. iom_use( 'uocetr_vsum_cumul' ) )  l_diaptr  = .TRUE.
-
- 
+      ! l_diaptr is defined with iom_use
+      !   --> dia_ptr_init must be done after the call to iom_init
+      !   --> cannot be .TRUE. without cpp key: key_iom -->  nbasin define by iom_init is initialized
+      l_diaptr = iom_use( 'zomsf'    ) .OR. iom_use( 'zotem'    ) .OR. iom_use( 'zosal'    ) .OR.  &
+         &       iom_use( 'zosrf'    ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR.  &
+         &       iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) .OR. iom_use( 'sophtadv' ) .OR.  &
+         &       iom_use( 'sopstadv' ) .OR. iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) .OR.  &
+         &       iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) .OR. iom_use( 'sopstvtr' ) .OR.  &
+         &       iom_use( 'sophtvtr' ) .OR. iom_use( 'uocetr_vsum_cumul' )
+      
       IF(lwp) THEN                     ! Control print
          WRITE(numout,*)
          WRITE(numout,*) 'dia_ptr_init : poleward transport and msf initialization'
          WRITE(numout,*) '~~~~~~~~~~~~'
-         WRITE(numout,*) '   Namelist namptr : set ptr parameters'
          WRITE(numout,*) '      Poleward heat & salt transport (T) or not (F)      l_diaptr  = ', l_diaptr
       ENDIF
@@ -347,5 +483,5 @@
          !
          IF( dia_ptr_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_ptr_init : unable to allocate arrays' )
-
+         !
          rc_pwatt = rc_pwatt * rho0_rcp          ! conversion from K.s-1 to PetaWatt
          rc_ggram = rc_ggram * rho0              ! conversion from m3/s to Gg/s
@@ -354,12 +490,14 @@
 
          btmsk(:,:,1) = tmask_i(:,:)                 
-         CALL iom_open( 'subbasins', inum,  ldstop = .FALSE.  )
-         CALL iom_get( inum, jpdom_global, 'atlmsk', btmsk(:,:,2) )   ! Atlantic basin
-         CALL iom_get( inum, jpdom_global, 'pacmsk', btmsk(:,:,3) )   ! Pacific  basin
-         CALL iom_get( inum, jpdom_global, 'indmsk', btmsk(:,:,4) )   ! Indian   basin
-         CALL iom_close( inum )
-         btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) )          ! Indo-Pacific basin
-         DO jn = 2, nptr
-            btmsk(:,:,jn) = btmsk(:,:,jn) * tmask_i(:,:)               ! interior domain only
+         IF( nbasin == 5 ) THEN   ! nbasin has been initialized in iom_init to define the axis "basin"
+            CALL iom_open( 'subbasins', inum )
+            CALL iom_get( inum, jpdom_global, 'atlmsk', btmsk(:,:,2) )   ! Atlantic basin
+            CALL iom_get( inum, jpdom_global, 'pacmsk', btmsk(:,:,3) )   ! Pacific  basin
+            CALL iom_get( inum, jpdom_global, 'indmsk', btmsk(:,:,4) )   ! Indian   basin
+            CALL iom_close( inum )
+            btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) )            ! Indo-Pacific basin
+         ENDIF
+         DO jn = 2, nbasin
+            btmsk(:,:,jn) = btmsk(:,:,jn) * tmask_i(:,:)                 ! interior domain only
          END DO
          ! JD : modification so that overturning streamfunction is available in Atlantic at 34S to compare with observations
@@ -370,6 +508,6 @@
          END WHERE
          btmsk34(:,:,1) = btmsk(:,:,1)                 
-         DO jn = 2, nptr
-            btmsk34(:,:,jn) = btmsk(:,:,jn) * zmsk(:,:)               ! interior domain only
+         DO jn = 2, nbasin
+            btmsk34(:,:,jn) = btmsk(:,:,jn) * zmsk(:,:)                  ! interior domain only
          ENDDO
 
@@ -382,4 +520,6 @@
          hstr_btr(:,:,:) = 0._wp           !
          hstr_vtr(:,:,:) = 0._wp           !
+         pvtr_int(:,:,:,:) = 0._wp
+         pzon_int(:,:,:,:) = 0._wp
          !
          ll_init = .FALSE.
@@ -399,61 +539,101 @@
       INTEGER                         , INTENT(in )  :: ktra  ! tracer index
       CHARACTER(len=3)                , INTENT(in)   :: cptr  ! transport type  'adv'/'ldf'/'eiv'
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: pvflx   ! 3D input array of advection/diffusion
+      REAL(wp), DIMENSION(A2D(nn_hls),jpk)    , INTENT(in)   :: pvflx ! 3D input array of advection/diffusion
+      REAL(wp), DIMENSION(A1Dj(nn_hls),nbasin)                 :: zsj   !
       INTEGER                                        :: jn    !
 
+      DO jn = 1, nbasin
+         zsj(:,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
+      ENDDO
       !
       IF( cptr == 'adv' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_adv(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_adv(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_adv(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_adv(:,jp_sal,:), zsj(:,:) )
+      ELSE IF( cptr == 'ldf' ) THEN
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_ldf(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_ldf(:,jp_sal,:), zsj(:,:) )
+      ELSE IF( cptr == 'eiv' ) THEN
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_eiv(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_eiv(:,jp_sal,:), zsj(:,:) )
+      ELSE IF( cptr == 'vtr' ) THEN
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_vtr(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_vtr(:,jp_sal,:), zsj(:,:) )
       ENDIF
       !
-      IF( cptr == 'ldf' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_ldf(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_ldf(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
+   END SUBROUTINE dia_ptr_hst
+
+
+   SUBROUTINE ptr_sum_2d( phstr, pva )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE ptr_sum_2d ***
+      !!----------------------------------------------------------------------
+      !! ** Purpose : Add two 2D arrays with (j,nbasin) dimensions
+      !!
+      !! ** Method  : - phstr = phstr + pva
+      !!              - Call mpp_sum if the final tile
+      !!
+      !! ** Action  : phstr
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(jpj,nbasin) , INTENT(inout)         ::  phstr  !
+      REAL(wp), DIMENSION(A1Dj(nn_hls),nbasin), INTENT(in)            ::  pva    !
+      INTEGER                                               ::  jj
+#if defined key_mpp_mpi
+      INTEGER, DIMENSION(1)           ::  ish1d
+      INTEGER, DIMENSION(2)           ::  ish2d
+      REAL(wp), DIMENSION(jpj*nbasin) ::  zwork
+#endif
+
+      DO jj = ntsj, ntej
+         phstr(jj,:) = phstr(jj,:)  + pva(jj,:)
+      END DO
+
+#if defined key_mpp_mpi
+      IF( ntile == 0 .OR. ntile == nijtile ) THEN
+         ish1d(1) = jpj*nbasin
+         ish2d(1) = jpj ; ish2d(2) = nbasin
+         zwork(:) = RESHAPE( phstr(:,:), ish1d )
+         CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl )
+         phstr(:,:) = RESHAPE( zwork, ish2d )
       ENDIF
-      !
-      IF( cptr == 'eiv' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_eiv(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_eiv(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
+#endif
+   END SUBROUTINE ptr_sum_2d
+
+
+   SUBROUTINE ptr_sum_3d( phstr, pva )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE ptr_sum_3d ***
+      !!----------------------------------------------------------------------
+      !! ** Purpose : Add two 3D arrays with (j,k,nbasin) dimensions
+      !!
+      !! ** Method  : - phstr = phstr + pva
+      !!              - Call mpp_sum if the final tile
+      !!
+      !! ** Action  : phstr
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(jpj,jpk,nbasin) , INTENT(inout)     ::  phstr  !
+      REAL(wp), DIMENSION(A1Dj(nn_hls),jpk,nbasin), INTENT(in)        ::  pva    !
+      INTEGER                                               ::  jj, jk
+#if defined key_mpp_mpi
+      INTEGER, DIMENSION(1)              ::  ish1d
+      INTEGER, DIMENSION(3)              ::  ish3d
+      REAL(wp), DIMENSION(jpj*jpk*nbasin)  ::  zwork
+#endif
+
+      DO jk = 1, jpk
+         DO jj = ntsj, ntej
+            phstr(jj,jk,:) = phstr(jj,jk,:)  + pva(jj,jk,:)
+         END DO
+      END DO
+
+#if defined key_mpp_mpi
+      IF( ntile == 0 .OR. ntile == nijtile ) THEN
+         ish1d(1) = jpj*jpk*nbasin
+         ish3d(1) = jpj ; ish3d(2) = jpk ; ish3d(3) = nbasin
+         zwork(:) = RESHAPE( phstr(:,:,:), ish1d )
+         CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl )
+         phstr(:,:,:) = RESHAPE( zwork, ish3d )
       ENDIF
-      !
-      IF( cptr == 'vtr' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_vtr(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_vtr(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-      ENDIF
-      !
-   END SUBROUTINE dia_ptr_hst
+#endif
+   END SUBROUTINE ptr_sum_3d
 
 
@@ -463,15 +643,18 @@
       !!----------------------------------------------------------------------
       INTEGER               ::   dia_ptr_alloc   ! return value
-      INTEGER, DIMENSION(3) ::   ierr
+      INTEGER, DIMENSION(2) ::   ierr
       !!----------------------------------------------------------------------
       ierr(:) = 0
       !
+      ! nbasin has been initialized in iom_init to define the axis "basin"
+      !
       IF( .NOT. ALLOCATED( btmsk ) ) THEN
-         ALLOCATE( btmsk(jpi,jpj,nptr)    , btmsk34(jpi,jpj,nptr),   &
-            &      hstr_adv(jpj,jpts,nptr), hstr_eiv(jpj,jpts,nptr), &
-            &      hstr_ove(jpj,jpts,nptr), hstr_btr(jpj,jpts,nptr), &
-            &      hstr_ldf(jpj,jpts,nptr), hstr_vtr(jpj,jpts,nptr), STAT=ierr(1)  )
-            !
-         ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2))
+         ALLOCATE( btmsk(jpi,jpj,nbasin)    , btmsk34(jpi,jpj,nbasin),   &
+            &      hstr_adv(jpj,jpts,nbasin), hstr_eiv(jpj,jpts,nbasin), &
+            &      hstr_ove(jpj,jpts,nbasin), hstr_btr(jpj,jpts,nbasin), &
+            &      hstr_ldf(jpj,jpts,nbasin), hstr_vtr(jpj,jpts,nbasin), STAT=ierr(1)  )
+            !
+         ALLOCATE( pvtr_int(jpj,jpk,jpts+2,nbasin), &
+            &      pzon_int(jpj,jpk,jpts+1,nbasin), STAT=ierr(2) )
          !
          dia_ptr_alloc = MAXVAL( ierr )
@@ -493,23 +676,15 @@
       !! ** Action  : - p_fval: i-k-mean poleward flux of pvflx
       !!----------------------------------------------------------------------
-      REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk)  ::   pvflx  ! mask flux array at V-point
-      REAL(wp), INTENT(in), DIMENSION(jpi,jpj)      ::   pmsk   ! Optional 2D basin mask
+      REAL(wp), INTENT(in), DIMENSION(A2D(nn_hls),jpk)  ::   pvflx  ! mask flux array at V-point
+      REAL(wp), INTENT(in), DIMENSION(jpi,jpj)  ::   pmsk   ! Optional 2D basin mask
       !
       INTEGER                  ::   ji, jj, jk   ! dummy loop arguments
-      INTEGER                  ::   ijpj         ! ???
-      REAL(wp), POINTER, DIMENSION(:) :: p_fval  ! function value
+      REAL(wp), DIMENSION(A1Dj(nn_hls)) :: p_fval  ! function value
       !!--------------------------------------------------------------------
       !
-      p_fval => p_fval1d
-
-      ijpj = jpj
       p_fval(:) = 0._wp
       DO_3D( 0, 0, 0, 0, 1, jpkm1 )
          p_fval(jj) = p_fval(jj) + pvflx(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj)
       END_3D
-#if defined key_mpp_mpi
-      CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl)
-#endif
-      !
    END FUNCTION ptr_sj_3d
 
@@ -526,23 +701,15 @@
       !! ** Action  : - p_fval: i-k-mean poleward flux of pvflx
       !!----------------------------------------------------------------------
-      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pvflx  ! mask flux array at V-point
+      REAL(wp) , INTENT(in), DIMENSION(A2D(nn_hls))     ::   pvflx  ! mask flux array at V-point
       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pmsk   ! Optional 2D basin mask
       !
       INTEGER                  ::   ji,jj       ! dummy loop arguments
-      INTEGER                  ::   ijpj        ! ???
-      REAL(wp), POINTER, DIMENSION(:) :: p_fval ! function value
+      REAL(wp), DIMENSION(A1Dj(nn_hls)) :: p_fval ! function value
       !!--------------------------------------------------------------------
-      ! 
-      p_fval => p_fval1d
-
-      ijpj = jpj
+      !
       p_fval(:) = 0._wp
       DO_2D( 0, 0, 0, 0 )
          p_fval(jj) = p_fval(jj) + pvflx(ji,jj) * pmsk(ji,jj) * tmask_i(ji,jj)
       END_2D
-#if defined key_mpp_mpi
-      CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl )
-#endif
-      ! 
    END FUNCTION ptr_sj_2d
 
@@ -570,5 +737,4 @@
             p_fval(ji,jj) = p_fval(ji,jj-1) + pva(ji,jj) * tmask_i(ji,jj)
          END_2D
-         CALL lbc_lnk( 'diaptr', p_fval, 'U', -1.0_wp )
       END DO
       ! 
@@ -589,19 +755,11 @@
       !!
       IMPLICIT none
-      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pta    ! mask flux array at V-point
-      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)     ::   pmsk   ! Optional 2D basin mask
+      REAL(wp) , INTENT(in), DIMENSION(A2D(nn_hls),jpk) ::   pta    ! mask flux array at V-point
+      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pmsk   ! Optional 2D basin mask
       !!
       INTEGER                           :: ji, jj, jk ! dummy loop arguments
-      REAL(wp), POINTER, DIMENSION(:,:) :: p_fval     ! return function value
-#if defined key_mpp_mpi
-      INTEGER, DIMENSION(1) ::   ish
-      INTEGER, DIMENSION(2) ::   ish2
-      INTEGER               ::   ijpjjpk
-      REAL(wp), DIMENSION(jpj*jpk) ::   zwork    ! mask flux array at V-point
-#endif
+      REAL(wp), DIMENSION(A1Dj(nn_hls),jpk) :: p_fval     ! return function value
       !!--------------------------------------------------------------------
       !
-      p_fval => p_fval2d
-
       p_fval(:,:) = 0._wp
       !
@@ -609,13 +767,4 @@
          p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj)
       END_3D
-      !
-#if defined key_mpp_mpi
-      ijpjjpk = jpj*jpk
-      ish(1) = ijpjjpk  ;   ish2(1) = jpj   ;   ish2(2) = jpk
-      zwork(1:ijpjjpk) = RESHAPE( p_fval, ish )
-      CALL mpp_sum( 'diaptr', zwork, ijpjjpk, ncomm_znl )
-      p_fval(:,:) = RESHAPE( zwork, ish2 )
-#endif
-      !
    END FUNCTION ptr_sjk
 
