Spaces:
Sleeping
Sleeping
Update app.py
Browse files
app.py
CHANGED
|
@@ -11,9 +11,11 @@
|
|
| 11 |
# - Guided mode computes a realistic turn between 18:22→18:40 and then southern path.
|
| 12 |
# - Guided validation timestamps: 18:22, 18:40, 19:40, 20:40, 21:40, 22:40, 00:10, 00:20
|
| 13 |
#
|
| 14 |
-
#
|
| 15 |
-
#
|
| 16 |
-
#
|
|
|
|
|
|
|
| 17 |
|
| 18 |
import os
|
| 19 |
import re
|
|
@@ -61,12 +63,15 @@ RING_TIMES = {
|
|
| 61 |
"Ping Ring 7": "24:10:00Z",
|
| 62 |
}
|
| 63 |
|
|
|
|
| 64 |
BFO_SCALE_HZ = 5.0
|
| 65 |
BTO_SCALE_US = 50.0
|
| 66 |
|
|
|
|
| 67 |
DEFAULT_WEIGHT_BFO = 1.0
|
| 68 |
DEFAULT_WEIGHT_BTO = 2.0
|
| 69 |
|
|
|
|
| 70 |
SATLOC_IDX_MEAS_BFO = 4
|
| 71 |
SATLOC_IDX_MEAS_BTO = 5
|
| 72 |
|
|
@@ -78,9 +83,11 @@ RADAR_FIX_TIME_UTC = datetime(2014, 3, 7, 18, 22, 0, tzinfo=timezone.utc)
|
|
| 78 |
TURN_COMPLETE_TIME_UTC = datetime(2014, 3, 7, 18, 40, 0, tzinfo=timezone.utc)
|
| 79 |
TURN_WINDOW_MIN = (TURN_COMPLETE_TIME_UTC - RADAR_FIX_TIME_UTC).total_seconds() / 60.0
|
| 80 |
|
|
|
|
| 81 |
RADAR_INITIAL_TRACK_DEG = 296.0
|
| 82 |
|
| 83 |
# ===== Guided validation timestamps of interest =====
|
|
|
|
| 84 |
GUIDED_VALIDATION_LABELS = [
|
| 85 |
"18:22:00",
|
| 86 |
"18:40:00",
|
|
@@ -91,6 +98,7 @@ GUIDED_VALIDATION_LABELS = [
|
|
| 91 |
"24:10:00", # 00:10Z
|
| 92 |
"24:20:00", # 00:20Z
|
| 93 |
]
|
|
|
|
| 94 |
|
| 95 |
# Turn modeling parameters
|
| 96 |
TURN_BANK_MIN_DEG = 15.0
|
|
@@ -98,13 +106,13 @@ TURN_BANK_MAX_DEG = 30.0
|
|
| 98 |
TURN_BANK_DEFAULT_DEG = 25.0
|
| 99 |
TURN_STEP_SEC = 30.0
|
| 100 |
|
| 101 |
-
# Guided heading search
|
| 102 |
GUIDED_HEADING_MIN = 150.0
|
| 103 |
GUIDED_HEADING_MAX = 210.0
|
| 104 |
STRAIGHT_STEP_DEG = 0.5
|
| 105 |
|
| 106 |
-
# Southern schedule
|
| 107 |
-
GUIDED_LEG_SCHEDULE_MIN = [12.0, 60.0, 60.0, 60.0, 60.0, 90.0, 10.0]
|
| 108 |
|
| 109 |
# Drift
|
| 110 |
DEFAULT_DRIFT_MAX_DEG = 2.0
|
|
@@ -411,9 +419,11 @@ def build_turn_legs(
|
|
| 411 |
|
| 412 |
legs = []
|
| 413 |
|
|
|
|
| 414 |
if straight_sec > 1.0:
|
| 415 |
legs.append((h0, float(speed_knots), straight_sec / 60.0))
|
| 416 |
|
|
|
|
| 417 |
if turn_sec > 1.0 and delta_abs > 1e-6:
|
| 418 |
step_sec = max(5.0, float(step_sec))
|
| 419 |
n = int(math.ceil(turn_sec / step_sec))
|
|
@@ -439,45 +449,6 @@ def build_turn_legs(
|
|
| 439 |
return legs, info
|
| 440 |
|
| 441 |
|
| 442 |
-
# ================= SIM CORE + STATS =================
|
| 443 |
-
def simulate_path(start_lat, start_lon, legs, method):
|
| 444 |
-
lat, lon = float(start_lat), float(start_lon)
|
| 445 |
-
pts = [(lat, lon)]
|
| 446 |
-
seg_dists_m = []
|
| 447 |
-
total_time_min = 0.0
|
| 448 |
-
|
| 449 |
-
for heading, speed_knots, minutes in legs:
|
| 450 |
-
dist_m = float(speed_knots) * 1852.0 * (float(minutes) / 60.0)
|
| 451 |
-
if method == "VINCENTY":
|
| 452 |
-
(lat, lon), _ = fp.vincenty_direct((lat, lon), float(heading), dist_m)
|
| 453 |
-
else:
|
| 454 |
-
lat, lon = fp.rhumb_direct((lat, lon), float(heading), dist_m)
|
| 455 |
-
|
| 456 |
-
pts.append((lat, lon))
|
| 457 |
-
seg_dists_m.append(dist_m)
|
| 458 |
-
total_time_min += float(minutes)
|
| 459 |
-
|
| 460 |
-
return pts, seg_dists_m, total_time_min
|
| 461 |
-
|
| 462 |
-
|
| 463 |
-
def compute_stats(path_pts, seg_dists_m, total_time_min, legs):
|
| 464 |
-
total_dist_m = float(np.sum(seg_dists_m)) if seg_dists_m else 0.0
|
| 465 |
-
total_dist_km = total_dist_m / 1000.0
|
| 466 |
-
total_dist_nm = total_dist_m / 1852.0
|
| 467 |
-
|
| 468 |
-
total_time_hr = total_time_min / 60.0
|
| 469 |
-
avg_gs_knots = (total_dist_nm / total_time_hr) if total_time_hr > 0 else 0.0
|
| 470 |
-
|
| 471 |
-
end_lat, end_lon = path_pts[-1]
|
| 472 |
-
return "\n\n".join([
|
| 473 |
-
f"**Legs:** {len(legs)}",
|
| 474 |
-
f"**Total time:** {total_time_min:.1f} min ({total_time_hr:.2f} hr)",
|
| 475 |
-
f"**Total distance:** {total_dist_nm:.1f} nm ({total_dist_km:.1f} km)",
|
| 476 |
-
f"**Average ground speed:** {avg_gs_knots:.1f} kt",
|
| 477 |
-
f"**Final position (simulation end):** lat {end_lat:.4f}, lon {end_lon:.4f}",
|
| 478 |
-
])
|
| 479 |
-
|
| 480 |
-
|
| 481 |
# ================= TIME HELPERS =================
|
| 482 |
def _parse_start_time(s: str) -> datetime:
|
| 483 |
s = (s or "").strip()
|
|
@@ -551,7 +522,46 @@ def build_altitude_plot(start_dt: datetime, total_time_min: float, alt0_ft: floa
|
|
| 551 |
return fig
|
| 552 |
|
| 553 |
|
| 554 |
-
# =================
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 555 |
def _predict_state_from_leg_start(method, leg_start_lat, leg_start_lon, hdg_deg, spd_knots, t_inside_leg_sec):
|
| 556 |
dist_m = float(spd_knots) * 1852.0 * (float(t_inside_leg_sec) / 3600.0)
|
| 557 |
if method == "VINCENTY":
|
|
@@ -627,6 +637,143 @@ def _map_marker_times(start_dt, legs):
|
|
| 627 |
return times
|
| 628 |
|
| 629 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 630 |
# ================= VALIDATION TABLE =================
|
| 631 |
def _lookup_measured_at_dt(dt_ping: datetime):
|
| 632 |
if hasattr(fp, "nearest_measured_for_time"):
|
|
@@ -752,7 +899,8 @@ def build_leaflet_map(path_pts, legs, method, start_dt):
|
|
| 752 |
- Ping rings always drawn
|
| 753 |
- Search areas always drawn (color-coded)
|
| 754 |
- Validation markers are the authoritative markers
|
| 755 |
-
-
|
|
|
|
| 756 |
"""
|
| 757 |
start_lat, start_lon = path_pts[0]
|
| 758 |
is_guided = (start_dt == RADAR_FIX_TIME_UTC)
|
|
@@ -835,7 +983,7 @@ def build_leaflet_map(path_pts, legs, method, start_dt):
|
|
| 835 |
tooltip = f"{ring_name} ({ring_time})" if ring_time else ring_name
|
| 836 |
folium.PolyLine(ring_pts, weight=2, opacity=0.6, color="yellow", tooltip=tooltip).add_to(m)
|
| 837 |
|
| 838 |
-
# Flight path
|
| 839 |
folium.PolyLine(path_pts, color="white", weight=8, opacity=0.3).add_to(m)
|
| 840 |
folium.PolyLine(path_pts, color="red", weight=4, opacity=1.0, tooltip="Flight path").add_to(m)
|
| 841 |
|
|
@@ -852,8 +1000,7 @@ def build_leaflet_map(path_pts, legs, method, start_dt):
|
|
| 852 |
),
|
| 853 |
).add_to(m)
|
| 854 |
|
| 855 |
-
#
|
| 856 |
-
# Manual mode can still show it as a convenience.
|
| 857 |
if not is_guided:
|
| 858 |
folium.CircleMarker(
|
| 859 |
location=path_pts[-1],
|
|
@@ -875,6 +1022,7 @@ def build_leaflet_map(path_pts, legs, method, start_dt):
|
|
| 875 |
st = _state_at_dt(path_pts[0][0], path_pts[0][1], legs, method, start_dt, dt_k)
|
| 876 |
if st is None:
|
| 877 |
continue
|
|
|
|
| 878 |
lat = st["lat"]
|
| 879 |
lon = st["lon"]
|
| 880 |
marker_pts_for_bounds.append((lat, lon))
|
|
@@ -885,8 +1033,10 @@ def build_leaflet_map(path_pts, legs, method, start_dt):
|
|
| 885 |
spd = st["speed_kt"]
|
| 886 |
dist_km = st["dist_km"]
|
| 887 |
|
|
|
|
|
|
|
| 888 |
popup_lines = [
|
| 889 |
-
f"<b>Validation point</b>",
|
| 890 |
f"UTC time: {tlabel}",
|
| 891 |
f"Lat: {lat:.5f}",
|
| 892 |
f"Lon: {lon:.5f}",
|
|
@@ -896,16 +1046,29 @@ def build_leaflet_map(path_pts, legs, method, start_dt):
|
|
| 896 |
f"Distance flown: {dist_km:.2f} km",
|
| 897 |
]
|
| 898 |
|
| 899 |
-
|
| 900 |
-
|
| 901 |
-
|
| 902 |
-
|
| 903 |
-
|
| 904 |
-
|
| 905 |
-
|
| 906 |
-
|
| 907 |
-
|
| 908 |
-
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 909 |
|
| 910 |
# Permanent labels
|
| 911 |
for name, lat, lon in WAYPOINTS:
|
|
@@ -987,125 +1150,30 @@ def make_legs_from_schedule(headings, speed_knots, schedule_minutes):
|
|
| 987 |
return legs
|
| 988 |
|
| 989 |
|
| 990 |
-
# ================= Minimal Scoring (kept from your previous working builds) =================
|
| 991 |
-
def _ping_points_after(start_dt: datetime):
|
| 992 |
-
pts = []
|
| 993 |
-
for row in fp.satLoc:
|
| 994 |
-
label = fp._extract_time_string(row[0])
|
| 995 |
-
meas_bfo = row[SATLOC_IDX_MEAS_BFO]
|
| 996 |
-
meas_bto = row[SATLOC_IDX_MEAS_BTO]
|
| 997 |
-
if meas_bfo is None:
|
| 998 |
-
continue
|
| 999 |
-
try:
|
| 1000 |
-
dt_ping = fp._time_label_to_dt(label, fp.SAT_EPOCH_UTC)
|
| 1001 |
-
except Exception:
|
| 1002 |
-
continue
|
| 1003 |
-
if dt_ping < start_dt:
|
| 1004 |
-
continue
|
| 1005 |
-
pts.append((label, dt_ping, float(meas_bfo), None if meas_bto is None else float(meas_bto)))
|
| 1006 |
-
return pts
|
| 1007 |
-
|
| 1008 |
-
|
| 1009 |
-
def score_legs_bfo_bto(
|
| 1010 |
-
start_lat, start_lon, method, start_dt, end_dt, legs,
|
| 1011 |
-
alt_ft, vs_fpm, descent_start_dt,
|
| 1012 |
-
weight_bfo: float, weight_bto: float
|
| 1013 |
-
):
|
| 1014 |
-
pts = _ping_points_after(start_dt)
|
| 1015 |
-
pts = [(lbl, dt, bfo, bto) for (lbl, dt, bfo, bto) in pts if start_dt <= dt <= end_dt]
|
| 1016 |
-
if len(pts) < 2:
|
| 1017 |
-
return None
|
| 1018 |
-
|
| 1019 |
-
leg_secs = [float(m) * 60.0 for _, _, m in legs]
|
| 1020 |
-
boundaries = np.cumsum([0.0] + leg_secs)
|
| 1021 |
-
total_end = boundaries[-1]
|
| 1022 |
-
|
| 1023 |
-
leg_start_positions = _precompute_leg_starts(start_lat, start_lon, legs, method)
|
| 1024 |
-
|
| 1025 |
-
bfo_preds, bfo_meas = [], []
|
| 1026 |
-
bto_preds, bto_meas = [], []
|
| 1027 |
-
|
| 1028 |
-
for _lbl, dt_ping, meas_bfo, meas_bto in pts:
|
| 1029 |
-
t_sec = (dt_ping - start_dt).total_seconds()
|
| 1030 |
-
if t_sec < 0 or t_sec > total_end:
|
| 1031 |
-
continue
|
| 1032 |
-
|
| 1033 |
-
i = int(np.searchsorted(boundaries, t_sec, side="right") - 1)
|
| 1034 |
-
i = max(0, min(i, len(legs) - 1))
|
| 1035 |
-
t_inside = t_sec - boundaries[i]
|
| 1036 |
-
|
| 1037 |
-
leg_start_lat, leg_start_lon = leg_start_positions[i]
|
| 1038 |
-
hdg, spd, _ = legs[i]
|
| 1039 |
-
lat, lon = _predict_state_from_leg_start(method, leg_start_lat, leg_start_lon, hdg, spd, t_inside)
|
| 1040 |
-
|
| 1041 |
-
alt_now_ft = altitude_ft_at_dt(dt_ping, start_dt, float(alt_ft), float(vs_fpm), descent_start_dt)
|
| 1042 |
-
alt_m = alt_now_ft * 0.3048
|
| 1043 |
-
|
| 1044 |
-
vs_mps = float(vs_fpm) * 0.00508
|
| 1045 |
-
if descent_start_dt is not None and dt_ping < descent_start_dt:
|
| 1046 |
-
vs_mps = 0.0
|
| 1047 |
-
|
| 1048 |
-
spd_mps = float(spd) * 1852.0 / 3600.0
|
| 1049 |
-
pred_bfo = float(fp.calculateBFO_at_time(spd_mps, float(hdg), float(lat), float(lon), alt_m, vs_mps, dt_ping))
|
| 1050 |
-
bfo_preds.append(pred_bfo)
|
| 1051 |
-
bfo_meas.append(float(meas_bfo))
|
| 1052 |
-
|
| 1053 |
-
if meas_bto is not None:
|
| 1054 |
-
pred_bto = float(fp.calculateBTO_at_time(float(lat), float(lon), alt_m, dt_ping))
|
| 1055 |
-
bto_preds.append(pred_bto)
|
| 1056 |
-
bto_meas.append(float(meas_bto))
|
| 1057 |
-
|
| 1058 |
-
if len(bfo_meas) < 2:
|
| 1059 |
-
return None
|
| 1060 |
-
|
| 1061 |
-
e_bfo = np.array(bfo_preds) - np.array(bfo_meas)
|
| 1062 |
-
rmse_bfo = float(np.sqrt(np.mean(e_bfo * e_bfo)))
|
| 1063 |
-
|
| 1064 |
-
rmse_bto = None
|
| 1065 |
-
if len(bto_meas) >= 2:
|
| 1066 |
-
e_bto = np.array(bto_preds) - np.array(bto_meas)
|
| 1067 |
-
rmse_bto = float(np.sqrt(np.mean(e_bto * e_bto)))
|
| 1068 |
-
|
| 1069 |
-
combined = float(weight_bfo) * (rmse_bfo / BFO_SCALE_HZ)
|
| 1070 |
-
if rmse_bto is not None:
|
| 1071 |
-
combined += float(weight_bto) * (rmse_bto / BTO_SCALE_US)
|
| 1072 |
-
|
| 1073 |
-
return {
|
| 1074 |
-
"n_bfo": len(bfo_meas),
|
| 1075 |
-
"n_bto": len(bto_meas),
|
| 1076 |
-
"rmse_bfo": rmse_bfo,
|
| 1077 |
-
"rmse_bto": rmse_bto,
|
| 1078 |
-
"combined": combined,
|
| 1079 |
-
}
|
| 1080 |
-
|
| 1081 |
-
|
| 1082 |
# ================= Guided Solvers =================
|
| 1083 |
def _guided_start_dt():
|
| 1084 |
return RADAR_FIX_TIME_UTC
|
| 1085 |
|
| 1086 |
|
| 1087 |
-
def _guided_end_dt():
|
| 1088 |
-
total_min = TURN_WINDOW_MIN + float(np.sum(GUIDED_LEG_SCHEDULE_MIN))
|
| 1089 |
-
return RADAR_FIX_TIME_UTC + timedelta(minutes=total_min)
|
| 1090 |
-
|
| 1091 |
-
|
| 1092 |
def solve_guided_straight(
|
| 1093 |
method, speed_knots, alt_ft, vs_fpm,
|
| 1094 |
descent_start_choice, weight_bfo, weight_bto
|
| 1095 |
):
|
|
|
|
| 1096 |
start_dt = _guided_start_dt()
|
| 1097 |
-
end_dt =
|
| 1098 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1099 |
descent_start_dt = None
|
| 1100 |
if (descent_start_choice or "").startswith("After Ping Ring 6"):
|
| 1101 |
descent_start_dt = dt_ping6_utc()
|
| 1102 |
elif (descent_start_choice or "").startswith("After Ping Ring 7"):
|
| 1103 |
descent_start_dt = dt_ping7_utc()
|
| 1104 |
|
| 1105 |
-
best = None
|
| 1106 |
-
best_h = None
|
| 1107 |
-
best_full_legs = None
|
| 1108 |
-
|
| 1109 |
h = GUIDED_HEADING_MIN
|
| 1110 |
while h <= GUIDED_HEADING_MAX + 1e-9:
|
| 1111 |
southern_legs = make_legs_from_schedule(h, float(speed_knots), GUIDED_LEG_SCHEDULE_MIN)
|
|
@@ -1116,27 +1184,32 @@ def solve_guided_straight(
|
|
| 1116 |
)
|
| 1117 |
full_legs = preturn_legs + southern_legs
|
| 1118 |
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1119 |
s = score_legs_bfo_bto(
|
| 1120 |
RADAR_FIX_LAT, RADAR_FIX_LON, method, start_dt, end_dt, full_legs,
|
| 1121 |
alt_ft, vs_fpm, descent_start_dt,
|
| 1122 |
weight_bfo, weight_bto
|
| 1123 |
)
|
| 1124 |
-
if s is not None and (
|
| 1125 |
-
|
| 1126 |
best_h = float(h)
|
| 1127 |
best_full_legs = full_legs
|
| 1128 |
|
| 1129 |
h += STRAIGHT_STEP_DEG
|
| 1130 |
|
| 1131 |
-
if
|
| 1132 |
return None, "Guided solve: not enough ping points to score.", None
|
| 1133 |
|
| 1134 |
report = "\n".join([
|
| 1135 |
-
"**Guided mode: Straight + computed turn (18:22→18:40)**",
|
| 1136 |
f"- Speed (constrained): **{float(speed_knots):.1f} kt**",
|
| 1137 |
f"- Best southern heading: **{best_h:.1f}°**",
|
| 1138 |
-
f"-
|
| 1139 |
-
f"-
|
|
|
|
| 1140 |
])
|
| 1141 |
|
| 1142 |
return best_full_legs, report, start_dt
|
|
@@ -1147,67 +1220,168 @@ def solve_guided_straight_with_drift(
|
|
| 1147 |
descent_start_choice, weight_bfo, weight_bto,
|
| 1148 |
drift_max_deg, drift_weight
|
| 1149 |
):
|
| 1150 |
-
#
|
| 1151 |
-
# - drift applies by perturbing the southern heading by +/- drift_max around nominal and picking best.
|
| 1152 |
start_dt = _guided_start_dt()
|
| 1153 |
-
end_dt =
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1154 |
|
|
|
|
| 1155 |
descent_start_dt = None
|
| 1156 |
if (descent_start_choice or "").startswith("After Ping Ring 6"):
|
| 1157 |
descent_start_dt = dt_ping6_utc()
|
| 1158 |
elif (descent_start_choice or "").startswith("After Ping Ring 7"):
|
| 1159 |
descent_start_dt = dt_ping7_utc()
|
| 1160 |
|
| 1161 |
-
|
| 1162 |
-
|
| 1163 |
-
|
|
|
|
| 1164 |
|
| 1165 |
h0 = GUIDED_HEADING_MIN
|
| 1166 |
while h0 <= GUIDED_HEADING_MAX + 1e-9:
|
| 1167 |
-
|
| 1168 |
-
|
| 1169 |
-
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1170 |
|
| 1171 |
-
|
| 1172 |
-
|
| 1173 |
-
speed_knots=float(speed_knots),
|
| 1174 |
-
method=method,
|
| 1175 |
-
heading_final_deg=float(h_eff),
|
| 1176 |
-
)
|
| 1177 |
-
full_legs = preturn_legs + southern_legs
|
| 1178 |
|
| 1179 |
-
|
| 1180 |
-
|
| 1181 |
-
|
| 1182 |
-
|
| 1183 |
-
|
| 1184 |
-
|
| 1185 |
-
best = s
|
| 1186 |
-
best_h0 = float(h0)
|
| 1187 |
-
best_full_legs = full_legs
|
| 1188 |
|
| 1189 |
h0 += 1.0
|
| 1190 |
|
| 1191 |
-
if
|
| 1192 |
return None, "Guided drift solve: not enough ping points to score.", None
|
| 1193 |
|
| 1194 |
report = "\n".join([
|
| 1195 |
-
"**Guided mode: Straight + drift + computed turn (18:22→18:40)**",
|
| 1196 |
f"- Speed (constrained): **{float(speed_knots):.1f} kt**",
|
| 1197 |
f"- Best nominal heading: **{best_h0:.1f}°**",
|
| 1198 |
-
f"-
|
| 1199 |
-
f"-
|
|
|
|
|
|
|
| 1200 |
])
|
| 1201 |
|
| 1202 |
return best_full_legs, report, start_dt
|
| 1203 |
|
| 1204 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1205 |
# ================= RUNNERS =================
|
| 1206 |
def _run_common(start_lat, start_lon, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice):
|
| 1207 |
-
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1208 |
|
| 1209 |
html_map = build_leaflet_map(path_pts, legs, method, start_dt)
|
| 1210 |
-
stats_md = compute_stats(path_pts,
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1211 |
|
| 1212 |
descent_start_dt = None
|
| 1213 |
if (descent_start_choice or "").startswith("After Ping Ring 6"):
|
|
@@ -1220,17 +1394,13 @@ def _run_common(start_lat, start_lon, legs, method, start_dt, alt_ft, vs_fpm, de
|
|
| 1220 |
start_lat, start_lon, method, start_dt, legs, alt_ft, vs_fpm, descent_start_choice
|
| 1221 |
)
|
| 1222 |
|
| 1223 |
-
# lightweight empty BFO plot placeholder (your earlier full version had full BFO graph)
|
| 1224 |
-
bfo_fig = go.Figure()
|
| 1225 |
-
bfo_fig.update_layout(title="BFO Analysis", xaxis_title="Time", yaxis_title="Hz", height=520)
|
| 1226 |
-
bfo_md = ""
|
| 1227 |
-
|
| 1228 |
return html_map, stats_md, bfo_fig, bfo_md, alt_fig, headers, rows
|
| 1229 |
|
| 1230 |
|
| 1231 |
def run_sim_manual(start_lat, start_lon, legs_df, method, start_time_str, alt_ft, vs_fpm, descent_start_choice):
|
| 1232 |
legs = parse_legs_df(legs_df)
|
| 1233 |
start_dt = _parse_start_time(start_time_str)
|
|
|
|
| 1234 |
html_map, stats_md, bfo_fig, bfo_md, alt_fig, headers, rows = _run_common(
|
| 1235 |
start_lat, start_lon, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice
|
| 1236 |
)
|
|
@@ -1271,6 +1441,7 @@ def run_guided(
|
|
| 1271 |
], value=[])
|
| 1272 |
return gr.update(), f"### ❌ {report_md}", "", "", fig, "", alt_fig, empty_table
|
| 1273 |
|
|
|
|
| 1274 |
html_map, stats_md, bfo_fig, bfo_md, alt_fig, headers, rows = _run_common(
|
| 1275 |
RADAR_FIX_LAT, RADAR_FIX_LON, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice
|
| 1276 |
)
|
|
@@ -1301,7 +1472,7 @@ with gr.Blocks() as demo:
|
|
| 1301 |
f"- **Last radar spot (fixed):** 18:22:00Z @ **{RADAR_FIX_LAT:.3f}N, {RADAR_FIX_LON:.3f}E**\n"
|
| 1302 |
f"- **Turn completion (fixed):** 18:40:00Z\n"
|
| 1303 |
f"- **Initial track (fixed):** {RADAR_INITIAL_TRACK_DEG:.1f}°\n"
|
| 1304 |
-
f"- **Final resting point
|
| 1305 |
)
|
| 1306 |
|
| 1307 |
with gr.Row():
|
|
|
|
| 11 |
# - Guided mode computes a realistic turn between 18:22→18:40 and then southern path.
|
| 12 |
# - Guided validation timestamps: 18:22, 18:40, 19:40, 20:40, 21:40, 22:40, 00:10, 00:20
|
| 13 |
#
|
| 14 |
+
# FIXES IMPLEMENTED IN THIS VERSION:
|
| 15 |
+
# 1) Remove trajectory (red line) AFTER 00:20
|
| 16 |
+
# - In guided mode, we CLIP the simulated legs so the polyline ends exactly at 00:20Z.
|
| 17 |
+
# 2) Restore BFO plot (it disappeared previously)
|
| 18 |
+
# - Re-added full BFO analysis computation + chart + MAE/RMSE markdown.
|
| 19 |
|
| 20 |
import os
|
| 21 |
import re
|
|
|
|
| 63 |
"Ping Ring 7": "24:10:00Z",
|
| 64 |
}
|
| 65 |
|
| 66 |
+
# Normalization (used for combining BFO + BTO into a single score in guided solve)
|
| 67 |
BFO_SCALE_HZ = 5.0
|
| 68 |
BTO_SCALE_US = 50.0
|
| 69 |
|
| 70 |
+
# EXPLICIT WEIGHTS
|
| 71 |
DEFAULT_WEIGHT_BFO = 1.0
|
| 72 |
DEFAULT_WEIGHT_BTO = 2.0
|
| 73 |
|
| 74 |
+
# satLoc indices
|
| 75 |
SATLOC_IDX_MEAS_BFO = 4
|
| 76 |
SATLOC_IDX_MEAS_BTO = 5
|
| 77 |
|
|
|
|
| 83 |
TURN_COMPLETE_TIME_UTC = datetime(2014, 3, 7, 18, 40, 0, tzinfo=timezone.utc)
|
| 84 |
TURN_WINDOW_MIN = (TURN_COMPLETE_TIME_UTC - RADAR_FIX_TIME_UTC).total_seconds() / 60.0
|
| 85 |
|
| 86 |
+
# Initial track shown in your radar annotation
|
| 87 |
RADAR_INITIAL_TRACK_DEG = 296.0
|
| 88 |
|
| 89 |
# ===== Guided validation timestamps of interest =====
|
| 90 |
+
# (00:10 and 00:20 are represented as 24:10 and 24:20 in your SAT epoch helper)
|
| 91 |
GUIDED_VALIDATION_LABELS = [
|
| 92 |
"18:22:00",
|
| 93 |
"18:40:00",
|
|
|
|
| 98 |
"24:10:00", # 00:10Z
|
| 99 |
"24:20:00", # 00:20Z
|
| 100 |
]
|
| 101 |
+
GUIDED_FINAL_LABEL = "24:20:00" # 00:20Z FINAL RESTING POINT
|
| 102 |
|
| 103 |
# Turn modeling parameters
|
| 104 |
TURN_BANK_MIN_DEG = 15.0
|
|
|
|
| 106 |
TURN_BANK_DEFAULT_DEG = 25.0
|
| 107 |
TURN_STEP_SEC = 30.0
|
| 108 |
|
| 109 |
+
# Guided heading search (south-ish)
|
| 110 |
GUIDED_HEADING_MIN = 150.0
|
| 111 |
GUIDED_HEADING_MAX = 210.0
|
| 112 |
STRAIGHT_STEP_DEG = 0.5
|
| 113 |
|
| 114 |
+
# Southern schedule (kept)
|
| 115 |
+
GUIDED_LEG_SCHEDULE_MIN = [12.0, 60.0, 60.0, 60.0, 60.0, 90.0, 10.0]
|
| 116 |
|
| 117 |
# Drift
|
| 118 |
DEFAULT_DRIFT_MAX_DEG = 2.0
|
|
|
|
| 419 |
|
| 420 |
legs = []
|
| 421 |
|
| 422 |
+
# initial straight on radar track
|
| 423 |
if straight_sec > 1.0:
|
| 424 |
legs.append((h0, float(speed_knots), straight_sec / 60.0))
|
| 425 |
|
| 426 |
+
# discretized turn segments
|
| 427 |
if turn_sec > 1.0 and delta_abs > 1e-6:
|
| 428 |
step_sec = max(5.0, float(step_sec))
|
| 429 |
n = int(math.ceil(turn_sec / step_sec))
|
|
|
|
| 449 |
return legs, info
|
| 450 |
|
| 451 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 452 |
# ================= TIME HELPERS =================
|
| 453 |
def _parse_start_time(s: str) -> datetime:
|
| 454 |
s = (s or "").strip()
|
|
|
|
| 522 |
return fig
|
| 523 |
|
| 524 |
|
| 525 |
+
# ================= SIM CORE =================
|
| 526 |
+
def simulate_path(start_lat, start_lon, legs, method):
|
| 527 |
+
lat, lon = float(start_lat), float(start_lon)
|
| 528 |
+
pts = [(lat, lon)]
|
| 529 |
+
seg_dists_m = []
|
| 530 |
+
total_time_min = 0.0
|
| 531 |
+
|
| 532 |
+
for heading, speed_knots, minutes in legs:
|
| 533 |
+
dist_m = float(speed_knots) * 1852.0 * (float(minutes) / 60.0)
|
| 534 |
+
if method == "VINCENTY":
|
| 535 |
+
(lat, lon), _ = fp.vincenty_direct((lat, lon), float(heading), dist_m)
|
| 536 |
+
else:
|
| 537 |
+
lat, lon = fp.rhumb_direct((lat, lon), float(heading), dist_m)
|
| 538 |
+
|
| 539 |
+
pts.append((lat, lon))
|
| 540 |
+
seg_dists_m.append(dist_m)
|
| 541 |
+
total_time_min += float(minutes)
|
| 542 |
+
|
| 543 |
+
return pts, seg_dists_m, total_time_min
|
| 544 |
+
|
| 545 |
+
|
| 546 |
+
def compute_stats(path_pts, seg_dists_m, total_time_min, legs):
|
| 547 |
+
total_dist_m = float(np.sum(seg_dists_m)) if seg_dists_m else 0.0
|
| 548 |
+
total_dist_km = total_dist_m / 1000.0
|
| 549 |
+
total_dist_nm = total_dist_m / 1852.0
|
| 550 |
+
|
| 551 |
+
total_time_hr = total_time_min / 60.0
|
| 552 |
+
avg_gs_knots = (total_dist_nm / total_time_hr) if total_time_hr > 0 else 0.0
|
| 553 |
+
|
| 554 |
+
end_lat, end_lon = path_pts[-1]
|
| 555 |
+
return "\n\n".join([
|
| 556 |
+
f"**Legs:** {len(legs)}",
|
| 557 |
+
f"**Total time:** {total_time_min:.1f} min ({total_time_hr:.2f} hr)",
|
| 558 |
+
f"**Total distance:** {total_dist_nm:.1f} nm ({total_dist_km:.1f} km)",
|
| 559 |
+
f"**Average ground speed:** {avg_gs_knots:.1f} kt",
|
| 560 |
+
f"**Final position (polyline end):** lat {end_lat:.4f}, lon {end_lon:.4f}",
|
| 561 |
+
])
|
| 562 |
+
|
| 563 |
+
|
| 564 |
+
# ================= EXACT STATE HELPERS =================
|
| 565 |
def _predict_state_from_leg_start(method, leg_start_lat, leg_start_lon, hdg_deg, spd_knots, t_inside_leg_sec):
|
| 566 |
dist_m = float(spd_knots) * 1852.0 * (float(t_inside_leg_sec) / 3600.0)
|
| 567 |
if method == "VINCENTY":
|
|
|
|
| 637 |
return times
|
| 638 |
|
| 639 |
|
| 640 |
+
# ================= CLIP TRAJECTORY TO 00:20 (GUIDED) =================
|
| 641 |
+
def clip_legs_to_end_dt(legs, start_dt: datetime, end_dt: datetime):
|
| 642 |
+
"""
|
| 643 |
+
Returns a new legs list truncated so total duration ends exactly at end_dt.
|
| 644 |
+
If legs already shorter than end_dt-start_dt, returns legs unchanged.
|
| 645 |
+
"""
|
| 646 |
+
total_target_sec = (end_dt - start_dt).total_seconds()
|
| 647 |
+
if total_target_sec <= 0:
|
| 648 |
+
return []
|
| 649 |
+
|
| 650 |
+
out = []
|
| 651 |
+
remaining = float(total_target_sec)
|
| 652 |
+
|
| 653 |
+
for (hdg, spd, minutes) in legs:
|
| 654 |
+
sec = float(minutes) * 60.0
|
| 655 |
+
if remaining <= 0:
|
| 656 |
+
break
|
| 657 |
+
if sec <= remaining + 1e-9:
|
| 658 |
+
out.append((float(hdg), float(spd), float(minutes)))
|
| 659 |
+
remaining -= sec
|
| 660 |
+
else:
|
| 661 |
+
# truncate this leg
|
| 662 |
+
trunc_min = max(0.0, remaining / 60.0)
|
| 663 |
+
if trunc_min > 1e-6:
|
| 664 |
+
out.append((float(hdg), float(spd), float(trunc_min)))
|
| 665 |
+
remaining = 0.0
|
| 666 |
+
break
|
| 667 |
+
|
| 668 |
+
return out
|
| 669 |
+
|
| 670 |
+
|
| 671 |
+
# ================= BFO ANALYSIS (RESTORED) =================
|
| 672 |
+
def _ping_points_after(start_dt: datetime):
|
| 673 |
+
pts = []
|
| 674 |
+
for row in fp.satLoc:
|
| 675 |
+
label = fp._extract_time_string(row[0])
|
| 676 |
+
meas_bfo = row[SATLOC_IDX_MEAS_BFO]
|
| 677 |
+
meas_bto = row[SATLOC_IDX_MEAS_BTO]
|
| 678 |
+
if meas_bfo is None:
|
| 679 |
+
continue
|
| 680 |
+
try:
|
| 681 |
+
dt_ping = fp._time_label_to_dt(label, fp.SAT_EPOCH_UTC)
|
| 682 |
+
except Exception:
|
| 683 |
+
continue
|
| 684 |
+
if dt_ping < start_dt:
|
| 685 |
+
continue
|
| 686 |
+
pts.append((label, dt_ping, float(meas_bfo), None if meas_bto is None else float(meas_bto)))
|
| 687 |
+
return pts
|
| 688 |
+
|
| 689 |
+
|
| 690 |
+
def bfo_analysis_with_legs(start_lat, start_lon, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice):
|
| 691 |
+
pts = _ping_points_after(start_dt)
|
| 692 |
+
|
| 693 |
+
descent_start_dt = None
|
| 694 |
+
if (descent_start_choice or "").startswith("After Ping Ring 6"):
|
| 695 |
+
descent_start_dt = dt_ping6_utc()
|
| 696 |
+
elif (descent_start_choice or "").startswith("After Ping Ring 7"):
|
| 697 |
+
descent_start_dt = dt_ping7_utc()
|
| 698 |
+
|
| 699 |
+
leg_secs = [float(m) * 60.0 for _, _, m in legs]
|
| 700 |
+
boundaries = np.cumsum([0.0] + leg_secs)
|
| 701 |
+
total_end = float(boundaries[-1])
|
| 702 |
+
|
| 703 |
+
leg_start_positions = _precompute_leg_starts(start_lat, start_lon, legs, method)
|
| 704 |
+
|
| 705 |
+
x_labels, measured, calculated = [], [], []
|
| 706 |
+
|
| 707 |
+
for label, dt_ping, meas_bfo, _meas_bto in pts:
|
| 708 |
+
t_sec = (dt_ping - start_dt).total_seconds()
|
| 709 |
+
if t_sec < 0 or t_sec > total_end:
|
| 710 |
+
continue
|
| 711 |
+
|
| 712 |
+
i = int(np.searchsorted(boundaries, t_sec, side="right") - 1)
|
| 713 |
+
i = max(0, min(i, len(legs) - 1))
|
| 714 |
+
t_inside = float(t_sec) - float(boundaries[i])
|
| 715 |
+
|
| 716 |
+
leg_start_lat, leg_start_lon = leg_start_positions[i]
|
| 717 |
+
hdg, spd, _m = legs[i]
|
| 718 |
+
|
| 719 |
+
lat, lon = _predict_state_from_leg_start(method, leg_start_lat, leg_start_lon, hdg, spd, t_inside)
|
| 720 |
+
|
| 721 |
+
alt_now_ft = altitude_ft_at_dt(dt_ping, start_dt, float(alt_ft), float(vs_fpm), descent_start_dt)
|
| 722 |
+
alt_m = alt_now_ft * 0.3048
|
| 723 |
+
|
| 724 |
+
vs_mps = float(vs_fpm) * 0.00508
|
| 725 |
+
if descent_start_dt is not None and dt_ping < descent_start_dt:
|
| 726 |
+
vs_mps = 0.0
|
| 727 |
+
|
| 728 |
+
spd_mps = float(spd) * 1852.0 / 3600.0
|
| 729 |
+
pred_bfo = float(fp.calculateBFO_at_time(spd_mps, float(hdg), float(lat), float(lon), alt_m, vs_mps, dt_ping))
|
| 730 |
+
|
| 731 |
+
x_labels.append(label)
|
| 732 |
+
measured.append(float(meas_bfo))
|
| 733 |
+
calculated.append(float(pred_bfo))
|
| 734 |
+
|
| 735 |
+
fig = go.Figure()
|
| 736 |
+
if not x_labels:
|
| 737 |
+
fig.update_layout(title="BFO Analysis", xaxis_title="Time", yaxis_title="Hz", height=520)
|
| 738 |
+
return fig, "**BFO scoring:** no valid ping points inside your simulated time window."
|
| 739 |
+
|
| 740 |
+
m_arr = np.array(measured, dtype=float)
|
| 741 |
+
c_arr = np.array(calculated, dtype=float)
|
| 742 |
+
err = c_arr - m_arr
|
| 743 |
+
mae = float(np.mean(np.abs(err)))
|
| 744 |
+
rmse = float(np.sqrt(np.mean(err * err)))
|
| 745 |
+
|
| 746 |
+
fig.add_trace(go.Scatter(
|
| 747 |
+
x=x_labels, y=measured,
|
| 748 |
+
mode="markers",
|
| 749 |
+
name="Measured BFO",
|
| 750 |
+
marker=dict(size=11, symbol="circle")
|
| 751 |
+
))
|
| 752 |
+
fig.add_trace(go.Scatter(
|
| 753 |
+
x=x_labels, y=calculated,
|
| 754 |
+
mode="lines+markers",
|
| 755 |
+
name="Calculated BFO",
|
| 756 |
+
line=dict(width=3),
|
| 757 |
+
marker=dict(size=10, symbol="diamond")
|
| 758 |
+
))
|
| 759 |
+
fig.update_layout(
|
| 760 |
+
title="BFO Analysis",
|
| 761 |
+
xaxis_title="Time",
|
| 762 |
+
yaxis_title="Hz",
|
| 763 |
+
legend=dict(x=1.02, y=1.0),
|
| 764 |
+
margin=dict(l=60, r=160, t=80, b=60),
|
| 765 |
+
height=520,
|
| 766 |
+
)
|
| 767 |
+
|
| 768 |
+
score_md = (
|
| 769 |
+
f"**BFO scoring (Hz)**\n\n"
|
| 770 |
+
f"- Points: **{len(m_arr)}**\n"
|
| 771 |
+
f"- MAE: **{mae:.2f} Hz**\n"
|
| 772 |
+
f"- RMSE: **{rmse:.2f} Hz**"
|
| 773 |
+
)
|
| 774 |
+
return fig, score_md
|
| 775 |
+
|
| 776 |
+
|
| 777 |
# ================= VALIDATION TABLE =================
|
| 778 |
def _lookup_measured_at_dt(dt_ping: datetime):
|
| 779 |
if hasattr(fp, "nearest_measured_for_time"):
|
|
|
|
| 899 |
- Ping rings always drawn
|
| 900 |
- Search areas always drawn (color-coded)
|
| 901 |
- Validation markers are the authoritative markers
|
| 902 |
+
- Guided mode does NOT draw a purple End marker
|
| 903 |
+
- Trajectory is clipped to 00:20 in guided mode by clipping legs before simulate_path
|
| 904 |
"""
|
| 905 |
start_lat, start_lon = path_pts[0]
|
| 906 |
is_guided = (start_dt == RADAR_FIX_TIME_UTC)
|
|
|
|
| 983 |
tooltip = f"{ring_name} ({ring_time})" if ring_time else ring_name
|
| 984 |
folium.PolyLine(ring_pts, weight=2, opacity=0.6, color="yellow", tooltip=tooltip).add_to(m)
|
| 985 |
|
| 986 |
+
# Flight path (polyline)
|
| 987 |
folium.PolyLine(path_pts, color="white", weight=8, opacity=0.3).add_to(m)
|
| 988 |
folium.PolyLine(path_pts, color="red", weight=4, opacity=1.0, tooltip="Flight path").add_to(m)
|
| 989 |
|
|
|
|
| 1000 |
),
|
| 1001 |
).add_to(m)
|
| 1002 |
|
| 1003 |
+
# Manual mode can show an "End" marker. Guided mode: no extra end marker.
|
|
|
|
| 1004 |
if not is_guided:
|
| 1005 |
folium.CircleMarker(
|
| 1006 |
location=path_pts[-1],
|
|
|
|
| 1022 |
st = _state_at_dt(path_pts[0][0], path_pts[0][1], legs, method, start_dt, dt_k)
|
| 1023 |
if st is None:
|
| 1024 |
continue
|
| 1025 |
+
|
| 1026 |
lat = st["lat"]
|
| 1027 |
lon = st["lon"]
|
| 1028 |
marker_pts_for_bounds.append((lat, lon))
|
|
|
|
| 1033 |
spd = st["speed_kt"]
|
| 1034 |
dist_km = st["dist_km"]
|
| 1035 |
|
| 1036 |
+
is_final = (start_dt == RADAR_FIX_TIME_UTC) and (dt_k == _dt_from_label(GUIDED_FINAL_LABEL))
|
| 1037 |
+
|
| 1038 |
popup_lines = [
|
| 1039 |
+
f"<b>{'FINAL RESTING POINT' if is_final else 'Validation point'}</b>",
|
| 1040 |
f"UTC time: {tlabel}",
|
| 1041 |
f"Lat: {lat:.5f}",
|
| 1042 |
f"Lon: {lon:.5f}",
|
|
|
|
| 1046 |
f"Distance flown: {dist_km:.2f} km",
|
| 1047 |
]
|
| 1048 |
|
| 1049 |
+
# Make final point visually distinct
|
| 1050 |
+
if is_final:
|
| 1051 |
+
folium.CircleMarker(
|
| 1052 |
+
location=(lat, lon),
|
| 1053 |
+
radius=9,
|
| 1054 |
+
fill=True,
|
| 1055 |
+
color="magenta",
|
| 1056 |
+
fill_color="magenta",
|
| 1057 |
+
fill_opacity=0.95,
|
| 1058 |
+
tooltip=f"{tlabel} (FINAL)",
|
| 1059 |
+
popup=folium.Popup("<br>".join(popup_lines), max_width=360),
|
| 1060 |
+
).add_to(m)
|
| 1061 |
+
else:
|
| 1062 |
+
folium.CircleMarker(
|
| 1063 |
+
location=(lat, lon),
|
| 1064 |
+
radius=6,
|
| 1065 |
+
fill=True,
|
| 1066 |
+
color="white",
|
| 1067 |
+
fill_color="white",
|
| 1068 |
+
fill_opacity=1.0,
|
| 1069 |
+
tooltip=tlabel,
|
| 1070 |
+
popup=folium.Popup("<br>".join(popup_lines), max_width=360),
|
| 1071 |
+
).add_to(m)
|
| 1072 |
|
| 1073 |
# Permanent labels
|
| 1074 |
for name, lat, lon in WAYPOINTS:
|
|
|
|
| 1150 |
return legs
|
| 1151 |
|
| 1152 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1153 |
# ================= Guided Solvers =================
|
| 1154 |
def _guided_start_dt():
|
| 1155 |
return RADAR_FIX_TIME_UTC
|
| 1156 |
|
| 1157 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1158 |
def solve_guided_straight(
|
| 1159 |
method, speed_knots, alt_ft, vs_fpm,
|
| 1160 |
descent_start_choice, weight_bfo, weight_bto
|
| 1161 |
):
|
| 1162 |
+
# We clip everything to final 00:20Z; end_dt is fixed by label.
|
| 1163 |
start_dt = _guided_start_dt()
|
| 1164 |
+
end_dt = _dt_from_label(GUIDED_FINAL_LABEL)
|
| 1165 |
|
| 1166 |
+
best_h = None
|
| 1167 |
+
best_full_legs = None
|
| 1168 |
+
best_score = None
|
| 1169 |
+
|
| 1170 |
+
# descent start
|
| 1171 |
descent_start_dt = None
|
| 1172 |
if (descent_start_choice or "").startswith("After Ping Ring 6"):
|
| 1173 |
descent_start_dt = dt_ping6_utc()
|
| 1174 |
elif (descent_start_choice or "").startswith("After Ping Ring 7"):
|
| 1175 |
descent_start_dt = dt_ping7_utc()
|
| 1176 |
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1177 |
h = GUIDED_HEADING_MIN
|
| 1178 |
while h <= GUIDED_HEADING_MAX + 1e-9:
|
| 1179 |
southern_legs = make_legs_from_schedule(h, float(speed_knots), GUIDED_LEG_SCHEDULE_MIN)
|
|
|
|
| 1184 |
)
|
| 1185 |
full_legs = preturn_legs + southern_legs
|
| 1186 |
|
| 1187 |
+
# clip to 00:20 so scoring window and plotting window are consistent
|
| 1188 |
+
full_legs = clip_legs_to_end_dt(full_legs, start_dt, end_dt)
|
| 1189 |
+
|
| 1190 |
+
# score using BFO/BTO inside start_dt→end_dt duration only
|
| 1191 |
s = score_legs_bfo_bto(
|
| 1192 |
RADAR_FIX_LAT, RADAR_FIX_LON, method, start_dt, end_dt, full_legs,
|
| 1193 |
alt_ft, vs_fpm, descent_start_dt,
|
| 1194 |
weight_bfo, weight_bto
|
| 1195 |
)
|
| 1196 |
+
if s is not None and (best_score is None or s["combined"] < best_score["combined"]):
|
| 1197 |
+
best_score = s
|
| 1198 |
best_h = float(h)
|
| 1199 |
best_full_legs = full_legs
|
| 1200 |
|
| 1201 |
h += STRAIGHT_STEP_DEG
|
| 1202 |
|
| 1203 |
+
if best_full_legs is None:
|
| 1204 |
return None, "Guided solve: not enough ping points to score.", None
|
| 1205 |
|
| 1206 |
report = "\n".join([
|
| 1207 |
+
"**Guided mode: Straight + computed turn (18:22→18:40), clipped to 00:20**",
|
| 1208 |
f"- Speed (constrained): **{float(speed_knots):.1f} kt**",
|
| 1209 |
f"- Best southern heading: **{best_h:.1f}°**",
|
| 1210 |
+
f"- BFO points: **{best_score['n_bfo']}** | RMSE_BFO: **{best_score['rmse_bfo']:.2f} Hz**",
|
| 1211 |
+
f"- BTO points: **{best_score['n_bto']}** | RMSE_BTO: **{best_score['rmse_bto']:.1f} µs**" if best_score.get("rmse_bto") is not None else
|
| 1212 |
+
f"- BTO points: **{best_score['n_bto']}** | RMSE_BTO: n/a",
|
| 1213 |
])
|
| 1214 |
|
| 1215 |
return best_full_legs, report, start_dt
|
|
|
|
| 1220 |
descent_start_choice, weight_bfo, weight_bto,
|
| 1221 |
drift_max_deg, drift_weight
|
| 1222 |
):
|
| 1223 |
+
# lightweight drift version: use nominal heading search, but add a small constant drift
|
|
|
|
| 1224 |
start_dt = _guided_start_dt()
|
| 1225 |
+
end_dt = _dt_from_label(GUIDED_FINAL_LABEL)
|
| 1226 |
+
|
| 1227 |
+
best_h0 = None
|
| 1228 |
+
best_full_legs = None
|
| 1229 |
+
best_score = None
|
| 1230 |
|
| 1231 |
+
# descent start
|
| 1232 |
descent_start_dt = None
|
| 1233 |
if (descent_start_choice or "").startswith("After Ping Ring 6"):
|
| 1234 |
descent_start_dt = dt_ping6_utc()
|
| 1235 |
elif (descent_start_choice or "").startswith("After Ping Ring 7"):
|
| 1236 |
descent_start_dt = dt_ping7_utc()
|
| 1237 |
|
| 1238 |
+
# try drift candidates within ±drift_max_deg
|
| 1239 |
+
drift_candidates = list(np.arange(-float(drift_max_deg), float(drift_max_deg) + 1e-9, DRIFT_STEP_DEG))
|
| 1240 |
+
if not drift_candidates:
|
| 1241 |
+
drift_candidates = [0.0]
|
| 1242 |
|
| 1243 |
h0 = GUIDED_HEADING_MIN
|
| 1244 |
while h0 <= GUIDED_HEADING_MAX + 1e-9:
|
| 1245 |
+
for d in drift_candidates:
|
| 1246 |
+
h_eff = _wrap_heading_deg(float(h0) + float(d))
|
| 1247 |
+
|
| 1248 |
+
southern_legs = make_legs_from_schedule(h_eff, float(speed_knots), GUIDED_LEG_SCHEDULE_MIN)
|
| 1249 |
+
preturn_legs, _turn_info = build_turn_legs(
|
| 1250 |
+
speed_knots=float(speed_knots),
|
| 1251 |
+
method=method,
|
| 1252 |
+
heading_final_deg=float(h_eff),
|
| 1253 |
+
)
|
| 1254 |
+
full_legs = preturn_legs + southern_legs
|
| 1255 |
+
|
| 1256 |
+
# clip to 00:20
|
| 1257 |
+
full_legs = clip_legs_to_end_dt(full_legs, start_dt, end_dt)
|
| 1258 |
+
|
| 1259 |
+
s = score_legs_bfo_bto(
|
| 1260 |
+
RADAR_FIX_LAT, RADAR_FIX_LON, method, start_dt, end_dt, full_legs,
|
| 1261 |
+
alt_ft, vs_fpm, descent_start_dt,
|
| 1262 |
+
weight_bfo, weight_bto
|
| 1263 |
+
)
|
| 1264 |
+
if s is None:
|
| 1265 |
+
continue
|
| 1266 |
|
| 1267 |
+
# add drift penalty (RMS here is just |d| for constant drift)
|
| 1268 |
+
combined = float(s["combined"]) + float(drift_weight) * (abs(float(d)) / DRIFT_SCALE_DEG)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1269 |
|
| 1270 |
+
if best_score is None or combined < best_score["combined_plus"]:
|
| 1271 |
+
best_score = dict(s)
|
| 1272 |
+
best_score["combined_plus"] = combined
|
| 1273 |
+
best_score["drift_used"] = float(d)
|
| 1274 |
+
best_h0 = float(h0)
|
| 1275 |
+
best_full_legs = full_legs
|
|
|
|
|
|
|
|
|
|
| 1276 |
|
| 1277 |
h0 += 1.0
|
| 1278 |
|
| 1279 |
+
if best_full_legs is None:
|
| 1280 |
return None, "Guided drift solve: not enough ping points to score.", None
|
| 1281 |
|
| 1282 |
report = "\n".join([
|
| 1283 |
+
"**Guided mode: Straight + drift + computed turn (18:22→18:40), clipped to 00:20**",
|
| 1284 |
f"- Speed (constrained): **{float(speed_knots):.1f} kt**",
|
| 1285 |
f"- Best nominal heading: **{best_h0:.1f}°**",
|
| 1286 |
+
f"- Drift used: **{best_score.get('drift_used', 0.0):+.1f}°** | Drift penalty weight: **{float(drift_weight):.2f}**",
|
| 1287 |
+
f"- BFO points: **{best_score['n_bfo']}** | RMSE_BFO: **{best_score['rmse_bfo']:.2f} Hz**",
|
| 1288 |
+
f"- BTO points: **{best_score['n_bto']}** | RMSE_BTO: **{best_score['rmse_bto']:.1f} µs**" if best_score.get("rmse_bto") is not None else
|
| 1289 |
+
f"- BTO points: **{best_score['n_bto']}** | RMSE_BTO: n/a",
|
| 1290 |
])
|
| 1291 |
|
| 1292 |
return best_full_legs, report, start_dt
|
| 1293 |
|
| 1294 |
|
| 1295 |
+
# ================= Scoring helper (used in guided solvers) =================
|
| 1296 |
+
def score_legs_bfo_bto(
|
| 1297 |
+
start_lat, start_lon, method, start_dt, end_dt, legs,
|
| 1298 |
+
alt_ft, vs_fpm, descent_start_dt,
|
| 1299 |
+
weight_bfo: float, weight_bto: float
|
| 1300 |
+
):
|
| 1301 |
+
pts = _ping_points_after(start_dt)
|
| 1302 |
+
pts = [(lbl, dt, bfo, bto) for (lbl, dt, bfo, bto) in pts if start_dt <= dt <= end_dt]
|
| 1303 |
+
if len(pts) < 2:
|
| 1304 |
+
return None
|
| 1305 |
+
|
| 1306 |
+
leg_secs = [float(m) * 60.0 for _, _, m in legs]
|
| 1307 |
+
boundaries = np.cumsum([0.0] + leg_secs)
|
| 1308 |
+
total_end = float(boundaries[-1])
|
| 1309 |
+
|
| 1310 |
+
leg_start_positions = _precompute_leg_starts(start_lat, start_lon, legs, method)
|
| 1311 |
+
|
| 1312 |
+
bfo_preds, bfo_meas = [], []
|
| 1313 |
+
bto_preds, bto_meas = [], []
|
| 1314 |
+
|
| 1315 |
+
for _lbl, dt_ping, meas_bfo, meas_bto in pts:
|
| 1316 |
+
t_sec = (dt_ping - start_dt).total_seconds()
|
| 1317 |
+
if t_sec < 0 or t_sec > total_end:
|
| 1318 |
+
continue
|
| 1319 |
+
|
| 1320 |
+
i = int(np.searchsorted(boundaries, t_sec, side="right") - 1)
|
| 1321 |
+
i = max(0, min(i, len(legs) - 1))
|
| 1322 |
+
t_inside = float(t_sec) - float(boundaries[i])
|
| 1323 |
+
|
| 1324 |
+
leg_start_lat, leg_start_lon = leg_start_positions[i]
|
| 1325 |
+
hdg, spd, _ = legs[i]
|
| 1326 |
+
lat, lon = _predict_state_from_leg_start(method, leg_start_lat, leg_start_lon, hdg, spd, t_inside)
|
| 1327 |
+
|
| 1328 |
+
alt_now_ft = altitude_ft_at_dt(dt_ping, start_dt, float(alt_ft), float(vs_fpm), descent_start_dt)
|
| 1329 |
+
alt_m = alt_now_ft * 0.3048
|
| 1330 |
+
|
| 1331 |
+
vs_mps = float(vs_fpm) * 0.00508
|
| 1332 |
+
if descent_start_dt is not None and dt_ping < descent_start_dt:
|
| 1333 |
+
vs_mps = 0.0
|
| 1334 |
+
|
| 1335 |
+
spd_mps = float(spd) * 1852.0 / 3600.0
|
| 1336 |
+
pred_bfo = float(fp.calculateBFO_at_time(spd_mps, float(hdg), float(lat), float(lon), alt_m, vs_mps, dt_ping))
|
| 1337 |
+
bfo_preds.append(pred_bfo)
|
| 1338 |
+
bfo_meas.append(float(meas_bfo))
|
| 1339 |
+
|
| 1340 |
+
if meas_bto is not None:
|
| 1341 |
+
pred_bto = float(fp.calculateBTO_at_time(float(lat), float(lon), alt_m, dt_ping))
|
| 1342 |
+
bto_preds.append(pred_bto)
|
| 1343 |
+
bto_meas.append(float(meas_bto))
|
| 1344 |
+
|
| 1345 |
+
if len(bfo_meas) < 2:
|
| 1346 |
+
return None
|
| 1347 |
+
|
| 1348 |
+
e_bfo = np.array(bfo_preds) - np.array(bfo_meas)
|
| 1349 |
+
rmse_bfo = float(np.sqrt(np.mean(e_bfo * e_bfo)))
|
| 1350 |
+
|
| 1351 |
+
rmse_bto = None
|
| 1352 |
+
if len(bto_meas) >= 2:
|
| 1353 |
+
e_bto = np.array(bto_preds) - np.array(bto_meas)
|
| 1354 |
+
rmse_bto = float(np.sqrt(np.mean(e_bto * e_bto)))
|
| 1355 |
+
|
| 1356 |
+
combined = float(weight_bfo) * (rmse_bfo / BFO_SCALE_HZ)
|
| 1357 |
+
if rmse_bto is not None:
|
| 1358 |
+
combined += float(weight_bto) * (rmse_bto / BTO_SCALE_US)
|
| 1359 |
+
|
| 1360 |
+
return {
|
| 1361 |
+
"n_bfo": len(bfo_meas),
|
| 1362 |
+
"n_bto": len(bto_meas),
|
| 1363 |
+
"rmse_bfo": rmse_bfo,
|
| 1364 |
+
"rmse_bto": rmse_bto,
|
| 1365 |
+
"combined": combined,
|
| 1366 |
+
}
|
| 1367 |
+
|
| 1368 |
+
|
| 1369 |
# ================= RUNNERS =================
|
| 1370 |
def _run_common(start_lat, start_lon, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice):
|
| 1371 |
+
# In guided mode, ensure the polyline ends at 00:20 by clipping legs to end_dt first.
|
| 1372 |
+
if start_dt == RADAR_FIX_TIME_UTC:
|
| 1373 |
+
end_dt = _dt_from_label(GUIDED_FINAL_LABEL)
|
| 1374 |
+
legs = clip_legs_to_end_dt(legs, start_dt, end_dt)
|
| 1375 |
+
|
| 1376 |
+
path_pts, seg_dists_m, total_time_min = simulate_path(start_lat, start_lon, legs, method)
|
| 1377 |
|
| 1378 |
html_map = build_leaflet_map(path_pts, legs, method, start_dt)
|
| 1379 |
+
stats_md = compute_stats(path_pts, seg_dists_m, total_time_min, legs)
|
| 1380 |
+
|
| 1381 |
+
# BFO plot (restored)
|
| 1382 |
+
bfo_fig, bfo_md = bfo_analysis_with_legs(
|
| 1383 |
+
start_lat, start_lon, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice
|
| 1384 |
+
)
|
| 1385 |
|
| 1386 |
descent_start_dt = None
|
| 1387 |
if (descent_start_choice or "").startswith("After Ping Ring 6"):
|
|
|
|
| 1394 |
start_lat, start_lon, method, start_dt, legs, alt_ft, vs_fpm, descent_start_choice
|
| 1395 |
)
|
| 1396 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
| 1397 |
return html_map, stats_md, bfo_fig, bfo_md, alt_fig, headers, rows
|
| 1398 |
|
| 1399 |
|
| 1400 |
def run_sim_manual(start_lat, start_lon, legs_df, method, start_time_str, alt_ft, vs_fpm, descent_start_choice):
|
| 1401 |
legs = parse_legs_df(legs_df)
|
| 1402 |
start_dt = _parse_start_time(start_time_str)
|
| 1403 |
+
|
| 1404 |
html_map, stats_md, bfo_fig, bfo_md, alt_fig, headers, rows = _run_common(
|
| 1405 |
start_lat, start_lon, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice
|
| 1406 |
)
|
|
|
|
| 1441 |
], value=[])
|
| 1442 |
return gr.update(), f"### ❌ {report_md}", "", "", fig, "", alt_fig, empty_table
|
| 1443 |
|
| 1444 |
+
# Guided always starts at radar fix
|
| 1445 |
html_map, stats_md, bfo_fig, bfo_md, alt_fig, headers, rows = _run_common(
|
| 1446 |
RADAR_FIX_LAT, RADAR_FIX_LON, legs, method, start_dt, alt_ft, vs_fpm, descent_start_choice
|
| 1447 |
)
|
|
|
|
| 1472 |
f"- **Last radar spot (fixed):** 18:22:00Z @ **{RADAR_FIX_LAT:.3f}N, {RADAR_FIX_LON:.3f}E**\n"
|
| 1473 |
f"- **Turn completion (fixed):** 18:40:00Z\n"
|
| 1474 |
f"- **Initial track (fixed):** {RADAR_INITIAL_TRACK_DEG:.1f}°\n"
|
| 1475 |
+
f"- **Final resting point:** **00:20:00Z** (trajectory clipped exactly at this time)"
|
| 1476 |
)
|
| 1477 |
|
| 1478 |
with gr.Row():
|