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

Last change on this file since 2903 was 2817, checked in by knoop, 7 years ago

Preliminary gust module interface implemented

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