ROMS
Loading...
Searching...
No Matches
def_std.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#ifdef STD_MODEL
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! This module creates standard deviation output file, which is to !
14! model the 4D-Var Bacground Error Covariance, B. It defines its !
15! dimensions, attributes, and variables using either the standard !
16! NetCDF library or the Parallel-IO (PIO) library. !
17! !
18!=======================================================================
19!
20 USE mod_param
21 USE mod_parallel
22# ifdef BIOLOGY
23 USE mod_biology
24# endif
25 USE mod_fourdvar
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_std
42 PRIVATE :: def_std_nf90
43# if defined PIO_LIB && defined DISTRIBUTE
44 PRIVATE :: def_std_pio
45# endif
46!
47 CONTAINS
48!
49!***********************************************************************
50 SUBROUTINE def_std (ng)
51!***********************************************************************
52!
53! Imported variable declarations.
54!
55 integer, intent(in) :: ng
56!
57! Local variable declarations.
58!
59 character (len=*), parameter :: MyFile = &
60 & __FILE__
61!
62!-----------------------------------------------------------------------
63! Create a new history file according to IO type.
64!-----------------------------------------------------------------------
65!
66 SELECT CASE (std(5,ng)%IOtype)
67 CASE (io_nf90)
68 CALL def_std_nf90 (ng)
69
70# if defined PIO_LIB && defined DISTRIBUTE
71 CASE (io_pio)
72 CALL def_std_pio (ng)
73# endif
74 CASE DEFAULT
75 IF (master) WRITE (stdout,10) std(5,ng)%IOtype
76 exit_flag=3
77 END SELECT
78 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
79!
80 10 FORMAT (' DEF_STD - Illegal output file type, io_type = ',i0, &
81 & /,11x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
82!
83 RETURN
84 END SUBROUTINE def_std
85!
86!***********************************************************************
87 SUBROUTINE def_std_nf90 (ng)
88!***********************************************************************
89!
90 USE mod_netcdf
91!
92! Imported variable declarations.
93!
94 integer, intent(in) :: ng
95!
96! Local variable declarations.
97!
98 logical :: got_var(NV)
99!
100 integer, parameter :: Natt = 25
101
102 integer :: i, j, ifield, itrc, nrec, nvd, nvd3, nvd4
103 integer :: recdim, status, varid
104 integer :: NinnerDim
105# ifdef ADJUST_BOUNDARY
106 integer :: IorJdim, brecdim
107# endif
108# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
109 integer :: frecdim
110# endif
111 integer :: DimIDs(nDimID)
112 integer :: t2dgrd(3), u2dgrd(3), v2dgrd(3)
113# ifdef ADJUST_BOUNDARY
114 integer :: t2dobc(4)
115# endif
116
117# ifdef SOLVE3D
118 integer :: t3dgrd(4), u3dgrd(4), v3dgrd(4), w3dgrd(4)
119# ifdef ADJUST_BOUNDARY
120 integer :: t3dobc(5)
121# endif
122# ifdef ADJUST_STFLUX
123 integer :: t3dfrc(4)
124# endif
125# endif
126# ifdef ADJUST_WSTRESS
127 integer :: u3dfrc(4), v3dfrc(4)
128# endif
129!
130 real(r8) :: Aval(6)
131!
132 character (len=256) :: ncname
133 character (len=MaxLen) :: Vinfo(Natt)
134
135 character (len=*), parameter :: MyFile = &
136 & __FILE__//", def_std_nf90"
137!
138 sourcefile=myfile
139!
140!-----------------------------------------------------------------------
141! Set and report file name.
142!-----------------------------------------------------------------------
143!
144 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
145 ncname=std(5,ng)%name
146!
147 IF (master) THEN
148 IF (ldefstd(ng)) THEN
149 WRITE (stdout,10) ng, trim(ncname)
150 ELSE
151 WRITE (stdout,20) ng, trim(ncname)
152 END IF
153 END IF
154!
155!=======================================================================
156! Create a new Background standard deviation file.
157!=======================================================================
158!
159 define : IF (ldefstd(ng)) THEN
160 CALL netcdf_create (ng, inlm, trim(ncname), std(5,ng)%ncid)
161 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
162 IF (master) WRITE (stdout,30) trim(ncname)
163 RETURN
164 END IF
165!
166!-----------------------------------------------------------------------
167! Define file dimensions.
168!-----------------------------------------------------------------------
169!
170 dimids=0
171!
172 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xi_rho', &
173 & iobounds(ng)%xi_rho, dimids( 1))
174 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
175
176 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xi_u', &
177 & iobounds(ng)%xi_u, dimids( 2))
178 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
179
180 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xi_v', &
181 & iobounds(ng)%xi_v, dimids( 3))
182 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
183
184 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xi_psi', &
185 & iobounds(ng)%xi_psi, dimids( 4))
186 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
187
188 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'eta_rho', &
189 & iobounds(ng)%eta_rho, dimids( 5))
190 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
191
192 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'eta_u', &
193 & iobounds(ng)%eta_u, dimids( 6))
194 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
195
196 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'eta_v', &
197 & iobounds(ng)%eta_v, dimids( 7))
198 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
199
200 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'eta_psi', &
201 & iobounds(ng)%eta_psi, dimids( 8))
202 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
203
204# ifdef ADJUST_BOUNDARY
205 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'IorJ', &
206 & iobounds(ng)%IorJ, iorjdim)
207 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
208# endif
209
210# if defined WRITE_WATER && defined MASKING
211 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xy_rho', &
212 & iobounds(ng)%xy_rho, dimids(17))
213 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
214
215 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xy_u', &
216 & iobounds(ng)%xy_u, dimids(18))
217 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
218
219 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xy_v', &
220 & iobounds(ng)%xy_v, dimids(19))
221 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
222# endif
223
224# ifdef SOLVE3D
225# if defined WRITE_WATER && defined MASKING
226 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xyz_rho', &
227 & iobounds(ng)%xy_rho*n(ng), dimids(20))
228 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
229
230 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xyz_u', &
231 & iobounds(ng)%xy_u*n(ng), dimids(21))
232 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
233
234 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xyz_v', &
235 & iobounds(ng)%xy_v*n(ng), dimids(22))
236 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
237
238 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xyz_w', &
239 & iobounds(ng)%xy_rho*(n(ng)+1), dimids(23))
240 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
241# endif
242
243 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'N', &
244 & n(ng), dimids( 9))
245 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
246
247 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 's_rho', &
248 & n(ng), dimids( 9))
249 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
250
251 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 's_w', &
252 & n(ng)+1, dimids(10))
253 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
254
255 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'tracer', &
256 & nt(ng), dimids(11))
257 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
258
259# ifdef SEDIMENT
260 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'NST', &
261 & nst, dimids(32))
262 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
263
264 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Nbed', &
265 & nbed, dimids(16))
266 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
267
268# if defined WRITE_WATER && defined MASKING
269 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'xybed', &
270 & iobounds(ng)%xy_rho*nbed, dimids(24))
271 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
272# endif
273# endif
274
275# ifdef ECOSIM
276 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Nbands', &
277 & nbands, dimids(33))
278 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
279
280 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Nphy', &
281 & nphy, dimids(25))
282 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
283
284 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Nbac', &
285 & nbac, dimids(26))
286 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
287
288 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Ndom', &
289 & ndom, dimids(27))
290 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
291
292 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Nfec', &
293 & nfec, dimids(28))
294 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
295# endif
296# endif
297
298 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'boundary', &
299 & 4, dimids(14))
300 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
301
302# ifdef FOUR_DVAR
303 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Nstate', &
304 & nstatevar(ng), dimids(29))
305 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
306# endif
307
308 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'Ninner', &
309 & ninner, ninnerdim)
310 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
311
312# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
313 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'frc_adjust', &
314 & nfrec(ng), dimids(30))
315 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
316# endif
317
318# ifdef ADJUST_BOUNDARY
319 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, 'obc_adjust', &
320 & nbrec(ng), dimids(31))
321 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
322# endif
323
324 status=def_dim(ng, inlm, std(5,ng)%ncid, ncname, &
325 & trim(adjustl(vname(5,idtime))), &
326 & nf90_unlimited, dimids(12))
327 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
328
329 recdim=dimids(12)
330# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
331 frecdim=dimids(30)
332# endif
333# ifdef ADJUST_BOUNDARY
334 brecdim=dimids(31)
335# endif
336!
337! Set number of dimensions for output variables.
338!
339# if defined WRITE_WATER && defined MASKING
340 nvd3=2
341 nvd4=2
342# else
343 nvd3=3
344 nvd4=4
345# endif
346!
347! Define dimension vectors for staggered tracer type variables.
348!
349# if defined WRITE_WATER && defined MASKING
350 t2dgrd(1)=dimids(17)
351 t2dgrd(2)=dimids(12)
352# ifdef SOLVE3D
353 t3dgrd(1)=dimids(20)
354 t3dgrd(2)=dimids(12)
355# endif
356# else
357 t2dgrd(1)=dimids( 1)
358 t2dgrd(2)=dimids( 5)
359 t2dgrd(3)=dimids(12)
360# ifdef SOLVE3D
361 t3dgrd(1)=dimids( 1)
362 t3dgrd(2)=dimids( 5)
363 t3dgrd(3)=dimids( 9)
364 t3dgrd(4)=dimids(12)
365# endif
366# ifdef ADJUST_STFLUX
367 t3dfrc(1)=dimids( 1)
368 t3dfrc(2)=dimids( 5)
369 t3dfrc(3)=frecdim
370 t3dfrc(4)=dimids(12)
371# endif
372# endif
373# ifdef ADJUST_BOUNDARY
374 t2dobc(1)=iorjdim
375 t2dobc(2)=dimids(14)
376 t2dobc(3)=brecdim
377 t2dobc(4)=dimids(12)
378# ifdef SOLVE3D
379 t3dobc(1)=iorjdim
380 t3dobc(2)=dimids( 9)
381 t3dobc(3)=dimids(14)
382 t3dobc(4)=brecdim
383 t3dobc(5)=dimids(12)
384# endif
385# endif
386!
387! Define dimension vectors for staggered u-momentum type variables.
388!
389# if defined WRITE_WATER && defined MASKING
390 u2dgrd(1)=dimids(18)
391 u2dgrd(2)=dimids(12)
392# ifdef SOLVE3D
393 u3dgrd(1)=dimids(21)
394 u3dgrd(2)=dimids(12)
395# endif
396# else
397 u2dgrd(1)=dimids( 2)
398 u2dgrd(2)=dimids( 6)
399 u2dgrd(3)=dimids(12)
400# ifdef SOLVE3D
401 u3dgrd(1)=dimids( 2)
402 u3dgrd(2)=dimids( 6)
403 u3dgrd(3)=dimids( 9)
404 u3dgrd(4)=dimids(12)
405# endif
406# ifdef ADJUST_WSTRESS
407 u3dfrc(1)=dimids( 2)
408 u3dfrc(2)=dimids( 6)
409 u3dfrc(3)=frecdim
410 u3dfrc(4)=dimids(12)
411# endif
412# endif
413!
414! Define dimension vectors for staggered v-momentum type variables.
415!
416# if defined WRITE_WATER && defined MASKING
417 v2dgrd(1)=dimids(19)
418 v2dgrd(2)=dimids(12)
419# ifdef SOLVE3D
420 v3dgrd(1)=dimids(22)
421 v3dgrd(2)=dimids(12)
422# endif
423# else
424 v2dgrd(1)=dimids( 3)
425 v2dgrd(2)=dimids( 7)
426 v2dgrd(3)=dimids(12)
427# ifdef SOLVE3D
428 v3dgrd(1)=dimids( 3)
429 v3dgrd(2)=dimids( 7)
430 v3dgrd(3)=dimids( 9)
431 v3dgrd(4)=dimids(12)
432# endif
433# ifdef ADJUST_WSTRESS
434 v3dfrc(1)=dimids( 3)
435 v3dfrc(2)=dimids( 7)
436 v3dfrc(3)=frecdim
437 v3dfrc(4)=dimids(12)
438# endif
439# endif
440# ifdef SOLVE3D
441!
442! Define dimension vector for staggered w-momentum type variables.
443!
444# if defined WRITE_WATER && defined MASKING
445 w3dgrd(1)=dimids(23)
446 w3dgrd(2)=dimids(12)
447# else
448 w3dgrd(1)=dimids( 1)
449 w3dgrd(2)=dimids( 5)
450 w3dgrd(3)=dimids(10)
451 w3dgrd(4)=dimids(12)
452# endif
453# endif
454!
455! Initialize unlimited time record dimension.
456!
457 std(5,ng)%Rindex=0
458!
459! Initialize local information variable arrays.
460!
461 DO i=1,natt
462 DO j=1,len(vinfo(1))
463 vinfo(i)(j:j)=' '
464 END DO
465 END DO
466 DO i=1,6
467 aval(i)=0.0_r8
468 END DO
469!
470!-----------------------------------------------------------------------
471! Define time-recordless information variables.
472!-----------------------------------------------------------------------
473!
474 CALL def_info (ng, inlm, std(5,ng)%ncid, ncname, dimids)
475 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
476!
477!-----------------------------------------------------------------------
478! Define time-varying variables.
479!-----------------------------------------------------------------------
480!
481! Define model time.
482!
483 vinfo( 1)=vname(1,idtime)
484 vinfo( 2)=vname(2,idtime)
485 WRITE (vinfo( 3),'(a,a)') 'seconds since ', trim(rclock%string)
486 vinfo( 4)=trim(rclock%calendar)
487 vinfo(14)=vname(4,idtime)
488 vinfo(21)=vname(6,idtime)
489 status=def_var(ng, inlm, std(5,ng)%ncid, std(5,ng)%Vid(idtime), &
490 & nf_tout, 1, (/recdim/), aval, vinfo, ncname, &
491 & setparaccess = .true.)
492 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
493!
494! Define free-surface standard deviation.
495!
496 vinfo( 1)=vname(1,idfsur)
497 WRITE (vinfo( 2),40) trim(vname(2,idfsur))
498 vinfo( 3)='meter2'
499 vinfo(14)=vname(4,idfsur)
500 vinfo(16)=vname(1,idtime)
501# if defined WRITE_WATER && defined MASKING
502 vinfo(20)='mask_rho'
503# endif
504 vinfo(21)=vname(6,idfsur)
505 vinfo(22)='coordinates'
506 aval(5)=real(iinfo(1,idfsur,ng),r8)
507 status=def_var(ng, inlm, std(5,ng)%ncid, std(5,ng)%Vid(idfsur), &
508 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
509 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
510!
511! Define 2D U-momentum component standard deviation.
512!
513 vinfo( 1)=vname(1,idubar)
514 WRITE (vinfo( 2),40) trim(vname(2,idubar))
515 vinfo( 3)='meter2 second-2'
516 vinfo(14)=vname(4,idubar)
517 vinfo(16)=vname(1,idtime)
518# if defined WRITE_WATER && defined MASKING
519 vinfo(20)='mask_u'
520# endif
521 vinfo(21)=vname(6,idubar)
522 vinfo(22)='coordinates'
523 aval(5)=real(iinfo(1,idubar,ng),r8)
524 status=def_var(ng, inlm, std(5,ng)%ncid, std(5,ng)%Vid(idubar), &
525 & nf_fout, nvd3, u2dgrd, aval, vinfo, ncname)
526 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
527!
528! Define 2D V-momentum component standard deviation.
529!
530 vinfo( 1)=vname(1,idvbar)
531 WRITE (vinfo( 2),40) trim(vname(2,idvbar))
532 vinfo( 3)='meter2 second-2'
533 vinfo(14)=vname(4,idvbar)
534 vinfo(16)=vname(1,idtime)
535# if defined WRITE_WATER && defined MASKING
536 vinfo(20)='mask_v'
537# endif
538 vinfo(21)=vname(6,idvbar)
539 vinfo(22)='coordinates'
540 aval(5)=real(iinfo(1,idvbar,ng),r8)
541 status=def_var(ng, inlm, std(5,ng)%ncid, std(5,ng)%Vid(idvbar), &
542 & nf_fout, nvd3, v2dgrd, aval, vinfo, ncname)
543 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
544
545# ifdef SOLVE3D
546!
547! Define 3D U-momentum component standard deviation.
548!
549 vinfo( 1)=vname(1,iduvel)
550 WRITE (vinfo( 2),40) trim(vname(2,iduvel))
551 vinfo( 3)='meter2 second-2'
552 vinfo(14)=vname(4,iduvel)
553 vinfo(16)=vname(1,idtime)
554# if defined WRITE_WATER && defined MASKING
555 vinfo(20)='mask_u'
556# endif
557 vinfo(21)=vname(6,iduvel)
558 vinfo(22)='coordinates'
559 aval(5)=real(iinfo(1,iduvel,ng),r8)
560 status=def_var(ng, inlm, std(5,ng)%ncid, std(5,ng)%Vid(iduvel), &
561 & nf_fout, nvd4, u3dgrd, aval, vinfo, ncname)
562 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
563!
564! Define 3D V-momentum component standard deviation.
565!
566 vinfo( 1)=vname(1,idvvel)
567 WRITE (vinfo( 2),40) trim(vname(2,idvvel))
568 vinfo( 3)='meter2 second-2'
569 vinfo(14)=vname(4,idvvel)
570 vinfo(16)=vname(1,idtime)
571# if defined WRITE_WATER && defined MASKING
572 vinfo(20)='mask_v'
573# endif
574 vinfo(21)=vname(6,idvvel)
575 vinfo(22)='coordinates'
576 aval(5)=real(iinfo(1,idvvel,ng),r8)
577 status=def_var(ng, inlm, std(5,ng)%ncid, std(5,ng)%Vid(idvvel), &
578 & nf_fout, nvd4, v3dgrd, aval, vinfo, ncname)
579 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
580!
581! Define tracer type variables standard deviation.
582!
583 DO itrc=1,nt(ng)
584 vinfo( 1)=vname(1,idtvar(itrc))
585 WRITE (vinfo( 2),40) trim(vname(2,idtvar(itrc)))
586 IF (itrc.eq.itemp) THEN
587 vinfo( 3)='Celsius2'
588 ELSE IF (itrc.eq.isalt) THEN
589 vinfo( 3)='nondimensional'
590 ELSE
591 WRITE (vinfo( 3),50) trim(vname(3,idtvar(itrc)))
592 END IF
593 vinfo(14)=vname(4,idtvar(itrc))
594 vinfo(16)=vname(1,idtime)
595# ifdef SEDIMENT
596 DO i=1,nst
597 IF (itrc.eq.idsed(i)) THEN
598 WRITE (vinfo(19),60) 1000.0_r8*sd50(i,ng)
599 END IF
600 END DO
601# endif
602# if defined WRITE_WATER && defined MASKING
603 vinfo(20)='mask_rho'
604# endif
605 vinfo(21)=vname(6,idtvar(itrc))
606 vinfo(22)='coordinates'
607 aval(5)=real(r3dvar,r8)
608 status=def_var(ng, inlm, std(5,ng)%ncid, std(5,ng)%Tid(itrc), &
609 & nf_fout, nvd4, t3dgrd, aval, vinfo, ncname)
610 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
611 END DO
612# endif
613!
614!-----------------------------------------------------------------------
615! Leave definition mode.
616!-----------------------------------------------------------------------
617!
618 CALL netcdf_enddef (ng, inlm, ncname, std(5,ng)%ncid)
619 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
620!
621!-----------------------------------------------------------------------
622! Write out time-recordless, information variables.
623!-----------------------------------------------------------------------
624!
625 CALL wrt_info (ng, inlm, std(5,ng)%ncid, ncname)
626 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
627 END IF define
628!
629!=======================================================================
630! Open an existing posterior error covariance matrix, check its
631! contents, and prepare for appending data.
632!=======================================================================
633!
634 query: IF (.not.ldefstd(ng)) THEN
635 ncname=std(5,ng)%name
636!
637! Open posterior error covariance file for read/write.
638!
639 CALL netcdf_open (ng, inlm, ncname, 1, std(5,ng)%ncid)
640 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
641 WRITE (stdout,70) trim(ncname)
642 RETURN
643 END IF
644!
645! Inquire about the dimensions and check for consistency.
646!
647 CALL netcdf_check_dim (ng, inlm, ncname, &
648 & ncid = std(5,ng)%ncid)
649 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
650!
651! Inquire about the variables.
652!
653 CALL netcdf_inq_var (ng, inlm, ncname, &
654 & ncid = std(5,ng)%ncid)
655 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
656!
657! Initialize logical switches.
658!
659 DO i=1,nv
660 got_var(i)=.false.
661 END DO
662!
663! Scan variable list from input NetCDF and activate switches for
664! posterior error covariance matrix variables. Get variable IDs.
665!
666 DO i=1,n_var
667 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
668 got_var(idtime)=.true.
669 std(5,ng)%Vid(idtime)=var_id(i)
670 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
671 got_var(idfsur)=.true.
672 std(5,ng)%Vid(idfsur)=var_id(i)
673# ifndef SOLVE3D
674 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
675 got_var(idubar)=.true.
676 std(5,ng)%Vid(idubar)=var_id(i)
677 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
678 got_var(idvbar)=.true.
679 std(5,ng)%Vid(idvbar)=var_id(i)
680# endif
681# ifdef SOLVE3D
682 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
683 got_var(iduvel)=.true.
684 std(5,ng)%Vid(iduvel)=var_id(i)
685 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
686 got_var(idvvel)=.true.
687 std(5,ng)%Vid(idvvel)=var_id(i)
688# endif
689 END IF
690# ifdef SOLVE3D
691 DO itrc=1,nt(ng)
692 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
693 got_var(idtvar(itrc))=.true.
694 std(5,ng)%Tid(itrc)=var_id(i)
695 END IF
696 END DO
697# endif
698 END DO
699!
700! Check if posterior error covariance matrix variables are available
701! in input NetCDF file.
702!
703 IF (.not.got_var(idtime)) THEN
704 IF (master) WRITE (stdout,80) trim(vname(1,idtime)), &
705 & trim(ncname)
706 exit_flag=3
707 RETURN
708 END IF
709 IF (.not.got_var(idfsur)) THEN
710 IF (master) WRITE (stdout,80) trim(vname(1,idfsur)), &
711 & trim(ncname)
712 exit_flag=3
713 RETURN
714 END IF
715# ifndef SOLVE3D
716 IF (.not.got_var(idubar)) THEN
717 IF (master) WRITE (stdout,80) trim(vname(1,idubar)), &
718 & trim(ncname)
719 exit_flag=3
720 RETURN
721 END IF
722 IF (.not.got_var(idvbar)) THEN
723 IF (master) WRITE (stdout,80) trim(vname(1,idvbar)), &
724 & trim(ncname)
725 exit_flag=3
726 RETURN
727 END IF
728# endif
729# ifdef SOLVE3D
730 IF (.not.got_var(iduvel)) THEN
731 IF (master) WRITE (stdout,80) trim(vname(1,iduvel)), &
732 & trim(ncname)
733 exit_flag=3
734 RETURN
735 END IF
736 IF (.not.got_var(idvvel)) THEN
737 IF (master) WRITE (stdout,80) trim(vname(1,idvvel)), &
738 & trim(ncname)
739 exit_flag=3
740 RETURN
741 END IF
742# endif
743# ifdef SOLVE3D
744 DO itrc=1,nt(ng)
745 IF (.not.got_var(idtvar(itrc))) THEN
746 IF (master) WRITE (stdout,80) trim(vname(1,idtvar(itrc))), &
747 & trim(ncname)
748 exit_flag=3
749 RETURN
750 END IF
751 END DO
752# endif
753!
754! Set unlimited time record dimension to the appropriate value.
755!
756 std(5,ng)%Rindex=rec_size
757 END IF query
758!
759 10 FORMAT (2x,'DEF_STD_NF90 - creating standard deviation', &
760 & ' file,',t56,'Grid ',i2.2,': ',a)
761 20 FORMAT (2x,'DEF_STD_NF90 - inquiring standard deviation', &
762 & ' file,',t56,'Grid ',i2.2,': ',a)
763 30 FORMAT (/,' DEF_STD_NF90 - unable to create standard deviation' &
764 & ' NetCDF file:',1x,a)
765 40 FORMAT (a,', standard deviation')
766 50 FORMAT ('(',a,')^2')
767 60 FORMAT (1pe11.4,1x,'millimeter')
768 70 FORMAT (/,' DEF_STD_NF90 - unable to open standard deviation', &
769 & ' NetCDF file: ',a)
770 80 FORMAT (/,' DEF_STD_NF90 - unable to find variable: ',a,2x, &
771 & ' in standard deviation NetCDF file: ',a)
772!
773 RETURN
774 END SUBROUTINE def_std_nf90
775
776# if defined PIO_LIB && defined DISTRIBUTE
777!
778!***********************************************************************
779 SUBROUTINE def_std_pio (ng)
780!***********************************************************************
781!
783!
784! Imported variable declarations.
785!
786 integer, intent(in) :: ng
787!
788! Local variable declarations.
789!
790 logical :: got_var(NV)
791!
792 integer, parameter :: Natt = 25
793
794 integer :: i, j, ifield, itrc, nrec, nvd, nvd3, nvd4
795 integer :: recdim, status
796 integer :: NinnerDim
797# ifdef ADJUST_BOUNDARY
798 integer :: IorJdim, brecdim
799# endif
800# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
801 integer :: frecdim
802# endif
803 integer :: DimIDs(nDimID)
804 integer :: t2dgrd(3), u2dgrd(3), v2dgrd(3)
805# ifdef ADJUST_BOUNDARY
806 integer :: t2dobc(4)
807# endif
808# ifdef SOLVE3D
809 integer :: t3dgrd(4), u3dgrd(4), v3dgrd(4), w3dgrd(4)
810# ifdef ADJUST_BOUNDARY
811 integer :: t3dobc(5)
812# endif
813# ifdef ADJUST_STFLUX
814 integer :: t3dfrc(4)
815# endif
816# endif
817# ifdef ADJUST_WSTRESS
818 integer :: u3dfrc(4), v3dfrc(4)
819# endif
820!
821 real(r8) :: Aval(6)
822!
823 character (len=256) :: ncname
824 character (len=MaxLen) :: Vinfo(Natt)
825
826 character (len=*), parameter :: MyFile = &
827 & __FILE__//", def_std_pio"
828!
829 TYPE (Var_desc_t) :: varDesc
830!
831 sourcefile=myfile
832!
833!-----------------------------------------------------------------------
834! Set and report file name.
835!-----------------------------------------------------------------------
836!
837 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
838 ncname=std(5,ng)%name
839!
840 IF (master) THEN
841 IF (ldefstd(ng)) THEN
842 WRITE (stdout,10) ng, trim(ncname)
843 ELSE
844 WRITE (stdout,20) ng, trim(ncname)
845 END IF
846 END IF
847!
848!=======================================================================
849! Create a new posterior Background standard deviation file.
850!=======================================================================
851!
852 define : IF (ldefstd(ng)) THEN
853 CALL pio_netcdf_create (ng, inlm, trim(ncname), &
854 & std(5,ng)%pioFile)
855 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
856 IF (master) WRITE (stdout,30) trim(ncname)
857 RETURN
858 END IF
859!
860!-----------------------------------------------------------------------
861! Define file dimensions.
862!-----------------------------------------------------------------------
863!
864 dimids=0
865!
866 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xi_rho', &
867 & iobounds(ng)%xi_rho, dimids( 1))
868 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
869
870 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xi_u', &
871 & iobounds(ng)%xi_u, dimids( 2))
872 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
873
874 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xi_v', &
875 & iobounds(ng)%xi_v, dimids( 3))
876 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
877
878 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xi_psi', &
879 & iobounds(ng)%xi_psi, dimids( 4))
880 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
881
882 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'eta_rho', &
883 & iobounds(ng)%eta_rho, dimids( 5))
884 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
885
886 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'eta_u', &
887 & iobounds(ng)%eta_u, dimids( 6))
888 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
889
890 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'eta_v', &
891 & iobounds(ng)%eta_v, dimids( 7))
892 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
893
894 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'eta_psi', &
895 & iobounds(ng)%eta_psi, dimids( 8))
896 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
897
898# ifdef ADJUST_BOUNDARY
899 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'IorJ', &
900 & iobounds(ng)%IorJ, iorjdim)
901 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
902# endif
903
904# if defined WRITE_WATER && defined MASKING
905 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xy_rho', &
906 & iobounds(ng)%xy_rho, dimids(17))
907 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
908
909 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xy_u', &
910 & iobounds(ng)%xy_u, dimids(18))
911 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
912
913 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xy_v', &
914 & iobounds(ng)%xy_v, dimids(19))
915 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
916# endif
917
918# ifdef SOLVE3D
919# if defined WRITE_WATER && defined MASKING
920 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xyz_rho', &
921 & iobounds(ng)%xy_rho*n(ng), dimids(20))
922 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
923
924 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xyz_u', &
925 & iobounds(ng)%xy_u*n(ng), dimids(21))
926 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
927
928 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xyz_v', &
929 & iobounds(ng)%xy_v*n(ng), dimids(22))
930 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
931
932 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xyz_w', &
933 & iobounds(ng)%xy_rho*(n(ng)+1), dimids(23))
934 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
935# endif
936
937 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'N', &
938 & n(ng), dimids( 9))
939 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
940
941 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 's_rho', &
942 & n(ng), dimids( 9))
943 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
944
945 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 's_w', &
946 & n(ng)+1, dimids(10))
947 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
948
949 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'tracer', &
950 & nt(ng), dimids(11))
951 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
952
953# ifdef SEDIMENT
954 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'NST', &
955 & nst, dimids(32))
956 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
957
958 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Nbed', &
959 & nbed, dimids(16))
960 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
961
962# if defined WRITE_WATER && defined MASKING
963 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'xybed', &
964 & iobounds(ng)%xy_rho*nbed, dimids(24))
965 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
966# endif
967# endif
968
969# ifdef ECOSIM
970 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Nbands', &
971 & nbands, dimids(33))
972 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
973
974 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Nphy', &
975 & nphy, dimids(25))
976 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
977
978 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Nbac', &
979 & nbac, dimids(26))
980 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
981
982 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Ndom', &
983 & ndom, dimids(27))
984 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
985
986 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Nfec', &
987 & nfec, dimids(28))
988 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
989# endif
990# endif
991
992 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'boundary', &
993 & 4, dimids(14))
994 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
995
996# ifdef FOUR_DVAR
997 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Nstate', &
998 & nstatevar(ng), dimids(29))
999 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1000# endif
1001
1002 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, 'Ninner', &
1003 & ninner, ninnerdim)
1004 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1005
1006# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1007 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, &
1008 & 'frc_adjust', nfrec(ng), dimids(30))
1009 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1010# endif
1011
1012# ifdef ADJUST_BOUNDARY
1013 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, &
1014 & 'obc_adjust', nbrec(ng), dimids(31))
1015 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1016# endif
1017
1018 status=def_dim(ng, inlm, std(5,ng)%pioFile, ncname, &
1019 & trim(adjustl(vname(5,idtime))), &
1020 & pio_unlimited, dimids(12))
1021 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1022
1023 recdim=dimids(12)
1024# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1025 frecdim=dimids(30)
1026# endif
1027# ifdef ADJUST_BOUNDARY
1028 brecdim=dimids(31)
1029# endif
1030!
1031! Set number of dimensions for output variables.
1032!
1033# if defined WRITE_WATER && defined MASKING
1034 nvd3=2
1035 nvd4=2
1036# else
1037 nvd3=3
1038 nvd4=4
1039# endif
1040!
1041! Define dimension vectors for staggered tracer type variables.
1042!
1043# if defined WRITE_WATER && defined MASKING
1044 t2dgrd(1)=dimids(17)
1045 t2dgrd(2)=dimids(12)
1046# ifdef SOLVE3D
1047 t3dgrd(1)=dimids(20)
1048 t3dgrd(2)=dimids(12)
1049# endif
1050# else
1051 t2dgrd(1)=dimids( 1)
1052 t2dgrd(2)=dimids( 5)
1053 t2dgrd(3)=dimids(12)
1054# ifdef SOLVE3D
1055 t3dgrd(1)=dimids( 1)
1056 t3dgrd(2)=dimids( 5)
1057 t3dgrd(3)=dimids( 9)
1058 t3dgrd(4)=dimids(12)
1059# endif
1060# ifdef ADJUST_STFLUX
1061 t3dfrc(1)=dimids( 1)
1062 t3dfrc(2)=dimids( 5)
1063 t3dfrc(3)=frecdim
1064 t3dfrc(4)=dimids(12)
1065# endif
1066# endif
1067# ifdef ADJUST_BOUNDARY
1068 t2dobc(1)=iorjdim
1069 t2dobc(2)=dimids(14)
1070 t2dobc(3)=brecdim
1071 t2dobc(4)=dimids(12)
1072# ifdef SOLVE3D
1073 t3dobc(1)=iorjdim
1074 t3dobc(2)=dimids( 9)
1075 t3dobc(3)=dimids(14)
1076 t3dobc(4)=brecdim
1077 t3dobc(5)=dimids(12)
1078# endif
1079# endif
1080!
1081! Define dimension vectors for staggered u-momentum type variables.
1082!
1083# if defined WRITE_WATER && defined MASKING
1084 u2dgrd(1)=dimids(18)
1085 u2dgrd(2)=dimids(12)
1086# ifdef SOLVE3D
1087 u3dgrd(1)=dimids(21)
1088 u3dgrd(2)=dimids(12)
1089# endif
1090# else
1091 u2dgrd(1)=dimids( 2)
1092 u2dgrd(2)=dimids( 6)
1093 u2dgrd(3)=dimids(12)
1094# ifdef SOLVE3D
1095 u3dgrd(1)=dimids( 2)
1096 u3dgrd(2)=dimids( 6)
1097 u3dgrd(3)=dimids( 9)
1098 u3dgrd(4)=dimids(12)
1099# endif
1100# ifdef ADJUST_WSTRESS
1101 u3dfrc(1)=dimids( 2)
1102 u3dfrc(2)=dimids( 6)
1103 u3dfrc(3)=frecdim
1104 u3dfrc(4)=dimids(12)
1105# endif
1106# endif
1107!
1108! Define dimension vectors for staggered v-momentum type variables.
1109!
1110# if defined WRITE_WATER && defined MASKING
1111 v2dgrd(1)=dimids(19)
1112 v2dgrd(2)=dimids(12)
1113# ifdef SOLVE3D
1114 v3dgrd(1)=dimids(22)
1115 v3dgrd(2)=dimids(12)
1116# endif
1117# else
1118 v2dgrd(1)=dimids( 3)
1119 v2dgrd(2)=dimids( 7)
1120 v2dgrd(3)=dimids(12)
1121# ifdef SOLVE3D
1122 v3dgrd(1)=dimids( 3)
1123 v3dgrd(2)=dimids( 7)
1124 v3dgrd(3)=dimids( 9)
1125 v3dgrd(4)=dimids(12)
1126# endif
1127# ifdef ADJUST_WSTRESS
1128 v3dfrc(1)=dimids( 3)
1129 v3dfrc(2)=dimids( 7)
1130 v3dfrc(3)=frecdim
1131 v3dfrc(4)=dimids(12)
1132# endif
1133# endif
1134# ifdef SOLVE3D
1135!
1136! Define dimension vector for staggered w-momentum type variables.
1137!
1138# if defined WRITE_WATER && defined MASKING
1139 w3dgrd(1)=dimids(23)
1140 w3dgrd(2)=dimids(12)
1141# else
1142 w3dgrd(1)=dimids( 1)
1143 w3dgrd(2)=dimids( 5)
1144 w3dgrd(3)=dimids(10)
1145 w3dgrd(4)=dimids(12)
1146# endif
1147# endif
1148!
1149! Initialize unlimited time record dimension.
1150!
1151 std(5,ng)%Rindex=0
1152!
1153! Initialize local information variable arrays.
1154!
1155 DO i=1,natt
1156 DO j=1,len(vinfo(1))
1157 vinfo(i)(j:j)=' '
1158 END DO
1159 END DO
1160 DO i=1,6
1161 aval(i)=0.0_r8
1162 END DO
1163!
1164!-----------------------------------------------------------------------
1165! Define time-recordless information variables.
1166!-----------------------------------------------------------------------
1167!
1168 CALL def_info (ng, inlm, std(5,ng)%pioFile, ncname, dimids)
1169 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1170!
1171!-----------------------------------------------------------------------
1172! Define time-varying variables.
1173!-----------------------------------------------------------------------
1174!
1175! Define model time.
1176!
1177 vinfo( 1)=vname(1,idtime)
1178 vinfo( 2)=vname(2,idtime)
1179 WRITE (vinfo( 3),'(a,a)') 'seconds since ', trim(rclock%string)
1180 vinfo( 4)=trim(rclock%calendar)
1181 vinfo(14)=vname(4,idtime)
1182 vinfo(21)=vname(6,idtime)
1183 std(5,ng)%pioVar(idtime)%dkind=pio_tout
1184 std(5,ng)%pioVar(idtime)%gtype=0
1185!
1186 status=def_var(ng, inlm, std(5,ng)%pioFile, &
1187 & std(5,ng)%pioVar(idtime)%vd, &
1188 & pio_tout, 1, (/recdim/), aval, vinfo, ncname, &
1189 & setparaccess = .true.)
1190 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1191!
1192! Define free-surface standard deviation.
1193!
1194 vinfo( 1)=vname(1,idfsur)
1195 WRITE (vinfo( 2),40) trim(vname(2,idfsur))
1196 vinfo( 3)='meter2'
1197 vinfo(14)=vname(4,idfsur)
1198 vinfo(16)=vname(1,idtime)
1199# if defined WRITE_WATER && defined MASKING
1200 vinfo(20)='mask_rho'
1201# endif
1202 vinfo(21)=vname(6,idfsur)
1203 vinfo(22)='coordinates'
1204 aval(5)=real(iinfo(1,idfsur,ng),r8)
1205 std(5,ng)%pioVar(idfsur)%dkind=pio_fout
1206 std(5,ng)%pioVar(idfsur)%gtype=r2dvar
1207!
1208 status=def_var(ng, inlm, std(5,ng)%pioFile, &
1209 & std(5,ng)%pioVar(idfsur)%vd, &
1210 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1211 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1212!
1213! Define 2D U-momentum component standard deviation.
1214!
1215 vinfo( 1)=vname(1,idubar)
1216 WRITE (vinfo( 2),40) trim(vname(2,idubar))
1217 vinfo( 3)='meter2 second-2'
1218 vinfo(14)=vname(4,idubar)
1219 vinfo(16)=vname(1,idtime)
1220# if defined WRITE_WATER && defined MASKING
1221 vinfo(20)='mask_u'
1222# endif
1223 vinfo(21)=vname(6,idubar)
1224 vinfo(22)='coordinates'
1225 aval(5)=real(iinfo(1,idubar,ng),r8)
1226 std(5,ng)%pioVar(idubar)%dkind=pio_fout
1227 std(5,ng)%pioVar(idubar)%gtype=u2dvar
1228!
1229 status=def_var(ng, inlm, std(5,ng)%pioFile, &
1230 & std(5,ng)%pioVar(idubar)%vd, &
1231 & pio_fout, nvd3, u2dgrd, aval, vinfo, ncname)
1232 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1233!
1234! Define 2D V-momentum component standard deviation.
1235!
1236 vinfo( 1)=vname(1,idvbar)
1237 WRITE (vinfo( 2),40) trim(vname(2,idvbar))
1238 vinfo( 3)='meter2 second-2'
1239 vinfo(14)=vname(4,idvbar)
1240 vinfo(16)=vname(1,idtime)
1241# if defined WRITE_WATER && defined MASKING
1242 vinfo(20)='mask_v'
1243# endif
1244 vinfo(21)=vname(6,idvbar)
1245 vinfo(22)='coordinates'
1246 aval(5)=real(iinfo(1,idvbar,ng),r8)
1247 std(5,ng)%pioVar(idvbar)%dkind=pio_fout
1248 std(5,ng)%pioVar(idvbar)%gtype=v2dvar
1249!
1250 status=def_var(ng, inlm, std(5,ng)%pioFile, &
1251 & std(5,ng)%pioVar(idvbar)%vd, &
1252 & pio_fout, nvd3, v2dgrd, aval, vinfo, ncname)
1253 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1254
1255# ifdef SOLVE3D
1256!
1257! Define 3D U-momentum component standard deviation.
1258!
1259 vinfo( 1)=vname(1,iduvel)
1260 WRITE (vinfo( 2),40) trim(vname(2,iduvel))
1261 vinfo( 3)='meter2 second-2'
1262 vinfo(14)=vname(4,iduvel)
1263 vinfo(16)=vname(1,idtime)
1264# if defined WRITE_WATER && defined MASKING
1265 vinfo(20)='mask_u'
1266# endif
1267 vinfo(21)=vname(6,iduvel)
1268 vinfo(22)='coordinates'
1269 aval(5)=real(iinfo(1,iduvel,ng),r8)
1270 std(5,ng)%pioVar(iduvel)%dkind=pio_fout
1271 std(5,ng)%pioVar(iduvel)%gtype=u3dvar
1272!
1273 status=def_var(ng, inlm, std(5,ng)%pioFile, &
1274 & std(5,ng)%pioVar(iduvel)%vd, &
1275 & pio_fout, nvd4, u3dgrd, aval, vinfo, ncname)
1276 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1277!
1278! Define 3D V-momentum component standard deviation.
1279!
1280 vinfo( 1)=vname(1,idvvel)
1281 WRITE (vinfo( 2),40) trim(vname(2,idvvel))
1282 vinfo( 3)='meter2 second-2'
1283 vinfo(14)=vname(4,idvvel)
1284 vinfo(16)=vname(1,idtime)
1285# if defined WRITE_WATER && defined MASKING
1286 vinfo(20)='mask_v'
1287# endif
1288 vinfo(21)=vname(6,idvvel)
1289 vinfo(22)='coordinates'
1290 aval(5)=real(iinfo(1,idvvel,ng),r8)
1291 std(5,ng)%pioVar(idvvel)%dkind=pio_fout
1292 std(5,ng)%pioVar(idvvel)%gtype=v3dvar
1293!
1294 status=def_var(ng, inlm, std(5,ng)%pioFile, &
1295 & std(5,ng)%pioVar(idvvel)%vd, &
1296 & pio_fout, nvd4, v3dgrd, aval, vinfo, ncname)
1297 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1298!
1299! Define tracer type variables standard deviation.
1300!
1301 DO itrc=1,nt(ng)
1302 vinfo( 1)=vname(1,idtvar(itrc))
1303 WRITE (vinfo( 2),40) trim(vname(2,idtvar(itrc)))
1304 IF (itrc.eq.itemp) THEN
1305 vinfo( 3)='Celsius2'
1306 ELSE IF (itrc.eq.isalt) THEN
1307 vinfo( 3)='nondimensional'
1308 ELSE
1309 WRITE (vinfo( 3),50) trim(vname(3,idtvar(itrc)))
1310 END IF
1311 vinfo(14)=vname(4,idtvar(itrc))
1312 vinfo(16)=vname(1,idtime)
1313# ifdef SEDIMENT
1314 DO i=1,nst
1315 IF (itrc.eq.idsed(i)) THEN
1316 WRITE (vinfo(19),60) 1000.0_r8*sd50(i,ng)
1317 END IF
1318 END DO
1319# endif
1320# if defined WRITE_WATER && defined MASKING
1321 vinfo(20)='mask_rho'
1322# endif
1323 vinfo(21)=vname(6,idtvar(itrc))
1324 vinfo(22)='coordinates'
1325 aval(5)=real(r3dvar,r8)
1326 std(5,ng)%pioTrc(itrc)%dkind=pio_fout
1327 std(5,ng)%pioTrc(itrc)%gtype=r3dvar
1328!
1329 status=def_var(ng, inlm, std(5,ng)%pioFile, &
1330 & std(5,ng)%pioTrc(itrc)%vd, &
1331 & pio_fout, nvd4, t3dgrd, aval, vinfo, ncname)
1332 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1333 END DO
1334# endif
1335!
1336!-----------------------------------------------------------------------
1337! Leave definition mode.
1338!-----------------------------------------------------------------------
1339!
1340 CALL pio_netcdf_enddef (ng, inlm, ncname, std(5,ng)%pioFile)
1341 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1342!
1343!-----------------------------------------------------------------------
1344! Write out time-recordless, information variables.
1345!-----------------------------------------------------------------------
1346!
1347 CALL wrt_info (ng, inlm, std(5,ng)%pioFile, ncname)
1348 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1349 END IF define
1350!
1351!=======================================================================
1352! Open an existing posterior error covariance matrix, check its
1353! contents, and prepare for appending data.
1354!=======================================================================
1355!
1356 query: IF (.not.ldefstd(ng)) THEN
1357 ncname=std(5,ng)%name
1358!
1359! Open posterior error covariance file for read/write.
1360!
1361 CALL pio_netcdf_open (ng, inlm, ncname, 1, std(5,ng)%pioFile)
1362 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1363 WRITE (stdout,70) trim(ncname)
1364 RETURN
1365 END IF
1366!
1367! Inquire about the dimensions and check for consistency.
1368!
1369 CALL pio_netcdf_check_dim (ng, inlm, ncname, &
1370 & piofile = std(5,ng)%pioFile)
1371 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1372!
1373! Inquire about the variables.
1374!
1375 CALL pio_netcdf_inq_var (ng, inlm, ncname, &
1376 & piofile = std(5,ng)%pioFile)
1377 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1378!
1379! Initialize logical switches.
1380!
1381 DO i=1,nv
1382 got_var(i)=.false.
1383 END DO
1384!
1385! Scan variable list from input NetCDF and activate switches for
1386! posterior error covariance matrix variables. Get variable IDs.
1387!
1388 DO i=1,n_var
1389 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
1390 got_var(idtime)=.true.
1391 std(5,ng)%pioVar(idtime)%vd=var_desc(i)
1392 std(5,ng)%pioVar(idtime)%dkind=pio_tout
1393 std(5,ng)%pioVar(idtime)%gtype=0
1394 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
1395 got_var(idfsur)=.true.
1396 std(5,ng)%pioVar(idfsur)%vd=var_desc(i)
1397 std(5,ng)%pioVar(idfsur)%dkind=pio_fout
1398 std(5,ng)%pioVar(idfsur)%gtype=r2dvar
1399# ifndef SOLVE3D
1400 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
1401 got_var(idubar)=.true.
1402 std(5,ng)%pioVar(idubar)%vd=var_desc(i)
1403 std(5,ng)%pioVar(idubar)%dkind=pio_fout
1404 std(5,ng)%pioVar(idubar)%gtype=u2dvar
1405 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
1406 got_var(idvbar)=.true.
1407 std(5,ng)%pioVar(idvbar)%vd=var_desc(i)
1408 std(5,ng)%pioVar(idvbar)%dkind=pio_fout
1409 std(5,ng)%pioVar(idvbar)%gtype=v2dvar
1410# endif
1411# ifdef SOLVE3D
1412 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
1413 got_var(iduvel)=.true.
1414 std(5,ng)%pioVar(iduvel)%vd=var_desc(i)
1415 std(5,ng)%pioVar(iduvel)%dkind=pio_fout
1416 std(5,ng)%pioVar(iduvel)%gtype=u3dvar
1417 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
1418 got_var(idvvel)=.true.
1419 std(5,ng)%pioVar(idvvel)%vd=var_desc(i)
1420 std(5,ng)%pioVar(idvvel)%dkind=pio_fout
1421 std(5,ng)%pioVar(idvvel)%gtype=v3dvar
1422# endif
1423 END IF
1424# ifdef SOLVE3D
1425 DO itrc=1,nt(ng)
1426 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
1427 got_var(idtvar(itrc))=.true.
1428 std(5,ng)%pioTrc(itrc)%vd=var_desc(i)
1429 std(5,ng)%pioTrc(itrc)%dkind=pio_fout
1430 std(5,ng)%pioTrc(itrc)%gtype=r3dvar
1431 END IF
1432 END DO
1433# endif
1434 END DO
1435!
1436! Check if posterior error covariance matrix variables are available
1437! in input NetCDF file.
1438!
1439 IF (.not.got_var(idtime)) THEN
1440 IF (master) WRITE (stdout,80) trim(vname(1,idtime)), &
1441 & trim(ncname)
1442 exit_flag=3
1443 RETURN
1444 END IF
1445 IF (.not.got_var(idfsur)) THEN
1446 IF (master) WRITE (stdout,80) trim(vname(1,idfsur)), &
1447 & trim(ncname)
1448 exit_flag=3
1449 RETURN
1450 END IF
1451 IF (.not.got_var(idubar)) THEN
1452 IF (master) WRITE (stdout,80) trim(vname(1,idubar)), &
1453 & trim(ncname)
1454 exit_flag=3
1455 RETURN
1456 END IF
1457 IF (.not.got_var(idvbar)) THEN
1458 IF (master) WRITE (stdout,80) trim(vname(1,idvbar)), &
1459 & trim(ncname)
1460 exit_flag=3
1461 RETURN
1462 END IF
1463# ifdef SOLVE3D
1464 IF (.not.got_var(iduvel)) THEN
1465 IF (master) WRITE (stdout,80) trim(vname(1,iduvel)), &
1466 & trim(ncname)
1467 exit_flag=3
1468 RETURN
1469 END IF
1470 IF (.not.got_var(idvvel)) THEN
1471 IF (master) WRITE (stdout,80) trim(vname(1,idvvel)), &
1472 & trim(ncname)
1473 exit_flag=3
1474 RETURN
1475 END IF
1476# endif
1477# ifdef SOLVE3D
1478 DO itrc=1,nt(ng)
1479 IF (.not.got_var(idtvar(itrc))) THEN
1480 IF (master) WRITE (stdout,80) trim(vname(1,idtvar(itrc))), &
1481 & trim(ncname)
1482 exit_flag=3
1483 RETURN
1484 END IF
1485 END DO
1486# endif
1487!
1488! Set unlimited time record dimension to the appropriate value.
1489!
1490 std(5,ng)%Rindex=rec_size
1491 END IF query
1492!
1493 10 FORMAT (2x,'DEF_STD_PIO - creating standard deviation', &
1494 & ' file,',t56,'Grid ',i2.2,': ',a)
1495 20 FORMAT (2x,'DEF_STD_PIO - inquiring standard deviation', &
1496 & ' file,',t56,'Grid ',i2.2,': ',a)
1497 30 FORMAT (/,' DEF_STD_PIO - unable to create standard deviation', &
1498 & ' NetCDF file:',1x,a)
1499 40 FORMAT (a,', standard deviation')
1500 50 FORMAT ('(',a,')^2')
1501 60 FORMAT (1pe11.4,1x,'millimeter')
1502 70 FORMAT (/,' DEF_STD_PIO - unable to open standard deviation', &
1503 & ' NetCDF file: ',a)
1504 80 FORMAT (/,' DEF_STD_PIO - unable to find variable: ',a,2x, &
1505 & ' in standard deviation NetCDF file: ',a)
1506!
1507 RETURN
1508 END SUBROUTINE def_std_pio
1509# endif
1510#endif
1511 END MODULE def_std_mod
integer, parameter nfec
Definition ecosim_mod.h:204
integer, parameter nbac
Definition ecosim_mod.h:202
integer, parameter ndom
Definition ecosim_mod.h:203
integer, parameter nbands
Definition ecosim_mod.h:201
integer, parameter nphy
Definition ecosim_mod.h:205
integer, dimension(:), allocatable nstatevar
type(t_io), dimension(:,:), allocatable std
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer idubar
integer idvvel
integer, parameter io_pio
Definition mod_ncparam.F:96
integer idfsur
integer, dimension(:), allocatable idtvar
integer iduvel
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer, dimension(:,:,:), allocatable iinfo
integer idvbar
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, 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
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, dimension(:), allocatable nt
Definition mod_param.F:489
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
integer nst
Definition mod_param.F:521
integer, parameter v3dvar
Definition mod_param.F:723
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)
integer ninner
integer, dimension(:), allocatable nfrec
type(t_clock) rclock
integer exit_flag
integer isalt
integer itemp
integer, dimension(:), allocatable nbrec
integer noerror
integer, dimension(:), allocatable idsed
real(r8), dimension(:,:), allocatable sd50
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52