Home

Context Navigation

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

Last change on this file since 3040 was 3037, checked in by gronemeier, 6 years ago
renamed input dimension "depth" to "zsoil"
Property svn:keywords set to `Id`
File size: 210.3 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |