Генерация гексагонов внутри полигона (OSM)
Код взять из статьи по ссылке: https://habr.com/ru/post/579838/.
Необходимо разделить область города Москвы на гексагоны. В статье выше аналогичный код для города Краснодара.
При изменении зоны в cities и в polygon_krd на 'Москва' возникают различные ошибки:
1.
m = folium.Map(location=[sum(lats)/len(lats), sum(lons)/len(lons)], zoom_start=13, tiles='cartodbpositron')
overlay = gpd.GeoSeries(geometry).to_json()
folium.GeoJson(overlay, name = 'boundary').add_to(m)
ZeroDivisionError: division by zero
geoJson = {'type':'Polygon','coordinates': [np.column_stack((np.array(geoJson['coordinates'][0])[:, 1], np.array(geoJson['coordinates'][0])
IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed
Скажите, как можно решить такие проблемы? Видела, что можно использовать h3pandas, но не удается найти аналог метода geo_to_h3 в этой библиотеке. Буду признательна!
Сам код:
import geopandas as gpd
import pandas as pd
import numpy as np
import json
from h3 import h3
import h3pandas
import folium
import osmnx as ox
from shapely import wkt
from folium.plugins import HeatMap
from shapely.geometry import Polygon
def visualize_polygons(geometry):
lats, lons = get_lat_lon(geometry)
# в строке ниже появляется ошибка
m = folium.Map(location=[sum(lats)/len(lats), sum(lons)/len(lons)], zoom_start=13, tiles='cartodbpositron')
overlay = gpd.GeoSeries(geometry).to_json()
folium.GeoJson(overlay, name = 'boundary').add_to(m)
return m
# выводим центроиды полигонов
def get_lat_lon(geometry):
lon = geometry.apply(lambda x: x.x if x.type == 'Point' else x.centroid.x)
lat = geometry.apply(lambda x: x.y if x.type == 'Point' else x.centroid.y)
return lat, lon
# выгрузим границы Москвы из OSM, после изменения Краснодар на Москву появляются ошибки
cities = ['Москва']
polygon_krd = ox.geometries_from_place(cities, {'boundary':'administrative'}).reset_index()
polygon_krd = polygon_krd[(polygon_krd['name'] == 'Москва')]
# посмотрим что получилось
visualize_polygons(polygon_krd['geometry'])
Создается функция для генерации гексагонов внутри полигона:
def create_hexagons(geoJson):
polyline = geoJson['coordinates'][0]
polyline.append(polyline[0])
lat = [p[0] for p in polyline]
lng = [p[1] for p in polyline]
m = folium.Map(location=[sum(lat)/len(lat), sum(lng)/len(lng)], zoom_start=13, tiles='cartodbpositron')
my_PolyLine=folium.PolyLine(locations=polyline,weight=8,color="green")
m.add_child(my_PolyLine)
hexagons = list(h3.polyfill(geoJson, 8))
polylines = []
lat = []
lng = []
for hex in hexagons:
polygons = h3.h3_set_to_multi_polygon([hex], geo_json=False)
# flatten polygons into loops.
outlines = [loop for polygon in polygons for loop in polygon]
polyline = [outline + [outline[0]] for outline in outlines][0]
lat.extend(map(lambda v:v[0],polyline))
lng.extend(map(lambda v:v[1],polyline))
polylines.append(polyline)
for polyline in polylines:
my_PolyLine=folium.PolyLine(locations=polyline,weight=3,color='red')
m.add_child(my_PolyLine)
polylines_x = []
for j in range(len(polylines)):
a = np.column_stack((np.array(polylines[j])[:,1],np.array(polylines[j])[:,0])).tolist()
polylines_x.append([(a[i][0], a[i][1]) for i in range(len(a))])
polygons_hex = pd.Series(polylines_x).apply(lambda x: Polygon(x))
return m, polygons_hex, polylines
# polygon_hex , polylines - геометрии гексагонов в разных форматах
# сгенерим гексагоны внутри полигона г. Москвы
geoJson = json.loads(gpd.GeoSeries(polygon_krd['geometry']).to_json())
geoJson = geoJson['features'][0]['geometry']
# в строках ниже появляется ошибка
geoJson = {'type':'Polygon','coordinates': [np.column_stack((np.array(geoJson['coordinates'][0])[:, 1],
а np.array(geoJson['coordinates'][0])[:, 0])).tolist()]}
m, polygons, polylines = create_hexagons(geoJson)
m
Ответы (1 шт):
Можно выбрать любой multypolygon, содержащий город или регион, который вам необходим:
polygon_krd[(polygon_krd['name'] == 'Центральный федеральный округ')]
А далее с помощью explode() можно разбить multipolygon на отдельные polygon(s):
polygon_krd = polygon_krd.explode()
polygon_krd