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

Functions/Subroutines

subroutine, public wrt_error (ng, tile, kout, nout)
 
subroutine, private wrt_error_nf90 (ng, tile, kout, nout, lbij, ubij, lbi, ubi, lbj, ubj)
 
subroutine, private wrt_error_pio (ng, tile, kout, nout, lbij, ubij, lbi, ubi, lbj, ubj)
 

Function/Subroutine Documentation

◆ wrt_error()

subroutine, public wrt_error_mod::wrt_error ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) kout,
integer, intent(in) nout )

Definition at line 61 of file wrt_error.F.

62!***********************************************************************
63!
64! Imported variable declarations.
65!
66 integer, intent(in) :: ng, tile, kout, nout
67!
68! Local variable declarations.
69!
70# ifdef ADJUST_BOUNDARY
71 integer :: LBij, UBij
72# endif
73 integer :: LBi, UBi, LBj, UBj
74!
75 character (len=*), parameter :: MyFile = &
76 & __FILE__
77!
78!-----------------------------------------------------------------------
79! Write out time-averaged fields according to IO type.
80!-----------------------------------------------------------------------
81!
82# ifdef ADJUST_BOUNDARY
83 lbij=bounds(ng)%LBij
84 ubij=bounds(ng)%UBij
85# endif
86 lbi=bounds(ng)%LBi(tile)
87 ubi=bounds(ng)%UBi(tile)
88 lbj=bounds(ng)%LBj(tile)
89 ubj=bounds(ng)%UBj(tile)
90!
91 SELECT CASE (err(ng)%IOtype)
92 CASE (io_nf90)
93 CALL wrt_error_nf90 (ng, tile, kout, nout, &
94# ifdef ADJUST_BOUNDARY
95 & lbij, ubij, &
96# endif
97 & lbi, ubi, lbj, ubj)
98
99# if defined PIO_LIB && defined DISTRIBUTE
100 CASE (io_pio)
101 CALL wrt_error_pio (ng, tile, kout, nout, &
102# ifdef ADJUST_BOUNDARY
103 & lbij, ubij, &
104# endif
105 & lbi, ubi, lbj, ubj)
106# endif
107 CASE DEFAULT
108 IF (master) WRITE (stdout,10) err(ng)%IOtype
109 exit_flag=3
110 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
111 END SELECT
112!
113 10 FORMAT (' WRT_ERROR - Illegal output file type, io_type = ',i0, &
114 & /,13x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
115!
116 RETURN

References mod_param::bounds, mod_iounits::err, mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::io_nf90, mod_ncparam::io_pio, mod_parallel::master, mod_scalars::noerror, mod_iounits::stdout, wrt_error_nf90(), and wrt_error_pio().

Referenced by rbl4dvar_mod::posterior_error().

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

◆ wrt_error_nf90()

subroutine, private wrt_error_mod::wrt_error_nf90 ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) kout,
integer, intent(in) nout,
integer, intent(in) lbij,
integer, intent(in) ubij,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj )
private

Definition at line 120 of file wrt_error.F.

125!***********************************************************************
126!
127 USE mod_netcdf
128!
129! Imported variable declarations.
130!
131 integer, intent(in) :: ng, tile, kout, nout
132# ifdef ADJUST_BOUNDARY
133 integer, intent(in) :: LBij, UBij
134# endif
135 integer, intent(in) :: LBi, UBi, LBj, UBj
136!
137! Local variable declarations.
138!
139 integer :: Fcount, i, j, gfactor, gtype, status
140# ifdef SOLVE3D
141 integer :: itrc, k
142# endif
143!
144 real(dp) :: scale
145!
146 character (len=*), parameter :: MyFile = &
147 & __FILE__//", wrt_error_nf90"
148!
149 sourcefile=myfile
150!
151!-----------------------------------------------------------------------
152! Write out full posterior error covariance (diagonal) matrix.
153!-----------------------------------------------------------------------
154!
155 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
156!
157! Set grid type factor to write full (gfactor=1) fields or water
158! points (gfactor=-1) fields only.
159!
160# if defined WRITE_WATER && defined MASKING
161 gfactor=-1
162# else
163 gfactor=1
164# endif
165!
166! Set time record index.
167!
168 err(ng)%Rindex=err(ng)%Rindex+1
169 fcount=err(ng)%Fcount
170 err(ng)%Nrec(fcount)=err(ng)%Nrec(fcount)+1
171!
172! Report.
173!
174# ifdef SOLVE3D
175# ifdef NESTING
176 IF (master) WRITE (stdout,10) kout, nout, err(ng)%Rindex, ng
177# else
178 IF (master) WRITE (stdout,10) kout, nout, err(ng)%Rindex
179# endif
180# else
181# ifdef NESTING
182 IF (master) WRITE (stdout,10) kout, err(ng)%Rindex, ng
183# else
184 IF (master) WRITE (stdout,10) kout, err(ng)%Rindex
185# endif
186# endif
187!
188! Write out model time (s).
189!
190 CALL netcdf_put_fvar (ng, itlm, err(ng)%name, &
191 & trim(vname(1,idtime)), time(ng:), &
192 & (/err(ng)%Rindex/), (/1/), &
193 & ncid = err(ng)%ncid, &
194 & varid = err(ng)%Vid(idtime))
195 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
196!
197! Write out inner-loop Lanczos vectors tridiagonal system.
198!
199 CALL netcdf_put_fvar (ng, itlm, err(ng)%name, &
200 & 'zLanczos_coef', zlanczos_coef, &
201 & (/1,1/), (/ninner,ninner/), &
202 & ncid = err(ng)%ncid)
203 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
204
205 CALL netcdf_put_fvar (ng, itlm, err(ng)%name, &
206 & 'zLanczos_inv', zlanczos_inv, &
207 & (/1,1/), (/ninner,ninner/), &
208 & ncid = err(ng)%ncid)
209 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
210
211 CALL netcdf_put_fvar (ng, itlm, err(ng)%name, &
212 & 'zLanczos_err', zlanczos_err, &
213 & (/1,1/), (/ninner,ninner/), &
214 & ncid = err(ng)%ncid)
215 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
216!
217! Write out free-surface error variance.
218!
219 scale=1.0_dp
220 gtype=gfactor*r2dvar
221 status=nf_fwrite2d(ng, itlm, err(ng)%ncid, idfsur, &
222 & err(ng)%Vid(idfsur), &
223 & err(ng)%Rindex, gtype, &
224 & lbi, ubi, lbj, ubj, scale, &
225# ifdef MASKING
226 & grid(ng) % rmask, &
227# endif
228 & ocean(ng)% tl_zeta(:,:,kout))
229 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
230 IF (master) THEN
231 WRITE (stdout,20) trim(vname(1,idfsur)), err(ng)%Rindex
232 END IF
233 exit_flag=3
234 ioerror=status
235 RETURN
236 END IF
237
238# ifdef ADJUST_BOUNDARY
239!
240! Write out free-surface open boundaries error variance.
241!
242 IF (any(lobc(:,isfsur,ng))) THEN
243 scale=1.0_dp
244 status=nf_fwrite2d_bry(ng, itlm, err(ng)%name, err(ng)%ncid, &
245 & vname(1,idsbry(isfsur)), &
246 & err(ng)%Vid(idsbry(isfsur)), &
247 & err(ng)%Rindex, r2dvar, &
248 & lbij, ubij, nbrec(ng), scale, &
249 & boundary(ng) % tl_zeta_obc(lbij:,:,:, &
250 & kout))
251 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
252 IF (master) THEN
253 WRITE (stdout,20) trim(vname(1,idsbry(isfsur))), &
254 & err(ng)%Rindex
255 END IF
256 exit_flag=3
257 ioerror=status
258 RETURN
259 END IF
260 END IF
261# endif
262!
263! Write out 2D U-momentum component error variance.
264!
265 scale=1.0_dp
266 gtype=gfactor*u2dvar
267 status=nf_fwrite2d(ng, itlm, err(ng)%ncid, idubar, &
268 & err(ng)%Vid(idubar), &
269 & err(ng)%Rindex, gtype, &
270 & lbi, ubi, lbj, ubj, scale, &
271# ifdef MASKING
272 & grid(ng) % umask_full, &
273# endif
274 & ocean(ng) % tl_ubar(:,:,kout))
275 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
276 IF (master) THEN
277 WRITE (stdout,20) trim(vname(1,idubar)), err(ng)%Rindex
278 END IF
279 exit_flag=3
280 ioerror=status
281 RETURN
282 END IF
283
284# ifdef ADJUST_BOUNDARY
285!
286! Write out 2D U-momentum component open boundaries error variance.
287!
288 IF (any(lobc(:,isubar,ng))) THEN
289 scale=1.0_dp
290 status=nf_fwrite2d_bry(ng, itlm, err(ng)%name, err(ng)%ncid, &
291 & vname(1,idsbry(isubar)), &
292 & err(ng)%Vid(idsbry(isubar)), &
293 & err(ng)%Rindex, u2dvar, &
294 & lbij, ubij, nbrec(ng), scale, &
295 & boundary(ng) % tl_ubar_obc(lbij:,:,:, &
296 & kout))
297 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
298 IF (master) THEN
299 WRITE (stdout,20) trim(vname(1,idsbry(isubar))), &
300 & err(ng)%Rindex
301 END IF
302 exit_flag=3
303 ioerror=status
304 RETURN
305 END IF
306 END IF
307# endif
308!
309! Write out 2D V-momentum component error variance.
310!
311 scale=1.0_dp
312 gtype=gfactor*v2dvar
313 status=nf_fwrite2d(ng, itlm, err(ng)%ncid, idvbar, &
314 & err(ng)%Vid(idvbar), &
315 & err(ng)%Rindex, gtype, &
316 & lbi, ubi, lbj, ubj, scale, &
317# ifdef MASKING
318 & grid(ng) % vmask_full, &
319# endif
320 & ocean(ng) % tl_vbar(:,:,kout))
321 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
322 IF (master) THEN
323 WRITE (stdout,20) trim(vname(1,idvbar)), err(ng)%Rindex
324 END IF
325 exit_flag=3
326 ioerror=status
327 RETURN
328 END IF
329
330# ifdef ADJUST_BOUNDARY
331!
332! Write out 2D V-momentum component open boundaries error variance.
333!
334 IF (any(lobc(:,isvbar,ng))) THEN
335 scale=1.0_dp
336 status=nf_fwrite2d_bry(ng, itlm, err(ng)%name, err(ng)%ncid, &
337 & vname(1,idsbry(isvbar)), &
338 & err(ng)%Vid(idsbry(isvbar)), &
339 & err(ng)%Rindex, v2dvar, &
340 & lbij, ubij, nbrec(ng), scale, &
341 & boundary(ng) % tl_vbar_obc(lbij:,:,:, &
342 & kout))
343 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
344 IF (master) THEN
345 WRITE (stdout,20) trim(vname(1,idsbry(isvbar))), &
346 & err(ng)%Rindex
347 END IF
348 exit_flag=3
349 ioerror=status
350 RETURN
351 END IF
352 END IF
353# endif
354
355# ifdef SOLVE3D
356!
357! Write out 3D U-momentum component error variance.
358!
359 scale=1.0_dp
360 gtype=gfactor*u3dvar
361 status=nf_fwrite3d(ng, itlm, err(ng)%ncid, iduvel, &
362 & err(ng)%Vid(iduvel), &
363 & err(ng)%Rindex, gtype, &
364 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
365# ifdef MASKING
366 & grid(ng) % umask_full, &
367# endif
368 & ocean(ng) % tl_u(:,:,:,nout))
369 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
370 IF (master) THEN
371 WRITE (stdout,20) trim(vname(1,iduvel)), err(ng)%Rindex
372 END IF
373 exit_flag=3
374 ioerror=status
375 RETURN
376 END IF
377
378# ifdef ADJUST_BOUNDARY
379!
380! Write out 3D U-momentum component open boundaries error variance.
381!
382 IF (any(lobc(:,isuvel,ng))) THEN
383 scale=1.0_dp
384 status=nf_fwrite3d_bry(ng, itlm, err(ng)%name, err(ng)%ncid, &
385 & vname(1,idsbry(isuvel)), &
386 & err(ng)%Vid(idsbry(isuvel)), &
387 & err(ng)%Rindex, u3dvar, &
388 & lbij, ubij, 1, n(ng), nbrec(ng), scale, &
389 & boundary(ng) % tl_u_obc(lbij:,:,:,:, &
390 & nout))
391 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
392 IF (master) THEN
393 WRITE (stdout,20) trim(vname(1,idsbry(isuvel))), &
394 & err(ng)%Rindex
395 END IF
396 exit_flag=3
397 ioerror=status
398 RETURN
399 END IF
400 END IF
401# endif
402!
403! Write out 3D V-momentum component error variance.
404!
405 scale=1.0_dp
406 gtype=gfactor*v3dvar
407 status=nf_fwrite3d(ng, itlm, err(ng)%ncid, idvvel, &
408 & err(ng)%Vid(idvvel), &
409 & err(ng)%Rindex, gtype, &
410 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
411# ifdef MASKING
412 & grid(ng) % vmask_full, &
413# endif
414 & ocean(ng) % tl_v(:,:,:,nout))
415 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
416 IF (master) THEN
417 WRITE (stdout,20) trim(vname(1,idvvel)), err(ng)%Rindex
418 END IF
419 exit_flag=3
420 ioerror=status
421 RETURN
422 END IF
423
424# ifdef ADJUST_BOUNDARY
425!
426! Write out 3D V-momentum component open boundaries error variance.
427!
428 IF (any(lobc(:,isvvel,ng))) THEN
429 scale=1.0_dp
430 status=nf_fwrite3d_bry(ng, itlm, err(ng)%name, err(ng)%ncid, &
431 & vname(1,idsbry(isvvel)), &
432 & err(ng)%Vid(idsbry(isvvel)), &
433 & err(ng)%Rindex, v3dvar, &
434 & lbij, ubij, 1, n(ng), nbrec(ng), scale, &
435 & boundary(ng) % tl_v_obc(lbij:,:,:,:, &
436 & nout))
437 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
438 IF (master) THEN
439 WRITE (stdout,20) trim(vname(1,idsbry(isvvel))), &
440 & err(ng)%Rindex
441 END IF
442 exit_flag=3
443 ioerror=status
444 RETURN
445 END IF
446 END IF
447# endif
448!
449! Write out tracer type variables error variance.
450!
451 DO itrc=1,nt(ng)
452 scale=1.0_dp
453 gtype=gfactor*r3dvar
454 status=nf_fwrite3d(ng, itlm, err(ng)%ncid, idtvar(itrc), &
455 & err(ng)%Tid(itrc), &
456 & err(ng)%Rindex, gtype, &
457 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
458# ifdef MASKING
459 & grid(ng) % rmask, &
460# endif
461 & ocean(ng) % tl_t(:,:,:,nout,itrc))
462 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
463 IF (master) THEN
464 WRITE (stdout,20) trim(vname(1,idtvar(itrc))), &
465 & err(ng)%Rindex
466 END IF
467 exit_flag=3
468 ioerror=status
469 RETURN
470 END IF
471 END DO
472
473# ifdef ADJUST_BOUNDARY
474!
475! Write out tracers open boundaries error variance.
476!
477 DO itrc=1,nt(ng)
478 IF (any(lobc(:,istvar(itrc),ng))) THEN
479 scale=1.0_dp
480 status=nf_fwrite3d_bry(ng, itlm, err(ng)%name, err(ng)%ncid, &
481 & vname(1,idsbry(istvar(itrc))), &
482 & err(ng)%Vid(idsbry(istvar(itrc))), &
483 & err(ng)%Rindex, r3dvar, &
484 & lbij, ubij, 1, n(ng), nbrec(ng), &
485 & scale, &
486 & boundary(ng) % tl_t_obc(lbij:,:,:,:, &
487 & nout,itrc))
488 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
489 IF (master) THEN
490 WRITE (stdout,20) trim(vname(1,idsbry(istvar(itrc)))), &
491 & err(ng)%Rindex
492 END IF
493 exit_flag=3
494 ioerror=status
495 RETURN
496 END IF
497 END IF
498 END DO
499# endif
500
501# ifdef ADJUST_STFLUX
502!
503! Write out surface net tracers fluxes error variance. Notice that
504! fluxes have their own fixed time-dimension (of size Nfrec) to allow
505! 4DVar adjustments at other times in addition to initialization time.
506!
507 DO itrc=1,nt(ng)
508 IF (lstflux(itrc,ng)) THEN
509 scale=1.0_dp
510 gtype=gfactor*r3dvar
511 status=nf_fwrite3d(ng, itlm, err(ng)%ncid, idtsur(itrc), &
512 & err(ng)%Vid(idtsur(itrc)), &
513 & err(ng)%Rindex, gtype, &
514 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
515# ifdef MASKING
516 & grid(ng) % rmask, &
517# endif
518 & forces(ng) % tl_tflux(:,:,:,kout,itrc))
519 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
520 IF (master) THEN
521 WRITE (stdout,20) trim(vname(1,idtsur(itrc))), &
522 & err(ng)%Rindex
523 END IF
524 exit_flag=3
525 ioerror=status
526 RETURN
527 END IF
528 END IF
529 END DO
530# endif
531# endif
532# ifdef ADJUST_WSTRESS
533!
534! Write out surface U-momentum stress error variance. Notice that the
535! stress has its own fixed time-dimension (of size Nfrec) to allow
536! 4DVar adjustments at other times in addition to initialization time.
537!
538 scale=1.0_dp
539 gtype=gfactor*u3dvar
540 status=nf_fwrite3d(ng, itlm, err(ng)%ncid, idusms, &
541 & err(ng)%Vid(idusms), &
542 & err(ng)%Rindex, gtype, &
543 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
544# ifdef MASKING
545 & grid(ng) % umask, &
546# endif
547 & forces(ng) % tl_ustr(:,:,:,kout))
548 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
549 IF (master) THEN
550 WRITE (stdout,20) trim(vname(1,idusms)), err(ng)%Rindex
551 END IF
552 exit_flag=3
553 ioerror=status
554 RETURN
555 END IF
556!
557! Write out surface V-momentum stress error variance.
558!
559 scale=1.0_dp
560 gtype=gfactor*v3dvar
561 status=nf_fwrite3d(ng, itlm, err(ng)%ncid, idvsms, &
562 & err(ng)%Vid(idvsms), &
563 & err(ng)%Rindex, gtype, &
564 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
565# ifdef MASKING
566 & grid(ng) % vmask, &
567# endif
568 & forces(ng) % tl_vstr(:,:,:,kout))
569 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
570 IF (master) THEN
571 WRITE (stdout,20) trim(vname(1,idvsms)), err(ng)%Rindex
572 END IF
573 exit_flag=3
574 ioerror=status
575 RETURN
576 END IF
577# endif
578!
579!-----------------------------------------------------------------------
580! Synchronize posterior error covariance NetCDF file to disk to allow
581! other processes to access data immediately after it is written.
582!-----------------------------------------------------------------------
583!
584 CALL netcdf_sync (ng, itlm, err(ng)%name, err(ng)%ncid)
585 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
586!
587 10 FORMAT (2x,'WRT_ERROR_NF90 - writing error', t42, &
588# ifdef SOLVE3D
589# ifdef NESTING
590 & 'fields (Index=',i1,',',i1,') in record = ',i0,t92,i2.2)
591# else
592 & 'fields (Index=',i1,',',i1,') in record = ',i0)
593# endif
594# else
595# ifdef NESTING
596 & 'fields (Index=',i1,') in record = ',i0,t92,i2.2)
597# else
598 & 'fields (Index=',i1,') in record = ',i0)
599# endif
600# endif
601 20 FORMAT (/,' WRT_ERROR_NF90 - error while writing variable: ',a, &
602 & /,18x,'into 4DVar error NetCDF file for time record: ',i0)
603!
604 RETURN
subroutine, public netcdf_sync(ng, model, ncname, ncid)

References mod_boundary::boundary, mod_iounits::err, mod_scalars::exit_flag, mod_forces::forces, strings_mod::founderror(), mod_grid::grid, mod_ncparam::idfsur, mod_ncparam::idsbry, 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_iounits::ioerror, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_param::itlm, mod_scalars::lobc, mod_scalars::lstflux, mod_parallel::master, mod_param::n, mod_scalars::nbrec, mod_netcdf::netcdf_sync(), mod_scalars::nfrec, mod_scalars::ninner, mod_scalars::noerror, mod_param::nt, mod_ocean::ocean, mod_param::r2dvar, mod_param::r3dvar, mod_iounits::sourcefile, mod_iounits::stdout, mod_scalars::time, mod_param::u2dvar, mod_param::u3dvar, mod_param::v2dvar, mod_param::v3dvar, mod_ncparam::vname, mod_fourdvar::zlanczos_coef, mod_fourdvar::zlanczos_err, and mod_fourdvar::zlanczos_inv.

Referenced by wrt_error().

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

◆ wrt_error_pio()

subroutine, private wrt_error_mod::wrt_error_pio ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) kout,
integer, intent(in) nout,
integer, intent(in) lbij,
integer, intent(in) ubij,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj )
private

Definition at line 610 of file wrt_error.F.

615!***********************************************************************
616!
618!
619! Imported variable declarations.
620!
621 integer, intent(in) :: ng, tile, kout, nout
622# ifdef ADJUST_BOUNDARY
623 integer, intent(in) :: LBij, UBij
624# endif
625 integer, intent(in) :: LBi, UBi, LBj, UBj
626!
627! Local variable declarations.
628!
629 integer :: Fcount, i, ifield, j, status
630# ifdef SOLVE3D
631 integer :: itrc, k
632# endif
633!
634 real(dp) :: scale
635!
636 character (len=*), parameter :: MyFile = &
637 & __FILE__//", wrt_error_nf90"
638!
639 TYPE (IO_desc_t), pointer :: ioDesc
640!
641 sourcefile=myfile
642!
643!-----------------------------------------------------------------------
644! Write out full posterior error covariance (diagonal) matrix.
645!-----------------------------------------------------------------------
646!
647 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
648!
649! Set time record index.
650!
651 err(ng)%Rindex=err(ng)%Rindex+1
652 fcount=err(ng)%Fcount
653 err(ng)%Nrec(fcount)=err(ng)%Nrec(fcount)+1
654!
655! Report.
656!
657# ifdef SOLVE3D
658# ifdef NESTING
659 IF (master) WRITE (stdout,20) kout, nout, err(ng)%Rindex, ng
660# else
661 IF (master) WRITE (stdout,20) kout, nout, err(ng)%Rindex
662# endif
663# else
664# ifdef NESTING
665 IF (master) WRITE (stdout,20) kout, err(ng)%Rindex, ng
666# else
667 IF (master) WRITE (stdout,20) kout, err(ng)%Rindex
668# endif
669# endif
670!
671! Write out model time (s).
672!
673 CALL pio_netcdf_put_fvar (ng, itlm, err(ng)%name, &
674 & trim(vname(1,idtime)), time(ng:), &
675 & (/err(ng)%Rindex/), (/1/), &
676 & piofile = err(ng)%pioFile, &
677 & piovar = err(ng)%pioVar(idtime)%vd)
678 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
679!
680! Write out inner-loop Lanczos vectors tridiagonal system.
681!
682 CALL pio_netcdf_put_fvar (ng, itlm, err(ng)%name, &
683 & 'zLanczos_coef', zlanczos_coef, &
684 & (/1,1/), (/ninner,ninner/), &
685 & piofile = err(ng)%pioFile)
686 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
687
688 CALL pio_netcdf_put_fvar (ng, itlm, err(ng)%name, &
689 & 'zLanczos_inv', zlanczos_inv, &
690 & (/1,1/), (/ninner,ninner/), &
691 & piofile = err(ng)%pioFile)
692 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
693
694 CALL pio_netcdf_put_fvar (ng, itlm, err(ng)%name, &
695 & 'zLanczos_err', zlanczos_err, &
696 & (/1,1/), (/ninner,ninner/), &
697 & piofile = err(ng)%pioFile)
698 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
699!
700! Write out free-surface error variance.
701!
702 scale=1.0_dp
703 IF (err(ng)%pioVar(idfsur)%dkind.eq.pio_double) THEN
704 iodesc => iodesc_dp_r2dvar(ng)
705 ELSE
706 iodesc => iodesc_sp_r2dvar(ng)
707 END IF
708!
709 status=nf_fwrite2d(ng, itlm, err(ng)%pioFile, idfsur, &
710 & err(ng)%pioVar(idfsur), &
711 & err(ng)%Rindex, &
712 & iodesc, &
713 & lbi, ubi, lbj, ubj, scale, &
714# ifdef MASKING
715 & grid(ng) % rmask, &
716# endif
717 & ocean(ng)% tl_zeta(:,:,kout))
718 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
719 IF (master) THEN
720 WRITE (stdout,20) trim(vname(1,idfsur)), err(ng)%Rindex
721 END IF
722 exit_flag=3
723 ioerror=status
724 RETURN
725 END IF
726
727# ifdef ADJUST_BOUNDARY
728!
729! Write out free-surface open boundaries error variance.
730!
731 IF (any(lobc(:,isfsur,ng))) THEN
732 scale=1.0_dp
733 IF (err(ng)%pioVar(idsbry(isfsur))%dkind.eq.pio_double) THEN
734 iodesc => iodesc_dp_r2dobc(ng)
735 ELSE
736 iodesc => iodesc_sp_r2dobc(ng)
737 END IF
738!
739 status=nf_fwrite2d_bry(ng, itlm, err(ng)%name, &
740 & err(ng)%pioFile, &
741 & vname(1,idsbry(isfsur)), &
742 & err(ng)%pioVar(idsbry(isfsur)), &
743 & err(ng)%Rindex, &
744 & iodesc, &
745 & lbij, ubij, nbrec(ng), scale, &
746 & boundary(ng) % tl_zeta_obc(lbij:,:,:, &
747 & kout))
748 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
749 IF (master) THEN
750 WRITE (stdout,20) trim(vname(1,idsbry(isfsur))), &
751 & err(ng)%Rindex
752 END IF
753 exit_flag=3
754 ioerror=status
755 RETURN
756 END IF
757 END IF
758# endif
759!
760! Write out 2D U-momentum component error variance.
761!
762 scale=1.0_dp
763 IF (err(ng)%pioVar(idubar)%dkind.eq.pio_double) THEN
764 iodesc => iodesc_dp_u2dvar(ng)
765 ELSE
766 iodesc => iodesc_sp_u2dvar(ng)
767 END IF
768!
769 status=nf_fwrite2d(ng, itlm, err(ng)%pioFile, idubar, &
770 & err(ng)%pioVar(idubar), &
771 & err(ng)%Rindex, &
772 & iodesc, &
773 & lbi, ubi, lbj, ubj, scale, &
774# ifdef MASKING
775 & grid(ng) % umask_full, &
776# endif
777 & ocean(ng) % tl_ubar(:,:,kout))
778 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
779 IF (master) THEN
780 WRITE (stdout,20) trim(vname(1,idubar)), err(ng)%Rindex
781 END IF
782 exit_flag=3
783 ioerror=status
784 RETURN
785 END IF
786
787# ifdef ADJUST_BOUNDARY
788!
789! Write out 2D U-momentum component open boundaries error variance.
790!
791 IF (any(lobc(:,isubar,ng))) THEN
792 scale=1.0_dp
793 IF (err(ng)%pioVar(idsbry(isubar))%dkind.eq.pio_double) THEN
794 iodesc => iodesc_dp_u2dobc(ng)
795 ELSE
796 iodesc => iodesc_sp_u2dobc(ng)
797 END IF
798!
799 status=nf_fwrite2d_bry(ng, itlm, err(ng)%name, &
800 & err(ng)%pioFile, &
801 & vname(1,idsbry(isubar)), &
802 & err(ng)%pioVar(idsbry(isubar)), &
803 & err(ng)%Rindex, &
804 & iodesc, &
805 & lbij, ubij, nbrec(ng), scale, &
806 & boundary(ng) % tl_ubar_obc(lbij:,:,:, &
807 & kout))
808 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
809 IF (master) THEN
810 WRITE (stdout,20) trim(vname(1,idsbry(isubar))), &
811 & err(ng)%Rindex
812 END IF
813 exit_flag=3
814 ioerror=status
815 RETURN
816 END IF
817 END IF
818# endif
819!
820! Write out 2D V-momentum component error variance.
821!
822 scale=1.0_dp
823 IF (err(ng)%pioVar(idvbar)%dkind.eq.pio_double) THEN
824 iodesc => iodesc_dp_v2dvar(ng)
825 ELSE
826 iodesc => iodesc_sp_v2dvar(ng)
827 END IF
828!
829 status=nf_fwrite2d(ng, itlm, err(ng)%pioFile, idvbar, &
830 & err(ng)%pioVar(idvbar), &
831 & err(ng)%Rindex, &
832 & iodesc, &
833 & lbi, ubi, lbj, ubj, scale, &
834# ifdef MASKING
835 & grid(ng) % vmask_full, &
836# endif
837 & ocean(ng) % tl_vbar(:,:,kout))
838 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
839 IF (master) THEN
840 WRITE (stdout,20) trim(vname(1,idvbar)), err(ng)%Rindex
841 END IF
842 exit_flag=3
843 ioerror=status
844 RETURN
845 END IF
846
847# ifdef ADJUST_BOUNDARY
848!
849! Write out 2D V-momentum component open boundaries error variance.
850!
851 IF (any(lobc(:,isvbar,ng))) THEN
852 scale=1.0_dp
853 IF (err(ng)%pioVar(idsbry(isvbar))%dkind.eq.pio_double) THEN
854 iodesc => iodesc_dp_v2dobc(ng)
855 ELSE
856 iodesc => iodesc_sp_v2dobc(ng)
857 END IF
858!
859 status=nf_fwrite2d_bry(ng, itlm, err(ng)%name, &
860 & err(ng)%pioFile, &
861 & vname(1,idsbry(isvbar)), &
862 & err(ng)%pioVar(idsbry(isvbar)), &
863 & err(ng)%Rindex, &
864 & iodesc, &
865 & lbij, ubij, nbrec(ng), scale, &
866 & boundary(ng) % tl_vbar_obc(lbij:,:,:, &
867 & kout))
868 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
869 IF (master) THEN
870 WRITE (stdout,20) trim(vname(1,idsbry(isvbar))), &
871 & err(ng)%Rindex
872 END IF
873 exit_flag=3
874 ioerror=status
875 RETURN
876 END IF
877 END IF
878# endif
879
880# ifdef SOLVE3D
881!
882! Write out 3D U-momentum component error variance.
883!
884 scale=1.0_dp
885 IF (err(ng)%pioVar(iduvel)%dkind.eq.pio_double) THEN
886 iodesc => iodesc_dp_u3dvar(ng)
887 ELSE
888 iodesc => iodesc_sp_u3dvar(ng)
889 END IF
890!
891 status=nf_fwrite3d(ng, itlm, err(ng)%pioFile, iduvel, &
892 & err(ng)%pioVar(iduvel), &
893 & err(ng)%Rindex, &
894 & iodesc, &
895 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
896# ifdef MASKING
897 & grid(ng) % umask_full, &
898# endif
899 & ocean(ng) % tl_u(:,:,:,nout))
900 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
901 IF (master) THEN
902 WRITE (stdout,20) trim(vname(1,iduvel)), err(ng)%Rindex
903 END IF
904 exit_flag=3
905 ioerror=status
906 RETURN
907 END IF
908
909# ifdef ADJUST_BOUNDARY
910!
911! Write out 3D U-momentum component open boundaries error variance.
912!
913 IF (any(lobc(:,isuvel,ng))) THEN
914 scale=1.0_dp
915 IF (err(ng)%pioVar(idsbry(isuvel))%dkind.eq.pio_double) THEN
916 iodesc => iodesc_dp_u3dobc(ng)
917 ELSE
918 iodesc => iodesc_sp_u3dobc(ng)
919 END IF
920!
921 status=nf_fwrite3d_bry(ng, itlm, err(ng)%name, &
922 & err(ng)%pioFile, &
923 & vname(1,idsbry(isuvel)), &
924 & err(ng)%pioVar(idsbry(isuvel)), &
925 & err(ng)%Rindex, &
926 & iodesc, &
927 & lbij, ubij, 1, n(ng), nbrec(ng), scale, &
928 & boundary(ng) % tl_u_obc(lbij:,:,:,:, &
929 & nout))
930 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
931 IF (master) THEN
932 WRITE (stdout,20) trim(vname(1,idsbry(isuvel))), &
933 & err(ng)%Rindex
934 END IF
935 exit_flag=3
936 ioerror=status
937 RETURN
938 END IF
939 END IF
940# endif
941!
942! Write out 3D V-momentum component error variance.
943!
944 scale=1.0_dp
945 IF (err(ng)%pioVar(idfsur)%dkind.eq.pio_double) THEN
946 iodesc => iodesc_dp_v3dvar(ng)
947 ELSE
948 iodesc => iodesc_sp_v3dvar(ng)
949 END IF
950!
951 status=nf_fwrite3d(ng, itlm, err(ng)%pioFile, idvvel, &
952 & err(ng)%pioVar(idvvel), &
953 & err(ng)%Rindex, &
954 & iodesc, &
955 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
956# ifdef MASKING
957 & grid(ng) % vmask_full, &
958# endif
959 & ocean(ng) % tl_v(:,:,:,nout))
960 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
961 IF (master) THEN
962 WRITE (stdout,20) trim(vname(1,idvvel)), err(ng)%Rindex
963 END IF
964 exit_flag=3
965 ioerror=status
966 RETURN
967 END IF
968
969# ifdef ADJUST_BOUNDARY
970!
971! Write out 3D V-momentum component open boundaries error variance.
972!
973 IF (any(lobc(:,isvvel,ng))) THEN
974 scale=1.0_dp
975 IF (err(ng)%pioVar(idsbry(isvvel))%dkind.eq.pio_double) THEN
976 iodesc => iodesc_dp_v3dobc(ng)
977 ELSE
978 iodesc => iodesc_sp_v3dobc(ng)
979 END IF
980!
981 status=nf_fwrite3d_bry(ng, itlm, err(ng)%name, &
982 & err(ng)%pioFile, &
983 & vname(1,idsbry(isvvel)), &
984 & err(ng)%pioVar(idsbry(isvvel)), &
985 & err(ng)%Rindex, &
986 & iodesc, &
987 & lbij, ubij, 1, n(ng), nbrec(ng), scale, &
988 & boundary(ng) % tl_v_obc(lbij:,:,:,:, &
989 & nout))
990 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
991 IF (master) THEN
992 WRITE (stdout,20) trim(vname(1,idsbry(isvvel))), &
993 & err(ng)%Rindex
994 END IF
995 exit_flag=3
996 ioerror=status
997 RETURN
998 END IF
999 END IF
1000# endif
1001!
1002! Write out tracer type variables error variance.
1003!
1004 DO itrc=1,nt(ng)
1005 scale=1.0_dp
1006 IF (err(ng)%pioTrc(itrc)%dkind.eq.pio_double) THEN
1007 iodesc => iodesc_dp_r3dvar(ng)
1008 ELSE
1009 iodesc => iodesc_sp_r3dvar(ng)
1010 END IF
1011!
1012 status=nf_fwrite3d(ng, itlm, err(ng)%pioFile, idtvar(itrc), &
1013 & err(ng)%pioTrc(itrc), &
1014 & err(ng)%Rindex, &
1015 & iodesc, &
1016 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
1017# ifdef MASKING
1018 & grid(ng) % rmask, &
1019# endif
1020 & ocean(ng) % tl_t(:,:,:,nout,itrc))
1021 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1022 IF (master) THEN
1023 WRITE (stdout,20) trim(vname(1,idtvar(itrc))), &
1024 & err(ng)%Rindex
1025 END IF
1026 exit_flag=3
1027 ioerror=status
1028 RETURN
1029 END IF
1030 END DO
1031
1032# ifdef ADJUST_BOUNDARY
1033!
1034! Write out tracers open boundaries error variance.
1035!
1036 DO itrc=1,nt(ng)
1037 IF (any(lobc(:,istvar(itrc),ng))) THEN
1038 scale=1.0_dp
1039 ifield=idsbry(istvar(itrc))
1040 IF (err(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1041 iodesc => iodesc_dp_r3dobc(ng)
1042 ELSE
1043 iodesc => iodesc_sp_r3dobc(ng)
1044 END IF
1045!
1046 status=nf_fwrite3d_bry(ng, itlm, err(ng)%name, &
1047 & err(ng)%pioFile, &
1048 & vname(1,ifield), &
1049 & err(ng)%pioVar(ifield), &
1050 & err(ng)%Rindex, &
1051 & iodesc, &
1052 & lbij, ubij, 1, n(ng), nbrec(ng), &
1053 & scale, &
1054 & boundary(ng) % tl_t_obc(lbij:,:,:,:, &
1055 & nout,itrc))
1056 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1057 IF (master) THEN
1058 WRITE (stdout,20) trim(vname(1,ifield), err(ng)%Rindex
1059 END IF
1060 exit_flag=3
1061 ioerror=status
1062 RETURN
1063 END IF
1064 END IF
1065 END DO
1066# endif
1067
1068# ifdef ADJUST_STFLUX
1069!
1070! Write out surface net tracers fluxes error variance. Notice that
1071! fluxes have their own fixed time-dimension (of size Nfrec) to allow
1072! 4DVar adjustments at other times in addition to initialization time.
1073!
1074 DO itrc=1,nt(ng)
1075 IF (lstflux(itrc,ng)) THEN
1076 scale=1.0_dp
1077 IF (err(ng)%pioVar(idtsur(itrc))%dkind.eq.pio_double) THEN
1078 iodesc => iodesc_dp_r2dfrc(ng)
1079 ELSE
1080 iodesc => iodesc_sp_r2dfrc(ng)
1081 END IF
1082!
1083 status=nf_fwrite3d(ng, itlm, err(ng)%pioFile, idtsur(itrc), &
1084 & err(ng)%pioVar(idtsur(itrc)), &
1085 & err(ng)%Rindex, &
1086 & iodesc, &
1087 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1088# ifdef MASKING
1089 & grid(ng) % rmask, &
1090# endif
1091 & forces(ng) % tl_tflux(:,:,:,kout,itrc))
1092 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1093 IF (master) THEN
1094 WRITE (stdout,20) trim(vname(1,idtsur(itrc))), &
1095 & err(ng)%Rindex
1096 END IF
1097 exit_flag=3
1098 ioerror=status
1099 RETURN
1100 END IF
1101 END IF
1102 END DO
1103# endif
1104# endif
1105# ifdef ADJUST_WSTRESS
1106!
1107! Write out surface U-momentum stress error variance. Notice that the
1108! stress has its own fixed time-dimension (of size Nfrec) to allow
1109! 4DVar adjustments at other times in addition to initialization time.
1110!
1111 scale=1.0_dp
1112 IF (err(ng)%pioVar(idusms)%dkind.eq.pio_double) THEN
1113 iodesc => iodesc_dp_u2dfrc(ng)
1114 ELSE
1115 iodesc => iodesc_sp_u2dfrc(ng)
1116 END IF
1117!
1118 status=nf_fwrite3d(ng, itlm, err(ng)%pioFile, idusms, &
1119 & err(ng)%pioVar(idusms), &
1120 & err(ng)%Rindex, &
1121 & iodesc, &
1122 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1123# ifdef MASKING
1124 & grid(ng) % umask, &
1125# endif
1126 & forces(ng) % tl_ustr(:,:,:,kout))
1127 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1128 IF (master) THEN
1129 WRITE (stdout,20) trim(vname(1,idusms)), err(ng)%Rindex
1130 END IF
1131 exit_flag=3
1132 ioerror=status
1133 RETURN
1134 END IF
1135!
1136! Write out surface V-momentum stress error variance.
1137!
1138 scale=1.0_dp
1139 IF (err(ng)%pioVar(idvsms)%dkind.eq.pio_double) THEN
1140 iodesc => iodesc_dp_v2dfrc(ng)
1141 ELSE
1142 iodesc => iodesc_sp_v2dfrc(ng)
1143 END IF
1144!
1145 status=nf_fwrite3d(ng, itlm, err(ng)%pioFile, idvsms, &
1146 & err(ng)%pioVar(idvsms), &
1147 & err(ng)%Rindex, &
1148 & iodesc, &
1149 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1150# ifdef MASKING
1151 & grid(ng) % vmask, &
1152# endif
1153 & forces(ng) % tl_vstr(:,:,:,kout))
1154 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1155 IF (master) THEN
1156 WRITE (stdout,20) trim(vname(1,idvsms)), err(ng)%Rindex
1157 END IF
1158 exit_flag=3
1159 ioerror=status
1160 RETURN
1161 END IF
1162# endif
1163!
1164!-----------------------------------------------------------------------
1165! Synchronize posterior error covariance NetCDF file to disk to allow
1166! other processes to access data immediately after it is written.
1167!-----------------------------------------------------------------------
1168!
1169 CALL pio_netcdf_sync (ng, itlm, err(ng)%name, err(ng)%pioFile)
1170 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1171!
1172 10 FORMAT (2x,'WRT_ERROR_PIO - writing error', t42, &
1173# ifdef SOLVE3D
1174# ifdef NESTING
1175 & 'fields (Index=',i1,',',i1,') in record = ',i0,t92,i2.2)
1176# else
1177 & 'fields (Index=',i1,',',i1,') in record = ',i0)
1178# endif
1179# else
1180# ifdef NESTING
1181 & 'fields (Index=',i1,') in record = ',i0,t92,i2.2)
1182# else
1183 & 'fields (Index=',i1,') in record = ',i0)
1184# endif
1185# endif
1186 20 FORMAT (/,' WRT_ERROR_PIO - error while writing variable: ',a, &
1187 & /,18x,'into 4DVar error NetCDF file for time record: ',i0)
1188!
1189 RETURN
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dobc
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_v3dobc
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_u3dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dobc
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dobc
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dobc
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dfrc
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_dp_r3dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_v3dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dfrc

