From 2d40ed7b415bcb6a358c9f221160c4a419f01f19 Mon Sep 17 00:00:00 2001 From: Artem Bityutskiy Date: Sat, 17 Jul 2010 09:43:08 +0300 Subject: libs: make crc32 and fec to be libraries Signed-off-by: Artem Bityutskiy --- lib/Makefile | 6 +- lib/crc32.c | 95 ------ lib/fec.c | 904 --------------------------------------------------------- lib/libcrc32.c | 104 +++++++ lib/libfec.c | 904 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 5 files changed, 1011 insertions(+), 1002 deletions(-) delete mode 100644 lib/crc32.c delete mode 100644 lib/fec.c create mode 100644 lib/libcrc32.c create mode 100644 lib/libfec.c (limited to 'lib') diff --git a/lib/Makefile b/lib/Makefile index 82194e4..621ed63 100644 --- a/lib/Makefile +++ b/lib/Makefile @@ -6,13 +6,13 @@ SUBDIRS = # CFLAGS += -Werror CPPFLAGS += -I../include -LIBS = libmtd -TARGETS = libmtd.a +LIBS = libmtd libcrc32 +TARGETS = libmtd.a libcrc32.a include ../common.mk $(BUILDDIR)/libmtd.a: $(addprefix $(BUILDDIR)/,\ - libmtd.o libmtd_legacy.o crc32.o fec.o) + libmtd.o libmtd_legacy.o libcrc32.o libfec.o) clean:: rm -f $(addsuffix .a, $(LIBS)) diff --git a/lib/crc32.c b/lib/crc32.c deleted file mode 100644 index 6b1e50c..0000000 --- a/lib/crc32.c +++ /dev/null @@ -1,95 +0,0 @@ -/* - * COPYRIGHT (C) 1986 Gary S. Brown. You may use this program, or - * code or tables extracted from it, as desired without restriction. - * - * First, the polynomial itself and its table of feedback terms. The - * polynomial is - * X^32+X^26+X^23+X^22+X^16+X^12+X^11+X^10+X^8+X^7+X^5+X^4+X^2+X^1+X^0 - * - * Note that we take it "backwards" and put the highest-order term in - * the lowest-order bit. The X^32 term is "implied"; the LSB is the - * X^31 term, etc. The X^0 term (usually shown as "+1") results in - * the MSB being 1 - * - * Note that the usual hardware shift register implementation, which - * is what we're using (we're merely optimizing it by doing eight-bit - * chunks at a time) shifts bits into the lowest-order term. In our - * implementation, that means shifting towards the right. Why do we - * do it this way? Because the calculated CRC must be transmitted in - * order from highest-order term to lowest-order term. UARTs transmit - * characters in order from LSB to MSB. By storing the CRC this way - * we hand it to the UART in the order low-byte to high-byte; the UART - * sends each low-bit to hight-bit; and the result is transmission bit - * by bit from highest- to lowest-order term without requiring any bit - * shuffling on our part. Reception works similarly - * - * The feedback terms table consists of 256, 32-bit entries. Notes - * - * The table can be generated at runtime if desired; code to do so - * is shown later. It might not be obvious, but the feedback - * terms simply represent the results of eight shift/xor opera - * tions for all combinations of data and CRC register values - * - * The values must be right-shifted by eight bits by the "updcrc - * logic; the shift must be unsigned (bring in zeroes). On some - * hardware you could probably optimize the shift in assembler by - * using byte-swap instructions - * polynomial $edb88320 - */ - -#include - -const uint32_t crc32_table[256] = { - 0x00000000L, 0x77073096L, 0xee0e612cL, 0x990951baL, 0x076dc419L, - 0x706af48fL, 0xe963a535L, 0x9e6495a3L, 0x0edb8832L, 0x79dcb8a4L, - 0xe0d5e91eL, 0x97d2d988L, 0x09b64c2bL, 0x7eb17cbdL, 0xe7b82d07L, - 0x90bf1d91L, 0x1db71064L, 0x6ab020f2L, 0xf3b97148L, 0x84be41deL, - 0x1adad47dL, 0x6ddde4ebL, 0xf4d4b551L, 0x83d385c7L, 0x136c9856L, - 0x646ba8c0L, 0xfd62f97aL, 0x8a65c9ecL, 0x14015c4fL, 0x63066cd9L, - 0xfa0f3d63L, 0x8d080df5L, 0x3b6e20c8L, 0x4c69105eL, 0xd56041e4L, - 0xa2677172L, 0x3c03e4d1L, 0x4b04d447L, 0xd20d85fdL, 0xa50ab56bL, - 0x35b5a8faL, 0x42b2986cL, 0xdbbbc9d6L, 0xacbcf940L, 0x32d86ce3L, - 0x45df5c75L, 0xdcd60dcfL, 0xabd13d59L, 0x26d930acL, 0x51de003aL, - 0xc8d75180L, 0xbfd06116L, 0x21b4f4b5L, 0x56b3c423L, 0xcfba9599L, - 0xb8bda50fL, 0x2802b89eL, 0x5f058808L, 0xc60cd9b2L, 0xb10be924L, - 0x2f6f7c87L, 0x58684c11L, 0xc1611dabL, 0xb6662d3dL, 0x76dc4190L, - 0x01db7106L, 0x98d220bcL, 0xefd5102aL, 0x71b18589L, 0x06b6b51fL, - 0x9fbfe4a5L, 0xe8b8d433L, 0x7807c9a2L, 0x0f00f934L, 0x9609a88eL, - 0xe10e9818L, 0x7f6a0dbbL, 0x086d3d2dL, 0x91646c97L, 0xe6635c01L, - 0x6b6b51f4L, 0x1c6c6162L, 0x856530d8L, 0xf262004eL, 0x6c0695edL, - 0x1b01a57bL, 0x8208f4c1L, 0xf50fc457L, 0x65b0d9c6L, 0x12b7e950L, - 0x8bbeb8eaL, 0xfcb9887cL, 0x62dd1ddfL, 0x15da2d49L, 0x8cd37cf3L, - 0xfbd44c65L, 0x4db26158L, 0x3ab551ceL, 0xa3bc0074L, 0xd4bb30e2L, - 0x4adfa541L, 0x3dd895d7L, 0xa4d1c46dL, 0xd3d6f4fbL, 0x4369e96aL, - 0x346ed9fcL, 0xad678846L, 0xda60b8d0L, 0x44042d73L, 0x33031de5L, - 0xaa0a4c5fL, 0xdd0d7cc9L, 0x5005713cL, 0x270241aaL, 0xbe0b1010L, - 0xc90c2086L, 0x5768b525L, 0x206f85b3L, 0xb966d409L, 0xce61e49fL, - 0x5edef90eL, 0x29d9c998L, 0xb0d09822L, 0xc7d7a8b4L, 0x59b33d17L, - 0x2eb40d81L, 0xb7bd5c3bL, 0xc0ba6cadL, 0xedb88320L, 0x9abfb3b6L, - 0x03b6e20cL, 0x74b1d29aL, 0xead54739L, 0x9dd277afL, 0x04db2615L, - 0x73dc1683L, 0xe3630b12L, 0x94643b84L, 0x0d6d6a3eL, 0x7a6a5aa8L, - 0xe40ecf0bL, 0x9309ff9dL, 0x0a00ae27L, 0x7d079eb1L, 0xf00f9344L, - 0x8708a3d2L, 0x1e01f268L, 0x6906c2feL, 0xf762575dL, 0x806567cbL, - 0x196c3671L, 0x6e6b06e7L, 0xfed41b76L, 0x89d32be0L, 0x10da7a5aL, - 0x67dd4accL, 0xf9b9df6fL, 0x8ebeeff9L, 0x17b7be43L, 0x60b08ed5L, - 0xd6d6a3e8L, 0xa1d1937eL, 0x38d8c2c4L, 0x4fdff252L, 0xd1bb67f1L, - 0xa6bc5767L, 0x3fb506ddL, 0x48b2364bL, 0xd80d2bdaL, 0xaf0a1b4cL, - 0x36034af6L, 0x41047a60L, 0xdf60efc3L, 0xa867df55L, 0x316e8eefL, - 0x4669be79L, 0xcb61b38cL, 0xbc66831aL, 0x256fd2a0L, 0x5268e236L, - 0xcc0c7795L, 0xbb0b4703L, 0x220216b9L, 0x5505262fL, 0xc5ba3bbeL, - 0xb2bd0b28L, 0x2bb45a92L, 0x5cb36a04L, 0xc2d7ffa7L, 0xb5d0cf31L, - 0x2cd99e8bL, 0x5bdeae1dL, 0x9b64c2b0L, 0xec63f226L, 0x756aa39cL, - 0x026d930aL, 0x9c0906a9L, 0xeb0e363fL, 0x72076785L, 0x05005713L, - 0x95bf4a82L, 0xe2b87a14L, 0x7bb12baeL, 0x0cb61b38L, 0x92d28e9bL, - 0xe5d5be0dL, 0x7cdcefb7L, 0x0bdbdf21L, 0x86d3d2d4L, 0xf1d4e242L, - 0x68ddb3f8L, 0x1fda836eL, 0x81be16cdL, 0xf6b9265bL, 0x6fb077e1L, - 0x18b74777L, 0x88085ae6L, 0xff0f6a70L, 0x66063bcaL, 0x11010b5cL, - 0x8f659effL, 0xf862ae69L, 0x616bffd3L, 0x166ccf45L, 0xa00ae278L, - 0xd70dd2eeL, 0x4e048354L, 0x3903b3c2L, 0xa7672661L, 0xd06016f7L, - 0x4969474dL, 0x3e6e77dbL, 0xaed16a4aL, 0xd9d65adcL, 0x40df0b66L, - 0x37d83bf0L, 0xa9bcae53L, 0xdebb9ec5L, 0x47b2cf7fL, 0x30b5ffe9L, - 0xbdbdf21cL, 0xcabac28aL, 0x53b39330L, 0x24b4a3a6L, 0xbad03605L, - 0xcdd70693L, 0x54de5729L, 0x23d967bfL, 0xb3667a2eL, 0xc4614ab8L, - 0x5d681b02L, 0x2a6f2b94L, 0xb40bbe37L, 0xc30c8ea1L, 0x5a05df1bL, - 0x2d02ef8dL -}; diff --git a/lib/fec.c b/lib/fec.c deleted file mode 100644 index 6d9020f..0000000 --- a/lib/fec.c +++ /dev/null @@ -1,904 +0,0 @@ -/* - * fec.c -- forward error correction based on Vandermonde matrices - * 980624 - * (C) 1997-98 Luigi Rizzo (luigi@iet.unipi.it) - * - * Portions derived from code by Phil Karn (karn@ka9q.ampr.org), - * Robert Morelos-Zaragoza (robert@spectra.eng.hawaii.edu) and Hari - * Thirumoorthy (harit@spectra.eng.hawaii.edu), Aug 1995 - * - * Redistribution and use in source and binary forms, with or without - * modification, are permitted provided that the following conditions - * are met: - * - * 1. Redistributions of source code must retain the above copyright - * notice, this list of conditions and the following disclaimer. - * 2. Redistributions in binary form must reproduce the above - * copyright notice, this list of conditions and the following - * disclaimer in the documentation and/or other materials - * provided with the distribution. - * - * THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND - * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, - * THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A - * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS - * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, - * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, - * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, - * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY - * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR - * TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT - * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY - * OF SUCH DAMAGE. - */ - -/* - * The following parameter defines how many bits are used for - * field elements. The code supports any value from 2 to 16 - * but fastest operation is achieved with 8 bit elements - * This is the only parameter you may want to change. - */ -#ifndef GF_BITS -#define GF_BITS 8 /* code over GF(2**GF_BITS) - change to suit */ -#endif - -#include -#include -#include - -/* - * stuff used for testing purposes only - */ - -#ifdef TEST -#define DEB(x) -#define DDB(x) x -#define DEBUG 0 /* minimal debugging */ -#ifdef MSDOS -#include -struct timeval { - unsigned long ticks; -}; -#define gettimeofday(x, dummy) { (x)->ticks = clock() ; } -#define DIFF_T(a,b) (1+ 1000000*(a.ticks - b.ticks) / CLOCKS_PER_SEC ) -typedef unsigned long u_long ; -typedef unsigned short u_short ; -#else /* typically, unix systems */ -#include -#define DIFF_T(a,b) \ - (1+ 1000000*(a.tv_sec - b.tv_sec) + (a.tv_usec - b.tv_usec) ) -#endif - -#define TICK(t) \ - {struct timeval x ; \ - gettimeofday(&x, NULL) ; \ - t = x.tv_usec + 1000000* (x.tv_sec & 0xff ) ; \ - } -#define TOCK(t) \ - { u_long t1 ; TICK(t1) ; \ - if (t1 < t) t = 256000000 + t1 - t ; \ - else t = t1 - t ; \ - if (t == 0) t = 1 ;} - -u_long ticks[10]; /* vars for timekeeping */ -#else -#define DEB(x) -#define DDB(x) -#define TICK(x) -#define TOCK(x) -#endif /* TEST */ - -/* - * You should not need to change anything beyond this point. - * The first part of the file implements linear algebra in GF. - * - * gf is the type used to store an element of the Galois Field. - * Must constain at least GF_BITS bits. - * - * Note: unsigned char will work up to GF(256) but int seems to run - * faster on the Pentium. We use int whenever have to deal with an - * index, since they are generally faster. - */ -#if (GF_BITS < 2 && GF_BITS >16) -#error "GF_BITS must be 2 .. 16" -#endif -#if (GF_BITS <= 8) -typedef unsigned char gf; -#else -typedef unsigned short gf; -#endif - -#define GF_SIZE ((1 << GF_BITS) - 1) /* powers of \alpha */ - -/* - * Primitive polynomials - see Lin & Costello, Appendix A, - * and Lee & Messerschmitt, p. 453. - */ -static char *allPp[] = { /* GF_BITS polynomial */ - NULL, /* 0 no code */ - NULL, /* 1 no code */ - "111", /* 2 1+x+x^2 */ - "1101", /* 3 1+x+x^3 */ - "11001", /* 4 1+x+x^4 */ - "101001", /* 5 1+x^2+x^5 */ - "1100001", /* 6 1+x+x^6 */ - "10010001", /* 7 1 + x^3 + x^7 */ - "101110001", /* 8 1+x^2+x^3+x^4+x^8 */ - "1000100001", /* 9 1+x^4+x^9 */ - "10010000001", /* 10 1+x^3+x^10 */ - "101000000001", /* 11 1+x^2+x^11 */ - "1100101000001", /* 12 1+x+x^4+x^6+x^12 */ - "11011000000001", /* 13 1+x+x^3+x^4+x^13 */ - "110000100010001", /* 14 1+x+x^6+x^10+x^14 */ - "1100000000000001", /* 15 1+x+x^15 */ - "11010000000010001" /* 16 1+x+x^3+x^12+x^16 */ -}; - - -/* - * To speed up computations, we have tables for logarithm, exponent - * and inverse of a number. If GF_BITS <= 8, we use a table for - * multiplication as well (it takes 64K, no big deal even on a PDA, - * especially because it can be pre-initialized an put into a ROM!), - * otherwhise we use a table of logarithms. - * In any case the macro gf_mul(x,y) takes care of multiplications. - */ - -static gf gf_exp[2*GF_SIZE]; /* index->poly form conversion table */ -static int gf_log[GF_SIZE + 1]; /* Poly->index form conversion table */ -static gf inverse[GF_SIZE+1]; /* inverse of field elem. */ - /* inv[\alpha**i]=\alpha**(GF_SIZE-i-1) */ - -/* - * modnn(x) computes x % GF_SIZE, where GF_SIZE is 2**GF_BITS - 1, - * without a slow divide. - */ -static inline gf -modnn(int x) -{ - while (x >= GF_SIZE) { - x -= GF_SIZE; - x = (x >> GF_BITS) + (x & GF_SIZE); - } - return x; -} - -#define SWAP(a,b,t) {t tmp; tmp=a; a=b; b=tmp;} - -/* - * gf_mul(x,y) multiplies two numbers. If GF_BITS<=8, it is much - * faster to use a multiplication table. - * - * USE_GF_MULC, GF_MULC0(c) and GF_ADDMULC(x) can be used when multiplying - * many numbers by the same constant. In this case the first - * call sets the constant, and others perform the multiplications. - * A value related to the multiplication is held in a local variable - * declared with USE_GF_MULC . See usage in addmul1(). - */ -#if (GF_BITS <= 8) -static gf gf_mul_table[GF_SIZE + 1][GF_SIZE + 1]; - -#define gf_mul(x,y) gf_mul_table[x][y] - -#define USE_GF_MULC register gf * __gf_mulc_ -#define GF_MULC0(c) __gf_mulc_ = gf_mul_table[c] -#define GF_ADDMULC(dst, x) dst ^= __gf_mulc_[x] - -static void -init_mul_table() -{ - int i, j; - for (i=0; i< GF_SIZE+1; i++) - for (j=0; j< GF_SIZE+1; j++) - gf_mul_table[i][j] = gf_exp[modnn(gf_log[i] + gf_log[j]) ] ; - - for (j=0; j< GF_SIZE+1; j++) - gf_mul_table[0][j] = gf_mul_table[j][0] = 0; -} -#else /* GF_BITS > 8 */ -static inline gf -gf_mul(x,y) -{ - if ( (x) == 0 || (y)==0 ) return 0; - - return gf_exp[gf_log[x] + gf_log[y] ] ; -} -#define init_mul_table() - -#define USE_GF_MULC register gf * __gf_mulc_ -#define GF_MULC0(c) __gf_mulc_ = &gf_exp[ gf_log[c] ] -#define GF_ADDMULC(dst, x) { if (x) dst ^= __gf_mulc_[ gf_log[x] ] ; } -#endif - -/* - * Generate GF(2**m) from the irreducible polynomial p(X) in p[0]..p[m] - * Lookup tables: - * index->polynomial form gf_exp[] contains j= \alpha^i; - * polynomial form -> index form gf_log[ j = \alpha^i ] = i - * \alpha=x is the primitive element of GF(2^m) - * - * For efficiency, gf_exp[] has size 2*GF_SIZE, so that a simple - * multiplication of two numbers can be resolved without calling modnn - */ - -/* - * i use malloc so many times, it is easier to put checks all in - * one place. - */ -static void * -my_malloc(int sz, char *err_string) -{ - void *p = malloc( sz ); - if (p == NULL) { - fprintf(stderr, "-- malloc failure allocating %s\n", err_string); - exit(1) ; - } - return p ; -} - -#define NEW_GF_MATRIX(rows, cols) \ - (gf *)my_malloc(rows * cols * sizeof(gf), " ## __LINE__ ## " ) - -/* - * initialize the data structures used for computations in GF. - */ -static void -generate_gf(void) -{ - int i; - gf mask; - char *Pp = allPp[GF_BITS] ; - - mask = 1; /* x ** 0 = 1 */ - gf_exp[GF_BITS] = 0; /* will be updated at the end of the 1st loop */ - /* - * first, generate the (polynomial representation of) powers of \alpha, - * which are stored in gf_exp[i] = \alpha ** i . - * At the same time build gf_log[gf_exp[i]] = i . - * The first GF_BITS powers are simply bits shifted to the left. - */ - for (i = 0; i < GF_BITS; i++, mask <<= 1 ) { - gf_exp[i] = mask; - gf_log[gf_exp[i]] = i; - /* - * If Pp[i] == 1 then \alpha ** i occurs in poly-repr - * gf_exp[GF_BITS] = \alpha ** GF_BITS - */ - if ( Pp[i] == '1' ) - gf_exp[GF_BITS] ^= mask; - } - /* - * now gf_exp[GF_BITS] = \alpha ** GF_BITS is complete, so can als - * compute its inverse. - */ - gf_log[gf_exp[GF_BITS]] = GF_BITS; - /* - * Poly-repr of \alpha ** (i+1) is given by poly-repr of - * \alpha ** i shifted left one-bit and accounting for any - * \alpha ** GF_BITS term that may occur when poly-repr of - * \alpha ** i is shifted. - */ - mask = 1 << (GF_BITS - 1 ) ; - for (i = GF_BITS + 1; i < GF_SIZE; i++) { - if (gf_exp[i - 1] >= mask) - gf_exp[i] = gf_exp[GF_BITS] ^ ((gf_exp[i - 1] ^ mask) << 1); - else - gf_exp[i] = gf_exp[i - 1] << 1; - gf_log[gf_exp[i]] = i; - } - /* - * log(0) is not defined, so use a special value - */ - gf_log[0] = GF_SIZE ; - /* set the extended gf_exp values for fast multiply */ - for (i = 0 ; i < GF_SIZE ; i++) - gf_exp[i + GF_SIZE] = gf_exp[i] ; - - /* - * again special cases. 0 has no inverse. This used to - * be initialized to GF_SIZE, but it should make no difference - * since noone is supposed to read from here. - */ - inverse[0] = 0 ; - inverse[1] = 1; - for (i=2; i<=GF_SIZE; i++) - inverse[i] = gf_exp[GF_SIZE-gf_log[i]]; -} - -/* - * Various linear algebra operations that i use often. - */ - -/* - * addmul() computes dst[] = dst[] + c * src[] - * This is used often, so better optimize it! Currently the loop is - * unrolled 16 times, a good value for 486 and pentium-class machines. - * The case c=0 is also optimized, whereas c=1 is not. These - * calls are unfrequent in my typical apps so I did not bother. - * - * Note that gcc on - */ -#define addmul(dst, src, c, sz) \ - if (c != 0) addmul1(dst, src, c, sz) - -#define UNROLL 16 /* 1, 4, 8, 16 */ -static void -addmul1(gf *dst1, gf *src1, gf c, int sz) -{ - USE_GF_MULC ; - register gf *dst = dst1, *src = src1 ; - gf *lim = &dst[sz - UNROLL + 1] ; - - GF_MULC0(c) ; - -#if (UNROLL > 1) /* unrolling by 8/16 is quite effective on the pentium */ - for (; dst < lim ; dst += UNROLL, src += UNROLL ) { - GF_ADDMULC( dst[0] , src[0] ); - GF_ADDMULC( dst[1] , src[1] ); - GF_ADDMULC( dst[2] , src[2] ); - GF_ADDMULC( dst[3] , src[3] ); -#if (UNROLL > 4) - GF_ADDMULC( dst[4] , src[4] ); - GF_ADDMULC( dst[5] , src[5] ); - GF_ADDMULC( dst[6] , src[6] ); - GF_ADDMULC( dst[7] , src[7] ); -#endif -#if (UNROLL > 8) - GF_ADDMULC( dst[8] , src[8] ); - GF_ADDMULC( dst[9] , src[9] ); - GF_ADDMULC( dst[10] , src[10] ); - GF_ADDMULC( dst[11] , src[11] ); - GF_ADDMULC( dst[12] , src[12] ); - GF_ADDMULC( dst[13] , src[13] ); - GF_ADDMULC( dst[14] , src[14] ); - GF_ADDMULC( dst[15] , src[15] ); -#endif - } -#endif - lim += UNROLL - 1 ; - for (; dst < lim; dst++, src++ ) /* final components */ - GF_ADDMULC( *dst , *src ); -} - -/* - * computes C = AB where A is n*k, B is k*m, C is n*m - */ -static void -matmul(gf *a, gf *b, gf *c, int n, int k, int m) -{ - int row, col, i ; - - for (row = 0; row < n ; row++) { - for (col = 0; col < m ; col++) { - gf *pa = &a[ row * k ]; - gf *pb = &b[ col ]; - gf acc = 0 ; - for (i = 0; i < k ; i++, pa++, pb += m ) - acc ^= gf_mul( *pa, *pb ) ; - c[ row * m + col ] = acc ; - } - } -} - -#ifdef DEBUG -/* - * returns 1 if the square matrix is identiy - * (only for test) - */ -static int -is_identity(gf *m, int k) -{ - int row, col ; - for (row=0; row 1) { - fprintf(stderr, "singular matrix\n"); - goto fail ; - } - } - } - } - if (icol == -1) { - fprintf(stderr, "XXX pivot not found!\n"); - goto fail ; - } -found_piv: - ++(ipiv[icol]) ; - /* - * swap rows irow and icol, so afterwards the diagonal - * element will be correct. Rarely done, not worth - * optimizing. - */ - if (irow != icol) { - for (ix = 0 ; ix < k ; ix++ ) { - SWAP( src[irow*k + ix], src[icol*k + ix], gf) ; - } - } - indxr[col] = irow ; - indxc[col] = icol ; - pivot_row = &src[icol*k] ; - c = pivot_row[icol] ; - if (c == 0) { - fprintf(stderr, "singular matrix 2\n"); - goto fail ; - } - if (c != 1 ) { /* otherwhise this is a NOP */ - /* - * this is done often , but optimizing is not so - * fruitful, at least in the obvious ways (unrolling) - */ - DEB( pivswaps++ ; ) - c = inverse[ c ] ; - pivot_row[icol] = 1 ; - for (ix = 0 ; ix < k ; ix++ ) - pivot_row[ix] = gf_mul(c, pivot_row[ix] ); - } - /* - * from all rows, remove multiples of the selected row - * to zero the relevant entry (in fact, the entry is not zero - * because we know it must be zero). - * (Here, if we know that the pivot_row is the identity, - * we can optimize the addmul). - */ - id_row[icol] = 1; - if (memcmp(pivot_row, id_row, k*sizeof(gf)) != 0) { - for (p = src, ix = 0 ; ix < k ; ix++, p += k ) { - if (ix != icol) { - c = p[icol] ; - p[icol] = 0 ; - addmul(p, pivot_row, c, k ); - } - } - } - id_row[icol] = 0; - } /* done all columns */ - for (col = k-1 ; col >= 0 ; col-- ) { - if (indxr[col] <0 || indxr[col] >= k) - fprintf(stderr, "AARGH, indxr[col] %d\n", indxr[col]); - else if (indxc[col] <0 || indxc[col] >= k) - fprintf(stderr, "AARGH, indxc[col] %d\n", indxc[col]); - else - if (indxr[col] != indxc[col] ) { - for (row = 0 ; row < k ; row++ ) { - SWAP( src[row*k + indxr[col]], src[row*k + indxc[col]], gf) ; - } - } - } - error = 0 ; -fail: - free(indxc); - free(indxr); - free(ipiv); - free(id_row); - free(temp_row); - return error ; -} - -/* - * fast code for inverting a vandermonde matrix. - * XXX NOTE: It assumes that the matrix - * is not singular and _IS_ a vandermonde matrix. Only uses - * the second column of the matrix, containing the p_i's. - * - * Algorithm borrowed from "Numerical recipes in C" -- sec.2.8, but - * largely revised for my purposes. - * p = coefficients of the matrix (p_i) - * q = values of the polynomial (known) - */ - -int -invert_vdm(gf *src, int k) -{ - int i, j, row, col ; - gf *b, *c, *p; - gf t, xx ; - - if (k == 1) /* degenerate case, matrix must be p^0 = 1 */ - return 0 ; - /* - * c holds the coefficient of P(x) = Prod (x - p_i), i=0..k-1 - * b holds the coefficient for the matrix inversion - */ - c = NEW_GF_MATRIX(1, k); - b = NEW_GF_MATRIX(1, k); - - p = NEW_GF_MATRIX(1, k); - - for ( j=1, i = 0 ; i < k ; i++, j+=k ) { - c[i] = 0 ; - p[i] = src[j] ; /* p[i] */ - } - /* - * construct coeffs. recursively. We know c[k] = 1 (implicit) - * and start P_0 = x - p_0, then at each stage multiply by - * x - p_i generating P_i = x P_{i-1} - p_i P_{i-1} - * After k steps we are done. - */ - c[k-1] = p[0] ; /* really -p(0), but x = -x in GF(2^m) */ - for (i = 1 ; i < k ; i++ ) { - gf p_i = p[i] ; /* see above comment */ - for (j = k-1 - ( i - 1 ) ; j < k-1 ; j++ ) - c[j] ^= gf_mul( p_i, c[j+1] ) ; - c[k-1] ^= p_i ; - } - - for (row = 0 ; row < k ; row++ ) { - /* - * synthetic division etc. - */ - xx = p[row] ; - t = 1 ; - b[k-1] = 1 ; /* this is in fact c[k] */ - for (i = k-2 ; i >= 0 ; i-- ) { - b[i] = c[i+1] ^ gf_mul(xx, b[i+1]) ; - t = gf_mul(xx, t) ^ b[i] ; - } - for (col = 0 ; col < k ; col++ ) - src[col*k + row] = gf_mul(inverse[t], b[col] ); - } - free(c) ; - free(b) ; - free(p) ; - return 0 ; -} - -static int fec_initialized = 0 ; -static void -init_fec() -{ - TICK(ticks[0]); - generate_gf(); - TOCK(ticks[0]); - DDB(fprintf(stderr, "generate_gf took %ldus\n", ticks[0]);) - TICK(ticks[0]); - init_mul_table(); - TOCK(ticks[0]); - DDB(fprintf(stderr, "init_mul_table took %ldus\n", ticks[0]);) - fec_initialized = 1 ; -} - -/* - * This section contains the proper FEC encoding/decoding routines. - * The encoding matrix is computed starting with a Vandermonde matrix, - * and then transforming it into a systematic matrix. - */ - -#define FEC_MAGIC 0xFECC0DEC - -struct fec_parms { - u_long magic ; - int k, n ; /* parameters of the code */ - gf *enc_matrix ; -} ; - -void -fec_free(struct fec_parms *p) -{ - if (p==NULL || - p->magic != ( ( (FEC_MAGIC ^ p->k) ^ p->n) ^ (int)(p->enc_matrix)) ) { - fprintf(stderr, "bad parameters to fec_free\n"); - return ; - } - free(p->enc_matrix); - free(p); -} - -/* - * create a new encoder, returning a descriptor. This contains k,n and - * the encoding matrix. - */ -struct fec_parms * -fec_new(int k, int n) -{ - int row, col ; - gf *p, *tmp_m ; - - struct fec_parms *retval ; - - if (fec_initialized == 0) - init_fec(); - - if (k > GF_SIZE + 1 || n > GF_SIZE + 1 || k > n ) { - fprintf(stderr, "Invalid parameters k %d n %d GF_SIZE %d\n", - k, n, GF_SIZE ); - return NULL ; - } - retval = my_malloc(sizeof(struct fec_parms), "new_code"); - retval->k = k ; - retval->n = n ; - retval->enc_matrix = NEW_GF_MATRIX(n, k); - retval->magic = ( ( FEC_MAGIC ^ k) ^ n) ^ (int)(retval->enc_matrix) ; - tmp_m = NEW_GF_MATRIX(n, k); - /* - * fill the matrix with powers of field elements, starting from 0. - * The first row is special, cannot be computed with exp. table. - */ - tmp_m[0] = 1 ; - for (col = 1; col < k ; col++) - tmp_m[col] = 0 ; - for (p = tmp_m + k, row = 0; row < n-1 ; row++, p += k) { - for ( col = 0 ; col < k ; col ++ ) - p[col] = gf_exp[modnn(row*col)]; - } - - /* - * quick code to build systematic matrix: invert the top - * k*k vandermonde matrix, multiply right the bottom n-k rows - * by the inverse, and construct the identity matrix at the top. - */ - TICK(ticks[3]); - invert_vdm(tmp_m, k); /* much faster than invert_mat */ - matmul(tmp_m + k*k, tmp_m, retval->enc_matrix + k*k, n - k, k, k); - /* - * the upper matrix is I so do not bother with a slow multiply - */ - memset(retval->enc_matrix, '\0', k*k*sizeof(gf) ); - for (p = retval->enc_matrix, col = 0 ; col < k ; col++, p += k+1 ) - *p = 1 ; - free(tmp_m); - TOCK(ticks[3]); - - DDB(fprintf(stderr, "--- %ld us to build encoding matrix\n", - ticks[3]);) - DEB(pr_matrix(retval->enc_matrix, n, k, "encoding_matrix");) - return retval ; -} - -/* - * fec_encode accepts as input pointers to n data packets of size sz, - * and produces as output a packet pointed to by fec, computed - * with index "index". - */ -void -fec_encode(struct fec_parms *code, gf *src[], gf *fec, int index, int sz) -{ - int i, k = code->k ; - gf *p ; - - if (GF_BITS > 8) - sz /= 2 ; - - if (index < k) - memcpy(fec, src[index], sz*sizeof(gf) ) ; - else if (index < code->n) { - p = &(code->enc_matrix[index*k] ); - memset(fec, '\0', sz*sizeof(gf)); - for (i = 0; i < k ; i++) - addmul(fec, src[i], p[i], sz ) ; - } else - fprintf(stderr, "Invalid index %d (max %d)\n", - index, code->n - 1 ); -} - -void fec_encode_linear(struct fec_parms *code, gf *src, gf *fec, int index, int sz) -{ - int i, k = code->k ; - gf *p ; - - if (GF_BITS > 8) - sz /= 2 ; - - if (index < k) - memcpy(fec, src + (index * sz), sz*sizeof(gf) ) ; - else if (index < code->n) { - p = &(code->enc_matrix[index*k] ); - memset(fec, '\0', sz*sizeof(gf)); - for (i = 0; i < k ; i++) - addmul(fec, src + (i * sz), p[i], sz ) ; - } else - fprintf(stderr, "Invalid index %d (max %d)\n", - index, code->n - 1 ); -} -/* - * shuffle move src packets in their position - */ -static int -shuffle(gf *pkt[], int index[], int k) -{ - int i; - - for ( i = 0 ; i < k ; ) { - if (index[i] >= k || index[i] == i) - i++ ; - else { - /* - * put pkt in the right position (first check for conflicts). - */ - int c = index[i] ; - - if (index[c] == c) { - DEB(fprintf(stderr, "\nshuffle, error at %d\n", i);) - return 1 ; - } - SWAP(index[i], index[c], int) ; - SWAP(pkt[i], pkt[c], gf *) ; - } - } - DEB( /* just test that it works... */ - for ( i = 0 ; i < k ; i++ ) { - if (index[i] < k && index[i] != i) { - fprintf(stderr, "shuffle: after\n"); - for (i=0; ik ; - gf *p, *matrix = NEW_GF_MATRIX(k, k); - - TICK(ticks[9]); - for (i = 0, p = matrix ; i < k ; i++, p += k ) { -#if 1 /* this is simply an optimization, not very useful indeed */ - if (index[i] < k) { - memset(p, '\0', k*sizeof(gf) ); - p[i] = 1 ; - } else -#endif - if (index[i] < code->n ) - memcpy(p, &(code->enc_matrix[index[i]*k]), k*sizeof(gf) ); - else { - fprintf(stderr, "decode: invalid index %d (max %d)\n", - index[i], code->n - 1 ); - free(matrix) ; - return NULL ; - } - } - TICK(ticks[9]); - if (invert_mat(matrix, k)) { - free(matrix); - matrix = NULL ; - } - TOCK(ticks[9]); - return matrix ; -} - -/* - * fec_decode receives as input a vector of packets, the indexes of - * packets, and produces the correct vector as output. - * - * Input: - * code: pointer to code descriptor - * pkt: pointers to received packets. They are modified - * to store the output packets (in place) - * index: pointer to packet indexes (modified) - * sz: size of each packet - */ -int -fec_decode(struct fec_parms *code, gf *pkt[], int index[], int sz) -{ - gf *m_dec ; - gf **new_pkt ; - int row, col , k = code->k ; - - if (GF_BITS > 8) - sz /= 2 ; - - if (shuffle(pkt, index, k)) /* error if true */ - return 1 ; - m_dec = build_decode_matrix(code, pkt, index); - - if (m_dec == NULL) - return 1 ; /* error */ - /* - * do the actual decoding - */ - new_pkt = my_malloc (k * sizeof (gf * ), "new pkt pointers" ); - for (row = 0 ; row < k ; row++ ) { - if (index[row] >= k) { - new_pkt[row] = my_malloc (sz * sizeof (gf), "new pkt buffer" ); - memset(new_pkt[row], '\0', sz * sizeof(gf) ) ; - for (col = 0 ; col < k ; col++ ) - addmul(new_pkt[row], pkt[col], m_dec[row*k + col], sz) ; - } - } - /* - * move pkts to their final destination - */ - for (row = 0 ; row < k ; row++ ) { - if (index[row] >= k) { - memcpy(pkt[row], new_pkt[row], sz*sizeof(gf)); - free(new_pkt[row]); - } - } - free(new_pkt); - free(m_dec); - - return 0; -} - -/*********** end of FEC code -- beginning of test code ************/ - -#if (TEST || DEBUG) -void -test_gf() -{ - int i ; - /* - * test gf tables. Sufficiently tested... - */ - for (i=0; i<= GF_SIZE; i++) { - if (gf_exp[gf_log[i]] != i) - fprintf(stderr, "bad exp/log i %d log %d exp(log) %d\n", - i, gf_log[i], gf_exp[gf_log[i]]); - - if (i != 0 && gf_mul(i, inverse[i]) != 1) - fprintf(stderr, "bad mul/inv i %d inv %d i*inv(i) %d\n", - i, inverse[i], gf_mul(i, inverse[i]) ); - if (gf_mul(0,i) != 0) - fprintf(stderr, "bad mul table 0,%d\n",i); - if (gf_mul(i,0) != 0) - fprintf(stderr, "bad mul table %d,0\n",i); - } -} -#endif /* TEST */ diff --git a/lib/libcrc32.c b/lib/libcrc32.c new file mode 100644 index 0000000..d556e91 --- /dev/null +++ b/lib/libcrc32.c @@ -0,0 +1,104 @@ +/* + * COPYRIGHT (C) 1986 Gary S. Brown. You may use this program, or + * code or tables extracted from it, as desired without restriction. + * + * First, the polynomial itself and its table of feedback terms. The + * polynomial is + * X^32+X^26+X^23+X^22+X^16+X^12+X^11+X^10+X^8+X^7+X^5+X^4+X^2+X^1+X^0 + * + * Note that we take it "backwards" and put the highest-order term in + * the lowest-order bit. The X^32 term is "implied"; the LSB is the + * X^31 term, etc. The X^0 term (usually shown as "+1") results in + * the MSB being 1 + * + * Note that the usual hardware shift register implementation, which + * is what we're using (we're merely optimizing it by doing eight-bit + * chunks at a time) shifts bits into the lowest-order term. In our + * implementation, that means shifting towards the right. Why do we + * do it this way? Because the calculated CRC must be transmitted in + * order from highest-order term to lowest-order term. UARTs transmit + * characters in order from LSB to MSB. By storing the CRC this way + * we hand it to the UART in the order low-byte to high-byte; the UART + * sends each low-bit to hight-bit; and the result is transmission bit + * by bit from highest- to lowest-order term without requiring any bit + * shuffling on our part. Reception works similarly + * + * The feedback terms table consists of 256, 32-bit entries. Notes + * + * The table can be generated at runtime if desired; code to do so + * is shown later. It might not be obvious, but the feedback + * terms simply represent the results of eight shift/xor opera + * tions for all combinations of data and CRC register values + * + * The values must be right-shifted by eight bits by the "updcrc + * logic; the shift must be unsigned (bring in zeroes). On some + * hardware you could probably optimize the shift in assembler by + * using byte-swap instructions + * polynomial $edb88320 + */ + +#include + +const static uint32_t crc32_table[256] = { + 0x00000000L, 0x77073096L, 0xee0e612cL, 0x990951baL, 0x076dc419L, + 0x706af48fL, 0xe963a535L, 0x9e6495a3L, 0x0edb8832L, 0x79dcb8a4L, + 0xe0d5e91eL, 0x97d2d988L, 0x09b64c2bL, 0x7eb17cbdL, 0xe7b82d07L, + 0x90bf1d91L, 0x1db71064L, 0x6ab020f2L, 0xf3b97148L, 0x84be41deL, + 0x1adad47dL, 0x6ddde4ebL, 0xf4d4b551L, 0x83d385c7L, 0x136c9856L, + 0x646ba8c0L, 0xfd62f97aL, 0x8a65c9ecL, 0x14015c4fL, 0x63066cd9L, + 0xfa0f3d63L, 0x8d080df5L, 0x3b6e20c8L, 0x4c69105eL, 0xd56041e4L, + 0xa2677172L, 0x3c03e4d1L, 0x4b04d447L, 0xd20d85fdL, 0xa50ab56bL, + 0x35b5a8faL, 0x42b2986cL, 0xdbbbc9d6L, 0xacbcf940L, 0x32d86ce3L, + 0x45df5c75L, 0xdcd60dcfL, 0xabd13d59L, 0x26d930acL, 0x51de003aL, + 0xc8d75180L, 0xbfd06116L, 0x21b4f4b5L, 0x56b3c423L, 0xcfba9599L, + 0xb8bda50fL, 0x2802b89eL, 0x5f058808L, 0xc60cd9b2L, 0xb10be924L, + 0x2f6f7c87L, 0x58684c11L, 0xc1611dabL, 0xb6662d3dL, 0x76dc4190L, + 0x01db7106L, 0x98d220bcL, 0xefd5102aL, 0x71b18589L, 0x06b6b51fL, + 0x9fbfe4a5L, 0xe8b8d433L, 0x7807c9a2L, 0x0f00f934L, 0x9609a88eL, + 0xe10e9818L, 0x7f6a0dbbL, 0x086d3d2dL, 0x91646c97L, 0xe6635c01L, + 0x6b6b51f4L, 0x1c6c6162L, 0x856530d8L, 0xf262004eL, 0x6c0695edL, + 0x1b01a57bL, 0x8208f4c1L, 0xf50fc457L, 0x65b0d9c6L, 0x12b7e950L, + 0x8bbeb8eaL, 0xfcb9887cL, 0x62dd1ddfL, 0x15da2d49L, 0x8cd37cf3L, + 0xfbd44c65L, 0x4db26158L, 0x3ab551ceL, 0xa3bc0074L, 0xd4bb30e2L, + 0x4adfa541L, 0x3dd895d7L, 0xa4d1c46dL, 0xd3d6f4fbL, 0x4369e96aL, + 0x346ed9fcL, 0xad678846L, 0xda60b8d0L, 0x44042d73L, 0x33031de5L, + 0xaa0a4c5fL, 0xdd0d7cc9L, 0x5005713cL, 0x270241aaL, 0xbe0b1010L, + 0xc90c2086L, 0x5768b525L, 0x206f85b3L, 0xb966d409L, 0xce61e49fL, + 0x5edef90eL, 0x29d9c998L, 0xb0d09822L, 0xc7d7a8b4L, 0x59b33d17L, + 0x2eb40d81L, 0xb7bd5c3bL, 0xc0ba6cadL, 0xedb88320L, 0x9abfb3b6L, + 0x03b6e20cL, 0x74b1d29aL, 0xead54739L, 0x9dd277afL, 0x04db2615L, + 0x73dc1683L, 0xe3630b12L, 0x94643b84L, 0x0d6d6a3eL, 0x7a6a5aa8L, + 0xe40ecf0bL, 0x9309ff9dL, 0x0a00ae27L, 0x7d079eb1L, 0xf00f9344L, + 0x8708a3d2L, 0x1e01f268L, 0x6906c2feL, 0xf762575dL, 0x806567cbL, + 0x196c3671L, 0x6e6b06e7L, 0xfed41b76L, 0x89d32be0L, 0x10da7a5aL, + 0x67dd4accL, 0xf9b9df6fL, 0x8ebeeff9L, 0x17b7be43L, 0x60b08ed5L, + 0xd6d6a3e8L, 0xa1d1937eL, 0x38d8c2c4L, 0x4fdff252L, 0xd1bb67f1L, + 0xa6bc5767L, 0x3fb506ddL, 0x48b2364bL, 0xd80d2bdaL, 0xaf0a1b4cL, + 0x36034af6L, 0x41047a60L, 0xdf60efc3L, 0xa867df55L, 0x316e8eefL, + 0x4669be79L, 0xcb61b38cL, 0xbc66831aL, 0x256fd2a0L, 0x5268e236L, + 0xcc0c7795L, 0xbb0b4703L, 0x220216b9L, 0x5505262fL, 0xc5ba3bbeL, + 0xb2bd0b28L, 0x2bb45a92L, 0x5cb36a04L, 0xc2d7ffa7L, 0xb5d0cf31L, + 0x2cd99e8bL, 0x5bdeae1dL, 0x9b64c2b0L, 0xec63f226L, 0x756aa39cL, + 0x026d930aL, 0x9c0906a9L, 0xeb0e363fL, 0x72076785L, 0x05005713L, + 0x95bf4a82L, 0xe2b87a14L, 0x7bb12baeL, 0x0cb61b38L, 0x92d28e9bL, + 0xe5d5be0dL, 0x7cdcefb7L, 0x0bdbdf21L, 0x86d3d2d4L, 0xf1d4e242L, + 0x68ddb3f8L, 0x1fda836eL, 0x81be16cdL, 0xf6b9265bL, 0x6fb077e1L, + 0x18b74777L, 0x88085ae6L, 0xff0f6a70L, 0x66063bcaL, 0x11010b5cL, + 0x8f659effL, 0xf862ae69L, 0x616bffd3L, 0x166ccf45L, 0xa00ae278L, + 0xd70dd2eeL, 0x4e048354L, 0x3903b3c2L, 0xa7672661L, 0xd06016f7L, + 0x4969474dL, 0x3e6e77dbL, 0xaed16a4aL, 0xd9d65adcL, 0x40df0b66L, + 0x37d83bf0L, 0xa9bcae53L, 0xdebb9ec5L, 0x47b2cf7fL, 0x30b5ffe9L, + 0xbdbdf21cL, 0xcabac28aL, 0x53b39330L, 0x24b4a3a6L, 0xbad03605L, + 0xcdd70693L, 0x54de5729L, 0x23d967bfL, 0xb3667a2eL, 0xc4614ab8L, + 0x5d681b02L, 0x2a6f2b94L, 0xb40bbe37L, 0xc30c8ea1L, 0x5a05df1bL, + 0x2d02ef8dL +}; + +uint32_t crc32(uint32_t val, const void *ss, int len) +{ + const unsigned char *s = ss; + + while (--len >= 0) + val = crc32_table[(val ^ *s++) & 0xff] ^ (val >> 8); + return val; +} diff --git a/lib/libfec.c b/lib/libfec.c new file mode 100644 index 0000000..6d9020f --- /dev/null +++ b/lib/libfec.c @@ -0,0 +1,904 @@ +/* + * fec.c -- forward error correction based on Vandermonde matrices + * 980624 + * (C) 1997-98 Luigi Rizzo (luigi@iet.unipi.it) + * + * Portions derived from code by Phil Karn (karn@ka9q.ampr.org), + * Robert Morelos-Zaragoza (robert@spectra.eng.hawaii.edu) and Hari + * Thirumoorthy (harit@spectra.eng.hawaii.edu), Aug 1995 + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above + * copyright notice, this list of conditions and the following + * disclaimer in the documentation and/or other materials + * provided with the distribution. + * + * THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND + * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, + * THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A + * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS + * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, + * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, + * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY + * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR + * TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT + * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY + * OF SUCH DAMAGE. + */ + +/* + * The following parameter defines how many bits are used for + * field elements. The code supports any value from 2 to 16 + * but fastest operation is achieved with 8 bit elements + * This is the only parameter you may want to change. + */ +#ifndef GF_BITS +#define GF_BITS 8 /* code over GF(2**GF_BITS) - change to suit */ +#endif + +#include +#include +#include + +/* + * stuff used for testing purposes only + */ + +#ifdef TEST +#define DEB(x) +#define DDB(x) x +#define DEBUG 0 /* minimal debugging */ +#ifdef MSDOS +#include +struct timeval { + unsigned long ticks; +}; +#define gettimeofday(x, dummy) { (x)->ticks = clock() ; } +#define DIFF_T(a,b) (1+ 1000000*(a.ticks - b.ticks) / CLOCKS_PER_SEC ) +typedef unsigned long u_long ; +typedef unsigned short u_short ; +#else /* typically, unix systems */ +#include +#define DIFF_T(a,b) \ + (1+ 1000000*(a.tv_sec - b.tv_sec) + (a.tv_usec - b.tv_usec) ) +#endif + +#define TICK(t) \ + {struct timeval x ; \ + gettimeofday(&x, NULL) ; \ + t = x.tv_usec + 1000000* (x.tv_sec & 0xff ) ; \ + } +#define TOCK(t) \ + { u_long t1 ; TICK(t1) ; \ + if (t1 < t) t = 256000000 + t1 - t ; \ + else t = t1 - t ; \ + if (t == 0) t = 1 ;} + +u_long ticks[10]; /* vars for timekeeping */ +#else +#define DEB(x) +#define DDB(x) +#define TICK(x) +#define TOCK(x) +#endif /* TEST */ + +/* + * You should not need to change anything beyond this point. + * The first part of the file implements linear algebra in GF. + * + * gf is the type used to store an element of the Galois Field. + * Must constain at least GF_BITS bits. + * + * Note: unsigned char will work up to GF(256) but int seems to run + * faster on the Pentium. We use int whenever have to deal with an + * index, since they are generally faster. + */ +#if (GF_BITS < 2 && GF_BITS >16) +#error "GF_BITS must be 2 .. 16" +#endif +#if (GF_BITS <= 8) +typedef unsigned char gf; +#else +typedef unsigned short gf; +#endif + +#define GF_SIZE ((1 << GF_BITS) - 1) /* powers of \alpha */ + +/* + * Primitive polynomials - see Lin & Costello, Appendix A, + * and Lee & Messerschmitt, p. 453. + */ +static char *allPp[] = { /* GF_BITS polynomial */ + NULL, /* 0 no code */ + NULL, /* 1 no code */ + "111", /* 2 1+x+x^2 */ + "1101", /* 3 1+x+x^3 */ + "11001", /* 4 1+x+x^4 */ + "101001", /* 5 1+x^2+x^5 */ + "1100001", /* 6 1+x+x^6 */ + "10010001", /* 7 1 + x^3 + x^7 */ + "101110001", /* 8 1+x^2+x^3+x^4+x^8 */ + "1000100001", /* 9 1+x^4+x^9 */ + "10010000001", /* 10 1+x^3+x^10 */ + "101000000001", /* 11 1+x^2+x^11 */ + "1100101000001", /* 12 1+x+x^4+x^6+x^12 */ + "11011000000001", /* 13 1+x+x^3+x^4+x^13 */ + "110000100010001", /* 14 1+x+x^6+x^10+x^14 */ + "1100000000000001", /* 15 1+x+x^15 */ + "11010000000010001" /* 16 1+x+x^3+x^12+x^16 */ +}; + + +/* + * To speed up computations, we have tables for logarithm, exponent + * and inverse of a number. If GF_BITS <= 8, we use a table for + * multiplication as well (it takes 64K, no big deal even on a PDA, + * especially because it can be pre-initialized an put into a ROM!), + * otherwhise we use a table of logarithms. + * In any case the macro gf_mul(x,y) takes care of multiplications. + */ + +static gf gf_exp[2*GF_SIZE]; /* index->poly form conversion table */ +static int gf_log[GF_SIZE + 1]; /* Poly->index form conversion table */ +static gf inverse[GF_SIZE+1]; /* inverse of field elem. */ + /* inv[\alpha**i]=\alpha**(GF_SIZE-i-1) */ + +/* + * modnn(x) computes x % GF_SIZE, where GF_SIZE is 2**GF_BITS - 1, + * without a slow divide. + */ +static inline gf +modnn(int x) +{ + while (x >= GF_SIZE) { + x -= GF_SIZE; + x = (x >> GF_BITS) + (x & GF_SIZE); + } + return x; +} + +#define SWAP(a,b,t) {t tmp; tmp=a; a=b; b=tmp;} + +/* + * gf_mul(x,y) multiplies two numbers. If GF_BITS<=8, it is much + * faster to use a multiplication table. + * + * USE_GF_MULC, GF_MULC0(c) and GF_ADDMULC(x) can be used when multiplying + * many numbers by the same constant. In this case the first + * call sets the constant, and others perform the multiplications. + * A value related to the multiplication is held in a local variable + * declared with USE_GF_MULC . See usage in addmul1(). + */ +#if (GF_BITS <= 8) +static gf gf_mul_table[GF_SIZE + 1][GF_SIZE + 1]; + +#define gf_mul(x,y) gf_mul_table[x][y] + +#define USE_GF_MULC register gf * __gf_mulc_ +#define GF_MULC0(c) __gf_mulc_ = gf_mul_table[c] +#define GF_ADDMULC(dst, x) dst ^= __gf_mulc_[x] + +static void +init_mul_table() +{ + int i, j; + for (i=0; i< GF_SIZE+1; i++) + for (j=0; j< GF_SIZE+1; j++) + gf_mul_table[i][j] = gf_exp[modnn(gf_log[i] + gf_log[j]) ] ; + + for (j=0; j< GF_SIZE+1; j++) + gf_mul_table[0][j] = gf_mul_table[j][0] = 0; +} +#else /* GF_BITS > 8 */ +static inline gf +gf_mul(x,y) +{ + if ( (x) == 0 || (y)==0 ) return 0; + + return gf_exp[gf_log[x] + gf_log[y] ] ; +} +#define init_mul_table() + +#define USE_GF_MULC register gf * __gf_mulc_ +#define GF_MULC0(c) __gf_mulc_ = &gf_exp[ gf_log[c] ] +#define GF_ADDMULC(dst, x) { if (x) dst ^= __gf_mulc_[ gf_log[x] ] ; } +#endif + +/* + * Generate GF(2**m) from the irreducible polynomial p(X) in p[0]..p[m] + * Lookup tables: + * index->polynomial form gf_exp[] contains j= \alpha^i; + * polynomial form -> index form gf_log[ j = \alpha^i ] = i + * \alpha=x is the primitive element of GF(2^m) + * + * For efficiency, gf_exp[] has size 2*GF_SIZE, so that a simple + * multiplication of two numbers can be resolved without calling modnn + */ + +/* + * i use malloc so many times, it is easier to put checks all in + * one place. + */ +static void * +my_malloc(int sz, char *err_string) +{ + void *p = malloc( sz ); + if (p == NULL) { + fprintf(stderr, "-- malloc failure allocating %s\n", err_string); + exit(1) ; + } + return p ; +} + +#define NEW_GF_MATRIX(rows, cols) \ + (gf *)my_malloc(rows * cols * sizeof(gf), " ## __LINE__ ## " ) + +/* + * initialize the data structures used for computations in GF. + */ +static void +generate_gf(void) +{ + int i; + gf mask; + char *Pp = allPp[GF_BITS] ; + + mask = 1; /* x ** 0 = 1 */ + gf_exp[GF_BITS] = 0; /* will be updated at the end of the 1st loop */ + /* + * first, generate the (polynomial representation of) powers of \alpha, + * which are stored in gf_exp[i] = \alpha ** i . + * At the same time build gf_log[gf_exp[i]] = i . + * The first GF_BITS powers are simply bits shifted to the left. + */ + for (i = 0; i < GF_BITS; i++, mask <<= 1 ) { + gf_exp[i] = mask; + gf_log[gf_exp[i]] = i; + /* + * If Pp[i] == 1 then \alpha ** i occurs in poly-repr + * gf_exp[GF_BITS] = \alpha ** GF_BITS + */ + if ( Pp[i] == '1' ) + gf_exp[GF_BITS] ^= mask; + } + /* + * now gf_exp[GF_BITS] = \alpha ** GF_BITS is complete, so can als + * compute its inverse. + */ + gf_log[gf_exp[GF_BITS]] = GF_BITS; + /* + * Poly-repr of \alpha ** (i+1) is given by poly-repr of + * \alpha ** i shifted left one-bit and accounting for any + * \alpha ** GF_BITS term that may occur when poly-repr of + * \alpha ** i is shifted. + */ + mask = 1 << (GF_BITS - 1 ) ; + for (i = GF_BITS + 1; i < GF_SIZE; i++) { + if (gf_exp[i - 1] >= mask) + gf_exp[i] = gf_exp[GF_BITS] ^ ((gf_exp[i - 1] ^ mask) << 1); + else + gf_exp[i] = gf_exp[i - 1] << 1; + gf_log[gf_exp[i]] = i; + } + /* + * log(0) is not defined, so use a special value + */ + gf_log[0] = GF_SIZE ; + /* set the extended gf_exp values for fast multiply */ + for (i = 0 ; i < GF_SIZE ; i++) + gf_exp[i + GF_SIZE] = gf_exp[i] ; + + /* + * again special cases. 0 has no inverse. This used to + * be initialized to GF_SIZE, but it should make no difference + * since noone is supposed to read from here. + */ + inverse[0] = 0 ; + inverse[1] = 1; + for (i=2; i<=GF_SIZE; i++) + inverse[i] = gf_exp[GF_SIZE-gf_log[i]]; +} + +/* + * Various linear algebra operations that i use often. + */ + +/* + * addmul() computes dst[] = dst[] + c * src[] + * This is used often, so better optimize it! Currently the loop is + * unrolled 16 times, a good value for 486 and pentium-class machines. + * The case c=0 is also optimized, whereas c=1 is not. These + * calls are unfrequent in my typical apps so I did not bother. + * + * Note that gcc on + */ +#define addmul(dst, src, c, sz) \ + if (c != 0) addmul1(dst, src, c, sz) + +#define UNROLL 16 /* 1, 4, 8, 16 */ +static void +addmul1(gf *dst1, gf *src1, gf c, int sz) +{ + USE_GF_MULC ; + register gf *dst = dst1, *src = src1 ; + gf *lim = &dst[sz - UNROLL + 1] ; + + GF_MULC0(c) ; + +#if (UNROLL > 1) /* unrolling by 8/16 is quite effective on the pentium */ + for (; dst < lim ; dst += UNROLL, src += UNROLL ) { + GF_ADDMULC( dst[0] , src[0] ); + GF_ADDMULC( dst[1] , src[1] ); + GF_ADDMULC( dst[2] , src[2] ); + GF_ADDMULC( dst[3] , src[3] ); +#if (UNROLL > 4) + GF_ADDMULC( dst[4] , src[4] ); + GF_ADDMULC( dst[5] , src[5] ); + GF_ADDMULC( dst[6] , src[6] ); + GF_ADDMULC( dst[7] , src[7] ); +#endif +#if (UNROLL > 8) + GF_ADDMULC( dst[8] , src[8] ); + GF_ADDMULC( dst[9] , src[9] ); + GF_ADDMULC( dst[10] , src[10] ); + GF_ADDMULC( dst[11] , src[11] ); + GF_ADDMULC( dst[12] , src[12] ); + GF_ADDMULC( dst[13] , src[13] ); + GF_ADDMULC( dst[14] , src[14] ); + GF_ADDMULC( dst[15] , src[15] ); +#endif + } +#endif + lim += UNROLL - 1 ; + for (; dst < lim; dst++, src++ ) /* final components */ + GF_ADDMULC( *dst , *src ); +} + +/* + * computes C = AB where A is n*k, B is k*m, C is n*m + */ +static void +matmul(gf *a, gf *b, gf *c, int n, int k, int m) +{ + int row, col, i ; + + for (row = 0; row < n ; row++) { + for (col = 0; col < m ; col++) { + gf *pa = &a[ row * k ]; + gf *pb = &b[ col ]; + gf acc = 0 ; + for (i = 0; i < k ; i++, pa++, pb += m ) + acc ^= gf_mul( *pa, *pb ) ; + c[ row * m + col ] = acc ; + } + } +} + +#ifdef DEBUG +/* + * returns 1 if the square matrix is identiy + * (only for test) + */ +static int +is_identity(gf *m, int k) +{ + int row, col ; + for (row=0; row 1) { + fprintf(stderr, "singular matrix\n"); + goto fail ; + } + } + } + } + if (icol == -1) { + fprintf(stderr, "XXX pivot not found!\n"); + goto fail ; + } +found_piv: + ++(ipiv[icol]) ; + /* + * swap rows irow and icol, so afterwards the diagonal + * element will be correct. Rarely done, not worth + * optimizing. + */ + if (irow != icol) { + for (ix = 0 ; ix < k ; ix++ ) { + SWAP( src[irow*k + ix], src[icol*k + ix], gf) ; + } + } + indxr[col] = irow ; + indxc[col] = icol ; + pivot_row = &src[icol*k] ; + c = pivot_row[icol] ; + if (c == 0) { + fprintf(stderr, "singular matrix 2\n"); + goto fail ; + } + if (c != 1 ) { /* otherwhise this is a NOP */ + /* + * this is done often , but optimizing is not so + * fruitful, at least in the obvious ways (unrolling) + */ + DEB( pivswaps++ ; ) + c = inverse[ c ] ; + pivot_row[icol] = 1 ; + for (ix = 0 ; ix < k ; ix++ ) + pivot_row[ix] = gf_mul(c, pivot_row[ix] ); + } + /* + * from all rows, remove multiples of the selected row + * to zero the relevant entry (in fact, the entry is not zero + * because we know it must be zero). + * (Here, if we know that the pivot_row is the identity, + * we can optimize the addmul). + */ + id_row[icol] = 1; + if (memcmp(pivot_row, id_row, k*sizeof(gf)) != 0) { + for (p = src, ix = 0 ; ix < k ; ix++, p += k ) { + if (ix != icol) { + c = p[icol] ; + p[icol] = 0 ; + addmul(p, pivot_row, c, k ); + } + } + } + id_row[icol] = 0; + } /* done all columns */ + for (col = k-1 ; col >= 0 ; col-- ) { + if (indxr[col] <0 || indxr[col] >= k) + fprintf(stderr, "AARGH, indxr[col] %d\n", indxr[col]); + else if (indxc[col] <0 || indxc[col] >= k) + fprintf(stderr, "AARGH, indxc[col] %d\n", indxc[col]); + else + if (indxr[col] != indxc[col] ) { + for (row = 0 ; row < k ; row++ ) { + SWAP( src[row*k + indxr[col]], src[row*k + indxc[col]], gf) ; + } + } + } + error = 0 ; +fail: + free(indxc); + free(indxr); + free(ipiv); + free(id_row); + free(temp_row); + return error ; +} + +/* + * fast code for inverting a vandermonde matrix. + * XXX NOTE: It assumes that the matrix + * is not singular and _IS_ a vandermonde matrix. Only uses + * the second column of the matrix, containing the p_i's. + * + * Algorithm borrowed from "Numerical recipes in C" -- sec.2.8, but + * largely revised for my purposes. + * p = coefficients of the matrix (p_i) + * q = values of the polynomial (known) + */ + +int +invert_vdm(gf *src, int k) +{ + int i, j, row, col ; + gf *b, *c, *p; + gf t, xx ; + + if (k == 1) /* degenerate case, matrix must be p^0 = 1 */ + return 0 ; + /* + * c holds the coefficient of P(x) = Prod (x - p_i), i=0..k-1 + * b holds the coefficient for the matrix inversion + */ + c = NEW_GF_MATRIX(1, k); + b = NEW_GF_MATRIX(1, k); + + p = NEW_GF_MATRIX(1, k); + + for ( j=1, i = 0 ; i < k ; i++, j+=k ) { + c[i] = 0 ; + p[i] = src[j] ; /* p[i] */ + } + /* + * construct coeffs. recursively. We know c[k] = 1 (implicit) + * and start P_0 = x - p_0, then at each stage multiply by + * x - p_i generating P_i = x P_{i-1} - p_i P_{i-1} + * After k steps we are done. + */ + c[k-1] = p[0] ; /* really -p(0), but x = -x in GF(2^m) */ + for (i = 1 ; i < k ; i++ ) { + gf p_i = p[i] ; /* see above comment */ + for (j = k-1 - ( i - 1 ) ; j < k-1 ; j++ ) + c[j] ^= gf_mul( p_i, c[j+1] ) ; + c[k-1] ^= p_i ; + } + + for (row = 0 ; row < k ; row++ ) { + /* + * synthetic division etc. + */ + xx = p[row] ; + t = 1 ; + b[k-1] = 1 ; /* this is in fact c[k] */ + for (i = k-2 ; i >= 0 ; i-- ) { + b[i] = c[i+1] ^ gf_mul(xx, b[i+1]) ; + t = gf_mul(xx, t) ^ b[i] ; + } + for (col = 0 ; col < k ; col++ ) + src[col*k + row] = gf_mul(inverse[t], b[col] ); + } + free(c) ; + free(b) ; + free(p) ; + return 0 ; +} + +static int fec_initialized = 0 ; +static void +init_fec() +{ + TICK(ticks[0]); + generate_gf(); + TOCK(ticks[0]); + DDB(fprintf(stderr, "generate_gf took %ldus\n", ticks[0]);) + TICK(ticks[0]); + init_mul_table(); + TOCK(ticks[0]); + DDB(fprintf(stderr, "init_mul_table took %ldus\n", ticks[0]);) + fec_initialized = 1 ; +} + +/* + * This section contains the proper FEC encoding/decoding routines. + * The encoding matrix is computed starting with a Vandermonde matrix, + * and then transforming it into a systematic matrix. + */ + +#define FEC_MAGIC 0xFECC0DEC + +struct fec_parms { + u_long magic ; + int k, n ; /* parameters of the code */ + gf *enc_matrix ; +} ; + +void +fec_free(struct fec_parms *p) +{ + if (p==NULL || + p->magic != ( ( (FEC_MAGIC ^ p->k) ^ p->n) ^ (int)(p->enc_matrix)) ) { + fprintf(stderr, "bad parameters to fec_free\n"); + return ; + } + free(p->enc_matrix); + free(p); +} + +/* + * create a new encoder, returning a descriptor. This contains k,n and + * the encoding matrix. + */ +struct fec_parms * +fec_new(int k, int n) +{ + int row, col ; + gf *p, *tmp_m ; + + struct fec_parms *retval ; + + if (fec_initialized == 0) + init_fec(); + + if (k > GF_SIZE + 1 || n > GF_SIZE + 1 || k > n ) { + fprintf(stderr, "Invalid parameters k %d n %d GF_SIZE %d\n", + k, n, GF_SIZE ); + return NULL ; + } + retval = my_malloc(sizeof(struct fec_parms), "new_code"); + retval->k = k ; + retval->n = n ; + retval->enc_matrix = NEW_GF_MATRIX(n, k); + retval->magic = ( ( FEC_MAGIC ^ k) ^ n) ^ (int)(retval->enc_matrix) ; + tmp_m = NEW_GF_MATRIX(n, k); + /* + * fill the matrix with powers of field elements, starting from 0. + * The first row is special, cannot be computed with exp. table. + */ + tmp_m[0] = 1 ; + for (col = 1; col < k ; col++) + tmp_m[col] = 0 ; + for (p = tmp_m + k, row = 0; row < n-1 ; row++, p += k) { + for ( col = 0 ; col < k ; col ++ ) + p[col] = gf_exp[modnn(row*col)]; + } + + /* + * quick code to build systematic matrix: invert the top + * k*k vandermonde matrix, multiply right the bottom n-k rows + * by the inverse, and construct the identity matrix at the top. + */ + TICK(ticks[3]); + invert_vdm(tmp_m, k); /* much faster than invert_mat */ + matmul(tmp_m + k*k, tmp_m, retval->enc_matrix + k*k, n - k, k, k); + /* + * the upper matrix is I so do not bother with a slow multiply + */ + memset(retval->enc_matrix, '\0', k*k*sizeof(gf) ); + for (p = retval->enc_matrix, col = 0 ; col < k ; col++, p += k+1 ) + *p = 1 ; + free(tmp_m); + TOCK(ticks[3]); + + DDB(fprintf(stderr, "--- %ld us to build encoding matrix\n", + ticks[3]);) + DEB(pr_matrix(retval->enc_matrix, n, k, "encoding_matrix");) + return retval ; +} + +/* + * fec_encode accepts as input pointers to n data packets of size sz, + * and produces as output a packet pointed to by fec, computed + * with index "index". + */ +void +fec_encode(struct fec_parms *code, gf *src[], gf *fec, int index, int sz) +{ + int i, k = code->k ; + gf *p ; + + if (GF_BITS > 8) + sz /= 2 ; + + if (index < k) + memcpy(fec, src[index], sz*sizeof(gf) ) ; + else if (index < code->n) { + p = &(code->enc_matrix[index*k] ); + memset(fec, '\0', sz*sizeof(gf)); + for (i = 0; i < k ; i++) + addmul(fec, src[i], p[i], sz ) ; + } else + fprintf(stderr, "Invalid index %d (max %d)\n", + index, code->n - 1 ); +} + +void fec_encode_linear(struct fec_parms *code, gf *src, gf *fec, int index, int sz) +{ + int i, k = code->k ; + gf *p ; + + if (GF_BITS > 8) + sz /= 2 ; + + if (index < k) + memcpy(fec, src + (index * sz), sz*sizeof(gf) ) ; + else if (index < code->n) { + p = &(code->enc_matrix[index*k] ); + memset(fec, '\0', sz*sizeof(gf)); + for (i = 0; i < k ; i++) + addmul(fec, src + (i * sz), p[i], sz ) ; + } else + fprintf(stderr, "Invalid index %d (max %d)\n", + index, code->n - 1 ); +} +/* + * shuffle move src packets in their position + */ +static int +shuffle(gf *pkt[], int index[], int k) +{ + int i; + + for ( i = 0 ; i < k ; ) { + if (index[i] >= k || index[i] == i) + i++ ; + else { + /* + * put pkt in the right position (first check for conflicts). + */ + int c = index[i] ; + + if (index[c] == c) { + DEB(fprintf(stderr, "\nshuffle, error at %d\n", i);) + return 1 ; + } + SWAP(index[i], index[c], int) ; + SWAP(pkt[i], pkt[c], gf *) ; + } + } + DEB( /* just test that it works... */ + for ( i = 0 ; i < k ; i++ ) { + if (index[i] < k && index[i] != i) { + fprintf(stderr, "shuffle: after\n"); + for (i=0; ik ; + gf *p, *matrix = NEW_GF_MATRIX(k, k); + + TICK(ticks[9]); + for (i = 0, p = matrix ; i < k ; i++, p += k ) { +#if 1 /* this is simply an optimization, not very useful indeed */ + if (index[i] < k) { + memset(p, '\0', k*sizeof(gf) ); + p[i] = 1 ; + } else +#endif + if (index[i] < code->n ) + memcpy(p, &(code->enc_matrix[index[i]*k]), k*sizeof(gf) ); + else { + fprintf(stderr, "decode: invalid index %d (max %d)\n", + index[i], code->n - 1 ); + free(matrix) ; + return NULL ; + } + } + TICK(ticks[9]); + if (invert_mat(matrix, k)) { + free(matrix); + matrix = NULL ; + } + TOCK(ticks[9]); + return matrix ; +} + +/* + * fec_decode receives as input a vector of packets, the indexes of + * packets, and produces the correct vector as output. + * + * Input: + * code: pointer to code descriptor + * pkt: pointers to received packets. They are modified + * to store the output packets (in place) + * index: pointer to packet indexes (modified) + * sz: size of each packet + */ +int +fec_decode(struct fec_parms *code, gf *pkt[], int index[], int sz) +{ + gf *m_dec ; + gf **new_pkt ; + int row, col , k = code->k ; + + if (GF_BITS > 8) + sz /= 2 ; + + if (shuffle(pkt, index, k)) /* error if true */ + return 1 ; + m_dec = build_decode_matrix(code, pkt, index); + + if (m_dec == NULL) + return 1 ; /* error */ + /* + * do the actual decoding + */ + new_pkt = my_malloc (k * sizeof (gf * ), "new pkt pointers" ); + for (row = 0 ; row < k ; row++ ) { + if (index[row] >= k) { + new_pkt[row] = my_malloc (sz * sizeof (gf), "new pkt buffer" ); + memset(new_pkt[row], '\0', sz * sizeof(gf) ) ; + for (col = 0 ; col < k ; col++ ) + addmul(new_pkt[row], pkt[col], m_dec[row*k + col], sz) ; + } + } + /* + * move pkts to their final destination + */ + for (row = 0 ; row < k ; row++ ) { + if (index[row] >= k) { + memcpy(pkt[row], new_pkt[row], sz*sizeof(gf)); + free(new_pkt[row]); + } + } + free(new_pkt); + free(m_dec); + + return 0; +} + +/*********** end of FEC code -- beginning of test code ************/ + +#if (TEST || DEBUG) +void +test_gf() +{ + int i ; + /* + * test gf tables. Sufficiently tested... + */ + for (i=0; i<= GF_SIZE; i++) { + if (gf_exp[gf_log[i]] != i) + fprintf(stderr, "bad exp/log i %d log %d exp(log) %d\n", + i, gf_log[i], gf_exp[gf_log[i]]); + + if (i != 0 && gf_mul(i, inverse[i]) != 1) + fprintf(stderr, "bad mul/inv i %d inv %d i*inv(i) %d\n", + i, inverse[i], gf_mul(i, inverse[i]) ); + if (gf_mul(0,i) != 0) + fprintf(stderr, "bad mul table 0,%d\n",i); + if (gf_mul(i,0) != 0) + fprintf(stderr, "bad mul table %d,0\n",i); + } +} +#endif /* TEST */ -- cgit v1.2.3