diff --git a/src/OCE/DIA/diahth.F90 b/src/OCE/DIA/diahth.F90
index c3d7f071aa77344564f91a13d970d882ff196cff..2cee430d390d650c650df5e2745bb32a5cfa2ce3 100644
--- a/src/OCE/DIA/diahth.F90
+++ b/src/OCE/DIA/diahth.F90
@@ -86,22 +86,22 @@ CONTAINS
INTEGER, INTENT( in ) :: kt ! ocean time-step index
INTEGER, INTENT( in ) :: Kmm ! ocean time level index
!!
- INTEGER :: ji, jj, jk ! dummy loop arguments
- REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth
- REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth
- REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth
- REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop
- REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace
- REAL(wp), DIMENSION(A2D(0)) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2
- REAL(wp), DIMENSION(A2D(0)) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2
- REAL(wp), DIMENSION(A2D(0)) :: zrho10_3 ! MLD: rho = rho10m + zrho3
- REAL(wp), DIMENSION(A2D(0)) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC)
- REAL(wp), DIMENSION(A2D(0)) :: ztinv ! max of temperature inversion
- REAL(wp), DIMENSION(A2D(0)) :: zdepinv ! depth of temperature inversion
- REAL(wp), DIMENSION(A2D(0)) :: zrho0_3 ! MLD rho = rho(surf) = 0.03
- REAL(wp), DIMENSION(A2D(0)) :: zrho0_1 ! MLD rho = rho(surf) = 0.01
- REAL(wp), DIMENSION(A2D(0)) :: zmaxdzT ! max of dT/dz
- REAL(wp), DIMENSION(A2D(0)) :: zdelr ! delta rho equivalent to deltaT = 0.2
+ INTEGER :: ji, jj, jk ! dummy loop arguments
+ REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth
+ REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth
+ REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth
+ REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop
+ REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace
+ REAL(wp), DIMENSION(jpi,jpj) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2
+ REAL(wp), DIMENSION(jpi,jpj) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2
+ REAL(wp), DIMENSION(jpi,jpj) :: zrho10_3 ! MLD: rho = rho10m + zrho3
+ REAL(wp), DIMENSION(jpi,jpj) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC)
+ REAL(wp), DIMENSION(jpi,jpj) :: ztinv ! max of temperature inversion
+ REAL(wp), DIMENSION(jpi,jpj) :: zdepinv ! depth of temperature inversion
+ REAL(wp), DIMENSION(jpi,jpj) :: zrho0_3 ! MLD rho = rho(surf) = 0.03
+ REAL(wp), DIMENSION(jpi,jpj) :: zrho0_1 ! MLD rho = rho(surf) = 0.01
+ REAL(wp), DIMENSION(jpi,jpj) :: zmaxdzT ! max of dT/dz
+ REAL(wp), DIMENSION(jpi,jpj) :: zdelr ! delta rho equivalent to deltaT = 0.2
!!----------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('dia_hth')
@@ -126,12 +126,10 @@ CONTAINS
!
IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN
! initialization
- DO_2D( 0, 0, 0, 0) ! loop from bottom to 2
- ztinv (ji,jj) = 0._wp
- zdepinv(ji,jj) = 0._wp
- zmaxdzT(ji,jj) = 0._wp
- END_2D
- DO_2D( 0, 0, 0, 0 )
+ ztinv (:,:) = 0._wp
+ zdepinv(:,:) = 0._wp
+ zmaxdzT(:,:) = 0._wp
+ DO_2D( 1, 1, 1, 1 )
zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
hth (ji,jj) = zztmp
zabs2 (ji,jj) = zztmp
@@ -140,7 +138,7 @@ CONTAINS
zpycn (ji,jj) = zztmp
END_2D
IF( nla10 > 1 ) THEN
- DO_2D( 0, 0, 0, 0 )
+ DO_2D( 1, 1, 1, 1 )
zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
zrho0_3(ji,jj) = zztmp
zrho0_1(ji,jj) = zztmp
@@ -149,7 +147,7 @@ CONTAINS
! Preliminary computation
! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC)
- DO_2D( 0, 0, 0, 0 )
+ DO_2D( 1, 1, 1, 1 )
IF( tmask(ji,jj,nla10) == 1. ) THEN
zu = 1779.50 + 11.250 * ts(ji,jj,nla10,jp_tem,Kmm) - 3.80 * ts(ji,jj,nla10,jp_sal,Kmm) &
& - 0.0745 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm) &
@@ -171,7 +169,7 @@ CONTAINS
! MLD: rho = rho(1) + zrho3 !
! MLD: rho = rho(1) + zrho1 !
! ------------------------------------------------------------- !
- DO_3DS( 0, 0, 0, 0, jpkm1, 2, -1 ) ! loop from bottom to 2
+ DO_3DS( 1, 1, 1, 1, jpkm1, 2, -1 ) ! loop from bottom to 2
!
zzdep = gdepw(ji,jj,jk,Kmm)
zztmp = ( ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) ) &
@@ -208,7 +206,7 @@ CONTAINS
! temperature inversion: max( 0, max of tn - tn(10m) ) !
! depth of temperature inversion !
! ------------------------------------------------------------- !
- DO_3DS( 0, 0, 0, 0, jpkm1, nlb10, -1 ) ! loop from bottom to nlb10
+ DO_3DS( 1, 1, 1, 1, jpkm1, nlb10, -1 ) ! loop from bottom to nlb10
!
zzdep = gdepw(ji,jj,jk,Kmm) * tmask(ji,jj,1)
!
@@ -300,16 +298,13 @@ CONTAINS
!
INTEGER :: ji, jj, jk, iid
REAL(wp) :: zztmp, zzdep
- INTEGER, DIMENSION(A2D(0)) :: iktem
+ INTEGER, DIMENSION(jpi,jpj) :: iktem
! --------------------------------------- !
! search deepest level above ptem !
! --------------------------------------- !
- DO_2D( 0, 0, 0, 0 )
- iktem(ji,jj) = 1
- END_2D
-
- DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! beware temperature is not always decreasing with depth => loop from top to bottom
+ iktem(:,:) = 1
+ DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! beware temperature is not always decreasing with depth => loop from top to bottom
zztmp = ts(ji,jj,jk,jp_tem,Kmm)
IF( zztmp >= ptem ) iktem(ji,jj) = jk
END_3D
@@ -317,7 +312,7 @@ CONTAINS
! ------------------------------- !
! Depth of ptem isotherm !
! ------------------------------- !
- DO_2D( 0, 0, 0, 0 )
+ DO_2D( 1, 1, 1, 1 )
!
zzdep = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) ! depth of the ocean bottom
!
@@ -344,29 +339,18 @@ CONTAINS
REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc
!
INTEGER :: ji, jj, jk, ik
- REAL(wp), DIMENSION(A2D(0)) :: zthick
- INTEGER , DIMENSION(A2D(0)) :: ilevel
+ REAL(wp), DIMENSION(jpi,jpj) :: zthick
+ INTEGER , DIMENSION(jpi,jpj) :: ilevel
! surface boundary condition
- IF( .NOT. ln_linssh ) THEN
- DO_2D( 0, 0, 0, 0 )
- zthick(ji,jj) = 0._wp
- phtc (ji,jj) = 0._wp
- END_2D
- ELSE
- DO_2D( 0, 0, 0, 0 )
- zthick(ji,jj) = ssh(ji,jj,Kmm)
- phtc (ji,jj) = pt(ji,jj,1) * ssh(ji,jj,Kmm) * tmask(ji,jj,1)
- END_2D
+ IF( .NOT. ln_linssh ) THEN ; zthick(:,:) = 0._wp ; phtc(:,:) = 0._wp
+ ELSE ; zthick(:,:) = ssh(:,:,Kmm) ; phtc(:,:) = pt(:,:,1) * ssh(:,:,Kmm) * tmask(:,:,1)
ENDIF
!
- DO_2D( 0, 0, 0, 0 )
- ilevel(ji,jj) = 1
- END_2D
- !
- DO_3D( 0, 0, 0, 0, 1, jpkm1 )
+ ilevel(:,:) = 1
+ DO_3D( 1, 1, 1, 1, 1, jpkm1 )
IF( ( gdepw(ji,jj,jk+1,Kmm) < pdep ) .AND. ( tmask(ji,jj,jk) == 1 ) ) THEN
ilevel(ji,jj) = jk+1
zthick(ji,jj) = zthick(ji,jj) + e3t(ji,jj,jk,Kmm)
@@ -374,7 +358,7 @@ CONTAINS
ENDIF
END_3D
!
- DO_2D( 0, 0, 0, 0 )
+ DO_2D( 1, 1, 1, 1 )
ik = ilevel(ji,jj)
IF( tmask(ji,jj,ik) == 1 ) THEN
zthick(ji,jj) = MIN ( gdepw(ji,jj,ik+1,Kmm), pdep ) - zthick(ji,jj) ! remaining thickness to reach dephw pdep