ewsdr/BeaconFEC.pas

unit BeaconFEC;

{
  QO-100 центральный маяк — FEC-бэкенд (Milestone 2).

  Декодирует AO-40 FEC поверх потока мягких символов (400/с, после Manchester+
  дифф-декода фронтенда). Цепочка (ровно как gr-satellites ao40_fec_deframer):

    мягкие символы 400/с
      → distributed sync (65 бит, шаг 80, окно 5200, порог 8, обе полярности)
      → matrix deinterleave 80×65: out[i]=data[80*(i%65)+i/65], skip 65, take 5132
      → Viterbi r=1/2 k=7, полиномы [0x4F, -0x6D] (CCSDS), terminated → 2560 бит/320 байт
      → CCSDS additive descramble (LFSR mask 0xA9, seed 0xFF, len 7, XOR по битам)
      → 2× RS(255,223) укороч. pad=95 → (160,128), conventional basis, 2-way интерлив
      → 256 байт кадра → OnFrame.

  Viterbi и RS — через libfec Карна (разделяемая libfec.so/.dll/.dylib, load-time
  линковка как WDSP/IIO). libfec ставится отдельно на каждой ОС.
}

{$IFDEF FPC}{$MODE Delphi}{$ENDIF}
{$PACKRECORDS C}

interface

uses
  Classes, SysUtils;

const
  AO40_SYNC_LEN    = 65;
  AO40_SYNC_STEP   = 80;
  AO40_FRAME_SYMS  = AO40_SYNC_LEN * AO40_SYNC_STEP;   // 5200
  AO40_IL_ROWS     = 80;
  AO40_IL_COLS     = 65;
  AO40_VITERBI_SYM = 5132;                             // после skip 65
  AO40_DATA_BITS   = 2560;                             // выход Viterbi (без 6 tail)
  AO40_CONV_BYTES  = AO40_DATA_BITS div 8;             // 320
  AO40_RS_NROOTS   = 32;
  AO40_RS_PAD      = 95;                               // 255-160
  AO40_RS_NN       = 160;                              // укороченная длина
  AO40_RS_KK       = 128;
  AO40_INTERLEAVE  = 2;
  AO40_FRAME_BYTES = AO40_RS_KK * AO40_INTERLEAVE;     // 256
  AO40_SYNC_THRESH = 8;                                // допуск ошибок синхры

type
  TBeaconFrame = array[0..AO40_FRAME_BYTES-1] of Byte;
  TBeaconFrameEvent = procedure(const Frame: TBeaconFrame; RSErrors: Integer) of object;

  TBeaconFEC = class
  private
    FSync:     array[0..AO40_SYNC_LEN-1] of Byte;
    FRing:     array[0..AO40_FRAME_SYMS-1] of Single;
    FRingPos:  Integer;                 // позиция следующей записи = старейший элемент
    FRingCnt:  Integer;
    FCooldown: Integer;                 // не перепроверять сразу после успешного кадра
    FVp:       Pointer;                 // viterbi27 instance
    FFrames:   Int64;
    FOnFrame:  TBeaconFrameEvent;

    function  CheckSyncAt: Integer;     // -1 нет; 0 норм. полярность; 1 инверсия
    procedure DecodeCurrentFrame(Invert: Boolean);
    procedure Descramble(var Buf: array of Byte; Count: Integer);
  public
    constructor Create;
    destructor Destroy; override;
    procedure Reset;
    procedure PushSoftSymbol(const S: Single);
    function  SelfTest(out Msg: string): Boolean;
    property  Frames: Int64 read FFrames;
    property  OnFrame: TBeaconFrameEvent read FOnFrame write FOnFrame;
  end;

implementation

// ---------------- libfec bindings (Karn) ----------------
// Динамическая (load-time) линковка против разделяемой libfec — как WDSP.
// libfec ставится отдельно на каждой ОС (Linux: make install из исходников;
// Windows/macOS — отдельная сборка).
//
// {$LINKLIB fec} обязателен: одни лишь `external FEC_LIB` на функциях НЕ передают
// библиотеку линкеру (ld) на FPC/Darwin — символы остаются undefined. Точно так
// же работает WDSP ({$linklib wdsp}). Это НЕ форсирует статику: при наличии и
// libfec.dylib, и libfec.a в пути поиска ld предпочитает .dylib (как для wdsp).
//
// Только Darwin: на Windows {$LINKLIB fec} заставляет ld искать импорт-либу
// (libfec.a/.dll.a) и падает с "Import library not found for fec" — там load-time
// линковка делается самим FPC из external 'libfec.dll'. На Linux аналогично (как WDSP).
{$IFDEF DARWIN}{$LINKLIB fec}{$ENDIF}
const
{$IFDEF WINDOWS}
  FEC_LIB = 'libfec.dll';
{$ENDIF}
{$IFDEF LINUX}
  FEC_LIB = 'libfec.so';
{$ENDIF}
{$IFDEF DARWIN}
  FEC_LIB = 'libfec.dylib';
{$ENDIF}

