ROMS
Loading...
Searching...
No Matches
interp_floats.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined NONLINEAR && defined FLOATS
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Mark Hadfield !
8! Licensed under a MIT/X style license John M. Klinck !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine interpolates requested field at the float trajectory !
13! locations. !
14! !
15! On Input: !
16! !
17! ng Nested grid number. !
18! LBi I-dimension Lower bound. !
19! UBi I-dimension Upper bound. !
20! LBj J-dimension Lower bound. !
21! UBj J-dimension Upper bound. !
22! LBk K-dimension Lower bound. !
23! UBk K-dimension Upper bound. !
24! Lstr Starting float index to process. !
25! Lend Ending float index to process. !
26! itime Floats time level to process. !
27! ifield ID of field to compute. !
28! isBval Lateral boundary variable index. !
29! gtype Grid type. If negative, interpolate floats slopes. !
30! maskit Should the field be masked? Ignored if Land/Sea !
31! masking is not active. !
32! Fspval Special values for unbounded float. !
33! nudg Vertical random walk term to be added to the field. !
34! pm Inverse grid spacing (1/m) in the XI-direction. !
35! pn Inverse grid spacing (1/m) in the ETA-direction. !
36! Hz Vertical thicknesses (m). !
37! Amask Field Land/Sea mask. !
38! A Field to interpolate from. !
39! my_thread Float parallel thread bounded switch. !
40! bounded Float grid bounded status switch. !
41! !
42! On Output: !
43! !
44! track Interpolated field: track(ifield,itime,:). !
45! !
46!=======================================================================
47!
48 implicit none
49!
50 PRIVATE
51 PUBLIC :: interp_floats
52!
53 CONTAINS
54!
55!***********************************************************************
56 SUBROUTINE interp_floats (ng, LBi, UBi, LBj, UBj, LBk, UBk, &
57 & Lstr, Lend, itime, ifield, isBval, &
58 & gtype, maskit, Fspval, nudg, &
59 & pm, pn, &
60# ifdef SOLVE3D
61 & Hz, &
62# endif
63# ifdef MASKING
64 & Amask, &
65# endif
66 & A, my_thread, bounded, track)
67!***********************************************************************
68!
69 USE mod_param
70 USE mod_ncparam
71 USE mod_floats
72 USE mod_scalars
73!
74 implicit none
75!
76! Imported variable declarations.
77!
78 integer, intent(in) :: ng, lbi, ubi, lbj, ubj, lbk, ubk
79 integer, intent(in) :: lstr, lend, itime, ifield, isbval, gtype
80
81 logical, intent(in) :: maskit
82 logical, intent(in) :: my_thread(lstr:lend)
83 logical, intent(in) :: bounded(nfloats(ng))
84
85 real(dp), intent(in) :: fspval
86
87 real(r8), intent(in) :: nudg(lstr:lend)
88
89 real(r8), intent(in) :: pm(lbi:ubi,lbj:ubj)
90 real(r8), intent(in) :: pn(lbi:ubi,lbj:ubj)
91# ifdef SOLVE3D
92 real(r8), intent(in) :: hz(lbi:ubi,lbj:ubj,ubk)
93# endif
94# ifdef MASKING
95 real(r8), intent(in) :: amask(lbi:ubi,lbj:ubj)
96# endif
97 real(r8), intent(in) :: a(lbi:ubi,lbj:ubj,lbk:ubk)
98
99 real(r8), intent(inout) :: track(nfv(ng),0:nft,nfloats(ng))
100!
101! Local variable declarations.
102!
103 logical :: irvar, iuvar, jrvar, jvvar, krvar, kwvar, lmask
104 logical :: halo
105
106 integer :: ir, iu, jr, jv, kr, kw, l
107 integer :: i1, i2, j1, j2, k1, k2, khm, khp, vtype
108
109 real(r8) :: p1, p2, q1, q2, r1, r2
110 real(r8) :: s111, s211, s121, s221, s112, s212, s122, s222
111 real(r8) :: t111, t211, t121, t221, t112, t212, t122, t222
112
113# ifdef MASKING
114 integer :: irn, irnm1, irnp1, jrn, jrnm1, jrnp1
115
116 real(r8) :: cff1, cff2, cff3
117# endif
118!
119!-----------------------------------------------------------------------
120! Initialize various internal variables.
121!-----------------------------------------------------------------------
122!
123! Determine variable type switches.
124!
125 vtype=abs(gtype)
126 irvar=(vtype.eq.r2dvar).or.(vtype.eq.r3dvar).or. &
127 & (vtype.eq.v2dvar).or.(vtype.eq.v3dvar).or. &
128 & (vtype.eq.w3dvar)
129 jrvar=(vtype.eq.r2dvar).or.(vtype.eq.r3dvar).or. &
130 & (vtype.eq.u2dvar).or.(vtype.eq.u3dvar).or. &
131 & (vtype.eq.w3dvar)
132 iuvar=(vtype.eq.u2dvar).or.(vtype.eq.u3dvar)
133 jvvar=(vtype.eq.v2dvar).or.(vtype.eq.v3dvar)
134 krvar=(vtype.eq.r3dvar).or.(vtype.eq.u3dvar).or. &
135 & (vtype.eq.v3dvar)
136 kwvar=(vtype.eq.w3dvar)
137!
138! Determine whether to allow for masking in horizontal interpolation.
139!
140# ifdef MASKING
141 lmask=maskit
142# else
143 lmask=.false.
144# endif
145!
146! If not interpolating slope, set multipliers to 1.
147!
148 IF (gtype.ge.0) THEN
149 s111=1.0_r8
150 s121=1.0_r8
151 s211=1.0_r8
152 s221=1.0_r8
153 s112=1.0_r8
154 s122=1.0_r8
155 s212=1.0_r8
156 s222=1.0_r8
157 t111=1.0_r8
158 t121=1.0_r8
159 t211=1.0_r8
160 t221=1.0_r8
161 t112=1.0_r8
162 t122=1.0_r8
163 t212=1.0_r8
164 t222=1.0_r8
165 END IF
166!
167!-----------------------------------------------------------------------
168! Loop through floats.
169!-----------------------------------------------------------------------
170!
171 DO l=lstr,lend
172 IF (my_thread(l)) THEN
173 IF (.not.bounded(l)) THEN
174 track(ifield,itime,l)=fspval
175 ELSE
176!
177! Calculate indices and weights for vertical interpolation, if any.
178!
179 IF (krvar) THEN
180 kr=int(track(izgrd,itime,l)+0.5_r8)
181 k1=min(max(kr ,1),n(ng))
182 k2=min(max(kr+1,1),n(ng))
183 r2=real(k2-k1,r8)*(track(izgrd,itime,l)+ &
184 & 0.5_r8-real(k1,r8))
185 ELSE IF (kwvar) THEN
186 kw=int(track(izgrd,itime,l))
187 k1=min(max(kw ,0),n(ng))
188 k2=min(max(kw+1,0),n(ng))
189 r2=real(k2-k1,r8)*(track(izgrd,itime,l)-real(k1,r8))
190 ELSE
191 k1=1
192 k2=1
193 r2=0.0_r8
194 END IF
195 r1=1.0_r8-r2
196!
197!-----------------------------------------------------------------------
198! Interpolation on RHO-columns.
199!-----------------------------------------------------------------------
200!
201 IF (irvar.and.jrvar) THEN
202!
203 ir=int(track(ixgrd,itime,l))
204 jr=int(track(iygrd,itime,l))
205!
206 i1=min(max(ir ,0),lm(ng)+1)
207 i2=min(max(ir+1,1),lm(ng)+1)
208 j1=min(max(jr ,0),mm(ng)+1)
209 j2=min(max(jr+1,1),mm(ng)+1)
210!
211 p2=real(i2-i1,r8)*(track(ixgrd,itime,l)-real(i1,r8))
212 q2=real(j2-j1,r8)*(track(iygrd,itime,l)-real(j1,r8))
213 p1=1.0_r8-p2
214 q1=1.0_r8-q2
215# ifdef SOLVE3D
216!
217 IF (gtype.eq.-w3dvar) THEN
218 khm=min(max(k1 ,1),n(ng))
219 khp=min(max(k1+1,1),n(ng))
220 s111=2.0_r8*pm(i1,j1)*pn(i1,j1)/ &
221 & (hz(i1,j1,khm)+hz(i1,j1,khp))
222 s211=2.0_r8*pm(i2,j1)*pn(i2,j1)/ &
223 & (hz(i2,j1,khm)+hz(i2,j1,khp))
224 s121=2.0_r8*pm(i1,j2)*pn(i1,j2)/ &
225 & (hz(i1,j2,khm)+hz(i1,j2,khp))
226 s221=2.0_r8*pm(i2,j2)*pn(i2,j2)/ &
227 & (hz(i2,j2,khm)+hz(i2,j2,khp))
228 t111=2.0_r8/(hz(i1,j1,khm)+hz(i1,j1,khp))
229 t211=2.0_r8/(hz(i2,j1,khm)+hz(i2,j1,khp))
230 t121=2.0_r8/(hz(i1,j2,khm)+hz(i1,j2,khp))
231 t221=2.0_r8/(hz(i2,j2,khm)+hz(i2,j2,khp))
232 khm=min(max(k2 ,1),n(ng))
233 khp=min(max(k2+1,1),n(ng))
234 s112=2.0_r8*pm(i1,j1)*pn(i1,j1)/ &
235 & (hz(i1,j1,khm)+hz(i1,j1,khp))
236 s212=2.0_r8*pm(i2,j1)*pn(i2,j1)/ &
237 & (hz(i2,j1,khm)+hz(i2,j1,khp))
238 s122=2.0_r8*pm(i1,j2)*pn(i1,j2)/ &
239 & (hz(i1,j2,khm)+hz(i1,j2,khp))
240 s222=2.0_r8*pm(i2,j2)*pn(i2,j2)/ &
241 & (hz(i2,j2,khm)+hz(i2,j2,khp))
242 t112=2.0_r8/(hz(i1,j1,khm)+hz(i1,j1,khp))
243 t212=2.0_r8/(hz(i2,j1,khm)+hz(i2,j1,khp))
244 t122=2.0_r8/(hz(i1,j2,khm)+hz(i1,j2,khp))
245 t222=2.0_r8/(hz(i2,j2,khm)+hz(i2,j2,khp))
246 END IF
247# endif
248!
249 IF (lmask) THEN
250# ifdef MASKING
251 cff1=p1*q1*r1*amask(i1,j1)+ &
252 & p2*q1*r1*amask(i2,j1)+ &
253 & p1*q2*r1*amask(i1,j2)+ &
254 & p2*q2*r1*amask(i2,j2)+ &
255 & p1*q1*r2*amask(i1,j1)+ &
256 & p2*q1*r2*amask(i2,j1)+ &
257 & p1*q2*r2*amask(i1,j2)+ &
258 & p2*q2*r2*amask(i2,j2)
259 IF (cff1.gt.0.0_r8) THEN
260 cff2=p1*q1*r1*amask(i1,j1)*s111*a(i1,j1,k1)+ &
261 & p2*q1*r1*amask(i2,j1)*s211*a(i2,j1,k1)+ &
262 & p1*q2*r1*amask(i1,j2)*s121*a(i1,j2,k1)+ &
263 & p2*q2*r1*amask(i2,j2)*s221*a(i2,j2,k1)+ &
264 & p1*q1*r2*amask(i1,j1)*s112*a(i1,j1,k2)+ &
265 & p2*q1*r2*amask(i2,j1)*s212*a(i2,j1,k2)+ &
266 & p1*q2*r2*amask(i1,j2)*s122*a(i1,j2,k2)+ &
267 & p2*q2*r2*amask(i2,j2)*s222*a(i2,j2,k2)
268 cff3=(p1*q1*r1*amask(i1,j1)*t111+ &
269 & p2*q1*r1*amask(i2,j1)*t211+ &
270 & p1*q2*r1*amask(i1,j2)*t121+ &
271 & p2*q2*r1*amask(i2,j2)*t221+ &
272 & p1*q1*r2*amask(i1,j1)*t112+ &
273 & p2*q1*r2*amask(i2,j1)*t212+ &
274 & p1*q2*r2*amask(i1,j2)*t122+ &
275 & p2*q2*r2*amask(i2,j2)*t222)*nudg(l)
276 track(ifield,itime,l)=cff2/cff1+cff3
277 ELSE
278 track(ifield,itime,l)=0.0_r8
279 END IF
280# endif
281 ELSE
282 track(ifield,itime,l)=p1*q1*r1*s111*a(i1,j1,k1)+ &
283 & p2*q1*r1*s211*a(i2,j1,k1)+ &
284 & p1*q2*r1*s121*a(i1,j2,k1)+ &
285 & p2*q2*r1*s221*a(i2,j2,k1)+ &
286 & p1*q1*r2*s112*a(i1,j1,k2)+ &
287 & p2*q1*r2*s212*a(i2,j1,k2)+ &
288 & p1*q2*r2*s122*a(i1,j2,k2)+ &
289 & p2*q2*r2*s222*a(i2,j2,k2)+ &
290 & (p1*q1*r1*t111+ &
291 & p2*q1*r1*t211+ &
292 & p1*q2*r1*t121+ &
293 & p2*q2*r1*t221+ &
294 & p1*q1*r2*t112+ &
295 & p2*q1*r2*t212+ &
296 & p1*q2*r2*t122+ &
297 & p2*q2*r2*t222)*nudg(l)
298 END IF
299!
300!-----------------------------------------------------------------------
301! Interpolation at horizontal velocity points.
302!-----------------------------------------------------------------------
303!
304 ELSE
305 ir=int(track(ixgrd,itime,l))
306 jr=int(track(iygrd,itime,l))
307 iu=int(track(ixgrd,itime,l)+0.5_r8)
308 jv=int(track(iygrd,itime,l)+0.5_r8)
309!
310 halo=.false.
311# ifdef MASKING
312!
313! Is the float inside the halo (i.e. inside a masked grid cell or
314! within 0.5*delta of one)? Depending on which portion of the cell
315! the float is in, we look for adjacent and diagonally adjacent
316! masked cells.
317!
318! Note that the halo-checking code is evaluated every time
319! interp_floats is called for a velocity variable with masking
320! activated. This is at least twice as often as it needs to be
321! called, as the float does not move between interpolating U and V.
322!
323! The special-case code for periodic boundary conditions may be
324! unnecessary
325!
326 IF (lmask) THEN
327 irn=nint(track(ixgrd,itime,l))
328 jrn=nint(track(iygrd,itime,l))
329 IF (ewperiodic(ng)) THEN
330 IF (irn.ge.lm(ng)) THEN
331 irnp1=irn+1-lm(ng)
332 ELSE
333 irnp1=irn+1
334 END IF
335 IF (irn.le.1) THEN
336 irnm1=irn-1+lm(ng)
337 ELSE
338 irnm1=irn-1
339 END IF
340 ELSE
341 irnm1=irn-1
342 irnp1=irn+1
343 END IF
344 IF (nsperiodic(ng)) THEN
345 IF (jrn.ge.mm(ng)) THEN
346 jrnp1=jrn+1-mm(ng)
347 ELSE
348 jrnp1=jrn+1
349 END IF
350 IF (jrn.le.1) THEN
351 jrnm1=jrn-1+mm(ng)
352 ELSE
353 jrnm1=jrn-1
354 END IF
355 ELSE
356 jrnm1=jrn-1
357 jrnp1=jrn+1
358 END IF
359 IF (amask(irn,jrn).lt.0.5_r8) THEN
360 halo=.true.
361 ELSE IF ((ir.lt.irn).and. &
362 & (amask(irn-1,jrn).lt.0.5_r8)) THEN
363 halo=.true.
364 ELSE IF ((ir.eq.irn).and. &
365 & (amask(irn+1,jrn).lt.0.5_r8)) THEN
366 halo=.true.
367 ELSE IF ((jr.lt.jrn).and. &
368 & (amask(irn,jrn-1).lt.0.5_r8)) THEN
369 halo=.true.
370 ELSE IF ((jr.eq.jrn).and. &
371 & (amask(irn,jrn+1).lt.0.5_r8)) THEN
372 halo=.true.
373 ELSE IF ((ir.lt.irn).and.(jr.lt.jrn).and. &
374 & (amask(irn-1,jrn-1).lt.0.5_r8)) THEN
375 halo=.true.
376 ELSE IF ((ir.eq.irn).and.(jr.lt.jrn).and. &
377 & (amask(irn+1,jrn-1).lt.0.5_r8)) THEN
378 halo=.true.
379 ELSE IF ((ir.lt.irn).and.(jr.eq.jrn).and. &
380 & (amask(irn-1,jrn+1).lt.0.5_r8)) THEN
381 halo=.true.
382 ELSE IF ((ir.eq.irn).and.(jr.eq.jrn).and. &
383 & (amask(irn+1,jrn+1).lt.0.5_r8)) THEN
384 halo=.true.
385 END IF
386 END IF
387# endif
388!
389!-----------------------------------------------------------------------
390! Interpolation at U-points.
391!-----------------------------------------------------------------------
392!
393 IF (iuvar) THEN
394 IF (halo) THEN
395# ifdef MASKING
396!
397! Velocity interpolation inside the halo is linear in the parallel
398! direction and nearest-neighbour in the perpendicular direction.
399! The latter ensures that the perpendicular velocity is zero
400! everywhere on the perimeter of a masked cell.
401!
402 i1=min(max(iu ,1),lm(ng)+1)
403 i2=min(max(iu+1,1),lm(ng)+1)
404 j1=jrn
405!
406 p2=real(i2-i1,r8)* &
407 & (track(ixgrd,itime,l)-real(i1,r8)+0.5_r8)
408 p1=1.0_r8-p2
409 q1=1.0_r8
410!
411 IF (gtype.lt.0) THEN
412 s111=0.5_r8*(pm(i1-1,j1)+pm(i1,j1))
413 s211=0.5_r8*(pm(i2-1,j1)+pm(i2,j1))
414 s112=s111
415 s212=s112
416 END IF
417!
418 track(ifield,itime,l)=p1*q1*r1*s111*a(i1,j1,k1)+ &
419 & p2*q1*r1*s211*a(i2,j1,k1)+ &
420 & p1*q1*r2*s112*a(i1,j1,k2)+ &
421 & p2*q1*r2*s212*a(i2,j1,k2)+ &
422 & nudg(l)
423# endif
424 ELSE
425!
426! Bilinear interpolation outside halo.
427!
428 i1=min(max(iu ,1),lm(ng)+1)
429 i2=min(max(iu+1,1),lm(ng)+1)
430 j1=min(max(jr ,0),mm(ng)+1)
431 j2=min(max(jr+1,0),mm(ng)+1)
432!
433 p2=real(i2-i1,r8)* &
434 & (track(ixgrd,itime,l)-real(i1,r8)+0.5_r8)
435 q2=real(j2-j1,r8)* &
436 & (track(iygrd,itime,l)-real(j1,r8))
437 p1=1.0_r8-p2
438 q1=1.0_r8-q2
439!
440 IF (gtype.lt.0) THEN
441 s111=0.5_r8*(pm(i1-1,j1)+pm(i1,j1))
442 s211=0.5_r8*(pm(i2-1,j1)+pm(i2,j1))
443 s121=0.5_r8*(pm(i1-1,j2)+pm(i1,j2))
444 s221=0.5_r8*(pm(i2-1,j2)+pm(i2,j2))
445 s112=s111
446 s212=s112
447 s122=s121
448 s222=s221
449 END IF
450!
451 track(ifield,itime,l)=p1*q1*r1*s111*a(i1,j1,k1)+ &
452 & p2*q1*r1*s211*a(i2,j1,k1)+ &
453 & p1*q2*r1*s121*a(i1,j2,k1)+ &
454 & p2*q2*r1*s221*a(i2,j2,k1)+ &
455 & p1*q1*r2*s112*a(i1,j1,k2)+ &
456 & p2*q1*r2*s212*a(i2,j1,k2)+ &
457 & p1*q2*r2*s122*a(i1,j2,k2)+ &
458 & p2*q2*r2*s222*a(i2,j2,k2)+ &
459 & nudg(l)
460 END IF
461!
462!-----------------------------------------------------------------------
463! Interpolation at V-points.
464!-----------------------------------------------------------------------
465!
466 ELSE IF (jvvar) THEN
467 IF (halo) THEN
468# ifdef MASKING
469!
470! Velocity interpolation inside the halo is linear in the parallel
471! direction and nearest-neighbour in the perpendicular direction.
472! The latter ensures that the perpendicular velocity is zero
473! everywhere on the perimeter of a masked cell.
474!
475 i1=irn
476 j1=min(max(jv ,1),mm(ng)+1)
477 j2=min(max(jv+1,1),mm(ng)+1)
478!
479 q2=real(j2-j1,r8)* &
480 & (track(iygrd,itime,l)-real(j1,r8)+0.5_r8)
481 p1=1.0_r8
482 q1=1.0_r8-q2
483!
484 IF (gtype.lt.0) THEN
485 s111=0.5_r8*(pn(i1,j1-1)+pn(i1,j1))
486 s121=0.5_r8*(pn(i1,j2-1)+pn(i1,j2))
487 s112=s111
488 s122=s121
489 END IF
490!
491 track(ifield,itime,l)=p1*q1*r1*s111*a(i1,j1,k1)+ &
492 & p1*q2*r1*s121*a(i1,j2,k1)+ &
493 & p1*q1*r2*s112*a(i1,j1,k2)+ &
494 & p1*q2*r2*s122*a(i1,j2,k2)+ &
495 & nudg(l)
496# endif
497 ELSE
498!
499! Bilinear interpolation outside halo.
500!
501 i1=min(max(ir ,0),lm(ng)+1)
502 i2=min(max(ir+1,1),lm(ng)+1)
503 j1=min(max(jv ,1),mm(ng)+1)
504 j2=min(max(jv+1,1),mm(ng)+1)
505!
506 p2=real(i2-i1,r8)* &
507 & (track(ixgrd,itime,l)-real(i1,r8))
508 q2=real(j2-j1,r8)* &
509 & (track(iygrd,itime,l)-real(j1,r8)+0.5_r8)
510 p1=1.0_r8-p2
511 q1=1.0_r8-q2
512!
513 IF (gtype.lt.0) THEN
514 s111=0.5_r8*(pn(i1,j1-1)+pn(i1,j1))
515 s211=0.5_r8*(pn(i2,j1-1)+pn(i2,j1))
516 s121=0.5_r8*(pn(i1,j2-1)+pn(i1,j2))
517 s221=0.5_r8*(pn(i2,j2-1)+pn(i2,j2))
518 s112=s111
519 s212=s112
520 s122=s121
521 s222=s221
522 END IF
523!
524 track(ifield,itime,l)=p1*q1*r1*s111*a(i1,j1,k1)+ &
525 & p2*q1*r1*s211*a(i2,j1,k1)+ &
526 & p1*q2*r1*s121*a(i1,j2,k1)+ &
527 & p2*q2*r1*s221*a(i2,j2,k1)+ &
528 & p1*q1*r2*s112*a(i1,j1,k2)+ &
529 & p2*q1*r2*s212*a(i2,j1,k2)+ &
530 & p1*q2*r2*s122*a(i1,j2,k2)+ &
531 & p2*q2*r2*s222*a(i2,j2,k2)+ &
532 & nudg(l)
533 END IF
534 END IF
535 END IF
536 END IF
537 END IF
538 END DO
539
540 RETURN
541 END SUBROUTINE interp_floats
542#endif
543 END MODULE interp_floats_mod
subroutine, public interp_floats(ng, lbi, ubi, lbj, ubj, lbk, ubk, lstr, lend, itime, ifield, isbval, gtype, maskit, fspval, nudg, pm, pn, hz, amask, a, my_thread, bounded, track)
integer, parameter iygrd
Definition mod_floats.F:82
integer, parameter ixgrd
Definition mod_floats.F:81
integer, parameter izgrd
Definition mod_floats.F:83
integer, dimension(:), allocatable nfv
Definition mod_param.F:547
integer, dimension(:), allocatable nfloats
Definition mod_param.F:543
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, parameter r3dvar
Definition mod_param.F:721
integer, parameter u3dvar
Definition mod_param.F:722
integer, parameter nft
Definition mod_param.F:539
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, parameter u2dvar
Definition mod_param.F:718
integer, parameter w3dvar
Definition mod_param.F:724
integer, dimension(:), allocatable mm
Definition mod_param.F:456
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
integer, parameter v3dvar
Definition mod_param.F:723
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic