121 lines
5.1 KiB
Python
121 lines
5.1 KiB
Python
"""Ad-hoc probe for the signal-chain float model (detection firing + Doppler sign)."""
|
|
import math
|
|
import sys
|
|
|
|
sys.path.insert(0, "src")
|
|
from defense_sim.sensors import signal_chain as sc
|
|
import numpy as np
|
|
|
|
cfg = sc.RadarSignalConfig()
|
|
|
|
|
|
def noiseless_frame(proc, v, range_m=1000.0):
|
|
ci = np.zeros((32, 64))
|
|
cq = np.zeros((32, 64))
|
|
for c in range(32):
|
|
pi_, pq_, rb = proc.target_row(range_m, v, cfg.rcs_dbsm_default, c)
|
|
for delta, wgt in ((0, 1.0), (-1, .2), (1, .2), (-2, .1), (2, .1)):
|
|
row = rb + delta
|
|
if 0 <= row < 64:
|
|
ci[c, row] += wgt * pi_
|
|
cq[c, row] += wgt * pq_
|
|
return ci, cq
|
|
|
|
|
|
# 1) single-frame behaviour (noiseless + noisy), Doppler sign sanity
|
|
# Velocities stay within the single-PRI unambiguous Doppler range (~±42 m/s long
|
|
# PRI) so the wrapped estimate must match truth without unwrapping.
|
|
for v in (30.0, -30.0, 0.0):
|
|
proc = sc.ReceiverProcessor(np.random.default_rng(7), cfg)
|
|
f = proc.generate_frame([{"range_m": 1000.0, "radial_velocity_m_s": v}])
|
|
dets = proc.process_frame(*f)
|
|
top = sorted(dets, key=lambda d: -d["snr_db"])[:4]
|
|
print(f"noiseless v={v:+.0f}: n_dets={len(dets)} top_bins={[t['doppler_bin'] for t in top]} "
|
|
f"vels={[round(t['radial_velocity_m_s'], 2) for t in top]}")
|
|
|
|
for v in (50.0, -50.0):
|
|
proc = sc.ReceiverProcessor(np.random.default_rng(11), cfg)
|
|
ci, cq = proc.generate_frame([dict(range_m=1000.0, radial_velocity_m_s=v)])
|
|
dets = proc.process_frame(ci, cq)
|
|
top = sorted(dets, key=lambda d: -d["snr_db"])[:4]
|
|
print(f"noisy v={v:+.0f}: n_dets={len(dets)} top_bins={[t['doppler_bin'] for t in top]} "
|
|
f"vels={[round(t['radial_velocity_m_s'], 2) for t in top]}")
|
|
|
|
# 2) three successive dwell frames with scene motion (range shrinking)
|
|
proc = sc.ReceiverProcessor(np.random.default_rng(42), cfg)
|
|
for step, (rng_m, v) in enumerate([(1500.0, 60.0), (1440.0, 60.0), (1380.0, 60.0)]):
|
|
ci, cq = proc.generate_frame([dict(range_m=rng_m, radial_velocity_m_s=v)])
|
|
dets = proc.process_frame(ci, cq)
|
|
top = sorted(dets, key=lambda x: -x["snr_db"])[:3]
|
|
print(f"frame{step}: n={len(dets)} top={[(t['range_bin'], t['doppler_bin'], round(t['radial_velocity_m_s'], 1)) for t in top]}")
|
|
|
|
# 3) intermediate map stats (single noisy frame, noise floor 0.5)
|
|
proc2 = sc.ReceiverProcessor(np.random.default_rng(7), cfg)
|
|
ci, cq = proc2.generate_frame([dict(range_m=1000.0, radial_velocity_m_s=50.0)])
|
|
chirp_stack = ci + 1j * cq
|
|
gw = sc.FFT_SIZE // sc.RANGE_BINS
|
|
padded = np.zeros((32, 64, gw), dtype=np.complex128)
|
|
padded[:, :, 0] = chirp_stack
|
|
range_fft = np.fft.fft(padded.reshape(32, sc.FFT_SIZE), axis=1)
|
|
rows = []
|
|
for chirp_row in range_fft:
|
|
groups = chirp_row.reshape(64, gw)
|
|
mags = np.abs(groups.real) + np.abs(groups.imag)
|
|
winners = np.argmax(mags, axis=1)
|
|
rows.append(groups[np.arange(64), winners])
|
|
decimated = np.stack(rows)
|
|
mti = np.zeros_like(decimated)
|
|
mti[1:] = decimated[1:] - decimated[:-1]
|
|
doppler_map = np.zeros((64, 32), dtype=np.complex128)
|
|
for sf in (0, 1):
|
|
block = mti[sf * 16:(sf + 1) * 16, :]
|
|
windowed = block.real * proc2.hamming[:, None] + 1j * (block.imag * proc2.hamming[:, None])
|
|
doppler_map[:, sf * 16:(sf + 1) * 16] = np.fft.fft(windowed, axis=0).T
|
|
for dbin in range(32):
|
|
within = dbin & 0xF
|
|
if within < 2 or within > (15 - 2 + 1):
|
|
doppler_map[:, dbin] = 0
|
|
mags_map = np.abs(doppler_map.real) + np.abs(doppler_map.imag)
|
|
print("mti decim std:", float(np.std(mti)))
|
|
print("doppler map mag stats: mean", mags_map.mean(), "max", mags_map.max(), "p99.5", np.percentile(mags_map, 99.5))
|
|
|
|
# noise-only background map for comparison
|
|
proc3 = sc.ReceiverProcessor(np.random.default_rng(7), cfg)
|
|
ci0, cq0 = proc3.generate_frame([])
|
|
chirp_stack0 = ci0 + 1j * cq0
|
|
padded0 = np.zeros((32, 64, gw), dtype=np.complex128)
|
|
padded0[:, :, 0] = chirp_stack0
|
|
rf0 = np.fft.fft(padded0.reshape(32, sc.FFT_SIZE), axis=1)
|
|
rows0 = []
|
|
for chirp_row in rf0:
|
|
groups = chirp_row.reshape(64, gw)
|
|
mags = np.abs(groups.real) + np.abs(groups.imag)
|
|
winners = np.argmax(mags, axis=1)
|
|
rows0.append(groups[np.arange(64), winners])
|
|
dec0 = np.stack(rows0)
|
|
mti0 = np.zeros_like(dec0)
|
|
mti0[1:] = dec0[1:] - dec0[:-1]
|
|
dop0 = np.zeros((64, 32), dtype=np.complex128)
|
|
for sf in (0, 1):
|
|
block = mti0[sf * 16:(sf + 1) * 16, :]
|
|
windowed = block.real * proc2.hamming[:, None] + 1j * (block.imag * proc2.hamming[:, None])
|
|
dop0[:, sf * 16:(sf + 1) * 16] = np.fft.fft(windowed, axis=0).T
|
|
for dbin in range(32):
|
|
within = dbin & 0xF
|
|
if within < 2 or within > (15 - 2 + 1):
|
|
dop0[:, dbin] = 0
|
|
mags0 = np.abs(dop0.real) + np.abs(dop0.imag)
|
|
print("noise-only map: mean", mags0.mean(), "max", mags0.max())
|
|
print("target+noise map at (42,9):", mags_map[42, 9], " noise-only:", mags0[42, 9])
|
|
|
|
# 4) sigma sanity + false-alarm budget
|
|
fa = 0
|
|
frames = 0
|
|
for seed in range(5):
|
|
p = sc.ReceiverProcessor(np.random.default_rng(100 + seed), cfg)
|
|
ci_o, cq_o = p.generate_frame([])
|
|
d = p.process_frame(ci_o, cq_o)
|
|
fa += len(d)
|
|
frames += 1
|
|
print("noise-only false alarms per frame:", fa / frames, "(budget: ~", cfg.false_alarm_probability, ")")
|