mdbartos / pysheds

:earth_americas: Simple and fast watershed delineation in python.
GNU General Public License v3.0
706 stars 191 forks source link

catchment crashes when dir points out of array boundaries #123

Closed AlienAtSystem closed 2 years ago

AlienAtSystem commented 4 years ago

The following code run over the attached image crashes via an IndexError. This appears to be the case because when it tries to calculate the catchment for the point on the edge, it tries to recurse into a cell outside of the image boundary.

Tamriel_Heightmap_detail2

from pysheds.grid import Grid
import numpy as np
import matplotlib.pyplot as plt
from affine import Affine

path="Tamriel_Heightmap_detail2.png"

img=plt.imread(path)
if len(img.shape)>2:
    img=img[:,:,0]

grid = Grid()
grid.add_gridded_data(img,'dem',affine=Affine.identity(), crs=grid.crs, nodata=0,mask=np.logical_not(np.isnan(img)))
grid.fill_depressions('dem', out_name='flooded_dem')
grid.resolve_flats('flooded_dem', out_name='inflated_dem')
#           [N, NE, E, SE, S, SW, W, NW]
dirmap = (64,  128,  1,   2,    4,   8,    16,  32)
grid.flowdir(data='inflated_dem', out_name='dir', dirmap=dirmap)
peak=(63,18)
grid.catchment(x=peak[1],y=peak[0],data='dir',xytype='index',out_name='c',dirmap=dirmap,recursionlimit=15000)
plt.imshow(grid.dir,alpha=0.5)
plt.show()
mdbartos commented 2 years ago

Greetings, this should be fixed in pysheds v0.3 (see #165), but needs to be tested for edge cases.