ROMS
Loading...
Searching...
No Matches
set_ngfld.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 SUBROUTINE set_ngfld (ng, model, ifield, &
3 & LBi, UBi, UBj, Istr, Iend, Jrec, &
4 & Finp, Fout, update)
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! This routine time-interpolates requested non-grided field from !
14! time snapshots of input data. !
15! !
16! On Input: !
17! !
18! ng Nested grid number. !
19! model Calling model identifier. !
20! ifield Field ID. !
21! LBi Finp/Fout 1st dimension lower-bound value. !
22! UBi Finp/Fout 1st dimension upper-bound value. !
23! UBj Finp/Fout 2nd dimension upper-bound value, if any. !
24! Otherwise, a value of one is expected. !
25! Istr Starting location to process in the 1st dimension. !
26! Iend Ending location to process in the 1st dimension. !
27! Jrec Number of records to process in the 2nd dimenision, !
28! if any, Otherwise, a value of one is expected. !
29! Finp Latest two-snapshopts of field to interpolate. !
30! !
31! On Output: !
32! !
33! Fout Interpolated field. !
34! update Switch indicating successful interpolation. !
35! !
36!=======================================================================
37!
38 USE mod_param
39 USE mod_parallel
40 USE mod_iounits
41 USE mod_ncparam
42 USE mod_scalars
43!
44 implicit none
45!
46! Imported variable declarations.
47!
48 logical, intent(out) :: update
49
50 integer, intent(in) :: ng, model, ifield
51 integer, intent(in) :: LBi, UBi, UBj, Istr, Iend, Jrec
52
53 real(r8), intent(in) :: Finp(LBi:UBi,UBj,2)
54
55 real(r8), intent(out) :: Fout(LBi:UBi,UBj)
56!
57! Local variable declarations.
58!
59 logical :: Lonerec
60
61 integer :: Tindex, i, it1, it2, j
62
63 real(dp) :: SecScale, fac, fac1, fac2
64!
65!----------------------------------------------------------------------
66! Set up requested field from data snapshots.
67!----------------------------------------------------------------------
68!
69! Get requested field information from global storage.
70!
71 lonerec=linfo(3,ifield,ng)
72 tindex=iinfo(8,ifield,ng)
73 update=.true.
74!
75! Set linear, time interpolation factors. Fractional seconds are
76! rounded to the nearest milliseconds integer towards zero in the
77! time interpolation weights.
78!
79 secscale=1000.0_dp ! seconds to milliseconds
80 it1=3-tindex
81 it2=tindex
82 fac1=anint((tintrp(it2,ifield,ng)-time(ng))*secscale,dp)
83 fac2=anint((time(ng)-tintrp(it1,ifield,ng))*secscale,dp)
84!
85! Load time-invariant data. Time interpolation is not necessary.
86!
87 IF (lonerec) THEN
88 DO j=1,jrec
89 DO i=istr,iend
90 fout(i,j)=finp(i,j,tindex)
91 END DO
92 END DO
93!
94! Time-interpolate.
95!
96 ELSE IF (((fac1*fac2).ge.0.0_dp).and.(fac1+fac2).gt.0.0_dp) THEN
97 fac=1.0_dp/(fac1+fac2)
98 fac1=fac*fac1 ! nondimensional
99 fac2=fac*fac2 ! nondimensional
100 DO j=1,jrec
101 DO i=istr,iend
102 fout(i,j)=fac1*finp(i,j,it1)+fac2*finp(i,j,it2)
103 END DO
104 END DO
105!
106! Activate synchronization flag if a new time record needs to be
107! read in at the next time step.
108!
109 IF ((time(ng)+dt(ng)).gt.tintrp(it2,ifield,ng)) THEN
110 synchro_flag(ng)=.true.
111 END IF
112!
113! Unable to interpolate field. Activate error flag to quit.
114!
115 ELSE
116 IF (master) THEN
117 WRITE (stdout,10) trim(vname(1,ifield)), tdays(ng), &
118 & finfo(1,ifield,ng), finfo(2,ifield,ng), &
119 & finfo(3,ifield,ng), finfo(4,ifield,ng), &
120 & tintrp(it1,ifield,ng)*sec2day, &
121 & tintrp(it2,ifield,ng)*sec2day, &
122 & fac1*sec2day/secscale, &
123 & fac2*sec2day/secscale
124 END IF
125 10 FORMAT (/,' SET_NGFLD - current model time', &
126 & ' exceeds ending value for variable: ',a, &
127 & /,14x,'TDAYS = ',f15.4, &
128 & /,14x,'Data Tmin = ',f15.4,2x,'Data Tmax = ',f15.4, &
129 & /,14x,'Data Tstr = ',f15.4,2x,'Data Tend = ',f15.4, &
130 & /,14x,'TINTRP1 = ',f15.4,2x,'TINTRP2 = ',f15.4, &
131 & /,14x,'FAC1 = ',f15.4,2x,'FAC2 = ',f15.4)
132 exit_flag=2
133 update=.false.
134 END IF
135 RETURN
136 END SUBROUTINE set_ngfld
integer stdout
logical, dimension(:,:,:), allocatable linfo
real(dp), dimension(:,:,:), allocatable tintrp
real(dp), dimension(:,:,:), allocatable finfo
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:,:,:), allocatable iinfo
logical master
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable synchro_flag
real(dp), dimension(:), allocatable tdays
real(dp), parameter sec2day
integer exit_flag
real(dp), dimension(:), allocatable time
subroutine set_ngfld(ng, model, ifield, lbi, ubi, ubj, istr, iend, jrec, finp, fout, update)
Definition set_ngfld.F:5