Home

Context Navigation

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

Last change on this file since 760 was 760, checked in by raasch, 13 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 65.1 KB

Rev	Line
[1]	1	SUBROUTINE data_output_2d( mode, av )
	2
	3	!------------------------------------------------------------------------------!
[254]	4	! Current revisions:
[1]	5	! -----------------
[760]	6	!
[674]	7	!
[392]	8	! Former revisions:
	9	! -----------------
	10	! $Id: data_output_2d.f90 760 2011-09-15 14:37:54Z raasch $
	11	!
[760]	12	! 759 2011-09-15 13:58:31Z raasch
	13	! Splitting of parallel I/O
	14	!
[730]	15	! 729 2011-05-26 10:33:34Z heinze
	16	! Exchange ghost layers for p regardless of used pressure solver (except SOR).
	17	!
[692]	18	! 691 2011-03-04 08:45:30Z maronga
	19	! Replaced simulated_time by time_since_reference_point
	20	!
[674]	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	!
[668]	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	!
[623]	30	! 622 2010-12-10 08:08:13Z raasch
	31	! optional barriers included in order to speed up collective operations
	32	!
[494]	33	! 493 2010-03-01 08:30:24Z raasch
	34	! NetCDF4 support (parallel output)
	35	!
[392]	36	! 367 2009-08-25 08:35:52Z maronga
[291]	37	! simulated_time in NetCDF output replaced by time_since_reference_point.
[263]	38	! Output of NetCDF messages with aid of message handling routine.
[336]	39	! Bugfix: averaging along z is not allowed for 2d quantities (e.g. u* and z0)
[254]	40	! Output of messages replaced by message handling routine.
[343]	41	! Output of user defined 2D (XY) arrays at z=nzb+1 is now possible
[355]	42	! Bugfix: to_be_resorted => s_av for time-averaged scalars
[367]	43	! Calculation of shf* and qsws* added.
[1]	44	!
[226]	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	!
[98]	49	! 96 2007-06-04 08:07:41Z raasch
	50	! Output of density and salinity
	51	!
[77]	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	!
[3]	56	! RCS Log replace by Id keyword, revision history cleaned up
	57	!
[1]	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
[493]	99	INTEGER :: av, ngp, file_id, i, if, is, iis, j, k, l, layer_xy, n, psi, &
[559]	100	s, sender, &
[1]	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
[667]	133	ALLOCATE( level_z(nzb:nzt+1), local_2d(nxlg:nxrg,nysg:nyng) )
[1]	134
[493]	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
[1]	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 )
[667]	149	ALLOCATE( total_2d(-nbgp:nx+nbgp,-nbgp:ny+nbgp) )
[1]	150	#endif
	151	ENDIF
	152	ENDIF
	153
	154	CASE ( 'xz' )
	155	s = 2
[667]	156	ALLOCATE( local_2d(nxlg:nxrg,nzb:nzt+1) )
[1]	157
[493]	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
[1]	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 )
[667]	172	ALLOCATE( total_2d(-nbgp:nx+nbgp,nzb:nzt+1) )
[1]	173	#endif
	174	ENDIF
	175	ENDIF
	176
	177	CASE ( 'yz' )
	178
	179	s = 3
[667]	180	ALLOCATE( local_2d(nysg:nyng,nzb:nzt+1) )
[1]	181
[493]	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
[1]	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 )
[667]	196	ALLOCATE( total_2d(-nbgp:ny+nbgp,nzb:nzt+1) )
[1]	197	#endif
	198	ENDIF
	199	ENDIF
	200
	201	CASE DEFAULT
	202
[254]	203	message_string = 'unknown cross-section: ' // TRIM( mode )
	204	CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
[1]	205
	206	END SELECT
	207
	208	!
	209	!-- Allocate a temporary array for resorting (kji -> ijk).
[667]	210	ALLOCATE( local_pf(nxlg:nxrg,nysg:nyng,nzb:nzt+1) )
[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
[667]	237	DO i = nxlg, nxrg
	238	DO j = nysg, nyng
[1]	239	local_pf(i,j,nzb+1) = SUM( ql(nzb:nzt,j,i) * &
	240	dzw(1:nzt+1) )
	241	ENDDO
	242	ENDDO
	243	ELSE
[667]	244	DO i = nxlg, nxrg
	245	DO j = nysg, nyng
[1]	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
[729]	256	IF ( psolver /= 'sor' ) CALL exchange_horiz( p, nbgp )
[1]	257	to_be_resorted => p
	258	ELSE
[729]	259	IF ( psolver /= 'sor' ) CALL exchange_horiz( p_av, nbgp )
[1]	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
[215]	266	IF ( simulated_time >= particle_advection_start ) THEN
	267	tend = prt_count
[667]	268	CALL exchange_horiz( tend, nbgp )
[215]	269	ELSE
	270	tend = 0.0
	271	ENDIF
[667]	272	DO i = nxlg, nxrg
	273	DO j = nysg, nyng
[1]	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
[667]	281	CALL exchange_horiz( pc_av, nbgp )
[1]	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
[215]	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
[1]	304	ENDDO
	305	ENDDO
	306	ENDDO
[667]	307	CALL exchange_horiz( tend, nbgp )
[215]	308	ELSE
	309	tend = 0.0
[667]	310	END IF
	311	DO i = nxlg, nxrg
	312	DO j = nysg, nyng
[1]	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
[667]	320	CALL exchange_horiz( pr_av, nbgp )
[1]	321	to_be_resorted => pr_av
	322	ENDIF
	323
[72]	324	CASE ( 'pra*_xy' ) ! 2d-array / integral quantity => no av
	325	CALL exchange_horiz_2d( precipitation_amount )
[667]	326	DO i = nxlg, nxrg
	327	DO j = nysg, nyng
[72]	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 )
[667]	339	DO i = nxlg, nxrg
	340	DO j = nysg, nyng
[72]	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 )
[667]	346	DO i = nxlg, nxrg
	347	DO j = nysg, nyng
[72]	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
[1]	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
[667]	361	DO i = nxlg, nxrg
	362	DO j = nysg, nyng
[1]	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
[354]	417	CASE ( 'qsws*_xy' ) ! 2d-array
	418	IF ( av == 0 ) THEN
[667]	419	DO i = nxlg, nxrg
	420	DO j = nysg, nyng
[354]	421	local_pf(i,j,nzb+1) = qsws(j,i)
	422	ENDDO
	423	ENDDO
	424	ELSE
[667]	425	DO i = nxlg, nxrg
	426	DO j = nysg, nyng
[354]	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
[1]	435	CASE ( 'qv_xy', 'qv_xz', 'qv_yz' )
	436	IF ( av == 0 ) THEN
[667]	437	DO i = nxlg, nxrg
	438	DO j = nysg, nyng
[1]	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
[96]	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
[1]	457	CASE ( 's_xy', 's_xz', 's_yz' )
	458	IF ( av == 0 ) THEN
	459	to_be_resorted => q
	460	ELSE
[355]	461	to_be_resorted => s_av
[1]	462	ENDIF
	463
[96]	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
[354]	471	CASE ( 'shf*_xy' ) ! 2d-array
	472	IF ( av == 0 ) THEN
[667]	473	DO i = nxlg, nxrg
	474	DO j = nysg, nyng
[354]	475	local_pf(i,j,nzb+1) = shf(j,i)
	476	ENDDO
	477	ENDDO
	478	ELSE
[667]	479	DO i = nxlg, nxrg
	480	DO j = nysg, nyng
[354]	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
[1]	489	CASE ( 't*_xy' ) ! 2d-array
	490	IF ( av == 0 ) THEN
[667]	491	DO i = nxlg, nxrg
	492	DO j = nysg, nyng
[1]	493	local_pf(i,j,nzb+1) = ts(j,i)
	494	ENDDO
	495	ENDDO
	496	ELSE
[667]	497	DO i = nxlg, nxrg
	498	DO j = nysg, nyng
[1]	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
[667]	523	DO i = nxlg, nxrg
	524	DO j = nysg, nyng
[1]	525	local_pf(i,j,nzb+1) = us(j,i)
	526	ENDDO
	527	ENDDO
	528	ELSE
[667]	529	DO i = nxlg, nxrg
	530	DO j = nysg, nyng
[1]	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
[72]	569	CASE ( 'z0*_xy' ) ! 2d-array
	570	IF ( av == 0 ) THEN
[667]	571	DO i = nxlg, nxrg
	572	DO j = nysg, nyng
[72]	573	local_pf(i,j,nzb+1) = z0(j,i)
	574	ENDDO
	575	ENDDO
	576	ELSE
[667]	577	DO i = nxlg, nxrg
	578	DO j = nysg, nyng
[72]	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
[1]	587	CASE DEFAULT
	588	!
	589	!-- User defined quantity
	590	CALL user_data_output_2d( av, do2d(av,if), found, grid, &
[343]	591	local_pf, two_d )
[1]	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
[343]	598	ELSEIF ( grid == 'zu1' ) THEN
	599	IF ( mode == 'xy' ) level_z(nzb+1) = zu(nzb+1)
[1]	600	ENDIF
	601
	602	IF ( .NOT. found ) THEN
[274]	603	message_string = 'no output provided for: ' // &
	604	TRIM( do2d(av,if) )
[254]	605	CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
[1]	606	ENDIF
	607
	608	END SELECT
	609
	610	!
	611	!-- Resort the array to be output, if not done above
	612	IF ( .NOT. resorted ) THEN
[667]	613	DO i = nxlg, nxrg
	614	DO j = nysg, nyng
[1]	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
[493]	641	IF ( myid == 0 .OR. netcdf_data_format > 2 ) THEN
[1]	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
[493]	645	IF ( ( .NOT. data_output_2d_on_each_pe .AND. &
	646	netcdf_output ) .OR. netcdf_data_format > 2 ) &
	647	THEN
[1]	648	#if defined( __netcdf )
	649	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	650	id_var_time_xy(av), &
[291]	651	(/ time_since_reference_point /), &
[1]	652	start = (/ do2d_xy_time_count(av) /), &
	653	count = (/ 1 /) )
[493]	654	CALL handle_netcdf_error( 'data_output_2d', 53 )
[1]	655	#endif
	656	ENDIF
	657	ENDIF
	658	ENDIF
	659	!
	660	!-- If required, carry out averaging along z
[336]	661	IF ( section(is,s) == -1 .AND. .NOT. two_d ) THEN
[1]	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]	667	DO j = nysg, nyng
	668	DO i = nxlg, nxrg
[1]	669	local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
	670	ENDDO
	671	ENDDO
	672	ENDDO
	673
[667]	674	local_2d = local_2d / ( nzt -nzb + 2.0)
[1]	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 )
[493]	684	IF ( netcdf_output .AND. netcdf_data_format > 2 ) THEN
[1]	685	!
[493]	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
[1]	695	#if defined( __netcdf )
[493]	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 /) )
[1]	728	ENDIF
[493]	729
	730	CALL handle_netcdf_error( 'data_output_2d', 55 )
[1]	731	#endif
[493]	732	ELSE
[1]	733
[493]	734	IF ( data_output_2d_on_each_pe ) THEN
[1]	735	!
[493]	736	!-- Output of partial arrays on each PE
	737	#if defined( __netcdf )
	738	IF ( netcdf_output .AND. myid == 0 ) THEN
[691]	739	WRITE ( 21 ) time_since_reference_point, &
[493]	740	do2d_xy_time_count(av), av
	741	ENDIF
	742	#endif
[759]	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
[559]	752
[493]	753	ELSE
[1]	754	!
[493]	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
[667]	761	ngp = ( nxrg-nxlg+1 ) * ( nyng-nysg+1 )
[493]	762	IF ( myid == 0 ) THEN
[1]	763	!
[493]	764	!-- Local array can be relocated directly.
[667]	765	total_2d(nxlg:nxrg,nysg:nyng) = local_2d
[1]	766	!
[493]	767	!-- Receive data from all other PEs.
	768	DO n = 1, numprocs-1
[1]	769	!
[493]	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
[1]	784	!
[493]	785	!-- Output of the total cross-section.
	786	IF ( iso2d_output ) THEN
[667]	787	WRITE (21) total_2d(-nbgp:nx+nbgp,-nbgp:ny+nbgp)
[493]	788	ENDIF
	789	!
	790	!-- Relocate the local array for the next loop increment
	791	DEALLOCATE( local_2d )
[667]	792	ALLOCATE( local_2d(nxlg:nxrg,nysg:nyng) )
[1]	793
	794	#if defined( __netcdf )
[493]	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), &
[1]	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 /) )
[493]	802	ELSE
	803	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	804	id_var_do2d(av,if), &
[1]	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 /) )
[493]	808	ENDIF
	809	CALL handle_netcdf_error( 'data_output_2d', 54 )
[1]	810	ENDIF
	811	#endif
	812
[493]	813	ELSE
[1]	814	!
[493]	815	!-- First send the local index limits to PE0
[667]	816	ind(1) = nxlg; ind(2) = nxrg
	817	ind(3) = nysg; ind(4) = nyng
[493]	818	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
	819	comm2d, ierr )
[1]	820	!
[493]	821	!-- Send data to PE0
[667]	822	CALL MPI_SEND( local_2d(nxlg,nysg), ngp, &
[493]	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 )
[1]	830	ENDIF
[493]	831
[1]	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
[493]	852	CALL handle_netcdf_error( 'data_output_2d', 447 )
[1]	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	!
[108]	888	!-- Update the NetCDF xz cross section time axis
[493]	889	IF ( myid == 0 .OR. netcdf_data_format > 2 ) THEN
	890
[1]	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
[493]	894	IF ( ( .NOT. data_output_2d_on_each_pe .AND. &
	895	netcdf_output ) .OR. netcdf_data_format > 2 ) &
	896	THEN
[1]	897	#if defined( __netcdf )
	898	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
	899	id_var_time_xz(av), &
[291]	900	(/ time_since_reference_point /), &
[1]	901	start = (/ do2d_xz_time_count(av) /), &
	902	count = (/ 1 /) )
[493]	903	CALL handle_netcdf_error( 'data_output_2d', 56 )
[1]	904	#endif
	905	ENDIF
	906	ENDIF
[493]	907
[1]	908	ENDIF
[667]	909
[1]	910	!
	911	!-- If required, carry out averaging along y
	912	IF ( section(is,s) == -1 ) THEN
	913
[667]	914	ALLOCATE( local_2d_l(nxlg:nxrg,nzb:nzt+1) )
[1]	915	local_2d_l = 0.0
[667]	916	ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
[1]	917	!
	918	!-- First local averaging on the PE
	919	DO k = nzb, nzt+1
	920	DO j = nys, nyn
[667]	921	DO i = nxlg, nxrg
[1]	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
[622]	930	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[667]	931	CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb), &
	932	local_2d(nxlg,nzb), ngp, MPI_REAL, &
[1]	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 )
[493]	953	IF ( netcdf_output .AND. netcdf_data_format > 2 ) THEN
[1]	954	!
[493]	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
[667]	964	ALLOCATE( local_2d_l(nxlg:nxrg,nzb:nzt+1) )
[493]	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
[667]	973	ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
[622]	974	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[667]	975	CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb), &
	976	local_2d(nxlg,nzb), ngp, &
[493]	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.
[1]	996	#if defined( __netcdf )
[493]	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
[1]	1030	#endif
	1031	ENDIF
	1032
	1033	ELSE
	1034
[493]	1035	IF ( data_output_2d_on_each_pe ) THEN
[1]	1036	!
[493]	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
[691]	1042	WRITE ( 22 ) time_since_reference_point, &
[493]	1043	do2d_xz_time_count(av), av
	1044	ENDIF
	1045	#endif
[759]	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
[493]	1063
	1064	ELSE
[1]	1065	!
[493]	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
[667]	1072	ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
[493]	1073	IF ( myid == 0 ) THEN
[1]	1074	!
[493]	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
[667]	1080	total_2d(nxlg:nxrg,nzb:nzt+1) = local_2d
[493]	1081	ENDIF
[1]	1082	!
[493]	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, &
[1]	1091	status, ierr )
[493]	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
[667]	1109	WRITE (22) total_2d(-nbgp:nx+nbgp,nzb:nzt+1)
[1]	1110	ENDIF
	1111	!
[493]	1112	!-- Relocate the local array for the next loop increment
	1113	DEALLOCATE( local_2d )
[667]	1114	ALLOCATE( local_2d(nxlg:nxrg,nzb:nzt+1) )
[1]	1115
	1116	#if defined( __netcdf )
[493]	1117	IF ( netcdf_output ) THEN
	1118	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
[1]	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 /) )
[493]	1123	CALL handle_netcdf_error( 'data_output_2d', 58 )
	1124	ENDIF
[1]	1125	#endif
	1126
[493]	1127	ELSE
[1]	1128	!
[493]	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
[667]	1135	ind(1) = nxlg; ind(2) = nxrg
[493]	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
[667]	1146	CALL MPI_SEND( local_2d(nxlg,nzb), ngp, &
[493]	1147	MPI_REAL, 0, 1, comm2d, ierr )
	1148	ENDIF
[1]	1149	ENDIF
	1150	!
[493]	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 )
[1]	1155	ENDIF
[493]	1156
[1]	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 /) )
[493]	1169	CALL handle_netcdf_error( 'data_output_2d', 451 )
[1]	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	!
[493]	1198	!-- Update the NetCDF yz cross section time axis
	1199	IF ( myid == 0 .OR. netcdf_data_format > 2 ) THEN
	1200
[1]	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
[493]	1204	IF ( ( .NOT. data_output_2d_on_each_pe .AND. &
	1205	netcdf_output ) .OR. netcdf_data_format > 2 ) &
	1206	THEN
[1]	1207	#if defined( __netcdf )
	1208	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
	1209	id_var_time_yz(av), &
[291]	1210	(/ time_since_reference_point /), &
[1]	1211	start = (/ do2d_yz_time_count(av) /), &
	1212	count = (/ 1 /) )
[263]	1213	CALL handle_netcdf_error( 'data_output_2d', 59 )
[1]	1214	#endif
	1215	ENDIF
	1216	ENDIF
[493]	1217
[1]	1218	ENDIF
	1219	!
	1220	!-- If required, carry out averaging along x
	1221	IF ( section(is,s) == -1 ) THEN
	1222
[667]	1223	ALLOCATE( local_2d_l(nysg:nyng,nzb:nzt+1) )
[1]	1224	local_2d_l = 0.0
[667]	1225	ngp = ( nyng-nysg+1 ) * ( nzt-nzb+2 )
[1]	1226	!
	1227	!-- First local averaging on the PE
	1228	DO k = nzb, nzt+1
[667]	1229	DO j = nysg, nyng
[1]	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
[622]	1239	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[667]	1240	CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb), &
	1241	local_2d(nysg,nzb), ngp, MPI_REAL, &
[1]	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 )
[493]	1262	IF ( netcdf_output .AND. netcdf_data_format > 2 ) THEN
[1]	1263	!
[493]	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
[667]	1273	ALLOCATE( local_2d_l(nysg:nyng,nzb:nzt+1) )
[493]	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
[667]	1282	ngp = ( nyng-nysg+1 ) * ( nzt-nzb + 2 )
[622]	1283	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[667]	1284	CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb), &
	1285	local_2d(nysg,nzb), ngp, &
[493]	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.
[1]	1305	#if defined( __netcdf )
[493]	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
[1]	1339	#endif
	1340	ENDIF
	1341
	1342	ELSE
	1343
[493]	1344	IF ( data_output_2d_on_each_pe ) THEN
[1]	1345	!
[493]	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
[691]	1351	WRITE ( 23 ) time_since_reference_point, &
[493]	1352	do2d_yz_time_count(av), av
	1353	ENDIF
	1354	#endif
[759]	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
[493]	1372
	1373	ELSE
[1]	1374	!
[493]	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
[667]	1381	ngp = ( nyng-nysg+1 ) * ( nzt-nzb+2 )
[493]	1382	IF ( myid == 0 ) THEN
[1]	1383	!
[493]	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
[667]	1389	total_2d(nysg:nyng,nzb:nzt+1) = local_2d
[493]	1390	ENDIF
[1]	1391	!
[493]	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, &
[1]	1400	status, ierr )
[493]	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)
[1]	1419	ENDIF
	1420	!
[493]	1421	!-- Relocate the local array for the next loop increment
	1422	DEALLOCATE( local_2d )
[667]	1423	ALLOCATE( local_2d(nysg:nyng,nzb:nzt+1) )
[1]	1424
	1425	#if defined( __netcdf )
[493]	1426	IF ( netcdf_output ) THEN
	1427	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
[1]	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 /) )
[493]	1432	CALL handle_netcdf_error( 'data_output_2d', 61 )
	1433	ENDIF
[1]	1434	#endif
	1435
[493]	1436	ELSE
[1]	1437	!
[493]	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
[667]	1444	ind(1) = nysg; ind(2) = nyng
[493]	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
[667]	1455	CALL MPI_SEND( local_2d(nysg,nzb), ngp, &
[493]	1456	MPI_REAL, 0, 1, comm2d, ierr )
	1457	ENDIF
[1]	1458	ENDIF
	1459	!
[493]	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 )
[1]	1464	ENDIF
[493]	1465
[1]	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 /) )
[493]	1478	CALL handle_netcdf_error( 'data_output_2d', 452 )
[1]	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
[759]	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
[1]	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.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |