ROMS
Loading...
Searching...
No Matches
def_state_mod Module Reference

Functions/Subroutines

subroutine, public def_state (ng, model, label, lcreate, s, lclose)
 
subroutine, private def_state_nf90 (ng, model, label, lcreate, s, lclose)
 
subroutine, private def_state_pio (ng, model, label, lcreate, s, lclose)
 

Function/Subroutine Documentation

◆ def_state()

subroutine, public def_state_mod::def_state ( integer, intent(in) ng,
integer, intent(in) model,
character (len=*), intent(in) label,
logical, intent(in) lcreate,
type(t_io), dimension(ngrids), intent(inout) s,
logical, intent(in) lclose )

Definition at line 58 of file def_state.F.

59!***********************************************************************
60!
61! Imported variable declarations.
62!
63 logical, intent(in) :: Lcreate, Lclose
64!
65 integer, intent(in) :: ng, model
66!
67 character (len=*), intent(in) :: label
68!
69 TYPE(T_IO), intent(inout) :: S(Ngrids)
70!
71! Local variable declarations.
72!
73 character (len=*), parameter :: MyFile = &
74 & __FILE__//", def_state_nf90"
75!
76!-----------------------------------------------------------------------
77! Create a new history file according to IO type.
78!-----------------------------------------------------------------------
79!
80 SELECT CASE (s(ng)%IOtype)
81 CASE (io_nf90)
82 CALL def_state_nf90 (ng, model, label, lcreate, s, lclose)
83
84#if defined PIO_LIB && defined DISTRIBUTE
85 CASE (io_pio)
86 CALL def_state_pio (ng, model, label, lcreate, s, lclose)
87#endif
88 CASE DEFAULT
89 IF (master) WRITE (stdout,10) avg(ng)%IOtype
90 exit_flag=3
91 END SELECT
92 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
93!
94 10 FORMAT (' DEF_STATE - Illegal output file type, io_type = ',i0, &
95 & /,13x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
96!
97 RETURN

References mod_iounits::avg, def_state_nf90(), def_state_pio(), mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::io_nf90, mod_ncparam::io_pio, mod_parallel::master, mod_scalars::noerror, and mod_iounits::stdout.

Here is the call graph for this function:

◆ def_state_nf90()

subroutine, private def_state_mod::def_state_nf90 ( integer, intent(in) ng,
integer, intent(in) model,
character (len=*), intent(in) label,
logical, intent(in) lcreate,
type(t_io), dimension(ngrids), intent(inout) s,
logical, intent(in) lclose )
private

Definition at line 101 of file def_state.F.

102!***********************************************************************
103!
104 USE mod_netcdf
105!
106! Imported variable declarations.
107!
108 logical, intent(in) :: Lcreate, Lclose
109!
110 integer, intent(in) :: ng, model
111!
112 character (len=*), intent(in) :: label
113!
114 TYPE(T_IO), intent(inout) :: S(Ngrids)
115!
116! Local variable declarations.
117!
118 logical :: got_var(NV)
119!
120 integer, parameter :: Natt = 25
121
122 integer :: i, j, ifield, itrc, nrec, nvd, nvd3, nvd4
123 integer :: recdim, status, varid
124#ifdef ADJUST_BOUNDARY
125 integer :: IorJdim, brecdim
126#endif
127#if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
128 integer :: frecdim
129#endif
130 integer :: DimIDs(nDimID)
131 integer :: t2dgrd(3), u2dgrd(3), v2dgrd(3)
132#ifdef ADJUST_BOUNDARY
133 integer :: t2dobc(4)
134#endif
135#ifdef SOLVE3D
136 integer :: t3dgrd(4), u3dgrd(4), v3dgrd(4), w3dgrd(4)
137# ifdef ADJUST_BOUNDARY
138 integer :: t3dobc(5)
139# endif
140# ifdef ADJUST_STFLUX
141 integer :: t3dfrc(4)
142# endif
143#endif
144#ifdef ADJUST_WSTRESS
145 integer :: u3dfrc(4), v3dfrc(4)
146#endif
147!
148 real(r8) :: Aval(6)
149!
150 character (len=256) :: ncname
151 character (len=MaxLen) :: Vinfo(Natt)
152
153 character (len=*), parameter :: MyFile = &
154 & __FILE__
155!
156 sourcefile=myfile
157!
158!-----------------------------------------------------------------------
159! Set and report file name.
160!-----------------------------------------------------------------------
161!
162 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
163 ncname=s(ng)%name
164!
165 IF (master) THEN
166 IF (lcreate) THEN
167 WRITE (stdout,10) trim(label), ng, trim(ncname)
168 ELSE
169 WRITE (stdout,20) trim(label), ng, trim(ncname)
170 END IF
171 END IF
172!
173!=======================================================================
174! Create a new ocean state NetCDF file.
175!=======================================================================
176!
177 define : IF (lcreate) THEN
178 CALL netcdf_create (ng, model, trim(ncname), s(ng)%ncid)
179 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
180 IF (master) WRITE (stdout,30) trim(label), trim(ncname)
181 RETURN
182 END IF
183!
184!-----------------------------------------------------------------------
185! Define file dimensions.
186!-----------------------------------------------------------------------
187!
188 dimids=0
189!
190 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xi_rho', &
191 & iobounds(ng)%xi_rho, dimids( 1))
192 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
193
194 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xi_u', &
195 & iobounds(ng)%xi_u, dimids( 2))
196 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
197
198 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xi_v', &
199 & iobounds(ng)%xi_v, dimids( 3))
200 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
201
202 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xi_psi', &
203 & iobounds(ng)%xi_psi, dimids( 4))
204 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
205
206 status=def_dim(ng, model, s(ng)%ncid, ncname, 'eta_rho', &
207 & iobounds(ng)%eta_rho, dimids( 5))
208 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
209
210 status=def_dim(ng, model, s(ng)%ncid, ncname, 'eta_u', &
211 & iobounds(ng)%eta_u, dimids( 6))
212 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
213
214 status=def_dim(ng, model, s(ng)%ncid, ncname, 'eta_v', &
215 & iobounds(ng)%eta_v, dimids( 7))
216 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
217
218 status=def_dim(ng, model, s(ng)%ncid, ncname, 'eta_psi', &
219 & iobounds(ng)%eta_psi, dimids( 8))
220 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
221
222#ifdef ADJUST_BOUNDARY
223 status=def_dim(ng, model, s(ng)%ncid, ncname, 'IorJ', &
224 & iobounds(ng)%IorJ, iorjdim)
225 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
226#endif
227
228#if defined WRITE_WATER && defined MASKING
229 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xy_rho', &
230 & iobounds(ng)%xy_rho, dimids(17))
231 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
232
233 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xy_u', &
234 & iobounds(ng)%xy_u, dimids(18))
235 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
236
237 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xy_v', &
238 & iobounds(ng)%xy_v, dimids(19))
239 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
240#endif
241#ifdef SOLVE3D
242# if defined WRITE_WATER && defined MASKING
243 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xyz_rho', &
244 & iobounds(ng)%xy_rho*n(ng), dimids(20))
245 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
246
247 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xyz_u', &
248 & iobounds(ng)%xy_u*n(ng), dimids(21))
249 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
250
251 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xyz_v', &
252 & iobounds(ng)%xy_v*n(ng), dimids(22))
253 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
254
255 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xyz_w', &
256 & iobounds(ng)%xy_rho*(n(ng)+1), dimids(23))
257 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
258# endif
259
260 status=def_dim(ng, model, s(ng)%ncid, ncname, 'N', &
261 & n(ng), dimids( 9))
262 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
263
264 status=def_dim(ng, model, s(ng)%ncid, ncname, 's_rho', &
265 & n(ng), dimids( 9))
266 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
267
268 status=def_dim(ng, model, s(ng)%ncid, ncname, 's_w', &
269 & n(ng)+1, dimids(10))
270 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
271
272 status=def_dim(ng, model, s(ng)%ncid, ncname, 'tracer', &
273 & nt(ng), dimids(11))
274 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
275
276# ifdef SEDIMENT
277 status=def_dim(ng, model, s(ng)%ncid, ncname, 'NST', &
278 & nst, dimids(32))
279 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
280
281 status=def_dim(ng, model, s(ng)%ncid, ncname, 'Nbed', &
282 & nbed, dimids(16))
283 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
284
285# if defined WRITE_WATER && defined MASKING
286 status=def_dim(ng, model, s(ng)%ncid, ncname, 'xybed', &
287 & iobounds(ng)%xy_rho*nbed, dimids(24))
288 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
289# endif
290# endif
291
292# ifdef ECOSIM
293 status=def_dim(ng, inlm, s(ng)%ncid, ncname, 'Nbands', &
294 & nbands, dimids(33))
295 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
296
297 status=def_dim(ng, model, s(ng)%ncid, ncname, 'Nphy', &
298 & nphy, dimids(25))
299 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
300
301 status=def_dim(ng, model, s(ng)%ncid, ncname, 'Nbac', &
302 & nbac, dimids(26))
303 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
304
305 status=def_dim(ng, model, s(ng)%ncid, ncname, 'Ndom', &
306 & ndom, dimids(27))
307 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
308
309 status=def_dim(ng, model, s(ng)%ncid, ncname, 'Nfec', &
310 & nfec, dimids(28))
311 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
312# endif
313#endif
314
315 status=def_dim(ng, model, s(ng)%ncid, ncname, 'boundary', &
316 & 4, dimids(14))
317 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
318
319#ifdef FOUR_DVAR
320 status=def_dim(ng, model, s(ng)%ncid, ncname, 'Nstate', &
321 & nstatevar(ng), dimids(29))
322 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
323#endif
324
325#if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
326 status=def_dim(ng, model, s(ng)%ncid, ncname, 'frc_adjust', &
327 & nfrec(ng), dimids(30))
328 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
329#endif
330
331#ifdef ADJUST_BOUNDARY
332 status=def_dim(ng, model, s(ng)%ncid, ncname, 'obc_adjust', &
333 & nbrec(ng), dimids(31))
334 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
335#endif
336
337 status=def_dim(ng, model, s(ng)%ncid, ncname, &
338 & trim(adjustl(vname(5,idtime))), &
339 & nf90_unlimited, dimids(12))
340 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
341
342 recdim=dimids(12)
343#if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
344 frecdim=dimids(30)
345#endif
346#ifdef ADJUST_BOUNDARY
347 brecdim=dimids(31)
348#endif
349!
350! Set number of dimensions for output variables.
351!
352#if defined WRITE_WATER && defined MASKING
353 nvd3=2
354 nvd4=2
355#else
356 nvd3=3
357 nvd4=4
358#endif
359!
360! Define dimension vectors for staggered tracer type variables.
361!
362#if defined WRITE_WATER && defined MASKING
363 t2dgrd(1)=dimids(17)
364 t2dgrd(2)=dimids(12)
365# ifdef SOLVE3D
366 t3dgrd(1)=dimids(20)
367 t3dgrd(2)=dimids(12)
368# endif
369#else
370 t2dgrd(1)=dimids( 1)
371 t2dgrd(2)=dimids( 5)
372 t2dgrd(3)=dimids(12)
373# ifdef SOLVE3D
374 t3dgrd(1)=dimids( 1)
375 t3dgrd(2)=dimids( 5)
376 t3dgrd(3)=dimids( 9)
377 t3dgrd(4)=dimids(12)
378# endif
379# ifdef ADJUST_STFLUX
380 t3dfrc(1)=dimids( 1)
381 t3dfrc(2)=dimids( 5)
382 t3dfrc(3)=frecdim
383 t3dfrc(4)=dimids(12)
384# endif
385#endif
386#ifdef ADJUST_BOUNDARY
387 t2dobc(1)=iorjdim
388 t2dobc(2)=dimids(14)
389 t2dobc(3)=brecdim
390 t2dobc(4)=dimids(12)
391# ifdef SOLVE3D
392 t3dobc(1)=iorjdim
393 t3dobc(2)=dimids( 9)
394 t3dobc(3)=dimids(14)
395 t3dobc(4)=brecdim
396 t3dobc(5)=dimids(12)
397# endif
398#endif
399!
400! Define dimension vectors for staggered u-momentum type variables.
401!
402#if defined WRITE_WATER && defined MASKING
403 u2dgrd(1)=dimids(18)
404 u2dgrd(2)=dimids(12)
405# ifdef SOLVE3D
406 u3dgrd(1)=dimids(21)
407 u3dgrd(2)=dimids(12)
408# endif
409#else
410 u2dgrd(1)=dimids( 2)
411 u2dgrd(2)=dimids( 6)
412 u2dgrd(3)=dimids(12)
413# ifdef SOLVE3D
414 u3dgrd(1)=dimids( 2)
415 u3dgrd(2)=dimids( 6)
416 u3dgrd(3)=dimids( 9)
417 u3dgrd(4)=dimids(12)
418# endif
419# ifdef ADJUST_WSTRESS
420 u3dfrc(1)=dimids( 2)
421 u3dfrc(2)=dimids( 6)
422 u3dfrc(3)=frecdim
423 u3dfrc(4)=dimids(12)
424# endif
425#endif
426!
427! Define dimension vectors for staggered v-momentum type variables.
428!
429#if defined WRITE_WATER && defined MASKING
430 v2dgrd(1)=dimids(19)
431 v2dgrd(2)=dimids(12)
432# ifdef SOLVE3D
433 v3dgrd(1)=dimids(22)
434 v3dgrd(2)=dimids(12)
435# endif
436#else
437 v2dgrd(1)=dimids( 3)
438 v2dgrd(2)=dimids( 7)
439 v2dgrd(3)=dimids(12)
440# ifdef SOLVE3D
441 v3dgrd(1)=dimids( 3)
442 v3dgrd(2)=dimids( 7)
443 v3dgrd(3)=dimids( 9)
444 v3dgrd(4)=dimids(12)
445# endif
446# ifdef ADJUST_WSTRESS
447 v3dfrc(1)=dimids( 3)
448 v3dfrc(2)=dimids( 7)
449 v3dfrc(3)=frecdim
450 v3dfrc(4)=dimids(12)
451# endif
452#endif
453#ifdef SOLVE3D
454!
455! Define dimension vector for staggered w-momentum type variables.
456!
457# if defined WRITE_WATER && defined MASKING
458 w3dgrd(1)=dimids(23)
459 w3dgrd(2)=dimids(12)
460# else
461 w3dgrd(1)=dimids( 1)
462 w3dgrd(2)=dimids( 5)
463 w3dgrd(3)=dimids(10)
464 w3dgrd(4)=dimids(12)
465# endif
466#endif
467!
468! Initialize unlimited time record dimension.
469!
470 s(ng)%Rindex=0
471!
472! Initialize local information variable arrays.
473!
474 DO i=1,natt
475 DO j=1,len(vinfo(1))
476 vinfo(i)(j:j)=' '
477 END DO
478 END DO
479 DO i=1,6
480 aval(i)=0.0_r8
481 END DO
482!
483!-----------------------------------------------------------------------
484! Define time-recordless information variables.
485!-----------------------------------------------------------------------
486!
487 CALL def_info (ng, model, s(ng)%ncid, ncname, dimids)
488 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
489!
490!-----------------------------------------------------------------------
491! Define time-varying variables.
492!-----------------------------------------------------------------------
493!
494! Define model time.
495!
496 vinfo( 1)=vname(1,idtime)
497 vinfo( 2)=vname(2,idtime)
498 WRITE (vinfo( 3),'(a,a)') 'seconds since ', trim(rclock%string)
499 vinfo( 4)=trim(rclock%calendar)
500 vinfo(14)=vname(4,idtime)
501 vinfo(21)=vname(6,idtime)
502 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idtime), &
503 & nf_tout, 1, (/recdim/), aval, vinfo,ncname, &
504 & setparaccess = .true.)
505 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
506!
507! Define free-surface.
508!
509 vinfo( 1)=vname(1,idfsur)
510 vinfo( 2)=vname(2,idfsur)
511 vinfo( 3)='nondimensional'
512 vinfo(14)=vname(4,idfsur)
513 vinfo(16)=vname(1,idtime)
514#if defined WRITE_WATER && defined MASKING
515 vinfo(20)='mask_rho'
516#endif
517 vinfo(21)=vname(6,idfsur)
518 vinfo(22)='coordinates'
519 aval(5)=real(iinfo(1,idfsur,ng),r8)
520 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idfsur), &
521 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
522 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
523
524#ifdef ADJUST_BOUNDARY
525!
526! Define free-surface open boundaries.
527!
528 IF (any(lobc(:,isfsur,ng))) THEN
529 ifield=idsbry(isfsur)
530 vinfo( 1)=vname(1,ifield)
531 vinfo( 2)=vname(2,ifield)
532 vinfo( 3)='nondimensional'
533 vinfo(14)=vname(4,ifield)
534 vinfo(16)=vname(1,idtime)
535 vinfo(21)=vname(6,ifield)
536 aval(5)=real(iinfo(1,ifield,ng),r8)
537 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(ifield), &
538 & nf_fout, 4, t2dobc, aval, vinfo, ncname, &
539 & setfillval = .false.)
540 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
541 END IF
542#endif
543!
544! Define 2D U-momentum component.
545!
546 vinfo( 1)=vname(1,idubar)
547 vinfo( 2)=vname(2,idubar)
548 vinfo( 3)='nondimensional'
549 vinfo(14)=vname(4,idubar)
550 vinfo(16)=vname(1,idtime)
551#if defined WRITE_WATER && defined MASKING
552 vinfo(20)='mask_u'
553#endif
554 vinfo(21)=vname(6,idubar)
555 vinfo(22)='coordinates'
556 aval(5)=real(iinfo(1,idubar,ng),r8)
557 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idubar), &
558 & nf_fout, nvd3, u2dgrd, aval, vinfo, ncname)
559 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
560
561#ifdef ADJUST_BOUNDARY
562!
563! Define 2D U-momentum component open boundaries.
564!
565 IF (any(lobc(:,isubar,ng))) THEN
566 ifield=idsbry(isubar)
567 vinfo( 1)=vname(1,ifield)
568 vinfo( 2)=vname(2,ifield)
569 vinfo( 3)='nondimensional'
570 vinfo(14)=vname(4,ifield)
571 vinfo(16)=vname(1,idtime)
572 vinfo(21)=vname(6,ifield)
573 aval(5)=real(iinfo(1,ifield,ng),r8)
574 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(ifield), &
575 & nf_fout, 4, t2dobc, aval, vinfo, ncname, &
576 & setfillval = .false.)
577 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
578 END IF
579#endif
580!
581! Define 2D V-momentum component.
582!
583 vinfo( 1)=vname(1,idvbar)
584 vinfo( 2)=vname(2,idvbar)
585 vinfo( 3)='nondimensional'
586 vinfo(14)=vname(4,idvbar)
587 vinfo(16)=vname(1,idtime)
588#if defined WRITE_WATER && defined MASKING
589 vinfo(20)='mask_v'
590#endif
591 vinfo(21)=vname(6,idvbar)
592 vinfo(22)='coordinates'
593 aval(5)=real(iinfo(1,idvbar,ng),r8)
594 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idvbar), &
595 & nf_fout, nvd3, v2dgrd, aval, vinfo, ncname)
596 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
597
598#ifdef ADJUST_BOUNDARY
599!
600! Define 2D V-momentum component open boundaries.
601!
602 IF (any(lobc(:,isvbar,ng))) THEN
603 ifield=idsbry(isvbar)
604 vinfo( 1)=vname(1,ifield)
605 vinfo( 2)=vname(2,ifield)
606 vinfo( 3)='nondimensional'
607 vinfo(14)=vname(4,ifield)
608 vinfo(16)=vname(1,idtime)
609 vinfo(21)=vname(6,ifield)
610 aval(5)=real(iinfo(1,ifield,ng),r8)
611 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(ifield), &
612 & nf_fout, 4, t2dobc, aval, vinfo, ncname, &
613 & setfillval = .false.)
614 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
615 END IF
616#endif
617#ifdef SOLVE3D
618!
619! Define 3D U-momentum component.
620!
621 vinfo( 1)=vname(1,iduvel)
622 vinfo( 2)=vname(2,iduvel)
623 vinfo( 3)='nondimensional'
624 vinfo(14)=vname(4,iduvel)
625 vinfo(16)=vname(1,idtime)
626# if defined WRITE_WATER && defined MASKING
627 vinfo(20)='mask_u'
628# endif
629 vinfo(21)=vname(6,iduvel)
630 vinfo(22)='coordinates'
631 aval(5)=real(iinfo(1,iduvel,ng),r8)
632 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(iduvel), &
633 & nf_fout, nvd4, u3dgrd, aval, vinfo, ncname)
634 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
635
636# ifdef ADJUST_BOUNDARY
637!
638! Define 3D U-momentum component open boundaries.
639!
640 IF (any(lobc(:,isuvel,ng))) THEN
641 ifield=idsbry(isuvel)
642 vinfo( 1)=vname(1,ifield)
643 vinfo( 2)=vname(2,ifield)
644 vinfo( 3)='nondimensional'
645 vinfo(14)=vname(4,ifield)
646 vinfo(16)=vname(1,idtime)
647 vinfo(21)=vname(6,ifield)
648 aval(5)=real(iinfo(1,ifield,ng),r8)
649 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(ifield), &
650 & nf_fout, 5, t3dobc, aval, vinfo, ncname, &
651 & setfillval = .false.)
652 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
653 END IF
654# endif
655!
656! Define 3D V-momentum component.
657!
658 vinfo( 1)=vname(1,idvvel)
659 vinfo( 2)=vname(2,idvvel)
660 vinfo( 3)='nondimensional'
661 vinfo(14)=vname(4,idvvel)
662 vinfo(16)=vname(1,idtime)
663# if defined WRITE_WATER && defined MASKING
664 vinfo(20)='mask_v'
665# endif
666 vinfo(21)=vname(6,idvvel)
667 vinfo(22)='coordinates'
668 aval(5)=real(iinfo(1,idvvel,ng),r8)
669 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idvvel), &
670 & nf_fout, nvd4, v3dgrd, aval, vinfo, ncname)
671 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
672
673# ifdef ADJUST_BOUNDARY
674!
675! Define 3D V-momentum component open boundaries.
676!
677 IF (any(lobc(:,isvvel,ng))) THEN
678 ifield=idsbry(isvvel)
679 vinfo( 1)=vname(1,ifield)
680 vinfo( 2)=vname(2,ifield)
681 vinfo( 3)='nondimensional'
682 vinfo(14)=vname(4,ifield)
683 vinfo(16)=vname(1,idtime)
684 vinfo(21)=vname(6,ifield)
685 aval(5)=real(iinfo(1,ifield,ng),r8)
686 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(ifield), &
687 & nf_fout, 5, t3dobc, aval, vinfo, ncname, &
688 & setfillval = .false.)
689 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
690 END IF
691# endif
692!
693! Define tracer type variables.
694!
695 DO itrc=1,nt(ng)
696 vinfo( 1)=vname(1,idtvar(itrc))
697 vinfo( 2)=vname(2,idtvar(itrc))
698 vinfo( 3)='nondimensional'
699 vinfo(14)=vname(4,idtvar(itrc))
700 vinfo(16)=vname(1,idtime)
701# ifdef SEDIMENT
702 DO i=1,nst
703 IF (itrc.eq.idsed(i)) THEN
704 WRITE (vinfo(19),40) 1000.0_r8*sd50(i,ng)
705 END IF
706 END DO
707# endif
708# if defined WRITE_WATER && defined MASKING
709 vinfo(20)='mask_rho'
710# endif
711 vinfo(21)=vname(6,idtvar(itrc))
712 vinfo(22)='coordinates'
713 aval(5)=real(r3dvar,r8)
714 status=def_var(ng, model, s(ng)%ncid, s(ng)%Tid(itrc), &
715 & nf_fout, nvd4, t3dgrd, aval, vinfo, ncname)
716 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
717 END DO
718
719# ifdef ADJUST_BOUNDARY
720!
721! Define tracer type variables open boundaries.
722!
723 DO itrc=1,nt(ng)
724 IF (any(lobc(:,istvar(itrc),ng))) THEN
725 ifield=idsbry(istvar(itrc))
726 vinfo( 1)=vname(1,ifield)
727 vinfo( 2)=vname(2,ifield)
728 vinfo( 3)='nondimensional'
729 vinfo(14)=vname(4,ifield)
730 vinfo(16)=vname(1,idtime)
731# ifdef SEDIMENT
732 DO i=1,nst
733 IF (itrc.eq.idsed(i)) THEN
734 WRITE (vinfo(19),40) 1000.0_r8*sd50(i,ng)
735 END IF
736 END DO
737# endif
738 vinfo(21)=vname(6,ifield)
739 aval(5)=real(iinfo(1,ifield,ng),r8)
740 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(ifield), &
741 & nf_fout, 5, t3dobc, aval, vinfo, ncname, &
742 & setfillval = .false.)
743 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
744 END IF
745 END DO
746# endif
747# ifdef ADJUST_STFLUX
748!
749! Define surface tracer fluxes.
750!
751 DO itrc=1,nt(ng)
752 IF (lstflux(itrc,ng)) THEN
753 vinfo( 1)=vname(1,idtsur(itrc))
754 vinfo( 2)=vname(2,idtsur(itrc))
755 vinfo( 3)='nondimensional'
756 IF (itrc.eq.itemp) THEN
757 vinfo(11)='upward flux, cooling'
758 vinfo(12)='downward flux, heating'
759 ELSE IF (itrc.eq.isalt) THEN
760 vinfo(11)='upward flux, freshening (net precipitation)'
761 vinfo(12)='downward flux, salting (net evaporation)'
762 END IF
763 vinfo(14)=vname(4,idtsur(itrc))
764 vinfo(16)=vname(1,idtime)
765# if defined WRITE_WATER && defined MASKING
766 vinfo(20)='mask_rho'
767# endif
768 vinfo(21)=vname(6,idtsur(itrc))
769 vinfo(22)='coordinates'
770 aval(5)=real(r2dvar,r8)
771 status=def_var(ng, model, s(ng)%ncid, &
772 & s(ng)%Vid(idtsur(itrc)), &
773 & nf_fout, nvd4, t3dfrc, aval, vinfo, ncname)
774 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
775 END IF
776 END DO
777# endif
778#endif
779#ifdef ADJUST_WSTRESS
780!
781! Define surface U-momentum stress.
782!
783 vinfo( 1)=vname(1,idusms)
784 vinfo( 2)=vname(2,idusms)
785 vinfo( 3)='nondimensional'
786 vinfo(14)=vname(4,idusms)
787 vinfo(16)=vname(1,idtime)
788# if defined WRITE_WATER && defined MASKING
789 vinfo(20)='mask_u'
790# endif
791 vinfo(21)=vname(6,idusms)
792 vinfo(22)='coordinates'
793 aval(5)=real(u2dvar,r8)
794 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idusms), &
795 & nf_fout, nvd4, u3dfrc, aval, vinfo, ncname)
796 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
797!
798! Define surface V-momentum stress.
799!
800 vinfo( 1)=vname(1,idvsms)
801 vinfo( 2)=vname(2,idvsms)
802 vinfo( 3)='nondimensional'
803 vinfo(14)=vname(4,idvsms)
804 vinfo(16)=vname(1,idtime)
805# if defined WRITE_WATER && defined MASKING
806 vinfo(20)='mask_v'
807# endif
808 vinfo(21)=vname(6,idvsms)
809 vinfo(22)='coordinates'
810 aval(5)=real(v2dvar,r8)
811 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idvsms), &
812 & nf_fout, nvd4, v3dfrc, aval, vinfo, ncname)
813 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
814#endif
815!
816!-----------------------------------------------------------------------
817! Leave definition mode.
818!-----------------------------------------------------------------------
819!
820 CALL netcdf_enddef (ng, model, ncname, s(ng)%ncid)
821 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
822!
823!-----------------------------------------------------------------------
824! Write out time-recordless, information variables.
825!-----------------------------------------------------------------------
826!
827 CALL wrt_info (ng, model, s(ng)%ncid, ncname)
828 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
829!
830!-----------------------------------------------------------------------
831! Write out time-recordless, information variables.
832!-----------------------------------------------------------------------
833!
834 IF (lclose) THEN
835 CALL netcdf_close (ng, model, s(ng)%ncid, ncname, .false.)
836 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
837 END IF
838
839 END IF define
840!
841!=======================================================================
842! Open an existing ocean state NetCDF file, check its contents, and
843! prepare for appending data.
844!=======================================================================
845!
846 query: IF (.not.lcreate) THEN
847 ncname=s(ng)%name
848!
849! Open ocean state file for read/write.
850!
851 CALL netcdf_open (ng, model, ncname, 1, s(ng)%ncid)
852 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
853 WRITE (stdout,50) trim(label), trim(ncname)
854 RETURN
855 END IF
856!
857! Inquire about the dimensions and check for consistency.
858!
859 CALL netcdf_check_dim (ng, model, ncname, &
860 & ncid = s(ng)%ncid)
861 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
862!
863! Inquire about the variables.
864!
865 CALL netcdf_inq_var (ng, model, ncname, &
866 & ncid = s(ng)%ncid)
867 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
868!
869! Initialize logical switches.
870!
871 DO i=1,nv
872 got_var(i)=.false.
873 END DO
874!
875! Scan variable list from input NetCDF and activate switches for
876! Hessian eigenvectors variables. Get variable IDs.
877!
878 DO i=1,n_var
879 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
880 got_var(idtime)=.true.
881 s(ng)%Vid(idtime)=var_id(i)
882 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
883 got_var(idfsur)=.true.
884 s(ng)%Vid(idfsur)=var_id(i)
885 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
886 got_var(idubar)=.true.
887 s(ng)%Vid(idubar)=var_id(i)
888 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
889 got_var(idvbar)=.true.
890 s(ng)%Vid(idvbar)=var_id(i)
891#ifdef ADJUST_BOUNDARY
892 ELSE IF (trim(var_name(i)).eq. &
893 & trim(vname(1,idsbry(isfsur)))) THEN
894 got_var(idsbry(isfsur))=.true.
895 s(ng)%Vid(idsbry(isfsur))=var_id(i)
896 ELSE IF (trim(var_name(i)).eq. &
897 & trim(vname(1,idsbry(isubar)))) THEN
898 got_var(idsbry(isubar))=.true.
899 s(ng)%Vid(idsbry(isubar))=var_id(i)
900 ELSE IF (trim(var_name(i)).eq. &
901 & trim(vname(1,idsbry(isvbar)))) THEN
902 got_var(idsbry(isvbar))=.true.
903 s(ng)%Vid(idsbry(isvbar))=var_id(i)
904#endif
905#ifdef ADJUST_WSTRESS
906 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idusms))) THEN
907 got_var(idusms)=.true.
908 s(ng)%Vid(idusms)=var_id(i)
909 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvsms))) THEN
910 got_var(idvsms)=.true.
911 s(ng)%Vid(idvsms)=var_id(i)
912#endif
913#ifdef SOLVE3D
914 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
915 got_var(iduvel)=.true.
916 s(ng)%Vid(iduvel)=var_id(i)
917 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
918 got_var(idvvel)=.true.
919 s(ng)%Vid(idvvel)=var_id(i)
920# ifdef ADJUST_BOUNDARY
921 ELSE IF (trim(var_name(i)).eq. &
922 & trim(vname(1,idsbry(isuvel)))) THEN
923 got_var(idsbry(isuvel))=.true.
924 s(ng)%Vid(idsbry(isuvel))=var_id(i)
925 ELSE IF (trim(var_name(i)).eq. &
926 & trim(vname(1,idsbry(isvvel)))) THEN
927 got_var(idsbry(isvvel))=.true.
928 s(ng)%Vid(idsbry(isvvel))=var_id(i)
929# endif
930#endif
931 END IF
932#ifdef SOLVE3D
933 DO itrc=1,nt(ng)
934 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
935 got_var(idtvar(itrc))=.true.
936 s(ng)%Tid(itrc)=var_id(i)
937# ifdef ADJUST_BOUNDARY
938 ELSE IF (trim(var_name(i)).eq. &
939 & trim(vname(1,idsbry(istvar(itrc))))) THEN
940 got_var(idsbry(istvar(itrc)))=.true.
941 s(ng)%Vid(idsbry(istvar(itrc)))=var_id(i)
942# endif
943# ifdef ADJUST_STFLUX
944 ELSE IF (trim(var_name(i)).eq. &
945 & trim(vname(1,idtsur(itrc)))) THEN
946 got_var(idtsur(itrc))=.true.
947 s(ng)%Vid(idtsur(itrc))=var_id(i)
948# endif
949 END IF
950 END DO
951#endif
952 END DO
953!
954! Check if ocean state variables are available in input NetCDF file.
955!
956 IF (.not.got_var(idtime)) THEN
957 IF (master) WRITE (stdout,60) trim(vname(1,idtime)), &
958 & trim(label), trim(ncname)
959 exit_flag=3
960 RETURN
961 END IF
962 IF (.not.got_var(idfsur)) THEN
963 IF (master) WRITE (stdout,60) trim(vname(1,idfsur)), &
964 & trim(label), trim(ncname)
965 exit_flag=3
966 RETURN
967 END IF
968 IF (.not.got_var(idubar)) THEN
969 IF (master) WRITE (stdout,60) trim(vname(1,idubar)), &
970 & trim(label), trim(ncname)
971 exit_flag=3
972 RETURN
973 END IF
974 IF (.not.got_var(idvbar)) THEN
975 IF (master) WRITE (stdout,60) trim(vname(1,idvbar)), &
976 & trim(label), trim(ncname)
977 exit_flag=3
978 RETURN
979 END IF
980#ifdef ADJUST_BOUNDARY
981 IF (.not.got_var(idsbry(isfsur))) THEN
982 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isfsur))), &
983 & trim(label), trim(ncname)
984 exit_flag=3
985 RETURN
986 END IF
987 IF (.not.got_var(idsbry(isubar))) THEN
988 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isubar))), &
989 & trim(label), trim(ncname)
990 exit_flag=3
991 RETURN
992 END IF
993 IF (.not.got_var(idsbry(isvbar))) THEN
994 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isvbar))), &
995 & trim(label), trim(ncname)
996 exit_flag=3
997 RETURN
998 END IF
999#endif
1000#ifdef ADJUST_WSTRESS
1001 IF (.not.got_var(idusms)) THEN
1002 IF (master) WRITE (stdout,60) trim(vname(1,idusms)), &
1003 & trim(label), trim(ncname)
1004 exit_flag=3
1005 RETURN
1006 END IF
1007 IF (.not.got_var(idvsms)) THEN
1008 IF (master) WRITE (stdout,60) trim(vname(1,idvsms)), &
1009 & trim(label), trim(ncname)
1010 exit_flag=3
1011 RETURN
1012 END IF
1013#endif
1014#ifdef SOLVE3D
1015 IF (.not.got_var(iduvel)) THEN
1016 IF (master) WRITE (stdout,60) trim(vname(1,iduvel)), &
1017 & trim(label), trim(ncname)
1018 exit_flag=3
1019 RETURN
1020 END IF
1021 IF (.not.got_var(idvvel)) THEN
1022 IF (master) WRITE (stdout,60) trim(vname(1,idvvel)), &
1023 & trim(label), trim(ncname)
1024 exit_flag=3
1025 RETURN
1026 END IF
1027# ifdef ADJUST_BOUNDARY
1028 IF (.not.got_var(idsbry(isuvel))) THEN
1029 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isuvel))), &
1030 & trim(label), trim(ncname)
1031 exit_flag=3
1032 RETURN
1033 END IF
1034 IF (.not.got_var(idsbry(isvvel))) THEN
1035 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isvvel))), &
1036 & trim(label), trim(ncname)
1037 RETURN
1038 END IF
1039# endif
1040#endif
1041#ifdef SOLVE3D
1042 DO itrc=1,nt(ng)
1043 IF (.not.got_var(idtvar(itrc))) THEN
1044 IF (master) WRITE (stdout,60) trim(vname(1,idtvar(itrc))), &
1045 & trim(label), trim(ncname)
1046 exit_flag=3
1047 RETURN
1048 END IF
1049# ifdef ADJUST_BOUNDARY
1050 IF (.not.got_var(idsbry(istvar(itrc)))) THEN
1051 IF (master) WRITE (stdout,60) &
1052 & trim(vname(1,idsbry(istvar(itrc)))), &
1053 & trim(label), trim(ncname)
1054 exit_flag=3
1055 RETURN
1056 END IF
1057# endif
1058# ifdef ADJUST_STFLUX
1059 IF (.not.got_var(idtsur(itrc)).and.lstflux(itrc,ng)) THEN
1060 IF (master) WRITE (stdout,60) trim(vname(1,idtsur(itrc))), &
1061 & trim(label), trim(ncname)
1062 exit_flag=3
1063 RETURN
1064 END IF
1065# endif
1066 END DO
1067#endif
1068!
1069! Set unlimited time record dimension to the appropriate value.
1070!
1071 s(ng)%Rindex=rec_size
1072 END IF query
1073!
1074 10 FORMAT (2x,'DEF_STATE_NF90 - creating ',a,' file,',t56, &
1075 & 'Grid ',i2.2,': ',a)
1076 20 FORMAT (2x,'DEF_STATE_NF90 - inquiring ',a,' file,',t56, &
1077 & 'Grid ',i2.2,': ',a)
1078 30 FORMAT (/,' DEF_STATE_NF90 - unable to create ',a, &
1079 & ' NetCDF file:',1x,a)
1080 40 FORMAT (1pe11.4,1x,'millimeter')
1081 50 FORMAT (/,' DEF_STATE_NF90 - unable to open ',a, &
1082 & ' NetCDF file: ',a)
1083 60 FORMAT (/,' DEF_STATE_NF90 - unable to find variable: ',a,2x, &
1084 & ' in ',a,' NetCDF file: ',a)
1085!
1086 RETURN
subroutine, public netcdf_close(ng, model, ncid, ncname, lupdate)
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)

References mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::idfsur, mod_ncparam::idsbry, mod_sediment::idsed, mod_ncparam::idtime, mod_ncparam::idtsur, mod_ncparam::idtvar, mod_ncparam::idubar, mod_ncparam::idusms, mod_ncparam::iduvel, mod_ncparam::idvbar, mod_ncparam::idvsms, mod_ncparam::idvvel, mod_ncparam::iinfo, mod_param::inlm, mod_param::iobounds, mod_scalars::isalt, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::itemp, mod_scalars::lobc, mod_scalars::lstflux, mod_parallel::master, mod_param::n, mod_netcdf::n_var, mod_biology::nbac, mod_biology::nbands, mod_param::nbed, mod_scalars::nbrec, mod_biology::ndom, mod_netcdf::netcdf_check_dim(), mod_netcdf::netcdf_close(), mod_netcdf::netcdf_create(), mod_netcdf::netcdf_enddef(), mod_netcdf::netcdf_inq_var(), mod_netcdf::netcdf_open(), mod_netcdf::nf_fout, mod_netcdf::nf_tout, mod_biology::nfec, mod_scalars::nfrec, mod_scalars::noerror, mod_biology::nphy, mod_param::nst, mod_fourdvar::nstatevar, mod_param::nt, mod_ncparam::nv, mod_param::r2dvar, mod_param::r3dvar, mod_scalars::rclock, mod_netcdf::rec_size, mod_sediment::sd50, mod_iounits::sourcefile, mod_iounits::stdout, mod_param::u2dvar, mod_param::v2dvar, mod_netcdf::var_id, mod_netcdf::var_name, and mod_ncparam::vname.

Referenced by def_state().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ def_state_pio()

subroutine, private def_state_mod::def_state_pio ( integer, intent(in) ng,
integer, intent(in) model,
character (len=*), intent(in) label,
logical, intent(in) lcreate,
type(t_io), dimension(ngrids), intent(inout) s,
logical, intent(in) lclose )
private

Definition at line 1092 of file def_state.F.

1093!***********************************************************************
1094!
1095 USE mod_pio_netcdf
1096!
1097! Imported variable declarations.
1098!
1099 logical, intent(in) :: Lcreate, Lclose
1100!
1101 integer, intent(in) :: ng, model
1102!
1103 character (len=*), intent(in) :: label
1104!
1105 TYPE(T_IO), intent(inout) :: S(Ngrids)
1106!
1107! Local variable declarations.
1108!
1109 logical :: got_var(NV)
1110!
1111 integer, parameter :: Natt = 25
1112
1113 integer :: i, j, ifield, itrc, nrec, nvd, nvd3, nvd4
1114 integer :: recdim, status, varid
1115# ifdef ADJUST_BOUNDARY
1116 integer :: IorJdim, brecdim
1117# endif
1118# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1119 integer :: frecdim
1120# endif
1121 integer :: DimIDs(nDimID)
1122 integer :: t2dgrd(3), u2dgrd(3), v2dgrd(3)
1123# ifdef ADJUST_BOUNDARY
1124 integer :: t2dobc(4)
1125# endif
1126# ifdef SOLVE3D
1127 integer :: t3dgrd(4), u3dgrd(4), v3dgrd(4), w3dgrd(4)
1128# ifdef ADJUST_BOUNDARY
1129 integer :: t3dobc(5)
1130# endif
1131# ifdef ADJUST_STFLUX
1132 integer :: t3dfrc(4)
1133# endif
1134# endif
1135# ifdef ADJUST_WSTRESS
1136 integer :: u3dfrc(4), v3dfrc(4)
1137# endif
1138!
1139 real(r8) :: Aval(6)
1140!
1141 character (len=256) :: ncname
1142 character (len=MaxLen) :: Vinfo(Natt)
1143
1144 character (len=*), parameter :: MyFile = &
1145 & __FILE__//", def_state_pio"
1146!
1147 sourcefile=myfile
1148!
1149!-----------------------------------------------------------------------
1150! Set and report file name.
1151!-----------------------------------------------------------------------
1152!
1153 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1154 ncname=s(ng)%name
1155!
1156 IF (master) THEN
1157 IF (lcreate) THEN
1158 WRITE (stdout,10) trim(label), ng, trim(ncname)
1159 ELSE
1160 WRITE (stdout,20) trim(label), ng, trim(ncname)
1161 END IF
1162 END IF
1163!
1164!=======================================================================
1165! Create a new ocean state NetCDF file.
1166!=======================================================================
1167!
1168 define : IF (lcreate) THEN
1169 CALL pio_netcdf_create (ng, model, trim(ncname), s(ng)%pioFile)
1170 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1171 IF (master) WRITE (stdout,30) trim(label), trim(ncname)
1172 RETURN
1173 END IF
1174!
1175!-----------------------------------------------------------------------
1176! Define file dimensions.
1177!-----------------------------------------------------------------------
1178!
1179 dimids=0
1180!
1181 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xi_rho', &
1182 & iobounds(ng)%xi_rho, dimids( 1))
1183 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1184
1185 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xi_u', &
1186 & iobounds(ng)%xi_u, dimids( 2))
1187 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1188
1189 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xi_v', &
1190 & iobounds(ng)%xi_v, dimids( 3))
1191 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1192
1193 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xi_psi', &
1194 & iobounds(ng)%xi_psi, dimids( 4))
1195 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1196
1197 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'eta_rho', &
1198 & iobounds(ng)%eta_rho, dimids( 5))
1199 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1200
1201 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'eta_u', &
1202 & iobounds(ng)%eta_u, dimids( 6))
1203 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1204
1205 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'eta_v', &
1206 & iobounds(ng)%eta_v, dimids( 7))
1207 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1208
1209 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'eta_psi', &
1210 & iobounds(ng)%eta_psi, dimids( 8))
1211 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1212
1213# ifdef ADJUST_BOUNDARY
1214 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'IorJ', &
1215 & iobounds(ng)%IorJ, iorjdim)
1216 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1217# endif
1218
1219# if defined WRITE_WATER && defined MASKING
1220 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xy_rho', &
1221 & iobounds(ng)%xy_rho, dimids(17))
1222 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1223
1224 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xy_u', &
1225 & iobounds(ng)%xy_u, dimids(18))
1226 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1227
1228 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xy_v', &
1229 & iobounds(ng)%xy_v, dimids(19))
1230 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1231# endif
1232# ifdef SOLVE3D
1233# if defined WRITE_WATER && defined MASKING
1234 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xyz_rho', &
1235 & iobounds(ng)%xy_rho*n(ng), dimids(20))
1236 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1237
1238 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xyz_u', &
1239 & iobounds(ng)%xy_u*n(ng), dimids(21))
1240 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1241
1242 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xyz_v', &
1243 & iobounds(ng)%xy_v*n(ng), dimids(22))
1244 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1245
1246 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xyz_w', &
1247 & iobounds(ng)%xy_rho*(n(ng)+1), dimids(23))
1248 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1249# endif
1250
1251 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'N', &
1252 & n(ng), dimids( 9))
1253 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1254
1255 status=def_dim(ng, model, s(ng)%pioFile, ncname, 's_rho', &
1256 & n(ng), dimids( 9))
1257 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1258
1259 status=def_dim(ng, model, s(ng)%pioFile, ncname, 's_w', &
1260 & n(ng)+1, dimids(10))
1261 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1262
1263 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'tracer', &
1264 & nt(ng), dimids(11))
1265 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1266
1267# ifdef SEDIMENT
1268 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'NST', &
1269 & nst, dimids(32))
1270 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1271
1272 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'Nbed', &
1273 & nbed, dimids(16))
1274 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1275
1276# if defined WRITE_WATER && defined MASKING
1277 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'xybed', &
1278 & iobounds(ng)%xy_rho*nbed, dimids(24))
1279 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1280# endif
1281# endif
1282
1283# ifdef ECOSIM
1284 status=def_dim(ng, inlm, s(ng)%pioFile, ncname, 'Nbands', &
1285 & nbands, dimids(33))
1286 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1287
1288 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'Nphy', &
1289 & nphy, dimids(25))
1290 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1291
1292 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'Nbac', &
1293 & nbac, dimids(26))
1294 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1295
1296 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'Ndom', &
1297 & ndom, dimids(27))
1298 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1299
1300 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'Nfec', &
1301 & nfec, dimids(28))
1302 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1303# endif
1304# endif
1305
1306 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'boundary', &
1307 & 4, dimids(14))
1308 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1309
1310# ifdef FOUR_DVAR
1311 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'Nstate', &
1312 & nstatevar(ng), dimids(29))
1313 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1314# endif
1315
1316# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1317 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'frc_adjust', &
1318 & nfrec(ng), dimids(30))
1319 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1320# endif
1321
1322# ifdef ADJUST_BOUNDARY
1323 status=def_dim(ng, model, s(ng)%pioFile, ncname, 'obc_adjust', &
1324 & nbrec(ng), dimids(31))
1325 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1326# endif
1327
1328 status=def_dim(ng, model, s(ng)%pioFile, ncname, &
1329 & trim(adjustl(vname(5,idtime))), &
1330 & pio_unlimited, dimids(12))
1331 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1332
1333 recdim=dimids(12)
1334# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1335 frecdim=dimids(30)
1336# endif
1337# ifdef ADJUST_BOUNDARY
1338 brecdim=dimids(31)
1339# endif
1340!
1341! Set number of dimensions for output variables.
1342!
1343# if defined WRITE_WATER && defined MASKING
1344 nvd3=2
1345 nvd4=2
1346# else
1347 nvd3=3
1348 nvd4=4
1349# endif
1350!
1351! Define dimension vectors for staggered tracer type variables.
1352!
1353# if defined WRITE_WATER && defined MASKING
1354 t2dgrd(1)=dimids(17)
1355 t2dgrd(2)=dimids(12)
1356# ifdef SOLVE3D
1357 t3dgrd(1)=dimids(20)
1358 t3dgrd(2)=dimids(12)
1359# endif
1360# else
1361 t2dgrd(1)=dimids( 1)
1362 t2dgrd(2)=dimids( 5)
1363 t2dgrd(3)=dimids(12)
1364# ifdef SOLVE3D
1365 t3dgrd(1)=dimids( 1)
1366 t3dgrd(2)=dimids( 5)
1367 t3dgrd(3)=dimids( 9)
1368 t3dgrd(4)=dimids(12)
1369# endif
1370# ifdef ADJUST_STFLUX
1371 t3dfrc(1)=dimids( 1)
1372 t3dfrc(2)=dimids( 5)
1373 t3dfrc(3)=frecdim
1374 t3dfrc(4)=dimids(12)
1375# endif
1376# endif
1377# ifdef ADJUST_BOUNDARY
1378 t2dobc(1)=iorjdim
1379 t2dobc(2)=dimids(14)
1380 t2dobc(3)=brecdim
1381 t2dobc(4)=dimids(12)
1382# ifdef SOLVE3D
1383 t3dobc(1)=iorjdim
1384 t3dobc(2)=dimids( 9)
1385 t3dobc(3)=dimids(14)
1386 t3dobc(4)=brecdim
1387 t3dobc(5)=dimids(12)
1388# endif
1389# endif
1390!
1391! Define dimension vectors for staggered u-momentum type variables.
1392!
1393# if defined WRITE_WATER && defined MASKING
1394 u2dgrd(1)=dimids(18)
1395 u2dgrd(2)=dimids(12)
1396# ifdef SOLVE3D
1397 u3dgrd(1)=dimids(21)
1398 u3dgrd(2)=dimids(12)
1399# endif
1400# else
1401 u2dgrd(1)=dimids( 2)
1402 u2dgrd(2)=dimids( 6)
1403 u2dgrd(3)=dimids(12)
1404# ifdef SOLVE3D
1405 u3dgrd(1)=dimids( 2)
1406 u3dgrd(2)=dimids( 6)
1407 u3dgrd(3)=dimids( 9)
1408 u3dgrd(4)=dimids(12)
1409# endif
1410# ifdef ADJUST_WSTRESS
1411 u3dfrc(1)=dimids( 2)
1412 u3dfrc(2)=dimids( 6)
1413 u3dfrc(3)=frecdim
1414 u3dfrc(4)=dimids(12)
1415# endif
1416# endif
1417!
1418! Define dimension vectors for staggered v-momentum type variables.
1419!
1420# if defined WRITE_WATER && defined MASKING
1421 v2dgrd(1)=dimids(19)
1422 v2dgrd(2)=dimids(12)
1423# ifdef SOLVE3D
1424 v3dgrd(1)=dimids(22)
1425 v3dgrd(2)=dimids(12)
1426# endif
1427# else
1428 v2dgrd(1)=dimids( 3)
1429 v2dgrd(2)=dimids( 7)
1430 v2dgrd(3)=dimids(12)
1431# ifdef SOLVE3D
1432 v3dgrd(1)=dimids( 3)
1433 v3dgrd(2)=dimids( 7)
1434 v3dgrd(3)=dimids( 9)
1435 v3dgrd(4)=dimids(12)
1436# endif
1437# ifdef ADJUST_WSTRESS
1438 v3dfrc(1)=dimids( 3)
1439 v3dfrc(2)=dimids( 7)
1440 v3dfrc(3)=frecdim
1441 v3dfrc(4)=dimids(12)
1442# endif
1443# endif
1444# ifdef SOLVE3D
1445!
1446! Define dimension vector for staggered w-momentum type variables.
1447!
1448# if defined WRITE_WATER && defined MASKING
1449 w3dgrd(1)=dimids(23)
1450 w3dgrd(2)=dimids(12)
1451# else
1452 w3dgrd(1)=dimids( 1)
1453 w3dgrd(2)=dimids( 5)
1454 w3dgrd(3)=dimids(10)
1455 w3dgrd(4)=dimids(12)
1456# endif
1457# endif
1458!
1459! Initialize unlimited time record dimension.
1460!
1461 s(ng)%Rindex=0
1462!
1463! Initialize local information variable arrays.
1464!
1465 DO i=1,natt
1466 DO j=1,len(vinfo(1))
1467 vinfo(i)(j:j)=' '
1468 END DO
1469 END DO
1470 DO i=1,6
1471 aval(i)=0.0_r8
1472 END DO
1473!
1474!-----------------------------------------------------------------------
1475! Define time-recordless information variables.
1476!-----------------------------------------------------------------------
1477!
1478 CALL def_info (ng, model, s(ng)%pioFile, ncname, dimids)
1479 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1480!
1481!-----------------------------------------------------------------------
1482! Define time-varying variables.
1483!-----------------------------------------------------------------------
1484!
1485! Define model time.
1486!
1487 vinfo( 1)=vname(1,idtime)
1488 vinfo( 2)=vname(2,idtime)
1489 WRITE (vinfo( 3),'(a,a)') 'seconds since ', trim(rclock%string)
1490 vinfo( 4)=trim(rclock%calendar)
1491 vinfo(14)=vname(4,idtime)
1492 vinfo(21)=vname(6,idtime)
1493 s(ng)%pioVar(idtime)%dkind=pio_tout
1494 s(ng)%pioVar(idtime)%gtype=0
1495!
1496 status=def_var(ng, model, s(ng)%pioFile, &
1497 & s(ng)%pioVar(idtime)%vd, &
1498 & pio_tout, 1, (/recdim/), aval, vinfo,ncname, &
1499 & setparaccess = .true.)
1500 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1501!
1502! Define free-surface.
1503!
1504 vinfo( 1)=vname(1,idfsur)
1505 vinfo( 2)=vname(2,idfsur)
1506 vinfo( 3)='nondimensional'
1507 vinfo(14)=vname(4,idfsur)
1508 vinfo(16)=vname(1,idtime)
1509# if defined WRITE_WATER && defined MASKING
1510 vinfo(20)='mask_rho'
1511# endif
1512 vinfo(21)=vname(6,idfsur)
1513 vinfo(22)='coordinates'
1514 aval(5)=real(iinfo(1,idfsur,ng),r8)
1515 s(ng)%pioVar(idfsur)%dkind=pio_fout
1516 s(ng)%pioVar(idfsur)%gtype=r2dvar
1517!
1518 status=def_var(ng, model, s(ng)%pioFile, &
1519 & s(ng)%pioVar(idfsur)%vd, &
1520 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1521 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1522
1523# ifdef ADJUST_BOUNDARY
1524!
1525! Define free-surface open boundaries.
1526!
1527 IF (any(lobc(:,isfsur,ng))) THEN
1528 ifield=idsbry(isfsur)
1529 vinfo( 1)=vname(1,ifield)
1530 vinfo( 2)=vname(2,ifield)
1531 vinfo( 3)='nondimensional'
1532 vinfo(14)=vname(4,ifield)
1533 vinfo(16)=vname(1,idtime)
1534 vinfo(21)=vname(6,ifield)
1535 aval(5)=real(iinfo(1,ifield,ng),r8)
1536 s(ng)%pioVar(ifield)%dkind=pio_fout
1537 s(ng)%pioVar(ifield)%gtype=r2dobc
1538!
1539 status=def_var(ng, model, s(ng)%pioFile, &
1540 & s(ng)%pioVar(ifield)%vd, &
1541 & pio_fout, 4, t2dobc, aval, vinfo, ncname, &
1542 & setfillval = .false.)
1543 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1544 END IF
1545# endif
1546!
1547! Define 2D U-momentum component.
1548!
1549 vinfo( 1)=vname(1,idubar)
1550 vinfo( 2)=vname(2,idubar)
1551 vinfo( 3)='nondimensional'
1552 vinfo(14)=vname(4,idubar)
1553 vinfo(16)=vname(1,idtime)
1554# if defined WRITE_WATER && defined MASKING
1555 vinfo(20)='mask_u'
1556# endif
1557 vinfo(21)=vname(6,idubar)
1558 vinfo(22)='coordinates'
1559 aval(5)=real(iinfo(1,idubar,ng),r8)
1560 s(ng)%pioVar(idubar)%dkind=pio_fout
1561 s(ng)%pioVar(idubar)%gtype=u2dvar
1562!
1563 status=def_var(ng, model, s(ng)%pioFile, &
1564 & s(ng)%pioVar(idubar)%vd, &
1565 & pio_fout, nvd3, u2dgrd, aval, vinfo, ncname)
1566 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1567
1568# ifdef ADJUST_BOUNDARY
1569!
1570! Define 2D U-momentum component open boundaries.
1571!
1572 IF (any(lobc(:,isubar,ng))) THEN
1573 ifield=idsbry(isubar)
1574 vinfo( 1)=vname(1,ifield)
1575 vinfo( 2)=vname(2,ifield)
1576 vinfo( 3)='nondimensional'
1577 vinfo(14)=vname(4,ifield)
1578 vinfo(16)=vname(1,idtime)
1579 vinfo(21)=vname(6,ifield)
1580 aval(5)=real(iinfo(1,ifield,ng),r8)
1581 s(ng)%pioVar(ifield)%dkind=pio_fout
1582 s(ng)%pioVar(ifield)%gtype=u2dobc
1583!
1584 status=def_var(ng, model, s(ng)%pioFile, &
1585 & s(ng)%pioVar(ifield)%vd, &
1586 & pio_fout, 4, t2dobc, aval, vinfo, ncname, &
1587 & setfillval = .false.)
1588 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1589 END IF
1590# endif
1591!
1592! Define 2D V-momentum component.
1593!
1594 vinfo( 1)=vname(1,idvbar)
1595 vinfo( 2)=vname(2,idvbar)
1596 vinfo( 3)='nondimensional'
1597 vinfo(14)=vname(4,idvbar)
1598 vinfo(16)=vname(1,idtime)
1599# if defined WRITE_WATER && defined MASKING
1600 vinfo(20)='mask_v'
1601# endif
1602 vinfo(21)=vname(6,idvbar)
1603 vinfo(22)='coordinates'
1604 aval(5)=real(iinfo(1,idvbar,ng),r8)
1605 s(ng)%pioVar(idvbar)%dkind=pio_fout
1606 s(ng)%pioVar(idvbar)%gtype=v2dvar
1607!
1608 status=def_var(ng, model, s(ng)%pioFile, &
1609 & s(ng)%pioVar(idvbar)%vd, &
1610 & pio_fout, nvd3, v2dgrd, aval, vinfo, ncname)
1611 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1612
1613# ifdef ADJUST_BOUNDARY
1614!
1615! Define 2D V-momentum component open boundaries.
1616!
1617 IF (any(lobc(:,isvbar,ng))) THEN
1618 ifield=idsbry(isvbar)
1619 vinfo( 1)=vname(1,ifield)
1620 vinfo( 2)=vname(2,ifield)
1621 vinfo( 3)='nondimensional'
1622 vinfo(14)=vname(4,ifield)
1623 vinfo(16)=vname(1,idtime)
1624 vinfo(21)=vname(6,ifield)
1625 aval(5)=real(iinfo(1,ifield,ng),r8)
1626 s(ng)%pioVar(ifield)%dkind=pio_fout
1627 s(ng)%pioVar(ifield)%gtype=v2dobc
1628!
1629 status=def_var(ng, model, s(ng)%pioFile, &
1630 & s(ng)%pioVar(ifield)%vd, &
1631 & pio_fout, 4, t2dobc, aval, vinfo, ncname, &
1632 & setfillval = .false.)
1633 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1634 END IF
1635# endif
1636# ifdef SOLVE3D
1637!
1638! Define 3D U-momentum component.
1639!
1640 vinfo( 1)=vname(1,iduvel)
1641 vinfo( 2)=vname(2,iduvel)
1642 vinfo( 3)='nondimensional'
1643 vinfo(14)=vname(4,iduvel)
1644 vinfo(16)=vname(1,idtime)
1645# if defined WRITE_WATER && defined MASKING
1646 vinfo(20)='mask_u'
1647# endif
1648 vinfo(21)=vname(6,iduvel)
1649 vinfo(22)='coordinates'
1650 aval(5)=real(iinfo(1,iduvel,ng),r8)
1651 s(ng)%pioVar(iduvel)%dkind=pio_fout
1652 s(ng)%pioVar(iduvel)%gtype=u3dvar
1653!
1654 status=def_var(ng, model, s(ng)%pioFile, &
1655 & s(ng)%pioVar(iduvel)%vd, &
1656 & pio_fout, nvd4, u3dgrd, aval, vinfo, ncname)
1657 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1658
1659# ifdef ADJUST_BOUNDARY
1660!
1661! Define 3D U-momentum component open boundaries.
1662!
1663 IF (any(lobc(:,isuvel,ng))) THEN
1664 ifield=idsbry(isuvel)
1665 vinfo( 1)=vname(1,ifield)
1666 vinfo( 2)=vname(2,ifield)
1667 vinfo( 3)='nondimensional'
1668 vinfo(14)=vname(4,ifield)
1669 vinfo(16)=vname(1,idtime)
1670 vinfo(21)=vname(6,ifield)
1671 aval(5)=real(iinfo(1,ifield,ng),r8)
1672 s(ng)%pioVar(ifield)%dkind=pio_fout
1673 s(ng)%pioVar(ifield)%gtype=u3dobc
1674!
1675 status=def_var(ng, model, s(ng)%pioFile, &
1676 & s(ng)%pioVar(ifield)%vd, &
1677 & pio_fout, 5, t3dobc, aval, vinfo, ncname, &
1678 & setfillval = .false.)
1679 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1680 END IF
1681# endif
1682!
1683! Define 3D V-momentum component.
1684!
1685 vinfo( 1)=vname(1,idvvel)
1686 vinfo( 2)=vname(2,idvvel)
1687 vinfo( 3)='nondimensional'
1688 vinfo(14)=vname(4,idvvel)
1689 vinfo(16)=vname(1,idtime)
1690# if defined WRITE_WATER && defined MASKING
1691 vinfo(20)='mask_v'
1692# endif
1693 vinfo(21)=vname(6,idvvel)
1694 vinfo(22)='coordinates'
1695 aval(5)=real(iinfo(1,idvvel,ng),r8)
1696 s(ng)%pioVar(idvvel)%dkind=pio_fout
1697 s(ng)%pioVar(idvvel)%gtype=v3dvar
1698!
1699 status=def_var(ng, model, s(ng)%pioFile, &
1700 & s(ng)%pioVar(idvvel)%vd, &
1701 & pio_fout, nvd4, v3dgrd, aval, vinfo, ncname)
1702 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1703
1704# ifdef ADJUST_BOUNDARY
1705!
1706! Define 3D V-momentum component open boundaries.
1707!
1708 IF (any(lobc(:,isvvel,ng))) THEN
1709 ifield=idsbry(isvvel)
1710 vinfo( 1)=vname(1,ifield)
1711 vinfo( 2)=vname(2,ifield)
1712 vinfo( 3)='nondimensional'
1713 vinfo(14)=vname(4,ifield)
1714 vinfo(16)=vname(1,idtime)
1715 vinfo(21)=vname(6,ifield)
1716 aval(5)=real(iinfo(1,ifield,ng),r8)
1717 s(ng)%pioVar(ifield)%dkind=pio_fout
1718 s(ng)%pioVar(ifield)%gtype=v3dobc
1719!
1720 status=def_var(ng, model, s(ng)%pioFile, &
1721 & s(ng)%pioVar(ifield)%vd, &
1722 & pio_fout, 5, t3dobc, aval, vinfo, ncname, &
1723 & setfillval = .false.)
1724 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1725 END IF
1726# endif
1727!
1728! Define tracer type variables.
1729!
1730 DO itrc=1,nt(ng)
1731 vinfo( 1)=vname(1,idtvar(itrc))
1732 vinfo( 2)=vname(2,idtvar(itrc))
1733 vinfo( 3)='nondimensional'
1734 vinfo(14)=vname(4,idtvar(itrc))
1735 vinfo(16)=vname(1,idtime)
1736# ifdef SEDIMENT
1737 DO i=1,nst
1738 IF (itrc.eq.idsed(i)) THEN
1739 WRITE (vinfo(19),40) 1000.0_r8*sd50(i,ng)
1740 END IF
1741 END DO
1742# endif
1743# if defined WRITE_WATER && defined MASKING
1744 vinfo(20)='mask_rho'
1745# endif
1746 vinfo(21)=vname(6,idtvar(itrc))
1747 vinfo(22)='coordinates'
1748 aval(5)=real(r3dvar,r8)
1749 s(ng)%pioTrc(itrc)%dkind=pio_fout
1750 s(ng)%pioTrc(itrc)%gtype=r3dvar
1751!
1752 status=def_var(ng, model, s(ng)%pioFile, &
1753 & s(ng)%pioTrc(itrc)%vd, &
1754 & pio_fout, nvd4, t3dgrd, aval, vinfo, ncname)
1755 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1756 END DO
1757
1758# ifdef ADJUST_BOUNDARY
1759!
1760! Define tracer type variables open boundaries.
1761!
1762 DO itrc=1,nt(ng)
1763 IF (any(lobc(:,istvar(itrc),ng))) THEN
1764 ifield=idsbry(istvar(itrc))
1765 vinfo( 1)=vname(1,ifield)
1766 vinfo( 2)=vname(2,ifield)
1767 vinfo( 3)='nondimensional'
1768 vinfo(14)=vname(4,ifield)
1769 vinfo(16)=vname(1,idtime)
1770# ifdef SEDIMENT
1771 DO i=1,nst
1772 IF (itrc.eq.idsed(i)) THEN
1773 WRITE (vinfo(19),40) 1000.0_r8*sd50(i,ng)
1774 END IF
1775 END DO
1776# endif
1777 vinfo(21)=vname(6,ifield)
1778 aval(5)=real(iinfo(1,ifield,ng),r8)
1779 s(ng)%pioVar(ifield)%dkind=pio_fout
1780 s(ng)%pioVar(ifield)%gtype=r3dobc
1781!
1782 status=def_var(ng, model, s(ng)%pioFile, &
1783 & s(ng)%pioVar(ifield)%vd, &
1784 & pio_fout, 5, t3dobc, aval, vinfo, ncname, &
1785 & setfillval = .false.)
1786 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1787 END IF
1788 END DO
1789# endif
1790# ifdef ADJUST_STFLUX
1791!
1792! Define surface tracer fluxes.
1793!
1794 DO itrc=1,nt(ng)
1795 IF (lstflux(itrc,ng)) THEN
1796 vinfo( 1)=vname(1,idtsur(itrc))
1797 vinfo( 2)=vname(2,idtsur(itrc))
1798 vinfo( 3)='nondimensional'
1799 IF (itrc.eq.itemp) THEN
1800 vinfo(11)='upward flux, cooling'
1801 vinfo(12)='downward flux, heating'
1802 ELSE IF (itrc.eq.isalt) THEN
1803 vinfo(11)='upward flux, freshening (net precipitation)'
1804 vinfo(12)='downward flux, salting (net evaporation)'
1805 END IF
1806 vinfo(14)=vname(4,idtsur(itrc))
1807 vinfo(16)=vname(1,idtime)
1808# if defined WRITE_WATER && defined MASKING
1809 vinfo(20)='mask_rho'
1810# endif
1811 vinfo(21)=vname(6,idtsur(itrc))
1812 vinfo(22)='coordinates'
1813 aval(5)=real(r2dvar,r8)
1814 s(ng)%pioVar(idtsur(itrc))%dkind=pio_fout
1815 s(ng)%pioVar(idtsur(itrc))%gtype=r2dvar
1816!
1817 status=def_var(ng, model, s(ng)%pioFile, &
1818 & s(ng)%pioVar(idtsur(itrc))%vd, &
1819 & pio_fout, nvd4, t3dfrc, aval, vinfo, ncname)
1820 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1821 END IF
1822 END DO
1823# endif
1824# endif
1825# ifdef ADJUST_WSTRESS
1826!
1827! Define surface U-momentum stress.
1828!
1829 vinfo( 1)=vname(1,idusms)
1830 vinfo( 2)=vname(2,idusms)
1831 vinfo( 3)='nondimensional'
1832 vinfo(14)=vname(4,idusms)
1833 vinfo(16)=vname(1,idtime)
1834# if defined WRITE_WATER && defined MASKING
1835 vinfo(20)='mask_u'
1836# endif
1837 vinfo(21)=vname(6,idusms)
1838 vinfo(22)='coordinates'
1839 aval(5)=real(u2dvar,r8)
1840 s(ng)%pioVar(idusms)%dkind=pio_fout
1841 s(ng)%pioVar(idusms)%gtype=u2dvar
1842!
1843 status=def_var(ng, model, s(ng)%pioFile, &
1844 & s(ng)%pioVar(idusms)%vd, &
1845 & pio_fout, nvd4, u3dfrc, aval, vinfo, ncname)
1846 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1847!
1848! Define surface V-momentum stress.
1849!
1850 vinfo( 1)=vname(1,idvsms)
1851 vinfo( 2)=vname(2,idvsms)
1852 vinfo( 3)='nondimensional'
1853 vinfo(14)=vname(4,idvsms)
1854 vinfo(16)=vname(1,idtime)
1855# if defined WRITE_WATER && defined MASKING
1856 vinfo(20)='mask_v'
1857# endif
1858 vinfo(21)=vname(6,idvsms)
1859 vinfo(22)='coordinates'
1860 aval(5)=real(v2dvar,r8)
1861 s(ng)%pioVar(idvsms)%dkind=pio_fout
1862 s(ng)%pioVar(idvsms)%gtype=v2dvar
1863!
1864 status=def_var(ng, model, s(ng)%pioFile, &
1865 & s(ng)%pioVar(idvsms)%vd, &
1866 & pio_fout, nvd4, v3dfrc, aval, vinfo, ncname)
1867 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1868# endif
1869!
1870!-----------------------------------------------------------------------
1871! Leave definition mode.
1872!-----------------------------------------------------------------------
1873!
1874 CALL pio_netcdf_enddef (ng, model, ncname, s(ng)%pioFile)
1875 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1876!
1877!-----------------------------------------------------------------------
1878! Write out time-recordless, information variables.
1879!-----------------------------------------------------------------------
1880!
1881 CALL wrt_info (ng, model, s(ng)%pioFile, ncname)
1882 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1883!
1884!-----------------------------------------------------------------------
1885! Write out time-recordless, information variables.
1886!-----------------------------------------------------------------------
1887!
1888 IF (lclose) THEN
1889 CALL pio_netcdf_close (ng, model, s(ng)%pioFile, ncname, &
1890 & .false.)
1891 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1892 END IF
1893
1894 END IF define
1895!
1896!=======================================================================
1897! Open an existing ocean state NetCDF file, check its contents, and
1898! prepare for appending data.
1899!=======================================================================
1900!
1901 query: IF (.not.lcreate) THEN
1902 ncname=s(ng)%name
1903!
1904! Open ocean state file for read/write.
1905!
1906 CALL pio_netcdf_open (ng, model, ncname, 1, s(ng)%pioFile)
1907 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1908 WRITE (stdout,50) trim(label), trim(ncname)
1909 RETURN
1910 END IF
1911!
1912! Inquire about the dimensions and check for consistency.
1913!
1914 CALL pio_netcdf_check_dim (ng, model, ncname, &
1915 & piofile = s(ng)%pioFile)
1916 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1917!
1918! Inquire about the variables.
1919!
1920 CALL pio_netcdf_inq_var (ng, model, ncname, &
1921 & piofile = s(ng)%pioFile)
1922 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1923!
1924! Initialize logical switches.
1925!
1926 DO i=1,nv
1927 got_var(i)=.false.
1928 END DO
1929!
1930! Scan variable list from input NetCDF and activate switches for
1931! Hessian eigenvectors variables. Get variable IDs.
1932!
1933 DO i=1,n_var
1934 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
1935 got_var(idtime)=.true.
1936 s(ng)%pioVar(idtime)%vd=var_desc(i)
1937 s(ng)%pioVar(idtime)%dkind=pio_tout
1938 s(ng)%pioVar(idtime)%gtype=0
1939 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
1940 got_var(idfsur)=.true.
1941 s(ng)%pioVar(idfsur)%vd=var_desc(i)
1942 s(ng)%pioVar(idfsur)%dkind=pio_tout
1943 s(ng)%pioVar(idfsur)%gtype=r2dvar
1944 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
1945 got_var(idubar)=.true.
1946 s(ng)%pioVar(idubar)%vd=var_desc(i)
1947 s(ng)%pioVar(idubar)%dkind=pio_fout
1948 s(ng)%pioVar(idubar)%gtype=u2dvar
1949 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
1950 got_var(idvbar)=.true.
1951 s(ng)%pioVar(idvbar)%vd=var_desc(i)
1952 s(ng)%pioVar(idvbar)%dkind=pio_fout
1953 s(ng)%pioVar(idvbar)%gtype=v2dvar
1954# ifdef ADJUST_BOUNDARY
1955 ELSE IF (trim(var_name(i)).eq. &
1956 & trim(vname(1,idsbry(isfsur)))) THEN
1957 got_var(idsbry(isfsur))=.true.
1958 s(ng)%pioVar(idsbry(isfsur))%vd=var_desc(i)
1959 s(ng)%pioVar(idsbry(isfsur))%dkind=pio_fout
1960 s(ng)%pioVar(idsbry(isfsur))%gtype=r2dobc
1961 ELSE IF (trim(var_name(i)).eq. &
1962 & trim(vname(1,idsbry(isubar)))) THEN
1963 got_var(idsbry(isubar))=.true.
1964 s(ng)%pioVar(idsbry(isubar))%vd=var_desc(i)
1965 s(ng)%pioVar(idsbry(isubar))%dkind=pio_fout
1966 s(ng)%pioVar(idsbry(isubar))%gtype=u2dobc
1967 ELSE IF (trim(var_name(i)).eq. &
1968 & trim(vname(1,idsbry(isvbar)))) THEN
1969 got_var(idsbry(isvbar))=.true.
1970 s(ng)%pioVar(idsbry(isvbar))%vd=var_desc(i)
1971 s(ng)%pioVar(idsbry(isvbar))%dkind=pio_fout
1972 s(ng)%pioVar(idsbry(isvbar))%gtype=v2dobc
1973# endif
1974# ifdef ADJUST_WSTRESS
1975 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idusms))) THEN
1976 got_var(idusms)=.true.
1977 s(ng)%pioVar(idusms)%vd=var_desc(i)
1978 s(ng)%pioVar(idusms)%dkind=pio_fout
1979 s(ng)%pioVar(idusms)%gtype=u2dvar
1980 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvsms))) THEN
1981 got_var(idvsms)=.true.
1982 s(ng)%pioVar(idvsms)%vd=var_desc(i)
1983 s(ng)%pioVar(idvsms)%dkind=pio_fout
1984 s(ng)%pioVar(idvsms)%gtype=v2dvar
1985# endif
1986# ifdef SOLVE3D
1987 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
1988 got_var(iduvel)=.true.
1989 s(ng)%pioVar(iduvel)%vd=var_desc(i)
1990 s(ng)%pioVar(iduvel)%dkind=pio_fout
1991 s(ng)%pioVar(iduvel)%gtype=u3dvar
1992 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
1993 got_var(idvvel)=.true.
1994 s(ng)%pioVar(idvvel)%vd=var_desc(i)
1995 s(ng)%pioVar(idvvel)%dkind=pio_fout
1996 s(ng)%pioVar(idvvel)%gtype=v3dvar
1997# ifdef ADJUST_BOUNDARY
1998 ELSE IF (trim(var_name(i)).eq. &
1999 & trim(vname(1,idsbry(isuvel)))) THEN
2000 got_var(idsbry(isuvel))=.true.
2001 s(ng)%pioVar(idsbry(isuvel))%vd=var_desc(i)
2002 s(ng)%pioVar(idsbry(isuvel))%dkind=pio_fout
2003 s(ng)%pioVar(idsbry(isuvel))%gtype=u3dobc
2004 ELSE IF (trim(var_name(i)).eq. &
2005 & trim(vname(1,idsbry(isvvel)))) THEN
2006 got_var(idsbry(isvvel))=.true.
2007 s(ng)%pioVar(idsbry(isvvel))%vd=var_desc(i)
2008 s(ng)%pioVar(idsbry(isvvel))%dkind=pio_fout
2009 s(ng)%pioVar(idsbry(isvvel))%gtype=v3dobc
2010# endif
2011# endif
2012 END IF
2013# ifdef SOLVE3D
2014 DO itrc=1,nt(ng)
2015 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
2016 got_var(idtvar(itrc))=.true.
2017 s(ng)%pioTrc(itrc)%vd=var_desc(i)
2018 s(ng)%pioTrc(itrc)%dkind=pio_fout
2019 s(ng)%pioTrc(itrc)%gtype=r3dvar
2020# ifdef ADJUST_BOUNDARY
2021 ELSE IF (trim(var_name(i)).eq. &
2022 & trim(vname(1,idsbry(istvar(itrc))))) THEN
2023 got_var(idsbry(istvar(itrc)))=.true.
2024 s(ng)%pioVar(idsbry(istvar(itrc)))%vd=var_desc(i)
2025 s(ng)%pioVar(idsbry(istvar(itrc)))%dkind=pio_fout
2026 s(ng)%pioVar(idsbry(istvar(itrc)))%gtype=r3dobc
2027# endif
2028# ifdef ADJUST_STFLUX
2029 ELSE IF (trim(var_name(i)).eq. &
2030 & trim(vname(1,idtsur(itrc)))) THEN
2031 got_var(idtsur(itrc))=.true.
2032 s(ng)%pioVar(idtsur(itrc))%vd=var_desc(i)
2033 s(ng)%pioVar(idtsur(itrc))%dkind=pio_fout
2034 s(ng)%pioVar(idtsur(itrc))%gtype=r2dvar
2035# endif
2036 END IF
2037 END DO
2038# endif
2039 END DO
2040!
2041! Check if ocean state variables are available in input NetCDF file.
2042!
2043 IF (.not.got_var(idtime)) THEN
2044 IF (master) WRITE (stdout,60) trim(vname(1,idtime)), &
2045 & trim(label), trim(ncname)
2046 exit_flag=3
2047 RETURN
2048 END IF
2049 IF (.not.got_var(idfsur)) THEN
2050 IF (master) WRITE (stdout,60) trim(vname(1,idfsur)), &
2051 & trim(label), trim(ncname)
2052 exit_flag=3
2053 RETURN
2054 END IF
2055 IF (.not.got_var(idubar)) THEN
2056 IF (master) WRITE (stdout,60) trim(vname(1,idubar)), &
2057 & trim(label), trim(ncname)
2058 exit_flag=3
2059 RETURN
2060 END IF
2061 IF (.not.got_var(idvbar)) THEN
2062 IF (master) WRITE (stdout,60) trim(vname(1,idvbar)), &
2063 & trim(label), trim(ncname)
2064 exit_flag=3
2065 RETURN
2066 END IF
2067# ifdef ADJUST_BOUNDARY
2068 IF (.not.got_var(idsbry(isfsur))) THEN
2069 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isfsur))), &
2070 & trim(label), trim(ncname)
2071 exit_flag=3
2072 RETURN
2073 END IF
2074 IF (.not.got_var(idsbry(isubar))) THEN
2075 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isubar))), &
2076 & trim(label), trim(ncname)
2077 exit_flag=3
2078 RETURN
2079 END IF
2080 IF (.not.got_var(idsbry(isvbar))) THEN
2081 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isvbar))), &
2082 & trim(label), trim(ncname)
2083 exit_flag=3
2084 RETURN
2085 END IF
2086# endif
2087# ifdef ADJUST_WSTRESS
2088 IF (.not.got_var(idusms)) THEN
2089 IF (master) WRITE (stdout,60) trim(vname(1,idusms)), &
2090 & trim(label), trim(ncname)
2091 exit_flag=3
2092 RETURN
2093 END IF
2094 IF (.not.got_var(idvsms)) THEN
2095 IF (master) WRITE (stdout,60) trim(vname(1,idvsms)), &
2096 & trim(label), trim(ncname)
2097 exit_flag=3
2098 RETURN
2099 END IF
2100# endif
2101# ifdef SOLVE3D
2102 IF (.not.got_var(iduvel)) THEN
2103 IF (master) WRITE (stdout,60) trim(vname(1,iduvel)), &
2104 & trim(label), trim(ncname)
2105 exit_flag=3
2106 RETURN
2107 END IF
2108 IF (.not.got_var(idvvel)) THEN
2109 IF (master) WRITE (stdout,60) trim(vname(1,idvvel)), &
2110 & trim(label), trim(ncname)
2111 exit_flag=3
2112 RETURN
2113 END IF
2114# ifdef ADJUST_BOUNDARY
2115 IF (.not.got_var(idsbry(isuvel))) THEN
2116 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isuvel))), &
2117 & trim(label), trim(ncname)
2118 exit_flag=3
2119 RETURN
2120 END IF
2121 IF (.not.got_var(idsbry(isvvel))) THEN
2122 IF (master) WRITE (stdout,60) trim(vname(1,idsbry(isvvel))), &
2123 & trim(label), trim(ncname)
2124 RETURN
2125 END IF
2126# endif
2127# endif
2128# ifdef SOLVE3D
2129 DO itrc=1,nt(ng)
2130 IF (.not.got_var(idtvar(itrc))) THEN
2131 IF (master) WRITE (stdout,60) trim(vname(1,idtvar(itrc))), &
2132 & trim(label), trim(ncname)
2133 exit_flag=3
2134 RETURN
2135 END IF
2136# ifdef ADJUST_BOUNDARY
2137 IF (.not.got_var(idsbry(istvar(itrc)))) THEN
2138 IF (master) WRITE (stdout,60) &
2139 & trim(vname(1,idsbry(istvar(itrc)))), &
2140 & trim(label), trim(ncname)
2141 exit_flag=3
2142 RETURN
2143 END IF
2144# endif
2145# ifdef ADJUST_STFLUX
2146 IF (.not.got_var(idtsur(itrc)).and.lstflux(itrc,ng)) THEN
2147 IF (master) WRITE (stdout,60) trim(vname(1,idtsur(itrc))), &
2148 & trim(label), trim(ncname)
2149 exit_flag=3
2150 RETURN
2151 END IF
2152# endif
2153 END DO
2154# endif
2155!
2156! Set unlimited time record dimension to the appropriate value.
2157!
2158 s(ng)%Rindex=rec_size
2159 END IF query
2160!
2161 10 FORMAT (2x,'DEF_STATE_PIO - creating ',a,' file,',t56, &
2162 & 'Grid ',i2.2,': ',a)
2163 20 FORMAT (2x,'DEF_STATE_PIO - inquiring ',a,' file,',t56, &
2164 & 'Grid ',i2.2,': ',a)
2165 30 FORMAT (/,' DEF_STATE_PIO - unable to create ',a, &
2166 & ' NetCDF file:',1x,a)
2167 40 FORMAT (1pe11.4,1x,'millimeter')
2168 50 FORMAT (/,' DEF_STATE_PIO - unable to open ',a, &
2169 & ' NetCDF file: ',a)
2170 60 FORMAT (/,' DEF_STATE_PIO - unable to find variable: ',a,2x, &
2171 & ' in ',a,' NetCDF file: ',a)
2172!
2173 RETURN
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_close(ng, model, piofile, ncname, lupdate)
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)

References mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::idfsur, mod_ncparam::idsbry, mod_sediment::idsed, mod_ncparam::idtime, mod_ncparam::idtsur, mod_ncparam::idtvar, mod_ncparam::idubar, mod_ncparam::idusms, mod_ncparam::iduvel, mod_ncparam::idvbar, mod_ncparam::idvsms, mod_ncparam::idvvel, mod_ncparam::iinfo, mod_param::inlm, mod_param::iobounds, mod_scalars::isalt, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::itemp, mod_scalars::lobc, mod_scalars::lstflux, mod_parallel::master, mod_param::n, mod_biology::nbac, mod_biology::nbands, mod_param::nbed, mod_scalars::nbrec, mod_biology::ndom, mod_biology::nfec, mod_scalars::nfrec, mod_scalars::noerror, mod_biology::nphy, mod_param::nst, mod_fourdvar::nstatevar, mod_param::nt, mod_ncparam::nv, mod_pio_netcdf::pio_fout, mod_pio_netcdf::pio_netcdf_check_dim(), mod_pio_netcdf::pio_netcdf_close(), mod_pio_netcdf::pio_netcdf_create(), mod_pio_netcdf::pio_netcdf_enddef(), mod_pio_netcdf::pio_netcdf_inq_var(), mod_pio_netcdf::pio_netcdf_open(), mod_pio_netcdf::pio_tout, mod_param::r2dobc, mod_param::r2dvar, mod_param::r3dobc, mod_param::r3dvar, mod_scalars::rclock, mod_sediment::sd50, mod_iounits::sourcefile, mod_iounits::stdout, mod_param::u2dobc, mod_param::u2dvar, mod_param::u3dobc, mod_param::u3dvar, mod_param::v2dobc, mod_param::v2dvar, mod_param::v3dobc, mod_param::v3dvar, mod_pio_netcdf::var_desc, and mod_ncparam::vname.

Referenced by def_state().

Here is the call graph for this function:
Here is the caller graph for this function: