Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |