Home

Context Navigation

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

Last change on this file since 254 was 226, checked in by raasch, 16 years ago
preparations for the next release
Property svn:keywords set to `Id`
File size: 143.9 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |