Home

Context Navigation

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

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