20#if defined PIO_LIB && defined DISTRIBUTE
27#if defined PIO_LIB && defined DISTRIBUTE
35#if defined PIO_LIB && defined DISTRIBUTE
71#if defined PIO_LIB && defined DISTRIBUTE
78#if defined PARALLEL_IO && defined DISTRIBUTE
83 & tindex, gtype, Vsize, &
84 & LBi, UBi, LBj, UBj, &
91 & Lregrid)
RESULT (status)
102 logical,
intent(out),
optional :: lregrid
104 integer,
intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
105 integer,
intent(in) :: lbi, ubi, lbj, ubj
106 integer,
intent(in) :: vsize(4)
108 integer(i8b),
intent(out),
optional :: checksum
110 real(dp),
intent(in) :: ascl
111 real(r8),
intent(out) :: amin
112 real(r8),
intent(out) :: amax
114 character (len=*),
intent(in) :: ncname
115 character (len=*),
intent(in) :: ncvname
119 real(r8),
intent(in) :: amask(lbi:,lbj:)
121 real(r8),
intent(out) :: adat(lbi:,lbj:)
124 real(r8),
intent(in) :: amask(lbi:ubi,lbj:ubj)
126 real(r8),
intent(out) :: adat(lbi:ubi,lbj:ubj)
131 logical :: lchecksum, interpolate
133 logical,
dimension(3) :: foundit
135 integer :: i, ic, ij, j, jc, np, mynpts, npts
136 integer :: imin, imax, isize, jmin, jmax, jsize, ijsize
137 integer :: istr, iend, jstr, jend
138 integer :: ioff, joff, ijoff
139 integer :: ilen, itile, jlen, jtile, ijlen
140 integer :: cgrid, mytype, ghost, status, wtype
142 integer,
dimension(3) :: start, total
144 real(r8) :: afactor, aoffset, aspval
146 real(r8),
parameter :: inival = 0.0_r8
148 real(r8),
dimension(2) :: rbuffer
149 real(r8),
dimension(3) :: attvalue
151 real(r8),
allocatable :: awrk(:,:)
152# if defined MASKING && defined READ_WATER
153 real(r8),
allocatable :: a2d(:)
155 real(r8),
allocatable :: wrk(:)
157 character (len= 3),
dimension(2) :: op_handle
158 character (len=12),
dimension(3) :: attname
160 character (len=*),
parameter :: myfile = &
161 & __FILE__//
", nf90_fread2d"
182 IF (model.eq.
iadm)
THEN
191 SELECT CASE (abs(mytype))
229 IF (((vsize(1).gt.0).and.(vsize(1).ne.isize)).or. &
230 & ((vsize(2).gt.0).and.(vsize(2).ne.jsize)))
THEN
235 IF (
PRESENT(lregrid))
THEN
252 attname(1)=
'scale_factor'
253 attname(2)=
'add_offset '
254 attname(3)=
'_FillValue '
257 & attvalue, foundit, &
264 IF (.not.foundit(1))
THEN
270 IF (.not.foundit(2))
THEN
276 IF (.not.foundit(3))
THEN
284 IF (
PRESENT(checksum))
THEN
301 SELECT CASE (abs(mytype))
323 IF (interpolate)
THEN
324 IF (.not.
allocated(awrk))
THEN
325 allocate ( awrk(ilen,jlen) )
328 IF (.not.
allocated(wrk))
THEN
329 allocate ( wrk(npts) )
333 IF (.not.
allocated(wrk))
THEN
341 IF (interpolate)
THEN
343 & istr, iend, jstr, jend)
354 mynpts=total(1)*total(2)
358 status=nf90_get_var(ncid, ncvarid, wrk, start, total)
363 IF (status.eq.nf90_noerr)
THEN
367 IF (abs(wrk(i)).ge.abs(aspval))
THEN
370 wrk(i)=ascl*(afactor*wrk(i)+aoffset)
371 amin=min(amin,wrk(i))
372 amax=max(amax,wrk(i))
375 IF ((abs(amin).ge.abs(aspval)).and. &
376 & (abs(amax).ge.abs(aspval)))
THEN
387 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
394 IF (interpolate)
THEN
395 CALL mp_collect (ng, model, npts, inival, wrk)
418# if defined MASKING && defined READ_WATER
431 SELECT CASE (abs(mytype))
462 IF (.not.
allocated(a2d))
THEN
463 allocate ( a2d(ijsize) )
466 IF (.not.
allocated(wrk))
THEN
467 allocate ( wrk(npts) )
482 status=nf90_get_var(ncid, ncvarid, wrk(istr:), start, total)
487 IF (status.eq.nf90_noerr)
THEN
488 CALL mp_collect (ng, model, npts, inival, wrk)
496 IF (abs(wrk(i)).ge.abs(aspval))
THEN
499 wrk(i)=ascl*(afactor*wrk(i)+aoffset)
500 amin=min(amin,wrk(i))
501 amax=max(amax,wrk(i))
504 IF ((abs(amin).ge.abs(aspval)).and. &
505 & (abs(amax).ge.abs(aspval)))
THEN
515 ij=
scalars(ng)%IJwater(np,wtype)
538 IF (interpolate.and.(status.eq.nf90_noerr))
THEN
539 SELECT CASE (abs(mytype))
544 & ilen, jlen, awrk, amin, amax, &
545 & lbi, ubi, lbj, ubj, &
546 & imin, imax, jmin, jmax, &
548 &
grid(ng) % pmask, &
550 &
grid(ng) % MyLon, &
557 & ilen, jlen, awrk, amin, amax, &
558 & lbi, ubi, lbj, ubj, &
559 & imin, imax, jmin, jmax, &
561 &
grid(ng) % pmask, &
563 &
grid(ng) % MyLon, &
572 & ilen, jlen, awrk, amin, amax, &
573 & lbi, ubi, lbj, ubj, &
574 & imin, imax, jmin, jmax, &
576 &
grid(ng) % rmask, &
578 &
grid(ng) % MyLon, &
585 & ilen, jlen, awrk, amin, amax, &
586 & lbi, ubi, lbj, ubj, &
587 & imin, imax, jmin, jmax, &
589 &
grid(ng) % rmask, &
591 &
grid(ng) % MyLon, &
600 & ilen, jlen, awrk, amin, amax, &
601 & lbi, ubi, lbj, ubj, &
602 & imin, imax, jmin, jmax, &
604 &
grid(ng) % umask, &
606 &
grid(ng) % MyLon, &
613 & ilen, jlen, awrk, amin, amax, &
614 & lbi, ubi, lbj, ubj, &
615 & imin, imax, jmin, jmax, &
617 &
grid(ng) % umask, &
619 &
grid(ng) % MyLon, &
628 & ilen, jlen, awrk, amin, amax, &
629 & lbi, ubi, lbj, ubj, &
630 & imin, imax, jmin, jmax, &
632 &
grid(ng) % vmask, &
634 &
grid(ng) % MyLon, &
641 & ilen, jlen, awrk, amin, amax, &
642 & lbi, ubi, lbj, ubj, &
643 & imin, imax, jmin, jmax, &
645 &
grid(ng) % vmask, &
647 &
grid(ng) % MyLon, &
656 & ilen, jlen, awrk, amin, amax, &
657 & lbi, ubi, lbj, ubj, &
658 & imin, imax, jmin, jmax, &
660 &
grid(ng) % rmask, &
662 &
grid(ng) % MyLon, &
669 & ilen, jlen, awrk, amin, amax, &
670 & lbi, ubi, lbj, ubj, &
671 & imin, imax, jmin, jmax, &
673 &
grid(ng) % rmask, &
675 &
grid(ng) % MyLon, &
687 IF (interpolate)
THEN
688 IF (
allocated(awrk))
THEN
693# if defined MASKING && defined READ_WATER
694 IF (
allocated(a2d))
THEN
699 IF (
allocated(wrk))
THEN
711 & ncvname, ncvarid, &
712 & tindex, gtype, Vsize, &
713 & LBi, UBi, LBj, UBj, &
714 & Ascl, Amin, Amax, &
720 & Lregrid)
RESULT (status)
733 logical,
intent(out),
optional :: lregrid
735 integer,
intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
736 integer,
intent(in) :: lbi, ubi, lbj, ubj
737 integer,
intent(in) :: vsize(4)
739 integer(i8b),
intent(out),
optional :: checksum
741 real(dp),
intent(in) :: ascl
742 real(r8),
intent(out) :: amin
743 real(r8),
intent(out) :: amax
745 character (len=*),
intent(in) :: ncname
746 character (len=*),
intent(in) :: ncvname
750 real(r8),
intent(in) :: amask(lbi:,lbj:)
752 real(r8),
intent(out) :: adat(lbi:,lbj:)
755 real(r8),
intent(in) :: amask(lbi:ubi,lbj:ubj)
757 real(r8),
intent(out) :: adat(lbi:ubi,lbj:ubj)
762 logical :: lchecksum, interpolate
764 logical,
dimension(3) :: foundit
766 integer :: i, j, ic, npts, nwpts, status, wtype
767 integer :: is, ie, js, je
768 integer :: imin, imax, jmin, jmax
769 integer :: ilen, jlen, ijlen
770 integer :: cgrid, mytype, ghost
774 integer,
dimension(3) :: start, total
776 real(r8) :: afactor, aoffset, aspval
778 real(r8),
dimension(3) :: attvalue
780 real(r8),
allocatable :: cwrk(:)
781 real(r8),
allocatable :: wrk(:)
783 character (len=12),
dimension(3) :: attname
785 character (len=*),
parameter :: myfile = &
786 & __FILE__//
", nf90_fread2d"
800 SELECT CASE (abs(mytype))
854 IF (model.eq.
iadm)
THEN
869 IF (((vsize(1).gt.0).and.(vsize(1).ne.ilen)).or. &
870 & ((vsize(2).gt.0).and.(vsize(2).ne.jlen)))
THEN
875 IF (
PRESENT(lregrid))
THEN
893 attname(1)=
'scale_factor'
894 attname(2)=
'add_offset '
895 attname(3)=
'_FillValue '
898 & attvalue, foundit, &
905 IF (.not.foundit(1))
THEN
911 IF (.not.foundit(2))
THEN
917 IF (.not.foundit(3))
THEN
923# if defined READ_WATER && defined MASKING
928 SELECT CASE (abs(mytype))
945 nwpts=(
lm(ng)+2)*(
mm(ng)+2)
950 IF (mytype.gt.0)
THEN
958# if defined READ_WATER && defined MASKING
973 IF (.not.
allocated(wrk))
THEN
974 IF (interpolate)
THEN
975 allocate ( wrk(npts) )
977 allocate ( wrk(npts+2) )
984 IF (
PRESENT(checksum))
THEN
997 status=nf90_get_var(ncid, ncvarid, wrk, start, total)
998 IF (status.eq.nf90_noerr)
THEN
1002 IF (abs(wrk(i)).ge.abs(aspval))
THEN
1005 wrk(i)=ascl*(afactor*wrk(i)+aoffset)
1006 amin=min(amin,wrk(i))
1007 amax=max(amax,wrk(i))
1010 IF ((abs(amin).ge.abs(aspval)).and. &
1011 & (abs(amax).ge.abs(aspval)))
THEN
1020 IF (
founderror(status, nf90_noerr, __line__, myfile))
THEN
1030 IF (.not.interpolate)
THEN
1036 & nghost, mytype, amin, amax, &
1037# if defined READ_WATER && defined MASKING
1038 & nwpts,
scalars(ng)%IJwater(:,wtype), &
1046 IF (mytype.gt.0)
THEN
1054# if defined MASKING || defined READ_WATER
1059 IF (amask(i,j).gt.0.0_r8)
THEN
1081 IF (interpolate)
THEN
1082 SELECT CASE (abs(mytype))
1087 & ilen, jlen, wrk, amin, amax, &
1088 & lbi, ubi, lbj, ubj, &
1093 &
grid(ng) % MyLon, &
1094 &
grid(ng) % lonp, &
1095 &
grid(ng) % latp, &
1100 & ilen, jlen, wrk, amin, amax, &
1101 & lbi, ubi, lbj, ubj, &
1106 &
grid(ng) % MyLon, &
1115 & ilen, jlen, wrk, amin, amax, &
1116 & lbi, ubi, lbj, ubj, &
1119 &
grid(ng) % rmask, &
1121 &
grid(ng) % MyLon, &
1122 &
grid(ng) % lonr, &
1123 &
grid(ng) % latr, &
1128 & ilen, jlen, wrk, amin, amax, &
1129 & lbi, ubi, lbj, ubj, &
1132 &
grid(ng) % rmask, &
1134 &
grid(ng) % MyLon, &
1143 & ilen, jlen, wrk, amin, amax, &
1144 & lbi, ubi, lbj, ubj, &
1147 &
grid(ng) % umask, &
1149 &
grid(ng) % MyLon, &
1150 &
grid(ng) % lonu, &
1151 &
grid(ng) % latu, &
1156 & ilen, jlen, wrk, amin, amax, &
1157 & lbi, ubi, lbj, ubj, &
1160 &
grid(ng) % umask, &
1162 &
grid(ng) % MyLon, &
1171 & ilen, jlen, wrk, amin, amax, &
1172 & lbi, ubi, lbj, ubj, &
1175 &
grid(ng) % vmask, &
1177 &
grid(ng) % MyLon, &
1178 &
grid(ng) % lonv, &
1179 &
grid(ng) % latv, &
1184 & ilen, jlen, wrk, amin, amax, &
1185 & lbi, ubi, lbj, ubj, &
1188 &
grid(ng) % vmask, &
1190 &
grid(ng) % MyLon, &
1199 & ilen, jlen, wrk, amin, amax, &
1200 & lbi, ubi, lbj, ubj, &
1203 &
grid(ng) % rmask, &
1205 &
grid(ng) % MyLon, &
1206 &
grid(ng) % lonr, &
1207 &
grid(ng) % latr, &
1212 & ilen, jlen, wrk, amin, amax, &
1213 & lbi, ubi, lbj, ubj, &
1216 &
grid(ng) % rmask, &
1218 &
grid(ng) % MyLon, &
1232 npts=(imax-imin+1)*(jmax-jmin+1)
1233 IF (.not.
allocated(cwrk))
allocate ( cwrk(npts) )
1234 cwrk = pack(adat(imin:imax, jmin:jmax), .true.)
1235 CALL get_hash (cwrk, npts, checksum, .true.)
1237 npts=(ie-is+1)*(je-js+1)
1238 IF (.not.
allocated(cwrk))
allocate ( cwrk(npts) )
1239 cwrk = pack(adat(is:ie, js:je), .true.)
1240 CALL get_hash (cwrk, npts, checksum)
1242 IF (
allocated(cwrk))
deallocate (cwrk)
1249 IF (
allocated(wrk))
THEN
1257#if defined PIO_LIB && defined DISTRIBUTE
1261 & ncvname, pioVar, &
1262 & tindex, pioDesc, Vsize, &
1263 & LBi, UBi, LBj, UBj, &
1264 & Ascl, Amin, Amax, &
1270 & Lregrid)
RESULT (status)
1280 logical,
intent(out),
optional :: lregrid
1282 integer,
intent(in) :: ng, model, tindex
1283 integer,
intent(in) :: lbi, ubi, lbj, ubj
1284 integer,
intent(in) :: vsize(4)
1286 integer(i8b),
intent(out),
optional :: checksum
1288 real(dp),
intent(in) :: ascl
1289 real(r8),
intent(out) :: amin
1290 real(r8),
intent(out) :: amax
1292 character (len=*),
intent(in) :: ncname
1293 character (len=*),
intent(in) :: ncvname
1295# ifdef ASSUMED_SHAPE
1297 real(r8),
intent(in) :: amask(lbi:,lbj:)
1299 real(r8),
intent(out) :: adat(lbi:,lbj:)
1302 real(r8),
intent(in) :: amask(lbi:ubi,lbj:ubj)
1304 real(r8),
intent(out) :: adat(lbi:ubi,lbj:ubj)
1307 TYPE (file_desc_t),
intent(inout) :: piofile
1308 TYPE (io_desc_t),
intent(inout) :: piodesc
1313 logical :: lchecksum, interpolate
1315 logical,
dimension(3) :: foundit
1317 integer :: i, j, npts, status
1318 integer :: is, ie, js, je
1319 integer :: imin, imax, jmin, jmax
1320 integer :: ilen, jlen, ijlen
1321 integer :: cgrid, ghost, dkind, gtype
1323 integer,
dimension(3) :: start, total
1325 real(r8) :: afactor, aoffset, aspval, avalue
1326 real(r8) :: my_amin, my_amax
1328 real(r8),
dimension(3) :: attvalue
1329 real(r8),
dimension(2) :: rbuffer
1331 real(r4),
pointer :: awrk4(:,:)
1332 real(r8),
pointer :: awrk8(:,:)
1333 real(r8),
allocatable :: cwrk(:)
1334 real(r8),
allocatable :: wrk(:,:)
1336 character (len=12),
dimension(3) :: attname
1337 character (len= 3),
dimension(2) :: op_handle
1339 character (len=*),
parameter :: myfile = &
1340 & __FILE__//
", pio_fread2d"
1359 SELECT CASE (abs(gtype))
1413 IF (((vsize(1).gt.0).and.(vsize(1).ne.ilen)).or. &
1414 & ((vsize(2).gt.0).and.(vsize(2).ne.jlen)))
THEN
1419 IF (
PRESENT(lregrid))
THEN
1437 attname(1)=
'scale_factor'
1438 attname(2)=
'add_offset '
1439 attname(3)=
'_FillValue '
1442 & attvalue, foundit, &
1443 & piofile = piofile)
1449 IF (.not.foundit(1))
THEN
1455 IF (.not.foundit(2))
THEN
1461 IF (.not.foundit(3))
THEN
1471 IF (interpolate)
THEN
1472 IF (.not.
allocated(wrk))
THEN
1473 allocate ( wrk(ilen,jlen) )
1490 & piofile = piofile, &
1493 & broadcast = .false., &
1501 wrk(i,j)=ascl*wrk(i,j)
1508 IF (
PRESENT(checksum))
THEN
1519 IF (.not.interpolate)
THEN
1525 IF (dkind.eq.pio_double)
THEN
1526 IF (.not.
associated(awrk8))
THEN
1527 allocate ( awrk8(imin:imax, jmin:jmax) )
1531 IF (.not.
associated(awrk4))
THEN
1532 allocate ( awrk4(imin:imax, jmin:jmax) )
1539 IF (tindex.gt.0)
THEN
1540 CALL pio_setframe (piofile, &
1542 & int(tindex, kind=pio_offset_kind))
1547 IF (dkind.eq.pio_double)
THEN
1548 CALL pio_read_darray (piofile, &
1551 & awrk8(imin:,jmin:), &
1556 IF (abs(awrk8(i,j)).ge.abs(aspval))
THEN
1559 avalue=ascl*(afactor*awrk8(i,j)+aoffset)
1561 my_amin=min(my_amin,avalue)
1562 my_amax=max(my_amax,avalue)
1566 IF (
associated(awrk8))
deallocate (awrk8)
1571 CALL pio_read_darray (piofile, &
1574 & awrk4(imin:,jmin:), &
1579 IF (abs(awrk4(i,j)).ge.abs(real(aspval,r4)))
THEN
1582 avalue=real(ascl*(afactor*awrk4(i,j)+aoffset),r8)
1584 my_amin=real(min(my_amin,avalue),r8)
1585 my_amax=real(max(my_amax,avalue),r8)
1589 IF (
associated(awrk4))
deallocate (awrk4)
1598 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
1602 IF ((abs(amin).ge.abs(
spval)).and. &
1603 & (abs(amax).ge.abs(
spval)))
THEN
1614 IF (interpolate)
THEN
1618 CALL regrid_pio (ng, model, ncname, piofile, &
1620 & ilen, jlen, wrk, amin, amax, &
1621 & lbi, ubi, lbj, ubj, &
1622 & imin, imax, jmin, jmax, &
1626 &
grid(ng) % MyLon, &
1627 &
grid(ng) % lonp, &
1628 &
grid(ng) % latp, &
1631 CALL regrid_pio (ng, model, ncname, piofile, &
1633 & ilen, jlen, wrk, amin, amax, &
1634 & lbi, ubi, lbj, ubj, &
1635 & imin, imax, jmin, jmax, &
1639 &
grid(ng) % MyLon, &
1646 CALL regrid_pio (ng, model, ncname, piofile, &
1648 & ilen, jlen, wrk, amin, amax, &
1649 & lbi, ubi, lbj, ubj, &
1650 & imin, imax, jmin, jmax, &
1652 &
grid(ng) % rmask, &
1654 &
grid(ng) % MyLon, &
1655 &
grid(ng) % lonr, &
1656 &
grid(ng) % latr, &
1659 CALL regrid_pio (ng, model, ncname, piofile, &
1661 & ilen, jlen, wrk, amin, amax, &
1662 & lbi, ubi, lbj, ubj, &
1663 & imin, imax, jmin, jmax, &
1665 &
grid(ng) % rmask, &
1667 &
grid(ng) % MyLon, &
1674 CALL regrid_pio (ng, model, ncname, piofile, &
1676 & ilen, jlen, wrk, amin, amax, &
1677 & lbi, ubi, lbj, ubj, &
1678 & imin, imax, jmin, jmax, &
1680 &
grid(ng) % umask, &
1682 &
grid(ng) % MyLon, &
1683 &
grid(ng) % lonu, &
1684 &
grid(ng) % latu, &
1687 CALL regrid_pio (ng, model, ncname, piofile, &
1689 & ilen, jlen, wrk, amin, amax, &
1690 & lbi, ubi, lbj, ubj, &
1691 & imin, imax, jmin, jmax, &
1693 &
grid(ng) % umask, &
1695 &
grid(ng) % MyLon, &
1702 CALL regrid_pio (ng, model, ncname, piofile, &
1704 & ilen, jlen, wrk, amin, amax, &
1705 & lbi, ubi, lbj, ubj, &
1706 & imin, imax, jmin, jmax, &
1708 &
grid(ng) % vmask, &
1710 &
grid(ng) % MyLon, &
1711 &
grid(ng) % lonv, &
1712 &
grid(ng) % latv, &
1715 CALL regrid_pio (ng, model, ncname, piofile, &
1717 & ilen, jlen, wrk, amin, amax, &
1718 & lbi, ubi, lbj, ubj, &
1719 & imin, imax, jmin, jmax, &
1721 &
grid(ng) % vmask, &
1723 &
grid(ng) % MyLon, &
1730 CALL regrid_pio (ng, model, ncname, piofile, &
1732 & ilen, jlen, wrk, amin, amax, &
1733 & lbi, ubi, lbj, ubj, &
1734 & imin, imax, jmin, jmax, &
1736 &
grid(ng) % rmask, &
1738 &
grid(ng) % MyLon, &
1739 &
grid(ng) % lonr, &
1740 &
grid(ng) % latr, &
1743 CALL regrid_pio (ng, model, ncname, piofile, &
1745 & ilen, jlen, wrk, amin, amax, &
1746 & lbi, ubi, lbj, ubj, &
1747 & imin, imax, jmin, jmax, &
1749 &
grid(ng) % rmask, &
1751 &
grid(ng) % MyLon, &
1760 IF (
allocated(wrk))
THEN
1768 npts=(imax-imin+1)*(jmax-jmin+1)
1769 IF (.not.
allocated(cwrk))
allocate ( cwrk(npts) )
1770 cwrk=pack(adat(imin:imax, jmin:jmax), .true.)
1771 CALL get_hash (cwrk, npts, checksum, .true.)
1772 IF (
allocated(cwrk))
deallocate (cwrk)
subroutine mp_scatter2d(ng, model, lbi, ubi, lbj, ubj, nghost, gtype, amin, amax, nwpts, ij_water, npts, a, awrk)
subroutine, public tile_bounds_2d(ng, tile, imax, jmax, itile, jtile, istr, iend, jstr, jend)
subroutine, public tile_bounds_1d(ng, tile, imax, istr, iend)
subroutine, public get_hash(a, asize, hash, lreduce)
type(t_grid), dimension(:), allocatable grid
integer, dimension(:), allocatable tilesize
type(t_bounds), dimension(:), allocatable bounds
integer, parameter r3dvar
type(t_iobounds), dimension(:), allocatable iobounds
integer, parameter u3dvar
integer, dimension(:), allocatable lm
integer, parameter u2dvar
integer, parameter p2dvar
integer, dimension(:), allocatable mm
integer, parameter r2dvar
integer, parameter v2dvar
integer, parameter p3dvar
integer, parameter v3dvar
real(dp), parameter spval
real(dp), parameter spval_check
type(t_scalars), dimension(:), allocatable scalars
integer function nf90_fread2d(ng, model, ncname, ncid, ncvname, ncvarid, tindex, gtype, vsize, lbi, ubi, lbj, ubj, ascl, amin, amax, amask, adat, checksum, lregrid)
integer function pio_fread2d(ng, model, ncname, piofile, ncvname, piovar, tindex, piodesc, vsize, lbi, ubi, lbj, ubj, ascl, amin, amax, amask, adat, checksum, lregrid)
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)
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)
logical function, public founderror(flag, noerr, line, routine)