source: palm/trunk/SOURCE/data_output_3d.f90 @ 2723

Last change on this file since 2723 was 2718, checked in by maronga, 7 years ago

deleting of deprecated files; headers updated where needed

  • Property svn:keywords set to Id
File size: 27.0 KB
Line 
1!> @file data_output_3d.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2018 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: data_output_3d.f90 2718 2018-01-02 08:49:38Z maronga $
27! Corrected "Former revisions" section
28!
29! 2696 2017-12-14 17:12:51Z kanani
30! Change in file header (GPL part)
31! Implementation of turbulence_closure_mod (TG)
32! Implementation of chemistry module (FK)
33! Set fill values at topography grid points or e.g. non-natural-type surface
34! in case of LSM output (MS)
35!
36! 2512 2017-10-04 08:26:59Z raasch
37! upper bounds of 3d output changed from nx+1,ny+1 to nx,ny
38! no output of ghost layer data
39!
40! 2292 2017-06-20 09:51:42Z schwenkel
41! Implementation of new microphysic scheme: cloud_scheme = 'morrison'
42! includes two more prognostic equations for cloud drop concentration (nc) 
43! and cloud water content (qc).
44!
45! 2233 2017-05-30 18:08:54Z suehring
46!
47! 2232 2017-05-30 17:47:52Z suehring
48! Adjustments to new topography concept
49!
50! 2209 2017-04-19 09:34:46Z kanani
51! Added plant canopy model output
52!
53! 2031 2016-10-21 15:11:58Z knoop
54! renamed variable rho to rho_ocean and rho_av to rho_ocean_av
55!
56! 2011 2016-09-19 17:29:57Z kanani
57! Flag urban_surface is now defined in module control_parameters,
58! changed prefix for urban surface model output to "usm_",
59! introduced control parameter varnamelength for LEN of trimvar.
60!
61! 2007 2016-08-24 15:47:17Z kanani
62! Added support for new urban surface model (temporary modifications of
63! SELECT CASE ( ) necessary, see variable trimvar)
64!
65! 2000 2016-08-20 18:09:15Z knoop
66! Forced header and separation lines into 80 columns
67!
68! 1980 2016-07-29 15:51:57Z suehring
69! Bugfix, in order to steer user-defined output, setting flag found explicitly
70! to .F.
71!
72! 1976 2016-07-27 13:28:04Z maronga
73! Output of radiation quantities is now done directly in the respective module
74!
75! 1972 2016-07-26 07:52:02Z maronga
76! Output of land surface quantities is now done directly in the respective module.
77! Unnecessary directive __parallel removed.
78!
79! 1960 2016-07-12 16:34:24Z suehring
80! Scalar surface flux added
81!
82! 1849 2016-04-08 11:33:18Z hoffmann
83! prr moved to arrays_3d
84!
85! 1822 2016-04-07 07:49:42Z hoffmann
86! prr vertical dimensions set to nzb_do to nzt_do. Unused variables deleted.
87!
88! 1808 2016-04-05 19:44:00Z raasch
89! test output removed
90!
91! 1783 2016-03-06 18:36:17Z raasch
92! name change of netcdf routines and module + related changes
93!
94! 1745 2016-02-05 13:06:51Z gronemeier
95! Bugfix: test if time axis limit exceeds moved to point after call of check_open
96!
97! 1691 2015-10-26 16:17:44Z maronga
98! Added output of radiative heating rates for RRTMG
99!
100! 1682 2015-10-07 23:56:08Z knoop
101! Code annotations made doxygen readable
102!
103! 1585 2015-04-30 07:05:52Z maronga
104! Added support for RRTMG
105!
106! 1551 2015-03-03 14:18:16Z maronga
107! Added suppport for land surface model and radiation model output. In the course
108! of this action, the limits for vertical loops have been changed (from nzb and
109! nzt+1 to nzb_do and nzt_do, respectively in order to allow soil model output).
110! Moreover, a new vertical grid zs was introduced.
111!
112! 1359 2014-04-11 17:15:14Z hoffmann
113! New particle structure integrated.
114!
115! 1353 2014-04-08 15:21:23Z heinze
116! REAL constants provided with KIND-attribute
117!
118! 1327 2014-03-21 11:00:16Z raasch
119! parts concerning avs output removed,
120! -netcdf output queries
121!
122! 1320 2014-03-20 08:40:49Z raasch
123! ONLY-attribute added to USE-statements,
124! kind-parameters added to all INTEGER and REAL declaration statements,
125! kinds are defined in new module kinds,
126! old module precision_kind is removed,
127! revision history before 2012 removed,
128! comment fields (!:) to be used for variable explanations added to
129! all variable declaration statements
130!
131! 1318 2014-03-17 13:35:16Z raasch
132! barrier argument removed from cpu_log,
133! module interfaces removed
134!
135! 1308 2014-03-13 14:58:42Z fricke
136! Check, if the limit of the time dimension is exceeded for parallel output
137! To increase the performance for parallel output, the following is done:
138! - Update of time axis is only done by PE0
139!
140! 1244 2013-10-31 08:16:56Z raasch
141! Bugfix for index bounds in case of 3d-parallel output
142!
143! 1115 2013-03-26 18:16:16Z hoffmann
144! ql is calculated by calc_liquid_water_content
145!
146! 1106 2013-03-04 05:31:38Z raasch
147! array_kind renamed precision_kind
148!
149! 1076 2012-12-05 08:30:18Z hoffmann
150! Bugfix in output of ql
151!
152! 1053 2012-11-13 17:11:03Z hoffmann
153! +nr, qr, prr, qc and averaged quantities
154!
155! 1036 2012-10-22 13:43:42Z raasch
156! code put under GPL (PALM 3.9)
157!
158! 1031 2012-10-19 14:35:30Z raasch
159! netCDF4 without parallel file support implemented
160!
161! 1007 2012-09-19 14:30:36Z franke
162! Bugfix: missing calculation of ql_vp added
163!
164! Revision 1.1  1997/09/03 06:29:36  raasch
165! Initial revision
166!
167!
168! Description:
169! ------------
170!> Output of the 3D-arrays in netCDF and/or AVS format.
171!------------------------------------------------------------------------------!
172 SUBROUTINE data_output_3d( av )
173 
174
175    USE arrays_3d,                                                             &
176        ONLY:  e, nc, nr, p, pt, prr, q, qc, ql, ql_c, ql_v, qr, rho_ocean, s, &
177               sa, tend, u, v, vpt, w
178       
179    USE averaging
180       
181#if defined( __chem )
182    USE chemistry_model_mod,                                                   &
183        ONLY:  chem_data_output_3d
184#endif       
185
186    USE cloud_parameters,                                                      &
187        ONLY:  l_d_cp, pt_d_t
188       
189    USE control_parameters,                                                    &
190        ONLY:  air_chemistry, cloud_physics, do3d, do3d_no, do3d_time_count,   &
191               io_blocks, io_group, land_surface, message_string,              &
192               ntdim_3d, nz_do3d,                                              &
193               psolver, simulated_time, time_since_reference_point,            &
194               urban_surface, varnamelength
195       
196    USE cpulog,                                                                &
197        ONLY:  log_point, cpu_log
198       
199    USE indices,                                                               &
200        ONLY:  nbgp, nx, nxl, nxr, ny, nyn, nys, nzb, wall_flags_0
201       
202    USE kinds
203   
204    USE land_surface_model_mod,                                                &
205        ONLY: lsm_data_output_3d, nzb_soil, nzt_soil
206
207#if defined( __netcdf )
208    USE NETCDF
209#endif
210
211    USE netcdf_interface,                                                      &
212        ONLY:  fill_value, id_set_3d, id_var_do3d, id_var_time_3d, nc_stat,    &
213               netcdf_data_format, netcdf_handle_error
214       
215    USE particle_attributes,                                                   &
216        ONLY:  grid_particles, number_of_particles, particles,                 &
217               particle_advection_start, prt_count
218       
219    USE pegrid
220
221    USE plant_canopy_model_mod,                                                &
222        ONLY:  pcm_data_output_3d, plant_canopy
223       
224    USE radiation_model_mod,                                                   &
225        ONLY:  nzub, nzut, radiation, radiation_data_output_3d
226
227    USE turbulence_closure_mod,                                                &
228        ONLY:  tcm_data_output_3d
229
230    USE urban_surface_mod,                                                     &
231        ONLY:  usm_data_output_3d
232
233
234    IMPLICIT NONE
235
236    INTEGER(iwp) ::  av        !<
237    INTEGER(iwp) ::  flag_nr   !< number of masking flag
238    INTEGER(iwp) ::  i         !<
239    INTEGER(iwp) ::  if        !<
240    INTEGER(iwp) ::  j         !<
241    INTEGER(iwp) ::  k         !<
242    INTEGER(iwp) ::  n         !<
243    INTEGER(iwp) ::  nzb_do    !< vertical lower limit for data output
244    INTEGER(iwp) ::  nzt_do    !< vertical upper limit for data output
245
246    LOGICAL      ::  found     !<
247    LOGICAL      ::  resorted  !<
248
249    REAL(wp)     ::  mean_r    !<
250    REAL(wp)     ::  s_r2      !<
251    REAL(wp)     ::  s_r3      !<
252
253    REAL(sp), DIMENSION(:,:,:), ALLOCATABLE ::  local_pf  !<
254
255    REAL(wp), DIMENSION(:,:,:), POINTER ::  to_be_resorted  !<
256
257    CHARACTER (LEN=varnamelength) ::  trimvar  !< TRIM of output-variable string
258
259!
260!-- Return, if nothing to output
261    IF ( do3d_no(av) == 0 )  RETURN
262
263    CALL cpu_log (log_point(14),'data_output_3d','start')
264
265!
266!-- Open output file.
267!-- For classic or 64bit netCDF output on more than one PE, each PE opens its
268!-- own file and writes the data of its subdomain in binary format. After the
269!-- run, these files are combined to one NetCDF file by combine_plot_fields.
270!-- For netCDF4/HDF5 output, data is written in parallel into one file.
271    IF ( netcdf_data_format < 5 )  THEN
272       CALL check_open( 30 )
273       IF ( myid == 0 )  CALL check_open( 106+av*10 )
274    ELSE
275       CALL check_open( 106+av*10 )
276    ENDIF
277
278!
279!-- For parallel netcdf output the time axis must be limited. Return, if this
280!-- limit is exceeded. This could be the case, if the simulated time exceeds
281!-- the given end time by the length of the given output interval.
282    IF ( netcdf_data_format > 4 )  THEN
283       IF ( do3d_time_count(av) + 1 > ntdim_3d(av) )  THEN
284          WRITE ( message_string, * ) 'Output of 3d data is not given at t=',  &
285                                      simulated_time, '&because the maximum ', & 
286                                      'number of output time levels is ',      &
287                                      'exceeded.'
288          CALL message( 'data_output_3d', 'PA0387', 0, 1, 0, 6, 0 )
289          CALL cpu_log( log_point(14), 'data_output_3d', 'stop' )
290          RETURN
291       ENDIF
292    ENDIF
293
294!
295!-- Update the netCDF time axis
296!-- In case of parallel output, this is only done by PE0 to increase the
297!-- performance.
298#if defined( __netcdf )
299    do3d_time_count(av) = do3d_time_count(av) + 1
300    IF ( myid == 0 )  THEN
301       nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_time_3d(av),           &
302                               (/ time_since_reference_point /),            &
303                               start = (/ do3d_time_count(av) /),           &
304                               count = (/ 1 /) )
305       CALL netcdf_handle_error( 'data_output_3d', 376 )
306    ENDIF
307#endif
308
309!
310!-- Loop over all variables to be written.
311    if = 1
312
313    DO  WHILE ( do3d(av,if)(1:1) /= ' ' )
314
315!
316!--    Temporary solution to account for data output within the new urban
317!--    surface model (urban_surface_mod.f90), see also SELECT CASE ( trimvar ).
318!--    Store the array chosen on the temporary array.
319       trimvar = TRIM( do3d(av,if) )
320       IF ( urban_surface  .AND.  trimvar(1:4) == 'usm_' )  THEN
321          trimvar = 'usm_output'
322          resorted = .TRUE.
323          nzb_do   = nzub
324          nzt_do   = nzut
325       ELSE
326          resorted = .FALSE.
327          nzb_do   = nzb
328          nzt_do   = nz_do3d
329       ENDIF
330!
331!--    Set flag to steer output of radiation, land-surface, or user-defined
332!--    quantities
333       found = .FALSE.
334!
335!--    Allocate a temporary array with the desired output dimensions.
336       ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )
337!
338!--    Before each output, set array local_pf to fill value
339       local_pf = fill_value
340!
341!--    Set masking flag for topography for not resorted arrays
342       flag_nr = 0
343
344       SELECT CASE ( trimvar )
345
346          CASE ( 'e' )
347             IF ( av == 0 )  THEN
348                to_be_resorted => e
349             ELSE
350                to_be_resorted => e_av
351             ENDIF
352
353          CASE ( 'lpt' )
354             IF ( av == 0 )  THEN
355                to_be_resorted => pt
356             ELSE
357                to_be_resorted => lpt_av
358             ENDIF
359
360          CASE ( 'nc' )
361             IF ( av == 0 )  THEN
362                to_be_resorted => nc
363             ELSE
364                to_be_resorted => nc_av
365             ENDIF
366
367          CASE ( 'nr' )
368             IF ( av == 0 )  THEN
369                to_be_resorted => nr
370             ELSE
371                to_be_resorted => nr_av
372             ENDIF
373
374          CASE ( 'p' )
375             IF ( av == 0 )  THEN
376                IF ( psolver /= 'sor' )  CALL exchange_horiz( p, nbgp )
377                to_be_resorted => p
378             ELSE
379                IF ( psolver /= 'sor' )  CALL exchange_horiz( p_av, nbgp )
380                to_be_resorted => p_av
381             ENDIF
382
383          CASE ( 'pc' )  ! particle concentration (requires ghostpoint exchange)
384             IF ( av == 0 )  THEN
385                IF ( simulated_time >= particle_advection_start )  THEN
386                   tend = prt_count
387                ELSE
388                   tend = 0.0_wp
389                ENDIF
390                DO  i = nxl, nxr
391                   DO  j = nys, nyn
392                      DO  k = nzb_do, nzt_do
393                         local_pf(i,j,k) = tend(k,j,i)
394                      ENDDO
395                   ENDDO
396                ENDDO
397                resorted = .TRUE.
398             ELSE
399                to_be_resorted => pc_av
400             ENDIF
401
402          CASE ( 'pr' )  ! mean particle radius (effective radius)
403             IF ( av == 0 )  THEN
404                IF ( simulated_time >= particle_advection_start )  THEN
405                   DO  i = nxl, nxr
406                      DO  j = nys, nyn
407                         DO  k = nzb_do, nzt_do
408                            number_of_particles = prt_count(k,j,i)
409                            IF (number_of_particles <= 0)  CYCLE
410                            particles => grid_particles(k,j,i)%particles(1:number_of_particles)
411                            s_r2 = 0.0_wp
412                            s_r3 = 0.0_wp
413                            DO  n = 1, number_of_particles
414                               IF ( particles(n)%particle_mask )  THEN
415                                  s_r2 = s_r2 + particles(n)%radius**2 * &
416                                         particles(n)%weight_factor
417                                  s_r3 = s_r3 + particles(n)%radius**3 * &
418                                         particles(n)%weight_factor
419                               ENDIF
420                            ENDDO
421                            IF ( s_r2 > 0.0_wp )  THEN
422                               mean_r = s_r3 / s_r2
423                            ELSE
424                               mean_r = 0.0_wp
425                            ENDIF
426                            tend(k,j,i) = mean_r
427                         ENDDO
428                      ENDDO
429                   ENDDO
430                ELSE
431                   tend = 0.0_wp
432                ENDIF
433                DO  i = nxl, nxr
434                   DO  j = nys, nyn
435                      DO  k = nzb_do, nzt_do
436                         local_pf(i,j,k) = tend(k,j,i)
437                      ENDDO
438                   ENDDO
439                ENDDO
440                resorted = .TRUE.
441             ELSE
442                to_be_resorted => pr_av
443             ENDIF
444
445          CASE ( 'prr' )
446             IF ( av == 0 )  THEN
447                DO  i = nxl, nxr
448                   DO  j = nys, nyn
449                      DO  k = nzb_do, nzt_do
450                         local_pf(i,j,k) = prr(k,j,i)
451                      ENDDO
452                   ENDDO
453                ENDDO
454             ELSE
455                DO  i = nxl, nxr
456                   DO  j = nys, nyn
457                      DO  k = nzb_do, nzt_do
458                         local_pf(i,j,k) = prr_av(k,j,i)
459                      ENDDO
460                   ENDDO
461                ENDDO
462             ENDIF
463             resorted = .TRUE.
464
465          CASE ( 'pt' )
466             IF ( av == 0 )  THEN
467                IF ( .NOT. cloud_physics ) THEN
468                   to_be_resorted => pt
469                ELSE
470                   DO  i = nxl, nxr
471                      DO  j = nys, nyn
472                         DO  k = nzb_do, nzt_do
473                            local_pf(i,j,k) = pt(k,j,i) + l_d_cp *             &
474                                                          pt_d_t(k) *          &
475                                                          ql(k,j,i)
476                         ENDDO
477                      ENDDO
478                   ENDDO
479                   resorted = .TRUE.
480                ENDIF
481             ELSE
482                to_be_resorted => pt_av
483             ENDIF
484
485          CASE ( 'q' )
486             IF ( av == 0 )  THEN
487                to_be_resorted => q
488             ELSE
489                to_be_resorted => q_av
490             ENDIF
491
492          CASE ( 'qc' )
493             IF ( av == 0 )  THEN
494                to_be_resorted => qc
495             ELSE
496                to_be_resorted => qc_av
497             ENDIF
498
499          CASE ( 'ql' )
500             IF ( av == 0 )  THEN
501                to_be_resorted => ql
502             ELSE
503                to_be_resorted => ql_av
504             ENDIF
505
506          CASE ( 'ql_c' )
507             IF ( av == 0 )  THEN
508                to_be_resorted => ql_c
509             ELSE
510                to_be_resorted => ql_c_av
511             ENDIF
512
513          CASE ( 'ql_v' )
514             IF ( av == 0 )  THEN
515                to_be_resorted => ql_v
516             ELSE
517                to_be_resorted => ql_v_av
518             ENDIF
519
520          CASE ( 'ql_vp' )
521             IF ( av == 0 )  THEN
522                IF ( simulated_time >= particle_advection_start )  THEN
523                   DO  i = nxl, nxr
524                      DO  j = nys, nyn
525                         DO  k = nzb_do, nzt_do
526                            number_of_particles = prt_count(k,j,i)
527                            IF (number_of_particles <= 0)  CYCLE
528                            particles => grid_particles(k,j,i)%particles(1:number_of_particles)
529                            DO  n = 1, number_of_particles
530                               IF ( particles(n)%particle_mask )  THEN
531                                  tend(k,j,i) =  tend(k,j,i) +                 &
532                                                 particles(n)%weight_factor /  &
533                                                 prt_count(k,j,i)
534                               ENDIF
535                            ENDDO
536                         ENDDO
537                      ENDDO
538                   ENDDO
539                ELSE
540                   tend = 0.0_wp
541                ENDIF
542                DO  i = nxl, nxr
543                   DO  j = nys, nyn
544                      DO  k = nzb_do, nzt_do
545                         local_pf(i,j,k) = tend(k,j,i)
546                      ENDDO
547                   ENDDO
548                ENDDO
549                resorted = .TRUE.
550             ELSE
551                to_be_resorted => ql_vp_av
552             ENDIF
553
554          CASE ( 'qr' )
555             IF ( av == 0 )  THEN
556                to_be_resorted => qr
557             ELSE
558                to_be_resorted => qr_av
559             ENDIF
560
561          CASE ( 'qv' )
562             IF ( av == 0 )  THEN
563                DO  i = nxl, nxr
564                   DO  j = nys, nyn
565                      DO  k = nzb_do, nzt_do
566                         local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
567                      ENDDO
568                   ENDDO
569                ENDDO
570                resorted = .TRUE.
571             ELSE
572                to_be_resorted => qv_av
573             ENDIF
574
575          CASE ( 'rho_ocean' )
576             IF ( av == 0 )  THEN
577                to_be_resorted => rho_ocean
578             ELSE
579                to_be_resorted => rho_ocean_av
580             ENDIF
581
582          CASE ( 's' )
583             IF ( av == 0 )  THEN
584                to_be_resorted => s
585             ELSE
586                to_be_resorted => s_av
587             ENDIF
588
589          CASE ( 'sa' )
590             IF ( av == 0 )  THEN
591                to_be_resorted => sa
592             ELSE
593                to_be_resorted => sa_av
594             ENDIF
595
596          CASE ( 'u' )
597             flag_nr = 1
598             IF ( av == 0 )  THEN
599                to_be_resorted => u
600             ELSE
601                to_be_resorted => u_av
602             ENDIF
603
604          CASE ( 'v' )
605             flag_nr = 2
606             IF ( av == 0 )  THEN
607                to_be_resorted => v
608             ELSE
609                to_be_resorted => v_av
610             ENDIF
611
612          CASE ( 'vpt' )
613             IF ( av == 0 )  THEN
614                to_be_resorted => vpt
615             ELSE
616                to_be_resorted => vpt_av
617             ENDIF
618
619          CASE ( 'w' )
620             flag_nr = 3
621             IF ( av == 0 )  THEN
622                to_be_resorted => w
623             ELSE
624                to_be_resorted => w_av
625             ENDIF
626!             
627!--       Block of urban surface model outputs   
628          CASE ( 'usm_output' )
629             CALL usm_data_output_3d( av, do3d(av,if), found, local_pf,     &
630                                         nzb_do, nzt_do )
631
632          CASE DEFAULT
633
634!
635!--          Land surface quantity
636             IF ( land_surface )  THEN
637!
638!--             For soil model quantities, it is required to re-allocate local_pf
639                nzb_do = nzb_soil
640                nzt_do = nzt_soil
641
642                DEALLOCATE ( local_pf )
643                ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )
644                local_pf = fill_value
645
646                CALL lsm_data_output_3d( av, do3d(av,if), found, local_pf )
647                resorted = .TRUE.
648
649!
650!--             If no soil model variable was found, re-allocate local_pf
651                IF ( .NOT. found )  THEN
652                   nzb_do = nzb
653                   nzt_do = nz_do3d
654
655                   DEALLOCATE ( local_pf )
656                   ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )                 
657                ENDIF
658
659             ENDIF
660
661             IF ( .NOT. found )  THEN
662                CALL tcm_data_output_3d( av, do3d(av,if), found, local_pf )
663                resorted = .TRUE.
664             ENDIF
665
666!
667!--          Radiation quantity
668             IF ( .NOT. found  .AND.  radiation )  THEN
669                CALL radiation_data_output_3d( av, do3d(av,if), found,         &
670                                               local_pf )
671                resorted = .TRUE.
672             ENDIF
673
674!
675!--          Chemistry model output
676#if defined( __chem )
677             IF ( .NOT. found  .AND.  air_chemistry )  THEN
678                CALL chem_data_output_3d( av, do3d(av,if), found,              &
679                                          local_pf )
680                resorted = .TRUE.
681             ENDIF
682#endif
683
684!
685!--          Plant canopy model output
686             IF ( .NOT. found  .AND.  plant_canopy )  THEN
687                CALL pcm_data_output_3d( av, do3d(av,if), found, local_pf,     &
688                                         fill_value )
689                resorted = .TRUE.
690             ENDIF
691
692!
693!--          User defined quantity
694             IF ( .NOT. found )  THEN
695                CALL user_data_output_3d( av, do3d(av,if), found, local_pf,    &
696                                          nzb_do, nzt_do )
697                resorted = .TRUE.
698             ENDIF
699
700             IF ( .NOT. found )  THEN
701                message_string =  'no output available for: ' //               &
702                                  TRIM( do3d(av,if) )
703                CALL message( 'data_output_3d', 'PA0182', 0, 0, 0, 6, 0 )
704             ENDIF
705
706       END SELECT
707
708!
709!--    Resort the array to be output, if not done above
710       IF ( .NOT. resorted )  THEN
711          DO  i = nxl, nxr
712             DO  j = nys, nyn
713                DO  k = nzb_do, nzt_do
714                   local_pf(i,j,k) = MERGE(                                    &
715                                      to_be_resorted(k,j,i),                   &
716                                      REAL( fill_value, KIND = wp ),           &
717                                      BTEST( wall_flags_0(k,j,i), flag_nr ) )
718                ENDDO
719             ENDDO
720          ENDDO
721       ENDIF
722
723!
724!--    Output of the 3D-array
725#if defined( __parallel )
726       IF ( netcdf_data_format < 5 )  THEN
727!
728!--       Non-parallel netCDF output. Data is output in parallel in
729!--       FORTRAN binary format here, and later collected into one file by
730!--       combine_plot_fields
731          IF ( myid == 0 )  THEN
732             WRITE ( 30 )  time_since_reference_point,                   &
733                           do3d_time_count(av), av
734          ENDIF
735          DO  i = 0, io_blocks-1
736             IF ( i == io_group )  THEN
737                WRITE ( 30 )  nxl, nxr, nys, nyn, nzb_do, nzt_do
738                WRITE ( 30 )  local_pf(:,:,nzb_do:nzt_do)
739             ENDIF
740
741             CALL MPI_BARRIER( comm2d, ierr )
742
743          ENDDO
744
745       ELSE
746#if defined( __netcdf )
747!
748!--       Parallel output in netCDF4/HDF5 format.
749!          IF ( nxr == nx  .AND.  nyn /= ny )  THEN
750!             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
751!                               local_pf(nxl:nxr+1,nys:nyn,nzb_do:nzt_do),    &
752!                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
753!                count = (/ nxr-nxl+2, nyn-nys+1, nzt_do-nzb_do+1, 1 /) )
754!          ELSEIF ( nxr /= nx  .AND.  nyn == ny )  THEN
755!             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
756!                               local_pf(nxl:nxr,nys:nyn+1,nzb_do:nzt_do),    &
757!                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
758!                count = (/ nxr-nxl+1, nyn-nys+2, nzt_do-nzb_do+1, 1 /) )
759!          ELSEIF ( nxr == nx  .AND.  nyn == ny )  THEN
760!             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
761!                             local_pf(nxl:nxr+1,nys:nyn+1,nzb_do:nzt_do  ),  &
762!                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
763!                count = (/ nxr-nxl+2, nyn-nys+2, nzt_do-nzb_do+1, 1 /) )
764!          ELSE
765             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
766                                 local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do),    &
767                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
768                count = (/ nxr-nxl+1, nyn-nys+1, nzt_do-nzb_do+1, 1 /) )
769!          ENDIF
770          CALL netcdf_handle_error( 'data_output_3d', 386 )
771#endif
772       ENDIF
773#else
774#if defined( __netcdf )
775       nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),        &
776                         local_pf(nxl:nxr+1,nys:nyn+1,nzb_do:nzt_do),        &
777                         start = (/ 1, 1, 1, do3d_time_count(av) /),     &
778                         count = (/ nx+1, ny+1, nzt_do-nzb_do+1, 1 /) )
779       CALL netcdf_handle_error( 'data_output_3d', 446 )
780#endif
781#endif
782
783       if = if + 1
784
785!
786!--    Deallocate temporary array
787       DEALLOCATE ( local_pf )
788
789    ENDDO
790
791    CALL cpu_log( log_point(14), 'data_output_3d', 'stop' )
792
793!
794!-- Formats.
7953300 FORMAT ('variable ',I4,'  file=',A,'  filetype=unformatted  skip=',I12/   &
796             'label = ',A,A)
797
798 END SUBROUTINE data_output_3d
Note: See TracBrowser for help on using the repository browser.