Ignore:
Timestamp:
Oct 22, 2018 5:32:49 PM (6 years ago)
Author:
eckhard
Message:

inifor: Added computation of geostrophic winds from COSMO input

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/UTIL/inifor/src/grid.f90

    r3183 r3395  
    2626! -----------------
    2727! $Id$
     28! Added computation of geostrophic winds form COSMO pressure fields
     29! Introduced averaging grids and internal 'output' variables for computation of
     30!     geostrophic and large-scale forcings
     31! Merged qv, uv and vg output variables and 'temperature' IO group into new
     32!     'thermodynamics' IP group
     33! Cleaned up terminal output, show some messages only with --debug option
     34! Improved variable naming and readability
     35!
     36!
     37! 3183 2018-07-27 14:25:55Z suehring
    2838! Introduced new PALM grid stretching
    2939! Updated variable names and metadata for PIDS v1.9 compatibility
     
    6171        ONLY:  DATE, EARTH_RADIUS, TO_RADIANS, TO_DEGREES, PI, dp, hp, sp,     &
    6272               SNAME, LNAME, PATH, FORCING_STEP, WATER_ID, FILL_ITERATIONS,    &
    63                BETA, P_SL, T_SL, BETA, RD, G, P_REF, RD_PALM, CP_PALM, RHO_L
     73               BETA, P_SL, T_SL, BETA, RD, RV, G, P_REF, RD_PALM, CP_PALM,     &
     74               RHO_L, OMEGA, HECTO
    6475    USE io,                                                                    &
    6576        ONLY:  get_netcdf_variable, get_netcdf_attribute,                      &
     
    6879        ONLY:  NF90_MAX_NAME, NF90_MAX_VAR_DIMS
    6980    USE transform,                                                             &
    70         ONLY:  rotate_to_cosmo, find_horizontal_neighbours,                    &
     81        ONLY:  average_2d, rotate_to_cosmo, find_horizontal_neighbours,&
    7182               compute_horizontal_interp_weights,                              &
    72                find_vertical_neighbours_and_weights, interpolate_2d,           &
    73                gamma_from_hemisphere, phic_to_phin, lamc_to_lamn, average_2d,  &
    74                project, centre_velocities, phi2phirot, rla2rlarot, uv2uvrot
     83               find_vertical_neighbours_and_weights_interp,                    &
     84               find_vertical_neighbours_and_weights_average, interpolate_2d,   &
     85               gamma_from_hemisphere, phic_to_phin, lamc_to_lamn,              &
     86               project, centre_velocities, phi2phirot, rla2rlarot, uv2uvrot,   &
     87               phirot2phi, rlarot2rla
    7588    USE types
    7689    USE util
     
    8093    SAVE
    8194   
     95    REAL(dp) ::  averaging_angle   = 0.0_dp       !< latitudal and longitudal width of averaging regions [rad]
     96    REAL(dp) ::  averaging_width_ns = 0.0_dp       !< longitudal width of averaging regions [m]
     97    REAL(dp) ::  averaging_width_ew = 0.0_dp       !< latitudal width of averaging regions [m]
    8298    REAL(dp) ::  phi_equat         = 0.0_dp       !< latitude of rotated equator of COSMO-DE grid [rad]
    8399    REAL(dp) ::  phi_n             = 0.0_dp       !< latitude of rotated pole of COSMO-DE grid [rad]
     
    87103    REAL(dp) ::  phi_cn            = 0.0_dp       !< latitude of the rotated pole relative to the COSMO-DE grid [rad]
    88104    REAL(dp) ::  lambda_cn         = 0.0_dp       !< longitude of the rotated pole relative to the COSMO-DE grid [rad]
     105    REAL(dp) ::  lam_centre        = 0.0_dp       !< longitude of the PLAM domain centre in the source (COSMO rotated-pole) system [rad]
     106    REAL(dp) ::  phi_centre        = 0.0_dp       !< latitude of the PLAM domain centre in the source (COSMO rotated-pole) system [rad]
     107    REAL(dp) ::  lam_east          = 0.0_dp       !< longitude of the east central-averaging-domain boundary in the source (COSMO rotated-pole) system [rad]
     108    REAL(dp) ::  lam_west          = 0.0_dp       !< longitude of the west central-averaging-domain boundary in the source (COSMO rotated-pole) system [rad]
     109    REAL(dp) ::  phi_north         = 0.0_dp       !< latitude of the north central-averaging-domain boundary in the source (COSMO rotated-pole) system [rad]
     110    REAL(dp) ::  phi_south         = 0.0_dp       !< latitude of the south central-averaging-domain boundary in the source (COSMO rotated-pole) system [rad]
    89111    REAL(dp) ::  gam               = 0.0_dp       !< angle for working around phirot2phi/rlarot2rla bug
    90112    REAL(dp) ::  dx                = 0.0_dp       !< PALM-4U grid spacing in x direction [m]
     
    100122    REAL(dp) ::  dyi               = 0.0_dp       !< inverse PALM-4U grid spacing in y direction [m^-1]
    101123    REAL(dp) ::  dzi               = 0.0_dp       !< inverse PALM-4U grid spacing in z direction [m^-1]
     124    REAL(dp) ::  f3                = 0.0_dp       !< Coriolis parameter
    102125    REAL(dp) ::  lx                = 0.0_dp       !< PALM-4U domain size in x direction [m]
    103126    REAL(dp) ::  ly                = 0.0_dp       !< PALM-4U domain size in y direction [m]
    104     REAL(dp) ::  ug                = 0.0_dp       !< geostrophic wind in x direction [m/s]
    105     REAL(dp) ::  vg                = 0.0_dp       !< geostrophic wind in y direction [m/s]
    106127    REAL(dp) ::  p0                = 0.0_dp       !< PALM-4U surface pressure, at z0 [Pa]
    107128    REAL(dp) ::  x0                = 0.0_dp       !< x coordinate of PALM-4U Earth tangent [m]
     
    140161    INTEGER ::  dz_stretch_level_end_index(9)               !< vertical grid level index until which the vertical grid spacing is stretched
    141162    INTEGER ::  dz_stretch_level_start_index(9)             !< vertical grid level index above which the vertical grid spacing is stretched
    142     INTEGER ::  i     !< indexing variable
    143     INTEGER ::  average_imin, average_imax, average_jmin, average_jmax !< index ranges for profile averaging
    144     INTEGER ::  k     !< indexing variable
    145163    INTEGER ::  nt    !< number of output time steps
    146164    INTEGER ::  nx    !< number of PALM-4U grid points in x direction
     
    164182    LOGICAL ::  profile_grids_required
    165183
    166     INTEGER ::  n_invar = 0 !< number of variables in the input variable table
    167     INTEGER ::  n_outvar = 0 !< number of variables in the output variable table
    168184    TYPE(nc_var), ALLOCATABLE, TARGET ::  input_var_table(:)  !< table of input variables
    169185    TYPE(nc_var), ALLOCATABLE, TARGET ::  output_var_table(:) !< table of input variables
     
    219235    TYPE(grid_definition), TARGET ::  w_south_intermediate            !<
    220236    TYPE(grid_definition), TARGET ::  w_top_intermediate              !<
    221     TYPE(grid_definition), TARGET ::  scalar_profile_grid             !<
    222     TYPE(grid_definition), TARGET ::  scalar_profile_intermediate     !<
    223     TYPE(grid_definition), TARGET ::  w_profile_grid                  !<
    224     TYPE(grid_definition), TARGET ::  w_profile_intermediate          !<
     237
     238    TYPE(grid_definition), TARGET ::  north_averaged_scalar_profile
     239    TYPE(grid_definition), TARGET ::  south_averaged_scalar_profile
     240    TYPE(grid_definition), TARGET ::  west_averaged_scalar_profile
     241    TYPE(grid_definition), TARGET ::  east_averaged_scalar_profile
     242    TYPE(grid_definition), TARGET ::  averaged_scalar_profile
     243    TYPE(grid_definition), TARGET ::  averaged_w_profile
    225244
    226245    TYPE(io_group), ALLOCATABLE, TARGET ::  io_group_list(:)  !< List of I/O groups, which group together output variables that share the same input variable
     
    237256    CHARACTER(LEN=PATH), ALLOCATABLE, DIMENSION(:) ::  radiation_files
    238257    CHARACTER(LEN=SNAME) ::  input_prefix         !< Prefix of input files, e.g. 'laf' for COSMO-DE analyses
     258    CHARACTER(LEN=SNAME) ::  flow_prefix          !< Prefix of flow input files, e.g. 'laf' for COSMO-DE analyses
     259    CHARACTER(LEN=SNAME) ::  soil_prefix          !< Prefix of soil input files, e.g. 'laf' for COSMO-DE analyses
     260    CHARACTER(LEN=SNAME) ::  radiation_prefix     !< Prefix of radiation input files, e.g. 'laf' for COSMO-DE analyses
     261    CHARACTER(LEN=SNAME) ::  soilmoisture_prefix  !< Prefix of input files for soil moisture spin-up, e.g. 'laf' for COSMO-DE analyses
    239262    CHARACTER(LEN=SNAME) ::  flow_suffix          !< Suffix of flow input files, e.g. 'flow'
    240263    CHARACTER(LEN=SNAME) ::  soil_suffix          !< Suffix of soil input files, e.g. 'soil'
     
    249272   
    250273    SUBROUTINE setup_parameters()
     274       INTEGER ::  k
    251275
    252276!
     
    282306
    283307       ! Default atmospheric parameters
    284        ug = 0.0_dp
    285        vg = 0.0_dp
     308       cfg % ug = 0.0_dp
     309       cfg % vg = 0.0_dp
    286310       cfg % p0 = P_SL
    287311
     
    296320       input_prefix = 'laf'  ! 'laf' for COSMO-DE analyses
    297321       cfg % flow_prefix = input_prefix
     322       cfg % input_prefix = input_prefix
    298323       cfg % soil_prefix = input_prefix
    299324       cfg % radiation_prefix = input_prefix
     
    304329       radiation_suffix = '-rad'
    305330       soilmoisture_suffix = '-soilmoisture'
     331
     332       cfg % debug = .FALSE.
     333       cfg % averaging_angle = 2.0_dp
    306334!
    307335!------------------------------------------------------------------------------
     
    318346       cfg % ic_mode = 'volume'
    319347       cfg % bc_mode = 'real'
    320 
    321        ! Update default file names and domain centre
     348       cfg % averaging_mode = 'level'
     349
     350       ! Overwrite defaults with user configuration
    322351       CALL parse_command_line_arguments( cfg )
     352
     353       flow_prefix = TRIM(cfg % input_prefix)
     354       radiation_prefix = TRIM(cfg % input_prefix)
     355       soil_prefix = TRIM(cfg % input_prefix)
     356       soilmoisture_prefix = TRIM(cfg % input_prefix)
     357       IF (cfg % flow_prefix_is_set)  flow_prefix = TRIM(cfg % flow_prefix)
     358       IF (cfg % radiation_prefix_is_set)  radiation_prefix = TRIM(cfg % radiation_prefix)
     359       IF (cfg % soil_prefix_is_set)  soil_prefix = TRIM(cfg % soil_prefix)
     360       IF (cfg % soilmoisture_prefix_is_set)  soilmoisture_prefix = TRIM(cfg % soilmoisture_prefix)
    323361
    324362       output_file % name = cfg % output_file
     
    336374       END IF
    337375
     376
     377       ! Set default file paths, if not specified by user.
    338378       CALL normalize_path(cfg % input_path)
    339379       IF (TRIM(cfg % hhl_file) == '')  cfg % hhl_file = TRIM(cfg % input_path) // 'hhl.nc'
     
    361401
    362402       ! Generate input file lists
    363        CALL get_input_file_list(cfg % start_date, start_hour_flow, end_hour, step_hour,  &
    364           cfg % input_path, cfg % flow_prefix, flow_suffix, flow_files)
    365        CALL get_input_file_list(cfg % start_date, start_hour_soil, end_hour, step_hour,  &
    366           cfg % input_path, cfg % soil_prefix, soil_suffix, soil_files)
    367        CALL get_input_file_list(cfg % start_date, start_hour_radiation, end_hour, step_hour, &
    368           cfg % input_path, cfg % radiation_prefix, radiation_suffix, radiation_files)
    369        CALL get_input_file_list(cfg % start_date, start_hour_soilmoisture, end_hour_soilmoisture, step_hour, &
    370           cfg % input_path, cfg % soilmoisture_prefix, soilmoisture_suffix, soil_moisture_files)
     403       CALL get_input_file_list(                                               &
     404          cfg % start_date, start_hour_flow, end_hour, step_hour,              &
     405          cfg % input_path, flow_prefix, flow_suffix, flow_files)
     406       CALL get_input_file_list(                                               &
     407          cfg % start_date, start_hour_soil, end_hour, step_hour,              &
     408          cfg % input_path, soil_prefix, soil_suffix, soil_files)
     409       CALL get_input_file_list(                                               &
     410          cfg % start_date, start_hour_radiation, end_hour, step_hour,         &
     411          cfg % input_path, radiation_prefix, radiation_suffix, radiation_files)
     412       CALL get_input_file_list(                                               &
     413          cfg % start_date, start_hour_soilmoisture, end_hour_soilmoisture, step_hour, &
     414          cfg % input_path, soilmoisture_prefix, soilmoisture_suffix, soil_moisture_files)
    371415
    372416!
     
    499543
    500544 CALL run_control('time', 'comp')
     545
     546!------------------------------------------------------------------------------
     547! Section 4.3: INIFOR averaging domains
     548!------------------------------------------------------------------------------
     549
     550    ! Compute coordiantes of the PALM centre in the source (COSMO) system
     551    phi_centre = phirot2phi(                                                   &
     552       phirot = project(0.5_dp*ly, y0, EARTH_RADIUS) * TO_DEGREES,             &
     553       rlarot = project(0.5_dp*lx, x0, EARTH_RADIUS) * TO_DEGREES,             &
     554       polphi = phi_cn * TO_DEGREES, pollam = lambda_cn * TO_DEGREES,          &
     555       polgam = gam * TO_DEGREES                                               &
     556    ) * TO_RADIANS
     557
     558    lam_centre = rlarot2rla(                                                   &
     559       phirot = project(0.5_dp*ly, y0, EARTH_RADIUS) * TO_DEGREES,             &
     560       rlarot = project(0.5_dp*lx, x0, EARTH_RADIUS) * TO_DEGREES,             &
     561       polphi = phi_cn * TO_DEGREES, pollam = lambda_cn * TO_DEGREES,          &
     562       polgam = gam * TO_DEGREES                                               &
     563    ) * TO_RADIANS
     564
     565    message = "PALM-4U centre:" // NEW_LINE('') // &
     566       "           lon (lambda) = " // &
     567       TRIM(real_to_str_f(lam_centre * TO_DEGREES)) // " deg" // NEW_LINE(' ') // &
     568       "           lat (phi   ) = " // &
     569       TRIM(real_to_str_f(phi_centre * TO_DEGREES)) // " deg (COSMO-DE rotated-pole)"
     570    CALL report( 'setup_parameters', message )
     571
     572
     573    ! Compute boundaries of the central averaging box
     574    averaging_angle = cfg % averaging_angle * TO_RADIANS
     575    lam_east = lam_centre + 0.5_dp * averaging_angle
     576    lam_west = lam_centre - 0.5_dp * averaging_angle
     577    phi_north = phi_centre + 0.5_dp * averaging_angle
     578    phi_south = phi_centre - 0.5_dp * averaging_angle
     579    averaging_width_ew = averaging_angle * COS(phi_centre) * EARTH_RADIUS
     580    averaging_width_ns = averaging_angle * EARTH_RADIUS
     581
     582    ! Coriolis parameter
     583    f3 = 2.0_dp * OMEGA * SIN(                                                 &
     584       TO_RADIANS*phirot2phi( phi_centre * TO_DEGREES, lam_centre * TO_DEGREES,&
     585                              phi_n * TO_DEGREES, lambda_n * TO_DEGREES,       &
     586                              gam * TO_DEGREES )                               &
     587       )
    501588
    502589    END SUBROUTINE setup_parameters
     
    603690               ymin = 0.5_dp * dy, ymax = ly - 0.5_dp * dy,                    &
    604691               x0=x0, y0=y0, z0 = z0,                                          &
    605                nx = nx, ny = ny, nz = nlev - 2)
     692               nx = nx, ny = ny, nz = nlev - 1)
    606693
    607694      IF (boundary_variables_required)  THEN
     
    824911                  ymin =  0.5_dp * dy, ymax = ly - 0.5_dp * dy,                   &
    825912                  x0=x0, y0=y0, z0 = z0,                                          &
    826                   nx = 0, ny = ny, nz = nlev - 2)
     913                  nx = 0, ny = ny, nz = nlev - 1)
    827914
    828915          CALL init_grid_definition('boundary intermediate', grid=w_west_intermediate,         &
     
    830917                  ymin =  0.5_dp * dy, ymax = ly - 0.5_dp * dy,                   &
    831918                  x0=x0, y0=y0, z0 = z0,                                          &
    832                   nx = 0, ny = ny, nz = nlev - 2)
     919                  nx = 0, ny = ny, nz = nlev - 1)
    833920
    834921          CALL init_grid_definition('boundary intermediate', grid=w_north_intermediate,        &
     
    836923                  ymin = ly + 0.5_dp * dy, ymax = ly + 0.5_dp * dy,               &
    837924                  x0=x0, y0=y0, z0 = z0,                                          &
    838                   nx = nx, ny = 0, nz = nlev - 2)
     925                  nx = nx, ny = 0, nz = nlev - 1)
    839926
    840927          CALL init_grid_definition('boundary intermediate', grid=w_south_intermediate,        &
     
    842929                  ymin = -0.5_dp * dy, ymax = -0.5_dp * dy,                       &
    843930                  x0=x0, y0=y0, z0 = z0,                                          &
    844                   nx = nx, ny = 0, nz = nlev - 2)
     931                  nx = nx, ny = 0, nz = nlev - 1)
    845932
    846933          CALL init_grid_definition('boundary intermediate', grid=w_top_intermediate,          &
     
    848935                  ymin =  0.5_dp * dy, ymax = ly - 0.5_dp * dy,                   &
    849936                  x0=x0, y0=y0, z0 = z0,                                          &
    850                   nx = nx, ny = ny, nz = nlev - 2)
     937                  nx = nx, ny = ny, nz = nlev - 1)
    851938       END IF
    852939
     
    856943!------------------------------------------------------------------------------
    857944
    858        profile_grids_required = ( TRIM(cfg % ic_mode) == 'profile' .OR.              &
     945       CALL init_averaging_grid(averaged_scalar_profile, cosmo_grid,           &
     946               x = 0.5_dp * lx, y = 0.5_dp * ly, z = z, z0 = z0,               &
     947               lonmin = lam_west, lonmax = lam_east,                           &
     948               latmin = phi_south, latmax = phi_north,                         &
     949               kind='scalar')
     950
     951       CALL init_averaging_grid(averaged_w_profile, cosmo_grid,                &
     952               x = 0.5_dp * lx, y = 0.5_dp * ly, z = zw, z0 = z0,              &
     953               lonmin = lam_west, lonmax = lam_east,                           &
     954               latmin = phi_south, latmax = phi_north,                         &
     955               kind='w')
     956
     957       CALL init_averaging_grid(south_averaged_scalar_profile, cosmo_grid,     &
     958               x = 0.5_dp * lx, y = 0.5_dp * ly, z = z, z0 = z0,               &
     959               lonmin = lam_west, lonmax = lam_east,                           &
     960               latmin = phi_centre - averaging_angle,                          &
     961               latmax = phi_centre,                                            &
     962               kind='scalar')
     963
     964       CALL init_averaging_grid(north_averaged_scalar_profile, cosmo_grid,     &
     965               x = 0.5_dp * lx, y = 0.5_dp * ly, z = z, z0 = z0,               &
     966               lonmin = lam_west, lonmax = lam_east,                           &
     967               latmin = phi_centre,                                            &
     968               latmax = phi_centre + averaging_angle,                          &
     969               kind='scalar')
     970
     971       CALL init_averaging_grid(west_averaged_scalar_profile, cosmo_grid,      &
     972               x = 0.5_dp * lx, y = 0.5_dp * ly, z = z, z0 = z0,               &
     973               lonmin = lam_centre - averaging_angle,                          &
     974               lonmax = lam_centre,                                            &
     975               latmin = phi_south, latmax = phi_north,                         &
     976               kind='scalar')
     977
     978       CALL init_averaging_grid(east_averaged_scalar_profile, cosmo_grid,      &
     979               x = 0.5_dp * lx, y = 0.5_dp * ly, z = z, z0 = z0,               &
     980               lonmin = lam_centre,                                            &
     981               lonmax = lam_centre + averaging_angle,                          &
     982               latmin = phi_south, latmax = phi_north,                         &
     983               kind='scalar')
     984
     985       profile_grids_required = ( TRIM(cfg % ic_mode) == 'profile' .OR.        &
    859986                                  TRIM(cfg % bc_mode) == 'ideal' )
    860987
    861        IF (profile_grids_required)  THEN
    862           CALL init_grid_definition('boundary', grid=scalar_profile_grid,      &
    863                   xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
    864                   ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
    865                   x0=x0, y0=y0, z0 = z0,                                       &
    866                   nx = 0, ny = 0, nz = nz, z=z)
    867 
    868           CALL init_grid_definition('boundary', grid=w_profile_grid,           &
    869                   xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
    870                   ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
    871                   x0=x0, y0=y0, z0 = z0,                                       &
    872                   nx = 0, ny = 0, nz = nz - 1, z=zw)
    873 
    874           CALL init_grid_definition('boundary', grid=scalar_profile_intermediate,&
    875                   xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
    876                   ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
    877                   x0=x0, y0=y0, z0 = z0,                                       &
    878                   nx = 0, ny = 0, nz = nlev - 2, z=z)
    879 
    880           CALL init_grid_definition('boundary', grid=w_profile_intermediate,   &
    881                   xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
    882                   ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
    883                   x0=x0, y0=y0, z0 = z0,                                       &
    884                   nx = 0, ny = 0, nz = nlev - 2, z=zw)
    885        END IF
     988       !IF (profile_grids_required)  THEN
     989       !   ! Here, I use the 'boundary' kind, so the vertical coordinate uses the
     990       !   ! PALM z coordinate.
     991       !   CALL init_grid_definition('boundary', grid=scalar_profile_grid,      &
     992       !           xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
     993       !           ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
     994       !           x0=x0, y0=y0, z0 = z0,                                       &
     995       !           nx = 0, ny = 0, nz = nz, z=z)
     996
     997       !   CALL init_grid_definition('boundary', grid=w_profile_grid,           &
     998       !           xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
     999       !           ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
     1000       !           x0=x0, y0=y0, z0 = z0,                                       &
     1001       !           nx = 0, ny = 0, nz = nz - 1, z=zw)
     1002
     1003       !   CALL init_grid_definition('boundary', grid=scalar_profile_intermediate,&
     1004       !           xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
     1005       !           ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
     1006       !           x0=x0, y0=y0, z0 = z0,                                       &
     1007       !           nx = 0, ny = 0, nz = nlev - 2, z=z)
     1008
     1009       !   CALL init_grid_definition('boundary', grid=w_profile_intermediate,   &
     1010       !           xmin = 0.5_dp * lx, xmax = 0.5_dp * lx,                      &
     1011       !           ymin = 0.5_dp * ly, ymax = 0.5_dp * ly,                      &
     1012       !           x0=x0, y0=y0, z0 = z0,                                       &
     1013       !           nx = 0, ny = 0, nz = nlev - 2, z=zw)
     1014       !END IF
    8861015
    8871016!                                                                             
     
    9301059
    9311060       IF (TRIM(cfg % ic_mode) == 'profile')  THEN
    932            CALL setup_averaging(cosmo_grid, palm_intermediate,                 &
    933                 average_imin, average_imax, average_jmin, average_jmax)
     1061           !TODO: remove this conditional if not needed.
    9341062       END IF
    9351063       
     
    9381066
    9391067
    940     SUBROUTINE setup_interpolation(cosmo_grid, palm_grid, palm_intermediate, kind, ic_mode)
     1068    SUBROUTINE setup_interpolation(cosmo_grid, grid, intermediate_grid, kind, ic_mode)
    9411069
    9421070       TYPE(grid_definition), INTENT(IN), TARGET    ::  cosmo_grid
    943        TYPE(grid_definition), INTENT(INOUT), TARGET ::  palm_grid, palm_intermediate
     1071       TYPE(grid_definition), INTENT(INOUT), TARGET ::  grid, intermediate_grid
    9441072       CHARACTER, INTENT(IN)                        ::  kind
    9451073       CHARACTER(LEN=*), INTENT(IN), OPTIONAL       ::  ic_mode
    9461074
    947        TYPE(grid_definition), POINTER      ::  grid
    948        REAL(dp), DIMENSION(:), POINTER     ::  lat, lon
    949        REAL(dp), DIMENSION(:,:,:), POINTER ::  h
    950 
    951        LOGICAL :: setup_vertical
    952 
    953        setup_vertical = .TRUE.
    954        IF (PRESENT(ic_mode))  THEN
    955           IF (TRIM(ic_mode) == 'profile')  setup_vertical = .FALSE.
    956        END IF
     1075       REAL(dp), DIMENSION(:), POINTER     ::  cosmo_lat, cosmo_lon
     1076       REAL(dp), DIMENSION(:,:,:), POINTER ::  cosmo_h
     1077
     1078       LOGICAL :: setup_volumetric
    9571079
    9581080!------------------------------------------------------------------------------
     
    9641086       CASE('s') ! scalars
    9651087
    966           lat => cosmo_grid % lat
    967           lon => cosmo_grid % lon
    968           h   => cosmo_grid % hfl
     1088          cosmo_lat => cosmo_grid % lat
     1089          cosmo_lon => cosmo_grid % lon
     1090          cosmo_h   => cosmo_grid % hfl
    9691091
    9701092       CASE('w') ! vertical velocity
    9711093
    972           lat => cosmo_grid % lat
    973           lon => cosmo_grid % lon
    974           h   => cosmo_grid % hhl
     1094          cosmo_lat => cosmo_grid % lat
     1095          cosmo_lon => cosmo_grid % lon
     1096          cosmo_h   => cosmo_grid % hhl
    9751097
    9761098       CASE('u') ! x velocity
    9771099
    978           lat => cosmo_grid % lat
    979           lon => cosmo_grid % lonu
    980           h   => cosmo_grid % hfl
     1100          cosmo_lat => cosmo_grid % lat
     1101          cosmo_lon => cosmo_grid % lonu
     1102          cosmo_h   => cosmo_grid % hfl
    9811103
    9821104       CASE('v') ! y velocity
    9831105
    984           lat => cosmo_grid % latv
    985           lon => cosmo_grid % lon
    986           h   => cosmo_grid % hfl
     1106          cosmo_lat => cosmo_grid % latv
     1107          cosmo_lon => cosmo_grid % lon
     1108          cosmo_h   => cosmo_grid % hfl
    9871109
    9881110       CASE DEFAULT
     
    9931115       END SELECT
    9941116
    995        grid => palm_intermediate
    996 
    997        CALL find_horizontal_neighbours(lat, lon,                               &
    998           grid % clat, grid % clon, grid % ii, grid % jj)
    999 
    1000        CALL compute_horizontal_interp_weights(lat, lon,                        &
    1001           grid % clat, grid % clon, grid % ii, grid % jj, grid % w_horiz)
     1117       CALL find_horizontal_neighbours(cosmo_lat, cosmo_lon,                   &
     1118          intermediate_grid % clat, intermediate_grid % clon,                  &
     1119          intermediate_grid % ii, intermediate_grid % jj)
     1120
     1121       CALL compute_horizontal_interp_weights(cosmo_lat, cosmo_lon,            &
     1122          intermediate_grid % clat, intermediate_grid % clon,                  &
     1123          intermediate_grid % ii, intermediate_grid % jj,                      &
     1124          intermediate_grid % w_horiz)
    10021125
    10031126!------------------------------------------------------------------------------
     
    10051128!------------------------------------------------------------------------------
    10061129
    1007        IF (setup_vertical)  THEN
    1008           ALLOCATE( grid % h(0:grid % nx, 0:grid % ny, 0:grid % nz) )
    1009           grid % h(:,:,:) = - EARTH_RADIUS
     1130       ! If profile initialization is chosen, we--somewhat counterintuitively--
     1131       ! don't need to compute vertical interpolation weights. At least, we
     1132       ! don't need them on the intermediate grid, which fills the entire PALM
     1133       ! domain volume. Instead we need vertical weights for the intermediate
     1134       ! profile grids, which get computed in setup_averaging().
     1135       setup_volumetric = .TRUE.
     1136       IF (PRESENT(ic_mode))  THEN
     1137          IF (TRIM(ic_mode) == 'profile')  setup_volumetric = .FALSE.
     1138       END IF
     1139
     1140       IF (setup_volumetric)  THEN
     1141          ALLOCATE( intermediate_grid % h(0:intermediate_grid % nx,            &
     1142                                          0:intermediate_grid % ny,            &
     1143                                          0:intermediate_grid % nz) )
     1144          intermediate_grid % h(:,:,:) = - EARTH_RADIUS
    10101145
    10111146          ! For w points, use hhl, for scalars use hfl
    10121147          ! compute the full heights for the intermediate grids
    1013           CALL interpolate_2d(cosmo_grid % hfl, grid % h, grid)
    1014           CALL find_vertical_neighbours_and_weights(palm_grid, grid)
     1148          CALL interpolate_2d(cosmo_h, intermediate_grid % h, intermediate_grid)
     1149          CALL find_vertical_neighbours_and_weights_interp(grid, intermediate_grid)
    10151150       END IF
    10161151       
    10171152    END SUBROUTINE setup_interpolation
    1018 
    1019 
    1020     SUBROUTINE setup_averaging(cosmo_grid, palm_intermediate, imin, imax, jmin, jmax)
    1021 
    1022        TYPE(grid_definition), INTENT(IN) ::  cosmo_grid, palm_intermediate
    1023        INTEGER, INTENT(INOUT)            ::  imin, imax, jmin, jmax
    1024 
    1025        TYPE(grid_definition), POINTER    ::  grid
    1026        REAL(dp)                          ::  lonmin_pos,lonmax_pos, latmin_pos, latmax_pos
    1027        REAL(dp)                          ::  cosmo_dxi, cosmo_dyi
    1028 
    1029        cosmo_dxi = 1.0_dp / (cosmo_grid % lon(1) - cosmo_grid % lon(0))
    1030        cosmo_dyi = 1.0_dp / (cosmo_grid % lat(1) - cosmo_grid % lat(0))
    1031 
    1032        ! find horizontal index ranges for profile averaging
    1033        lonmin_pos = (MINVAL(palm_intermediate % clon(:,:)) - cosmo_grid % lon(0)) * cosmo_dxi
    1034        lonmax_pos = (MAXVAL(palm_intermediate % clon(:,:)) - cosmo_grid % lon(0)) * cosmo_dxi
    1035        latmin_pos = (MINVAL(palm_intermediate % clat(:,:)) - cosmo_grid % lat(0)) * cosmo_dyi
    1036        latmax_pos = (MAXVAL(palm_intermediate % clat(:,:)) - cosmo_grid % lat(0)) * cosmo_dyi
    1037 
    1038        imin = FLOOR(lonmin_pos)
    1039        imax = CEILING(lonmax_pos)
    1040        jmin = FLOOR(latmin_pos)
    1041        jmax = CEILING(latmax_pos)
    1042        
    1043        ! average heights for intermediate scalar and w profile grids
    1044        grid => scalar_profile_intermediate
    1045        ALLOCATE( grid % h(0:grid % nx, 0:grid % ny, 0:grid % nz) )
    1046        grid % h(:,:,:) = - EARTH_RADIUS
    1047        CALL average_2d(cosmo_grid % hfl, grid % h(0,0,:), imin, imax, jmin, jmax)
    1048        CALL find_vertical_neighbours_and_weights(scalar_profile_grid, grid)
    1049 
    1050        grid => w_profile_intermediate
    1051        ALLOCATE( grid % h(0:grid % nx, 0:grid % ny, 0:grid % nz) )
    1052        grid % h(:,:,:) = - EARTH_RADIUS
    1053        CALL average_2d(cosmo_grid % hhl, grid % h(0,0,:), imin, imax, jmin, jmax)
    1054        CALL find_vertical_neighbours_and_weights(w_profile_grid, grid)
    1055 
    1056     END SUBROUTINE setup_averaging
    1057 
    10581153
    10591154!------------------------------------------------------------------------------!
     
    12911386
    12921387    END SUBROUTINE init_grid_definition
     1388
     1389
     1390!------------------------------------------------------------------------------!
     1391! Description:
     1392! ------------
     1393!> Initializes averagin_grid-type variables.
     1394!>
     1395!> Input parameters:
     1396!> -----------------
     1397!>
     1398!> cosmo_grid : grid_definition-type varieble initialized as 'cosmo-de' grid
     1399!>    providing longitudes and latitudes of the parent grid
     1400!>
     1401!> x, y : location of the profile to be averaged in the PALM system [m]
     1402!>
     1403!> z : vertical grid coordinates of the profile in the PALM system [m]
     1404!>
     1405!> <lat/lon><min/max>: boundaries of the averaging region in the parent system,
     1406!>    i.e. the COSMO-DE rotated-pole system. [rad]
     1407!>
     1408!> kind : kind of quantity to be averaged using this averaging grid.
     1409!>    Destinguishes COSMO-DE scalar and w-velocity levels. Note that finding the
     1410!>    parent/COSMO columns for the region in get_cosmo_averaging_region() is
     1411!>    independent of 'kind' b/ecause INIFOR uses column-centred u and v velocity
     1412!>    components, which are computed in the preprocessing step.
     1413!>
     1414!> Output parameters:
     1415!> ------------------
     1416!> avg_grid : averaging_grid variable to be initialized.
     1417!------------------------------------------------------------------------------!
     1418    SUBROUTINE init_averaging_grid(avg_grid, cosmo_grid, x, y, z, z0,          &
     1419       lonmin, lonmax, latmin, latmax, kind)
     1420
     1421       TYPE(grid_definition), INTENT(INOUT) ::  avg_grid
     1422       TYPE(grid_definition), INTENT(IN)    ::  cosmo_grid
     1423       REAL(dp), INTENT(IN)                 ::  x, y, z0
     1424       REAL(dp), INTENT(IN), TARGET         ::  z(:)
     1425       REAL(dp), INTENT(IN)                 ::  lonmin, lonmax, latmin, latmax
     1426
     1427       CHARACTER(LEN=*), INTENT(IN)         ::  kind
     1428
     1429       LOGICAL                              :: level_based_averaging
     1430
     1431       INTEGER :: i, j
     1432
     1433       ALLOCATE( avg_grid % x(1) )
     1434       ALLOCATE( avg_grid % y(1) )
     1435       avg_grid % x(1) = x
     1436       avg_grid % y(1) = y
     1437       avg_grid % z => z
     1438       avg_grid % z0 = z0
     1439
     1440       avg_grid % nz = SIZE(z, 1)
     1441
     1442       ALLOCATE( avg_grid % lon(2) )
     1443       ALLOCATE( avg_grid % lat(2) )
     1444       avg_grid % lon(1:2) = (/lonmin, lonmax/)
     1445       avg_grid % lat(1:2) = (/latmin, latmax/)
     1446
     1447       avg_grid % kind = TRIM(kind)
     1448
     1449       ! Find and store COSMO columns that fall into the coordinate range
     1450       ! given by avg_grid % clon, % clat
     1451       CALL get_cosmo_averaging_region(avg_grid, cosmo_grid)
     1452
     1453       ALLOCATE (avg_grid % kkk(avg_grid % n_columns, avg_grid % nz, 2) )
     1454       ALLOCATE (avg_grid % w(avg_grid % n_columns, avg_grid % nz, 2) )
     1455       ! Compute average COSMO levels in the averaging region
     1456       SELECT CASE(avg_grid % kind)
     1457
     1458       CASE('scalar', 'u', 'v')
     1459          avg_grid % cosmo_h => cosmo_grid % hfl
     1460
     1461       CASE('w')
     1462          avg_grid % cosmo_h => cosmo_grid % hhl
     1463
     1464       CASE DEFAULT
     1465          message = "Averaging grid kind '" // TRIM(avg_grid % kind) // &
     1466                    "' is not supported. Use 'scalar', 'u', or 'v'."
     1467          CALL abort('get_cosmo_averaging_region', message)
     1468
     1469       END SELECT
     1470
     1471       ! For level-besed averaging, compute average heights
     1472       level_based_averaging = .TRUE.
     1473       IF (level_based_averaging)  THEN
     1474          ALLOCATE(avg_grid % h(1,1,SIZE(avg_grid % cosmo_h, 3)) )
     1475 
     1476          CALL average_2d(avg_grid % cosmo_h, avg_grid % h(1,1,:),                            &
     1477                          avg_grid % iii, avg_grid % jjj)
     1478       END IF
     1479
     1480       ! Copy averaged grid to all COSMO columnts, leads to computing the same
     1481       ! vertical interpolation weights for all columns and to COSMO grid level
     1482       ! based averaging onto the averaged COSMO heights.
     1483       IF ( TRIM(cfg % averaging_mode) == 'level' )  THEN
     1484          FORALL(                                                              &
     1485             i=1:SIZE(avg_grid % cosmo_h, 1),                                  &
     1486             j=1:SIZE(avg_grid % cosmo_h, 2)                                   &
     1487          )   avg_grid % cosmo_h(i,j,:) = avg_grid % h(1,1,:)
     1488       END IF
     1489
     1490       ! Compute vertical weights and neighbours
     1491       CALL find_vertical_neighbours_and_weights_average(avg_grid)
     1492
     1493    END SUBROUTINE init_averaging_grid
     1494
     1495
     1496    SUBROUTINE get_cosmo_averaging_region(avg_grid, cosmo_grid)
     1497       TYPE(grid_definition), INTENT(INOUT)         ::  avg_grid
     1498       TYPE(grid_definition), TARGET, INTENT(IN)    ::  cosmo_grid
     1499
     1500       REAL(dp), DIMENSION(:), POINTER              ::  cosmo_lon, cosmo_lat
     1501       REAL(dp)                                     ::  dlon, dlat
     1502
     1503       INTEGER ::  i, j, imin, imax, jmin, jmax, l, nx, ny
     1504
     1505
     1506       SELECT CASE( TRIM(avg_grid % kind) )
     1507
     1508       CASE('scalar', 'w')
     1509          cosmo_lon => cosmo_grid % lon
     1510          cosmo_lat => cosmo_grid % lat
     1511
     1512       CASE('u')
     1513          cosmo_lon => cosmo_grid % lonu
     1514          cosmo_lat => cosmo_grid % lat
     1515
     1516       CASE('v')
     1517          cosmo_lon => cosmo_grid % lon
     1518          cosmo_lat => cosmo_grid % latv
     1519
     1520       CASE DEFAULT
     1521          message = "Averaging grid kind '" // TRIM(avg_grid % kind) // &
     1522                    "' is not supported. Use 'scalar', 'u', or 'v'."
     1523          CALL abort('get_cosmo_averaging_region', message)
     1524
     1525       END SELECT
     1526
     1527
     1528       ! FIXME: make dlon, dlat parameters of the grid_defintion type
     1529       dlon = cosmo_lon(1) - cosmo_lon(0)
     1530       dlat = cosmo_lat(1) - cosmo_lat(0)
     1531
     1532       imin = CEILING( (avg_grid % lon(1) - cosmo_lon(0)) / dlon )
     1533       imax = FLOOR(   (avg_grid % lon(2) - cosmo_lon(0)) / dlon )
     1534
     1535       jmin = CEILING( (avg_grid % lat(1) - cosmo_lat(0)) / dlat )
     1536       jmax = FLOOR(   (avg_grid % lat(2) - cosmo_lat(0)) / dlat )
     1537       
     1538       message = "Averaging over "//                                           &
     1539                 TRIM(str(imin)) // " < i < " // TRIM(str(imax)) //            &
     1540                 " and " //                                                    &
     1541                 TRIM(str(jmin)) // " < j < " // TRIM(str(jmax))
     1542       CALL report( 'get_cosmo_averaging_region', message )
     1543
     1544       nx = imax - imin + 1
     1545       ny = jmax - jmin + 1
     1546       avg_grid % n_columns = nx * ny
     1547
     1548       ALLOCATE( avg_grid % iii(avg_grid % n_columns),                                 &
     1549                 avg_grid % jjj(avg_grid % n_columns) )
     1550
     1551       l = 0
     1552       DO j = jmin, jmax
     1553       DO i = imin, imax
     1554          l = l + 1
     1555          avg_grid % iii(l) = i
     1556          avg_grid % jjj(l) = j
     1557       END DO
     1558       END DO
     1559
     1560    END SUBROUTINE get_cosmo_averaging_region
    12931561
    12941562
     
    17181986       )
    17191987
    1720        !potential temperature, surface pressure, including nudging and subsidence
     1988       !potential temperature, surface pressure, specific humidity including
     1989       !nudging and subsidence, and geostrophic winds ug, vg
    17211990       io_group_list(3) = init_io_group(                                       &
    17221991          in_files = flow_files,                                               &
    1723           out_vars = [output_var_table(3:8), output_var_table(42:42),          &
    1724                       output_var_table(49:51)],                                &
    1725           in_var_list = (/input_var_table(3), input_var_table(17)/),           &
    1726           kind = 'temperature'                                                 &
    1727        )
    1728 
    1729        !specific humidity including nudging and subsidence
    1730        io_group_list(4) = init_io_group(                                       &
    1731           in_files = flow_files,                                               &
    1732           out_vars = [output_var_table(9:14), output_var_table(52:54)],        &
    1733           in_var_list = input_var_table(4:4),                                  &
    1734           kind = 'scalar'                                                      &
    1735        )
    1736 
    1737        !u and v velocity and geostrophic winds ug, vg
     1992          out_vars = [output_var_table(56:64),                                 & ! internal averaged density and pressure profiles
     1993                      output_var_table(3:8), output_var_table(9:14),           &
     1994                      output_var_table(42:42), output_var_table(43:44),        &
     1995                      output_var_table(49:51), output_var_table(52:54)],       &
     1996          in_var_list = (/input_var_table(3), input_var_table(17),             & ! T, P, QV
     1997                          input_var_table(4) /),                               &
     1998          kind = 'thermodynamics',                                             &
     1999          n_output_quantities = 4                                              & ! P, Theta, Rho, qv
     2000       )
     2001
     2002       !Moved to therodynamic io_group
     2003       !io_group_list(4) = init_io_group(                                       &
     2004       !   in_files = flow_files,                                               &
     2005       !   out_vars = [output_var_table(9:14), output_var_table(52:54)],        &
     2006       !   in_var_list = input_var_table(4:4),                                  &
     2007       !   kind = 'scalar'                                                      &
     2008       !)
     2009
     2010       !u and v velocity, ug anv vg moved to thermodynamic io group.
    17382011       io_group_list(5) = init_io_group(                                       &
    17392012          in_files = flow_files,                                               &
    1740           out_vars = [output_var_table(15:26), output_var_table(43:46)],       &
     2013          out_vars = [output_var_table(15:26), output_var_table(45:46)],       &
    17412014          !out_vars = output_var_table(15:20),                                  &
    17422015          in_var_list = input_var_table(5:6),                                  &
     
    18452118
    18462119
    1847     FUNCTION init_io_group(in_files, out_vars, in_var_list, kind) RESULT(group)
     2120    FUNCTION init_io_group(in_files, out_vars, in_var_list, kind,              &
     2121                           n_output_quantities) RESULT(group)
    18482122       CHARACTER(LEN=PATH), INTENT(IN) ::  in_files(:)
    18492123       CHARACTER(LEN=*), INTENT(IN)    ::  kind
    18502124       TYPE(nc_var), INTENT(IN)        ::  out_vars(:)
    18512125       TYPE(nc_var), INTENT(IN)        ::  in_var_list(:)
     2126       INTEGER, OPTIONAL               ::  n_output_quantities
    18522127
    18532128       TYPE(io_group)                  ::  group
     
    18552130       group % nt = SIZE(in_files)
    18562131       group % nv = SIZE(out_vars)
     2132       group % n_inputs = SIZE(in_var_list)
    18572133       group % kind = TRIM(kind)
    18582134
    1859        ALLOCATE(group % in_var_list( SIZE(in_var_list) ))
     2135       IF ( PRESENT(n_output_quantities) )  THEN
     2136          group % n_output_quantities = n_output_quantities
     2137       ELSE
     2138          group % n_output_quantities = group % n_inputs
     2139       END IF
     2140
     2141       ALLOCATE(group % in_var_list(group % n_inputs))
    18602142       ALLOCATE(group % in_files(group % nt))
    18612143       ALLOCATE(group % out_vars(group % nv))
     
    18762158    SUBROUTINE fini_grids()
    18772159
    1878        CALL report('fini_grids', 'Deallocating grids')
     2160       CALL report('fini_grids', 'Deallocating grids', cfg % debug)
    18792161       
    18802162       DEALLOCATE(x, y, z, xu, yv, zw, z_column, zw_column)
     
    19042186    SUBROUTINE setup_variable_tables(ic_mode)
    19052187       CHARACTER(LEN=*), INTENT(IN) ::  ic_mode
     2188       INTEGER                      ::  n_invar = 0  !< number of variables in the input variable table
     2189       INTEGER                      ::  n_outvar = 0 !< number of variables in the output variable table
    19062190       TYPE(nc_var), POINTER        ::  var
    19072191
     
    19142198
    19152199       n_invar = 17
    1916        n_outvar = 55
     2200       n_outvar = 64
    19172201       ALLOCATE( input_var_table(n_invar) )
    19182202       ALLOCATE( output_var_table(n_outvar) )
     
    20522336       )
    20532337       IF (TRIM(ic_mode) == 'profile')  THEN
    2054           output_var_table(3) % grid => scalar_profile_grid
    2055           output_var_table(3) % intermediate_grid => scalar_profile_intermediate
     2338          output_var_table(3) % averaging_grid => averaged_scalar_profile
    20562339       END IF
    20572340
     
    21222405          units             = "kg/kg",                                         &
    21232406          kind              = "init scalar",                                   &
    2124           input_id          = 1,                                               &
     2407          input_id          = 3,                                               &
    21252408          output_file       = output_file,                                     &
    21262409          grid              = palm_grid,                                       &
     
    21292412       )
    21302413       IF (TRIM(ic_mode) == 'profile')  THEN
    2131           output_var_table(9) % grid => scalar_profile_grid
    2132           output_var_table(9) % intermediate_grid => scalar_profile_intermediate
     2414          output_var_table(9) % averaging_grid => averaged_scalar_profile
    21332415       END IF
    21342416
     
    21392421          units             = "kg/kg",                                         &
    21402422          kind              = "left scalar",                                   &
    2141           input_id          = 1,                                               &
     2423          input_id          = 3,                                               &
    21422424          output_file       = output_file,                                     &
    21432425          grid              = scalars_west_grid,                               &
     
    21512433          units             = "kg/kg",                                         &
    21522434          kind              = "right scalar",                                  &
    2153           input_id          = 1,                                               &
     2435          input_id          = 3,                                               &
    21542436          output_file       = output_file,                                     &
    21552437          grid              = scalars_east_grid,                               &
     
    21632445          units             = "kg/kg",                                         &
    21642446          kind              = "north scalar",                                  &
    2165           input_id          = 1,                                               &
     2447          input_id          = 3,                                               &
    21662448          output_file       = output_file,                                     &
    21672449          grid              = scalars_north_grid,                              &
     
    21752457          units             = "kg/kg",                                         &
    21762458          kind              = "south scalar",                                  &
    2177           input_id          = 1,                                               &
     2459          input_id          = 3,                                               &
    21782460          output_file       = output_file,                                     &
    21792461          grid              = scalars_south_grid,                              &
     
    21872469          units             = "kg/kg",                                         &
    21882470          kind              = "top scalar",                                    &
    2189           input_id          = 1,                                               &
     2471          input_id          = 3,                                               &
    21902472          output_file       = output_file,                                     &
    21912473          grid              = scalars_top_grid,                                &
     
    22062488       )
    22072489       IF (TRIM(ic_mode) == 'profile')  THEN
    2208           output_var_table(15) % grid => scalar_profile_grid
    2209           output_var_table(15) % intermediate_grid => scalar_profile_intermediate
     2490          output_var_table(15) % averaging_grid => averaged_scalar_profile
    22102491       END IF
    22112492
     
    22832564       )
    22842565       IF (TRIM(ic_mode) == 'profile')  THEN
    2285           output_var_table(21) % grid => scalar_profile_grid
    2286           output_var_table(21) % intermediate_grid => scalar_profile_intermediate
     2566          output_var_table(21) % averaging_grid => averaged_scalar_profile
    22872567       END IF
    22882568
     
    23602640       )
    23612641       IF (TRIM(ic_mode) == 'profile')  THEN
    2362           output_var_table(27) % grid => w_profile_grid
    2363           output_var_table(27) % intermediate_grid => w_profile_intermediate
     2642          output_var_table(27) % averaging_grid => averaged_w_profile
    23642643       END IF
    23652644
     
    25462825          intermediate_grid = palm_intermediate                                &
    25472826       )
     2827       output_var_table(42) % averaging_grid => averaged_scalar_profile
    25482828
    25492829       output_var_table(43) = init_nc_var(                                     &
     
    25522832          long_name         = "geostrophic wind (u component)",                &
    25532833          units             = "m/s",                                           &
    2554           kind              = "constant scalar profile",                       &
     2834          kind              = "geostrophic",                                   &
    25552835          input_id          = 1,                                               &
    25562836          output_file       = output_file,                                     &
    2557           grid              = scalar_profile_grid,                             &
    2558           intermediate_grid = scalar_profile_intermediate                      &
     2837          grid              = averaged_scalar_profile,                         &
     2838          intermediate_grid = averaged_scalar_profile                          &
    25592839       )
    25602840
     
    25642844          long_name         = "geostrophic wind (v component)",                &
    25652845          units             = "m/s",                                           &
    2566           kind              = "constant scalar profile",                       &
     2846          kind              = "geostrophic",                                   &
    25672847          input_id          = 1,                                               &
    25682848          output_file       = output_file,                                     &
    2569           grid              = scalar_profile_grid,                             &
    2570           intermediate_grid = scalar_profile_intermediate                      &
     2849          grid              = averaged_scalar_profile,                         &
     2850          intermediate_grid = averaged_scalar_profile                          &
    25712851       )
    25722852
     
    25762856          long_name         = "wind component in x direction",                 &
    25772857          units             = "m/s",                                           &
    2578           kind              = "large-scale scalar forcing",                    &
     2858          kind              = "geostrophic",                                   &
    25792859          input_id          = 1,                                               &
    25802860          output_file       = output_file,                                     &
    2581           grid              = scalar_profile_grid,                             &
    2582           intermediate_grid = scalar_profile_intermediate                      &
     2861          grid              = averaged_scalar_profile,                         &
     2862          intermediate_grid = averaged_scalar_profile                          &
    25832863       )
    25842864       output_var_table(45) % to_be_processed = ls_forcing_variables_required
     
    25922872          input_id          = 1,                                               &
    25932873          output_file       = output_file,                                     &
    2594           grid              = scalar_profile_grid,                             &
    2595           intermediate_grid = scalar_profile_intermediate                      &
     2874          grid              = averaged_scalar_profile,                         &
     2875          intermediate_grid = averaged_scalar_profile                          &
    25962876       )
    25972877       output_var_table(46) % to_be_processed = ls_forcing_variables_required
     
    26052885          input_id          = 1,                                               &
    26062886          output_file       = output_file,                                     &
    2607           grid              = w_profile_grid,                                  &
    2608           intermediate_grid = w_profile_intermediate                           &
     2887          grid              = averaged_scalar_profile,                         &
     2888          intermediate_grid = averaged_scalar_profile                          &
    26092889       )
    26102890       output_var_table(47) % to_be_processed = ls_forcing_variables_required
     
    26182898          input_id          = 1,                                               &
    26192899          output_file       = output_file,                                     &
    2620           grid              = w_profile_grid,                                  &
    2621           intermediate_grid = w_profile_intermediate                           &
     2900          grid              = averaged_w_profile,                              &
     2901          intermediate_grid = averaged_w_profile                               &
    26222902       )
    26232903       output_var_table(48) % to_be_processed = ls_forcing_variables_required
     
    26322912          input_id          = 1,                                               &
    26332913          output_file       = output_file,                                     &
    2634           grid              = scalar_profile_grid,                             &
    2635           intermediate_grid = scalar_profile_intermediate                      &
     2914          grid              = averaged_scalar_profile,                         &
     2915          intermediate_grid = averaged_scalar_profile                          &
    26362916       )
    26372917       output_var_table(49) % to_be_processed = ls_forcing_variables_required
     
    26452925          input_id          = 1,                                               &
    26462926          output_file       = output_file,                                     &
    2647           grid              = scalar_profile_grid,                             &
    2648           intermediate_grid = scalar_profile_intermediate                      &
     2927          grid              = averaged_scalar_profile,                         &
     2928          intermediate_grid = averaged_scalar_profile                          &
    26492929       )
    26502930       output_var_table(50) % to_be_processed = ls_forcing_variables_required
     
    26582938          input_id          = 1,                                               &
    26592939          output_file       = output_file,                                     &
    2660           grid              = scalar_profile_grid,                             &
    2661           intermediate_grid = scalar_profile_intermediate                      &
     2940          grid              = averaged_scalar_profile,                         &
     2941          intermediate_grid = averaged_scalar_profile                          &
    26622942       )
    26632943       output_var_table(51) % to_be_processed = ls_forcing_variables_required
     
    26692949          units             = "kg/kg/s",                                       &
    26702950          kind              = "large-scale scalar forcing",                    &
    2671           input_id          = 1,                                               &
    2672           output_file       = output_file,                                     &
    2673           grid              = scalar_profile_grid,                             &
    2674           intermediate_grid = scalar_profile_intermediate                      &
     2951          input_id          = 3,                                               &
     2952          output_file       = output_file,                                     &
     2953          grid              = averaged_scalar_profile,                         &
     2954          intermediate_grid = averaged_scalar_profile                          &
    26752955       )
    26762956       output_var_table(52) % to_be_processed = ls_forcing_variables_required
     
    26832963          units             = "kg/kg/s",                                       &
    26842964          kind              = "large-scale scalar forcing",                    &
    2685           input_id          = 1,                                               &
    2686           output_file       = output_file,                                     &
    2687           grid              = scalar_profile_grid,                             &
    2688           intermediate_grid = scalar_profile_intermediate                      &
     2965          input_id          = 3,                                               &
     2966          output_file       = output_file,                                     &
     2967          grid              = averaged_scalar_profile,                         &
     2968          intermediate_grid = averaged_scalar_profile                          &
    26892969       )
    26902970       output_var_table(53) % to_be_processed = ls_forcing_variables_required
     
    26962976          units             = "kg/kg",                                         &
    26972977          kind              = "large-scale scalar forcing",                    &
    2698           input_id          = 1,                                               &
    2699           output_file       = output_file,                                     &
    2700           grid              = scalar_profile_grid,                             &
    2701           intermediate_grid = scalar_profile_intermediate                      &
     2978          input_id          = 3,                                               &
     2979          output_file       = output_file,                                     &
     2980          grid              = averaged_scalar_profile,                         &
     2981          intermediate_grid = averaged_scalar_profile                          &
    27022982       )
    27032983       output_var_table(54) % to_be_processed = ls_forcing_variables_required
     
    27112991          input_id          = 1,                                               &
    27122992          output_file       = output_file,                                     &
    2713           grid              = scalar_profile_grid,                             &
    2714           intermediate_grid = scalar_profile_intermediate                      &
     2993          grid              = averaged_scalar_profile,                         &
     2994          intermediate_grid = averaged_scalar_profile                          &
    27152995       )
    27162996       output_var_table(55) % to_be_processed = ls_forcing_variables_required
    27172997
     2998
     2999       output_var_table(56) = init_nc_var(                                     &
     3000          name              = 'internal_density_centre',                              &
     3001          std_name          = "",                                              &
     3002          long_name         = "",                                              &
     3003          units             = "",                                              &
     3004          kind              = "internal profile",                              &
     3005          input_id          = 4,                                               &
     3006          output_file       = output_file,                                     &
     3007          grid              = averaged_scalar_profile,                         &
     3008          intermediate_grid = averaged_scalar_profile                          &
     3009       )
     3010       output_var_table(56) % averaging_grid => averaged_scalar_profile
     3011
     3012
     3013       output_var_table(57) = init_nc_var(                                     &
     3014          name              = 'internal_density_north',                       &
     3015          std_name          = "",                                              &
     3016          long_name         = "",                                              &
     3017          units             = "",                                              &
     3018          kind              = "internal profile",                              &
     3019          input_id          = 4,                                               &
     3020          output_file       = output_file,                                     &
     3021          grid              = north_averaged_scalar_profile,                   &
     3022          intermediate_grid = north_averaged_scalar_profile                    &
     3023       )
     3024       output_var_table(57) % averaging_grid => north_averaged_scalar_profile
     3025       output_var_table(57) % to_be_processed = .NOT. cfg % ug_is_set
     3026
     3027
     3028       output_var_table(58) = init_nc_var(                                     &
     3029          name              = 'internal_density_south',                       &
     3030          std_name          = "",                                              &
     3031          long_name         = "",                                              &
     3032          units             = "",                                              &
     3033          kind              = "internal profile",                              &
     3034          input_id          = 4,                                               &
     3035          output_file       = output_file,                                     &
     3036          grid              = south_averaged_scalar_profile,                   &
     3037          intermediate_grid = south_averaged_scalar_profile                    &
     3038       )
     3039       output_var_table(58) % averaging_grid => south_averaged_scalar_profile
     3040       output_var_table(58) % to_be_processed = .NOT. cfg % ug_is_set
     3041
     3042
     3043       output_var_table(59) = init_nc_var(                                     &
     3044          name              = 'internal_density_east',                        &
     3045          std_name          = "",                                              &
     3046          long_name         = "",                                              &
     3047          units             = "",                                              &
     3048          kind              = "internal profile",                              &
     3049          input_id          = 4,                                               &
     3050          output_file       = output_file,                                     &
     3051          grid              = east_averaged_scalar_profile,                    &
     3052          intermediate_grid = east_averaged_scalar_profile                     &
     3053       )
     3054       output_var_table(59) % averaging_grid => east_averaged_scalar_profile
     3055       output_var_table(59) % to_be_processed = .NOT. cfg % ug_is_set
     3056
     3057
     3058       output_var_table(60) = init_nc_var(                                     &
     3059          name              = 'internal_density_west',                        &
     3060          std_name          = "",                                              &
     3061          long_name         = "",                                              &
     3062          units             = "",                                              &
     3063          kind              = "internal profile",                              &
     3064          input_id          = 4,                                               &
     3065          output_file       = output_file,                                     &
     3066          grid              = west_averaged_scalar_profile,                    &
     3067          intermediate_grid = west_averaged_scalar_profile                     &
     3068       )
     3069       output_var_table(60) % averaging_grid => west_averaged_scalar_profile
     3070       output_var_table(60) % to_be_processed = .NOT. cfg % ug_is_set
     3071
     3072       output_var_table(61) = init_nc_var(                                     &
     3073          name              = 'internal_pressure_north',                       &
     3074          std_name          = "",                                              &
     3075          long_name         = "",                                              &
     3076          units             = "",                                              &
     3077          kind              = "internal profile",                              &
     3078          input_id          = 2,                                               &
     3079          output_file       = output_file,                                     &
     3080          grid              = north_averaged_scalar_profile,                   &
     3081          intermediate_grid = north_averaged_scalar_profile                    &
     3082       )
     3083       output_var_table(61) % averaging_grid => north_averaged_scalar_profile
     3084       output_var_table(61) % to_be_processed = .NOT. cfg % ug_is_set
     3085
     3086
     3087       output_var_table(62) = init_nc_var(                                     &
     3088          name              = 'internal_pressure_south',                       &
     3089          std_name          = "",                                              &
     3090          long_name         = "",                                              &
     3091          units             = "",                                              &
     3092          kind              = "internal profile",                              &
     3093          input_id          = 2,                                               &
     3094          output_file       = output_file,                                     &
     3095          grid              = south_averaged_scalar_profile,                   &
     3096          intermediate_grid = south_averaged_scalar_profile                    &
     3097       )
     3098       output_var_table(62) % averaging_grid => south_averaged_scalar_profile
     3099       output_var_table(62) % to_be_processed = .NOT. cfg % ug_is_set
     3100
     3101
     3102       output_var_table(63) = init_nc_var(                                     &
     3103          name              = 'internal_pressure_east',                        &
     3104          std_name          = "",                                              &
     3105          long_name         = "",                                              &
     3106          units             = "",                                              &
     3107          kind              = "internal profile",                              &
     3108          input_id          = 2,                                               &
     3109          output_file       = output_file,                                     &
     3110          grid              = east_averaged_scalar_profile,                    &
     3111          intermediate_grid = east_averaged_scalar_profile                     &
     3112       )
     3113       output_var_table(63) % averaging_grid => east_averaged_scalar_profile
     3114       output_var_table(63) % to_be_processed = .NOT. cfg % ug_is_set
     3115
     3116
     3117       output_var_table(64) = init_nc_var(                                     &
     3118          name              = 'internal_pressure_west',                        &
     3119          std_name          = "",                                              &
     3120          long_name         = "",                                              &
     3121          units             = "",                                              &
     3122          kind              = "internal profile",                              &
     3123          input_id          = 2,                                               &
     3124          output_file       = output_file,                                     &
     3125          grid              = west_averaged_scalar_profile,                    &
     3126          intermediate_grid = west_averaged_scalar_profile                     &
     3127       )
     3128       output_var_table(64) % averaging_grid => west_averaged_scalar_profile
     3129       output_var_table(64) % to_be_processed = .NOT. cfg % ug_is_set
    27183130
    27193131       ! Attributes shared among all variables
     
    27323144!------------------------------------------------------------------------------!
    27333145    FUNCTION init_nc_var(name, std_name, long_name, units, kind, input_id,     &
    2734                          grid, intermediate_grid, output_file, is_profile      &
    2735              ) RESULT(var)
     3146                         grid, intermediate_grid, output_file, is_profile)     &
     3147       RESULT(var)
    27363148
    27373149       CHARACTER(LEN=*), INTENT(IN)      ::  name, std_name, long_name, units, kind
     
    27723184          var % dimvarids(1:3)  = output_file % dimvarids_soil
    27733185          var % to_be_processed = init_variables_required
     3186          var % is_internal     = .FALSE.
    27743187          var % task            = "interpolate_2d"
    27753188
     
    27813194          var % dimvarids(1:3)  = output_file % dimvarids_scl
    27823195          var % to_be_processed = init_variables_required
     3196          var % is_internal     = .FALSE.
    27833197          var % task            = "interpolate_3d"
    27843198
     
    27943208          var % dimvarids(3)    = output_file % dimvarids_scl(3)
    27953209          var % to_be_processed = init_variables_required
     3210          var % is_internal     = .FALSE.
    27963211          var % task            = "interpolate_3d"
    27973212
     
    28073222          var % dimvarids(3)    = output_file % dimvarids_scl(3)
    28083223          var % to_be_processed = init_variables_required
     3224          var % is_internal     = .FALSE.
    28093225          var % task            = "interpolate_3d"
    28103226
     
    28203236          var % dimvarids(3)    = output_file % dimvarids_vel(3)
    28213237          var % to_be_processed = init_variables_required
     3238          var % is_internal     = .FALSE.
    28223239          var % task            = "interpolate_3d"
    28233240
     
    28293246          var % dimvarids(1)    = output_file % dimvarids_scl(3) !z
    28303247          var % to_be_processed = init_variables_required
     3248          var % is_internal     = .FALSE.
    28313249          var % task            = "average profile"
    28323250
     
    28383256          var % dimvarids(1)    = output_file % dimvarids_vel(3) !z
    28393257          var % to_be_processed = init_variables_required
     3258          var % is_internal     = .FALSE.
    28403259          var % task            = "average profile"
    28413260
     
    28483267          var % dimvarids(1:2)  = output_file % dimvarids_soil(1:2)
    28493268          var % to_be_processed = boundary_variables_required
     3269          var % is_internal     = .FALSE.
    28503270          var % task            = "interpolate_2d"
    28513271
     
    28603280          var % dimvarids(2)    = output_file % dimvarids_scl(3)
    28613281          var % to_be_processed = boundary_variables_required
     3282          var % is_internal     = .FALSE.
    28623283          var % task            = "interpolate_3d"
    28633284
     
    28723293          var % dimvarids(2)    = output_file % dimvarids_scl(3)
    28733294          var % to_be_processed = boundary_variables_required
     3295          var % is_internal     = .FALSE.
    28743296          var % task            = "interpolate_3d"
    28753297
     
    28843306          var % dimvarids(2)    = output_file % dimvarids_scl(2)
    28853307          var % to_be_processed = boundary_variables_required
     3308          var % is_internal     = .FALSE.
    28863309          var % task            = "interpolate_3d"
    28873310
     
    28963319          var % dimvarids(2)    = output_file % dimvarids_scl(3)
    28973320          var % to_be_processed = boundary_variables_required
     3321          var % is_internal     = .FALSE.
    28983322          var % task            = "interpolate_3d"
    28993323
     
    29083332          var % dimvarids(2)    = output_file % dimvarids_scl(3)
    29093333          var % to_be_processed = boundary_variables_required
     3334          var % is_internal     = .FALSE.
    29103335          var % task            = "interpolate_3d"
    29113336
     
    29203345          var % dimvarids(2)    = output_file % dimvarids_scl(2)
    29213346          var % to_be_processed = boundary_variables_required
     3347          var % is_internal     = .FALSE.
    29223348          var % task            = "interpolate_3d"
    29233349
     
    29323358          var % dimvarids(2)    = output_file % dimvarids_scl(3)
    29333359          var % to_be_processed = boundary_variables_required
     3360          var % is_internal     = .FALSE.
    29343361          var % task            = "interpolate_3d"
    29353362
     
    29443371          var % dimvarids(2)    = output_file % dimvarids_scl(3)
    29453372          var % to_be_processed = boundary_variables_required
     3373          var % is_internal     = .FALSE.
    29463374          var % task            = "interpolate_3d"
    29473375
     
    29563384          var % dimvarids(2)    = output_file % dimvarids_vel(2)
    29573385          var % to_be_processed = boundary_variables_required
     3386          var % is_internal     = .FALSE.
    29583387          var % task            = "interpolate_3d"
    29593388
     
    29683397          var % dimvarids(2)    = output_file % dimvarids_vel(3)
    29693398          var % to_be_processed = boundary_variables_required
     3399          var % is_internal     = .FALSE.
    29703400          var % task            = "interpolate_3d"
    29713401
     
    29803410          var % dimvarids(2)    = output_file % dimvarids_vel(3)
    29813411          var % to_be_processed = boundary_variables_required
     3412          var % is_internal     = .FALSE.
    29823413          var % task            = "interpolate_3d"
    29833414
     
    29883419          var % dimvarids(1)    = output_file % dimvarid_time
    29893420          var % to_be_processed = .TRUE.
    2990           var % task            = "average scalar"
     3421          var % is_internal     = .FALSE.
     3422          var % task            = "average profile"
    29913423
    29923424       CASE( 'constant scalar profile' )
     
    29983430          var % dimvarids(1)    = output_file % dimvarids_scl(3)
    29993431          var % to_be_processed = .TRUE.
     3432          var % is_internal     = .FALSE.
    30003433          var % task            = "set profile"
    30013434
     
    30083441          var % dimvarids(1)    = output_file % dimvarids_scl(3)
    30093442          var % to_be_processed = ls_forcing_variables_required
     3443          var % is_internal     = .FALSE.
    30103444          var % task            = "average large-scale profile"
     3445
     3446       CASE( 'geostrophic' )
     3447          var % lod             = -1
     3448          var % ndim            = 2
     3449          var % dimids(2)       = output_file % dimid_time    !t
     3450          var % dimids(1)       = output_file % dimids_scl(3) !z
     3451          var % dimvarids(2)    = output_file % dimvarid_time
     3452          var % dimvarids(1)    = output_file % dimvarids_scl(3)
     3453          var % to_be_processed = .TRUE.
     3454          var % is_internal     = .FALSE.
     3455          var % task            = "geostrophic winds"
    30113456
    30123457       CASE( 'large-scale w forcing' )
     
    30183463          var % dimvarids(1)    = output_file % dimvarids_vel(3)
    30193464          var % to_be_processed = ls_forcing_variables_required
     3465          var % is_internal     = .FALSE.
    30203466          var % task            = "average large-scale profile"
     3467
     3468       CASE( 'internal profile' )
     3469          var % lod             = -1
     3470          var % ndim            = 2
     3471          var % dimids(2)       = output_file % dimid_time    !t
     3472          var % dimids(1)       = output_file % dimids_scl(3) !z
     3473          var % dimvarids(2)    = output_file % dimvarid_time
     3474          var % dimvarids(1)    = output_file % dimvarids_scl(3)
     3475          var % to_be_processed = .TRUE.
     3476          var % is_internal     = .TRUE.
     3477          var % task            = "internal profile"
    30213478
    30223479       CASE DEFAULT
     
    30313488    SUBROUTINE fini_variables()
    30323489
    3033        CALL report('fini_variables', 'Deallocating variable table')
     3490       CALL report('fini_variables', 'Deallocating variable table', cfg % debug)
    30343491       DEALLOCATE( input_var_table )
    30353492
     
    30393496    SUBROUTINE fini_io_groups()
    30403497
    3041        CALL report('fini_io_groups', 'Deallocating IO groups')
     3498       CALL report('fini_io_groups', 'Deallocating IO groups', cfg % debug)
    30423499       DEALLOCATE( io_group_list )
    30433500
     
    30473504    SUBROUTINE fini_file_lists()
    30483505       
    3049        CALL report('fini_file_lists', 'Deallocating file lists')
     3506       CALL report('fini_file_lists', 'Deallocating file lists', cfg % debug)
    30503507       DEALLOCATE( flow_files, soil_files, radiation_files, soil_moisture_files )
    30513508
     
    31343591             ! rotated velocities
    31353592             DO k = 1, nz
    3136              DO j = 2, ny
    3137              DO i = 2, nx
     3593             DO j = 1, ny
     3594             DO i = 1, nx
    31383595                CALL uv2uvrot( urot = preprocess_buffer(1) % array(i,j,k),     &
    31393596                               vrot = preprocess_buffer(2) % array(i,j,k),     &
     
    31713628          CALL report('preprocess', message)
    31723629       
    3173        CASE( 'temperature' ) ! P and T
     3630       CASE( 'thermodynamics' ) ! T, P, QV
    31743631          nx = SIZE(input_buffer(1) % array, 1)
    31753632          ny = SIZE(input_buffer(1) % array, 2)
     
    31923649
    31933650                input_buffer (2) % array(i,j,:) =                              &
    3194                    input_buffer (2) % array(i,j,:) + basic_state_pressure(:)
     3651                   HECTO * input_buffer (2) % array(i,j,:) +                   &
     3652                   basic_state_pressure(:)
    31953653
    31963654             END DO
     
    32043662
    32053663          END IF
     3664          ! pressure
     3665          input_buffer(2) % is_preprocessed = .TRUE.
     3666
     3667          ! Copy temperature to last input buffer array
     3668          ALLOCATE(                                                            &
     3669              input_buffer( group % n_output_quantities ) % array (nx, ny, nz) &
     3670          )
     3671          input_buffer(group % n_output_quantities) % array(:,:,:) =           &
     3672              input_buffer(1) % array(:,:,:)
    32063673
    32073674          ! Convert absolute to potential temperature
    3208           input_buffer(1) % array(:,:,:) = input_buffer(1) % array(:,:,:) *    &
    3209              EXP( RD_PALM/CP_PALM * LOG( P_REF / input_buffer(2) % array(:,:,:) ))
    3210 
    3211           input_buffer(1:2) % is_preprocessed = .TRUE.
     3675          CALL potential_temperature(                                          &
     3676             t = input_buffer(1) % array(:,:,:),                               &
     3677             p = input_buffer(2) % array(:,:,:),                               &
     3678             p_ref = P_REF,                                                    &
     3679             r = RD_PALM,                                                      &
     3680             cp = CP_PALM                                                      &
     3681          )
     3682
     3683          ! potential temperature
     3684          input_buffer(1) % is_preprocessed = .TRUE.
     3685
     3686          ! Convert temperature copy to density
     3687          CALL moist_density(                                                  &
     3688             t_rho = input_buffer(group % n_output_quantities) % array(:,:,:), &
     3689             p = input_buffer(2) % array(:,:,:),                               &
     3690             qv = input_buffer(3) % array(:,:,:),                              &
     3691             rd = RD,                                                          &
     3692             rv = RV                                                           &
     3693          )
     3694
     3695          ! qv
     3696          input_buffer(3) % is_preprocessed = .TRUE.
     3697
     3698          ! density
     3699          input_buffer(group % n_output_quantities) % is_preprocessed = .TRUE.
     3700
     3701
    32123702          message = "Input buffers for group '" // TRIM(group % kind) // "'"//&
    32133703             " preprocessed sucessfully."
     
    33363826
    33373827       CASE DEFAULT
    3338           message = "Buffer kind '" // TRIM(group % kind) // "' is not supported."
     3828          message = "IO group kind '" // TRIM(group % kind) // "' is not supported."
    33393829          CALL abort('prerpocess', message)
    33403830
Note: See TracChangeset for help on using the changeset viewer.