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/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

Context Navigation

source: branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 @ 8406

Last change on this file since 8406 was 8400, checked in by timgraham, 7 years ago

GMED ticket 335:

Merge dev_r5518_GO6_package_inc_asm into package branch to make everything easier for data assimilation
No effect on configs without data assimilation

File size: 58.9 KB

Rev	Line
[3]	1	MODULE diawri
	2	!!======================================================================
	3	!! * MODULE diawri *
	4	!! Ocean diagnostics : write ocean output files
	5	!!=====================================================================
[2528]	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	!!----------------------------------------------------------------------
[3]	20
	21	!!----------------------------------------------------------------------
[2528]	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	!!----------------------------------------------------------------------
[3]	25	USE oce ! ocean dynamics and tracers
	26	USE dom_oce ! ocean space and time domain
[4292]	27	USE dynadv, ONLY: ln_dynadv_vec
[3]	28	USE zdf_oce ! ocean vertical physics
	29	USE ldftra_oce ! ocean active tracers: lateral physics
	30	USE ldfdyn_oce ! ocean dynamics: lateral physics
[3294]	31	USE traldf_iso_grif, ONLY : psix_eiv, psiy_eiv
[3]	32	USE sol_oce ! solver variables
[888]	33	USE sbc_oce ! Surface boundary condition: ocean fields
	34	USE sbc_ice ! Surface boundary condition: ice fields
[3609]	35	USE icb_oce ! Icebergs
	36	USE icbdia ! Iceberg budgets
[888]	37	USE sbcssr ! restoring term toward SST/SSS climatology
[3]	38	USE phycst ! physical constants
	39	USE zdfmxl ! mixed layer
	40	USE dianam ! build name of file (routine)
[6491]	41	USE zdftke ! vertical physics: one-equation scheme
[3]	42	USE zdfddm ! vertical physics: double diffusion
	43	USE diahth ! thermocline diagnostics
	44	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	45	USE in_out_manager ! I/O manager
[216]	46	USE diadimg ! dimg direct access file format output
[1359]	47	USE iom
[460]	48	USE ioipsl
[8400]	49	USE dynspg_oce, ONLY: un_adv, vn_adv ! barotropic velocities
	50	USE insitu_tem, ONLY: insitu_t, theta2t
[1482]	51	#if defined key_lim2
	52	USE limwri_2
[4161]	53	#elif defined key_lim3
	54	USE limwri
[1482]	55	#endif
[2715]	56	USE lib_mpp ! MPP library
[3294]	57	USE timing ! preformance summary
	58	USE wrk_nemo ! working array
[2528]	59
[3]	60	IMPLICIT NONE
	61	PRIVATE
	62
[2528]	63	PUBLIC dia_wri ! routines called by step.F90
	64	PUBLIC dia_wri_state
[2715]	65	PUBLIC dia_wri_alloc ! Called by nemogcm module
[3]	66
[2528]	67	INTEGER :: nid_T, nz_T, nh_T, ndim_T, ndim_hT ! grid_T file
[3609]	68	INTEGER :: nb_T , ndim_bT ! grid_T file
[2528]	69	INTEGER :: nid_U, nz_U, nh_U, ndim_U, ndim_hU ! grid_U file
	70	INTEGER :: nid_V, nz_V, nh_V, ndim_V, ndim_hV ! grid_V file
	71	INTEGER :: nid_W, nz_W, nh_W ! grid_W file
	72	INTEGER :: ndex(1) ! ???
[2715]	73	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
	74	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V
[3609]	75	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_bT
[3]	76
	77	!! * Substitutions
	78	# include "zdfddm_substitute.h90"
[1756]	79	# include "domzgr_substitute.h90"
	80	# include "vectopt_loop_substitute.h90"
[3]	81	!!----------------------------------------------------------------------
[2528]	82	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[5217]	83	!! $Id$
[2528]	84	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	85	!!----------------------------------------------------------------------
	86	CONTAINS
	87
[2715]	88	INTEGER FUNCTION dia_wri_alloc()
	89	!!----------------------------------------------------------------------
	90	INTEGER, DIMENSION(2) :: ierr
	91	!!----------------------------------------------------------------------
	92	ierr = 0
	93	ALLOCATE( ndex_hT(jpijpj) , ndex_T(jpijpj*jpk) , &
	94	& ndex_hU(jpijpj) , ndex_U(jpijpj*jpk) , &
	95	& ndex_hV(jpijpj) , ndex_V(jpijpj*jpk) , STAT=ierr(1) )
	96	!
	97	dia_wri_alloc = MAXVAL(ierr)
	98	IF( lk_mpp ) CALL mpp_sum( dia_wri_alloc )
	99	!
	100	END FUNCTION dia_wri_alloc
	101
[216]	102	#if defined key_dimgout
[3]	103	!!----------------------------------------------------------------------
[2528]	104	!! 'key_dimgout' DIMG output file
[3]	105	!!----------------------------------------------------------------------
	106	# include "diawri_dimg.h90"
	107
	108	#else
	109	!!----------------------------------------------------------------------
	110	!! Default option NetCDF output file
	111	!!----------------------------------------------------------------------
[2528]	112	# if defined key_iomput
[3]	113	!!----------------------------------------------------------------------
[2528]	114	!! 'key_iomput' use IOM library
	115	!!----------------------------------------------------------------------
[2715]	116
[1561]	117	SUBROUTINE dia_wri( kt )
[1482]	118	!!---------------------------------------------------------------------
	119	!! * ROUTINE dia_wri *
	120	!!
	121	!! ** Purpose : Standard output of opa: dynamics and tracer fields
	122	!! NETCDF format is used by default
	123	!!
	124	!! ** Method : use iom_put
	125	!!----------------------------------------------------------------------
[1756]	126	!!
[1482]	127	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[1756]	128	!!
	129	INTEGER :: ji, jj, jk ! dummy loop indices
[5460]	130	INTEGER :: jkbot !
[1756]	131	REAL(wp) :: zztmp, zztmpx, zztmpy !
[3294]	132	!!
[4761]	133	REAL(wp), POINTER, DIMENSION(:,:) :: z2d ! 2D workspace
[3294]	134	REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d ! 3D workspace
[1482]	135	!!----------------------------------------------------------------------
	136	!
[3294]	137	IF( nn_timing == 1 ) CALL timing_start('dia_wri')
	138	!
[4990]	139	CALL wrk_alloc( jpi , jpj , z2d )
[3294]	140	CALL wrk_alloc( jpi , jpj, jpk , z3d )
[2715]	141	!
[1482]	142	! Output the initial state and forcings
	143	IF( ninist == 1 ) THEN
	144	CALL dia_wri_state( 'output.init', kt )
	145	ninist = 0
	146	ENDIF
[3]	147
[6498]	148	! Output of initial vertical scale factor
	149	CALL iom_put("e3t_0", e3t_0(:,:,:) )
	150	CALL iom_put("e3u_0", e3t_0(:,:,:) )
	151	CALL iom_put("e3v_0", e3t_0(:,:,:) )
	152	!
	153	CALL iom_put( "e3t" , fse3t_n(:,:,:) )
	154	CALL iom_put( "e3u" , fse3u_n(:,:,:) )
	155	CALL iom_put( "e3v" , fse3v_n(:,:,:) )
	156	CALL iom_put( "e3w" , fse3w_n(:,:,:) )
	157	IF( iom_use("e3tdef") ) &
	158	CALL iom_put( "e3tdef" , ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 )
[7335]	159	CALL iom_put("tpt_dep", fsdept_n(:,:,:) )
[7651]	160	CALL iom_put("wpt_dep", fsdepw_n(:,:,:) )
[5461]	161
[6498]	162
[5461]	163	CALL iom_put( "ssh" , sshn ) ! sea surface height
[5107]	164
	165	CALL iom_put( "toce", tsn(:,:,:,jp_tem) ) ! 3D temperature
[8400]	166	CALL theta2t ! in-situ temperature conversion
	167	CALL iom_put( "tinsitu", insitu_t(:,:,:)) ! in-situ temperature
[5107]	168	CALL iom_put( "sst", tsn(:,:,1,jp_tem) ) ! surface temperature
	169	IF ( iom_use("sbt") ) THEN
[4990]	170	DO jj = 1, jpj
	171	DO ji = 1, jpi
[5460]	172	jkbot = mbkt(ji,jj)
	173	z2d(ji,jj) = tsn(ji,jj,jkbot,jp_tem)
[4990]	174	END DO
[5107]	175	END DO
	176	CALL iom_put( "sbt", z2d ) ! bottom temperature
	177	ENDIF
	178
	179	CALL iom_put( "soce", tsn(:,:,:,jp_sal) ) ! 3D salinity
	180	CALL iom_put( "sss", tsn(:,:,1,jp_sal) ) ! surface salinity
	181	IF ( iom_use("sbs") ) THEN
[4990]	182	DO jj = 1, jpj
	183	DO ji = 1, jpi
[5460]	184	jkbot = mbkt(ji,jj)
	185	z2d(ji,jj) = tsn(ji,jj,jkbot,jp_sal)
[4990]	186	END DO
[5107]	187	END DO
	188	CALL iom_put( "sbs", z2d ) ! bottom salinity
	189	ENDIF
[5463]	190
	191	IF ( iom_use("taubot") ) THEN ! bottom stress
	192	z2d(:,:) = 0._wp
	193	DO jj = 2, jpjm1
	194	DO ji = fs_2, fs_jpim1 ! vector opt.
	195	zztmpx = ( bfrua(ji ,jj) * un(ji ,jj,mbku(ji ,jj)) &
	196	& + bfrua(ji-1,jj) * un(ji-1,jj,mbku(ji-1,jj)) )
	197	zztmpy = ( bfrva(ji, jj) * vn(ji,jj ,mbkv(ji,jj )) &
	198	& + bfrva(ji,jj-1) * vn(ji,jj-1,mbkv(ji,jj-1)) )
	199	z2d(ji,jj) = rau0 * SQRT( zztmpx * zztmpx + zztmpy * zztmpy ) * tmask(ji,jj,1)
	200	!
	201	ENDDO
	202	ENDDO
	203	CALL lbc_lnk( z2d, 'T', 1. )
	204	CALL iom_put( "taubot", z2d )
	205	ENDIF
[8400]	206
[5107]	207	CALL iom_put( "uoce", un(:,:,:) ) ! 3D i-current
	208	CALL iom_put( "ssu", un(:,:,1) ) ! surface i-current
	209	IF ( iom_use("sbu") ) THEN
[4990]	210	DO jj = 1, jpj
	211	DO ji = 1, jpi
[5460]	212	jkbot = mbku(ji,jj)
	213	z2d(ji,jj) = un(ji,jj,jkbot)
[4990]	214	END DO
[5107]	215	END DO
	216	CALL iom_put( "sbu", z2d ) ! bottom i-current
	217	ENDIF
[5463]	218	#if defined key_dynspg_ts
	219	CALL iom_put( "ubar", un_adv(:,:) ) ! barotropic i-current
	220	#else
	221	CALL iom_put( "ubar", un_b(:,:) ) ! barotropic i-current
	222	#endif
[5107]	223
	224	CALL iom_put( "voce", vn(:,:,:) ) ! 3D j-current
	225	CALL iom_put( "ssv", vn(:,:,1) ) ! surface j-current
	226	IF ( iom_use("sbv") ) THEN
[4990]	227	DO jj = 1, jpj
	228	DO ji = 1, jpi
[5460]	229	jkbot = mbkv(ji,jj)
	230	z2d(ji,jj) = vn(ji,jj,jkbot)
[4990]	231	END DO
[5107]	232	END DO
	233	CALL iom_put( "sbv", z2d ) ! bottom j-current
[4990]	234	ENDIF
[5463]	235	#if defined key_dynspg_ts
	236	CALL iom_put( "vbar", vn_adv(:,:) ) ! barotropic j-current
	237	#else
	238	CALL iom_put( "vbar", vn_b(:,:) ) ! barotropic j-current
	239	#endif
[1482]	240
[5461]	241	CALL iom_put( "woce", wn ) ! vertical velocity
	242	IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN ! vertical mass transport & its square value
	243	! Caution: in the VVL case, it only correponds to the baroclinic mass transport.
	244	z2d(:,:) = rau0 * e12t(:,:)
	245	DO jk = 1, jpk
	246	z3d(:,:,jk) = wn(:,:,jk) * z2d(:,:)
	247	END DO
	248	CALL iom_put( "w_masstr" , z3d )
	249	IF( iom_use('w_masstr2') ) CALL iom_put( "w_masstr2", z3d(:,:,:) * z3d(:,:,:) )
	250	ENDIF
	251
[5107]	252	CALL iom_put( "avt" , avt ) ! T vert. eddy diff. coef.
	253	CALL iom_put( "avm" , avmu ) ! T vert. eddy visc. coef.
[6491]	254	IF( lk_zdftke ) THEN
	255	CALL iom_put( "tke" , en ) ! TKE budget: Turbulent Kinetic Energy
	256	CALL iom_put( "tke_niw" , e_niw ) ! TKE budget: Near-inertial waves
	257	ENDIF
[5107]	258	CALL iom_put( "avs" , fsavs(:,:,:) ) ! S vert. eddy diff. coef. (useful only with key_zdfddm)
[6498]	259	! Log of eddy diff coef
	260	IF( iom_use('logavt') ) CALL iom_put( "logavt", LOG( MAX( 1.e-20_wp, avt (:,:,:) ) ) )
	261	IF( iom_use('logavs') ) CALL iom_put( "logavs", LOG( MAX( 1.e-20_wp, fsavs(:,:,:) ) ) )
[5107]	262
	263	IF ( iom_use("sstgrad") .OR. iom_use("sstgrad2") ) THEN
[4990]	264	DO jj = 2, jpjm1 ! sst gradient
	265	DO ji = fs_2, fs_jpim1 ! vector opt.
	266	zztmp = tsn(ji,jj,1,jp_tem)
	267	zztmpx = ( tsn(ji+1,jj ,1,jp_tem) - zztmp ) / e1u(ji,jj) + ( zztmp - tsn(ji-1,jj ,1,jp_tem) ) / e1u(ji-1,jj )
	268	zztmpy = ( tsn(ji ,jj+1,1,jp_tem) - zztmp ) / e2v(ji,jj) + ( zztmp - tsn(ji ,jj-1,1,jp_tem) ) / e2v(ji ,jj-1)
	269	z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy ) &
	270	& * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1)
	271	END DO
[1756]	272	END DO
[4990]	273	CALL lbc_lnk( z2d, 'T', 1. )
	274	CALL iom_put( "sstgrad2", z2d ) ! square of module of sst gradient
	275	z2d(:,:) = SQRT( z2d(:,:) )
	276	CALL iom_put( "sstgrad" , z2d ) ! module of sst gradient
	277	ENDIF
	278
[4761]	279	! clem: heat and salt content
[4990]	280	IF( iom_use("heatc") ) THEN
	281	z2d(:,:) = 0._wp
	282	DO jk = 1, jpkm1
[5107]	283	DO jj = 1, jpj
	284	DO ji = 1, jpi
[4990]	285	z2d(ji,jj) = z2d(ji,jj) + fse3t(ji,jj,jk) * tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk)
	286	END DO
[4761]	287	END DO
	288	END DO
[4990]	289	CALL iom_put( "heatc", (rau0 * rcp) * z2d ) ! vertically integrated heat content (J/m2)
	290	ENDIF
	291
	292	IF( iom_use("saltc") ) THEN
	293	z2d(:,:) = 0._wp
	294	DO jk = 1, jpkm1
[5107]	295	DO jj = 1, jpj
	296	DO ji = 1, jpi
