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

Last change on this file since 359 was 359, checked in by opalod, 18 years ago

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