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