Home

Context Navigation

source: palm/trunk/SOURCE/netcdf_interface_mod.f90 @ 1941

Last change on this file since 1941 was 1851, checked in by maronga, 9 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 219.0 KB

Line
1	!> @file netcdf_interface_mod.f90
2	!--------------------------------------------------------------------------------!
3	! This file is part of PALM.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the terms
6	! of the GNU General Public License as published by the Free Software Foundation,
7	! either version 3 of the License, or (at your option) any later version.
8	!
9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
10	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12	!
13	! You should have received a copy of the GNU General Public License along with
14	! PALM. If not, see <http://www.gnu.org/licenses/>.
15	!
16	! Copyright 1997-2016 Leibniz Universitaet Hannover
17	!--------------------------------------------------------------------------------!
18	!
19	! Current revisions:
20	! ------------------
21	!
22	!
23	! Former revisions:
24	! -----------------
25	! $Id: netcdf_interface_mod.f90 1851 2016-04-08 13:32:50Z raasch $
26	!
27	! 1850 2016-04-08 13:29:27Z maronga
28	! Module renamed
29	!
30	!
31	! 1833 2016-04-07 14:23:03Z raasch
32	! spectrum renamed spectra_mod
33	!
34	! 1786 2016-03-08 05:49:27Z raasch
35	! Bugfix: id_var_time_sp made public
36	!
37	! 1783 2016-03-06 18:36:17Z raasch
38	! netcdf interface has been modularized, former file netcdf renamed to
39	! netcdf_interface, creation of netcdf-dimensions and -variables moved to
40	! specific new subroutines create_netcdf_dim and create_netcdf_var,
41	! compression (deflation) of variables implemented,
42	! ibmy special cpp directive removed
43	!
44	! 1745 2016-02-05 13:06:51Z gronemeier
45	! Bugfix: recalculating ntdim_3d, ntdim_2d_xy/xz/yz when checking the
46	! extensibility of an existing file (only when using parallel NetCDF).
47	!
48	! 1691 2015-10-26 16:17:44Z maronga
49	! Added output of radiative heating rates for RRTMG. Corrected output of
50	! radiative fluxes
51	!
52	! 1682 2015-10-07 23:56:08Z knoop
53	! Code annotations made doxygen readable
54	!
55	! 1596 2015-05-21 09:34:28Z gronemeier
56	! Bugfix in masked data output. Read 'zu_3d' when trying to extend masked data
57	!
58	! 1551 2015-03-03 14:18:16Z maronga
59	! Added support for land surface model and radiation model output. In the course
60	! of this action a new vertical grid zs (soil) was introduced.
61	!
62	! 1353 2014-04-08 15:21:23Z heinze
63	! REAL constants provided with KIND-attribute
64	!
65	! 1322 2014-03-20 16:38:49Z raasch
66	! Forgotten ONLY-attribute added to USE-statements
67	!
68	! 1320 2014-03-20 08:40:49Z raasch
69	! ONLY-attribute added to USE-statements,
70	! kind-parameters added to all INTEGER and REAL declaration statements,
71	! kinds are defined in new module kinds,
72	! revision history before 2012 removed,
73	! comment fields (!:) to be used for variable explanations added to
74	! all variable declaration statements
75	!
76	! 1308 2014-03-13 14:58:42Z fricke
77	! +ntime_count, oldmode
78	! Adjust NF90_CREATE and NF90_OPEN statement for parallel output
79	! To increase the performance for parallel output, the following is done:
80	! - Limit time dimension
81	! - Values of axis data are only written by PE0
82	! - No fill is set for all variables
83	! Check the number of output time levels for restart jobs
84	!
85	! 1206 2013-07-18 12:49:16Z witha
86	! Bugfix: typo in preprocessor directive in subroutine open_write_netcdf_file
87	!
88	! 1092 2013-02-02 11:24:22Z raasch
89	! unused variables removed
90	!
91	! 1053 2012-11-13 17:11:03Z hoffmann
92	! +qr, nr, prr
93	!
94	! 1036 2012-10-22 13:43:42Z raasch
95	! code put under GPL (PALM 3.9)
96	!
97	! 1031 2012-10-19 14:35:30Z raasch
98	! netCDF4 without parallel file support implemented, new routines
99	! create_netcdf_file and open_write_netcdf_file at end
100	!
101	! 992 2012-09-05 15:08:26Z hoffmann
102	! Removal of the informative messages PA0352 and PA0353.
103	!
104	! 983 2012-08-21 14:17:57Z hoffmann
105	! Bugfix in cross_profiles.
106	!
107	! 964 2012-07-26 09:14:24Z raasch
108	! rev 951 and 959 reformatted
109	!
110	! 959 2012-07-24 13:13:41Z hoffmann
111	! Bugfix in cross_profiles. It is not allowed to arrange more than 100
112	! profiles with cross_profiles.
113	!
114	! 951 2012-07-19 14:22:52Z hoffmann
115	! cross_profiles, profile_rows, profile_columns are written to netCDF header
116	!
117	! Revision 1.1 2005/05/18 15:37:16 raasch
118	! Initial revision
119	!
120	!
121	! Description:
122	! ------------
123	!> In case of extend = .FALSE.:
124	!> Define all necessary dimensions, axes and variables for the different
125	!> netCDF datasets. This subroutine is called from check_open after a new
126	!> dataset is created. It leaves the open netCDF files ready to write.
127	!>
128	!> In case of extend = .TRUE.:
129	!> Find out if dimensions and variables of an existing file match the values
130	!> of the actual run. If so, get all necessary information (ids, etc.) from
131	!> this file.
132	!>
133	!> Parameter av can assume values 0 (non-averaged data) and 1 (time averaged
134	!> data)
135	!>
136	!> @todo calculation of output time levels for parallel NetCDF still does not
137	!> cover every exception (change of dt_do, end_time in restart)
138	!------------------------------------------------------------------------------!
139	MODULE netcdf_interface
140
141	USE control_parameters, ONLY: max_masks
142	USE kinds
143	#if defined( __netcdf )
144	USE NETCDF
145	#endif
146
147	PRIVATE
148
149	INTEGER(iwp), PARAMETER :: dopr_norm_num = 7, dopts_num = 29, dots_max = 100
150
151	CHARACTER (LEN=6), DIMENSION(dopr_norm_num) :: dopr_norm_names = &
152	(/ 'wpt0 ', 'ws2 ', 'tsw2 ', 'ws3 ', 'ws2tsw', 'wstsw2', &
153	'z_i ' /)
154
155	CHARACTER (LEN=6), DIMENSION(dopr_norm_num) :: dopr_norm_longnames = &
156	(/ 'wpt0 ', 'w2 ', 'tw2 ', 'w3 ', 'w2tw', 'wt*w2', &
157	'z_i ' /)
158
159	CHARACTER (LEN=7), DIMENSION(dopts_num) :: dopts_label = &
160	(/ 'tnpt ', 'x_ ', 'y_ ', 'z_ ', 'z_abs ', 'u ', &
161	'v ', 'w ', 'u" ', 'v" ', 'w" ', 'npt_up ', &
162	'w_up ', 'w_down ', 'radius ', 'r_min ', 'r_max ', 'npt_max', &
163	'npt_min', 'x2 ', 'y2 ', 'z2 ', 'u2 ', 'v*2 ', &
164	'w2 ', 'u"2 ', 'v"2 ', 'w"2 ', 'npt2 ' /)
165
166	CHARACTER (LEN=7), DIMENSION(dopts_num) :: dopts_unit = &
167	(/ 'number ', 'm ', 'm ', 'm ', 'm ', 'm/s ', &
168	'm/s ', 'm/s ', 'm/s ', 'm/s ', 'm/s ', 'number ', &
169	'm/s ', 'm/s ', 'm ', 'm ', 'm ', 'number ', &
170	'number ', 'm2 ', 'm2 ', 'm2 ', 'm2/s2 ', 'm2/s2 ', &
171	'm2/s2 ', 'm2/s2 ', 'm2/s2 ', 'm2/s2 ', 'number2' /)
172
173	INTEGER(iwp) :: dots_num = 31 !< number of timeseries defined by default
174	INTEGER(iwp) :: dots_soil = 24 !< starting index for soil-timeseries
175	INTEGER(iwp) :: dots_rad = 32 !< starting index for radiation-timeseries
176
177	CHARACTER (LEN=13), DIMENSION(dots_max) :: dots_label = &
178	(/ 'E ', 'E* ', 'dt ', &
179	'u* ', 'th* ', 'umax ', &
180	'vmax ', 'wmax ', 'div_new ', &
181	'div_old ', 'z_i_wpt ', 'z_i_pt ', &
182	'w* ', 'w"pt"0 ', 'w"pt" ', &
183	'wpt ', 'pt(0) ', 'pt(z_mo) ', &
184	'w"u"0 ', 'w"v"0 ', 'w"q"0 ', &
185	'ol ', 'q* ', 'ghf_eb ', &
186	'shf_eb ', 'qsws_eb ', 'qsws_liq_eb ', &
187	'qsws_soil_eb ', 'qsws_veg_eb ', 'r_a ', &
188	'r_s ', 'rad_net ', 'rad_lw_in ', &
189	'rad_lw_out ', 'rad_sw_in ', 'rad_sw_out ', &
190	'rrtm_aldif ', 'rrtm_aldir ', 'rrtm_asdif ', &
191	'rrtm_asdir ', &
192	( 'unknown ', i9 = 1, dots_max-40 ) /)
193
194	CHARACTER (LEN=13), DIMENSION(dots_max) :: dots_unit = &
195	(/ 'm2/s2 ', 'm2/s2 ', 's ', &
196	'm/s ', 'K ', 'm/s ', &
197	'm/s ', 'm/s ', 's-1 ', &
198	's-1 ', 'm ', 'm ', &
199	'm/s ', 'K m/s ', 'K m/s ', &
200	'K m/s ', 'K ', 'K ', &
201	'm2/s2 ', 'm2/s2 ', 'kg m/s ', &
202	'm ', 'kg/kg ', ' ', &
203	' ', ' ', ' ', &
204	' ', 'W/m2 ', 's/m ', &
205	' ', 'W/m2 ', 'W/m2 ', &
206	'W/m2 ', 'W/m2 ', 'W/m2 ', &
207	' ', ' ', ' ', &
208	' ', &
209	( 'unknown ', i9 = 1, dots_max-40 ) /)
210
211	CHARACTER (LEN=9), DIMENSION(300) :: dopr_unit = 'unknown'
212
213	CHARACTER (LEN=7), DIMENSION(0:1,100) :: do2d_unit, do3d_unit
214
215	CHARACTER (LEN=16), DIMENSION(25) :: prt_var_names = &
216	(/ 'pt_age ', 'pt_dvrp_size ', 'pt_origin_x ', &
217	'pt_origin_y ', 'pt_origin_z ', 'pt_radius ', &
218	'pt_speed_x ', 'pt_speed_y ', 'pt_speed_z ', &
219	'pt_weight_factor', 'pt_x ', 'pt_y ', &
220	'pt_z ', 'pt_color ', 'pt_group ', &
221	'pt_tailpoints ', 'pt_tail_id ', 'pt_density_ratio', &
222	'pt_exp_arg ', 'pt_exp_term ', 'not_used ', &
223	'not_used ', 'not_used ', 'not_used ', &
224	'not_used ' /)
225
226	CHARACTER (LEN=16), DIMENSION(25) :: prt_var_units = &
227	(/ 'seconds ', 'meters ', 'meters ', &
228	'meters ', 'meters ', 'meters ', &
229	'm/s ', 'm/s ', 'm/s ', &
230	'factor ', 'meters ', 'meters ', &
231	'meters ', 'none ', 'none ', &
232	'none ', 'none ', 'ratio ', &
233	'none ', 'none ', 'not_used ', &
234	'not_used ', 'not_used ', 'not_used ', &
235	'not_used ' /)
236
237	CHARACTER(LEN=20), DIMENSION(11) :: netcdf_precision = ' '
238	CHARACTER(LEN=40) :: netcdf_data_format_string
239
240	INTEGER(iwp) :: id_dim_prtnum, id_dim_time_pr, id_dim_time_prt, &
241	id_dim_time_pts, id_dim_time_sp, id_dim_time_ts, id_dim_x_sp, &
242	id_dim_y_sp, id_dim_zu_sp, id_dim_zw_sp, id_set_pr, &
243	id_set_prt, id_set_pts, id_set_sp, id_set_ts, id_var_prtnum, &
244	id_var_rnop_prt, id_var_time_pr, id_var_time_prt, &
245	id_var_time_pts, id_var_time_sp, id_var_time_ts, id_var_x_sp, &
246	id_var_y_sp, id_var_zu_sp, id_var_zw_sp, nc_stat
247
248	INTEGER(iwp), DIMENSION(0:1) :: id_dim_time_xy, id_dim_time_xz, &
249	id_dim_time_yz, id_dim_time_3d, id_dim_x_xy, id_dim_xu_xy, &
250	id_dim_x_xz, id_dim_xu_xz, id_dim_x_yz, id_dim_xu_yz, &
251	id_dim_x_3d, id_dim_xu_3d, id_dim_y_xy, id_dim_yv_xy, &
252	id_dim_y_xz, id_dim_yv_xz, id_dim_y_yz, id_dim_yv_yz, &
253	id_dim_y_3d, id_dim_yv_3d, id_dim_zs_xy, id_dim_zs_xz, &
254	id_dim_zs_yz, id_dim_zs_3d, id_dim_zu_xy, id_dim_zu1_xy, &
255	id_dim_zu_xz, id_dim_zu_yz, id_dim_zu_3d, id_dim_zw_xy, &
256	id_dim_zw_xz, id_dim_zw_yz, id_dim_zw_3d, id_set_xy, &
257	id_set_xz, id_set_yz, id_set_3d, id_var_ind_x_yz, &
258	id_var_ind_y_xz, id_var_ind_z_xy, id_var_time_xy, &
259	id_var_time_xz, id_var_time_yz, id_var_time_3d, id_var_x_xy, &
260	id_var_xu_xy, id_var_x_xz, id_var_xu_xz, id_var_x_yz, &
261	id_var_xu_yz, id_var_x_3d, id_var_xu_3d, id_var_y_xy, &
262	id_var_yv_xy, id_var_y_xz, id_var_yv_xz, id_var_y_yz, &
263	id_var_yv_yz, id_var_y_3d, id_var_yv_3d, id_var_zs_xy, &
264	id_var_zs_xz, id_var_zs_yz, id_var_zs_3d, id_var_zusi_xy, &
265	id_var_zusi_3d, id_var_zu_xy, id_var_zu1_xy, id_var_zu_xz, &
266	id_var_zu_yz, id_var_zu_3d, id_var_zwwi_xy, id_var_zwwi_3d, &
267	id_var_zw_xy, id_var_zw_xz, id_var_zw_yz, id_var_zw_3d
268
269	INTEGER :: netcdf_data_format = 2 !< NetCDF3 64bit offset format
270	INTEGER :: netcdf_deflate = 0 !< NetCDF compression, default: no
271	!< compression
272
273	INTEGER(iwp), DIMENSION(10) :: id_var_dospx, id_var_dospy
274	INTEGER(iwp), DIMENSION(20) :: id_var_prt
275	INTEGER(iwp), DIMENSION(11) :: nc_precision
276	INTEGER(iwp), DIMENSION(dopr_norm_num) :: id_var_norm_dopr
277
278	INTEGER(iwp), DIMENSION(dopts_num,0:10) :: id_var_dopts
279	INTEGER(iwp), DIMENSION(0:1,100) :: id_var_do2d, id_var_do3d
280	INTEGER(iwp), DIMENSION(100,0:9) :: id_dim_z_pr, id_var_dopr, &
281	id_var_z_pr
282	INTEGER(iwp), DIMENSION(dots_max,0:9) :: id_var_dots
283
284	!
285	!-- Masked output
286	CHARACTER (LEN=7), DIMENSION(max_masks,0:1,100) :: domask_unit
287
288	LOGICAL :: output_for_t0 = .FALSE.
289
290	INTEGER(iwp), DIMENSION(1:max_masks,0:1) :: id_dim_time_mask, id_dim_x_mask, &
291	id_dim_xu_mask, id_dim_y_mask, id_dim_yv_mask, id_dim_zs_mask, &
292	id_dim_zu_mask, id_dim_zw_mask, &
293	id_set_mask, &
294	id_var_time_mask, id_var_x_mask, id_var_xu_mask, &
295	id_var_y_mask, id_var_yv_mask, id_var_zs_mask, &
296	id_var_zu_mask, id_var_zw_mask, &
297	id_var_zusi_mask, id_var_zwwi_mask
298
299	INTEGER(iwp), DIMENSION(1:max_masks,0:1,100) :: id_var_domask
300
301
302	PUBLIC domask_unit, dopr_unit, dopts_num, dots_label, dots_max, dots_num, &
303	dots_rad, dots_soil, dots_unit, do2d_unit, do3d_unit, id_set_mask, &
304	id_set_pr, id_set_prt, id_set_pts, id_set_sp, id_set_ts, &
305	id_set_xy, id_set_xz, id_set_yz, id_set_3d, id_var_domask, &
306	id_var_dopr, id_var_dopts, id_var_dospx, id_var_dospy, &
307	id_var_dots, id_var_do2d, id_var_do3d, id_var_norm_dopr, &
308	id_var_time_mask, id_var_time_pr, id_var_time_pts, id_var_time_sp, &
309	id_var_time_ts, id_var_time_xy, id_var_time_xz, id_var_time_yz, &
310	id_var_time_3d, nc_stat, netcdf_data_format, &
311	netcdf_data_format_string, netcdf_deflate, netcdf_precision, &
312	output_for_t0
313
314	SAVE
315
316	INTERFACE netcdf_create_dim
317	MODULE PROCEDURE netcdf_create_dim
318	END INTERFACE netcdf_create_dim
319
320	INTERFACE netcdf_create_file
321	MODULE PROCEDURE netcdf_create_file
322	END INTERFACE netcdf_create_file
323
324	INTERFACE netcdf_create_var
325	MODULE PROCEDURE netcdf_create_var
326	END INTERFACE netcdf_create_var
327
328	INTERFACE netcdf_define_header
329	MODULE PROCEDURE netcdf_define_header
330	END INTERFACE netcdf_define_header
331
332	INTERFACE netcdf_handle_error
333	MODULE PROCEDURE netcdf_handle_error
334	END INTERFACE netcdf_handle_error
335
336	INTERFACE netcdf_open_write_file
337	MODULE PROCEDURE netcdf_open_write_file
338	END INTERFACE netcdf_open_write_file
339
340	PUBLIC netcdf_create_file, netcdf_define_header, netcdf_handle_error, &
341	netcdf_open_write_file
342
343	CONTAINS
344
345	SUBROUTINE netcdf_define_header( callmode, extend, av )
346
347	#if defined( __netcdf )
348
349	USE arrays_3d, &
350	ONLY: zu, zw
351
352	USE constants, &
353	ONLY: pi
354
355	USE control_parameters, &
356	ONLY: averaging_interval, averaging_interval_pr, &
357	data_output_pr, domask, dopr_n, &
358	dopr_time_count, dopts_time_count, dots_time_count, &
359	do2d, do2d_xz_time_count, do3d, &
360	do2d_yz_time_count, dt_data_output_av, dt_do2d_xy, dt_do2d_xz, &
361	dt_do2d_yz, dt_do3d, mask_size, do2d_xy_time_count, &
362	do3d_time_count, domask_time_count, end_time, mask_i_global, &
363	mask_j_global, mask_k_global, message_string, mid, ntdim_2d_xy, &
364	ntdim_2d_xz, ntdim_2d_yz, ntdim_3d, nz_do3d, prt_time_count, &
365	run_description_header, section, simulated_time, &
366	simulated_time_at_begin, skip_time_data_output_av, &
367	skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz, &
368	skip_time_do3d, topography
369
370	USE grid_variables, &
371	ONLY: dx, dy, zu_s_inner, zw_w_inner
372
373	USE indices, &
374	ONLY: nx, ny, nz ,nzb, nzt
375
376	USE kinds
377
378	USE land_surface_model_mod, &
379	ONLY: land_surface, nzb_soil, nzt_soil, nzs, zs
380
381	USE pegrid
382
383	USE particle_attributes, &
384	ONLY: maximum_number_of_particles, number_of_particle_groups
385
386	USE profil_parameter, &
387	ONLY: crmax, cross_profiles, dopr_index, profile_columns, profile_rows
388
389	USE spectra_mod, &
390	ONLY: averaging_interval_sp, comp_spectra_level, data_output_sp, dosp_time_count, spectra_direction
391
392	USE statistics, &
393	ONLY: hom, statistic_regions
394
395
396	IMPLICIT NONE
397
398	CHARACTER (LEN=2) :: suffix !<
399	CHARACTER (LEN=2), INTENT (IN) :: callmode !<
400	CHARACTER (LEN=3) :: suffix1 !<
401	CHARACTER (LEN=4) :: grid_x !<
402	CHARACTER (LEN=4) :: grid_y !<
403	CHARACTER (LEN=4) :: grid_z !<
404	CHARACTER (LEN=6) :: mode !<
405	CHARACTER (LEN=10) :: netcdf_var_name !<
406	CHARACTER (LEN=10) :: precision !<
407	CHARACTER (LEN=10) :: var !<
408	CHARACTER (LEN=80) :: time_average_text !<
409	CHARACTER (LEN=2000) :: char_cross_profiles !<
410	CHARACTER (LEN=2000) :: var_list !<
411	CHARACTER (LEN=2000) :: var_list_old !<
412
413	CHARACTER (LEN=100), DIMENSION(1:crmax) :: cross_profiles_adj !<
414	CHARACTER (LEN=100), DIMENSION(1:crmax) :: cross_profiles_char !<
415
416	INTEGER(iwp) :: av !<
417	INTEGER(iwp) :: cross_profiles_count !<
418	INTEGER(iwp) :: cross_profiles_maxi !<
419	INTEGER(iwp) :: delim !<
420	INTEGER(iwp) :: delim_old !<
421	INTEGER(iwp) :: file_id !<
422	INTEGER(iwp) :: i !<
423	INTEGER(iwp) :: id_last !<
424	INTEGER(iwp) :: id_x !<
425	INTEGER(iwp) :: id_y !<
426	INTEGER(iwp) :: id_z !<
427	INTEGER(iwp) :: j !<
428	INTEGER(iwp) :: k !<
429	INTEGER(iwp) :: kk !<
430	INTEGER(iwp) :: ns !<
431	INTEGER(iwp) :: ns_do !< actual value of ns for soil model data
432	INTEGER(iwp) :: ns_old !<
433	INTEGER(iwp) :: ntime_count !< number of time levels found in file
434	INTEGER(iwp) :: nz_old !<
435
436	INTEGER(iwp), SAVE :: oldmode !<
437
438	INTEGER(iwp), DIMENSION(1) :: id_dim_time_old !<
439	INTEGER(iwp), DIMENSION(1) :: id_dim_x_yz_old !<
440	INTEGER(iwp), DIMENSION(1) :: id_dim_y_xz_old !<
441	INTEGER(iwp), DIMENSION(1) :: id_dim_zu_sp_old !<
442	INTEGER(iwp), DIMENSION(1) :: id_dim_zu_xy_old !<
443	INTEGER(iwp), DIMENSION(1) :: id_dim_zu_3d_old !<
444	INTEGER(iwp), DIMENSION(1) :: id_dim_zu_mask_old !<
445
446
447	INTEGER(iwp), DIMENSION(1:crmax) :: cross_profiles_numb !<
448
449	LOGICAL :: found !<
450
451	LOGICAL, INTENT (INOUT) :: extend !<
452
453	LOGICAL, SAVE :: init_netcdf = .FALSE. !<
454
455	REAL(wp), DIMENSION(1) :: last_time_coordinate !< last time value in file
456
457	REAL(wp), DIMENSION(:), ALLOCATABLE :: netcdf_data !<
458	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: netcdf_data_2d !<
459
460	!
461	!-- Initializing actions
462	IF ( .NOT. init_netcdf ) THEN
463	!
464	!-- Check and set accuracy for netCDF output. First set default value
465	nc_precision = NF90_REAL4
466
467	i = 1
468	DO WHILE ( netcdf_precision(i) /= ' ' )
469	j = INDEX( netcdf_precision(i), '_' )
470	IF ( j == 0 ) THEN
471	WRITE ( message_string, * ) 'netcdf_precision must contain a ', &
472	'"_"netcdf_precision(', i, ')="', &
473	TRIM( netcdf_precision(i) ),'"'
474	CALL message( 'netcdf_define_header', 'PA0241', 2, 2, 0, 6, 0 )
475	ENDIF
476
477	var = netcdf_precision(i)(1:j-1)
478	precision = netcdf_precision(i)(j+1:)
479
480	IF ( precision == 'NF90_REAL4' ) THEN
481	j = NF90_REAL4
482	ELSEIF ( precision == 'NF90_REAL8' ) THEN
483	j = NF90_REAL8
484	ELSE
485	WRITE ( message_string, * ) 'illegal netcdf precision: ', &
486	'netcdf_precision(', i, ')="', &
487	TRIM( netcdf_precision(i) ),'"'
488	CALL message( 'netcdf_define_header', 'PA0242', 1, 2, 0, 6, 0 )
489	ENDIF
490
491	SELECT CASE ( var )
492	CASE ( 'xy' )
493	nc_precision(1) = j
494	CASE ( 'xz' )
495	nc_precision(2) = j
496	CASE ( 'yz' )
497	nc_precision(3) = j
498	CASE ( '2d' )
499	nc_precision(1:3) = j
500	CASE ( '3d' )
501	nc_precision(4) = j
502	CASE ( 'pr' )
503	nc_precision(5) = j
504	CASE ( 'ts' )
505	nc_precision(6) = j
506	CASE ( 'sp' )
507	nc_precision(7) = j
508	CASE ( 'prt' )
509	nc_precision(8) = j
510	CASE ( 'masks' )
511	nc_precision(11) = j
512	CASE ( 'all' )
513	nc_precision = j
514
515	CASE DEFAULT
516	WRITE ( message_string, * ) 'unknown variable in inipar ', &
517	'assignment: netcdf_precision(', i, ')="', &
518	TRIM( netcdf_precision(i) ),'"'
519	CALL message( 'netcdf_define_header', 'PA0243', 1, 2, 0, 6, 0 )
520
521	END SELECT
522
523	i = i + 1
524	IF ( i > 50 ) EXIT
525	ENDDO
526
527	!
528	!-- Check for allowed parameter range
529	IF ( netcdf_deflate < 0 .OR. netcdf_deflate > 9 ) THEN
530	WRITE ( message_string, '(A,I3,A)' ) 'netcdf_deflate out of ' // &
531	'range & given value: ', netcdf_deflate, &
532	', allowed range: 0-9'
533	CALL message( 'netcdf_define_header', 'PA0355', 2, 2, 0, 6, 0 )
534	ENDIF
535	!
536	!-- Data compression does not work with parallel NetCDF/HDF5
537	IF ( netcdf_deflate > 0 .AND. netcdf_data_format /= 3 ) THEN
538	message_string = 'netcdf_deflate reset to 0'
539	CALL message( 'netcdf_define_header', 'PA0356', 0, 1, 0, 6, 0 )
540
541	netcdf_deflate = 0
542	ENDIF
543
544	init_netcdf = .TRUE.
545
546	ENDIF
547
548	!
549	!-- Determine the mode to be processed
550	IF ( extend ) THEN
551	mode = callmode // '_ext'
552	ELSE
553	mode = callmode // '_new'
554	ENDIF
555
556	!
557	!-- Select the mode to be processed. Possibilities are 3d, mask, xy, xz, yz,
558	!-- pr and ts.
559	SELECT CASE ( mode )
560
561	CASE ( 'ma_new' )
562
563	!
564	!-- decompose actual parameter file_id (=formal parameter av) into
565	!-- mid and av
566	file_id = av
567	IF ( file_id <= 200+max_masks ) THEN
568	mid = file_id - 200
569	av = 0
570	ELSE
571	mid = file_id - (200+max_masks)
572	av = 1
573	ENDIF
574
575	!
576	!-- Define some global attributes of the dataset
577	nc_stat = NF90_PUT_ATT( id_set_mask(mid,av), NF90_GLOBAL, &
578	'Conventions', 'COARDS' )
579	CALL netcdf_handle_error( 'netcdf_define_header', 464 )
580
581	IF ( av == 0 ) THEN
582	time_average_text = ' '
583	ELSE
584	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
585	averaging_interval
586	ENDIF
587	nc_stat = NF90_PUT_ATT( id_set_mask(mid,av), NF90_GLOBAL, 'title', &
588	TRIM( run_description_header ) // &
589	TRIM( time_average_text ) )
590	CALL netcdf_handle_error( 'netcdf_define_header', 465 )
591	IF ( av == 1 ) THEN
592	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
593	nc_stat = NF90_PUT_ATT( id_set_mask(mid,av), NF90_GLOBAL, &
594	'time_avg', TRIM( time_average_text ) )
595	CALL netcdf_handle_error( 'netcdf_define_header', 466 )
596	ENDIF
597
598	!
599	!-- Define time coordinate for volume data (unlimited dimension)
600	CALL netcdf_create_dim( id_set_mask(mid,av), 'time', NF90_UNLIMITED, &
601	id_dim_time_mask(mid,av), 467 )
602	CALL netcdf_create_var( id_set_mask(mid,av), &
603	(/ id_dim_time_mask(mid,av) /), 'time', &
604	NF90_DOUBLE, id_var_time_mask(mid,av), &
605	'seconds', '', 468, 469, 000 )
606	!
607	!-- Define spatial dimensions and coordinates:
608	!-- Define vertical coordinate grid (zu grid)
609	CALL netcdf_create_dim( id_set_mask(mid,av), 'zu_3d', &
610	mask_size(mid,3), id_dim_zu_mask(mid,av), &
611	470 )
612	CALL netcdf_create_var( id_set_mask(mid,av), &
613	(/ id_dim_zu_mask(mid,av) /), 'zu_3d', &
614	NF90_DOUBLE, id_var_zu_mask(mid,av), &
615	'meters', '', 471, 472, 000 )
616	!
617	!-- Define vertical coordinate grid (zw grid)
618	CALL netcdf_create_dim( id_set_mask(mid,av), 'zw_3d', &
619	mask_size(mid,3), id_dim_zw_mask(mid,av), &
620	473 )
621	CALL netcdf_create_var( id_set_mask(mid,av), &
622	(/ id_dim_zw_mask(mid,av) /), 'zw_3d', &
623	NF90_DOUBLE, id_var_zw_mask(mid,av), &
624	'meters', '', 474, 475, 000 )
625	!
626	!-- Define x-axis (for scalar position)
627	CALL netcdf_create_dim( id_set_mask(mid,av), 'x', mask_size(mid,1), &
628	id_dim_x_mask(mid,av), 476 )
629	CALL netcdf_create_var( id_set_mask(mid,av), &
630	(/ id_dim_x_mask(mid,av) /), 'x', &
631	NF90_DOUBLE, id_var_x_mask(mid,av), &
632	'meters', '', 477, 478, 000 )
633	!
634	!-- Define x-axis (for u position)
635	CALL netcdf_create_dim( id_set_mask(mid,av), 'xu', mask_size(mid,1), &
636	id_dim_xu_mask(mid,av), 479 )
637	CALL netcdf_create_var( id_set_mask(mid,av), &
638	(/ id_dim_xu_mask(mid,av) /), 'xu', &
639	NF90_DOUBLE, id_var_xu_mask(mid,av), &
640	'meters', '', 480, 481, 000 )
641	!
642	!-- Define y-axis (for scalar position)
643	CALL netcdf_create_dim( id_set_mask(mid,av), 'y', mask_size(mid,2), &
644	id_dim_y_mask(mid,av), 482 )
645	CALL netcdf_create_var( id_set_mask(mid,av), &
646	(/ id_dim_y_mask(mid,av) /), 'y', &
647	NF90_DOUBLE, id_var_y_mask(mid,av), &
648	'meters', '', 483, 484, 000 )
649	!
650	!-- Define y-axis (for v position)
651	CALL netcdf_create_dim( id_set_mask(mid,av), 'yv', mask_size(mid,2), &
652	id_dim_yv_mask(mid,av), 485 )
653	CALL netcdf_create_var( id_set_mask(mid,av), &
654	(/ id_dim_yv_mask(mid,av) /), &
655	'yv', NF90_DOUBLE, id_var_yv_mask(mid,av), &
656	'meters', '', 486, 487, 000 )
657	!
658	!-- In case of non-flat topography define 2d-arrays containing the height
659	!-- information
660	IF ( TRIM( topography ) /= 'flat' ) THEN
661	!
662	!-- Define zusi = zu(nzb_s_inner)
663	CALL netcdf_create_var( id_set_mask(mid,av), &
664	(/ id_dim_x_mask(mid,av), &
665	id_dim_y_mask(mid,av) /), 'zusi', &
666	NF90_DOUBLE, id_var_zusi_mask(mid,av), &
667	'meters', 'zu(nzb_s_inner)', 488, 489, &
668	490 )
669	!
670	!-- Define zwwi = zw(nzb_w_inner)
671	CALL netcdf_create_var( id_set_mask(mid,av), &
672	(/ id_dim_x_mask(mid,av), &
673	id_dim_y_mask(mid,av) /), 'zwwi', &
674	NF90_DOUBLE, id_var_zwwi_mask(mid,av), &
675	'meters', 'zw(nzb_w_inner)', 491, 492, &
676	493 )
677	ENDIF
678
679	IF ( land_surface ) THEN
680	!
681	!-- Define vertical coordinate grid (zw grid)
682	CALL netcdf_create_dim( id_set_mask(mid,av), 'zs_3d', &
683	mask_size(mid,3), id_dim_zs_mask(mid,av), &
684	536 )
685	CALL netcdf_create_var( id_set_mask(mid,av), &
686	(/ id_dim_zs_mask(mid,av) /), 'zs_3d', &
687	NF90_DOUBLE, id_var_zs_mask(mid,av), &
688	'meters', '', 537, 555, 000 )
689	ENDIF
690
691	!
692	!-- Define the variables
693	var_list = ';'
694	i = 1
695
696	DO WHILE ( domask(mid,av,i)(1:1) /= ' ' )
697
698	!
699	!-- Check for the grid
700	found = .TRUE.
701	SELECT CASE ( domask(mid,av,i) )
702	!
703	!-- Most variables are defined on the scalar grid
704	CASE ( 'e', 'lpt', 'nr', 'p', 'pc', 'pr', 'prr', 'pt', 'q', &
705	'qc', 'ql', 'ql_c', 'ql_v', 'ql_vp', 'qr', 'qv', &
706	'rad_lw_cs_hr', 'rad_lw_hr', 'rad_sw_cs_hr', &
707	'rad_sw_hr', 'rho', 's', 'sa', 'vpt' )
708
709	grid_x = 'x'
710	grid_y = 'y'
711	grid_z = 'zu'
712	!
713	!-- u grid
714	CASE ( 'u' )
715
716	grid_x = 'xu'
717	grid_y = 'y'
718	grid_z = 'zu'
719	!
720	!-- v grid
721	CASE ( 'v' )
722
723	grid_x = 'x'
724	grid_y = 'yv'
725	grid_z = 'zu'
726	!
727	!-- w grid
728	CASE ( 'rad_lw_in', 'rad_lw_out', 'rad_sw_in', 'rad_sw_out', &
729	'w' )
730
731	grid_x = 'x'
732	grid_y = 'y'
733	grid_z = 'zw'
734	!
735	!-- soil grid
736	CASE ( 'm_soil', 't_soil' )
737
738	grid_x = 'x'
739	grid_y = 'y'
740	grid_z = 'zs'
741
742	CASE DEFAULT
743	!
744	!-- Check for user-defined quantities
745	CALL user_define_netcdf_grid( domask(mid,av,i), found, &
746	grid_x, grid_y, grid_z )
747
748	IF ( .NOT. found ) THEN
749	WRITE ( message_string, * ) 'no grid defined for', &
750	' variable ', TRIM( domask(mid,av,i) )
751	CALL message( 'define_netcdf_header', 'PA0244', 0, 1, 0, &
752	6, 0 )
753	ENDIF
754
755	END SELECT
756
757	!
758	!-- Select the respective dimension ids
759	IF ( grid_x == 'x' ) THEN
760	id_x = id_dim_x_mask(mid,av)
761	ELSEIF ( grid_x == 'xu' ) THEN
762	id_x = id_dim_xu_mask(mid,av)
763	ENDIF
764
765	IF ( grid_y == 'y' ) THEN
766	id_y = id_dim_y_mask(mid,av)
767	ELSEIF ( grid_y == 'yv' ) THEN
768	id_y = id_dim_yv_mask(mid,av)
769	ENDIF
770
771	IF ( grid_z == 'zu' ) THEN
772	id_z = id_dim_zu_mask(mid,av)
773	ELSEIF ( grid_z == 'zw' ) THEN
774	id_z = id_dim_zw_mask(mid,av)
775	ELSEIF ( grid_z == "zs" ) THEN
776	id_z = id_dim_zs_mask(mid,av)
777	ENDIF
778
779	!
780	!-- Define the grid
781	CALL netcdf_create_var( id_set_mask(mid,av), (/ id_x, id_y, id_z, &
782	id_dim_time_mask(mid,av) /), &
783	domask(mid,av,i), nc_precision(11), &
784	id_var_domask(mid,av,i), &
785	TRIM( domask_unit(mid,av,i) ), &
786	domask(mid,av,i), 494, 495, 496 )
787
788	var_list = TRIM( var_list ) // TRIM( domask(mid,av,i) ) // ';'
789
790	i = i + 1
791
792	ENDDO
793
794	!
795	!-- No arrays to output
796	IF ( i == 1 ) RETURN
797
798	!
799	!-- Write the list of variables as global attribute (this is used by
800	!-- restart runs and by combine_plot_fields)
801	nc_stat = NF90_PUT_ATT( id_set_mask(mid,av), NF90_GLOBAL, &
802	'VAR_LIST', var_list )
803	CALL netcdf_handle_error( 'netcdf_define_header', 497 )
804
805	!
806	!-- Leave netCDF define mode
807	nc_stat = NF90_ENDDEF( id_set_mask(mid,av) )
808	CALL netcdf_handle_error( 'netcdf_define_header', 498 )
809
810	!
811	!-- Write data for x (shifted by +dx/2) and xu axis
812	ALLOCATE( netcdf_data(mask_size(mid,1)) )
813
814	netcdf_data = ( mask_i_global(mid,:mask_size(mid,1)) + 0.5_wp ) * dx
815
816	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_x_mask(mid,av), &
817	netcdf_data, start = (/ 1 /), &
818	count = (/ mask_size(mid,1) /) )
819	CALL netcdf_handle_error( 'netcdf_define_header', 499 )
820
821	netcdf_data = mask_i_global(mid,:mask_size(mid,1)) * dx
822
823	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_xu_mask(mid,av),&
824	netcdf_data, start = (/ 1 /), &
825	count = (/ mask_size(mid,1) /) )
826	CALL netcdf_handle_error( 'netcdf_define_header', 500 )
827
828	DEALLOCATE( netcdf_data )
829
830	!
831	!-- Write data for y (shifted by +dy/2) and yv axis
832	ALLOCATE( netcdf_data(mask_size(mid,2)) )
833
834	netcdf_data = ( mask_j_global(mid,:mask_size(mid,2)) + 0.5_wp ) * dy
835
836	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_y_mask(mid,av), &
837	netcdf_data, start = (/ 1 /), &
838	count = (/ mask_size(mid,2) /))
839	CALL netcdf_handle_error( 'netcdf_define_header', 501 )
840
841	netcdf_data = mask_j_global(mid,:mask_size(mid,2)) * dy
842
843	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_yv_mask(mid,av), &
844	netcdf_data, start = (/ 1 /), &
845	count = (/ mask_size(mid,2) /))
846	CALL netcdf_handle_error( 'netcdf_define_header', 502 )
847
848	DEALLOCATE( netcdf_data )
849
850	!
851	!-- Write zu and zw data (vertical axes)
852	ALLOCATE( netcdf_data(mask_size(mid,3)) )
853
854	netcdf_data = zu( mask_k_global(mid,:mask_size(mid,3)) )
855
856	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_zu_mask(mid,av), &
857	netcdf_data, start = (/ 1 /), &
858	count = (/ mask_size(mid,3) /) )
859	CALL netcdf_handle_error( 'netcdf_define_header', 503 )
860
861	netcdf_data = zw( mask_k_global(mid,:mask_size(mid,3)) )
862
863	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_zw_mask(mid,av), &
864	netcdf_data, start = (/ 1 /), &
865	count = (/ mask_size(mid,3) /) )
866	CALL netcdf_handle_error( 'netcdf_define_header', 504 )
867
868	DEALLOCATE( netcdf_data )
869
870	!
871	!-- In case of non-flat topography write height information
872	IF ( TRIM( topography ) /= 'flat' ) THEN
873
874	ALLOCATE( netcdf_data_2d(mask_size(mid,1),mask_size(mid,2)) )
875	netcdf_data_2d = zu_s_inner( mask_i_global(mid,:mask_size(mid,1)),&
876	mask_j_global(mid,:mask_size(mid,2)) )
877
878	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), &
879	id_var_zusi_mask(mid,av), &
880	netcdf_data_2d, &
881	start = (/ 1, 1 /), &
882	count = (/ mask_size(mid,1), &
883	mask_size(mid,2) /) )
884	CALL netcdf_handle_error( 'netcdf_define_header', 505 )
885
886	netcdf_data_2d = zw_w_inner( mask_i_global(mid,:mask_size(mid,1)),&
887	mask_j_global(mid,:mask_size(mid,2)) )
888
889	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), &
890	id_var_zwwi_mask(mid,av), &
891	netcdf_data_2d, &
892	start = (/ 1, 1 /), &
893	count = (/ mask_size(mid,1), &
894	mask_size(mid,2) /) )
895	CALL netcdf_handle_error( 'netcdf_define_header', 506 )
896
897	DEALLOCATE( netcdf_data_2d )
898
899	ENDIF
900
901	IF ( land_surface ) THEN
902	!
903	!-- Write zs data (vertical axes for soil model), use negative values
904	!-- to indicate soil depth
905	ALLOCATE( netcdf_data(mask_size(mid,3)) )
906
907	netcdf_data = zs( mask_k_global(mid,:mask_size(mid,3)) )
908
909	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_zs_mask(mid,av), &
910	netcdf_data, start = (/ 1 /), &
911	count = (/ mask_size(mid,3) /) )
912	CALL netcdf_handle_error( 'netcdf_define_header', 538 )
913
914	DEALLOCATE( netcdf_data )
915
916	ENDIF
917
918	!
919	!-- restore original parameter file_id (=formal parameter av) into av
920	av = file_id
921
922
923	CASE ( 'ma_ext' )
924
925	!
926	!-- decompose actual parameter file_id (=formal parameter av) into
927	!-- mid and av
928	file_id = av
929	IF ( file_id <= 200+max_masks ) THEN
930	mid = file_id - 200
931	av = 0
932	ELSE
933	mid = file_id - (200+max_masks)
934	av = 1
935	ENDIF
936
937	!
938	!-- Get the list of variables and compare with the actual run.
939	!-- First var_list_old has to be reset, since GET_ATT does not assign
940	!-- trailing blanks.
941	var_list_old = ' '
942	nc_stat = NF90_GET_ATT( id_set_mask(mid,av), NF90_GLOBAL, 'VAR_LIST',&
943	var_list_old )
944	CALL netcdf_handle_error( 'netcdf_define_header', 507 )
945
946	var_list = ';'
947	i = 1
948	DO WHILE ( domask(mid,av,i)(1:1) /= ' ' )
949	var_list = TRIM(var_list) // TRIM( domask(mid,av,i) ) // ';'
950	i = i + 1
951	ENDDO
952
953	IF ( av == 0 ) THEN
954	var = '(mask)'
955	ELSE
956	var = '(mask_av)'
957	ENDIF
958
959	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
960	WRITE ( message_string, * ) 'netCDF file for ', TRIM( var ), &
961	' data for mask', mid, ' from previous run found,', &
962	'&but this file cannot be extended due to variable ', &
963	'mismatch.&New file is created instead.'
964	CALL message( 'define_netcdf_header', 'PA0335', 0, 1, 0, 6, 0 )
965	extend = .FALSE.
966	RETURN
967	ENDIF
968
969	!
970	!-- Get and compare the number of vertical gridpoints
971	nc_stat = NF90_INQ_VARID( id_set_mask(mid,av), 'zu_3d', &
972	id_var_zu_mask(mid,av) )
973	CALL netcdf_handle_error( 'netcdf_define_header', 508 )
974
975	nc_stat = NF90_INQUIRE_VARIABLE( id_set_mask(mid,av), &
976	id_var_zu_mask(mid,av), &
977	dimids = id_dim_zu_mask_old )
978	CALL netcdf_handle_error( 'netcdf_define_header', 509 )
979	id_dim_zu_mask(mid,av) = id_dim_zu_mask_old(1)
980
981	nc_stat = NF90_INQUIRE_DIMENSION( id_set_mask(mid,av), &
982	id_dim_zu_mask(mid,av), &
983	len = nz_old )
984	CALL netcdf_handle_error( 'netcdf_define_header', 510 )
985
986	IF ( mask_size(mid,3) /= nz_old ) THEN
987	WRITE ( message_string, * ) 'netCDF file for ', TRIM( var ), &
988	' data for mask', mid, ' from previous run found,', &
989	'&but this file cannot be extended due to mismatch in ', &
990	' number of&vertical grid points.', &
991	'&New file is created instead.'
992	CALL message( 'define_netcdf_header', 'PA0336', 0, 1, 0, 6, 0 )
993	extend = .FALSE.
994	RETURN
995	ENDIF
996
997	!
998	!-- Get the id of the time coordinate (unlimited coordinate) and its
999	!-- last index on the file. The next time level is plmask..count+1.
1000	!-- The current time must be larger than the last output time
1001	!-- on the file.
1002	nc_stat = NF90_INQ_VARID( id_set_mask(mid,av), 'time', &
1003	id_var_time_mask(mid,av) )
1004	CALL netcdf_handle_error( 'netcdf_define_header', 511 )
1005
1006	nc_stat = NF90_INQUIRE_VARIABLE( id_set_mask(mid,av), &
1007	id_var_time_mask(mid,av), &
1008	dimids = id_dim_time_old )
1009	CALL netcdf_handle_error( 'netcdf_define_header', 512 )
1010	id_dim_time_mask(mid,av) = id_dim_time_old(1)
1011
1012	nc_stat = NF90_INQUIRE_DIMENSION( id_set_mask(mid,av), &
1013	id_dim_time_mask(mid,av), &
1014	len = domask_time_count(mid,av) )
1015	CALL netcdf_handle_error( 'netcdf_define_header', 513 )
1016
1017	nc_stat = NF90_GET_VAR( id_set_mask(mid,av), &
1018	id_var_time_mask(mid,av), &
1019	last_time_coordinate, &
1020	start = (/ domask_time_count(mid,av) /), &
1021	count = (/ 1 /) )
1022	CALL netcdf_handle_error( 'netcdf_define_header', 514 )
1023
1024	IF ( last_time_coordinate(1) >= simulated_time ) THEN
1025	WRITE ( message_string, * ) 'netCDF file for ', TRIM( var ), &
1026	' data for mask', mid, ' from previous run found,', &
1027	'&but this file cannot be extended because the current ', &
1028	'output time&is less or equal than the last output time ', &
1029	'on this file.&New file is created instead.'
1030	CALL message( 'define_netcdf_header', 'PA0337', 0, 1, 0, 6, 0 )
1031	domask_time_count(mid,av) = 0
1032	extend = .FALSE.
1033	RETURN
1034	ENDIF
1035
1036	!
1037	!-- Dataset seems to be extendable.
1038	!-- Now get the variable ids.
1039	i = 1
1040	DO WHILE ( domask(mid,av,i)(1:1) /= ' ' )
1041	nc_stat = NF90_INQ_VARID( id_set_mask(mid,av), &
1042	TRIM( domask(mid,av,i) ), &
1043	id_var_domask(mid,av,i) )
1044	CALL netcdf_handle_error( 'netcdf_define_header', 515 )
1045	i = i + 1
1046	ENDDO
1047
1048	!
1049	!-- Update the title attribute on file
1050	!-- In order to avoid 'data mode' errors if updated attributes are larger
1051	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
1052	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
1053	!-- performance loss due to data copying; an alternative strategy would be
1054	!-- to ensure equal attribute size in a job chain. Maybe revise later.
1055	IF ( av == 0 ) THEN
1056	time_average_text = ' '
1057	ELSE
1058	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
1059	averaging_interval
1060	ENDIF
1061	nc_stat = NF90_REDEF( id_set_mask(mid,av) )
1062	CALL netcdf_handle_error( 'netcdf_define_header', 516 )
1063	nc_stat = NF90_PUT_ATT( id_set_mask(mid,av), NF90_GLOBAL, 'title', &
1064	TRIM( run_description_header ) // &
1065	TRIM( time_average_text ) )
1066	CALL netcdf_handle_error( 'netcdf_define_header', 517 )
1067	nc_stat = NF90_ENDDEF( id_set_mask(mid,av) )
1068	CALL netcdf_handle_error( 'netcdf_define_header', 518 )
1069	WRITE ( message_string, * ) 'netCDF file for ', TRIM( var ), &
1070	' data for mask', mid, ' from previous run found.', &
1071	'&This file will be extended.'
1072	CALL message( 'define_netcdf_header', 'PA0338', 0, 0, 0, 6, 0 )
1073	!
1074	!-- restore original parameter file_id (=formal parameter av) into av
1075	av = file_id
1076
1077
1078	CASE ( '3d_new' )
1079
1080	!
1081	!-- Define some global attributes of the dataset
1082	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'Conventions', &
1083	'COARDS' )
1084	CALL netcdf_handle_error( 'netcdf_define_header', 62 )
1085
1086	IF ( av == 0 ) THEN
1087	time_average_text = ' '
1088	ELSE
1089	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
1090	averaging_interval
1091	ENDIF
1092	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'title', &
1093	TRIM( run_description_header ) // &
1094	TRIM( time_average_text ) )
1095	CALL netcdf_handle_error( 'netcdf_define_header', 63 )
1096	IF ( av == 1 ) THEN
1097	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
1098	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'time_avg', &
1099	TRIM( time_average_text ) )
1100	CALL netcdf_handle_error( 'netcdf_define_header', 63 )
1101	ENDIF
1102
1103	!
1104	!-- Define time coordinate for volume data.
1105	!-- For parallel output the time dimensions has to be limited, otherwise
1106	!-- the performance drops significantly.
1107	IF ( netcdf_data_format < 5 ) THEN
1108	CALL netcdf_create_dim( id_set_3d(av), 'time', NF90_UNLIMITED, &
1109	id_dim_time_3d(av), 64 )
1110	ELSE
1111	CALL netcdf_create_dim( id_set_3d(av), 'time', ntdim_3d(av), &
1112	id_dim_time_3d(av), 523 )
1113	ENDIF
1114
1115	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_time_3d(av) /), &
1116	'time', NF90_DOUBLE, id_var_time_3d(av), &
1117	'seconds', '', 65, 66, 00 )
1118	!
1119	!-- Define spatial dimensions and coordinates:
1120	!-- Define vertical coordinate grid (zu grid)
1121	CALL netcdf_create_dim( id_set_3d(av), 'zu_3d', nz_do3d-nzb+1, &
1122	id_dim_zu_3d(av), 67 )
1123	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_zu_3d(av) /), &
1124	'zu_3d', NF90_DOUBLE, id_var_zu_3d(av), &
1125	'meters', '', 68, 69, 00 )
1126	!
1127	!-- Define vertical coordinate grid (zw grid)
1128	CALL netcdf_create_dim( id_set_3d(av), 'zw_3d', nz_do3d-nzb+1, &
1129	id_dim_zw_3d(av), 70 )
1130	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_zw_3d(av) /), &
1131	'zw_3d', NF90_DOUBLE, id_var_zw_3d(av), &
1132	'meters', '', 71, 72, 00 )
1133	!
1134	!-- Define x-axis (for scalar position)
1135	CALL netcdf_create_dim( id_set_3d(av), 'x', nx+2, id_dim_x_3d(av), &
1136	73 )
1137	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_x_3d(av) /), 'x', &
1138	NF90_DOUBLE, id_var_x_3d(av), 'meters', '', &
1139	74, 75, 00 )
1140	!
1141	!-- Define x-axis (for u position)
1142	CALL netcdf_create_dim( id_set_3d(av), 'xu', nx+2, id_dim_xu_3d(av), &
1143	358 )
1144	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_xu_3d(av) /), 'xu', &
1145	NF90_DOUBLE, id_var_xu_3d(av), 'meters', '', &
1146	359, 360, 000 )
1147	!
1148	!-- Define y-axis (for scalar position)
1149	CALL netcdf_create_dim( id_set_3d(av), 'y', ny+2, id_dim_y_3d(av), &
1150	76 )
1151	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_y_3d(av) /), 'y', &
1152	NF90_DOUBLE, id_var_y_3d(av), 'meters', '', &
1153	77, 78, 00 )
1154	!
1155	!-- Define y-axis (for v position)
1156	CALL netcdf_create_dim( id_set_3d(av), 'yv', ny+2, id_dim_yv_3d(av), &
1157	361 )
1158	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_yv_3d(av) /), 'yv', &
1159	NF90_DOUBLE, id_var_yv_3d(av), 'meters', '', &
1160	362, 363, 000 )
1161	!
1162	!-- In case of non-flat topography define 2d-arrays containing the height
1163	!-- information
1164	IF ( TRIM( topography ) /= 'flat' ) THEN
1165	!
1166	!-- Define zusi = zu(nzb_s_inner)
1167	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_x_3d(av), &
1168	id_dim_y_3d(av) /), 'zusi', NF90_DOUBLE, &
1169	id_var_zusi_3d(av), 'meters', &
1170	'zu(nzb_s_inner)', 413, 414, 415 )
1171	!
1172	!-- Define zwwi = zw(nzb_w_inner)
1173	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_x_3d(av), &
1174	id_dim_y_3d(av) /), 'zwwi', NF90_DOUBLE, &
1175	id_var_zwwi_3d(av), 'meters', &
1176	'zw(nzb_w_inner)', 416, 417, 418 )
1177
1178	ENDIF
1179
1180	IF ( land_surface ) THEN
1181	!
1182	!-- Define vertical coordinate grid (zs grid)
1183	CALL netcdf_create_dim( id_set_3d(av), 'zs_3d', &
1184	nzt_soil-nzb_soil+1, id_dim_zs_3d(av), 70 )
1185	CALL netcdf_create_var( id_set_3d(av), (/ id_dim_zs_3d(av) /), &
1186	'zs_3d', NF90_DOUBLE, id_var_zs_3d(av), &
1187	'meters', '', 71, 72, 00 )
1188
1189	ENDIF
1190
1191	!
1192	!-- Define the variables
1193	var_list = ';'
1194	i = 1
1195
1196	DO WHILE ( do3d(av,i)(1:1) /= ' ' )
1197
1198	!
1199	!-- Check for the grid
1200	found = .TRUE.
1201	SELECT CASE ( do3d(av,i) )
1202	!
1203	!-- Most variables are defined on the scalar grid
1204	CASE ( 'e', 'lpt', 'nr', 'p', 'pc', 'pr', 'prr', 'pt', 'q', &
1205	'qc', 'ql', 'ql_c', 'ql_v', 'ql_vp', 'qr', 'qv', 'rho', &
1206	's', 'sa', 'vpt' , 'rad_lw_cs_hr', 'rad_lw_hr', &
1207	'rad_sw_cs_hr', 'rad_sw_hr' )
1208
1209	grid_x = 'x'
1210	grid_y = 'y'
1211	grid_z = 'zu'
1212	!
1213	!-- u grid
1214	CASE ( 'u' )
1215
1216	grid_x = 'xu'
1217	grid_y = 'y'
1218	grid_z = 'zu'
1219	!
1220	!-- v grid
1221	CASE ( 'v' )
1222
1223	grid_x = 'x'
1224	grid_y = 'yv'
1225	grid_z = 'zu'
1226	!
1227	!-- w grid
1228	CASE ( 'rad_lw_in', 'rad_lw_out', 'rad_sw_in', 'rad_sw_out', &
1229	'w' )
1230
1231	grid_x = 'x'
1232	grid_y = 'y'
1233	grid_z = 'zw'
1234	!
1235	!-- soil grid
1236	CASE ( 'm_soil', 't_soil' )
1237
1238	grid_x = 'x'
1239	grid_y = 'y'
1240	grid_z = 'zs'
1241
1242	CASE DEFAULT
1243	!
1244	!-- Check for user-defined quantities
1245	CALL user_define_netcdf_grid( do3d(av,i), found, grid_x, &
1246	grid_y, grid_z )
1247
1248	IF ( .NOT. found ) THEN
1249	WRITE ( message_string, * ) 'no grid defined for varia', &
1250	'ble ', TRIM( do3d(av,i) )
1251	CALL message( 'define_netcdf_header', 'PA0244', 0, 1, 0, &
1252	6, 0 )
1253	ENDIF
1254
1255	END SELECT
1256
1257	!
1258	!-- Select the respective dimension ids
1259	IF ( grid_x == 'x' ) THEN
1260	id_x = id_dim_x_3d(av)
1261	ELSEIF ( grid_x == 'xu' ) THEN
1262	id_x = id_dim_xu_3d(av)
1263	ENDIF
1264
1265	IF ( grid_y == 'y' ) THEN
1266	id_y = id_dim_y_3d(av)
1267	ELSEIF ( grid_y == 'yv' ) THEN
1268	id_y = id_dim_yv_3d(av)
1269	ENDIF
1270
1271	IF ( grid_z == 'zu' ) THEN
1272	id_z = id_dim_zu_3d(av)
1273	ELSEIF ( grid_z == 'zw' ) THEN
1274	id_z = id_dim_zw_3d(av)
1275	ELSEIF ( grid_z == 'zs' ) THEN
1276	id_z = id_dim_zs_3d(av)
1277	ENDIF
1278
1279	!
1280	!-- Define the grid
1281	CALL netcdf_create_var( id_set_3d(av),(/ id_x, id_y, id_z, &
1282	id_dim_time_3d(av) /), do3d(av,i), &
1283	nc_precision(4), id_var_do3d(av,i), &
1284	TRIM( do3d_unit(av,i) ), do3d(av,i), 79, &
1285	80, 357 )
1286	#if defined( __netcdf4_parallel )
1287	IF ( netcdf_data_format > 4 ) THEN
1288	!
1289	!-- Set no fill for every variable to increase performance.
1290	nc_stat = NF90_DEF_VAR_FILL( id_set_3d(av), &
1291	id_var_do3d(av,i), &
1292	1, 0 )
1293	CALL netcdf_handle_error( 'netcdf_define_header', 532 )
1294	!
1295	!-- Set collective io operations for parallel io
1296	nc_stat = NF90_VAR_PAR_ACCESS( id_set_3d(av), &
1297	id_var_do3d(av,i), &
1298	NF90_COLLECTIVE )
1299	CALL netcdf_handle_error( 'netcdf_define_header', 445 )
1300	ENDIF
1301	#endif
1302	var_list = TRIM( var_list ) // TRIM( do3d(av,i) ) // ';'
1303
1304	i = i + 1
1305
1306	ENDDO
1307
1308	!
1309	!-- No arrays to output
1310	IF ( i == 1 ) RETURN
1311
1312	!
1313	!-- Write the list of variables as global attribute (this is used by
1314	!-- restart runs and by combine_plot_fields)
1315	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'VAR_LIST', &
1316	var_list )
1317	CALL netcdf_handle_error( 'netcdf_define_header', 81 )
1318
1319	!
1320	!-- Set general no fill, otherwise the performance drops significantly for
1321	!-- parallel output.
1322	nc_stat = NF90_SET_FILL( id_set_3d(av), NF90_NOFILL, oldmode )
1323	CALL netcdf_handle_error( 'netcdf_define_header', 528 )
1324
1325	!
1326	!-- Leave netCDF define mode
1327	nc_stat = NF90_ENDDEF( id_set_3d(av) )
1328	CALL netcdf_handle_error( 'netcdf_define_header', 82 )
1329
1330	!
1331	!-- These data are only written by PE0 for parallel output to increase
1332	!-- the performance.
1333	IF ( myid == 0 .OR. netcdf_data_format < 5 ) THEN
1334	!
1335	!-- Write data for x (shifted by +dx/2) and xu axis
1336	ALLOCATE( netcdf_data(0:nx+1) )
1337
1338	DO i = 0, nx+1
1339	netcdf_data(i) = ( i + 0.5 ) * dx
1340	ENDDO
1341
1342	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_x_3d(av), &
1343	netcdf_data, start = (/ 1 /), &
1344	count = (/ nx+2 /) )
1345	CALL netcdf_handle_error( 'netcdf_define_header', 83 )
1346
1347	DO i = 0, nx+1
1348	netcdf_data(i) = i * dx
1349	ENDDO
1350
1351	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_xu_3d(av), &
1352	netcdf_data, start = (/ 1 /), &
1353	count = (/ nx+2 /) )
1354	CALL netcdf_handle_error( 'netcdf_define_header', 385 )
1355
1356	DEALLOCATE( netcdf_data )
1357
1358	!
1359	!-- Write data for y (shifted by +dy/2) and yv axis
1360	ALLOCATE( netcdf_data(0:ny+1) )
1361
1362	DO i = 0, ny+1
1363	netcdf_data(i) = ( i + 0.5_wp ) * dy
1364	ENDDO
1365
1366	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_y_3d(av), &
1367	netcdf_data, start = (/ 1 /), &
1368	count = (/ ny+2 /) )
1369	CALL netcdf_handle_error( 'netcdf_define_header', 84 )
1370
1371	DO i = 0, ny+1
1372	netcdf_data(i) = i * dy
1373	ENDDO
1374
1375	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_yv_3d(av), &
1376	netcdf_data, start = (/ 1 /), &
1377	count = (/ ny+2 /))
1378	CALL netcdf_handle_error( 'netcdf_define_header', 387 )
1379
1380	DEALLOCATE( netcdf_data )
1381
1382	!
1383	!-- Write zu and zw data (vertical axes)
1384	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_zu_3d(av), &
1385	zu(nzb:nz_do3d), start = (/ 1 /), &
1386	count = (/ nz_do3d-nzb+1 /) )
1387	CALL netcdf_handle_error( 'netcdf_define_header', 85 )
1388
1389
1390	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_zw_3d(av), &
1391	zw(nzb:nz_do3d), start = (/ 1 /), &
1392	count = (/ nz_do3d-nzb+1 /) )
1393	CALL netcdf_handle_error( 'netcdf_define_header', 86 )
1394
1395	!
1396	!-- In case of non-flat topography write height information
1397	IF ( TRIM( topography ) /= 'flat' ) THEN
1398
1399	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_zusi_3d(av), &
1400	zu_s_inner(0:nx+1,0:ny+1), &
1401	start = (/ 1, 1 /), &
1402	count = (/ nx+2, ny+2 /) )
1403	CALL netcdf_handle_error( 'netcdf_define_header', 419 )
1404
1405	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_zwwi_3d(av), &
1406	zw_w_inner(0:nx+1,0:ny+1), &
1407	start = (/ 1, 1 /), &
1408	count = (/ nx+2, ny+2 /) )
1409	CALL netcdf_handle_error( 'netcdf_define_header', 420 )
1410
1411	ENDIF
1412
1413	IF ( land_surface ) THEN
1414	!
1415	!-- Write zs grid
1416	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_zs_3d(av), &
1417	- zs(nzb_soil:nzt_soil), start = (/ 1 /), &
1418	count = (/ nzt_soil-nzb_soil+1 /) )
1419	CALL netcdf_handle_error( 'netcdf_define_header', 86 )
1420	ENDIF
1421
1422	ENDIF
1423
1424	CASE ( '3d_ext' )
1425
1426	!
1427	!-- Get the list of variables and compare with the actual run.
1428	!-- First var_list_old has to be reset, since GET_ATT does not assign
1429	!-- trailing blanks.
1430	var_list_old = ' '
1431	nc_stat = NF90_GET_ATT( id_set_3d(av), NF90_GLOBAL, 'VAR_LIST', &
1432	var_list_old )
1433	CALL netcdf_handle_error( 'netcdf_define_header', 87 )
1434
1435	var_list = ';'
1436	i = 1
1437	DO WHILE ( do3d(av,i)(1:1) /= ' ' )
1438	var_list = TRIM(var_list) // TRIM( do3d(av,i) ) // ';'
1439	i = i + 1
1440	ENDDO
1441
1442	IF ( av == 0 ) THEN
1443	var = '(3d)'
1444	ELSE
1445	var = '(3d_av)'
1446	ENDIF
1447
1448	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
1449	message_string = 'netCDF file for volume data ' // &
1450	TRIM( var ) // ' from previous run found,' // &
1451	'&but this file cannot be extended due to' // &
1452	' variable mismatch.' // &
1453	'&New file is created instead.'
1454	CALL message( 'define_netcdf_header', 'PA0245', 0, 1, 0, 6, 0 )
1455	extend = .FALSE.
1456	RETURN
1457	ENDIF
1458
1459	!
1460	!-- Get and compare the number of vertical gridpoints
1461	nc_stat = NF90_INQ_VARID( id_set_3d(av), 'zu_3d', id_var_zu_3d(av) )
1462	CALL netcdf_handle_error( 'netcdf_define_header', 88 )
1463
1464	nc_stat = NF90_INQUIRE_VARIABLE( id_set_3d(av), id_var_zu_3d(av), &
1465	dimids = id_dim_zu_3d_old )
1466	CALL netcdf_handle_error( 'netcdf_define_header', 89 )
1467	id_dim_zu_3d(av) = id_dim_zu_3d_old(1)
1468
1469	nc_stat = NF90_INQUIRE_DIMENSION( id_set_3d(av), id_dim_zu_3d(av), &
1470	len = nz_old )
1471	CALL netcdf_handle_error( 'netcdf_define_header', 90 )
1472
1473	IF ( nz_do3d-nzb+1 /= nz_old ) THEN
1474	message_string = 'netCDF file for volume data ' // &
1475	TRIM( var ) // ' from previous run found,' // &
1476	'&but this file cannot be extended due to' // &
1477	' mismatch in number of' // &
1478	'&vertical grid points (nz_do3d).' // &
1479	'&New file is created instead.'
1480	CALL message( 'define_netcdf_header', 'PA0246', 0, 1, 0, 6, 0 )
1481	extend = .FALSE.
1482	RETURN
1483	ENDIF
1484
1485	!
1486	!-- Get the id of the time coordinate (unlimited coordinate) and its
1487	!-- last index on the file. The next time level is pl3d..count+1.
1488	!-- The current time must be larger than the last output time
1489	!-- on the file.
1490	nc_stat = NF90_INQ_VARID( id_set_3d(av), 'time', id_var_time_3d(av) )
1491	CALL netcdf_handle_error( 'netcdf_define_header', 91 )
1492
1493	nc_stat = NF90_INQUIRE_VARIABLE( id_set_3d(av), id_var_time_3d(av), &
1494	dimids = id_dim_time_old )
1495	CALL netcdf_handle_error( 'netcdf_define_header', 92 )
1496
1497	id_dim_time_3d(av) = id_dim_time_old(1)
1498
1499	nc_stat = NF90_INQUIRE_DIMENSION( id_set_3d(av), id_dim_time_3d(av), &
1500	len = ntime_count )
1501	CALL netcdf_handle_error( 'netcdf_define_header', 93 )
1502
1503	!
1504	!-- For non-parallel output use the last output time level of the netcdf
1505	!-- file because the time dimension is unlimited. In case of parallel
1506	!-- output the variable ntime_count could get the value of 9*10E36 because
1507	!-- the time dimension is limited.
1508	IF ( netcdf_data_format < 5 ) do3d_time_count(av) = ntime_count
1509
1510	nc_stat = NF90_GET_VAR( id_set_3d(av), id_var_time_3d(av), &
1511	last_time_coordinate, &
1512	start = (/ do3d_time_count(av) /), &
1513	count = (/ 1 /) )
1514	CALL netcdf_handle_error( 'netcdf_define_header', 94 )
1515
1516	IF ( last_time_coordinate(1) >= simulated_time ) THEN
1517	message_string = 'netCDF file for volume data ' // &
1518	TRIM( var ) // ' from previous run found,' // &
1519	'&but this file cannot be extended becaus' // &
1520	'e the current output time' // &
1521	'&is less or equal than the last output t' // &
1522	'ime on this file.' // &
1523	'&New file is created instead.'
1524	CALL message( 'define_netcdf_header', 'PA0247', 0, 1, 0, 6, 0 )
1525	do3d_time_count(av) = 0
1526	extend = .FALSE.
1527	RETURN
1528	ENDIF
1529
1530	IF ( netcdf_data_format > 4 ) THEN
1531	!
1532	!-- Set needed time levels (ntdim) to
1533	!-- saved time levels + to be saved time levels.
1534	IF ( av == 0 ) THEN
1535	ntdim_3d(0) = do3d_time_count(0) + CEILING( &
1536	( end_time - MAX( skip_time_do3d, &
1537	simulated_time_at_begin ) &
1538	) / dt_do3d )
1539	ELSE
1540	ntdim_3d(1) = do3d_time_count(1) + CEILING( &
1541	( end_time - MAX( skip_time_data_output_av, &
1542	simulated_time_at_begin ) &
1543	) / dt_data_output_av )
1544	ENDIF
1545
1546	!
1547	!-- Check if the needed number of output time levels is increased
1548	!-- compared to the number of time levels in the existing file.
1549	IF ( ntdim_3d(av) > ntime_count ) THEN
1550	message_string = 'netCDF file for volume data ' // &
1551	TRIM( var ) // ' from previous run found,' // &
1552	'&but this file cannot be extended becaus' // &
1553	'e the number of output time levels&has b' // &
1554	'een increased compared to the previous s' // &
1555	'imulation.' // &
1556	'&New file is created instead.'
1557	CALL message( 'define_netcdf_header', 'PA0388', 0, 1, 0, 6, 0 )
1558	do3d_time_count(av) = 0
1559	extend = .FALSE.
1560	!
1561	!-- Recalculate the needed time levels for the new file.
1562	IF ( av == 0 ) THEN
1563	ntdim_3d(0) = CEILING( &
1564	( end_time - MAX( skip_time_do3d, &
1565	simulated_time_at_begin ) &
1566	) / dt_do3d )
1567	ELSE
1568	ntdim_3d(1) = CEILING( &
1569	( end_time - MAX( skip_time_data_output_av, &
1570	simulated_time_at_begin ) &
1571	) / dt_data_output_av )
1572	ENDIF
1573	RETURN
1574	ENDIF
1575	ENDIF
1576
1577	!
1578	!-- Dataset seems to be extendable.
1579	!-- Now get the variable ids.
1580	i = 1
1581	DO WHILE ( do3d(av,i)(1:1) /= ' ' )
1582	nc_stat = NF90_INQ_VARID( id_set_3d(av), TRIM( do3d(av,i) ), &
1583	id_var_do3d(av,i) )
1584	CALL netcdf_handle_error( 'netcdf_define_header', 95 )
1585	#if defined( __netcdf4_parallel )
1586	!
1587	!-- Set collective io operations for parallel io
1588	IF ( netcdf_data_format > 4 ) THEN
1589	nc_stat = NF90_VAR_PAR_ACCESS( id_set_3d(av), &
1590	id_var_do3d(av,i), &
1591	NF90_COLLECTIVE )
1592	CALL netcdf_handle_error( 'netcdf_define_header', 453 )
1593	ENDIF
1594	#endif
1595	i = i + 1
1596	ENDDO
1597
1598	!
1599	!-- Update the title attribute on file
1600	!-- In order to avoid 'data mode' errors if updated attributes are larger
1601	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
1602	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
1603	!-- performance loss due to data copying; an alternative strategy would be
1604	!-- to ensure equal attribute size. Maybe revise later.
1605	IF ( av == 0 ) THEN
1606	time_average_text = ' '
1607	ELSE
1608	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
1609	averaging_interval
1610	ENDIF
1611	nc_stat = NF90_REDEF( id_set_3d(av) )
1612	CALL netcdf_handle_error( 'netcdf_define_header', 429 )
1613	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'title', &
1614	TRIM( run_description_header ) // &
1615	TRIM( time_average_text ) )
1616	CALL netcdf_handle_error( 'netcdf_define_header', 96 )
1617	nc_stat = NF90_ENDDEF( id_set_3d(av) )
1618	CALL netcdf_handle_error( 'netcdf_define_header', 430 )
1619	message_string = 'netCDF file for volume data ' // &
1620	TRIM( var ) // ' from previous run found.' // &
1621	'&This file will be extended.'
1622	CALL message( 'define_netcdf_header', 'PA0248', 0, 0, 0, 6, 0 )
1623
1624	CASE ( 'xy_new' )
1625
1626	!
1627	!-- Define some global attributes of the dataset
1628	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'Conventions', &
1629	'COARDS' )
1630	CALL netcdf_handle_error( 'netcdf_define_header', 97 )
1631
1632	IF ( av == 0 ) THEN
1633	time_average_text = ' '
1634	ELSE
1635	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
1636	averaging_interval
1637	ENDIF
1638	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'title', &
1639	TRIM( run_description_header ) // &
1640	TRIM( time_average_text ) )
1641	CALL netcdf_handle_error( 'netcdf_define_header', 98 )
1642	IF ( av == 1 ) THEN
1643	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
1644	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'time_avg', &
1645	TRIM( time_average_text ) )
1646	CALL netcdf_handle_error( 'netcdf_define_header', 98 )
1647	ENDIF
1648
1649	!
1650	!-- Define time coordinate for xy sections.
1651	!-- For parallel output the time dimensions has to be limited, otherwise
1652	!-- the performance drops significantly.
1653	IF ( netcdf_data_format < 5 ) THEN
1654	CALL netcdf_create_dim( id_set_xy(av), 'time', NF90_UNLIMITED, &
1655	id_dim_time_xy(av), 99 )
1656	ELSE
1657	CALL netcdf_create_dim( id_set_xy(av), 'time', ntdim_2d_xy(av), &
1658	id_dim_time_xy(av), 524 )
1659	ENDIF
1660
1661	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_time_xy(av) /), &
1662	'time', NF90_DOUBLE, id_var_time_xy(av), &
1663	'seconds', '', 100, 101, 000 )
1664	!
1665	!-- Define the spatial dimensions and coordinates for xy-sections.
1666	!-- First, determine the number of horizontal sections to be written.
1667	IF ( section(1,1) == -9999 ) THEN
1668	RETURN
1669	ELSE
1670	ns = 1
1671	DO WHILE ( section(ns,1) /= -9999 .AND. ns <= 100 )
1672	ns = ns + 1
1673	ENDDO
1674	ns = ns - 1
1675	ENDIF
1676
1677	!
1678	!-- Define vertical coordinate grid (zu grid)
1679	CALL netcdf_create_dim( id_set_xy(av), 'zu_xy', ns, &
1680	id_dim_zu_xy(av), 102 )
1681	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_zu_xy(av) /), &
1682	'zu_xy', NF90_DOUBLE, id_var_zu_xy(av), &
1683	'meters', '', 103, 104, 000 )
1684	!
1685	!-- Define vertical coordinate grid (zw grid)
1686	CALL netcdf_create_dim( id_set_xy(av), 'zw_xy', ns, &
1687	id_dim_zw_xy(av), 105 )
1688	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_zw_xy(av) /), &
1689	'zw_xy', NF90_DOUBLE, id_var_zw_xy(av), &
1690	'meters', '', 106, 107, 000 )
1691
1692	IF ( land_surface ) THEN
1693
1694	ns_do = 0
1695	DO WHILE ( section(ns_do+1,1) < nzs )
1696	ns_do = ns_do + 1
1697	ENDDO
1698	!
1699	!-- Define vertical coordinate grid (zs grid)
1700	CALL netcdf_create_dim( id_set_xy(av), 'zs_xy', ns_do, &
1701	id_dim_zs_xy(av), 539 )
1702	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_zs_xy(av) /), &
1703	'zs_xy', NF90_DOUBLE, id_var_zs_xy(av), &
1704	'meters', '', 540, 541, 000 )
1705
1706	ENDIF
1707
1708	!
1709	!-- Define a pseudo vertical coordinate grid for the surface variables
1710	!-- u* and t* to store their height level
1711	CALL netcdf_create_dim( id_set_xy(av), 'zu1_xy', 1, &
1712	id_dim_zu1_xy(av), 108 )
1713	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_zu1_xy(av) /), &
1714	'zu1_xy', NF90_DOUBLE, id_var_zu1_xy(av), &
1715	'meters', '', 109, 110, 000 )
1716	!
1717	!-- Define a variable to store the layer indices of the horizontal cross
1718	!-- sections, too
1719	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_zu_xy(av) /), &
1720	'ind_z_xy', NF90_DOUBLE, &
1721	id_var_ind_z_xy(av), 'gridpoints', '', 111, &
1722	112, 000 )
1723	!
1724	!-- Define x-axis (for scalar position)
1725	CALL netcdf_create_dim( id_set_xy(av), 'x', nx+2, id_dim_x_xy(av), &
1726	113 )
1727	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_x_xy(av) /), 'x', &
1728	NF90_DOUBLE, id_var_x_xy(av), 'meters', '', &
1729	114, 115, 000 )
1730	!
1731	!-- Define x-axis (for u position)
1732	CALL netcdf_create_dim( id_set_xy(av), 'xu', nx+2, &
1733	id_dim_xu_xy(av), 388 )
1734	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_xu_xy(av) /), 'xu', &
1735	NF90_DOUBLE, id_var_xu_xy(av), 'meters', '', &
1736	389, 390, 000 )
1737	!
1738	!-- Define y-axis (for scalar position)
1739	CALL netcdf_create_dim( id_set_xy(av), 'y', ny+2, id_dim_y_xy(av), &
1740	116 )
1741	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_y_xy(av) /), 'y', &
1742	NF90_DOUBLE, id_var_y_xy(av), 'meters', '', &
1743	117, 118, 000 )
1744	!
1745	!-- Define y-axis (for scalar position)
1746	CALL netcdf_create_dim( id_set_xy(av), 'yv', ny+2, &
1747	id_dim_yv_xy(av), 364 )
1748	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_yv_xy(av) /), 'yv', &
1749	NF90_DOUBLE, id_var_yv_xy(av), 'meters', '', &
1750	365, 366, 000 )
1751	!
1752	!-- In case of non-flat topography define 2d-arrays containing the height
1753	!-- information
1754	IF ( TRIM( topography ) /= 'flat' ) THEN
1755	!
1756	!-- Define zusi = zu(nzb_s_inner)
1757	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_x_xy(av), &
1758	id_dim_y_xy(av) /), 'zusi', NF90_DOUBLE, &
1759	id_var_zusi_xy(av), 'meters', &
1760	'zu(nzb_s_inner)', 421, 422, 423 )
1761	!
1762	!-- Define zwwi = zw(nzb_w_inner)
1763	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_x_xy(av), &
1764	id_dim_y_xy(av) /), 'zwwi', NF90_DOUBLE, &
1765	id_var_zwwi_xy(av), 'meters', &
1766	'zw(nzb_w_inner)', 424, 425, 426 )
1767
1768	ENDIF
1769
1770	!
1771	!-- Define the variables
1772	var_list = ';'
1773	i = 1
1774
1775	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
1776
1777	IF ( INDEX( do2d(av,i), 'xy' ) /= 0 ) THEN
1778	!
1779	!-- If there is a star in the variable name (u* or t*), it is a
1780	!-- surface variable. Define it with id_dim_zu1_xy.
1781	IF ( INDEX( do2d(av,i), '*' ) /= 0 ) THEN
1782
1783	CALL netcdf_create_var( id_set_xy(av), (/ id_dim_x_xy(av), &
1784	id_dim_y_xy(av), id_dim_zu1_xy(av), &
1785	id_dim_time_xy(av) /), do2d(av,i), &
1786	nc_precision(1), id_var_do2d(av,i), &
1787	TRIM( do2d_unit(av,i) ), &
1788	do2d(av,i), 119, 120, 354 )
1789
1790	ELSE
1791
1792	!
1793	!-- Check for the grid
1794	found = .TRUE.
1795	SELECT CASE ( do2d(av,i) )
1796	!
1797	!-- Most variables are defined on the zu grid
1798	CASE ( 'e_xy', 'lpt_xy', 'nr_xy', 'p_xy', 'pc_xy', &
1799	'pr_xy', 'prr_xy', 'pt_xy', 'q_xy', 'qc_xy', &
1800	'ql_xy', 'ql_c_xy', 'ql_v_xy', 'ql_vp_xy', &
1801	'qr_xy', 'qv_xy', 'rad_lw_cs_hr_xy', &
1802	'rad_lw_hr_xy', 'rad_sw_cs_hr_xy', 'rad_sw_hr_xy',&
1803	'rho_xy', 's_xy', 'sa_xy', 'vpt_xy' )
1804
1805	grid_x = 'x'
1806	grid_y = 'y'
1807	grid_z = 'zu'
1808	!
1809	!-- u grid
1810	CASE ( 'u_xy' )
1811
1812	grid_x = 'xu'
1813	grid_y = 'y'
1814	grid_z = 'zu'
1815	!
1816	!-- v grid
1817	CASE ( 'v_xy' )
1818
1819	grid_x = 'x'
1820	grid_y = 'yv'
1821	grid_z = 'zu'
1822	!
1823	!-- w grid
1824	CASE ( 'rad_lw_in_xy', 'rad_lw_out_xy', 'rad_sw_in_xy', &
1825	'rad_sw_out_xy' , 'w_xy' )
1826
1827	grid_x = 'x'
1828	grid_y = 'y'
1829	grid_z = 'zw'
1830	!
1831	!-- soil grid
1832	CASE ( 'm_soil_xy', 't_soil_xy' )
1833	grid_x = 'x'
1834	grid_y = 'y'
1835	grid_z = 'zs'
1836
1837	CASE DEFAULT
1838	!
1839	!-- Check for user-defined quantities
1840	CALL user_define_netcdf_grid( do2d(av,i), found, &
1841	grid_x, grid_y, grid_z )
1842
1843	IF ( .NOT. found ) THEN
1844	WRITE ( message_string, * ) 'no grid defined for', &
1845	' variable ', TRIM( do2d(av,i) )
1846	CALL message( 'define_netcdf_header', 'PA0244', &
1847	0, 1, 0, 6, 0 )
1848	ENDIF
1849
1850	END SELECT
1851
1852	!
1853	!-- Select the respective dimension ids
1854	IF ( grid_x == 'x' ) THEN
1855	id_x = id_dim_x_xy(av)
1856	ELSEIF ( grid_x == 'xu' ) THEN
1857	id_x = id_dim_xu_xy(av)
1858	ENDIF
1859
1860	IF ( grid_y == 'y' ) THEN
1861	id_y = id_dim_y_xy(av)
1862	ELSEIF ( grid_y == 'yv' ) THEN
1863	id_y = id_dim_yv_xy(av)
1864	ENDIF
1865
1866	IF ( grid_z == 'zu' ) THEN
1867	id_z = id_dim_zu_xy(av)
1868	ELSEIF ( grid_z == 'zw' ) THEN
1869	id_z = id_dim_zw_xy(av)
1870	ELSEIF ( grid_z == 'zs' ) THEN
1871	id_z = id_dim_zs_xy(av)
1872	ENDIF
1873
1874	!
1875	!-- Define the grid
1876	CALL netcdf_create_var( id_set_xy(av), (/ id_x, id_y, id_z, &
1877	id_dim_time_xy(av) /), do2d(av,i), &
1878	nc_precision(1), id_var_do2d(av,i), &
1879	TRIM( do2d_unit(av,i) ), &
1880	do2d(av,i), 119, 120, 354 )
1881
1882	ENDIF
1883
1884	#if defined( __netcdf4_parallel )
1885	IF ( netcdf_data_format > 4 ) THEN
1886	!
1887	!-- Set no fill for every variable to increase performance.
1888	nc_stat = NF90_DEF_VAR_FILL( id_set_xy(av), &
1889	id_var_do2d(av,i), &
1890	1, 0 )
1891	CALL netcdf_handle_error( 'netcdf_define_header', 533 )
1892	!
1893	!-- Set collective io operations for parallel io
1894	nc_stat = NF90_VAR_PAR_ACCESS( id_set_xy(av), &
1895	id_var_do2d(av,i), &
1896	NF90_COLLECTIVE )
1897	CALL netcdf_handle_error( 'netcdf_define_header', 448 )
1898	ENDIF
1899	#endif
1900	var_list = TRIM( var_list) // TRIM( do2d(av,i) ) // ';'
1901
1902	ENDIF
1903
1904	i = i + 1
1905
1906	ENDDO
1907
1908	!
1909	!-- No arrays to output. Close the netcdf file and return.
1910	IF ( i == 1 ) RETURN
1911
1912	!
1913	!-- Write the list of variables as global attribute (this is used by
1914	!-- restart runs and by combine_plot_fields)
1915	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'VAR_LIST', &
1916	var_list )
1917	CALL netcdf_handle_error( 'netcdf_define_header', 121 )
1918
1919	!
1920	!-- Set general no fill, otherwise the performance drops significantly for
1921	!-- parallel output.
1922	nc_stat = NF90_SET_FILL( id_set_xy(av), NF90_NOFILL, oldmode )
1923	CALL netcdf_handle_error( 'netcdf_define_header', 529 )
1924
1925	!
1926	!-- Leave netCDF define mode
1927	nc_stat = NF90_ENDDEF( id_set_xy(av) )
1928	CALL netcdf_handle_error( 'netcdf_define_header', 122 )
1929
1930	!
1931	!-- These data are only written by PE0 for parallel output to increase
1932	!-- the performance.
1933	IF ( myid == 0 .OR. netcdf_data_format < 5 ) THEN
1934
1935	!
1936	!-- Write axis data: z_xy, x, y
1937	ALLOCATE( netcdf_data(1:ns) )
1938
1939	!
1940	!-- Write zu data
1941	DO i = 1, ns
1942	IF( section(i,1) == -1 ) THEN
1943	netcdf_data(i) = -1.0_wp ! section averaged along z
1944	ELSE
1945	netcdf_data(i) = zu( section(i,1) )
1946	ENDIF
1947	ENDDO
1948	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zu_xy(av), &
1949	netcdf_data, start = (/ 1 /), &
1950	count = (/ ns /) )
1951	CALL netcdf_handle_error( 'netcdf_define_header', 123 )
1952
1953	!
1954	!-- Write zw data
1955	DO i = 1, ns
1956	IF( section(i,1) == -1 ) THEN
1957	netcdf_data(i) = -1.0_wp ! section averaged along z
1958	ELSE
1959	netcdf_data(i) = zw( section(i,1) )
1960	ENDIF
1961	ENDDO
1962	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zw_xy(av), &
1963	netcdf_data, start = (/ 1 /), &
1964	count = (/ ns /) )
1965	CALL netcdf_handle_error( 'netcdf_define_header', 124 )
1966
1967	!
1968	!-- Write zs data
1969	IF ( land_surface ) THEN
1970	ns_do = 0
1971	DO i = 1, ns
1972	IF( section(i,1) == -1 ) THEN
1973	netcdf_data(i) = 1.0_wp ! section averaged along z
1974	ns_do = ns_do + 1
1975	ELSEIF ( section(i,1) < nzs ) THEN
1976	netcdf_data(i) = - zs( section(i,1) )
1977	ns_do = ns_do + 1
1978	ENDIF
1979	ENDDO
1980
1981	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zs_xy(av), &
1982	netcdf_data(1:ns_do), start = (/ 1 /), &
1983	count = (/ ns_do /) )
1984	CALL netcdf_handle_error( 'netcdf_define_header', 124 )
1985
1986	ENDIF
1987
1988	!
1989	!-- Write gridpoint number data
1990	netcdf_data(1:ns) = section(1:ns,1)
1991	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_ind_z_xy(av), &
1992	netcdf_data, start = (/ 1 /), &
1993	count = (/ ns /) )
1994	CALL netcdf_handle_error( 'netcdf_define_header', 125 )
1995
1996	DEALLOCATE( netcdf_data )
1997
1998	!
1999	!-- Write the cross section height u, t
2000	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zu1_xy(av), &
2001	(/ zu(nzb+1) /), start = (/ 1 /), &
2002	count = (/ 1 /) )
2003	CALL netcdf_handle_error( 'netcdf_define_header', 126 )
2004
2005	!
2006	!-- Write data for x (shifted by +dx/2) and xu axis
2007	ALLOCATE( netcdf_data(0:nx+1) )
2008
2009	DO i = 0, nx+1
2010	netcdf_data(i) = ( i + 0.5_wp ) * dx
2011	ENDDO
2012
2013	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_x_xy(av), &
2014	netcdf_data, start = (/ 1 /), &
2015	count = (/ nx+2 /) )
2016	CALL netcdf_handle_error( 'netcdf_define_header', 127 )
2017
2018	DO i = 0, nx+1
2019	netcdf_data(i) = i * dx
2020	ENDDO
2021
2022	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_xu_xy(av), &
2023	netcdf_data, start = (/ 1 /), &
2024	count = (/ nx+2 /) )
2025	CALL netcdf_handle_error( 'netcdf_define_header', 367 )
2026
2027	DEALLOCATE( netcdf_data )
2028
2029	!
2030	!-- Write data for y (shifted by +dy/2) and yv axis
2031	ALLOCATE( netcdf_data(0:ny+1) )
2032
2033	DO i = 0, ny+1
2034	netcdf_data(i) = ( i + 0.5_wp ) * dy
2035	ENDDO
2036
2037	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_y_xy(av), &
2038	netcdf_data, start = (/ 1 /), &
2039	count = (/ ny+2 /))
2040	CALL netcdf_handle_error( 'netcdf_define_header', 128 )
2041
2042	DO i = 0, ny+1
2043	netcdf_data(i) = i * dy
2044	ENDDO
2045
2046	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_yv_xy(av), &
2047	netcdf_data, start = (/ 1 /), &
2048	count = (/ ny+2 /))
2049	CALL netcdf_handle_error( 'netcdf_define_header', 368 )
2050
2051	DEALLOCATE( netcdf_data )
2052
2053	!
2054	!-- In case of non-flat topography write height information
2055	IF ( TRIM( topography ) /= 'flat' ) THEN
2056
2057	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zusi_xy(av), &
2058	zu_s_inner(0:nx+1,0:ny+1), &
2059	start = (/ 1, 1 /), &
2060	count = (/ nx+2, ny+2 /) )
2061	CALL netcdf_handle_error( 'netcdf_define_header', 427 )
2062
2063	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zwwi_xy(av), &
2064	zw_w_inner(0:nx+1,0:ny+1), &
2065	start = (/ 1, 1 /), &
2066	count = (/ nx+2, ny+2 /) )
2067	CALL netcdf_handle_error( 'netcdf_define_header', 428 )
2068
2069	ENDIF
2070
2071
2072
2073	ENDIF
2074
2075	CASE ( 'xy_ext' )
2076
2077	!
2078	!-- Get the list of variables and compare with the actual run.
2079	!-- First var_list_old has to be reset, since GET_ATT does not assign
2080	!-- trailing blanks.
2081	var_list_old = ' '
2082	nc_stat = NF90_GET_ATT( id_set_xy(av), NF90_GLOBAL, 'VAR_LIST', &
2083	var_list_old )
2084	CALL netcdf_handle_error( 'netcdf_define_header', 129 )
2085
2086	var_list = ';'
2087	i = 1
2088	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
2089	IF ( INDEX( do2d(av,i), 'xy' ) /= 0 ) THEN
2090	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
2091	ENDIF
2092	i = i + 1
2093	ENDDO
2094
2095	IF ( av == 0 ) THEN
2096	var = '(xy)'
2097	ELSE
2098	var = '(xy_av)'
2099	ENDIF
2100
2101	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
2102	message_string = 'netCDF file for cross-sections ' // &
2103	TRIM( var ) // ' from previous run found,' // &
2104	'& but this file cannot be extended due to' // &
2105	' variable mismatch.' // &
2106	'&New file is created instead.'
2107	CALL message( 'define_netcdf_header', 'PA0249', 0, 1, 0, 6, 0 )
2108	extend = .FALSE.
2109	RETURN
2110	ENDIF
2111
2112	!
2113	!-- Calculate the number of current sections
2114	ns = 1
2115	DO WHILE ( section(ns,1) /= -9999 .AND. ns <= 100 )
2116	ns = ns + 1
2117	ENDDO
2118	ns = ns - 1
2119
2120	!
2121	!-- Get and compare the number of horizontal cross sections
2122	nc_stat = NF90_INQ_VARID( id_set_xy(av), 'zu_xy', id_var_zu_xy(av) )
2123	CALL netcdf_handle_error( 'netcdf_define_header', 130 )
2124
2125	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xy(av), id_var_zu_xy(av), &
2126	dimids = id_dim_zu_xy_old )
2127	CALL netcdf_handle_error( 'netcdf_define_header', 131 )
2128	id_dim_zu_xy(av) = id_dim_zu_xy_old(1)
2129
2130	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xy(av), id_dim_zu_xy(av), &
2131	len = ns_old )
2132	CALL netcdf_handle_error( 'netcdf_define_header', 132 )
2133
2134	IF ( ns /= ns_old ) THEN
2135	message_string = 'netCDF file for cross-sections ' // &
2136	TRIM( var ) // ' from previous run found,' // &
2137	'&but this file cannot be extended due to' // &
2138	' mismatch in number of' // &
2139	'&cross sections.' // &
2140	'&New file is created instead.'
2141	CALL message( 'define_netcdf_header', 'PA0250', 0, 1, 0, 6, 0 )
2142	extend = .FALSE.
2143	RETURN
2144	ENDIF
2145
2146	!
2147	!-- Get and compare the heights of the cross sections
2148	ALLOCATE( netcdf_data(1:ns_old) )
2149
2150	nc_stat = NF90_GET_VAR( id_set_xy(av), id_var_zu_xy(av), netcdf_data )
2151	CALL netcdf_handle_error( 'netcdf_define_header', 133 )
2152
2153	DO i = 1, ns
2154	IF ( section(i,1) /= -1 ) THEN
2155	IF ( zu(section(i,1)) /= netcdf_data(i) ) THEN
2156	message_string = 'netCDF file for cross-sections ' // &
2157	TRIM( var ) // ' from previous run found,' // &
2158	'&but this file cannot be extended' // &
2159	' due to mismatch in cross' // &
2160	'&section levels.' // &
2161	'&New file is created instead.'
2162	CALL message( 'define_netcdf_header', 'PA0251', &
2163	0, 1, 0, 6, 0 )
2164	extend = .FALSE.
2165	RETURN
2166	ENDIF
2167	ELSE
2168	IF ( -1.0_wp /= netcdf_data(i) ) THEN
2169	message_string = 'netCDF file for cross-sections ' // &
2170	TRIM( var ) // ' from previous run found,' // &
2171	'&but this file cannot be extended' // &
2172	' due to mismatch in cross' // &
2173	'&section levels.' // &
2174	'&New file is created instead.'
2175	CALL message( 'define_netcdf_header', 'PA0251',&
2176	0, 1, 0, 6, 0 )
2177	extend = .FALSE.
2178	RETURN
2179	ENDIF
2180	ENDIF
2181	ENDDO
2182
2183	DEALLOCATE( netcdf_data )
2184
2185	!
2186	!-- Get the id of the time coordinate (unlimited coordinate) and its
2187	!-- last index on the file. The next time level is do2d..count+1.
2188	!-- The current time must be larger than the last output time
2189	!-- on the file.
2190	nc_stat = NF90_INQ_VARID( id_set_xy(av), 'time', id_var_time_xy(av) )
2191	CALL netcdf_handle_error( 'netcdf_define_header', 134 )
2192
2193	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xy(av), id_var_time_xy(av), &
2194	dimids = id_dim_time_old )
2195	CALL netcdf_handle_error( 'netcdf_define_header', 135 )
2196	id_dim_time_xy(av) = id_dim_time_old(1)
2197
2198	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xy(av), id_dim_time_xy(av), &
2199	len = ntime_count )
2200	CALL netcdf_handle_error( 'netcdf_define_header', 136 )
2201
2202	!
2203	!-- For non-parallel output use the last output time level of the netcdf
2204	!-- file because the time dimension is unlimited. In case of parallel
2205	!-- output the variable ntime_count could get the value of 9*10E36 because
2206	!-- the time dimension is limited.
2207	IF ( netcdf_data_format < 5 ) do2d_xy_time_count(av) = ntime_count
2208
2209	nc_stat = NF90_GET_VAR( id_set_xy(av), id_var_time_xy(av), &
2210	last_time_coordinate, &
2211	start = (/ do2d_xy_time_count(av) /), &
2212	count = (/ 1 /) )
2213	CALL netcdf_handle_error( 'netcdf_define_header', 137 )
2214
2215	IF ( last_time_coordinate(1) >= simulated_time ) THEN
2216	message_string = 'netCDF file for cross sections ' // &
2217	TRIM( var ) // ' from previous run found,' // &
2218	'&but this file cannot be extended becaus' // &
2219	'e the current output time' // &
2220	'&is less or equal than the last output t' // &
2221	'ime on this file.' // &
2222	'&New file is created instead.'
2223	CALL message( 'define_netcdf_header', 'PA0252', 0, 1, 0, 6, 0 )
2224	do2d_xy_time_count(av) = 0
2225	extend = .FALSE.
2226	RETURN
2227	ENDIF
2228
2229	IF ( netcdf_data_format > 4 ) THEN
2230	!
2231	!-- Set needed time levels (ntdim) to
2232	!-- saved time levels + to be saved time levels.
2233	IF ( av == 0 ) THEN
2234	ntdim_2d_xy(0) = do2d_xy_time_count(0) + CEILING( &
2235	( end_time - MAX( skip_time_do2d_xy, &
2236	simulated_time_at_begin ) &
2237	) / dt_do2d_xy )
2238	ELSE
2239	ntdim_2d_xy(1) = do2d_xy_time_count(1) + CEILING( &
2240	( end_time - MAX( skip_time_data_output_av, &
2241	simulated_time_at_begin ) &
2242	) / dt_data_output_av )
2243	ENDIF
2244
2245	!
2246	!-- Check if the needed number of output time levels is increased
2247	!-- compared to the number of time levels in the existing file.
2248	IF ( ntdim_2d_xy(av) > ntime_count ) THEN
2249	message_string = 'netCDF file for cross sections ' // &
2250	TRIM( var ) // ' from previous run found,' // &
2251	'&but this file cannot be extended becaus' // &
2252	'e the number of output time levels&has b' // &
2253	'een increased compared to the previous s' // &
2254	'imulation.' // &
2255	'&New file is created instead.'
2256	CALL message( 'define_netcdf_header', 'PA0389', 0, 1, 0, 6, 0 )
2257	do2d_xy_time_count(av) = 0
2258	extend = .FALSE.
2259	!
2260	!-- Recalculate the needed time levels for the new file.
2261	IF ( av == 0 ) THEN
2262	ntdim_2d_xy(0) = CEILING( &
2263	( end_time - MAX( skip_time_do2d_xy, &
2264	simulated_time_at_begin ) &
2265	) / dt_do2d_xy )
2266	ELSE
2267	ntdim_2d_xy(1) = CEILING( &
2268	( end_time - MAX( skip_time_data_output_av, &
2269	simulated_time_at_begin ) &
2270	) / dt_data_output_av )
2271	ENDIF
2272	RETURN
2273	ENDIF
2274	ENDIF
2275
2276	!
2277	!-- Dataset seems to be extendable.
2278	!-- Now get the variable ids.
2279	i = 1
2280	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
2281	IF ( INDEX( do2d(av,i), 'xy' ) /= 0 ) THEN
2282	nc_stat = NF90_INQ_VARID( id_set_xy(av), do2d(av,i), &
2283	id_var_do2d(av,i) )
2284	CALL netcdf_handle_error( 'netcdf_define_header', 138 )
2285	#if defined( __netcdf4_parallel )
2286	!
2287	!-- Set collective io operations for parallel io
2288	IF ( netcdf_data_format > 4 ) THEN
2289	nc_stat = NF90_VAR_PAR_ACCESS( id_set_xy(av), &
2290	id_var_do2d(av,i), &
2291	NF90_COLLECTIVE )
2292	CALL netcdf_handle_error( 'netcdf_define_header', 454 )
2293	ENDIF
2294	#endif
2295	ENDIF
2296	i = i + 1
2297	ENDDO
2298
2299	!
2300	!-- Update the title attribute on file
2301	!-- In order to avoid 'data mode' errors if updated attributes are larger
2302	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
2303	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
2304	!-- performance loss due to data copying; an alternative strategy would be
2305	!-- to ensure equal attribute size in a job chain. Maybe revise later.
2306	IF ( av == 0 ) THEN
2307	time_average_text = ' '
2308	ELSE
2309	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
2310	averaging_interval
2311	ENDIF
2312	nc_stat = NF90_REDEF( id_set_xy(av) )
2313	CALL netcdf_handle_error( 'netcdf_define_header', 431 )
2314	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'title', &
2315	TRIM( run_description_header ) // &
2316	TRIM( time_average_text ) )
2317	CALL netcdf_handle_error( 'netcdf_define_header', 139 )
2318	nc_stat = NF90_ENDDEF( id_set_xy(av) )
2319	CALL netcdf_handle_error( 'netcdf_define_header', 432 )
2320	message_string = 'netCDF file for cross-sections ' // &
2321	TRIM( var ) // ' from previous run found.' // &
2322	'&This file will be extended.'
2323	CALL message( 'define_netcdf_header', 'PA0253', 0, 0, 0, 6, 0 )
2324
2325
2326	CASE ( 'xz_new' )
2327
2328	!
2329	!-- Define some global attributes of the dataset
2330	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'Conventions', &
2331	'COARDS' )
2332	CALL netcdf_handle_error( 'netcdf_define_header', 140 )
2333
2334	IF ( av == 0 ) THEN
2335	time_average_text = ' '
2336	ELSE
2337	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
2338	averaging_interval
2339	ENDIF
2340	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'title', &
2341	TRIM( run_description_header ) // &
2342	TRIM( time_average_text ) )
2343	CALL netcdf_handle_error( 'netcdf_define_header', 141 )
2344	IF ( av == 1 ) THEN
2345	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
2346	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'time_avg', &
2347	TRIM( time_average_text ) )
2348	CALL netcdf_handle_error( 'netcdf_define_header', 141 )
2349	ENDIF
2350
2351	!
2352	!-- Define time coordinate for xz sections.
2353	!-- For parallel output the time dimensions has to be limited, otherwise
2354	!-- the performance drops significantly.
2355	IF ( netcdf_data_format < 5 ) THEN
2356	CALL netcdf_create_dim( id_set_xz(av), 'time', NF90_UNLIMITED, &
2357	id_dim_time_xz(av), 142 )
2358	ELSE
2359	CALL netcdf_create_dim( id_set_xz(av), 'time', ntdim_2d_xz(av), &
2360	id_dim_time_xz(av), 525 )
2361	ENDIF
2362
2363	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_time_xz(av) /), &
2364	'time', NF90_DOUBLE, id_var_time_xz(av), &
2365	'seconds', '', 143, 144, 000 )
2366	!
2367	!-- Define the spatial dimensions and coordinates for xz-sections.
2368	!-- First, determine the number of vertical sections to be written.
2369	IF ( section(1,2) == -9999 ) THEN
2370	RETURN
2371	ELSE
2372	ns = 1
2373	DO WHILE ( section(ns,2) /= -9999 .AND. ns <= 100 )
2374	ns = ns + 1
2375	ENDDO
2376	ns = ns - 1
2377	ENDIF
2378
2379	!
2380	!-- Define y-axis (for scalar position)
2381	CALL netcdf_create_dim( id_set_xz(av), 'y_xz', ns, id_dim_y_xz(av), &
2382	145 )
2383	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_y_xz(av) /), &
2384	'y_xz', NF90_DOUBLE, id_var_y_xz(av), &
2385	'meters', '', 146, 147, 000 )
2386	!
2387	!-- Define y-axis (for v position)
2388	CALL netcdf_create_dim( id_set_xz(av), 'yv_xz', ns, &
2389	id_dim_yv_xz(av), 369 )
2390	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_yv_xz(av) /), &
2391	'yv_xz', NF90_DOUBLE, id_var_yv_xz(av), &
2392	'meters', '', 370, 371, 000 )
2393	!
2394	!-- Define a variable to store the layer indices of the vertical cross
2395	!-- sections
2396	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_y_xz(av) /), &
2397	'ind_y_xz', NF90_DOUBLE, &
2398	id_var_ind_y_xz(av), 'gridpoints', '', 148, &
2399	149, 000 )
2400	!
2401	!-- Define x-axis (for scalar position)
2402	CALL netcdf_create_dim( id_set_xz(av), 'x', nx+2, id_dim_x_xz(av), &
2403	150 )
2404	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_x_xz(av) /), 'x', &
2405	NF90_DOUBLE, id_var_x_xz(av), 'meters', '', &
2406	151, 152, 000 )
2407	!
2408	!-- Define x-axis (for u position)
2409	CALL netcdf_create_dim( id_set_xz(av), 'xu', nx+2, id_dim_xu_xz(av), &
2410	372 )
2411	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_xu_xz(av) /), 'xu', &
2412	NF90_DOUBLE, id_var_xu_xz(av), 'meters', '', &
2413	373, 374, 000 )
2414	!
2415	!-- Define the three z-axes (zu, zw, and zs)
2416	CALL netcdf_create_dim( id_set_xz(av), 'zu', nz+2, id_dim_zu_xz(av), &
2417	153 )
2418	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_zu_xz(av) /), 'zu', &
2419	NF90_DOUBLE, id_var_zu_xz(av), 'meters', '', &
2420	154, 155, 000 )
2421	CALL netcdf_create_dim( id_set_xz(av), 'zw', nz+2, id_dim_zw_xz(av), &
2422	156 )
2423	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_zw_xz(av) /), 'zw', &
2424	NF90_DOUBLE, id_var_zw_xz(av), 'meters', '', &
2425	157, 158, 000 )
2426
2427	IF ( land_surface ) THEN
2428
2429	CALL netcdf_create_dim( id_set_xz(av), 'zs', nzs, &
2430	id_dim_zs_xz(av), 542 )
2431	CALL netcdf_create_var( id_set_xz(av), (/ id_dim_zs_xz(av) /), &
2432	'zs', NF90_DOUBLE, id_var_zs_xz(av), &
2433	'meters', '', 543, 544, 000 )
2434
2435	ENDIF
2436
2437	!
2438	!-- Define the variables
2439	var_list = ';'
2440	i = 1
2441
2442	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
2443
2444	IF ( INDEX( do2d(av,i), 'xz' ) /= 0 ) THEN
2445
2446	!
2447	!-- Check for the grid
2448	found = .TRUE.
2449	SELECT CASE ( do2d(av,i) )
2450	!
2451	!-- Most variables are defined on the zu grid
2452	CASE ( 'e_xz', 'lpt_xz', 'nr_xz', 'p_xz', 'pc_xz', 'pr_xz', &
2453	'prr_xz', 'pt_xz', 'q_xz', 'qc_xz', 'ql_xz', &
2454	'ql_c_xz', 'ql_v_xz', 'ql_vp_xz', 'qr_xz', 'qv_xz', &
2455	'rad_lw_cs_hr_xz', 'rad_lw_hr_xz', &
2456	'rad_sw_cs_hr_xz', 'rad_sw_hr_xz''rho_xz', 's_xz', &
2457	'sa_xz', 'vpt_xz' )
2458
2459	grid_x = 'x'
2460	grid_y = 'y'
2461	grid_z = 'zu'
2462	!
2463	!-- u grid
2464	CASE ( 'u_xz' )
2465
2466	grid_x = 'xu'
2467	grid_y = 'y'
2468	grid_z = 'zu'
2469	!
2470	!-- v grid
2471	CASE ( 'v_xz' )
2472
2473	grid_x = 'x'
2474	grid_y = 'yv'
2475	grid_z = 'zu'
2476	!
2477	!-- w grid
2478	CASE ( 'rad_lw_in_xz', 'rad_lw_out_xz', 'rad_sw_in_xz', &
2479	'rad_sw_out_xz', 'w_xz' )
2480
2481	grid_x = 'x'
2482	grid_y = 'y'
2483	grid_z = 'zw'
2484
2485	!
2486	!-- soil grid
2487	CASE ( 'm_soil_xz', 't_soil_xz' )
2488
2489	grid_x = 'x'
2490	grid_y = 'y'
2491	grid_z = 'zs'
2492
2493	CASE DEFAULT
2494	!
2495	!-- Check for user-defined quantities
2496	CALL user_define_netcdf_grid( do2d(av,i), found, &
2497	grid_x, grid_y, grid_z )
2498	IF ( .NOT. found ) THEN
2499	WRITE ( message_string, * ) 'no grid defined for', &
2500	' variable ', TRIM( do2d(av,i) )
2501	CALL message( 'define_netcdf_header', 'PA0244', &
2502	0, 1, 0, 6, 0 )
2503	ENDIF
2504
2505	END SELECT
2506
2507	!
2508	!-- Select the respective dimension ids
2509	IF ( grid_x == 'x' ) THEN
2510	id_x = id_dim_x_xz(av)
2511	ELSEIF ( grid_x == 'xu' ) THEN
2512	id_x = id_dim_xu_xz(av)
2513	ENDIF
2514
2515	IF ( grid_y == 'y' ) THEN
2516	id_y = id_dim_y_xz(av)
2517	ELSEIF ( grid_y == 'yv' ) THEN
2518	id_y = id_dim_yv_xz(av)
2519	ENDIF
2520
2521	IF ( grid_z == 'zu' ) THEN
2522	id_z = id_dim_zu_xz(av)
2523	ELSEIF ( grid_z == 'zw' ) THEN
2524	id_z = id_dim_zw_xz(av)
2525	ELSEIF ( grid_z == 'zs' ) THEN
2526	id_z = id_dim_zs_xz(av)
2527	ENDIF
2528
2529	!
2530	!-- Define the grid
2531	CALL netcdf_create_var( id_set_xz(av), (/ id_x, id_y, id_z, &
2532	id_dim_time_xz(av) /), do2d(av,i), &
2533	nc_precision(2), id_var_do2d(av,i), &
2534	TRIM( do2d_unit(av,i) ), do2d(av,i), &
2535	159, 160, 355 )
2536
2537	#if defined( __netcdf4_parallel )
2538
2539	IF ( netcdf_data_format > 4 ) THEN
2540	!
2541	!-- Set no fill for every variable to increase performance.
2542	nc_stat = NF90_DEF_VAR_FILL( id_set_xz(av), &
2543	id_var_do2d(av,i), &
2544	1, 0 )
2545	CALL netcdf_handle_error( 'netcdf_define_header', 534 )
2546	!
2547	!-- Set independent io operations for parallel io. Collective io
2548	!-- is only allowed in case of a 1d-decomposition along x,
2549	!-- because otherwise, not all PEs have output data.
2550	IF ( npey == 1 ) THEN
2551	nc_stat = NF90_VAR_PAR_ACCESS( id_set_xz(av), &
2552	id_var_do2d(av,i), &
2553	NF90_COLLECTIVE )
2554	ELSE
2555	!
2556	!-- Test simulations showed that the output of cross sections
2557	!-- by all PEs in data_output_2d using NF90_COLLECTIVE is
2558	!-- faster than the output by the first row of PEs in
2559	!-- x-direction using NF90_INDEPENDENT.
2560	nc_stat = NF90_VAR_PAR_ACCESS( id_set_xz(av), &
2561	id_var_do2d(av,i), &
2562	NF90_COLLECTIVE )
2563	! nc_stat = NF90_VAR_PAR_ACCESS( id_set_xz(av), &
2564	! id_var_do2d(av,i), &
2565	! NF90_INDEPENDENT )
2566	ENDIF
2567	CALL netcdf_handle_error( 'netcdf_define_header', 449 )
2568	ENDIF
2569	#endif
2570	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
2571
2572	ENDIF
2573
2574	i = i + 1
2575
2576	ENDDO
2577
2578	!
2579	!-- No arrays to output. Close the netcdf file and return.
2580	IF ( i == 1 ) RETURN
2581
2582	!
2583	!-- Write the list of variables as global attribute (this is used by
2584	!-- restart runs and by combine_plot_fields)
2585	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'VAR_LIST', &
2586	var_list )
2587	CALL netcdf_handle_error( 'netcdf_define_header', 161 )
2588
2589	!
2590	!-- Set general no fill, otherwise the performance drops significantly for
2591	!-- parallel output.
2592	nc_stat = NF90_SET_FILL( id_set_xz(av), NF90_NOFILL, oldmode )
2593	CALL netcdf_handle_error( 'netcdf_define_header', 530 )
2594
2595	!
2596	!-- Leave netCDF define mode
2597	nc_stat = NF90_ENDDEF( id_set_xz(av) )
2598	CALL netcdf_handle_error( 'netcdf_define_header', 162 )
2599
2600	!
2601	!-- These data are only written by PE0 for parallel output to increase
2602	!-- the performance.
2603	IF ( myid == 0 .OR. netcdf_data_format < 5 ) THEN
2604
2605	!
2606	!-- Write axis data: y_xz, x, zu, zw
2607	ALLOCATE( netcdf_data(1:ns) )
2608
2609	!
2610	!-- Write y_xz data (shifted by +dy/2)
2611	DO i = 1, ns
2612	IF( section(i,2) == -1 ) THEN
2613	netcdf_data(i) = -1.0_wp ! section averaged along y
2614	ELSE
2615	netcdf_data(i) = ( section(i,2) + 0.5_wp ) * dy
2616	ENDIF
2617	ENDDO
2618	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_y_xz(av), &
2619	netcdf_data, start = (/ 1 /), &
2620	count = (/ ns /) )
2621	CALL netcdf_handle_error( 'netcdf_define_header', 163 )
2622
2623	!
2624	!-- Write yv_xz data
2625	DO i = 1, ns
2626	IF( section(i,2) == -1 ) THEN
2627	netcdf_data(i) = -1.0_wp ! section averaged along y
2628	ELSE
2629	netcdf_data(i) = section(i,2) * dy
2630	ENDIF
2631	ENDDO
2632	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_yv_xz(av), &
2633	netcdf_data, start = (/ 1 /), &
2634	count = (/ ns /) )
2635	CALL netcdf_handle_error( 'netcdf_define_header', 375 )
2636
2637	!
2638	!-- Write gridpoint number data
2639	netcdf_data(1:ns) = section(1:ns,2)
2640	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_ind_y_xz(av), &
2641	netcdf_data, start = (/ 1 /), &
2642	count = (/ ns /) )
2643	CALL netcdf_handle_error( 'netcdf_define_header', 164 )
2644
2645
2646	DEALLOCATE( netcdf_data )
2647
2648	!
2649	!-- Write data for x (shifted by +dx/2) and xu axis
2650	ALLOCATE( netcdf_data(0:nx+1) )
2651
2652	DO i = 0, nx+1
2653	netcdf_data(i) = ( i + 0.5_wp ) * dx
2654	ENDDO
2655
2656	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_x_xz(av), &
2657	netcdf_data, start = (/ 1 /), &
2658	count = (/ nx+2 /) )
2659	CALL netcdf_handle_error( 'netcdf_define_header', 165 )
2660
2661	DO i = 0, nx+1
2662	netcdf_data(i) = i * dx
2663	ENDDO
2664
2665	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_xu_xz(av), &
2666	netcdf_data, start = (/ 1 /), &
2667	count = (/ nx+2 /) )
2668	CALL netcdf_handle_error( 'netcdf_define_header', 377 )
2669
2670	DEALLOCATE( netcdf_data )
2671
2672	!
2673	!-- Write zu and zw data (vertical axes)
2674	ALLOCATE( netcdf_data(0:nz+1) )
2675
2676	netcdf_data(0:nz+1) = zu(nzb:nzt+1)
2677	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_zu_xz(av), &
2678	netcdf_data, start = (/ 1 /), &
2679	count = (/ nz+2 /) )
2680	CALL netcdf_handle_error( 'netcdf_define_header', 166 )
2681
2682	netcdf_data(0:nz+1) = zw(nzb:nzt+1)
2683	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_zw_xz(av), &
2684	netcdf_data, start = (/ 1 /), &
2685	count = (/ nz+2 /) )
2686	CALL netcdf_handle_error( 'netcdf_define_header', 167 )
2687
2688	!
2689	!-- Write zs data
2690	IF ( land_surface ) THEN
2691	netcdf_data(0:nzs-1) = - zs(nzb_soil:nzt_soil)
2692	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_zs_xz(av), &
2693	netcdf_data(0:nzs), start = (/ 1 /), &
2694	count = (/ nzt_soil-nzb_soil+1 /) )
2695	CALL netcdf_handle_error( 'netcdf_define_header', 548 )
2696	ENDIF
2697
2698
2699	DEALLOCATE( netcdf_data )
2700
2701	ENDIF
2702
2703
2704	CASE ( 'xz_ext' )
2705
2706	!
2707	!-- Get the list of variables and compare with the actual run.
2708	!-- First var_list_old has to be reset, since GET_ATT does not assign
2709	!-- trailing blanks.
2710	var_list_old = ' '
2711	nc_stat = NF90_GET_ATT( id_set_xz(av), NF90_GLOBAL, 'VAR_LIST', &
2712	var_list_old )
2713	CALL netcdf_handle_error( 'netcdf_define_header', 168 )
2714
2715	var_list = ';'
2716	i = 1
2717	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
2718	IF ( INDEX( do2d(av,i), 'xz' ) /= 0 ) THEN
2719	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
2720	ENDIF
2721	i = i + 1
2722	ENDDO
2723
2724	IF ( av == 0 ) THEN
2725	var = '(xz)'
2726	ELSE
2727	var = '(xz_av)'
2728	ENDIF
2729
2730	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
2731	message_string = 'netCDF file for cross-sections ' // &
2732	TRIM( var ) // ' from previous run found,' // &
2733	'& but this file cannot be extended due to' // &
2734	' variable mismatch.' // &
2735	'&New file is created instead.'
2736	CALL message( 'define_netcdf_header', 'PA0249', 0, 1, 0, 6, 0 )
2737	extend = .FALSE.
2738	RETURN
2739	ENDIF
2740
2741	!
2742	!-- Calculate the number of current sections
2743	ns = 1
2744	DO WHILE ( section(ns,2) /= -9999 .AND. ns <= 100 )
2745	ns = ns + 1
2746	ENDDO
2747	ns = ns - 1
2748
2749	!
2750	!-- Get and compare the number of vertical cross sections
2751	nc_stat = NF90_INQ_VARID( id_set_xz(av), 'y_xz', id_var_y_xz(av) )
2752	CALL netcdf_handle_error( 'netcdf_define_header', 169 )
2753
2754	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xz(av), id_var_y_xz(av), &
2755	dimids = id_dim_y_xz_old )
2756	CALL netcdf_handle_error( 'netcdf_define_header', 170 )
2757	id_dim_y_xz(av) = id_dim_y_xz_old(1)
2758
2759	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xz(av), id_dim_y_xz(av), &
2760	len = ns_old )
2761	CALL netcdf_handle_error( 'netcdf_define_header', 171 )
2762
2763	IF ( ns /= ns_old ) THEN
2764	message_string = 'netCDF file for cross-sections ' // &
2765	TRIM( var ) // ' from previous run found,' // &
2766	'&but this file cannot be extended due to' // &
2767	' mismatch in number of' // &
2768	'&cross sections.' // &
2769	'&New file is created instead.'
2770	CALL message( 'define_netcdf_header', 'PA0250', 0, 1, 0, 6, 0 )
2771	extend = .FALSE.
2772	RETURN
2773	ENDIF
2774
2775	!
2776	!-- Get and compare the heights of the cross sections
2777	ALLOCATE( netcdf_data(1:ns_old) )
2778
2779	nc_stat = NF90_GET_VAR( id_set_xz(av), id_var_y_xz(av), netcdf_data )
2780	CALL netcdf_handle_error( 'netcdf_define_header', 172 )
2781
2782	DO i = 1, ns
2783	IF ( section(i,2) /= -1 ) THEN
2784	IF ( ( ( section(i,2) + 0.5 ) * dy ) /= netcdf_data(i) ) THEN
2785	message_string = 'netCDF file for cross-sections ' // &
2786	TRIM( var ) // ' from previous run found,' // &
2787	'&but this file cannot be extended' // &
2788	' due to mismatch in cross' // &
2789	'&section levels.' // &
2790	'&New file is created instead.'
2791	CALL message( 'define_netcdf_header', 'PA0251', &
2792	0, 1, 0, 6, 0 )
2793	extend = .FALSE.
2794	RETURN
2795	ENDIF
2796	ELSE
2797	IF ( -1.0_wp /= netcdf_data(i) ) THEN
2798	message_string = 'netCDF file for cross-sections ' // &
2799	TRIM( var ) // ' from previous run found,' // &
2800	'&but this file cannot be extended' // &
2801	' due to mismatch in cross' // &
2802	'&section levels.' // &
2803	'&New file is created instead.'
2804	CALL message( 'define_netcdf_header', 'PA0251', &
2805	0, 1, 0, 6, 0 )
2806	extend = .FALSE.
2807	RETURN
2808	ENDIF
2809	ENDIF
2810	ENDDO
2811
2812	DEALLOCATE( netcdf_data )
2813
2814	!
2815	!-- Get the id of the time coordinate (unlimited coordinate) and its
2816	!-- last index on the file. The next time level is do2d..count+1.
2817	!-- The current time must be larger than the last output time
2818	!-- on the file.
2819	nc_stat = NF90_INQ_VARID( id_set_xz(av), 'time', id_var_time_xz(av) )
2820	CALL netcdf_handle_error( 'netcdf_define_header', 173 )
2821
2822	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xz(av), id_var_time_xz(av), &
2823	dimids = id_dim_time_old )
2824	CALL netcdf_handle_error( 'netcdf_define_header', 174 )
2825	id_dim_time_xz(av) = id_dim_time_old(1)
2826
2827	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xz(av), id_dim_time_xz(av), &
2828	len = ntime_count )
2829	CALL netcdf_handle_error( 'netcdf_define_header', 175 )
2830
2831	!
2832	!-- For non-parallel output use the last output time level of the netcdf
2833	!-- file because the time dimension is unlimited. In case of parallel
2834	!-- output the variable ntime_count could get the value of 9*10E36 because
2835	!-- the time dimension is limited.
2836	IF ( netcdf_data_format < 5 ) do2d_xz_time_count(av) = ntime_count
2837
2838	nc_stat = NF90_GET_VAR( id_set_xz(av), id_var_time_xz(av), &
2839	last_time_coordinate, &
2840	start = (/ do2d_xz_time_count(av) /), &
2841	count = (/ 1 /) )
2842	CALL netcdf_handle_error( 'netcdf_define_header', 176 )
2843
2844	IF ( last_time_coordinate(1) >= simulated_time ) THEN
2845	message_string = 'netCDF file for cross sections ' // &
2846	TRIM( var ) // ' from previous run found,' // &
2847	'&but this file cannot be extended becaus' // &
2848	'e the current output time' // &
2849	'&is less or equal than the last output t' // &
2850	'ime on this file.' // &
2851	'&New file is created instead.'
2852	CALL message( 'define_netcdf_header', 'PA0252', 0, 1, 0, 6, 0 )
2853	do2d_xz_time_count(av) = 0
2854	extend = .FALSE.
2855	RETURN
2856	ENDIF
2857
2858	IF ( netcdf_data_format > 4 ) THEN
2859	!
2860	!-- Set needed time levels (ntdim) to
2861	!-- saved time levels + to be saved time levels.
2862	IF ( av == 0 ) THEN
2863	ntdim_2d_xz(0) = do2d_xz_time_count(0) + CEILING( &
2864	( end_time - MAX( skip_time_do2d_xz, &
2865	simulated_time_at_begin ) &
2866	) / dt_do2d_xz )
2867	ELSE
2868	ntdim_2d_xz(1) = do2d_xz_time_count(1) + CEILING( &
2869	( end_time - MAX( skip_time_data_output_av, &
2870	simulated_time_at_begin ) &
2871	) / dt_data_output_av )
2872	ENDIF
2873
2874	!
2875	!-- Check if the needed number of output time levels is increased
2876	!-- compared to the number of time levels in the existing file.
2877	IF ( ntdim_2d_xz(av) > ntime_count ) THEN
2878	message_string = 'netCDF file for cross sections ' // &
2879	TRIM( var ) // ' from previous run found,' // &
2880	'&but this file cannot be extended becaus' // &
2881	'e the number of output time levels&has b' // &
2882	'een increased compared to the previous s' // &
2883	'imulation.' // &
2884	'&New file is created instead.'
2885	CALL message( 'define_netcdf_header', 'PA0390', 0, 1, 0, 6, 0 )
2886	do2d_xz_time_count(av) = 0
2887	extend = .FALSE.
2888	!
2889	!-- Recalculate the needed time levels for the new file.
2890	IF ( av == 0 ) THEN
2891	ntdim_2d_xz(0) = CEILING( &
2892	( end_time - MAX( skip_time_do2d_xz, &
2893	simulated_time_at_begin ) &
2894	) / dt_do2d_xz )
2895	ELSE
2896	ntdim_2d_xz(1) = CEILING( &
2897	( end_time - MAX( skip_time_data_output_av, &
2898	simulated_time_at_begin ) &
2899	) / dt_data_output_av )
2900	ENDIF
2901	RETURN
2902	ENDIF
2903	ENDIF
2904
2905	!
2906	!-- Dataset seems to be extendable.
2907	!-- Now get the variable ids.
2908	i = 1
2909	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
2910	IF ( INDEX( do2d(av,i), 'xz' ) /= 0 ) THEN
2911	nc_stat = NF90_INQ_VARID( id_set_xz(av), do2d(av,i), &
2912	id_var_do2d(av,i) )
2913	CALL netcdf_handle_error( 'netcdf_define_header', 177 )
2914	#if defined( __netcdf4_parallel )
2915	!
2916	!-- Set independent io operations for parallel io. Collective io
2917	!-- is only allowed in case of a 1d-decomposition along x, because
2918	!-- otherwise, not all PEs have output data.
2919	IF ( netcdf_data_format > 4 ) THEN
2920	IF ( npey == 1 ) THEN
2921	nc_stat = NF90_VAR_PAR_ACCESS( id_set_xz(av), &
2922	id_var_do2d(av,i), &
2923	NF90_COLLECTIVE )
2924	ELSE
2925	!
2926	!-- Test simulations showed that the output of cross sections
2927	!-- by all PEs in data_output_2d using NF90_COLLECTIVE is
2928	!-- faster than the output by the first row of PEs in
2929	!-- x-direction using NF90_INDEPENDENT.
2930	nc_stat = NF90_VAR_PAR_ACCESS( id_set_xz(av), &
2931	id_var_do2d(av,i), &
2932	NF90_COLLECTIVE )
2933	! nc_stat = NF90_VAR_PAR_ACCESS( id_set_xz(av), &
2934	! id_var_do2d(av,i), &
2935	! NF90_INDEPENDENT )
2936	ENDIF
2937	CALL netcdf_handle_error( 'netcdf_define_header', 455 )
2938	ENDIF
2939	#endif
2940	ENDIF
2941	i = i + 1
2942	ENDDO
2943
2944	!
2945	!-- Update the title attribute on file
2946	!-- In order to avoid 'data mode' errors if updated attributes are larger
2947	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
2948	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
2949	!-- performance loss due to data copying; an alternative strategy would be
2950	!-- to ensure equal attribute size in a job chain. Maybe revise later.
2951	IF ( av == 0 ) THEN
2952	time_average_text = ' '
2953	ELSE
2954	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
2955	averaging_interval
2956	ENDIF
2957	nc_stat = NF90_REDEF( id_set_xz(av) )
2958	CALL netcdf_handle_error( 'netcdf_define_header', 433 )
2959	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'title', &
2960	TRIM( run_description_header ) // &
2961	TRIM( time_average_text ) )
2962	CALL netcdf_handle_error( 'netcdf_define_header', 178 )
2963	nc_stat = NF90_ENDDEF( id_set_xz(av) )
2964	CALL netcdf_handle_error( 'netcdf_define_header', 434 )
2965	message_string = 'netCDF file for cross-sections ' // &
2966	TRIM( var ) // ' from previous run found.' // &
2967	'&This file will be extended.'
2968	CALL message( 'define_netcdf_header', 'PA0253', 0, 0, 0, 6, 0 )
2969
2970
2971	CASE ( 'yz_new' )
2972
2973	!
2974	!-- Define some global attributes of the dataset
2975	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'Conventions', &
2976	'COARDS' )
2977	CALL netcdf_handle_error( 'netcdf_define_header', 179 )
2978
2979	IF ( av == 0 ) THEN
2980	time_average_text = ' '
2981	ELSE
2982	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
2983	averaging_interval
2984	ENDIF
2985	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'title', &
2986	TRIM( run_description_header ) // &
2987	TRIM( time_average_text ) )
2988	CALL netcdf_handle_error( 'netcdf_define_header', 180 )
2989	IF ( av == 1 ) THEN
2990	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
2991	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'time_avg', &
2992	TRIM( time_average_text ) )
2993	CALL netcdf_handle_error( 'netcdf_define_header', 180 )
2994	ENDIF
2995
2996	!
2997	!-- Define time coordinate for yz sections.
2998	!-- For parallel output the time dimensions has to be limited, otherwise
2999	!-- the performance drops significantly.
3000	IF ( netcdf_data_format < 5 ) THEN
3001	CALL netcdf_create_dim( id_set_yz(av), 'time', NF90_UNLIMITED, &
3002	id_dim_time_yz(av), 181 )
3003	ELSE
3004	CALL netcdf_create_dim( id_set_yz(av), 'time', ntdim_2d_yz(av), &
3005	id_dim_time_yz(av), 526 )
3006	ENDIF
3007
3008	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_time_yz(av) /), &
3009	'time', NF90_DOUBLE, id_var_time_yz(av), &
3010	'seconds', '', 182, 183, 000 )
3011	!
3012	!-- Define the spatial dimensions and coordinates for yz-sections.
3013	!-- First, determine the number of vertical sections to be written.
3014	IF ( section(1,3) == -9999 ) THEN
3015	RETURN
3016	ELSE
3017	ns = 1
3018	DO WHILE ( section(ns,3) /= -9999 .AND. ns <= 100 )
3019	ns = ns + 1
3020	ENDDO
3021	ns = ns - 1
3022	ENDIF
3023
3024	!
3025	!-- Define x axis (for scalar position)
3026	CALL netcdf_create_dim( id_set_yz(av), 'x_yz', ns, id_dim_x_yz(av), &
3027	184 )
3028	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_x_yz(av) /), &
3029	'x_yz', NF90_DOUBLE, id_var_x_yz(av), &
3030	'meters', '', 185, 186, 000 )
3031	!
3032	!-- Define x axis (for u position)
3033	CALL netcdf_create_dim( id_set_yz(av), 'xu_yz', ns, &
3034	id_dim_xu_yz(av), 377 )
3035	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_xu_yz(av) /), &
3036	'xu_yz', NF90_DOUBLE, id_var_xu_yz(av), &
3037	'meters', '', 378, 379, 000 )
3038	!
3039	!-- Define a variable to store the layer indices of the vertical cross
3040	!-- sections
3041	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_x_yz(av) /), &
3042	'ind_x_yz', NF90_DOUBLE, &
3043	id_var_ind_x_yz(av), 'gridpoints', '', 187, &
3044	188, 000 )
3045	!
3046	!-- Define y-axis (for scalar position)
3047	CALL netcdf_create_dim( id_set_yz(av), 'y', ny+2, id_dim_y_yz(av), &
3048	189 )
3049	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_y_yz(av) /), 'y', &
3050	NF90_DOUBLE, id_var_y_yz(av), 'meters', '', &
3051	190, 191, 000 )
3052	!
3053	!-- Define y-axis (for v position)
3054	CALL netcdf_create_dim( id_set_yz(av), 'yv', ny+2, id_dim_yv_yz(av), &
3055	380 )
3056	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_yv_yz(av) /), 'yv', &
3057	NF90_DOUBLE, id_var_yv_yz(av), 'meters', '', &
3058	381, 382, 000 )
3059	!
3060	!-- Define the two z-axes (zu and zw)
3061	CALL netcdf_create_dim( id_set_yz(av), 'zu', nz+2, id_dim_zu_yz(av), &
3062	192 )
3063	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_zu_yz(av) /), 'zu', &
3064	NF90_DOUBLE, id_var_zu_yz(av), 'meters', '', &
3065	193, 194, 000 )
3066
3067	CALL netcdf_create_dim( id_set_yz(av), 'zw', nz+2, id_dim_zw_yz(av), &
3068	195 )
3069	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_zw_yz(av) /), 'zw', &
3070	NF90_DOUBLE, id_var_zw_yz(av), 'meters', '', &
3071	196, 197, 000 )
3072
3073	IF ( land_surface ) THEN
3074
3075	CALL netcdf_create_dim( id_set_yz(av), 'zs', nzs, &
3076	id_dim_zs_yz(av), 545 )
3077	CALL netcdf_create_var( id_set_yz(av), (/ id_dim_zs_yz(av) /), &
3078	'zs', NF90_DOUBLE, id_var_zs_yz(av), &
3079	'meters', '', 546, 547, 000 )
3080
3081	ENDIF
3082
3083	!
3084	!-- Define the variables
3085	var_list = ';'
3086	i = 1
3087
3088	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
3089
3090	IF ( INDEX( do2d(av,i), 'yz' ) /= 0 ) THEN
3091
3092	!
3093	!-- Check for the grid
3094	found = .TRUE.
3095	SELECT CASE ( do2d(av,i) )
3096	!
3097	!-- Most variables are defined on the zu grid
3098	CASE ( 'e_yz', 'lpt_yz', 'nr_yz', 'p_yz', 'pc_yz', 'pr_yz', &
3099	'prr_yz', 'pt_yz', 'q_yz', 'qc_yz', 'ql_yz', &
3100	'ql_c_yz', 'ql_v_yz', 'ql_vp_yz', 'qr_yz', 'qv_yz', &
3101	'rad_lw_cs_hr_yz', 'rad_lw_hr_yz', &
3102	'rad_sw_cs_hr_yz', 'rad_sw_hr_yz''rho_yz', 's_yz', &
3103	'sa_yz', 'vpt_yz' )
3104
3105	grid_x = 'x'
3106	grid_y = 'y'
3107	grid_z = 'zu'
3108	!
3109	!-- u grid
3110	CASE ( 'u_yz' )
3111
3112	grid_x = 'xu'
3113	grid_y = 'y'
3114	grid_z = 'zu'
3115	!
3116	!-- v grid
3117	CASE ( 'v_yz' )
3118
3119	grid_x = 'x'
3120	grid_y = 'yv'
3121	grid_z = 'zu'
3122	!
3123	!-- w grid
3124	CASE ( 'rad_lw_in_yz', 'rad_lw_out_yz', 'rad_sw_in_yz', &
3125	'rad_sw_out_yz', 'w_yz' )
3126
3127	grid_x = 'x'
3128	grid_y = 'y'
3129	grid_z = 'zw'
3130	!
3131	!-- soil grid
3132	CASE ( 'm_soil_yz', 't_soil_yz' )
3133
3134	grid_x = 'x'
3135	grid_y = 'y'
3136	grid_z = 'zs'
3137
3138	CASE DEFAULT
3139	!
3140	!-- Check for user-defined quantities
3141	CALL user_define_netcdf_grid( do2d(av,i), found, &
3142	grid_x, grid_y, grid_z )
3143
3144	IF ( .NOT. found ) THEN
3145	WRITE ( message_string, * ) 'no grid defined for', &
3146	' variable ', TRIM( do2d(av,i) )
3147	CALL message( 'define_netcdf_header', 'PA0244', &
3148	0, 1, 0, 6, 0 )
3149	ENDIF
3150
3151	END SELECT
3152
3153	!
3154	!-- Select the respective dimension ids
3155	IF ( grid_x == 'x' ) THEN
3156	id_x = id_dim_x_yz(av)
3157	ELSEIF ( grid_x == 'xu' ) THEN
3158	id_x = id_dim_xu_yz(av)
3159	ENDIF
3160
3161	IF ( grid_y == 'y' ) THEN
3162	id_y = id_dim_y_yz(av)
3163	ELSEIF ( grid_y == 'yv' ) THEN
3164	id_y = id_dim_yv_yz(av)
3165	ENDIF
3166
3167	IF ( grid_z == 'zu' ) THEN
3168	id_z = id_dim_zu_yz(av)
3169	ELSEIF ( grid_z == 'zw' ) THEN
3170	id_z = id_dim_zw_yz(av)
3171	ELSEIF ( grid_z == 'zs' ) THEN
3172	id_z = id_dim_zs_yz(av)
3173	ENDIF
3174
3175	!
3176	!-- Define the grid
3177	CALL netcdf_create_var( id_set_yz(av), (/ id_x, id_y, id_z, &
3178	id_dim_time_yz(av) /), do2d(av,i), &
3179	nc_precision(3), id_var_do2d(av,i), &
3180	TRIM( do2d_unit(av,i) ), do2d(av,i), &
3181	198, 199, 356 )
3182
3183	#if defined( __netcdf4_parallel )
3184	IF ( netcdf_data_format > 4 ) THEN
3185	!
3186	!-- Set no fill for every variable to increase performance.
3187	nc_stat = NF90_DEF_VAR_FILL( id_set_yz(av), &
3188	id_var_do2d(av,i), &
3189	1, 0 )
3190	CALL netcdf_handle_error( 'netcdf_define_header', 535 )
3191	!
3192	!-- Set independent io operations for parallel io. Collective io
3193	!-- is only allowed in case of a 1d-decomposition along y,
3194	!-- because otherwise, not all PEs have output data.
3195	IF ( npex == 1 ) THEN
3196	nc_stat = NF90_VAR_PAR_ACCESS( id_set_yz(av), &
3197	id_var_do2d(av,i), &
3198	NF90_COLLECTIVE )
3199	ELSE
3200	!
3201	!-- Test simulations showed that the output of cross sections
3202	!-- by all PEs in data_output_2d using NF90_COLLECTIVE is
3203	!-- faster than the output by the first row of PEs in
3204	!-- y-direction using NF90_INDEPENDENT.
3205	nc_stat = NF90_VAR_PAR_ACCESS( id_set_yz(av), &
3206	id_var_do2d(av,i), &
3207	NF90_COLLECTIVE )
3208	! nc_stat = NF90_VAR_PAR_ACCESS( id_set_yz(av), &
3209	! id_var_do2d(av,i), &
3210	! NF90_INDEPENDENT )
3211	ENDIF
3212	CALL netcdf_handle_error( 'netcdf_define_header', 450 )
3213	ENDIF
3214	#endif
3215	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
3216
3217	ENDIF
3218
3219	i = i + 1
3220
3221	ENDDO
3222
3223	!
3224	!-- No arrays to output. Close the netcdf file and return.
3225	IF ( i == 1 ) RETURN
3226
3227	!
3228	!-- Write the list of variables as global attribute (this is used by
3229	!-- restart runs and by combine_plot_fields)
3230	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'VAR_LIST', &
3231	var_list )
3232	CALL netcdf_handle_error( 'netcdf_define_header', 200 )
3233
3234	!
3235	!-- Set general no fill, otherwise the performance drops significantly for
3236	!-- parallel output.
3237	nc_stat = NF90_SET_FILL( id_set_yz(av), NF90_NOFILL, oldmode )
3238	CALL netcdf_handle_error( 'netcdf_define_header', 531 )
3239
3240	!
3241	!-- Leave netCDF define mode
3242	nc_stat = NF90_ENDDEF( id_set_yz(av) )
3243	CALL netcdf_handle_error( 'netcdf_define_header', 201 )
3244
3245	!
3246	!-- These data are only written by PE0 for parallel output to increase
3247	!-- the performance.
3248	IF ( myid == 0 .OR. netcdf_data_format < 5 ) THEN
3249
3250	!
3251	!-- Write axis data: x_yz, y, zu, zw
3252	ALLOCATE( netcdf_data(1:ns) )
3253
3254	!
3255	!-- Write x_yz data (shifted by +dx/2)
3256	DO i = 1, ns
3257	IF( section(i,3) == -1 ) THEN
3258	netcdf_data(i) = -1.0_wp ! section averaged along x
3259	ELSE
3260	netcdf_data(i) = ( section(i,3) + 0.5_wp ) * dx
3261	ENDIF
3262	ENDDO
3263	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_x_yz(av), &
3264	netcdf_data, start = (/ 1 /), &
3265	count = (/ ns /) )
3266	CALL netcdf_handle_error( 'netcdf_define_header', 202 )
3267
3268	!
3269	!-- Write x_yz data (xu grid)
3270	DO i = 1, ns
3271	IF( section(i,3) == -1 ) THEN
3272	netcdf_data(i) = -1.0_wp ! section averaged along x
3273	ELSE
3274	netcdf_data(i) = section(i,3) * dx
3275	ENDIF
3276	ENDDO
3277	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_xu_yz(av), &
3278	netcdf_data, start = (/ 1 /), &
3279	count = (/ ns /) )
3280	CALL netcdf_handle_error( 'netcdf_define_header', 383 )
3281
3282	!
3283	!-- Write gridpoint number data
3284	netcdf_data(1:ns) = section(1:ns,3)
3285	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_ind_x_yz(av), &
3286	netcdf_data, start = (/ 1 /), &
3287	count = (/ ns /) )
3288	CALL netcdf_handle_error( 'netcdf_define_header', 203 )
3289
3290	DEALLOCATE( netcdf_data )
3291
3292	!
3293	!-- Write data for y (shifted by +dy/2) and yv axis
3294	ALLOCATE( netcdf_data(0:ny+1) )
3295
3296	DO j = 0, ny+1
3297	netcdf_data(j) = ( j + 0.5_wp ) * dy
3298	ENDDO
3299
3300	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_y_yz(av), &
3301	netcdf_data, start = (/ 1 /), &
3302	count = (/ ny+2 /) )
3303	CALL netcdf_handle_error( 'netcdf_define_header', 204 )
3304
3305	DO j = 0, ny+1
3306	netcdf_data(j) = j * dy
3307	ENDDO
3308
3309	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_yv_yz(av), &
3310	netcdf_data, start = (/ 1 /), &
3311	count = (/ ny+2 /) )
3312	CALL netcdf_handle_error( 'netcdf_define_header', 384 )
3313
3314	DEALLOCATE( netcdf_data )
3315
3316	!
3317	!-- Write zu and zw data (vertical axes)
3318	ALLOCATE( netcdf_data(0:nz+1) )
3319
3320	netcdf_data(0:nz+1) = zu(nzb:nzt+1)
3321	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_zu_yz(av), &
3322	netcdf_data, start = (/ 1 /), &
3323	count = (/ nz+2 /) )
3324	CALL netcdf_handle_error( 'netcdf_define_header', 205 )
3325
3326	netcdf_data(0:nz+1) = zw(nzb:nzt+1)
3327	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_zw_yz(av), &
3328	netcdf_data, start = (/ 1 /), &
3329	count = (/ nz+2 /) )
3330	CALL netcdf_handle_error( 'netcdf_define_header', 206 )
3331
3332	DEALLOCATE( netcdf_data )
3333
3334	ENDIF
3335
3336
3337	CASE ( 'yz_ext' )
3338
3339	!
3340	!-- Get the list of variables and compare with the actual run.
3341	!-- First var_list_old has to be reset, since GET_ATT does not assign
3342	!-- trailing blanks.
3343	var_list_old = ' '
3344	nc_stat = NF90_GET_ATT( id_set_yz(av), NF90_GLOBAL, 'VAR_LIST', &
3345	var_list_old )
3346	CALL netcdf_handle_error( 'netcdf_define_header', 207 )
3347
3348	var_list = ';'
3349	i = 1
3350	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
3351	IF ( INDEX( do2d(av,i), 'yz' ) /= 0 ) THEN
3352	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
3353	ENDIF
3354	i = i + 1
3355	ENDDO
3356
3357	IF ( av == 0 ) THEN
3358	var = '(yz)'
3359	ELSE
3360	var = '(yz_av)'
3361	ENDIF
3362
3363	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
3364	message_string = 'netCDF file for cross-sections ' // &
3365	TRIM( var ) // ' from previous run found,' // &
3366	'& but this file cannot be extended due to' // &
3367	' variable mismatch.' // &
3368	'&New file is created instead.'
3369	CALL message( 'define_netcdf_header', 'PA0249', 0, 1, 0, 6, 0 )
3370	extend = .FALSE.
3371	RETURN
3372	ENDIF
3373
3374	!
3375	!-- Calculate the number of current sections
3376	ns = 1
3377	DO WHILE ( section(ns,3) /= -9999 .AND. ns <= 100 )
3378	ns = ns + 1
3379	ENDDO
3380	ns = ns - 1
3381
3382	!
3383	!-- Get and compare the number of vertical cross sections
3384	nc_stat = NF90_INQ_VARID( id_set_yz(av), 'x_yz', id_var_x_yz(av) )
3385	CALL netcdf_handle_error( 'netcdf_define_header', 208 )
3386
3387	nc_stat = NF90_INQUIRE_VARIABLE( id_set_yz(av), id_var_x_yz(av), &
3388	dimids = id_dim_x_yz_old )
3389	CALL netcdf_handle_error( 'netcdf_define_header', 209 )
3390	id_dim_x_yz(av) = id_dim_x_yz_old(1)
3391
3392	nc_stat = NF90_INQUIRE_DIMENSION( id_set_yz(av), id_dim_x_yz(av), &
3393	len = ns_old )
3394	CALL netcdf_handle_error( 'netcdf_define_header', 210 )
3395
3396	IF ( ns /= ns_old ) THEN
3397	message_string = 'netCDF file for cross-sections ' // &
3398	TRIM( var ) // ' from previous run found,' // &
3399	'&but this file cannot be extended due to' // &
3400	' mismatch in number of' // &
3401	'&cross sections.' // &
3402	'&New file is created instead.'
3403	CALL message( 'define_netcdf_header', 'PA0250', 0, 1, 0, 6, 0 )
3404	extend = .FALSE.
3405	RETURN
3406	ENDIF
3407
3408	!
3409	!-- Get and compare the heights of the cross sections
3410	ALLOCATE( netcdf_data(1:ns_old) )
3411
3412	nc_stat = NF90_GET_VAR( id_set_yz(av), id_var_x_yz(av), netcdf_data )
3413	CALL netcdf_handle_error( 'netcdf_define_header', 211 )
3414
3415	DO i = 1, ns
3416	IF ( section(i,3) /= -1 ) THEN
3417	IF ( ( ( section(i,3) + 0.5 ) * dx ) /= netcdf_data(i) ) THEN
3418	message_string = 'netCDF file for cross-sections ' // &
3419	TRIM( var ) // ' from previous run found,' // &
3420	'&but this file cannot be extended' // &
3421	' due to mismatch in cross' // &
3422	'&section levels.' // &
3423	'&New file is created instead.'
3424	CALL message( 'define_netcdf_header', 'PA0251', &
3425	0, 1, 0, 6, 0 )
3426	extend = .FALSE.
3427	RETURN
3428	ENDIF
3429	ELSE
3430	IF ( -1.0_wp /= netcdf_data(i) ) THEN
3431	message_string = 'netCDF file for cross-sections ' // &
3432	TRIM( var ) // ' from previous run found,' // &
3433	'&but this file cannot be extended' // &
3434	' due to mismatch in cross' // &
3435	'&section levels.' // &
3436	'&New file is created instead.'
3437	CALL message( 'define_netcdf_header', 'PA0251', &
3438	0, 1, 0, 6, 0 )
3439	extend = .FALSE.
3440	RETURN
3441	ENDIF
3442	ENDIF
3443	ENDDO
3444
3445	DEALLOCATE( netcdf_data )
3446
3447	!
3448	!-- Get the id of the time coordinate (unlimited coordinate) and its
3449	!-- last index on the file. The next time level is pl2d..count+1.
3450	!-- The current time must be larger than the last output time
3451	!-- on the file.
3452	nc_stat = NF90_INQ_VARID( id_set_yz(av), 'time', id_var_time_yz(av) )
3453	CALL netcdf_handle_error( 'netcdf_define_header', 212 )
3454
3455	nc_stat = NF90_INQUIRE_VARIABLE( id_set_yz(av), id_var_time_yz(av), &
3456	dimids = id_dim_time_old )
3457	CALL netcdf_handle_error( 'netcdf_define_header', 213 )
3458	id_dim_time_yz(av) = id_dim_time_old(1)
3459
3460	nc_stat = NF90_INQUIRE_DIMENSION( id_set_yz(av), id_dim_time_yz(av), &
3461	len = ntime_count )
3462	CALL netcdf_handle_error( 'netcdf_define_header', 214 )
3463
3464	!
3465	!-- For non-parallel output use the last output time level of the netcdf
3466	!-- file because the time dimension is unlimited. In case of parallel
3467	!-- output the variable ntime_count could get the value of 9*10E36 because
3468	!-- the time dimension is limited.
3469	IF ( netcdf_data_format < 5 ) do2d_yz_time_count(av) = ntime_count
3470
3471	nc_stat = NF90_GET_VAR( id_set_yz(av), id_var_time_yz(av), &
3472	last_time_coordinate, &
3473	start = (/ do2d_yz_time_count(av) /), &
3474	count = (/ 1 /) )
3475	CALL netcdf_handle_error( 'netcdf_define_header', 215 )
3476
3477	IF ( last_time_coordinate(1) >= simulated_time ) THEN
3478	message_string = 'netCDF file for cross sections ' // &
3479	TRIM( var ) // ' from previous run found,' // &
3480	'&but this file cannot be extended becaus' // &
3481	'e the current output time' // &
3482	'&is less or equal than the last output t' // &
3483	'ime on this file.' // &
3484	'&New file is created instead.'
3485	CALL message( 'define_netcdf_header', 'PA0252', 0, 1, 0, 6, 0 )
3486	do2d_yz_time_count(av) = 0
3487	extend = .FALSE.
3488	RETURN
3489	ENDIF
3490
3491	IF ( netcdf_data_format > 4 ) THEN
3492	!
3493	!-- Set needed time levels (ntdim) to
3494	!-- saved time levels + to be saved time levels.
3495	IF ( av == 0 ) THEN
3496	ntdim_2d_yz(0) = do2d_yz_time_count(0) + CEILING( &
3497	( end_time - MAX( skip_time_do2d_yz, &
3498	simulated_time_at_begin ) &
3499	) / dt_do2d_yz )
3500	ELSE
3501	ntdim_2d_yz(1) = do2d_yz_time_count(1) + CEILING( &
3502	( end_time - MAX( skip_time_data_output_av, &
3503	simulated_time_at_begin ) &
3504	) / dt_data_output_av )
3505	ENDIF
3506
3507	!
3508	!-- Check if the needed number of output time levels is increased
3509	!-- compared to the number of time levels in the existing file.
3510	IF ( ntdim_2d_yz(av) > ntime_count ) THEN
3511	message_string = 'netCDF file for cross sections ' // &
3512	TRIM( var ) // ' from previous run found,' // &
3513	'&but this file cannot be extended becaus' // &
3514	'e the number of output time levels&has b' // &
3515	'een increased compared to the previous s' // &
3516	'imulation.' // &
3517	'&New file is created instead.'
3518	CALL message( 'define_netcdf_header', 'PA0391', 0, 1, 0, 6, 0 )
3519	do2d_yz_time_count(av) = 0
3520	extend = .FALSE.
3521	!
3522	!-- Recalculate the needed time levels for the new file.
3523	IF ( av == 0 ) THEN
3524	ntdim_2d_yz(0) = CEILING( &
3525	( end_time - MAX( skip_time_do2d_yz, &
3526	simulated_time_at_begin ) &
3527	) / dt_do2d_yz )
3528	ELSE
3529	ntdim_2d_yz(1) = CEILING( &
3530	( end_time - MAX( skip_time_data_output_av, &
3531	simulated_time_at_begin ) &
3532	) / dt_data_output_av )
3533	ENDIF
3534	RETURN
3535	ENDIF
3536	ENDIF
3537
3538	!
3539	!-- Dataset seems to be extendable.
3540	!-- Now get the variable ids.
3541	i = 1
3542	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
3543	IF ( INDEX( do2d(av,i), 'yz' ) /= 0 ) THEN
3544	nc_stat = NF90_INQ_VARID( id_set_yz(av), do2d(av,i), &
3545	id_var_do2d(av,i) )
3546	CALL netcdf_handle_error( 'netcdf_define_header', 216 )
3547	#if defined( __netcdf4_parallel )
3548	!
3549	!-- Set independent io operations for parallel io. Collective io
3550	!-- is only allowed in case of a 1d-decomposition along y, because
3551	!-- otherwise, not all PEs have output data.
3552	IF ( netcdf_data_format > 4 ) THEN
3553	IF ( npex == 1 ) THEN
3554	nc_stat = NF90_VAR_PAR_ACCESS( id_set_yz(av), &
3555	id_var_do2d(av,i), &
3556	NF90_COLLECTIVE )
3557	ELSE
3558	!
3559	!-- Test simulations showed that the output of cross sections
3560	!-- by all PEs in data_output_2d using NF90_COLLECTIVE is
3561	!-- faster than the output by the first row of PEs in
3562	!-- y-direction using NF90_INDEPENDENT.
3563	nc_stat = NF90_VAR_PAR_ACCESS( id_set_yz(av), &
3564	id_var_do2d(av,i), &
3565	NF90_COLLECTIVE )
3566	! nc_stat = NF90_VAR_PAR_ACCESS( id_set_yz(av), &
3567	! id_var_do2d(av,i), &
3568	! NF90_INDEPENDENT )
3569	ENDIF
3570	CALL netcdf_handle_error( 'netcdf_define_header', 450 )
3571	ENDIF
3572	#endif
3573	ENDIF
3574	i = i + 1
3575	ENDDO
3576
3577	!
3578	!-- Update the title attribute on file
3579	!-- In order to avoid 'data mode' errors if updated attributes are larger
3580	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
3581	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
3582	!-- performance loss due to data copying; an alternative strategy would be
3583	!-- to ensure equal attribute size in a job chain. Maybe revise later.
3584	IF ( av == 0 ) THEN
3585	time_average_text = ' '
3586	ELSE
3587	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
3588	averaging_interval
3589	ENDIF
3590	nc_stat = NF90_REDEF( id_set_yz(av) )
3591	CALL netcdf_handle_error( 'netcdf_define_header', 435 )
3592	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'title', &
3593	TRIM( run_description_header ) // &
3594	TRIM( time_average_text ) )
3595	CALL netcdf_handle_error( 'netcdf_define_header', 217 )
3596	nc_stat = NF90_ENDDEF( id_set_yz(av) )
3597	CALL netcdf_handle_error( 'netcdf_define_header', 436 )
3598	message_string = 'netCDF file for cross-sections ' // &
3599	TRIM( var ) // ' from previous run found.' // &
3600	'&This file will be extended.'
3601	CALL message( 'define_netcdf_header', 'PA0253', 0, 0, 0, 6, 0 )
3602
3603
3604	CASE ( 'pr_new' )
3605
3606	!
3607	!-- Define some global attributes of the dataset
3608	IF ( averaging_interval_pr /= 0.0_wp ) THEN
3609	WRITE (time_average_text,'('', '',F7.1,'' s average'')') &
3610	averaging_interval_pr
3611	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'title', &
3612	TRIM( run_description_header ) // &
3613	TRIM( time_average_text ) )
3614	CALL netcdf_handle_error( 'netcdf_define_header', 218 )
3615
3616	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval_pr
3617	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'time_avg', &
3618	TRIM( time_average_text ) )
3619	ELSE
3620	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'title', &
3621	TRIM( run_description_header ) )
3622	ENDIF
3623	CALL netcdf_handle_error( 'netcdf_define_header', 219 )
3624
3625	!
3626	!-- Write number of columns and rows of coordinate systems to be plotted
3627	!-- on one page to the netcdf header.
3628	!-- This information can be used by palmplot.
3629	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, &
3630	'no_rows', &
3631	profile_rows )
3632	CALL netcdf_handle_error( 'netcdf_define_header', 519 )
3633
3634	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, &
3635	'no_columns', &
3636	profile_columns )
3637	CALL netcdf_handle_error( 'netcdf_define_header', 520 )
3638
3639
3640	cross_profiles_adj = ADJUSTL( cross_profiles )
3641	cross_profiles_numb = 999999
3642	cross_profiles_char = ''
3643
3644	!
3645	!-- Each profile defined in cross_profiles is written to an array
3646	!-- (cross_profiles_char). The number of the respective coordinate
3647	!-- system is assigned in a second array (cross_profiles_numb).
3648	k = 1
3649
3650	DO i = 1, crmax
3651
3652	IF ( TRIM( cross_profiles_adj(i) ) == ' ' ) EXIT
3653	delim_old = 0
3654
3655	DO j = 1, crmax
3656	delim = INDEX( cross_profiles_adj(i)(delim_old+1:), ' ' )
3657	IF ( delim == 1 ) EXIT
3658	kk = MIN( crmax, k )
3659	cross_profiles_char(kk) = cross_profiles_adj(i)(delim_old+1: &
3660	delim_old+delim-1)
3661	cross_profiles_numb(kk) = i
3662	k = k + 1
3663	cross_profiles_maxi = i
3664	delim_old = delim_old + delim
3665	ENDDO
3666
3667	ENDDO
3668
3669	cross_profiles_count = MIN( crmax, k-1 )
3670	!
3671	!-- Check if all profiles defined in cross_profiles are defined in
3672	!-- data_output_pr. If not, they will be skipped.
3673	DO i = 1, cross_profiles_count
3674	DO j = 1, dopr_n
3675
3676	IF ( TRIM(cross_profiles_char(i)) == TRIM(data_output_pr(j)) ) &
3677	THEN
3678	EXIT
3679	ENDIF
3680
3681	IF ( j == dopr_n ) THEN
3682	cross_profiles_numb(i) = 999999
3683	ENDIF
3684
3685	ENDDO
3686	ENDDO
3687
3688	DO i = 1, crmax
3689	IF ( cross_profiles_numb(i) == 999999 ) THEN
3690	DO j = i + 1, crmax
3691	IF ( cross_profiles_numb(j) /= 999999 ) THEN
3692	cross_profiles_char(i) = cross_profiles_char(j)
3693	cross_profiles_numb(i) = cross_profiles_numb(j)
3694	cross_profiles_numb(j) = 999999
3695	EXIT
3696	ENDIF
3697	ENDDO
3698	ENDIF
3699	ENDDO
3700
3701	DO i = 1, crmax-1
3702	IF ( cross_profiles_numb(i + 1) == 999999 ) THEN
3703	cross_profiles_count = i
3704	EXIT
3705	ENDIF
3706	ENDDO
3707	!
3708	!-- Check if all profiles defined in data_output_pr are defined in
3709	!-- cross_profiles. If not, they will be added to cross_profiles.
3710	DO i = 1, dopr_n
3711	DO j = 1, cross_profiles_count
3712
3713	IF ( TRIM(cross_profiles_char(j)) == TRIM(data_output_pr(i))) &
3714	THEN
3715	EXIT
3716	ENDIF
3717
3718	IF (( j == cross_profiles_count ) .AND. &
3719	( cross_profiles_count <= crmax - 1)) THEN
3720	cross_profiles_count = cross_profiles_count + 1
3721	cross_profiles_maxi = cross_profiles_maxi + 1
3722	cross_profiles_char(MIN( crmax, cross_profiles_count )) = &
3723	TRIM( data_output_pr(i) )
3724	cross_profiles_numb(MIN( crmax, cross_profiles_count )) = &
3725	cross_profiles_maxi
3726	ENDIF
3727
3728	ENDDO
3729	ENDDO
3730
3731	IF ( cross_profiles_count >= crmax ) THEN
3732	message_string = 'It is not allowed to arrange more than ' &
3733	// '100 profiles with & cross_profiles. Apart'&
3734	// ' from that, all profiles are saved & to ' &
3735	// 'the netCDF file.'
3736	CALL message( 'define_netcdf_header', 'PA0354', 0, 0, 0, 6, 0 )
3737	ENDIF
3738
3739	!
3740	!-- Writing cross_profiles to netcdf header. This information can be
3741	!-- used by palmplot. Each profile is separated by ",", each cross is
3742	!-- separated by ";".
3743	char_cross_profiles = ';'
3744	id_last = 1
3745	cross_profiles_count = MIN( cross_profiles_count, crmax )
3746
3747	DO i = 1, cross_profiles_count
3748
3749	IF ( cross_profiles_numb(i) /= 999999 ) THEN
3750	IF ( TRIM( char_cross_profiles ) == ';' ) THEN
3751	char_cross_profiles = TRIM( char_cross_profiles ) // &
3752	TRIM( cross_profiles_char(i) )
3753	ELSEIF ( id_last == cross_profiles_numb(i) ) THEN
3754	char_cross_profiles = TRIM( char_cross_profiles ) // &
3755	',' // TRIM( cross_profiles_char(i) )
3756	ELSE
3757	char_cross_profiles = TRIM( char_cross_profiles ) // &
3758	';' // TRIM( cross_profiles_char(i) )
3759	ENDIF
3760	id_last = cross_profiles_numb(i)
3761	ENDIF
3762
3763	ENDDO
3764
3765	char_cross_profiles = TRIM( char_cross_profiles ) // ';'
3766
3767	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'cross_profiles', &
3768	TRIM( char_cross_profiles ) )
3769	CALL netcdf_handle_error( 'netcdf_define_header', 521 )
3770
3771	!
3772	!-- Define time coordinate for profiles (unlimited dimension)
3773	CALL netcdf_create_dim( id_set_pr, 'time', NF90_UNLIMITED, &
3774	id_dim_time_pr, 220 )
3775	CALL netcdf_create_var( id_set_pr, (/ id_dim_time_pr /), 'time', &
3776	NF90_DOUBLE, id_var_time_pr, 'seconds', '', &
3777	221, 222, 000 )
3778	!
3779	!-- Define the variables
3780	var_list = ';'
3781	DO i = 1, dopr_n
3782
3783	IF ( statistic_regions == 0 ) THEN
3784
3785	!
3786	!-- Define the z-axes (each variable gets its own z-axis)
3787	CALL netcdf_create_dim( id_set_pr, &
3788	'z' // TRIM( data_output_pr(i) ), &
3789	nzt+2-nzb, id_dim_z_pr(i,0), 223 )
3790	CALL netcdf_create_var( id_set_pr, (/ id_dim_z_pr(i,0) /), &
3791	'z' // TRIM( data_output_pr(i) ), &
3792	NF90_DOUBLE, id_var_z_pr(i,0), &
3793	'meters', '', 224, 225, 000 )
3794	!
3795	!-- Define the variable
3796	CALL netcdf_create_var( id_set_pr, (/ id_dim_z_pr(i,0), &
3797	id_dim_time_pr /), data_output_pr(i), &
3798	nc_precision(5), id_var_dopr(i,0), &
3799	TRIM( dopr_unit(i) ), &
3800	TRIM( data_output_pr(i) ), 226, 227, &
3801	228 )
3802
3803	var_list = TRIM( var_list ) // TRIM( data_output_pr(i) ) // ';'
3804
3805	ELSE
3806	!
3807	!-- If statistic regions are defined, add suffix _SR+#SR to the
3808	!-- names
3809	DO j = 0, statistic_regions
3810	WRITE ( suffix, '(''_'',I1)' ) j
3811
3812	!
3813	!-- Define the z-axes (each variable gets it own z-axis)
3814	CALL netcdf_create_dim( id_set_pr, 'z' // &
3815	TRIM(data_output_pr(i)) // suffix, &
3816	nzt+2-nzb, id_dim_z_pr(i,j), 229 )
3817	CALL netcdf_create_var( id_set_pr, (/ id_dim_z_pr(i,j) /), &
3818	'z' // TRIM(data_output_pr(i)) // &
3819	suffix, NF90_DOUBLE, &
3820	id_var_z_pr(i,j), 'meters', '', &
3821	230, 231, 000 )
3822	!
3823	!-- Define the variable
3824	CALL netcdf_create_var( id_set_pr, (/ id_dim_z_pr(i,j), &
3825	id_dim_time_pr /), &
3826	TRIM(data_output_pr(i)) // suffix, &
3827	nc_precision(5), id_var_dopr(i,j), &
3828	TRIM( dopr_unit(i) ), &
3829	TRIM( data_output_pr(i) ) // &
3830	' SR ', 232, 233, 234 )
3831
3832	var_list = TRIM( var_list ) // TRIM( data_output_pr(i) ) // &
3833	suffix // ';'
3834
3835	ENDDO
3836
3837	ENDIF
3838
3839	ENDDO
3840
3841	!
3842	!-- Write the list of variables as global attribute (this is used by
3843	!-- restart runs)
3844	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'VAR_LIST', var_list )
3845	CALL netcdf_handle_error( 'netcdf_define_header', 235 )
3846
3847	!
3848	!-- Define normalization variables (as time series)
3849	DO i = 1, dopr_norm_num
3850
3851	CALL netcdf_create_var( id_set_pr, (/ id_dim_time_pr /), &
3852	'NORM_' // TRIM( dopr_norm_names(i) ), &
3853	nc_precision(5), id_var_norm_dopr(i), &
3854	'', TRIM( dopr_norm_longnames(i) ), 236, &
3855	000, 237 )
3856
3857	ENDDO
3858
3859	!
3860	!-- Leave netCDF define mode
3861	nc_stat = NF90_ENDDEF( id_set_pr )
3862	CALL netcdf_handle_error( 'netcdf_define_header', 238 )
3863
3864	!
3865	!-- Write z-axes data
3866	DO i = 1, dopr_n
3867	DO j = 0, statistic_regions
3868
3869	nc_stat = NF90_PUT_VAR( id_set_pr, id_var_z_pr(i,j), &
3870	hom(nzb:nzt+1,2,dopr_index(i),0), &
3871	start = (/ 1 /), &
3872	count = (/ nzt-nzb+2 /) )
3873	CALL netcdf_handle_error( 'netcdf_define_header', 239 )
3874
3875	ENDDO
3876	ENDDO
3877
3878
3879	CASE ( 'pr_ext' )
3880
3881	!
3882	!-- Get the list of variables and compare with the actual run.
3883	!-- First var_list_old has to be reset, since GET_ATT does not assign
3884	!-- trailing blanks.
3885	var_list_old = ' '
3886	nc_stat = NF90_GET_ATT( id_set_pr, NF90_GLOBAL, 'VAR_LIST', &
3887	var_list_old )
3888	CALL netcdf_handle_error( 'netcdf_define_header', 240 )
3889
3890	var_list = ';'
3891	DO i = 1, dopr_n
3892
3893	IF ( statistic_regions == 0 ) THEN
3894	var_list = TRIM( var_list ) // TRIM( data_output_pr(i) ) // ';'
3895	ELSE
3896	DO j = 0, statistic_regions
3897	WRITE ( suffix, '(''_'',I1)' ) j
3898	var_list = TRIM( var_list ) // TRIM( data_output_pr(i) ) // &
3899	suffix // ';'
3900	ENDDO
3901	ENDIF
3902
3903	ENDDO
3904
3905	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
3906	message_string = 'netCDF file for vertical profiles ' // &
3907	'from previous run found,' // &
3908	'& but this file cannot be extended due to' // &
3909	' variable mismatch.' // &
3910	'&New file is created instead.'
3911	CALL message( 'define_netcdf_header', 'PA0254', 0, 1, 0, 6, 0 )
3912	extend = .FALSE.
3913	RETURN
3914	ENDIF
3915
3916	!
3917	!-- Get the id of the time coordinate (unlimited coordinate) and its
3918	!-- last index on the file. The next time level is dopr..count+1.
3919	!-- The current time must be larger than the last output time
3920	!-- on the file.
3921	nc_stat = NF90_INQ_VARID( id_set_pr, 'time', id_var_time_pr )
3922	CALL netcdf_handle_error( 'netcdf_define_header', 241 )
3923
3924	nc_stat = NF90_INQUIRE_VARIABLE( id_set_pr, id_var_time_pr, &
3925	dimids = id_dim_time_old )
3926	CALL netcdf_handle_error( 'netcdf_define_header', 242 )
3927	id_dim_time_pr = id_dim_time_old(1)
3928
3929	nc_stat = NF90_INQUIRE_DIMENSION( id_set_pr, id_dim_time_pr, &
3930	len = dopr_time_count )
3931	CALL netcdf_handle_error( 'netcdf_define_header', 243 )
3932
3933	nc_stat = NF90_GET_VAR( id_set_pr, id_var_time_pr, &
3934	last_time_coordinate, &
3935	start = (/ dopr_time_count /), &
3936	count = (/ 1 /) )
3937	CALL netcdf_handle_error( 'netcdf_define_header', 244 )
3938
3939	IF ( last_time_coordinate(1) >= simulated_time ) THEN
3940	message_string = 'netCDF file for vertical profiles ' // &
3941	'from previous run found,' // &
3942	'&but this file cannot be extended becaus' // &
3943	'e the current output time' // &
3944	'&is less or equal than the last output t' // &
3945	'ime on this file.' // &
3946	'&New file is created instead.'
3947	CALL message( 'define_netcdf_header', 'PA0255', 0, 1, 0, 6, 0 )
3948	dopr_time_count = 0
3949	extend = .FALSE.
3950	RETURN
3951	ENDIF
3952
3953	!
3954	!-- Dataset seems to be extendable.
3955	!-- Now get the variable ids.
3956	i = 1
3957	DO i = 1, dopr_n
3958
3959	IF ( statistic_regions == 0 ) THEN
3960	nc_stat = NF90_INQ_VARID( id_set_pr, data_output_pr(i), &
3961	id_var_dopr(i,0) )
3962	CALL netcdf_handle_error( 'netcdf_define_header', 245 )
3963	ELSE
3964	DO j = 0, statistic_regions
3965	WRITE ( suffix, '(''_'',I1)' ) j
3966	netcdf_var_name = TRIM( data_output_pr(i) ) // suffix
3967	nc_stat = NF90_INQ_VARID( id_set_pr, netcdf_var_name, &
3968	id_var_dopr(i,j) )
3969	CALL netcdf_handle_error( 'netcdf_define_header', 246 )
3970	ENDDO
3971	ENDIF
3972
3973	ENDDO
3974
3975	!
3976	!-- Get ids of the normalization variables
3977	DO i = 1, dopr_norm_num
3978	nc_stat = NF90_INQ_VARID( id_set_pr, &
3979	'NORM_' // TRIM( dopr_norm_names(i) ), &
3980	id_var_norm_dopr(i) )
3981	CALL netcdf_handle_error( 'netcdf_define_header', 247 )
3982	ENDDO
3983
3984	!
3985	!-- Update the title attribute on file
3986	!-- In order to avoid 'data mode' errors if updated attributes are larger
3987	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
3988	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
3989	!-- performance loss due to data copying; an alternative strategy would be
3990	!-- to ensure equal attribute size in a job chain. Maybe revise later.
3991	IF ( averaging_interval_pr == 0.0_wp ) THEN
3992	time_average_text = ' '
3993	ELSE
3994	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
3995	averaging_interval_pr
3996	ENDIF
3997	nc_stat = NF90_REDEF( id_set_pr )
3998	CALL netcdf_handle_error( 'netcdf_define_header', 437 )
3999	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'title', &
4000	TRIM( run_description_header ) // &
4001	TRIM( time_average_text ) )
4002	CALL netcdf_handle_error( 'netcdf_define_header', 248 )
4003
4004	nc_stat = NF90_ENDDEF( id_set_pr )
4005	CALL netcdf_handle_error( 'netcdf_define_header', 438 )
4006	message_string = 'netCDF file for vertical profiles ' // &
4007	'from previous run found.' // &
4008	'&This file will be extended.'
4009	CALL message( 'define_netcdf_header', 'PA0256', 0, 0, 0, 6, 0 )
4010
4011
4012	CASE ( 'ts_new' )
4013
4014	!
4015	!-- Define some global attributes of the dataset
4016	nc_stat = NF90_PUT_ATT( id_set_ts, NF90_GLOBAL, 'title', &
4017	TRIM( run_description_header ) )
4018	CALL netcdf_handle_error( 'netcdf_define_header', 249 )
4019
4020	!
4021	!-- Define time coordinate for time series (unlimited dimension)
4022	CALL netcdf_create_dim( id_set_ts, 'time', NF90_UNLIMITED, &
4023	id_dim_time_ts, 250 )
4024	CALL netcdf_create_var( id_set_ts, (/ id_dim_time_ts /), 'time', &
4025	NF90_DOUBLE, id_var_time_ts, 'seconds', '', &
4026	251, 252, 000 )
4027	!
4028	!-- Define the variables
4029	var_list = ';'
4030	DO i = 1, dots_num
4031
4032	IF ( statistic_regions == 0 ) THEN
4033
4034	CALL netcdf_create_var( id_set_ts, (/ id_dim_time_ts /), &
4035	dots_label(i), nc_precision(6), &
4036	id_var_dots(i,0), &
4037	TRIM( dots_unit(i) ), &
4038	TRIM( dots_label(i) ), 253, 254, 255 )
4039
4040	ELSE
4041	!
4042	!-- If statistic regions are defined, add suffix _SR+#SR to the
4043	!-- names
4044	DO j = 0, statistic_regions
4045	WRITE ( suffix, '(''_'',I1)' ) j
4046
4047	CALL netcdf_create_var( id_set_ts, (/ id_dim_time_ts /), &
4048	TRIM( dots_label(i) ) // suffix, &
4049	nc_precision(6), id_var_dots(i,j), &
4050	TRIM( dots_unit(i) ), &
4051	TRIM( dots_label(i) ) // ' SR ' // &
4052	suffix(2:2), 256, 257, 347)
4053
4054	ENDDO
4055
4056	ENDIF
4057
4058	ENDDO
4059
4060	!
4061	!-- Write the list of variables as global attribute (this is used by
4062	!-- restart runs)
4063	nc_stat = NF90_PUT_ATT( id_set_ts, NF90_GLOBAL, 'VAR_LIST', var_list )
4064	CALL netcdf_handle_error( 'netcdf_define_header', 258 )
4065
4066	!
4067	!-- Leave netCDF define mode
4068	nc_stat = NF90_ENDDEF( id_set_ts )
4069	CALL netcdf_handle_error( 'netcdf_define_header', 259 )
4070
4071
4072	CASE ( 'ts_ext' )
4073
4074	!
4075	!-- Get the list of variables and compare with the actual run.
4076	!-- First var_list_old has to be reset, since GET_ATT does not assign
4077	!-- trailing blanks.
4078	var_list_old = ' '
4079	nc_stat = NF90_GET_ATT( id_set_ts, NF90_GLOBAL, 'VAR_LIST', &
4080	var_list_old )
4081	CALL netcdf_handle_error( 'netcdf_define_header', 260 )
4082
4083	var_list = ';'
4084	i = 1
4085	DO i = 1, dots_num
4086
4087	IF ( statistic_regions == 0 ) THEN
4088	var_list = TRIM( var_list ) // TRIM( dots_label(i) ) // ';'
4089	ELSE
4090	DO j = 0, statistic_regions
4091	WRITE ( suffix, '(''_'',I1)' ) j
4092	var_list = TRIM( var_list ) // TRIM( dots_label(i) ) // &
4093	suffix // ';'
4094	ENDDO
4095	ENDIF
4096
4097	ENDDO
4098
4099	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
4100	message_string = 'netCDF file for time series ' // &
4101	'from previous run found,' // &
4102	'& but this file cannot be extended due to' // &
4103	' variable mismatch.' // &
4104	'&New file is created instead.'
4105	CALL message( 'define_netcdf_header', 'PA0257', 0, 1, 0, 6, 0 )
4106	extend = .FALSE.
4107	RETURN
4108	ENDIF
4109
4110	!
4111	!-- Get the id of the time coordinate (unlimited coordinate) and its
4112	!-- last index on the file. The next time level is dots..count+1.
4113	!-- The current time must be larger than the last output time
4114	!-- on the file.
4115	nc_stat = NF90_INQ_VARID( id_set_ts, 'time', id_var_time_ts )
4116	CALL netcdf_handle_error( 'netcdf_define_header', 261 )
4117
4118	nc_stat = NF90_INQUIRE_VARIABLE( id_set_ts, id_var_time_ts, &
4119	dimids = id_dim_time_old )
4120	CALL netcdf_handle_error( 'netcdf_define_header', 262 )
4121	id_dim_time_ts = id_dim_time_old(1)
4122
4123	nc_stat = NF90_INQUIRE_DIMENSION( id_set_ts, id_dim_time_ts, &
4124	len = dots_time_count )
4125	CALL netcdf_handle_error( 'netcdf_define_header', 263 )
4126
4127	nc_stat = NF90_GET_VAR( id_set_ts, id_var_time_ts, &
4128	last_time_coordinate, &
4129	start = (/ dots_time_count /), &
4130	count = (/ 1 /) )
4131	CALL netcdf_handle_error( 'netcdf_define_header', 264 )
4132
4133	IF ( last_time_coordinate(1) >= simulated_time ) THEN
4134	message_string = 'netCDF file for time series ' // &
4135	'from previous run found,' // &
4136	'&but this file cannot be extended becaus' // &
4137	'e the current output time' // &
4138	'&is less or equal than the last output t' // &
4139	'ime on this file.' // &
4140	'&New file is created instead.'
4141	CALL message( 'define_netcdf_header', 'PA0258', 0, 1, 0, 6, 0 )
4142	dots_time_count = 0
4143	extend = .FALSE.
4144	RETURN
4145	ENDIF
4146
4147	!
4148	!-- Dataset seems to be extendable.
4149	!-- Now get the variable ids
4150	i = 1
4151	DO i = 1, dots_num
4152
4153	IF ( statistic_regions == 0 ) THEN
4154	nc_stat = NF90_INQ_VARID( id_set_ts, dots_label(i), &
4155	id_var_dots(i,0) )
4156	CALL netcdf_handle_error( 'netcdf_define_header', 265 )
4157	ELSE
4158	DO j = 0, statistic_regions
4159	WRITE ( suffix, '(''_'',I1)' ) j
4160	netcdf_var_name = TRIM( dots_label(i) ) // suffix
4161	nc_stat = NF90_INQ_VARID( id_set_ts, netcdf_var_name, &
4162	id_var_dots(i,j) )
4163	CALL netcdf_handle_error( 'netcdf_define_header', 266 )
4164	ENDDO
4165	ENDIF
4166
4167	ENDDO
4168
4169	!
4170	!-- Update the title attribute on file
4171	!-- In order to avoid 'data mode' errors if updated attributes are larger
4172	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
4173	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
4174	!-- performance loss due to data copying; an alternative strategy would be
4175	!-- to ensure equal attribute size in a job chain. Maybe revise later.
4176	nc_stat = NF90_REDEF( id_set_ts )
4177	CALL netcdf_handle_error( 'netcdf_define_header', 439 )
4178	nc_stat = NF90_PUT_ATT( id_set_ts, NF90_GLOBAL, 'title', &
4179	TRIM( run_description_header ) )
4180	CALL netcdf_handle_error( 'netcdf_define_header', 267 )
4181	nc_stat = NF90_ENDDEF( id_set_ts )
4182	CALL netcdf_handle_error( 'netcdf_define_header', 440 )
4183	message_string = 'netCDF file for time series ' // &
4184	'from previous run found.' // &
4185	'&This file will be extended.'
4186	CALL message( 'define_netcdf_header', 'PA0259', 0, 0, 0, 6, 0 )
4187
4188
4189	CASE ( 'sp_new' )
4190
4191	!
4192	!-- Define some global attributes of the dataset
4193	IF ( averaging_interval_sp /= 0.0_wp ) THEN
4194	WRITE (time_average_text,'('', '',F7.1,'' s average'')') &
4195	averaging_interval_sp
4196	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
4197	TRIM( run_description_header ) // &
4198	TRIM( time_average_text ) )
4199	CALL netcdf_handle_error( 'netcdf_define_header', 268 )
4200
4201	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval_sp
4202	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'time_avg', &
4203	TRIM( time_average_text ) )
4204	ELSE
4205	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
4206	TRIM( run_description_header ) )
4207	ENDIF
4208	CALL netcdf_handle_error( 'netcdf_define_header', 269 )
4209
4210	!
4211	!-- Define time coordinate for spectra (unlimited dimension)
4212	CALL netcdf_create_dim( id_set_sp, 'time', NF90_UNLIMITED, &
4213	id_dim_time_sp, 270 )
4214	CALL netcdf_create_var( id_set_sp, (/ id_dim_time_sp /), 'time', &
4215	NF90_DOUBLE, id_var_time_sp, 'seconds', '', &
4216	271, 272, 000 )
4217	!
4218	!-- Define the spatial dimensions and coordinates for spectra.
4219	!-- First, determine the number of vertical levels for which spectra
4220	!-- are to be output.
4221	ns = 1
4222	DO WHILE ( comp_spectra_level(ns) /= 999999 .AND. ns <= 100 )
4223	ns = ns + 1
4224	ENDDO
4225	ns = ns - 1
4226
4227	!
4228	!-- Define vertical coordinate grid (zu grid)
4229	CALL netcdf_create_dim( id_set_sp, 'zu_sp', ns, id_dim_zu_sp, 273 )
4230	CALL netcdf_create_var( id_set_sp, (/ id_dim_zu_sp /), 'zu_sp', &
4231	NF90_DOUBLE, id_var_zu_sp, 'meters', '', &
4232	274, 275, 000 )
4233	!
4234	!-- Define vertical coordinate grid (zw grid)
4235	CALL netcdf_create_dim( id_set_sp, 'zw_sp', ns, id_dim_zw_sp, 276 )
4236	CALL netcdf_create_var( id_set_sp, (/ id_dim_zw_sp /), 'zw_sp', &
4237	NF90_DOUBLE, id_var_zw_sp, 'meters', '', &
4238	277, 278, 000 )
4239	!
4240	!-- Define x-axis
4241	CALL netcdf_create_dim( id_set_sp, 'k_x', nx/2, id_dim_x_sp, 279 )
4242	CALL netcdf_create_var( id_set_sp, (/ id_dim_x_sp /), 'k_x', &
4243	NF90_DOUBLE, id_var_x_sp, 'm-1', '', 280, &
4244	281, 000 )
4245	!
4246	!-- Define y-axis
4247	CALL netcdf_create_dim(id_set_sp, 'k_y', ny/2, id_dim_y_sp, 282 )
4248	CALL netcdf_create_var( id_set_sp, (/ id_dim_y_sp /), 'k_y', &
4249	NF90_DOUBLE, id_var_y_sp, 'm-1', '', 283, &
4250	284, 000 )
4251	!
4252	!-- Define the variables
4253	var_list = ';'
4254	i = 1
4255	DO WHILE ( data_output_sp(i) /= ' ' .AND. i <= 10 )
4256	!
4257	!-- First check for the vertical grid
4258	found = .TRUE.
4259	SELECT CASE ( data_output_sp(i) )
4260	!
4261	!-- Most variables are defined on the zu levels
4262	CASE ( 'e', 'lpt', 'nr', 'p', 'pc', 'pr', 'prr', 'pt', 'q', &
4263	'qc', 'ql', 'ql_c', 'ql_v', 'ql_vp', 'qr', 'qv', 'rho',&
4264	's', 'sa', 'u', 'v', 'vpt' )
4265
4266	grid_z = 'zu'
4267	!
4268	!-- zw levels
4269	CASE ( 'w' )
4270
4271	grid_z = 'zw'
4272
4273	CASE DEFAULT
4274	!
4275	!-- Check for user-defined quantities (found, grid_x and grid_y
4276	!-- are dummies)
4277	CALL user_define_netcdf_grid( data_output_sp(i), found, &
4278	grid_x, grid_y, grid_z )
4279
4280	END SELECT
4281
4282	IF ( INDEX( spectra_direction(i), 'x' ) /= 0 ) THEN
4283
4284	!
4285	!-- Define the variable
4286	netcdf_var_name = TRIM( data_output_sp(i) ) // '_x'
4287	IF ( TRIM( grid_z ) == 'zw' ) THEN
4288	CALL netcdf_create_var( id_set_sp, (/ id_dim_x_sp, &
4289	id_dim_zw_sp, id_dim_time_sp /), &
4290	netcdf_var_name, nc_precision(7), &
4291	id_var_dospx(i), 'unknown', &
4292	netcdf_var_name, 285, 286, 287 )
4293	ELSE
4294	CALL netcdf_create_var( id_set_sp, (/ id_dim_x_sp, &
4295	id_dim_zu_sp, id_dim_time_sp /), &
4296	netcdf_var_name, nc_precision(7), &
4297	id_var_dospx(i), 'unknown', &
4298	netcdf_var_name, 285, 286, 287 )
4299	ENDIF
4300
4301	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
4302
4303	ENDIF
4304
4305	IF ( INDEX( spectra_direction(i), 'y' ) /= 0 ) THEN
4306
4307	!
4308	!-- Define the variable
4309	netcdf_var_name = TRIM( data_output_sp(i) ) // '_y'
4310	IF ( TRIM( grid_z ) == 'zw' ) THEN
4311	CALL netcdf_create_var( id_set_sp, (/ id_dim_y_sp, &
4312	id_dim_zw_sp, id_dim_time_sp /), &
4313	netcdf_var_name, nc_precision(7), &
4314	id_var_dospy(i), 'unknown', &
4315	netcdf_var_name, 288, 289, 290 )
4316	ELSE
4317	CALL netcdf_create_var( id_set_sp, (/ id_dim_y_sp, &
4318	id_dim_zu_sp, id_dim_time_sp /), &
4319	netcdf_var_name, nc_precision(7), &
4320	id_var_dospy(i), 'unknown', &
4321	netcdf_var_name, 288, 289, 290 )
4322	ENDIF
4323
4324	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
4325
4326	ENDIF
4327
4328	i = i + 1
4329
4330	ENDDO
4331
4332	!
4333	!-- Write the list of variables as global attribute (this is used by
4334	!-- restart runs)
4335	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'VAR_LIST', var_list )
4336	CALL netcdf_handle_error( 'netcdf_define_header', 291 )
4337
4338	!
4339	!-- Leave netCDF define mode
4340	nc_stat = NF90_ENDDEF( id_set_sp )
4341	CALL netcdf_handle_error( 'netcdf_define_header', 292 )
4342
4343	!
4344	!-- Write axis data: zu_sp, zw_sp, k_x, k_y
4345	ALLOCATE( netcdf_data(1:ns) )
4346
4347	!
4348	!-- Write zu data
4349	netcdf_data(1:ns) = zu( comp_spectra_level(1:ns) )
4350	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_zu_sp, netcdf_data, &
4351	start = (/ 1 /), count = (/ ns /) )
4352	CALL netcdf_handle_error( 'netcdf_define_header', 293 )
4353
4354	!
4355	!-- Write zw data
4356	netcdf_data(1:ns) = zw( comp_spectra_level(1:ns) )
4357	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_zw_sp, netcdf_data, &
4358	start = (/ 1 /), count = (/ ns /) )
4359	CALL netcdf_handle_error( 'netcdf_define_header', 294 )
4360
4361	DEALLOCATE( netcdf_data )
4362
4363	!
4364	!-- Write data for x and y axis (wavenumbers)
4365	ALLOCATE( netcdf_data(nx/2) )
4366	DO i = 1, nx/2
4367	netcdf_data(i) = 2.0_wp * pi * i / ( dx * ( nx + 1 ) )
4368	ENDDO
4369
4370	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_x_sp, netcdf_data, &
4371	start = (/ 1 /), count = (/ nx/2 /) )
4372	CALL netcdf_handle_error( 'netcdf_define_header', 295 )
4373
4374	DEALLOCATE( netcdf_data )
4375
4376	ALLOCATE( netcdf_data(ny/2) )
4377	DO i = 1, ny/2
4378	netcdf_data(i) = 2.0_wp * pi * i / ( dy * ( ny + 1 ) )
4379	ENDDO
4380
4381	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_y_sp, netcdf_data, &
4382	start = (/ 1 /), count = (/ ny/2 /) )
4383	CALL netcdf_handle_error( 'netcdf_define_header', 296 )
4384
4385	DEALLOCATE( netcdf_data )
4386
4387
4388	CASE ( 'sp_ext' )
4389
4390	!
4391	!-- Get the list of variables and compare with the actual run.
4392	!-- First var_list_old has to be reset, since GET_ATT does not assign
4393	!-- trailing blanks.
4394	var_list_old = ' '
4395	nc_stat = NF90_GET_ATT( id_set_sp, NF90_GLOBAL, 'VAR_LIST', &
4396	var_list_old )
4397	CALL netcdf_handle_error( 'netcdf_define_header', 297 )
4398	var_list = ';'
4399	i = 1
4400	DO WHILE ( data_output_sp(i) /= ' ' .AND. i <= 10 )
4401
4402	IF ( INDEX( spectra_direction(i), 'x' ) /= 0 ) THEN
4403	netcdf_var_name = TRIM( data_output_sp(i) ) // '_x'
4404	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
4405	ENDIF
4406
4407	IF ( INDEX( spectra_direction(i), 'y' ) /= 0 ) THEN
4408	netcdf_var_name = TRIM( data_output_sp(i) ) // '_y'
4409	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
4410	ENDIF
4411
4412	i = i + 1
4413
4414	ENDDO
4415
4416	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
4417	message_string = 'netCDF file for spectra ' // &
4418	'from previous run found,' // &
4419	'& but this file cannot be extended due to' // &
4420	' variable mismatch.' // &
4421	'&New file is created instead.'
4422	CALL message( 'define_netcdf_header', 'PA0260', 0, 1, 0, 6, 0 )
4423	extend = .FALSE.
4424	RETURN
4425	ENDIF
4426
4427	!
4428	!-- Determine the number of current vertical levels for which spectra
4429	!-- shall be output
4430	ns = 1
4431	DO WHILE ( comp_spectra_level(ns) /= 999999 .AND. ns <= 100 )
4432	ns = ns + 1
4433	ENDDO
4434	ns = ns - 1
4435
4436	!
4437	!-- Get and compare the number of vertical levels
4438	nc_stat = NF90_INQ_VARID( id_set_sp, 'zu_sp', id_var_zu_sp )
4439	CALL netcdf_handle_error( 'netcdf_define_header', 298 )
4440
4441	nc_stat = NF90_INQUIRE_VARIABLE( id_set_sp, id_var_zu_sp, &
4442	dimids = id_dim_zu_sp_old )
4443	CALL netcdf_handle_error( 'netcdf_define_header', 299 )
4444	id_dim_zu_sp = id_dim_zu_sp_old(1)
4445
4446	nc_stat = NF90_INQUIRE_DIMENSION( id_set_sp, id_dim_zu_sp, &
4447	len = ns_old )
4448	CALL netcdf_handle_error( 'netcdf_define_header', 300 )
4449
4450	IF ( ns /= ns_old ) THEN
4451	message_string = 'netCDF file for spectra ' // &
4452	' from previous run found,' // &
4453	'&but this file cannot be extended due to' // &
4454	' mismatch in number of' // &
4455	'&vertical levels.' // &
4456	'&New file is created instead.'
4457	CALL message( 'define_netcdf_header', 'PA0261', 0, 1, 0, 6, 0 )
4458	extend = .FALSE.
4459	RETURN
4460	ENDIF
4461
4462	!
4463	!-- Get and compare the heights of the cross sections
4464	ALLOCATE( netcdf_data(1:ns_old) )
4465
4466	nc_stat = NF90_GET_VAR( id_set_sp, id_var_zu_sp, netcdf_data )
4467	CALL netcdf_handle_error( 'netcdf_define_header', 301 )
4468
4469	DO i = 1, ns
4470	IF ( zu(comp_spectra_level(i)) /= netcdf_data(i) ) THEN
4471	message_string = 'netCDF file for spectra ' // &
4472	' from previous run found,' // &
4473	'&but this file cannot be extended due to' // &
4474	' mismatch in heights of' // &
4475	'&vertical levels.' // &
4476	'&New file is created instead.'
4477	CALL message( 'define_netcdf_header', 'PA0262', 0, 1, 0, 6, 0 )
4478	extend = .FALSE.
4479	RETURN
4480	ENDIF
4481	ENDDO
4482
4483	DEALLOCATE( netcdf_data )
4484
4485	!
4486	!-- Get the id of the time coordinate (unlimited coordinate) and its
4487	!-- last index on the file. The next time level is plsp..count+1.
4488	!-- The current time must be larger than the last output time
4489	!-- on the file.
4490	nc_stat = NF90_INQ_VARID( id_set_sp, 'time', id_var_time_sp )
4491	CALL netcdf_handle_error( 'netcdf_define_header', 302 )
4492
4493	nc_stat = NF90_INQUIRE_VARIABLE( id_set_sp, id_var_time_sp, &
4494	dimids = id_dim_time_old )
4495	CALL netcdf_handle_error( 'netcdf_define_header', 303 )
4496	id_dim_time_sp = id_dim_time_old(1)
4497
4498	nc_stat = NF90_INQUIRE_DIMENSION( id_set_sp, id_dim_time_sp, &
4499	len = dosp_time_count )
4500	CALL netcdf_handle_error( 'netcdf_define_header', 304 )
4501
4502	nc_stat = NF90_GET_VAR( id_set_sp, id_var_time_sp, &
4503	last_time_coordinate, &
4504	start = (/ dosp_time_count /), &
4505	count = (/ 1 /) )
4506	CALL netcdf_handle_error( 'netcdf_define_header', 305 )
4507
4508	IF ( last_time_coordinate(1) >= simulated_time ) THEN
4509	message_string = 'netCDF file for spectra ' // &
4510	'from previous run found,' // &
4511	'&but this file cannot be extended becaus' // &
4512	'e the current output time' // &
4513	'&is less or equal than the last output t' // &
4514	'ime on this file.' // &
4515	'&New file is created instead.'
4516	CALL message( 'define_netcdf_header', 'PA0263', 0, 1, 0, 6, 0 )
4517	dosp_time_count = 0
4518	extend = .FALSE.
4519	RETURN
4520	ENDIF
4521
4522	!
4523	!-- Dataset seems to be extendable.
4524	!-- Now get the variable ids.
4525	i = 1
4526	DO WHILE ( data_output_sp(i) /= ' ' .AND. i <= 10 )
4527
4528	IF ( INDEX( spectra_direction(i), 'x' ) /= 0 ) THEN
4529	netcdf_var_name = TRIM( data_output_sp(i) ) // '_x'
4530	nc_stat = NF90_INQ_VARID( id_set_sp, netcdf_var_name, &
4531	id_var_dospx(i) )
4532	CALL netcdf_handle_error( 'netcdf_define_header', 306 )
4533	ENDIF
4534
4535	IF ( INDEX( spectra_direction(i), 'y' ) /= 0 ) THEN
4536	netcdf_var_name = TRIM( data_output_sp(i) ) // '_y'
4537	nc_stat = NF90_INQ_VARID( id_set_sp, netcdf_var_name, &
4538	id_var_dospy(i) )
4539	CALL netcdf_handle_error( 'netcdf_define_header', 307 )
4540	ENDIF
4541
4542	i = i + 1
4543
4544	ENDDO
4545
4546	!
4547	!-- Update the title attribute on file
4548	!-- In order to avoid 'data mode' errors if updated attributes are larger
4549	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
4550	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
4551	!-- performance loss due to data copying; an alternative strategy would be
4552	!-- to ensure equal attribute size in a job chain. Maybe revise later.
4553	nc_stat = NF90_REDEF( id_set_sp )
4554	CALL netcdf_handle_error( 'netcdf_define_header', 441 )
4555	IF ( averaging_interval_sp /= 0.0_wp ) THEN
4556	WRITE (time_average_text,'('', '',F7.1,'' s average'')') &
4557	averaging_interval_sp
4558	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
4559	TRIM( run_description_header ) // &
4560	TRIM( time_average_text ) )
4561	CALL netcdf_handle_error( 'netcdf_define_header', 308 )
4562
4563	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval_sp
4564	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'time_avg', &
4565	TRIM( time_average_text ) )
4566	ELSE
4567	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
4568	TRIM( run_description_header ) )
4569	ENDIF
4570	CALL netcdf_handle_error( 'netcdf_define_header', 309 )
4571	nc_stat = NF90_ENDDEF( id_set_sp )
4572	CALL netcdf_handle_error( 'netcdf_define_header', 442 )
4573	message_string = 'netCDF file for spectra ' // &
4574	'from previous run found.' // &
4575	'&This file will be extended.'
4576	CALL message( 'define_netcdf_header', 'PA0264', 0, 0, 0, 6, 0 )
4577
4578
4579	CASE ( 'pt_new' )
4580
4581	!
4582	!-- Define some global attributes of the dataset
4583	nc_stat = NF90_PUT_ATT( id_set_prt, NF90_GLOBAL, 'title', &
4584	TRIM( run_description_header ) )
4585	CALL netcdf_handle_error( 'netcdf_define_header', 310 )
4586
4587	!
4588	!-- Define time coordinate for particles (unlimited dimension)
4589	CALL netcdf_create_dim( id_set_prt, 'time', NF90_UNLIMITED, &
4590	id_dim_time_prt, 311 )
4591	CALL netcdf_create_var( id_set_prt, (/ id_dim_time_prt /), 'time', &
4592	NF90_DOUBLE, id_var_time_prt, 'seconds', '', &
4593	312, 313, 000 )
4594	!
4595	!-- Define particle coordinate (maximum particle number)
4596	IF ( netcdf_data_format < 3 ) THEN
4597	CALL netcdf_create_dim( id_set_prt, 'particle_number', &
4598	maximum_number_of_particles, &
4599	id_dim_prtnum, 314 )
4600	ELSE
4601	!
4602	!-- netCDF4 allows more than one unlimited dimension
4603	CALL netcdf_create_dim( id_set_prt, 'particle_number', &
4604	NF90_UNLIMITED, id_dim_prtnum, 314 )
4605	ENDIF
4606
4607	CALL netcdf_create_var( id_set_prt, (/ id_dim_prtnum /), &
4608	'particle_number', NF90_DOUBLE, &
4609	id_var_prtnum, 'particle number', '', 315, &
4610	316, 000 )
4611	!
4612	!-- Define variable which contains the real number of particles in use
4613	CALL netcdf_create_var( id_set_prt, (/ id_dim_time_prt /), &
4614	'real_num_of_prt', NF90_DOUBLE, &
4615	id_var_rnop_prt, 'particle number', '', 317, &
4616	318, 000 )
4617	!
4618	!-- Define the variables
4619	DO i = 1, 17
4620
4621	CALL netcdf_create_var( id_set_prt, (/ id_dim_prtnum, &
4622	id_dim_time_prt /), prt_var_names(i), &
4623	nc_precision(8), id_var_prt(i), &
4624	TRIM( prt_var_units(i) ), &
4625	TRIM( prt_var_names(i) ), 319, 320, 321 )
4626
4627	ENDDO
4628
4629	!
4630	!-- Leave netCDF define mode
4631	nc_stat = NF90_ENDDEF( id_set_prt )
4632	CALL netcdf_handle_error( 'netcdf_define_header', 322 )
4633
4634
4635	CASE ( 'pt_ext' )
4636
4637	!
4638	!-- Get the id of the time coordinate (unlimited coordinate) and its
4639	!-- last index on the file. The next time level is prt..count+1.
4640	!-- The current time must be larger than the last output time
4641	!-- on the file.
4642	nc_stat = NF90_INQ_VARID( id_set_prt, 'time', id_var_time_prt )
4643	CALL netcdf_handle_error( 'netcdf_define_header', 323 )
4644
4645	nc_stat = NF90_INQUIRE_VARIABLE( id_set_prt, id_var_time_prt, &
4646	dimids = id_dim_time_old )
4647	CALL netcdf_handle_error( 'netcdf_define_header', 324 )
4648	id_dim_time_prt = id_dim_time_old(1)
4649
4650	nc_stat = NF90_INQUIRE_DIMENSION( id_set_prt, id_dim_time_prt, &
4651	len = prt_time_count )
4652	CALL netcdf_handle_error( 'netcdf_define_header', 325 )
4653
4654	nc_stat = NF90_GET_VAR( id_set_prt, id_var_time_prt, &
4655	last_time_coordinate, &
4656	start = (/ prt_time_count /), &
4657	count = (/ 1 /) )
4658	CALL netcdf_handle_error( 'netcdf_define_header', 326 )
4659
4660	IF ( last_time_coordinate(1) >= simulated_time ) THEN
4661	message_string = 'netCDF file for particles ' // &
4662	'from previous run found,' // &
4663	'&but this file cannot be extended becaus' // &
4664	'e the current output time' // &
4665	'&is less or equal than the last output t' // &
4666	'ime on this file.' // &
4667	'&New file is created instead.'
4668	CALL message( 'define_netcdf_header', 'PA0265', 0, 1, 0, 6, 0 )
4669	prt_time_count = 0
4670	extend = .FALSE.
4671	RETURN
4672	ENDIF
4673
4674	!
4675	!-- Dataset seems to be extendable.
4676	!-- Now get the variable ids.
4677	nc_stat = NF90_INQ_VARID( id_set_prt, 'real_num_of_prt', &
4678	id_var_rnop_prt )
4679	CALL netcdf_handle_error( 'netcdf_define_header', 327 )
4680
4681	DO i = 1, 17
4682
4683	nc_stat = NF90_INQ_VARID( id_set_prt, prt_var_names(i), &
4684	id_var_prt(i) )
4685	CALL netcdf_handle_error( 'netcdf_define_header', 328 )
4686
4687	ENDDO
4688
4689	message_string = 'netCDF file for particles ' // &
4690	'from previous run found.' // &
4691	'&This file will be extended.'
4692	CALL message( 'define_netcdf_header', 'PA0266', 0, 0, 0, 6, 0 )
4693
4694
4695
4696	CASE ( 'ps_new' )
4697
4698	!
4699	!-- Define some global attributes of the dataset
4700	nc_stat = NF90_PUT_ATT( id_set_pts, NF90_GLOBAL, 'title', &
4701	TRIM( run_description_header ) )
4702	CALL netcdf_handle_error( 'netcdf_define_header', 396 )
4703
4704	!
4705	!-- Define time coordinate for particle time series (unlimited dimension)
4706	CALL netcdf_create_dim( id_set_pts, 'time', NF90_UNLIMITED, &
4707	id_dim_time_pts, 397 )
4708	CALL netcdf_create_var( id_set_pts, (/ id_dim_time_pts /), 'time', &
4709	NF90_DOUBLE, id_var_time_pts, 'seconds', '', &
4710	398, 399, 000 )
4711	!
4712	!-- Define the variables. If more than one particle group is defined,
4713	!-- define seperate variables for each group
4714	var_list = ';'
4715	DO i = 1, dopts_num
4716
4717	DO j = 0, number_of_particle_groups
4718
4719	IF ( j == 0 ) THEN
4720	suffix1 = ''
4721	ELSE
4722	WRITE ( suffix1, '(''_'',I2.2)' ) j
4723	ENDIF
4724
4725	IF ( j == 0 ) THEN
4726	CALL netcdf_create_var( id_set_pts, (/ id_dim_time_pts /), &
4727	TRIM( dopts_label(i) ) // suffix1, &
4728	nc_precision(6), id_var_dopts(i,j), &
4729	TRIM( dopts_unit(i) ), &
4730	TRIM( dopts_label(i) ), 400, 401, &
4731	402 )
4732	ELSE
4733	CALL netcdf_create_var( id_set_pts, (/ id_dim_time_pts /), &
4734	TRIM( dopts_label(i) ) // suffix1, &
4735	nc_precision(6), id_var_dopts(i,j), &
4736	TRIM( dopts_unit(i) ), &
4737	TRIM( dopts_label(i) ) // ' PG ' // &
4738	suffix1(2:3), 400, 401, 402 )
4739	ENDIF
4740
4741	var_list = TRIM( var_list ) // TRIM( dopts_label(i) ) // &
4742	suffix1 // ';'
4743
4744	IF ( number_of_particle_groups == 1 ) EXIT
4745
4746	ENDDO
4747
4748	ENDDO
4749
4750	!
4751	!-- Write the list of variables as global attribute (this is used by
4752	!-- restart runs)
4753	nc_stat = NF90_PUT_ATT( id_set_pts, NF90_GLOBAL, 'VAR_LIST', &
4754	var_list )
4755	CALL netcdf_handle_error( 'netcdf_define_header', 403 )
4756
4757
4758	!
4759	!-- Leave netCDF define mode
4760	nc_stat = NF90_ENDDEF( id_set_pts )
4761	CALL netcdf_handle_error( 'netcdf_define_header', 404 )
4762
4763
4764	CASE ( 'ps_ext' )
4765
4766	!
4767	!-- Get the list of variables and compare with the actual run.
4768	!-- First var_list_old has to be reset, since GET_ATT does not assign
4769	!-- trailing blanks.
4770	var_list_old = ' '
4771	nc_stat = NF90_GET_ATT( id_set_pts, NF90_GLOBAL, 'VAR_LIST', &
4772	var_list_old )
4773	CALL netcdf_handle_error( 'netcdf_define_header', 405 )
4774
4775	var_list = ';'
4776	i = 1
4777	DO i = 1, dopts_num
4778
4779	DO j = 0, number_of_particle_groups
4780
4781	IF ( j == 0 ) THEN
4782	suffix1 = ''
4783	ELSE
4784	WRITE ( suffix1, '(''_'',I2.2)' ) j
4785	ENDIF
4786
4787	var_list = TRIM( var_list ) // TRIM( dopts_label(i) ) // &
4788	suffix1 // ';'
4789
4790	IF ( number_of_particle_groups == 1 ) EXIT
4791
4792	ENDDO
4793
4794	ENDDO
4795
4796	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
4797	message_string = 'netCDF file for particle time series ' // &
4798	'from previous run found,' // &
4799	'& but this file cannot be extended due to' // &
4800	' variable mismatch.' // &
4801	'&New file is created instead.'
4802	CALL message( 'define_netcdf_header', 'PA0267', 0, 1, 0, 6, 0 )
4803	extend = .FALSE.
4804	RETURN
4805	ENDIF
4806
4807	!
4808	!-- Get the id of the time coordinate (unlimited coordinate) and its
4809	!-- last index on the file. The next time level is dots..count+1.
4810	!-- The current time must be larger than the last output time
4811	!-- on the file.
4812	nc_stat = NF90_INQ_VARID( id_set_pts, 'time', id_var_time_pts )
4813	CALL netcdf_handle_error( 'netcdf_define_header', 406 )
4814
4815	nc_stat = NF90_INQUIRE_VARIABLE( id_set_pts, id_var_time_pts, &
4816	dimids = id_dim_time_old )
4817	CALL netcdf_handle_error( 'netcdf_define_header', 407 )
4818	id_dim_time_pts = id_dim_time_old(1)
4819
4820	nc_stat = NF90_INQUIRE_DIMENSION( id_set_pts, id_dim_time_pts, &
4821	len = dopts_time_count )
4822	CALL netcdf_handle_error( 'netcdf_define_header', 408 )
4823
4824	nc_stat = NF90_GET_VAR( id_set_pts, id_var_time_pts, &
4825	last_time_coordinate, &
4826	start = (/ dopts_time_count /), &
4827	count = (/ 1 /) )
4828	CALL netcdf_handle_error( 'netcdf_define_header', 409 )
4829
4830	IF ( last_time_coordinate(1) >= simulated_time ) THEN
4831	message_string = 'netCDF file for particle time series ' // &
4832	'from previous run found,' // &
4833	'&but this file cannot be extended becaus' // &
4834	'e the current output time' // &
4835	'&is less or equal than the last output t' // &
4836	'ime on this file.' // &
4837	'&New file is created instead.'
4838	CALL message( 'define_netcdf_header', 'PA0268', 0, 1, 0, 6, 0 )
4839	dopts_time_count = 0
4840	extend = .FALSE.
4841	RETURN
4842	ENDIF
4843
4844	!
4845	!-- Dataset seems to be extendable.
4846	!-- Now get the variable ids
4847	i = 1
4848	DO i = 1, dopts_num
4849
4850	DO j = 0, number_of_particle_groups
4851
4852	IF ( j == 0 ) THEN
4853	suffix1 = ''
4854	ELSE
4855	WRITE ( suffix1, '(''_'',I2.2)' ) j
4856	ENDIF
4857
4858	netcdf_var_name = TRIM( dopts_label(i) ) // suffix1
4859
4860	nc_stat = NF90_INQ_VARID( id_set_pts, netcdf_var_name, &
4861	id_var_dopts(i,j) )
4862	CALL netcdf_handle_error( 'netcdf_define_header', 410 )
4863
4864	IF ( number_of_particle_groups == 1 ) EXIT
4865
4866	ENDDO
4867
4868	ENDDO
4869
4870	!
4871	!-- Update the title attribute on file
4872	!-- In order to avoid 'data mode' errors if updated attributes are larger
4873	!-- than their original size, NF90_PUT_ATT is called in 'define mode'
4874	!-- enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible
4875	!-- performance loss due to data copying; an alternative strategy would be
4876	!-- to ensure equal attribute size in a job chain. Maybe revise later.
4877	nc_stat = NF90_REDEF( id_set_pts )
4878	CALL netcdf_handle_error( 'netcdf_define_header', 443 )
4879	nc_stat = NF90_PUT_ATT( id_set_pts, NF90_GLOBAL, 'title', &
4880	TRIM( run_description_header ) )
4881	CALL netcdf_handle_error( 'netcdf_define_header', 411 )
4882	nc_stat = NF90_ENDDEF( id_set_pts )
4883	CALL netcdf_handle_error( 'netcdf_define_header', 444 )
4884	message_string = 'netCDF file for particle time series ' // &
4885	'from previous run found.' // &
4886	'&This file will be extended.'
4887	CALL message( 'netcdf_define_header', 'PA0269', 0, 0, 0, 6, 0 )
4888
4889
4890	CASE DEFAULT
4891
4892	message_string = 'mode "' // TRIM( mode) // '" not supported'
4893	CALL message( 'netcdf_define_header', 'PA0270', 0, 0, 0, 6, 0 )
4894
4895	END SELECT
4896
4897	#endif
4898	END SUBROUTINE netcdf_define_header
4899
4900
4901	!------------------------------------------------------------------------------!
4902	! Description:
4903	! ------------
4904	!> Creates a netCDF file and give back the id. The parallel flag has to be TRUE
4905	!> for parallel netCDF output support.
4906	!------------------------------------------------------------------------------!
4907
4908	SUBROUTINE netcdf_create_file( filename , id, parallel, errno )
4909	#if defined( __netcdf )
4910
4911	USE pegrid
4912
4913	IMPLICIT NONE
4914
4915	CHARACTER (LEN=*), INTENT(IN) :: filename
4916	INTEGER, INTENT(IN) :: errno
4917	INTEGER, INTENT(OUT) :: id
4918	LOGICAL, INTENT(IN) :: parallel
4919
4920
4921	!
4922	!-- Create a new netCDF output file with requested netCDF format
4923	IF ( netcdf_data_format == 1 ) THEN
4924	!
4925	!-- Classic netCDF format
4926	nc_stat = NF90_CREATE( filename, NF90_NOCLOBBER, id )
4927
4928	ELSEIF ( netcdf_data_format == 2 ) THEN
4929	!
4930	!-- 64bit-offset format
4931	nc_stat = NF90_CREATE( filename, &
4932	IOR( NF90_NOCLOBBER, NF90_64BIT_OFFSET ), id )
4933
4934	#if defined( __netcdf4 )
4935	ELSEIF ( netcdf_data_format == 3 .OR. &
4936	( .NOT. parallel .AND. netcdf_data_format == 5 ) ) THEN
4937	!
4938	!-- netCDF4/HDF5 format
4939	nc_stat = NF90_CREATE( filename, IOR( NF90_NOCLOBBER, NF90_NETCDF4 ), id )
4940
4941	ELSEIF ( netcdf_data_format == 4 .OR. &
4942	( .NOT. parallel .AND. netcdf_data_format == 6 ) ) THEN
4943	!
4944	!-- netCDF4/HDF5 format with classic model flag
4945	nc_stat = NF90_CREATE( filename, &
4946	IOR( NF90_NOCLOBBER, &
4947	IOR( NF90_CLASSIC_MODEL, NF90_HDF5 ) ), id )
4948
4949	#if defined( __netcdf4_parallel )
4950	ELSEIF ( netcdf_data_format == 5 .AND. parallel ) THEN
4951	!
4952	!-- netCDF4/HDF5 format, parallel
4953	nc_stat = NF90_CREATE( filename, &
4954	IOR( NF90_NOCLOBBER, &
4955	IOR( NF90_NETCDF4, NF90_MPIIO ) ), &
4956	id, COMM = comm2d, INFO = MPI_INFO_NULL )
4957
4958	ELSEIF ( netcdf_data_format == 6 .AND. parallel ) THEN
4959	!
4960	!-- netCDF4/HDF5 format with classic model flag, parallel
4961	nc_stat = NF90_CREATE( filename, &
4962	IOR( NF90_NOCLOBBER, &
4963	IOR( NF90_MPIIO, &
4964	IOR( NF90_CLASSIC_MODEL, NF90_HDF5 ) ) ), &
4965	id, COMM = comm2d, INFO = MPI_INFO_NULL )
4966
4967	#endif
4968	#endif
4969	ENDIF
4970
4971	CALL netcdf_handle_error( 'netcdf_create_file', errno )
4972	#endif
4973	END SUBROUTINE netcdf_create_file
4974
4975
4976	!------------------------------------------------------------------------------!
4977	! Description:
4978	! ------------
4979	!> Opens an existing netCDF file for writing and gives back the id.
4980	!> The parallel flag has to be TRUE for parallel netCDF output support.
4981	!------------------------------------------------------------------------------!
4982
4983	SUBROUTINE netcdf_open_write_file( filename, id, parallel, errno )
4984	#if defined( __netcdf )
4985
4986	USE pegrid
4987
4988	IMPLICIT NONE
4989
4990	CHARACTER (LEN=*), INTENT(IN) :: filename
4991	INTEGER, INTENT(IN) :: errno
4992	INTEGER, INTENT(OUT) :: id
4993	LOGICAL, INTENT(IN) :: parallel
4994
4995
4996	IF ( netcdf_data_format < 5 .OR. .NOT. parallel ) THEN
4997	nc_stat = NF90_OPEN( filename, NF90_WRITE, id )
4998	#if defined( __netcdf4 )
4999	#if defined( __netcdf4_parallel )
5000	ELSEIF ( netcdf_data_format > 4 .AND. parallel ) THEN
5001	nc_stat = NF90_OPEN( filename, IOR( NF90_WRITE, NF90_MPIIO ), id, &
5002	COMM = comm2d, INFO = MPI_INFO_NULL )
5003	#endif
5004	#endif
5005	ENDIF
5006
5007	CALL netcdf_handle_error( 'netcdf_open_write_file', errno )
5008	#endif
5009	END SUBROUTINE netcdf_open_write_file
5010
5011
5012	!------------------------------------------------------------------------------!
5013	! Description:
5014	! ------------
5015	!> Prints out a text message corresponding to the current status.
5016	!------------------------------------------------------------------------------!
5017
5018	SUBROUTINE netcdf_handle_error( routine_name, errno )
5019	#if defined( __netcdf )
5020
5021
5022	USE control_parameters, &
5023	ONLY: message_string
5024
5025	IMPLICIT NONE
5026
5027	CHARACTER(LEN=6) :: message_identifier
5028	CHARACTER(LEN=*) :: routine_name
5029
5030	INTEGER(iwp) :: errno
5031
5032	IF ( nc_stat /= NF90_NOERR ) THEN
5033
5034	WRITE( message_identifier, '(''NC'',I4.4)' ) errno
5035	message_string = TRIM( NF90_STRERROR( nc_stat ) )
5036
5037	CALL message( routine_name, message_identifier, 2, 2, 0, 6, 1 )
5038
5039	ENDIF
5040
5041	#endif
5042	END SUBROUTINE netcdf_handle_error
5043
5044
5045	!------------------------------------------------------------------------------!
5046	! Description:
5047	! ------------
5048	!> Create a dimension in NetCDF file
5049	!------------------------------------------------------------------------------!
5050
5051	SUBROUTINE netcdf_create_dim(ncid, dim_name, ncdim_type, ncdim_id, error_no)
5052
5053	#if defined( __netcdf )
5054
5055	USE kinds
5056
5057	IMPLICIT NONE
5058
5059	CHARACTER(LEN=*), INTENT(IN) :: dim_name
5060
5061	INTEGER, INTENT(IN) :: error_no
5062	INTEGER, INTENT(IN) :: ncid
5063	INTEGER, INTENT(OUT) :: ncdim_id
5064	INTEGER, INTENT(IN) :: ncdim_type
5065
5066	!
5067	!-- Define time coordinate for volume data (unlimited dimension)
5068	nc_stat = NF90_DEF_DIM( ncid, dim_name, ncdim_type, ncdim_id )
5069	CALL netcdf_handle_error( 'netcdf_create_dim', error_no )
5070
5071	#endif
5072
5073	END SUBROUTINE netcdf_create_dim
5074
5075
5076	!------------------------------------------------------------------------------!
5077	! Description:
5078	! ------------
5079	!> Create a one dimensional variable in specific units in NetCDF file
5080	!------------------------------------------------------------------------------!
5081
5082	SUBROUTINE netcdf_create_var( ncid, dim_id, var_name, var_type, var_id, &
5083	unit_name, long_name, error_no1, error_no2, &
5084	error_no3 )
5085
5086	#if defined( __netcdf )
5087	IMPLICIT NONE
5088
5089	CHARACTER(LEN=*), INTENT(IN) :: long_name
5090	CHARACTER(LEN=*), INTENT(IN) :: unit_name
5091	CHARACTER(LEN=*), INTENT(IN) :: var_name
5092
5093	INTEGER, INTENT(IN) :: error_no1
5094	INTEGER, INTENT(IN) :: error_no2
5095	INTEGER, INTENT(IN) :: error_no3
5096	INTEGER, INTENT(IN) :: ncid
5097	INTEGER, INTENT(OUT) :: var_id
5098	INTEGER, INTENT(IN) :: var_type
5099
5100	INTEGER, DIMENSION(:), INTENT(IN) :: dim_id
5101
5102	!
5103	!-- Define variable
5104	nc_stat = NF90_DEF_VAR( ncid, var_name, var_type, dim_id, var_id )
5105	CALL netcdf_handle_error( 'netcdf_create_var', error_no1 )
5106
5107	#if defined( __netcdf4 )
5108	!
5109	!-- Check if variable should be deflate (including shuffling)
5110	!-- and if it is possible (only NetCDF4 with HDF5 supports compression)
5111	IF ( netcdf_data_format > 2 .AND. netcdf_deflate > 0 ) THEN
5112	nc_stat = NF90_DEF_VAR_DEFLATE( ncid, var_id, 1, 1, netcdf_deflate )
5113	CALL netcdf_handle_error( 'netcdf_create_var_deflate', error_no1 )
5114	ENDIF
5115	#endif
5116	!
5117	!-- Set unit name if set
5118	IF ( unit_name /= '' ) THEN
5119	nc_stat = NF90_PUT_ATT( ncid, var_id, 'units', unit_name )
5120	CALL netcdf_handle_error( 'netcdf_create_var', error_no2 )
5121	ENDIF
5122
5123	!
5124	!-- Set long name if set
5125	IF ( long_name /= '' ) THEN
5126	nc_stat = NF90_PUT_ATT( ncid, var_id, 'long_name', long_name )
5127	CALL netcdf_handle_error( 'netcdf_create_var', error_no3 )
5128	ENDIF
5129
5130	#endif
5131	END SUBROUTINE netcdf_create_var
5132
5133	END MODULE netcdf_interface

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |