Author: nextgens
Date: 2009-03-30 13:57:14 +0000 (Mon, 30 Mar 2009)
New Revision: 26259

Added:
   trunk/contrib/fec/src/
   trunk/contrib/fec/src/csrc/
   trunk/contrib/fec/src/csrc/Makefile
   trunk/contrib/fec/src/csrc/README
   trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native16Code.h
   trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native8Code.h
   trunk/contrib/fec/src/csrc/fec.3
   trunk/contrib/fec/src/csrc/fec.c
   trunk/contrib/fec/src/csrc/fec.h
   trunk/contrib/fec/src/csrc/fec16-jinterf.c
   trunk/contrib/fec/src/csrc/fec8-jinterf.c
   trunk/contrib/fec/src/csrc/fec_win32.c
   trunk/contrib/fec/src/csrc/test.c
   trunk/contrib/fec/src/csrc/w32
Log:
contrib: unzip the FEC C files so that people can work on it

Added: trunk/contrib/fec/src/csrc/Makefile
===================================================================
--- trunk/contrib/fec/src/csrc/Makefile                         (rev 0)
+++ trunk/contrib/fec/src/csrc/Makefile 2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,62 @@
+#
+# makefile for vdm.
+#
+# fec.S.980624a is an optimized version for use on 486 and old pentium
+# machines. It is only for GF_BITS=8 and generally does not work
+# fast on systems with multiple instruction pipelines (PentiumPro,
+# Pentium2)
+# same for fec.S16.980624a , for use with GF_BITS=16
+#
+# gcc does something strange, so check the various opt. levels for
+# best performance (or write addmul1 in assembly code).
+#
+# Standard compilation with -O9 works well for PentiumPro and Pentium2
+# machines.
+#
+
+CC=gcc
+# COPT= -O9 -funroll-loops
+COPT= -O1 
+CFLAGS=$(COPT) -Wall # -DTEST
+SRCS= fec.c Makefile test.c fec.s.980621e \
+       fec.S.980624a \
+       fec.S16.980624a
+DOCS= README fec.3
+ALLSRCS= $(SRCS) $(DOCS) fec.h
+
+fec: libfec8.so libfec16.so test.o
+       $(CC) $(CFLAGS) -o fec fec8.o test.o
+
+libfec8.so: fec8.o
+       $(CC) $(CFLAGS) -DGF_BITS=8 -c -I$(JAVA_HOME)/include/ \
+               -I$(JAVA_HOME)/include/linux fec8-jinterf.c \
+               -o fec8-jinterf.o ; \
+       mkdir -p ../../lib/fec-linux-x86/lib/linux/x86 ; \
+        ld -shared fec8-jinterf.o fec8.o -o \
+               ../../lib/fec-linux-x86/lib/linux/x86/libfec8.so
+
+fec8.o: fec.h fec8.S
+       $(CC) $(CFLAGS) -DGF_BITS=8 -c -o fec8.o fec8.S
+
+fec8.S: fec.c Makefile
+       $(CC) $(CFLAGS) -DGF_BITS=8 -S -o fec8.S fec.c
+
+libfec16.so: fec16.o
+       $(CC) $(CFLAGS) -DGF_BITS=16 -c -I$(JAVA_HOME)/include/ \
+               -I$(JAVA_HOME)/include/linux fec16-jinterf.c \
+               -o fec16-jinterf.o ; \
+       mkdir -p ../../lib/fec-linux-x86/lib/linux/x86 ; \
+        ld -shared fec16-jinterf.o fec16.o -o \
+               ../../lib/fec-linux-x86/lib/linux/x86/libfec16.so
+
+fec16.o: fec.h fec16.S
+       $(CC) $(CFLAGS) -DGF_BITS=16 -c -o fec16.o fec16.S
+
+fec16.S: fec.c Makefile
+       $(CC) $(CFLAGS) -DGF_BITS=16 -S -o fec16.S fec.c
+
+clean:
+       - rm *.o *.S fec
+
+tgz: $(ALLSRCS)
+       tar cvzf vdm`date +%y%m%d`.tgz $(ALLSRCS)

Added: trunk/contrib/fec/src/csrc/README
===================================================================
--- trunk/contrib/fec/src/csrc/README                           (rev 0)
+++ trunk/contrib/fec/src/csrc/README   2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,65 @@
+--------    Erasure codes based on Vandermonde matrices    ---------
+--------  (C) 1996-1998  Luigi Rizzo (luigi at iet.unipi.it)  ---------
+
+See fec.c for other contributors.
+
+fec.c contains an implementation of an encoder/decoder for an erasure
+code based on Vandermonde matrices computed over GF(2^m), m=2..16
+
+PRINCIPLE OF OPERATION
+
+The encoded data is computer as
+ 
+       y = E x
+
+where x is a k-vector with source data, y is an n-vector with the
+redundant info, and E is an n*k matrix derived from a Vandermonde
+matrix. The code is systematic.
+
+At the receiver, any subset y' of k elements from y allows the
+reconstruction of the whole x by solving the system
+
+       y' = E' x
+
+where E' is made of rows from E corresponding to the received
+elements.
+
+The complexity of matrix inversion is O(k*l^2) where l is the number
+of elements not in x available at the receiver. This might seem
+large, but data elements are in fact be packets of large size, so
+the inversion cost can be amortized over the size of the packet.
+For practical applications (k and l as large as 30, packet sizes
+of 1KB) the cost can be neglected.
+
+In addition, each of the l lost data packets has a reconstruction
+cost O(k), (obviously) similar to the cost of the encoding phase.
+Being the code systematic, you can express encoding and decoding costs
+roughly as
+
+    Encoding speed:    (c_e/l) MB/s
+    Decoding speed:    (c_d/l) MB/s
+       
+
+PERFORMANCE
+
+The values of c_d and c_e (similar) are very sensitive to issues
+such as cache hit rate, memory speed, field size (8 or 16 bits),
+etc.  Also some machines are better than others in working with
+bytes or 16-bit words.  With the June'98 implementation I have
+measured the following values for c_e and c_d (8-bit version; 16-bit
+version has a penalty between 3 and 4).
+
+    Hardware           C version       Optimized version(*)
+
+    PentiumII 266      62              33
+    PentiumPRO 200     56              30
+    Pentium133         14.5            18
+    486dx2/66           4.05            4.55
+
+(*) The 'optimized' version has some manual optimizations of the assembly
+    code generated by the compiler. It is generally faster for machines
+    with a single instruction pipeline, and generally slower for
+    machines with multiple pipelines.
+
+See the manpage for detailed usage information.
+

Added: trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native16Code.h
===================================================================
--- trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native16Code.h             
                (rev 0)
+++ trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native16Code.h     
2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,46 @@
+/* DO NOT EDIT THIS FILE - it is machine generated */
+#include <jni.h>
+/* Header for class com_onionnetworks_fec_Native16Code */
+
+#ifndef _Included_com_onionnetworks_fec_Native16Code
+#define _Included_com_onionnetworks_fec_Native16Code
+#ifdef __cplusplus
+extern "C" {
+#endif
+/* Inaccessible static: initialized */
+/*
+ * Class:     com_onionnetworks_fec_Native16Code
+ * Method:    nativeEncode
+ * Signature: (I[[B[I[I[[B[III)V
+ */
+JNIEXPORT void JNICALL Java_com_onionnetworks_fec_Native16Code_nativeEncode
+  (JNIEnv *, jobject, jint, jobjectArray, jintArray, jintArray, jobjectArray, 
jintArray, jint, jint);
+
+/*
+ * Class:     com_onionnetworks_fec_Native16Code
+ * Method:    nativeDecode
+ * Signature: (I[[B[I[III)V
+ */
+JNIEXPORT void JNICALL Java_com_onionnetworks_fec_Native16Code_nativeDecode
+  (JNIEnv *, jobject, jint, jobjectArray, jintArray, jintArray, jint, jint);
+
+/*
+ * Class:     com_onionnetworks_fec_Native16Code
+ * Method:    nativeNewFEC
+ * Signature: (II)I
+ */
+JNIEXPORT jint JNICALL Java_com_onionnetworks_fec_Native16Code_nativeNewFEC
+  (JNIEnv *, jobject, jint, jint);
+
+/*
+ * Class:     com_onionnetworks_fec_Native16Code
+ * Method:    nativeFreeFEC
+ * Signature: (I)V
+ */
+JNIEXPORT void JNICALL Java_com_onionnetworks_fec_Native16Code_nativeFreeFEC
+  (JNIEnv *, jobject, jint);
+
+#ifdef __cplusplus
+}
+#endif
+#endif

Added: trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native8Code.h
===================================================================
--- trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native8Code.h              
                (rev 0)
+++ trunk/contrib/fec/src/csrc/com_onionnetworks_fec_Native8Code.h      
2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,46 @@
+/* DO NOT EDIT THIS FILE - it is machine generated */
+#include <jni.h>
+/* Header for class com_onionnetworks_fec_Native8Code */
+
+#ifndef _Included_com_onionnetworks_fec_Native8Code
+#define _Included_com_onionnetworks_fec_Native8Code
+#ifdef __cplusplus
+extern "C" {
+#endif
+/* Inaccessible static: initialized */
+/*
+ * Class:     com_onionnetworks_fec_Native8Code
+ * Method:    nativeEncode
+ * Signature: (I[[B[I[I[[B[III)V
+ */
+JNIEXPORT void JNICALL Java_com_onionnetworks_fec_Native8Code_nativeEncode
+  (JNIEnv *, jobject, jint, jobjectArray, jintArray, jintArray, jobjectArray, 
jintArray, jint, jint);
+
+/*
+ * Class:     com_onionnetworks_fec_Native8Code
+ * Method:    nativeDecode
+ * Signature: (I[[B[I[III)V
+ */
+JNIEXPORT void JNICALL Java_com_onionnetworks_fec_Native8Code_nativeDecode
+  (JNIEnv *, jobject, jint, jobjectArray, jintArray, jintArray, jint, jint);
+
+/*
+ * Class:     com_onionnetworks_fec_Native8Code
+ * Method:    nativeNewFEC
+ * Signature: (II)I
+ */
+JNIEXPORT jint JNICALL Java_com_onionnetworks_fec_Native8Code_nativeNewFEC
+  (JNIEnv *, jobject, jint, jint);
+
+/*
+ * Class:     com_onionnetworks_fec_Native8Code
+ * Method:    nativeFreeFEC
+ * Signature: (I)V
+ */
+JNIEXPORT void JNICALL Java_com_onionnetworks_fec_Native8Code_nativeFreeFEC
+  (JNIEnv *, jobject, jint);
+
+#ifdef __cplusplus
+}
+#endif
+#endif

Added: trunk/contrib/fec/src/csrc/fec.3
===================================================================
--- trunk/contrib/fec/src/csrc/fec.3                            (rev 0)
+++ trunk/contrib/fec/src/csrc/fec.3    2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,112 @@
+.\"
+.\" Copyright 1998 by Luigi Rizzo, Dip. Ingegneria dell'Informazione,
+.\" Universitaet Berlin.  See the source code for copyright details.
+.\" THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE.
+.\"
+.Dd July 15, 1998
+.Dt FEC 3
+.Os
+.Sh NAME
+.Nm fec_new, fec_encode, fec_encode, fec_free
+.Nd An erasure code in GF(2^m)
+.Sh SYNOPSIS
+.Fd #include <fec.h>
+.Ft void *
+.Fn fec_new "int k" "int n"
+.Ft void
+.Fn fec_encode "void *code" "void *data[]" "void *dst" "int i" "int sz"
+.Ft int
+.Fn fec_decode "void *code" "void *data[]" "int i[]" "int sz"
+.Ft void *
+.Fn fec_free "void *code"
+.Sh "DESCRIPTION"
+This library implements a simple (n,k)
+erasure code based on Vandermonde matrices.
+The encoder takes 
+.Fa k
+packets of size
+.Fa sz
+each, and is able to produce up to
+.Fa n
+different encoded packets, numbered from 0 to n-1,
+such that any subset of
+.Fa k
+of them permits reconstruction of the original data.
+.Pp
+The data structures necessary for the encoding/decoding must
+first be created using calling
+.Fn fec_new
+with the desired parameters. The code descriptor returned by the function
+must be passed to other functions, and destroyed calling
+.Fn fec_free
+.Pp
+Allowed values for k and n depend on a compile-time value
+of
+.Fa GF_BITS
+and must be k <= n <= 2^GF_BITS.
+Best performance is achieved with GF_BITS=8, although the code supports
+also GF_BITS=16.
+.Pp
+Encoding is done by calling
+.Fn fec_encode
+and passing it pointers to the code descriptor, the source and
+destination data packets, the index of the packet to be produced,
+and the size of the packet.
+
+.Pp Decoding is done calling
+.Fn fec_decode
+with pointers to the code, received packets, indexes of received
+packets, and packet size. Decoding is done in place, possibly
+shuffling the arrays passed as parameters.  Decoding is deterministic
+as long as the received packets are different. The decoding procedure
+does some limited testing on this and returns if parameters are
+invalid.
+
+.Sh EXAMPLE
+.nf
+#include <fec.h>
+
+/*
+ * example of sender code
+ */
+void *code ;
+int n, k ;
+
+void *src[] ;
+void *pkt ;
+
+code = new_code (k, n );
+
+for (i = 0 ; i < k ; i++ )
+    src[i] = .. pointer to i-th source packet ..
+for (each packet to transmit) {
+   i = ... index of the packet ;
+   fec_encode(code, src, pkt, i, size) ;
+   .. use packet in pkt
+}
+fec_free(code) ;
+
+/*
+ * example of receiver code
+ */
+void *code ;
+int n, k ;
+
+void *data[] ;
+int *ix[] ;
+
+code = new_code (k, n );
+
+for (i = 0 ; i < k ; i++ ) {
+    ... receive a new packet ...
+    data[i] = .. pointer to i-th source packet ..
+    ix[i] = .. index of i-th source packet ..
+}
+fec_decode(code, data, ix, size) ;
+/*
+ * now data[] has pointers to the source packets
+ */
+   
+.SH BUGS
+Please direct bug reports to luigi at iet.unipi.it .
+.Sh "SEE ALSO"

Added: trunk/contrib/fec/src/csrc/fec.c
===================================================================
--- trunk/contrib/fec/src/csrc/fec.c                            (rev 0)
+++ trunk/contrib/fec/src/csrc/fec.c    2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,899 @@
+/*
+ * fec.c -- forward error correction based on Vandermonde matrices
+ * 980624
+ * (C) 1997-98 Luigi Rizzo (luigi at iet.unipi.it)
+ *
+ * Portions derived from code by Phil Karn (karn at ka9q.ampr.org),
+ * Robert Morelos-Zaragoza (robert at spectra.eng.hawaii.edu) and Hari
+ * Thirumoorthy (harit at 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 <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+
+/*
+ * compatibility stuff
+ */
+#ifdef MSDOS   /* but also for others, e.g. sun... */
+#define NEED_BCOPY
+#define bcmp(a,b,n) memcmp(a,b,n)
+#endif
+
+#ifdef NEED_BCOPY
+#define bcopy(s, d, siz)        memcpy((d), (s), (siz))
+#define bzero(d, siz)   memset((d), '\0', (siz))
+#endif
+
+/*
+ * stuff used for testing purposes only
+ */
+
+#ifdef TEST
+#define DEB(x)
+#define DDB(x) x
+#define        DEBUG   0       /* minimal debugging */
+#ifdef MSDOS
+#include <time.h>
+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 <sys/time.h>
+#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<k; row++)
+       for (col=0; col<k; col++)
+           if ( (row==col && *m != 1) ||
+                (row!=col && *m != 0) )
+                return 0 ;
+           else
+               m++ ;
+    return 1 ;
+}
+#endif /* debug */
+
+/*
+ * invert_mat() takes a matrix and produces its inverse
+ * k is the size of the matrix.
+ * (Gauss-Jordan, adapted from Numerical Recipes in C)
+ * Return non-zero if singular.
+ */
+DEB( int pivloops=0; int pivswaps=0 ; /* diagnostic */)
+static int
+invert_mat(gf *src, int k)
+{
+    gf c, *p ;
+    int irow, icol, row, col, i, ix ;
+
+    int error = 1 ;
+    int *indxc = my_malloc(k*sizeof(int), "indxc");
+    int *indxr = my_malloc(k*sizeof(int), "indxr");
+    int *ipiv = my_malloc(k*sizeof(int), "ipiv");
+    gf *id_row = NEW_GF_MATRIX(1, k);
+    gf *temp_row = NEW_GF_MATRIX(1, k);
+
+    bzero(id_row, k*sizeof(gf));
+    DEB( pivloops=0; pivswaps=0 ; /* diagnostic */ )
+    /*
+     * ipiv marks elements already used as pivots.
+     */
+    for (i = 0; i < k ; i++)
+       ipiv[i] = 0 ;
+
+    for (col = 0; col < k ; col++) {
+       gf *pivot_row ;
+       /*
+        * Zeroing column 'col', look for a non-zero element.
+        * First try on the diagonal, if it fails, look elsewhere.
+        */
+       irow = icol = -1 ;
+       if (ipiv[col] != 1 && src[col*k + col] != 0) {
+           irow = col ;
+           icol = col ;
+           goto found_piv ;
+       }
+       for (row = 0 ; row < k ; row++) {
+           if (ipiv[row] != 1) {
+               for (ix = 0 ; ix < k ; ix++) {
+                   DEB( pivloops++ ; )
+                   if (ipiv[ix] == 0) {
+                       if (src[row*k + ix] != 0) {
+                           irow = row ;
+                           icol = ix ;
+                           goto found_piv ;
+                       }
+                   } else if (ipiv[ix] > 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 (bcmp(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
+     */
+    bzero(retval->enc_matrix, 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)
+         bcopy(src[index], fec, sz*sizeof(gf) ) ;
+    else if (index < code->n) {
+       p = &(code->enc_matrix[index*k] );
+        bzero(fec, 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 );
+}
+
+/*
+ * 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; i<k ; i++) fprintf(stderr, "%3d ", index[i]);
+           fprintf(stderr, "\n");
+           return 1 ;
+       }
+    }
+    )
+    return 0 ;
+}
+
+/*
+ * build_decode_matrix constructs the encoding matrix given the
+ * indexes. The matrix must be already allocated as
+ * a vector of k*k elements, in row-major order
+ */
+static gf *
+build_decode_matrix(struct fec_parms *code, gf *pkt[], int index[])
+{
+    int i , k = code->k ;
+    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) {
+           bzero(p, k*sizeof(gf) );
+           p[i] = 1 ;
+       } else
+#endif
+       if (index[i] < code->n )
+           bcopy( &(code->enc_matrix[index[i]*k]), p, 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" );
+           bzero(new_pkt[row], 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) {
+           bcopy(new_pkt[row], pkt[row], sz*sizeof(gf));
+           free(new_pkt[row]);
+            index[row] = 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 */

Added: trunk/contrib/fec/src/csrc/fec.h
===================================================================
--- trunk/contrib/fec/src/csrc/fec.h                            (rev 0)
+++ trunk/contrib/fec/src/csrc/fec.h    2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,52 @@
+/*
+ * fec.c -- forward error correction based on Vandermonde matrices
+ * 980614
+ * (C) 1997-98 Luigi Rizzo (luigi at iet.unipi.it)
+ *
+ * Portions derived from code by Phil Karn (karn at ka9q.ampr.org),
+ * Robert Morelos-Zaragoza (robert at spectra.eng.hawaii.edu) and Hari
+ * Thirumoorthy (harit at 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
+
+#define        GF_SIZE ((1 << GF_BITS) - 1)    /* powers of \alpha */
+void fec_free(void *p) ;
+void * fec_new(int k, int n) ;
+void init_fec() ;
+void fec_encode(void *code, void *src[], void *dst, int index, int sz) ;
+int fec_decode(void *code, void *pkt[], int index[], int sz) ;
+
+/* end of file */

Added: trunk/contrib/fec/src/csrc/fec16-jinterf.c
===================================================================
--- trunk/contrib/fec/src/csrc/fec16-jinterf.c                          (rev 0)
+++ trunk/contrib/fec/src/csrc/fec16-jinterf.c  2009-03-30 13:57:14 UTC (rev 
26259)
@@ -0,0 +1,201 @@
+#include <fcntl.h>
+#include <string.h>
+#include "jni.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <malloc.h>
+#include "com_onionnetworks_fec_Native16Code.h"
+#include "fec.h"
+
+/*
+ * encode
+ *
+ * @param code This int is actually stores a memory address that points to
+ * an fec_parms struct.
+ */
+JNIEXPORT void JNICALL
+    Java_com_onionnetworks_fec_Native16Code_nativeEncode
+    (JNIEnv *env, jobject obj, jint code, jobjectArray src, jintArray srcOff,
+     jintArray index, jobjectArray ret, jintArray retOff, jint k, 
+     jint packetLength) {
+    
+    jint *localSrcOff, *localIndex, *localRetOff;
+       jbyteArray *inArr, *retArr;
+    jbyte **inarr, **retarr;
+       jobject result = NULL;
+
+       int i, numRet;
+
+       /* allocate memory for the arrays */
+    inArr  = (jbyteArray *) malloc(sizeof(jbyteArray) * k);
+    retArr = (jbyteArray *) malloc(sizeof(jbyteArray) * k);        
+
+    inarr  = (jbyte **) malloc(sizeof(jbyte *) * k);
+    retarr = (jbyte **) malloc(sizeof(jbyte *) * k);
+
+    numRet = (*env)->GetArrayLength(env, ret);
+
+       /* PushLocalFrame reserves enough space for local variable references */
+       if ((*env)->PushLocalFrame(env, k*2+numRet+3) < 0) {
+               return; /* exception OutOfMemoryError */
+       }
+
+    localSrcOff = (*env)->GetIntArrayElements(env, srcOff, NULL);
+    if (localSrcOff == NULL) {
+        return; /* exception occured */
+    }
+
+    localIndex = (*env)->GetIntArrayElements(env, index, NULL);
+    if (localIndex == NULL) {
+        return; /* exception occured */
+    }
+
+    localRetOff = (*env)->GetIntArrayElements(env, retOff, NULL);
+    if (localRetOff == NULL) {
+        return; /* exception occured */
+    }
+
+    for (i=0;i<k;i++) {
+               inArr[i] = ((*env)->GetObjectArrayElement(env, src, i));
+                       if (inArr[i] == NULL) {
+                               return; /* exception occured */
+                       }
+
+               inarr[i] = (*env)->GetPrimitiveArrayCritical(env, inArr[i], 0); 
+        if (inarr[i] == NULL) {
+            return; /* exception occured */
+        }
+        inarr[i] += localSrcOff[i];
+    }
+
+    for (i=0;i<numRet;i++) {
+               retArr[i] = ((*env)->GetObjectArrayElement(env, ret, i));
+        if (retArr[i] == NULL) {
+            return; /* exception occured */
+        }
+
+               retarr[i] = (*env)->GetPrimitiveArrayCritical(env, retArr[i], 
0); 
+        if (retarr[i] == NULL) {
+            return; /* exception occured */
+        }
+        retarr[i] += localRetOff[i];
+    }
+
+    for (i=0;i<numRet;i++) {
+        fec_encode((void *)code, (void **)inarr, (void *)retarr[i], 
+                   (int)localIndex[i], (int)packetLength); 
+    }
+
+    for (i=0;i<k;i++) {
+        inarr[i] -= localSrcOff[i]; 
+               (*env)->ReleasePrimitiveArrayCritical(env, inArr[i], inarr[i], 
0);
+    } 
+ 
+    for (i=0;i<numRet;i++) {
+        retarr[i] -= localRetOff[i];
+               (*env)->ReleasePrimitiveArrayCritical(env, retArr[i], 
retarr[i], 0); 
+    }
+
+    (*env)->ReleaseIntArrayElements(env, srcOff, localSrcOff, 0);
+    (*env)->ReleaseIntArrayElements(env, index, localIndex, 0);
+    (*env)->ReleaseIntArrayElements(env, retOff, localRetOff, 0);
+
+       /* free the memory reserved by PushLocalFrame() */
+       result = (*env)->PopLocalFrame(env, result);
+
+       /* free() may not be necessary. complements malloc() */
+       free(inArr);
+       free(retArr);
+       free(inarr);
+       free(retarr);
+
+}
+
+
+/*
+ * The data[] MUST be preshuffled before this call is made or it WILL NOT
+ * WORK!  It is very difficult to make Java aware that the pointers have
+ * been shuffled in the encode() call, so we must pre-shuffle the data
+ * so that encode doesn't move any pointers around.
+ */
+JNIEXPORT void JNICALL
+    Java_com_onionnetworks_fec_Native16Code_nativeDecode
+    (JNIEnv *env, jobject obj, jint code, jobjectArray data, jintArray dataOff,
+     jintArray whichdata, jint k, jint packetLength) {
+
+    jint *localWhich, *localDataOff;
+       jbyteArray *inArr;
+       jbyte **inarr;
+       jobject result = NULL;
+       
+       int i;
+
+       /* allocate memory for the arrays */
+       inArr = (jbyteArray *) malloc(sizeof(jbyteArray) * k);
+       inarr = (jbyte **) malloc(sizeof(jbyte *) * k);
+
+    localDataOff = (*env)->GetIntArrayElements(env, dataOff, NULL);
+    if (localDataOff == NULL) {
+        return;  /* exception occured */
+    }
+
+    localWhich = (*env)->GetIntArrayElements(env, whichdata, NULL);
+    if (localWhich == NULL) {
+        return;  /* exception occured */
+    }
+
+       /* PushLocalFrame reserves enough space for local variable references */
+       if ((*env)->PushLocalFrame(env, k) < 0) {
+               return; /* exception: OutOfMemoryError */
+       }
+
+    for (i=0;i<k;i++) {
+       inArr[i] = ((*env)->GetObjectArrayElement(env, data, i));
+        if (inArr[i] == NULL) {
+            return;  /* exception occured */
+        }
+       inarr[i] = (*env)->GetPrimitiveArrayCritical(env, inArr[i], 0); 
+        if (inarr[i] == NULL) {
+            return;  /* exception occured */
+        }
+        inarr[i] += localDataOff[i];
+    }
+
+    fec_decode((void *)code, (void **)inarr, (int *)localWhich, 
(int)packetLength);
+
+    for (i = 0; i < k; i++) {
+        inarr[i] -= localDataOff[i];
+        (*env)->SetObjectArrayElement(env, data, i, inArr[i]);
+    }
+
+    for (i = 0; i < k; i++) {
+               (*env)->ReleasePrimitiveArrayCritical(env, inArr[i], inarr[i], 
0); 
+    }
+
+    (*env)->ReleaseIntArrayElements(env, whichdata, localWhich, 0);
+    (*env)->ReleaseIntArrayElements(env, dataOff, localDataOff, 0);
+
+       /* free the memory reserved by PushLocalFrame() */
+       result = (*env)->PopLocalFrame(env, result);
+
+       /* free() may not be necessary. complements malloc() */
+       free(inArr);
+       free(inarr);
+
+}
+
+JNIEXPORT jint JNICALL
+    Java_com_onionnetworks_fec_Native16Code_nativeNewFEC
+    (JNIEnv * env, jobject obj, jint k, jint n) {
+    
+    return (int)fec_new(k,n); 
+
+}
+
+JNIEXPORT void JNICALL
+    Java_com_onionnetworks_fec_Native16Code_nativeFreeFEC
+    (JNIEnv * env, jobject obj, jint code) {
+    
+    fec_free((void *)code); 
+
+}

Added: trunk/contrib/fec/src/csrc/fec8-jinterf.c
===================================================================
--- trunk/contrib/fec/src/csrc/fec8-jinterf.c                           (rev 0)
+++ trunk/contrib/fec/src/csrc/fec8-jinterf.c   2009-03-30 13:57:14 UTC (rev 
26259)
@@ -0,0 +1,201 @@
+#include <fcntl.h>
+#include <string.h>
+#include "jni.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <malloc.h>
+#include "com_onionnetworks_fec_Native8Code.h"
+#include "fec.h"
+
+/*
+ * encode
+ *
+ * @param code This int is actually stores a memory address that points to
+ * an fec_parms struct.
+ */
+JNIEXPORT void JNICALL
+    Java_com_onionnetworks_fec_Native8Code_nativeEncode
+    (JNIEnv *env, jobject obj, jint code, jobjectArray src, jintArray srcOff,
+     jintArray index, jobjectArray ret, jintArray retOff, jint k, 
+     jint packetLength) {
+    
+    jint *localSrcOff, *localIndex, *localRetOff;
+       jbyteArray *inArr, *retArr;
+    jbyte **inarr, **retarr;
+       jobject result = NULL;
+
+       int i, numRet;
+
+       /* allocate memory for the arrays */
+    inArr  = (jbyteArray *) malloc(sizeof(jbyteArray) * k);
+    retArr = (jbyteArray *) malloc(sizeof(jbyteArray) * k);        
+
+    inarr  = (jbyte **) malloc(sizeof(jbyte *) * k);
+    retarr = (jbyte **) malloc(sizeof(jbyte *) * k);
+
+    numRet = (*env)->GetArrayLength(env, ret);
+
+       /* PushLocalFrame reserves enough space for local variable references */
+       if ((*env)->PushLocalFrame(env, k*2+numRet+3) < 0) {
+               return; /* exception OutOfMemoryError */
+       }
+
+    localSrcOff = (*env)->GetIntArrayElements(env, srcOff, NULL);
+    if (localSrcOff == NULL) {
+        return; /* exception occured */
+    }
+
+    localIndex = (*env)->GetIntArrayElements(env, index, NULL);
+    if (localIndex == NULL) {
+        return; /* exception occured */
+    }
+
+    localRetOff = (*env)->GetIntArrayElements(env, retOff, NULL);
+    if (localRetOff == NULL) {
+        return; /* exception occured */
+    }
+
+    for (i=0;i<k;i++) {
+               inArr[i] = ((*env)->GetObjectArrayElement(env, src, i));
+                       if (inArr[i] == NULL) {
+                               return; /* exception occured */
+                       }
+
+               inarr[i] = (*env)->GetPrimitiveArrayCritical(env, inArr[i], 0); 
+        if (inarr[i] == NULL) {
+            return; /* exception occured */
+        }
+        inarr[i] += localSrcOff[i]; 
+    }
+
+    for (i=0;i<numRet;i++) {
+               retArr[i] = ((*env)->GetObjectArrayElement(env, ret, i));
+        if (retArr[i] == NULL) {
+            return; /* exception occured */
+        }
+
+               retarr[i] = (*env)->GetPrimitiveArrayCritical(env, retArr[i], 
0); 
+        if (retarr[i] == NULL) {
+            return; /* exception occured */
+        }
+        retarr[i] += localRetOff[i];
+    }
+
+    for (i=0;i<numRet;i++) {
+        fec_encode((void *)code, (void **)inarr, (void *)retarr[i], 
+                   (int)localIndex[i], (int)packetLength); 
+    }
+
+    for (i=0;i<k;i++) {
+        inarr[i] -= localSrcOff[i]; 
+               (*env)->ReleasePrimitiveArrayCritical(env, inArr[i], inarr[i], 
0);
+    } 
+ 
+    for (i=0;i<numRet;i++) {
+        retarr[i] -= localRetOff[i];
+               (*env)->ReleasePrimitiveArrayCritical(env, retArr[i], 
retarr[i], 0); 
+    }
+
+    (*env)->ReleaseIntArrayElements(env, srcOff, localSrcOff, 0);
+    (*env)->ReleaseIntArrayElements(env, index, localIndex, 0);
+    (*env)->ReleaseIntArrayElements(env, retOff, localRetOff, 0);
+
+       /* free the memory reserved by PushLocalFrame() */
+       result = (*env)->PopLocalFrame(env, result);
+
+       /* free() complements malloc() */
+       free(inArr);
+       free(retArr);
+       free(inarr);
+       free(retarr);
+
+}
+
+
+/*
+ * The data[] MUST be preshuffled before this call is made or it WILL NOT
+ * WORK!  It is very difficult to make Java aware that the pointers have
+ * been shuffled in the encode() call, so we must pre-shuffle the data
+ * so that encode doesn't move any pointers around.
+ */
+JNIEXPORT void JNICALL
+    Java_com_onionnetworks_fec_Native8Code_nativeDecode
+    (JNIEnv *env, jobject obj, jint code, jobjectArray data, jintArray dataOff,
+     jintArray whichdata, jint k, jint packetLength) {
+
+    jint *localWhich, *localDataOff;
+       jbyteArray *inArr;
+       jbyte **inarr;
+       jobject result = NULL;
+       
+       int i;
+
+       /* allocate memory for the arrays */
+       inArr = (jbyteArray *) malloc(sizeof(jbyteArray) * k);
+       inarr = (jbyte **) malloc(sizeof(jbyte *) * k);
+
+    localDataOff = (*env)->GetIntArrayElements(env, dataOff, NULL);
+    if (localDataOff == NULL) {
+        return;  /* exception occured */
+    }
+
+    localWhich = (*env)->GetIntArrayElements(env, whichdata, NULL);
+    if (localWhich == NULL) {
+        return;  /* exception occured */
+    }
+
+       /* PushLocalFrame reserves enough space for local variable references */
+       if ((*env)->PushLocalFrame(env, k) < 0) {
+               return; /* exception: OutOfMemoryError */
+       }
+
+    for (i=0;i<k;i++) {
+       inArr[i] = ((*env)->GetObjectArrayElement(env, data, i));
+        if (inArr[i] == NULL) {
+            return;  /* exception occured */
+        }
+       inarr[i] = (*env)->GetPrimitiveArrayCritical(env, inArr[i], 0); 
+        if (inarr[i] == NULL) {
+            return;  /* exception occured */
+        }
+        inarr[i] += localDataOff[i];
+    }
+
+    fec_decode((void *)code, (void **)inarr, (int *)localWhich, 
(int)packetLength);
+
+    for (i = 0; i < k; i++) {
+        inarr[i] -= localDataOff[i];
+        (*env)->SetObjectArrayElement(env, data, i, inArr[i]);
+    }
+
+    for (i = 0; i < k; i++) {
+               (*env)->ReleasePrimitiveArrayCritical(env, inArr[i], inarr[i], 
0); 
+    }
+
+    (*env)->ReleaseIntArrayElements(env, whichdata, localWhich, 0);
+    (*env)->ReleaseIntArrayElements(env, dataOff, localDataOff, 0);
+
+       /* free the memory reserved by PushLocalFrame() */
+       result = (*env)->PopLocalFrame(env, result);
+
+       /* free() may not be necessary. complements malloc() */
+       free(inArr);
+       free(inarr);
+
+}
+
+JNIEXPORT jint JNICALL
+    Java_com_onionnetworks_fec_Native8Code_nativeNewFEC
+    (JNIEnv * env, jobject obj, jint k, jint n) {
+    
+    return (int)fec_new(k,n); 
+
+}
+
+JNIEXPORT void JNICALL
+    Java_com_onionnetworks_fec_Native8Code_nativeFreeFEC
+    (JNIEnv * env, jobject obj, jint code) {
+    
+    fec_free((void *)code); 
+
+}

Added: trunk/contrib/fec/src/csrc/fec_win32.c
===================================================================
--- trunk/contrib/fec/src/csrc/fec_win32.c                              (rev 0)
+++ trunk/contrib/fec/src/csrc/fec_win32.c      2009-03-30 13:57:14 UTC (rev 
26259)
@@ -0,0 +1,882 @@
+/*
+ * fec.c -- forward error correction based on Vandermonde matrices
+ * 980624
+ * (C) 1997-98 Luigi Rizzo (luigi at iet.unipi.it)
+ *
+ * Portions derived from code by Phil Karn (karn at ka9q.ampr.org),
+ * Robert Morelos-Zaragoza (robert at spectra.eng.hawaii.edu) and Hari
+ * Thirumoorthy (harit at 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 <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+
+
+/*
+ * compatibility stuff
+ */
+
+#define NEED_BCOPY 0
+#define bcmp(a,b,n) memcmp(a,b,n)
+
+
+#ifdef NEED_BCOPY
+#define bcopy(s, d, siz)        memcpy((d), (s), (siz))
+#define bzero(d, siz)                  memset((d), '\0', (siz))
+#endif
+
+/*
+ * stuff used for testing purposes only
+ */
+
+
+#define DEB(x)
+#define DDB(x)
+#define TICK(x)
+#define TOCK(x)
+
+typedef unsigned long u_long;
+
+/*
+ * 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 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 8 /* 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) ;
+       //printf("ADDMUL");
+#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 );
+
+       //printf("ADDMUL-out");
+}
+
+/*
+ * 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<k; row++)
+       for (col=0; col<k; col++)
+           if ( (row==col && *m != 1) ||
+                (row!=col && *m != 0) )
+                return 0 ;
+           else
+               m++ ;
+    return 1 ;
+}
+#endif /* debug */
+
+/*
+ * invert_mat() takes a matrix and produces its inverse
+ * k is the size of the matrix.
+ * (Gauss-Jordan, adapted from Numerical Recipes in C)
+ * Return non-zero if singular.
+ */
+DEB( int pivloops=0; int pivswaps=0 ; /* diagnostic */)
+static int
+invert_mat(gf *src, int k)
+{
+    gf c, *p ;
+    int irow, icol, row, col, i, ix ;
+
+    int error = 1 ;
+    int *indxc = my_malloc(k*sizeof(int), "indxc");
+    int *indxr = my_malloc(k*sizeof(int), "indxr");
+    int *ipiv = my_malloc(k*sizeof(int), "ipiv");
+    gf *id_row = NEW_GF_MATRIX(1, k);
+    gf *temp_row = NEW_GF_MATRIX(1, k);
+
+    bzero(id_row, k*sizeof(gf));
+    DEB( pivloops=0; pivswaps=0 ; /* diagnostic */ )
+    /*
+     * ipiv marks elements already used as pivots.
+     */
+    for (i = 0; i < k ; i++)
+       ipiv[i] = 0 ;
+
+    for (col = 0; col < k ; col++) {
+       gf *pivot_row ;
+       /*
+        * Zeroing column 'col', look for a non-zero element.
+        * First try on the diagonal, if it fails, look elsewhere.
+        */
+       irow = icol = -1 ;
+       if (ipiv[col] != 1 && src[col*k + col] != 0) {
+           irow = col ;
+           icol = col ;
+           goto found_piv ;
+       }
+       for (row = 0 ; row < k ; row++) {
+           if (ipiv[row] != 1) {
+               for (ix = 0 ; ix < k ; ix++) {
+                   DEB( pivloops++ ; )
+                   if (ipiv[ix] == 0) {
+                       if (src[row*k + ix] != 0) {
+                           irow = row ;
+                           icol = ix ;
+                           goto found_piv ;
+                       }
+                   } else if (ipiv[ix] > 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 (bcmp(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 ) {
+       //printf(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
+     */
+    bzero(retval->enc_matrix, 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) {
+
+         bcopy(src[index], fec, sz*sizeof(gf) ) ;
+
+       } else if (index < code->n) {
+
+               p = &(code->enc_matrix[index*k] );
+
+               //printf("sz:"+sz+" gf:"+gf);
+        bzero(fec, sz*sizeof(gf));
+
+               for (i = 0; i < k ; i++) {
+
+                       //printf("src[i]:"+src[i]);
+            addmul(fec, src[i], p[i], sz ) ;
+
+               }
+    } //else
+       //printf(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; i<k ; i++) fprintf(stderr, "%3d ", index[i]);
+           //fprintf(stderr, "\n");
+           return 1 ;
+       }
+    }
+    )
+    return 0 ;
+}
+
+/*
+ * build_decode_matrix constructs the encoding matrix given the
+ * indexes. The matrix must be already allocated as
+ * a vector of k*k elements, in row-major order
+ */
+static gf *
+build_decode_matrix(struct fec_parms *code, gf *pkt[], int index[])
+{
+    int i , k = code->k ;
+    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) {
+           bzero(p, k*sizeof(gf) );
+           p[i] = 1 ;
+       } else
+#endif
+       if (index[i] < code->n )
+           bcopy( &(code->enc_matrix[index[i]*k]), p, 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" );
+           bzero(new_pkt[row], 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) {
+           bcopy(new_pkt[row], pkt[row], sz*sizeof(gf));
+           free(new_pkt[row]);
+            index[row] = 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 */

Added: trunk/contrib/fec/src/csrc/test.c
===================================================================
--- trunk/contrib/fec/src/csrc/test.c                           (rev 0)
+++ trunk/contrib/fec/src/csrc/test.c   2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,273 @@
+/*
+ * test.c -- test code for FEC library
+ *
+ * (C) 1997-98 Luigi Rizzo (luigi at iet.unipi.it)
+ */
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include "fec.h"
+
+/*
+ * compatibility stuff
+ */
+#ifdef MSDOS   /* but also for others, e.g. sun... */
+#define NEED_BCOPY
+#define bcmp(a,b,n) memcmp(a,b,n)
+#endif
+
+#ifdef NEED_BCOPY
+#define bcopy(s, d, siz)        memcpy((d), (s), (siz))
+#define bzero(d, siz)   memset((d), '\0', (siz))
+#endif
+
+/*
+ * stuff used for testing purposes only
+ */
+
+#define DEB(x)
+#define DDB(x) x
+#define        DEBUG   0       /* minimal debugging */
+#ifdef MSDOS
+#include <time.h>
+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 <sys/time.h>
+#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 */
+
+void *
+my_malloc(int sz, char *s)
+{
+    void *p = malloc(sz) ;
+    if (p != NULL)
+       return p ;
+    fprintf(stderr, "test: malloc failure for %d bytes in <%s>\n",
+       sz, s);
+    exit(1);
+}
+int
+pr_matrix(void *m1, int rows, int cols, char *s)
+{
+    int r, c;
+#if GF_BITS >8
+    u_char *m = m1 ;
+#else
+    u_short *m = m1 ;
+#endif
+    fprintf(stderr,"%s\n", s);
+    for (r=0; r<rows; r++) {
+       for (c=0; c<cols; c++)
+           fprintf(stderr,"%3d  ",m[r*cols+c]);
+       fprintf(stderr,"\n");
+    }
+    fprintf(stderr,"\n");
+    return 0;
+}
+
+/*
+ * the following is only test code and can be safely omitted
+ * in applications.
+ * Creates k packets of size sz of random data, encodes them,
+ * and tries to decode.
+ * Index contains the permutation entry.
+ */
+
+int
+test_decode(void *code, int k, int index[], int sz, char *s)
+{
+    int errors;
+    int reconstruct = 0 ;
+    int item, i ;
+
+    static int prev_k = 0, prev_sz = 0;
+    static u_char **d_original = NULL, **d_src = NULL ;
+
+    if (sz < 1 || sz > 8192) {
+       fprintf(stderr, "test_decode: size %d invalid, must be 1..8K\n",
+               sz);
+       return 1 ;
+    }
+    if (k < 1 || k > GF_SIZE + 1) {
+       fprintf(stderr, "test_decode: k %d invalid, must be 1..%d\n",
+               k, GF_SIZE + 1 );
+       return 2 ;
+    }
+    if (prev_k != k || prev_sz != sz) {
+       if (d_original != NULL) {
+           for (i = 0 ; i < prev_k ; i++ ) {
+               free(d_original[i]);
+               free(d_src[i]);
+           }
+           free(d_original);
+           free(d_src);
+           d_original = NULL ;
+           d_src = NULL ;
+       }
+    }
+    prev_k = k ;
+    prev_sz = sz ;
+    if (d_original == NULL) {
+       d_original = my_malloc(k * sizeof(void *), "d_original ptr");
+       d_src = my_malloc(k * sizeof(void *), "d_src ptr");
+
+       for (i = 0 ; i < k ; i++ ) {
+           d_original[i] = my_malloc(sz, "d_original data");
+           d_src[i] = my_malloc(sz, "d_src data");
+       }
+       /*
+        * build sample data
+        */
+       for (i = 0 ; i < k ; i++ ) {
+           for (item=0; item < sz; item++)
+               d_original[i][item] = ((item ^ i) + 3) & GF_SIZE;
+       }
+    }
+
+    errors = 0 ;
+
+    for( i = 0 ; i < k ; i++ )
+       if (index[i] >= k ) reconstruct ++ ;
+
+    TICK(ticks[2]);
+    for( i = 0 ; i < k ; i++ )
+       fec_encode(code, d_original, d_src[i], index[i], sz );
+    TOCK(ticks[2]);
+
+    TICK(ticks[1]);
+    if (fec_decode(code, d_src, index, sz)) {
+       fprintf(stderr, "detected singular matrix for %s  \n", s);
+       return 1 ;
+    }
+    TOCK(ticks[1]);
+
+    for (i=0; i<k; i++)
+       if (bcmp(d_original[i], d_src[i], sz )) {
+           errors++;
+           fprintf(stderr, "error reconstructing block %d\n", i);
+       }
+    if (errors)
+       fprintf(stderr, "Errors reconstructing %d blocks out of %d\n",
+           errors, k);
+
+    fprintf(stderr,
+       "  k %3d, l %3d  c_enc %10.6f MB/s c_dec %10.6f MB/s     \r",
+       k, reconstruct,
+       (double)(k * sz * reconstruct)/(double)ticks[2],
+       (double)(k * sz * reconstruct)/(double)ticks[1]);
+    return errors ;
+}
+
+#if 0
+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
+
+#define KK 64 /* 255 */
+#define SZ 1024
+int
+main(int argc, char *argv[])
+{
+    char buf[256];
+    void *code ;
+
+    int kk ;
+    int i ;
+
+    int *ixs ;
+
+    int lim = GF_SIZE + 1 ;
+
+    if (lim > 1024) lim = 1024 ;
+
+#if 0
+    test_gf();
+#endif
+    for ( kk = KK ; kk > 2 ; kk-- ) {
+       code = fec_new(kk, lim);
+       ixs = my_malloc(kk * sizeof(int), "ixs" );
+
+       for (i=0; i<kk; i++) ixs[i] = kk - i ;
+       sprintf(buf, "kk=%d, kk - i", kk); 
+       test_decode(code, kk, ixs, SZ, buf);
+
+       for (i=0; i<kk; i++) ixs[i] = i ;
+       test_decode(code, kk, ixs, SZ, "i");
+
+if (0) {
+       for (i=0; i<kk; i++) ixs[i] = i ;
+       ixs[0] = ixs[kk/2] ;
+       test_decode(code, kk, ixs, SZ, "0 = 1 (error expected)");
+       }
+
+if (0)
+       for (i= lim-1 ; i >= kk ; i--) {
+           int j ;
+           for (j=0; j<KK; j++) ixs[j] = kk - j ;
+           ixs[0] = i ;
+           test_decode(code, kk, ixs, SZ, "0 = big");
+       }
+
+if (0)
+       for (i= lim - kk ; i >= 0 && i>= lim - kk - 4 ; i--) {
+           int j ;
+           for (j=0; j<kk; j++)
+               ixs[j] = kk -1 - j + i ;
+           test_decode(code, kk, ixs, SZ, "shifted j");
+       }
+if (1)  {
+       int j, max_i0 = KK/2 ;
+       if (max_i0 + KK > lim)
+           max_i0 = lim - KK ;
+       for (i= 0 ; i <= max_i0 ; i++) {
+           for (j=0; j<kk; j++)
+               ixs[j] = j + i ;
+           test_decode(code, kk, ixs, SZ, "shifted j");
+       }
+       }
+       fprintf(stderr, "\n");
+       free(ixs);
+       fec_free(code);
+    }
+    return 0;
+}

Added: trunk/contrib/fec/src/csrc/w32
===================================================================
--- trunk/contrib/fec/src/csrc/w32                              (rev 0)
+++ trunk/contrib/fec/src/csrc/w32      2009-03-30 13:57:14 UTC (rev 26259)
@@ -0,0 +1,65 @@
+MAKE := make -f w32
+
+CPP := cl.exe
+
+MODE ?= Debug
+
+CPP_OPTS := /nologo /I c:/jdk1.3/include /I c:/jdk1.3/include/win32 \
+       /D WIN32 /D _WINDOWS /D _MBCS /D _USRDLL /D FEC_EXPORTS 
+
+ifeq ($(MODE),Debug)
+Mode := Debug
+CPP_OPTS := /nologo /MTd /W3 /GZ /ZI /Od /D DEBUG $(CPP_OPTS)
+else
+MODE := Release
+CPP_OPTS := /nologo /MT /W3 /O1 /D NDEBUG $(CPP_OPTS)
+endif
+
+LIBS := kernel32 user32 #advapi32 shell32 gdi32 winspool cmdlg32 ole32 
oleaut32 uuid odbc32 odbccp32 
+
+LDFLAGS=$(patsubt %,%.lib,$(LIBS)) /nologo /dll /incremental:no /machine:I386 \
+       /out:$(MODE)/fec$(BITS).dll /implib:$(MODE)/fec$(BITS).lib \
+       /OPT:REF /MAP
+
+LD=link.exe
+
+LDOBJS= $(MODE)/fec_win32.obj $(MODE)/fec$(BITS)-jinterf.obj
+
+
+.PHONY: all feclib debug-all release-all bits8-all bits16-all clean
+
+all: debug-all release-all
+
+feclib: $(MODE)/fec$(BITS).dll
+
+debug-all: Debug 
+       $(MAKE) BITS=8 MODE=Debug feclib
+       $(MAKE) BITS=16 MODE=Debug feclib
+
+release-all: Release
+       $(MAKE) BITS=8 MODE=Release feclib
+       $(MAKE) BITS=16 MODE=Release feclib
+
+bits8-all: Debug Release
+       $(MAKE) BITS=8 MODE=Debug feclib
+       $(MAKE) BITS=8 MODE=Release feclib
+
+bits16-all: Debug Release
+       $(MAKE) BITS=16 MODE=Debug feclib
+       $(MAKE) BITS=16 MODE=Release feclib
+
+clean: 
+       rm -f Debug Release
+
+Debug:
+       mkdir -p Debug
+
+Release:
+       mkdir -p Release
+
+$(MODE)/fec$(BITS).dll : $(MODE) $(DEF_FILE) $(LDOBJS)
+       $(LD) $(LDFLAGS) $(LDOBJS)
+
+
+$(MODE)/%.obj: %.c
+       $(CPP) $(CPP_OPTS) /c /Fo"$@" $<


Reply via email to