Home

Context Navigation

source: palm/trunk/SOURCE/netcdf_data_input_mod.f90 @ 3036

Last change on this file since 3036 was 3036, checked in by gronemeier, 6 years ago
renamed input variables according to UC2 data standard
Property svn:keywords set to `Id`
File size: 210.3 KB

Line
1	!> @file netcdf_data_input_mod.f90
2	!------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the
6	! terms of the GNU General Public License as published by the Free Software
7	! Foundation, either version 3 of the License, or (at your option) any later
8	! version.
9	!
10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
13	!
14	! You should have received a copy of the GNU General Public License along with
15	! PALM. If not, see <http://www.gnu.org/licenses/>.
16	!
17	! Copyright 1997-2018 Leibniz Universitaet Hannover
18	!------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! -----------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: netcdf_data_input_mod.f90 3036 2018-05-24 10:18:26Z gronemeier $
27	! Revision of input vars according to UC2 data standard
28	! - renamed 'orography_2D' to 'zt'
29	! - renamed 'buildings_2D' to 'buildings_2d'
30	! - renamed 'buildings_3D' to 'buildings_3d'
31	! - renamed 'leaf_are_density' to 'lad'
32	! - renamed 'basal_are_density' to 'bad'
33	! - renamed 'root_are_density_lad' to 'root_area_dens_r'
34	! - renamed 'root_are_density_lsm' to 'root_area_dens_s'
35	! - renamed 'ls_forcing_ug' to 'tend_ug'
36	! - renamed 'ls_forcing_vg' to 'tend_vg'
37	!
38	! 3019 2018-05-13 07:05:43Z maronga
39	! Improved reading speed of large NetCDF files
40	!
41	! 2963 2018-04-12 14:47:44Z suehring
42	! - Revise checks for static input variables.
43	! - Introduce index for vegetation/wall, pavement/green-wall and water/window
44	! surfaces, for clearer access of surface fraction, albedo, emissivity, etc. .
45	!
46	! 2958 2018-04-11 15:38:13Z suehring
47	! Synchronize longitude and latitude between nested model domains, values are
48	! taken from the root model.
49	!
50	! 2955 2018-04-09 15:14:01Z suehring
51	! Extend checks for consistent setting of buildings, its ID and type.
52	! Add log-points to measure CPU time of NetCDF data input.
53	!
54	! 2953 2018-04-09 11:26:02Z suehring
55	! Bugfix in checks for initialization data
56	!
57	! 2947 2018-04-04 18:01:41Z suehring
58	! Checks for dynamic input revised
59	!
60	! 2946 2018-04-04 17:01:23Z suehring
61	! Bugfix for revision 2945, perform checks only if dynamic input file is
62	! available.
63	!
64	! 2945 2018-04-04 16:27:14Z suehring
65	! - Mimic for topography input slightly revised, in order to enable consistency
66	! checks
67	! - Add checks for dimensions in dynamic input file and move already existing
68	! checks
69	!
70	! 2938 2018-03-27 15:52:42Z suehring
71	! Initial read of geostrophic wind components from dynamic driver.
72	!
73	! 2773 2018-01-30 14:12:54Z suehring
74	! Revise checks for surface_fraction.
75	!
76	! 2925 2018-03-23 14:54:11Z suehring
77	! Check for further inconsistent settings of surface_fractions.
78	! Some messages slightly rephrased and error numbers renamed.
79	!
80	! 2898 2018-03-15 13:03:01Z suehring
81	! Check if each building has a type. Further, check if dimensions in static
82	! input file match the model dimensions.
83	!
84	! 2897 2018-03-15 11:47:16Z suehring
85	! Relax restrictions for topography input, terrain and building heights can be
86	! input separately and are not mandatory any more.
87	!
88	! 2874 2018-03-13 10:55:42Z knoop
89	! Bugfix: wrong placement of netcdf cpp-macros fixed
90	!
91	! 2794 2018-02-07 14:09:43Z knoop
92	! Check if 3D building input is consistent to numeric grid.
93	!
94	! 2773 2018-01-30 14:12:54Z suehring
95	! - Enable initialization with 3D topography.
96	! - Move check for correct initialization in nesting mode to check_parameters.
97	!
98	! 2772 2018-01-29 13:10:35Z suehring
99	! Initialization of simulation independent on land-surface model.
100	!
101	! 2746 2018-01-15 12:06:04Z suehring
102	! Read plant-canopy variables independently on land-surface model usage
103	!
104	! 2718 2018-01-02 08:49:38Z maronga
105	! Corrected "Former revisions" section
106	!
107	! 2711 2017-12-20 17:04:49Z suehring
108	! Rename subroutine close_file to avoid double-naming.
109	!
110	! 2700 2017-12-15 14:12:35Z suehring
111	!
112	! 2696 2017-12-14 17:12:51Z kanani
113	! Initial revision (suehring)
114	!
115	!
116	!
117	!
118	! Authors:
119	! --------
120	! @author Matthias Suehring
121	!
122	! Description:
123	! ------------
124	!> Modulue contains routines to input data according to Palm input data
125	!> standart using dynamic and static input files.
126	!>
127	!> @todo - Order input alphabetically
128	!> @todo - Revise error messages and error numbers
129	!> @todo - Input of missing quantities (chemical species, emission rates)
130	!> @todo - Defninition and input of still missing variable attributes
131	!> @todo - Input of initial geostrophic wind profiles with cyclic conditions.
132	!------------------------------------------------------------------------------!
133	MODULE netcdf_data_input_mod
134
135	USE control_parameters, &
136	ONLY: coupling_char, io_blocks, io_group
137
138	USE cpulog, &
139	ONLY: cpu_log, log_point_s
140
141	USE kinds
142
143	#if defined ( __netcdf )
144	USE NETCDF
145	#endif
146
147	USE pegrid
148
149	USE surface_mod, &
150	ONLY: ind_pav_green, ind_veg_wall, ind_wat_win
151	!
152	!-- Define type for dimensions.
153	TYPE dims_xy
154	INTEGER(iwp) :: nx !< dimension length in x
155	INTEGER(iwp) :: ny !< dimension length in y
156	INTEGER(iwp) :: nz !< dimension length in z
157	REAL(wp), DIMENSION(:), ALLOCATABLE :: x !< dimension array in x
158	REAL(wp), DIMENSION(:), ALLOCATABLE :: y !< dimension array in y
159	REAL(wp), DIMENSION(:), ALLOCATABLE :: z !< dimension array in z
160	END TYPE dims_xy
161	!
162	!-- Define data type for nesting in larger-scale models like COSMO.
163	!-- Data type comprises u, v, w, pt, and q at lateral and top boundaries.
164	TYPE force_type
165
166	CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE :: var_names
167
168	INTEGER(iwp) :: nt !< number of time levels in dynamic input file
169	INTEGER(iwp) :: nzu !< number of vertical levels on scalar grid in dynamic input file
170	INTEGER(iwp) :: nzw !< number of vertical levels on w grid in dynamic input file
171	INTEGER(iwp) :: tind !< time index for reference time in large-scale forcing data
172	INTEGER(iwp) :: tind_p !< time index for following time in large-scale forcing data
173
174	LOGICAL :: init = .FALSE.
175	LOGICAL :: interpolated = .FALSE.
176	LOGICAL :: from_file = .FALSE.
177
178	REAL(wp), DIMENSION(:), ALLOCATABLE :: surface_pressure !< time dependent surface pressure
179	REAL(wp), DIMENSION(:), ALLOCATABLE :: time !< time levels in dynamic input file
180	REAL(wp), DIMENSION(:), ALLOCATABLE :: zu_atmos !< vertical levels at scalar grid in dynamic input file
181	REAL(wp), DIMENSION(:), ALLOCATABLE :: zw_atmos !< vertical levels at w grid in dynamic input file
182
183	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ug !< domain-averaged geostrophic component
184	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vg !< domain-averaged geostrophic component
185
186	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: u_left !< u-component at left boundary
187	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: v_left !< v-component at left boundary
188	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: w_left !< w-component at left boundary
189	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: q_left !< mixing ratio at left boundary
190	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: pt_left !< potentital temperautre at left boundary
191
192	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: u_north !< u-component at north boundary
193	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: v_north !< v-component at north boundary
194	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: w_north !< w-component at north boundary
195	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: q_north !< mixing ratio at north boundary
196	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: pt_north !< potentital temperautre at north boundary
197
198	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: u_right !< u-component at right boundary
199	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: v_right !< v-component at right boundary
200	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: w_right !< w-component at right boundary
201	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: q_right !< mixing ratio at right boundary
202	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: pt_right !< potentital temperautre at right boundary
203
204	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: u_south !< u-component at south boundary
205	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: v_south !< v-component at south boundary
206	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: w_south !< w-component at south boundary
207	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: q_south !< mixing ratio at south boundary
208	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: pt_south !< potentital temperautre at south boundary
209
210	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: u_top !< u-component at top boundary
211	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: v_top !< v-component at top boundary
212	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: w_top !< w-component at top boundary
213	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: q_top !< mixing ratio at top boundary
214	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: pt_top !< potentital temperautre at top boundary
215
216	END TYPE force_type
217
218	TYPE init_type
219
220	INTEGER(iwp) :: lod_msoil !< level of detail - soil moisture
221	INTEGER(iwp) :: lod_pt !< level of detail - pt
222	INTEGER(iwp) :: lod_q !< level of detail - q
223	INTEGER(iwp) :: lod_tsoil !< level of detail - soil temperature
224	INTEGER(iwp) :: lod_u !< level of detail - u-component
225	INTEGER(iwp) :: lod_v !< level of detail - v-component
226	INTEGER(iwp) :: lod_w !< level of detail - w-component
227	INTEGER(iwp) :: nx !< number of scalar grid points along x in dynamic input file
228	INTEGER(iwp) :: nxu !< number of u grid points along x in dynamic input file
229	INTEGER(iwp) :: ny !< number of scalar grid points along y in dynamic input file
230	INTEGER(iwp) :: nyv !< number of v grid points along y in dynamic input file
231	INTEGER(iwp) :: nzs !< number of vertical soil levels in dynamic input file
232	INTEGER(iwp) :: nzu !< number of vertical levels on scalar grid in dynamic input file
233	INTEGER(iwp) :: nzw !< number of vertical levels on w grid in dynamic input file
234
235	LOGICAL :: from_file_msoil = .FALSE. !< flag indicating whether soil moisture is already initialized from file
236	LOGICAL :: from_file_pt = .FALSE. !< flag indicating whether pt is already initialized from file
237	LOGICAL :: from_file_q = .FALSE. !< flag indicating whether q is already initialized from file
238	LOGICAL :: from_file_tsoil = .FALSE. !< flag indicating whether soil temperature is already initialized from file
239	LOGICAL :: from_file_u = .FALSE. !< flag indicating whether u is already initialized from file
240	LOGICAL :: from_file_ug = .FALSE. !< flag indicating whether ug is already initialized from file
241	LOGICAL :: from_file_v = .FALSE. !< flag indicating whether v is already initialized from file
242	LOGICAL :: from_file_vg = .FALSE. !< flag indicating whether ug is already initialized from file
243	LOGICAL :: from_file_w = .FALSE. !< flag indicating whether w is already initialized from file
244
245
246	REAL(wp) :: fill_msoil !< fill value for soil moisture
247	REAL(wp) :: fill_pt !< fill value for pt
248	REAL(wp) :: fill_q !< fill value for q
249	REAL(wp) :: fill_tsoil !< fill value for soil temperature
250	REAL(wp) :: fill_u !< fill value for u
251	REAL(wp) :: fill_v !< fill value for v
252	REAL(wp) :: fill_w !< fill value for w
253	REAL(wp) :: latitude !< latitude of the southern model boundary
254	REAL(wp) :: longitude !< longitude of the western model boundary
255
256	REAL(wp), DIMENSION(:), ALLOCATABLE :: msoil_init !< initial vertical profile of soil moisture
257	REAL(wp), DIMENSION(:), ALLOCATABLE :: pt_init !< initial vertical profile of pt
258	REAL(wp), DIMENSION(:), ALLOCATABLE :: q_init !< initial vertical profile of q
259	REAL(wp), DIMENSION(:), ALLOCATABLE :: tsoil_init !< initial vertical profile of soil temperature
260	REAL(wp), DIMENSION(:), ALLOCATABLE :: u_init !< initial vertical profile of u
261	REAL(wp), DIMENSION(:), ALLOCATABLE :: ug_init !< initial vertical profile of ug
262	REAL(wp), DIMENSION(:), ALLOCATABLE :: v_init !< initial vertical profile of v
263	REAL(wp), DIMENSION(:), ALLOCATABLE :: vg_init !< initial vertical profile of ug
264	REAL(wp), DIMENSION(:), ALLOCATABLE :: w_init !< initial vertical profile of w
265	REAL(wp), DIMENSION(:), ALLOCATABLE :: z_soil !< vertical levels in soil in dynamic input file, used for interpolation
266	REAL(wp), DIMENSION(:), ALLOCATABLE :: zu_atmos !< vertical levels at scalar grid in dynamic input file, used for interpolation
267	REAL(wp), DIMENSION(:), ALLOCATABLE :: zw_atmos !< vertical levels at w grid in dynamic input file, used for interpolation
268
269
270	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: msoil !< initial 3d soil moisture provide by Inifor and interpolated onto soil grid
271	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: tsoil !< initial 3d soil temperature provide by Inifor and interpolated onto soil grid
272
273	END TYPE init_type
274
275	!
276	!-- Define data structures for different input data types.
277	!-- 8-bit Integer 2D
278	TYPE int_2d_8bit
279	INTEGER(KIND=1) :: fill = -127 !< fill value
280	INTEGER(KIND=1), DIMENSION(:,:), ALLOCATABLE :: var !< respective variable
281
282	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
283	END TYPE int_2d_8bit
284	!
285	!-- 32-bit Integer 2D
286	TYPE int_2d_32bit
287	INTEGER(iwp) :: fill = -9999 !< fill value
288	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: var !< respective variable
289
290	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
291	END TYPE int_2d_32bit
292
293	!
294	!-- Define data type to read 2D real variables
295	TYPE real_2d
296	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
297
298	REAL(wp) :: fill = -9999.9_wp !< fill value
299	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: var !< respective variable
300	END TYPE real_2d
301
302	!
303	!-- Define data type to read 2D real variables
304	TYPE real_3d
305	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
306
307	INTEGER(iwp) :: nz !< number of grid points along vertical dimension
308
309	REAL(wp) :: fill = -9999.9_wp !< fill value
310	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: var !< respective variable
311	END TYPE real_3d
312	!
313	!-- Define data structure where the dimension and type of the input depends
314	!-- on the given level of detail.
315	!-- For buildings, the input is either 2D float, or 3d byte.
316	TYPE build_in
317	INTEGER(iwp) :: lod = 1 !< level of detail
318	INTEGER(KIND=1) :: fill2 = -127 !< fill value for lod = 2
319	INTEGER(iwp) :: nz !< number of vertical layers in file
320	INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE :: var_3d !< 3d variable (lod = 2)
321
322	REAL(wp), DIMENSION(:), ALLOCATABLE :: z !< vertical coordinate for 3D building, used for consistency check
323
324	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
325
326	REAL(wp) :: fill1 = -9999.9_wp !< fill values for lod = 1
327	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: var_2d !< 2d variable (lod = 1)
328	END TYPE build_in
329
330	!
331	!-- For soil_type, the input is either 2D or 3D one-byte integer.
332	TYPE soil_in
333	INTEGER(iwp) :: lod = 1 !< level of detail
334	INTEGER(KIND=1) :: fill = -127 !< fill value for lod = 2
335	INTEGER(KIND=1), DIMENSION(:,:), ALLOCATABLE :: var_2d !< 2d variable (lod = 1)
336	INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE :: var_3d !< 3d variable (lod = 2)
337
338	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
339	END TYPE soil_in
340
341	!
342	!-- Define data type for fractions between surface types
343	TYPE fracs
344	INTEGER(iwp) :: nf !< total number of fractions
345	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nfracs !< dimension array for fraction
346
347	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
348
349	REAL(wp) :: fill = -9999.9_wp !< fill value
350	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: frac !< respective fraction between different surface types
351	END TYPE fracs
352	!
353	!-- Data type for parameter lists, Depending on the given level of detail,
354	!-- the input is 3D or 4D
355	TYPE pars
356	INTEGER(iwp) :: lod = 1 !< level of detail
357	INTEGER(iwp) :: np !< total number of parameters
358	INTEGER(iwp) :: nz !< vertical dimension - number of soil layers
359	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: layers !< dimension array for soil layers
360	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: pars !< dimension array for parameters
361
362	LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used
363
364	REAL(wp) :: fill = -9999.9_wp !< fill value
365	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: pars_xy !< respective parameters, level of detail = 1
366	REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE :: pars_xyz !< respective parameters, level of detail = 2
367	END TYPE pars
368	!
369	!-- Define variables
370	TYPE(dims_xy) :: dim_static !< data structure for x, y-dimension in static input file
371
372	TYPE(force_type) :: force !< data structure for data input at lateral and top boundaries (provided by Inifor)
373
374	TYPE(init_type) :: init_3d !< data structure for the initialization of the 3D flow and soil fields
375	TYPE(init_type) :: init_model !< data structure for the initialization of the model
376
377	!
378	!-- Define 2D variables of type NC_BYTE
379	TYPE(int_2d_8bit) :: albedo_type_f !< input variable for albedo type
380	TYPE(int_2d_8bit) :: building_type_f !< input variable for building type
381	TYPE(int_2d_8bit) :: pavement_type_f !< input variable for pavenment type
382	TYPE(int_2d_8bit) :: street_crossing_f !< input variable for water type
383	TYPE(int_2d_8bit) :: street_type_f !< input variable for water type
384	TYPE(int_2d_8bit) :: vegetation_type_f !< input variable for vegetation type
385	TYPE(int_2d_8bit) :: water_type_f !< input variable for water type
386
387	!
388	!-- Define 2D variables of type NC_INT
389	TYPE(int_2d_32bit) :: building_id_f !< input variable for building ID
390	!
391	!-- Define 2D variables of type NC_FLOAT
392	TYPE(real_2d) :: terrain_height_f !< input variable for terrain height
393	!
394	!-- Define 3D variables of type NC_FLOAT
395	TYPE(real_3d) :: basal_area_density_f !< input variable for basal area density - resolved vegetation
396	TYPE(real_3d) :: leaf_area_density_f !< input variable for leaf area density - resolved vegetation
397	TYPE(real_3d) :: root_area_density_lad_f !< input variable for root area density - resolved vegetation
398	TYPE(real_3d) :: root_area_density_lsm_f !< input variable for root area density - parametrized vegetation
399
400	!
401	!-- Define input variable for buildings
402	TYPE(build_in) :: buildings_f !< input variable for buildings
403	!
404	!-- Define input variables for soil_type
405	TYPE(soil_in) :: soil_type_f !< input variable for soil type
406
407	TYPE(fracs) :: surface_fraction_f !< input variable for surface fraction
408
409	TYPE(pars) :: albedo_pars_f !< input variable for albedo parameters
410	TYPE(pars) :: building_pars_f !< input variable for building parameters
411	TYPE(pars) :: pavement_pars_f !< input variable for pavement parameters
412	TYPE(pars) :: pavement_subsurface_pars_f !< input variable for pavement parameters
413	TYPE(pars) :: soil_pars_f !< input variable for soil parameters
414	TYPE(pars) :: vegetation_pars_f !< input variable for vegetation parameters
415	TYPE(pars) :: water_pars_f !< input variable for water parameters
416
417
418	CHARACTER(LEN=3) :: char_lod = 'lod' !< name of level-of-detail attribute in NetCDF file
419
420	CHARACTER(LEN=10) :: char_fill = '_FillValue' !< name of fill value attribute in NetCDF file
421	CHARACTER(LEN=10) :: char_lon = 'origin_lon' !< name of global attribute for longitude in NetCDF file
422	CHARACTER(LEN=10) :: char_lat = 'origin_lat' !< name of global attribute for latitude in NetCDF file
423
424	CHARACTER(LEN=100) :: input_file_static = 'PIDS_STATIC' !< Name of file which comprises static input data
425	CHARACTER(LEN=100) :: input_file_dynamic = 'PIDS_DYNAMIC' !< Name of file which comprises dynamic input data
426
427	INTEGER(iwp) :: nc_stat !< return value of nf90 function call
428
429	LOGICAL :: input_pids_static = .FALSE. !< Flag indicating whether Palm-input-data-standard file containing static information exists
430	LOGICAL :: input_pids_dynamic = .FALSE. !< Flag indicating whether Palm-input-data-standard file containing dynamic information exists
431
432	LOGICAL :: collective_read = .FALSE. !< Enable NetCDF collective read
433
434	SAVE
435
436	PRIVATE
437
438	INTERFACE netcdf_data_input_interpolate
439	MODULE PROCEDURE netcdf_data_input_interpolate_1d
440	MODULE PROCEDURE netcdf_data_input_interpolate_1d_soil
441	MODULE PROCEDURE netcdf_data_input_interpolate_2d
442	MODULE PROCEDURE netcdf_data_input_interpolate_3d
443	END INTERFACE netcdf_data_input_interpolate
444
445	INTERFACE netcdf_data_input_check_dynamic
446	MODULE PROCEDURE netcdf_data_input_check_dynamic
447	END INTERFACE netcdf_data_input_check_dynamic
448
449	INTERFACE netcdf_data_input_check_static
450	MODULE PROCEDURE netcdf_data_input_check_static
451	END INTERFACE netcdf_data_input_check_static
452
453	INTERFACE netcdf_data_input_inquire_file
454	MODULE PROCEDURE netcdf_data_input_inquire_file
455	END INTERFACE netcdf_data_input_inquire_file
456
457	INTERFACE netcdf_data_input_init
458	MODULE PROCEDURE netcdf_data_input_init
459	END INTERFACE netcdf_data_input_init
460
461	INTERFACE netcdf_data_input_init_3d
462	MODULE PROCEDURE netcdf_data_input_init_3d
463	END INTERFACE netcdf_data_input_init_3d
464
465	INTERFACE netcdf_data_input_lsf
466	MODULE PROCEDURE netcdf_data_input_lsf
467	END INTERFACE netcdf_data_input_lsf
468
469	INTERFACE netcdf_data_input_surface_data
470	MODULE PROCEDURE netcdf_data_input_surface_data
471	END INTERFACE netcdf_data_input_surface_data
472
473	INTERFACE netcdf_data_input_topo
474	MODULE PROCEDURE netcdf_data_input_topo
475	END INTERFACE netcdf_data_input_topo
476
477	INTERFACE get_variable
478	MODULE PROCEDURE get_variable_1d_int
479	MODULE PROCEDURE get_variable_1d_real
480	MODULE PROCEDURE get_variable_2d_int8
481	MODULE PROCEDURE get_variable_2d_int32
482	MODULE PROCEDURE get_variable_2d_real
483	MODULE PROCEDURE get_variable_3d_int8
484	MODULE PROCEDURE get_variable_3d_real
485	MODULE PROCEDURE get_variable_3d_real_v
486	MODULE PROCEDURE get_variable_4d_real
487	END INTERFACE get_variable
488
489	INTERFACE get_variable_pr
490	MODULE PROCEDURE get_variable_pr
491	END INTERFACE get_variable_pr
492
493	INTERFACE get_attribute
494	MODULE PROCEDURE get_attribute_real
495	MODULE PROCEDURE get_attribute_int8
496	MODULE PROCEDURE get_attribute_int32
497	MODULE PROCEDURE get_attribute_string
498	END INTERFACE get_attribute
499
500	!
501	!-- Public variables
502	PUBLIC albedo_pars_f, albedo_type_f, basal_area_density_f, buildings_f, &
503	building_id_f, building_pars_f, building_type_f, force, init_3d, &
504	init_model, input_file_static, input_pids_static, &
505	input_pids_dynamic, leaf_area_density_f, &
506	pavement_pars_f, pavement_subsurface_pars_f, pavement_type_f, &
507	root_area_density_lad_f, root_area_density_lsm_f, soil_pars_f, &
508	soil_type_f, street_crossing_f, street_type_f, surface_fraction_f, &
509	terrain_height_f, vegetation_pars_f, vegetation_type_f, &
510	water_pars_f, water_type_f
511
512	!
513	!-- Public subroutines
514	PUBLIC netcdf_data_input_check_dynamic, netcdf_data_input_check_static, &
515	netcdf_data_input_inquire_file, &
516	netcdf_data_input_init, netcdf_data_input_init_3d, &
517	netcdf_data_input_interpolate, netcdf_data_input_lsf, &
518	netcdf_data_input_surface_data, netcdf_data_input_topo
519
520	CONTAINS
521
522	!------------------------------------------------------------------------------!
523	! Description:
524	! ------------
525	!> Inquires whether NetCDF input files according to Palm-input-data standard
526	!> exist. Moreover, basic checks are performed.
527	!------------------------------------------------------------------------------!
528	SUBROUTINE netcdf_data_input_inquire_file
529
530	USE control_parameters, &
531	ONLY: land_surface, message_string, topo_no_distinct, urban_surface
532
533	IMPLICIT NONE
534
535	LOGICAL :: check_nest !< flag indicating whether a check passed or not
536
537	#if defined ( __netcdf )
538	INQUIRE( FILE = TRIM( input_file_static ) // TRIM( coupling_char ), &
539	EXIST = input_pids_static )
540	INQUIRE( FILE = TRIM( input_file_dynamic ) // TRIM( coupling_char ), &
541	EXIST = input_pids_dynamic )
542	#endif
543
544	!
545	!-- As long as topography can be input via ASCII format, no distinction
546	!-- between building and terrain can be made. This case, classify all
547	!-- surfaces as default type. Same in case land-surface and urban-surface
548	!-- model are not applied.
549	IF ( .NOT. input_pids_static ) THEN
550	topo_no_distinct = .TRUE.
551	ENDIF
552
553	END SUBROUTINE netcdf_data_input_inquire_file
554
555	!------------------------------------------------------------------------------!
556	! Description:
557	! ------------
558	!> Reads global attributes required for initialization of the model.
559	!------------------------------------------------------------------------------!
560	SUBROUTINE netcdf_data_input_init
561
562	IMPLICIT NONE
563
564	INTEGER(iwp) :: id_mod !< NetCDF id of input file
565	INTEGER(iwp) :: ii !< running index for IO blocks
566
567	IF ( .NOT. input_pids_static ) RETURN
568
569	DO ii = 0, io_blocks-1
570	IF ( ii == io_group ) THEN
571	#if defined ( __netcdf )
572	!
573	!-- Open file in read-only mode
574	CALL open_read_file( TRIM( input_file_static ) // &
575	TRIM( coupling_char ), id_mod )
576	!
577	!-- Read global attribute for latitude and longitude
578	CALL get_attribute( id_mod, char_lat, &
579	init_model%latitude, .TRUE. )
580
581	CALL get_attribute( id_mod, char_lon, &
582	init_model%longitude, .TRUE. )
583	!
584	!-- Finally, close input file
585	CALL close_input_file( id_mod )
586	#endif
587	ENDIF
588	#if defined( __parallel )
589	CALL MPI_BARRIER( comm2d, ierr )
590	#endif
591	ENDDO
592	!
593	!-- In case of nested runs, each model domain might have different longitude
594	!-- and latitude, which would result in different Coriolis parameters and
595	!-- sun-zenith angles. To avoid this, longitude and latitude in each model
596	!-- domain will be set to the values of the root model. Please note, this
597	!-- synchronization is required already here.
598	#if defined( __parallel )
599	CALL MPI_BCAST( init_model%latitude, 1, MPI_REAL, 0, &
600	MPI_COMM_WORLD, ierr )
601	CALL MPI_BCAST( init_model%longitude, 1, MPI_REAL, 0, &
602	MPI_COMM_WORLD, ierr )
603	#endif
604
605
606	END SUBROUTINE netcdf_data_input_init
607
608	!------------------------------------------------------------------------------!
609	! Description:
610	! ------------
611	!> Reads surface classification data, such as vegetation and soil type, etc. .
612	!------------------------------------------------------------------------------!
613	SUBROUTINE netcdf_data_input_surface_data
614
615	USE control_parameters, &
616	ONLY: bc_lr_cyc, bc_ns_cyc, land_surface, message_string, &
617	plant_canopy, urban_surface
618
619	USE indices, &
620	ONLY: nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg
621
622
623	IMPLICIT NONE
624
625	CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE :: var_names !< variable names in static input file
626
627	INTEGER(iwp) :: i !< running index along x-direction
628	INTEGER(iwp) :: ii !< running index for IO blocks
629	INTEGER(iwp) :: id_surf !< NetCDF id of input file
630	INTEGER(iwp) :: j !< running index along y-direction
631	INTEGER(iwp) :: k !< running index along z-direction
632	INTEGER(iwp) :: k2 !< running index
633	INTEGER(iwp) :: num_vars !< number of variables in input file
634	INTEGER(iwp) :: nz_soil !< number of soil layers in file
635
636	INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg) :: var_exchange_int !< dummy variables used to exchange 32-bit Integer arrays
637	INTEGER(iwp), DIMENSION(:,:,:), ALLOCATABLE :: var_dum_int_3d !< dummy variables used to exchange real arrays
638
639	REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) :: var_exchange_real !< dummy variables used to exchange real arrays
640
641	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: var_dum_real_3d !< dummy variables used to exchange real arrays
642	REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE :: var_dum_real_4d !< dummy variables used to exchange real arrays
643
644	!
645	!-- If not static input file is available, skip this routine
646	IF ( .NOT. input_pids_static ) RETURN
647	!
648	!-- Measure CPU time
649	CALL cpu_log( log_point_s(82), 'NetCDF input', 'start' )
650	!
651	!-- Read plant canopy variables.
652	IF ( plant_canopy ) THEN
653	DO ii = 0, io_blocks-1
654	IF ( ii == io_group ) THEN
655	#if defined ( __netcdf )
656	!
657	!-- Open file in read-only mode
658	CALL open_read_file( TRIM( input_file_static ) // &
659	TRIM( coupling_char ) , id_surf )
660	!
661	!-- At first, inquire all variable names.
662	!-- This will be used to check whether an optional input variable
663	!-- exist or not.
664	CALL inquire_num_variables( id_surf, num_vars )
665
666	ALLOCATE( var_names(1:num_vars) )
667	CALL inquire_variable_names( id_surf, var_names )
668
669	!
670	!-- Read leaf area density - resolved vegetation
671	IF ( check_existence( var_names, 'lad' ) ) THEN
672	leaf_area_density_f%from_file = .TRUE.
673	CALL get_attribute( id_surf, char_fill, &
674	leaf_area_density_f%fill, &
675	.FALSE., 'lad' )
676	!
677	!-- Inquire number of vertical vegetation layer
678	CALL get_dimension_length( id_surf, leaf_area_density_f%nz, &
679	'zlad' )
680	!
681	!-- Allocate variable for leaf-area density
682	ALLOCATE( leaf_area_density_f%var( &
683	0:leaf_area_density_f%nz-1, &
684	nys:nyn,nxl:nxr) )
685
686	CALL get_variable( id_surf, 'lad', &
687	nxl, nxr, nys, nyn, &
688	leaf_area_density_f%var(:,nys:nyn, nxl:nxr) )
689
690	ELSE
691	leaf_area_density_f%from_file = .FALSE.
692	ENDIF
693
694	!
695	!-- Read basal area density - resolved vegetation
696	IF ( check_existence( var_names, 'bad' ) ) THEN
697	basal_area_density_f%from_file = .TRUE.
698	CALL get_attribute( id_surf, char_fill, &
699	basal_area_density_f%fill, &
700	.FALSE., 'bad' )
701	!
702	!-- Inquire number of vertical vegetation layer
703	CALL get_dimension_length( id_surf, &
704	basal_area_density_f%nz, &
705	'zlad' )
706	!
707	!-- Allocate variable
708	ALLOCATE( basal_area_density_f%var( &
709	0:basal_area_density_f%nz-1, &
710	nys:nyn,nxl:nxr) )
711
712	CALL get_variable( id_surf, 'bad', &
713	nxl, nxr, nys, nyn, &
714	basal_area_density_f%var(:,nys:nyn, nxl:nxr) )
715	ELSE
716	basal_area_density_f%from_file = .FALSE.
717	ENDIF
718
719	!
720	!-- Read root area density - resolved vegetation
721	IF ( check_existence( var_names, 'root_area_dens_r' ) ) THEN
722	root_area_density_lad_f%from_file = .TRUE.
723	CALL get_attribute( id_surf, char_fill, &
724	root_area_density_lad_f%fill, &
725	.FALSE., 'root_area_dens_r' )
726	!
727	!-- Inquire number of vertical soil layers
728	CALL get_dimension_length( id_surf, &
729	root_area_density_lad_f%nz, &
730	'zsoil' )
731	!
732	!-- Allocate variable
733	ALLOCATE( root_area_density_lad_f%var &
734	(0:root_area_density_lad_f%nz-1,&
735	nys:nyn,nxl:nxr) )
736
737	CALL get_variable( id_surf, 'root_area_dens_r', &
738	nxl, nxr, nys, nyn, &
739	root_area_density_lad_f%var(:,nys:nyn, nxl:nxr) )
740	ELSE
741	root_area_density_lad_f%from_file = .FALSE.
742	ENDIF
743	!
744	!-- Finally, close input file
745	CALL close_input_file( id_surf )
746	#endif
747	ENDIF
748	#if defined( __parallel )
749	CALL MPI_BARRIER( comm2d, ierr )
750	#endif
751	ENDDO
752	!
753	!-- Deallocate variable list. Will be re-allocated in case further
754	!-- variables are read from file.
755	IF ( ALLOCATED( var_names ) ) DEALLOCATE( var_names )
756
757	ENDIF
758	!
759	!-- Skip the following if no land-surface or urban-surface module are
760	!-- applied. This case, no one of the following variables is used anyway.
761	IF ( .NOT. land_surface .OR. .NOT. urban_surface ) RETURN
762	!
763	!-- Initialize dummy arrays used for ghost-point exchange
764	var_exchange_int = 0
765	var_exchange_real = 0.0_wp
766
767	DO ii = 0, io_blocks-1
768	IF ( ii == io_group ) THEN
769	#if defined ( __netcdf )
770	!
771	!-- Open file in read-only mode
772	CALL open_read_file( TRIM( input_file_static ) // &
773	TRIM( coupling_char ) , id_surf )
774
775	!
776	!-- Inquire all variable names.
777	!-- This will be used to check whether an optional input variable exist
778	!-- or not.
779	CALL inquire_num_variables( id_surf, num_vars )
780
781	ALLOCATE( var_names(1:num_vars) )
782	CALL inquire_variable_names( id_surf, var_names )
783
784	!
785	!-- Read vegetation type and required attributes
786	IF ( check_existence( var_names, 'vegetation_type' ) ) THEN
787	vegetation_type_f%from_file = .TRUE.
788	CALL get_attribute( id_surf, char_fill, &
789	vegetation_type_f%fill, &
790	.FALSE., 'vegetation_type' )
791	!
792	!-- PE-wise reading of 2D vegetation type.
793	ALLOCATE ( vegetation_type_f%var(nys:nyn,nxl:nxr) )
794
795	DO i = nxl, nxr
796	CALL get_variable( id_surf, 'vegetation_type', &
797	i, vegetation_type_f%var(:,i) )
798	ENDDO
799	ELSE
800	vegetation_type_f%from_file = .FALSE.
801	ENDIF
802
803	!
804	!-- Read soil type and required attributes
805	IF ( check_existence( var_names, 'soil_type' ) ) THEN
806	soil_type_f%from_file = .TRUE.
807	!
808	!-- Note, lod is currently not on file; skip for the moment
809	! CALL get_attribute( id_surf, char_lod, &
810	! soil_type_f%lod, &
811	! .FALSE., 'soil_type' )
812	CALL get_attribute( id_surf, char_fill, &
813	soil_type_f%fill, &
814	.FALSE., 'soil_type' )
815
816	IF ( soil_type_f%lod == 1 ) THEN
817	!
818	!-- PE-wise reading of 2D soil type.
819	ALLOCATE ( soil_type_f%var_2d(nys:nyn,nxl:nxr) )
820	DO i = nxl, nxr
821	CALL get_variable( id_surf, 'soil_type', &
822	i, soil_type_f%var_2d(:,i) )
823	ENDDO
824	ELSEIF ( soil_type_f%lod == 2 ) THEN
825	!
826	!-- Obtain number of soil layers from file.
827	CALL get_dimension_length( id_surf, nz_soil, 'zsoil' )
828	!
829	!-- PE-wise reading of 3D soil type.
830	ALLOCATE ( soil_type_f%var_3d(0:nz_soil,nys:nyn,nxl:nxr) )
831	DO i = nxl, nxr
832	DO j = nys, nyn
833	CALL get_variable( id_surf, 'soil_type', i, j, &
834	soil_type_f%var_3d(:,j,i) )
835	ENDDO
836	ENDDO
837	ENDIF
838	ELSE
839	soil_type_f%from_file = .FALSE.
840	ENDIF
841
842	!
843	!-- Read pavement type and required attributes
844	IF ( check_existence( var_names, 'pavement_type' ) ) THEN
845	pavement_type_f%from_file = .TRUE.
846	CALL get_attribute( id_surf, char_fill, &
847	pavement_type_f%fill, .FALSE., &
848	'pavement_type' )
849	!
850	!-- PE-wise reading of 2D pavement type.
851	ALLOCATE ( pavement_type_f%var(nys:nyn,nxl:nxr) )
852	DO i = nxl, nxr
853	CALL get_variable( id_surf, 'pavement_type', &
854	i, pavement_type_f%var(:,i) )
855	ENDDO
856	ELSE
857	pavement_type_f%from_file = .FALSE.
858	ENDIF
859
860	!
861	!-- Read water type and required attributes
862	IF ( check_existence( var_names, 'water_type' ) ) THEN
863	water_type_f%from_file = .TRUE.
864	CALL get_attribute( id_surf, char_fill, water_type_f%fill, &
865	.FALSE., 'water_type' )
866	!
867	!-- PE-wise reading of 2D water type.
868	ALLOCATE ( water_type_f%var(nys:nyn,nxl:nxr) )
869	DO i = nxl, nxr
870	CALL get_variable( id_surf, 'water_type', i, &
871	water_type_f%var(:,i) )
872	ENDDO
873	ELSE
874	water_type_f%from_file = .FALSE.
875	ENDIF
876	!
877	!-- Read surface fractions and related information
878	IF ( check_existence( var_names, 'surface_fraction' ) ) THEN
879	surface_fraction_f%from_file = .TRUE.
880	CALL get_attribute( id_surf, char_fill, &
881	surface_fraction_f%fill, &
882	.FALSE., 'surface_fraction' )
883	!
884	!-- Inquire number of surface fractions
885	CALL get_dimension_length( id_surf, &
886	surface_fraction_f%nf, &
887	'nsurface_fraction' )
888	!
889	!-- Allocate dimension array and input array for surface fractions
890	ALLOCATE( surface_fraction_f%nfracs(0:surface_fraction_f%nf-1) )
891	ALLOCATE( surface_fraction_f%frac(0:surface_fraction_f%nf-1, &
892	nys:nyn,nxl:nxr) )
893	!
894	!-- Get dimension of surface fractions
895	CALL get_variable( id_surf, 'nsurface_fraction', &
896	surface_fraction_f%nfracs )
897	!
898	!-- Read surface fractions
899	CALL get_variable( id_surf, 'surface_fraction', &
900	nxl, nxr, nys, nyn, &
901	surface_fraction_f%frac(:,nys:nyn, nxl:nxr))
902	ELSE
903	surface_fraction_f%from_file = .FALSE.
904	ENDIF
905	!
906	!-- Read building parameters and related information
907	IF ( check_existence( var_names, 'building_pars' ) ) THEN
908	building_pars_f%from_file = .TRUE.
909	CALL get_attribute( id_surf, char_fill, &
910	building_pars_f%fill, &
911	.FALSE., 'building_pars' )
912	!
913	!-- Inquire number of building parameters
914	CALL get_dimension_length( id_surf, &
915	building_pars_f%np, &
916	'nbuilding_pars' )
917	!
918	!-- Allocate dimension array and input array for building parameters
919	ALLOCATE( building_pars_f%pars(0:building_pars_f%np-1) )
920	ALLOCATE( building_pars_f%pars_xy(0:building_pars_f%np-1, &
921	nys:nyn,nxl:nxr) )
922	!
923	!-- Get dimension of building parameters
924	CALL get_variable( id_surf, 'nbuilding_pars', &
925	building_pars_f%pars )
926	!
927	!-- Read building_pars
928	CALL get_variable( id_surf, 'building_pars', &
929	nxl, nxr, nys, nyn, &
930	building_pars_f%pars_xy(:,nys:nyn, nxl:nxr) )
931	ELSE
932	building_pars_f%from_file = .FALSE.
933	ENDIF
934
935	!
936	!-- Read albedo type and required attributes
937	IF ( check_existence( var_names, 'albedo_type' ) ) THEN
938	albedo_type_f%from_file = .TRUE.
939	CALL get_attribute( id_surf, char_fill, albedo_type_f%fill, &
940	.FALSE., 'albedo_type' )
941	!
942	!-- PE-wise reading of 2D water type.
943	ALLOCATE ( albedo_type_f%var(nys:nyn,nxl:nxr) )
944	DO i = nxl, nxr
945	CALL get_variable( id_surf, 'albedo_type', &
946	i, albedo_type_f%var(:,i) )
947	ENDDO
948	ELSE
949	albedo_type_f%from_file = .FALSE.
950	ENDIF
951	!
952	!-- Read albedo parameters and related information
953	IF ( check_existence( var_names, 'albedo_pars' ) ) THEN
954	albedo_pars_f%from_file = .TRUE.
955	CALL get_attribute( id_surf, char_fill, albedo_pars_f%fill, &
956	.FALSE., 'albedo_pars' )
957	!
958	!-- Inquire number of albedo parameters
959	CALL get_dimension_length( id_surf, albedo_pars_f%np, &
960	'nalbedo_pars' )
961	!
962	!-- Allocate dimension array and input array for albedo parameters
963	ALLOCATE( albedo_pars_f%pars(0:albedo_pars_f%np-1) )
964	ALLOCATE( albedo_pars_f%pars_xy(0:albedo_pars_f%np-1, &
965	nys:nyn,nxl:nxr) )
966	!
967	!-- Get dimension of albedo parameters
968	CALL get_variable( id_surf, 'nalbedo_pars', albedo_pars_f%pars )
969
970	DO i = nxl, nxr
971	DO j = nys, nyn
972	CALL get_variable( id_surf, 'albedo_pars', i, j, &
973	albedo_pars_f%pars_xy(:,j,i) )
974	ENDDO
975	ENDDO
976	ELSE
977	albedo_pars_f%from_file = .FALSE.
978	ENDIF
979
980	!
981	!-- Read pavement parameters and related information
982	IF ( check_existence( var_names, 'pavement_pars' ) ) THEN
983	pavement_pars_f%from_file = .TRUE.
984	CALL get_attribute( id_surf, char_fill, &
985	pavement_pars_f%fill, &
986	.FALSE., 'pavement_pars' )
987	!
988	!-- Inquire number of pavement parameters
989	CALL get_dimension_length( id_surf, pavement_pars_f%np, &
990	'npavement_pars' )
991	!
992	!-- Allocate dimension array and input array for pavement parameters
993	ALLOCATE( pavement_pars_f%pars(0:pavement_pars_f%np-1) )
994	ALLOCATE( pavement_pars_f%pars_xy(0:pavement_pars_f%np-1, &
995	nys:nyn,nxl:nxr) )
996	!
997	!-- Get dimension of pavement parameters
998	CALL get_variable( id_surf, 'npavement_pars', &
999	pavement_pars_f%pars )
1000
1001	DO i = nxl, nxr
1002	DO j = nys, nyn
1003	CALL get_variable( id_surf, 'pavement_pars', i, j, &
1004	pavement_pars_f%pars_xy(:,j,i) )
1005	ENDDO
1006	ENDDO
1007	ELSE
1008	pavement_pars_f%from_file = .FALSE.
1009	ENDIF
1010
1011	!
1012	!-- Read pavement subsurface parameters and related information
1013	IF ( check_existence( var_names, 'pavement_subsurface_pars' ) ) &
1014	THEN
1015	pavement_subsurface_pars_f%from_file = .TRUE.
1016	CALL get_attribute( id_surf, char_fill, &
1017	pavement_subsurface_pars_f%fill, &
1018	.FALSE., 'pavement_subsurface_pars' )
1019	!
1020	!-- Inquire number of parameters
1021	CALL get_dimension_length( id_surf, &
1022	pavement_subsurface_pars_f%np, &
1023	'npavement_subsurface_pars' )
1024	!
1025	!-- Inquire number of soil layers
1026	CALL get_dimension_length( id_surf, &
1027	pavement_subsurface_pars_f%nz, &
1028	'zsoil' )
1029	!
1030	!-- Allocate dimension array and input array for pavement parameters
1031	ALLOCATE( pavement_subsurface_pars_f%pars &
1032	(0:pavement_subsurface_pars_f%np-1) )
1033	ALLOCATE( pavement_subsurface_pars_f%pars_xyz &
1034	(0:pavement_subsurface_pars_f%np-1, &
1035	0:pavement_subsurface_pars_f%nz-1, &
1036	nys:nyn,nxl:nxr) )
1037	!
1038	!-- Get dimension of pavement parameters
1039	CALL get_variable( id_surf, 'npavement_subsurface_pars', &
1040	pavement_subsurface_pars_f%pars )
1041
1042	DO i = nxl, nxr
1043	DO j = nys, nyn
1044	CALL get_variable( &
1045	id_surf, 'pavement_subsurface_pars', &
1046	i, j, &
1047	pavement_subsurface_pars_f%pars_xyz(:,:,j,i),&
1048	pavement_subsurface_pars_f%nz, &
1049	pavement_subsurface_pars_f%np )
1050	ENDDO
1051	ENDDO
1052	ELSE
1053	pavement_subsurface_pars_f%from_file = .FALSE.
1054	ENDIF
1055
1056
1057	!
1058	!-- Read vegetation parameters and related information
1059	IF ( check_existence( var_names, 'vegetation_pars' ) ) THEN
1060	vegetation_pars_f%from_file = .TRUE.
1061	CALL get_attribute( id_surf, char_fill, &
1062	vegetation_pars_f%fill, &
1063	.FALSE., 'vegetation_pars' )
1064	!
1065	!-- Inquire number of vegetation parameters
1066	CALL get_dimension_length( id_surf, vegetation_pars_f%np, &
1067	'nvegetation_pars' )
1068	!
1069	!-- Allocate dimension array and input array for surface fractions
1070	ALLOCATE( vegetation_pars_f%pars(0:vegetation_pars_f%np-1) )
1071	ALLOCATE( vegetation_pars_f%pars_xy(0:vegetation_pars_f%np-1, &
1072	nys:nyn,nxl:nxr) )
1073	!
1074	!-- Get dimension of the parameters
1075	CALL get_variable( id_surf, 'nvegetation_pars', &
1076	vegetation_pars_f%pars )
1077
1078	CALL get_variable( id_surf, 'vegetation_pars', &
1079	nxl, nxr, nys, nyn, &
1080	vegetation_pars_f%pars_xy(:,nys:nyn,nxl:nxr ) )
1081	ELSE
1082	vegetation_pars_f%from_file = .FALSE.
1083	ENDIF
1084
1085	!
1086	!-- Read root parameters/distribution and related information
1087	IF ( check_existence( var_names, 'soil_pars' ) ) THEN
1088	soil_pars_f%from_file = .TRUE.
1089	CALL get_attribute( id_surf, char_fill, &
1090	soil_pars_f%fill, &
1091	.FALSE., 'soil_pars' )
1092
1093	CALL get_attribute( id_surf, char_lod, &
1094	soil_pars_f%lod, &
1095	.FALSE., 'soil_pars' )
1096
1097	!
1098	!-- Inquire number of soil parameters
1099	CALL get_dimension_length( id_surf, &
1100	soil_pars_f%np, &
1101	'nsoil_pars' )
1102	!
1103	!-- Read parameters array
1104	ALLOCATE( soil_pars_f%pars(0:soil_pars_f%np-1) )
1105	CALL get_variable( id_surf, 'nsoil_pars', soil_pars_f%pars )
1106
1107	!
1108	!-- In case of level of detail 2, also inquire number of vertical
1109	!-- soil layers, allocate memory and read the respective dimension
1110	IF ( soil_pars_f%lod == 2 ) THEN
1111	CALL get_dimension_length( id_surf, soil_pars_f%nz, 'zsoil' )
1112
1113	ALLOCATE( soil_pars_f%layers(0:soil_pars_f%nz-1) )
1114	CALL get_variable( id_surf, 'zsoil', soil_pars_f%layers )
1115
1116	ENDIF
1117
1118	!
1119	!-- Read soil parameters, depending on level of detail
1120	IF ( soil_pars_f%lod == 1 ) THEN
1121	ALLOCATE( soil_pars_f%pars_xy(0:soil_pars_f%np-1, &
1122	nys:nyn,nxl:nxr) )
1123	DO i = nxl, nxr
1124	DO j = nys, nyn
1125	CALL get_variable( id_surf, 'soil_pars', i, j, &
1126	soil_pars_f%pars_xy(:,j,i) )
1127	ENDDO
1128	ENDDO
1129
1130	ELSEIF ( soil_pars_f%lod == 2 ) THEN
1131	ALLOCATE( soil_pars_f%pars_xyz(0:soil_pars_f%np-1, &
1132	0:soil_pars_f%nz-1, &
1133	nys:nyn,nxl:nxr) )
1134	DO i = nxl, nxr
1135	DO j = nys, nyn
1136	CALL get_variable( id_surf, 'soil_pars', i, j, &
1137	soil_pars_f%pars_xyz(:,:,j,i), &
1138	soil_pars_f%nz, soil_pars_f%np )
1139	ENDDO
1140	ENDDO
1141	ENDIF
1142	ELSE
1143	soil_pars_f%from_file = .FALSE.
1144	ENDIF
1145
1146	!
1147	!-- Read water parameters and related information
1148	IF ( check_existence( var_names, 'water_pars' ) ) THEN
1149	water_pars_f%from_file = .TRUE.
1150	CALL get_attribute( id_surf, char_fill, &
1151	water_pars_f%fill, &
1152	.FALSE., 'water_pars' )
1153	!
1154	!-- Inquire number of water parameters
1155	CALL get_dimension_length( id_surf, &
1156	water_pars_f%np, &
1157	'nwater_pars' )
1158	!
1159	!-- Allocate dimension array and input array for water parameters
1160	ALLOCATE( water_pars_f%pars(0:water_pars_f%np-1) )
1161	ALLOCATE( water_pars_f%pars_xy(0:water_pars_f%np-1, &
1162	nys:nyn,nxl:nxr) )
1163	!
1164	!-- Get dimension of water parameters
1165	CALL get_variable( id_surf, 'nwater_pars', water_pars_f%pars )
1166
1167	DO i = nxl, nxr
1168	DO j = nys, nyn
1169	CALL get_variable( id_surf, 'water_pars', i, j, &
1170	water_pars_f%pars_xy(:,j,i) )
1171	ENDDO
1172	ENDDO
1173	ELSE
1174	water_pars_f%from_file = .FALSE.
1175	ENDIF
1176	!
1177	!-- Read root area density - parametrized vegetation
1178	IF ( check_existence( var_names, 'root_area_dens_s' ) ) THEN
1179	root_area_density_lsm_f%from_file = .TRUE.
1180	CALL get_attribute( id_surf, char_fill, &
1181	root_area_density_lsm_f%fill, &
1182	.FALSE., 'root_area_dens_s' )
1183	!
1184	!-- Obtain number of soil layers from file and allocate variable
1185	CALL get_dimension_length( id_surf, root_area_density_lsm_f%nz,&
1186	'zsoil' )
1187	ALLOCATE( root_area_density_lsm_f%var &
1188	(0:root_area_density_lsm_f%nz-1, &
1189	nys:nyn,nxl:nxr) )
1190
1191	!
1192	!-- Read root-area density
1193	DO i = nxl, nxr
1194	DO j = nys, nyn
1195	CALL get_variable( id_surf, 'root_area_dens_s', &
1196	i, j, &
1197	root_area_density_lsm_f%var(:,j,i) )
1198	ENDDO
1199	ENDDO
1200
1201	ELSE
1202	root_area_density_lsm_f%from_file = .FALSE.
1203	ENDIF
1204	!
1205	!-- Read street type and street crossing
1206	IF ( check_existence( var_names, 'street_type' ) ) THEN
1207	street_type_f%from_file = .TRUE.
1208	CALL get_attribute( id_surf, char_fill, &
1209	street_type_f%fill, .FALSE., &
1210	'street_type' )
1211	!
1212	!-- PE-wise reading of 2D pavement type.
1213	ALLOCATE ( street_type_f%var(nys:nyn,nxl:nxr) )
1214	DO i = nxl, nxr
1215	CALL get_variable( id_surf, 'street_type', &
1216	i, street_type_f%var(:,i) )
1217	ENDDO
1218	ELSE
1219	street_type_f%from_file = .FALSE.
1220	ENDIF
1221
1222	IF ( check_existence( var_names, 'street_crossing' ) ) THEN
1223	street_crossing_f%from_file = .TRUE.
1224	CALL get_attribute( id_surf, char_fill, &
1225	street_crossing_f%fill, .FALSE., &
1226	'street_crossing' )
1227	!
1228	!-- PE-wise reading of 2D pavement type.
1229	ALLOCATE ( street_crossing_f%var(nys:nyn,nxl:nxr) )
1230	DO i = nxl, nxr
1231	CALL get_variable( id_surf, 'street_crossing', &
1232	i, street_crossing_f%var(:,i) )
1233	ENDDO
1234	ELSE
1235	street_crossing_f%from_file = .FALSE.
1236	ENDIF
1237	!
1238	!-- Still missing: root_resolved and building_surface_pars.
1239	!-- Will be implemented as soon as they are available.
1240
1241	!
1242	!-- Finally, close input file
1243	CALL close_input_file( id_surf )
1244	#endif
1245	ENDIF
1246	#if defined( __parallel )
1247	CALL MPI_BARRIER( comm2d, ierr )
1248	#endif
1249	ENDDO
1250	!
1251	!-- End of CPU measurement
1252	CALL cpu_log( log_point_s(82), 'NetCDF input', 'stop' )
1253	!
1254	!-- Exchange 1 ghost points for surface variables. Please note, ghost point
1255	!-- exchange for 3D parameter lists should be revised by using additional
1256	!-- MPI datatypes or rewriting exchange_horiz.
1257	!-- Moreover, varialbes will be resized in the following, including ghost
1258	!-- points.
1259	!-- Start with 2D Integer variables. Please note, for 8-bit integer
1260	!-- variables must be swapt to 32-bit integer before calling exchange_horiz.
1261	IF ( albedo_type_f%from_file ) THEN
1262	var_exchange_int = INT( albedo_type_f%fill, KIND = 1 )
1263	var_exchange_int(nys:nyn,nxl:nxr) = &
1264	INT( albedo_type_f%var(nys:nyn,nxl:nxr), KIND = 4 )
1265	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1266	DEALLOCATE( albedo_type_f%var )
1267	ALLOCATE( albedo_type_f%var(nysg:nyng,nxlg:nxrg) )
1268	albedo_type_f%var = INT( var_exchange_int, KIND = 1 )
1269	ENDIF
1270	IF ( pavement_type_f%from_file ) THEN
1271	var_exchange_int = INT( pavement_type_f%fill, KIND = 1 )
1272	var_exchange_int(nys:nyn,nxl:nxr) = &
1273	INT( pavement_type_f%var(nys:nyn,nxl:nxr), KIND = 4 )
1274	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1275	DEALLOCATE( pavement_type_f%var )
1276	ALLOCATE( pavement_type_f%var(nysg:nyng,nxlg:nxrg) )
1277	pavement_type_f%var = INT( var_exchange_int, KIND = 1 )
1278	ENDIF
1279	IF ( soil_type_f%from_file .AND. ALLOCATED( soil_type_f%var_2d ) ) THEN
1280	var_exchange_int = INT( soil_type_f%fill, KIND = 1 )
1281	var_exchange_int(nys:nyn,nxl:nxr) = &
1282	INT( soil_type_f%var_2d(nys:nyn,nxl:nxr), KIND = 4 )
1283	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1284	DEALLOCATE( soil_type_f%var_2d )
1285	ALLOCATE( soil_type_f%var_2d(nysg:nyng,nxlg:nxrg) )
1286	soil_type_f%var_2d = INT( var_exchange_int, KIND = 1 )
1287	ENDIF
1288	IF ( vegetation_type_f%from_file ) THEN
1289	var_exchange_int = INT( vegetation_type_f%fill, KIND = 1 )
1290	var_exchange_int(nys:nyn,nxl:nxr) = &
1291	INT( vegetation_type_f%var(nys:nyn,nxl:nxr), KIND = 4 )
1292	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1293	DEALLOCATE( vegetation_type_f%var )
1294	ALLOCATE( vegetation_type_f%var(nysg:nyng,nxlg:nxrg) )
1295	vegetation_type_f%var = INT( var_exchange_int, KIND = 1 )
1296	ENDIF
1297	IF ( water_type_f%from_file ) THEN
1298	var_exchange_int = INT( water_type_f%fill, KIND = 1 )
1299	var_exchange_int(nys:nyn,nxl:nxr) = &
1300	INT( water_type_f%var(nys:nyn,nxl:nxr), KIND = 4 )
1301	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1302	DEALLOCATE( water_type_f%var )
1303	ALLOCATE( water_type_f%var(nysg:nyng,nxlg:nxrg) )
1304	water_type_f%var = INT( var_exchange_int, KIND = 1 )
1305	ENDIF
1306	!
1307	!-- Exchange 1 ghost point for 3/4-D variables. For the sake of simplicity,
1308	!-- loop further dimensions to use 2D exchange routines.
1309	!-- This should be revised later by introducing new MPI datatypes.
1310	IF ( soil_type_f%from_file .AND. ALLOCATED( soil_type_f%var_3d ) ) &
1311	THEN
1312	ALLOCATE( var_dum_int_3d(0:nz_soil,nys:nyn,nxl:nxr) )
1313	var_dum_int_3d = soil_type_f%var_3d
1314	DEALLOCATE( soil_type_f%var_3d )
1315	ALLOCATE( soil_type_f%var_3d(0:nz_soil,nysg:nyng,nxlg:nxrg) )
1316	soil_type_f%var_3d = soil_type_f%fill
1317
1318	DO k = 0, nz_soil
1319	var_exchange_int(nys:nyn,nxl:nxr) = var_dum_int_3d(k,nys:nyn,nxl:nxr)
1320	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1321	soil_type_f%var_3d(k,:,:) = INT( var_exchange_int(:,:), KIND = 1 )
1322	ENDDO
1323	DEALLOCATE( var_dum_int_3d )
1324	ENDIF
1325
1326	IF ( surface_fraction_f%from_file ) THEN
1327	ALLOCATE( var_dum_real_3d(0:surface_fraction_f%nf-1,nys:nyn,nxl:nxr) )
1328	var_dum_real_3d = surface_fraction_f%frac
1329	DEALLOCATE( surface_fraction_f%frac )
1330	ALLOCATE( surface_fraction_f%frac(0:surface_fraction_f%nf-1, &
1331	nysg:nyng,nxlg:nxrg) )
1332	surface_fraction_f%frac = surface_fraction_f%fill
1333
1334	DO k = 0, surface_fraction_f%nf-1
1335	var_exchange_real(nys:nyn,nxl:nxr) = var_dum_real_3d(k,nys:nyn,nxl:nxr)
1336	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1337	surface_fraction_f%frac(k,:,:) = var_exchange_real(:,:)
1338	ENDDO
1339	DEALLOCATE( var_dum_real_3d )
1340	ENDIF
1341
1342	IF ( building_pars_f%from_file ) THEN
1343	ALLOCATE( var_dum_real_3d(0:building_pars_f%np-1,nys:nyn,nxl:nxr) )
1344	var_dum_real_3d = building_pars_f%pars_xy
1345	DEALLOCATE( building_pars_f%pars_xy )
1346	ALLOCATE( building_pars_f%pars_xy(0:building_pars_f%np-1, &
1347	nysg:nyng,nxlg:nxrg) )
1348	building_pars_f%pars_xy = building_pars_f%fill
1349	DO k = 0, building_pars_f%np-1
1350	var_exchange_real(nys:nyn,nxl:nxr) = &
1351	var_dum_real_3d(k,nys:nyn,nxl:nxr)
1352	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1353	building_pars_f%pars_xy(k,:,:) = var_exchange_real(:,:)
1354	ENDDO
1355	DEALLOCATE( var_dum_real_3d )
1356	ENDIF
1357
1358	IF ( albedo_pars_f%from_file ) THEN
1359	ALLOCATE( var_dum_real_3d(0:albedo_pars_f%np-1,nys:nyn,nxl:nxr) )
1360	var_dum_real_3d = albedo_pars_f%pars_xy
1361	DEALLOCATE( albedo_pars_f%pars_xy )
1362	ALLOCATE( albedo_pars_f%pars_xy(0:albedo_pars_f%np-1, &
1363	nysg:nyng,nxlg:nxrg) )
1364	albedo_pars_f%pars_xy = albedo_pars_f%fill
1365	DO k = 0, albedo_pars_f%np-1
1366	var_exchange_real(nys:nyn,nxl:nxr) = &
1367	var_dum_real_3d(k,nys:nyn,nxl:nxr)
1368	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1369	albedo_pars_f%pars_xy(k,:,:) = var_exchange_real(:,:)
1370	ENDDO
1371	DEALLOCATE( var_dum_real_3d )
1372	ENDIF
1373
1374	IF ( pavement_pars_f%from_file ) THEN
1375	ALLOCATE( var_dum_real_3d(0:pavement_pars_f%np-1,nys:nyn,nxl:nxr) )
1376	var_dum_real_3d = pavement_pars_f%pars_xy
1377	DEALLOCATE( pavement_pars_f%pars_xy )
1378	ALLOCATE( pavement_pars_f%pars_xy(0:pavement_pars_f%np-1, &
1379	nysg:nyng,nxlg:nxrg) )
1380	pavement_pars_f%pars_xy = pavement_pars_f%fill
1381	DO k = 0, pavement_pars_f%np-1
1382	var_exchange_real(nys:nyn,nxl:nxr) = &
1383	var_dum_real_3d(k,nys:nyn,nxl:nxr)
1384	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1385	pavement_pars_f%pars_xy(k,:,:) = var_exchange_real(:,:)
1386	ENDDO
1387	DEALLOCATE( var_dum_real_3d )
1388	ENDIF
1389
1390	IF ( vegetation_pars_f%from_file ) THEN
1391	ALLOCATE( var_dum_real_3d(0:vegetation_pars_f%np-1,nys:nyn,nxl:nxr) )
1392	var_dum_real_3d = vegetation_pars_f%pars_xy
1393	DEALLOCATE( vegetation_pars_f%pars_xy )
1394	ALLOCATE( vegetation_pars_f%pars_xy(0:vegetation_pars_f%np-1, &
1395	nysg:nyng,nxlg:nxrg) )
1396	vegetation_pars_f%pars_xy = vegetation_pars_f%fill
1397	DO k = 0, vegetation_pars_f%np-1
1398	var_exchange_real(nys:nyn,nxl:nxr) = &
1399	var_dum_real_3d(k,nys:nyn,nxl:nxr)
1400	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1401	vegetation_pars_f%pars_xy(k,:,:) = var_exchange_real(:,:)
1402	ENDDO
1403	DEALLOCATE( var_dum_real_3d )
1404	ENDIF
1405
1406	IF ( water_pars_f%from_file ) THEN
1407	ALLOCATE( var_dum_real_3d(0:water_pars_f%np-1,nys:nyn,nxl:nxr) )
1408	var_dum_real_3d = water_pars_f%pars_xy
1409	DEALLOCATE( water_pars_f%pars_xy )
1410	ALLOCATE( water_pars_f%pars_xy(0:water_pars_f%np-1, &
1411	nysg:nyng,nxlg:nxrg) )
1412	water_pars_f%pars_xy = water_pars_f%fill
1413	DO k = 0, water_pars_f%np-1
1414	var_exchange_real(nys:nyn,nxl:nxr) = &
1415	var_dum_real_3d(k,nys:nyn,nxl:nxr)
1416	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1417	water_pars_f%pars_xy(k,:,:) = var_exchange_real(:,:)
1418	ENDDO
1419	DEALLOCATE( var_dum_real_3d )
1420	ENDIF
1421
1422	IF ( root_area_density_lsm_f%from_file ) THEN
1423	ALLOCATE( var_dum_real_3d(0:root_area_density_lsm_f%nz-1,nys:nyn,nxl:nxr) )
1424	var_dum_real_3d = root_area_density_lsm_f%var
1425	DEALLOCATE( root_area_density_lsm_f%var )
1426	ALLOCATE( root_area_density_lsm_f%var(0:root_area_density_lsm_f%nz-1,&
1427	nysg:nyng,nxlg:nxrg) )
1428	root_area_density_lsm_f%var = root_area_density_lsm_f%fill
1429
1430	DO k = 0, root_area_density_lsm_f%nz-1
1431	var_exchange_real(nys:nyn,nxl:nxr) = &
1432	var_dum_real_3d(k,nys:nyn,nxl:nxr)
1433	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1434	root_area_density_lsm_f%var(k,:,:) = var_exchange_real(:,:)
1435	ENDDO
1436	DEALLOCATE( var_dum_real_3d )
1437	ENDIF
1438
1439	IF ( soil_pars_f%from_file ) THEN
1440	IF ( soil_pars_f%lod == 1 ) THEN
1441
1442	ALLOCATE( var_dum_real_3d(0:soil_pars_f%np-1,nys:nyn,nxl:nxr) )
1443	var_dum_real_3d = soil_pars_f%pars_xy
1444	DEALLOCATE( soil_pars_f%pars_xy )
1445	ALLOCATE( soil_pars_f%pars_xy(0:soil_pars_f%np-1, &
1446	nysg:nyng,nxlg:nxrg) )
1447	soil_pars_f%pars_xy = soil_pars_f%fill
1448
1449	DO k = 0, soil_pars_f%np-1
1450	var_exchange_real(nys:nyn,nxl:nxr) = &
1451	var_dum_real_3d(k,nys:nyn,nxl:nxr)
1452	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1453	soil_pars_f%pars_xy(k,:,:) = var_exchange_real(:,:)
1454	ENDDO
1455	DEALLOCATE( var_dum_real_3d )
1456	ELSEIF ( soil_pars_f%lod == 2 ) THEN
1457	ALLOCATE( var_dum_real_4d(0:soil_pars_f%np-1, &
1458	0:soil_pars_f%nz-1, &
1459	nys:nyn,nxl:nxr) )
1460	var_dum_real_4d = soil_pars_f%pars_xyz
1461	DEALLOCATE( soil_pars_f%pars_xyz )
1462	ALLOCATE( soil_pars_f%pars_xyz(0:soil_pars_f%np-1, &
1463	0:soil_pars_f%nz-1, &
1464	nysg:nyng,nxlg:nxrg) )
1465	soil_pars_f%pars_xyz = soil_pars_f%fill
1466
1467	DO k2 = 0, soil_pars_f%nz-1
1468	DO k = 0, soil_pars_f%np-1
1469	var_exchange_real(nys:nyn,nxl:nxr) = &
1470	var_dum_real_4d(k,k2,nys:nyn,nxl:nxr)
1471	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1472
1473	soil_pars_f%pars_xyz(k,k2,:,:) = var_exchange_real(:,:)
1474	ENDDO
1475	ENDDO
1476	DEALLOCATE( var_dum_real_4d )
1477	ENDIF
1478	ENDIF
1479
1480	IF ( pavement_subsurface_pars_f%from_file ) THEN
1481	ALLOCATE( var_dum_real_4d(0:pavement_subsurface_pars_f%np-1, &
1482	0:pavement_subsurface_pars_f%nz-1, &
1483	nys:nyn,nxl:nxr) )
1484	var_dum_real_4d = pavement_subsurface_pars_f%pars_xyz
1485	DEALLOCATE( pavement_subsurface_pars_f%pars_xyz )
1486	ALLOCATE( pavement_subsurface_pars_f%pars_xyz &
1487	(0:pavement_subsurface_pars_f%np-1, &
1488	0:pavement_subsurface_pars_f%nz-1, &
1489	nysg:nyng,nxlg:nxrg) )
1490	pavement_subsurface_pars_f%pars_xyz = pavement_subsurface_pars_f%fill
1491
1492	DO k2 = 0, pavement_subsurface_pars_f%nz-1
1493	DO k = 0, pavement_subsurface_pars_f%np-1
1494	var_exchange_real(nys:nyn,nxl:nxr) = &
1495	var_dum_real_4d(k,k2,nys:nyn,nxl:nxr)
1496	CALL exchange_horiz_2d( var_exchange_real, nbgp )
1497	pavement_subsurface_pars_f%pars_xyz(k,k2,:,:) = &
1498	var_exchange_real(:,:)
1499	ENDDO
1500	ENDDO
1501	DEALLOCATE( var_dum_real_4d )
1502	ENDIF
1503
1504	!
1505	!-- In case of non-cyclic boundary conditions, set Neumann conditions at the
1506	!-- lateral boundaries.
1507	IF ( .NOT. bc_ns_cyc ) THEN
1508	IF ( nys == 0 ) THEN
1509	IF ( albedo_type_f%from_file ) &
1510	albedo_type_f%var(-1,:) = albedo_type_f%var(0,:)
1511	IF ( pavement_type_f%from_file ) &
1512	pavement_type_f%var(-1,:) = pavement_type_f%var(0,:)
1513	IF ( soil_type_f%from_file ) THEN
1514	IF ( ALLOCATED( soil_type_f%var_2d ) ) THEN
1515	soil_type_f%var_2d(-1,:) = soil_type_f%var_2d(0,:)
1516	ELSE
1517	soil_type_f%var_3d(:,-1,:) = soil_type_f%var_3d(:,0,:)
1518	ENDIF
1519	ENDIF
1520	IF ( vegetation_type_f%from_file ) &
1521	vegetation_type_f%var(-1,:) = vegetation_type_f%var(0,:)
1522	IF ( water_type_f%from_file ) &
1523	water_type_f%var(-1,:) = water_type_f%var(0,:)
1524	IF ( surface_fraction_f%from_file ) &
1525	surface_fraction_f%frac(:,-1,:) = surface_fraction_f%frac(:,0,:)
1526	IF ( building_pars_f%from_file ) &
1527	building_pars_f%pars_xy(:,-1,:) = building_pars_f%pars_xy(:,0,:)
1528	IF ( albedo_pars_f%from_file ) &
1529	albedo_pars_f%pars_xy(:,-1,:) = albedo_pars_f%pars_xy(:,0,:)
1530	IF ( pavement_pars_f%from_file ) &
1531	pavement_pars_f%pars_xy(:,-1,:) = pavement_pars_f%pars_xy(:,0,:)
1532	IF ( vegetation_pars_f%from_file ) &
1533	vegetation_pars_f%pars_xy(:,-1,:) = &
1534	vegetation_pars_f%pars_xy(:,0,:)
1535	IF ( water_pars_f%from_file ) &
1536	water_pars_f%pars_xy(:,-1,:) = water_pars_f%pars_xy(:,0,:)
1537	IF ( root_area_density_lsm_f%from_file ) &
1538	root_area_density_lsm_f%var(:,-1,:) = &
1539	root_area_density_lsm_f%var(:,0,:)
1540	IF ( soil_pars_f%from_file ) THEN
1541	IF ( soil_pars_f%lod == 1 ) THEN
1542	soil_pars_f%pars_xy(:,-1,:) = soil_pars_f%pars_xy(:,0,:)
1543	ELSE
1544	soil_pars_f%pars_xyz(:,:,-1,:) = soil_pars_f%pars_xyz(:,:,0,:)
1545	ENDIF
1546	ENDIF
1547	IF ( pavement_subsurface_pars_f%from_file ) &
1548	pavement_subsurface_pars_f%pars_xyz(:,:,-1,:) = &
1549	pavement_subsurface_pars_f%pars_xyz(:,:,0,:)
1550	ENDIF
1551
1552	IF ( nyn == ny ) THEN
1553	IF ( albedo_type_f%from_file ) &
1554	albedo_type_f%var(ny+1,:) = albedo_type_f%var(ny,:)
1555	IF ( pavement_type_f%from_file ) &
1556	pavement_type_f%var(ny+1,:) = pavement_type_f%var(ny,:)
1557	IF ( soil_type_f%from_file ) THEN
1558	IF ( ALLOCATED( soil_type_f%var_2d ) ) THEN
1559	soil_type_f%var_2d(ny+1,:) = soil_type_f%var_2d(ny,:)
1560	ELSE
1561	soil_type_f%var_3d(:,ny+1,:) = soil_type_f%var_3d(:,ny,:)
1562	ENDIF
1563	ENDIF
1564	IF ( vegetation_type_f%from_file ) &
1565	vegetation_type_f%var(ny+1,:) = vegetation_type_f%var(ny,:)
1566	IF ( water_type_f%from_file ) &
1567	water_type_f%var(ny+1,:) = water_type_f%var(ny,:)
1568	IF ( surface_fraction_f%from_file ) &
1569	surface_fraction_f%frac(:,ny+1,:) = &
1570	surface_fraction_f%frac(:,ny,:)
1571	IF ( building_pars_f%from_file ) &
1572	building_pars_f%pars_xy(:,ny+1,:) = &
1573	building_pars_f%pars_xy(:,ny,:)
1574	IF ( albedo_pars_f%from_file ) &
1575	albedo_pars_f%pars_xy(:,ny+1,:) = albedo_pars_f%pars_xy(:,ny,:)
1576	IF ( pavement_pars_f%from_file ) &
1577	pavement_pars_f%pars_xy(:,ny+1,:) = &
1578	pavement_pars_f%pars_xy(:,ny,:)
1579	IF ( vegetation_pars_f%from_file ) &
1580	vegetation_pars_f%pars_xy(:,ny+1,:) = &
1581	vegetation_pars_f%pars_xy(:,ny,:)
1582	IF ( water_pars_f%from_file ) &
1583	water_pars_f%pars_xy(:,ny+1,:) = water_pars_f%pars_xy(:,ny,:)
1584	IF ( root_area_density_lsm_f%from_file ) &
1585	root_area_density_lsm_f%var(:,ny+1,:) = &
1586	root_area_density_lsm_f%var(:,ny,:)
1587	IF ( soil_pars_f%from_file ) THEN
1588	IF ( soil_pars_f%lod == 1 ) THEN
1589	soil_pars_f%pars_xy(:,ny+1,:) = soil_pars_f%pars_xy(:,ny,:)
1590	ELSE
1591	soil_pars_f%pars_xyz(:,:,ny+1,:) = &
1592	soil_pars_f%pars_xyz(:,:,ny,:)
1593	ENDIF
1594	ENDIF
1595	IF ( pavement_subsurface_pars_f%from_file ) &
1596	pavement_subsurface_pars_f%pars_xyz(:,:,ny+1,:) = &
1597	pavement_subsurface_pars_f%pars_xyz(:,:,ny,:)
1598	ENDIF
1599	ENDIF
1600
1601	IF ( .NOT. bc_lr_cyc ) THEN
1602	IF ( nxl == 0 ) THEN
1603	IF ( albedo_type_f%from_file ) &
1604	albedo_type_f%var(:,-1) = albedo_type_f%var(:,0)
1605	IF ( pavement_type_f%from_file ) &
1606	pavement_type_f%var(:,-1) = pavement_type_f%var(:,0)
1607	IF ( soil_type_f%from_file ) THEN
1608	IF ( ALLOCATED( soil_type_f%var_2d ) ) THEN
1609	soil_type_f%var_2d(:,-1) = soil_type_f%var_2d(:,0)
1610	ELSE
1611	soil_type_f%var_3d(:,:,-1) = soil_type_f%var_3d(:,:,0)
1612	ENDIF
1613	ENDIF
1614	IF ( vegetation_type_f%from_file ) &
1615	vegetation_type_f%var(:,-1) = vegetation_type_f%var(:,0)
1616	IF ( water_type_f%from_file ) &
1617	water_type_f%var(:,-1) = water_type_f%var(:,0)
1618	IF ( surface_fraction_f%from_file ) &
1619	surface_fraction_f%frac(:,:,-1) = surface_fraction_f%frac(:,:,0)
1620	IF ( building_pars_f%from_file ) &
1621	building_pars_f%pars_xy(:,:,-1) = building_pars_f%pars_xy(:,:,0)
1622	IF ( albedo_pars_f%from_file ) &
1623	albedo_pars_f%pars_xy(:,:,-1) = albedo_pars_f%pars_xy(:,:,0)
1624	IF ( pavement_pars_f%from_file ) &
1625	pavement_pars_f%pars_xy(:,:,-1) = pavement_pars_f%pars_xy(:,:,0)
1626	IF ( vegetation_pars_f%from_file ) &
1627	vegetation_pars_f%pars_xy(:,:,-1) = &
1628	vegetation_pars_f%pars_xy(:,:,0)
1629	IF ( water_pars_f%from_file ) &
1630	water_pars_f%pars_xy(:,:,-1) = water_pars_f%pars_xy(:,:,0)
1631	IF ( root_area_density_lsm_f%from_file ) &
1632	root_area_density_lsm_f%var(:,:,-1) = &
1633	root_area_density_lsm_f%var(:,:,0)
1634	IF ( soil_pars_f%from_file ) THEN
1635	IF ( soil_pars_f%lod == 1 ) THEN
1636	soil_pars_f%pars_xy(:,:,-1) = soil_pars_f%pars_xy(:,:,0)
1637	ELSE
1638	soil_pars_f%pars_xyz(:,:,:,-1) = soil_pars_f%pars_xyz(:,:,:,0)
1639	ENDIF
1640	ENDIF
1641	IF ( pavement_subsurface_pars_f%from_file ) &
1642	pavement_subsurface_pars_f%pars_xyz(:,:,:,-1) = &
1643	pavement_subsurface_pars_f%pars_xyz(:,:,:,0)
1644	ENDIF
1645
1646	IF ( nxr == nx ) THEN
1647	IF ( albedo_type_f%from_file ) &
1648	albedo_type_f%var(:,nx+1) = albedo_type_f%var(:,nx)
1649	IF ( pavement_type_f%from_file ) &
1650	pavement_type_f%var(:,nx+1) = pavement_type_f%var(:,nx)
1651	IF ( soil_type_f%from_file ) THEN
1652	IF ( ALLOCATED( soil_type_f%var_2d ) ) THEN
1653	soil_type_f%var_2d(:,nx+1) = soil_type_f%var_2d(:,nx)
1654	ELSE
1655	soil_type_f%var_3d(:,:,nx+1) = soil_type_f%var_3d(:,:,nx)
1656	ENDIF
1657	ENDIF
1658	IF ( vegetation_type_f%from_file ) &
1659	vegetation_type_f%var(:,nx+1) = vegetation_type_f%var(:,nx)
1660	IF ( water_type_f%from_file ) &
1661	water_type_f%var(:,nx+1) = water_type_f%var(:,nx)
1662	IF ( surface_fraction_f%from_file ) &
1663	surface_fraction_f%frac(:,:,nx+1) = &
1664	surface_fraction_f%frac(:,:,nx)
1665	IF ( building_pars_f%from_file ) &
1666	building_pars_f%pars_xy(:,:,nx+1) = &
1667	building_pars_f%pars_xy(:,:,nx)
1668	IF ( albedo_pars_f%from_file ) &
1669	albedo_pars_f%pars_xy(:,:,nx+1) = albedo_pars_f%pars_xy(:,:,nx)
1670	IF ( pavement_pars_f%from_file ) &
1671	pavement_pars_f%pars_xy(:,:,nx+1) = &
1672	pavement_pars_f%pars_xy(:,:,nx)
1673	IF ( vegetation_pars_f%from_file ) &
1674	vegetation_pars_f%pars_xy(:,:,nx+1) = &
1675	vegetation_pars_f%pars_xy(:,:,nx)
1676	IF ( water_pars_f%from_file ) &
1677	water_pars_f%pars_xy(:,:,nx+1) = water_pars_f%pars_xy(:,:,nx)
1678	IF ( root_area_density_lsm_f%from_file ) &
1679	root_area_density_lsm_f%var(:,:,nx+1) = &
1680	root_area_density_lsm_f%var(:,:,nx)
1681	IF ( soil_pars_f%from_file ) THEN
1682	IF ( soil_pars_f%lod == 1 ) THEN
1683	soil_pars_f%pars_xy(:,:,nx+1) = soil_pars_f%pars_xy(:,:,nx)
1684	ELSE
1685	soil_pars_f%pars_xyz(:,:,:,nx+1) = &
1686	soil_pars_f%pars_xyz(:,:,:,nx)
1687	ENDIF
1688	ENDIF
1689	IF ( pavement_subsurface_pars_f%from_file ) &
1690	pavement_subsurface_pars_f%pars_xyz(:,:,:,nx+1) = &
1691	pavement_subsurface_pars_f%pars_xyz(:,:,:,nx)
1692	ENDIF
1693	ENDIF
1694
1695	END SUBROUTINE netcdf_data_input_surface_data
1696
1697	!------------------------------------------------------------------------------!
1698	! Description:
1699	! ------------
1700	!> Reads orography and building information.
1701	!------------------------------------------------------------------------------!
1702	SUBROUTINE netcdf_data_input_topo
1703
1704	USE control_parameters, &
1705	ONLY: bc_lr_cyc, bc_ns_cyc, message_string, topography
1706
1707	USE indices, &
1708	ONLY: nbgp, nx, nxl, nxr, ny, nyn, nys, nzb, nzt
1709
1710
1711	IMPLICIT NONE
1712
1713	CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE :: var_names !< variable names in static input file
1714
1715
1716	INTEGER(iwp) :: i !< running index along x-direction
1717	INTEGER(iwp) :: ii !< running index for IO blocks
1718	INTEGER(iwp) :: id_topo !< NetCDF id of topograhy input file
1719	INTEGER(iwp) :: j !< running index along y-direction
1720	INTEGER(iwp) :: k !< running index along z-direction
1721	INTEGER(iwp) :: num_vars !< number of variables in netcdf input file
1722	INTEGER(iwp) :: skip_n_rows !< counting variable to skip rows while reading topography file
1723
1724	INTEGER(iwp), DIMENSION(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) :: var_exchange_int !< dummy variables used to exchange 32-bit Integer arrays
1725
1726	REAL(wp) :: dum !< dummy variable to skip columns while reading topography file
1727	!
1728	!-- CPU measurement
1729	CALL cpu_log( log_point_s(83), 'NetCDF/ASCII input topo', 'start' )
1730
1731	DO ii = 0, io_blocks-1
1732	IF ( ii == io_group ) THEN
1733	!
1734	!-- Input via palm-input data standard
1735	IF ( input_pids_static ) THEN
1736	#if defined ( __netcdf )
1737	!
1738	!-- Open file in read-only mode
1739	CALL open_read_file( TRIM( input_file_static ) // &
1740	TRIM( coupling_char ), id_topo )
1741
1742	!
1743	!-- At first, inquire all variable names.
1744	!-- This will be used to check whether an input variable exist
1745	!-- or not.
1746	CALL inquire_num_variables( id_topo, num_vars )
1747	!
1748	!-- Allocate memory to store variable names and inquire them.
1749	ALLOCATE( var_names(1:num_vars) )
1750	CALL inquire_variable_names( id_topo, var_names )
1751	!
1752	!-- Read x, y - dimensions. Only required for consistency checks.
1753	CALL get_dimension_length( id_topo, dim_static%nx, 'x' )
1754	CALL get_dimension_length( id_topo, dim_static%ny, 'y' )
1755	ALLOCATE( dim_static%x(0:dim_static%nx-1) )
1756	ALLOCATE( dim_static%y(0:dim_static%ny-1) )
1757	CALL get_variable( id_topo, 'x', dim_static%x )
1758	CALL get_variable( id_topo, 'y', dim_static%y )
1759	!
1760	!-- Terrain height. First, get variable-related _FillValue attribute
1761	IF ( check_existence( var_names, 'zt' ) ) THEN
1762	terrain_height_f%from_file = .TRUE.
1763	CALL get_attribute( id_topo, char_fill, &
1764	terrain_height_f%fill, &
1765	.FALSE., 'zt' )
1766	!
1767	!-- PE-wise reading of 2D terrain height.
1768	ALLOCATE ( terrain_height_f%var(nys:nyn,nxl:nxr) )
1769	DO i = nxl, nxr
1770	CALL get_variable( id_topo, 'zt', &
1771	i, terrain_height_f%var(:,i) )
1772	ENDDO
1773	ELSE
1774	terrain_height_f%from_file = .FALSE.
1775	ENDIF
1776
1777	!
1778	!-- Read building height. First, read its _FillValue attribute,
1779	!-- as well as lod attribute
1780	buildings_f%from_file = .FALSE.
1781	IF ( check_existence( var_names, 'buildings_2d' ) ) THEN
1782	buildings_f%from_file = .TRUE.
1783	CALL get_attribute( id_topo, char_lod, buildings_f%lod, &
1784	.FALSE., 'buildings_2d' )
1785
1786	CALL get_attribute( id_topo, char_fill, &
1787	buildings_f%fill1, &
1788	.FALSE., 'buildings_2d' )
1789
1790	!
1791	!-- Read 2D topography
1792	IF ( buildings_f%lod == 1 ) THEN
1793	ALLOCATE ( buildings_f%var_2d(nys:nyn,nxl:nxr) )
1794	DO i = nxl, nxr
1795	CALL get_variable( id_topo, 'buildings_2d', &
1796	i, buildings_f%var_2d(:,i) )
1797	ENDDO
1798	ELSE
1799	message_string = 'NetCDF attribute lod ' // &
1800	'(level of detail) is not set ' // &
1801	'properly for buildings_2d.'
1802	CALL message( 'netcdf_data_input_mod', 'NDI000', &
1803	1, 2, 0, 6, 0 )
1804	ENDIF
1805	ENDIF
1806	!
1807	!-- If available, also read 3D building information. If both are
1808	!-- available, use 3D information.
1809	IF ( check_existence( var_names, 'buildings_3d' ) ) THEN
1810	buildings_f%from_file = .TRUE.
1811	CALL get_attribute( id_topo, char_lod, buildings_f%lod, &
1812	.FALSE., 'buildings_3d' )
1813
1814	CALL get_attribute( id_topo, char_fill, &
1815	buildings_f%fill2, &
1816	.FALSE., 'buildings_3d' )
1817
1818	CALL get_dimension_length( id_topo, buildings_f%nz, 'z' )
1819
1820	IF ( buildings_f%lod == 2 ) THEN
1821	ALLOCATE( buildings_f%z(nzb:buildings_f%nz-1) )
1822	CALL get_variable( id_topo, 'z', buildings_f%z )
1823
1824	ALLOCATE( buildings_f%var_3d(nzb:buildings_f%nz-1, &
1825	nys:nyn,nxl:nxr) )
1826	buildings_f%var_3d = 0
1827	!
1828	!-- Read data PE-wise. Read yz-slices.
1829	DO i = nxl, nxr
1830	DO j = nys, nyn
1831	CALL get_variable( id_topo, 'buildings_3d', &
1832	i, j, &
1833	buildings_f%var_3d(:,j,i) )
1834	ENDDO
1835	ENDDO
1836	ELSE
1837	message_string = 'NetCDF attribute lod ' // &
1838	'(level of detail) is not set ' // &
1839	'properly for buildings_3d.'
1840	CALL message( 'netcdf_data_input_mod', 'NDI001', &
1841	1, 2, 0, 6, 0 )
1842	ENDIF
1843	ENDIF
1844	!
1845	!-- Read building IDs and its FillValue attribute. Further required
1846	!-- for mapping buildings on top of orography.
1847	IF ( check_existence( var_names, 'building_id' ) ) THEN
1848	building_id_f%from_file = .TRUE.
1849	CALL get_attribute( id_topo, char_fill, &
1850	building_id_f%fill, .FALSE., &
1851	'building_id' )
1852
1853
1854	ALLOCATE ( building_id_f%var(nys:nyn,nxl:nxr) )
1855	DO i = nxl, nxr
1856	CALL get_variable( id_topo, 'building_id', &
1857	i, building_id_f%var(:,i) )
1858	ENDDO
1859	ELSE
1860	building_id_f%from_file = .FALSE.
1861	ENDIF
1862	!
1863	!-- Read building_type and required attributes.
1864	IF ( check_existence( var_names, 'building_type' ) ) THEN
1865	building_type_f%from_file = .TRUE.
1866	CALL get_attribute( id_topo, char_fill, &
1867	building_type_f%fill, .FALSE., &
1868	'building_type' )
1869
1870	ALLOCATE ( building_type_f%var(nys:nyn,nxl:nxr) )
1871	DO i = nxl, nxr
1872	CALL get_variable( id_topo, 'building_type', &
1873	i, building_type_f%var(:,i) )
1874	ENDDO
1875	ELSE
1876	building_type_f%from_file = .FALSE.
1877	ENDIF
1878
1879	!
1880	!-- Close topography input file
1881	CALL close_input_file( id_topo )
1882	#else
1883	CONTINUE
1884	#endif
1885	!
1886	!-- ASCII input
1887	ELSEIF ( TRIM( topography ) == 'read_from_file' ) THEN
1888
1889	OPEN( 90, FILE='TOPOGRAPHY_DATA'//TRIM( coupling_char ), &
1890	STATUS='OLD', FORM='FORMATTED', ERR=10 )
1891	!
1892	!-- Read topography PE-wise. Rows are read from nyn to nys, columns
1893	!-- are read from nxl to nxr. At first, ny-nyn rows need to be skipped.
1894	skip_n_rows = 0
1895	DO WHILE ( skip_n_rows < ny - nyn )
1896	READ( 90, * )
1897	skip_n_rows = skip_n_rows + 1
1898	ENDDO
1899	!
1900	!-- Read data from nyn to nys and nxl to nxr. Therefore, skip
1901	!-- column until nxl-1 is reached
1902	ALLOCATE ( buildings_f%var_2d(nys:nyn,nxl:nxr) )
1903	DO j = nyn, nys, -1
1904	READ( 90, *, ERR=11, END=11 ) &
1905	( dum, i = 0, nxl-1 ), &
1906	( buildings_f%var_2d(j,i), i = nxl, nxr )
1907	ENDDO
1908
1909	GOTO 12
1910
1911	10 message_string = 'file TOPOGRAPHY'//TRIM( coupling_char )// &
1912	' does not exist'
1913	CALL message( 'netcdf_data_input_mod', 'PA0208', 1, 2, 0, 6, 0 )
1914
1915	11 message_string = 'errors in file TOPOGRAPHY_DATA'// &
1916	TRIM( coupling_char )
1917	CALL message( 'netcdf_data_input_mod', 'PA0209', 1, 2, 0, 6, 0 )
1918
1919	12 CLOSE( 90 )
1920	buildings_f%from_file = .TRUE.
1921
1922	ENDIF
1923
1924	ENDIF
1925	#if defined( __parallel )
1926	CALL MPI_BARRIER( comm2d, ierr )
1927	#endif
1928	ENDDO
1929	!
1930	!-- End of CPU measurement
1931	CALL cpu_log( log_point_s(83), 'NetCDF/ASCII input topo', 'stop' )
1932	!
1933	!-- Check for minimum requirement to setup building topography. If buildings
1934	!-- are provided, also an ID and a type are required.
1935	!-- Note, doing this check in check_parameters
1936	!-- will be too late (data will be used for grid inititialization before).
1937	IF ( input_pids_static ) THEN
1938	IF ( buildings_f%from_file .AND. &
1939	.NOT. building_id_f%from_file ) THEN
1940	message_string = 'If building heigths are prescribed in ' // &
1941	'static input file, also an ID is required.'
1942	CALL message( 'netcdf_data_input_mod', 'NDI002', 1, 2, 0, 6, 0 )
1943	ENDIF
1944	ENDIF
1945	!
1946	!-- In case no terrain height is provided by static input file, allocate
1947	!-- array nevertheless and set terrain height to 0, which simplifies
1948	!-- topography initialization.
1949	IF ( .NOT. terrain_height_f%from_file ) THEN
1950	ALLOCATE ( terrain_height_f%var(nys:nyn,nxl:nxr) )
1951	terrain_height_f%var = 0.0_wp
1952	ENDIF
1953	!
1954	!-- Finally, exchange 1 ghost point for building ID and type.
1955	!-- In case of non-cyclic boundary conditions set Neumann conditions at the
1956	!-- lateral boundaries.
1957	IF ( building_id_f%from_file ) THEN
1958	var_exchange_int = building_id_f%fill
1959	var_exchange_int(nys:nyn,nxl:nxr) = building_id_f%var(nys:nyn,nxl:nxr)
1960	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1961	DEALLOCATE( building_id_f%var )
1962	ALLOCATE( building_id_f%var(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
1963	building_id_f%var = var_exchange_int
1964
1965	IF ( .NOT. bc_ns_cyc ) THEN
1966	IF ( nys == 0 ) building_id_f%var(-1,:) = building_id_f%var(0,:)
1967	IF ( nyn == ny ) building_id_f%var(ny+1,:) = building_id_f%var(ny,:)
1968	ENDIF
1969	IF ( .NOT. bc_lr_cyc ) THEN
1970	IF ( nxl == 0 ) building_id_f%var(:,-1) = building_id_f%var(:,0)
1971	IF ( nxr == nx ) building_id_f%var(:,nx+1) = building_id_f%var(:,nx)
1972	ENDIF
1973	ENDIF
1974
1975	IF ( building_type_f%from_file ) THEN
1976	var_exchange_int = INT( building_type_f%fill, KIND = 4 )
1977	var_exchange_int(nys:nyn,nxl:nxr) = &
1978	INT( building_type_f%var(nys:nyn,nxl:nxr), KIND = 4 )
1979	CALL exchange_horiz_2d_int( var_exchange_int, nys, nyn, nxl, nxr, nbgp )
1980	DEALLOCATE( building_type_f%var )
1981	ALLOCATE( building_type_f%var(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
1982	building_type_f%var = INT( var_exchange_int, KIND = 1 )
1983
1984	IF ( .NOT. bc_ns_cyc ) THEN
1985	IF ( nys == 0 ) building_type_f%var(-1,:) = building_type_f%var(0,:)
1986	IF ( nyn == ny ) building_type_f%var(ny+1,:) = building_type_f%var(ny,:)
1987	ENDIF
1988	IF ( .NOT. bc_lr_cyc ) THEN
1989	IF ( nxl == 0 ) building_type_f%var(:,-1) = building_type_f%var(:,0)
1990	IF ( nxr == nx ) building_type_f%var(:,nx+1) = building_type_f%var(:,nx)
1991	ENDIF
1992	ENDIF
1993
1994	END SUBROUTINE netcdf_data_input_topo
1995
1996	!------------------------------------------------------------------------------!
1997	! Description:
1998	! ------------
1999	!> Reads initialization data of u, v, w, pt, q, geostrophic wind components,
2000	!> as well as soil moisture and soil temperature, derived from larger-scale
2001	!> model (COSMO) by Inifor.
2002	!------------------------------------------------------------------------------!
2003	SUBROUTINE netcdf_data_input_init_3d
2004
2005	USE arrays_3d, &
2006	ONLY: q, pt, u, v, w
2007
2008	USE control_parameters, &
2009	ONLY: bc_lr_cyc, bc_ns_cyc, forcing, humidity, land_surface, &
2010	message_string, neutral, surface_pressure
2011
2012	USE indices, &
2013	ONLY: nx, nxl, nxlu, nxr, ny, nyn, nys, nysv, nzb, nz, nzt
2014
2015	IMPLICIT NONE
2016
2017	CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE :: var_names
2018
2019	INTEGER(iwp) :: i !< running index along x-direction
2020	INTEGER(iwp) :: ii !< running index for IO blocks
2021	INTEGER(iwp) :: id_dynamic !< NetCDF id of dynamic input file
2022	INTEGER(iwp) :: j !< running index along y-direction
2023	INTEGER(iwp) :: k !< running index along z-direction
2024	INTEGER(iwp) :: num_vars !< number of variables in netcdf input file
2025	INTEGER(iwp) :: off_i !< offset in x-direction used for reading the u-component
2026	INTEGER(iwp) :: off_j !< offset in y-direction used for reading the v-component
2027
2028	LOGICAL :: check_passed !< flag indicating if a check passed
2029
2030	!
2031	!-- Skip routine if no input file with dynamic input data is available.
2032	IF ( .NOT. input_pids_dynamic ) RETURN
2033	!
2034	!-- Please note, Inifor is designed to provide initial data for u and v for
2035	!-- the prognostic grid points in case of lateral Dirichlet conditions.
2036	!-- This means that Inifor provides data from nxlu:nxr (for u) and
2037	!-- from nysv:nyn (for v) at the left and south domain boundary, respectively.
2038	!-- However, as work-around for the moment, PALM will run with cyclic
2039	!-- conditions and will be initialized with data provided by Inifor
2040	!-- boundaries in case of Dirichlet.
2041	!-- Hence, simply set set nxlu/nysv to 1 (will be reset to its original value
2042	!-- at the end of this routine.
2043	IF ( bc_lr_cyc .AND. nxl == 0 ) nxlu = 1
2044	IF ( bc_ns_cyc .AND. nys == 0 ) nysv = 1
2045
2046	!
2047	!-- CPU measurement
2048	CALL cpu_log( log_point_s(85), 'NetCDF input init', 'start' )
2049
2050	DO ii = 0, io_blocks-1
2051	IF ( ii == io_group ) THEN
2052	#if defined ( __netcdf )
2053	!
2054	!-- Open file in read-only mode
2055	CALL open_read_file( TRIM( input_file_dynamic ) // &
2056	TRIM( coupling_char ), id_dynamic )
2057
2058	!
2059	!-- At first, inquire all variable names.
2060	CALL inquire_num_variables( id_dynamic, num_vars )
2061	!
2062	!-- Allocate memory to store variable names.
2063	ALLOCATE( var_names(1:num_vars) )
2064	CALL inquire_variable_names( id_dynamic, var_names )
2065	!
2066	!-- Read vertical dimension of scalar und w grid. Will be used for
2067	!-- inter- and extrapolation in case of stretched numeric grid.
2068	!-- This will be removed when Inifor is able to handle stretched grids.
2069	CALL get_dimension_length( id_dynamic, init_3d%nzu, 'z' )
2070	CALL get_dimension_length( id_dynamic, init_3d%nzw, 'zw' )
2071	CALL get_dimension_length( id_dynamic, init_3d%nzs, 'depth' )
2072	!
2073	!-- Read also the horizontal dimensions. These are used just used fo
2074	!-- checking the compatibility with the PALM grid before reading.
2075	CALL get_dimension_length( id_dynamic, init_3d%nx, 'x' )
2076	CALL get_dimension_length( id_dynamic, init_3d%nxu, 'xu' )
2077	CALL get_dimension_length( id_dynamic, init_3d%ny, 'y' )
2078	CALL get_dimension_length( id_dynamic, init_3d%nyv, 'yv' )
2079
2080	!
2081	!-- Check for correct horizontal and vertical dimension. Please note,
2082	!-- checks are performed directly here and not called from
2083	!-- check_parameters as some varialbes are still not allocated there.
2084	!-- Moreover, please note, u- and v-grid has 1 grid point less on
2085	!-- Inifor grid.
2086	IF ( init_3d%nx-1 /= nx .OR. init_3d%nxu-1 /= nx - 1 .OR. &
2087	init_3d%ny-1 /= ny .OR. init_3d%nyv-1 /= ny - 1 ) THEN
2088	message_string = 'Number of inifor horizontal grid points '// &
2089	'does not match the number of numeric grid '//&
2090	'points.'
2091	CALL message( 'netcdf_data_input_mod', 'NDI003', 1, 2, 0, 6, 0 )
2092	ENDIF
2093
2094	IF ( init_3d%nzu-1 /= nz ) THEN
2095	message_string = 'Number of inifor vertical grid points ' // &
2096	'does not match the number of numeric grid '//&
2097	'points.'
2098	CALL message( 'netcdf_data_input_mod', 'NDI003', 1, 2, 0, 6, 0 )
2099	ENDIF
2100	!
2101	!-- Read vertical dimensions. Later, these are required for eventual
2102	!-- inter- and extrapolations of the initialization data.
2103	IF ( check_existence( var_names, 'z' ) ) THEN
2104	ALLOCATE( init_3d%zu_atmos(1:init_3d%nzu) )
2105	CALL get_variable( id_dynamic, 'z', init_3d%zu_atmos )
2106	ENDIF
2107	IF ( check_existence( var_names, 'zw' ) ) THEN
2108	ALLOCATE( init_3d%zw_atmos(1:init_3d%nzw) )
2109	CALL get_variable( id_dynamic, 'zw', init_3d%zw_atmos )
2110	ENDIF
2111	IF ( check_existence( var_names, 'depth' ) ) THEN
2112	ALLOCATE( init_3d%z_soil(1:init_3d%nzs) )
2113	CALL get_variable( id_dynamic, 'depth', init_3d%z_soil )
2114	ENDIF
2115	!
2116	!-- Read initial geostrophic wind components at
2117	!-- t = 0 (index 1 in file).
2118	IF ( check_existence( var_names, 'tend_ug' ) ) THEN
2119	ALLOCATE( init_3d%ug_init(nzb:nzt+1) )
2120	CALL get_variable_pr( id_dynamic, 'tend_ug', 1, &
2121	init_3d%ug_init )
2122	init_3d%from_file_ug = .TRUE.
2123	ELSE
2124	init_3d%from_file_ug = .FALSE.
2125	ENDIF
2126	IF ( check_existence( var_names, 'tend_vg' ) ) THEN
2127	ALLOCATE( init_3d%vg_init(nzb:nzt+1) )
2128	CALL get_variable_pr( id_dynamic, 'tend_vg', 1, &
2129	init_3d%vg_init )
2130	init_3d%from_file_vg = .TRUE.
2131	ELSE
2132	init_3d%from_file_vg = .FALSE.
2133	ENDIF
2134	!
2135	!-- Read inital 3D data of u, v, w, pt and q,
2136	!-- derived from COSMO model. Read PE-wise yz-slices.
2137	!-- Please note, the u-, v- and w-component are defined on different
2138	!-- grids with one element less in the x-, y-,
2139	!-- and z-direction, respectively. Hence, reading is subdivided
2140	!-- into separate loops. Moreover, i and j are used
2141	!-- as start index in the NF90 interface.
2142	!-- The passed arguments for u, and v are (i,j)-1, respectively,
2143	!-- in contrast to the remaining quantities. This is because in case
2144	!-- of forcing is applied, the input data for u and v has one
2145	!-- element less along the x- and y-direction respectively.
2146	!-- Read u-component
2147	IF ( check_existence( var_names, 'init_u' ) ) THEN
2148	!
2149	!-- Read attributes for the fill value and level-of-detail
2150	CALL get_attribute( id_dynamic, char_fill, init_3d%fill_u, &
2151	.FALSE., 'init_u' )
2152	CALL get_attribute( id_dynamic, char_lod, init_3d%lod_u, &
2153	.FALSE., 'init_u' )
2154	!
2155	!-- level-of-detail 1 - read initialization profile
2156	IF ( init_3d%lod_u == 1 ) THEN
2157	ALLOCATE( init_3d%u_init(nzb:nzt+1) )
2158	init_3d%u_init = 0.0_wp
2159
2160	CALL get_variable( id_dynamic, 'init_u', &
2161	init_3d%u_init(nzb+1:nzt+1) )
2162	!
2163	!-- level-of-detail 2 - read 3D initialization data
2164	ELSEIF ( init_3d%lod_u == 2 ) THEN
2165	!
2166	!-- Set offset value. In case of Dirichlet conditions at the left
2167	!-- domain boundary, the u component starts at nxl+1. This case,
2168	!-- the passed start-index for reading the NetCDF data is shifted
2169	!-- by -1.
2170	off_i = 1 !MERGE( 1, 0, forcing )
2171
2172	DO i = nxlu, nxr
2173	DO j = nys, nyn
2174	CALL get_variable( id_dynamic, 'init_u', i-off_i, j, &
2175	u(nzb+1:nzt+1,j,i) )
2176	ENDDO
2177	ENDDO
2178
2179	ENDIF
2180	init_3d%from_file_u = .TRUE.
2181	ENDIF
2182	!
2183	!-- Read v-component
2184	IF ( check_existence( var_names, 'init_v' ) ) THEN
2185	!
2186	!-- Read attributes for the fill value and level-of-detail
2187	CALL get_attribute( id_dynamic, char_fill, init_3d%fill_v, &
2188	.FALSE., 'init_v' )
2189	CALL get_attribute( id_dynamic, char_lod, init_3d%lod_v, &
2190	.FALSE., 'init_v' )
2191	!
2192	!-- level-of-detail 1 - read initialization profile
2193	IF ( init_3d%lod_v == 1 ) THEN
2194	ALLOCATE( init_3d%v_init(nzb:nzt+1) )
2195	init_3d%v_init = 0.0_wp
2196
2197	CALL get_variable( id_dynamic, 'init_v', &
2198	init_3d%v_init(nzb+1:nzt+1) )
2199
2200	!
2201	!-- level-of-detail 2 - read 3D initialization data
2202	ELSEIF ( init_3d%lod_v == 2 ) THEN
2203	!
2204	!-- Set offset value. In case of Dirichlet conditions at the south
2205	!-- domain boundary, the v component starts at nys+1. This case,
2206	!-- the passed start-index for reading the NetCDF data is shifted
2207	!-- by -1.
2208	off_j = 1 !MERGE( 1, 0, forcing )
2209
2210	DO i = nxl, nxr
2211	DO j = nysv, nyn
2212	CALL get_variable( id_dynamic, 'init_v', i, j-off_j, &
2213	v(nzb+1:nzt+1,j,i) )
2214	ENDDO
2215	ENDDO
2216
2217	ENDIF
2218	init_3d%from_file_v = .TRUE.
2219	ENDIF
2220	!
2221	!-- Read w-component
2222	IF ( check_existence( var_names, 'init_w' ) ) THEN
2223	!
2224	!-- Read attributes for the fill value and level-of-detail
2225	CALL get_attribute( id_dynamic, char_fill, init_3d%fill_w, &
2226	.FALSE., 'init_w' )
2227	CALL get_attribute( id_dynamic, char_lod, init_3d%lod_w, &
2228	.FALSE., 'init_w' )
2229	!
2230	!-- level-of-detail 1 - read initialization profile
2231	IF ( init_3d%lod_w == 1 ) THEN
2232	ALLOCATE( init_3d%w_init(nzb:nzt+1) )
2233	init_3d%w_init = 0.0_wp
2234
2235	CALL get_variable( id_dynamic, 'init_w', &
2236	init_3d%w_init(nzb+1:nzt) )
2237
2238	!
2239	!-- level-of-detail 2 - read 3D initialization data
2240	ELSEIF ( init_3d%lod_w == 2 ) THEN
2241	DO i = nxl, nxr
2242	DO j = nys, nyn
2243	CALL get_variable( id_dynamic, 'init_w', i, j, &
2244	w(nzb+1:nzt,j,i) )
2245	ENDDO
2246	ENDDO
2247
2248	ENDIF
2249	init_3d%from_file_w = .TRUE.
2250	ENDIF
2251	!
2252	!-- Read potential temperature
2253	IF ( .NOT. neutral ) THEN
2254	IF ( check_existence( var_names, 'init_pt' ) ) THEN
2255	!
2256	!-- Read attributes for the fill value and level-of-detail
2257	CALL get_attribute( id_dynamic, char_fill, init_3d%fill_pt, &
2258	.FALSE., 'init_pt' )
2259	CALL get_attribute( id_dynamic, char_lod, init_3d%lod_pt, &
2260	.FALSE., 'init_pt' )
2261	!
2262	!-- level-of-detail 1 - read initialization profile
2263	IF ( init_3d%lod_pt == 1 ) THEN
2264	ALLOCATE( init_3d%pt_init(nzb:nzt+1) )
2265
2266	CALL get_variable( id_dynamic, 'init_pt', &
2267	init_3d%pt_init(nzb+1:nzt+1) )
2268	!
2269	!-- Set Neumann surface boundary condition for initial profil
2270	init_3d%pt_init(nzb) = init_3d%pt_init(nzb+1)
2271	!
2272	!-- level-of-detail 2 - read 3D initialization data
2273	ELSEIF ( init_3d%lod_pt == 2 ) THEN
2274	DO i = nxl, nxr
2275	DO j = nys, nyn
2276	CALL get_variable( id_dynamic, 'init_pt', i, j, &
2277	pt(nzb+1:nzt+1,j,i) )
2278	ENDDO
2279	ENDDO
2280
2281	ENDIF
2282	init_3d%from_file_pt = .TRUE.
2283	ENDIF
2284	ENDIF
2285	!
2286	!-- Read mixing ratio
2287	IF ( humidity ) THEN
2288	IF ( check_existence( var_names, 'init_qv' ) ) THEN
2289	!
2290	!-- Read attributes for the fill value and level-of-detail
2291	CALL get_attribute( id_dynamic, char_fill, init_3d%fill_q, &
2292	.FALSE., 'init_qv' )
2293	CALL get_attribute( id_dynamic, char_lod, init_3d%lod_q, &
2294	.FALSE., 'init_qv' )
2295	!
2296	!-- level-of-detail 1 - read initialization profile
2297	IF ( init_3d%lod_q == 1 ) THEN
2298	ALLOCATE( init_3d%q_init(nzb:nzt+1) )
2299
2300	CALL get_variable( id_dynamic, 'init_qv', &
2301	init_3d%q_init(nzb+1:nzt+1) )
2302	!
2303	!-- Set Neumann surface boundary condition for initial profil
2304	init_3d%q_init(nzb) = init_3d%q_init(nzb+1)
2305
2306	!
2307	!-- level-of-detail 2 - read 3D initialization data
2308	ELSEIF ( init_3d%lod_q == 2 ) THEN
2309	DO i = nxl, nxr
2310	DO j = nys, nyn
2311	CALL get_variable( id_dynamic, 'init_qv', i, j, &
2312	q(nzb+1:nzt+1,j,i) )
2313	ENDDO
2314	ENDDO
2315
2316	ENDIF
2317	init_3d%from_file_q = .TRUE.
2318	ENDIF
2319	ENDIF
2320	!
2321	!-- Read soil moisture
2322	IF ( land_surface ) THEN
2323
2324	IF ( check_existence( var_names, 'init_soil_m' ) ) THEN
2325	!
2326	!-- Read attributes for the fill value and level-of-detail
2327	CALL get_attribute( id_dynamic, char_fill, &
2328	init_3d%fill_msoil, &
2329	.FALSE., 'init_soil_m' )
2330	CALL get_attribute( id_dynamic, char_lod, &
2331	init_3d%lod_msoil, &
2332	.FALSE., 'init_soil_m' )
2333	!
2334	!-- level-of-detail 1 - read initialization profile
2335	IF ( init_3d%lod_msoil == 1 ) THEN
2336	ALLOCATE( init_3d%msoil_init(0:init_3d%nzs-1) )
2337
2338	CALL get_variable( id_dynamic, 'init_soil_m', &
2339	init_3d%msoil_init(0:init_3d%nzs-1) )
2340	!
2341	!-- level-of-detail 2 - read 3D initialization data
2342	ELSEIF ( init_3d%lod_msoil == 2 ) THEN
2343	ALLOCATE ( init_3d%msoil(0:init_3d%nzs-1,nys:nyn,nxl:nxr) )
2344	DO i = nxl, nxr
2345	DO j = nys, nyn
2346	CALL get_variable( id_dynamic, 'init_soil_m', i, j,&
2347	init_3d%msoil(0:init_3d%nzs-1,j,i) )
2348	ENDDO
2349	ENDDO
2350	ENDIF
2351	init_3d%from_file_msoil = .TRUE.
2352	ENDIF
2353	!
2354	!-- Read soil temperature
2355	IF ( check_existence( var_names, 'init_soil_t' ) ) THEN
2356	!
2357	!-- Read attributes for the fill value and level-of-detail
2358	CALL get_attribute( id_dynamic, char_fill, &
2359	init_3d%fill_tsoil, &
2360	.FALSE., 'init_soil_t' )
2361	CALL get_attribute( id_dynamic, char_lod, &
2362	init_3d%lod_tsoil, &
2363	.FALSE., 'init_soil_t' )
2364	!
2365	!-- level-of-detail 1 - read initialization profile
2366	IF ( init_3d%lod_tsoil == 1 ) THEN
2367	ALLOCATE( init_3d%tsoil_init(0:init_3d%nzs-1) )
2368
2369	CALL get_variable( id_dynamic, 'init_soil_t', &
2370	init_3d%tsoil_init(0:init_3d%nzs-1) )
2371
2372	!
2373	!-- level-of-detail 2 - read 3D initialization data
2374	ELSEIF ( init_3d%lod_tsoil == 2 ) THEN
2375	ALLOCATE ( init_3d%tsoil(0:init_3d%nzs-1,nys:nyn,nxl:nxr) )
2376	DO i = nxl, nxr
2377	DO j = nys, nyn
2378	CALL get_variable( id_dynamic, 'init_soil_t', i, j,&
2379	init_3d%tsoil(0:init_3d%nzs-1,j,i) )
2380	ENDDO
2381	ENDDO
2382	ENDIF
2383	init_3d%from_file_tsoil = .TRUE.
2384	ENDIF
2385	ENDIF
2386	!
2387	!-- Close input file
2388	CALL close_input_file( id_dynamic )
2389	#endif
2390	ENDIF
2391	#if defined( __parallel )
2392	CALL MPI_BARRIER( comm2d, ierr )
2393	#endif
2394	ENDDO
2395	!
2396	!-- End of CPU measurement
2397	CALL cpu_log( log_point_s(85), 'NetCDF input init', 'stop' )
2398	!
2399	!-- Finally, check if the input data has any fill values. Please note,
2400	!-- checks depend on the LOD of the input data.
2401	IF ( init_3d%from_file_u ) THEN
2402	check_passed = .TRUE.
2403	IF ( init_3d%lod_u == 1 ) THEN
2404	IF ( ANY( init_3d%u_init(nzb+1:nzt+1) == init_3d%fill_u ) ) &
2405	check_passed = .FALSE.
2406	ELSEIF ( init_3d%lod_u == 2 ) THEN
2407	IF ( ANY( u(nzb+1:nzt+1,nys:nyn,nxlu:nxr) == init_3d%fill_u ) ) &
2408	check_passed = .FALSE.
2409	ENDIF
2410	IF ( .NOT. check_passed ) THEN
2411	message_string = 'NetCDF input for u_init must not contain ' // &
2412	'any _FillValues'
2413	CALL message( 'netcdf_data_input_mod', 'NDI004', 2, 2, 0, 6, 0 )
2414	ENDIF
2415	ENDIF
2416
2417	IF ( init_3d%from_file_v ) THEN
2418	check_passed = .TRUE.
2419	IF ( init_3d%lod_v == 1 ) THEN
2420	IF ( ANY( init_3d%v_init(nzb+1:nzt+1) == init_3d%fill_v ) ) &
2421	check_passed = .FALSE.
2422	ELSEIF ( init_3d%lod_v == 2 ) THEN
2423	IF ( ANY( v(nzb+1:nzt+1,nysv:nyn,nxl:nxr) == init_3d%fill_v ) ) &
2424	check_passed = .FALSE.
2425	ENDIF
2426	IF ( .NOT. check_passed ) THEN
2427	message_string = 'NetCDF input for v_init must not contain ' // &
2428	'any _FillValues'
2429	CALL message( 'netcdf_data_input_mod', 'NDI005', 2, 2, 0, 6, 0 )
2430	ENDIF
2431	ENDIF
2432
2433	IF ( init_3d%from_file_w ) THEN
2434	check_passed = .TRUE.
2435	IF ( init_3d%lod_w == 1 ) THEN
2436	IF ( ANY( init_3d%w_init(nzb+1:nzt) == init_3d%fill_w ) ) &
2437	check_passed = .FALSE.
2438	ELSEIF ( init_3d%lod_w == 2 ) THEN
2439	IF ( ANY( w(nzb+1:nzt,nys:nyn,nxl:nxr) == init_3d%fill_w ) ) &
2440	check_passed = .FALSE.
2441	ENDIF
2442	IF ( .NOT. check_passed ) THEN
2443	message_string = 'NetCDF input for w_init must not contain ' // &
2444	'any _FillValues'
2445	CALL message( 'netcdf_data_input_mod', 'NDI006', 2, 2, 0, 6, 0 )
2446	ENDIF
2447	ENDIF
2448
2449	IF ( init_3d%from_file_pt ) THEN
2450	check_passed = .TRUE.
2451	IF ( init_3d%lod_pt == 1 ) THEN
2452	IF ( ANY( init_3d%pt_init(nzb+1:nzt+1) == init_3d%fill_pt ) ) &
2453	check_passed = .FALSE.
2454	ELSEIF ( init_3d%lod_pt == 2 ) THEN
2455	IF ( ANY( pt(nzb+1:nzt+1,nys:nyn,nxl:nxr) == init_3d%fill_pt ) ) &
2456	check_passed = .FALSE.
2457	ENDIF
2458	IF ( .NOT. check_passed ) THEN
2459	message_string = 'NetCDF input for pt_init must not contain ' // &
2460	'any _FillValues'
2461	CALL message( 'netcdf_data_input_mod', 'NDI007', 2, 2, 0, 6, 0 )
2462	ENDIF
2463	ENDIF
2464
2465	IF ( init_3d%from_file_q ) THEN
2466	check_passed = .TRUE.
2467	IF ( init_3d%lod_q == 1 ) THEN
2468	IF ( ANY( init_3d%q_init(nzb+1:nzt+1) == init_3d%fill_q ) ) &
2469	check_passed = .FALSE.
2470	ELSEIF ( init_3d%lod_q == 2 ) THEN
2471	IF ( ANY( q(nzb+1:nzt+1,nys:nyn,nxl:nxr) == init_3d%fill_q ) ) &
2472	check_passed = .FALSE.
2473	ENDIF
2474	IF ( .NOT. check_passed ) THEN
2475	message_string = 'NetCDF input for q_init must not contain ' // &
2476	'any _FillValues'
2477	CALL message( 'netcdf_data_input_mod', 'NDI008', 2, 2, 0, 6, 0 )
2478	ENDIF
2479	ENDIF
2480	!
2481	!-- Workaround for cyclic conditions. Please see above for further explanation.
2482	IF ( bc_lr_cyc .AND. nxl == 0 ) nxlu = nxl
2483	IF ( bc_ns_cyc .AND. nys == 0 ) nysv = nys
2484
2485	END SUBROUTINE netcdf_data_input_init_3d
2486
2487	!------------------------------------------------------------------------------!
2488	! Description:
2489	! ------------
2490	!> Reads data at lateral and top boundaries derived from larger-scale model
2491	!> (COSMO) by Inifor.
2492	!------------------------------------------------------------------------------!
2493	SUBROUTINE netcdf_data_input_lsf
2494
2495	USE control_parameters, &
2496	ONLY: bc_lr_cyc, bc_ns_cyc, force_bound_l, force_bound_n, &
2497	force_bound_r, force_bound_s, &
2498	forcing, humidity, message_string, neutral, simulated_time
2499
2500
2501	USE indices, &
2502	ONLY: nxl, nxlu, nxr, nyn, nys, nysv, nzb, nzt
2503
2504	IMPLICIT NONE
2505
2506
2507	INTEGER(iwp) :: i !< running index along x-direction
2508	INTEGER(iwp) :: ii !< running index for IO blocks
2509	INTEGER(iwp) :: id_dynamic !< NetCDF id of dynamic input file
2510	INTEGER(iwp) :: j !< running index along y-direction
2511	INTEGER(iwp) :: k !< running index along z-direction
2512	INTEGER(iwp) :: num_vars !< number of variables in netcdf input file
2513	INTEGER(iwp) :: t !< running index time dimension
2514
2515	REAL(wp) :: dum !< dummy variable to skip columns while reading topography file
2516
2517	force%from_file = MERGE( .TRUE., .FALSE., input_pids_dynamic )
2518	!
2519	!-- Skip input if no forcing from larger-scale models is applied.
2520	IF ( .NOT. forcing ) RETURN
2521
2522	!
2523	!-- CPU measurement
2524	CALL cpu_log( log_point_s(86), 'NetCDF input forcing', 'start' )
2525
2526	DO ii = 0, io_blocks-1
2527	IF ( ii == io_group ) THEN
2528	#if defined ( __netcdf )
2529	!
2530	!-- Open file in read-only mode
2531	CALL open_read_file( TRIM( input_file_dynamic ) // &
2532	TRIM( coupling_char ), id_dynamic )
2533	!
2534	!-- Initialize INIFOR forcing.
2535	IF ( .NOT. force%init ) THEN
2536	!
2537	!-- At first, inquire all variable names.
2538	CALL inquire_num_variables( id_dynamic, num_vars )
2539	!
2540	!-- Allocate memory to store variable names.
2541	ALLOCATE( force%var_names(1:num_vars) )
2542	CALL inquire_variable_names( id_dynamic, force%var_names )
2543	!
2544	!-- Read time dimension, allocate memory and finally read time array
2545	CALL get_dimension_length( id_dynamic, force%nt, 'time' )
2546
2547	IF ( check_existence( force%var_names, 'time' ) ) THEN
2548	ALLOCATE( force%time(0:force%nt-1) )
2549	CALL get_variable( id_dynamic, 'time', force%time )
2550	ENDIF
2551	!
2552	!-- Read vertical dimension of scalar und w grid
2553	CALL get_dimension_length( id_dynamic, force%nzu, 'z' )
2554	CALL get_dimension_length( id_dynamic, force%nzw, 'zw' )
2555
2556	IF ( check_existence( force%var_names, 'z' ) ) THEN
2557	ALLOCATE( force%zu_atmos(1:force%nzu) )
2558	CALL get_variable( id_dynamic, 'z', force%zu_atmos )
2559	ENDIF
2560	IF ( check_existence( force%var_names, 'zw' ) ) THEN
2561	ALLOCATE( force%zw_atmos(1:force%nzw) )
2562	CALL get_variable( id_dynamic, 'zw', force%zw_atmos )
2563	ENDIF
2564
2565	!
2566	!-- Read surface pressure
2567	IF ( check_existence( force%var_names, &
2568	'surface_forcing_surface_pressure' ) ) THEN
2569	ALLOCATE( force%surface_pressure(0:force%nt-1) )
2570	CALL get_variable( id_dynamic, &
2571	'surface_forcing_surface_pressure', &
2572	force%surface_pressure )
2573	ENDIF
2574	!
2575	!-- Set control flag to indicate that initialization is already done
2576	force%init = .TRUE.
2577
2578	ENDIF
2579
2580	!
2581	!-- Obtain time index for current input starting at 0.
2582	!-- @todo: At the moment time, in INIFOR and simulated time correspond
2583	!-- to each other. If required, adjust to daytime.
2584	force%tind = MINLOC( ABS( force%time - simulated_time ), DIM = 1 )&
2585	- 1
2586	force%tind_p = force%tind + 1
2587	!
2588	!-- Read geostrophic wind components. In case of forcing, this is only
2589	!-- required if cyclic boundary conditions are applied.
2590	IF ( bc_lr_cyc .AND. bc_ns_cyc ) THEN
2591	DO t = force%tind, force%tind_p
2592	CALL get_variable_pr( id_dynamic, 'tend_ug', t+1, &
2593	force%ug(t-force%tind,:) )
2594	CALL get_variable_pr( id_dynamic, 'tend_vg', t+1, &
2595	force%ug(t-force%tind,:) )
2596	ENDDO
2597	ENDIF
2598	!
2599	!-- Read data at lateral and top boundaries. Please note, at left and
2600	!-- right domain boundary, yz-layers are read for u, v, w, pt and q.
2601	!-- For the v-component, the data starts at nysv, while for the other
2602	!-- quantities the data starts at nys. This is equivalent at the north
2603	!-- and south domain boundary for the u-component.
2604	!-- The function get_variable_bc assumes the start indices with respect
2605	!-- to the netcdf file convention (data starts at index 1). For this
2606	!-- reason, nys+1 / nxl+1 are passed instead of nys / nxl. For the
2607	!-- the u- and v-component at the north/south, and left/right boundary,
2608	!-- nxlu and nysv are passed, respectively, since these always starts
2609	!-- at index 1 in case of forcing.
2610
2611	IF ( force_bound_l ) THEN
2612	DO j = nys, nyn
2613	DO t = force%tind, force%tind_p
2614	CALL get_variable_bc( id_dynamic, 'ls_forcing_left_u', &
2615	t+1, &
2616	nzb+1, nzt+1-(nzb+1)+1, &
2617	j+1, 1, &
2618	force%u_left(t-force%tind,nzb+1:nzt+1,j) )
2619	ENDDO
2620	ENDDO
2621
2622	DO j = nysv, nyn
2623	DO t = force%tind, force%tind_p
2624	CALL get_variable_bc( id_dynamic, 'ls_forcing_left_v', &
2625	t+1, &
2626	nzb+1, nzt+1-(nzb+1)+1, &
2627	j, 1, &
2628	force%v_left(t-force%tind,nzb+1:nzt+1,j) )
2629	ENDDO
2630	ENDDO
2631	DO j = nys, nyn
2632	DO t = force%tind, force%tind_p
2633	CALL get_variable_bc( id_dynamic, &
2634	'ls_forcing_left_w', &
2635	t+1, &
2636	nzb+1, nzt-(nzb+1) + 1, &
2637	j+1, 1, &
2638	force%w_left(t-force%tind,nzb+1:nzt,j) )
2639	ENDDO
2640	ENDDO
2641	IF ( .NOT. neutral ) THEN
2642	DO j = nys, nyn
2643	DO t = force%tind, force%tind_p
2644	CALL get_variable_bc( id_dynamic, &
2645	'ls_forcing_left_pt', &
2646	t+1, &
2647	nzb+1, nzt+1-(nzb+1)+1, &
2648	j+1, 1, &
2649	force%pt_left(t-force%tind,nzb+1:nzt+1,j) )
2650	ENDDO
2651	ENDDO
2652	ENDIF
2653	IF ( humidity ) THEN
2654	DO j = nys, nyn
2655	DO t = force%tind, force%tind_p
2656	CALL get_variable_bc( id_dynamic, &
2657	'ls_forcing_left_qv', &
2658	t+1, &
2659	nzb+1, nzt+1-(nzb+1)+1, &
2660	j+1, 1, &
2661	force%q_left(t-force%tind,nzb+1:nzt+1,j) )
2662	ENDDO
2663	ENDDO
2664	ENDIF
2665	ENDIF
2666
2667	IF ( force_bound_r ) THEN
2668	DO j = nys, nyn
2669	DO t = force%tind, force%tind_p
2670	CALL get_variable_bc( id_dynamic, 'ls_forcing_right_u', &
2671	t+1, &
2672	nzb+1, nzt+1-(nzb+1)+1, &
2673	j+1, 1, &
2674	force%u_right(t-force%tind,nzb+1:nzt+1,j) )
2675	ENDDO
2676	ENDDO
2677	DO j = nysv, nyn
2678	DO t = force%tind, force%tind_p
2679	CALL get_variable_bc( id_dynamic, 'ls_forcing_right_v', &
2680	t+1, &
2681	nzb+1, nzt+1-(nzb+1)+1, &
2682	j, 1, &
2683	force%v_right(t-force%tind,nzb+1:nzt+1,j) )
2684	ENDDO
2685	ENDDO
2686	DO j = nys, nyn
2687	DO t = force%tind, force%tind_p
2688	CALL get_variable_bc( id_dynamic, 'ls_forcing_right_w', &
2689	t+1, &
2690	nzb+1, nzt-(nzb+1)+1, &
2691	j+1, 1, &
2692	force%w_right(t-force%tind,nzb+1:nzt,j) )
2693	ENDDO
2694	ENDDO
2695	IF ( .NOT. neutral ) THEN
2696	DO j = nys, nyn
2697	DO t = force%tind, force%tind_p
2698	CALL get_variable_bc( id_dynamic, &
2699	'ls_forcing_right_pt', &
2700	t+1, &
2701	nzb+1, nzt+1-(nzb+1)+1, &
2702	j+1, 1, &
2703	force%pt_right(t-force%tind,nzb+1:nzt+1,j) )
2704	ENDDO
2705	ENDDO
2706	ENDIF
2707	IF ( humidity ) THEN
2708	DO j = nys, nyn
2709	DO t = force%tind, force%tind_p
2710	CALL get_variable_bc( id_dynamic, &
2711	'ls_forcing_right_qv', &
2712	t+1, &
2713	nzb+1, nzt+1-(nzb+1)+1, &
2714	j+1, 1, &
2715	force%q_right(t-force%tind,nzb+1:nzt+1,j) )
2716	ENDDO
2717	ENDDO
2718	ENDIF
2719	ENDIF
2720
2721	IF ( force_bound_n ) THEN
2722	DO i = nxlu, nxr
2723	DO t = force%tind, force%tind_p
2724	CALL get_variable_bc( id_dynamic, 'ls_forcing_north_u', &
2725	t+1, &
2726	nzb+1, nzt+1-(nzb+1)+1, &
2727	i, 1, &
2728	force%u_north(t-force%tind,nzb+1:nzt+1,i) )
2729	ENDDO
2730	ENDDO
2731	DO i = nxl, nxr
2732	DO t = force%tind, force%tind_p
2733	CALL get_variable_bc( id_dynamic, 'ls_forcing_north_v', &
2734	t+1, &
2735	nzb+1, nzt+1-(nzb+1)+1, &
2736	i+1, 1, &
2737	force%v_north(t-force%tind,nzb+1:nzt+1,i) )
2738	ENDDO
2739	ENDDO
2740	DO i = nxl, nxr
2741	DO t = force%tind, force%tind_p
2742	CALL get_variable_bc( id_dynamic, 'ls_forcing_north_w', &
2743	t+1, &
2744	nzb+1, nzt-(nzb+1)+1, &
2745	i+1, 1, &
2746	force%w_north(t-force%tind,nzb+1:nzt,i) )
2747	ENDDO
2748	ENDDO
2749	IF ( .NOT. neutral ) THEN
2750	DO i = nxl, nxr
2751	DO t = force%tind, force%tind_p
2752	CALL get_variable_bc( id_dynamic, &
2753	'ls_forcing_north_pt', &
2754	t+1, &
2755	nzb+1, nzt+1-(nzb+1)+1, &
2756	i+1, 1, &
2757	force%pt_north(t-force%tind,nzb+1:nzt+1,i) )
2758	ENDDO
2759	ENDDO
2760	ENDIF
2761	IF ( humidity ) THEN
2762	DO i = nxl, nxr
2763	DO t = force%tind, force%tind_p
2764	CALL get_variable_bc( id_dynamic, &
2765	'ls_forcing_north_qv', &
2766	t+1, &
2767	nzb+1, nzt+1-(nzb+1)+1, &
2768	i+1, 1, &
2769	force%q_north(t-force%tind,nzb+1:nzt+1,i) )
2770	ENDDO
2771	ENDDO
2772	ENDIF
2773	ENDIF
2774
2775	IF ( force_bound_s ) THEN
2776	DO i = nxlu, nxr
2777	DO t = force%tind, force%tind_p
2778	CALL get_variable_bc( id_dynamic, 'ls_forcing_south_u', &
2779	t+1, &
2780	nzb+1, nzt+1-(nzb+1)+1, &
2781	i, 1, &
2782	force%u_south(t-force%tind,nzb+1:nzt+1,i) )
2783	ENDDO
2784	ENDDO
2785	DO i = nxl, nxr
2786	DO t = force%tind, force%tind_p
2787	CALL get_variable_bc( id_dynamic, 'ls_forcing_south_v', &
2788	t+1, &
2789	nzb+1, nzt+1-(nzb+1)+1, &
2790	i+1, 1, &
2791	force%v_south(t-force%tind,nzb+1:nzt+1,i) )
2792	ENDDO
2793	ENDDO
2794	DO i = nxl, nxr
2795	DO t = force%tind, force%tind_p
2796	CALL get_variable_bc( id_dynamic, 'ls_forcing_south_w', &
2797	t+1, &
2798	nzb+1, nzt-(nzb+1)+1, &
2799	i+1, 1, &
2800	force%w_south(t-force%tind,nzb+1:nzt,i) )
2801	ENDDO
2802	ENDDO
2803	IF ( .NOT. neutral ) THEN
2804	DO i = nxl, nxr
2805	DO t = force%tind, force%tind_p
2806	CALL get_variable_bc( id_dynamic, &
2807	'ls_forcing_south_pt', &
2808	t+1, &
2809	nzb+1, nzt+1-(nzb+1)+1, &
2810	i+1, 1, &
2811	force%pt_south(t-force%tind,nzb+1:nzt+1,i) )
2812	ENDDO
2813	ENDDO
2814	ENDIF
2815	IF ( humidity ) THEN
2816	DO i = nxl, nxr
2817	DO t = force%tind, force%tind_p
2818	CALL get_variable_bc( id_dynamic, &
2819	'ls_forcing_south_qv', &
2820	t+1, &
2821	nzb+1, nzt+1-(nzb+1)+1, &
2822	i+1, 1, &
2823	force%q_south(t-force%tind,nzb+1:nzt+1,i) )
2824	ENDDO
2825	ENDDO
2826	ENDIF
2827	ENDIF
2828	!
2829	!-- Top boundary
2830	DO i = nxlu, nxr
2831	DO t = force%tind, force%tind_p
2832	CALL get_variable_bc( id_dynamic, 'ls_forcing_top_u', &
2833	t+1, &
2834	nys+1, nyn-nys+1, &
2835	i, 1, &
2836	force%u_top(t-force%tind,nys:nyn,i) )
2837	ENDDO
2838	ENDDO
2839	DO i = nxl, nxr
2840	DO t = force%tind, force%tind_p
2841	CALL get_variable_bc( id_dynamic, 'ls_forcing_top_v', &
2842	t+1, &
2843	nysv, nyn-nysv+1, &
2844	i+1, 1, &
2845	force%v_top(t-force%tind,nysv:nyn,i) )
2846	ENDDO
2847	ENDDO
2848	DO i = nxl, nxr
2849	DO t = force%tind, force%tind_p
2850	CALL get_variable_bc( id_dynamic, 'ls_forcing_top_w', &
2851	t+1, &
2852	nys+1, nyn-nys+1, &
2853	i+1, 1, &
2854	force%w_top(t-force%tind,nys:nyn,i) )
2855	ENDDO
2856	ENDDO
2857	IF ( .NOT. neutral ) THEN
2858	DO i = nxl, nxr
2859	DO t = force%tind, force%tind_p
2860	CALL get_variable_bc( id_dynamic, 'ls_forcing_top_pt', &
2861	t+1, &
2862	nys+1, nyn-nys+1, &
2863	i+1, 1, &
2864	force%pt_top(t-force%tind,nys:nyn,i) )
2865	ENDDO
2866	ENDDO
2867	ENDIF
2868	IF ( humidity ) THEN
2869	DO i = nxl, nxr
2870	DO t = force%tind, force%tind_p
2871	CALL get_variable_bc( id_dynamic, 'ls_forcing_top_qv', &
2872	t+1, &
2873	nys+1, nyn-nys+1, &
2874	i+1, 1, &
2875	force%q_top(t-force%tind,nys:nyn,i) )
2876	ENDDO
2877	ENDDO
2878	ENDIF
2879
2880	!
2881	!-- Close input file
2882	CALL close_input_file( id_dynamic )
2883	#endif
2884	ENDIF
2885	#if defined( __parallel )
2886	CALL MPI_BARRIER( comm2d, ierr )
2887	#endif
2888	ENDDO
2889
2890	!
2891	!-- End of CPU measurement
2892	CALL cpu_log( log_point_s(86), 'NetCDF input forcing', 'stop' )
2893
2894	!
2895	!-- Finally, after data input set control flag indicating that vertical
2896	!-- inter- and/or extrapolation is required.
2897	!-- Please note, inter/extrapolation of INIFOR data is only a workaroud,
2898	!-- as long as INIFOR delivers vertically equidistant data.
2899	force%interpolated = .FALSE.
2900
2901	END SUBROUTINE netcdf_data_input_lsf
2902
2903
2904	!------------------------------------------------------------------------------!
2905	! Description:
2906	! ------------
2907	!> Checks input file for consistency and minimum requirements.
2908	!------------------------------------------------------------------------------!
2909	SUBROUTINE netcdf_data_input_check_dynamic
2910
2911	USE control_parameters, &
2912	ONLY: initializing_actions, forcing, message_string
2913
2914	IMPLICIT NONE
2915
2916	!
2917	!-- In case of forcing, check whether dynamic input file is present
2918	IF ( .NOT. input_pids_dynamic .AND. forcing ) THEN
2919	message_string = 'forcing = .TRUE. requires dynamic input file ' // &
2920	TRIM( input_file_dynamic ) // TRIM( coupling_char )
2921	CALL message( 'netcdf_data_input_mod', 'NDI009', 1, 2, 0, 6, 0 )
2922	ENDIF
2923	!
2924	!-- Dynamic input file must also be present if initialization via inifor is
2925	!-- prescribed.
2926	IF ( .NOT. input_pids_dynamic .AND. &
2927	TRIM( initializing_actions ) == 'inifor' ) THEN
2928	message_string = 'initializing_actions = inifor requires dynamic ' //&
2929	'input file ' // TRIM( input_file_dynamic ) // &
2930	TRIM( coupling_char )
2931	CALL message( 'netcdf_data_input_mod', 'NDI010', 1, 2, 0, 6, 0 )
2932	ENDIF
2933
2934	END SUBROUTINE netcdf_data_input_check_dynamic
2935
2936	!------------------------------------------------------------------------------!
2937	! Description:
2938	! ------------
2939	!> Checks input file for consistency and minimum requirements.
2940	!------------------------------------------------------------------------------!
2941	SUBROUTINE netcdf_data_input_check_static
2942
2943	USE arrays_3d, &
2944	ONLY: zu
2945
2946	USE control_parameters, &
2947	ONLY: land_surface, message_string, urban_surface
2948
2949	USE grid_variables, &
2950	ONLY: dx, dy
2951
2952	USE indices, &
2953	ONLY: nx, nxl, nxr, ny, nyn, nys
2954
2955	IMPLICIT NONE
2956
2957	INTEGER(iwp) :: i !< loop index along x-direction
2958	INTEGER(iwp) :: j !< loop index along y-direction
2959	INTEGER(iwp) :: n_surf !< number of different surface types at given location
2960
2961	LOGICAL :: check_passed !< flag indicating if a check passed
2962
2963	!
2964	!-- Return if no static input file is available
2965	IF ( .NOT. input_pids_static ) RETURN
2966	!
2967	!-- Check whether dimension size in input file matches the model dimensions
2968	IF ( dim_static%nx-1 /= nx .OR. dim_static%ny-1 /= ny ) THEN
2969	message_string = 'Static input file: horizontal dimension in ' // &
2970	'x- and/or y-direction ' // &
2971	'do not match the respective model dimension'
2972	CALL message( 'netcdf_data_input_mod', 'NDI011', 1, 2, 0, 6, 0 )
2973	ENDIF
2974	!
2975	!-- Check if grid spacing of provided input data matches the respective
2976	!-- grid spacing in the model.
2977	IF ( dim_static%x(1) - dim_static%x(0) /= dx .OR. &
2978	dim_static%y(1) - dim_static%y(0) /= dy ) THEN
2979	message_string = 'Static input file: horizontal grid spacing ' // &
2980	'in x- and/or y-direction ' // &
2981	'do not match the respective model grid spacing.'
2982	CALL message( 'netcdf_data_input_mod', 'NDI012', 1, 2, 0, 6, 0 )
2983	ENDIF
2984	!
2985	!-- Check orography for fill-values. For the moment, give an error message.
2986	!-- More advanced methods, e.g. a nearest neighbor algorithm as used in GIS
2987	!-- systems might be implemented later.
2988	!-- Please note, if no terrain height is provided, it is set to 0.
2989	IF ( ANY( terrain_height_f%var == terrain_height_f%fill ) ) THEN
2990	message_string = 'NetCDF variable zt is not ' // &
2991	'allowed to have missing data'
2992	CALL message( 'netcdf_data_input_mod', 'NDI013', 2, 2, 0, 6, 0 )
2993	ENDIF
2994	!
2995	!-- If 3D buildings are read, check if building information is consistent
2996	!-- to numeric grid.
2997	IF ( buildings_f%from_file ) THEN
2998	IF ( buildings_f%lod == 2 ) THEN
2999	IF ( buildings_f%nz > SIZE( zu ) ) THEN
3000	message_string = 'Reading 3D building data - too much ' // &
3001	'data points along the vertical coordinate.'
3002	CALL message( 'netcdf_data_input_mod', 'NDI014', 2, 2, 0, 6, 0 )
3003	ENDIF
3004
3005	IF ( ANY( buildings_f%z(0:buildings_f%nz-1) /= &
3006	zu(0:buildings_f%nz-1) ) ) THEN
3007	message_string = 'Reading 3D building data - vertical ' // &
3008	'coordinate do not match numeric grid.'
3009	CALL message( 'netcdf_data_input_mod', 'NDI015', 2, 2, 0, 6, 0 )
3010	ENDIF
3011	ENDIF
3012	ENDIF
3013
3014	!
3015	!-- Skip further checks concerning buildings and natural surface properties
3016	!-- if no urban surface and land surface model are applied.
3017	IF ( .NOT. land_surface .OR. .NOT. urban_surface ) RETURN
3018	!
3019	!-- Check for minimum requirement of surface-classification data in case
3020	!-- static input file is used.
3021	IF ( ( .NOT. vegetation_type_f%from_file .OR. &
3022	.NOT. pavement_type_f%from_file .OR. &
3023	.NOT. water_type_f%from_file .OR. &
3024	.NOT. soil_type_f%from_file ) .OR. &
3025	( urban_surface .AND. .NOT. building_type_f%from_file ) ) THEN
3026	message_string = 'Minimum requirement for surface classification ' //&
3027	'is not fulfilled. At least ' // &
3028	'vegetation_type, pavement_type, ' // &
3029	'soil_type and water_type are '// &
3030	'required. If urban-surface model is applied, ' // &
3031	'also building_type ist required'
3032	CALL message( 'netcdf_data_input_mod', 'NDI016', 1, 2, 0, 6, 0 )
3033	ENDIF
3034	!
3035	!-- Check for general availability of input variables.
3036	!-- If vegetation_type is 0 at any location, vegetation_pars as well as
3037	!-- root_area_dens_s are required.
3038	IF ( vegetation_type_f%from_file ) THEN
3039	IF ( ANY( vegetation_type_f%var == 0 ) ) THEN
3040	IF ( .NOT. vegetation_pars_f%from_file ) THEN
3041	message_string = 'If vegegation_type = 0 at any location, ' // &
3042	'vegetation_pars is required'
3043	CALL message( 'netcdf_data_input_mod', 'NDI017', 2, 2, 0, 6, 0 )
3044	ENDIF
3045	IF ( .NOT. root_area_density_lsm_f%from_file ) THEN
3046	message_string = 'If vegegation_type = 0 at any location, ' // &
3047	'root_area_dens_s is required'
3048	CALL message( 'netcdf_data_input_mod', 'NDI018', 2, 2, 0, 6, 0 )
3049	ENDIF
3050	ENDIF
3051	ENDIF
3052	!
3053	!-- If soil_type is zero at any location, soil_pars is required.
3054	IF ( soil_type_f%from_file ) THEN
3055	check_passed = .TRUE.
3056	IF ( ALLOCATED( soil_type_f%var_2d ) ) THEN
3057	IF ( ANY( soil_type_f%var_2d == 0 ) ) THEN
3058	IF ( .NOT. soil_pars_f%from_file ) check_passed = .FALSE.
3059	ENDIF
3060	ELSE
3061	IF ( ANY( soil_type_f%var_3d == 0 ) ) THEN
3062	IF ( .NOT. soil_pars_f%from_file ) check_passed = .FALSE.
3063	ENDIF
3064	ENDIF
3065	IF ( .NOT. check_passed ) THEN
3066	message_string = 'If soil_type = 0 at any location, ' // &
3067	'soil_pars is required'
3068	CALL message( 'netcdf_data_input_mod', 'NDI019', 2, 2, 0, 6, 0 )
3069	ENDIF
3070	ENDIF
3071	!
3072	!-- If building_type is zero at any location, building_pars is required.
3073	IF ( building_type_f%from_file ) THEN
3074	IF ( ANY( building_type_f%var == 0 ) ) THEN
3075	IF ( .NOT. building_pars_f%from_file ) THEN
3076	message_string = 'If building_type = 0 at any location, ' // &
3077	'building_pars is required'
3078	CALL message( 'netcdf_data_input_mod', 'NDI020', 2, 2, 0, 6, 0 )
3079	ENDIF
3080	ENDIF
3081	ENDIF
3082	!
3083	!-- If albedo_type is zero at any location, albedo_pars is required.
3084	IF ( albedo_type_f%from_file ) THEN
3085	IF ( ANY( albedo_type_f%var == 0 ) ) THEN
3086	IF ( .NOT. albedo_pars_f%from_file ) THEN
3087	message_string = 'If albedo_type = 0 at any location, ' // &
3088	'albedo_pars is required'
3089	CALL message( 'netcdf_data_input_mod', 'NDI021', 2, 2, 0, 6, 0 )
3090	ENDIF
3091	ENDIF
3092	ENDIF
3093	!
3094	!-- If pavement_type is zero at any location, pavement_pars is required.
3095	IF ( pavement_type_f%from_file ) THEN
3096	IF ( ANY( pavement_type_f%var == 0 ) ) THEN
3097	IF ( .NOT. pavement_pars_f%from_file ) THEN
3098	message_string = 'If pavement_type = 0 at any location, ' // &
3099	'pavement_pars is required'
3100	CALL message( 'netcdf_data_input_mod', 'NDI022', 2, 2, 0, 6, 0 )
3101	ENDIF
3102	ENDIF
3103	ENDIF
3104	!
3105	!-- If pavement_type is zero at any location, also pavement_subsurface_pars
3106	!-- is required.
3107	IF ( pavement_type_f%from_file ) THEN
3108	IF ( ANY( pavement_type_f%var == 0 ) ) THEN
3109	IF ( .NOT. pavement_subsurface_pars_f%from_file ) THEN
3110	message_string = 'If pavement_type = 0 at any location, ' // &
3111	'pavement_subsurface_pars is required'
3112	CALL message( 'netcdf_data_input_mod', 'NDI023', 2, 2, 0, 6, 0 )
3113	ENDIF
3114	ENDIF
3115	ENDIF
3116	!
3117	!-- If water_type is zero at any location, water_pars is required.
3118	IF ( water_type_f%from_file ) THEN
3119	IF ( ANY( water_type_f%var == 0 ) ) THEN
3120	IF ( .NOT. water_pars_f%from_file ) THEN
3121	message_string = 'If water_type = 0 at any location, ' // &
3122	'water_pars is required'
3123	CALL message( 'netcdf_data_input_mod', 'NDI024', 2, 2, 0, 6, 0 )
3124	ENDIF
3125	ENDIF
3126	ENDIF
3127	!
3128	!-- Check for local consistency of the input data.
3129	DO i = nxl, nxr
3130	DO j = nys, nyn
3131	!
3132	!-- For each (y,x)-location at least one of the parameters
3133	!-- vegetation_type, pavement_type, building_type, or water_type
3134	!-- must be set to a nonÂmissing value.
3135	IF ( vegetation_type_f%var(j,i) == vegetation_type_f%fill .AND. &
3136	pavement_type_f%var(j,i) == pavement_type_f%fill .AND. &
3137	building_type_f%var(j,i) == building_type_f%fill .AND. &
3138	water_type_f%var(j,i) == water_type_f%fill ) THEN
3139	WRITE( message_string, * ) 'At least one of the parameters '// &
3140	'vegetation_type, pavement_type, ' // &
3141	'building_type, or water_type must be set '// &
3142	'to a non-missing value. Grid point: ', j, i
3143	CALL message( 'netcdf_data_input_mod', 'NDI025', 2, 2, 0, 6, 0 )
3144	ENDIF
3145	!
3146	!-- Note that a soil_type is required for each location (y,x) where
3147	!-- either vegetation_type or pavement_type is a nonÂmissing value.
3148	IF ( ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill .OR. &
3149	pavement_type_f%var(j,i) /= pavement_type_f%fill ) ) THEN
3150	check_passed = .TRUE.
3151	IF ( ALLOCATED( soil_type_f%var_2d ) ) THEN
3152	IF ( soil_type_f%var_2d(j,i) == soil_type_f%fill ) &
3153	check_passed = .FALSE.
3154	ELSE
3155	IF ( ANY( soil_type_f%var_3d(:,j,i) == soil_type_f%fill) ) &
3156	check_passed = .FALSE.
3157	ENDIF
3158
3159	IF ( .NOT. check_passed ) THEN
3160	message_string = 'soil_type is required for each '// &
3161	'location (y,x) where vegetation_type or ' // &
3162	'pavement_type is a non-missing value.'
3163	CALL message( 'netcdf_data_input_mod', 'NDI026', &
3164	2, 2, 0, 6, 0 )
3165	ENDIF
3166	ENDIF
3167	!
3168	!-- Check for consistency of surface fraction. If more than one type
3169	!-- is set, surface fraction need to be given and the sum must not
3170	!-- be larger than 1.
3171	n_surf = 0
3172	IF ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill ) &
3173	n_surf = n_surf + 1
3174	IF ( water_type_f%var(j,i) /= water_type_f%fill ) &
3175	n_surf = n_surf + 1
3176	IF ( pavement_type_f%var(j,i) /= pavement_type_f%fill ) &
3177	n_surf = n_surf + 1
3178
3179	IF ( n_surf > 1 ) THEN
3180	IF ( ANY ( surface_fraction_f%frac(:,j,i) == &
3181	surface_fraction_f%fill ) ) THEN
3182	message_string = 'If more than one surface type is ' // &
3183	'given at a location, surface_fraction ' // &
3184	'must be provided.'
3185	CALL message( 'netcdf_data_input_mod', 'NDI027', &
3186	2, 2, 0, 6, 0 )
3187	ENDIF
3188	IF ( SUM ( surface_fraction_f%frac(:,j,i) ) > 1.0_wp ) THEN
3189	message_string = 'surface_fraction must not exceed 1'
3190	CALL message( 'netcdf_data_input_mod', 'NDI028', &
3191	2, 2, 0, 6, 0 )
3192	ENDIF
3193	ENDIF
3194	!
3195	!-- Check for further mismatches, e.g. vegetation_type is set but
3196	!-- surface vegetation fraction is zero.
3197	IF ( ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill .AND.&
3198	( surface_fraction_f%frac(ind_veg_wall,j,i) == 0.0_wp .OR. &
3199	surface_fraction_f%frac(ind_veg_wall,j,i) == &
3200	surface_fraction_f%fill ) &
3201	) .OR. &
3202	( pavement_type_f%var(j,i) /= pavement_type_f%fill .AND. &
3203	( surface_fraction_f%frac(ind_pav_green,j,i) == 0.0_wp .OR. &
3204	surface_fraction_f%frac(ind_pav_green,j,i) == &
3205	surface_fraction_f%fill ) &
3206	) .OR. &
3207	( water_type_f%var(j,i) /= water_type_f%fill .AND. &
3208	( surface_fraction_f%frac(ind_wat_win,j,i) == 0.0_wp .OR. &
3209	surface_fraction_f%frac(ind_wat_win,j,i) == &
3210	surface_fraction_f%fill ) &
3211	) ) THEN
3212	WRITE( message_string, * ) 'Mismatch in setting of ' // &
3213	'surface_fraction. Vegetation-, pavement-, or '// &
3214	'water surface is given at (i,j) = ( ', i, j, &
3215	' ), but surface fraction is 0 for the given type.'
3216	CALL message( 'netcdf_data_input_mod', 'NDI029', &
3217	2, 2, 0, 6, 0 )
3218	ENDIF
3219	!
3220	!-- Check for further mismatches, e.g. vegetation_type is not set
3221	!-- surface vegetation fraction is non-zero.
3222	IF ( ( vegetation_type_f%var(j,i) == vegetation_type_f%fill .AND.&
3223	( surface_fraction_f%frac(ind_veg_wall,j,i) /= 0.0_wp .AND. &
3224	surface_fraction_f%frac(ind_veg_wall,j,i) /= &
3225	surface_fraction_f%fill ) &
3226	) .OR. &
3227	( pavement_type_f%var(j,i) == pavement_type_f%fill .AND. &
3228	( surface_fraction_f%frac(ind_pav_green,j,i) /= 0.0_wp .AND. &
3229	surface_fraction_f%frac(ind_pav_green,j,i) /= &
3230	surface_fraction_f%fill ) &
3231	) .OR. &
3232	( water_type_f%var(j,i) == water_type_f%fill .AND. &
3233	( surface_fraction_f%frac(ind_wat_win,j,i) /= 0.0_wp .AND. &
3234	surface_fraction_f%frac(ind_wat_win,j,i) /= &
3235	surface_fraction_f%fill ) &
3236	) ) THEN
3237	WRITE( message_string, * ) 'Mismatch in setting of ' // &
3238	'surface_fraction. Vegetation-, pavement-, or '// &
3239	'water surface is not given at (i,j) = ( ', i, j, &
3240	' ), but surface fraction is not 0 for the ' // &
3241	'given type.'
3242	CALL message( 'netcdf_data_input_mod', 'NDI030', &
3243	2, 2, 0, 6, 0 )
3244	ENDIF
3245	!
3246	!-- Check vegetation_pars. If vegetation_type is 0, all parameters
3247	!-- need to be set, otherwise, single parameters set by
3248	!-- vegetation_type can be overwritten.
3249	IF ( vegetation_type_f%from_file ) THEN
3250	IF ( vegetation_type_f%var(j,i) == 0 ) THEN
3251	IF ( ANY( vegetation_pars_f%pars_xy(:,j,i) == &
3252	vegetation_pars_f%fill ) ) THEN
3253	message_string = 'If vegetation_type(y,x) = 0, all ' // &
3254	'parameters of vegetation_pars at '// &
3255	'this location must be set.'
3256	CALL message( 'netcdf_data_input_mod', 'NDI031', &
3257	2, 2, 0, 6, 0 )
3258	ENDIF
3259	ENDIF
3260	ENDIF
3261	!
3262	!-- Check root distribution. If vegetation_type is 0, all levels must
3263	!-- be set.
3264	IF ( vegetation_type_f%from_file ) THEN
3265	IF ( vegetation_type_f%var(j,i) == 0 ) THEN
3266	IF ( ANY( root_area_density_lsm_f%var(:,j,i) == &
3267	root_area_density_lsm_f%fill ) ) THEN
3268	message_string = 'If vegetation_type(y,x) = 0, all ' // &
3269	'levels of root_area_dens_s ' // &
3270	'must be set at this location.'
3271	CALL message( 'netcdf_data_input_mod', 'NDI032', &
3272	2, 2, 0, 6, 0 )
3273	ENDIF
3274	ENDIF
3275	ENDIF
3276	!
3277	!-- Check soil parameters. If soil_type is 0, all parameters
3278	!-- must be set.
3279	IF ( soil_type_f%from_file ) THEN
3280	check_passed = .TRUE.
3281	IF ( ALLOCATED( soil_type_f%var_2d ) ) THEN
3282	IF ( soil_type_f%var_2d(j,i) == 0 ) THEN
3283	IF ( ANY( soil_pars_f%pars_xy(:,j,i) == &
3284	soil_pars_f%fill ) ) check_passed = .FALSE.
3285	ENDIF
3286	ELSE
3287	IF ( ANY( soil_type_f%var_3d(:,j,i) == 0 ) ) THEN
3288	IF ( ANY( soil_pars_f%pars_xy(:,j,i) == &
3289	soil_pars_f%fill ) ) check_passed = .FALSE.
3290	ENDIF
3291	ENDIF
3292	IF ( .NOT. check_passed ) THEN
3293	message_string = 'If soil_type(y,x) = 0, all levels of ' //&
3294	'soil_pars at this location must be set.'
3295	CALL message( 'netcdf_data_input_mod', 'NDI033', &
3296	2, 2, 0, 6, 0 )
3297	ENDIF
3298	ENDIF
3299
3300	!
3301	!-- Check building parameters. If building_type is 0, all parameters
3302	!-- must be set.
3303	IF ( building_type_f%from_file ) THEN
3304	IF ( building_type_f%var(j,i) == 0 ) THEN
3305	IF ( ANY( building_pars_f%pars_xy(:,j,i) == &
3306	building_pars_f%fill ) ) THEN
3307	message_string = 'If building_type(y,x) = 0, all ' // &
3308	'parameters of building_pars at this '//&
3309	'location must be set.'
3310	CALL message( 'netcdf_data_input_mod', 'NDI034', &
3311	2, 2, 0, 6, 0 )
3312	ENDIF
3313	ENDIF
3314	ENDIF
3315	!
3316	!-- Check if building_type is set at each building and vice versa.
3317	IF ( building_type_f%from_file .AND. buildings_f%from_file ) THEN
3318	IF ( buildings_f%lod == 1 ) THEN
3319	IF ( buildings_f%var_2d(j,i) /= buildings_f%fill1 .AND. &
3320	building_type_f%var(j,i) == building_type_f%fill ) THEN
3321
3322	WRITE( message_string, * ) 'Each location where a ' // &
3323	'building is set requires a type ' // &
3324	'( and vice versa ) in case the ' // &
3325	'urban-surface model is applied. ' // &
3326	'i, j = ', i, j
3327	CALL message( 'netcdf_data_input_mod', 'NDI035', &
3328	2, 2, 0, 6, 0 )
3329	ENDIF
3330	ENDIF
3331	IF ( buildings_f%lod == 2 ) THEN
3332	IF ( ANY( buildings_f%var_3d(:,j,i) == 1 ) .AND. &
3333	building_type_f%var(j,i) == building_type_f%fill ) THEN
3334	WRITE( message_string, * ) 'Each location where a ' // &
3335	'building is set requires a type ' // &
3336	'( and vice versa ) in case the ' // &
3337	'urban-surface model is applied. ' // &
3338	'i, j = ', i, j
3339	CALL message( 'netcdf_data_input_mod', 'NDI035', &
3340	2, 2, 0, 6, 0 )
3341	ENDIF
3342	ENDIF
3343	ENDIF
3344	!
3345	!-- Check if at each location where a building is present also an ID
3346	!-- is set and vice versa.
3347	IF ( buildings_f%from_file ) THEN
3348	IF ( buildings_f%lod == 1 ) THEN
3349	IF ( buildings_f%var_2d(j,i) /= buildings_f%fill1 .AND. &
3350	building_id_f%var(j,i) == building_id_f%fill ) THEN
3351	WRITE( message_string, * ) 'Each location where a ' // &
3352	'building is set requires an ID ' // &
3353	'( and vice versa ). i, j = ', i, j
3354	CALL message( 'netcdf_data_input_mod', 'NDI036', &
3355	2, 2, 0, 6, 0 )
3356	ENDIF
3357	ELSEIF ( buildings_f%lod == 2 ) THEN
3358	IF ( ANY( buildings_f%var_3d(:,j,i) == 1 ) .AND. &
3359	building_id_f%var(j,i) == building_id_f%fill ) THEN
3360	WRITE( message_string, * ) 'Each location where a ' // &
3361	'building is set requires an ID ' // &
3362	'( and vice versa ). i, j = ', i, j
3363	CALL message( 'netcdf_data_input_mod', 'NDI036', &
3364	2, 2, 0, 6, 0 )
3365	ENDIF
3366	ENDIF
3367	ENDIF
3368	!
3369	!-- Check if at each location where a building ID or a -type is set
3370	!-- also a bulding is defined.
3371	IF ( buildings_f%from_file ) THEN
3372	IF ( buildings_f%lod == 1 ) THEN
3373	IF ( buildings_f%var_2d(j,i) /= buildings_f%fill1 .AND. &
3374	building_id_f%var(j,i) == building_id_f%fill ) THEN
3375	WRITE( message_string, * ) 'Each building grid point '// &
3376	'requires an ID.', i, j
3377	CALL message( 'netcdf_data_input_mod', 'NDI036', &
3378	2, 2, 0, 6, 0 )
3379	ENDIF
3380	ELSEIF ( buildings_f%lod == 2 ) THEN
3381	IF ( ANY( buildings_f%var_3d(:,j,i) == 1 ) &
3382	.AND. building_id_f%var(j,i) == building_id_f%fill ) THEN
3383	WRITE( message_string, * ) 'Each building grid point '// &
3384	'requires an ID.', i, j
3385	CALL message( 'netcdf_data_input_mod', 'NDI036', &
3386	2, 2, 0, 6, 0 )
3387	ENDIF
3388	ENDIF
3389	ENDIF
3390	!
3391	!-- Check albedo parameters. If albedo_type is 0, all parameters
3392	!-- must be set.
3393	IF ( albedo_type_f%from_file ) THEN
3394	IF ( albedo_type_f%var(j,i) == 0 ) THEN
3395	IF ( ANY( albedo_pars_f%pars_xy(:,j,i) == &
3396	albedo_pars_f%fill ) ) THEN
3397	message_string = 'If albedo_type(y,x) = 0, all ' // &
3398	'parameters of albedo_pars at this ' // &
3399	'location must be set.'
3400	CALL message( 'netcdf_data_input_mod', 'NDI037', &
3401	2, 2, 0, 6, 0 )
3402	ENDIF
3403	ENDIF
3404	ENDIF
3405
3406	!
3407	!-- Check pavement parameters. If pavement_type is 0, all parameters
3408	!-- of pavement_pars must be set at this location.
3409	IF ( pavement_type_f%from_file ) THEN
3410	IF ( pavement_type_f%var(j,i) == 0 ) THEN
3411	IF ( ANY( pavement_pars_f%pars_xy(:,j,i) == &
3412	pavement_pars_f%fill ) ) THEN
3413	message_string = 'If pavement_type(y,x) = 0, all ' // &
3414	'parameters of pavement_pars at this '//&
3415	'location must be set.'
3416	CALL message( 'netcdf_data_input_mod', 'NDI038', &
3417	2, 2, 0, 6, 0 )
3418	ENDIF
3419	ENDIF
3420	ENDIF
3421	!
3422	!-- Check pavement-subsurface parameters. If pavement_type is 0,
3423	!-- all parameters of pavement_subsurface_pars must be set at this
3424	!-- location.
3425	IF ( pavement_type_f%from_file ) THEN
3426	IF ( pavement_type_f%var(j,i) == 0 ) THEN
3427	IF ( ANY( pavement_subsurface_pars_f%pars_xyz(:,:,j,i) == &
3428	pavement_subsurface_pars_f%fill ) ) THEN
3429	message_string = 'If pavement_type(y,x) = 0, all ' // &
3430	'parameters of ' // &
3431	'pavement_subsurface_pars at this '// &
3432	'location must be set.'
3433	CALL message( 'netcdf_data_input_mod', 'NDI039', &
3434	2, 2, 0, 6, 0 )
3435	ENDIF
3436	ENDIF
3437	ENDIF
3438
3439	!
3440	!-- Check water parameters. If water_type is 0, all parameters
3441	!-- must be set at this location.
3442	IF ( water_type_f%from_file ) THEN
3443	IF ( water_type_f%var(j,i) == 0 ) THEN
3444	IF ( ANY( water_pars_f%pars_xy(:,j,i) == &
3445	water_pars_f%fill ) ) THEN
3446	message_string = 'If water_type(y,x) = 0, all ' // &
3447	'parameters of water_pars at this ' // &
3448	'location must be set.'
3449	CALL message( 'netcdf_data_input_mod', 'NDI040', &
3450	2, 2, 0, 6, 0 )
3451	ENDIF
3452	ENDIF
3453	ENDIF
3454
3455	ENDDO
3456	ENDDO
3457
3458	END SUBROUTINE netcdf_data_input_check_static
3459
3460	!------------------------------------------------------------------------------!
3461	! Description:
3462	! ------------
3463	!> Vertical interpolation and extrapolation of 1D variables.
3464	!------------------------------------------------------------------------------!
3465	SUBROUTINE netcdf_data_input_interpolate_1d( var, z_grid, z_file)
3466
3467	IMPLICIT NONE
3468
3469	LOGICAL :: top !< flag indicating extrapolation at model top
3470
3471	INTEGER(iwp) :: k !< running index z-direction file
3472	INTEGER(iwp) :: kk !< running index z-direction stretched model grid
3473	INTEGER(iwp) :: kl !< lower index bound along z-direction
3474	INTEGER(iwp) :: ku !< upper index bound along z-direction
3475	INTEGER(iwp) :: nz_file !< number of vertical levels on file
3476
3477
3478	REAL(wp), DIMENSION(:) :: z_grid !< grid levels on numeric grid
3479	REAL(wp), DIMENSION(:) :: z_file !< grid levels on file grid
3480	REAL(wp), DIMENSION(:), INTENT(INOUT) :: var !< treated variable
3481	REAL(wp), DIMENSION(:), ALLOCATABLE :: var_tmp !< temporary variable
3482
3483
3484	kl = LBOUND(var,1)
3485	ku = UBOUND(var,1)
3486	ALLOCATE( var_tmp(kl:ku) )
3487
3488	DO k = kl, ku
3489
3490	kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )
3491
3492	IF ( kk < ku ) THEN
3493	IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN
3494	var_tmp(k) = var(kk) + &
3495	( var(kk+1) - var(kk) ) / &
3496	( z_file(kk+1) - z_file(kk) ) * &
3497	( z_grid(k) - z_file(kk) )
3498
3499	ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN
3500	var_tmp(k) = var(kk-1) + &
3501	( var(kk) - var(kk-1) ) / &
3502	( z_file(kk) - z_file(kk-1) ) * &
3503	( z_grid(k) - z_file(kk-1) )
3504	ENDIF
3505	!
3506	!-- Extrapolate
3507	ELSE
3508
3509	var_tmp(k) = var(ku) + ( var(ku) - var(ku-1) ) / &
3510	( z_file(ku) - z_file(ku-1) ) * &
3511	( z_grid(k) - z_file(ku) )
3512
3513	ENDIF
3514
3515	ENDDO
3516	var(:) = var_tmp(:)
3517
3518	DEALLOCATE( var_tmp )
3519
3520
3521	END SUBROUTINE netcdf_data_input_interpolate_1d
3522
3523
3524	!------------------------------------------------------------------------------!
3525	! Description:
3526	! ------------
3527	!> Vertical interpolation and extrapolation of 1D variables from Inifor grid
3528	!> onto Palm grid, where both have same dimension. Please note, the passed
3529	!> paramter list in 1D version is different compared to 2D version.
3530	!------------------------------------------------------------------------------!
3531	SUBROUTINE netcdf_data_input_interpolate_1d_soil( var, var_file, &
3532	z_grid, z_file, &
3533	nzb_var, nzt_var, &
3534	nzb_file, nzt_file )
3535
3536	IMPLICIT NONE
3537
3538	INTEGER(iwp) :: i !< running index x-direction
3539	INTEGER(iwp) :: j !< running index y-direction
3540	INTEGER(iwp) :: k !< running index z-direction file
3541	INTEGER(iwp) :: kk !< running index z-direction stretched model grid
3542	INTEGER(iwp) :: ku !< upper index bound along z-direction for varialbe from file
3543	INTEGER(iwp) :: nzb_var !< lower bound of final array
3544	INTEGER(iwp) :: nzt_var !< upper bound of final array
3545	INTEGER(iwp) :: nzb_file !< lower bound of file array
3546	INTEGER(iwp) :: nzt_file !< upper bound of file array
3547
3548	! LOGICAL, OPTIONAL :: depth !< flag indicating reverse z-axis, i.e. depth instead of height, e.g. in case of ocean or soil
3549
3550	REAL(wp), DIMENSION(nzb_var:nzt_var) :: z_grid !< grid levels on numeric grid
3551	REAL(wp), DIMENSION(nzb_file:nzt_file) :: z_file !< grid levels on file grid
3552	REAL(wp), DIMENSION(nzb_var:nzt_var) :: var !< treated variable
3553	REAL(wp), DIMENSION(nzb_file:nzt_file) :: var_file !< temporary variable
3554
3555	ku = nzt_file
3556
3557	DO k = nzb_var, nzt_var
3558	!
3559	!-- Determine index on Inifor grid which is closest to the actual height
3560	kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )
3561	!
3562	!-- If closest index on Inifor grid is smaller than top index,
3563	!-- interpolate the data
3564	IF ( kk < nzt_file ) THEN
3565	IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN
3566	var(k) = var_file(kk) + ( var_file(kk+1) - var_file(kk) ) / &
3567	( z_file(kk+1) - z_file(kk) ) * &
3568	( z_grid(k) - z_file(kk) )
3569
3570	ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN
3571	var(k) = var_file(kk-1) + ( var_file(kk) - var_file(kk-1) ) / &
3572	( z_file(kk) - z_file(kk-1) ) * &
3573	( z_grid(k) - z_file(kk-1) )
3574	ENDIF
3575	!
3576	!-- Extrapolate if actual height is above the highest Inifor level
3577	ELSE
3578	var(k) = var_file(ku) + ( var_file(ku) - var_file(ku-1) ) / &
3579	( z_file(ku) - z_file(ku-1) ) * &
3580	( z_grid(k) - z_file(ku) )
3581
3582	ENDIF
3583
3584	ENDDO
3585
3586	END SUBROUTINE netcdf_data_input_interpolate_1d_soil
3587
3588	!------------------------------------------------------------------------------!
3589	! Description:
3590	! ------------
3591	!> Vertical interpolation and extrapolation of 2D variables at lateral boundaries.
3592	!------------------------------------------------------------------------------!
3593	SUBROUTINE netcdf_data_input_interpolate_2d( var, z_grid, z_file)
3594
3595	IMPLICIT NONE
3596
3597	LOGICAL :: top !< flag indicating extrapolation at model top
3598
3599	INTEGER(iwp) :: i !< running index x- or y -direction
3600	INTEGER(iwp) :: il !< lower index bound along x- or y-direction
3601	INTEGER(iwp) :: iu !< upper index bound along x- or y-direction
3602	INTEGER(iwp) :: k !< running index z-direction file
3603	INTEGER(iwp) :: kk !< running index z-direction stretched model grid
3604	INTEGER(iwp) :: kl !< lower index bound along z-direction
3605	INTEGER(iwp) :: ku !< upper index bound along z-direction
3606	INTEGER(iwp) :: nz_file !< number of vertical levels on file
3607
3608
3609	REAL(wp), DIMENSION(:) :: z_grid !< grid levels on numeric grid
3610	REAL(wp), DIMENSION(:) :: z_file !< grid levels on file grid
3611	REAL(wp), DIMENSION(:,:), INTENT(INOUT) :: var !< treated variable
3612	REAL(wp), DIMENSION(:), ALLOCATABLE :: var_tmp !< temporary variable
3613
3614
3615	il = LBOUND(var,2)
3616	iu = UBOUND(var,2)
3617	kl = LBOUND(var,1)
3618	ku = UBOUND(var,1)
3619	ALLOCATE( var_tmp(kl:ku) )
3620
3621	DO i = il, iu
3622	DO k = kl, ku
3623
3624	kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )
3625
3626	IF ( kk < ku ) THEN
3627	IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN
3628	var_tmp(k) = var(kk,i) + &
3629	( var(kk+1,i) - var(kk,i) ) / &
3630	( z_file(kk+1) - z_file(kk) ) * &
3631	( z_grid(k) - z_file(kk) )
3632
3633	ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN
3634	var_tmp(k) = var(kk-1,i) + &
3635	( var(kk,i) - var(kk-1,i) ) / &
3636	( z_file(kk) - z_file(kk-1) ) * &
3637	( z_grid(k) - z_file(kk-1) )
3638	ENDIF
3639	!
3640	!-- Extrapolate
3641	ELSE
3642
3643	var_tmp(k) = var(ku,i) + ( var(ku,i) - var(ku-1,i) ) / &
3644	( z_file(ku) - z_file(ku-1) ) * &
3645	( z_grid(k) - z_file(ku) )
3646
3647	ENDIF
3648
3649	ENDDO
3650	var(:,i) = var_tmp(:)
3651
3652	ENDDO
3653
3654	DEALLOCATE( var_tmp )
3655
3656
3657	END SUBROUTINE netcdf_data_input_interpolate_2d
3658
3659	!------------------------------------------------------------------------------!
3660	! Description:
3661	! ------------
3662	!> Vertical interpolation and extrapolation of 3D variables.
3663	!------------------------------------------------------------------------------!
3664	SUBROUTINE netcdf_data_input_interpolate_3d( var, z_grid, z_file )
3665
3666	IMPLICIT NONE
3667
3668	INTEGER(iwp) :: i !< running index x-direction
3669	INTEGER(iwp) :: il !< lower index bound along x-direction
3670	INTEGER(iwp) :: iu !< upper index bound along x-direction
3671	INTEGER(iwp) :: j !< running index y-direction
3672	INTEGER(iwp) :: jl !< lower index bound along x-direction
3673	INTEGER(iwp) :: ju !< upper index bound along x-direction
3674	INTEGER(iwp) :: k !< running index z-direction file
3675	INTEGER(iwp) :: kk !< running index z-direction stretched model grid
3676	INTEGER(iwp) :: kl !< lower index bound along z-direction
3677	INTEGER(iwp) :: ku !< upper index bound along z-direction
3678	INTEGER(iwp) :: nz_file !< number of vertical levels on file
3679
3680	REAL(wp), DIMENSION(:) :: z_grid !< grid levels on numeric grid
3681	REAL(wp), DIMENSION(:) :: z_file !< grid levels on file grid
3682	REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) :: var !< treated variable
3683	REAL(wp), DIMENSION(:), ALLOCATABLE :: var_tmp !< temporary variable
3684
3685	il = LBOUND(var,3)
3686	iu = UBOUND(var,3)
3687	jl = LBOUND(var,2)
3688	ju = UBOUND(var,2)
3689	kl = LBOUND(var,1)
3690	ku = UBOUND(var,1)
3691
3692	ALLOCATE( var_tmp(kl:ku) )
3693
3694	DO i = il, iu
3695	DO j = jl, ju
3696	DO k = kl, ku
3697
3698	kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )
3699
3700	IF ( kk < ku ) THEN
3701	IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN
3702	var_tmp(k) = var(kk,j,i) + &
3703	( var(kk+1,j,i) - var(kk,j,i) ) / &
3704	( z_file(kk+1) - z_file(kk) ) * &
3705	( z_grid(k) - z_file(kk) )
3706
3707	ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN
3708	var_tmp(k) = var(kk-1,j,i) + &
3709	( var(kk,j,i) - var(kk-1,j,i) ) / &
3710	( z_file(kk) - z_file(kk-1) ) * &
3711	( z_grid(k) - z_file(kk-1) )
3712	ENDIF
3713	!
3714	!-- Extrapolate
3715	ELSE
3716	var_tmp(k) = var(ku,j,i) + &
3717	( var(ku,j,i) - var(ku-1,j,i) ) / &
3718	( z_file(ku) - z_file(ku-1) ) * &
3719	( z_grid(k) - z_file(ku) )
3720
3721	ENDIF
3722	ENDDO
3723	var(:,j,i) = var_tmp(:)
3724	ENDDO
3725	ENDDO
3726
3727	DEALLOCATE( var_tmp )
3728
3729
3730	END SUBROUTINE netcdf_data_input_interpolate_3d
3731
3732	!------------------------------------------------------------------------------!
3733	! Description:
3734	! ------------
3735	!> Checks if a given variables is on file
3736	!------------------------------------------------------------------------------!
3737	FUNCTION check_existence( vars_in_file, var_name )
3738
3739	IMPLICIT NONE
3740
3741	CHARACTER(LEN=*) :: var_name !< variable to be checked
3742	CHARACTER(LEN=*), DIMENSION(:) :: vars_in_file !< list of variables in file
3743
3744	INTEGER(iwp) :: i !< loop variable
3745
3746	LOGICAL :: check_existence !< flag indicating whether a variable exist or not - actual return value
3747
3748	i = 1
3749	check_existence = .FALSE.
3750	DO WHILE ( i <= SIZE( vars_in_file ) )
3751	check_existence = TRIM( vars_in_file(i) ) == TRIM( var_name ) .OR. &
3752	check_existence
3753	i = i + 1
3754	ENDDO
3755
3756	RETURN
3757
3758	END FUNCTION check_existence
3759
3760
3761	!------------------------------------------------------------------------------!
3762	! Description:
3763	! ------------
3764	!> Closes an existing netCDF file.
3765	!------------------------------------------------------------------------------!
3766	SUBROUTINE close_input_file( id )
3767	#if defined( __netcdf )
3768
3769	USE pegrid
3770
3771	IMPLICIT NONE
3772
3773	INTEGER(iwp), INTENT(INOUT) :: id !< file id
3774
3775	nc_stat = NF90_CLOSE( id )
3776	CALL handle_error( 'close', 537 )
3777	#endif
3778	END SUBROUTINE close_input_file
3779
3780	!------------------------------------------------------------------------------!
3781	! Description:
3782	! ------------
3783	!> Opens an existing netCDF file for reading only and returns its id.
3784	!------------------------------------------------------------------------------!
3785	SUBROUTINE open_read_file( filename, id )
3786	#if defined( __netcdf )
3787
3788	USE pegrid
3789
3790	IMPLICIT NONE
3791
3792	CHARACTER (LEN=*), INTENT(IN) :: filename !< filename
3793	INTEGER(iwp), INTENT(INOUT) :: id !< file id
3794	LOGICAL :: file_open = .FALSE.
3795
3796	#if defined( __netcdf4_parallel )
3797	! if __netcdf4_parallel is defined, parrallel NetCDF will be used unconditionally
3798	nc_stat = NF90_OPEN( filename, IOR( NF90_WRITE, NF90_MPIIO ), id, &
3799	COMM = comm2d, INFO = MPI_INFO_NULL )
3800	IF(nc_stat /= NF90_NOERR ) THEN !possible NetCDF 3 file
3801	nc_stat = NF90_OPEN( filename, NF90_NOWRITE, id )
3802	collective_read = .FALSE.
3803	write(9,*) 'open ',TRIM(filename),' as serial NetCDF file on every PE'
3804	ELSE
3805	collective_read = .TRUE.
3806	write(9,*) 'open ',TRIM(filename),' using parallel NetCDF'
3807	END IF
3808	#else
3809	! All MPI processes open und read
3810	nc_stat = NF90_OPEN( filename, NF90_NOWRITE, id )
3811	#endif
3812
3813	CALL handle_error( 'open_read_file', 536 )
3814
3815	#endif
3816	END SUBROUTINE open_read_file
3817
3818	!------------------------------------------------------------------------------!
3819	! Description:
3820	! ------------
3821	!> Reads global or variable-related attributes of type INTEGER (32-bit)
3822	!------------------------------------------------------------------------------!
3823	SUBROUTINE get_attribute_int32( id, attribute_name, value, global, &
3824	variable_name )
3825
3826	USE pegrid
3827
3828	IMPLICIT NONE
3829
3830	CHARACTER(LEN=*) :: attribute_name !< attribute name
3831	CHARACTER(LEN=*), OPTIONAL :: variable_name !< variable name
3832
3833	INTEGER(iwp), INTENT(IN) :: id !< file id
3834	INTEGER(iwp) :: id_var !< variable id
3835	INTEGER(iwp), INTENT(INOUT) :: value !< read value
3836
3837	LOGICAL, INTENT(IN) :: global !< flag indicating global attribute
3838	#if defined( __netcdf )
3839
3840	!
3841	!-- Read global attribute
3842	IF ( global ) THEN
3843	nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3844	CALL handle_error( 'get_attribute_int32 global', 522 )
3845	!
3846	!-- Read attributes referring to a single variable. Therefore, first inquire
3847	!-- variable id
3848	ELSE
3849	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3850	CALL handle_error( 'get_attribute_int32', 522 )
3851	nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3852	CALL handle_error( 'get_attribute_int32', 522 )
3853	ENDIF
3854	#endif
3855	END SUBROUTINE get_attribute_int32
3856
3857	!------------------------------------------------------------------------------!
3858	! Description:
3859	! ------------
3860	!> Reads global or variable-related attributes of type INTEGER (8-bit)
3861	!------------------------------------------------------------------------------!
3862	SUBROUTINE get_attribute_int8( id, attribute_name, value, global, &
3863	variable_name )
3864
3865	USE pegrid
3866
3867	IMPLICIT NONE
3868
3869	CHARACTER(LEN=*) :: attribute_name !< attribute name
3870	CHARACTER(LEN=*), OPTIONAL :: variable_name !< variable name
3871
3872	INTEGER(iwp), INTENT(IN) :: id !< file id
3873	INTEGER(iwp) :: id_var !< variable id
3874	INTEGER(KIND=1), INTENT(INOUT) :: value !< read value
3875
3876	LOGICAL, INTENT(IN) :: global !< flag indicating global attribute
3877	#if defined( __netcdf )
3878
3879	!
3880	!-- Read global attribute
3881	IF ( global ) THEN
3882	nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3883	CALL handle_error( 'get_attribute_int8 global', 523 )
3884	!
3885	!-- Read attributes referring to a single variable. Therefore, first inquire
3886	!-- variable id
3887	ELSE
3888	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3889	CALL handle_error( 'get_attribute_int8', 523 )
3890	nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3891	CALL handle_error( 'get_attribute_int8', 523 )
3892	ENDIF
3893	#endif
3894	END SUBROUTINE get_attribute_int8
3895
3896	!------------------------------------------------------------------------------!
3897	! Description:
3898	! ------------
3899	!> Reads global or variable-related attributes of type REAL
3900	!------------------------------------------------------------------------------!
3901	SUBROUTINE get_attribute_real( id, attribute_name, value, global, &
3902	variable_name )
3903
3904	USE pegrid
3905
3906	IMPLICIT NONE
3907
3908	CHARACTER(LEN=*) :: attribute_name !< attribute name
3909	CHARACTER(LEN=*), OPTIONAL :: variable_name !< variable name
3910
3911	INTEGER(iwp), INTENT(IN) :: id !< file id
3912	INTEGER(iwp) :: id_var !< variable id
3913
3914	LOGICAL, INTENT(IN) :: global !< flag indicating global attribute
3915
3916	REAL(wp), INTENT(INOUT) :: value !< read value
3917	#if defined( __netcdf )
3918
3919
3920	!
3921	!-- Read global attribute
3922	IF ( global ) THEN
3923	nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3924	CALL handle_error( 'get_attribute_real global', 524 )
3925	!
3926	!-- Read attributes referring to a single variable. Therefore, first inquire
3927	!-- variable id
3928	ELSE
3929	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3930	CALL handle_error( 'get_attribute_real', 524 )
3931	nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3932	CALL handle_error( 'get_attribute_real', 524 )
3933	ENDIF
3934	#endif
3935	END SUBROUTINE get_attribute_real
3936
3937	!------------------------------------------------------------------------------!
3938	! Description:
3939	! ------------
3940	!> Reads global or variable-related attributes of type CHARACTER
3941	!> Remark: reading attributes of type NF_STRING return an error code 56 -
3942	!> Attempt to convert between text & numbers.
3943	!------------------------------------------------------------------------------!
3944	SUBROUTINE get_attribute_string( id, attribute_name, value, global, &
3945	variable_name )
3946
3947	USE pegrid
3948
3949	IMPLICIT NONE
3950
3951	CHARACTER(LEN=*) :: attribute_name !< attribute name
3952	CHARACTER(LEN=*), OPTIONAL :: variable_name !< variable name
3953	CHARACTER(LEN=*), INTENT(INOUT) :: value !< read value
3954
3955	INTEGER(iwp), INTENT(IN) :: id !< file id
3956	INTEGER(iwp) :: id_var !< variable id
3957
3958	LOGICAL, INTENT(IN) :: global !< flag indicating global attribute
3959	#if defined( __netcdf )
3960
3961	!
3962	!-- Read global attribute
3963	IF ( global ) THEN
3964	nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3965	CALL handle_error( 'get_attribute_string global', 525 )
3966	!
3967	!-- Read attributes referring to a single variable. Therefore, first inquire
3968	!-- variable id
3969	ELSE
3970	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3971	CALL handle_error( 'get_attribute_string', 525 )
3972
3973	nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3974	CALL handle_error( 'get_attribute_string',525 )
3975
3976	ENDIF
3977	#endif
3978	END SUBROUTINE get_attribute_string
3979
3980
3981
3982	!------------------------------------------------------------------------------!
3983	! Description:
3984	! ------------
3985	!> Get dimension array for a given dimension
3986	!------------------------------------------------------------------------------!
3987	SUBROUTINE get_dimension_length( id, dim_len, variable_name )
3988	#if defined( __netcdf )
3989
3990	USE pegrid
3991
3992	IMPLICIT NONE
3993
3994	CHARACTER(LEN=*) :: variable_name !< dimension name
3995	CHARACTER(LEN=100) :: dum !< dummy variable to receive return character
3996
3997	INTEGER(iwp) :: dim_len !< dimension size
3998	INTEGER(iwp), INTENT(IN) :: id !< file id
3999	INTEGER(iwp) :: id_dim !< dimension id
4000
4001	!
4002	!-- First, inquire dimension ID
4003	nc_stat = NF90_INQ_DIMID( id, TRIM( variable_name ), id_dim )
4004	CALL handle_error( 'get_dimension_length', 526 )
4005	!
4006	!-- Inquire dimension length
4007	nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim, dum, LEN = dim_len )
4008	CALL handle_error( 'get_dimension_length', 526 )
4009
4010	#endif
4011	END SUBROUTINE get_dimension_length
4012
4013	!------------------------------------------------------------------------------!
4014	! Description:
4015	! ------------
4016	!> Reads a 1D integer variable from file.
4017	!------------------------------------------------------------------------------!
4018	SUBROUTINE get_variable_1d_int( id, variable_name, var )
4019
4020	USE pegrid
4021
4022	IMPLICIT NONE
4023
4024	CHARACTER(LEN=*) :: variable_name !< variable name
4025
4026	INTEGER(iwp), INTENT(IN) :: id !< file id
4027	INTEGER(iwp) :: id_var !< dimension id
4028
4029	INTEGER(iwp), DIMENSION(:), INTENT(INOUT) :: var !< variable to be read
4030	#if defined( __netcdf )
4031
4032	!
4033	!-- First, inquire variable ID
4034	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4035	CALL handle_error( 'get_variable_1d_int', 527 )
4036	!
4037	!-- Inquire dimension length
4038	nc_stat = NF90_GET_VAR( id, id_var, var )
4039	CALL handle_error( 'get_variable_1d_int', 527 )
4040
4041	#endif
4042	END SUBROUTINE get_variable_1d_int
4043
4044	!------------------------------------------------------------------------------!
4045	! Description:
4046	! ------------
4047	!> Reads a 1D float variable from file.
4048	!------------------------------------------------------------------------------!
4049	SUBROUTINE get_variable_1d_real( id, variable_name, var )
4050
4051	USE pegrid
4052
4053	IMPLICIT NONE
4054
4055	CHARACTER(LEN=*) :: variable_name !< variable name
4056
4057	INTEGER(iwp), INTENT(IN) :: id !< file id
4058	INTEGER(iwp) :: id_var !< dimension id
4059
4060	REAL(wp), DIMENSION(:), INTENT(INOUT) :: var !< variable to be read
4061	#if defined( __netcdf )
4062
4063	!
4064	!-- First, inquire variable ID
4065	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4066	CALL handle_error( 'get_variable_1d_real', 527 )
4067	!
4068	!-- Inquire dimension length
4069	nc_stat = NF90_GET_VAR( id, id_var, var )
4070	CALL handle_error( 'get_variable_1d_real', 527 )
4071
4072	#endif
4073	END SUBROUTINE get_variable_1d_real
4074
4075
4076	!------------------------------------------------------------------------------!
4077	! Description:
4078	! ------------
4079	!> Reads a time-dependent 1D float variable from file.
4080	!------------------------------------------------------------------------------!
4081	SUBROUTINE get_variable_pr( id, variable_name, t, var )
4082	#if defined( __netcdf )
4083
4084	USE pegrid
4085
4086	IMPLICIT NONE
4087
4088	CHARACTER(LEN=*) :: variable_name !< variable name
4089
4090	INTEGER(iwp), INTENT(IN) :: id !< file id
4091	INTEGER(iwp), DIMENSION(1:2) :: id_dim !< dimension ids
4092	INTEGER(iwp) :: id_var !< dimension id
4093	INTEGER(iwp) :: n_file !< number of data-points in file along z dimension
4094	INTEGER(iwp), INTENT(IN) :: t !< timestep number
4095
4096	REAL(wp), DIMENSION(:), INTENT(INOUT) :: var !< variable to be read
4097
4098	!
4099	!-- First, inquire variable ID
4100	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4101	!
4102	!-- Inquire dimension size of vertical dimension
4103	nc_stat = NF90_INQUIRE_VARIABLE( id, id_var, DIMIDS = id_dim )
4104	nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim(1), LEN = n_file )
4105	!
4106	!-- Read variable.
4107	nc_stat = NF90_GET_VAR( id, id_var, var, &
4108	start = (/ 1, t /), &
4109	count = (/ n_file, 1 /) )
4110	CALL handle_error( 'get_variable_pr', 527 )
4111
4112	#endif
4113	END SUBROUTINE get_variable_pr
4114
4115
4116	!------------------------------------------------------------------------------!
4117	! Description:
4118	! ------------
4119	!> Reads a 2D REAL variable from a file. Reading is done processor-wise,
4120	!> i.e. each core reads its own domain in slices along x.
4121	!------------------------------------------------------------------------------!
4122	SUBROUTINE get_variable_2d_real( id, variable_name, i, var )
4123
4124	USE indices
4125	USE pegrid
4126
4127	IMPLICIT NONE
4128
4129	CHARACTER(LEN=*) :: variable_name !< variable name
4130
4131	INTEGER(iwp), INTENT(IN) :: i !< index along x direction
4132	INTEGER(iwp), INTENT(IN) :: id !< file id
4133	INTEGER(iwp) :: id_var !< variable id
4134
4135	REAL(wp), DIMENSION(nys:nyn), INTENT(INOUT) :: var !< variable to be read
4136	#if defined( __netcdf )
4137	!
4138	!-- Inquire variable id
4139	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4140	!
4141	!-- Get variable
4142	nc_stat = NF90_GET_VAR( id, id_var, var(nys:nyn), &
4143	start = (/ i+1, nys+1 /), &
4144	count = (/ 1, nyn - nys + 1 /) )
4145
4146	CALL handle_error( 'get_variable_2d_real', 528 )
4147	#endif
4148	END SUBROUTINE get_variable_2d_real
4149
4150	!------------------------------------------------------------------------------!
4151	! Description:
4152	! ------------
4153	!> Reads a 2D 32-bit INTEGER variable from file. Reading is done processor-wise,
4154	!> i.e. each core reads its own domain in slices along x.
4155	!------------------------------------------------------------------------------!
4156	SUBROUTINE get_variable_2d_int32( id, variable_name, i, var )
4157
4158	USE indices
4159	USE pegrid
4160
4161	IMPLICIT NONE
4162
4163	CHARACTER(LEN=*) :: variable_name !< variable name
4164
4165	INTEGER(iwp), INTENT(IN) :: i !< index along x direction
4166	INTEGER(iwp), INTENT(IN) :: id !< file id
4167	INTEGER(iwp) :: id_var !< variable id
4168	INTEGER(iwp), DIMENSION(nys:nyn), INTENT(INOUT) :: var !< variable to be read
4169	#if defined( __netcdf )
4170	!
4171	!-- Inquire variable id
4172	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4173	!
4174	!-- Get variable
4175	nc_stat = NF90_GET_VAR( id, id_var, var(nys:nyn), &
4176	start = (/ i+1, nys+1 /), &
4177	count = (/ 1, nyn - nys + 1 /) )
4178
4179	CALL handle_error( 'get_variable_2d_int32', 529 )
4180	#endif
4181	END SUBROUTINE get_variable_2d_int32
4182
4183	!------------------------------------------------------------------------------!
4184	! Description:
4185	! ------------
4186	!> Reads a 2D 8-bit INTEGER variable from file. Reading is done processor-wise,
4187	!> i.e. each core reads its own domain in slices along x.
4188	!------------------------------------------------------------------------------!
4189	SUBROUTINE get_variable_2d_int8( id, variable_name, i, var )
4190
4191	USE indices
4192	USE pegrid
4193
4194	IMPLICIT NONE
4195
4196	CHARACTER(LEN=*) :: variable_name !< variable name
4197
4198	INTEGER(iwp), INTENT(IN) :: i !< index along x direction
4199	INTEGER(iwp), INTENT(IN) :: id !< file id
4200	INTEGER(iwp) :: id_var !< variable id
4201	INTEGER(KIND=1), DIMENSION(nys:nyn), INTENT(INOUT) :: var !< variable to be read
4202	#if defined( __netcdf )
4203	!
4204	!-- Inquire variable id
4205	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4206	!
4207	!-- Get variable
4208	nc_stat = NF90_GET_VAR( id, id_var, var(nys:nyn), &
4209	start = (/ i+1, nys+1 /), &
4210	count = (/ 1, nyn - nys + 1 /) )
4211
4212	CALL handle_error( 'get_variable_2d_int8', 530 )
4213	#endif
4214	END SUBROUTINE get_variable_2d_int8
4215
4216	!------------------------------------------------------------------------------!
4217	! Description:
4218	! ------------
4219	!> Reads a 3D 8-bit INTEGER variable from file.
4220	!------------------------------------------------------------------------------!
4221	SUBROUTINE get_variable_3d_int8( id, variable_name, i, j, var )
4222
4223	USE indices
4224	USE pegrid
4225
4226	IMPLICIT NONE
4227
4228	CHARACTER(LEN=*) :: variable_name !< variable name
4229
4230	INTEGER(iwp), INTENT(IN) :: i !< index along x direction
4231	INTEGER(iwp), INTENT(IN) :: id !< file id
4232	INTEGER(iwp) :: id_var !< variable id
4233	INTEGER(iwp), INTENT(IN) :: j !< index along y direction
4234	INTEGER(iwp) :: n_file !< number of data-points along 3rd dimension
4235
4236	INTEGER(iwp), DIMENSION(1:3) :: id_dim
4237
4238	INTEGER( KIND = 1 ), DIMENSION(nzb:nzt+1), INTENT(INOUT) :: var !< variable to be read
4239	#if defined( __netcdf )
4240
4241	!
4242	!-- Inquire variable id
4243	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4244	!
4245	!-- Get length of first dimension, required for the count parameter.
4246	!-- Therefore, first inquired dimension ids
4247	nc_stat = NF90_INQUIRE_VARIABLE( id, id_var, DIMIDS = id_dim )
4248	nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim(3), LEN = n_file )
4249	!
4250	!-- Get variable
4251	nc_stat = NF90_GET_VAR( id, id_var, var, &
4252	start = (/ i+1, j+1, 1 /), &
4253	count = (/ 1, 1, n_file /) )
4254
4255	CALL handle_error( 'get_variable_3d_int8', 531 )
4256	#endif
4257	END SUBROUTINE get_variable_3d_int8
4258
4259
4260	!------------------------------------------------------------------------------!
4261	! Description:
4262	! ------------
4263	!> Reads a 3D float variable from file.
4264	!------------------------------------------------------------------------------!
4265	SUBROUTINE get_variable_3d_real( id, variable_name, i, j, var )
4266
4267	USE indices
4268	USE pegrid
4269
4270	IMPLICIT NONE
4271
4272	CHARACTER(LEN=*) :: variable_name !< variable name
4273
4274	INTEGER(iwp), INTENT(IN) :: i !< index along x direction
4275	INTEGER(iwp), INTENT(IN) :: id !< file id
4276	INTEGER(iwp) :: id_var !< variable id
4277	INTEGER(iwp), INTENT(IN) :: j !< index along y direction
4278	INTEGER(iwp) :: n3 !< number of data-points along 3rd dimension
4279
4280	INTEGER(iwp), DIMENSION(3) :: id_dim
4281
4282	REAL(wp), DIMENSION(:), INTENT(INOUT) :: var !< variable to be read
4283	#if defined( __netcdf )
4284
4285	!
4286	!-- Inquire variable id
4287	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4288	!
4289	!-- Get length of first dimension, required for the count parameter.
4290	!-- Therefore, first inquired dimension ids
4291	nc_stat = NF90_INQUIRE_VARIABLE( id, id_var, DIMIDS = id_dim )
4292	nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim(3), LEN = n3 )
4293	!
4294	!-- Get variable
4295	nc_stat = NF90_GET_VAR( id, id_var, var, &
4296	start = (/ i+1, j+1, 1 /), &
4297	count = (/ 1, 1, n3 /) )
4298
4299	CALL handle_error( 'get_variable_3d_real', 532 )
4300
4301	#endif
4302	END SUBROUTINE get_variable_3d_real
4303
4304	!------------------------------------------------------------------------------!
4305	! Description:
4306	! ------------
4307	!> Reads a 3D float array from file.
4308	!------------------------------------------------------------------------------!
4309	SUBROUTINE get_variable_3d_real_v( id, variable_name, is, ie, js, je, var )
4310
4311	USE indices
4312	USE pegrid
4313
4314	IMPLICIT NONE
4315
4316	CHARACTER(LEN=*) :: variable_name !< variable name
4317
4318	INTEGER(iwp), INTENT(IN) :: is,ie !< index range along x direction
4319	INTEGER(iwp), INTENT(IN) :: id !< file id
4320	INTEGER(iwp) :: id_var !< variable id
4321	INTEGER(iwp), INTENT(IN) :: js,je !< index range along y direction
4322	INTEGER(iwp) :: n3 !< number of data-points along 3rd dimension
4323
4324	INTEGER(iwp) :: i,j,k
4325	INTEGER(iwp), DIMENSION(3) :: id_dim
4326
4327	REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) :: var !< variable to be read
4328	#if defined( __netcdf )
4329	!
4330	!-- Inside the ...static NetCDF files, the array is stored as float.
4331	!-- Therefore single precision is sufficiant for the temporary array
4332
4333	REAL(sp), DIMENSION(:,:,:), ALLOCATABLE :: tmp_var !< temporary array to read NetCDF data in i,j,k direction
4334
4335	!kk Please check, if it is time consuming to do the inquire every time
4336	!
4337	!-- Inquire variable id
4338	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4339	!
4340	!-- Get length of first dimension, required for the count parameter.
4341	!-- Therefore, first inquired dimension ids
4342	nc_stat = NF90_INQUIRE_VARIABLE( id, id_var, DIMIDS = id_dim )
4343	nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim(3), LEN = n3 )
4344
4345	IF(collective_read) THEN
4346	nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4347	ENDIF
4348
4349	!
4350	!-- Allocate temporary array ro read NetCDF data in i,j,k direction
4351
4352	ALLOCATE(tmp_var(is:ie,js:je,n3))
4353	!
4354	!-- Get variable
4355	!-- Read complete local 3-D array in oone call
4356
4357	nc_stat = NF90_GET_VAR( id, id_var, tmp_var, &
4358	start = (/ is+1, js+1, 1 /), &
4359	count = (/ ie-is+1, je-js+1, n3 /) )
4360
4361	CALL handle_error( 'get_variable_3d_real', 532 )
4362
4363	!
4364	!-- Resort data in k,j,i direction
4365
4366	DO i=is,ie
4367	DO j=js,je
4368	DO K=1,n3
4369	var (k,j-js+1,i-is+1) = tmp_var(i,j,k)
4370	END DO
4371	END DO
4372	END DO
4373
4374	DEALLOCATE(tmp_var)
4375
4376	#endif
4377	END SUBROUTINE get_variable_3d_real_v
4378
4379
4380	!------------------------------------------------------------------------------!
4381	! Description:
4382	! ------------
4383	!> Reads a 4D float variable from file. Note, in constrast to 3D versions,
4384	!> dimensions are already inquired and passed so that they are known here.
4385	!------------------------------------------------------------------------------!
4386	SUBROUTINE get_variable_4d_real( id, variable_name, i, j, var, n3, n4 )
4387
4388	USE indices
4389	USE pegrid
4390
4391	IMPLICIT NONE
4392
4393	CHARACTER(LEN=*) :: variable_name !< variable name
4394
4395	INTEGER(iwp), INTENT(IN) :: i !< index along x direction
4396	INTEGER(iwp), INTENT(IN) :: id !< file id
4397	INTEGER(iwp) :: id_var !< variable id
4398	INTEGER(iwp), INTENT(IN) :: j !< index along y direction
4399	INTEGER(iwp), INTENT(IN) :: n3 !< number of data-points along 3rd dimension
4400	INTEGER(iwp), INTENT(IN) :: n4 !< number of data-points along 4th dimension
4401
4402	INTEGER(iwp), DIMENSION(3) :: id_dim
4403
4404	REAL(wp), DIMENSION(:,:), INTENT(INOUT) :: var !< variable to be read
4405	#if defined( __netcdf )
4406
4407	!
4408	!-- Inquire variable id
4409	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4410	!
4411	!-- Get variable
4412	nc_stat = NF90_GET_VAR( id, id_var, var, &
4413	start = (/ i+1, j+1, 1, 1 /), &
4414	count = (/ 1, 1, n3, n4 /) )
4415
4416	CALL handle_error( 'get_variable_4d_real', 533 )
4417	#endif
4418	END SUBROUTINE get_variable_4d_real
4419
4420
4421
4422	!------------------------------------------------------------------------------!
4423	! Description:
4424	! ------------
4425	!> Reads a 3D float variable at left, right, north, south and top boundaries.
4426	!------------------------------------------------------------------------------!
4427	SUBROUTINE get_variable_bc( id, variable_name, t_start, &
4428	i2_s, count_2, i3_s, count_3, var )
4429
4430	USE indices
4431	USE pegrid
4432
4433	IMPLICIT NONE
4434
4435	CHARACTER(LEN=*) :: variable_name !< variable name
4436
4437	INTEGER(iwp) :: count_2 !< number of elements in second dimension
4438	INTEGER(iwp) :: count_3 !< number of elements in third dimension (usually 1)
4439	INTEGER(iwp) :: i2_s !< start index of second dimension
4440	INTEGER(iwp) :: i3_s !< start index of third dimension
4441	INTEGER(iwp), INTENT(IN) :: id !< file id
4442	INTEGER(iwp) :: id_var !< variable id
4443	INTEGER(iwp) :: t_start !< start index at time dimension with respect to netcdf convention
4444
4445	REAL(wp), DIMENSION(:), INTENT(INOUT) :: var !< input variable
4446	#if defined( __netcdf )
4447
4448	!
4449	!-- Inquire variable id
4450	nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4451	!
4452	!-- Get variable
4453	nc_stat = NF90_GET_VAR( id, id_var, var, &
4454	start = (/ i3_s, i2_s, t_start /), &
4455	count = (/ count_3, count_2, 1 /) )
4456
4457	CALL handle_error( 'get_variable_bc', 532 )
4458	#endif
4459	END SUBROUTINE get_variable_bc
4460
4461
4462
4463	!------------------------------------------------------------------------------!
4464	! Description:
4465	! ------------
4466	!> Inquires the number of variables in a file
4467	!------------------------------------------------------------------------------!
4468	SUBROUTINE inquire_num_variables( id, num_vars )
4469
4470	USE indices
4471	USE pegrid
4472
4473	IMPLICIT NONE
4474
4475	INTEGER(iwp), INTENT(IN) :: id !< file id
4476	INTEGER(iwp), INTENT(INOUT) :: num_vars !< number of variables in a file
4477	#if defined( __netcdf )
4478
4479	nc_stat = NF90_INQUIRE( id, NVARIABLES = num_vars )
4480	CALL handle_error( 'inquire_num_variables', 534 )
4481
4482	#endif
4483	END SUBROUTINE inquire_num_variables
4484
4485
4486	!------------------------------------------------------------------------------!
4487	! Description:
4488	! ------------
4489	!> Inquires the variable names belonging to a file.
4490	!------------------------------------------------------------------------------!
4491	SUBROUTINE inquire_variable_names( id, var_names )
4492
4493	USE indices
4494	USE pegrid
4495
4496	IMPLICIT NONE
4497
4498	CHARACTER(LEN=*), DIMENSION(:), INTENT(INOUT) :: var_names !< return variable - variable names
4499	INTEGER(iwp) :: i !< loop variable
4500	INTEGER(iwp), INTENT(IN) :: id !< file id
4501	INTEGER(iwp) :: num_vars !< number of variables (unused return parameter)
4502	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: varids !< dummy array to strore variable ids temporarily
4503	#if defined( __netcdf )
4504
4505	ALLOCATE( varids(1:SIZE(var_names)) )
4506	nc_stat = NF90_INQ_VARIDS( id, NVARS = num_vars, VARIDS = varids )
4507	CALL handle_error( 'inquire_variable_names', 535 )
4508
4509	DO i = 1, SIZE(var_names)
4510	nc_stat = NF90_INQUIRE_VARIABLE( id, varids(i), NAME = var_names(i) )
4511	CALL handle_error( 'inquire_variable_names', 535 )
4512	ENDDO
4513
4514	DEALLOCATE( varids )
4515	#endif
4516	END SUBROUTINE inquire_variable_names
4517
4518	!------------------------------------------------------------------------------!
4519	! Description:
4520	! ------------
4521	!> Prints out a text message corresponding to the current status.
4522	!------------------------------------------------------------------------------!
4523	SUBROUTINE handle_error( routine_name, errno )
4524
4525	USE control_parameters, &
4526	ONLY: message_string
4527
4528	IMPLICIT NONE
4529
4530	CHARACTER(LEN=6) :: message_identifier
4531	CHARACTER(LEN=*) :: routine_name
4532
4533	INTEGER(iwp) :: errno
4534	#if defined( __netcdf )
4535
4536	IF ( nc_stat /= NF90_NOERR ) THEN
4537
4538	WRITE( message_identifier, '(''NC'',I4.4)' ) errno
4539	message_string = TRIM( NF90_STRERROR( nc_stat ) )
4540
4541	CALL message( routine_name, message_identifier, 2, 2, 0, 6, 1 )
4542
4543	ENDIF
4544
4545	#endif
4546	END SUBROUTINE handle_error
4547
4548
4549	END MODULE netcdf_data_input_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |