Plotting GeoDataFrames with holes

Despite this

and this

and this

and this

I had no clue for how to properly create a GIS plot when my input data is a GeoDataFrame that contains holes.

Thankfully, ChatGPT came to the rescue. It was able to create a function gdf_to_cds for converting GeoDataFrames to ColumnDataSource.

I am dropping this here in the hope that google may show it to other people.
If developmennt time allowed it a function like this should be a part of bokeh itself, in my humble opinion.

from bokeh.models import ColumnDataSource, MultiPolygons
from bokeh.plotting import figure, show
import geopandas as gpd
from shapely.geometry import Polygon, MultiPolygon

# Create a GeoDataFrame with Polygon and MultiPolygon geometries
gdf = gpd.GeoDataFrame({
    'name': ['Polygon 1', 'Polygon 2', 'MultiPolygon'],
    'geometry': [
        Polygon([(0, 0), (1, 0), (1, 1), (0, 1)], holes=[[(0.2, 0.2), (0.2, 0.8), (0.8, 0.8), (0.8, 0.2)]]),
        Polygon([(2, 2), (3, 2), (3, 3), (2, 3)], holes=[[(2.2, 2.2), (2.2, 2.8), (2.8, 2.8), (2.8, 2.2)]]),
        MultiPolygon([
            Polygon([(4, 4), (5, 4), (5, 5), (4, 5)], holes=[[(4.2, 4.2), (4.2, 4.8), (4.8, 4.8), (4.8, 4.2)]]),
            Polygon([(6, 6), (7, 6), (7, 7), (6, 7)], holes=[[(6.2, 6.2), (6.2, 6.8), (6.8, 6.8), (6.8, 6.2)]])
        ])
    ]
})

def gdf_to_cds(gdf):
    from shapely.geometry import MultiPolygon
    gdf['geometry'] = gdf['geometry'].apply(lambda x: MultiPolygon([x]) if x.geom_type == 'Polygon' else x)

    xs = []
    ys = []

    for mp in gdf['geometry']:
        xs.append([[p.exterior.xy[0].tolist(), *[hole.xy[0].tolist() for hole in p.interiors]] for p in mp.geoms])
        ys.append([[p.exterior.xy[1].tolist(), *[hole.xy[1].tolist() for hole in p.interiors]] for p in mp.geoms])

    cds_dict = gdf.drop(columns='geometry').to_dict(orient='list')
    cds_dict['xs'] = xs
    cds_dict['ys'] = ys

    return ColumnDataSource(cds_dict)

# Convert the GeoDataFrame to a ColumnDataSource
cds = gdf_to_cds(gdf)

# Create a Bokeh plot with a MultiPolygons glyph
p = figure(match_aspect=True)
p.multi_polygons(xs='xs', ys='ys', source=cds)

show(p)

Output:

1 Like

FWIW, as someone that heavily uses geopandas–>bokeh for my work, having functions provided to us like this would be nice (i’ve written my own), but we can’t expect the devs to write “translator” functions for whatever library/package we might be coming from, even more so when apparently chatGPT can do it for us :joy:

Also, note that this function only works for geodataframes containing polygons and/or multipolygon geometries. You should ask chatGPT for linestring/multilinestring and point/multipoint and geometry collection translation too…

Another major pitfall to watch out for is projection/crs woes → to get your geospatial data on a figure with a tile server, you’ll need to reproject to epsg:3857… but once you’ve done that, you can’t calculate areas/lengths off the geometry correctly anymore without doing harder math (e.g. converting to lat long and using haversine’s formula for lengths etc.).

1 Like

This topic was automatically closed 90 days after the last reply. New replies are no longer allowed.