ROMS
Loading...
Searching...
No Matches
nf_fread2d_bry_mod Module Reference

Data Types

interface  nf_fread2d_bry
 

Functions/Subroutines

integer function nf90_fread2d_bry (ng, model, ncname, ncid, ncvname, ncvarid, tindex, gtype, lbij, ubij, nrec, ascl, amin, amax, abry, checksum)
 
integer function pio_fread2d_bry (ng, model, ncname, piofile, ncvname, piovar, tindex, piodesc, lbij, ubij, nrec, ascl, amin, amax, abry, checksum)
 

Function/Subroutine Documentation

◆ nf90_fread2d_bry()

integer function nf_fread2d_bry_mod::nf90_fread2d_bry ( integer, intent(in) ng,
integer, intent(in) model,
character (len=*), intent(in) ncname,
integer, intent(in) ncid,
character (len=*), intent(in) ncvname,
integer, intent(in) ncvarid,
integer, intent(in) tindex,
integer, intent(in) gtype,
integer, intent(in) lbij,
integer, intent(in) ubij,
integer, intent(in) nrec,
real(dp), intent(in) ascl,
real(r8), intent(out) amin,
real(r8), intent(out) amax,
real(r8), dimension(lbij:,:,:), intent(out) abry,
integer(i8b), intent(out), optional checksum )

Definition at line 80 of file nf_fread2d_bry.F.

