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