Home

Context Navigation

netcdf_interface_mod.f90 @ 1977

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