86!***********************************************************************
87!
88 USE mod_netcdf
89!
90! Imported variable declarations.
91!
92 integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
93 integer, intent(in) :: LBij, UBij, Nrec
94!
95 integer(i8b), intent(out), optional :: checksum
96!
97 real(dp), intent(in) :: Ascl
98 real(r8), intent(out) :: Amin
99 real(r8), intent(out) :: Amax
100!
101 character (len=*), intent(in) :: ncname
102 character (len=*), intent(in) :: ncvname
103
104#ifdef ASSUMED_SHAPE
105 real(r8), intent(out) :: Abry(LBij:,:,:)
106#else
107 real(r8), intent(out) :: Abry(LBij:UBij,4,Nrec)
108#endif
109!
110! Local variable declarations.
111!
112 logical :: Lchecksum
113 logical, dimension(3) :: foundit
114 logical, dimension(4) :: bounded
115!
116 integer :: bc, ghost, i, ib, ic, ij, ir, j, tile
117 integer :: Cgrid, IorJ, Imin, Imax, Jmin, Jmax, Npts
118 integer :: Istr, Iend, Jstr, Jend
119 integer, dimension(4) :: start, total
120
121 integer :: status
122!
123 real(r8) :: Afactor, Aoffset, Aspval
124!
125 real(r8), allocatable :: Cwrk(:) ! used for checksum
126
127 real(r8), dimension(3) :: AttValue
128
129#if !defined PARALLEL_IO && defined DISTRIBUTE
130 real(r8), dimension(3) :: rbuffer
131#endif
132 real(r8), dimension(LBij:UBij,4,Nrec) :: wrk
133!
134 character (len=12), dimension(3) :: AttName
135
136 character (len=*), parameter :: MyFile = &
137 & __FILE__//", pio_fread2d_bry"
138!
139!-----------------------------------------------------------------------
140! Set starting and ending indices to process.
141!-----------------------------------------------------------------------
142!
143 status=nf90_noerr
144!
145! Set first and last grid point according to staggered C-grid
146! classification.
147!
148! Notice that (Imin,Jmin) and (Imax,Jmax) are the corner of the
149! computational tile. If ghost=0, ghost points are not processed.
150! They will be processed elsewhere by the appropriate call to any
151! of the routines in "mp_exchange.F". If ghost=1, the ghost points
152! are read.
153!
154 IF (model.eq.iadm) THEN
155 ghost=0 ! non-overlapping, no ghost points
156 ELSE
157 ghost=1 ! overlapping, read ghost points
158 END IF
159!
160 SELECT CASE (gtype)
161 CASE (p2dvar, p3dvar)
162 cgrid=1
163 CASE (r2dvar, r3dvar)
164 cgrid=2
165 CASE (u2dvar, u3dvar)
166 cgrid=3
167 CASE (v2dvar, v3dvar)
168 cgrid=4
169 CASE DEFAULT
170 cgrid=2
171 END SELECT
172!
173#ifdef DISTRIBUTE
174 tile=myrank
175 imin=bounds(ng)%Imin(cgrid,ghost,tile)
176 imax=bounds(ng)%Imax(cgrid,ghost,tile)
177 jmin=bounds(ng)%Jmin(cgrid,ghost,tile)
178 jmax=bounds(ng)%Jmax(cgrid,ghost,tile)
179#else
180 tile=-1
181 imin=lbij
182 imax=ubij
183 jmin=lbij
184 jmax=ubij
185#endif
186!
187 iorj=iobounds(ng)%IorJ
188 npts=iorj*4*nrec
189!
190! Get tile bounds.
191!
192 istr=bounds(ng)%Istr (tile)
193 iend=bounds(ng)%Iend (tile)
194 jstr=bounds(ng)%Jstr (tile)
195 jend=bounds(ng)%Jend (tile)
196!
197! Set switch to process boundary data by their associated tiles.
198!
199 bounded(iwest )=domain(ng)%Western_Edge(tile)
200 bounded(ieast )=domain(ng)%Eastern_Edge(tile)
201 bounded(isouth)=domain(ng)%Southern_Edge(tile)
202 bounded(inorth)=domain(ng)%Northern_Edge(tile)
203!
204! Set NetCDF dimension counters for processing requested field.
205!
206 start(1)=1
207 total(1)=iorj
208 start(2)=1
209 total(2)=4
210 start(3)=1
211 total(3)=nrec
212 start(4)=tindex
213 total(4)=1
214!
215! Check if the following attributes: "scale_factor", "add_offset", and
216! "_FillValue" are present in the input NetCDF variable:
217!
218! If the "scale_value" attribute is present, the data is multiplied by
219! this factor after reading.
220! If the "add_offset" attribute is present, this value is added to the
221! data after reading.
222! If both "scale_factor" and "add_offset" attributes are present, the
223! data are first scaled before the offset is added.
224! If the "_FillValue" attribute is present, the data having this value
225! is treated as missing and it is replaced with zero. This feature it
226! is usually related with the land/sea masking.
227!
228 attname(1)='scale_factor'
229 attname(2)='add_offset '
230 attname(3)='_FillValue '
231
232 CALL netcdf_get_fatt (ng, model, ncname, ncvarid, attname, &
233 & attvalue, foundit, &
234 & ncid = ncid)
235 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
236 status=ioerror
237 RETURN
238 END IF
239
240 IF (.not.foundit(1)) THEN
241 afactor=1.0_r8
242 ELSE
243 afactor=attvalue(1)
244 END IF
245
246 IF (.not.foundit(2)) THEN
247 aoffset=0.0_r8
248 ELSE
249 aoffset=attvalue(2)
250 END IF
251
252 IF (.not.foundit(3)) THEN
253 aspval=spval_check
254 ELSE
255 aspval=attvalue(3)
256 END IF
257!
258! Initialize checsum value.
259!
260 IF (PRESENT(checksum)) THEN
261 lchecksum=.true.
262 checksum=0_i8b
263 ELSE
264 lchecksum=.false.
265 END IF
266!
267!-----------------------------------------------------------------------
268! Read in requested data and scale it.
269!-----------------------------------------------------------------------
270!
271 wrk=0.0_r8
272
273 IF (inpthread) THEN
274 status=nf90_get_var(ncid, ncvarid, wrk(lbij:,:,:), start, total)
275 IF (status.eq.nf90_noerr) THEN
276 amin=spval
277 amax=-spval
278 DO ir=1,nrec
279 DO ib=1,4
280 DO ij=lbij,ubij
281 IF (abs(wrk(ij,ib,ir)).ge.abs(aspval)) THEN
282 wrk(ij,ib,ir)=0.0_r8
283 ELSE
284 wrk(ij,ib,ir)=ascl*(afactor*wrk(ij,ib,ir)+aoffset)
285 amin=min(amin,wrk(ij,ib,ir))
286 amax=max(amax,wrk(ij,ib,ir))
287 END IF
288 END DO
289 END DO
290 END DO
291 IF ((abs(amin).ge.abs(aspval)).and. &
292 & (abs(amax).ge.abs(aspval))) THEN
293 amin=0.0_r8 ! the entire data is all
294 amax=0.0_r8 ! field value, _FillValue
295 END IF
296!
297 IF (lchecksum) THEN
298 npts=(ubij-lbij+1)*nrec*4
299 IF (.not.allocated(cwrk)) allocate ( cwrk(npts) )
300 cwrk=pack(wrk(lbij:ubij, 1:4, 1:nrec), .true.)
301 CALL get_hash (cwrk, npts, checksum)
302 IF (allocated(cwrk)) deallocate (cwrk)
303 END IF
304 END IF
305 END IF
306
307#if !defined PARALLEL_IO && defined DISTRIBUTE
308!
309 rbuffer(1)=real(status,r8)
310 rbuffer(2)=amin
311 rbuffer(3)=amax
312 CALL mp_bcastf (ng, model, rbuffer)
313 status=int(rbuffer(1))
314 amin=rbuffer(2)
315 amax=rbuffer(3)
316#endif
317!
318 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
319 exit_flag=2
320 ioerror=status
321 RETURN
322 END IF
323
324#if !defined PARALLEL_IO && defined DISTRIBUTE
325!
326! Broadcast data to all spawned nodes.
327!
328 CALL mp_bcastf (ng, model, wrk)
329#endif
330!
331!-----------------------------------------------------------------------
332! Unpack read data.
333!-----------------------------------------------------------------------
334!
335 abry=0.0_r8
336
337 DO ir=1,nrec
338 DO ib=1,4
339 IF (bounded(ib)) THEN
340 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
341 DO j=jmin,jmax
342 abry(j,ib,ir)=wrk(j,ib,ir)
343 END DO
344 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
345 DO i=imin,imax
346 abry(i,ib,ir)=wrk(i,ib,ir)
347 END DO
348 END IF
349 END IF
350 END DO
351 END DO
352!
353 RETURN

◆ pio_fread2d_bry()

integer function nf_fread2d_bry_mod::pio_fread2d_bry ( integer, intent(in) ng,
integer, intent(in) model,
character (len=*), intent(in) ncname,
type (file_desc_t), intent(inout) piofile,
character (len=*), intent(in) ncvname,
type (my_vardesc), intent(inout) piovar,
integer, intent(in) tindex,
type (io_desc_t), intent(inout) piodesc,
integer, intent(in) lbij,
integer, intent(in) ubij,
integer, intent(in) nrec,
real(dp), intent(in) ascl,
real(r8), intent(out) amin,
real(r8), intent(out) amax,
real(r8), dimension(lbij:,:,:), intent(out) abry,
integer(i8b), intent(out), optional checksum )

Definition at line 359 of file nf_fread2d_bry.F.

365!***********************************************************************
366!
368!
369! Imported variable declarations.
370!
371 integer, intent(in) :: ng, model, tindex
372 integer, intent(in) :: LBij, UBij, Nrec
373!
374 integer(i8b), intent(out), optional :: checksum
375!
376 real(dp), intent(in) :: Ascl
377 real(r8), intent(out) :: Amin
378 real(r8), intent(out) :: Amax
379!
380 character (len=*), intent(in) :: ncname
381 character (len=*), intent(in) :: ncvname
382
383# ifdef ASSUMED_SHAPE
384 real(r8), intent(out) :: Abry(LBij:,:,:)
385# else
386 real(r8), intent(out) :: Abry(LBij:UBij,4,Nrec)
387# endif
388!
389 TYPE (File_desc_t), intent(inout) :: pioFile
390 TYPE (IO_Desc_t), intent(inout) :: pioDesc
391 TYPE (My_VarDesc), intent(inout) :: pioVar
392!
393! Local variable declarations.
394!
395 logical :: Lchecksum
396 logical, dimension(3) :: foundit
397 logical, dimension(4) :: bounded
398!
399 integer :: i, ib, ic, ij, ir, j
400 integer :: Cgrid, dkind, ghost, gtype, tile
401 integer :: IorJ, Imin, Imax, Jmin, Jmax, Npts
402 integer :: Istr, Iend, Jstr, Jend
403 integer, dimension(4) :: start, total
404
405 integer :: status
406!
407 real(r8) :: Afactor, Aoffset, Aspval
408!
409 real(r8), allocatable :: Cwrk(:) ! used for checksum
410
411 real(r8), dimension(3) :: AttValue
412
413 real(r8), dimension(LBij:UBij,4,Nrec) :: wrk
414!
415 character (len=12), dimension(3) :: AttName
416
417 character (len=*), parameter :: MyFile = &
418 & __FILE__//", pio_fread2d_bry"
419!
420!-----------------------------------------------------------------------
421! Set starting and ending indices to process.
422!-----------------------------------------------------------------------
423!
424 status=pio_noerr
425
426! Set first and last grid point according to staggered C-grid
427! classification.
428!
429! Notice that (Imin,Jmin) and (Imax,Jmax) are the corner of the
430! computational tile. If ghost=0, ghost points are not processed.
431! They will be processed elsewhere by the appropriate call to any
432! of the routines in "mp_exchange.F". If ghost=1, the ghost points
433! are read.
434!
435 IF (model.eq.iadm) THEN
436 ghost=0 ! non-overlapping, no ghost points
437 ELSE
438 ghost=1 ! overlapping, read ghost points
439 END IF
440 dkind=piovar%dkind
441 gtype=piovar%gtype
442!
443 SELECT CASE (gtype)
444 CASE (p2dvar, p3dvar)
445 cgrid=1
446 CASE (r2dvar, r3dvar)
447 cgrid=2
448 CASE (u2dvar, u3dvar)
449 cgrid=3
450 CASE (v2dvar, v3dvar)
451 cgrid=4
452 CASE DEFAULT
453 cgrid=2
454 END SELECT
455!
456 tile=myrank
457 imin=bounds(ng)%Imin(cgrid,ghost,tile)
458 imax=bounds(ng)%Imax(cgrid,ghost,tile)
459 jmin=bounds(ng)%Jmin(cgrid,ghost,tile)
460 jmax=bounds(ng)%Jmax(cgrid,ghost,tile)
461!
462 iorj=iobounds(ng)%IorJ
463 npts=iorj*4*nrec
464!
465! Get tile bounds.
466!
467 istr=bounds(ng)%Istr(tile)
468 iend=bounds(ng)%Iend(tile)
469 jstr=bounds(ng)%Jstr(tile)
470 jend=bounds(ng)%Jend(tile)
471!
472! Set switch to process boundary data by their associated tiles.
473!
474 bounded(iwest )=domain(ng)%Western_Edge(tile)
475 bounded(ieast )=domain(ng)%Eastern_Edge(tile)
476 bounded(isouth)=domain(ng)%Southern_Edge(tile)
477 bounded(inorth)=domain(ng)%Northern_Edge(tile)
478!
479! Set NetCDF dimension counters for processing requested field.
480!
481 start(1)=1
482 total(1)=iorj
483 start(2)=1
484 total(2)=4
485 start(3)=1
486 total(3)=nrec
487 start(4)=tindex
488 total(4)=1
489!
490! Check if the following attributes: "scale_factor", "add_offset", and
491! "_FillValue" are present in the input NetCDF variable:
492!
493! If the "scale_value" attribute is present, the data is multiplied by
494! this factor after reading.
495! If the "add_offset" attribute is present, this value is added to the
496! data after reading.
497! If both "scale_factor" and "add_offset" attributes are present, the
498! data are first scaled before the offset is added.
499! If the "_FillValue" attribute is present, the data having this value
500! is treated as missing and it is replaced with zero. This feature it
501! is usually related with the land/sea masking.
502!
503 attname(1)='scale_factor'
504 attname(2)='add_offset '
505 attname(3)='_FillValue '
506
507 CALL pio_netcdf_get_fatt (ng, model, ncname, piovar%vd, attname, &
508 & attvalue, foundit, &
509 & piofile = piofile)
510 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
511 status=ioerror
512 RETURN
513 END IF
514
515 IF (.not.foundit(1)) THEN
516 afactor=1.0_r8
517 ELSE
518 afactor=attvalue(1)
519 END IF
520
521 IF (.not.foundit(2)) THEN
522 aoffset=0.0_r8
523 ELSE
524 aoffset=attvalue(2)
525 END IF
526
527 IF (.not.foundit(3)) THEN
528 aspval=spval_check
529 ELSE
530 aspval=attvalue(3)
531 END IF
532!
533! Initialize checsum value.
534!
535 IF (PRESENT(checksum)) THEN
536 lchecksum=.true.
537 checksum=0_i8b
538 ELSE
539 lchecksum=.false.
540 END IF
541!
542!-----------------------------------------------------------------------
543! Read in 2D boundary data and scale it.
544!-----------------------------------------------------------------------
545!
546 wrk=0.0_r8
547
548 status=pio_get_var(piofile, piovar%vd, start, total, &
549 & wrk(lbij:,:,:))
550 IF (status.eq.pio_noerr) THEN
551 amin=spval
552 amax=-spval
553 DO ir=1,nrec
554 DO ib=1,4
555 DO ij=lbij,ubij
556 IF (abs(wrk(ij,ib,ir)).ge.abs(aspval)) THEN
557 wrk(ij,ib,ir)=0.0_r8
558 ELSE
559 wrk(ij,ib,ir)=ascl*(afactor*wrk(ij,ib,ir)+aoffset)
560 amin=min(amin,wrk(ij,ib,ir))
561 amax=max(amax,wrk(ij,ib,ir))
562 END IF
563 END DO
564 END DO
565 END DO
566 IF ((abs(amin).ge.abs(aspval)).and. &
567 & (abs(amax).ge.abs(aspval))) THEN
568 amin=0.0_r8 ! the entire data is all
569 amax=0.0_r8 ! field value, _FillValue
570 END IF
571!
572 IF (lchecksum) THEN
573 npts=(ubij-lbij+1)*nrec*4
574 IF (.not.allocated(cwrk)) allocate ( cwrk(npts) )
575 cwrk=pack(wrk(lbij:ubij, 1:4, 1:nrec), .true.)
576 CALL get_hash (cwrk, npts, checksum)
577 IF (allocated(cwrk)) deallocate (cwrk)
578 END IF
579 END IF
580!
581 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
582 exit_flag=2
583 ioerror=status
584 RETURN
585 END IF
586!
587!-----------------------------------------------------------------------
588! Unpack read data.
589!-----------------------------------------------------------------------
590!
591 abry=0.0_r8
592
593 DO ir=1,nrec
594 DO ib=1,4
595 IF (bounded(ib)) THEN
596 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
597 DO j=jmin,jmax
598 abry(j,ib,ir)=wrk(j,ib,ir)
599 END DO
600 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
601 DO i=imin,imax
602 abry(i,ib,ir)=wrk(i,ib,ir)
603 END DO
604 END IF
605 END IF
606 END DO
607 END DO
608!
609 RETURN