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.

diawri.F90 in branches/2015/dev_CMCC_merge_2015/NEMOGCM/NEMO/SAS_SRC – NEMO

Context Navigation

source: branches/2015/dev_CMCC_merge_2015/NEMOGCM/NEMO/SAS_SRC/diawri.F90 @ 6051

Last change on this file since 6051 was 6051, checked in by lovato, 8 years ago
Merge branches/2015/dev_r5056_CMCC4_simplification (see ticket #1456)
Property svn:keywords set to `Id`
File size: 21.5 KB

Line
1	MODULE diawri
2	!!======================================================================
3	!! * MODULE diawri *
4	!! Ocean diagnostics : write ocean output files
5	!!=====================================================================
6	!! History : OPA ! 1991-03 (M.-A. Foujols) Original code
7	!! 4.0 ! 1991-11 (G. Madec)
8	!! ! 1992-06 (M. Imbard) correction restart file
9	!! ! 1992-07 (M. Imbard) split into diawri and rstwri
10	!! ! 1993-03 (M. Imbard) suppress writibm
11	!! ! 1998-01 (C. Levy) NETCDF format using ioipsl INTERFACE
12	!! ! 1999-02 (E. Guilyardi) name of netCDF files + variables
13	!! 8.2 ! 2000-06 (M. Imbard) Original code (diabort.F)
14	!! NEMO 1.0 ! 2002-06 (A.Bozec, E. Durand) Original code (diainit.F)
15	!! - ! 2002-09 (G. Madec) F90: Free form and module
16	!! - ! 2002-12 (G. Madec) merge of diabort and diainit, F90
17	!! ! 2005-11 (V. Garnier) Surface pressure gradient organization
18	!! 3.2 ! 2008-11 (B. Lemaire) creation from old diawri
19	!!----------------------------------------------------------------------
20
21	!!----------------------------------------------------------------------
22	!! dia_wri : create the standart output files
23	!! dia_wri_state : create an output NetCDF file for a single instantaeous ocean state and forcing fields
24	!!----------------------------------------------------------------------
25	USE oce ! ocean dynamics and tracers
26	USE dom_oce ! ocean space and time domain
27	USE zdf_oce ! ocean vertical physics
28	USE sbc_oce ! Surface boundary condition: ocean fields
29	USE sbc_ice ! Surface boundary condition: ice fields
30	USE sbcssr ! restoring term toward SST/SSS climatology
31	USE phycst ! physical constants
32	USE zdfmxl ! mixed layer
33	USE dianam ! build name of file (routine)
34	USE zdfddm ! vertical physics: double diffusion
35	USE diahth ! thermocline diagnostics
36	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
37	USE in_out_manager ! I/O manager
38	USE diaar5, ONLY : lk_diaar5
39	USE iom
40	USE ioipsl
41	#if defined key_lim2
42	USE limwri_2
43	#endif
44	USE lib_mpp ! MPP library
45	USE timing ! preformance summary
46	USE wrk_nemo ! working array
47
48	IMPLICIT NONE
49	PRIVATE
50
51	PUBLIC dia_wri ! routines called by step.F90
52	PUBLIC dia_wri_state
53	PUBLIC dia_wri_alloc ! Called by nemogcm module
54
55	INTEGER :: nid_T, nz_T, nh_T, ndim_T, ndim_hT ! grid_T file
56	INTEGER :: nid_U, nz_U, nh_U, ndim_U, ndim_hU ! grid_U file
57	INTEGER :: nid_V, nz_V, nh_V, ndim_V, ndim_hV ! grid_V file
58	INTEGER :: ndex(1) ! ???
59	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
60
61	!! * Substitutions
62	# include "zdfddm_substitute.h90"
63	# include "domzgr_substitute.h90"
64	# include "vectopt_loop_substitute.h90"
65	!!----------------------------------------------------------------------
66	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
67	!! $Id$
68	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
69	!!----------------------------------------------------------------------
70	CONTAINS
71
72	INTEGER FUNCTION dia_wri_alloc()
73	!!----------------------------------------------------------------------
74	INTEGER :: ierr
75	!!----------------------------------------------------------------------
76	!
77	ALLOCATE( ndex_hT(jpijpj), ndex_hU(jpijpj), ndex_hV(jpi*jpj), STAT=dia_wri_alloc )
78	IF( lk_mpp ) CALL mpp_sum( dia_wri_alloc )
79	!
80	END FUNCTION dia_wri_alloc
81
82	!!----------------------------------------------------------------------
83	!! Default option NetCDF output file
84	!!----------------------------------------------------------------------
85	#if defined key_iomput
86	!!----------------------------------------------------------------------
87	!! 'key_iomput' use IOM library
88	!!----------------------------------------------------------------------
89
90	SUBROUTINE dia_wri( kt )
91	!!---------------------------------------------------------------------
92	!! * ROUTINE dia_wri *
93	!!
94	!! ** Purpose : Standard output of opa: dynamics and tracer fields
95	!! NETCDF format is used by default
96	!! Standalone surface scheme
97	!!
98	!! ** Method : use iom_put
99	!!----------------------------------------------------------------------
100	!!
101	INTEGER, INTENT( in ) :: kt ! ocean time-step index
102	!!----------------------------------------------------------------------
103	!
104	!! no relevant 2D arrays to write in iomput case
105	!
106	END SUBROUTINE dia_wri
107
108	#else
109	!!----------------------------------------------------------------------
110	!! Default option use IOIPSL library
111	!!----------------------------------------------------------------------
112
113	SUBROUTINE dia_wri( kt )
114	!!---------------------------------------------------------------------
115	!! * ROUTINE dia_wri *
116	!!
117	!! ** Purpose : Standard output of opa: dynamics and tracer fields
118	!! NETCDF format is used by default
119	!!
120	!! ** Method : At the beginning of the first time step (nit000),
121	!! define all the NETCDF files and fields
122	!! At each time step call histdef to compute the mean if ncessary
123	!! Each nwrite time step, output the instantaneous or mean fields
124	!!----------------------------------------------------------------------
125	!!
126	INTEGER, INTENT( in ) :: kt ! ocean time-step index
127	!!
128	LOGICAL :: ll_print = .FALSE. ! =T print and flush numout
129	CHARACTER (len=40) :: clhstnam, clop, clmx ! local names
130	INTEGER :: inum = 11 ! temporary logical unit
131	INTEGER :: ji, jj, jk ! dummy loop indices
132	INTEGER :: ierr ! error code return from allocation
133	INTEGER :: iimi, iima, ipk, it, itmod, ijmi, ijma ! local integers
134	REAL(wp) :: zsto, zout, zmax, zjulian ! local scalars
135	!!----------------------------------------------------------------------
136	!
137	IF( nn_timing == 1 ) CALL timing_start('dia_wri')
138	!
139	! Output the initial state and forcings
140	IF( ninist == 1 ) THEN
141	CALL dia_wri_state( 'output.init', kt )
142	ninist = 0
143	ENDIF
144	!
145	! 0. Initialisation
146	! -----------------
147
148	! local variable for debugging
149	ll_print = .FALSE.
150	ll_print = ll_print .AND. lwp
151
152	! Define frequency of output and means
153	IF( ln_mskland ) THEN ; clop = "only(x)" ! put 1.e+20 on land (very expensive!!)
154	ELSE ; clop = "x" ! no use of the mask value (require less cpu time)
155	ENDIF
156	#if defined key_diainstant
157	zsto = nwrite * rdt
158	clop = "inst("//TRIM(clop)//")"
159	#else
160	zsto=rdt
161	clop = "ave("//TRIM(clop)//")"
162	#endif
163	zout = nwrite * rdt
164	zmax = ( nitend - nit000 + 1 ) * rdt
165
166	! Define indices of the horizontal output zoom and vertical limit storage
167	iimi = 1 ; iima = jpi
168	ijmi = 1 ; ijma = jpj
169	ipk = jpk
170
171	! define time axis
172	it = kt
173	itmod = kt - nit000 + 1
174
175
176	! 1. Define NETCDF files and fields at beginning of first time step
177	! -----------------------------------------------------------------
178
179	IF( kt == nit000 ) THEN
180
181	! Define the NETCDF files (one per grid)
182
183	! Compute julian date from starting date of the run
184	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
185	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
186	IF(lwp)WRITE(numout,*)
187	IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear, &
188	& ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
189	IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma, &
190	' limit storage in depth = ', ipk
191
192	! WRITE root name in date.file for use by postpro
193	IF(lwp) THEN
194	CALL dia_nam( clhstnam, nwrite,' ' )
195	CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
196	WRITE(inum,*) clhstnam
197	CLOSE(inum)
198	ENDIF
199
200	! Define the T grid FILE ( nid_T )
201
202	CALL dia_nam( clhstnam, nwrite, 'grid_T' )
203	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
204	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
205	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
206	& nit000-1, zjulian, rdt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )
207	CALL histvert( nid_T, "deptht", "Vertical T levels", & ! Vertical grid: gdept
208	& "m", ipk, gdept_1d, nz_T, "down" )
209	! ! Index of ocean points
210	CALL wheneq( jpi*jpj , tmask, 1, 1., ndex_hT, ndim_hT ) ! surface
211
212	! Define the U grid FILE ( nid_U )
213
214	CALL dia_nam( clhstnam, nwrite, 'grid_U' )
215	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
216	CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu, & ! Horizontal grid: glamu and gphiu
217	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
218	& nit000-1, zjulian, rdt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )
219	CALL histvert( nid_U, "depthu", "Vertical U levels", & ! Vertical grid: gdept
220	& "m", ipk, gdept_1d, nz_U, "down" )
221	! ! Index of ocean points
222	CALL wheneq( jpi*jpj , umask, 1, 1., ndex_hU, ndim_hU ) ! surface
223
224	! Define the V grid FILE ( nid_V )
225
226	CALL dia_nam( clhstnam, nwrite, 'grid_V' ) ! filename
227	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
228	CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv, & ! Horizontal grid: glamv and gphiv
229	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
230	& nit000-1, zjulian, rdt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )
231	CALL histvert( nid_V, "depthv", "Vertical V levels", & ! Vertical grid : gdept
232	& "m", ipk, gdept_1d, nz_V, "down" )
233	! ! Index of ocean points
234	CALL wheneq( jpi*jpj , vmask, 1, 1., ndex_hV, ndim_hV ) ! surface
235
236	! No W grid FILE
237
238	! Declare all the output fields as NETCDF variables
239
240	! !!! nid_T : 3D
241	CALL histdef( nid_T, "sst_m", "Sea Surface temperature" , "C" , & ! sst
242	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
243	CALL histdef( nid_T, "sss_m", "Sea Surface Salinity" , "PSU" , & ! sss
244	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
245	CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux" , "Kg/m2/s", & ! (emp-rnf)
246	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
247	CALL histdef( nid_T, "sosfldow", "downward salt flux" , "PSU/m2/s", & ! (sfx)
248	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
249	CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! qns + qsr
250	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
251	CALL histdef( nid_T, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! qsr
252	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
253	CALL histdef( nid_T, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
254	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
255	CALL histdef( nid_T, "sowindsp", "wind speed at 10m" , "m/s" , & ! wndm
256	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
257
258	CALL histend( nid_T, snc4chunks=snc4set )
259
260	! !!! nid_U : 3D
261	CALL histdef( nid_U, "ssu_m", "Velocity component in x-direction", "m/s" , & ! ssu
262	& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
263	CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau
264	& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
265
266	CALL histend( nid_U, snc4chunks=snc4set )
267
268	! !!! nid_V : 3D
269	CALL histdef( nid_V, "ssv_m", "Velocity component in y-direction", "m/s", & ! ssv_m
270	& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
271	CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau
272	& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
273
274	CALL histend( nid_V, snc4chunks=snc4set )
275
276	IF(lwp) WRITE(numout,*)
277	IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
278	IF(ll_print) CALL FLUSH(numout )
279
280	ENDIF
281
282	! 2. Start writing data
283	! ---------------------
284
285	! ndex(1) est utilise ssi l'avant dernier argument est diffferent de
286	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
287	! donne le nombre d'elements, et ndex la liste des indices a sortir
288
289	IF( lwp .AND. MOD( itmod, nwrite ) == 0 ) THEN
290	WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
291	WRITE(numout,*) '~~~~~~ '
292	ENDIF
293
294	! Write fields on T grid
295	CALL histwrite( nid_T, "sst_m", it, sst_m, ndim_hT, ndex_hT ) ! sea surface temperature
296	CALL histwrite( nid_T, "sss_m", it, sss_m, ndim_hT, ndex_hT ) ! sea surface salinity
297	CALL histwrite( nid_T, "sowaflup", it, (emp - rnf ) , ndim_hT, ndex_hT ) ! upward water flux
298	CALL histwrite( nid_T, "sosfldow", it, sfx , ndim_hT, ndex_hT ) ! downward salt flux
299	! (includes virtual salt flux beneath ice
300	! in linear free surface case)
301
302	CALL histwrite( nid_T, "sohefldo", it, qns + qsr , ndim_hT, ndex_hT ) ! total heat flux
303	CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux
304	CALL histwrite( nid_T, "soicecov", it, fr_i , ndim_hT, ndex_hT ) ! ice fraction
305	CALL histwrite( nid_T, "sowindsp", it, wndm , ndim_hT, ndex_hT ) ! wind speed
306
307	! Write fields on U grid
308	CALL histwrite( nid_U, "ssu_m" , it, ssu_m , ndim_hU, ndex_hU ) ! i-current speed
309	CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress
310
311	! Write fields on V grid
312	CALL histwrite( nid_V, "ssv_m" , it, ssv_m , ndim_hV, ndex_hV ) ! j-current speed
313	CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress
314
315	! 3. Close all files
316	! ---------------------------------------
317	IF( kt == nitend ) THEN
318	CALL histclo( nid_T )
319	CALL histclo( nid_U )
320	CALL histclo( nid_V )
321	ENDIF
322	!
323	IF( nn_timing == 1 ) CALL timing_stop('dia_wri')
324	!
325	END SUBROUTINE dia_wri
326	#endif
327
328	SUBROUTINE dia_wri_state( cdfile_name, kt )
329	!!---------------------------------------------------------------------
330	!! * ROUTINE dia_wri_state *
331	!!
332	!! ** Purpose : create a NetCDF file named cdfile_name which contains
333	!! the instantaneous ocean state and forcing fields.
334	!! Used to find errors in the initial state or save the last
335	!! ocean state in case of abnormal end of a simulation
336	!!
337	!! ** Method : NetCDF files using ioipsl
338	!! File 'output.init.nc' is created if ninist = 1 (namelist)
339	!! File 'output.abort.nc' is created in case of abnormal job end
340	!!----------------------------------------------------------------------
341	CHARACTER (len=* ), INTENT( in ) :: cdfile_name ! name of the file created
342	INTEGER , INTENT( in ) :: kt ! ocean time-step index
343	!!
344	CHARACTER (len=32) :: clname
345	CHARACTER (len=40) :: clop
346	INTEGER :: id_i , nz_i, nh_i
347	INTEGER, DIMENSION(1) :: idex ! local workspace
348	REAL(wp) :: zsto, zout, zmax, zjulian
349	!!----------------------------------------------------------------------
350	!
351	IF( nn_timing == 1 ) CALL timing_start('dia_wri_state')
352
353	! 0. Initialisation
354	! -----------------
355
356	! Define name, frequency of output and means
357	clname = cdfile_name
358	IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
359	zsto = rdt
360	clop = "inst(x)" ! no use of the mask value (require less cpu time)
361	zout = rdt
362	zmax = ( nitend - nit000 + 1 ) * rdt
363
364	IF(lwp) WRITE(numout,*)
365	IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
366	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ and forcing fields file created '
367	IF(lwp) WRITE(numout,*) ' and named :', clname, '.nc'
368
369
370	! 1. Define NETCDF files and fields at beginning of first time step
371	! -----------------------------------------------------------------
372
373	! Compute julian date from starting date of the run
374	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian ) ! time axis
375	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
376	CALL histbeg( clname, jpi, glamt, jpj, gphit, &
377	1, jpi, 1, jpj, nit000-1, zjulian, rdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
378	CALL histvert( id_i, "deptht", "Vertical T levels", & ! Vertical grid : gdept
379	"m", jpk, gdept_1d, nz_i, "down")
380
381	! Declare all the output fields as NetCDF variables
382
383	CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S", & ! net freshwater
384	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
385	CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2" , & ! net heat flux
386	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
387	CALL histdef( id_i, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! solar flux
388	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
389	CALL histdef( id_i, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
390	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
391	CALL histdef( id_i, "sozotaux", "Zonal Wind Stress" , "N/m2" , & ! i-wind stress
392	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
393	CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2" , & ! j-wind stress
394	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
395
396	#if defined key_lim2
397	CALL lim_wri_state_2( kt, id_i, nh_i )
398	#else
399	CALL histend( id_i, snc4chunks=snc4set )
400	#endif
401
402	! 2. Start writing data
403	! ---------------------
404	! idex(1) est utilise ssi l'avant dernier argument est diffferent de
405	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
406	! donne le nombre d'elements, et idex la liste des indices a sortir
407	idex(1) = 1 ! init to avoid compil warning
408
409	! Write all fields on T grid
410	CALL histwrite( id_i, "sowaflup", kt, emp , jpi*jpj , idex ) ! freshwater budget
411	CALL histwrite( id_i, "sohefldo", kt, qsr + qns , jpi*jpj , idex ) ! total heat flux
412	CALL histwrite( id_i, "soshfldo", kt, qsr , jpi*jpj , idex ) ! solar heat flux
413	CALL histwrite( id_i, "soicecov", kt, fr_i , jpi*jpj , idex ) ! ice fraction
414	CALL histwrite( id_i, "sozotaux", kt, utau , jpi*jpj , idex ) ! i-wind stress
415	CALL histwrite( id_i, "sometauy", kt, vtau , jpi*jpj , idex ) ! j-wind stress
416
417	! 3. Close the file
418	! -----------------
419	CALL histclo( id_i )
420	#if ! defined key_iomput
421	IF( ninist /= 1 ) THEN
422	CALL histclo( nid_T )
423	CALL histclo( nid_U )
424	CALL histclo( nid_V )
425	ENDIF
426	#endif
427
428	IF( nn_timing == 1 ) CALL timing_stop('dia_wri_state')
429	!
430
431	END SUBROUTINE dia_wri_state
432	!!======================================================================
433	END MODULE diawri

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

Download in other formats: