Home

Context Navigation

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

Last change on this file since 1552 was 1552, checked in by maronga, 9 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 210.4 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |