diff options
Diffstat (limited to 'ARToolKitPlus/src/extra/BCH.cxx')
| -rw-r--r-- | ARToolKitPlus/src/extra/BCH.cxx | 648 |
1 files changed, 0 insertions, 648 deletions
diff --git a/ARToolKitPlus/src/extra/BCH.cxx b/ARToolKitPlus/src/extra/BCH.cxx deleted file mode 100644 index 294eab1..0000000 --- a/ARToolKitPlus/src/extra/BCH.cxx +++ /dev/null @@ -1,648 +0,0 @@ -/* ========================================================================
- * PROJECT: ARToolKitPlus
- * ========================================================================
- *
- * See src/extra/BCH_original.txt for details/orig. copyright notice
- *
- * Copyright of the derived and new portions of this work
- * (C) 2006 Graz University of Technology
- *
- * This framework is free software; you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation; either version 2 of the License, or
- * (at your option) any later version.
- *
- * This framework is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this framework; if not, write to the Free Software
- * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
- *
- * For further information please contact
- * Dieter Schmalstieg
- * <schmalstieg@icg.tu-graz.ac.at>
- * Graz University of Technology,
- * Institut for Computer Graphics and Vision,
- * Inffeldgasse 16a, 8010 Graz, Austria.
- * ========================================================================
- ** @author Thomas Pintaric
- *
- * $Id: BCH.cxx 162 2006-04-19 21:28:10Z grabner $
- * @file
- * ======================================================================== */
-
-
-#include <ARToolKitPlus/extra/BCH.h>
-#include "assert.h"
-
-#include <vector>
-
-namespace ARToolKitPlus {
-
-
-static bool
-_isBitSet(_64bits n, int which_bit)
-{
- return ((n>>which_bit)&1) != 0;
-}
-
-
-static void
-_setBit(_64bits &n, int which_bit)
-{
- const _64bits one = 1;
- n |= (one<<which_bit);
-}
-
-/*
-static void
-_clearBit(_64bits &n, int which_bit)
-{
- const _64bits one = 1;
- n &= ~(one<<which_bit);
-}
-*/
-
-/*
-static void
-_copyBit(_64bits &dest_n, const int dest_bit, const _64bits src_n, const int src_bit)
-{
- const _64bits one = 1;
- if(((src_n>>src_bit)&1) != 0) dest_n |= (one<<dest_bit);
- else dest_n &= ~(one<<dest_bit);
-}
-*/
-
-/*
-static int
-_countOnes(const _64bits src_n)
-{
- int cnt = 0;
- for(int i=0; i<64; i++)
- cnt += (int)((src_n>>i)&1);
- return(cnt);
-}
-*/
-
-static int*
-toBitPattern(int b[], _64bits n, int n_bits)
-{
- for(int i=0; i<n_bits; i++)
- b[i] = (_isBitSet(n,i) ? 1:0);
- return b;
-}
-
-static _64bits
-fromBitPattern(int b[], int n_bits)
-{
- _64bits n = 0;
- for(int i=0; i<n_bits; i++)
- if(b[i] == 1) _setBit(n,i);
- return(n);
-}
-
-/*
-static void
-printBitPattern(_64bits n, int n_bits)
-{
- for(int i=0; i<n_bits; i++)
- {
- if(_isBitSet(n,i)) printf("1");
- else printf("0");
- }
-}
-*/
-
-/*
-static void
-printBitArray(int bb[], int bb_length)
-{
- for(int i=0; i<bb_length; i++) printf("%i",bb[i]);
-}
-*/
-
-
-BCH::BCH()
-{
- initialize(BCH_DEFAULT_M, BCH_DEFAULT_LENGTH, BCH_DEFAULT_T);
-}
-
-BCH::BCH(int _m, int _length, int _t)
-{
- initialize(_m,_length,_t);
-}
-
-
-void BCH::initialize(int _m, int _length, int _t)
-{
- int i, ninf;
-
- m = _m;
- length = _length;
- t = _t;
-
- p.resize(BCH_MAX_P);
- alpha_to.resize(BCH_MAX_LUT);
- index_of.resize(BCH_MAX_LUT);
- g.resize(BCH_MAX_LUT);
-
- _elp.resize(BCH_MAX_LUT);
- for(i=0; i<BCH_MAX_LUT; i++) _elp[i].resize(BCH_MAX_LUT);
- _d.resize(BCH_MAX_LUT);
- _l.resize(BCH_MAX_LUT);
- _u_lu.resize(BCH_MAX_LUT);
- _s.resize(BCH_MAX_LUT);
- _root.resize(BCH_MAX_LUT);
- _loc.resize(BCH_MAX_LUT);
- _reg.resize(BCH_MAX_LUT);
-
-
- for (i=1; i<m; i++)
- p[i] = 0;
- p[0] = p[m] = 1;
- if (m == 2) p[1] = 1;
- else if (m == 3) p[1] = 1;
- else if (m == 4) p[1] = 1;
- else if (m == 5) p[2] = 1;
- else if (m == 6) p[1] = 1;
- else if (m == 7) p[1] = 1;
- else if (m == 8) p[4] = p[5] = p[6] = 1;
- else if (m == 9) p[4] = 1;
- else if (m == 10) p[3] = 1;
- else if (m == 11) p[2] = 1;
- else if (m == 12) p[3] = p[4] = p[7] = 1;
- else if (m == 13) p[1] = p[3] = p[4] = 1;
- else if (m == 14) p[1] = p[11] = p[12] = 1;
- else if (m == 15) p[1] = 1;
- else if (m == 16) p[2] = p[3] = p[5] = 1;
- else if (m == 17) p[3] = 1;
- else if (m == 18) p[7] = 1;
- else if (m == 19) p[1] = p[5] = p[6] = 1;
- else if (m == 20) p[3] = 1;
-
- n = 1;
- for (i = 0; i <= m; i++)
- {
- n *= 2;
- }
- n = n / 2 - 1;
- ninf = (n + 1) / 2 - 1;
-
- generate_gf(); /* Construct the Galois Field GF(2**m) */
- gen_poly(t); /* Compute the generator polynomial of BCH code */
-}
-
-
-void BCH::generate_gf()
-/*
- * Generate field GF(2**m) from the irreducible polynomial p(X) with
- * coefficients in p[0]..p[m].
- *
- * Lookup tables:
- * index->polynomial form: alpha_to[] contains j=alpha^i;
- * polynomial form -> index form: index_of[j=alpha^i] = i
- *
- * alpha=2 is the primitive element of GF(2**m)
- */
-{
- int i, mask;
-
- mask = 1;
- alpha_to[m] = 0;
- for (i = 0; i < m; i++) {
- alpha_to[i] = mask;
- index_of[alpha_to[i]] = i;
- if (p[i] != 0)
- alpha_to[m] ^= mask;
- mask <<= 1;
- }
- index_of[alpha_to[m]] = m;
- mask >>= 1;
- for (i = m + 1; i < n; i++) {
- if (alpha_to[i - 1] >= mask)
- alpha_to[i] = alpha_to[m] ^ ((alpha_to[i - 1] ^ mask) << 1);
- else
- alpha_to[i] = alpha_to[i - 1] << 1;
- index_of[alpha_to[i]] = i;
- }
- index_of[0] = -1;
-}
-
-
-bool BCH::gen_poly(int _t)
-/*
- * Compute the generator polynomial of a binary BCH code. Fist generate the
- * cycle sets modulo 2**m - 1, cycle[][] = (i, 2*i, 4*i, ..., 2^l*i). Then
- * determine those cycle sets that contain integers in the set of (d-1)
- * consecutive integers {1..(d-1)}. The generator polynomial is calculated
- * as the product of linear factors of the form (x+alpha^i), for every i in
- * the above cycle sets.
- */
-{
- int ii, jj, ll, kaux;
- int test, aux, nocycles, root, noterms, rdncy;
- int cycle[1024][21], size[1024], min[1024], zeros[1024];
-
- /* Generate cycle sets modulo n, n = 2**m - 1 */
- cycle[0][0] = 0;
- size[0] = 1;
- cycle[1][0] = 1;
- size[1] = 1;
- jj = 1; /* cycle set index */
-
- #ifdef _BCH_VERBOSE_
- if (m > 9) {
- printf("Computing cycle sets modulo %d\n", n);
- printf("(This may take some time)...\n");
- }
- #endif
-
- do {
- /* Generate the jj-th cycle set */
- ii = 0;
- do {
- ii++;
- cycle[jj][ii] = (cycle[jj][ii - 1] * 2) % n;
- size[jj]++;
- aux = (cycle[jj][ii] * 2) % n;
- } while (aux != cycle[jj][0]);
- /* Next cycle set representative */
- ll = 0;
- do {
- ll++;
- test = 0;
- for (ii = 1; ((ii <= jj) && (!test)); ii++)
- /* Examine previous cycle sets */
- for (kaux = 0; ((kaux < size[ii]) && (!test)); kaux++)
- if (ll == cycle[ii][kaux])
- test = 1;
- } while ((test) && (ll < (n - 1)));
- if (!(test)) {
- jj++; /* next cycle set index */
- cycle[jj][0] = ll;
- size[jj] = 1;
- }
- } while (ll < (n - 1));
- nocycles = jj; /* number of cycle sets modulo n */
-
- //printf("Enter the error correcting capability, t: ");
- //scanf("%d", &t);
- t = _t;
-
- d = 2 * t + 1;
-
- /* Search for roots 1, 2, ..., d-1 in cycle sets */
- kaux = 0;
- rdncy = 0;
- for (ii = 1; ii <= nocycles; ii++) {
- min[kaux] = 0;
- test = 0;
- for (jj = 0; ((jj < size[ii]) && (!test)); jj++)
- for (root = 1; ((root < d) && (!test)); root++)
- if (root == cycle[ii][jj]) {
- test = 1;
- min[kaux] = ii;
- }
- if (min[kaux]) {
- rdncy += size[min[kaux]];
- kaux++;
- }
- }
- noterms = kaux;
- kaux = 1;
- for (ii = 0; ii < noterms; ii++)
- for (jj = 0; jj < size[min[ii]]; jj++) {
- zeros[kaux] = cycle[min[ii]][jj];
- kaux++;
- }
-
- k = length - rdncy;
-
- if (k<0)
- {
- #ifdef _BCH_VERBOSE_
- printf("Parameters invalid!\n");
- #endif
- return(false);
- }
-
- #ifdef _BCH_VERBOSE_
- printf("This is a (%d, %d, %d) binary BCH code\n", length, k, d);
- #endif
-
- /* Compute the generator polynomial */
- g[0] = alpha_to[zeros[1]];
- g[1] = 1; /* g(x) = (X + zeros[1]) initially */
- for (ii = 2; ii <= rdncy; ii++) {
- g[ii] = 1;
- for (jj = ii - 1; jj > 0; jj--)
- if (g[jj] != 0)
- g[jj] = g[jj - 1] ^ alpha_to[(index_of[g[jj]] + zeros[ii]) % n];
- else
- g[jj] = g[jj - 1];
- g[0] = alpha_to[(index_of[g[0]] + zeros[ii]) % n];
- }
-
- #ifdef _BCH_VERBOSE_
- printf("Generator polynomial:\ng(x) = ");
- for (ii = 0; ii <= rdncy; ii++) {
- printf("%d", g[ii]);
- if (ii && ((ii % 50) == 0))
- printf("\n");
- }
- printf("\n");
- #endif
- return(true);
-}
-
-
-
-void BCH::encode_bch(int *bb, const int *data)
-/*
- * Compute redundacy bb[], the coefficients of b(x). The redundancy
- * polynomial b(x) is the remainder after dividing x^(length-k)*data(x)
- * by the generator polynomial g(x).
- */
-{
- int i, j;
- int feedback;
-
- for (i = 0; i < length - k; i++)
- bb[i] = 0;
- for (i = k - 1; i >= 0; i--) {
- feedback = data[i] ^ bb[length - k - 1];
- if (feedback != 0) {
- for (j = length - k - 1; j > 0; j--)
- if (g[j] != 0)
- bb[j] = bb[j - 1] ^ feedback;
- else
- bb[j] = bb[j - 1];
- bb[0] = g[0] && feedback;
- } else {
- for (j = length - k - 1; j > 0; j--)
- bb[j] = bb[j - 1];
- bb[0] = 0;
- }
- }
-}
-
-
-int BCH::decode_bch(int *recd)
-/*
- * Simon Rockliff's implementation of Berlekamp's algorithm.
- *
- * Assume we have received bits in recd[i], i=0..(n-1).
- *
- * Compute the 2*t syndromes by substituting alpha^i into rec(X) and
- * evaluating, storing the syndromes in s[i], i=1..2t (leave s[0] zero) .
- * Then we use the Berlekamp algorithm to find the error location polynomial
- * elp[i].
- *
- * If the degree of the elp is >t, then we cannot correct all the errors, and
- * we have detected an uncorrectable error pattern. We output the information
- * bits uncorrected.
- *
- * If the degree of elp is <=t, we substitute alpha^i , i=1..n into the elp
- * to get the roots, hence the inverse roots, the error location numbers.
- * This step is usually called "Chien's search".
- *
- * If the number of errors located is not equal the degree of the elp, then
- * the decoder assumes that there are more than t errors and cannot correct
- * them, only detect them. We output the information bits uncorrected.
- */
-{
- int i, j, u, q, t2, count = 0, syn_error = 0;
- bool too_many_errors = false;
- //int elp[BCH_MAX_LUT][BCH_MAX_LUT], d[BCH_MAX_LUT], l[BCH_MAX_LUT], u_lu[BCH_MAX_LUT], s[BCH_MAX_LUT];
- //int root[BCH_MAX_SQ], loc[BCH_MAX_SQ], reg[BCH_MAX_SQ];
-
- t2 = 2 * t;
-
- /* first form the syndromes */
- #ifdef _BCH_VERBOSE_
- printf("S(x) = ");
- #endif
- for (i = 1; i <= t2; i++) {
- _s[i] = 0;
- for (j = 0; j < length; j++)
- if (recd[j] != 0)
- _s[i] ^= alpha_to[(i * j) % n];
- if (_s[i] != 0)
- syn_error = 1; /* set error flag if non-zero syndrome */
-/*
- * Note: If the code is used only for ERROR DETECTION, then
- * exit program here indicating the presence of errors.
- */
- /* convert syndrome from polynomial form to index form */
- _s[i] = index_of[_s[i]];
- #ifdef _BCH_VERBOSE_
- printf("%3d ", _s[i]);
- #endif
- }
-
- #ifdef _BCH_VERBOSE_
- printf("\n");
- #endif
-
- if (syn_error) { /* if there are errors, try to correct them */
- /*
- * Compute the error location polynomial via the Berlekamp
- * iterative algorithm. Following the terminology of Lin and
- * Costello's book : d[u] is the 'mu'th discrepancy, where
- * u='mu'+1 and 'mu' (the Greek letter!) is the step number
- * ranging from -1 to 2*t (see L&C), l[u] is the degree of
- * the elp at that step, and u_l[u] is the difference between
- * the step number and the degree of the elp.
- */
- /* initialise table entries */
- _d[0] = 0; /* index form */
- _d[1] = _s[1]; /* index form */
- _elp[0][0] = 0; /* index form */
- _elp[1][0] = 1; /* polynomial form */
- for (i = 1; i < t2; i++) {
- _elp[0][i] = -1; /* index form */
- _elp[1][i] = 0; /* polynomial form */
- }
- _l[0] = 0;
- _l[1] = 0;
- _u_lu[0] = -1;
- _u_lu[1] = 0;
- u = 0;
-
- do {
- u++;
- if (_d[u] == -1) {
- _l[u + 1] = _l[u];
- for (i = 0; i <= _l[u]; i++) {
- _elp[u + 1][i] = _elp[u][i];
- _elp[u][i] = index_of[_elp[u][i]];
- }
- } else
- /*
- * search for words with greatest _u_lu[q] for
- * which _d[q]!=0
- */
- {
- q = u - 1;
- while ((_d[q] == -1) && (q > 0))
- q--;
- /* have found first non-zero _d[q] */
- if (q > 0) {
- j = q;
- do {
- j--;
- if ((_d[j] != -1) && (_u_lu[q] < _u_lu[j]))
- q = j;
- } while (j > 0);
- }
-
- /*
- * have now found q such that _d[u]!=0 and
- * _u_lu[q] is maximum
- */
- /* store degree of new elp polynomial */
- if (_l[u] > _l[q] + u - q)
- _l[u + 1] = _l[u];
- else
- _l[u + 1] = _l[q] + u - q;
-
- /* form new elp(x) */
- for (i = 0; i < t2; i++)
- _elp[u + 1][i] = 0;
- for (i = 0; i <= _l[q]; i++)
- if (_elp[q][i] != -1)
- _elp[u + 1][i + u - q] =
- alpha_to[(_d[u] + n - _d[q] + _elp[q][i]) % n];
- for (i = 0; i <= _l[u]; i++) {
- _elp[u + 1][i] ^= _elp[u][i];
- _elp[u][i] = index_of[_elp[u][i]];
- }
- }
- _u_lu[u + 1] = u - _l[u + 1];
-
- /* form (u+1)th discrepancy */
- if (u < t2) {
- /* no discrepancy computed on last iteration */
- if (_s[u + 1] != -1)
- _d[u + 1] = alpha_to[_s[u + 1]];
- else
- _d[u + 1] = 0;
- for (i = 1; i <= _l[u + 1]; i++)
- if ((_s[u + 1 - i] != -1) && (_elp[u + 1][i] != 0))
- _d[u + 1] ^= alpha_to[(_s[u + 1 - i]
- + index_of[_elp[u + 1][i]]) % n];
- /* put _d[u+1] into index form */
- _d[u + 1] = index_of[_d[u + 1]];
- }
- } while ((u < t2) && (_l[u + 1] <= t));
-
- u++;
- if (_l[u] <= t) {/* Can correct errors */
- /* put elp into index form */
- for (i = 0; i <= _l[u]; i++)
- _elp[u][i] = index_of[_elp[u][i]];
-
-
-#ifdef _BCH_VERBOSE_
- printf("sigma(x) = ");
- for (i = 0; i <= _l[u]; i++)
- printf("%3d ", _elp[u][i]);
- printf("\n");
- printf("Roots: ");
-#endif
-
- /* Chien search: find roots of the error location polynomial */
- for (i = 1; i <= _l[u]; i++)
- _reg[i] = _elp[u][i];
- count = 0;
- for (i = 1; i <= n; i++) {
- q = 1;
- for (j = 1; j <= _l[u]; j++)
- if (_reg[j] != -1) {
- _reg[j] = (_reg[j] + j) % n;
- q ^= alpha_to[_reg[j]];
- }
- if (!q) { /* store root and error
- * location number indices */
- _root[count] = i;
- _loc[count] = n - i;
- count++;
- #ifdef _BCH_VERBOSE_
- printf("%3d ", n - i);
- #endif
- }
- }
- #ifdef _BCH_VERBOSE_
- printf("\n");
- #endif
- if (count == _l[u])
- {
- /* no. roots = degree of elp hence <= t errors */
- for (i = 0; i < _l[u]; i++)
- {
- if(_loc[i]<BCH_DEFAULT_LENGTH) recd[_loc[i]] ^= 1;
- else too_many_errors = true;
- }
- }
- else /* elp has degree >t hence cannot solve */
- {
- return(_l[u]); // number of errors
- #ifdef _BCH_VERBOSE_
- printf("Incomplete decoding: errors detected\n");
- #endif
- }
- }
- }
- if(too_many_errors) return(BCH_DEFAULT_T+1);
- else return(syn_error == 0 ? 0 : _l[u]); // number of errors
-}
-
-
-void BCH::encode(int encoded_bits[BCH_DEFAULT_LENGTH], const _64bits orig_n)
-{
- assert(k == BCH_DEFAULT_K && length == BCH_DEFAULT_LENGTH);
- int orig_bits[BCH_DEFAULT_K];
- toBitPattern(orig_bits,orig_n,k);
- encode_bch(encoded_bits,orig_bits);
-
- for (int i = 0; i < k; i++)
- encoded_bits[i + length - k] = orig_bits[i];
-
- //printf("[x] encoded_bits: ");
- //printBitArray(encoded_bits, BCH_DEFAULT_LENGTH);
- //printf("\n");
-
-}
-
-bool BCH::decode(int &err_n, _64bits &orig_n, const int encoded_bits[BCH_DEFAULT_LENGTH])
-{
- assert(k == BCH_DEFAULT_K && length == BCH_DEFAULT_LENGTH);
- int temp_bits[BCH_DEFAULT_LENGTH];
- for(int i=0; i<BCH_DEFAULT_LENGTH; i++) temp_bits[i] = encoded_bits[i];
- err_n = decode_bch(temp_bits);
- if(err_n > t) return(false);
- orig_n = fromBitPattern(&temp_bits[length - k],BCH_DEFAULT_K);
- return(true);
-}
-
-void BCH::encode(_64bits &encoded_n, const _64bits orig_n)
-{
- assert(k == BCH_DEFAULT_K && length == BCH_DEFAULT_LENGTH);
- int encoded_bits[BCH_DEFAULT_LENGTH];
- encode(encoded_bits,orig_n);
- encoded_n = fromBitPattern(encoded_bits,BCH_DEFAULT_LENGTH);
-}
-
-bool BCH::decode(int &err_n, _64bits &orig_n, const _64bits encoded_n)
-{
- assert(k == BCH_DEFAULT_K && length == BCH_DEFAULT_LENGTH);
- int encoded_bits[BCH_DEFAULT_LENGTH];
- toBitPattern(encoded_bits,encoded_n,BCH_DEFAULT_LENGTH);
- return(decode(err_n, orig_n,encoded_bits));
-}
-
-
-} // namespace ARToolKitPlus
|