Index: trunk/NEMO/OPA_SRC/DOM/domhgr.F90
===================================================================
--- trunk/NEMO/OPA_SRC/DOM/domhgr.F90	(revision 434)
+++ trunk/NEMO/OPA_SRC/DOM/domhgr.F90	(revision 473)
@@ -4,4 +4,17 @@
    !! Ocean initialization : domain initialization
    !!==============================================================================
+   !! History :       !  88-03  (G. Madec)
+   !!                 !  91-11  (G. Madec)
+   !!                 !  92-06  (M. Imbard)
+   !!                 !  96-01  (G. Madec)  terrain following coordinates
+   !!                 !  97-02  (G. Madec)  print mesh informations
+   !!                 !  99-11  (M. Imbard) NetCDF format with IO-IPSL
+   !!                 !  00-08  (D. Ludicone) Reduced section at Bab el Mandeb
+   !!                 !  01-09  (M. Levy)  eel config: grid in km, beta-plane
+   !!            8.5  !  02-08  (G. Madec)  F90: Free form and module, namelist
+   !!            9.0  !  04-01  (A.M. Treguier, J.M. Molines) Case 4 (Mercator mesh)
+   !!                           use of parameters in par_CONFIG-Rxx.h90, not in namelist
+   !!                 !  04-05  (A. Koch-Larrouy) Add Gyre configuration 
+   !!----------------------------------------------------------------------
 
    !!----------------------------------------------------------------------
@@ -27,5 +40,5 @@
    !!   OPA 9.0 , LOCEAN-IPSL (2005) 
    !! $Header$ 
-   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
+   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
 
@@ -84,20 +97,7 @@
       !!        define ff: coriolis factor at f-point
       !!
-      !! References :
-      !!      Marti, Madec and Delecluse, 1992, j. geophys. res., in press.
-      !!
-      !! History :
-      !!        !  88-03  (G. Madec)
-      !!        !  91-11  (G. Madec)
-      !!        !  92-06  (M. Imbard)
-      !!        !  96-01  (G. Madec)  terrain following coordinates
-      !!        !  97-02  (G. Madec)  print mesh informations
-      !!        !  01-09  (M. Levy)  eel config: grid in km, beta-plane
-      !!   8.5  !  02-08  (G. Madec)  F90: Free form and module, namelist
-      !!   9.0  !  04-01  (A.M. Treguier, J.M. Molines) Case 4 (Mercator mesh)
-      !!                  use of parameters in par_CONFIG-Rxx.h90, not in namelist
-      !!        !  04-05  (A. Koch-Larrouy) Add Gyre configuration 
+      !! References :   Marti, Madec and Delecluse, 1992, JGR
+      !!                Madec, Imbard, 1996, Clim. Dyn.
       !!----------------------------------------------------------------------
-      !! * local declarations
       INTEGER  ::   ji, jj              ! dummy loop indices
       INTEGER  ::   ii0, ii1, ij0, ij1  ! temporary integers
@@ -164,4 +164,35 @@
             IF(lwp) WRITE(numout,*)
             IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Gibraltar Strait'
+            !
+            ii0 = 627   ;   ii1 = 628        ! Bosphore Strait (e2u = 10 km)
+            ij0 = 343   ;   ij1 = 343   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3
+            IF(lwp) WRITE(numout,*)
+            IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Bosphore Strait'
+            !
+            ii0 =  93   ;   ii1 =  94        ! Sumba Strait (e2u = 40 km)
+            ij0 = 232   ;   ij1 = 232   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  40.e3
+            IF(lwp) WRITE(numout,*)
+            IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Sumba Strait'
+            !
+            ii0 = 103   ;   ii1 = 103        ! Ombai Strait (e2u = 15 km)
+            ij0 = 232   ;   ij1 = 232   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  15.e3
+            IF(lwp) WRITE(numout,*)
+            IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Ombai Strait'
+            !
+            ii0 =  15   ;   ii1 =  15        ! Palk Strait (e2u = 10 km)
+            ij0 = 270   ;   ij1 = 270   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3
+            IF(lwp) WRITE(numout,*)
+            IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e2u at the Palk Strait'
+            !
+            ii0 =  87   ;   ii1 =  87        ! Lombok Strait (e1v = 10 km)
+            ij0 = 232   ;   ij1 = 233   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3
+            IF(lwp) WRITE(numout,*)
+            IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e1v at the Lombok Strait'
+            !
+            !
+            ii0 = 662   ;   ii1 = 662        ! Bab el Mandeb (e1v = 25 km)
+            ij0 = 276   ;   ij1 = 276   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  25.e3
+            IF(lwp) WRITE(numout,*)
+            IF(lwp) WRITE(numout,*) '             orca_r05: Reduced e1v at the Bab el Mandeb'
             !
          ENDIF
@@ -269,9 +300,5 @@
          IF(lwp) WRITE(numout,*) '          geographical mesh on the sphere, MERCATOR type'
          IF(lwp) WRITE(numout,*) '          longitudinal/latitudinal spacing given by ppe1_deg'
-         IF ( ppgphi0 == -90 ) THEN
-                IF(lwp) WRITE(numout,*) ' Mercator grid cannot start at south pole !!!! '
-                IF(lwp) WRITE(numout,*) ' We stop '
-                STOP
-         ENDIF
+         IF ( ppgphi0 == -90 ) CALL ctl_stop( ' Mercator grid cannot start at south pole !!!! ' )
 
          !  Find index corresponding to the equator, given the grid spacing e1_deg
@@ -368,7 +395,6 @@
 
       CASE DEFAULT
-         IF(lwp) WRITE(numout,cform_err)
-         IF(lwp) WRITE(numout,*) '          bad flag value for jphgr_msh = ', jphgr_msh
-         nstop = nstop + 1
+         WRITE(ctmp1,*) '          bad flag value for jphgr_msh = ', jphgr_msh
+         CALL ctl_stop( ctmp1 )
 
       END SELECT
@@ -480,9 +506,5 @@
       IF( nperio == 2 ) THEN
          znorme = SQRT( SUM( gphiu(:,2) * gphiu(:,2) ) ) / FLOAT( jpi )
-         IF( znorme > 1.e-13 ) THEN
-            IF(lwp) WRITE(numout,cform_err)
-            IF(lwp) WRITE(numout,*) ' ===>>>> : symmetrical condition: rerun with good equator line'
-            nstop = nstop + 1
-         ENDIF
+         IF( znorme > 1.e-13 ) CALL ctl_stop( ' ===>>>> : symmetrical condition: rerun with good equator line' )
       ENDIF
 
@@ -499,183 +521,42 @@
       !!      or semi-analytical method. It is read in a NetCDF file. 
       !!     
-      !! References :
-      !!      Marti, Madec and Delecluse, 1992, JGR, 97, 12,763-12,766.
-      !!      Madec, Imbard, 1996, Clim. Dyn., 12, 381-388.
-      !!
-      !! History :
-      !!        !         (O. Marti)  Original code
-      !!        !  91-03  (G. Madec)
-      !!        !  92-07  (M. Imbard)
-      !!        !  99-11  (M. Imbard) NetCDF format with IOIPSL
-      !!        !  00-08  (D. Ludicone) Reduced section at Bab el Mandeb
-      !!   8.5  !  02-06  (G. Madec)  F90: Free form
       !!----------------------------------------------------------------------
-      !! * Modules used
-      USE ioipsl
-
-      !! * Local declarations
-      LOGICAL ::   llog = .FALSE.
-      CHARACTER(len=21) ::   clname
-      INTEGER  ::   ji, jj              ! dummy loop indices
-      INTEGER  ::   inum                ! temporary logical unit
-      INTEGER  ::   ilev, itime         ! temporary integers
-      REAL(wp) ::   zdt, zdate0         ! temporary scalars
-      REAL(wp) ::   zdept(1)            ! temporary workspace
-      REAL(wp), DIMENSION(jpidta,jpjdta) ::   &
-         zlamt, zphit, zdta             ! temporary workspace (NetCDF read)
+      USE iom
+
+      INTEGER ::   inum   ! temporary logical unit
       !!----------------------------------------------------------------------
-      clname = 'coordinates'
-#if defined key_agrif
-      if ( .NOT. Agrif_Root() ) then
-         clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
-      endif
-#endif         
-
-
-      ! 1. Read of the grid coordinates and scale factors
-      ! -------------------------------------------------
 
       IF(lwp) THEN
          WRITE(numout,*)
          WRITE(numout,*) 'hgr_read : read the horizontal coordinates'
-         WRITE(numout,*) '~~~~~~~~~~~      jpiglo = ', jpiglo, ' jpjglo = ', jpjglo, ' jpk = ', jpk
+         WRITE(numout,*) '~~~~~~~~      jpiglo = ', jpiglo, ' jpjglo = ', jpjglo, ' jpk = ', jpk
       ENDIF
-
-      ! read the file
-      itime = 0
-      ilev = 1   
-      zlamt(:,:) = 0.e0
-      zphit(:,:) = 0.e0
-      CALL restini( clname, jpidta, jpjdta, zlamt , zphit,   &
-         &                  ilev  , zdept , 'NONE',          &
-         &                  itime , zdate0, zdt   , inum, domain_id=nidom )
-
-      CALL restget( inum, 'glamt', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            glamt(ji,jj) = zdta(mig(ji),mjg(jj))
-         END DO
-      END DO
-      CALL restget( inum, 'glamu', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            glamu(ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'glamv', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            glamv(ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'glamf', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            glamf(ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'gphit', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            gphit(ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'gphiu', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            gphiu(ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'gphiv', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            gphiv(ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'gphif', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            gphif(ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e1t', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e1t  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e1u', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e1u  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e1v', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e1v  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e1f', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e1f  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e2t', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e2t  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e2u', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e2u  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e2v', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e2v  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-      CALL restget( inum, 'e2f', jpidta, jpjdta, 1, itime, llog, zdta )
-      DO jj = 1, nlcj
-         DO ji = 1, nlci
-            e2f  (ji,jj) = zdta(mig(ji),mjg(jj))                    
-         END DO
-      END DO
-
-      CALL restclo( inum )
-
-      ! set extra rows add in mpp to none zero values
-      DO jj = nlcj+1, jpj
-         DO ji = 1, nlci
-            glamt(ji,jj) = glamt(ji,1)   ;   gphit(ji,jj) = gphit(ji,1)
-            glamu(ji,jj) = glamu(ji,1)   ;   gphiu(ji,jj) = gphiu(ji,1)
-            glamv(ji,jj) = glamv(ji,1)   ;   gphiv(ji,jj) = gphiv(ji,1)
-            glamf(ji,jj) = glamf(ji,1)   ;   gphif(ji,jj) = gphif(ji,1)
-            e1t  (ji,jj) = e1t  (ji,1)   ;   e2t  (ji,jj) = e2t  (ji,1)
-            e1u  (ji,jj) = e1u  (ji,1)   ;   e2u  (ji,jj) = e2u  (ji,1)
-            e1v  (ji,jj) = e1v  (ji,1)   ;   e2v  (ji,jj) = e2v  (ji,1)
-            e1f  (ji,jj) = e1f  (ji,1)   ;   e2f  (ji,jj) = e2f  (ji,1)
-         END DO
-      END DO
-
-      ! set extra columns add in mpp to none zero values
-      DO ji = nlci+1, jpi
-         glamt(ji,:) = glamt(1,:)   ;   gphit(ji,:) = gphit(1,:)
-         glamu(ji,:) = glamu(1,:)   ;   gphiu(ji,:) = gphiu(1,:)
-         glamv(ji,:) = glamv(1,:)   ;   gphiv(ji,:) = gphiv(1,:)
-         glamf(ji,:) = glamf(1,:)   ;   gphif(ji,:) = gphif(1,:)
-         e1t  (ji,:) = e1t  (1,:)   ;   e2t  (ji,:) = e2t  (1,:)
-         e1u  (ji,:) = e1u  (1,:)   ;   e2u  (ji,:) = e2u  (1,:)
-         e1v  (ji,:) = e1v  (1,:)   ;   e2v  (ji,:) = e2v  (1,:)
-         e1f  (ji,:) = e1f  (1,:)   ;   e2f  (ji,:) = e2f  (1,:)
-      END DO
-
-   END SUBROUTINE hgr_read
-
+      
+      CALL iom_open( 'coordinates', inum )
+      
+      CALL iom_get( inum, jpdom_data, 'glamt', glamt )
+      CALL iom_get( inum, jpdom_data, 'glamu', glamu )
+      CALL iom_get( inum, jpdom_data, 'glamv', glamv )
+      CALL iom_get( inum, jpdom_data, 'glamf', glamf )
+      
+      CALL iom_get( inum, jpdom_data, 'gphit', gphit )
+      CALL iom_get( inum, jpdom_data, 'gphiu', gphiu )
+      CALL iom_get( inum, jpdom_data, 'gphiv', gphiv )
+      CALL iom_get( inum, jpdom_data, 'gphif', gphif )
+      
+      CALL iom_get( inum, jpdom_data, 'e1t', e1t )
+      CALL iom_get( inum, jpdom_data, 'e1u', e1u )
+      CALL iom_get( inum, jpdom_data, 'e1v', e1v )
+      CALL iom_get( inum, jpdom_data, 'e1f', e1f )
+      
+      CALL iom_get( inum, jpdom_data, 'e2t', e2t )
+      CALL iom_get( inum, jpdom_data, 'e2u', e2u )
+      CALL iom_get( inum, jpdom_data, 'e2v', e2v )
+      CALL iom_get( inum, jpdom_data, 'e2f', e2f )
+      
+      CALL iom_close( inum )
+      
+    END SUBROUTINE hgr_read
+    
    !!======================================================================
 END MODULE domhgr
