git.gag.com Git - debian/gnuradio/blob - gnuradio-core/src/lib/viterbi/metrics.c

   1 /*
   2  * Copyright 1995 Phil Karn, KA9Q
   3  * Copyright 2008 Free Software Foundation, Inc.
   4  *
   5  * This file is part of GNU Radio
   6  *
   7  * GNU Radio is free software; you can redistribute it and/or modify
   8  * it under the terms of the GNU General Public License as published by
   9  * the Free Software Foundation; either version 3, or (at your option)
  10  * any later version.
  11  *
  12  * GNU Radio is distributed in the hope that it will be useful,
  13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15  * GNU General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU General Public License
  18  * along with GNU Radio; see the file COPYING.  If not, write to
  19  * the Free Software Foundation, Inc., 51 Franklin Street,
  20  * Boston, MA 02110-1301, USA.
  21  */
  22
  23 /*
  24  * Generate metric tables for a soft-decision convolutional decoder
  25  * assuming gaussian noise on a PSK channel.
  26  *
  27  * Works from "first principles" by evaluating the normal probability
  28  * function and then computing the log-likelihood function
  29  * for every possible received symbol value
  30  *
  31  */
  32
  33 /* Symbols are offset-binary, with 128 corresponding to an erased (no
  34  * information) symbol
  35  */
  36 #define OFFSET  128
  37
  38 #include <stdlib.h>
  39 #include <math.h>
  40
  41 /* Normal function integrated from -Inf to x. Range: 0-1 */
  42 #define normal(x)       (0.5 + 0.5*erf((x)/M_SQRT2))
  43
  44 /* Logarithm base 2 */
  45 #define log2(x) (log(x)*M_LOG2E)
  46
  47 /* Generate log-likelihood metrics for 8-bit soft quantized channel
  48  * assuming AWGN and BPSK
  49  */
  50 void
  51 gen_met(int mettab[2][256],     /* Metric table, [sent sym][rx symbol] */
  52         int amp,                /* Signal amplitude, units */
  53         double esn0,            /* Es/N0 ratio in dB */
  54         double bias,            /* Metric bias; 0 for viterbi, rate for sequential */
  55         int scale)              /* Scale factor */
  56 {
  57   double noise;
  58   int s,bit;
  59   double metrics[2][256];
  60   double p0,p1;
  61
  62   /* Es/N0 as power ratio */
  63   esn0 = pow(10.,esn0/10);
  64
  65   noise = 0.5/esn0;     /* only half the noise for BPSK */
  66   noise = sqrt(noise);  /* noise/signal Voltage ratio */
  67
  68   /* Zero is a special value, since this sample includes all
  69    * lower samples that were clipped to this value, i.e., it
  70    * takes the whole lower tail of the curve
  71    */
  72   p1 = normal(((0-OFFSET+0.5)/amp - 1)/noise);  /* P(s|1) */
  73
  74   /* Prob of this value occurring for a 0-bit */        /* P(s|0) */
  75   p0 = normal(((0-OFFSET+0.5)/amp + 1)/noise);
  76   metrics[0][0] = log2(2*p0/(p1+p0)) - bias;
  77   metrics[1][0] = log2(2*p1/(p1+p0)) - bias;
  78
  79   for(s=1;s<255;s++){
  80     /* P(s|1), prob of receiving s given 1 transmitted */
  81     p1 = normal(((s-OFFSET+0.5)/amp - 1)/noise) -
  82       normal(((s-OFFSET-0.5)/amp - 1)/noise);
  83
  84     /* P(s|0), prob of receiving s given 0 transmitted */
  85     p0 = normal(((s-OFFSET+0.5)/amp + 1)/noise) -
  86       normal(((s-OFFSET-0.5)/amp + 1)/noise);
  87
  88 #ifdef notdef
  89     printf("P(%d|1) = %lg, P(%d|0) = %lg\n",s,p1,s,p0);
  90 #endif
  91     metrics[0][s] = log2(2*p0/(p1+p0)) - bias;
  92     metrics[1][s] = log2(2*p1/(p1+p0)) - bias;
  93   }
  94   /* 255 is also a special value */
  95   /* P(s|1) */
  96   p1 = 1 - normal(((255-OFFSET-0.5)/amp - 1)/noise);
  97   /* P(s|0) */
  98   p0 = 1 - normal(((255-OFFSET-0.5)/amp + 1)/noise);
  99
 100   metrics[0][255] = log2(2*p0/(p1+p0)) - bias;
 101   metrics[1][255] = log2(2*p1/(p1+p0)) - bias;
 102 #ifdef  notdef
 103   /* The probability of a raw symbol error is the probability
 104    * that a 1-bit would be received as a sample with value
 105    * 0-128. This is the offset normal curve integrated from -Inf to 0.
 106    */
 107   printf("symbol Pe = %lg\n",normal(-1/noise));
 108 #endif
 109   for(bit=0;bit<2;bit++){
 110     for(s=0;s<256;s++){
 111       /* Scale and round to nearest integer */
 112       mettab[bit][s] = floor(metrics[bit][s] * scale + 0.5);
 113 #ifdef  notdef
 114       printf("metrics[%d][%d] = %lg, mettab = %d\n",
 115              bit,s,metrics[bit][s],mettab[bit][s]);
 116 #endif
 117     }
 118   }
 119 }