diff --git a/retina_simulation/generator.py b/retina_simulation/generator.py index ce3bd31..3a54d4c 100644 --- a/retina_simulation/generator.py +++ b/retina_simulation/generator.py @@ -20,6 +20,26 @@ # Each tower: (lat, lon, alt_ft, freq_hz, callsign) # Modeled after real VHF/UHF broadcast transmitters suitable for passive radar +# Effective radiated power, dBm, keyed by callsign. Kept beside the tower +# tuples rather than inside them so the 5-element unpacking used throughout +# this module stays valid. Only populated for sites taken from real FCC +# records; illuminator selection falls back to _DEFAULT_EIRP_DBM otherwise. +# +# The spread here is 65 dB — Caesars Head at 92.2 dBm against Spartanburg at +# 27.0 — so a fleet that ignores EIRP treats a 0.5 W transmitter as the equal +# of a megawatt one. +_TOWER_EIRP_DBM = { + "WYFF": 92.2, # Caesars Head + "WMYA-TV": 90.9, # Fountain Inn + "WNTV": 84.7, # Paris Mountain — Tower Finder's top pick for the metro + "WLOS": 83.7, # Mt Pisgah + "WSPA-TV": 77.4, # Hogback Mtn + "BLP00776": 57.0, # near the core, low power + "W07DT-D": 50.0, # Tryon NC + "BLP01065": 27.0, # Spartanburg — geometrically valuable, radiologically weak +} +_DEFAULT_EIRP_DBM = 80.0 + _TOWERS_US = [ # East Coast (33.75667, -84.33184, 1600, 195_000_000, "WSB-TV"), # Atlanta @@ -34,6 +54,22 @@ (30.33270, -81.65560, 1200, 575_000_000, "WJXT"), # Jacksonville (36.85260, -75.97820, 1300, 539_000_000, "WAVY"), # Norfolk (35.78700, -78.78170, 1500, 563_000_000, "WRAL"), # Raleigh + # ── Greenville SC ──────────────────────────────────────────────────────── + # Real FCC facilities from the Tower Finder illuminator search, one entry + # per *distinct site*. Twenty stations serve this market but they share + # only eight masts — Paris Mountain alone carries WNTV, WRET-TV, WGGS-TV, + # W10AJ-D and five LPTVs. Co-sited transmitters are worthless as a + # bistatic pair (identical geometry, identical ellipse), so the table lists + # sites and the strongest station at each. + # alt is the radiating centre AMSL (ground + antenna height), in feet. + (34.941222, -82.410278, 3315, 183_000_000, "WNTV"), # Paris Mountain + (34.647500, -82.270000, 1847, 599_000_000, "WMYA-TV"), # Fountain Inn — south + (35.170194, -82.290500, 5437, 201_000_000, "WSPA-TV"), # Hogback Mtn + (35.111944, -82.606389, 5058, 569_000_000, "WYFF"), # Caesars Head + (35.222222, -82.549444, 5220, 213_000_000, "WLOS"), # Mt Pisgah + (34.970111, -81.948391, 794, 195_000_000, "BLP01065"), # Spartanburg — east + (35.266278, -82.244111, 3186, 177_000_000, "W07DT-D"), # Tryon NC + (34.875111, -82.338211, 984, 183_000_000, "BLP00776"), # near the core # Midwest (41.87150, -87.62440, 1650, 191_000_000, "WBBM-TV"), # Chicago (42.33140, -83.04580, 1200, 551_000_000, "WXYZ-TV"), # Detroit @@ -154,9 +190,34 @@ (33.74900, -84.38800, 1050, 199_000_000, "WSB-RING", 33.6407, -84.4277), # ATL (39.73920, -104.99030, 5300, 201_000_000, "KCNC-RING", 39.8561, -104.6737), # DEN (39.09970, -94.57860, 900, 203_000_000, "KMBC-RING", 39.2976, -94.7139), # MCI Kansas City + # WSPA-TV, RF ch 11 (201 MHz) on Hogback Mtn, 31 km NNW of GSP. It shares + # 201 MHz with the DEN ring above, which is harmless: rings never overlap + # geographically. (An earlier version of this note said the market had no + # other VHF station — it has several; see _TOWERS_US. The ring uses this + # one because it is the strongest VHF site with a clear line to the core.) + (35.170194, -82.290500, 5437, 201_000_000, "WSPA-RING", 34.8957, -82.2189), # GSP Greenville SC ] +# ── Metro areas ─────────────────────────────────────────────────────────────── +# Shared by the generator's --metro scoping and the orchestrator's --metros +# post-generation filter, so the two can never disagree about where a metro is. +_KNOWN_METROS = { + "atl": {"name": "Atlanta", "lat": 33.749, "lon": -84.388, "radius_nm": 80}, + "gvl": {"name": "Greenville", "lat": 34.852, "lon": -82.394, "radius_nm": 60}, + "clt": {"name": "Charlotte", "lat": 35.227, "lon": -80.843, "radius_nm": 70}, + "nyc": {"name": "New York", "lat": 40.748, "lon": -73.986, "radius_nm": 80}, + "dca": {"name": "Washington DC", "lat": 38.935, "lon": -77.079, "radius_nm": 70}, + "chi": {"name": "Chicago", "lat": 41.872, "lon": -87.624, "radius_nm": 80}, + "den": {"name": "Denver", "lat": 39.739, "lon": -104.990, "radius_nm": 80}, + "lax": {"name": "Los Angeles", "lat": 34.052, "lon": -118.244, "radius_nm": 80}, + "dfw": {"name": "Dallas-Fort Worth", "lat": 32.897, "lon": -97.038, "radius_nm": 80}, + "kc": {"name": "Kansas City", "lat": 39.298, "lon": -94.714, "radius_nm": 70}, +} + +_NM_TO_KM = 1.852 + + def _haversine_km(lat1: float, lon1: float, lat2: float, lon2: float) -> float: """Return great-circle distance in km between two lat/lon points.""" R = 6371.0 @@ -252,23 +313,31 @@ class GeneratedNodeConfig: tx_alt_ft: float fc_hz: float fs_hz: float = 2_000_000.0 - beam_width_deg: float = 40.0 + beam_width_deg: float = 42.0 max_range_km: float = 50.0 region: str = "us" tx_callsign: str = "" beam_azimuth_deg: float | None = None # explicit Yagi aim; None → broadside + # Bistatic range limit: (RX→target) + (target→TX) − baseline. That sum is + # what the delay measures and what sets received power, so it is the + # physical detection limit; max_range_km is a monostatic approximation + # retained for hardware nodes. None → omitted from the wire format. + max_bistatic_range_km: float | None = None def _node_dict(node: GeneratedNodeConfig) -> dict: - """Serialize a node, omitting an unset beam_azimuth_deg. + """Serialize a node, omitting unset optional geometry keys. The backend solver does `float(node_cfg["beam_azimuth_deg"])` whenever the key is present, so a null on the wire would crash it. Dropping the key for broadside nodes makes the backend fall back to its own broadside auto-aim. + max_bistatic_range_km is dropped for the same reason and so that its + absence unambiguously means "use the monostatic rule". """ d = asdict(node) - if d.get("beam_azimuth_deg") is None: - d.pop("beam_azimuth_deg", None) + for _optional in ("beam_azimuth_deg", "max_bistatic_range_km"): + if d.get(_optional) is None: + d.pop(_optional, None) return d @@ -623,6 +692,266 @@ def _place_rx_on_land( return (round(tx_lat, 6), round(tx_lon, 6)) +def _subtended_deg(from_lat, from_lon, a, b) -> float: + """Angle between two towers as seen from a point, in degrees. + + This — not the towers' bearing separation from the metro centre — is what + conditions a two-illuminator fix. The bistatic range gradient is + b = u_tx + u_rx, so two measurements are independent to the extent their + transmitters lie in different directions *from the target*. Two towers on + a similar bearing from the metro core but far apart in range still subtend + a usable angle across most of the coverage area. + """ + ba = _bearing_between(from_lat, from_lon, a[0], a[1]) + bb = _bearing_between(from_lat, from_lon, b[0], b[1]) + return abs((ba - bb + 180.0) % 360.0 - 180.0) + + +def _beam_footprint(rx_lat, rx_lon, beam_azimuth_deg, beam_width_deg, max_range_km): + """Sample points across a receiver's beam, for evaluating pair geometry. + + Deliberately samples the *edges* as well as the centre: a pair can condition + well at boresight and collapse at the beam edge, and selecting on a single + representative point would bake that blind spot in. + """ + R = 6371.0 + pts = [] + half = beam_width_deg / 2.0 + for frac in (0.35, 0.7, 1.0): + for off in (-half, -half / 2, 0.0, half / 2, half): + br = math.radians((beam_azimuth_deg + off) % 360.0) + d = max_range_km * frac + pts.append( + ( + rx_lat + math.degrees((d * math.cos(br)) / R), + rx_lon + math.degrees((d * math.sin(br)) / (R * math.cos(math.radians(rx_lat)))), + ) + ) + return pts + + +def _pick_illuminator_pair(rx_lat, rx_lon, beam_azimuth_deg, beam_width_deg, max_range_km, towers, min_eirp_dbm): + """Choose the two towers giving the best-conditioned pair for this receiver. + + Scored on the *worst* subtended angle across the beam footprint rather than + the mean, because a pair that degenerates anywhere in the footprint is + unreliable there. Towers are already one-per-site in _TOWERS_US; co-sited + transmitters would share an ellipse and be worthless as a pair. + + Returns (tower_a, tower_b) or None when nothing clears min_eirp_dbm. + """ + usable = [t for t in towers if _TOWER_EIRP_DBM.get(t[4], _DEFAULT_EIRP_DBM) >= min_eirp_dbm] + if len(usable) < 2: + return None + probes = _beam_footprint(rx_lat, rx_lon, beam_azimuth_deg, beam_width_deg, max_range_km) + best, best_score = None, -1.0 + for i in range(len(usable)): + for j in range(i + 1, len(usable)): + a, b = usable[i], usable[j] + if _haversine_km(a[0], a[1], b[0], b[1]) < 1.0: + continue # same mast + worst = min(_subtended_deg(p[0], p[1], a, b) for p in probes) + if worst > best_score: + best, best_score = (a, b), worst + return best + + +def _generate_dual_sites( + n_sites: int, + core_lat: float, + core_lon: float, + towers: list[tuple], + metro_radius_km: float, + prefix: str = "synth-GVL", + beam_width_deg: float = 42.0, + max_bistatic_range_km: float = 60.0, + min_eirp_dbm: float = 40.0, + aim: str = "core", + aim_jitter_deg: float = 30.0, +) -> list[dict]: + """Generate n_sites receivers, each running two nodes on two illuminators. + + This is how a real passive-radar site is built: one antenna, one RX + position, several receiver chains tuned to different transmitters. The two + nodes therefore share rx position, altitude, beam azimuth, beam width and + range — they differ only in which tower they listen to. + + The geometric payoff is that the two bistatic ellipses share a focus (the + common RX), so they intersect in at most two points and the beam almost + always excludes one. A single site localises on its own, with the residual + ambiguity bounded by the antenna pattern rather than by a second receiver + tens of km away. + + Position is well determined this way; velocity is not. One pair gives two + Doppler projections for three velocity components, so it stays + under-determined unless the level-flight assumption is applied (which the + association stage does). Two pairs — four nodes — give eight residuals + against five unknowns, and only then do the solver's residual gates regain + the discriminating power they lack at n=2. + + Sites are placed at random across the metro rather than ringed around a + core, so their beams overlap each other far less than a ring's do. + """ + if n_sites <= 0 or len(towers) < 2: + return [] + + nodes = [] + for i in range(n_sites): + node_id = f"{prefix}-{i + 1:04d}" + rx_lat, rx_lon = _place_rx_on_land( + core_lat, + core_lon, + dist_min_km=5.0, + dist_max_km=max(10.0, metro_radius_km * 0.85), + display_node_id=node_id, + ) + # Aim. "random" spreads sectors and minimises inter-site overlap, but + # a beam pointed away from the traffic sees nothing: with 85% of + # aircraft routed through the metro core, random aiming left most sites + # idle and collapsed the solve rate by an order of magnitude. + # + # "core" aims at the core with jitter, which is also what a real + # operator would do — receivers are sited to cover the airspace of + # interest. Inter-site overlap is not the enemy here the way it is for + # a ring: each dual site already self-solves from its own two + # illuminators, so a second site overlapping it upgrades the fix to + # four nodes rather than manufacturing a two-node ambiguity. + if aim == "random": + beam_azimuth = random.uniform(0.0, 360.0) + else: + beam_azimuth = ( + _bearing_between(rx_lat, rx_lon, core_lat, core_lon) + random.uniform(-aim_jitter_deg, aim_jitter_deg) + ) % 360.0 + pair = _pick_illuminator_pair( + rx_lat, + rx_lon, + beam_azimuth, + beam_width_deg, + max_bistatic_range_km, + towers, + min_eirp_dbm, + ) + if pair is None: + continue + rx_alt_ft = round(random.uniform(100, 1500), 1) + for suffix, tower in zip("ab", pair): + tx_lat, tx_lon, tx_alt_ft, fc_hz, callsign = tower + node = GeneratedNodeConfig( + node_id=f"{node_id}{suffix}", + rx_lat=round(rx_lat, 6), + rx_lon=round(rx_lon, 6), + rx_alt_ft=rx_alt_ft, + tx_lat=tx_lat, + tx_lon=tx_lon, + tx_alt_ft=tx_alt_ft, + fc_hz=fc_hz, + fs_hz=2_000_000, + beam_width_deg=round(beam_width_deg, 1), + max_range_km=round(max_bistatic_range_km, 1), + region="us", + tx_callsign=callsign, + beam_azimuth_deg=round(beam_azimuth, 2), + max_bistatic_range_km=round(max_bistatic_range_km, 1), + ) + nodes.append(_node_dict(node)) + return nodes + + +def _generate_metro_solo( + n: int, + core_lat: float, + core_lon: float, + towers: list[tuple], + metro_radius_km: float, + prefix: str = "synth-SOLO", + beam_width_deg: float = 42.0, + max_bistatic_range_km: float = 60.0, + ring_radius_km: float = 18.0, + start_bearing_deg: float = 30.0, +) -> list[dict]: + """Generate n isolated receivers on the metro rim, each aimed outward. + + These exist to keep the single-node ellipse-arc path exercised. Every + receiver in the coverage ring aims *inward* at the core, so their beams all + intersect and essentially every detection associates into a multinode + solve — leaving the single-node arc code with almost no live coverage. + + Isolation here is by beam geometry, not distance: the metro is far too + small for the nationwide pool's 400 km separation. Each solo RX sits near + the rim and points away from the core, so its sector cannot intersect the + inward-aimed ring beams no matter how the range circles overlap. The + association overlap zone is computed from beam sectors + (compute_overlap_zone), so this is the property that actually decides + whether detections stay single-node. + + Each also takes its own illuminator rather than the shared ring TX, which + puts its Doppler outside the association gate — a second, independent + reason not to pair. + """ + if n <= 0 or not towers: + return [] + + R = 6371.0 + # Sit between the ring envelope and the metro edge. Far enough out that + # the outward beam looks away from ring airspace; inside the metro radius + # so the node stays on the map with the rest of the fleet. + rim_km = max(ring_radius_km * 2.0, metro_radius_km * 0.80) + + nodes = [] + for i in range(n): + # Offset from the ring's own start bearing so a solo node never lands + # on top of a ring receiver. + bearing_deg = (start_bearing_deg + 360.0 * i / max(n, 1)) % 360.0 + bearing_rad = math.radians(bearing_deg) + dlat = (rim_km * math.cos(bearing_rad)) / R + dlon = (rim_km * math.sin(bearing_rad)) / (R * math.cos(math.radians(core_lat))) + rx_lat = core_lat + math.degrees(dlat) + rx_lon = core_lon + math.degrees(dlon) + + node_id = f"{prefix}-{i + 1:04d}" + if not _candidate_is_safe(rx_lat, rx_lon, node_id): + rx_lat, rx_lon = _place_rx_on_land( + core_lat, + core_lon, + dist_min_km=rim_km - 10, + dist_max_km=rim_km + 10, + display_node_id=node_id, + ) + + # Pick the illuminator *furthest* from the core among the metro towers: + # its baseline points away from ring airspace, so the bistatic ellipse + # opens outward too rather than folding back over the mesh. + tower = max( + towers, + key=lambda t: _haversine_km(t[0], t[1], rx_lat, rx_lon), + ) + tx_lat, tx_lon, tx_alt_ft, fc_hz, callsign = tower + + # Aim directly away from the core — the isolation property. + beam_azimuth = (_bearing_between(rx_lat, rx_lon, core_lat, core_lon) + 180.0) % 360.0 + + node = GeneratedNodeConfig( + node_id=node_id, + rx_lat=round(rx_lat, 6), + rx_lon=round(rx_lon, 6), + rx_alt_ft=round(random.uniform(100, 1500), 1), + tx_lat=tx_lat, + tx_lon=tx_lon, + tx_alt_ft=tx_alt_ft, + fc_hz=fc_hz, + fs_hz=2_000_000, + beam_width_deg=round(beam_width_deg, 1), + max_range_km=round(max_bistatic_range_km, 1), + region="us", + tx_callsign=callsign, + beam_azimuth_deg=round(beam_azimuth, 2), + max_bistatic_range_km=round(max_bistatic_range_km, 1), + ) + nodes.append(_node_dict(node)) + + return nodes + + def _generate_coverage_ring( n: int, core_lat: float, @@ -630,7 +959,7 @@ def _generate_coverage_ring( tx_tower: tuple, prefix: str = "synth-RING", radius_km: float = 18.0, - beam_width_deg: float = 50.0, + beam_width_deg: float = 42.0, max_range_km: float = 60.0, aim: str = "core", start_bearing_deg: float = 0.0, @@ -685,12 +1014,183 @@ def _generate_coverage_ring( region="us", tx_callsign=callsign, beam_azimuth_deg=round(beam_azimuth, 2), + # Ring receivers model real passive-radar reach, so they are limited + # on bistatic range. max_range_km stays for consumers that have not + # been taught the bistatic rule. + max_bistatic_range_km=round(max_range_km, 1), + ) + nodes.append(_node_dict(node)) + + return nodes + + +def _generate_scatter_sites( + n: int, + core_lat: float, + core_lon: float, + towers: list[tuple], + metro_radius_km: float, + prefix: str = "synth-SCAT", + beam_width_deg: float = 42.0, + max_bistatic_range_km: float = 60.0, + aim_sigma_deg: float = 25.0, + frac_off_core: float = 0.25, + clump_sigma_km: float = 7.0, +) -> list[dict]: + """Generate n receivers scattered the way a real deployment actually lands. + + The ring and dual layouts are *designs*: someone chose where every receiver + goes and what it points at, to buy geometry. A community fleet is not + designed. Receivers appear where operators live and point where operators + care, and the resulting geometry is whatever falls out. This layout models + that, so the solver is measured against the fleet it will really get rather + than the one we would have built. + + Four properties, each with a reason: + + - **Clumped, not uniform.** Sites concentrate around the metro core and + around the towns the metro's real broadcast towers serve (an FM/TV tower + is sited for population, so the tower list doubles as a population proxy). + Clumping matters because co-located receivers see near-parallel bistatic + range gradients — the high-GDOP case a ring is specifically built to avoid. + - **Own illuminator per site.** Each picks a nearby strong tower rather than + sharing one, weighted toward short baselines the way an operator picks the + station that comes in best. Diverse illuminators mean per-node Doppler no + longer sits inside one shared association gate. + - **Aimed by hand.** Most point roughly at the core airspace with real + pointing error; a minority point somewhere else entirely. + - **Heterogeneous hardware.** Beamwidth and bistatic reach vary per site. + 60 km is what a good setup achieves, not what an average one does, so the + reach distribution is skewed below it with a tail reaching it. + """ + if n <= 0 or not towers: + return [] + + R = 6371.0 + + # Population anchors: the core, weighted heavily, plus each real tower. + anchors = [(core_lat, core_lon)] * max(2, len(towers) // 2) + anchors += [(t[0], t[1]) for t in towers] + + def _scatter_around(anchor_lat, anchor_lon): + d_north = random.gauss(0, clump_sigma_km) + d_east = random.gauss(0, clump_sigma_km) + lat = anchor_lat + math.degrees(d_north / R) + lon = anchor_lon + math.degrees(d_east / (R * math.cos(math.radians(anchor_lat)))) + return lat, lon + + nodes = [] + for i in range(n): + node_id = f"{prefix}-{i + 1:04d}" + + rx_lat = rx_lon = None + for _ in range(40): + a_lat, a_lon = random.choice(anchors) + cand_lat, cand_lon = _scatter_around(a_lat, a_lon) + # Stay inside the metro the rest of the fleet and the traffic model + # live in; a site outside it just never sees an aircraft. + if _haversine_km(cand_lat, cand_lon, core_lat, core_lon) > metro_radius_km: + continue + if _candidate_is_safe(cand_lat, cand_lon, node_id): + rx_lat, rx_lon = cand_lat, cand_lon + break + if rx_lat is None: + rx_lat, rx_lon = _place_rx_on_land( + core_lat, + core_lon, + dist_min_km=2.0, + dist_max_km=max(5.0, metro_radius_km * 0.8), + display_node_id=node_id, + ) + + # Illuminator: prefer a short baseline, but not deterministically — + # 1/d² weighting reproduces "whichever strong station comes in best" + # without every site in a clump converging on the same tower. + weighted = [] + for t in towers: + d = _haversine_km(rx_lat, rx_lon, t[0], t[1]) + if d < 4.0 or d > 75.0: + continue + weighted.append((t, 1.0 / (d * d))) + if weighted: + total = sum(w for _t, w in weighted) + pick = random.uniform(0, total) + acc = 0.0 + tower = weighted[-1][0] + for t, w in weighted: + acc += w + if acc >= pick: + tower = t + break + else: + tower = min(towers, key=lambda t: _haversine_km(rx_lat, rx_lon, t[0], t[1])) + tx_lat, tx_lon, tx_alt_ft, fc_hz, callsign = tower + + if random.random() < frac_off_core: + beam_azimuth = random.uniform(0.0, 360.0) + else: + beam_azimuth = ( + _bearing_between(rx_lat, rx_lon, core_lat, core_lon) + random.gauss(0, aim_sigma_deg) + ) % 360.0 + + # All fleet antennas are identical 42-degree Yagis — no width jitter. + width = beam_width_deg + # Mode at 0.7 of the ceiling: most setups fall short of the best case. + # NOTE: `reach` is a DIFFERENTIAL range (Δ = R_tx + R_rx − L), and it + # is assigned to both max_bistatic_range_km (its true meaning) and + # max_range_km (the RX-circle approximation of it). Consumers that + # know the bistatic rule read the first; the second is only a + # fallback and reads up to 2x too large away from the transmitter — + # same convention as the ring path. + reach = max_bistatic_range_km * random.triangular(0.40, 1.0, 0.70) + + node = GeneratedNodeConfig( + node_id=node_id, + rx_lat=round(rx_lat, 6), + rx_lon=round(rx_lon, 6), + rx_alt_ft=round(random.uniform(100, 1500), 1), + tx_lat=tx_lat, + tx_lon=tx_lon, + tx_alt_ft=tx_alt_ft, + fc_hz=fc_hz, + fs_hz=2_000_000, + beam_width_deg=round(width, 1), + max_range_km=round(reach, 1), + region="us", + tx_callsign=callsign, + beam_azimuth_deg=round(beam_azimuth, 2), + max_bistatic_range_km=round(reach, 1), ) nodes.append(_node_dict(node)) return nodes +def _resolve_metro(metro: str) -> dict: + """Resolve a metro code (e.g. "gvl") to its descriptor in _KNOWN_METROS.""" + key = metro.strip().lower() + if key not in _KNOWN_METROS: + raise ValueError(f"Unknown metro {metro!r} (available: {', '.join(sorted(_KNOWN_METROS))})") + return _KNOWN_METROS[key] + + +def _towers_in_metro(towers: list, metro: dict) -> list: + """Towers within the metro's own radius of its centre.""" + radius_km = metro["radius_nm"] * _NM_TO_KM + return [t for t in towers if _haversine_km(t[0], t[1], metro["lat"], metro["lon"]) <= radius_km] + + +def _rings_in_metro(ring_spec: list, metro: dict) -> list: + """Ring specs whose airspace core lies inside the metro. + + Cores are matched (not the illuminators) because the core is what the ring + and its coverage cell are built around — a TX can legitimately sit outside + the metro radius while lighting an airspace inside it. + """ + radius_km = metro["radius_nm"] * _NM_TO_KM + return [s for s in ring_spec if _haversine_km(s[5], s[6], metro["lat"], metro["lon"]) <= radius_km] + + def _active_rings(n_cluster: int, n_clusters: int, ring_spec: list = _RING_TXS): """Yield (ring_id, spec, size) for each ring the budget actually produces. @@ -715,6 +1215,8 @@ def coverage_cells( n_clusters: int = 1, ring_spec: list = _RING_TXS, traffic_radius_km: float = 70.0, + metro: str | None = None, + layout: str = "ring", ) -> list[dict]: """First-class metro-cell descriptors for the active rings. @@ -722,7 +1224,33 @@ def coverage_cells( never reconstructed from receiver positions — so water-displaced receivers cannot drift the hub-radial aim point. ops_weight defaults to ring size; a caller with real traffic figures can override the spec to inject ops/yr. + + When *metro* is set the spec is narrowed to that metro's rings with the same + filter generate_fleet uses, so the cells always describe the nodes that were + actually generated. """ + if layout in ("dual", "scatter"): + # These layouts replace the ring; emitting ring cells for them wrote a + # fleet_config.json describing an airspace no generated node was + # placed around — exactly the disagreement the docstring above rules + # out. Both layouts orbit the metro core, so one core cell carries + # the traffic weighting. + if not metro: + return [] + m = _resolve_metro(metro) + return [ + { + "ring_id": f"synth-{layout.upper()}", + "core_lat": m["lat"], + "core_lon": m["lon"], + "radius_km": traffic_radius_km, + "ops_weight": float(max(1, n_cluster)), + "illuminator": "", + } + ] + if metro: + ring_spec = _rings_in_metro(ring_spec, _resolve_metro(metro)) + n_clusters = min(n_clusters, len(ring_spec)) cells = [] for ring_id, spec, size in _active_rings(n_cluster, n_clusters, ring_spec): tx_lat, tx_lon, tx_alt_ft, fc_hz, callsign, core_lat, core_lon = spec @@ -748,10 +1276,16 @@ def generate_fleet( n_cluster: int = 8, n_clusters: int = 1, ring_radius_km: float = 18.0, - ring_beam_width_deg: float = 50.0, + ring_beam_width_deg: float = 42.0, ring_max_range_km: float = 60.0, ring_aim: str = "core", ring_spec: list = _RING_TXS, + metro: str | None = None, + layout: str = "ring", + illuminator_band: str = "any", + dual_min_eirp_dbm: float = 40.0, + dual_aim: str = "core", + dual_fraction: float = 0.0, ) -> list[dict]: """Generate a fleet of synthetic node configurations. @@ -773,6 +1307,10 @@ def generate_fleet( illuminator. Diverse look angles give low-GDOP, velocity-observable multinode fixes. Ring slots are carved out of the metro allocation so total stays n_nodes. + layout="dual" and layout="scatter" each spend that same n_cluster budget on a + different arrangement instead of the ring — see _generate_dual_sites and + _generate_scatter_sites. + Args: n_nodes: Total nodes to generate (100-1000). regions: List of regions to distribute across ["us", "eu", "au"]. @@ -787,6 +1325,16 @@ def generate_fleet( ring_aim: "core" (aim at metro core) or "broadside" (perp to TX). ring_spec: Metro ring table (defaults to _RING_TXS); inject to add metros or change illuminators without editing library source. + metro: Restrict the whole fleet to one metro area (a _KNOWN_METROS code + such as "gvl"). Towers and rings outside that metro's radius are + dropped and solo/rural placement is disabled, so every node lands in + the one metro. None (default) keeps the continent-wide behaviour. + dual_fraction: Fraction (0.0-1.0) of n_nodes to additionally run as + dual-illuminator sites (see _generate_dual_sites), carved out of + the ring/metro budget and appended after it. Requires metro + (ValueError otherwise, same as layout="dual"/"scatter"). Ignored + when layout == "dual" — the whole cluster budget is already dual + sites there. Returns: List of node config dicts ready for fleet_config.json. @@ -794,27 +1342,53 @@ def generate_fleet( if regions is None: regions = ["us"] + # Reproducibility caveat: --seed fixes the RNG stream (it is re-seeded + # again after the network tower lookup below, so lookup retries cannot + # shift it), but the *content* of the Tower API response and the + # availability of shapely for the land check both feed placement + # decisions. Same seed + same tower cache + same optional deps ⇒ same + # fleet; a flaky API or a machine without shapely will differ. random.seed(seed) + metro_area = _resolve_metro(metro) if metro else None + tower_db = { "us": _TOWERS_US, "eu": _TOWERS_EU, "au": _TOWERS_AU, } - # Solo towers — only available for US region (where rural towers are defined) - solo_towers = _TOWERS_SOLO_US if "us" in regions else [] + # Solo towers — only available for US region (where rural towers are defined). + # + # The nationwide pool separates receivers by 400 km, which cannot apply + # inside a metro, so --metro uses metro-scoped solo placement instead (see + # _metro_solo_nodes below). Both exist for the same reason: solo receivers + # are the only way to exercise the single-node ellipse-arc path. Without + # them every detection lands in an overlap zone and associates into a + # multinode solve — measured on the Greenville fleet as 15 of 16 nodes + # overlapping 1-11 neighbours, and single-node arcs nearly absent. + solo_towers = _TOWERS_SOLO_US if ("us" in regions and not metro_area) else [] # Distribute nodes across regions proportionally to tower count available_towers = [] for region in regions: towers = tower_db.get(region, []) + if metro_area: + towers = _towers_in_metro(towers, metro_area) for t in towers: available_towers.append((region, t)) if not available_towers: + if metro_area: + raise ValueError( + f"No towers within {metro_area['radius_nm']} nm of {metro_area['name']} for regions: {regions}" + ) raise ValueError(f"No towers available for regions: {regions}") + if metro_area: + ring_spec = _rings_in_metro(ring_spec, metro_area) + n_clusters = min(n_clusters, len(ring_spec)) + # ── Pre-fetch real towers from Tower API for metro areas ────────────────── # Each metro area gets multiple real towers so nodes in the same city # use DIFFERENT transmitters instead of all sharing the same one. @@ -829,10 +1403,9 @@ def _cache_key(lat, lon): from retina_simulation.tower_resolver import lookup_metro_towers except ImportError: from tower_resolver import lookup_metro_towers - all_metro_centers = [] - for region in regions: - for t in tower_db.get(region, []): - all_metro_centers.append(t) + # Only the towers actually in play — under --metro this is a handful + # of centres instead of every metro on the continent. + all_metro_centers = [t for _region, t in available_towers] metro_api_towers_raw = lookup_metro_towers(all_metro_centers, radius_km=80, limit=50) # Map back to (lat, lon) → tower list for t in all_metro_centers: @@ -843,12 +1416,76 @@ def _cache_key(lat, lon): import logging logging.warning("Tower API lookup failed, using hardcoded towers: %s", exc) + # Re-seed after the lookup: any RNG the HTTP/cache path consumed (or will + # consume differently on retry) must not shift the placement stream. + random.seed(seed) # Allocate solo and cluster node counts, carving both from metro allocation - n_solo = max(1, round(n_nodes * solo_fraction)) if solo_towers else 0 - n_cluster = max(0, n_cluster) + if solo_towers: + n_solo = max(1, round(n_nodes * solo_fraction)) + elif metro_area: + # Metro-scoped solo receivers, placed on the rim and aimed outward. + n_solo = max(1, round(n_nodes * solo_fraction)) + else: + n_solo = 0 + # No rings survived the metro filter → give their budget back to metro nodes + # instead of silently generating fewer nodes than asked for. + if layout in ("dual", "scatter"): + # These layouts replace the ring and take their budget straight from + # n_cluster. Gating it on the ring table surviving the metro filter + # (n_clusters > 0) zeroed the whole layout for any metro without a + # _RING_TXS entry — the ring table is irrelevant to them. + n_cluster = max(0, n_cluster) + else: + n_cluster = max(0, n_cluster) if n_clusters > 0 else 0 n_metro = max(0, n_nodes - n_solo - n_cluster) + # dual_fraction carve: a slice of the SAME layout's ring/metro budget + # additionally runs as dual-illuminator sites (see _generate_dual_sites), + # independent of layout="dual" (already all-dual there — the whole + # cluster budget went to dual sites, so dual_fraction is a no-op). + # Solo carve is untouched: solo answers a different question + # (single-node ellipse-arc coverage) than dual does. + n_dual_sites = 0 + if dual_fraction > 0 and layout != "dual": + if not metro_area: + raise ValueError("dual_fraction > 0 requires --metro: dual sites are placed around a metro core") + n_dual_nodes = min(round(n_nodes * dual_fraction / 2) * 2, n_cluster + n_metro) + _from_cluster = min(n_dual_nodes, n_cluster) + n_cluster -= _from_cluster + n_metro -= n_dual_nodes - _from_cluster + n_dual_sites = n_dual_nodes // 2 + + def _dual_fraction_sites(n_sites: int) -> list[dict]: + """Extra dual-illuminator sites for the dual_fraction carve above. + + Mirrors the layout="dual" branch's own _generate_dual_sites call + (same prefix/towers/aim construction) so a fraction-carved site is + indistinguishable from a full dual-layout one. Called last — after + every other node in this fleet has consumed its RNG draws — so + dual_fraction=0.0 (n_sites=0, no call) reproduces today's scene + byte-for-byte at the same seed. + """ + if n_sites <= 0: + return [] + _dual_towers = [t for _r, t in available_towers] or _TOWERS_US + if illuminator_band == "vhf": + _vhf = [t for t in _dual_towers if t[3] < 300e6] + if len(_vhf) >= 2: + _dual_towers = _vhf + return _generate_dual_sites( + n_sites=n_sites, + core_lat=metro_area["lat"], + core_lon=metro_area["lon"], + towers=_dual_towers, + metro_radius_km=metro_area["radius_nm"] * _NM_TO_KM, + prefix=f"synth-{metro.strip().upper()}-DUAL" if metro else "synth-DUAL", + beam_width_deg=ring_beam_width_deg, + max_bistatic_range_km=ring_max_range_km, + min_eirp_dbm=dual_min_eirp_dbm, + aim=dual_aim, + ) + # Track how many times each API tower has been used (per metro) for # round-robin distribution — avoids all nodes sharing one tower. _metro_tower_idx: dict[str, int] = {} @@ -872,7 +1509,9 @@ def _cache_key(lat, lon): fc_hz = t["fc_hz"] callsign = t["tx_callsign"] - region_prefix = region.upper() + # Metro-scoped fleets are named for the metro (synth-GVL-0001) rather + # than the continent, so node IDs say where they actually are. + region_prefix = metro.strip().upper() if metro else region.upper() node_id = f"synth-{region_prefix}-{i + 1:04d}" rx_lat, rx_lon = _place_rx_on_land( tx_lat, @@ -884,7 +1523,7 @@ def _cache_key(lat, lon): rx_alt_ft = random.uniform(100, 2000) node_fc = fc_hz + random.choice([-500000, 0, 0, 0, 500000]) - beam_width = random.uniform(35, 45) + beam_width = 42.0 # identical 42-degree Yagis fleet-wide max_range = random.uniform(35, 55) node = GeneratedNodeConfig( @@ -901,11 +1540,32 @@ def _cache_key(lat, lon): max_range_km=round(max_range, 1), region=region, tx_callsign=callsign, + # Every bistatic receiver is bounded by differential range; a circle + # on the RX is never the true footprint. These base nodes were the + # last path still declaring only a monostatic limit, which left them + # gating and rendering as circles while the ring, solo and dual + # paths all used the ellipse. The randomised value carries over + # unchanged — it is the same number, read correctly. + max_bistatic_range_km=round(max_range, 1), ) nodes.append(_node_dict(node)) # --- Solo nodes (isolated — strictly one node per unique tower position) --- - if n_solo > 0: + if n_solo > 0 and metro_area: + # Metro-scoped: isolation comes from aiming away from the core, not + # from the nationwide pool's 400 km separation (impossible in a metro). + nodes.extend( + _generate_metro_solo( + n=n_solo, + core_lat=metro_area["lat"], + core_lon=metro_area["lon"], + towers=[t for _region, t in available_towers] or _TOWERS_US, + metro_radius_km=metro_area["radius_nm"] * _NM_TO_KM, + ring_radius_km=ring_radius_km, + max_bistatic_range_km=ring_max_range_km, + ) + ) + elif n_solo > 0: # All US positions that must be avoided when extending the pool # Avoid positions: metro towers only. Named solo towers are gated # inside _extend_solo_pool with the same min_sep check so they are @@ -935,7 +1595,7 @@ def _cache_key(lat, lon): ) rx_alt_ft = random.uniform(100, 1500) - beam_width = random.uniform(35, 45) + beam_width = 42.0 # identical 42-degree Yagis fleet-wide max_range = random.uniform(35, 55) node = GeneratedNodeConfig( @@ -952,6 +1612,11 @@ def _cache_key(lat, lon): max_range_km=round(max_range, 1), region="us", tx_callsign=callsign, + # Same bistatic bound the base/ring/dual paths declare — these + # were the last nodes gating as monostatic circles, and the + # one path specifically meant to exercise single-node ellipse + # arcs. Same number, read correctly (see the base-node note). + max_bistatic_range_km=round(max_range, 1), ) nodes.append(_node_dict(node)) @@ -962,6 +1627,59 @@ def _cache_key(lat, lon): # shared VHF TX keeps Doppler inside the association gate. Spreading the # budget across metros puts overlap coverage where traffic actually flies. ring_nodes = [] + if layout == "dual": + # Dual-illuminator sites replace the coverage ring entirely: they are + # two different answers to the same problem. A ring buys geometry by + # surrounding the airspace with receivers that all overlap; a dual site + # buys it at the receiver, from two illuminators sharing one antenna. + if not metro_area: + # Without a core to place sites around, the layout silently + # returned a fleet ~n_cluster nodes short of --nodes. + raise ValueError("--layout dual requires --metro: dual sites are placed around a metro core") + if metro_area: + _dual_towers = [t for _r, t in available_towers] or _TOWERS_US + if illuminator_band == "vhf": + # All-VHF is a legitimate configuration to test on its own: + # VHF is the better illuminator on physics, and restricting to + # one band removes the cross-band question from the result. + _vhf = [t for t in _dual_towers if t[3] < 300e6] + if len(_vhf) >= 2: + _dual_towers = _vhf + ring_nodes.extend( + _generate_dual_sites( + n_sites=max(0, n_cluster) // 2, + core_lat=metro_area["lat"], + core_lon=metro_area["lon"], + towers=_dual_towers, + metro_radius_km=metro_area["radius_nm"] * _NM_TO_KM, + prefix=f"synth-{metro.strip().upper()}-DUAL" if metro else "synth-DUAL", + beam_width_deg=ring_beam_width_deg, + max_bistatic_range_km=ring_max_range_km, + min_eirp_dbm=dual_min_eirp_dbm, + aim=dual_aim, + ) + ) + return nodes + ring_nodes + if layout == "scatter": + # Like "dual", this replaces the ring rather than adding to it — the + # point is a fleet nobody placed, and a designed ring alongside it + # would carry the geometry the layout exists to do without. + if not metro_area: + raise ValueError("--layout scatter requires --metro: scatter sites are placed around a metro core") + if metro_area: + ring_nodes.extend( + _generate_scatter_sites( + n=max(0, n_cluster), + core_lat=metro_area["lat"], + core_lon=metro_area["lon"], + towers=[t for _r, t in available_towers] or _TOWERS_US, + metro_radius_km=metro_area["radius_nm"] * _NM_TO_KM, + prefix=f"synth-{metro.strip().upper()}-SCAT" if metro else "synth-SCAT", + beam_width_deg=ring_beam_width_deg, + max_bistatic_range_km=ring_max_range_km, + ) + ) + return nodes + ring_nodes + _dual_fraction_sites(n_dual_sites) for ring_id, spec, size in _active_rings(n_cluster, n_clusters, ring_spec): tx_lat, tx_lon, tx_alt_ft, fc_hz, callsign, core_lat, core_lon = spec ring_nodes.extend( @@ -979,13 +1697,26 @@ def _cache_key(lat, lon): ) nodes = ring_nodes + nodes # prepend so ring IDs are first - return nodes + # Appended last: every other node above has already consumed its RNG + # draws, so dual_fraction=0.0 (n_dual_sites=0) leaves that stream + # untouched and reproduces today's scene byte-for-byte at the same seed. + return nodes + _dual_fraction_sites(n_dual_sites) def fleet_summary(nodes: list[dict]) -> dict: """Compute a summary of the fleet configuration.""" from collections import Counter + if not nodes: + # min()/max() below raise on an empty fleet — report it instead. + return { + "total_nodes": 0, + "regions": {}, + "unique_towers": 0, + "towers_by_usage": {}, + "lat_range": None, + "lon_range": None, + } regions = Counter(n["region"] for n in nodes) towers = Counter(n["tx_callsign"] for n in nodes) return { @@ -1008,8 +1739,44 @@ def main(): parser = argparse.ArgumentParser(description="Generate fleet of synthetic node configs") parser.add_argument("--nodes", type=int, default=200, help="Number of nodes (100-1000)") parser.add_argument("--regions", type=str, default="us", help="Comma-separated regions: us,eu,au") + parser.add_argument( + "--metro", + type=str, + default=None, + choices=sorted(_KNOWN_METROS), + help="Restrict the whole fleet to one metro area (drops solo/rural nodes)", + ) parser.add_argument("--output", type=str, default="fleet_config.json", help="Output file path") parser.add_argument("--seed", type=int, default=42, help="Random seed") + parser.add_argument( + "--layout", + choices=("ring", "dual", "scatter"), + default="ring", + help="ring: receivers circling a core, all aimed inward. " + "dual: receivers scattered across the metro, each " + "running two nodes on two illuminators from one " + "antenna (n-cluster is the node budget, so half " + "that many sites). " + "scatter: an undesigned community fleet — sites " + "clumped where people live, each on its own nearby " + "illuminator, hand-aimed, heterogeneous hardware.", + ) + parser.add_argument( + "--illuminator-band", choices=("any", "vhf"), default="any", help="restrict dual-site illuminators to VHF" + ) + parser.add_argument( + "--dual-min-eirp-dbm", type=float, default=40.0, help="floor on the weaker illuminator of a dual pair" + ) + parser.add_argument( + "--dual-aim", + choices=("core", "random"), + default="core", + help="core: aim dual sites at the metro core with jitter. " + "random: scatter sectors to minimise inter-site " + "overlap — measured to starve the layout, because " + "85%% of traffic runs through the core and the solve " + "rate collapsed from 105 to 9.", + ) parser.add_argument( "--n-cluster", "--n-ring", @@ -1024,13 +1791,14 @@ def main(): dest="n_clusters", type=int, default=5, - help="Number of metro coverage rings (more = multinode spread across the map)", + help="Number of metro coverage rings (more = multinode spread across " + "the map). Capped at the number of rings that survive --metro.", ) parser.add_argument( "--ring-radius-km", type=float, default=18.0, help="Receiver ring radius around each metro core" ) parser.add_argument( - "--ring-beam-width-deg", type=float, default=50.0, help="Yagi half-power beamwidth for ring receivers" + "--ring-beam-width-deg", type=float, default=42.0, help="Yagi half-power beamwidth for ring receivers" ) parser.add_argument("--ring-max-range-km", type=float, default=60.0, help="Detection range for ring receivers") parser.add_argument( @@ -1040,6 +1808,15 @@ def main(): choices=["core", "broadside"], help="Aim ring beams at the metro core or broadside to TX", ) + parser.add_argument( + "--dual-fraction", + type=float, + default=0.0, + help="Fraction (0.0-1.0) of --nodes to additionally run as " + "dual-illuminator sites, carved out of the ring/metro " + "budget and appended after it. Requires --metro. " + "Ignored for --layout dual (already all-dual).", + ) args = parser.parse_args() regions = [r.strip().lower() for r in args.regions.split(",")] @@ -1049,18 +1826,33 @@ def main(): seed=args.seed, n_cluster=args.n_cluster, n_clusters=args.n_clusters, + layout=args.layout, + illuminator_band=args.illuminator_band, + dual_min_eirp_dbm=args.dual_min_eirp_dbm, + dual_aim=args.dual_aim, + dual_fraction=args.dual_fraction, ring_radius_km=args.ring_radius_km, ring_beam_width_deg=args.ring_beam_width_deg, ring_max_range_km=args.ring_max_range_km, ring_aim=args.ring_aim, + metro=args.metro, ) - cells = coverage_cells(n_cluster=args.n_cluster, n_clusters=args.n_clusters) + cells = coverage_cells(n_cluster=args.n_cluster, n_clusters=args.n_clusters, metro=args.metro, layout=args.layout) summary = fleet_summary(nodes) config = { "fleet": { "generated_at": __import__("datetime").datetime.now(__import__("datetime").timezone.utc).isoformat(), "summary": summary, + # Stamps the scene actually generated so the orchestrator can + # detect a drift between this and a polled backend config and + # self-restart for regeneration (see orchestrator._poll_simulation_config). + "scene": { + "n_nodes": args.nodes, + "dual_fraction": args.dual_fraction, + "layout": args.layout, + "seed": args.seed, + }, }, "nodes": nodes, "cells": cells, diff --git a/retina_simulation/node.py b/retina_simulation/node.py index 5356ac6..dacdc38 100644 --- a/retina_simulation/node.py +++ b/retina_simulation/node.py @@ -251,7 +251,14 @@ def _target_detection(self, target: SyntheticTarget) -> dict: delay += random.gauss(0, 0.1) # ~0.1 μs noise (GPS-disciplined SDR) doppler += random.gauss(0, 2.0) # ~2 Hz noise - # SNR depends on distance (closer = stronger) + # SNR depends on distance (closer = stronger). + # KNOWN SIMPLIFICATION: this is a one-way 10 dB/decade falloff on the + # RX-relative distance. Real bistatic received power goes as + # 1/(R_tx² · R_rx²) — ~40 dB/decade split across both legs — so + # synthetic SNR falls off far more gently than hardware will. Any + # SNR-derived gate tuned on this model needs re-tuning on real + # captures; changing the model is a measured follow-up, not a + # drive-by (it reshapes every detection threshold downstream). dist = _norm(pos) base_snr = 25 - 10 * math.log10(max(dist, 1)) snr = max(base_snr + random.gauss(0, 2), 4.0) @@ -871,6 +878,8 @@ def _stream_multi_node_tcp( fs_hz=nd.get("fs_hz", 2_000_000.0), beam_width_deg=nd.get("beam_width_deg", 41.0), max_range_km=nd.get("max_range_km", 50.0), + # Absent → monostatic range rule, so hardware nodes are unaffected. + max_bistatic_range_km=nd.get("max_bistatic_range_km"), ) world.add_node(wc) node_configs.append(wc) diff --git a/retina_simulation/orchestrator.py b/retina_simulation/orchestrator.py index d6790e7..26ab0b2 100644 --- a/retina_simulation/orchestrator.py +++ b/retina_simulation/orchestrator.py @@ -32,11 +32,20 @@ import time from datetime import datetime, timezone -from retina_simulation.generator import coverage_cells, fleet_summary, generate_fleet -from retina_simulation.tower_resolver import apply_tower_assignments, resolve_towers - # Add parent dir so we can import simulation packages -from retina_simulation.world import MetroCell, NodeConfig, SimulationWorld +from retina_simulation.generator import ( + _KNOWN_METROS, + coverage_cells, + fleet_summary, + generate_fleet, +) +from retina_simulation.tower_resolver import apply_tower_assignments, resolve_towers +from retina_simulation.world import ( + MetroCell, + NodeConfig, + SimulationWorld, + waypoints_for_metro, +) logging.basicConfig( level=logging.INFO, @@ -237,9 +246,11 @@ def __init__( hub_radial: bool = True, metro_traffic_frac: float = 0.6, cells: list[dict] | None = None, + metro: str | None = None, ): self.node_configs = node_configs self.cells = cells or [] + self.metro = metro self.host = host self.port = port self.mode = mode @@ -278,13 +289,22 @@ def _build_world(self): center_lat = sum(lats) / len(lats) center_lon = sum(lons) / len(lons) - self.world = SimulationWorld(center_lat=center_lat, center_lon=center_lon) - # Aircraft population is a nationwide field (alive map), independent of - # sensor count — floor it high so the whole US stays populated, not just - # the metros. --min-aircraft/--max-aircraft override. + self.world = SimulationWorld( + center_lat=center_lat, + center_lon=center_lon, + waypoints=waypoints_for_metro(self.metro), + ) + # Aircraft population is a field, not a per-sensor count. Nationwide it is + # floored high so the whole US stays populated; a metro-scoped fleet covers + # ~1/40th of that area, so the same floor would pack hundreds of aircraft + # into one terminal area. --min-aircraft/--max-aircraft override either way. n = len(self.node_configs) - auto_min_aircraft = max(150, n // 2) - auto_max_aircraft = max(300, n) + if self.metro: + auto_min_aircraft = max(15, n // 2) + auto_max_aircraft = max(30, n) + else: + auto_min_aircraft = max(150, n // 2) + auto_max_aircraft = max(300, n) self.world.min_aircraft = self.min_aircraft or auto_min_aircraft self.world.max_aircraft = self.max_aircraft or auto_max_aircraft if self.world.max_aircraft < self.world.min_aircraft: @@ -309,9 +329,13 @@ def _build_world(self): tx_alt_ft=cfg["tx_alt_ft"], fc_hz=cfg["fc_hz"], fs_hz=cfg.get("fs_hz", 2_000_000), - beam_width_deg=self.beam_width_deg or cfg.get("beam_width_deg", 40), + beam_width_deg=self.beam_width_deg or cfg.get("beam_width_deg", 42), max_range_km=self.max_range_km or cfg.get("max_range_km", 50), beam_azimuth_deg=cfg.get("beam_azimuth_deg"), # None → broadside in add_node + # Without this the world falls back to the monostatic RX-radius + # rule while the handshake tells the server the bistatic limit — + # nodes then "detect" 1.6x beyond what the server will accept. + max_bistatic_range_km=cfg.get("max_bistatic_range_km"), ) self.world.add_node(node) @@ -320,13 +344,15 @@ def _build_world(self): self.world.frac_metro_traffic = self.metro_traffic_frac log.info( - "SimulationWorld: center=(%.2f, %.2f), %d nodes, %d-%d aircraft, %d metro cells", + "SimulationWorld: center=(%.2f, %.2f), %d nodes, %d-%d aircraft, %d metro cells, waypoint net=%s (%d)", center_lat, center_lon, len(self.node_configs), self.world.min_aircraft, self.world.max_aircraft, len(self.world.metro_cells), + self.metro or "nationwide", + len(self.world.waypoints), ) async def _connect_batch(self, configs: list[dict]) -> list[dict]: @@ -642,6 +668,35 @@ def save_ground_truth(self, path: str): log.info("Ground truth saved: %s (%d snapshots)", path, len(self.ground_truth)) +def build_ground_truth_payload(aircraft_summaries: list[dict]) -> list[dict]: + """Remap world aircraft summaries to the server ground-truth push schema. + + Dark objects have no ADS-B hex, so their stable object id doubles as the + ground-truth key; entries with neither are unidentifiable and dropped. + """ + payload_aircraft = [] + for ac in aircraft_summaries: + hex_code = ac.get("adsb_hex") or ac.get("id", "") + if not hex_code: + continue + payload_aircraft.append( + { + "hex": hex_code, + "lat": ac["lat"], + "lon": ac["lon"], + "alt_m": ac["alt_km"] * 1000, + "heading": ac.get("heading", 0), + "speed_ms": ac.get("speed_ms", 0), + "object_type": ac.get("object_type", "aircraft"), + "is_anomalous": ac.get("is_anomalous", False), + "has_adsb": ac.get("has_adsb", False), + "adsb_callsign": ac.get("adsb_callsign") or None, + "anomaly_event": ac.get("anomaly_event") or None, + } + ) + return payload_aircraft + + async def _push_ground_truth_live( orchestrator: FleetOrchestrator, base_url: str, @@ -669,24 +724,7 @@ async def _push_ground_truth_live( try: aircraft = orchestrator.world.get_aircraft_summary() - # Remap field names to what the server endpoint expects - payload_aircraft = [] - for ac in aircraft: - hex_code = ac.get("adsb_hex") or ac.get("id", "") - if not hex_code: - continue - payload_aircraft.append( - { - "hex": hex_code, - "lat": ac["lat"], - "lon": ac["lon"], - "alt_m": ac["alt_km"] * 1000, - "heading": ac.get("heading", 0), - "speed_ms": ac.get("speed_ms", 0), - "object_type": ac.get("object_type", "aircraft"), - "is_anomalous": ac.get("is_anomalous", False), - } - ) + payload_aircraft = build_ground_truth_payload(aircraft) if payload_aircraft: body = json.dumps( @@ -802,8 +840,19 @@ async def _poll_simulation_config( orchestrator: FleetOrchestrator, base_url: str, interval_s: float = 5.0, + scene: dict | None = None, ): - """Poll /api/simulation/config every interval_s and apply updated spawn fractions to world.""" + """Poll /api/simulation/config every interval_s and apply updated spawn fractions to world. + + Also watches for a scene change (n_nodes / dual_fraction): those are + baked into the fleet at container boot (see generator.py), so unlike the + spawn fractions above they cannot be applied in-process. When ``scene`` + (this container's own stamped {n_nodes, dual_fraction, layout, seed}) is + present and differs from what the backend now reports, this logs a loud + WARN and calls ``orchestrator.stop()`` — every loop then exits, the + process exits 0, and `restart: unless-stopped` relaunches into + fleet-entrypoint.sh, which fetches the desired scene before regenerating. + """ import urllib.request log.info("Simulation config polling started (url=%s, interval=%.1fs)", base_url, interval_s) @@ -831,8 +880,10 @@ def _fetch(): cfg = await loop.run_in_executor(None, _fetch) updated_at = cfg.get("_updated_at", 0.0) if updated_at > last_updated_at: - orchestrator.world.frac_anomalous = float(cfg.get("frac_anomalous", 0.05)) - orchestrator.world.frac_drone = float(cfg.get("frac_drone", 0.10)) + # Fallback matches SimulationWorld's default (anomalies off), so a + # payload missing the key cannot silently switch them back on. + orchestrator.world.frac_anomalous = float(cfg.get("frac_anomalous", 0.0)) + orchestrator.world.frac_drone = float(cfg.get("frac_drone", 0.0)) orchestrator.world.frac_dark = float(cfg.get("frac_dark", 0.15)) if "min_aircraft" in cfg: orchestrator.world.min_aircraft = int(cfg["min_aircraft"]) @@ -847,6 +898,44 @@ def _fetch(): orchestrator.world.min_aircraft, orchestrator.world.max_aircraft, ) + + # Scene-change detection. Absent scene stamp (stale volume, + # in-process generation) or absent config keys (never PUT) → + # no comparison, no restart. + scene_n_nodes = scene.get("n_nodes") if scene else None + scene_dual_fraction = scene.get("dual_fraction") if scene else None + scene_diff = False + if scene: + if "n_nodes" in cfg and scene_n_nodes is not None and int(cfg["n_nodes"]) != int(scene_n_nodes): + scene_diff = True + if ( + "dual_fraction" in cfg + and scene_dual_fraction is not None + and abs(float(cfg["dual_fraction"]) - float(scene_dual_fraction)) > 1e-6 + ): + scene_diff = True + # max_range_km needs no stamp: the orchestrator itself holds + # the running value, so this runs even when `scene` is None. + # Applying a range change requires regenerating node configs + # (every node's cfg is built from it at construction), which + # is exactly the restart path — the poll loop deliberately + # does NOT live-apply it. + if "max_range_km" in cfg and abs(float(cfg["max_range_km"]) - float(orchestrator.max_range_km)) > 1e-6: + scene_diff = True + if scene_diff: + log.warning( + "Scene change requested (n_nodes=%s dual_fraction=%s " + "max_range_km=%s, running n_nodes=%s dual_fraction=%s " + "max_range_km=%s) — shutting down for regeneration", + cfg.get("n_nodes"), + cfg.get("dual_fraction"), + cfg.get("max_range_km"), + scene_n_nodes, + scene_dual_fraction, + orchestrator.max_range_km, + ) + await orchestrator.stop() + return except Exception as e: log.debug("Config poll failed: %s", e) @@ -984,21 +1073,6 @@ def _get_json(endpoint_url): log.debug("Validation check failed: %s", e) -# ── Predefined metro areas ────────────────────────────────────────────────── -_KNOWN_METROS = { - "atl": {"name": "Atlanta", "lat": 33.749, "lon": -84.388, "radius_nm": 80}, - "gvl": {"name": "Greenville", "lat": 34.852, "lon": -82.394, "radius_nm": 60}, - "clt": {"name": "Charlotte", "lat": 35.227, "lon": -80.843, "radius_nm": 70}, - "nyc": {"name": "New York", "lat": 40.748, "lon": -73.986, "radius_nm": 80}, - "dca": {"name": "Washington DC", "lat": 38.935, "lon": -77.079, "radius_nm": 70}, - "chi": {"name": "Chicago", "lat": 41.872, "lon": -87.624, "radius_nm": 80}, - "den": {"name": "Denver", "lat": 39.739, "lon": -104.990, "radius_nm": 80}, - "lax": {"name": "Los Angeles", "lat": 34.052, "lon": -118.244, "radius_nm": 80}, - "dfw": {"name": "Dallas-Fort Worth", "lat": 32.897, "lon": -97.038, "radius_nm": 80}, - "kc": {"name": "Kansas City", "lat": 39.298, "lon": -94.714, "radius_nm": 70}, -} - - def _parse_metro_areas(metros_str: str) -> list[dict]: """Parse comma-separated metro codes into area dicts for AdsbLolClient.""" result = [] @@ -1014,6 +1088,15 @@ def _parse_metro_areas(metros_str: str) -> list[dict]: async def main_async(args): """Main async entry point.""" # Load or generate fleet config + # Scene stamp (n_nodes/dual_fraction/layout/seed) the generator wrote into + # the config it produced — read here so _poll_simulation_config can detect + # a backend-requested scene change and self-restart for regeneration. + # Only ever populated via the loaded-config-file path: fleet-entrypoint.sh + # always runs the generator CLI (which stamps it) before starting us with + # --config. None here means "no comparison" — a stale volume with no + # stamp, or the in-process generation fallback below, never triggers a + # restart loop. + scene = None if args.config and os.path.exists(args.config): with open(args.config) as f: data = json.load(f) @@ -1022,6 +1105,7 @@ async def main_async(args): # Fallback: maybe it's the old nodes_config.json format all_nodes = data.get("nodes", []) cells = data.get("cells", []) + scene = data.get("fleet", {}).get("scene") else: log.info("No config file, generating %d nodes...", args.nodes) regions = [r.strip() for r in args.regions.split(",")] @@ -1031,8 +1115,14 @@ async def main_async(args): seed=args.seed, n_cluster=args.n_cluster, n_clusters=args.n_clusters, + metro=getattr(args, "metro", None), + ) + cells = coverage_cells( + n_cluster=args.n_cluster, + n_clusters=args.n_clusters, + metro=getattr(args, "metro", None), + layout=getattr(args, "layout", "ring"), ) - cells = coverage_cells(n_cluster=args.n_cluster, n_clusters=args.n_clusters) # When --metros is specified, filter nodes to only those near selected metros if getattr(args, "metros", "") and args.metros: @@ -1091,6 +1181,7 @@ def _near_any_metro(node): hub_radial=not args.no_hub_radial, metro_traffic_frac=args.metro_traffic_frac, cells=cells, + metro=getattr(args, "metro", None), ) # Build shared simulation world @@ -1143,12 +1234,16 @@ def _near_any_metro(node): orchestrator, args.validation_url, interval_s=5.0, + scene=scene, ) ) # Real ADS-B from adsb.lol — inject real air traffic when metro areas are configured. - if args.validation_url and hasattr(args, "metros") and args.metros: - metro_areas = _parse_metro_areas(args.metros) + # --metro (generation-time scoping) implies the same area for real traffic, so a + # Greenville-only fleet gets Greenville ADS-B without also passing --metros. + adsb_metros = getattr(args, "metros", "") or getattr(args, "metro", "") or "" + if args.validation_url and adsb_metros: + metro_areas = _parse_metro_areas(adsb_metros) if metro_areas: tasks.append( _push_real_adsb( @@ -1210,7 +1305,8 @@ def main(): default=5, help="Number of distinct metro rings to fan the --n-cluster budget " "across (5 = Dallas, Chicago, Atlanta, Denver, Kansas City; " - "1 = single Dallas ring). Matches the generator default.", + "1 = single Dallas ring). Capped at the number of rings that " + "survive --metro. Matches the generator default.", ) parser.add_argument("--host", type=str, default="localhost", help="Server hostname") parser.add_argument("--port", type=int, default=3012, help="Server TCP port") @@ -1246,12 +1342,22 @@ def main(): parser.add_argument( "--ground-truth-path", type=str, default="ground_truth.json", help="Path to save ground truth data" ) + parser.add_argument( + "--metro", + type=str, + default=None, + choices=sorted(_KNOWN_METROS), + help="Generate the whole fleet inside one metro area (drops " + "solo/rural nodes and non-local rings). Applies at " + "generation time, so it needs no --config.", + ) parser.add_argument( "--metros", type=str, default="", help="Comma-separated metro codes to focus on (e.g. atl,gvl). " - "Filters fleet to these metros and injects real ADS-B from adsb.lol. " + "Filters an already-generated fleet to these metros and " + "injects real ADS-B from adsb.lol. " f"Available: {','.join(_KNOWN_METROS.keys())}", ) parser.add_argument( diff --git a/retina_simulation/world.py b/retina_simulation/world.py index 10a6ec0..5366e55 100644 --- a/retina_simulation/world.py +++ b/retina_simulation/world.py @@ -73,6 +73,40 @@ (45.5898, -122.5951), # PDX Portland ] +# ── Regional waypoint nets ──────────────────────────────────────────────────── +# A metro-scoped fleet has no receivers outside its own metro, so cross-country +# en-route traffic is pure waste: it burns simulation budget on aircraft no node +# can ever see, and puts ground-truth tracks on the map thousands of km from the +# only coverage that exists. Selecting a regional net keeps the background +# traffic inside the region the fleet actually covers. +# +# Keyed by the same metro codes as generator._KNOWN_METROS. +_REGIONAL_WAYPOINTS: dict[str, list[tuple[float, float]]] = { + "gvl": [ + (34.8957, -82.2189), # GSP Greenville-Spartanburg + (34.8479, -82.3499), # GMU Greenville Downtown + (34.9157, -81.9565), # SPA Spartanburg Downtown Memorial + (34.4946, -82.7093), # AND Anderson Regional + (35.4362, -82.5418), # AVL Asheville + (35.2144, -80.9473), # CLT Charlotte + (34.8964, -81.0572), # RKH Rock Hill + (33.9388, -81.1195), # CAE Columbia + (34.4984, -81.9573), # GRD Greenwood + (35.7565, -81.6790), # HKY Hickory + ], +} + + +def waypoints_for_metro(metro: str | None) -> list[tuple[float, float]]: + """Waypoint net for a metro code, falling back to the nationwide net. + + An unknown or absent code returns _US_WAYPOINTS, so callers that don't know + about regional scoping keep the original coast-to-coast behaviour. + """ + if not metro: + return _US_WAYPOINTS + return _REGIONAL_WAYPOINTS.get(metro.strip().lower(), _US_WAYPOINTS) + @dataclass class SimulatedAircraft: @@ -90,6 +124,12 @@ class SimulatedAircraft: # Heading (degrees from north, clockwise) heading_deg: float speed_km_s: float + # Cruise speed the in-trail separation modulation recovers toward once a + # conflict clears (speed_km_s is the CURRENT speed, which separation and + # anomaly events both mutate). 0.0 = "not set" for pre-existing callers + # that construct aircraft directly; separation then treats the current + # speed as cruise. + base_speed_km_s: float = 0.0 # Type and classification has_adsb: bool = False is_anomalous: bool = False @@ -99,6 +139,9 @@ class SimulatedAircraft: # Lifecycle created_at: float = 0.0 lifetime_s: float = 600.0 + # Expired and rerouted toward the region edge for retirement (see + # SimulationWorld._route_out); guards the reroute from re-firing. + departing: bool = False # Waypoint navigation waypoints: list = field(default_factory=list) waypoint_idx: int = 0 @@ -128,8 +171,16 @@ class NodeConfig: min_doppler: float = 15.0 # Detection geometry beam_azimuth_deg: float | None = None # None → auto broadside in add_node - beam_width_deg: float = 41.0 # Yagi half-power beamwidth (40-42° spec) - max_range_km: float = 50.0 # maximum detection range + beam_width_deg: float = 42.0 # Yagi half-power beamwidth (fleet spec) + max_range_km: float = 50.0 # maximum RX→target range (monostatic) + # Maximum *bistatic* range: (RX→target) + (target→TX) − baseline, i.e. the + # differential range the delay measurement actually represents, and what + # sets received power via the bistatic radar equation. Physically the + # correct limit — it makes the footprint an ellipse with foci at RX and TX + # rather than a circle around the RX. + # None keeps the older monostatic rule, so real hardware nodes carrying + # only max_range_km are unaffected. + max_bistatic_range_km: float | None = None def config_hash(config: NodeConfig) -> str: @@ -225,12 +276,37 @@ def _bistatic_doppler(target_enu, vel_enu, tx_enu, rx_enu, freq_hz): # ── Flight corridor route generation ───────────────────────────────────────── +_NATIONWIDE_MIN_LEG_KM = 400.0 + -def _pick_route(center_lat: float, center_lon: float, max_dist_km: float = 300) -> list[tuple[float, float]]: +def _min_leg_km(waypoints: list[tuple[float, float]]) -> float: + """A "this is a real en-route leg" threshold scaled to a waypoint net. + + Scaled from the net's own extent, capped at the nationwide 400 km so the + continent-wide net keeps its original behaviour exactly. A regional net has + no pair 400 km apart, so without the scaling every destination would be + rejected and the fallback would hand back coast-to-coast routes — the exact + thing regional scoping exists to prevent. + """ + if len(waypoints) < 2: + return 0.0 + lats = [wp[0] for wp in waypoints] + lons = [wp[1] for wp in waypoints] + span_km = _haversine_km(min(lats), min(lons), max(lats), max(lons)) + return min(_NATIONWIDE_MIN_LEG_KM, 0.35 * span_km) + + +def _pick_route( + center_lat: float, + center_lon: float, + max_dist_km: float = 300, + waypoints: list[tuple[float, float]] | None = None, +) -> list[tuple[float, float]]: """Pick a sequence of 2-4 waypoints near center forming a realistic route.""" - nearby = [wp for wp in _US_WAYPOINTS if _haversine_km(center_lat, center_lon, wp[0], wp[1]) < max_dist_km] + waypoints = waypoints if waypoints is not None else _US_WAYPOINTS + nearby = [wp for wp in waypoints if _haversine_km(center_lat, center_lon, wp[0], wp[1]) < max_dist_km] if len(nearby) < 2: - nearby = sorted(_US_WAYPOINTS, key=lambda wp: _haversine_km(center_lat, center_lon, wp[0], wp[1]))[:6] + nearby = sorted(waypoints, key=lambda wp: _haversine_km(center_lat, center_lon, wp[0], wp[1]))[:6] n_waypoints = random.randint(2, min(4, len(nearby))) start = random.choice(nearby) @@ -256,9 +332,18 @@ def _pick_route(center_lat: float, center_lon: float, max_dist_km: float = 300) class SimulationWorld: """Shared simulation world with aircraft and multiple observer nodes.""" - def __init__(self, center_lat: float = 34.0, center_lon: float = -84.0): + def __init__( + self, + center_lat: float = 34.85, + center_lon: float = -82.39, + waypoints: list[tuple[float, float]] | None = None, + ): self.center_lat = center_lat self.center_lon = center_lon + # En-route waypoint net for background traffic. Defaults to the + # nationwide list; set a regional net (waypoints_for_metro) to keep + # background aircraft inside a metro-scoped fleet's coverage. + self.waypoints = waypoints if waypoints is not None else _US_WAYPOINTS self.aircraft: list[SimulatedAircraft] = [] self.nodes: dict[str, NodeConfig] = {} self._next_id = 1 @@ -266,9 +351,21 @@ def __init__(self, center_lat: float = 34.0, center_lon: float = -84.0): # Target count range self.min_aircraft = 5 self.max_aircraft = 15 - # Object type spawn fractions (adjustable at runtime) - self.frac_anomalous: float = 0.05 - self.frac_drone: float = 0.10 + # Object type spawn fractions (adjustable at runtime). + # frac_anomalous doubles as the master gate for ALL anomaly generation: + # at 0 it also suppresses the mid-flight anomaly scheduler + # (_maybe_schedule_anomaly), which would otherwise keep turning normal + # commercial aircraft anomalous at its own hardcoded rate. One switch, + # so "anomalies off" means off — and raising it turns everything back on. + self.frac_anomalous: float = 0.0 + # Drones off by default, matching the backend's simulation_config + # default (user call, 2026-08: fixed-wing scene only). The old 0.10 + # here spawned a handful of drones in the window between world + # construction and the first config poll on EVERY fleet restart — + # visible as "a few drones exist even with the slider at 0" until + # they aged out. The orchestrator poll's absent-key fallback must + # stay matched to this value. + self.frac_drone: float = 0.0 self.frac_dark: float = 0.15 # remaining fraction = commercial aircraft with ADS-B # Hub-radial flight planning: when metro_cells is non-empty, this @@ -278,6 +375,19 @@ def __init__(self, center_lat: float = 34.0, center_lon: float = -84.0): self.metro_cells: list[MetroCell] = [] self.frac_metro_traffic: float = 0.6 self.arrival_departure_overflight_weights = (0.45, 0.35, 0.20) + # Traffic separation. Real traffic never stacks: terminal in-trail + # minima are ~3 NM and crossing flows are altitude-split, but the + # hub-radial planner converges every metro spawn on the same ~2 km + # core with no deconfliction at all, so the fleet routinely flew + # pairs inside the solver's association gates (delay gate ≈ 1-3 km) + # — making association ambiguity a property of the SIMULATOR, not of + # realistic traffic. Two mechanisms, both gated on these knobs: + # spawn poses resample away from live traffic, and in-flight + # conflicts resolve by slowing the later-created aircraft + # (_enforce_separation). A pair needs BOTH horizontal and vertical + # proximity to count as a conflict. + self.min_separation_km: float = 5.0 + self.min_vertical_sep_km: float = 0.6 # ~2000 ft def add_node(self, config: NodeConfig): """Register a synthetic node in the simulation. @@ -316,12 +426,14 @@ def _choose_spawn_pose(self) -> tuple[float, float, list]: return self._nationwide_pose() def _nationwide_pose(self) -> tuple[float, float, list]: - """Cross-country en-route traffic on the national waypoint net, spawned - anywhere along the leg (not just at airports) and independent of node - placement — the nationwide background that keeps the map alive.""" - start = random.choice(_US_WAYPOINTS) - far = [wp for wp in _US_WAYPOINTS if _haversine_km(start[0], start[1], wp[0], wp[1]) > 400] - dest = random.choice(far or _US_WAYPOINTS) + """En-route traffic on self.waypoints, spawned anywhere along the leg (not + just at airports) and independent of node placement — the background that + keeps the map alive. Nationwide by default; regional under metro scoping.""" + net = self.waypoints + start = random.choice(net) + min_leg = _min_leg_km(net) + far = [wp for wp in net if _haversine_km(start[0], start[1], wp[0], wp[1]) > min_leg] + dest = random.choice(far or net) t = random.uniform(0.0, 1.0) lat = start[0] + t * (dest[0] - start[0]) + random.gauss(0, 0.3) lon = start[1] + t * (dest[1] - start[1]) + random.gauss(0, 0.3) @@ -376,25 +488,35 @@ def _fallback_pose(self) -> tuple[float, float, list]: baseline_bearing = _bearing_deg(anchor.rx_lat, anchor.rx_lon, anchor.tx_lat, anchor.tx_lon) perp_rad = math.radians((baseline_bearing + 90.0) % 360.0) dist_km = random.uniform(5.0, anchor.max_range_km * 0.7) - anchor_lat = anchor.rx_lat + (dist_km * math.cos(perp_rad)) / 111.32 - cos_lat = math.cos(math.radians(anchor.rx_lat)) - anchor_lon = anchor.rx_lon + (dist_km * math.sin(perp_rad)) / (111.32 * max(cos_lat, 1e-6)) + # R_EARTH-derived like every other conversion in this file — these + # were the last 111.32 literals, 0.11% off the rest of the sim. + anchor_lat, anchor_lon, _ = _enu_to_lla( + dist_km * math.sin(perp_rad), + dist_km * math.cos(perp_rad), + 0.0, + anchor.rx_lat, + anchor.rx_lon, + 0.0, + ) else: anchor_lat, anchor_lon = self.center_lat, self.center_lon lat = anchor_lat + random.gauss(0, 0.03) lon = anchor_lon + random.gauss(0, 0.03) - route = [(lat, lon)] + _pick_route(anchor_lat, anchor_lon, max_dist_km=50) + route = [(lat, lon)] + _pick_route(anchor_lat, anchor_lon, max_dist_km=50, waypoints=self.waypoints) return lat, lon, route def _spawn_aircraft(self, mode: str = "detection") -> SimulatedAircraft: """Spawn a new aircraft along a realistic flight corridor. - Object types are selected probabilistically: - - 70% commercial aircraft (with ADS-B in adsb/anomalous modes) - - 15% dark aircraft (no ADS-B transponder) - - 10% drones (low/slow) - - 5% anomalous objects (erratic behavior) + Object types are selected probabilistically from the frac_* instance + attributes: dark aircraft (no transponder), drones (low/slow), anomalous + objects (erratic), and commercial aircraft with ADS-B as the remainder. + frac_anomalous defaults to 0 — see __init__. + + mode="anomalous" is an explicit opt-in that injects anomalies regardless + of frac_anomalous; it is a testing mode and is not used by any deployment + (every compose profile sets FLEET_MODE=adsb). """ oid = f"obj-{self._next_id:05d}" self._next_id += 1 @@ -417,7 +539,21 @@ def _spawn_aircraft(self, mode: str = "detection") -> SimulatedAircraft: else: object_type = "aircraft" # commercial — will get ADS-B in adsb modes + # Resample the spawn pose away from live traffic — best-effort, never + # fails: after the attempt budget the farthest candidate wins, so a + # saturated core degrades to "as separated as the planner could get" + # rather than blocking spawns (step() spawns in a while-loop up to + # min_aircraft; a hard reject there would spin forever). lat, lon, route = self._choose_spawn_pose() + best = (self._nearest_traffic_km(lat, lon), lat, lon, route) + for _ in range(9): + if best[0] >= self.min_separation_km: + break + lat, lon, route = self._choose_spawn_pose() + d = self._nearest_traffic_km(lat, lon) + if d > best[0]: + best = (d, lat, lon, route) + _, lat, lon, route = best if mode == "anomalous" and random.random() < 0.2: is_anomalous = True @@ -436,7 +572,15 @@ def _spawn_aircraft(self, mode: str = "detection") -> SimulatedAircraft: # Anomalous objects also get ADS-B — anomalous means unusual flight # behaviour (speed/altitude/heading changes), NOT transponder absence. - if mode in ("adsb", "anomalous") and object_type != "drone" and roll >= 0.30 or is_anomalous: + # + # The floor must be the SAME cumulative boundary the type roll used for + # "commercial" above. It was previously hardcoded to 0.30, which silently + # assumed the original defaults (0.05 + 0.10 + 0.15). Any other fractions + # — including anything set through the Physics Settings slider at + # runtime — pushed part of the commercial band below the literal, so + # those aircraft were spawned with no transponder and showed up as dark. + adsb_roll_floor = self.frac_anomalous + self.frac_drone + self.frac_dark + if (mode in ("adsb", "anomalous") and object_type != "drone" and roll >= adsb_roll_floor) or is_anomalous: has_adsb = True adsb_hex = f"{random.randint(0x100000, 0xFFFFFF):06x}" letters = "ABCDEFGHIJKLMNOPQRSTUVWXYZ" @@ -461,6 +605,7 @@ def _spawn_aircraft(self, mode: str = "detection") -> SimulatedAircraft: vel_up=vel_up, heading_deg=heading, speed_km_s=speed_km_s, + base_speed_km_s=speed_km_s, has_adsb=has_adsb, is_anomalous=is_anomalous, object_type=object_type, @@ -473,12 +618,70 @@ def _spawn_aircraft(self, mode: str = "detection") -> SimulatedAircraft: **self._maybe_schedule_anomaly(is_anomalous, object_type), ) + # Expired aircraft are retired only once they are at least this far from + # the world center — a target vanishing overhead reads as a tracking bug + # on the map, so retirement happens off at the edges, beyond the ~60 km + # node coverage of a metro-scoped fleet. (A nationwide fleet spawns most + # aircraft beyond this radius, which degrades to the old expire-anywhere + # behaviour — acceptable, nothing deployed runs unscoped.) + retire_edge_km: float = 70.0 + # Extra flight time an expired aircraft gets to actually REACH the edge. + # Metro-routed traffic (85% of spawns) circulates inside the metro and + # rarely crosses retire_edge_km on its own, so under the old 2x-lifetime + # hard cap most of it still vanished mid-view — the edge gate only ever + # covered planes that happened to be leaving anyway. On expiry the + # aircraft is now rerouted outward (_route_out) and 70 km at the slowest + # commercial speed (~0.1 km/s) is ~700 s, so 900 s means the backstop + # below fires only for genuinely stuck aircraft. + exit_grace_s: float = 900.0 + + def _route_out(self, ac: "SimulatedAircraft") -> None: + """Point an expired aircraft at the region edge and let it fly out. + + Replaces the remaining route with one waypoint past retire_edge_km on + the bearing away from the world center (current heading when the + aircraft sits at the center itself), so retirement is something the + viewer watches happen at the edge instead of a dot blinking out. + """ + if _haversine_km(self.center_lat, self.center_lon, ac.lat, ac.lon) > 1.0: + brg = _bearing_deg(self.center_lat, self.center_lon, ac.lat, ac.lon) + else: + brg = ac.heading_deg + dist_km = self.retire_edge_km + 15.0 + brg_rad = math.radians(brg) + wp_lat = self.center_lat + dist_km * math.cos(brg_rad) / 111.32 + wp_lon = self.center_lon + dist_km * math.sin(brg_rad) / (111.32 * math.cos(math.radians(self.center_lat))) + ac.waypoints = [(wp_lat, wp_lon)] + ac.waypoint_idx = 0 + ac.departing = True + + def _should_retire(self, ac: "SimulatedAircraft") -> bool: + age = self._time - ac.created_at + if age < ac.lifetime_s: + return False + if ac.object_type == "drone": + # Drones loop low and slow and are expected to churn; an amber + # X-frame vanishing reads as turnover, not a tracking bug. + if age > ac.lifetime_s * 2: + return True + elif age > ac.lifetime_s + self.exit_grace_s: + # Leak backstop only — a departing aircraft normally crosses the + # edge well inside the grace window. + return True + return _haversine_km(self.center_lat, self.center_lon, ac.lat, ac.lon) > self.retire_edge_km + def step(self, dt: float, mode: str = "detection"): """Advance simulation by dt seconds.""" self._time += dt - # Remove expired aircraft - self.aircraft = [ac for ac in self.aircraft if (self._time - ac.created_at) < ac.lifetime_s] + # Expired non-drones head for the edge before the retire filter sees + # them past the edge (see _route_out / _should_retire). + for ac in self.aircraft: + if not ac.departing and ac.object_type != "drone" and self._time - ac.created_at >= ac.lifetime_s: + self._route_out(ac) + + # Retire expired aircraft (edge-gated — see _should_retire) + self.aircraft = [ac for ac in self.aircraft if not self._should_retire(ac)] # Spawn to maintain target count while len(self.aircraft) < self.min_aircraft: @@ -490,13 +693,64 @@ def step(self, dt: float, mode: str = "detection"): for ac in self.aircraft: self._update_aircraft(ac, dt) + self._enforce_separation(dt) + + def _nearest_traffic_km(self, lat: float, lon: float) -> float: + """Horizontal distance to the closest live aircraft (inf when none). + + Horizontal-only on purpose: spawn altitude is rolled after the pose, + so a vertical allowance here would be checked against a value that + does not exist yet — and a stricter horizontal-only bubble at spawn + costs nothing (the resample budget absorbs it).""" + if not self.aircraft: + return float("inf") + return min(_haversine_km(lat, lon, ac.lat, ac.lon) for ac in self.aircraft) + + def _enforce_separation(self, dt: float): + """In-trail speed modulation: the later-created aircraft of any + conflicting pair slows toward 70% of cruise until the conflict clears, + then recovers toward cruise. + + Speed-only, never heading: heading belongs to the waypoint router, and + a lateral dodge here would fight it every tick. Slowing the trailer + resolves both in-trail stacking (the leader pulls away) and crossing + conflicts (the leader crosses first) — the same sequencing terminal + control actually applies. Exemptions are the point, not an + optimisation: anomalous aircraft keep their erratic close approaches + (that proximity IS the signature the radar network exists to catch), + and drones live below min_vertical_sep_km of the jet flow anyway. + The blend is dt-based (~3 s time constant) so speed changes read as + throttle on the Doppler channel, not steps.""" + flow = [ac for ac in self.aircraft if not ac.is_anomalous and ac.object_type == "aircraft"] + flow.sort(key=lambda ac: ac.created_at) + slowed: set[str] = set() + for i, trail in enumerate(flow): + for lead in flow[:i]: + if abs(lead.alt_km - trail.alt_km) >= self.min_vertical_sep_km: + continue + if _haversine_km(lead.lat, lead.lon, trail.lat, trail.lon) < self.min_separation_km: + slowed.add(trail.object_id) + break + blend = 1.0 - math.exp(-dt / 3.0) + for ac in flow: + base = ac.base_speed_km_s or ac.speed_km_s + target = base * 0.7 if ac.object_id in slowed else base + ac.speed_km_s += (target - ac.speed_km_s) * blend + # ── Mid-flight anomaly scheduling ──────────────────────────────────────── _ANOMALY_EVENTS = ["hijack", "spoof", "orbit", "altitude_jump", "id_swap"] def _maybe_schedule_anomaly(self, is_anomalous: bool, object_type: str) -> dict: """Return kwargs to schedule a mid-flight anomaly event on ~8% of - normal commercial aircraft. Already-anomalous or drones are skipped.""" + normal commercial aircraft. Already-anomalous or drones are skipped. + + Gated on frac_anomalous: this rate is larger than the spawn-time + fraction, so without the gate, zeroing frac_anomalous would still leave + the majority of anomalies running (they just appear 30-120s late). + """ + if self.frac_anomalous <= 0: + return {} if is_anomalous or object_type == "drone": return {} if random.random() > 0.08: @@ -524,8 +778,17 @@ def _update_aircraft(self, ac: SimulatedAircraft, dt: float): ac.lat += dlat ac.lon += dlon ac.alt_km += ac.vel_up * dt + # Level off: vel_up was set once at spawn and integrated forever, so + # every aircraft eventually saturated against the altitude clamps — + # the whole fleet ended up pinned to the floor or the ceiling. Unlike + # vel_east/north (recomputed from heading each tick), there is no + # vertical navigation, so decay toward level flight with a ~5 min + # time constant and stop climbing at the clamp. + ac.vel_up *= math.exp(-dt / 300.0) # Clamp altitude + if ac.alt_km <= 0.1 or ac.alt_km >= 15.0: + ac.vel_up = 0.0 ac.alt_km = max(0.1, min(ac.alt_km, 15.0)) # ── Orbit anomaly: circle in place instead of following waypoints ──── @@ -624,10 +887,45 @@ def _fire_anomaly_event(self, ac: SimulatedAircraft): ac.speed_km_s = random.uniform(0.35, 0.50) def _aircraft_in_detection_cone(self, ac: SimulatedAircraft, node: NodeConfig) -> bool: - """Check if aircraft is within the node's detection cone.""" - dist = _haversine_km(node.rx_lat, node.rx_lon, ac.lat, ac.lon) - if dist > node.max_range_km: - return False + """Check if aircraft is within the node's detection cone. + + Range is limited on *bistatic* range when the node declares one — the + sum of both legs minus the baseline, which is what the delay actually + measures and what sets received power. A monostatic RX-distance limit + ignores the TX leg entirely, so it accepts targets far behind the + transmitter and rejects near ones on a long baseline. Nodes without + max_bistatic_range_km keep the monostatic rule so real hardware is + unaffected. + """ + if node.max_bistatic_range_km is not None: + rx_alt_km = node.rx_alt_ft * 0.3048 / 1000.0 + tx_alt_km = node.tx_alt_ft * 0.3048 / 1000.0 + target_enu = _lla_to_enu( + ac.lat, + ac.lon, + ac.alt_km, + node.rx_lat, + node.rx_lon, + rx_alt_km, + ) + tx_enu = _lla_to_enu( + node.tx_lat, + node.tx_lon, + tx_alt_km, + node.rx_lat, + node.rx_lon, + rx_alt_km, + ) + # _bistatic_delay returns the differential range in µs; multiply + # back by c to compare in km. Reused rather than recomputing the + # two legs so the gate and the emitted delay can never disagree. + diff_range_km = _bistatic_delay(target_enu, tx_enu, (0.0, 0.0, 0.0)) * C_KM_US + if diff_range_km > node.max_bistatic_range_km: + return False + else: + dist = _haversine_km(node.rx_lat, node.rx_lon, ac.lat, ac.lon) + if dist > node.max_range_km: + return False bearing = _bearing_deg(node.rx_lat, node.rx_lon, ac.lat, ac.lon) angle_diff = abs((bearing - node.beam_azimuth_deg + 180) % 360 - 180) @@ -770,6 +1068,8 @@ def get_aircraft_summary(self) -> list[dict]: "is_anomalous": ac.is_anomalous, "object_type": ac.object_type, "adsb_hex": ac.adsb_hex, + "adsb_callsign": ac.adsb_callsign, + "anomaly_event": ac.anomaly_event, } for ac in self.aircraft ] diff --git a/tests/test_bistatic_range.py b/tests/test_bistatic_range.py new file mode 100644 index 0000000..b55061e --- /dev/null +++ b/tests/test_bistatic_range.py @@ -0,0 +1,130 @@ +"""Detection range limited on bistatic range rather than RX distance. + +A monostatic limit compares only the RX->target leg against max_range_km, +ignoring the transmitter entirely. That is not what a passive radar is +limited by: the delay measures (RX->target) + (target->TX) - baseline, and +that sum is what sets received power via the bistatic radar equation. The +difference is not cosmetic — at a fixed RX distance the bistatic range varies +by more than the whole budget depending on which way the target lies relative +to the transmitter. + +Nodes without max_bistatic_range_km keep the monostatic rule so real hardware +(which only ever carries max_range_km) is unaffected. +""" + +import math + +from retina_simulation.world import ( + NodeConfig, + SimulatedAircraft, + SimulationWorld, + _haversine_km, +) + +_RX_LAT, _RX_LON = 34.85, -82.39 +_BASELINE_KM = 40.0 +_TX_LON = _RX_LON + _BASELINE_KM / (111.32 * math.cos(math.radians(_RX_LAT))) + + +def _node(max_bistatic_range_km): + # Beam aimed west, wide open, so only the range rule can reject. + return NodeConfig( + node_id="bistatic-node", + rx_lat=_RX_LAT, + rx_lon=_RX_LON, + rx_alt_ft=0.0, + tx_lat=_RX_LAT, + tx_lon=_TX_LON, + tx_alt_ft=0.0, + beam_azimuth_deg=270.0, + beam_width_deg=200.0, + max_range_km=50.0, + max_bistatic_range_km=max_bistatic_range_km, + ) + + +def _aircraft(bearing_deg, range_km): + br = math.radians(bearing_deg) + return SimulatedAircraft( + object_id="probe", + lat=_RX_LAT + math.degrees((range_km * math.cos(br)) / 6371.0), + lon=_RX_LON + math.degrees((range_km * math.sin(br)) / (6371.0 * math.cos(math.radians(_RX_LAT)))), + alt_km=0.0, + vel_east=0.0, + vel_north=0.0, + vel_up=0.0, + heading_deg=0.0, + speed_km_s=0.2, + ) + + +def _bistatic_km(ac): + r_rx = _haversine_km(_RX_LAT, _RX_LON, ac.lat, ac.lon) + r_tx = _haversine_km(_RX_LAT, _TX_LON, ac.lat, ac.lon) + return r_rx + r_tx - _BASELINE_KM + + +class TestBistaticRangeGate: + def test_rejects_on_tx_leg_where_monostatic_accepts(self): + """The whole point: 35 km from the RX, inside max_range_km=50, but + directly away from the TX so the bistatic range is 70 km.""" + world = SimulationWorld() + ac = _aircraft(270.0, 35.0) + assert _bistatic_km(ac) > 60.0 + assert world._aircraft_in_detection_cone(ac, _node(None)) is True + assert world._aircraft_in_detection_cone(ac, _node(60.0)) is False + + def test_accepts_inside_the_bistatic_budget(self): + world = SimulationWorld() + ac = _aircraft(270.0, 20.0) + assert _bistatic_km(ac) < 60.0 + assert world._aircraft_in_detection_cone(ac, _node(60.0)) is True + + def test_same_rx_range_different_tx_leg(self): + """Bistatic range must vary with bearing at a fixed RX distance — + otherwise the gate has silently stayed monostatic.""" + near_tx = _aircraft(330.0, 28.0) + away_tx = _aircraft(270.0, 28.0) + r_near = _haversine_km(_RX_LAT, _RX_LON, near_tx.lat, near_tx.lon) + r_away = _haversine_km(_RX_LAT, _RX_LON, away_tx.lat, away_tx.lon) + assert abs(r_near - r_away) < 0.1, "probes must be at equal RX range" + assert _bistatic_km(near_tx) < _bistatic_km(away_tx) + + def test_absent_key_keeps_monostatic_behaviour(self): + """Real hardware carries only max_range_km and must be untouched.""" + world = SimulationWorld() + mono = _node(None) + assert mono.max_bistatic_range_km is None + assert world._aircraft_in_detection_cone(_aircraft(270.0, 45.0), mono) is True + assert world._aircraft_in_detection_cone(_aircraft(270.0, 55.0), mono) is False + + def test_beam_still_applies_under_the_bistatic_rule(self): + """Range is not the only gate — an out-of-beam target stays rejected.""" + world = SimulationWorld() + node = _node(60.0) + node.beam_width_deg = 40.0 + assert world._aircraft_in_detection_cone(_aircraft(90.0, 10.0), node) is False + + +class TestEveryNodeDeclaresABistaticLimit: + """A circle on the receiver is never a bistatic node's true footprint. + + The ring, solo and dual paths all declared max_bistatic_range_km; the + generic region-node path did not, so those nodes alone kept gating and + rendering as circles — visible on staging as three synthetic nodes + reporting bistatic=None while every other node reported 60.0. + """ + + def test_metro_fleet_is_uniformly_bistatic(self): + from retina_simulation.generator import generate_fleet + + fleet = generate_fleet(n_nodes=15, metro="gvl", n_cluster=10, n_clusters=1, use_tower_api=False, seed=42) + missing = [n["node_id"] for n in fleet if n.get("max_bistatic_range_km") is None] + assert not missing, f"nodes still monostatic: {missing}" + + def test_the_limit_matches_the_declared_range(self): + from retina_simulation.generator import generate_fleet + + fleet = generate_fleet(n_nodes=15, metro="gvl", n_cluster=10, n_clusters=1, use_tower_api=False, seed=7) + for n in fleet: + assert n["max_bistatic_range_km"] == n["max_range_km"], n["node_id"] diff --git a/tests/test_dual_fraction.py b/tests/test_dual_fraction.py new file mode 100644 index 0000000..f01c8f6 --- /dev/null +++ b/tests/test_dual_fraction.py @@ -0,0 +1,299 @@ +"""dual_fraction: a slice of the ring/scatter budget additionally run as +dual-illuminator sites (see generator._generate_dual_sites), independent of +--layout dual (already all-dual there). + +Also covers the orchestrator's scene-change detection: n_nodes/dual_fraction +are baked into the fleet at container boot, so the only way to apply a +change is a full regeneration — _poll_simulation_config exits the process +via orchestrator.stop() when the backend-reported scene drifts from the +container's own stamped scene. +""" + +import asyncio +import json +import re + +import pytest + +from retina_simulation.generator import generate_fleet + +_DUAL_ID_RE = re.compile(r"-DUAL-\d{4}[ab]$") + + +def _fleet(**kw): + kw.setdefault("use_tower_api", False) + kw.setdefault("seed", 42) + return generate_fleet(**kw) + + +class TestDualFractionCarve: + def test_dual_ids_appear_in_rx_sharing_pairs_scatter_gvl(self): + fleet = _fleet( + n_nodes=30, + metro="gvl", + n_cluster=30, + n_clusters=1, + layout="scatter", + dual_fraction=0.4, + ) + dual_nodes = [n for n in fleet if _DUAL_ID_RE.search(n["node_id"])] + # n_dual_sites = round(30 * 0.4 / 2) = 6 sites -> up to 12 nodes, + # clamped to (n_cluster=30 after solo carve) + n_metro budget. + assert dual_nodes + assert len(dual_nodes) == 12 + + sites = {} + for n in dual_nodes: + sites.setdefault(n["node_id"][:-1], []).append(n) + assert len(sites) == 6 + for a, b in sites.values(): + # Same antenna, same mast — only the transmitter differs (mirrors + # test_dual_sites.py's rx-sharing assertion for layout="dual"). + assert (a["rx_lat"], a["rx_lon"], a["rx_alt_ft"]) == (b["rx_lat"], b["rx_lon"], b["rx_alt_ft"]) + assert (a["tx_lat"], a["tx_lon"]) != (b["tx_lat"], b["tx_lon"]) + + def test_dual_ids_appear_in_ring_layout_too(self): + # Spec: "append dual sites after the existing layout output (scatter + # AND ring branches)" — ring is the default layout. + fleet = _fleet( + n_nodes=30, + metro="gvl", + n_cluster=30, + n_clusters=1, + layout="ring", + dual_fraction=0.4, + ) + dual_nodes = [n for n in fleet if _DUAL_ID_RE.search(n["node_id"])] + assert dual_nodes + + def test_dual_fraction_ignored_for_layout_dual(self): + # layout="dual" is already all-dual; dual_fraction must be a no-op, + # not an additional carve on top of an already-fully-dual cluster. + without = _fleet( + n_nodes=16, + metro="gvl", + n_cluster=16, + n_clusters=1, + layout="dual", + ) + with_frac = _fleet( + n_nodes=16, + metro="gvl", + n_cluster=16, + n_clusters=1, + layout="dual", + dual_fraction=0.5, + ) + assert without == with_frac + + +class TestDeterminismRegression: + def test_dual_fraction_zero_reproduces_todays_scene(self): + """dual_fraction=0.0 must not perturb the RNG stream of pre-existing + nodes — sites are appended AFTER the existing layout output, so a + zero carve means zero extra RNG draws and an identical fleet.""" + baseline = _fleet(n_nodes=50) + explicit_zero = _fleet(n_nodes=50, dual_fraction=0.0) + assert explicit_zero == baseline + + def test_dual_fraction_zero_with_metro_layout_ring(self): + baseline = _fleet(n_nodes=30, metro="gvl", n_cluster=30, n_clusters=1, layout="ring") + explicit_zero = _fleet( + n_nodes=30, + metro="gvl", + n_cluster=30, + n_clusters=1, + layout="ring", + dual_fraction=0.0, + ) + assert explicit_zero == baseline + + def test_same_seed_same_dual_fraction_is_deterministic(self): + a = _fleet(n_nodes=30, metro="gvl", n_cluster=30, n_clusters=1, layout="scatter", dual_fraction=0.4) + b = _fleet(n_nodes=30, metro="gvl", n_cluster=30, n_clusters=1, layout="scatter", dual_fraction=0.4) + assert a == b + + +class TestDualFractionClampAndGuards: + def test_dual_fraction_one_clamps_to_available_budget(self): + # n_dual_nodes is clamped to n_cluster + n_metro — it can never + # exceed the whole layout's node budget even at dual_fraction=1.0. + fleet = _fleet( + n_nodes=30, + metro="gvl", + n_cluster=30, + n_clusters=1, + layout="scatter", + dual_fraction=1.0, + ) + dual_nodes = [n for n in fleet if _DUAL_ID_RE.search(n["node_id"])] + non_solo_non_dual = [n for n in fleet if not _DUAL_ID_RE.search(n["node_id"]) and "SOLO" not in n["node_id"]] + # Whole cluster budget (30, after the metro-scoped solo carve leaves + # it untouched) went to dual sites; nothing left for scatter nodes. + assert len(dual_nodes) <= 30 + assert non_solo_non_dual == [] + + def test_dual_fraction_without_metro_raises(self): + with pytest.raises(ValueError): + _fleet(n_nodes=30, layout="ring", dual_fraction=0.3) + + +class _StubWorld: + frac_anomalous = 0.0 + frac_drone = 0.0 + frac_dark = 0.0 + min_aircraft = 1 + max_aircraft = 1 + + +class _StubOrchestrator: + def __init__(self, max_range_km=0.0): + self._running = True + self.world = _StubWorld() + self.stop_calls = 0 + # The orchestrator itself holds the running value — max_range_km + # needs no scene stamp, unlike n_nodes/dual_fraction below. + self.max_range_km = max_range_km + + async def stop(self): + self.stop_calls += 1 + self._running = False + + +class _FakeResponse: + def __init__(self, payload): + self._body = json.dumps(payload).encode() + + def read(self): + return self._body + + def __enter__(self): + return self + + def __exit__(self, *exc): + return False + + +def _run_one_poll(monkeypatch, orchestrator, scene, cfg): + """Run _poll_simulation_config for exactly one fetch cycle. + + Patches urllib.request.urlopen (imported locally inside the polled + function) to return `cfg` once, then flips orchestrator._running off so + the `while orchestrator._running` loop exits on its next check — needed + for the "no diff" cases where the function itself never returns early. + """ + import urllib.request + + from retina_simulation.orchestrator import _poll_simulation_config + + def _fake_urlopen(url, context=None, timeout=None): + orchestrator._running = False + return _FakeResponse(cfg) + + monkeypatch.setattr(urllib.request, "urlopen", _fake_urlopen) + asyncio.run( + _poll_simulation_config( + orchestrator, + "http://validation.example", + interval_s=0.0, + scene=scene, + ) + ) + + +class TestScenePollDetection: + def test_stop_called_on_scene_diff(self, monkeypatch): + orch = _StubOrchestrator() + scene = {"n_nodes": 30, "dual_fraction": 0.0} + cfg = {"_updated_at": 1.0, "n_nodes": 32, "dual_fraction": 0.2} + _run_one_poll(monkeypatch, orch, scene, cfg) + assert orch.stop_calls == 1 + + def test_stop_called_on_dual_fraction_diff_alone(self, monkeypatch): + orch = _StubOrchestrator() + scene = {"n_nodes": 30, "dual_fraction": 0.0} + cfg = {"_updated_at": 1.0, "n_nodes": 30, "dual_fraction": 0.2} + _run_one_poll(monkeypatch, orch, scene, cfg) + assert orch.stop_calls == 1 + + def test_stop_not_called_when_scene_matches(self, monkeypatch): + orch = _StubOrchestrator() + scene = {"n_nodes": 30, "dual_fraction": 0.0} + cfg = {"_updated_at": 1.0, "n_nodes": 30, "dual_fraction": 0.0} + _run_one_poll(monkeypatch, orch, scene, cfg) + assert orch.stop_calls == 0 + + def test_stop_not_called_within_float_tolerance(self, monkeypatch): + orch = _StubOrchestrator() + scene = {"n_nodes": 30, "dual_fraction": 0.2} + cfg = {"_updated_at": 1.0, "n_nodes": 30, "dual_fraction": 0.2 + 1e-9} + _run_one_poll(monkeypatch, orch, scene, cfg) + assert orch.stop_calls == 0 + + def test_stop_not_called_when_scene_is_empty(self, monkeypatch): + # Absent scene stamp (stale volume, in-process generation fallback) + # -> no comparison, no restart, regardless of how different cfg is. + orch = _StubOrchestrator() + cfg = {"_updated_at": 1.0, "n_nodes": 99, "dual_fraction": 0.9} + _run_one_poll(monkeypatch, orch, None, cfg) + assert orch.stop_calls == 0 + + def test_stop_not_called_when_config_omits_scene_keys(self, monkeypatch): + # Backend never had a PUT for n_nodes/dual_fraction (only-if-set + # pattern) -> absent keys -> no comparison. + orch = _StubOrchestrator() + scene = {"n_nodes": 30, "dual_fraction": 0.0} + cfg = {"_updated_at": 1.0} + _run_one_poll(monkeypatch, orch, scene, cfg) + assert orch.stop_calls == 0 + + +class TestScenePollMaxRangeKm: + """max_range_km is a scene key too — applying it requires regenerating + every node's config, so it goes through the same restart path as + n_nodes/dual_fraction. Unlike those two it needs no scene stamp: the + orchestrator itself holds the running value, so the comparison runs even + when `scene` is None.""" + + def test_stop_called_on_max_range_km_diff(self, monkeypatch): + orch = _StubOrchestrator(max_range_km=0.0) + cfg = {"_updated_at": 1.0, "max_range_km": 150.0} + _run_one_poll(monkeypatch, orch, None, cfg) + assert orch.stop_calls == 1 + + def test_stop_called_on_max_range_km_diff_with_scene_present(self, monkeypatch): + # Runs independent of the scene stamp -- present-but-matching scene + # keys must not suppress the max_range_km comparison. + orch = _StubOrchestrator(max_range_km=0.0) + scene = {"n_nodes": 30, "dual_fraction": 0.0} + cfg = {"_updated_at": 1.0, "n_nodes": 30, "dual_fraction": 0.0, "max_range_km": 200.0} + _run_one_poll(monkeypatch, orch, scene, cfg) + assert orch.stop_calls == 1 + + def test_stop_not_called_when_max_range_km_matches(self, monkeypatch): + orch = _StubOrchestrator(max_range_km=100.0) + cfg = {"_updated_at": 1.0, "max_range_km": 100.0} + _run_one_poll(monkeypatch, orch, None, cfg) + assert orch.stop_calls == 0 + + def test_stop_not_called_within_float_tolerance(self, monkeypatch): + orch = _StubOrchestrator(max_range_km=100.0) + cfg = {"_updated_at": 1.0, "max_range_km": 100.0 + 1e-9} + _run_one_poll(monkeypatch, orch, None, cfg) + assert orch.stop_calls == 0 + + def test_stop_not_called_when_max_range_km_absent(self, monkeypatch): + # Never PUT (only-if-set pattern) -> absent key -> no comparison, + # even though the running value is nonzero. + orch = _StubOrchestrator(max_range_km=75.0) + cfg = {"_updated_at": 1.0} + _run_one_poll(monkeypatch, orch, None, cfg) + assert orch.stop_calls == 0 + + def test_works_with_scene_stamp_none_and_no_scene_diff(self, monkeypatch): + # scene=None must not itself raise or skip the max_range_km check -- + # only n_nodes/dual_fraction are guarded by the stamp. + orch = _StubOrchestrator(max_range_km=50.0) + cfg = {"_updated_at": 1.0, "max_range_km": 50.0} + _run_one_poll(monkeypatch, orch, None, cfg) + assert orch.stop_calls == 0 diff --git a/tests/test_dual_sites.py b/tests/test_dual_sites.py new file mode 100644 index 0000000..0bf2e91 --- /dev/null +++ b/tests/test_dual_sites.py @@ -0,0 +1,187 @@ +"""Dual-illuminator sites: one antenna, one RX, two transmitters. + +This is how a real passive-radar site is built, and the geometric payoff is +that the two bistatic ellipses share a focus (the common RX). They therefore +intersect in at most two points and the beam almost always excludes one, so a +single site localises on its own with the residual ambiguity bounded by the +antenna pattern rather than by a second receiver tens of km away. +""" + +import math + +from retina_simulation.generator import ( + _subtended_deg, + generate_fleet, +) + + +def _dual_sites(**kw): + kw.setdefault("n_nodes", 16) + kw.setdefault("metro", "gvl") + kw.setdefault("n_cluster", 16) + kw.setdefault("n_clusters", 1) + kw.setdefault("use_tower_api", False) + kw.setdefault("seed", 42) + kw.setdefault("layout", "dual") + fleet = generate_fleet(**kw) + sites = {} + for n in fleet: + if "DUAL" in n["node_id"]: + sites.setdefault(n["node_id"][:-1], []).append(n) + return fleet, sites + + +class TestDualSiteStructure: + def test_every_site_emits_exactly_two_nodes(self): + _fleet, sites = _dual_sites() + assert sites + assert all(len(v) == 2 for v in sites.values()) + + def test_the_pair_shares_one_receiver(self): + """Same antenna, same mast — only the transmitter differs.""" + _fleet, sites = _dual_sites() + for a, b in sites.values(): + assert (a["rx_lat"], a["rx_lon"], a["rx_alt_ft"]) == (b["rx_lat"], b["rx_lon"], b["rx_alt_ft"]) + + def test_the_pair_shares_one_detection_area(self): + _fleet, sites = _dual_sites() + for a, b in sites.values(): + assert a["beam_azimuth_deg"] == b["beam_azimuth_deg"] + assert a["beam_width_deg"] == b["beam_width_deg"] + assert a["max_bistatic_range_km"] == b["max_bistatic_range_km"] + + def test_the_pair_uses_two_different_transmitters(self): + """Two receivers on one mast would share an ellipse and be useless.""" + _fleet, sites = _dual_sites() + for a, b in sites.values(): + assert (a["tx_lat"], a["tx_lon"]) != (b["tx_lat"], b["tx_lon"]) + sep = math.hypot((a["tx_lat"] - b["tx_lat"]) * 111.32, (a["tx_lon"] - b["tx_lon"]) * 91.0) + assert sep > 1.0, "co-sited transmitters give an identical ellipse" + + def test_vhf_band_restriction_is_honoured(self): + _fleet, sites = _dual_sites(illuminator_band="vhf") + for pair in sites.values(): + for n in pair: + assert n["fc_hz"] < 300e6 + + def test_sites_are_scattered_not_ringed(self): + """A ring puts every receiver the same distance from the core; these + are drawn across the metro, which is what keeps inter-site overlap + low.""" + _fleet, sites = _dual_sites() + core = (34.852, -82.394) + d = [math.hypot((p[0]["rx_lat"] - core[0]) * 111.32, (p[0]["rx_lon"] - core[1]) * 91.0) for p in sites.values()] + assert max(d) - min(d) > 20.0 + + +class TestIlluminatorSelection: + def test_pairs_are_geometrically_usable(self): + """Selection maximises the worst subtended angle over the beam + footprint; a pair whose transmitters lie in nearly the same direction + gives two near-parallel ellipses and a degenerate intersection.""" + _fleet, sites = _dual_sites() + angles = [ + _subtended_deg(a["rx_lat"], a["rx_lon"], (a["tx_lat"], a["tx_lon"]), (b["tx_lat"], b["tx_lon"])) + for a, b in sites.values() + ] + assert sum(angles) / len(angles) > 25.0 + + def test_eirp_floor_excludes_weak_illuminators(self): + """Spartanburg (27 dBm) is geometrically valuable but radiologically + weak; a high floor must drop it.""" + _fleet, sites = _dual_sites(dual_min_eirp_dbm=75.0) + used = {n["tx_callsign"] for pair in sites.values() for n in pair} + assert "BLP01065" not in used + assert "W07DT-D" not in used + + +class TestAim: + def test_core_aim_points_beams_at_the_traffic(self): + """Random aiming starved the layout — with most aircraft routed through + the metro core, a randomly-pointed beam sees nothing and the solve rate + collapsed by an order of magnitude.""" + _fleet, sites = _dual_sites(dual_aim="core") + core = (34.852, -82.394) + off = [] + for pair in sites.values(): + n = pair[0] + want = ( + math.degrees( + math.atan2((core[1] - n["rx_lon"]) * math.cos(math.radians(core[0])), core[0] - n["rx_lat"]) + ) + % 360.0 + ) + off.append(abs((n["beam_azimuth_deg"] - want + 180) % 360 - 180)) + assert sum(off) / len(off) < 45.0 + + def test_random_aim_is_still_available(self): + _fleet, sites = _dual_sites(dual_aim="random") + assert sites + + +class TestGeneratorCLI: + """main() must be able to supply every argument generate_fleet takes. + + --dual-aim was added to generate_fleet and to the offline bench but never + registered on the generator's own parser, while main() passed + args.dual_aim regardless. Nothing caught it: every test and the bench call + generate_fleet() directly, and the CLI runs only in the fleet container's + entrypoint — so the break surfaced as a staging deploy coming up with zero + synthetic nodes. + """ + + def _parser_dests(self): + import argparse + from unittest import mock + + import retina_simulation.generator as gen + + captured = {} + real_parse = argparse.ArgumentParser.parse_args + + def _capture(self, *a, **kw): + captured["dests"] = {act.dest for act in self._actions} + raise SystemExit(0) # stop before generating a fleet + + with mock.patch.object(argparse.ArgumentParser, "parse_args", _capture): + try: + gen.main() + except SystemExit: + pass + assert real_parse is argparse.ArgumentParser.parse_args + return captured["dests"] + + def test_every_arg_main_reads_is_registered(self): + import ast + import inspect + + import retina_simulation.generator as gen + + src = inspect.getsource(gen.main) + tree = ast.parse(src.lstrip()) + used = { + node.attr + for node in ast.walk(tree) + if isinstance(node, ast.Attribute) and isinstance(node.value, ast.Name) and node.value.id == "args" + } + missing = used - self._parser_dests() + assert not missing, f"main() reads unregistered args: {sorted(missing)}" + + def test_dual_aim_is_registered_and_defaults_to_core(self): + import argparse + from unittest import mock + + import retina_simulation.generator as gen + + seen = {} + + def _capture(self, *a, **kw): + seen["ns"] = argparse.Namespace(**{act.dest: act.default for act in self._actions}) + raise SystemExit(0) + + with mock.patch.object(argparse.ArgumentParser, "parse_args", _capture): + try: + gen.main() + except SystemExit: + pass + assert seen["ns"].dual_aim == "core" diff --git a/tests/test_ground_truth_push.py b/tests/test_ground_truth_push.py new file mode 100644 index 0000000..66b1a2d --- /dev/null +++ b/tests/test_ground_truth_push.py @@ -0,0 +1,60 @@ +"""Ground-truth push payload schema (build_ground_truth_payload). + +The server's debug map shows simulated parameters per ground-truth object, so +the push payload must carry the ADS-B/anomaly attributes the world already +tracks — a silent field drop here is invisible until someone clicks a dot. +""" + +from retina_simulation.orchestrator import build_ground_truth_payload + + +def _summary(**overrides) -> dict: + base = { + "id": "obj-0001", + "lat": 34.85, + "lon": -82.4, + "alt_km": 9.5, + "heading": 270.0, + "speed_ms": 230.0, + "has_adsb": True, + "is_anomalous": False, + "object_type": "aircraft", + "adsb_hex": "a1b2c3", + "adsb_callsign": "ABC1234", + "anomaly_event": None, + } + base.update(overrides) + return base + + +class TestBuildGroundTruthPayload: + def test_passes_through_adsb_and_anomaly_fields(self): + out = build_ground_truth_payload([_summary(anomaly_event="hijack", is_anomalous=True)]) + assert len(out) == 1 + entry = out[0] + assert entry["hex"] == "a1b2c3" + assert entry["has_adsb"] is True + assert entry["adsb_callsign"] == "ABC1234" + assert entry["anomaly_event"] == "hijack" + assert entry["is_anomalous"] is True + assert entry["alt_m"] == 9500.0 + + def test_dark_object_falls_back_to_object_id(self): + out = build_ground_truth_payload([_summary(adsb_hex=None, adsb_callsign=None, has_adsb=False)]) + assert out[0]["hex"] == "obj-0001" + assert out[0]["has_adsb"] is False + assert out[0]["adsb_callsign"] is None + + def test_entry_without_hex_or_id_is_dropped(self): + out = build_ground_truth_payload([_summary(adsb_hex=None, id="")]) + assert out == [] + + def test_missing_optional_keys_default(self): + # A summary from an older world build without the new keys must not + # crash and must default sanely. + minimal = {"id": "obj-2", "lat": 1.0, "lon": 2.0, "alt_km": 3.0} + out = build_ground_truth_payload([minimal]) + assert out[0]["has_adsb"] is False + assert out[0]["adsb_callsign"] is None + assert out[0]["anomaly_event"] is None + assert out[0]["object_type"] == "aircraft" diff --git a/tests/test_retirement.py b/tests/test_retirement.py new file mode 100644 index 0000000..c9212ca --- /dev/null +++ b/tests/test_retirement.py @@ -0,0 +1,113 @@ +"""Aircraft retirement must happen at the region edge, never mid-scene. + +Staging measurement (2 min, 1 Hz feed sampling): every ground-truth vanish +was a lifetime expiry, and lifetimes uniform(180, 900) s expire aircraft +wherever they happen to be — including directly over the metro core, where a +dot blinking out reads as a tracking bug. Expiry marks the aircraft and +reroutes it toward the edge (_route_out); it keeps flying until it clears +retire_edge_km. The old 2x-lifetime hard cap fired regardless of position, +so metro-routed traffic — which rarely crosses the edge on its own — still +vanished mid-view; non-drones now get exit_grace_s of outbound flight and +only a genuinely stuck aircraft hits the backstop. Drones keep the 2x cap: +they loop low and slow and are expected to churn. +""" + +from retina_simulation.world import SimulationWorld, _haversine_km + + +def _world(): + w = SimulationWorld(center_lat=35.0, center_lon=-82.4) + w.min_aircraft = 0 + w.max_aircraft = 0 + return w + + +def _plant(w, lat, lon, created_at=0.0, lifetime_s=100.0): + ac = w._spawn_aircraft("adsb") + ac.lat, ac.lon = lat, lon + ac.created_at = created_at + ac.lifetime_s = lifetime_s + w.aircraft = [ac] + return ac + + +class TestEdgeGatedRetirement: + def test_unexpired_aircraft_is_kept(self): + w = _world() + _plant(w, 35.0, -82.4, lifetime_s=1000.0) + w._time = 100.0 + w.step(0.0, mode="adsb") + assert len(w.aircraft) == 1 + + def test_expired_but_central_aircraft_keeps_flying(self): + w = _world() + _plant(w, 35.05, -82.35) # a few km from center + w._time = 150.0 # expired (100 s lifetime), inside the exit grace + w.step(0.0, mode="adsb") + assert len(w.aircraft) == 1 + + def test_expired_aircraft_at_the_edge_is_retired(self): + w = _world() + _plant(w, 35.0, -81.5) # ~82 km east of center — beyond retire_edge_km + w._time = 150.0 + w.step(0.0, mode="adsb") + assert w.aircraft == [] + + def test_exit_grace_backstop_retires_a_stuck_aircraft_anywhere(self): + w = _world() + _plant(w, 35.0, -82.4) # dead center + w._time = 100.0 + w.exit_grace_s + 1.0 + w.step(0.0, mode="adsb") + assert w.aircraft == [] + + def test_expired_central_aircraft_survives_the_old_2x_cap(self): + w = _world() + _plant(w, 35.0, -82.4) + w._time = 201.0 # past the old 2x cap, well inside lifetime + grace + w.step(0.0, mode="adsb") + assert len(w.aircraft) == 1 + + +class TestFlyOutRetirement: + def test_expiry_reroutes_toward_the_edge(self): + w = _world() + ac = _plant(w, 35.05, -82.35) # northeast of center + ac.waypoints = [(35.0, -82.4), (35.1, -82.3)] # looping metro route + ac.waypoint_idx = 0 + w._time = 150.0 + w.step(0.0, mode="adsb") + + assert ac.departing is True + assert len(ac.waypoints) == 1 + wp_lat, wp_lon = ac.waypoints[0] + # The single exit waypoint sits beyond the retire edge... + assert _haversine_km(w.center_lat, w.center_lon, wp_lat, wp_lon) > w.retire_edge_km + # ...on the bearing away from center (northeast, like the aircraft). + assert wp_lat > w.center_lat and wp_lon > w.center_lon + + def test_reroute_fires_once(self): + w = _world() + ac = _plant(w, 35.05, -82.35) + w._time = 150.0 + w.step(0.0, mode="adsb") + first_route = list(ac.waypoints) + w.step(1.0, mode="adsb") + assert ac.waypoints == first_route # not re-planned every step + + def test_an_aircraft_at_the_exact_center_departs_on_its_heading(self): + w = _world() + ac = _plant(w, 35.0, -82.4) # bearing from center undefined + ac.heading_deg = 90.0 + w._time = 150.0 + w.step(0.0, mode="adsb") + wp_lat, wp_lon = ac.waypoints[0] + assert wp_lon > w.center_lon # due east, per the heading + assert abs(wp_lat - w.center_lat) < 0.05 + + def test_a_drone_keeps_the_2x_cap(self): + w = _world() + ac = _plant(w, 35.0, -82.4) + ac.object_type = "drone" + w._time = 201.0 # past 2x lifetime + w.step(0.0, mode="adsb") + assert w.aircraft == [] diff --git a/tests/test_scatter_sites.py b/tests/test_scatter_sites.py new file mode 100644 index 0000000..1828fb2 --- /dev/null +++ b/tests/test_scatter_sites.py @@ -0,0 +1,222 @@ +"""Scatter layout: the fleet nobody designed. + +Ring and dual both buy geometry on purpose — receivers placed to surround the +airspace, or paired on one mast to share a focus. A community deployment does +neither: sites land where operators live, on whatever illuminator comes in +best, pointed by hand. These tests pin the properties that make the layout a +fair test of the solver rather than a second ring under another name. +""" + +import math +from collections import Counter + +from retina_simulation.generator import ( + _KNOWN_METROS, + _haversine_km, + generate_fleet, +) + + +def _scatter(**kw): + kw.setdefault("n_nodes", 16) + kw.setdefault("metro", "gvl") + kw.setdefault("n_cluster", 12) + kw.setdefault("n_clusters", 1) + kw.setdefault("use_tower_api", False) + kw.setdefault("seed", 42) + kw.setdefault("layout", "scatter") + fleet = generate_fleet(**kw) + return fleet, [n for n in fleet if "SCAT" in n["node_id"]] + + +GVL = _KNOWN_METROS["gvl"] + + +class TestScatterStructure: + def test_the_layout_produces_the_requested_node_budget(self): + _fleet, scat = _scatter(n_cluster=12) + assert len(scat) == 12 + + def test_no_ring_nodes_survive(self): + """Scatter replaces the ring; a ring alongside it would supply the + designed geometry the layout exists to do without.""" + fleet, _scat = _scatter() + assert not [n for n in fleet if "RING" in n["node_id"]] + + def test_every_site_stays_inside_the_metro(self): + """A site outside the metro never sees an aircraft — the traffic model + is metro-scoped, so it would be a permanently dark node.""" + _fleet, scat = _scatter() + radius_km = GVL["radius_nm"] * 1.852 + for n in scat: + assert _haversine_km(n["rx_lat"], n["rx_lon"], GVL["lat"], GVL["lon"]) <= radius_km + + def test_sites_carry_an_explicit_hand_aimed_azimuth(self): + _fleet, scat = _scatter() + assert all("beam_azimuth_deg" in n for n in scat) + assert all(0.0 <= n["beam_azimuth_deg"] < 360.0 for n in scat) + + +class TestUndesignedGeometry: + def test_sites_do_not_share_one_illuminator(self): + """The ring's shared TX is the thing that keeps every node's Doppler + inside one association gate. A real fleet has no such luxury.""" + _fleet, scat = _scatter() + towers = Counter(n["tx_callsign"] for n in scat) + assert len(towers) >= 3 + # And no single tower carries the fleet the way the ring TX does. + assert towers.most_common(1)[0][1] < len(scat) * 0.75 + + def test_illuminators_are_nearby_ones(self): + """Weighting is 1/d^2 — an operator uses the station that comes in + best, not one 200 km away.""" + _fleet, scat = _scatter() + for n in scat: + baseline = _haversine_km(n["rx_lat"], n["rx_lon"], n["tx_lat"], n["tx_lon"]) + assert baseline <= 75.0 + + def test_aim_is_not_uniformly_at_the_core(self): + """The ring aims every beam at the core exactly. Here most point + roughly there and some point elsewhere, so beams do not all intersect.""" + _fleet, scat = _scatter() + errs = [] + for n in scat: + to_core = ( + math.degrees( + math.atan2( + math.sin(math.radians(GVL["lon"] - n["rx_lon"])) * math.cos(math.radians(GVL["lat"])), + math.cos(math.radians(n["rx_lat"])) * math.sin(math.radians(GVL["lat"])) + - math.sin(math.radians(n["rx_lat"])) + * math.cos(math.radians(GVL["lat"])) + * math.cos(math.radians(GVL["lon"] - n["rx_lon"])), + ) + ) + % 360.0 + ) + errs.append(abs((n["beam_azimuth_deg"] - to_core + 180.0) % 360.0 - 180.0)) + assert max(errs) > 45.0 # somebody points well off-core + assert sum(e < 45.0 for e in errs) >= len(errs) // 2 # most do not + + def test_reach_varies_but_beamwidth_is_uniform(self): + """Reach varies per site (60 km is what a good setup achieves, not an + average one), but every antenna is the same 42-degree Yagi — width + jitter was removed deliberately.""" + _fleet, scat = _scatter() + reaches = {n["max_bistatic_range_km"] for n in scat} + widths = {n["beam_width_deg"] for n in scat} + assert len(reaches) > 1 + assert widths == {42.0} + assert max(n["max_bistatic_range_km"] for n in scat) <= 60.0 + + def test_max_range_and_bistatic_limit_agree(self): + """Consumers that have not been taught the bistatic rule read + max_range_km; they must not see a different number.""" + _fleet, scat = _scatter() + for n in scat: + assert n["max_range_km"] == n["max_bistatic_range_km"] + + def test_sites_clump_rather_than_spacing_evenly(self): + """A ring's nearest-neighbour distances are all equal by construction. + Clumping is what puts near-parallel range gradients in the fleet, which + is the geometry the solver actually has to cope with.""" + _fleet, scat = _scatter() + nearest = [] + for a in scat: + nearest.append( + min(_haversine_km(a["rx_lat"], a["rx_lon"], b["rx_lat"], b["rx_lon"]) for b in scat if b is not a) + ) + assert max(nearest) / max(min(nearest), 0.1) > 3.0 + + +class TestDeterminism: + def test_same_seed_same_fleet(self): + a, _ = _scatter(seed=7) + b, _ = _scatter(seed=7) + assert a == b + + def test_different_seed_different_placement(self): + _fa, sa = _scatter(seed=7) + _fb, sb = _scatter(seed=8) + assert [(n["rx_lat"], n["rx_lon"]) for n in sa] != [(n["rx_lat"], n["rx_lon"]) for n in sb] + + +class TestLayoutFleetFaults: + """Stage-1 fixes: short/empty fleets and ring cells for ringless layouts.""" + + def test_scatter_without_metro_raises_instead_of_a_short_fleet(self): + import pytest + + with pytest.raises(ValueError, match="requires --metro"): + generate_fleet(n_nodes=16, n_cluster=12, layout="scatter", use_tower_api=False, seed=42) + + def test_dual_without_metro_raises_instead_of_a_short_fleet(self): + import pytest + + with pytest.raises(ValueError, match="requires --metro"): + generate_fleet(n_nodes=16, n_cluster=12, layout="dual", use_tower_api=False, seed=42) + + def test_scatter_budget_is_not_gated_on_the_ring_table(self): + """A metro with no _RING_TXS entry used to zero n_cluster for scatter, + which has no rings at all.""" + fleet, scat = _scatter(n_clusters=0) + assert len(scat) == 12 + assert len(fleet) == 16 + + def test_coverage_cells_do_not_describe_rings_for_ringless_layouts(self): + from retina_simulation.generator import coverage_cells + + cells = coverage_cells(n_cluster=12, n_clusters=1, metro="gvl", layout="scatter") + assert len(cells) == 1 + assert cells[0]["ring_id"] == "synth-SCATTER" + assert abs(cells[0]["core_lat"] - GVL["lat"]) < 1e-9 + # And without a metro there is nothing to describe. + assert coverage_cells(n_cluster=12, layout="scatter") == [] + + def test_ring_layout_cells_are_unchanged(self): + from retina_simulation.generator import coverage_cells + + ring = coverage_cells(n_cluster=12, n_clusters=1, metro="gvl", layout="ring") + default = coverage_cells(n_cluster=12, n_clusters=1, metro="gvl") + assert ring == default + + def test_solo_nodes_declare_the_bistatic_limit(self): + fleet = generate_fleet(n_nodes=30, n_cluster=8, n_clusters=1, use_tower_api=False, seed=42) + solos = [n for n in fleet if "SOLO" in n["node_id"]] + assert solos, "expected solo nodes in a nationwide fleet" + for n in solos: + assert n.get("max_bistatic_range_km") == n["max_range_km"] + + def test_empty_fleet_summary_does_not_crash(self): + from retina_simulation.generator import fleet_summary + + s = fleet_summary([]) + assert s["total_nodes"] == 0 + + +class TestVerticalRateDecays: + def test_vel_up_decays_toward_level_flight(self): + import random + + from retina_simulation.world import SimulationWorld + + random.seed(3) + w = SimulationWorld() + w.step(1.0) # spawn traffic + tagged = list(w.aircraft) + assert tagged + for ac in tagged: + ac.vel_up = 0.003 # force the worst-case spawn climb rate + ac.alt_km = 8.0 + # Outlive the 600 s window: the subject is vel_up decay, and + # whether seed 3's lifetime rolls happen to retire the tagged + # aircraft first is RNG-stream trivia (spawn-separation resampling + # legitimately consumes extra draws). + ac.lifetime_s = 10_000.0 + for _ in range(600): # 600 s of simulation + w.step(1.0) + survivors = [ac for ac in w.aircraft if ac in tagged] + assert survivors + # vel_up used to integrate forever; every aircraft saturated at the + # 15 km ceiling. With the decay nobody is pinned there. + assert all(abs(ac.vel_up) < 0.001 for ac in survivors) + assert all(ac.alt_km < 15.0 for ac in survivors) diff --git a/tests/test_separation.py b/tests/test_separation.py new file mode 100644 index 0000000..78d349d --- /dev/null +++ b/tests/test_separation.py @@ -0,0 +1,153 @@ +"""Traffic separation: spawn-pose resampling and in-trail speed modulation. + +The hub-radial planner converges metro spawns on the same ~2 km core, so +without separation the fleet routinely flew pairs inside the solver's +association gates. Spawn poses now resample away from live traffic +(best-effort), and in-flight conflicts slow the later-created aircraft +toward 70% of cruise until clear. Anomalous aircraft are exempt — erratic +close approaches are the anomaly signature, not a bug. +""" + +import random + +from retina_simulation.world import ( + MetroCell, + NodeConfig, + SimulatedAircraft, + SimulationWorld, + _haversine_km, +) + +_CORE_LAT, _CORE_LON = 34.85, -82.39 + + +def _world(): + world = SimulationWorld() + world.add_node(NodeConfig(node_id="anchor")) + world.metro_cells = [MetroCell(core_lat=_CORE_LAT, core_lon=_CORE_LON, radius_km=70.0)] + world.frac_metro_traffic = 1.0 + return world + + +def _plane(oid, lat, lon, alt_km=8.0, speed=0.2, created_at=0.0, **kw): + return SimulatedAircraft( + object_id=oid, + lat=lat, + lon=lon, + alt_km=alt_km, + vel_east=0.0, + vel_north=speed, + vel_up=0.0, + heading_deg=0.0, + speed_km_s=speed, + base_speed_km_s=speed, + created_at=created_at, + waypoints=[(lat, lon), (lat + 2.0, lon)], + waypoint_idx=1, + **kw, + ) + + +class TestNearestTraffic: + def test_empty_world_is_unconstrained(self): + assert _world()._nearest_traffic_km(_CORE_LAT, _CORE_LON) == float("inf") + + def test_reports_closest_of_several(self): + world = _world() + world.aircraft = [ + _plane("a", _CORE_LAT + 0.5, _CORE_LON), + _plane("b", _CORE_LAT + 0.1, _CORE_LON), + ] + d = world._nearest_traffic_km(_CORE_LAT, _CORE_LON) + assert abs(d - _haversine_km(_CORE_LAT, _CORE_LON, _CORE_LAT + 0.1, _CORE_LON)) < 1e-6 + + +class TestSpawnSeparation: + def test_spawns_avoid_live_traffic(self): + """With a 70 km cell and one blocker parked on the core, resampled + spawns should clear min_separation_km essentially always — the pose + space is enormous relative to one 5 km bubble.""" + random.seed(42) + world = _world() + world.aircraft = [_plane("blocker", _CORE_LAT, _CORE_LON)] + clear = 0 + for _ in range(30): + ac = world._spawn_aircraft(mode="adsb") + if _haversine_km(ac.lat, ac.lon, _CORE_LAT, _CORE_LON) >= world.min_separation_km: + clear += 1 + assert clear >= 28 # best-effort: allow the odd saturated roll + + def test_best_effort_never_blocks_spawn(self): + """A blocker on every candidate pose still yields an aircraft — + separation degrades, spawning never deadlocks (step() spawns in a + while-loop up to min_aircraft).""" + random.seed(7) + world = _world() + world.metro_cells[0].radius_km = 1.0 # pose space smaller than the bubble + world.aircraft = [_plane("blocker", _CORE_LAT, _CORE_LON)] + ac = world._spawn_aircraft(mode="adsb") + assert ac is not None + + def test_spawned_aircraft_carry_cruise_speed(self): + random.seed(3) + ac = _world()._spawn_aircraft(mode="adsb") + assert ac.base_speed_km_s == ac.speed_km_s > 0 + + +class TestInTrailModulation: + def test_trailing_conflict_slows_toward_seventy_percent(self): + world = _world() + lead = _plane("lead", _CORE_LAT, _CORE_LON, created_at=0.0) + trail = _plane("trail", _CORE_LAT + 0.01, _CORE_LON, created_at=10.0) + world.aircraft = [lead, trail] + for _ in range(30): + world._enforce_separation(1.0) + assert trail.speed_km_s < 0.72 * trail.base_speed_km_s + assert lead.speed_km_s > 0.95 * lead.base_speed_km_s + + def test_vertical_separation_is_no_conflict(self): + world = _world() + lead = _plane("lead", _CORE_LAT, _CORE_LON, alt_km=6.0, created_at=0.0) + trail = _plane("trail", _CORE_LAT + 0.01, _CORE_LON, alt_km=9.0, created_at=10.0) + world.aircraft = [lead, trail] + world._enforce_separation(1.0) + assert trail.speed_km_s == trail.base_speed_km_s + + def test_recovers_to_cruise_when_clear(self): + world = _world() + trail = _plane("trail", _CORE_LAT, _CORE_LON, created_at=10.0) + trail.speed_km_s = 0.7 * trail.base_speed_km_s + world.aircraft = [trail] + for _ in range(30): + world._enforce_separation(1.0) + assert trail.speed_km_s > 0.99 * trail.base_speed_km_s + + def test_anomalous_and_drones_exempt(self): + world = _world() + lead = _plane("lead", _CORE_LAT, _CORE_LON, created_at=0.0) + anom = _plane("anom", _CORE_LAT + 0.005, _CORE_LON, created_at=5.0, is_anomalous=True, object_type="anomalous") + drone = _plane("drone", _CORE_LAT + 0.005, _CORE_LON, created_at=6.0, object_type="drone") + world.aircraft = [lead, anom, drone] + world._enforce_separation(1.0) + assert anom.speed_km_s == anom.base_speed_km_s + assert drone.speed_km_s == drone.base_speed_km_s + + def test_step_reduces_close_pairs_over_time(self): + """End-to-end: a running world holds materially fewer sub-5 km + same-level pairs than the pre-separation planner produced.""" + random.seed(1234) + world = _world() + world.min_aircraft, world.max_aircraft = 20, 25 + for _ in range(600): + world.step(1.0, mode="adsb") + flow = [ac for ac in world.aircraft if not ac.is_anomalous and ac.object_type == "aircraft"] + conflicts = sum( + 1 + for i, a in enumerate(flow) + for b in flow[i + 1 :] + if abs(a.alt_km - b.alt_km) < world.min_vertical_sep_km + and _haversine_km(a.lat, a.lon, b.lat, b.lon) < world.min_separation_km + ) + # Not zero — crossing flows transit each other's bubbles while the + # modulation sequences them — but stacking is gone. + assert conflicts <= 2