Home

Context Navigation

source: palm/trunk/SOURCE/netcdf.f90 @ 1586

Last change on this file since 1586 was 1585, checked in by maronga, 9 years ago
Added support for RRTMG radiation code
Property svn:keywords set to `Id`
File size: 211.0 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |