Home

Context Navigation

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

Last change on this file since 2295 was 2292, checked in by schwenkel, 8 years ago
implementation of new bulk microphysics scheme
Property svn:keywords set to `Id`
File size: 247.7 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |