Home

Context Navigation

netcdf_data_input_mod.f90 @ 2958

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