{$IFDEF UNIX}{$LINKLIB m}{$ENDIF}   // libfec использует log() из libm

function  create_viterbi27(len: LongInt): Pointer; cdecl; external FEC_LIB;
function  init_viterbi27(vp: Pointer; starting_state: LongInt): LongInt; cdecl; external FEC_LIB;
function  update_viterbi27_blk(vp: Pointer; syms: PByte; npairs: LongInt): LongInt; cdecl; external FEC_LIB;
function  chainback_viterbi27(vp: Pointer; data: PByte; nbits: LongWord; endstate: LongWord): LongInt; cdecl; external FEC_LIB;
procedure delete_viterbi27(vp: Pointer); cdecl; external FEC_LIB;

function  decode_rs_8(data: PByte; eras_pos: PLongInt; no_eras: LongInt; pad: LongInt): LongInt; cdecl; external FEC_LIB;
procedure encode_rs_8(data: PByte; parity: PByte; pad: LongInt); cdecl; external FEC_LIB;

// ---------------- helpers ----------------

function ParityB(x: LongWord): Byte; inline;
// чётность числа единичных битов: 1=нечётно, 0=чётно
begin
  x := x xor (x shr 16);
  x := x xor (x shr 8);
  x := x xor (x shr 4);
  x := x xor (x shr 2);
  x := x xor (x shr 1);
  Result := x and 1;
end;

const
  // 65-битная распределённая синхра AO-40 (gr-satellites _syncword)
  AO40_SYNCWORD: array[0..AO40_SYNC_LEN-1] of Byte = (
    1,1,1,1,1,1,1,0,0,0,0,1,1,1,0,1,1,1,1,0,0,1,0,1,1,0,0,1,0,0,1,0,
    0,0,0,0,0,1,0,0,0,1,0,0,1,1,0,0,0,1,0,1,1,1,0,1,0,1,1,0,1,1,0,0,0);

// ---------------- TBeaconFEC ----------------

constructor TBeaconFEC.Create;
var i: Integer;
begin
  inherited Create;
  for i := 0 to AO40_SYNC_LEN-1 do FSync[i] := AO40_SYNCWORD[i];
  FVp := create_viterbi27(AO40_DATA_BITS);
  Reset;
end;

destructor TBeaconFEC.Destroy;
begin
  if FVp <> nil then delete_viterbi27(FVp);
  inherited Destroy;
end;

procedure TBeaconFEC.Reset;
var i: Integer;
begin
  for i := 0 to AO40_FRAME_SYMS-1 do FRing[i] := 0;
  FRingPos := 0; FRingCnt := 0; FCooldown := 0;
end;

procedure TBeaconFEC.PushSoftSymbol(const S: Single);
var r: Integer;
begin
  FRing[FRingPos] := S;
  FRingPos := FRingPos + 1;
  if FRingPos >= AO40_FRAME_SYMS then FRingPos := 0;
  if FRingCnt < AO40_FRAME_SYMS then Inc(FRingCnt);

  if FCooldown > 0 then begin Dec(FCooldown); Exit; end;
  if FRingCnt < AO40_FRAME_SYMS then Exit;

  r := CheckSyncAt;
  if r >= 0 then
  begin
    DecodeCurrentFrame(r = 1);
    FCooldown := AO40_FRAME_SYMS;   // следующий кадр — через полный фрейм
  end;
end;

function TBeaconFEC.CheckSyncAt: Integer;
// Кандидат-кадр = текущее содержимое кольца (старейший→новейший), старт = FRingPos.
// Синхра на относительных позициях j*step. Возвращает полярность или -1.
var
  j, idx, hard, matches: Integer;
begin
  matches := 0;
  for j := 0 to AO40_SYNC_LEN-1 do
  begin
    idx := FRingPos + j * AO40_SYNC_STEP;
    if idx >= AO40_FRAME_SYMS then Dec(idx, AO40_FRAME_SYMS);
    if FRing[idx] > 0.0 then hard := 1 else hard := 0;   // gr-конвенция: +soft = «1»
    if hard = FSync[j] then Inc(matches);
  end;
  if matches >= AO40_SYNC_LEN - AO40_SYNC_THRESH then Exit(0);   // норм. полярность
  if matches <= AO40_SYNC_THRESH then Exit(1);                   // инверсия
  Result := -1;
end;

procedure TBeaconFEC.Descramble(var Buf: array of Byte; Count: Integer);
// CCSDS additive descrambler = GNU Radio additive_scrambler_bb(0xA9, 0xFF, 7).
// Fibonacci LFSR, seed 0xFF, XOR с каждым битом (MSB-first), reset на старте кадра.
var
  state: Byte;
  i, b, outp, newbit, bit: Integer;
  v: Byte;
begin
  state := $FF;
  for i := 0 to Count-1 do
  begin
    v := Buf[i];
    for b := 7 downto 0 do
    begin
      outp := state and 1;
      newbit := ParityB(state and $A9);
      state := (state shr 1) or (newbit shl 7);
      bit := (v shr b) and 1;
      bit := bit xor outp;
      if bit <> 0 then v := v or (1 shl b) else v := v and not (1 shl b);
    end;
    Buf[i] := v;
  end;
end;

procedure TBeaconFEC.DecodeCurrentFrame(Invert: Boolean);
var
  win:    array[0..AO40_FRAME_SYMS-1] of Single;
  deint:  array[0..AO40_FRAME_SYMS-1] of Single;
  syms:   array[0..AO40_VITERBI_SYM-1] of Byte;
  conv:   array[0..AO40_CONV_BYTES-1] of Byte;
  cw:     array[0..AO40_RS_NN-1] of Byte;
  frame:  TBeaconFrame;
  i, idx, k, j, sgn, totErr, res: Integer;
  v: Double;
begin
  // 1. снимок кольца в хронологическом порядке (+ инверсия полярности)
  if Invert then sgn := -1 else sgn := 1;
  for i := 0 to AO40_FRAME_SYMS-1 do
  begin
    idx := FRingPos + i;
    if idx >= AO40_FRAME_SYMS then Dec(idx, AO40_FRAME_SYMS);
    win[i] := sgn * FRing[idx];
  end;

  // 2. матричный деинтерливер 80×65, out[i]=data[80*(i%65)+i/65]
  for i := 0 to AO40_FRAME_SYMS-1 do
    deint[i] := win[AO40_IL_ROWS * (i mod AO40_IL_COLS) + (i div AO40_IL_COLS)];

  // 3. skip 65 (синхра) → 5132 мягких символа → байты для libfec (255=«1»).
  // РЕМАП КОНВЕНЦИИ Viterbi: AO-40/gr-satellites кодирует полиномами [0x4F, -0x6D]
  // (cc_decoder), а libfec viterbi27 ЗАШИВАЕТ [0x6D, 0x4F] БЕЗ инверсии. Эмпирически
  // (брут-форс против gr post_viterbi_reference, см. историю) точное соответствие:
  //   libfec[2p]   = -soft[2p+1]   (своп пары + инверсия первого символа)
  //   libfec[2p+1] =  soft[2p]
  // Без ремапа Viterbi выдаёт мусор → RS не корректируется.
  for i := 0 to (AO40_VITERBI_SYM div 2) - 1 do
  begin
    v := 128.0 - deint[AO40_SYNC_LEN + 2*i + 1] * 100.0;     // -soft[2p+1]
    if v > 255.0 then v := 255.0 else if v < 0.0 then v := 0.0;
    syms[2*i]   := Round(v);
    v := 128.0 + deint[AO40_SYNC_LEN + 2*i] * 100.0;         //  soft[2p]
    if v > 255.0 then v := 255.0 else if v < 0.0 then v := 0.0;
    syms[2*i+1] := Round(v);
  end;

  // 4. Viterbi r=1/2 k=7 terminated → 2560 бит = 320 байт (MSB-first)
  init_viterbi27(FVp, 0);
  update_viterbi27_blk(FVp, @syms[0], (AO40_DATA_BITS + 6));   // 2566 пар
  chainback_viterbi27(FVp, @conv[0], AO40_DATA_BITS, 0);

  // 5. CCSDS дескремблер
  Descramble(conv, AO40_CONV_BYTES);

  // 6. RS: 2-way байт-интерлив, decode_rs_8 укороч. pad=95
  totErr := 0;
  for j := 0 to AO40_INTERLEAVE-1 do
  begin
    for k := 0 to AO40_RS_NN-1 do
      cw[k] := conv[j + k * AO40_INTERLEAVE];
    res := decode_rs_8(@cw[0], nil, 0, AO40_RS_PAD);
    if res < 0 then Exit;                   // некорректируемо → кадр невалиден
    Inc(totErr, res);
    for k := 0 to AO40_RS_KK-1 do
      frame[j + k * AO40_INTERLEAVE] := cw[k];
  end;

  Inc(FFrames);
  if Assigned(FOnFrame) then FOnFrame(frame, totErr);
end;

// ---------------- self-test ----------------

function TBeaconFEC.SelfTest(out Msg: string): Boolean;
// Проверяет линковку libfec и базовую корректность: RS encode/decode round-trip
// с внесёнными ошибками (укорочение pad=95) + Viterbi encode/decode round-trip.
var
  data: array[0..AO40_RS_KK-1] of Byte;
  par:  array[0..AO40_RS_NROOTS-1] of Byte;
  cw:   array[0..AO40_RS_NN-1] of Byte;
  i, res, encst, sym0, sym1, errs: Integer;
  bits: array[0..AO40_DATA_BITS-1] of Byte;
  vsyms: array[0..AO40_VITERBI_SYM-1] of Byte;
  outb: array[0..AO40_CONV_BYTES-1] of Byte;
  b: Integer;
begin
  Result := False;

  // --- RS round-trip ---
  for i := 0 to AO40_RS_KK-1 do data[i] := (i * 37 + 11) and $FF;
  encode_rs_8(@data[0], @par[0], AO40_RS_PAD);
  for i := 0 to AO40_RS_KK-1 do cw[i] := data[i];
  for i := 0 to AO40_RS_NROOTS-1 do cw[AO40_RS_KK + i] := par[i];
  // вносим 16 байтовых ошибок (= предел исправления для nroots=32)
  for i := 0 to 15 do cw[i * 9] := cw[i * 9] xor $A5;
  res := decode_rs_8(@cw[0], nil, 0, AO40_RS_PAD);
  if res < 0 then begin Msg := 'RS decode FAILED (link/params?)'; Exit; end;
  for i := 0 to AO40_RS_KK-1 do
    if cw[i] <> data[i] then begin Msg := Format('RS mismatch at %d', [i]); Exit; end;

  // --- Viterbi round-trip ---
  // libfec viterbi27 зашивает полиномы V27POLYA=$6D, V27POLYB=$4F (CCSDS) БЕЗ
  // инверсии и без API смены. Энкодер строго как у Карна (vtest27.c):
  //   sr=(sr shl 1)|bit; sym0=parity(sr & $6D); sym1=parity(sr & $4F).
  // Мягкий «1» = 255, «0» = 0. + 6 хвостовых нулей (terminated).
  Randomize;
  for i := 0 to AO40_DATA_BITS-1 do bits[i] := Random(2);
  encst := 0;
  for i := 0 to AO40_DATA_BITS + 6 - 1 do
  begin
    if i < AO40_DATA_BITS then b := bits[i] else b := 0;
    encst := ((encst shl 1) or b) and $7F;
    sym0 := ParityB(encst and $6D);
    sym1 := ParityB(encst and $4F);
    vsyms[2*i]   := sym0 * 255;
    vsyms[2*i+1] := sym1 * 255;
  end;
  init_viterbi27(FVp, 0);
  update_viterbi27_blk(FVp, @vsyms[0], AO40_DATA_BITS + 6);
  chainback_viterbi27(FVp, @outb[0], AO40_DATA_BITS, 0);
  errs := 0;
  for i := 0 to AO40_DATA_BITS-1 do
  begin
    b := (outb[i div 8] shr (7 - (i mod 8))) and 1;   // MSB-first
    if b <> bits[i] then Inc(errs);
  end;
  if errs <> 0 then
  begin
    Msg := Format('Viterbi round-trip: %d/%d bit errors (convention mismatch)',
                  [errs, AO40_DATA_BITS]);
    Exit;
  end;

  Msg := Format('OK: RS corrected 16 errs; Viterbi 0/%d bit errors', [AO40_DATA_BITS]);
  Result := True;
end;

end.