ROMS
Loading...
Searching...
No Matches
nf_fwrite2d.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 module writes out a generic floating point 2D array into !
12! an output file using either the standard NetCDF library or the !
13! Parallel-IO (PIO) library. !
14! !
15! On Input: !
16! !
17! ng Nested grid number !
18! model Calling model identifier !
19! ncid NetCDF file ID !
20! ifield Field metadata index (integer) !
21#if defined PIO_LIB && defined DISTRIBUTE
22!or pioFile PIO file descriptor structure, TYPE(file_desc_t) !
23! pioFile%fh file handler !
24! pioFile%iosystem IO system descriptor (struct) !
25#endif
26! ncvarid NetCDF variable ID !
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 write !
34! gtype Grid type. If negative, only write water points !
35#if defined PIO_LIB && defined DISTRIBUTE
36!or pioDesc IO data decomposition descriptor, TYPE(IO_desc_t) !
37#endif
38! LBi I-dimension Lower bound !
39! UBi I-dimension Upper bound !
40! LBj J-dimension Lower bound !
41! UBj J-dimension Upper bound !
42! Amask land/Sea mask, if any (real) !
43! Ascl Factor to scale field before writing (real) !
44! Adat Field to write out (real) !
45! SetFillVal Logical switch to set fill value in land areas !
46! (OPTIONAL) !
47! ExtractField Field extraction flag (integer, OPTIONAL) !
48! ExtractField = 0 no extraction !
49! ExtractField = 1 extraction by interpolation !
50! ExtractField > 1 extraction by decimation !
51! !
52! On Output: !
53! !
54! status Error flag (integer) !
55! MinValue Minimum value (real, OPTIONAL) !
56! MaxValue Maximum value (real, OPTIONAL) !
57! !
58#ifdef POSITIVE_ZERO
59! Starting F95 zero values can be signed (-0 or +0) following the !
60! IEEE 754 floating point standard. This may produce different !
61! output data in serial and parallel applications. Since comparing !
62! serial and parallel output is essential for tracking parallel !
63! partition bugs, "positive zero" is enforced. !
64! !
65#endif
66!=======================================================================
67!
68 USE mod_param
69 USE mod_parallel
70 USE mod_ncparam
71 USE mod_scalars
72!
74!
75 implicit none
76!
77 INTERFACE nf_fwrite2d
78 MODULE PROCEDURE nf90_fwrite2d
79#if defined PIO_LIB && defined DISTRIBUTE
80 MODULE PROCEDURE pio_fwrite2d
81#endif
82 END INTERFACE nf_fwrite2d
83!
84 CONTAINS
85
86#if defined PARALLEL_IO && defined DISTRIBUTE
87!
88!***********************************************************************
89 FUNCTION nf90_fwrite2d (ng, model, ncid, ifield, &
90 & ncvarid, tindex, gtype, &
91 & LBi, UBi, LBj, UBj, Ascl, &
92# ifdef MASKING
93 & Amask, &
94# endif
95 & Adat, &
96 & SetFillVal, &
97 & ExtractField, &
98 & MinValue, MaxValue) RESULT (status)
99!***********************************************************************
100!
101 USE mod_netcdf
102
103# if defined WRITE_WATER && defined MASKING
104!
105 USE distribute_mod, ONLY : mp_collect
106# endif
108!
109! Imported variable declarations.
110!
111 logical, intent(in), optional :: setfillval
112!
113 integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
114 integer, intent(in) :: ifield
115 integer, intent(in) :: lbi, ubi, lbj, ubj
116!
117 integer, intent(in), optional :: extractfield
118!
119 real(dp), intent(in) :: ascl
120!
121# ifdef ASSUMED_SHAPE
122# ifdef MASKING
123 real(r8), intent(in) :: amask(lbi:,lbj:)
124# endif
125 real(r8), intent(in) :: adat(lbi:,lbj:)
126# else
127# ifdef MASKING
128 real(r8), intent(in) :: amask(lbi:ubi,lbj:ubj)
129# endif
130 real(r8), intent(in) :: adat(lbi:ubi,lbj:ubj)
131# endif
132 real(r8), intent(out), optional :: minvalue
133 real(r8), intent(out), optional :: maxvalue
134!
135! Local variable declarations.
136!
137 logical :: landfill
138!
139 integer :: extract_flag
140 integer :: i, ic, j, jc, npts
141 integer :: imin, imax, jmin, jmax
142 integer :: ioff, joff
143 integer :: istr, iend
144 integer :: ilen, isize, jlen, jsize, ijlen
145 integer :: mytype
146 integer :: status
147
148 integer, dimension(3) :: start, total
149!
150 real(r8), parameter :: inival = 0.0_r8
151
152 real(r8), allocatable :: awrk(:)
153!
154!-----------------------------------------------------------------------
155! Set starting and ending indices to process.
156!-----------------------------------------------------------------------
157!
158 status=nf90_noerr
159!
160! Set first and last grid point according to staggered C-grid
161! classification.
162!
163 mytype=gtype
164!
165 SELECT CASE (abs(mytype))
166 CASE (p2dvar, p3dvar)
167 imin=bounds(ng)%Istr (myrank)
168 imax=bounds(ng)%IendR(myrank)
169 jmin=bounds(ng)%Jstr (myrank)
170 jmax=bounds(ng)%JendR(myrank)
171 isize=iobounds(ng)%xi_psi
172 jsize=iobounds(ng)%eta_psi
173 CASE (r2dvar, r3dvar)
174 imin=bounds(ng)%IstrR(myrank)
175 imax=bounds(ng)%IendR(myrank)
176 jmin=bounds(ng)%JstrR(myrank)
177 jmax=bounds(ng)%JendR(myrank)
178 isize=iobounds(ng)%xi_rho
179 jsize=iobounds(ng)%eta_rho
180 CASE (u2dvar, u3dvar)
181 imin=bounds(ng)%Istr (myrank)
182 imax=bounds(ng)%IendR(myrank)
183 jmin=bounds(ng)%JstrR(myrank)
184 jmax=bounds(ng)%JendR(myrank)
185 isize=iobounds(ng)%xi_u
186 jsize=iobounds(ng)%eta_u
187 CASE (v2dvar, v3dvar)
188 imin=bounds(ng)%IstrR(myrank)
189 imax=bounds(ng)%IendR(myrank)
190 jmin=bounds(ng)%Jstr (myrank)
191 jmax=bounds(ng)%JendR(myrank)
192 isize=iobounds(ng)%xi_v
193 jsize=iobounds(ng)%eta_v
194 CASE DEFAULT
195 imin=bounds(ng)%IstrR(myrank)
196 imax=bounds(ng)%IendR(myrank)
197 jmin=bounds(ng)%JstrR(myrank)
198 jmax=bounds(ng)%JendR(myrank)
199 isize=iobounds(ng)%xi_rho
200 jsize=iobounds(ng)%eta_rho
201 END SELECT
202!
203 ilen=imax-imin+1
204 jlen=jmax-jmin+1
205 ijlen=ilen*jlen
206!
207! Set switch to replace land areas with fill value, spval.
208!
209# ifdef MASKING
210 IF (PRESENT(setfillval)) THEN
211 landfill=setfillval
212 ELSE
213 landfill=tindex.gt.0
214 END IF
215# else
216 landfill=.false.
217# endif
218!
219! If appropriate, set the field extraction flag to the provided grid
220! geometry through interpolation or decimation.
221!
222 IF (PRESENT(extractfield)) THEN
223 extract_flag=extractfield
224 ELSE
225 extract_flag=0
226 END IF
227!
228!-----------------------------------------------------------------------
229! Parallel I/O: Pack tile data into 1D array in column-major order.
230# ifdef MASKING
231! Overwrite masked points with special value.
232# endif
233!-----------------------------------------------------------------------
234!
235 IF (gtype.gt.0) THEN
236!
237! Set offsets due the NetCDF dimensions. Recall that some output
238! variables not always start at one.
239!
240 SELECT CASE (abs(mytype))
241 CASE (p2dvar, p3dvar)
242 ioff=0
243 joff=0
244 CASE (r2dvar, r3dvar)
245 ioff=1
246 joff=1
247 CASE (u2dvar, u3dvar)
248 ioff=0
249 joff=1
250 CASE (v2dvar, v3dvar)
251 ioff=1
252 joff=0
253 CASE DEFAULT
254 ioff=1
255 joff=1
256 END SELECT
257!
258! Allocate and initialize scratch work array.
259!
260 npts=ijlen
261 IF (.not.allocated(awrk)) THEN
262 allocate ( awrk(npts) )
263 awrk=inival
264 END IF
265!
266! Pack and scale tile data.
267!
268 ic=0
269 DO j=jmin,jmax
270 DO i=imin,imax
271 ic=ic+1
272 awrk(ic)=adat(i,j)*ascl
273#ifdef POSITIVE_ZERO
274 IF (abs(awrk(ic)).eq.0.0_r8) THEN
275 awrk(ic)=0.0_r8 ! impose positive zero
276 END IF
277#endif
278# ifdef MASKING
279 IF ((amask(i,j).eq.0.0_r8).and.landfill) THEN
280 awrk(ic)=spval
281 END IF
282# endif
283 END DO
284 END DO
285!
286! Write out data: all parallel nodes write their own packed tile data.
287!
288 start(1)=imin+ioff
289 total(1)=ilen
290 start(2)=jmin+joff
291 total(2)=jlen
292 start(3)=tindex
293 total(3)=1
294
295 status=nf90_put_var(ncid, ncvarid, awrk, start, total)
296 END IF
297
298# if defined WRITE_WATER && defined MASKING
299!
300!-----------------------------------------------------------------------
301! Parallel I/O: Remove land points and pack tile data into 1D array.
302!-----------------------------------------------------------------------
303!
304 IF (gtype.lt.0) THEN
305!
306! Set offsets due array packing into 1D array in column-major order.
307!
308 SELECT CASE (abs(mytype))
309 CASE (p2dvar, p3dvar)
310 ioff=0
311 joff=1
312 CASE (r2dvar, r3dvar)
313 ioff=1
314 joff=0
315 CASE (u2dvar, u3dvar)
316 ioff=0
317 joff=0
318 CASE (v2dvar, v3dvar)
319 ioff=1
320 joff=1
321 CASE DEFAULT
322 ioff=1
323 joff=0
324 END SELECT
325!
326! Allocate and initialize scratch work array.
327!
328 npts=isize*jsize
329 IF (.not.allocated(awrk)) THEN
330 allocate ( awrk(npts) )
331 awrk=inival
332 END IF
333!
334! Scale and gather data from all spawned nodes. Store data into a 1D
335! global array, packed in column-major order. Flag land point with
336! special value.
337!
338 DO j=jmin,jmax
339 jc=(j-joff)*isize
340 DO i=imin,imax
341 ic=i+ioff+jc
342 awrk(ic)=adat(i,j)*ascl
343#ifdef POSITIVE_ZERO
344 IF (abs(awrk(ic)).eq.0.0_r8) THEN
345 awrk(ic)=0.0_r8 ! impose positive zero
346 END IF
347#endif
348 IF (amask(i,j).eq.0.0_r8) THEN
349 awrk(ic)=spval
350 END IF
351 END DO
352 END DO
353!
354! Global reduction of work array.
355!
356 CALL mp_collect (ng, model, npts, inival, awrk)
357!
358! Remove land points.
359!
360 ic=0
361 DO i=1,npts
362 IF (awrk(i).lt.spval) THEN
363 ic=ic+1
364 awrk(ic)=awrk(i)
365 END IF
366 END DO
367 npts=ic
368!
369! Write out data: all parallel nodes write a section of the packed
370! data.
371!
372 CALL tile_bounds_1d (ng, myrank, npts, istr, iend)
373
374 start(1)=istr
375 total(1)=iend-istr+1
376 start(2)=tindex
377 total(2)=1
378
379 status=nf90_put_var(ncid, ncvarid, awrk(istr:), start, total)
380 END IF
381# endif
382!
383!-----------------------------------------------------------------------
384! Deallocate scratch work array.
385!-----------------------------------------------------------------------
386!
387 IF (allocated(awrk)) THEN
388 deallocate (awrk)
389 END IF
390
391 RETURN
392 END FUNCTION nf90_fwrite2d
393
394#else
395
396!
397!***********************************************************************
398 FUNCTION nf90_fwrite2d (ng, model, ncid, ifield, &
399 & ncvarid, tindex, gtype, &
400 & LBi, UBi, LBj, UBj, Ascl, &
401# ifdef MASKING
402 & Amask, &
403# endif
404 & Adat, &
405 & SetFillVal, &
406 & ExtractField, &
407 & MinValue, MaxValue) RESULT (status)
408!***********************************************************************
409!
410 USE mod_netcdf
411
412# ifdef DISTRIBUTE
413!
415# endif
416# ifdef OUTPUT_STATS
417 USE stats_mod, ONLY : stats_2dfld
418# endif
419!
420! Imported variable declarations.
421!
422 logical, intent(in), optional :: setfillval
423!
424 integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
425 integer, intent(in) :: ifield
426 integer, intent(in) :: lbi, ubi, lbj, ubj
427!
428 integer, intent(in), optional :: extractfield
429!
430 real(dp), intent(in) :: ascl
431!
432# ifdef ASSUMED_SHAPE
433# ifdef MASKING
434 real(r8), intent(in) :: amask(lbi:,lbj:)
435# endif
436 real(r8), intent(in) :: adat(lbi:,lbj:)
437# else
438# ifdef MASKING
439 real(r8), intent(in) :: amask(lbi:ubi,lbj:ubj)
440# endif
441 real(r8), intent(in) :: adat(lbi:ubi,lbj:ubj)
442# endif
443 real(r8), intent(out), optional :: minvalue
444 real(r8), intent(out), optional :: maxvalue
445!
446! Local variable declarations.
447!
448 logical :: landfill
449!
450 integer :: extract_flag
451 integer :: npts, i, tile
452 integer :: status
453
454 integer, dimension(3) :: start, total
455!
456 real(r8), dimension((Lm(ng)+2)*(Mm(ng)+2)) :: awrk
457
458# ifdef OUTPUT_STATS
459!
460 TYPE (t_stats) :: stats
461# endif
462!
463!-----------------------------------------------------------------------
464! Initialize local variables.
465!-----------------------------------------------------------------------
466!
467 status=nf90_noerr
468!
469! Set parallel tile.
470!
471# ifdef DISTRIBUTE
472 tile=myrank
473# else
474 tile=-1
475# endif
476!
477! Set switch to replace land areas with fill value, spval.
478!
479# ifdef MASKING
480 IF (PRESENT(setfillval)) THEN
481 landfill=setfillval
482 ELSE
483 landfill=tindex.gt.0
484 END IF
485# else
486 landfill=.false.
487# endif
488!
489! If appropriate, set the field extraction flag to the provided grid
490! geometry through interpolation or decimation.
491!
492 IF (PRESENT(extractfield)) THEN
493 extract_flag=extractfield
494 ELSE
495 extract_flag=0
496 END IF
497
498# ifdef OUTPUT_STATS
499!
500! Initialize output variable statistics structure.
501!
502 stats % checksum=0_i8b
503 stats % count=0
504 stats % min=spval
505 stats % max=-spval
506 stats % avg=0.0_r8
507 stats % rms=0.0_r8
508# endif
509!
510! Initialize local array to avoid denormalized numbers. This
511! facilitates processing and debugging.
512!
513 awrk=0.0_r8
514!
515!-----------------------------------------------------------------------
516! Pack 2D field data into 1D array.
517!-----------------------------------------------------------------------
518!
519 CALL pack_field2d (ng, model, tile, &
520 & gtype, ifield, tindex, &
521 & landfill, extract_flag, &
522 & lbi, ubi, lbj, ubj, &
523# ifdef MASKING
524 & amask, &
525# endif
526 & ascl, adat, &
527 & start, total, npts, awrk)
528!
529!-----------------------------------------------------------------------
530! If applicable, compute output field minimum and maximum values.
531!-----------------------------------------------------------------------
532!
533 IF (PRESENT(minvalue)) THEN
534 IF (outthread) THEN
535 minvalue=spval
536 maxvalue=-spval
537 DO i=1,npts
538 IF (abs(awrk(i)).lt.spval) THEN
539 minvalue=min(minvalue,awrk(i))
540 maxvalue=max(maxvalue,awrk(i))
541 END IF
542 END DO
543 END IF
544 END IF
545!
546!-----------------------------------------------------------------------
547! Write output buffer into NetCDF file.
548!-----------------------------------------------------------------------
549!
550 IF (outthread) THEN
551 status=nf90_put_var(ncid, ncvarid, awrk, start, total)
552 END IF
553
554# ifdef OUTPUT_STATS
555!
556!-----------------------------------------------------------------------
557! Compute and report output field statistics.
558!-----------------------------------------------------------------------
559!
560 CALL stats_2dfld (ng, tile, model, gtype, stats, &
561 & extract_flag, &
562 & lbi, ubi, lbj, ubj, &
563 & adat, &
564# ifdef MASKING
565 & fmask = amask, &
566# endif
567 & debug = .false.)
568 IF (outthread) THEN
569 WRITE (stdout,10) trim(vname(1,ifield)), stats%min, stats%max, &
570 & stats%checksum
571 10 FORMAT (4x,'- ',a,':',t35,'Min = ',1p,e15.8,', Max = ', &
572 & 1p,e15.8,', CheckSum = ',i0)
573 END IF
574# endif
575# ifdef DISTRIBUTE
576!
577!-----------------------------------------------------------------------
578! Broadcast IO error flag to all nodes.
579!-----------------------------------------------------------------------
580!
581 CALL mp_bcasti (ng, model, status)
582# endif
583!
584 RETURN
585 END FUNCTION nf90_fwrite2d
586#endif
587
588#if defined PIO_LIB && defined DISTRIBUTE
589!
590!***********************************************************************
591 FUNCTION pio_fwrite2d (ng, model, pioFile, ifield, &
592 & pioVar, tindex, pioDesc, &
593 & LBi, UBi, LBj, UBj, Ascl, &
594# ifdef MASKING
595 & Amask, &
596# endif
597 & Adat, &
598 & SetFillVal, &
599 & ExtractField, &
600 & MinValue, MaxValue) RESULT (status)
601!***********************************************************************
602!
604!
605 USE distribute_mod, ONLY : mp_reduce
606# ifdef OUTPUT_STATS
607 USE stats_mod, ONLY : stats_2dfld
608# endif
609!
610! Imported variable declarations.
611!
612 logical, intent(in), optional :: setfillval
613!
614 integer, intent(in) :: ng, model, tindex
615 integer, intent(in) :: ifield
616 integer, intent(in) :: lbi, ubi, lbj, ubj
617!
618 integer, intent(in), optional :: extractfield
619!
620 real(dp), intent(in) :: ascl
621!
622# ifdef ASSUMED_SHAPE
623# ifdef MASKING
624 real(r8), intent(in) :: amask(lbi:,lbj:)
625# endif
626 real(r8), intent(in) :: adat(lbi:,lbj:)
627# else
628# ifdef MASKING
629 real(r8), intent(in) :: amask(lbi:ubi,lbj:ubj)
630# endif
631 real(r8), intent(in) :: adat(lbi:ubi,lbj:ubj)
632# endif
633 real(r8), intent(out), optional :: minvalue
634 real(r8), intent(out), optional :: maxvalue
635!
636 TYPE (file_desc_t), intent(inout) :: piofile
637 TYPE (io_desc_t), intent(inout) :: piodesc
638 TYPE (my_vardesc), intent(inout) :: piovar
639!
640! Local variable declarations.
641!
642 logical :: landfill, lminmax
643
644 logical, pointer :: lwater(:,:)
645!
646 integer :: extract_flag
647 integer :: i, j, tile
648 integer :: imin, imax, jmin, jmax
649 integer :: cgrid, dkind, ghost, gtype
650 integer :: status
651
652 integer, dimension(3) :: start, total
653!
654 real(r8), dimension(2) :: rbuffer
655
656 real(r4), pointer :: awrk4(:,:)
657 real(r8), pointer :: awrk8(:,:)
658
659# ifdef OUTPUT_STATS
660!
661 TYPE (t_stats) :: stats
662# endif
663!
664 character (len= 3), dimension(2) :: op_handle
665!
666!-----------------------------------------------------------------------
667! Set starting and ending indices to process.
668!-----------------------------------------------------------------------
669!
670 status=pio_noerr
671!
672! Set first and last tile computational grid point according to the
673! staggered C-grid location. Ghost points are not included.
674!
675 ghost=0
676 dkind=piovar%dkind
677 gtype=piovar%gtype
678!
679 SELECT CASE (gtype)
680 CASE (p2dvar, p3dvar)
681 cgrid=1 ! PSI-points
682 CASE (r2dvar, r3dvar)
683 cgrid=2 ! RHO-points
684 CASE (u2dvar, u3dvar)
685 cgrid=3 ! U-points
686 CASE (v2dvar, v3dvar)
687 cgrid=4 ! V-points
688 CASE DEFAULT
689 cgrid=2 ! RHO-points
690 END SELECT
691!
692 imin=bounds(ng)%Imin(cgrid,ghost,myrank)
693 imax=bounds(ng)%Imax(cgrid,ghost,myrank)
694 jmin=bounds(ng)%Jmin(cgrid,ghost,myrank)
695 jmax=bounds(ng)%Jmax(cgrid,ghost,myrank)
696!
697! Set switch to compute minimum and maximum values.
698!
699 IF (PRESENT(minvalue)) THEN
700 lminmax=.true.
701 IF (.not.associated(lwater)) THEN
702 allocate ( lwater(lbi:ubi,lbj:ubj) )
703 lwater=.true.
704 END IF
705 ELSE
706 lminmax=.false.
707 END IF
708!
709! Set switch to replace land areas with fill value, spval.
710!
711# ifdef MASKING
712 IF (PRESENT(setfillval)) THEN
713 landfill=setfillval
714 ELSE
715 landfill=tindex.gt.0
716 END IF
717# else
718 landfill=.false.
719# endif
720!
721! If appropriate, set the field extraction flag to the provided grid
722! geometry through interpolation or decimation.
723!
724 IF (PRESENT(extractfield)) THEN
725 extract_flag=extractfield
726 ELSE
727 extract_flag=0
728 END IF
729
730# ifdef OUTPUT_STATS
731!
732! Initialize output variable statistics structure.
733!
734 stats % checksum=0_i8b
735 stats % count=0
736 stats % min=spval
737 stats % max=-spval
738 stats % avg=0.0_r8
739 stats % rms=0.0_r8
740# endif
741!
742!-----------------------------------------------------------------------
743! Write out data into NetCDF file.
744!-----------------------------------------------------------------------
745!
746! Allocate, initialize and load data into local array used for
747! writing. Overwrite masked points with special value.
748!
749 IF (dkind.eq.pio_double) THEN ! double precision
750 IF (.not.associated(awrk8)) THEN
751 allocate ( awrk8(lbi:ubi,lbj:ubj) )
752 awrk8=0.0_r8
753 END IF
754!
755 DO j=jmin,jmax
756 DO i=imin,imax
757 awrk8(i,j)=adat(i,j)*ascl
758# ifdef MASKING
759 IF((amask(i,j).eq.0.0_r8).and.landfill) THEN
760 awrk8(i,j)=spval
761 IF (lminmax) lwater(i,j)=.false.
762 END iF
763# endif
764 END DO
765 END DO
766 IF (lminmax) THEN
767 rbuffer(1)=minval(awrk8, mask=lwater)
768 rbuffer(2)=maxval(awrk8, mask=lwater)
769 END IF
770 ELSE ! single precision
771 IF (.not.associated(awrk4)) THEN
772 allocate ( awrk4(lbi:ubi,lbj:ubj) )
773 awrk4=0.0_r4
774 END IF
775!
776 DO j=jmin,jmax
777 DO i=imin,imax
778 awrk4(i,j)=real(adat(i,j)*ascl, r4)
779
780# ifdef MASKING
781 IF((amask(i,j).eq.0.0_r8).and.landfill) THEN
782 awrk4(i,j)=real(spval, r4)
783 IF (lminmax) lwater(i,j)=.false.
784 END iF
785# endif
786 END DO
787 END DO
788 IF (lminmax) THEN
789 rbuffer(1)=real(minval(awrk4, mask=lwater),r8)
790 rbuffer(2)=real(maxval(awrk4, mask=lwater),r8)
791 END IF
792 END IF
793!
794! Set unlimited time dimension record to write, if any.
795!
796 IF (tindex.gt.0) THEN
797 CALL pio_setframe (piofile, &
798 & piovar%vd, &
799 & int(tindex, kind=pio_offset_kind))
800 END IF
801!
802! Write out data into NetCDF.
803!
804 IF (dkind.eq.pio_double) THEN ! double precision
805 CALL pio_write_darray (piofile, &
806 & piovar%vd, &
807 & piodesc, &
808 & awrk8(imin:imax,jmin:jmax), &
809 & status)
810 ELSE ! single precision
811 CALL pio_write_darray (piofile, &
812 & piovar%vd, &
813 & piodesc, &
814 & awrk4(imin:imax,jmin:jmax), &
815 & status)
816 END IF
817
818# ifdef OUTPUT_STATS
819!
820!-----------------------------------------------------------------------
821! Compute and report output field statistics.
822!-----------------------------------------------------------------------
823!
824# ifdef DISTRIBUTE
825 tile=myrank
826# else
827 tile=-1
828# endif
829 CALL stats_2dfld (ng, tile, model, gtype, stats, &
830 & extract_flag, &
831 & lbi, ubi, lbj, ubj, &
832 & adat, &
833# ifdef MASKING
834 & fmask = amask, &
835# endif
836 & debug = .false.)
837 IF (outthread) THEN
838 WRITE (stdout,10) trim(vname(1,ifield)), stats%min, stats%max, &
839 & stats%checksum
840 10 FORMAT (4x,'- ',a,':',t35,'Min = ',1p,e15.8,', Max = ', &
841 & 1p,e15.8,', CheckSum = ',i0)
842 END IF
843# endif
844!
845!-----------------------------------------------------------------------
846! If applicable, compute global minimum and maximum values.
847!-----------------------------------------------------------------------
848!
849 IF (lminmax) THEN
850 op_handle(1)='MIN'
851 op_handle(2)='MAX'
852 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
853 minvalue=rbuffer(1)
854 maxvalue=rbuffer(2)
855 IF (associated(lwater)) deallocate (lwater)
856 END IF
857!
858! Deallocate local array.
859!
860 IF (dkind.eq.pio_double) THEN
861 IF (associated(awrk8)) deallocate (awrk8)
862 ELSE
863 IF (associated(awrk4)) deallocate (awrk4)
864 END IF
865!
866 RETURN
867 END FUNCTION pio_fwrite2d
868#endif
869!
870 END MODULE nf_fwrite2d_mod
subroutine mp_gather2d(ng, model, lbi, ubi, lbj, ubj, tindex, gtype, ascl, amask, a, npts, awrk, setfillval)
subroutine, public tile_bounds_1d(ng, tile, imax, istr, iend)
Definition get_bounds.F:921
character(len=maxlen), dimension(6, 0:nv) vname
logical outthread
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, 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(dp), parameter spval
integer function pio_fwrite2d(ng, model, piofile, ifield, piovar, tindex, piodesc, lbi, ubi, lbj, ubj, ascl, amask, adat, setfillval, extractfield, minvalue, maxvalue)
integer function nf90_fwrite2d(ng, model, ncid, ifield, ncvarid, tindex, gtype, lbi, ubi, lbj, ubj, ascl, amask, adat, setfillval, extractfield, minvalue, maxvalue)
Definition nf_fwrite2d.F:99
subroutine, public pack_field2d(ng, model, tile, gtype, ifield, tindex, landfill, extract_flag, lbi, ubi, lbj, ubj, amask, ascl, adat, start, total, npts, awrk)
Definition pack_field.F:513
subroutine, public stats_2dfld(ng, tile, model, gtype, s, extract_flag, lbi, ubi, lbj, ubj, f, fmask, debug)
Definition stats.F:47