source: palm/trunk/SOURCE/data_output_2d.f90 @ 760

Last change on this file since 760 was 760, checked in by raasch, 11 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 65.1 KB
Line 
1 SUBROUTINE data_output_2d( mode, av )
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: data_output_2d.f90 760 2011-09-15 14:37:54Z raasch $
11!
12! 759 2011-09-15 13:58:31Z raasch
13! Splitting of parallel I/O
14!
15! 729 2011-05-26 10:33:34Z heinze
16! Exchange ghost layers for p regardless of used pressure solver (except SOR).
17!
18! 691 2011-03-04 08:45:30Z maronga
19! Replaced simulated_time by time_since_reference_point
20!
21! 673 2011-01-18 16:19:48Z suehring
22! When using Multigrid or SOR solver an additional CALL exchange_horiz is
23! is needed for pressure output.
24!
25! 667 2010-12-23 12:06:00Z suehring/gryschka
26! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng in loops and
27! allocation of arrays local_2d and total_2d.
28! Calls of exchange_horiz are modiefied.
29!
30! 622 2010-12-10 08:08:13Z raasch
31! optional barriers included in order to speed up collective operations
32!
33! 493 2010-03-01 08:30:24Z raasch
34! NetCDF4 support (parallel output)
35!
36! 367 2009-08-25 08:35:52Z maronga
37! simulated_time in NetCDF output replaced by time_since_reference_point.
38! Output of NetCDF messages with aid of message handling routine.
39! Bugfix: averaging along z is not allowed for 2d quantities (e.g. u* and z0)
40! Output of messages replaced by message handling routine.
41! Output of user defined 2D (XY) arrays at z=nzb+1 is now possible
42! Bugfix: to_be_resorted => s_av for time-averaged scalars
43! Calculation of shf* and qsws* added.
44!
45! 215 2008-11-18 09:54:31Z raasch
46! Bugfix: no output of particle concentration and radius unless particles
47! have been started
48!
49! 96 2007-06-04 08:07:41Z raasch
50! Output of density and salinity
51!
52! 75 2007-03-22 09:54:05Z raasch
53! Output of precipitation amount/rate and roughness length,
54! 2nd+3rd argument removed from exchange horiz
55!
56! RCS Log replace by Id keyword, revision history cleaned up
57!
58! Revision 1.5  2006/08/22 13:50:29  raasch
59! xz and yz cross sections now up to nzt+1
60!
61! Revision 1.2  2006/02/23 10:19:22  raasch
62! Output of time-averaged data, output of averages along x, y, or z,
63! output of user-defined quantities,
64! section data are copied from local_pf to local_2d before they are output,
65! output of particle concentration and mean radius,
66! Former subroutine plot_2d renamed data_output_2d, pl2d.. renamed do2d..,
67! anz renamed ngp, ebene renamed section, pl2d_.._anz renamed do2d_.._n
68!
69! Revision 1.1  1997/08/11 06:24:09  raasch
70! Initial revision
71!
72!
73! Description:
74! ------------
75! Data output of horizontal cross-sections in NetCDF format or binary format
76! compatible to old graphic software iso2d.
77! Attention: The position of the sectional planes is still not always computed
78! ---------  correctly. (zu is used always)!
79!------------------------------------------------------------------------------!
80
81    USE arrays_3d
82    USE averaging
83    USE cloud_parameters
84    USE control_parameters
85    USE cpulog
86    USE grid_variables
87    USE indices
88    USE interfaces
89    USE netcdf_control
90    USE particle_attributes
91    USE pegrid
92
93    IMPLICIT NONE
94
95    CHARACTER (LEN=2)  ::  do2d_mode, mode
96    CHARACTER (LEN=4)  ::  grid
97    CHARACTER (LEN=25) ::  section_chr
98    CHARACTER (LEN=50) ::  rtext
99    INTEGER ::  av, ngp, file_id, i, if, is, iis, j, k, l, layer_xy, n, psi, &
100                s, sender, &
101                ind(4)
102    LOGICAL ::  found, resorted, two_d
103    REAL    ::  mean_r, s_r3, s_r4
104    REAL, DIMENSION(:), ALLOCATABLE ::      level_z
105    REAL, DIMENSION(:,:), ALLOCATABLE ::    local_2d, local_2d_l
106    REAL, DIMENSION(:,:,:), ALLOCATABLE ::  local_pf
107#if defined( __parallel )
108    REAL, DIMENSION(:,:),   ALLOCATABLE ::  total_2d
109#endif
110    REAL, DIMENSION(:,:,:), POINTER ::  to_be_resorted
111
112    NAMELIST /LOCAL/  rtext
113
114    CALL cpu_log (log_point(3),'data_output_2d','start')
115
116!
117!-- Immediate return, if no output is requested (no respective sections
118!-- found in parameter data_output)
119    IF ( mode == 'xy'  .AND.  .NOT. data_output_xy(av) )  RETURN
120    IF ( mode == 'xz'  .AND.  .NOT. data_output_xz(av) )  RETURN
121    IF ( mode == 'yz'  .AND.  .NOT. data_output_yz(av) )  RETURN
122
123    two_d = .FALSE.    ! local variable to distinguish between output of pure 2D
124                       ! arrays and cross-sections of 3D arrays.
125
126!
127!-- Depending on the orientation of the cross-section, the respective output
128!-- files have to be opened.
129    SELECT CASE ( mode )
130
131       CASE ( 'xy' )
132          s = 1
133          ALLOCATE( level_z(nzb:nzt+1), local_2d(nxlg:nxrg,nysg:nyng) )
134
135!
136!--       Classic and 64bit offset NetCDF output is done only on PE0.
137!--       netCDF4/HDF5 output is done in parallel on all PEs.
138          IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
139          THEN
140             CALL check_open( 101+av*10 )
141          ENDIF
142
143          IF ( data_output_2d_on_each_pe )  THEN
144             CALL check_open( 21 )
145          ELSE
146             IF ( myid == 0 )  THEN
147                IF ( iso2d_output )  CALL check_open( 21 )
148#if defined( __parallel )
149                ALLOCATE( total_2d(-nbgp:nx+nbgp,-nbgp:ny+nbgp) )
150#endif
151             ENDIF
152          ENDIF
153
154       CASE ( 'xz' )
155          s = 2
156          ALLOCATE( local_2d(nxlg:nxrg,nzb:nzt+1) )
157
158!
159!--       Classic and 64bit offset NetCDF output is done only on PE0.
160!--       netCDF4/HDF5 output may be done in parallel on all PEs.
161          IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
162          THEN
163             CALL check_open( 102+av*10 )
164          ENDIF
165
166          IF ( data_output_2d_on_each_pe )  THEN
167             CALL check_open( 22 )
168          ELSE
169             IF ( myid == 0 )  THEN
170                IF ( iso2d_output )  CALL check_open( 22 )
171#if defined( __parallel )
172                ALLOCATE( total_2d(-nbgp:nx+nbgp,nzb:nzt+1) )
173#endif
174             ENDIF
175          ENDIF
176
177       CASE ( 'yz' )
178
179          s = 3
180          ALLOCATE( local_2d(nysg:nyng,nzb:nzt+1) )
181
182!
183!--       Classic and 64bit offset NetCDF output is done only on PE0.
184!--       netCDF4/HDF5 output may be done in parallel on all PEs.
185          IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
186          THEN
187             CALL check_open( 103+av*10 )
188          ENDIF
189
190          IF ( data_output_2d_on_each_pe )  THEN
191             CALL check_open( 23 )
192          ELSE
193             IF ( myid == 0 )  THEN
194                IF ( iso2d_output )  CALL check_open( 23 )
195#if defined( __parallel )
196                ALLOCATE( total_2d(-nbgp:ny+nbgp,nzb:nzt+1) )
197#endif
198             ENDIF
199          ENDIF
200
201       CASE DEFAULT
202
203          message_string = 'unknown cross-section: ' // TRIM( mode )
204          CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
205
206    END SELECT
207
208!
209!-- Allocate a temporary array for resorting (kji -> ijk).
210    ALLOCATE( local_pf(nxlg:nxrg,nysg:nyng,nzb:nzt+1) )
211
212!
213!-- Loop of all variables to be written.
214!-- Output dimensions chosen
215    if = 1
216    l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
217    do2d_mode = do2d(av,if)(l-1:l)
218
219    DO  WHILE ( do2d(av,if)(1:1) /= ' ' )
220
221       IF ( do2d_mode == mode )  THEN
222!
223!--       Store the array chosen on the temporary array.
224          resorted = .FALSE.
225          SELECT CASE ( TRIM( do2d(av,if) ) )
226
227             CASE ( 'e_xy', 'e_xz', 'e_yz' )
228                IF ( av == 0 )  THEN
229                   to_be_resorted => e
230                ELSE
231                   to_be_resorted => e_av
232                ENDIF
233                IF ( mode == 'xy' )  level_z = zu
234
235             CASE ( 'lwp*_xy' )        ! 2d-array
236                IF ( av == 0 )  THEN
237                   DO  i = nxlg, nxrg
238                      DO  j = nysg, nyng
239                         local_pf(i,j,nzb+1) = SUM( ql(nzb:nzt,j,i) * &
240                                                    dzw(1:nzt+1) )
241                      ENDDO
242                   ENDDO
243                ELSE
244                   DO  i = nxlg, nxrg
245                      DO  j = nysg, nyng
246                         local_pf(i,j,nzb+1) = lwp_av(j,i)
247                      ENDDO
248                   ENDDO
249                ENDIF
250                resorted = .TRUE.
251                two_d = .TRUE.
252                level_z(nzb+1) = zu(nzb+1)
253
254             CASE ( 'p_xy', 'p_xz', 'p_yz' )
255                IF ( av == 0 )  THEN
256                   IF ( psolver /= 'sor' )  CALL exchange_horiz( p, nbgp )
257                   to_be_resorted => p
258                ELSE
259                   IF ( psolver /= 'sor' )  CALL exchange_horiz( p_av, nbgp )
260                   to_be_resorted => p_av
261                ENDIF
262                IF ( mode == 'xy' )  level_z = zu
263
264             CASE ( 'pc_xy', 'pc_xz', 'pc_yz' )  ! particle concentration
265                IF ( av == 0 )  THEN
266                   IF ( simulated_time >= particle_advection_start )  THEN
267                      tend = prt_count
268                      CALL exchange_horiz( tend, nbgp )
269                   ELSE
270                      tend = 0.0
271                   ENDIF
272                   DO  i = nxlg, nxrg
273                      DO  j = nysg, nyng
274                         DO  k = nzb, nzt+1
275                            local_pf(i,j,k) = tend(k,j,i)
276                         ENDDO
277                      ENDDO
278                   ENDDO
279                   resorted = .TRUE.
280                ELSE
281                   CALL exchange_horiz( pc_av, nbgp )
282                   to_be_resorted => pc_av
283                ENDIF
284
285             CASE ( 'pr_xy', 'pr_xz', 'pr_yz' )  ! mean particle radius
286                IF ( av == 0 )  THEN
287                   IF ( simulated_time >= particle_advection_start )  THEN
288                      DO  i = nxl, nxr
289                         DO  j = nys, nyn
290                            DO  k = nzb, nzt+1
291                               psi = prt_start_index(k,j,i)
292                               s_r3 = 0.0
293                               s_r4 = 0.0
294                               DO  n = psi, psi+prt_count(k,j,i)-1
295                                  s_r3 = s_r3 + particles(n)%radius**3
296                                  s_r4 = s_r4 + particles(n)%radius**4
297                               ENDDO
298                               IF ( s_r3 /= 0.0 )  THEN
299                                  mean_r = s_r4 / s_r3
300                               ELSE
301                                  mean_r = 0.0
302                               ENDIF
303                               tend(k,j,i) = mean_r
304                            ENDDO
305                         ENDDO
306                      ENDDO
307                      CALL exchange_horiz( tend, nbgp )
308                   ELSE
309                      tend = 0.0
310                   END IF
311                   DO  i = nxlg, nxrg
312                      DO  j = nysg, nyng
313                         DO  k = nzb, nzt+1
314                            local_pf(i,j,k) = tend(k,j,i)
315                         ENDDO
316                      ENDDO
317                   ENDDO
318                   resorted = .TRUE.
319                ELSE
320                   CALL exchange_horiz( pr_av, nbgp )
321                   to_be_resorted => pr_av
322                ENDIF
323
324             CASE ( 'pra*_xy' )        ! 2d-array / integral quantity => no av
325                CALL exchange_horiz_2d( precipitation_amount )
326                   DO  i = nxlg, nxrg
327                      DO  j = nysg, nyng
328                      local_pf(i,j,nzb+1) =  precipitation_amount(j,i)
329                   ENDDO
330                ENDDO
331                precipitation_amount = 0.0   ! reset for next integ. interval
332                resorted = .TRUE.
333                two_d = .TRUE.
334                level_z(nzb+1) = zu(nzb+1)
335
336             CASE ( 'prr*_xy' )        ! 2d-array
337                IF ( av == 0 )  THEN
338                   CALL exchange_horiz_2d( precipitation_rate )
339                   DO  i = nxlg, nxrg
340                      DO  j = nysg, nyng
341                         local_pf(i,j,nzb+1) =  precipitation_rate(j,i)
342                      ENDDO
343                   ENDDO
344                ELSE
345                   CALL exchange_horiz_2d( precipitation_rate_av )
346                   DO  i = nxlg, nxrg
347                      DO  j = nysg, nyng
348                         local_pf(i,j,nzb+1) =  precipitation_rate_av(j,i)
349                      ENDDO
350                   ENDDO
351                ENDIF
352                resorted = .TRUE.
353                two_d = .TRUE.
354                level_z(nzb+1) = zu(nzb+1)
355
356             CASE ( 'pt_xy', 'pt_xz', 'pt_yz' )
357                IF ( av == 0 )  THEN
358                   IF ( .NOT. cloud_physics ) THEN
359                      to_be_resorted => pt
360                   ELSE
361                   DO  i = nxlg, nxrg
362                      DO  j = nysg, nyng
363                            DO  k = nzb, nzt+1
364                               local_pf(i,j,k) = pt(k,j,i) + l_d_cp *    &
365                                                             pt_d_t(k) * &
366                                                             ql(k,j,i)
367                            ENDDO
368                         ENDDO
369                      ENDDO
370                      resorted = .TRUE.
371                   ENDIF
372                ELSE
373                   to_be_resorted => pt_av
374                ENDIF
375                IF ( mode == 'xy' )  level_z = zu
376
377             CASE ( 'q_xy', 'q_xz', 'q_yz' )
378                IF ( av == 0 )  THEN
379                   to_be_resorted => q
380                ELSE
381                   to_be_resorted => q_av
382                ENDIF
383                IF ( mode == 'xy' )  level_z = zu
384
385             CASE ( 'ql_xy', 'ql_xz', 'ql_yz' )
386                IF ( av == 0 )  THEN
387                   to_be_resorted => ql
388                ELSE
389                   to_be_resorted => ql_av
390                ENDIF
391                IF ( mode == 'xy' )  level_z = zu
392
393             CASE ( 'ql_c_xy', 'ql_c_xz', 'ql_c_yz' )
394                IF ( av == 0 )  THEN
395                   to_be_resorted => ql_c
396                ELSE
397                   to_be_resorted => ql_c_av
398                ENDIF
399                IF ( mode == 'xy' )  level_z = zu
400
401             CASE ( 'ql_v_xy', 'ql_v_xz', 'ql_v_yz' )
402                IF ( av == 0 )  THEN
403                   to_be_resorted => ql_v
404                ELSE
405                   to_be_resorted => ql_v_av
406                ENDIF
407                IF ( mode == 'xy' )  level_z = zu
408
409             CASE ( 'ql_vp_xy', 'ql_vp_xz', 'ql_vp_yz' )
410                IF ( av == 0 )  THEN
411                   to_be_resorted => ql_vp
412                ELSE
413                   to_be_resorted => ql_vp_av
414                ENDIF
415                IF ( mode == 'xy' )  level_z = zu
416
417             CASE ( 'qsws*_xy' )        ! 2d-array
418                IF ( av == 0 ) THEN
419                   DO  i = nxlg, nxrg
420                      DO  j = nysg, nyng
421                         local_pf(i,j,nzb+1) =  qsws(j,i)
422                      ENDDO
423                   ENDDO
424                ELSE
425                   DO  i = nxlg, nxrg
426                      DO  j = nysg, nyng 
427                         local_pf(i,j,nzb+1) =  qsws_av(j,i)
428                      ENDDO
429                   ENDDO
430                ENDIF
431                resorted = .TRUE.
432                two_d = .TRUE.
433                level_z(nzb+1) = zu(nzb+1)
434
435             CASE ( 'qv_xy', 'qv_xz', 'qv_yz' )
436                IF ( av == 0 )  THEN
437                   DO  i = nxlg, nxrg
438                      DO  j = nysg, nyng
439                         DO  k = nzb, nzt+1
440                            local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
441                         ENDDO
442                      ENDDO
443                   ENDDO
444                   resorted = .TRUE.
445                ELSE
446                   to_be_resorted => qv_av
447                ENDIF
448                IF ( mode == 'xy' )  level_z = zu
449
450             CASE ( 'rho_xy', 'rho_xz', 'rho_yz' )
451                IF ( av == 0 )  THEN
452                   to_be_resorted => rho
453                ELSE
454                   to_be_resorted => rho_av
455                ENDIF
456
457             CASE ( 's_xy', 's_xz', 's_yz' )
458                IF ( av == 0 )  THEN
459                   to_be_resorted => q
460                ELSE
461                   to_be_resorted => s_av
462                ENDIF
463
464             CASE ( 'sa_xy', 'sa_xz', 'sa_yz' )
465                IF ( av == 0 )  THEN
466                   to_be_resorted => sa
467                ELSE
468                   to_be_resorted => sa_av
469                ENDIF
470
471             CASE ( 'shf*_xy' )        ! 2d-array
472                IF ( av == 0 ) THEN
473                   DO  i = nxlg, nxrg
474                      DO  j = nysg, nyng
475                         local_pf(i,j,nzb+1) =  shf(j,i)
476                      ENDDO
477                   ENDDO
478                ELSE
479                   DO  i = nxlg, nxrg
480                      DO  j = nysg, nyng
481                         local_pf(i,j,nzb+1) =  shf_av(j,i)
482                      ENDDO
483                   ENDDO
484                ENDIF
485                resorted = .TRUE.
486                two_d = .TRUE.
487                level_z(nzb+1) = zu(nzb+1)
488
489             CASE ( 't*_xy' )        ! 2d-array
490                IF ( av == 0 )  THEN
491                   DO  i = nxlg, nxrg
492                      DO  j = nysg, nyng
493                         local_pf(i,j,nzb+1) = ts(j,i)
494                      ENDDO
495                   ENDDO
496                ELSE
497                   DO  i = nxlg, nxrg
498                      DO  j = nysg, nyng
499                         local_pf(i,j,nzb+1) = ts_av(j,i)
500                      ENDDO
501                   ENDDO
502                ENDIF
503                resorted = .TRUE.
504                two_d = .TRUE.
505                level_z(nzb+1) = zu(nzb+1)
506
507             CASE ( 'u_xy', 'u_xz', 'u_yz' )
508                IF ( av == 0 )  THEN
509                   to_be_resorted => u
510                ELSE
511                   to_be_resorted => u_av
512                ENDIF
513                IF ( mode == 'xy' )  level_z = zu
514!
515!--             Substitute the values generated by "mirror" boundary condition
516!--             at the bottom boundary by the real surface values.
517                IF ( do2d(av,if) == 'u_xz'  .OR.  do2d(av,if) == 'u_yz' )  THEN
518                   IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0
519                ENDIF
520
521             CASE ( 'u*_xy' )        ! 2d-array
522                IF ( av == 0 )  THEN
523                   DO  i = nxlg, nxrg
524                      DO  j = nysg, nyng
525                         local_pf(i,j,nzb+1) = us(j,i)
526                      ENDDO
527                   ENDDO
528                ELSE
529                   DO  i = nxlg, nxrg
530                      DO  j = nysg, nyng
531                         local_pf(i,j,nzb+1) = us_av(j,i)
532                      ENDDO
533                   ENDDO
534                ENDIF
535                resorted = .TRUE.
536                two_d = .TRUE.
537                level_z(nzb+1) = zu(nzb+1)
538
539             CASE ( 'v_xy', 'v_xz', 'v_yz' )
540                IF ( av == 0 )  THEN
541                   to_be_resorted => v
542                ELSE
543                   to_be_resorted => v_av
544                ENDIF
545                IF ( mode == 'xy' )  level_z = zu
546!
547!--             Substitute the values generated by "mirror" boundary condition
548!--             at the bottom boundary by the real surface values.
549                IF ( do2d(av,if) == 'v_xz'  .OR.  do2d(av,if) == 'v_yz' )  THEN
550                   IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0
551                ENDIF
552
553             CASE ( 'vpt_xy', 'vpt_xz', 'vpt_yz' )
554                IF ( av == 0 )  THEN
555                   to_be_resorted => vpt
556                ELSE
557                   to_be_resorted => vpt_av
558                ENDIF
559                IF ( mode == 'xy' )  level_z = zu
560
561             CASE ( 'w_xy', 'w_xz', 'w_yz' )
562                IF ( av == 0 )  THEN
563                   to_be_resorted => w
564                ELSE
565                   to_be_resorted => w_av
566                ENDIF
567                IF ( mode == 'xy' )  level_z = zw
568
569             CASE ( 'z0*_xy' )        ! 2d-array
570                IF ( av == 0 ) THEN
571                   DO  i = nxlg, nxrg
572                      DO  j = nysg, nyng
573                         local_pf(i,j,nzb+1) =  z0(j,i)
574                      ENDDO
575                   ENDDO
576                ELSE
577                   DO  i = nxlg, nxrg
578                      DO  j = nysg, nyng
579                         local_pf(i,j,nzb+1) =  z0_av(j,i)
580                      ENDDO
581                   ENDDO
582                ENDIF
583                resorted = .TRUE.
584                two_d = .TRUE.
585                level_z(nzb+1) = zu(nzb+1)
586
587             CASE DEFAULT
588!
589!--             User defined quantity
590                CALL user_data_output_2d( av, do2d(av,if), found, grid, &
591                                          local_pf, two_d )
592                resorted = .TRUE.
593
594                IF ( grid == 'zu' )  THEN
595                   IF ( mode == 'xy' )  level_z = zu
596                ELSEIF ( grid == 'zw' )  THEN
597                   IF ( mode == 'xy' )  level_z = zw
598                ELSEIF ( grid == 'zu1' ) THEN
599                   IF ( mode == 'xy' )  level_z(nzb+1) = zu(nzb+1)
600                ENDIF
601
602                IF ( .NOT. found )  THEN
603                   message_string = 'no output provided for: ' //    &
604                                    TRIM( do2d(av,if) )
605                   CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
606                ENDIF
607
608          END SELECT
609
610!
611!--       Resort the array to be output, if not done above
612          IF ( .NOT. resorted )  THEN
613             DO  i = nxlg, nxrg
614                DO  j = nysg, nyng
615                   DO  k = nzb, nzt+1
616                      local_pf(i,j,k) = to_be_resorted(k,j,i)
617                   ENDDO
618                ENDDO
619             ENDDO
620          ENDIF
621
622!
623!--       Output of the individual cross-sections, depending on the cross-
624!--       section mode chosen.
625          is = 1
626   loop1: DO  WHILE ( section(is,s) /= -9999  .OR.  two_d )
627
628             SELECT CASE ( mode )
629
630                CASE ( 'xy' )
631!
632!--                Determine the cross section index
633                   IF ( two_d )  THEN
634                      layer_xy = nzb+1
635                   ELSE
636                      layer_xy = section(is,s)
637                   ENDIF
638
639!
640!--                Update the NetCDF xy cross section time axis
641                   IF ( myid == 0  .OR.  netcdf_data_format > 2 )  THEN
642                      IF ( simulated_time /= do2d_xy_last_time(av) )  THEN
643                         do2d_xy_time_count(av) = do2d_xy_time_count(av) + 1
644                         do2d_xy_last_time(av)  = simulated_time
645                         IF ( ( .NOT. data_output_2d_on_each_pe  .AND. &
646                              netcdf_output )  .OR.  netcdf_data_format > 2 ) &
647                         THEN
648#if defined( __netcdf )
649                            nc_stat = NF90_PUT_VAR( id_set_xy(av),             &
650                                                    id_var_time_xy(av),        &
651                                             (/ time_since_reference_point /), &
652                                         start = (/ do2d_xy_time_count(av) /), &
653                                                    count = (/ 1 /) )
654                            CALL handle_netcdf_error( 'data_output_2d', 53 )
655#endif
656                         ENDIF
657                      ENDIF
658                   ENDIF
659!
660!--                If required, carry out averaging along z
661                   IF ( section(is,s) == -1  .AND.  .NOT. two_d )  THEN
662
663                      local_2d = 0.0
664!
665!--                   Carry out the averaging (all data are on the PE)
666                      DO  k = nzb, nzt+1
667                         DO  j = nysg, nyng
668                            DO  i = nxlg, nxrg
669                               local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
670                            ENDDO
671                         ENDDO
672                      ENDDO
673
674                      local_2d = local_2d / ( nzt -nzb + 2.0)
675
676                   ELSE
677!
678!--                   Just store the respective section on the local array
679                      local_2d = local_pf(:,:,layer_xy)
680
681                   ENDIF
682
683#if defined( __parallel )
684                   IF ( netcdf_output  .AND.  netcdf_data_format > 2 )  THEN
685!
686!--                   Output in NetCDF4/HDF5 format.
687!--                   Do not output redundant ghost point data except for the
688!--                   boundaries of the total domain.
689                      IF ( two_d ) THEN
690                         iis = 1
691                      ELSE
692                         iis = is
693                      ENDIF
694
695#if defined( __netcdf )
696                      IF ( nxr == nx  .AND.  nyn /= ny )  THEN
697                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
698                                                 id_var_do2d(av,if),           &
699                                                 local_2d(nxl:nxr+1,nys:nyn),  &
700                                                 start = (/ nxl+1, nys+1, iis, &
701                                                    do2d_xy_time_count(av) /), &
702                                                 count = (/ nxr-nxl+2,         &
703                                                            nyn-nys+1, 1, 1 /) )
704                      ELSEIF ( nxr /= nx  .AND.  nyn == ny )  THEN
705                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
706                                                 id_var_do2d(av,if),           &
707                                                 local_2d(nxl:nxr,nys:nyn+1),  &
708                                                 start = (/ nxl+1, nys+1, iis, &
709                                                    do2d_xy_time_count(av) /), &
710                                                 count = (/ nxr-nxl+1,         &
711                                                            nyn-nys+2, 1, 1 /) )
712                      ELSEIF ( nxr == nx  .AND.  nyn == ny )  THEN
713                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
714                                                 id_var_do2d(av,if),           &
715                                                 local_2d(nxl:nxr+1,nys:nyn+1),&
716                                                 start = (/ nxl+1, nys+1, iis, &
717                                                    do2d_xy_time_count(av) /), &
718                                                 count = (/ nxr-nxl+2,          &
719                                                            nyn-nys+2, 1, 1 /) )
720                      ELSE
721                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
722                                                 id_var_do2d(av,if),           &
723                                                 local_2d(nxl:nxr,nys:nyn),    &
724                                                 start = (/ nxl+1, nys+1, iis, &
725                                                    do2d_xy_time_count(av) /), &
726                                                 count = (/ nxr-nxl+1,         &
727                                                            nyn-nys+1, 1, 1 /) )
728                      ENDIF
729
730                      CALL handle_netcdf_error( 'data_output_2d', 55 )
731#endif
732                   ELSE
733
734                      IF ( data_output_2d_on_each_pe )  THEN
735!
736!--                      Output of partial arrays on each PE
737#if defined( __netcdf )
738                         IF ( netcdf_output  .AND.  myid == 0 )  THEN
739                            WRITE ( 21 )  time_since_reference_point, &
740                                          do2d_xy_time_count(av), av
741                         ENDIF
742#endif
743                         DO  i = 0, io_blocks-1
744                            IF ( i == io_group )  THEN
745                               WRITE ( 21 )  nxlg, nxrg, nysg, nyng
746                               WRITE ( 21 )  local_2d
747                            ENDIF
748#if defined( __parallel )
749                            CALL MPI_BARRIER( comm2d, ierr )
750#endif
751                         ENDDO
752
753                      ELSE
754!
755!--                      PE0 receives partial arrays from all processors and
756!--                      then outputs them. Here a barrier has to be set,
757!--                      because otherwise "-MPI- FATAL: Remote protocol queue
758!--                      full" may occur.
759                         CALL MPI_BARRIER( comm2d, ierr )
760
761                         ngp = ( nxrg-nxlg+1 ) * ( nyng-nysg+1 )
762                         IF ( myid == 0 )  THEN
763!
764!--                         Local array can be relocated directly.
765                            total_2d(nxlg:nxrg,nysg:nyng) = local_2d
766!
767!--                         Receive data from all other PEs.
768                            DO  n = 1, numprocs-1
769!
770!--                            Receive index limits first, then array.
771!--                            Index limits are received in arbitrary order from
772!--                            the PEs.
773                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,    &
774                                              MPI_ANY_SOURCE, 0, comm2d, &
775                                              status, ierr )
776                               sender = status(MPI_SOURCE)
777                               DEALLOCATE( local_2d )
778                               ALLOCATE( local_2d(ind(1):ind(2),ind(3):ind(4)) )
779                               CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp,  &
780                                              MPI_REAL, sender, 1, comm2d,   &
781                                              status, ierr )
782                               total_2d(ind(1):ind(2),ind(3):ind(4)) = local_2d
783                            ENDDO
784!
785!--                         Output of the total cross-section.
786                            IF ( iso2d_output )  THEN
787                               WRITE (21)  total_2d(-nbgp:nx+nbgp,-nbgp:ny+nbgp)
788                            ENDIF
789!
790!--                         Relocate the local array for the next loop increment
791                            DEALLOCATE( local_2d )
792                            ALLOCATE( local_2d(nxlg:nxrg,nysg:nyng) )
793
794#if defined( __netcdf )
795                            IF ( netcdf_output )  THEN
796                               IF ( two_d ) THEN
797                                  nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
798                                                          id_var_do2d(av,if),  &
799                                                      total_2d(0:nx+1,0:ny+1), &
800                                start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
801                                                count = (/ nx+2, ny+2, 1, 1 /) )
802                               ELSE
803                                  nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
804                                                          id_var_do2d(av,if),  &
805                                                      total_2d(0:nx+1,0:ny+1), &
806                               start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
807                                                count = (/ nx+2, ny+2, 1, 1 /) )
808                               ENDIF
809                               CALL handle_netcdf_error( 'data_output_2d', 54 )
810                            ENDIF
811#endif
812
813                         ELSE
814!
815!--                         First send the local index limits to PE0
816                            ind(1) = nxlg; ind(2) = nxrg
817                            ind(3) = nysg; ind(4) = nyng
818                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
819                                           comm2d, ierr )
820!
821!--                         Send data to PE0
822                            CALL MPI_SEND( local_2d(nxlg,nysg), ngp, &
823                                           MPI_REAL, 0, 1, comm2d, ierr )
824                         ENDIF
825!
826!--                      A barrier has to be set, because otherwise some PEs may
827!--                      proceed too fast so that PE0 may receive wrong data on
828!--                      tag 0
829                         CALL MPI_BARRIER( comm2d, ierr )
830                      ENDIF
831
832                   ENDIF
833#else
834                   IF ( iso2d_output )  THEN
835                      WRITE (21)  local_2d(nxl:nxr+1,nys:nyn+1)
836                   ENDIF
837#if defined( __netcdf )
838                   IF ( netcdf_output )  THEN
839                      IF ( two_d ) THEN
840                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
841                                                 id_var_do2d(av,if),           &
842                                                local_2d(nxl:nxr+1,nys:nyn+1), &
843                                start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
844                                              count = (/ nx+2, ny+2, 1, 1 /) )
845                      ELSE
846                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
847                                                 id_var_do2d(av,if),           &
848                                                local_2d(nxl:nxr+1,nys:nyn+1), &
849                               start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
850                                              count = (/ nx+2, ny+2, 1, 1 /) )
851                      ENDIF
852                      CALL handle_netcdf_error( 'data_output_2d', 447 )
853                   ENDIF
854#endif
855#endif
856                   do2d_xy_n = do2d_xy_n + 1
857!
858!--                Write LOCAL parameter set for ISO2D.
859                   IF ( myid == 0  .AND.  iso2d_output )  THEN
860                      IF ( section(is,s) /= -1 )  THEN
861                         WRITE ( section_chr, '(''z = '',F7.2,'' m  (GP '',I3, &
862                                               &'')'')'                        &
863                               )  level_z(layer_xy), layer_xy
864                      ELSE
865                         section_chr = 'averaged along z'
866                      ENDIF
867                      IF ( av == 0 )  THEN
868                         rtext = TRIM( do2d(av,if) ) // '  t = ' //    &
869                                 TRIM( simulated_time_chr ) // '  ' // &
870                                 TRIM( section_chr )
871                      ELSE
872                         rtext = TRIM( do2d(av,if) ) // '  averaged t = ' // &
873                                 TRIM( simulated_time_chr ) // '  ' //       &
874                                 TRIM( section_chr )
875                      ENDIF
876                      WRITE (27,LOCAL)
877                   ENDIF
878!
879!--                For 2D-arrays (e.g. u*) only one cross-section is available.
880!--                Hence exit loop of output levels.
881                   IF ( two_d )  THEN
882                      two_d = .FALSE.
883                      EXIT loop1
884                   ENDIF
885
886                CASE ( 'xz' )
887!
888!--                Update the NetCDF xz cross section time axis
889                   IF ( myid == 0  .OR.  netcdf_data_format > 2 )  THEN
890
891                      IF ( simulated_time /= do2d_xz_last_time(av) )  THEN
892                         do2d_xz_time_count(av) = do2d_xz_time_count(av) + 1
893                         do2d_xz_last_time(av)  = simulated_time
894                         IF ( ( .NOT. data_output_2d_on_each_pe  .AND.        &
895                              netcdf_output )  .OR.  netcdf_data_format > 2 ) &
896                         THEN
897#if defined( __netcdf )
898                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
899                                                    id_var_time_xz(av),        &
900                                             (/ time_since_reference_point /), &
901                                         start = (/ do2d_xz_time_count(av) /), &
902                                                    count = (/ 1 /) )
903                            CALL handle_netcdf_error( 'data_output_2d', 56 )
904#endif
905                         ENDIF
906                      ENDIF
907
908                   ENDIF
909
910!
911!--                If required, carry out averaging along y
912                   IF ( section(is,s) == -1 )  THEN
913
914                      ALLOCATE( local_2d_l(nxlg:nxrg,nzb:nzt+1) )
915                      local_2d_l = 0.0
916                      ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
917!
918!--                   First local averaging on the PE
919                      DO  k = nzb, nzt+1
920                         DO  j = nys, nyn
921                            DO  i = nxlg, nxrg
922                               local_2d_l(i,k) = local_2d_l(i,k) + &
923                                                 local_pf(i,j,k)
924                            ENDDO
925                         ENDDO
926                      ENDDO
927#if defined( __parallel )
928!
929!--                   Now do the averaging over all PEs along y
930                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
931                      CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb),              &
932                                          local_2d(nxlg,nzb), ngp, MPI_REAL, &
933                                          MPI_SUM, comm1dy, ierr )
934#else
935                      local_2d = local_2d_l
936#endif
937                      local_2d = local_2d / ( ny + 1.0 )
938
939                      DEALLOCATE( local_2d_l )
940
941                   ELSE
942!
943!--                   Just store the respective section on the local array
944!--                   (but only if it is available on this PE!)
945                      IF ( section(is,s) >= nys  .AND.  section(is,s) <= nyn ) &
946                      THEN
947                         local_2d = local_pf(:,section(is,s),nzb:nzt+1)
948                      ENDIF
949
950                   ENDIF
951
952#if defined( __parallel )
953                   IF ( netcdf_output  .AND.  netcdf_data_format > 2 )  THEN
954!
955!--                   ATTENTION: The following lines are a workaround, because
956!--                              independet output does not work with the
957!--                              current NetCDF4 installation. Therefore, data
958!--                              are transferred from PEs having the cross
959!--                              sections to other PEs along y having no cross
960!--                              section data. Some of these data are the
961!--                              output.
962!--                   BEGIN WORKAROUND---------------------------------------
963                      IF ( npey /= 1  .AND.  section(is,s) /= -1)  THEN
964                         ALLOCATE( local_2d_l(nxlg:nxrg,nzb:nzt+1) )
965                         local_2d_l = 0.0
966                         IF ( section(is,s) >= nys .AND. section(is,s) <= nyn )&
967                         THEN
968                            local_2d_l = local_2d
969                         ENDIF
970#if defined( __parallel )
971!
972!--                      Distribute data over all PEs along y
973                         ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
974                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
975                         CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb),            &
976                                             local_2d(nxlg,nzb), ngp,         &
977                                             MPI_REAL, MPI_SUM, comm1dy, ierr )
978#else
979                         local_2d = local_2d_l
980#endif
981                         DEALLOCATE( local_2d_l )
982                      ENDIF
983!--                   END WORKAROUND-----------------------------------------
984
985!
986!--                   Output in NetCDF4/HDF5 format.
987!--                   Output only on those PEs where the respective cross
988!--                   sections reside. Cross sections averaged along y are
989!--                   output on the respective first PE along y (myidy=0).
990                      IF ( ( section(is,s) >= nys  .AND.  &
991                             section(is,s) <= nyn )  .OR.  &
992                           ( section(is,s) == -1  .AND.  myidy == 0 ) )  THEN
993!
994!--                      Do not output redundant ghost point data except for the
995!--                      boundaries of the total domain.
996#if defined( __netcdf )
997                         IF ( nxr == nx )  THEN
998                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
999                                                id_var_do2d(av,if),            &
1000                                                local_2d(nxl:nxr+1,nzb:nzt+1), &
1001                                                start = (/ nxl+1, is, 1,       &
1002                                                    do2d_xz_time_count(av) /), &
1003                                                count = (/ nxr-nxl+2, 1,       &
1004                                                           nzt+2, 1 /) )
1005                         ELSE
1006                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
1007                                                id_var_do2d(av,if),            &
1008                                                local_2d(nxl:nxr,nzb:nzt+1),   &
1009                                                start = (/ nxl+1, is, 1,       &
1010                                                    do2d_xz_time_count(av) /), &
1011                                                count = (/ nxr-nxl+1, 1,       &
1012                                                           nzt+2, 1 /) )
1013                         ENDIF
1014
1015                         CALL handle_netcdf_error( 'data_output_2d', 57 )
1016
1017                      ELSE
1018!
1019!--                      Output on other PEs. Only one point is output!!
1020!--                      ATTENTION: This is a workaround (see above)!!
1021                         IF ( npey /= 1 )  THEN
1022                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
1023                                                    id_var_do2d(av,if),        &
1024                                                    local_2d(nxl:nxl,nzb:nzb), &
1025                                                    start = (/ nxl+1, is, 1,   &
1026                                                    do2d_xz_time_count(av) /), &
1027                                                    count = (/ 1, 1, 1, 1 /) )
1028                            CALL handle_netcdf_error( 'data_output_2d', 451 )
1029                         ENDIF
1030#endif
1031                      ENDIF
1032
1033                   ELSE
1034
1035                      IF ( data_output_2d_on_each_pe )  THEN
1036!
1037!--                      Output of partial arrays on each PE. If the cross
1038!--                      section does not reside on the PE, output special
1039!--                      index values.
1040#if defined( __netcdf )
1041                         IF ( netcdf_output  .AND.  myid == 0 )  THEN
1042                            WRITE ( 22 )  time_since_reference_point, &
1043                                          do2d_xz_time_count(av), av
1044                         ENDIF
1045#endif
1046                         DO  i = 0, io_blocks-1
1047                            IF ( i == io_group )  THEN
1048                               IF ( ( section(is,s) >= nys  .AND.   &
1049                                      section(is,s) <= nyn )  .OR.  &
1050                                    ( section(is,s) == -1  .AND.    &
1051                                      nys-1 == -1 ) )               &
1052                               THEN
1053                                  WRITE (22)  nxlg, nxrg, nzb, nzt+1
1054                                  WRITE (22)  local_2d
1055                               ELSE
1056                                  WRITE (22)  -1, -1, -1, -1
1057                               ENDIF
1058                            ENDIF
1059#if defined( __parallel )
1060                            CALL MPI_BARRIER( comm2d, ierr )
1061#endif
1062                         ENDDO
1063
1064                      ELSE
1065!
1066!--                      PE0 receives partial arrays from all processors of the
1067!--                      respective cross section and outputs them. Here a
1068!--                      barrier has to be set, because otherwise
1069!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1070                         CALL MPI_BARRIER( comm2d, ierr )
1071
1072                         ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
1073                         IF ( myid == 0 )  THEN
1074!
1075!--                         Local array can be relocated directly.
1076                            IF ( ( section(is,s) >= nys  .AND.                 &
1077                                   section(is,s) <= nyn )  .OR.                &
1078                                 ( section(is,s) == -1  .AND.  nys-1 == -1 ) ) &
1079                            THEN
1080                               total_2d(nxlg:nxrg,nzb:nzt+1) = local_2d
1081                            ENDIF
1082!
1083!--                         Receive data from all other PEs.
1084                            DO  n = 1, numprocs-1
1085!
1086!--                            Receive index limits first, then array.
1087!--                            Index limits are received in arbitrary order from
1088!--                            the PEs.
1089                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,     &
1090                                              MPI_ANY_SOURCE, 0, comm2d,  &
1091                                              status, ierr )
1092!
1093!--                            Not all PEs have data for XZ-cross-section.
1094                               IF ( ind(1) /= -9999 )  THEN
1095                                  sender = status(MPI_SOURCE)
1096                                  DEALLOCATE( local_2d )
1097                                  ALLOCATE( local_2d(ind(1):ind(2), &
1098                                                     ind(3):ind(4)) )
1099                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1100                                                 MPI_REAL, sender, 1, comm2d,  &
1101                                                 status, ierr )
1102                                  total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1103                                                                        local_2d
1104                               ENDIF
1105                            ENDDO
1106!
1107!--                         Output of the total cross-section.
1108                            IF ( iso2d_output )  THEN
1109                               WRITE (22)  total_2d(-nbgp:nx+nbgp,nzb:nzt+1)
1110                            ENDIF
1111!
1112!--                         Relocate the local array for the next loop increment
1113                            DEALLOCATE( local_2d )
1114                            ALLOCATE( local_2d(nxlg:nxrg,nzb:nzt+1) )
1115
1116#if defined( __netcdf )
1117                            IF ( netcdf_output )  THEN
1118                               nc_stat = NF90_PUT_VAR( id_set_xz(av),          &
1119                                                    id_var_do2d(av,if),        &
1120                                                    total_2d(0:nx+1,nzb:nzt+1),&
1121                               start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1122                                                count = (/ nx+2, 1, nz+2, 1 /) )
1123                               CALL handle_netcdf_error( 'data_output_2d', 58 )
1124                            ENDIF
1125#endif
1126
1127                         ELSE
1128!
1129!--                         If the cross section resides on the PE, send the
1130!--                         local index limits, otherwise send -9999 to PE0.
1131                            IF ( ( section(is,s) >= nys  .AND.                 &
1132                                   section(is,s) <= nyn )  .OR.                &
1133                                 ( section(is,s) == -1  .AND.  nys-1 == -1 ) ) &
1134                            THEN
1135                               ind(1) = nxlg; ind(2) = nxrg
1136                               ind(3) = nzb;   ind(4) = nzt+1
1137                            ELSE
1138                               ind(1) = -9999; ind(2) = -9999
1139                               ind(3) = -9999; ind(4) = -9999
1140                            ENDIF
1141                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1142                                           comm2d, ierr )
1143!
1144!--                         If applicable, send data to PE0.
1145                            IF ( ind(1) /= -9999 )  THEN
1146                               CALL MPI_SEND( local_2d(nxlg,nzb), ngp, &
1147                                              MPI_REAL, 0, 1, comm2d, ierr )
1148                            ENDIF
1149                         ENDIF
1150!
1151!--                      A barrier has to be set, because otherwise some PEs may
1152!--                      proceed too fast so that PE0 may receive wrong data on
1153!--                      tag 0
1154                         CALL MPI_BARRIER( comm2d, ierr )
1155                      ENDIF
1156
1157                   ENDIF
1158#else
1159                   IF ( iso2d_output )  THEN
1160                      WRITE (22)  local_2d(nxl:nxr+1,nzb:nzt+1)
1161                   ENDIF
1162#if defined( __netcdf )
1163                   IF ( netcdf_output )  THEN
1164                      nc_stat = NF90_PUT_VAR( id_set_xz(av),                   &
1165                                              id_var_do2d(av,if),              &
1166                                              local_2d(nxl:nxr+1,nzb:nzt+1),   &
1167                               start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1168                                              count = (/ nx+2, 1, nz+2, 1 /) )
1169                      CALL handle_netcdf_error( 'data_output_2d', 451 )
1170                   ENDIF
1171#endif
1172#endif
1173                   do2d_xz_n = do2d_xz_n + 1
1174!
1175!--                Write LOCAL-parameter set for ISO2D.
1176                   IF ( myid == 0  .AND.  iso2d_output )  THEN
1177                      IF ( section(is,s) /= -1 )  THEN
1178                         WRITE ( section_chr, '(''y = '',F8.2,'' m  (GP '',I3, &
1179                                               &'')'')'                        &
1180                               )  section(is,s)*dy, section(is,s)
1181                      ELSE
1182                         section_chr = 'averaged along y'
1183                      ENDIF
1184                      IF ( av == 0 )  THEN
1185                         rtext = TRIM( do2d(av,if) ) // '  t = ' //    &
1186                                 TRIM( simulated_time_chr ) // '  ' // &
1187                                 TRIM( section_chr )
1188                      ELSE
1189                         rtext = TRIM( do2d(av,if) ) // '  averaged t = ' // &
1190                                 TRIM( simulated_time_chr ) // '  ' //       &
1191                                 TRIM( section_chr )
1192                      ENDIF
1193                      WRITE (28,LOCAL)
1194                   ENDIF
1195
1196                CASE ( 'yz' )
1197!
1198!--                Update the NetCDF yz cross section time axis
1199                   IF ( myid == 0  .OR.  netcdf_data_format > 2 )  THEN
1200
1201                      IF ( simulated_time /= do2d_yz_last_time(av) )  THEN
1202                         do2d_yz_time_count(av) = do2d_yz_time_count(av) + 1
1203                         do2d_yz_last_time(av)  = simulated_time
1204                         IF ( ( .NOT. data_output_2d_on_each_pe  .AND.        &
1205                              netcdf_output )  .OR.  netcdf_data_format > 2 ) &
1206                         THEN
1207#if defined( __netcdf )
1208                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1209                                                    id_var_time_yz(av),        &
1210                                             (/ time_since_reference_point /), &
1211                                         start = (/ do2d_yz_time_count(av) /), &
1212                                                    count = (/ 1 /) )
1213                            CALL handle_netcdf_error( 'data_output_2d', 59 )
1214#endif
1215                         ENDIF
1216                      ENDIF
1217
1218                   ENDIF
1219!
1220!--                If required, carry out averaging along x
1221                   IF ( section(is,s) == -1 )  THEN
1222
1223                      ALLOCATE( local_2d_l(nysg:nyng,nzb:nzt+1) )
1224                      local_2d_l = 0.0
1225                      ngp = ( nyng-nysg+1 ) * ( nzt-nzb+2 )
1226!
1227!--                   First local averaging on the PE
1228                      DO  k = nzb, nzt+1
1229                         DO  j = nysg, nyng
1230                            DO  i = nxl, nxr
1231                               local_2d_l(j,k) = local_2d_l(j,k) + &
1232                                                 local_pf(i,j,k)
1233                            ENDDO
1234                         ENDDO
1235                      ENDDO
1236#if defined( __parallel )
1237!
1238!--                   Now do the averaging over all PEs along x
1239                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
1240                      CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb),              &
1241                                          local_2d(nysg,nzb), ngp, MPI_REAL, &
1242                                          MPI_SUM, comm1dx, ierr )
1243#else
1244                      local_2d = local_2d_l
1245#endif
1246                      local_2d = local_2d / ( nx + 1.0 )
1247
1248                      DEALLOCATE( local_2d_l )
1249
1250                   ELSE
1251!
1252!--                   Just store the respective section on the local array
1253!--                   (but only if it is available on this PE!)
1254                      IF ( section(is,s) >= nxl  .AND.  section(is,s) <= nxr ) &
1255                      THEN
1256                         local_2d = local_pf(section(is,s),:,nzb:nzt+1)
1257                      ENDIF
1258
1259                   ENDIF
1260
1261#if defined( __parallel )
1262                   IF ( netcdf_output  .AND.  netcdf_data_format > 2 )  THEN
1263!
1264!--                   ATTENTION: The following lines are a workaround, because
1265!--                              independet output does not work with the
1266!--                              current NetCDF4 installation. Therefore, data
1267!--                              are transferred from PEs having the cross
1268!--                              sections to other PEs along y having no cross
1269!--                              section data. Some of these data are the
1270!--                              output.
1271!--                   BEGIN WORKAROUND---------------------------------------
1272                      IF ( npex /= 1  .AND.  section(is,s) /= -1)  THEN
1273                         ALLOCATE( local_2d_l(nysg:nyng,nzb:nzt+1) )
1274                         local_2d_l = 0.0
1275                         IF ( section(is,s) >= nxl .AND. section(is,s) <= nxr )&
1276                         THEN
1277                            local_2d_l = local_2d
1278                         ENDIF
1279#if defined( __parallel )
1280!
1281!--                      Distribute data over all PEs along x
1282                         ngp = ( nyng-nysg+1 ) * ( nzt-nzb + 2 )
1283                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1284                         CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb),            &
1285                                             local_2d(nysg,nzb), ngp,         &
1286                                             MPI_REAL, MPI_SUM, comm1dx, ierr )
1287#else
1288                         local_2d = local_2d_l
1289#endif
1290                         DEALLOCATE( local_2d_l )
1291                      ENDIF
1292!--                   END WORKAROUND-----------------------------------------
1293
1294!
1295!--                   Output in NetCDF4/HDF5 format.
1296!--                   Output only on those PEs where the respective cross
1297!--                   sections reside. Cross sections averaged along x are
1298!--                   output on the respective first PE along x (myidx=0).
1299                      IF ( ( section(is,s) >= nxl  .AND.  &
1300                             section(is,s) <= nxr )  .OR.  &
1301                           ( section(is,s) == -1  .AND.  myidx == 0 ) )  THEN
1302!
1303!--                      Do not output redundant ghost point data except for the
1304!--                      boundaries of the total domain.
1305#if defined( __netcdf )
1306                         IF ( nyn == ny )  THEN
1307                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1308                                                id_var_do2d(av,if),            &
1309                                                local_2d(nys:nyn+1,nzb:nzt+1), &
1310                                                start = (/ is, nys+1, 1,       &
1311                                                    do2d_yz_time_count(av) /), &
1312                                                count = (/ 1, nyn-nys+2,       &
1313                                                           nzt+2, 1 /) )
1314                         ELSE
1315                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1316                                                id_var_do2d(av,if),            &
1317                                                local_2d(nys:nyn,nzb:nzt+1),   &
1318                                                start = (/ is, nys+1, 1,       &
1319                                                    do2d_yz_time_count(av) /), &
1320                                                count = (/ 1, nyn-nys+1,       &
1321                                                           nzt+2, 1 /) )
1322                         ENDIF
1323
1324                         CALL handle_netcdf_error( 'data_output_2d', 60 )
1325
1326                      ELSE
1327!
1328!--                      Output on other PEs. Only one point is output!!
1329!--                      ATTENTION: This is a workaround (see above)!!
1330                         IF ( npex /= 1 )  THEN
1331                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1332                                                    id_var_do2d(av,if),        &
1333                                                    local_2d(nys:nys,nzb:nzb), &
1334                                                    start = (/ is, nys+1, 1,   &
1335                                                    do2d_yz_time_count(av) /), &
1336                                                    count = (/ 1, 1, 1, 1 /) )
1337                            CALL handle_netcdf_error( 'data_output_2d', 452 )
1338                         ENDIF
1339#endif
1340                      ENDIF
1341
1342                   ELSE
1343
1344                      IF ( data_output_2d_on_each_pe )  THEN
1345!
1346!--                      Output of partial arrays on each PE. If the cross
1347!--                      section does not reside on the PE, output special
1348!--                      index values.
1349#if defined( __netcdf )
1350                         IF ( netcdf_output  .AND.  myid == 0 )  THEN
1351                            WRITE ( 23 )  time_since_reference_point, &
1352                                          do2d_yz_time_count(av), av
1353                         ENDIF
1354#endif
1355                         DO  i = 0, io_blocks-1
1356                            IF ( i == io_group )  THEN
1357                               IF ( ( section(is,s) >= nxl  .AND.   &
1358                                      section(is,s) <= nxr )  .OR.  &
1359                                    ( section(is,s) == -1  .AND.    &
1360                                      nxl-1 == -1 ) )               &
1361                               THEN
1362                                  WRITE (23)  nysg, nyng, nzb, nzt+1
1363                                  WRITE (23)  local_2d
1364                               ELSE
1365                                  WRITE (23)  -1, -1, -1, -1
1366                               ENDIF
1367                            ENDIF
1368#if defined( __parallel )
1369                            CALL MPI_BARRIER( comm2d, ierr )
1370#endif
1371                         ENDDO
1372
1373                      ELSE
1374!
1375!--                      PE0 receives partial arrays from all processors of the
1376!--                      respective cross section and outputs them. Here a
1377!--                      barrier has to be set, because otherwise
1378!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1379                         CALL MPI_BARRIER( comm2d, ierr )
1380
1381                         ngp = ( nyng-nysg+1 ) * ( nzt-nzb+2 )
1382                         IF ( myid == 0 )  THEN
1383!
1384!--                         Local array can be relocated directly.
1385                            IF ( ( section(is,s) >= nxl  .AND.                 &
1386                                   section(is,s) <= nxr )   .OR.               &
1387                                 ( section(is,s) == -1  .AND.  nxl-1 == -1 ) ) &
1388                            THEN
1389                               total_2d(nysg:nyng,nzb:nzt+1) = local_2d
1390                            ENDIF
1391!
1392!--                         Receive data from all other PEs.
1393                            DO  n = 1, numprocs-1
1394!
1395!--                            Receive index limits first, then array.
1396!--                            Index limits are received in arbitrary order from
1397!--                            the PEs.
1398                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,     &
1399                                              MPI_ANY_SOURCE, 0, comm2d,  &
1400                                              status, ierr )
1401!
1402!--                            Not all PEs have data for YZ-cross-section.
1403                               IF ( ind(1) /= -9999 )  THEN
1404                                  sender = status(MPI_SOURCE)
1405                                  DEALLOCATE( local_2d )
1406                                  ALLOCATE( local_2d(ind(1):ind(2), &
1407                                                     ind(3):ind(4)) )
1408                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1409                                                 MPI_REAL, sender, 1, comm2d,  &
1410                                                 status, ierr )
1411                                  total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1412                                                                        local_2d
1413                               ENDIF
1414                            ENDDO
1415!
1416!--                         Output of the total cross-section.
1417                            IF ( iso2d_output )  THEN
1418                               WRITE (23)  total_2d(0:ny+1,nzb:nzt+1)
1419                            ENDIF
1420!
1421!--                         Relocate the local array for the next loop increment
1422                            DEALLOCATE( local_2d )
1423                            ALLOCATE( local_2d(nysg:nyng,nzb:nzt+1) )
1424
1425#if defined( __netcdf )
1426                            IF ( netcdf_output )  THEN
1427                               nc_stat = NF90_PUT_VAR( id_set_yz(av),          &
1428                                                    id_var_do2d(av,if),        &
1429                                                    total_2d(0:ny+1,nzb:nzt+1),&
1430                               start = (/ is, 1, 1, do2d_yz_time_count(av) /), &
1431                                                count = (/ 1, ny+2, nz+2, 1 /) )
1432                               CALL handle_netcdf_error( 'data_output_2d', 61 )
1433                            ENDIF
1434#endif
1435
1436                         ELSE
1437!
1438!--                         If the cross section resides on the PE, send the
1439!--                         local index limits, otherwise send -9999 to PE0.
1440                            IF ( ( section(is,s) >= nxl  .AND.                 &
1441                                   section(is,s) <= nxr )  .OR.                &
1442                                 ( section(is,s) == -1  .AND.  nxl-1 == -1 ) ) &
1443                            THEN
1444                               ind(1) = nysg; ind(2) = nyng
1445                               ind(3) = nzb;   ind(4) = nzt+1
1446                            ELSE
1447                               ind(1) = -9999; ind(2) = -9999
1448                               ind(3) = -9999; ind(4) = -9999
1449                            ENDIF
1450                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1451                                           comm2d, ierr )
1452!
1453!--                         If applicable, send data to PE0.
1454                            IF ( ind(1) /= -9999 )  THEN
1455                               CALL MPI_SEND( local_2d(nysg,nzb), ngp, &
1456                                              MPI_REAL, 0, 1, comm2d, ierr )
1457                            ENDIF
1458                         ENDIF
1459!
1460!--                      A barrier has to be set, because otherwise some PEs may
1461!--                      proceed too fast so that PE0 may receive wrong data on
1462!--                      tag 0
1463                         CALL MPI_BARRIER( comm2d, ierr )
1464                      ENDIF
1465
1466                   ENDIF
1467#else
1468                   IF ( iso2d_output )  THEN
1469                      WRITE (23)  local_2d(nys:nyn+1,nzb:nzt+1)
1470                   ENDIF
1471#if defined( __netcdf )
1472                   IF ( netcdf_output )  THEN
1473                      nc_stat = NF90_PUT_VAR( id_set_yz(av),                   &
1474                                              id_var_do2d(av,if),              &
1475                                              local_2d(nys:nyn+1,nzb:nzt+1),   &
1476                               start = (/ is, 1, 1, do2d_xz_time_count(av) /), &
1477                                              count = (/ 1, ny+2, nz+2, 1 /) )
1478                      CALL handle_netcdf_error( 'data_output_2d', 452 )
1479                   ENDIF
1480#endif
1481#endif
1482                   do2d_yz_n = do2d_yz_n + 1
1483!
1484!--                Write LOCAL-parameter set for ISO2D.
1485                   IF ( myid == 0  .AND.  iso2d_output )  THEN
1486                      IF ( section(is,s) /= -1 )  THEN
1487                         WRITE ( section_chr, '(''x = '',F8.2,'' m  (GP '',I3, &
1488                                               &'')'')'                        &
1489                               )  section(is,s)*dx, section(is,s)
1490                      ELSE
1491                         section_chr = 'averaged along x'
1492                      ENDIF
1493                      IF ( av == 0 )  THEN
1494                         rtext = TRIM( do2d(av,if) ) // '  t = ' //    &
1495                                 TRIM( simulated_time_chr ) // '  ' // &
1496                                 TRIM( section_chr )
1497                      ELSE
1498                         rtext = TRIM( do2d(av,if) ) // '  averaged t = ' // &
1499                                 TRIM( simulated_time_chr ) // '  ' //       &
1500                                 TRIM( section_chr )
1501                      ENDIF
1502                      WRITE (29,LOCAL)
1503                   ENDIF
1504
1505             END SELECT
1506
1507             is = is + 1
1508          ENDDO loop1
1509
1510       ENDIF
1511
1512       if = if + 1
1513       l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
1514       do2d_mode = do2d(av,if)(l-1:l)
1515
1516    ENDDO
1517
1518!
1519!-- Deallocate temporary arrays.
1520    IF ( ALLOCATED( level_z ) )  DEALLOCATE( level_z )
1521    DEALLOCATE( local_pf, local_2d )
1522#if defined( __parallel )
1523    IF ( .NOT.  data_output_2d_on_each_pe  .AND.  myid == 0 )  THEN
1524       DEALLOCATE( total_2d )
1525    ENDIF
1526#endif
1527
1528!
1529!-- Close plot output file.
1530    file_id = 20 + s
1531
1532    IF ( data_output_2d_on_each_pe )  THEN
1533       DO  i = 0, io_blocks-1
1534          IF ( i == io_group )  THEN
1535             CALL close_file( file_id )
1536          ENDIF
1537#if defined( __parallel )
1538          CALL MPI_BARRIER( comm2d, ierr )
1539#endif
1540       ENDDO
1541    ELSE
1542       IF ( myid == 0 )  CALL close_file( file_id )
1543    ENDIF
1544
1545
1546    CALL cpu_log (log_point(3),'data_output_2d','stop','nobarrier')
1547
1548 END SUBROUTINE data_output_2d
Note: See TracBrowser for help on using the repository browser.