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

Data Types

interface  get_cycle
 

Functions/Subroutines

subroutine get_cycle_nf90 (ng, model, ifield, job, lmulti, ncname, ncid, tvarname, ntime, smday, tid, lcycle, clength, tindex, tstr, tend, tmin, tmax, tscale)
 
subroutine get_cycle_pio (ng, model, ifield, job, lmulti, ncname, piofile, tvarname, ntime, smday, tpiovar, lcycle, clength, tindex, tstr, tend, tmin, tmax, tscale)
 

Function/Subroutine Documentation

◆ get_cycle_nf90()

subroutine get_cycle_mod::get_cycle_nf90 ( integer, intent(in) ng,
integer, intent(in) model,
integer, intent(in) ifield,
integer, intent(in) job,
logical, intent(in) lmulti,
character (len=*), intent(in) ncname,
integer, intent(in) ncid,
character (len=*), intent(in) tvarname,
integer, intent(in) ntime,
real(dp), intent(in) smday,
integer, intent(out) tid,
logical, intent(out) lcycle,
real(dp), intent(out) clength,
integer, intent(out) tindex,
real(dp), intent(out) tstr,
real(dp), intent(out) tend,
real(dp), intent(out) tmin,
real(dp), intent(out) tmax,
real(dp), intent(out) tscale )

Definition at line 92 of file get_cycle.F.

