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.
domwri.F90 in trunk/NEMO/OPA_SRC/DOM – NEMO

source: trunk/NEMO/OPA_SRC/DOM/domwri.F90 @ 1161

Last change on this file since 1161 was 1161, checked in by smasson, 16 years ago

improve mesh_mask file, see ticket:228
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 13.2 KB
Line 
1MODULE domwri
2   !!======================================================================
3   !!                       ***  MODULE domwri  ***
4   !! Ocean initialization : write the ocean domain mesh ask file(s)
5   !!======================================================================
6
7   !!----------------------------------------------------------------------
8   !!   dom_wri        : create and write 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   !! * Modules used
14   USE dom_oce         ! ocean space and time domain
15   USE in_out_manager
16   USE iom
17   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
18   USE lib_mpp
19
20   IMPLICIT NONE
21   PRIVATE
22
23   !! * Accessibility
24   PUBLIC dom_wri        ! routine called by inidom.F90
25   !!----------------------------------------------------------------------
26   !!   OPA 9.0 , LOCEAN-IPSL (2005)
27   !! $Id$
28   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
29   !!----------------------------------------------------------------------
30
31CONTAINS
32
33   SUBROUTINE dom_wri
34      !!----------------------------------------------------------------------
35      !!                  ***  ROUTINE dom_wri  ***
36      !!                   
37      !! ** Purpose :   Create the NetCDF file(s) which contain(s) all the
38      !!      ocean domain informations (mesh and mask arrays). This (these)
39      !!      file(s) is (are) used for visualisation (SAXO software) and
40      !!      diagnostic computation.
41      !!
42      !! ** Method  :   Write in a file all the arrays generated in routines
43      !!      domhgr, domzgr, and dommsk. Note: the file contain depends on
44      !!      the vertical coord. used (z-coord, partial steps, s-coord)
45      !!                    nmsh = 1  :   'mesh_mask.nc' file
46      !!                         = 2  :   'mesh.nc' and mask.nc' files
47      !!                         = 3  :   'mesh_hgr.nc', 'mesh_zgr.nc' and
48      !!                                  'mask.nc' files
49      !!      For huge size domain, use option 2 or 3 depending on your
50      !!      vertical coordinate.
51      !!
52      !! ** output file :
53      !!      meshmask.nc  : domain size, horizontal grid-point position,
54      !!                     masks, depth and vertical scale factors
55      !!
56      !! History :
57      !!        !  97-02  (G. Madec)  Original code
58      !!        !  99-11  (M. Imbard)  NetCDF FORMAT with IOIPSL
59      !!   9.0  !  02-08  (G. Madec)  F90 and several file
60      !!----------------------------------------------------------------------
61      INTEGER  ::   inum0   ! temprary units for 'mesh_mask.nc' file
62      INTEGER  ::   inum1   ! temprary units for 'mesh.nc'      file
63      INTEGER  ::   inum2   ! temprary units for 'mask.nc'      file
64      INTEGER  ::   inum3   ! temprary units for 'mesh_hgr.nc'  file
65      INTEGER  ::   inum4   ! temprary units for 'mesh_zgr.nc'  file
66      REAL(wp), DIMENSION(jpi,jpj) ::    zprt   ! temporary array for bathymetry
67      CHARACTER (len=21) ::   clnam0   ! filename (mesh and mask informations)
68      CHARACTER (len=21) ::   clnam1   ! filename (mesh informations)
69      CHARACTER (len=21) ::   clnam2   ! filename (mask informations)
70      CHARACTER (len=21) ::   clnam3   ! filename (horizontal mesh informations)
71      CHARACTER (len=21) ::   clnam4   ! filename (vertical   mesh informations)
72      !!----------------------------------------------------------------------
73
74       IF(lwp) WRITE(numout,*)
75       IF(lwp) WRITE(numout,*) 'dom_wri : create NetCDF mesh and mask information file(s)'
76       IF(lwp) WRITE(numout,*) '~~~~~~~'
77
78       clnam0 = 'mesh_mask'  ! filename (mesh and mask informations)
79       clnam1 = 'mesh'       ! filename (mesh informations)
80       clnam2 = 'mask'       ! filename (mask informations)
81       clnam3 = 'mesh_hgr'   ! filename (horizontal mesh informations)
82       clnam4 = 'mesh_zgr'   ! filename (vertical   mesh informations)
83
84      SELECT CASE (nmsh)
85         !                                  ! ============================
86      CASE ( 1 )                            !  create 'mesh_mask.nc' file
87         !                                  ! ============================
88         CALL iom_open( TRIM(clnam0), inum0, ldwrt = .TRUE., kiolib = jprstlib )
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         CALL iom_open( TRIM(clnam1), inum1, ldwrt = .TRUE., kiolib = jprstlib )
98         CALL iom_open( TRIM(clnam2), inum2, ldwrt = .TRUE., kiolib = jprstlib )
99         inum3 = inum1                                            ! put mesh informations
100         inum4 = inum1                                            ! in unit inum1
101         !                                  ! ============================
102      CASE ( 3 )                            !  create 'mesh_hgr.nc'
103         !                                  !         'mesh_zgr.nc' and
104         !                                  !         'mask.nc'     files
105         !                                  ! ============================
106         CALL iom_open( TRIM(clnam2), inum2, ldwrt = .TRUE., kiolib = jprstlib )
107         CALL iom_open( TRIM(clnam3), inum3, ldwrt = .TRUE., kiolib = jprstlib )
108         CALL iom_open( TRIM(clnam4), inum4, ldwrt = .TRUE., kiolib = jprstlib )
109         
110      END SELECT
111     
112      !                                                         ! masks (inum2)
113      CALL iom_rstput( 0, 0, inum2, 'tmask', tmask, ktype = jp_i1 )     !    ! land-sea mask
114      CALL iom_rstput( 0, 0, inum2, 'umask', umask, ktype = jp_i1 )
115      CALL iom_rstput( 0, 0, inum2, 'vmask', vmask, ktype = jp_i1 )
116      CALL iom_rstput( 0, 0, inum2, 'fmask', fmask, ktype = jp_i1 )
117     
118     
119      zprt = tmask(:,:,1) * dom_uniq('T')                               !    ! unique point mask
120      CALL iom_rstput( 0, 0, inum2, 'tmaskutil', zprt, ktype = jp_i1 ) 
121      zprt = umask(:,:,1) * dom_uniq('U')
122      CALL iom_rstput( 0, 0, inum2, 'umaskutil', zprt, ktype = jp_i1 ) 
123      zprt = vmask(:,:,1) * dom_uniq('V')
124      CALL iom_rstput( 0, 0, inum2, 'vmaskutil', zprt, ktype = jp_i1 ) 
125      zprt = fmask(:,:,1) * dom_uniq('F')
126      CALL iom_rstput( 0, 0, inum2, 'fmaskutil', zprt, ktype = jp_i1 ) 
127
128      !                                                         ! horizontal mesh (inum3)
129      CALL iom_rstput( 0, 0, inum3, 'glamt', glamt, ktype = jp_r4 )     !    ! latitude
130      CALL iom_rstput( 0, 0, inum3, 'glamu', glamu, ktype = jp_r4 )
131      CALL iom_rstput( 0, 0, inum3, 'glamv', glamv, ktype = jp_r4 )
132      CALL iom_rstput( 0, 0, inum3, 'glamf', glamf, ktype = jp_r4 )
133     
134      CALL iom_rstput( 0, 0, inum3, 'gphit', gphit, ktype = jp_r4 )     !    ! longitude
135      CALL iom_rstput( 0, 0, inum3, 'gphiu', gphiu, ktype = jp_r4 )
136      CALL iom_rstput( 0, 0, inum3, 'gphiv', gphiv, ktype = jp_r4 )
137      CALL iom_rstput( 0, 0, inum3, 'gphif', gphif, ktype = jp_r4 )
138     
139      CALL iom_rstput( 0, 0, inum3, 'e1t', e1t, ktype = jp_r8 )         !    ! e1 scale factors
140      CALL iom_rstput( 0, 0, inum3, 'e1u', e1u, ktype = jp_r8 )
141      CALL iom_rstput( 0, 0, inum3, 'e1v', e1v, ktype = jp_r8 )
142      CALL iom_rstput( 0, 0, inum3, 'e1f', e1f, ktype = jp_r8 )
143     
144      CALL iom_rstput( 0, 0, inum3, 'e2t', e2t, ktype = jp_r8 )         !    ! e2 scale factors
145      CALL iom_rstput( 0, 0, inum3, 'e2u', e2u, ktype = jp_r8 )
146      CALL iom_rstput( 0, 0, inum3, 'e2v', e2v, ktype = jp_r8 )
147      CALL iom_rstput( 0, 0, inum3, 'e2f', e2f, ktype = jp_r8 )
148     
149      CALL iom_rstput( 0, 0, inum3, 'ff', ff, ktype = jp_r8 )           !    ! coriolis factor
150     
151!       note that mbathy has been modified in dommsk or in solver.
152!       it is the number of non-zero "w" levels in the water, and the minimum
153!       value (on land) is 2. We define zprt as the number of "T" points in the ocean
154!       at any location, and zero on land.
155!
156      zprt = tmask(:,:,1)*(mbathy-1)
157      CALL iom_rstput( 0, 0, inum4, 'mbathy', zprt, ktype = jp_i2 )
158           
159#if ! defined key_zco
160      IF( ln_sco ) THEN                                         ! s-coordinate
161         CALL iom_rstput( 0, 0, inum4, 'hbatt', hbatt )         !    ! depth
162         CALL iom_rstput( 0, 0, inum4, 'hbatu', hbatu ) 
163         CALL iom_rstput( 0, 0, inum4, 'hbatv', hbatv )
164         CALL iom_rstput( 0, 0, inum4, 'hbatf', hbatf )
165         
166         CALL iom_rstput( 0, 0, inum4, 'gsigt', gsigt )         !    ! scaling coef.
167         CALL iom_rstput( 0, 0, inum4, 'gsigw', gsigw ) 
168         CALL iom_rstput( 0, 0, inum4, 'gsi3w', gsi3w )
169         CALL iom_rstput( 0, 0, inum4, 'esigt', esigt )
170         CALL iom_rstput( 0, 0, inum4, 'esigw', esigw )
171         
172         CALL iom_rstput( 0, 0, inum4, 'e3t', e3t )             !    ! scale factors
173         CALL iom_rstput( 0, 0, inum4, 'e3u', e3u )
174         CALL iom_rstput( 0, 0, inum4, 'e3v', e3v )
175         CALL iom_rstput( 0, 0, inum4, 'e3w', e3w )
176         
177         CALL iom_rstput( 0, 0, inum4, 'gdept_0' , gdept_0 )    !    ! stretched system
178         CALL iom_rstput( 0, 0, inum4, 'gdepw_0' , gdepw_0 )
179      ENDIF
180     
181      IF( ln_zps ) THEN                                         ! z-coordinate - partial steps
182         CALL iom_rstput( 0, 0, inum4, 'hdept' , hdept  )       !    ! bottom depth
183         CALL iom_rstput( 0, 0, inum4, 'hdepw' , hdepw  ) 
184         
185         CALL iom_rstput( 0, 0, inum4, 'e3t_ps', e3tp   )       !    ! bottom scale factors
186         CALL iom_rstput( 0, 0, inum4, 'e3w_ps', e3wp   ) 
187   
188         CALL iom_rstput( 0, 0, inum4, 'gdept_0', gdept_0 )     !    ! reference z-coord.
189         CALL iom_rstput( 0, 0, inum4, 'gdepw_0', gdepw_0 )
190         CALL iom_rstput( 0, 0, inum4, 'e3t_0'  , e3t_0   )
191         CALL iom_rstput( 0, 0, inum4, 'e3w_0'  , e3w_0   )
192      ENDIF
193     
194#endif
195     
196      IF( ln_zco ) THEN
197         !                                                      ! z-coordinate - full steps
198         CALL iom_rstput( 0, 0, inum4, 'gdept_0', gdept_0 )     !    ! depth
199         CALL iom_rstput( 0, 0, inum4, 'gdepw_0', gdepw_0 )
200         CALL iom_rstput( 0, 0, inum4, 'e3t_0'  , e3t_0   )     !    ! scale factors
201         CALL iom_rstput( 0, 0, inum4, 'e3w_0'  , e3w_0   )
202      ENDIF
203      !                                     ! ============================
204      !                                     !        close the files
205      !                                     ! ============================
206      SELECT CASE ( nmsh )
207      CASE ( 1 )               
208         CALL iom_close( inum0 )
209      CASE ( 2 )
210         CALL iom_close( inum1 )
211         CALL iom_close( inum2 )
212      CASE ( 3 )
213         CALL iom_close( inum2 )
214         CALL iom_close( inum3 )
215         CALL iom_close( inum4 )
216      END SELECT
217     
218   END SUBROUTINE dom_wri
219
220
221   FUNCTION dom_uniq( cdgrd )   RESULT( puniq )
222      !!----------------------------------------------------------------------
223      !!                  ***  ROUTINE dom_uniq  ***
224      !!                   
225      !! ** Purpose :   identify unique point of a grid (TUVF)
226      !!
227      !! ** Method  :   1) aplly lbc_lnk on an array with different values for each element
228      !!                2) check which elements have been changed
229      !!
230      !!----------------------------------------------------------------------
231      CHARACTER(len=1)            , INTENT(in   ) ::  cdgrd   !
232      REAL(wp), DIMENSION(jpi,jpj)                ::  puniq   !
233
234      REAL(wp), DIMENSION(jpi,jpj  ) ::  ztstref   ! array with different values for each element
235      REAL(wp)                       ::  zshift    ! shift value link to the process number
236      LOGICAL , DIMENSION(jpi,jpj,1) ::  lldbl     ! is the point unique or not?
237      INTEGER                        ::  ji        ! dummy loop indices
238      !!----------------------------------------------------------------------
239
240      ! build an array with different values for each element
241      ! in mpp: make sure that these values are different even between process
242      ! -> apply a shift value according to the process number
243      zshift = jpi * jpj * ( narea - 1 )
244      ztstref(:,:) = RESHAPE( (/ (zshift + REAL(ji, wp), ji = 1, jpi*jpj) /), (/ jpi, jpj /) )
245   
246      puniq(:,:) = ztstref(:,:)                   ! default definition
247      CALL lbc_lnk( puniq, cdgrd, 1. )            ! apply boundary conditions
248      lldbl(:,:,1) = puniq(:,:) == ztstref(:,:)   ! check which values have been changed
249     
250      puniq(:,:) = 1.                             ! default definition
251      ! fill only the inner part of the cpu with llbl converted into real
252      puniq(nldi:nlei,nldj:nlej) = REAL(COUNT( lldbl(nldi:nlei,nldj:nlej,:), dim = 3 ), wp)
253
254   END FUNCTION dom_uniq
255
256
257   !!======================================================================
258END MODULE domwri
Note: See TracBrowser for help on using the repository browser.