Home

Context Navigation

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

Last change on this file since 2232 was 2232, checked in by suehring, 7 years ago
Adjustments according new topography and surface-modelling concept implemented
Property svn:keywords set to `Id`
File size: 84.3 KB

Rev	Line
[1682]	1	!> @file data_output_2d.f90
[2000]	2	!------------------------------------------------------------------------------!
[1036]	3	! This file is part of PALM.
	4	!
[2000]	5	! PALM is free software: you can redistribute it and/or modify it under the
	6	! terms of the GNU General Public License as published by the Free Software
	7	! Foundation, either version 3 of the License, or (at your option) any later
	8	! version.
[1036]	9	!
	10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
	11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
	12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
	13	!
	14	! You should have received a copy of the GNU General Public License along with
	15	! PALM. If not, see <http://www.gnu.org/licenses/>.
	16	!
[2101]	17	! Copyright 1997-2017 Leibniz Universitaet Hannover
[2000]	18	!------------------------------------------------------------------------------!
[1036]	19	!
[254]	20	! Current revisions:
[1]	21	! -----------------
[2232]	22	! Adjustments to new surface concept
[1961]	23	!
[2191]	24	!
[1552]	25	! Former revisions:
	26	! -----------------
	27	! $Id: data_output_2d.f90 2232 2017-05-30 17:47:52Z suehring $
	28	!
[2191]	29	! 2190 2017-03-21 12:16:43Z raasch
	30	! bugfix for r2031: string rho replaced by rho_ocean
	31	!
[2032]	32	! 2031 2016-10-21 15:11:58Z knoop
	33	! renamed variable rho to rho_ocean and rho_av to rho_ocean_av
	34	!
[2001]	35	! 2000 2016-08-20 18:09:15Z knoop
	36	! Forced header and separation lines into 80 columns
	37	!
[1981]	38	! 1980 2016-07-29 15:51:57Z suehring
	39	! Bugfix, in order to steer user-defined output, setting flag found explicitly
	40	! to .F.
	41	!
[1977]	42	! 1976 2016-07-27 13:28:04Z maronga
	43	! Output of radiation quantities is now done directly in the respective module
	44	!
[1973]	45	! 1972 2016-07-26 07:52:02Z maronga
[1976]	46	! Output of land surface quantities is now done directly in the respective
	47	! module
[1973]	48	!
[1961]	49	! 1960 2016-07-12 16:34:24Z suehring
	50	! Scalar surface flux added
	51	! Rename INTEGER variable s into s_ind, as s is already assigned to scalar
	52	!
[1851]	53	! 1849 2016-04-08 11:33:18Z hoffmann
	54	! precipitation_amount, precipitation_rate, prr moved to arrays_3d
	55	!
[1823]	56	! 1822 2016-04-07 07:49:42Z hoffmann
	57	! Output of bulk cloud physics simplified.
	58	!
[1789]	59	! 1788 2016-03-10 11:01:04Z maronga
	60	! Added output of z0q
	61	!
[1784]	62	! 1783 2016-03-06 18:36:17Z raasch
	63	! name change of netcdf routines and module + related changes
	64	!
[1746]	65	! 1745 2016-02-05 13:06:51Z gronemeier
	66	! Bugfix: test if time axis limit exceeds moved to point after call of check_open
	67	!
[1704]	68	! 1703 2015-11-02 12:38:44Z raasch
	69	! bugfix for output of single (*) xy-sections in case of parallel netcdf I/O
	70	!
[1702]	71	! 1701 2015-11-02 07:43:04Z maronga
	72	! Bugfix in output of RRTGM data
	73	!
[1692]	74	! 1691 2015-10-26 16:17:44Z maronga
	75	! Added output of Obukhov length (ol) and radiative heating rates for RRTMG.
	76	! Formatting corrections.
	77	!
[1683]	78	! 1682 2015-10-07 23:56:08Z knoop
	79	! Code annotations made doxygen readable
	80	!
[1586]	81	! 1585 2015-04-30 07:05:52Z maronga
	82	! Added support for RRTMG
	83	!
[1556]	84	! 1555 2015-03-04 17:44:27Z maronga
	85	! Added output of r_a and r_s
	86	!
[1552]	87	! 1551 2015-03-03 14:18:16Z maronga
[1551]	88	! Added suppport for land surface model and radiation model output. In the course
	89	! of this action, the limits for vertical loops have been changed (from nzb and
	90	! nzt+1 to nzb_do and nzt_do, respectively in order to allow soil model output).
	91	! Moreover, a new vertical grid zs was introduced.
[1329]	92	!
[1360]	93	! 1359 2014-04-11 17:15:14Z hoffmann
	94	! New particle structure integrated.
	95	!
[1354]	96	! 1353 2014-04-08 15:21:23Z heinze
	97	! REAL constants provided with KIND-attribute
	98	!
[1329]	99	! 1327 2014-03-21 11:00:16Z raasch
	100	! parts concerning iso2d output removed,
	101	! -netcdf output queries
	102	!
[1321]	103	! 1320 2014-03-20 08:40:49Z raasch
[1320]	104	! ONLY-attribute added to USE-statements,
	105	! kind-parameters added to all INTEGER and REAL declaration statements,
	106	! kinds are defined in new module kinds,
	107	! revision history before 2012 removed,
	108	! comment fields (!:) to be used for variable explanations added to
	109	! all variable declaration statements
[1309]	110	!
[1319]	111	! 1318 2014-03-17 13:35:16Z raasch
	112	! barrier argument removed from cpu_log.
	113	! module interfaces removed
	114	!
[1312]	115	! 1311 2014-03-14 12:13:39Z heinze
	116	! bugfix: close #if defined( __netcdf )
	117	!
[1309]	118	! 1308 2014-03-13 14:58:42Z fricke
[1308]	119	! +local_2d_sections, local_2d_sections_l, ns
	120	! Check, if the limit of the time dimension is exceeded for parallel output
	121	! To increase the performance for parallel output, the following is done:
	122	! - Update of time axis is only done by PE0
	123	! - Cross sections are first stored on a local array and are written
	124	! collectively to the output file by all PEs.
[674]	125	!
[1116]	126	! 1115 2013-03-26 18:16:16Z hoffmann
	127	! ql is calculated by calc_liquid_water_content
	128	!
[1077]	129	! 1076 2012-12-05 08:30:18Z hoffmann
	130	! Bugfix in output of ql
	131	!
[1066]	132	! 1065 2012-11-22 17:42:36Z hoffmann
	133	! Bugfix: Output of cross sections of ql
	134	!
[1054]	135	! 1053 2012-11-13 17:11:03Z hoffmann
	136	! +qr, nr, qc and cross sections
	137	!
[1037]	138	! 1036 2012-10-22 13:43:42Z raasch
	139	! code put under GPL (PALM 3.9)
	140	!
[1035]	141	! 1031 2012-10-19 14:35:30Z raasch
	142	! netCDF4 without parallel file support implemented
	143	!
[1008]	144	! 1007 2012-09-19 14:30:36Z franke
	145	! Bugfix: missing calculation of ql_vp added
	146	!
[979]	147	! 978 2012-08-09 08:28:32Z fricke
	148	! +z0h
	149	!
[1]	150	! Revision 1.1 1997/08/11 06:24:09 raasch
	151	! Initial revision
	152	!
	153	!
	154	! Description:
	155	! ------------
[1682]	156	!> Data output of horizontal cross-sections in netCDF format or binary format
	157	!> compatible to old graphic software iso2d.
	158	!> Attention: The position of the sectional planes is still not always computed
	159	!> --------- correctly. (zu is used always)!
[1]	160	!------------------------------------------------------------------------------!
[1682]	161	SUBROUTINE data_output_2d( mode, av )
	162
[1]	163
[1320]	164	USE arrays_3d, &
[2232]	165	ONLY: dzw, e, nr, p, pt, precipitation_amount, precipitation_rate, &
	166	prr, q, qc, ql, ql_c, ql_v, ql_vp, qr, rho_ocean, s, sa, &
	167	tend, u, v, vpt, w, zu, zw
[1320]	168
[1]	169	USE averaging
[1320]	170
	171	USE cloud_parameters, &
[1849]	172	ONLY: hyrho, l_d_cp, pt_d_t
[1320]	173
	174	USE control_parameters, &
	175	ONLY: cloud_physics, data_output_2d_on_each_pe, data_output_xy, &
	176	data_output_xz, data_output_yz, do2d, &
	177	do2d_xy_last_time, do2d_xy_n, do2d_xy_time_count, &
	178	do2d_xz_last_time, do2d_xz_n, do2d_xz_time_count, &
	179	do2d_yz_last_time, do2d_yz_n, do2d_yz_time_count, &
[2232]	180	ibc_uv_b, io_blocks, io_group, land_surface, message_string, &
[1822]	181	ntdim_2d_xy, ntdim_2d_xz, ntdim_2d_yz, &
	182	psolver, section, simulated_time, simulated_time_chr, &
	183	time_since_reference_point
[1320]	184
	185	USE cpulog, &
	186	ONLY: cpu_log, log_point
	187
	188	USE grid_variables, &
	189	ONLY: dx, dy
	190
	191	USE indices, &
	192	ONLY: nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, &
	193	nz, nzb, nzt
	194
	195	USE kinds
[1551]	196
	197	USE land_surface_model_mod, &
[2232]	198	ONLY: lsm_data_output_2d, zs
[1551]	199
[1783]	200	#if defined( __netcdf )
	201	USE NETCDF
	202	#endif
[1320]	203
[1783]	204	USE netcdf_interface, &
	205	ONLY: id_set_xy, id_set_xz, id_set_yz, id_var_do2d, id_var_time_xy, &
	206	id_var_time_xz, id_var_time_yz, nc_stat, netcdf_data_format, &
	207	netcdf_handle_error
	208
[1320]	209	USE particle_attributes, &
[1359]	210	ONLY: grid_particles, number_of_particles, particle_advection_start, &
	211	particles, prt_count
[1320]	212
[1]	213	USE pegrid
	214
[1551]	215	USE radiation_model_mod, &
[1976]	216	ONLY: radiation, radiation_data_output_2d
[1551]	217
[2232]	218	USE surface_mod, &
	219	ONLY: surf_def_h, surf_lsm_h, surf_usm_h
	220
[1]	221	IMPLICIT NONE
	222
[1682]	223	CHARACTER (LEN=2) :: do2d_mode !<
	224	CHARACTER (LEN=2) :: mode !<
	225	CHARACTER (LEN=4) :: grid !<
	226	CHARACTER (LEN=25) :: section_chr !<
	227	CHARACTER (LEN=50) :: rtext !<
[1320]	228
[1682]	229	INTEGER(iwp) :: av !<
	230	INTEGER(iwp) :: ngp !<
	231	INTEGER(iwp) :: file_id !<
	232	INTEGER(iwp) :: i !<
	233	INTEGER(iwp) :: if !<
	234	INTEGER(iwp) :: is !<
	235	INTEGER(iwp) :: iis !<
	236	INTEGER(iwp) :: j !<
	237	INTEGER(iwp) :: k !<
	238	INTEGER(iwp) :: l !<
	239	INTEGER(iwp) :: layer_xy !<
[2232]	240	INTEGER(iwp) :: m !<
[1682]	241	INTEGER(iwp) :: n !<
[1703]	242	INTEGER(iwp) :: nis !<
[1682]	243	INTEGER(iwp) :: ns !<
	244	INTEGER(iwp) :: nzb_do !< lower limit of the data field (usually nzb)
	245	INTEGER(iwp) :: nzt_do !< upper limit of the data field (usually nzt+1)
	246	INTEGER(iwp) :: psi !<
[1960]	247	INTEGER(iwp) :: s_ind !<
[1682]	248	INTEGER(iwp) :: sender !<
	249	INTEGER(iwp) :: ind(4) !<
[1320]	250
[1682]	251	LOGICAL :: found !<
	252	LOGICAL :: resorted !<
	253	LOGICAL :: two_d !<
[1320]	254
[1682]	255	REAL(wp) :: mean_r !<
	256	REAL(wp) :: s_r2 !<
	257	REAL(wp) :: s_r3 !<
[1320]	258
[1682]	259	REAL(wp), DIMENSION(:), ALLOCATABLE :: level_z !<
	260	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: local_2d !<
	261	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: local_2d_l !<
[2232]	262	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tmp_2d !< temporary field used to exchange surface-related quantities
	263
[1682]	264	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: local_pf !<
	265	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: local_2d_sections !<
	266	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: local_2d_sections_l !<
[1359]	267
[1]	268	#if defined( __parallel )
[1682]	269	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: total_2d !<
[1]	270	#endif
[1682]	271	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !<
[1]	272
	273	NAMELIST /LOCAL/ rtext
	274
	275	!
	276	!-- Immediate return, if no output is requested (no respective sections
	277	!-- found in parameter data_output)
	278	IF ( mode == 'xy' .AND. .NOT. data_output_xy(av) ) RETURN
	279	IF ( mode == 'xz' .AND. .NOT. data_output_xz(av) ) RETURN
	280	IF ( mode == 'yz' .AND. .NOT. data_output_yz(av) ) RETURN
	281
[1308]	282	CALL cpu_log (log_point(3),'data_output_2d','start')
	283
[1]	284	two_d = .FALSE. ! local variable to distinguish between output of pure 2D
	285	! arrays and cross-sections of 3D arrays.
	286
	287	!
	288	!-- Depending on the orientation of the cross-section, the respective output
	289	!-- files have to be opened.
	290	SELECT CASE ( mode )
	291
	292	CASE ( 'xy' )
[1960]	293	s_ind = 1
[667]	294	ALLOCATE( level_z(nzb:nzt+1), local_2d(nxlg:nxrg,nysg:nyng) )
[1]	295
[1308]	296	IF ( netcdf_data_format > 4 ) THEN
	297	ns = 1
[1960]	298	DO WHILE ( section(ns,s_ind) /= -9999 .AND. ns <= 100 )
[1308]	299	ns = ns + 1
	300	ENDDO
	301	ns = ns - 1
	302	ALLOCATE( local_2d_sections(nxlg:nxrg,nysg:nyng,1:ns) )
[1353]	303	local_2d_sections = 0.0_wp
[1308]	304	ENDIF
	305
[493]	306	!
[1031]	307	!-- Parallel netCDF4/HDF5 output is done on all PEs, all other on PE0 only
[1327]	308	IF ( myid == 0 .OR. netcdf_data_format > 4 ) THEN
[493]	309	CALL check_open( 101+av*10 )
	310	ENDIF
[1]	311	IF ( data_output_2d_on_each_pe ) THEN
	312	CALL check_open( 21 )
	313	ELSE
	314	IF ( myid == 0 ) THEN
	315	#if defined( __parallel )
[667]	316	ALLOCATE( total_2d(-nbgp:nx+nbgp,-nbgp:ny+nbgp) )
[1]	317	#endif
	318	ENDIF
	319	ENDIF
	320
	321	CASE ( 'xz' )
[1960]	322	s_ind = 2
[667]	323	ALLOCATE( local_2d(nxlg:nxrg,nzb:nzt+1) )
[1]	324
[1308]	325	IF ( netcdf_data_format > 4 ) THEN
	326	ns = 1
[1960]	327	DO WHILE ( section(ns,s_ind) /= -9999 .AND. ns <= 100 )
[1308]	328	ns = ns + 1
	329	ENDDO
	330	ns = ns - 1
	331	ALLOCATE( local_2d_sections(nxlg:nxrg,1:ns,nzb:nzt+1) )
	332	ALLOCATE( local_2d_sections_l(nxlg:nxrg,1:ns,nzb:nzt+1) )
[1353]	333	local_2d_sections = 0.0_wp; local_2d_sections_l = 0.0_wp
[1308]	334	ENDIF
	335
[493]	336	!
[1031]	337	!-- Parallel netCDF4/HDF5 output is done on all PEs, all other on PE0 only
[1327]	338	IF ( myid == 0 .OR. netcdf_data_format > 4 ) THEN
[493]	339	CALL check_open( 102+av*10 )
	340	ENDIF
[1]	341
	342	IF ( data_output_2d_on_each_pe ) THEN
	343	CALL check_open( 22 )
	344	ELSE
	345	IF ( myid == 0 ) THEN
	346	#if defined( __parallel )
[667]	347	ALLOCATE( total_2d(-nbgp:nx+nbgp,nzb:nzt+1) )
[1]	348	#endif
	349	ENDIF
	350	ENDIF
	351
	352	CASE ( 'yz' )
[1960]	353	s_ind = 3
[667]	354	ALLOCATE( local_2d(nysg:nyng,nzb:nzt+1) )
[1]	355
[1308]	356	IF ( netcdf_data_format > 4 ) THEN
	357	ns = 1
[1960]	358	DO WHILE ( section(ns,s_ind) /= -9999 .AND. ns <= 100 )
[1308]	359	ns = ns + 1
	360	ENDDO
	361	ns = ns - 1
	362	ALLOCATE( local_2d_sections(1:ns,nysg:nyng,nzb:nzt+1) )
	363	ALLOCATE( local_2d_sections_l(1:ns,nysg:nyng,nzb:nzt+1) )
[1353]	364	local_2d_sections = 0.0_wp; local_2d_sections_l = 0.0_wp
[1308]	365	ENDIF
	366
[493]	367	!
[1031]	368	!-- Parallel netCDF4/HDF5 output is done on all PEs, all other on PE0 only
[1327]	369	IF ( myid == 0 .OR. netcdf_data_format > 4 ) THEN
[493]	370	CALL check_open( 103+av*10 )
	371	ENDIF
[1]	372
	373	IF ( data_output_2d_on_each_pe ) THEN
	374	CALL check_open( 23 )
	375	ELSE
	376	IF ( myid == 0 ) THEN
	377	#if defined( __parallel )
[667]	378	ALLOCATE( total_2d(-nbgp:ny+nbgp,nzb:nzt+1) )
[1]	379	#endif
	380	ENDIF
	381	ENDIF
	382
	383	CASE DEFAULT
[254]	384	message_string = 'unknown cross-section: ' // TRIM( mode )
	385	CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
[1]	386
	387	END SELECT
	388
	389	!
[1745]	390	!-- For parallel netcdf output the time axis must be limited. Return, if this
	391	!-- limit is exceeded. This could be the case, if the simulated time exceeds
	392	!-- the given end time by the length of the given output interval.
	393	IF ( netcdf_data_format > 4 ) THEN
	394	IF ( mode == 'xy' .AND. do2d_xy_time_count(av) + 1 > &
	395	ntdim_2d_xy(av) ) THEN
	396	WRITE ( message_string, * ) 'Output of xy cross-sections is not ', &
	397	'given at t=', simulated_time, '&because the', &
	398	' maximum number of output time levels is exceeded.'
	399	CALL message( 'data_output_2d', 'PA0384', 0, 1, 0, 6, 0 )
	400	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
	401	RETURN
	402	ENDIF
	403	IF ( mode == 'xz' .AND. do2d_xz_time_count(av) + 1 > &
	404	ntdim_2d_xz(av) ) THEN
	405	WRITE ( message_string, * ) 'Output of xz cross-sections is not ', &
	406	'given at t=', simulated_time, '&because the', &
	407	' maximum number of output time levels is exceeded.'
	408	CALL message( 'data_output_2d', 'PA0385', 0, 1, 0, 6, 0 )
	409	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
	410	RETURN
	411	ENDIF
	412	IF ( mode == 'yz' .AND. do2d_yz_time_count(av) + 1 > &
	413	ntdim_2d_yz(av) ) THEN
	414	WRITE ( message_string, * ) 'Output of yz cross-sections is not ', &
	415	'given at t=', simulated_time, '&because the', &
	416	' maximum number of output time levels is exceeded.'
	417	CALL message( 'data_output_2d', 'PA0386', 0, 1, 0, 6, 0 )
	418	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
	419	RETURN
	420	ENDIF
	421	ENDIF
	422
	423	!
[1]	424	!-- Allocate a temporary array for resorting (kji -> ijk).
[667]	425	ALLOCATE( local_pf(nxlg:nxrg,nysg:nyng,nzb:nzt+1) )
[2232]	426	local_pf = 0.0
	427	!
	428	!-- Allocate temporary array used for exchanging ghoist points of surface-data
	429	ALLOCATE( tmp_2d(nysg:nyng,nxlg:nxrg) )
	430	tmp_2d = 0.0
[1]	431
	432	!
	433	!-- Loop of all variables to be written.
	434	!-- Output dimensions chosen
	435	if = 1
	436	l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
	437	do2d_mode = do2d(av,if)(l-1:l)
	438
	439	DO WHILE ( do2d(av,if)(1:1) /= ' ' )
	440
	441	IF ( do2d_mode == mode ) THEN
[1980]	442	!
	443	!-- Set flag to steer output of radiation, land-surface, or user-defined
	444	!-- quantities
	445	found = .FALSE.
[1551]	446
	447	nzb_do = nzb
	448	nzt_do = nzt+1
[1]	449	!
	450	!-- Store the array chosen on the temporary array.
	451	resorted = .FALSE.
	452	SELECT CASE ( TRIM( do2d(av,if) ) )
	453
	454	CASE ( 'e_xy', 'e_xz', 'e_yz' )
	455	IF ( av == 0 ) THEN
	456	to_be_resorted => e
	457	ELSE
	458	to_be_resorted => e_av
	459	ENDIF
	460	IF ( mode == 'xy' ) level_z = zu
	461
[771]	462	CASE ( 'lpt_xy', 'lpt_xz', 'lpt_yz' )
	463	IF ( av == 0 ) THEN
	464	to_be_resorted => pt
	465	ELSE
	466	to_be_resorted => lpt_av
	467	ENDIF
	468	IF ( mode == 'xy' ) level_z = zu
	469
[1]	470	CASE ( 'lwp*_xy' ) ! 2d-array
	471	IF ( av == 0 ) THEN
[667]	472	DO i = nxlg, nxrg
	473	DO j = nysg, nyng
[1320]	474	local_pf(i,j,nzb+1) = SUM( ql(nzb:nzt,j,i) * &
[1]	475	dzw(1:nzt+1) )
	476	ENDDO
	477	ENDDO
	478	ELSE
[667]	479	DO i = nxlg, nxrg
	480	DO j = nysg, nyng
[1]	481	local_pf(i,j,nzb+1) = lwp_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
[1053]	489	CASE ( 'nr_xy', 'nr_xz', 'nr_yz' )
	490	IF ( av == 0 ) THEN
	491	to_be_resorted => nr
	492	ELSE
	493	to_be_resorted => nr_av
	494	ENDIF
	495	IF ( mode == 'xy' ) level_z = zu
	496
[1691]	497	CASE ( 'ol*_xy' ) ! 2d-array
	498	IF ( av == 0 ) THEN
[2232]	499	DO m = 1, surf_def_h(0)%ns
	500	i = surf_def_h(0)%i(m)
	501	j = surf_def_h(0)%j(m)
	502	tmp_2d(j,i) = surf_def_h(0)%ol(m)
	503	ENDDO
	504	DO m = 1, surf_lsm_h%ns
	505	i = surf_lsm_h%i(m)
	506	j = surf_lsm_h%j(m)
	507	tmp_2d(j,i) = surf_lsm_h%ol(m)
	508	ENDDO
	509	DO m = 1, surf_usm_h%ns
	510	i = surf_usm_h%i(m)
	511	j = surf_usm_h%j(m)
	512	tmp_2d(j,i) = surf_usm_h%ol(m)
	513	ENDDO
	514
	515	CALL exchange_horiz_2d( tmp_2d, nbgp )
	516
[1691]	517	DO i = nxlg, nxrg
	518	DO j = nysg, nyng
[2232]	519	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[1691]	520	ENDDO
	521	ENDDO
[2232]	522
[1691]	523	ELSE
	524	DO i = nxlg, nxrg
	525	DO j = nysg, nyng
	526	local_pf(i,j,nzb+1) = ol_av(j,i)
	527	ENDDO
	528	ENDDO
	529	ENDIF
	530	resorted = .TRUE.
	531	two_d = .TRUE.
	532	level_z(nzb+1) = zu(nzb+1)
	533
[1]	534	CASE ( 'p_xy', 'p_xz', 'p_yz' )
	535	IF ( av == 0 ) THEN
[729]	536	IF ( psolver /= 'sor' ) CALL exchange_horiz( p, nbgp )
[1]	537	to_be_resorted => p
	538	ELSE
[729]	539	IF ( psolver /= 'sor' ) CALL exchange_horiz( p_av, nbgp )
[1]	540	to_be_resorted => p_av
	541	ENDIF
	542	IF ( mode == 'xy' ) level_z = zu
	543
	544	CASE ( 'pc_xy', 'pc_xz', 'pc_yz' ) ! particle concentration
	545	IF ( av == 0 ) THEN
[215]	546	IF ( simulated_time >= particle_advection_start ) THEN
	547	tend = prt_count
[667]	548	CALL exchange_horiz( tend, nbgp )
[215]	549	ELSE
[1353]	550	tend = 0.0_wp
[215]	551	ENDIF
[667]	552	DO i = nxlg, nxrg
	553	DO j = nysg, nyng
[1]	554	DO k = nzb, nzt+1
	555	local_pf(i,j,k) = tend(k,j,i)
	556	ENDDO
	557	ENDDO
	558	ENDDO
	559	resorted = .TRUE.
	560	ELSE
[667]	561	CALL exchange_horiz( pc_av, nbgp )
[1]	562	to_be_resorted => pc_av
	563	ENDIF
	564
[1359]	565	CASE ( 'pr_xy', 'pr_xz', 'pr_yz' ) ! mean particle radius (effective radius)
[1]	566	IF ( av == 0 ) THEN
[215]	567	IF ( simulated_time >= particle_advection_start ) THEN
	568	DO i = nxl, nxr
	569	DO j = nys, nyn
	570	DO k = nzb, nzt+1
[1359]	571	number_of_particles = prt_count(k,j,i)
	572	IF (number_of_particles <= 0) CYCLE
	573	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
	574	s_r2 = 0.0_wp
[1353]	575	s_r3 = 0.0_wp
[1359]	576	DO n = 1, number_of_particles
	577	IF ( particles(n)%particle_mask ) THEN
	578	s_r2 = s_r2 + particles(n)%radius*2 &
	579	particles(n)%weight_factor
	580	s_r3 = s_r3 + particles(n)%radius*3 &
	581	particles(n)%weight_factor
	582	ENDIF
[215]	583	ENDDO
[1359]	584	IF ( s_r2 > 0.0_wp ) THEN
	585	mean_r = s_r3 / s_r2
[215]	586	ELSE
[1353]	587	mean_r = 0.0_wp
[215]	588	ENDIF
	589	tend(k,j,i) = mean_r
[1]	590	ENDDO
	591	ENDDO
	592	ENDDO
[667]	593	CALL exchange_horiz( tend, nbgp )
[215]	594	ELSE
[1353]	595	tend = 0.0_wp
[1359]	596	ENDIF
[667]	597	DO i = nxlg, nxrg
	598	DO j = nysg, nyng
[1]	599	DO k = nzb, nzt+1
	600	local_pf(i,j,k) = tend(k,j,i)
	601	ENDDO
	602	ENDDO
	603	ENDDO
	604	resorted = .TRUE.
	605	ELSE
[667]	606	CALL exchange_horiz( pr_av, nbgp )
[1]	607	to_be_resorted => pr_av
	608	ENDIF
	609
[72]	610	CASE ( 'pra*_xy' ) ! 2d-array / integral quantity => no av
	611	CALL exchange_horiz_2d( precipitation_amount )
[667]	612	DO i = nxlg, nxrg
	613	DO j = nysg, nyng
[72]	614	local_pf(i,j,nzb+1) = precipitation_amount(j,i)
	615	ENDDO
	616	ENDDO
[1353]	617	precipitation_amount = 0.0_wp ! reset for next integ. interval
[72]	618	resorted = .TRUE.
	619	two_d = .TRUE.
	620	level_z(nzb+1) = zu(nzb+1)
	621
	622	CASE ( 'prr*_xy' ) ! 2d-array
[1822]	623	IF ( av == 0 ) THEN
	624	CALL exchange_horiz_2d( prr(nzb+1,:,:) )
	625	DO i = nxlg, nxrg
	626	DO j = nysg, nyng
	627	local_pf(i,j,nzb+1) = prr(nzb+1,j,i) * hyrho(nzb+1)
[1053]	628	ENDDO
[1822]	629	ENDDO
[1053]	630	ELSE
[1822]	631	CALL exchange_horiz_2d( prr_av(nzb+1,:,:) )
	632	DO i = nxlg, nxrg
	633	DO j = nysg, nyng
	634	local_pf(i,j,nzb+1) = prr_av(nzb+1,j,i) * hyrho(nzb+1)
[1053]	635	ENDDO
[1822]	636	ENDDO
[1053]	637	ENDIF
	638	resorted = .TRUE.
	639	two_d = .TRUE.
	640	level_z(nzb+1) = zu(nzb+1)
	641
	642	CASE ( 'prr_xy', 'prr_xz', 'prr_yz' )
[72]	643	IF ( av == 0 ) THEN
[1053]	644	CALL exchange_horiz( prr, nbgp )
[667]	645	DO i = nxlg, nxrg
	646	DO j = nysg, nyng
[1053]	647	DO k = nzb, nzt+1
[1822]	648	local_pf(i,j,k) = prr(k,j,i) * hyrho(nzb+1)
[1053]	649	ENDDO
[72]	650	ENDDO
	651	ENDDO
	652	ELSE
[1053]	653	CALL exchange_horiz( prr_av, nbgp )
[667]	654	DO i = nxlg, nxrg
	655	DO j = nysg, nyng
[1053]	656	DO k = nzb, nzt+1
[1822]	657	local_pf(i,j,k) = prr_av(k,j,i) * hyrho(nzb+1)
[1053]	658	ENDDO
[72]	659	ENDDO
	660	ENDDO
	661	ENDIF
	662	resorted = .TRUE.
[1053]	663	IF ( mode == 'xy' ) level_z = zu
[72]	664
[1]	665	CASE ( 'pt_xy', 'pt_xz', 'pt_yz' )
	666	IF ( av == 0 ) THEN
	667	IF ( .NOT. cloud_physics ) THEN
	668	to_be_resorted => pt
	669	ELSE
[667]	670	DO i = nxlg, nxrg
	671	DO j = nysg, nyng
[1]	672	DO k = nzb, nzt+1
[1320]	673	local_pf(i,j,k) = pt(k,j,i) + l_d_cp * &
	674	pt_d_t(k) * &
[1]	675	ql(k,j,i)
	676	ENDDO
	677	ENDDO
	678	ENDDO
	679	resorted = .TRUE.
	680	ENDIF
	681	ELSE
	682	to_be_resorted => pt_av
	683	ENDIF
	684	IF ( mode == 'xy' ) level_z = zu
	685
	686	CASE ( 'q_xy', 'q_xz', 'q_yz' )
	687	IF ( av == 0 ) THEN
	688	to_be_resorted => q
	689	ELSE
	690	to_be_resorted => q_av
	691	ENDIF
	692	IF ( mode == 'xy' ) level_z = zu
	693
[1053]	694	CASE ( 'qc_xy', 'qc_xz', 'qc_yz' )
[1]	695	IF ( av == 0 ) THEN
[1115]	696	to_be_resorted => qc
[1]	697	ELSE
[1115]	698	to_be_resorted => qc_av
[1]	699	ENDIF
	700	IF ( mode == 'xy' ) level_z = zu
	701
[1053]	702	CASE ( 'ql_xy', 'ql_xz', 'ql_yz' )
	703	IF ( av == 0 ) THEN
[1115]	704	to_be_resorted => ql
[1053]	705	ELSE
[1115]	706	to_be_resorted => ql_av
[1053]	707	ENDIF
	708	IF ( mode == 'xy' ) level_z = zu
	709
[1]	710	CASE ( 'ql_c_xy', 'ql_c_xz', 'ql_c_yz' )
	711	IF ( av == 0 ) THEN
	712	to_be_resorted => ql_c
	713	ELSE
	714	to_be_resorted => ql_c_av
	715	ENDIF
	716	IF ( mode == 'xy' ) level_z = zu
	717
	718	CASE ( 'ql_v_xy', 'ql_v_xz', 'ql_v_yz' )
	719	IF ( av == 0 ) THEN
	720	to_be_resorted => ql_v
	721	ELSE
	722	to_be_resorted => ql_v_av
	723	ENDIF
	724	IF ( mode == 'xy' ) level_z = zu
	725
	726	CASE ( 'ql_vp_xy', 'ql_vp_xz', 'ql_vp_yz' )
	727	IF ( av == 0 ) THEN
[1007]	728	IF ( simulated_time >= particle_advection_start ) THEN
	729	DO i = nxl, nxr
	730	DO j = nys, nyn
	731	DO k = nzb, nzt+1
[1359]	732	number_of_particles = prt_count(k,j,i)
	733	IF (number_of_particles <= 0) CYCLE
	734	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
	735	DO n = 1, number_of_particles
	736	IF ( particles(n)%particle_mask ) THEN
	737	tend(k,j,i) = tend(k,j,i) + &
	738	particles(n)%weight_factor / &
	739	prt_count(k,j,i)
	740	ENDIF
[1007]	741	ENDDO
	742	ENDDO
	743	ENDDO
	744	ENDDO
	745	CALL exchange_horiz( tend, nbgp )
	746	ELSE
[1353]	747	tend = 0.0_wp
[1359]	748	ENDIF
[1007]	749	DO i = nxlg, nxrg
	750	DO j = nysg, nyng
	751	DO k = nzb, nzt+1
	752	local_pf(i,j,k) = tend(k,j,i)
	753	ENDDO
	754	ENDDO
	755	ENDDO
	756	resorted = .TRUE.
	757	ELSE
	758	CALL exchange_horiz( ql_vp_av, nbgp )
[1]	759	to_be_resorted => ql_vp
	760	ENDIF
	761	IF ( mode == 'xy' ) level_z = zu
	762
[1053]	763	CASE ( 'qr_xy', 'qr_xz', 'qr_yz' )
	764	IF ( av == 0 ) THEN
	765	to_be_resorted => qr
	766	ELSE
	767	to_be_resorted => qr_av
	768	ENDIF
	769	IF ( mode == 'xy' ) level_z = zu
	770
[354]	771	CASE ( 'qsws*_xy' ) ! 2d-array
	772	IF ( av == 0 ) THEN
[2232]	773	DO m = 1, surf_def_h(0)%ns
	774	i = surf_def_h(0)%i(m)
	775	j = surf_def_h(0)%j(m)
	776	tmp_2d(j,i) = surf_def_h(0)%qsws(m)
	777	ENDDO
	778	DO m = 1, surf_lsm_h%ns
	779	i = surf_lsm_h%i(m)
	780	j = surf_lsm_h%j(m)
	781	tmp_2d(j,i) = surf_lsm_h%qsws(m)
	782	ENDDO
	783	DO m = 1, surf_usm_h%ns
	784	i = surf_usm_h%i(m)
	785	j = surf_usm_h%j(m)
	786	tmp_2d(j,i) = surf_usm_h%qsws(m)
	787	ENDDO
	788
	789	CALL exchange_horiz_2d( tmp_2d, nbgp )
	790
[667]	791	DO i = nxlg, nxrg
	792	DO j = nysg, nyng
[2232]	793	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[354]	794	ENDDO
	795	ENDDO
	796	ELSE
[667]	797	DO i = nxlg, nxrg
	798	DO j = nysg, nyng
[354]	799	local_pf(i,j,nzb+1) = qsws_av(j,i)
	800	ENDDO
	801	ENDDO
	802	ENDIF
	803	resorted = .TRUE.
	804	two_d = .TRUE.
	805	level_z(nzb+1) = zu(nzb+1)
	806
[1]	807	CASE ( 'qv_xy', 'qv_xz', 'qv_yz' )
	808	IF ( av == 0 ) THEN
[667]	809	DO i = nxlg, nxrg
	810	DO j = nysg, nyng
[1]	811	DO k = nzb, nzt+1
	812	local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
	813	ENDDO
	814	ENDDO
	815	ENDDO
	816	resorted = .TRUE.
	817	ELSE
	818	to_be_resorted => qv_av
	819	ENDIF
	820	IF ( mode == 'xy' ) level_z = zu
	821
[1551]	822
[1585]	823
[2190]	824	CASE ( 'rho_ocean_xy', 'rho_ocean_xz', 'rho_ocean_yz' )
[96]	825	IF ( av == 0 ) THEN
[2031]	826	to_be_resorted => rho_ocean
[96]	827	ELSE
[2031]	828	to_be_resorted => rho_ocean_av
[96]	829	ENDIF
	830
[1]	831	CASE ( 's_xy', 's_xz', 's_yz' )
	832	IF ( av == 0 ) THEN
[1960]	833	to_be_resorted => s
[1]	834	ELSE
[355]	835	to_be_resorted => s_av
[1]	836	ENDIF
	837
[96]	838	CASE ( 'sa_xy', 'sa_xz', 'sa_yz' )
	839	IF ( av == 0 ) THEN
	840	to_be_resorted => sa
	841	ELSE
	842	to_be_resorted => sa_av
	843	ENDIF
	844
[354]	845	CASE ( 'shf*_xy' ) ! 2d-array
	846	IF ( av == 0 ) THEN
[2232]	847	DO m = 1, surf_def_h(0)%ns
	848	i = surf_def_h(0)%i(m)
	849	j = surf_def_h(0)%j(m)
	850	tmp_2d(j,i) = surf_def_h(0)%shf(m)
	851	ENDDO
	852	DO m = 1, surf_lsm_h%ns
	853	i = surf_lsm_h%i(m)
	854	j = surf_lsm_h%j(m)
	855	tmp_2d(j,i) = surf_lsm_h%shf(m)
	856	ENDDO
	857	DO m = 1, surf_usm_h%ns
	858	i = surf_usm_h%i(m)
	859	j = surf_usm_h%j(m)
	860	tmp_2d(j,i) = surf_usm_h%shf(m)
	861	ENDDO
	862
	863	CALL exchange_horiz_2d( tmp_2d, nbgp )
	864
[667]	865	DO i = nxlg, nxrg
	866	DO j = nysg, nyng
[2232]	867	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[354]	868	ENDDO
	869	ENDDO
	870	ELSE
[667]	871	DO i = nxlg, nxrg
	872	DO j = nysg, nyng
[354]	873	local_pf(i,j,nzb+1) = shf_av(j,i)
	874	ENDDO
	875	ENDDO
	876	ENDIF
	877	resorted = .TRUE.
	878	two_d = .TRUE.
	879	level_z(nzb+1) = zu(nzb+1)
[1960]	880
	881	CASE ( 'ssws*_xy' ) ! 2d-array
	882	IF ( av == 0 ) THEN
[2232]	883	DO m = 1, surf_def_h(0)%ns
	884	i = surf_def_h(0)%i(m)
	885	j = surf_def_h(0)%j(m)
	886	tmp_2d(j,i) = surf_def_h(0)%ssws(m)
	887	ENDDO
	888	DO m = 1, surf_lsm_h%ns
	889	i = surf_lsm_h%i(m)
	890	j = surf_lsm_h%j(m)
	891	tmp_2d(j,i) = surf_lsm_h%ssws(m)
	892	ENDDO
	893	DO m = 1, surf_usm_h%ns
	894	i = surf_usm_h%i(m)
	895	j = surf_usm_h%j(m)
	896	tmp_2d(j,i) = surf_usm_h%ssws(m)
	897	ENDDO
	898
	899	CALL exchange_horiz_2d( tmp_2d, nbgp )
	900
[1960]	901	DO i = nxlg, nxrg
	902	DO j = nysg, nyng
[2232]	903	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[1960]	904	ENDDO
	905	ENDDO
	906	ELSE
	907	DO i = nxlg, nxrg
	908	DO j = nysg, nyng
	909	local_pf(i,j,nzb+1) = ssws_av(j,i)
	910	ENDDO
	911	ENDDO
	912	ENDIF
	913	resorted = .TRUE.
	914	two_d = .TRUE.
	915	level_z(nzb+1) = zu(nzb+1)
[1551]	916
[1]	917	CASE ( 't*_xy' ) ! 2d-array
	918	IF ( av == 0 ) THEN
[2232]	919	DO m = 1, surf_def_h(0)%ns
	920	i = surf_def_h(0)%i(m)
	921	j = surf_def_h(0)%j(m)
	922	tmp_2d(j,i) = surf_def_h(0)%ts(m)
	923	ENDDO
	924	DO m = 1, surf_lsm_h%ns
	925	i = surf_lsm_h%i(m)
	926	j = surf_lsm_h%j(m)
	927	tmp_2d(j,i) = surf_lsm_h%ts(m)
	928	ENDDO
	929	DO m = 1, surf_usm_h%ns
	930	i = surf_usm_h%i(m)
	931	j = surf_usm_h%j(m)
	932	tmp_2d(j,i) = surf_usm_h%ts(m)
	933	ENDDO
	934
	935	CALL exchange_horiz_2d( tmp_2d, nbgp )
	936
[667]	937	DO i = nxlg, nxrg
	938	DO j = nysg, nyng
[2232]	939	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[1]	940	ENDDO
	941	ENDDO
[2232]	942
[1]	943	ELSE
[667]	944	DO i = nxlg, nxrg
	945	DO j = nysg, nyng
[1]	946	local_pf(i,j,nzb+1) = ts_av(j,i)
	947	ENDDO
	948	ENDDO
	949	ENDIF
	950	resorted = .TRUE.
	951	two_d = .TRUE.
	952	level_z(nzb+1) = zu(nzb+1)
	953
	954	CASE ( 'u_xy', 'u_xz', 'u_yz' )
	955	IF ( av == 0 ) THEN
	956	to_be_resorted => u
	957	ELSE
	958	to_be_resorted => u_av
	959	ENDIF
	960	IF ( mode == 'xy' ) level_z = zu
	961	!
	962	!-- Substitute the values generated by "mirror" boundary condition
	963	!-- at the bottom boundary by the real surface values.
	964	IF ( do2d(av,if) == 'u_xz' .OR. do2d(av,if) == 'u_yz' ) THEN
[1353]	965	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0_wp
[1]	966	ENDIF
	967
	968	CASE ( 'u*_xy' ) ! 2d-array
	969	IF ( av == 0 ) THEN
[2232]	970	DO m = 1, surf_def_h(0)%ns
	971	i = surf_def_h(0)%i(m)
	972	j = surf_def_h(0)%j(m)
	973	tmp_2d(j,i) = surf_def_h(0)%us(m)
	974	ENDDO
	975	DO m = 1, surf_lsm_h%ns
	976	i = surf_lsm_h%i(m)
	977	j = surf_lsm_h%j(m)
	978	tmp_2d(j,i) = surf_lsm_h%us(m)
	979	ENDDO
	980	DO m = 1, surf_usm_h%ns
	981	i = surf_usm_h%i(m)
	982	j = surf_usm_h%j(m)
	983	tmp_2d(j,i) = surf_usm_h%us(m)
	984	ENDDO
	985
	986	CALL exchange_horiz_2d( tmp_2d, nbgp )
	987
[667]	988	DO i = nxlg, nxrg
	989	DO j = nysg, nyng
[2232]	990	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[1]	991	ENDDO
	992	ENDDO
	993	ELSE
[667]	994	DO i = nxlg, nxrg
	995	DO j = nysg, nyng
[1]	996	local_pf(i,j,nzb+1) = us_av(j,i)
	997	ENDDO
	998	ENDDO
	999	ENDIF
	1000	resorted = .TRUE.
	1001	two_d = .TRUE.
	1002	level_z(nzb+1) = zu(nzb+1)
	1003
	1004	CASE ( 'v_xy', 'v_xz', 'v_yz' )
	1005	IF ( av == 0 ) THEN
	1006	to_be_resorted => v
	1007	ELSE
	1008	to_be_resorted => v_av
	1009	ENDIF
	1010	IF ( mode == 'xy' ) level_z = zu
	1011	!
	1012	!-- Substitute the values generated by "mirror" boundary condition
	1013	!-- at the bottom boundary by the real surface values.
	1014	IF ( do2d(av,if) == 'v_xz' .OR. do2d(av,if) == 'v_yz' ) THEN
[1353]	1015	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0_wp
[1]	1016	ENDIF
	1017
	1018	CASE ( 'vpt_xy', 'vpt_xz', 'vpt_yz' )
	1019	IF ( av == 0 ) THEN
	1020	to_be_resorted => vpt
	1021	ELSE
	1022	to_be_resorted => vpt_av
	1023	ENDIF
	1024	IF ( mode == 'xy' ) level_z = zu
	1025
	1026	CASE ( 'w_xy', 'w_xz', 'w_yz' )
	1027	IF ( av == 0 ) THEN
	1028	to_be_resorted => w
	1029	ELSE
	1030	to_be_resorted => w_av
	1031	ENDIF
	1032	IF ( mode == 'xy' ) level_z = zw
	1033
[72]	1034	CASE ( 'z0*_xy' ) ! 2d-array
	1035	IF ( av == 0 ) THEN
[2232]	1036	DO m = 1, surf_def_h(0)%ns
	1037	i = surf_def_h(0)%i(m)
	1038	j = surf_def_h(0)%j(m)
	1039	tmp_2d(j,i) = surf_def_h(0)%z0(m)
	1040	ENDDO
	1041	DO m = 1, surf_lsm_h%ns
	1042	i = surf_lsm_h%i(m)
	1043	j = surf_lsm_h%j(m)
	1044	tmp_2d(j,i) = surf_lsm_h%z0(m)
	1045	ENDDO
	1046	DO m = 1, surf_usm_h%ns
	1047	i = surf_usm_h%i(m)
	1048	j = surf_usm_h%j(m)
	1049	tmp_2d(j,i) = surf_usm_h%z0(m)
	1050	ENDDO
	1051
	1052	CALL exchange_horiz_2d( tmp_2d, nbgp )
	1053
