Spaces:
Sleeping
Sleeping
Update app.py
Browse files
app.py
CHANGED
|
@@ -10,24 +10,13 @@
|
|
| 10 |
# - Guided mode computes a realistic turn between 18:22→18:40 and then southern path.
|
| 11 |
# - Guided validation timestamps: 18:22, 18:40, 19:40, 20:40, 21:40, 22:40, 00:10, 00:20
|
| 12 |
#
|
| 13 |
-
#
|
| 14 |
-
#
|
| 15 |
-
#
|
| 16 |
-
#
|
| 17 |
-
#
|
| 18 |
-
#
|
| 19 |
-
#
|
| 20 |
-
# - Southern trajectory uses PER-SEGMENT speeds (one per ping interval) in [430, 510] kt
|
| 21 |
-
# to match Ping Rings 3..7 and 00:20 as tightly as possible.
|
| 22 |
-
# - This is much faster than brute-forcing a huge grid.
|
| 23 |
-
# 4) HARD BTO constraint applies at:
|
| 24 |
-
# 18:40, 19:40, 20:40, 21:40, 22:40, 00:10, 00:20
|
| 25 |
-
# The solver tries to make these as close as possible (and will usually get extremely tight),
|
| 26 |
-
# while keeping runtime low.
|
| 27 |
-
#
|
| 28 |
-
# NOTE:
|
| 29 |
-
# - “Exactly at ping rings” is interpreted as: at each listed waypoint time, the computed BTO
|
| 30 |
-
# matches the measured BTO as closely as numerically possible within the allowed speed ranges.
|
| 31 |
|
| 32 |
import os
|
| 33 |
import re
|
|
@@ -79,7 +68,6 @@ RING_TIMES = {
|
|
| 79 |
BFO_SCALE_HZ = 5.0
|
| 80 |
BTO_SCALE_US = 50.0
|
| 81 |
|
| 82 |
-
# EXPLICIT WEIGHTS (used only for combined reporting; solver is BTO-first)
|
| 83 |
DEFAULT_WEIGHT_BFO = 1.0
|
| 84 |
DEFAULT_WEIGHT_BTO = 2.0
|
| 85 |
|
|
@@ -94,7 +82,6 @@ RADAR_FIX_TIME_UTC = datetime(2014, 3, 7, 18, 22, 0, tzinfo=timezone.utc)
|
|
| 94 |
|
| 95 |
TURN_COMPLETE_TIME_UTC = datetime(2014, 3, 7, 18, 40, 0, tzinfo=timezone.utc)
|
| 96 |
|
| 97 |
-
# Initial track shown in your radar annotation
|
| 98 |
RADAR_INITIAL_TRACK_DEG = 296.0
|
| 99 |
|
| 100 |
# ===== Guided validation timestamps of interest =====
|
|
@@ -116,7 +103,7 @@ TURN_BANK_MAX_DEG = 30.0
|
|
| 116 |
TURN_BANK_DEFAULT_DEG = 25.0
|
| 117 |
TURN_STEP_SEC = 30.0
|
| 118 |
|
| 119 |
-
# =====
|
| 120 |
TURN_SPEED_MIN_KT = 470.0
|
| 121 |
TURN_SPEED_MAX_KT = 520.0
|
| 122 |
|
|
@@ -131,21 +118,13 @@ REFINE_HEADING_HALF_WINDOW_DEG = 2.0
|
|
| 131 |
REFINE_HEADING_STEP_DEG = 0.5
|
| 132 |
|
| 133 |
# Southern ping-interval schedule (18:40→00:20 = 340 min)
|
| 134 |
-
# segments:
|
| 135 |
-
# 18:40→19:40 (60)
|
| 136 |
-
# 19:40→20:40 (60)
|
| 137 |
-
# 20:40→21:40 (60)
|
| 138 |
-
# 21:40→22:40 (60)
|
| 139 |
-
# 22:40→00:10 (90)
|
| 140 |
-
# 00:10→00:20 (10)
|
| 141 |
GUIDED_SOUTH_SEG_MIN = [60.0, 60.0, 60.0, 60.0, 90.0, 10.0]
|
| 142 |
GUIDED_SOUTH_SEG_LABELS_END = ["19:40:00", "20:40:00", "21:40:00", "22:40:00", "24:10:00", "24:20:00"]
|
| 143 |
|
| 144 |
-
# ===== BTO
|
| 145 |
-
|
| 146 |
-
BTO_HARD_TOL_US = 10.0
|
| 147 |
|
| 148 |
-
# HARD BTO
|
| 149 |
HARD_BTO_LABELS = ["18:40:00", "19:40:00", "20:40:00", "21:40:00", "22:40:00", "24:10:00", "24:20:00"]
|
| 150 |
|
| 151 |
# Search area colors (optional)
|
|
@@ -1130,8 +1109,22 @@ def _refine_range(center: float, half_window: float, step: float, lo: float, hi:
|
|
| 1130 |
return _frange(lo2, hi2, step)
|
| 1131 |
|
| 1132 |
|
| 1133 |
-
def
|
| 1134 |
-
return
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1135 |
|
| 1136 |
|
| 1137 |
def _simulate_one_leg_end(method: str, lat0: float, lon0: float, hdg: float, spd_kt: float, minutes: float):
|
|
@@ -1143,20 +1136,14 @@ def _simulate_one_leg_end(method: str, lat0: float, lon0: float, hdg: float, spd
|
|
| 1143 |
return float(lat1), float(lon1)
|
| 1144 |
|
| 1145 |
|
| 1146 |
-
def
|
| 1147 |
alt_now_ft = altitude_ft_at_dt(dt_ping, start_dt, float(alt_ft), float(vs_fpm), descent_start_dt)
|
| 1148 |
alt_m = alt_now_ft * 0.3048
|
| 1149 |
return float(fp.calculateBTO_at_time(float(lat), float(lon), alt_m, dt_ping))
|
| 1150 |
|
| 1151 |
|
| 1152 |
-
def
|
| 1153 |
-
|
| 1154 |
-
return None if bto is None else float(bto)
|
| 1155 |
-
|
| 1156 |
-
|
| 1157 |
-
def _golden_min_abs(fn, lo: float, hi: float, iters: int = 26):
|
| 1158 |
-
# minimize |fn(x)| over [lo, hi] using golden section on abs value
|
| 1159 |
-
grc = (math.sqrt(5.0) - 1.0) / 2.0 # ~0.618
|
| 1160 |
a, b = float(lo), float(hi)
|
| 1161 |
c = b - grc * (b - a)
|
| 1162 |
d = a + grc * (b - a)
|
|
@@ -1189,40 +1176,29 @@ def _fit_segment_speed(
|
|
| 1189 |
vmin: float,
|
| 1190 |
vmax: float,
|
| 1191 |
):
|
| 1192 |
-
meas =
|
| 1193 |
if meas is None:
|
| 1194 |
-
# if no measurement, just pick mid
|
| 1195 |
v = 0.5 * (vmin + vmax)
|
| 1196 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v, minutes)
|
| 1197 |
-
|
| 1198 |
-
return float(v), float(bto - (bto if bto is None else bto)), (lat1, lon1)
|
| 1199 |
|
| 1200 |
def err(v):
|
| 1201 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v, minutes)
|
| 1202 |
-
bto =
|
| 1203 |
return float(bto - meas)
|
| 1204 |
|
| 1205 |
e_lo = float(err(vmin))
|
| 1206 |
e_hi = float(err(vmax))
|
| 1207 |
|
| 1208 |
-
# If sign change
|
| 1209 |
-
if e_lo == 0.0:
|
| 1210 |
-
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, vmin, minutes)
|
| 1211 |
-
return float(vmin), 0.0, (lat1, lon1)
|
| 1212 |
-
if e_hi == 0.0:
|
| 1213 |
-
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, vmax, minutes)
|
| 1214 |
-
return float(vmax), 0.0, (lat1, lon1)
|
| 1215 |
-
|
| 1216 |
if (e_lo < 0.0 and e_hi > 0.0) or (e_lo > 0.0 and e_hi < 0.0):
|
| 1217 |
a, b = float(vmin), float(vmax)
|
| 1218 |
fa, fb = e_lo, e_hi
|
| 1219 |
-
# 28 iters is plenty for sub-1e-6 kt resolution; we’ll stop early if close
|
| 1220 |
for _ in range(28):
|
| 1221 |
m = 0.5 * (a + b)
|
| 1222 |
fm = float(err(m))
|
| 1223 |
-
if abs(fm) < 0.2:
|
| 1224 |
a = b = m
|
| 1225 |
-
fa = fb = fm
|
| 1226 |
break
|
| 1227 |
if (fa < 0 and fm > 0) or (fa > 0 and fm < 0):
|
| 1228 |
b, fb = m, fm
|
|
@@ -1230,14 +1206,12 @@ def _fit_segment_speed(
|
|
| 1230 |
a, fa = m, fm
|
| 1231 |
v_best = 0.5 * (a + b)
|
| 1232 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v_best, minutes)
|
| 1233 |
-
|
| 1234 |
-
return float(v_best), float(final_err), (lat1, lon1)
|
| 1235 |
|
| 1236 |
-
# No sign change:
|
| 1237 |
-
v_best,
|
| 1238 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v_best, minutes)
|
| 1239 |
-
|
| 1240 |
-
return float(v_best), float(final_err), (lat1, lon1)
|
| 1241 |
|
| 1242 |
|
| 1243 |
def _build_south_legs_piecewise(h_south: float, speeds: list[float], minutes: list[float]):
|
|
@@ -1250,10 +1224,6 @@ def _build_south_legs_piecewise(h_south: float, speeds: list[float], minutes: li
|
|
| 1250 |
|
| 1251 |
|
| 1252 |
# ================= Guided Solver =================
|
| 1253 |
-
def _guided_start_dt():
|
| 1254 |
-
return RADAR_FIX_TIME_UTC
|
| 1255 |
-
|
| 1256 |
-
|
| 1257 |
def solve_guided_straight(
|
| 1258 |
method, alt_ft, vs_fpm,
|
| 1259 |
descent_start_choice, weight_bfo, weight_bto
|
|
@@ -1267,69 +1237,37 @@ def solve_guided_straight(
|
|
| 1267 |
elif (descent_start_choice or "").startswith("After Ping Ring 7"):
|
| 1268 |
descent_start_dt = dt_ping7_utc()
|
| 1269 |
|
| 1270 |
-
# Precompute ping points window for later reporting/scoring (BFO plot still uses full set)
|
| 1271 |
pts_all = _ping_points_after(start_dt)
|
| 1272 |
pts_window = [(lbl, dt, bfo, bto) for (lbl, dt, bfo, bto) in pts_all if start_dt <= dt <= end_dt]
|
| 1273 |
if len(pts_window) < 2:
|
| 1274 |
return None, "Guided solve: not enough ping points to score.", None
|
| 1275 |
|
| 1276 |
-
|
| 1277 |
-
|
| 1278 |
-
|
| 1279 |
-
meas_bto_1840 = _measured_bto_at(dt_1840)
|
| 1280 |
|
| 1281 |
-
|
| 1282 |
-
return None, "Guided solve: missing measured BTO at 18:40 in satLoc.", None
|
| 1283 |
-
|
| 1284 |
-
def evaluate_heading(h_south: float):
|
| 1285 |
-
"""
|
| 1286 |
-
For a given south heading:
|
| 1287 |
-
1) Choose a turn speed in [470,520] to match BTO at 18:40 (Ping Ring 2).
|
| 1288 |
-
2) From the 18:40 endpoint, choose per-segment speeds in [430,510] to match
|
| 1289 |
-
BTO at 19:40..00:20 endpoints sequentially.
|
| 1290 |
-
Returns: (total_abs_err_sum, max_abs_err, legs, details)
|
| 1291 |
-
"""
|
| 1292 |
-
h_south = float(h_south)
|
| 1293 |
-
|
| 1294 |
-
# --- Turn speed fit: minimize |BTO(18:40)-meas|
|
| 1295 |
-
def turn_err(v_turn):
|
| 1296 |
-
turn_legs, _ = build_turn_legs(
|
| 1297 |
-
speed_knots=float(v_turn),
|
| 1298 |
-
method=method,
|
| 1299 |
-
heading_final_deg=h_south,
|
| 1300 |
-
)
|
| 1301 |
-
# endpoint at 18:40 is the end of turn legs
|
| 1302 |
-
pts, _seg, _tmin = simulate_path(RADAR_FIX_LAT, RADAR_FIX_LON, turn_legs, method)
|
| 1303 |
-
lat_1840, lon_1840 = pts[-1]
|
| 1304 |
-
bto = _bto_at(method, lat_1840, lon_1840, dt_1840, alt_ft, vs_fpm, start_dt, descent_start_dt)
|
| 1305 |
-
return float(bto - meas_bto_1840)
|
| 1306 |
-
|
| 1307 |
-
v_turn_best, abs_turn_err = _golden_min_abs(turn_err, TURN_SPEED_MIN_KT, TURN_SPEED_MAX_KT, iters=28)
|
| 1308 |
|
|
|
|
|
|
|
| 1309 |
turn_legs, turn_info = build_turn_legs(
|
| 1310 |
-
speed_knots=float(
|
| 1311 |
method=method,
|
| 1312 |
-
heading_final_deg=h_south,
|
| 1313 |
)
|
| 1314 |
pts_turn, _, _ = simulate_path(RADAR_FIX_LAT, RADAR_FIX_LON, turn_legs, method)
|
| 1315 |
lat_cur, lon_cur = pts_turn[-1]
|
| 1316 |
|
| 1317 |
-
#
|
| 1318 |
south_speeds = []
|
| 1319 |
south_errs = []
|
| 1320 |
-
seg_end_positions = []
|
| 1321 |
-
dt_prev = dt_1840
|
| 1322 |
-
|
| 1323 |
-
total_abs = float(abs_turn_err)
|
| 1324 |
-
max_abs = float(abs_turn_err)
|
| 1325 |
-
|
| 1326 |
for mins, end_lbl in zip(GUIDED_SOUTH_SEG_MIN, GUIDED_SOUTH_SEG_LABELS_END):
|
| 1327 |
dt_end = _dt_from_label(end_lbl)
|
| 1328 |
v_seg, err_seg, (lat1, lon1) = _fit_segment_speed(
|
| 1329 |
method=method,
|
| 1330 |
lat0=lat_cur,
|
| 1331 |
lon0=lon_cur,
|
| 1332 |
-
hdg=h_south,
|
| 1333 |
minutes=float(mins),
|
| 1334 |
dt_end=dt_end,
|
| 1335 |
alt_ft=alt_ft,
|
|
@@ -1341,98 +1279,138 @@ def solve_guided_straight(
|
|
| 1341 |
)
|
| 1342 |
south_speeds.append(float(v_seg))
|
| 1343 |
south_errs.append(float(err_seg))
|
| 1344 |
-
seg_end_positions.append((float(lat1), float(lon1)))
|
| 1345 |
lat_cur, lon_cur = float(lat1), float(lon1)
|
| 1346 |
-
dt_prev = dt_end
|
| 1347 |
|
| 1348 |
-
|
| 1349 |
-
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1350 |
|
| 1351 |
-
|
| 1352 |
-
|
| 1353 |
-
|
|
|
|
|
|
|
|
|
|
| 1354 |
|
| 1355 |
details = {
|
| 1356 |
"h_south": float(h_south),
|
| 1357 |
-
"v_turn": float(
|
| 1358 |
-
"turn_abs_err_us": float(abs_turn_err),
|
| 1359 |
"south_speeds": south_speeds,
|
| 1360 |
"south_errs_us": south_errs,
|
| 1361 |
-
"max_abs_err_us": float(max_abs),
|
| 1362 |
-
"sum_abs_err_us": float(total_abs),
|
| 1363 |
"turn_info": turn_info,
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1364 |
}
|
| 1365 |
-
return
|
|
|
|
|
|
|
|
|
|
| 1366 |
|
| 1367 |
-
|
| 1368 |
-
|
| 1369 |
-
|
| 1370 |
-
|
|
|
|
| 1371 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1372 |
coarse_headings = _frange(GUIDED_HEADING_MIN, GUIDED_HEADING_MAX, COARSE_HEADING_STEP_DEG)
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1373 |
for h in coarse_headings:
|
| 1374 |
-
|
| 1375 |
-
|
| 1376 |
-
|
| 1377 |
-
|
| 1378 |
-
|
| 1379 |
-
|
| 1380 |
-
|
| 1381 |
-
|
| 1382 |
-
|
| 1383 |
-
|
|
|
|
|
|
|
| 1384 |
refine_headings = _refine_range(h0, REFINE_HEADING_HALF_WINDOW_DEG, REFINE_HEADING_STEP_DEG, GUIDED_HEADING_MIN, GUIDED_HEADING_MAX)
|
| 1385 |
|
| 1386 |
for h in refine_headings:
|
| 1387 |
-
|
| 1388 |
-
|
| 1389 |
-
if best is None or key < best:
|
| 1390 |
-
best = key
|
| 1391 |
-
best_details = details
|
| 1392 |
-
best_legs = legs
|
| 1393 |
-
|
| 1394 |
-
# Build a compact report
|
| 1395 |
-
south_speeds = best_details["south_speeds"]
|
| 1396 |
-
south_errs = best_details["south_errs_us"]
|
| 1397 |
-
|
| 1398 |
-
# Also compute RMSE_BTO on the enforced timestamps for reporting
|
| 1399 |
-
enforced_dts = [_dt_from_label(lbl) for lbl in HARD_BTO_LABELS]
|
| 1400 |
-
bto_errors = []
|
| 1401 |
-
for dt_k in enforced_dts:
|
| 1402 |
-
st = _state_at_dt(RADAR_FIX_LAT, RADAR_FIX_LON, best_legs, method, start_dt, dt_k)
|
| 1403 |
-
if st is None:
|
| 1404 |
-
continue
|
| 1405 |
-
meas = _measured_bto_at(dt_k)
|
| 1406 |
-
if meas is None:
|
| 1407 |
continue
|
| 1408 |
-
|
| 1409 |
-
|
|
|
|
| 1410 |
|
| 1411 |
-
|
|
|
|
|
|
|
| 1412 |
|
| 1413 |
report_lines = [
|
| 1414 |
"**Guided mode: Straight + computed turn (18:22→18:40), clipped to 00:20**",
|
| 1415 |
"- Solver mode: **BTO-first ring-fit (fast)**",
|
| 1416 |
-
f"- Best southern heading: **{
|
| 1417 |
-
f"- Turn speed
|
| 1418 |
-
f"- Turn BTO error at 18:40: **{best_details['turn_abs_err_us']:.2f} µs**",
|
| 1419 |
"",
|
| 1420 |
"**Southern per-segment speeds (each constrained to 430–510 kt):**",
|
| 1421 |
-
f"- 18:40→19:40: **{south_speeds[0]:.2f} kt** (
|
| 1422 |
-
f"- 19:40→20:40: **{south_speeds[1]:.2f} kt** (
|
| 1423 |
-
f"- 20:40→21:40: **{south_speeds[2]:.2f} kt** (
|
| 1424 |
-
f"- 21:40→22:40: **{south_speeds[3]:.2f} kt** (
|
| 1425 |
-
f"- 22:40→00:10: **{south_speeds[4]:.2f} kt** (
|
| 1426 |
-
f"- 00:10→00:20: **{south_speeds[5]:.2f} kt** (
|
| 1427 |
"",
|
| 1428 |
-
f"- Enforced BTO timestamps: **{', '.join(
|
| 1429 |
-
f"- Max |BTO error|
|
| 1430 |
-
f"- RMSE_BTO
|
| 1431 |
"",
|
| 1432 |
-
|
| 1433 |
]
|
| 1434 |
|
| 1435 |
-
return
|
| 1436 |
|
| 1437 |
|
| 1438 |
# ================= RUNNERS =================
|
|
@@ -1529,7 +1507,7 @@ with gr.Blocks() as demo:
|
|
| 1529 |
f"- **Final resting point:** **00:20:00Z** (trajectory clipped exactly at this time)\n\n"
|
| 1530 |
f"**Turn speed bounds:** {TURN_SPEED_MIN_KT:.0f}–{TURN_SPEED_MAX_KT:.0f} kt\n"
|
| 1531 |
f"**Southern speed bounds:** {SOUTH_SPEED_MIN_KT:.0f}–{SOUTH_SPEED_MAX_KT:.0f} kt (per segment)\n"
|
| 1532 |
-
f"**Ring-fit goal:** match BTO at
|
| 1533 |
)
|
| 1534 |
|
| 1535 |
with gr.Row():
|
|
@@ -1552,10 +1530,10 @@ with gr.Blocks() as demo:
|
|
| 1552 |
gr.Markdown(
|
| 1553 |
"- Guided mode has **no speed inputs**.\n"
|
| 1554 |
"- Solver uses **fast BTO-first ring-fit** (turn speed + per-segment southern speeds).\n"
|
| 1555 |
-
"-
|
| 1556 |
)
|
| 1557 |
|
| 1558 |
-
gr.Markdown("### Guided scoring weights (for reporting
|
| 1559 |
with gr.Row():
|
| 1560 |
weight_bfo = gr.Slider(0.0, 5.0, value=DEFAULT_WEIGHT_BFO, step=0.1, label="Weight: BFO")
|
| 1561 |
weight_bto = gr.Slider(0.0, 10.0, value=DEFAULT_WEIGHT_BTO, step=0.1, label="Weight: BTO")
|
|
|
|
| 10 |
# - Guided mode computes a realistic turn between 18:22→18:40 and then southern path.
|
| 11 |
# - Guided validation timestamps: 18:22, 18:40, 19:40, 20:40, 21:40, 22:40, 00:10, 00:20
|
| 12 |
#
|
| 13 |
+
# THIS VERSION:
|
| 14 |
+
# - Guided solver is BTO-first "ring-fit", with fast speed optimization
|
| 15 |
+
# - IMPORTANT: If a timestamp has no measured BTO (e.g. 18:40 in your satLoc),
|
| 16 |
+
# the solver automatically does NOT enforce BTO there.
|
| 17 |
+
# - Turn speed range: 470–520 kt
|
| 18 |
+
# - Southern segment speed range: 430–510 kt (per segment)
|
| 19 |
+
# - Waypoints at ping times are driven to ping rings for those timestamps with measured BTO.
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 20 |
|
| 21 |
import os
|
| 22 |
import re
|
|
|
|
| 68 |
BFO_SCALE_HZ = 5.0
|
| 69 |
BTO_SCALE_US = 50.0
|
| 70 |
|
|
|
|
| 71 |
DEFAULT_WEIGHT_BFO = 1.0
|
| 72 |
DEFAULT_WEIGHT_BTO = 2.0
|
| 73 |
|
|
|
|
| 82 |
|
| 83 |
TURN_COMPLETE_TIME_UTC = datetime(2014, 3, 7, 18, 40, 0, tzinfo=timezone.utc)
|
| 84 |
|
|
|
|
| 85 |
RADAR_INITIAL_TRACK_DEG = 296.0
|
| 86 |
|
| 87 |
# ===== Guided validation timestamps of interest =====
|
|
|
|
| 103 |
TURN_BANK_DEFAULT_DEG = 25.0
|
| 104 |
TURN_STEP_SEC = 30.0
|
| 105 |
|
| 106 |
+
# ===== Speed bounds per your request =====
|
| 107 |
TURN_SPEED_MIN_KT = 470.0
|
| 108 |
TURN_SPEED_MAX_KT = 520.0
|
| 109 |
|
|
|
|
| 118 |
REFINE_HEADING_STEP_DEG = 0.5
|
| 119 |
|
| 120 |
# Southern ping-interval schedule (18:40→00:20 = 340 min)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 121 |
GUIDED_SOUTH_SEG_MIN = [60.0, 60.0, 60.0, 60.0, 90.0, 10.0]
|
| 122 |
GUIDED_SOUTH_SEG_LABELS_END = ["19:40:00", "20:40:00", "21:40:00", "22:40:00", "24:10:00", "24:20:00"]
|
| 123 |
|
| 124 |
+
# ===== BTO "tight" target (solver tries tighter; used only for reporting) =====
|
| 125 |
+
BTO_TARGET_TIGHT_US = 1.0 # purely informational
|
|
|
|
| 126 |
|
| 127 |
+
# HARD BTO labels list (will be AUTO-FILTERED to those with measured BTO)
|
| 128 |
HARD_BTO_LABELS = ["18:40:00", "19:40:00", "20:40:00", "21:40:00", "22:40:00", "24:10:00", "24:20:00"]
|
| 129 |
|
| 130 |
# Search area colors (optional)
|
|
|
|
| 1109 |
return _frange(lo2, hi2, step)
|
| 1110 |
|
| 1111 |
|
| 1112 |
+
def _guided_start_dt():
|
| 1113 |
+
return RADAR_FIX_TIME_UTC
|
| 1114 |
+
|
| 1115 |
+
|
| 1116 |
+
def _lookup_measured_bto(dt_ping: datetime):
|
| 1117 |
+
_bfo, bto = _lookup_measured_at_dt(dt_ping)
|
| 1118 |
+
return None if bto is None else float(bto)
|
| 1119 |
+
|
| 1120 |
+
|
| 1121 |
+
def _available_hard_bto_labels():
|
| 1122 |
+
out = []
|
| 1123 |
+
for lbl in HARD_BTO_LABELS:
|
| 1124 |
+
dt = _dt_from_label(lbl)
|
| 1125 |
+
if _lookup_measured_bto(dt) is not None:
|
| 1126 |
+
out.append(lbl)
|
| 1127 |
+
return out
|
| 1128 |
|
| 1129 |
|
| 1130 |
def _simulate_one_leg_end(method: str, lat0: float, lon0: float, hdg: float, spd_kt: float, minutes: float):
|
|
|
|
| 1136 |
return float(lat1), float(lon1)
|
| 1137 |
|
| 1138 |
|
| 1139 |
+
def _bto_calc(method: str, lat: float, lon: float, dt_ping: datetime, alt_ft: float, vs_fpm: float, start_dt: datetime, descent_start_dt: datetime | None):
|
| 1140 |
alt_now_ft = altitude_ft_at_dt(dt_ping, start_dt, float(alt_ft), float(vs_fpm), descent_start_dt)
|
| 1141 |
alt_m = alt_now_ft * 0.3048
|
| 1142 |
return float(fp.calculateBTO_at_time(float(lat), float(lon), alt_m, dt_ping))
|
| 1143 |
|
| 1144 |
|
| 1145 |
+
def _golden_min_abs(fn, lo: float, hi: float, iters: int = 24):
|
| 1146 |
+
grc = (math.sqrt(5.0) - 1.0) / 2.0
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1147 |
a, b = float(lo), float(hi)
|
| 1148 |
c = b - grc * (b - a)
|
| 1149 |
d = a + grc * (b - a)
|
|
|
|
| 1176 |
vmin: float,
|
| 1177 |
vmax: float,
|
| 1178 |
):
|
| 1179 |
+
meas = _lookup_measured_bto(dt_end)
|
| 1180 |
if meas is None:
|
|
|
|
| 1181 |
v = 0.5 * (vmin + vmax)
|
| 1182 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v, minutes)
|
| 1183 |
+
return float(v), float("nan"), (lat1, lon1)
|
|
|
|
| 1184 |
|
| 1185 |
def err(v):
|
| 1186 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v, minutes)
|
| 1187 |
+
bto = _bto_calc(method, lat1, lon1, dt_end, alt_ft, vs_fpm, start_dt, descent_start_dt)
|
| 1188 |
return float(bto - meas)
|
| 1189 |
|
| 1190 |
e_lo = float(err(vmin))
|
| 1191 |
e_hi = float(err(vmax))
|
| 1192 |
|
| 1193 |
+
# If sign change: bisection root (fast + tight)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1194 |
if (e_lo < 0.0 and e_hi > 0.0) or (e_lo > 0.0 and e_hi < 0.0):
|
| 1195 |
a, b = float(vmin), float(vmax)
|
| 1196 |
fa, fb = e_lo, e_hi
|
|
|
|
| 1197 |
for _ in range(28):
|
| 1198 |
m = 0.5 * (a + b)
|
| 1199 |
fm = float(err(m))
|
| 1200 |
+
if abs(fm) < 0.2:
|
| 1201 |
a = b = m
|
|
|
|
| 1202 |
break
|
| 1203 |
if (fa < 0 and fm > 0) or (fa > 0 and fm < 0):
|
| 1204 |
b, fb = m, fm
|
|
|
|
| 1206 |
a, fa = m, fm
|
| 1207 |
v_best = 0.5 * (a + b)
|
| 1208 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v_best, minutes)
|
| 1209 |
+
return float(v_best), float(err(v_best)), (lat1, lon1)
|
|
|
|
| 1210 |
|
| 1211 |
+
# No sign change: minimize |err| on bounds
|
| 1212 |
+
v_best, _ = _golden_min_abs(err, vmin, vmax, iters=28)
|
| 1213 |
lat1, lon1 = _simulate_one_leg_end(method, lat0, lon0, hdg, v_best, minutes)
|
| 1214 |
+
return float(v_best), float(err(v_best)), (lat1, lon1)
|
|
|
|
| 1215 |
|
| 1216 |
|
| 1217 |
def _build_south_legs_piecewise(h_south: float, speeds: list[float], minutes: list[float]):
|
|
|
|
| 1224 |
|
| 1225 |
|
| 1226 |
# ================= Guided Solver =================
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1227 |
def solve_guided_straight(
|
| 1228 |
method, alt_ft, vs_fpm,
|
| 1229 |
descent_start_choice, weight_bfo, weight_bto
|
|
|
|
| 1237 |
elif (descent_start_choice or "").startswith("After Ping Ring 7"):
|
| 1238 |
descent_start_dt = dt_ping7_utc()
|
| 1239 |
|
|
|
|
| 1240 |
pts_all = _ping_points_after(start_dt)
|
| 1241 |
pts_window = [(lbl, dt, bfo, bto) for (lbl, dt, bfo, bto) in pts_all if start_dt <= dt <= end_dt]
|
| 1242 |
if len(pts_window) < 2:
|
| 1243 |
return None, "Guided solve: not enough ping points to score.", None
|
| 1244 |
|
| 1245 |
+
enforced_labels = _available_hard_bto_labels()
|
| 1246 |
+
if not enforced_labels:
|
| 1247 |
+
return None, "Guided solve: no measured BTO available at any enforced ring times.", None
|
|
|
|
| 1248 |
|
| 1249 |
+
enforced_dts = [_dt_from_label(lbl) for lbl in enforced_labels]
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1250 |
|
| 1251 |
+
def evaluate_heading_and_turnspeed(h_south: float, v_turn: float):
|
| 1252 |
+
# Build turn legs
|
| 1253 |
turn_legs, turn_info = build_turn_legs(
|
| 1254 |
+
speed_knots=float(v_turn),
|
| 1255 |
method=method,
|
| 1256 |
+
heading_final_deg=float(h_south),
|
| 1257 |
)
|
| 1258 |
pts_turn, _, _ = simulate_path(RADAR_FIX_LAT, RADAR_FIX_LON, turn_legs, method)
|
| 1259 |
lat_cur, lon_cur = pts_turn[-1]
|
| 1260 |
|
| 1261 |
+
# Fit south segments sequentially (per segment speed)
|
| 1262 |
south_speeds = []
|
| 1263 |
south_errs = []
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1264 |
for mins, end_lbl in zip(GUIDED_SOUTH_SEG_MIN, GUIDED_SOUTH_SEG_LABELS_END):
|
| 1265 |
dt_end = _dt_from_label(end_lbl)
|
| 1266 |
v_seg, err_seg, (lat1, lon1) = _fit_segment_speed(
|
| 1267 |
method=method,
|
| 1268 |
lat0=lat_cur,
|
| 1269 |
lon0=lon_cur,
|
| 1270 |
+
hdg=float(h_south),
|
| 1271 |
minutes=float(mins),
|
| 1272 |
dt_end=dt_end,
|
| 1273 |
alt_ft=alt_ft,
|
|
|
|
| 1279 |
)
|
| 1280 |
south_speeds.append(float(v_seg))
|
| 1281 |
south_errs.append(float(err_seg))
|
|
|
|
| 1282 |
lat_cur, lon_cur = float(lat1), float(lon1)
|
|
|
|
| 1283 |
|
| 1284 |
+
south_legs = _build_south_legs_piecewise(float(h_south), south_speeds, GUIDED_SOUTH_SEG_MIN)
|
| 1285 |
+
full_legs = clip_legs_to_end_dt(turn_legs + south_legs, start_dt, end_dt)
|
| 1286 |
+
|
| 1287 |
+
# Compute BTO errors at enforced timestamps
|
| 1288 |
+
bto_errs = []
|
| 1289 |
+
for dt_k in enforced_dts:
|
| 1290 |
+
meas = _lookup_measured_bto(dt_k)
|
| 1291 |
+
if meas is None:
|
| 1292 |
+
continue
|
| 1293 |
+
st = _state_at_dt(RADAR_FIX_LAT, RADAR_FIX_LON, full_legs, method, start_dt, dt_k)
|
| 1294 |
+
if st is None:
|
| 1295 |
+
continue
|
| 1296 |
+
bto_calc = _bto_calc(method, st["lat"], st["lon"], dt_k, alt_ft, vs_fpm, start_dt, descent_start_dt)
|
| 1297 |
+
bto_errs.append(float(bto_calc - meas))
|
| 1298 |
|
| 1299 |
+
if not bto_errs:
|
| 1300 |
+
return None
|
| 1301 |
+
|
| 1302 |
+
max_abs = float(np.max(np.abs(bto_errs)))
|
| 1303 |
+
sum_abs = float(np.sum(np.abs(bto_errs)))
|
| 1304 |
+
rmse = float(np.sqrt(np.mean(np.array(bto_errs) ** 2)))
|
| 1305 |
|
| 1306 |
details = {
|
| 1307 |
"h_south": float(h_south),
|
| 1308 |
+
"v_turn": float(v_turn),
|
|
|
|
| 1309 |
"south_speeds": south_speeds,
|
| 1310 |
"south_errs_us": south_errs,
|
|
|
|
|
|
|
| 1311 |
"turn_info": turn_info,
|
| 1312 |
+
"enforced_labels": enforced_labels,
|
| 1313 |
+
"max_abs_err_us": max_abs,
|
| 1314 |
+
"sum_abs_err_us": sum_abs,
|
| 1315 |
+
"rmse_bto_us": rmse,
|
| 1316 |
}
|
| 1317 |
+
return (max_abs, sum_abs, rmse), full_legs, details
|
| 1318 |
+
|
| 1319 |
+
def best_for_heading(h_south: float):
|
| 1320 |
+
h_south = float(h_south)
|
| 1321 |
|
| 1322 |
+
# Coarse scan turn speed (fast, avoids expensive nested golden everywhere)
|
| 1323 |
+
coarse_turns = list(np.arange(TURN_SPEED_MIN_KT, TURN_SPEED_MAX_KT + 1e-9, 5.0))
|
| 1324 |
+
best_key = None
|
| 1325 |
+
best_legs = None
|
| 1326 |
+
best_details = None
|
| 1327 |
|
| 1328 |
+
for v in coarse_turns:
|
| 1329 |
+
out = evaluate_heading_and_turnspeed(h_south, v)
|
| 1330 |
+
if out is None:
|
| 1331 |
+
continue
|
| 1332 |
+
key, legs, det = out
|
| 1333 |
+
if best_key is None or key < best_key:
|
| 1334 |
+
best_key, best_legs, best_details = key, legs, det
|
| 1335 |
+
|
| 1336 |
+
if best_key is None:
|
| 1337 |
+
return None
|
| 1338 |
+
|
| 1339 |
+
# Refine around best coarse speed via golden minimize max_abs_err
|
| 1340 |
+
v0 = float(best_details["v_turn"])
|
| 1341 |
+
v_lo = _clamp(v0 - 7.5, TURN_SPEED_MIN_KT, TURN_SPEED_MAX_KT)
|
| 1342 |
+
v_hi = _clamp(v0 + 7.5, TURN_SPEED_MIN_KT, TURN_SPEED_MAX_KT)
|
| 1343 |
+
|
| 1344 |
+
def obj(v):
|
| 1345 |
+
out2 = evaluate_heading_and_turnspeed(h_south, v)
|
| 1346 |
+
if out2 is None:
|
| 1347 |
+
return 1e9
|
| 1348 |
+
(max_abs, sum_abs, _rmse), _legs, _det = out2
|
| 1349 |
+
# primary max abs, secondary sum abs
|
| 1350 |
+
return float(max_abs) + 1e-3 * float(sum_abs)
|
| 1351 |
+
|
| 1352 |
+
v_best, _ = _golden_min_abs(obj, v_lo, v_hi, iters=22)
|
| 1353 |
+
out3 = evaluate_heading_and_turnspeed(h_south, v_best)
|
| 1354 |
+
if out3 is None:
|
| 1355 |
+
return best_key, best_legs, best_details
|
| 1356 |
+
key3, legs3, det3 = out3
|
| 1357 |
+
return key3, legs3, det3
|
| 1358 |
+
|
| 1359 |
+
# Stage 1: coarse headings
|
| 1360 |
coarse_headings = _frange(GUIDED_HEADING_MIN, GUIDED_HEADING_MAX, COARSE_HEADING_STEP_DEG)
|
| 1361 |
+
global_best = None
|
| 1362 |
+
global_legs = None
|
| 1363 |
+
global_det = None
|
| 1364 |
+
|
| 1365 |
for h in coarse_headings:
|
| 1366 |
+
out = best_for_heading(h)
|
| 1367 |
+
if out is None:
|
| 1368 |
+
continue
|
| 1369 |
+
key, legs, det = out
|
| 1370 |
+
if global_best is None or key < global_best:
|
| 1371 |
+
global_best, global_legs, global_det = key, legs, det
|
| 1372 |
+
|
| 1373 |
+
if global_legs is None:
|
| 1374 |
+
return None, "Guided solve: no solution found within speed bounds.", None
|
| 1375 |
+
|
| 1376 |
+
# Stage 2: refine headings around best
|
| 1377 |
+
h0 = float(global_det["h_south"])
|
| 1378 |
refine_headings = _refine_range(h0, REFINE_HEADING_HALF_WINDOW_DEG, REFINE_HEADING_STEP_DEG, GUIDED_HEADING_MIN, GUIDED_HEADING_MAX)
|
| 1379 |
|
| 1380 |
for h in refine_headings:
|
| 1381 |
+
out = best_for_heading(h)
|
| 1382 |
+
if out is None:
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1383 |
continue
|
| 1384 |
+
key, legs, det = out
|
| 1385 |
+
if global_best is None or key < global_best:
|
| 1386 |
+
global_best, global_legs, global_det = key, legs, det
|
| 1387 |
|
| 1388 |
+
# Reporting
|
| 1389 |
+
south_speeds = global_det["south_speeds"]
|
| 1390 |
+
south_errs = global_det["south_errs_us"]
|
| 1391 |
|
| 1392 |
report_lines = [
|
| 1393 |
"**Guided mode: Straight + computed turn (18:22→18:40), clipped to 00:20**",
|
| 1394 |
"- Solver mode: **BTO-first ring-fit (fast)**",
|
| 1395 |
+
f"- Best southern heading: **{global_det['h_south']:.2f}°**",
|
| 1396 |
+
f"- Turn speed range: {TURN_SPEED_MIN_KT:.0f}–{TURN_SPEED_MAX_KT:.0f} kt | Best: **{global_det['v_turn']:.2f} kt**",
|
|
|
|
| 1397 |
"",
|
| 1398 |
"**Southern per-segment speeds (each constrained to 430–510 kt):**",
|
| 1399 |
+
f"- 18:40→19:40: **{south_speeds[0]:.2f} kt** (segment end err {south_errs[0]:+.2f} µs)",
|
| 1400 |
+
f"- 19:40→20:40: **{south_speeds[1]:.2f} kt** (segment end err {south_errs[1]:+.2f} µs)",
|
| 1401 |
+
f"- 20:40→21:40: **{south_speeds[2]:.2f} kt** (segment end err {south_errs[2]:+.2f} µs)",
|
| 1402 |
+
f"- 21:40→22:40: **{south_speeds[3]:.2f} kt** (segment end err {south_errs[3]:+.2f} µs)",
|
| 1403 |
+
f"- 22:40→00:10: **{south_speeds[4]:.2f} kt** (segment end err {south_errs[4]:+.2f} µs)",
|
| 1404 |
+
f"- 00:10→00:20: **{south_speeds[5]:.2f} kt** (segment end err {south_errs[5]:+.2f} µs)",
|
| 1405 |
"",
|
| 1406 |
+
f"- Enforced BTO timestamps (measured available): **{', '.join(global_det['enforced_labels'])}**",
|
| 1407 |
+
f"- Max |BTO error| at enforced times: **{global_det['max_abs_err_us']:.2f} µs**",
|
| 1408 |
+
f"- RMSE_BTO at enforced times: **{global_det['rmse_bto_us']:.2f} µs**",
|
| 1409 |
"",
|
| 1410 |
+
"*If 18:40 has no measured BTO in your satLoc, it is automatically not enforced.*",
|
| 1411 |
]
|
| 1412 |
|
| 1413 |
+
return global_legs, "\n".join(report_lines), start_dt
|
| 1414 |
|
| 1415 |
|
| 1416 |
# ================= RUNNERS =================
|
|
|
|
| 1507 |
f"- **Final resting point:** **00:20:00Z** (trajectory clipped exactly at this time)\n\n"
|
| 1508 |
f"**Turn speed bounds:** {TURN_SPEED_MIN_KT:.0f}–{TURN_SPEED_MAX_KT:.0f} kt\n"
|
| 1509 |
f"**Southern speed bounds:** {SOUTH_SPEED_MIN_KT:.0f}–{SOUTH_SPEED_MAX_KT:.0f} kt (per segment)\n"
|
| 1510 |
+
f"**Ring-fit goal:** match BTO at ping times where measured BTO exists (auto-detected from satLoc)."
|
| 1511 |
)
|
| 1512 |
|
| 1513 |
with gr.Row():
|
|
|
|
| 1530 |
gr.Markdown(
|
| 1531 |
"- Guided mode has **no speed inputs**.\n"
|
| 1532 |
"- Solver uses **fast BTO-first ring-fit** (turn speed + per-segment southern speeds).\n"
|
| 1533 |
+
"- If 18:40 has no measured BTO in your satLoc, it is **automatically skipped** for enforcement."
|
| 1534 |
)
|
| 1535 |
|
| 1536 |
+
gr.Markdown("### Guided scoring weights (for reporting only)")
|
| 1537 |
with gr.Row():
|
| 1538 |
weight_bfo = gr.Slider(0.0, 5.0, value=DEFAULT_WEIGHT_BFO, step=0.1, label="Weight: BFO")
|
| 1539 |
weight_bto = gr.Slider(0.0, 10.0, value=DEFAULT_WEIGHT_BTO, step=0.1, label="Weight: BTO")
|