run UMEP.ShadowGenerator from python

[s-martiz]

Dear UMEP-community,

I'm trying to run UMEP from python to avoid using QGIS and facilitate debugging.
My question is quite simple > Is there a simple example of the shadow_generator.py that I can run from python to better understand how it works ?
Let’s say something like :
ShadowGenerator(geojson_object, position, date, hour) > .tiff object representing the shadow of the zone.

Simon

[biglimp]

Hi Simon,
Here is some code to debug the shadow algorithm. This code requires gdal but actually, any 2D numpy array would do as the dsm variable. All the functions are available in the UMEP tool.

import numpy as np
from osgeo import gdal
import matplotlib.pyplot as plt
from Utilities.SEBESOLWEIGCommonFiles import shadowingfunctions as sh
import time
# import wallalgorithms as wa # if you want shadow on walls. Otherwise, not needed. Just create matrices like np.zeros([rows, cols]). Same # goes for cdsm and tdsm (vegetation) 

filepath_dsm = './data/DSM_KRbig.asc'
filepath_cdsm = './data/CDSM_KRbig.asc'
usevegdem = 0 
trans = 0.05 # transmissivity of shortwave radiation through vegetation
trunkratio = 0.23 # hieght of trunk zone based on cdsm
altitude = 54.8  # sun altitude in degrees
azimuth = 187.7 # sun azimuth in degrees

gdal_dsm = gdal.Open(filepath_dsm)
dsm = gdal_dsm.ReadAsArray().astype(float)
sizex = dsm.shape[0]  # rows
sizey = dsm.shape[1]  # cols
rows = dsm.shape[0]
cols = dsm.shape[1]
geotransform = gdal_dsm.GetGeoTransform()
scale = 1 / geotransform[1]

gdal_cdsm = gdal.Open(filepath_cdsm)
vegdsm = gdal_cdsm.ReadAsArray().astype(float)
gdal_bh = gdal.Open('./data/heightnew_KRbig.tif')
walls = gdal_bh.ReadAsArray().astype(float)
gdal_ba = gdal.Open('./data/aspectnew_KRbig.tif')
wallaspect = gdal_ba.ReadAsArray().astype(float)

vegdsm2 = vegdsm * trunkratio

# walls = wa.findwalls(dsm, 3)



if usevegdem == 1:
    # Vegetation transmittivity of shortwave radiation
    psi = trans
    # amaxvalue
    vegmax = vegdsm.max()
    amaxvalue = dsm.max() - dsm.min()
    amaxvalue = np.maximum(amaxvalue, vegmax)

    # Elevation vegdsms if buildingDEM includes ground heights
    vegdem = vegdsm + dsm
    vegdem[vegdem == dsm] = 0
    vegdem2 = vegdsm2 + dsm
    vegdem2[vegdem2 == dsm] = 0

    # % Bush separation
    bush = np.logical_not((vegdem2 * vegdem)) * vegdem
else:
    psi = 1.
    bush = np.zeros([rows, cols])
    amaxvalue = 0

#f = sh.shadowingfunctionglobalradiation(dsm, azimuth, altitude, scale) # without wall shadows
t0 = time.time()
sh1, wallsh, wallsun, facesh, facesun = sh.shadowingfunction_wallheight_13(dsm, azimuth, altitude, scale, walls, wallaspect*np.pi/180.)
print('Time: ' + str(time.time() - t0))

#  Vegetation included
# vegsh, sh, _, wallsh, wallsun, wallshve, _, facesun = sh.shadowingfunction_wallheight_23(dsm, vegdem, vegdem2, azimuth,
#                                                 altitude, scale, amaxvalue, bush, walls, wallaspect * np.pi / 180.)

plt.imshow(f)
plt.show()

Hopefully you get the code to work. /Fredrik

[bbrangeo]

filepath_dsm = './data/DSM_KRbig.tif'

[bbrangeo]

from Utilities.SEBESOLWEIGCommonFiles import shadowingfunction_wallheight_13

...

sh1, wallsh, wallsun, facesh, facesun = shadowingfunction_wallheight_13(dsm, azimuth, altitude, scale, walls, wallaspect*np.pi/180.)

[s-martiz]

Thanks @bbrangeo and @biglimp
Still, I can't find this files in the Goteborg example

heightnew_KRbig.tif
aspectnew_KRbig.tif

[biglimp]

They are not included. Create them in Wall and height aspect tool in UMEP.

[s-martiz]

I finally manage to make it works
thank you @biglimp