Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

domwri.F90 @ 2240

Last change on this file since 2240 was 2240, checked in by cetlod, 14 years ago
Suppression of key_zco everywhere in the code
Property svn:eol-style set to `native` Property svn:keywords set to `Id`
File size: 16.0 KB

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

Note: See TracBrowser for help on using the repository browser.

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

Context Navigation

source: branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/DOM/domwri.F90 @ 2240

Download in other formats: