Home

Context Navigation

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

Last change on this file since 30 was 4, checked in by raasch, 18 years ago
Id keyword set as property for all *.f90 files
Property svn:keywords set to `Id`
File size: 138.0 KB

Line
1	#if defined( __ibmy_special )
2	@PROCESS NOOPTimize
3	#endif
4	SUBROUTINE define_netcdf_header( callmode, extend, av )
5
6	!------------------------------------------------------------------------------!
7	! Current revisions:
8	! ------------------
9	!
10	!
11	! Former revisions:
12	! -----------------
13	! $Id: netcdf.f90 4 2007-02-13 11:33:16Z raasch $
14	! RCS Log replace by Id keyword, revision history cleaned up
15	!
16	! Revision 1.12 2006/09/26 19:35:16 raasch
17	! Bugfix yv coordinates for yz cross sections
18	!
19	! Revision 1.1 2005/05/18 15:37:16 raasch
20	! Initial revision
21	!
22	!
23	! Description:
24	! ------------
25	! In case of extend = .FALSE.:
26	! Define all necessary dimensions, axes and variables for the different
27	! NetCDF datasets. This subroutine is called from check_open after a new
28	! dataset is created. It leaves the open NetCDF files ready to write.
29	!
30	! In case of extend = .TRUE.:
31	! Find out if dimensions and variables of an existing file match the values
32	! of the actual run. If so, get all necessary informations (ids, etc.) from
33	! this file.
34	!
35	! Parameter av can assume values 0 (non-averaged data) and 1 (time averaged
36	! data)
37	!------------------------------------------------------------------------------!
38	#if defined( __netcdf )
39
40	USE arrays_3d
41	USE constants
42	USE control_parameters
43	USE grid_variables
44	USE indices
45	USE netcdf_control
46	USE pegrid
47	USE particle_attributes
48	USE profil_parameter
49	USE spectrum
50	USE statistics
51
52
53	IMPLICIT NONE
54
55	CHARACTER (LEN=2) :: suffix
56	CHARACTER (LEN=2), INTENT (IN) :: callmode
57	CHARACTER (LEN=3) :: suffix1
58	CHARACTER (LEN=4) :: grid_x, grid_y, grid_z
59	CHARACTER (LEN=6) :: mode
60	CHARACTER (LEN=10) :: netcdf_var_name, netcdf_var_name_base, &
61	precision, var
62	CHARACTER (LEN=80) :: time_average_text
63	CHARACTER (LEN=2000) :: var_list, var_list_old
64
65	INTEGER :: av, i, id_x, id_y, id_z, j, ns, ns_old, nz_old
66
67	INTEGER, DIMENSION(1) :: id_dim_time_old, id_dim_x_yz_old, &
68	id_dim_y_xz_old, id_dim_zu_sp_old, &
69	id_dim_zu_xy_old, id_dim_zu_3d_old
70
71	LOGICAL :: found
72
73	LOGICAL, INTENT (INOUT) :: extend
74
75	LOGICAL, SAVE :: init_netcdf = .FALSE.
76
77	REAL, DIMENSION(1) :: last_time_coordinate
78
79	REAL, DIMENSION(:), ALLOCATABLE :: netcdf_data
80
81
82	!
83	!-- Initializing actions (return to calling routine check_parameters afterwards)
84	IF ( .NOT. init_netcdf ) THEN
85	!
86	!-- Check and set accuracy for NetCDF output. First set default value
87	nc_precision = NF90_REAL4
88
89	i = 1
90	DO WHILE ( netcdf_precision(i) /= ' ' )
91	j = INDEX( netcdf_precision(i), '_' )
92	IF ( j == 0 ) THEN
93	IF ( myid == 0 ) THEN
94	PRINT*, '+++ define_netcdf_header: netcdf_precision must ', &
95	'contain a "_" netcdf_precision(', i, ')="', &
96	TRIM( netcdf_precision(i) ),'"'
97	ENDIF
98	CALL local_stop
99	ENDIF
100
101	var = netcdf_precision(i)(1:j-1)
102	precision = netcdf_precision(i)(j+1:)
103
104	IF ( precision == 'NF90_REAL4' ) THEN
105	j = NF90_REAL4
106	ELSEIF ( precision == 'NF90_REAL8' ) THEN
107	j = NF90_REAL8
108	ELSE
109	IF ( myid == 0 ) THEN
110	PRINT*, '+++ define_netcdf_header: illegal netcdf precision: ',&
111	'netcdf_precision(', i, ')="', &
112	TRIM( netcdf_precision(i) ),'"'
113	ENDIF
114	CALL local_stop
115	ENDIF
116
117	SELECT CASE ( var )
118	CASE ( 'xy' )
119	nc_precision(1) = j
120	CASE ( 'xz' )
121	nc_precision(2) = j
122	CASE ( 'yz' )
123	nc_precision(3) = j
124	CASE ( '2d' )
125	nc_precision(1:3) = j
126	CASE ( '3d' )
127	nc_precision(4) = j
128	CASE ( 'pr' )
129	nc_precision(5) = j
130	CASE ( 'ts' )
131	nc_precision(6) = j
132	CASE ( 'sp' )
133	nc_precision(7) = j
134	CASE ( 'prt' )
135	nc_precision(8) = j
136	CASE ( 'all' )
137	nc_precision = j
138
139	CASE DEFAULT
140	IF ( myid == 0 ) THEN
141	PRINT*, '+++ define_netcdf_header: unknown variable in ', &
142	'inipar assignment: netcdf_precision(', i, ')="',&
143	TRIM( netcdf_precision(i) ),'"'
144	ENDIF
145	CALL local_stop
146
147	END SELECT
148
149	i = i + 1
150	IF ( i > 10 ) EXIT
151	ENDDO
152
153	init_netcdf = .TRUE.
154
155	RETURN
156
157	ENDIF
158
159	!
160	!-- Determine the mode to be processed
161	IF ( extend ) THEN
162	mode = callmode // '_ext'
163	ELSE
164	mode = callmode // '_new'
165	ENDIF
166
167	!
168	!-- Select the mode to be processed. Possibilities are xy, xz, yz, pr and ts.
169	SELECT CASE ( mode )
170
171	CASE ( '3d_new' )
172
173	!
174	!-- Define some global attributes of the dataset
175	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'Conventions', &
176	'COARDS' )
177	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 62 )
178
179	IF ( av == 0 ) THEN
180	time_average_text = ' '
181	ELSE
182	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
183	averaging_interval
184	ENDIF
185	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'title', &
186	TRIM( run_description_header ) // &
187	TRIM( time_average_text ) )
188	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 63 )
189	IF ( av == 1 ) THEN
190	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
191	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'time_avg', &
192	TRIM( time_average_text ) )
193	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 63 )
194	ENDIF
195
196	!
197	!-- Define time coordinate for volume data (unlimited dimension)
198	nc_stat = NF90_DEF_DIM( id_set_3d(av), 'time', NF90_UNLIMITED, &
199	id_dim_time_3d(av) )
200	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 64 )
201
202	nc_stat = NF90_DEF_VAR( id_set_3d(av), 'time', NF90_DOUBLE, &
203	id_dim_time_3d(av), id_var_time_3d(av) )
204	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 65 )
205
206	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_time_3d(av), 'units', &
207	'seconds')
208	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 66 )
209
210	!
211	!-- Define spatial dimensions and coordinates:
212	!-- Define vertical coordinate grid (zu grid)
213	nc_stat = NF90_DEF_DIM( id_set_3d(av), 'zu_3d', nz_do3d-nzb+1, &
214	id_dim_zu_3d(av) )
215	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 67 )
216
217	nc_stat = NF90_DEF_VAR( id_set_3d(av), 'zu_3d', NF90_DOUBLE, &
218	id_dim_zu_3d(av), id_var_zu_3d(av) )
219	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 68 )
220
221	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_zu_3d(av), 'units', &
222	'meters' )
223	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 69 )
224
225	!
226	!-- Define vertical coordinate grid (zw grid)
227	nc_stat = NF90_DEF_DIM( id_set_3d(av), 'zw_3d', nz_do3d-nzb+1, &
228	id_dim_zw_3d(av) )
229	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 70 )
230
231	nc_stat = NF90_DEF_VAR( id_set_3d(av), 'zw_3d', NF90_DOUBLE, &
232	id_dim_zw_3d(av), id_var_zw_3d(av) )
233	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 71 )
234
235	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_zw_3d(av), 'units', &
236	'meters' )
237	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 72 )
238
239	!
240	!-- Define x-axis (for scalar position)
241	nc_stat = NF90_DEF_DIM( id_set_3d(av), 'x', nx+2, id_dim_x_3d(av) )
242	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 73 )
243
244	nc_stat = NF90_DEF_VAR( id_set_3d(av), 'x', NF90_DOUBLE, &
245	id_dim_x_3d(av), id_var_x_3d(av) )
246	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 74 )
247
248	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_x_3d(av), 'units', &
249	'meters' )
250	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 75 )
251
252	!
253	!-- Define x-axis (for u position)
254	nc_stat = NF90_DEF_DIM( id_set_3d(av), 'xu', nx+2, id_dim_xu_3d(av) )
255	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 358 )
256
257	nc_stat = NF90_DEF_VAR( id_set_3d(av), 'xu', NF90_DOUBLE, &
258	id_dim_xu_3d(av), id_var_xu_3d(av) )
259	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 359 )
260
261	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_xu_3d(av), 'units', &
262	'meters' )
263	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 360 )
264
265	!
266	!-- Define y-axis (for scalar position)
267	nc_stat = NF90_DEF_DIM( id_set_3d(av), 'y', ny+2, id_dim_y_3d(av) )
268	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 76 )
269
270	nc_stat = NF90_DEF_VAR( id_set_3d(av), 'y', NF90_DOUBLE, &
271	id_dim_y_3d(av), id_var_y_3d(av) )
272	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 77 )
273
274	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_y_3d(av), 'units', &
275	'meters' )
276	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 78 )
277
278	!
279	!-- Define y-axis (for v position)
280	nc_stat = NF90_DEF_DIM( id_set_3d(av), 'yv', ny+2, id_dim_yv_3d(av) )
281	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 361 )
282
283	nc_stat = NF90_DEF_VAR( id_set_3d(av), 'yv', NF90_DOUBLE, &
284	id_dim_yv_3d(av), id_var_yv_3d(av) )
285	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 362 )
286
287	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_yv_3d(av), 'units', &
288	'meters' )
289	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 363 )
290
291	!
292	!-- Define the variables
293	var_list = ';'
294	i = 1
295
296	DO WHILE ( do3d(av,i)(1:1) /= ' ' )
297
298	!
299	!-- Check for the grid
300	found = .TRUE.
301	SELECT CASE ( do3d(av,i) )
302	!
303	!-- Most variables are defined on the scalar grid
304	CASE ( 'e', 'p', 'pc', 'pr', 'pt', 'q', 'ql', 'ql_c', 'ql_v', &
305	'ql_vp', 'qv', 's', 'vpt' )
306
307	grid_x = 'x'
308	grid_y = 'y'
309	grid_z = 'zu'
310	!
311	!-- u grid
312	CASE ( 'u' )
313
314	grid_x = 'xu'
315	grid_y = 'y'
316	grid_z = 'zu'
317	!
318	!-- v grid
319	CASE ( 'v' )
320
321	grid_x = 'x'
322	grid_y = 'yv'
323	grid_z = 'zu'
324	!
325	!-- w grid
326	CASE ( 'w' )
327
328	grid_x = 'x'
329	grid_y = 'y'
330	grid_z = 'zw'
331
332	CASE DEFAULT
333	!
334	!-- Check for user-defined quantities
335	CALL user_define_netcdf_grid( do3d(av,i), found, grid_x, &
336	grid_y, grid_z )
337
338	END SELECT
339
340	!
341	!-- Select the respective dimension ids
342	IF ( grid_x == 'x' ) THEN
343	id_x = id_dim_x_3d(av)
344	ELSEIF ( grid_x == 'xu' ) THEN
345	id_x = id_dim_xu_3d(av)
346	ENDIF
347
348	IF ( grid_y == 'y' ) THEN
349	id_y = id_dim_y_3d(av)
350	ELSEIF ( grid_y == 'yv' ) THEN
351	id_y = id_dim_yv_3d(av)
352	ENDIF
353
354	IF ( grid_z == 'zu' ) THEN
355	id_z = id_dim_zu_3d(av)
356	ELSEIF ( grid_z == 'zw' ) THEN
357	id_z = id_dim_zw_3d(av)
358	ENDIF
359
360	!
361	!-- Define the grid
362	nc_stat = NF90_DEF_VAR( id_set_3d(av), do3d(av,i), &
363	nc_precision(4), &
364	(/ id_x, id_y, id_z, id_dim_time_3d(av) /), &
365	id_var_do3d(av,i) )
366
367	IF ( .NOT. found ) THEN
368	PRINT*, '+++ define_netcdf_header: no grid defined for', &
369	' variable ', do3d(av,i)
370	ENDIF
371
372	var_list = TRIM( var_list ) // TRIM( do3d(av,i) ) // ';'
373
374	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 79 )
375	!
376	!-- Store the 'real' name of the variable (with *, for example)
377	!-- in the long_name attribute. This is evaluated by Ferret,
378	!-- for example.
379	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_do3d(av,i), &
380	'long_name', do3d(av,i) )
381	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 80 )
382	!
383	!-- Define the variable's unit
384	nc_stat = NF90_PUT_ATT( id_set_3d(av), id_var_do3d(av,i), &
385	'units', TRIM( do3d_unit(av,i) ) )
386	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 357 )
387
388	i = i + 1
389
390	ENDDO
391
392	!
393	!-- No arrays to output
394	IF ( i == 1 ) RETURN
395	! nc_stat = NF90_CLOSE( id_set_3d(av) )
396	! IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 9999 )
397	! RETURN
398	! ENDIF
399
400	!
401	!-- Write the list of variables as global attribute (this is used by
402	!-- restart runs and by combine_plot_fields)
403	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'VAR_LIST', &
404	var_list )
405	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 81 )
406
407	!
408	!-- Leave NetCDF define mode
409	nc_stat = NF90_ENDDEF( id_set_3d(av) )
410	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 82 )
411
412	!
413	!-- Write data for x and xu axis (shifted by -dx/2)
414	ALLOCATE( netcdf_data(0:nx+1) )
415
416	DO i = 0, nx+1
417	netcdf_data(i) = i * dx
418	ENDDO
419
420	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_x_3d(av), netcdf_data, &
421	start = (/ 1 /), count = (/ nx+2 /) )
422	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 83 )
423
424	DO i = 0, nx+1
425	netcdf_data(i) = ( i - 0.5 ) * dx
426	ENDDO
427
428	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_xu_3d(av), &
429	netcdf_data, start = (/ 1 /), &
430	count = (/ nx+2 /) )
431	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 385 )
432
433	DEALLOCATE( netcdf_data )
434
435	!
436	!-- Write data for y and yv axis (shifted by -dy/2)
437	ALLOCATE( netcdf_data(0:ny+1) )
438
439	DO i = 0, ny+1
440	netcdf_data(i) = i * dy
441	ENDDO
442
443	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_y_3d(av), netcdf_data, &
444	start = (/ 1 /), count = (/ ny+2 /))
445	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 84 )
446
447	DO i = 0, ny+1
448	netcdf_data(i) = ( i - 0.5 ) * dy
449	ENDDO
450
451	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_yv_3d(av), &
452	netcdf_data, start = (/ 1 /), &
453	count = (/ ny+2 /))
454	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 387 )
455
456	DEALLOCATE( netcdf_data )
457
458	!
459	!-- Write zu and zw data (vertical axes)
460	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_zu_3d(av), &
461	zu(nzb:nz_do3d), start = (/ 1 /), &
462	count = (/ nz_do3d-nzb+1 /) )
463	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 85 )
464
465	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_zw_3d(av), &
466	zw(nzb:nz_do3d), start = (/ 1 /), &
467	count = (/ nz_do3d-nzb+1 /) )
468	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 86 )
469
470
471	CASE ( '3d_ext' )
472
473	!
474	!-- Get the list of variables and compare with the actual run.
475	!-- First var_list_old has to be reset, since GET_ATT does not assign
476	!-- trailing blanks.
477	var_list_old = ' '
478	nc_stat = NF90_GET_ATT( id_set_3d(av), NF90_GLOBAL, 'VAR_LIST', &
479	var_list_old )
480	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 87 )
481
482	var_list = ';'
483	i = 1
484	DO WHILE ( do3d(av,i)(1:1) /= ' ' )
485	var_list = TRIM(var_list) // TRIM( do3d(av,i) ) // ';'
486	i = i + 1
487	ENDDO
488
489	IF ( av == 0 ) THEN
490	var = '(3d)'
491	ELSE
492	var = '(3d_av)'
493	ENDIF
494
495	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
496	PRINT*, '+++ WARNING: NetCDF file for volume data ' // &
497	TRIM( var ) // ' from previuos run found,'
498	PRINT*, ' but this file cannot be extended due to' // &
499	' variable mismatch.'
500	PRINT*, ' New file is created instead.'
501	extend = .FALSE.
502	RETURN
503	ENDIF
504
505	!
506	!-- Get and compare the number of vertical gridpoints
507	nc_stat = NF90_INQ_VARID( id_set_3d(av), 'zu_3d', id_var_zu_3d(av) )
508	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 88 )
509
510	nc_stat = NF90_INQUIRE_VARIABLE( id_set_3d(av), id_var_zu_3d(av), &
511	dimids = id_dim_zu_3d_old )
512	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 89 )
513	id_dim_zu_3d(av) = id_dim_zu_3d_old(1)
514
515	nc_stat = NF90_INQUIRE_DIMENSION( id_set_3d(av), id_dim_zu_3d(av), &
516	len = nz_old )
517	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 90 )
518
519	IF ( nz_do3d-nzb+1 /= nz_old ) THEN
520	PRINT*, '+++ WARNING: NetCDF file for volume data ' // &
521	TRIM( var ) // ' from previuos run found,'
522	PRINT*, ' but this file cannot be extended due to' // &
523	' mismatch in number of'
524	PRINT*, ' vertical grid points (nz_do3d).'
525	PRINT*, ' New file is created instead.'
526	extend = .FALSE.
527	RETURN
528	ENDIF
529
530	!
531	!-- Get the id of the time coordinate (unlimited coordinate) and its
532	!-- last index on the file. The next time level is pl3d..count+1.
533	!-- The current time must be larger than the last output time
534	!-- on the file.
535	nc_stat = NF90_INQ_VARID( id_set_3d(av), 'time', id_var_time_3d(av) )
536	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 91 )
537
538	nc_stat = NF90_INQUIRE_VARIABLE( id_set_3d(av), id_var_time_3d(av), &
539	dimids = id_dim_time_old )
540	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 92 )
541	id_dim_time_3d(av) = id_dim_time_old(1)
542
543	nc_stat = NF90_INQUIRE_DIMENSION( id_set_3d(av), id_dim_time_3d(av), &
544	len = do3d_time_count(av) )
545	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 93 )
546
547	nc_stat = NF90_GET_VAR( id_set_3d(av), id_var_time_3d(av), &
548	last_time_coordinate, &
549	start = (/ do3d_time_count(av) /), &
550	count = (/ 1 /) )
551	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 94 )
552
553	IF ( last_time_coordinate(1) >= simulated_time ) THEN
554	PRINT*, '+++ WARNING: NetCDF file for volume data ' // &
555	TRIM( var ) // ' from previuos run found,'
556	PRINT*, ' but this file cannot be extended becaus' // &
557	'e the current output time'
558	PRINT*, ' is less or equal than the last output t' // &
559	'ime on this file.'
560	PRINT*, ' New file is created instead.'
561	do3d_time_count(av) = 0
562	extend = .FALSE.
563	RETURN
564	ENDIF
565
566	!
567	!-- Dataset seems to be extendable.
568	!-- Now get the variable ids.
569	i = 1
570	DO WHILE ( do3d(av,i)(1:1) /= ' ' )
571	nc_stat = NF90_INQ_VARID( id_set_3d(av), TRIM( do3d(av,i) ), &
572	id_var_do3d(av,i) )
573	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 95 )
574	i = i + 1
575	ENDDO
576
577	!
578	!-- Change the titel attribute on file
579	nc_stat = NF90_PUT_ATT( id_set_3d(av), NF90_GLOBAL, 'title', &
580	TRIM( run_description_header ) )
581	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 96 )
582	PRINT, '** NetCDF file for volume data ' // TRIM( var ) // &
583	' from previous run found.'
584	PRINT*, ' This file will be extended.'
585
586
587	CASE ( 'xy_new' )
588
589	!
590	!-- Define some global attributes of the dataset
591	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'Conventions', &
592	'COARDS' )
593	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 97 )
594
595	IF ( av == 0 ) THEN
596	time_average_text = ' '
597	ELSE
598	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
599	averaging_interval
600	ENDIF
601	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'title', &
602	TRIM( run_description_header ) // &
603	TRIM( time_average_text ) )
604	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 98 )
605	IF ( av == 1 ) THEN
606	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
607	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'time_avg', &
608	TRIM( time_average_text ) )
609	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 98 )
610	ENDIF
611
612	!
613	!-- Define time coordinate for xy sections (unlimited dimension)
614	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'time', NF90_UNLIMITED, &
615	id_dim_time_xy(av) )
616	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 99 )
617
618	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'time', NF90_DOUBLE, &
619	id_dim_time_xy(av), id_var_time_xy(av) )
620	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 100 )
621
622	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_time_xy(av), 'units', &
623	'seconds')
624	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 101 )
625
626	!
627	!-- Define the spatial dimensions and coordinates for xy-sections.
628	!-- First, determine the number of horizontal sections to be written.
629	IF ( section(1,1) == -9999 ) THEN
630	RETURN
631	ELSE
632	ns = 1
633	DO WHILE ( section(ns,1) /= -9999 .AND. ns <= 100 )
634	ns = ns + 1
635	ENDDO
636	ns = ns - 1
637	ENDIF
638
639	!
640	!-- Define vertical coordinate grid (zu grid)
641	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'zu_xy', ns, id_dim_zu_xy(av) )
642	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 102 )
643
644	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'zu_xy', NF90_DOUBLE, &
645	id_dim_zu_xy(av), id_var_zu_xy(av) )
646	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 103 )
647
648	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_zu_xy(av), 'units', &
649	'meters' )
650	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 104 )
651
652	!
653	!-- Define vertical coordinate grid (zw grid)
654	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'zw_xy', ns, id_dim_zw_xy(av) )
655	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 105 )
656
657	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'zw_xy', NF90_DOUBLE, &
658	id_dim_zw_xy(av), id_var_zw_xy(av) )
659	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 106 )
660
661	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_zw_xy(av), 'units', &
662	'meters' )
663	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 107 )
664
665	!
666	!-- Define a pseudo vertical coordinate grid for the surface variables
667	!-- u* and t* to store their height level
668	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'zu1_xy', 1, &
669	id_dim_zu1_xy(av) )
670	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 108 )
671
672	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'zu1_xy', NF90_DOUBLE, &
673	id_dim_zu1_xy(av), id_var_zu1_xy(av) )
674	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 109 )
675
676	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_zu1_xy(av), 'units', &
677	'meters' )
678	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 110 )
679
680	!
681	!-- Define a variable to store the layer indices of the horizontal cross
682	!-- sections, too
683	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'ind_z_xy', NF90_DOUBLE, &
684	id_dim_zu_xy(av), id_var_ind_z_xy(av) )
685	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 111 )
686
687	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_ind_z_xy(av), 'units', &
688	'gridpoints')
689	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 112 )
690
691	!
692	!-- Define x-axis (for scalar position)
693	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'x', nx+2, id_dim_x_xy(av) )
694	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 113 )
695
696	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'x', NF90_DOUBLE, &
697	id_dim_x_xy(av), id_var_x_xy(av) )
698	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 114 )
699
700	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_x_xy(av), 'units', &
701	'meters' )
702	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 115 )
703
704	!
705	!-- Define x-axis (for u position)
706	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'xu', nx+2, id_dim_xu_xy(av) )
707	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 388 )
708
709	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'xu', NF90_DOUBLE, &
710	id_dim_xu_xy(av), id_var_xu_xy(av) )
711	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 389 )
712
713	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_xu_xy(av), 'units', &
714	'meters' )
715	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 390 )
716
717	!
718	!-- Define y-axis (for scalar position)
719	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'y', ny+2, id_dim_y_xy(av) )
720	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 116 )
721
722	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'y', NF90_DOUBLE, &
723	id_dim_y_xy(av), id_var_y_xy(av) )
724	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 117 )
725
726	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_y_xy(av), 'units', &
727	'meters' )
728	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 118 )
729
730	!
731	!-- Define y-axis (for scalar position)
732	nc_stat = NF90_DEF_DIM( id_set_xy(av), 'yv', ny+2, id_dim_yv_xy(av) )
733	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 364 )
734
735	nc_stat = NF90_DEF_VAR( id_set_xy(av), 'yv', NF90_DOUBLE, &
736	id_dim_yv_xy(av), id_var_yv_xy(av) )
737	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 365 )
738
739	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_yv_xy(av), 'units', &
740	'meters' )
741	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 366 )
742
743	!
744	!-- Define the variables
745	var_list = ';'
746	i = 1
747
748	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
749
750	IF ( INDEX( do2d(av,i), 'xy' ) /= 0 ) THEN
751	!
752	!-- If there is a star in the variable name (u* or t*), it is a
753	!-- surface variable. Define it with id_dim_zu1_xy.
754	IF ( INDEX( do2d(av,i), '*' ) /= 0 ) THEN
755	!
756	!-- First, remove those characters not allowed by NetCDF
757	netcdf_var_name = do2d(av,i)
758	CALL clean_netcdf_varname( netcdf_var_name )
759
760	nc_stat = NF90_DEF_VAR( id_set_xy(av), netcdf_var_name, &
761	nc_precision(1), &
762	(/ id_dim_x_xy(av), id_dim_y_xy(av),&
763	id_dim_zu1_xy(av), &
764	id_dim_time_xy(av) /), &
765	id_var_do2d(av,i) )
766
767	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
768
769	ELSE
770
771	!
772	!-- Check for the grid
773	found = .TRUE.
774	SELECT CASE ( do2d(av,i) )
775	!
776	!-- Most variables are defined on the zu grid
777	CASE ( 'e_xy', 'p_xy', 'pc_xy', 'pr_xy', 'pt_xy', 'q_xy',&
778	'ql_xy', 'ql_c_xy', 'ql_v_xy', 'ql_vp_xy', &
779	'qv_xy', 's_xy', 'vpt_xy' )
780
781	grid_x = 'x'
782	grid_y = 'y'
783	grid_z = 'zu'
784	!
785	!-- u grid
786	CASE ( 'u_xy' )
787
788	grid_x = 'xu'
789	grid_y = 'y'
790	grid_z = 'zu'
791	!
792	!-- v grid
793	CASE ( 'v_xy' )
794
795	grid_x = 'x'
796	grid_y = 'yv'
797	grid_z = 'zu'
798	!
799	!-- w grid
800	CASE ( 'w_xy' )
801
802	grid_x = 'x'
803	grid_y = 'y'
804	grid_z = 'zw'
805
806	CASE DEFAULT
807	!
808	!-- Check for user-defined quantities
809	CALL user_define_netcdf_grid( do2d(av,i), found, &
810	grid_x, grid_y, grid_z )
811
812	END SELECT
813
814	!
815	!-- Select the respective dimension ids
816	IF ( grid_x == 'x' ) THEN
817	id_x = id_dim_x_xy(av)
818	ELSEIF ( grid_x == 'xu' ) THEN
819	id_x = id_dim_xu_xy(av)
820	ENDIF
821
822	IF ( grid_y == 'y' ) THEN
823	id_y = id_dim_y_xy(av)
824	ELSEIF ( grid_y == 'yv' ) THEN
825	id_y = id_dim_yv_xy(av)
826	ENDIF
827
828	IF ( grid_z == 'zu' ) THEN
829	id_z = id_dim_zu_xy(av)
830	ELSEIF ( grid_z == 'zw' ) THEN
831	id_z = id_dim_zw_xy(av)
832	ENDIF
833
834	!
835	!-- Define the grid
836	nc_stat = NF90_DEF_VAR( id_set_xy(av), do2d(av,i), &
837	nc_precision(1), &
838	(/ id_x, id_y, id_z, id_dim_time_xy(av) /), &
839	id_var_do2d(av,i) )
840
841	IF ( .NOT. found ) THEN
842	PRINT*, '+++ define_netcdf_header: no grid defined ', &
843	'for variable ', do2d(av,i)
844	ENDIF
845
846	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
847
848	ENDIF
849
850	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 119 )
851	!
852	!-- Store the 'real' name of the variable (with *, for example)
853	!-- in the long_name attribute. This is evaluated by Ferret,
854	!-- for example.
855	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_do2d(av,i), &
856	'long_name', do2d(av,i) )
857	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 120 )
858	!
859	!-- Define the variable's unit
860	nc_stat = NF90_PUT_ATT( id_set_xy(av), id_var_do2d(av,i), &
861	'units', TRIM( do2d_unit(av,i) ) )
862	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 354 )
863	ENDIF
864
865	i = i + 1
866
867	ENDDO
868
869	!
870	!-- No arrays to output. Close the netcdf file and return.
871	IF ( i == 1 ) RETURN
872	! nc_stat = NF90_CLOSE( id_set_xy(av) )
873	! IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 9999 )
874	! RETURN
875	! ENDIF
876
877	!
878	!-- Write the list of variables as global attribute (this is used by
879	!-- restart runs and by combine_plot_fields)
880	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'VAR_LIST', &
881	var_list )
882	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 121 )
883
884	!
885	!-- Leave NetCDF define mode
886	nc_stat = NF90_ENDDEF( id_set_xy(av) )
887	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 122 )
888
889	!
890	!-- Write axis data: z_xy, x, y
891	ALLOCATE( netcdf_data(1:ns) )
892
893	!
894	!-- Write zu data
895	DO i = 1, ns
896	IF( section(i,1) == -1 ) THEN
897	netcdf_data(i) = -1.0 ! section averaged along z
898	ELSE
899	netcdf_data(i) = zu( section(i,1) )
900	ENDIF
901	ENDDO
902	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zu_xy(av), &
903	netcdf_data, start = (/ 1 /), &
904	count = (/ ns /) )
905	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 123 )
906
907	!
908	!-- Write zw data
909	DO i = 1, ns
910	IF( section(i,1) == -1 ) THEN
911	netcdf_data(i) = -1.0 ! section averaged along z
912	ELSE
913	netcdf_data(i) = zw( section(i,1) )
914	ENDIF
915	ENDDO
916	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zw_xy(av), &
917	netcdf_data, start = (/ 1 /), &
918	count = (/ ns /) )
919	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 124 )
920
921	!
922	!-- Write gridpoint number data
923	netcdf_data(1:ns) = section(1:ns,1)
924	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_ind_z_xy(av), &
925	netcdf_data, start = (/ 1 /), &
926	count = (/ ns /) )
927	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 125 )
928
929	DEALLOCATE( netcdf_data )
930
931	!
932	!-- Write the cross section height u, t
933	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_zu1_xy(av), &
934	(/ zu(nzb+1) /), start = (/ 1 /), &
935	count = (/ 1 /) )
936	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 126 )
937
938	!
939	!-- Write data for x and xu axis (shifted by -dx/2)
940	ALLOCATE( netcdf_data(0:nx+1) )
941
942	DO i = 0, nx+1
943	netcdf_data(i) = i * dx
944	ENDDO
945
946	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_x_xy(av), netcdf_data, &
947	start = (/ 1 /), count = (/ nx+2 /) )
948	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 127 )
949
950	DO i = 0, nx+1
951	netcdf_data(i) = ( i - 0.5 ) * dx
952	ENDDO
953
954	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_xu_xy(av), &
955	netcdf_data, start = (/ 1 /), &
956	count = (/ nx+2 /) )
957	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 367 )
958
959	DEALLOCATE( netcdf_data )
960
961	!
962	!-- Write data for y and yv axis (shifted by -dy/2)
963	ALLOCATE( netcdf_data(0:ny+1) )
964
965	DO i = 0, ny+1
966	netcdf_data(i) = i * dy
967	ENDDO
968
969	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_y_xy(av), netcdf_data, &
970	start = (/ 1 /), count = (/ ny+2 /))
971	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 128 )
972
973	DO i = 0, ny+1
974	netcdf_data(i) = ( i - 0.5 ) * dy
975	ENDDO
976
977	nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_yv_xy(av), &
978	netcdf_data, start = (/ 1 /), &
979	count = (/ ny+2 /))
980	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 368 )
981
982	DEALLOCATE( netcdf_data )
983
984
985	CASE ( 'xy_ext' )
986
987	!
988	!-- Get the list of variables and compare with the actual run.
989	!-- First var_list_old has to be reset, since GET_ATT does not assign
990	!-- trailing blanks.
991	var_list_old = ' '
992	nc_stat = NF90_GET_ATT( id_set_xy(av), NF90_GLOBAL, 'VAR_LIST', &
993	var_list_old )
994	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 129 )
995
996	var_list = ';'
997	i = 1
998	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
999	IF ( INDEX( do2d(av,i), 'xy' ) /= 0 ) THEN
1000	netcdf_var_name = do2d(av,i)
1001	CALL clean_netcdf_varname( netcdf_var_name )
1002	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
1003	ENDIF
1004	i = i + 1
1005	ENDDO
1006
1007	IF ( av == 0 ) THEN
1008	var = '(xy)'
1009	ELSE
1010	var = '(xy_av)'
1011	ENDIF
1012
1013	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
1014	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1015	TRIM( var ) // ' from previuos run found,'
1016	PRINT*, ' but this file cannot be extended due to' // &
1017	' variable mismatch.'
1018	PRINT*, ' New file is created instead.'
1019	extend = .FALSE.
1020	RETURN
1021	ENDIF
1022
1023	!
1024	!-- Calculate the number of current sections
1025	ns = 1
1026	DO WHILE ( section(ns,1) /= -9999 .AND. ns <= 100 )
1027	ns = ns + 1
1028	ENDDO
1029	ns = ns - 1
1030
1031	!
1032	!-- Get and compare the number of horizontal cross sections
1033	nc_stat = NF90_INQ_VARID( id_set_xy(av), 'zu_xy', id_var_zu_xy(av) )
1034	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 130 )
1035
1036	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xy(av), id_var_zu_xy(av), &
1037	dimids = id_dim_zu_xy_old )
1038	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 131 )
1039	id_dim_zu_xy(av) = id_dim_zu_xy_old(1)
1040
1041	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xy(av), id_dim_zu_xy(av), &
1042	len = ns_old )
1043	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 132 )
1044
1045	IF ( ns /= ns_old ) THEN
1046	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1047	TRIM( var ) // ' from previuos run found,'
1048	PRINT*, ' but this file cannot be extended due to' // &
1049	' mismatch in number of'
1050	PRINT*, ' cross sections.'
1051	PRINT*, ' New file is created instead.'
1052	extend = .FALSE.
1053	RETURN
1054	ENDIF
1055
1056	!
1057	!-- Get and compare the heights of the cross sections
1058	ALLOCATE( netcdf_data(1:ns_old) )
1059
1060	nc_stat = NF90_GET_VAR( id_set_xy(av), id_var_zu_xy(av), netcdf_data )
1061	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 133 )
1062
1063	DO i = 1, ns
1064	IF ( section(i,1) /= -1 ) THEN
1065	IF ( zu(section(i,1)) /= netcdf_data(i) ) THEN
1066	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1067	TRIM( var ) // ' from previuos run found,'
1068	PRINT*, ' but this file cannot be extended' // &
1069	' due to mismatch in cross'
1070	PRINT*, ' section levels.'
1071	PRINT*, ' New file is created instead.'
1072	extend = .FALSE.
1073	RETURN
1074	ENDIF
1075	ELSE
1076	IF ( -1.0 /= netcdf_data(i) ) THEN
1077	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1078	TRIM( var ) // ' from previuos run found,'
1079	PRINT*, ' but this file cannot be extended' // &
1080	' due to mismatch in cross'
1081	PRINT*, ' section levels.'
1082	PRINT*, ' New file is created instead.'
1083	extend = .FALSE.
1084	RETURN
1085	ENDIF
1086	ENDIF
1087	ENDDO
1088
1089	DEALLOCATE( netcdf_data )
1090
1091	!
1092	!-- Get the id of the time coordinate (unlimited coordinate) and its
1093	!-- last index on the file. The next time level is do2d..count+1.
1094	!-- The current time must be larger than the last output time
1095	!-- on the file.
1096	nc_stat = NF90_INQ_VARID( id_set_xy(av), 'time', id_var_time_xy(av) )
1097	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 134 )
1098
1099	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xy(av), id_var_time_xy(av), &
1100	dimids = id_dim_time_old )
1101	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 135 )
1102	id_dim_time_xy(av) = id_dim_time_old(1)
1103
1104	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xy(av), id_dim_time_xy(av), &
1105	len = do2d_xy_time_count(av) )
1106	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 136 )
1107
1108	nc_stat = NF90_GET_VAR( id_set_xy(av), id_var_time_xy(av), &
1109	last_time_coordinate, &
1110	start = (/ do2d_xy_time_count(av) /), &
1111	count = (/ 1 /) )
1112	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 137 )
1113
1114	IF ( last_time_coordinate(1) >= simulated_time ) THEN
1115	PRINT*, '+++ WARNING: NetCDF file for cross sections ' // &
1116	TRIM( var ) // ' from previuos run found,'
1117	PRINT*, ' but this file cannot be extended becaus' // &
1118	'e the current output time'
1119	PRINT*, ' is less or equal than the last output t' // &
1120	'ime on this file.'
1121	PRINT*, ' New file is created instead.'
1122	do2d_xy_time_count(av) = 0
1123	extend = .FALSE.
1124	RETURN
1125	ENDIF
1126
1127	!
1128	!-- Dataset seems to be extendable.
1129	!-- Now get the variable ids.
1130	i = 1
1131	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
1132	IF ( INDEX( do2d(av,i), 'xy' ) /= 0 ) THEN
1133	netcdf_var_name = do2d(av,i)
1134	CALL clean_netcdf_varname( netcdf_var_name )
1135	nc_stat = NF90_INQ_VARID( id_set_xy(av), netcdf_var_name, &
1136	id_var_do2d(av,i) )
1137	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 138 )
1138	ENDIF
1139	i = i + 1
1140	ENDDO
1141
1142	!
1143	!-- Change the titel attribute on file
1144	nc_stat = NF90_PUT_ATT( id_set_xy(av), NF90_GLOBAL, 'title', &
1145	TRIM( run_description_header ) )
1146	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 139 )
1147	PRINT, '** NetCDF file for cross-sections ' // TRIM( var ) // &
1148	' from previous run found.'
1149	PRINT*, ' This file will be extended.'
1150
1151
1152	CASE ( 'xz_new' )
1153
1154	!
1155	!-- Define some global attributes of the dataset
1156	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'Conventions', &
1157	'COARDS' )
1158	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 140 )
1159
1160	IF ( av == 0 ) THEN
1161	time_average_text = ' '
1162	ELSE
1163	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
1164	averaging_interval
1165	ENDIF
1166	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'title', &
1167	TRIM( run_description_header ) // &
1168	TRIM( time_average_text ) )
1169	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 141 )
1170	IF ( av == 1 ) THEN
1171	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
1172	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'time_avg', &
1173	TRIM( time_average_text ) )
1174	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 141 )
1175	ENDIF
1176
1177	!
1178	!-- Define time coordinate for xz sections (unlimited dimension)
1179	nc_stat = NF90_DEF_DIM( id_set_xz(av), 'time', NF90_UNLIMITED, &
1180	id_dim_time_xz(av) )
1181	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 142 )
1182
1183	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'time', NF90_DOUBLE, &
1184	id_dim_time_xz(av), id_var_time_xz(av) )
1185	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 143 )
1186
1187	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_time_xz(av), 'units', &
1188	'seconds')
1189	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 144 )
1190
1191	!
1192	!-- Define the spatial dimensions and coordinates for xz-sections.
1193	!-- First, determine the number of vertical sections to be written.
1194	IF ( section(1,2) == -9999 ) THEN
1195	RETURN
1196	ELSE
1197	ns = 1
1198	DO WHILE ( section(ns,2) /= -9999 .AND. ns <= 100 )
1199	ns = ns + 1
1200	ENDDO
1201	ns = ns - 1
1202	ENDIF
1203
1204	!
1205	!-- Define y-axis (for scalar position)
1206	nc_stat = NF90_DEF_DIM( id_set_xz(av), 'y_xz', ns, id_dim_y_xz(av) )
1207	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 145 )
1208
1209	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'y_xz', NF90_DOUBLE, &
1210	id_dim_y_xz(av), id_var_y_xz(av) )
1211	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 146 )
1212
1213	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_y_xz(av), 'units', &
1214	'meters' )
1215	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 147 )
1216
1217	!
1218	!-- Define y-axis (for v position)
1219	nc_stat = NF90_DEF_DIM( id_set_xz(av), 'yv_xz', ns, id_dim_yv_xz(av) )
1220	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 369 )
1221
1222	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'yv_xz', NF90_DOUBLE, &
1223	id_dim_yv_xz(av), id_var_yv_xz(av) )
1224	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 370 )
1225
1226	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_yv_xz(av), 'units', &
1227	'meters' )
1228	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 371 )
1229
1230	!
1231	!-- Define a variable to store the layer indices of the vertical cross
1232	!-- sections
1233	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'ind_y_xz', NF90_DOUBLE, &
1234	id_dim_y_xz(av), id_var_ind_y_xz(av) )
1235	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 148 )
1236
1237	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_ind_y_xz(av), 'units', &
1238	'gridpoints')
1239	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 149 )
1240
1241	!
1242	!-- Define x-axis (for scalar position)
1243	nc_stat = NF90_DEF_DIM( id_set_xz(av), 'x', nx+2, id_dim_x_xz(av) )
1244	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 150 )
1245
1246	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'x', NF90_DOUBLE, &
1247	id_dim_x_xz(av), id_var_x_xz(av) )
1248	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 151 )
1249
1250	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_x_xz(av), 'units', &
1251	'meters' )
1252	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 152 )
1253
1254	!
1255	!-- Define x-axis (for u position)
1256	nc_stat = NF90_DEF_DIM( id_set_xz(av), 'xu', nx+2, id_dim_xu_xz(av) )
1257	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 372 )
1258
1259	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'xu', NF90_DOUBLE, &
1260	id_dim_xu_xz(av), id_var_xu_xz(av) )
1261	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 373 )
1262
1263	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_xu_xz(av), 'units', &
1264	'meters' )
1265	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 374 )
1266
1267	!
1268	!-- Define the two z-axes (zu and zw)
1269	nc_stat = NF90_DEF_DIM( id_set_xz(av), 'zu', nz+2, id_dim_zu_xz(av) )
1270	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 153 )
1271
1272	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'zu', NF90_DOUBLE, &
1273	id_dim_zu_xz(av), id_var_zu_xz(av) )
1274	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 154 )
1275
1276	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_zu_xz(av), 'units', &
1277	'meters' )
1278	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 155 )
1279
1280	nc_stat = NF90_DEF_DIM( id_set_xz(av), 'zw', nz+2, id_dim_zw_xz(av) )
1281	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 156 )
1282
1283	nc_stat = NF90_DEF_VAR( id_set_xz(av), 'zw', NF90_DOUBLE, &
1284	id_dim_zw_xz(av), id_var_zw_xz(av) )
1285	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 157 )
1286
1287	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_zw_xz(av), 'units', &
1288	'meters' )
1289	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 158 )
1290
1291
1292	!
1293	!-- Define the variables
1294	var_list = ';'
1295	i = 1
1296
1297	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
1298
1299	IF ( INDEX( do2d(av,i), 'xz' ) /= 0 ) THEN
1300
1301	!
1302	!-- Check for the grid
1303	found = .TRUE.
1304	SELECT CASE ( do2d(av,i) )
1305	!
1306	!-- Most variables are defined on the zu grid
1307	CASE ( 'e_xz', 'p_xz', 'pc_xz', 'pr_xz', 'pt_xz', 'q_xz', &
1308	'ql_xz', 'ql_c_xz', 'ql_v_xz', 'ql_vp_xz', 'qv_xz', &
1309	's_xz', 'vpt_xz' )
1310
1311	grid_x = 'x'
1312	grid_y = 'y'
1313	grid_z = 'zu'
1314	!
1315	!-- u grid
1316	CASE ( 'u_xz' )
1317
1318	grid_x = 'xu'
1319	grid_y = 'y'
1320	grid_z = 'zu'
1321	!
1322	!-- v grid
1323	CASE ( 'v_xz' )
1324
1325	grid_x = 'x'
1326	grid_y = 'yv'
1327	grid_z = 'zu'
1328	!
1329	!-- w grid
1330	CASE ( 'w_xz' )
1331
1332	grid_x = 'x'
1333	grid_y = 'y'
1334	grid_z = 'zw'
1335
1336	CASE DEFAULT
1337	!
1338	!-- Check for user-defined quantities
1339	CALL user_define_netcdf_grid( do2d(av,i), found, &
1340	grid_x, grid_y, grid_z )
1341
1342	END SELECT
1343
1344	!
1345	!-- Select the respective dimension ids
1346	IF ( grid_x == 'x' ) THEN
1347	id_x = id_dim_x_xz(av)
1348	ELSEIF ( grid_x == 'xu' ) THEN
1349	id_x = id_dim_xu_xz(av)
1350	ENDIF
1351
1352	IF ( grid_y == 'y' ) THEN
1353	id_y = id_dim_y_xz(av)
1354	ELSEIF ( grid_y == 'yv' ) THEN
1355	id_y = id_dim_yv_xz(av)
1356	ENDIF
1357
1358	IF ( grid_z == 'zu' ) THEN
1359	id_z = id_dim_zu_xz(av)
1360	ELSEIF ( grid_z == 'zw' ) THEN
1361	id_z = id_dim_zw_xz(av)
1362	ENDIF
1363
1364	!
1365	!-- Define the grid
1366	nc_stat = NF90_DEF_VAR( id_set_xz(av), do2d(av,i), &
1367	nc_precision(2), &
1368	(/ id_x, id_y, id_z, id_dim_time_xz(av) /), &
1369	id_var_do2d(av,i) )
1370
1371	IF ( .NOT. found ) THEN
1372	PRINT*, '+++ define_netcdf_header: no grid defined for', &
1373	' variable ', do2d(av,i)
1374	ENDIF
1375
1376	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
1377
1378	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 159 )
1379	!
1380	!-- Store the 'real' name of the variable (with *, for example)
1381	!-- in the long_name attribute. This is evaluated by Ferret,
1382	!-- for example.
1383	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_do2d(av,i), &
1384	'long_name', do2d(av,i) )
1385	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 160 )
1386	!
1387	!-- Define the variable's unit
1388	nc_stat = NF90_PUT_ATT( id_set_xz(av), id_var_do2d(av,i), &
1389	'units', TRIM( do2d_unit(av,i) ) )
1390	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 355 )
1391	ENDIF
1392
1393	i = i + 1
1394
1395	ENDDO
1396
1397	!
1398	!-- No arrays to output. Close the netcdf file and return.
1399	IF ( i == 1 ) RETURN
1400	! nc_stat = NF90_CLOSE( id_set_xz(av) )
1401	! IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 9999 )
1402	! RETURN
1403	! ENDIF
1404
1405	!
1406	!-- Write the list of variables as global attribute (this is used by
1407	!-- restart runs and by combine_plot_fields)
1408	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'VAR_LIST', &
1409	var_list )
1410	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 161 )
1411
1412	!
1413	!-- Leave NetCDF define mode
1414	nc_stat = NF90_ENDDEF( id_set_xz(av) )
1415	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 162 )
1416
1417	!
1418	!-- Write axis data: y_xz, x, zu, zw
1419	ALLOCATE( netcdf_data(1:ns) )
1420
1421	!
1422	!-- Write y_xz data
1423	DO i = 1, ns
1424	IF( section(i,2) == -1 ) THEN
1425	netcdf_data(i) = -1.0 ! section averaged along y
1426	ELSE
1427	netcdf_data(i) = section(i,2) * dy
1428	ENDIF
1429	ENDDO
1430	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_y_xz(av), netcdf_data, &
1431	start = (/ 1 /), count = (/ ns /) )
1432	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 163 )
1433
1434	!
1435	!-- Write yv_xz data
1436	DO i = 1, ns
1437	IF( section(i,2) == -1 ) THEN
1438	netcdf_data(i) = -1.0 ! section averaged along y
1439	ELSE
1440	netcdf_data(i) = ( section(i,2) - 0.5 ) * dy
1441	ENDIF
1442	ENDDO
1443	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_yv_xz(av), &
1444	netcdf_data, start = (/ 1 /), &
1445	count = (/ ns /) )
1446	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 375 )
1447
1448	!
1449	!-- Write gridpoint number data
1450	netcdf_data(1:ns) = section(1:ns,2)
1451	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_ind_y_xz(av), &
1452	netcdf_data, start = (/ 1 /), &
1453	count = (/ ns /) )
1454	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 164 )
1455
1456	DEALLOCATE( netcdf_data )
1457
1458	!
1459	!-- Write data for x and xu axis (shifted by -dx/2)
1460	ALLOCATE( netcdf_data(0:nx+1) )
1461
1462	DO i = 0, nx+1
1463	netcdf_data(i) = i * dx
1464	ENDDO
1465
1466	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_x_xz(av), netcdf_data, &
1467	start = (/ 1 /), count = (/ nx+2 /) )
1468	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 165 )
1469
1470	DO i = 0, nx+1
1471	netcdf_data(i) = ( i - 0.5 ) * dx
1472	ENDDO
1473
1474	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_xu_xz(av), &
1475	netcdf_data, start = (/ 1 /), &
1476	count = (/ nx+2 /) )
1477	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 376 )
1478
1479	DEALLOCATE( netcdf_data )
1480
1481	!
1482	!-- Write zu and zw data (vertical axes)
1483	ALLOCATE( netcdf_data(0:nz+1) )
1484
1485	netcdf_data(0:nz+1) = zu(nzb:nzt+1)
1486	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_zu_xz(av), &
1487	netcdf_data, start = (/ 1 /), &
1488	count = (/ nz+2 /) )
1489	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 166 )
1490
1491	netcdf_data(0:nz+1) = zw(nzb:nzt+1)
1492	nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_zw_xz(av), &
1493	netcdf_data, start = (/ 1 /), &
1494	count = (/ nz+2 /) )
1495	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 167 )
1496
1497	DEALLOCATE( netcdf_data )
1498
1499
1500	CASE ( 'xz_ext' )
1501
1502	!
1503	!-- Get the list of variables and compare with the actual run.
1504	!-- First var_list_old has to be reset, since GET_ATT does not assign
1505	!-- trailing blanks.
1506	var_list_old = ' '
1507	nc_stat = NF90_GET_ATT( id_set_xz(av), NF90_GLOBAL, 'VAR_LIST', &
1508	var_list_old )
1509	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 168 )
1510
1511	var_list = ';'
1512	i = 1
1513	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
1514	IF ( INDEX( do2d(av,i), 'xz' ) /= 0 ) THEN
1515	netcdf_var_name = do2d(av,i)
1516	CALL clean_netcdf_varname( netcdf_var_name )
1517	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
1518	ENDIF
1519	i = i + 1
1520	ENDDO
1521
1522	IF ( av == 0 ) THEN
1523	var = '(xz)'
1524	ELSE
1525	var = '(xz_av)'
1526	ENDIF
1527
1528	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
1529	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1530	TRIM( var ) // ' from previuos run found,'
1531	PRINT*, ' but this file cannot be extended due to' // &
1532	' variable mismatch.'
1533	PRINT*, ' New file is created instead.'
1534	extend = .FALSE.
1535	RETURN
1536	ENDIF
1537
1538	!
1539	!-- Calculate the number of current sections
1540	ns = 1
1541	DO WHILE ( section(ns,2) /= -9999 .AND. ns <= 100 )
1542	ns = ns + 1
1543	ENDDO
1544	ns = ns - 1
1545
1546	!
1547	!-- Get and compare the number of vertical cross sections
1548	nc_stat = NF90_INQ_VARID( id_set_xz(av), 'y_xz', id_var_y_xz(av) )
1549	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 169 )
1550
1551	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xz(av), id_var_y_xz(av), &
1552	dimids = id_dim_y_xz_old )
1553	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 170 )
1554	id_dim_y_xz(av) = id_dim_y_xz_old(1)
1555
1556	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xz(av), id_dim_y_xz(av), &
1557	len = ns_old )
1558	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 171 )
1559
1560	IF ( ns /= ns_old ) THEN
1561	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1562	TRIM( var ) // ' from previuos run found,'
1563	PRINT*, ' but this file cannot be extended due to' // &
1564	' mismatch in number of'
1565	PRINT*, ' cross sections.'
1566	PRINT*, ' New file is created instead.'
1567	extend = .FALSE.
1568	RETURN
1569	ENDIF
1570
1571	!
1572	!-- Get and compare the heights of the cross sections
1573	ALLOCATE( netcdf_data(1:ns_old) )
1574
1575	nc_stat = NF90_GET_VAR( id_set_xz(av), id_var_y_xz(av), netcdf_data )
1576	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 172 )
1577
1578	DO i = 1, ns
1579	IF ( section(i,2) /= -1 ) THEN
1580	IF ( ( section(i,2) * dy ) /= netcdf_data(i) ) THEN
1581	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1582	TRIM( var ) // ' from previuos run found,'
1583	PRINT*, ' but this file cannot be extended' // &
1584	' due to mismatch in cross'
1585	PRINT*, ' section indices.'
1586	PRINT*, ' New file is created instead.'
1587	extend = .FALSE.
1588	RETURN
1589	ENDIF
1590	ELSE
1591	IF ( -1.0 /= netcdf_data(i) ) THEN
1592	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
1593	TRIM( var ) // ' from previuos run found,'
1594	PRINT*, ' but this file cannot be extended' // &
1595	' due to mismatch in cross'
1596	PRINT*, ' section indices.'
1597	PRINT*, ' New file is created instead.'
1598	extend = .FALSE.
1599	RETURN
1600	ENDIF
1601	ENDIF
1602	ENDDO
1603
1604	DEALLOCATE( netcdf_data )
1605
1606	!
1607	!-- Get the id of the time coordinate (unlimited coordinate) and its
1608	!-- last index on the file. The next time level is do2d..count+1.
1609	!-- The current time must be larger than the last output time
1610	!-- on the file.
1611	nc_stat = NF90_INQ_VARID( id_set_xz(av), 'time', id_var_time_xz(av) )
1612	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 173 )
1613
1614	nc_stat = NF90_INQUIRE_VARIABLE( id_set_xz(av), id_var_time_xz(av), &
1615	dimids = id_dim_time_old )
1616	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 174 )
1617	id_dim_time_xz(av) = id_dim_time_old(1)
1618
1619	nc_stat = NF90_INQUIRE_DIMENSION( id_set_xz(av), id_dim_time_xz(av), &
1620	len = do2d_xz_time_count(av) )
1621	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 175 )
1622
1623	nc_stat = NF90_GET_VAR( id_set_xz(av), id_var_time_xz(av), &
1624	last_time_coordinate, &
1625	start = (/ do2d_xz_time_count(av) /), &
1626	count = (/ 1 /) )
1627	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 176 )
1628
1629	IF ( last_time_coordinate(1) >= simulated_time ) THEN
1630	PRINT*, '+++ WARNING: NetCDF file for cross sections ' // &
1631	TRIM( var ) // ' from previuos run found,'
1632	PRINT*, ' but this file cannot be extended becaus' // &
1633	'e the current output time'
1634	PRINT*, ' is less or equal than the last output t' // &
1635	'ime on this file.'
1636	PRINT*, ' New file is created instead.'
1637	do2d_xz_time_count(av) = 0
1638	extend = .FALSE.
1639	RETURN
1640	ENDIF
1641
1642	!
1643	!-- Dataset seems to be extendable.
1644	!-- Now get the variable ids.
1645	i = 1
1646	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
1647	IF ( INDEX( do2d(av,i), 'xz' ) /= 0 ) THEN
1648	netcdf_var_name = do2d(av,i)
1649	CALL clean_netcdf_varname( netcdf_var_name )
1650	nc_stat = NF90_INQ_VARID( id_set_xz(av), netcdf_var_name, &
1651	id_var_do2d(av,i) )
1652	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 177 )
1653	ENDIF
1654	i = i + 1
1655	ENDDO
1656
1657	!
1658	!-- Change the titel attribute on file
1659	nc_stat = NF90_PUT_ATT( id_set_xz(av), NF90_GLOBAL, 'title', &
1660	TRIM( run_description_header ) )
1661	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 178 )
1662	PRINT, '** NetCDF file for cross-sections ' // TRIM( var ) // &
1663	' from previous run found.'
1664	PRINT*, ' This file will be extended.'
1665
1666
1667	CASE ( 'yz_new' )
1668
1669	!
1670	!-- Define some global attributes of the dataset
1671	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'Conventions', &
1672	'COARDS' )
1673	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 179 )
1674
1675	IF ( av == 0 ) THEN
1676	time_average_text = ' '
1677	ELSE
1678	WRITE (time_average_text, '('', '',F7.1,'' s average'')') &
1679	averaging_interval
1680	ENDIF
1681	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'title', &
1682	TRIM( run_description_header ) // &
1683	TRIM( time_average_text ) )
1684	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 180 )
1685	IF ( av == 1 ) THEN
1686	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval
1687	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'time_avg', &
1688	TRIM( time_average_text ) )
1689	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 180 )
1690	ENDIF
1691
1692	!
1693	!-- Define time coordinate for yz sections (unlimited dimension)
1694	nc_stat = NF90_DEF_DIM( id_set_yz(av), 'time', NF90_UNLIMITED, &
1695	id_dim_time_yz(av) )
1696	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 181 )
1697
1698	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'time', NF90_DOUBLE, &
1699	id_dim_time_yz(av), id_var_time_yz(av) )
1700	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 182 )
1701
1702	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_time_yz(av), 'units', &
1703	'seconds')
1704	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 183 )
1705
1706	!
1707	!-- Define the spatial dimensions and coordinates for yz-sections.
1708	!-- First, determine the number of vertical sections to be written.
1709	IF ( section(1,3) == -9999 ) THEN
1710	RETURN
1711	ELSE
1712	ns = 1
1713	DO WHILE ( section(ns,3) /= -9999 .AND. ns <= 100 )
1714	ns = ns + 1
1715	ENDDO
1716	ns = ns - 1
1717	ENDIF
1718
1719	!
1720	!-- Define x axis (for scalar position)
1721	nc_stat = NF90_DEF_DIM( id_set_yz(av), 'x_yz', ns, id_dim_x_yz(av) )
1722	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 184 )
1723
1724	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'x_yz', NF90_DOUBLE, &
1725	id_dim_x_yz(av), id_var_x_yz(av) )
1726	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 185 )
1727
1728	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_x_yz(av), 'units', &
1729	'meters' )
1730	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 186 )
1731
1732	!
1733	!-- Define x axis (for u position)
1734	nc_stat = NF90_DEF_DIM( id_set_yz(av), 'xu_yz', ns, id_dim_xu_yz(av) )
1735	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 377 )
1736
1737	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'xu_yz', NF90_DOUBLE, &
1738	id_dim_xu_yz(av), id_var_xu_yz(av) )
1739	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 378 )
1740
1741	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_xu_yz(av), 'units', &
1742	'meters' )
1743	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 379 )
1744
1745	!
1746	!-- Define a variable to store the layer indices of the vertical cross
1747	!-- sections
1748	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'ind_x_yz', NF90_DOUBLE, &
1749	id_dim_x_yz(av), id_var_ind_x_yz(av) )
1750	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 187 )
1751
1752	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_ind_x_yz(av), 'units', &
1753	'gridpoints')
1754	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 188 )
1755
1756	!
1757	!-- Define y-axis (for scalar position)
1758	nc_stat = NF90_DEF_DIM( id_set_yz(av), 'y', ny+2, id_dim_y_yz(av) )
1759	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 189 )
1760
1761	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'y', NF90_DOUBLE, &
1762	id_dim_y_yz(av), id_var_y_yz(av) )
1763	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 190 )
1764
1765	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_y_yz(av), 'units', &
1766	'meters' )
1767	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 191 )
1768
1769	!
1770	!-- Define y-axis (for v position)
1771	nc_stat = NF90_DEF_DIM( id_set_yz(av), 'yv', ny+2, id_dim_yv_yz(av) )
1772	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 380 )
1773
1774	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'yv', NF90_DOUBLE, &
1775	id_dim_yv_yz(av), id_var_yv_yz(av) )
1776	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 381 )
1777
1778	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_yv_yz(av), 'units', &
1779	'meters' )
1780	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 382 )
1781
1782	!
1783	!-- Define the two z-axes (zu and zw)
1784	nc_stat = NF90_DEF_DIM( id_set_yz(av), 'zu', nz+2, id_dim_zu_yz(av) )
1785	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 192 )
1786
1787	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'zu', NF90_DOUBLE, &
1788	id_dim_zu_yz(av), id_var_zu_yz(av) )
1789	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 193 )
1790
1791	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_zu_yz(av), 'units', &
1792	'meters' )
1793	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 194 )
1794
1795	nc_stat = NF90_DEF_DIM( id_set_yz(av), 'zw', nz+2, id_dim_zw_yz(av) )
1796	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 195 )
1797
1798	nc_stat = NF90_DEF_VAR( id_set_yz(av), 'zw', NF90_DOUBLE, &
1799	id_dim_zw_yz(av), id_var_zw_yz(av) )
1800	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 196 )
1801
1802	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_zw_yz(av), 'units', &
1803	'meters' )
1804	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 197 )
1805
1806
1807	!
1808	!-- Define the variables
1809	var_list = ';'
1810	i = 1
1811
1812	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
1813
1814	IF ( INDEX( do2d(av,i), 'yz' ) /= 0 ) THEN
1815
1816	!
1817	!-- Check for the grid
1818	found = .TRUE.
1819	SELECT CASE ( do2d(av,i) )
1820	!
1821	!-- Most variables are defined on the zu grid
1822	CASE ( 'e_yz', 'p_yz', 'pc_yz', 'pr_yz', 'pt_yz', 'q_yz', &
1823	'ql_yz', 'ql_c_yz', 'ql_v_yz', 'ql_vp_yz', 'qv_yz', &
1824	's_yz', 'vpt_yz' )
1825
1826	grid_x = 'x'
1827	grid_y = 'y'
1828	grid_z = 'zu'
1829	!
1830	!-- u grid
1831	CASE ( 'u_yz' )
1832
1833	grid_x = 'xu'
1834	grid_y = 'y'
1835	grid_z = 'zu'
1836	!
1837	!-- v grid
1838	CASE ( 'v_yz' )
1839
1840	grid_x = 'x'
1841	grid_y = 'yv'
1842	grid_z = 'zu'
1843	!
1844	!-- w grid
1845	CASE ( 'w_yz' )
1846
1847	grid_x = 'x'
1848	grid_y = 'y'
1849	grid_z = 'zw'
1850
1851	CASE DEFAULT
1852	!
1853	!-- Check for user-defined quantities
1854	CALL user_define_netcdf_grid( do2d(av,i), found, &
1855	grid_x, grid_y, grid_z )
1856
1857	END SELECT
1858
1859	!
1860	!-- Select the respective dimension ids
1861	IF ( grid_x == 'x' ) THEN
1862	id_x = id_dim_x_yz(av)
1863	ELSEIF ( grid_x == 'xu' ) THEN
1864	id_x = id_dim_xu_yz(av)
1865	ENDIF
1866
1867	IF ( grid_y == 'y' ) THEN
1868	id_y = id_dim_y_yz(av)
1869	ELSEIF ( grid_y == 'yv' ) THEN
1870	id_y = id_dim_yv_yz(av)
1871	ENDIF
1872
1873	IF ( grid_z == 'zu' ) THEN
1874	id_z = id_dim_zu_yz(av)
1875	ELSEIF ( grid_z == 'zw' ) THEN
1876	id_z = id_dim_zw_yz(av)
1877	ENDIF
1878
1879	!
1880	!-- Define the grid
1881	nc_stat = NF90_DEF_VAR( id_set_yz(av), do2d(av,i), &
1882	nc_precision(3), &
1883	(/ id_x, id_y, id_z, id_dim_time_yz(av) /), &
1884	id_var_do2d(av,i) )
1885
1886	IF ( .NOT. found ) THEN
1887	PRINT*, '+++ define_netcdf_header: no grid defined for', &
1888	' variable ', do2d(av,i)
1889	ENDIF
1890
1891	var_list = TRIM( var_list ) // TRIM( do2d(av,i) ) // ';'
1892
1893	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 198 )
1894	!
1895	!-- Store the 'real' name of the variable (with *, for example)
1896	!-- in the long_name attribute. This is evaluated by Ferret,
1897	!-- for example.
1898	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_do2d(av,i), &
1899	'long_name', do2d(av,i) )
1900	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 199 )
1901	!
1902	!-- Define the variable's unit
1903	nc_stat = NF90_PUT_ATT( id_set_yz(av), id_var_do2d(av,i), &
1904	'units', TRIM( do2d_unit(av,i) ) )
1905	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 356 )
1906	ENDIF
1907
1908	i = i + 1
1909
1910	ENDDO
1911
1912	!
1913	!-- No arrays to output. Close the netcdf file and return.
1914	IF ( i == 1 ) RETURN
1915	! nc_stat = NF90_CLOSE( id_set_yz(av) )
1916	! IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 9999 )
1917	! RETURN
1918	! ENDIF
1919
1920	!
1921	!-- Write the list of variables as global attribute (this is used by
1922	!-- restart runs and by combine_plot_fields)
1923	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'VAR_LIST', &
1924	var_list )
1925	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 200 )
1926
1927	!
1928	!-- Leave NetCDF define mode
1929	nc_stat = NF90_ENDDEF( id_set_yz(av) )
1930	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 201 )
1931
1932	!
1933	!-- Write axis data: x_yz, y, zu, zw
1934	ALLOCATE( netcdf_data(1:ns) )
1935
1936	!
1937	!-- Write x_yz data
1938	DO i = 1, ns
1939	IF( section(i,3) == -1 ) THEN
1940	netcdf_data(i) = -1.0 ! section averaged along x
1941	ELSE
1942	netcdf_data(i) = section(i,3) * dx
1943	ENDIF
1944	ENDDO
1945	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_x_yz(av), netcdf_data, &
1946	start = (/ 1 /), count = (/ ns /) )
1947	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 202 )
1948
1949	!
1950	!-- Write x_yz data (xu grid)
1951	DO i = 1, ns
1952	IF( section(i,3) == -1 ) THEN
1953	netcdf_data(i) = -1.0 ! section averaged along x
1954	ELSE
1955	netcdf_data(i) = (section(i,3)-0.5) * dx
1956	ENDIF
1957	ENDDO
1958	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_xu_yz(av), netcdf_data, &
1959	start = (/ 1 /), count = (/ ns /) )
1960	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 383 )
1961
1962	!
1963	!-- Write gridpoint number data
1964	netcdf_data(1:ns) = section(1:ns,3)
1965	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_ind_x_yz(av), &
1966	netcdf_data, start = (/ 1 /), &
1967	count = (/ ns /) )
1968	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 203 )
1969
1970	DEALLOCATE( netcdf_data )
1971
1972	!
1973	!-- Write data for y and yv axis (shifted by -dy/2)
1974	ALLOCATE( netcdf_data(0:ny+1) )
1975
1976	DO j = 0, ny+1
1977	netcdf_data(j) = j * dy
1978	ENDDO
1979
1980	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_y_yz(av), netcdf_data, &
1981	start = (/ 1 /), count = (/ ny+2 /) )
1982	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 204 )
1983
1984	DO j = 0, ny+1
1985	netcdf_data(j) = ( j - 0.5 ) * dy
1986	ENDDO
1987
1988	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_yv_yz(av), &
1989	netcdf_data, start = (/ 1 /), &
1990	count = (/ ny+2 /) )
1991	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 384 )
1992
1993	DEALLOCATE( netcdf_data )
1994
1995	!
1996	!-- Write zu and zw data (vertical axes)
1997	ALLOCATE( netcdf_data(0:nz+1) )
1998
1999	netcdf_data(0:nz+1) = zu(nzb:nzt+1)
2000	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_zu_yz(av), &
2001	netcdf_data, start = (/ 1 /), &
2002	count = (/ nz+2 /) )
2003	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 205 )
2004
2005	netcdf_data(0:nz+1) = zw(nzb:nzt+1)
2006	nc_stat = NF90_PUT_VAR( id_set_yz(av), id_var_zw_yz(av), &
2007	netcdf_data, start = (/ 1 /), &
2008	count = (/ nz+2 /) )
2009	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 206 )
2010
2011	DEALLOCATE( netcdf_data )
2012
2013
2014	CASE ( 'yz_ext' )
2015
2016	!
2017	!-- Get the list of variables and compare with the actual run.
2018	!-- First var_list_old has to be reset, since GET_ATT does not assign
2019	!-- trailing blanks.
2020	var_list_old = ' '
2021	nc_stat = NF90_GET_ATT( id_set_yz(av), NF90_GLOBAL, 'VAR_LIST', &
2022	var_list_old )
2023	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 207 )
2024
2025	var_list = ';'
2026	i = 1
2027	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
2028	IF ( INDEX( do2d(av,i), 'yz' ) /= 0 ) THEN
2029	netcdf_var_name = do2d(av,i)
2030	CALL clean_netcdf_varname( netcdf_var_name )
2031	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
2032	ENDIF
2033	i = i + 1
2034	ENDDO
2035
2036	IF ( av == 0 ) THEN
2037	var = '(yz)'
2038	ELSE
2039	var = '(yz_av)'
2040	ENDIF
2041
2042	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
2043	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
2044	TRIM( var ) // ' from previuos run found,'
2045	PRINT*, ' but this file cannot be extended due to' // &
2046	' variable mismatch.'
2047	PRINT*, ' New file is created instead.'
2048	extend = .FALSE.
2049	RETURN
2050	ENDIF
2051
2052	!
2053	!-- Calculate the number of current sections
2054	ns = 1
2055	DO WHILE ( section(ns,3) /= -9999 .AND. ns <= 100 )
2056	ns = ns + 1
2057	ENDDO
2058	ns = ns - 1
2059
2060	!
2061	!-- Get and compare the number of vertical cross sections
2062	nc_stat = NF90_INQ_VARID( id_set_yz(av), 'x_yz', id_var_x_yz(av) )
2063	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 208 )
2064
2065	nc_stat = NF90_INQUIRE_VARIABLE( id_set_yz(av), id_var_x_yz(av), &
2066	dimids = id_dim_x_yz_old )
2067	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 209 )
2068	id_dim_x_yz(av) = id_dim_x_yz_old(1)
2069
2070	nc_stat = NF90_INQUIRE_DIMENSION( id_set_yz(av), id_dim_x_yz(av), &
2071	len = ns_old )
2072	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 210 )
2073
2074	IF ( ns /= ns_old ) THEN
2075	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
2076	TRIM( var ) // ' from previuos run found,'
2077	PRINT*, ' but this file cannot be extended due to' // &
2078	' mismatch in number of'
2079	PRINT*, ' cross sections.'
2080	PRINT*, ' New file is created instead.'
2081	extend = .FALSE.
2082	RETURN
2083	ENDIF
2084
2085	!
2086	!-- Get and compare the heights of the cross sections
2087	ALLOCATE( netcdf_data(1:ns_old) )
2088
2089	nc_stat = NF90_GET_VAR( id_set_yz(av), id_var_x_yz(av), netcdf_data )
2090	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 211 )
2091
2092	DO i = 1, ns
2093	IF ( section(i,3) /= -1 ) THEN
2094	IF ( ( section(i,3) * dx ) /= netcdf_data(i) ) THEN
2095	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
2096	TRIM( var ) // ' from previuos run found,'
2097	PRINT*, ' but this file cannot be extended' // &
2098	' due to mismatch in cross'
2099	PRINT*, ' section indices.'
2100	PRINT*, ' New file is created instead.'
2101	extend = .FALSE.
2102	RETURN
2103	ENDIF
2104	ELSE
2105	IF ( -1.0 /= netcdf_data(i) ) THEN
2106	PRINT*, '+++ WARNING: NetCDF file for cross-sections ' // &
2107	TRIM( var ) // ' from previuos run found,'
2108	PRINT*, ' but this file cannot be extended' // &
2109	' due to mismatch in cross'
2110	PRINT*, ' section indices.'
2111	PRINT*, ' New file is created instead.'
2112	extend = .FALSE.
2113	RETURN
2114	ENDIF
2115	ENDIF
2116	ENDDO
2117
2118	DEALLOCATE( netcdf_data )
2119
2120	!
2121	!-- Get the id of the time coordinate (unlimited coordinate) and its
2122	!-- last index on the file. The next time level is pl2d..count+1.
2123	!-- The current time must be larger than the last output time
2124	!-- on the file.
2125	nc_stat = NF90_INQ_VARID( id_set_yz(av), 'time', id_var_time_yz(av) )
2126	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 212 )
2127
2128	nc_stat = NF90_INQUIRE_VARIABLE( id_set_yz(av), id_var_time_yz(av), &
2129	dimids = id_dim_time_old )
2130	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 213 )
2131	id_dim_time_yz(av) = id_dim_time_old(1)
2132
2133	nc_stat = NF90_INQUIRE_DIMENSION( id_set_yz(av), id_dim_time_yz(av), &
2134	len = do2d_yz_time_count(av) )
2135	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 214 )
2136
2137	nc_stat = NF90_GET_VAR( id_set_yz(av), id_var_time_yz(av), &
2138	last_time_coordinate, &
2139	start = (/ do2d_yz_time_count(av) /), &
2140	count = (/ 1 /) )
2141	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 215 )
2142
2143	IF ( last_time_coordinate(1) >= simulated_time ) THEN
2144	PRINT*, '+++ WARNING: NetCDF file for cross sections ' // &
2145	TRIM( var ) // ' from previuos run found,'
2146	PRINT*, ' but this file cannot be extended becaus' // &
2147	'e the current output time'
2148	PRINT*, ' is less or equal than the last output t' // &
2149	'ime on this file.'
2150	PRINT*, ' New file is created instead.'
2151	do2d_yz_time_count(av) = 0
2152	extend = .FALSE.
2153	RETURN
2154	ENDIF
2155
2156	!
2157	!-- Dataset seems to be extendable.
2158	!-- Now get the variable ids.
2159	i = 1
2160	DO WHILE ( do2d(av,i)(1:1) /= ' ' )
2161	IF ( INDEX( do2d(av,i), 'yz' ) /= 0 ) THEN
2162	netcdf_var_name = do2d(av,i)
2163	CALL clean_netcdf_varname( netcdf_var_name )
2164	nc_stat = NF90_INQ_VARID( id_set_yz(av), netcdf_var_name, &
2165	id_var_do2d(av,i) )
2166	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 216 )
2167	ENDIF
2168	i = i + 1
2169	ENDDO
2170
2171	!
2172	!-- Change the titel attribute on file
2173	nc_stat = NF90_PUT_ATT( id_set_yz(av), NF90_GLOBAL, 'title', &
2174	TRIM( run_description_header ) )
2175	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 217 )
2176	PRINT, '** NetCDF file for cross-sections ' // TRIM( var ) // &
2177	' from previous run found.'
2178	PRINT*, ' This file will be extended.'
2179
2180
2181	CASE ( 'pr_new' )
2182
2183	!
2184	!-- Define some global attributes of the dataset
2185	IF ( averaging_interval_pr /= 0.0 ) THEN
2186	WRITE (time_average_text,'('', '',F7.1,'' s average'')') &
2187	averaging_interval_pr
2188	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'title', &
2189	TRIM( run_description_header ) // &
2190	TRIM( time_average_text ) )
2191	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 218 )
2192
2193	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval_pr
2194	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'time_avg', &
2195	TRIM( time_average_text ) )
2196	ELSE
2197	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'title', &
2198	TRIM( run_description_header ) )
2199	ENDIF
2200	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 219 )
2201
2202	!
2203	!-- Define time coordinate for profiles (unlimited dimension)
2204	nc_stat = NF90_DEF_DIM( id_set_pr, 'time', NF90_UNLIMITED, &
2205	id_dim_time_pr )
2206	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 220 )
2207
2208	nc_stat = NF90_DEF_VAR( id_set_pr, 'time', NF90_DOUBLE, &
2209	id_dim_time_pr, id_var_time_pr )
2210	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 221 )
2211
2212	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_time_pr, 'units', 'seconds')
2213	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 222 )
2214
2215	!
2216	!-- Define the variables
2217	var_list = ';'
2218	DO i = 1, dopr_n
2219	!
2220	!-- First, remove those characters not allowed by NetCDF
2221	netcdf_var_name = data_output_pr(i)
2222	CALL clean_netcdf_varname( netcdf_var_name )
2223
2224	IF ( statistic_regions == 0 ) THEN
2225
2226	!
2227	!-- Define the z-axes (each variable gets its own z-axis)
2228	nc_stat = NF90_DEF_DIM( id_set_pr, 'z'//TRIM(netcdf_var_name), &
2229	nzt+2-nzb, id_dim_z_pr(i,0) )
2230	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 223 )
2231
2232	nc_stat = NF90_DEF_VAR( id_set_pr, 'z'//TRIM(netcdf_var_name), &
2233	NF90_DOUBLE, id_dim_z_pr(i,0), &
2234	id_var_z_pr(i,0) )
2235	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 224 )
2236
2237	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_z_pr(i,0), 'units', &
2238	'meters' )
2239	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 225 )
2240
2241	!
2242	!-- Define the variable
2243	nc_stat = NF90_DEF_VAR( id_set_pr, netcdf_var_name, &
2244	nc_precision(5), (/ id_dim_z_pr(i,0), &
2245	id_dim_time_pr /), id_var_dopr(i,0) )
2246	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 226 )
2247
2248	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_dopr(i,0), &
2249	'long_name', TRIM( data_output_pr(i) ) )
2250	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 227 )
2251
2252	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_dopr(i,0), &
2253	'units', TRIM( dopr_unit(i) ) )
2254	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 228 )
2255
2256	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
2257
2258	ELSE
2259	!
2260	!-- If statistic regions are defined, add suffix _SR+#SR to the
2261	!-- names
2262	DO j = 0, statistic_regions
2263	WRITE ( suffix, '(''_'',I1)' ) j
2264
2265	!
2266	!-- Define the z-axes (each variable gets it own z-axis)
2267	nc_stat = NF90_DEF_DIM( id_set_pr, &
2268	'z'//TRIM(netcdf_var_name)//suffix, &
2269	nzt+2-nzb, id_dim_z_pr(i,j) )
2270	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 229 )
2271
2272	nc_stat = NF90_DEF_VAR( id_set_pr, &
2273	'z'//TRIM(netcdf_var_name)//suffix, &
2274	nc_precision(5), id_dim_z_pr(i,j), &
2275	id_var_z_pr(i,j) )
2276	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 230 )
2277
2278	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_z_pr(i,j), &
2279	'units', 'meters' )
2280	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 231 )
2281
2282	!
2283	!-- Define the variable
2284	nc_stat = NF90_DEF_VAR( id_set_pr, &
2285	TRIM( netcdf_var_name ) // suffix, &
2286	nc_precision(5), &
2287	(/ id_dim_z_pr(i,j), &
2288	id_dim_time_pr /), id_var_dopr(i,j) )
2289	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 232 )
2290
2291	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_dopr(i,j), &
2292	'long_name', &
2293	TRIM( data_output_pr(i) ) // ' SR ' &
2294	// suffix(2:2) )
2295	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 233 )
2296
2297	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_dopr(i,j), &
2298	'units', TRIM( dopr_unit(i) ) )
2299	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 234 )
2300
2301	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // &
2302	suffix // ';'
2303
2304	ENDDO
2305
2306	ENDIF
2307
2308	ENDDO
2309
2310	!
2311	!-- Write the list of variables as global attribute (this is used by
2312	!-- restart runs)
2313	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'VAR_LIST', var_list )
2314	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 235 )
2315
2316	!
2317	!-- Define normalization variables (as time series)
2318	DO i = 1, dopr_norm_num
2319
2320	nc_stat = NF90_DEF_VAR( id_set_pr, 'NORM_' // &
2321	TRIM( dopr_norm_names(i) ), &
2322	nc_precision(5), (/ id_dim_time_pr /), &
2323	id_var_norm_dopr(i) )
2324	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 236 )
2325
2326	nc_stat = NF90_PUT_ATT( id_set_pr, id_var_norm_dopr(i), &
2327	'long_name', &
2328	TRIM( dopr_norm_longnames(i) ) )
2329	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 237 )
2330
2331	ENDDO
2332
2333	!
2334	!-- Leave NetCDF define mode
2335	nc_stat = NF90_ENDDEF( id_set_pr )
2336	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 238 )
2337
2338	!
2339	!-- Write z-axes data
2340	DO i = 1, dopr_n
2341	DO j = 0, statistic_regions
2342
2343	nc_stat = NF90_PUT_VAR( id_set_pr, id_var_z_pr(i,j), &
2344	hom(nzb:nzt+1,2,dopr_index(i),0), &
2345	start = (/ 1 /), &
2346	count = (/ nzt-nzb+2 /) )
2347	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 239 )
2348
2349	ENDDO
2350	ENDDO
2351
2352
2353	CASE ( 'pr_ext' )
2354
2355	!
2356	!-- Get the list of variables and compare with the actual run.
2357	!-- First var_list_old has to be reset, since GET_ATT does not assign
2358	!-- trailing blanks.
2359	var_list_old = ' '
2360	nc_stat = NF90_GET_ATT( id_set_pr, NF90_GLOBAL, 'VAR_LIST', &
2361	var_list_old )
2362	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 240 )
2363
2364	var_list = ';'
2365	DO i = 1, dopr_n
2366
2367	netcdf_var_name = data_output_pr(i)
2368	CALL clean_netcdf_varname( netcdf_var_name )
2369
2370	IF ( statistic_regions == 0 ) THEN
2371	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
2372	ELSE
2373	DO j = 0, statistic_regions
2374	WRITE ( suffix, '(''_'',I1)' ) j
2375	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // &
2376	suffix // ';'
2377	ENDDO
2378	ENDIF
2379
2380	ENDDO
2381
2382	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
2383	PRINT*, '+++ WARNING: NetCDF file for vertical profiles from' // &
2384	' previuos run found,'
2385	PRINT*, ' but this file cannot be extended due to' // &
2386	' variable mismatch.'
2387	PRINT*, ' New file is created instead.'
2388	extend = .FALSE.
2389	RETURN
2390	ENDIF
2391
2392	!
2393	!-- Get the id of the time coordinate (unlimited coordinate) and its
2394	!-- last index on the file. The next time level is dopr..count+1.
2395	!-- The current time must be larger than the last output time
2396	!-- on the file.
2397	nc_stat = NF90_INQ_VARID( id_set_pr, 'time', id_var_time_pr )
2398	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 241 )
2399
2400	nc_stat = NF90_INQUIRE_VARIABLE( id_set_pr, id_var_time_pr, &
2401	dimids = id_dim_time_old )
2402	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 242 )
2403	id_dim_time_pr = id_dim_time_old(1)
2404
2405	nc_stat = NF90_INQUIRE_DIMENSION( id_set_pr, id_dim_time_pr, &
2406	len = dopr_time_count )
2407	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 243 )
2408
2409	nc_stat = NF90_GET_VAR( id_set_pr, id_var_time_pr, &
2410	last_time_coordinate, &
2411	start = (/ dopr_time_count /), &
2412	count = (/ 1 /) )
2413	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 244 )
2414
2415	IF ( last_time_coordinate(1) >= simulated_time ) THEN
2416	PRINT*, '+++ WARNING: NetCDF file for vertical profiles from' // &
2417	' previuos run found,'
2418	PRINT*, ' but this file cannot be extended becaus' // &
2419	'e the current output time'
2420	PRINT*, ' is less or equal than the last output t' // &
2421	'ime on this file.'
2422	PRINT*, ' New file is created instead.'
2423	dopr_time_count = 0
2424	extend = .FALSE.
2425	RETURN
2426	ENDIF
2427
2428	!
2429	!-- Dataset seems to be extendable.
2430	!-- Now get the variable ids.
2431	i = 1
2432	DO i = 1, dopr_n
2433
2434	netcdf_var_name_base = data_output_pr(i)
2435	CALL clean_netcdf_varname( netcdf_var_name_base )
2436
2437	IF ( statistic_regions == 0 ) THEN
2438	nc_stat = NF90_INQ_VARID( id_set_pr, netcdf_var_name_base, &
2439	id_var_dopr(i,0) )
2440	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 245 )
2441	ELSE
2442	DO j = 0, statistic_regions
2443	WRITE ( suffix, '(''_'',I1)' ) j
2444	netcdf_var_name = TRIM( netcdf_var_name_base ) // suffix
2445	nc_stat = NF90_INQ_VARID( id_set_pr, netcdf_var_name, &
2446	id_var_dopr(i,j) )
2447	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 246 )
2448	ENDDO
2449	ENDIF
2450
2451	ENDDO
2452
2453	!
2454	!-- Get ids of the normalization variables
2455	DO i = 1, dopr_norm_num
2456	nc_stat = NF90_INQ_VARID( id_set_pr, &
2457	'NORM_' // TRIM( dopr_norm_names(i) ), &
2458	id_var_norm_dopr(i) )
2459	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 247 )
2460	ENDDO
2461
2462	!
2463	!-- Change the title attribute on file
2464	nc_stat = NF90_PUT_ATT( id_set_pr, NF90_GLOBAL, 'title', &
2465	TRIM( run_description_header ) )
2466	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 248 )
2467	PRINT, '** NetCDF file for vertical profiles from previous run' // &
2468	' found.'
2469	PRINT*, ' This file will be extended.'
2470
2471
2472	CASE ( 'ts_new' )
2473
2474	!
2475	!-- Define some global attributes of the dataset
2476	nc_stat = NF90_PUT_ATT( id_set_ts, NF90_GLOBAL, 'title', &
2477	TRIM( run_description_header ) )
2478	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 249 )
2479
2480	!
2481	!-- Define time coordinate for time series (unlimited dimension)
2482	nc_stat = NF90_DEF_DIM( id_set_ts, 'time', NF90_UNLIMITED, &
2483	id_dim_time_ts )
2484	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 250 )
2485
2486	nc_stat = NF90_DEF_VAR( id_set_ts, 'time', NF90_DOUBLE, &
2487	id_dim_time_ts, id_var_time_ts )
2488	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 251 )
2489
2490	nc_stat = NF90_PUT_ATT( id_set_ts, id_var_time_ts, 'units', 'seconds')
2491	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 252 )
2492
2493	!
2494	!-- Define the variables
2495	var_list = ';'
2496	DO i = 1, dots_num
2497	!
2498	!-- First, remove those characters not allowed by NetCDF
2499	netcdf_var_name = dots_label(i)
2500	CALL clean_netcdf_varname( netcdf_var_name )
2501
2502	IF ( statistic_regions == 0 ) THEN
2503
2504	nc_stat = NF90_DEF_VAR( id_set_ts, netcdf_var_name, &
2505	nc_precision(6), (/ id_dim_time_ts /), &
2506	id_var_dots(i,0) )
2507	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 253 )
2508
2509	nc_stat = NF90_PUT_ATT( id_set_ts, id_var_dots(i,0), &
2510	'long_name', TRIM( dots_label(i) ) )
2511	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 254 )
2512
2513	nc_stat = NF90_PUT_ATT( id_set_ts, id_var_dots(i,0), &
2514	'units', TRIM( dots_unit(i) ) )
2515	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 255 )
2516
2517	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
2518
2519	ELSE
2520	!
2521	!-- If statistic regions are defined, add suffix _SR+#SR to the
2522	!-- names
2523	DO j = 0, statistic_regions
2524	WRITE ( suffix, '(''_'',I1)' ) j
2525
2526	nc_stat = NF90_DEF_VAR( id_set_ts, &
2527	TRIM( netcdf_var_name ) // suffix, &
2528	nc_precision(6), &
2529	(/ id_dim_time_ts /), &
2530	id_var_dots(i,j) )
2531	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 256 )
2532
2533	nc_stat = NF90_PUT_ATT( id_set_ts, id_var_dots(i,j), &
2534	'long_name', &
2535	TRIM( dots_label(i) ) // ' SR ' // &
2536	suffix(2:2) )
2537	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 257 )
2538
2539	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // &
2540	suffix // ';'
2541
2542	ENDDO
2543
2544	ENDIF
2545
2546	ENDDO
2547
2548	!
2549	!-- Write the list of variables as global attribute (this is used by
2550	!-- restart runs)
2551	nc_stat = NF90_PUT_ATT( id_set_ts, NF90_GLOBAL, 'VAR_LIST', var_list )
2552	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 258 )
2553
2554	!
2555	!-- Leave NetCDF define mode
2556	nc_stat = NF90_ENDDEF( id_set_ts )
2557	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 259 )
2558
2559
2560	CASE ( 'ts_ext' )
2561
2562	!
2563	!-- Get the list of variables and compare with the actual run.
2564	!-- First var_list_old has to be reset, since GET_ATT does not assign
2565	!-- trailing blanks.
2566	var_list_old = ' '
2567	nc_stat = NF90_GET_ATT( id_set_ts, NF90_GLOBAL, 'VAR_LIST', &
2568	var_list_old )
2569	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 260 )
2570
2571	var_list = ';'
2572	i = 1
2573	DO i = 1, dots_num
2574
2575	netcdf_var_name = dots_label(i)
2576	CALL clean_netcdf_varname( netcdf_var_name )
2577
2578	IF ( statistic_regions == 0 ) THEN
2579	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // ';'
2580	ELSE
2581	DO j = 0, statistic_regions
2582	WRITE ( suffix, '(''_'',I1)' ) j
2583	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // &
2584	suffix // ';'
2585	ENDDO
2586	ENDIF
2587
2588	ENDDO
2589
2590	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
2591	PRINT*, '+++ WARNING: NetCDF file for time series from' //&
2592	' previuos run found,'
2593	PRINT*, ' but this file cannot be extended due to' // &
2594	' variable mismatch.'
2595	PRINT*, ' New file is created instead.'
2596	extend = .FALSE.
2597	RETURN
2598	ENDIF
2599
2600	!
2601	!-- Get the id of the time coordinate (unlimited coordinate) and its
2602	!-- last index on the file. The next time level is dots..count+1.
2603	!-- The current time must be larger than the last output time
2604	!-- on the file.
2605	nc_stat = NF90_INQ_VARID( id_set_ts, 'time', id_var_time_ts )
2606	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 261 )
2607
2608	nc_stat = NF90_INQUIRE_VARIABLE( id_set_ts, id_var_time_ts, &
2609	dimids = id_dim_time_old )
2610	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 262 )
2611	id_dim_time_ts = id_dim_time_old(1)
2612
2613	nc_stat = NF90_INQUIRE_DIMENSION( id_set_ts, id_dim_time_ts, &
2614	len = dots_time_count )
2615	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 263 )
2616
2617	nc_stat = NF90_GET_VAR( id_set_ts, id_var_time_ts, &
2618	last_time_coordinate, &
2619	start = (/ dots_time_count /), &
2620	count = (/ 1 /) )
2621	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 264 )
2622
2623	IF ( last_time_coordinate(1) >= simulated_time ) THEN
2624	PRINT*, '+++ WARNING: NetCDF file for time series from' // &
2625	' previuos run found,'
2626	PRINT*, ' but this file cannot be extended becaus' // &
2627	'e the current output time'
2628	PRINT*, ' is less or equal than the last output t' // &
2629	'ime on this file.'
2630	PRINT*, ' New file is created instead.'
2631	dots_time_count = 0
2632	extend = .FALSE.
2633	RETURN
2634	ENDIF
2635
2636	!
2637	!-- Dataset seems to be extendable.
2638	!-- Now get the variable ids
2639	i = 1
2640	DO i = 1, dots_num
2641
2642	netcdf_var_name_base = dots_label(i)
2643	CALL clean_netcdf_varname( netcdf_var_name_base )
2644
2645	IF ( statistic_regions == 0 ) THEN
2646	nc_stat = NF90_INQ_VARID( id_set_ts, netcdf_var_name_base, &
2647	id_var_dots(i,0) )
2648	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 265 )
2649	ELSE
2650	DO j = 0, statistic_regions
2651	WRITE ( suffix, '(''_'',I1)' ) j
2652	netcdf_var_name = TRIM( netcdf_var_name_base ) // suffix
2653	nc_stat = NF90_INQ_VARID( id_set_ts, netcdf_var_name, &
2654	id_var_dots(i,j) )
2655	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 266 )
2656	ENDDO
2657	ENDIF
2658
2659	ENDDO
2660
2661	!
2662	!-- Change the title attribute on file
2663	nc_stat = NF90_PUT_ATT( id_set_ts, NF90_GLOBAL, 'title', &
2664	TRIM( run_description_header ) )
2665	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 267 )
2666	PRINT, '** NetCDF file for time series from previous run found.'
2667	PRINT*, ' This file will be extended.'
2668
2669
2670	CASE ( 'sp_new' )
2671
2672	!
2673	!-- Define some global attributes of the dataset
2674	IF ( averaging_interval_sp /= 0.0 ) THEN
2675	WRITE (time_average_text,'('', '',F7.1,'' s average'')') &
2676	averaging_interval_sp
2677	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
2678	TRIM( run_description_header ) // &
2679	TRIM( time_average_text ) )
2680	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 268 )
2681
2682	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval_sp
2683	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'time_avg', &
2684	TRIM( time_average_text ) )
2685	ELSE
2686	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
2687	TRIM( run_description_header ) )
2688	ENDIF
2689	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 269 )
2690
2691	!
2692	!-- Define time coordinate for spectra (unlimited dimension)
2693	nc_stat = NF90_DEF_DIM( id_set_sp, 'time', NF90_UNLIMITED, &
2694	id_dim_time_sp )
2695	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 270 )
2696
2697	nc_stat = NF90_DEF_VAR( id_set_sp, 'time', NF90_DOUBLE, &
2698	id_dim_time_sp, id_var_time_sp )
2699	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 271 )
2700
2701	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_time_sp, 'units', 'seconds')
2702	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 272 )
2703
2704	!
2705	!-- Define the spatial dimensions and coordinates for spectra.
2706	!-- First, determine the number of vertical levels for which spectra
2707	!-- are to be output.
2708	ns = 1
2709	DO WHILE ( comp_spectra_level(ns) /= 999999 .AND. ns <= 100 )
2710	ns = ns + 1
2711	ENDDO
2712	ns = ns - 1
2713
2714	!
2715	!-- Define vertical coordinate grid (zu grid)
2716	nc_stat = NF90_DEF_DIM( id_set_sp, 'zu_sp', ns, id_dim_zu_sp )
2717	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 273 )
2718
2719	nc_stat = NF90_DEF_VAR( id_set_sp, 'zu_sp', NF90_DOUBLE, &
2720	id_dim_zu_sp, id_var_zu_sp )
2721	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 274 )
2722
2723	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_zu_sp, 'units', 'meters' )
2724	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 275 )
2725
2726	!
2727	!-- Define vertical coordinate grid (zw grid)
2728	nc_stat = NF90_DEF_DIM( id_set_sp, 'zw_sp', ns, id_dim_zw_sp )
2729	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 276 )
2730
2731	nc_stat = NF90_DEF_VAR( id_set_sp, 'zw_sp', NF90_DOUBLE, &
2732	id_dim_zw_sp, id_var_zw_sp )
2733	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 277 )
2734
2735	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_zw_sp, 'units', 'meters' )
2736	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 278 )
2737
2738	!
2739	!-- Define x-axis
2740	nc_stat = NF90_DEF_DIM( id_set_sp, 'k_x', nx/2, id_dim_x_sp )
2741	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 279 )
2742
2743	nc_stat = NF90_DEF_VAR( id_set_sp, 'k_x', NF90_DOUBLE, id_dim_x_sp, &
2744	id_var_x_sp )
2745	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 280 )
2746
2747	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_x_sp, 'units', 'm-1' )
2748	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 281 )
2749
2750	!
2751	!-- Define y-axis
2752	nc_stat = NF90_DEF_DIM( id_set_sp, 'k_y', ny/2, id_dim_y_sp )
2753	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 282 )
2754
2755	nc_stat = NF90_DEF_VAR( id_set_sp, 'k_y', NF90_DOUBLE, id_dim_y_sp, &
2756	id_var_y_sp )
2757	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 283 )
2758
2759	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_y_sp, 'units', 'm-1' )
2760	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 284 )
2761
2762	!
2763	!-- Define the variables
2764	var_list = ';'
2765	i = 1
2766	DO WHILE ( data_output_sp(i) /= ' ' .AND. i <= 10 )
2767
2768	IF ( INDEX( spectra_direction(i), 'x' ) /= 0 ) THEN
2769
2770	!
2771	!-- Define the variable
2772	netcdf_var_name = TRIM( data_output_sp(i) ) // '_x'
2773	IF ( data_output_sp(i) == 'w' ) THEN
2774	nc_stat = NF90_DEF_VAR( id_set_sp, netcdf_var_name, &
2775	nc_precision(7), (/ id_dim_x_sp, &
2776	id_dim_zw_sp, id_dim_time_sp /), &
2777	id_var_dospx(i) )
2778	ELSE
2779	nc_stat = NF90_DEF_VAR( id_set_sp, netcdf_var_name, &
2780	nc_precision(7), (/ id_dim_x_sp, &
2781	id_dim_zu_sp, id_dim_time_sp /), &
2782	id_var_dospx(i) )
2783	ENDIF
2784	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 285 )
2785
2786	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_dospx(i), &
2787	'long_name', netcdf_var_name )
2788	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 286 )
2789
2790	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_dospx(i), &
2791	'units', 'unknown' )
2792	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 287 )
2793
2794	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
2795
2796	ENDIF
2797
2798	IF ( INDEX( spectra_direction(i), 'y' ) /= 0 ) THEN
2799
2800	!
2801	!-- Define the variable
2802	netcdf_var_name = TRIM( data_output_sp(i) ) // '_y'
2803	IF ( data_output_sp(i) == 'w' ) THEN
2804	nc_stat = NF90_DEF_VAR( id_set_sp, netcdf_var_name, &
2805	nc_precision(7), (/ id_dim_y_sp, &
2806	id_dim_zw_sp, id_dim_time_sp /), &
2807	id_var_dospy(i) )
2808	ELSE
2809	nc_stat = NF90_DEF_VAR( id_set_sp, netcdf_var_name, &
2810	nc_precision(7), (/ id_dim_y_sp, &
2811	id_dim_zu_sp, id_dim_time_sp /), &
2812	id_var_dospy(i) )
2813	ENDIF
2814	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 288 )
2815
2816	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_dospy(i), &
2817	'long_name', netcdf_var_name )
2818	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 289 )
2819
2820	nc_stat = NF90_PUT_ATT( id_set_sp, id_var_dospy(i), &
2821	'units', 'unknown' )
2822	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 290 )
2823
2824	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
2825
2826	ENDIF
2827
2828	i = i + 1
2829
2830	ENDDO
2831
2832	!
2833	!-- Write the list of variables as global attribute (this is used by
2834	!-- restart runs)
2835	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'VAR_LIST', var_list )
2836	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 291 )
2837
2838	!
2839	!-- Leave NetCDF define mode
2840	nc_stat = NF90_ENDDEF( id_set_sp )
2841	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 292 )
2842
2843	!
2844	!-- Write axis data: zu_sp, zw_sp, k_x, k_y
2845	ALLOCATE( netcdf_data(1:ns) )
2846
2847	!
2848	!-- Write zu data
2849	netcdf_data(1:ns) = zu( comp_spectra_level(1:ns) )
2850	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_zu_sp, netcdf_data, &
2851	start = (/ 1 /), count = (/ ns /) )
2852	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 293 )
2853
2854	!
2855	!-- Write zw data
2856	netcdf_data(1:ns) = zw( comp_spectra_level(1:ns) )
2857	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_zw_sp, netcdf_data, &
2858	start = (/ 1 /), count = (/ ns /) )
2859	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 294 )
2860
2861	DEALLOCATE( netcdf_data )
2862
2863	!
2864	!-- Write data for x and y axis (wavenumbers)
2865	ALLOCATE( netcdf_data(nx/2) )
2866	DO i = 1, nx/2
2867	netcdf_data(i) = 2.0 * pi * i / ( dx * ( nx + 1 ) )
2868	ENDDO
2869
2870	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_x_sp, netcdf_data, &
2871	start = (/ 1 /), count = (/ nx/2 /) )
2872	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 295 )
2873
2874	DEALLOCATE( netcdf_data )
2875
2876	ALLOCATE( netcdf_data(ny/2) )
2877	DO i = 1, ny/2
2878	netcdf_data(i) = 2.0 * pi * i / ( dy * ( ny + 1 ) )
2879	ENDDO
2880
2881	nc_stat = NF90_PUT_VAR( id_set_sp, id_var_y_sp, netcdf_data, &
2882	start = (/ 1 /), count = (/ ny/2 /) )
2883	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 296 )
2884
2885	DEALLOCATE( netcdf_data )
2886
2887
2888	CASE ( 'sp_ext' )
2889
2890	!
2891	!-- Get the list of variables and compare with the actual run.
2892	!-- First var_list_old has to be reset, since GET_ATT does not assign
2893	!-- trailing blanks.
2894	var_list_old = ' '
2895	nc_stat = NF90_GET_ATT( id_set_sp, NF90_GLOBAL, 'VAR_LIST', &
2896	var_list_old )
2897	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 297 )
2898
2899	var_list = ';'
2900	i = 1
2901	DO WHILE ( data_output_sp(i) /= ' ' .AND. i <= 10 )
2902
2903	IF ( INDEX( spectra_direction(i), 'x' ) /= 0 ) THEN
2904	netcdf_var_name = TRIM( data_output_sp(i) ) // '_x'
2905	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
2906	ENDIF
2907
2908	IF ( INDEX( spectra_direction(i), 'y' ) /= 0 ) THEN
2909	netcdf_var_name = TRIM( data_output_sp(i) ) // '_y'
2910	var_list = TRIM( var_list ) // TRIM( netcdf_var_name ) // ';'
2911	ENDIF
2912
2913	i = i + 1
2914
2915	ENDDO
2916
2917	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
2918	PRINT*, '+++ WARNING: NetCDF file for spectra from previous ' //&
2919	'run found,'
2920	PRINT*, ' but this file cannot be extended due to' // &
2921	' variable mismatch.'
2922	PRINT*, ' New file is created instead.'
2923	extend = .FALSE.
2924	RETURN
2925	ENDIF
2926
2927	!
2928	!-- Determine the number of current vertical levels for which spectra
2929	!-- shall be output
2930	ns = 1
2931	DO WHILE ( comp_spectra_level(ns) /= 999999 .AND. ns <= 100 )
2932	ns = ns + 1
2933	ENDDO
2934	ns = ns - 1
2935
2936	!
2937	!-- Get and compare the number of vertical levels
2938	nc_stat = NF90_INQ_VARID( id_set_sp, 'zu_sp', id_var_zu_sp )
2939	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 298 )
2940
2941	nc_stat = NF90_INQUIRE_VARIABLE( id_set_sp, id_var_zu_sp, &
2942	dimids = id_dim_zu_sp_old )
2943	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 299 )
2944	id_dim_zu_sp = id_dim_zu_sp_old(1)
2945
2946	nc_stat = NF90_INQUIRE_DIMENSION( id_set_sp, id_dim_zu_sp, &
2947	len = ns_old )
2948	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 300 )
2949
2950	IF ( ns /= ns_old ) THEN
2951	PRINT*, '+++ WARNING: NetCDF file for spectra from previous ' //&
2952	'run found,'
2953	PRINT*, ' but this file cannot be extended due to' // &
2954	' mismatch in number of'
2955	PRINT*, ' vertical levels.'
2956	PRINT*, ' New file is created instead.'
2957	extend = .FALSE.
2958	RETURN
2959	ENDIF
2960
2961	!
2962	!-- Get and compare the heights of the cross sections
2963	ALLOCATE( netcdf_data(1:ns_old) )
2964
2965	nc_stat = NF90_GET_VAR( id_set_sp, id_var_zu_sp, netcdf_data )
2966	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 301 )
2967
2968	DO i = 1, ns
2969	IF ( zu(comp_spectra_level(i)) /= netcdf_data(i) ) THEN
2970	PRINT*, '+++ WARNING: NetCDF file for spectra from previou' // &
2971	's run found,'
2972	PRINT*, ' but this file cannot be extended due' // &
2973	' to mismatch in heights'
2974	PRINT*, ' of vertical levels.'
2975	PRINT*, ' New file is created instead.'
2976	extend = .FALSE.
2977	RETURN
2978	ENDIF
2979	ENDDO
2980
2981	DEALLOCATE( netcdf_data )
2982
2983	!
2984	!-- Get the id of the time coordinate (unlimited coordinate) and its
2985	!-- last index on the file. The next time level is plsp..count+1.
2986	!-- The current time must be larger than the last output time
2987	!-- on the file.
2988	nc_stat = NF90_INQ_VARID( id_set_sp, 'time', id_var_time_sp )
2989	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 302 )
2990
2991	nc_stat = NF90_INQUIRE_VARIABLE( id_set_sp, id_var_time_sp, &
2992	dimids = id_dim_time_old )
2993	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 303 )
2994	id_dim_time_sp = id_dim_time_old(1)
2995
2996	nc_stat = NF90_INQUIRE_DIMENSION( id_set_sp, id_dim_time_sp, &
2997	len = dosp_time_count )
2998	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 304 )
2999
3000	nc_stat = NF90_GET_VAR( id_set_sp, id_var_time_sp, &
3001	last_time_coordinate, &
3002	start = (/ dosp_time_count /), &
3003	count = (/ 1 /) )
3004	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 305 )
3005
3006	IF ( last_time_coordinate(1) >= simulated_time ) THEN
3007	PRINT*, '+++ WARNING: NetCDF file for spectra from previous ' // &
3008	'run found,'
3009	PRINT*, ' but this file cannot be extended becaus' // &
3010	'e the current output time'
3011	PRINT*, ' is less or equal than the last output t' // &
3012	'ime on this file.'
3013	PRINT*, ' New file is created instead.'
3014	dosp_time_count = 0
3015	extend = .FALSE.
3016	RETURN
3017	ENDIF
3018
3019	!
3020	!-- Dataset seems to be extendable.
3021	!-- Now get the variable ids.
3022	i = 1
3023	DO WHILE ( data_output_sp(i) /= ' ' .AND. i <= 10 )
3024
3025	IF ( INDEX( spectra_direction(i), 'x' ) /= 0 ) THEN
3026	netcdf_var_name = TRIM( data_output_sp(i) ) // '_x'
3027	nc_stat = NF90_INQ_VARID( id_set_sp, netcdf_var_name, &
3028	id_var_dospx(i) )
3029	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 306 )
3030	ENDIF
3031
3032	IF ( INDEX( spectra_direction(i), 'y' ) /= 0 ) THEN
3033	netcdf_var_name = TRIM( data_output_sp(i) ) // '_y'
3034	nc_stat = NF90_INQ_VARID( id_set_sp, netcdf_var_name, &
3035	id_var_dospy(i) )
3036	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 307 )
3037	ENDIF
3038
3039	i = i + 1
3040
3041	ENDDO
3042
3043	!
3044	!-- Change the titel attribute on file
3045	IF ( averaging_interval_sp /= 0.0 ) THEN
3046	WRITE (time_average_text,'('', '',F7.1,'' s average'')') &
3047	averaging_interval_sp
3048	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
3049	TRIM( run_description_header ) // &
3050	TRIM( time_average_text ) )
3051	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 308 )
3052
3053	WRITE ( time_average_text,'(F7.1,'' s avg'')' ) averaging_interval_sp
3054	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'time_avg', &
3055	TRIM( time_average_text ) )
3056	ELSE
3057	nc_stat = NF90_PUT_ATT( id_set_sp, NF90_GLOBAL, 'title', &
3058	TRIM( run_description_header ) )
3059	ENDIF
3060	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 309 )
3061
3062	PRINT, '** NetCDF file for spectra from previous run found.'
3063	PRINT*, ' This file will be extended.'
3064
3065
3066	CASE ( 'pt_new' )
3067
3068	!
3069	!-- Define some global attributes of the dataset
3070	nc_stat = NF90_PUT_ATT( id_set_prt, NF90_GLOBAL, 'title', &
3071	TRIM( run_description_header ) )
3072	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 310 )
3073
3074	!
3075	!-- Define time coordinate for particles (unlimited dimension)
3076	nc_stat = NF90_DEF_DIM( id_set_prt, 'time', NF90_UNLIMITED, &
3077	id_dim_time_prt )
3078	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 311 )
3079
3080	nc_stat = NF90_DEF_VAR( id_set_prt, 'time', NF90_DOUBLE, &
3081	id_dim_time_prt, id_var_time_prt )
3082	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 312 )
3083
3084	nc_stat = NF90_PUT_ATT( id_set_prt, id_var_time_prt, 'units', &
3085	'seconds' )
3086	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 313 )
3087
3088	!
3089	!-- Define particle coordinate (maximum particle number)
3090	nc_stat = NF90_DEF_DIM( id_set_prt, 'particle_number', &
3091	maximum_number_of_particles, id_dim_prtnum )
3092	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 314 )
3093
3094	nc_stat = NF90_DEF_VAR( id_set_prt, 'particle_number', NF90_DOUBLE, &
3095	id_dim_prtnum, id_var_prtnum )
3096	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 315 )
3097
3098	nc_stat = NF90_PUT_ATT( id_set_prt, id_var_prtnum, 'units', &
3099	'particle number' )
3100	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 316 )
3101
3102	!
3103	!-- Define variable which contains the real number of particles in use
3104	nc_stat = NF90_DEF_VAR( id_set_prt, 'real_num_of_prt', NF90_INT, &
3105	id_dim_time_prt, id_var_rnop_prt )
3106	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 317 )
3107
3108	nc_stat = NF90_PUT_ATT( id_set_prt, id_var_rnop_prt, 'units', &
3109	'particle number' )
3110	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 318 )
3111
3112	!
3113	!-- Define the variables
3114	DO i = 1, 17
3115
3116	nc_stat = NF90_DEF_VAR( id_set_prt, prt_var_names(i), &
3117	nc_precision(8), &
3118	(/ id_dim_prtnum, id_dim_time_prt /), &
3119	id_var_prt(i) )
3120	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 319 )
3121
3122	nc_stat = NF90_PUT_ATT( id_set_prt, id_var_prt(i), &
3123	'long_name', TRIM( prt_var_names(i) ) )
3124	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 320 )
3125
3126	nc_stat = NF90_PUT_ATT( id_set_prt, id_var_prt(i), &
3127	'units', TRIM( prt_var_units(i) ) )
3128	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 321 )
3129
3130	ENDDO
3131
3132	!
3133	!-- Leave NetCDF define mode
3134	nc_stat = NF90_ENDDEF( id_set_prt )
3135	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 322 )
3136
3137
3138	CASE ( 'pt_ext' )
3139
3140	!
3141	!-- Get the id of the time coordinate (unlimited coordinate) and its
3142	!-- last index on the file. The next time level is prt..count+1.
3143	!-- The current time must be larger than the last output time
3144	!-- on the file.
3145	nc_stat = NF90_INQ_VARID( id_set_prt, 'time', id_var_time_prt )
3146	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 323 )
3147
3148	nc_stat = NF90_INQUIRE_VARIABLE( id_set_prt, id_var_time_prt, &
3149	dimids = id_dim_time_old )
3150	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 324 )
3151	id_dim_time_prt = id_dim_time_old(1)
3152
3153	nc_stat = NF90_INQUIRE_DIMENSION( id_set_prt, id_dim_time_prt, &
3154	len = prt_time_count )
3155	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 325 )
3156
3157	nc_stat = NF90_GET_VAR( id_set_prt, id_var_time_prt, &
3158	last_time_coordinate, &
3159	start = (/ prt_time_count /), &
3160	count = (/ 1 /) )
3161	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 326 )
3162
3163	IF ( last_time_coordinate(1) >= simulated_time ) THEN
3164	PRINT*, '+++ WARNING: NetCDF file for particles from previous ' //&
3165	'run found,'
3166	PRINT*, ' but this file cannot be extended becaus' // &
3167	'e the current output time'
3168	PRINT*, ' is less or equal than the last output t' // &
3169	'ime on this file.'
3170	PRINT*, ' New file is created instead.'
3171	prt_time_count = 0
3172	extend = .FALSE.
3173	RETURN
3174	ENDIF
3175
3176	!
3177	!-- Dataset seems to be extendable.
3178	!-- Now get the variable ids.
3179	nc_stat = NF90_INQ_VARID( id_set_prt, 'real_num_of_prt', &
3180	id_var_rnop_prt )
3181	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 327 )
3182
3183	DO i = 1, 17
3184
3185	nc_stat = NF90_INQ_VARID( id_set_prt, prt_var_names(i), &
3186	id_var_prt(i) )
3187	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 328 )
3188
3189	ENDDO
3190
3191	IF ( myid == 0 ) THEN
3192	PRINT, '** NetCDF file for particles from previous run found.'
3193	PRINT*, ' This file will be extended.'
3194	ENDIF
3195
3196
3197	CASE ( 'ps_new' )
3198
3199	!
3200	!-- Define some global attributes of the dataset
3201	nc_stat = NF90_PUT_ATT( id_set_pts, NF90_GLOBAL, 'title', &
3202	TRIM( run_description_header ) )
3203	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 396 )
3204
3205	!
3206	!-- Define time coordinate for particle time series (unlimited dimension)
3207	nc_stat = NF90_DEF_DIM( id_set_pts, 'time', NF90_UNLIMITED, &
3208	id_dim_time_pts )
3209	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 397 )
3210
3211	nc_stat = NF90_DEF_VAR( id_set_pts, 'time', NF90_DOUBLE, &
3212	id_dim_time_pts, id_var_time_pts )
3213	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 398 )
3214
3215	nc_stat = NF90_PUT_ATT( id_set_pts, id_var_time_pts, 'units', &
3216	'seconds')
3217	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 399 )
3218
3219	!
3220	!-- Define the variables. If more than one particle group is defined,
3221	!-- define seperate variables for each group
3222	var_list = ';'
3223	DO i = 1, dopts_num
3224
3225	!
3226	!-- First, remove those characters not allowed by NetCDF
3227	netcdf_var_name = dopts_label(i)
3228	CALL clean_netcdf_varname( netcdf_var_name )
3229
3230	DO j = 0, number_of_particle_groups
3231
3232	IF ( j == 0 ) THEN
3233	suffix1 = ''
3234	ELSE
3235	WRITE ( suffix1, '(''_'',I2.2)' ) j
3236	ENDIF
3237
3238	nc_stat = NF90_DEF_VAR( id_set_pts, &
3239	TRIM( netcdf_var_name ) // suffix1, &
3240	nc_precision(6), &
3241	(/ id_dim_time_pts /), &
3242	id_var_dopts(i,j) )
3243	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 400 )
3244
3245	IF ( j == 0 ) THEN
3246	nc_stat = NF90_PUT_ATT( id_set_pts, id_var_dopts(i,j), &
3247	'long_name', &
3248	TRIM( dopts_label(i) ) )
3249	ELSE
3250	nc_stat = NF90_PUT_ATT( id_set_pts, id_var_dopts(i,j), &
3251	'long_name', &
3252	TRIM( dopts_label(i) ) // ' PG ' // &
3253	suffix1(2:3) )
3254	ENDIF
3255	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 401 )
3256
3257	nc_stat = NF90_PUT_ATT( id_set_pts, id_var_dopts(i,j), &
3258	'units', TRIM( dopts_unit(i) ) )
3259	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 402 )
3260
3261	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // &
3262	suffix1 // ';'
3263
3264	IF ( number_of_particle_groups == 1 ) EXIT
3265
3266	ENDDO
3267
3268	ENDDO
3269
3270	!
3271	!-- Write the list of variables as global attribute (this is used by
3272	!-- restart runs)
3273	nc_stat = NF90_PUT_ATT( id_set_pts, NF90_GLOBAL, 'VAR_LIST', &
3274	var_list )
3275	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 403 )
3276
3277	!
3278	!-- Leave NetCDF define mode
3279	nc_stat = NF90_ENDDEF( id_set_pts )
3280	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 404 )
3281
3282
3283	CASE ( 'ps_ext' )
3284
3285	!
3286	!-- Get the list of variables and compare with the actual run.
3287	!-- First var_list_old has to be reset, since GET_ATT does not assign
3288	!-- trailing blanks.
3289	var_list_old = ' '
3290	nc_stat = NF90_GET_ATT( id_set_pts, NF90_GLOBAL, 'VAR_LIST', &
3291	var_list_old )
3292	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 405 )
3293
3294	var_list = ';'
3295	i = 1
3296	DO i = 1, dopts_num
3297
3298	netcdf_var_name = dopts_label(i)
3299	CALL clean_netcdf_varname( netcdf_var_name )
3300
3301	DO j = 0, number_of_particle_groups
3302
3303	IF ( j == 0 ) THEN
3304	suffix1 = ''
3305	ELSE
3306	WRITE ( suffix1, '(''_'',I2.2)' ) j
3307	ENDIF
3308
3309	var_list = TRIM(var_list) // TRIM(netcdf_var_name) // &
3310	suffix1 // ';'
3311
3312	IF ( number_of_particle_groups == 1 ) EXIT
3313
3314	ENDDO
3315
3316	ENDDO
3317
3318	IF ( TRIM( var_list ) /= TRIM( var_list_old ) ) THEN
3319	PRINT*, '+++ WARNING: NetCDF file for particle time series ' //&
3320	'from previuos run found,'
3321	PRINT*, ' but this file cannot be extended due to' // &
3322	' variable mismatch.'
3323	PRINT*, ' New file is created instead.'
3324	extend = .FALSE.
3325	RETURN
3326	ENDIF
3327
3328	!
3329	!-- Get the id of the time coordinate (unlimited coordinate) and its
3330	!-- last index on the file. The next time level is dots..count+1.
3331	!-- The current time must be larger than the last output time
3332	!-- on the file.
3333	nc_stat = NF90_INQ_VARID( id_set_pts, 'time', id_var_time_pts )
3334	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 406 )
3335
3336	nc_stat = NF90_INQUIRE_VARIABLE( id_set_pts, id_var_time_pts, &
3337	dimids = id_dim_time_old )
3338	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 407 )
3339	id_dim_time_pts = id_dim_time_old(1)
3340
3341	nc_stat = NF90_INQUIRE_DIMENSION( id_set_pts, id_dim_time_pts, &
3342	len = dopts_time_count )
3343	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 408 )
3344
3345	nc_stat = NF90_GET_VAR( id_set_pts, id_var_time_pts, &
3346	last_time_coordinate, &
3347	start = (/ dopts_time_count /), &
3348	count = (/ 1 /) )
3349	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 409 )
3350
3351	IF ( last_time_coordinate(1) >= simulated_time ) THEN
3352	PRINT*, '+++ WARNING: NetCDF file for time series from' // &
3353	' previuos run found,'
3354	PRINT*, ' but this file cannot be extended becaus' // &
3355	'e the current output time'
3356	PRINT*, ' is less or equal than the last output t' // &
3357	'ime on this file.'
3358	PRINT*, ' New file is created instead.'
3359	dopts_time_count = 0
3360	extend = .FALSE.
3361	RETURN
3362	ENDIF
3363
3364	!
3365	!-- Dataset seems to be extendable.
3366	!-- Now get the variable ids
3367	i = 1
3368	DO i = 1, dopts_num
3369
3370	netcdf_var_name_base = dopts_label(i)
3371	CALL clean_netcdf_varname( netcdf_var_name_base )
3372
3373	DO j = 0, number_of_particle_groups
3374
3375	IF ( j == 0 ) THEN
3376	suffix1 = ''
3377	ELSE
3378	WRITE ( suffix1, '(''_'',I2.2)' ) j
3379	ENDIF
3380
3381	netcdf_var_name = TRIM( netcdf_var_name_base ) // suffix1
3382
3383	nc_stat = NF90_INQ_VARID( id_set_pts, netcdf_var_name, &
3384	id_var_dopts(i,j) )
3385	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 410 )
3386
3387	IF ( number_of_particle_groups == 1 ) EXIT
3388
3389	ENDDO
3390
3391	ENDDO
3392
3393	!
3394	!-- Change the title attribute on file
3395	nc_stat = NF90_PUT_ATT( id_set_pts, NF90_GLOBAL, 'title', &
3396	TRIM( run_description_header ) )
3397	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 411 )
3398	PRINT, '** NetCDF file for particle time series from previous ', &
3399	'run found.'
3400	PRINT*, ' This file will be extended.'
3401
3402
3403	CASE DEFAULT
3404
3405	PRINT*, '+++ define_netcdf_header: mode "', mode, '" not supported'
3406
3407	END SELECT
3408
3409	#endif
3410	END SUBROUTINE define_netcdf_header
3411
3412
3413
3414	SUBROUTINE handle_netcdf_error( errno )
3415	#if defined( __netcdf )
3416
3417	!------------------------------------------------------------------------------!
3418	!
3419	! Description:
3420	! ------------
3421	! Prints out a text message corresponding to the current status.
3422	!------------------------------------------------------------------------------!
3423
3424	USE netcdf
3425	USE netcdf_control
3426	USE pegrid
3427
3428	IMPLICIT NONE
3429
3430	INTEGER :: errno
3431
3432	IF ( nc_stat /= NF90_NOERR ) THEN
3433	PRINT*, '+++ netcdf error ', errno,': ', TRIM( NF90_STRERROR( nc_stat ) )
3434	#if defined( __parallel )
3435	CALL MPI_ABORT( comm2d, 9999, ierr )
3436	#else
3437	CALL local_stop
3438	#endif
3439	ENDIF
3440
3441	#endif
3442	END SUBROUTINE handle_netcdf_error
3443
3444
3445
3446	SUBROUTINE clean_netcdf_varname( string )
3447	#if defined( __netcdf )
3448
3449	!------------------------------------------------------------------------------!
3450	!
3451	! Description:
3452	! ------------
3453	! Replace those characters in string which are not allowed by NetCDF.
3454	!------------------------------------------------------------------------------!
3455
3456	USE netcdf_control
3457
3458	IMPLICIT NONE
3459
3460	CHARACTER (LEN=10), INTENT(INOUT) :: string
3461
3462	INTEGER :: i, ic
3463
3464	DO i = 1, replace_num
3465	DO
3466	ic = INDEX( string, replace_char(i) )
3467	IF ( ic /= 0 ) THEN
3468	string(ic:ic) = replace_by(i)
3469	ELSE
3470	EXIT
3471	ENDIF
3472	ENDDO
3473	ENDDO
3474
3475	#endif
3476	END SUBROUTINE clean_netcdf_varname

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |