WSJT-X/lib/ldpc/decode.c
Steven Franke 5ac886855d Add ldpc sandbox folder.
git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@6437 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
2016-01-25 00:04:21 +00:00

385 lines
8.4 KiB
C
Executable File

/* DECODE.C - Decode blocks of received data. */
/* Copyright (c) 1995-2012 by Radford M. Neal.
*
* Permission is granted for anyone to copy, use, modify, and distribute
* these programs and accompanying documents for any purpose, provided
* this copyright notice is retained and prominently displayed, and note
* is made of any changes made to these programs. These programs and
* documents are distributed without any warranty, express or implied.
* As the programs were written for research purposes only, they have not
* been tested to the degree that would be advisable in any important
* application. All use of these programs is entirely at the user's own
* risk.
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "rand.h"
#include "alloc.h"
#include "blockio.h"
#include "open.h"
#include "mod2sparse.h"
#include "mod2dense.h"
#include "mod2convert.h"
#include "channel.h"
#include "rcode.h"
#include "check.h"
#include "dec.h"
void usage(void);
/* MAIN PROGRAM. */
int main
( int argc,
char **argv
)
{
char *pchk_file, *rfile, *dfile, *pfile;
char **meth;
FILE *rf, *df, *pf;
char *dblk, *pchk;
double *lratio;
double *bitpr;
double *awn_data; /* Places to store channel data */
int *bsc_data;
unsigned iters; /* Unsigned because can be huge for enum */
double tot_iter; /* Double because can be huge for enum */
double chngd, tot_changed; /* Double because can be fraction if lratio==1*/
int tot_valid;
char junk;
int valid;
int i, j, k;
/* Look at initial flag arguments. */
table = 0;
blockio_flush = 0;
while (argc>1)
{
if (strcmp(argv[1],"-t")==0)
{ if (table!=0) usage();
table = 1;
}
else if (strcmp(argv[1],"-T")==0)
{ if (table!=0) usage();
table = 2;
}
else if (strcmp(argv[1],"-f")==0)
{ if (blockio_flush!=0) usage();
blockio_flush = 1;
}
else
{ break;
}
argc -= 1;
argv += 1;
}
/* Look at arguments up to the decoding method specification. */
if (!(pchk_file = argv[1])
|| !(rfile = argv[2])
|| !(dfile = argv[3]))
{ usage();
}
if (argv[4]==0 || argv[5]==0) usage();
k = channel_parse(argv+4,argc-4);
if (k<=0)
{ pfile = argv[4];
k = channel_parse(argv+5,argc-5);
if (k<=0) usage();
meth = argv+5+k;
}
else
{ pfile = 0;
meth = argv+4+k;
}
/* Look at the specification of the decoding method, which starts at meth and
continues to the end of the command line (marked by a zero pointer). */
if (!meth[0]) usage();
if (strcmp(meth[0],"prprp")==0)
{ dec_method = Prprp;
if (!meth[1] || sscanf(meth[1],"%d%c",&max_iter,&junk)!=1 || meth[2])
{ usage();
}
}
else if (strcmp(meth[0],"enum-block")==0)
{ dec_method = Enum_block;
if (!(gen_file = meth[1]) || meth[2]) usage();
}
else if (strcmp(meth[0],"enum-bit")==0)
{ dec_method = Enum_bit;
if (!(gen_file = meth[1]) || meth[2]) usage();
}
else
{ usage();
}
/* Check that we aren't overusing standard input or output. */
if ((strcmp(pchk_file,"-")==0)
+ (strcmp(rfile,"-")==0) > 1)
{ fprintf(stderr,"Can't read more than one stream from standard input\n");
exit(1);
}
if ((table>0)
+ (strcmp(dfile,"-")==0)
+ (pfile!=0 && strcmp(pfile,"-")==0) > 1)
{ fprintf(stderr,"Can't send more than one stream to standard output\n");
exit(1);
}
/* Read parity check file. */
read_pchk(pchk_file);
if (N<=M)
{ fprintf(stderr,
"Number of bits (%d) should be greater than number of checks (%d)\n",N,M);
exit(1);
}
/* Open file of received data. */
rf = open_file_std(rfile,"r");
if (rf==NULL)
{ fprintf(stderr,"Can't open file of received data: %s\n",rfile);
exit(1);
}
/* Create file for decoded data. */
df = open_file_std(dfile,"w");
if (df==NULL)
{ fprintf(stderr,"Can't create file for decoded data: %s\n",dfile);
exit(1);
}
/* Create file for bit probabilities, if specified. */
if (pfile)
{ pf = open_file_std(pfile,"w");
if (pf==NULL)
{ fprintf(stderr,"Can't create file for bit probabilities: %s\n",pfile);
exit(1);
}
}
/* Allocate space for data from channel. */
switch (channel)
{ case BSC:
{ bsc_data = chk_alloc (N, sizeof *bsc_data);
break;
}
case AWGN: case AWLN:
{ awn_data = chk_alloc (N, sizeof *awn_data);
break;
}
default:
{ abort();
}
}
/* Allocate other space. */
dblk = chk_alloc (N, sizeof *dblk);
lratio = chk_alloc (N, sizeof *lratio);
pchk = chk_alloc (M, sizeof *pchk);
bitpr = chk_alloc (N, sizeof *bitpr);
/* Print header for summary table. */
if (table==1)
{ printf(" block iterations valid changed\n");
}
/* Do the setup for the decoding method. */
switch (dec_method)
{ case Prprp:
{ prprp_decode_setup();
break;
}
case Enum_block: case Enum_bit:
{ enum_decode_setup();
break;
}
default: abort();
}
/* Read received blocks, decode, and write decoded blocks. */
tot_iter = 0;
tot_valid = 0;
tot_changed = 0;
for (block_no = 0; ; block_no++)
{
/* Read block from received file, exit if end-of-file encountered. */
for (i = 0; i<N; i++)
{ int c;
switch (channel)
{ case BSC:
{ c = fscanf(rf,"%1d",&bsc_data[i]);
break;
}
case AWGN: case AWLN:
{ c = fscanf(rf,"%lf",&awn_data[i]);
break;
}
default: abort();
}
if (c==EOF)
{ if (i>0)
{ fprintf(stderr,
"Warning: Short block (%d long) at end of received file ignored\n",i);
}
goto done;
}
if (c<1 || channel==BSC && bsc_data[i]!=0 && bsc_data[i]!=1)
{ fprintf(stderr,"File of received data is garbled\n");
exit(1);
}
}
/* Find likelihood ratio for each bit. */
switch (channel)
{ case BSC:
{ for (i = 0; i<N; i++)
{ lratio[i] = bsc_data[i]==1 ? (1-error_prob) / error_prob
: error_prob / (1-error_prob);
}
break;
}
case AWGN:
{ for (i = 0; i<N; i++)
{ lratio[i] = exp(2*awn_data[i]/(std_dev*std_dev));
}
break;
}
case AWLN:
{ for (i = 0; i<N; i++)
{ double e, d1, d0;
e = exp(-(awn_data[i]-1)/lwidth);
d1 = 1 / ((1+e)*(1+1/e));
e = exp(-(awn_data[i]+1)/lwidth);
d0 = 1 / ((1+e)*(1+1/e));
lratio[i] = d1/d0;
}
break;
}
default: abort();
}
/* Try to decode using the specified method. */
switch (dec_method)
{ case Prprp:
{ iters = prprp_decode (H, lratio, dblk, pchk, bitpr);
break;
}
case Enum_block: case Enum_bit:
{ iters = enum_decode (lratio, dblk, bitpr, dec_method==Enum_block);
break;
}
default: abort();
}
/* See if it worked, and how many bits were changed. */
valid = check(H,dblk,pchk)==0;
chngd = changed(lratio,dblk,N);
tot_iter += iters;
tot_valid += valid;
tot_changed += chngd;
/* Print summary table entry. */
if (table==1)
{ printf ("%7d %10f %d %8.1f\n",
block_no, (double)iters, valid, (double)chngd);
/* iters is printed as a double to avoid problems if it's >= 2^31 */
fflush(stdout);
}
/* Write decoded block. */
blockio_write(df,dblk,N);
/* Write bit probabilities, if asked to. */
if (pfile)
{ for (j = 0; j<N; j++)
{ fprintf(pf," %.5f",bitpr[j]);
}
fprintf(pf,"\n");
}
/* Check for errors when writing. */
if (ferror(df) || pfile && ferror(pf))
{ break;
}
}
/* Finish up. */
done:
fprintf(stderr,
"Decoded %d blocks, %d valid. Average %.1f iterations, %.0f%% bit changes\n",
block_no, tot_valid, (double)tot_iter/block_no,
100.0*(double)tot_changed/(N*block_no));
if (ferror(df) || fclose(df)!=0)
{ fprintf(stderr,"Error writing decoded blocks to %s\n",dfile);
exit(1);
}
if (pfile)
{ if (ferror(pf) || fclose(pf)!=0)
{ fprintf(stderr,"Error writing bit probabilities to %s\n",dfile);
exit(1);
}
}
exit(0);
}
/* PRINT USAGE MESSAGE AND EXIT. */
void usage(void)
{ fprintf(stderr,"Usage:\n");
fprintf(stderr,
" decode [ -f ] [ -t | -T ] pchk-file received-file decoded-file [ bp-file ] channel method\n");
channel_usage();
fprintf(stderr,
"Method: enum-block gen-file | enum-bit gen-file | prprp [-]max-iterations\n");
exit(1);
}