ROMS
Loading...
Searching...
No Matches
nf_fread2d.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 function reads in a generic floating point 2D array from an !
12! input NetCDF file. !
13! !
14! On Input: !
15! !
16! ng Nested grid number (integer) !
17! model Calling model identifier (integer) !
18! ncname NetCDF file name (string) !
19! ncid NetCDF file ID (integer) !
20#if defined PIO_LIB && defined DISTRIBUTE
21!or pioFile PIO file descriptor structure, TYPE(file_desc_t) !
22! pioFile%fh file handler !
23! pioFile%iosystem IO system descriptor (struct) !
24#endif
25! ncvname NetCDF variable name (string) !
26! ncvarid NetCDF variable ID (integer) !
27#if defined PIO_LIB && defined DISTRIBUTE
28!or pioVar PIO variable descriptor structure, TYPE(My_VarDesc) !
29! pioVar%vd variable descriptor TYPE(Var_Desc_t) !
30! pioVar%dkind variable data kind !
31! pioVar%gtype variable C-gridtype !
32#endif
33! tindex NetCDF time record index to read (integer) !
34! gtype C-grid type (integer) !
35#if defined PIO_LIB && defined DISTRIBUTE
36!or pioDesc IO data decomposition descriptor, TYPE(IO_desc_t) !
37#endif
38! Vsize Variable dimensions in NetCDF file (integer 1D array) !
39! LBi I-dimension Lower bound (integer) !
40! UBi I-dimension Upper bound (integer) !
41! LBj J-dimension Lower bound (integer) !
42! UBj J-dimension Upper bound (integer) !
43! Ascl Factor to scale field after reading (real). !
44! Amask Land/Sea mask, if any (real 2D array) !
45! !
46! On Output: !
47! !
48! Amin Field minimum value (real) !
49! Amax Field maximum value (real) !
50! Adat Field to read in (real 2D array) !
51! checksum Field checksum value (integer*8; OPTIONAL) !
52! Lregrid Field regrid interpolation switch (logical; OPTIONAL) !
53! status Result error flag (integer) !
54! !
55!=======================================================================
56!
57 USE mod_param
58 USE mod_parallel
59 USE mod_grid
60 USE mod_iounits
61 USE mod_ncparam
62 USE mod_scalars
63!
64 USE get_hash_mod, ONLY : get_hash
65 USE strings_mod, ONLY : founderror
66!
67 implicit none
68!
69 INTERFACE nf_fread2d
70 MODULE PROCEDURE nf90_fread2d
71#if defined PIO_LIB && defined DISTRIBUTE
72 MODULE PROCEDURE pio_fread2d
73#endif
74 END INTERFACE nf_fread2d
75!
76 CONTAINS
77
78#if defined PARALLEL_IO && defined DISTRIBUTE
79!
80!***********************************************************************
81 FUNCTION nf90_fread2d (ng, model, ncname, ncid, &
82 & ncvname, ncvarid, &
83 & tindex, gtype, Vsize, &
84 & LBi, UBi, LBj, UBj, &
85 & Ascl, Amin, Amax, &
86# ifdef MASKING
87 & Amask, &
88# endif
89 & Adat, &
90 & checksum, &
91 & Lregrid) RESULT (status)
92!***********************************************************************
93!
94 USE mod_netcdf
95!
98 USE regrid_mod, ONLY : regrid_nf90
99!
100! Imported variable declarations.
101!
102 logical, intent(out), optional :: lregrid
103!
104 integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
105 integer, intent(in) :: lbi, ubi, lbj, ubj
106 integer, intent(in) :: vsize(4)
107!
108 integer(i8b), intent(out), optional :: checksum
109!
110 real(dp), intent(in) :: ascl
111 real(r8), intent(out) :: amin
112 real(r8), intent(out) :: amax
113!
114 character (len=*), intent(in) :: ncname
115 character (len=*), intent(in) :: ncvname
116!
117# ifdef ASSUMED_SHAPE
118# ifdef MASKING
119 real(r8), intent(in) :: amask(lbi:,lbj:)
120# endif
121 real(r8), intent(out) :: adat(lbi:,lbj:)
122# else
123# ifdef MASKING
124 real(r8), intent(in) :: amask(lbi:ubi,lbj:ubj)
125# endif
126 real(r8), intent(out) :: adat(lbi:ubi,lbj:ubj)
127# endif
128!
129! Local variable declarations.
130!
131 logical :: lchecksum, interpolate
132!
133 logical, dimension(3) :: foundit
134!
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
141
142 integer, dimension(3) :: start, total
143!
144 real(r8) :: afactor, aoffset, aspval
145
146 real(r8), parameter :: inival = 0.0_r8
147
148 real(r8), dimension(2) :: rbuffer
149 real(r8), dimension(3) :: attvalue
150
151 real(r8), allocatable :: awrk(:,:)
152# if defined MASKING && defined READ_WATER
153 real(r8), allocatable :: a2d(:)
154# endif
155 real(r8), allocatable :: wrk(:)
156!
157 character (len= 3), dimension(2) :: op_handle
158 character (len=12), dimension(3) :: attname
159
160 character (len=*), parameter :: myfile = &
161 & __FILE__//", nf90_fread2d"
162!
163!-----------------------------------------------------------------------
164! Set starting and ending indices to process.
165!-----------------------------------------------------------------------
166!
167 status=nf90_noerr
168!
169! Set first and last grid point according to staggered C-grid
170! classification. Set the offsets for variables with starting
171! zero-index. Recall the NetCDF does not support a zero-index.
172!
173! Notice that (Imin,Jmin) and (Imax,Jmax) are the corner of the
174! computational tile. If ghost=0, ghost points are not processed.
175! They will be processed elsewhere by the appropriate call to any
176! of the routines in "mp_exchange.F". If ghost=1, the ghost points
177! are read.
178!
179# ifdef NO_READ_GHOST
180 ghost=0 ! non-overlapping, no ghost points
181# else
182 IF (model.eq.iadm) THEN
183 ghost=0 ! non-overlapping, no ghost points
184 ELSE
185 ghost=1 ! overlapping, read ghost points
186 END IF
187# endif
188
189 mytype=gtype
190
191 SELECT CASE (abs(mytype))
192 CASE (p2dvar, p3dvar)
193 cgrid=1
194 isize=iobounds(ng)%xi_psi
195 jsize=iobounds(ng)%eta_psi
196 CASE (r2dvar, r3dvar)
197 cgrid=2
198 isize=iobounds(ng)%xi_rho
199 jsize=iobounds(ng)%eta_rho
200 CASE (u2dvar, u3dvar)
201 cgrid=3
202 isize=iobounds(ng)%xi_u
203 jsize=iobounds(ng)%eta_u
204 CASE (v2dvar, v3dvar)
205 cgrid=4
206 isize=iobounds(ng)%xi_v
207 jsize=iobounds(ng)%eta_v
208 CASE DEFAULT
209 cgrid=2
210 isize=iobounds(ng)%xi_rho
211 jsize=iobounds(ng)%eta_rho
212 END SELECT
213
214 imin=bounds(ng)%Imin(cgrid,ghost,myrank)
215 imax=bounds(ng)%Imax(cgrid,ghost,myrank)
216 jmin=bounds(ng)%Jmin(cgrid,ghost,myrank)
217 jmax=bounds(ng)%Jmax(cgrid,ghost,myrank)
218
219 ilen=imax-imin+1
220 jlen=jmax-jmin+1
221!
222! Determine if interpolating from coarse gridded data to model grid
223! is required. This is only allowed for gridded 2D fields. This is
224! convenient for atmospheric forcing datasets that are usually on
225! coarser grids. The user can provide coarser gridded data to avoid
226! very large input files.
227!
228 interpolate=.false.
229 IF (((vsize(1).gt.0).and.(vsize(1).ne.isize)).or. &
230 & ((vsize(2).gt.0).and.(vsize(2).ne.jsize))) THEN
231 interpolate=.true.
232 ilen=vsize(1)
233 jlen=vsize(2)
234 END IF
235 IF (PRESENT(lregrid)) THEN
236 lregrid=interpolate
237 END IF
238!
239! Check if the following attributes: "scale_factor", "add_offset", and
240! "_FillValue" are present in the input NetCDF variable:
241!
242! If the "scale_value" attribute is present, the data is multiplied by
243! this factor after reading.
244! If the "add_offset" attribute is present, this value is added to the
245! data after reading.
246! If both "scale_factor" and "add_offset" attributes are present, the
247! data are first scaled before the offset is added.
248! If the "_FillValue" attribute is present, the data having this value
249! is treated as missing and it is replaced with zero. This feature it
250! is usually related with the land/sea masking.
251!
252 attname(1)='scale_factor'
253 attname(2)='add_offset '
254 attname(3)='_FillValue '
255
256 CALL netcdf_get_fatt (ng, model, ncname, ncvarid, attname, &
257 & attvalue, foundit, &
258 & ncid = ncid)
259 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
260 status=ioerror
261 RETURN
262 END IF
263
264 IF (.not.foundit(1)) THEN
265 afactor=1.0_r8
266 ELSE
267 afactor=attvalue(1)
268 END IF
269
270 IF (.not.foundit(2)) THEN
271 aoffset=0.0_r8
272 ELSE
273 aoffset=attvalue(2)
274 END IF
275
276 IF (.not.foundit(3)) THEN
277 aspval=spval_check
278 ELSE
279 aspval=attvalue(3)
280 END IF
281!
282! Initialize checsum value.
283!
284 IF (PRESENT(checksum)) THEN
285 lchecksum=.true.
286 checksum=0_i8b
287 ELSE
288 lchecksum=.false.
289 END IF
290!
291!-----------------------------------------------------------------------
292! Parallel I/O: Read in tile data from requested field and scale it.
293! Processing both water and land points.
294!-----------------------------------------------------------------------
295!
296 IF (gtype.gt.0) THEN
297!
298! Set offsets due the NetCDF dimensions. Recall that some output
299! variables not always start at one.
300!
301 SELECT CASE (abs(mytype))
302 CASE (p2dvar, p3dvar)
303 ioff=0
304 joff=0
305 CASE (r2dvar, r3dvar)
306 ioff=1
307 joff=1
308 CASE (u2dvar, u3dvar)
309 ioff=0
310 joff=1
311 CASE (v2dvar, v3dvar)
312 ioff=1
313 joff=0
314 CASE DEFAULT
315 ioff=1
316 joff=1
317 END SELECT
318
319 npts=ilen*jlen
320!
321! Allocate scratch work arrays.
322!
323 IF (interpolate) THEN
324 IF (.not.allocated(awrk)) THEN
325 allocate ( awrk(ilen,jlen) )
326 awrk=inival
327 END IF
328 IF (.not.allocated(wrk)) THEN
329 allocate ( wrk(npts) )
330 wrk=inival
331 END IF
332 ELSE
333 IF (.not.allocated(wrk)) THEN
334 allocate ( wrk(tilesize(ng)) )
335 wrk=0.0_r8
336 END IF
337 END IF
338!
339! Read in data: all parallel nodes read their own tile data.
340!
341 IF (interpolate) THEN
342 CALL tile_bounds_2d (ng, myrank, ilen, jlen, itile, jtile, &
343 & istr, iend, jstr, jend)
344 start(1)=istr
345 total(1)=iend-istr+1
346 start(2)=jstr
347 total(2)=jend-jstr+1
348 ELSE
349 start(1)=imin+ioff
350 total(1)=ilen
351 start(2)=jmin+joff
352 total(2)=jlen
353 END IF
354 mynpts=total(1)*total(2)
355 start(3)=tindex
356 total(3)=1
357
358 status=nf90_get_var(ncid, ncvarid, wrk, start, total)
359!
360! Scale read data and process fill values, if any. Compute minimum
361! and maximum values.
362!
363 IF (status.eq.nf90_noerr) THEN
364 amin=spval
365 amax=-spval
366 DO i=1,mynpts
367 IF (abs(wrk(i)).ge.abs(aspval)) THEN
368 wrk(i)=0.0_r8 ! masked with _FillValue
369 ELSE
370 wrk(i)=ascl*(afactor*wrk(i)+aoffset)
371 amin=min(amin,wrk(i))
372 amax=max(amax,wrk(i))
373 END IF
374 END DO
375 IF ((abs(amin).ge.abs(aspval)).and. &
376 & (abs(amax).ge.abs(aspval))) THEN
377 amin=0.0_r8 ! the entire data is all
378 amax=0.0_r8 ! field value, _FillValue
379 END IF
380!
381! Set minimum and maximum values: global reduction.
382!
383 rbuffer(1)=amin
384 op_handle(1)='MIN'
385 rbuffer(2)=amax
386 op_handle(2)='MAX'
387 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
388 amin=rbuffer(1)
389 amax=rbuffer(2)
390!
391! Unpack read data. If not interpolating, the data is loaded into
392! model array.
393!
394 IF (interpolate) THEN
395 CALL mp_collect (ng, model, npts, inival, wrk)
396 ic=0
397 DO j=jstr,jend
398 DO i=istr,iend
399 ic=ic+1
400 awrk(i,j)=wrk(ic)
401 END DO
402 END DO
403 ELSE
404 ic=0
405 DO j=jmin,jmax
406 DO i=imin,imax
407 ic=ic+1
408 adat(i,j)=wrk(ic)
409 END DO
410 END DO
411 END IF
412 ELSE
413 exit_flag=2
414 ioerror=status
415 END IF
416 END IF
417
418# if defined MASKING && defined READ_WATER
419!
420!-----------------------------------------------------------------------
421! Parallel I/O: Read in tile data from requested field and scale it.
422! Processing water points only. Interpolation is not
423! allowed.
424!-----------------------------------------------------------------------
425!
426 IF (gtype.lt.0) THEN
427!
428! Set number of points to process, grid type switch, and offsets due
429! array packing into 1D array in column-major order.
430!
431 SELECT CASE (abs(mytype))
432 CASE (p2dvar)
433 npts=iobounds(ng)%xy_psi
434 wtype=p2dvar
435 ioff=0
436 joff=1
437 CASE (r2dvar)
438 npts=iobounds(ng)%xy_rho
439 wtype=r2dvar
440 ioff=1
441 joff=0
442 CASE (u2dvar)
443 npts=iobounds(ng)%xy_u
444 wtype=u2dvar
445 ioff=0
446 joff=0
447 CASE (v2dvar)
448 npts=iobounds(ng)%xy_v
449 wtype=v2dvar
450 ioff=1
451 joff=1
452 CASE DEFAULT
453 npts=iobounds(ng)%xy_rho
454 wtype=r2dvar
455 ioff=1
456 joff=0
457 END SELECT
458 ijsize=isize*jsize
459!
460! Allocate scratch work arrays.
461!
462 IF (.not.allocated(a2d)) THEN
463 allocate ( a2d(ijsize) )
464 a2d=inival
465 END IF
466 IF (.not.allocated(wrk)) THEN
467 allocate ( wrk(npts) )
468 wrk=inival
469 END IF
470!
471! Read in data: all parallel nodes read a segment of the 1D data.
472! Recall that water points are pack in the NetCDF file in a single
473! dimension.
474!
475 CALL tile_bounds_1d (ng, myrank, npts, istr, iend)
476
477 start(1)=istr
478 total(1)=iend-istr+1
479 start(2)=1
480 total(2)=tindex
481
482 status=nf90_get_var(ncid, ncvarid, wrk(istr:), start, total)
483!
484! Global reduction of work array. We need this because the packing
485! of the water point only affects the model tile partition.
486!
487 IF (status.eq.nf90_noerr) THEN
488 CALL mp_collect (ng, model, npts, inival, wrk)
489!
490! Scale read data and process fill values, if any. Compute minimum
491! and maximum values.
492!
493 amin=spval
494 amax=-spval
495 DO i=1,npts
496 IF (abs(wrk(i)).ge.abs(aspval)) THEN
497 wrk(i)=0.0_r8 ! masked with _FillValue
498 ELSE
499 wrk(i)=ascl*(afactor*wrk(i)+aoffset)
500 amin=min(amin,wrk(i))
501 amax=max(amax,wrk(i))
502 END IF
503 END DO
504 IF ((abs(amin).ge.abs(aspval)).and. &
505 & (abs(amax).ge.abs(aspval))) THEN
506 amin=0.0_r8 ! the entire data is all
507 amax=0.0_r8 ! field value, _FillValue
508 END IF
509!
510! Unpack read data. This is tricky in parallel I/O. The cheapeast
511! thing to do is reconstruct a packed global array and then select
512! the appropriate values for the tile.
513!
514 DO np=1,npts
515 ij=scalars(ng)%IJwater(np,wtype)
516 a2d(ij)=wrk(np)
517 END DO
518
519 DO j=jmin,jmax
520 jc=(j-joff)*isize
521 DO i=imin,imax
522 ij=i+ioff+jc
523 adat(i,j)=a2d(ij)
524 END DO
525 END DO
526 ELSE
527 exit_flag=2
528 ioerror=status
529 END IF
530 END IF
531# endif
532!
533!-----------------------------------------------------------------------
534! Parallel I/O: If interpolating from gridded data, read its associated
535! locations and interpolate.
536!-----------------------------------------------------------------------
537!
538 IF (interpolate.and.(status.eq.nf90_noerr)) THEN
539 SELECT CASE (abs(mytype))
540 CASE (p2dvar, p3dvar)
541 IF (spherical) THEN
542 CALL regrid_nf90 (ng, model, ncname, ncid, &
543 & ncvname, ncvarid, mytype, interpflag, &
544 & ilen, jlen, awrk, amin, amax, &
545 & lbi, ubi, lbj, ubj, &
546 & imin, imax, jmin, jmax, &
547# ifdef MASKING
548 & grid(ng) % pmask, &
549# endif
550 & grid(ng) % MyLon, &
551 & grid(ng) % lonp, &
552 & grid(ng) % latp, &
553 & adat)
554 ELSE
555 CALL regrid_nf90 (ng, model, ncname, ncid, &
556 & ncvname, ncvarid, mytype, interpflag, &
557 & ilen, jlen, awrk, amin, amax, &
558 & lbi, ubi, lbj, ubj, &
559 & imin, imax, jmin, jmax, &
560# ifdef MASKING
561 & grid(ng) % pmask, &
562# endif
563 & grid(ng) % MyLon, &
564 & grid(ng) % xp, &
565 & grid(ng) % yp, &
566 & adat)
567 END IF
568 CASE (r2dvar, r3dvar)
569 IF (spherical) THEN
570 CALL regrid_nf90 (ng, model, ncname, ncid, &
571 & ncvname, ncvarid, mytype, interpflag, &
572 & ilen, jlen, awrk, amin, amax, &
573 & lbi, ubi, lbj, ubj, &
574 & imin, imax, jmin, jmax, &
575# ifdef MASKING
576 & grid(ng) % rmask, &
577# endif
578 & grid(ng) % MyLon, &
579 & grid(ng) % lonr, &
580 & grid(ng) % latr, &
581 & adat)
582 ELSE
583 CALL regrid_nf90 (ng, model, ncname, ncid, &
584 & ncvname, ncvarid, mytype, interpflag, &
585 & ilen, jlen, awrk, amin, amax, &
586 & lbi, ubi, lbj, ubj, &
587 & imin, imax, jmin, jmax, &
588# ifdef MASKING
589 & grid(ng) % rmask, &
590# endif
591 & grid(ng) % MyLon, &
592 & grid(ng) % xr, &
593 & grid(ng) % yr, &
594 & adat)
595 END IF
596 CASE (u2dvar, u3dvar)
597 IF (spherical) THEN
598 CALL regrid_nf90 (ng, model, ncname, ncid, &
599 & ncvname, ncvarid, mytype, interpflag, &
600 & ilen, jlen, awrk, amin, amax, &
601 & lbi, ubi, lbj, ubj, &
602 & imin, imax, jmin, jmax, &
603# ifdef MASKING
604 & grid(ng) % umask, &
605# endif
606 & grid(ng) % MyLon, &
607 & grid(ng) % lonu, &
608 & grid(ng) % latu, &
609 & adat)
610 ELSE
611 CALL regrid_nf90 (ng, model, ncname, ncid, &
612 & ncvname, ncvarid, mytype, interpflag, &
613 & ilen, jlen, awrk, amin, amax, &
614 & lbi, ubi, lbj, ubj, &
615 & imin, imax, jmin, jmax, &
616# ifdef MASKING
617 & grid(ng) % umask, &
618# endif
619 & grid(ng) % MyLon, &
620 & grid(ng) % xu, &
621 & grid(ng) % yu, &
622 & adat)
623 END IF
624 CASE (v2dvar, v3dvar)
625 IF (spherical) THEN
626 CALL regrid_nf90 (ng, model, ncname, ncid, &
627 & ncvname, ncvarid, mytype, interpflag, &
628 & ilen, jlen, awrk, amin, amax, &
629 & lbi, ubi, lbj, ubj, &
630 & imin, imax, jmin, jmax, &
631# ifdef MASKING
632 & grid(ng) % vmask, &
633# endif
634 & grid(ng) % MyLon, &
635 & grid(ng) % lonv, &
636 & grid(ng) % latv, &
637 & adat)
638 ELSE
639 CALL regrid_nf90 (ng, model, ncname, ncid, &
640 & ncvname, ncvarid, mytype, interpflag, &
641 & ilen, jlen, awrk, amin, amax, &
642 & lbi, ubi, lbj, ubj, &
643 & imin, imax, jmin, jmax, &
644# ifdef MASKING
645 & grid(ng) % vmask, &
646# endif
647 & grid(ng) % MyLon, &
648 & grid(ng) % xv, &
649 & grid(ng) % yv, &
650 & adat)
651 END IF
652 CASE DEFAULT
653 IF (spherical) THEN
654 CALL regrid_nf90 (ng, model, ncname, ncid, &
655 & ncvname, ncvarid, mytype, interpflag, &
656 & ilen, jlen, awrk, amin, amax, &
657 & lbi, ubi, lbj, ubj, &
658 & imin, imax, jmin, jmax, &
659# ifdef MASKING
660 & grid(ng) % rmask, &
661# endif
662 & grid(ng) % MyLon, &
663 & grid(ng) % lonr, &
664 & grid(ng) % latr, &
665 & adat)
666 ELSE
667 CALL regrid_nf90 (ng, model, ncname, ncid, &
668 & ncvname, ncvarid, mytype, interpflag, &
669 & ilen, jlen, awrk, amin, amax, &
670 & lbi, ubi, lbj, ubj, &
671 & imin, imax, jmin, jmax, &
672# ifdef MASKING
673 & grid(ng) % rmask, &
674# endif
675 & grid(ng) % MyLon, &
676 & grid(ng) % xr, &
677 & grid(ng) % yr, &
678 & adat)
679 END IF
680 END SELECT
681 END IF
682!
683!-----------------------------------------------------------------------
684! Deallocate scratch work arrays.
685!-----------------------------------------------------------------------
686!
687 IF (interpolate) THEN
688 IF (allocated(awrk)) THEN
689 deallocate ( awrk )
690 END IF
691 END IF
692
693# if defined MASKING && defined READ_WATER
694 IF (allocated(a2d)) THEN
695 deallocate (a2d)
696 END IF
697# endif
698
699 IF (allocated(wrk)) THEN
700 deallocate (wrk)
701 END IF
702!
703 RETURN
704 END FUNCTION nf90_fread2d
705
706#else
707
708!
709!***********************************************************************
710 FUNCTION nf90_fread2d (ng, model, ncname, ncid, &
711 & ncvname, ncvarid, &
712 & tindex, gtype, Vsize, &
713 & LBi, UBi, LBj, UBj, &
714 & Ascl, Amin, Amax, &
715# ifdef MASKING
716 & Amask, &
717# endif
718 & Adat, &
719 & checksum, &
720 & Lregrid) RESULT (status)
721!***********************************************************************
722!
723 USE mod_netcdf
724
725# ifdef DISTRIBUTE
726!
728# endif
729 USE regrid_mod, ONLY : regrid_nf90
730!
731! Imported variable declarations.
732!
733 logical, intent(out), optional :: lregrid
734!
735 integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
736 integer, intent(in) :: lbi, ubi, lbj, ubj
737 integer, intent(in) :: vsize(4)
738!
739 integer(i8b), intent(out), optional :: checksum
740!
741 real(dp), intent(in) :: ascl
742 real(r8), intent(out) :: amin
743 real(r8), intent(out) :: amax
744!
745 character (len=*), intent(in) :: ncname
746 character (len=*), intent(in) :: ncvname
747!
748# ifdef ASSUMED_SHAPE
749# ifdef MASKING
750 real(r8), intent(in) :: amask(lbi:,lbj:)
751# endif
752 real(r8), intent(out) :: adat(lbi:,lbj:)
753# else
754# ifdef MASKING
755 real(r8), intent(in) :: amask(lbi:ubi,lbj:ubj)
756# endif
757 real(r8), intent(out) :: adat(lbi:ubi,lbj:ubj)
758# endif
759!
760! Local variable declarations.
761!
762 logical :: lchecksum, interpolate
763
764 logical, dimension(3) :: foundit
765!
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
771# ifdef DISTRIBUTE
772 integer :: nghost
773# endif
774 integer, dimension(3) :: start, total
775!
776 real(r8) :: afactor, aoffset, aspval
777
778 real(r8), dimension(3) :: attvalue
779
780 real(r8), allocatable :: cwrk(:) ! used for checksum
781 real(r8), allocatable :: wrk(:)
782!
783 character (len=12), dimension(3) :: attname
784
785 character (len=*), parameter :: myfile = &
786 & __FILE__//", nf90_fread2d"
787!
788!-----------------------------------------------------------------------
789! Set starting and ending indices to process.
790!-----------------------------------------------------------------------
791!
792 status=nf90_noerr
793!
794! Set global (interior plus boundary) starting and ending grid cell
795! indices in the I- and J-directions according to staggered C-grid
796! classification.
797!
798 mytype=gtype
799
800 SELECT CASE (abs(mytype))
801 CASE (p2dvar)
802 cgrid=1 ! PSI-points
803 is=iobounds(ng)%ILB_psi
804 ie=iobounds(ng)%IUB_psi
805 js=iobounds(ng)%JLB_psi
806 je=iobounds(ng)%JUB_psi
807 CASE (r2dvar)
808 cgrid=2 ! RHO-points
809 is=iobounds(ng)%ILB_rho
810 ie=iobounds(ng)%IUB_rho
811 js=iobounds(ng)%JLB_rho
812 je=iobounds(ng)%JUB_rho
813 CASE (u2dvar)
814 cgrid=3 ! U-points
815 is=iobounds(ng)%ILB_u
816 ie=iobounds(ng)%IUB_u
817 js=iobounds(ng)%JLB_u
818 je=iobounds(ng)%JUB_u
819 CASE (v2dvar)
820 cgrid=4 ! V-points
821 is=iobounds(ng)%ILB_v
822 ie=iobounds(ng)%IUB_v
823 js=iobounds(ng)%JLB_v
824 je=iobounds(ng)%JUB_v
825 CASE DEFAULT
826 cgrid=2 ! RHO-points
827 is=iobounds(ng)%ILB_rho
828 ie=iobounds(ng)%IUB_rho
829 js=iobounds(ng)%JLB_rho
830 je=iobounds(ng)%JUB_rho
831 END SELECT
832
833 ilen=ie-is+1
834 jlen=je-js+1
835!
836! Set the tile computational I- and J-bounds (no ghost points).
837!
838 ghost=0
839 imin=bounds(ng)%Imin(cgrid,ghost,myrank)
840 imax=bounds(ng)%Imax(cgrid,ghost,myrank)
841 jmin=bounds(ng)%Jmin(cgrid,ghost,myrank)
842 jmax=bounds(ng)%Jmax(cgrid,ghost,myrank)
843
844# ifdef DISTRIBUTE
845!
846! Set the number of tile ghost points, Nghost, to scatter in
847! distributed-memory applications. If Nghost=0, the ghost points
848! are not processed. They will be processed elsewhere by the
849! appropriate call to any of the routines in "mp_exchange.F".
850!
851# ifdef NO_READ_GHOST
852 nghost=0 ! no ghost points exchange
853# else
854 IF (model.eq.iadm) THEN
855 nghost=0 ! no ghost points exchange
856 ELSE
857 nghost=nghostpoints ! do ghost points exchange
858 END IF
859# endif
860# endif
861!
862! Determine if interpolating from coarse gridded data to model grid
863! is required. This is only allowed for gridded 2D fields. This is
864! convenient for atmospheric forcing datasets that are usually on
865! coarser grids. The user can provide coarser gridded data to avoid
866! very large input files.
867!
868 interpolate=.false.
869 IF (((vsize(1).gt.0).and.(vsize(1).ne.ilen)).or. &
870 & ((vsize(2).gt.0).and.(vsize(2).ne.jlen))) THEN
871 interpolate=.true.
872 ilen=vsize(1)
873 jlen=vsize(2)
874 END IF
875 IF (PRESENT(lregrid)) THEN
876 lregrid=interpolate
877 END IF
878 ijlen=ilen*jlen
879!
880! Check if the following attributes: "scale_factor", "add_offset", and
881! "_FillValue" are present in the input NetCDF variable:
882!
883! If the "scale_value" attribute is present, the data is multiplied by
884! this factor after reading.
885! If the "add_offset" attribute is present, this value is added to the
886! data after reading.
887! If both "scale_factor" and "add_offset" attributes are present, the
888! data are first scaled before the offset is added.
889! If the "_FillValue" attribute is present, the data having this value
890! is treated as missing and it is replaced with zero. This feature it
891! is usually related with the land/sea masking.
892!
893 attname(1)='scale_factor'
894 attname(2)='add_offset '
895 attname(3)='_FillValue '
896
897 CALL netcdf_get_fatt (ng, model, ncname, ncvarid, attname, &
898 & attvalue, foundit, &
899 & ncid = ncid)
900 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
901 status=ioerror
902 RETURN
903 END IF
904
905 IF (.not.foundit(1)) THEN
906 afactor=1.0_r8
907 ELSE
908 afactor=attvalue(1)
909 END IF
910
911 IF (.not.foundit(2)) THEN
912 aoffset=0.0_r8
913 ELSE
914 aoffset=attvalue(2)
915 END IF
916
917 IF (.not.foundit(3)) THEN
918 aspval=spval_check
919 ELSE
920 aspval=attvalue(3)
921 END IF
922
923# if defined READ_WATER && defined MASKING
924!
925! If processing water points only, set number of points and type
926! switch.
927!
928 SELECT CASE (abs(mytype))
929 CASE (p2dvar)
930 npts=iobounds(ng)%xy_psi
931 wtype=p2dvar
932 CASE (r2dvar)
933 npts=iobounds(ng)%xy_rho
934 wtype=r2dvar
935 CASE (u2dvar)
936 npts=iobounds(ng)%xy_u
937 wtype=u2dvar
938 CASE (v2dvar)
939 npts=iobounds(ng)%xy_v
940 wtype=v2dvar
941 CASE DEFAULT
942 npts=iobounds(ng)%xy_rho
943 wtype=r2dvar
944 END SELECT
945 nwpts=(lm(ng)+2)*(mm(ng)+2)
946# endif
947!
948! Set NetCDF dimension counters for processing requested field.
949!
950 IF (mytype.gt.0) THEN
951 npts=ijlen
952 start(1)=1
953 total(1)=ilen
954 start(2)=1
955 total(2)=jlen
956 start(3)=tindex
957 total(3)=1
958# if defined READ_WATER && defined MASKING
959 ELSE
960 start(1)=1
961 total(1)=npts
962 start(2)=1
963 total(2)=tindex
964# endif
965 END IF
966!
967! Allocate scratch work vector. The dimension of this vector is
968! unknown when interpolating input data to model grid. Notice
969! that the array length is increased by two because the minimum
970! and maximum values are appended in distributed-memory
971! communications.
972!
973 IF (.not.allocated(wrk)) THEN
974 IF (interpolate) THEN
975 allocate ( wrk(npts) )
976 ELSE
977 allocate ( wrk(npts+2) )
978 END IF
979 wrk=0.0_r8
980 END IF
981!
982! Initialize checsum value.
983!
984 IF (PRESENT(checksum)) THEN
985 lchecksum=.true.
986 checksum=0_i8b
987 ELSE
988 lchecksum=.false.
989 END IF
990!
991!-----------------------------------------------------------------------
992! Serial I/O: Read in requested field and scale it.
993!-----------------------------------------------------------------------
994!
995 status=nf90_noerr
996 IF (inpthread) THEN
997 status=nf90_get_var(ncid, ncvarid, wrk, start, total)
998 IF (status.eq.nf90_noerr) THEN
999 amin=spval
1000 amax=-spval
1001 DO i=1,npts
1002 IF (abs(wrk(i)).ge.abs(aspval)) THEN
1003 wrk(i)=0.0_r8 ! masked with _FillValue
1004 ELSE
1005 wrk(i)=ascl*(afactor*wrk(i)+aoffset)
1006 amin=min(amin,wrk(i))
1007 amax=max(amax,wrk(i))
1008 END IF
1009 END DO
1010 IF ((abs(amin).ge.abs(aspval)).and. &
1011 & (abs(amax).ge.abs(aspval))) THEN
1012 amin=0.0_r8 ! the entire data is all
1013 amax=0.0_r8 ! field value, _FillValue
1014 END IF
1015 END IF
1016 END IF
1017# ifdef DISTRIBUTE
1018 CALL mp_bcasti (ng, model, status)
1019# endif
1020 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1021 exit_flag=2
1022 ioerror=status
1023 RETURN
1024 END IF
1025!
1026!-----------------------------------------------------------------------
1027! Serial I/O: If not interpolating, unpack read field.
1028!-----------------------------------------------------------------------
1029!
1030 IF (.not.interpolate) THEN
1031# ifdef DISTRIBUTE
1032!
1033! Scatter read data over the distributed memory tiles.
1034!
1035 CALL mp_scatter2d (ng, model, lbi, ubi, lbj, ubj, &
1036 & nghost, mytype, amin, amax, &
1037# if defined READ_WATER && defined MASKING
1038 & nwpts, scalars(ng)%IJwater(:,wtype), &
1039# endif
1040 & npts, wrk, adat)
1041# else
1042!
1043! Unpack data into the global array: serial, serial with partitions,
1044! and shared-memory applications.
1045!
1046 IF (mytype.gt.0) THEN
1047 ic=0
1048 DO j=js,je
1049 DO i=is,ie
1050 ic=ic+1
1051 adat(i,j)=wrk(ic)
1052 END DO
1053 END DO
1054# if defined MASKING || defined READ_WATER
1055 ELSE
1056 ic=0
1057 DO j=js,je
1058 DO i=is,ie
1059 IF (amask(i,j).gt.0.0_r8) THEN
1060 ic=ic+1
1061 adat(i,j)=wrk(ic)
1062 ELSE
1063 adat(i,j)=0.0_r8
1064 END IF
1065 END DO
1066 END DO
1067# endif
1068 END IF
1069# endif
1070# ifdef DISTRIBUTE
1071 ELSE
1072 CALL mp_bcastf (ng, model, wrk)
1073# endif
1074 END IF
1075!
1076!-----------------------------------------------------------------------
1077! Serial I/O: If interpolating from gridded data, read its associated
1078! locations and interpolate.
1079!-----------------------------------------------------------------------
1080!
1081 IF (interpolate) THEN
1082 SELECT CASE (abs(mytype))
1083 CASE (p2dvar, p3dvar)
1084 IF (spherical) THEN
1085 CALL regrid_nf90 (ng, model, ncname, ncid, &
1086 & ncvname, ncvarid, mytype, interpflag, &
1087 & ilen, jlen, wrk, amin, amax, &
1088 & lbi, ubi, lbj, ubj, &
1089 & is, ie, js, je, &
1090# ifdef MASKING
1091 & amask, &
1092# endif
1093 & grid(ng) % MyLon, &
1094 & grid(ng) % lonp, &
1095 & grid(ng) % latp, &
1096 & adat)
1097 ELSE
1098 CALL regrid_nf90 (ng, model, ncname, ncid, &
1099 & ncvname, ncvarid, mytype, interpflag, &
1100 & ilen, jlen, wrk, amin, amax, &
1101 & lbi, ubi, lbj, ubj, &
1102 & is, ie, js, je, &
1103# ifdef MASKING
1104 & amask, &
1105# endif
1106 & grid(ng) % MyLon, &
1107 & grid(ng) % xp, &
1108 & grid(ng) % yp, &
1109 & adat)
1110 END IF
1111 CASE (r2dvar, r3dvar)
1112 IF (spherical) THEN
1113 CALL regrid_nf90 (ng, model, ncname, ncid, &
1114 & ncvname, ncvarid, mytype, interpflag, &
1115 & ilen, jlen, wrk, amin, amax, &
1116 & lbi, ubi, lbj, ubj, &
1117 & is, ie, js, je, &
1118# ifdef MASKING
1119 & grid(ng) % rmask, &
1120# endif
1121 & grid(ng) % MyLon, &
1122 & grid(ng) % lonr, &
1123 & grid(ng) % latr, &
1124 & adat)
1125 ELSE
1126 CALL regrid_nf90 (ng, model, ncname, ncid, &
1127 & ncvname, ncvarid, mytype, interpflag, &
1128 & ilen, jlen, wrk, amin, amax, &
1129 & lbi, ubi, lbj, ubj, &
1130 & is, ie, js, je, &
1131# ifdef MASKING
1132 & grid(ng) % rmask, &
1133# endif
1134 & grid(ng) % MyLon, &
1135 & grid(ng) % xr, &
1136 & grid(ng) % yr, &
1137 & adat)
1138 END IF
1139 CASE (u2dvar, u3dvar)
1140 IF (spherical) THEN
1141 CALL regrid_nf90 (ng, model, ncname, ncid, &
1142 & ncvname, ncvarid, mytype, interpflag, &
1143 & ilen, jlen, wrk, amin, amax, &
1144 & lbi, ubi, lbj, ubj, &
1145 & is, ie, js, je, &
1146# ifdef MASKING
1147 & grid(ng) % umask, &
1148# endif
1149 & grid(ng) % MyLon, &
1150 & grid(ng) % lonu, &
1151 & grid(ng) % latu, &
1152 & adat)
1153 ELSE
1154 CALL regrid_nf90 (ng, model, ncname, ncid, &
1155 & ncvname, ncvarid, mytype, interpflag, &
1156 & ilen, jlen, wrk, amin, amax, &
1157 & lbi, ubi, lbj, ubj, &
1158 & is, ie, js, je, &
1159# ifdef MASKING
1160 & grid(ng) % umask, &
1161# endif
1162 & grid(ng) % MyLon, &
1163 & grid(ng) % xu, &
1164 & grid(ng) % yu, &
1165 & adat)
1166 END IF
1167 CASE (v2dvar, v3dvar)
1168 IF (spherical) THEN
1169 CALL regrid_nf90 (ng, model, ncname, ncid, &
1170 & ncvname, ncvarid, mytype, interpflag, &
1171 & ilen, jlen, wrk, amin, amax, &
1172 & lbi, ubi, lbj, ubj, &
1173 & is, ie, js, je, &
1174# ifdef MASKING
1175 & grid(ng) % vmask, &
1176# endif
1177 & grid(ng) % MyLon, &
1178 & grid(ng) % lonv, &
1179 & grid(ng) % latv, &
1180 & adat)
1181 ELSE
1182 CALL regrid_nf90 (ng, model, ncname, ncid, &
1183 & ncvname, ncvarid, mytype, interpflag, &
1184 & ilen, jlen, wrk, amin, amax, &
1185 & lbi, ubi, lbj, ubj, &
1186 & is, ie, js, je, &
1187# ifdef MASKING
1188 & grid(ng) % vmask, &
1189# endif
1190 & grid(ng) % MyLon, &
1191 & grid(ng) % xv, &
1192 & grid(ng) % yv, &
1193 & adat)
1194 END IF
1195 CASE DEFAULT
1196 IF (spherical) THEN
1197 CALL regrid_nf90 (ng, model, ncname, ncid, &
1198 & ncvname, ncvarid, mytype, interpflag, &
1199 & ilen, jlen, wrk, amin, amax, &
1200 & lbi, ubi, lbj, ubj, &
1201 & is, ie, js, je, &
1202# ifdef MASKING
1203 & grid(ng) % rmask, &
1204# endif
1205 & grid(ng) % MyLon, &
1206 & grid(ng) % lonr, &
1207 & grid(ng) % latr, &
1208 & adat)
1209 ELSE
1210 CALL regrid_nf90 (ng, model, ncname, ncid, &
1211 & ncvname, ncvarid, mytype, interpflag, &
1212 & ilen, jlen, wrk, amin, amax, &
1213 & lbi, ubi, lbj, ubj, &
1214 & is, ie, js, je, &
1215# ifdef MASKING
1216 & grid(ng) % rmask, &
1217# endif
1218 & grid(ng) % MyLon, &
1219 & grid(ng) % xr, &
1220 & grid(ng) % yr, &
1221 & adat)
1222 END IF
1223 END SELECT
1224 END IF
1225!
1226!-----------------------------------------------------------------------
1227! If requested, compute data checksum value.
1228!-----------------------------------------------------------------------
1229!
1230 IF (lchecksum) THEN
1231# ifdef DISTRIBUTE
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.)
1236# else
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)
1241# endif
1242 IF (allocated(cwrk)) deallocate (cwrk)
1243 END IF
1244!
1245!-----------------------------------------------------------------------
1246! Deallocate scratch work vector.
1247!-----------------------------------------------------------------------
1248!
1249 IF (allocated(wrk)) THEN
1250 deallocate (wrk)
1251 END IF
1252!
1253 RETURN
1254 END FUNCTION nf90_fread2d
1255#endif
1256
1257#if defined PIO_LIB && defined DISTRIBUTE
1258!
1259!***********************************************************************
1260 FUNCTION pio_fread2d (ng, model, ncname, pioFile, &
1261 & ncvname, pioVar, &
1262 & tindex, pioDesc, Vsize, &
1263 & LBi, UBi, LBj, UBj, &
1264 & Ascl, Amin, Amax, &
1265# ifdef MASKING
1266 & Amask, &
1267# endif
1268 & Adat, &
1269 & checksum, &
1270 & Lregrid) RESULT (status)
1271!***********************************************************************
1272!
1273 USE mod_pio_netcdf
1274!
1275 USE distribute_mod, ONLY : mp_reduce
1276 USE regrid_mod, ONLY : regrid_pio
1277!
1278! Imported variable declarations.
1279!
1280 logical, intent(out), optional :: lregrid
1281!
1282 integer, intent(in) :: ng, model, tindex
1283 integer, intent(in) :: lbi, ubi, lbj, ubj
1284 integer, intent(in) :: vsize(4)
1285!
1286 integer(i8b), intent(out), optional :: checksum
1287!
1288 real(dp), intent(in) :: ascl
1289 real(r8), intent(out) :: amin
1290 real(r8), intent(out) :: amax
1291!
1292 character (len=*), intent(in) :: ncname
1293 character (len=*), intent(in) :: ncvname
1294!
1295# ifdef ASSUMED_SHAPE
1296# ifdef MASKING
1297 real(r8), intent(in) :: amask(lbi:,lbj:)
1298# endif
1299 real(r8), intent(out) :: adat(lbi:,lbj:)
1300# else
1301# ifdef MASKING
1302 real(r8), intent(in) :: amask(lbi:ubi,lbj:ubj)
1303# endif
1304 real(r8), intent(out) :: adat(lbi:ubi,lbj:ubj)
1305# endif
1306!
1307 TYPE (file_desc_t), intent(inout) :: piofile
1308 TYPE (io_desc_t), intent(inout) :: piodesc
1309 TYPE (my_vardesc), intent(inout) :: piovar
1310!
1311! Local variable declarations.
1312!
1313 logical :: lchecksum, interpolate
1314
1315 logical, dimension(3) :: foundit
1316!
1317 integer :: i, j, npts, status
1318 integer :: is, ie, js, je ! global bounds
1319 integer :: imin, imax, jmin, jmax ! tile bounds
1320 integer :: ilen, jlen, ijlen
1321 integer :: cgrid, ghost, dkind, gtype
1322
1323 integer, dimension(3) :: start, total
1324!
1325 real(r8) :: afactor, aoffset, aspval, avalue
1326 real(r8) :: my_amin, my_amax
1327
1328 real(r8), dimension(3) :: attvalue
1329 real(r8), dimension(2) :: rbuffer
1330!
1331 real(r4), pointer :: awrk4(:,:) ! single precision
1332 real(r8), pointer :: awrk8(:,:) ! double precision
1333 real(r8), allocatable :: cwrk(:) ! used for checksum
1334 real(r8), allocatable :: wrk(:,:) ! interpolating, regrid
1335!
1336 character (len=12), dimension(3) :: attname
1337 character (len= 3), dimension(2) :: op_handle
1338
1339 character (len=*), parameter :: myfile = &
1340 & __FILE__//", pio_fread2d"
1341!
1342!-----------------------------------------------------------------------
1343! Set starting and ending indices to process.
1344!-----------------------------------------------------------------------
1345!
1346 status=pio_noerr
1347 amin=spval
1348 amax=-spval
1349 my_amin=spval
1350 my_amax=-spval
1351!
1352! Set global (interior plus boundary) starting and ending grid cell
1353! indices in the I- and J-directions according to staggered C-grid
1354! classification.
1355!
1356 dkind=piovar%dkind
1357 gtype=piovar%gtype
1358!
1359 SELECT CASE (abs(gtype))
1360 CASE (p2dvar, p3dvar)
1361 cgrid=1 ! PSI-points
1362 is=iobounds(ng)%ILB_psi
1363 ie=iobounds(ng)%IUB_psi
1364 js=iobounds(ng)%JLB_psi
1365 je=iobounds(ng)%JUB_psi
1366 CASE (r2dvar, r3dvar)
1367 cgrid=2 ! RHO-points
1368 is=iobounds(ng)%ILB_rho
1369 ie=iobounds(ng)%IUB_rho
1370 js=iobounds(ng)%JLB_rho
1371 je=iobounds(ng)%JUB_rho
1372 CASE (u2dvar, u3dvar)
1373 cgrid=3 ! U-points
1374 is=iobounds(ng)%ILB_u
1375 ie=iobounds(ng)%IUB_u
1376 js=iobounds(ng)%JLB_u
1377 je=iobounds(ng)%JUB_u
1378 CASE (v2dvar, v3dvar)
1379 cgrid=4 ! V-points
1380 is=iobounds(ng)%ILB_v
1381 ie=iobounds(ng)%IUB_v
1382 js=iobounds(ng)%JLB_v
1383 je=iobounds(ng)%JUB_v
1384 CASE DEFAULT
1385 cgrid=2 ! RHO-points
1386 is=iobounds(ng)%ILB_rho
1387 ie=iobounds(ng)%IUB_rho
1388 js=iobounds(ng)%JLB_rho
1389 je=iobounds(ng)%JUB_rho
1390 END SELECT
1391!
1392! Compute the global lengths of the I- and J-directions. They are
1393! needed to check if regridding is required.
1394!
1395 ilen=ie-is+1
1396 jlen=je-js+1
1397!
1398! Set the tile computational I- and J-bounds (no ghost points).
1399!
1400 ghost=0
1401 imin=bounds(ng)%Imin(cgrid,ghost,myrank)
1402 imax=bounds(ng)%Imax(cgrid,ghost,myrank)
1403 jmin=bounds(ng)%Jmin(cgrid,ghost,myrank)
1404 jmax=bounds(ng)%Jmax(cgrid,ghost,myrank)
1405!
1406! Determine if interpolating from coarse gridded data to model grid
1407! is required. This is only allowed for gridded 2D fields. This is
1408! convenient for atmospheric forcing datasets that are usually on
1409! coarser grids. The user can provide coarser gridded data to avoid
1410! very large input files.
1411!
1412 interpolate=.false.
1413 IF (((vsize(1).gt.0).and.(vsize(1).ne.ilen)).or. &
1414 & ((vsize(2).gt.0).and.(vsize(2).ne.jlen))) THEN
1415 interpolate=.true.
1416 ilen=vsize(1) ! data and state variable are incongruent
1417 jlen=vsize(2) ! horizontal interpolation is required
1418 END IF
1419 IF (PRESENT(lregrid)) THEN
1420 lregrid=interpolate
1421 END IF
1422 ijlen=ilen*jlen
1423!
1424! Check if the following attributes: "scale_factor", "add_offset", and
1425! "_FillValue" are present in the input NetCDF variable:
1426!
1427! If the "scale_value" attribute is present, the data is multiplied by
1428! this factor after reading.
1429! If the "add_offset" attribute is present, this value is added to the
1430! data after reading.
1431! If both "scale_factor" and "add_offset" attributes are present, the
1432! data are first scaled before the offset is added.
1433! If the "_FillValue" attribute is present, the data having this value
1434! is treated as missing and it is replaced with zero. This feature it
1435! is usually related with the land/sea masking.
1436!
1437 attname(1)='scale_factor'
1438 attname(2)='add_offset '
1439 attname(3)='_FillValue '
1440
1441 CALL pio_netcdf_get_fatt (ng, model, ncname, piovar%vd, attname, &
1442 & attvalue, foundit, &
1443 & piofile = piofile)
1444 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1445 status=ioerror
1446 RETURN
1447 END IF
1448
1449 IF (.not.foundit(1)) THEN
1450 afactor=1.0_r8
1451 ELSE
1452 afactor=attvalue(1)
1453 END IF
1454
1455 IF (.not.foundit(2)) THEN
1456 aoffset=0.0_r8
1457 ELSE
1458 aoffset=attvalue(2)
1459 END IF
1460
1461 IF (.not.foundit(3)) THEN
1462 aspval=spval_check
1463 ELSE
1464 aspval=attvalue(3)
1465 END IF
1466!
1467!-----------------------------------------------------------------------
1468! If interpolting, read in input data with the non-tiled PIO interface.
1469!-----------------------------------------------------------------------
1470!
1471 IF (interpolate) THEN
1472 IF (.not.allocated(wrk)) THEN
1473 allocate ( wrk(ilen,jlen) )
1474 END IF
1475 wrk=0.0_r8
1476!
1477 start(1)=1
1478 total(1)=ilen
1479 start(2)=1
1480 total(2)=jlen
1481 start(3)=tindex
1482 total(3)=1
1483!
1484! The generic routine "pio_netcdf_get_fvar" checks if the attributes
1485! "scale_factor", "add_offset", and "_FillValue" are present in the
1486! input NetCDF variable and it applies its operators.
1487!
1488 CALL pio_netcdf_get_fvar (ng, model, ncname, &
1489 & ncvname, wrk, &
1490 & piofile = piofile, &
1491 & start = start, &
1492 & total = total, &
1493 & broadcast = .false., &
1494 & min_val = amin, &
1495 & max_val = amax)
1496!
1497! Multipy by user metadata scale.
1498!
1499 DO j=1,jlen
1500 DO i=1,ilen
1501 wrk(i,j)=ascl*wrk(i,j)
1502 END DO
1503 END DO
1504 END IF
1505!
1506! Initialize checsum value.
1507!
1508 IF (PRESENT(checksum)) THEN
1509 lchecksum=.true.
1510 checksum=0_i8b
1511 ELSE
1512 lchecksum=.false.
1513 END IF
1514!
1515!-----------------------------------------------------------------------
1516! If not interpolated, read in requested tiled field and scale it.
1517!-----------------------------------------------------------------------
1518!
1519 IF (.not.interpolate) THEN
1520!
1521! Allocate and initialize local array used for reading. The local array
1522! needs to be of the same precision as "A" and its IO decomposition
1523! descriptor "pioDesc".
1524!
1525 IF (dkind.eq.pio_double) THEN ! double precision
1526 IF (.not.associated(awrk8)) THEN
1527 allocate ( awrk8(imin:imax, jmin:jmax) )
1528 END IF
1529 awrk8=0.0_r8
1530 ELSE ! single precision
1531 IF (.not.associated(awrk4)) THEN
1532 allocate ( awrk4(imin:imax, jmin:jmax) )
1533 END IF
1534 awrk4=0.0_r4
1535 END IF
1536!
1537! Set unlimited time dimension record to write, if any.
1538!
1539 IF (tindex.gt.0) THEN
1540 CALL pio_setframe (piofile, &
1541 & piovar%vd, &
1542 & int(tindex, kind=pio_offset_kind))
1543 END IF
1544!
1545! Read in and load double precision data from NetCDF file.
1546!
1547 IF (dkind.eq.pio_double) THEN
1548 CALL pio_read_darray (piofile, &
1549 & piovar%vd, &
1550 & piodesc, &
1551 & awrk8(imin:,jmin:), &
1552 & status)
1553!
1554 DO j=jmin,jmax
1555 DO i=imin,imax
1556 IF (abs(awrk8(i,j)).ge.abs(aspval)) THEN
1557 adat(i,j)=0.0_r8 ! masked with _FillValue
1558 ELSE
1559 avalue=ascl*(afactor*awrk8(i,j)+aoffset)
1560 adat(i,j)=avalue
1561 my_amin=min(my_amin,avalue)
1562 my_amax=max(my_amax,avalue)
1563 END IF
1564 END DO
1565 END DO
1566 IF (associated(awrk8)) deallocate (awrk8)
1567!
1568! Read in and load single precision data from NetCDF file.
1569!
1570 ELSE
1571 CALL pio_read_darray (piofile, &
1572 & piovar%vd, &
1573 & piodesc, &
1574 & awrk4(imin:,jmin:), &
1575 & status)
1576!
1577 DO j=jmin,jmax
1578 DO i=imin,imax
1579 IF (abs(awrk4(i,j)).ge.abs(real(aspval,r4))) THEN
1580 adat(i,j)=0.0_r8 ! masked with _FillValue
1581 ELSE
1582 avalue=real(ascl*(afactor*awrk4(i,j)+aoffset),r8)
1583 adat(i,j)=avalue
1584 my_amin=real(min(my_amin,avalue),r8)
1585 my_amax=real(max(my_amax,avalue),r8)
1586 END IF
1587 END DO
1588 END DO
1589 IF (associated(awrk4)) deallocate (awrk4)
1590 END IF
1591!
1592! Compute global minimum and maximum values.
1593!
1594 rbuffer(1)=my_amin
1595 rbuffer(2)=my_amax
1596 op_handle(1)='MIN'
1597 op_handle(2)='MAX'
1598 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
1599 amin=rbuffer(1)
1600 amax=rbuffer(2)
1601!
1602 IF ((abs(amin).ge.abs(spval)).and. &
1603 & (abs(amax).ge.abs(spval))) THEN
1604 amin=0.0_r8 ! the entire data is all
1605 amax=0.0_r8 ! field value, _FillValue
1606 END IF ! and was zeroth out
1607 END IF
1608!
1609!-----------------------------------------------------------------------
1610! If interpolating from gridded data, read its associated locations
1611! and interpolate.
1612!-----------------------------------------------------------------------
1613!
1614 IF (interpolate) THEN
1615 SELECT CASE (gtype)
1616 CASE (p2dvar, p3dvar)
1617 IF (spherical) THEN
1618 CALL regrid_pio (ng, model, ncname, piofile, &
1619 & ncvname, piovar, gtype, interpflag, &
1620 & ilen, jlen, wrk, amin, amax, &
1621 & lbi, ubi, lbj, ubj, &
1622 & imin, imax, jmin, jmax, &
1623# ifdef MASKING
1624 & amask, &
1625# endif
1626 & grid(ng) % MyLon, &
1627 & grid(ng) % lonp, &
1628 & grid(ng) % latp, &
1629 & adat)
1630 ELSE
1631 CALL regrid_pio (ng, model, ncname, piofile, &
1632 & ncvname, piovar, gtype, interpflag, &
1633 & ilen, jlen, wrk, amin, amax, &
1634 & lbi, ubi, lbj, ubj, &
1635 & imin, imax, jmin, jmax, &
1636# ifdef MASKING
1637 & amask, &
1638# endif
1639 & grid(ng) % MyLon, &
1640 & grid(ng) % xp, &
1641 & grid(ng) % yp, &
1642 & adat)
1643 END IF
1644 CASE (r2dvar, r3dvar)
1645 IF (spherical) THEN
1646 CALL regrid_pio (ng, model, ncname, piofile, &
1647 & ncvname, piovar, gtype, interpflag, &
1648 & ilen, jlen, wrk, amin, amax, &
1649 & lbi, ubi, lbj, ubj, &
1650 & imin, imax, jmin, jmax, &
1651# ifdef MASKING
1652 & grid(ng) % rmask, &
1653# endif
1654 & grid(ng) % MyLon, &
1655 & grid(ng) % lonr, &
1656 & grid(ng) % latr, &
1657 & adat)
1658 ELSE
1659 CALL regrid_pio (ng, model, ncname, piofile, &
1660 & ncvname, piovar, gtype, interpflag, &
1661 & ilen, jlen, wrk, amin, amax, &
1662 & lbi, ubi, lbj, ubj, &
1663 & imin, imax, jmin, jmax, &
1664# ifdef MASKING
1665 & grid(ng) % rmask, &
1666# endif
1667 & grid(ng) % MyLon, &
1668 & grid(ng) % xr, &
1669 & grid(ng) % yr, &
1670 & adat)
1671 END IF
1672 CASE (u2dvar, u3dvar)
1673 IF (spherical) THEN
1674 CALL regrid_pio (ng, model, ncname, piofile, &
1675 & ncvname, piovar, gtype, interpflag, &
1676 & ilen, jlen, wrk, amin, amax, &
1677 & lbi, ubi, lbj, ubj, &
1678 & imin, imax, jmin, jmax, &
1679# ifdef MASKING
1680 & grid(ng) % umask, &
1681# endif
1682 & grid(ng) % MyLon, &
1683 & grid(ng) % lonu, &
1684 & grid(ng) % latu, &
1685 & adat)
1686 ELSE
1687 CALL regrid_pio (ng, model, ncname, piofile, &
1688 & ncvname, piovar, gtype, interpflag, &
1689 & ilen, jlen, wrk, amin, amax, &
1690 & lbi, ubi, lbj, ubj, &
1691 & imin, imax, jmin, jmax, &
1692# ifdef MASKING
1693 & grid(ng) % umask, &
1694# endif
1695 & grid(ng) % MyLon, &
1696 & grid(ng) % xu, &
1697 & grid(ng) % yu, &
1698 & adat)
1699 END IF
1700 CASE (v2dvar, v3dvar)
1701 IF (spherical) THEN
1702 CALL regrid_pio (ng, model, ncname, piofile, &
1703 & ncvname, piovar, gtype, interpflag, &
1704 & ilen, jlen, wrk, amin, amax, &
1705 & lbi, ubi, lbj, ubj, &
1706 & imin, imax, jmin, jmax, &
1707# ifdef MASKING
1708 & grid(ng) % vmask, &
1709# endif
1710 & grid(ng) % MyLon, &
1711 & grid(ng) % lonv, &
1712 & grid(ng) % latv, &
1713 & adat)
1714 ELSE
1715 CALL regrid_pio (ng, model, ncname, piofile, &
1716 & ncvname, piovar, gtype, interpflag, &
1717 & ilen, jlen, wrk, amin, amax, &
1718 & lbi, ubi, lbj, ubj, &
1719 & imin, imax, jmin, jmax, &
1720# ifdef MASKING
1721 & grid(ng) % vmask, &
1722# endif
1723 & grid(ng) % MyLon, &
1724 & grid(ng) % xv, &
1725 & grid(ng) % yv, &
1726 & adat)
1727 END IF
1728 CASE DEFAULT
1729 IF (spherical) THEN
1730 CALL regrid_pio (ng, model, ncname, piofile, &
1731 & ncvname, piovar, gtype, interpflag, &
1732 & ilen, jlen, wrk, amin, amax, &
1733 & lbi, ubi, lbj, ubj, &
1734 & imin, imax, jmin, jmax, &
1735# ifdef MASKING
1736 & grid(ng) % rmask, &
1737# endif
1738 & grid(ng) % MyLon, &
1739 & grid(ng) % lonr, &
1740 & grid(ng) % latr, &
1741 & adat)
1742 ELSE
1743 CALL regrid_pio (ng, model, ncname, piofile, &
1744 & ncvname, piovar, gtype, interpflag, &
1745 & ilen, jlen, wrk, amin, amax, &
1746 & lbi, ubi, lbj, ubj, &
1747 & imin, imax, jmin, jmax, &
1748# ifdef MASKING
1749 & grid(ng) % rmask, &
1750# endif
1751 & grid(ng) % MyLon, &
1752 & grid(ng) % xr, &
1753 & grid(ng) % yr, &
1754 & adat)
1755 END IF
1756 END SELECT
1757!
1758! Deallocate regridding work array.
1759!
1760 IF (allocated(wrk)) THEN
1761 deallocate (wrk)
1762 END IF
1763 END IF
1764!
1765! If requested, compute input data checksum value.
1766!
1767 IF (lchecksum) 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)
1773 END IF
1774!
1775 RETURN
1776 END FUNCTION pio_fread2d
1777#endif
1778 END MODULE nf_fread2d_mod
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)
Definition get_bounds.F:974
subroutine, public tile_bounds_1d(ng, tile, imax, istr, iend)
Definition get_bounds.F:921
subroutine, public get_hash(a, asize, hash, lreduce)
Definition get_hash.F:72
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer ioerror
logical inpthread
integer, dimension(:), allocatable tilesize
Definition mod_param.F:705
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter r3dvar
Definition mod_param.F:721
type(t_iobounds), dimension(:), allocatable iobounds
Definition mod_param.F:282
integer nghostpoints
Definition mod_param.F:710
integer, parameter iadm
Definition mod_param.F:665
integer, parameter u3dvar
Definition mod_param.F:722
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, parameter u2dvar
Definition mod_param.F:718
integer, parameter p2dvar
Definition mod_param.F:716
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 p3dvar
Definition mod_param.F:720
integer, parameter v3dvar
Definition mod_param.F:723
logical spherical
real(dp), parameter spval
real(dp), parameter spval_check
integer interpflag
integer exit_flag
type(t_scalars), dimension(:), allocatable scalars
Definition mod_scalars.F:65
integer noerror
integer function nf90_fread2d(ng, model, ncname, ncid, ncvname, ncvarid, tindex, gtype, vsize, lbi, ubi, lbj, ubj, ascl, amin, amax, amask, adat, checksum, lregrid)
Definition nf_fread2d.F:92
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)
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
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52