Home

Context Navigation

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

Last change on this file since 712 was 692, checked in by maronga, 14 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 63.8 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |