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

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

diawri.F90 in branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

Context Navigation

source: branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 @ 3419

Last change on this file since 3419 was 3419, checked in by acc, 12 years ago
Branch: dev_r3385_NOCS04_HAMF; #665. Stage 4 of 2012 development: Changes to get LIM3 working with embedded sea-ice. Working with a reduced (halved) timestep but exhibiting the same stability problems as LIM2_EVP with standard ORCA2 settings.
Property svn:keywords set to `Id`
File size: 42.7 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
	27	USE zdf_oce ! ocean vertical physics
	28	USE ldftra_oce ! ocean active tracers: lateral physics
	29	USE ldfdyn_oce ! ocean dynamics: lateral physics
[3294]	30	USE traldf_iso_grif, ONLY : psix_eiv, psiy_eiv
[3]	31	USE sol_oce ! solver variables
[888]	32	USE sbc_oce ! Surface boundary condition: ocean fields
	33	USE sbc_ice ! Surface boundary condition: ice fields
	34	USE sbcssr ! restoring term toward SST/SSS climatology
[3]	35	USE phycst ! physical constants
	36	USE zdfmxl ! mixed layer
	37	USE dianam ! build name of file (routine)
	38	USE zdfddm ! vertical physics: double diffusion
	39	USE diahth ! thermocline diagnostics
	40	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	41	USE in_out_manager ! I/O manager
[216]	42	USE diadimg ! dimg direct access file format output
[1756]	43	USE diaar5, ONLY : lk_diaar5
[1359]	44	USE iom
[460]	45	USE ioipsl
[1482]	46	#if defined key_lim2
	47	USE limwri_2
	48	#endif
[2715]	49	USE lib_mpp ! MPP library
[3294]	50	USE timing ! preformance summary
	51	USE wrk_nemo ! working array
[2528]	52
[3]	53	IMPLICIT NONE
	54	PRIVATE
	55
[2528]	56	PUBLIC dia_wri ! routines called by step.F90
	57	PUBLIC dia_wri_state
[2715]	58	PUBLIC dia_wri_alloc ! Called by nemogcm module
[3]	59
[2528]	60	INTEGER :: nid_T, nz_T, nh_T, ndim_T, ndim_hT ! grid_T file
	61	INTEGER :: nid_U, nz_U, nh_U, ndim_U, ndim_hU ! grid_U file
	62	INTEGER :: nid_V, nz_V, nh_V, ndim_V, ndim_hV ! grid_V file
	63	INTEGER :: nid_W, nz_W, nh_W ! grid_W file
	64	INTEGER :: ndex(1) ! ???
[2715]	65	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
	66	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V
[3]	67
	68	!! * Substitutions
	69	# include "zdfddm_substitute.h90"
[1756]	70	# include "domzgr_substitute.h90"
	71	# include "vectopt_loop_substitute.h90"
[3]	72	!!----------------------------------------------------------------------
[2528]	73	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
	74	!! $Id $
	75	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	76	!!----------------------------------------------------------------------
	77	CONTAINS
	78
[2715]	79	INTEGER FUNCTION dia_wri_alloc()
	80	!!----------------------------------------------------------------------
	81	INTEGER, DIMENSION(2) :: ierr
	82	!!----------------------------------------------------------------------
	83	!
	84	ierr = 0
	85	!
	86	ALLOCATE( ndex_hT(jpijpj) , ndex_T(jpijpj*jpk) , &
	87	& ndex_hU(jpijpj) , ndex_U(jpijpj*jpk) , &
	88	& ndex_hV(jpijpj) , ndex_V(jpijpj*jpk) , STAT=ierr(1) )
	89	!
	90	dia_wri_alloc = MAXVAL(ierr)
	91	IF( lk_mpp ) CALL mpp_sum( dia_wri_alloc )
	92	!
	93	END FUNCTION dia_wri_alloc
	94
[216]	95	#if defined key_dimgout
[3]	96	!!----------------------------------------------------------------------
[2528]	97	!! 'key_dimgout' DIMG output file
[3]	98	!!----------------------------------------------------------------------
	99	# include "diawri_dimg.h90"
	100
	101	#else
	102	!!----------------------------------------------------------------------
	103	!! Default option NetCDF output file
	104	!!----------------------------------------------------------------------
[2528]	105	# if defined key_iomput
[3]	106	!!----------------------------------------------------------------------
[2528]	107	!! 'key_iomput' use IOM library
	108	!!----------------------------------------------------------------------
[2715]	109
[1561]	110	SUBROUTINE dia_wri( kt )
[1482]	111	!!---------------------------------------------------------------------
	112	!! * ROUTINE dia_wri *
	113	!!
	114	!! ** Purpose : Standard output of opa: dynamics and tracer fields
	115	!! NETCDF format is used by default
	116	!!
	117	!! ** Method : use iom_put
	118	!!----------------------------------------------------------------------
[1756]	119	!!
[1482]	120	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[1756]	121	!!
	122	INTEGER :: ji, jj, jk ! dummy loop indices
	123	REAL(wp) :: zztmp, zztmpx, zztmpy !
[3294]	124	!!
	125	REAL(wp), POINTER, DIMENSION(:,:) :: z2d ! 2D workspace
	126	REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d ! 3D workspace
[1482]	127	!!----------------------------------------------------------------------
	128	!
[3294]	129	IF( nn_timing == 1 ) CALL timing_start('dia_wri')
	130	!
	131	CALL wrk_alloc( jpi , jpj , z2d )
	132	CALL wrk_alloc( jpi , jpj, jpk , z3d )
