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

Last change on this file since 3 was 3, checked in by opalod, 20 years ago

Initial revision
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 33.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          ! ???
18   USE ocesbc          ! surface thermohaline fluxes
19   USE taumod          ! surface stress
20   USE flxrnf          ! ???
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
30   IMPLICIT NONE
31   PRIVATE
32
33   !! * Accessibility
34   PUBLIC dia_wri                 ! routines called by step.F90
35   PUBLIC dia_wri_state
36   PUBLIC dia_wri_dimg            ! called by trd_mld (eg)
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 , LODYC-IPSL  (2003)
54   !!----------------------------------------------------------------------
55
56CONTAINS
57
58#if defined key_fdir 
59   !!----------------------------------------------------------------------
60   !!   'key_fdir' :                              direct access output file
61   !!----------------------------------------------------------------------
62#   include "diawri_fdir.h90"
63
64#elif defined key_dimgout
65   !!----------------------------------------------------------------------
66   !!   dia_wri       : create the dimg direct access output file (mpp)
67   !!----------------------------------------------------------------------
68#   include "diawri_dimg.h90"
69
70#else
71   !!----------------------------------------------------------------------
72   !!   Default option                                   NetCDF output file
73   !!----------------------------------------------------------------------
74   !!   dia_wri       : create the standart NetCDF output files
75   !!   dia_wri_state : create an output NetCDF file for a single
76   !!                   instantaeous ocean state and forcing fields
77   !!----------------------------------------------------------------------
78
79   SUBROUTINE dia_wri( kt, kindic )
80      !!---------------------------------------------------------------------
81      !!                  ***  ROUTINE dia_wri  ***
82      !!                   
83      !! ** Purpose :   Standard output of opa: dynamics and tracer fields
84      !!      NETCDF format is used by default
85      !!
86      !! ** Method  :   At the beginning of the first time step (nit000),
87      !!      define all the NETCDF files and fields
88      !!      At each time step call histdef to compute the mean if ncessary
89      !!      Each nwrite time step, output the instantaneous or mean fields
90      !!      IF kindic <0, output of fields before the model interruption.
91      !!      IF kindic =0, time step loop
92      !!      IF kindic >0, output of fields before the time step loop
93      !!
94      !! History :
95      !!        !  91-03  (M.-A. Foujols)  Original code
96      !!        !  91-11  (G. Madec)
97      !!        !  92-06  (M. Imbard)  correction restart file
98      !!        !  92-07  (M. Imbard)  split into diawri and rstwri
99      !!        !  93-03  (M. Imbard)  suppress writibm
100      !!        !  98-01  (C. Levy)  NETCDF format using ioipsl INTERFACE
101      !!        !  99-02  (E. Guilyardi)  name of netCDF files + variables
102      !!   8.5  !  02-09  (G. Madec)  F90: Free form and module
103      !!----------------------------------------------------------------------
104      !! * Modules used
105      USE ioipsl
106
107      !! * Arguments
108      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
109      INTEGER, INTENT( in ) ::   kindic  !
110
111      !! * Local declarations
112      LOGICAL ::   ll_print = .FALSE.    ! =T print and flush numout
113      CHARACTER (len=40) ::           &
114         clhstnam, clop, clmx            ! temporary names
115      INTEGER ::   inum = 11             ! temporary logical unit
116      INTEGER ::   &
117         iimi, iima, ipk, it,         &  ! temporary integers
118         ijmi, ijma                      !    "          "
119      REAL(wp) ::   &
120         zsto, zout, zmax,            &  ! temporary scalars
121         zjulian, zdt                    !    "         "
122      REAL(wp), DIMENSION(jpi,jpj) :: &
123         zw2d                            ! temporary workspace
124      !!----------------------------------------------------------------------
125     
126      ! 0. Initialisation
127      ! -----------------
128     
129      ! local variable for debugging
130      ll_print = .FALSE.
131      ll_print = ll_print .AND. lwp
132
133      ! Define frequency of output and means
134      zdt = rdt
135      IF( nacc == 1 ) zdt = rdtmin
136#if defined key_diainstant
137      zsto = nwrite * zdt
138      clop = "inst(x)"           ! no use of the mask value (require less cpu time)
139      !!! clop="inst(only(x))"   ! put 1.e+20 on land (very expensive!!)
140#else
141      zsto=zdt
142      clop="ave(x)"              ! no use of the mask value (require less cpu time)
143      !!! clop="ave(only(x))"    ! put 1.e+20 on land (very expensive!!)
144#endif
145      zout = nwrite * zdt
146      zmax = ( nitend - nit000 + 1 ) * zdt
147
148      ! Define indices of the horizontal output zoom and vertical limit storage
149      iimi = 1      ;      iima = jpi
150      ijmi = 1      ;      ijma = jpj
151      ipk = jpk
152
153      ! define time axis
154      it = kt - nit000 + 1
155
156
157      ! 1. Define NETCDF files and fields at beginning of first time step
158      ! -----------------------------------------------------------------
159
160      IF(ll_print) WRITE(numout,*) 'dia_wri kt = ', kt, ' kindic ', kindic
161
162      IF( kt == nit000 ) THEN
163
164         ! Define the NETCDF files (one per grid)
165         
166         ! Compute julian date from starting date of the run
167         CALL ymds2ju( nyear, nmonth, nday, 0.e0, zjulian )
168         IF(lwp)WRITE(numout,*)
169         IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear,   &
170            &                    ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
171         IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma,   &
172                                 ' limit storage in depth = ', ipk
173
174         ! WRITE root name in date.file for use by postpro
175         CALL dia_nam( clhstnam, nwrite,' ' )
176         CALL ctlopn( inum, 'date.file', 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', 1, numout, lwp, 1 )
177         WRITE(inum,*) clhstnam
178         CLOSE(inum)
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 )
187         CALL histvert( nid_T, "deptht", "Vertical T levels",      &  ! Vertical grid: gdept
188            &           "m", ipk, gdept, 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 )
200         CALL histvert( nid_U, "depthu", "Vertical U levels",      &  ! Vertical grid: gdept
201            &           "m", ipk, gdept, 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 )
213         CALL histvert( nid_V, "depthv", "Vertical V levels",      &  ! Vertical grid : gdept
214            &          "m", ipk, gdept, 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 )
226         CALL histvert( nid_W, "depthw", "Vertical W levels",      &  ! Vertical grid: gdepw
227            &          "m", ipk, gdepw, 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_fsc
243         CALL histdef( nid_T, "sossheig", "Sea Surface Height"                 , "m"      ,   &  ! ssh
244            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
245#else
246         CALL histdef( nid_T, "sobarstf","Barotropic StreamFunction"           , "m3/s2"  ,   &  ! bsf
247            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
248#endif
249#if defined key_dynspg_fsc && defined key_ice_lim
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"   ,   &  ! qt
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 Cover"                          , "[0,1]"  ,   &  ! freeze
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#if ( defined key_coupled && ! defined key_ice_lim )
286         CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp
287            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
288         CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping"        , "Kg/m2/s",   &  ! erp
289            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
290         CALL histdef( nid_T, "sosafldp", "Surface salt flux: Damping"         , "Kg/m2/s",   &  ! erp * sn
291            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
292#endif
293#if defined key_diaspr
294         CALL histdef( nid_T, "sosurfps", "Surface Pressure"                   , "cm"     ,   &  ! sp
295            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
296#endif
297         clmx ="l_max(only(x))"    ! max index on a period
298         CALL histdef( nid_T, "sobowlin", "Bowl Index"                         , "W-point",   &  ! bowl INDEX
299            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clmx, zsto, zout )
300#if defined key_diahth
301         CALL histdef( nid_T, "sothedep", "Thermocline Depth"                  , "m"      ,   & ! hth
302            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
303         CALL histdef( nid_T, "so20chgt", "Depth of 20C isotherm"              , "m"      ,   & ! hd20
304            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
305         CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm"              , "m"      ,   & ! hd28
306            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
307         CALL histdef( nid_T, "sohtc300", "Heat content 300 m"                 , "W"      ,   & ! htc3
308            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
309#endif
310
311#if defined key_ice_lim && defined key_coupled
312         CALL histdef( nid_T,"soicetem" , "Ice Surface Temperature"            , "K"      ,   &  ! tn_ice
313            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
314         CALL histdef( nid_T,"soicealb" , "Ice Albedo"                         , "[0,1]"  ,   &  ! alb_ice
315            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
316#endif
317
318         CALL histend( nid_T )
319
320         !                                                                                      !!! nid_U : 3D
321         CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! un
322            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
323#if defined key_diaeiv
324         CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current"                  , "m/s"    ,   &  ! u_eiv
325            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
326#endif
327         !                                                                                      !!! nid_U : 2D
328         CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! taux
329            &          jpi, jpj, nh_U, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
330#if ! defined key_dynspg_fsc
331         CALL histdef( nid_U, "sozospgx", "Zonal Surface Pressure Gradient"    , "N/kg"   ,   &  ! spgu
332            &          jpi, jpj, nh_U, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
333#endif
334
335         CALL histend( nid_U )
336
337         !                                                                                      !!! nid_V : 3D
338         CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vn
339            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
340#if defined key_diaeiv
341         CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current"             , "m/s"    ,   &  ! v_eiv
342            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
343#endif
344         !                                                                                      !!! nid_V : 2D
345         CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! tauy
346            &          jpi, jpj, nh_V, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
347#if ! defined key_dynspg_fsc
348         CALL histdef( nid_V, "somespgy", "Meridional Surface Pressure Grad."  , "N/kg"   ,   &  ! spgv
349            &          jpi, jpj, nh_V, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
350#endif
351
352         CALL histend( nid_V )
353
354         !                                                                                      !!! nid_W : 3D
355         CALL histdef( nid_W, "vovecrtz", "Vertical Velocity"                  , "m/s"    ,   &  ! wn
356            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
357#if defined key_diaeiv
358         CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity"              , "m/s"    ,   &  ! w_eiv
359            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
360#endif
361         CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity"          , "m2/s"   ,   &  ! avt
362            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
363         IF( lk_zdfddm ) THEN
364            CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity"    , "m2/s"   ,   &  ! avs
365               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
366            CALL histdef( nid_W,"voddmrra","Heat/Salt buoyancy Ratio"          , "-"      ,   &  ! rrau
367               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )
368         ENDIF
369         !                                                                                      !!! nid_W : 2D
370#if defined key_traldf_c2d
371         CALL histdef( nid_W, "soleahtw", "lateral eddy diffusivity"           , "m2/s"   ,   &  ! ahtw
372            &          jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
373         IF( lk_traldf_eiv ) THEN
374            CALL histdef( nid_W, "soleaeiw", "eddy induced vel. coeff. at w-point", "m2/s",   &  ! aeiw
375               &       jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
376         ENDIF
377#endif
378
379         CALL histend( nid_W )
380
381         IF(lwp) WRITE(numout,*)
382         IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
383         IF(ll_print) CALL FLUSH(numout )
384
385      ENDIF
386
387      ! 2. Start writing data
388      ! ---------------------
389
390      ! ndex(1) est utilise ssi l'avant dernier argument est diffferent de
391      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
392      ! donne le nombre d'elements, et ndex la liste des indices a sortir
393
394      IF( lwp .AND. MOD( kt, nwrite ) == 0 ) THEN
395         WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
396         WRITE(numout,*) '~~~~~~ '
397      ENDIF
398
399      ! Write fields on T grid
400      CALL histwrite( nid_T, "votemper", it, tn            , ndim_T , ndex_T  )   ! temperature
401      CALL histwrite( nid_T, "vosaline", it, sn            , ndim_T , ndex_T  )   ! salinity
402      CALL histwrite( nid_T, "sosstsst", it, tn(:,:,1)     , ndim_hT, ndex_hT )   ! sea surface temperature
403      CALL histwrite( nid_T, "sosaline", it, sn(:,:,1)     , ndim_hT, ndex_hT )   ! sea surface salinity
404#if defined key_dynspg_fsc
405      CALL histwrite( nid_T, "sossheig", it, sshn          , ndim_hT, ndex_hT )   ! sea surface height
406#else
407      CALL histwrite( nid_T, "sobarstf", it, bsfn          , ndim_hT, ndex_hT )   ! barotropic streamfunction
408#endif
409#if defined key_dynspg_fsc && defined key_ice_lim
410      CALL histwrite( nid_T, "iowaflup", it, fsalt(:,:)    , ndim_hT, ndex_hT )   ! ice=>ocean water flux
411      CALL histwrite( nid_T, "sowaflep", it, fmass(:,:)    , ndim_hT, ndex_hT )   ! atmos=>ocean water flux
412#endif
413      CALL histwrite( nid_T, "sowaflup", it, emp           , ndim_hT, ndex_hT )   ! upward water flux
414      CALL histwrite( nid_T, "sorunoff", it, runoff        , ndim_hT, ndex_hT )   ! runoff
415      CALL histwrite( nid_T, "sowaflcd", it, emps          , ndim_hT, ndex_hT )   ! c/d water flux
416      zw2d(:,:) = emps(:,:) * sn(:,:,1) * tmask(:,:,1)
417      CALL histwrite( nid_T, "sosalflx", it, zw2d          , ndim_hT, ndex_hT )   ! c/d salt flux
418      CALL histwrite( nid_T, "sohefldo", it, qt            , ndim_hT, ndex_hT )   ! total heat flux
419      CALL histwrite( nid_T, "soshfldo", it, qsr           , ndim_hT, ndex_hT )   ! solar heat flux
420      CALL histwrite( nid_T, "somxl010", it, hmlp          , ndim_hT, ndex_hT )   ! mixed layer depth
421      CALL histwrite( nid_T, "somixhgt", it, hmld          , ndim_hT, ndex_hT )   ! turbocline depth
422      CALL histwrite( nid_T, "soicecov", it, freeze        , ndim_hT, ndex_hT )   ! ice cover
423#if ! defined key_coupled
424      CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping
425      CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping
426      zw2d(:,:) = erp(:,:) * sn(:,:,1) * tmask(:,:,1)
427      CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping
428#endif
429#if ( defined key_coupled && ! defined key_ice_lim )
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_diaspr
436      CALL histwrite( nid_T, "sosurfps", it, gps           , ndim_hT, ndex_hT )   ! surface pressure
437#endif
438         zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)
439      CALL histwrite( nid_T, "sobowlin", it, zw2d          , ndim_hT, ndex_hT )   ! ???
440
441#if defined key_diahth
442      CALL histwrite( nid_T, "sothedep", it, hth           , ndim_hT, ndex_hT )   ! depth of the thermocline
443      CALL histwrite( nid_T, "so20chgt", it, hd20          , ndim_hT, ndex_hT )   ! depth of the 20 isotherm
444      CALL histwrite( nid_T, "so28chgt", it, hd28          , ndim_hT, ndex_hT )   ! depth of the 28 isotherm
445      CALL histwrite( nid_T, "sohtc300", it, htc3          , ndim_hT, ndex_hT )   ! first 300m heaat content
446#endif
447#if defined key_ice_lim &&  defined key_coupled 
448      CALL histwrite( nid_T, "soicetem", it, tn_ice        , ndim_hT, ndex_hT )   ! surf. ice temperature
449      CALL histwrite( nid_T, "soicealb", it, alb_ice       , ndim_hT, ndex_hT )   ! ice albedo
450#endif
451         ! Write fields on U grid
452      CALL histwrite( nid_U, "vozocrtx", it, un            , ndim_U , ndex_U )    ! i-current
453#if defined key_diaeiv
454      CALL histwrite( nid_U, "vozoeivu", it, u_eiv         , ndim_U , ndex_U )    ! i-eiv current
455#endif
456      CALL histwrite( nid_U, "sozotaux", it, taux          , ndim_hU, ndex_hU )   ! i-wind stress
457#if ! defined key_dynspg_fsc
458      CALL lbc_lnk( spgu, 'U', -1. )
459      CALL histwrite( nid_U, "sozospgx", it, spgu          , ndim_hU, ndex_hU )   ! i-surf. press. grad.
460#endif
461
462         ! Write fields on V grid
463      CALL histwrite( nid_V, "vomecrty", it, vn            , ndim_V , ndex_V  )   ! j-current
464#if defined key_diaeiv
465      CALL histwrite( nid_V, "vomeeivv", it, v_eiv         , ndim_V , ndex_V  )   ! j-eiv current
466#endif
467      CALL histwrite( nid_V, "sometauy", it, tauy          , ndim_hV, ndex_hV )   ! j-wind stress
468#if ! defined key_dynspg_fsc
469      CALL lbc_lnk( spgv, 'V', -1. )
470      CALL histwrite( nid_V, "somespgy", it, spgv          , ndim_hV, ndex_hV )   ! j-surf. pressure grad.
471#endif
472
473         ! Write fields on W grid
474      CALL histwrite( nid_W, "vovecrtz", it, wn             , ndim_T, ndex_T )    ! vert. current
475#   if defined key_diaeiv
476      CALL histwrite( nid_W, "voveeivw", it, w_eiv          , ndim_T, ndex_T )    ! vert. eiv current
477#   endif
478      CALL histwrite( nid_W, "votkeavt", it, avt            , ndim_T, ndex_T )    ! T vert. eddy diff. coef.
479      IF( lk_zdfddm ) THEN
480         CALL histwrite( nid_W, "voddmavs", it, fsavs(:,:,:), ndim_T, ndex_T )    ! S vert. eddy diff. coef.
481      ENDIF
482#if defined key_traldf_c2d
483      CALL histwrite( nid_W, "soleahtw", it, ahtw          , ndim_hT, ndex_hT )   ! lateral eddy diff. coef.
484      IF( lk_traldf_eiv ) THEN
485         CALL histwrite( nid_W, "soleaeiw", it, aeiw       , ndim_hT, ndex_hT )   ! EIV coefficient at w-point
486      ENDIF
487#endif
488
489      ! 3. Synchronise and close all files
490      ! ---------------------------------------
491      IF( MOD( kt, nwrite ) == 0 .OR. kindic < 0 ) THEN
492         CALL histsync( nid_T )
493         CALL histsync( nid_U )
494         CALL histsync( nid_V )
495         CALL histsync( nid_W )
496      ENDIF
497
498      !  Create an output files (output.abort.nc) if S < 0 or u > 20 m/s
499      IF( kindic < 0 )   CALL dia_wri_state( 'output.abort' )
500
501      IF( kt == nitend .OR. kindic < 0 ) THEN
502         CALL histclo( nid_T )
503         CALL histclo( nid_U )
504         CALL histclo( nid_V )
505         CALL histclo( nid_W )
506      ENDIF
507
508   END SUBROUTINE dia_wri
509
510
511   SUBROUTINE dia_wri_state( cdfile_name )
512      !!---------------------------------------------------------------------
513      !!                 ***  ROUTINE dia_wri_state  ***
514      !!       
515      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
516      !!      the instantaneous ocean state and forcing fields.
517      !!        Used to find errors in the initial state or save the last
518      !!      ocean state in case of abnormal end of a simulation
519      !!
520      !! ** Method  :   NetCDF files using ioipsl
521      !!      File 'output.init.nc'  is created if ninist = 1 (namelist)
522      !!      File 'output.abort.nc' is created in case of abnormal job end
523      !!
524      !! History :
525      !!   8.2  !  00-06  (M. Imbard)  Original code (diabort.F)
526      !!   8.5  !  02-06  (A.Bozec, E. Durand)  Original code (diainit.F)
527      !!   9.0  !  02-12  (G. Madec)  merge of diabort and diainit, F90
528      !!----------------------------------------------------------------------
529      !! * Modules used
530      USE ioipsl
531
532      !! * Arguments
533      CHARACTER (len=* ), INTENT( in ) ::   &
534         cdfile_name      ! name of the file created
535
536      !! * Local declarations
537      CHARACTER (len=40) :: clop
538      INTEGER  ::   &
539         id_i , nz_i, nh_i       
540      INTEGER, DIMENSION(1) ::   &
541         idex             ! temprary workspace
542      REAL(wp) ::   &
543         zsto, zout, zmax,   &
544         zjulian, zdt
545      !!----------------------------------------------------------------------
546
547      IF(lwp) WRITE(numout,*)
548      IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
549      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created '
550      IF(lwp) WRITE(numout,*) '                and named :', cdfile_name, '.nc'
551     
552      ! 0. Initialisation
553      ! -----------------
554     
555      ! Define frequency of output and means
556      zdt  = rdt
557      zsto = rdt
558      clop = "inst(x)"           ! no use of the mask value (require less cpu time)
559      zout = rdt
560      zmax = ( nitend - nit000 + 1 ) * zdt
561
562      ! 1. Define NETCDF files and fields at beginning of first time step
563      ! -----------------------------------------------------------------
564
565      ! Compute julian date from starting date of the run
566      CALL ymds2ju( nyear, nmonth, nday, 0.e0, zjulian )         ! time axis
567      CALL histbeg( cdfile_name, jpi, glamt, jpj, gphit,   &
568          1, jpi, 1, jpj, 0, zjulian, zdt, nh_i, id_i )       ! Horizontal grid : glamt and gphit
569      CALL histvert( id_i, "deptht", "Vertical T levels",   &    ! Vertical grid : gdept
570          "m", jpk, gdept, nz_i)
571
572      ! Declare all the output fields as NetCDF variables
573
574      CALL histdef( id_i, "vosaline", "Salinity"              , "PSU"    ,   &   ! salinity
575         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
576      CALL histdef( id_i, "votemper", "Temperature"           , "C"      ,   &   ! temperature
577         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
578      CALL histdef( id_i, "vozocrtx", "Zonal Current"         , "m/s"    ,   &   ! zonal current
579         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
580      CALL histdef( id_i, "vomecrty", "Meridional Current"    , "m/s"    ,   &   ! meridonal current
581         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
582      CALL histdef( id_i, "vovecrtz", "Vertical Velocity"     , "m/s"    ,   &   ! vertical current
583         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
584      CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S",   &   ! net freshwater
585         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
586      CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2"   ,   &   ! net heat flux
587         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
588      CALL histdef( id_i, "soshfldo", "Shortwave Radiation"   , "W/m2"   ,   &   ! solar flux
589         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
590      CALL histdef( id_i, "soicecov", "Ice fraction"          , "[0,1]"  ,   &   ! freeze
591         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
592      CALL histdef( id_i, "sozotaux", "Zonal Wind Stress"     , "N/m2"   ,   &   ! i-wind stress
593         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
594      CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2"   ,   &   ! j-wind stress
595         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
596
597      CALL histend( id_i )
598
599      ! 2. Start writing data
600      ! ---------------------
601      ! idex(1) est utilise ssi l'avant dernier argument est diffferent de
602      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
603      ! donne le nombre d'elements, et idex la liste des indices a sortir
604      idex(1) = 1   ! init to avoid compil warning
605     
606      ! Write all fields on T grid
607      CALL histwrite( id_i, "votemper", 1, tn    , jpi*jpj*jpk, idex )    ! now temperature
608      CALL histwrite( id_i, "vosaline", 1, sn    , jpi*jpj*jpk, idex )    ! now salinity
609      CALL histwrite( id_i, "vozocrtx", 1, un    , jpi*jpj*jpk, idex )    ! now i-velocity
610      CALL histwrite( id_i, "vomecrty", 1, vn    , jpi*jpj*jpk, idex )    ! now j-velocity
611      CALL histwrite( id_i, "vovecrtz", 1, wn    , jpi*jpj*jpk, idex )    ! now k-velocity
612      CALL histwrite( id_i, "sowaflup", 1, emp   , jpi*jpj    , idex )    ! freshwater budget
613      CALL histwrite( id_i, "sohefldo", 1, qt    , jpi*jpj    , idex )    ! total heat flux
614      CALL histwrite( id_i, "soshfldo", 1, qsr   , jpi*jpj    , idex )    ! total heat flux
615      CALL histwrite( id_i, "soicecov", 1, freeze, jpi*jpj    , idex )    ! ice cover
616      CALL histwrite( id_i, "sozotaux", 1, taux  , jpi*jpj    , idex )    ! i-wind stress
617      CALL histwrite( id_i, "sometauy", 1, tauy  , jpi*jpj    , idex )    ! j-wind stress
618
619      ! 3. Close the file
620      ! -----------------
621      CALL histclo( id_i )
622
623   END SUBROUTINE dia_wri_state
624 
625   SUBROUTINE dia_wri_dimg( cd_name, cd_exper, ptab, klev, cd_type )
626   REAL(wp), DIMENSION(:,:,:) :: ptab
627   INTEGER :: klev
628   CHARACTER(LEN=80) :: cd_name, cd_exper,cd_type
629   print *, ' This print must never occur ', cd_name, cd_exper, ptab, klev, cd_type
630   PRINT *, ' this routine is here just for compilation '
631   
632   END SUBROUTINE dia_wri_dimg
633#endif
634   !!======================================================================
635END MODULE diawri
Note: See TracBrowser for help on using the repository browser.