ROMS
Loading...
Searching...
No Matches
state_join.F
Go to the documentation of this file.
1#include "cppdefs.h"
3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine reads and writes ROMS NLM background state solution !
12! for specified input and output records. Primarily, it is used to !
13! concatenate the NLM solution from concurrent interval window !
14! history NetCDF files. The resulting solution is used to linearize !
15! the tangent linear and adjoint kernels. !
16! !
17! On Input: !
18! !
19! ng Nested grid number (integer) !
20! model Calling model identifier (integer) !
21! InpName Input history filename (string) !
22! InpName Output history filename (string) !
23! InpStrRec Starting Input history time record (integer) !
24! InpEndRec Ending Input history time record (integer) !
25! OutStrRec Starting Output history time record (integer) !
26! !
27!=======================================================================
28!
29 USE mod_param
30 USE mod_parallel
31 USE mod_coupling
32 USE mod_grid
33 USE mod_iounits
34 USE mod_mixing
35 USE mod_ncparam
36 USE mod_netcdf
37 USE mod_ocean
38 USE mod_scalars
39!
40 USE dateclock_mod, ONLY : time_string
41 USE nf_fread2d_mod, ONLY : nf_fread2d
43#ifdef SOLVE3D
44 USE nf_fread3d_mod, ONLY : nf_fread3d
46#endif
47 USE strings_mod, ONLY : founderror
48!
49 implicit none
50!
51 PUBLIC :: state_join
52 PRIVATE :: state_join_nf90
53# if defined PIO_LIB && defined DISTRIBUTE
54 PRIVATE :: state_join_pio
55# endif
56!
57 CONTAINS
58!
59!***********************************************************************
60 SUBROUTINE state_join (ng, tile, model, InpName, OutName, IOtype, &
61 & InpStrRec, InpEndRec, OutStrRec)
62!***********************************************************************
63!
64! Imported variable declarations.
65!
66 integer, intent(in) :: ng, tile, model, iotype
67 integer, intent(in) :: inpstrrec, inpendrec
68
69 integer, intent(inout) :: outstrrec
70!
71 character (len=*) :: inpname, outname
72!
73! Local variable declarations.
74!
75 integer :: lbi, ubi, lbj, ubj
76!
77 character (len=*), parameter :: myfile = &
78 & __FILE__
79!
80!-----------------------------------------------------------------------
81! Create a new history file according to IO type.
82!-----------------------------------------------------------------------
83!
84 lbi=bounds(ng)%LBi(tile)
85 ubi=bounds(ng)%UBi(tile)
86 lbj=bounds(ng)%LBj(tile)
87 ubj=bounds(ng)%UBj(tile)
88!
89 SELECT CASE (iotype)
90 CASE (io_nf90)
91 CALL state_join_nf90 (ng, tile, model, inpname, outname, &
92 & inpstrrec, inpendrec, outstrrec, &
93 & lbi, ubi, lbj, ubj)
94
95#if defined PIO_LIB && defined DISTRIBUTE
96 CASE (io_pio)
97 CALL state_join_pio (ng, tile, model, inpname, outname, &
98 & inpstrrec, inpendrec, outstrrec, &
99 & lbi, ubi, lbj, ubj)
100#endif
101 CASE DEFAULT
102 IF (master) WRITE (stdout,10) iotype
103 exit_flag=2
104 END SELECT
105 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
106!
107 10 FORMAT (' STATE_JOIN - Illegal input file type, io_type = ',i0, &
108 & /,14x,'Check KeyWord ''INP_LIB'' in ''roms.in''.')
109!
110 RETURN
111 END SUBROUTINE state_join
112!
113!***********************************************************************
114 SUBROUTINE state_join_nf90 (ng, tile, model, InpName, OutName, &
115 & InpStrRec, InpEndRec, OutStrRec, &
116 & LBi, UBi, LBj, UBj)
117!***********************************************************************
118!
119! Imported variable declarations.
120!
121 integer, intent(in) :: ng, tile, model
122 integer, intent(in) :: inpstrrec, inpendrec
123 integer, intent(in) :: lbi, ubi, lbj, ubj
124
125 integer, intent(inout) :: outstrrec
126!
127 character (len=*) :: inpname, outname
128!
129! Local variable declarations.
130!
131 integer :: inpid, inprec, outid, outrec, tindex
132 integer :: gtype, i, ic, itrc, status
133
134 integer :: vsize(4)
135!
136 real(r8) :: fmin, fmax
137 real(dp) :: fscl, stime
138!
139 character (len=15) :: tstring
140 character (len=22) :: t_code
141
142 character (len=*), parameter :: myfile = &
143 & __FILE__//", state_join_nf90"
144!
145 sourcefile=myfile
146!
147!-----------------------------------------------------------------------
148! Process forward background solution.
149!-----------------------------------------------------------------------
150!
151! Open Input NetCDF file for reading.
152!
153 CALL netcdf_open (ng, model, inpname, 0, inpid)
154 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
155 WRITE (stdout,10) trim(inpname)
156 RETURN
157 END IF
158!
159! Open Output NetCDF file for writing.
160!
161 CALL netcdf_open (ng, model, outname, 1, outid)
162 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
163 WRITE (stdout,10) trim(inpname)
164 RETURN
165 END IF
166!
167! Inquire about the input variables.
168!
169 CALL netcdf_inq_var (ng, model, inpname, &
170 & ncid = inpid)
171 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
172!
173! Set Vsize to zero to deativate interpolation of input data to model
174! grid in "nf_fread2d".
175!
176 DO i=1,4
177 vsize(i)=0
178 END DO
179!
180! Scan the variable list and read in needed variables.
181!
182 outrec=outstrrec-1
183 tindex=1
184 fscl=1.0_dp
185!
186 rec_loop : DO inprec=inpstrrec,inpendrec
187 outrec=outrec+1
188 var_loop : DO i=1,n_var
189!
190! Time.
191!
192 check1 : IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
193 CALL netcdf_get_time (ng, model, inpname, &
194 & trim(var_name(i)), &
195 & rclock%DateNumber, stime, &
196 & ncid = inpid, &
197 & start = (/inprec/), &
198 & total = (/1/))
199 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
200 IF (master) THEN
201 WRITE (stdout,30) trim(vname(1,idtime)), inprec, &
202 & trim(inpname)
203 END IF
204 exit_flag=2
205 RETURN
206 END IF
207!
208 CALL netcdf_put_fvar (ng, model, outname, &
209 & trim(var_name(i)), stime, &
210 & (/outrec/), (/1/), &
211 & ncid = outid, &
212 & varid = var_id(i))
213 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
214 IF (master) THEN
215 WRITE (stdout,40) trim(vname(1,idtime)), outrec, &
216 & trim(outname)
217 END IF
218 exit_flag=3
219 RETURN
220 ELSE
221 CALL time_string (stime, t_code)
222 WRITE (tstring,'(f15.4)') stime*sec2day
223 IF (master) THEN
224 WRITE (stdout,20) t_code, &
225 & ng, trim(adjustl(tstring)), inprec, &
226 & tindex, trim(inpname), &
227 & ng, trim(adjustl(tstring)), outrec, &
228 & tindex, trim(outname)
229 END IF
230 END IF
231!
232! Free-surface.
233!
234 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
235 gtype=r2dvar
236 status=nf_fread2d(ng, model, inpname, inpid, &
237 & var_name(i), var_id(i), &
238 & inprec, gtype, vsize, &
239 & lbi, ubi, lbj, ubj, &
240 & fscl, fmin, fmax, &
241#ifdef MASKING
242 & grid(ng) % rmask, &
243#endif
244 & ocean(ng) % zeta(:,:,tindex))
245 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
246 IF (master) THEN
247 WRITE (stdout,30) trim(vname(1,idfsur)), inprec, &
248 & trim(inpname)
249 END IF
250 exit_flag=2
251 ioerror=status
252 RETURN
253 END IF
254!
255 status=nf_fwrite2d(ng, model, outid, idfsur, &
256 & var_id(i), outrec, gtype, &
257 & lbi, ubi, lbj, ubj, fscl, &
258#ifdef MASKING
259 & grid(ng) % rmask, &
260#endif
261 & ocean(ng) % zeta(:,:,tindex), &
262 & setfillval = .false.)
263 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
264 IF (master) THEN
265 WRITE (stdout,40) trim(vname(1,idfsur)), outrec, &
266 & trim(outname)
267 END IF
268 exit_flag=3
269 ioerror=status
270 RETURN
271 ELSE
272 IF (master) THEN
273 WRITE (stdout,50) trim(vname(2,idfsur)), fmin, fmax
274 END IF
275 END IF
276
277#ifdef FORWARD_RHS
278!
279! Free-surface equation right-hand-side term.
280!
281 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrzet))) THEN
282 gtype=r2dvar
283 status=nf_fread2d(ng, model, inpname, inpid, &
284 & var_name(i), var_id(i), &
285 & inprec, gtype, vsize, &
286 & lbi, ubi, lbj, ubj, &
287 & fscl, fmin, fmax, &
288# ifdef MASKING
289 & grid(ng) % rmask, &
290# endif
291 & ocean(ng) % rzeta(:,:,tindex))
292 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
293 IF (master) THEN
294 WRITE (stdout,30) trim(vname(1,idrzet)), inprec, &
295 & trim(inpname)
296 END IF
297 exit_flag=2
298 ioerror=status
299 RETURN
300 END IF
301!
302 status=nf_fwrite2d(ng, model, outid, idrzet, &
303 & var_id(i), outrec, gtype, &
304 & lbi, ubi, lbj, ubj, fscl, &
305# ifdef MASKING
306 & grid(ng) % rmask, &
307# endif
308 & ocean(ng) % rzeta(:,:,tindex))
309 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
310 IF (master) THEN
311 WRITE (stdout,40) trim(vname(1,idrzet)), outrec, &
312 & trim(outname)
313 END IF
314 exit_flag=3
315 ioerror=status
316 RETURN
317 ELSE
318 IF (master) THEN
319 WRITE (stdout,50) trim(vname(2,idrzet)), fmin, fmax
320 END IF
321 END IF
322#endif
323!
324! 2D U-momentum component.
325!
326 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
327 gtype=u2dvar
328 status=nf_fread2d(ng, model, inpname, inpid, &
329 & var_name(i), var_id(i), &
330 & inprec, gtype, vsize, &
331 & lbi, ubi, lbj, ubj, &
332 & fscl, fmin, fmax, &
333#ifdef MASKING
334 & grid(ng) % umask_full, &
335#endif
336 & ocean(ng) % ubar(:,:,tindex))
337 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
338 IF (master) THEN
339 WRITE (stdout,30) trim(vname(1,idubar)), inprec, &
340 & trim(inpname)
341 END IF
342 exit_flag=2
343 ioerror=status
344 RETURN
345 END IF
346!
347 status=nf_fwrite2d(ng, model, outid, idubar, &
348 & var_id(i), outrec, gtype, &
349 & lbi, ubi, lbj, ubj, fscl, &
350#ifdef MASKING
351 & grid(ng) % umask_full, &
352#endif
353 & ocean(ng) % ubar(:,:,tindex))
354 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
355 IF (master) THEN
356 WRITE (stdout,40) trim(vname(1,idubar)), outrec, &
357 & trim(outname)
358 END IF
359 exit_flag=3
360 ioerror=status
361 RETURN
362 ELSE
363 IF (master) THEN
364 WRITE (stdout,50) trim(vname(2,idubar)), fmin, fmax
365 END IF
366 END IF
367
368#ifdef FORWARD_RHS
369!
370! 2D U-equation right-hand-side term.
371!
372 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idru2d))) THEN
373 gtype=u2dvar
374 status=nf_fread2d(ng, model, inpname, inpid, &
375 & var_name(i), var_id(i), &
376 & inprec, gtype, vsize, &
377 & lbi, ubi, lbj, ubj, &
378 & fscl, fmin, fmax, &
379# ifdef MASKING
380 & grid(ng) % umask_full, &
381# endif
382 & ocean(ng) % rubar(:,:,tindex))
383 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
384 IF (master) THEN
385 WRITE (stdout,30) trim(vname(1,idru2d)), inprec, &
386 & trim(inpname)
387 END IF
388 exit_flag=2
389 ioerror=status
390 RETURN
391 END IF
392!
393 status=nf_fwrite2d(ng, model, outid, idru2d, &
394 & var_id(i), outrec, gtype, &
395 & lbi, ubi, lbj, ubj, fscl, &
396# ifdef MASKING
397 & grid(ng) % umask_full, &
398# endif
399 & ocean(ng) % rubar(:,:,tindex))
400 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
401 IF (master) THEN
402 WRITE (stdout,40) trim(vname(1,idru2d)), outrec, &
403 & trim(outname)
404 END IF
405 exit_flag=3
406 ioerror=status
407 RETURN
408 ELSE
409 IF (master) THEN
410 WRITE (stdout,50) trim(vname(2,idru2d)), fmin, fmax
411 END IF
412 END IF
413#endif
414!
415! 2D V-momentum component.
416!
417 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
418 gtype=v2dvar
419 status=nf_fread2d(ng, model, inpname, inpid, &
420 & var_name(i), var_id(i), &
421 & inprec, gtype, vsize, &
422 & lbi, ubi, lbj, ubj, &
423 & fscl, fmin, fmax, &
424#ifdef MASKING
425 & grid(ng) % vmask_full, &
426#endif
427 & ocean(ng) % vbar(:,:,tindex))
428 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
429 IF (master) THEN
430 WRITE (stdout,30) trim(vname(1,idvbar)), inprec, &
431 & trim(inpname)
432 END IF
433 exit_flag=2
434 ioerror=status
435 RETURN
436 END IF
437!
438 status=nf_fwrite2d(ng, model, outid, idvbar, &
439 & var_id(i), outrec, gtype, &
440 & lbi, ubi, lbj, ubj, fscl, &
441#ifdef MASKING
442 & grid(ng) % vmask_full, &
443#endif
444 & ocean(ng) % vbar(:,:,tindex))
445 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
446 IF (master) THEN
447 WRITE (stdout,40) trim(vname(1,idvbar)), outrec, &
448 & trim(outname)
449 END IF
450 exit_flag=3
451 ioerror=status
452 RETURN
453 ELSE
454 IF (master) THEN
455 WRITE (stdout,50) trim(vname(2,idvbar)), fmin, fmax
456 END IF
457 END IF
458
459#ifdef FORWARD_RHS
460!
461! 2D V-equation right-hand-side term.
462!
463 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrv2d))) THEN
464 gtype=v2dvar
465 status=nf_fread2d(ng, model, inpname, inpid, &
466 & var_name(i), var_id(i), &
467 & inprec, gtype, vsize, &
468 & lbi, ubi, lbj, ubj, &
469 & fscl, fmin, fmax, &
470# ifdef MASKING
471 & grid(ng) % vmask_full, &
472# endif
473 & ocean(ng) % rvbar(:,:,tindex))
474 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
475 IF (master) THEN
476 WRITE (stdout,30) trim(vname(1,idrv2d)), inprec, &
477 & trim(inpname)
478 END IF
479 exit_flag=2
480 ioerror=status
481 RETURN
482 END IF
483!
484 status=nf_fwrite2d(ng, model, outid, idrv2d, &
485 & var_id(i), outrec, gtype, &
486 & lbi, ubi, lbj, ubj, fscl, &
487# ifdef MASKING
488 & grid(ng) % vmask_full, &
489# endif
490 & ocean(ng) % rvbar(:,:,tindex))
491 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
492 IF (master) THEN
493 WRITE (stdout,40) trim(vname(1,idrv2d)), outrec, &
494 & trim(outname)
495 END IF
496 exit_flag=3
497 ioerror=status
498 RETURN
499 ELSE
500 IF (master) THEN
501 WRITE (stdout,50) trim(vname(2,idrv2d)), fmin, fmax
502 END IF
503 END IF
504#endif
505#ifdef SOLVE3D
506# ifdef FORWARD_RHS
507!
508! 2D U-equation right-hand-side forcing term.
509!
510 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idruct))) THEN
511 gtype=u2dvar
512 status=nf_fread2d(ng, model, inpname, inpid, &
513 & var_name(i), var_id(i), &
514 & inprec, gtype, vsize, &
515 & lbi, ubi, lbj, ubj, &
516 & fscl, fmin, fmax, &
517# ifdef MASKING
518 & grid(ng) % umask_full, &
519# endif
520 & ocean(ng) % rufrc)
521 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
522 IF (master) THEN
523 WRITE (stdout,30) trim(vname(1,idruct)), inprec, &
524 & trim(inpname)
525 END IF
526 exit_flag=2
527 ioerror=status
528 RETURN
529 END IF
530!
531 status=nf_fwrite2d(ng, model, outid, idruct, &
532 & var_id(i), outrec, gtype, &
533 & lbi, ubi, lbj, ubj, fscl, &
534# ifdef MASKING
535 & grid(ng) % umask_full, &
536# endif
537 & ocean(ng) % rufrc)
538 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
539 IF (master) THEN
540 WRITE (stdout,40) trim(vname(1,idruct)), outrec, &
541 & trim(outname)
542 END IF
543 exit_flag=3
544 ioerror=status
545 RETURN
546 ELSE
547 IF (master) THEN
548 WRITE (stdout,50) trim(vname(2,idruct)), fmin, fmax
549 END IF
550 END IF
551# endif
552!
553! Time-averaged U-flux component for 2D equations.
554!
555 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idufx1))) THEN
556 gtype=u2dvar
557 status=nf_fread2d(ng, model, inpname, inpid, &
558 & var_name(i), var_id(i), &
559 & inprec, gtype, vsize, &
560 & lbi, ubi, lbj, ubj, &
561 & fscl, fmin, fmax, &
562# ifdef MASKING
563 & grid(ng) % umask, &
564# endif
565 & coupling(ng) % DU_avg1)
566 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
567 IF (master) THEN
568 WRITE (stdout,30) trim(vname(1,idufx1)), inprec, &
569 & trim(inpname)
570 END IF
571 exit_flag=2
572 ioerror=status
573 RETURN
574 END IF
575!
576 status=nf_fwrite2d(ng, model, outid, idufx1, &
577 & var_id(i), outrec, gtype, &
578 & lbi, ubi, lbj, ubj, fscl, &
579# ifdef MASKING
580 & grid(ng) % umask, &
581# endif
582 & coupling(ng) % DU_avg1)
583 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
584 IF (master) THEN
585 WRITE (stdout,40) trim(vname(1,idufx1)), outrec, &
586 & trim(outname)
587 END IF
588 exit_flag=3
589 ioerror=status
590 RETURN
591 ELSE
592 IF (master) THEN
593 WRITE (stdout,50) trim(vname(2,idufx1)), fmin, fmax
594 END IF
595 END IF
596!
597! Time-averaged U-flux component for 3D equations.
598!
599 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idufx2))) THEN
600 gtype=u2dvar
601 status=nf_fread2d(ng, model, inpname, inpid, &
602 & var_name(i), var_id(i), &
603 & inprec, gtype, vsize, &
604 & lbi, ubi, lbj, ubj, &
605 & fscl, fmin, fmax, &
606# ifdef MASKING
607 & grid(ng) % umask, &
608# endif
609 & coupling(ng) % DU_avg2)
610 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
611 IF (master) THEN
612 WRITE (stdout,30) trim(vname(1,idufx2)), inprec, &
613 & trim(inpname)
614 END IF
615 exit_flag=2
616 ioerror=status
617 RETURN
618 END IF
619!
620 status=nf_fwrite2d(ng, model, outid, idufx2, &
621 & var_id(i), outrec, gtype, &
622 & lbi, ubi, lbj, ubj, fscl, &
623# ifdef MASKING
624 & grid(ng) % umask, &
625# endif
626 & coupling(ng) % DU_avg2)
627 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
628 IF (master) THEN
629 WRITE (stdout,40) trim(vname(1,idufx2)), outrec, &
630 & trim(outname)
631 END IF
632 exit_flag=3
633 ioerror=status
634 RETURN
635 ELSE
636 IF (master) THEN
637 WRITE (stdout,50) trim(vname(2,idufx2)), fmin, fmax
638 END IF
639 END IF
640
641# ifdef FORWARD_RHS
642!
643! 2D V-equation right-hand-side forcing term.
644!
645 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrvct))) THEN
646 gtype=v2dvar
647 status=nf_fread2d(ng, model, inpname, inpid, &
648 & var_name(i), var_id(i), &
649 & inprec, gtype, vsize, &
650 & lbi, ubi, lbj, ubj, &
651 & fscl, fmin, fmax, &
652# ifdef MASKING
653 & grid(ng) % vmask_full, &
654# endif
655 & ocean(ng) % rvfrc)
656 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
657 IF (master) THEN
658 WRITE (stdout,30) trim(vname(1,idrvct)), inprec, &
659 & trim(inpname)
660 END IF
661 exit_flag=2
662 ioerror=status
663 RETURN
664 END IF
665!
666 status=nf_fwrite2d(ng, model, outid, idrvct, &
667 & var_id(i), outrec, gtype, &
668 & lbi, ubi, lbj, ubj, fscl, &
669# ifdef MASKING
670 & grid(ng) % vmask_full, &
671# endif
672 & ocean(ng) % rvfrc)
673 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
674 IF (master) THEN
675 WRITE (stdout,40) trim(vname(1,idrvct)), outrec, &
676 & trim(outname)
677 END IF
678 exit_flag=3
679 ioerror=status
680 RETURN
681 ELSE
682 IF (master) THEN
683 WRITE (stdout,50) trim(vname(2,idrvct)), fmin, fmax
684 END IF
685 END IF
686# endif
687!
688! Time-averaged V-flux component for 2D equations.
689!
690 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvfx1))) THEN
691 gtype=v2dvar
692 status=nf_fread2d(ng, model, inpname, inpid, &
693 & var_name(i), var_id(i), &
694 & inprec, gtype, vsize, &
695 & lbi, ubi, lbj, ubj, &
696 & fscl, fmin, fmax, &
697# ifdef MASKING
698 & grid(ng) % vmask_full, &
699# endif
700 & coupling(ng) % DV_avg1)
701 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
702 IF (master) THEN
703 WRITE (stdout,30) trim(vname(1,idvfx1)), inprec, &
704 & trim(inpname)
705 END IF
706 exit_flag=2
707 ioerror=status
708 RETURN
709 END IF
710!
711 status=nf_fwrite2d(ng, model, outid, idvfx1, &
712 & var_id(i), outrec, gtype, &
713 & lbi, ubi, lbj, ubj, fscl, &
714# ifdef MASKING
715 & grid(ng) % vmask_full, &
716# endif
717 & coupling(ng) % DV_avg1)
718 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
719 IF (master) THEN
720 WRITE (stdout,40) trim(vname(1,idvfx1)), outrec, &
721 & trim(outname)
722 END IF
723 exit_flag=3
724 ioerror=status
725 RETURN
726 ELSE
727 IF (master) THEN
728 WRITE (stdout,50) trim(vname(2,idvfx1)), fmin, fmax
729 END IF
730 END IF
731!
732! Time-averaged U-flux component for 3D equations.
733!
734 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvfx2))) THEN
735 gtype=v2dvar
736 status=nf_fread2d(ng, model, inpname, inpid, &
737 & var_name(i), var_id(i), &
738 & inprec, gtype, vsize, &
739 & lbi, ubi, lbj, ubj, &
740 & fscl, fmin, fmax, &
741# ifdef MASKING
742 & grid(ng) % vmask_full, &
743# endif
744 & coupling(ng) % DV_avg2)
745 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
746 IF (master) THEN
747 WRITE (stdout,30) trim(vname(1,idvfx2)), inprec, &
748 & trim(inpname)
749 END IF
750 exit_flag=2
751 ioerror=status
752 RETURN
753 END IF
754!
755 status=nf_fwrite2d(ng, model, outid, idvfx2, &
756 & var_id(i), outrec, gtype, &
757 & lbi, ubi, lbj, ubj, fscl, &
758# ifdef MASKING
759 & grid(ng) % vmask_full, &
760# endif
761 & coupling(ng) % DV_avg2)
762 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
763 IF (master) THEN
764 WRITE (stdout,40) trim(vname(1,idvfx2)), outrec, &
765 & trim(outname)
766 END IF
767 exit_flag=3
768 ioerror=status
769 RETURN
770 ELSE
771 IF (master) THEN
772 WRITE (stdout,50) trim(vname(2,idvfx2)), fmin, fmax
773 END IF
774 END IF
775!
776! 3D U-momentum component.
777!
778 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
779 gtype=u3dvar
780 status=nf_fread3d(ng, model, inpname, inpid, &
781 & var_name(i), var_id(i), &
782 & inprec, gtype, vsize, &
783 & lbi, ubi, lbj, ubj, 1, n(ng), &
784 & fscl, fmin, fmax, &
785# ifdef MASKING
786 & grid(ng) % umask_full, &
787# endif
788 & ocean(ng) % u(:,:,:,tindex))
789 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
790 IF (master) THEN
791 WRITE (stdout,30) trim(vname(1,iduvel)), inprec, &
792 & trim(inpname)
793 END IF
794 exit_flag=2
795 ioerror=status
796 RETURN
797 END IF
798!
799 status=nf_fwrite3d(ng, model, outid, iduvel, &
800 & var_id(i), outrec, gtype, &
801 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
802# ifdef MASKING
803 & grid(ng) % umask_full, &
804# endif
805 & ocean(ng) % u(:,:,:,tindex))
806 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
807 IF (master) THEN
808 WRITE (stdout,40) trim(vname(1,iduvel)), outrec, &
809 & trim(outname)
810 END IF
811 exit_flag=3
812 ioerror=status
813 RETURN
814 ELSE
815 IF (master) THEN
816 WRITE (stdout,50) trim(vname(2,iduvel)), fmin, fmax
817 END IF
818 END IF
819
820# ifdef FORWARD_RHS
821!
822! 3D U-momentum right-hand-side term.
823!
824 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idru3d))) THEN
825 gtype=u3dvar
826 status=nf_fread3d(ng, model, inpname, inpid, &
827 & var_name(i), var_id(i), &
828 & inprec, gtype, vsize, &
829 & lbi, ubi, lbj, ubj, 1, n(ng), &
830 & fscl, fmin, fmax, &
831# ifdef MASKING
832 & grid(ng) % umask_full, &
833# endif
834 & ocean(ng) % ru(:,:,:,tindex))
835 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
836 IF (master) THEN
837 WRITE (stdout,30) trim(vname(1,idru3d)), inprec, &
838 & trim(inpname)
839 END IF
840 exit_flag=2
841 ioerror=status
842 RETURN
843 END IF
844!
845 status=nf_fwrite3d(ng, model, outid, idru3d, &
846 & var_id(i), outrec, gtype, &
847 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
848# ifdef MASKING
849 & grid(ng) % umask_full, &
850# endif
851 & ocean(ng) % ru(:,:,:,tindex))
852 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
853 IF (master) THEN
854 WRITE (stdout,40) trim(vname(1,idru3d)), outrec, &
855 & trim(outname)
856 END IF
857 exit_flag=3
858 ioerror=status
859 RETURN
860 ELSE
861 IF (master) THEN
862 WRITE (stdout,50) trim(vname(2,idru3d)), fmin, fmax
863 END IF
864 END IF
865# endif
866!
867! 3D V-momentum component.
868!
869 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
870 gtype=v3dvar
871 status=nf_fread3d(ng, model, inpname, inpid, &
872 & var_name(i), var_id(i), &
873 & inprec, gtype, vsize, &
874 & lbi, ubi, lbj, ubj, 1, n(ng), &
875 & fscl, fmin, fmax, &
876# ifdef MASKING
877 & grid(ng) % vmask_full, &
878# endif
879 & ocean(ng) % v(:,:,:,tindex))
880 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
881 IF (master) THEN
882 WRITE (stdout,30) trim(vname(1,idvvel)), inprec, &
883 & trim(inpname)
884 END IF
885 exit_flag=2
886 ioerror=status
887 RETURN
888 END IF
889!
890 status=nf_fwrite3d(ng, model, outid, idvvel, &
891 & var_id(i), outrec, gtype, &
892 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
893# ifdef MASKING
894 & grid(ng) % vmask_full, &
895# endif
896 & ocean(ng) % v(:,:,:,tindex))
897 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
898 IF (master) THEN
899 WRITE (stdout,40) trim(vname(1,idvvel)), outrec, &
900 & trim(outname)
901 END IF
902 exit_flag=3
903 ioerror=status
904 RETURN
905 ELSE
906 IF (master) THEN
907 WRITE (stdout,50) trim(vname(2,idvvel)), fmin, fmax
908 END IF
909 END IF
910
911# ifdef FORWARD_RHS
912!
913! 3D U-momentum right-hand-side term.
914!
915 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrv3d))) THEN
916 gtype=v3dvar
917 status=nf_fread3d(ng, model, inpname, inpid, &
918 & var_name(i), var_id(i), &
919 & inprec, gtype, vsize, &
920 & lbi, ubi, lbj, ubj, 1, n(ng), &
921 & fscl, fmin, fmax, &
922# ifdef MASKING
923 & grid(ng) % vmask_full, &
924# endif
925 & ocean(ng) % rv(:,:,:,tindex))
926 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
927 IF (master) THEN
928 WRITE (stdout,30) trim(vname(1,idrv3d)), inprec, &
929 & trim(inpname)
930 END IF
931 exit_flag=2
932 ioerror=status
933 RETURN
934 END IF
935!
936 status=nf_fwrite3d(ng, model, outid, idrv3d, &
937 & var_id(i), outrec, gtype, &
938 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
939# ifdef MASKING
940 & grid(ng) % vmask_full, &
941# endif
942 & ocean(ng) % rv(:,:,:,tindex))
943 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
944 IF (master) THEN
945 WRITE (stdout,40) trim(vname(1,idrv3d)), outrec, &
946 & trim(outname)
947 END IF
948 exit_flag=3
949 ioerror=status
950 RETURN
951 ELSE
952 IF (master) THEN
953 WRITE (stdout,50) trim(vname(2,idrv3d)), fmin, fmax
954 END IF
955 END IF
956# endif
957# ifdef FORWARD_MIXING
958!
959! Vertical viscosity.
960!
961 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvis))) THEN
962 gtype=w3dvar
963 status=nf_fread3d(ng, model, inpname, inpid, &
964 & var_name(i), var_id(i), &
965 & inprec, gtype, vsize, &
966 & lbi, ubi, lbj, ubj, 0, n(ng), &
967 & fscl, fmin, fmax, &
968# ifdef MASKING
969 & grid(ng) % rmask, &
970# endif
971 & mixing(ng) % Akv)
972 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
973 IF (master) THEN
974 WRITE (stdout,30) trim(vname(1,idvvis)), inprec, &
975 & trim(inpname)
976 END IF
977 exit_flag=2
978 ioerror=status
979 RETURN
980 END IF
981!
982 status=nf_fwrite3d(ng, model, outid, idvvis, &
983 & var_id(i), outrec, gtype, &
984 & lbi, ubi, lbj, ubj, 0, n(ng), fscl, &
985# ifdef MASKING
986 & grid(ng) % rmask, &
987# endif
988 & mixing(ng) % Akv)
989 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
990 IF (master) THEN
991 WRITE (stdout,40) trim(vname(1,idvvis)), outrec, &
992 & trim(outname)
993 END IF
994 exit_flag=3
995 ioerror=status
996 RETURN
997 ELSE
998 IF (master) THEN
999 WRITE (stdout,50) trim(vname(2,idvvis)), fmin, fmax
1000 END IF
1001 END IF
1002# endif
1003#endif
1004 END IF check1
1005
1006#ifdef SOLVE3D
1007!
1008! Tracer type variables.
1009!
1010 tracer1_loop : DO itrc=1,nt(ng)
1011 gtype=r3dvar
1012 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
1013 status=nf_fread3d(ng, model, inpname, inpid, &
1014 & var_name(i), var_id(i), &
1015 & inprec, gtype, vsize, &
1016 & lbi, ubi, lbj, ubj, 1, n(ng), &
1017 & fscl, fmin, fmax, &
1018# ifdef MASKING
1019 & grid(ng) % rmask, &
1020# endif
1021 & ocean(ng) % t(:,:,:,tindex,itrc))
1022 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1023 IF (master) THEN
1024 WRITE (stdout,30) trim(vname(1,idtvar(itrc))), &
1025 & inprec, trim(inpname)
1026 END IF
1027 exit_flag=2
1028 ioerror=status
1029 RETURN
1030 END IF
1031!
1032 status=nf_fwrite3d(ng, model, outid, idtvar(itrc), &
1033 & var_id(i), outrec, gtype, &
1034 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
1035# ifdef MASKING
1036 & grid(ng) % rmask, &
1037# endif
1038 & ocean(ng) % t(:,:,:,tindex,itrc))
1039 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1040 IF (master) THEN
1041 WRITE (stdout,40) trim(vname(1,idtvar(itrc))), &
1042 & outrec, trim(outname)
1043 END IF
1044 exit_flag=3
1045 ioerror=status
1046 RETURN
1047 ELSE
1048 IF (master) THEN
1049 WRITE (stdout,50) trim(vname(2,idtvar(itrc))), &
1050 & fmin, fmax
1051 END IF
1052 END IF
1053 END IF
1054 END DO tracer1_loop
1055!
1056! Tracer type variables vertical diffusion.
1057!
1058 tracer2_loop : DO itrc=1,nat
1059 gtype=w3dvar
1060 IF (trim(var_name(i)).eq.trim(vname(1,iddiff(itrc)))) THEN
1061 status=nf_fread3d(ng, model, inpname, inpid, &
1062 & var_name(i), var_id(i), &
1063 & inprec, gtype, vsize, &
1064 & lbi, ubi, lbj, ubj, 0, n(ng), &
1065 & fscl, fmin, fmax, &
1066# ifdef MASKING
1067 & grid(ng) % rmask, &
1068# endif
1069 & mixing(ng) % Akt(:,:,:,itrc))
1070 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1071 IF (master) THEN
1072 WRITE (stdout,30) trim(vname(1,iddiff(itrc))), &
1073 & inprec, trim(inpname)
1074 END IF
1075 exit_flag=2
1076 ioerror=status
1077 RETURN
1078 END IF
1079!
1080 status=nf_fwrite3d(ng, model, outid, iddiff(itrc), &
1081 & var_id(i), outrec, gtype, &
1082 & lbi, ubi, lbj, ubj, 0, n(ng), fscl, &
1083# ifdef MASKING
1084 & grid(ng) % rmask, &
1085# endif
1086 & mixing(ng) % Akt(:,:,:,itrc))
1087 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1088 IF (master) THEN
1089 WRITE (stdout,40) trim(vname(1,iddiff(itrc))), &
1090 & outrec, trim(outname)
1091 END IF
1092 exit_flag=3
1093 ioerror=status
1094 RETURN
1095 ELSE
1096 IF (master) THEN
1097 WRITE (stdout,50) trim(vname(2,iddiff(itrc))), &
1098 & fmin, fmax
1099 END IF
1100 END IF
1101 END IF
1102 END DO tracer2_loop
1103#endif
1104
1105 END DO var_loop
1106
1107 END DO rec_loop
1108!
1109! Update output record value.
1110!
1111 outstrrec=outrec
1112!
1113! Close input and output files for compleate synchronization.
1114!
1115 CALL netcdf_close (ng, inlm, inpid, inpname, .false.)
1116 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1117!
1118 CALL netcdf_close (ng, inlm, outid, outname, .false.)
1119 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1120!
1121 10 FORMAT (/,' STATE_JOIN_NF90 - unable to open grid NetCDF file:', &
1122 & 1x,a)
1123 20 FORMAT ('NLM: STATE_JOIN_NF90 - Concatenating state fields,', &
1124 & t75,a,/,23x,'(Grid ',i2.2,', t = ',a,', InpRec=',i4.4, &
1125 & ', Index=',i1,', InpFile: ',a,')', &
1126 & /,19x,'(Grid ',i2.2,', t = ',a,', OutRec=',i4.4, &
1127 & ', Index=',i1,', OutFile: ',a,')')
1128 30 FORMAT (/,' STATE_JOIN_NF90 - error while reading variable: ', &
1129 & a,2x,'at time record = ',i0, &
1130 & /,19x,'in input NetCDF file: ',a)
1131 40 FORMAT (/,' STATE_JOIN_NF90 - error while writing variable: ', &
1132 & a,2x,'at time record = ',i0, &
1133 & /,19x,'in output NetCDF file: ',a)
1134 50 FORMAT (21x,'- ',a,/,23x,'(Min = ',1p,e15.8, &
1135 & ' Max = ',1p,e15.8,')')
1136!
1137 RETURN
1138 END SUBROUTINE state_join_nf90
1139
1140#if defined PIO_LIB && defined DISTRIBUTE
1141!
1142!***********************************************************************
1143 SUBROUTINE state_join_pio (ng, tile, model, InpName, OutName, &
1144 & InpStrRec, InpEndRec, OutStrRec, &
1145 & LBi, UBi, LBj, UBj)
1146!***********************************************************************
1147!
1148 USE mod_pio_netcdf
1149!
1150! Imported variable declarations.
1151!
1152 integer, intent(in) :: ng, tile, model
1153 integer, intent(in) :: inpstrrec, inpendrec
1154 integer, intent(in) :: lbi, ubi, lbj, ubj
1155
1156 integer, intent(inout) :: outstrrec
1157!
1158 character (len=*) :: inpname, outname
1159!
1160! Local variable declarations.
1161!
1162 integer :: inprec, outrec, tindex
1163 integer :: i, ic, itrc, status
1164
1165 integer :: vsize(4)
1166!
1167 real(r8) :: fmin, fmax
1168 real(dp) :: fscl, stime
1169!
1170 character (len=15) :: tstring
1171 character (len=22) :: t_code
1172
1173 character (len=*), parameter :: myfile = &
1174 & __FILE__//", state_join_pio"
1175!
1176 TYPE (io_desc_t), pointer :: iodesc
1177 TYPE (my_vardesc), pointer :: piovar
1178
1179 TYPE (file_desc_t) :: piofileinp, piofileout
1180!
1181 sourcefile=myfile
1182!
1183!-----------------------------------------------------------------------
1184! Process forward background solution.
1185!-----------------------------------------------------------------------
1186!
1187! Open Input NetCDF file for reading.
1188!
1189 CALL pio_netcdf_open (ng, model, inpname, 0, piofileinp)
1190 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1191 WRITE (stdout,10) trim(inpname)
1192 RETURN
1193 END IF
1194!
1195! Open Output NetCDF file for writing.
1196!
1197 CALL pio_netcdf_open (ng, model, outname, 1, piofileout)
1198 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1199 WRITE (stdout,10) trim(inpname)
1200 RETURN
1201 END IF
1202!
1203! Inquire about the input variables.
1204!
1205 CALL pio_netcdf_inq_var (ng, model, inpname, &
1206 & piofile = piofileinp)
1207 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1208!
1209! Set Vsize to zero to deativate interpolation of input data to model
1210! grid in "nf_fread2d".
1211!
1212 DO i=1,4
1213 vsize(i)=0
1214 END DO
1215!
1216! Scan the variable list and read in needed variables.
1217!
1218 outrec=outstrrec-1
1219 tindex=1
1220 fscl=1.0_dp
1221!
1222 rec_loop : DO inprec=inpstrrec,inpendrec
1223 outrec=outrec+1
1224 var_loop : DO i=1,n_var
1225!
1226! Time.
1227!
1228 check1 : IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
1229 CALL pio_netcdf_get_time (ng, model, inpname, &
1230 & trim(var_name(i)), &
1231 & rclock%DateNumber, stime, &
1232 & piofile = piofileinp, &
1233 & start = (/inprec/), &
1234 & total = (/1/))
1235 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1236 IF (master) THEN
1237 WRITE (stdout,30) trim(vname(1,idtime)), inprec, &
1238 & trim(inpname)
1239 END IF
1240 exit_flag=2
1241 RETURN
1242 END IF
1243!
1244 CALL pio_netcdf_put_fvar (ng, model, outname, &
1245 & trim(var_name(i)), stime, &
1246 & (/outrec/), (/1/), &
1247 & piofile = piofileout, &
1248 & piovar = var_desc(i))
1249 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1250 IF (master) THEN
1251 WRITE (stdout,40) trim(vname(1,idtime)), outrec, &
1252 & trim(outname)
1253 END IF
1254 exit_flag=3
1255 RETURN
1256 ELSE
1257 CALL time_string (stime, t_code)
1258 WRITE (tstring,'(f15.4)') stime*sec2day
1259 IF (master) THEN
1260 WRITE (stdout,20) t_code, &
1261 & ng, trim(adjustl(tstring)), inprec, &
1262 & tindex, trim(inpname), &
1263 & ng, trim(adjustl(tstring)), outrec, &
1264 & tindex, trim(outname)
1265 END IF
1266 END IF
1267!
1268! Free-surface.
1269!
1270 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
1271
1272 piovar%vd=var_desc(i)
1273 IF (kind(ocean(ng)%zeta).eq.8) THEN
1274 piovar%dkind=pio_double
1275 iodesc => iodesc_dp_r2dvar(ng)
1276 ELSE
1277 piovar%dkind=pio_real
1278 iodesc => iodesc_sp_r2dvar(ng)
1279 END IF
1280 piovar%gtype=r2dvar
1281!
1282 status=nf_fread2d(ng, model, inpname, piofileinp, &
1283 & var_name(i), piovar, &
1284 & inprec, iodesc, vsize, &
1285 & lbi, ubi, lbj, ubj, &
1286 & fscl, fmin, fmax, &
1287# ifdef MASKING
1288 & grid(ng) % rmask, &
1289# endif
1290 & ocean(ng) % zeta(:,:,tindex))
1291 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1292 IF (master) THEN
1293 WRITE (stdout,30) trim(vname(1,idfsur)), inprec, &
1294 & trim(inpname)
1295 END IF
1296 exit_flag=2
1297 ioerror=status
1298 RETURN
1299 END IF
1300!
1301 IF (var_type(i).eq.pio_double) THEN
1302 piovar%dkind=pio_double
1303 iodesc => iodesc_dp_r2dvar(ng)
1304 ELSE
1305 piovar%dkind=pio_real
1306 iodesc => iodesc_sp_r2dvar(ng)
1307 END IF
1308!
1309 status=nf_fwrite2d(ng, model, piofileout, idfsur, &
1310 & piovar, outrec, iodesc, &
1311 & lbi, ubi, lbj, ubj, fscl, &
1312# ifdef MASKING
1313 & grid(ng) % rmask, &
1314# endif
1315 & ocean(ng) % zeta(:,:,tindex), &
1316 & setfillval = .false.)
1317 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1318 IF (master) THEN
1319 WRITE (stdout,40) trim(vname(1,idfsur)), outrec, &
1320 & trim(outname)
1321 END IF
1322 exit_flag=3
1323 ioerror=status
1324 RETURN
1325 ELSE
1326 IF (master) THEN
1327 WRITE (stdout,50) trim(vname(2,idfsur)), fmin, fmax
1328 END IF
1329 END IF
1330
1331# ifdef FORWARD_RHS
1332!
1333! Free-surface equation right-hand-side term.
1334!
1335 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrzet))) THEN
1336
1337 piovar%vd=var_desc(i)
1338 IF (kind(ocean(ng)%rzeta).eq.8) THEN
1339 piovar%dkind=pio_double
1340 iodesc => iodesc_dp_r2dvar(ng)
1341 ELSE
1342 piovar%dkind=pio_real
1343 iodesc => iodesc_sp_r2dvar(ng)
1344 END IF
1345 piovar%gtype=r2dvar
1346!
1347 status=nf_fread2d(ng, model, inpname, piofileinp, &
1348 & var_name(i), piovar, &
1349 & inprec, iodesc, vsize, &
1350 & lbi, ubi, lbj, ubj, &
1351 & fscl, fmin, fmax, &
1352# ifdef MASKING
1353 & grid(ng) % rmask, &
1354# endif
1355 & ocean(ng) % rzeta(:,:,tindex))
1356 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1357 IF (master) THEN
1358 WRITE (stdout,30) trim(vname(1,idrzet)), inprec, &
1359 & trim(inpname)
1360 END IF
1361 exit_flag=2
1362 ioerror=status
1363 RETURN
1364 END IF
1365!
1366 IF (var_type(i).eq.pio_double) THEN
1367 piovar%dkind=pio_double
1368 iodesc => iodesc_dp_r2dvar(ng)
1369 ELSE
1370 piovar%dkind=pio_real
1371 iodesc => iodesc_sp_r2dvar(ng)
1372 END IF
1373!
1374 status=nf_fwrite2d(ng, model, piofileout, idrzet, &
1375 & piovar, outrec, iodesc, &
1376 & lbi, ubi, lbj, ubj, fscl, &
1377# ifdef MASKING
1378 & grid(ng) % rmask, &
1379# endif
1380 & ocean(ng) % rzeta(:,:,tindex))
1381 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1382 IF (master) THEN
1383 WRITE (stdout,40) trim(vname(1,idrzet)), outrec, &
1384 & trim(outname)
1385 END IF
1386 exit_flag=3
1387 ioerror=status
1388 RETURN
1389 ELSE
1390 IF (master) THEN
1391 WRITE (stdout,50) trim(vname(2,idrzet)), fmin, fmax
1392 END IF
1393 END IF
1394# endif
1395!
1396! 2D U-momentum component.
1397!
1398 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
1399
1400 piovar%vd=var_desc(i)
1401 IF (kind(ocean(ng)%ubar).eq.8) THEN
1402 piovar%dkind=pio_double
1403 iodesc => iodesc_dp_u2dvar(ng)
1404 ELSE
1405 piovar%dkind=pio_real
1406 iodesc => iodesc_sp_u2dvar(ng)
1407 END IF
1408 piovar%gtype=u2dvar
1409!
1410 status=nf_fread2d(ng, model, inpname, piofileinp, &
1411 & var_name(i), piovar, &
1412 & inprec, iodesc, vsize, &
1413 & lbi, ubi, lbj, ubj, &
1414 & fscl, fmin, fmax, &
1415# ifdef MASKING
1416 & grid(ng) % umask_full, &
1417# endif
1418 & ocean(ng) % ubar(:,:,tindex))
1419 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1420 IF (master) THEN
1421 WRITE (stdout,30) trim(vname(1,idubar)), inprec, &
1422 & trim(inpname)
1423 END IF
1424 exit_flag=2
1425 ioerror=status
1426 RETURN
1427 END IF
1428!
1429 IF (var_type(i).eq.pio_double) THEN
1430 piovar%dkind=pio_double
1431 iodesc => iodesc_dp_u2dvar(ng)
1432 ELSE
1433 piovar%dkind=pio_real
1434 iodesc => iodesc_sp_u2dvar(ng)
1435 END IF
1436!
1437 status=nf_fwrite2d(ng, model, piofileout, idubar, &
1438 & piovar, outrec, iodesc, &
1439 & lbi, ubi, lbj, ubj, fscl, &
1440# ifdef MASKING
1441 & grid(ng) % umask_full, &
1442# endif
1443 & ocean(ng) % ubar(:,:,tindex))
1444 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1445 IF (master) THEN
1446 WRITE (stdout,40) trim(vname(1,idubar)), outrec, &
1447 & trim(outname)
1448 END IF
1449 exit_flag=3
1450 ioerror=status
1451 RETURN
1452 ELSE
1453 IF (master) THEN
1454 WRITE (stdout,50) trim(vname(2,idubar)), fmin, fmax
1455 END IF
1456 END IF
1457
1458# ifdef FORWARD_RHS
1459!
1460! 2D U-equation right-hand-side term.
1461!
1462 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idru2d))) THEN
1463
1464 piovar%vd=var_desc(i)
1465 IF (kind(ocean(ng)%rubar).eq.8) THEN
1466 piovar%dkind=pio_double
1467 iodesc => iodesc_dp_u2dvar(ng)
1468 ELSE
1469 piovar%dkind=pio_real
1470 iodesc => iodesc_sp_u2dvar(ng)
1471 END IF
1472 piovar%gtype=u2dvar
1473!
1474 status=nf_fread2d(ng, model, inpname, piofileinp, &
1475 & var_name(i), piovar, &
1476 & inprec, iodesc, vsize, &
1477 & lbi, ubi, lbj, ubj, &
1478 & fscl, fmin, fmax, &
1479# ifdef MASKING
1480 & grid(ng) % umask_full, &
1481# endif
1482 & ocean(ng) % rubar(:,:,tindex))
1483 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1484 IF (master) THEN
1485 WRITE (stdout,30) trim(vname(1,idru2d)), inprec, &
1486 & trim(inpname)
1487 END IF
1488 exit_flag=2
1489 ioerror=status
1490 RETURN
1491 END IF
1492!
1493 IF (var_type(i).eq.pio_double) THEN
1494 piovar%dkind=pio_double
1495 iodesc => iodesc_dp_u2dvar(ng)
1496 ELSE
1497 piovar%dkind=pio_real
1498 iodesc => iodesc_sp_u2dvar(ng)
1499 END IF
1500!
1501 status=nf_fwrite2d(ng, model, piofileout, idru2d, &
1502 & piovar, outrec, iodesc, &
1503 & lbi, ubi, lbj, ubj, fscl, &
1504# ifdef MASKING
1505 & grid(ng) % umask_full, &
1506# endif
1507 & ocean(ng) % rubar(:,:,tindex))
1508 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1509 IF (master) THEN
1510 WRITE (stdout,40) trim(vname(1,idru2d)), outrec, &
1511 & trim(outname)
1512 END IF
1513 exit_flag=3
1514 ioerror=status
1515 RETURN
1516 ELSE
1517 IF (master) THEN
1518 WRITE (stdout,50) trim(vname(2,idru2d)), fmin, fmax
1519 END IF
1520 END IF
1521# endif
1522!
1523! 2D V-momentum component.
1524!
1525 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
1526
1527 piovar%vd=var_desc(i)
1528 IF (kind(ocean(ng)%vbar).eq.8) THEN
1529 piovar%dkind=pio_double
1530 iodesc => iodesc_dp_v2dvar(ng)
1531 ELSE
1532 piovar%dkind=pio_real
1533 iodesc => iodesc_sp_v2dvar(ng)
1534 END IF
1535 piovar%gtype=v2dvar
1536!
1537 status=nf_fread2d(ng, model, inpname, piofileinp, &
1538 & var_name(i), piovar, &
1539 & inprec, iodesc, vsize, &
1540 & lbi, ubi, lbj, ubj, &
1541 & fscl, fmin, fmax, &
1542# ifdef MASKING
1543 & grid(ng) % vmask_full, &
1544# endif
1545 & ocean(ng) % vbar(:,:,tindex))
1546 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1547 IF (master) THEN
1548 WRITE (stdout,30) trim(vname(1,idvbar)), inprec, &
1549 & trim(inpname)
1550 END IF
1551 exit_flag=2
1552 ioerror=status
1553 RETURN
1554 END IF
1555!
1556 IF (var_type(i).eq.pio_double) THEN
1557 piovar%dkind=pio_double
1558 iodesc => iodesc_dp_v2dvar(ng)
1559 ELSE
1560 piovar%dkind=pio_real
1561 iodesc => iodesc_sp_v2dvar(ng)
1562 END IF
1563!
1564 status=nf_fwrite2d(ng, model, piofileout, idvbar, &
1565 & piovar, outrec, iodesc, &
1566 & lbi, ubi, lbj, ubj, fscl, &
1567# ifdef MASKING
1568 & grid(ng) % vmask_full, &
1569# endif
1570 & ocean(ng) % vbar(:,:,tindex))
1571 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1572 IF (master) THEN
1573 WRITE (stdout,40) trim(vname(1,idvbar)), outrec, &
1574 & trim(outname)
1575 END IF
1576 exit_flag=3
1577 ioerror=status
1578 RETURN
1579 ELSE
1580 IF (master) THEN
1581 WRITE (stdout,50) trim(vname(2,idvbar)), fmin, fmax
1582 END IF
1583 END IF
1584
1585# ifdef FORWARD_RHS
1586!
1587! 2D V-equation right-hand-side term.
1588!
1589 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrv2d))) THEN
1590
1591 piovar%vd=var_desc(i)
1592 IF (kind(ocean(ng)%rvbar).eq.8) THEN
1593 piovar%dkind=pio_double
1594 iodesc => iodesc_dp_v2dvar(ng)
1595 ELSE
1596 piovar%dkind=pio_real
1597 iodesc => iodesc_sp_v2dvar(ng)
1598 END IF
1599 piovar%gtype=v2dvar
1600!
1601 status=nf_fread2d(ng, model, inpname, piofileinp, &
1602 & var_name(i), piovar, &
1603 & inprec, iodesc, vsize, &
1604 & lbi, ubi, lbj, ubj, &
1605 & fscl, fmin, fmax, &
1606# ifdef MASKING
1607 & grid(ng) % vmask_full, &
1608# endif
1609 & ocean(ng) % rvbar(:,:,tindex))
1610 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1611 IF (master) THEN
1612 WRITE (stdout,30) trim(vname(1,idrv2d)), inprec, &
1613 & trim(inpname)
1614 END IF
1615 exit_flag=2
1616 ioerror=status
1617 RETURN
1618 END IF
1619!
1620 IF (var_type(i).eq.pio_double) THEN
1621 piovar%dkind=pio_double
1622 iodesc => iodesc_dp_v2dvar(ng)
1623 ELSE
1624 piovar%dkind=pio_real
1625 iodesc => iodesc_sp_v2dvar(ng)
1626 END IF
1627!
1628 status=nf_fwrite2d(ng, model, piofileout, idrv2d &
1629 & piovar, outrec, iodesc, &
1630 & lbi, ubi, lbj, ubj, fscl, &
1631# ifdef MASKING
1632 & grid(ng) % vmask_full, &
1633# endif
1634 & ocean(ng) % rvbar(:,:,tindex))
1635 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1636 IF (master) THEN
1637 WRITE (stdout,40) trim(vname(1,idrv2d)), outrec, &
1638 & trim(outname)
1639 END IF
1640 exit_flag=3
1641 ioerror=status
1642 RETURN
1643 ELSE
1644 IF (master) THEN
1645 WRITE (stdout,50) trim(vname(2,idrv2d)), fmin, fmax
1646 END IF
1647 END IF
1648# endif
1649# ifdef SOLVE3D
1650# ifdef FORWARD_RHS
1651!
1652! 2D U-equation right-hand-side forcing term.
1653!
1654 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idruct))) THEN
1655
1656 piovar%vd=var_desc(i)
1657 IF (kind(ocean(ng)%rufrc).eq.8) THEN
1658 piovar%dkind=pio_double
1659 iodesc => iodesc_dp_u2dvar(ng)
1660 ELSE
1661 piovar%dkind=pio_real
1662 iodesc => iodesc_sp_u2dvar(ng)
1663 END IF
1664 piovar%gtype=u2dvar
1665!
1666 status=nf_fread2d(ng, model, inpname, piofileinp, &
1667 & var_name(i), piovar, &
1668 & inprec, iodesc, vsize, &
1669 & lbi, ubi, lbj, ubj, &
1670 & fscl, fmin, fmax, &
1671# ifdef MASKING
1672 & grid(ng) % umask_full, &
1673# endif
1674 & ocean(ng) % rufrc)
1675 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1676 IF (master) THEN
1677 WRITE (stdout,30) trim(vname(1,idruct)), inprec, &
1678 & trim(inpname)
1679 END IF
1680 exit_flag=2
1681 ioerror=status
1682 RETURN
1683 END IF
1684!
1685 IF (var_type(i).eq.pio_double) THEN
1686 piovar%dkind=pio_double
1687 iodesc => iodesc_dp_u2dvar(ng)
1688 ELSE
1689 piovar%dkind=pio_real
1690 iodesc => iodesc_sp_u2dvar(ng)
1691 END IF
1692!
1693 status=nf_fwrite2d(ng, model, piofileout, idruct, &
1694 & piovar, outrec, iodesc, &
1695 & lbi, ubi, lbj, ubj, fscl, &
1696# ifdef MASKING
1697 & grid(ng) % umask_full, &
1698# endif
1699 & ocean(ng) % rufrc)
1700 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1701 IF (master) THEN
1702 WRITE (stdout,40) trim(vname(1,idruct)), outrec, &
1703 & trim(outname)
1704 END IF
1705 exit_flag=3
1706 ioerror=status
1707 RETURN
1708 ELSE
1709 IF (master) THEN
1710 WRITE (stdout,50) trim(vname(2,idruct)), fmin, fmax
1711 END IF
1712 END IF
1713# endif
1714!
1715! Time-averaged U-flux component for 2D equations.
1716!
1717 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idufx1))) THEN
1718
1719 piovar%vd=var_desc(i)
1720 IF (kind(coupling(ng)%DU_avg1).eq.8) THEN
1721 piovar%dkind=pio_double
1722 iodesc => iodesc_dp_u2dvar(ng)
1723 ELSE
1724 piovar%dkind=pio_real
1725 iodesc => iodesc_sp_u2dvar(ng)
1726 END IF
1727 piovar%gtype=u2dvar
1728!
1729 status=nf_fread2d(ng, model, inpname, piofileinp, &
1730 & var_name(i), piovar, &
1731 & inprec, iodesc, vsize, &
1732 & lbi, ubi, lbj, ubj, &
1733 & fscl, fmin, fmax, &
1734# ifdef MASKING
1735 & grid(ng) % umask, &
1736# endif
1737 & coupling(ng) % DU_avg1)
1738 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1739 IF (master) THEN
1740 WRITE (stdout,30) trim(vname(1,idufx1)), inprec, &
1741 & trim(inpname)
1742 END IF
1743 exit_flag=2
1744 ioerror=status
1745 RETURN
1746 END IF
1747!
1748 IF (var_type(i).eq.pio_double) THEN
1749 piovar%dkind=pio_double
1750 iodesc => iodesc_dp_u2dvar(ng)
1751 ELSE
1752 piovar%dkind=pio_real
1753 iodesc => iodesc_sp_u2dvar(ng)
1754 END IF
1755!
1756 status=nf_fwrite2d(ng, model, piofileout, idufx1, &
1757 & piovar, outrec, iodesc, &
1758 & lbi, ubi, lbj, ubj, fscl, &
1759# ifdef MASKING
1760 & grid(ng) % umask, &
1761# endif
1762 & coupling(ng) % DU_avg1)
1763 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1764 IF (master) THEN
1765 WRITE (stdout,40) trim(vname(1,idufx1)), outrec, &
1766 & trim(outname)
1767 END IF
1768 exit_flag=3
1769 ioerror=status
1770 RETURN
1771 ELSE
1772 IF (master) THEN
1773 WRITE (stdout,50) trim(vname(2,idufx1)), fmin, fmax
1774 END IF
1775 END IF
1776!
1777! Time-averaged U-flux component for 3D equations.
1778!
1779 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idufx2))) THEN
1780
1781 piovar%vd=var_desc(i)
1782 IF (kind(coupling(ng)%DU_avg2).eq.8) THEN
1783 piovar%dkind=pio_double
1784 iodesc => iodesc_dp_u2dvar(ng)
1785 ELSE
1786 piovar%dkind=pio_real
1787 iodesc => iodesc_sp_u2dvar(ng)
1788 END IF
1789 piovar%gtype=u2dvar
1790!
1791 status=nf_fread2d(ng, model, inpname, piofileinp, &
1792 & var_name(i), piovar, &
1793 & inprec, iodesc, vsize, &
1794 & lbi, ubi, lbj, ubj, &
1795 & fscl, fmin, fmax, &
1796# ifdef MASKING
1797 & grid(ng) % umask, &
1798# endif
1799 & coupling(ng) % DU_avg2)
1800 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1801 IF (master) THEN
1802 WRITE (stdout,30) trim(vname(1,idufx2)), inprec, &
1803 & trim(inpname)
1804 END IF
1805 exit_flag=2
1806 ioerror=status
1807 RETURN
1808 END IF
1809!
1810 IF (var_type(i).eq.pio_double) THEN
1811 piovar%dkind=pio_double
1812 iodesc => iodesc_dp_u2dvar(ng)
1813 ELSE
1814 piovar%dkind=pio_real
1815 iodesc => iodesc_sp_u2dvar(ng)
1816 END IF
1817!
1818 status=nf_fwrite2d(ng, model, piofileout, idufx2, &
1819 & piovar, outrec, iodesc, &
1820 & lbi, ubi, lbj, ubj, fscl, &
1821# ifdef MASKING
1822 & grid(ng) % umask, &
1823# endif
1824 & coupling(ng) % DU_avg2)
1825 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1826 IF (master) THEN
1827 WRITE (stdout,40) trim(vname(1,idufx2)), outrec, &
1828 & trim(outname)
1829 END IF
1830 exit_flag=3
1831 ioerror=status
1832 RETURN
1833 ELSE
1834 IF (master) THEN
1835 WRITE (stdout,50) trim(vname(2,idufx2)), fmin, fmax
1836 END IF
1837 END IF
1838
1839# ifdef FORWARD_RHS
1840!
1841! 2D V-equation right-hand-side forcing term.
1842!
1843 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrvct))) THEN
1844
1845 piovar%vd=var_desc(i)
1846 IF (kind(ocean(ng)%rvfrc).eq.8) THEN
1847 piovar%dkind=pio_double
1848 iodesc => iodesc_dp_v2dvar(ng)
1849 ELSE
1850 piovar%dkind=pio_real
1851 iodesc => iodesc_sp_v2dvar(ng)
1852 END IF
1853 piovar%gtype=v2dvar
1854!
1855 status=nf_fread2d(ng, model, inpname, piofileinp, &
1856 & var_name(i), piovar, &
1857 & inprec, iodesc, vsize, &
1858 & lbi, ubi, lbj, ubj, &
1859 & fscl, fmin, fmax, &
1860# ifdef MASKING
1861 & grid(ng) % vmask_full, &
1862# endif
1863 & ocean(ng) % rvfrc)
1864 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1865 IF (master) THEN
1866 WRITE (stdout,30) trim(vname(1,idrvct)), inprec, &
1867 & trim(inpname)
1868 END IF
1869 exit_flag=2
1870 ioerror=status
1871 RETURN
1872 END IF
1873!
1874 IF (var_type(i).eq.pio_double) THEN
1875 piovar%dkind=pio_double
1876 iodesc => iodesc_dp_v2dvar(ng)
1877 ELSE
1878 piovar%dkind=pio_real
1879 iodesc => iodesc_sp_v2dvar(ng)
1880 END IF
1881!
1882 status=nf_fwrite2d(ng, model, piofileout, idrvct, &
1883 & piovar, outrec, iodesc, &
1884 & lbi, ubi, lbj, ubj, fscl, &
1885# ifdef MASKING
1886 & grid(ng) % vmask_full, &
1887# endif
1888 & ocean(ng) % rvfrc)
1889 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1890 IF (master) THEN
1891 WRITE (stdout,40) trim(vname(1,idrvct)), outrec, &
1892 & trim(outname)
1893 END IF
1894 exit_flag=3
1895 ioerror=status
1896 RETURN
1897 ELSE
1898 IF (master) THEN
1899 WRITE (stdout,50) trim(vname(2,idrvct)), fmin, fmax
1900 END IF
1901 END IF
1902# endif
1903!
1904! Time-averaged V-flux component for 2D equations.
1905!
1906 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvfx1))) THEN
1907
1908 piovar%vd=var_desc(i)
1909 IF (kind(coupling(ng)%DV_avg1).eq.8) THEN
1910 piovar%dkind=pio_double
1911 iodesc => iodesc_dp_v2dvar(ng)
1912 ELSE
1913 piovar%dkind=pio_real
1914 iodesc => iodesc_sp_v2dvar(ng)
1915 END IF
1916 piovar%gtype=v2dvar
1917!
1918 status=nf_fread2d(ng, model, inpname, piofileinp, &
1919 & var_name(i), piovar, &
1920 & inprec, iodesc, vsize, &
1921 & lbi, ubi, lbj, ubj, &
1922 & fscl, fmin, fmax, &
1923# ifdef MASKING
1924 & grid(ng) % vmask_full, &
1925# endif
1926 & coupling(ng) % DV_avg1)
1927 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1928 IF (master) THEN
1929 WRITE (stdout,30) trim(vname(1,idvfx1)), inprec, &
1930 & trim(inpname)
1931 END IF
1932 exit_flag=2
1933 ioerror=status
1934 RETURN
1935 END IF
1936!
1937 IF (var_type(i).eq.pio_double) THEN
1938 piovar%dkind=pio_double
1939 iodesc => iodesc_dp_v2dvar(ng)
1940 ELSE
1941 piovar%dkind=pio_real
1942 iodesc => iodesc_sp_v2dvar(ng)
1943 END IF
1944!
1945 status=nf_fwrite2d(ng, model, piofileout, idvfx1, &
1946 & piovar, outrec, iodesc, &
1947 & lbi, ubi, lbj, ubj, fscl, &
1948# ifdef MASKING
1949 & grid(ng) % vmask_full, &
1950# endif
1951 & coupling(ng) % DV_avg1)
1952 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1953 IF (master) THEN
1954 WRITE (stdout,40) trim(vname(1,idvfx1)), outrec, &
1955 & trim(outname)
1956 END IF
1957 exit_flag=3
1958 ioerror=status
1959 RETURN
1960 ELSE
1961 IF (master) THEN
1962 WRITE (stdout,50) trim(vname(2,idvfx1)), fmin, fmax
1963 END IF
1964 END IF
1965!
1966! Time-averaged U-flux component for 3D equations.
1967!
1968 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvfx2))) THEN
1969
1970 piovar%vd=var_desc(i)
1971 IF (kind(coupling(ng)%DV_avg2).eq.8) THEN
1972 piovar%dkind=pio_double
1973 iodesc => iodesc_dp_v2dvar(ng)
1974 ELSE
1975 piovar%dkind=pio_real
1976 iodesc => iodesc_sp_v2dvar(ng)
1977 END IF
1978 piovar%gtype=v2dvar
1979!
1980 status=nf_fread2d(ng, model, inpname, piofileinp, &
1981 & var_name(i), piovar, &
1982 & inprec, iodesc, vsize, &
1983 & lbi, ubi, lbj, ubj, &
1984 & fscl, fmin, fmax, &
1985# ifdef MASKING
1986 & grid(ng) % vmask_full, &
1987# endif
1988 & coupling(ng) % DV_avg2)
1989 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1990 IF (master) THEN
1991 WRITE (stdout,30) trim(vname(1,idvfx2)), inprec, &
1992 & trim(inpname)
1993 END IF
1994 exit_flag=2
1995 ioerror=status
1996 RETURN
1997 END IF
1998!
1999 IF (var_type(i).eq.pio_double) THEN
2000 piovar%dkind=pio_double
2001 iodesc => iodesc_dp_v2dvar(ng)
2002 ELSE
2003 piovar%dkind=pio_real
2004 iodesc => iodesc_sp_v2dvar(ng)
2005 END IF
2006!
2007 status=nf_fwrite2d(ng, model, piofileout, idvfx2, &
2008 & piovar, outrec, iodesc, &
2009 & lbi, ubi, lbj, ubj, fscl, &
2010# ifdef MASKING
2011 & grid(ng) % vmask_full, &
2012# endif
2013 & coupling(ng) % DV_avg2)
2014 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2015 IF (master) THEN
2016 WRITE (stdout,40) trim(vname(1,idvfx2)), outrec, &
2017 & trim(outname)
2018 END IF
2019 exit_flag=3
2020 ioerror=status
2021 RETURN
2022 ELSE
2023 IF (master) THEN
2024 WRITE (stdout,50) trim(vname(2,idvfx2)), fmin, fmax
2025 END IF
2026 END IF
2027!
2028! 3D U-momentum component.
2029!
2030 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
2031
2032 piovar%vd=var_desc(i)
2033 IF (kind(ocean(ng)%u).eq.8) THEN
2034 piovar%dkind=pio_double
2035 iodesc => iodesc_dp_u3dvar(ng)
2036 ELSE
2037 piovar%dkind=pio_real
2038 iodesc => iodesc_sp_u3dvar(ng)
2039 END IF
2040 piovar%gtype=u3dvar
2041!
2042 status=nf_fread3d(ng, model, inpname, piofileinp, &
2043 & var_name(i), piovar, &
2044 & inprec, iodesc, vsize, &
2045 & lbi, ubi, lbj, ubj, 1, n(ng), &
2046 & fscl, fmin, fmax, &
2047# ifdef MASKING
2048 & grid(ng) % umask_full, &
2049# endif
2050 & ocean(ng) % u(:,:,:,tindex))
2051 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2052 IF (master) THEN
2053 WRITE (stdout,30) trim(vname(1,iduvel)), inprec, &
2054 & trim(inpname)
2055 END IF
2056 exit_flag=2
2057 ioerror=status
2058 RETURN
2059 END IF
2060!
2061 IF (var_type(i).eq.pio_double) THEN
2062 piovar%dkind=pio_double
2063 iodesc => iodesc_dp_u3dvar(ng)
2064 ELSE
2065 piovar%dkind=pio_real
2066 iodesc => iodesc_sp_u3dvar(ng)
2067 END IF
2068!
2069 status=nf_fwrite3d(ng, model, piofileout, iduvel, &
2070 & piovar, outrec, iodesc, &
2071 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
2072# ifdef MASKING
2073 & grid(ng) % umask_full, &
2074# endif
2075 & ocean(ng) % u(:,:,:,tindex))
2076 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2077 IF (master) THEN
2078 WRITE (stdout,40) trim(vname(1,iduvel)), outrec, &
2079 & trim(outname)
2080 END IF
2081 exit_flag=3
2082 ioerror=status
2083 RETURN
2084 ELSE
2085 IF (master) THEN
2086 WRITE (stdout,50) trim(vname(2,iduvel)), fmin, fmax
2087 END IF
2088 END IF
2089
2090# ifdef FORWARD_RHS
2091!
2092! 3D U-momentum right-hand-side term.
2093!
2094 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idru3d))) THEN
2095
2096 piovar%vd=var_desc(i)
2097 IF (kind(ocean(ng)%ru).eq.8) THEN
2098 piovar%dkind=pio_double
2099 iodesc => iodesc_dp_u3dvar(ng)
2100 ELSE
2101 piovar%dkind=pio_real
2102 iodesc => iodesc_sp_u3dvar(ng)
2103 END IF
2104 piovar%gtype=u3dvar
2105!
2106 status=nf_fread3d(ng, model, inpname, piofileinp, &
2107 & var_name(i), piovar, &
2108 & inprec, iodesc, vsize, &
2109 & lbi, ubi, lbj, ubj, 1, n(ng), &
2110 & fscl, fmin, fmax, &
2111# ifdef MASKING
2112 & grid(ng) % umask_full, &
2113# endif
2114 & ocean(ng) % ru(:,:,:,tindex))
2115 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2116 IF (master) THEN
2117 WRITE (stdout,30) trim(vname(1,idru3d)), inprec, &
2118 & trim(inpname)
2119 END IF
2120 exit_flag=2
2121 ioerror=status
2122 RETURN
2123 END IF
2124!
2125 IF (var_type(i).eq.pio_double) THEN
2126 piovar%dkind=pio_double
2127 iodesc => iodesc_dp_u3dvar(ng)
2128 ELSE
2129 piovar%dkind=pio_real
2130 iodesc => iodesc_sp_u3dvar(ng)
2131 END IF
2132!
2133 status=nf_fwrite3d(ng, model, piofileout, idru3d, &
2134 & piovar, outrec, iodesc, &
2135 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
2136# ifdef MASKING
2137 & grid(ng) % umask_full, &
2138# endif
2139 & ocean(ng) % ru(:,:,:,tindex))
2140 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2141 IF (master) THEN
2142 WRITE (stdout,40) trim(vname(1,idru3d)), outrec, &
2143 & trim(outname)
2144 END IF
2145 exit_flag=3
2146 ioerror=status
2147 RETURN
2148 ELSE
2149 IF (master) THEN
2150 WRITE (stdout,50) trim(vname(2,idru3d)), fmin, fmax
2151 END IF
2152 END IF
2153# endif
2154!
2155! 3D V-momentum component.
2156!
2157 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
2158
2159 piovar%vd=var_desc(i)
2160 IF (kind(ocean(ng)%v).eq.8) THEN
2161 piovar%dkind=pio_double
2162 iodesc => iodesc_dp_v3dvar(ng)
2163 ELSE
2164 piovar%dkind=pio_real
2165 iodesc => iodesc_sp_v3dvar(ng)
2166 END IF
2167 piovar%gtype=v3dvar
2168!
2169 status=nf_fread3d(ng, model, inpname, piofileinp, &
2170 & var_name(i), piovar, &
2171 & inprec, iodesc, vsize, &
2172 & lbi, ubi, lbj, ubj, 1, n(ng), &
2173 & fscl, fmin, fmax, &
2174# ifdef MASKING
2175 & grid(ng) % vmask_full, &
2176# endif
2177 & ocean(ng) % v(:,:,:,tindex))
2178 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2179 IF (master) THEN
2180 WRITE (stdout,30) trim(vname(1,idvvel)), inprec, &
2181 & trim(inpname)
2182 END IF
2183 exit_flag=2
2184 ioerror=status
2185 RETURN
2186 END IF
2187!
2188 IF (var_type(i).eq.pio_double) THEN
2189 piovar%dkind=pio_double
2190 iodesc => iodesc_dp_v3dvar(ng)
2191 ELSE
2192 piovar%dkind=pio_real
2193 iodesc => iodesc_sp_v3dvar(ng)
2194 END IF
2195!
2196 status=nf_fwrite3d(ng, model, piofileout, idvvel, &
2197 & piovar, outrec, iodesc, &
2198 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
2199# ifdef MASKING
2200 & grid(ng) % vmask_full, &
2201# endif
2202 & ocean(ng) % v(:,:,:,tindex))
2203 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2204 IF (master) THEN
2205 WRITE (stdout,40) trim(vname(1,idvvel)), outrec, &
2206 & trim(outname)
2207 END IF
2208 exit_flag=3
2209 ioerror=status
2210 RETURN
2211 ELSE
2212 IF (master) THEN
2213 WRITE (stdout,50) trim(vname(2,idvvel)), fmin, fmax
2214 END IF
2215 END IF
2216
2217# ifdef FORWARD_RHS
2218!
2219! 3D V-momentum right-hand-side term.
2220!
2221 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idrv3d))) THEN
2222
2223 piovar%vd=var_desc(i)
2224 IF (kind(ocean(ng)%rv).eq.8) THEN
2225 piovar%dkind=pio_double
2226 iodesc => iodesc_dp_v3dvar(ng)
2227 ELSE
2228 piovar%dkind=pio_real
2229 iodesc => iodesc_sp_v3dvar(ng)
2230 END IF
2231 piovar%gtype=v3dvar
2232!
2233 status=nf_fread3d(ng, model, inpname, piofileinp, &
2234 & var_name(i), piovar, &
2235 & inprec, iodesc, vsize, &
2236 & lbi, ubi, lbj, ubj, 1, n(ng), &
2237 & fscl, fmin, fmax, &
2238# ifdef MASKING
2239 & grid(ng) % vmask_full, &
2240# endif
2241 & ocean(ng) % rv(:,:,:,tindex))
2242 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2243 IF (master) THEN
2244 WRITE (stdout,30) trim(vname(1,idrv3d)), inprec, &
2245 & trim(inpname)
2246 END IF
2247 exit_flag=2
2248 ioerror=status
2249 RETURN
2250 END IF
2251!
2252 IF (var_type(i).eq.pio_double) THEN
2253 piovar%dkind=pio_double
2254 iodesc => iodesc_dp_v3dvar(ng)
2255 ELSE
2256 piovar%dkind=pio_real
2257 iodesc => iodesc_sp_v3dvar(ng)
2258 END IF
2259!
2260 status=nf_fwrite3d(ng, model, piofileout, idrv3d, &
2261 & piovar, outrec, iodesc, &
2262 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
2263# ifdef MASKING
2264 & grid(ng) % vmask_full, &
2265# endif
2266 & ocean(ng) % rv(:,:,:,tindex))
2267 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2268 IF (master) THEN
2269 WRITE (stdout,40) trim(vname(1,idrv3d)), outrec, &
2270 & trim(outname)
2271 END IF
2272 exit_flag=3
2273 ioerror=status
2274 RETURN
2275 ELSE
2276 IF (master) THEN
2277 WRITE (stdout,50) trim(vname(2,idrv3d)), fmin, fmax
2278 END IF
2279 END IF
2280# endif
2281# ifdef FORWARD_MIXING
2282!
2283! Vertical viscosity.
2284!
2285 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvis))) THEN
2286
2287 piovar%vd=var_desc(i)
2288 IF (kind(mixing(ng)%Akv).eq.8) THEN
2289 piovar%dkind=pio_double
2290 iodesc => iodesc_dp_w3dvar(ng)
2291 ELSE
2292 piovar%dkind=pio_real
2293 iodesc => iodesc_sp_w3dvar(ng)
2294 END IF
2295 piovar%gtype=w3dvar
2296!
2297 status=nf_fread3d(ng, model, inpname, piofileinp, &
2298 & var_name(i), piovar, &
2299 & inprec, iodesc, vsize, &
2300 & lbi, ubi, lbj, ubj, 0, n(ng), &
2301 & fscl, fmin, fmax, &
2302# ifdef MASKING
2303 & grid(ng) % rmask, &
2304# endif
2305 & mixing(ng) % Akv)
2306 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2307 IF (master) THEN
2308 WRITE (stdout,30) trim(vname(1,idvvis)), inprec, &
2309 & trim(inpname)
2310 END IF
2311 exit_flag=2
2312 ioerror=status
2313 RETURN
2314 END IF
2315!
2316 IF (var_type(i).eq.pio_double) THEN
2317 piovar%dkind=pio_double
2318 iodesc => iodesc_dp_w3dvar(ng)
2319 ELSE
2320 piovar%dkind=pio_real
2321 iodesc => iodesc_sp_w3dvar(ng)
2322 END IF
2323!
2324 status=nf_fwrite3d(ng, model, piofileout, idvvis, &
2325 & piovar, outrec, iodesc, &
2326 & lbi, ubi, lbj, ubj, 0, n(ng), fscl, &
2327# ifdef MASKING
2328 & grid(ng) % rmask, &
2329# endif
2330 & mixing(ng) % Akv)
2331 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2332 IF (master) THEN
2333 WRITE (stdout,40) trim(vname(1,idvvis)), outrec, &
2334 & trim(outname)
2335 END IF
2336 exit_flag=3
2337 ioerror=status
2338 RETURN
2339 ELSE
2340 IF (master) THEN
2341 WRITE (stdout,50) trim(vname(2,idvvis)), fmin, fmax
2342 END IF
2343 END IF
2344# endif
2345# endif
2346 END IF check1
2347
2348# ifdef SOLVE3D
2349!
2350! Tracer type variables.
2351!
2352 tracer1_loop : DO itrc=1,nt(ng)
2353!
2354 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
2355
2356 piovar%vd=var_desc(i)
2357 IF (kind(ocean(ng)%t).eq.8) THEN
2358 piovar%dkind=pio_double
2359 iodesc => iodesc_dp_r3dvar(ng)
2360 ELSE
2361 piovar%dkind=pio_real
2362 iodesc => iodesc_sp_r3dvar(ng)
2363 END IF
2364 piovar%gtype=r3dvar
2365!
2366 status=nf_fread3d(ng, model, inpname, piofileinp, &
2367 & var_name(i), piovar, &
2368 & inprec, iodesc, vsize, &
2369 & lbi, ubi, lbj, ubj, 1, n(ng), &
2370 & fscl, fmin, fmax, &
2371# ifdef MASKING
2372 & grid(ng) % rmask, &
2373# endif
2374 & ocean(ng) % t(:,:,:,tindex,itrc))
2375 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2376 IF (master) THEN
2377 WRITE (stdout,30) trim(vname(1,idtvar(itrc))), &
2378 & inprec, trim(inpname)
2379 END IF
2380 exit_flag=2
2381 ioerror=status
2382 RETURN
2383 END IF
2384!
2385 IF (var_type(i).eq.pio_double) THEN
2386 piovar%dkind=pio_double
2387 iodesc => iodesc_dp_r3dvar(ng)
2388 ELSE
2389 piovar%dkind=pio_real
2390 iodesc => iodesc_sp_r3dvar(ng)
2391 END IF
2392!
2393 status=nf_fwrite3d(ng, model, piofileout, idtvar(itrc), &
2394 & piovar, outrec, iodesc, &
2395 & lbi, ubi, lbj, ubj, 1, n(ng), fscl, &
2396# ifdef MASKING
2397 & grid(ng) % rmask, &
2398# endif
2399 & ocean(ng) % t(:,:,:,tindex,itrc))
2400 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2401 IF (master) THEN
2402 WRITE (stdout,40) trim(vname(1,idtvar(itrc))), &
2403 & outrec, trim(outname)
2404 END IF
2405 exit_flag=3
2406 ioerror=status
2407 RETURN
2408 ELSE
2409 IF (master) THEN
2410 WRITE (stdout,50) trim(vname(2,idtvar(itrc))), &
2411 & fmin, fmax
2412 END IF
2413 END IF
2414 END IF
2415 END DO tracer1_loop
2416!
2417! Tracer type variables vertical diffusion.
2418!
2419 tracer2_loop : DO itrc=1,nat
2420 piovar%vd=var_desc(i)
2421 piovar%dkind=pio_type
2422 piovar%gtype=w3dvar
2423!
2424 IF (trim(var_name(i)).eq.trim(vname(1,iddiff(itrc)))) THEN
2425
2426 piovar%vd=var_desc(i)
2427 IF (kind(mixing(ng)%Akt).eq.8) THEN
2428 piovar%dkind=pio_double
2429 iodesc => iodesc_dp_w3dvar(ng)
2430 ELSE
2431 piovar%dkind=pio_real
2432 iodesc => iodesc_sp_w3dvar(ng)
2433 END IF
2434 piovar%gtype=w3dvar
2435!
2436 status=nf_fread3d(ng, model, inpname, piofileinp, &
2437 & var_name(i), piovar, &
2438 & inprec, iodesc, vsize, &
2439 & lbi, ubi, lbj, ubj, 0, n(ng), &
2440 & fscl, fmin, fmax, &
2441# ifdef MASKING
2442 & grid(ng) % rmask, &
2443# endif
2444 & mixing(ng) % Akt(:,:,:,itrc))
2445 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2446 IF (master) THEN
2447 WRITE (stdout,30) trim(vname(1,iddiff(itrc))), &
2448 & inprec, trim(inpname)
2449 END IF
2450 exit_flag=2
2451 ioerror=status
2452 RETURN
2453 END IF
2454!
2455 IF (var_type(i).eq.pio_double) THEN
2456 piovar%dkind=pio_double
2457 iodesc => iodesc_dp_w3dvar(ng)
2458 ELSE
2459 piovar%dkind=pio_real
2460 iodesc => iodesc_sp_w3dvar(ng)
2461 END IF
2462!
2463 status=nf_fwrite3d(ng, model, piofileout, iddiff(itrc), &
2464 & piovar, outrec, iodesc, &
2465 & lbi, ubi, lbj, ubj, 0, n(ng), fscl, &
2466# ifdef MASKING
2467 & grid(ng) % rmask, &
2468# endif
2469 & mixing(ng) % Akt(:,:,:,itrc))
2470 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
2471 IF (master) THEN
2472 WRITE (stdout,40) trim(vname(1,iddiff(itrc))), &
2473 & outrec, trim(outname)
2474 END IF
2475 exit_flag=3
2476 ioerror=status
2477 RETURN
2478 ELSE
2479 IF (master) THEN
2480 WRITE (stdout,50) trim(vname(2,iddiff(itrc))), &
2481 & fmin, fmax
2482 END IF
2483 END IF
2484 END IF
2485 END DO tracer2_loop
2486# endif
2487
2488 END DO var_loop
2489
2490 END DO rec_loop
2491!
2492! Update output record value.
2493!
2494 outstrrec=outrec
2495!
2496! Close input and output files for compleate synchronization.
2497!
2498 CALL pio_netcdf_close (ng, inlm, piofileinp, inpname, .false.)
2499 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2500!
2501 CALL pio_netcdf_close (ng, inlm, piofileout, outname, .false.)
2502 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2503!
2504 10 FORMAT (/,' STATE_JOIN_PIO - unable to open grid NetCDF file:', &
2505 & 1x,a)
2506 20 FORMAT ('NLM: STATE_JOIN_PIO - Concatenating state fields,', &
2507 & t75,a,/,23x,'(Grid ',i2.2,', t = ',a,', InpRec=',i4.4, &
2508 & ', Index=',i1,', InpFile: ',a,')', &
2509 & /,19x,'(Grid ',i2.2,', t = ',a,', OutRec=',i4.4, &
2510 & ', Index=',i1,', OutFile: ',a,')')
2511 30 FORMAT (/,' STATE_JOIN_PIO - error while reading variable: ', &
2512 & a,2x,'at time record = ',i0, &
2513 & /,18x,'in input NetCDF file: ',a)
2514 40 FORMAT (/,' STATE_JOIN_PIO - error while writing variable: ', &
2515 & a,2x,'at time record = ',i0, &
2516 & /,18x,'in output NetCDF file: ',a)
2517 50 FORMAT (21x,'- ',a,/,23x,'(Min = ',1p,e15.8, &
2518 & ' Max = ',1p,e15.8,')')
2519!
2520 RETURN
2521 END SUBROUTINE state_join_pio
2522#endif
2523
2524 END MODULE state_join_mod
subroutine, public time_string(mytime, date_string)
Definition dateclock.F:1272
type(t_coupling), dimension(:), allocatable coupling
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer ioerror
integer stdout
character(len=256) sourcefile
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
integer, dimension(2) iddiff
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer idrv3d
integer idubar
integer idvvel
integer, parameter io_pio
Definition mod_ncparam.F:96
integer idvfx2
integer idru2d
integer idfsur
integer, dimension(:), allocatable idtvar
integer idvfx1
integer idufx2
integer iduvel
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer idru3d
integer idvvis
integer idrzet
integer idrvct
integer idufx1
integer idruct
integer idrv2d
integer idvbar
subroutine, public netcdf_close(ng, model, ncid, ncname, lupdate)
subroutine, public netcdf_open(ng, model, ncname, omode, 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
subroutine, public netcdf_inq_var(ng, model, ncname, ncid, myvarname, searchvar, varid, nvardim, nvaratt)
integer, dimension(mvars) var_type
Definition mod_netcdf.F:163
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
logical master
integer nat
Definition mod_param.F:499
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, parameter w3dvar
Definition mod_param.F:724
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_sp_w3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dvar
integer, parameter pio_type
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dvar
type(var_desc_t), dimension(:), pointer var_desc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
subroutine, public pio_netcdf_inq_var(ng, model, ncname, piofile, myvarname, searchvar, piovar, nvardim, nvaratt)
subroutine, public pio_netcdf_close(ng, model, piofile, ncname, lupdate)
subroutine, public pio_netcdf_open(ng, model, ncname, omode, piofile)
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_w3dvar
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
type(t_clock) rclock
real(dp), parameter sec2day
integer exit_flag
integer noerror
subroutine, private state_join_pio(ng, tile, model, inpname, outname, inpstrrec, inpendrec, outstrrec, lbi, ubi, lbj, ubj)
subroutine, public state_join(ng, tile, model, inpname, outname, iotype, inpstrrec, inpendrec, outstrrec)
Definition state_join.F:62
subroutine, private state_join_nf90(ng, tile, model, inpname, outname, inpstrrec, inpendrec, outstrrec, lbi, ubi, lbj, ubj)
Definition state_join.F:117
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52