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

source: branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 @ 3952

Last change on this file since 3952 was 3938, checked in by flavoni, 11 years ago

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