[2715]	133	!
[1482]	134	! Output the initial state and forcings
	135	IF( ninist == 1 ) THEN
	136	CALL dia_wri_state( 'output.init', kt )
	137	ninist = 0
	138	ENDIF
[3]	139
[3294]	140	CALL iom_put( "toce" , tsn(:,:,:,jp_tem) ) ! temperature
	141	CALL iom_put( "soce" , tsn(:,:,:,jp_sal) ) ! salinity
	142	CALL iom_put( "sst" , tsn(:,:,1,jp_tem) ) ! sea surface temperature
	143	CALL iom_put( "sst2" , tsn(:,:,1,jp_tem) * tsn(:,:,1,jp_tem) ) ! square of sea surface temperature
	144	CALL iom_put( "sss" , tsn(:,:,1,jp_sal) ) ! sea surface salinity
	145	CALL iom_put( "sss2" , tsn(:,:,1,jp_sal) * tsn(:,:,1,jp_sal) ) ! square of sea surface salinity
	146	CALL iom_put( "uoce" , un ) ! i-current
	147	CALL iom_put( "voce" , vn ) ! j-current
[1482]	148
[3294]	149	CALL iom_put( "avt" , avt ) ! T vert. eddy diff. coef.
	150	CALL iom_put( "avm" , avmu ) ! T vert. eddy visc. coef.
[1482]	151	IF( lk_zdfddm ) THEN
[3294]	152	CALL iom_put( "avs" , fsavs(:,:,:) ) ! S vert. eddy diff. coef.
[1482]	153	ENDIF
	154
[1756]	155	DO jj = 2, jpjm1 ! sst gradient
	156	DO ji = fs_2, fs_jpim1 ! vector opt.
[3294]	157	zztmp = tsn(ji,jj,1,jp_tem)
	158	zztmpx = ( tsn(ji+1,jj ,1,jp_tem) - zztmp ) / e1u(ji,jj) + ( zztmp - tsn(ji-1,jj ,1,jp_tem) ) / e1u(ji-1,jj )
	159	zztmpy = ( tsn(ji ,jj+1,1,jp_tem) - zztmp ) / e2v(ji,jj) + ( zztmp - tsn(ji ,jj-1,1,jp_tem) ) / e2v(ji ,jj-1)
[1756]	160	z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy ) &
	161	& * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1)
	162	END DO
	163	END DO
	164	CALL lbc_lnk( z2d, 'T', 1. )
[2561]	165	CALL iom_put( "sstgrad2", z2d ) ! square of module of sst gradient
[1756]	166	!CDIR NOVERRCHK
	167	z2d(:,:) = SQRT( z2d(:,:) )
[2561]	168	CALL iom_put( "sstgrad" , z2d ) ! module of sst gradient
[1756]	169
	170	IF( lk_diaar5 ) THEN
	171	z3d(:,:,jpk) = 0.e0
	172	DO jk = 1, jpkm1
	173	z3d(:,:,jk) = rau0 * un(:,:,jk) * e1u(:,:) * fse3u(:,:,jk)
	174	END DO
	175	CALL iom_put( "u_masstr", z3d ) ! mass transport in i-direction
	176	zztmp = 0.5 * rcp
	177	z2d(:,:) = 0.e0
	178	DO jk = 1, jpkm1
	179	DO jj = 2, jpjm1
	180	DO ji = fs_2, fs_jpim1 ! vector opt.
[3294]	181	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * zztmp * ( tsn(ji,jj,jk,jp_tem) + tsn(ji+1,jj,jk,jp_tem) )
[1756]	182	END DO
	183	END DO
	184	END DO
	185	CALL lbc_lnk( z2d, 'U', -1. )
	186	CALL iom_put( "u_heattr", z2d ) ! heat transport in i-direction
	187	DO jk = 1, jpkm1
	188	z3d(:,:,jk) = rau0 * vn(:,:,jk) * e2v(:,:) * fse3v(:,:,jk)
	189	END DO
	190	CALL iom_put( "v_masstr", z3d ) ! mass transport in j-direction
	191	z2d(:,:) = 0.e0
	192	DO jk = 1, jpkm1
	193	DO jj = 2, jpjm1
	194	DO ji = fs_2, fs_jpim1 ! vector opt.
[3294]	195	z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * zztmp * ( tsn(ji,jj,jk,jp_tem) + tsn(ji,jj+1,jk,jp_tem) )
[1756]	196	END DO
	197	END DO
	198	END DO
	199	CALL lbc_lnk( z2d, 'V', -1. )
	200	CALL iom_put( "v_heattr", z2d ) ! heat transport in i-direction
	201	ENDIF
[2528]	202	!
[3294]	203	CALL wrk_dealloc( jpi , jpj , z2d )
	204	CALL wrk_dealloc( jpi , jpj, jpk , z3d )
[2715]	205	!
[3294]	206	IF( nn_timing == 1 ) CALL timing_stop('dia_wri')
	207	!
[1482]	208	END SUBROUTINE dia_wri
	209
	210	#else
[2528]	211	!!----------------------------------------------------------------------
	212	!! Default option use IOIPSL library
	213	!!----------------------------------------------------------------------
	214
[1561]	215	SUBROUTINE dia_wri( kt )
[3]	216	!!---------------------------------------------------------------------
	217	!! * ROUTINE dia_wri *
	218	!!
	219	!! ** Purpose : Standard output of opa: dynamics and tracer fields
	220	!! NETCDF format is used by default
	221	!!
	222	!! ** Method : At the beginning of the first time step (nit000),
	223	!! define all the NETCDF files and fields
	224	!! At each time step call histdef to compute the mean if ncessary
	225	!! Each nwrite time step, output the instantaneous or mean fields
	226	!!----------------------------------------------------------------------
[2715]	227	!!
[3]	228	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[2528]	229	!!
	230	LOGICAL :: ll_print = .FALSE. ! =T print and flush numout
	231	CHARACTER (len=40) :: clhstnam, clop, clmx ! local names
	232	INTEGER :: inum = 11 ! temporary logical unit
[3294]	233	INTEGER :: ji, jj, jk ! dummy loop indices
	234	INTEGER :: ierr ! error code return from allocation
[2528]	235	INTEGER :: iimi, iima, ipk, it, itmod, ijmi, ijma ! local integers
	236	REAL(wp) :: zsto, zout, zmax, zjulian, zdt ! local scalars
[3294]	237	!!
	238	REAL(wp), POINTER, DIMENSION(:,:) :: zw2d ! 2D workspace
	239	REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d ! 3D workspace
[3]	240	!!----------------------------------------------------------------------
[3294]	241	!
	242	IF( nn_timing == 1 ) CALL timing_start('dia_wri')
[1482]	243	!
[3294]	244	CALL wrk_alloc( jpi , jpj , zw2d )
	245	IF ( ln_traldf_gdia ) call wrk_alloc( jpi , jpj , jpk , zw3d )
[2715]	246	!
[1482]	247	! Output the initial state and forcings
	248	IF( ninist == 1 ) THEN
	249	CALL dia_wri_state( 'output.init', kt )
	250	ninist = 0
	251	ENDIF
	252	!
[3]	253	! 0. Initialisation
	254	! -----------------
[632]	255
[3]	256	! local variable for debugging
	257	ll_print = .FALSE.
	258	ll_print = ll_print .AND. lwp
	259
	260	! Define frequency of output and means
	261	zdt = rdt
	262	IF( nacc == 1 ) zdt = rdtmin
[1312]	263	IF( ln_mskland ) THEN ; clop = "only(x)" ! put 1.e+20 on land (very expensive!!)
	264	ELSE ; clop = "x" ! no use of the mask value (require less cpu time)
	265	ENDIF
[3]	266	#if defined key_diainstant
	267	zsto = nwrite * zdt
[1312]	268	clop = "inst("//TRIM(clop)//")"
[3]	269	#else
	270	zsto=zdt
[1312]	271	clop = "ave("//TRIM(clop)//")"
[3]	272	#endif
	273	zout = nwrite * zdt
	274	zmax = ( nitend - nit000 + 1 ) * zdt
	275
	276	! Define indices of the horizontal output zoom and vertical limit storage
	277	iimi = 1 ; iima = jpi
	278	ijmi = 1 ; ijma = jpj
	279	ipk = jpk
	280
	281	! define time axis
[1334]	282	it = kt
	283	itmod = kt - nit000 + 1
[3]	284
	285
	286	! 1. Define NETCDF files and fields at beginning of first time step
	287	! -----------------------------------------------------------------
	288
	289	IF( kt == nit000 ) THEN
	290
	291	! Define the NETCDF files (one per grid)
[632]	292
[3]	293	! Compute julian date from starting date of the run
[1309]	294	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
	295	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
[3]	296	IF(lwp)WRITE(numout,*)
	297	IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear, &
	298	& ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
	299	IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma, &
	300	' limit storage in depth = ', ipk
	301
	302	! WRITE root name in date.file for use by postpro
[1581]	303	IF(lwp) THEN
[895]	304	CALL dia_nam( clhstnam, nwrite,' ' )
[1581]	305	CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
[895]	306	WRITE(inum,*) clhstnam
	307	CLOSE(inum)
	308	ENDIF
[632]	309
[3]	310	! Define the T grid FILE ( nid_T )
[632]	311
[3]	312	CALL dia_nam( clhstnam, nwrite, 'grid_T' )
	313	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
	314	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
	315	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	316	& nit000-1, zjulian, zdt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )
[3]	317	CALL histvert( nid_T, "deptht", "Vertical T levels", & ! Vertical grid: gdept
[1334]	318	& "m", ipk, gdept_0, nz_T, "down" )
[3]	319	! ! Index of ocean points
	320	CALL wheneq( jpijpjipk, tmask, 1, 1., ndex_T , ndim_T ) ! volume
	321	CALL wheneq( jpi*jpj , tmask, 1, 1., ndex_hT, ndim_hT ) ! surface
	322
	323	! Define the U grid FILE ( nid_U )
	324
	325	CALL dia_nam( clhstnam, nwrite, 'grid_U' )
	326	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
	327	CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu, & ! Horizontal grid: glamu and gphiu
	328	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	329	& nit000-1, zjulian, zdt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )
[3]	330	CALL histvert( nid_U, "depthu", "Vertical U levels", & ! Vertical grid: gdept
[1334]	331	& "m", ipk, gdept_0, nz_U, "down" )
[3]	332	! ! Index of ocean points
	333	CALL wheneq( jpijpjipk, umask, 1, 1., ndex_U , ndim_U ) ! volume
	334	CALL wheneq( jpi*jpj , umask, 1, 1., ndex_hU, ndim_hU ) ! surface
	335
	336	! Define the V grid FILE ( nid_V )
	337
	338	CALL dia_nam( clhstnam, nwrite, 'grid_V' ) ! filename
	339	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
	340	CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv, & ! Horizontal grid: glamv and gphiv
	341	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	342	& nit000-1, zjulian, zdt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )
[3]	343	CALL histvert( nid_V, "depthv", "Vertical V levels", & ! Vertical grid : gdept
[1334]	344	& "m", ipk, gdept_0, nz_V, "down" )
[3]	345	! ! Index of ocean points
	346	CALL wheneq( jpijpjipk, vmask, 1, 1., ndex_V , ndim_V ) ! volume
	347	CALL wheneq( jpi*jpj , vmask, 1, 1., ndex_hV, ndim_hV ) ! surface
	348
	349	! Define the W grid FILE ( nid_W )
	350
	351	CALL dia_nam( clhstnam, nwrite, 'grid_W' ) ! filename
	352	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
	353	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
	354	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
[2528]	355	& nit000-1, zjulian, zdt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set )
[3]	356	CALL histvert( nid_W, "depthw", "Vertical W levels", & ! Vertical grid: gdepw
[1334]	357	& "m", ipk, gdepw_0, nz_W, "down" )
[3]	358
[632]	359
[3]	360	! Declare all the output fields as NETCDF variables
	361
	362	! !!! nid_T : 3D
	363	CALL histdef( nid_T, "votemper", "Temperature" , "C" , & ! tn
	364	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	365	CALL histdef( nid_T, "vosaline", "Salinity" , "PSU" , & ! sn
	366	& jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
	367	! !!! nid_T : 2D
	368	CALL histdef( nid_T, "sosstsst", "Sea Surface temperature" , "C" , & ! sst
	369	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	370	CALL histdef( nid_T, "sosaline", "Sea Surface Salinity" , "PSU" , & ! sss
	371	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[359]	372	CALL histdef( nid_T, "sossheig", "Sea Surface Height" , "m" , & ! ssh
[3]	373	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1528]	374	!!$#if defined key_lim3 \|\| defined key_lim2
[888]	375	!!$ ! sowaflup = sowaflep + sorunoff + sowafldp + a term associated to
	376	!!$ ! internal damping to Levitus that can be diagnosed from others
	377	!!$ ! sowaflcd = sowaflep + sorunoff + sowafldp + iowaflup
	378	!!$ CALL histdef( nid_T, "iowaflup", "Ice=>ocean net freshwater" , "kg/m2/s", & ! fsalt
	379	!!$ & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	380	!!$ CALL histdef( nid_T, "sowaflep", "atmos=>ocean net freshwater" , "kg/m2/s", & ! fmass
	381	!!$ & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	382	!!$#endif
[2528]	383	CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux" , "Kg/m2/s", & ! (emp-rnf)
[3]	384	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[888]	385	!!$ CALL histdef( nid_T, "sorunoff", "Runoffs" , "Kg/m2/s", & ! runoffs
	386	!!$ & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[2528]	387	CALL histdef( nid_T, "sowaflcd", "concentration/dilution water flux" , "kg/m2/s", & ! (emps-rnf)
[3]	388	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[2528]	389	CALL histdef( nid_T, "sosalflx", "Surface Salt Flux" , "Kg/m2/s", & ! (emps-rnf) * sn
[3]	390	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[888]	391	CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! qns + qsr
[3]	392	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	393	CALL histdef( nid_T, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! qsr
	394	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1585]	395	CALL histdef( nid_T, "somixhgt", "Turbocline Depth" , "m" , & ! hmld
	396	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[3]	397	CALL histdef( nid_T, "somxl010", "Mixed Layer Depth 0.01" , "m" , & ! hmlp
	398	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1037]	399	CALL histdef( nid_T, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
[3]	400	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1649]	401	CALL histdef( nid_T, "sowindsp", "wind speed at 10m" , "m/s" , & ! wndm
	402	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[3]	403	#if ! defined key_coupled
	404	CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping" , "W/m2" , & ! qrp
	405	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	406	CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping" , "Kg/m2/s", & ! erp
	407	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	408	CALL histdef( nid_T, "sosafldp", "Surface salt flux: damping" , "Kg/m2/s", & ! erp * sn
	409	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	410	#endif
	411
[632]	412
	413
[888]	414	#if ( defined key_coupled && ! defined key_lim3 && ! defined key_lim2 )
[3]	415	CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping" , "W/m2" , & ! qrp
	416	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	417	CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping" , "Kg/m2/s", & ! erp
	418	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	419	CALL histdef( nid_T, "sosafldp", "Surface salt flux: Damping" , "Kg/m2/s", & ! erp * sn
	420	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	421	#endif
	422	clmx ="l_max(only(x))" ! max index on a period
	423	CALL histdef( nid_T, "sobowlin", "Bowl Index" , "W-point", & ! bowl INDEX
	424	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clmx, zsto, zout )
	425	#if defined key_diahth
	426	CALL histdef( nid_T, "sothedep", "Thermocline Depth" , "m" , & ! hth
	427	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	428	CALL histdef( nid_T, "so20chgt", "Depth of 20C isotherm" , "m" , & ! hd20
	429	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	430	CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm" , "m" , & ! hd28
	431	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	432	CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , "W" , & ! htc3
	433	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	434	#endif
	435
[888]	436	#if defined key_coupled
	437	# if defined key_lim3
	438	Must be adapted to LIM3
	439	# else
[3]	440	CALL histdef( nid_T,"soicetem" , "Ice Surface Temperature" , "K" , & ! tn_ice
	441	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	442	CALL histdef( nid_T,"soicealb" , "Ice Albedo" , "[0,1]" , & ! alb_ice
	443	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[888]	444	# endif
[3]	445	#endif
	446
[2528]	447	CALL histend( nid_T, snc4chunks=snc4set )
[3]	448
	449	! !!! nid_U : 3D
	450	CALL histdef( nid_U, "vozocrtx", "Zonal Current" , "m/s" , & ! un
	451	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
[3294]	452	IF( ln_traldf_gdia ) THEN
	453	CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current" , "m/s" , & ! u_eiv
	454	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
	455	ELSE
[3]	456	#if defined key_diaeiv
[3294]	457	CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current" , "m/s" , & ! u_eiv
[3]	458	& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
	459	#endif
[3294]	460	END IF
[3]	461	! !!! nid_U : 2D
[888]	462	CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau
[3]	463	& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
	464
[2528]	465	CALL histend( nid_U, snc4chunks=snc4set )
[3]	466
	467	! !!! nid_V : 3D
	468	CALL histdef( nid_V, "vomecrty", "Meridional Current" , "m/s" , & ! vn
	469	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
[3294]	470	IF( ln_traldf_gdia ) THEN
	471	CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current" , "m/s" , & ! v_eiv
	472	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
	473	ELSE
[3]	474	#if defined key_diaeiv
[3294]	475	CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current" , "m/s" , & ! v_eiv
[3]	476	& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
	477	#endif
[3294]	478	END IF
[3]	479	! !!! nid_V : 2D
[888]	480	CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau
[3]	481	& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
	482
[2528]	483	CALL histend( nid_V, snc4chunks=snc4set )
[3]	484
	485	! !!! nid_W : 3D
	486	CALL histdef( nid_W, "vovecrtz", "Vertical Velocity" , "m/s" , & ! wn
	487	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[3294]	488	IF( ln_traldf_gdia ) THEN
	489	CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity" , "m/s" , & ! w_eiv
	490	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	491	ELSE
[3]	492	#if defined key_diaeiv
[3294]	493	CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity" , "m/s" , & ! w_eiv
	494	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[3]	495	#endif
[3294]	496	END IF
[3]	497	CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity" , "m2/s" , & ! avt
	498	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[255]	499	CALL histdef( nid_W, "votkeavm", "Vertical Eddy Viscosity" , "m2/s" , & ! avmu
	500	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	501
[3]	502	IF( lk_zdfddm ) THEN
	503	CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity" , "m2/s" , & ! avs
	504	& jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
	505	ENDIF
	506	! !!! nid_W : 2D
	507	#if defined key_traldf_c2d
	508	CALL histdef( nid_W, "soleahtw", "lateral eddy diffusivity" , "m2/s" , & ! ahtw
	509	& jpi, jpj, nh_W, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
[23]	510	# if defined key_traldf_eiv
[3]	511	CALL histdef( nid_W, "soleaeiw", "eddy induced vel. coeff. at w-point", "m2/s", & ! aeiw
	512	& jpi, jpj, nh_W, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
[23]	513	# endif
[3]	514	#endif
	515
[2528]	516	CALL histend( nid_W, snc4chunks=snc4set )
[3]	517
	518	IF(lwp) WRITE(numout,*)
	519	IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
	520	IF(ll_print) CALL FLUSH(numout )
	521
	522	ENDIF
	523
	524	! 2. Start writing data
	525	! ---------------------
	526
	527	! ndex(1) est utilise ssi l'avant dernier argument est diffferent de
	528	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
	529	! donne le nombre d'elements, et ndex la liste des indices a sortir
	530
[1334]	531	IF( lwp .AND. MOD( itmod, nwrite ) == 0 ) THEN
[3]	532	WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
	533	WRITE(numout,*) '~~~~~~ '
	534	ENDIF
	535
	536	! Write fields on T grid
[3294]	537	CALL histwrite( nid_T, "votemper", it, tsn(:,:,:,jp_tem), ndim_T , ndex_T ) ! temperature
	538	CALL histwrite( nid_T, "vosaline", it, tsn(:,:,:,jp_sal), ndim_T , ndex_T ) ! salinity
	539	CALL histwrite( nid_T, "sosstsst", it, tsn(:,:,1,jp_tem), ndim_hT, ndex_hT ) ! sea surface temperature
	540	CALL histwrite( nid_T, "sosaline", it, tsn(:,:,1,jp_sal), ndim_hT, ndex_hT ) ! sea surface salinity
[3]	541	CALL histwrite( nid_T, "sossheig", it, sshn , ndim_hT, ndex_hT ) ! sea surface height
[1528]	542	!!$#if defined key_lim3 \|\| defined key_lim2
[888]	543	!!$ CALL histwrite( nid_T, "iowaflup", it, fsalt(:,:) , ndim_hT, ndex_hT ) ! ice=>ocean water flux
	544	!!$ CALL histwrite( nid_T, "sowaflep", it, fmass(:,:) , ndim_hT, ndex_hT ) ! atmos=>ocean water flux
	545	!!$#endif
[2528]	546	CALL histwrite( nid_T, "sowaflup", it, ( emp-rnf ) , ndim_hT, ndex_hT ) ! upward water flux
[888]	547	!!$ CALL histwrite( nid_T, "sorunoff", it, runoff , ndim_hT, ndex_hT ) ! runoff
[3402]	548	CALL histwrite( nid_T, "sowaflcd", it, ( sfx -rnf ) , ndim_hT, ndex_hT ) ! c/d water flux
	549	zw2d(:,:) = ( sfx (:,:) - rnf(:,:) ) * tsn(:,:,1,jp_sal) * tmask(:,:,1)
[3]	550	CALL histwrite( nid_T, "sosalflx", it, zw2d , ndim_hT, ndex_hT ) ! c/d salt flux
[888]	551	CALL histwrite( nid_T, "sohefldo", it, qns + qsr , ndim_hT, ndex_hT ) ! total heat flux
[3]	552	CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux
[1585]	553	CALL histwrite( nid_T, "somixhgt", it, hmld , ndim_hT, ndex_hT ) ! turbocline depth
[3]	554	CALL histwrite( nid_T, "somxl010", it, hmlp , ndim_hT, ndex_hT ) ! mixed layer depth
[1037]	555	CALL histwrite( nid_T, "soicecov", it, fr_i , ndim_hT, ndex_hT ) ! ice fraction
[1649]	556	CALL histwrite( nid_T, "sowindsp", it, wndm , ndim_hT, ndex_hT ) ! wind speed
[3]	557	#if ! defined key_coupled
	558	CALL histwrite( nid_T, "sohefldp", it, qrp , ndim_hT, ndex_hT ) ! heat flux damping
	559	CALL histwrite( nid_T, "sowafldp", it, erp , ndim_hT, ndex_hT ) ! freshwater flux damping
[3294]	560	IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1)
[3]	561	CALL histwrite( nid_T, "sosafldp", it, zw2d , ndim_hT, ndex_hT ) ! salt flux damping
	562	#endif
[888]	563	#if ( defined key_coupled && ! defined key_lim3 && ! defined key_lim2 )
[3]	564	CALL histwrite( nid_T, "sohefldp", it, qrp , ndim_hT, ndex_hT ) ! heat flux damping
	565	CALL histwrite( nid_T, "sowafldp", it, erp , ndim_hT, ndex_hT ) ! freshwater flux damping
[3294]	566	IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1)
[3]	567	CALL histwrite( nid_T, "sosafldp", it, zw2d , ndim_hT, ndex_hT ) ! salt flux damping
	568	#endif
[2528]	569	zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)
[3]	570	CALL histwrite( nid_T, "sobowlin", it, zw2d , ndim_hT, ndex_hT ) ! ???
	571
	572	#if defined key_diahth
	573	CALL histwrite( nid_T, "sothedep", it, hth , ndim_hT, ndex_hT ) ! depth of the thermocline
	574	CALL histwrite( nid_T, "so20chgt", it, hd20 , ndim_hT, ndex_hT ) ! depth of the 20 isotherm
	575	CALL histwrite( nid_T, "so28chgt", it, hd28 , ndim_hT, ndex_hT ) ! depth of the 28 isotherm
	576	CALL histwrite( nid_T, "sohtc300", it, htc3 , ndim_hT, ndex_hT ) ! first 300m heaat content
	577	#endif
[888]	578
	579	#if defined key_coupled
	580	# if defined key_lim3
	581	Must be adapted for LIM3
[3]	582	CALL histwrite( nid_T, "soicetem", it, tn_ice , ndim_hT, ndex_hT ) ! surf. ice temperature
	583	CALL histwrite( nid_T, "soicealb", it, alb_ice , ndim_hT, ndex_hT ) ! ice albedo
[888]	584	# else
[1463]	585	CALL histwrite( nid_T, "soicetem", it, tn_ice(:,:,1) , ndim_hT, ndex_hT ) ! surf. ice temperature
	586	CALL histwrite( nid_T, "soicealb", it, alb_ice(:,:,1), ndim_hT, ndex_hT ) ! ice albedo
[888]	587	# endif
[3]	588	#endif
	589	! Write fields on U grid
	590	CALL histwrite( nid_U, "vozocrtx", it, un , ndim_U , ndex_U ) ! i-current
[3294]	591	IF( ln_traldf_gdia ) THEN
	592	IF (.not. ALLOCATED(psix_eiv))THEN
	593	ALLOCATE( psix_eiv(jpi,jpj,jpk) , psiy_eiv(jpi,jpj,jpk) , STAT=ierr )
	594	IF( lk_mpp ) CALL mpp_sum ( ierr )
	595	IF( ierr > 0 ) CALL ctl_stop('STOP', 'diawri: unable to allocate psi{x,y}_eiv')
	596	psix_eiv(:,:,:) = 0.0_wp
	597	psiy_eiv(:,:,:) = 0.0_wp
	598	ENDIF
	599	DO jk=1,jpkm1
	600	zw3d(:,:,jk) = (psix_eiv(:,:,jk+1) - psix_eiv(:,:,jk))/fse3u(:,:,jk) ! u_eiv = -dpsix/dz
	601	END DO
	602	zw3d(:,:,jpk) = 0._wp
	603	CALL histwrite( nid_U, "vozoeivu", it, zw3d, ndim_U , ndex_U ) ! i-eiv current
	604	ELSE
[3]	605	#if defined key_diaeiv
[3294]	606	CALL histwrite( nid_U, "vozoeivu", it, u_eiv, ndim_U , ndex_U ) ! i-eiv current
[3]	607	#endif
[3294]	608	ENDIF
[888]	609	CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress
[3]	610
	611	! Write fields on V grid
	612	CALL histwrite( nid_V, "vomecrty", it, vn , ndim_V , ndex_V ) ! j-current
[3294]	613	IF( ln_traldf_gdia ) THEN
	614	DO jk=1,jpk-1
	615	zw3d(:,:,jk) = (psiy_eiv(:,:,jk+1) - psiy_eiv(:,:,jk))/fse3v(:,:,jk) ! v_eiv = -dpsiy/dz
	616	END DO
	617	zw3d(:,:,jpk) = 0._wp
	618	CALL histwrite( nid_V, "vomeeivv", it, zw3d, ndim_V , ndex_V ) ! j-eiv current
	619	ELSE
[3]	620	#if defined key_diaeiv
[3294]	621	CALL histwrite( nid_V, "vomeeivv", it, v_eiv, ndim_V , ndex_V ) ! j-eiv current
[3]	622	#endif
[3294]	623	ENDIF
[888]	624	CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress
[3]	625
	626	! Write fields on W grid
	627	CALL histwrite( nid_W, "vovecrtz", it, wn , ndim_T, ndex_T ) ! vert. current
[3294]	628	IF( ln_traldf_gdia ) THEN
	629	DO jk=1,jpk-1
	630	DO jj = 2, jpjm1
	631	DO ji = fs_2, fs_jpim1 ! vector opt.
	632	zw3d(ji,jj,jk) = (psiy_eiv(ji,jj,jk) - psiy_eiv(ji,jj-1,jk))/e2v(ji,jj) + &
	633	& (psix_eiv(ji,jj,jk) - psix_eiv(ji-1,jj,jk))/e1u(ji,jj) ! w_eiv = dpsiy/dy + dpsiy/dx
	634	END DO
	635	END DO
	636	END DO
	637	zw3d(:,:,jpk) = 0._wp
	638	CALL histwrite( nid_W, "voveeivw", it, zw3d , ndim_T, ndex_T ) ! vert. eiv current
	639	ELSE
[3]	640	# if defined key_diaeiv
[3294]	641	CALL histwrite( nid_W, "voveeivw", it, w_eiv , ndim_T, ndex_T ) ! vert. eiv current
[3]	642	# endif
[3294]	643	ENDIF
[3]	644	CALL histwrite( nid_W, "votkeavt", it, avt , ndim_T, ndex_T ) ! T vert. eddy diff. coef.
[255]	645	CALL histwrite( nid_W, "votkeavm", it, avmu , ndim_T, ndex_T ) ! T vert. eddy visc. coef.
[3]	646	IF( lk_zdfddm ) THEN
	647	CALL histwrite( nid_W, "voddmavs", it, fsavs(:,:,:), ndim_T, ndex_T ) ! S vert. eddy diff. coef.
	648	ENDIF
	649	#if defined key_traldf_c2d
	650	CALL histwrite( nid_W, "soleahtw", it, ahtw , ndim_hT, ndex_hT ) ! lateral eddy diff. coef.
[23]	651	# if defined key_traldf_eiv
[3]	652	CALL histwrite( nid_W, "soleaeiw", it, aeiw , ndim_hT, ndex_hT ) ! EIV coefficient at w-point
[23]	653	# endif
[3]	654	#endif
	655
[1318]	656	! 3. Close all files
	657	! ---------------------------------------
[1561]	658	IF( kt == nitend ) THEN
[3]	659	CALL histclo( nid_T )
	660	CALL histclo( nid_U )
	661	CALL histclo( nid_V )
	662	CALL histclo( nid_W )
	663	ENDIF
[2528]	664	!
[3294]	665	CALL wrk_dealloc( jpi , jpj , zw2d )
	666	IF ( ln_traldf_gdia ) call wrk_dealloc( jpi , jpj , jpk , zw3d )
[2715]	667	!
[3294]	668	IF( nn_timing == 1 ) CALL timing_stop('dia_wri')
	669	!
[3]	670	END SUBROUTINE dia_wri
[1482]	671	# endif
[3]	672
[1567]	673	#endif
	674
[1334]	675	SUBROUTINE dia_wri_state( cdfile_name, kt )
[3]	676	!!---------------------------------------------------------------------
	677	!! * ROUTINE dia_wri_state *
	678	!!
	679	!! ** Purpose : create a NetCDF file named cdfile_name which contains
	680	!! the instantaneous ocean state and forcing fields.
	681	!! Used to find errors in the initial state or save the last
	682	!! ocean state in case of abnormal end of a simulation
	683	!!
	684	!! ** Method : NetCDF files using ioipsl
	685	!! File 'output.init.nc' is created if ninist = 1 (namelist)
	686	!! File 'output.abort.nc' is created in case of abnormal job end
	687	!!----------------------------------------------------------------------
[1334]	688	CHARACTER (len=* ), INTENT( in ) :: cdfile_name ! name of the file created
	689	INTEGER , INTENT( in ) :: kt ! ocean time-step index
[2528]	690	!!
[648]	691	CHARACTER (len=32) :: clname
[3]	692	CHARACTER (len=40) :: clop
[2528]	693	INTEGER :: id_i , nz_i, nh_i
	694	INTEGER, DIMENSION(1) :: idex ! local workspace
	695	REAL(wp) :: zsto, zout, zmax, zjulian, zdt
[3]	696	!!----------------------------------------------------------------------
[3294]	697	!
[3419]	698	! IF( nn_timing == 1 ) CALL timing_start('dia_wri_state') ! not sure this works for routines not called in first timestep
[3]	699
	700	! 0. Initialisation
	701	! -----------------
[632]	702
[648]	703	! Define name, frequency of output and means
	704	clname = cdfile_name
[1792]	705	IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
[3]	706	zdt = rdt
	707	zsto = rdt
	708	clop = "inst(x)" ! no use of the mask value (require less cpu time)
	709	zout = rdt
	710	zmax = ( nitend - nit000 + 1 ) * zdt
	711
[648]	712	IF(lwp) WRITE(numout,*)
	713	IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
	714	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ and forcing fields file created '
	715	IF(lwp) WRITE(numout,*) ' and named :', clname, '.nc'
	716
	717
[3]	718	! 1. Define NETCDF files and fields at beginning of first time step
	719	! -----------------------------------------------------------------
	720
	721	! Compute julian date from starting date of the run
[1310]	722	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian ) ! time axis
	723	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
[648]	724	CALL histbeg( clname, jpi, glamt, jpj, gphit, &
[2528]	725	1, jpi, 1, jpj, nit000-1, zjulian, zdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
[3]	726	CALL histvert( id_i, "deptht", "Vertical T levels", & ! Vertical grid : gdept
[1334]	727	"m", jpk, gdept_0, nz_i, "down")
[3]	728
	729	! Declare all the output fields as NetCDF variables
	730
	731	CALL histdef( id_i, "vosaline", "Salinity" , "PSU" , & ! salinity
	732	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	733	CALL histdef( id_i, "votemper", "Temperature" , "C" , & ! temperature
	734	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
[359]	735	CALL histdef( id_i, "sossheig", "Sea Surface Height" , "m" , & ! ssh
	736	& jpi, jpj, nh_i, 1 , 1, 1 , nz_i, 32, clop, zsto, zout )
[3]	737	CALL histdef( id_i, "vozocrtx", "Zonal Current" , "m/s" , & ! zonal current
	738	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	739	CALL histdef( id_i, "vomecrty", "Meridional Current" , "m/s" , & ! meridonal current
	740	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	741	CALL histdef( id_i, "vovecrtz", "Vertical Velocity" , "m/s" , & ! vertical current
	742	& jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
	743	CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S", & ! net freshwater
	744	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	745	CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2" , & ! net heat flux
	746	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	747	CALL histdef( id_i, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! solar flux
	748	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
[1037]	749	CALL histdef( id_i, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
[3]	750	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	751	CALL histdef( id_i, "sozotaux", "Zonal Wind Stress" , "N/m2" , & ! i-wind stress
	752	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	753	CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2" , & ! j-wind stress
	754	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
	755
[1482]	756	#if defined key_lim2
	757	CALL lim_wri_state_2( kt, id_i, nh_i )
	758	#else
[2528]	759	CALL histend( id_i, snc4chunks=snc4set )
[1482]	760	#endif
[3]	761
	762	! 2. Start writing data
	763	! ---------------------
	764	! idex(1) est utilise ssi l'avant dernier argument est diffferent de
	765	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
	766	! donne le nombre d'elements, et idex la liste des indices a sortir
	767	idex(1) = 1 ! init to avoid compil warning
[632]	768
[3]	769	! Write all fields on T grid
[3294]	770	CALL histwrite( id_i, "votemper", kt, tsn(:,:,:,jp_tem), jpijpjjpk, idex ) ! now temperature
	771	CALL histwrite( id_i, "vosaline", kt, tsn(:,:,:,jp_sal), jpijpjjpk, idex ) ! now salinity
	772	CALL histwrite( id_i, "sossheig", kt, sshn , jpi*jpj , idex ) ! sea surface height
	773	CALL histwrite( id_i, "vozocrtx", kt, un , jpijpjjpk, idex ) ! now i-velocity
	774	CALL histwrite( id_i, "vomecrty", kt, vn , jpijpjjpk, idex ) ! now j-velocity
	775	CALL histwrite( id_i, "vovecrtz", kt, wn , jpijpjjpk, idex ) ! now k-velocity
	776	CALL histwrite( id_i, "sowaflup", kt, (emp-rnf ) , jpi*jpj , idex ) ! freshwater budget
	777	CALL histwrite( id_i, "sohefldo", kt, qsr + qns , jpi*jpj , idex ) ! total heat flux
	778	CALL histwrite( id_i, "soshfldo", kt, qsr , jpi*jpj , idex ) ! solar heat flux
	779	CALL histwrite( id_i, "soicecov", kt, fr_i , jpi*jpj , idex ) ! ice fraction
	780	CALL histwrite( id_i, "sozotaux", kt, utau , jpi*jpj , idex ) ! i-wind stress
	781	CALL histwrite( id_i, "sometauy", kt, vtau , jpi*jpj , idex ) ! j-wind stress
[3]	782
	783	! 3. Close the file
	784	! -----------------
	785	CALL histclo( id_i )
[1567]	786	#if ! defined key_iomput && ! defined key_dimgout
[1561]	787	IF( ninist /= 1 ) THEN
	788	CALL histclo( nid_T )
	789	CALL histclo( nid_U )
	790	CALL histclo( nid_V )
	791	CALL histclo( nid_W )
	792	ENDIF
	793	#endif
[3294]	794
[3419]	795	! IF( nn_timing == 1 ) CALL timing_stop('dia_wri_state') ! not sure this works for routines not called in first timestep
[3294]	796	!
[3]	797
	798	END SUBROUTINE dia_wri_state
	799	!!======================================================================
	800	END MODULE diawri

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

Download in other formats: