Home

Context Navigation

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

Last change on this file since 984 was 984, checked in by hoffmann, 12 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 185.8 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |