ROMS
Loading...
Searching...
No Matches
get_3dfld.F
Go to the documentation of this file.
1#include "cppdefs.h"
2#ifdef SOLVE3D
3 SUBROUTINE get_3dfld (ng, model, ifield, ncid, &
4# if defined PIO_LIB && defined DISTRIBUTE
5 & pioFile, &
6# endif
7 nfiles, S, update, &
8 & LBi, UBi, LBj, UBj, LBk, UBk, &
9 & Iout, Irec, &
10# ifdef MASKING
11 & Fmask, &
12# endif
13 & Fout)
14!
15!git $Id$
16!================================================== Hernan G. Arango ===
17! Copyright (c) 2002-2025 The ROMS Group !
18! Licensed under a MIT/X style license !
19! See License_ROMS.md !
20!=======================================================================
21! !
22! This routine reads in requested 3D field (point or grided) from !
23! specified NetCDF file. Forward time processing. !
24! !
25! On Input: !
26! !
27! ng Nested grid number. !
28! model Calling model identifier. !
29! ifield Field ID. !
30! ncid NetCDF file ID. !
31# if defined PIO_LIB && defined DISTRIBUTE
32! pioFile PIO file descriptor structure, TYPE(file_desc_t) !
33! pioFile%fh file handler !
34! pioFile%iosystem IO system descriptor (struct) !
35# endif
36! nfiles Number of input NetCDF files. !
37! S I/O derived type structure, TYPE(T_IO). !
38! LBi I-dimension Lower bound. !
39! UBi I-dimension Upper bound. !
40! LBj J-dimension Lower bound. !
41! UBj J-dimension Upper bound. !
42! LBk K-dimension Lower bound. !
43! UBk K-dimension Upper bound. !
44! Iout Size of the outer dimension, if any. Otherwise, !
45! Iout must be set to one by the calling program. !
46! Irec Number of 3D field records to read. !
47! Fmask Land/Sea mask, if any. !
48! !
49! On Output: !
50! !
51! Fout Read field. !
52! update Switch indicating reading of the requested field !
53! the current time step. !
54! !
55!=======================================================================
56!
57 USE mod_param
58 USE mod_parallel
59 USE mod_iounits
60 USE mod_ncparam
61 USE mod_scalars
62!
63 USE strings_mod, ONLY : founderror
64!
65 implicit none
66!
67! Imported variable declarations.
68!
69 logical, intent(out) :: update
70!
71 integer, intent(in) :: ng, model, ifield, nfiles
72 integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk
73 integer, intent(in) :: Iout, Irec
74 integer, intent(inout) :: ncid
75!
76# if defined PIO_LIB && defined DISTRIBUTE
77 TYPE (File_desc_t), intent(inout) :: pioFile
78# endif
79 TYPE(t_io), intent(inout) :: S(nfiles)
80!
81# ifdef MASKING
82 real(r8), intent(in) :: Fmask(LBi:UBi,LBj:UBj)
83# endif
84 real(r8), intent(inout) :: Fout(LBi:UBi,LBj:UBj,LBk:UBk,Iout)
85!
86! Local variable declarations.
87!
88 character (len=*), parameter :: MyFile = &
89 & __FILE__
90!
91!-----------------------------------------------------------------------
92! Read in requested 3D field according to IO type.
93!-----------------------------------------------------------------------
94!
95 SELECT CASE (s(1)%IOtype)
96 CASE (io_nf90)
97 CALL get_3dfld_nf90 (ng, model, ifield, ncid, &
98 & nfiles, s, update, &
99 & lbi, ubi, lbj, ubj, lbk, ubk, &
100 & iout, irec, &
101# ifdef MASKING
102 & fmask, &
103# endif
104 & fout)
105
106# if defined PIO_LIB && defined DISTRIBUTE
107 CASE (io_pio)
108 CALL get_3dfld_pio (ng, model, ifield, piofile, &
109 & nfiles, s, update, &
110 & lbi, ubi, lbj, ubj, lbk, ubk, &
111 & iout, irec, &
112# ifdef MASKING
113 & fmask, &
114# endif
115 & fout)
116
117# endif
118 CASE DEFAULT
119 IF (master) WRITE (stdout,10) s(1)%IOtype
120 exit_flag=2
121 END SELECT
122 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
123!
124 10 FORMAT (' GET_3DFLD - Illegal input file type, io_type = ',i0, &
125 & /,13x,'Check KeyWord ''INP_LIB'' in ''roms.in''.')
126!
127 RETURN
128 END SUBROUTINE get_3dfld
129!
130!***********************************************************************
131 SUBROUTINE get_3dfld_nf90 (ng, model, ifield, ncid, &
132 & nfiles, S, update, &
133 & LBi, UBi, LBj, UBj, LBk, UBk, &
134 & Iout, Irec, &
135# ifdef MASKING
136 & Fmask, &
137# endif
138 & Fout)
139!***********************************************************************
140!
141 USE mod_param
142 USE mod_parallel
143 USE mod_iounits
144 USE mod_ncparam
145 USE mod_netcdf
146 USE mod_scalars
147!
148 USE dateclock_mod, ONLY : time_string
149 USE inquiry_mod, ONLY : inquiry
150 USE nf_fread3d_mod, ONLY : nf_fread3d
151 USE strings_mod, ONLY : founderror
152!
153 implicit none
154!
155! Imported variable declarations.
156!
157 logical, intent(out) :: update
158!
159 integer, intent(in) :: ng, model, ifield, nfiles
160 integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk
161 integer, intent(in) :: Iout, Irec
162 integer, intent(inout) :: ncid
163!
164 TYPE(t_io), intent(inout) :: S(nfiles)
165!
166# ifdef MASKING
167 real(r8), intent(in) :: Fmask(LBi:UBi,LBj:UBj)
168# endif
169 real(r8), intent(inout) :: Fout(LBi:UBi,LBj:UBj,LBk:UBk,Iout)
170!
171! Local variable declarations.
172!
173 logical :: Lgridded, Linquire, Liocycle, Lmulti, Lonerec
174!
175 integer :: Nrec, Tid, Tindex, Trec, Vid, Vtype
176 integer :: i, job, lend, lstr, lvar, status
177 integer :: Vsize(4)
178# ifdef CHECKSUM
179 integer(i8b) :: Fhash
180# endif
181!
182 real(r8) :: Fmax, Fmin, Fval, my_Fmin, my_Fmax
183
184 real(dp) :: Clength, Tdelta, Tend
185 real(dp) :: Tmax, Tmin, Tmono, Tscale, Tstr
186 real(dp) :: Tsec, Tval
187!
188 character (len=22) :: t_code
189
190 character (len=*), parameter :: MyFile = &
191 & __FILE__//", get_3dfld_nf90"
192!
193 sourcefile=myfile
194!
195!-----------------------------------------------------------------------
196! Initialize.
197!-----------------------------------------------------------------------
198!
199 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
200!
201! Determine if inquiring about the field to process in input NetCDF
202! file(s). This usually happens on first call or when the field
203! time records are split (saved) in several multi-files.
204!
205 linquire=.false.
206 lmulti=.false.
207 IF (iic(ng).eq.0) linquire=.true.
208 IF (.not.linquire.and. &
209 & ((iinfo(10,ifield,ng).gt.1).and. &
210 & (linfo( 6,ifield,ng).or. &
211 & (finfo( 2,ifield,ng)*day2sec.lt.time(ng))))) THEN
212 linquire=.true.
213 lmulti=.true.
214 END IF
215!
216! If appropriate, inquire about the contents of input NetCDF file and
217! fill information arrays.
218!
219! Also, if appropriate, deactivate the Linfo(6,ifield,ng) switch after
220! the query for the UPPER snapshot interpolant from previous multifile
221! in the list. The switch was activated previously to indicate the
222! processing of the FIRST record of the file for the LOWER snapshot
223! interpolant.
224!
225 IF (linquire) THEN
226 job=3
227 CALL inquiry (ng, model, job, iout, irec, 1, ifield, ncid, &
228 & lmulti, nfiles, s)
229 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
230 IF (linfo(6,ifield,ng)) THEN
231 linfo(6,ifield,ng)=.false.
232 END IF
233 END IF
234!
235!-----------------------------------------------------------------------
236! If appropriate, read in new data.
237!-----------------------------------------------------------------------
238!
239# ifdef CHECKSUM
240 fhash=0_i8b
241# endif
242 update=.false.
243 tmono=finfo(7,ifield,ng)
244!
245 IF ((tmono.lt.time(ng)).or.(iic(ng).eq.0).or. &
246 & (iic(ng).eq.ntstart(ng))) THEN
247!
248! Load properties for requested field from information arrays.
249!
250 lgridded=linfo(1,ifield,ng)
251 liocycle=linfo(2,ifield,ng)
252 lonerec =linfo(3,ifield,ng)
253 vtype =iinfo(1,ifield,ng)
254 vid =iinfo(2,ifield,ng)
255 tid =iinfo(3,ifield,ng)
256 nrec =iinfo(4,ifield,ng)
257 vsize(1)=iinfo(5,ifield,ng)
258 vsize(2)=iinfo(6,ifield,ng)
259 vsize(3)=iinfo(7,ifield,ng)
260 tindex =iinfo(8,ifield,ng)
261 trec =iinfo(9,ifield,ng)
262 tmin =finfo(1,ifield,ng)
263 tmax =finfo(2,ifield,ng)
264 clength =finfo(5,ifield,ng)
265 tscale =finfo(6,ifield,ng)
266 ncfile =cinfo(ifield,ng)
267!
268 IF (liocycle) THEN
269 trec=mod(trec,nrec)+1
270 ELSE
271 trec=trec+1
272 END IF
273 iinfo(9,ifield,ng)=trec
274!
275 IF (trec.le.nrec) THEN
276!
277! Set rolling index for two-time record storage of input data. If
278! "Iout" is unity, input data is stored in timeless array by the
279! calling program. If Irec > 1, this routine is used to read a 3D
280! field varying in another non-time dimension.
281!
282 IF (irec.eq.1) THEN
283 IF (iout.eq.1) THEN
284 tindex=1
285 ELSE
286 tindex=3-tindex
287 END IF
288 iinfo(8,ifield,ng)=tindex
289 END IF
290!
291! Read in time coordinate and scale it to day units.
292!
293 IF (tid.ge.0) THEN
294 CALL netcdf_get_time (ng, model, ncfile, tname(ifield), &
295 & rclock%DateNumber, tval, &
296 & ncid = ncid, &
297 & start = (/trec/), &
298 & total = (/1/))
299 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
300 IF (master) WRITE (stdout,40) trim(tname(ifield)), trec
301 RETURN
302 END IF
303 END IF
304 tval=tval*tscale
305 vtime(tindex,ifield,ng)=tval
306!
307! Activate switch Linfo(6,ifield,ng) if processing the LAST record of
308! the file for the LOWER time snapshot. We need to get the UPPER time
309! snapshot from NEXT multifile.
310!
311 IF ((trec.eq.nrec).and.(tval*day2sec.le.time(ng))) THEN
312 linfo(6,ifield,ng)=.true.
313 END IF
314!
315! Read in 3D-grided or point data.
316!
317 IF (vid.ge.0) THEN
318 IF (lgridded) THEN
319 fmin=spval
320 fmax=-spval
321 IF (irec.gt.1) THEN
322 DO i=1,irec
323 status=nf_fread3d(ng, model, ncfile, ncid, &
324 & vname(1,ifield), vid, &
325 & i, vtype, vsize, &
326 & lbi, ubi, lbj, ubj, lbk, ubk, &
327 & fscale(ifield,ng), &
328 & my_fmin, my_fmax, &
329# ifdef MASKING
330 & fmask, &
331# endif
332# ifdef CHECKSUM
333 & fout(:,:,:,i), &
334 & checksum = fhash)
335# else
336 & fout(:,:,:,i))
337# endif
338 fmin=min(fmin,my_fmin)
339 fmax=max(fmax,my_fmax)
340 END DO
341 finfo(8,ifield,ng)=fmin
342 finfo(9,ifield,ng)=fmax
343 ELSE
344 status=nf_fread3d(ng, model, ncfile, ncid, &
345 & vname(1,ifield), vid, &
346 & trec, vtype, vsize, &
347 & lbi, ubi, lbj, ubj, lbk, ubk, &
348 & fscale(ifield,ng), fmin, fmax, &
349# ifdef MASKING
350 & fmask, &
351# endif
352# ifdef CHECKSUM
353 & fout(:,:,:,tindex), &
354 & checksum = fhash)
355# else
356 & fout(:,:,:,tindex))
357# endif
358 finfo(8,ifield,ng)=fmin
359 finfo(9,ifield,ng)=fmax
360 END IF
361 ELSE
362 CALL netcdf_get_fvar (ng, model, ncfile, &
363 & vname(1,ifield), fval, &
364 & ncid = ncid, &
365 & start = (/trec/), &
366 & total = (/1/))
367 fval=fval*fscale(ifield,ng)
368 fpoint(tindex,ifield,ng)=fval
369 fmin=fval
370 fmax=fval
371 END IF
372 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
373 IF (master) THEN
374 WRITE (stdout,40) trim(vname(1,ifield)), trec
375 END IF
376 RETURN
377 END IF
378 IF (master) THEN
379 IF (irec.gt.1) THEN
380 WRITE (stdout,50) trim(vname(2,ifield)), ng, fmin, fmax
381 ELSE
382 lstr=scan(ncfile,'/',back=.true.)+1
383 lend=len_trim(ncfile)
384 lvar=min(46,len_trim(vname(2,ifield)))
385 tsec=tval*day2sec
386 CALL time_string (tsec, t_code)
387 WRITE (stdout,60) vname(2,ifield)(1:lvar), t_code, &
388 & ng, trec, tindex, ncfile(lstr:lend), &
389 & tmin, tmax, tval, fmin, fmax
390 END IF
391# ifdef CHECKSUM
392 WRITE (stdout,70) fhash
393# endif
394 END IF
395 update=.true.
396 END IF
397 END IF
398!
399! Increment the local time variable "Tmono" by the interval between
400! snapshots. If the interval is negative, indicating cycling, add in
401! a cycle length. Load time value (sec) into "Tintrp" which used
402! during interpolation between snapshots.
403!
404 IF (.not.lonerec) THEN
405 tdelta=vtime(tindex,ifield,ng)-vtime(3-tindex,ifield,ng)
406 IF (liocycle.and.(tdelta.lt.0.0_r8)) THEN
407 tdelta=tdelta+clength
408 END IF
409 tmono=tmono+tdelta*day2sec
410 finfo(7,ifield,ng)=tmono
411 tintrp(tindex,ifield,ng)=tmono
412 END IF
413 END IF
414!
415 10 FORMAT (/,' GET_3DFLD_NF90 - unable to find dimension ',a, &
416 & /,18x,'for variable: ',a,/,18x,'in file: ',a, &
417 & /,18x,'file is not CF compliant...')
418 20 FORMAT (/,' GET_3DFLD_NF90 - unable to find requested variable:', &
419 & 1x,a,/,18x,'in input NetCDF file: ',a)
420 30 FORMAT (/,' GET_3DFLD_NF90 - unable to open input NetCDF', &
421 & ' file: ',a)
422 40 FORMAT (/,' GET_3DFLD_NF90 - error while reading variable: ',a, &
423 & 2x,' at TIME index = ',i0)
424 50 FORMAT (2x,'GET_3DFLD_NF90 - ',a,/,22x,'(Grid = ',i2.2, &
425 & ', Min = ',1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
426 60 FORMAT (2x,'GET_3DFLD_NF90 - ',a,',',t75,a,/,22x, &
427 & '(Grid=',i2.2,', Rec=',i0,', Index=',i1, &
428 & ', File: ',a,')',/,22x, &
429 & '(Tmin= ', f15.4, ' Tmax= ', f15.4,')', &
430 & t71, 't = ', f15.4 ,/,22x, &
431 & '(Min = ', 1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
432# ifdef CHECKSUM
433 70 FORMAT (22x,'(CheckSum = ',i0,')')
434# endif
435!
436 RETURN
437 END SUBROUTINE get_3dfld_nf90
438
439# if defined PIO_LIB && defined DISTRIBUTE
440!
441!***********************************************************************
442 SUBROUTINE get_3dfld_pio (ng, model, ifield, pioFile, &
443 & nfiles, S, update, &
444 & LBi, UBi, LBj, UBj, LBk, UBk, &
445 & Iout, Irec, &
446# ifdef MASKING
447 & Fmask, &
448# endif
449 & Fout)
450!***********************************************************************
451!
452 USE mod_param
453 USE mod_parallel
454 USE mod_iounits
455 USE mod_ncparam
457 USE mod_scalars
458!
459 USE dateclock_mod, ONLY : time_string
460 USE inquiry_mod, ONLY : inquiry
461 USE nf_fread3d_mod, ONLY : nf_fread3d
462 USE strings_mod, ONLY : founderror
463!
464 implicit none
465!
466! Imported variable declarations.
467!
468 logical, intent(out) :: update
469!
470 integer, intent(in) :: ng, model, ifield, nfiles
471 integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk
472 integer, intent(in) :: Iout, Irec
473!
474 TYPE (File_desc_t), intent(inout) :: pioFile
475 TYPE(t_io), intent(inout) :: S(nfiles)
476!
477# ifdef MASKING
478 real(r8), intent(in) :: Fmask(LBi:UBi,LBj:UBj)
479# endif
480 real(r8), intent(inout) :: Fout(LBi:UBi,LBj:UBj,LBk:UBk,Iout)
481!
482! Local variable declarations.
483!
484 logical :: Lgridded, Linquire, Liocycle, Lmulti, Lonerec
485!
486 integer :: Nrec, Tindex, Trec, Vtype
487 integer :: i, job, lend, lstr, lvar, status
488 integer :: Vsize(4)
489# ifdef CHECKSUM
490 integer(i8b) :: Fhash
491# endif
492!
493 real(r8) :: Fmax, Fmin, Fval, my_Fmin, my_Fmax
494
495 real(dp) :: Clength, Tdelta, Tend
496 real(dp) :: Tmax, Tmin, Tmono, Tscale, Tstr
497 real(dp) :: Tsec, Tval
498!
499 character (len=22) :: t_code
500
501 character (len=*), parameter :: MyFile = &
502 & __FILE__//", get_3dfld_pio"
503!
504 TYPE (IO_Desc_t), pointer :: ioDesc
505 TYPE (My_VarDesc) :: TpioVar, VpioVar
506!
507 sourcefile=myfile
508!
509!-----------------------------------------------------------------------
510! Initialize.
511!-----------------------------------------------------------------------
512!
513 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
514!
515! Determine if inquiring about the field to process in input NetCDF
516! file(s). This usually happens on first call or when the field
517! time records are split (saved) in several multi-files.
518!
519 linquire=.false.
520 lmulti=.false.
521 IF (iic(ng).eq.0) linquire=.true.
522 IF (.not.linquire.and. &
523 & ((iinfo(10,ifield,ng).gt.1).and. &
524 & (linfo( 6,ifield,ng).or. &
525 & (finfo( 2,ifield,ng)*day2sec.lt.time(ng))))) THEN
526 linquire=.true.
527 lmulti=.true.
528 END IF
529!
530! If appropriate, inquire about the contents of input NetCDF file and
531! fill information arrays.
532!
533! Also, if appropriate, deactivate the Linfo(6,ifield,ng) switch after
534! the query for the UPPER snapshot interpolant from previous multifile
535! in the list. The switch was activated previously to indicate the
536! processing of the FIRST record of the file for the LOWER snapshot
537! interpolant.
538!
539 IF (linquire) THEN
540 job=3
541 CALL inquiry (ng, model, job, iout, irec, 1, ifield, piofile, &
542 & lmulti, nfiles, s)
543 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
544 IF (linfo(6,ifield,ng)) THEN
545 linfo(6,ifield,ng)=.false.
546 END IF
547 END IF
548!
549!-----------------------------------------------------------------------
550! If appropriate, read in new data.
551!-----------------------------------------------------------------------
552!
553# ifdef CHECKSUM
554 fhash=0_i8b
555# endif
556 update=.false.
557 tmono=finfo(7,ifield,ng)
558!
559 IF ((tmono.lt.time(ng)).or.(iic(ng).eq.0).or. &
560 & (iic(ng).eq.ntstart(ng))) THEN
561!
562! Load properties for requested field from information arrays.
563!
564 lgridded=linfo(1,ifield,ng)
565 liocycle=linfo(2,ifield,ng)
566 lonerec =linfo(3,ifield,ng)
567 vtype =iinfo(1,ifield,ng)
568 vpiovar =dinfo(1,ifield,ng)
569 tpiovar =dinfo(2,ifield,ng)
570 nrec =iinfo(4,ifield,ng)
571 vsize(1)=iinfo(5,ifield,ng)
572 vsize(2)=iinfo(6,ifield,ng)
573 vsize(3)=iinfo(7,ifield,ng)
574 tindex =iinfo(8,ifield,ng)
575 trec =iinfo(9,ifield,ng)
576 tmin =finfo(1,ifield,ng)
577 tmax =finfo(2,ifield,ng)
578 clength =finfo(5,ifield,ng)
579 tscale =finfo(6,ifield,ng)
580 ncfile =cinfo(ifield,ng)
581!
582 IF (liocycle) THEN
583 trec=mod(trec,nrec)+1
584 ELSE
585 trec=trec+1
586 END IF
587 iinfo(9,ifield,ng)=trec
588!
589 IF (trec.le.nrec) THEN
590!
591! Set rolling index for two-time record storage of input data. If
592! "Iout" is unity, input data is stored in timeless array by the
593! calling program. If Irec > 1, this routine is used to read a 3D
594! field varying in another non-time dimension.
595!
596 IF (irec.eq.1) THEN
597 IF (iout.eq.1) THEN
598 tindex=1
599 ELSE
600 tindex=3-tindex
601 END IF
602 iinfo(8,ifield,ng)=tindex
603 END IF
604!
605! Read in time coordinate and scale it to day units.
606!
607 IF (tpiovar%vd%varID.ge.0) THEN
608 CALL pio_netcdf_get_time (ng, model, ncfile, tname(ifield), &
609 & rclock%DateNumber, tval, &
610 & piofile = piofile, &
611 & start = (/trec/), &
612 & total = (/1/))
613 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
614 IF (master) WRITE (stdout,40) trim(tname(ifield)), trec
615 RETURN
616 END IF
617 END IF
618 tval=tval*tscale
619 vtime(tindex,ifield,ng)=tval
620!
621! Activate switch Linfo(6,ifield,ng) if processing the LAST record of
622! the file for the LOWER time snapshot. We need to get the UPPER time
623! snapshot from NEXT multifile.
624!
625 IF ((trec.eq.nrec).and.(tval*day2sec.le.time(ng))) THEN
626 linfo(6,ifield,ng)=.true.
627 END IF
628!
629! Read in 3D-grided or point data.
630!
631 IF (vpiovar%vd%varID.ge.0) THEN
632 SELECT CASE (vtype)
633# ifdef SEDIMENT
634 CASE (b3dvar)
635 IF (kind(fout).eq.8) THEN
636 iodesc => iodesc_dp_b3dvar(ng)
637 ELSE
638 iodesc => iodesc_sp_b3dvar(ng)
639 END IF
640# endif
641# if defined DIAGNOSTICS_BIO && defined ECOSIM
642 CASE (l3dvar)
643 IF (kind(fout).eq.8) THEN
644 iodesc => iodesc_dp_l3dvar(ng)
645 ELSE
646 iodesc => iodesc_sp_l3dvar(ng)
647 END IF
648# endif
649 CASE (p3dvar)
650 IF (kind(fout).eq.8) THEN
651 iodesc => iodesc_dp_p3dvar(ng)
652 ELSE
653 iodesc => iodesc_sp_p3dvar(ng)
654 END IF
655 CASE (r3dvar)
656 IF (kind(fout).eq.8) THEN
657 iodesc => iodesc_dp_r3dvar(ng)
658 ELSE
659 iodesc => iodesc_sp_r3dvar(ng)
660 END IF
661 CASE (u3dvar)
662 IF (kind(fout).eq.8) THEN
663 iodesc => iodesc_dp_u3dvar(ng)
664 ELSE
665 iodesc => iodesc_sp_u3dvar(ng)
666 END IF
667 CASE (v3dvar)
668 IF (kind(fout).eq.8) THEN
669 iodesc => iodesc_dp_v3dvar(ng)
670 ELSE
671 iodesc => iodesc_sp_v3dvar(ng)
672 END IF
673 CASE (w3dvar)
674 IF (kind(fout).eq.8) THEN
675 iodesc => iodesc_dp_w3dvar(ng)
676 ELSE
677 iodesc => iodesc_sp_w3dvar(ng)
678 END IF
679 END SELECT
680!
681 IF (lgridded) THEN
682 fmin=spval
683 fmax=-spval
684 IF (irec.gt.1) THEN
685 DO i=1,irec
686 status=nf_fread3d(ng, model, ncfile, piofile, &
687 & vname(1,ifield), vpiovar, &
688 & i, iodesc, vsize, &
689 & lbi, ubi, lbj, ubj, lbk, ubk, &
690 & fscale(ifield,ng), &
691 & my_fmin, my_fmax, &
692# ifdef MASKING
693 & fmask, &
694# endif
695# ifdef CHECKSUM
696 & fout(:,:,:,i), &
697 & checksum = fhash)
698# else
699 & fout(:,:,:,i))
700# endif
701 fmin=min(fmin,my_fmin)
702 fmax=max(fmax,my_fmax)
703 END DO
704 finfo(8,ifield,ng)=fmin
705 finfo(9,ifield,ng)=fmax
706 ELSE
707 status=nf_fread3d(ng, model, ncfile, piofile, &
708 & vname(1,ifield), vpiovar, &
709 & trec, iodesc, vsize, &
710 & lbi, ubi, lbj, ubj, lbk, ubk, &
711 & fscale(ifield,ng), fmin, fmax, &
712# ifdef MASKING
713 & fmask, &
714# endif
715# ifdef CHECKSUM
716 & fout(:,:,:,tindex), &
717 & checksum = fhash)
718# else
719 & fout(:,:,:,tindex))
720# endif
721 finfo(8,ifield,ng)=fmin
722 finfo(9,ifield,ng)=fmax
723 END IF
724 ELSE
725 CALL pio_netcdf_get_fvar (ng, model, ncfile, &
726 & vname(1,ifield), fval, &
727 & piofile = piofile, &
728 & start = (/trec/), &
729 & total = (/1/))
730 fval=fval*fscale(ifield,ng)
731 fpoint(tindex,ifield,ng)=fval
732 fmin=fval
733 fmax=fval
734 END IF
735 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
736 IF (master) THEN
737 WRITE (stdout,40) trim(vname(1,ifield)), trec
738 END IF
739 RETURN
740 END IF
741 IF (master) THEN
742 IF (irec.gt.1) THEN
743 WRITE (stdout,50) trim(vname(2,ifield)), ng, fmin, fmax
744 ELSE
745 lstr=scan(ncfile,'/',back=.true.)+1
746 lend=len_trim(ncfile)
747 lvar=min(46,len_trim(vname(2,ifield)))
748 tsec=tval*day2sec
749 CALL time_string (tsec, t_code)
750 WRITE (stdout,60) vname(2,ifield)(1:lvar), t_code, &
751 & ng, trec, tindex, ncfile(lstr:lend), &
752 & tmin, tmax, tval, fmin, fmax
753 END IF
754# ifdef CHECKSUM
755 WRITE (stdout,70) fhash
756# endif
757 END IF
758 update=.true.
759 END IF
760 END IF
761!
762! Increment the local time variable "Tmono" by the interval between
763! snapshots. If the interval is negative, indicating cycling, add in
764! a cycle length. Load time value (sec) into "Tintrp" which used
765! during interpolation between snapshots.
766!
767 IF (.not.lonerec) THEN
768 tdelta=vtime(tindex,ifield,ng)-vtime(3-tindex,ifield,ng)
769 IF (liocycle.and.(tdelta.lt.0.0_r8)) THEN
770 tdelta=tdelta+clength
771 END IF
772 tmono=tmono+tdelta*day2sec
773 finfo(7,ifield,ng)=tmono
774 tintrp(tindex,ifield,ng)=tmono
775 END IF
776 END IF
777!
778 10 FORMAT (/,' GET_3DFLD_PIO - unable to find dimension ',a, &
779 & /,17x,'for variable: ',a,/,17x,'in file: ',a, &
780 & /,17x,'file is not CF compliant...')
781 20 FORMAT (/,' GET_3DFLD_PIO - unable to find requested variable:', &
782 & 1x,a,/,17x,'in input NetCDF file: ',a)
783 30 FORMAT (/,' GET_3DFLD_PIO - unable to open input NetCDF', &
784 & ' file: ',a)
785 40 FORMAT (/,' GET_3DFLD_PIO - error while reading variable: ',a, &
786 & 2x,' at TIME index = ',i0)
787 50 FORMAT (2x,'GET_3DFLD_PIO - ',a,/,22x,'(Grid = ',i2.2, &
788 & ', Min = ',1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
789 60 FORMAT (2x,'GET_3DFLD_PIO - ',a,',',t75,a,/,22x, &
790 & '(Grid=',i2.2,', Rec=',i0,', Index=',i1, &
791 & ', File: ',a,')',/,22x, &
792 & '(Tmin= ', f15.4, ' Tmax= ', f15.4,')', &
793 & t71, 't = ', f15.4 ,/,22x, &
794 & '(Min = ', 1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
795# ifdef CHECKSUM
796 70 FORMAT (22x,'(CheckSum = ',i0,')')
797# endif
798!
799 RETURN
800 END SUBROUTINE get_3dfld_pio
801# endif
802#else
803 SUBROUTINE get_3dfld
804 RETURN
805 END SUBROUTINE get_3dfld
806#endif
subroutine get_3dfld_nf90(ng, model, ifield, ncid, nfiles, s, update, lbi, ubi, lbj, ubj, lbk, ubk, iout, irec, fmask, fout)
Definition get_3dfld.F:139
subroutine get_3dfld(ng, model, ifield, ncid, piofile, nfiles, s, update, lbi, ubi, lbj, ubj, lbk, ubk, iout, irec, fmask, fout)
Definition get_3dfld.F:14
subroutine get_3dfld_pio(ng, model, ifield, piofile, nfiles, s, update, lbi, ubi, lbj, ubj, lbk, ubk, iout, irec, fmask, fout)
Definition get_3dfld.F:450
subroutine, public time_string(mytime, date_string)
Definition dateclock.F:1272
character(len=256) ncfile
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
character(len=256), dimension(:,:), allocatable cinfo
logical, dimension(:,:,:), allocatable linfo
integer, parameter io_pio
Definition mod_ncparam.F:96
type(my_vardesc), dimension(:,:,:), pointer dinfo
real(dp), dimension(:,:,:), allocatable vtime
real(dp), dimension(:,:,:), allocatable tintrp
real(dp), dimension(:,:,:), allocatable fpoint
real(dp), dimension(:,:), allocatable fscale
character(len=46), dimension(0:nv) tname
real(dp), dimension(:,:,:), allocatable finfo
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:,:,:), allocatable iinfo
logical master
integer, parameter b3dvar
Definition mod_param.F:725
integer, parameter r3dvar
Definition mod_param.F:721
integer, parameter u3dvar
Definition mod_param.F:722
integer, parameter w3dvar
Definition mod_param.F:724
integer, parameter l3dvar
Definition mod_param.F:726
integer, parameter p3dvar
Definition mod_param.F:720
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
type(io_desc_t), dimension(:), pointer iodesc_dp_p3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_b3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_l3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_l3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_b3dvar
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_sp_p3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dvar
real(dp), parameter day2sec
real(dp), parameter spval
integer, dimension(:), allocatable iic
type(t_clock) rclock
integer exit_flag
real(dp), dimension(:), allocatable time
integer, dimension(:), allocatable ntstart
integer noerror
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52