[4990]	297	z2d(ji,jj) = z2d(ji,jj) + fse3t(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk)
	298	END DO
	299	END DO
	300	END DO
	301	CALL iom_put( "saltc", rau0 * z2d ) ! vertically integrated salt content (PSU*kg/m2)
	302	ENDIF
[4840]	303	!
[4990]	304	IF ( iom_use("eken") ) THEN
	305	rke(:,:,jk) = 0._wp ! kinetic energy
	306	DO jk = 1, jpkm1
	307	DO jj = 2, jpjm1
	308	DO ji = fs_2, fs_jpim1 ! vector opt.
	309	zztmp = 1._wp / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) )
	310	zztmpx = 0.5 * ( un(ji-1,jj,jk) * un(ji-1,jj,jk) * e2u(ji-1,jj) * fse3u(ji-1,jj,jk) &
	311	& + un(ji ,jj,jk) * un(ji ,jj,jk) * e2u(ji ,jj) * fse3u(ji ,jj,jk) ) &
	312	& * zztmp
	313	!
	314	zztmpy = 0.5 * ( vn(ji,jj-1,jk) * vn(ji,jj-1,jk) * e1v(ji,jj-1) * fse3v(ji,jj-1,jk) &
	315	& + vn(ji,jj ,jk) * vn(ji,jj ,jk) * e1v(ji,jj ) * fse3v(ji,jj ,jk) ) &
	316	& * zztmp
	317	!
	318	rke(ji,jj,jk) = 0.5_wp * ( zztmpx + zztmpy )
	319	!
	320	ENDDO
[4840]	321	ENDDO
	322	ENDDO
[4990]	323	CALL lbc_lnk( rke, 'T', 1. )
	324	CALL iom_put( "eken", rke )
	325	ENDIF
[6498]	326	!
	327	CALL iom_put( "hdiv", hdivn ) ! Horizontal divergence
	328	!
[7179]	329	IF( iom_use("u_masstr") .OR. iom_use("u_masstr_vint") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN
[1756]	330	z3d(:,:,jpk) = 0.e0
[7179]	331	z2d(:,:) = 0.e0
[1756]	332	DO jk = 1, jpkm1
[4761]	333	z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * fse3u(:,:,jk) * umask(:,:,jk)
[7179]	334	z2d(:,:) = z2d(:,:) + z3d(:,:,jk)
[1756]	335	END DO
	336	CALL iom_put( "u_masstr", z3d ) ! mass transport in i-direction
[7179]	337	CALL iom_put( "u_masstr_vint", z2d ) ! mass transport in i-direction vertical sum
[4990]	338	ENDIF
	339
	340	IF( iom_use("u_heattr") ) THEN
	341	z2d(:,:) = 0.e0
	342	DO jk = 1, jpkm1
	343	DO jj = 2, jpjm1
	344	DO ji = fs_2, fs_jpim1 ! vector opt.
	345	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji+1,jj,jk,jp_tem) )
	346	END DO
	347	END DO
	348	END DO
	349	CALL lbc_lnk( z2d, 'U', -1. )
	350	CALL iom_put( "u_heattr", (0.5 * rcp) * z2d ) ! heat transport in i-direction
	351	ENDIF
[4761]	352
[4990]	353	IF( iom_use("u_salttr") ) THEN
[1756]	354	z2d(:,:) = 0.e0
	355	DO jk = 1, jpkm1
	356	DO jj = 2, jpjm1
	357	DO ji = fs_2, fs_jpim1 ! vector opt.
[4990]	358	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji+1,jj,jk,jp_sal) )
[1756]	359	END DO
	360	END DO
	361	END DO
	362	CALL lbc_lnk( z2d, 'U', -1. )
[4990]	363	CALL iom_put( "u_salttr", 0.5 * z2d ) ! heat transport in i-direction
	364	ENDIF
[4761]	365
[4990]	366
	367	IF( iom_use("v_masstr") .OR. iom_use("v_heattr") .OR. iom_use("v_salttr") ) THEN
[4761]	368	z3d(:,:,jpk) = 0.e0
[1756]	369	DO jk = 1, jpkm1
[4761]	370	z3d(:,:,jk) = rau0 * vn(:,:,jk) * e1v(:,:) * fse3v(:,:,jk) * vmask(:,:,jk)
[1756]	371	END DO
	372	CALL iom_put( "v_masstr", z3d ) ! mass transport in j-direction
[4990]	373	ENDIF
	374
	375	IF( iom_use("v_heattr") ) THEN
	376	z2d(:,:) = 0.e0
	377	DO jk = 1, jpkm1
	378	DO jj = 2, jpjm1
	379	DO ji = fs_2, fs_jpim1 ! vector opt.
	380	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji,jj+1,jk,jp_tem) )
	381	END DO
	382	END DO
	383	END DO
	384	CALL lbc_lnk( z2d, 'V', -1. )
	385	CALL iom_put( "v_heattr", (0.5 * rcp) * z2d ) ! heat transport in j-direction
	386	ENDIF
[4761]	387
[4990]	388	IF( iom_use("v_salttr") ) THEN
[1756]	389	z2d(:,:) = 0.e0
	390	DO jk = 1, jpkm1
	391	DO jj = 2, jpjm1
	392	DO ji = fs_2, fs_jpim1 ! vector opt.
[4990]	393	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji,jj+1,jk,jp_sal) )
[1756]	394	END DO
	395	END DO
	396	END DO
	397	CALL lbc_lnk( z2d, 'V', -1. )
[4990]	398	CALL iom_put( "v_salttr", 0.5 * z2d ) ! heat transport in j-direction
[1756]	399	ENDIF
[7179]	400
	401	! Vertical integral of temperature
	402	IF( iom_use("tosmint") ) THEN
	403	z2d(:,:)=0._wp
	404	DO jk = 1, jpkm1
	405	DO jj = 2, jpjm1
	406	DO ji = fs_2, fs_jpim1 ! vector opt.
	407	z2d(ji,jj) = z2d(ji,jj) + rau0 * fse3t(ji,jj,jk) * tsn(ji,jj,jk,jp_tem)
	408	END DO
	409	END DO
	410	END DO
	411	CALL lbc_lnk( z2d, 'T', -1. )
	412	CALL iom_put( "tosmint", z2d )
	413	ENDIF
	414
	415	! Vertical integral of salinity
	416	IF( iom_use("somint") ) THEN
	417	z2d(:,:)=0._wp
	418	DO jk = 1, jpkm1
	419	DO jj = 2, jpjm1
	420	DO ji = fs_2, fs_jpim1 ! vector opt.
	421	z2d(ji,jj) = z2d(ji,jj) + rau0 * fse3t(ji,jj,jk) * tsn(ji,jj,jk,jp_sal)
	422	END DO
	423	END DO
	424	END DO
	425	CALL lbc_lnk( z2d, 'T', -1. )
	426	CALL iom_put( "somint", z2d )
	427	ENDIF
	428
	429	CALL iom_put( "bn2", rn2 ) !Brunt-Vaisala buoyancy frequency (N^2)
[2528]	430	!
[4990]	431	CALL wrk_dealloc( jpi , jpj , z2d )
[3294]	432	CALL wrk_dealloc( jpi , jpj, jpk , z3d )
[2715]	433	!
[3294]	434	IF( nn_timing == 1 ) CALL timing_stop('dia_wri')
	435	!
[1482]	436	END SUBROUTINE dia_wri
	437
	438	#else
[2528]	439	!!----------------------------------------------------------------------
	440	!! Default option use IOIPSL library
	441	!!----------------------------------------------------------------------
	442
[1561]	443	SUBROUTINE dia_wri( kt )
[3]	444	!!---------------------------------------------------------------------
	445	!! * ROUTINE dia_wri *
	446	!!
	447	!! ** Purpose : Standard output of opa: dynamics and tracer fields
	448	!! NETCDF format is used by default
	449	!!
	450	!! ** Method : At the beginning of the first time step (nit000),
	451	!! define all the NETCDF files and fields
	452	!! At each time step call histdef to compute the mean if ncessary
	453	!! Each nwrite time step, output the instantaneous or mean fields
	454	!!----------------------------------------------------------------------