97!***********************************************************************
98!
99 USE mod_netcdf
100!
101! Imported variable declarations.
102!
103 logical, intent(in) :: Lmulti
104 logical, intent(out) :: Lcycle
105
106 integer, intent(in) :: ng, model, ifield, job, ncid, ntime
107
108 integer, intent(out) :: Tindex, Tid
109
110 real(dp), intent(in) :: smday
111
112 real(dp), intent(out) :: Tmax, Tmin, Tend, Tscale, Tstr, clength
113
114 character (len=*), intent(in) :: ncname
115 character (len=*), intent(in) :: TvarName
116!
117! Local variable declarations.
118!
119 logical :: TimeLess
120 logical :: Linside, LowerBound, Upperbound
121!
122 integer :: i, nvatt, nvdim
123!
124 real(dp) :: mday, tstart
125 real(dp) :: Tval(ntime)
126!
127 character (len=40) :: tunits
128
129 character (len=*), parameter :: MyFile = &
130 & __FILE__
131!
132 sourcefile=myfile
133!
134!-----------------------------------------------------------------------
135! Find time cycling parameters, if any.
136!-----------------------------------------------------------------------
137!
138! Initialize.
139!
140 lcycle=.false.
141 linside=.false.
142 lowerbound=.false.
143 timeless=.false.
144 upperbound=.false.
145 tindex=0
146 clength=0.0_r8
147 tstr=0.0_r8
148 tscale=1.0_r8
149!
150! Determine if the associated time variable is actually a time
151! variable.
152!
153 IF ((index(trim(tvarname),'period').gt.0).or. &
154 & (trim(tvarname).eq.'river')) THEN
155 timeless=.true.
156 END IF
157!
158! Inquire input NetCDF file about the time variable.
159!
160 IF (ntime.ge.1) THEN
161 CALL netcdf_inq_var (ng, model, ncname, &
162 & ncid = ncid, &
163 & myvarname = trim(tvarname), &
164 & varid = tid, &
165 & nvardim = nvdim, &
166 & nvaratt = nvatt)
167 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
168!
169! Check if time cycling attribute is present and then read in time
170! cycle length. Check time coordinate units and determine time
171! scale. The internal processing of all fields requires time in
172! day units. Check if more than one time record is available.
173!
174 DO i=1,nvatt
175 IF (trim(var_aname(i)).eq.'cycle_length') THEN
176 lcycle=.true.
177 IF (var_afloat(i).gt.0.0_r8) THEN
178 clength=var_afloat(i)
179 ELSE IF (var_aint(i).gt.0) THEN ! no CF compliance
180 clength=real(var_aint(i),r8) ! attribute is an integer
181 ELSE
182 IF (master) WRITE (stdout,10) trim(var_aname(i)), &
183 & trim(tvarname)
184 exit_flag=2
185 RETURN
186 END IF
187 ELSE IF (trim(var_aname(i)).eq.'units') THEN
188 tunits=trim(var_achar(i))
189 IF (tunits(1:6).eq.'second') THEN
190 tscale=sec2day
191 END IF
192 END IF
193 END DO
194 END IF
195!
196! Read in time variable values.
197!
198 CALL netcdf_get_time (ng, model, ncname, tvarname, &
199 & rclock%DateNumber, tval, &
200 & ncid = ncid, &
201 & start = (/1/), &
202 & total = (/ntime/), &
203 & min_val = tmin, &
204 & max_val = tmax)
205 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
206!
207! Scale time variable to days. Determine the minimum and maximum time
208! values available
209!
210 DO i=1,ntime
211 tval(i)=tval(i)*tscale
212 END DO
213 tmin=tmin*tscale
214 tmax=tmax*tscale
215 IF (lcycle) THEN
216 mday=mod(smday,clength)
217 ELSE
218 mday=smday
219 END IF
220!
221! Is the model time inside the data time range? If not, check if the
222! data just has the LOWER- or the UPPER-snapshot interpolant.
223!
224 IF ((tmin.le.mday).and.(mday.le.tmax)) THEN
225 linside=.true.
226 ELSE IF (mday.ge.tmax) THEN
227 lowerbound=.true.
228 ELSE IF (mday.le.tmin) THEN
229 upperbound=.true.
230 END IF
231!
232! If processing timeless variable, set processing parameter.
233!
234 IF (timeless) THEN
235 tstr=tmin
236 tend=tmax
237 tindex=ntime
238 RETURN
239 END IF
240!
241! If processing field split in several files, find UPPER time-snapshot
242! and its associated time record (Tindex).
243!
244 IF (lmulti) THEN
245 IF (job.gt.0) THEN ! forward time-stepping logic
246 DO i=1,ntime
247 IF (tval(i).gt.mday) THEN
248 tindex=i-1
249 tend=tval(i)
250 EXIT
251 END IF
252 END DO
253 ELSE ! backward time-stepping logic
254 DO i=ntime,1,-1
255 IF (tval(i).le.mday) THEN
256 tindex=i+1
257 tstr=tval(i)
258 EXIT
259 END IF
260 END DO
261 END IF
262!
263! If not processing a multi-file field, find LOWER time-snapshot
264! and its associated time record (Tindex).
265!
266 ELSE
267 IF (job.gt.0) THEN ! forward: Tval(i) =< mday =< Tval(i+1)
268 IF (linside) THEN
269 tstart=tmin
270 IF (ntime.eq.1) THEN
271 tstr=tmin
272 tindex=1
273 ELSE
274 DO i=2,ntime
275 IF ((tstart.le.mday).and.(mday.le.tval(i))) THEN
276 tstr=tstart
277 tindex=i-1 ! one is added when processing
278 EXIT
279 END IF
280 tstart=tval(i)
281 END DO
282 END IF
283 ELSE
284 tstr=tmax ! LowerBound for next multifile or
285 tindex=ntime ! time cycling
286 END IF
287 ELSE ! backward: Tval(i) =< mday =< Tval(i+1)
288 IF (linside) THEN
289 tstart=tmin
290 IF (ntime.eq.1) THEN
291 tstr=tmin
292 tindex=1
293 ELSE
294 DO i=2,ntime
295 IF ((tstart.le.mday).and.(mday.le.tval(i))) THEN
296 tstr=tval(i)
297 tindex=i
298 EXIT
299 END IF
300 tstart=tval(i)
301 END DO
302 END IF
303 ELSE IF (lowerbound) THEN
304 tstr=tmax
305 tindex=ntime
306 ELSE IF (upperbound) THEN
307 tstr=tmin
308 tindex=1
309 END IF
310 END IF
311 END IF
312!
313! If processing a multi-file field, set LOWER time-snapshot. It
314! is the last value from previous file. Otherwise, set UPPER
315! time-snapshot.
316!
317 IF (lmulti) THEN
318 IF (job.gt.0) THEN
319 tstr=finfo(2,ifield,ng) ! Tmax from previous file
320 ELSE
321 tend=finfo(1,ifield,ng) ! Tmin from previous file
322 END IF
323 ELSE
324 IF (lcycle.and.(tindex.eq.ntime)) THEN
325 tend=tmin
326 ELSE
327 IF (job.gt.0) THEN
328 i=min(ntime,tindex+1)
329 tend=tval(i)
330 ELSE
331 i=max(1,tindex-1)
332 tend=tval(i)
333 END IF
334 END IF
335 END IF
336!
337! If not cycling, stop execution if there is not field data
338! available for current model time. Avoid check on tidal data
339! since time is in terms of frequencies.
340!
341 IF (.not.timeless) THEN
342 IF (.not.lcycle.and.(ntime.gt.1)) THEN
343 IF (lmulti) THEN
344 IF (job.gt.0) THEN
345 IF (smday.gt.tmax) THEN
346 IF (master) WRITE (stdout,20) trim(tvarname), &
347 & tmax, smday
348 exit_flag=2
349 RETURN
350 ELSE IF (linfo(6,ifield,ng)) THEN
351 IF (tmin.lt.smday) THEN
352 IF (master) THEN
353 WRITE (stdout,30) &
354 & 'Upper snapshot time for multi-file variable:', &
355 & trim(tvarname), &
356 & trim(vname(1,ifield)), &
357 & 'is less than current model time.', &
358 & 'Tmin = ', tmin, smday
359 END IF
360 exit_flag=2
361 RETURN
362 END IF
363 END IF
364 ELSE
365 IF (linfo(5,ifield,ng)) THEN
366 IF (tmin.gt.smday) THEN
367 IF (master) THEN
368 WRITE (stdout,30) &
369 & 'Lower snapshot time for multi-file variable:', &
370 & trim(tvarname), &
371 & trim(vname(1,ifield)), &
372 & 'is greater than current model time.', &
373 & 'Tmin = ', tmin, smday
374 END IF
375 exit_flag=2
376 RETURN
377 END IF
378 ELSE
379 IF (smday.lt.tmax) THEN
380 IF (master) THEN
381 WRITE (stdout,30) &
382 & 'starting time for multi-file variable:', &
383 & trim(tvarname), &
384 & trim(vname(1,ifield)), &
385 & 'is greater than current model time.', &
386 & 'Tmax = ', tmax, smday
387 END IF
388 exit_flag=2
389 RETURN
390 END IF
391 END IF
392 END IF
393 ELSE
394 IF (.not.upperbound.and.(smday.lt.tmin)) THEN
395 IF (master) WRITE (stdout,40) trim(tvarname), &
396 & trim(vname(1,ifield)), &
397 & tmin, smday
398 exit_flag=2
399 RETURN
400 END IF
401 END IF
402 END IF
403 END IF
404!
405 10 FORMAT (/,' GET_CYCLE_NF90 - unable to get value for attribute:', &
406 & 1x,a,/,18x,'in variable: ',a, &
407 & /,18x,'This attribute value is expected to be of', &
408 & /,18x,'the same external type as the variable.')
409 20 FORMAT (/,' GET_CYCLE_NF90 - ending time for multi-file', &
410 & ' variable: ',a,/,18x,'is less than current model time.', &
411 & /,18x,'TMAX = ',f15.4,2x,'TDAYS = ',f15.4)
412 30 FORMAT (/,' GET_CYCLE_NF90 - ',a,1x,a,2x,'(',a,')',/,18x,a, &
413 & /,18x,a,f15.4,2x,'TDAYS = ',f15.4)
414 40 FORMAT (/,' GET_CYCLE_Nf90 - starting time for variable: ',a,2x, &
415 & '(',a,')',/,18x,'is greater than current model time.', &
416 & /,18x,'TMIN = ',f15.4,2x,'TDAYS = ',f15.4)
417!
418 RETURN
integer, dimension(nvara) var_aint
Definition mod_netcdf.F:178
character(len=1024), dimension(nvara) var_achar
Definition mod_netcdf.F:183
real(r8), dimension(nvara) var_afloat
Definition mod_netcdf.F:179
character(len=100), dimension(nvara) var_aname
Definition mod_netcdf.F:181
subroutine, public netcdf_inq_var(ng, model, ncname, ncid, myvarname, searchvar, varid, nvardim, nvaratt)

