Context Navigation

diawri.F90 @ 9830

Last change on this file since 9830 was 9830, checked in by frrh, 6 years ago

Merge revisions 9607:9721 of/branches/UKMO/dev_r5518_GO6_diag_bitcomp
into GO6 package branch.

This change ensures most 2D and 3D diagnostics produced by NEMO and MEDUSA
are bit reproducible on different PE decompositions.

Command used:
svn merge -r 9607:9721 svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/UKMO/dev_r5518_GO6_diag_bitcomp

Met Office GMED ticket 389 refers.

This change applies a mask to all duplicate grid points
output on diagnistic grids for T, U and V points. i.e. it masks
any wrap columns and duplicated grid points across the N-fold.
Fields affected are all "standard" NEMO diagnostics (scalar and
diaptr diagnostics are not on "normal" grids).

It also introduces some corrections/initialisations to achieve
PE decomposition bit comparison.

Most 2D or 3D fields are now bit comparable on different PE decompositions.
Only diaptr diagnostics can not be guaranteed bit reproducible
(due to their method of computation).

This change does nothing to CICE output.

Model evolution is unaffected.

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

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

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

Context Navigation

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

Download in other formats: