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 @ 1312

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

add a namelist logical to mask land points in NetCDF outputs, see ticket:322

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