ROMS
Loading...
Searching...
No Matches
get_3dfldr.F
Go to the documentation of this file.
1#include "cppdefs.h"
2#if defined ADJOINT && defined SOLVE3D
3 SUBROUTINE get_3dfldr (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. Backward 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 for backward time processing according
93! to IO type.
94!-----------------------------------------------------------------------
95!
96 SELECT CASE (s(1)%IOtype)
97 CASE (io_nf90)
98 CALL get_3dfldr_nf90 (ng, model, ifield, ncid, &
99 & nfiles, s, update, &
100 & lbi, ubi, lbj, ubj, lbk, ubk, &
101 & iout, irec, &
102# ifdef MASKING
103 & fmask, &
104# endif
105 & fout)
106
107# if defined PIO_LIB && defined DISTRIBUTE
108 CASE (io_pio)
109 CALL get_3dfldr_pio (ng, model, ifield, piofile, &
110 & nfiles, s, update, &
111 & lbi, ubi, lbj, ubj, lbk, ubk, &
112 & iout, irec, &
113# ifdef MASKING
114 & fmask, &
115# endif
116 & fout)
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_3DFLDR - Illegal inpput type, io_type = ',i0, &
125 & /,14x,'Check KeyWord ''INP_LIB'' in ''roms.in''.')
126!
127 RETURN
128 END SUBROUTINE get_3dfldr
129!
130!***********************************************************************
131 SUBROUTINE get_3dfldr_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_3dfldr_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( 5,ifield,ng).or. &
211 & (finfo( 1,ifield,ng)*day2sec.gt.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(5,ifield,ng) switch after
220! the query for the LOWER 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 UPPER 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(5,ifield,ng)) THEN
231 linfo(5,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.gt.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 IF (trec.le.0) trec=nrec+trec
271 ELSE
272 trec=trec-1
273 END IF
274 iinfo(9,ifield,ng)=trec
275!
276 IF ((1.le.trec).and.(trec.le.nrec)) THEN
277!
278! Set rolling index for two-time record storage of input data. If
279! "Iout" is unity, input data is stored in timeless array by the
280! calling program. If Irec > 1, this routine is used to read a 3D
281! field varying in another non-time dimension.
282!
283 IF (irec.eq.1) THEN
284 IF (iout.eq.1) THEN
285 tindex=1
286 ELSE
287 tindex=3-tindex
288 END IF
289 iinfo(8,ifield,ng)=tindex
290 END IF
291!
292! Read in time coordinate and scale it to day units.
293!
294 IF (tid.ge.0) THEN
295 CALL netcdf_get_time (ng, model, ncfile, tname(ifield), &
296 & rclock%DateNumber, tval, &
297 & ncid = ncid, &
298 & start = (/trec/), &
299 & total = (/1/))
300 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
301 IF (master) WRITE (stdout,40) trim(tname(ifield)), trec
302 RETURN
303 END IF
304 END IF
305 tval=tval*tscale
306 vtime(tindex,ifield,ng)=tval
307!
308! Activate switch Linfo(5,ifield,ng) if processing the FIRST record of
309! the file for the UPPER time snapshot. We need to get the LOWER time
310! snapshot from previous multifile.
311!
312 IF ((trec.eq.1).and.(tval*day2sec.ge.time(ng))) THEN
313 linfo(5,ifield,ng)=.true.
314 END IF
315!
316! Read in 3D-grided or point data.
317!
318 IF (vid.ge.0) THEN
319 IF (lgridded) THEN
320 fmin=spval
321 fmax=-spval
322 IF (irec.gt.1) THEN
323 DO i=1,irec
324 status=nf_fread3d(ng, model, ncfile, ncid, &
325 & vname(1,ifield), vid, &
326 & i, vtype, vsize, &
327 & lbi, ubi, lbj, ubj, lbk, ubk, &
328 & fscale(ifield,ng), &
329 & my_fmin, my_fmax, &
330# ifdef MASKING
331 & fmask, &
332# endif
333# ifdef CHECKSUM
334 & fout(:,:,:,i), &
335 & checksum = fhash)
336# else
337 & fout(:,:,:,i))
338# endif
339 fmin=min(fmin,my_fmin)
340 fmax=max(fmax,my_fmax)
341 END DO
342 finfo(8,ifield,ng)=fmin
343 finfo(9,ifield,ng)=fmax
344 ELSE
345 status=nf_fread3d(ng, model, ncfile, ncid, &
346 & vname(1,ifield), vid, &
347 & trec, vtype, vsize, &
348 & lbi, ubi, lbj, ubj, lbk, ubk, &
349 & fscale(ifield,ng), fmin, fmax, &
350# ifdef MASKING
351 & fmask, &
352# endif
353# ifdef CHECKSUM
354 & fout(:,:,:,tindex), &
355 & checksum = fhash)
356# else
357 & fout(:,:,:,tindex))
358# endif
359 finfo(8,ifield,ng)=fmin
360 finfo(9,ifield,ng)=fmax
361 END IF
362 ELSE
363 CALL netcdf_get_fvar (ng, model, ncfile, vname(1,ifield), &
364 & fval, &
365 & ncid = ncid, &
366 & start = (/trec/), &
367 & total = (/1/))
368 fval=fval*fscale(ifield,ng)
369 fpoint(tindex,ifield,ng)=fval
370 fmin=fval
371 fmax=fval
372 END IF
373 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
374 IF (master) THEN
375 WRITE (stdout,40) trim(vname(1,ifield)), trec
376 END IF
377 RETURN
378 END IF
379 IF (master) THEN
380 IF (irec.gt.1) THEN
381 WRITE (stdout,50) trim(vname(2,ifield)), ng, fmin, fmax
382 ELSE
383 lstr=scan(ncfile,'/',back=.true.)+1
384 lend=len_trim(ncfile)
385 lvar=min(46,len_trim(vname(2,ifield)))
386 tsec=tval*day2sec
387 CALL time_string (tsec, t_code)
388 WRITE (stdout,60) vname(2,ifield)(1:lvar), t_code, &
389 & ng, trec, tindex, ncfile(lstr:lend), &
390 & tmin, tmax, tval, fmin, fmax
391 END IF
392#ifdef CHECKSUM
393 WRITE (stdout,70) fhash
394#endif
395 END IF
396 update=.true.
397 END IF
398 END IF
399!
400! Decrease the local time variable "Tmono" by the interval between
401! snapshots. If the interval is negative, indicating cycling, add in
402! a cycle length. Load time value (sec) into "Tintrp" which used
403! during interpolation between snapshots.
404!
405 IF (.not.lonerec) THEN
406 tdelta=vtime(3-tindex,ifield,ng)-vtime(tindex,ifield,ng)
407 IF (liocycle.and.(tdelta.lt.0.0_r8)) THEN
408 tdelta=tdelta+clength
409 END IF
410 tmono=tmono-tdelta*day2sec
411 finfo(7,ifield,ng)=tmono
412 tintrp(tindex,ifield,ng)=tmono
413 END IF
414 END IF
415!
416 10 FORMAT (/,' GET_3DFLDR_NF90 - unable to find dimension ',a, &
417 & /,19x,'for variable: ',a,/,19x,'in file: ',a, &
418 & /,19x,'file is not CF compliant...')
419 20 FORMAT (/,' GET_3DFLDR_NF90 - unable to find requested', &
420 & ' variable: ',a,/,19x,'in input NetCDF file: ',a)
421 30 FORMAT (/,' GET_3DFLDR_NF90 - unable to open input NetCDF', &
422 & ' file: ',a)
423 40 FORMAT (/,' GET_3DFLDR_NF90 - error while reading variable: ',a, &
424 & 2x,' at TIME index = ',i0)
425 50 FORMAT (2x,'GET_3DFLDR_NF90 - ',a,/,22x,'(Grid = ',i2.2, &
426 & ', Min = ',1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
427 60 FORMAT (2x,'GET_3DFLDR_NF90 - ',a,',',t75,a,/,22x, &
428 & '(Grid=',i2.2,', Rec=',i0,', Index=',i1, &
429 & ', File: ',a,')',/,22x, &
430 & '(Tmin= ', f15.4, ' Tmax= ', f15.4,')', &
431 & t71, 't = ', f15.4 ,/,22x, &
432 & '(Min = ', 1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
433# ifdef CHECKSUM
434 70 FORMAT (22x,'(CheckSum = ',i0,')')
435# endif
436!
437 RETURN
438 END SUBROUTINE get_3dfldr_nf90
439
440# if defined PIO_LIB && defined DISTRIBUTE
441!
442!***********************************************************************
443 SUBROUTINE get_3dfldr_pio (ng, model, ifield, pioFile, &
444 & nfiles, S, update, &
445 & LBi, UBi, LBj, UBj, LBk, UBk, &
446 & Iout, Irec, &
447# ifdef MASKING
448 & Fmask, &
449# endif
450 & Fout)
451!***********************************************************************
452!
453 USE mod_param
454 USE mod_parallel
455 USE mod_iounits
456 USE mod_ncparam
458 USE mod_scalars
459!
460 USE dateclock_mod, ONLY : time_string
461 USE inquiry_mod, ONLY : inquiry
462 USE nf_fread3d_mod, ONLY : nf_fread3d
463 USE strings_mod, ONLY : founderror
464!
465 implicit none
466!
467! Imported variable declarations.
468!
469 logical, intent(out) :: update
470!
471 integer, intent(in) :: ng, model, ifield, nfiles
472 integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk
473 integer, intent(in) :: Iout, Irec
474!
475 TYPE (File_desc_t), intent(inout) :: pioFile
476 TYPE(t_io), intent(inout) :: S(nfiles)
477!
478# ifdef MASKING
479 real(r8), intent(in) :: Fmask(LBi:UBi,LBj:UBj)
480# endif
481 real(r8), intent(inout) :: Fout(LBi:UBi,LBj:UBj,LBk:UBk,Iout)
482!
483! Local variable declarations.
484!
485 logical :: Lgridded, Linquire, Liocycle, Lmulti, Lonerec
486!
487 integer :: Nrec, Tindex, Trec, Vtype
488 integer :: i, job, lend, lstr, lvar, status
489 integer :: Vsize(4)
490# ifdef CHECKSUM
491 integer(i8b) :: Fhash
492# endif
493!
494 real(r8) :: Fmax, Fmin, Fval, my_Fmin, my_Fmax
495
496 real(dp) :: Clength, Tdelta, Tend
497 real(dp) :: Tmax, Tmin, Tmono, Tscale, Tstr
498 real(dp) :: Tsec, Tval
499!
500 character (len=22) :: t_code
501
502 character (len=*), parameter :: MyFile = &
503 & __FILE__//", get_3dfldr_pio"
504!
505 TYPE (IO_Desc_t), pointer :: ioDesc
506 TYPE (My_VarDesc) :: TpioVar, VpioVar
507!
508 sourcefile=myfile
509!
510!-----------------------------------------------------------------------
511! Initialize.
512!-----------------------------------------------------------------------
513!
514 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
515!
516! Determine if inquiring about the field to process in input NetCDF
517! file(s). This usually happens on first call or when the field
518! time records are split (saved) in several multi-files.
519!
520 linquire=.false.
521 lmulti=.false.
522 IF (iic(ng).eq.0) linquire=.true.
523 IF (.not.linquire.and. &
524 & ((iinfo(10,ifield,ng).gt.1).and. &
525 & (linfo( 5,ifield,ng).or. &
526 & (finfo( 1,ifield,ng)*day2sec.gt.time(ng))))) THEN
527 linquire=.true.
528 lmulti=.true.
529 END IF
530!
531! If appropriate, inquire about the contents of input NetCDF file and
532! fill information arrays.
533!
534! Also, if appropriate, deactivate the Linfo(5,ifield,ng) switch after
535! the query for the LOWER snapshot interpolant from previous multifile
536! in the list. The switch was activated previously to indicate the
537! processing of the FIRST record of the file for the UPPER snapshot
538! interpolant.
539!
540 IF (linquire) THEN
541 job=-3
542 CALL inquiry (ng, model, job, iout, irec, 1, ifield, piofile, &
543 & lmulti, nfiles, s)
544 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
545 IF (linfo(5,ifield,ng)) THEN
546 linfo(5,ifield,ng)=.false.
547 END IF
548 END IF
549!
550!-----------------------------------------------------------------------
551! If appropriate, read in new data.
552!-----------------------------------------------------------------------
553!
554# ifdef CHECKSUM
555 fhash=0_i8b
556# endif
557 update=.false.
558 tmono=finfo(7,ifield,ng)
559!
560 IF ((tmono.gt.time(ng)).or.(iic(ng).eq.0).or. &
561 & (iic(ng).eq.ntstart(ng))) THEN
562!
563! Load properties for requested field from information arrays.
564!
565 lgridded=linfo(1,ifield,ng)
566 liocycle=linfo(2,ifield,ng)
567 lonerec =linfo(3,ifield,ng)
568 vtype =iinfo(1,ifield,ng)
569 vpiovar =dinfo(1,ifield,ng)
570 tpiovar =dinfo(2,ifield,ng)
571 nrec =iinfo(4,ifield,ng)
572 vsize(1)=iinfo(5,ifield,ng)
573 vsize(2)=iinfo(6,ifield,ng)
574 vsize(3)=iinfo(7,ifield,ng)
575 tindex =iinfo(8,ifield,ng)
576 trec =iinfo(9,ifield,ng)
577 tmin =finfo(1,ifield,ng)
578 tmax =finfo(2,ifield,ng)
579 clength =finfo(5,ifield,ng)
580 tscale =finfo(6,ifield,ng)
581 ncfile =cinfo(ifield,ng)
582!
583 IF (liocycle) THEN
584 trec=mod(trec,nrec)-1
585 IF (trec.le.0) trec=nrec+trec
586 ELSE
587 trec=trec-1
588 END IF
589 iinfo(9,ifield,ng)=trec
590!
591 IF ((1.le.trec).and.(trec.le.nrec)) THEN
592!
593! Set rolling index for two-time record storage of input data. If
594! "Iout" is unity, input data is stored in timeless array by the
595! calling program. If Irec > 1, this routine is used to read a 3D
596! field varying in another non-time dimension.
597!
598 IF (irec.eq.1) THEN
599 IF (iout.eq.1) THEN
600 tindex=1
601 ELSE
602 tindex=3-tindex
603 END IF
604 iinfo(8,ifield,ng)=tindex
605 END IF
606!
607! Read in time coordinate and scale it to day units.
608!
609 IF (tpiovar%vd%varID.ge.0) THEN
610 CALL pio_netcdf_get_time (ng, model, ncfile, tname(ifield), &
611 & rclock%DateNumber, tval, &
612 & piofile = piofile, &
613 & start = (/trec/), &
614 & total = (/1/))
615 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
616 IF (master) WRITE (stdout,40) trim(tname(ifield)), trec
617 RETURN
618 END IF
619 END IF
620 tval=tval*tscale
621 vtime(tindex,ifield,ng)=tval
622!
623! Activate switch Linfo(5,ifield,ng) if processing the FIRST record of
624! the file for the UPPER time snapshot. We need to get the LOWER time
625! snapshot from previous multifile.
626!
627 IF ((trec.eq.1).and.(tval*day2sec.ge.time(ng))) THEN
628 linfo(5,ifield,ng)=.true.
629 END IF
630!
631! Read in 3D-grided or point data.
632!
633 IF (vpiovar%vd%varID.ge.0) THEN
634 SELECT CASE (vtype)
635# ifdef SEDIMENT
636 CASE (b3dvar)
637 IF (kind(fout).eq.8) THEN
638 iodesc => iodesc_dp_b3dvar(ng)
639 ELSE
640 iodesc => iodesc_sp_b3dvar(ng)
641 END IF
642# endif
643# if defined DIAGNOSTICS_BIO && defined ECOSIM
644 CASE (l3dvar)
645 IF (kind(fout).eq.8) THEN
646 iodesc => iodesc_dp_l3dvar(ng)
647 ELSE
648 iodesc => iodesc_sp_l3dvar(ng)
649 END IF
650# endif
651 CASE (p3dvar)
652 IF (kind(fout).eq.8) THEN
653 iodesc => iodesc_dp_p3dvar(ng)
654 ELSE
655 iodesc => iodesc_sp_p3dvar(ng)
656 END IF
657 CASE (r3dvar)
658 IF (kind(fout).eq.8) THEN
659 iodesc => iodesc_dp_r3dvar(ng)
660 ELSE
661 iodesc => iodesc_sp_r3dvar(ng)
662 END IF
663 CASE (u3dvar)
664 IF (kind(fout).eq.8) THEN
665 iodesc => iodesc_dp_u3dvar(ng)
666 ELSE
667 iodesc => iodesc_sp_u3dvar(ng)
668 END IF
669 CASE (v3dvar)
670 IF (kind(fout).eq.8) THEN
671 iodesc => iodesc_dp_v3dvar(ng)
672 ELSE
673 iodesc => iodesc_sp_v3dvar(ng)
674 END IF
675 CASE (w3dvar)
676 IF (kind(fout).eq.8) THEN
677 iodesc => iodesc_dp_w3dvar(ng)
678 ELSE
679 iodesc => iodesc_sp_w3dvar(ng)
680 END IF
681 END SELECT
682!
683 IF (lgridded) THEN
684 fmin=spval
685 fmax=-spval
686 IF (irec.gt.1) THEN
687 DO i=1,irec
688 status=nf_fread3d(ng, model, ncfile, piofile, &
689 & vname(1,ifield), vpiovar, &
690 & i, iodesc, vsize, &
691 & lbi, ubi, lbj, ubj, lbk, ubk, &
692 & fscale(ifield,ng), &
693 & my_fmin, my_fmax, &
694# ifdef MASKING
695 & fmask, &
696# endif
697# ifdef CHECKSUM
698 & fout(:,:,:,i), &
699 & checksum = fhash)
700# else
701 & fout(:,:,:,i))
702# endif
703 fmin=min(fmin,my_fmin)
704 fmax=min(fmax,my_fmax)
705 END DO
706 finfo(8,ifield,ng)=fmin
707 finfo(9,ifield,ng)=fmax
708 ELSE
709 status=nf_fread3d(ng, model, ncfile, piofile, &
710 & vname(1,ifield), vpiovar, &
711 & trec, iodesc, vsize, &
712 & lbi, ubi, lbj, ubj, lbk, ubk, &
713 & fscale(ifield,ng), fmin, fmax, &
714# ifdef MASKING
715 & fmask, &
716# endif
717# ifdef CHECKSUM
718 & fout(:,:,:,tindex), &
719 & checksum = fhash)
720# else
721 & fout(:,:,:,tindex))
722# endif
723 finfo(8,ifield,ng)=fmin
724 finfo(9,ifield,ng)=fmax
725 END IF
726 ELSE
727 CALL pio_netcdf_get_fvar (ng, model, ncfile, &
728 & vname(1,ifield), fval, &
729 & piofile = piofile, &
730 & start = (/trec/), &
731 & total = (/1/))
732 fval=fval*fscale(ifield,ng)
733 fpoint(tindex,ifield,ng)=fval
734 fmin=fval
735 fmax=fval
736 END IF
737 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
738 IF (master) THEN
739 WRITE (stdout,40) trim(vname(1,ifield)), trec
740 END IF
741 RETURN
742 END IF
743 IF (master) THEN
744 IF (irec.gt.1) THEN
745 WRITE (stdout,50) trim(vname(2,ifield)), ng, fmin, fmax
746 ELSE
747 lstr=scan(ncfile,'/',back=.true.)+1
748 lend=len_trim(ncfile)
749 lvar=min(46,len_trim(vname(2,ifield)))
750 tsec=tval*day2sec
751 CALL time_string (tsec, t_code)
752 WRITE (stdout,60) vname(2,ifield)(1:lvar), t_code, &
753 & ng, trec, tindex, ncfile(lstr:lend), &
754 & tmin, tmax, tval, fmin, fmax
755 END IF
756# ifdef CHECKSUM
757 WRITE (stdout,70) fhash
758# endif
759 END IF
760 update=.true.
761 END IF
762 END IF
763!
764! Decrease the local time variable "Tmono" by the interval between
765! snapshots. If the interval is negative, indicating cycling, add in
766! a cycle length. Load time value (sec) into "Tintrp" which used
767! during interpolation between snapshots.
768!
769 IF (.not.lonerec) THEN
770 tdelta=vtime(3-tindex,ifield,ng)-vtime(tindex,ifield,ng)
771 IF (liocycle.and.(tdelta.lt.0.0_r8)) THEN
772 tdelta=tdelta+clength
773 END IF
774 tmono=tmono-tdelta*day2sec
775 finfo(7,ifield,ng)=tmono
776 tintrp(tindex,ifield,ng)=tmono
777 END IF
778 END IF
779!
780 10 FORMAT (/,' GET_3DFLDR_PIO - unable to find dimension ',a, &
781 & /,18x,'for variable: ',a,/,18x,'in file: ',a, &
782 & /,18x,'file is not CF compliant...')
783 20 FORMAT (/,' GET_3DFLDR_PIO - unable to find requested', &
784 & ' variable: ',a,/,18x,'in input NetCDF file: ',a)
785 30 FORMAT (/,' GET_3DFLDR_PIO - unable to open input NetCDF', &
786 & ' file: ',a)
787 40 FORMAT (/,' GET_3DFLDR_PIO - error while reading variable: ',a, &
788 & 2x,' at TIME index = ',i0)
789 50 FORMAT (2x,'GET_3DFLDR_PIO - ',a,/,22x,'(Grid = ',i2.2, &
790 & ', Min = ',1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
791 60 FORMAT (2x,'GET_3DFLDR_PIO - ',a,',',t75,a,/,22x, &
792 & '(Grid=',i2.2,', Rec=',i0,', Index=',i1, &
793 & ', File: ',a,')',/,22x, &
794 & '(Tmin= ', f15.4, ' Tmax= ', f15.4,')', &
795 & t71, 't = ', f15.4 ,/,22x, &
796 & '(Min = ', 1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
797# ifdef CHECKSUM
798 70 FORMAT (22x,'(CheckSum = ',i0,')')
799# endif
800!
801 RETURN
802 END SUBROUTINE get_3dfldr_pio
803# endif
804#else
805 SUBROUTINE get_3dfldr
806 RETURN
807 END SUBROUTINE get_3dfldr
808#endif
subroutine get_3dfldr_pio(ng, model, ifield, piofile, nfiles, s, update, lbi, ubi, lbj, ubj, lbk, ubk, iout, irec, fmask, fout)
Definition get_3dfldr.F:451
subroutine get_3dfldr(ng, model, ifield, ncid, piofile, nfiles, s, update, lbi, ubi, lbj, ubj, lbk, ubk, iout, irec, fmask, fout)
Definition get_3dfldr.F:14
subroutine get_3dfldr_nf90(ng, model, ifield, ncid, nfiles, s, update, lbi, ubi, lbj, ubj, lbk, ubk, iout, irec, fmask, fout)
Definition get_3dfldr.F:139
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