sync.py
Rama main · réplica de lectura del repositorio original.
"""
Stage 1: FT8 signal candidate detection (sync search).
Method: compute a spectrogram (STFT) of the whole audio window with time
oversampling (e.g. every 1/2 symbol) and frequency oversampling (every 1/2
tone), then slide the 7-symbol Costas pattern across the (time, frequency)
grid and look for local correlation maxima (sum of power on the expected
Costas tones minus the average power on the other tones in that window).
An FT8 signal has THREE occurrences of the Costas pattern (at symbol
positions 0, 36, 72), so we sum the correlation across all three for more
reliable detection in noise.
"""
import numpy as np
from params import (SAMPLE_RATE, SAMPLES_PER_SYMBOL, TONE_SPACING, N_TONES,
COSTAS, COSTAS_POS, N_SYM)
def compute_magnitude_spectrogram(audio, freq_osr=2, time_osr=2, f_min=200, f_max=3000):
"""
Computes an oversampled power spectrogram.
freq_osr: frequency oversampling relative to TONE_SPACING (2 = 3.125Hz
steps), achieved via ZERO-PADDING (not lengthening the window), so
each FFT window still covers exactly 1 symbol and doesn't mix
energy from adjacent symbols.
time_osr: time oversampling relative to SAMPLES_PER_SYMBOL (2 = half-symbol steps)
Returns: (mag, n_blocks, n_bins, freq_step, time_step_samples, bin_min)
mag[time_block, freq_bin] = power
"""
n = SAMPLES_PER_SYMBOL # analysis window length = ALWAYS 1 symbol
nfft = n * freq_osr # zero-pad to this length (frequency resolution)
freq_step = SAMPLE_RATE / nfft # Hz per FFT bin
time_step = SAMPLES_PER_SYMBOL // time_osr # samples between consecutive STFT windows
n_samples = len(audio)
n_blocks = (n_samples - n) // time_step + 1
if n_blocks < 1:
return None
window = np.hanning(n)
bin_min = int(f_min / freq_step)
bin_max = int(f_max / freq_step)
mag = np.zeros((n_blocks, bin_max - bin_min), dtype=np.float64)
for i in range(n_blocks):
start = i * time_step
seg = audio[start:start + n]
if len(seg) < n:
seg = np.pad(seg, (0, n - len(seg)))
spec = np.fft.rfft(seg * window, n=nfft) # zero-padding here, window is still length n
power = np.abs(spec) ** 2
mag[i, :] = power[bin_min:bin_max]
return mag, n_blocks, mag.shape[1], freq_step, time_step, bin_min
def find_candidates(audio, freq_osr=2, time_osr=2, f_min=200, f_max=3000,
max_candidates=30, min_score=0.4):
"""
Searches audio for FT8 signal candidates.
Returns a list of dicts: {freq_hz, time_offset_s, score, block0, bin0}
Scoring robust to overlapping stations: for each Costas position we
compute (power_of_correct_tone - average_of_8_tones) / average_of_8_tones
on LINEAR power. The old method (log(expected) - log(max_of_the_other_7))
penalized neighbors: a strong signal nearby raised the max and pushed the
margin below threshold, losing even strong stations (e.g. E20LXN -3dB
with a neighbor at 2009Hz). The "relative to the average" method is
robust to a single neighbor.
Vectorized version: 3D arrays, one numpy operation per Costas term.
"""
result = compute_magnitude_spectrogram(audio, freq_osr, time_osr, f_min, f_max)
if result is None:
return []
mag, n_blocks, n_bins, freq_step, time_step, bin_min = result
bins_per_tone = freq_osr
sym_blocks = time_osr
max_time_block = n_blocks - N_SYM * sym_blocks
max_freq_bin = n_bins - (N_TONES - 1) * bins_per_tone
if max_time_block < 1 or max_freq_bin < 1:
return []
t0_step = max(1, sym_blocks // 2)
f0_step = max(1, bins_per_tone // 2)
t0_vals = np.arange(0, max_time_block, t0_step)
f0_vals = np.arange(0, max_freq_bin, f0_step)
n_t0 = len(t0_vals)
n_f0 = len(f0_vals)
score_accum = np.zeros((n_t0, n_f0), dtype=np.float64)
count_accum = 0
for costas_sym_offset in COSTAS_POS:
for k, tone in enumerate(COSTAS):
sym_idx = costas_sym_offset + k
tblocks = t0_vals + sym_idx * sym_blocks # (n_t0,)
valid_t = tblocks < n_blocks
if not np.any(valid_t):
continue
# For each f0 (n_f0,) and each tone (8,), compute the bin index
tone_offsets = np.arange(N_TONES) * bins_per_tone # (8,)
fbins = f0_vals[:, None] + tone_offsets[None, :] # (n_f0, 8)
valid_f = np.all((fbins >= 0) & (fbins < n_bins), axis=1) # (n_f0,)
if not np.any(valid_f):
continue
# Slice out the power (LINEAR, not log) for (tblock, fbin): (n_t0, n_f0, 8)
tb_idx = np.clip(tblocks, 0, n_blocks - 1)
fb_idx = np.clip(fbins, 0, n_bins - 1)
block_rows = mag[tb_idx] # (n_t0, n_bins)
powers = block_rows[:, fb_idx] # (n_t0, n_f0, 8)
expected = powers[:, :, tone] # (n_t0, n_f0)
mean_all = np.mean(powers, axis=2) # (n_t0, n_f0)
# contribution: how much the correct tone exceeds the AVERAGE (normalized).
# Robust to 1 neighbor (a neighbor only raises the average by ~1/8).
with np.errstate(divide='ignore', invalid='ignore'):
contrib = np.where(mean_all > 1e-12,
(expected - mean_all) / mean_all, 0.0)
valid_mask = valid_t[:, None] & valid_f[None, :]
score_accum += np.where(valid_mask, contrib, 0.0)
count_accum += 1
if count_accum == 0:
return []
score_map = score_accum / (len(COSTAS_POS) * len(COSTAS))
candidates = []
# Non-max suppression: find local maxima in score_map
flat_idx = np.argsort(score_map, axis=None)[::-1]
taken = np.zeros_like(score_map, dtype=bool)
# suppression radius expressed in score_map INDEX UNITS (not t0_step/f0_step
# directly from samples), so we convert it against the actual grid step
suppress_t = max(1, (sym_blocks * 2) // t0_step)
suppress_f = max(1, (bins_per_tone * 2) // f0_step)
for idx in flat_idx:
i_t, i_f = np.unravel_index(idx, score_map.shape)
sc = score_map[i_t, i_f]
if sc < min_score or not np.isfinite(sc):
break
if taken[i_t, i_f]:
continue
t_lo, t_hi = max(0, i_t - suppress_t), min(n_t0, i_t + suppress_t + 1)
f_lo, f_hi = max(0, i_f - suppress_f), min(n_f0, i_f + suppress_f + 1)
if taken[t_lo:t_hi, f_lo:f_hi].any():
continue
taken[t_lo:t_hi, f_lo:f_hi] = True
f0 = f0_vals[i_f]
t0 = t0_vals[i_t]
freq_hz = (f0 + bin_min) * freq_step
time_offset_s = t0 * time_step / SAMPLE_RATE
candidates.append({
'freq_hz': freq_hz,
'time_offset_s': time_offset_s,
'score': float(sc),
'block0': t0,
'bin0': f0 + bin_min,
})
if len(candidates) >= max_candidates:
break
candidates.sort(key=lambda c: -c['score'])
return candidates