[2715]	455	!!
[3]	456	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[2528]	457	!!
	458	LOGICAL :: ll_print = .FALSE. ! =T print and flush numout
	459	CHARACTER (len=40) :: clhstnam, clop, clmx ! local names
	460	INTEGER :: inum = 11 ! temporary logical unit
[3294]	461	INTEGER :: ji, jj, jk ! dummy loop indices
	462	INTEGER :: ierr ! error code return from allocation
[2528]	463	INTEGER :: iimi, iima, ipk, it, itmod, ijmi, ijma ! local integers
[3609]	464	INTEGER :: jn, ierror ! local integers
[2528]	465	REAL(wp) :: zsto, zout, zmax, zjulian, zdt ! local scalars
[3294]	466	!!
	467	REAL(wp), POINTER, DIMENSION(:,:) :: zw2d ! 2D workspace
	468	REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d ! 3D workspace
[3]	469	!!----------------------------------------------------------------------
[3294]	470	!
	471	IF( nn_timing == 1 ) CALL timing_start('dia_wri')
[1482]	472	!
[3294]	473	CALL wrk_alloc( jpi , jpj , zw2d )
[4292]	474	IF ( ln_traldf_gdia .OR. lk_vvl ) call wrk_alloc( jpi , jpj , jpk , zw3d )
[2715]	475	!
[1482]	476	! Output the initial state and forcings
	477	IF( ninist == 1 ) THEN
	478	CALL dia_wri_state( 'output.init', kt )
	479	ninist = 0
	480	ENDIF
	481	!
[3]	482	! 0. Initialisation
	483	! -----------------
[632]	484
[3]	485	! local variable for debugging
	486	ll_print = .FALSE.
	487	ll_print = ll_print .AND. lwp
	488
	489	! Define frequency of output and means
	490	zdt = rdt
	491	IF( nacc == 1 ) zdt = rdtmin
[6487]	492	clop = "x" ! no use of the mask value (require less cpu time, and otherwise the model crashes)
[3]	493	#if defined key_diainstant
	494	zsto = nwrite * zdt
[1312]	495	clop = "inst("//TRIM(clop)//")"
[3]	496	#else
	497	zsto=zdt
[1312]	498	clop = "ave("//TRIM(clop)//")"
[3]	499	#endif
	500	zout = nwrite * zdt
	501	zmax = ( nitend - nit000 + 1 ) * zdt
	502
	503	! Define indices of the horizontal output zoom and vertical limit storage
	504	iimi = 1 ; iima = jpi
	505	ijmi = 1 ; ijma = jpj
	506	ipk = jpk
	507
	508	! define time axis
[1334]	509	it = kt
	510	itmod = kt - nit000 + 1
[3]	511
	512
	513	! 1. Define NETCDF files and fields at beginning of first time step
	514	! -----------------------------------------------------------------
	515
	516	IF( kt == nit000 ) THEN
	517
	518	! Define the NETCDF files (one per grid)
[632]	519
[3]	520	! Compute julian date from starting date of the run
[1309]	521	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
	522	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
[3]	523	IF(lwp)WRITE(numout,*)
	524	IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear, &
	525	& ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
	526	IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma, &
	527	' limit storage in depth = ', ipk
	528
	529	! WRITE root name in date.file for use by postpro
[1581]	530	IF(lwp) THEN
[895]	531	CALL dia_nam( clhstnam, nwrite,' ' )
[1581]	532	CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
[895]	533	WRITE(inum,*) clhstnam
	534	CLOSE(inum)
	535	ENDIF
[632]	536
[3]	537	! Define the T grid FILE ( nid_T )
[632]	538
[3]	539	CALL dia_nam( clhstnam, nwrite, 'grid_T' )
	540	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
	541	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
	542	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	543	& nit000-1, zjulian, zdt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )
[3]	544	CALL histvert( nid_T, "deptht", "Vertical T levels", & ! Vertical grid: gdept
[4292]	545	& "m", ipk, gdept_1d, nz_T, "down" )
[3]	546	! ! Index of ocean points
	547	CALL wheneq( jpijpjipk, tmask, 1, 1., ndex_T , ndim_T ) ! volume
	548	CALL wheneq( jpi*jpj , tmask, 1, 1., ndex_hT, ndim_hT ) ! surface
[3609]	549	!
	550	IF( ln_icebergs ) THEN
	551	!
	552	!! allocation cant go in dia_wri_alloc because ln_icebergs is only set after
	553	!! that routine is called from nemogcm, so do it here immediately before its needed
	554	ALLOCATE( ndex_bT(jpijpjnclasses), STAT=ierror )
	555	IF( lk_mpp ) CALL mpp_sum( ierror )
	556	IF( ierror /= 0 ) THEN
	557	CALL ctl_stop('dia_wri: failed to allocate iceberg diagnostic array')
	558	RETURN
	559	ENDIF
	560	!
	561	!! iceberg vertical coordinate is class number
	562	CALL histvert( nid_T, "class", "Iceberg class", & ! Vertical grid: class
	563	& "number", nclasses, class_num, nb_T )
	564	!
	565	!! each class just needs the surface index pattern
	566	ndim_bT = 3
	567	DO jn = 1,nclasses
	568	ndex_bT((jn-1)jpijpj+1:jnjpijpj) = ndex_hT(1:jpi*jpj)
	569	ENDDO
	570	!
	571	ENDIF
[3]	572
	573	! Define the U grid FILE ( nid_U )
	574
	575	CALL dia_nam( clhstnam, nwrite, 'grid_U' )
	576	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
	577	CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu, & ! Horizontal grid: glamu and gphiu
	578	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	579	& nit000-1, zjulian, zdt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )
[3]	580	CALL histvert( nid_U, "depthu", "Vertical U levels", & ! Vertical grid: gdept
[4292]	581	& "m", ipk, gdept_1d, nz_U, "down" )
[3]	582	! ! Index of ocean points
	583	CALL wheneq( jpijpjipk, umask, 1, 1., ndex_U , ndim_U ) ! volume
	584	CALL wheneq( jpi*jpj , umask, 1, 1., ndex_hU, ndim_hU ) ! surface
	585
	586	! Define the V grid FILE ( nid_V )
	587
	588	CALL dia_nam( clhstnam, nwrite, 'grid_V' ) ! filename
	589	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
	590	CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv, & ! Horizontal grid: glamv and gphiv
	591	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	592	& nit000-1, zjulian, zdt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )
[3]	593	CALL histvert( nid_V, "depthv", "Vertical V levels", & ! Vertical grid : gdept
[4292]	594	& "m", ipk, gdept_1d, nz_V, "down" )
[3]	595	! ! Index of ocean points
	596	CALL wheneq( jpijpjipk, vmask, 1, 1., ndex_V , ndim_V ) ! volume
	597	CALL wheneq( jpi*jpj , vmask, 1, 1., ndex_hV, ndim_hV ) ! surface
	598
	599	! Define the W grid FILE ( nid_W )
	600
	601	CALL dia_nam( clhstnam, nwrite, 'grid_W' ) ! filename
	602	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
	603	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
	604	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	605	& nit000-1, zjulian, zdt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set )
[3]	606	CALL histvert( nid_W, "depthw", "Vertical W levels", & ! Vertical grid: gdepw
[4292]	607	& "m", ipk, gdepw_1d, nz_W, "down" )
[3]	608
[632]	609
[3]	610	! Declare all the output fields as NETCDF variables
	611
	612	! !!! nid_T : 3D
	613	CALL histdef( nid_T, "votemper", "Temperature" , "C" , & ! tn
	614	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	615	CALL histdef( nid_T, "vosaline", "Salinity" , "PSU" , & ! sn
	616	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
[4292]	617	IF( lk_vvl ) THEN
	618	CALL histdef( nid_T, "vovvle3t", "Level thickness" , "m" ,& ! e3t_n
	619	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	620	CALL histdef( nid_T, "vovvldep", "T point depth" , "m" ,& ! e3t_n
	621	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	622	CALL histdef( nid_T, "vovvldef", "Squared level deformation" , "%^2" ,& ! e3t_n
	623	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	624	ENDIF
[3]	625	! !!! nid_T : 2D
	626	CALL histdef( nid_T, "sosstsst", "Sea Surface temperature" , "C" , & ! sst
	627	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	628	CALL histdef( nid_T, "sosaline", "Sea Surface Salinity" , "PSU" , & ! sss
	629	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[359]	630	CALL histdef( nid_T, "sossheig", "Sea Surface Height" , "m" , & ! ssh
[3]	631	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[2528]	632	CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux" , "Kg/m2/s", & ! (emp-rnf)
[3]	633	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[4570]	634	CALL histdef( nid_T, "sorunoff", "River runoffs" , "Kg/m2/s", & ! runoffs
	635	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[3625]	636	CALL histdef( nid_T, "sosfldow", "downward salt flux" , "PSU/m2/s", & ! sfx
[3]	637	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[4292]	638	IF( .NOT. lk_vvl ) THEN
	639	CALL histdef( nid_T, "sosst_cd", "Concentration/Dilution term on temperature" & ! emp * tsn(:,:,1,jp_tem)
[3625]	640	& , "KgC/m2/s", & ! sosst_cd
[4292]	641	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	642	CALL histdef( nid_T, "sosss_cd", "Concentration/Dilution term on salinity" & ! emp * tsn(:,:,1,jp_sal)
[3625]	643	& , "KgPSU/m2/s",& ! sosss_cd
[4292]	644	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	645	ENDIF
[888]	646	CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! qns + qsr
[3]	647	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	648	CALL histdef( nid_T, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! qsr
	649	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1585]	650	CALL histdef( nid_T, "somixhgt", "Turbocline Depth" , "m" , & ! hmld
	651	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[3]	652	CALL histdef( nid_T, "somxl010", "Mixed Layer Depth 0.01" , "m" , & ! hmlp
	653	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1037]	654	CALL histdef( nid_T, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
[3]	655	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1649]	656	CALL histdef( nid_T, "sowindsp", "wind speed at 10m" , "m/s" , & ! wndm
	657	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[3609]	658	!
	659	IF( ln_icebergs ) THEN
	660	CALL histdef( nid_T, "calving" , "calving mass input" , "kg/s" , &
	661	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	662	CALL histdef( nid_T, "calving_heat" , "calving heat flux" , "XXXX" , &
	663	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	664	CALL histdef( nid_T, "berg_floating_melt" , "Melt rate of icebergs + bits" , "kg/m2/s", &
	665	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	666	CALL histdef( nid_T, "berg_stored_ice" , "Accumulated ice mass by class" , "kg" , &
	667	& jpi, jpj, nh_T, nclasses , 1, nclasses , nb_T , 32, clop, zsto, zout )
	668	IF( ln_bergdia ) THEN
	669	CALL histdef( nid_T, "berg_melt" , "Melt rate of icebergs" , "kg/m2/s", &
	670	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	671	CALL histdef( nid_T, "berg_buoy_melt" , "Buoyancy component of iceberg melt rate" , "kg/m2/s", &
	672	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	673	CALL histdef( nid_T, "berg_eros_melt" , "Erosion component of iceberg melt rate" , "kg/m2/s", &
	674	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	675	CALL histdef( nid_T, "berg_conv_melt" , "Convective component of iceberg melt rate", "kg/m2/s", &
	676	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	677	CALL histdef( nid_T, "berg_virtual_area" , "Virtual coverage by icebergs" , "m2" , &
	678	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	679	CALL histdef( nid_T, "bits_src" , "Mass source of bergy bits" , "kg/m2/s", &
	680	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	681	CALL histdef( nid_T, "bits_melt" , "Melt rate of bergy bits" , "kg/m2/s", &
	682	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	683	CALL histdef( nid_T, "bits_mass" , "Bergy bit density field" , "kg/m2" , &
	684	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	685	CALL histdef( nid_T, "berg_mass" , "Iceberg density field" , "kg/m2" , &
	686	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	687	CALL histdef( nid_T, "berg_real_calving" , "Calving into iceberg class" , "kg/s" , &
	688	& jpi, jpj, nh_T, nclasses , 1, nclasses , nb_T , 32, clop, zsto, zout )
	689	ENDIF
	690	ENDIF
	691
[5407]	692	IF( .NOT. ln_cpl ) THEN
[4990]	693	CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping" , "W/m2" , & ! qrp
	694	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	695	CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping" , "Kg/m2/s", & ! erp
	696	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	697	CALL histdef( nid_T, "sosafldp", "Surface salt flux: damping" , "Kg/m2/s", & ! erp * sn
	698	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	699	ENDIF
[3]	700
[5407]	701	IF( ln_cpl .AND. nn_ice <= 1 ) THEN
[4990]	702	CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping" , "W/m2" , & ! qrp
	703	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	704	CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping" , "Kg/m2/s", & ! erp
	705	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	706	CALL histdef( nid_T, "sosafldp", "Surface salt flux: Damping" , "Kg/m2/s", & ! erp * sn
	707	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	708	ENDIF
	709
[3]	710	clmx ="l_max(only(x))" ! max index on a period
	711	CALL histdef( nid_T, "sobowlin", "Bowl Index" , "W-point", & ! bowl INDEX
	712	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clmx, zsto, zout )
	713	#if defined key_diahth
	714	CALL histdef( nid_T, "sothedep", "Thermocline Depth" , "m" , & ! hth
	715	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	716	CALL histdef( nid_T, "so20chgt", "Depth of 20C isotherm" , "m" , & ! hd20
	717	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	718	CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm" , "m" , & ! hd28
	719	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	720	CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , "W" , & ! htc3
	721	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	722	#endif
	723
[5407]	724	IF( ln_cpl .AND. nn_ice == 2 ) THEN
[4990]	725	CALL histdef( nid_T,"soicetem" , "Ice Surface Temperature" , "K" , & ! tn_ice
	726	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	727	CALL histdef( nid_T,"soicealb" , "Ice Albedo" , "[0,1]" , & ! alb_ice
	728	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	729	ENDIF
[3]	730
[2528]	731	CALL histend( nid_T, snc4chunks=snc4set )
[3]	732
	733	! !!! nid_U : 3D
	734	CALL histdef( nid_U, "vozocrtx", "Zonal Current" , "m/s" , & ! un
	735	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
[3294]	736	IF( ln_traldf_gdia ) THEN
	737	CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current" , "m/s" , & ! u_eiv
	738	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
	739	ELSE
[3]	740	#if defined key_diaeiv
[3294]	741	CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current" , "m/s" , & ! u_eiv
[3]	742	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
	743	#endif
[3294]	744	END IF
[3]	745	! !!! nid_U : 2D
[888]	746	CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau
[3]	747	& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
	748
[2528]	749	CALL histend( nid_U, snc4chunks=snc4set )
[3]	750
	751	! !!! nid_V : 3D
	752	CALL histdef( nid_V, "vomecrty", "Meridional Current" , "m/s" , & ! vn
	753	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
[3294]	754	IF( ln_traldf_gdia ) THEN
	755	CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current" , "m/s" , & ! v_eiv
	756	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
	757	ELSE
[3]	758	#if defined key_diaeiv
[3294]	759	CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current" , "m/s" , & ! v_eiv
[3]	760	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
	761	#endif
[3294]	762	END IF
[3]	763	! !!! nid_V : 2D
[888]	764	CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau
[3]	765	& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
	766
[2528]	767	CALL histend( nid_V, snc4chunks=snc4set )
[3]	768
	769	! !!! nid_W : 3D
	770	CALL histdef( nid_W, "vovecrtz", "Vertical Velocity" , "m/s" , & ! wn
	771	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[3294]	772	IF( ln_traldf_gdia ) THEN
	773	CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity" , "m/s" , & ! w_eiv
	774	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	775	ELSE
[3]	776	#if defined key_diaeiv
[3294]	777	CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity" , "m/s" , & ! w_eiv
	778	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[3]	779	#endif
[3294]	780	END IF
[3]	781	CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity" , "m2/s" , & ! avt
	782	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[255]	783	CALL histdef( nid_W, "votkeavm", "Vertical Eddy Viscosity" , "m2/s" , & ! avmu
	784	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	785
[3]	786	IF( lk_zdfddm ) THEN
	787	CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity" , "m2/s" , & ! avs
	788	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	789	ENDIF
	790	! !!! nid_W : 2D
	791	#if defined key_traldf_c2d
	792	CALL histdef( nid_W, "soleahtw", "lateral eddy diffusivity" , "m2/s" , & ! ahtw
	793	& jpi, jpj, nh_W, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
[23]	794	# if defined key_traldf_eiv
[3]	795	CALL histdef( nid_W, "soleaeiw", "eddy induced vel. coeff. at w-point", "m2/s", & ! aeiw
	796	& jpi, jpj, nh_W, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
[23]	797	# endif
[3]	798	#endif
	799
[2528]	800	CALL histend( nid_W, snc4chunks=snc4set )
[3]	801
	802	IF(lwp) WRITE(numout,*)
	803	IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
	804	IF(ll_print) CALL FLUSH(numout )
	805
	806	ENDIF
	807
	808	! 2. Start writing data
	809	! ---------------------
	810
[4292]	811	! ndex(1) est utilise ssi l'avant dernier argument est different de
[3]	812	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
	813	! donne le nombre d'elements, et ndex la liste des indices a sortir
	814
[1334]	815	IF( lwp .AND. MOD( itmod, nwrite ) == 0 ) THEN
[3]	816	WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
	817	WRITE(numout,*) '~~~~~~ '
	818	ENDIF
	819
[4292]	820	IF( lk_vvl ) THEN
[4492]	821	CALL histwrite( nid_T, "votemper", it, tsn(:,:,:,jp_tem) * fse3t_n(:,:,:) , ndim_T , ndex_T ) ! heat content
	822	CALL histwrite( nid_T, "vosaline", it, tsn(:,:,:,jp_sal) * fse3t_n(:,:,:) , ndim_T , ndex_T ) ! salt content
	823	CALL histwrite( nid_T, "sosstsst", it, tsn(:,:,1,jp_tem) * fse3t_n(:,:,1) , ndim_hT, ndex_hT ) ! sea surface heat content
	824	CALL histwrite( nid_T, "sosaline", it, tsn(:,:,1,jp_sal) * fse3t_n(:,:,1) , ndim_hT, ndex_hT ) ! sea surface salinity content
[4292]	825	ELSE
	826	CALL histwrite( nid_T, "votemper", it, tsn(:,:,:,jp_tem) , ndim_T , ndex_T ) ! temperature
	827	CALL histwrite( nid_T, "vosaline", it, tsn(:,:,:,jp_sal) , ndim_T , ndex_T ) ! salinity
	828	CALL histwrite( nid_T, "sosstsst", it, tsn(:,:,1,jp_tem) , ndim_hT, ndex_hT ) ! sea surface temperature
	829	CALL histwrite( nid_T, "sosaline", it, tsn(:,:,1,jp_sal) , ndim_hT, ndex_hT ) ! sea surface salinity
	830	ENDIF
	831	IF( lk_vvl ) THEN
	832	zw3d(:,:,:) = ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2
	833	CALL histwrite( nid_T, "vovvle3t", it, fse3t_n (:,:,:) , ndim_T , ndex_T ) ! level thickness
	834	CALL histwrite( nid_T, "vovvldep", it, fsdept_n(:,:,:) , ndim_T , ndex_T ) ! t-point depth
	835	CALL histwrite( nid_T, "vovvldef", it, zw3d , ndim_T , ndex_T ) ! level thickness deformation
	836	ENDIF
[3]	837	CALL histwrite( nid_T, "sossheig", it, sshn , ndim_hT, ndex_hT ) ! sea surface height
[2528]	838	CALL histwrite( nid_T, "sowaflup", it, ( emp-rnf ) , ndim_hT, ndex_hT ) ! upward water flux
[4570]	839	CALL histwrite( nid_T, "sorunoff", it, rnf , ndim_hT, ndex_hT ) ! river runoffs
[3625]	840	CALL histwrite( nid_T, "sosfldow", it, sfx , ndim_hT, ndex_hT ) ! downward salt flux
	841	! (includes virtual salt flux beneath ice
	842	! in linear free surface case)
[4292]	843	IF( .NOT. lk_vvl ) THEN
	844	zw2d(:,:) = emp (:,:) * tsn(:,:,1,jp_tem)
	845	CALL histwrite( nid_T, "sosst_cd", it, zw2d, ndim_hT, ndex_hT ) ! c/d term on sst
	846	zw2d(:,:) = emp (:,:) * tsn(:,:,1,jp_sal)
	847	CALL histwrite( nid_T, "sosss_cd", it, zw2d, ndim_hT, ndex_hT ) ! c/d term on sss
	848	ENDIF
[888]	849	CALL histwrite( nid_T, "sohefldo", it, qns + qsr , ndim_hT, ndex_hT ) ! total heat flux
[3]	850	CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux
[1585]	851	CALL histwrite( nid_T, "somixhgt", it, hmld , ndim_hT, ndex_hT ) ! turbocline depth
[3]	852	CALL histwrite( nid_T, "somxl010", it, hmlp , ndim_hT, ndex_hT ) ! mixed layer depth
[1037]	853	CALL histwrite( nid_T, "soicecov", it, fr_i , ndim_hT, ndex_hT ) ! ice fraction
[1649]	854	CALL histwrite( nid_T, "sowindsp", it, wndm , ndim_hT, ndex_hT ) ! wind speed
[3609]	855	!
	856	IF( ln_icebergs ) THEN
	857	!
	858	CALL histwrite( nid_T, "calving" , it, berg_grid%calving , ndim_hT, ndex_hT )
	859	CALL histwrite( nid_T, "calving_heat" , it, berg_grid%calving_hflx , ndim_hT, ndex_hT )
	860	CALL histwrite( nid_T, "berg_floating_melt" , it, berg_grid%floating_melt, ndim_hT, ndex_hT )
	861	!
	862	CALL histwrite( nid_T, "berg_stored_ice" , it, berg_grid%stored_ice , ndim_bT, ndex_bT )
	863	!
	864	IF( ln_bergdia ) THEN
	865	CALL histwrite( nid_T, "berg_melt" , it, berg_melt , ndim_hT, ndex_hT )
	866	CALL histwrite( nid_T, "berg_buoy_melt" , it, buoy_melt , ndim_hT, ndex_hT )
	867	CALL histwrite( nid_T, "berg_eros_melt" , it, eros_melt , ndim_hT, ndex_hT )
	868	CALL histwrite( nid_T, "berg_conv_melt" , it, conv_melt , ndim_hT, ndex_hT )
	869	CALL histwrite( nid_T, "berg_virtual_area" , it, virtual_area , ndim_hT, ndex_hT )
	870	CALL histwrite( nid_T, "bits_src" , it, bits_src , ndim_hT, ndex_hT )
	871	CALL histwrite( nid_T, "bits_melt" , it, bits_melt , ndim_hT, ndex_hT )
	872	CALL histwrite( nid_T, "bits_mass" , it, bits_mass , ndim_hT, ndex_hT )
	873	CALL histwrite( nid_T, "berg_mass" , it, berg_mass , ndim_hT, ndex_hT )
	874	!
	875	CALL histwrite( nid_T, "berg_real_calving" , it, real_calving , ndim_bT, ndex_bT )
	876	ENDIF
	877	ENDIF
	878
[5407]	879	IF( .NOT. ln_cpl ) THEN
[4990]	880	CALL histwrite( nid_T, "sohefldp", it, qrp , ndim_hT, ndex_hT ) ! heat flux damping
	881	CALL histwrite( nid_T, "sowafldp", it, erp , ndim_hT, ndex_hT ) ! freshwater flux damping
[3294]	882	IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1)
[4990]	883	CALL histwrite( nid_T, "sosafldp", it, zw2d , ndim_hT, ndex_hT ) ! salt flux damping
	884	ENDIF
[5407]	885	IF( ln_cpl .AND. nn_ice <= 1 ) THEN
[4990]	886	CALL histwrite( nid_T, "sohefldp", it, qrp , ndim_hT, ndex_hT ) ! heat flux damping
	887	CALL histwrite( nid_T, "sowafldp", it, erp , ndim_hT, ndex_hT ) ! freshwater flux damping
	888	IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1)
	889	CALL histwrite( nid_T, "sosafldp", it, zw2d , ndim_hT, ndex_hT ) ! salt flux damping
	890	ENDIF
	891	! zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)
	892	! CALL histwrite( nid_T, "sobowlin", it, zw2d , ndim_hT, ndex_hT ) ! ???
[3]	893
	894	#if defined key_diahth
	895	CALL histwrite( nid_T, "sothedep", it, hth , ndim_hT, ndex_hT ) ! depth of the thermocline
	896	CALL histwrite( nid_T, "so20chgt", it, hd20 , ndim_hT, ndex_hT ) ! depth of the 20 isotherm
	897	CALL histwrite( nid_T, "so28chgt", it, hd28 , ndim_hT, ndex_hT ) ! depth of the 28 isotherm
	898	CALL histwrite( nid_T, "sohtc300", it, htc3 , ndim_hT, ndex_hT ) ! first 300m heaat content
	899	#endif
[888]	900
[5407]	901	IF( ln_cpl .AND. nn_ice == 2 ) THEN
[4990]	902	CALL histwrite( nid_T, "soicetem", it, tn_ice(:,:,1) , ndim_hT, ndex_hT ) ! surf. ice temperature
	903	CALL histwrite( nid_T, "soicealb", it, alb_ice(:,:,1), ndim_hT, ndex_hT ) ! ice albedo
	904	ENDIF
	905
[3]	906	CALL histwrite( nid_U, "vozocrtx", it, un , ndim_U , ndex_U ) ! i-current
[3294]	907	IF( ln_traldf_gdia ) THEN
	908	IF (.not. ALLOCATED(psix_eiv))THEN
	909	ALLOCATE( psix_eiv(jpi,jpj,jpk) , psiy_eiv(jpi,jpj,jpk) , STAT=ierr )
	910	IF( lk_mpp ) CALL mpp_sum ( ierr )
	911	IF( ierr > 0 ) CALL ctl_stop('STOP', 'diawri: unable to allocate psi{x,y}_eiv')
	912	psix_eiv(:,:,:) = 0.0_wp
	913	psiy_eiv(:,:,:) = 0.0_wp
	914	ENDIF
	915	DO jk=1,jpkm1
	916	zw3d(:,:,jk) = (psix_eiv(:,:,jk+1) - psix_eiv(:,:,jk))/fse3u(:,:,jk) ! u_eiv = -dpsix/dz
	917	END DO
	918	zw3d(:,:,jpk) = 0._wp
	919	CALL histwrite( nid_U, "vozoeivu", it, zw3d, ndim_U , ndex_U ) ! i-eiv current
	920	ELSE
[3]	921	#if defined key_diaeiv
[3294]	922	CALL histwrite( nid_U, "vozoeivu", it, u_eiv, ndim_U , ndex_U ) ! i-eiv current
[3]	923	#endif
[3294]	924	ENDIF
[888]	925	CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress
[3]	926
	927	CALL histwrite( nid_V, "vomecrty", it, vn , ndim_V , ndex_V ) ! j-current
[3294]	928	IF( ln_traldf_gdia ) THEN
	929	DO jk=1,jpk-1
	930	zw3d(:,:,jk) = (psiy_eiv(:,:,jk+1) - psiy_eiv(:,:,jk))/fse3v(:,:,jk) ! v_eiv = -dpsiy/dz
	931	END DO
	932	zw3d(:,:,jpk) = 0._wp
	933	CALL histwrite( nid_V, "vomeeivv", it, zw3d, ndim_V , ndex_V ) ! j-eiv current
	934	ELSE
[3]	935	#if defined key_diaeiv
[3294]	936	CALL histwrite( nid_V, "vomeeivv", it, v_eiv, ndim_V , ndex_V ) ! j-eiv current
[3]	937	#endif
[3294]	938	ENDIF
[888]	939	CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress
[3]	940
	941	CALL histwrite( nid_W, "vovecrtz", it, wn , ndim_T, ndex_T ) ! vert. current
[3294]	942	IF( ln_traldf_gdia ) THEN
	943	DO jk=1,jpk-1
	944	DO jj = 2, jpjm1
	945	DO ji = fs_2, fs_jpim1 ! vector opt.
	946	zw3d(ji,jj,jk) = (psiy_eiv(ji,jj,jk) - psiy_eiv(ji,jj-1,jk))/e2v(ji,jj) + &
	947	& (psix_eiv(ji,jj,jk) - psix_eiv(ji-1,jj,jk))/e1u(ji,jj) ! w_eiv = dpsiy/dy + dpsiy/dx
	948	END DO
	949	END DO
	950	END DO
	951	zw3d(:,:,jpk) = 0._wp
	952	CALL histwrite( nid_W, "voveeivw", it, zw3d , ndim_T, ndex_T ) ! vert. eiv current
	953	ELSE
[3]	954	# if defined key_diaeiv
[3294]	955	CALL histwrite( nid_W, "voveeivw", it, w_eiv , ndim_T, ndex_T ) ! vert. eiv current
[3]	956	# endif
[3294]	957	ENDIF
[3]	958	CALL histwrite( nid_W, "votkeavt", it, avt , ndim_T, ndex_T ) ! T vert. eddy diff. coef.
[255]	959	CALL histwrite( nid_W, "votkeavm", it, avmu , ndim_T, ndex_T ) ! T vert. eddy visc. coef.
[3]	960	IF( lk_zdfddm ) THEN
	961	CALL histwrite( nid_W, "voddmavs", it, fsavs(:,:,:), ndim_T, ndex_T ) ! S vert. eddy diff. coef.
	962	ENDIF
	963	#if defined key_traldf_c2d
	964	CALL histwrite( nid_W, "soleahtw", it, ahtw , ndim_hT, ndex_hT ) ! lateral eddy diff. coef.
[23]	965	# if defined key_traldf_eiv
[3]	966	CALL histwrite( nid_W, "soleaeiw", it, aeiw , ndim_hT, ndex_hT ) ! EIV coefficient at w-point
[23]	967	# endif
[3]	968	#endif
	969
[1318]	970	! 3. Close all files
	971	! ---------------------------------------
[1561]	972	IF( kt == nitend ) THEN
[3]	973	CALL histclo( nid_T )
	974	CALL histclo( nid_U )
	975	CALL histclo( nid_V )
	976	CALL histclo( nid_W )
	977	ENDIF
[2528]	978	!
[3294]	979	CALL wrk_dealloc( jpi , jpj , zw2d )
[4292]	980	IF ( ln_traldf_gdia .OR. lk_vvl ) call wrk_dealloc( jpi , jpj , jpk , zw3d )
[2715]	981	!
[3294]	982	IF( nn_timing == 1 ) CALL timing_stop('dia_wri')
	983	!
[3]	984	END SUBROUTINE dia_wri
[1482]	985	# endif
[3]	986
[1567]	987	#endif
	988
[1334]	989	SUBROUTINE dia_wri_state( cdfile_name, kt )
[3]	990	!!---------------------------------------------------------------------
	991	!! * ROUTINE dia_wri_state *
	992	!!
	993	!! ** Purpose : create a NetCDF file named cdfile_name which contains
	994	!! the instantaneous ocean state and forcing fields.
	995	!! Used to find errors in the initial state or save the last
	996	!! ocean state in case of abnormal end of a simulation
	997	!!
	998	!! ** Method : NetCDF files using ioipsl
	999	!! File 'output.init.nc' is created if ninist = 1 (namelist)
	1000	!! File 'output.abort.nc' is created in case of abnormal job end
	1001	!!----------------------------------------------------------------------
[1334]	1002	CHARACTER (len=* ), INTENT( in ) :: cdfile_name ! name of the file created
	1003	INTEGER , INTENT( in ) :: kt ! ocean time-step index
[2528]	1004	!!
[648]	1005	CHARACTER (len=32) :: clname
[3]	1006	CHARACTER (len=40) :: clop
[2528]	1007	INTEGER :: id_i , nz_i, nh_i
	1008	INTEGER, DIMENSION(1) :: idex ! local workspace
	1009	REAL(wp) :: zsto, zout, zmax, zjulian, zdt
[3]	1010	!!----------------------------------------------------------------------
[3294]	1011	!
[3625]	1012	! IF( nn_timing == 1 ) CALL timing_start('dia_wri_state') ! not sure this works for routines not called in first timestep
[3]	1013
	1014	! 0. Initialisation
	1015	! -----------------
[632]	1016
[648]	1017	! Define name, frequency of output and means
	1018	clname = cdfile_name
[1792]	1019	IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
[3]	1020	zdt = rdt
	1021	zsto = rdt
	1022	clop = "inst(x)" ! no use of the mask value (require less cpu time)
	1023	zout = rdt
	1024	zmax = ( nitend - nit000 + 1 ) * zdt
	1025
[648]	1026	IF(lwp) WRITE(numout,*)
	1027	IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
	1028	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ and forcing fields file created '
	1029	IF(lwp) WRITE(numout,*) ' and named :', clname, '.nc'
	1030
	1031
[3]	1032	! 1. Define NETCDF files and fields at beginning of first time step
	1033	! -----------------------------------------------------------------
	1034
	1035	! Compute julian date from starting date of the run
[1310]	1036	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian ) ! time axis
	1037	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
[648]	1038	CALL histbeg( clname, jpi, glamt, jpj, gphit, &
[2528]	1039	1, jpi, 1, jpj, nit000-1, zjulian, zdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
[3]	1040	CALL histvert( id_i, "deptht", "Vertical T levels", & ! Vertical grid : gdept
[4292]	1041	"m", jpk, gdept_1d, nz_i, "down")
[3]	1042
	1043	! Declare all the output fields as NetCDF variables
	1044
	1045	CALL histdef( id_i, "vosaline", "Salinity" , "PSU" , & ! salinity
	1046	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	1047	CALL histdef( id_i, "votemper", "Temperature" , "C" , & ! temperature
	1048	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
[359]	1049	CALL histdef( id_i, "sossheig", "Sea Surface Height" , "m" , & ! ssh
	1050	& jpi, jpj, nh_i, 1 , 1, 1 , nz_i, 32, clop, zsto, zout )
[3]	1051	CALL histdef( id_i, "vozocrtx", "Zonal Current" , "m/s" , & ! zonal current
	1052	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	1053	CALL histdef( id_i, "vomecrty", "Meridional Current" , "m/s" , & ! meridonal current
	1054	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	1055	CALL histdef( id_i, "vovecrtz", "Vertical Velocity" , "m/s" , & ! vertical current
	1056	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	1057	CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S", & ! net freshwater
	1058	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	1059	CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2" , & ! net heat flux
	1060	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	1061	CALL histdef( id_i, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! solar flux
	1062	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1037]	1063	CALL histdef( id_i, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
[3]	1064	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	1065	CALL histdef( id_i, "sozotaux", "Zonal Wind Stress" , "N/m2" , & ! i-wind stress
	1066	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	1067	CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2" , & ! j-wind stress
	1068	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[4292]	1069	IF( lk_vvl ) THEN
	1070	CALL histdef( id_i, "vovvldep", "T point depth" , "m" , & ! t-point depth
	1071	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
[6487]	1072	CALL histdef( id_i, "vovvle3t", "T point thickness" , "m" , & ! t-point depth
	1073	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
[4292]	1074	END IF
[3]	1075
[1482]	1076	#if defined key_lim2
	1077	CALL lim_wri_state_2( kt, id_i, nh_i )
[4161]	1078	#elif defined key_lim3
	1079	CALL lim_wri_state( kt, id_i, nh_i )
[1482]	1080	#else
[2528]	1081	CALL histend( id_i, snc4chunks=snc4set )
[1482]	1082	#endif
[3]	1083
	1084	! 2. Start writing data
	1085	! ---------------------
	1086	! idex(1) est utilise ssi l'avant dernier argument est diffferent de
	1087	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
	1088	! donne le nombre d'elements, et idex la liste des indices a sortir
	1089	idex(1) = 1 ! init to avoid compil warning
[632]	1090
[3]	1091	! Write all fields on T grid
[3294]	1092	CALL histwrite( id_i, "votemper", kt, tsn(:,:,:,jp_tem), jpijpjjpk, idex ) ! now temperature
	1093	CALL histwrite( id_i, "vosaline", kt, tsn(:,:,:,jp_sal), jpijpjjpk, idex ) ! now salinity
	1094	CALL histwrite( id_i, "sossheig", kt, sshn , jpi*jpj , idex ) ! sea surface height
	1095	CALL histwrite( id_i, "vozocrtx", kt, un , jpijpjjpk, idex ) ! now i-velocity
	1096	CALL histwrite( id_i, "vomecrty", kt, vn , jpijpjjpk, idex ) ! now j-velocity
	1097	CALL histwrite( id_i, "vovecrtz", kt, wn , jpijpjjpk, idex ) ! now k-velocity
	1098	CALL histwrite( id_i, "sowaflup", kt, (emp-rnf ) , jpi*jpj , idex ) ! freshwater budget
	1099	CALL histwrite( id_i, "sohefldo", kt, qsr + qns , jpi*jpj , idex ) ! total heat flux
	1100	CALL histwrite( id_i, "soshfldo", kt, qsr , jpi*jpj , idex ) ! solar heat flux
	1101	CALL histwrite( id_i, "soicecov", kt, fr_i , jpi*jpj , idex ) ! ice fraction
	1102	CALL histwrite( id_i, "sozotaux", kt, utau , jpi*jpj , idex ) ! i-wind stress
	1103	CALL histwrite( id_i, "sometauy", kt, vtau , jpi*jpj , idex ) ! j-wind stress
[6487]	1104	IF( lk_vvl ) THEN
	1105	CALL histwrite( id_i, "vovvldep", kt, fsdept_n(:,:,:), jpijpjjpk, idex )! T-cell depth
	1106	CALL histwrite( id_i, "vovvle3t", kt, fse3t_n (:,:,:), jpijpjjpk, idex )! T-cell thickness
	1107	END IF
[3]	1108
	1109	! 3. Close the file
	1110	! -----------------
	1111	CALL histclo( id_i )
[1567]	1112	#if ! defined key_iomput && ! defined key_dimgout
[1561]	1113	IF( ninist /= 1 ) THEN
	1114	CALL histclo( nid_T )
	1115	CALL histclo( nid_U )
	1116	CALL histclo( nid_V )
	1117	CALL histclo( nid_W )
	1118	ENDIF
	1119	#endif
[8280]	1120
	1121	IF (cdfile_name == "output.abort") THEN
	1122	CALL ctl_stop('MPPSTOP', 'NEMO abort from dia_wri_state')
	1123	END IF
[3294]	1124
[3625]	1125	! IF( nn_timing == 1 ) CALL timing_stop('dia_wri_state') ! not sure this works for routines not called in first timestep
[3294]	1126	!
[3]	1127
	1128	END SUBROUTINE dia_wri_state
	1129	!!======================================================================
	1130	END MODULE diawri

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

Download in other formats: