mirror of
https://github.com/ruvnet/RuView
synced 2026-08-11 20:41:44 +00:00
chore: organise examples/research-sota/ into 9 thematic folders with READMEs (#744)
User request: organise examples/research-sota/ into folders with READMEs and main overview. Moved 46 files into 9 thematic folders by thread family + research category: 01-physics-floor/ (R1, R6, R6.1) — bedrock primitives 02-placement/ (R6.2 family, 7 sub-ticks) — antenna placement 03-spatial-intelligence/ (R5, R7) — saliency + mincut 04-rssi/ (R8, R9) — RSSI-only sensing 05-cross-room-reid/ (R3 arc, 3 ticks) — cross-room identity 06-structure-detection/ (R12 arc, 3 ticks) — PABS + closed loop 07-negative-results/ (R13) — productive failure 08-verticals/ (R10, R11) — wildlife + maritime physics 09-quantum-fusion/ (R20.1) — ADR-114 quantum-classical demo Each folder has its own README.md documenting: - Scripts + headlines table - Why this folder bounds / composes with others - Sample output / honest scope - Cross-references to related loop notes + ADRs Main README.md at the top covers: - Folder map with thread numbers - Cross-folder dependency graph - Headline findings table (8 entries) - Reading order for newcomers (4 scripts in suggested order) - Honest scope (synthetic-physics caveats) All git mv operations preserve file history. Total: 46 files moved, 10 new READMEs (main + 9 sub) totalling ~1300 lines of organising documentation.
This commit is contained in:
rUv
@@ -0,0 +1,208 @@
|
||||
#!/usr/bin/env python3
|
||||
"""R12.1 — Pose-PABS closed loop.
|
||||
|
||||
See docs/research/sota-2026-05-22/R12_1-pose-pabs-closed-loop.md.
|
||||
|
||||
R12 PABS (tick 19) had a false-alarm problem: subject moving 10 cm gave
|
||||
PABS = 22,000x natural drift floor. R12 PABS noted: 'Real production
|
||||
PABS needs a pose-aware forward model updating from pose_tracker.rs in
|
||||
real-time. The actual structure-detection signal is PABS-after-pose-
|
||||
update.'
|
||||
|
||||
This tick implements the closed loop in synthetic form:
|
||||
1. Subject moves on a continuous trajectory
|
||||
2. 'Pose tracker' estimates the subject position (with noise)
|
||||
3. Forward model uses the ESTIMATED position to predict expected CSI
|
||||
4. PABS = |observed - expected| using the pose-updated expected
|
||||
5. At tick T_intrude, insert an unexpected second subject
|
||||
6. Measure: does PABS-after-pose-update spike at T_intrude vs being
|
||||
noisy during subject motion?
|
||||
|
||||
Pure NumPy.
|
||||
"""
|
||||
|
||||
from __future__ import annotations
|
||||
|
||||
import argparse
|
||||
import json
|
||||
from pathlib import Path
|
||||
import numpy as np
|
||||
|
||||
C = 2.998e8
|
||||
|
||||
|
||||
def wavelength_m(freq_ghz: float) -> float:
|
||||
return C / (freq_ghz * 1e9)
|
||||
|
||||
|
||||
def csi_contribution(pos, refl, tx, rx, sub_freqs_hz):
|
||||
d_tx = np.linalg.norm(pos - tx)
|
||||
d_rx = np.linalg.norm(pos - rx)
|
||||
d_direct = np.linalg.norm(tx - rx)
|
||||
delta_l = d_tx + d_rx - d_direct
|
||||
amp = refl / max(d_tx * d_rx, 1e-3)
|
||||
phase = 2 * np.pi * sub_freqs_hz * delta_l / C
|
||||
return amp * np.exp(1j * phase)
|
||||
|
||||
|
||||
def simulate(scatterers, tx, rx, freq_ghz, n_sub=52, sub_spacing_khz=312.5):
|
||||
sub_offsets = (np.arange(n_sub) - n_sub // 2) * sub_spacing_khz * 1e3
|
||||
sub_freqs = freq_ghz * 1e9 + sub_offsets
|
||||
total = np.zeros(n_sub, dtype=complex)
|
||||
for s in scatterers:
|
||||
total += csi_contribution(np.asarray(s["pos"]), s["refl"],
|
||||
np.asarray(tx), np.asarray(rx), sub_freqs)
|
||||
return total
|
||||
|
||||
|
||||
def human_body(cx, cy):
|
||||
return [
|
||||
{"pos": [cx, cy ], "refl": 0.10}, # head
|
||||
{"pos": [cx, cy ], "refl": 0.50}, # chest
|
||||
{"pos": [cx - 0.20, cy ], "refl": 0.10}, # arms
|
||||
{"pos": [cx + 0.20, cy ], "refl": 0.10},
|
||||
{"pos": [cx - 0.10, cy - 0.40], "refl": 0.10}, # legs
|
||||
{"pos": [cx + 0.10, cy - 0.40], "refl": 0.10},
|
||||
]
|
||||
|
||||
|
||||
def walls():
|
||||
return [
|
||||
{"pos": [0.5, 4.5], "refl": 0.30},
|
||||
{"pos": [4.5, 4.5], "refl": 0.25},
|
||||
{"pos": [0.5, 0.5], "refl": 0.20},
|
||||
{"pos": [4.5, 0.5], "refl": 0.15},
|
||||
]
|
||||
|
||||
|
||||
def pabs(observed, predicted):
|
||||
res = observed - predicted
|
||||
e_obs = np.linalg.norm(observed) ** 2
|
||||
return float(np.linalg.norm(res) ** 2 / max(e_obs, 1e-12))
|
||||
|
||||
|
||||
def pose_tracker_estimate(true_pos, std_noise=0.05, rng=None):
|
||||
"""Simulate a pose tracker with ~5 cm position noise.
|
||||
Real pose_tracker.rs achieves this at ~95% PCK@20."""
|
||||
rng = rng or np.random.default_rng(0)
|
||||
return true_pos + rng.standard_normal(2) * std_noise
|
||||
|
||||
|
||||
def main():
|
||||
parser = argparse.ArgumentParser()
|
||||
parser.add_argument("--out", default="examples/research-sota/r12_1_pose_pabs_results.json")
|
||||
args = parser.parse_args()
|
||||
|
||||
tx = np.array([0.0, 2.5])
|
||||
rx = np.array([5.0, 2.5])
|
||||
freq = 2.4
|
||||
rng = np.random.default_rng(7)
|
||||
|
||||
# Subject walks from (2.0, 2.0) to (3.0, 3.5) over 50 frames
|
||||
n_frames = 50
|
||||
trajectory = np.linspace([2.0, 2.0], [3.0, 3.5], n_frames)
|
||||
walls_static = walls()
|
||||
|
||||
# Intruder enters at frame T_intrude
|
||||
T_intrude = 25
|
||||
intruder_pos = (1.5, 1.5)
|
||||
|
||||
# Two PABS pipelines:
|
||||
# (a) FIXED expected scene (R12 PABS naive — expects subject at start position)
|
||||
# (b) POSE-UPDATED expected scene (R12.1 — uses pose-tracker estimate)
|
||||
fixed_subject_pos = trajectory[0] # never updated
|
||||
fixed_expected = human_body(*fixed_subject_pos) + walls_static
|
||||
y_fixed = simulate(fixed_expected, tx, rx, freq)
|
||||
|
||||
pabs_fixed = []
|
||||
pabs_pose_updated = []
|
||||
pose_estimates = []
|
||||
|
||||
for t in range(n_frames):
|
||||
true_pos = trajectory[t]
|
||||
# Build the observed scene
|
||||
scene_obs = human_body(*true_pos) + walls_static
|
||||
if t >= T_intrude:
|
||||
scene_obs = scene_obs + human_body(*intruder_pos)
|
||||
y_obs = simulate(scene_obs, tx, rx, freq)
|
||||
|
||||
# (a) Fixed expected
|
||||
pabs_fixed.append(pabs(y_obs, y_fixed))
|
||||
|
||||
# (b) Pose-updated expected
|
||||
est_pos = pose_tracker_estimate(true_pos, std_noise=0.05, rng=rng)
|
||||
pose_estimates.append(est_pos.tolist())
|
||||
expected_pose = human_body(*est_pos) + walls_static
|
||||
y_pose = simulate(expected_pose, tx, rx, freq)
|
||||
pabs_pose_updated.append(pabs(y_obs, y_pose))
|
||||
|
||||
pabs_fixed = np.array(pabs_fixed)
|
||||
pabs_pose_updated = np.array(pabs_pose_updated)
|
||||
|
||||
# Analysis:
|
||||
# During T<T_intrude: pose-updated should be LOW (pose tracker explains subject)
|
||||
# During T>=T_intrude: pose-updated should SPIKE (intruder unexplained)
|
||||
# Fixed should be HIGH throughout (subject motion always unexplained)
|
||||
|
||||
pre_intrude_fixed_mean = pabs_fixed[:T_intrude].mean()
|
||||
post_intrude_fixed_mean = pabs_fixed[T_intrude:].mean()
|
||||
pre_intrude_pose_mean = pabs_pose_updated[:T_intrude].mean()
|
||||
post_intrude_pose_mean = pabs_pose_updated[T_intrude:].mean()
|
||||
|
||||
pose_intruder_lift = post_intrude_pose_mean / max(pre_intrude_pose_mean, 1e-9)
|
||||
fixed_intruder_lift = post_intrude_fixed_mean / max(pre_intrude_fixed_mean, 1e-9)
|
||||
|
||||
out = {
|
||||
"config": {
|
||||
"n_frames": n_frames,
|
||||
"trajectory_start": trajectory[0].tolist(),
|
||||
"trajectory_end": trajectory[-1].tolist(),
|
||||
"T_intrude": T_intrude,
|
||||
"intruder_pos": list(intruder_pos),
|
||||
"pose_tracker_std_m": 0.05,
|
||||
},
|
||||
"pabs_fixed": pabs_fixed.tolist(),
|
||||
"pabs_pose_updated": pabs_pose_updated.tolist(),
|
||||
"pre_intrude_means": {
|
||||
"fixed": float(pre_intrude_fixed_mean),
|
||||
"pose": float(pre_intrude_pose_mean),
|
||||
},
|
||||
"post_intrude_means": {
|
||||
"fixed": float(post_intrude_fixed_mean),
|
||||
"pose": float(post_intrude_pose_mean),
|
||||
},
|
||||
"intruder_detection_lift": {
|
||||
"fixed": fixed_intruder_lift,
|
||||
"pose": pose_intruder_lift,
|
||||
},
|
||||
}
|
||||
Path(args.out).parent.mkdir(parents=True, exist_ok=True)
|
||||
Path(args.out).write_text(json.dumps(out, indent=2))
|
||||
|
||||
print("=== R12.1 pose-PABS closed loop ===")
|
||||
print(f" Subject walks {n_frames} frames from {trajectory[0]} to {trajectory[-1]}")
|
||||
print(f" Intruder enters at frame {T_intrude} at position {intruder_pos}")
|
||||
print(f" Pose tracker noise: 5 cm std (ADR-079 ~95% PCK@20 quality)")
|
||||
print()
|
||||
print(f"=== Mean PABS by phase ===")
|
||||
print(f" Phase Fixed-expected Pose-updated")
|
||||
print(f" Pre-intruder (T<25): {pre_intrude_fixed_mean:>14.4f} {pre_intrude_pose_mean:>13.4f}")
|
||||
print(f" Post-intruder (T>=25): {post_intrude_fixed_mean:>14.4f} {post_intrude_pose_mean:>13.4f}")
|
||||
print()
|
||||
print(f"=== Intruder detection lift ===")
|
||||
print(f" FIXED-expected pipeline: {fixed_intruder_lift:>7.2f}x (R12 naive)")
|
||||
print(f" POSE-UPDATED pipeline: {pose_intruder_lift:>7.2f}x (R12.1 closed loop)")
|
||||
print()
|
||||
if pose_intruder_lift > fixed_intruder_lift * 3:
|
||||
verdict = "CLOSED LOOP WORKS: pose-PABS lift > 3x the naive baseline. False-alarm problem from R12 PABS resolved."
|
||||
elif pose_intruder_lift > 2.0:
|
||||
verdict = "CLOSED LOOP WORKS: pose-PABS lift > 2x baseline. Intruder detection clean."
|
||||
else:
|
||||
verdict = "MARGINAL: pose-PABS lift not decisive vs baseline. May need temporal averaging."
|
||||
print(f"VERDICT: {verdict}")
|
||||
print()
|
||||
print(f"Wrote {args.out}")
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
main()
|
||||
@@ -0,0 +1,135 @@
|
||||
{
|
||||
"config": {
|
||||
"n_frames": 50,
|
||||
"trajectory_start": [
|
||||
2.0,
|
||||
2.0
|
||||
],
|
||||
"trajectory_end": [
|
||||
3.0,
|
||||
3.5
|
||||
],
|
||||
"T_intrude": 25,
|
||||
"intruder_pos": [
|
||||
1.5,
|
||||
1.5
|
||||
],
|
||||
"pose_tracker_std_m": 0.05
|
||||
},
|
||||
"pabs_fixed": [
|
||||
0.0,
|
||||
0.23976021993699137,
|
||||
1.333289923835776,
|
||||
4.7449972298645005,
|
||||
16.132302954344752,
|
||||
57.31864185847987,
|
||||
34.59671192160786,
|
||||
11.19613115945127,
|
||||
5.077096413694479,
|
||||
2.8125145174844848,
|
||||
1.7357497400150317,
|
||||
1.1422331113156927,
|
||||
0.7902984026449109,
|
||||
0.5844695055883886,
|
||||
0.48864852071817233,
|
||||
0.49495019610807023,
|
||||
0.5992183799548572,
|
||||
0.7707784100562064,
|
||||
0.9509356710764513,
|
||||
1.1010310944881865,
|
||||
1.2286767924050106,
|
||||
1.3666209606880533,
|
||||
1.5555622650632148,
|
||||
1.8511220775066175,
|
||||
2.3569113678968043,
|
||||
23.64420568922056,
|
||||
24.766708919894374,
|
||||
12.440097343342567,
|
||||
5.835505088452743,
|
||||
3.016239220001779,
|
||||
1.6368370866065183,
|
||||
0.8521752953170693,
|
||||
0.35830915433305105,
|
||||
0.06898386583751527,
|
||||
0.11286933302231912,
|
||||
1.49823836553597,
|
||||
11.73405853896596,
|
||||
15.012383585890914,
|
||||
5.44051226107576,
|
||||
2.450306678228625,
|
||||
1.144765319492743,
|
||||
0.43860379597713645,
|
||||
0.6217089528021075,
|
||||
40.28090119216048,
|
||||
9.742961313951346,
|
||||
2.4076884969330483,
|
||||
0.8288916761760434,
|
||||
0.12070720537158618,
|
||||
0.66996511955866,
|
||||
28.778255288508806
|
||||
],
|
||||
"pabs_pose_updated": [
|
||||
0.0397808142334705,
|
||||
0.5104513448136311,
|
||||
0.8158108392380339,
|
||||
0.9465194410415606,
|
||||
0.5508926517254545,
|
||||
0.6594979498306511,
|
||||
2.0582347819010445,
|
||||
0.6060528733141695,
|
||||
0.12736172431501477,
|
||||
0.5159119899356763,
|
||||
0.01556708655354054,
|
||||
0.007342537186192009,
|
||||
0.002804857511672747,
|
||||
0.020407791283141442,
|
||||
0.00023421796933611544,
|
||||
0.004093746595234462,
|
||||
0.008881014219198688,
|
||||
0.012739000996667617,
|
||||
0.028360834638721005,
|
||||
0.0004098514050666686,
|
||||
0.00010859128727197401,
|
||||
0.00016902339492389355,
|
||||
0.054732157887574226,
|
||||
0.0006514193522454603,
|
||||
0.6018761650863446,
|
||||
10.405708813283992,
|
||||
1.6307427510614485,
|
||||
0.7535171230661254,
|
||||
0.6341883054891835,
|
||||
1.1494872301598305,
|
||||
0.4973417823824021,
|
||||
0.5908828843636849,
|
||||
0.19423577429400954,
|
||||
1.4642997355851366,
|
||||
0.08691356242442586,
|
||||
1.3298358192934818,
|
||||
3.4730881799534568,
|
||||
0.11532793333150544,
|
||||
1.7292922842852005,
|
||||
2.527226823962975,
|
||||
0.26166589945633334,
|
||||
0.27967362635220994,
|
||||
0.13730251197140705,
|
||||
22.685535567483463,
|
||||
0.8599415629887098,
|
||||
1.0779487716387626,
|
||||
1.9983295809816795,
|
||||
1.2202290817498453,
|
||||
1.0205655174952935,
|
||||
14.910181149340993
|
||||
],
|
||||
"pre_intrude_means": {
|
||||
"fixed": 6.018746107769026,
|
||||
"pose": 0.30355570822863354
|
||||
},
|
||||
"post_intrude_means": {
|
||||
"fixed": 7.756075151466307,
|
||||
"pose": 2.841338490895822
|
||||
},
|
||||
"intruder_detection_lift": {
|
||||
"fixed": 1.2886529872816415,
|
||||
"pose": 9.360187978266477
|
||||
}
|
||||
}
|
||||
@@ -0,0 +1,203 @@
|
||||
#!/usr/bin/env python3
|
||||
"""R12 PABS — Physics-Anchored Background Subtraction structure detection.
|
||||
|
||||
See docs/research/sota-2026-05-22/R12-pabs-implementation.md.
|
||||
|
||||
R12 NEGATIVE concluded that naive SVD-spectrum-cosine-distance failed
|
||||
because the eigenshift was indistinguishable from natural drift. The
|
||||
deferred revision: 'PABS over Fresnel basis'. R6.1 just shipped the
|
||||
multi-scatterer Fresnel forward operator, so PABS is now implementable.
|
||||
|
||||
PABS = norm(y_observed - y_predicted)
|
||||
where y_predicted is computed from R6.1's multi-scatterer model
|
||||
using a population-prior body assumption.
|
||||
|
||||
Scenarios tested:
|
||||
A. Empty room (no occupant) — baseline PABS
|
||||
B. Subject standing (expected) — small PABS (expected occupant)
|
||||
C. Subject + added furniture (1 new piece) — large PABS (new structure)
|
||||
D. Subject + 2nd subject (unexpected person) — large PABS
|
||||
E. Subject + wall reflector moved (drift) — comparison vs natural drift
|
||||
|
||||
This is the experiment R12 wanted but couldn't run without R6.1. Pure NumPy.
|
||||
"""
|
||||
|
||||
from __future__ import annotations
|
||||
|
||||
import argparse
|
||||
import json
|
||||
from pathlib import Path
|
||||
import numpy as np
|
||||
|
||||
C = 2.998e8
|
||||
|
||||
|
||||
def wavelength_m(freq_ghz: float) -> float:
|
||||
return C / (freq_ghz * 1e9)
|
||||
|
||||
|
||||
def path_delta_m(scatterer_pos, tx_pos, rx_pos):
|
||||
d_tx = np.linalg.norm(scatterer_pos - tx_pos)
|
||||
d_rx = np.linalg.norm(scatterer_pos - rx_pos)
|
||||
d_direct = np.linalg.norm(tx_pos - rx_pos)
|
||||
return d_tx + d_rx - d_direct
|
||||
|
||||
|
||||
def csi_contribution(scatterer_pos, reflectivity, tx_pos, rx_pos, sub_freqs_hz):
|
||||
delta_l = path_delta_m(scatterer_pos, tx_pos, rx_pos)
|
||||
d_tx = np.linalg.norm(scatterer_pos - tx_pos)
|
||||
d_rx = np.linalg.norm(scatterer_pos - rx_pos)
|
||||
amp = reflectivity / max(d_tx * d_rx, 1e-3)
|
||||
phase = 2 * np.pi * sub_freqs_hz * delta_l / C
|
||||
return amp * np.exp(1j * phase)
|
||||
|
||||
|
||||
def simulate(scatterers, tx_pos, rx_pos, freq_ghz, n_sub=52, sub_spacing_khz=312.5):
|
||||
sub_offsets = (np.arange(n_sub) - n_sub // 2) * sub_spacing_khz * 1e3
|
||||
sub_freqs = freq_ghz * 1e9 + sub_offsets
|
||||
total = np.zeros(n_sub, dtype=complex)
|
||||
for s in scatterers:
|
||||
total += csi_contribution(np.asarray(s["pos"]), s["refl"],
|
||||
np.asarray(tx_pos), np.asarray(rx_pos), sub_freqs)
|
||||
return total
|
||||
|
||||
|
||||
def human_body(center_x, center_y):
|
||||
return [
|
||||
{"pos": [center_x, center_y ], "refl": 0.10, "name": "head"},
|
||||
{"pos": [center_x, center_y ], "refl": 0.50, "name": "chest"},
|
||||
{"pos": [center_x - 0.20, center_y ], "refl": 0.10, "name": "left_arm"},
|
||||
{"pos": [center_x + 0.20, center_y ], "refl": 0.10, "name": "right_arm"},
|
||||
{"pos": [center_x - 0.10, center_y - 0.40], "refl": 0.10, "name": "left_leg"},
|
||||
{"pos": [center_x + 0.10, center_y - 0.40], "refl": 0.10, "name": "right_leg"},
|
||||
]
|
||||
|
||||
|
||||
def static_wall_reflectors(amplitudes=(0.3, 0.2, 0.15, 0.1)):
|
||||
"""Four wall reflectors at fixed positions -- typical bedroom multipath."""
|
||||
return [
|
||||
{"pos": [0.5, 4.5], "refl": amplitudes[0], "name": "wall_NW"},
|
||||
{"pos": [4.5, 4.5], "refl": amplitudes[1], "name": "wall_NE"},
|
||||
{"pos": [0.5, 0.5], "refl": amplitudes[2], "name": "wall_SW"},
|
||||
{"pos": [4.5, 0.5], "refl": amplitudes[3], "name": "wall_SE"},
|
||||
]
|
||||
|
||||
|
||||
def pabs(y_observed, y_predicted):
|
||||
"""L2 norm of the residual, normalised by signal energy."""
|
||||
residual = y_observed - y_predicted
|
||||
energy = np.linalg.norm(y_observed) ** 2
|
||||
return float(np.linalg.norm(residual) ** 2 / max(energy, 1e-12))
|
||||
|
||||
|
||||
def main():
|
||||
parser = argparse.ArgumentParser()
|
||||
parser.add_argument("--out", default="examples/research-sota/r12_pabs_results.json")
|
||||
args = parser.parse_args()
|
||||
|
||||
tx = np.array([0.0, 2.5])
|
||||
rx = np.array([5.0, 2.5])
|
||||
freq_ghz = 2.4
|
||||
walls = static_wall_reflectors()
|
||||
|
||||
# ===== Build the "expected" scene model (subject + walls) =====
|
||||
# This is what PABS predicts as the baseline.
|
||||
subject_expected = human_body(2.5, 2.75)
|
||||
expected_scene = subject_expected + walls
|
||||
y_expected = simulate(expected_scene, tx, rx, freq_ghz)
|
||||
|
||||
# ===== Scenario A: empty room (no occupant) =====
|
||||
y_empty = simulate(walls, tx, rx, freq_ghz)
|
||||
pabs_A = pabs(y_empty, y_expected)
|
||||
|
||||
# ===== Scenario B: subject standing where expected =====
|
||||
y_B = simulate(subject_expected + walls, tx, rx, freq_ghz)
|
||||
pabs_B = pabs(y_B, y_expected)
|
||||
|
||||
# ===== Scenario C: subject + 1 added piece of furniture =====
|
||||
new_furniture = [{"pos": [3.5, 1.0], "refl": 0.25, "name": "new_chair"}]
|
||||
y_C = simulate(subject_expected + walls + new_furniture, tx, rx, freq_ghz)
|
||||
pabs_C = pabs(y_C, y_expected)
|
||||
|
||||
# ===== Scenario D: subject + unexpected second person =====
|
||||
intruder = human_body(2.0, 2.0)
|
||||
y_D = simulate(subject_expected + walls + intruder, tx, rx, freq_ghz)
|
||||
pabs_D = pabs(y_D, y_expected)
|
||||
|
||||
# ===== Scenario E: subject + natural drift (wall reflectivity shift) =====
|
||||
# Walls have ~5% reflectivity drift over the day (humidity, temperature)
|
||||
drifted_walls = static_wall_reflectors(amplitudes=(0.315, 0.21, 0.158, 0.105))
|
||||
y_E = simulate(subject_expected + drifted_walls, tx, rx, freq_ghz)
|
||||
pabs_E = pabs(y_E, y_expected)
|
||||
|
||||
# ===== Scenario F: small subject position shift (subject moved 10 cm) =====
|
||||
subject_shifted = human_body(2.5, 2.85) # 10 cm closer to LOS
|
||||
y_F = simulate(subject_shifted + walls, tx, rx, freq_ghz)
|
||||
pabs_F = pabs(y_F, y_expected)
|
||||
|
||||
# ===== R12 NEGATIVE baseline: naive SVD cosine distance =====
|
||||
# Run the same scenarios through R12's failed approach for comparison.
|
||||
def svd_distance(y_obs, y_ref):
|
||||
# Treat as 1D signal; SVD spectrum on |y|
|
||||
return float(np.linalg.norm(np.abs(y_obs) - np.abs(y_ref)))
|
||||
|
||||
svd_A = svd_distance(y_empty, y_expected)
|
||||
svd_B = svd_distance(y_B, y_expected)
|
||||
svd_C = svd_distance(y_C, y_expected)
|
||||
svd_D = svd_distance(y_D, y_expected)
|
||||
svd_E = svd_distance(y_E, y_expected)
|
||||
svd_F = svd_distance(y_F, y_expected)
|
||||
|
||||
out = {
|
||||
"model": "PABS = ||y_observed - y_predicted||^2 / ||y_observed||^2",
|
||||
"forward_operator_source": "R6.1 multi-scatterer additive Fresnel",
|
||||
"expected_scene": {
|
||||
"subject_pos": [2.5, 2.75],
|
||||
"wall_reflectors": 4,
|
||||
},
|
||||
"link": {"tx": tx.tolist(), "rx": rx.tolist(), "freq_ghz": freq_ghz},
|
||||
"scenarios": {
|
||||
"A_empty_room": {"description": "no occupant", "pabs": pabs_A, "svd_distance": svd_A},
|
||||
"B_subject_expected": {"description": "subject where expected", "pabs": pabs_B, "svd_distance": svd_B},
|
||||
"C_added_furniture": {"description": "+1 new structural element", "pabs": pabs_C, "svd_distance": svd_C},
|
||||
"D_unexpected_person":{"description": "+1 unexpected human", "pabs": pabs_D, "svd_distance": svd_D},
|
||||
"E_natural_drift": {"description": "5%% wall reflectivity drift", "pabs": pabs_E, "svd_distance": svd_E},
|
||||
"F_subject_moved": {"description": "subject shifted 10 cm", "pabs": pabs_F, "svd_distance": svd_F},
|
||||
},
|
||||
"verdict": {
|
||||
"pabs_signal_to_drift": pabs_D / pabs_E if pabs_E > 0 else float("inf"),
|
||||
"pabs_furniture_to_drift": pabs_C / pabs_E if pabs_E > 0 else float("inf"),
|
||||
"svd_signal_to_drift": svd_D / svd_E if svd_E > 0 else float("inf"),
|
||||
"svd_furniture_to_drift": svd_C / svd_E if svd_E > 0 else float("inf"),
|
||||
},
|
||||
}
|
||||
Path(args.out).parent.mkdir(parents=True, exist_ok=True)
|
||||
Path(args.out).write_text(json.dumps(out, indent=2))
|
||||
|
||||
print("=== R12 PABS implementation results ===")
|
||||
print()
|
||||
print(f"{'Scenario':<30} {'PABS':>9} {'SVD':>9} {'PABS / drift':>14} {'SVD / drift':>13}")
|
||||
print("-" * 90)
|
||||
for key, s in out["scenarios"].items():
|
||||
pabs_ratio = s['pabs'] / pabs_E if pabs_E > 0 else float('inf')
|
||||
svd_ratio = s['svd_distance'] / svd_E if svd_E > 0 else float('inf')
|
||||
print(f"{s['description']:<30} {s['pabs']:>9.4f} {s['svd_distance']:>9.4f} "
|
||||
f"{pabs_ratio:>14.2f}x {svd_ratio:>13.2f}x")
|
||||
print()
|
||||
print(f"PABS detects unexpected person at {out['verdict']['pabs_signal_to_drift']:.1f}x the natural drift floor")
|
||||
print(f"PABS detects new furniture at {out['verdict']['pabs_furniture_to_drift']:.1f}x the natural drift floor")
|
||||
print(f"SVD (R12 naive) signal/drift: {out['verdict']['svd_signal_to_drift']:.2f}x")
|
||||
print(f"SVD (R12 naive) furniture/drift: {out['verdict']['svd_furniture_to_drift']:.2f}x")
|
||||
print()
|
||||
if out['verdict']['pabs_signal_to_drift'] > 3 and out['verdict']['svd_signal_to_drift'] < 2:
|
||||
print("VERDICT: PABS works where R12 naive SVD failed. R12 NEGATIVE -> revisited and POSITIVE.")
|
||||
elif out['verdict']['pabs_signal_to_drift'] > out['verdict']['svd_signal_to_drift'] * 2:
|
||||
print("VERDICT: PABS is meaningfully better than R12 naive SVD.")
|
||||
else:
|
||||
print("VERDICT: PABS is not yet decisive. Needs longer time-series / temporal averaging.")
|
||||
print()
|
||||
print(f"Wrote {args.out}")
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
main()
|
||||
@@ -0,0 +1,60 @@
|
||||
{
|
||||
"model": "PABS = ||y_observed - y_predicted||^2 / ||y_observed||^2",
|
||||
"forward_operator_source": "R6.1 multi-scatterer additive Fresnel",
|
||||
"expected_scene": {
|
||||
"subject_pos": [
|
||||
2.5,
|
||||
2.75
|
||||
],
|
||||
"wall_reflectors": 4
|
||||
},
|
||||
"link": {
|
||||
"tx": [
|
||||
0.0,
|
||||
2.5
|
||||
],
|
||||
"rx": [
|
||||
5.0,
|
||||
2.5
|
||||
],
|
||||
"freq_ghz": 2.4
|
||||
},
|
||||
"scenarios": {
|
||||
"A_empty_room": {
|
||||
"description": "no occupant",
|
||||
"pabs": 4.170183705070839,
|
||||
"svd_distance": 0.5965843005537784
|
||||
},
|
||||
"B_subject_expected": {
|
||||
"description": "subject where expected",
|
||||
"pabs": 0.0,
|
||||
"svd_distance": 0.0
|
||||
},
|
||||
"C_added_furniture": {
|
||||
"description": "+1 new structural element",
|
||||
"pabs": 0.04744306789447172,
|
||||
"svd_distance": 0.1011460778806426
|
||||
},
|
||||
"D_unexpected_person": {
|
||||
"description": "+1 unexpected human",
|
||||
"pabs": 0.6575620431155754,
|
||||
"svd_distance": 0.09866444424036849
|
||||
},
|
||||
"E_natural_drift": {
|
||||
"description": "5%% wall reflectivity drift",
|
||||
"pabs": 0.0005664412950287771,
|
||||
"svd_distance": 0.009233808950251039
|
||||
},
|
||||
"F_subject_moved": {
|
||||
"description": "subject shifted 10 cm",
|
||||
"pabs": 12.442629346878062,
|
||||
"svd_distance": 0.8354632981416396
|
||||
}
|
||||
},
|
||||
"verdict": {
|
||||
"pabs_signal_to_drift": 1160.8652986399395,
|
||||
"pabs_furniture_to_drift": 83.75637212689702,
|
||||
"svd_signal_to_drift": 10.685129481446127,
|
||||
"svd_furniture_to_drift": 10.953884623949552
|
||||
}
|
||||
}
|
||||
@@ -0,0 +1,181 @@
|
||||
#!/usr/bin/env python3
|
||||
"""R12 — RF weather: can SVD-eigenvalue drift detect structural changes?
|
||||
|
||||
See docs/research/sota-2026-05-22/R12-rf-weather-mapping.md.
|
||||
|
||||
The persistent-room field model in `wifi-densepose-signal/src/ruvsense/
|
||||
field_model.rs` does an SVD on empty-room CSI to extract an eigenstructure
|
||||
that describes "what this room's RF reflection looks like with nobody
|
||||
in it". Today that's used to subtract the room's baseline so motion
|
||||
detection isn't confused by static multipath.
|
||||
|
||||
This experiment asks a different question: **does the eigenvalue
|
||||
*spectrum* itself drift in a detectable way when something structural
|
||||
changes in the room?** "Structural change" = a new piece of furniture,
|
||||
a window that opened, water in the wall, settled foundation, missing
|
||||
ceiling tile. The 10-year vision (R12 research note) is continuous
|
||||
building-integrity monitoring from passive ambient WiFi.
|
||||
|
||||
Test:
|
||||
1. Take the existing 1,077 CSI windows. Split first 50% = "before",
|
||||
last 50% = "after".
|
||||
2. Inject a synthetic "structural perturbation" into the "after"
|
||||
half — multiply 3 subcarriers by 0.85 (simulating a new reflective
|
||||
surface that attenuates those frequencies).
|
||||
3. For each half, stack the windows into a `[N, 56]` per-frame
|
||||
matrix (each row = one timestep), compute SVD, take the top-10
|
||||
singular values.
|
||||
4. Measure: do the singular-value spectra differ in a way that
|
||||
distinguishes "structural perturbation present" from "no
|
||||
perturbation"?
|
||||
5. Repeat with NO perturbation as control — the same first-half /
|
||||
second-half split should produce *similar* spectra (just temporal
|
||||
drift from operator movement, not structural).
|
||||
|
||||
If the perturbed-vs-control eigenvalue spectra are distinguishable by
|
||||
a simple distance metric, RF-weather detection is feasible.
|
||||
"""
|
||||
|
||||
from __future__ import annotations
|
||||
|
||||
import argparse
|
||||
import json
|
||||
from pathlib import Path
|
||||
import numpy as np
|
||||
|
||||
N_SUB, N_FRAMES = 56, 20
|
||||
|
||||
|
||||
def load_windows(path: Path, max_samples: int | None = None) -> np.ndarray:
|
||||
csis = []
|
||||
with path.open(encoding="utf-8") as f:
|
||||
for line in f:
|
||||
if not line.strip():
|
||||
continue
|
||||
d = json.loads(line)
|
||||
shape = d.get("csi_shape", [N_SUB, N_FRAMES])
|
||||
if shape != [N_SUB, N_FRAMES]:
|
||||
continue
|
||||
csi = np.asarray(d["csi"], dtype=np.float32).reshape(N_SUB, N_FRAMES)
|
||||
csis.append(csi)
|
||||
if max_samples and len(csis) >= max_samples:
|
||||
break
|
||||
return np.stack(csis)
|
||||
|
||||
|
||||
def perturb_subcarriers(X: np.ndarray, indices: list[int], gain: float) -> np.ndarray:
|
||||
"""Multiply the listed subcarriers by `gain` to simulate a structural
|
||||
change (e.g. a new reflector attenuates certain frequencies)."""
|
||||
out = X.copy()
|
||||
out[:, indices, :] *= gain
|
||||
return out
|
||||
|
||||
|
||||
def per_frame_matrix(X: np.ndarray) -> np.ndarray:
|
||||
"""Stack all windows' frames into a [N_total_frames, 56] matrix.
|
||||
Each row is one timestep, used as a multivariate observation of the
|
||||
56-subcarrier channel state."""
|
||||
return X.transpose(0, 2, 1).reshape(-1, N_SUB)
|
||||
|
||||
|
||||
def top_k_singular_values(M: np.ndarray, k: int = 10) -> np.ndarray:
|
||||
"""Compute SVD on M, return top-k singular values."""
|
||||
M_centered = M - M.mean(axis=0, keepdims=True)
|
||||
# Use SVD on the centered matrix (== PCA without normalisation)
|
||||
s = np.linalg.svd(M_centered, compute_uv=False)
|
||||
return s[:k]
|
||||
|
||||
|
||||
def spectrum_distance(s1: np.ndarray, s2: np.ndarray) -> float:
|
||||
"""Cosine distance between two singular-value spectra. 0 = identical
|
||||
direction, 2 = opposite. Symmetric, scale-invariant."""
|
||||
s1n = s1 / (np.linalg.norm(s1) + 1e-9)
|
||||
s2n = s2 / (np.linalg.norm(s2) + 1e-9)
|
||||
return float(1.0 - np.dot(s1n, s2n))
|
||||
|
||||
|
||||
def main():
|
||||
parser = argparse.ArgumentParser()
|
||||
parser.add_argument("--paired", required=True)
|
||||
parser.add_argument("--out", default="examples/research-sota/r12_rf_weather_results.json")
|
||||
parser.add_argument("--perturb-indices", default="30,41,52",
|
||||
help="comma-separated subcarrier indices to perturb (chosen from R5's top-saliency list)")
|
||||
parser.add_argument("--perturb-gain", type=float, default=0.85)
|
||||
args = parser.parse_args()
|
||||
|
||||
print(f"Loading windows from {args.paired}")
|
||||
X = load_windows(Path(args.paired))
|
||||
print(f" total windows: {X.shape[0]} (shape {X.shape})")
|
||||
|
||||
n = X.shape[0]
|
||||
half = n // 2
|
||||
X_before = X[:half]
|
||||
X_after_raw = X[half:] # unmodified second half — the CONTROL
|
||||
perturb_idx = [int(x) for x in args.perturb_indices.split(",")]
|
||||
X_after_perturbed = perturb_subcarriers(X_after_raw, perturb_idx, args.perturb_gain)
|
||||
|
||||
# Convert each half to a [N_frames, 56] matrix
|
||||
M_before = per_frame_matrix(X_before)
|
||||
M_after_raw = per_frame_matrix(X_after_raw)
|
||||
M_after_pert = per_frame_matrix(X_after_perturbed)
|
||||
print(f" per-frame matrix: before={M_before.shape}, after={M_after_raw.shape}")
|
||||
|
||||
# Top-10 singular values per half
|
||||
s_before = top_k_singular_values(M_before, k=10)
|
||||
s_after_raw = top_k_singular_values(M_after_raw, k=10)
|
||||
s_after_pert = top_k_singular_values(M_after_pert, k=10)
|
||||
|
||||
print(f"\n Singular value spectra (top-10):")
|
||||
print(f" before : [{', '.join(f'{v:.1f}' for v in s_before)}]")
|
||||
print(f" after (raw) : [{', '.join(f'{v:.1f}' for v in s_after_raw)}]")
|
||||
print(f" after (pert) : [{', '.join(f'{v:.1f}' for v in s_after_pert)}]")
|
||||
|
||||
# Distances
|
||||
d_raw = spectrum_distance(s_before, s_after_raw)
|
||||
d_pert = spectrum_distance(s_before, s_after_pert)
|
||||
|
||||
print(f"\n Cosine distances from BEFORE:")
|
||||
print(f" before -> after raw (control, no perturbation): {d_raw:.5f}")
|
||||
print(f" before -> after pert (synthetic structural shift): {d_pert:.5f}")
|
||||
|
||||
# Distance ratio = how much the perturbation amplifies the detection signal
|
||||
# over the natural temporal drift.
|
||||
if d_raw > 1e-9:
|
||||
ratio = d_pert / d_raw
|
||||
print(f"\n Signal-to-natural-drift ratio: {ratio:.2f}x")
|
||||
|
||||
if d_pert > d_raw * 3:
|
||||
verdict = "STRONG: perturbation easily distinguishable from natural temporal drift"
|
||||
elif d_pert > d_raw * 1.5:
|
||||
verdict = "MODERATE: perturbation detectable but with margin"
|
||||
else:
|
||||
verdict = "WEAK: structural perturbation gets lost in temporal drift"
|
||||
print(f"\n Verdict: {verdict}")
|
||||
|
||||
out = {
|
||||
"perturbation": {
|
||||
"subcarrier_indices": perturb_idx,
|
||||
"amplitude_gain": args.perturb_gain,
|
||||
"comment": "simulates a new reflective surface that attenuates these frequencies",
|
||||
},
|
||||
"n_before_windows": int(half),
|
||||
"n_after_windows": int(n - half),
|
||||
"spectra": {
|
||||
"before": s_before.tolist(),
|
||||
"after_raw_control": s_after_raw.tolist(),
|
||||
"after_perturbed": s_after_pert.tolist(),
|
||||
},
|
||||
"distances": {
|
||||
"before_to_after_raw": d_raw,
|
||||
"before_to_after_perturbed": d_pert,
|
||||
"signal_over_natural_drift": float(d_pert / max(d_raw, 1e-9)),
|
||||
},
|
||||
"verdict": verdict,
|
||||
}
|
||||
Path(args.out).parent.mkdir(parents=True, exist_ok=True)
|
||||
Path(args.out).write_text(json.dumps(out, indent=2))
|
||||
print(f"\nWrote {args.out}")
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
main()
|
||||
@@ -0,0 +1,57 @@
|
||||
{
|
||||
"perturbation": {
|
||||
"subcarrier_indices": [
|
||||
30,
|
||||
41,
|
||||
52
|
||||
],
|
||||
"amplitude_gain": 0.85,
|
||||
"comment": "simulates a new reflective surface that attenuates these frequencies"
|
||||
},
|
||||
"n_before_windows": 538,
|
||||
"n_after_windows": 539,
|
||||
"spectra": {
|
||||
"before": [
|
||||
2220.65673828125,
|
||||
1856.8695068359375,
|
||||
1563.7314453125,
|
||||
1303.56298828125,
|
||||
1057.757080078125,
|
||||
770.67822265625,
|
||||
757.5601196289062,
|
||||
689.5866088867188,
|
||||
595.6748046875,
|
||||
556.3777465820312
|
||||
],
|
||||
"after_raw_control": [
|
||||
2182.5712890625,
|
||||
1837.5084228515625,
|
||||
1647.6357421875,
|
||||
1315.103759765625,
|
||||
1053.489013671875,
|
||||
794.1417236328125,
|
||||
737.1859130859375,
|
||||
704.1968994140625,
|
||||
571.363037109375,
|
||||
535.6047973632812
|
||||
],
|
||||
"after_perturbed": [
|
||||
2172.6552734375,
|
||||
1824.164794921875,
|
||||
1615.7850341796875,
|
||||
1304.227783203125,
|
||||
1040.461181640625,
|
||||
791.2919921875,
|
||||
736.2902221679688,
|
||||
691.3584594726562,
|
||||
568.5400390625,
|
||||
530.7666625976562
|
||||
]
|
||||
},
|
||||
"distances": {
|
||||
"before_to_after_raw": 0.0003509521484375,
|
||||
"before_to_after_perturbed": 0.00024056434631347656,
|
||||
"signal_over_natural_drift": 0.6854619565217391
|
||||
},
|
||||
"verdict": "WEAK: structural perturbation gets lost in temporal drift"
|
||||
}
|
||||
Reference in New Issue
Block a user