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

Last change on this file since 2696 was 2696, checked in by kanani, 4 years ago

Merge of branch palm4u into trunk

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