Changeset 4796


Ignore:
Timestamp:
Nov 25, 2020 10:54:32 PM (3 years ago)
Author:
gronemeier
Message:

add an example to create_basic_static_driver.py which creates static driver of urban_environment test case; update static driver of urban_environment test case

Location:
palm/trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SCRIPTS/create_basic_static_driver.py

    r4370 r4796  
    1 #------------------------------------------------------------------------------#
     1#!/usr/bin/env python3
     2# -------------------------------------------------------------------- #
    23# This file is part of the PALM model system.
    34#
    4 # PALM is free software: you can redistribute it and/or modify it under the
    5 # terms of the GNU General Public License as published by the Free Software
    6 # Foundation, either version 3 of the License, or (at your option) any later
    7 # version.
     5# PALM is free software: you can redistribute it and/or modify it under
     6# the terms of the GNU General Public License as published by the Free
     7# Software Foundation, either version 3 of the License, or (at your
     8# option) any later version.
    89#
    9 # PALM is distributed in the hope that it will be useful, but WITHOUT ANY
    10 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
    11 # A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
     10# PALM is distributed in the hope that it will be useful, but WITHOUT
     11# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
     12# FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
     13# for more details.
    1214#
    13 # You should have received a copy of the GNU General Public License along with
    14 # PALM. If not, see <http://www.gnu.org/licenses/>.
     15# You should have received a copy of the GNU General Public License
     16# along with PALM. If not, see <http://www.gnu.org/licenses/>.
    1517#
    1618# Copyright 1997-2020 Leibniz Universitaet Hannover
    17 #------------------------------------------------------------------------------#
    18 
     19# -------------------------------------------------------------------- #
     20
     21import datetime
    1922import math
    20 from   netCDF4 import Dataset
     23import pytz
     24
     25from netCDF4 import Dataset
    2126import numpy as np
    22 import datetime
    2327
    2428
    2529class StaticDriver:
    26     """ This is an example script to generate static drivers for PALM.
    27    
     30    """This is an example script to generate static drivers for PALM.
     31
    2832    You can use it as a starting point for creating your setup specific
    29     driver. 
     33    driver.
    3034    """
    31    
    32    
     35
    3336    def __init__(self):
    34         """ Open the static driver as NetCDF4 file. Here, you have to give the
    35         full path to the static driver that shall be created. Existing file
    36         with same name is deleted.
     37        """Open the static driver as netCDF4 file. Here, you have to
     38        give the full path to the static driver that shall be created.
     39        Existing file with same name is deleted.
    3740        """
    3841        print('Opening file...')
    39         self.nc_file = Dataset('path/to/file.nc', 'w', format='NETCDF4')
    40 
     42        self.nc_file = Dataset('example_static_file.nc', 'w', format='NETCDF4')
    4143
    4244    def write_global_attributes(self):
    43         """ Write global attributes to static driver. """
     45        """Write global attributes to static driver."""
    4446        print("Writing global attributes...")
    45        
    46         # mandatory global attributes
    47         self.nc_file.origin_lon = 55.0 # used to initialize coriolis parameter
    48         self.nc_file.origin_lat = 0.0 # (overwrite initialization_parameters)
    49         self.nc_file.origin_time = '2000-06-21 12:00:00 +00'
    50         self.nc_file.origin_x = 308124
    51         self.nc_file.origin_y = 6098908
     47
     48        # Optional global attributes
     49        # --------------------------
     50        self.nc_file.title = 'Example PALM static driver'
     51        self.nc_file.author = 'PALM user'
     52        self.nc_file.institution = 'Institut of Meteorology and Climatology,' \
     53            'Leibniz University Hannover'
     54        self.nc_file.comment = 'Generic crossing example'
     55        self.nc_file.creation_date = \
     56            pytz.utc.localize(datetime.datetime.utcnow()).strftime('%y-%m-%d %H:%M:%S %z')[:-2]
     57        self.nc_file.history = ''
     58        self.nc_file.keywords = 'example, PALM-4U'
     59        self.nc_file.license = ''
     60        self.nc_file.palm_version = ''
     61        self.nc_file.references = ''
     62        self.nc_file.source = 'PALM trunk'
     63        self.nc_file.version = '1'
     64
     65        # Mandatory global attributes
     66        # ---------------------------
     67        self.nc_file.Conventions = 'CF-1.7'
     68        self.nc_file.origin_lat = 52.50965  # (overwrite initialization_parameters)
     69        self.nc_file.origin_lon = 13.3139  # Used to initialize Coriolis parameter
     70        self.nc_file.origin_time = '2019-03-06 10:00:00 +00'
     71        self.nc_file.origin_x = 3455249.0
     72        self.nc_file.origin_y = 5424815.0
    5273        self.nc_file.origin_z = 0.0
    5374        self.nc_file.rotation_angle = 0.0
    54        
    55         # optional global attributes
    56         self.nc_file.author = 'Your Name'
    57         self.nc_file.comment = 'Miscellaneous information about the data ' \
    58                                'or methods to produce it.'
    59         self.nc_file.creation_date = str(datetime.datetime.now())
    60         self.nc_file.institution = 'Institut of Meteorology and Climatology,' \
    61                                    'Leibniz University Hannover'
    62         self.nc_file.history = ''
    63         self.nc_file.palm_revision = ''
    64         self.nc_file.title = 'Static driver for some arbitrary PALM setup'
    65 
    6675
    6776    def define_dimensions(self):
    68         """ Set dimensions on which variables are defined. """
     77        """Set dimensions on which variables are defined."""
    6978        print("Writing dimensions...")
    70        
    71         # specify general grid parameters
    72         # these values must equal those set in the initialization_parameters
    73         self.nx = 39
    74         self.ny = 39
    75         self.nz = 40
    76         dx = 50
    77         dy = 50
    78         dz = 50
    79        
    80         # create soil grid (only relevant if land surface module is used)
     79
     80        # Specify general grid parameters
     81        # These values must equal to those set in the initialization_parameters
     82        self.nx = 19
     83        self.ny = 19
     84        self.nz = 20
     85        dx = 2
     86        dy = 2
     87        dz = 2
     88
     89        # Create soil grid (only relevant if land surface module is used)
    8190        dz_soil = np.array((0.01, 0.02, 0.04, 0.06, 0.14, 0.26, 0.54, 1.86))
    8291        zsoil_fullLayers = np.zeros_like(dz_soil)
    83         zsoil_fullLayers = np.around([np.sum(dz_soil[:zs]) for zs in np.arange(1,len(dz_soil)+1)],2)
     92        zsoil_fullLayers = np.around(
     93            [np.sum(dz_soil[:zs]) for zs in np.arange(1, len(dz_soil)+1)], 2)
    8494        zsoil_array = zsoil_fullLayers - dz_soil/2.
    85        
    86         # create plant canopy grid (only relevant of plant canopy module is used)
    87         zlad_array = np.arange(0,10,1)
    88        
    89         # mandatory dimensions
    90         self.nc_file.createDimension('x' ,self.nx+1)
    91         self.x = self.nc_file.createVariable('x', 'f8', ('x',))
     95
     96        # Coordinates
     97        # -----------
     98        self.nc_file.createDimension('x', self.nx+1)
     99        self.x = self.nc_file.createVariable('x', 'f4', ('x',))
     100        self.x.long_name = 'distance to origin in x-direction'
    92101        self.x.units = 'm'
    93         self.x.standard_name = 'x coordinate of cell centers'
    94         self.x[:] = np.arange(0,(self.nx+1)*dx,dx)
    95        
    96         self.nc_file.createDimension('xu', self.nx)
    97         self.xu = self.nc_file.createVariable('xu', 'f8', ('xu',))
    98         self.xu.units = 'm'
    99         self.xu.standard_name = 'x coordinate of cell edges'
    100         self.xu[:] = np.arange(dx/2.,self.nx*dx,dx)
    101        
     102        self.x.axis = 'X'
     103        self.x[:] = np.arange(0, (self.nx+1)*dx, dx) + 0.5 * dx
     104
    102105        self.nc_file.createDimension('y', self.ny+1)
    103         self.y = self.nc_file.createVariable('y', 'f8', ('y',))
     106        self.y = self.nc_file.createVariable('y', 'f4', ('y',))
     107        self.y.long_name = 'distance to origin in y-direction'
    104108        self.y.units = 'm'
    105         self.y.standard_name = 'y coordinate of cell centers'
    106         self.y[:] = np.arange(0,(self.ny+1)*dy,dy)
    107        
    108         self.nc_file.createDimension('yv', self.ny)
    109         self.yv = self.nc_file.createVariable('yv', 'f8', ('yv',))
    110         self.yv.units = 'm'
    111         self.yv.standard_name = 'y coordinate of cell edges'
    112         self.yv[:] = np.arange(dy/2.,self.ny*dy,dy)
    113        
    114         # if your simulation uses a stretched vertical grid, you need to
    115         # modify the z and zw coordinates,
    116         # e.g. z_array = (...) and zw_array = (...)
    117         z_array = np.append(0, np.arange(dz/2,(self.nz+1)*dz,dz))
    118         self.nc_file.createDimension('z',self.nz+2)
    119         self.z = self.nc_file.createVariable('z', 'f8', ('z',))
     109        self.y.axis = 'Y'
     110        self.y[:] = np.arange(0, (self.ny+1)*dy, dy) + 0.5 * dy
     111
     112        # NOTE if your simulation uses a stretched vertical grid, you
     113        # need to modify the z coordinates, e.g. z_array = (...)
     114        z_array = np.append(0, np.arange(dz/2, (self.nz)*dz, dz))
     115        self.nc_file.createDimension('z', self.nz+1)
     116        self.z = self.nc_file.createVariable('z', 'f4', ('z',))
     117        self.z.long_name = 'height above origin'
    120118        self.z.units = 'm'
    121         self.z.standard_name = 'z coordinate of cell centers'
     119        self.z.axis = 'Z'
     120        self.z.positive = 'up'
    122121        self.z[:] = z_array
    123122
    124         zw_array = np.arange(0,(self.nz+2)*dz,dz)
    125         self.nc_file.createDimension('zw', self.nz+2)
    126         self.zw = self.nc_file.createVariable('zw', 'f8', ('zw',))
    127         self.zw.units = 'm'
    128         self.zw.standard_name = 'z coordinate of cell edges'
    129         self.zw[:] = zw_array
    130        
    131         # optional dimensions, uncomment if needed
    132         #self.nc_file.createDimension('zsoil', len(zsoil_array))
    133         #self.zsoil = self.nc_file.createVariable('zsoil', 'f8', ('zsoil',))
    134         #self.zsoil.positive = 'down'
    135         #self.zsoil.units = 'm'
    136         #self.zsoil.standard_name = 'depth below land'
    137         #self.zsoil[:] = zsoil_array
    138        
    139         #self.nc_file.createDimension('zlad', len(zlad_array))
    140         #self.zlad = self.nc_file.createVariable('zlad', 'f8', ('zlad',))
    141         #self.zlad.units = 'm'
    142         #self.zlad.standard_name = 'z coordinate of resolved plant canopy'
    143         #self.zlad[:] = zlad_array
    144        
    145         #self.nc_file.createDimension('nalbedo_pars', 3)
    146         #self.nalbedo_pars = self.nc_file.createVariable(
    147             #'nalbedo_pars', 'i1', ('nalbedo_pars',))
    148         #self.nalbedo_pars[:] = np.arange(0,3)
    149        
    150         #self.nc_file.createDimension('nbuilding_pars', 46)
    151         #self.nbuilding_pars = self.nc_file.createVariable(
    152             #'nbuilding_pars', 'i4', ('nbuilding_pars',))
    153         #self.nbuilding_pars[:] = np.arange(0,46)
    154        
    155         #self.nc_file.createDimension('nsoil_pars', 8)
    156         #self.nsoil_pars = self.nc_file.createVariable(
    157             #'nsoil_pars', 'i1', ('nsoil_pars',))
    158         #self.nsoil_pars[:] = np.arange(0,8)
    159        
    160         #self.nc_file.createDimension('nsurface_fraction', 3)
    161         #self.nsurface_fraction = self.nc_file.createVariable(
    162             #'nsurface_fraction', 'i1', ('nsurface_fraction',))
    163         #self.nsurface_fraction[:] = np.arange(0,3)
    164        
    165         #self.nc_file.createDimension('nvegetation_pars', 12)
    166         #self.nvegetation_pars = self.nc_file.createVariable(
    167             #'nvegetation_pars', 'i1', ('nvegetation_pars',))
    168         #self.nvegetation_pars[:] = np.arange(0,12)
    169        
    170         #self.nc_file.createDimension('npavement_pars', 4)
    171         #self.npavement_pars = self.nc_file.createVariable(
    172             #'npavement_pars', 'i1', ('npavement_pars',))
    173         #self.npavement_pars[:] = np.arange(0,4)
    174 
    175         #self.nc_file.createDimension('npavement_subsurface_pars', 2)
    176         #self.npavement_subsurface_pars = self.nc_file.createVariable(
    177             #'npavement_subsurface_pars', 'i1', ('npavement_subsurface_pars',))
    178         #self.npavement_subsurface_pars[:] = np.arange(0,2)
    179        
    180         #self.nc_file.createDimension('nwater_pars', 6)
    181         #self.nwater_pars = self.nc_file.createVariable(
    182             #'nwater_pars', 'i1', ('nwater_pars',))
    183         #self.nwater_pars[:] = np.arange(0,6)
    184        
    185 
    186     def add_variables(self):
    187         """ Uncomment variables below as you like.
    188        
    189         Be aware that some variables depend on others. For a description of
    190         each variable, please have a look at the documentation at
     123        zlad_array = self.z[:6]
     124        self.nc_file.createDimension('zlad', len(zlad_array))
     125        self.zlad = self.nc_file.createVariable('zlad', 'f4', ('zlad',))
     126        self.zlad.long_name = 'height above ground'
     127        self.zlad.units = 'm'
     128        self.zlad.axis = 'Z'
     129        self.zlad.positive = 'up'
     130        self.zlad[:] = zlad_array
     131
     132        # self.nc_file.createDimension('zsoil', len(zsoil_array))
     133        # self.zsoil = self.nc_file.createVariable('zsoil', 'f4', ('zsoil',))
     134        # self.zsoil.long_name = 'depth in the soil'
     135        # self.zsoil.units = 'm'
     136        # self.zsoil.axis = 'Z'
     137        # self.zsoil.positive = 'down'
     138        # self.zsoil[:] = zsoil_array
     139
     140        # Other dimensions
     141        # ----------------
     142        # self.nc_file.createDimension('nalbedo_pars', 3)
     143        # self.nalbedo_pars = self.nc_file.createVariable(
     144        #     'nalbedo_pars', 'i4', ('nalbedo_pars',))
     145        # self.nalbedo_pars[:] = np.arange(8)
     146        #
     147        # self.nc_file.createDimension('nbuilding_pars', 46)
     148        # self.nbuilding_pars = self.nc_file.createVariable(
     149        #     'nbuilding_pars', 'i4', ('nbuilding_pars',))
     150        # self.nbuilding_pars[:] = np.arange(46)
     151        #
     152        # self.nc_file.createDimension('npavement_pars', 4)
     153        # self.npavement_pars = self.nc_file.createVariable(
     154        #     'npavement_pars', 'i4', ('npavement_pars',))
     155        # self.npavement_pars[:] = np.arange(4)
     156        #
     157        # self.nc_file.createDimension('npavement_subsurface_pars', 2)
     158        # self.npavement_subsurface_pars = self.nc_file.createVariable(
     159        #     'npavement_subsurface_pars', 'i4', ('npavement_subsurface_pars',))
     160        # self.npavement_subsurface_pars[:] = np.arange(2)
     161        #
     162        # self.nc_file.createDimension('nsoil_pars', 8)
     163        # self.nsoil_pars = self.nc_file.createVariable(
     164        #     'nsoil_pars', 'i4', ('nsoil_pars',))
     165        # self.nsoil_pars[:] = np.arange(8)
     166        #
     167        self.nc_file.createDimension('nsurface_fraction', 3)
     168        self.nsurface_fraction = self.nc_file.createVariable(
     169            'nsurface_fraction', 'i4', ('nsurface_fraction',))
     170        self.nsurface_fraction[:] = np.arange(3)
     171
     172        # self.nc_file.createDimension('nvegetation_pars', 12)
     173        # self.nvegetation_pars = self.nc_file.createVariable(
     174        #     'nvegetation_pars', 'i4', ('nvegetation_pars',))
     175        # self.nvegetation_pars[:] = np.arange(12)
     176        #
     177        # self.nc_file.createDimension('nwater_pars', 6)
     178        # self.nwater_pars = self.nc_file.createVariable(
     179        #     'nwater_pars', 'i4', ('nwater_pars',))
     180        # self.nwater_pars[:] = np.arange(7)
     181
     182    def define_variables(self):
     183        """Define variables for the static driver.
     184
     185        Be aware that some variables depend on others. For a description
     186        of each variable, please have a look at the documentation at
    191187        palm.muk.uni-hannover.de/trac/wiki/doc/app/iofiles/pids/static
    192        
     188
    193189        An example of how you modify the variables is given below:
    194        
    195         building_2d_array = np.ones((self.ny+1,self.nx+1)) * -9999.9
     190
     191        building_2d_array = np.ones((self.ny+1, self.nx+1)) * -9999.0
    196192        south_wall, north_wall, left_wall, right_wall = 20, 25, 20, 25
    197         building_2d_array[south_wall:north_wall,left_wall:right_wall] = 50
     193        building_2d_array[
     194            south_wall:north_wall,
     195            left_wall:right_wall
     196            ] = 50
    198197        nc_buildings_2d = self.nc_file.createVariable(
    199             'buildings_2d', 'f4', ('y','x'),fill_value=-9999.9)
     198            'buildings_2d', 'f4', ('y', 'x'), fill_value=-9999.0)
    200199        nc_buildings_2d.lod = 1
    201         nc_buildings_2d[:,:] = building_2d_array
     200        nc_buildings_2d[:, :] = building_2d_array
    202201        """
    203202        print("Writing variables...")
    204        
    205         ## topography variables
    206         #nc_buildings_2d = self.nc_file.createVariable(
    207             #'buildings_2d', 'f4', ('y','x'), fill_value=-9999.9)
    208         #nc_buildings_2d.lod = 1
    209         #nc_buildings_2d[:,:] = -9999.9
    210        
    211         #nc_buildings_3d = self.nc_file.createVariable(
    212             #'buildings_3d', 'i1', ('z','y','x'),fill_value=-127)
    213         #nc_buildings_3d.lod = 2
    214         #nc_buildings_3d[:,:,:] = -127
    215        
    216         #nc_building_id = self.nc_file.createVariable(
    217             #'building_id', 'i1', ('y','x'), fill_value=-127)
    218         #nc_building_id[:,:] = -127
    219        
    220         #nc_zt = self.nc_file.createVariable(
    221             #'zt', 'f4', ('y','x'), fill_value=-9999.9)
    222         #nc_zt.long_name = 'orography'
    223         #nc_zt.units = 'm'
    224         #nc_zt[:,:] = -9999.9
    225        
    226         ## surface variables
    227         #nc_albedo_pars = self.nc_file.createVariable(
    228             #'albedo_pars', 'f4', ('nalbedo_pars','y','x'), fill_value=-9999.9)
    229         #nc_albedo_pars.long_name = 'albedo parameters'
    230         #nc_albedo_pars[:,:,:] = -9999.9
    231        
    232         #nc_albedo_type = self.nc_file.createVariable(
    233             #'albedo_type', 'f4', ('y','x'), fill_value=-9999.9)
    234         #nc_albedo_type[:,:] = -9999.9
    235        
    236         #nc_building_pars = self.nc_file.createVariable(
    237             #'building_pars', 'f4', ('nbuilding_pars','y','x'), fill_value=-9999.9)
    238         #nc_building_pars.long_name = 'building parameters'
    239         #nc_building_pars[:,:,:] = -9999.9
    240        
    241         #nc_building_type = self.nc_file.createVariable(
    242             #'building_type', 'i1', ('y','x'), fill_value=-127)
    243         #nc_building_type[:,:] = -127
    244        
    245         #nc_pavement_pars = self.nc_file.createVariable(
    246             #'pavement_pars', 'f4', ('npavement_pars','y','x'), fill_value=-9999.9)
    247         #nc_pavement_pars[:,:,:] = -9999.9
    248        
    249         #nc_pavement_subsurface_pars = self.nc_file.createVariable(
    250             #'pavement_subsurface_pars', 'i1',
    251             #('npavement_subsurface_pars','zsoil','y','x'), fill_value=-9999.9)
    252         #nc_pavement_subsurface_pars[:,:,:,:] = -9999.9
    253        
    254         #nc_pavement_type = self.nc_file.createVariable(
    255             #'pavement_type', 'i1', ('y','x'), fill_value=-127)
    256         #nc_pavement_type[:,:] = -127
    257        
    258         #nc_soil_pars = self.nc_file.createVariable(
    259             #'soil_pars', 'f', ('nsoil_pars','zsoil','y','x'), fill_value=-9999.9)
    260         #nc_soil_pars.long_name = 'soil parameters'
    261         #nc_soil_pars.lod = 2 # if set to 1, adjust dimensions of soil_pars
    262         #nc_soil_pars[:,:,:,:] = -9999.9
    263        
    264         #nc_soil_type = self.nc_file.createVariable(
    265             #'soil_type', 'i1', ('zsoil','y','x'), fill_value=-127)
    266         #nc_soil_type.lod = 2 # if set to 1, adjust dimensions of soil_type
    267         #nc_soil_type[:,:,:] = -127
    268        
    269         ## required, if several surface types are defined at the same location
    270         #nc_surface_fraction = self.nc_file.createVariable(
    271             #'surface_fraction', 'f', ('nsurface_fraction','y','x'), fill_value=-9999.9)
    272         #nc_surface_fraction[0,:,:] = -9999.9 # vegetation fraction
    273         #nc_surface_fraction[1,:,:] = -9999.9 # pavement fraction
    274         #nc_surface_fraction[2,:,:] = -9999.9 # water fraction
    275        
    276         #nc_vegetation_pars = self.nc_file.createVariable(
    277             #'vegetation_pars', 'f', ('nvegetation_pars','y','x'), fill_value=-9999.9)
    278         #nc_vegetation_pars.long_name = 'vegetation parameters'
    279         #nc_vegetation_pars[:,:,:] = -9999.9
    280        
    281         #nc_vegetation_type = self.nc_file.createVariable(
    282             #'vegetation_type', 'i1', ('y','x'), fill_value=-127)
    283         #nc_vegetation_type[:,:] = -127
    284        
    285         #nc_street_type = self.nc_file.createVariable(
    286             #'street_type', 'i1', ('y','x'),
    287             #fill_value=-127)
    288         #nc_street_type[:,:] = -127
    289        
    290         #nc_water_pars = self.nc_file.createVariable(
    291             #'water_pars', 'i1', ('nwater_pars','y','x'), fill_value=-9999.9)
    292         #nc_water_pars[:,:,:] = -9999.9
    293        
    294         #nc_water_type = self.nc_file.createVariable(
    295             #'water_type', 'i1',('y','x'), fill_value=-127)
    296         #nc_water_type[:,:] = -127
    297        
    298         ## resolved plant canopy
    299         #nc_lad = self.nc_file.createVariable(
    300             #'lad', 'f4', ('zlad','y','x'), fill_value=-9999.9)
    301         #nc_lad.long_name = 'leaf area density'
    302         #nc_lad[:,:,:] = -9999.9
    303        
    304         #nc_root_area_dens_s = self.nc_file.createVariable(
    305             #'root_area_dens_s', 'f4', ('zsoil','y','x'), fill_value=-9999.9)
    306         #nc_root_area_dens_s.long_name = 'parameterized root area density'
    307         #nc_root_area_dens_s.units = '1/m'
    308         #nc_root_area_dens_s[:,:,:] = -9999.9
    309        
    310        
     203
     204        # Topography set-up
     205        # -----------------
     206        nc_building_id = self.nc_file.createVariable(
     207            'building_id', 'i4', ('y', 'x'), fill_value=-9999)
     208        nc_building_id.long_name = "building id number"
     209        nc_building_id.units = "1"
     210        nc_building_id[:, :] = nc_building_id._FillValue
     211
     212        nc_buildings_2d = self.nc_file.createVariable(
     213            'buildings_2d', 'f4', ('y', 'x'), fill_value=-9999.0)
     214        nc_buildings_2d.long_name = "building height"
     215        nc_buildings_2d.units = "m"
     216        nc_buildings_2d.lod = np.int32(1)
     217        nc_buildings_2d[:, :] = nc_buildings_2d._FillValue
     218
     219        nc_buildings_3d = self.nc_file.createVariable(
     220            'buildings_3d', 'i1', ('z', 'y', 'x'), fill_value=-127)
     221        nc_buildings_3d.long_name = "building flag"
     222        nc_buildings_3d.units = "1"
     223        nc_buildings_3d.lod = np.int32(2)
     224        nc_buildings_3d[:, :, :] = nc_buildings_3d._FillValue
     225
     226        nc_zt = self.nc_file.createVariable(
     227            'zt', 'f4', ('y', 'x'), fill_value=-9999.0)
     228        nc_zt.long_name = 'terrain height'
     229        nc_zt.units = 'm'
     230        nc_zt[:, :] = nc_zt._FillValue
     231
     232        # nc_z0 = self.nc_file.createVariable(
     233        #     'z0', 'f4', ('y', 'x'), fill_value=-9999.0)
     234        # nc_z0.long_name = 'roughness length for momentum'
     235        # nc_z0.units = 'm'
     236        # nc_z0[:, :] = nc_z0._FillValue
     237
     238        # Main surface clasification
     239        # --------------------------
     240        # nc_albedo_type = self.nc_file.createVariable(
     241        #     'albedo_type', 'f4', ('y', 'x'), fill_value=-9999.0)
     242        # nc_albedo_type.long_name = "albedo type"
     243        # nc_albedo_type.units = "1"
     244        # nc_albedo_type[:, :] = nc_albedo_type._FillValue
     245        #
     246        nc_building_type = self.nc_file.createVariable(
     247            'building_type', 'i1', ('y', 'x'), fill_value=-127)
     248        nc_building_type.long_name = "building type classification"
     249        nc_building_type.units = "1"
     250        nc_building_type[:, :] = nc_building_type._FillValue
     251
     252        nc_pavement_type = self.nc_file.createVariable(
     253            'pavement_type', 'i1', ('y', 'x'), fill_value=-127)
     254        nc_pavement_type.long_name = "pavement type classification"
     255        nc_pavement_type.units = "1"
     256        nc_pavement_type[:, :] = nc_pavement_type._FillValue
     257
     258        nc_soil_type = self.nc_file.createVariable(
     259            'soil_type', 'i1', ('y', 'x'), fill_value=-127)
     260        nc_soil_type.long_name = "soil type classification"
     261        nc_soil_type.units = "1"
     262        nc_soil_type.lod = np.int32(1)
     263        nc_soil_type[:, :] = nc_soil_type._FillValue
     264
     265        # NOTE alternatively, define different soil types for different
     266        # soil layers
     267        # nc_soil_type = self.nc_file.createVariable(
     268        #     'soil_type', 'i1', ('zsoil', 'y', 'x'), fill_value=-127)
     269        # nc_soil_type.long_name = "soil type"
     270        # nc_soil_type.units = "1"
     271        # nc_soil_type.lod = np.int32(2)
     272        # nc_soil_type[:, :, :] = nc_soil_type._FillValue
     273        #
     274        # nc_street_crossing = self.nc_file.createVariable(
     275        #     'street_crossing', 'i1', ('y', 'x'),
     276        #     fill_value=-127)
     277        # nc_street_crossing.long_name = "street crossing"
     278        # nc_street_crossing.units = "1"
     279        # nc_street_crossing.valid_range = np.byte(1)
     280        # nc_street_crossing[:, :] = nc_street_crossing._FillValue
     281        #
     282        nc_street_type = self.nc_file.createVariable(
     283            'street_type', 'i1', ('y', 'x'),
     284            fill_value=-127)
     285        nc_street_type.long_name = "street type classification"
     286        nc_street_type.units = "1"
     287        nc_street_type[:, :] = nc_street_type._FillValue
     288
     289        nc_surface_fraction = self.nc_file.createVariable(
     290            'surface_fraction', 'f4', ('nsurface_fraction', 'y', 'x'), fill_value=-9999.0)
     291        nc_surface_fraction.long_name = "surface fraction"
     292        nc_surface_fraction.units = "1"
     293        nc_surface_fraction[0, :, :] = nc_surface_fraction._FillValue  # vegetation fraction
     294        nc_surface_fraction[1, :, :] = nc_surface_fraction._FillValue  # pavement fraction
     295        nc_surface_fraction[2, :, :] = nc_surface_fraction._FillValue  # water fraction
     296
     297        nc_vegetation_type = self.nc_file.createVariable(
     298            'vegetation_type', 'i1', ('y', 'x'), fill_value=-127)
     299        nc_vegetation_type.long_name = "vegetation type classification"
     300        nc_vegetation_type.units = "1"
     301        nc_vegetation_type[:, :] = nc_vegetation_type._FillValue
     302
     303        nc_water_type = self.nc_file.createVariable(
     304            'water_type', 'i1', ('y', 'x'), fill_value=-127)
     305        nc_water_type.long_name = "water type classification"
     306        nc_water_type.units = "1"
     307        nc_water_type[:, :] = nc_water_type._FillValue
     308
     309        # Pixel-based surface and sub-surface parameters
     310        # ----------------------------------------------
     311        # nc_albedo_pars = self.nc_file.createVariable(
     312        #     'albedo_pars', 'f4', ('nalbedo_pars', 'y', 'x'), fill_value=-9999.0)
     313        # nc_albedo_pars.long_name = 'albedo parameters'
     314        # nc_albedo_pars[:, :, :] = nc_albedo_pars._FillValue
     315        #
     316        # nc_building_pars = self.nc_file.createVariable(
     317        #     'building_pars', 'f4', ('nbuilding_pars', 'y', 'x'), fill_value=-9999.0)
     318        # nc_building_pars.long_name = 'building parameters'
     319        # nc_building_pars[:, :, :] = nc_building_pars._FillValue
     320        #
     321        # nc_pavement_pars = self.nc_file.createVariable(
     322        #     'pavement_pars', 'f4', ('npavement_pars', 'y', 'x'), fill_value=-9999.0)
     323        # nc_pavement_pars.long_name = "pavement parameters"
     324        # nc_pavement_pars[:, :, :] = nc_pavement_pars._FillValue
     325        #
     326        # nc_pavement_subsurface_pars = self.nc_file.createVariable(
     327        #     'pavement_subsurface_pars', 'i1',
     328        #     ('npavement_subsurface_pars','zsoil', 'y', 'x'), fill_value=-9999.0)
     329        # nc_pavement_subsurface_pars.long_name = "pavement subsurface parameters"
     330        # nc_pavement_subsurface_pars[:, :, :, :] = nc_pavement_subsurface_pars._FillValue
     331        #
     332        # nc_soil_pars = self.nc_file.createVariable(
     333        #     'soil_pars', 'f', ('nsoil_pars','zsoil', 'y', 'x'), fill_value=-9999.0)
     334        # nc_soil_pars.long_name = "soil parameters"
     335        # nc_soil_pars.lod = np.int32(2) # if set to 1, adjust dimensions of soil_pars
     336        # nc_soil_pars[:, :, :, :] = nc_soil_pars._FillValue
     337        #
     338        # nc_vegetation_pars = self.nc_file.createVariable(
     339        #     'vegetation_pars', 'f', ('nvegetation_pars', 'y', 'x'), fill_value=-9999.0)
     340        # nc_vegetation_pars.long_name = 'vegetation parameters'
     341        # nc_vegetation_pars[:, :, :] = nc_vegetation_pars._FillValue
     342        #
     343        # nc_water_pars = self.nc_file.createVariable(
     344        #     'water_pars', 'i1', ('nwater_pars', 'y', 'x'), fill_value=-9999.0)
     345        # nc_water_pars.long_name = "water parameters"
     346        # nc_water_pars[:, :, :] = nc_water_pars._FillValue
     347
     348        # Vegetation parameters
     349        # ---------------------
     350        nc_lad = self.nc_file.createVariable(
     351            'lad', 'f4', ('zlad', 'y', 'x'), fill_value=-9999.0)
     352        nc_lad.long_name = "leaf area density"
     353        nc_lad.units = "m2 m-3"
     354        nc_lad[:, :, :] = nc_lad._FillValue
     355
     356        # nc_bad = self.nc_file.createVariable(
     357        #     'bad', 'f4', ('zlad', 'y', 'x'), fill_value=-9999.0)
     358        # nc_bad.long_name = "basal area density"
     359        # nc_bad.units = "m2 m-3"
     360        # nc_bad[:, :, :] = nc_bad._FillValue
     361        #
     362        # nc_root_area_dens_r = self.nc_file.createVariable(
     363        #     'root_area_dens_s', 'f4', ('zsoil', 'y', 'x'), fill_value=-9999.0)
     364        # nc_root_area_dens_r.long_name = 'root area density of resolved vegetation'
     365        # nc_root_area_dens_r.units = '1'
     366        # nc_root_area_dens_r[:, :, :] = nc_root_area_dens_r._FillValue
     367        #
     368        # nc_root_area_dens_s = self.nc_file.createVariable(
     369        #     'root_area_dens_s', 'f4', ('zsoil', 'y', 'x'), fill_value=-9999.0)
     370        # nc_root_area_dens_s.long_name = 'root area density of parameterized vegetation'
     371        # nc_root_area_dens_s.units = '1'
     372        # nc_root_area_dens_s[:, :, :] = nc_root_area_dens_s._FillValue
     373        #
     374        # nc_tree_id = self.nc_file.createVariable(
     375        #     'tree_id', 'i4', ('y', 'x'), fill_value=-9999.0)
     376        # nc_tree_id.long_name = "tree id"
     377        # nc_tree_id.units = "1"
     378        # nc_tree_id[:, :, :] = nc_tree_id._FillValue
     379
     380        # Set topography
     381        # --------------
     382        nc_building_id[:5, :5] = 1
     383        nc_building_id[-5:, :5] = 2
     384        nc_building_id[:5, -5:] = 3
     385        nc_building_id[-5:, -5:] = 4
     386
     387        nc_buildings_2d[:, :] = np.where(
     388            nc_building_id[:, :] > nc_building_id._FillValue,
     389            nc_building_id[:, :] * 10.0,
     390            nc_buildings_2d[:, :])
     391
     392        for k in range(self.nz+1):
     393            nc_buildings_3d[k, :, :] = np.where(nc_buildings_2d[:, :] > self.z[k], 1, 0)
     394
     395        nc_zt[:, :] = 4.0
     396
     397        # Set surface types
     398        # -----------------
     399        nc_building_type[:, :] = np.where(
     400            nc_building_id[:, :] > nc_building_id._FillValue,
     401            nc_building_id[:, :] * 1,
     402            nc_building_type[:, :])
     403
     404        nc_pavement_type[9:11, :5] = 1
     405        nc_pavement_type[:, 7:13] = 2
     406
     407        nc_vegetation_type[:, 5:7] = 3
     408        nc_vegetation_type[:, 13:15] = 3
     409        nc_vegetation_type[5:15, 15:] = 3
     410        nc_vegetation_type[7:9, 15:17] = 1
     411        nc_vegetation_type[11:13, 15:17] = 1
     412
     413        nc_water_type[5:9, :5] = 2
     414        nc_water_type[11:15, :5] = 2
     415
     416        nc_soil_type[:, :] = np.where(
     417            nc_vegetation_type[:, :] > nc_vegetation_type._FillValue,
     418            2,
     419            nc_soil_type[:, :])
     420        nc_soil_type[:, :] = np.where(
     421            nc_pavement_type[:, :] > nc_pavement_type._FillValue,
     422            2,
     423            nc_soil_type[:, :])
     424
     425        nc_street_type[:, :] = np.where(nc_pavement_type[:, :] == 1, 11, nc_street_type[:, :])
     426        nc_street_type[:, :] = np.where(nc_pavement_type[:, :] == 2, 13, nc_street_type[:, :])
     427
     428        nc_surface_fraction[0, :, :] = np.where(
     429            nc_building_id[:, :] > nc_building_id._FillValue,
     430            nc_surface_fraction[0, :, :],
     431            0)
     432        nc_surface_fraction[2, :, :] = nc_surface_fraction[1, :, :] = nc_surface_fraction[0, :, :]
     433        nc_surface_fraction[0, :, :] = np.where(
     434            nc_vegetation_type[:, :] > nc_vegetation_type._FillValue,
     435            1,
     436            nc_surface_fraction[0, :, :])
     437        nc_surface_fraction[1, :, :] = np.where(
     438            nc_pavement_type[:, :] > nc_pavement_type._FillValue,
     439            1,
     440            nc_surface_fraction[1, :, :])
     441        nc_surface_fraction[2, :, :] = np.where(
     442            nc_water_type[:, :] > nc_water_type._FillValue,
     443            1,
     444            nc_surface_fraction[2, :, :])
     445
     446        # Set vegetation
     447        # --------------
     448        lad_profile = [0.0, 0.01070122, 0.1070122, 0.3130108, 0.3879193, 0.1712195]
     449        for k in range(len(self.zlad)):
     450            nc_lad[k, 7:9, 15:17] = lad_profile[k]
     451            nc_lad[k, 11:13, 15:17] = lad_profile[k]
     452
    311453    def finalize(self):
    312         """ Close file """
     454        """Close file."""
    313455        print("Closing file...")
    314        
     456
    315457        self.nc_file.close()
    316458
     
    320462    driver.write_global_attributes()
    321463    driver.define_dimensions()
    322     driver.add_variables()
     464    driver.define_variables()
    323465    driver.finalize()
    324    
Note: See TracChangeset for help on using the changeset viewer.