[667]	1054	DO i = nxlg, nxrg
	1055	DO j = nysg, nyng
[2232]	1056	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[72]	1057	ENDDO
	1058	ENDDO
[2232]	1059
[72]	1060	ELSE
[667]	1061	DO i = nxlg, nxrg
	1062	DO j = nysg, nyng
[72]	1063	local_pf(i,j,nzb+1) = z0_av(j,i)
	1064	ENDDO
	1065	ENDDO
	1066	ENDIF
	1067	resorted = .TRUE.
	1068	two_d = .TRUE.
	1069	level_z(nzb+1) = zu(nzb+1)
	1070
[978]	1071	CASE ( 'z0h*_xy' ) ! 2d-array
	1072	IF ( av == 0 ) THEN
[2232]	1073	DO m = 1, surf_def_h(0)%ns
	1074	i = surf_def_h(0)%i(m)
	1075	j = surf_def_h(0)%j(m)
	1076	tmp_2d(j,i) = surf_def_h(0)%z0h(m)
	1077	ENDDO
	1078	DO m = 1, surf_lsm_h%ns
	1079	i = surf_lsm_h%i(m)
	1080	j = surf_lsm_h%j(m)
	1081	tmp_2d(j,i) = surf_lsm_h%z0h(m)
	1082	ENDDO
	1083	DO m = 1, surf_usm_h%ns
	1084	i = surf_usm_h%i(m)
	1085	j = surf_usm_h%j(m)
	1086	tmp_2d(j,i) = surf_usm_h%z0h(m)
	1087	ENDDO
	1088
	1089	CALL exchange_horiz_2d( tmp_2d, nbgp )
	1090
[978]	1091	DO i = nxlg, nxrg
	1092	DO j = nysg, nyng
[2232]	1093	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[978]	1094	ENDDO
	1095	ENDDO
	1096	ELSE
	1097	DO i = nxlg, nxrg
	1098	DO j = nysg, nyng
	1099	local_pf(i,j,nzb+1) = z0h_av(j,i)
	1100	ENDDO
	1101	ENDDO
	1102	ENDIF
	1103	resorted = .TRUE.
	1104	two_d = .TRUE.
	1105	level_z(nzb+1) = zu(nzb+1)
	1106
[1788]	1107	CASE ( 'z0q*_xy' ) ! 2d-array
	1108	IF ( av == 0 ) THEN
[2232]	1109	DO m = 1, surf_def_h(0)%ns
	1110	i = surf_def_h(0)%i(m)
	1111	j = surf_def_h(0)%j(m)
	1112	tmp_2d(j,i) = surf_def_h(0)%z0q(m)
	1113	ENDDO
	1114	DO m = 1, surf_lsm_h%ns
	1115	i = surf_lsm_h%i(m)
	1116	j = surf_lsm_h%j(m)
	1117	tmp_2d(j,i) = surf_lsm_h%z0q(m)
	1118	ENDDO
	1119	DO m = 1, surf_usm_h%ns
	1120	i = surf_usm_h%i(m)
	1121	j = surf_usm_h%j(m)
	1122	tmp_2d(j,i) = surf_usm_h%z0q(m)
	1123	ENDDO
	1124
	1125	CALL exchange_horiz_2d( tmp_2d, nbgp )
	1126
[1788]	1127	DO i = nxlg, nxrg
	1128	DO j = nysg, nyng
[2232]	1129	local_pf(i,j,nzb+1) = tmp_2d(j,i)
[1788]	1130	ENDDO
	1131	ENDDO
	1132	ELSE
	1133	DO i = nxlg, nxrg
	1134	DO j = nysg, nyng
	1135	local_pf(i,j,nzb+1) = z0q_av(j,i)
	1136	ENDDO
	1137	ENDDO
	1138	ENDIF
	1139	resorted = .TRUE.
	1140	two_d = .TRUE.
	1141	level_z(nzb+1) = zu(nzb+1)
	1142
[1]	1143	CASE DEFAULT
[1972]	1144
[1]	1145	!
[1972]	1146	!-- Land surface model quantity
	1147	IF ( land_surface ) THEN
	1148	CALL lsm_data_output_2d( av, do2d(av,if), found, grid, mode,&
	1149	local_pf, two_d, nzb_do, nzt_do )
	1150	ENDIF
	1151
	1152	!
[1976]	1153	!-- Radiation quantity
	1154	IF ( .NOT. found .AND. radiation ) THEN
	1155	CALL radiation_data_output_2d( av, do2d(av,if), found, grid,&
	1156	mode, local_pf, two_d )
	1157	ENDIF
	1158
	1159	!
[1]	1160	!-- User defined quantity
[1972]	1161	IF ( .NOT. found ) THEN
	1162	CALL user_data_output_2d( av, do2d(av,if), found, grid, &
	1163	local_pf, two_d, nzb_do, nzt_do )
	1164	ENDIF
	1165
[1]	1166	resorted = .TRUE.
	1167
	1168	IF ( grid == 'zu' ) THEN
	1169	IF ( mode == 'xy' ) level_z = zu
	1170	ELSEIF ( grid == 'zw' ) THEN
	1171	IF ( mode == 'xy' ) level_z = zw
[343]	1172	ELSEIF ( grid == 'zu1' ) THEN
	1173	IF ( mode == 'xy' ) level_z(nzb+1) = zu(nzb+1)
[1551]	1174	ELSEIF ( grid == 'zs' ) THEN
	1175	IF ( mode == 'xy' ) level_z = zs
[1]	1176	ENDIF
	1177
	1178	IF ( .NOT. found ) THEN
[1320]	1179	message_string = 'no output provided for: ' // &
[274]	1180	TRIM( do2d(av,if) )
[254]	1181	CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
[1]	1182	ENDIF
	1183
	1184	END SELECT
	1185
	1186	!
	1187	!-- Resort the array to be output, if not done above
	1188	IF ( .NOT. resorted ) THEN
[667]	1189	DO i = nxlg, nxrg
	1190	DO j = nysg, nyng
[1551]	1191	DO k = nzb_do, nzt_do
[1]	1192	local_pf(i,j,k) = to_be_resorted(k,j,i)
	1193	ENDDO
	1194	ENDDO
	1195	ENDDO
	1196	ENDIF
	1197
	1198	!
	1199	!-- Output of the individual cross-sections, depending on the cross-
	1200	!-- section mode chosen.
	1201	is = 1
[1960]	1202	loop1: DO WHILE ( section(is,s_ind) /= -9999 .OR. two_d )
[1]	1203
	1204	SELECT CASE ( mode )
	1205
	1206	CASE ( 'xy' )
	1207	!
	1208	!-- Determine the cross section index
	1209	IF ( two_d ) THEN
	1210	layer_xy = nzb+1
	1211	ELSE
[1960]	1212	layer_xy = section(is,s_ind)
[1]	1213	ENDIF
	1214
	1215	!
[1551]	1216	!-- Exit the loop for layers beyond the data output domain
	1217	!-- (used for soil model)
[1691]	1218	IF ( layer_xy > nzt_do ) THEN
[1551]	1219	EXIT loop1
	1220	ENDIF
	1221
	1222	!
[1308]	1223	!-- Update the netCDF xy cross section time axis.
	1224	!-- In case of parallel output, this is only done by PE0
	1225	!-- to increase the performance.
	1226	IF ( simulated_time /= do2d_xy_last_time(av) ) THEN
	1227	do2d_xy_time_count(av) = do2d_xy_time_count(av) + 1
	1228	do2d_xy_last_time(av) = simulated_time
	1229	IF ( myid == 0 ) THEN
[1327]	1230	IF ( .NOT. data_output_2d_on_each_pe &
	1231	.OR. netcdf_data_format > 4 ) &
[493]	1232	THEN
[1]	1233	#if defined( __netcdf )
	1234	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1235	id_var_time_xy(av), &
[291]	1236	(/ time_since_reference_point /), &
[1]	1237	start = (/ do2d_xy_time_count(av) /), &
	1238	count = (/ 1 /) )
[1783]	1239	CALL netcdf_handle_error( 'data_output_2d', 53 )
[1]	1240	#endif
	1241	ENDIF
	1242	ENDIF
	1243	ENDIF
	1244	!
	1245	!-- If required, carry out averaging along z
[1960]	1246	IF ( section(is,s_ind) == -1 .AND. .NOT. two_d ) THEN
[1]	1247
[1353]	1248	local_2d = 0.0_wp
[1]	1249	!
	1250	!-- Carry out the averaging (all data are on the PE)
[1551]	1251	DO k = nzb_do, nzt_do
[667]	1252	DO j = nysg, nyng
	1253	DO i = nxlg, nxrg
[1]	1254	local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
	1255	ENDDO
	1256	ENDDO
	1257	ENDDO
	1258
[1551]	1259	local_2d = local_2d / ( nzt_do - nzb_do + 1.0_wp)
[1]	1260
	1261	ELSE
	1262	!
	1263	!-- Just store the respective section on the local array
	1264	local_2d = local_pf(:,:,layer_xy)
	1265
	1266	ENDIF
	1267
	1268	#if defined( __parallel )
[1327]	1269	IF ( netcdf_data_format > 4 ) THEN
[1]	1270	!
[1031]	1271	!-- Parallel output in netCDF4/HDF5 format.
[493]	1272	IF ( two_d ) THEN
	1273	iis = 1
	1274	ELSE
	1275	iis = is
	1276	ENDIF
	1277
[1]	1278	#if defined( __netcdf )
[1308]	1279	!
	1280	!-- For parallel output, all cross sections are first stored
	1281	!-- here on a local array and will be written to the output
	1282	!-- file afterwards to increase the performance.
	1283	DO i = nxlg, nxrg
	1284	DO j = nysg, nyng
	1285	local_2d_sections(i,j,iis) = local_2d(i,j)
	1286	ENDDO
	1287	ENDDO
[1]	1288	#endif
[493]	1289	ELSE
[1]	1290
[493]	1291	IF ( data_output_2d_on_each_pe ) THEN
[1]	1292	!
[493]	1293	!-- Output of partial arrays on each PE
	1294	#if defined( __netcdf )
[1327]	1295	IF ( myid == 0 ) THEN
[1320]	1296	WRITE ( 21 ) time_since_reference_point, &
[493]	1297	do2d_xy_time_count(av), av
	1298	ENDIF
	1299	#endif
[759]	1300	DO i = 0, io_blocks-1
	1301	IF ( i == io_group ) THEN
[1551]	1302	WRITE ( 21 ) nxlg, nxrg, nysg, nyng, nysg, nyng
[759]	1303	WRITE ( 21 ) local_2d
	1304	ENDIF
	1305	#if defined( __parallel )
	1306	CALL MPI_BARRIER( comm2d, ierr )
	1307	#endif
	1308	ENDDO
[559]	1309
[493]	1310	ELSE
[1]	1311	!
[493]	1312	!-- PE0 receives partial arrays from all processors and
	1313	!-- then outputs them. Here a barrier has to be set,
	1314	!-- because otherwise "-MPI- FATAL: Remote protocol queue
	1315	!-- full" may occur.
	1316	CALL MPI_BARRIER( comm2d, ierr )
	1317
[667]	1318	ngp = ( nxrg-nxlg+1 ) * ( nyng-nysg+1 )
[493]	1319	IF ( myid == 0 ) THEN
[1]	1320	!
[493]	1321	!-- Local array can be relocated directly.
[667]	1322	total_2d(nxlg:nxrg,nysg:nyng) = local_2d
[1]	1323	!
[493]	1324	!-- Receive data from all other PEs.
	1325	DO n = 1, numprocs-1
[1]	1326	!
[493]	1327	!-- Receive index limits first, then array.
	1328	!-- Index limits are received in arbitrary order from
	1329	!-- the PEs.
[1320]	1330	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
	1331	MPI_ANY_SOURCE, 0, comm2d, &
[493]	1332	status, ierr )
	1333	sender = status(MPI_SOURCE)
	1334	DEALLOCATE( local_2d )
	1335	ALLOCATE( local_2d(ind(1):ind(2),ind(3):ind(4)) )
[1320]	1336	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
	1337	MPI_REAL, sender, 1, comm2d, &
[493]	1338	status, ierr )
	1339	total_2d(ind(1):ind(2),ind(3):ind(4)) = local_2d
	1340	ENDDO
[1]	1341	!
[493]	1342	!-- Relocate the local array for the next loop increment
	1343	DEALLOCATE( local_2d )
[667]	1344	ALLOCATE( local_2d(nxlg:nxrg,nysg:nyng) )
[1]	1345
	1346	#if defined( __netcdf )
[1327]	1347	IF ( two_d ) THEN
	1348	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1349	id_var_do2d(av,if), &
	1350	total_2d(0:nx+1,0:ny+1), &
	1351	start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
	1352	count = (/ nx+2, ny+2, 1, 1 /) )
	1353	ELSE
	1354	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1355	id_var_do2d(av,if), &
	1356	total_2d(0:nx+1,0:ny+1), &
	1357	start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
	1358	count = (/ nx+2, ny+2, 1, 1 /) )
[1]	1359	ENDIF
[1783]	1360	CALL netcdf_handle_error( 'data_output_2d', 54 )
[1]	1361	#endif
	1362
[493]	1363	ELSE
[1]	1364	!
[493]	1365	!-- First send the local index limits to PE0
[667]	1366	ind(1) = nxlg; ind(2) = nxrg
	1367	ind(3) = nysg; ind(4) = nyng
[1320]	1368	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
[493]	1369	comm2d, ierr )
[1]	1370	!
[493]	1371	!-- Send data to PE0
[1320]	1372	CALL MPI_SEND( local_2d(nxlg,nysg), ngp, &
[493]	1373	MPI_REAL, 0, 1, comm2d, ierr )
	1374	ENDIF
	1375	!
	1376	!-- A barrier has to be set, because otherwise some PEs may
	1377	!-- proceed too fast so that PE0 may receive wrong data on
	1378	!-- tag 0
	1379	CALL MPI_BARRIER( comm2d, ierr )
[1]	1380	ENDIF
[493]	1381
[1]	1382	ENDIF
	1383	#else
	1384	#if defined( __netcdf )
[1327]	1385	IF ( two_d ) THEN
	1386	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1387	id_var_do2d(av,if), &
	1388	local_2d(nxl:nxr+1,nys:nyn+1), &
	1389	start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
	1390	count = (/ nx+2, ny+2, 1, 1 /) )
	1391	ELSE
	1392	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1393	id_var_do2d(av,if), &
	1394	local_2d(nxl:nxr+1,nys:nyn+1), &
	1395	start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
	1396	count = (/ nx+2, ny+2, 1, 1 /) )
[1]	1397	ENDIF
[1783]	1398	CALL netcdf_handle_error( 'data_output_2d', 447 )
[1]	1399	#endif
	1400	#endif
	1401	do2d_xy_n = do2d_xy_n + 1
	1402	!
	1403	!-- For 2D-arrays (e.g. u*) only one cross-section is available.
	1404	!-- Hence exit loop of output levels.
	1405	IF ( two_d ) THEN
[1703]	1406	IF ( netcdf_data_format < 5 ) two_d = .FALSE.
[1]	1407	EXIT loop1
	1408	ENDIF
	1409
	1410	CASE ( 'xz' )
	1411	!
[1308]	1412	!-- Update the netCDF xz cross section time axis.
	1413	!-- In case of parallel output, this is only done by PE0
	1414	!-- to increase the performance.
	1415	IF ( simulated_time /= do2d_xz_last_time(av) ) THEN
	1416	do2d_xz_time_count(av) = do2d_xz_time_count(av) + 1
	1417	do2d_xz_last_time(av) = simulated_time
	1418	IF ( myid == 0 ) THEN
[1327]	1419	IF ( .NOT. data_output_2d_on_each_pe &
	1420	.OR. netcdf_data_format > 4 ) &
[493]	1421	THEN
[1]	1422	#if defined( __netcdf )
	1423	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
	1424	id_var_time_xz(av), &
[291]	1425	(/ time_since_reference_point /), &
[1]	1426	start = (/ do2d_xz_time_count(av) /), &
	1427	count = (/ 1 /) )
[1783]	1428	CALL netcdf_handle_error( 'data_output_2d', 56 )
[1]	1429	#endif
	1430	ENDIF
	1431	ENDIF
	1432	ENDIF
[667]	1433
[1]	1434	!
	1435	!-- If required, carry out averaging along y
[1960]	1436	IF ( section(is,s_ind) == -1 ) THEN
[1]	1437
[1551]	1438	ALLOCATE( local_2d_l(nxlg:nxrg,nzb_do:nzt_do) )
[1353]	1439	local_2d_l = 0.0_wp
[1551]	1440	ngp = ( nxrg-nxlg + 1 ) * ( nzt_do-nzb_do + 1 )
[1]	1441	!
	1442	!-- First local averaging on the PE
[1551]	1443	DO k = nzb_do, nzt_do
[1]	1444	DO j = nys, nyn
[667]	1445	DO i = nxlg, nxrg
[1320]	1446	local_2d_l(i,k) = local_2d_l(i,k) + &
[1]	1447	local_pf(i,j,k)
	1448	ENDDO
	1449	ENDDO
	1450	ENDDO
	1451	#if defined( __parallel )
	1452	!
	1453	!-- Now do the averaging over all PEs along y
[622]	1454	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[1551]	1455	CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb_do), &
	1456	local_2d(nxlg,nzb_do), ngp, MPI_REAL, &
[1]	1457	MPI_SUM, comm1dy, ierr )
	1458	#else
	1459	local_2d = local_2d_l
	1460	#endif
[1353]	1461	local_2d = local_2d / ( ny + 1.0_wp )
[1]	1462
	1463	DEALLOCATE( local_2d_l )
	1464
	1465	ELSE
	1466	!
	1467	!-- Just store the respective section on the local array
	1468	!-- (but only if it is available on this PE!)
[1960]	1469	IF ( section(is,s_ind) >= nys .AND. section(is,s_ind) <= nyn ) &
[1]	1470	THEN
[1960]	1471	local_2d = local_pf(:,section(is,s_ind),nzb_do:nzt_do)
[1]	1472	ENDIF
	1473
	1474	ENDIF
	1475
	1476	#if defined( __parallel )
[1327]	1477	IF ( netcdf_data_format > 4 ) THEN
[1]	1478	!
[1031]	1479	!-- Output in netCDF4/HDF5 format.
[493]	1480	!-- Output only on those PEs where the respective cross
	1481	!-- sections reside. Cross sections averaged along y are
	1482	!-- output on the respective first PE along y (myidy=0).
[1960]	1483	IF ( ( section(is,s_ind) >= nys .AND. &
	1484	section(is,s_ind) <= nyn ) .OR. &
	1485	( section(is,s_ind) == -1 .AND. myidy == 0 ) ) THEN
[1]	1486	#if defined( __netcdf )
[493]	1487	!
[1308]	1488	!-- For parallel output, all cross sections are first
	1489	!-- stored here on a local array and will be written to the
	1490	!-- output file afterwards to increase the performance.
	1491	DO i = nxlg, nxrg
[1551]	1492	DO k = nzb_do, nzt_do
[1308]	1493	local_2d_sections_l(i,is,k) = local_2d(i,k)
	1494	ENDDO
	1495	ENDDO
[1]	1496	#endif
	1497	ENDIF
	1498
	1499	ELSE
	1500
[493]	1501	IF ( data_output_2d_on_each_pe ) THEN
[1]	1502	!
[493]	1503	!-- Output of partial arrays on each PE. If the cross
	1504	!-- section does not reside on the PE, output special
	1505	!-- index values.
	1506	#if defined( __netcdf )
[1327]	1507	IF ( myid == 0 ) THEN
[1320]	1508	WRITE ( 22 ) time_since_reference_point, &
[493]	1509	do2d_xz_time_count(av), av
	1510	ENDIF
	1511	#endif
[759]	1512	DO i = 0, io_blocks-1
	1513	IF ( i == io_group ) THEN
[1960]	1514	IF ( ( section(is,s_ind) >= nys .AND. &
	1515	section(is,s_ind) <= nyn ) .OR. &
	1516	( section(is,s_ind) == -1 .AND. &
[1320]	1517	nys-1 == -1 ) ) &
[759]	1518	THEN
[1551]	1519	WRITE (22) nxlg, nxrg, nzb_do, nzt_do, nzb, nzt+1
[759]	1520	WRITE (22) local_2d
	1521	ELSE
[1551]	1522	WRITE (22) -1, -1, -1, -1, -1, -1
[759]	1523	ENDIF
	1524	ENDIF
	1525	#if defined( __parallel )
	1526	CALL MPI_BARRIER( comm2d, ierr )
	1527	#endif
	1528	ENDDO
[493]	1529
	1530	ELSE
[1]	1531	!
[493]	1532	!-- PE0 receives partial arrays from all processors of the
	1533	!-- respective cross section and outputs them. Here a
	1534	!-- barrier has to be set, because otherwise
	1535	!-- "-MPI- FATAL: Remote protocol queue full" may occur.
	1536	CALL MPI_BARRIER( comm2d, ierr )
	1537
[1551]	1538	ngp = ( nxrg-nxlg + 1 ) * ( nzt_do-nzb_do + 1 )
[493]	1539	IF ( myid == 0 ) THEN
[1]	1540	!
[493]	1541	!-- Local array can be relocated directly.
[1960]	1542	IF ( ( section(is,s_ind) >= nys .AND. &
	1543	section(is,s_ind) <= nyn ) .OR. &
	1544	( section(is,s_ind) == -1 .AND. &
	1545	nys-1 == -1 ) ) THEN
[1551]	1546	total_2d(nxlg:nxrg,nzb_do:nzt_do) = local_2d
[493]	1547	ENDIF
[1]	1548	!
[493]	1549	!-- Receive data from all other PEs.
	1550	DO n = 1, numprocs-1
	1551	!
	1552	!-- Receive index limits first, then array.
	1553	!-- Index limits are received in arbitrary order from
	1554	!-- the PEs.
[1320]	1555	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
	1556	MPI_ANY_SOURCE, 0, comm2d, &
[1]	1557	status, ierr )
[493]	1558	!
	1559	!-- Not all PEs have data for XZ-cross-section.
	1560	IF ( ind(1) /= -9999 ) THEN
	1561	sender = status(MPI_SOURCE)
	1562	DEALLOCATE( local_2d )
[1320]	1563	ALLOCATE( local_2d(ind(1):ind(2), &
[493]	1564	ind(3):ind(4)) )
	1565	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
	1566	MPI_REAL, sender, 1, comm2d, &
	1567	status, ierr )
[1320]	1568	total_2d(ind(1):ind(2),ind(3):ind(4)) = &
[493]	1569	local_2d
	1570	ENDIF
	1571	ENDDO
	1572	!
	1573	!-- Relocate the local array for the next loop increment
	1574	DEALLOCATE( local_2d )
[1551]	1575	ALLOCATE( local_2d(nxlg:nxrg,nzb_do:nzt_do) )
[1]	1576
	1577	#if defined( __netcdf )
[1327]	1578	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
	1579	id_var_do2d(av,if), &
[1551]	1580	total_2d(0:nx+1,nzb_do:nzt_do),&
[1327]	1581	start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
[1551]	1582	count = (/ nx+2, 1, nzt_do-nzb_do+1, 1 /) )
[1783]	1583	CALL netcdf_handle_error( 'data_output_2d', 58 )
[1]	1584	#endif
	1585
[493]	1586	ELSE
[1]	1587	!
[493]	1588	!-- If the cross section resides on the PE, send the
	1589	!-- local index limits, otherwise send -9999 to PE0.
[1960]	1590	IF ( ( section(is,s_ind) >= nys .AND. &
	1591	section(is,s_ind) <= nyn ) .OR. &
	1592	( section(is,s_ind) == -1 .AND. nys-1 == -1 ) ) &
[493]	1593	THEN
[667]	1594	ind(1) = nxlg; ind(2) = nxrg
[1551]	1595	ind(3) = nzb_do; ind(4) = nzt_do
[493]	1596	ELSE
	1597	ind(1) = -9999; ind(2) = -9999
	1598	ind(3) = -9999; ind(4) = -9999
	1599	ENDIF
[1320]	1600	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
[493]	1601	comm2d, ierr )
	1602	!
	1603	!-- If applicable, send data to PE0.
	1604	IF ( ind(1) /= -9999 ) THEN
[1551]	1605	CALL MPI_SEND( local_2d(nxlg,nzb_do), ngp, &
[493]	1606	MPI_REAL, 0, 1, comm2d, ierr )
	1607	ENDIF
[1]	1608	ENDIF
	1609	!
[493]	1610	!-- A barrier has to be set, because otherwise some PEs may
	1611	!-- proceed too fast so that PE0 may receive wrong data on
	1612	!-- tag 0
	1613	CALL MPI_BARRIER( comm2d, ierr )
[1]	1614	ENDIF
[493]	1615
[1]	1616	ENDIF
	1617	#else
	1618	#if defined( __netcdf )
[1327]	1619	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
	1620	id_var_do2d(av,if), &
[1551]	1621	local_2d(nxl:nxr+1,nzb_do:nzt_do), &
[1327]	1622	start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
[1551]	1623	count = (/ nx+2, 1, nzt_do-nzb_do+1, 1 /) )
[1783]	1624	CALL netcdf_handle_error( 'data_output_2d', 451 )
[1]	1625	#endif
	1626	#endif
	1627	do2d_xz_n = do2d_xz_n + 1
	1628
	1629	CASE ( 'yz' )
	1630	!
[1308]	1631	!-- Update the netCDF yz cross section time axis.
	1632	!-- In case of parallel output, this is only done by PE0
	1633	!-- to increase the performance.
	1634	IF ( simulated_time /= do2d_yz_last_time(av) ) THEN
	1635	do2d_yz_time_count(av) = do2d_yz_time_count(av) + 1
	1636	do2d_yz_last_time(av) = simulated_time
	1637	IF ( myid == 0 ) THEN
[1327]	1638	IF ( .NOT. data_output_2d_on_each_pe &
	1639	.OR. netcdf_data_format > 4 ) &
[493]	1640	THEN
[1]	1641	#if defined( __netcdf )
	1642	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
	1643	id_var_time_yz(av), &
[291]	1644	(/ time_since_reference_point /), &
[1]	1645	start = (/ do2d_yz_time_count(av) /), &
	1646	count = (/ 1 /) )
[1783]	1647	CALL netcdf_handle_error( 'data_output_2d', 59 )
[1]	1648	#endif
	1649	ENDIF
	1650	ENDIF
[1308]	1651	ENDIF
[493]	1652
[1]	1653	!
	1654	!-- If required, carry out averaging along x
[1960]	1655	IF ( section(is,s_ind) == -1 ) THEN
[1]	1656
[1551]	1657	ALLOCATE( local_2d_l(nysg:nyng,nzb_do:nzt_do) )
[1353]	1658	local_2d_l = 0.0_wp
[1551]	1659	ngp = ( nyng-nysg+1 ) * ( nzt_do-nzb_do+1 )
[1]	1660	!
	1661	!-- First local averaging on the PE
[1551]	1662	DO k = nzb_do, nzt_do
[667]	1663	DO j = nysg, nyng
[1]	1664	DO i = nxl, nxr
[1320]	1665	local_2d_l(j,k) = local_2d_l(j,k) + &
[1]	1666	local_pf(i,j,k)
	1667	ENDDO
	1668	ENDDO
	1669	ENDDO
	1670	#if defined( __parallel )
	1671	!
	1672	!-- Now do the averaging over all PEs along x
[622]	1673	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[1551]	1674	CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb_do), &
	1675	local_2d(nysg,nzb_do), ngp, MPI_REAL, &
[1]	1676	MPI_SUM, comm1dx, ierr )
	1677	#else
	1678	local_2d = local_2d_l
	1679	#endif
[1353]	1680	local_2d = local_2d / ( nx + 1.0_wp )
[1]	1681
	1682	DEALLOCATE( local_2d_l )
	1683
	1684	ELSE
	1685	!
	1686	!-- Just store the respective section on the local array
	1687	!-- (but only if it is available on this PE!)
[1960]	1688	IF ( section(is,s_ind) >= nxl .AND. section(is,s_ind) <= nxr ) &
[1]	1689	THEN
[1960]	1690	local_2d = local_pf(section(is,s_ind),:,nzb_do:nzt_do)
[1]	1691	ENDIF
	1692
	1693	ENDIF
	1694
	1695	#if defined( __parallel )
[1327]	1696	IF ( netcdf_data_format > 4 ) THEN
[1]	1697	!
[1031]	1698	!-- Output in netCDF4/HDF5 format.
[493]	1699	!-- Output only on those PEs where the respective cross
	1700	!-- sections reside. Cross sections averaged along x are
	1701	!-- output on the respective first PE along x (myidx=0).
[1960]	1702	IF ( ( section(is,s_ind) >= nxl .AND. &
	1703	section(is,s_ind) <= nxr ) .OR. &
	1704	( section(is,s_ind) == -1 .AND. myidx == 0 ) ) THEN
[1]	1705	#if defined( __netcdf )
[493]	1706	!
[1308]	1707	!-- For parallel output, all cross sections are first
	1708	!-- stored here on a local array and will be written to the
	1709	!-- output file afterwards to increase the performance.
	1710	DO j = nysg, nyng
[1551]	1711	DO k = nzb_do, nzt_do
[1308]	1712	local_2d_sections_l(is,j,k) = local_2d(j,k)
	1713	ENDDO
	1714	ENDDO
[1]	1715	#endif
	1716	ENDIF
	1717
	1718	ELSE
	1719
[493]	1720	IF ( data_output_2d_on_each_pe ) THEN
[1]	1721	!
[493]	1722	!-- Output of partial arrays on each PE. If the cross
	1723	!-- section does not reside on the PE, output special
	1724	!-- index values.
	1725	#if defined( __netcdf )
[1327]	1726	IF ( myid == 0 ) THEN
[1320]	1727	WRITE ( 23 ) time_since_reference_point, &
[493]	1728	do2d_yz_time_count(av), av
	1729	ENDIF
	1730	#endif
[759]	1731	DO i = 0, io_blocks-1
	1732	IF ( i == io_group ) THEN
[1960]	1733	IF ( ( section(is,s_ind) >= nxl .AND. &
	1734	section(is,s_ind) <= nxr ) .OR. &
	1735	( section(is,s_ind) == -1 .AND. &
[1320]	1736	nxl-1 == -1 ) ) &
[759]	1737	THEN
[1551]	1738	WRITE (23) nysg, nyng, nzb_do, nzt_do, nzb, nzt+1
[759]	1739	WRITE (23) local_2d
	1740	ELSE
[1551]	1741	WRITE (23) -1, -1, -1, -1, -1, -1
[759]	1742	ENDIF
	1743	ENDIF
	1744	#if defined( __parallel )
	1745	CALL MPI_BARRIER( comm2d, ierr )
	1746	#endif
	1747	ENDDO
[493]	1748
	1749	ELSE
[1]	1750	!
[493]	1751	!-- PE0 receives partial arrays from all processors of the
	1752	!-- respective cross section and outputs them. Here a
	1753	!-- barrier has to be set, because otherwise
	1754	!-- "-MPI- FATAL: Remote protocol queue full" may occur.
	1755	CALL MPI_BARRIER( comm2d, ierr )
	1756
[1551]	1757	ngp = ( nyng-nysg+1 ) * ( nzt_do-nzb_do+1 )
[493]	1758	IF ( myid == 0 ) THEN
[1]	1759	!
[493]	1760	!-- Local array can be relocated directly.
[1960]	1761	IF ( ( section(is,s_ind) >= nxl .AND. &
	1762	section(is,s_ind) <= nxr ) .OR. &
	1763	( section(is,s_ind) == -1 .AND. nxl-1 == -1 ) ) &
[493]	1764	THEN
[1551]	1765	total_2d(nysg:nyng,nzb_do:nzt_do) = local_2d
[493]	1766	ENDIF
[1]	1767	!
[493]	1768	!-- Receive data from all other PEs.
	1769	DO n = 1, numprocs-1
	1770	!
	1771	!-- Receive index limits first, then array.
	1772	!-- Index limits are received in arbitrary order from
	1773	!-- the PEs.
[1320]	1774	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
	1775	MPI_ANY_SOURCE, 0, comm2d, &
[1]	1776	status, ierr )
[493]	1777	!
	1778	!-- Not all PEs have data for YZ-cross-section.
	1779	IF ( ind(1) /= -9999 ) THEN
	1780	sender = status(MPI_SOURCE)
	1781	DEALLOCATE( local_2d )
[1320]	1782	ALLOCATE( local_2d(ind(1):ind(2), &
[493]	1783	ind(3):ind(4)) )
	1784	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
	1785	MPI_REAL, sender, 1, comm2d, &
	1786	status, ierr )
[1320]	1787	total_2d(ind(1):ind(2),ind(3):ind(4)) = &
[493]	1788	local_2d
	1789	ENDIF
	1790	ENDDO
	1791	!
	1792	!-- Relocate the local array for the next loop increment
	1793	DEALLOCATE( local_2d )
[1551]	1794	ALLOCATE( local_2d(nysg:nyng,nzb_do:nzt_do) )
[1]	1795
	1796	#if defined( __netcdf )
[1327]	1797	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
	1798	id_var_do2d(av,if), &
[1551]	1799	total_2d(0:ny+1,nzb_do:nzt_do),&
[1327]	1800	start = (/ is, 1, 1, do2d_yz_time_count(av) /), &
[1551]	1801	count = (/ 1, ny+2, nzt_do-nzb_do+1, 1 /) )
[1783]	1802	CALL netcdf_handle_error( 'data_output_2d', 61 )
[1]	1803	#endif
	1804
[493]	1805	ELSE
[1]	1806	!
[493]	1807	!-- If the cross section resides on the PE, send the
	1808	!-- local index limits, otherwise send -9999 to PE0.
[1960]	1809	IF ( ( section(is,s_ind) >= nxl .AND. &
	1810	section(is,s_ind) <= nxr ) .OR. &
	1811	( section(is,s_ind) == -1 .AND. nxl-1 == -1 ) ) &
[493]	1812	THEN
[667]	1813	ind(1) = nysg; ind(2) = nyng
[1551]	1814	ind(3) = nzb_do; ind(4) = nzt_do
[493]	1815	ELSE
	1816	ind(1) = -9999; ind(2) = -9999
	1817	ind(3) = -9999; ind(4) = -9999
	1818	ENDIF
[1320]	1819	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
[493]	1820	comm2d, ierr )
	1821	!
	1822	!-- If applicable, send data to PE0.
	1823	IF ( ind(1) /= -9999 ) THEN
[1551]	1824	CALL MPI_SEND( local_2d(nysg,nzb_do), ngp, &
[493]	1825	MPI_REAL, 0, 1, comm2d, ierr )
	1826	ENDIF
[1]	1827	ENDIF
	1828	!
[493]	1829	!-- A barrier has to be set, because otherwise some PEs may
	1830	!-- proceed too fast so that PE0 may receive wrong data on
	1831	!-- tag 0
	1832	CALL MPI_BARRIER( comm2d, ierr )
[1]	1833	ENDIF
[493]	1834
[1]	1835	ENDIF
	1836	#else
	1837	#if defined( __netcdf )
[1327]	1838	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
	1839	id_var_do2d(av,if), &
[1551]	1840	local_2d(nys:nyn+1,nzb_do:nzt_do), &
[1327]	1841	start = (/ is, 1, 1, do2d_xz_time_count(av) /), &
[1551]	1842	count = (/ 1, ny+2, nzt_do-nzb_do+1, 1 /) )
[1783]	1843	CALL netcdf_handle_error( 'data_output_2d', 452 )
[1]	1844	#endif
	1845	#endif
	1846	do2d_yz_n = do2d_yz_n + 1
	1847
	1848	END SELECT
	1849
	1850	is = is + 1
	1851	ENDDO loop1
	1852
[1308]	1853	!
	1854	!-- For parallel output, all data were collected before on a local array
	1855	!-- and are written now to the netcdf file. This must be done to increase
	1856	!-- the performance of the parallel output.
	1857	#if defined( __netcdf )
[1327]	1858	IF ( netcdf_data_format > 4 ) THEN
[1308]	1859
	1860	SELECT CASE ( mode )
	1861
	1862	CASE ( 'xy' )
	1863	IF ( two_d ) THEN
[1703]	1864	nis = 1
	1865	two_d = .FALSE.
[1308]	1866	ELSE
[1703]	1867	nis = ns
[1308]	1868	ENDIF
	1869	!
	1870	!-- Do not output redundant ghost point data except for the
	1871	!-- boundaries of the total domain.
	1872	IF ( nxr == nx .AND. nyn /= ny ) THEN
	1873	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1874	id_var_do2d(av,if), &
	1875	local_2d_sections(nxl:nxr+1, &
[1703]	1876	nys:nyn,1:nis), &
[1308]	1877	start = (/ nxl+1, nys+1, 1, &
	1878	do2d_xy_time_count(av) /), &
	1879	count = (/ nxr-nxl+2, &
[1703]	1880	nyn-nys+1, nis, 1 &
[1308]	1881	/) )
	1882	ELSEIF ( nxr /= nx .AND. nyn == ny ) THEN
	1883	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1884	id_var_do2d(av,if), &
	1885	local_2d_sections(nxl:nxr, &
[1703]	1886	nys:nyn+1,1:nis), &
[1308]	1887	start = (/ nxl+1, nys+1, 1, &
	1888	do2d_xy_time_count(av) /), &
	1889	count = (/ nxr-nxl+1, &
[1703]	1890	nyn-nys+2, nis, 1 &
[1308]	1891	/) )
	1892	ELSEIF ( nxr == nx .AND. nyn == ny ) THEN
	1893	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1894	id_var_do2d(av,if), &
	1895	local_2d_sections(nxl:nxr+1, &
[1703]	1896	nys:nyn+1,1:nis), &
[1308]	1897	start = (/ nxl+1, nys+1, 1, &
	1898	do2d_xy_time_count(av) /), &
	1899	count = (/ nxr-nxl+2, &
[1703]	1900	nyn-nys+2, nis, 1 &
[1308]	1901	/) )
	1902	ELSE
	1903	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
	1904	id_var_do2d(av,if), &
	1905	local_2d_sections(nxl:nxr, &
[1703]	1906	nys:nyn,1:nis), &
[1308]	1907	start = (/ nxl+1, nys+1, 1, &
	1908	do2d_xy_time_count(av) /), &
	1909	count = (/ nxr-nxl+1, &
[1703]	1910	nyn-nys+1, nis, 1 &
[1308]	1911	/) )
	1912	ENDIF
	1913
[1783]	1914	CALL netcdf_handle_error( 'data_output_2d', 55 )
[1308]	1915
	1916	CASE ( 'xz' )
	1917	!
	1918	!-- First, all PEs get the information of all cross-sections.
	1919	!-- Then the data are written to the output file by all PEs
	1920	!-- while NF90_COLLECTIVE is set in subroutine
	1921	!-- define_netcdf_header. Although redundant information are
	1922	!-- written to the output file in that case, the performance
	1923	!-- is significantly better compared to the case where only
	1924	!-- the first row of PEs in x-direction (myidx = 0) is given
	1925	!-- the output while NF90_INDEPENDENT is set.
	1926	IF ( npey /= 1 ) THEN
	1927
	1928	#if defined( __parallel )
	1929	!
	1930	!-- Distribute data over all PEs along y
[1551]	1931	ngp = ( nxrg-nxlg+1 ) * ( nzt_do-nzb_do+1 ) * ns
[1308]	1932	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[1551]	1933	CALL MPI_ALLREDUCE( local_2d_sections_l(nxlg,1,nzb_do), &
	1934	local_2d_sections(nxlg,1,nzb_do), &
[1308]	1935	ngp, MPI_REAL, MPI_SUM, comm1dy, &
	1936	ierr )
	1937	#else
	1938	local_2d_sections = local_2d_sections_l
	1939	#endif
	1940	ENDIF
	1941	!
	1942	!-- Do not output redundant ghost point data except for the
	1943	!-- boundaries of the total domain.
	1944	IF ( nxr == nx ) THEN
	1945	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
	1946	id_var_do2d(av,if), &
	1947	local_2d_sections(nxl:nxr+1,1:ns, &
[1551]	1948	nzb_do:nzt_do), &
[1308]	1949	start = (/ nxl+1, 1, 1, &
	1950	do2d_xz_time_count(av) /), &
[1551]	1951	count = (/ nxr-nxl+2, ns, nzt_do-nzb_do+1, &
[1308]	1952	1 /) )
	1953	ELSE
	1954	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
	1955	id_var_do2d(av,if), &
	1956	local_2d_sections(nxl:nxr,1:ns, &
[1551]	1957	nzb_do:nzt_do), &
[1308]	1958	start = (/ nxl+1, 1, 1, &
	1959	do2d_xz_time_count(av) /), &
[1551]	1960	count = (/ nxr-nxl+1, ns, nzt_do-nzb_do+1, &
[1308]	1961	1 /) )
	1962	ENDIF
	1963
[1783]	1964	CALL netcdf_handle_error( 'data_output_2d', 57 )
[1308]	1965
	1966	CASE ( 'yz' )
	1967	!
	1968	!-- First, all PEs get the information of all cross-sections.
	1969	!-- Then the data are written to the output file by all PEs
	1970	!-- while NF90_COLLECTIVE is set in subroutine
	1971	!-- define_netcdf_header. Although redundant information are
	1972	!-- written to the output file in that case, the performance
	1973	!-- is significantly better compared to the case where only
	1974	!-- the first row of PEs in y-direction (myidy = 0) is given
	1975	!-- the output while NF90_INDEPENDENT is set.
	1976	IF ( npex /= 1 ) THEN
	1977
	1978	#if defined( __parallel )
	1979	!
	1980	!-- Distribute data over all PEs along x
	1981	ngp = ( nyng-nysg+1 ) * ( nzt-nzb + 2 ) * ns
	1982	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[1551]	1983	CALL MPI_ALLREDUCE( local_2d_sections_l(1,nysg,nzb_do), &
	1984	local_2d_sections(1,nysg,nzb_do), &
[1308]	1985	ngp, MPI_REAL, MPI_SUM, comm1dx, &
	1986	ierr )
	1987	#else
	1988	local_2d_sections = local_2d_sections_l
	1989	#endif
	1990	ENDIF
	1991	!
	1992	!-- Do not output redundant ghost point data except for the
	1993	!-- boundaries of the total domain.
	1994	IF ( nyn == ny ) THEN
	1995	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
	1996	id_var_do2d(av,if), &
	1997	local_2d_sections(1:ns, &
[1551]	1998	nys:nyn+1,nzb_do:nzt_do), &
[1308]	1999	start = (/ 1, nys+1, 1, &
	2000	do2d_yz_time_count(av) /), &
	2001	count = (/ ns, nyn-nys+2, &
[1551]	2002	nzt_do-nzb_do+1, 1 /) )
[1308]	2003	ELSE
	2004	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
	2005	id_var_do2d(av,if), &
	2006	local_2d_sections(1:ns,nys:nyn, &
[1551]	2007	nzb_do:nzt_do), &
[1308]	2008	start = (/ 1, nys+1, 1, &
	2009	do2d_yz_time_count(av) /), &
	2010	count = (/ ns, nyn-nys+1, &
[1551]	2011	nzt_do-nzb_do+1, 1 /) )
[1308]	2012	ENDIF
	2013
[1783]	2014	CALL netcdf_handle_error( 'data_output_2d', 60 )
[1308]	2015
	2016	CASE DEFAULT
	2017	message_string = 'unknown cross-section: ' // TRIM( mode )
	2018	CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
	2019
	2020	END SELECT
	2021
	2022	ENDIF
[1311]	2023	#endif
[1]	2024	ENDIF
	2025
	2026	if = if + 1
	2027	l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
	2028	do2d_mode = do2d(av,if)(l-1:l)
	2029
	2030	ENDDO
	2031
	2032	!
	2033	!-- Deallocate temporary arrays.
	2034	IF ( ALLOCATED( level_z ) ) DEALLOCATE( level_z )
[1308]	2035	IF ( netcdf_data_format > 4 ) THEN
	2036	DEALLOCATE( local_pf, local_2d, local_2d_sections )
	2037	IF( mode == 'xz' .OR. mode == 'yz' ) DEALLOCATE( local_2d_sections_l )
	2038	ENDIF
[1]	2039	#if defined( __parallel )
	2040	IF ( .NOT. data_output_2d_on_each_pe .AND. myid == 0 ) THEN
	2041	DEALLOCATE( total_2d )
	2042	ENDIF
	2043	#endif
	2044
	2045	!
	2046	!-- Close plot output file.
[1960]	2047	file_id = 20 + s_ind
[1]	2048
	2049	IF ( data_output_2d_on_each_pe ) THEN
[759]	2050	DO i = 0, io_blocks-1
	2051	IF ( i == io_group ) THEN
	2052	CALL close_file( file_id )
	2053	ENDIF
	2054	#if defined( __parallel )
	2055	CALL MPI_BARRIER( comm2d, ierr )
	2056	#endif
	2057	ENDDO
[1]	2058	ELSE
	2059	IF ( myid == 0 ) CALL close_file( file_id )
	2060	ENDIF
	2061
[1318]	2062	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
[1]	2063
	2064	END SUBROUTINE data_output_2d

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |