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

Context Navigation

source: branches/UKMO/AMM15_v3_6_STABLE_package/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 @ 8059

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

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

Download in other formats: