Open barneydobson opened 2 months ago
Example code to do this with @cheginit 's new wrapper:
import shapely
bbox = shapely.box(*dem.rio.bounds())
u, x, y = zip(
*[
(u, float(p["x"]), float(p["y"]))
for u, p in g.nodes(data=True)
if shapely.Point(p["x"], p["y"]).within(bbox)
]
)
gpd.GeoSeries(gpd.points_from_xy(x, y), crs=g.graph["crs"]).to_file("pour_pts.shp")
whitebox_tools("BreachDepressionsLeastCost", ["-i=dem.tif", "-o=dem_corr.tif"])
whitebox_tools("D8Pointer", ["-i=dem_corr.tif", "-o=fdir.tif"])
whitebox_tools("D8FlowAccumulation", ["-i=fdir.tif", "--pntr", "-o=streams.tif"])
whitebox_tools("JensonSnapPourPoints", ["--pour_pts=pour_pts.shp", "--streams=streams.tif", "--snap_dist=15.0", "-o=pour_pts_snapped.shp"])
whitebox_tools("Watershed", ["--d8_pntr=fdir.tif", "--pour_pts=pour_pts_snapped.shp", "-o=watersheds.tif"])
with rasterio.open("watersheds.tif") as src:
gdf_bas = vectorize(
src.read(1).astype(np.int32), 0, src.transform, src.crs, name="basin"
)
gdf_bas = gdf_bas[gdf_bas.basin != src.nodata]
gdf_bas["id"] = [u[x - 1] for x in gdf_bas["basin"]]
A more elegant way of calling whitebox_tools
is using a dict
:
wbt_args = {
"BreachDepressionsLeastCost": ["-i=dem.tif", "-o=dem_corr.tif"],
"D8Pointer": ["-i=dem_corr.tif", "-o=fdir.tif"],
"D8FlowAccumulation": ["-i=fdir.tif", "--pntr", "-o=streams.tif"],
"JensonSnapPourPoints": ["--pour_pts=pour_pts.shp", "--streams=streams.tif", "--snap_dist=15.0", "-o=pour_pts_snapped.shp"],
"Watershed": ["--d8_pntr=fdir.tif", "--pour_pts=pour_pts_snapped.shp", "-o=watersheds.tif"]
}
for tool, args in wbt_args.items():
whitebox_tools(tool, args)
With #266 we have
whitebox
which introduces a range of powerful tools for catchment delineation.A summary of obvious tasks:
pyflwdir
with this in the long run.