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

Last change on this file since 4392 was 4392, checked in by pavelkrc, 20 months ago

Merge branch resler (updated documentation, optional flux tracing, bugfix)

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