ROMS
Loading...
Searching...
No Matches
regrid.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE regrid_mod
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 interpolates gridded data, Finp, to model locations !
12! Xout and Yout. Input data are allow to have longitudes ranging !
13! from [-180 180] or [0 360]. Usually, the ROMS application have !
14! longitudes in the range [-180 180]. !
15! !
16! The Imin, Imax, Jmin, and Jmax indices are the global values of !
17! application domain in serial I/O. However, in parallel I/O they !
18! are the tiled values. !
19! !
20! On Input: !
21! !
22! ng Nested grid number (integer) !
23! model Calling model identifier (integer) !
24! ncname NetCDF file name (string) !
25! ncid NetCDF file ID (integer) !
26#if defined PIO_LIB && defined DISTRIBUTE
27!or pioFile PIO file descriptor structure, TYPE(file_desc_t) !
28! pioFile%fh file handler !
29! pioFile%iosystem IO system descriptor (struct) !
30#endif
31! ncvname NetCDF variable name (string) !
32! ncvarid NetCDF variable ID (integer) !
33#if defined PIO_LIB && defined DISTRIBUTE
34!or pioDesc PIO variable descriptor structure, TYPE(My_VarDesc) !
35! pioVar%vd variable descriptor TYPE(Var_Desc_t) !
36! pioVar%dkind variable data kind !
37! pioVar%gtype variable C-gridtype !
38#endif
39! gtype C-grid type (integer) !
40! iflag Interpolation flag (integer, 0: linear, 1: cubic) !
41! Nx X-dimension size for gridded data, Finp (integer) !
42! Ny Y-dimension size for gridded data, Finp (integer) !
43! Finp Gridded data to interpolate from (real) !
44! Amin Gridded data minimum value (integer) !
45! Amax Gridded data maximum value (integer) !
46! LBi Fout I-dimension Lower bound (integer) !
47! UBi Fout I-dimension Upper bound (integer) !
48! LBj Fout J-dimension Lower bound (integer) !
49! UBj Fout J-dimension Upper bound (integer) !
50! Imin Fout starting data I-index (integer) !
51! Imax Fout ending data I-index (integer) !
52! Jmin Fout starting data J-index (integer) !
53! Jmax Fout ending data J-index (integer) !
54# ifdef MASKING
55! MyMask Fout Land/Sea mask (real) !
56# endif
57! MyXout Work X-locations (longitude) for regridding (real) !
58! Xout X-locations (longitude) to interpolate (real) !
59! Yout Y-locations (latitude) to interpolate (real) !
60! !
61! On Output: !
62! !
63! Fout Interpolated field (real) !
64! !
65!=======================================================================
66!
67 USE mod_param
68 USE mod_parallel
69 USE mod_iounits
70 USE mod_ncparam
71 USE mod_scalars
72 USE roms_interpolate_mod
73!
74#ifdef DISTRIBUTE
75 USE distribute_mod, ONLY : mp_reduce
76#endif
78#if defined DISTRIBUTE && defined REGRID_SHAPIRO
80#endif
81#ifdef REGRID_SHAPIRO
82 USE shapiro_mod, ONLY : shapiro2d_tile
83#endif
84 USE strings_mod, ONLY : founderror
85!
86 implicit none
87!
88! Interface for same name routine overloading.
89!
90 PUBLIC :: regrid_nf90
91#if defined PIO_LIB && defined DISTRIBUTE
92 PUBLIC :: regrid_pio
93#endif
94!
95 CONTAINS
96!
97!***********************************************************************
98 SUBROUTINE regrid_nf90 (ng, model, ncname, ncid, &
99 & ncvname, ncvarid, gtype, iflag, &
100 & Nx, Ny, Finp, Amin, Amax, &
101 & LBi, UBi, LBj, UBj, &
102 & Imin, Imax, Jmin, Jmax, &
103#ifdef MASKING
104 & MyMask, &
105#endif
106 & MyXout, Xout, Yout, Fout)
107!***********************************************************************
108!
109! Imported variable declarations.
110!
111 integer, intent(in) :: ng, model, ncid, ncvarid, gtype, iflag
112 integer, intent(in) :: nx, ny
113 integer, intent(in) :: lbi, ubi, lbj, ubj
114 integer, intent(in) :: imin, imax, jmin, jmax
115!
116 real(r8), intent(inout) :: amin, amax
117!
118 real(r8), intent(inout) :: finp(nx,ny)
119
120 real(r8), intent(inout) :: myxout(lbi:ubi,lbj:ubj)
121
122#ifdef MASKING
123 real(r8), intent(in) :: mymask(lbi:ubi,lbj:ubj)
124#endif
125 real(r8), intent(in) :: xout(lbi:ubi,lbj:ubj)
126 real(r8), intent(in) :: yout(lbi:ubi,lbj:ubj)
127
128 real(r8), intent(out) :: fout(lbi:ubi,lbj:ubj)
129!
130 character (len=*), intent(in) :: ncname
131 character (len=*), intent(in) :: ncvname
132!
133! Local variable declarations
134!
135 logical :: eastlon, rectangular
136!
137 integer :: i, j
138 integer :: istr, iend, jstr, jend
139#ifdef REGRID_SHAPIRO
140 integer :: imins, imaxs, jmins, jmaxs, tile
141#endif
142#ifdef DISTRIBUTE
143 integer :: cgrid, ghost
144#endif
145!
146 real(r8), parameter :: ijspv = 0.0_r8
147
148 real(r8) :: my_min, my_max, xmin, xmax, ymin, ymax
149 real(r8) :: mylonmin, mylonmax
150
151 real(r8), dimension(Nx,Ny) :: angle
152 real(r8), dimension(Nx,Ny) :: xinp
153 real(r8), dimension(Nx,Ny) :: yinp
154
155 real(r8), dimension(LBi:UBi,LBj:UBj) :: iout
156 real(r8), dimension(LBi:UBi,LBj:UBj) :: jout
157
158#ifdef DISTRIBUTE
159 real(r8), dimension(2) :: rbuffer
160!
161 character (len=3), dimension(2) :: op_handle
162#endif
163!
164 character (len=*), parameter :: myfile = &
165 & __FILE__//", regrid_nf90"
166!
167!-----------------------------------------------------------------------
168! Get input variable coordinates.
169!-----------------------------------------------------------------------
170!
171 CALL get_varcoords (ng, model, ncname, ncid, &
172 & ncvname, ncvarid, nx, ny, &
173 & xmin, xmax, xinp, ymin, ymax, yinp, &
174 & rectangular)
175 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
176!
177! Set input gridded data rotation angle.
178!
179 DO i=1,nx
180 DO j=1,ny
181 angle(i,j)=0.0_r8
182 END DO
183 END DO
184!
185! Initialize local fractional coordinates arrays to avoid
186! deframentation.
187!
188 iout=0.0_r8
189 jout=0.0_r8
190!
191! Copy longitude coordinate Xout to MyXout. If the longitude of the
192! data is from a global grid [0-360] or in degrees_east, convert Xout
193! to east longitudes (MyXout) to facilitate regridding. In such case,
194! positive multiples of 360 map to 360 and negative multiples of 360
195! map to zero using the MODULO intrinsic Fortran function.
196!
197 IF ((xmin.ge.0.0_r8).and.(xmax.gt.0.0_r8).and. &
198 & ((xmax-xmin).gt.315.0_r8)) THEN
199 eastlon=.true.
200 mylonmin=modulo(lonmin(ng), 360.0_r8)
201 IF ((mylonmin.eq.0.0_r8).and. &
202 & (lonmin(ng).gt.0.0_r8)) mylonmin=360.0_r8
203 mylonmax=modulo(lonmax(ng), 360.0_r8)
204 IF ((mylonmax.eq.0.0_r8).and. &
205 & (lonmax(ng).gt.0.0_r8)) mylonmax=360.0_r8
206 ELSE
207 eastlon=.false.
208 mylonmin=lonmin(ng)
209 mylonmax=lonmax(ng)
210 END IF
211 IF (eastlon) THEN
212 DO j=lbj,ubj
213 DO i=lbi,ubi
214 myxout(i,j)=modulo(xout(i,j), 360.0_r8) ! range [0 360]
215 IF ((myxout(i,j).eq.0.0_r8).and. &
216 & (xout(i,j).gt.0.0_r8)) myxout(i,j)=360.0_r8
217 END DO
218 END DO
219 ELSE
220 DO j=lbj,ubj
221 DO i=lbi,ubi
222 myxout(i,j)=xout(i,j) ! range [-180 180]
223 END DO
224 END DO
225 END IF
226!
227!-----------------------------------------------------------------------
228! Check if gridded data contains model grid.
229!-----------------------------------------------------------------------
230!
231 IF ((mylonmin .lt.xmin).or. &
232 & (mylonmax .gt.xmax).or. &
233 & (latmin(ng).lt.ymin).or. &
234 & (latmax(ng).gt.ymax)) THEN
235 IF (master) THEN
236 WRITE (stdout,10) xmin, xmax, ymin, ymax, &
237 & mylonmin , mylonmax, &
238 & latmin(ng), latmax(ng)
239 10 FORMAT (/, ' REGRID - input gridded data does not contain', &
240 & ' model grid:', /, &
241 & /,10x,'Gridded: LonMin = ',f9.4,' LonMax = ',f9.4, &
242 & /,10x,' LatMin = ',f9.4,' LatMax = ',f9.4, &
243 & /,10x,'Model: LonMin = ',f9.4,' LonMax = ',f9.4, &
244 & /,10x,' LatMin = ',f9.4,' LatMax = ',f9.4)
245 END IF
246 exit_flag=4
247 RETURN
248 END IF
249!
250!-----------------------------------------------------------------------
251! Interpolate (bilinear or bicubic) to requested positions.
252!-----------------------------------------------------------------------
253!
254! Set tile starting and ending indices.
255!
256#ifdef DISTRIBUTE
257 ghost=0 ! non-overlapping, no ghost points
258
259 SELECT CASE (abs(gtype))
260 CASE (p2dvar, p3dvar)
261 cgrid=1
262 CASE (r2dvar, r3dvar)
263 cgrid=2
264 CASE (u2dvar, u3dvar)
265 cgrid=3
266 CASE (v2dvar, v3dvar)
267 cgrid=4
268 CASE DEFAULT
269 cgrid=1
270 END SELECT
271
272 istr=bounds(ng)%Imin(cgrid,ghost,myrank)
273 iend=bounds(ng)%Imax(cgrid,ghost,myrank)
274 jstr=bounds(ng)%Jmin(cgrid,ghost,myrank)
275 jend=bounds(ng)%Jmax(cgrid,ghost,myrank)
276#else
277 istr=imin
278 iend=imax
279 jstr=jmin
280 jend=jmax
281#endif
282!
283! Find fractional indices (Iout,Jout) of the grid cells in Finp
284! containing positions to intepolate.
285!
286 CALL hindices (ng, 1, nx, 1, ny, 1, nx, 1, ny, &
287 & angle, xinp, yinp, &
288 & lbi, ubi, lbj, ubj, &
289 & istr, iend, jstr, jend, &
290 & myxout, yout, &
291 & iout, jout, &
292 & ijspv, rectangular)
293
294 IF (iflag.eq.linear) THEN
295 CALL linterp2d (ng, 1, nx, 1, ny, &
296 & xinp, yinp, finp, &
297 & lbi, ubi, lbj, ubj, &
298 & istr, iend, jstr, jend, &
299 & iout, jout, myxout, yout, &
300 & fout, my_min, my_max)
301 ELSE IF (iflag.eq.cubic) THEN
302 CALL cinterp2d (ng, 1, nx, 1, ny, &
303 & xinp, yinp, finp, &
304 & lbi, ubi, lbj, ubj, &
305 & istr, iend, jstr, jend, &
306 & iout, jout, myxout, yout, &
307 & fout, my_min, my_max)
308 END IF
309
310#ifdef REGRID_SHAPIRO
311!
312! Apply Shapiro filter to interpolated date. Often, input data is
313! coarse and field have small scale noise. So, it is desirable to
314! the 2-delta noise.
315!
316# ifdef DISTRIBUTE
317 tile=myrank
318 CALL mp_exchange2d (ng, tile, 1, 1, &
319 & lbi, ubi, lbj, ubj, &
320 & nghostpoints, &
321 & ewperiodic(ng), nsperiodic(ng), &
322 & fout)
323# else
324 tile=-1
325# endif
326# ifdef NESTING
327 imins=bounds(ng)%Istr(tile)-4
328 imaxs=bounds(ng)%Iend(tile)+3
329 jmins=bounds(ng)%Jstr(tile)-4
330 jmaxs=bounds(ng)%Jend(tile)+3
331# else
332 imins=bounds(ng)%Istr(tile)-3
333 imaxs=bounds(ng)%Iend(tile)+3
334 jmins=bounds(ng)%Jstr(tile)-3
335 jmaxs=bounds(ng)%Jend(tile)+3
336# endif
337!
338 CALL shapiro2d_tile (ng, tile, 1, &
339 & lbi, ubi, lbj, ubj, &
340 & imins, imaxs, jmins, jmaxs, &
341# ifdef MASKING
342 & mymask, &
343# endif
344 & fout)
345#endif
346!
347! Compute global interpolated field minimum and maximum values.
348! Notice that gridded data values are overwritten.
349!
350#ifdef DISTRIBUTE
351 rbuffer(1)=my_min
352 rbuffer(2)=my_max
353 op_handle(1)='MIN'
354 op_handle(2)='MAX'
355 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
356 amin=rbuffer(1)
357 amax=rbuffer(2)
358#else
359 amin=my_min
360 amax=my_max
361#endif
362!
363 RETURN
364 END SUBROUTINE regrid_nf90
365
366#if defined PIO_LIB && defined DISTRIBUTE
367!
368!***********************************************************************
369 SUBROUTINE regrid_pio (ng, model, ncname, pioFile, &
370 & ncvname, pioVar, gtype, iflag, &
371 & Nx, Ny, Finp, Amin, Amax, &
372 & LBi, UBi, LBj, UBj, &
373 & Imin, Imax, Jmin, Jmax, &
374# ifdef MASKING
375 & MyMask, &
376# endif
377 & MyXout, Xout, Yout, Fout)
378!***********************************************************************
379!
381!
382! Imported variable declarations.
383!
384 integer, intent(in) :: ng, model, gtype, iflag
385 integer, intent(in) :: nx, ny
386 integer, intent(in) :: lbi, ubi, lbj, ubj
387 integer, intent(in) :: imin, imax, jmin, jmax
388!
389 real(r8), intent(inout) :: amin, amax
390!
391 real(r8), intent(inout) :: finp(nx,ny)
392
393 real(r8), intent(inout) :: myxout(lbi:ubi,lbj:ubj)
394
395# ifdef MASKING
396 real(r8), intent(in) :: mymask(lbi:ubi,lbj:ubj)
397# endif
398 real(r8), intent(in) :: xout(lbi:ubi,lbj:ubj)
399 real(r8), intent(in) :: yout(lbi:ubi,lbj:ubj)
400
401 real(r8), intent(out) :: fout(lbi:ubi,lbj:ubj)
402!
403 character (len=*), intent(in) :: ncname
404 character (len=*), intent(in) :: ncvname
405!
406 TYPE (file_desc_t), intent(inout) :: piofile
407 TYPE (my_vardesc), intent(inout) :: piovar
408!
409! Local variable declarations
410!
411 logical :: eastlon, rectangular
412!
413 integer :: i, j
414 integer :: istr, iend, jstr, jend
415# ifdef REGRID_SHAPIRO
416 integer :: imins, imaxs, jmins, jmaxs, tile
417# endif
418# ifdef DISTRIBUTE
419 integer :: cgrid, ghost
420# endif
421!
422 real(r8), parameter :: ijspv = 0.0_r8
423
424 real(r8) :: my_min, my_max, xmin, xmax, ymin, ymax
425 real(r8) :: mylonmin, mylonmax
426
427 real(r8), dimension(Nx,Ny) :: angle
428 real(r8), dimension(Nx,Ny) :: xinp
429 real(r8), dimension(Nx,Ny) :: yinp
430
431 real(r8), dimension(LBi:UBi,LBj:UBj) :: iout
432 real(r8), dimension(LBi:UBi,LBj:UBj) :: jout
433
434# ifdef DISTRIBUTE
435 real(r8), dimension(2) :: rbuffer
436!
437 character (len=3), dimension(2) :: op_handle
438# endif
439!
440 character (len=*), parameter :: myfile = &
441 & __FILE__//", regrid_pio"
442!
443!-----------------------------------------------------------------------
444! Get input variable coordinates.
445!-----------------------------------------------------------------------
446!
447 CALL get_varcoords (ng, model, ncname, piofile, &
448 & ncvname, piovar, nx, ny, &
449 & xmin, xmax, xinp, ymin, ymax, yinp, &
450 & rectangular)
451 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
452!
453! Set input gridded data rotation angle.
454!
455 DO i=1,nx
456 DO j=1,ny
457 angle(i,j)=0.0_r8
458 END DO
459 END DO
460!
461! Initialize local fractional coordinates arrays to avoid
462! deframentation.
463!
464 iout=0.0_r8
465 jout=0.0_r8
466!
467! Copy longitude coordinate Xout to MyXout. If the longitude of the
468! data is from a global grid [0-360] or in degrees_east, convert Xout
469! to east longitudes (MyXout) to facilitate regridding. In such case,
470! positive multiples of 360 map to 360 and negative multiples of 360
471! map to zero using the MODULO intrinsic Fortran function.
472!
473 IF ((xmin.ge.0.0_r8).and.(xmax.gt.0.0_r8).and. &
474 & ((xmax-xmin).gt.315.0_r8)) THEN
475 eastlon=.true.
476 mylonmin=modulo(lonmin(ng), 360.0_r8)
477 IF ((mylonmin.eq.0.0_r8).and. &
478 & (lonmin(ng).gt.0.0_r8)) mylonmin=360.0_r8
479 mylonmax=modulo(lonmax(ng), 360.0_r8)
480 IF ((mylonmax.eq.0.0_r8).and. &
481 & (lonmax(ng).gt.0.0_r8)) mylonmax=360.0_r8
482 ELSE
483 eastlon=.false.
484 mylonmin=lonmin(ng)
485 mylonmax=lonmax(ng)
486 END IF
487 IF (eastlon) THEN
488 DO j=lbj,ubj
489 DO i=lbi,ubi
490 myxout(i,j)=modulo(xout(i,j), 360.0_r8) ! range [0 360]
491 IF ((myxout(i,j).eq.0.0_r8).and. &
492 & (xout(i,j).gt.0.0_r8)) myxout(i,j)=360.0_r8
493 END DO
494 END DO
495 ELSE
496 DO j=lbj,ubj
497 DO i=lbi,ubi
498 myxout(i,j)=xout(i,j) ! range [-180 180]
499 END DO
500 END DO
501 END IF
502!
503!-----------------------------------------------------------------------
504! Check if gridded data contains model grid.
505!-----------------------------------------------------------------------
506!
507 IF ((mylonmin .lt.xmin).or. &
508 & (mylonmax .gt.xmax).or. &
509 & (latmin(ng).lt.ymin).or. &
510 & (latmax(ng).gt.ymax)) THEN
511 IF (master) THEN
512 WRITE (stdout,10) xmin, xmax, ymin, ymax, &
513 & mylonmin , mylonmax, &
514 & latmin(ng), latmax(ng)
515 10 FORMAT (/, ' REGRID - input gridded data does not contain', &
516 & ' model grid:', /, &
517 & /,10x,'Gridded: LonMin = ',f9.4,' LonMax = ',f9.4, &
518 & /,10x,' LatMin = ',f9.4,' LatMax = ',f9.4, &
519 & /,10x,'Model: LonMin = ',f9.4,' LonMax = ',f9.4, &
520 & /,10x,' LatMin = ',f9.4,' LatMax = ',f9.4)
521 END IF
522 exit_flag=4
523 RETURN
524 END IF
525!
526!-----------------------------------------------------------------------
527! Interpolate (bilinear or bicubic) to requested positions.
528!-----------------------------------------------------------------------
529!
530! Set tile starting and ending indices.
531!
532# ifdef DISTRIBUTE
533 ghost=0 ! non-overlapping, no ghost points
534
535 SELECT CASE (abs(gtype))
536 CASE (p2dvar, p3dvar)
537 cgrid=1
538 CASE (r2dvar, r3dvar)
539 cgrid=2
540 CASE (u2dvar, u3dvar)
541 cgrid=3
542 CASE (v2dvar, v3dvar)
543 cgrid=4
544 CASE DEFAULT
545 cgrid=1
546 END SELECT
547
548 istr=bounds(ng)%Imin(cgrid,ghost,myrank)
549 iend=bounds(ng)%Imax(cgrid,ghost,myrank)
550 jstr=bounds(ng)%Jmin(cgrid,ghost,myrank)
551 jend=bounds(ng)%Jmax(cgrid,ghost,myrank)
552# else
553 istr=imin
554 iend=imax
555 jstr=jmin
556 jend=jmax
557# endif
558!
559! Find fractional indices (Iout,Jout) of the grid cells in Finp
560! containing positions to intepolate.
561!
562 CALL hindices (ng, 1, nx, 1, ny, 1, nx, 1, ny, &
563 & angle, xinp, yinp, &
564 & lbi, ubi, lbj, ubj, &
565 & istr, iend, jstr, jend, &
566 & myxout, yout, &
567 & iout, jout, &
568 & ijspv, rectangular)
569
570 IF (iflag.eq.linear) THEN
571 CALL linterp2d (ng, 1, nx, 1, ny, &
572 & xinp, yinp, finp, &
573 & lbi, ubi, lbj, ubj, &
574 & istr, iend, jstr, jend, &
575 & iout, jout, myxout, yout, &
576 & fout, my_min, my_max)
577 ELSE IF (iflag.eq.cubic) THEN
578 CALL cinterp2d (ng, 1, nx, 1, ny, &
579 & xinp, yinp, finp, &
580 & lbi, ubi, lbj, ubj, &
581 & istr, iend, jstr, jend, &
582 & iout, jout, myxout, yout, &
583 & fout, my_min, my_max)
584 END IF
585
586# ifdef REGRID_SHAPIRO
587!
588! Apply Shapiro filter to interpolated date. Often, input data is
589! coarse and field have small scale noise. So, it is desirable to
590! the 2-delta noise.
591!
592# ifdef DISTRIBUTE
593 tile=myrank
594 CALL mp_exchange2d (ng, tile, 1, 1, &
595 & lbi, ubi, lbj, ubj, &
596 & nghostpoints, &
597 & ewperiodic(ng), nsperiodic(ng), &
598 & fout)
599# else
600 tile=-1
601# endif
602# ifdef NESTING
603 imins=bounds(ng)%Istr(tile)-4
604 imaxs=bounds(ng)%Iend(tile)+3
605 jmins=bounds(ng)%Jstr(tile)-4
606 jmaxs=bounds(ng)%Jend(tile)+3
607# else
608 imins=bounds(ng)%Istr(tile)-3
609 imaxs=bounds(ng)%Iend(tile)+3
610 jmins=bounds(ng)%Jstr(tile)-3
611 jmaxs=bounds(ng)%Jend(tile)+3
612# endif
613!
614 CALL shapiro2d_tile (ng, tile, 1, &
615 & lbi, ubi, lbj, ubj, &
616 & imins, imaxs, jmins, jmaxs, &
617# ifdef MASKING
618 & mymask, &
619# endif
620 & fout)
621# endif
622!
623! Compute global interpolated field minimum and maximum values.
624! Notice that gridded data values are overwritten.
625!
626# ifdef DISTRIBUTE
627 rbuffer(1)=my_min
628 rbuffer(2)=my_max
629 op_handle(1)='MIN'
630 op_handle(2)='MAX'
631 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
632 amin=rbuffer(1)
633 amax=rbuffer(2)
634# else
635 amin=my_min
636 amax=my_max
637# endif
638!
639 RETURN
640 END SUBROUTINE regrid_pio
641#endif
642 END MODULE regrid_mod
subroutine hindices(ng, lbi, ubi, lbj, ubj, is, ie, js, je, angler, xgrd, ygrd, lbm, ubm, lbn, ubn, ms, me, ns, ne, xpos, ypos, ipos, jpos, ijspv, rectangular)
subroutine cinterp2d(ng, lbx, ubx, lby, uby, xinp, yinp, finp, lbi, ubi, lbj, ubj, istr, iend, jstr, jend, iout, jout, xout, yout, fout, minval, maxval)
subroutine linterp2d(ng, lbx, ubx, lby, uby, xinp, yinp, finp, lbi, ubi, lbj, ubj, istr, iend, jstr, jend, iout, jout, xout, yout, fout, minval, maxval)
integer stdout
logical master
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter r3dvar
Definition mod_param.F:721
integer nghostpoints
Definition mod_param.F:710
integer, parameter u3dvar
Definition mod_param.F:722
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
integer, parameter p3dvar
Definition mod_param.F:720
integer, parameter v3dvar
Definition mod_param.F:723
real(r8), dimension(:), allocatable latmax
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
integer, parameter linear
real(r8), dimension(:), allocatable latmin
real(r8), dimension(:), allocatable lonmax
integer, parameter cubic
integer exit_flag
real(r8), dimension(:), allocatable lonmin
integer noerror
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public regrid_nf90(ng, model, ncname, ncid, ncvname, ncvarid, gtype, iflag, nx, ny, finp, amin, amax, lbi, ubi, lbj, ubj, imin, imax, jmin, jmax, mymask, myxout, xout, yout, fout)
Definition regrid.F:107
subroutine, public regrid_pio(ng, model, ncname, piofile, ncvname, piovar, gtype, iflag, nx, ny, finp, amin, amax, lbi, ubi, lbj, ubj, imin, imax, jmin, jmax, mymask, myxout, xout, yout, fout)
Definition regrid.F:378
subroutine shapiro2d_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, amask, a)
Definition shapiro.F:33
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52