Home

Context Navigation

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

Last change on this file since 3157 was 3049, checked in by Giersch, 7 years ago
Revision history corrected
Property svn:keywords set to `Id`
File size: 248.0 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |