Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |