Home

Context Navigation

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

Last change on this file since 269 was 263, checked in by heinze, 16 years ago
Output of NetCDF messages with aid of message handling routine.
Property svn:keywords set to `Id`
File size: 142.5 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |