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/UKMO/r6232_collate_bgc_diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

source: branches/UKMO/r6232_collate_bgc_diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 @ 11132

Last change on this file since 11132 was 11132, checked in by jcastill, 5 years ago

Remove svn keywords

File size: 56.6 KB
RevLine 
[3]1MODULE diawri
2   !!======================================================================
3   !!                     ***  MODULE  diawri  ***
4   !! Ocean diagnostics :  write ocean output files
5   !!=====================================================================
[2528]6   !! History :  OPA  ! 1991-03  (M.-A. Foujols)  Original code
7   !!            4.0  ! 1991-11  (G. Madec)
8   !!                 ! 1992-06  (M. Imbard)  correction restart file
9   !!                 ! 1992-07  (M. Imbard)  split into diawri and rstwri
10   !!                 ! 1993-03  (M. Imbard)  suppress writibm
11   !!                 ! 1998-01  (C. Levy)  NETCDF format using ioipsl INTERFACE
12   !!                 ! 1999-02  (E. Guilyardi)  name of netCDF files + variables
13   !!            8.2  ! 2000-06  (M. Imbard)  Original code (diabort.F)
14   !!   NEMO     1.0  ! 2002-06  (A.Bozec, E. Durand)  Original code (diainit.F)
15   !!             -   ! 2002-09  (G. Madec)  F90: Free form and module
16   !!             -   ! 2002-12  (G. Madec)  merge of diabort and diainit, F90
17   !!                 ! 2005-11  (V. Garnier) Surface pressure gradient organization
18   !!            3.2  ! 2008-11  (B. Lemaire) creation from old diawri
19   !!----------------------------------------------------------------------
[3]20
21   !!----------------------------------------------------------------------
[2528]22   !!   dia_wri       : create the standart output files
23   !!   dia_wri_state : create an output NetCDF file for a single instantaeous ocean state and forcing fields
24   !!----------------------------------------------------------------------
[3]25   USE oce             ! ocean dynamics and tracers
26   USE dom_oce         ! ocean space and time domain
[4292]27   USE dynadv, ONLY: ln_dynadv_vec
[3]28   USE zdf_oce         ! ocean vertical physics
29   USE ldftra_oce      ! ocean active tracers: lateral physics
30   USE ldfdyn_oce      ! ocean dynamics: lateral physics
[3294]31   USE traldf_iso_grif, ONLY : psix_eiv, psiy_eiv
[3]32   USE sol_oce         ! solver variables
[888]33   USE sbc_oce         ! Surface boundary condition: ocean fields
34   USE sbc_ice         ! Surface boundary condition: ice fields
[3609]35   USE icb_oce         ! Icebergs
36   USE icbdia          ! Iceberg budgets
[888]37   USE sbcssr          ! restoring term toward SST/SSS climatology
[3]38   USE phycst          ! physical constants
39   USE zdfmxl          ! mixed layer
40   USE dianam          ! build name of file (routine)
41   USE zdfddm          ! vertical  physics: double diffusion
42   USE diahth          ! thermocline diagnostics
43   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
44   USE in_out_manager  ! I/O manager
[216]45   USE diadimg         ! dimg direct access file format output
[1359]46   USE iom
[460]47   USE ioipsl
[5463]48   USE dynspg_oce, ONLY: un_adv, vn_adv ! barotropic velocities     
49
[1482]50#if defined key_lim2
51   USE limwri_2 
[4161]52#elif defined key_lim3
53   USE limwri 
[1482]54#endif
[2715]55   USE lib_mpp         ! MPP library
[3294]56   USE timing          ! preformance summary
57   USE wrk_nemo        ! working array
[2528]58
[3]59   IMPLICIT NONE
60   PRIVATE
61
[2528]62   PUBLIC   dia_wri                 ! routines called by step.F90
63   PUBLIC   dia_wri_state
[2715]64   PUBLIC   dia_wri_alloc           ! Called by nemogcm module
[3]65
[2528]66   INTEGER ::   nid_T, nz_T, nh_T, ndim_T, ndim_hT   ! grid_T file
[3609]67   INTEGER ::          nb_T              , ndim_bT   ! grid_T file
[2528]68   INTEGER ::   nid_U, nz_U, nh_U, ndim_U, ndim_hU   ! grid_U file
69   INTEGER ::   nid_V, nz_V, nh_V, ndim_V, ndim_hV   ! grid_V file
70   INTEGER ::   nid_W, nz_W, nh_W                    ! grid_W file
71   INTEGER ::   ndex(1)                              ! ???
[2715]72   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
73   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V
[3609]74   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_bT
[3]75
76   !! * Substitutions
77#  include "zdfddm_substitute.h90"
[1756]78#  include "domzgr_substitute.h90"
79#  include "vectopt_loop_substitute.h90"
[3]80   !!----------------------------------------------------------------------
[2528]81   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
[5217]82   !! $Id$
[2528]83   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
[3]84   !!----------------------------------------------------------------------
85CONTAINS
86
[2715]87   INTEGER FUNCTION dia_wri_alloc()
88      !!----------------------------------------------------------------------
89      INTEGER, DIMENSION(2) :: ierr
90      !!----------------------------------------------------------------------
91      ierr = 0
92      ALLOCATE( ndex_hT(jpi*jpj) , ndex_T(jpi*jpj*jpk) ,     &
93         &      ndex_hU(jpi*jpj) , ndex_U(jpi*jpj*jpk) ,     &
94         &      ndex_hV(jpi*jpj) , ndex_V(jpi*jpj*jpk) , STAT=ierr(1) )
95         !
96      dia_wri_alloc = MAXVAL(ierr)
97      IF( lk_mpp )   CALL mpp_sum( dia_wri_alloc )
98      !
99  END FUNCTION dia_wri_alloc
100
[216]101#if defined key_dimgout
[3]102   !!----------------------------------------------------------------------
[2528]103   !!   'key_dimgout'                                      DIMG output file
[3]104   !!----------------------------------------------------------------------
105#   include "diawri_dimg.h90"
106
107#else
108   !!----------------------------------------------------------------------
109   !!   Default option                                   NetCDF output file
110   !!----------------------------------------------------------------------
[2528]111# if defined key_iomput
[3]112   !!----------------------------------------------------------------------
[2528]113   !!   'key_iomput'                                        use IOM library
114   !!----------------------------------------------------------------------
[2715]115
[1561]116   SUBROUTINE dia_wri( kt )
[1482]117      !!---------------------------------------------------------------------
118      !!                  ***  ROUTINE dia_wri  ***
119      !!                   
120      !! ** Purpose :   Standard output of opa: dynamics and tracer fields
121      !!      NETCDF format is used by default
122      !!
123      !! ** Method  :  use iom_put
124      !!----------------------------------------------------------------------
[1756]125      !!
[1482]126      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
[1756]127      !!
128      INTEGER                      ::   ji, jj, jk              ! dummy loop indices
[5460]129      INTEGER                      ::   jkbot                   !
[1756]130      REAL(wp)                     ::   zztmp, zztmpx, zztmpy   !
[3294]131      !!
[4761]132      REAL(wp), POINTER, DIMENSION(:,:)   :: z2d      ! 2D workspace
[3294]133      REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d      ! 3D workspace
[1482]134      !!----------------------------------------------------------------------
135      !
[3294]136      IF( nn_timing == 1 )   CALL timing_start('dia_wri')
137      !
[4990]138      CALL wrk_alloc( jpi , jpj      , z2d )
[3294]139      CALL wrk_alloc( jpi , jpj, jpk , z3d )
[2715]140      !
[1482]141      ! Output the initial state and forcings
142      IF( ninist == 1 ) THEN                       
143         CALL dia_wri_state( 'output.init', kt )
144         ninist = 0
145      ENDIF
[3]146
[5107]147      IF( .NOT.lk_vvl ) THEN
148         CALL iom_put( "e3t" , fse3t_n(:,:,:) )
149         CALL iom_put( "e3u" , fse3u_n(:,:,:) )
150         CALL iom_put( "e3v" , fse3v_n(:,:,:) )
151         CALL iom_put( "e3w" , fse3w_n(:,:,:) )
152      ENDIF
[5461]153
154      CALL iom_put( "ssh" , sshn )                 ! sea surface height
155      if( iom_use('ssh2') )   CALL iom_put( "ssh2", sshn(:,:) * sshn(:,:) )   ! square of sea surface height
[5107]156     
157      CALL iom_put( "toce", tsn(:,:,:,jp_tem) )    ! 3D temperature
158      CALL iom_put(  "sst", tsn(:,:,1,jp_tem) )    ! surface temperature
159      IF ( iom_use("sbt") ) THEN
[4990]160         DO jj = 1, jpj
161            DO ji = 1, jpi
[5460]162               jkbot = mbkt(ji,jj)
163               z2d(ji,jj) = tsn(ji,jj,jkbot,jp_tem)
[4990]164            END DO
[5107]165         END DO
166         CALL iom_put( "sbt", z2d )                ! bottom temperature
167      ENDIF
168     
169      CALL iom_put( "soce", tsn(:,:,:,jp_sal) )    ! 3D salinity
170      CALL iom_put(  "sss", tsn(:,:,1,jp_sal) )    ! surface salinity
171      IF ( iom_use("sbs") ) THEN
[4990]172         DO jj = 1, jpj
173            DO ji = 1, jpi
[5460]174               jkbot = mbkt(ji,jj)
175               z2d(ji,jj) = tsn(ji,jj,jkbot,jp_sal)
[4990]176            END DO
[5107]177         END DO
178         CALL iom_put( "sbs", z2d )                ! bottom salinity
179      ENDIF
[5463]180
181      IF ( iom_use("taubot") ) THEN                ! bottom stress
182         z2d(:,:) = 0._wp
183         DO jj = 2, jpjm1
184            DO ji = fs_2, fs_jpim1   ! vector opt.
185               zztmpx = (  bfrua(ji  ,jj) * un(ji  ,jj,mbku(ji  ,jj))  &
186                      &  + bfrua(ji-1,jj) * un(ji-1,jj,mbku(ji-1,jj))  )     
187               zztmpy = (  bfrva(ji,  jj) * vn(ji,jj  ,mbkv(ji,jj  ))  &
188                      &  + bfrva(ji,jj-1) * vn(ji,jj-1,mbkv(ji,jj-1))  ) 
189               z2d(ji,jj) = rau0 * SQRT( zztmpx * zztmpx + zztmpy * zztmpy ) * tmask(ji,jj,1) 
190               !
191            ENDDO
192         ENDDO
193         CALL lbc_lnk( z2d, 'T', 1. )
194         CALL iom_put( "taubot", z2d )           
195      ENDIF
[5107]196         
197      CALL iom_put( "uoce", un(:,:,:)         )    ! 3D i-current
198      CALL iom_put(  "ssu", un(:,:,1)         )    ! surface i-current
199      IF ( iom_use("sbu") ) THEN
[4990]200         DO jj = 1, jpj
201            DO ji = 1, jpi
[5460]202               jkbot = mbku(ji,jj)
203               z2d(ji,jj) = un(ji,jj,jkbot)
[4990]204            END DO
[5107]205         END DO
206         CALL iom_put( "sbu", z2d )                ! bottom i-current
207      ENDIF
[5463]208#if defined key_dynspg_ts
209      CALL iom_put(  "ubar", un_adv(:,:)      )    ! barotropic i-current
210#else
211      CALL iom_put(  "ubar", un_b(:,:)        )    ! barotropic i-current
212#endif
[5107]213     
214      CALL iom_put( "voce", vn(:,:,:)         )    ! 3D j-current
215      CALL iom_put(  "ssv", vn(:,:,1)         )    ! surface j-current
216      IF ( iom_use("sbv") ) THEN
[4990]217         DO jj = 1, jpj
218            DO ji = 1, jpi
[5460]219               jkbot = mbkv(ji,jj)
220               z2d(ji,jj) = vn(ji,jj,jkbot)
[4990]221            END DO
[5107]222         END DO
223         CALL iom_put( "sbv", z2d )                ! bottom j-current
[4990]224      ENDIF
[5463]225#if defined key_dynspg_ts
226      CALL iom_put(  "vbar", vn_adv(:,:)      )    ! barotropic j-current
227#else
228      CALL iom_put(  "vbar", vn_b(:,:)        )    ! barotropic j-current
229#endif
[1482]230
[5461]231      CALL iom_put( "woce", wn )                   ! vertical velocity
232      IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN   ! vertical mass transport & its square value
233         ! Caution: in the VVL case, it only correponds to the baroclinic mass transport.
234         z2d(:,:) = rau0 * e12t(:,:)
235         DO jk = 1, jpk
236            z3d(:,:,jk) = wn(:,:,jk) * z2d(:,:)
237         END DO
238         CALL iom_put( "w_masstr" , z3d ) 
239         IF( iom_use('w_masstr2') )   CALL iom_put( "w_masstr2", z3d(:,:,:) * z3d(:,:,:) )
240      ENDIF
241
[5107]242      CALL iom_put( "avt" , avt                        )    ! T vert. eddy diff. coef.
243      CALL iom_put( "avm" , avmu                       )    ! T vert. eddy visc. coef.
244      CALL iom_put( "avs" , fsavs(:,:,:)               )    ! S vert. eddy diff. coef. (useful only with key_zdfddm)
245
246      IF ( iom_use("sstgrad") .OR. iom_use("sstgrad2") ) THEN
[4990]247         DO jj = 2, jpjm1                                    ! sst gradient
248            DO ji = fs_2, fs_jpim1   ! vector opt.
249               zztmp      = tsn(ji,jj,1,jp_tem)
250               zztmpx     = ( tsn(ji+1,jj  ,1,jp_tem) - zztmp ) / e1u(ji,jj) + ( zztmp - tsn(ji-1,jj  ,1,jp_tem) ) / e1u(ji-1,jj  )
251               zztmpy     = ( tsn(ji  ,jj+1,1,jp_tem) - zztmp ) / e2v(ji,jj) + ( zztmp - tsn(ji  ,jj-1,1,jp_tem) ) / e2v(ji  ,jj-1)
252               z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy )   &
253                  &              * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1)
254            END DO
[1756]255         END DO
[4990]256         CALL lbc_lnk( z2d, 'T', 1. )
257         CALL iom_put( "sstgrad2",  z2d               )    ! square of module of sst gradient
258         z2d(:,:) = SQRT( z2d(:,:) )
259         CALL iom_put( "sstgrad" ,  z2d               )    ! module of sst gradient
260      ENDIF
261         
[4761]262      ! clem: heat and salt content
[4990]263      IF( iom_use("heatc") ) THEN
264         z2d(:,:)  = 0._wp 
265         DO jk = 1, jpkm1
[5107]266            DO jj = 1, jpj
267               DO ji = 1, jpi
[4990]268                  z2d(ji,jj) = z2d(ji,jj) + fse3t(ji,jj,jk) * tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk)
269               END DO
[4761]270            END DO
271         END DO
[4990]272         CALL iom_put( "heatc", (rau0 * rcp) * z2d )    ! vertically integrated heat content (J/m2)
273      ENDIF
274
275      IF( iom_use("saltc") ) THEN
276         z2d(:,:)  = 0._wp 
277         DO jk = 1, jpkm1
[5107]278            DO jj = 1, jpj
279               DO ji = 1, jpi
[4990]280                  z2d(ji,jj) = z2d(ji,jj) + fse3t(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk)
281               END DO
282            END DO
283         END DO
284         CALL iom_put( "saltc", rau0 * z2d )   ! vertically integrated salt content (PSU*kg/m2)
285      ENDIF
[4840]286      !
[4990]287      IF ( iom_use("eken") ) THEN
288         rke(:,:,jk) = 0._wp                               !      kinetic energy
289         DO jk = 1, jpkm1
290            DO jj = 2, jpjm1
291               DO ji = fs_2, fs_jpim1   ! vector opt.
292                  zztmp   = 1._wp / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) )
293                  zztmpx  = 0.5 * (  un(ji-1,jj,jk) * un(ji-1,jj,jk) * e2u(ji-1,jj) * fse3u(ji-1,jj,jk)    &
294                     &             + un(ji  ,jj,jk) * un(ji  ,jj,jk) * e2u(ji  ,jj) * fse3u(ji  ,jj,jk) )  &
295                     &          *  zztmp 
296                  !
297                  zztmpy  = 0.5 * (  vn(ji,jj-1,jk) * vn(ji,jj-1,jk) * e1v(ji,jj-1) * fse3v(ji,jj-1,jk)    &
298                     &             + vn(ji,jj  ,jk) * vn(ji,jj  ,jk) * e1v(ji,jj  ) * fse3v(ji,jj  ,jk) )  &
299                     &          *  zztmp 
300                  !
301                  rke(ji,jj,jk) = 0.5_wp * ( zztmpx + zztmpy )
302                  !
303               ENDDO
[4840]304            ENDDO
305         ENDDO
[4990]306         CALL lbc_lnk( rke, 'T', 1. )
307         CALL iom_put( "eken", rke )           
308      ENDIF
309         
310      IF( iom_use("u_masstr") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN
[1756]311         z3d(:,:,jpk) = 0.e0
312         DO jk = 1, jpkm1
[4761]313            z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * fse3u(:,:,jk) * umask(:,:,jk)
[1756]314         END DO
315         CALL iom_put( "u_masstr", z3d )                  ! mass transport in i-direction
[4990]316      ENDIF
317     
318      IF( iom_use("u_heattr") ) THEN
319         z2d(:,:) = 0.e0 
320         DO jk = 1, jpkm1
321            DO jj = 2, jpjm1
322               DO ji = fs_2, fs_jpim1   ! vector opt.
323                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji+1,jj,jk,jp_tem) )
324               END DO
325            END DO
326         END DO
327         CALL lbc_lnk( z2d, 'U', -1. )
328         CALL iom_put( "u_heattr", (0.5 * rcp) * z2d )    ! heat transport in i-direction
329      ENDIF
[4761]330
[4990]331      IF( iom_use("u_salttr") ) THEN
[1756]332         z2d(:,:) = 0.e0 
333         DO jk = 1, jpkm1
334            DO jj = 2, jpjm1
335               DO ji = fs_2, fs_jpim1   ! vector opt.
[4990]336                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji+1,jj,jk,jp_sal) )
[1756]337               END DO
338            END DO
339         END DO
340         CALL lbc_lnk( z2d, 'U', -1. )
[4990]341         CALL iom_put( "u_salttr", 0.5 * z2d )            ! heat transport in i-direction
342      ENDIF
[4761]343
[4990]344     
345      IF( iom_use("v_masstr") .OR. iom_use("v_heattr") .OR. iom_use("v_salttr") ) THEN
[4761]346         z3d(:,:,jpk) = 0.e0
[1756]347         DO jk = 1, jpkm1
[4761]348            z3d(:,:,jk) = rau0 * vn(:,:,jk) * e1v(:,:) * fse3v(:,:,jk) * vmask(:,:,jk)
[1756]349         END DO
350         CALL iom_put( "v_masstr", z3d )                  ! mass transport in j-direction
[4990]351      ENDIF
352     
353      IF( iom_use("v_heattr") ) THEN
354         z2d(:,:) = 0.e0 
355         DO jk = 1, jpkm1
356            DO jj = 2, jpjm1
357               DO ji = fs_2, fs_jpim1   ! vector opt.
358                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji,jj+1,jk,jp_tem) )
359               END DO
360            END DO
361         END DO
362         CALL lbc_lnk( z2d, 'V', -1. )
363         CALL iom_put( "v_heattr", (0.5 * rcp) * z2d )    !  heat transport in j-direction
364      ENDIF
[4761]365
[4990]366      IF( iom_use("v_salttr") ) THEN
[1756]367         z2d(:,:) = 0.e0 
368         DO jk = 1, jpkm1
369            DO jj = 2, jpjm1
370               DO ji = fs_2, fs_jpim1   ! vector opt.
[4990]371                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji,jj+1,jk,jp_sal) )
[1756]372               END DO
373            END DO
374         END DO
375         CALL lbc_lnk( z2d, 'V', -1. )
[4990]376         CALL iom_put( "v_salttr", 0.5 * z2d )            !  heat transport in j-direction
[1756]377      ENDIF
[2528]378      !
[4990]379      CALL wrk_dealloc( jpi , jpj      , z2d )
[3294]380      CALL wrk_dealloc( jpi , jpj, jpk , z3d )
[2715]381      !
[3294]382      IF( nn_timing == 1 )   CALL timing_stop('dia_wri')
383      !
[1482]384   END SUBROUTINE dia_wri
385
386#else
[2528]387   !!----------------------------------------------------------------------
388   !!   Default option                                  use IOIPSL  library
389   !!----------------------------------------------------------------------
390
[1561]391   SUBROUTINE dia_wri( kt )
[3]392      !!---------------------------------------------------------------------
393      !!                  ***  ROUTINE dia_wri  ***
394      !!                   
395      !! ** Purpose :   Standard output of opa: dynamics and tracer fields
396      !!      NETCDF format is used by default
397      !!
398      !! ** Method  :   At the beginning of the first time step (nit000),
399      !!      define all the NETCDF files and fields
400      !!      At each time step call histdef to compute the mean if ncessary
401      !!      Each nwrite time step, output the instantaneous or mean fields
402      !!----------------------------------------------------------------------
[2715]403      !!
[3]404      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
[2528]405      !!
406      LOGICAL ::   ll_print = .FALSE.                        ! =T print and flush numout
407      CHARACTER (len=40) ::   clhstnam, clop, clmx           ! local names
408      INTEGER  ::   inum = 11                                ! temporary logical unit
[3294]409      INTEGER  ::   ji, jj, jk                               ! dummy loop indices
410      INTEGER  ::   ierr                                     ! error code return from allocation
[2528]411      INTEGER  ::   iimi, iima, ipk, it, itmod, ijmi, ijma   ! local integers
[3609]412      INTEGER  ::   jn, ierror                               ! local integers
[2528]413      REAL(wp) ::   zsto, zout, zmax, zjulian, zdt           ! local scalars
[3294]414      !!
415      REAL(wp), POINTER, DIMENSION(:,:)   :: zw2d       ! 2D workspace
416      REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d       ! 3D workspace
[3]417      !!----------------------------------------------------------------------
[3294]418      !
419      IF( nn_timing == 1 )   CALL timing_start('dia_wri')
[1482]420      !
[3294]421      CALL wrk_alloc( jpi , jpj      , zw2d )
[4292]422      IF ( ln_traldf_gdia .OR. lk_vvl )  call wrk_alloc( jpi , jpj , jpk  , zw3d )
[2715]423      !
[1482]424      ! Output the initial state and forcings
425      IF( ninist == 1 ) THEN                       
426         CALL dia_wri_state( 'output.init', kt )
427         ninist = 0
428      ENDIF
429      !
[3]430      ! 0. Initialisation
431      ! -----------------
[632]432
[3]433      ! local variable for debugging
434      ll_print = .FALSE.
435      ll_print = ll_print .AND. lwp
436
437      ! Define frequency of output and means
438      zdt = rdt
439      IF( nacc == 1 ) zdt = rdtmin
[5565]440      clop = "x"         ! no use of the mask value (require less cpu time, and otherwise the model crashes)
[3]441#if defined key_diainstant
442      zsto = nwrite * zdt
[1312]443      clop = "inst("//TRIM(clop)//")"
[3]444#else
445      zsto=zdt
[1312]446      clop = "ave("//TRIM(clop)//")"
[3]447#endif
448      zout = nwrite * zdt
449      zmax = ( nitend - nit000 + 1 ) * zdt
450
451      ! Define indices of the horizontal output zoom and vertical limit storage
452      iimi = 1      ;      iima = jpi
453      ijmi = 1      ;      ijma = jpj
454      ipk = jpk
455
456      ! define time axis
[1334]457      it = kt
458      itmod = kt - nit000 + 1
[3]459
460
461      ! 1. Define NETCDF files and fields at beginning of first time step
462      ! -----------------------------------------------------------------
463
464      IF( kt == nit000 ) THEN
465
466         ! Define the NETCDF files (one per grid)
[632]467
[3]468         ! Compute julian date from starting date of the run
[1309]469         CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
470         zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
[3]471         IF(lwp)WRITE(numout,*)
472         IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear,   &
473            &                    ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
474         IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma,   &
475                                 ' limit storage in depth = ', ipk
476
477         ! WRITE root name in date.file for use by postpro
[1581]478         IF(lwp) THEN
[895]479            CALL dia_nam( clhstnam, nwrite,' ' )
[1581]480            CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
[895]481            WRITE(inum,*) clhstnam
482            CLOSE(inum)
483         ENDIF
[632]484
[3]485         ! Define the T grid FILE ( nid_T )
[632]486
[3]487         CALL dia_nam( clhstnam, nwrite, 'grid_T' )
488         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename
489         CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit
490            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
[2528]491            &          nit000-1, zjulian, zdt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )
[3]492         CALL histvert( nid_T, "deptht", "Vertical T levels",      &  ! Vertical grid: gdept
[4292]493            &           "m", ipk, gdept_1d, nz_T, "down" )
[3]494         !                                                            ! Index of ocean points
495         CALL wheneq( jpi*jpj*ipk, tmask, 1, 1., ndex_T , ndim_T  )      ! volume
496         CALL wheneq( jpi*jpj    , tmask, 1, 1., ndex_hT, ndim_hT )      ! surface
[3609]497         !
498         IF( ln_icebergs ) THEN
499            !
500            !! allocation cant go in dia_wri_alloc because ln_icebergs is only set after
501            !! that routine is called from nemogcm, so do it here immediately before its needed
502            ALLOCATE( ndex_bT(jpi*jpj*nclasses), STAT=ierror )
503            IF( lk_mpp )   CALL mpp_sum( ierror )
504            IF( ierror /= 0 ) THEN
505               CALL ctl_stop('dia_wri: failed to allocate iceberg diagnostic array')
506               RETURN
507            ENDIF
508            !
509            !! iceberg vertical coordinate is class number
510            CALL histvert( nid_T, "class", "Iceberg class",      &  ! Vertical grid: class
511               &           "number", nclasses, class_num, nb_T )
512            !
513            !! each class just needs the surface index pattern
514            ndim_bT = 3
515            DO jn = 1,nclasses
516               ndex_bT((jn-1)*jpi*jpj+1:jn*jpi*jpj) = ndex_hT(1:jpi*jpj)
517            ENDDO
518            !
519         ENDIF
[3]520
521         ! Define the U grid FILE ( nid_U )
522
523         CALL dia_nam( clhstnam, nwrite, 'grid_U' )
524         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename
525         CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu,           &  ! Horizontal grid: glamu and gphiu
526            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
[2528]527            &          nit000-1, zjulian, zdt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )
[3]528         CALL histvert( nid_U, "depthu", "Vertical U levels",      &  ! Vertical grid: gdept
[4292]529            &           "m", ipk, gdept_1d, nz_U, "down" )
[3]530         !                                                            ! Index of ocean points
531         CALL wheneq( jpi*jpj*ipk, umask, 1, 1., ndex_U , ndim_U  )      ! volume
532         CALL wheneq( jpi*jpj    , umask, 1, 1., ndex_hU, ndim_hU )      ! surface
533
534         ! Define the V grid FILE ( nid_V )
535
536         CALL dia_nam( clhstnam, nwrite, 'grid_V' )                   ! filename
537         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
538         CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv,           &  ! Horizontal grid: glamv and gphiv
539            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
[2528]540            &          nit000-1, zjulian, zdt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )
[3]541         CALL histvert( nid_V, "depthv", "Vertical V levels",      &  ! Vertical grid : gdept
[4292]542            &          "m", ipk, gdept_1d, nz_V, "down" )
[3]543         !                                                            ! Index of ocean points
544         CALL wheneq( jpi*jpj*ipk, vmask, 1, 1., ndex_V , ndim_V  )      ! volume
545         CALL wheneq( jpi*jpj    , vmask, 1, 1., ndex_hV, ndim_hV )      ! surface
546
547         ! Define the W grid FILE ( nid_W )
548
549         CALL dia_nam( clhstnam, nwrite, 'grid_W' )                   ! filename
550         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
551         CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit
552            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
[2528]553            &          nit000-1, zjulian, zdt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set )
[3]554         CALL histvert( nid_W, "depthw", "Vertical W levels",      &  ! Vertical grid: gdepw
[4292]555            &          "m", ipk, gdepw_1d, nz_W, "down" )
[3]556
[632]557
[3]558         ! Declare all the output fields as NETCDF variables
559
560         !                                                                                      !!! nid_T : 3D
561         CALL histdef( nid_T, "votemper", "Temperature"                        , "C"      ,   &  ! tn
562            &          jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
563         CALL histdef( nid_T, "vosaline", "Salinity"                           , "PSU"    ,   &  ! sn
564            &          jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
[4292]565         IF(  lk_vvl  ) THEN
566            CALL histdef( nid_T, "vovvle3t", "Level thickness"                    , "m"      ,&  ! e3t_n
567            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
568            CALL histdef( nid_T, "vovvldep", "T point depth"                      , "m"      ,&  ! e3t_n
569            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
570            CALL histdef( nid_T, "vovvldef", "Squared level deformation"          , "%^2"    ,&  ! e3t_n
571            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
572         ENDIF
[3]573         !                                                                                      !!! nid_T : 2D
574         CALL histdef( nid_T, "sosstsst", "Sea Surface temperature"            , "C"      ,   &  ! sst
575            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
576         CALL histdef( nid_T, "sosaline", "Sea Surface Salinity"               , "PSU"    ,   &  ! sss
577            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[359]578         CALL histdef( nid_T, "sossheig", "Sea Surface Height"                 , "m"      ,   &  ! ssh
[3]579            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[2528]580         CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux"              , "Kg/m2/s",   &  ! (emp-rnf)
[3]581            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[4570]582         CALL histdef( nid_T, "sorunoff", "River runoffs"                      , "Kg/m2/s",   &  ! runoffs
583            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[3625]584         CALL histdef( nid_T, "sosfldow", "downward salt flux"                 , "PSU/m2/s",  &  ! sfx
[3]585            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[4292]586         IF(  .NOT. lk_vvl  ) THEN
587            CALL histdef( nid_T, "sosst_cd", "Concentration/Dilution term on temperature"     &  ! emp * tsn(:,:,1,jp_tem)
[3625]588            &                                                                  , "KgC/m2/s",  &  ! sosst_cd
[4292]589            &             jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
590            CALL histdef( nid_T, "sosss_cd", "Concentration/Dilution term on salinity"        &  ! emp * tsn(:,:,1,jp_sal)
[3625]591            &                                                                  , "KgPSU/m2/s",&  ! sosss_cd
[4292]592            &             jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
593         ENDIF
[888]594         CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux"             , "W/m2"   ,   &  ! qns + qsr
[3]595            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
596         CALL histdef( nid_T, "soshfldo", "Shortwave Radiation"                , "W/m2"   ,   &  ! qsr
597            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[1585]598         CALL histdef( nid_T, "somixhgt", "Turbocline Depth"                   , "m"      ,   &  ! hmld
599            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[3]600         CALL histdef( nid_T, "somxl010", "Mixed Layer Depth 0.01"             , "m"      ,   &  ! hmlp
601            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[1037]602         CALL histdef( nid_T, "soicecov", "Ice fraction"                       , "[0,1]"  ,   &  ! fr_i
[3]603            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[1649]604         CALL histdef( nid_T, "sowindsp", "wind speed at 10m"                  , "m/s"    ,   &  ! wndm
605            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[3609]606!
607         IF( ln_icebergs ) THEN
608            CALL histdef( nid_T, "calving"             , "calving mass input"                       , "kg/s"   , &
609               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
610            CALL histdef( nid_T, "calving_heat"        , "calving heat flux"                        , "XXXX"   , &
611               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
612            CALL histdef( nid_T, "berg_floating_melt"  , "Melt rate of icebergs + bits"             , "kg/m2/s", &
613               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
614            CALL histdef( nid_T, "berg_stored_ice"     , "Accumulated ice mass by class"            , "kg"     , &
615               &          jpi, jpj, nh_T, nclasses  , 1, nclasses  , nb_T , 32, clop, zsto, zout )
616            IF( ln_bergdia ) THEN
617               CALL histdef( nid_T, "berg_melt"           , "Melt rate of icebergs"                    , "kg/m2/s", &
618                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
619               CALL histdef( nid_T, "berg_buoy_melt"      , "Buoyancy component of iceberg melt rate"  , "kg/m2/s", &
620                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
621               CALL histdef( nid_T, "berg_eros_melt"      , "Erosion component of iceberg melt rate"   , "kg/m2/s", &
622                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
623               CALL histdef( nid_T, "berg_conv_melt"      , "Convective component of iceberg melt rate", "kg/m2/s", &
624                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
625               CALL histdef( nid_T, "berg_virtual_area"   , "Virtual coverage by icebergs"             , "m2"     , &
626                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
627               CALL histdef( nid_T, "bits_src"           , "Mass source of bergy bits"                , "kg/m2/s", &
628                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
629               CALL histdef( nid_T, "bits_melt"          , "Melt rate of bergy bits"                  , "kg/m2/s", &
630                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
631               CALL histdef( nid_T, "bits_mass"          , "Bergy bit density field"                  , "kg/m2"  , &
632                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
633               CALL histdef( nid_T, "berg_mass"           , "Iceberg density field"                    , "kg/m2"  , &
634                  &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
635               CALL histdef( nid_T, "berg_real_calving"   , "Calving into iceberg class"               , "kg/s"   , &
636                  &          jpi, jpj, nh_T, nclasses  , 1, nclasses  , nb_T , 32, clop, zsto, zout )
637            ENDIF
638         ENDIF
639
[5407]640         IF( .NOT. ln_cpl ) THEN
[4990]641            CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp
642               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
643            CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping"        , "Kg/m2/s",   &  ! erp
644               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
645            CALL histdef( nid_T, "sosafldp", "Surface salt flux: damping"         , "Kg/m2/s",   &  ! erp * sn
646               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
647         ENDIF
[3]648
[5407]649         IF( ln_cpl .AND. nn_ice <= 1 ) THEN
[4990]650            CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp
651               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
652            CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping"        , "Kg/m2/s",   &  ! erp
653               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
654            CALL histdef( nid_T, "sosafldp", "Surface salt flux: Damping"         , "Kg/m2/s",   &  ! erp * sn
655               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
656         ENDIF
657         
[3]658         clmx ="l_max(only(x))"    ! max index on a period
659         CALL histdef( nid_T, "sobowlin", "Bowl Index"                         , "W-point",   &  ! bowl INDEX
660            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clmx, zsto, zout )
661#if defined key_diahth
662         CALL histdef( nid_T, "sothedep", "Thermocline Depth"                  , "m"      ,   & ! hth
663            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
664         CALL histdef( nid_T, "so20chgt", "Depth of 20C isotherm"              , "m"      ,   & ! hd20
665            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
666         CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm"              , "m"      ,   & ! hd28
667            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
668         CALL histdef( nid_T, "sohtc300", "Heat content 300 m"                 , "W"      ,   & ! htc3
669            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
670#endif
671
[5407]672         IF( ln_cpl .AND. nn_ice == 2 ) THEN
[4990]673            CALL histdef( nid_T,"soicetem" , "Ice Surface Temperature"            , "K"      ,   &  ! tn_ice
674               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
675            CALL histdef( nid_T,"soicealb" , "Ice Albedo"                         , "[0,1]"  ,   &  ! alb_ice
676               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
677         ENDIF
[3]678
[2528]679         CALL histend( nid_T, snc4chunks=snc4set )
[3]680
681         !                                                                                      !!! nid_U : 3D
682         CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! un
683            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
[3294]684         IF( ln_traldf_gdia ) THEN
685            CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current"                  , "m/s"    ,   &  ! u_eiv
686                 &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
687         ELSE
[3]688#if defined key_diaeiv
[3294]689            CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current"                  , "m/s"    ,   &  ! u_eiv
[3]690            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
691#endif
[3294]692         END IF
[3]693         !                                                                                      !!! nid_U : 2D
[888]694         CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! utau
[3]695            &          jpi, jpj, nh_U, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
696
[2528]697         CALL histend( nid_U, snc4chunks=snc4set )
[3]698
699         !                                                                                      !!! nid_V : 3D
700         CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vn
701            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
[3294]702         IF( ln_traldf_gdia ) THEN
703            CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current"             , "m/s"    ,   &  ! v_eiv
704                 &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
705         ELSE 
[3]706#if defined key_diaeiv
[3294]707            CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current"             , "m/s"    ,   &  ! v_eiv
[3]708            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
709#endif
[3294]710         END IF
[3]711         !                                                                                      !!! nid_V : 2D
[888]712         CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! vtau
[3]713            &          jpi, jpj, nh_V, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
714
[2528]715         CALL histend( nid_V, snc4chunks=snc4set )
[3]716
717         !                                                                                      !!! nid_W : 3D
718         CALL histdef( nid_W, "vovecrtz", "Vertical Velocity"                  , "m/s"    ,   &  ! wn
719            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[3294]720         IF( ln_traldf_gdia ) THEN
721            CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity"              , "m/s"    ,   &  ! w_eiv
722                 &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
723         ELSE
[3]724#if defined key_diaeiv
[3294]725            CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity"              , "m/s"    ,   &  ! w_eiv
726                 &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[3]727#endif
[3294]728         END IF
[3]729         CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity"          , "m2/s"   ,   &  ! avt
730            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
[255]731         CALL histdef( nid_W, "votkeavm", "Vertical Eddy Viscosity"             , "m2/s"  ,   &  ! avmu
732            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
733
[3]734         IF( lk_zdfddm ) THEN
735            CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity"    , "m2/s"   ,   &  ! avs
736               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
737         ENDIF
738         !                                                                                      !!! nid_W : 2D
739#if defined key_traldf_c2d
740         CALL histdef( nid_W, "soleahtw", "lateral eddy diffusivity"           , "m2/s"   ,   &  ! ahtw
741            &          jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
[23]742# if defined key_traldf_eiv 
[3]743            CALL histdef( nid_W, "soleaeiw", "eddy induced vel. coeff. at w-point", "m2/s",   &  ! aeiw
744               &       jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
[23]745# endif
[3]746#endif
747
[2528]748         CALL histend( nid_W, snc4chunks=snc4set )
[3]749
750         IF(lwp) WRITE(numout,*)
751         IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
752         IF(ll_print) CALL FLUSH(numout )
753
754      ENDIF
755
756      ! 2. Start writing data
757      ! ---------------------
758
[4292]759      ! ndex(1) est utilise ssi l'avant dernier argument est different de
[3]760      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
761      ! donne le nombre d'elements, et ndex la liste des indices a sortir
762
[1334]763      IF( lwp .AND. MOD( itmod, nwrite ) == 0 ) THEN
[3]764         WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
765         WRITE(numout,*) '~~~~~~ '
766      ENDIF
767
[4292]768      IF( lk_vvl ) THEN
[4492]769         CALL histwrite( nid_T, "votemper", it, tsn(:,:,:,jp_tem) * fse3t_n(:,:,:) , ndim_T , ndex_T  )   ! heat content
770         CALL histwrite( nid_T, "vosaline", it, tsn(:,:,:,jp_sal) * fse3t_n(:,:,:) , ndim_T , ndex_T  )   ! salt content
771         CALL histwrite( nid_T, "sosstsst", it, tsn(:,:,1,jp_tem) * fse3t_n(:,:,1) , ndim_hT, ndex_hT )   ! sea surface heat content
772         CALL histwrite( nid_T, "sosaline", it, tsn(:,:,1,jp_sal) * fse3t_n(:,:,1) , ndim_hT, ndex_hT )   ! sea surface salinity content
[4292]773      ELSE
774         CALL histwrite( nid_T, "votemper", it, tsn(:,:,:,jp_tem) , ndim_T , ndex_T  )   ! temperature
775         CALL histwrite( nid_T, "vosaline", it, tsn(:,:,:,jp_sal) , ndim_T , ndex_T  )   ! salinity
776         CALL histwrite( nid_T, "sosstsst", it, tsn(:,:,1,jp_tem) , ndim_hT, ndex_hT )   ! sea surface temperature
777         CALL histwrite( nid_T, "sosaline", it, tsn(:,:,1,jp_sal) , ndim_hT, ndex_hT )   ! sea surface salinity
778      ENDIF
779      IF( lk_vvl ) THEN
780         zw3d(:,:,:) = ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2
781         CALL histwrite( nid_T, "vovvle3t", it, fse3t_n (:,:,:) , ndim_T , ndex_T  )   ! level thickness
782         CALL histwrite( nid_T, "vovvldep", it, fsdept_n(:,:,:) , ndim_T , ndex_T  )   ! t-point depth
783         CALL histwrite( nid_T, "vovvldef", it, zw3d             , ndim_T , ndex_T  )   ! level thickness deformation
784      ENDIF
[3]785      CALL histwrite( nid_T, "sossheig", it, sshn          , ndim_hT, ndex_hT )   ! sea surface height
[2528]786      CALL histwrite( nid_T, "sowaflup", it, ( emp-rnf )   , ndim_hT, ndex_hT )   ! upward water flux
[4570]787      CALL histwrite( nid_T, "sorunoff", it, rnf           , ndim_hT, ndex_hT )   ! river runoffs
[3625]788      CALL histwrite( nid_T, "sosfldow", it, sfx           , ndim_hT, ndex_hT )   ! downward salt flux
789                                                                                  ! (includes virtual salt flux beneath ice
790                                                                                  ! in linear free surface case)
[4292]791      IF( .NOT. lk_vvl ) THEN
792         zw2d(:,:) = emp (:,:) * tsn(:,:,1,jp_tem)
793         CALL histwrite( nid_T, "sosst_cd", it, zw2d, ndim_hT, ndex_hT )          ! c/d term on sst
794         zw2d(:,:) = emp (:,:) * tsn(:,:,1,jp_sal)
795         CALL histwrite( nid_T, "sosss_cd", it, zw2d, ndim_hT, ndex_hT )          ! c/d term on sss
796      ENDIF
[888]797      CALL histwrite( nid_T, "sohefldo", it, qns + qsr     , ndim_hT, ndex_hT )   ! total heat flux
[3]798      CALL histwrite( nid_T, "soshfldo", it, qsr           , ndim_hT, ndex_hT )   ! solar heat flux
[1585]799      CALL histwrite( nid_T, "somixhgt", it, hmld          , ndim_hT, ndex_hT )   ! turbocline depth
[3]800      CALL histwrite( nid_T, "somxl010", it, hmlp          , ndim_hT, ndex_hT )   ! mixed layer depth
[1037]801      CALL histwrite( nid_T, "soicecov", it, fr_i          , ndim_hT, ndex_hT )   ! ice fraction   
[1649]802      CALL histwrite( nid_T, "sowindsp", it, wndm          , ndim_hT, ndex_hT )   ! wind speed   
[3609]803!
804      IF( ln_icebergs ) THEN
805         !
806         CALL histwrite( nid_T, "calving"             , it, berg_grid%calving      , ndim_hT, ndex_hT ) 
807         CALL histwrite( nid_T, "calving_heat"        , it, berg_grid%calving_hflx , ndim_hT, ndex_hT )         
808         CALL histwrite( nid_T, "berg_floating_melt"  , it, berg_grid%floating_melt, ndim_hT, ndex_hT ) 
809         !
810         CALL histwrite( nid_T, "berg_stored_ice"     , it, berg_grid%stored_ice   , ndim_bT, ndex_bT )
811         !
812         IF( ln_bergdia ) THEN
813            CALL histwrite( nid_T, "berg_melt"           , it, berg_melt        , ndim_hT, ndex_hT   ) 
814            CALL histwrite( nid_T, "berg_buoy_melt"      , it, buoy_melt        , ndim_hT, ndex_hT   ) 
815            CALL histwrite( nid_T, "berg_eros_melt"      , it, eros_melt        , ndim_hT, ndex_hT   ) 
816            CALL histwrite( nid_T, "berg_conv_melt"      , it, conv_melt        , ndim_hT, ndex_hT   ) 
817            CALL histwrite( nid_T, "berg_virtual_area"   , it, virtual_area     , ndim_hT, ndex_hT   ) 
818            CALL histwrite( nid_T, "bits_src"            , it, bits_src         , ndim_hT, ndex_hT   ) 
819            CALL histwrite( nid_T, "bits_melt"           , it, bits_melt        , ndim_hT, ndex_hT   ) 
820            CALL histwrite( nid_T, "bits_mass"           , it, bits_mass        , ndim_hT, ndex_hT   ) 
821            CALL histwrite( nid_T, "berg_mass"           , it, berg_mass        , ndim_hT, ndex_hT   ) 
822            !
823            CALL histwrite( nid_T, "berg_real_calving"   , it, real_calving     , ndim_bT, ndex_bT   )
824         ENDIF
825      ENDIF
826
[5407]827      IF( .NOT. ln_cpl ) THEN
[4990]828         CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping
829         CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping
[3294]830         IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1)
[4990]831         CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping
832      ENDIF
[5407]833      IF( ln_cpl .AND. nn_ice <= 1 ) THEN
[4990]834         CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping
835         CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping
836         IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1)
837         CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping
838      ENDIF
839!      zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)
840!      CALL histwrite( nid_T, "sobowlin", it, zw2d          , ndim_hT, ndex_hT )   ! ???
[3]841
842#if defined key_diahth
843      CALL histwrite( nid_T, "sothedep", it, hth           , ndim_hT, ndex_hT )   ! depth of the thermocline
844      CALL histwrite( nid_T, "so20chgt", it, hd20          , ndim_hT, ndex_hT )   ! depth of the 20 isotherm
845      CALL histwrite( nid_T, "so28chgt", it, hd28          , ndim_hT, ndex_hT )   ! depth of the 28 isotherm
846      CALL histwrite( nid_T, "sohtc300", it, htc3          , ndim_hT, ndex_hT )   ! first 300m heaat content
847#endif
[888]848
[5407]849      IF( ln_cpl .AND. nn_ice == 2 ) THEN
[4990]850         CALL histwrite( nid_T, "soicetem", it, tn_ice(:,:,1) , ndim_hT, ndex_hT )   ! surf. ice temperature
851         CALL histwrite( nid_T, "soicealb", it, alb_ice(:,:,1), ndim_hT, ndex_hT )   ! ice albedo
852      ENDIF
853
[3]854      CALL histwrite( nid_U, "vozocrtx", it, un            , ndim_U , ndex_U )    ! i-current
[3294]855      IF( ln_traldf_gdia ) THEN
856         IF (.not. ALLOCATED(psix_eiv))THEN
857            ALLOCATE( psix_eiv(jpi,jpj,jpk) , psiy_eiv(jpi,jpj,jpk) , STAT=ierr )
858            IF( lk_mpp   )   CALL mpp_sum ( ierr )
859            IF( ierr > 0 )   CALL ctl_stop('STOP', 'diawri: unable to allocate psi{x,y}_eiv')
860            psix_eiv(:,:,:) = 0.0_wp
861            psiy_eiv(:,:,:) = 0.0_wp
862         ENDIF
863         DO jk=1,jpkm1
864            zw3d(:,:,jk) = (psix_eiv(:,:,jk+1) - psix_eiv(:,:,jk))/fse3u(:,:,jk)  ! u_eiv = -dpsix/dz
865         END DO
866         zw3d(:,:,jpk) = 0._wp
867         CALL histwrite( nid_U, "vozoeivu", it, zw3d, ndim_U , ndex_U )           ! i-eiv current
868      ELSE
[3]869#if defined key_diaeiv
[3294]870         CALL histwrite( nid_U, "vozoeivu", it, u_eiv, ndim_U , ndex_U )          ! i-eiv current
[3]871#endif
[3294]872      ENDIF
[888]873      CALL histwrite( nid_U, "sozotaux", it, utau          , ndim_hU, ndex_hU )   ! i-wind stress
[3]874
875      CALL histwrite( nid_V, "vomecrty", it, vn            , ndim_V , ndex_V  )   ! j-current
[3294]876      IF( ln_traldf_gdia ) THEN
877         DO jk=1,jpk-1
878            zw3d(:,:,jk) = (psiy_eiv(:,:,jk+1) - psiy_eiv(:,:,jk))/fse3v(:,:,jk)  ! v_eiv = -dpsiy/dz
879         END DO
880         zw3d(:,:,jpk) = 0._wp
881         CALL histwrite( nid_V, "vomeeivv", it, zw3d, ndim_V , ndex_V )           ! j-eiv current
882      ELSE
[3]883#if defined key_diaeiv
[3294]884         CALL histwrite( nid_V, "vomeeivv", it, v_eiv, ndim_V , ndex_V )          ! j-eiv current
[3]885#endif
[3294]886      ENDIF
[888]887      CALL histwrite( nid_V, "sometauy", it, vtau          , ndim_hV, ndex_hV )   ! j-wind stress
[3]888
889      CALL histwrite( nid_W, "vovecrtz", it, wn             , ndim_T, ndex_T )    ! vert. current
[3294]890      IF( ln_traldf_gdia ) THEN
891         DO jk=1,jpk-1
892            DO jj = 2, jpjm1
893               DO ji = fs_2, fs_jpim1  ! vector opt.
894                  zw3d(ji,jj,jk) = (psiy_eiv(ji,jj,jk) - psiy_eiv(ji,jj-1,jk))/e2v(ji,jj) + &
895                       &    (psix_eiv(ji,jj,jk) - psix_eiv(ji-1,jj,jk))/e1u(ji,jj) ! w_eiv = dpsiy/dy + dpsiy/dx
896               END DO
897            END DO
898         END DO
899         zw3d(:,:,jpk) = 0._wp
900         CALL histwrite( nid_W, "voveeivw", it, zw3d          , ndim_T, ndex_T )    ! vert. eiv current
901      ELSE
[3]902#   if defined key_diaeiv
[3294]903         CALL histwrite( nid_W, "voveeivw", it, w_eiv          , ndim_T, ndex_T )    ! vert. eiv current
[3]904#   endif
[3294]905      ENDIF
[3]906      CALL histwrite( nid_W, "votkeavt", it, avt            , ndim_T, ndex_T )    ! T vert. eddy diff. coef.
[255]907      CALL histwrite( nid_W, "votkeavm", it, avmu           , ndim_T, ndex_T )    ! T vert. eddy visc. coef.
[3]908      IF( lk_zdfddm ) THEN
909         CALL histwrite( nid_W, "voddmavs", it, fsavs(:,:,:), ndim_T, ndex_T )    ! S vert. eddy diff. coef.
910      ENDIF
911#if defined key_traldf_c2d
912      CALL histwrite( nid_W, "soleahtw", it, ahtw          , ndim_hT, ndex_hT )   ! lateral eddy diff. coef.
[23]913# if defined key_traldf_eiv
[3]914         CALL histwrite( nid_W, "soleaeiw", it, aeiw       , ndim_hT, ndex_hT )   ! EIV coefficient at w-point
[23]915# endif
[3]916#endif
917
[1318]918      ! 3. Close all files
919      ! ---------------------------------------
[1561]920      IF( kt == nitend ) THEN
[3]921         CALL histclo( nid_T )
922         CALL histclo( nid_U )
923         CALL histclo( nid_V )
924         CALL histclo( nid_W )
925      ENDIF
[2528]926      !
[3294]927      CALL wrk_dealloc( jpi , jpj      , zw2d )
[4292]928      IF ( ln_traldf_gdia .OR. lk_vvl )  call wrk_dealloc( jpi , jpj , jpk  , zw3d )
[2715]929      !
[3294]930      IF( nn_timing == 1 )   CALL timing_stop('dia_wri')
931      !
[3]932   END SUBROUTINE dia_wri
[1482]933# endif
[3]934
[1567]935#endif
936
[1334]937   SUBROUTINE dia_wri_state( cdfile_name, kt )
[3]938      !!---------------------------------------------------------------------
939      !!                 ***  ROUTINE dia_wri_state  ***
940      !!       
941      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
942      !!      the instantaneous ocean state and forcing fields.
943      !!        Used to find errors in the initial state or save the last
944      !!      ocean state in case of abnormal end of a simulation
945      !!
946      !! ** Method  :   NetCDF files using ioipsl
947      !!      File 'output.init.nc'  is created if ninist = 1 (namelist)
948      !!      File 'output.abort.nc' is created in case of abnormal job end
949      !!----------------------------------------------------------------------
[1334]950      CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created
951      INTEGER           , INTENT( in ) ::   kt               ! ocean time-step index
[2528]952      !!
[648]953      CHARACTER (len=32) :: clname
[3]954      CHARACTER (len=40) :: clop
[2528]955      INTEGER  ::   id_i , nz_i, nh_i       
956      INTEGER, DIMENSION(1) ::   idex             ! local workspace
957      REAL(wp) ::   zsto, zout, zmax, zjulian, zdt
[3]958      !!----------------------------------------------------------------------
[3294]959      !
[3625]960!     IF( nn_timing == 1 )   CALL timing_start('dia_wri_state') ! not sure this works for routines not called in first timestep
[3]961
962      ! 0. Initialisation
963      ! -----------------
[632]964
[648]965      ! Define name, frequency of output and means
966      clname = cdfile_name
[1792]967      IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
[3]968      zdt  = rdt
969      zsto = rdt
970      clop = "inst(x)"           ! no use of the mask value (require less cpu time)
971      zout = rdt
972      zmax = ( nitend - nit000 + 1 ) * zdt
973
[648]974      IF(lwp) WRITE(numout,*)
975      IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
976      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created '
977      IF(lwp) WRITE(numout,*) '                and named :', clname, '.nc'
978
979
[3]980      ! 1. Define NETCDF files and fields at beginning of first time step
981      ! -----------------------------------------------------------------
982
983      ! Compute julian date from starting date of the run
[1310]984      CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )         ! time axis
985      zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
[648]986      CALL histbeg( clname, jpi, glamt, jpj, gphit,   &
[2528]987          1, jpi, 1, jpj, nit000-1, zjulian, zdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
[3]988      CALL histvert( id_i, "deptht", "Vertical T levels",   &    ! Vertical grid : gdept
[4292]989          "m", jpk, gdept_1d, nz_i, "down")
[3]990
991      ! Declare all the output fields as NetCDF variables
992
993      CALL histdef( id_i, "vosaline", "Salinity"              , "PSU"    ,   &   ! salinity
994         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
995      CALL histdef( id_i, "votemper", "Temperature"           , "C"      ,   &   ! temperature
996         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
[359]997      CALL histdef( id_i, "sossheig", "Sea Surface Height"    , "m"      ,   &  ! ssh
998         &          jpi, jpj, nh_i, 1  , 1, 1  , nz_i, 32, clop, zsto, zout )
[3]999      CALL histdef( id_i, "vozocrtx", "Zonal Current"         , "m/s"    ,   &   ! zonal current
1000         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
1001      CALL histdef( id_i, "vomecrty", "Meridional Current"    , "m/s"    ,   &   ! meridonal current
1002         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
1003      CALL histdef( id_i, "vovecrtz", "Vertical Velocity"     , "m/s"    ,   &   ! vertical current
1004         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
1005      CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S",   &   ! net freshwater
1006         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
1007      CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2"   ,   &   ! net heat flux
1008         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
1009      CALL histdef( id_i, "soshfldo", "Shortwave Radiation"   , "W/m2"   ,   &   ! solar flux
1010         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[1037]1011      CALL histdef( id_i, "soicecov", "Ice fraction"          , "[0,1]"  ,   &   ! fr_i
[3]1012         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
1013      CALL histdef( id_i, "sozotaux", "Zonal Wind Stress"     , "N/m2"   ,   &   ! i-wind stress
1014         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
1015      CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2"   ,   &   ! j-wind stress
1016         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
[4292]1017      IF( lk_vvl ) THEN
1018         CALL histdef( id_i, "vovvldep", "T point depth"         , "m"      ,   &   ! t-point depth
1019            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
[6204]1020         CALL histdef( id_i, "vovvle3t", "T point thickness"         , "m"      ,   &   ! t-point depth
1021            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
[4292]1022      END IF
[3]1023
[1482]1024#if defined key_lim2
1025      CALL lim_wri_state_2( kt, id_i, nh_i )
[4161]1026#elif defined key_lim3
1027      CALL lim_wri_state( kt, id_i, nh_i )
[1482]1028#else
[2528]1029      CALL histend( id_i, snc4chunks=snc4set )
[1482]1030#endif
[3]1031
1032      ! 2. Start writing data
1033      ! ---------------------
1034      ! idex(1) est utilise ssi l'avant dernier argument est diffferent de
1035      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
1036      ! donne le nombre d'elements, et idex la liste des indices a sortir
1037      idex(1) = 1   ! init to avoid compil warning
[632]1038
[3]1039      ! Write all fields on T grid
[3294]1040      CALL histwrite( id_i, "votemper", kt, tsn(:,:,:,jp_tem), jpi*jpj*jpk, idex )    ! now temperature
1041      CALL histwrite( id_i, "vosaline", kt, tsn(:,:,:,jp_sal), jpi*jpj*jpk, idex )    ! now salinity
1042      CALL histwrite( id_i, "sossheig", kt, sshn             , jpi*jpj    , idex )    ! sea surface height
1043      CALL histwrite( id_i, "vozocrtx", kt, un               , jpi*jpj*jpk, idex )    ! now i-velocity
1044      CALL histwrite( id_i, "vomecrty", kt, vn               , jpi*jpj*jpk, idex )    ! now j-velocity
1045      CALL histwrite( id_i, "vovecrtz", kt, wn               , jpi*jpj*jpk, idex )    ! now k-velocity
1046      CALL histwrite( id_i, "sowaflup", kt, (emp-rnf )       , jpi*jpj    , idex )    ! freshwater budget
1047      CALL histwrite( id_i, "sohefldo", kt, qsr + qns        , jpi*jpj    , idex )    ! total heat flux
1048      CALL histwrite( id_i, "soshfldo", kt, qsr              , jpi*jpj    , idex )    ! solar heat flux
1049      CALL histwrite( id_i, "soicecov", kt, fr_i             , jpi*jpj    , idex )    ! ice fraction
1050      CALL histwrite( id_i, "sozotaux", kt, utau             , jpi*jpj    , idex )    ! i-wind stress
1051      CALL histwrite( id_i, "sometauy", kt, vtau             , jpi*jpj    , idex )    ! j-wind stress
[6204]1052      IF( lk_vvl ) THEN
1053         CALL histwrite( id_i, "vovvldep", kt, fsdept_n(:,:,:), jpi*jpj*jpk, idex )!  T-cell depth       
1054         CALL histwrite( id_i, "vovvle3t", kt, fse3t_n (:,:,:), jpi*jpj*jpk, idex )!  T-cell thickness 
1055      END IF
[3]1056
1057      ! 3. Close the file
1058      ! -----------------
1059      CALL histclo( id_i )
[1567]1060#if ! defined key_iomput && ! defined key_dimgout
[1561]1061      IF( ninist /= 1  ) THEN
1062         CALL histclo( nid_T )
1063         CALL histclo( nid_U )
1064         CALL histclo( nid_V )
1065         CALL histclo( nid_W )
1066      ENDIF
1067#endif
[3294]1068       
[3625]1069!     IF( nn_timing == 1 )   CALL timing_stop('dia_wri_state') ! not sure this works for routines not called in first timestep
[3294]1070      !
[3]1071
1072   END SUBROUTINE dia_wri_state
1073   !!======================================================================
1074END MODULE diawri
Note: See TracBrowser for help on using the repository browser.