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 trunk/NEMO/OPA_SRC/DIA – NEMO

source: trunk/NEMO/OPA_SRC/DIA/diawri.F90 @ 1397

Last change on this file since 1397 was 1359, checked in by smasson, 15 years ago

first implementation of iom_put, see ticket:387

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 39.4 KB
Line 
1MODULE diawri
2   !!======================================================================
3   !!                     ***  MODULE  diawri  ***
4   !! Ocean diagnostics :  write ocean output files
5   !!=====================================================================
6
7   !!----------------------------------------------------------------------
8   !! * Modules used
9   USE oce             ! ocean dynamics and tracers
10   USE dom_oce         ! ocean space and time domain
11   USE zdf_oce         ! ocean vertical physics
12   USE ldftra_oce      ! ocean active tracers: lateral physics
13   USE ldfdyn_oce      ! ocean dynamics: lateral physics
14   USE sol_oce         ! solver variables
15   USE ice_oce         ! ice variables
16   USE sbc_oce         ! Surface boundary condition: ocean fields
17   USE sbc_ice         ! Surface boundary condition: ice fields
18   USE sbcssr          ! restoring term toward SST/SSS climatology
19   USE phycst          ! physical constants
20   USE zdfmxl          ! mixed layer
21   USE daymod          ! calendar
22   USE dianam          ! build name of file (routine)
23   USE zdfddm          ! vertical  physics: double diffusion
24   USE diahth          ! thermocline diagnostics
25   USE diaspr          ! surface pressure diagnostics (rigid lid case)
26   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
27   USE in_out_manager  ! I/O manager
28   USE diadimg         ! dimg direct access file format output
29   USE iom
30   USE ioipsl
31
32   IMPLICIT NONE
33   PRIVATE
34
35   !! * Accessibility
36   PUBLIC dia_wri                 ! routines called by step.F90
37   PUBLIC dia_wri_state
38
39   !! * Module variables
40   INTEGER ::   &
41      nid_T, nz_T, nh_T, ndim_T, ndim_hT,      &   ! grid_T file
42      nid_U, nz_U, nh_U, ndim_U, ndim_hU,      &   ! grid_U file
43      nid_V, nz_V, nh_V, ndim_V, ndim_hV,      &   ! grid_V file
44      nid_W, nz_W, nh_W,                       &   ! grid_W file
45      ndex(1)                                      ! ???
46   INTEGER, DIMENSION(jpi*jpj) ::   &
47      ndex_hT, ndex_hU, ndex_hV
48   INTEGER, DIMENSION(jpi*jpj*jpk) ::   &
49      ndex_T, ndex_U, ndex_V
50
51   !! * Substitutions
52#  include "zdfddm_substitute.h90"
53   !!----------------------------------------------------------------------
54   !!   OPA 9.0 , LOCEAN-IPSL (2005)
55   !! $Id$
56   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
57   !!----------------------------------------------------------------------
58
59CONTAINS
60
61#if defined key_dimgout
62   !!----------------------------------------------------------------------
63   !!   dia_wri       : create the dimg direct access output file (mpp)
64   !!----------------------------------------------------------------------
65#   include "diawri_dimg.h90"
66
67#else
68   !!----------------------------------------------------------------------
69   !!   Default option                                   NetCDF output file
70   !!----------------------------------------------------------------------
71   !!   dia_wri       : create the standart NetCDF output files
72   !!   dia_wri_state : create an output NetCDF file for a single
73   !!                   instantaeous ocean state and forcing fields
74   !!----------------------------------------------------------------------
75
76   SUBROUTINE dia_wri( kt, kindic )
77      !!---------------------------------------------------------------------
78      !!                  ***  ROUTINE dia_wri  ***
79      !!                   
80      !! ** Purpose :   Standard output of opa: dynamics and tracer fields
81      !!      NETCDF format is used by default
82      !!
83      !! ** Method  :   At the beginning of the first time step (nit000),
84      !!      define all the NETCDF files and fields
85      !!      At each time step call histdef to compute the mean if ncessary
86      !!      Each nwrite time step, output the instantaneous or mean fields
87      !!      IF kindic <0, output of fields before the model interruption.
88      !!      IF kindic =0, time step loop
89      !!      IF kindic >0, output of fields before the time step loop
90      !!
91      !! History :
92      !!        !  91-03  (M.-A. Foujols)  Original code
93      !!        !  91-11  (G. Madec)
94      !!        !  92-06  (M. Imbard)  correction restart file
95      !!        !  92-07  (M. Imbard)  split into diawri and rstwri
96      !!        !  93-03  (M. Imbard)  suppress writibm
97      !!        !  98-01  (C. Levy)  NETCDF format using ioipsl INTERFACE
98      !!        !  99-02  (E. Guilyardi)  name of netCDF files + variables
99      !!   8.5  !  02-09  (G. Madec)  F90: Free form and module
100      !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization
101      !!----------------------------------------------------------------------
102      !! * Arguments
103      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
104      INTEGER, INTENT( in ) ::   kindic  !
105
106      !! * Local declarations
107      LOGICAL ::   ll_print = .FALSE.    ! =T print and flush numout
108      CHARACTER (len=40) ::           &
109         clhstnam, clop, clmx            ! temporary names
110      INTEGER ::   inum = 11             ! temporary logical unit
111      INTEGER ::   &
112         iimi, iima, ipk, it, itmod,  &  ! temporary integers
113         ijmi, ijma                      !    "          "
114      REAL(wp) ::   &
115         zsto, zout, zmax,            &  ! temporary scalars
116         zjulian, zdt                    !    "         "
117      REAL(wp), DIMENSION(jpi,jpj) :: &
118         zw2d                            ! temporary workspace
119      CHARACTER (len=80) :: clname
120      !!----------------------------------------------------------------------
121
122      ! 0. Initialisation
123      ! -----------------
124
125      ! local variable for debugging
126      ll_print = .FALSE.
127      ll_print = ll_print .AND. lwp
128
129      ! Define frequency of output and means
130      zdt = rdt
131      IF( nacc == 1 ) zdt = rdtmin
132      IF( ln_mskland )   THEN   ;   clop = "only(x)"   ! put 1.e+20 on land (very expensive!!)
133      ELSE                      ;   clop = "x"         ! no use of the mask value (require less cpu time)
134      ENDIF
135#if defined key_diainstant
136      zsto = nwrite * zdt
137      clop = "inst("//TRIM(clop)//")"
138#else
139      zsto=zdt
140      clop = "ave("//TRIM(clop)//")"
141#endif
142      zout = nwrite * zdt
143      zmax = ( nitend - nit000 + 1 ) * zdt
144
145      ! Define indices of the horizontal output zoom and vertical limit storage
146      iimi = 1      ;      iima = jpi
147      ijmi = 1      ;      ijma = jpj
148      ipk = jpk
149
150      ! define time axis
151      it = kt
152      itmod = kt - nit000 + 1
153
154
155      ! 1. Define NETCDF files and fields at beginning of first time step
156      ! -----------------------------------------------------------------
157
158      IF(ll_print) WRITE(numout,*) 'dia_wri kt = ', kt, ' kindic ', kindic
159
160      IF( kt == nit000 ) THEN
161
162         ! Define the NETCDF files (one per grid)
163
164         ! Compute julian date from starting date of the run
165         CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
166         zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
167         IF(lwp)WRITE(numout,*)
168         IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear,   &
169            &                    ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
170         IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma,   &
171                                 ' limit storage in depth = ', ipk
172
173         ! WRITE root name in date.file for use by postpro
174         IF( lwp) THEN
175            CALL dia_nam( clhstnam, nwrite,' ' )
176            clname = 'date.file'
177            CALL ctlopn( inum, clname,  'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', 1, numout, lwp, 1 )
178            WRITE(inum,*) clhstnam
179            CLOSE(inum)
180         ENDIF
181
182         ! Define the T grid FILE ( nid_T )
183
184         CALL dia_nam( clhstnam, nwrite, 'grid_T' )
185         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename
186         CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit
187            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
188            &          nit000-1, zjulian, zdt, nh_T, nid_T, domain_id=nidom )
189         CALL histvert( nid_T, "deptht", "Vertical T levels",      &  ! Vertical grid: gdept
190            &           "m", ipk, gdept_0, nz_T, "down" )
191         !                                                            ! Index of ocean points
192         CALL wheneq( jpi*jpj*ipk, tmask, 1, 1., ndex_T , ndim_T  )      ! volume
193         CALL wheneq( jpi*jpj    , tmask, 1, 1., ndex_hT, ndim_hT )      ! surface
194
195         ! Define the U grid FILE ( nid_U )
196
197         CALL dia_nam( clhstnam, nwrite, 'grid_U' )
198         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename
199         CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu,           &  ! Horizontal grid: glamu and gphiu
200            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
201            &          nit000-1, zjulian, zdt, nh_U, nid_U, domain_id=nidom )
202         CALL histvert( nid_U, "depthu", "Vertical U levels",      &  ! Vertical grid: gdept
203            &           "m", ipk, gdept_0, nz_U, "down" )
204         !                                                            ! Index of ocean points
205         CALL wheneq( jpi*jpj*ipk, umask, 1, 1., ndex_U , ndim_U  )      ! volume
206         CALL wheneq( jpi*jpj    , umask, 1, 1., ndex_hU, ndim_hU )      ! surface
207
208         ! Define the V grid FILE ( nid_V )
209
210         CALL dia_nam( clhstnam, nwrite, 'grid_V' )                   ! filename
211         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
212         CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv,           &  ! Horizontal grid: glamv and gphiv
213            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
214            &          nit000-1, zjulian, zdt, nh_V, nid_V, domain_id=nidom )
215         CALL histvert( nid_V, "depthv", "Vertical V levels",      &  ! Vertical grid : gdept
216            &          "m", ipk, gdept_0, nz_V, "down" )
217         !                                                            ! Index of ocean points
218         CALL wheneq( jpi*jpj*ipk, vmask, 1, 1., ndex_V , ndim_V  )      ! volume
219         CALL wheneq( jpi*jpj    , vmask, 1, 1., ndex_hV, ndim_hV )      ! surface
220
221         ! Define the W grid FILE ( nid_W )
222
223         CALL dia_nam( clhstnam, nwrite, 'grid_W' )                   ! filename
224         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
225         CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit
226            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
227            &          nit000-1, zjulian, zdt, nh_W, nid_W, domain_id=nidom )
228         CALL histvert( nid_W, "depthw", "Vertical W levels",      &  ! Vertical grid: gdepw
229            &          "m", ipk, gdepw_0, nz_W, "down" )
230
231
232         ! Declare all the output fields as NETCDF variables
233
234         !                                                                                      !!! nid_T : 3D
235         CALL histdef( nid_T, "votemper", "Temperature"                        , "C"      ,   &  ! tn
236            &          jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
237         CALL histdef( nid_T, "vosaline", "Salinity"                           , "PSU"    ,   &  ! sn
238            &          jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout )
239         !                                                                                      !!! nid_T : 2D
240         CALL histdef( nid_T, "sosstsst", "Sea Surface temperature"            , "C"      ,   &  ! sst
241            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
242         CALL histdef( nid_T, "sosaline", "Sea Surface Salinity"               , "PSU"    ,   &  ! sss
243            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
244#if defined key_dynspg_rl
245         CALL histdef( nid_T, "sobarstf","Barotropic StreamFunction"           , "m3/s2"  ,   &  ! bsf
246            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
247#else
248         CALL histdef( nid_T, "sossheig", "Sea Surface Height"                 , "m"      ,   &  ! ssh
249            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
250#endif
251!!$#if ! defined key_dynspg_rl && ( defined key_lim3 || defined key_lim2 )
252!!$         ! sowaflup = sowaflep + sorunoff + sowafldp + a term associated to
253!!$         !    internal damping to Levitus that can be diagnosed from others
254!!$         ! sowaflcd = sowaflep + sorunoff + sowafldp + iowaflup
255!!$         CALL histdef( nid_T, "iowaflup", "Ice=>ocean net freshwater"          , "kg/m2/s",   &  ! fsalt
256!!$            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
257!!$         CALL histdef( nid_T, "sowaflep", "atmos=>ocean net freshwater"        , "kg/m2/s",   &  ! fmass
258!!$            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
259!!$#endif
260         CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux"              , "Kg/m2/s",   &  ! emp
261            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
262!!$         CALL histdef( nid_T, "sorunoff", "Runoffs"                            , "Kg/m2/s",   &  ! runoffs
263!!$            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
264         CALL histdef( nid_T, "sowaflcd", "concentration/dilution water flux"  , "kg/m2/s",   &  ! emps
265            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
266         CALL histdef( nid_T, "sosalflx", "Surface Salt Flux"                  , "Kg/m2/s",   &  ! emps * sn
267            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
268         CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux"             , "W/m2"   ,   &  ! qns + qsr
269            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
270         CALL histdef( nid_T, "soshfldo", "Shortwave Radiation"                , "W/m2"   ,   &  ! qsr
271            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
272         CALL histdef( nid_T, "somxl010", "Mixed Layer Depth 0.01"             , "m"      ,   &  ! hmlp
273            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
274         CALL histdef( nid_T, "somixhgt", "Turbocline Depth"                   , "m"      ,   &  ! hmld
275            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
276         CALL histdef( nid_T, "soicecov", "Ice fraction"                       , "[0,1]"  ,   &  ! fr_i
277            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
278#if ! defined key_coupled
279         CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp
280            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
281         CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping"        , "Kg/m2/s",   &  ! erp
282            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
283         CALL histdef( nid_T, "sosafldp", "Surface salt flux: damping"         , "Kg/m2/s",   &  ! erp * sn
284            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
285#endif
286
287
288
289#if ( defined key_coupled && ! defined key_lim3 && ! defined key_lim2 )
290         CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp
291            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
292         CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping"        , "Kg/m2/s",   &  ! erp
293            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
294         CALL histdef( nid_T, "sosafldp", "Surface salt flux: Damping"         , "Kg/m2/s",   &  ! erp * sn
295            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
296#endif
297#if defined key_diaspr
298         CALL histdef( nid_T, "sosurfps", "Surface Pressure"                   , "cm"     ,   &  ! sp
299            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
300#endif
301         clmx ="l_max(only(x))"    ! max index on a period
302         CALL histdef( nid_T, "sobowlin", "Bowl Index"                         , "W-point",   &  ! bowl INDEX
303            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clmx, zsto, zout )
304#if defined key_diahth
305         CALL histdef( nid_T, "sothedep", "Thermocline Depth"                  , "m"      ,   & ! hth
306            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
307         CALL histdef( nid_T, "so20chgt", "Depth of 20C isotherm"              , "m"      ,   & ! hd20
308            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
309         CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm"              , "m"      ,   & ! hd28
310            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
311         CALL histdef( nid_T, "sohtc300", "Heat content 300 m"                 , "W"      ,   & ! htc3
312            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
313#endif
314
315#if defined key_coupled 
316# if defined key_lim3
317         Must be adapted to LIM3
318# else
319         CALL histdef( nid_T,"soicetem" , "Ice Surface Temperature"            , "K"      ,   &  ! tn_ice
320            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
321         CALL histdef( nid_T,"soicealb" , "Ice Albedo"                         , "[0,1]"  ,   &  ! alb_ice
322            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
323# endif 
324#endif
325
326         CALL histend( nid_T )
327
328         !                                                                                      !!! nid_U : 3D
329         CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! un
330            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
331#if defined key_diaeiv
332         CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current"                  , "m/s"    ,   &  ! u_eiv
333            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
334#endif
335         !                                                                                      !!! nid_U : 2D
336         CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! utau
337            &          jpi, jpj, nh_U, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
338#if defined key_dynspg_rl
339         CALL histdef( nid_U, "sozospgx", "Zonal Surface Pressure Gradient"    , "N/kg"   ,   &  ! spgu
340            &          jpi, jpj, nh_U, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
341#endif
342
343         CALL histend( nid_U )
344
345         !                                                                                      !!! nid_V : 3D
346         CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vn
347            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
348#if defined key_diaeiv
349         CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current"             , "m/s"    ,   &  ! v_eiv
350            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
351#endif
352         !                                                                                      !!! nid_V : 2D
353         CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! vtau
354            &          jpi, jpj, nh_V, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
355#if defined key_dynspg_rl
356         CALL histdef( nid_V, "somespgy", "Meridional Surface Pressure Grad."  , "N/kg"   ,   &  ! spgv
357            &          jpi, jpj, nh_V, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
358#endif
359
360         CALL histend( nid_V )
361
362         !                                                                                      !!! nid_W : 3D
363         CALL histdef( nid_W, "vovecrtz", "Vertical Velocity"                  , "m/s"    ,   &  ! wn
364            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
365#if defined key_diaeiv
366         CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity"              , "m/s"    ,   &  ! w_eiv
367            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
368#endif
369         CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity"          , "m2/s"   ,   &  ! avt
370            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
371
372         CALL histdef( nid_W, "votkeevd", "Enhanced Vertical Diffusivity",       "m2/s"   ,   &  ! avt_evd
373            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
374         !
375         CALL histdef( nid_W, "votkeavm", "Vertical Eddy Viscosity"             , "m2/s"  ,   &  ! avmu
376            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
377
378         CALL histdef( nid_W, "votkeevm", "Enhanced Vertical Viscosity",         "m2/s"   ,   &  ! avmu_evd
379            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
380
381         IF( lk_zdfddm ) THEN
382            CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity"    , "m2/s"   ,   &  ! avs
383               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
384         ENDIF
385         !                                                                                      !!! nid_W : 2D
386#if defined key_traldf_c2d
387         CALL histdef( nid_W, "soleahtw", "lateral eddy diffusivity"           , "m2/s"   ,   &  ! ahtw
388            &          jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
389# if defined key_traldf_eiv 
390            CALL histdef( nid_W, "soleaeiw", "eddy induced vel. coeff. at w-point", "m2/s",   &  ! aeiw
391               &       jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
392# endif
393#endif
394
395         CALL histend( nid_W )
396
397         IF(lwp) WRITE(numout,*)
398         IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
399         IF(ll_print) CALL FLUSH(numout )
400
401      ENDIF
402
403      ! 2. Start writing data
404      ! ---------------------
405
406      ! ndex(1) est utilise ssi l'avant dernier argument est diffferent de
407      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
408      ! donne le nombre d'elements, et ndex la liste des indices a sortir
409
410      IF( lwp .AND. MOD( itmod, nwrite ) == 0 ) THEN
411         WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
412         WRITE(numout,*) '~~~~~~ '
413      ENDIF
414     
415!---------------------------------------------------------------------------
416
417      CALL iom_put("votemper",tn)
418      CALL iom_put("vosaline",sn)
419      CALL iom_put("sosstsst",tn(:,:,1))   ! sea surface temperature
420      CALL iom_put("sosaline",sn(:,:,1))   ! sea surface salinity
421#if defined key_dynspg_rl
422      CALL iom_put("sobarstf",bsfn)   ! barotropic streamfunction
423#else
424      CALL iom_put("sossheig",sshn)   ! sea surface height
425#endif
426
427      CALL iom_put("sowaflup",emp )   ! upward water flux
428      CALL iom_put("sowaflcd",emps)   ! c/d water flux
429      zw2d(:,:) = emps(:,:) * sn(:,:,1) * tmask(:,:,1)
430      CALL iom_put("sosalflx",zw2d)   ! c/d salt flux
431      CALL iom_put("sohefldo",qns + qsr )   ! total heat flux
432      CALL iom_put("soshfldo",qsr)   ! solar heat flux
433      CALL iom_put("somxl010",hmlp)   ! mixed layer depth
434      CALL iom_put("somixhgt",hmld)   ! turbocline depth
435      CALL iom_put("soicecov",fr_i)   ! ice fraction
436#if ! defined key_coupled
437      CALL iom_put("sohefldp",qrp)   ! heat flux damping
438      CALL iom_put("sowafldp",erp)   ! freshwater flux damping
439      zw2d(:,:) = erp(:,:) * sn(:,:,1) * tmask(:,:,1)
440      CALL iom_put("sosafldp", zw2d)   ! salt flux damping
441#endif
442
443#if ( defined key_coupled && ! defined key_lim3 && ! defined key_lim2 )
444      CALL iom_put("sohefldp",qrp)   ! heat flux damping
445      CALL iom_put("sowafldp",erp)   ! freshwater flux damping
446         zw2d(:,:) = erp(:,:) * sn(:,:,1) * tmask(:,:,1)
447      CALL iom_put("sosafldp",zw2d)   ! salt flux damping
448#endif
449#if defined key_diaspr
450      CALL iom_put("sosurfps",gps)   ! surface pressure
451#endif
452         zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)
453      CALL iom_put("sobowlin",zw2d)   ! ???
454
455#if defined key_diahth
456      CALL iom_put("sothedep",hth)   ! depth of the thermocline
457      CALL iom_put("so20chgt",hd20)   ! depth of the 20 isotherm
458      CALL iom_put("so28chgt",hd28)   ! depth of the 28 isotherm
459      CALL iom_put("sohtc300",htc3)   ! first 300m heaat content
460#endif
461
462#if defined key_coupled 
463#  if defined key_lim3
464      Must be adapted for LIM3
465#  else
466      CALL iom_put("soicetem",tn_ice)   ! surf. ice temperature
467      CALL iom_put("soicealb",alb_ice)   ! ice albedo
468#  endif
469#endif
470         ! Write fields on U grid
471      CALL iom_put("vozocrtx",un)    ! i-current
472#if defined key_diaeiv
473      CALL iom_put("vozoeivu",u_eiv)    ! i-eiv current
474#endif
475      CALL iom_put("sozotaux",utau)   ! i-wind stress
476#if defined key_dynspg_rl
477      CALL lbc_lnk( spgu, 'U', -1. )
478      CALL iom_put("sozospgx",spgu)   ! i-surf. press. grad.
479#endif
480
481         ! Write fields on V grid
482      CALL iom_put("vomecrty",vn)   ! j-current
483#if defined key_diaeiv
484      CALL iom_put("vomeeivv",v_eiv)   ! j-eiv current
485#endif
486      CALL iom_put("sometauy", vtau)   ! j-wind stress
487#if defined key_dynspg_rl
488      CALL lbc_lnk( spgv, 'V', -1. )
489      CALL iom_put("somespgy",spgv)   ! j-surf. pressure grad.
490#endif
491
492         ! Write fields on W grid
493      CALL iom_put("vovecrtz",wn)    ! vert. current
494#   if defined key_diaeiv
495      CALL iom_put("voveeivw",w_eiv)    ! vert. eiv current
496#   endif
497      CALL iom_put("votkeavt",avt)    ! T vert. eddy diff. coef.
498      CALL iom_put("votkeevd",avt_evd)    ! T enhan. vert. eddy diff. coef.
499      CALL iom_put("votkeavm",avmu)    ! T vert. eddy visc. coef.
500      CALL iom_put("votkeevm",avmu_evd)    ! T enhan. vert. eddy visc. coef.
501      IF( lk_zdfddm ) THEN
502         CALL iom_put("voddmavs",fsavs(:,:,:) )    ! S vert. eddy diff. coef.
503      ENDIF
504#if defined key_traldf_c2d
505      CALL iom_put("soleahtw",ahtw)   ! lateral eddy diff. coef.
506# if defined key_traldf_eiv
507      CALL iom_put("soleaeiw",aeiw)   ! EIV coefficient at w-point
508# endif
509#endif
510!---------------------------------------------------------------------------
511
512      ! Write fields on T grid
513      CALL histwrite( nid_T, "votemper", it, tn            , ndim_T , ndex_T  )   ! temperature
514      CALL histwrite( nid_T, "vosaline", it, sn            , ndim_T , ndex_T  )   ! salinity
515      CALL histwrite( nid_T, "sosstsst", it, tn(:,:,1)     , ndim_hT, ndex_hT )   ! sea surface temperature
516      CALL histwrite( nid_T, "sosaline", it, sn(:,:,1)     , ndim_hT, ndex_hT )   ! sea surface salinity
517#if defined key_dynspg_rl
518      CALL histwrite( nid_T, "sobarstf", it, bsfn          , ndim_hT, ndex_hT )   ! barotropic streamfunction
519#else
520      CALL histwrite( nid_T, "sossheig", it, sshn          , ndim_hT, ndex_hT )   ! sea surface height
521#endif
522!!$#if ! defined key_dynspg_rl && ( defined key_lim3 || defined key_lim2 )
523!!$      CALL histwrite( nid_T, "iowaflup", it, fsalt(:,:)    , ndim_hT, ndex_hT )   ! ice=>ocean water flux
524!!$      CALL histwrite( nid_T, "sowaflep", it, fmass(:,:)    , ndim_hT, ndex_hT )   ! atmos=>ocean water flux
525!!$#endif
526      CALL histwrite( nid_T, "sowaflup", it, emp           , ndim_hT, ndex_hT )   ! upward water flux
527!!$      CALL histwrite( nid_T, "sorunoff", it, runoff        , ndim_hT, ndex_hT )   ! runoff
528      CALL histwrite( nid_T, "sowaflcd", it, emps          , ndim_hT, ndex_hT )   ! c/d water flux
529      zw2d(:,:) = emps(:,:) * sn(:,:,1) * tmask(:,:,1)
530      CALL histwrite( nid_T, "sosalflx", it, zw2d          , ndim_hT, ndex_hT )   ! c/d salt flux
531      CALL histwrite( nid_T, "sohefldo", it, qns + qsr     , ndim_hT, ndex_hT )   ! total heat flux
532      CALL histwrite( nid_T, "soshfldo", it, qsr           , ndim_hT, ndex_hT )   ! solar heat flux
533      CALL histwrite( nid_T, "somxl010", it, hmlp          , ndim_hT, ndex_hT )   ! mixed layer depth
534      CALL histwrite( nid_T, "somixhgt", it, hmld          , ndim_hT, ndex_hT )   ! turbocline depth
535      CALL histwrite( nid_T, "soicecov", it, fr_i          , ndim_hT, ndex_hT )   ! ice fraction   
536#if ! defined key_coupled
537      CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping
538      CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping
539      zw2d(:,:) = erp(:,:) * sn(:,:,1) * tmask(:,:,1)
540      CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping
541#endif
542#if ( defined key_coupled && ! defined key_lim3 && ! defined key_lim2 )
543      CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping
544      CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping
545         zw2d(:,:) = erp(:,:) * sn(:,:,1) * tmask(:,:,1)
546      CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping
547#endif
548#if defined key_diaspr
549      CALL histwrite( nid_T, "sosurfps", it, gps           , ndim_hT, ndex_hT )   ! surface pressure
550#endif
551         zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)
552      CALL histwrite( nid_T, "sobowlin", it, zw2d          , ndim_hT, ndex_hT )   ! ???
553
554#if defined key_diahth
555      CALL histwrite( nid_T, "sothedep", it, hth           , ndim_hT, ndex_hT )   ! depth of the thermocline
556      CALL histwrite( nid_T, "so20chgt", it, hd20          , ndim_hT, ndex_hT )   ! depth of the 20 isotherm
557      CALL histwrite( nid_T, "so28chgt", it, hd28          , ndim_hT, ndex_hT )   ! depth of the 28 isotherm
558      CALL histwrite( nid_T, "sohtc300", it, htc3          , ndim_hT, ndex_hT )   ! first 300m heaat content
559#endif
560
561#if defined key_coupled 
562# if defined key_lim3
563      Must be adapted for LIM3
564      CALL histwrite( nid_T, "soicetem", it, tn_ice        , ndim_hT, ndex_hT )   ! surf. ice temperature
565      CALL histwrite( nid_T, "soicealb", it, alb_ice       , ndim_hT, ndex_hT )   ! ice albedo
566# else
567      CALL histwrite( nid_T, "soicetem", it, tn_ice        , ndim_hT, ndex_hT )   ! surf. ice temperature
568      CALL histwrite( nid_T, "soicealb", it, alb_ice       , ndim_hT, ndex_hT )   ! ice albedo
569# endif
570#endif
571         ! Write fields on U grid
572      CALL histwrite( nid_U, "vozocrtx", it, un            , ndim_U , ndex_U )    ! i-current
573#if defined key_diaeiv
574      CALL histwrite( nid_U, "vozoeivu", it, u_eiv         , ndim_U , ndex_U )    ! i-eiv current
575#endif
576      CALL histwrite( nid_U, "sozotaux", it, utau          , ndim_hU, ndex_hU )   ! i-wind stress
577#if defined key_dynspg_rl
578      CALL lbc_lnk( spgu, 'U', -1. )
579      CALL histwrite( nid_U, "sozospgx", it, spgu          , ndim_hU, ndex_hU )   ! i-surf. press. grad.
580#endif
581
582         ! Write fields on V grid
583      CALL histwrite( nid_V, "vomecrty", it, vn            , ndim_V , ndex_V  )   ! j-current
584#if defined key_diaeiv
585      CALL histwrite( nid_V, "vomeeivv", it, v_eiv         , ndim_V , ndex_V  )   ! j-eiv current
586#endif
587      CALL histwrite( nid_V, "sometauy", it, vtau          , ndim_hV, ndex_hV )   ! j-wind stress
588#if defined key_dynspg_rl
589      CALL lbc_lnk( spgv, 'V', -1. )
590      CALL histwrite( nid_V, "somespgy", it, spgv          , ndim_hV, ndex_hV )   ! j-surf. pressure grad.
591#endif
592
593         ! Write fields on W grid
594      CALL histwrite( nid_W, "vovecrtz", it, wn             , ndim_T, ndex_T )    ! vert. current
595#   if defined key_diaeiv
596      CALL histwrite( nid_W, "voveeivw", it, w_eiv          , ndim_T, ndex_T )    ! vert. eiv current
597#   endif
598      CALL histwrite( nid_W, "votkeavt", it, avt            , ndim_T, ndex_T )    ! T vert. eddy diff. coef.
599      CALL histwrite( nid_W, "votkeevd", it, avt_evd        , ndim_T, ndex_T )    ! T enhan. vert. eddy diff. coef.
600      CALL histwrite( nid_W, "votkeavm", it, avmu           , ndim_T, ndex_T )    ! T vert. eddy visc. coef.
601      CALL histwrite( nid_W, "votkeevm", it, avmu_evd       , ndim_T, ndex_T )    ! T enhan. vert. eddy visc. coef.
602      IF( lk_zdfddm ) THEN
603         CALL histwrite( nid_W, "voddmavs", it, fsavs(:,:,:), ndim_T, ndex_T )    ! S vert. eddy diff. coef.
604      ENDIF
605#if defined key_traldf_c2d
606      CALL histwrite( nid_W, "soleahtw", it, ahtw          , ndim_hT, ndex_hT )   ! lateral eddy diff. coef.
607# if defined key_traldf_eiv
608         CALL histwrite( nid_W, "soleaeiw", it, aeiw       , ndim_hT, ndex_hT )   ! EIV coefficient at w-point
609# endif
610#endif
611
612      !  Create an output files (output.abort.nc) if S < 0 or u > 20 m/s
613      IF( kindic < 0 )   CALL dia_wri_state( 'output.abort', kt )
614
615      ! 3. Close all files
616      ! ---------------------------------------
617      IF( kt == nitend .OR. kindic < 0 ) THEN
618         CALL histclo( nid_T )
619         CALL histclo( nid_U )
620         CALL histclo( nid_V )
621         CALL histclo( nid_W )
622      ENDIF
623
624   END SUBROUTINE dia_wri
625
626
627   SUBROUTINE dia_wri_state( cdfile_name, kt )
628      !!---------------------------------------------------------------------
629      !!                 ***  ROUTINE dia_wri_state  ***
630      !!       
631      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
632      !!      the instantaneous ocean state and forcing fields.
633      !!        Used to find errors in the initial state or save the last
634      !!      ocean state in case of abnormal end of a simulation
635      !!
636      !! ** Method  :   NetCDF files using ioipsl
637      !!      File 'output.init.nc'  is created if ninist = 1 (namelist)
638      !!      File 'output.abort.nc' is created in case of abnormal job end
639      !!
640      !! History :
641      !!   8.2  !  00-06  (M. Imbard)  Original code (diabort.F)
642      !!   8.5  !  02-06  (A.Bozec, E. Durand)  Original code (diainit.F)
643      !!   9.0  !  02-12  (G. Madec)  merge of diabort and diainit, F90
644      !!    "   !  05-11  (V. Garnier) Surface pressure gradient organization
645      !!----------------------------------------------------------------------
646      !! * Arguments
647      CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created
648      INTEGER           , INTENT( in ) ::   kt               ! ocean time-step index
649
650      !! * Local declarations
651      CHARACTER (len=32) :: clname
652      CHARACTER (len=40) :: clop
653      INTEGER  ::   &
654         id_i , nz_i, nh_i       
655      INTEGER, DIMENSION(1) ::   &
656         idex             ! temprary workspace
657      REAL(wp) ::   &
658         zsto, zout, zmax,   &
659         zjulian, zdt
660      !!----------------------------------------------------------------------
661
662      ! 0. Initialisation
663      ! -----------------
664
665      ! Define name, frequency of output and means
666      clname = cdfile_name
667#if defined key_agrif
668      if ( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
669#endif
670      zdt  = rdt
671      zsto = rdt
672      clop = "inst(x)"           ! no use of the mask value (require less cpu time)
673      zout = rdt
674      zmax = ( nitend - nit000 + 1 ) * zdt
675
676      IF(lwp) WRITE(numout,*)
677      IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
678      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created '
679      IF(lwp) WRITE(numout,*) '                and named :', clname, '.nc'
680
681
682      ! 1. Define NETCDF files and fields at beginning of first time step
683      ! -----------------------------------------------------------------
684
685      ! Compute julian date from starting date of the run
686      CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )         ! time axis
687      zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
688      CALL histbeg( clname, jpi, glamt, jpj, gphit,   &
689          1, jpi, 1, jpj, nit000-1, zjulian, zdt, nh_i, id_i, domain_id=nidom )          ! Horizontal grid : glamt and gphit
690      CALL histvert( id_i, "deptht", "Vertical T levels",   &    ! Vertical grid : gdept
691          "m", jpk, gdept_0, nz_i, "down")
692
693      ! Declare all the output fields as NetCDF variables
694
695      CALL histdef( id_i, "vosaline", "Salinity"              , "PSU"    ,   &   ! salinity
696         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
697      CALL histdef( id_i, "votemper", "Temperature"           , "C"      ,   &   ! temperature
698         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
699#if defined key_dynspg_rl
700      CALL histdef( id_i, "sobarstf","Barotropic StreamFunction", "m3/s2"  ,   &  ! bsf
701         &          jpi, jpj, nh_i, 1  , 1, 1  , nz_i, 32, clop, zsto, zout )
702#else
703      CALL histdef( id_i, "sossheig", "Sea Surface Height"    , "m"      ,   &  ! ssh
704         &          jpi, jpj, nh_i, 1  , 1, 1  , nz_i, 32, clop, zsto, zout )
705#endif
706      CALL histdef( id_i, "vozocrtx", "Zonal Current"         , "m/s"    ,   &   ! zonal current
707         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
708      CALL histdef( id_i, "vomecrty", "Meridional Current"    , "m/s"    ,   &   ! meridonal current
709         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
710      CALL histdef( id_i, "vovecrtz", "Vertical Velocity"     , "m/s"    ,   &   ! vertical current
711         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
712      CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S",   &   ! net freshwater
713         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
714      CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2"   ,   &   ! net heat flux
715         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
716      CALL histdef( id_i, "soshfldo", "Shortwave Radiation"   , "W/m2"   ,   &   ! solar flux
717         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
718      CALL histdef( id_i, "soicecov", "Ice fraction"          , "[0,1]"  ,   &   ! fr_i
719         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
720      CALL histdef( id_i, "sozotaux", "Zonal Wind Stress"     , "N/m2"   ,   &   ! i-wind stress
721         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
722      CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2"   ,   &   ! j-wind stress
723         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
724
725      CALL histend( id_i )
726
727      ! 2. Start writing data
728      ! ---------------------
729      ! idex(1) est utilise ssi l'avant dernier argument est diffferent de
730      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
731      ! donne le nombre d'elements, et idex la liste des indices a sortir
732      idex(1) = 1   ! init to avoid compil warning
733
734      ! Write all fields on T grid
735      CALL histwrite( id_i, "votemper", kt, tn      , jpi*jpj*jpk, idex )    ! now temperature
736      CALL histwrite( id_i, "vosaline", kt, sn      , jpi*jpj*jpk, idex )    ! now salinity
737#if defined key_dynspg_rl
738      CALL histwrite( id_i, "sobarstf", kt, bsfn     , jpi*jpj    , idex )    ! barotropic streamfunction
739#else
740      CALL histwrite( id_i, "sossheig", kt, sshn     , jpi*jpj    , idex )    ! sea surface height
741#endif
742      CALL histwrite( id_i, "vozocrtx", kt, un       , jpi*jpj*jpk, idex )    ! now i-velocity
743      CALL histwrite( id_i, "vomecrty", kt, vn       , jpi*jpj*jpk, idex )    ! now j-velocity
744      CALL histwrite( id_i, "vovecrtz", kt, wn       , jpi*jpj*jpk, idex )    ! now k-velocity
745      CALL histwrite( id_i, "sowaflup", kt, emp      , jpi*jpj    , idex )    ! freshwater budget
746      CALL histwrite( id_i, "sohefldo", kt, qsr + qns, jpi*jpj    , idex )    ! total heat flux
747      CALL histwrite( id_i, "soshfldo", kt, qsr      , jpi*jpj    , idex )    ! solar heat flux
748      CALL histwrite( id_i, "soicecov", kt, fr_i     , jpi*jpj    , idex )    ! ice fraction
749      CALL histwrite( id_i, "sozotaux", kt, utau     , jpi*jpj    , idex )    ! i-wind stress
750      CALL histwrite( id_i, "sometauy", kt, vtau     , jpi*jpj    , idex )    ! j-wind stress
751
752      ! 3. Close the file
753      ! -----------------
754      CALL histclo( id_i )
755
756   END SUBROUTINE dia_wri_state
757
758#endif
759   !!======================================================================
760END MODULE diawri
Note: See TracBrowser for help on using the repository browser.