Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |