BUG in `workflows.gis.ground_elevation_from_dem`

[DirkEilander]

bug

The workflows.gis.ground_elevation_from_dem assumes only a single building can exist within a grid cell and uses a random filling of missing data values. Note that this is really problematic for impact assessments as the building elevation makes a huge!! difference on the results.

First all building ids are rasterized into a single grid at the resolution of the DEM, thus only one building can exist per grid cell.

hydromt_fiat/hydromt_fiat/workflows/gis.py

Line 298 in 8d37a58

rasterized = rasterize(

Then, and this is really problematic, a bffil is used to fill in gaps of buildings which don't have data because of the previous assumption. This uses data from random (i.e. the order of buildings is random!) other buildings.

hydromt_fiat/hydromt_fiat/workflows/gis.py

Line 323 in 8d37a58

exposure_db["ground_elevtn"].bfill(inplace=True)

possible solutions

• use another filling algorithm by e.g. sampling from the rasterized building id map at the centroid of buildings with missing values.
• use hydromt.raster.zonal_stats. much slower but it does the job right. There might also be other packages (a quick google search brought me to exactextract) that do this correctly.
• use hydromt.raster.sample to get fast results based on the centroid.

[panosatha]

@DirkEilander I have noticed this last week as well. I have made an alternative method that is using either the standar .sel method from xarray to sample at point locations, or the rasterstats package (https://pythonhosted.org/rasterstats/rasterstats.html#submodules) to get the average values based on zonal stats. I played around with the exactextract package you mentioned and speed is similar (but probably results more accurate in some case) so I could replace the rasterstats with this.

I can try to make a PR next week.

[DirkEilander]

@panosatha note that I already suggested a fix in #505