Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |