Home

Context Navigation

netcdf.f90 @ 964

Last change on this file since 964 was 964, checked in by raasch, 13 years ago

old profil-parameters (cross_xtext, cross_normalized_x, etc. ) and respective code removed
(check_open, check_parameters, close_file, data_output_profiles, data_output_spectra, header, modules, parin)

reformatting (netcdf)

append feature removed from unit 14 (check_open)

Property svn:keywords set to Id

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

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

Context Navigation

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

Download in other formats: