Home

Context Navigation

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

Last change on this file since 992 was 992, checked in by hoffmann, 12 years ago
removal of two informative messages in netcdf.f90
Property svn:keywords set to `Id`
File size: 184.8 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |