ROMS
Loading...
Searching...
No Matches
def_diags.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef DIAGNOSTICS
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This module creates output TIME-AVERAGED DIAGNOSTIC file using !
13! either the standard NetCDF library or the Parallel-IO (PIO) !
14! library. It defines its dimensions, attributes, and variables. !
15! !
16!=======================================================================
17!
18 USE mod_param
19 USE mod_parallel
20# ifdef BIOLOGY
21 USE mod_biology
22# endif
23# ifdef FOUR_DVAR
24 USE mod_fourdvar
25# endif
26 USE mod_iounits
27 USE mod_ncparam
28 USE mod_scalars
29# ifdef SEDIMENT
30 USE mod_sediment
31# endif
32!
33 USE def_dim_mod, ONLY : def_dim
34 USE def_info_mod, ONLY : def_info
35 USE def_var_mod, ONLY : def_var
36 USE strings_mod, ONLY : founderror
37 USE wrt_info_mod, ONLY : wrt_info
38!
39 implicit none
40!
41 PUBLIC :: def_diags
42 PRIVATE :: def_diags_nf90
43# if defined PIO_LIB && defined DISTRIBUTE
44 PRIVATE :: def_diags_pio
45# endif
46!
47 CONTAINS
48!
49!***********************************************************************
50 SUBROUTINE def_diags (ng, ldef)
51!***********************************************************************
52!
53! Imported variable declarations.
54!
55 logical, intent(in) :: ldef
56!
57 integer, intent(in) :: ng
58!
59! Local variable declarations.
60!
61 character (len=*), parameter :: myfile = &
62 & __FILE__
63!
64!-----------------------------------------------------------------------
65! Create a new history file according to IO type.
66!-----------------------------------------------------------------------
67!
68 SELECT CASE (dia(ng)%IOtype)
69 CASE (io_nf90)
70 CALL def_diags_nf90 (ng, ldef)
71
72# if defined PIO_LIB && defined DISTRIBUTE
73 CASE (io_pio)
74 CALL def_diags_pio (ng, ldef)
75# endif
76 CASE DEFAULT
77 IF (master) WRITE (stdout,10) dia(ng)%IOtype
78 exit_flag=3
79 END SELECT
80 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
81!
82 10 FORMAT (' DEF_DIAGS - Illegal output file type, io_type = ',i0, &
83 & /,13x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
84!
85 RETURN
86 END SUBROUTINE def_diags
87!
88!***********************************************************************
89 SUBROUTINE def_diags_nf90 (ng, ldef)
90!***********************************************************************
91!
92 USE mod_netcdf
93!
94! Imported variable declarations.
95!
96 integer, intent(in) :: ng
97
98 logical, intent(in) :: ldef
99!
100! Local variable declarations.
101!
102 logical :: got_var(nv)
103!
104 integer, parameter :: natt = 25
105
106 integer :: i, ifield, itrc, ivar, j, nvd3, nvd4, nvd5
107 integer :: recdim, status
108# if defined WRITE_WATER && defined MASKING
109 integer :: xy_pdim, xyz_pdim
110# endif
111 integer :: dimids(ndimid)
112 integer :: t2dgrd(3), u2dgrd(3), v2dgrd(3)
113# if defined ECOSIM && defined DIAGNOSTICS_BIO
114 integer :: l3dgrd(4), l4dgrd(5)
115# endif
116
117# ifdef SOLVE3D
118# ifdef SEDIMENT
119 integer :: b3dgrd(4)
120# endif
121 integer :: t3dgrd(4), u3dgrd(4), v3dgrd(4), w3dgrd(4)
122# endif
123!
124 real(r8) :: aval(6)
125!
126 character (len= 13) :: prefix
127 character (len=256) :: ncname
128 character (len=MaxLen) :: vinfo(natt)
129
130 character (len=*), parameter :: myfile = &
131 & __FILE__//", def_diags_nf90"
132!
133 sourcefile=myfile
134!
135!-----------------------------------------------------------------------
136! Set and report file name.
137!-----------------------------------------------------------------------
138!
139 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
140 ncname=dia(ng)%name
141!
142 IF (master) THEN
143 IF (ldef) THEN
144 WRITE (stdout,10) ng, trim(ncname)
145 ELSE
146 WRITE (stdout,20) ng, trim(ncname)
147 END IF
148 END IF
149!
150!=======================================================================
151! Create a new diagnostics NetCDF file.
152!=======================================================================
153!
154 define : IF (ldef) THEN
155 CALL netcdf_create (ng, inlm, trim(ncname), dia(ng)%ncid)
156 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
157 IF (master) WRITE (stdout,30) trim(ncname)
158 RETURN
159 END IF
160!
161!-----------------------------------------------------------------------
162! Define file dimensions.
163!-----------------------------------------------------------------------
164!
165 dimids=0
166!
167 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xi_rho', &
168 & iobounds(ng)%xi_rho, dimids( 1))
169 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
170
171 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xi_u', &
172 & iobounds(ng)%xi_u, dimids( 2))
173 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
174
175 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xi_v', &
176 & iobounds(ng)%xi_v, dimids( 3))
177 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
178
179 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xi_psi', &
180 & iobounds(ng)%xi_psi, dimids( 4))
181 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
182
183 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'eta_rho', &
184 & iobounds(ng)%eta_rho, dimids( 5))
185 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
186
187 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'eta_u', &
188 & iobounds(ng)%eta_u, dimids( 6))
189 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
190
191 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'eta_v', &
192 & iobounds(ng)%eta_v, dimids( 7))
193 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
194
195 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'eta_psi', &
196 & iobounds(ng)%eta_psi, dimids( 8))
197 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
198
199# if defined WRITE_WATER && defined MASKING
200 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xy_psi', &
201 & iobounds(ng)%xy_psi, xy_pdim)
202 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
203
204 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xy_rho', &
205 & iobounds(ng)%xy_rho, dimids(17))
206 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
207
208 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xy_u', &
209 & iobounds(ng)%xy_u, dimids(18))
210 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
211
212 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xy_v', &
213 & iobounds(ng)%xy_v, dimids(19))
214 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
215# endif
216
217# ifdef SOLVE3D
218# if defined WRITE_WATER && defined MASKING
219 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xyz_psi', &
220 & iobounds(ng)%xy_psi*n(ng), xyz_pdim)
221 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
222
223 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xyz_rho', &
224 & iobounds(ng)%xy_rho*n(ng), dimids(20))
225 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
226
227 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xyz_u', &
228 & iobounds(ng)%xy_u*n(ng), dimids(21))
229 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
230
231 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xyz_v', &
232 & iobounds(ng)%xy_v*n(ng), dimids(22))
233 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
234
235 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xyz_w', &
236 & iobounds(ng)%xy_rho*(n(ng)+1), dimids(23))
237 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
238# endif
239
240 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 's_rho', &
241 & n(ng), dimids( 9))
242 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
243
244 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 's_w', &
245 & n(ng)+1, dimids(10))
246 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
247
248 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'tracer', &
249 & nt(ng), dimids(11))
250 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
251
252# ifdef SEDIMENT
253 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'NST', &
254 & nst, dimids(32))
255 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
256
257 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'Nbed', &
258 & nbed, dimids(16))
259 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
260
261# if defined WRITE_WATER && defined MASKING
262 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'xybed', &
263 & iobounds(ng)%xy_rho*nbed, dimids(24))
264 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
265# endif
266# endif
267# ifdef ECOSIM
268 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'Nbands', &
269 & ndbands, dimids(33))
270 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
271
272 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'Nphy', &
273 & nphy, dimids(25))
274 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
275
276 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'Nbac', &
277 & nbac, dimids(26))
278 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
279
280 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'Ndom', &
281 & ndom, dimids(27))
282 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
283
284 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'Nfec', &
285 & nfec, dimids(28))
286 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
287# endif
288# endif
289
290 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'boundary', &
291 & 4, dimids(14))
292 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
293
294# ifdef FOUR_DVAR
295 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, 'Nstate', &
296 & nstatevar(ng), dimids(29))
297 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
298# endif
299
300 status=def_dim(ng, inlm, dia(ng)%ncid, ncname, &
301 & trim(adjustl(vname(5,idtime))), &
302 & nf90_unlimited, dimids(12))
303 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
304
305 recdim=dimids(12)
306!
307! Set number of dimensions for output variables.
308!
309#if defined WRITE_WATER && defined MASKING
310 nvd3=2
311 nvd4=2
312 nvd5=2
313#else
314 nvd3=3
315 nvd4=4
316 nvd5=5
317#endif
318!
319! Define dimension vectors for staggered tracer type variables.
320!
321# if defined WRITE_WATER && defined MASKING
322 t2dgrd(1)=dimids(17)
323 t2dgrd(2)=dimids(12)
324# ifdef SOLVE3D
325 t3dgrd(1)=dimids(20)
326 t3dgrd(2)=dimids(12)
327# endif
328# else
329 t2dgrd(1)=dimids( 1)
330 t2dgrd(2)=dimids( 5)
331 t2dgrd(3)=dimids(12)
332# ifdef SOLVE3D
333 t3dgrd(1)=dimids( 1)
334 t3dgrd(2)=dimids( 5)
335 t3dgrd(3)=dimids( 9)
336 t3dgrd(4)=dimids(12)
337# endif
338# endif
339!
340! Define dimension vectors for staggered u-momentum type variables.
341!
342# if defined WRITE_WATER && defined MASKING
343 u2dgrd(1)=dimids(18)
344 u2dgrd(2)=dimids(12)
345# ifdef SOLVE3D
346 u3dgrd(1)=dimids(21)
347 u3dgrd(2)=dimids(12)
348# endif
349# else
350 u2dgrd(1)=dimids( 2)
351 u2dgrd(2)=dimids( 6)
352 u2dgrd(3)=dimids(12)
353# ifdef SOLVE3D
354 u3dgrd(1)=dimids( 2)
355 u3dgrd(2)=dimids( 6)
356 u3dgrd(3)=dimids( 9)
357 u3dgrd(4)=dimids(12)
358# endif
359# endif
360!
361! Define dimension vectors for staggered v-momentum type variables.
362!
363# if defined WRITE_WATER && defined MASKING
364 v2dgrd(1)=dimids(19)
365 v2dgrd(2)=dimids(12)
366# ifdef SOLVE3D
367 v3dgrd(1)=dimids(22)
368 v3dgrd(2)=dimids(12)
369# endif
370# else
371 v2dgrd(1)=dimids( 3)
372 v2dgrd(2)=dimids( 7)
373 v2dgrd(3)=dimids(12)
374# ifdef SOLVE3D
375 v3dgrd(1)=dimids( 3)
376 v3dgrd(2)=dimids( 7)
377 v3dgrd(3)=dimids( 9)
378 v3dgrd(4)=dimids(12)
379# endif
380# endif
381# ifdef SOLVE3D
382!
383! Define dimension vector for staggered w-momentum type variables.
384!
385# if defined WRITE_WATER && defined MASKING
386 w3dgrd(1)=dimids(23)
387 w3dgrd(2)=dimids(12)
388# else
389 w3dgrd(1)=dimids( 1)
390 w3dgrd(2)=dimids( 5)
391 w3dgrd(3)=dimids(10)
392 w3dgrd(4)=dimids(12)
393# endif
394# ifdef SEDIMENT
395!
396! Define dimension vector for sediment bed layer type variables.
397!
398# if defined WRITE_WATER && defined MASKING
399 b3dgrd(1)=dimids(24)
400 b3dgrd(2)=dimids(12)
401# else
402 b3dgrd(1)=dimids( 1)
403 b3dgrd(2)=dimids( 5)
404 b3dgrd(3)=dimids(16)
405 b3dgrd(4)=dimids(12)
406# endif
407# endif
408# if defined ECOSIM && defined DIAGNOSTICS_BIO
409!
410! Define dimension vector for spectral light type variables.
411!
412 l3dgrd(1)=dimids( 1)
413 l3dgrd(2)=dimids( 5)
414 l3dgrd(3)=dimids(33)
415 l3dgrd(4)=dimids(12)
416!
417 l4dgrd(1)=dimids( 1)
418 l4dgrd(2)=dimids( 5)
419 l4dgrd(3)=dimids( 9)
420 l4dgrd(4)=dimids(33)
421 l4dgrd(5)=dimids(12)
422# endif
423# endif
424!
425! Initialize unlimited time record dimension.
426!
427 dia(ng)%Rindex=0
428!
429! Initialize local information variable arrays.
430!
431 DO i=1,natt
432 DO j=1,len(vinfo(1))
433 vinfo(i)(j:j)=' '
434 END DO
435 END DO
436 DO i=1,6
437 aval(i)=0.0_r8
438 END DO
439!
440! Set long name prefix string.
441!
442!! Prefix='time-averaged'
443 prefix=char(32) ! blank space
444!
445!-----------------------------------------------------------------------
446! Define time-recordless information variables.
447!-----------------------------------------------------------------------
448!
449 CALL def_info (ng, inlm, dia(ng)%ncid, ncname, dimids)
450 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
451!
452!-----------------------------------------------------------------------
453! Define variables and their attributes.
454!-----------------------------------------------------------------------
455!
456! Define model time.
457!
458 vinfo( 1)=vname(1,idtime)
459 WRITE (vinfo( 2),'(a,a)') 'averaged ', trim(vname(2,idtime))
460 WRITE (vinfo( 3),'(a,a)') 'seconds since ', trim(rclock%string)
461 vinfo( 4)=trim(rclock%calendar)
462 vinfo(14)=vname(4,idtime)
463 vinfo(21)=vname(6,idtime)
464 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(idtime), &
465 & nf_tout, 1, (/recdim/), aval, vinfo, ncname, &
466 & setparaccess = .true.)
467 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
468!
469! Define time-averaged free-surface. It is needed for CF compliance for
470! vertical coordinates formula.
471!
472 vinfo( 1)=vname(1,idfsur)
473 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idfsur))
474 vinfo( 3)=vname(3,idfsur)
475 vinfo(14)=vname(4,idfsur)
476 vinfo(16)=vname(1,idtime)
477# if !defined WET_DRY && (defined WRITE_WATER && defined MASKING)
478 vinfo(20)='mask_rho'
479# endif
480 vinfo(21)=vname(6,idfsur)
481 vinfo(22)='coordinates'
482 aval(5)=real(iinfo(1,idfsur,ng),r8)
483 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(idfsur), &
484 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
485 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
486
487# ifdef DIAGNOSTICS_UV
488!
489! Define 2D momentum diagnostic fields.
490!
491 DO ivar=1,ndm2d
492 ifield=iddu2d(ivar)
493 IF (dout(ifield,ng)) THEN
494 vinfo( 1)=vname(1,ifield)
495 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
496 vinfo( 3)=vname(3,ifield)
497 vinfo(14)=vname(4,ifield)
498 vinfo(16)=vname(1,idtime)
499# if defined WRITE_WATER && defined MASKING
500 vinfo(20)='mask_u'
501# endif
502 vinfo(21)=vname(6,ifield)
503 vinfo(22)='coordinates'
504 aval(5)=real(u2dvar,r8)
505 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
506 & nf_fout, nvd3, u2dgrd, aval, vinfo, ncname)
507 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
508 END IF
509 ifield=iddv2d(ivar)
510 IF (dout(ifield,ng)) THEN
511 vinfo( 1)=vname(1,ifield)
512 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
513 vinfo( 3)=vname(3,ifield)
514 vinfo(14)=vname(4,ifield)
515 vinfo(16)=vname(1,idtime)
516# if defined WRITE_WATER && defined MASKING
517 vinfo(20)='mask_v'
518# endif
519 vinfo(21)=vname(6,ifield)
520 vinfo(22)='coordinates'
521 aval(5)=real(v2dvar,r8)
522 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
523 & nf_fout, nvd3, v2dgrd, aval, vinfo, ncname)
524 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
525 END IF
526 END DO
527
528# ifdef SOLVE3D
529!
530! Define 3D momentum diagnostic fields.
531!
532 DO ivar=1,ndm3d
533 ifield=iddu3d(ivar)
534 IF (dout(ifield,ng)) THEN
535 vinfo( 1)=vname(1,ifield)
536 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
537 vinfo( 3)=vname(3,ifield)
538 vinfo(14)=vname(4,ifield)
539 vinfo(16)=vname(1,idtime)
540# if defined WRITE_WATER && defined MASKING
541 vinfo(20)='mask_u'
542# endif
543 vinfo(21)=vname(6,ifield)
544 vinfo(22)='coordinates'
545 aval(5)=real(u3dvar,r8)
546 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
547 & nf_fout, nvd4, u3dgrd, aval, vinfo, ncname)
548 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
549 END IF
550 ifield=iddv3d(ivar)
551 IF (dout(ifield,ng)) THEN
552 vinfo( 1)=vname(1,ifield)
553 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
554 vinfo( 3)=vname(3,ifield)
555 vinfo(14)=vname(4,ifield)
556 vinfo(16)=vname(1,idtime)
557# if defined WRITE_WATER && defined MASKING
558 vinfo(20)='mask_v'
559# endif
560 vinfo(21)=vname(6,ifield)
561 vinfo(22)='coordinates'
562 aval(5)=real(v3dvar,r8)
563 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
564 & nf_fout, nvd4, v3dgrd, aval, vinfo, ncname)
565 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
566 END IF
567 END DO
568# endif
569# endif
570# ifdef DIAGNOSTICS_TS
571!
572! Define tracer diagnostic fields.
573!
574 DO itrc=1,nt(ng)
575 DO ivar=1,ndt
576 ifield=iddtrc(itrc,ivar)
577 IF (dout(ifield,ng)) THEN
578 vinfo( 1)=vname(1,ifield)
579 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
580 vinfo( 3)=vname(3,ifield)
581 vinfo(14)=vname(4,ifield)
582 vinfo(16)=vname(1,idtime)
583# if defined WRITE_WATER && defined MASKING
584 vinfo(20)='mask_rho'
585# endif
586 vinfo(21)=vname(6,ifield)
587 vinfo(22)='coordinates'
588 aval(5)=real(r3dvar,r8)
589 status=def_var(ng, inlm, dia(ng)%ncid, &
590 & dia(ng)%Vid(ifield), nf_fout, &
591 & nvd4, t3dgrd, aval, vinfo, ncname)
593 & __line__, myfile)) RETURN
594 END IF
595 END DO
596 END DO
597# endif
598
599# ifdef DIAGNOSTICS_BIO
600# if defined BIO_FENNEL || defined HYPOXIA_SRM
601!
602! Define 2D biological diagnostic fields.
603!
604 DO ivar=1,ndbio2d
605 ifield=idbio2(ivar)
606 IF (dout(ifield,ng)) THEN
607 vinfo( 1)=vname(1,ifield)
608 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
609 vinfo( 3)=vname(3,ifield)
610 vinfo(14)=vname(4,ifield)
611 vinfo(16)=vname(1,idtime)
612# if defined WRITE_WATER && defined MASKING
613 vinfo(20)='mask_rho'
614# endif
615 vinfo(21)=vname(6,ifield)
616 vinfo(22)='coordinates'
617 aval(5)=real(r2dvar,r8)
618 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
619 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
620 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
621 END IF
622 END DO
623# endif
624# if defined BIO_FENNEL
625!
626! Define 3D biological diagnostic fields.
627!
628 DO ivar=1,ndbio3d
629 ifield=idbio3(ivar)
630 IF (dout(ifield,ng)) THEN
631 vinfo( 1)=vname(1,ifield)
632 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
633 vinfo( 3)=vname(3,ifield)
634 vinfo(14)=vname(4,ifield)
635 vinfo(16)=vname(1,idtime)
636# if defined WRITE_WATER && defined MASKING
637 vinfo(20)='mask_rho'
638# endif
639 vinfo(21)=vname(6,ifield)
640 vinfo(22)='coordinates'
641 aval(5)=real(r3dvar,r8)
642 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
643 & nf_fout, nvd4, t3dgrd, aval, vinfo, ncname)
644 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
645 END IF
646 END DO
647# elif defined ECOSIM
648!
649! Define 3D biological diagnostic fields.
650!
651 DO ivar=1,ndbio3d
652 ifield=idbio3(ivar)
653 IF (dout(ifield,ng)) THEN
654 vinfo( 1)=vname(1,ifield)
655 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
656 vinfo( 3)=vname(3,ifield)
657 vinfo(14)=vname(4,ifield)
658 vinfo(16)=vname(1,idtime)
659# if defined WRITE_WATER && defined MASKING
660 vinfo(20)='mask_rho'
661# endif
662 vinfo(21)=vname(6,ifield)
663 vinfo(22)='coordinates'
664 aval(5)=real(iinfo(1,ifield,ng),r8)
665 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
666 & nf_fout, nvd4, l3dgrd, aval, vinfo, ncname)
667 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
668 END IF
669 END DO
670!
671! Define 4D biological diagnostic fields.
672!
673 DO ivar=1,ndbio4d
674 ifield=idbio4(ivar)
675 IF (dout(ifield,ng)) THEN
676 vinfo( 1)=vname(1,ifield)
677 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
678 vinfo( 3)=vname(3,ifield)
679 vinfo(14)=vname(4,ifield)
680 vinfo(16)=vname(1,idtime)
681# if defined WRITE_WATER && defined MASKING
682 vinfo(20)='mask_rho'
683# endif
684 vinfo(21)=vname(6,ifield)
685 vinfo(22)='coordinates'
686 aval(5)=real(iinfo(1,ifield,ng),r8)
687 status=def_var(ng, inlm, dia(ng)%ncid, dia(ng)%Vid(ifield), &
688 & nf_fout, nvd5, l4dgrd, aval, vinfo, ncname)
689 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
690 END IF
691 END DO
692# endif
693# endif
694!
695!-----------------------------------------------------------------------
696! Leave definition mode.
697!-----------------------------------------------------------------------
698!
699 CALL netcdf_enddef (ng, inlm, ncname, dia(ng)%ncid)
700 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
701!
702!-----------------------------------------------------------------------
703! Write out time-recordless, information variables.
704!-----------------------------------------------------------------------
705!
706 CALL wrt_info (ng, inlm, dia(ng)%ncid, ncname)
707 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
708 END IF define
709!
710!=======================================================================
711! Open an existing diagnostics file, check its contents, and prepare
712! for appending data.
713!=======================================================================
714!
715 query : IF (.not.ldef) THEN
716 ncname=dia(ng)%name
717!
718! Open diagnostics file for read/write.
719!
720 CALL netcdf_open (ng, inlm, ncname, 1, dia(ng)%ncid)
721 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
722 WRITE (stdout,50) trim(ncname)
723 RETURN
724 END IF
725!
726! Inquire about the dimensions and check for consistency.
727!
728 CALL netcdf_check_dim (ng, inlm, ncname, &
729 & ncid = dia(ng)%ncid)
730 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
731!
732! Inquire about the variables.
733!
734 CALL netcdf_inq_var (ng, inlm, ncname, &
735 & ncid = dia(ng)%ncid)
736 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
737!
738! Initialize logical switches.
739!
740 DO i=1,nv
741 got_var(i)=.false.
742 END DO
743!
744! Scan variable list from input NetCDF and activate switches for
745! diagnostics variables. Get variable IDs.
746!
747 DO i=1,n_var
748 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
749 got_var(idtime)=.true.
750 dia(ng)%Vid(idtime)=var_id(i)
751 END IF
752# ifdef DIAGNOSTICS_UV
753 DO ivar=1,ndm2d
754 IF (trim(var_name(i)).eq. &
755 & trim(vname(1,iddu2d(ivar)))) THEN
756 got_var(iddu2d(ivar))=.true.
757 dia(ng)%Vid(iddu2d(ivar))=var_id(i)
758 ELSE IF (trim(var_name(i)).eq. &
759 & trim(vname(1,iddv2d(ivar)))) THEN
760 got_var(iddv2d(ivar))=.true.
761 dia(ng)%Vid(iddv2d(ivar))=var_id(i)
762 END IF
763 END DO
764# ifdef SOLVE3D
765 DO ivar=1,ndm3d
766 IF (trim(var_name(i)).eq. &
767 & trim(vname(1,iddu3d(ivar)))) THEN
768 got_var(iddu3d(ivar))=.true.
769 dia(ng)%Vid(iddu3d(ivar))=var_id(i)
770 ELSE IF (trim(var_name(i)).eq. &
771 & trim(vname(1,iddv3d(ivar)))) THEN
772 got_var(iddv3d(ivar))=.true.
773 dia(ng)%Vid(iddv3d(ivar))=var_id(i)
774 END IF
775 END DO
776# endif
777# endif
778# ifdef DIAGNOSTICS_TS
779 DO itrc=1,nt(ng)
780 DO ivar=1,ndt
781 ifield=iddtrc(itrc,ivar)
782 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
783 got_var(ifield)=.true.
784 dia(ng)%Vid(ifield)=var_id(i)
785 END IF
786 END DO
787 END DO
788# endif
789# ifdef DIAGNOSTICS_BIO
790# if defined BIO_FENNEL || defined HYPOXIA_SRM
791 DO ivar=1,ndbio2d
792 ifield=idbio2(ivar)
793 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
794 got_var(ifield)=.true.
795 dia(ng)%Vid(ifield)=var_id(i)
796 END IF
797 END DO
798# endif
799# if defined BIO_FENNEL
800 DO ivar=1,ndbio3d
801 ifield=idbio3(ivar)
802 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
803 got_var(ifield)=.true.
804 dia(ng)%Vid(ifield)=var_id(i)
805 END IF
806 END DO
807# elif defined ECOSIM
808 DO ivar=1,ndbio3d
809 ifield=idbio3(ivar)
810 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
811 got_var(ifield)=.true.
812 dia(ng)%Vid(ifield)=var_id(i)
813 END IF
814 END DO
815 DO ivar=1,ndbio4d
816 ifield=idbio4(ivar)
817 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
818 got_var(ifield)=.true.
819 dia(ng)%Vid(ifield)=var_id(i)
820 END IF
821 END DO
822# endif
823# endif
824 END DO
825!
826! Check if diagnostics variables are available in input NetCDF file.
827!
828 IF (.not.got_var(idtime)) THEN
829 IF (master) WRITE (stdout,60) trim(vname(1,idtime)), &
830 & trim(ncname)
831 exit_flag=3
832 RETURN
833 END IF
834# ifdef DIAGNOSTICS_UV
835 DO ivar=1,ndm2d
836 ifield=iddu2d(ivar)
837 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
838 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
839 & trim(ncname)
840 exit_flag=3
841 RETURN
842 END IF
843 ifield=iddv2d(ivar)
844 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
845 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
846 & trim(ncname)
847 exit_flag=3
848 RETURN
849 END IF
850 END DO
851# ifdef SOLVE3D
852 DO ivar=1,ndm3d
853 ifield=iddu3d(ivar)
854 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
855 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
856 & trim(ncname)
857 exit_flag=3
858 RETURN
859 END IF
860 ifield=iddv3d(ivar)
861 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
862 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
863 & trim(ncname)
864 exit_flag=3
865 RETURN
866 END IF
867 END DO
868# endif
869# endif
870# ifdef DIAGNOSTICS_TS
871 DO itrc=1,nt(ng)
872 DO ivar=1,ndt
873 ifield=iddtrc(itrc,ivar)
874 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
875 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
876 & trim(ncname)
877 exit_flag=3
878 RETURN
879 END IF
880 END DO
881 END DO
882# endif
883# ifdef DIAGNOSTICS_BIO
884# if defined BIO_FENNEL || defined HYPOXIA_SRM
885 DO ivar=1,ndbio2d
886 ifield=idbio2(ivar)
887 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
888 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
889 & trim(ncname)
890 exit_flag=3
891 RETURN
892 END IF
893 END DO
894# endif
895# if defined BIO_FENNEL
896 DO ivar=1,ndbio3d
897 ifield=idbio3(ivar)
898 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
899 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
900 & trim(ncname)
901 exit_flag=3
902 RETURN
903 END IF
904 END DO
905# elif defined ECOSIM
906 DO ivar=1,ndbio3d
907 ifield=idbio3(ivar)
908 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
909 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
910 & trim(ncname)
911 exit_flag=3
912 RETURN
913 END IF
914 END DO
915 DO ivar=1,ndbio4d
916 ifield=idbio4(ivar)
917 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
918 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
919 & trim(ncname)
920 exit_flag=3
921 RETURN
922 END IF
923 END DO
924# endif
925# endif
926!
927! Set unlimited time record dimension to the appropriate value.
928!
929 IF (nrst(ng).eq.ndia(ng)) THEN
930 IF (ndefdia(ng).gt.0) THEN
931 dia(ng)%Rindex=((ntstart(ng)-1)- &
932 & ndefdia(ng)*((ntstart(ng)-1)/ndefdia(ng)))/ &
933 & ndia(ng)
934 ELSE
935 dia(ng)%Rindex=(ntstart(ng)-1)/ndia(ng)
936 END IF
937 ELSE
938 dia(ng)%Rindex=rec_size
939 END IF
940 END IF query
941!
942! Set initial averaged time.
943!
944 IF (ntsdia(ng).eq.1) THEN
945 diatime(ng)=time(ng)-0.5_r8*real(ndia(ng),r8)*dt(ng)
946 ELSE
947 diatime(ng)=time(ng)+real(ntsdia(ng),r8)*dt(ng)- &
948 & 0.5_r8*real(ndia(ng),r8)*dt(ng)
949 END IF
950!
951 10 FORMAT (2x,'DEF_DIAGS_NF90 - creating diagnostics file,',t56, &
952 & 'Grid ',i2.2,': ',a)
953 20 FORMAT (2x,'DEF_DIAGS_NF90 - inquiring diagnostics file,',t56, &
954 & 'Grid ',i2.2,': ', a)
955 30 FORMAT (/,' DEF_DIAGS_NF90 - unable to create diagnostics NetCDF',&
956 & ' file: ',a)
957 40 FORMAT (1pe11.4,1x,'millimeter')
958 50 FORMAT (/,' DEF_DIAGS_NF90 - unable to open diagnostics NetCDF', &
959 & ' file: ',a)
960 60 FORMAT (/,' DEF_DIAGS_NF90 - unable to find variable: ',a,2x, &
961 & ' in diagnostics NetCDF file: ',a)
962!
963 RETURN
964 END SUBROUTINE def_diags_nf90
965
966# if defined PIO_LIB && defined DISTRIBUTE
967!
968!***********************************************************************
969 SUBROUTINE def_diags_pio (ng, ldef)
970!***********************************************************************
971!
973!
974! Imported variable declarations.
975!
976 integer, intent(in) :: ng
977
978 logical, intent(in) :: ldef
979!
980! Local variable declarations.
981!
982 logical :: got_var(nv)
983!
984 integer, parameter :: natt = 25
985
986 integer :: i, ifield, itrc, ivar, j, nvd3, nvd4, nvd5
987 integer :: recdim, status
988# if defined WRITE_WATER && defined MASKING
989 integer :: xy_pdim, xyz_pdim
990# endif
991 integer :: dimids(ndimid)
992 integer :: t2dgrd(3), u2dgrd(3), v2dgrd(3)
993# if defined ECOSIM && defined DIAGNOSTICS_BIO
994 integer :: l3dgrd(4), l4dgrd(5)
995# endif
996
997# ifdef SOLVE3D
998# ifdef SEDIMENT
999 integer :: b3dgrd(4)
1000# endif
1001 integer :: t3dgrd(4), u3dgrd(4), v3dgrd(4), w3dgrd(4)
1002# endif
1003!
1004 real(r8) :: aval(6)
1005!
1006 character (len= 13) :: prefix
1007 character (len=256) :: ncname
1008 character (len=MaxLen) :: vinfo(natt)
1009
1010 character (len=*), parameter :: myfile = &
1011 & __FILE__//", def_diags_pio"
1012!
1013 sourcefile=myfile
1014!
1015!-----------------------------------------------------------------------
1016! Set and report file name.
1017!-----------------------------------------------------------------------
1018!
1019 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1020 ncname=dia(ng)%name
1021!
1022 IF (master) THEN
1023 IF (ldef) THEN
1024 WRITE (stdout,10) ng, trim(ncname)
1025 ELSE
1026 WRITE (stdout,20) ng, trim(ncname)
1027 END IF
1028 END IF
1029!
1030!=======================================================================
1031! Create a new diagnostics NetCDF file.
1032!=======================================================================
1033!
1034 define : IF (ldef) THEN
1035 CALL pio_netcdf_create (ng, inlm, trim(ncname), dia(ng)%pioFile)
1036 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1037 IF (master) WRITE (stdout,30) trim(ncname)
1038 RETURN
1039 END IF
1040!
1041!-----------------------------------------------------------------------
1042! Define file dimensions.
1043!-----------------------------------------------------------------------
1044!
1045 dimids=0
1046!
1047 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xi_rho', &
1048 & iobounds(ng)%xi_rho, dimids( 1))
1049 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1050
1051 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xi_u', &
1052 & iobounds(ng)%xi_u, dimids( 2))
1053 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1054
1055 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xi_v', &
1056 & iobounds(ng)%xi_v, dimids( 3))
1057 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1058
1059 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xi_psi', &
1060 & iobounds(ng)%xi_psi, dimids( 4))
1061 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1062
1063 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'eta_rho', &
1064 & iobounds(ng)%eta_rho, dimids( 5))
1065 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1066
1067 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'eta_u', &
1068 & iobounds(ng)%eta_u, dimids( 6))
1069 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1070
1071 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'eta_v', &
1072 & iobounds(ng)%eta_v, dimids( 7))
1073 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1074
1075 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'eta_psi', &
1076 & iobounds(ng)%eta_psi, dimids( 8))
1077 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1078
1079# if defined WRITE_WATER && defined MASKING
1080 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xy_psi', &
1081 & iobounds(ng)%xy_psi, xy_pdim)
1082 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1083
1084 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xy_rho', &
1085 & iobounds(ng)%xy_rho, dimids(17))
1086 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1087
1088 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xy_u', &
1089 & iobounds(ng)%xy_u, dimids(18))
1090 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1091
1092 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xy_v', &
1093 & iobounds(ng)%xy_v, dimids(19))
1094 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1095# endif
1096
1097# ifdef SOLVE3D
1098# if defined WRITE_WATER && defined MASKING
1099 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xyz_psi', &
1100 & iobounds(ng)%xy_psi*n(ng), xyz_pdim)
1101 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1102
1103 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xyz_rho', &
1104 & iobounds(ng)%xy_rho*n(ng), dimids(20))
1105 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1106
1107 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xyz_u', &
1108 & iobounds(ng)%xy_u*n(ng), dimids(21))
1109 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1110
1111 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xyz_v', &
1112 & iobounds(ng)%xy_v*n(ng), dimids(22))
1113 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1114
1115 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xyz_w', &
1116 & iobounds(ng)%xy_rho*(n(ng)+1), dimids(23))
1117 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1118# endif
1119
1120 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 's_rho', &
1121 & n(ng), dimids( 9))
1122 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1123
1124 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 's_w', &
1125 & n(ng)+1, dimids(10))
1126 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1127
1128 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'tracer', &
1129 & nt(ng), dimids(11))
1130 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1131
1132# ifdef SEDIMENT
1133 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'NST', &
1134 & nst, dimids(32))
1135 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1136
1137 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'Nbed', &
1138 & nbed, dimids(16))
1139 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1140
1141# if defined WRITE_WATER && defined MASKING
1142 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'xybed', &
1143 & iobounds(ng)%xy_rho*nbed, dimids(24))
1144 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1145# endif
1146# endif
1147# ifdef ECOSIM
1148 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'Nbands', &
1149 & ndbands, dimids(33))
1150 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1151
1152 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'Nphy', &
1153 & nphy, dimids(25))
1154 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1155
1156 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'Nbac', &
1157 & nbac, dimids(26))
1158 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1159
1160 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'Ndom', &
1161 & ndom, dimids(27))
1162 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1163
1164 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'Nfec', &
1165 & nfec, dimids(28))
1166 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1167# endif
1168# endif
1169
1170 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'boundary', &
1171 & 4, dimids(14))
1172 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1173
1174# ifdef FOUR_DVAR
1175 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, 'Nstate', &
1176 & nstatevar(ng), dimids(29))
1177 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1178# endif
1179
1180 status=def_dim(ng, inlm, dia(ng)%pioFile, ncname, &
1181 & trim(adjustl(vname(5,idtime))), &
1182 & pio_unlimited, dimids(12))
1183 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1184
1185 recdim=dimids(12)
1186!
1187! Set number of dimensions for output variables.
1188!
1189# if defined WRITE_WATER && defined MASKING
1190 nvd3=2
1191 nvd4=2
1192 nvd5=2
1193# else
1194 nvd3=3
1195 nvd4=4
1196 nvd5=5
1197# endif
1198!
1199! Define dimension vectors for staggered tracer type variables.
1200!
1201# if defined WRITE_WATER && defined MASKING
1202 t2dgrd(1)=dimids(17)
1203 t2dgrd(2)=dimids(12)
1204# ifdef SOLVE3D
1205 t3dgrd(1)=dimids(20)
1206 t3dgrd(2)=dimids(12)
1207# endif
1208# else
1209 t2dgrd(1)=dimids( 1)
1210 t2dgrd(2)=dimids( 5)
1211 t2dgrd(3)=dimids(12)
1212# ifdef SOLVE3D
1213 t3dgrd(1)=dimids( 1)
1214 t3dgrd(2)=dimids( 5)
1215 t3dgrd(3)=dimids( 9)
1216 t3dgrd(4)=dimids(12)
1217# endif
1218# endif
1219!
1220! Define dimension vectors for staggered u-momentum type variables.
1221!
1222# if defined WRITE_WATER && defined MASKING
1223 u2dgrd(1)=dimids(18)
1224 u2dgrd(2)=dimids(12)
1225# ifdef SOLVE3D
1226 u3dgrd(1)=dimids(21)
1227 u3dgrd(2)=dimids(12)
1228# endif
1229# else
1230 u2dgrd(1)=dimids( 2)
1231 u2dgrd(2)=dimids( 6)
1232 u2dgrd(3)=dimids(12)
1233# ifdef SOLVE3D
1234 u3dgrd(1)=dimids( 2)
1235 u3dgrd(2)=dimids( 6)
1236 u3dgrd(3)=dimids( 9)
1237 u3dgrd(4)=dimids(12)
1238# endif
1239# endif
1240!
1241! Define dimension vectors for staggered v-momentum type variables.
1242!
1243# if defined WRITE_WATER && defined MASKING
1244 v2dgrd(1)=dimids(19)
1245 v2dgrd(2)=dimids(12)
1246# ifdef SOLVE3D
1247 v3dgrd(1)=dimids(22)
1248 v3dgrd(2)=dimids(12)
1249# endif
1250# else
1251 v2dgrd(1)=dimids( 3)
1252 v2dgrd(2)=dimids( 7)
1253 v2dgrd(3)=dimids(12)
1254# ifdef SOLVE3D
1255 v3dgrd(1)=dimids( 3)
1256 v3dgrd(2)=dimids( 7)
1257 v3dgrd(3)=dimids( 9)
1258 v3dgrd(4)=dimids(12)
1259# endif
1260# endif
1261# ifdef SOLVE3D
1262!
1263! Define dimension vector for staggered w-momentum type variables.
1264!
1265# if defined WRITE_WATER && defined MASKING
1266 w3dgrd(1)=dimids(23)
1267 w3dgrd(2)=dimids(12)
1268# else
1269 w3dgrd(1)=dimids( 1)
1270 w3dgrd(2)=dimids( 5)
1271 w3dgrd(3)=dimids(10)
1272 w3dgrd(4)=dimids(12)
1273# endif
1274# ifdef SEDIMENT
1275!
1276! Define dimension vector for sediment bed layer type variables.
1277!
1278# if defined WRITE_WATER && defined MASKING
1279 b3dgrd(1)=dimids(24)
1280 b3dgrd(2)=dimids(12)
1281# else
1282 b3dgrd(1)=dimids( 1)
1283 b3dgrd(2)=dimids( 5)
1284 b3dgrd(3)=dimids(16)
1285 b3dgrd(4)=dimids(12)
1286# endif
1287# endif
1288# if defined ECOSIM && defined DIAGNOSTICS_BIO
1289!
1290! Define dimension vector for spectral light type variables.
1291!
1292 l3dgrd(1)=dimids( 1)
1293 l3dgrd(2)=dimids( 5)
1294 l3dgrd(3)=dimids(33)
1295 l3dgrd(4)=dimids(12)
1296!
1297 l4dgrd(1)=dimids( 1)
1298 l4dgrd(2)=dimids( 5)
1299 l4dgrd(3)=dimids( 9)
1300 l4dgrd(4)=dimids(33)
1301 l4dgrd(5)=dimids(12)
1302# endif
1303# endif
1304!
1305! Initialize unlimited time record dimension.
1306!
1307 dia(ng)%Rindex=0
1308!
1309! Initialize local information variable arrays.
1310!
1311 DO i=1,natt
1312 DO j=1,len(vinfo(1))
1313 vinfo(i)(j:j)=' '
1314 END DO
1315 END DO
1316 DO i=1,6
1317 aval(i)=0.0_r8
1318 END DO
1319!
1320! Set long name prefix string.
1321!
1322!! Prefix='time-averaged'
1323 prefix=char(32) ! blank space
1324!
1325!-----------------------------------------------------------------------
1326! Define time-recordless information variables.
1327!-----------------------------------------------------------------------
1328!
1329 CALL def_info (ng, inlm, dia(ng)%pioFile, ncname, dimids)
1330 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1331!
1332!-----------------------------------------------------------------------
1333! Define variables and their attributes.
1334!-----------------------------------------------------------------------
1335!
1336! Define model time.
1337!
1338 vinfo( 1)=vname(1,idtime)
1339 WRITE (vinfo( 2),'(a,a)') 'averaged ', trim(vname(2,idtime))
1340 WRITE (vinfo( 3),'(a,a)') 'seconds since ', trim(rclock%string)
1341 vinfo( 4)=trim(rclock%calendar)
1342 vinfo(14)=vname(4,idtime)
1343 vinfo(21)=vname(6,idtime)
1344 dia(ng)%pioVar(idtime)%dkind=pio_tout
1345 dia(ng)%pioVar(idtime)%gtype=0
1346!
1347 status=def_var(ng, inlm, dia(ng)%pioFile, &
1348 & dia(ng)%pioVar(idtime)%vd, &
1349 & pio_tout, 1, (/recdim/), aval, vinfo, ncname, &
1350 & setparaccess = .true.)
1351 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1352!
1353! Define time-averaged free-surface. It is needed for CF compliance for
1354! vertical coordinates formula.
1355!
1356 vinfo( 1)=vname(1,idfsur)
1357 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idfsur))
1358 vinfo( 3)=vname(3,idfsur)
1359 vinfo(14)=vname(4,idfsur)
1360 vinfo(16)=vname(1,idtime)
1361# if !defined WET_DRY && (defined WRITE_WATER && defined MASKING)
1362 vinfo(20)='mask_rho'
1363# endif
1364 vinfo(21)=vname(6,idfsur)
1365 vinfo(22)='coordinates'
1366 aval(5)=real(iinfo(1,idfsur,ng),r8)
1367 avg(ng)%pioVar(idfsur)%dkind=pio_fout
1368 avg(ng)%pioVar(idfsur)%gtype=r2dvar
1369!
1370 status=def_var(ng, inlm, dia(ng)%pioFile, &
1371 & dia(ng)%pioVar(idfsur)%vd, &
1372 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1373 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1374
1375# ifdef DIAGNOSTICS_UV
1376!
1377! Define 2D momentum diagnostic fields.
1378!
1379 DO ivar=1,ndm2d
1380 ifield=iddu2d(ivar)
1381 IF (dout(ifield,ng)) THEN
1382 vinfo( 1)=vname(1,ifield)
1383 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1384 vinfo( 3)=vname(3,ifield)
1385 vinfo(14)=vname(4,ifield)
1386 vinfo(16)=vname(1,idtime)
1387# if defined WRITE_WATER && defined MASKING
1388 vinfo(20)='mask_u'
1389# endif
1390 vinfo(21)=vname(6,ifield)
1391 vinfo(22)='coordinates'
1392 aval(5)=real(u2dvar,r8)
1393 dia(ng)%pioVar(ifield)%dkind=pio_fout
1394 dia(ng)%pioVar(ifield)%gtype=u2dvar
1395!
1396 status=def_var(ng, inlm, dia(ng)%pioFile, &
1397 & dia(ng)%pioVar(ifield)%vd, &
1398 & pio_fout, nvd3, u2dgrd, aval, vinfo, ncname)
1399 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1400 END IF
1401!
1402 ifield=iddv2d(ivar)
1403 IF (dout(ifield,ng)) THEN
1404 vinfo( 1)=vname(1,ifield)
1405 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1406 vinfo( 3)=vname(3,ifield)
1407 vinfo(14)=vname(4,ifield)
1408 vinfo(16)=vname(1,idtime)
1409# if defined WRITE_WATER && defined MASKING
1410 vinfo(20)='mask_v'
1411# endif
1412 vinfo(21)=vname(6,ifield)
1413 vinfo(22)='coordinates'
1414 aval(5)=real(v2dvar,r8)
1415 dia(ng)%pioVar(ifield)%dkind=pio_fout
1416 dia(ng)%pioVar(ifield)%gtype=v2dvar
1417!
1418 status=def_var(ng, inlm, dia(ng)%pioFile, &
1419 & dia(ng)%pioVar(ifield)%vd, &
1420 & pio_fout, nvd3, v2dgrd, aval, vinfo, ncname)
1421 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1422 END IF
1423 END DO
1424
1425# ifdef SOLVE3D
1426!
1427! Define 3D momentum diagnostic fields.
1428!
1429 DO ivar=1,ndm3d
1430 ifield=iddu3d(ivar)
1431 IF (dout(ifield,ng)) THEN
1432 vinfo( 1)=vname(1,ifield)
1433 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1434 vinfo( 3)=vname(3,ifield)
1435 vinfo(14)=vname(4,ifield)
1436 vinfo(16)=vname(1,idtime)
1437# if defined WRITE_WATER && defined MASKING
1438 vinfo(20)='mask_u'
1439# endif
1440 vinfo(21)=vname(6,ifield)
1441 vinfo(22)='coordinates'
1442 aval(5)=real(u3dvar,r8)
1443 dia(ng)%pioVar(ifield)%dkind=pio_fout
1444 dia(ng)%pioVar(ifield)%gtype=u3dvar
1445!
1446 status=def_var(ng, inlm, dia(ng)%pioFile, &
1447 & dia(ng)%pioVar(ifield)%vd, &
1448 & pio_fout, nvd4, u3dgrd, aval, vinfo, ncname)
1449 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1450 END IF
1451!
1452 ifield=iddv3d(ivar)
1453 IF (dout(ifield,ng)) THEN
1454 vinfo( 1)=vname(1,ifield)
1455 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1456 vinfo( 3)=vname(3,ifield)
1457 vinfo(14)=vname(4,ifield)
1458 vinfo(16)=vname(1,idtime)
1459# if defined WRITE_WATER && defined MASKING
1460 vinfo(20)='mask_v'
1461# endif
1462 vinfo(21)=vname(6,ifield)
1463 vinfo(22)='coordinates'
1464 aval(5)=real(v3dvar,r8)
1465 dia(ng)%pioVar(ifield)%dkind=pio_fout
1466 dia(ng)%pioVar(ifield)%gtype=v3dvar
1467!
1468 status=def_var(ng, inlm, dia(ng)%pioFile, &
1469 & dia(ng)%pioVar(ifield)%vd, &
1470 & pio_fout, nvd4, v3dgrd, aval, vinfo, ncname)
1471 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1472 END IF
1473 END DO
1474# endif
1475# endif
1476# ifdef DIAGNOSTICS_TS
1477!
1478! Define tracer diagnostic fields.
1479!
1480 DO itrc=1,nt(ng)
1481 DO ivar=1,ndt
1482 ifield=iddtrc(itrc,ivar)
1483 IF (dout(ifield,ng)) THEN
1484 vinfo( 1)=vname(1,ifield)
1485 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1486 vinfo( 3)=vname(3,ifield)
1487 vinfo(14)=vname(4,ifield)
1488 vinfo(16)=vname(1,idtime)
1489# if defined WRITE_WATER && defined MASKING
1490 vinfo(20)='mask_rho'
1491# endif
1492 vinfo(21)=vname(6,ifield)
1493 vinfo(22)='coordinates'
1494 aval(5)=real(r3dvar,r8)
1495 dia(ng)%pioVar(ifield)%dkind=pio_fout
1496 dia(ng)%pioVar(ifield)%gtype=r3dvar
1497!
1498 status=def_var(ng, inlm, dia(ng)%pioFile, &
1499 & dia(ng)%pioVar(ifield)%vd, &
1500 & pio_fout, nvd4, t3dgrd, aval, vinfo, &
1501 & ncname)
1503 & __line__, myfile)) RETURN
1504 END IF
1505 END DO
1506 END DO
1507# endif
1508
1509# ifdef DIAGNOSTICS_BIO
1510# if defined BIO_FENNEL || defined HYPOXIA_SRM
1511!
1512! Define 2D biological diagnostic fields.
1513!
1514 DO ivar=1,ndbio2d
1515 ifield=idbio2(ivar)
1516 IF (dout(ifield,ng)) THEN
1517 vinfo( 1)=vname(1,ifield)
1518 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1519 vinfo( 3)=vname(3,ifield)
1520 vinfo(14)=vname(4,ifield)
1521 vinfo(16)=vname(1,idtime)
1522# if defined WRITE_WATER && defined MASKING
1523 vinfo(20)='mask_rho'
1524# endif
1525 vinfo(21)=vname(6,ifield)
1526 vinfo(22)='coordinates'
1527 aval(5)=real(r2dvar,r8)
1528 dia(ng)%pioVar(ifield)%dkind=pio_fout
1529 dia(ng)%pioVar(ifield)%gtype=r2dvar
1530!
1531 status=def_var(ng, inlm, dia(ng)%pioFile, &
1532 & dia(ng)%pioVar(ifield)%vd, &
1533 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1534 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1535 END IF
1536 END DO
1537# endif
1538# if defined BIO_FENNEL
1539!
1540! Define 3D biological diagnostic fields.
1541!
1542 DO ivar=1,ndbio3d
1543 ifield=idbio3(ivar)
1544 IF (dout(ifield,ng)) THEN
1545 vinfo( 1)=vname(1,ifield)
1546 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1547 vinfo( 3)=vname(3,ifield)
1548 vinfo(14)=vname(4,ifield)
1549 vinfo(16)=vname(1,idtime)
1550# if defined WRITE_WATER && defined MASKING
1551 vinfo(20)='mask_rho'
1552# endif
1553 vinfo(21)=vname(6,ifield)
1554 vinfo(22)='coordinates'
1555 aval(5)=real(r3dvar,r8)
1556 dia(ng)%pioVar(ifield)%dkind=pio_fout
1557 dia(ng)%pioVar(ifield)%gtype=r3dvar
1558!
1559 status=def_var(ng, inlm, dia(ng)%pioFile, &
1560 & dia(ng)%pioVar(ifield)%vd, &
1561 & pio_fout, nvd4, t3dgrd, aval, vinfo, ncname)
1562 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1563 END IF
1564 END DO
1565
1566# elif defined ECOSIM
1567!
1568! Define 3D biological diagnostic fields.
1569!
1570 DO ivar=1,ndbio3d
1571 ifield=idbio3(ivar)
1572 IF (dout(ifield,ng)) THEN
1573 vinfo( 1)=vname(1,ifield)
1574 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1575 vinfo( 3)=vname(3,ifield)
1576 vinfo(14)=vname(4,ifield)
1577 vinfo(16)=vname(1,idtime)
1578# if defined WRITE_WATER && defined MASKING
1579 vinfo(20)='mask_rho'
1580# endif
1581 vinfo(21)=vname(6,ifield)
1582 vinfo(22)='coordinates'
1583 aval(5)=real(iinfo(1,ifield,ng),r8)
1584 dia(ng)%pioVar(ifield)%dkind=pio_fout
1585 dia(ng)%pioVar(ifield)%gtype=l3dvar
1586!
1587 status=def_var(ng, inlm, dia(ng)%pioFile, &
1588 & dia(ng)%pioVar(ifield)%vd, &
1589 & pio_fout, nvd4, l3dgrd, aval, vinfo, ncname)
1590 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1591 END IF
1592 END DO
1593!
1594! Define 4D biological diagnostic fields.
1595!
1596 DO ivar=1,ndbio4d
1597 ifield=idbio4(ivar)
1598 IF (dout(ifield,ng)) THEN
1599 vinfo( 1)=vname(1,ifield)
1600 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,ifield))
1601 vinfo( 3)=vname(3,ifield)
1602 vinfo(14)=vname(4,ifield)
1603 vinfo(16)=vname(1,idtime)
1604# if defined WRITE_WATER && defined MASKING
1605 vinfo(20)='mask_rho'
1606# endif
1607 vinfo(21)=vname(6,ifield)
1608 vinfo(22)='coordinates'
1609 aval(5)=real(iinfo(1,ifield,ng),r8)
1610 dia(ng)%pioVar(ifield)%dkind=pio_fout
1611 dia(ng)%pioVar(ifield)%gtype=l4dvar
1612!
1613 status=def_var(ng, inlm, dia(ng)%pioFile, &
1614 & dia(ng)%pioVar(ifield)%vd, &
1615 & pio_fout, nvd5, l4dgrd, aval, vinfo, ncname)
1616 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1617 END IF
1618 END DO
1619# endif
1620# endif
1621!
1622!-----------------------------------------------------------------------
1623! Leave definition mode.
1624!-----------------------------------------------------------------------
1625!
1626 CALL pio_netcdf_enddef (ng, inlm, ncname, dia(ng)%pioFile)
1627 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1628!
1629!-----------------------------------------------------------------------
1630! Write out time-recordless, information variables.
1631!-----------------------------------------------------------------------
1632!
1633 CALL wrt_info (ng, inlm, dia(ng)%pioFile, ncname)
1634 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1635 END IF define
1636!
1637!=======================================================================
1638! Open an existing diagnostics file, check its contents, and prepare
1639! for appending data.
1640!=======================================================================
1641!
1642 query : IF (.not.ldef) THEN
1643 ncname=dia(ng)%name
1644!
1645! Open diagnostics file for read/write.
1646!
1647 CALL pio_netcdf_open (ng, inlm, ncname, 1, dia(ng)%pioFile)
1648 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1649 WRITE (stdout,50) trim(ncname)
1650 RETURN
1651 END IF
1652!
1653! Inquire about the dimensions and check for consistency.
1654!
1655 CALL pio_netcdf_check_dim (ng, inlm, ncname, &
1656 & piofile = dia(ng)%pioFile)
1657 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1658!
1659! Inquire about the variables.
1660!
1661 CALL pio_netcdf_inq_var (ng, inlm, ncname, &
1662 & piofile = dia(ng)%pioFile)
1663 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1664!
1665! Initialize logical switches.
1666!
1667 DO i=1,nv
1668 got_var(i)=.false.
1669 END DO
1670!
1671! Scan variable list from input NetCDF and activate switches for
1672! diagnostics variables. Get variable IDs.
1673!
1674 DO i=1,n_var
1675 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
1676 got_var(idtime)=.true.
1677 dia(ng)%pioVar(idtime)%vd=var_desc(i)
1678 dia(ng)%pioVar(idtime)%dkind=pio_tout
1679 dia(ng)%pioVar(idtime)%gtype=0
1680 END IF
1681# ifdef DIAGNOSTICS_UV
1682 DO ivar=1,ndm2d
1683 IF (trim(var_name(i)).eq. &
1684 & trim(vname(1,iddu2d(ivar)))) THEN
1685 got_var(iddu2d(ivar))=.true.
1686 dia(ng)%pioVar(iddu2d(ivar))%vd=var_desc(i)
1687 dia(ng)%pioVar(iddu2d(ivar))%dkind=pio_fout
1688 dia(ng)%pioVar(iddu2d(ivar))%gtype=u2dvar
1689 ELSE IF (trim(var_name(i)).eq. &
1690 & trim(vname(1,iddv2d(ivar)))) THEN
1691 got_var(iddv2d(ivar))=.true.
1692 dia(ng)%pioVar(iddv2d(ivar))%vd=var_desc(i)
1693 dia(ng)%pioVar(iddv2d(ivar))%dkind=pio_fout
1694 dia(ng)%pioVar(iddv2d(ivar))%gtype=v2dvar
1695 END IF
1696 END DO
1697# ifdef SOLVE3D
1698 DO ivar=1,ndm3d
1699 IF (trim(var_name(i)).eq. &
1700 & trim(vname(1,iddu3d(ivar)))) THEN
1701 got_var(iddu3d(ivar))=.true.
1702 dia(ng)%pioVar(iddu3d(ivar))%vd=var_desc(i)
1703 dia(ng)%pioVar(iddu3d(ivar))%dkind=pio_fout
1704 dia(ng)%pioVar(iddu3d(ivar))%gtype=u3dvar
1705 ELSE IF (trim(var_name(i)).eq. &
1706 & trim(vname(1,iddv3d(ivar)))) THEN
1707 got_var(iddv3d(ivar))=.true.
1708 dia(ng)%pioVar(iddv3d(ivar))%vd=var_desc(i)
1709 dia(ng)%pioVar(iddv3d(ivar))%dkind=pio_fout
1710 dia(ng)%pioVar(iddv3d(ivar))%gtype=v3dvar
1711 END IF
1712 END DO
1713# endif
1714# endif
1715# ifdef DIAGNOSTICS_TS
1716 DO itrc=1,nt(ng)
1717 DO ivar=1,ndt
1718 ifield=iddtrc(itrc,ivar)
1719 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
1720 got_var(ifield)=.true.
1721 dia(ng)%pioVar(ifield)%vd=var_desc(i)
1722 dia(ng)%pioVar(ifield)%dkind=pio_fout
1723 dia(ng)%pioVar(ifield)%gtype=r3dvar
1724 END IF
1725 END DO
1726 END DO
1727# endif
1728# ifdef DIAGNOSTICS_BIO
1729# if defined BIO_FENNEL || defined HYPOXIA_SRM
1730 DO ivar=1,ndbio2d
1731 ifield=idbio2(ivar)
1732 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
1733 got_var(ifield)=.true.
1734 dia(ng)%pioVar(ifield)%vd=var_desc(i)
1735 dia(ng)%pioVar(ifield)%dkind=pio_fout
1736 dia(ng)%pioVar(ifield)%gtype=r2dvar
1737 END IF
1738 END DO
1739# endif
1740# if defined BIO_FENNEL
1741 DO ivar=1,ndbio3d
1742 ifield=idbio3(ivar)
1743 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
1744 got_var(ifield)=.true.
1745 dia(ng)%pioVar(ifield)%vd=var_desc(i)
1746 dia(ng)%pioVar(ifield)%dkind=pio_fout
1747 dia(ng)%pioVar(ifield)%gtype=r3dvar
1748 END IF
1749 END DO
1750# elif defined ECOSIM
1751 DO ivar=1,ndbio3d
1752 ifield=idbio3(ivar)
1753 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
1754 got_var(ifield)=.true.
1755 dia(ng)%pioVar(ifield)%vd=var_desc(i)
1756 dia(ng)%pioVar(ifield)%dkind=pio_fout
1757 dia(ng)%pioVar(ifield)%gtype=l3dvar
1758 END IF
1759 END DO
1760 DO ivar=1,ndbio4d
1761 ifield=idbio4(ivar)
1762 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
1763 got_var(ifield)=.true.
1764 dia(ng)%pioVar(ifield)%vd=var_desc(i)
1765 dia(ng)%pioVar(ifield)%dkind=pio_fout
1766 dia(ng)%pioVar(ifield)%gtype=l4dvar
1767 END IF
1768 END DO
1769# endif
1770# endif
1771 END DO
1772!
1773! Check if diagnostics variables are available in input NetCDF file.
1774!
1775 IF (.not.got_var(idtime)) THEN
1776 IF (master) WRITE (stdout,60) trim(vname(1,idtime)), &
1777 & trim(ncname)
1778 exit_flag=3
1779 RETURN
1780 END IF
1781# ifdef DIAGNOSTICS_UV
1782 DO ivar=1,ndm2d
1783 ifield=iddu2d(ivar)
1784 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1785 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1786 & trim(ncname)
1787 exit_flag=3
1788 RETURN
1789 END IF
1790 ifield=iddv2d(ivar)
1791 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1792 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1793 & trim(ncname)
1794 exit_flag=3
1795 RETURN
1796 END IF
1797 END DO
1798# ifdef SOLVE3D
1799 DO ivar=1,ndm3d
1800 ifield=iddu3d(ivar)
1801 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1802 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1803 & trim(ncname)
1804 exit_flag=3
1805 RETURN
1806 END IF
1807 ifield=iddv3d(ivar)
1808 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1809 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1810 & trim(ncname)
1811 exit_flag=3
1812 RETURN
1813 END IF
1814 END DO
1815# endif
1816# endif
1817# ifdef DIAGNOSTICS_TS
1818 DO itrc=1,nt(ng)
1819 DO ivar=1,ndt
1820 ifield=iddtrc(itrc,ivar)
1821 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1822 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1823 & trim(ncname)
1824 exit_flag=3
1825 RETURN
1826 END IF
1827 END DO
1828 END DO
1829# endif
1830# ifdef DIAGNOSTICS_BIO
1831# if defined BIO_FENNEL || defined HYPOXIA_SRM
1832 DO ivar=1,ndbio2d
1833 ifield=idbio2(ivar)
1834 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1835 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1836 & trim(ncname)
1837 exit_flag=3
1838 RETURN
1839 END IF
1840 END DO
1841# endif
1842# if defined BIO_FENNEL
1843 DO ivar=1,ndbio3d
1844 ifield=idbio3(ivar)
1845 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1846 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1847 & trim(ncname)
1848 exit_flag=3
1849 RETURN
1850 END IF
1851 END DO
1852# elif defined ECOSIM
1853 DO ivar=1,ndbio3d
1854 ifield=idbio3(ivar)
1855 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1856 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1857 & trim(ncname)
1858 exit_flag=3
1859 RETURN
1860 END IF
1861 END DO
1862 DO ivar=1,ndbio4d
1863 ifield=idbio4(ivar)
1864 IF (.not.got_var(ifield).and.dout(ifield,ng)) THEN
1865 IF (master) WRITE (stdout,60) trim(vname(1,ifield)), &
1866 & trim(ncname)
1867 exit_flag=3
1868 RETURN
1869 END IF
1870 END DO
1871# endif
1872# endif
1873!
1874! Set unlimited time record dimension to the appropriate value.
1875!
1876 IF (nrst(ng).eq.ndia(ng)) THEN
1877 IF (ndefdia(ng).gt.0) THEN
1878 dia(ng)%Rindex=((ntstart(ng)-1)- &
1879 & ndefdia(ng)*((ntstart(ng)-1)/ndefdia(ng)))/ &
1880 & ndia(ng)
1881 ELSE
1882 dia(ng)%Rindex=(ntstart(ng)-1)/ndia(ng)
1883 END IF
1884 ELSE
1885 dia(ng)%Rindex=rec_size
1886 END IF
1887 END IF query
1888!
1889! Set initial averaged time.
1890!
1891 IF (ntsdia(ng).eq.1) THEN
1892 diatime(ng)=time(ng)-0.5_r8*real(ndia(ng),r8)*dt(ng)
1893 ELSE
1894 diatime(ng)=time(ng)+real(ntsdia(ng),r8)*dt(ng)- &
1895 & 0.5_r8*real(ndia(ng),r8)*dt(ng)
1896 END IF
1897!
1898 10 FORMAT (2x,'DEF_DIAGS_PIO - creating diagnostics file,',t56, &
1899 & 'Grid ',i2.2,': ',a)
1900 20 FORMAT (2x,'DEF_DIAGS_PIO - inquiring diagnostics file,',t56, &
1901 & 'Grid ',i2.2,': ',a)
1902 30 FORMAT (/,' DEF_DIAGS_PIO - unable to create diagnostics NetCDF', &
1903 & ' file: ',a)
1904 40 FORMAT (1pe11.4,1x,'millimeter')
1905 50 FORMAT (/,' DEF_DIAGS_PIO - unable to open diagnostics NetCDF', &
1906 & ' file: ',a)
1907 60 FORMAT (/,' DEF_DIAGS_PIO - unable to find variable: ',a,2x, &
1908 & ' in diagnostics NetCDF file: ',a)
1909!
1910 RETURN
1911 END SUBROUTINE def_diags_pio
1912# endif
1913#endif
1914 END MODULE def_diags_mod
subroutine, public def_diags(ng, ldef)
Definition def_diags.F:51
subroutine, private def_diags_nf90(ng, ldef)
Definition def_diags.F:90
subroutine, private def_diags_pio(ng, ldef)
Definition def_diags.F:970
integer, parameter nfec
Definition ecosim_mod.h:204
integer, dimension(:), allocatable idbio4
Definition ecosim_mod.h:304
integer, parameter nbac
Definition ecosim_mod.h:202
integer, parameter ndom
Definition ecosim_mod.h:203
integer ndbands
Definition ecosim_mod.h:212
integer, dimension(:), allocatable idbio3
Definition ecosim_mod.h:298
integer, dimension(:), allocatable idbio2
Definition fennel_mod.h:105
integer, parameter nphy
Definition ecosim_mod.h:205
integer, dimension(:), allocatable nstatevar
type(t_io), dimension(:), allocatable avg
integer stdout
character(len=256) sourcefile
type(t_io), dimension(:), allocatable dia
integer, dimension(:), allocatable iddu2d
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer, dimension(:), allocatable iddv2d
integer, parameter nv
integer, parameter io_pio
Definition mod_ncparam.F:96
integer, dimension(:), allocatable iddu3d
integer idfsur
integer, dimension(:), allocatable iddv3d
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer, dimension(:,:,:), allocatable iinfo
integer, parameter ndimid
logical, dimension(:,:), allocatable dout
integer, dimension(:,:), allocatable iddtrc
integer, parameter nf_tout
Definition mod_netcdf.F:207
subroutine, public netcdf_check_dim(ng, model, ncname, ncid)
Definition mod_netcdf.F:538
subroutine, public netcdf_open(ng, model, ncname, omode, ncid)
integer, parameter nf_fout
Definition mod_netcdf.F:188
subroutine, public netcdf_enddef(ng, model, ncname, ncid)
character(len=100), dimension(mvars) var_name
Definition mod_netcdf.F:169
integer, dimension(mvars) var_id
Definition mod_netcdf.F:160
integer n_var
Definition mod_netcdf.F:152
integer rec_size
Definition mod_netcdf.F:156
subroutine, public netcdf_create(ng, model, ncname, ncid)
subroutine, public netcdf_inq_var(ng, model, ncname, ncid, myvarname, searchvar, varid, nvardim, nvaratt)
logical master
integer ndbio2d
Definition mod_param.F:584
integer, parameter inlm
Definition mod_param.F:662
integer nbed
Definition mod_param.F:517
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, parameter r3dvar
Definition mod_param.F:721
integer, parameter l4dvar
Definition mod_param.F:727
integer ndbio4d
Definition mod_param.F:586
type(t_iobounds), dimension(:), allocatable iobounds
Definition mod_param.F:282
integer, parameter u3dvar
Definition mod_param.F:722
integer, parameter u2dvar
Definition mod_param.F:718
integer ndm3d
Definition mod_param.F:579
integer ndt
Definition mod_param.F:574
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer ndbio3d
Definition mod_param.F:585
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter l3dvar
Definition mod_param.F:726
integer, parameter v2dvar
Definition mod_param.F:719
integer nst
Definition mod_param.F:521
integer, parameter v3dvar
Definition mod_param.F:723
integer ndm2d
Definition mod_param.F:578
integer, parameter pio_fout
type(var_desc_t), dimension(:), pointer var_desc
subroutine, public pio_netcdf_create(ng, model, ncname, piofile)
subroutine, public pio_netcdf_inq_var(ng, model, ncname, piofile, myvarname, searchvar, piovar, nvardim, nvaratt)
subroutine, public pio_netcdf_open(ng, model, ncname, omode, piofile)
subroutine, public pio_netcdf_check_dim(ng, model, ncname, piofile)
integer, parameter pio_tout
subroutine, public pio_netcdf_enddef(ng, model, ncname, piofile)
real(dp), dimension(:), allocatable dt
integer, dimension(:), allocatable nrst
type(t_clock) rclock
integer exit_flag
real(dp), dimension(:), allocatable diatime
real(dp), dimension(:), allocatable time
integer, dimension(:), allocatable ntstart
integer, dimension(:), allocatable ndia
integer, dimension(:), allocatable ntsdia
integer noerror
integer, dimension(:), allocatable ndefdia
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52