◆ get_cycle_pio()

subroutine get_cycle_mod::get_cycle_pio ( integer, intent(in) ng,
integer, intent(in) model,
integer, intent(in) ifield,
integer, intent(in) job,
logical, intent(in) lmulti,
character (len=*), intent(in) ncname,
type (file_desc_t), intent(inout) piofile,
character (len=*), intent(in) tvarname,
integer, intent(in) ntime,
real(dp), intent(in) smday,
type (var_desc_t), intent(out) tpiovar,
logical, intent(out) lcycle,
real(dp), intent(out) clength,
integer, intent(out) tindex,
real(dp), intent(out) tstr,
real(dp), intent(out) tend,
real(dp), intent(out) tmin,
real(dp), intent(out) tmax,
real(dp), intent(out) tscale )

Definition at line 424 of file get_cycle.F.

429!***********************************************************************
430!
432!
433! Imported variable declarations.
434!
435 logical, intent(in) :: Lmulti
436 logical, intent(out) :: Lcycle
437!
438 integer, intent(in) :: ng, model, ifield, job, ntime
439
440 integer, intent(out) :: Tindex
441!
442 real(dp), intent(in) :: smday
443
444 real(dp), intent(out) :: Tmax, Tmin, Tend, Tscale, Tstr, clength
445!
446 character (len=*), intent(in) :: ncname
447 character (len=*), intent(in) :: TvarName
448!
449 TYPE (File_desc_t), intent(inout) :: pioFile
450 TYPE (Var_desc_t), intent(out) :: TpioVar
451!
452! Local variable declarations.
453!
454 logical :: TimeLess
455 logical :: Linside, LowerBound, Upperbound
456!
457 integer :: i, nvatt, nvdim
458!
459 real(dp) :: mday, tstart
460 real(dp) :: Tval(ntime)
461!
462 character (len=40) :: tunits
463
464 character (len=*), parameter :: MyFile = &
465 & __FILE__
466!
467 sourcefile=myfile
468!
469!-----------------------------------------------------------------------
470! Find time cycling parameters, if any.
471!-----------------------------------------------------------------------
472!
473! Initialize.
474!
475 lcycle=.false.
476 linside=.false.
477 lowerbound=.false.
478 timeless=.false.
479 upperbound=.false.
480 tindex=0
481 clength=0.0_r8
482 tstr=0.0_r8
483 tscale=1.0_r8
484!
485! Determine if the associated time variable is actually a time
486! variable.
487!
488 IF ((index(trim(tvarname),'period').gt.0).or. &
489 & (trim(tvarname).eq.'river')) THEN
490 timeless=.true.
491 END IF
492!
493! Inquire input NetCDF file about the time variable.
494!
495 IF (ntime.ge.1) THEN
496 CALL pio_netcdf_inq_var (ng, model, ncname, &
497 & piofile = piofile, &
498 & myvarname = trim(tvarname), &
499 & piovar = tpiovar, &
500 & nvardim = nvdim, &
501 & nvaratt = nvatt)
502 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
503!
504! Check if time cycling attribute is present and then read in time
505! cycle length. Check time coordinate units and determine time
506! scale. The internal processing of all fields requires time in
507! day units. Check if more than one time record is available.
508!
509 DO i=1,nvatt
510 IF (trim(var_aname(i)).eq.'cycle_length') THEN
511 lcycle=.true.
512 IF (var_afloat(i).gt.0.0_r8) THEN
513 clength=var_afloat(i)
514 ELSE IF (var_aint(i).gt.0) THEN ! no CF compliance
515 clength=real(var_aint(i),r8) ! attribute is an integer
516 ELSE
517 IF (master) WRITE (stdout,10) trim(var_aname(i)), &
518 & trim(tvarname)
519 exit_flag=2
520 RETURN
521 END IF
522 ELSE IF (trim(var_aname(i)).eq.'units') THEN
523 tunits=trim(var_achar(i))
524 IF (tunits(1:6).eq.'second') THEN
525 tscale=sec2day
526 END IF
527 END IF
528 END DO
529 END IF
530!
531! Read in time variable values.
532!
533 CALL pio_netcdf_get_time (ng, model, ncname, tvarname, &
534 & rclock%DateNumber, tval, &
535 & piofile = piofile, &
536 & start = (/1/), &
537 & total = (/ntime/), &
538 & min_val = tmin, &
539 & max_val = tmax)
540 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
541!
542! Scale time variable to days. Determine the minimum and maximum time
543! values available
544!
545 DO i=1,ntime
546 tval(i)=tval(i)*tscale
547 END DO
548 tmin=tmin*tscale
549 tmax=tmax*tscale
550 IF (lcycle) THEN
551 mday=mod(smday,clength)
552 ELSE
553 mday=smday
554 END IF
555!
556! Is the model time inside the data time range? If not, check if the
557! data just has the LOWER- or the UPPER-snapshot interpolant.
558!
559 IF ((tmin.le.mday).and.(mday.le.tmax)) THEN
560 linside=.true.
561 ELSE IF (mday.ge.tmax) THEN
562 lowerbound=.true.
563 ELSE IF (mday.le.tmin) THEN
564 upperbound=.true.
565 END IF
566!
567! If processing timeless variable, set processing parameter.
568!
569 IF (timeless) THEN
570 tstr=tmin
571 tend=tmax
572 tindex=ntime
573 RETURN
574 END IF
575!
576! If processing field split in several files, find UPPER time-snapshot
577! and its associated time record (Tindex).
578!
579 IF (lmulti) THEN
580 IF (job.gt.0) THEN ! forward time-stepping logic
581 DO i=1,ntime
582 IF (tval(i).gt.mday) THEN
583 tindex=i-1
584 tend=tval(i)
585 EXIT
586 END IF
587 END DO
588 ELSE ! backward time-stepping logic
589 DO i=ntime,1,-1
590 IF (tval(i).le.mday) THEN
591 tindex=i+1
592 tstr=tval(i)
593 EXIT
594 END IF
595 END DO
596 END IF
597!
598! If not processing a multi-file field, find LOWER time-snapshot
599! and its associated time record (Tindex).
600!
601 ELSE
602 IF (job.gt.0) THEN ! forward: Tval(i) =< mday =< Tval(i+1)
603 IF (linside) THEN
604 tstart=tmin
605 IF (ntime.eq.1) THEN
606 tstr=tmin
607 tindex=1
608 ELSE
609 DO i=2,ntime
610 IF ((tstart.le.mday).and.(mday.le.tval(i))) THEN
611 tstr=tstart
612 tindex=i-1 ! one is added when processing
613 EXIT
614 END IF
615 tstart=tval(i)
616 END DO
617 END IF
618 ELSE
619 tstr=tmax ! LowerBound for next multifile or
620 tindex=ntime ! time cycling
621 END IF
622 ELSE ! backward: Tval(i) =< mday =< Tval(i+1)
623 IF (linside) THEN
624 tstart=tmin
625 IF (ntime.eq.1) THEN
626 tstr=tmin
627 tindex=1
628 ELSE
629 DO i=2,ntime
630 IF ((tstart.le.mday).and.(mday.le.tval(i))) THEN
631 tstr=tval(i)
632 tindex=i
633 EXIT
634 END IF
635 tstart=tval(i)
636 END DO
637 END IF
638 ELSE IF (lowerbound) THEN
639 tstr=tmax
640 tindex=ntime
641 ELSE IF (upperbound) THEN
642 tstr=tmin
643 tindex=1
644 END IF
645 END IF
646 END IF
647!
648! If processing a multi-file field, set LOWER time-snapshot. It
649! is the last value from previous file. Otherwise, set UPPER
650! time-snapshot.
651!
652 IF (lmulti) THEN
653 IF (job.gt.0) THEN
654 tstr=finfo(2,ifield,ng) ! Tmax from previous file
655 ELSE
656 tend=finfo(1,ifield,ng) ! Tmin from previous file
657 END IF
658 ELSE
659 IF (lcycle.and.(tindex.eq.ntime)) THEN
660 tend=tmin
661 ELSE
662 IF (job.gt.0) THEN
663 i=min(ntime,tindex+1)
664 tend=tval(i)
665 ELSE
666 i=max(1,tindex-1)
667 tend=tval(i)
668 END IF
669 END IF
670 END IF
671!
672! If not cycling, stop execution if there is not field data
673! available for current model time. Avoid check on tidal data
674! since time is in terms of frequencies.
675!
676 IF (.not.timeless) THEN
677 IF (.not.lcycle.and.(ntime.gt.1)) THEN
678 IF (lmulti) THEN
679 IF (job.gt.0) THEN
680 IF (smday.gt.tmax) THEN
681 IF (master) WRITE (stdout,20) trim(tvarname), &
682 & tmax, smday
683 exit_flag=2
684 RETURN
685 ELSE IF (linfo(6,ifield,ng)) THEN
686 IF (tmin.lt.smday) THEN
687 IF (master) THEN
688 WRITE (stdout,30) &
689 & 'Upper snapshot time for multi-file variable:', &
690 & trim(tvarname), &
691 & trim(vname(1,ifield)), &
692 & 'is less than current model time.', &
693 & 'Tmin = ', tmin, smday
694 END IF
695 exit_flag=2
696 RETURN
697 END IF
698 END IF
699 ELSE
700 IF (linfo(5,ifield,ng)) THEN
701 IF (tmin.gt.smday) THEN
702 IF (master) THEN
703 WRITE (stdout,30) &
704 & 'Lower snapshot time for multi-file variable:', &
705 & trim(tvarname), &
706 & trim(vname(1,ifield)), &
707 & 'is greater than current model time.', &
708 & 'Tmin = ', tmin, smday
709 END IF
710 exit_flag=2
711 RETURN
712 END IF
713 ELSE
714 IF (smday.lt.tmax) THEN
715 IF (master) THEN
716 WRITE (stdout,30) &
717 & 'starting time for multi-file variable:', &
718 & trim(tvarname), &
719 & trim(vname(1,ifield)), &
720 & 'is greater than current model time.', &
721 & 'Tmax = ', tmax, smday
722 END IF
723 exit_flag=2
724 RETURN
725 END IF
726 END IF
727 END IF
728 ELSE
729 IF (.not.upperbound.and.(smday.lt.tmin)) THEN
730 IF (master) WRITE (stdout,40) trim(tvarname), &
731 & trim(vname(1,ifield)), &
732 & tmin, smday
733 exit_flag=2
734 RETURN
735 END IF
736 END IF
737 END IF
738 END IF
739!
740 10 FORMAT (/,' GET_CYCLE_PIO - unable to get value for attribute:', &
741 & 1x,a,/,17x,'in variable: ',a, &
742 & /,17x,'This attribute value is expected to be of', &
743 & /,17x,'the same external type as the variable.')
744 20 FORMAT (/,' GET_CYCLE_PIO - ending time for multi-file', &
745 & ' variable: ',a,/,17x,'is less than current model time.', &
746 & /,17x,'TMAX = ',f15.4,2x,'TDAYS = ',f15.4)
747 30 FORMAT (/,' GET_CYCLE_PIO - ',a,1x,a,2x,'(',a,')',/,17x,a, &
748 & /,17x,a,f15.4,2x,'TDAYS = ',f15.4)
749 40 FORMAT (/,' GET_CYCLE_PIO - starting time for variable: ',a,2x, &
750 & '(',a,')',/,17x,'is greater than current model time.', &
751 & /,17x,'TMIN = ',f15.4,2x,'TDAYS = ',f15.4)
752!
753 RETURN
subroutine, public pio_netcdf_inq_var(ng, model, ncname, piofile, myvarname, searchvar, piovar, nvardim, nvaratt)