References mod_boundary::boundary, mod_iounits::err, mod_scalars::exit_flag, mod_forces::forces, strings_mod::founderror(), mod_grid::grid, mod_ncparam::idfsur, mod_ncparam::idsbry, 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_pio_netcdf::iodesc_dp_r2dfrc, mod_pio_netcdf::iodesc_dp_r2dobc, mod_pio_netcdf::iodesc_dp_r2dvar, mod_pio_netcdf::iodesc_dp_r3dobc, mod_pio_netcdf::iodesc_dp_r3dvar, mod_pio_netcdf::iodesc_dp_u2dfrc, mod_pio_netcdf::iodesc_dp_u2dobc, mod_pio_netcdf::iodesc_dp_u2dvar, mod_pio_netcdf::iodesc_dp_u3dobc, mod_pio_netcdf::iodesc_dp_u3dvar, mod_pio_netcdf::iodesc_dp_v2dfrc, mod_pio_netcdf::iodesc_dp_v2dobc, mod_pio_netcdf::iodesc_dp_v2dvar, mod_pio_netcdf::iodesc_dp_v3dobc, mod_pio_netcdf::iodesc_dp_v3dvar, mod_pio_netcdf::iodesc_sp_r2dfrc, mod_pio_netcdf::iodesc_sp_r2dobc, mod_pio_netcdf::iodesc_sp_r2dvar, mod_pio_netcdf::iodesc_sp_r3dobc, mod_pio_netcdf::iodesc_sp_r3dvar, mod_pio_netcdf::iodesc_sp_u2dfrc, mod_pio_netcdf::iodesc_sp_u2dobc, mod_pio_netcdf::iodesc_sp_u2dvar, mod_pio_netcdf::iodesc_sp_u3dobc, mod_pio_netcdf::iodesc_sp_u3dvar, mod_pio_netcdf::iodesc_sp_v2dfrc, mod_pio_netcdf::iodesc_sp_v2dobc, mod_pio_netcdf::iodesc_sp_v2dvar, mod_pio_netcdf::iodesc_sp_v3dobc, mod_pio_netcdf::iodesc_sp_v3dvar, mod_iounits::ioerror, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_param::itlm, mod_scalars::lobc, mod_scalars::lstflux, mod_parallel::master, mod_param::n, mod_scalars::nbrec, mod_scalars::nfrec, mod_scalars::ninner, mod_scalars::noerror, mod_param::nt, mod_ocean::ocean, mod_pio_netcdf::pio_netcdf_sync(), mod_iounits::sourcefile, mod_iounits::stdout, mod_scalars::time, mod_ncparam::vname, mod_fourdvar::zlanczos_coef, mod_fourdvar::zlanczos_err, and mod_fourdvar::zlanczos_inv.

Referenced by wrt_error().

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