Как удалить океан с помощью OpenStreetMap? ⇐ Python

[Anonymous]

Я работаю с данными библиотеки OpenStreetMap (osmnx, в Python), чтобы извлечь границу Сальвадора, Бразилия. Код успешно удаляет водоемы и получает материк и острова (я их опускаю для простоты), но есть проблема с искусственной треугольной границей, которая простирается в океан.
Этот треугольник образован частью административной границы, которая почему-то сливается с естественной береговой линией и не вычитается океаном. Как бы вы попытались удалить эту искусственную треугольную границу в OpenStreetMap? Есть ли какой-нибудь другой фильтр, который мне следует попробовать? Я хочу найти более общее решение, чем делать «дебафф, бафф» вручную. Я пробовал использовать береговую линию, но результаты оказались не очень хорошими.

Код: Выделить всё

import osmnx as ox
import shapely

def get_city_boundary(place_name):
city_gdf = ox.geocode_to_gdf(place_name)
city_polygon = city_gdf.loc[0, "geometry"]

expand = 0
minx, miny, maxx, maxy = city_polygon.bounds
bounding_box = shapely.geometry.box(minx - expand, miny - expand, maxx + expand, maxy + expand)

try:
water_tags = {"natural": ["bay"], "maritime": "yes"}
water_gdf = ox.features_from_polygon(bounding_box, tags=water_tags)
if not water_gdf.empty:
water_union = water_gdf.geometry.union_all()
city_polygon = city_polygon.difference(water_union)
except Exception as e:
print(f"Error processing water areas: {e}")

city_polygon = city_polygon.buffer(0)
### change around here a bit if you want the islands of the city to appear
if isinstance(city_polygon, shapely.geometry.MultiPolygon):
polygons = list(city_polygon.geoms)
city_polygon = max(polygons, key=lambda p: p.area)

return city_polygon

place_name = "Salvador, Bahia, Brazil"
processed_polygon = get_city_boundary(place_name)

processed_gdf = ox.geocode_to_gdf(place_name).copy()
processed_gdf.loc[0, "geometry"] = processed_polygon

if processed_gdf.crs is None:
processed_gdf.set_crs(epsg=4326, inplace=True)

output_file = "salvador_processed_boundary.shp"
processed_gdf.to_file(output_file)

Подробнее здесь: https://stackoverflow.com/questions/793 ... nstreetmap