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

Last change on this file since 4434 was 4434, checked in by oliver.maas, 4 years ago

added optional netcdf data input for wtm array input parameters

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