ROMS
Loading...
Searching...
No Matches
wrt_std.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#if defined 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 writes 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 USE mod_forces
23 USE mod_fourdvar
24 USE mod_grid
25 USE mod_iounits
26 USE mod_mixing
27 USE mod_ncparam
28 USE mod_ocean
29 USE mod_scalars
30 USE mod_stepping
31!
33# ifdef SOLVE3D
35# endif
36 USE strings_mod, ONLY : founderror
37!
38 implicit none
39!
40 PUBLIC :: wrt_std
41 PRIVATE :: wrt_std_nf90
42# if defined PIO_LIB && defined DISTRIBUTE
43 PRIVATE :: wrt_std_pio
44# endif
45!
46 CONTAINS
47!
48!***********************************************************************
49 SUBROUTINE wrt_std (ng, tile, Lstd)
50!***********************************************************************
51!
52! Imported variable declarations.
53!
54 integer, intent(in) :: ng, tile, Lstd
55!
56! Local variable declarations.
57!
58 integer :: LBi, UBi, LBj, UBj
59!
60 character (len=*), parameter :: MyFile = &
61 & __FILE__
62!
63!-----------------------------------------------------------------------
64! Write out Background standard deviation fields according to IO type.
65!-----------------------------------------------------------------------
66!
67 lbi=bounds(ng)%LBi(tile)
68 ubi=bounds(ng)%UBi(tile)
69 lbj=bounds(ng)%LBj(tile)
70 ubj=bounds(ng)%UBj(tile)
71!
72 SELECT CASE (std(5,ng)%IOtype)
73 CASE (io_nf90)
74 CALL wrt_std_nf90 (ng, tile, lstd, &
75 & lbi, ubi, lbj, ubj)
76
77# if defined PIO_LIB && defined DISTRIBUTE
78 CASE (io_pio)
79 CALL wrt_std_pio (ng, tile, lstd, &
80 & lbi, ubi, lbj, ubj)
81# endif
82 CASE DEFAULT
83 IF (master) WRITE (stdout,10) std(5,ng)%IOtype
84 exit_flag=3
85 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
86 END SELECT
87 lwrtstd(ng)=.false.
88!
89 10 FORMAT (' WRT_STD - Illegal output file type, io_type = ',i0, &
90 & /,13x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
91!
92 RETURN
93 END SUBROUTINE wrt_std
94!
95!***********************************************************************
96 SUBROUTINE wrt_std_nf90 (ng, tile, Lstd, &
97 & LBi, UBi, LBj, UBj)
98!***********************************************************************
99!
100 USE mod_netcdf
101!
102! Imported variable declarations.
103!
104 integer, intent(in) :: ng, tile, Lstd
105 integer, intent(in) :: LBi, UBi, LBj, UBj
106!
107! Local variable declarations.
108!
109 integer :: Fcount, i, j, gfactor, gtype, status
110# ifdef SOLVE3D
111 integer :: itrc, k
112# endif
113!
114 real(dp) :: scale
115!
116 character (len=*), parameter :: MyFile = &
117 & __FILE__//", wrt_std_nf90"
118!
119 sourcefile=myfile
120!
121!-----------------------------------------------------------------------
122! Write out full posterior error covariance (diagonal) matrix.
123!-----------------------------------------------------------------------
124!
125 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
126!
127! Set grid type factor to write full (gfactor=1) fields or water
128! points (gfactor=-1) fields only.
129!
130# if defined WRITE_WATER && defined MASKING
131 gfactor=-1
132# else
133 gfactor=1
134# endif
135!
136! Set time record index.
137!
138 std(5,ng)%Rindex=std(5,ng)%Rindex+1
139 fcount=std(5,ng)%Fcount
140 std(5,ng)%Nrec(fcount)=std(5,ng)%Nrec(fcount)+1
141!
142! Report.
143!
144# ifdef SOLVE3D
145# ifdef NESTING
146 IF (master) WRITE (stdout,10) lstd, lstd, std(5,ng)%Rindex, ng
147# else
148 IF (master) WRITE (stdout,10) lstd, lstd, std(5,ng)%Rindex
149# endif
150# else
151# ifdef NESTING
152 IF (master) WRITE (stdout,10) lstd, std(5,ng)%Rindex, ng
153# else
154 IF (master) WRITE (stdout,10) lstd, std(5,ng)%Rindex
155# endif
156# endif
157!
158! Write out model time (s).
159!
160 CALL netcdf_put_fvar (ng, inlm, std(5,ng)%name, &
161 & trim(vname(1,idtime)), time(ng:), &
162 & (/std(5,ng)%Rindex/), (/1/), &
163 & ncid = std(5,ng)%ncid, &
164 & varid = std(5,ng)%Vid(idtime))
165 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
166!
167! Write out free-surface standard deviation.
168!
169 scale=1.0_dp
170 gtype=gfactor*r2dvar
171 status=nf_fwrite2d(ng, inlm, std(5,ng)%ncid, idfsur, &
172 & std(5,ng)%Vid(idfsur), &
173 & std(5,ng)%Rindex, gtype, &
174 & lbi, ubi, lbj, ubj, scale, &
175# ifdef MASKING
176 & grid(ng) % rmask, &
177# endif
178 & ocean(ng)% e_zeta(:,:,lstd))
179 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
180 IF (master) THEN
181 WRITE (stdout,20) trim(vname(1,idfsur)), std(5,ng)%Rindex
182 END IF
183 exit_flag=3
184 ioerror=status
185 RETURN
186 END IF
187!
188! Write out 2D U-momentum component standard deviation.
189!
190 scale=1.0_dp
191 gtype=gfactor*u2dvar
192 status=nf_fwrite2d(ng, inlm, std(5,ng)%ncid, idubar, &
193 & std(5,ng)%Vid(idubar), &
194 & std(5,ng)%Rindex, gtype, &
195 & lbi, ubi, lbj, ubj, scale, &
196# ifdef MASKING
197 & grid(ng) % umask_full, &
198# endif
199 & ocean(ng) % e_ubar(:,:,lstd))
200 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
201 IF (master) THEN
202 WRITE (stdout,20) trim(vname(1,idubar)), std(5,ng)%Rindex
203 END IF
204 exit_flag=3
205 ioerror=status
206 RETURN
207 END IF
208!
209! Write out 2D V-momentum component standard deviation.
210!
211 scale=1.0_dp
212 gtype=gfactor*v2dvar
213 status=nf_fwrite2d(ng, inlm, std(5,ng)%ncid, idvbar, &
214 & std(5,ng)%Vid(idvbar), &
215 & std(5,ng)%Rindex, gtype, &
216 & lbi, ubi, lbj, ubj, scale, &
217# ifdef MASKING
218 & grid(ng) % vmask_full, &
219# endif
220 & ocean(ng) % e_vbar(:,:,lstd))
221 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
222 IF (master) THEN
223 WRITE (stdout,20) trim(vname(1,idvbar)), std(5,ng)%Rindex
224 END IF
225 exit_flag=3
226 ioerror=status
227 RETURN
228 END IF
229
230# ifdef SOLVE3D
231!
232! Write out 3D U-momentum component standard deviation.
233!
234 scale=1.0_dp
235 gtype=gfactor*u3dvar
236 status=nf_fwrite3d(ng, inlm, std(5,ng)%ncid, iduvel, &
237 & std(5,ng)%Vid(iduvel), &
238 & std(5,ng)%Rindex, gtype, &
239 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
240# ifdef MASKING
241 & grid(ng) % umask_full, &
242# endif
243 & ocean(ng) % e_u(:,:,:,lstd))
244 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
245 IF (master) THEN
246 WRITE (stdout,20) trim(vname(1,iduvel)), std(5,ng)%Rindex
247 END IF
248 exit_flag=3
249 ioerror=status
250 RETURN
251 END IF
252!
253! Write out 3D V-momentum component standard deviation.
254!
255 scale=1.0_dp
256 gtype=gfactor*v3dvar
257 status=nf_fwrite3d(ng, inlm, std(5,ng)%ncid, idvvel, &
258 & std(5,ng)%Vid(idvvel), &
259 & std(5,ng)%Rindex, gtype, &
260 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
261# ifdef MASKING
262 & grid(ng) % vmask_full, &
263# endif
264 & ocean(ng) % e_v(:,:,:,lstd))
265 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
266 IF (master) THEN
267 WRITE (stdout,20) trim(vname(1,idvvel)), std(5,ng)%Rindex
268 END IF
269 exit_flag=3
270 ioerror=status
271 RETURN
272 END IF
273!
274! Write out tracer type variables standard deviation.
275!
276 DO itrc=1,nt(ng)
277 scale=1.0_dp
278 gtype=gfactor*r3dvar
279 status=nf_fwrite3d(ng, inlm, std(5,ng)%ncid, idtvar(itrc), &
280 & std(5,ng)%Tid(itrc), &
281 & std(5,ng)%Rindex, gtype, &
282 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
283# ifdef MASKING
284 & grid(ng) % rmask, &
285# endif
286 & ocean(ng) % e_t(:,:,:,lstd,itrc))
287 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
288 IF (master) THEN
289 WRITE (stdout,20) trim(vname(1,idtvar(itrc))), &
290 & std(5,ng)%Rindex
291 END IF
292 exit_flag=3
293 ioerror=status
294 RETURN
295 END IF
296 END DO
297# endif
298!
299!-----------------------------------------------------------------------
300! Synchronize Background standard deviation NetCDF file to disk to
301! allow other processes to access data immediately after it is written.
302!-----------------------------------------------------------------------
303!
304 CALL netcdf_sync (ng, inlm, std(5,ng)%name, std(5,ng)%ncid)
305 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
306!
307 10 FORMAT (2x,'WRT_STD_NF90 - writing standard deviation', t42, &
308# ifdef SOLVE3D
309# ifdef NESTING
310 & 'fields (Index=',i1,',',i1,') in record = ',i0,t92,i2.2)
311# else
312 & 'fields (Index=',i1,',',i1,') in record = ',i0)
313# endif
314# else
315# ifdef NESTING
316 & 'fields (Index=',i1,') in record = ',i0,t92,i2.2)
317# else
318 & 'fields (Index=',i1,') in record = ',i0)
319# endif
320# endif
321 20 FORMAT (/,' WRT_STD_NF90 - error while writing variable: ',a, &
322 & /,16x,'into 4DVar error NetCDF file for time record: ',i0)
323!
324 RETURN
325 END SUBROUTINE wrt_std_nf90
326
327# if defined PIO_LIB && defined DISTRIBUTE
328!
329!***********************************************************************
330 SUBROUTINE wrt_std_pio (ng, tile, Lstd, &
331 & LBi, UBi, LBj, UBj)
332!***********************************************************************
333!
335!
336! Imported variable declarations.
337!
338 integer, intent(in) :: ng, tile, Lstd
339 integer, intent(in) :: LBi, UBi, LBj, UBj
340!
341! Local variable declarations.
342!
343 integer :: Fcount, i, ifield, j, status
344# ifdef SOLVE3D
345 integer :: itrc, k
346# endif
347!
348 real(dp) :: scale
349!
350 character (len=*), parameter :: MyFile = &
351 & __FILE__//", wrt_std_nf90"
352!
353 TYPE (IO_desc_t), pointer :: ioDesc
354!
355 sourcefile=myfile
356!
357!-----------------------------------------------------------------------
358! Write out full posterior error covariance (diagonal) matrix.
359!-----------------------------------------------------------------------
360!
361 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
362!
363! Set time record index.
364!
365 std(5,ng)%Rindex=std(5,ng)%Rindex+1
366 fcount=std(5,ng)%Fcount
367 std(5,ng)%Nrec(fcount)=std(5,ng)%Nrec(fcount)+1
368!
369! Report.
370!
371# ifdef SOLVE3D
372# ifdef NESTING
373 IF (master) WRITE (stdout,20) lstd, lstd, std(5,ng)%Rindex, ng
374# else
375 IF (master) WRITE (stdout,20) lstd, lstd, std(5,ng)%Rindex
376# endif
377# else
378# ifdef NESTING
379 IF (master) WRITE (stdout,20) lstd, std(5,ng)%Rindex, ng
380# else
381 IF (master) WRITE (stdout,20) lstd, std(5,ng)%Rindex
382# endif
383# endif
384!
385! Write out model time (s).
386!
387 CALL pio_netcdf_put_fvar (ng, inlm, std(5,ng)%name, &
388 & trim(vname(1,idtime)), time(ng:), &
389 & (/std(5,ng)%Rindex/), (/1/), &
390 & piofile = std(5,ng)%pioFile, &
391 & piovar = std(5,ng)%pioVar(idtime)%vd)
392 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
393!
394! Write out free-surface standard deviation.
395!
396 scale=1.0_dp
397 IF (std(5,ng)%pioVar(idfsur)%dkind.eq.pio_double) THEN
398 iodesc => iodesc_dp_r2dvar(ng)
399 ELSE
400 iodesc => iodesc_sp_r2dvar(ng)
401 END IF
402!
403 status=nf_fwrite2d(ng, inlm, std(5,ng)%pioFile, idfsur, &
404 & std(5,ng)%pioVar(idfsur), &
405 & std(5,ng)%Rindex, &
406 & iodesc, &
407 & lbi, ubi, lbj, ubj, scale, &
408# ifdef MASKING
409 & grid(ng) % rmask, &
410# endif
411 & ocean(ng)% e_zeta(:,:,lstd))
412 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
413 IF (master) THEN
414 WRITE (stdout,20) trim(vname(1,idfsur)), std(5,ng)%Rindex
415 END IF
416 exit_flag=3
417 ioerror=status
418 RETURN
419 END IF
420
421# ifndef SOLVE3D
422!
423! Write out 2D U-momentum component standard deviation.
424!
425 scale=1.0_dp
426 IF (std(5,ng)%pioVar(idubar)%dkind.eq.pio_double) THEN
427 iodesc => iodesc_dp_u2dvar(ng)
428 ELSE
429 iodesc => iodesc_sp_u2dvar(ng)
430 END IF
431!
432 status=nf_fwrite2d(ng, inlm, std(5,ng)%pioFile, idubar, &
433 & std(5,ng)%pioVar(idubar), &
434 & std(5,ng)%Rindex, &
435 & iodesc, &
436 & lbi, ubi, lbj, ubj, scale, &
437# ifdef MASKING
438 & grid(ng) % umask_full, &
439# endif
440 & ocean(ng) % e_ubar(:,:,lstd))
441 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
442 IF (master) THEN
443 WRITE (stdout,20) trim(vname(1,idubar)), std(5,ng)%Rindex
444 END IF
445 exit_flag=3
446 ioerror=status
447 RETURN
448 END IF
449!
450! Write out 2D V-momentum component standard deviation.
451!
452 scale=1.0_dp
453 IF (std(5,ng)%pioVar(idvbar)%dkind.eq.pio_double) THEN
454 iodesc => iodesc_dp_v2dvar(ng)
455 ELSE
456 iodesc => iodesc_sp_v2dvar(ng)
457 END IF
458!
459 status=nf_fwrite2d(ng, inlm, std(5,ng)%pioFile, idvbar, &
460 & std(5,ng)%pioVar(idvbar), &
461 & std(5,ng)%Rindex, &
462 & iodesc, &
463 & lbi, ubi, lbj, ubj, scale, &
464# ifdef MASKING
465 & grid(ng) % vmask_full, &
466# endif
467 & ocean(ng) % e_vbar(:,:,lstd))
468 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
469 IF (master) THEN
470 WRITE (stdout,20) trim(vname(1,idvbar)), std(5,ng)%Rindex
471 END IF
472 exit_flag=3
473 ioerror=status
474 RETURN
475 END IF
476# endif
477
478# ifdef SOLVE3D
479!
480! Write out 3D U-momentum component standard deviation.
481!
482 scale=1.0_dp
483 IF (std(5,ng)%pioVar(iduvel)%dkind.eq.pio_double) THEN
484 iodesc => iodesc_dp_u3dvar(ng)
485 ELSE
486 iodesc => iodesc_sp_u3dvar(ng)
487 END IF
488!
489 status=nf_fwrite3d(ng, inlm, std(5,ng)%pioFile, iduvel, &
490 & std(5,ng)%pioVar(iduvel), &
491 & std(5,ng)%Rindex, &
492 & iodesc, &
493 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
494# ifdef MASKING
495 & grid(ng) % umask_full, &
496# endif
497 & ocean(ng) % e_u(:,:,:,lstd))
498 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
499 IF (master) THEN
500 WRITE (stdout,20) trim(vname(1,iduvel)), std(5,ng)%Rindex
501 END IF
502 exit_flag=3
503 ioerror=status
504 RETURN
505 END IF
506!
507! Write out 3D V-momentum component error variance.
508!
509 scale=1.0_dp
510 IF (std(5,ng)%pioVar(idfsur)%dkind.eq.pio_double) THEN
511 iodesc => iodesc_dp_v3dvar(ng)
512 ELSE
513 iodesc => iodesc_sp_v3dvar(ng)
514 END IF
515!
516 status=nf_fwrite3d(ng, inlm, std(5,ng)%pioFile, idvvel, &
517 & std(5,ng)%pioVar(idvvel), &
518 & std(5,ng)%Rindex, &
519 & iodesc, &
520 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
521# ifdef MASKING
522 & grid(ng) % vmask_full, &
523# endif
524 & ocean(ng) % e_v(:,:,:,lstd))
525 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
526 IF (master) THEN
527 WRITE (stdout,20) trim(vname(1,idvvel)), std(5,ng)%Rindex
528 END IF
529 exit_flag=3
530 ioerror=status
531 RETURN
532 END IF
533!
534! Write out tracer type variables standard deviation.
535!
536 DO itrc=1,nt(ng)
537 scale=1.0_dp
538 IF (std(5,ng)%pioTrc(itrc)%dkind.eq.pio_double) THEN
539 iodesc => iodesc_dp_r3dvar(ng)
540 ELSE
541 iodesc => iodesc_sp_r3dvar(ng)
542 END IF
543!
544 status=nf_fwrite3d(ng, inlm, std(5,ng)%pioFile, idtvar(itrc), &
545 & std(5,ng)%pioTrc(itrc), &
546 & std(5,ng)%Rindex, &
547 & iodesc, &
548 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
549# ifdef MASKING
550 & grid(ng) % rmask, &
551# endif
552 & ocean(ng) % e_t(:,:,:,lstd,itrc))
553 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
554 IF (master) THEN
555 WRITE (stdout,20) trim(vname(1,idtvar(itrc))), &
556 & std(5,ng)%Rindex
557 END IF
558 exit_flag=3
559 ioerror=status
560 RETURN
561 END IF
562 END DO
563# endif
564!
565!-----------------------------------------------------------------------
566! Synchronize Background standard deviation NetCDF file to disk to
567! allow other processes to access data immediately after it is written.
568!-----------------------------------------------------------------------
569!
570 CALL pio_netcdf_sync (ng, inlm, std(5,ng)%name, &
571 & std(5,ng)%pioFile)
572 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
573!
574 10 FORMAT (2x,'WRT_STD_PIO - writing standard deviation', t42, &
575# ifdef SOLVE3D
576# ifdef NESTING
577 & 'fields (Index=',i1,',',i1,') in record = ',i0,t92,i2.2)
578# else
579 & 'fields (Index=',i1,',',i1,') in record = ',i0)
580# endif
581# else
582# ifdef NESTING
583 & 'fields (Index=',i1,') in record = ',i0,t92,i2.2)
584# else
585 & 'fields (Index=',i1,') in record = ',i0)
586# endif
587# endif
588 20 FORMAT (/,' WRT_STD_PIO - error while writing variable: ',a, &
589 & /,15x,'into 4DVar error NetCDF file for time record: ',i0)
590!
591 RETURN
592 END SUBROUTINE wrt_std_pio
593# endif
594#endif
595 END MODULE wrt_std_mod
596
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_io), dimension(:,:), allocatable std
integer ioerror
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 idvbar
subroutine, public netcdf_sync(ng, model, ncname, ncid)
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
logical master
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter r3dvar
Definition mod_param.F:721
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, parameter v3dvar
Definition mod_param.F:723
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dvar
subroutine, public pio_netcdf_sync(ng, model, ncname, piofile)
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dvar
integer exit_flag
real(dp), dimension(:), allocatable time
integer noerror
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52