Index: branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c2d.h90
===================================================================
--- branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c2d.h90	(revision 4146)
+++ branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c2d.h90	(revision 4147)
@@ -160,107 +160,172 @@
       IF(lwp) WRITE(numout,*) '        orca ocean configuration'
 
-#if defined key_antarctic
-#     include "ldfdyn_antarctic.h90"
-#elif defined key_arctic
-#     include "ldfdyn_arctic.h90"
-#else
-      ! Read 2d integer array to specify western boundary increase in the
-      ! ===================== equatorial strip (20N-20S) defined at t-points
-
-      CALL ctl_opn( inum, 'ahmcoef', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )
-      READ(inum,9101) clexp, iim, ijm
-      READ(inum,'(/)')
-      ifreq = 40
-      il1 = 1
-      DO jn = 1, jpidta/ifreq+1
+      IF( cp_cfg == "orca" .AND. cp_cfz == "antarctic" ) THEN
+!
+! 1.2 Modify ahm
+! --------------
+         IF(lwp)WRITE(numout,*) ' inildf: Antarctic ocean'
+         IF(lwp)WRITE(numout,*) '         no tropics, no reduction of ahm'
+         IF(lwp)WRITE(numout,*) '         north boundary increase'
+
+         ahm1(:,:) = ahm0
+         ahm2(:,:) = ahm0
+
+         ijpt0=max(1,min(49 -njmpp+1,jpj))
+         ijpt1=max(0,min(49-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*2.
+            ahm1(:,jj)=ahm0*2.
+         END DO
+         ijpt0=max(1,min(48 -njmpp+1,jpj))
+         ijpt1=max(0,min(48-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.9
+            ahm1(:,jj)=ahm0*1.75
+         END DO
+         ijpt0=max(1,min(47 -njmpp+1,jpj))
+         ijpt1=max(0,min(47-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.5
+            ahm1(:,jj)=ahm0*1.25
+         END DO
+         ijpt0=max(1,min(46 -njmpp+1,jpj))
+         ijpt1=max(0,min(46-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.1
+         END DO
+
+      ELSE IF( cp_cfg == "orca" .AND. cp_cfz == "arctic" ) THEN
+! 1.2 Modify ahm 
+! --------------
+         IF(lwp)WRITE(numout,*) ' inildf: Arctic ocean'
+         IF(lwp)WRITE(numout,*) '         no tropics, no reduction of ahm'
+         IF(lwp)WRITE(numout,*) '         south and west boundary increase'
+
+
+         ahm1(:,:) = ahm0
+         ahm2(:,:) = ahm0
+
+         ijpt0=max(1,min(98-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(98-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*2.
+            ahm1(:,jj)=ahm0*2.
+         END DO
+         ijpt0=max(1,min(99-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(99-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.9
+            ahm1(:,jj)=ahm0*1.75
+         END DO
+         ijpt0=max(1,min(100-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(100-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.5
+            ahm1(:,jj)=ahm0*1.25
+         END DO
+         ijpt0=max(1,min(101-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(101-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.1
+         END DO
+      ELSE
+         ! Read 2d integer array to specify western boundary increase in the
+         ! ===================== equatorial strip (20N-20S) defined at t-points
+         
+         CALL ctl_opn( inum, 'ahmcoef', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )
+         READ(inum,9101) clexp, iim, ijm
          READ(inum,'(/)')
-         il2 = MIN( jpidta, il1+ifreq-1 )
-         READ(inum,9201) ( ii, ji = il1, il2, 5 )
-         READ(inum,'(/)')
-         DO jj = jpjdta, 1, -1
-            READ(inum,9202) ij, ( idata(ji,jj), ji = il1, il2 )
-         END DO
-         il1 = il1 + ifreq
-      END DO
-      
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            icof(ji,jj) = idata( mig(ji), mjg(jj) )
-         END DO
-      END DO
-      DO jj = nlcj+1, jpj
-         DO ji = 1, nlci
-            icof(ji,jj) = icof(ji,nlcj)
-         END DO
-      END DO
-      DO jj = 1, jpj
-         DO ji = nlci+1, jpi
-            icof(ji,jj) = icof(nlci,jj)
-         END DO
-      END DO
-
- 9101 FORMAT(1x,a15,2i8)
- 9201 FORMAT(3x,13(i3,12x))
- 9202 FORMAT(i3,41i3)
-
-
-      ! Set ahm1 and ahm2  ( T- and F- points) (used for laplacian operator)
-      ! =================
-      ! define ahm1 and ahm2 at the right grid point position
-      ! (USER: modify ahm1 and ahm2 following your desiderata)
-      
-      
-      ! Decrease ahm to zahmeq m2/s in the tropics
-      ! (from 90 to 20 degre: ahm = constant
-      ! from 20 to  2.5 degre: ahm = decrease in (1-cos)/2
-      ! from  2.5 to  0 degre: ahm = constant
-      ! symmetric in the south hemisphere)
-
-      zahmeq = aht0
-      
-      DO jj = 1, jpj
-         DO ji = 1, jpi
-            IF( ABS( gphif(ji,jj) ) >= 20. ) THEN
-               ahm2(ji,jj) =  ahm0
-            ELSEIF( ABS( gphif(ji,jj) ) <= 2.5 ) THEN
-               ahm2(ji,jj) =  zahmeq
-            ELSE
-               ahm2(ji,jj) = zahmeq + (ahm0-zahmeq)/2.   &
-                  * ( 1. - COS( rad * ( ABS(gphif(ji,jj))-2.5 ) * 180. / 17.5 ) )
-            ENDIF
-            IF( ABS( gphit(ji,jj) ) >= 20. ) THEN
-               ahm1(ji,jj) =  ahm0
-            ELSEIF( ABS( gphit(ji,jj) ) <= 2.5 ) THEN
-               ahm1(ji,jj) =  zahmeq
-            ELSE
-               ahm1(ji,jj) = zahmeq + (ahm0-zahmeq)/2.   &
-                  * ( 1. - COS( rad * ( ABS(gphit(ji,jj))-2.5 ) * 180. / 17.5 ) )
-            ENDIF
-         END DO
-      END DO
-
-      ! increase along western boundaries of equatorial strip
-      ! t-point
-      DO jj = 1, jpjm1
-         DO ji = 1, jpim1
-            zcoft = FLOAT( icof(ji,jj) ) / 100.
-            ahm1(ji,jj) = zcoft * ahm0 + (1.-zcoft) * ahm1(ji,jj) 
-         END DO
-      END DO
-      ! f-point
-      icof(:,:) = icof(:,:) * tmask(:,:,1)
-      DO jj = 1, jpjm1
-         DO ji = 1, jpim1   ! NO vector opt.
-            zmsk = tmask(ji,jj+1,1) + tmask(ji+1,jj+1,1) + tmask(ji,jj,1) + tmask(ji,jj+1,1)
-            IF( zmsk == 0. ) THEN
-               zcoff = 1.
-            ELSE
-               zcoff = FLOAT( icof(ji,jj+1) + icof(ji+1,jj+1) + icof(ji,jj) + icof(ji,jj+1) )   &
+         ifreq = 40
+         il1 = 1
+         DO jn = 1, jpidta/ifreq+1
+            READ(inum,'(/)')
+            il2 = MIN( jpidta, il1+ifreq-1 )
+            READ(inum,9201) ( ii, ji = il1, il2, 5 )
+            READ(inum,'(/)')
+            DO jj = jpjdta, 1, -1
+               READ(inum,9202) ij, ( idata(ji,jj), ji = il1, il2 )
+            END DO
+            il1 = il1 + ifreq
+         END DO
+
+         DO jj = 1, nlcj
+            DO ji = 1, nlci
+               icof(ji,jj) = idata( mig(ji), mjg(jj) )
+            END DO
+         END DO
+         DO jj = nlcj+1, jpj
+            DO ji = 1, nlci
+               icof(ji,jj) = icof(ji,nlcj)
+            END DO
+         END DO
+         DO jj = 1, jpj
+            DO ji = nlci+1, jpi
+               icof(ji,jj) = icof(nlci,jj)
+            END DO
+         END DO
+
+9101     FORMAT(1x,a15,2i8)
+9201     FORMAT(3x,13(i3,12x))
+9202     FORMAT(i3,41i3)
+
+
+         ! Set ahm1 and ahm2  ( T- and F- points) (used for laplacian operator)
+         ! =================
+         ! define ahm1 and ahm2 at the right grid point position
+         ! (USER: modify ahm1 and ahm2 following your desiderata)
+
+
+         ! Decrease ahm to zahmeq m2/s in the tropics
+         ! (from 90 to 20 degre: ahm = constant
+         ! from 20 to  2.5 degre: ahm = decrease in (1-cos)/2
+         ! from  2.5 to  0 degre: ahm = constant
+         ! symmetric in the south hemisphere)
+
+         zahmeq = aht0
+
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF( ABS( gphif(ji,jj) ) >= 20. ) THEN
+                  ahm2(ji,jj) =  ahm0
+               ELSEIF( ABS( gphif(ji,jj) ) <= 2.5 ) THEN
+                  ahm2(ji,jj) =  zahmeq
+               ELSE
+                  ahm2(ji,jj) = zahmeq + (ahm0-zahmeq)/2.   &
+                     * ( 1. - COS( rad * ( ABS(gphif(ji,jj))-2.5 ) * 180. / 17.5 ) )
+               ENDIF
+               IF( ABS( gphit(ji,jj) ) >= 20. ) THEN
+                  ahm1(ji,jj) =  ahm0
+               ELSEIF( ABS( gphit(ji,jj) ) <= 2.5 ) THEN
+                  ahm1(ji,jj) =  zahmeq
+               ELSE
+                  ahm1(ji,jj) = zahmeq + (ahm0-zahmeq)/2.   &
+                     * ( 1. - COS( rad * ( ABS(gphit(ji,jj))-2.5 ) * 180. / 17.5 ) )
+               ENDIF
+            END DO
+         END DO
+
+         ! increase along western boundaries of equatorial strip
+         ! t-point
+         DO jj = 1, jpjm1
+            DO ji = 1, jpim1
+               zcoft = FLOAT( icof(ji,jj) ) / 100.
+               ahm1(ji,jj) = zcoft * ahm0 + (1.-zcoft) * ahm1(ji,jj) 
+            END DO
+         END DO
+         ! f-point
+         icof(:,:) = icof(:,:) * tmask(:,:,1)
+         DO jj = 1, jpjm1
+            DO ji = 1, jpim1   ! NO vector opt.
+               zmsk = tmask(ji,jj+1,1) + tmask(ji+1,jj+1,1) + tmask(ji,jj,1) + tmask(ji,jj+1,1)
+               IF( zmsk == 0. ) THEN
+                  zcoff = 1.
+               ELSE
+                  zcoff = FLOAT( icof(ji,jj+1) + icof(ji+1,jj+1) + icof(ji,jj) + icof(ji,jj+1) )   &
                      / (zmsk * 100.)
-            ENDIF
-            ahm2(ji,jj) = zcoff * ahm0 + (1.-zcoff) * ahm2(ji,jj)
-         END DO
-      END DO
-#endif
+               ENDIF
+               ahm2(ji,jj) = zcoff * ahm0 + (1.-zcoff) * ahm2(ji,jj)
+            END DO
+         END DO
+      ENDIF
       
       ! Lateral boundary conditions on ( ahm1, ahm2 )
@@ -323,114 +388,180 @@
       IF(lwp) WRITE(numout,*) '        orca_r1 configuration'
 
-#if defined key_antarctic
-#     include "ldfdyn_antarctic.h90"
-#elif defined key_arctic
-#     include "ldfdyn_arctic.h90"
-#else
-      ! Read 2d integer array to specify western boundary increase in the
-      ! ===================== equatorial strip (20N-20S) defined at t-points
-
-      CALL ctl_opn( inum, 'ahmcoef', 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL',   &
-         &           1, numout, lwp )
-      REWIND inum
-      READ(inum,9101) clexp, iim, ijm
-      READ(inum,'(/)')
-      ifreq = 40
-      il1 = 1
-      DO jn = 1, jpidta/ifreq+1
+      IF( cp_cfg == "orca" .AND. cp_cfz == "antarctic" ) THEN
+!
+! 1.2 Modify ahm
+! --------------
+         IF(lwp)WRITE(numout,*) ' inildf: Antarctic ocean'
+         IF(lwp)WRITE(numout,*) '         no tropics, no reduction of ahm'
+         IF(lwp)WRITE(numout,*) '         north boundary increase'
+
+         ahm1(:,:) = ahm0
+         ahm2(:,:) = ahm0
+
+         ijpt0=max(1,min(49 -njmpp+1,jpj))
+         ijpt1=max(0,min(49-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*2.
+            ahm1(:,jj)=ahm0*2.
+         END DO
+         ijpt0=max(1,min(48 -njmpp+1,jpj))
+         ijpt1=max(0,min(48-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.9
+            ahm1(:,jj)=ahm0*1.75
+         END DO
+         ijpt0=max(1,min(47 -njmpp+1,jpj))
+         ijpt1=max(0,min(47-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.5
+            ahm1(:,jj)=ahm0*1.25
+         END DO
+         ijpt0=max(1,min(46 -njmpp+1,jpj))
+         ijpt1=max(0,min(46-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.1
+         END DO
+
+      ELSE IF( cp_cfg == "orca" .AND. cp_cfz == "arctic" ) THEN
+! 1.2 Modify ahm 
+! --------------
+         IF(lwp)WRITE(numout,*) ' inildf: Arctic ocean'
+         IF(lwp)WRITE(numout,*) '         no tropics, no reduction of ahm'
+         IF(lwp)WRITE(numout,*) '         south and west boundary increase'
+
+
+         ahm1(:,:) = ahm0
+         ahm2(:,:) = ahm0
+
+         ijpt0=max(1,min(98-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(98-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*2.
+            ahm1(:,jj)=ahm0*2.
+         END DO
+         ijpt0=max(1,min(99-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(99-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.9
+            ahm1(:,jj)=ahm0*1.75
+         END DO
+         ijpt0=max(1,min(100-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(100-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.5
+            ahm1(:,jj)=ahm0*1.25
+         END DO
+         ijpt0=max(1,min(101-jpjzoom+1-njmpp+1,jpj))
+         ijpt1=max(0,min(101-jpjzoom+1-njmpp+1,jpj-1))
+         DO jj=ijpt0,ijpt1
+            ahm2(:,jj)=ahm0*1.1
+         END DO
+      ELSE
+         
+         ! Read 2d integer array to specify western boundary increase in the
+         ! ===================== equatorial strip (20N-20S) defined at t-points
+         
+         CALL ctl_opn( inum, 'ahmcoef', 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL',   &
+            &           1, numout, lwp )
+         REWIND inum
+         READ(inum,9101) clexp, iim, ijm
          READ(inum,'(/)')
-         il2 = MIN( jpidta, il1+ifreq-1 )
-         READ(inum,9201) ( ii, ji = il1, il2, 5 )
-         READ(inum,'(/)')
-         DO jj = jpjdta, 1, -1
-            READ(inum,9202) ij, ( idata(ji,jj), ji = il1, il2 )
-         END DO
-         il1 = il1 + ifreq
-      END DO
-      
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            icof(ji,jj) = idata( mig(ji), mjg(jj) )
-         END DO
-      END DO
-      DO jj = nlcj+1, jpj
-         DO ji = 1, nlci
-            icof(ji,jj) = icof(ji,nlcj)
-         END DO
-      END DO
-      DO jj = 1, jpj
-         DO ji = nlci+1, jpi
-            icof(ji,jj) = icof(nlci,jj)
-         END DO
-      END DO
-
- 9101 FORMAT(1x,a15,2i8)
- 9201 FORMAT(3x,13(i3,12x))
- 9202 FORMAT(i3,41i3)
-
-
-      ! Set ahm1 and ahm2  ( T- and F- points) (used for laplacian operator)
-      ! =================
-      ! define ahm1 and ahm2 at the right grid point position
-      ! (USER: modify ahm1 and ahm2 following your desiderata)
-      
-      
-      ! Decrease ahm to zahmeq m2/s in the tropics
-      ! (from 90   to 20   degrees: ahm = scaled by local metrics
-      !  from 20   to  2.5 degrees: ahm = decrease in (1-cos)/2
-      !  from  2.5 to  0   degrees: ahm = constant
-      ! symmetric in the south hemisphere)
-
-      zahmeq = aht0
-      zam20s = ahm0*COS( rad * 20. )
-      
-      DO jj = 1, jpj
-         DO ji = 1, jpi
-            IF( ABS( gphif(ji,jj) ) >= 20. ) THEN
-!              leave as set in ldf_dyn_c2d
-            ELSEIF( ABS( gphif(ji,jj) ) <= 2.5 ) THEN
-               ahm2(ji,jj) =  zahmeq
-            ELSE
-               ahm2(ji,jj) =  zahmeq + (zam20s-zahmeq)/2.   &
-                  * ( 1. - COS( rad * ( ABS(gphif(ji,jj))-2.5 ) * 180. / 17.5 ) )
-            ENDIF
-            IF( ABS( gphit(ji,jj) ) >= 20. ) THEN
-!             leave as set in ldf_dyn_c2d
-            ELSEIF( ABS( gphit(ji,jj) ) <= 2.5 ) THEN
-               ahm1(ji,jj) =  zahmeq
-            ELSE
-               ahm1(ji,jj) =  zahmeq + (zam20s-zahmeq)/2.   &
-                  * ( 1. - COS( rad * ( ABS(gphit(ji,jj))-2.5 ) * 180. / 17.5 ) )
-            ENDIF
-         END DO
-      END DO
-
-      ! increase along western boundaries of equatorial strip
-      ! t-point
-      DO jj = 1, jpjm1
-         DO ji = 1, jpim1
-          IF( ABS( gphit(ji,jj) ) < 20. ) THEN
-            zcoft = FLOAT( icof(ji,jj) ) / 100.
-            ahm1(ji,jj) = zcoft * ahm0 + (1.-zcoft) * ahm1(ji,jj) 
-          ENDIF
-         END DO
-      END DO
-      ! f-point
-      icof(:,:) = icof(:,:) * tmask(:,:,1)
-      DO jj = 1, jpjm1
-         DO ji = 1, jpim1
-          IF( ABS( gphif(ji,jj) ) < 20. ) THEN
-            zmsk = tmask(ji,jj+1,1) + tmask(ji+1,jj+1,1) + tmask(ji,jj,1) + tmask(ji,jj+1,1)
-            IF( zmsk == 0. ) THEN
-               zcoff = 1.
-            ELSE
-               zcoff = FLOAT( icof(ji,jj+1) + icof(ji+1,jj+1) + icof(ji,jj) + icof(ji,jj+1) )   &
-                     / (zmsk * 100.)
-            ENDIF
-            ahm2(ji,jj) = zcoff * ahm0 + (1.-zcoff) * ahm2(ji,jj)
-          ENDIF
-         END DO
-      END DO
-#endif
+         ifreq = 40
+         il1 = 1
+         DO jn = 1, jpidta/ifreq+1
+            READ(inum,'(/)')
+            il2 = MIN( jpidta, il1+ifreq-1 )
+            READ(inum,9201) ( ii, ji = il1, il2, 5 )
+            READ(inum,'(/)')
+            DO jj = jpjdta, 1, -1
+               READ(inum,9202) ij, ( idata(ji,jj), ji = il1, il2 )
+            END DO
+            il1 = il1 + ifreq
+         END DO
+
+         DO jj = 1, nlcj
+            DO ji = 1, nlci
+               icof(ji,jj) = idata( mig(ji), mjg(jj) )
+            END DO
+         END DO
+         DO jj = nlcj+1, jpj
+            DO ji = 1, nlci
+               icof(ji,jj) = icof(ji,nlcj)
+            END DO
+         END DO
+         DO jj = 1, jpj
+            DO ji = nlci+1, jpi
+               icof(ji,jj) = icof(nlci,jj)
+            END DO
+         END DO
+
+9101     FORMAT(1x,a15,2i8)
+9201     FORMAT(3x,13(i3,12x))
+9202     FORMAT(i3,41i3)
+
+
+         ! Set ahm1 and ahm2  ( T- and F- points) (used for laplacian operator)
+         ! =================
+         ! define ahm1 and ahm2 at the right grid point position
+         ! (USER: modify ahm1 and ahm2 following your desiderata)
+
+
+         ! Decrease ahm to zahmeq m2/s in the tropics
+         ! (from 90   to 20   degrees: ahm = scaled by local metrics
+         !  from 20   to  2.5 degrees: ahm = decrease in (1-cos)/2
+         !  from  2.5 to  0   degrees: ahm = constant
+         ! symmetric in the south hemisphere)
+
+         zahmeq = aht0
+         zam20s = ahm0*COS( rad * 20. )
+
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF( ABS( gphif(ji,jj) ) >= 20. ) THEN
+                  !              leave as set in ldf_dyn_c2d
+               ELSEIF( ABS( gphif(ji,jj) ) <= 2.5 ) THEN
+                  ahm2(ji,jj) =  zahmeq
+               ELSE
+                  ahm2(ji,jj) =  zahmeq + (zam20s-zahmeq)/2.   &
+                     * ( 1. - COS( rad * ( ABS(gphif(ji,jj))-2.5 ) * 180. / 17.5 ) )
+               ENDIF
+               IF( ABS( gphit(ji,jj) ) >= 20. ) THEN
+                  !             leave as set in ldf_dyn_c2d
+               ELSEIF( ABS( gphit(ji,jj) ) <= 2.5 ) THEN
+                  ahm1(ji,jj) =  zahmeq
+               ELSE
+                  ahm1(ji,jj) =  zahmeq + (zam20s-zahmeq)/2.   &
+                     * ( 1. - COS( rad * ( ABS(gphit(ji,jj))-2.5 ) * 180. / 17.5 ) )
+               ENDIF
+            END DO
+         END DO
+
+         ! increase along western boundaries of equatorial strip
+         ! t-point
+         DO jj = 1, jpjm1
+            DO ji = 1, jpim1
+               IF( ABS( gphit(ji,jj) ) < 20. ) THEN
+                  zcoft = FLOAT( icof(ji,jj) ) / 100.
+                  ahm1(ji,jj) = zcoft * ahm0 + (1.-zcoft) * ahm1(ji,jj) 
+               ENDIF
+            END DO
+         END DO
+         ! f-point
+         icof(:,:) = icof(:,:) * tmask(:,:,1)
+         DO jj = 1, jpjm1
+            DO ji = 1, jpim1
+               IF( ABS( gphif(ji,jj) ) < 20. ) THEN
+                  zmsk = tmask(ji,jj+1,1) + tmask(ji+1,jj+1,1) + tmask(ji,jj,1) + tmask(ji,jj+1,1)
+                  IF( zmsk == 0. ) THEN
+                     zcoff = 1.
+                  ELSE
+                     zcoff = FLOAT( icof(ji,jj+1) + icof(ji+1,jj+1) + icof(ji,jj) + icof(ji,jj+1) )   &
+                        / (zmsk * 100.)
+                  ENDIF
+                  ahm2(ji,jj) = zcoff * ahm0 + (1.-zcoff) * ahm2(ji,jj)
+               ENDIF
+            END DO
+         END DO
+      ENDIF
       
       ! Lateral boundary conditions on ( ahm1, ahm2 )
