New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
domrea.F90 in branches/DEV_r2460_v3_3beta_NOL/NEMOGCM/NEMO/OFF_SRC – NEMO

source: branches/DEV_r2460_v3_3beta_NOL/NEMOGCM/NEMO/OFF_SRC/domrea.F90 @ 3381

Last change on this file since 3381 was 2457, checked in by cetlod, 14 years ago

Improve TOP & OFF components in v3.3beta, see ticket #774

  • Property svn:keywords set to Id
File size: 14.8 KB
Line 
1MODULE domrea
2   !!======================================================================
3   !!                       ***  MODULE domrea  ***
4   !! Ocean initialization : read the ocean domain meshmask file(s)
5   !!======================================================================
6
7   !!----------------------------------------------------------------------
8   !!   dom_rea        : read mesh and mask file(s)
9   !!                    nmsh = 1  :   mesh_mask file
10   !!                         = 2  :   mesh and mask file
11   !!                         = 3  :   mesh_hgr, mesh_zgr and mask
12   !!----------------------------------------------------------------------
13   USE dom_oce         ! ocean space and time domain
14   USE dommsk          ! domain: masks
15   USE lbclnk          !
16   USE in_out_manager  ! I/O manager
17
18   IMPLICIT NONE
19   PRIVATE
20
21   PUBLIC   dom_rea    ! routine called by inidom.F90
22   !!----------------------------------------------------------------------
23   !! NEMO/OFF 3.3 , NEMO Consortium (2010)
24   !! $Id$
25   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
26   !!----------------------------------------------------------------------
27
28CONTAINS
29
30#if   defined key_mpp_mpi   &&   defined key_dimgout
31   !!----------------------------------------------------------------------
32   !!   'key_mpp_mpi'     OR
33   !!   'key_dimgout' :         each processor makes its own direct access file
34   !!                      use build_nc_meshmask off line to retrieve
35   !!                      a ioipsl compliant meshmask file
36   !!----------------------------------------------------------------------
37#  include "domrea_dimg.h90"
38
39#else
40   !!----------------------------------------------------------------------
41   !!   Default option :                                        NetCDF file
42   !!----------------------------------------------------------------------
43
44   SUBROUTINE dom_rea
45      !!----------------------------------------------------------------------
46      !!                  ***  ROUTINE dom_rea  ***
47      !!                   
48      !! ** Purpose :  Read the NetCDF file(s) which contain(s) all the
49      !!      ocean domain informations (mesh and mask arrays). This (these)
50      !!      file(s) is (are) used for visualisation (SAXO software) and
51      !!      diagnostic computation.
52      !!
53      !! ** Method  :   Read in a file all the arrays generated in routines
54      !!      domhgr, domzgr, and dommsk. Note: the file contain depends on
55      !!      the vertical coord. used (z-coord, partial steps, s-coord)
56      !!                    nmsh = 1  :   'mesh_mask.nc' file
57      !!                         = 2  :   'mesh.nc' and mask.nc' files
58      !!                         = 3  :   'mesh_hgr.nc', 'mesh_zgr.nc' and
59      !!                                  'mask.nc' files
60      !!      For huge size domain, use option 2 or 3 depending on your
61      !!      vertical coordinate.
62      !!
63      !! ** input file :
64      !!      meshmask.nc  : domain size, horizontal grid-point position,
65      !!                     masks, depth and vertical scale factors
66      !!----------------------------------------------------------------------
67      USE iom
68      !!
69      INTEGER  ::   ji, jj, jk   ! dummy loop indices
70      INTEGER  ::   ik, inum0 , inum1 , inum2 , inum3 , inum4   ! local integers
71      REAL(wp) ::   zrefdep         ! local real
72      REAL(wp), DIMENSION(jpi,jpj) ::   zmbk   ! 2D workspace
73      !!----------------------------------------------------------------------
74
75      IF(lwp) WRITE(numout,*)
76      IF(lwp) WRITE(numout,*) 'dom_rea : read NetCDF mesh and mask information file(s)'
77      IF(lwp) WRITE(numout,*) '~~~~~~~'
78
79      zmbk(:,:) = 0._wp
80
81      SELECT CASE (nmsh)
82         !                                     ! ============================
83         CASE ( 1 )                            !  create 'mesh_mask.nc' file
84            !                                  ! ============================
85
86            IF(lwp) WRITE(numout,*) '          one file in "mesh_mask.nc" '
87            CALL iom_open( 'mesh_mask', inum0 )
88
89            inum2 = inum0                                            ! put all the informations
90            inum3 = inum0                                            ! in unit inum0
91            inum4 = inum0
92
93            !                                  ! ============================
94         CASE ( 2 )                            !  create 'mesh.nc' and
95            !                                  !         'mask.nc' files
96            !                                  ! ============================
97
98            IF(lwp) WRITE(numout,*) '          two files in "mesh.nc" and "mask.nc" '
99            CALL iom_open( 'mesh', inum1 )
100            CALL iom_open( 'mask', inum2 )
101
102            inum3 = inum1                                            ! put mesh informations
103            inum4 = inum1                                            ! in unit inum1
104
105            !                                  ! ============================
106         CASE ( 3 )                            !  create 'mesh_hgr.nc'
107            !                                  !         'mesh_zgr.nc' and
108            !                                  !         'mask.nc'     files
109            !                                  ! ============================
110
111            IF(lwp) WRITE(numout,*) '          three files in "mesh_hgr.nc" , "mesh_zgr.nc" and "mask.nc" '
112            CALL iom_open( 'mesh_hgr', inum3 ) ! create 'mesh_hgr.nc'
113            CALL iom_open( 'mesh_zgr', inum4 ) ! create 'mesh_zgr.nc'
114            CALL iom_open( 'mask'    , inum2 ) ! create 'mask.nc'
115
116         END SELECT
117
118         !                                                         ! masks (inum2)
119         CALL iom_get( inum2, jpdom_data, 'tmask', tmask )
120         CALL iom_get( inum2, jpdom_data, 'umask', umask )
121         CALL iom_get( inum2, jpdom_data, 'vmask', vmask )
122         CALL iom_get( inum2, jpdom_data, 'fmask', fmask )
123
124#if defined key_c1d
125         ! set umask and vmask equal tmask in 1D configuration
126         IF(lwp) WRITE(numout,*)
127         IF(lwp) WRITE(numout,*) '**********  1D configuration : set umask and vmask equal tmask ********'
128         IF(lwp) WRITE(numout,*) '**********                                                     ********'
129
130         umask(:,:,:) = tmask(:,:,:)
131         vmask(:,:,:) = tmask(:,:,:)
132#endif
133
134#if defined key_degrad
135         CALL iom_get( inum2, jpdom_data, 'facvolt', facvol )
136#endif
137
138         !                                                         ! horizontal mesh (inum3)
139         CALL iom_get( inum3, jpdom_data, 'glamt', glamt )
140         CALL iom_get( inum3, jpdom_data, 'glamu', glamu )
141         CALL iom_get( inum3, jpdom_data, 'glamv', glamv )
142         CALL iom_get( inum3, jpdom_data, 'glamf', glamf )
143
144         CALL iom_get( inum3, jpdom_data, 'gphit', gphit )
145         CALL iom_get( inum3, jpdom_data, 'gphiu', gphiu )
146         CALL iom_get( inum3, jpdom_data, 'gphiv', gphiv )
147         CALL iom_get( inum3, jpdom_data, 'gphif', gphif )
148
149         CALL iom_get( inum3, jpdom_data, 'e1t', e1t )
150         CALL iom_get( inum3, jpdom_data, 'e1u', e1u )
151         CALL iom_get( inum3, jpdom_data, 'e1v', e1v )
152         
153         CALL iom_get( inum3, jpdom_data, 'e2t', e2t )
154         CALL iom_get( inum3, jpdom_data, 'e2u', e2u )
155         CALL iom_get( inum3, jpdom_data, 'e2v', e2v )
156
157         CALL iom_get( inum3, jpdom_data, 'ff', ff )
158
159         CALL iom_get( inum4, jpdom_data, 'mbathy', zmbk )
160     
161         !
162         mbkt(:,:) = MAX( zmbk(:,:) * tmask(:,:,1), 1._wp )    ! bottom k-index of T-level (=1 over land)
163         
164         DO jj = 1, jpjm1                      ! bottom k-index of u- (v-) level
165            DO ji = 1, jpim1
166               mbku(ji,jj) = MIN(  mbkt(ji+1,jj  ) , mbkt(ji,jj)  )
167               mbkv(ji,jj) = MIN(  mbkt(ji  ,jj+1) , mbkt(ji,jj)  )
168            END DO
169         END DO
170         ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk
171         zmbk(:,:) = REAL( mbku(:,:), wp )   ;   CALL lbc_lnk(zmbk,'U',1.)   ;   mbku  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
172         zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
173         !                                                     ! bottom k-index of W-level = mbkt+1
174         !
175         IF( ln_sco ) THEN                                         ! s-coordinate
176            CALL iom_get( inum4, jpdom_data, 'hbatt', hbatt )
177            CALL iom_get( inum4, jpdom_data, 'hbatu', hbatu )
178            CALL iom_get( inum4, jpdom_data, 'hbatv', hbatv )
179            CALL iom_get( inum4, jpdom_data, 'hbatf', hbatf )
180           
181            CALL iom_get( inum4, jpdom_unknown, 'gsigt', gsigt ) ! scaling coef.
182            CALL iom_get( inum4, jpdom_unknown, 'gsigw', gsigw )
183            CALL iom_get( inum4, jpdom_unknown, 'gsi3w', gsi3w ) 
184            CALL iom_get( inum4, jpdom_unknown, 'esigt', esigt )
185            CALL iom_get( inum4, jpdom_unknown, 'esigw', esigw )
186
187            CALL iom_get( inum4, jpdom_data, 'e3t', e3t ) ! scale factors
188            CALL iom_get( inum4, jpdom_data, 'e3u', e3u )
189            CALL iom_get( inum4, jpdom_data, 'e3v', e3v )
190            CALL iom_get( inum4, jpdom_data, 'e3w', e3w )
191
192            CALL iom_get( inum4, jpdom_unknown, 'gdept_0', gdept_0 ) ! depth
193            CALL iom_get( inum4, jpdom_unknown, 'gdepw_0', gdepw_0 )
194         ENDIF
195
196 
197         IF( ln_zps ) THEN   
198            ! Vertical coordinates and scales factors
199            CALL iom_get( inum4, jpdom_unknown, 'gdept_0', gdept_0 ) ! depth
200            CALL iom_get( inum4, jpdom_unknown, 'gdepw_0', gdepw_0 )
201            CALL iom_get( inum4, jpdom_unknown, 'e3t_0'  , e3t_0   )
202            CALL iom_get( inum4, jpdom_unknown, 'e3w_0'  , e3w_0   )
203                                      ! z-coordinate - partial steps
204            IF( nmsh <= 6 ) THEN                                   !    ! 3D vertical scale factors
205              CALL iom_get( inum4, jpdom_data, 'e3t', e3t ) ! scale factors
206              CALL iom_get( inum4, jpdom_data, 'e3u', e3u )
207              CALL iom_get( inum4, jpdom_data, 'e3v', e3v )
208              CALL iom_get( inum4, jpdom_data, 'e3w', e3w )
209            ELSE                                                   !    ! 2D bottom scale factors
210              CALL iom_get( inum4, jpdom_data, 'e3t_ps', e3tp )
211              CALL iom_get( inum4, jpdom_data, 'e3w_ps', e3wp )
212            END IF
213
214            IF( iom_varid( inum4, 'gdept', ldstop = .FALSE. ) > 0 ) THEN
215              CALL iom_get( inum4, jpdom_data, 'gdept', gdept ) ! scale factors
216              CALL iom_get( inum4, jpdom_data, 'gdepw', gdepw )
217            ELSE                                                   !    ! 2D bottom depth
218              CALL iom_get( inum4, jpdom_data, 'hdept', hdept )   ! depth
219              CALL iom_get( inum4, jpdom_data, 'hdepw', hdepw )
220         
221              DO jk = 1,jpk
222                gdept(:,:,jk) = gdept_0(jk)
223                gdepw(:,:,jk) = gdepw_0(jk)
224              ENDDO
225
226              DO jj = 1, jpj
227                DO ji = 1, jpi
228                  ik = mbkt(ji,jj)
229                  ! ocean point only
230                  IF( ik > 0 ) THEN
231                     ! max ocean level case
232                     gdepw(ji,jj,ik+1) = hdepw(ji,jj)
233                     gdept(ji,jj,ik  ) = hdept(ji,jj)
234                     gdept(ji,jj,ik+1) = gdept(ji,jj,ik) + e3t(ji,jj,ik)
235                  ENDIF
236                END DO
237              END DO
238
239            ENDIF
240
241         ENDIF
242
243         IF( ln_zco ) THEN
244           ! Vertical coordinates and scales factors
245           CALL iom_get( inum4, jpdom_unknown, 'gdept_0', gdept_0 ) ! depth
246           CALL iom_get( inum4, jpdom_unknown, 'gdepw_0', gdepw_0 )
247           CALL iom_get( inum4, jpdom_unknown, 'e3t_0'  , e3t_0   )
248           CALL iom_get( inum4, jpdom_unknown, 'e3w_0'  , e3w_0   )
249         ENDIF
250
251!!gm BUG in s-coordinate this does not work!
252      ! deepest/shallowest W level Above/Below ~10m
253      zrefdep = 10._wp - ( 0.1_wp * MINVAL(e3w_0) )                  ! ref. depth with tolerance (10% of minimum layer thickness)
254      nlb10 = MINLOC( gdepw_0, mask = gdepw_0 > zrefdep, dim = 1 )   ! shallowest W level Below ~10m
255      nla10 = nlb10 - 1                                              ! deepest    W level Above ~10m
256!!gm end bug
257
258      ! Control printing : Grid informations (if not restart)
259      ! ----------------
260
261      IF(lwp .AND. .NOT.ln_rstart ) THEN
262         WRITE(numout,*)
263         WRITE(numout,*) '          longitude and e1 scale factors'
264         WRITE(numout,*) '          ------------------------------'
265         WRITE(numout,9300) ( ji, glamt(ji,1), glamu(ji,1),   &
266            glamv(ji,1), glamf(ji,1),   &
267            e1t(ji,1), e1u(ji,1),   &
268            e1v(ji,1), ji = 1, jpi,10)
2699300     FORMAT( 1x, i4, f8.2,1x, f8.2,1x, f8.2,1x, f8.2, 1x,    &
270            f19.10, 1x, f19.10, 1x, f19.10 )
271
272         WRITE(numout,*)
273         WRITE(numout,*) '          latitude and e2 scale factors'
274         WRITE(numout,*) '          -----------------------------'
275         WRITE(numout,9300) ( jj, gphit(1,jj), gphiu(1,jj),   &
276            &                     gphiv(1,jj), gphif(1,jj),   &
277            &                     e2t  (1,jj), e2u  (1,jj),   &
278            &                     e2v  (1,jj), jj = 1, jpj, 10 )
279      ENDIF
280
281
282      IF( nprint == 1 .AND. lwp ) THEN
283         WRITE(numout,*) '          e1u e2u '
284         CALL prihre( e1u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )
285         CALL prihre( e2u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )
286         WRITE(numout,*) '          e1v e2v  '
287         CALL prihre( e1v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )
288         CALL prihre( e2v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )
289      ENDIF
290
291      IF(lwp) THEN
292         WRITE(numout,*)
293         WRITE(numout,*) '              Reference z-coordinate depth and scale factors:'
294         WRITE(numout, "(9x,' level   gdept    gdepw     e3t      e3w  ')" )
295         WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, gdept_0(jk), gdepw_0(jk), e3t_0(jk), e3w_0(jk), jk = 1, jpk )
296      ENDIF
297
298      DO jk = 1, jpk
299         IF( e3w_0  (jk) <= 0._wp .OR. e3t_0  (jk) <= 0._wp )   CALL ctl_stop( ' e3w_0 or e3t_0 =< 0 ' )
300         IF( gdepw_0(jk) <  0._wp .OR. gdept_0(jk) <  0._wp )   CALL ctl_stop( ' gdepw_0 or gdept_0 < 0 ' )
301      END DO
302      !                                     ! ============================
303      !                                     !        close the files
304      !                                     ! ============================
305      SELECT CASE ( nmsh )
306         CASE ( 1 )               
307            CALL iom_close( inum0 )
308         CASE ( 2 )
309            CALL iom_close( inum1 )
310            CALL iom_close( inum2 )
311         CASE ( 3 )
312            CALL iom_close( inum2 )
313            CALL iom_close( inum3 )
314            CALL iom_close( inum4 )
315      END SELECT
316      !
317   END SUBROUTINE dom_rea
318
319#endif
320
321   !!======================================================================
322END MODULE domrea
Note: See TracBrowser for help on using the repository browser.