git.gag.com Git - debian/gnuradio/blob - gnuradio-core/src/lib/reed-solomon/init_rs.c

   1 /* Initialize a RS codec
   2  *
   3  * Copyright 2002 Phil Karn, KA9Q
   4  * May be used under the terms of the GNU General Public License (GPL)
   5  */
   6 #include <stdlib.h>
   7
   8 #ifdef CCSDS
   9 #include "ccsds.h"
  10 #elif defined(BIGSYM)
  11 #include "int.h"
  12 #else
  13 #include "char.h"
  14 #endif
  15
  16 #define NULL ((void *)0)
  17
  18 void FREE_RS(void *p){
  19   struct rs *rs = (struct rs *)p;
  20
  21   free(rs->alpha_to);
  22   free(rs->index_of);
  23   free(rs->genpoly);
  24   free(rs);
  25 }
  26
  27 /* Initialize a Reed-Solomon codec
  28  * symsize = symbol size, bits (1-8)
  29  * gfpoly = Field generator polynomial coefficients
  30  * fcr = first root of RS code generator polynomial, index form
  31  * prim = primitive element to generate polynomial roots
  32  * nroots = RS code generator polynomial degree (number of roots)
  33  */
  34 void *INIT_RS(unsigned int symsize,unsigned int gfpoly,unsigned fcr,unsigned prim,
  35                 unsigned int nroots){
  36   struct rs *rs;
  37   int i, j, sr,root,iprim;
  38
  39   if(symsize > 8*sizeof(DTYPE))
  40     return NULL; /* Need version with ints rather than chars */
  41
  42   if(fcr >= (1<<symsize))
  43     return NULL;
  44   if(prim == 0 || prim >= (1<<symsize))
  45     return NULL;
  46   if(nroots >= (1<<symsize))
  47     return NULL; /* Can't have more roots than symbol values! */
  48
  49   rs = (struct rs *)calloc(1,sizeof(struct rs));
  50   rs->mm = symsize;
  51   rs->nn = (1<<symsize)-1;
  52
  53   rs->alpha_to = (DTYPE *)malloc(sizeof(DTYPE)*(rs->nn+1));
  54   if(rs->alpha_to == NULL){
  55     free(rs);
  56     return NULL;
  57   }
  58   rs->index_of = (DTYPE *)malloc(sizeof(DTYPE)*(rs->nn+1));
  59   if(rs->index_of == NULL){
  60     free(rs->alpha_to);
  61     free(rs);
  62     return NULL;
  63   }
  64
  65   /* Generate Galois field lookup tables */
  66   rs->index_of[0] = A0; /* log(zero) = -inf */
  67   rs->alpha_to[A0] = 0; /* alpha**-inf = 0 */
  68   sr = 1;
  69   for(i=0;i<rs->nn;i++){
  70     rs->index_of[sr] = i;
  71     rs->alpha_to[i] = sr;
  72     sr <<= 1;
  73     if(sr & (1<<symsize))
  74       sr ^= gfpoly;
  75     sr &= rs->nn;
  76   }
  77   if(sr != 1){
  78     /* field generator polynomial is not primitive! */
  79     free(rs->alpha_to);
  80     free(rs->index_of);
  81     free(rs);
  82     return NULL;
  83   }
  84
  85   /* Form RS code generator polynomial from its roots */
  86   rs->genpoly = (DTYPE *)malloc(sizeof(DTYPE)*(nroots+1));
  87   if(rs->genpoly == NULL){
  88     free(rs->alpha_to);
  89     free(rs->index_of);
  90     free(rs);
  91     return NULL;
  92   }
  93   rs->fcr = fcr;
  94   rs->prim = prim;
  95   rs->nroots = nroots;
  96
  97   /* Find prim-th root of 1, used in decoding */
  98   for(iprim=1;(iprim % prim) != 0;iprim += rs->nn)
  99     ;
 100   rs->iprim = iprim / prim;
 101
 102   rs->genpoly[0] = 1;
 103   for (i = 0,root=fcr*prim; i < nroots; i++,root += prim) {
 104     rs->genpoly[i+1] = 1;
 105
 106     /* Multiply rs->genpoly[] by  @**(root + x) */
 107     for (j = i; j > 0; j--){
 108       if (rs->genpoly[j] != 0)
 109         rs->genpoly[j] = rs->genpoly[j-1] ^ rs->alpha_to[modnn(rs,rs->index_of[rs->genpoly[j]] + root)];
 110       else
 111         rs->genpoly[j] = rs->genpoly[j-1];
 112     }
 113     /* rs->genpoly[0] can never be zero */
 114     rs->genpoly[0] = rs->alpha_to[modnn(rs,rs->index_of[rs->genpoly[0]] + root)];
 115   }
 116   /* convert rs->genpoly[] to index form for quicker encoding */
 117   for (i = 0; i <= nroots; i++)
 118     rs->genpoly[i] = rs->index_of[rs->genpoly[i]];
 119
 120 #if 0
 121   printf ("genpoly:\n");
 122   for (i = nroots; i >= 0; i--){
 123     printf ("  %3d*X^%d\n", rs->alpha_to[rs->genpoly[i]], i);
 124   }
 125 #endif
 126
 127   return rs;
 128 }