ROMS
Loading...
Searching...
No Matches
set_2dfld.F
Go to the documentation of this file.
1#include "cppdefs.h"
3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine time-interpolates requested 2D field from snapshots !
12! of input data. !
13! !
14!=======================================================================
15!
16 implicit none
17
18 CONTAINS
19!
20!***********************************************************************
21 SUBROUTINE set_2dfld_tile (ng, tile, model, ifield, &
22 & LBi, UBi, LBj, UBj, &
23 & Finp, Fout, update, &
24 & SetBC)
25!***********************************************************************
26!
27 USE mod_param
28 USE mod_parallel
29 USE mod_iounits
30 USE mod_ncparam
31 USE mod_scalars
32!
34#ifdef DISTRIBUTE
36#endif
37!
38! Imported variable declarations.
39!
40 logical, intent(in), optional :: SetBC
41
42 logical, intent(out) :: update
43
44 integer, intent(in) :: ng, tile, model, ifield
45 integer, intent(in) :: LBi, UBi, LBj, UBj
46
47#ifdef ASSUMED_SHAPE
48 real(r8), intent(in) :: Finp(LBi:,LBj:,:)
49 real(r8), intent(out) :: Fout(LBi:,LBj:)
50#else
51 real(r8), intent(in) :: Finp(LBi:UBi,LBj:UBj,2)
52 real(r8), intent(out) :: Fout(LBi:UBi,LBj:UBj)
53#endif
54!
55! Local variable declarations.
56!
57 logical :: LapplyBC, Lgrided, Lonerec
58
59 integer :: Tindex, gtype, i, it1, it2, j
60
61 real(dp) :: SecScale, fac, fac1, fac2
62 real(r8) :: Fval
63
64#include "set_bounds.h"
65!
66!----------------------------------------------------------------------
67! Set-up requested field for current tile.
68!----------------------------------------------------------------------
69!
70! Set switch to apply boundary conditions.
71!
72 IF (PRESENT(setbc)) THEN
73 lapplybc=setbc
74 ELSE
75 lapplybc=.true.
76 END IF
77!
78! Get requested field information from global storage.
79!
80 lgrided=linfo(1,ifield,ng)
81 lonerec=linfo(3,ifield,ng)
82 gtype =iinfo(1,ifield,ng)
83 tindex =iinfo(8,ifield,ng)
84 update=.true.
85!
86! Set linear, time interpolation factors. Fractional seconds are
87! rounded to the nearest milliseconds integer towards zero in the
88! time interpolation weights.
89!
90 secscale=1000.0_dp ! seconds to milliseconds
91 it1=3-tindex
92 it2=tindex
93 fac1=anint((tintrp(it2,ifield,ng)-time(ng))*secscale,dp)
94 fac2=anint((time(ng)-tintrp(it1,ifield,ng))*secscale,dp)
95!
96! Load time-invariant data. Time interpolation is not necessary.
97!
98 IF (lonerec) THEN
99 IF (lgrided) THEN
100 DO j=jstrr,jendr
101 DO i=istrr,iendr
102 fout(i,j)=finp(i,j,tindex)
103 END DO
104 END DO
105 ELSE
106 fval=fpoint(tindex,ifield,ng)
107 DO j=jstrr,jendr
108 DO i=istrr,iendr
109 fout(i,j)=fval
110 END DO
111 END DO
112 END IF
113!
114! Time-interpolate from gridded or point data.
115!
116 ELSE IF (((fac1*fac2).ge.0.0_dp).and. &
117 & ((fac1+fac2).gt.0.0_dp)) THEN
118 fac=1.0_dp/(fac1+fac2)
119 fac1=fac*fac1 ! nondimensional
120 fac2=fac*fac2 ! nondimensional
121 IF (lgrided) THEN
122 DO j=jstrr,jendr
123 DO i=istrr,iendr
124 fout(i,j)=fac1*finp(i,j,it1)+fac2*finp(i,j,it2)
125 END DO
126 END DO
127 ELSE
128 fval=fac1*fpoint(it1,ifield,ng)+fac2*fpoint(it2,ifield,ng)
129 DO j=jstrr,jendr
130 DO i=istrr,iendr
131 fout(i,j)=fval
132 END DO
133 END DO
134 END IF
135!
136! Activate synchronization flag if a new time record needs to be
137! read in at the next time step.
138!
139 IF ((time(ng)+dt(ng)).gt.tintrp(it2,ifield,ng)) THEN
140 IF (domain(ng)%SouthWest_Test(tile)) synchro_flag(ng)=.true.
141 END IF
142!
143! Unable to set-up requested field. Activate error flag to quit.
144!
145 ELSE
146 IF (domain(ng)%SouthWest_Test(tile)) THEN
147 IF (master) THEN
148 WRITE (stdout,10) trim(vname(1,ifield)), tdays(ng), &
149 & finfo(1,ifield,ng), finfo(2,ifield,ng), &
150 & finfo(3,ifield,ng), finfo(4,ifield,ng), &
151 & tintrp(it1,ifield,ng)*sec2day, &
152 & tintrp(it2,ifield,ng)*sec2day, &
153 & fac1*sec2day/secscale, &
154 & fac2*sec2day/secscale
155 END IF
156 10 FORMAT (/,' SET_2DFLD - current model time', &
157 & ' exceeds ending value for variable: ',a, &
158 & /,14x,'TDAYS = ',f15.4, &
159 & /,14x,'Data Tmin = ',f15.4,2x,'Data Tmax = ',f15.4, &
160 & /,14x,'Data Tstr = ',f15.4,2x,'Data Tend = ',f15.4, &
161 & /,14x,'TINTRP1 = ',f15.4,2x,'TINTRP2 = ',f15.4, &
162 & /,14x,'FAC1 = ',f15.4,2x,'FAC2 = ',f15.4)
163 exit_flag=2
164 update=.false.
165 END IF
166 END IF
167!
168! Exchange boundary data.
169!
170 IF (update) THEN
171 IF (lapplybc.and.(ewperiodic(ng).or.nsperiodic(ng))) THEN
172 IF (gtype.eq.r2dvar) THEN
173 CALL exchange_r2d_tile (ng, tile, &
174 & lbi, ubi, lbj, ubj, &
175 & fout)
176 ELSE IF (gtype.eq.u2dvar) THEN
177 CALL exchange_u2d_tile (ng, tile, &
178 & lbi, ubi, lbj, ubj, &
179 & fout)
180 ELSE IF (gtype.eq.v2dvar) THEN
181 CALL exchange_v2d_tile (ng, tile, &
182 & lbi, ubi, lbj, ubj, &
183 & fout)
184 END IF
185 END IF
186
187#ifdef DISTRIBUTE
188 IF (.not.lapplybc) THEN
189 CALL mp_exchange2d (ng, tile, model, 1, &
190 & lbi, ubi, lbj, ubj, &
191 & nghostpoints, &
192 & .false., .false., &
193 & fout)
194 ELSE
195 CALL mp_exchange2d (ng, tile, model, 1, &
196 & lbi, ubi, lbj, ubj, &
197 & nghostpoints, &
198 & ewperiodic(ng), nsperiodic(ng), &
199 & fout)
200 END IF
201#endif
202 END IF
203
204 RETURN
205 END SUBROUTINE set_2dfld_tile
206
207 END MODULE set_2dfld_mod
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
integer stdout
logical, dimension(:,:,:), allocatable linfo
real(dp), dimension(:,:,:), allocatable tintrp
real(dp), dimension(:,:,:), allocatable fpoint
real(dp), dimension(:,:,:), allocatable finfo
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:,:,:), allocatable iinfo
logical master
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, parameter u2dvar
Definition mod_param.F:718
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable synchro_flag
real(dp), dimension(:), allocatable tdays
real(dp), parameter sec2day
integer exit_flag
real(dp), dimension(:), allocatable time
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine set_2dfld_tile(ng, tile, model, ifield, lbi, ubi, lbj, ubj, finp, fout, update, setbc)
Definition set_2dfld.F:25