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 NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/OCE/DIA – NEMO

Context Navigation

source: NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/OCE/DIA/diawri.F90 @ 13942

Last change on this file since 13942 was 13942, checked in by jchanut, 4 years ago
Merged with trunk @13787
Property svn:keywords set to `Id`
File size: 59.0 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
[5836]	19	!! 3.7 ! 2014-01 (G. Madec) remove eddy induced velocity from no-IOM output
	20	!! ! change name of output variables in dia_wri_state
[2528]	21	!!----------------------------------------------------------------------
[3]	22
	23	!!----------------------------------------------------------------------
[2528]	24	!! dia_wri : create the standart output files
	25	!! dia_wri_state : create an output NetCDF file for a single instantaeous ocean state and forcing fields
	26	!!----------------------------------------------------------------------
[9019]	27	USE oce ! ocean dynamics and tracers
[12377]	28	USE isf_oce
	29	USE isfcpl
	30	USE abl ! abl variables in case ln_abl = .true.
[9019]	31	USE dom_oce ! ocean space and time domain
	32	USE phycst ! physical constants
	33	USE dianam ! build name of file (routine)
	34	USE diahth ! thermocline diagnostics
	35	USE dynadv , ONLY: ln_dynadv_vec
	36	USE icb_oce ! Icebergs
	37	USE icbdia ! Iceberg budgets
	38	USE ldftra ! lateral physics: eddy diffusivity coef.
	39	USE ldfdyn ! lateral physics: eddy viscosity coef.
	40	USE sbc_oce ! Surface boundary condition: ocean fields
	41	USE sbc_ice ! Surface boundary condition: ice fields
	42	USE sbcssr ! restoring term toward SST/SSS climatology
	43	USE sbcwave ! wave parameters
	44	USE wet_dry ! wetting and drying
	45	USE zdf_oce ! ocean vertical physics
	46	USE zdfdrg ! ocean vertical physics: top/bottom friction
	47	USE zdfmxl ! mixed layer
[6140]	48	!
[9019]	49	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	50	USE in_out_manager ! I/O manager
	51	USE dia25h ! 25h Mean output
	52	USE iom !
	53	USE ioipsl !
[5463]	54
[9570]	55	#if defined key_si3
[10425]	56	USE ice
[9019]	57	USE icewri
[1482]	58	#endif
[2715]	59	USE lib_mpp ! MPP library
[3294]	60	USE timing ! preformance summary
[12377]	61	USE diu_bulk ! diurnal warm layer
	62	USE diu_coolskin ! Cool skin
[2528]	63
[3]	64	IMPLICIT NONE
	65	PRIVATE
	66
[2528]	67	PUBLIC dia_wri ! routines called by step.F90
	68	PUBLIC dia_wri_state
[2715]	69	PUBLIC dia_wri_alloc ! Called by nemogcm module
[12377]	70	#if ! defined key_iomput
	71	PUBLIC dia_wri_alloc_abl ! Called by sbcabl module (if ln_abl = .true.)
	72	#endif
[2528]	73	INTEGER :: nid_T, nz_T, nh_T, ndim_T, ndim_hT ! grid_T file
[3609]	74	INTEGER :: nb_T , ndim_bT ! grid_T file
[2528]	75	INTEGER :: nid_U, nz_U, nh_U, ndim_U, ndim_hU ! grid_U file
	76	INTEGER :: nid_V, nz_V, nh_V, ndim_V, ndim_hV ! grid_V file
	77	INTEGER :: nid_W, nz_W, nh_W ! grid_W file
[12377]	78	INTEGER :: nid_A, nz_A, nh_A, ndim_A, ndim_hA ! grid_ABL file
[2528]	79	INTEGER :: ndex(1) ! ???
[2715]	80	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
[12377]	81	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hA, ndex_A ! ABL
[2715]	82	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V
[3609]	83	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_bT
[3]	84
	85	!! * Substitutions
[12377]	86	# include "do_loop_substitute.h90"
[13237]	87	# include "domzgr_substitute.h90"
[3]	88	!!----------------------------------------------------------------------
[9598]	89	!! NEMO/OCE 4.0 , NEMO Consortium (2018)
[5217]	90	!! $Id$
[10068]	91	!! Software governed by the CeCILL license (see ./LICENSE)
[3]	92	!!----------------------------------------------------------------------
	93	CONTAINS
	94
[6140]	95	#if defined key_iomput
[3]	96	!!----------------------------------------------------------------------
[2528]	97	!! 'key_iomput' use IOM library
	98	!!----------------------------------------------------------------------
[9652]	99	INTEGER FUNCTION dia_wri_alloc()
	100	!
	101	dia_wri_alloc = 0
	102	!
	103	END FUNCTION dia_wri_alloc
[2715]	104
[9652]	105
[12377]	106	SUBROUTINE dia_wri( kt, Kmm )
[1482]	107	!!---------------------------------------------------------------------
	108	!! * ROUTINE dia_wri *
	109	!!
	110	!! ** Purpose : Standard output of opa: dynamics and tracer fields
	111	!! NETCDF format is used by default
	112	!!
	113	!! ** Method : use iom_put
	114	!!----------------------------------------------------------------------
	115	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[12377]	116	INTEGER, INTENT( in ) :: Kmm ! ocean time level index
[1756]	117	!!
[9019]	118	INTEGER :: ji, jj, jk ! dummy loop indices
	119	INTEGER :: ikbot ! local integer
[13942]	120	REAL(wp):: ze3
[9019]	121	REAL(wp):: zztmp , zztmpx ! local scalar
	122	REAL(wp):: zztmp2, zztmpy ! - -
	123	REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace
	124	REAL(wp), DIMENSION(jpi,jpj,jpk) :: z3d ! 3D workspace
[1482]	125	!!----------------------------------------------------------------------
	126	!
[9124]	127	IF( ln_timing ) CALL timing_start('dia_wri')
[3294]	128	!
[1482]	129	! Output the initial state and forcings
	130	IF( ninist == 1 ) THEN
[12377]	131	CALL dia_wri_state( Kmm, 'output.init' )
[1482]	132	ninist = 0
	133	ENDIF
[3]	134
[6351]	135	! Output of initial vertical scale factor
	136	CALL iom_put("e3t_0", e3t_0(:,:,:) )
[10114]	137	CALL iom_put("e3u_0", e3u_0(:,:,:) )
	138	CALL iom_put("e3v_0", e3v_0(:,:,:) )
[6351]	139	!
[13565]	140	IF ( iom_use("tpt_dep") ) THEN
	141	DO jk = 1, jpk
	142	z3d(:,:,jk) = gdept(:,:,jk,Kmm)
	143	END DO
	144	CALL iom_put( "tpt_dep", z3d(:,:,:) )
	145	ENDIF
	146
[13237]	147	IF ( iom_use("e3t") .OR. iom_use("e3tdef") ) THEN ! time-varying e3t
	148	DO jk = 1, jpk
	149	z3d(:,:,jk) = e3t(:,:,jk,Kmm)
	150	END DO
	151	CALL iom_put( "e3t" , z3d(:,:,:) )
	152	CALL iom_put( "e3tdef" , ( ( z3d(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 )
	153	ENDIF
	154	IF ( iom_use("e3u") ) THEN ! time-varying e3u
	155	DO jk = 1, jpk
	156	z3d(:,:,jk) = e3u(:,:,jk,Kmm)
	157	END DO
	158	CALL iom_put( "e3u" , z3d(:,:,:) )
	159	ENDIF
	160	IF ( iom_use("e3v") ) THEN ! time-varying e3v
	161	DO jk = 1, jpk
	162	z3d(:,:,jk) = e3v(:,:,jk,Kmm)
	163	END DO
	164	CALL iom_put( "e3v" , z3d(:,:,:) )
	165	ENDIF
	166	IF ( iom_use("e3w") ) THEN ! time-varying e3w
	167	DO jk = 1, jpk
	168	z3d(:,:,jk) = e3w(:,:,jk,Kmm)
	169	END DO
	170	CALL iom_put( "e3w" , z3d(:,:,:) )
	171	ENDIF
[5461]	172
[13237]	173	IF( ll_wd ) THEN ! sea surface height (brought back to the reference used for wetting and drying)
	174	CALL iom_put( "ssh" , (ssh(:,:,Kmm)+ssh_ref)*tmask(:,:,1) )
[9023]	175	ELSE
[12377]	176	CALL iom_put( "ssh" , ssh(:,:,Kmm) ) ! sea surface height
[9023]	177	ENDIF
	178
[7646]	179	IF( iom_use("wetdep") ) & ! wet depth
[12377]	180	CALL iom_put( "wetdep" , ht_0(:,:) + ssh(:,:,Kmm) )
[5107]	181
[12377]	182	CALL iom_put( "toce", ts(:,:,:,jp_tem,Kmm) ) ! 3D temperature
	183	CALL iom_put( "sst", ts(:,:,1,jp_tem,Kmm) ) ! surface temperature
[5107]	184	IF ( iom_use("sbt") ) THEN
[13942]	185	DO_2D( 0, 0, 0, 0 )
[12377]	186	ikbot = mbkt(ji,jj)
	187	z2d(ji,jj) = ts(ji,jj,ikbot,jp_tem,Kmm)
	188	END_2D
[5107]	189	CALL iom_put( "sbt", z2d ) ! bottom temperature
	190	ENDIF
	191
[12377]	192	CALL iom_put( "soce", ts(:,:,:,jp_sal,Kmm) ) ! 3D salinity
	193	CALL iom_put( "sss", ts(:,:,1,jp_sal,Kmm) ) ! surface salinity
[5107]	194	IF ( iom_use("sbs") ) THEN
[13942]	195	DO_2D( 0, 0, 0, 0 )
[12377]	196	ikbot = mbkt(ji,jj)
	197	z2d(ji,jj) = ts(ji,jj,ikbot,jp_sal,Kmm)
	198	END_2D
[5107]	199	CALL iom_put( "sbs", z2d ) ! bottom salinity
	200	ENDIF
[5463]	201
[13237]	202	#if ! defined key_qco
[13089]	203	CALL iom_put( "rhop", rhop(:,:,:) ) ! 3D potential density (sigma0)
[13237]	204	#endif
[13089]	205
[5463]	206	IF ( iom_use("taubot") ) THEN ! bottom stress
[12489]	207	zztmp = rho0 * 0.25
[7753]	208	z2d(:,:) = 0._wp
[13295]	209	DO_2D( 0, 0, 0, 0 )
[12377]	210	zztmp2 = ( ( rCdU_bot(ji+1,jj)+rCdU_bot(ji ,jj) ) * uu(ji ,jj,mbku(ji ,jj),Kmm) )**2 &
	211	& + ( ( rCdU_bot(ji ,jj)+rCdU_bot(ji-1,jj) ) * uu(ji-1,jj,mbku(ji-1,jj),Kmm) )**2 &
	212	& + ( ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj ) ) * vv(ji,jj ,mbkv(ji,jj ),Kmm) )**2 &
	213	& + ( ( rCdU_bot(ji,jj )+rCdU_bot(ji,jj-1) ) * vv(ji,jj-1,mbkv(ji,jj-1),Kmm) )**2
	214	z2d(ji,jj) = zztmp * SQRT( zztmp2 ) * tmask(ji,jj,1)
	215	!
	216	END_2D
[5463]	217	CALL iom_put( "taubot", z2d )
	218	ENDIF
[5107]	219
[12377]	220	CALL iom_put( "uoce", uu(:,:,:,Kmm) ) ! 3D i-current
	221	CALL iom_put( "ssu", uu(:,:,1,Kmm) ) ! surface i-current
[5107]	222	IF ( iom_use("sbu") ) THEN
[13942]	223	DO_2D( 0, 0, 0, 0 )
[12377]	224	ikbot = mbku(ji,jj)
	225	z2d(ji,jj) = uu(ji,jj,ikbot,Kmm)
	226	END_2D
[5107]	227	CALL iom_put( "sbu", z2d ) ! bottom i-current
	228	ENDIF
	229
[12377]	230	CALL iom_put( "voce", vv(:,:,:,Kmm) ) ! 3D j-current
	231	CALL iom_put( "ssv", vv(:,:,1,Kmm) ) ! surface j-current
[5107]	232	IF ( iom_use("sbv") ) THEN
[13942]	233	DO_2D( 0, 0, 0, 0 )
[12377]	234	ikbot = mbkv(ji,jj)
	235	z2d(ji,jj) = vv(ji,jj,ikbot,Kmm)
	236	END_2D
[5107]	237	CALL iom_put( "sbv", z2d ) ! bottom j-current
[4990]	238	ENDIF
[1482]	239
[12377]	240	IF( ln_zad_Aimp ) ww = ww + wi ! Recombine explicit and implicit parts of vertical velocity for diagnostic output
	241	CALL iom_put( "woce", ww ) ! vertical velocity
[13237]	242
[5461]	243	IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN ! vertical mass transport & its square value
	244	! Caution: in the VVL case, it only correponds to the baroclinic mass transport.
[12489]	245	z2d(:,:) = rho0 * e1e2t(:,:)
[5461]	246	DO jk = 1, jpk
[12377]	247	z3d(:,:,jk) = ww(:,:,jk) * z2d(:,:)
[5461]	248	END DO
	249	CALL iom_put( "w_masstr" , z3d )
	250	IF( iom_use('w_masstr2') ) CALL iom_put( "w_masstr2", z3d(:,:,:) * z3d(:,:,:) )
	251	ENDIF
[11418]	252	!
[12377]	253	IF( ln_zad_Aimp ) ww = ww - wi ! Remove implicit part of vertical velocity that was added for diagnostic output
[5461]	254
[9019]	255	CALL iom_put( "avt" , avt ) ! T vert. eddy diff. coef.
	256	CALL iom_put( "avs" , avs ) ! S vert. eddy diff. coef.
	257	CALL iom_put( "avm" , avm ) ! T vert. eddy visc. coef.
[5107]	258
[9019]	259	IF( iom_use('logavt') ) CALL iom_put( "logavt", LOG( MAX( 1.e-20_wp, avt(:,:,:) ) ) )
	260	IF( iom_use('logavs') ) CALL iom_put( "logavs", LOG( MAX( 1.e-20_wp, avs(:,:,:) ) ) )
[6351]	261
[13942]	262	IF ( iom_use("socegrad") .OR. iom_use("socegrad2") ) THEN
	263	z3d(:,:,jpk) = 0.
	264	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
	265	zztmp = ts(ji,jj,jk,jp_sal,Kmm)
	266	zztmpx = (ts(ji+1,jj,jk,jp_sal,Kmm) - zztmp) * r1_e1u(ji,jj) + (zztmp - ts(ji-1,jj ,jk,jp_sal,Kmm)) * r1_e1u(ji-1,jj)
	267	zztmpy = (ts(ji,jj+1,jk,jp_sal,Kmm) - zztmp) * r1_e2v(ji,jj) + (zztmp - ts(ji ,jj-1,jk,jp_sal,Kmm)) * r1_e2v(ji,jj-1)
	268	z3d(ji,jj,jk) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy ) &
	269	& * umask(ji,jj,jk) * umask(ji-1,jj,jk) * vmask(ji,jj,jk) * umask(ji,jj-1,jk)
	270	END_3D
	271	CALL iom_put( "socegrad2", z3d ) ! square of module of sal gradient
	272	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
	273	z3d(ji,jj,jk) = SQRT( z3d(ji,jj,jk) )
	274	END_3D
	275	CALL iom_put( "socegrad" , z3d ) ! module of sal gradient
	276	ENDIF
	277
[5107]	278	IF ( iom_use("sstgrad") .OR. iom_use("sstgrad2") ) THEN
[13942]	279	DO_2D( 0, 0, 0, 0 ) ! sst gradient
[12377]	280	zztmp = ts(ji,jj,1,jp_tem,Kmm)
	281	zztmpx = ( ts(ji+1,jj,1,jp_tem,Kmm) - zztmp ) * r1_e1u(ji,jj) + ( zztmp - ts(ji-1,jj ,1,jp_tem,Kmm) ) * r1_e1u(ji-1,jj)
	282	zztmpy = ( ts(ji,jj+1,1,jp_tem,Kmm) - zztmp ) * r1_e2v(ji,jj) + ( zztmp - ts(ji ,jj-1,1,jp_tem,Kmm) ) * r1_e2v(ji,jj-1)
	283	z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy ) &
	284	& * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1)
	285	END_2D
[9019]	286	CALL iom_put( "sstgrad2", z2d ) ! square of module of sst gradient
[13942]	287	DO_2D( 0, 0, 0, 0 )
	288	z2d(ji,jj) = SQRT( z2d(ji,jj) )
	289	END_2D
[9019]	290	CALL iom_put( "sstgrad" , z2d ) ! module of sst gradient
[4990]	291	ENDIF
	292
[9019]	293	! heat and salt contents
[4990]	294	IF( iom_use("heatc") ) THEN
[7753]	295	z2d(:,:) = 0._wp
[13942]	296	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	297	z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_tem,Kmm) * tmask(ji,jj,jk)
	298	END_3D
[12489]	299	CALL iom_put( "heatc", rho0_rcp * z2d ) ! vertically integrated heat content (J/m2)
[4990]	300	ENDIF
	301
	302	IF( iom_use("saltc") ) THEN
[7753]	303	z2d(:,:) = 0._wp
[13942]	304	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	305	z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) * tmask(ji,jj,jk)
	306	END_3D
[12489]	307	CALL iom_put( "saltc", rho0 * z2d ) ! vertically integrated salt content (PSU*kg/m2)
[4990]	308	ENDIF
[4840]	309	!
[13942]	310	IF( iom_use("salt2c") ) THEN
	311	z2d(:,:) = 0._wp
	312	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
	313	z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) * tmask(ji,jj,jk)
	314	END_3D
	315	CALL iom_put( "salt2c", rho0 * z2d ) ! vertically integrated salt content (PSU*kg/m2)
	316	ENDIF
	317	!
	318	IF ( iom_use("ke") .OR. iom_use("ke_int") ) THEN
[9399]	319	z3d(:,:,jpk) = 0._wp
[13295]	320	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[13942]	321	zztmpx = 0.5 * ( uu(ji-1,jj ,jk,Kmm) + uu(ji,jj,jk,Kmm) )
	322	zztmpy = 0.5 * ( vv(ji ,jj-1,jk,Kmm) + vv(ji,jj,jk,Kmm) )
	323	z3d(ji,jj,jk) = 0.5 * ( zztmpxzztmpx + zztmpyzztmpy )
[12377]	324	END_3D
[13942]	325	CALL iom_put( "ke", z3d ) ! kinetic energy
	326
	327	z2d(:,:) = 0._wp
	328	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
	329	z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * z3d(ji,jj,jk) * e1e2t(ji,jj) * tmask(ji,jj,jk)
	330	END_3D
	331	CALL iom_put( "ke_int", z2d ) ! vertically integrated kinetic energy
[4990]	332	ENDIF
[6351]	333	!
[12377]	334	CALL iom_put( "hdiv", hdiv ) ! Horizontal divergence
[13942]	335
	336	IF ( iom_use("relvor") .OR. iom_use("absvor") .OR. iom_use("potvor") ) THEN
	337
	338	z3d(:,:,jpk) = 0._wp
	339	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
	340	z3d(ji,jj,jk) = ( e2v(ji+1,jj ) * vv(ji+1,jj ,jk,Kmm) - e2v(ji,jj) * vv(ji,jj,jk,Kmm) &
	341	& - e1u(ji ,jj+1) * uu(ji ,jj+1,jk,Kmm) + e1u(ji,jj) * uu(ji,jj,jk,Kmm) ) * r1_e1e2f(ji,jj)
	342	END_3D
	343	CALL iom_put( "relvor", z3d ) ! relative vorticity
	344
	345	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
	346	z3d(ji,jj,jk) = ff_f(ji,jj) + z3d(ji,jj,jk)
	347	END_3D
	348	CALL iom_put( "absvor", z3d ) ! absolute vorticity
	349
	350	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
	351	ze3 = ( e3t(ji,jj+1,jk,Kmm)tmask(ji,jj+1,jk) + e3t(ji+1,jj+1,jk,Kmm)tmask(ji+1,jj+1,jk) &
	352	& + e3t(ji,jj ,jk,Kmm)tmask(ji,jj ,jk) + e3t(ji+1,jj ,jk,Kmm)tmask(ji+1,jj ,jk) )
	353	IF( ze3 /= 0._wp ) THEN ; ze3 = 4._wp / ze3
	354	ELSE ; ze3 = 0._wp
	355	ENDIF
	356	z3d(ji,jj,jk) = ze3 * z3d(ji,jj,jk)
	357	END_3D
	358	CALL iom_put( "potvor", z3d ) ! potential vorticity
	359
	360	ENDIF
[6351]	361	!
[7646]	362	IF( iom_use("u_masstr") .OR. iom_use("u_masstr_vint") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN
[7753]	363	z3d(:,:,jpk) = 0.e0
	364	z2d(:,:) = 0.e0
[1756]	365	DO jk = 1, jpkm1
[12489]	366	z3d(:,:,jk) = rho0 * uu(:,:,jk,Kmm) * e2u(:,:) * e3u(:,:,jk,Kmm) * umask(:,:,jk)
[7753]	367	z2d(:,:) = z2d(:,:) + z3d(:,:,jk)
[1756]	368	END DO
[9019]	369	CALL iom_put( "u_masstr" , z3d ) ! mass transport in i-direction
	370	CALL iom_put( "u_masstr_vint", z2d ) ! mass transport in i-direction vertical sum
[4990]	371	ENDIF
	372
	373	IF( iom_use("u_heattr") ) THEN
[9019]	374	z2d(:,:) = 0._wp
[13295]	375	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	376	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_tem,Kmm) + ts(ji+1,jj,jk,jp_tem,Kmm) )
	377	END_3D
[9019]	378	CALL iom_put( "u_heattr", 0.5rcp z2d ) ! heat transport in i-direction
[4990]	379	ENDIF
[4761]	380
[4990]	381	IF( iom_use("u_salttr") ) THEN
[7753]	382	z2d(:,:) = 0.e0
[13295]	383	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	384	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_sal,Kmm) + ts(ji+1,jj,jk,jp_sal,Kmm) )
	385	END_3D
[9019]	386	CALL iom_put( "u_salttr", 0.5 * z2d ) ! heat transport in i-direction
[4990]	387	ENDIF
[4761]	388
[4990]	389
	390	IF( iom_use("v_masstr") .OR. iom_use("v_heattr") .OR. iom_use("v_salttr") ) THEN
[7753]	391	z3d(:,:,jpk) = 0.e0
[1756]	392	DO jk = 1, jpkm1
[12489]	393	z3d(:,:,jk) = rho0 * vv(:,:,jk,Kmm) * e1v(:,:) * e3v(:,:,jk,Kmm) * vmask(:,:,jk)
[1756]	394	END DO
[9019]	395	CALL iom_put( "v_masstr", z3d ) ! mass transport in j-direction
[4990]	396	ENDIF
	397
	398	IF( iom_use("v_heattr") ) THEN
[7753]	399	z2d(:,:) = 0.e0
[13295]	400	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	401	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_tem,Kmm) + ts(ji,jj+1,jk,jp_tem,Kmm) )
	402	END_3D
[9019]	403	CALL iom_put( "v_heattr", 0.5rcp z2d ) ! heat transport in j-direction
[4990]	404	ENDIF
[4761]	405
[4990]	406	IF( iom_use("v_salttr") ) THEN
[9019]	407	z2d(:,:) = 0._wp
[13295]	408	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	409	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_sal,Kmm) + ts(ji,jj+1,jk,jp_sal,Kmm) )
	410	END_3D
[9019]	411	CALL iom_put( "v_salttr", 0.5 * z2d ) ! heat transport in j-direction
[1756]	412	ENDIF
[7646]	413
	414	IF( iom_use("tosmint") ) THEN
[9019]	415	z2d(:,:) = 0._wp
[13295]	416	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	417	z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_tem,Kmm)
	418	END_3D
[12489]	419	CALL iom_put( "tosmint", rho0 * z2d ) ! Vertical integral of temperature
[7646]	420	ENDIF
	421	IF( iom_use("somint") ) THEN
[7753]	422	z2d(:,:)=0._wp
[13295]	423	DO_3D( 0, 0, 0, 0, 1, jpkm1 )
[12377]	424	z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm)
	425	END_3D
[12489]	426	CALL iom_put( "somint", rho0 * z2d ) ! Vertical integral of salinity
[7646]	427	ENDIF
	428
[9019]	429	CALL iom_put( "bn2", rn2 ) ! Brunt-Vaisala buoyancy frequency (N^2)
[2528]	430	!
[12377]	431
	432	IF (ln_dia25h) CALL dia_25h( kt, Kmm ) ! 25h averaging
[6140]	433
[9124]	434	IF( ln_timing ) CALL timing_stop('dia_wri')
[3294]	435	!
[1482]	436	END SUBROUTINE dia_wri
	437
	438	#else
[2528]	439	!!----------------------------------------------------------------------
	440	!! Default option use IOIPSL library
	441	!!----------------------------------------------------------------------
	442
[9652]	443	INTEGER FUNCTION dia_wri_alloc()
	444	!!----------------------------------------------------------------------
	445	INTEGER, DIMENSION(2) :: ierr
	446	!!----------------------------------------------------------------------
[12493]	447	IF( nn_write == -1 ) THEN
	448	dia_wri_alloc = 0
	449	ELSE
	450	ierr = 0
	451	ALLOCATE( ndex_hT(jpijpj) , ndex_T(jpijpj*jpk) , &
	452	& ndex_hU(jpijpj) , ndex_U(jpijpj*jpk) , &
	453	& ndex_hV(jpijpj) , ndex_V(jpijpj*jpk) , STAT=ierr(1) )
[9652]	454	!
[12493]	455	dia_wri_alloc = MAXVAL(ierr)
	456	CALL mpp_sum( 'diawri', dia_wri_alloc )
	457	!
	458	ENDIF
[9652]	459	!
	460	END FUNCTION dia_wri_alloc
[12377]	461
	462	INTEGER FUNCTION dia_wri_alloc_abl()
	463	!!----------------------------------------------------------------------
	464	ALLOCATE( ndex_hA(jpijpj), ndex_A (jpijpj*jpkam1), STAT=dia_wri_alloc_abl)
	465	CALL mpp_sum( 'diawri', dia_wri_alloc_abl )
	466	!
	467	END FUNCTION dia_wri_alloc_abl
[9652]	468
	469
[12377]	470	SUBROUTINE dia_wri( kt, Kmm )
[3]	471	!!---------------------------------------------------------------------
	472	!! * ROUTINE dia_wri *
	473	!!
	474	!! ** Purpose : Standard output of opa: dynamics and tracer fields
	475	!! NETCDF format is used by default
	476	!!
	477	!! ** Method : At the beginning of the first time step (nit000),
	478	!! define all the NETCDF files and fields
	479	!! At each time step call histdef to compute the mean if ncessary
[11536]	480	!! Each nn_write time step, output the instantaneous or mean fields
[3]	481	!!----------------------------------------------------------------------
[5836]	482	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[12377]	483	INTEGER, INTENT( in ) :: Kmm ! ocean time level index
[5836]	484	!
[2528]	485	LOGICAL :: ll_print = .FALSE. ! =T print and flush numout
	486	CHARACTER (len=40) :: clhstnam, clop, clmx ! local names
	487	INTEGER :: inum = 11 ! temporary logical unit
[3294]	488	INTEGER :: ji, jj, jk ! dummy loop indices
	489	INTEGER :: ierr ! error code return from allocation
[2528]	490	INTEGER :: iimi, iima, ipk, it, itmod, ijmi, ijma ! local integers
[12377]	491	INTEGER :: ipka ! ABL
[3609]	492	INTEGER :: jn, ierror ! local integers
[6140]	493	REAL(wp) :: zsto, zout, zmax, zjulian ! local scalars
[5836]	494	!
[9019]	495	REAL(wp), DIMENSION(jpi,jpj) :: zw2d ! 2D workspace
[13237]	496	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zw3d, ze3t, zgdept ! 3D workspace
[12377]	497	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d_abl ! ABL 3D workspace
[3]	498	!!----------------------------------------------------------------------
[1482]	499	!
[9019]	500	IF( ninist == 1 ) THEN !== Output the initial state and forcings ==!
[12377]	501	CALL dia_wri_state( Kmm, 'output.init' )
[1482]	502	ninist = 0
	503	ENDIF
	504	!
[11536]	505	IF( nn_write == -1 ) RETURN ! we will never do any output
	506	!
	507	IF( ln_timing ) CALL timing_start('dia_wri')
	508	!
[3]	509	! 0. Initialisation
	510	! -----------------
[632]	511
[9019]	512	ll_print = .FALSE. ! local variable for debugging
[3]	513	ll_print = ll_print .AND. lwp
	514
	515	! Define frequency of output and means
[5566]	516	clop = "x" ! no use of the mask value (require less cpu time and otherwise the model crashes)
[3]	517	#if defined key_diainstant
[12489]	518	zsto = nn_write * rn_Dt
[1312]	519	clop = "inst("//TRIM(clop)//")"
[3]	520	#else
[12489]	521	zsto=rn_Dt
[1312]	522	clop = "ave("//TRIM(clop)//")"
[3]	523	#endif
[12489]	524	zout = nn_write * rn_Dt
	525	zmax = ( nitend - nit000 + 1 ) * rn_Dt
[3]	526
	527	! Define indices of the horizontal output zoom and vertical limit storage
[13286]	528	iimi = Nis0 ; iima = Nie0
	529	ijmi = Njs0 ; ijma = Nje0
[3]	530	ipk = jpk
[12377]	531	IF(ln_abl) ipka = jpkam1
[3]	532
	533	! define time axis
[1334]	534	it = kt
	535	itmod = kt - nit000 + 1
[3]	536
[13237]	537	! store e3t for subsitute
	538	DO jk = 1, jpk
	539	ze3t (:,:,jk) = e3t (:,:,jk,Kmm)
	540	zgdept(:,:,jk) = gdept(:,:,jk,Kmm)
	541	END DO
[3]	542
[13237]	543
[3]	544	! 1. Define NETCDF files and fields at beginning of first time step
	545	! -----------------------------------------------------------------
	546
	547	IF( kt == nit000 ) THEN
	548
	549	! Define the NETCDF files (one per grid)
[632]	550
[3]	551	! Compute julian date from starting date of the run
[12489]	552	CALL ymds2ju( nyear, nmonth, nday, rn_Dt, zjulian )
[1309]	553	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
[3]	554	IF(lwp)WRITE(numout,*)
	555	IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear, &
	556	& ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
	557	IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma, &
	558	' limit storage in depth = ', ipk
	559
	560	! WRITE root name in date.file for use by postpro
[1581]	561	IF(lwp) THEN
[11536]	562	CALL dia_nam( clhstnam, nn_write,' ' )
[1581]	563	CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
[895]	564	WRITE(inum,*) clhstnam
	565	CLOSE(inum)
	566	ENDIF
[632]	567
[3]	568	! Define the T grid FILE ( nid_T )
[632]	569
[11536]	570	CALL dia_nam( clhstnam, nn_write, 'grid_T' )
[3]	571	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
	572	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
	573	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[12489]	574	& nit000-1, zjulian, rn_Dt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )
[3]	575	CALL histvert( nid_T, "deptht", "Vertical T levels", & ! Vertical grid: gdept
[4292]	576	& "m", ipk, gdept_1d, nz_T, "down" )
[3]	577	! ! Index of ocean points
	578	CALL wheneq( jpijpjipk, tmask, 1, 1., ndex_T , ndim_T ) ! volume
	579	CALL wheneq( jpi*jpj , tmask, 1, 1., ndex_hT, ndim_hT ) ! surface
[3609]	580	!
	581	IF( ln_icebergs ) THEN
	582	!
	583	!! allocation cant go in dia_wri_alloc because ln_icebergs is only set after
	584	!! that routine is called from nemogcm, so do it here immediately before its needed
	585	ALLOCATE( ndex_bT(jpijpjnclasses), STAT=ierror )
[10425]	586	CALL mpp_sum( 'diawri', ierror )
[3609]	587	IF( ierror /= 0 ) THEN
	588	CALL ctl_stop('dia_wri: failed to allocate iceberg diagnostic array')
	589	RETURN
	590	ENDIF
	591	!
	592	!! iceberg vertical coordinate is class number
	593	CALL histvert( nid_T, "class", "Iceberg class", & ! Vertical grid: class
	594	& "number", nclasses, class_num, nb_T )
	595	!
	596	!! each class just needs the surface index pattern
	597	ndim_bT = 3
	598	DO jn = 1,nclasses
	599	ndex_bT((jn-1)jpijpj+1:jnjpijpj) = ndex_hT(1:jpi*jpj)
	600	ENDDO
	601	!
	602	ENDIF
[3]	603
	604	! Define the U grid FILE ( nid_U )
	605
[11536]	606	CALL dia_nam( clhstnam, nn_write, 'grid_U' )
[3]	607	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
	608	CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu, & ! Horizontal grid: glamu and gphiu
	609	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[12489]	610	& nit000-1, zjulian, rn_Dt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )
[3]	611	CALL histvert( nid_U, "depthu", "Vertical U levels", & ! Vertical grid: gdept
[4292]	612	& "m", ipk, gdept_1d, nz_U, "down" )
[3]	613	! ! Index of ocean points
	614	CALL wheneq( jpijpjipk, umask, 1, 1., ndex_U , ndim_U ) ! volume
	615	CALL wheneq( jpi*jpj , umask, 1, 1., ndex_hU, ndim_hU ) ! surface
	616
	617	! Define the V grid FILE ( nid_V )
	618
[11536]	619	CALL dia_nam( clhstnam, nn_write, 'grid_V' ) ! filename
[3]	620	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
	621	CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv, & ! Horizontal grid: glamv and gphiv
	622	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[12489]	623	& nit000-1, zjulian, rn_Dt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )
[3]	624	CALL histvert( nid_V, "depthv", "Vertical V levels", & ! Vertical grid : gdept
[4292]	625	& "m", ipk, gdept_1d, nz_V, "down" )
[3]	626	! ! Index of ocean points
	627	CALL wheneq( jpijpjipk, vmask, 1, 1., ndex_V , ndim_V ) ! volume
	628	CALL wheneq( jpi*jpj , vmask, 1, 1., ndex_hV, ndim_hV ) ! surface
	629
	630	! Define the W grid FILE ( nid_W )
	631
[11536]	632	CALL dia_nam( clhstnam, nn_write, 'grid_W' ) ! filename
[3]	633	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
	634	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
	635	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[12489]	636	& nit000-1, zjulian, rn_Dt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set )
[3]	637	CALL histvert( nid_W, "depthw", "Vertical W levels", & ! Vertical grid: gdepw
[4292]	638	& "m", ipk, gdepw_1d, nz_W, "down" )
[3]	639
[12377]	640	IF( ln_abl ) THEN
	641	! Define the ABL grid FILE ( nid_A )
	642	CALL dia_nam( clhstnam, nn_write, 'grid_ABL' )
	643	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
	644	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
	645	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[12489]	646	& nit000-1, zjulian, rn_Dt, nh_A, nid_A, domain_id=nidom, snc4chunks=snc4set )
[12377]	647	CALL histvert( nid_A, "ght_abl", "Vertical T levels", & ! Vertical grid: gdept
	648	& "m", ipka, ght_abl(2:jpka), nz_A, "up" )
	649	! ! Index of ocean points
	650	ALLOCATE( zw3d_abl(jpi,jpj,ipka) )
	651	zw3d_abl(:,:,:) = 1._wp
	652	CALL wheneq( jpijpjipka, zw3d_abl, 1, 1., ndex_A , ndim_A ) ! volume
	653	CALL wheneq( jpi*jpj , zw3d_abl, 1, 1., ndex_hA, ndim_hA ) ! surface
	654	DEALLOCATE(zw3d_abl)
	655	ENDIF
[13237]	656	!
[632]	657
[3]	658	! Declare all the output fields as NETCDF variables
	659
	660	! !!! nid_T : 3D
	661	CALL histdef( nid_T, "votemper", "Temperature" , "C" , & ! tn
	662	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	663	CALL histdef( nid_T, "vosaline", "Salinity" , "PSU" , & ! sn
	664	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
[6140]	665	IF( .NOT.ln_linssh ) THEN
[13237]	666	CALL histdef( nid_T, "vovvle3t", "Level thickness" , "m" ,& ! e3t n
[4292]	667	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
[13237]	668	CALL histdef( nid_T, "vovvldep", "T point depth" , "m" ,& ! e3t n
[4292]	669	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
[13237]	670	CALL histdef( nid_T, "vovvldef", "Squared level deformation" , "%^2" ,& ! e3t n
[4292]	671	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	672	ENDIF
[3]	673	! !!! nid_T : 2D
	674	CALL histdef( nid_T, "sosstsst", "Sea Surface temperature" , "C" , & ! sst
	675	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	676	CALL histdef( nid_T, "sosaline", "Sea Surface Salinity" , "PSU" , & ! sss
	677	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[359]	678	CALL histdef( nid_T, "sossheig", "Sea Surface Height" , "m" , & ! ssh
[3]	679	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[2528]	680	CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux" , "Kg/m2/s", & ! (emp-rnf)
[3]	681	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[4570]	682	CALL histdef( nid_T, "sorunoff", "River runoffs" , "Kg/m2/s", & ! runoffs
	683	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[3625]	684	CALL histdef( nid_T, "sosfldow", "downward salt flux" , "PSU/m2/s", & ! sfx
[3]	685	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[6140]	686	IF( ln_linssh ) THEN
[12377]	687	CALL histdef( nid_T, "sosst_cd", "Concentration/Dilution term on temperature" & ! emp * ts(:,:,1,jp_tem,Kmm)
[3625]	688	& , "KgC/m2/s", & ! sosst_cd
[4292]	689	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[12377]	690	CALL histdef( nid_T, "sosss_cd", "Concentration/Dilution term on salinity" & ! emp * ts(:,:,1,jp_sal,Kmm)
[3625]	691	& , "KgPSU/m2/s",& ! sosss_cd
[4292]	692	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	693	ENDIF
[888]	694	CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! qns + qsr
[3]	695	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	696	CALL histdef( nid_T, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! qsr
	697	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1585]	698	CALL histdef( nid_T, "somixhgt", "Turbocline Depth" , "m" , & ! hmld
	699	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[3]	700	CALL histdef( nid_T, "somxl010", "Mixed Layer Depth 0.01" , "m" , & ! hmlp
	701	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1037]	702	CALL histdef( nid_T, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
[3]	703	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1649]	704	CALL histdef( nid_T, "sowindsp", "wind speed at 10m" , "m/s" , & ! wndm
	705	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[12377]	706	!
	707	IF( ln_abl ) THEN
	708	CALL histdef( nid_A, "t_abl", "Potential Temperature" , "K" , & ! t_abl
	709	& jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
	710	CALL histdef( nid_A, "q_abl", "Humidity" , "kg/kg" , & ! q_abl
	711	& jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
	712	CALL histdef( nid_A, "u_abl", "Atmospheric U-wind " , "m/s" , & ! u_abl
	713	& jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
	714	CALL histdef( nid_A, "v_abl", "Atmospheric V-wind " , "m/s" , & ! v_abl
	715	& jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
	716	CALL histdef( nid_A, "tke_abl", "Atmospheric TKE " , "m2/s2" , & ! tke_abl
	717	& jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
	718	CALL histdef( nid_A, "avm_abl", "Atmospheric turbulent viscosity", "m2/s" , & ! avm_abl
	719	& jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
	720	CALL histdef( nid_A, "avt_abl", "Atmospheric turbulent diffusivity", "m2/s2", & ! avt_abl
	721	& jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
	722	CALL histdef( nid_A, "pblh", "Atmospheric boundary layer height " , "m", & ! pblh
	723	& jpi, jpj, nh_A, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	724	#if defined key_si3
	725	CALL histdef( nid_A, "oce_frac", "Fraction of open ocean" , " ", & ! ato_i
	726	& jpi, jpj, nh_A, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	727	#endif
	728	CALL histend( nid_A, snc4chunks=snc4set )
	729	ENDIF
	730	!
[3609]	731	IF( ln_icebergs ) THEN
	732	CALL histdef( nid_T, "calving" , "calving mass input" , "kg/s" , &
	733	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	734	CALL histdef( nid_T, "calving_heat" , "calving heat flux" , "XXXX" , &
	735	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	736	CALL histdef( nid_T, "berg_floating_melt" , "Melt rate of icebergs + bits" , "kg/m2/s", &
	737	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	738	CALL histdef( nid_T, "berg_stored_ice" , "Accumulated ice mass by class" , "kg" , &
	739	& jpi, jpj, nh_T, nclasses , 1, nclasses , nb_T , 32, clop, zsto, zout )
	740	IF( ln_bergdia ) THEN
	741	CALL histdef( nid_T, "berg_melt" , "Melt rate of icebergs" , "kg/m2/s", &
	742	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	743	CALL histdef( nid_T, "berg_buoy_melt" , "Buoyancy component of iceberg melt rate" , "kg/m2/s", &
	744	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	745	CALL histdef( nid_T, "berg_eros_melt" , "Erosion component of iceberg melt rate" , "kg/m2/s", &
	746	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	747	CALL histdef( nid_T, "berg_conv_melt" , "Convective component of iceberg melt rate", "kg/m2/s", &
	748	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	749	CALL histdef( nid_T, "berg_virtual_area" , "Virtual coverage by icebergs" , "m2" , &
	750	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	751	CALL histdef( nid_T, "bits_src" , "Mass source of bergy bits" , "kg/m2/s", &
	752	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	753	CALL histdef( nid_T, "bits_melt" , "Melt rate of bergy bits" , "kg/m2/s", &
	754	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	755	CALL histdef( nid_T, "bits_mass" , "Bergy bit density field" , "kg/m2" , &
	756	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	757	CALL histdef( nid_T, "berg_mass" , "Iceberg density field" , "kg/m2" , &
	758	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	759	CALL histdef( nid_T, "berg_real_calving" , "Calving into iceberg class" , "kg/s" , &
	760	& jpi, jpj, nh_T, nclasses , 1, nclasses , nb_T , 32, clop, zsto, zout )
	761	ENDIF
	762	ENDIF
	763
[11536]	764	IF( ln_ssr ) THEN
[4990]	765	CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping" , "W/m2" , & ! qrp
	766	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	767	CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping" , "Kg/m2/s", & ! erp
	768	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	769	CALL histdef( nid_T, "sosafldp", "Surface salt flux: damping" , "Kg/m2/s", & ! erp * sn
	770	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	771	ENDIF
[11536]	772
[3]	773	clmx ="l_max(only(x))" ! max index on a period
[5836]	774	! CALL histdef( nid_T, "sobowlin", "Bowl Index" , "W-point", & ! bowl INDEX
	775	! & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clmx, zsto, zout )
[3]	776	#if defined key_diahth
	777	CALL histdef( nid_T, "sothedep", "Thermocline Depth" , "m" , & ! hth
	778	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	779	CALL histdef( nid_T, "so20chgt", "Depth of 20C isotherm" , "m" , & ! hd20
	780	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	781	CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm" , "m" , & ! hd28
	782	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[7646]	783	CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , "J/m2" , & ! htc3
[3]	784	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	785	#endif
	786
[2528]	787	CALL histend( nid_T, snc4chunks=snc4set )
[3]	788
	789	! !!! nid_U : 3D
[12377]	790	CALL histdef( nid_U, "vozocrtx", "Zonal Current" , "m/s" , & ! uu(:,:,:,Kmm)
[3]	791	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
[7646]	792	IF( ln_wave .AND. ln_sdw) THEN
	793	CALL histdef( nid_U, "sdzocrtx", "Stokes Drift Zonal Current" , "m/s" , & ! usd
	794	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
	795	ENDIF
[3]	796	! !!! nid_U : 2D
[888]	797	CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau
[3]	798	& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
	799
[2528]	800	CALL histend( nid_U, snc4chunks=snc4set )
[3]	801
	802	! !!! nid_V : 3D
[12377]	803	CALL histdef( nid_V, "vomecrty", "Meridional Current" , "m/s" , & ! vv(:,:,:,Kmm)
[3]	804	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
[7646]	805	IF( ln_wave .AND. ln_sdw) THEN
	806	CALL histdef( nid_V, "sdmecrty", "Stokes Drift Meridional Current" , "m/s" , & ! vsd
	807	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
	808	ENDIF
[3]	809	! !!! nid_V : 2D
[888]	810	CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau
[3]	811	& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
	812
[2528]	813	CALL histend( nid_V, snc4chunks=snc4set )
[3]	814
	815	! !!! nid_W : 3D
[12377]	816	CALL histdef( nid_W, "vovecrtz", "Vertical Velocity" , "m/s" , & ! ww
[3]	817	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	818	CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity" , "m2/s" , & ! avt
	819	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[9019]	820	CALL histdef( nid_W, "votkeavm", "Vertical Eddy Viscosity" , "m2/s" , & ! avm
[255]	821	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	822
[9019]	823	IF( ln_zdfddm ) THEN
[3]	824	CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity" , "m2/s" , & ! avs
	825	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	826	ENDIF
[7646]	827
	828	IF( ln_wave .AND. ln_sdw) THEN
	829	CALL histdef( nid_W, "sdvecrtz", "Stokes Drift Vertical Current" , "m/s" , & ! wsd
	830	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	831	ENDIF
[3]	832	! !!! nid_W : 2D
[2528]	833	CALL histend( nid_W, snc4chunks=snc4set )
[3]	834
	835	IF(lwp) WRITE(numout,*)
	836	IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
	837	IF(ll_print) CALL FLUSH(numout )
	838
	839	ENDIF
	840
	841	! 2. Start writing data
	842	! ---------------------
	843
[4292]	844	! ndex(1) est utilise ssi l'avant dernier argument est different de
[3]	845	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
	846	! donne le nombre d'elements, et ndex la liste des indices a sortir
	847
[11536]	848	IF( lwp .AND. MOD( itmod, nn_write ) == 0 ) THEN
[3]	849	WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
	850	WRITE(numout,*) '~~~~~~ '
	851	ENDIF
	852
[6140]	853	IF( .NOT.ln_linssh ) THEN
[13237]	854	CALL histwrite( nid_T, "votemper", it, ts(:,:,:,jp_tem,Kmm) * ze3t(:,:,:) , ndim_T , ndex_T ) ! heat content
	855	CALL histwrite( nid_T, "vosaline", it, ts(:,:,:,jp_sal,Kmm) * ze3t(:,:,:) , ndim_T , ndex_T ) ! salt content
	856	CALL histwrite( nid_T, "sosstsst", it, ts(:,:,1,jp_tem,Kmm) * ze3t(:,:,1) , ndim_hT, ndex_hT ) ! sea surface heat content
	857	CALL histwrite( nid_T, "sosaline", it, ts(:,:,1,jp_sal,Kmm) * ze3t(:,:,1) , ndim_hT, ndex_hT ) ! sea surface salinity content
[4292]	858	ELSE
[12377]	859	CALL histwrite( nid_T, "votemper", it, ts(:,:,:,jp_tem,Kmm) , ndim_T , ndex_T ) ! temperature
	860	CALL histwrite( nid_T, "vosaline", it, ts(:,:,:,jp_sal,Kmm) , ndim_T , ndex_T ) ! salinity
	861	CALL histwrite( nid_T, "sosstsst", it, ts(:,:,1,jp_tem,Kmm) , ndim_hT, ndex_hT ) ! sea surface temperature
	862	CALL histwrite( nid_T, "sosaline", it, ts(:,:,1,jp_sal,Kmm) , ndim_hT, ndex_hT ) ! sea surface salinity
[4292]	863	ENDIF
[6140]	864	IF( .NOT.ln_linssh ) THEN
[13237]	865	zw3d(:,:,:) = ( ( ze3t(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2
	866	CALL histwrite( nid_T, "vovvle3t", it, ze3t (:,:,:) , ndim_T , ndex_T ) ! level thickness
	867	CALL histwrite( nid_T, "vovvldep", it, zgdept , ndim_T , ndex_T ) ! t-point depth
[4292]	868	CALL histwrite( nid_T, "vovvldef", it, zw3d , ndim_T , ndex_T ) ! level thickness deformation
	869	ENDIF
[12377]	870	CALL histwrite( nid_T, "sossheig", it, ssh(:,:,Kmm) , ndim_hT, ndex_hT ) ! sea surface height
[2528]	871	CALL histwrite( nid_T, "sowaflup", it, ( emp-rnf ) , ndim_hT, ndex_hT ) ! upward water flux
[4570]	872	CALL histwrite( nid_T, "sorunoff", it, rnf , ndim_hT, ndex_hT ) ! river runoffs
[3625]	873	CALL histwrite( nid_T, "sosfldow", it, sfx , ndim_hT, ndex_hT ) ! downward salt flux
	874	! (includes virtual salt flux beneath ice
	875	! in linear free surface case)
[6140]	876	IF( ln_linssh ) THEN
[12377]	877	zw2d(:,:) = emp (:,:) * ts(:,:,1,jp_tem,Kmm)
[4292]	878	CALL histwrite( nid_T, "sosst_cd", it, zw2d, ndim_hT, ndex_hT ) ! c/d term on sst
[12377]	879	zw2d(:,:) = emp (:,:) * ts(:,:,1,jp_sal,Kmm)
[4292]	880	CALL histwrite( nid_T, "sosss_cd", it, zw2d, ndim_hT, ndex_hT ) ! c/d term on sss
	881	ENDIF
[888]	882	CALL histwrite( nid_T, "sohefldo", it, qns + qsr , ndim_hT, ndex_hT ) ! total heat flux
[3]	883	CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux
[1585]	884	CALL histwrite( nid_T, "somixhgt", it, hmld , ndim_hT, ndex_hT ) ! turbocline depth
[3]	885	CALL histwrite( nid_T, "somxl010", it, hmlp , ndim_hT, ndex_hT ) ! mixed layer depth
[1037]	886	CALL histwrite( nid_T, "soicecov", it, fr_i , ndim_hT, ndex_hT ) ! ice fraction
[1649]	887	CALL histwrite( nid_T, "sowindsp", it, wndm , ndim_hT, ndex_hT ) ! wind speed
[12377]	888	!
	889	IF( ln_abl ) THEN
	890	ALLOCATE( zw3d_abl(jpi,jpj,jpka) )
	891	IF( ln_mskland ) THEN
	892	DO jk=1,jpka
	893	zw3d_abl(:,:,jk) = tmask(:,:,1)
	894	END DO
	895	ELSE
	896	zw3d_abl(:,:,:) = 1._wp
	897	ENDIF
	898	CALL histwrite( nid_A, "pblh" , it, pblh(:,:) *zw3d_abl(:,:,1 ), ndim_hA, ndex_hA ) ! pblh
	899	CALL histwrite( nid_A, "u_abl" , it, u_abl (:,:,2:jpka,nt_n )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! u_abl
	900	CALL histwrite( nid_A, "v_abl" , it, v_abl (:,:,2:jpka,nt_n )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! v_abl
	901	CALL histwrite( nid_A, "t_abl" , it, tq_abl (:,:,2:jpka,nt_n,1)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! t_abl
	902	CALL histwrite( nid_A, "q_abl" , it, tq_abl (:,:,2:jpka,nt_n,2)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! q_abl
	903	CALL histwrite( nid_A, "tke_abl", it, tke_abl (:,:,2:jpka,nt_n )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! tke_abl
	904	CALL histwrite( nid_A, "avm_abl", it, avm_abl (:,:,2:jpka )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! avm_abl
	905	CALL histwrite( nid_A, "avt_abl", it, avt_abl (:,:,2:jpka )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! avt_abl
	906	#if defined key_si3
	907	CALL histwrite( nid_A, "oce_frac" , it, ato_i(:,:) , ndim_hA, ndex_hA ) ! ato_i
	908	#endif
	909	DEALLOCATE(zw3d_abl)
	910	ENDIF
	911	!
[3609]	912	IF( ln_icebergs ) THEN
	913	!
	914	CALL histwrite( nid_T, "calving" , it, berg_grid%calving , ndim_hT, ndex_hT )
	915	CALL histwrite( nid_T, "calving_heat" , it, berg_grid%calving_hflx , ndim_hT, ndex_hT )
	916	CALL histwrite( nid_T, "berg_floating_melt" , it, berg_grid%floating_melt, ndim_hT, ndex_hT )
	917	!
	918	CALL histwrite( nid_T, "berg_stored_ice" , it, berg_grid%stored_ice , ndim_bT, ndex_bT )
	919	!
	920	IF( ln_bergdia ) THEN
	921	CALL histwrite( nid_T, "berg_melt" , it, berg_melt , ndim_hT, ndex_hT )
	922	CALL histwrite( nid_T, "berg_buoy_melt" , it, buoy_melt , ndim_hT, ndex_hT )
	923	CALL histwrite( nid_T, "berg_eros_melt" , it, eros_melt , ndim_hT, ndex_hT )
	924	CALL histwrite( nid_T, "berg_conv_melt" , it, conv_melt , ndim_hT, ndex_hT )
	925	CALL histwrite( nid_T, "berg_virtual_area" , it, virtual_area , ndim_hT, ndex_hT )
	926	CALL histwrite( nid_T, "bits_src" , it, bits_src , ndim_hT, ndex_hT )
	927	CALL histwrite( nid_T, "bits_melt" , it, bits_melt , ndim_hT, ndex_hT )
	928	CALL histwrite( nid_T, "bits_mass" , it, bits_mass , ndim_hT, ndex_hT )
	929	CALL histwrite( nid_T, "berg_mass" , it, berg_mass , ndim_hT, ndex_hT )
	930	!
	931	CALL histwrite( nid_T, "berg_real_calving" , it, real_calving , ndim_bT, ndex_bT )
	932	ENDIF
	933	ENDIF
	934
[11536]	935	IF( ln_ssr ) THEN
[4990]	936	CALL histwrite( nid_T, "sohefldp", it, qrp , ndim_hT, ndex_hT ) ! heat flux damping
	937	CALL histwrite( nid_T, "sowafldp", it, erp , ndim_hT, ndex_hT ) ! freshwater flux damping
[12377]	938	zw2d(:,:) = erp(:,:) * ts(:,:,1,jp_sal,Kmm) * tmask(:,:,1)
[4990]	939	CALL histwrite( nid_T, "sosafldp", it, zw2d , ndim_hT, ndex_hT ) ! salt flux damping
	940	ENDIF
	941	! zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)
	942	! CALL histwrite( nid_T, "sobowlin", it, zw2d , ndim_hT, ndex_hT ) ! ???
[3]	943
	944	#if defined key_diahth
	945	CALL histwrite( nid_T, "sothedep", it, hth , ndim_hT, ndex_hT ) ! depth of the thermocline
	946	CALL histwrite( nid_T, "so20chgt", it, hd20 , ndim_hT, ndex_hT ) ! depth of the 20 isotherm
	947	CALL histwrite( nid_T, "so28chgt", it, hd28 , ndim_hT, ndex_hT ) ! depth of the 28 isotherm
	948	CALL histwrite( nid_T, "sohtc300", it, htc3 , ndim_hT, ndex_hT ) ! first 300m heaat content
	949	#endif
[888]	950
[12377]	951	CALL histwrite( nid_U, "vozocrtx", it, uu(:,:,:,Kmm) , ndim_U , ndex_U ) ! i-current
[888]	952	CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress
[3]	953
[12377]	954	CALL histwrite( nid_V, "vomecrty", it, vv(:,:,:,Kmm) , ndim_V , ndex_V ) ! j-current
[888]	955	CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress
[3]	956
[11418]	957	IF( ln_zad_Aimp ) THEN
[12377]	958	CALL histwrite( nid_W, "vovecrtz", it, ww + wi , ndim_T, ndex_T ) ! vert. current
[11418]	959	ELSE
[12377]	960	CALL histwrite( nid_W, "vovecrtz", it, ww , ndim_T, ndex_T ) ! vert. current
[11418]	961	ENDIF
[3]	962	CALL histwrite( nid_W, "votkeavt", it, avt , ndim_T, ndex_T ) ! T vert. eddy diff. coef.
[9019]	963	CALL histwrite( nid_W, "votkeavm", it, avm , ndim_T, ndex_T ) ! T vert. eddy visc. coef.
	964	IF( ln_zdfddm ) THEN
	965	CALL histwrite( nid_W, "voddmavs", it, avs , ndim_T, ndex_T ) ! S vert. eddy diff. coef.
[3]	966	ENDIF
	967
[7646]	968	IF( ln_wave .AND. ln_sdw ) THEN
[9019]	969	CALL histwrite( nid_U, "sdzocrtx", it, usd , ndim_U , ndex_U ) ! i-StokesDrift-current
	970	CALL histwrite( nid_V, "sdmecrty", it, vsd , ndim_V , ndex_V ) ! j-StokesDrift-current
	971	CALL histwrite( nid_W, "sdvecrtz", it, wsd , ndim_T , ndex_T ) ! StokesDrift vert. current
[7646]	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 )
[12377]	981	IF(ln_abl) CALL histclo( nid_A )
[3]	982	ENDIF
[2528]	983	!
[9124]	984	IF( ln_timing ) CALL timing_stop('dia_wri')
[3294]	985	!
[3]	986	END SUBROUTINE dia_wri
[1567]	987	#endif
	988
[12377]	989	SUBROUTINE dia_wri_state( Kmm, cdfile_name )
[3]	990	!!---------------------------------------------------------------------
	991	!! * ROUTINE dia_wri_state *
	992	!!
	993	!! ** Purpose : create a NetCDF file named cdfile_name which contains
	994	!! the instantaneous ocean state and forcing fields.
	995	!! Used to find errors in the initial state or save the last
	996	!! ocean state in case of abnormal end of a simulation
	997	!!
	998	!! ** Method : NetCDF files using ioipsl
	999	!! File 'output.init.nc' is created if ninist = 1 (namelist)
	1000	!! File 'output.abort.nc' is created in case of abnormal job end
	1001	!!----------------------------------------------------------------------
[12377]	1002	INTEGER , INTENT( in ) :: Kmm ! time level index
[1334]	1003	CHARACTER (len=* ), INTENT( in ) :: cdfile_name ! name of the file created
[10425]	1004	!!
[12377]	1005	INTEGER :: inum, jk
[13237]	1006	REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t, zgdept ! 3D workspace !!st patch to use substitution
[3]	1007	!!----------------------------------------------------------------------
[3294]	1008	!
[13237]	1009	IF(lwp) THEN
	1010	WRITE(numout,*)
	1011	WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
	1012	WRITE(numout,*) '~~~~~~~~~~~~~ and forcing fields file created '
	1013	WRITE(numout,*) ' and named :', cdfile_name, '...nc'
	1014	ENDIF
[12649]	1015	!
[13237]	1016	DO jk = 1, jpk
	1017	ze3t(:,:,jk) = e3t(:,:,jk,Kmm)
	1018	zgdept(:,:,jk) = gdept(:,:,jk,Kmm)
	1019	END DO
	1020	!
[12649]	1021	CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE. )
	1022	!
[12377]	1023	CALL iom_rstput( 0, 0, inum, 'votemper', ts(:,:,:,jp_tem,Kmm) ) ! now temperature
	1024	CALL iom_rstput( 0, 0, inum, 'vosaline', ts(:,:,:,jp_sal,Kmm) ) ! now salinity
	1025	CALL iom_rstput( 0, 0, inum, 'sossheig', ssh(:,:,Kmm) ) ! sea surface height
	1026	CALL iom_rstput( 0, 0, inum, 'vozocrtx', uu(:,:,:,Kmm) ) ! now i-velocity
	1027	CALL iom_rstput( 0, 0, inum, 'vomecrty', vv(:,:,:,Kmm) ) ! now j-velocity
[11418]	1028	IF( ln_zad_Aimp ) THEN
[12377]	1029	CALL iom_rstput( 0, 0, inum, 'vovecrtz', ww + wi ) ! now k-velocity
[11418]	1030	ELSE
[12377]	1031	CALL iom_rstput( 0, 0, inum, 'vovecrtz', ww ) ! now k-velocity
[11418]	1032	ENDIF
[12377]	1033	CALL iom_rstput( 0, 0, inum, 'risfdep', risfdep ) ! now k-velocity
[13237]	1034	CALL iom_rstput( 0, 0, inum, 'ht' , ht(:,:) ) ! now water column height
[12649]	1035	!
[12377]	1036	IF ( ln_isf ) THEN
	1037	IF (ln_isfcav_mlt) THEN
	1038	CALL iom_rstput( 0, 0, inum, 'fwfisf_cav', fwfisf_cav ) ! now k-velocity
	1039	CALL iom_rstput( 0, 0, inum, 'rhisf_cav_tbl', rhisf_tbl_cav ) ! now k-velocity
	1040	CALL iom_rstput( 0, 0, inum, 'rfrac_cav_tbl', rfrac_tbl_cav ) ! now k-velocity
[12649]	1041	CALL iom_rstput( 0, 0, inum, 'misfkb_cav', REAL(misfkb_cav,wp) ) ! now k-velocity
	1042	CALL iom_rstput( 0, 0, inum, 'misfkt_cav', REAL(misfkt_cav,wp) ) ! now k-velocity
	1043	CALL iom_rstput( 0, 0, inum, 'mskisf_cav', REAL(mskisf_cav,wp), ktype = jp_i1 )
[12377]	1044	END IF
	1045	IF (ln_isfpar_mlt) THEN
[12649]	1046	CALL iom_rstput( 0, 0, inum, 'isfmsk_par', REAL(mskisf_par,wp) ) ! now k-velocity
[12377]	1047	CALL iom_rstput( 0, 0, inum, 'fwfisf_par', fwfisf_par ) ! now k-velocity
	1048	CALL iom_rstput( 0, 0, inum, 'rhisf_par_tbl', rhisf_tbl_par ) ! now k-velocity
	1049	CALL iom_rstput( 0, 0, inum, 'rfrac_par_tbl', rfrac_tbl_par ) ! now k-velocity
[12649]	1050	CALL iom_rstput( 0, 0, inum, 'misfkb_par', REAL(misfkb_par,wp) ) ! now k-velocity
	1051	CALL iom_rstput( 0, 0, inum, 'misfkt_par', REAL(misfkt_par,wp) ) ! now k-velocity
	1052	CALL iom_rstput( 0, 0, inum, 'mskisf_par', REAL(mskisf_par,wp), ktype = jp_i1 )
[12377]	1053	END IF
	1054	END IF
[12649]	1055	!
[7646]	1056	IF( ALLOCATED(ahtu) ) THEN
[10425]	1057	CALL iom_rstput( 0, 0, inum, 'ahtu', ahtu ) ! aht at u-point
	1058	CALL iom_rstput( 0, 0, inum, 'ahtv', ahtv ) ! aht at v-point
[7646]	1059	ENDIF
	1060	IF( ALLOCATED(ahmt) ) THEN
[10425]	1061	CALL iom_rstput( 0, 0, inum, 'ahmt', ahmt ) ! ahmt at u-point
	1062	CALL iom_rstput( 0, 0, inum, 'ahmf', ahmf ) ! ahmf at v-point
[7646]	1063	ENDIF
[10425]	1064	CALL iom_rstput( 0, 0, inum, 'sowaflup', emp - rnf ) ! freshwater budget
	1065	CALL iom_rstput( 0, 0, inum, 'sohefldo', qsr + qns ) ! total heat flux
	1066	CALL iom_rstput( 0, 0, inum, 'soshfldo', qsr ) ! solar heat flux
	1067	CALL iom_rstput( 0, 0, inum, 'soicecov', fr_i ) ! ice fraction
	1068	CALL iom_rstput( 0, 0, inum, 'sozotaux', utau ) ! i-wind stress
	1069	CALL iom_rstput( 0, 0, inum, 'sometauy', vtau ) ! j-wind stress
[6140]	1070	IF( .NOT.ln_linssh ) THEN
[13237]	1071	CALL iom_rstput( 0, 0, inum, 'vovvldep', zgdept ) ! T-cell depth
	1072	CALL iom_rstput( 0, 0, inum, 'vovvle3t', ze3t ) ! T-cell thickness
[10425]	1073	END IF
[7646]	1074	IF( ln_wave .AND. ln_sdw ) THEN
[10425]	1075	CALL iom_rstput( 0, 0, inum, 'sdzocrtx', usd ) ! now StokesDrift i-velocity
	1076	CALL iom_rstput( 0, 0, inum, 'sdmecrty', vsd ) ! now StokesDrift j-velocity
	1077	CALL iom_rstput( 0, 0, inum, 'sdvecrtz', wsd ) ! now StokesDrift k-velocity
[7646]	1078	ENDIF
[12377]	1079	IF ( ln_abl ) THEN
	1080	CALL iom_rstput ( 0, 0, inum, "uz1_abl", u_abl(:,:,2,nt_a ) ) ! now first level i-wind
	1081	CALL iom_rstput ( 0, 0, inum, "vz1_abl", v_abl(:,:,2,nt_a ) ) ! now first level j-wind
	1082	CALL iom_rstput ( 0, 0, inum, "tz1_abl", tq_abl(:,:,2,nt_a,1) ) ! now first level temperature
	1083	CALL iom_rstput ( 0, 0, inum, "qz1_abl", tq_abl(:,:,2,nt_a,2) ) ! now first level humidity
	1084	ENDIF
[12649]	1085	!
	1086	CALL iom_close( inum )
	1087	!
[10425]	1088	#if defined key_si3
	1089	IF( nn_ice == 2 ) THEN ! condition needed in case agrif + ice-model but no-ice in child grid
[12649]	1090	CALL iom_open( TRIM(cdfile_name)//'_ice', inum, ldwrt = .TRUE., kdlev = jpl, cdcomp = 'ICE' )
[10425]	1091	CALL ice_wri_state( inum )
[12649]	1092	CALL iom_close( inum )
[1561]	1093	ENDIF
[12933]	1094	!
[1561]	1095	#endif
[3]	1096	END SUBROUTINE dia_wri_state
[6140]	1097
[3]	1098	!!======================================================================
	1099	END MODULE diawri

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

Download in other formats: