Index: /branches/nemo_v3_3_beta/NEMOGCM/NEMO/OFF_SRC/domrea.F90
===================================================================
--- /branches/nemo_v3_3_beta/NEMOGCM/NEMO/OFF_SRC/domrea.F90	(revision 2463)
+++ /branches/nemo_v3_3_beta/NEMOGCM/NEMO/OFF_SRC/domrea.F90	(revision 2464)
@@ -4,4 +4,6 @@
    !! Ocean initialization : read the ocean domain meshmask file(s)
    !!======================================================================
+   !! History :  3.3  ! 2010-05  (C. Ethe)  Full reorganization of the off-line
+   !!----------------------------------------------------------------------
 
    !!----------------------------------------------------------------------
@@ -13,5 +15,5 @@
    USE dom_oce         ! ocean space and time domain
    USE dommsk          ! domain: masks
-   USE lbclnk          ! 
+   USE lbclnk          ! lateral boundary condition - MPP exchanges
    USE in_out_manager  ! I/O manager
 
@@ -23,22 +25,7 @@
    !! NEMO/OFF 3.3 , NEMO Consortium (2010)
    !! $Id$
-   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
-   !!----------------------------------------------------------------------
-
+   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
+   !!----------------------------------------------------------------------
 CONTAINS
-
-#if   defined key_mpp_mpi   &&   defined key_dimgout
-   !!----------------------------------------------------------------------
-   !!   'key_mpp_mpi'     OR
-   !!   'key_dimgout' :         each processor makes its own direct access file 
-   !!                      use build_nc_meshmask off line to retrieve 
-   !!                      a ioipsl compliant meshmask file
-   !!----------------------------------------------------------------------
-#  include "domrea_dimg.h90"
-
-#else
-   !!----------------------------------------------------------------------
-   !!   Default option :                                        NetCDF file
-   !!----------------------------------------------------------------------
 
    SUBROUTINE dom_rea
@@ -157,19 +144,8 @@
          CALL iom_get( inum3, jpdom_data, 'ff', ff )
 
-         CALL iom_get( inum4, jpdom_data, 'mbathy', zmbk )
-     
-         !
-         mbkt(:,:) = MAX( zmbk(:,:) * tmask(:,:,1), 1._wp )    ! bottom k-index of T-level (=1 over land)
-         !  
-         DO jj = 1, jpjm1                      ! bottom k-index of u- (v-) level
-            DO ji = 1, jpim1
-               mbku(ji,jj) = MIN(  mbkt(ji+1,jj  ) , mbkt(ji,jj)  )
-               mbkv(ji,jj) = MIN(  mbkt(ji  ,jj+1) , mbkt(ji,jj)  )
-            END DO
-         END DO
-         ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk 
-         zmbk(:,:) = REAL( mbku(:,:), wp )   ;   CALL lbc_lnk(zmbk,'U',1.)   ;   mbku  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
-         zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
-         !                                                     ! bottom k-index of W-level = mbkt+1
+         CALL iom_get( inum4, jpdom_data, 'mbathy', zmbk )              ! number of ocean t-points
+         mbathy(:,:) = INT( zmbk(:,:) )
+         
+         CALL zgr_bot_level                                             ! mbk. arrays (deepest ocean t-, u- & v-points
          !
          IF( ln_sco ) THEN                                         ! s-coordinate
@@ -195,56 +171,79 @@
 
  
-         IF( ln_zps ) THEN   
-            ! Vertical coordinates and scales factors
+         IF( ln_zps ) THEN                                           ! z-coordinate - partial steps
+            CALL iom_get( inum4, jpdom_unknown, 'gdept_0', gdept_0 )    ! reference depth
+            CALL iom_get( inum4, jpdom_unknown, 'gdepw_0', gdepw_0 )
+            CALL iom_get( inum4, jpdom_unknown, 'e3t_0'  , e3t_0   )    ! reference scale factors
+            CALL iom_get( inum4, jpdom_unknown, 'e3w_0'  , e3w_0   )
+            !
+            IF( nmsh <= 6 ) THEN                                        ! 3D vertical scale factors
+               CALL iom_get( inum4, jpdom_data, 'e3t', e3t )
+               CALL iom_get( inum4, jpdom_data, 'e3u', e3u )
+               CALL iom_get( inum4, jpdom_data, 'e3v', e3v )
+               CALL iom_get( inum4, jpdom_data, 'e3w', e3w )
+            ELSE                                                        ! 2D bottom scale factors
+               CALL iom_get( inum4, jpdom_data, 'e3t_ps', e3tp )
+               CALL iom_get( inum4, jpdom_data, 'e3w_ps', e3wp )
+               !                                                        ! deduces the 3D scale factors
+               DO jk = 1, jpk
+                  e3t(:,:,jk) = e3t_0(jk)                                     ! set to the ref. factors
+                  e3u(:,:,jk) = e3t_0(jk)
+                  e3v(:,:,jk) = e3t_0(jk)
+                  e3w(:,:,jk) = e3w_0(jk)
+               END DO
+               DO jj = 1,jpj                                                  ! adjust the deepest values
+                  DO ji = 1,jpi
+                     ik = mbkt(ji,jj)
+                     e3t(ji,jj,ik) = e3tp(ji,jj) * tmask(ji,jj,1) + e3t_0(1) * ( 1._wp - tmask(ji,jj,1) )
+                     e3w(ji,jj,ik) = e3wp(ji,jj) * tmask(ji,jj,1) + e3w_0(1) * ( 1._wp - tmask(ji,jj,1) )
+                  END DO
+               END DO
+               DO jk = 1,jpk                         ! Computed as the minimum of neighbooring scale factors
+                  DO jj = 1, jpjm1
+                     DO ji = 1, jpim1
+                        e3u(ji,jj,jk) = MIN( e3t(ji,jj,jk), e3t(ji+1,jj,jk) )
+                        e3v(ji,jj,jk) = MIN( e3t(ji,jj,jk), e3t(ji,jj+1,jk) )
+                     END DO
+                  END DO
+               END DO
+               CALL lbc_lnk( e3u , 'U', 1._wp )   ;   CALL lbc_lnk( e3uw, 'U', 1._wp )   ! lateral boundary conditions
+               CALL lbc_lnk( e3v , 'V', 1._wp )   ;   CALL lbc_lnk( e3vw, 'V', 1._wp )
+               !
+               DO jk = 1, jpk                        ! set to z-scale factor if zero (i.e. along closed boundaries)
+                  WHERE( e3u(:,:,jk) == 0._wp )   e3u(:,:,jk) = e3t_0(jk)
+                  WHERE( e3v(:,:,jk) == 0._wp )   e3v(:,:,jk) = e3t_0(jk)
+               END DO
+            END IF
+
+            IF( iom_varid( inum4, 'gdept', ldstop = .FALSE. ) > 0 ) THEN   ! 3D depth of t- and w-level
+               CALL iom_get( inum4, jpdom_data, 'gdept', gdept )
+               CALL iom_get( inum4, jpdom_data, 'gdepw', gdepw )
+            ELSE                                                           ! 2D bottom depth
+               CALL iom_get( inum4, jpdom_data, 'hdept', hdept )
+               CALL iom_get( inum4, jpdom_data, 'hdepw', hdepw )
+               !
+               DO jk = 1, jpk                                              ! deduces the 3D depth
+                  gdept(:,:,jk) = gdept_0(jk)
+                  gdepw(:,:,jk) = gdepw_0(jk)
+               END DO
+               DO jj = 1, jpj
+                  DO ji = 1, jpi
+                     ik = mbkt(ji,jj)
+                     IF( ik > 0 ) THEN
+                        gdepw(ji,jj,ik+1) = hdepw(ji,jj)
+                        gdept(ji,jj,ik  ) = hdept(ji,jj)
+                        gdept(ji,jj,ik+1) = gdept(ji,jj,ik) + e3t(ji,jj,ik)
+                     ENDIF
+                  END DO
+               END DO
+            ENDIF
+            !
+         ENDIF
+
+         IF( ln_zco ) THEN           ! Vertical coordinates and scales factors
             CALL iom_get( inum4, jpdom_unknown, 'gdept_0', gdept_0 ) ! depth
             CALL iom_get( inum4, jpdom_unknown, 'gdepw_0', gdepw_0 )
             CALL iom_get( inum4, jpdom_unknown, 'e3t_0'  , e3t_0   )
             CALL iom_get( inum4, jpdom_unknown, 'e3w_0'  , e3w_0   )
-                                      ! z-coordinate - partial steps
-            IF( nmsh <= 6 ) THEN                                   !    ! 3D vertical scale factors
-              CALL iom_get( inum4, jpdom_data, 'e3t', e3t ) ! scale factors
-              CALL iom_get( inum4, jpdom_data, 'e3u', e3u )
-              CALL iom_get( inum4, jpdom_data, 'e3v', e3v )
-              CALL iom_get( inum4, jpdom_data, 'e3w', e3w )
-            ELSE                                                   !    ! 2D bottom scale factors
-              CALL iom_get( inum4, jpdom_data, 'e3t_ps', e3tp )
-              CALL iom_get( inum4, jpdom_data, 'e3w_ps', e3wp )
-            END IF
-
-            IF( iom_varid( inum4, 'gdept', ldstop = .FALSE. ) > 0 ) THEN
-              CALL iom_get( inum4, jpdom_data, 'gdept', gdept ) ! scale factors
-              CALL iom_get( inum4, jpdom_data, 'gdepw', gdepw )
-            ELSE                                                   !    ! 2D bottom depth
-              CALL iom_get( inum4, jpdom_data, 'hdept', hdept )   ! depth
-              CALL iom_get( inum4, jpdom_data, 'hdepw', hdepw )
-         
-              DO jk = 1,jpk
-                gdept(:,:,jk) = gdept_0(jk)
-                gdepw(:,:,jk) = gdepw_0(jk)
-              ENDDO
-
-              DO jj = 1, jpj
-                DO ji = 1, jpi
-                  ik = mbkt(ji,jj)
-                  ! ocean point only 
-                  IF( ik > 0 ) THEN
-                     ! max ocean level case
-                     gdepw(ji,jj,ik+1) = hdepw(ji,jj)
-                     gdept(ji,jj,ik  ) = hdept(ji,jj)
-                     gdept(ji,jj,ik+1) = gdept(ji,jj,ik) + e3t(ji,jj,ik)
-                  ENDIF
-                END DO
-              END DO
-
-            ENDIF
-
-         ENDIF
-
-         IF( ln_zco ) THEN
-           ! Vertical coordinates and scales factors
-           CALL iom_get( inum4, jpdom_unknown, 'gdept_0', gdept_0 ) ! depth
-           CALL iom_get( inum4, jpdom_unknown, 'gdepw_0', gdepw_0 )
-           CALL iom_get( inum4, jpdom_unknown, 'e3t_0'  , e3t_0   )
-           CALL iom_get( inum4, jpdom_unknown, 'e3w_0'  , e3w_0   )
          ENDIF
 
@@ -317,5 +316,36 @@
    END SUBROUTINE dom_rea
 
-#endif
+
+   SUBROUTINE zgr_bot_level
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE zgr_bot_level  ***
+      !!
+      !! ** Purpose :   defines the vertical index of ocean bottom (mbk. arrays)
+      !!
+      !! ** Method  :   computes from mbathy with a minimum value of 1 over land
+      !!
+      !! ** Action  : - update mbathy: level bathymetry (in level index)
+      !!----------------------------------------------------------------------
+      INTEGER ::   ji, jj   ! dummy loop indices
+      REAL(wp), DIMENSION(jpi,jpj) ::   zmbk   ! 2D workspace 
+      !!----------------------------------------------------------------------
+      !
+      IF(lwp) WRITE(numout,*)
+      IF(lwp) WRITE(numout,*) '    zgr_bot_level : ocean bottom k-index of T-, U-, V- and W-levels '
+      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~'
+      !
+      mbkt(:,:) = MAX( mbathy(:,:) , 1 )    ! bottom k-index of T-level (=1 over land)
+      !
+      DO jj = 1, jpjm1                      ! bottom k-index of u- (v-) level
+         DO ji = 1, jpim1
+            mbku(ji,jj) = MIN(  mbkt(ji+1,jj  ) , mbkt(ji,jj)  )
+            mbkv(ji,jj) = MIN(  mbkt(ji  ,jj+1) , mbkt(ji,jj)  )
+         END DO
+      END DO
+      ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk 
+      zmbk(:,:) = REAL( mbku(:,:), wp )   ;   CALL lbc_lnk(zmbk,'U',1.)   ;   mbku  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
+      zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
+      !
+   END SUBROUTINE zgr_bot_level
 
    !!======================================================================
Index: /branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/sedini.F90
===================================================================
--- /branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/sedini.F90	(revision 2463)
+++ /branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/sedini.F90	(revision 2464)
@@ -123,5 +123,5 @@
       DO jj = 1, jpj
          DO ji = 1, jpi
-            ikt = MAX( INT( sbathy(ji,jj) )  - 1, 1 )
+            ikt = MAX( INT( sbathy(ji,jj) ), 1 )
             IF( tmask(ji,jj,ikt) == 1 ) epkbot(ji,jj) = zdta(ji,jj)
          ENDDO
