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

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 @ 3432

Last change on this file since 3432 was 3432, checked in by trackstand2, 12 years ago

Merge branch 'ksection_partition'

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