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

Last change on this file since 4389 was 4389, checked in by raasch, 20 months ago

Error messages refined for reading ASCII topo file, also reading of topo file revised so that statement labels and goto statements are not required any more

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