mdbartos / pysheds

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

How to exclude NoData values when calculating the flow direction of DEM? #256

Open mht2953658596 opened 2 months ago

mht2953658596 commented 2 months ago

When using PyShed to calculate the flow direction of DEM, I found that the flow direction was also calculated for areas where DEM values were NoData. How can I avoid this situation?

When I set NoData using the following code:

nodata = gdal.Open(raster_path)tGetRasterBand(1).ReadAsArray().GetNoDataValue()
grid = Grid.from_raster(raster_path, nodata=nodata)
dem = grid.read_raster(self.raster_path, band=1, nodata=grid.nodata)
flooded_dem = grid.fill_depressions(dem)
inflated_dem = grid.resolve_flats(flooded_dem)
fdir = grid.flowdir(inflated_dem)

it reports the following error:

image

mht2953658596 commented 2 months ago

If I don't set the NoData value, it will report the same error.

mht2953658596 commented 2 months ago

@cheginit

phillipaarestrup commented 1 month ago

I have found the same issue and also when loading some rasters. I think some of this is related to a new update in the Numpy package from 1.2x -> 2.0.0 (released in june). I only have these issues with the new version of Numpy and not when I downgrade it to 1.2x.

Your screenshow shows that you get a typeerror in sview.py which is a result of this function call np.can_cast(viewfinder.nodata, obj.dtype, casting='safe') in line 85

See also the discussion here on fx the function can_cast, which is used for checks in the pysheds package: https://github.com/numpy/numpy/issues/24856

mht2953658596 commented 1 month ago

@phillipaarestrup But my NumPy version is 1.24.3, and I am still encountering this issue.

JaweedNazary commented 1 month ago

I tried NumPy 2.0 and the nodata representation error pops up when running flodir. However, with NumPy 1.26, it is working pretty good,

JaweedNazary commented 1 month ago

@mht2953658596 try converting your raster to float before using the code.