WSJT-X/lib/ftrsd/ftrsd2.c
Joe Taylor 0e15a4cd5c Default parameter settings for JT65 decoding set for HF conditions.
More aggressive parameters are selected with higher values of 
"Aggressive decoding level" on the "Advanced" tab.  With settings
greater than 0, be sure to set Ftol=1000 if you want to decode over 
the whole passband.


git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@6330 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
2015-12-31 16:22:27 +00:00

236 lines
6.5 KiB
C

/*
ftrsd2.c
A soft-decision decoder for the JT65 (63,12) Reed-Solomon code.
This decoding scheme is built around Phil Karn's Berlekamp-Massey
errors and erasures decoder. The approach is inspired by a number of
publications, including the stochastic Chase decoder described
in "Stochastic Chase Decoding of Reed-Solomon Codes", by Leroux et al.,
IEEE Communications Letters, Vol. 14, No. 9, September 2010 and
"Soft-Decision Decoding of Reed-Solomon Codes Using Successive Error-
and-Erasure Decoding," by Soo-Woong Lee and B. V. K. Vijaya Kumar.
Steve Franke K9AN and Joe Taylor K1JT
*/
#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>
#include <time.h>
#include <string.h>
#include "rs2.h"
static void *rs;
void getpp_(int workdat[], float *pp);
void ftrsd2_(int mrsym[], int mrprob[], int mr2sym[], int mr2prob[],
int* ntrials0, int correct[], int param[], int ntry[])
{
int rxdat[63], rxprob[63], rxdat2[63], rxprob2[63];
int workdat[63];
int indexes[63];
int era_pos[51];
int i, j, numera, nerr, nn=63;
FILE *logfile = NULL;
int ntrials = *ntrials0;
int verbose = 0;
int nhard=0,nhard_min=32768,nsoft=0,nsoft_min=32768;
int nsofter=0,nsofter_min=32768,ntotal=0,ntotal_min=32768,ncandidates;
int nera_best=0;
float pp,pp1,pp2;
static unsigned int nseed;
// Power-percentage symbol metrics - composite gnnf/hf
int perr[8][8] = {
{ 4, 9, 11, 13, 14, 14, 15, 15},
{ 2, 20, 20, 30, 40, 50, 50, 50},
{ 7, 24, 27, 40, 50, 50, 50, 50},
{13, 25, 35, 46, 52, 70, 50, 50},
{17, 30, 42, 54, 55, 64, 71, 70},
{25, 39, 48, 57, 64, 66, 77, 77},
{32, 45, 54, 63, 66, 75, 78, 83},
{51, 58, 57, 66, 72, 77, 82, 86}};
if(verbose) {
logfile=fopen("/tmp/ftrsd.log","a");
if( !logfile ) {
printf("Unable to open ftrsd.log\n");
exit(1);
}
}
// Initialize the KA9Q Reed-Solomon encoder/decoder
unsigned int symsize=6, gfpoly=0x43, fcr=3, prim=1, nroots=51;
rs=init_rs_int(symsize, gfpoly, fcr, prim, nroots, 0);
// Reverse the received symbol vectors for BM decoder
for (i=0; i<63; i++) {
rxdat[i]=mrsym[62-i];
rxprob[i]=mrprob[62-i];
rxdat2[i]=mr2sym[62-i];
rxprob2[i]=mr2prob[62-i];
}
// Sort rxprob to find indexes of the least reliable symbols
int k, pass, tmp, nsym=63;
int probs[63];
for (i=0; i<63; i++) {
indexes[i]=i;
probs[i]=rxprob[i];
}
for (pass = 1; pass <= nsym-1; pass++) {
for (k = 0; k < nsym - pass; k++) {
if( probs[k] < probs[k+1] ) {
tmp = probs[k];
probs[k] = probs[k+1];
probs[k+1] = tmp;
tmp = indexes[k];
indexes[k] = indexes[k+1];
indexes[k+1] = tmp;
}
}
}
// See if we can decode using BM HDD, and calculate the syndrome vector.
memset(era_pos,0,51*sizeof(int));
numera=0;
memcpy(workdat,rxdat,sizeof(rxdat));
nerr=decode_rs_int(rs,workdat,era_pos,numera,1);
if( nerr >= 0 ) {
// Hard-decision decoding succeeded. Save codeword and some parameters.
nhard=0;
for (i=0; i<63; i++) {
if( workdat[i] != rxdat[i] ) nhard=nhard+1;
}
memcpy(correct,workdat,63*sizeof(int));
param[0]=0;
param[1]=nhard;
param[2]=0;
param[3]=0;
param[4]=0;
param[5]=0;
param[7]=1000*1000;
ntry[0]=0;
return;
}
/*
Hard-decision decoding failed. Try the FT soft-decision method.
Generate random erasure-locator vectors and see if any of them
decode. This will generate a list of "candidate" codewords. The
soft distance between each candidate codeword and the received
word is estimated by finding the largest (pp1) and second-largest
(pp2) outputs from a synchronized filter-bank operating on the
symbol spectra, and using these to decide which candidate
codeword is "best".
*/
nseed=1; //Seed for random numbers
float ratio;
int thresh, nsum;
int thresh0[63];
ncandidates=0;
nsum=0;
int ii,jj;
for (i=0; i<nn; i++) {
nsum=nsum+rxprob[i];
j = indexes[62-i];
ratio = (float)rxprob2[j]/((float)rxprob[j]+0.01);
ii = 7.999*ratio;
jj = (62-i)/8;
thresh0[i] = 1.3*perr[ii][jj];
}
if(nsum<=0) return;
pp1=0.0;
pp2=0.0;
for (k=1; k<=ntrials; k++) {
memset(era_pos,0,51*sizeof(int));
memcpy(workdat,rxdat,sizeof(rxdat));
/*
Mark a subset of the symbols as erasures.
Run through the ranked symbols, starting with the worst, i=0.
NB: j is the symbol-vector index of the symbol with rank i.
*/
numera=0;
for (i=0; i<nn; i++) {
j = indexes[62-i];
thresh=thresh0[i];
long int ir;
// Generate a random number ir, 0 <= ir < 100 (see POSIX.1-2001 example).
nseed = nseed * 1103515245 + 12345;
ir = (unsigned)(nseed/65536) % 32768;
ir = (100*ir)/32768;
if((ir < thresh ) && numera < 51) {
era_pos[numera]=j;
numera=numera+1;
}
}
nerr=decode_rs_int(rs,workdat,era_pos,numera,0);
if( nerr >= 0 ) {
// We have a candidate coderowd. Find its hard and soft distance from
// the received word. Also find pp1 and pp2 from the full array
// s3(64,63) of synchronized symbol spectra.
ncandidates=ncandidates+1;
nhard=0;
nsoft=0;
nsofter=0;
for (i=0; i<63; i++) {
if(workdat[i] != rxdat[i]) {
nhard=nhard+1;
nsofter=nsofter+rxprob[i];
if(workdat[i] != rxdat2[i]) {
nsoft=nsoft+rxprob[i];
}
} else {
nsofter=nsofter-rxprob[i];
}
}
nsoft=63*nsoft/nsum;
nsofter=63*nsofter/nsum;
ntotal=nsoft+nhard;
getpp_(workdat,&pp);
if(pp>pp1) {
pp2=pp1;
pp1=pp;
nsoft_min=nsoft;
nhard_min=nhard;
nsofter_min=nsofter;
ntotal_min=ntotal;
memcpy(correct,workdat,63*sizeof(int));
nera_best=numera;
ntry[0]=k;
} else {
if(pp>pp2 && pp!=pp1) pp2=pp;
}
// if(ntotal_min <= 81 && pp2/pp1 <= 0.87) break;
if(nhard_min <= 41 && ntotal_min <= 71) break;
}
if(k == ntrials) ntry[0]=k;
}
if(logfile) {
fprintf(logfile,"ncand %4d nhard %4d nsoft %4d nhard+nsoft %4d nsum %8d\n",
ncandidates,nhard_min,nsoft_min,ntotal_min,nsum);
fclose(logfile);
}
param[0]=ncandidates;
param[1]=nhard_min;
param[2]=nsoft_min;
param[3]=nera_best;
param[4]=1000.0*pp2/pp1;
param[5]=ntotal_min;
param[6]=ntry[0];
param[7]=1000.0*pp2;
param[8]=1000.0*pp1;
if(param[0]==0) param[2]=-1;
return;
}