pythonopenstreetmapstreamlitfoliumh3

Number of POIs per H3 hexagon for cities in streamlit app - hexagons not showing


I create a streamlit application, which is supposed to let the user select an area and then calculate the number of POIs per H3 Hexagon unit in this area. Somehow the hexagons are created but do not show up on the map. Here is the code:

import folium
from streamlit_folium import st_folium
import requests
import h3

# Set layout
st.set_page_config(layout="wide", page_title="POIs per Hexagon")

# Define colors for each amenity (example)
amenity_colors = {
    "restaurant": "blue",
    "bank": "green",
    "post_box": "red",
    "atm": "purple",
    "post_office": "orange",
    "marketplace": "pink",
    "pharmacy": "brown",
    "hospital": "gray",
    "clinic": "yellow",
    "doctor": "cyan",
    "dentist": "olive",
    "kindergarten": "lime",
    "childcare": "teal",
    "school": "magenta",
    "fast_food": "navy",
    "cafe": "maroon",
    "bar": "gold",
    "pub": "indigo",
    "ice_cream": "tan",
    "nightclub": "coral",
    "biergarten": "lavender",
    "community_centre": "slateblue",
    "library": "orchid",
    "theatre": "skyblue",
    "cinema": "seagreen",
    "arts_centre": "khaki",
    "events_venue": "darkred"
}

# Main Content
st.title("Calculate the POIs per Hexagon")
st.write("blablabla explaining text")

st.header("Choose a City")
st.write("Choose an area on the map below. The visible area of the map will be used to create the count.")

lat, lon = 49.40768, 8.69079  # Default coordinates

# Initialize session state for map markers and bounds
if "map_initialized" not in st.session_state:
    st.session_state.map_initialized = False
if "map_bounds" not in st.session_state:
    st.session_state.map_bounds = None
if "map_markers" not in st.session_state:
    st.session_state.map_markers = []


def get_center_from_bounds(bounds):
    """Calculate the center of the map from its bounds."""
    center_lat = (bounds['_southWest']['lat'] + bounds['_northEast']['lat']) / 2
    center_lon = (bounds['_southWest']['lng'] + bounds['_northEast']['lng']) / 2
    return center_lat, center_lon


def fetch_pois(sw_lat, sw_lng, ne_lat, ne_lng):
    """Fetch POIs (Points of Interest) from OSM using the Overpass API."""
    overpass_url = "http://overpass-api.de/api/interpreter"
    amenities = [
        "restaurant", "bank", "post_box", "atm", "post_office", "marketplace",
        "pharmacy", "hospital", "clinic", "doctor", "dentist", "kindergarten",
        "childcare", "school", "fast_food", "cafe", "bar", "pub", "ice_cream",
        "nightclub", "biergarten", "community_centre", "library", "theatre",
        "cinema", "arts_centre", "events_venue"
    ]
    
    nodes = "\n".join(
        f'node["amenity"="{amenity}"]({sw_lat},{sw_lng},{ne_lat},{ne_lng});' 
        for amenity in amenities
    )
    
    overpass_query = f"""
    [out:json];
    (
    {nodes}
    );
    out body;
    """
    
    response = requests.get(overpass_url, params={"data": overpass_query})
    return response


def add_markers_to_map(m, markers):
    """Add markers to the map with custom colors based on amenity type."""
    for marker in markers:
        amenity_type = marker['tooltip'].lower()  # assuming tooltip contains amenity type
        color = amenity_colors.get(amenity_type, 'gray')  # default to gray if type not found
        folium.Marker(location=marker['location'], popup=marker['popup'], tooltip=marker['tooltip'],
                      icon=folium.Icon(color=color)).add_to(m)


def display_map_with_bounds(lat, lon, zoom_start=12):
    """Display the map with given latitude, longitude, and zoom level."""
    m = folium.Map(location=[lat, lon], zoom_start=zoom_start)
    bounds_script = """
    <script>
    function updateBounds() {
        const bounds = map.getBounds();
        const data = {
            _northEast: bounds.getNorthEast(),
            _southWest: bounds.getSouthWest()
        };
        Streamlit.setComponentValue(data);
    }
    map.on('moveend', updateBounds);
    updateBounds();
    </script>
    """
    m.get_root().html.add_child(folium.Element(bounds_script))
    add_markers_to_map(m, st.session_state.map_markers)

    # Add legend
    legend_html = """
    <div style="position: fixed; bottom: 50px; left: 50px; z-index:1000; font-size:14px; background:white; padding: 10px; border: 2px solid black;">
      <p><strong>Legend</strong></p>
      """
    for amenity, color in amenity_colors.items():
        legend_html += f'<p><i class="fa fa-circle" style="color:{color}"></i> {amenity}</p>'
    legend_html += "</div>"

    m.get_root().html.add_child(folium.Element(legend_html))

    return st_folium(m, width=700, height=500)


def create_hexagons(pois, resolution):
    hex_counts = {}
    for poi in pois:
        hex_addr = h3.geo_to_h3(poi['location'][0], poi['location'][1], resolution)
        if hex_addr in hex_counts:
            hex_counts[hex_addr] += 1
        else:
            hex_counts[hex_addr] = 1
    return hex_counts


def add_hexagons_to_map(m, hex_counts, resolution):
    for hex_addr, count in hex_counts.items():
        hex_boundary = h3.h3_to_geo_boundary(hex_addr, geo_json=True)
        folium.Polygon(locations=hex_boundary, color='black', fill=True, fill_opacity=0.9,
                       tooltip=f'POIs: {count}').add_to(m)


# Display initial map or map with markers
if st.session_state.map_bounds:
    center_lat, center_lon = get_center_from_bounds(st.session_state.map_bounds)
    output = display_map_with_bounds(center_lat, center_lon)
else:
    output = display_map_with_bounds(lat, lon)

# Handle button click for "Create Score"
if st.button("Create Score"):
    if 'bounds' in output:
        bounds = output['bounds']
        st.session_state.map_bounds = bounds  # Save bounds in session state

        # create variables for coordinates
        if '_southWest' in bounds and '_northEast' in bounds:
            sw_lat = bounds['_southWest']['lat']
            sw_lng = bounds['_southWest']['lng']
            ne_lat = bounds['_northEast']['lat']
            ne_lng = bounds['_northEast']['lng']

            # fetch pois and store to variable
            response = fetch_pois(sw_lat, sw_lng, ne_lat, ne_lng)

            # if response is correct, then show in map as markers
            if response.status_code == 200:
                try:
                    data = response.json()
                    st.session_state.map_markers.clear()
                    pois = []
                    for element in data['elements']:
                        if element['type'] == 'node':
                            marker = {
                                'location': [element['lat'], element['lon']],
                                'popup': element.get('tags', {}).get('name', 'Unnamed'),
                                'tooltip': element.get('tags', {}).get('amenity', 'Restaurant')
                            }
                            pois.append(marker)
                            st.session_state.map_markers.append(marker)

                    # create hexagons
                    resolution = 8  # Adjust resolution as needed
                    hex_counts = create_hexagons(pois, resolution)

                    # Debug: Output hex_counts to verify
                    st.write(f"Hex Counts: {hex_counts}")

                    # get center coordinates of bbox
                    center_lat, center_lon = get_center_from_bounds(st.session_state.map_bounds)
                    # display map according to bbox
                    m = folium.Map(location=[center_lat, center_lon], zoom_start=12)
                    add_markers_to_map(m, st.session_state.map_markers)
                    add_hexagons_to_map(m, hex_counts, resolution)
                    st_folium(m, width=700, height=500)
                except ValueError:
                    st.error("Error decoding JSON response from Overpass API")
                    st.text(response.text)
            else:
                st.error(f"Error fetching data from Overpass API: {response.status_code}")
                st.text(response.text)
        else:
            st.write("Bounding box keys are missing.")
    else:
        st.write("No bounding box information available.")

For a few seconds, while generating, I can see that the hexagons exist with this: st.write(f"Hex Counts: {hex_counts}"). But they are not drawn on the map. This is the result:

enter image description here

Thankful for any hints!


Solution

  • There is a bunch of issues in your code but the main one (the invisible hexagons) is that for some reason h3.h3_to_geo_boundary is returning lon/lat coordinates instead of lat/lon (which folium expects) :

    def add_hexagons_to_map(hex_counts, resolution):
        for hex_addr, count in hex_counts.items():
            hex_boundary = h3.h3_to_geo_boundary(hex_addr, geo_json=True)
            poly = folium.Polygon(
                locations=[lonlat[::-1] for lonlat in hex_boundary], # << flip the lonlat
                color="black",
                fill=True,
                fill_opacity=0.3,
                tooltip=f"POIs: {count}",
            )
    

    Output (streamlit run .\app.py):