Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |