mdbartos / pysheds

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

Catchment contour to postgis #203

Open neocaos opened 1 year ago

neocaos commented 1 year ago

I'm pretty sure there is no exactly a pysheds issue, but:

I'm trying to infer the watersheds basins from a given point, using your library. Then get the watershed polygon and store into PostGIS database.

From your Readme.md @mdbartos

I've reached this point without problem:

catch = grid.catchment(x=x_snap, y=y_snap, fdir=fdir, dirmap=dirmap, 
                       xytype='coordinate')

as you probably know, it returns a NumPy ndim array

How can extract the basin contour into a poligyon? I've tested parsing the arrays into ogr.ring, then into ogr.polygon, but doesn't seems the proprety way.

Yes, I've close the polygon by inserting the first point at the end, but PostGis tells me st_isValid == false

Could you bring me some help,

maxschmi commented 1 year ago

Hello, I had the same issue. My solution was to use contourpy to create the catchment area as Shapely-Polygon:

import contourpy
from shapely.geometry import Polygon

cont_gen = contourpy.contour_generator(
    z=np.array(catch.astype(int))
)

Polygon([catch.affine * xy + np.array((catch.affine.a/2, catch.affine.e/2)) 
               for xy in cont_gen.lines(0.5)[0] ])

As contourpy only works with numeric values, I first convert the catch array to integer values. Then I draw the contour lines for the value 0.5, so the result will draw a line in the middle between all the 0 and 1 raster cells.