snlite scripts for geopandas and shapely

[zeffii]

maybe it makes sense to add a folder for scripts specifically for digging through these geometry layers. But for now i'll add them here..

"""
>in geom s
out verts v
out edges s
"""

try:
    import pandas as pd
    import shapely.geometry as sg
    from sverchok.data_structure import get_edge_loop
    from sverchok_gis.utils.shapely_parsing import get_edge_loops_from_multipoly

    for geo in geom[0]:
        coords = geo.__geo_interface__['coordinates']

        if isinstance(geo, sg.Polygon):
            array_2d = np.array(coords[0])
            edge_list = get_edge_loop(array_2d.shape[0])
        elif isinstance(geo, sg.MultiPolygon):
            array_2d = np.array(list(pd.core.common.flatten(coords))).reshape((-1, 2))
            edge_list = get_edge_loops_from_multipoly(coords)
        else:
            continue
        edges.append(edge_list)
        verts.append(np.c_[ array_2d, np.zeros(array_2d.shape[0])])

except Exception as err:
    print(err)

this appears to get all multipolygon edgenets.

[zeffii]

i;d like to point out that using >in geom s here (inside the directive), the > allows me to drop code like if geom and geom[0] (like you will have seen in earlier versions of this script..)

> makes sure the socket is connected and contains data, before it executes any part of the code after the directive.

[zeffii]

now ... i suppose generating ngons for the 2d map is "trivial", and worth starting with

[zeffii]

"""
>in geom s
out verts v
out edges s
out faces s
"""

try:
    import mathutils
    import pandas as pd
    import shapely.geometry as sg
    from sverchok.data_structure import get_edge_loop
    from sverchok_gis.utils.shapely_parsing import get_edge_loops_from_multipoly
    from sverchok_gis.utils.shapely_parsing import get_face_indices_from_multipoly

    for geo in geom[0]:
        coords = geo.__geo_interface__['coordinates']
        indices = []
        if isinstance(geo, sg.Polygon):
            array_2d = np.array(coords[0])
            edge_list = get_edge_loop(array_2d.shape[0])
            indices = mathutils.geometry.tessellate_polygon(coords)

        elif isinstance(geo, sg.MultiPolygon):
            array_2d = np.array(list(pd.core.common.flatten(coords))).reshape((-1, 2))
            edge_list = get_edge_loops_from_multipoly(coords)
            indices = get_face_indices_from_multipoly(coords)
        else:
            continue

        verts.append(np.c_[ array_2d, np.zeros(array_2d.shape[0])])
        edges.append(edge_list)
        faces.append(indices)

except Exception as err:
    print(err)

[zeffii]

• for the Flat geometry this is solved (until otherwise encountered).
• making country MultiPolygons representable in 3d/spherical conditions will involve doing a more complex tessellation
  • will require even subdivision of each MultiPolygon ( possible to do the subdivsion in 2d space using a Division intensity field which is a function of the Projection stretching going on in the 2d coordinate <----- i have no idea how any of that works. yet )

[Marcus-Richmond]

Yea, multipolygon was a beast I tried to avoid for while starting haha.

Here's another topic for a later date - correctly drawing polygons with holes in them. For instance, a (hypothetical) detailed polygon of any country will have holes in it for lake boundaries. I remember looking into this the other week but it seems like blender has to have edges connecting holes in polygons, but GIS data will just have a list of vertices for the outer edge and entire polygon, then another list for each inner hole (what QGIS calls rings). So I'm not sure if blender can draw geometry the same way QGIS can or what. Anyway, thinking aloud here, definitely a topic for another time.

[zeffii]

yeah lakes and rivers complicate things, but a polygon can't have holes. it just means a shape is filled using multiple polygons (the underlying edges of which are not drawn for the user to see) . Blender does have a Curve data type which accepts "rings" of vertices and tends do perform quite well.

[zeffii]

in 2D it is possible to convert a list of rings (or polylines) into a single "Curve" object, and you'll get something like this.

rather than writing code that will be slower in python than the native Blender Curve tesselation code, it's not too painful to write a small converter for shapely.geometry.Polygons / MultiPolygons -> Curve2D -> Mesh

[zeffii]

if you have a geodataframe which includes such ring based shapes. dump it here : ) and we'll experiment

[zeffii]

hahah, 7 years ago I came up with a similar GIS answer: https://blender.stackexchange.com/questions/39975/how-to-draw-a-postgis-polygon-with-holes-into-blender-using-python

somebody came up with another answer, which i think is more suitable now.

scanline script

"""
>in verts v
>in edges s
out verts_out v
out faces_out s
"""
import bmesh

bm = bmesh_from_pydata(verts[0], edges[0], [])
bmesh.ops.triangle_fill(bm, use_beauty=True, use_dissolve=False, edges=bm.edges[:])
new_verts, _, new_faces = pydata_from_bmesh(bm)
verts_out.append(new_verts)
faces_out.append(new_faces)

[zeffii]

i have been poking around in Shapely for the past week, this seems to be an appropriate place maybe to have a few real shapely nodes. Else I will add them directly to Sverchok eventually.

[Marcus-Richmond]

if you have a geodataframe which includes such ring based shapes. dump it here : ) and we'll experiment

sorry, not sure how I missed that a while back, but check this out: https://github.com/Marcus-Richmond/example-datasets/blob/main/Polygon_Rings.zip

This is an example geopackage with two layers: _Canadasimple and _Canadacomplex. The simple layer is a single polygon layer of canada, with a couple holes representing the larger lakes. The complex layer is the same boundary with the same holes/lakes, but has additional polygons representing islands within those lakes. The complex layer is also a multipolygon layer.