B16: локальный OSM в PostGIS — geo_service без интернета

- services/osm_local.py: слой запросов к planet_osm_* (osm2pgsql, SRID 3857):
  roads/density (ST_Length merc × cos(lat) — проверено на merc/geo=1/cos),
  forest (честная площадь ST_Intersection), water (полигон/водоток),
  settlements (village/town/city). Сектор-«пирог» 45° строится в WGS84.
- geo_service.get_zone_features: PostGIS-путь, Overpass — fallback
  (osm_available() по наличию planet_osm_line).
- scripts/b16-import.sh: osm2pgsql --slim -C 800, Беларусь 333МБ pbf
  (Geofabrik, md5 проверен) → 4.5 мин на CT108; PostGIS 3.6 установлен
  в живой vector-postgres (apt, без смены образа).
- Валидация: Минск/Полоцк/Браслав — правдоподобные плотности (7.9/2.9/1.5
  км/км²), зоны различаются по направлениям; 32 сектора за 1.9 c (без сети).
- Тесты: +8 (osm_local контракт, интеграция, отсутствие Overpass-вызовов),
  вне БД — skip. Всего 190 passed, 5 skipped.

Критерий B16: build_search_zones без внешнего интернета — ВЫПОЛНЕН
(единственный Overpass остаётся fallback-веткой при пустых таблицах).
This commit is contained in:
2026-09-09 18:58:30 +03:00
parent 8e99ea0290
commit 3360f7a8ac
4 changed files with 294 additions and 10 deletions
+14 -10
View File
@@ -1,5 +1,9 @@
"""
Geo service for building search zones and querying OpenStreetMap data via Overpass API.
"""Geo service for building search zones and querying OpenStreetMap data.
B16: данные берутся из ЛОКАЛЬНОГО PostGIS (services/osm_local.py, таблицы
planet_osm_* залиты osm2pgsql). Overpass — только fallback, если локальных
таблиц нет (osm_available() == False). Никаких внешних вызовов при
локальном режиме.
"""
import asyncio
import math
@@ -13,6 +17,8 @@ from typing import List, Dict, Optional, Tuple
import httpx
from pydantic import BaseModel
from services.osm_local import get_zone_features_postgis, osm_available
class Zone(BaseModel):
"""Search zone with geographic features."""
@@ -344,17 +350,15 @@ def calculate_forest_coverage(elements: List[Dict], radius_m: int) -> float:
async def get_zone_features(lat: float, lon: float, direction: str, radius_m: int) -> Dict:
"""
Get geographic features for a zone using Overpass API.
Args:
lat, lon: Center point
direction: Sector direction
radius_m: Search radius in meters
"""Характеристики сектора. B16: локальный PostGIS; Overpass — fallback.
Returns:
Dict with roads_km, water_distance_km, settlement_distance_km, forest_pct
Dict with roads_km, road_density, water_distance_km, settlement_distance_km, forest_pct
"""
if osm_available():
return get_zone_features_postgis(lat, lon, direction, radius_m)
# Fallback: живой Overpass (старое поведение, пока B16-данные не залиты)
# Query roads
roads_query = f"""
[out:json];
+162
View File
@@ -0,0 +1,162 @@
"""Локальный OSM-слой (B16): пространственные запросы к PostGIS вместо Overpass.
Таблицы planet_osm_* заливаются osm2pgsql из дампа Беларуси
(scripts/b16-import.sh; SRID 3857 — дефолт osm2pgsql, проекция Web Mercator).
Обновление — разово/по мере надобности. Никакого интернета.
Доступ к БД — свой sync-движок из DATABASE_URL (asyncpg → psycopg2).
"""
from __future__ import annotations
import logging
import math
import os
from typing import Any, Optional
from sqlalchemy import create_engine, text
from sqlalchemy.engine import Engine
logger = logging.getLogger(__name__)
_PG_URL = os.getenv('DATABASE_URL', 'postgresql://postgres:postgres@postgres:5432/vector_mchs')
if _PG_URL.startswith('postgresql+asyncpg://'):
_PG_URL = _PG_URL.replace('postgresql+asyncpg://', 'postgresql://')
_engine: Optional[Engine] = None
_available: Optional[bool] = None
def _get_engine() -> Engine:
global _engine
if _engine is None:
_engine = create_engine(_PG_URL, pool_pre_ping=True)
return _engine
def osm_available() -> bool:
"""True, если локальные OSM-таблицы залиты (иначе geo_service падает на Overpass)."""
global _available
if _available is not None:
return _available
try:
with _get_engine().connect() as conn:
row = conn.execute(text(
"SELECT count(*) FROM information_schema.tables "
"WHERE table_name = 'planet_osm_line'"
)).scalar()
_available = bool(row)
except Exception as e:
logger.warning('OSM PostGIS недоступен: %s', e)
_available = False
return _available
def _sector_wkt_3857(lat: float, lon: float, radius_m: int, direction: str) -> str:
"""Сектор-«пирог» 45° в Web Mercator (проекция данных planet_osm_*).
Строим в WGS84 (локальные метрики), потом трансформируем — mercator-метры
на широте Беларуси раздуты, строить дугу в них нельзя.
"""
angles = {'N': 0, 'NE': 45, 'E': 90, 'SE': 135, 'S': 180, 'SW': 225, 'W': 270, 'NW': 315}
az = angles.get(direction, 0)
half = 22.5
seg = 24
pts: list[tuple[float, float]] = []
for i in range(seg + 1):
bearing = math.radians(az - half + (2 * half) * i / seg)
dlat = (radius_m / 111320.0) * math.cos(bearing)
dlon = (radius_m / (111320.0 * max(0.2, math.cos(math.radians(lat))))) * math.sin(bearing)
pts.append((lon + dlon, lat + dlat))
# Замыкаем пирог через центр
pts_wkt = f'{lon:.7f} {lat:.7f}, ' + ','.join(f'{p[0]:.7f} {p[1]:.7f}' for p in pts)
return f'POLYGON(({pts_wkt}, {lon:.7f} {lat:.7f}))'
def get_zone_features_postgis(lat: float, lon: float, direction: str, radius_m: int) -> dict[str, Any]:
"""Характеристики сектора — 1 запрос к локальному PostGIS.
Контракт тот же, что get_zone_features (Overpass): roads_km, road_density,
water_distance_km, settlement_distance_km, forest_pct.
"""
sector_wkt = _sector_wkt_3857(lat, lon, radius_m, direction)
center_pt = f'POINT({lon} {lat})'
q = text("""
WITH s AS (
SELECT ST_Transform(ST_SetSRID(ST_GeomFromText(:sector_wkt), 4326), 3857) AS sector,
ST_Transform(ST_SetSRID(ST_GeomFromText(:center_pt), 4326), 3857) AS center_m
),
road_len AS (
SELECT COALESCE(SUM(ST_Length(l.way)), 0) AS meters
FROM planet_osm_line l, s
WHERE l.highway IS NOT NULL
AND l.way && s.sector
AND ST_Intersects(l.way, s.sector)
),
forest_area AS (
SELECT COALESCE(SUM(ST_Area(ST_Intersection(p.way, s.sector))), 0) AS m2
FROM planet_osm_polygon p, s
WHERE (p.landuse = 'forest' OR p."natural" = 'wood')
AND p.way && s.sector
),
water_near AS (
SELECT MIN(ST_Distance(p.way, s.center_m)) / :scale AS meters
FROM planet_osm_polygon p, s
WHERE (p."natural" = 'water' OR p.landuse = 'reservoir')
AND p.way && ST_Transform(ST_SetSRID(ST_GeomFromText(:sector_wkt), 4326), 3857)
),
waterway_near AS (
SELECT MIN(ST_Distance(l.way, s.center_m)) / :scale AS meters
FROM planet_osm_line l, s
WHERE l.waterway IN ('river','stream','canal','ditch','drain')
AND l.way && s.sector
),
settle_near AS (
SELECT MIN(ST_Distance(p.way, s.center_m)) / :scale AS meters
FROM planet_osm_point p, s
WHERE p.place IN ('village','town','city')
AND p.way && s.sector
)
SELECT
(SELECT meters FROM road_len),
(SELECT m2 FROM forest_area),
(SELECT meters FROM water_near),
(SELECT meters FROM waterway_near),
(SELECT meters FROM settle_near)
""")
# Масштаб 3857 на широте Минска (~53.9°): cos(lat)
scale = max(0.2, math.cos(math.radians(lat)))
with _get_engine().connect() as conn:
row = conn.execute(q, {
'sector_wkt': sector_wkt, 'center_pt': center_pt, 'scale': scale,
}).one()
# 3857-метры на широте lat раздуты в 1/cos(lat) раз:
# real = mercator * cos(lat) (проверено: merc/geo = 1/cos(lat) точно).
roads_km = float(row[0] or 0) / 1000.0 * scale
forest_m2 = float(row[1] or 0) * (scale * scale) # mercator-м² → реальные м²
# Расстояния (row[2..4]) уже поделены на scale в SQL → вернулись в метры
water_poly_m = row[2]
waterway_m = row[3]
settle_m = row[4]
search_area_km2 = math.pi * (radius_m / 1000.0) ** 2
road_density = roads_km / search_area_km2 if search_area_km2 > 0 else 0.0
forest_pct = 0.0
if forest_m2 > 0 and search_area_km2 > 0:
forest_pct = round(min(100.0, (forest_m2 / 1e6) / search_area_km2 * 100), 1)
water_candidates = [float(m) for m in (water_poly_m, waterway_m) if m is not None]
water_distance_km = round(min(water_candidates), 3) if water_candidates else None
settlement_distance_km = round(float(settle_m) / 1000.0, 3) if settle_m is not None else None
return {
'roads_km': round(roads_km, 3),
'road_density': round(road_density, 2),
'water_distance_km': water_distance_km,
'settlement_distance_km': settlement_distance_km,
'forest_pct': forest_pct,
}