source: palm/trunk/SOURCE/netcdf_data_input_mod.f90

Last change on this file was 4880, checked in by suehring, 7 months ago

Minor formatting adjustments and some comments added in get_variable_surf

  • Property svn:keywords set to Id
File size: 251.4 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 terms of the GNU General
6! Public License as published by the Free Software Foundation, either version 3 of the License, or
7! (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
10! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11! Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with PALM. If not, see
14! <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2021 Leibniz Universitaet Hannover
17!--------------------------------------------------------------------------------------------------!
18!
19! Current revisions:
20! -----------------
21!
22!
23! Former revisions:
24! -----------------
25! $Id: netcdf_data_input_mod.f90 4880 2021-02-18 12:08:41Z banzhafs $
26! Minor formatting adjustments and some comments added in get_variable_surf
27!
28! 4878 2021-02-18 09:47:49Z suehring
29! Rename resize_array into add_ghost_layers; remove number of passed indices; replace subroutine
30! calls with interface name
31!
32! 4877 2021-02-17 16:17:35Z suehring
33! Add interface for resize_array and add subroutine to resize 2d-real arrays
34!
35! 4828 2021-01-05 11:21:41Z Giersch
36! Check if netCDF file actually exists before opening it for reading.
37!
38! 4807 2020-12-02 21:02:28Z gronemeier
39! Deactivated reading of building_obstruction_full due to improper implementation (conflicts with
40! building_obstruction_f and it is also not used anywhere else in the code).
41!
42! 4806 2020-12-02 21:00:32Z gronemeier
43! file re-formatted to follow the PALM coding standard
44!
45! 4724 2020-10-06 17:20:39Z suehring
46! - New routines to read LOD=1 variables from dynamic input file
47! - add no_abort option to all get_attribute routines
48!
49! 4641 2020-08-13 09:57:07Z suehring
50! To follow (UC)2 standard, change default of attribute data_content
51!
52! 4507 2020-04-22 18:21:45Z gronemeier
53! - bugfix: check terrain height for fill values directly after reading
54! - changes:
55!   - remove check for negative zt
56!   - add reference height from input file upon PALM reference height (origin_z)
57!
58! 4457 2020-03-11 14:20:43Z raasch
59! use statement for exchange horiz added,
60! bugfixes for calls of exchange horiz 2d
61!
62! 4435 2020-03-03 10:38:41Z raasch
63! temporary bugfix to avoid compile problems with older NetCDFD libraries on IMUK machines
64!
65! 4434 2020-03-03 10:02:18Z oliver.maas
66! added optional netcdf data input for wtm array input parameters
67!
68! 4404 2020-02-12 17:01:53Z suehring
69! Fix misplaced preprocessor directives.
70!
71! 4401 2020-02-11 16:19:09Z suehring
72! Define a default list of coordinate reference system variables used when no static driver input is
73! available
74!
75! 4400 2020-02-10 20:32:41Z suehring
76! - Routine to inquire default fill values added
77! - netcdf_data_input_att and netcdf_data_input_var routines removed
78!
79! 4392 2020-01-31 16:14:57Z pavelkrc
80! (resler) Decrease length of reading buffer (fix problem of ifort/icc compilers)
81!
82! 4389 2020-01-29 08:22:42Z raasch
83! Error messages refined for reading ASCII topo file, also reading of topo file revised so that
84! statement labels and goto statements are not required any more
85!
86! 4388 2020-01-28 16:36:55Z raasch
87! bugfix for error messages while reading ASCII topo file
88!
89! 4387 2020-01-28 11:44:20Z banzhafs
90! Added subroutine get_variable_string_generic ( ) and added to interface get_variable to circumvent
91! unknown application-specific restrictions in existing function get_variable_string ( ), which is
92! retained for backward compatibility (ECC)
93!
94! 4370 2020-01-10 14:00:44Z raasch
95! collective read switched off on NEC Aurora to avoid hang situations
96!
97! 4362 2020-01-07 17:15:02Z suehring
98! Input of plant canopy variables from static driver moved to plant-canopy model
99!
100! 4360 2020-01-07 11:25:50Z suehring
101! Correct single message calls, local checks must be given by the respective mpi rank.
102!
103! 4346 2019-12-18 11:55:56Z motisi
104! Introduction of wall_flags_total_0, which currently sets bits based on static topography
105! information used in wall_flags_static_0
106!
107! 4329 2019-12-10 15:46:36Z motisi
108! Renamed wall_flags_0 to wall_flags_static_0
109!
110! 4321 2019-12-04 10:26:38Z pavelkrc
111! Further revise check for surface fractions
112!
113! 4313 2019-11-27 14:07:00Z suehring
114! Checks for surface fractions revised
115!
116! 4312 2019-11-27 14:06:25Z suehring
117! Open input files with read-only attribute instead of write attribute.
118!
119! 4280 2019-10-29 14:34:15Z monakurppa
120! Remove id_emis flags from get_variable_4d_to_3d_real and get_variable_5d_to_4d_real
121!
122! 4258 2019-10-07 13:29:08Z suehring
123! - Migrate input of soil temperature and moisture to land-surface model.
124! - Remove interpolate routines and move the only required subroutine to land-surface model.
125!
126! 4247 2019-09-30 10:18:24Z pavelkrc
127! Add reading and processing of building_surface_pars
128!
129! 4226 2019-09-10 17:03:24Z suehring
130! - Netcdf input routine for dimension length renamed
131! - Move offline-nesting-specific checks to nesting_offl_mod
132! - Module-specific input of boundary data for offline nesting moved to nesting_offl_mod
133! - Define module specific data type for offline nesting in nesting_offl_mod
134!
135! 4190 2019-08-27 15:42:37Z suehring
136! type real_1d changed to real_1d_3d
137!
138! 4186 2019-08-23 16:06:14Z suehring
139! Minor formatting adjustments
140!
141! 4182 2019-08-22 15:20:23Z scharf
142! Corrected "Former revisions" section
143!
144! 4178 2019-08-21 11:13:06Z suehring
145! Implement input of external radiation forcing. Therefore, provide public subroutines and
146! variables.
147!
148! 4150 2019-08-08 20:00:47Z suehring
149! Some variables are given the public attribute, in order to call netcdf input from single routines
150!
151! 4125 2019-07-29 13:31:44Z suehring
152! To enable netcdf-parallel access for lateral boundary data (dynamic input), zero number of
153! elements are passed to the respective get_variable routine for non-boundary cores.
154!
155! 4100 2019-07-17 08:11:29Z forkel
156! Made check for input_pids_dynamic and 'inifor' more general
157!
158! 4012 2019-05-31 15:19:05Z monakurppa
159!
160! 3994 2019-05-22 18:08:09Z suehring
161! Remove single location message
162!
163! 3976 2019-05-15 11:02:34Z hellstea
164! Remove unused variables from last commit
165!
166! 3969 2019-05-13 12:14:33Z suehring
167! - clean-up index notations for emission_values to eliminate magic numbers
168! - introduce temporary variable dum_var_5d as well as subroutines get_var_5d_real and
169!   get_var_5d_real_dynamic
170! - remove emission-specific code in generic get_variable routines
171! - in subroutine netcdf_data_input_chemistry_data change netCDF LOD 1 (default) emission_values to
172!   the following index order: z, y, x, species, category
173! - in subroutine netcdf_data_input_chemistry_data changed netCDF LOD 2 pre-processed
174!   emission_values to the following index order: time, z, y, x, species
175! - in type chem_emis_att_type replace nspec with n_emiss_species but retained nspec for backward
176!   compatibility with salsa_mod. (E.C. Chan)
177!
178! 3961 2019-05-08 16:12:31Z suehring
179! Revise checks for building IDs and types
180!
181! 3943 2019-05-02 09:50:41Z maronga
182! Temporarily disabled some (faulty) checks for static driver.
183!
184! 3942 2019-04-30 13:08:30Z kanani
185! Fix: increase LEN of all NetCDF attribute values (caused crash in netcdf_create_global_atts due to
186!      insufficient length)
187!
188! 3941 2019-04-30 09:48:33Z suehring
189! Move check for grid dimension to an earlier point in time when first array is read.
190! Improve checks for building types / IDs with respect to 2D/3D buildings.
191!
192! 3885 2019-04-11 11:29:34Z kanani
193! Changes related to global restructuring of location messages and introduction of additional debug
194! messages
195!
196! 3864 2019-04-05 09:01:56Z monakurppa
197! get_variable_4d_to_3d_real modified to enable read in data of type data(t,y,x,n) one timestep at
198! a time + some routines made public
199!
200! 3855 2019-04-03 10:00:59Z suehring
201! Typo removed
202!
203! 3854 2019-04-02 16:59:33Z suehring
204! Bugfix in one of the checks. Typo removed.
205!
206! 3744 2019-02-15 18:38:58Z suehring
207! Enable mesoscale offline nesting for chemistry variables as well as initialization of chemistry
208! via dynamic input file.
209!
210! 3705 2019-01-29 19:56:39Z suehring
211! Interface for attribute input of 8-bit and 32-bit integer
212!
213! 3704 2019-01-29 19:51:41Z suehring
214! unused variables removed
215!
216! 2696 2017-12-14 17:12:51Z kanani
217! Initial revision (suehring)
218!
219! Authors:
220! --------
221! @author Matthias Suehring
222! @author Edward C. Chan
223! @author Emanuele Russo
224!
225! Description:
226! ------------
227!> Modulue contains routines to input data according to Palm input data
228!> standart using dynamic and static input files.
229!> @todo - Chemistry: revise reading of netcdf file and ajdust formatting according to standard!!!
230!>         (ecc/done)
231!> @todo - Order input alphabetically
232!> @todo - Revise error messages and error numbers
233!> @todo - Input of missing quantities (chemical species, emission rates)
234!> @todo - Definition and input of still missing variable attributes (ecc/what are they?)
235!> @todo - Input of initial geostrophic wind profiles with cyclic conditions.
236!> @todo - remove z dimension from default_emission_data nad preproc_emission_data and correpsonding
237!>         subroutines get_var_5d_real and get_var_5d_dynamic (ecc)
238!> @todo - decpreciate chem_emis_att_type@nspec (ecc)
239!> @todo - depreciate subroutines get_variable_4d_to_3d_real and get_variable_5d_to_4d_real (ecc)
240!> @todo - introduce useful debug_message(s)
241!--------------------------------------------------------------------------------------------------!
242 MODULE netcdf_data_input_mod
243
244    USE control_parameters,                                                                        &
245        ONLY:  coupling_char, io_blocks, io_group
246
247    USE cpulog,                                                                                    &
248        ONLY:  cpu_log, log_point_s
249
250    USE indices,                                                                                   &
251        ONLY:  nbgp,                                                                               &
252               nxl,                                                                                &
253               nxr,                                                                                &
254               nys,                                                                                &
255               nyn
256
257    USE kinds
258
259#if defined ( __netcdf )
260    USE NETCDF
261#endif
262
263    USE pegrid
264
265    USE surface_mod,                                                                               &
266        ONLY:  ind_pav_green, ind_veg_wall, ind_wat_win
267!
268!-- Define type for dimensions.
269    TYPE dims_xy
270       INTEGER(iwp) :: nx                             !< dimension length in x
271       INTEGER(iwp) :: ny                             !< dimension length in y
272       INTEGER(iwp) :: nz                             !< dimension length in z
273       REAL(wp), DIMENSION(:), ALLOCATABLE :: x       !< dimension array in x
274       REAL(wp), DIMENSION(:), ALLOCATABLE :: y       !< dimension array in y
275       REAL(wp), DIMENSION(:), ALLOCATABLE :: z       !< dimension array in z
276    END TYPE dims_xy
277
278    TYPE init_type
279
280       CHARACTER(LEN=16) ::  init_char = 'init_atmosphere_'           !< leading substring for init variables
281       CHARACTER(LEN=23) ::  origin_time = '2000-01-01 00:00:00 +00'  !< reference time of input data
282
283       CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE ::  var_names_chem  !< list of chemistry variable names that can potentially be
284                                                                         !< on file
285
286       INTEGER(iwp) ::  lod_msoil !< level of detail - soil moisture
287       INTEGER(iwp) ::  lod_pt    !< level of detail - pt
288       INTEGER(iwp) ::  lod_q     !< level of detail - q
289       INTEGER(iwp) ::  lod_tsoil !< level of detail - soil temperature
290       INTEGER(iwp) ::  lod_u     !< level of detail - u-component
291       INTEGER(iwp) ::  lod_v     !< level of detail - v-component
292       INTEGER(iwp) ::  lod_w     !< level of detail - w-component
293       INTEGER(iwp) ::  nx        !< number of scalar grid points along x in dynamic input file
294       INTEGER(iwp) ::  nxu       !< number of u grid points along x in dynamic input file
295       INTEGER(iwp) ::  ny        !< number of scalar grid points along y in dynamic input file
296       INTEGER(iwp) ::  nyv       !< number of v grid points along y in dynamic input file
297       INTEGER(iwp) ::  nzs       !< number of vertical soil levels in dynamic input file
298       INTEGER(iwp) ::  nzu       !< number of vertical levels on scalar grid in dynamic input file
299       INTEGER(iwp) ::  nzw       !< number of vertical levels on w grid in dynamic input file
300
301       INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  lod_chem !< level of detail - chemistry variables
302
303       LOGICAL ::  from_file_msoil  = .FALSE. !< flag indicating whether soil moisture is already initialized from file
304       LOGICAL ::  from_file_pt     = .FALSE. !< flag indicating whether pt is already initialized from file
305       LOGICAL ::  from_file_q      = .FALSE. !< flag indicating whether q is already initialized from file
306       LOGICAL ::  from_file_tsoil  = .FALSE. !< flag indicating whether soil temperature is already initialized from file
307       LOGICAL ::  from_file_u      = .FALSE. !< flag indicating whether u is already initialized from file
308       LOGICAL ::  from_file_ug     = .FALSE. !< flag indicating whether ug is already initialized from file
309       LOGICAL ::  from_file_v      = .FALSE. !< flag indicating whether v is already initialized from file
310       LOGICAL ::  from_file_vg     = .FALSE. !< flag indicating whether ug is already initialized from file
311       LOGICAL ::  from_file_w      = .FALSE. !< flag indicating whether w is already initialized from file
312
313       LOGICAL, DIMENSION(:), ALLOCATABLE ::  from_file_chem !< flag indicating whether chemistry variable is read from file
314
315       REAL(wp) ::  fill_msoil              !< fill value for soil moisture
316       REAL(wp) ::  fill_pt                 !< fill value for pt
317       REAL(wp) ::  fill_q                  !< fill value for q
318       REAL(wp) ::  fill_tsoil              !< fill value for soil temperature
319       REAL(wp) ::  fill_u                  !< fill value for u
320       REAL(wp) ::  fill_v                  !< fill value for v
321       REAL(wp) ::  fill_w                  !< fill value for w
322       REAL(wp) ::  latitude = 0.0_wp       !< latitude of the lower left corner
323       REAL(wp) ::  longitude = 0.0_wp      !< longitude of the lower left corner
324       REAL(wp) ::  origin_x = 500000.0_wp  !< UTM easting of the lower left corner
325       REAL(wp) ::  origin_y = 0.0_wp       !< UTM northing of the lower left corner
326       REAL(wp) ::  origin_z = 0.0_wp       !< reference height of input data
327       REAL(wp) ::  rotation_angle = 0.0_wp !< rotation angle of input data
328
329       REAL(wp), DIMENSION(:), ALLOCATABLE ::  fill_chem    !< fill value - chemistry variables
330       REAL(wp), DIMENSION(:), ALLOCATABLE ::  msoil_1d     !< initial vertical profile of soil moisture
331       REAL(wp), DIMENSION(:), ALLOCATABLE ::  pt_init      !< initial vertical profile of pt
332       REAL(wp), DIMENSION(:), ALLOCATABLE ::  q_init       !< initial vertical profile of q
333       REAL(wp), DIMENSION(:), ALLOCATABLE ::  tsoil_1d     !< initial vertical profile of soil temperature
334       REAL(wp), DIMENSION(:), ALLOCATABLE ::  u_init       !< initial vertical profile of u
335       REAL(wp), DIMENSION(:), ALLOCATABLE ::  ug_init      !< initial vertical profile of ug
336       REAL(wp), DIMENSION(:), ALLOCATABLE ::  v_init       !< initial vertical profile of v
337       REAL(wp), DIMENSION(:), ALLOCATABLE ::  vg_init      !< initial vertical profile of ug
338       REAL(wp), DIMENSION(:), ALLOCATABLE ::  w_init       !< initial vertical profile of w
339       REAL(wp), DIMENSION(:), ALLOCATABLE ::  z_soil       !< vertical levels in soil in dynamic input file, used for interpolation
340       REAL(wp), DIMENSION(:), ALLOCATABLE ::  zu_atmos     !< vertical levels at scalar grid in dynamic input file, used for
341                                                            !< interpolation
342       REAL(wp), DIMENSION(:), ALLOCATABLE ::  zw_atmos     !< vertical levels at w grid in dynamic input file, used for
343                                                            !< interpolation
344
345       REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  chem_init  !< initial vertical profiles of chemistry variables
346
347
348       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  msoil_3d  !< initial 3d soil moisture provide by Inifor and interpolated onto
349                                                             !< soil grid
350       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  tsoil_3d  !< initial 3d soil temperature provide by Inifor and interpolated onto
351                                                             !< soil grid
352
353    END TYPE init_type
354!
355!-- Data type for the general information of chemistry emissions, do not dependent on the particular chemical species
356    TYPE chem_emis_att_type
357!
358!--    DIMENSIONS
359       INTEGER(iwp)                                 :: nspec=0            !< no of chem species provided in emission_values
360       INTEGER(iwp)                                 :: n_emiss_species=0  !< no of chem species provided in emission_values
361                                                                          !< same function as nspec, which will be depreciated (ecc)
362       INTEGER(iwp)                                 :: ncat=0             !< number of emission categories
363       INTEGER(iwp)                                 :: nvoc=0             !< number of VOC components
364       INTEGER(iwp)                                 :: npm=0              !< number of PM components
365       INTEGER(iwp)                                 :: nnox=2             !< number of NOx components: NO and NO2
366       INTEGER(iwp)                                 :: nsox=2             !< number of SOX components: SO and SO4
367       INTEGER(iwp)                                 :: nhoursyear         !< number of hours of a specific year in the HOURLY mode
368                                                                          !< of the default mode
369       INTEGER(iwp)                                 :: nmonthdayhour      !< number of month days and hours in the MDH mode
370                                                                          !< of the default mode
371       INTEGER(iwp)                                 :: dt_emission        !< Number of emissions timesteps for one year
372                                                                          !< in the pre-processed emissions case
373!
374!--    1d emission input variables
375       CHARACTER (LEN=25),ALLOCATABLE, DIMENSION(:) :: pm_name       !< Names of PM components
376       CHARACTER (LEN=25),ALLOCATABLE, DIMENSION(:) :: cat_name      !< Emission category names
377       CHARACTER (LEN=25),ALLOCATABLE, DIMENSION(:) :: species_name  !< Names of emission chemical species
378       CHARACTER (LEN=25),ALLOCATABLE, DIMENSION(:) :: voc_name      !< Names of VOCs components
379       CHARACTER (LEN=25)                           :: units         !< Units
380
381       INTEGER(iwp)                                 :: i_hour         !< indices for assigning emission values at different
382                                                                      !< timesteps
383       INTEGER(iwp),ALLOCATABLE, DIMENSION(:)       :: cat_index      !< Indices for emission categories
384       INTEGER(iwp),ALLOCATABLE, DIMENSION(:)       :: species_index  !< Indices for emission chem species
385
386       REAL(wp),ALLOCATABLE, DIMENSION(:)           :: xm             !< Molecular masses of emission chem species
387
388!--    2d emission input variables
389       REAL(wp),ALLOCATABLE, DIMENSION(:,:)         :: hourly_emis_time_factor  !< Time factors for HOURLY emissions (DEFAULT mode)
390       REAL(wp),ALLOCATABLE, DIMENSION(:,:)         :: mdh_emis_time_factor     !< Time factors for MDH emissions (DEFAULT mode)
391       REAL(wp),ALLOCATABLE, DIMENSION(:,:)         :: nox_comp                 !< Composition of NO and NO2
392       REAL(wp),ALLOCATABLE, DIMENSION(:,:)         :: sox_comp                 !< Composition of SO2 and SO4
393       REAL(wp),ALLOCATABLE, DIMENSION(:,:)         :: voc_comp                 !< Composition of VOC components (not fixed)
394
395!--    3d emission input variables
396       REAL(wp),ALLOCATABLE, DIMENSION(:,:,:)       :: pm_comp                  !< Composition of PM components (not fixed)
397
398    END TYPE chem_emis_att_type
399
400
401!-- Data type for the values of chemistry emissions
402    TYPE chem_emis_val_type
403
404       !REAL(wp),ALLOCATABLE, DIMENSION(:,:)     :: stack_height           !< stack height (ecc / to be implemented)
405       REAL(wp),ALLOCATABLE, DIMENSION(:,:,:)    :: default_emission_data  !< Emission input values for LOD1 (DEFAULT mode)
406       REAL(wp),ALLOCATABLE, DIMENSION(:,:,:,:)  :: preproc_emission_data  !< Emission input values for LOD2 (PRE-PROCESSED mode)
407
408    END TYPE chem_emis_val_type
409
410!
411!-- Define data structures for different input data types.
412!-- 8-bit Integer 2D
413    TYPE int_2d_8bit
414       INTEGER(KIND=1) ::  fill = -127                      !< fill value
415       INTEGER(KIND=1), DIMENSION(:,:), ALLOCATABLE ::  var !< respective variable
416
417       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
418                                        !< values are used
419    END TYPE int_2d_8bit
420!
421!-- 8-bit Integer 3D
422    TYPE int_3d_8bit
423       INTEGER(KIND=1) ::  fill = -127                           !< fill value
424       INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE ::  var_3d !< respective variable
425
426       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
427                                        !< values are used
428    END TYPE int_3d_8bit
429!
430!-- 32-bit Integer 2D
431    TYPE int_2d_32bit
432       INTEGER(iwp) ::  fill = -9999                      !< fill value
433       INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  var  !< respective variable
434
435       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
436                                        !< values are used
437    END TYPE int_2d_32bit
438!
439!-- Define data type to read 1D or 3D real variables.
440    TYPE real_1d_3d
441       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
442                                        !< values are used
443
444       INTEGER(iwp) ::  lod = -1        !< level-of-detail
445
446       REAL(wp) ::  fill = -9999.9_wp                  !< fill value
447
448       REAL(wp), DIMENSION(:),     ALLOCATABLE ::  var1d     !< respective 1D variable
449       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  var3d     !< respective 3D variable
450    END TYPE real_1d_3d
451!
452!-- Define data type to read 2D real variables
453    TYPE real_2d
454       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
455                                        !< values are used
456
457       INTEGER(iwp) ::  lod             !< level-of-detail
458
459       REAL(wp) ::  fill = -9999.9_wp                 !< fill value
460       REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  var  !< respective variable
461    END TYPE real_2d
462
463!
464!-- Define data type to read 3D real variables
465    TYPE real_3d
466       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
467                                        !< values are used
468
469       INTEGER(iwp) ::  nz   !< number of grid points along vertical dimension
470
471       REAL(wp) ::  fill = -9999.9_wp                   !< fill value
472       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  var  !< respective variable
473    END TYPE real_3d
474!
475!-- Define data structure where the dimension and type of the input depends on the given level of
476!-- detail.
477!-- For buildings, the input is either 2D float, or 3d byte.
478    TYPE build_in
479       INTEGER(iwp)    ::  lod = 1                               !< level of detail
480       INTEGER(KIND=1) ::  fill2 = -127                          !< fill value for lod = 2
481       INTEGER(iwp)    ::  nz                                    !< number of vertical layers in file
482       INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE ::  var_3d !< 3d variable (lod = 2)
483
484       REAL(wp), DIMENSION(:), ALLOCATABLE ::  z                 !< vertical coordinate for 3D building, used for consistency check
485
486       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
487                                        !< values are used
488
489       REAL(wp)                              ::  fill1 = -9999.9_wp !< fill values for lod = 1
490       REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  var_2d             !< 2d variable (lod = 1)
491       REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  oro_max            !< terraing height under particular buildings
492    END TYPE build_in
493
494!
495!-- For soil_type, the input is either 2D or 3D one-byte integer.
496    TYPE soil_in
497       INTEGER(iwp)                                   ::  lod = 1      !< level of detail
498       INTEGER(KIND=1)                                ::  fill = -127  !< fill value for lod = 2
499       INTEGER(KIND=1), DIMENSION(:,:), ALLOCATABLE   ::  var_2d       !< 2d variable (lod = 1)
500       INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE ::  var_3d       !< 3d variable (lod = 2)
501
502       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
503                                        !< values are used
504    END TYPE soil_in
505
506!
507!-- Define data type for fractions between surface types
508    TYPE fracs
509       INTEGER(iwp)                            ::  nf             !< total number of fractions
510       INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  nfracs         !< dimension array for fraction
511
512       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
513                                        !< values are used
514
515       REAL(wp)                                ::  fill = -9999.9_wp  !< fill value
516       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  frac               !< respective fraction between different surface types
517    END TYPE fracs
518!
519!-- Data type for parameter lists, Depending on the given level of detail, the input is 3D or 4D
520    TYPE pars
521       INTEGER(iwp)                            ::  lod = 1         !< level of detail
522       INTEGER(iwp)                            ::  np              !< total number of parameters
523       INTEGER(iwp)                            ::  nz              !< vertical dimension - number of soil layers
524       INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  layers          !< dimension array for soil layers
525       INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  pars            !< dimension array for parameters
526
527       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
528                                        !< values are used
529
530       REAL(wp)                                  ::  fill = -9999.9_wp !< fill value
531       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE   ::  pars_xy           !< respective parameters, level of detail = 1
532       REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::  pars_xyz          !< respective parameters, level of detail = 2
533    END TYPE pars
534!
535!-- Data type for surface parameter lists
536    TYPE pars_surf
537       INTEGER(iwp)                                ::  np          !< total number of parameters
538       INTEGER(iwp)                                ::  nsurf       !< number of local surfaces
539       INTEGER(iwp), DIMENSION(:,:,:), ALLOCATABLE ::  index_ji    !< index for beginning and end of surfaces at (j,i)
540       INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE   ::  coords      !< (k,j,i,norm_z,norm_y,norm_x)
541                                                                   !< k,j,i:                surface position
542                                                                   !< norm_z,norm_y,norm_x: surface normal vector
543
544       LOGICAL ::  from_file = .FALSE.  !< flag indicating whether an input variable is available and read from file or default
545                                        !< values are used
546
547       REAL(wp)                              ::  fill = -9999.9_wp !< fill value
548       REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pars              !< respective parameters per surface
549    END TYPE pars_surf
550!
551!-- Define type for global file attributes
552!-- Please refer to the PALM data standard for a detailed description of each attribute.
553    TYPE global_atts_type
554       CHARACTER(LEN=200) ::  acronym = ' '                           !< acronym of institution
555       CHARACTER(LEN=7)   ::  acronym_char = 'acronym'                !< name of attribute
556       CHARACTER(LEN=200) ::  author  = ' '                           !< first name, last name, email adress
557       CHARACTER(LEN=6)   ::  author_char = 'author'                  !< name of attribute
558       CHARACTER(LEN=200) ::  campaign = 'PALM-4U'                    !< name of campaign
559       CHARACTER(LEN=8)   ::  campaign_char = 'campaign'              !< name of attribute
560       CHARACTER(LEN=200) ::  comment = ' '                           !< comment to data
561       CHARACTER(LEN=7)   ::  comment_char = 'comment'                !< name of attribute
562       CHARACTER(LEN=200) ::  contact_person = ' '                    !< first name, last name, email adress
563       CHARACTER(LEN=14)  ::  contact_person_char = 'contact_person'  !< name of attribute
564       CHARACTER(LEN=200) ::  conventions = 'CF-1.7'                  !< netCDF convention
565       CHARACTER(LEN=11)  ::  conventions_char = 'Conventions'        !< name of attribute
566       CHARACTER(LEN=23 ) ::  creation_time = ' '                     !< creation time of data set
567       CHARACTER(LEN=13)  ::  creation_time_char = 'creation_time'    !< name of attribute
568       CHARACTER(LEN=200) ::  data_content = 'airmeteo'               !< content of data set
569       CHARACTER(LEN=12)  ::  data_content_char = 'data_content'      !< name of attribute
570       CHARACTER(LEN=200) ::  dependencies = ' '                      !< dependencies of data set
571       CHARACTER(LEN=12)  ::  dependencies_char = 'dependencies'      !< name of attribute
572       CHARACTER(LEN=200) ::  history = ' '                           !< information about data processing
573       CHARACTER(LEN=7)   ::  history_char = 'history'                !< name of attribute
574       CHARACTER(LEN=200) ::  institution = ' '                       !< name of responsible institution
575       CHARACTER(LEN=11)  ::  institution_char = 'institution'        !< name of attribute
576       CHARACTER(LEN=200) ::  keywords = ' '                          !< keywords of data set
577       CHARACTER(LEN=8)   ::  keywords_char = 'keywords'              !< name of attribute
578       CHARACTER(LEN=200) ::  licence = ' '                           !< licence of data set
579       CHARACTER(LEN=7)   ::  licence_char = 'licence'                !< name of attribute
580       CHARACTER(LEN=200) ::  location = ' '                          !< place which refers to data set
581       CHARACTER(LEN=8)   ::  location_char = 'location'              !< name of attribute
582       CHARACTER(LEN=10)  ::  origin_lat_char = 'origin_lat'          !< name of attribute
583       CHARACTER(LEN=10)  ::  origin_lon_char = 'origin_lon'          !< name of attribute
584       CHARACTER(LEN=23 ) ::  origin_time = '2000-01-01 00:00:00 +00' !< reference time
585       CHARACTER(LEN=11)  ::  origin_time_char = 'origin_time'        !< name of attribute
586       CHARACTER(LEN=8)   ::  origin_x_char = 'origin_x'              !< name of attribute
587       CHARACTER(LEN=8)   ::  origin_y_char = 'origin_y'              !< name of attribute
588       CHARACTER(LEN=8)   ::  origin_z_char = 'origin_z'              !< name of attribute
589       CHARACTER(LEN=12)  ::  palm_version_char = 'palm_version'      !< name of attribute
590       CHARACTER(LEN=200) ::  references = ' '                        !< literature referring to data set
591       CHARACTER(LEN=10)  ::  references_char = 'references'          !< name of attribute
592       CHARACTER(LEN=14)  ::  rotation_angle_char = 'rotation_angle'  !< name of attribute
593       CHARACTER(LEN=200) ::  site = ' '                              !< name of model domain
594       CHARACTER(LEN=4)   ::  site_char = 'site'                      !< name of attribute
595       CHARACTER(LEN=200) ::  source = ' '                            !< source of data set
596       CHARACTER(LEN=6)   ::  source_char = 'source'                  !< name of attribute
597       CHARACTER(LEN=200) ::  title = ' '                             !< title of data set
598       CHARACTER(LEN=5)   ::  title_char = 'title'                    !< name of attribute
599       CHARACTER(LEN=7)   ::  version_char = 'version'                !< name of attribute
600
601       INTEGER(iwp) ::  version              !< version of data set
602
603       REAL(wp) ::  fillvalue = -9999.0      !< default fill value
604       REAL(wp) ::  origin_lat               !< latitude of lower left corner
605       REAL(wp) ::  origin_lon               !< longitude of lower left corner
606       REAL(wp) ::  origin_x                 !< easting (UTM coordinate) of lower left corner
607       REAL(wp) ::  origin_y                 !< northing (UTM coordinate) of lower left corner
608       REAL(wp) ::  origin_z                 !< reference height
609       REAL(wp) ::  palm_version             !< PALM version of data set
610       REAL(wp) ::  rotation_angle           !< rotation angle of coordinate system of data set
611    END TYPE global_atts_type
612!
613!-- Define type for coordinate reference system (crs)
614    TYPE crs_type
615       CHARACTER(LEN=200) ::  epsg_code = 'EPSG:25831'                   !< EPSG code
616       CHARACTER(LEN=200) ::  grid_mapping_name = 'transverse_mercator'  !< name of grid mapping
617       CHARACTER(LEN=200) ::  long_name = 'coordinate reference system'  !< name of variable crs
618       CHARACTER(LEN=200) ::  units = 'm'                                !< unit of crs
619
620       REAL(wp) ::  false_easting = 500000.0_wp                   !< false easting
621       REAL(wp) ::  false_northing = 0.0_wp                       !< false northing
622       REAL(wp) ::  inverse_flattening = 298.257223563_wp         !< 1/f (default for WGS84)
623       REAL(wp) ::  latitude_of_projection_origin = 0.0_wp        !< latitude of projection origin
624       REAL(wp) ::  longitude_of_central_meridian = 3.0_wp        !< longitude of central meridian of UTM zone (default: zone 31)
625       REAL(wp) ::  longitude_of_prime_meridian = 0.0_wp          !< longitude of prime meridian
626       REAL(wp) ::  scale_factor_at_central_meridian = 0.9996_wp  !< scale factor of UTM coordinates
627       REAL(wp) ::  semi_major_axis = 6378137.0_wp                !< length of semi major axis (default for WGS84)
628    END TYPE crs_type
629
630!
631!-- Define variables
632    TYPE(crs_type)  ::  coord_ref_sys  !< coordinate reference system
633
634    TYPE(dims_xy)   ::  dim_static     !< data structure for x, y-dimension in static input file
635
636    TYPE(init_type) ::  init_3d     !< data structure for the initialization of the 3D flow and soil fields
637    TYPE(init_type) ::  init_model  !< data structure for the initialization of the model
638
639!
640!-- Define 2D variables of type NC_BYTE
641    TYPE(int_2d_8bit)  ::  albedo_type_f     !< input variable for albedo type
642    TYPE(int_2d_8bit)  ::  building_type_f   !< input variable for building type
643    TYPE(int_2d_8bit)  ::  pavement_type_f   !< input variable for pavenment type
644    TYPE(int_2d_8bit)  ::  street_crossing_f !< input variable for water type
645    TYPE(int_2d_8bit)  ::  street_type_f     !< input variable for water type
646    TYPE(int_2d_8bit)  ::  vegetation_type_f !< input variable for vegetation type
647    TYPE(int_2d_8bit)  ::  water_type_f      !< input variable for water type
648!
649!-- Define 3D variables of type NC_BYTE
650    TYPE(int_3d_8bit)  ::  building_obstruction_f    !< input variable for building obstruction
651    ! TYPE(int_3d_8bit)  ::  building_obstruction_full !< input variable for building obstruction
652!
653!-- Define 2D variables of type NC_INT
654    TYPE(int_2d_32bit) ::  building_id_f     !< input variable for building ID
655!
656!-- Define 2D variables of type NC_FLOAT
657    TYPE(real_2d) ::  terrain_height_f       !< input variable for terrain height
658    TYPE(real_2d) ::  uvem_irradiance_f      !< input variable for uvem irradiance lookup table
659    TYPE(real_2d) ::  uvem_integration_f     !< input variable for uvem integration
660!
661!-- Define 3D variables of type NC_FLOAT
662    TYPE(real_3d) ::  root_area_density_lsm_f !< input variable for root area density - parametrized vegetation
663    TYPE(real_3d) ::  uvem_radiance_f         !< input variable for uvem radiance lookup table
664    TYPE(real_3d) ::  uvem_projarea_f         !< input variable for uvem projection area lookup table
665!
666!-- Define input variable for buildings
667    TYPE(build_in) ::  buildings_f           !< input variable for buildings
668!
669!-- Define input variables for soil_type
670    TYPE(soil_in)  ::  soil_type_f           !< input variable for soil type
671
672    TYPE(fracs) ::  surface_fraction_f       !< input variable for surface fraction
673
674    TYPE(pars)  ::  albedo_pars_f              !< input variable for albedo parameters
675    TYPE(pars)  ::  building_pars_f            !< input variable for building parameters
676    TYPE(pars)  ::  pavement_pars_f            !< input variable for pavement parameters
677    TYPE(pars)  ::  pavement_subsurface_pars_f !< input variable for pavement parameters
678    TYPE(pars)  ::  soil_pars_f                !< input variable for soil parameters
679    TYPE(pars)  ::  vegetation_pars_f          !< input variable for vegetation parameters
680    TYPE(pars)  ::  water_pars_f               !< input variable for water parameters
681
682    TYPE(pars_surf)  ::  building_surface_pars_f  !< input variable for building surface parameters
683
684    TYPE(chem_emis_att_type)                             ::  chem_emis_att    !< Input Information of Chemistry Emission Data from
685                                                                              !< netcdf
686    TYPE(chem_emis_val_type), ALLOCATABLE, DIMENSION(:)  ::  chem_emis        !< Input Chemistry Emission Data from netcdf
687
688    CHARACTER(LEN=3)  ::  char_lod  = 'lod'         !< name of level-of-detail attribute in NetCDF file
689
690    CHARACTER(LEN=10) ::  char_fill = '_FillValue'        !< name of fill value attribute in NetCDF file
691
692    CHARACTER(LEN=100) ::  input_file_static  = 'PIDS_STATIC'  !< Name of file which comprises static input data
693    CHARACTER(LEN=100) ::  input_file_dynamic = 'PIDS_DYNAMIC' !< Name of file which comprises dynamic input data
694    CHARACTER(LEN=100) ::  input_file_chem    = 'PIDS_CHEM'    !< Name of file which comprises chemistry input data
695    CHARACTER(LEN=100) ::  input_file_uvem    = 'PIDS_UVEM'    !< Name of file which comprises static uv_exposure model input data
696    CHARACTER(LEN=100) ::  input_file_vm      = 'PIDS_VM'      !< Name of file which comprises virtual measurement data
697    CHARACTER(LEN=100) ::  input_file_wtm     = 'PIDS_WTM'     !< Name of file which comprises wind turbine model input data
698
699    CHARACTER(LEN=25), DIMENSION(:), ALLOCATABLE ::  string_values  !< output of string variables read from netcdf input files
700    CHARACTER(LEN=50), DIMENSION(:), ALLOCATABLE ::  vars_pids      !< variable in input file
701
702    INTEGER(iwp)                                     ::  id_emis        !< NetCDF id of input file for chemistry emissions: TBD: It
703                                                                        !< has to be removed
704
705    INTEGER(iwp) ::  nc_stat         !< return value of nf90 function call
706    INTEGER(iwp) ::  num_var_pids    !< number of variables in file
707    INTEGER(iwp) ::  pids_id         !< file id
708
709    LOGICAL ::  input_pids_static  = .FALSE.   !< Flag indicating whether Palm-input-data-standard file containing static
710                                               !< information exists
711    LOGICAL ::  input_pids_dynamic = .FALSE.   !< Flag indicating whether Palm-input-data-standard file containing dynamic
712                                               !< information exists
713    LOGICAL ::  input_pids_chem    = .FALSE.   !< Flag indicating whether Palm-input-data-standard file containing chemistry
714                                               !< information exists
715    LOGICAL ::  input_pids_uvem    = .FALSE.   !< Flag indicating whether uv-expoure-model input file containing static information
716                                               !< exists
717    LOGICAL ::  input_pids_vm      = .FALSE.   !< Flag indicating whether input file for virtual measurements exist
718    LOGICAL ::  input_pids_wtm     = .FALSE.   !< Flag indicating whether input file for wind turbine model exists
719
720    LOGICAL ::  collective_read = .FALSE.      !< Enable NetCDF collective read
721
722    REAL(wp), DIMENSION(8) ::  crs_list        !< list of coord_ref_sys values
723
724    TYPE(global_atts_type) ::  input_file_atts !< global attributes of input file
725
726    SAVE
727
728    PRIVATE
729
730    INTERFACE netcdf_data_input_check_dynamic
731       MODULE PROCEDURE netcdf_data_input_check_dynamic
732    END INTERFACE netcdf_data_input_check_dynamic
733
734    INTERFACE netcdf_data_input_check_static
735       MODULE PROCEDURE netcdf_data_input_check_static
736    END INTERFACE netcdf_data_input_check_static
737
738    INTERFACE netcdf_data_input_chemistry_data
739       MODULE PROCEDURE netcdf_data_input_chemistry_data
740    END INTERFACE netcdf_data_input_chemistry_data
741
742    INTERFACE get_dimension_length
743       MODULE PROCEDURE get_dimension_length
744    END INTERFACE get_dimension_length
745
746    INTERFACE inquire_fill_value
747       MODULE PROCEDURE inquire_fill_value_int
748       MODULE PROCEDURE inquire_fill_value_real
749    END INTERFACE inquire_fill_value
750
751    INTERFACE netcdf_data_input_inquire_file
752       MODULE PROCEDURE netcdf_data_input_inquire_file
753    END INTERFACE netcdf_data_input_inquire_file
754
755    INTERFACE netcdf_data_input_init
756       MODULE PROCEDURE netcdf_data_input_init
757    END INTERFACE netcdf_data_input_init
758
759    INTERFACE netcdf_data_input_init_3d
760       MODULE PROCEDURE netcdf_data_input_init_3d
761    END INTERFACE netcdf_data_input_init_3d
762
763    INTERFACE netcdf_data_input_surface_data
764       MODULE PROCEDURE netcdf_data_input_surface_data
765    END INTERFACE netcdf_data_input_surface_data
766
767    INTERFACE netcdf_data_input_uvem
768       MODULE PROCEDURE netcdf_data_input_uvem
769    END INTERFACE netcdf_data_input_uvem
770
771    INTERFACE get_variable
772       MODULE PROCEDURE get_variable_1d_char
773       MODULE PROCEDURE get_variable_1d_int
774       MODULE PROCEDURE get_variable_1d_real
775       MODULE PROCEDURE get_variable_2d_int8
776       MODULE PROCEDURE get_variable_2d_int32
777       MODULE PROCEDURE get_variable_2d_real
778       MODULE PROCEDURE get_variable_2d_real_dynamic
779       MODULE PROCEDURE get_variable_3d_int8
780       MODULE PROCEDURE get_variable_3d_real
781       MODULE PROCEDURE get_variable_3d_real_dynamic
782       MODULE PROCEDURE get_variable_4d_to_3d_real
783       MODULE PROCEDURE get_variable_4d_real
784       MODULE PROCEDURE get_variable_5d_to_4d_real
785       MODULE PROCEDURE get_variable_5d_real           ! (ecc) temp subroutine 4 reading 5D NC arrays
786       MODULE PROCEDURE get_variable_5d_real_dynamic   ! 2B removed as z is out of emission_values
787       MODULE PROCEDURE get_variable_string
788       MODULE PROCEDURE get_variable_string_generic    ! (ecc) generic string function
789
790    END INTERFACE get_variable
791
792    INTERFACE get_variable_pr
793       MODULE PROCEDURE get_variable_pr
794    END INTERFACE get_variable_pr
795
796    INTERFACE get_attribute
797       MODULE PROCEDURE get_attribute_real
798       MODULE PROCEDURE get_attribute_int8
799       MODULE PROCEDURE get_attribute_int32
800       MODULE PROCEDURE get_attribute_string
801    END INTERFACE get_attribute
802
803    INTERFACE add_ghost_layers
804       MODULE PROCEDURE add_ghost_layers_2d_int8
805       MODULE PROCEDURE add_ghost_layers_2d_int32
806       MODULE PROCEDURE add_ghost_layers_2d_real
807       MODULE PROCEDURE add_ghost_layers_3d_int8
808       MODULE PROCEDURE add_ghost_layers_3d_real
809       MODULE PROCEDURE add_ghost_layers_4d_real
810    END INTERFACE add_ghost_layers
811
812!
813!-- Public data structures
814    PUBLIC real_1d_3d,                                                                             &
815           real_2d,                                                                                &
816           real_3d
817!
818!-- Public variables
819    PUBLIC albedo_pars_f,                                                                          &
820           albedo_type_f,                                                                          &
821           buildings_f,                                                                            &
822           building_id_f,                                                                          &
823           building_pars_f,                                                                        &
824           building_surface_pars_f,                                                                &
825           building_type_f,                                                                        &
826           char_fill,                                                                              &
827           char_lod,                                                                               &
828           chem_emis,                                                                              &
829           chem_emis_att,                                                                          &
830           chem_emis_att_type,                                                                     &
831           chem_emis_val_type,                                                                     &
832           coord_ref_sys,                                                                          &
833           crs_list,                                                                               &
834           init_3d,                                                                                &
835           init_model,                                                                             &
836           input_file_atts,                                                                        &
837           input_file_dynamic,                                                                     &
838           input_file_static,                                                                      &
839           input_file_vm,                                                                          &
840           input_file_wtm,                                                                         &
841           input_pids_static,                                                                      &
842           input_pids_dynamic,                                                                     &
843           input_pids_vm,                                                                          &
844           input_pids_wtm,                                                                         &
845           num_var_pids,                                                                           &
846           pavement_pars_f,                                                                        &
847           pavement_subsurface_pars_f,                                                             &
848           pavement_type_f,                                                                        &
849           pids_id,                                                                                &
850           root_area_density_lsm_f,                                                                &
851           soil_pars_f,                                                                            &
852           soil_type_f,                                                                            &
853           street_crossing_f,                                                                      &
854           street_type_f,                                                                          &
855           surface_fraction_f,                                                                     &
856           terrain_height_f,                                                                       &
857           vars_pids,                                                                              &
858           vegetation_pars_f,                                                                      &
859           vegetation_type_f,                                                                      &
860           water_pars_f,                                                                           &
861           water_type_f
862!
863!-- Public uv exposure variables
864    PUBLIC building_obstruction_f,                                                                 &
865           input_file_uvem,                                                                        &
866           input_pids_uvem,                                                                        &
867           netcdf_data_input_uvem,                                                                 &
868           uvem_integration_f,                                                                     &
869           uvem_irradiance_f,                                                                      &
870           uvem_projarea_f,                                                                        &
871           uvem_radiance_f
872
873!
874!-- Public subroutines
875    PUBLIC add_ghost_layers,                                                                       &
876           check_existence,                                                                        &
877           close_input_file,                                                                       &
878           get_attribute,                                                                          &
879           get_dimension_length,                                                                   &
880           get_variable,                                                                           &
881           get_variable_pr,                                                                        &
882           inquire_fill_value,                                                                     &
883           inquire_num_variables,                                                                  &
884           inquire_variable_names,                                                                 &
885           netcdf_data_input_check_dynamic,                                                        &
886           netcdf_data_input_check_static,                                                         &
887           netcdf_data_input_chemistry_data,                                                       &
888           netcdf_data_input_inquire_file,                                                         &
889           netcdf_data_input_init,                                                                 &
890           netcdf_data_input_init_3d,                                                              &
891           netcdf_data_input_surface_data,                                                         &
892           netcdf_data_input_topo,                                                                 &
893           open_read_file
894
895
896 CONTAINS
897
898!--------------------------------------------------------------------------------------------------!
899! Description:
900! ------------
901!> Inquires whether NetCDF input files according to Palm-input-data standard exist. Moreover, basic
902!> checks are performed.
903!--------------------------------------------------------------------------------------------------!
904 SUBROUTINE netcdf_data_input_inquire_file
905
906    USE control_parameters,                                                                        &
907        ONLY:  topo_no_distinct
908
909    IMPLICIT NONE
910
911#if defined ( __netcdf )
912    INQUIRE( FILE = TRIM( input_file_static )   // TRIM( coupling_char ),                          &
913             EXIST = input_pids_static  )
914    INQUIRE( FILE = TRIM( input_file_dynamic ) // TRIM( coupling_char ),                           &
915             EXIST = input_pids_dynamic )
916    INQUIRE( FILE = TRIM( input_file_chem )    // TRIM( coupling_char ),                           &
917             EXIST = input_pids_chem )
918    INQUIRE( FILE = TRIM( input_file_uvem )    // TRIM( coupling_char ),                           &
919             EXIST = input_pids_uvem )
920    INQUIRE( FILE = TRIM( input_file_vm )      // TRIM( coupling_char ),                           &
921             EXIST = input_pids_vm )
922    INQUIRE( FILE = TRIM( input_file_wtm )     // TRIM( coupling_char ),                           &
923             EXIST = input_pids_wtm )
924#endif
925
926!
927!-- As long as topography can be input via ASCII format, no distinction between building and terrain
928!-- can be done. This case, classify all surfaces as default type. Same in case land-surface and
929!-- urban-surface model are not applied.
930    IF ( .NOT. input_pids_static )  THEN
931       topo_no_distinct = .TRUE.
932    ENDIF
933
934 END SUBROUTINE netcdf_data_input_inquire_file
935
936!--------------------------------------------------------------------------------------------------!
937! Description:
938! ------------
939!> Reads global attributes and coordinate reference system required for initialization of the model.
940!--------------------------------------------------------------------------------------------------!
941 SUBROUTINE netcdf_data_input_init
942
943    IMPLICIT NONE
944
945    INTEGER(iwp) ::  id_mod     !< NetCDF id of input file
946    INTEGER(iwp) ::  var_id_crs !< NetCDF id of variable crs
947
948!
949!-- Define default list of crs attributes. This is required for coordinate transformation.
950    crs_list = (/ coord_ref_sys%semi_major_axis,                                                   &
951                  coord_ref_sys%inverse_flattening,                                                &
952                  coord_ref_sys%longitude_of_prime_meridian,                                       &
953                  coord_ref_sys%longitude_of_central_meridian,                                     &
954                  coord_ref_sys%scale_factor_at_central_meridian,                                  &
955                  coord_ref_sys%latitude_of_projection_origin,                                     &
956                  coord_ref_sys%false_easting,                                                     &
957                  coord_ref_sys%false_northing /)
958
959    IF ( .NOT. input_pids_static )  RETURN
960
961#if defined ( __netcdf )
962!
963!-- Open file in read-only mode
964    CALL open_read_file( TRIM( input_file_static ) // TRIM( coupling_char ), id_mod )
965!
966!-- Read global attributes
967    CALL get_attribute( id_mod, input_file_atts%origin_lat_char,                                   &
968                        input_file_atts%origin_lat, .TRUE. )
969
970    CALL get_attribute( id_mod, input_file_atts%origin_lon_char,                                   &
971                        input_file_atts%origin_lon, .TRUE. )
972
973    CALL get_attribute( id_mod, input_file_atts%origin_time_char,                                  &
974                        input_file_atts%origin_time, .TRUE. )
975
976    CALL get_attribute( id_mod, input_file_atts%origin_x_char,                                     &
977                        input_file_atts%origin_x, .TRUE. )
978
979    CALL get_attribute( id_mod, input_file_atts%origin_y_char,                                     &
980                        input_file_atts%origin_y, .TRUE. )
981
982    CALL get_attribute( id_mod, input_file_atts%origin_z_char,                                     &
983                        input_file_atts%origin_z, .TRUE. )
984
985    CALL get_attribute( id_mod, input_file_atts%rotation_angle_char,                               &
986                        input_file_atts%rotation_angle, .TRUE. )
987
988    CALL get_attribute( id_mod, input_file_atts%author_char,                                       &
989                        input_file_atts%author, .TRUE., no_abort=.FALSE. )
990
991    CALL get_attribute( id_mod, input_file_atts%contact_person_char,                               &
992                        input_file_atts%contact_person, .TRUE., no_abort=.FALSE. )
993    CALL get_attribute( id_mod, input_file_atts%institution_char,                                  &
994                        input_file_atts%institution,    .TRUE., no_abort=.FALSE. )
995    CALL get_attribute( id_mod, input_file_atts%acronym_char,                                      &
996                        input_file_atts%acronym,        .TRUE., no_abort=.FALSE. )
997
998    CALL get_attribute( id_mod, input_file_atts%campaign_char,                                     &
999                        input_file_atts%campaign, .TRUE., no_abort=.FALSE. )
1000    CALL get_attribute( id_mod, input_file_atts%location_char,                                     &
1001                        input_file_atts%location, .TRUE., no_abort=.FALSE. )
1002    CALL get_attribute( id_mod, input_file_atts%site_char,                                         &
1003                        input_file_atts%site,     .TRUE., no_abort=.FALSE. )
1004
1005    CALL get_attribute( id_mod, input_file_atts%source_char,                                       &
1006                        input_file_atts%source,     .TRUE., no_abort=.FALSE. )
1007    CALL get_attribute( id_mod, input_file_atts%references_char,                                   &
1008                        input_file_atts%references, .TRUE., no_abort=.FALSE. )
1009    CALL get_attribute( id_mod, input_file_atts%keywords_char,                                     &
1010                        input_file_atts%keywords,   .TRUE., no_abort=.FALSE. )
1011    CALL get_attribute( id_mod, input_file_atts%licence_char,                                      &
1012                        input_file_atts%licence,    .TRUE., no_abort=.FALSE. )
1013    CALL get_attribute( id_mod, input_file_atts%comment_char,                                      &
1014                        input_file_atts%comment,    .TRUE., no_abort=.FALSE. )
1015!
1016!-- Read coordinate reference system if available
1017    nc_stat = NF90_INQ_VARID( id_mod, 'crs', var_id_crs )
1018    IF ( nc_stat == NF90_NOERR )  THEN
1019       CALL get_attribute( id_mod, 'epsg_code', coord_ref_sys%epsg_code, .FALSE., 'crs' )
1020       CALL get_attribute( id_mod, 'false_easting', coord_ref_sys%false_easting, .FALSE., 'crs' )
1021       CALL get_attribute( id_mod, 'false_northing', coord_ref_sys%false_northing, .FALSE., 'crs' )
1022       CALL get_attribute( id_mod, 'grid_mapping_name', coord_ref_sys%grid_mapping_name, .FALSE.,  &
1023                           'crs' )
1024       CALL get_attribute( id_mod, 'inverse_flattening', coord_ref_sys%inverse_flattening, .FALSE.,&
1025                           'crs' )
1026       CALL get_attribute( id_mod, 'latitude_of_projection_origin',                                &
1027                           coord_ref_sys%latitude_of_projection_origin, .FALSE., 'crs' )
1028       CALL get_attribute( id_mod, 'long_name', coord_ref_sys%long_name, .FALSE., 'crs' )
1029       CALL get_attribute( id_mod, 'longitude_of_central_meridian',                                &
1030                           coord_ref_sys%longitude_of_central_meridian, .FALSE., 'crs' )
1031       CALL get_attribute( id_mod, 'longitude_of_prime_meridian',                                  &
1032                           coord_ref_sys%longitude_of_prime_meridian, .FALSE., 'crs' )
1033       CALL get_attribute( id_mod, 'scale_factor_at_central_meridian',                             &
1034                           coord_ref_sys%scale_factor_at_central_meridian, .FALSE., 'crs' )
1035       CALL get_attribute( id_mod, 'semi_major_axis', coord_ref_sys%semi_major_axis, .FALSE., 'crs')
1036       CALL get_attribute( id_mod, 'units', coord_ref_sys%units, .FALSE., 'crs' )
1037    ELSE
1038!
1039!--    Calculate central meridian from origin_lon
1040       coord_ref_sys%longitude_of_central_meridian =                                               &
1041                                   CEILING( input_file_atts%origin_lon / 6.0_wp ) * 6.0_wp - 3.0_wp
1042    ENDIF
1043!
1044!-- Finally, close input file
1045    CALL close_input_file( id_mod )
1046#endif
1047!
1048!-- Copy latitude, longitude, origin_z, rotation angle on init type.
1049!-- Note: A shifting height might have already been saved to orgin_z in init_grid; therefore, do not
1050!--       override but add the reference height from the input file.
1051    init_model%latitude       = input_file_atts%origin_lat
1052    init_model%longitude      = input_file_atts%origin_lon
1053    init_model%origin_time    = input_file_atts%origin_time
1054    init_model%origin_x       = input_file_atts%origin_x
1055    init_model%origin_y       = input_file_atts%origin_y
1056    init_model%origin_z       = init_model%origin_z + input_file_atts%origin_z
1057    init_model%rotation_angle = input_file_atts%rotation_angle
1058
1059!
1060!-- Update list of crs attributes. This is required for coordinate transformation.
1061    crs_list = (/ coord_ref_sys%semi_major_axis,                                                   &
1062                  coord_ref_sys%inverse_flattening,                                                &
1063                  coord_ref_sys%longitude_of_prime_meridian,                                       &
1064                  coord_ref_sys%longitude_of_central_meridian,                                     &
1065                  coord_ref_sys%scale_factor_at_central_meridian,                                  &
1066                  coord_ref_sys%latitude_of_projection_origin,                                     &
1067                  coord_ref_sys%false_easting,                                                     &
1068                  coord_ref_sys%false_northing /)
1069!
1070!-- In case of nested runs, each model domain might have different longitude and latitude, which
1071!-- would result in different Coriolis parameters and sun-zenith angles. To avoid this, longitude
1072!-- and latitude in each model domain will be set to the values of the root model. Please note, this
1073!-- synchronization is required already here.
1074#if defined( __parallel )
1075    CALL MPI_BCAST( init_model%latitude,  1, MPI_REAL, 0, MPI_COMM_WORLD, ierr )
1076    CALL MPI_BCAST( init_model%longitude, 1, MPI_REAL, 0, MPI_COMM_WORLD, ierr )
1077#endif
1078
1079 END SUBROUTINE netcdf_data_input_init
1080
1081
1082!--------------------------------------------------------------------------------------------------!
1083! Description:
1084! ------------
1085!> Reads Chemistry NETCDF Input data, such as emission values, emission species, etc.
1086!--------------------------------------------------------------------------------------------------!
1087 SUBROUTINE netcdf_data_input_chemistry_data( emt_att, emt )
1088
1089    USE chem_modules,                                                                              &
1090        ONLY:  emiss_lod, surface_csflux_name, time_fac_type
1091
1092    USE control_parameters,                                                                        &
1093        ONLY:  message_string
1094
1095    USE indices,                                                                                   &
1096        ONLY:  nxl, nxr, nys, nyn
1097
1098    IMPLICIT NONE
1099
1100    TYPE(chem_emis_att_type), INTENT(INOUT)                            ::  emt_att
1101    TYPE(chem_emis_val_type), ALLOCATABLE, DIMENSION(:), INTENT(INOUT) ::  emt
1102
1103    INTEGER(iwp) ::  i, j, k      !< generic counters
1104    INTEGER(iwp) ::  ispec        !< index for number of emission species in input
1105    INTEGER(iwp) ::  len_dims     !< Length of dimension
1106    INTEGER(iwp) ::  num_vars     !< number of variables in netcdf input file
1107
1108!
1109!-- dum_var_4d are designed to read in emission_values from the chemistry netCDF file.
1110!-- Currently the vestigial "z" dimension in emission_values makes it a 5D array, hence the
1111!-- corresponding dum_var_5d array. When the "z" dimension is removed completely, dum_var_4d will be
1112!-- used instead (ecc 20190425)
1113
1114!    REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:)    ::  dum_var_4d  !< temp array 4 4D chem emission data
1115    REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:,:) ::  dum_var_5d  !< temp array 4 5D chem emission data
1116
1117!
1118!-- Start processing data
1119!
1120!-- Emission LOD 0 (Parameterized mode)
1121    IF ( emiss_lod == 0 )  THEN
1122
1123!
1124!--    For reference (ecc)
1125!       IF (TRIM(mode_emis) == "PARAMETERIZED" .OR. TRIM(mode_emis) == "parameterized") THEN
1126
1127       ispec=1
1128       emt_att%n_emiss_species = 0
1129
1130!
1131!--    Number of species
1132       DO  WHILE ( TRIM( surface_csflux_name( ispec ) ) /= 'novalue' )
1133
1134          emt_att%n_emiss_species = emt_att%n_emiss_species + 1
1135          ispec = ispec + 1
1136!
1137!--       Followling line retained for compatibility with salsa_mod which still uses emt_att%nspec
1138!--       heavily (ecc)
1139          emt_att%nspec = emt_att%nspec + 1
1140
1141       ENDDO
1142
1143!
1144!--    Allocate emission values data type arrays
1145       ALLOCATE( emt(emt_att%n_emiss_species) )
1146
1147!
1148!--    Read emission species names
1149
1150!
1151!--    Allocate space for strings
1152       ALLOCATE( emt_att%species_name(emt_att%n_emiss_species) )
1153
1154       DO  ispec = 1, emt_att%n_emiss_species
1155          emt_att%species_name(ispec) = TRIM( surface_csflux_name(ispec) )
1156       ENDDO
1157
1158!
1159!-- LOD 1 (default mode) and LOD 2 (pre-processed mode)
1160    ELSE
1161
1162#if defined ( __netcdf )
1163
1164       IF ( .NOT. input_pids_chem )  RETURN
1165
1166!
1167!--    First we allocate memory space for the emission species and then we differentiate between
1168!--    LOD 1 (default mode) and LOD 2 (pre-processed mode)
1169
1170!
1171!--    Open emission data file ( {palmcase}_chemistry )
1172       CALL open_read_file( TRIM( input_file_chem ) // TRIM( coupling_char ), id_emis )
1173
1174!
1175!--    Inquire number of variables
1176       CALL inquire_num_variables( id_emis, num_vars )
1177
1178!
1179!--    Get general dimension lengths: only # species and # categories.
1180!--    Other dimensions depend on the emission mode or specific components.
1181       CALL get_dimension_length( id_emis, emt_att%n_emiss_species, 'nspecies' )
1182
1183!
1184!--    Backward compatibility for salsa_mod (ecc)
1185       emt_att%nspec = emt_att%n_emiss_species
1186
1187!
1188!--    Allocate emission values data type arrays
1189       ALLOCATE( emt(emt_att%n_emiss_species) )
1190
1191!
1192!--    Reading in species names
1193
1194!
1195!--    Allocate memory for species names
1196       ALLOCATE( emt_att%species_name(emt_att%n_emiss_species) )
1197
1198!
1199!--    Retrieve variable name (again, should use n_emiss_strlen)
1200       CALL get_variable( id_emis, 'emission_name', string_values, emt_att%n_emiss_species )
1201       emt_att%species_name=string_values
1202
1203!
1204!--    Dealocate string_values previously allocated in get_variable call
1205       IF ( ALLOCATED( string_values ) )  DEALLOCATE( string_values )
1206
1207!
1208!--    Reading in species indices
1209
1210!
1211!--    Allocate memory for species indices
1212       ALLOCATE( emt_att%species_index(emt_att%n_emiss_species) )
1213
1214!
1215!--    Retrieve variable data
1216       CALL get_variable( id_emis, 'emission_index', emt_att%species_index )
1217!
1218!--    Now the routine has to distinguish between chemistry emission LOD 1 (default mode) and LOD 2
1219!--    (pre-processed mode)
1220
1221!
1222!--    Start of emission LOD 1 (default mode)
1223       IF ( emiss_lod == 1 )  THEN
1224!
1225!--       For reference (ecc)
1226!          IF (TRIM(mode_emis) == "DEFAULT" .OR. TRIM(mode_emis) == "default") THEN
1227
1228!
1229!--       Get number of emission categories
1230          CALL get_dimension_length( id_emis, emt_att%ncat, 'ncat' )
1231
1232!--       Reading in emission categories indices
1233          ALLOCATE( emt_att%cat_index(emt_att%ncat) )
1234
1235!
1236!--       Retrieve variable data
1237          CALL get_variable( id_emis, 'emission_cat_index', emt_att%cat_index )
1238
1239
1240!
1241!-- Loop through individual species to get basic information on VOC/PM/NOX/SOX
1242
1243!---------------------------------------------------------------------------------------------------
1244!-- Note - Check array indices for reading in names and species in LOD1 (default mode) for the
1245!--        various mode splits as all id_emis conditionals have been removed from get_var_functions.
1246!--        In theory this would mean all arrays should be read from 0 to n-1 (C convention) as
1247!--        opposed to 1 to n (FORTRAN convention). Keep this in mind !!
1248!--        (ecc 20190424)
1249!--------------------------------------------------------------------------------------------------
1250
1251          DO  ispec = 1, emt_att%n_emiss_species
1252
1253!
1254!--          VOC data (name and composition)
1255             IF ( TRIM( emt_att%species_name(ispec) ) == "VOC"  .OR.                               &
1256                  TRIM( emt_att%species_name(ispec) ) == "voc" )  THEN
1257
1258!
1259!--             VOC name
1260                CALL get_dimension_length( id_emis, emt_att%nvoc, 'nvoc' )
1261                ALLOCATE( emt_att%voc_name(emt_att%nvoc) )
1262                CALL get_variable ( id_emis,"emission_voc_name", string_values, emt_att%nvoc )
1263                emt_att%voc_name = string_values
1264                IF ( ALLOCATED( string_values ) )  DEALLOCATE( string_values )
1265
1266!
1267!--             VOC composition
1268                ALLOCATE( emt_att%voc_comp(emt_att%ncat,emt_att%nvoc) )
1269                CALL get_variable( id_emis, "composition_voc", emt_att%voc_comp, 1, emt_att%ncat,  &
1270                                    1, emt_att%nvoc )
1271
1272             ENDIF  ! VOC
1273
1274!
1275!--          PM data (name and composition)
1276             IF ( TRIM( emt_att%species_name(ispec) ) == "PM"  .OR.                                &
1277                  TRIM( emt_att%species_name(ispec) ) == "pm")  THEN
1278
1279!
1280!--             PM name
1281                CALL get_dimension_length( id_emis, emt_att%npm, 'npm' )
1282                ALLOCATE( emt_att%pm_name(emt_att%npm) )
1283                CALL get_variable ( id_emis, "pm_name", string_values, emt_att%npm )
1284                emt_att%pm_name = string_values
1285                IF ( ALLOCATED( string_values ) )  DEALLOCATE( string_values )
1286
1287!
1288!--             PM composition (PM1, PM2.5 and PM10)
1289                len_dims = 3  ! PM1, PM2.5, PM10
1290                ALLOCATE( emt_att%pm_comp(emt_att%ncat,emt_att%npm,len_dims) )
1291                CALL get_variable( id_emis, "composition_pm", emt_att%pm_comp, 1, emt_att%ncat, 1, &
1292                                   emt_att%npm, 1, len_dims )
1293
1294             ENDIF  ! PM
1295
1296!
1297!--          NOX (NO and NO2)
1298             IF ( TRIM( emt_att%species_name(ispec) ) == "NOX"  .OR.                               &
1299                  TRIM( emt_att%species_name(ispec) ) == "nox" )  THEN
1300
1301                ALLOCATE( emt_att%nox_comp(emt_att%ncat,emt_att%nnox) )
1302                CALL get_variable( id_emis, "composition_nox", emt_att%nox_comp, 1, emt_att%ncat,  &
1303                                   1, emt_att%nnox )
1304
1305             ENDIF  ! NOX
1306
1307!
1308!--          SOX (SO2 and SO4)
1309             IF ( TRIM( emt_att%species_name(ispec) ) == "SOX"  .OR.                               &
1310                  TRIM( emt_att%species_name(ispec) ) == "sox" )  THEN
1311
1312                ALLOCATE( emt_att%sox_comp(emt_att%ncat,emt_att%nsox) )
1313                CALL get_variable( id_emis, "composition_sox", emt_att%sox_comp, 1, emt_att%ncat,  &
1314                                   1, emt_att%nsox )
1315
1316             ENDIF  ! SOX
1317
1318          ENDDO  ! do ispec
1319
1320!
1321!--       Emission time scaling factors (hourly and MDH data)
1322
1323!
1324!--       Hour
1325          IF ( TRIM( time_fac_type ) == "HOUR"  .OR.  TRIM( time_fac_type ) == "hour" )  THEN
1326
1327             CALL get_dimension_length( id_emis, emt_att%nhoursyear, 'nhoursyear' )
1328             ALLOCATE( emt_att%hourly_emis_time_factor(emt_att%ncat,emt_att%nhoursyear) )
1329             CALL get_variable( id_emis, "emission_time_factors", emt_att%hourly_emis_time_factor, &
1330                                1, emt_att%ncat, 1, emt_att%nhoursyear )
1331
1332!
1333!--       MDH
1334          ELSEIF ( TRIM( time_fac_type )  ==  "MDH"  .OR.  TRIM( time_fac_type )  ==  "mdh" )  THEN
1335
1336             CALL get_dimension_length( id_emis, emt_att%nmonthdayhour, 'nmonthdayhour' )
1337             ALLOCATE( emt_att%mdh_emis_time_factor(emt_att%ncat,emt_att%nmonthdayhour) )
1338             CALL get_variable( id_emis, "emission_time_factors", emt_att%mdh_emis_time_factor,    &
1339                                1, emt_att%ncat, 1, emt_att%nmonthdayhour )
1340
1341!
1342!--       Error (time factor undefined)
1343          ELSE
1344
1345             message_string = 'We are in the DEFAULT chemistry emissions mode: '  //               &
1346                              '     !no time-factor type specified!'              //               &
1347                              'Please specify the value of time_fac_type:'        //               &
1348                              '         either "MDH" or "HOUR"'
1349             CALL message( 'netcdf_data_input_chemistry_data', 'CM0200', 2, 2, 0, 6, 0 )
1350
1351
1352          ENDIF  ! time_fac_type
1353
1354!
1355!--       Read in default (LOD1) emissions from chemisty netCDF file per species
1356
1357!
1358!--       Note - At the moment the data is read in per species, but in the future it would be much
1359!--              more sensible to read in per species per time step to reduce memory consumption
1360!--              and, to a smaller degree, dimensionality of data exchange (I expect this will be
1361!--              necessary when the problem size is large).
1362          DO  ispec = 1, emt_att%n_emiss_species
1363
1364!
1365!--          Allocate space for species specific emission values.
1366!--          Note - This array is extended by 1 cell in each horizontal direction to compensate for
1367!--                 an apparent linear offset. The reason of this offset is not known but it has
1368!--                 been determined to take place beyond the scope of this module, and has little to
1369!--                 do with index conventions.
1370!--                 That is, setting the array horizontal limit from nx0:nx1 to 1:(nx1-nx0+1) or
1371!--                 nx0+1:nx1+1 did not result in correct or definite behavior.
1372!--                 This must be looked at at some point by the Hannover team but for now this
1373!--                 workaround is deemed reasonable (ecc 20190417).
1374             IF ( .NOT. ALLOCATED( emt(ispec)%default_emission_data ) )  THEN
1375                ALLOCATE( emt(ispec)%default_emission_data(emt_att%ncat,nys:nyn+1,nxl:nxr+1) )
1376             ENDIF
1377!
1378!--          Allocate dummy variable w/ index order identical to that shown in the netCDF header
1379             ALLOCATE( dum_var_5d(1,nys:nyn,nxl:nxr,1,emt_att%ncat) )
1380!
1381!--          Get variable.  Be very careful
1382!--          I am using get_variable_5d_real_dynamic (note logical argument at the end)
1383!--          1) use Fortran index convention (i.e., 1 to N)
1384!--          2) index order must be in reverse order from above allocation order
1385             CALL get_variable( id_emis, "emission_values", dum_var_5d, 1, ispec,                  &
1386                                nxl+1,     nys+1,     1, emt_att%ncat, 1,                          &
1387                                nxr-nxl+1, nyn-nys+1, emt_att%dt_emission, .FALSE. )
1388!
1389!--          Assign temp array to data structure then deallocate temp array
1390!--          Note - Indices are shifted from nx0:nx1 to nx0+1:nx1+1 to offset the emission data
1391!--                 array to counter said domain offset.
1392!--                 (ecc 20190417)
1393             DO  k = 1, emt_att%ncat
1394                DO  j = nys+1, nyn+1
1395                   DO  i = nxl+1, nxr+1
1396                      emt(ispec)%default_emission_data(k,j,i) = dum_var_5d(1,j-1,i-1,1,k)
1397                   ENDDO
1398                ENDDO
1399             ENDDO
1400
1401             DEALLOCATE( dum_var_5d )
1402
1403          ENDDO  ! ispec
1404!
1405!--       Units
1406          CALL get_attribute( id_emis, "units", emt_att%units, .FALSE., "emission_values" )
1407
1408!
1409!--       End default mode
1410
1411
1412!
1413!--    Start LOD 2 (pre-processed mode)
1414       ELSEIF ( emiss_lod == 2 )  THEN
1415
1416!      For reference (ecc)
1417!       ELSEIF( TRIM( mode_emis ) == "PRE-PROCESSED" .OR. TRIM( mode_emis ) == "pre-processed")  THEN
1418
1419!
1420!--       For LOD 2 only VOC and emission data need to be read
1421!--       Note - Check array indices for reading in names and species in LOD2 (pre-processed mode)
1422!--       for the various mode splits as all id_emis conditionals have been removed from
1423!--       get_var_functions. In theory this would mean all arrays should be read from 0 to n-1
1424!--       (C convention) as opposed to 1 to n (FORTRAN convention). Keep this in mind!!
1425!--       (ecc 20190424)
1426          DO  ispec = 1, emt_att%n_emiss_species
1427
1428!
1429!--          VOC data (name and composition)
1430             IF ( TRIM( emt_att%species_name(ispec) ) == "VOC"  .OR.                               &
1431                  TRIM( emt_att%species_name(ispec) ) == "voc" )  THEN
1432
1433!
1434!--             VOC name
1435                CALL get_dimension_length( id_emis, emt_att%nvoc, 'nvoc' )
1436                ALLOCATE( emt_att%voc_name(emt_att%nvoc) )
1437                CALL get_variable( id_emis, "emission_voc_name", string_values, emt_att%nvoc )
1438                emt_att%voc_name = string_values
1439                IF ( ALLOCATED( string_values ) )  DEALLOCATE( string_values )
1440
1441!
1442!--             VOC composition
1443                ALLOCATE( emt_att%voc_comp(emt_att%ncat,emt_att%nvoc) )
1444                CALL get_variable( id_emis, "composition_voc", emt_att%voc_comp, 1, emt_att%ncat,  &
1445                                   1, emt_att%nvoc )
1446             ENDIF  ! VOC
1447
1448          ENDDO  ! ispec
1449
1450!
1451!--       Emission data
1452          CALL get_dimension_length( id_emis, emt_att%dt_emission, 'time' )
1453
1454!
1455!--       Read in pre-processed (LOD2) emissions from chemisty netCDF file per species
1456
1457!
1458!--       Note - At the moment the data is read in per species, but in the future it would be much
1459!--              more sensible to read in per species per time step to reduce memory consumption
1460!--              and, to a smaller degree, dimensionality of data exchange (I expect this will be
1461!--              necessary when the problem size is large).
1462          DO  ispec = 1, emt_att%n_emiss_species
1463
1464!
1465!--          Allocate space for species specific emission values.
1466!--          Note - This array is extended by 1 cell in each horizontal direction to compensate for
1467!--                 an apparent linear offset. The reason of this offset is not known but it has
1468!--                 been determined to take place beyond the scope of this module, and has little to
1469!--                 do with index conventions.
1470!--                 That is, setting the array horizontal limit from nx0:nx1 to 1:(nx1-nx0+1) or
1471!--                 nx0+1:nx1+1 did not result in correct or definite behavior.
1472!--                 This must be looked at at some point by the Hannover team but for now this
1473!--                 workaround is deemed reasonable (ecc 20190417).
1474
1475             IF ( .NOT. ALLOCATED( emt(ispec)%preproc_emission_data ) )  THEN
1476                ALLOCATE( emt(ispec)%preproc_emission_data(                                        &
1477                             emt_att%dt_emission, 1, nys:nyn+1, nxl:nxr+1) )
1478             ENDIF
1479!
1480!--          Allocate dummy variable w/ index order identical to that shown in the netCDF header
1481             ALLOCATE( dum_var_5d(emt_att%dt_emission,1,nys:nyn,nxl:nxr,1) )
1482!
1483!--          Get variable.  Be very careful.
1484!--          I am using get_variable_5d_real_dynamic (note logical argument at the end)
1485!--          1) use Fortran index convention (i.e., 1 to N)
1486!--          2) index order must be in reverse order from above allocation order
1487             CALL get_variable( id_emis, "emission_values", dum_var_5d, ispec,                     &
1488                                nxl+1,     nys+1,     1, 1, 1,                                     &
1489                                nxr-nxl+1, nyn-nys+1, 1, emt_att%dt_emission, .FALSE. )
1490!
1491!--          Assign temp array to data structure then deallocate temp array.
1492!--          Note - Indices are shifted from nx0:nx1 to nx0+1:nx1+1 to offset the emission data
1493!--                 array to counter mentioned unkonwn offset.
1494!--                 (ecc 20190417)
1495
1496             DO  k = 1, emt_att%dt_emission
1497                DO  j = nys+1, nyn+1
1498                   DO  i = nxl+1, nxr+1
1499                      emt(ispec)%preproc_emission_data(k,1,j,i) = dum_var_5d(k,1,j-1,i-1,1)
1500                   ENDDO
1501                ENDDO
1502             ENDDO
1503
1504             DEALLOCATE( dum_var_5d )
1505
1506          ENDDO  ! ispec
1507!
1508!--       Units
1509          CALL get_attribute( id_emis, "units", emt_att%units, .FALSE. , "emission_values" )
1510
1511       ENDIF  ! LOD1 & LOD2 (default and pre-processed mode)
1512
1513       CALL close_input_file( id_emis )
1514
1515#endif
1516
1517    ENDIF ! LOD0 (parameterized mode)
1518
1519 END SUBROUTINE netcdf_data_input_chemistry_data
1520
1521
1522!--------------------------------------------------------------------------------------------------!
1523! Description:
1524! ------------
1525!> Reads surface classification data, such as vegetation and soil type, etc. .
1526!--------------------------------------------------------------------------------------------------!
1527 SUBROUTINE netcdf_data_input_surface_data
1528
1529    USE control_parameters,                                                                        &
1530        ONLY:  land_surface, urban_surface
1531
1532    USE exchange_horiz_mod,                                                                        &
1533        ONLY:  exchange_horiz_2d, exchange_horiz_2d_byte, exchange_horiz_2d_int
1534
1535    USE indices,                                                                                   &
1536        ONLY:  nbgp, nxl, nxr, nyn, nys
1537
1538
1539    IMPLICIT NONE
1540
1541    CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE ::  var_names  !< variable names in static input file
1542
1543    INTEGER(iwp) ::  id_surf   !< NetCDF id of input file
1544    INTEGER(iwp) ::  k         !< running index along z-direction
1545    INTEGER(iwp) ::  k2        !< running index
1546    INTEGER(iwp) ::  num_vars  !< number of variables in input file
1547    INTEGER(iwp) ::  nz_soil   !< number of soil layers in file
1548
1549!
1550!-- If not static input file is available, skip this routine
1551    IF ( .NOT. input_pids_static )  RETURN
1552!
1553!-- Measure CPU time
1554    CALL cpu_log( log_point_s(82), 'NetCDF input', 'start' )
1555!
1556!-- Skip the following if no land-surface or urban-surface module are applied. This case, no one of
1557!-- the following variables is used anyway.
1558    IF ( .NOT. land_surface  .AND.  .NOT. urban_surface )  RETURN
1559
1560#if defined ( __netcdf )
1561!
1562!-- Open file in read-only mode
1563    CALL open_read_file( TRIM( input_file_static ) // TRIM( coupling_char ) , id_surf )
1564!
1565!-- Inquire all variable names.
1566!-- This will be used to check whether an optional input variable exists or not.
1567    CALL inquire_num_variables( id_surf, num_vars )
1568
1569    ALLOCATE( var_names(1:num_vars) )
1570    CALL inquire_variable_names( id_surf, var_names )
1571!
1572!-- Read vegetation type and required attributes
1573    IF ( check_existence( var_names, 'vegetation_type' ) )  THEN
1574       vegetation_type_f%from_file = .TRUE.
1575       CALL get_attribute( id_surf, char_fill, vegetation_type_f%fill, .FALSE., 'vegetation_type' )
1576
1577       ALLOCATE( vegetation_type_f%var(nys:nyn,nxl:nxr) )
1578
1579       CALL get_variable( id_surf, 'vegetation_type', vegetation_type_f%var, nxl, nxr, nys, nyn )
1580    ELSE
1581       vegetation_type_f%from_file = .FALSE.
1582    ENDIF
1583
1584!
1585!-- Read soil type and required attributes
1586    IF ( check_existence( var_names, 'soil_type' ) )  THEN
1587       soil_type_f%from_file = .TRUE.
1588!
1589!--    Note, lod is currently not on file; skip for the moment
1590!       CALL get_attribute( id_surf, char_lod,                                                      &
1591!                           soil_type_f%lod,                                                        &
1592!                           .FALSE., 'soil_type' )
1593       CALL get_attribute( id_surf, char_fill, soil_type_f%fill, .FALSE., 'soil_type' )
1594
1595       IF ( soil_type_f%lod == 1 )  THEN
1596
1597          ALLOCATE( soil_type_f%var_2d(nys:nyn,nxl:nxr)  )
1598
1599          CALL get_variable( id_surf, 'soil_type', soil_type_f%var_2d, nxl, nxr, nys, nyn )
1600
1601       ELSEIF ( soil_type_f%lod == 2 )  THEN
1602!
1603!--       Obtain number of soil layers from file.
1604          CALL get_dimension_length( id_surf, nz_soil, 'zsoil' )
1605
1606          ALLOCATE( soil_type_f%var_3d(0:nz_soil,nys:nyn,nxl:nxr) )
1607
1608          CALL get_variable( id_surf, 'soil_type', soil_type_f%var_3d, nxl, nxr, nys, nyn, 0,      &
1609                             nz_soil )
1610
1611       ENDIF
1612    ELSE
1613       soil_type_f%from_file = .FALSE.
1614    ENDIF
1615
1616!
1617!-- Read pavement type and required attributes
1618    IF ( check_existence( var_names, 'pavement_type' ) )  THEN
1619       pavement_type_f%from_file = .TRUE.
1620       CALL get_attribute( id_surf, char_fill, pavement_type_f%fill, .FALSE., 'pavement_type' )
1621
1622       ALLOCATE( pavement_type_f%var(nys:nyn,nxl:nxr)  )
1623
1624       CALL get_variable( id_surf, 'pavement_type', pavement_type_f%var, nxl, nxr, nys, nyn )
1625    ELSE
1626       pavement_type_f%from_file = .FALSE.
1627    ENDIF
1628
1629!
1630!-- Read water type and required attributes
1631    IF ( check_existence( var_names, 'water_type' ) )  THEN
1632       water_type_f%from_file = .TRUE.
1633       CALL get_attribute( id_surf, char_fill, water_type_f%fill, .FALSE., 'water_type' )
1634
1635       ALLOCATE( water_type_f%var(nys:nyn,nxl:nxr)  )
1636
1637       CALL get_variable( id_surf, 'water_type', water_type_f%var, nxl, nxr, nys, nyn )
1638
1639    ELSE
1640       water_type_f%from_file = .FALSE.
1641    ENDIF
1642!
1643!-- Read relative surface fractions of vegetation, pavement and water.
1644    IF ( check_existence( var_names, 'surface_fraction' ) )  THEN
1645       surface_fraction_f%from_file = .TRUE.
1646       CALL get_attribute( id_surf, char_fill, surface_fraction_f%fill, .FALSE., 'surface_fraction')
1647!
1648!--    Inquire number of surface fractions
1649       CALL get_dimension_length( id_surf, surface_fraction_f%nf, 'nsurface_fraction' )
1650!
1651!--    Allocate dimension array and input array for surface fractions
1652       ALLOCATE( surface_fraction_f%nfracs(0:surface_fraction_f%nf-1) )
1653       ALLOCATE( surface_fraction_f%frac(0:surface_fraction_f%nf-1,nys:nyn,nxl:nxr) )
1654!
1655!--    Get dimension of surface fractions
1656       CALL get_variable( id_surf, 'nsurface_fraction', surface_fraction_f%nfracs )
1657!
1658!--    Read surface fractions
1659       CALL get_variable( id_surf, 'surface_fraction', surface_fraction_f%frac,                    &
1660                          nxl, nxr, nys, nyn, 0, surface_fraction_f%nf-1 )
1661    ELSE
1662       surface_fraction_f%from_file = .FALSE.
1663    ENDIF
1664!
1665!-- Read building parameters and related information
1666    IF ( check_existence( var_names, 'building_pars' ) )  THEN
1667       building_pars_f%from_file = .TRUE.
1668       CALL get_attribute( id_surf, char_fill, building_pars_f%fill, .FALSE., 'building_pars' )
1669!
1670!--    Inquire number of building parameters
1671       CALL get_dimension_length( id_surf, building_pars_f%np, 'nbuilding_pars' )
1672!
1673!--    Allocate dimension array and input array for building parameters
1674       ALLOCATE( building_pars_f%pars(0:building_pars_f%np-1) )
1675       ALLOCATE( building_pars_f%pars_xy(0:building_pars_f%np-1,nys:nyn,nxl:nxr) )
1676!
1677!--    Get dimension of building parameters
1678       CALL get_variable( id_surf, 'nbuilding_pars', building_pars_f%pars )
1679!
1680!--    Read building_pars
1681       CALL get_variable( id_surf, 'building_pars', building_pars_f%pars_xy,                       &
1682                          nxl, nxr, nys, nyn, 0, building_pars_f%np-1 )
1683    ELSE
1684       building_pars_f%from_file = .FALSE.
1685    ENDIF
1686
1687!
1688!-- Read building surface parameters
1689    IF ( check_existence( var_names, 'building_surface_pars' ) )  THEN
1690       building_surface_pars_f%from_file = .TRUE.
1691       CALL get_attribute( id_surf, char_fill, building_surface_pars_f%fill, .FALSE.,              &
1692                           'building_surface_pars' )
1693!
1694!--    Read building_surface_pars
1695       CALL get_variable_surf( id_surf, 'building_surface_pars', building_surface_pars_f )
1696    ELSE
1697       building_surface_pars_f%from_file = .FALSE.
1698    ENDIF
1699
1700!
1701!-- Read albedo type and required attributes
1702    IF ( check_existence( var_names, 'albedo_type' ) )  THEN
1703       albedo_type_f%from_file = .TRUE.
1704       CALL get_attribute( id_surf, char_fill, albedo_type_f%fill, .FALSE.,  'albedo_type' )
1705
1706       ALLOCATE( albedo_type_f%var(nys:nyn,nxl:nxr)  )
1707
1708       CALL get_variable( id_surf, 'albedo_type', albedo_type_f%var, nxl, nxr, nys, nyn )
1709    ELSE
1710       albedo_type_f%from_file = .FALSE.
1711    ENDIF
1712
1713!
1714!-- Read albedo parameters and related information
1715    IF ( check_existence( var_names, 'albedo_pars' ) )  THEN
1716       albedo_pars_f%from_file = .TRUE.
1717       CALL get_attribute( id_surf, char_fill, albedo_pars_f%fill, .FALSE., 'albedo_pars' )
1718!
1719!--    Inquire number of albedo parameters
1720       CALL get_dimension_length( id_surf, albedo_pars_f%np, 'nalbedo_pars' )
1721!
1722!--    Allocate dimension array and input array for albedo parameters
1723       ALLOCATE( albedo_pars_f%pars(0:albedo_pars_f%np-1) )
1724       ALLOCATE( albedo_pars_f%pars_xy(0:albedo_pars_f%np-1,nys:nyn,nxl:nxr) )
1725!
1726!--    Get dimension of albedo parameters
1727       CALL get_variable( id_surf, 'nalbedo_pars', albedo_pars_f%pars )
1728
1729       CALL get_variable( id_surf, 'albedo_pars', albedo_pars_f%pars_xy,                           &
1730                          nxl, nxr, nys, nyn, 0, albedo_pars_f%np-1 )
1731    ELSE
1732       albedo_pars_f%from_file = .FALSE.
1733    ENDIF
1734
1735!
1736!-- Read pavement parameters and related information
1737    IF ( check_existence( var_names, 'pavement_pars' ) )  THEN
1738       pavement_pars_f%from_file = .TRUE.
1739       CALL get_attribute( id_surf, char_fill, pavement_pars_f%fill, .FALSE., 'pavement_pars' )
1740!
1741!--    Inquire number of pavement parameters
1742       CALL get_dimension_length( id_surf, pavement_pars_f%np, 'npavement_pars' )
1743!
1744!--    Allocate dimension array and input array for pavement parameters
1745       ALLOCATE( pavement_pars_f%pars(0:pavement_pars_f%np-1) )
1746       ALLOCATE( pavement_pars_f%pars_xy(0:pavement_pars_f%np-1,nys:nyn,nxl:nxr) )
1747!
1748!--    Get dimension of pavement parameters
1749       CALL get_variable( id_surf, 'npavement_pars', pavement_pars_f%pars )
1750
1751       CALL get_variable( id_surf, 'pavement_pars', pavement_pars_f%pars_xy,                       &
1752                          nxl, nxr, nys, nyn, 0, pavement_pars_f%np-1 )
1753    ELSE
1754       pavement_pars_f%from_file = .FALSE.
1755    ENDIF
1756
1757!
1758!-- Read pavement subsurface parameters and related information
1759    IF ( check_existence( var_names, 'pavement_subsurface_pars' ) )  THEN
1760       pavement_subsurface_pars_f%from_file = .TRUE.
1761       CALL get_attribute( id_surf, char_fill, pavement_subsurface_pars_f%fill, .FALSE.,           &
1762                           'pavement_subsurface_pars' )
1763!
1764!--    Inquire number of parameters
1765       CALL get_dimension_length( id_surf, pavement_subsurface_pars_f%np,                          &
1766                                  'npavement_subsurface_pars' )
1767!
1768!--    Inquire number of soil layers
1769       CALL get_dimension_length( id_surf, pavement_subsurface_pars_f%nz, 'zsoil' )
1770!
1771!--    Allocate dimension array and input array for pavement parameters
1772       ALLOCATE( pavement_subsurface_pars_f%pars(0:pavement_subsurface_pars_f%np-1) )
1773       ALLOCATE( pavement_subsurface_pars_f%pars_xyz(0:pavement_subsurface_pars_f%np-1,            &
1774                                                     0:pavement_subsurface_pars_f%nz-1,            &
1775                                                     nys:nyn,nxl:nxr) )
1776!
1777!--    Get dimension of pavement parameters
1778       CALL get_variable( id_surf, 'npavement_subsurface_pars', pavement_subsurface_pars_f%pars )
1779
1780       CALL get_variable( id_surf, 'pavement_subsurface_pars', pavement_subsurface_pars_f%pars_xyz,&
1781                          nxl, nxr, nys, nyn, 0, pavement_subsurface_pars_f%nz-1,                  &
1782                          0, pavement_subsurface_pars_f%np-1 )
1783    ELSE
1784       pavement_subsurface_pars_f%from_file = .FALSE.
1785    ENDIF
1786
1787
1788!
1789!-- Read vegetation parameters and related information
1790    IF ( check_existence( var_names, 'vegetation_pars' ) )  THEN
1791       vegetation_pars_f%from_file = .TRUE.
1792       CALL get_attribute( id_surf, char_fill, vegetation_pars_f%fill, .FALSE., 'vegetation_pars' )
1793!
1794!--    Inquire number of vegetation parameters
1795       CALL get_dimension_length( id_surf, vegetation_pars_f%np, 'nvegetation_pars' )
1796!
1797!--    Allocate dimension array and input array for surface fractions
1798       ALLOCATE( vegetation_pars_f%pars(0:vegetation_pars_f%np-1) )
1799       ALLOCATE( vegetation_pars_f%pars_xy(0:vegetation_pars_f%np-1,nys:nyn,nxl:nxr) )
1800!
1801!--    Get dimension of the parameters
1802       CALL get_variable( id_surf, 'nvegetation_pars', vegetation_pars_f%pars )
1803       CALL get_variable( id_surf, 'vegetation_pars', vegetation_pars_f%pars_xy,                   &
1804                          nxl, nxr, nys, nyn, 0, vegetation_pars_f%np-1 )
1805    ELSE
1806       vegetation_pars_f%from_file = .FALSE.
1807    ENDIF
1808
1809!
1810!-- Read root parameters/distribution and related information
1811    IF ( check_existence( var_names, 'soil_pars' ) )  THEN
1812       soil_pars_f%from_file = .TRUE.
1813       CALL get_attribute( id_surf, char_fill, soil_pars_f%fill, .FALSE., 'soil_pars' )
1814       CALL get_attribute( id_surf, char_lod, soil_pars_f%lod, .FALSE., 'soil_pars' )
1815
1816!
1817!--    Inquire number of soil parameters
1818       CALL get_dimension_length( id_surf, soil_pars_f%np, 'nsoil_pars' )
1819!
1820!--    Read parameters array
1821       ALLOCATE( soil_pars_f%pars(0:soil_pars_f%np-1) )
1822       CALL get_variable( id_surf, 'nsoil_pars', soil_pars_f%pars )
1823
1824!
1825!--    In case of level of detail 2, also inquire number of vertical soil layers, allocate memory
1826!--    and read the respective dimension.
1827       IF ( soil_pars_f%lod == 2 )  THEN
1828          CALL get_dimension_length( id_surf, soil_pars_f%nz, 'zsoil' )
1829
1830          ALLOCATE( soil_pars_f%layers(0:soil_pars_f%nz-1) )
1831          CALL get_variable( id_surf, 'zsoil', soil_pars_f%layers )
1832
1833       ENDIF
1834
1835!
1836!--    Read soil parameters, depending on level of detail
1837       IF ( soil_pars_f%lod == 1 )  THEN
1838          ALLOCATE( soil_pars_f%pars_xy(0:soil_pars_f%np-1,nys:nyn,nxl:nxr) )
1839
1840          CALL get_variable( id_surf, 'soil_pars', soil_pars_f%pars_xy,                            &
1841                             nxl, nxr, nys, nyn, 0, soil_pars_f%np-1 )
1842
1843       ELSEIF ( soil_pars_f%lod == 2 )  THEN
1844          ALLOCATE( soil_pars_f%pars_xyz(0:soil_pars_f%np-1,0:soil_pars_f%nz-1,nys:nyn,nxl:nxr) )
1845          CALL get_variable( id_surf, 'soil_pars', soil_pars_f%pars_xyz,                           &
1846                             nxl, nxr, nys, nyn, 0, soil_pars_f%nz-1, 0, soil_pars_f%np-1 )
1847       ENDIF
1848    ELSE
1849       soil_pars_f%from_file = .FALSE.
1850    ENDIF
1851
1852!
1853!-- Read water parameters and related information
1854    IF ( check_existence( var_names, 'water_pars' ) )  THEN
1855       water_pars_f%from_file = .TRUE.
1856       CALL get_attribute( id_surf, char_fill, water_pars_f%fill, .FALSE., 'water_pars' )
1857!
1858!--    Inquire number of water parameters
1859       CALL get_dimension_length( id_surf, water_pars_f%np, 'nwater_pars' )
1860!
1861!--    Allocate dimension array and input array for water parameters
1862       ALLOCATE( water_pars_f%pars(0:water_pars_f%np-1) )
1863       ALLOCATE( water_pars_f%pars_xy(0:water_pars_f%np-1,nys:nyn,nxl:nxr) )
1864!
1865!--    Get dimension of water parameters
1866       CALL get_variable( id_surf, 'nwater_pars', water_pars_f%pars )
1867       CALL get_variable( id_surf, 'water_pars', water_pars_f%pars_xy,                             &
1868                          nxl, nxr, nys, nyn, 0, water_pars_f%np-1 )
1869    ELSE
1870       water_pars_f%from_file = .FALSE.
1871    ENDIF
1872!
1873!-- Read root area density - parametrized vegetation
1874    IF ( check_existence( var_names, 'root_area_dens_s' ) )  THEN
1875       root_area_density_lsm_f%from_file = .TRUE.
1876       CALL get_attribute( id_surf, char_fill, root_area_density_lsm_f%fill, .FALSE.,              &
1877                           'root_area_dens_s' )
1878!
1879!--    Obtain number of soil layers from file and allocate variable
1880       CALL get_dimension_length( id_surf, root_area_density_lsm_f%nz, 'zsoil' )
1881       ALLOCATE( root_area_density_lsm_f%var(0:root_area_density_lsm_f%nz-1,nys:nyn,nxl:nxr) )
1882
1883!
1884!--    Read root-area density
1885       CALL get_variable( id_surf, 'root_area_dens_s', root_area_density_lsm_f%var,                &
1886                          nxl, nxr, nys, nyn, 0, root_area_density_lsm_f%nz-1 )
1887
1888    ELSE
1889       root_area_density_lsm_f%from_file = .FALSE.
1890    ENDIF
1891!
1892!-- Read street type and street crossing
1893    IF ( check_existence( var_names, 'street_type' ) )  THEN
1894       street_type_f%from_file = .TRUE.
1895       CALL get_attribute( id_surf, char_fill, street_type_f%fill, .FALSE., 'street_type' )
1896       ALLOCATE( street_type_f%var(nys:nyn,nxl:nxr) )
1897       CALL get_variable( id_surf, 'street_type', street_type_f%var, nxl, nxr, nys, nyn )
1898    ELSE
1899       street_type_f%from_file = .FALSE.
1900    ENDIF
1901
1902    IF ( check_existence( var_names, 'street_crossing' ) )  THEN
1903       street_crossing_f%from_file = .TRUE.
1904       CALL get_attribute( id_surf, char_fill, street_crossing_f%fill, .FALSE., 'street_crossing' )
1905       ALLOCATE( street_crossing_f%var(nys:nyn,nxl:nxr) )
1906       CALL get_variable( id_surf, 'street_crossing', street_crossing_f%var, nxl, nxr, nys, nyn )
1907    ELSE
1908       street_crossing_f%from_file = .FALSE.
1909    ENDIF
1910
1911!
1912!-- Still missing: root_resolved and building_surface_pars.
1913!-- Will be implemented as soon as they are available.
1914
1915!
1916!-- Finally, close input file
1917    CALL close_input_file( id_surf )
1918#endif
1919!
1920!-- End of CPU measurement
1921    CALL cpu_log( log_point_s(82), 'NetCDF input', 'stop' )
1922
1923!
1924!-- Exchange ghost points for surface variables. Therefore, resize variables.
1925    IF ( albedo_type_f%from_file )  THEN
1926       CALL add_ghost_layers( albedo_type_f%var )
1927       CALL exchange_horiz_2d_byte( albedo_type_f%var, nys, nyn, nxl, nxr, nbgp )
1928    ENDIF
1929    IF ( pavement_type_f%from_file )  THEN
1930       CALL add_ghost_layers( pavement_type_f%var )
1931       CALL exchange_horiz_2d_byte( pavement_type_f%var, nys, nyn, nxl, nxr, nbgp )
1932    ENDIF
1933    IF ( soil_type_f%from_file  .AND.  ALLOCATED( soil_type_f%var_2d ) )  THEN
1934       CALL add_ghost_layers( soil_type_f%var_2d )
1935       CALL exchange_horiz_2d_byte( soil_type_f%var_2d, nys, nyn, nxl, nxr, nbgp )
1936    ENDIF
1937    IF ( vegetation_type_f%from_file )  THEN
1938       CALL add_ghost_layers( vegetation_type_f%var )
1939       CALL exchange_horiz_2d_byte( vegetation_type_f%var, nys, nyn, nxl, nxr, nbgp )
1940    ENDIF
1941    IF ( water_type_f%from_file )  THEN
1942       CALL add_ghost_layers( water_type_f%var )
1943       CALL exchange_horiz_2d_byte( water_type_f%var, nys, nyn, nxl, nxr, nbgp )
1944    ENDIF
1945
1946!
1947!-- Exchange ghost points for 3/4-D variables. For the sake of simplicity, loop further dimensions
1948!-- to use 2D exchange routines. Unfortunately this is necessary, else new MPI-data types need to be
1949!-- introduced just for 2 variables.
1950    IF ( soil_type_f%from_file  .AND.  ALLOCATED( soil_type_f%var_3d ) )  THEN
1951       CALL add_ghost_layers( soil_type_f%var_3d, 0, nz_soil )
1952       DO  k = 0, nz_soil
1953          CALL exchange_horiz_2d_byte( soil_type_f%var_3d(k,:,:), nys, nyn, nxl, nxr, nbgp )
1954       ENDDO
1955    ENDIF
1956
1957    IF ( surface_fraction_f%from_file )  THEN
1958       CALL add_ghost_layers( surface_fraction_f%frac, 0, surface_fraction_f%nf-1 )
1959       DO  k = 0, surface_fraction_f%nf-1
1960          CALL exchange_horiz_2d( surface_fraction_f%frac(k,:,:) )
1961       ENDDO
1962    ENDIF
1963
1964    IF ( building_pars_f%from_file )  THEN
1965       CALL add_ghost_layers( building_pars_f%pars_xy, 0, building_pars_f%np-1 )
1966       DO  k = 0, building_pars_f%np-1
1967          CALL exchange_horiz_2d( building_pars_f%pars_xy(k,:,:) )
1968       ENDDO
1969    ENDIF
1970
1971    IF ( albedo_pars_f%from_file )  THEN
1972       CALL add_ghost_layers( albedo_pars_f%pars_xy, 0, albedo_pars_f%np-1 )
1973       DO  k = 0, albedo_pars_f%np-1
1974          CALL exchange_horiz_2d( albedo_pars_f%pars_xy(k,:,:) )
1975       ENDDO
1976    ENDIF
1977
1978    IF ( pavement_pars_f%from_file )  THEN
1979       CALL add_ghost_layers( pavement_pars_f%pars_xy, 0, pavement_pars_f%np-1 )
1980       DO  k = 0, pavement_pars_f%np-1
1981          CALL exchange_horiz_2d( pavement_pars_f%pars_xy(k,:,:) )
1982       ENDDO
1983    ENDIF
1984
1985    IF ( vegetation_pars_f%from_file )  THEN
1986       CALL add_ghost_layers( vegetation_pars_f%pars_xy, 0, vegetation_pars_f%np-1 )
1987       DO  k = 0, vegetation_pars_f%np-1
1988          CALL exchange_horiz_2d( vegetation_pars_f%pars_xy(k,:,:) )
1989       ENDDO
1990    ENDIF
1991
1992    IF ( water_pars_f%from_file )  THEN
1993       CALL add_ghost_layers( water_pars_f%pars_xy, 0, water_pars_f%np-1 )
1994       DO  k = 0, water_pars_f%np-1
1995          CALL exchange_horiz_2d( water_pars_f%pars_xy(k,:,:) )
1996       ENDDO
1997    ENDIF
1998
1999    IF ( root_area_density_lsm_f%from_file )  THEN
2000       CALL add_ghost_layers( root_area_density_lsm_f%var, 0, root_area_density_lsm_f%nz-1 )
2001       DO  k = 0, root_area_density_lsm_f%nz-1
2002          CALL exchange_horiz_2d( root_area_density_lsm_f%var(k,:,:) )
2003       ENDDO
2004    ENDIF
2005
2006    IF ( soil_pars_f%from_file )  THEN
2007       IF ( soil_pars_f%lod == 1 )  THEN
2008          CALL add_ghost_layers( soil_pars_f%pars_xy, 0, soil_pars_f%np-1 )
2009          DO  k = 0, soil_pars_f%np-1
2010             CALL exchange_horiz_2d( soil_pars_f%pars_xy(k,:,:) )
2011          ENDDO
2012
2013       ELSEIF ( soil_pars_f%lod == 2 )  THEN
2014          CALL add_ghost_layers( soil_pars_f%pars_xyz, 0, soil_pars_f%np-1, 0, soil_pars_f%nz-1 )
2015
2016          DO  k2 = 0, soil_pars_f%nz-1
2017             DO  k = 0, soil_pars_f%np-1
2018                CALL exchange_horiz_2d( soil_pars_f%pars_xyz(k,k2,:,:) )
2019             ENDDO
2020          ENDDO
2021       ENDIF
2022    ENDIF
2023
2024    IF ( pavement_subsurface_pars_f%from_file )  THEN
2025       CALL add_ghost_layers( pavement_subsurface_pars_f%pars_xyz,                                 &
2026                              0, pavement_subsurface_pars_f%np-1,                                  &
2027                              0, pavement_subsurface_pars_f%nz-1 )
2028
2029       DO  k2 = 0, pavement_subsurface_pars_f%nz-1
2030          DO  k = 0, pavement_subsurface_pars_f%np-1
2031             CALL exchange_horiz_2d( pavement_subsurface_pars_f%pars_xyz(k,k2,:,:) )
2032          ENDDO
2033       ENDDO
2034    ENDIF
2035
2036 END SUBROUTINE netcdf_data_input_surface_data
2037
2038
2039!--------------------------------------------------------------------------------------------------!
2040! Description:
2041! ------------
2042!> Reads uvem lookup table information.
2043!--------------------------------------------------------------------------------------------------!
2044 SUBROUTINE netcdf_data_input_uvem
2045
2046    USE indices,                                                                                   &
2047        ONLY:  nxl, nxr, nyn, nys
2048
2049    IMPLICIT NONE
2050
2051    CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE ::  var_names  !< variable names in static input file
2052
2053
2054    INTEGER(iwp) ::  id_uvem       !< NetCDF id of uvem lookup table input file
2055    INTEGER(iwp) ::  nli = 35      !< dimension length of lookup table in x
2056    INTEGER(iwp) ::  nlj =  9      !< dimension length of lookup table in y
2057    INTEGER(iwp) ::  nlk = 90      !< dimension length of lookup table in z
2058    INTEGER(iwp) ::  num_vars      !< number of variables in netcdf input file
2059!
2060!-- Input via uv exposure model lookup table input
2061    IF ( input_pids_uvem )  THEN
2062
2063#if defined ( __netcdf )
2064!
2065!--    Open file in read-only mode
2066       CALL open_read_file( TRIM( input_file_uvem ) // TRIM( coupling_char ), id_uvem )
2067!
2068!--    At first, inquire all variable names.
2069!--    This will be used to check whether an input variable exist or not.
2070       CALL inquire_num_variables( id_uvem, num_vars )
2071!
2072!--    Allocate memory to store variable names and inquire them.
2073       ALLOCATE( var_names(1:num_vars) )
2074       CALL inquire_variable_names( id_uvem, var_names )
2075!
2076!--    uvem integration
2077       IF ( check_existence( var_names, 'int_factors' ) )  THEN
2078          uvem_integration_f%from_file = .TRUE.
2079!
2080!--       Input 2D uvem integration.
2081          ALLOCATE( uvem_integration_f%var(0:nlj,0:nli)  )
2082          CALL get_variable( id_uvem, 'int_factors', uvem_integration_f%var, 0, nli, 0, nlj )
2083       ELSE
2084          uvem_integration_f%from_file = .FALSE.
2085       ENDIF
2086!
2087!--    uvem irradiance
2088       IF ( check_existence( var_names, 'irradiance' ) )  THEN
2089          uvem_irradiance_f%from_file = .TRUE.
2090!
2091!--       Input 2D uvem irradiance.
2092          ALLOCATE( uvem_irradiance_f%var(0:nlk, 0:2)  )
2093          CALL get_variable( id_uvem, 'irradiance', uvem_irradiance_f%var, 0, 2, 0, nlk )
2094       ELSE
2095          uvem_irradiance_f%from_file = .FALSE.
2096       ENDIF
2097!
2098!--    uvem porjection areas
2099       IF ( check_existence( var_names, 'projarea' ) )  THEN
2100          uvem_projarea_f%from_file = .TRUE.
2101!
2102!--       Input 3D uvem projection area (human geometgry)
2103          ALLOCATE( uvem_projarea_f%var(0:2,0:nlj,0:nli)  )
2104          CALL get_variable( id_uvem, 'projarea', uvem_projarea_f%var, 0, nli, 0, nlj, 0, 2 )
2105       ELSE
2106          uvem_projarea_f%from_file = .FALSE.
2107       ENDIF
2108!
2109!--    uvem radiance
2110       IF ( check_existence( var_names, 'radiance' ) )  THEN
2111          uvem_radiance_f%from_file = .TRUE.
2112!
2113!--       Input 3D uvem radiance
2114          ALLOCATE( uvem_radiance_f%var(0:nlk,0:nlj,0:nli)  )
2115          CALL get_variable( id_uvem, 'radiance', uvem_radiance_f%var, 0, nli, 0, nlj, 0, nlk )
2116       ELSE
2117          uvem_radiance_f%from_file = .FALSE.
2118       ENDIF
2119!
2120!--    Read building obstruction
2121!      @bug This part is deactivated due to improper implementation: conflicts with
2122!           building_obstruction_f and it is also not used anywhere else in the code
2123!        IF ( check_existence( var_names, 'obstruction' ) )  THEN
2124!           building_obstruction_full%from_file = .TRUE.
2125! !
2126! !--       Input 3D uvem building obstruction
2127!           ALLOCATE( building_obstruction_full%var_3d(0:44,0:2,0:2) )
2128!           CALL get_variable( id_uvem, 'obstruction', building_obstruction_full%var_3d, &
2129!                              0, 2, 0, 2, 0, 44 )
2130!        ELSE
2131!           building_obstruction_full%from_file = .FALSE.
2132!        ENDIF
2133!
2134       IF ( check_existence( var_names, 'obstruction' ) )  THEN
2135          building_obstruction_f%from_file = .TRUE.
2136!
2137!--       Input 3D uvem building obstruction
2138          ALLOCATE( building_obstruction_f%var_3d(0:44,nys:nyn,nxl:nxr) )
2139          CALL get_variable( id_uvem, 'obstruction', building_obstruction_f%var_3d,                &
2140                             nxl, nxr, nys, nyn, 0, 44 )
2141       ELSE
2142          building_obstruction_f%from_file = .FALSE.
2143       ENDIF
2144!
2145!--    Close uvem lookup table input file
2146       CALL close_input_file( id_uvem )
2147#else
2148       CONTINUE
2149#endif
2150    ENDIF
2151
2152 END SUBROUTINE netcdf_data_input_uvem
2153
2154
2155!--------------------------------------------------------------------------------------------------!
2156! Description:
2157! ------------
2158!> Reads orography and building information.
2159!--------------------------------------------------------------------------------------------------!
2160 SUBROUTINE netcdf_data_input_topo
2161
2162    USE control_parameters,                                                                        &
2163        ONLY:  message_string, topography
2164
2165    USE exchange_horiz_mod,                                                                        &
2166        ONLY:  exchange_horiz_2d_byte, exchange_horiz_2d_int
2167
2168    USE grid_variables,                                                                            &
2169        ONLY:  dx, dy
2170
2171    USE indices,                                                                                   &
2172        ONLY:  nbgp, nx, nxl, nxr, ny, nyn, nys, nzb
2173
2174
2175    IMPLICIT NONE
2176
2177    CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE ::  var_names  !< variable names in static input file
2178
2179
2180    INTEGER(iwp) ::  i             !< running index along x-direction
2181    INTEGER(iwp) ::  id_topo       !< NetCDF id of topograhy input file
2182    INTEGER(iwp) ::  ii            !< running index for IO blocks
2183    INTEGER(iwp) ::  io_status     !< status after reading the ascii topo file
2184    INTEGER(iwp) ::  j             !< running index along y-direction
2185    INTEGER(iwp) ::  num_vars      !< number of variables in netcdf input file
2186    INTEGER(iwp) ::  skip_n_rows   !< counting variable to skip rows while reading topography file
2187
2188    REAL(wp) ::  dum           !< dummy variable to skip columns while reading topography file
2189
2190
2191!
2192!-- CPU measurement
2193    CALL cpu_log( log_point_s(83), 'NetCDF/ASCII input topo', 'start' )
2194!
2195!-- Input via palm-input data standard
2196    IF ( input_pids_static )  THEN
2197#if defined ( __netcdf )
2198!
2199!--    Open file in read-only mode
2200       CALL open_read_file( TRIM( input_file_static ) // TRIM( coupling_char ), id_topo )
2201!
2202!--    At first, inquire all variable names.
2203!--    This will be used to check whether an input variable exists or not.
2204       CALL inquire_num_variables( id_topo, num_vars )
2205!
2206!--    Allocate memory to store variable names and inquire them.
2207       ALLOCATE( var_names(1:num_vars) )
2208       CALL inquire_variable_names( id_topo, var_names )
2209!
2210!--    Read x, y - dimensions. Only required for consistency checks.
2211       CALL get_dimension_length( id_topo, dim_static%nx, 'x' )
2212       CALL get_dimension_length( id_topo, dim_static%ny, 'y' )
2213       ALLOCATE( dim_static%x(0:dim_static%nx-1) )
2214       ALLOCATE( dim_static%y(0:dim_static%ny-1) )
2215       CALL get_variable( id_topo, 'x', dim_static%x )
2216       CALL get_variable( id_topo, 'y', dim_static%y )
2217!
2218!--    Check whether dimension size in input file matches the model dimensions
2219       IF ( dim_static%nx-1 /= nx  .OR.  dim_static%ny-1 /= ny )  THEN
2220          message_string = 'Static input file: horizontal dimension in ' //                        &
2221                           'x- and/or y-direction ' //                                             &
2222                           'do not match the respective model dimension'
2223          CALL message( 'netcdf_data_input_mod', 'PA0548', 1, 2, 0, 6, 0 )
2224       ENDIF
2225!
2226!--    Check if grid spacing of provided input data matches the respective grid spacing in the
2227!--    model.
2228       IF ( ABS( dim_static%x(1) - dim_static%x(0) - dx ) > 10E-6_wp  .OR.                         &
2229            ABS( dim_static%y(1) - dim_static%y(0) - dy ) > 10E-6_wp )  THEN
2230          message_string = 'Static input file: horizontal grid spacing ' //                        &
2231                           'in x- and/or y-direction ' //                                          &
2232                           'do not match the respective model grid spacing.'
2233          CALL message( 'netcdf_data_input_mod', 'PA0549', 1, 2, 0, 6, 0 )
2234       ENDIF
2235!
2236!--    Terrain height. First, get variable-related _FillValue attribute
2237       IF ( check_existence( var_names, 'zt' ) )  THEN
2238          terrain_height_f%from_file = .TRUE.
2239          CALL get_attribute( id_topo, char_fill, terrain_height_f%fill, .FALSE., 'zt' )
2240!
2241!--       Input 2D terrain height.
2242          ALLOCATE( terrain_height_f%var(nys:nyn,nxl:nxr)  )
2243
2244          CALL get_variable( id_topo, 'zt', terrain_height_f%var, nxl, nxr, nys, nyn )
2245
2246       ELSE
2247          terrain_height_f%from_file = .FALSE.
2248       ENDIF
2249
2250!
2251!--    Read building height. First, read its _FillValue attribute, as well as lod attribute
2252       buildings_f%from_file = .FALSE.
2253       IF ( check_existence( var_names, 'buildings_2d' ) )  THEN
2254          buildings_f%from_file = .TRUE.
2255          CALL get_attribute( id_topo, char_lod, buildings_f%lod, .FALSE., 'buildings_2d' )
2256          CALL get_attribute( id_topo, char_fill, buildings_f%fill1, .FALSE., 'buildings_2d' )
2257
2258!
2259!--       Read 2D buildings
2260          IF ( buildings_f%lod == 1 )  THEN
2261             ALLOCATE( buildings_f%var_2d(nys:nyn,nxl:nxr) )
2262             CALL get_variable( id_topo, 'buildings_2d',                                           &
2263                                buildings_f%var_2d,                                                &
2264                                nxl, nxr, nys, nyn )
2265          ELSE
2266             message_string = 'NetCDF attribute lod ' //                                           &
2267                              '(level of detail) is not set ' //                                   &
2268                              'properly for buildings_2d.'
2269             CALL message( 'netcdf_data_input_mod', 'PA0540', 1, 2, 0, 6, 0 )
2270          ENDIF
2271       ENDIF
2272!
2273!--    If available, also read 3D building information. If both are available, use 3D information.
2274       IF ( check_existence( var_names, 'buildings_3d' ) )  THEN
2275          buildings_f%from_file = .TRUE.
2276          CALL get_attribute( id_topo, char_lod, buildings_f%lod, .FALSE., 'buildings_3d' )
2277          CALL get_attribute( id_topo, char_fill, buildings_f%fill2, .FALSE., 'buildings_3d' )
2278          CALL get_dimension_length( id_topo, buildings_f%nz, 'z' )
2279!
2280!--       Read 3D buildings
2281          IF ( buildings_f%lod == 2 )  THEN
2282             ALLOCATE( buildings_f%z(nzb:buildings_f%nz-1) )
2283             CALL get_variable( id_topo, 'z', buildings_f%z )
2284             ALLOCATE( buildings_f%var_3d(nzb:buildings_f%nz-1, nys:nyn,nxl:nxr) )
2285             buildings_f%var_3d = 0
2286             CALL get_variable( id_topo, 'buildings_3d', buildings_f%var_3d, nxl, nxr, nys, nyn, 0,&
2287                                buildings_f%nz-1 )
2288          ELSE
2289             message_string = 'NetCDF attribute lod (level of detail) is not set ' //              &
2290                              'properly for buildings_3d.'
2291             CALL message( 'netcdf_data_input_mod', 'PA0541', 1, 2, 0, 6, 0 )
2292          ENDIF
2293       ENDIF
2294!
2295!--    Read building IDs and its FillValue attribute. Further required for mapping buildings on top
2296!--    of orography.
2297       IF ( check_existence( var_names, 'building_id' ) )  THEN
2298          building_id_f%from_file = .TRUE.
2299          CALL get_attribute( id_topo, char_fill, building_id_f%fill, .FALSE., 'building_id' )
2300          ALLOCATE( building_id_f%var(nys:nyn,nxl:nxr) )
2301          CALL get_variable( id_topo, 'building_id', building_id_f%var, nxl, nxr, nys, nyn )
2302       ELSE
2303          building_id_f%from_file = .FALSE.
2304       ENDIF
2305!
2306!--    Read building_type and required attributes.
2307       IF ( check_existence( var_names, 'building_type' ) )  THEN
2308          building_type_f%from_file = .TRUE.
2309          CALL get_attribute( id_topo, char_fill, building_type_f%fill, .FALSE., 'building_type' )
2310          ALLOCATE( building_type_f%var(nys:nyn,nxl:nxr) )
2311          CALL get_variable( id_topo, 'building_type', building_type_f%var, nxl, nxr, nys, nyn )
2312       ELSE
2313          building_type_f%from_file = .FALSE.
2314       ENDIF
2315!
2316!--    Close topography input file
2317       CALL close_input_file( id_topo )
2318#else
2319       CONTINUE
2320#endif
2321!
2322!-- ASCII input
2323    ELSEIF ( TRIM( topography ) == 'read_from_file' )  THEN
2324
2325       DO  ii = 0, io_blocks-1
2326          IF ( ii == io_group )  THEN
2327
2328             OPEN( 90, FILE='TOPOGRAPHY_DATA' // TRIM( coupling_char ), STATUS='OLD',              &
2329                       FORM='FORMATTED', IOSTAT=io_status )
2330
2331             IF ( io_status > 0 )  THEN
2332                message_string = 'file TOPOGRAPHY_DATA' //                                         &
2333                                 TRIM( coupling_char ) // ' does not exist'
2334                CALL message( 'netcdf_data_input_mod', 'PA0208', 1, 2, 0, 6, 0 )
2335             ENDIF
2336
2337!
2338!--          Read topography PE-wise. Rows are read from nyn to nys, columns are read from nxl to
2339!--          nxr. At first, ny-nyn rows need to be skipped.
2340             skip_n_rows = 0
2341             DO  WHILE ( skip_n_rows < ny - nyn )
2342                READ( 90, * )
2343                skip_n_rows = skip_n_rows + 1
2344             ENDDO
2345!
2346!--          Read data from nyn to nys and nxl to nxr. Therefore, skip column until nxl-1 is reached
2347             ALLOCATE( buildings_f%var_2d(nys:nyn,nxl:nxr) )
2348             DO  j = nyn, nys, -1
2349
2350                READ( 90, *, IOSTAT=io_status )  ( dum, i = 0, nxl-1 ),                            &
2351                                                 ( buildings_f%var_2d(j,i), i = nxl, nxr )
2352
2353                IF ( io_status > 0 )  THEN
2354                   WRITE( message_string, '(A,1X,I5,1X,A)' ) 'error reading line', ny-j+1,         &
2355                                          'of file TOPOGRAPHY_DATA' // TRIM( coupling_char )
2356                   CALL message( 'netcdf_data_input_mod', 'PA0209', 2, 2, myid, 6, 0 )
2357                ELSEIF ( io_status < 0 )  THEN
2358                   WRITE( message_string, '(A,1X,I5)' ) 'end of line or file detected for '//      &
2359                               'file TOPOGRAPHY_DATA' // TRIM( coupling_char ) // ' at line', ny-j+1
2360                   CALL message( 'netcdf_data_input_mod', 'PA0704', 2, 2, myid, 6, 0 )
2361                ENDIF
2362
2363             ENDDO
2364
2365             CLOSE( 90 )
2366             buildings_f%from_file = .TRUE.
2367
2368          ENDIF
2369#if defined( __parallel )
2370          CALL MPI_BARRIER( comm2d, ierr )
2371#endif
2372       ENDDO
2373
2374    ENDIF
2375!
2376!-- End of CPU measurement
2377    CALL cpu_log( log_point_s(83), 'NetCDF/ASCII input topo', 'stop' )
2378
2379!
2380!-- Check for minimum requirement to setup building topography. If buildings are provided, also an
2381!-- ID and a type are required.
2382!-- Note, doing this check in check_parameters will be too late (data will be used for grid
2383!-- inititialization before).
2384    IF ( input_pids_static )  THEN
2385       IF ( buildings_f%from_file  .AND.  .NOT. building_id_f%from_file )  THEN
2386          message_string = 'If building heights are prescribed in ' //                             &
2387                           'static input file, also an ID is required.'
2388          CALL message( 'netcdf_data_input_mod', 'PA0542', 1, 2, 0, 6, 0 )
2389       ENDIF
2390    ENDIF
2391
2392    IF ( terrain_height_f%from_file )  THEN
2393!
2394!--    Check orography for fill-values.
2395!--    For the moment, give an error message. More advanced methods, e.g. a nearest neighbor
2396!--    algorithm as used in GIS systems might be implemented later.
2397!--    Note: This check must be placed here as terrain_height_f is altered within init_grid which is
2398!--    called before netcdf_data_input_check_static
2399       IF ( ANY( terrain_height_f%var == terrain_height_f%fill ) )  THEN
2400          message_string = 'NetCDF variable zt is not ' // 'allowed to have missing data'
2401          CALL message( 'netcdf_data_input_mod', 'PA0550', 2, 2, myid, 6, 0 )
2402       ENDIF
2403    ELSE
2404!
2405!--    In case no terrain height is provided by static input file, allocate array nevertheless and
2406!--    set terrain height to 0, which simplifies topography initialization.
2407       ALLOCATE( terrain_height_f%var(nys:nyn,nxl:nxr) )
2408       terrain_height_f%var = 0.0_wp
2409    ENDIF
2410!
2411!-- Finally, exchange 1 ghost point for building ID and type.
2412!-- In case of non-cyclic boundary conditions set Neumann conditions at the lateral boundaries.
2413    IF ( building_id_f%from_file )  THEN
2414       CALL add_ghost_layers( building_id_f%var )
2415       CALL exchange_horiz_2d_int( building_id_f%var, nys, nyn, nxl, nxr, nbgp )
2416    ENDIF
2417
2418    IF ( building_type_f%from_file )  THEN
2419       CALL add_ghost_layers( building_type_f%var )
2420       CALL exchange_horiz_2d_byte( building_type_f%var, nys, nyn, nxl, nxr, nbgp )
2421    ENDIF
2422
2423 END SUBROUTINE netcdf_data_input_topo
2424
2425
2426!--------------------------------------------------------------------------------------------------!
2427! Description:
2428! ------------
2429!> Reads initialization data of u, v, w, pt, q, geostrophic wind components, as well as soil
2430!> moisture and soil temperature, derived from larger-scale model (COSMO) by Inifor.
2431!--------------------------------------------------------------------------------------------------!
2432 SUBROUTINE netcdf_data_input_init_3d
2433
2434    USE arrays_3d,                                                                                 &
2435        ONLY:  pt, q, u, v, w, zu, zw
2436
2437    USE control_parameters,                                                                        &
2438        ONLY:  air_chemistry, bc_lr_cyc, bc_ns_cyc, humidity, message_string, neutral
2439
2440    USE indices,                                                                                   &
2441        ONLY:  nx, nxl, nxlu, nxr, ny, nyn, nys, nysv, nzb, nz, nzt
2442
2443    IMPLICIT NONE
2444
2445    CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE ::  var_names
2446
2447    INTEGER(iwp) ::  id_dynamic  !< NetCDF id of dynamic input file
2448    INTEGER(iwp) ::  n           !< running index for chemistry variables
2449    INTEGER(iwp) ::  num_vars    !< number of variables in netcdf input file
2450
2451    LOGICAL      ::  check_passed         !< flag indicating if a check passed
2452    LOGICAL      ::  dynamic_3d = .TRUE.  !< flag indicating that 3D data is read from dynamic file
2453
2454!
2455!-- Skip routine if no input file with dynamic input data is available.
2456    IF ( .NOT. input_pids_dynamic )  RETURN
2457!
2458!-- Please note, Inifor is designed to provide initial data for u and v for the prognostic grid
2459!-- points in case of lateral Dirichlet conditions.
2460!-- This means that Inifor provides data from nxlu:nxr (for u) and from nysv:nyn (for v) at the left
2461!-- and south domain boundary, respectively.
2462!-- However, as work-around for the moment, PALM will run with cyclic conditions and will be
2463!-- initialized with data provided by Inifor boundaries in case of Dirichlet.
2464!-- Hence, simply set set nxlu/nysv to 1 (will be reset to its original value at the end of this
2465!-- routine).
2466    IF ( bc_lr_cyc  .AND.  nxl == 0 )  nxlu = 1
2467    IF ( bc_ns_cyc  .AND.  nys == 0 )  nysv = 1
2468
2469!
2470!-- CPU measurement
2471    CALL cpu_log( log_point_s(85), 'NetCDF input init', 'start' )
2472
2473#if defined ( __netcdf )
2474!
2475!-- Open file in read-only mode
2476    CALL open_read_file( TRIM( input_file_dynamic ) // TRIM( coupling_char ), id_dynamic )
2477
2478!
2479!-- At first, inquire all variable names.
2480    CALL inquire_num_variables( id_dynamic, num_vars )
2481!
2482!-- Allocate memory to store variable names.
2483    ALLOCATE( var_names(1:num_vars) )
2484    CALL inquire_variable_names( id_dynamic, var_names )
2485!
2486!-- Read vertical dimension of scalar und w grid.
2487    CALL get_dimension_length( id_dynamic, init_3d%nzu, 'z'     )
2488    CALL get_dimension_length( id_dynamic, init_3d%nzw, 'zw'    )
2489!
2490!-- Read also the horizontal dimensions. These are used just used for checking the compatibility
2491!-- with the PALM grid before reading.
2492    CALL get_dimension_length( id_dynamic, init_3d%nx,  'x'  )
2493    CALL get_dimension_length( id_dynamic, init_3d%nxu, 'xu' )
2494    CALL get_dimension_length( id_dynamic, init_3d%ny,  'y'  )
2495    CALL get_dimension_length( id_dynamic, init_3d%nyv, 'yv' )
2496
2497!
2498!-- Check for correct horizontal and vertical dimension. Please note, checks are performed directly
2499!-- here and not called from check_parameters as some varialbes are still not allocated there.
2500!-- Moreover, please note, u- and v-grid has 1 grid point less on Inifor grid.
2501    IF ( init_3d%nx-1 /= nx  .OR.  init_3d%nxu-1 /= nx - 1  .OR.                                   &
2502         init_3d%ny-1 /= ny  .OR.  init_3d%nyv-1 /= ny - 1 )  THEN
2503       message_string = 'Number of horizontal grid points in '//                                   &
2504                        'dynamic input file does not match ' //                                    &
2505                        'the number of numeric grid points.'
2506       CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 )
2507    ENDIF
2508
2509    IF ( init_3d%nzu /= nz )  THEN
2510       message_string = 'Number of vertical grid points in '//                                     &
2511                        'dynamic input file does not match ' //                                    &
2512                        'the number of numeric grid points.'
2513       CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 )
2514    ENDIF
2515!
2516!-- Read vertical dimensions. Later, these are required for eventual inter- and extrapolations of
2517!-- the initialization data.
2518    IF ( check_existence( var_names, 'z' ) )  THEN
2519       ALLOCATE( init_3d%zu_atmos(1:init_3d%nzu) )
2520       CALL get_variable( id_dynamic, 'z', init_3d%zu_atmos )
2521    ENDIF
2522    IF ( check_existence( var_names, 'zw' ) )  THEN
2523       ALLOCATE( init_3d%zw_atmos(1:init_3d%nzw) )
2524       CALL get_variable( id_dynamic, 'zw', init_3d%zw_atmos )
2525    ENDIF
2526!
2527!-- Check for consistency between vertical coordinates in dynamic driver and numeric grid.
2528!-- Please note, depending on compiler options both may be  equal up to a certain threshold, and
2529!-- differences between the numeric grid and vertical coordinate in the driver can built-up to
2530!-- 10E-1-10E-0 m. For this reason, the check is performed not for exactly matching values.
2531    IF ( ANY( ABS( zu(1:nzt)   - init_3d%zu_atmos(1:init_3d%nzu) ) > 10E-1 )  .OR.                 &
2532         ANY( ABS( zw(1:nzt-1) - init_3d%zw_atmos(1:init_3d%nzw) ) > 10E-1 ) )  THEN
2533       message_string = 'Vertical grid in dynamic driver does not ' // 'match the numeric grid.'
2534       CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 )
2535    ENDIF
2536!
2537!-- Read initial geostrophic wind components at t = 0 (index 1 in file).
2538    IF ( check_existence( var_names, 'ls_forcing_ug' ) )  THEN
2539       ALLOCATE( init_3d%ug_init(nzb:nzt+1) )
2540       init_3d%ug_init = 0.0_wp
2541       CALL get_variable_pr( id_dynamic, 'ls_forcing_ug', 1, init_3d%ug_init(1:nzt) )
2542!
2543!--    Set top-boundary condition (Neumann)
2544       init_3d%ug_init(nzt+1) = init_3d%ug_init(nzt)
2545       init_3d%from_file_ug = .TRUE.
2546    ELSE
2547       init_3d%from_file_ug = .FALSE.
2548    ENDIF
2549
2550    IF ( check_existence( var_names, 'ls_forcing_vg' ) )  THEN
2551       ALLOCATE( init_3d%vg_init(nzb:nzt+1) )
2552       init_3d%vg_init = 0.0_wp
2553       CALL get_variable_pr( id_dynamic, 'ls_forcing_vg', 1, init_3d%vg_init(1:nzt) )
2554!
2555!--    Set top-boundary condition (Neumann)
2556       init_3d%vg_init(nzt+1) = init_3d%vg_init(nzt)
2557       init_3d%from_file_vg = .TRUE.
2558    ELSE
2559       init_3d%from_file_vg = .FALSE.
2560    ENDIF
2561!
2562!-- Read inital 3D data of u, v, w, pt and q, derived from COSMO model. Read PE-wise yz-slices.
2563!-- Please note, the u-, v- and w-component are defined on different grids with one element less in
2564!-- the x-, y-, and z-direction, respectively. Hence, reading is subdivided into separate loops.
2565!-- Read u-component
2566    IF ( check_existence( var_names, 'init_atmosphere_u' ) )  THEN
2567!
2568!--    Read attributes for the fill value and level-of-detail
2569       CALL get_attribute( id_dynamic, char_fill, init_3d%fill_u, .FALSE., 'init_atmosphere_u' )
2570       CALL get_attribute( id_dynamic, char_lod, init_3d%lod_u, .FALSE., 'init_atmosphere_u' )
2571!
2572!--    level-of-detail 1 - read initialization profile
2573       IF ( init_3d%lod_u == 1 )  THEN
2574          ALLOCATE( init_3d%u_init(nzb:nzt+1) )
2575          init_3d%u_init = 0.0_wp
2576          CALL get_variable( id_dynamic, 'init_atmosphere_u', init_3d%u_init(nzb+1:nzt) )
2577!
2578!--       Set top-boundary condition (Neumann)
2579          init_3d%u_init(nzt+1) = init_3d%u_init(nzt)
2580!
2581!--    level-of-detail 2 - read 3D initialization data
2582       ELSEIF ( init_3d%lod_u == 2 )  THEN
2583          CALL get_variable( id_dynamic, 'init_atmosphere_u', u(nzb+1:nzt,nys:nyn,nxlu:nxr),       &
2584                             nxlu, nys+1, nzb+1, nxr-nxlu+1, nyn-nys+1, init_3d%nzu, dynamic_3d )
2585!
2586!--       Set value at leftmost model grid point nxl = 0. This is because Inifor provides data only
2587!--       from 1:nx-1 since it assumes non-cyclic conditions.
2588          IF ( nxl == 0 )  u(nzb+1:nzt,nys:nyn,nxl) = u(nzb+1:nzt,nys:nyn,nxlu)
2589!
2590!--       Set bottom and top-boundary
2591          u(nzb,:,:)   = u(nzb+1,:,:)
2592          u(nzt+1,:,:) = u(nzt,:,:)
2593       ENDIF
2594       init_3d%from_file_u = .TRUE.
2595    ELSE
2596       message_string = 'Missing initial data for u-component'
2597       CALL message( 'netcdf_data_input_mod', 'PA0544', 1, 2, 0, 6, 0 )
2598    ENDIF
2599!
2600!-- Read v-component
2601    IF ( check_existence( var_names, 'init_atmosphere_v' ) )  THEN
2602!
2603!--    Read attributes for the fill value and level-of-detail
2604       CALL get_attribute( id_dynamic, char_fill, init_3d%fill_v, .FALSE., 'init_atmosphere_v' )
2605       CALL get_attribute( id_dynamic, char_lod, init_3d%lod_v, .FALSE., 'init_atmosphere_v' )
2606!
2607!--    level-of-detail 1 - read initialization profile
2608       IF ( init_3d%lod_v == 1 )  THEN
2609          ALLOCATE( init_3d%v_init(nzb:nzt+1) )
2610          init_3d%v_init = 0.0_wp
2611          CALL get_variable( id_dynamic, 'init_atmosphere_v', init_3d%v_init(nzb+1:nzt) )
2612!
2613!--       Set top-boundary condition (Neumann)
2614          init_3d%v_init(nzt+1) = init_3d%v_init(nzt)
2615!
2616!--    level-of-detail 2 - read 3D initialization data
2617       ELSEIF ( init_3d%lod_v == 2 )  THEN
2618
2619          CALL get_variable( id_dynamic, 'init_atmosphere_v', v(nzb+1:nzt,nysv:nyn,nxl:nxr),       &
2620                             nxl+1, nysv, nzb+1, nxr-nxl+1, nyn-nysv+1, init_3d%nzu, dynamic_3d )
2621!
2622!--       Set value at southmost model grid point nys = 0. This is because Inifor provides data only
2623!--       from 1:ny-1 since it assumes non-cyclic conditions.
2624          IF ( nys == 0 )  v(nzb+1:nzt,nys,nxl:nxr) = v(nzb+1:nzt,nysv,nxl:nxr)
2625!
2626!--       Set bottom and top-boundary
2627          v(nzb,:,:)   = v(nzb+1,:,:)
2628          v(nzt+1,:,:) = v(nzt,:,:)
2629       ENDIF
2630       init_3d%from_file_v = .TRUE.
2631    ELSE
2632       message_string = 'Missing initial data for v-component'
2633       CALL message( 'netcdf_data_input_mod', 'PA0544', 1, 2, 0, 6, 0 )
2634    ENDIF
2635!
2636!-- Read w-component
2637    IF ( check_existence( var_names, 'init_atmosphere_w' ) )  THEN
2638!
2639!--    Read attributes for the fill value and level-of-detail
2640       CALL get_attribute( id_dynamic, char_fill, init_3d%fill_w, .FALSE., 'init_atmosphere_w' )
2641       CALL get_attribute( id_dynamic, char_lod, init_3d%lod_w, .FALSE., 'init_atmosphere_w' )
2642!
2643!--    level-of-detail 1 - read initialization profile
2644       IF ( init_3d%lod_w == 1 )  THEN
2645          ALLOCATE( init_3d%w_init(nzb:nzt+1) )
2646          init_3d%w_init = 0.0_wp
2647          CALL get_variable( id_dynamic, 'init_atmosphere_w', init_3d%w_init(nzb+1:nzt-1) )
2648!
2649!--       Set top-boundary condition (Neumann)
2650          init_3d%w_init(nzt:nzt+1) = init_3d%w_init(nzt-1)
2651!
2652!--    level-of-detail 2 - read 3D initialization data
2653       ELSEIF ( init_3d%lod_w == 2 )  THEN
2654
2655          CALL get_variable( id_dynamic, 'init_atmosphere_w', w(nzb+1:nzt-1,nys:nyn,nxl:nxr),      &
2656                             nxl+1, nys+1, nzb+1, nxr-nxl+1, nyn-nys+1, init_3d%nzw, dynamic_3d )
2657!
2658!--       Set bottom and top-boundary
2659          w(nzb,:,:)   = 0.0_wp
2660          w(nzt,:,:)   = w(nzt-1,:,:)
2661          w(nzt+1,:,:) = w(nzt-1,:,:)
2662       ENDIF
2663       init_3d%from_file_w = .TRUE.
2664    ELSE
2665       message_string = 'Missing initial data for w-component'
2666       CALL message( 'netcdf_data_input_mod', 'PA0544', 1, 2, 0, 6, 0 )
2667    ENDIF
2668!
2669!-- Read potential temperature
2670    IF ( .NOT. neutral )  THEN
2671       IF ( check_existence( var_names, 'init_atmosphere_pt' ) )  THEN
2672!
2673!--       Read attributes for the fill value and level-of-detail
2674          CALL get_attribute( id_dynamic, char_fill, init_3d%fill_pt, .FALSE., 'init_atmosphere_pt')
2675          CALL get_attribute( id_dynamic, char_lod, init_3d%lod_pt, .FALSE., 'init_atmosphere_pt' )
2676!
2677!--       level-of-detail 1 - read initialization profile
2678          IF ( init_3d%lod_pt == 1 )  THEN
2679             ALLOCATE( init_3d%pt_init(nzb:nzt+1) )
2680             CALL get_variable( id_dynamic, 'init_atmosphere_pt', init_3d%pt_init(nzb+1:nzt) )
2681!
2682!--          Set Neumann top and surface boundary condition for initial profil
2683             init_3d%pt_init(nzb)   = init_3d%pt_init(nzb+1)
2684             init_3d%pt_init(nzt+1) = init_3d%pt_init(nzt)
2685!
2686!--       level-of-detail 2 - read 3D initialization data
2687          ELSEIF ( init_3d%lod_pt == 2 )  THEN
2688
2689             CALL get_variable( id_dynamic, 'init_atmosphere_pt', pt(nzb+1:nzt,nys:nyn,nxl:nxr),   &
2690                                nxl+1, nys+1, nzb+1, nxr-nxl+1, nyn-nys+1, init_3d%nzu, dynamic_3d )
2691!
2692!--          Set bottom and top-boundary
2693             pt(nzb,:,:)   = pt(nzb+1,:,:)
2694             pt(nzt+1,:,:) = pt(nzt,:,:)
2695          ENDIF
2696          init_3d%from_file_pt = .TRUE.
2697       ELSE
2698          message_string = 'Missing initial data for ' // 'potential temperature'
2699          CALL message( 'netcdf_data_input_mod', 'PA0544', 1, 2, 0, 6, 0 )
2700       ENDIF
2701    ENDIF
2702
2703!
2704!-- Read mixing ratio
2705    IF ( humidity )  THEN
2706       IF ( check_existence( var_names, 'init_atmosphere_qv' ) )  THEN
2707!
2708!--       Read attributes for the fill value and level-of-detail
2709          CALL get_attribute( id_dynamic, char_fill, init_3d%fill_q, .FALSE., 'init_atmosphere_qv' )
2710          CALL get_attribute( id_dynamic, char_lod, init_3d%lod_q, .FALSE., 'init_atmosphere_qv' )
2711!
2712!--       level-of-detail 1 - read initialization profile
2713          IF ( init_3d%lod_q == 1 )  THEN
2714             ALLOCATE( init_3d%q_init(nzb:nzt+1) )
2715             CALL get_variable( id_dynamic, 'init_atmosphere_qv', init_3d%q_init(nzb+1:nzt) )
2716!
2717!--          Set bottom and top boundary condition (Neumann)
2718             init_3d%q_init(nzb)   = init_3d%q_init(nzb+1)
2719             init_3d%q_init(nzt+1) = init_3d%q_init(nzt)
2720!
2721!--       level-of-detail 2 - read 3D initialization data
2722          ELSEIF ( init_3d%lod_q == 2 )  THEN
2723
2724             CALL get_variable( id_dynamic, 'init_atmosphere_qv', q(nzb+1:nzt,nys:nyn,nxl:nxr),    &
2725                                nxl+1, nys+1, nzb+1, nxr-nxl+1, nyn-nys+1, init_3d%nzu, dynamic_3d )
2726!
2727!--          Set bottom and top-boundary
2728             q(nzb,:,:)   = q(nzb+1,:,:)
2729             q(nzt+1,:,:) = q(nzt,:,:)
2730          ENDIF
2731          init_3d%from_file_q = .TRUE.
2732       ELSE
2733          message_string = 'Missing initial data for ' // 'mixing ratio'
2734          CALL message( 'netcdf_data_input_mod', 'PA0544', 1, 2, 0, 6, 0 )
2735       ENDIF
2736    ENDIF
2737
2738!
2739!-- Read chemistry variables.
2740!-- Please note, for the moment, only LOD=1 is allowed
2741    IF ( air_chemistry )  THEN
2742!
2743!--    Allocate chemistry input profiles, as well as arrays for fill values and LOD's.
2744       ALLOCATE( init_3d%chem_init(nzb:nzt+1,1:UBOUND( init_3d%var_names_chem, 1 )) )
2745       ALLOCATE( init_3d%fill_chem(1:UBOUND( init_3d%var_names_chem, 1 )) )
2746       ALLOCATE( init_3d%lod_chem(1:UBOUND( init_3d%var_names_chem, 1 ))  )
2747
2748       DO  n = 1, UBOUND( init_3d%var_names_chem, 1 )
2749          IF ( check_existence( var_names, TRIM( init_3d%var_names_chem(n) ) ) )  THEN
2750!
2751!--          Read attributes for the fill value and level-of-detail
2752             CALL get_attribute( id_dynamic, char_fill, init_3d%fill_chem(n), .FALSE.,             &
2753                                 TRIM( init_3d%var_names_chem(n) ) )
2754             CALL get_attribute( id_dynamic, char_lod, init_3d%lod_chem(n), .FALSE.,               &
2755                                 TRIM( init_3d%var_names_chem(n) ) )
2756!
2757!--          Give message that only LOD=1 is allowed.
2758             IF ( init_3d%lod_chem(n) /= 1 )  THEN
2759                message_string = 'For chemistry variables only LOD=1 is ' // 'allowed.'
2760                CALL message( 'netcdf_data_input_mod', 'PA0586', 1, 2, 0, 6, 0 )
2761             ENDIF
2762!
2763!--          level-of-detail 1 - read initialization profile
2764             CALL get_variable( id_dynamic, TRIM( init_3d%var_names_chem(n) ),                     &
2765                                init_3d%chem_init(nzb+1:nzt,n) )
2766!
2767!--          Set bottom and top boundary condition (Neumann)
2768             init_3d%chem_init(nzb,n)   = init_3d%chem_init(nzb+1,n)
2769             init_3d%chem_init(nzt+1,n) = init_3d%chem_init(nzt,n)
2770             init_3d%from_file_chem(n)  = .TRUE.
2771          ENDIF
2772       ENDDO
2773    ENDIF
2774!
2775!-- Close input file
2776    CALL close_input_file( id_dynamic )
2777#endif
2778!
2779!-- End of CPU measurement
2780    CALL cpu_log( log_point_s(85), 'NetCDF input init', 'stop' )
2781!
2782!-- Finally, check if the input data has any fill values. Please note, checks depend on the LOD of
2783!-- the input data.
2784    IF ( init_3d%from_file_u )  THEN
2785       check_passed = .TRUE.
2786       IF ( init_3d%lod_u == 1 )  THEN
2787          IF ( ANY( init_3d%u_init(nzb+1:nzt+1) == init_3d%fill_u ) )  check_passed = .FALSE.
2788       ELSEIF ( init_3d%lod_u == 2 )  THEN
2789          IF ( ANY( u(nzb+1:nzt+1,nys:nyn,nxlu:nxr) == init_3d%fill_u ) )  check_passed = .FALSE.
2790       ENDIF
2791       IF ( .NOT. check_passed )  THEN
2792          message_string = 'NetCDF input for init_atmosphere_u must ' //                           &
2793                           'not contain any _FillValues'
2794          CALL message( 'netcdf_data_input_mod', 'PA0545', 2, 2, 0, 6, 0 )
2795       ENDIF
2796    ENDIF
2797
2798    IF ( init_3d%from_file_v )  THEN
2799       check_passed = .TRUE.
2800       IF ( init_3d%lod_v == 1 )  THEN
2801          IF ( ANY( init_3d%v_init(nzb+1:nzt+1) == init_3d%fill_v ) )  check_passed = .FALSE.
2802       ELSEIF ( init_3d%lod_v == 2 )  THEN
2803          IF ( ANY( v(nzb+1:nzt+1,nysv:nyn,nxl:nxr) == init_3d%fill_v ) )  check_passed = .FALSE.
2804       ENDIF
2805       IF ( .NOT. check_passed )  THEN
2806          message_string = 'NetCDF input for init_atmosphere_v must ' //                           &
2807                           'not contain any _FillValues'
2808          CALL message( 'netcdf_data_input_mod', 'PA0545', 2, 2, 0, 6, 0 )
2809       ENDIF
2810    ENDIF
2811
2812    IF ( init_3d%from_file_w )  THEN
2813       check_passed = .TRUE.
2814       IF ( init_3d%lod_w == 1 )  THEN
2815          IF ( ANY( init_3d%w_init(nzb+1:nzt) == init_3d%fill_w ) )  check_passed = .FALSE.
2816       ELSEIF ( init_3d%lod_w == 2 )  THEN
2817          IF ( ANY( w(nzb+1:nzt,nys:nyn,nxl:nxr) == init_3d%fill_w ) )  check_passed = .FALSE.
2818       ENDIF
2819       IF ( .NOT. check_passed )  THEN
2820          message_string = 'NetCDF input for init_atmosphere_w must ' //                           &
2821                           'not contain any _FillValues'
2822          CALL message( 'netcdf_data_input_mod', 'PA0545', 2, 2, 0, 6, 0 )
2823       ENDIF
2824    ENDIF
2825
2826    IF ( init_3d%from_file_pt )  THEN
2827       check_passed = .TRUE.
2828       IF ( init_3d%lod_pt == 1 )  THEN
2829          IF ( ANY( init_3d%pt_init(nzb+1:nzt+1) == init_3d%fill_pt ) )  check_passed = .FALSE.
2830       ELSEIF ( init_3d%lod_pt == 2 )  THEN
2831          IF ( ANY( pt(nzb+1:nzt+1,nys:nyn,nxl:nxr) == init_3d%fill_pt ) )  check_passed = .FALSE.
2832       ENDIF
2833       IF ( .NOT. check_passed )  THEN
2834          message_string = 'NetCDF input for init_atmosphere_pt must ' //                          &
2835                           'not contain any _FillValues'
2836          CALL message( 'netcdf_data_input_mod', 'PA0545', 2, 2, 0, 6, 0 )
2837       ENDIF
2838    ENDIF
2839
2840    IF ( init_3d%from_file_q )  THEN
2841       check_passed = .TRUE.
2842       IF ( init_3d%lod_q == 1 )  THEN
2843          IF ( ANY( init_3d%q_init(nzb+1:nzt+1) == init_3d%fill_q ) )  check_passed = .FALSE.
2844       ELSEIF ( init_3d%lod_q == 2 )  THEN
2845          IF ( ANY( q(nzb+1:nzt+1,nys:nyn,nxl:nxr) == init_3d%fill_q ) )  check_passed = .FALSE.
2846       ENDIF
2847       IF ( .NOT. check_passed )  THEN
2848          message_string = 'NetCDF input for init_atmosphere_q must ' //                           &
2849                           'not contain any _FillValues'
2850          CALL message( 'netcdf_data_input_mod', 'PA0545', 2, 2, 0, 6, 0 )
2851       ENDIF
2852    ENDIF
2853!
2854!-- Workaround for cyclic conditions. Please see above for further explanation.
2855    IF ( bc_lr_cyc  .AND.  nxl == 0 )  nxlu = nxl
2856    IF ( bc_ns_cyc  .AND.  nys == 0 )  nysv = nys
2857
2858 END SUBROUTINE netcdf_data_input_init_3d
2859
2860
2861!--------------------------------------------------------------------------------------------------!
2862! Description:
2863! ------------
2864!> Checks input file for consistency and minimum requirements.
2865!--------------------------------------------------------------------------------------------------!
2866 SUBROUTINE netcdf_data_input_check_dynamic
2867
2868    USE control_parameters,                                                                        &
2869        ONLY:  initializing_actions, message_string
2870
2871    IMPLICIT NONE
2872!
2873!-- Dynamic input file must also be present if initialization via inifor is prescribed.
2874    IF ( .NOT. input_pids_dynamic  .AND.  INDEX( initializing_actions, 'inifor' ) /= 0 )  THEN
2875       message_string = 'initializing_actions = inifor requires dynamic ' //                       &
2876                        'input file ' // TRIM( input_file_dynamic ) // TRIM( coupling_char )
2877       CALL message( 'netcdf_data_input_mod', 'PA0547', 1, 2, 0, 6, 0 )
2878    ENDIF
2879
2880 END SUBROUTINE netcdf_data_input_check_dynamic
2881
2882
2883!--------------------------------------------------------------------------------------------------!
2884! Description:
2885! ------------
2886!> Checks input file for consistency and minimum requirements.
2887!--------------------------------------------------------------------------------------------------!
2888 SUBROUTINE netcdf_data_input_check_static
2889
2890    USE arrays_3d,                                                                                 &
2891        ONLY:  zu
2892
2893    USE control_parameters,                                                                        &
2894        ONLY:  land_surface, message_string, urban_surface
2895
2896    USE indices,                                                                                   &
2897        ONLY:  nxl, nxr, nyn, nys, wall_flags_total_0
2898
2899    IMPLICIT NONE
2900
2901    INTEGER(iwp) ::  i       !< loop index along x-direction
2902    INTEGER(iwp) ::  j       !< loop index along y-direction
2903    INTEGER(iwp) ::  n_surf  !< number of different surface types at given location
2904
2905    LOGICAL      ::  check_passed  !< flag indicating if a check passed
2906
2907!
2908!-- Return if no static input file is available
2909    IF ( .NOT. input_pids_static )  RETURN
2910!
2911!-- Check for correct dimension of surface_fractions, should run from 0-2.
2912    IF ( surface_fraction_f%from_file )  THEN
2913       IF ( surface_fraction_f%nf-1 > 2 )  THEN
2914          message_string = 'nsurface_fraction must not be larger than 3.'
2915          CALL message( 'netcdf_data_input_mod', 'PA0580', 1, 2, 0, 6, 0 )
2916       ENDIF
2917    ENDIF
2918!
2919!-- If 3D buildings are read, check if building information is consistent to numeric grid.
2920    IF ( buildings_f%from_file )  THEN
2921       IF ( buildings_f%lod == 2 )  THEN
2922          IF ( buildings_f%nz > SIZE( zu ) )  THEN
2923             message_string = 'Reading 3D building data - too much ' //                            &
2924                              'data points along the vertical coordinate.'
2925             CALL message( 'netcdf_data_input_mod', 'PA0552', 2, 2, 0, 6, 0 )
2926          ENDIF
2927
2928          IF ( ANY( ABS( buildings_f%z(0:buildings_f%nz-1) - zu(0:buildings_f%nz-1) ) > 1E-6_wp ) )&
2929          THEN
2930             message_string = 'Reading 3D building data - vertical ' //                            &
2931                              'coordinate do not match numeric grid.'
2932             CALL message( 'netcdf_data_input_mod', 'PA0553', 2, 2, 0, 6, 0 )
2933          ENDIF
2934       ENDIF
2935    ENDIF
2936
2937!
2938!-- Skip further checks concerning buildings and natural surface properties if no urban surface and
2939!-- land surface model are applied.
2940    IF (  .NOT. land_surface  .AND.  .NOT. urban_surface )  RETURN
2941!
2942!-- Check for minimum requirement of surface-classification data in case static input file is used.
2943    IF ( ( .NOT. vegetation_type_f%from_file .OR.                                                  &
2944           .NOT. pavement_type_f%from_file   .OR.                                                  &
2945           .NOT. water_type_f%from_file      .OR.                                                  &
2946           .NOT. soil_type_f%from_file            )  .OR.                                          &
2947          ( urban_surface  .AND.  .NOT. building_type_f%from_file ) )  THEN
2948       message_string = 'Minimum requirement for surface classification is not fulfilled. At ' //  &
2949                        'least vegetation_type, pavement_type, soil_type and water_type are ' //   &
2950                        'required. If urban-surface model is applied, also building_type is required'
2951       CALL message( 'netcdf_data_input_mod', 'PA0554', 1, 2, 0, 6, 0 )
2952    ENDIF
2953!
2954!-- Check for general availability of input variables.
2955!-- If vegetation_type is 0 at any location, vegetation_pars as well as root_area_dens_s are
2956!-- required.
2957    IF ( vegetation_type_f%from_file )  THEN
2958       IF ( ANY( vegetation_type_f%var == 0 ) )  THEN
2959          IF ( .NOT. vegetation_pars_f%from_file )  THEN
2960             message_string = 'If vegetation_type = 0 at any location, ' //                        &
2961                              'vegetation_pars is required'
2962             CALL message( 'netcdf_data_input_mod', 'PA0555', 2, 2, myid, 6, 0 )
2963          ENDIF
2964          IF ( .NOT. root_area_density_lsm_f%from_file )  THEN
2965             message_string = 'If vegetation_type = 0 at any location, ' //                        &
2966                              'root_area_dens_s is required'
2967             CALL message( 'netcdf_data_input_mod', 'PA0556', 2, 2, myid, 6, 0 )
2968          ENDIF
2969       ENDIF
2970    ENDIF
2971!
2972!-- If soil_type is zero at any location, soil_pars is required.
2973    IF ( soil_type_f%from_file )  THEN
2974       check_passed = .TRUE.
2975       IF ( ALLOCATED( soil_type_f%var_2d ) )  THEN
2976          IF ( ANY( soil_type_f%var_2d == 0 ) )  THEN
2977             IF ( .NOT. soil_pars_f%from_file )  check_passed = .FALSE.
2978          ENDIF
2979       ELSE
2980          IF ( ANY( soil_type_f%var_3d == 0 ) )  THEN
2981             IF ( .NOT. soil_pars_f%from_file )  check_passed = .FALSE.
2982          ENDIF
2983       ENDIF
2984       IF ( .NOT. check_passed )  THEN
2985          message_string = 'If soil_type = 0 at any location, soil_pars is required'
2986          CALL message( 'netcdf_data_input_mod', 'PA0557', 2, 2, myid, 6, 0 )
2987       ENDIF
2988    ENDIF
2989!
2990!-- Buildings require a type in case of urban-surface model.
2991    IF ( buildings_f%from_file  .AND.  .NOT. building_type_f%from_file  )  THEN
2992       message_string = 'If buildings are provided, also building_type is required'
2993       CALL message( 'netcdf_data_input_mod', 'PA0581', 2, 2, 0, 6, 0 )
2994    ENDIF
2995!
2996!-- Buildings require an ID.
2997    IF ( buildings_f%from_file  .AND.  .NOT. building_id_f%from_file  )  THEN
2998       message_string = 'If buildings are provided, also building_id is required'
2999       CALL message( 'netcdf_data_input_mod', 'PA0582', 2, 2, 0, 6, 0 )
3000    ENDIF
3001!
3002!-- If building_type is zero at any location, building_pars is required.
3003    IF ( building_type_f%from_file )  THEN
3004       IF ( ANY( building_type_f%var == 0 ) )  THEN
3005          IF ( .NOT. building_pars_f%from_file )  THEN
3006             message_string = 'If building_type = 0 at any location, ' //                          &
3007                              'building_pars is required'
3008             CALL message( 'netcdf_data_input_mod', 'PA0558', 2, 2, myid, 6, 0 )
3009          ENDIF
3010       ENDIF
3011    ENDIF
3012!
3013!-- If building_type is provided, also building_id is needed (due to the filtering algorithm).
3014    IF ( building_type_f%from_file  .AND.  .NOT. building_id_f%from_file )  THEN
3015       message_string = 'If building_type is provided, also building_id is required'
3016       CALL message( 'netcdf_data_input_mod', 'PA0519', 2, 2, 0, 6, 0 )
3017    ENDIF
3018!
3019!-- If albedo_type is zero at any location, albedo_pars is required.
3020    IF ( albedo_type_f%from_file )  THEN
3021       IF ( ANY( albedo_type_f%var == 0 ) )  THEN
3022          IF ( .NOT. albedo_pars_f%from_file )  THEN
3023             message_string = 'If albedo_type = 0 at any location, albedo_pars is required'
3024             CALL message( 'netcdf_data_input_mod', 'PA0559', 2, 2, myid, 6, 0 )
3025          ENDIF
3026       ENDIF
3027    ENDIF
3028!
3029!-- If pavement_type is zero at any location, pavement_pars is required.
3030    IF ( pavement_type_f%from_file )  THEN
3031       IF ( ANY( pavement_type_f%var == 0 ) )  THEN
3032          IF ( .NOT. pavement_pars_f%from_file )  THEN
3033             message_string = 'If pavement_type = 0 at any location, pavement_pars is required'
3034             CALL message( 'netcdf_data_input_mod', 'PA0560', 2, 2, myid, 6, 0 )
3035          ENDIF
3036       ENDIF
3037    ENDIF
3038!
3039!-- If pavement_type is zero at any location, also pavement_subsurface_pars is required.
3040    IF ( pavement_type_f%from_file )  THEN
3041       IF ( ANY( pavement_type_f%var == 0 ) )  THEN
3042          IF ( .NOT. pavement_subsurface_pars_f%from_file )  THEN
3043             message_string = 'If pavement_type = 0 at any location, ' //                          &
3044                              'pavement_subsurface_pars is required'
3045             CALL message( 'netcdf_data_input_mod', 'PA0561', 2, 2, myid, 6, 0 )
3046          ENDIF
3047       ENDIF
3048    ENDIF
3049!
3050!-- If water_type is zero at any location, water_pars is required.
3051    IF ( water_type_f%from_file )  THEN
3052       IF ( ANY( water_type_f%var == 0 ) )  THEN
3053          IF ( .NOT. water_pars_f%from_file )  THEN
3054             message_string = 'If water_type = 0 at any location, water_pars is required'
3055             CALL message( 'netcdf_data_input_mod', 'PA0562', 2, 2, myid, 6, 0 )
3056          ENDIF
3057       ENDIF
3058    ENDIF
3059!
3060!-- Check for local consistency of the input data.
3061    DO  i = nxl, nxr
3062       DO  j = nys, nyn
3063!
3064!--       For each (y,x)-location at least one of the parameters vegetation_type, pavement_type,
3065!--       building_type, or water_type must be set to a non­missing value.
3066          IF ( land_surface  .AND.  .NOT. urban_surface )  THEN
3067             IF ( vegetation_type_f%var(j,i) == vegetation_type_f%fill  .AND.                      &
3068                  pavement_type_f%var(j,i)   == pavement_type_f%fill    .AND.                      &
3069                  water_type_f%var(j,i)      == water_type_f%fill )  THEN
3070                WRITE( message_string, * ) 'At least one of the parameters '//                     &
3071                                           'vegetation_type, pavement_type, ' //                   &
3072                                           'or water_type must be set '//                          &
3073                                           'to a non-missing value. Grid point: ', j, i
3074                CALL message( 'netcdf_data_input_mod', 'PA0563', 2, 2, myid, 6, 0 )
3075             ENDIF
3076          ELSEIF ( land_surface  .AND.  urban_surface )  THEN
3077             IF ( vegetation_type_f%var(j,i) == vegetation_type_f%fill  .AND.                      &
3078                  pavement_type_f%var(j,i)   == pavement_type_f%fill    .AND.                      &
3079                  building_type_f%var(j,i)   == building_type_f%fill    .AND.                      &
3080                  water_type_f%var(j,i)      == water_type_f%fill )  THEN
3081                WRITE( message_string, * ) 'At least one of the parameters '//                     &
3082                                           'vegetation_type, pavement_type, '  //                  &
3083                                           'building_type, or water_type must be set '//           &
3084                                           'to a non-missing value. Grid point: ', j, i
3085                CALL message( 'netcdf_data_input_mod', 'PA0563', 2, 2, myid, 6, 0 )
3086             ENDIF
3087          ENDIF
3088
3089!
3090!--       Note that a soil_type is required for each location (y,x) where either vegetation_type or
3091!--       pavement_type is a non­missing value.
3092          IF ( ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill  .OR.                        &
3093                 pavement_type_f%var(j,i)   /= pavement_type_f%fill ) )  THEN
3094             check_passed = .TRUE.
3095             IF ( ALLOCATED( soil_type_f%var_2d ) )  THEN
3096                IF ( soil_type_f%var_2d(j,i) == soil_type_f%fill )  check_passed = .FALSE.
3097             ELSE
3098                IF ( ANY( soil_type_f%var_3d(:,j,i) == soil_type_f%fill) )  check_passed = .FALSE.
3099             ENDIF
3100
3101             IF ( .NOT. check_passed )  THEN
3102                message_string = 'soil_type is required for each '//                               &
3103                                 'location (y,x) where vegetation_type or ' //                     &
3104                                 'pavement_type is a non-missing value.'
3105                CALL message( 'netcdf_data_input_mod', 'PA0564', 2, 2, myid, 6, 0 )
3106             ENDIF
3107          ENDIF
3108!
3109!--       Check for consistency of given types. At the moment, only one of vegetation, pavement, or
3110!--       water-type can be set. This is because no tile approach is yet implemented in the
3111!--       land-surface model. Later, when this is possible, surface fraction needs to be given and
3112!--       the sum must not be larger than 1. Please note, in case more than one type is given at a
3113!--       pixel, an error message will be given.
3114          n_surf = 0
3115          IF ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill )  n_surf = n_surf + 1
3116          IF ( water_type_f%var(j,i)      /= water_type_f%fill      )  n_surf = n_surf + 1
3117          IF ( pavement_type_f%var(j,i)   /= pavement_type_f%fill   )  n_surf = n_surf + 1
3118
3119          IF ( n_surf > 1 )  THEN
3120             WRITE( message_string, * ) 'More than one surface type (vegetation, ' //              &
3121                                        'pavement, water) is given at a location. ' //             &
3122                                        'Please note, this is not possible at ' //                 &
3123                                        'the moment as no tile approach has been ' // &
3124                                        'yet implemented. (i,j) = ', i, j
3125             CALL message( 'netcdf_data_input_mod', 'PA0565',               &
3126                            2, 2, myid, 6, 0 )
3127
3128!              IF ( .NOT. surface_fraction_f%from_file )  THEN
3129!                 message_string = 'More than one surface type (vegetation '//                     &
3130!                                  'pavement, water) is given at a location. '//                   &
3131!                                  'Please note, this is not possible at ' //                      &
3132!                                  'the moment as no tile approach is yet ' //                     &
3133!                                  'implemented.'
3134!                 message_string = 'If more than one surface type is ' //                          &
3135!                                  'given at a location, surface_fraction ' //                     &
3136!                                  'must be provided.'
3137!                 CALL message( 'netcdf_data_input_mod', 'PA0565', 2, 2, myid, 6, 0 )
3138!              ELSEIF ( ANY ( surface_fraction_f%frac(:,j,i) ==  surface_fraction_f%fill ) )  THEN
3139!                 message_string = 'If more than one surface type is ' //                          &
3140!                                  'given at a location, surface_fraction ' //                     &
3141!                                  'must be provided.'
3142!                 CALL message( 'netcdf_data_input_mod', 'PA0565', 2, 2, myid, 6, 0 )
3143!              ENDIF
3144          ENDIF
3145!
3146!--       Check for further mismatches. e.g. relative fractions exceed 1 or vegetation_type is set
3147!--       but surface vegetation fraction is zero, etc..
3148          IF ( surface_fraction_f%from_file )  THEN
3149!
3150!--          If surface fractions is given, also check that only one type is given.
3151             IF ( SUM( MERGE( 1, 0, surface_fraction_f%frac(:,j,i) /= 0.0_wp  .AND.                &
3152                                    surface_fraction_f%frac(:,j,i) /= surface_fraction_f%fill )    &
3153                     ) > 1 )                                                                       &
3154             THEN
3155                WRITE( message_string, * ) 'surface_fraction is given for more ' //                &
3156                                           'than one type. ' //                                    &
3157                                           'Please note, this is not possible at ' //              &
3158                                           'the moment as no tile approach has ' //                &
3159                                           'yet been implemented. (i, j) = ', i, j
3160                CALL message( 'netcdf_data_input_mod', 'PA0676', 2, 2, myid, 6, 0 )
3161             ENDIF
3162!
3163!--          Sum of relative fractions must be 1. Note, attributed to type conversions due to
3164!--          reading, the sum of surface fractions might be not exactly 1. Hence, the sum is check
3165!--          with a tolerance. Later, in the land-surface model, the relative fractions are
3166!--          normalized to one. Actually, surface fractions shall be _FillValue at building grid
3167!--          points, however, in order to relax this requirement and allow that surface-fraction can
3168!--          also be zero at these grid points, only perform this check at locations where some
3169!--          vegetation, pavement or water is defined.
3170             IF ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill  .OR.                       &
3171                  pavement_type_f%var(j,i)   /= pavement_type_f%fill    .OR.                       &
3172                  water_type_f%var(j,i)      /= water_type_f%fill )  THEN
3173                IF ( SUM ( surface_fraction_f%frac(0:2,j,i) ) > 1.0_wp + 1E-8_wp  .OR.             &
3174                     SUM ( surface_fraction_f%frac(0:2,j,i) ) < 1.0_wp - 1E-8_wp )  THEN
3175                   WRITE( message_string, * ) 'The sum of all land-surface fractions ' //          &
3176                                              'must equal 1. (i, j) = ', i, j
3177                   CALL message( 'netcdf_data_input_mod', 'PA0566', 2, 2, myid, 6, 0 )
3178                ENDIF
3179             ENDIF
3180!
3181!--          Relative fraction for a type must not be zero at locations where this type is set.
3182             IF ( ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill  .AND.                    &
3183                    ( surface_fraction_f%frac(ind_veg_wall,j,i) == 0.0_wp .OR.                     &
3184                      surface_fraction_f%frac(ind_veg_wall,j,i) == surface_fraction_f%fill )       &
3185                  )  .OR.                                                                          &
3186                  ( pavement_type_f%var(j,i) /= pavement_type_f%fill       .AND.                   &
3187                    ( surface_fraction_f%frac(ind_pav_green,j,i) == 0.0_wp .OR.                    &
3188                      surface_fraction_f%frac(ind_pav_green,j,i) == surface_fraction_f%fill )      &
3189                  )  .OR.                                                                          &
3190                  ( water_type_f%var(j,i) /= water_type_f%fill             .AND.                   &
3191                   ( surface_fraction_f%frac(ind_wat_win,j,i) == 0.0_wp      .OR.                  &
3192                     surface_fraction_f%frac(ind_wat_win,j,i) == surface_fraction_f%fill )         &
3193                  ) )  THEN
3194                WRITE( message_string, * ) 'Mismatch in setting of '     //                        &
3195                                           'surface_fraction. Vegetation-, pavement-, or ' //      &
3196                                           'water surface is given at (i,j) = ( ', i, j,           &
3197                                           ' ), but surface fraction is 0 for the given type.'
3198                CALL message( 'netcdf_data_input_mod', 'PA0567', 2, 2, myid, 6, 0 )
3199             ENDIF
3200!
3201!--          Relative fraction for a type must not contain non-zero values if this type is not set.
3202             IF ( ( vegetation_type_f%var(j,i) == vegetation_type_f%fill    .AND.                  &
3203                    ( surface_fraction_f%frac(ind_veg_wall,j,i) /= 0.0_wp   .AND.                  &
3204                      surface_fraction_f%frac(ind_veg_wall,j,i) /= surface_fraction_f%fill )       &
3205                  )  .OR.                                                                          &
3206                  ( pavement_type_f%var(j,i) == pavement_type_f%fill        .AND.                  &
3207                    ( surface_fraction_f%frac(ind_pav_green,j,i) /= 0.0_wp .AND.                   &
3208                      surface_fraction_f%frac(ind_pav_green,j,i) /= surface_fraction_f%fill )      &
3209                  )  .OR.                                                                          &
3210                  ( water_type_f%var(j,i) == water_type_f%fill           .AND.                     &
3211                    ( surface_fraction_f%frac(ind_wat_win,j,i) /= 0.0_wp .AND.                     &
3212                      surface_fraction_f%frac(ind_wat_win,j,i) /= surface_fraction_f%fill )        &
3213                  ) )  THEN
3214                WRITE( message_string, * ) 'Mismatch in setting of '     //                        &
3215                                           'surface_fraction. Vegetation-, pavement-, or '//       &
3216                                           'water surface is not given at (i,j) = ( ', i, j,       &
3217                                           ' ), but surface fraction is not 0 for the ' //         &
3218                                           'given type.'
3219                CALL message( 'netcdf_data_input_mod', 'PA0568', 2, 2, myid, 6, 0 )
3220             ENDIF
3221          ENDIF
3222!
3223!--       Check vegetation_pars. If vegetation_type is 0, all parameters need to be set, otherwise,
3224!--       single parameters set by vegetation_type can be overwritten.
3225          IF ( vegetation_type_f%from_file )  THEN
3226             IF ( vegetation_type_f%var(j,i) == 0 )  THEN
3227                IF ( ANY( vegetation_pars_f%pars_xy(:,j,i) == vegetation_pars_f%fill ) )  THEN
3228                   message_string = 'If vegetation_type(y,x) = 0, all '  //                        &
3229                                    'parameters of vegetation_pars at ' //                         &
3230                                    'this location must be set.'
3231                   CALL message( 'netcdf_data_input_mod', 'PA0569', 2, 2, myid, 6, 0 )
3232                ENDIF
3233             ENDIF
3234          ENDIF
3235!
3236!--       Check root distribution. If vegetation_type is 0, all levels must be set.
3237          IF ( vegetation_type_f%from_file )  THEN
3238             IF ( vegetation_type_f%var(j,i) == 0 )  THEN
3239                IF ( ANY( root_area_density_lsm_f%var(:,j,i) == root_area_density_lsm_f%fill ) )   &
3240                THEN
3241                   message_string = 'If vegetation_type(y,x) = 0, all ' //                         &
3242                                    'levels of root_area_dens_s ' //                               &
3243                                    'must be set at this location.'
3244                   CALL message( 'netcdf_data_input_mod', 'PA0570', 2, 2, myid, 6, 0 )
3245                ENDIF
3246             ENDIF
3247          ENDIF
3248!
3249!--       Check soil parameters. If soil_type is 0, all parameters must be set.
3250          IF ( soil_type_f%from_file )  THEN
3251             check_passed = .TRUE.
3252             IF ( ALLOCATED( soil_type_f%var_2d ) )  THEN
3253                IF ( soil_type_f%var_2d(j,i) == 0 )  THEN
3254                   IF ( ANY( soil_pars_f%pars_xy(:,j,i) == soil_pars_f%fill ) )                    &
3255                      check_passed = .FALSE.
3256                ENDIF
3257             ELSE
3258                IF ( ANY( soil_type_f%var_3d(:,j,i) == 0 ) )  THEN
3259                   IF ( ANY( soil_pars_f%pars_xy(:,j,i) == soil_pars_f%fill ) )                    &
3260                      check_passed = .FALSE.
3261                ENDIF
3262             ENDIF
3263             IF ( .NOT. check_passed )  THEN
3264                message_string = 'If soil_type(y,x) = 0, all levels of '  //                       &
3265                                 'soil_pars at this location must be set.'
3266                CALL message( 'netcdf_data_input_mod', 'PA0571', 2, 2, myid, 6, 0 )
3267             ENDIF
3268          ENDIF
3269
3270!
3271!--       Check building parameters. If building_type is 0, all parameters 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) == building_pars_f%fill ) )  THEN
3275                   message_string = 'If building_type(y,x) = 0, all ' //                           &
3276                                    'parameters of building_pars at this ' //                      &
3277                                    'location must be set.'
3278                   CALL message( 'netcdf_data_input_mod', 'PA0572', 2, 2, myid, 6, 0 )
3279                ENDIF
3280             ENDIF
3281          ENDIF
3282!
3283!--       Check if building_type is set at each building and vice versa.
3284!--       Please note, buildings are already processed and filtered.
3285!--       For this reason, consistency checks are based on wall_flags_total_0 rather than
3286!--       buildings_f (buildings are represented by bit 6 in wall_flags_total_0).
3287          IF ( building_type_f%from_file  .AND.  buildings_f%from_file )  THEN
3288             IF ( ANY( BTEST ( wall_flags_total_0(:,j,i), 6 ) )  .AND.                             &
3289                  building_type_f%var(j,i) == building_type_f%fill  .OR.                           &
3290                  .NOT. ANY( BTEST ( wall_flags_total_0(:,j,i), 6 ) )  .AND.                       &
3291                  building_type_f%var(j,i) /= building_type_f%fill )  THEN
3292                WRITE( message_string, * ) 'Each location where a ' //                             &
3293                                           'building is set requires a type ' //                   &
3294                                           '( and vice versa ) in case the ' //                    &
3295                                           'urban-surface model is applied. ' //                   &
3296                                           'i, j = ', i, j
3297                CALL message( 'netcdf_data_input_mod', 'PA0573', 2, 2, myid, 6, 0 )
3298             ENDIF
3299          ENDIF
3300!
3301!--       Check if at each location where a building is present also an ID is set and vice versa.
3302          IF ( buildings_f%from_file )  THEN
3303             IF ( ANY( BTEST ( wall_flags_total_0(:,j,i), 6 ) )  .AND.                             &
3304                  building_id_f%var(j,i) == building_id_f%fill  .OR.                               &
3305                  .NOT. ANY( BTEST ( wall_flags_total_0(:,j,i), 6 ) )  .AND.                       &
3306                  building_id_f%var(j,i) /= building_id_f%fill )  THEN
3307                WRITE( message_string, * ) 'Each location where a ' //                             &
3308                                           'building is set requires an ID ' //                    &
3309                                           '( and vice versa ). i, j = ', i, j
3310                CALL message( 'netcdf_data_input_mod', 'PA0574', 2, 2, myid, 6, 0 )
3311             ENDIF
3312          ENDIF
3313!
3314!--       Check if building ID is set where a bulding is defined.
3315          IF ( buildings_f%from_file )  THEN
3316             IF ( ANY( BTEST ( wall_flags_total_0(:,j,i), 6 ) )  .AND.                             &
3317                  building_id_f%var(j,i) == building_id_f%fill )  THEN
3318                WRITE( message_string, * ) 'Each building grid point ' // 'requires an ID.', i, j
3319                CALL message( 'netcdf_data_input_mod', 'PA0575', 2, 2, myid, 6, 0 )
3320             ENDIF
3321          ENDIF
3322!
3323!--       Check albedo parameters. If albedo_type is 0, all parameters must be set.
3324          IF ( albedo_type_f%from_file )  THEN
3325             IF ( albedo_type_f%var(j,i) == 0 )  THEN
3326                IF ( ANY( albedo_pars_f%pars_xy(:,j,i) == albedo_pars_f%fill ) )  THEN
3327                   message_string = 'If albedo_type(y,x) = 0, all ' //                             &
3328                                    'parameters of albedo_pars at this ' //                        &
3329                                    'location must be set.'
3330                   CALL message( 'netcdf_data_input_mod', 'PA0576', 2, 2, myid, 6, 0 )
3331                ENDIF
3332             ENDIF
3333          ENDIF
3334
3335!
3336!--       Check pavement parameters. If pavement_type is 0, all parameters of pavement_pars must be
3337!--       set at this location.
3338          IF ( pavement_type_f%from_file )  THEN
3339             IF ( pavement_type_f%var(j,i) == 0 )  THEN
3340                IF ( ANY( pavement_pars_f%pars_xy(:,j,i) == pavement_pars_f%fill ) )  THEN
3341                   message_string = 'If pavement_type(y,x) = 0, all ' //                           &
3342                                    'parameters of pavement_pars at this '//                       &
3343                                    'location must be set.'
3344                   CALL message( 'netcdf_data_input_mod', 'PA0577', 2, 2, myid, 6, 0 )
3345                ENDIF
3346             ENDIF
3347          ENDIF
3348!
3349!--       Check pavement-subsurface parameters. If pavement_type is 0, all parameters of
3350!--       pavement_subsurface_pars must be set at this location.
3351          IF ( pavement_type_f%from_file )  THEN
3352             IF ( pavement_type_f%var(j,i) == 0 )  THEN
3353                IF ( ANY( pavement_subsurface_pars_f%pars_xyz(:,:,j,i) ==                          &
3354                          pavement_subsurface_pars_f%fill ) )  THEN
3355                   message_string = 'If pavement_type(y,x) = 0, all '   //                         &
3356                                    'parameters of '                    //                         &
3357                                    'pavement_subsurface_pars at this ' //                         &
3358                                    'location must be set.'
3359                   CALL message( 'netcdf_data_input_mod', 'PA0578', 2, 2, myid, 6, 0 )
3360                ENDIF
3361             ENDIF
3362          ENDIF
3363
3364!
3365!--       Check water parameters. If water_type is 0, all parameters must be set  at this location.
3366          IF ( water_type_f%from_file )  THEN
3367             IF ( water_type_f%var(j,i) == 0 )  THEN
3368                IF ( ANY( water_pars_f%pars_xy(:,j,i) == water_pars_f%fill ) )  THEN
3369                   message_string = 'If water_type(y,x) = 0, all ' //                              &
3370                                    'parameters of water_pars at this ' //                         &
3371                                    'location must be set.'
3372                   CALL message( 'netcdf_data_input_mod', 'PA0579', 2, 2, myid, 6, 0 )
3373                ENDIF
3374             ENDIF
3375          ENDIF
3376
3377       ENDDO
3378    ENDDO
3379
3380 END SUBROUTINE netcdf_data_input_check_static
3381
3382
3383!--------------------------------------------------------------------------------------------------!
3384! Description:
3385! ------------
3386!> Resize 8-bit 2D Integer array: (nys:nyn,nxl:nxr) -> (nysg:nyng,nxlg:nxrg)
3387!--------------------------------------------------------------------------------------------------!
3388 SUBROUTINE add_ghost_layers_2d_int8( var )
3389
3390    IMPLICIT NONE
3391
3392    INTEGER(KIND=1), DIMENSION(:,:), ALLOCATABLE ::  var     !< treated variable
3393    INTEGER(KIND=1), DIMENSION(:,:), ALLOCATABLE ::  var_tmp !< temporary copy
3394!
3395!-- Allocate temporary variable
3396    ALLOCATE( var_tmp(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3397!
3398!-- Temporary copy of the variable
3399    var_tmp(nys:nyn,nxl:nxr) = var(nys:nyn,nxl:nxr)
3400!
3401!-- Resize the array
3402    DEALLOCATE( var )
3403    ALLOCATE( var(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3404!
3405!-- Transfer temporary copy back to original array
3406    var(nys:nyn,nxl:nxr) = var_tmp(nys:nyn,nxl:nxr)
3407
3408 END SUBROUTINE add_ghost_layers_2d_int8
3409
3410
3411!--------------------------------------------------------------------------------------------------!
3412! Description:
3413! ------------
3414!> Resize 32-bit 2D Integer array: (nys:nyn,nxl:nxr) -> (nysg:nyng,nxlg:nxrg)
3415!--------------------------------------------------------------------------------------------------!
3416 SUBROUTINE add_ghost_layers_2d_int32( var )
3417
3418    IMPLICIT NONE
3419
3420    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  var     !< treated variable
3421    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  var_tmp !< temporary copy
3422!
3423!-- Allocate temporary variable
3424    ALLOCATE( var_tmp(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3425!
3426!-- Temporary copy of the variable
3427    var_tmp(nys:nyn,nxl:nxr) = var(nys:nyn,nxl:nxr)
3428!
3429!-- Resize the array
3430    DEALLOCATE( var )
3431    ALLOCATE( var(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3432!
3433!-- Transfer temporary copy back to original array
3434    var(nys:nyn,nxl:nxr) = var_tmp(nys:nyn,nxl:nxr)
3435
3436 END SUBROUTINE add_ghost_layers_2d_int32
3437
3438!--------------------------------------------------------------------------------------------------!
3439! Description:
3440! ------------
3441!> Resize 2D float array: (nys:nyn,nxl:nxr) -> (nysg:nyng,nxlg:nxrg)
3442!--------------------------------------------------------------------------------------------------!
3443 SUBROUTINE add_ghost_layers_2d_real( var )
3444
3445    IMPLICIT NONE
3446
3447    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  var     !< treated variable
3448    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  var_tmp !< temporary copy
3449!
3450!-- Allocate temporary variable
3451    ALLOCATE( var_tmp(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3452!
3453!-- Temporary copy of the variable
3454    var_tmp(nys:nyn,nxl:nxr) = var(nys:nyn,nxl:nxr)
3455!
3456!-- Resize the array
3457    DEALLOCATE( var )
3458    ALLOCATE( var(nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3459!
3460!-- Transfer temporary copy back to original array
3461    var(nys:nyn,nxl:nxr) = var_tmp(nys:nyn,nxl:nxr)
3462
3463 END SUBROUTINE add_ghost_layers_2d_real
3464
3465
3466!--------------------------------------------------------------------------------------------------!
3467! Description:
3468! ------------
3469!> Resize 8-bit 3D Integer array: (:,nys:nyn,nxl:nxr) -> (:,nysg:nyng,nxlg:nxrg)
3470!--------------------------------------------------------------------------------------------------!
3471 SUBROUTINE add_ghost_layers_3d_int8( var, ks, ke )
3472
3473    IMPLICIT NONE
3474
3475    INTEGER(iwp) ::  ke  !< upper bound of treated array in z-direction
3476    INTEGER(iwp) ::  ks  !< lower bound of treated array in z-direction
3477
3478    INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE ::  var     !< treated variable
3479    INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE ::  var_tmp !< temporary copy
3480!
3481!-- Allocate temporary variable
3482    ALLOCATE( var_tmp(ks:ke,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3483!
3484!-- Temporary copy of the variable
3485    var_tmp(ks:ke,nys:nyn,nxl:nxr) = var(ks:ke,nys:nyn,nxl:nxr)
3486!
3487!-- Resize the array
3488    DEALLOCATE( var )
3489    ALLOCATE( var(ks:ke,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3490!
3491!-- Transfer temporary copy back to original array
3492    var(ks:ke,nys:nyn,nxl:nxr) = var_tmp(ks:ke,nys:nyn,nxl:nxr)
3493
3494 END SUBROUTINE add_ghost_layers_3d_int8
3495
3496
3497!--------------------------------------------------------------------------------------------------!
3498! Description:
3499! ------------
3500!> Resize 3D Real array: (:,nys:nyn,nxl:nxr) -> (:,nysg:nyng,nxlg:nxrg)
3501!--------------------------------------------------------------------------------------------------!
3502 SUBROUTINE add_ghost_layers_3d_real( var, ks, ke )
3503
3504    IMPLICIT NONE
3505
3506    INTEGER(iwp) ::  ke  !< upper bound of treated array in z-direction
3507    INTEGER(iwp) ::  ks  !< lower bound of treated array in z-direction
3508
3509    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  var     !< treated variable
3510    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  var_tmp !< temporary copy
3511!
3512!-- Allocate temporary variable
3513    ALLOCATE( var_tmp(ks:ke,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3514!
3515!-- Temporary copy of the variable
3516    var_tmp(ks:ke,nys:nyn,nxl:nxr) = var(ks:ke,nys:nyn,nxl:nxr)
3517!
3518!-- Resize the array
3519    DEALLOCATE( var )
3520    ALLOCATE( var(ks:ke,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3521!
3522!-- Transfer temporary copy back to original array
3523    var(ks:ke,nys:nyn,nxl:nxr) = var_tmp(ks:ke,nys:nyn,nxl:nxr)
3524
3525 END SUBROUTINE add_ghost_layers_3d_real
3526
3527
3528!--------------------------------------------------------------------------------------------------!
3529! Description:
3530! ------------
3531!> Resize 4D Real array: (:,:,nys:nyn,nxl:nxr) -> (:,nysg:nyng,nxlg:nxrg)
3532!--------------------------------------------------------------------------------------------------!
3533 SUBROUTINE add_ghost_layers_4d_real( var, k1s, k1e, k2s, k2e )
3534
3535    IMPLICIT NONE
3536
3537    INTEGER(iwp) ::  k1e  !< upper bound of treated array in z-direction
3538    INTEGER(iwp) ::  k1s  !< lower bound of treated array in z-direction
3539    INTEGER(iwp) ::  k2e  !< upper bound of treated array along parameter space
3540    INTEGER(iwp) ::  k2s  !< lower bound of treated array along parameter space
3541
3542    REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::  var      !< treated variable
3543    REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::  var_tmp  !< temporary copy
3544!
3545!-- Allocate temporary variable
3546    ALLOCATE( var_tmp(k1s:k1e,k2s:k2e,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3547!
3548!-- Temporary copy of the variable
3549    var_tmp(k1s:k1e,k2s:k2e,nys:nyn,nxl:nxr) = var(k1s:k1e,k2s:k2e,nys:nyn,nxl:nxr)
3550!
3551!-- Resize the array
3552    DEALLOCATE( var )
3553    ALLOCATE( var(k1s:k1e,k2s:k2e,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp) )
3554!
3555!-- Transfer temporary copy back to original array
3556    var(k1s:k1e,k2s:k2e,nys:nyn,nxl:nxr) = var_tmp(k1s:k1e,k2s:k2e,nys:nyn,nxl:nxr)
3557
3558 END SUBROUTINE add_ghost_layers_4d_real
3559
3560
3561!--------------------------------------------------------------------------------------------------!
3562! Description:
3563! ------------
3564!> Checks if a given variables is on file
3565!--------------------------------------------------------------------------------------------------!
3566 FUNCTION check_existence( vars_in_file, var_name )
3567
3568    IMPLICIT NONE
3569
3570    CHARACTER(LEN=*) ::  var_name                    !< variable to be checked
3571    CHARACTER(LEN=*), DIMENSION(:) ::  vars_in_file  !< list of variables in file
3572
3573    INTEGER(iwp) ::  i                              !< loop variable
3574
3575    LOGICAL ::  check_existence                     !< flag indicating whether a variable exist or not - actual return value
3576
3577    i = 1
3578    check_existence = .FALSE.
3579    DO  WHILE ( i <= SIZE( vars_in_file ) )
3580       check_existence = TRIM( vars_in_file(i) ) == TRIM( var_name )  .OR.  check_existence
3581       i = i + 1
3582    ENDDO
3583
3584    RETURN
3585
3586 END FUNCTION check_existence
3587
3588
3589!--------------------------------------------------------------------------------------------------!
3590! Description:
3591! ------------
3592!> Closes an existing netCDF file.
3593!--------------------------------------------------------------------------------------------------!
3594 SUBROUTINE close_input_file( id )
3595
3596    USE pegrid
3597
3598    IMPLICIT NONE
3599
3600    INTEGER(iwp), INTENT(INOUT) ::  id  !< file id
3601
3602#if defined( __netcdf )
3603    nc_stat = NF90_CLOSE( id )
3604    CALL handle_error( 'close', 540 )
3605#endif
3606
3607 END SUBROUTINE close_input_file
3608
3609
3610!--------------------------------------------------------------------------------------------------!
3611! Description:
3612! ------------
3613!> Opens an existing netCDF file for reading only and returns its id.
3614!--------------------------------------------------------------------------------------------------!
3615 SUBROUTINE open_read_file( filename, id )
3616
3617    USE control_parameters,                                                                        &
3618        ONLY: message_string
3619
3620    USE pegrid
3621
3622    IMPLICIT NONE
3623
3624    CHARACTER (LEN=*), INTENT(IN) ::  filename    !< filename
3625    INTEGER(iwp), INTENT(INOUT)   ::  id          !< file id
3626    LOGICAL                       ::  file_exists !< true if file exists
3627
3628#if defined( __netcdf )
3629!
3630!-- Check if requested file exists
3631    INQUIRE( FILE=filename, EXIST=file_exists )
3632
3633    IF ( .NOT. file_exists )  THEN
3634       WRITE( message_string, * ) 'Required input file "' // filename // '" does not exist!'
3635       CALL message( 'open_read_file', 'PA0516', 2, 2, 0, 6, 1 )
3636    ENDIF
3637
3638#if defined( __netcdf4_parallel )
3639!
3640!-- If __netcdf4_parallel is defined, parrallel NetCDF will be used.
3641    nc_stat = NF90_OPEN( filename, IOR( NF90_NOWRITE, NF90_MPIIO ), id, COMM = comm2d,             &
3642                         INFO = MPI_INFO_NULL )
3643!
3644!-- In case the previous open call fails, check for possible Netcdf 3 file, and open it. However,
3645!-- this case, disable parallel access.
3646    IF( nc_stat /= NF90_NOERR )  THEN
3647       nc_stat = NF90_OPEN( filename, NF90_NOWRITE, id )
3648       collective_read = .FALSE.
3649    ELSE
3650#if defined( __nec )
3651       collective_read = .FALSE.   ! collective read causes hang situations on NEC Aurora
3652#else
3653       collective_read = .TRUE.
3654#endif
3655    ENDIF
3656#else
3657!
3658!-- All MPI processes open the file and read it (but not in parallel).
3659    nc_stat = NF90_OPEN( filename, NF90_NOWRITE, id )
3660#endif
3661
3662    CALL handle_error( 'open_read_file', 539 )
3663
3664#endif
3665 END SUBROUTINE open_read_file
3666
3667
3668!--------------------------------------------------------------------------------------------------!
3669! Description:
3670! ------------
3671!> Reads global or variable-related attributes of type INTEGER (32-bit)
3672!--------------------------------------------------------------------------------------------------!
3673  SUBROUTINE get_attribute_int32( id, attribute_name, value, global, variable_name, no_abort )
3674
3675    USE pegrid
3676
3677    IMPLICIT NONE
3678
3679    CHARACTER(LEN=*)            ::  attribute_name   !< attribute name
3680    CHARACTER(LEN=*), OPTIONAL  ::  variable_name    !< variable name
3681
3682    INTEGER(iwp), INTENT(IN)    ::  id               !< file id
3683    INTEGER(iwp)                ::  id_var           !< variable id
3684    INTEGER(iwp), INTENT(INOUT) ::  value            !< read value
3685
3686    LOGICAL                       ::  check_error    !< flag indicating if handle_error shall be checked
3687    LOGICAL, INTENT(IN)           ::  global         !< flag indicating global attribute
3688    LOGICAL, INTENT(IN), OPTIONAL ::  no_abort       !< flag indicating if errors should be checked
3689#if defined( __netcdf )
3690
3691    IF ( PRESENT( no_abort ) )  THEN
3692       check_error = no_abort
3693    ELSE
3694       check_error = .TRUE.
3695    ENDIF
3696!
3697!-- Read global attribute
3698    IF ( global )  THEN
3699       nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3700       IF ( check_error)  CALL handle_error( 'get_attribute_int32 global', 522, attribute_name )
3701!
3702!-- Read attributes referring to a single variable. Therefore, first inquire variable id
3703    ELSE
3704       nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3705       IF ( check_error)  CALL handle_error( 'get_attribute_int32', 522, attribute_name )
3706       nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3707       IF ( check_error)  CALL handle_error( 'get_attribute_int32', 522, attribute_name )
3708    ENDIF
3709#endif
3710
3711 END SUBROUTINE get_attribute_int32
3712
3713
3714!--------------------------------------------------------------------------------------------------!
3715! Description:
3716! ------------
3717!> Reads global or variable-related attributes of type INTEGER (8-bit)
3718!--------------------------------------------------------------------------------------------------!
3719 SUBROUTINE get_attribute_int8( id, attribute_name, value, global, variable_name, no_abort )
3720
3721    USE pegrid
3722
3723    IMPLICIT NONE
3724
3725    CHARACTER(LEN=*)            ::  attribute_name   !< attribute name
3726    CHARACTER(LEN=*), OPTIONAL  ::  variable_name    !< variable name
3727
3728    INTEGER(iwp), INTENT(IN)    ::  id               !< file id
3729    INTEGER(iwp)                ::  id_var           !< variable id
3730    INTEGER(KIND=1), INTENT(INOUT) ::  value         !< read value
3731
3732    LOGICAL                       ::  check_error    !< flag indicating if handle_error shall be checked
3733    LOGICAL, INTENT(IN), OPTIONAL ::  no_abort       !< flag indicating if errors should be checked
3734    LOGICAL, INTENT(IN)           ::  global         !< flag indicating global attribute
3735#if defined( __netcdf )
3736
3737    IF ( PRESENT( no_abort ) )  THEN
3738       check_error = no_abort
3739    ELSE
3740       check_error = .TRUE.
3741    ENDIF
3742!
3743!-- Read global attribute
3744    IF ( global )  THEN
3745       nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3746       IF ( check_error)  CALL handle_error( 'get_attribute_int8 global', 523, attribute_name )
3747!
3748!-- Read attributes referring to a single variable. Therefore, first inquire variable id
3749    ELSE
3750       nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3751       IF ( check_error)  CALL handle_error( 'get_attribute_int8', 523, attribute_name )
3752       nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3753       IF ( check_error)  CALL handle_error( 'get_attribute_int8', 523, attribute_name )
3754    ENDIF
3755#endif
3756
3757 END SUBROUTINE get_attribute_int8
3758
3759
3760!--------------------------------------------------------------------------------------------------!
3761! Description:
3762! ------------
3763!> Reads global or variable-related attributes of type REAL
3764!--------------------------------------------------------------------------------------------------!
3765 SUBROUTINE get_attribute_real( id, attribute_name, value, global, variable_name, no_abort )
3766
3767    USE pegrid
3768
3769    IMPLICIT NONE
3770
3771    CHARACTER(LEN=*)            ::  attribute_name   !< attribute name
3772    CHARACTER(LEN=*), OPTIONAL  ::  variable_name    !< variable name
3773
3774    INTEGER(iwp), INTENT(IN)    ::  id               !< file id
3775    INTEGER(iwp)                ::  id_var           !< variable id
3776
3777    LOGICAL                       ::  check_error    !< flag indicating if handle_error shall be checked
3778    LOGICAL, INTENT(IN)           ::  global         !< flag indicating global attribute
3779    LOGICAL, INTENT(IN), OPTIONAL ::  no_abort       !< flag indicating if errors should be checked
3780
3781    REAL(wp), INTENT(INOUT)     ::  value            !< read value
3782#if defined( __netcdf )
3783
3784    IF ( PRESENT( no_abort ) )  THEN
3785       check_error = no_abort
3786    ELSE
3787       check_error = .TRUE.
3788    ENDIF
3789!
3790!-- Read global attribute
3791    IF ( global )  THEN
3792       nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3793       IF ( check_error)  CALL handle_error( 'get_attribute_real global', 524, attribute_name )
3794!
3795!-- Read attributes referring to a single variable. Therefore, first inquire variable id
3796    ELSE
3797       nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3798       IF ( check_error)  CALL handle_error( 'get_attribute_real', 524, attribute_name )
3799       nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3800       IF ( check_error)  CALL handle_error( 'get_attribute_real', 524, attribute_name )
3801    ENDIF
3802#endif
3803
3804 END SUBROUTINE get_attribute_real
3805
3806
3807!--------------------------------------------------------------------------------------------------!
3808! Description:
3809! ------------
3810!> Reads global or variable-related attributes of type CHARACTER
3811!> Remark: reading attributes of type NF_STRING return an error code 56 -
3812!> Attempt to convert between text & numbers.
3813!--------------------------------------------------------------------------------------------------!
3814 SUBROUTINE get_attribute_string( id, attribute_name, value, global, variable_name, no_abort )
3815
3816    USE pegrid
3817
3818    IMPLICIT NONE
3819
3820    CHARACTER(LEN=*)                ::  attribute_name   !< attribute name
3821    CHARACTER(LEN=*), INTENT(INOUT) ::  value            !< read value
3822    CHARACTER(LEN=*), OPTIONAL      ::  variable_name    !< variable name
3823
3824    INTEGER(iwp), INTENT(IN)    ::  id               !< file id
3825    INTEGER(iwp)                ::  id_var           !< variable id
3826
3827    LOGICAL ::  check_error                          !< flag indicating if handle_error shall be checked
3828    LOGICAL, INTENT(IN) ::  global                   !< flag indicating global attribute
3829    LOGICAL, INTENT(IN), OPTIONAL ::  no_abort       !< flag indicating if errors should be checked
3830#if defined( __netcdf )
3831
3832    IF ( PRESENT( no_abort ) )  THEN
3833       check_error = no_abort
3834    ELSE
3835       check_error = .TRUE.
3836    ENDIF
3837!
3838!-- Read global attribute
3839    IF ( global )  THEN
3840       nc_stat = NF90_GET_ATT( id, NF90_GLOBAL, TRIM( attribute_name ), value )
3841       IF ( check_error)  CALL handle_error( 'get_attribute_string global', 525, attribute_name )
3842!
3843!-- Read attributes referring to a single variable. Therefore, first inquire variable id
3844    ELSE
3845       nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3846       IF ( check_error)  CALL handle_error( 'get_attribute_string', 525, attribute_name )
3847
3848       nc_stat = NF90_GET_ATT( id, id_var, TRIM( attribute_name ), value )
3849       IF ( check_error)  CALL handle_error( 'get_attribute_string',525, attribute_name )
3850
3851    ENDIF
3852#endif
3853
3854 END SUBROUTINE get_attribute_string
3855
3856
3857!--------------------------------------------------------------------------------------------------!
3858! Description:
3859! ------------
3860!> Get dimension array for a given dimension
3861!--------------------------------------------------------------------------------------------------!
3862  SUBROUTINE get_dimension_length( id, dim_len, variable_name )
3863    USE pegrid
3864
3865    IMPLICIT NONE
3866
3867    CHARACTER(LEN=100)          ::  dum              !< dummy variable to receive return character
3868    CHARACTER(LEN=*)            ::  variable_name    !< dimension name
3869
3870    INTEGER(iwp)                ::  dim_len          !< dimension size
3871    INTEGER(iwp), INTENT(IN)    ::  id               !< file id
3872    INTEGER(iwp)                ::  id_dim           !< dimension id
3873
3874#if defined( __netcdf )
3875!
3876!-- First, inquire dimension ID
3877    nc_stat = NF90_INQ_DIMID( id, TRIM( variable_name ), id_dim )
3878    CALL handle_error( 'get_dimension_length', 526, variable_name )
3879!
3880!-- Inquire dimension length
3881    nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim, dum, LEN = dim_len )
3882    CALL handle_error( 'get_dimension_length', 526, variable_name )
3883
3884#endif
3885
3886 END SUBROUTINE get_dimension_length
3887
3888
3889!--------------------------------------------------------------------------------------------------!
3890! Description:
3891! ------------
3892!> Routine for reading-in a character string from the chem emissions netcdf input file.
3893!--------------------------------------------------------------------------------------------------!
3894 SUBROUTINE get_variable_string( id, variable_name, var_string, names_number )
3895
3896    USE indices
3897    USE pegrid
3898
3899    IMPLICIT NONE
3900
3901    CHARACTER (LEN=*)                                             :: variable_name    !< variable name
3902    CHARACTER (LEN=25), ALLOCATABLE, DIMENSION(:), INTENT(INOUT)  :: var_string
3903    CHARACTER (LEN=1), ALLOCATABLE, DIMENSION(:,:)                :: tmp_var_string   !< variable to be read
3904
3905    INTEGER(iwp)              :: i,j                    !< index to go through the length of the dimensions
3906    INTEGER(iwp), INTENT(IN)  :: id                     !< file id
3907    INTEGER(iwp)              :: id_var                 !< variable id
3908    INTEGER(iwp), INTENT(IN)  :: names_number           !< number of names
3909    INTEGER(iwp)              :: max_string_length=25   !< this is both the maximum length of a name, but also the number of the
3910                                                        !< components of the first dimensions (rows)
3911#if defined( __netcdf )
3912
3913    ALLOCATE( tmp_var_string(max_string_length,names_number) )
3914
3915    ALLOCATE( var_string(names_number) )
3916
3917!-- Inquire variable id
3918    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
3919
3920!
3921!-- Get variable
3922!-- Start cycle over the emission species
3923    DO  i = 1, names_number
3924!
3925!--    Read the first letter of each component
3926       nc_stat = NF90_GET_VAR( id, id_var, var_string(i), start = (/ 1,i /), count = (/ 1,1 /) )
3927       CALL handle_error( 'get_variable_string', 701 )
3928!
3929!--    Start cycle over charachters
3930       DO  j = 1, max_string_length
3931!
3932!--       Read the rest of the components of the name
3933          nc_stat = NF90_GET_VAR( id, id_var, tmp_var_string(j,i), start = (/ j,i /),              &
3934                                  count = (/ 1,1 /) )
3935          CALL handle_error( 'get_variable_string', 702 )
3936
3937          IF ( IACHAR( tmp_var_string(j,i) ) == 0 )  THEN
3938               tmp_var_string(j,i)=''
3939          ENDIF
3940
3941          IF ( j>1 )  THEN
3942!
3943!--          Concatenate first letter of the name and the others
3944             var_string(i)=TRIM( var_string(i) ) // TRIM( tmp_var_string(j,i) )
3945          ENDIF
3946       ENDDO
3947    ENDDO
3948
3949#endif
3950
3951 END SUBROUTINE get_variable_string
3952
3953
3954!
3955!--------------------------------------------------------------------------------------------------!
3956! Description:
3957! ------------
3958!> Generalized routine for reading strings from a netCDF variable to replace existing
3959!> get_variable_string ( )
3960!>
3961!> Improvements:
3962!>   - Expanded string size limit from 25 to 512
3963!>   - No longer removes spaces between text magically (this seems to have been aimed at a very
3964!>     specific application, but I don't know what)
3965!>   - Simplified implementation
3966!>
3967!> Caveats:
3968!>   - Somehow I could not get the subroutine to work with str_array(:,:) so I reverted to a
3969!>     hard-coded str_array(:,512), hopefully large enough for most general applications. This also
3970!>     means the character variable used for str_array must be of size (:,512)
3971!>     (ecc 20200128)
3972!--------------------------------------------------------------------------------------------------!
3973 SUBROUTINE get_variable_string_generic ( id, var_name, str_array, num_str, str_len )
3974
3975    IMPLICIT NONE
3976
3977    CHARACTER(LEN=*),                INTENT(IN)    :: var_name       !< netCDF variable name
3978    CHARACTER(LEN=512), ALLOCATABLE, INTENT(INOUT) :: str_array(:)   !< output string array
3979
3980    INTEGER(iwp)              :: buf_len   !< string buffer size
3981    INTEGER(iwp)              :: id_var    !< netCDF variable ID
3982    INTEGER(iwp), INTENT(IN)  :: id        !< netCDF file ID
3983    INTEGER(iwp)              :: k         !< generic counter
3984    INTEGER(iwp), INTENT(IN)  :: num_str   !< number of string elements in array
3985    INTEGER(iwp), INTENT(IN)  :: str_len   !< size of each string element
3986
3987#if defined( __netcdf )
3988
3989!
3990!-- Set buffer length to up to hard-coded string size
3991    buf_len = MIN( ABS(str_len), 512 )
3992
3993!
3994!-- Allocate necessary memories for string variables
3995    ALLOCATE( str_array(num_str) )
3996!
3997!-- Get variable id
3998    nc_stat = NF90_INQ_VARID( id, TRIM(var_name), id_var )
3999!
4000!-- Extract string variables
4001    DO  k = 1, num_str
4002       str_array(k) = ''
4003       nc_stat = NF90_GET_VAR( id, id_var, str_array(k), start = (/ 1, k /),                       &
4004                               count = (/ buf_len, 1 /)  )
4005       CALL handle_error ( 'get_variable_string_generic', 702 )
4006    ENDDO
4007
4008#endif
4009
4010 END SUBROUTINE get_variable_string_generic
4011
4012
4013!--------------------------------------------------------------------------------------------------!
4014!    Description:
4015!    ------------
4016!>   Reads a character variable in a 1D array
4017!--------------------------------------------------------------------------------------------------!
4018  SUBROUTINE get_variable_1d_char( id, variable_name, var )
4019
4020    USE pegrid
4021
4022    IMPLICIT NONE
4023
4024    CHARACTER(LEN=*)            ::  variable_name          !< variable name
4025
4026    CHARACTER(LEN=*), DIMENSION(:), INTENT(INOUT) ::  var  !< variable to be read
4027
4028    INTEGER(iwp)                ::  i                !< running index over variable dimension
4029    INTEGER(iwp), INTENT(IN)    ::  id               !< file id
4030    INTEGER(iwp)                ::  id_var           !< dimension id
4031
4032    INTEGER(iwp), DIMENSION(2)  ::  dimid            !< dimension IDs
4033    INTEGER(iwp), DIMENSION(2)  ::  dimsize          !< dimension size
4034
4035#if defined( __netcdf )
4036
4037!
4038!-- First, inquire variable ID
4039    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4040    CALL handle_error( 'get_variable_1d_int', 527, variable_name )
4041!
4042!-- Inquire dimension IDs
4043    nc_stat = NF90_INQUIRE_VARIABLE( id, id_var, dimids = dimid(1:2) )
4044    CALL handle_error( 'get_variable_1d_char', 527, variable_name )
4045!
4046!-- Read dimesnion length
4047    nc_stat = NF90_INQUIRE_DIMENSION( id, dimid(1), LEN = dimsize(1) )
4048    nc_stat = NF90_INQUIRE_DIMENSION( id, dimid(2), LEN = dimsize(2) )
4049
4050!
4051!-- Read character array. Note, each element is read individually, in order to better separate
4052!-- single strings.
4053    DO  i = 1, dimsize(2)
4054       nc_stat = NF90_GET_VAR( id, id_var, var(i),                                                 &
4055                               start = (/ 1, i /),                                                 &
4056                               count = (/ dimsize(1), 1 /) )
4057       CALL handle_error( 'get_variable_1d_char', 527, variable_name )
4058    ENDDO
4059
4060#endif
4061
4062 END SUBROUTINE get_variable_1d_char
4063
4064
4065!--------------------------------------------------------------------------------------------------!
4066! Description:
4067! ------------
4068!> Reads a 1D integer variable from file.
4069!--------------------------------------------------------------------------------------------------!
4070  SUBROUTINE get_variable_1d_int( id, variable_name, var )
4071
4072    USE pegrid
4073
4074    IMPLICIT NONE
4075
4076    CHARACTER(LEN=*)            ::  variable_name    !< variable name
4077
4078    INTEGER(iwp), INTENT(IN)    ::  id               !< file id
4079    INTEGER(iwp)                ::  id_var           !< dimension id
4080
4081    INTEGER(iwp), DIMENSION(:), INTENT(INOUT) ::  var  !< variable to be read
4082#if defined( __netcdf )
4083
4084!
4085!-- First, inquire variable ID
4086    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4087    CALL handle_error( 'get_variable_1d_int', 527, variable_name )
4088!
4089!-- Read variable
4090    nc_stat = NF90_GET_VAR( id, id_var, var )
4091    CALL handle_error( 'get_variable_1d_int', 527, variable_name )
4092
4093#endif
4094
4095 END SUBROUTINE get_variable_1d_int
4096
4097
4098!--------------------------------------------------------------------------------------------------!
4099! Description:
4100! ------------
4101!> Reads a 1D float variable from file.
4102!--------------------------------------------------------------------------------------------------!
4103  SUBROUTINE get_variable_1d_real( id, variable_name, var, is, count_elements )
4104
4105    USE pegrid
4106
4107    IMPLICIT NONE
4108
4109    CHARACTER(LEN=*)            ::  variable_name    !< variable name
4110
4111    INTEGER(iwp), INTENT(IN), OPTIONAL ::  count_elements   !< number of elements to be read
4112    INTEGER(iwp), INTENT(IN)           ::  id               !< file id
4113    INTEGER(iwp)                       ::  id_var           !< dimension id
4114    INTEGER(iwp), INTENT(IN), OPTIONAL ::  is               !< start index
4115
4116    REAL(wp), DIMENSION(:), INTENT(INOUT) ::  var    !< variable to be read
4117#if defined( __netcdf )
4118!
4119!-- First, inquire variable ID
4120    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4121    CALL handle_error( 'get_variable_1d_real', 528, variable_name )
4122!
4123!-- Read variable
4124    IF ( PRESENT( is ) )  THEN
4125       nc_stat = NF90_GET_VAR( id, id_var, var, start = (/ is /), count = (/ count_elements /) )
4126       CALL handle_error( 'get_variable_1d_real', 528, variable_name )
4127    ELSE
4128       nc_stat = NF90_GET_VAR( id, id_var, var )
4129       CALL handle_error( 'get_variable_1d_real', 528, variable_name )
4130    ENDIF
4131
4132#endif
4133
4134 END SUBROUTINE get_variable_1d_real
4135
4136
4137!--------------------------------------------------------------------------------------------------!
4138! Description:
4139! ------------
4140!> Reads a time-dependent 1D float variable from file.
4141!--------------------------------------------------------------------------------------------------!
4142 SUBROUTINE get_variable_pr( id, variable_name, t, var )
4143
4144    USE pegrid
4145
4146    IMPLICIT NONE
4147
4148    CHARACTER(LEN=*)                      ::  variable_name    !< variable name
4149
4150    INTEGER(iwp), INTENT(IN)              ::  id               !< file id
4151    INTEGER(iwp), DIMENSION(1:2)          ::  id_dim           !< dimension ids
4152    INTEGER(iwp)                          ::  id_var           !< dimension id
4153    INTEGER(iwp)                          ::  n_file           !< number of data-points in file along z dimension
4154    INTEGER(iwp), INTENT(IN)              ::  t                !< timestep number
4155
4156    REAL(wp), DIMENSION(:), INTENT(INOUT) ::  var              !< variable to be read
4157
4158#if defined( __netcdf )
4159!
4160!-- First, inquire variable ID
4161    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4162!
4163!-- Inquire dimension size of vertical dimension
4164    nc_stat = NF90_INQUIRE_VARIABLE( id, id_var, DIMIDS = id_dim )
4165    nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim(1), LEN = n_file )
4166!
4167!-- Read variable.
4168    nc_stat = NF90_GET_VAR( id, id_var, var,                                                       &
4169                            start = (/ 1,      t     /),                                           &
4170                            count = (/ n_file, 1     /) )
4171    CALL handle_error( 'get_variable_pr', 529, variable_name )
4172#endif
4173
4174 END SUBROUTINE get_variable_pr
4175
4176
4177!--------------------------------------------------------------------------------------------------!
4178! Description:
4179! ------------
4180!> Reads a per-surface pars variable from file. Because all surfaces are stored as flat 1-D array,
4181!> each PE has to scan the data and find the surface indices belonging to its subdomain. During this
4182!> scan, it also builds a necessary (j,i) index.
4183!--------------------------------------------------------------------------------------------------!
4184 SUBROUTINE get_variable_surf( id, variable_name, surf )
4185
4186    USE pegrid
4187
4188    USE indices,                                                                                   &
4189        ONLY:  nxl, nxr, nys, nyn
4190
4191    USE control_parameters,                                                                        &
4192        ONLY: dz, message_string
4193
4194    USE grid_variables,                                                                            &
4195        ONLY: ddx, ddy
4196
4197    USE basic_constants_and_equations_mod,                                                         &
4198        ONLY: pi
4199
4200    IMPLICIT NONE
4201
4202    INTEGER(iwp), PARAMETER                   ::  nsurf_pars_read = 2**15  !< read buffer size (value > 10^15 makes problem with
4203                                                                           !< ifort)
4204
4205    CHARACTER(LEN=*)                          ::  variable_name !< variable name
4206
4207    INTEGER(iwp)                              ::  i             !< grid index in x-direction
4208    INTEGER(iwp)                              ::  j             !< grid index in y-direction
4209    INTEGER(iwp), INTENT(IN)                  ::  id            !< file id
4210    INTEGER(iwp)                              ::  id_azimuth    !< azimuth variable id
4211    INTEGER(iwp)                              ::  id_var        !< variable id
4212    INTEGER(iwp)                              ::  id_xs         !< xs variable id
4213    INTEGER(iwp)                              ::  id_ys         !< ys variable id
4214    INTEGER(iwp)                              ::  id_zs         !< zs variable id
4215    INTEGER(iwp)                              ::  id_zenith     !< zeith variable id
4216    INTEGER(iwp)                              ::  is            !< local surface index
4217    INTEGER(iwp)                              ::  is0, isc      !< read surface start and count
4218    INTEGER(iwp)                              ::  isurf         !< netcdf surface index
4219    INTEGER(iwp)                              ::  nsurf         !< total number of surfaces in file
4220
4221    INTEGER(iwp), DIMENSION(6)                ::  coords        !< integer coordinates of surface location
4222    INTEGER(iwp), DIMENSION(2)                ::  id_dim        !< dimension ids
4223
4224    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  nsurf_ji      !< numbers of surfaces by coords
4225
4226    REAL(wp)                                  ::  oro_max_l     !< maximum terrain height under building
4227
4228    REAL(wp), DIMENSION(:), ALLOCATABLE       ::  azimuth       !< read buffer for azimuth(s)
4229    REAL(wp), DIMENSION(:), ALLOCATABLE       ::  zenith        !< read buffer for zenith(s)
4230    REAL(wp), DIMENSION(:), ALLOCATABLE       ::  xs            !< surface coordinate array of x-dimension
4231    REAL(wp), DIMENSION(:), ALLOCATABLE       ::  ys            !< surface coordinate array of y-dimension
4232    REAL(wp), DIMENSION(:), ALLOCATABLE       ::  zs            !< surface coordinate array of z-dimension
4233
4234    REAL(wp), DIMENSION(:,:), ALLOCATABLE     ::  pars_read     !< read buffer for the building parameters
4235
4236    TYPE(pars_surf)                           ::  surf          !< parameters variable to be loaded
4237
4238
4239#if defined( __netcdf )
4240!
4241!-- First, inquire variable ID's
4242    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4243    nc_stat = NF90_INQ_VARID( id, 'zs',                  id_zs )
4244    nc_stat = NF90_INQ_VARID( id, 'ys',                  id_ys )
4245    nc_stat = NF90_INQ_VARID( id, 'xs',                  id_xs )
4246    nc_stat = NF90_INQ_VARID( id, 'zenith',              id_zenith )
4247    nc_stat = NF90_INQ_VARID( id, 'azimuth',             id_azimuth )
4248!
4249!-- Inquire dimension sizes for the number of surfaces and parameters given
4250    nc_stat = NF90_INQUIRE_VARIABLE( id, id_var, DIMIDS = id_dim )
4251    nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim(1), LEN = nsurf )
4252    nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim(2), LEN = surf%np )
4253
4254    ALLOCATE( pars_read( nsurf_pars_read, surf%np ),                                               &
4255              zs(nsurf_pars_read), ys(nsurf_pars_read),                                            &
4256              xs(nsurf_pars_read), zenith(nsurf_pars_read),                                        &
4257              azimuth(nsurf_pars_read),                                                            &
4258              nsurf_ji(nys:nyn,nxl:nxr) )
4259
4260    nsurf_ji(:,:) = 0
4261!
4262!-- Scan surface coordinates, count locally
4263    is0 = 1
4264    DO
4265       isc = MIN(nsurf_pars_read, nsurf - is0 + 1)
4266       IF ( isc <= 0 )  EXIT
4267
4268       nc_stat = NF90_GET_VAR( id, id_ys, ys,                                                      &
4269                               start = (/ is0 /),                                                  &
4270                               count = (/ isc /) )
4271       nc_stat = NF90_GET_VAR( id, id_xs, xs,                                                      &
4272                               start = (/ is0 /),                                                  &
4273                               count = (/ isc /) )
4274       nc_stat = NF90_GET_VAR( id, id_zenith, zenith,                                              &
4275                               start = (/ is0 /), &
4276                               count = (/ isc /) )
4277       nc_stat = NF90_GET_VAR( id, id_azimuth, azimuth,                                            &
4278                               start = (/ is0 /),                                                  &
4279                               count = (/ isc /) )
4280       CALL handle_error( 'get_variable_surf', 682, 'azimuth' )
4281
4282       DO  isurf = 1, isc
4283!
4284!--       Parse coordinates, detect if belongs to subdomain
4285          coords = transform_coords( xs(isurf), ys(isurf), zenith(isurf), azimuth(isurf) )
4286          IF ( coords(2) < nys  .OR.  coords(2) > nyn  .OR.                                        &
4287               coords(3) < nxl  .OR.  coords(3) > nxr )  CYCLE
4288
4289          nsurf_ji(coords(2),coords(3)) = nsurf_ji(coords(2),coords(3)) + 1
4290       ENDDO
4291       is0 = is0 + isc
4292    ENDDO
4293!
4294!-- Populate reverse index from surface counts
4295    ALLOCATE( surf%index_ji( 2,nys:nyn,nxl:nxr ) )
4296    isurf = 1
4297    DO  j = nys, nyn
4298       DO  i = nxl, nxr
4299          surf%index_ji(:,j,i) = (/ isurf, isurf + nsurf_ji(j,i) - 1 /)
4300          isurf = isurf + nsurf_ji(j,i)
4301       ENDDO
4302    ENDDO
4303
4304    surf%nsurf = isurf - 1
4305    ALLOCATE( surf%pars(0:surf%np-1,surf%nsurf),                                                   &
4306              surf%coords(6,surf%nsurf) )
4307!
4308!-- Scan surfaces again, saving pars into allocated structures
4309    nsurf_ji(:,:) = 0
4310    is0 = 1
4311    DO
4312       isc = MIN(nsurf_pars_read, nsurf - is0 + 1)
4313       IF ( isc <= 0 )  EXIT
4314
4315       nc_stat = NF90_GET_VAR( id, id_var, pars_read(1:isc, 1:surf%np),                            &
4316                               start = (/ is0, 1       /),                                         &
4317                               count = (/ isc, surf%np /) )
4318       CALL handle_error( 'get_variable_surf', 683, variable_name )
4319
4320       nc_stat = NF90_GET_VAR( id, id_zs, zs,                                                      &
4321                               start = (/ is0 /),                                                  &
4322                               count = (/ isc /) )
4323       nc_stat = NF90_GET_VAR( id, id_ys, ys,                                                      &
4324                               start = (/ is0 /),                                                  &
4325                               count = (/ isc /) )
4326       nc_stat = NF90_GET_VAR( id, id_xs, xs,                                                      &
4327                               start = (/ is0 /),                                                  &
4328                               count = (/ isc /) )
4329       nc_stat = NF90_GET_VAR( id, id_zenith, zenith,                                              &
4330                               start = (/ is0 /),                                                  &
4331                               count = (/ isc /) )
4332       nc_stat = NF90_GET_VAR( id, id_azimuth, azimuth,                                            &
4333                               start = (/ is0 /),                                                  &
4334                               count = (/ isc /) )
4335
4336       DO  isurf = 1, isc
4337!
4338!--       Parse coordinates, detect if belongs to subdomain
4339          coords = transform_coords( xs(isurf), ys(isurf), zenith(isurf), azimuth(isurf) )
4340          IF ( coords(2) < nys  .OR.  coords(2) > nyn  .OR.                                        &
4341               coords(3) < nxl  .OR.  coords(3) > nxr )  CYCLE
4342!
4343!--       Determine maximum terrain under building (base z-coordinate). Using normal vector to
4344!--       locate building inner coordinates.
4345          oro_max_l = buildings_f%oro_max(coords(2)-coords(5),coords(3)-coords(6))
4346          IF ( oro_max_l == buildings_f%fill1 )  THEN
4347             WRITE( message_string, * ) 'Found building surface on '   //                          &
4348                                        'non-building coordinates (xs, ys, zenith, azimuth): ',    &
4349                                        xs(isurf), ys(isurf), zenith(isurf), azimuth(isurf)
4350             CALL message( 'get_variable_surf', 'PA0684', 2, 2, myid, 6, 0 )
4351          ENDIF
4352!
4353!--       Urban layer has no stretching, therefore using dz(1) instead of linear searching through
4354!--       zu/zw.
4355          coords(1) = NINT( ( zs(isurf) + oro_max_l ) / dz(1) + 0.5_wp + 0.5_wp * coords(4),       &
4356                            KIND=iwp )
4357!
4358!--       Save surface entry
4359          is = surf%index_ji(1,coords(2),coords(3)) + nsurf_ji(coords(2),coords(3))
4360          surf%pars(:,is) = pars_read(isurf,:)
4361          surf%coords(:,is) = coords(:)
4362
4363          nsurf_ji(coords(2),coords(3)) = nsurf_ji(coords(2),coords(3)) + 1
4364       ENDDO
4365
4366       is0 = is0 + isc
4367    ENDDO
4368
4369    DEALLOCATE( pars_read, zs, ys, xs, zenith, azimuth, nsurf_ji )
4370
4371 CONTAINS
4372
4373    PURE FUNCTION transform_coords( x, y, zenith, azimuth )
4374
4375       INTEGER(iwp), DIMENSION(6) ::  transform_coords !< (k,j,i,norm_z,norm_y,norm_x)
4376
4377       REAL(wp), INTENT(in) ::  azimuth  !< surface normal azimuth angle in degrees
4378       REAL(wp), INTENT(in) ::  x        !< surface centre coordinate in x in metres from origin
4379       REAL(wp), INTENT(in) ::  y        !< surface centre coordinate in y in metres from origin
4380       REAL(wp), INTENT(in) ::  zenith   !< surface normal zenith angle in degrees
4381
4382       transform_coords(4) = NINT( COS( zenith * pi / 180.0_wp ), KIND=iwp )
4383       IF ( transform_coords(4) == 0 )  THEN
4384          transform_coords(5) = NINT( COS( azimuth * pi / 180.0_wp ), KIND=iwp)
4385          transform_coords(6) = NINT( SIN( azimuth * pi / 180.0_wp ), KIND=iwp)
4386       ELSE
4387          transform_coords(5) = 0.0_wp
4388          transform_coords(6) = 0.0_wp
4389       ENDIF
4390
4391       transform_coords(1) = -999.0_wp ! not calculated here
4392       transform_coords(2) = NINT( y * ddy - 0.5_wp + 0.5_wp * transform_coords(5), KIND=iwp )
4393       transform_coords(3) = NINT( x * ddx - 0.5_wp + 0.5_wp * transform_coords(6), KIND=iwp )
4394
4395    END FUNCTION transform_coords
4396#endif
4397
4398 END SUBROUTINE get_variable_surf
4399
4400
4401!--------------------------------------------------------------------------------------------------!
4402! Description:
4403! ------------
4404!> Reads a 2D REAL variable from a file. Reading is done processor-wise, i.e. each core reads its
4405!> own domain in slices along x.
4406!--------------------------------------------------------------------------------------------------!
4407 SUBROUTINE get_variable_2d_real( id, variable_name, var, is, ie, js, je )
4408
4409    USE indices
4410    USE pegrid
4411
4412    IMPLICIT NONE
4413
4414    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4415
4416    INTEGER(iwp)                  ::  i               !< running index along x direction
4417    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4418    INTEGER(iwp)                  ::  id_var          !< variable id
4419    INTEGER(iwp)                  ::  ie              !< start index for subdomain input along x direction
4420    INTEGER(iwp)                  ::  is              !< end index for subdomain input along x direction
4421    INTEGER(iwp)                  ::  j               !< running index along y direction
4422    INTEGER(iwp)                  ::  je              !< start index for subdomain input along y direction
4423    INTEGER(iwp)                  ::  js              !< end index for subdomain input along y direction
4424
4425    REAL(wp), DIMENSION(:,:), ALLOCATABLE   ::  tmp   !< temporary variable to read data from file according
4426                                                      !< to its reverse memory access
4427    REAL(wp), DIMENSION(:,:), INTENT(INOUT) ::  var   !< variable to be read
4428
4429
4430#if defined( __netcdf )
4431!
4432!-- Inquire variable id
4433    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4434!
4435!-- Check for collective read-operation and set respective NetCDF flags if required.
4436    IF ( collective_read )  THEN
4437#if defined( __netcdf4_parallel )
4438       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4439#endif
4440    ENDIF
4441
4442!
4443!-- Allocate temporary variable according to memory access on file.
4444    ALLOCATE( tmp(is:ie,js:je) )
4445!
4446!-- Get variable
4447    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ is+1, js+1 /),                             &
4448                            count = (/ ie-is+1, je-js+1 /) )
4449    CALL handle_error( 'get_variable_2d_real', 530, variable_name )
4450!
4451!-- Resort data. Please note, dimension subscripts of var all start at 1.
4452    DO  i = is, ie
4453       DO  j = js, je
4454          var(j-js+1,i-is+1) = tmp(i,j)
4455       ENDDO
4456    ENDDO
4457
4458    DEALLOCATE( tmp )
4459#endif
4460
4461 END SUBROUTINE get_variable_2d_real
4462
4463
4464!--------------------------------------------------------------------------------------------------!
4465! Description:
4466! ------------
4467!> Reads a 2D 32-bit INTEGER variable from file. Reading is done processor-wise, i.e. each core
4468!> reads its own domain in slices along x.
4469!--------------------------------------------------------------------------------------------------!
4470 SUBROUTINE get_variable_2d_int32( id, variable_name, var, is, ie, js, je )
4471
4472    USE indices
4473    USE pegrid
4474
4475    IMPLICIT NONE
4476
4477    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4478
4479    INTEGER(iwp)                  ::  i               !< running index along x direction
4480    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4481    INTEGER(iwp)                  ::  id_var          !< variable id
4482    INTEGER(iwp)                  ::  ie              !< start index for subdomain input along x direction
4483    INTEGER(iwp)                  ::  is              !< end index for subdomain input along x direction
4484    INTEGER(iwp)                  ::  j               !< running index along y direction
4485    INTEGER(iwp)                  ::  je              !< start index for subdomain input along y direction
4486    INTEGER(iwp)                  ::  js              !< end index for subdomain input along y direction
4487
4488    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE   ::  tmp  !< temporary variable to read data from file according
4489                                                         !< to its reverse memory access
4490    INTEGER(iwp), DIMENSION(:,:), INTENT(INOUT) ::  var  !< variable to be read
4491
4492
4493#if defined( __netcdf )
4494!
4495!-- Inquire variable id
4496    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4497!
4498!-- Check for collective read-operation and set respective NetCDF flags if required.
4499    IF ( collective_read )  THEN
4500#if defined( __netcdf4_parallel )
4501       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4502#endif
4503    ENDIF
4504!
4505!-- Allocate temporary variable according to memory access on file.
4506    ALLOCATE( tmp(is:ie,js:je) )
4507!
4508!-- Get variable
4509    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ is+1, js+1 /),                             &
4510                            count = (/ ie-is+1, je-js+1 /) )
4511
4512    CALL handle_error( 'get_variable_2d_int32', 531, variable_name )
4513!
4514!-- Resort data. Please note, dimension subscripts of var all start at 1.
4515    DO  i = is, ie
4516       DO  j = js, je
4517          var(j-js+1,i-is+1) = tmp(i,j)
4518       ENDDO
4519    ENDDO
4520
4521    DEALLOCATE( tmp )
4522#endif
4523
4524 END SUBROUTINE get_variable_2d_int32
4525
4526
4527!--------------------------------------------------------------------------------------------------!
4528! Description:
4529! ------------
4530!> Reads a 2D 8-bit INTEGER variable from file. Reading is done processor-wise, i.e. each core reads
4531!> its own domain in slices along x.
4532!--------------------------------------------------------------------------------------------------!
4533 SUBROUTINE get_variable_2d_int8( id, variable_name, var, is, ie, js, je )
4534
4535    USE indices
4536    USE pegrid
4537
4538    IMPLICIT NONE
4539
4540    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4541
4542    INTEGER(iwp)                  ::  i               !< running index along x direction
4543    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4544    INTEGER(iwp)                  ::  id_var          !< variable id
4545    INTEGER(iwp)                  ::  ie              !< start index for subdomain input along x direction
4546    INTEGER(iwp)                  ::  is              !< end index for subdomain input along x direction
4547    INTEGER(iwp)                  ::  j               !< running index along y direction
4548    INTEGER(iwp)                  ::  je              !< start index for subdomain input along y direction
4549    INTEGER(iwp)                  ::  js              !< end index for subdomain input along y direction
4550
4551    INTEGER(KIND=1), DIMENSION(:,:), ALLOCATABLE   ::  tmp  !< temporary variable to read data from file according
4552                                                            !< to its reverse memory access
4553    INTEGER(KIND=1), DIMENSION(:,:), INTENT(INOUT) ::  var  !< variable to be read
4554
4555
4556#if defined( __netcdf )
4557!
4558!-- Inquire variable id
4559    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4560!
4561!-- Check for collective read-operation and set respective NetCDF flags if required.
4562    IF ( collective_read )  THEN
4563#if defined( __netcdf4_parallel )
4564       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4565#endif
4566    ENDIF
4567!
4568!-- Allocate temporary variable according to memory access on file.
4569    ALLOCATE( tmp(is:ie,js:je) )
4570!
4571!-- Get variable
4572    nc_stat = NF90_GET_VAR( id, id_var, tmp,                                                       &
4573                            start = (/ is+1,      js+1 /),                                         &
4574                            count = (/ ie-is + 1, je-js+1 /) )
4575
4576    CALL handle_error( 'get_variable_2d_int8', 532, variable_name )
4577!
4578!-- Resort data. Please note, dimension subscripts of var all start at 1.
4579    DO  i = is, ie
4580       DO  j = js, je
4581          var(j-js+1,i-is+1) = tmp(i,j)
4582       ENDDO
4583    ENDDO
4584
4585    DEALLOCATE( tmp )
4586#endif
4587
4588 END SUBROUTINE get_variable_2d_int8
4589
4590
4591!--------------------------------------------------------------------------------------------------!
4592! Description:
4593! ------------
4594!> Reads a 3D 8-bit INTEGER variable from file.
4595!--------------------------------------------------------------------------------------------------!
4596 SUBROUTINE get_variable_3d_int8( id, variable_name, var, is, ie, js, je, ks, ke )
4597
4598    USE indices
4599    USE pegrid
4600
4601    IMPLICIT NONE
4602
4603    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4604
4605    INTEGER(iwp)                  ::  i               !< index along x direction
4606    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4607    INTEGER(iwp)                  ::  id_var          !< variable id
4608    INTEGER(iwp)                  ::  ie              !< start index for subdomain input along x direction
4609    INTEGER(iwp)                  ::  is              !< end index for subdomain input along x direction
4610    INTEGER(iwp)                  ::  j               !< index along y direction
4611    INTEGER(iwp)                  ::  je              !< start index for subdomain input along y direction
4612    INTEGER(iwp)                  ::  js              !< end index for subdomain input along y direction
4613    INTEGER(iwp)                  ::  k               !< index along any 3rd dimension
4614    INTEGER(iwp)                  ::  ke              !< start index of 3rd dimension
4615    INTEGER(iwp)                  ::  ks              !< end index of 3rd dimension
4616
4617    INTEGER(KIND=1), DIMENSION(:,:,:), ALLOCATABLE   ::  tmp  !< temporary variable to read data from file according
4618                                                              !< to its reverse memory access
4619
4620    INTEGER(KIND=1), DIMENSION(:,:,:), INTENT(INOUT) ::  var  !< variable to be read
4621#if defined( __netcdf )
4622
4623!
4624!-- Inquire variable id
4625    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4626!
4627!-- Check for collective read-operation and set respective NetCDF flags if required.
4628    IF ( collective_read )  THEN
4629#if defined( __netcdf4_parallel )
4630       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4631#endif
4632    ENDIF
4633!
4634!-- Allocate temporary variable according to memory access on file.
4635    ALLOCATE( tmp(is:ie,js:je,ks:ke) )
4636!
4637!-- Get variable
4638    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ is+1, js+1, ks+1 /),                       &
4639                            count = (/ ie-is+1, je-js+1, ke-ks+1 /) )
4640
4641    CALL handle_error( 'get_variable_3d_int8', 533, variable_name )
4642!
4643!-- Resort data. Please note, dimension subscripts of var all start at 1.
4644    DO  i = is, ie
4645       DO  j = js, je
4646          DO  k = ks, ke
4647             var(k-ks+1,j-js+1,i-is+1) = tmp(i,j,k)
4648          ENDDO
4649       ENDDO
4650    ENDDO
4651
4652    DEALLOCATE( tmp )
4653#endif
4654
4655 END SUBROUTINE get_variable_3d_int8
4656
4657
4658!--------------------------------------------------------------------------------------------------!
4659! Description:
4660! ------------
4661!> Reads a 3D float variable from file.
4662!--------------------------------------------------------------------------------------------------!
4663 SUBROUTINE get_variable_3d_real( id, variable_name, var, is, ie, js, je, ks, ke )
4664
4665    USE indices
4666    USE pegrid
4667
4668    IMPLICIT NONE
4669
4670    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4671
4672    INTEGER(iwp)                  ::  i               !< index along x direction
4673    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4674    INTEGER(iwp)                  ::  id_var          !< variable id
4675    INTEGER(iwp)                  ::  ie              !< start index for subdomain input along x direction
4676    INTEGER(iwp)                  ::  is              !< end index for subdomain input along x direction
4677    INTEGER(iwp)                  ::  j               !< index along y direction
4678    INTEGER(iwp)                  ::  je              !< start index for subdomain input along y direction
4679    INTEGER(iwp)                  ::  js              !< end index for subdomain input along y direction
4680    INTEGER(iwp)                  ::  k               !< index along any 3rd dimension
4681    INTEGER(iwp)                  ::  ke              !< start index of 3rd dimension
4682    INTEGER(iwp)                  ::  ks              !< end index of 3rd dimension
4683
4684    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE   ::  tmp !< temporary variable to read data from file according
4685                                                      !< to its reverse memory access
4686
4687    REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) ::  var !< variable to be read
4688#if defined( __netcdf )
4689
4690!
4691!-- Inquire variable id
4692    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4693!
4694!-- Check for collective read-operation and set respective NetCDF flags if required.
4695    IF ( collective_read )  THEN
4696#if defined( __netcdf4_parallel )
4697       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4698#endif
4699    ENDIF
4700!
4701!-- Allocate temporary variable according to memory access on file.
4702    ALLOCATE( tmp(is:ie,js:je,ks:ke) )
4703!
4704!-- Get variable
4705    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ is+1, js+1, ks+1 /),                       &
4706                            count = (/ ie-is+1, je-js+1, ke-ks+1 /) )
4707
4708    CALL handle_error( 'get_variable_3d_real', 534, variable_name )
4709!
4710!-- Resort data. Please note, dimension subscripts of var all start at 1.
4711    DO  i = is, ie
4712       DO  j = js, je
4713          DO  k = ks, ke
4714             var(k-ks+1,j-js+1,i-is+1) = tmp(i,j,k)
4715          ENDDO
4716       ENDDO
4717    ENDDO
4718
4719    DEALLOCATE( tmp )
4720#endif
4721
4722 END SUBROUTINE get_variable_3d_real
4723
4724
4725!--------------------------------------------------------------------------------------------------!
4726! Description:
4727! ------------
4728!> Reads a 4D float variable from file.
4729!--------------------------------------------------------------------------------------------------!
4730 SUBROUTINE get_variable_4d_real( id, variable_name, var, is, ie, js, je, k1s, k1e, k2s, k2e )
4731
4732    USE indices
4733    USE pegrid
4734
4735    IMPLICIT NONE
4736
4737    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4738
4739    INTEGER(iwp)                  ::  i               !< index along x direction
4740    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4741    INTEGER(iwp)                  ::  id_var          !< variable id
4742    INTEGER(iwp)                  ::  ie              !< start index for subdomain input along x direction
4743    INTEGER(iwp)                  ::  is              !< end index for subdomain input along x direction
4744    INTEGER(iwp)                  ::  j               !< index along y direction
4745    INTEGER(iwp)                  ::  je              !< start index for subdomain input along y direction
4746    INTEGER(iwp)                  ::  js              !< end index for subdomain input along y direction
4747    INTEGER(iwp)                  ::  k1              !< index along 3rd direction
4748    INTEGER(iwp)                  ::  k1e             !< start index for 3rd dimension
4749    INTEGER(iwp)                  ::  k1s             !< end index for 3rd dimension
4750    INTEGER(iwp)                  ::  k2              !< index along 4th direction
4751    INTEGER(iwp)                  ::  k2e             !< start index for 4th dimension
4752    INTEGER(iwp)                  ::  k2s             !< end index for 4th dimension
4753
4754    REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE   ::  tmp  !< temporary variable to read data from file according
4755                                                         !< to its reverse memory access
4756    REAL(wp), DIMENSION(:,:,:,:), INTENT(INOUT) ::  var  !< variable to be read
4757
4758
4759#if defined( __netcdf )
4760
4761!
4762!-- Inquire variable id
4763    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4764!
4765!-- Check for collective read-operation and set respective NetCDF flags if required.
4766    IF ( collective_read )  THEN
4767#if defined( __netcdf4_parallel )
4768       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4769#endif
4770    ENDIF
4771
4772!
4773!-- Allocate temporary variable according to memory access on file.
4774    ALLOCATE( tmp(is:ie,js:je,k1s:k1e,k2s:k2e) )
4775!
4776!-- Get variable
4777    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ is+1, js+1, k1s+1, k2s+1 /),               &
4778                            count = (/ ie-is+1, je-js+1, k1e-k1s+1, k2e-k2s+1 /) )
4779
4780       CALL handle_error( 'get_variable_4d_real', 535, variable_name )
4781!
4782!-- Resort data. Please note, dimension subscripts of var all start at 1.
4783    DO  i = is, ie
4784       DO  j = js, je
4785          DO  k1 = k1s, k1e
4786             DO  k2 = k2s, k2e
4787                var(k2-k2s+1,k1-k1s+1,j-js+1,i-is+1) = tmp(i,j,k1,k2)
4788             ENDDO
4789          ENDDO
4790       ENDDO
4791    ENDDO
4792
4793    DEALLOCATE( tmp )
4794#endif
4795
4796 END SUBROUTINE get_variable_4d_real
4797
4798
4799!--------------------------------------------------------------------------------------------------!
4800! Description:
4801! ------------
4802!> Reads a 4D float variable from file and store it to a 3-d variable.
4803!--------------------------------------------------------------------------------------------------!
4804 SUBROUTINE get_variable_4d_to_3d_real( id, variable_name, var, ns, is, ie, js, je, ks, ke )
4805
4806    USE indices
4807    USE pegrid
4808
4809    IMPLICIT NONE
4810
4811    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4812
4813    INTEGER(iwp)                  ::  i               !< index along x direction
4814    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4815    INTEGER(iwp)                  ::  id_var          !< variable id
4816    INTEGER(iwp)                  ::  ie              !< end index for subdomain input along x direction
4817    INTEGER(iwp)                  ::  is              !< start index for subdomain input along x direction
4818    INTEGER(iwp)                  ::  j               !< index along y direction
4819    INTEGER(iwp)                  ::  je              !< end index for subdomain input along y direction
4820    INTEGER(iwp)                  ::  js              !< start index for subdomain input along y direction
4821    INTEGER(iwp)                  ::  k               !< index along any 4th dimension
4822    INTEGER(iwp)                  ::  ke              !< end index of 4th dimension
4823    INTEGER(iwp)                  ::  ks              !< start index of 4th dimension
4824    INTEGER(iwp)                  ::  ns              !< start index for subdomain input along n dimension
4825
4826    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE   ::  tmp !< temporary variable to read data from file according
4827                                                      !< to its reverse memory access
4828
4829    REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) ::  var  !< variable where the read data have to be stored:
4830                                                       !< one dimension is reduced in the process
4831
4832
4833#if defined( __netcdf )
4834!
4835!-- Inquire variable id
4836    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4837!
4838!-- Check for collective read-operation and set respective NetCDF flags if required.
4839    IF ( collective_read )  THEN
4840       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4841    ENDIF
4842!
4843!-- Allocate temporary variable according to memory access on file.
4844    ALLOCATE( tmp(is:ie,js:je,ks:ke) )
4845!
4846!-- Get variable
4847    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ is+1, js+1, ks+1, ns+1 /),                 &
4848                            count = (/ ie-is+1, je-js+1, ke-ks+1, 1   /) )
4849
4850    CALL handle_error( 'get_variable_4d_to_3d_real', 536, variable_name )
4851!
4852!-- Resort data. Please note, dimension subscripts of var all start at 1.
4853    DO  i = is, ie
4854       DO  j = js, je
4855          DO  k = ks, ke
4856             var(k-ks+1,j-js+1,i-is+1) = tmp(i,j,k)
4857          ENDDO
4858       ENDDO
4859    ENDDO
4860
4861   DEALLOCATE( tmp )
4862#endif
4863
4864 END SUBROUTINE get_variable_4d_to_3d_real
4865
4866
4867!--------------------------------------------------------------------------------------------------!
4868! Description:
4869! ------------
4870!> Reads a 3D float variables from dynamic driver with the last dimension only having 1 entry
4871!> (time,z). Please note, the passed arguments are start indices and number of elements in each
4872!> dimension, which is in contrast to the other 3d versions where start- and end indices are passed.
4873!> The different handling compared to get_variable_2d_real is due to its different start-index
4874!> treatment.
4875!--------------------------------------------------------------------------------------------------!
4876 SUBROUTINE get_variable_2d_real_dynamic( id, variable_name, var, i1s, i2s, count_1, count_2 )
4877
4878    USE indices
4879    USE pegrid
4880
4881    IMPLICIT NONE
4882
4883    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4884
4885    INTEGER(iwp)                  ::  count_1         !< number of elements to be read along 1st dimension (with respect to file)
4886    INTEGER(iwp)                  ::  count_2         !< number of elements to be read along 2nd dimension (with respect to file)
4887    INTEGER(iwp)                  ::  i1              !< running index along 1st dimension on file
4888    INTEGER(iwp)                  ::  i1s             !< start index for subdomain input along 1st dimension (with respect to file)
4889    INTEGER(iwp)                  ::  i2              !< running index along 2nd dimension on file
4890    INTEGER(iwp)                  ::  i2s             !< start index for subdomain input along 2nd dimension (with respect to file)
4891    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4892    INTEGER(iwp)                  ::  id_var          !< variable id
4893    INTEGER(iwp)                  ::  lb1             !< lower bound of 1st dimension (with respect to file)
4894    INTEGER(iwp)                  ::  lb2             !< lower bound of 2nd dimension (with respect to file)
4895    INTEGER(iwp)                  ::  ub1             !< upper bound of 1st dimension (with respect to file)
4896    INTEGER(iwp)                  ::  ub2             !< upper bound of 2nd dimension (with respect to file)
4897
4898    REAL(wp), DIMENSION(:,:), ALLOCATABLE   ::  tmp   !< temporary variable to read data from file according to its reverse memory
4899                                                      !< access
4900
4901    REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) ::  var !< input variable
4902
4903
4904#if defined( __netcdf )
4905!
4906!-- Inquire variable id.
4907    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4908!
4909!-- Allocate temporary variable according to memory access on file.
4910!-- Therefore, determine dimension bounds of input array.
4911    lb1 = LBOUND( var,2 )
4912    ub1 = UBOUND( var,2 )
4913    lb2 = LBOUND( var,1 )
4914    ub2 = UBOUND( var,1 )
4915
4916    ALLOCATE( tmp(lb1:ub1,lb2:ub2) )
4917!
4918!-- Get variable
4919    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ i1s, i2s /),                               &
4920                            count = (/ count_1, count_2 /) )
4921
4922    CALL handle_error( 'get_variable_2d_real_dynamic', 537, variable_name )
4923!
4924!-- Resort data. Please note, dimension subscripts of var all start at 1.
4925    DO i2 = lb2, ub2
4926       DO  i1 = lb1, ub1
4927          var(i2,i1,1) = tmp(i1,i2)
4928       ENDDO
4929    ENDDO
4930
4931    DEALLOCATE( tmp )
4932#endif
4933
4934 END SUBROUTINE get_variable_2d_real_dynamic
4935
4936
4937!--------------------------------------------------------------------------------------------------!
4938! Description:
4939! ------------
4940!> Reads a 3D float variables from dynamic driver, such as time-dependent xy-, xz- or yz-boundary
4941!> data as well as 3D initialization data. Please note, the passed arguments are start indices and
4942!> number of elements in each dimension, which is in contrast to the other 3d versions where start-
4943!> and end indices are passed. The different handling of 3D dynamic variables is due to its
4944!> asymmetry for the u- and v component.
4945!--------------------------------------------------------------------------------------------------!
4946 SUBROUTINE get_variable_3d_real_dynamic( id, variable_name, var, i1s, i2s, i3s,                   &
4947                                          count_1, count_2, count_3, par_access )
4948
4949    USE indices
4950    USE pegrid
4951
4952    IMPLICIT NONE
4953
4954    CHARACTER(LEN=*)              ::  variable_name   !< variable name
4955
4956    INTEGER(iwp)                  ::  count_1         !< number of elements to be read along 1st dimension (with respect to file)
4957    INTEGER(iwp)                  ::  count_2         !< number of elements to be read along 2nd dimension (with respect to file)
4958    INTEGER(iwp)                  ::  count_3         !< number of elements to be read along 3rd dimension (with respect to file)
4959    INTEGER(iwp)                  ::  i1              !< running index along 1st dimension on file
4960    INTEGER(iwp)                  ::  i1s             !< start index for subdomain input along 1st dimension (with respect to file)
4961    INTEGER(iwp)                  ::  i2              !< running index along 2nd dimension on file
4962    INTEGER(iwp)                  ::  i2s             !< start index for subdomain input along 2nd dimension (with respect to file)
4963    INTEGER(iwp)                  ::  i3              !< running index along 3rd dimension on file
4964    INTEGER(iwp)                  ::  i3s             !< start index of 3rd dimension, in dynamic file this is either time
4965                                                      !<(2D boundary) or z (3D)
4966    INTEGER(iwp), INTENT(IN)      ::  id              !< file id
4967    INTEGER(iwp)                  ::  id_var          !< variable id
4968    INTEGER(iwp)                  ::  lb1             !< lower bound of 1st dimension (with respect to file)
4969    INTEGER(iwp)                  ::  lb2             !< lower bound of 2nd dimension (with respect to file)
4970    INTEGER(iwp)                  ::  lb3             !< lower bound of 3rd dimension (with respect to file)
4971    INTEGER(iwp)                  ::  ub1             !< upper bound of 1st dimension (with respect to file)
4972    INTEGER(iwp)                  ::  ub2             !< upper bound of 2nd dimension (with respect to file)
4973    INTEGER(iwp)                  ::  ub3             !< upper bound of 3rd dimension (with respect to file)
4974
4975    LOGICAL                       ::  par_access      !< additional flag indicating whether parallel read operations should be
4976                                                      !< performed or not
4977
4978    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE   ::  tmp !< temporary variable to read data from file according
4979                                                      !< to its reverse memory access
4980
4981    REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) ::  var !< input variable
4982
4983
4984#if defined( __netcdf )
4985!
4986!-- Inquire variable id.
4987    nc_stat = NF90_INQ_VARID( id, TRIM( variable_name ), id_var )
4988!
4989!-- Check for collective read-operation and set respective NetCDF flags if required.
4990!-- Please note, in contrast to the other input routines where each PEs reads its subdomain data,
4991!-- dynamic input data not by all PEs, only by those which encompass lateral model boundaries.
4992!-- Hence, collective read operations are only enabled for top-boundary data.
4993    IF ( collective_read  .AND.  par_access )  THEN
4994#if defined( __netcdf4_parallel )
4995       nc_stat = NF90_VAR_PAR_ACCESS (id, id_var, NF90_COLLECTIVE)
4996#endif
4997    ENDIF
4998!
4999!-- Allocate temporary variable according to memory access on file.
5000!-- Therefore, determine dimension bounds of input array.
5001    lb1 = LBOUND( var,3 )
5002    ub1 = UBOUND( var,3 )
5003    lb2 = LBOUND( var,2 )
5004    ub2 = UBOUND( var,2 )
5005    lb3 = LBOUND( var,1 )
5006    ub3 = UBOUND( var,1 )
5007    ALLOCATE( tmp(lb1:ub1,lb2:ub2,lb3:ub3) )
5008!
5009!-- Get variable
5010    nc_stat = NF90_GET_VAR( id, id_var, tmp, start = (/ i1s, i2s, i3s /),                          &
5011                            count = (/ count_1, count_2, count_3 /) )
5012
5013    CALL handle_error( 'get_variable_3d_real_dynamic', 537, variable_name )
5014!
5015!-- Resort data. Please note, dimension subscripts of var all start at 1.
5016    DO  i3 = lb3, ub3
5017       DO i2 = lb2, ub2
5018          DO  i1 = lb1, ub1
5019             var(i3,i2,i1) = tmp(i1,i2,i3)
5020          ENDDO
5021       ENDDO
5022    ENDDO
5023
5024    DEALLOCATE( tmp )
5025#endif
5026
5027 END SUBROUTINE get_variable_3d_real_dynamic
5028
5029
5030!--------------------------------------------------------------------------------------------------!
5031! Description:
5032! ------------
5033!> Reads a 5D float variable from file and store it to a 4-d variable.
5034!--------------------------------------------------------------------------------------------------!
5035 SUBROUTINE