Patrick Lidstone
Self-hosted

WSPR: a native Weak-Signal Propagation Reporter codec

A clean-room, pure-NumPy implementation of the WSPR beacon mode — 50-bit message packing, a K=32 rate-½ convolutional code, bit-reversal interleaving, and continuous-phase 4-FSK at 1.4648 baud over a 2-minute slot — with an FFT-per-symbol demodulator and a stack (Zigangirov–Jelinek) sequential FEC decoder.

Rafe project · app/radio/wsprx/{constants,encode,decode,codec}.py · in-repo replacement for wsprd


Abstract

WSPR (pronounced "whisper", Weak Signal Propagation Reporter) is Joe Taylor K1JT's beacon protocol for probing HF/MF propagation with a handful of milliwatts. A fixed low-rate message — callsign, 4-character Maidenhead locator, transmit power in dBm — is packed to exactly 50 bits, protected by a constraint-length-32, rate-½ convolutional code (81 information+tail bits → 162 code bits), scrambled by a bit-reversal interleaver, merged bit-by-bit with a fixed 162-element synchronisation vector, and transmitted as continuous-phase 4-FSK at a symbol rate of 1.4648 baud with 1.4648 Hz tone spacing over one 2-minute transmit/receive slot. This document specifies the exact implementation in app/radio/wsprx/, which is a clean-room build from the published WSPR coding notes (K1JT, G4JNT) — not a port of the GPLv3 WSJT-X / wsprd source. It is a drop-in replacement for the external wsprd binary: encode_pcm(call, grid, dbm) renders a 120-second s16le mono PCM slot, and decode_slot(pcm) recovers the message. The encoder is symbol-exact against WSJT-X's own wsprcode/wsprd oracles (6/6 on prod); the decoder is a correct FFT-demod + stack-decoder baseline that reads down to roughly −28 dB, a few dB shy of wsprd's deep-search sensitivity. Everything runs in pure NumPy on a Raspberry-Pi-class box with no build dependencies.


1. Background

1.1 What WSPR is for

WSPR is a propagation-report mode, not a QSO mode. Stations transmit a beacon on a 2-minute schedule (aligned to even UTC minutes); receivers decode every beacon they can hear and upload the spots (who heard whom, at what SNR and frequency offset) to a central database. Because the message content is fixed and known-format — CALL GRID dBm — the protocol can spend its entire information budget on FEC and integration gain, which is why WSPR decodes far below the noise floor (published sensitivity around −28 to −34 dB in a 2.5 kHz reference bandwidth). The design trade is stark: ~50 bits of payload spread across ~111 seconds of signal in ~6 Hz of occupied bandwidth.

1.2 The signal at a glance

Every WSPR transmission is the same length and structure: 162 four-level FSK symbols, each 8192 samples long at a 12 kHz working sample rate — 0.6827 s per symbol, 110.592 s of signal, launched ~1 s into the slot and framed by trailing silence to fill 120 s. The four tones are spaced 12000/8192 = 1.46484375 Hz apart (the tone spacing equals the symbol rate — a minimum-shift-like orthogonal 4-FSK), so the total occupied bandwidth is about 4 × 1.46 ≈ 6 Hz. Nominal audio sub-carrier is 1500 Hz.

1.3 Why this structure decodes so well

Two facts make the DSP clean and are exploited directly by the implementation:

  • The tone spacing is exactly one FFT bin. A symbol is 8192 samples; a 8192-point FFT at 12 kHz has bin width 12000/8192 = 1.46484375 Hz — identical to the tone spacing. Each of the four tones therefore lands exactly on its own FFT bin, and because a tone fills the whole symbol, a rectangular window is the matched filter (no window loss, maximum SNR). Per-symbol tone powers are just four adjacent FFT bin magnitudes-squared.
  • The sync vector carries no data but is fully known. Half of every symbol's identity (its least-significant bit) is a fixed pseudo-random sequence, giving the receiver a strong, dense timing/frequency reference across the whole transmission — a matched sync template 162 symbols long.

1.4 Provenance and licence note

WSJT-X (which contains wsprd) is GPLv3; this project is PolyForm Noncommercial. The codec here is therefore implemented clean-room from published protocol descriptions, and validated behaviourally against WSJT-X's reference utilities (wsprcode, wsprd) as black-box oracles — observing input→output, never copying source. See docs/native-digimodes.md.


2. Signal, frame, and symbol structure (exact numbers)

All constants below are from app/radio/wsprx/constants.py.

Quantity Value Source constant
Message payload 50 bits (28-bit call + 22-bit grid/power) MSG_BITS = 50
FEC K=32, rate ½ convolutional CONV_K = 32
Generator polynomials 0xF2D05351, 0xE4613C47 CONV_POLY
Tail (flush) bits 31 zeros = K−1 TAIL_BITS = CONV_K - 1
Coded/channel symbols 162 = (50+31)×2 NSYM = 162
Modulation 4-FSK (tones 0…3) channel_symbols (encode.py)
Working sample rate 12000 Hz SAMPLE_RATE = 12000
Samples per symbol 8192 SPS = 8192
Tone spacing 1.46484375 Hz = 12000/8192 TONE_SPACING
Symbol rate (baud) 1.46484375 = 12000/8192 (spacing = baud)
Symbol duration 0.68266… s = 8192/12000
Signal duration 110.592 s = 162×8192/12000
Slot length 120 s SLOT_SECONDS = 120
Nominal audio sub-carrier 1500 Hz (default f0) synth (encode.py)
Occupied bandwidth ≈ 5.86 Hz (4 tones × spacing)

2.1 Symbol construction — the sync merge

Each channel symbol is a 2-bit value 0…3 formed by merging one interleaved data bit with one fixed sync bit (channel_symbols, encode.py):

Symbol[n] = SYNC[n] + 2 · Data[n]          n = 0 … 161
  • SYNC[n] ∈ {0,1} is the fixed 162-element synchronisation vector (§6.2).
  • Data[n] ∈ {0,1} is the interleaved coded bit.
  • Therefore the low bit of every symbol is SYNC[n] (tones 0/2 ⇒ sync 0, tones 1/3 ⇒ sync 1) and the high bit is the data bit (tones 0/1 ⇒ data 0, tones 2/3 ⇒ data 1). The demodulator relies on exactly this split (§5).

The four symbol values map to tone frequencies f0 + Symbol·TONE_SPACING, i.e. for f0 = 1500 Hz: 1500.000, 1501.465, 1502.930, 1504.395 Hz.

2.2 Slot timing

synth_slot (encode.py) places the 110.592 s signal into a 120 s slot:

[ 1.000 s lead-in silence | 110.592 s of 4-FSK signal | 8.408 s trailing silence ]
     12000 samples              1 327 104 samples            100 896 samples
                          total = 1 440 000 samples = 120 s

The 1 s lead-in (lead_s=1.0) matches WSPR's convention of starting transmission ~1 s after the even-minute boundary. The decoder's coarse time search (§5.1) spans 0–2.5 s to find it.


3. Message packing (50 bits)

The Type-1 (standard) WSPR message is CALLSIGN GRID POWER — a standard callsign, a 4-character Maidenhead locator, and power in dBm. It packs into a 28-bit callsign field and a 22-bit grid+power field.

3.1 Character value map

_cv(ch) (encode.py): space → 36; '0'..'9' → 0…9; 'A'..'Z' → 10…35.

3.2 Callsign normalisation

_normalize_call (encode.py) forces a digit into the third character position (index 2), which is what the packing arithmetic assumes:

  1. Uppercase and strip.
  2. If call[2] is not a digit but call[1] is, prepend a space (shifts a 2×2-prefix call like M0SUP" M0SUP", W1AW" W1AW").
  3. Right-pad with spaces and truncate to exactly 6 characters.
  4. If call[2] is still not a digit, raise ValueError (not WSPR-encodable).

Examples (from test_wspr.py): "G4JNT"→" G4JNT", "M0SUP"→" M0SUP", "PA0XYZ"→"PA0XYZ", "W1AW"→" W1AW ".

3.3 Callsign packing → 28 bits

pack_call (encode.py), on the 6-char normalised call c[0..5]:

n = cv(c0)
n = n·36 + cv(c1)
n = n·10 + cv(c2)
n = n·27 + (cv(c3) − 10)
n = n·27 + (cv(c4) − 10)
n = n·27 + (cv(c5) − 10)
return n & 0x0FFFFFFF          # 28 bits

The radices reflect the allowed alphabet per position: c0 ∈ {0-9,A-Z,space} (37), c1 ∈ {0-9,A-Z} (36), c2 ∈ {0-9} (10), c3..c5 ∈ {A-Z,space} (27 each; the −10 maps A…Z,space → 0…26). Maximum 37·36·10·27³ = 262 177 560 < 2²⁸. Worked example: pack_call("M0SUP") = 259 435 617 = 0x0F76AC61 (test_wspr.py).

3.4 Grid + power packing → 22 bits

pack_grid_power (encode.py), grid = 4-char Maidenhead G0 G1 G2 G3, dbm integer:

m1 = (179 − 10·(G0−'A') − (G2−'0'))·180 + 10·(G1−'A') + (G3−'0')
m  = m1·128 + dbm + 64
return m & 0x3FFFFF             # 22 bits

The 179 − … term flips the longitude field so grids run west-to-east; dbm + 64 offsets power into an unsigned field. Worked example: IO90, 30 dBm → m1 = (179−80−9)·180 + 140 + 0 = 16340; m = 16340·128 + 30 + 64 = 2 091 614 = 0x1FEA5E (test_wspr.py).

3.5 Bit assembly

_message_bits (encode.py) emits the 50 payload bits MSB-first: 28 bits of n (bit 27 → bit 0) followed by 22 bits of m (bit 21 → bit 0).


4. Encode — step by step

The full chain is channel_symbolssynth_slot (encode.py), wrapped by encode_pcm (codec.py).

4.1 Convolutional encoding (K=32, rate ½)

_convolutional (encode.py) appends 31 zero tail bits and runs the 81-bit stream through a 32-bit shift register, emitting two code bits per input bit:

reg = 0
for b in (50 message bits) + (31 zeros):
    reg = ((reg << 1) | b) & 0xFFFFFFFF          # 32-bit register, K=32
    out.append( parity( reg & 0xF2D05351 ) )     # generator 0
    out.append( parity( reg & 0xE4613C47 ) )     # generator 1
# 81 inputs × 2 = 162 output bits

parity(x) is popcount(x) & 1 (bin(...).count("1") & 1). The new input bit enters at the LSB; the two generator polynomials tap the 32-bit register. Result: 162 coded bits in encoder order.

4.2 Interleaving (bit-reversal permutation)

_interleave (encode.py) spreads the 162 coded bits by a bit-reversal address permutation. Counter i runs 0…255; each i is reversed as an 8-bit value by _bitrev8 (encode.py); addresses j = bitrev8(i) that fall in [0,162) are used as destination slots, in ascending-i order, consuming the coded bits sequentially:

p = 0
for i in 0 … 255:
    j = bitrev8(i)                 # reverse the 8 bits of i
    if j < 162:
        d[j] = coded[p]; p += 1    # p-th coded bit → position j
        if p == 162: break

Bit-reversal of an 8-bit index reaches 0…255; the 94 addresses ≥ 162 are skipped (the loop consumes i up to 254, skipping 93 of them). The forward map (coded-bit index p → symbol position j) begins and ends:

p:  0    1   2   3    4   5   6    7   8   9  10  11 …
j:  0  128  64  32  160  96  16  144  80  48 112   8 …
… 154 155 156 157 158 159 160 161
…  79  47 111  31 159  95  63 127

(Full table regenerated exactly by the bitrev8 recipe above; representative excerpt only.) The decoder inverts this with _deinterleave (decode.py), which reads out[p] = vals[j] over the identical i/j walk.

4.3 Sync merge → 162 channel symbols

channel_symbols (encode.py) merges the interleaved data with the sync vector: Symbol[n] = SYNC[n] + 2·data[n] (§2.1). Regression fixture: channel_symbols("M0SUP", "IO90", 30) gives the 162-symbol vector M0SUP_SYMS in test_wspr.py, verified by decoding the synthesised signal with wsprd (WSJT-X 2.7.0) on prod.

4.4 4-FSK synthesis

synth (encode.py) builds a continuous-phase waveform (no per-symbol phase resets, so the spectrum is clean):

freqs[k] = f0 + Symbol[⌊k/8192⌋] · 1.46484375        # per-sample instantaneous freq
phase[k] = 2π · cumsum(freqs)[k] / 12000             # running phase accumulator
x[k]     = sin(phase[k])                              # float32 in [-1, 1]

synth_slot (encode.py) inserts 1 s of lead-in silence and pads to 120 s (§2.2). encode_pcm (codec.py) clips to [−1,1], scales by 32767, rounds, and emits little-endian int16 (<i2) bytes — a 2 880 000-byte (120 s × 12000 × 2) mono PCM slot.

Note — TX purity. This is plain continuous-phase FSK. WSJT-X additionally applies Gaussian (GFSK) smoothing to the frequency steps for a narrower transmitted spectrum. The signals decode identically; GFSK is a cleanliness refinement, not a protocol requirement.


5. Decode — step by step

decode_slot (decode.py) is the top level: parse PCM → find_syncsoft_symbols_deinterleavesequential_decodeunpack.

Input s16le bytes are read as <i2 and scaled to float by /32768.0 (decode.py).

5.1 Per-symbol tone powers (matched-filter DSP)

_tone_power_matrix (decode.py) builds a 162 × 4096 power matrix. For each symbol n, it takes the 8192-sample block starting at t0 + n·8192 (zero-padded if the tail runs short), computes an 8192-point real FFT (np.fft.rfft), and stores real² + imag² for the 4096 positive-frequency bins:

block = sig[t0 + n·8192 : t0 + (n+1)·8192]      # rectangular window = matched filter
spec  = rfft(block)
mat[n, :] = spec.real² + spec.imag²             # bins 0 … 4095

Because a tone fills the symbol and lands exactly on a bin (§1.3), the rectangular window is optimal — no windowing is applied deliberately.

5.2 Sync search (time offset + base frequency bin)

find_sync (decode.py) jointly searches coarse start time and base tone bin:

  • Frequency bins: bin_lo = int(1400/1.46484375) = 955 to bin_hi = int(1600/1.46484375) = 1092 (the standard 1400–1600 Hz WSPR audio window; a 1500 Hz sub-carrier is bin 1024).
  • Time offsets: t0 in range(0, 30000, 2048) — 15 steps of ¼-symbol (0 … 28672 samples, i.e. 0 … ~2.39 s), covering the ~1 s nominal start.
  • For each t0 it computes one power matrix, then scores every candidate bin0 with _sync_metric, keeping the best (t0, bin0, metric).

_sync_metric (decode.py) correlates the known sync vector against tone parities at base bin bin0 (bounds-checked: 0 ≤ bin0 and bin0+3 < 4096, else −1e30). With p = mat[:, bin0:bin0+4] the four tone powers per symbol:

good = Σ_n  (p[n,0]+p[n,2]  if SYNC[n]==0  else  p[n,1]+p[n,3])   # correct-parity tones
bad  = Σ_n  (p[n,1]+p[n,3]  if SYNC[n]==0  else  p[n,0]+p[n,2])   # wrong-parity tones
tot  = p.sum(axis=1) + 1e-12
metric = Σ_n (good_n − bad_n) / tot_n

Tones {0,2} have sync-bit 0; tones {1,3} have sync-bit 1. The metric rewards energy in the correct-parity tones normalised per symbol; for a clean on-bin signal it approaches NSYM (≈ 162). This value is reported back as the (misnamed) "snr" field.

5.3 Soft symbol demodulation

soft_symbols (decode.py), at the winning (t0, bin0), converts tone powers to a soft data value per symbol. Given the sync bit, only two of the four tones are candidates (data=0 vs data=1):

p0[n] = p[n,0] if SYNC[n]==0 else p[n,1]     # power of the "data = 0" tone
p1[n] = p[n,2] if SYNC[n]==0 else p[n,3]     # power of the "data = 1" tone
soft[n]      = p1[n] − p0[n]                  # >0 ⇒ data bit 1
hard_data[n] = 1 if soft[n] > 0 else 0

soft[n] is a signed confidence for the interleaved data bit (positive ⇒ 1).

5.4 De-interleave

_deinterleave (decode.py) applies the inverse of §4.2 to the 162 soft values, producing conv_llr in convolutional-code order (positive value ⇒ code bit 1).

5.5 FEC decode — stack (Zigangirov–Jelinek) sequential decoder

sequential_decode (decode.py) is a best-first stack sequential decoder for the K=32 rate-½ code — the practical choice because a K=32 Viterbi trellis has 2³¹ states and is intractable, whereas a sequential decoder explores only the promising paths.

Metric normalisation and the Fano bit metric. The soft inputs are scaled to unit mean magnitude (r = conv_llr / mean(|conv_llr|)), and a calibration gain alpha = 2.0 maps them to log-likelihoods. The per-code-bit Fano metric is:

bit_metric(val, e) = 0.5 + log2 σ(alpha · val · (+1 if e==1 else −1))
     where σ(x) = 1/(1+e^−x),   log2 σ(x) = −log2(1 + e^−x)
     (guarded: x>30 ⇒ 0;  x<−30 ⇒ x/ln2)

e is the code bit predicted by a candidate branch. A matching bit contributes ≈ +0.5; a contradicted bit contributes a large negative penalty. The +0.5 bias per code bit makes the correct path's cumulative metric drift upward while wrong paths drift down, so the best-first search locks onto the true path and terminates quickly.

The search (min-heap keyed on −metric, so highest-metric node pops first):

node = (metric, depth, reg, bits)   starting at (0, 0, 0, ())
while heap not empty and pops < 200000:
    pop the highest-metric node; pops += 1
    if depth == 81:  return bits[:50]              # success: 50 message bits
    choices = (0,)      if depth ≥ 50   else (0,1)  # tail bits are known zeros
    for b in choices:
        reg2 = ((reg << 1) | b) & 0xFFFFFFFF
        e0 = popcount(reg2 & 0xF2D05351) & 1
        e1 = popcount(reg2 & 0xE4613C47) & 1
        bm = bit_metric(r[2·depth], e0) + bit_metric(r[2·depth+1], e1)
        push (metric + bm, depth+1, reg2, bits+(b,))
return None                                        # budget exhausted → decode fail

Key details:

  • Depth 81 = 50 information + 31 tail bits; the decoder walks the full 81-bit trellis and returns the first 50 bits.
  • Known-zero tail: for depth ≥ 50 the only branch is b=0, halving the branching factor over the final third of the trellis (a strong constraint that helps the search converge).
  • Budget: max_pops = 200000. Exhaustion returns None, which surfaces as an empty decode result.

5.6 Unpack → message string

unpack (decode.py) reverses §3. It reads the first 28 bits as n, the next 22 as m (both MSB-first), then:

c6 = ch(n%27+10); n//=27     # A-Z/space positions (inverse of ·27, −10)
c5 = ch(n%27+10); n//=27
c4 = ch(n%27+10); n//=27
c3 = ch(n%10);    n//=10     # digit
c2 = ch(n%36);    n//=36
c1 = ch(n)
call = (c1..c6).strip()

pwr = m%128 − 64;  m1 = m//128
lo = m1 % 180;  l2, l4 = lo//10, lo%10
v  = 179 − m1//180;  l1, l3 = v//10, v%10
grid = chr('A'+l1) + chr('A'+l2) + str(l3) + str(l4)

ch(v): v<10 → digit, v<36 → letter, else space. Validity gates return None (reject the decode) unless 0≤l1<18, 0≤l2<18, call non-empty, and 0 ≤ pwr ≤ 60. Output string: "CALL GRID PWR".

5.7 Result record

decode_slot returns a one-element list (or [] on failure):

{ "snr":  "%.0f" % metric,              # raw sync-correlation metric (NOT a dB SNR)
  "dt":   "%.1f" % (t0 / 12000),        # detected start time, seconds
  "freq": "%.4f" % (bin0 · 1.46484375), # base (lowest) tone frequency, Hz
  "msg":  "CALL GRID PWR" }

For the clean M0SUP IO90 30 self-test the record is {'snr':'162','dt':'1.0','freq':'1500.0000','msg':'M0SUP IO90 30'} (t0=12000, bin0=1024).

Caveat — snr is not an SNR. The "snr" field carries the sync-correlation metric (≈162 for a clean signal), not a calibrated dB figure like wsprd reports. Treat it as a detection-confidence score.


6. Constants & tables (every value or exact source)

6.1 Scalar constants (constants.py)

CONV_POLY    = (0xF2D05351, 0xE4613C47)   # generator polynomials
CONV_K       = 32                          # constraint length
MSG_BITS     = 50                          # payload bits
TAIL_BITS    = 31  (= CONV_K − 1)          # zero flush bits
NSYM         = 162  (= (50+31)·2)          # channel symbols
SAMPLE_RATE  = 12000                       # Hz
SPS          = 8192                         # samples/symbol
TONE_SPACING = 1.46484375  (= 12000/8192)  # Hz = baud
SLOT_SECONDS = 120                          # slot length

6.2 The 162-element sync vector (constants.py, verbatim)

Symbol[n] = SYNC[n] + 2·Data[n]; asserted length 162 (constants.py).

1,1,0,0,0,0,0,0,1,0,0,0,1,1,1,0, 0,0,1,0,0,1,0,1,1,1,1,0,0,0,0,0,
0,0,1,0,0,1,0,1,0,0,0,0,0,0,1,0, 1,1,0,0,1,1,0,1,0,0,0,1,1,0,1,0,
0,0,0,1,1,0,1,0,1,0,1,0,1,0,0,1, 0,0,1,0,1,1,0,0,0,1,1,0,1,0,1,0,
0,0,1,0,0,0,0,0,1,0,0,1,0,0,1,1, 1,0,1,1,0,0,1,1,0,1,0,0,0,1,1,1,
0,0,0,0,0,1,0,1,0,0,1,1,0,0,0,0, 0,0,0,1,1,0,1,0,1,1,0,0,0,1,1,0,
0,0

This is the standard WSPR synchronisation sequence; wsprcode independently confirms it (its sync + 2·data symbols match ours — see docs/native-digimodes.md).

6.3 Interleaver permutation

Not a stored table — generated on the fly by bit-reversal. The exact recipe (§4.2): walk i = 0…255; j = bitrev8(i) (reverse the low 8 bits); take j < 162 in ascending-i order as destination positions for the coded bits p = 0…161. The forward map (coded index → symbol position), first 24 and last 8 entries:

p :  0   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15
j :  0 128  64  32 160  96  16 144  80  48 112   8 136  72  40 104
p : 16  17  18  19  20  21  22  23  …  154 155 156 157 158 159 160 161
j : 24 152  88  56 120   4 132  68  …   79  47 111  31 159  95  63 127

_deinterleave uses the identical walk to invert it.

6.4 Convolutional code

  • Constraint length K = 32, rate ½, non-systematic feed-forward, no puncturing.
  • Register: 32-bit, new bit enters at LSB, masked & 0xFFFFFFFF.
  • Generators: g0 = 0xF2D05351, g1 = 0xE4613C47 (bit-reversal-free tap masks; output = popcount(reg & g) & 1).
  • Tail: 31 zeros flush the register to a known terminal state; the decoder exploits this (§5.5).

6.5 Message-field radices (§3)

Callsign: 37·36·10·27·27·27 positional radix (mask 0x0FFFFFFF, 28 bits). Grid/power: m1·128 + dbm + 64 with m1 = 180·(179 − 10·A − N1) + (10·B + N2) (mask 0x3FFFFF, 22 bits).

6.6 Known-answer vectors (test_wspr.py)

pack_call("M0SUP")            = 259435617   (0x0F76AC61)
pack_grid_power("IO90", 30)   = 2091614     (0x001FEA5E)
_normalize_call("G4JNT")      = " G4JNT"
_normalize_call("W1AW")       = " W1AW "
channel_symbols("M0SUP","IO90",30) = M0SUP_SYMS  (162 symbols, test_wspr.py)
encode_pcm(...) length        = 120·12000·2 = 2 880 000 bytes

7. Interoperability & validation

7.1 What it is exact against

The reference is WSJT-X's own wsprcode / wsprd utilities, used as behavioural oracles (black-box input→output), which keeps the port clean-room (no GPL source copied) while giving exact reference vectors.

  • Encoder — symbol-exact. wsprcode prints the packed bits and sync + 2·data channel symbols for a message; ours match. The M0SUP_SYMS fixture (test_wspr.py) was captured by synthesising the signal with encode_pcm and decoding it with wsprd (WSJT-X 2.7.0) on prod, recovering "M0SUP IO90 30". Validated across standard, space-prepended, and short callsigns — reported 6/6.

7.2 Oracle methodology

  1. Pick a message CALL GRID dBm.
  2. Generate channel symbols with wsprcode; compare byte-for-byte to channel_symbols(...).
  3. Render encode_pcm(...) to a WAV and feed it to prod's wsprd; confirm it decodes back to the exact input string, and that the reported sync symbols equal ours.
  4. For the native decoder, round-trip encode_pcm → decode_slot (clean and noise-added) and, for sensitivity, compare native vs wsprd decode rates on noisy signals.

7.3 In-repo regression tests (test_wspr.py)

  • test_pack — the two known packing values above.
  • test_normalize_call — the four normalisation cases.
  • test_symbols — all 162 symbols equal M0SUP_SYMS, each in 0…3.
  • test_encode_pcm_slot — slot is exactly 2 880 000 bytes.
  • test_decode_roundtrip_cleanM0SUP/IO90/30, G4JNT/IO90/37, PA0XYZ/JO22/23 each encode→decode to the exact string.
  • test_decode_roundtrip_noisyM0SUP IO90 30 recovered under additive Gaussian noise (σ = 15 on a unit-tone signal, seed 7).

All pass (python3 test_wspr.py → "all WSPR encoder+decoder tests passed").

7.4 Sensitivity vs wsprd

The native decoder (FFT-per-symbol demod + stack sequential decoder) decodes cleanly to about −28 dB; wsprd reaches roughly 2–3 dB deeper thanks to its deep-search machinery (multi-seed at threshold: native 0/6 where wsprd 5/6). The native decoder is a correct, faithful baseline — not yet a sensitivity match.


8. Limitations and future work

Honest gaps, roughly in order of impact:

  • Decode sensitivity (~2–3 dB shy of wsprd). Avenues, in expected-gain order: (a) a calibrated LLR soft metric derived from the tone-power statistics rather than the raw power difference p1 − p0; (b) sub-bin / fractional-frequency fitting and drift compensation (matters on real off-air signals — the synthetic tests are exactly on-bin); (c) a properly tuned Fano threshold or a larger stack budget; (d) coherent multi-symbol processing and decoded-signal subtraction (WSJT-X's "deep search").
  • Type-1 messages only. Standard callsign + 4-char grid + power. No compound/ add-on callsign forms (Type-2/3, hashed calls, 6-char grids). _normalize_call rejects calls without a digit in position 3.
  • snr is a sync-correlation score, not a dB SNR (§5.7). No wsprd-style dB estimate or frequency-drift report.
  • Single-signal decode. decode_slot finds and returns the single strongest signal in the 1400–1600 Hz window; it does not iterate to pull out multiple overlapping WSPR signals in one slot (real bands carry dozens).
  • Plain FSK on TX — no GFSK pulse-shaping, so the transmitted spectrum is slightly wider than WSJT-X's (§4.4). Decodes fine; a cleanliness refinement.
  • Runtime cost. find_sync recomputes a 162×8192 FFT set for each of 15 time offsets, and the stack decoder is pure-Python heap work — CPU-heavy, though the 120 s slot gives ample headroom on the prod box.

None of these are architectural; each is a bounded addition on a working base.


9. Implementation & reproduction

file role
app/radio/wsprx/constants.py all scalar constants + the 162-element sync vector
app/radio/wsprx/encode.py packing, K=32 conv encode, interleave, sync merge, 4-FSK synth
app/radio/wsprx/decode.py FFT demod, sync search, soft symbols, stack FEC decode, unpack
app/radio/wsprx/codec.py public encode_pcm / decode_slot
app/radio/wsprx/__init__.py re-exports encode_pcm, decode_slot
test_wspr.py known-answer + round-trip regression tests

Runtime API:

from app.radio.wsprx import encode_pcm, decode_slot

pcm  = encode_pcm("M0SUP", "IO90", 30)   # 120 s s16le mono @ 12 kHz (bytes)
spots = decode_slot(pcm)                  # [{'snr','dt','freq','msg'}] or []
# -> [{'snr':'162','dt':'1.0','freq':'1500.0000','msg':'M0SUP IO90 30'}]

Everything at runtime is NumPy; there are no learned tables or downloaded weights — the sync vector and generator polynomials are transcribed protocol facts, exactly as the other native modes' constants were.


10. References

  1. J. Taylor, K1JT — WSPR, part of the WSJT-X project. Protocol overview and the 2-minute beacon format. (Behavioural oracle: wsprcode, wsprd.)
  2. A. Talbot, G4JNT — "WSPR Coding Process" — the public description of the 50-bit packing, K=32 rate-½ convolutional code, interleaver, and 162-symbol sync merge that this clean-room implementation follows.
  3. K. S. Zigangirov / F. Jelinek — the stack (ZJ) sequential decoding algorithm used for the K=32 code (§5.5).
  4. R. M. Fano — the sequential-decoding bit metric whose calibrated form drives the stack search.
  5. docs/native-digimodes.md — Rafe's native WSJT-X-family effort, validation discipline, and the *code/*sim oracle methodology.
  6. docs/rvqvoice.md — the depth/structure/tone benchmark for this document.