Files
vector/services/geo_service.py
T
Sadmin 98e8e58023 feat: водоёмы в зоне поиска на карте (рекомендации профилей → геослой)
- GET /api/v1/water/{case_id}: полигоны водоёмов/болот OSM вокруг ТНП
  (Overpass, natural=water+wetland, кэш 72ч, circuit breaker от geo_service);
  бейдж-предупреждение по профилям случая (РАС/эпилепсия/не умеет плавать)
- geo_service: цепочка зеркал Overpass (maps.mail.ru → overpass-api.de):
  основной сервер имеет AAAA, на хостах без IPv6-маршрута httpx падал
  ConnectError'ом до IPv4-фолбэка; на CT108 фапало стабильно
- SearchMap: слой водоёмов под зонами (синие полигоны, болота пунктиром),
  попапы с площадью/дистанцией от ТНП; AnalysisResult грузит слой,
  показывает предупреждение профиля и счётчик водоёмов
- проверено на CT108: РАС-кейс (Минск) — 43 объекта, витя-кейс
  (Каменец-борисовская обл.) — 7 объектов, включая водохранилище Загацце
2026-09-09 12:48:44 +03:00

459 lines
14 KiB
Python

"""
Geo service for building search zones and querying OpenStreetMap data via Overpass API.
"""
import asyncio
import math
import json
import hashlib
import logging
import time
from datetime import datetime, timedelta
from pathlib import Path
from typing import List, Dict, Optional, Tuple
import httpx
from pydantic import BaseModel
class Zone(BaseModel):
"""Search zone with geographic features."""
direction: str # N, NE, E, SE, S, SW, W, NW
distance_km: float
forest_pct: float
road_density: float # km of roads per km²
water_distance_km: Optional[float]
settlement_distance_km: Optional[float]
# Cache configuration
CACHE_DIR = Path("/tmp/overpass_cache")
CACHE_TTL_HOURS = 24
# Порядок зеркал: основной сервер имеет AAAA, но на хостах без IPv6-маршрута
# (наш LXC в VPN-контуре) httpx падает на ConnectError ещё до IPv4-фолбэка.
# Поэтому основной — зеркало без AAAA (maps.mail.ru), overpass-api.de — резерв.
OVERPASS_URLS = [
"https://maps.mail.ru/osm/tools/overpass/api/interpreter",
"https://overpass-api.de/api/interpreter",
]
OVERPASS_URL = OVERPASS_URLS[0] # обратная совместимость для импортов
logger = logging.getLogger(__name__)
# --- Overpass circuit breaker ---------------------------------------------
# build_search_zones issues ~128 Overpass calls per analysis. When the host has
# no outbound connectivity every one of them burns the full connect timeout,
# which turns a single analysis into several minutes of waiting for results
# that are empty anyway. After OVERPASS_FAILURE_THRESHOLD consecutive
# transport failures we stop calling out until OVERPASS_COOLDOWN_SECONDS have
# passed. Callers get the same {'elements': []} they already got on error.
OVERPASS_FAILURE_THRESHOLD = 3
OVERPASS_COOLDOWN_SECONDS = 60.0
OVERPASS_CONNECT_TIMEOUT = 5.0
_overpass_failures = 0
_overpass_open_until = 0.0
def _overpass_circuit_open() -> bool:
"""True while the breaker is tripped (skip network, return empty fast)."""
if _overpass_failures < OVERPASS_FAILURE_THRESHOLD:
return False
if time.monotonic() >= _overpass_open_until:
_reset_overpass_circuit()
return False
return True
def _record_overpass_failure() -> None:
global _overpass_failures, _overpass_open_until
_overpass_failures += 1
if _overpass_failures == OVERPASS_FAILURE_THRESHOLD:
_overpass_open_until = time.monotonic() + OVERPASS_COOLDOWN_SECONDS
logger.warning(
"Overpass API недоступен (%d подряд неудачных запросов). "
"Геоданные отключены на %.0f c, анализ продолжается без них.",
_overpass_failures,
OVERPASS_COOLDOWN_SECONDS,
)
def _reset_overpass_circuit() -> None:
global _overpass_failures, _overpass_open_until
_overpass_failures = 0
_overpass_open_until = 0.0
# Direction mappings
SEARCH_DISTANCES = [500, 1000, 2000, 5000]
DIRECTIONS = ["N", "NE", "E", "SE", "S", "SW", "W", "NW"]
DIRECTION_ANGLES = {
"N": 0,
"NE": 45,
"E": 90,
"SE": 135,
"S": 180,
"SW": 225,
"W": 270,
"NW": 315
}
# Search distances computed dynamically from max_distance_km
def _build_search_distances(max_distance_km: float) -> list:
max_m = int(max_distance_km * 1000)
raw = [int(max_m * 0.25), int(max_m * 0.50), int(max_m * 0.75), max_m]
clamped = [max(200, min(5000, d)) for d in raw]
return sorted(set(clamped))
def haversine(lat1: float, lon1: float, lat2: float, lon2: float) -> float:
"""
Calculate distance between two points on Earth using Haversine formula.
Args:
lat1, lon1: First point coordinates
lat2, lon2: Second point coordinates
Returns:
Distance in kilometers
"""
R = 6371 # Earth radius in km
lat1_rad = math.radians(lat1)
lat2_rad = math.radians(lat2)
dlat = math.radians(lat2 - lat1)
dlon = math.radians(lon2 - lon1)
a = (math.sin(dlat / 2) ** 2 +
math.cos(lat1_rad) * math.cos(lat2_rad) * math.sin(dlon / 2) ** 2)
c = 2 * math.asin(math.sqrt(a))
return R * c
def get_sector_bounds(lat: float, lon: float, direction: str, radius_m: int) -> Tuple[float, float, float, float]:
"""
Calculate bounding box for a sector.
Args:
lat, lon: Center point
direction: Sector direction (N, NE, E, etc.)
radius_m: Radius in meters
Returns:
(min_lat, min_lon, max_lat, max_lon)
"""
# Convert radius to degrees (approximate)
radius_deg = radius_m / 111320 # 1 degree ≈ 111.32 km at equator
angle = DIRECTION_ANGLES[direction]
angle_rad = math.radians(angle)
# Calculate sector boundaries (45° sectors)
angle_start = angle - 22.5
angle_end = angle + 22.5
# Simple bounding box (can be optimized for actual sector shape)
lat_offset = radius_deg * math.cos(angle_rad)
lon_offset = radius_deg * math.sin(angle_rad) / math.cos(math.radians(lat))
min_lat = min(lat, lat + lat_offset) - radius_deg * 0.5
max_lat = max(lat, lat + lat_offset) + radius_deg * 0.5
min_lon = min(lon, lon + lon_offset) - radius_deg * 0.5
max_lon = max(lon, lon + lon_offset) + radius_deg * 0.5
return (min_lat, min_lon, max_lat, max_lon)
def get_cache_key(query: str) -> str:
"""Generate cache key from query."""
return hashlib.md5(query.encode()).hexdigest()
def get_cached_result(cache_key: str) -> Optional[Dict]:
"""Get cached Overpass API result if not expired."""
CACHE_DIR.mkdir(exist_ok=True)
cache_file = CACHE_DIR / f"{cache_key}.json"
if not cache_file.exists():
return None
try:
with open(cache_file, 'r') as f:
cached = json.load(f)
cached_time = datetime.fromisoformat(cached['timestamp'])
if datetime.now() - cached_time > timedelta(hours=CACHE_TTL_HOURS):
cache_file.unlink()
return None
return cached['data']
except Exception:
return None
def save_to_cache(cache_key: str, data: Dict):
"""Save Overpass API result to cache."""
CACHE_DIR.mkdir(exist_ok=True)
cache_file = CACHE_DIR / f"{cache_key}.json"
try:
with open(cache_file, 'w') as f:
json.dump({
'timestamp': datetime.now().isoformat(),
'data': data
}, f)
except Exception:
pass
async def query_overpass(query: str) -> Dict:
"""
Query Overpass API with caching.
Args:
query: Overpass QL query
Returns:
API response as dict
"""
cache_key = get_cache_key(query)
# Check cache
cached = get_cached_result(cache_key)
if cached is not None:
return cached
# Skip the network entirely while the breaker is tripped.
if _overpass_circuit_open():
return {'elements': []}
# Query API — по цепочке зеркал (первый успешный ответ кэшируется)
timeout = httpx.Timeout(30.0, connect=OVERPASS_CONNECT_TIMEOUT)
async with httpx.AsyncClient(timeout=timeout) as client:
for overpass_url in OVERPASS_URLS:
try:
response = await client.post(
overpass_url,
data={'data': query},
headers={'Content-Type': 'application/x-www-form-urlencoded'}
)
response.raise_for_status()
data = response.json()
# Save to cache
save_to_cache(cache_key, data)
_reset_overpass_circuit()
return data
except Exception as e:
logger.info("Overpass %s недоступен: %s", overpass_url, e)
continue
# Все зеркала недоступны
_record_overpass_failure()
# Return empty result on error
return {'elements': []}
def calculate_road_length(elements: List[Dict]) -> float:
"""
Calculate total road length from Overpass way elements.
Args:
elements: List of way elements from Overpass
Returns:
Total length in kilometers
"""
total_length = 0.0
for element in elements:
if element.get('type') != 'way':
continue
nodes = element.get('geometry', [])
if len(nodes) < 2:
continue
# Calculate length by summing distances between consecutive nodes
for i in range(len(nodes) - 1):
lat1, lon1 = nodes[i]['lat'], nodes[i]['lon']
lat2, lon2 = nodes[i + 1]['lat'], nodes[i + 1]['lon']
total_length += haversine(lat1, lon1, lat2, lon2)
return total_length
def find_nearest_distance(lat: float, lon: float, elements: List[Dict]) -> Optional[float]:
"""
Find distance to nearest element.
Args:
lat, lon: Reference point
elements: List of node elements from Overpass
Returns:
Distance in kilometers, or None if no elements
"""
if not elements:
return None
min_distance = float('inf')
for element in elements:
if element.get('type') != 'node':
continue
elem_lat = element.get('lat')
elem_lon = element.get('lon')
if elem_lat is None or elem_lon is None:
continue
distance = haversine(lat, lon, elem_lat, elem_lon)
min_distance = min(min_distance, distance)
return min_distance if min_distance != float('inf') else None
def calculate_forest_coverage(elements: List[Dict], radius_m: int) -> float:
"""
Estimate forest coverage percentage.
Args:
elements: List of way elements from Overpass
radius_m: Search radius in meters
Returns:
Forest coverage as percentage (0-100)
"""
if not elements:
return 0.0
# Approximate: count forest ways and estimate coverage
# This is a simplified calculation
forest_ways = len([e for e in elements if e.get('type') == 'way'])
# Rough heuristic: each forest way covers ~0.1 km²
# Total search area = π * r²
search_area_km2 = math.pi * (radius_m / 1000) ** 2
estimated_forest_km2 = forest_ways * 0.1
coverage_pct = min(100.0, (estimated_forest_km2 / search_area_km2) * 100)
return round(coverage_pct, 1)
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
Returns:
Dict with roads_km, water_distance_km, settlement_distance_km, forest_pct
"""
# Query roads
roads_query = f"""
[out:json];
(
way[highway](around:{radius_m},{lat},{lon});
);
out geom;
"""
# Query water bodies
water_query = f"""
[out:json];
(
node[natural=water](around:{radius_m},{lat},{lon});
way[natural=water](around:{radius_m},{lat},{lon});
);
out center;
"""
# Query settlements
settlement_query = f"""
[out:json];
(
node[place~"village|town|city"](around:{radius_m},{lat},{lon});
);
out;
"""
# Query forests
forest_query = f"""
[out:json];
(
way[landuse=forest](around:{radius_m},{lat},{lon});
way[natural=wood](around:{radius_m},{lat},{lon});
);
out geom;
"""
# The four queries are independent - issue them concurrently.
roads_data, water_data, settlement_data, forest_data = await asyncio.gather(
query_overpass(roads_query),
query_overpass(water_query),
query_overpass(settlement_query),
query_overpass(forest_query),
)
roads_km = calculate_road_length(roads_data.get('elements', []))
water_distance = find_nearest_distance(lat, lon, water_data.get('elements', []))
settlement_distance = find_nearest_distance(lat, lon, settlement_data.get('elements', []))
forest_pct = calculate_forest_coverage(forest_data.get('elements', []), radius_m)
# Calculate road density (km of roads per km²)
search_area_km2 = math.pi * (radius_m / 1000) ** 2
road_density = roads_km / search_area_km2 if search_area_km2 > 0 else 0.0
return {
'roads_km': roads_km,
'road_density': round(road_density, 2),
'water_distance_km': water_distance,
'settlement_distance_km': settlement_distance,
'forest_pct': forest_pct
}
async def build_search_zones(lat: float, lon: float, case_data: dict, max_distance_km: float = 3.0) -> List[Zone]:
"""
Build search zones around a point.
Creates 8 directional sectors (N, NE, E, SE, S, SW, W, NW) at multiple distances
(500m, 1000m, 2000m, 5000m) and queries geographic features for each.
Args:
lat: Latitude of search origin
lon: Longitude of search origin
case_data: Case information (for future enhancements)
Returns:
List of Zone objects with geographic features
"""
zones = []
for distance_m in _build_search_distances(max_distance_km):
for direction in DIRECTIONS:
# Get features for this zone
features = await get_zone_features(lat, lon, direction, distance_m)
zone = Zone(
direction=direction,
distance_km=distance_m / 1000,
forest_pct=features['forest_pct'],
road_density=features['road_density'],
water_distance_km=features['water_distance_km'],
settlement_distance_km=features['settlement_distance_km']
)
zones.append(zone)
return zones