Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |