Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |