Home

Context Navigation

netcdf_interface.f90 @ 1842

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

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

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

Context Navigation

source: palm/trunk/SOURCE/netcdf_interface.f90 @ 1842

Download in other formats: