Add running route planner with quality comparison and GPX export
This commit is contained in:
Peter Stockings
@@ -0,0 +1,271 @@
|
||||
"""Bounded pedestrian route generation and explainable ranking, using ORS GeoJSON."""
|
||||
|
||||
import itertools
|
||||
import json
|
||||
import math
|
||||
from collections import Counter
|
||||
from urllib.error import HTTPError, URLError
|
||||
from urllib.request import Request, urlopen
|
||||
|
||||
|
||||
class PlannerError(ValueError):
|
||||
pass
|
||||
|
||||
|
||||
def number(value, low, high, label):
|
||||
if isinstance(value, bool):
|
||||
raise PlannerError(f"Choose a valid {label}.")
|
||||
try:
|
||||
value = float(value)
|
||||
except (ValueError, TypeError):
|
||||
raise PlannerError(f"Choose a valid {label}.") from None
|
||||
if not math.isfinite(value) or not low <= value <= high:
|
||||
raise PlannerError(f"{label.capitalize()} must be between {low} and {high}.")
|
||||
return value
|
||||
|
||||
|
||||
def point(value):
|
||||
if not isinstance(value, list) or len(value) != 2:
|
||||
raise PlannerError("Select a location on the map or enter latitude, longitude.")
|
||||
return [number(value[0], -180, 180, "longitude"), number(value[1], -85, 85, "latitude")]
|
||||
|
||||
|
||||
def validate(data):
|
||||
if not isinstance(data, dict):
|
||||
raise PlannerError("Supply route preferences.")
|
||||
mode = data.get("mode", "loop")
|
||||
hill = data.get("hill", "rolling")
|
||||
surface = data.get("surface", "any")
|
||||
priority = data.get("priority", "quality")
|
||||
if not all(isinstance(value, str) for value in (mode, hill, surface, priority)):
|
||||
raise PlannerError("Choose valid route preferences.")
|
||||
if mode not in {"loop", "point"} or hill not in {"flat", "rolling", "hilly"}:
|
||||
raise PlannerError("Choose a valid route type and hill difficulty.")
|
||||
if surface not in {"any", "paved", "trails"} or priority not in {
|
||||
"quality",
|
||||
"distance",
|
||||
"hills",
|
||||
}:
|
||||
raise PlannerError("Choose valid route preferences.")
|
||||
start = point(data.get("start"))
|
||||
end = start if mode == "loop" else point(data.get("end"))
|
||||
distance = number(data.get("distance"), 1, 50, "distance (km)") * 1000
|
||||
if not isinstance(data.get("avoid_repeats", True), bool):
|
||||
raise PlannerError("Choose whether to minimise repeated sections.")
|
||||
if mode == "point" and metres(start, end) > distance * 1.05:
|
||||
raise PlannerError("Your target distance is shorter than the straight-line distance.")
|
||||
return {
|
||||
"mode": mode,
|
||||
"hill": hill,
|
||||
"surface": surface,
|
||||
"priority": priority,
|
||||
"start": start,
|
||||
"end": end,
|
||||
"distance": distance,
|
||||
"avoid_repeats": bool(data.get("avoid_repeats", True)),
|
||||
}
|
||||
|
||||
|
||||
def metres(a, b):
|
||||
lat1, lat2 = math.radians(a[1]), math.radians(b[1])
|
||||
dlat, dlon = lat2 - lat1, math.radians(b[0] - a[0])
|
||||
value = math.sin(dlat / 2) ** 2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon / 2) ** 2
|
||||
return 12742000 * math.asin(math.sqrt(min(1, value)))
|
||||
|
||||
|
||||
def destination(origin, distance, bearing):
|
||||
lat, lon = map(math.radians, [origin[1], origin[0]])
|
||||
angle, bearing = distance / 6371000, math.radians(bearing)
|
||||
target = math.asin(
|
||||
math.sin(lat) * math.cos(angle) + math.cos(lat) * math.sin(angle) * math.cos(bearing)
|
||||
)
|
||||
longitude = lon + math.atan2(
|
||||
math.sin(bearing) * math.sin(angle) * math.cos(lat),
|
||||
math.cos(angle) - math.sin(lat) * math.sin(target),
|
||||
)
|
||||
return [(math.degrees(longitude) + 180) % 360 - 180, math.degrees(target)]
|
||||
|
||||
|
||||
def provider_request(body, key):
|
||||
request = Request(
|
||||
"https://api.openrouteservice.org/v2/directions/foot-walking/geojson",
|
||||
data=json.dumps(body).encode(),
|
||||
headers={"Authorization": key, "Content-Type": "application/json"},
|
||||
method="POST",
|
||||
)
|
||||
try:
|
||||
with urlopen(request, timeout=12) as response:
|
||||
raw = response.read(4_000_001)
|
||||
if len(raw) > 4_000_000:
|
||||
raise PlannerError("Routing returned too much data. Try a shorter route.")
|
||||
result = json.loads(raw)
|
||||
if not isinstance(result, dict) or not isinstance(result.get("features"), list):
|
||||
raise PlannerError("The routing service returned an incomplete response.")
|
||||
return result
|
||||
except HTTPError as error:
|
||||
if error.code in (401, 403):
|
||||
raise PlannerError(
|
||||
"The routing service key was rejected. Contact the site operator."
|
||||
) from None
|
||||
if error.code == 429:
|
||||
raise PlannerError("The routing service is busy. Please try again later.") from None
|
||||
raise PlannerError("No route was found for these locations. Try nearby paths.") from None
|
||||
except (URLError, TimeoutError, OSError, json.JSONDecodeError):
|
||||
raise PlannerError("The routing service could not be reached. Please try again.") from None
|
||||
|
||||
|
||||
def coverage(extras, name, codes, total):
|
||||
if not isinstance(extras, dict) or not isinstance(extras.get(name, {}), dict):
|
||||
raise PlannerError("The routing service returned invalid quality data.")
|
||||
summary = extras.get(name, {}).get("summary", [])
|
||||
if not summary:
|
||||
return None
|
||||
if not isinstance(summary, list) or any(not isinstance(row, dict) for row in summary):
|
||||
raise PlannerError("The routing service returned invalid quality data.")
|
||||
distance = sum(
|
||||
number(row.get("distance", 0), 0, 200000, "quality distance")
|
||||
for row in summary
|
||||
if isinstance(row.get("value"), (int, float))
|
||||
and not isinstance(row.get("value"), bool)
|
||||
and row.get("value") in codes
|
||||
)
|
||||
return round(min(100, distance / total * 100), 1)
|
||||
|
||||
|
||||
def analyse(feature, preferences):
|
||||
if not isinstance(feature, dict) or not isinstance(feature.get("geometry"), dict):
|
||||
raise PlannerError("The routing service returned invalid route data.")
|
||||
geometry = feature.get("geometry", {})
|
||||
if geometry.get("type") != "LineString":
|
||||
raise PlannerError("The routing service returned an unsupported route.")
|
||||
coordinates = geometry.get("coordinates", [])
|
||||
if not isinstance(coordinates, list) or not 2 <= len(coordinates) <= 50000:
|
||||
raise PlannerError("The routing service returned an incomplete route.")
|
||||
clean = []
|
||||
for coordinate in coordinates:
|
||||
if not isinstance(coordinate, list) or len(coordinate) not in (2, 3):
|
||||
raise PlannerError("The routing service returned invalid coordinates.")
|
||||
xy = point(coordinate[:2])
|
||||
if len(coordinate) == 3:
|
||||
xy.append(number(coordinate[2], -12000, 12000, "elevation"))
|
||||
clean.append(xy)
|
||||
props = feature.get("properties", {})
|
||||
if not isinstance(props, dict) or not isinstance(props.get("summary"), dict):
|
||||
raise PlannerError("The routing service returned an incomplete summary.")
|
||||
summary = props.get("summary", {})
|
||||
distance = number(summary.get("distance"), 1, 200000, "route distance")
|
||||
ascent = props.get("ascent", summary.get("ascent"))
|
||||
ascent = number(ascent, 0, 20000, "ascent") if ascent is not None else None
|
||||
extras = props.get("extras", {})
|
||||
paths = coverage(extras, "waytype", {4, 5, 6, 7}, distance)
|
||||
roads = coverage(extras, "waytype", {1, 2, 3}, distance)
|
||||
paved = coverage(extras, "surface", {1, 3, 4, 5, 6, 14}, distance)
|
||||
unknown = coverage(extras, "surface", {0}, distance)
|
||||
if unknown is not None:
|
||||
known = coverage(extras, "surface", set(range(1, 19)), distance)
|
||||
unknown = round(max(unknown, 100 - (known or 0)), 1)
|
||||
# Approximate repeated edges, irrespective of travel direction; not a traffic/safety score.
|
||||
seen = Counter()
|
||||
repeated = 0
|
||||
travelled = 0
|
||||
profile = [[0, clean[0][2]]] if len(clean[0]) == 3 else []
|
||||
for a, b in itertools.pairwise(clean):
|
||||
length = metres(a, b)
|
||||
edge = tuple(sorted((tuple(round(x, 5) for x in a[:2]), tuple(round(x, 5) for x in b[:2]))))
|
||||
if seen[edge]:
|
||||
repeated += length
|
||||
seen[edge] += 1
|
||||
travelled += length
|
||||
if len(b) == 3:
|
||||
profile.append([travelled, b[2]])
|
||||
if len(profile) != len(clean):
|
||||
profile = []
|
||||
elif travelled:
|
||||
profile = [[round(p[0] / travelled * distance, 1), p[1]] for p in profile]
|
||||
sampled_profile = profile[:: max(1, math.ceil(len(profile) / 600))]
|
||||
if profile and sampled_profile[-1] != profile[-1]:
|
||||
sampled_profile.append(profile[-1])
|
||||
repeats = round(100 * repeated / travelled, 1) if travelled else 0
|
||||
deviation = (distance - preferences["distance"]) / preferences["distance"] * 100
|
||||
climb_target = {"flat": 3, "rolling": 12, "hilly": 25}[preferences["hill"]]
|
||||
hill_penalty = abs(ascent / (distance / 1000) - climb_target) if ascent is not None else 25
|
||||
quality_penalty = (roads if roads is not None else 50) / 5
|
||||
quality_penalty += (unknown if unknown is not None else 100) / 10
|
||||
if preferences["avoid_repeats"]:
|
||||
quality_penalty += repeats / 3
|
||||
if preferences["surface"] == "paved":
|
||||
quality_penalty += (100 - paved if paved is not None else 100) / 5
|
||||
if preferences["surface"] == "trails":
|
||||
quality_penalty += (100 - paths if paths is not None else 100) / 5
|
||||
score = abs(deviation) * (3 if preferences["priority"] == "distance" else 1)
|
||||
score += hill_penalty * (2 if preferences["priority"] == "hills" else 0.5)
|
||||
score += quality_penalty * (2 if preferences["priority"] == "quality" else 0.5)
|
||||
return {
|
||||
"coordinates": clean,
|
||||
"distance_m": distance,
|
||||
"ascent_m": ascent,
|
||||
"profile": sampled_profile,
|
||||
"quality": {
|
||||
"paths_pct": paths,
|
||||
"roads_pct": roads,
|
||||
"paved_pct": paved,
|
||||
"unknown_surface_pct": unknown,
|
||||
"repeats_pct": repeats,
|
||||
},
|
||||
"deviation_pct": round(deviation, 1),
|
||||
"within_target": abs(deviation) <= 5,
|
||||
"preferences": preferences,
|
||||
"rank_cost": score,
|
||||
}
|
||||
|
||||
|
||||
def generate(preferences, key):
|
||||
if not key:
|
||||
raise PlannerError("Route generation needs ORS_API_KEY configured by the site operator.")
|
||||
routes, errors, signatures = [], [], set()
|
||||
for index in range(3):
|
||||
body = {
|
||||
"elevation": True,
|
||||
"instructions": False,
|
||||
"extra_info": ["surface", "waytype"],
|
||||
"coordinates": [preferences["start"]],
|
||||
}
|
||||
if preferences["mode"] == "loop":
|
||||
body["options"] = {
|
||||
"round_trip": {"length": preferences["distance"], "points": 3, "seed": index + 1}
|
||||
}
|
||||
else:
|
||||
start, end = preferences["start"], preferences["end"]
|
||||
body["coordinates"] = [start, end]
|
||||
if index:
|
||||
separation = metres(start, end)
|
||||
midpoint = [(start[0] + end[0]) / 2, (start[1] + end[1]) / 2]
|
||||
offset = (
|
||||
math.sqrt(max(0, (preferences["distance"] / 2) ** 2 - (separation / 2) ** 2))
|
||||
* 0.65
|
||||
)
|
||||
bearing = math.degrees(
|
||||
math.atan2(
|
||||
(end[0] - start[0]) * math.cos(math.radians(midpoint[1])), end[1] - start[1]
|
||||
)
|
||||
)
|
||||
waypoint = destination(
|
||||
midpoint, max(100, offset), bearing + (90 if index == 1 else -90)
|
||||
)
|
||||
body["coordinates"] = [start, waypoint, end]
|
||||
try:
|
||||
result = provider_request(body, key)
|
||||
for feature in result.get("features", [])[:1]:
|
||||
route = analyse(feature, preferences)
|
||||
signature = json.dumps(route["coordinates"])
|
||||
if signature not in signatures:
|
||||
signatures.add(signature)
|
||||
routes.append(route)
|
||||
except PlannerError as error:
|
||||
errors.append(str(error))
|
||||
if not routes:
|
||||
raise PlannerError(errors[0] if errors else "No routes were found. Try a different start.")
|
||||
routes.sort(key=lambda route: route["rank_cost"])
|
||||
for index, route in enumerate(routes):
|
||||
route["label"] = "Best match" if index == 0 else f"Alternative {index + 1}"
|
||||
return routes, errors
|
||||
Reference in New Issue
Block a user