Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |