Home

Context Navigation

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

Last change on this file since 1036 was 1036, checked in by raasch, 11 years ago
code has been put under the GNU General Public License (v3)
Property svn:keywords set to `Id`
File size: 68.2 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |