[MN-dev] [mndiff]: r82 - in trunk/noe: Makefile galois.c galois.h rs.c rs.h test.c

michael subversion at mplayerhq.hu
Tue Jul 10 00:53:52 CEST 2007


Author: michael
Date: Tue Jul 10 00:53:52 2007
New Revision: 82

Log:
version from unknown date but definitly before 2006-06-22 13:51


Modified:
   trunk/noe/Makefile
   trunk/noe/galois.c
   trunk/noe/galois.h
   trunk/noe/rs.c
   trunk/noe/rs.h
   trunk/noe/test.c

Modified: trunk/noe/Makefile
==============================================================================
--- trunk/noe/Makefile	(original)
+++ trunk/noe/Makefile	Tue Jul 10 00:53:52 2007
@@ -17,12 +17,15 @@ CFLAGS  = -g -Wall -O4 $(OPTFLAGS) -I. $
 %_10001.o: %.c
 	$(CC) $(CFLAGS) -c -DSIZE=0x10001 -o $@ $<
 
-all: test_100 test_101 test_10001 libnoe_101.a libnoe_10001.a
+all: test_100 test_101 test_10001
 
 libnoe_100.a: $(LIBOBJS:.o=_100.o)
 	$(AR) rc $@ $^
 
-libnoe_1%1.a: $(LIBOBJS:.o=_1%1.o) gfft_1%1.o
+libnoe_101.a: $(LIBOBJS:.o=_101.o) gfft_101.o
+	$(AR) rc $@ $^
+
+libnoe_10001.a: $(LIBOBJS:.o=_10001.o) gfft_10001.o
 	$(AR) rc $@ $^
 
 test: all

Modified: trunk/noe/galois.c
==============================================================================
--- trunk/noe/galois.c	(original)
+++ trunk/noe/galois.c	Tue Jul 10 00:53:52 2007
@@ -27,24 +27,25 @@
 #include "galois.h"
 #include "gfft.h"
 
-GFF4Element EXT(exp)[SIZE];
+GFF4Element EXT(exp)[4*SIZE];
 GFF4Element EXT(log)[SIZE];
 
 void EXT(init)(){
     GFF4Element ge= 1;
     int i;
     
-    for(i=0; i<SIZE; i++){
+    for(i=0; i<=(SIZE-2)*2; i++){
         EXT(exp)[i]= ge;
         EXT(log)[ge]= i%(SIZE-1);
 #if M!=0
-        assert(PRIMITIVE_ELEMENT==2);
         ge+=ge;
         if(ge&SIZE) ge ^= M;
 #else
-        ge= prod(ge, PRIMITIVE_ELEMENT);
+        ge= reduce(ge * PRIMITIVE_ELEMENT);
 #endif
     }
+    EXT(log)[0]= (SIZE-2)*2 + 1;
+
 #if M==0
     noe_gfft_init();
 #endif
@@ -52,13 +53,7 @@ void EXT(init)(){
 
 int noe_log2(unsigned int v) // from ffmpeg
 {
-    int n;
-
-    n = 0;
-    if (v & 0xffff0000) {
-        v >>= 16;
-        n += 16;
-    }
+    int n= 0;
     if (v & 0xff00) {
         v >>= 8;
         n += 8;
@@ -98,10 +93,8 @@ void EXT(prodPoly)(GFF4Element *dst, GFF
             j  = NOE_MAX(0, i - order2);
             end= NOE_MIN(i, order1);
 
-            for(; j<=end; j++){
-                if(M) acc^= prod(src1[j], src2[i-j]);
-                else  acc+= prod(src1[j], src2[i-j]);
-            }
+            for(; j<=end; j++)
+                acc ADD_OP prod(src1[j], src2[i-j]);
             dst[i]= reduce(acc);
         }
     }else{
@@ -128,6 +121,7 @@ void EXT(prodPoly)(GFF4Element *dst, GFF
     }
 }
 
+#ifdef ALTERNATIVE_SOLVERS
 void EXT(prodMatrix)(GFF4Element *d1[2], GFF4Element *d2[2], GFF4Element *s1[2], GFF4Element *s2[2]){
 //FIXME improve threshold decission
     if((s1[0][0] + s1[1][0] + s2[0][0] + s2[1][0]) < 64 || M/*|| 
@@ -195,6 +189,7 @@ void EXT(prodMatrix)(GFF4Element *d1[2],
         }
     }
 }
+#endif
 
 /**
  * returns the first order2 coefficients of the product of the 2 polynoms.
@@ -218,10 +213,8 @@ void EXT(partialProdPoly)(GFF4Element *d
             j  = NOE_MAX(0, i - order2);
             end= NOE_MIN(i, order1);
 
-            for(; j<=end; j++){
-                if(M) acc^= prod(src1[j], src2[i-j]);
-                else  acc+= prod(src1[j], src2[i-j]);
-            }
+            for(; j<=end; j++)
+                acc ADD_OP prod(src1[j], src2[i-j]);
             dst[i]= reduce(acc);
         }
     }else{
@@ -262,13 +255,11 @@ void EXT(printPoly)(GFF4Element *src){
  * Evaluates the src polynom at x.
  */
 GFF4Element EXT(evalPoly)(GFF4Element *src, GFF4Element x){
-    unsigned int acc=0;
+    unsigned int acc=src[ src[0] + 1 ];
     int j;
 
-    for(j=src[0]+1; j>0; j--){
-//FIXME skip reduce
+    for(j=src[0]; j>0; j--)
         acc = sum(prod(acc, x), src[j]);
-    }
 
     return acc;
 }
@@ -283,8 +274,7 @@ void EXT(getDerivative)(GFF4Element *dst
         dst[i  ]= src[i+1];
     }
     dst[0]= i-3;
-    assert(dst[0] >= 0);
-    assert(dst[ dst[0]+1 ]); //FIXME this isnt guranteed but seems to be always true
+    EXT(getOrder)(dst);
 #else
     for(i=1; i<=src[0]; i++){
         dst[i  ]= prod(src[i+1], i);
@@ -296,62 +286,37 @@ void EXT(getDerivative)(GFF4Element *dst
 /**
  * ...
  */
-void EXT(synthPoly)(GFF4Element *dst, GFF4Element *src, int count){
-    if(count<20){
-        int i;
-
-        SET_POLY0(dst, 1);
-
-        for(i=0; i<count; i++){
-            assert(src[3*i] == 1);
-            EXT(prodPoly)(dst, dst, src+3*i);
-        }
-    }else{
-        int pass, i, countLeft;
-        const int passCount= noe_log2(count-1)+1;
-        GFF4Element temp[2][4<<noe_log2(count-1)];
-
-        assert(passCount>0);
-
-        countLeft= count;
-
-    //  EXT(printPoly)(src, 1);
-    //  EXT(printPoly)(src+2, 1);
+void EXT(synthPoly)(GFF4Element *dst, int *src, int count){
+        int pass, i;
+        GFF4Element temp[2][4*count];
+        GFF4Element *temp1= count<=2 ? dst : temp[0];
 
-        for(i=0; i<(count>>1); i++){
-            EXT(prodPoly)(temp[0]+4*i, src+6*i, src+6*i+3);
+        for(i=0; i+1<count; i+=2){
+            GFF4Element a= neg(EXT(exp)[src[i+0]]);
+            GFF4Element b= neg(EXT(exp)[src[i+1]]);
+            SET_POLY2(temp1+2*i, 1, sum(a,b), prod(a,b));
         }
-        if(count&1){
-            memcpy(temp[0]+4*i, src+6*i, 3*sizeof(GFF4Element));
-            countLeft++;
+        if(i<count){
+            SET_POLY1(temp1+2*i, 1, neg(EXT(exp)[src[i+0]]));
         }
-        countLeft>>=1;
-
-    //  EXT(printPoly)(temp[0], 2);
+        count= (count+1)>>1;
 
-        for(pass=1; countLeft>1; pass++){
+        for(pass=1; count>1; pass++){
             const int step= 2<<pass;
-            GFF4Element *temp1= temp[ pass&1   ];
+            GFF4Element *temp1= count<=2 ? dst : temp[ pass&1   ];
             GFF4Element *temp2= temp[(pass&1)^1];
 
-            for(i=0; i<(countLeft>>1); i++){
-                EXT(prodPoly)(temp1+2*step*i, temp2+2*step*i, temp2+2*step*i+step);
+            for(i=0; i+1<count; i+=2){
+                EXT(prodPoly)(temp1, temp2, temp2+step);
+                temp1 += 2*step;
+                temp2 += 2*step;
             }
 
-            if(countLeft&1){
-                memcpy(temp1+2*step*i, temp2+2*step*i, step*sizeof(GFF4Element));
-                countLeft++;
+            if(i<count){
+                memcpy(temp1, temp2, step*sizeof(GFF4Element));
             }
-            countLeft>>=1;
+            count= (count+1)>>1;
         }
-
-        assert(pass==passCount);
-
-    //  EXT(printPoly)(temp[0], count);
-
-        memcpy(dst, temp[(pass&1)^1], (count+2)*sizeof(GFF4Element));
-    //  EXT(printPoly)(dst, count);
-    }
 }
 
 void EXT(getOrder)(GFF4Element *src){
@@ -359,6 +324,7 @@ void EXT(getOrder)(GFF4Element *src){
         src[0]--;
 }
 
+#ifdef ALTERNATIVE_SOLVERS
 /**
  * ...
  * Note: rem can be identical to nom
@@ -475,3 +441,4 @@ void EXT(scaledSumPoly)(GFF4Element *dst
         EXT(getOrder)(dst-1); //FIXME we need this due to src2= 0
     }
 }
+#endif

Modified: trunk/noe/galois.h
==============================================================================
--- trunk/noe/galois.h	(original)
+++ trunk/noe/galois.h	Tue Jul 10 00:53:52 2007
@@ -21,16 +21,19 @@
 #define SHIFT 8
 #define M 0x11D
 #define EXT(name) noe_ ## name ## _100
+#define ADD_OP ^=
 #elif SIZE == 0x10001
 #define PRIMITIVE_ELEMENT 3
 #define SHIFT 16
 #define M 0
 #define EXT(name) name ## _10001
+#define ADD_OP +=
 #elif SIZE == 0x101
 #define PRIMITIVE_ELEMENT 3
 #define SHIFT 8
 #define M 0
 #define EXT(name) noe_ ## name ## _101
+#define ADD_OP +=
 #else
 #error wrong SIZE
 #endif
@@ -40,7 +43,7 @@
 
 typedef unsigned int GFF4Element;
 
-extern GFF4Element EXT(exp)[SIZE];
+extern GFF4Element EXT(exp)[4*SIZE];
 extern GFF4Element EXT(log)[SIZE];
 
 void EXT(prodPoly)(GFF4Element *dst, GFF4Element *src1, GFF4Element *src2);
@@ -49,13 +52,17 @@ void EXT(init)();
 void EXT(printPoly)(GFF4Element *src);
 void EXT(getDerivative)(GFF4Element *dst, GFF4Element *src);
 GFF4Element EXT(evalPoly)(GFF4Element *src, GFF4Element x);
-void EXT(synthPoly)(GFF4Element *dst, GFF4Element *src, int count);
+void EXT(synthPoly)(GFF4Element *dst, int *src, int count);
+#ifdef ALTERNATIVE_SOLVERS
 void EXT(divPoly)(GFF4Element *quot, GFF4Element *rem, GFF4Element *nom, GFF4Element *denom);
 void EXT(diffPoly)(GFF4Element *dst, GFF4Element *src1, GFF4Element *src2);
 void EXT(sumPoly)(GFF4Element *dst, GFF4Element *src1, GFF4Element *src2);
+#endif
 void EXT(getOrder)(GFF4Element *src);
+#ifdef ALTERNATIVE_SOLVERS
 void EXT(scaledSumPoly)(GFF4Element *dst, GFF4Element *src1, GFF4Element *src2, GFF4Element f, int s);
 void EXT(prodMatrix)(GFF4Element *d1[2], GFF4Element *d2[2], GFF4Element *s1[2], GFF4Element *s2[2]);
+#endif
 int noe_log2(unsigned int v);
 
 static inline GFF4Element reduce(GFF4Element a){
@@ -67,30 +74,12 @@ static inline GFF4Element reduce(GFF4Ele
 }
 
 static inline GFF4Element prod(GFF4Element a, GFF4Element b){
-#if M!=0
-    if(a==0 || b==0) return 0;
-    a= EXT(log)[a] + EXT(log)[b];
-    if(a>=SIZE) a-= SIZE-1;
-    return EXT(exp)[a];
-#else
-//printf("%d %d\n", a,b);
-#if 0
-    uint64_t x= ((uint64_t)a)*b;
-    return reduce(x + (x>>32));
+#if M!=0 || SIZE < 0x500
+    return EXT(exp)[EXT(log)[a] + EXT(log)[b]];
 #else
     if(a==MINUS1)  return -b + ((((int)-b)>>31)&SIZE);
     else           return reduce(a*b);
 #endif
-#endif
-}
-
-static inline GFF4Element sum(GFF4Element a, GFF4Element b){
-#if M!=0
-    return a^b;
-#else
-    a += b - SIZE;
-    return a + ((((int)a)>>31)&SIZE);
-#endif
 }
 
 static inline GFF4Element diff(GFF4Element a, GFF4Element b){
@@ -108,11 +97,14 @@ static inline GFF4Element diff(GFF4Eleme
 #define neg(a) (SIZE - (a))
 #endif
 
+#define sum(a,b) diff(a,neg(b))
+
 static inline GFF4Element inv(GFF4Element a){
     assert(a!=0);
 
-    return EXT(exp)[MINUS1 - EXT(log)[a]];
+    return (EXT(exp)+MINUS1)[- EXT(log)[a]];
 }
 
 #define SET_POLY0(dst, coeff0) (dst)[0]=0; (dst)[1]=coeff0;
 #define SET_POLY1(dst, coeff0, coeff1) (dst)[0]=1; (dst)[1]=coeff0; (dst)[2]=coeff1;
+#define SET_POLY2(dst, coeff0, coeff1, coeff2) (dst)[0]=2; (dst)[1]=coeff0; (dst)[2]=coeff1; (dst)[3]=coeff2;

Modified: trunk/noe/rs.c
==============================================================================
--- trunk/noe/rs.c	(original)
+++ trunk/noe/rs.c	Tue Jul 10 00:53:52 2007
@@ -31,90 +31,49 @@
 #undef NDEBUG
 #include <assert.h>
 
-/**
- * 
- */
-void EXT(getRsGenerator)(GFF4Element *dst, int first, int order){
-    int i;
-    GFF4Element factor[3*order];
-
-    for(i=0; i<order; i++){
-        SET_POLY1(factor+3*i, neg(EXT(exp)[ first + i ]), 1);
-    }
-    
-    EXT(synthPoly)(dst, factor, order);
-}
-
 void EXT(getSyndrom)(GFF4Element *syn, GFF4Element *src, int order){
     *syn++= order;
-#if M==0
-    EXT(gfft)(syn, src, SHIFT);
-    EXT(permute)(syn, syn, order+1);
-#else
-{
-    int i;
-    int bak= src[-1];
-    src[-1]= SIZE-2;
-    for(i=0; i<order; i++){
-        syn[i+1]= EXT(evalPoly)(src-1, EXT(exp)[1 + i]);
+    if(order>6 && !M){
+        EXT(gfft)(syn, src, SHIFT);
+        EXT(permute)(syn, syn, order+1);
+    }else{
+        int i;
+        for(i=1; i<=order; i++){
+            int bak= src[-1];
+            src[-1]= SIZE-2;
+            syn[i]= EXT(evalPoly)(src-1, EXT(exp)[i]);
+            src[-1]= bak;
+        }
     }
-    src[-1]= bak;
-}
-#endif
     syn[0]= 1;
 }
-#if 1
-noe_RsEncoder *EXT(getRsEncoder)(int parityCount){
+
+noe_RsEncoder *EXT(getRsEncoder)(int parityCount, int tLocation){
     const int dataCount= SIZE - parityCount - 1;
-    int i;  
-    GFF4Element *locator, *factor;
-    GFF4Element locatorDerivative[parityCount+1];
     noe_RsEncoder *encoder= memalign(16, sizeof(noe_RsEncoder));
-    GFF4Element locationFactor[3*parityCount];
     
-    locator= encoder->parityLocator= memalign(16, sizeof(GFF4Element)*(parityCount+2));
-    factor=  encoder->parityFactor = memalign(16, sizeof(GFF4Element)* parityCount   );
+    encoder->tPoly=NULL;
 
     encoder->parityCount= parityCount;
-    
-    for(i=0; i<parityCount; i++){
-        SET_POLY1(locationFactor+3*i, 1, neg(EXT(exp)[dataCount + i]));
-    }
-    EXT(synthPoly)(locator, locationFactor, parityCount);
-#if 0
-for(i=0; i<parityCount; i++){
-    if(EXT(evalPoly)(locator, inv(EXT(exp)[dataCount + i]))){
-        printf("Internal Error 1\n");
-    }
-}
-if(!EXT(evalPoly)(locator, inv(EXT(exp)[dataCount - 1])))
-    printf("Internal Error 2\n");
-#endif
-    EXT(getDerivative)(locatorDerivative, locator);
-    assert(locatorDerivative[0] == parityCount-1);
 
-    for(i=0; i<parityCount; i++){
-        GFF4Element X= EXT(exp)[dataCount + i];
-        GFF4Element invX= EXT(exp)[MINUS1 - (dataCount + i)];
-        GFF4Element denom;
-        
-        assert(X == inv(invX));
-
-        denom= EXT(evalPoly)(locatorDerivative, invX); //FIXME do in freq domain if parityCount>1000
-        factor[i]= prod(X, inv(denom));
+    if(tLocation>=0 && tLocation<dataCount && !M){
+        GFF4Element temp[SIZE - 1];
+        memset(temp, 0, sizeof(temp));
+        temp[tLocation]= 1;
+        EXT(rsEncode)(encoder, temp);
+        encoder->tPoly= memalign(16, sizeof(GFF4Element)*(parityCount));
+        encoder->tLocation= tLocation;
+        memcpy(encoder->tPoly, temp + dataCount, parityCount*sizeof(GFF4Element));
     }
     
     return encoder;
 }
 
+#define free(p) if(p) free(p); p=NULL;
 void EXT(freeRsEncoder)(noe_RsEncoder *encoder){
     if(!encoder) return;
     
-    if(encoder->parityLocator) free(encoder->parityLocator);
-    encoder->parityLocator= NULL;
-    
-    if(encoder->parityFactor) free(encoder->parityFactor);
-    encoder->parityFactor= NULL;
+    free(encoder->tPoly);
 
     free(encoder);
 }
@@ -126,39 +85,59 @@ void EXT(rsEncode)(noe_RsEncoder *encode
     int i;
     const int parityCount= encoder->parityCount;
     const int dataCount= SIZE - parityCount - 1;
-    GFF4Element syn[SIZE];
-    GFF4Element omega[SIZE];
-    GFF4Element *locator= encoder->parityLocator;
-    GFF4Element *factor= encoder->parityFactor;
+    int erasureLocations[parityCount];
 
-    memset(data + dataCount, 0, parityCount*sizeof(GFF4Element));
-    
-    EXT(getSyndrom)(syn, data, parityCount);
+    for(i=0; i<parityCount; i++)
+        erasureLocations[i]= dataCount + i;
 
-    EXT(partialProdPoly)(omega, locator, syn, syn[0]);
+    EXT(rsDecode)(data, erasureLocations, parityCount, parityCount);
+}
 
-    if(parityCount < 1000 || M){
-        for(i=0; i<parityCount; i++){
-            GFF4Element invX= EXT(exp)[MINUS1 - (dataCount + i)];
-            
-            GFF4Element parity= prod(EXT(evalPoly)(omega, invX), factor[i]);
-            
-            data[dataCount + i]= parity;
+#if M==0
+//FIXME avoid the noe_RsEncoder dependancy
+int EXT(rsTransform)(noe_RsEncoder *e, GFF4Element *data, int encode){
+    int i, j;
+    const int parityCount= e->parityCount;
+    const int dataCount= SIZE - parityCount - 1;
+    char stats[SIZE];
+    GFF4Element t= data[ e->tLocation ];
+
+    memset(stats, 0, sizeof(stats));
+
+    for(i=0; i<parityCount; i++){
+        GFF4Element p= diff(data[i+dataCount], prod(e->tPoly[i], t));
+
+        stats[ prod(MINUS1 - p, inv(e->tPoly[i])) ]|=1; //FIXME get rid of the neg
+    }
+
+    if(encode){
+        for(j=i=0; i<SIZE; i++){
+            if(stats[i]) 
+                continue;
+            if(j == t)
+                break;
+            j++;
         }
+        if(i>=SIZE)
+            return -1;
     }else{
-        memset(omega + parityCount + 1 + 1, 0, (SIZE - 2 - parityCount)*sizeof(GFF4Element));
-        EXT(gfft)(omega+1, omega+1, SHIFT);
-        for(i=0; i<parityCount; i++){
-            int index= SIZE - dataCount - i - 1;
-            
-            GFF4Element parity= prod(omega[bitReverse(index)+1], factor[i]);
-            
-            data[dataCount + i]= parity;
+        for(j=i=0; j<t; j++){
+            if(stats[j]) 
+                continue;
+            i++;
         }
     }
+    data[ e->tLocation ]= i;
+    t= diff(i, t);
+
+    for(i=0; i<parityCount; i++){
+        data[i+dataCount]= sum(data[i+dataCount], prod(e->tPoly[i], t));
+    }
+    return 0;
 }
 #endif
 
+#ifdef ALTERNATIVE_SOLVERS
 #define XCHG(a,b) {void *tmp= a; a= b; b= tmp;}
 
 static void hgcd(GFF4Element *ma[2], GFF4Element *mb[2], GFF4Element *a, GFF4Element *b, int len){
@@ -221,7 +200,7 @@ static void hgcd(GFF4Element *ma[2], GFF
         SET_POLY0(mb[1], 1);
 
         for(;;){
-            GFF4Element q= SIZE - prod(a[ a[0]+1 ], inv(b[ b[0]+1 ]));
+            GFF4Element q= neg(prod(a[ a[0]+1 ], inv(b[ b[0]+1 ])));
             int s= a[0] - b[0];
 
             assert(s >= 0);
@@ -333,6 +312,7 @@ static int rsEuclid(GFF4Element *locator
         EXT(diffPoly)(locator[di], locator[di], quot);
     }
 }
+#endif
 
 /**
  * ...
@@ -354,13 +334,9 @@ static int rsBerlekampMassey(GFF4Element
     for(k=1+erasureCount; k<=parityCount; k++){
         const int k2= k - erasureCount;
         GFF4Element delta= syn[k];
-
         for(j=1; j<=Lix[-1]; j++)
-#if M!=0
-            delta^= prod(Lix[j], syn[k-j]);
-#else
-            delta+= prod(Lix[j], syn[k-j]);
-#endif
+            delta ADD_OP prod(Lix[j], syn[k-j]);
+
         delta= reduce(delta);
 
         if(delta){
@@ -369,18 +345,22 @@ static int rsBerlekampMassey(GFF4Element
                 for(j=Tshift; j<=Tix[-1]; j++)
                     Lix[j]= sum(Lix[j], prod(factor, Tix[j-Tshift]));
             }else{
-                //Lix[-1] >= Tshift-1 has always been true during tests
                 //Tshift == 2 true in >99% of the cases
                 if(Tshift == 2 && Lix[-1]>0){
                     for(j=Tix[-1]; j>Lix[-1]; j--)
                         Tix[j]=             prod(factor, Tix[j-2]);
                     for(; j>=2; j--)
                         Tix[j]= sum(Lix[j], prod(factor, Tix[j-2]));
-                }else{
-                    for(j=Tix[-1]; j>NOE_MAX(Lix[-1], Tshift-1); j--)
+                }else if(Lix[-1] >= Tshift-1){
+                    for(j=Tix[-1]; j>Lix[-1]; j--)
                         Tix[j]=             prod(factor, Tix[j-Tshift]);
                     for(; j>=Tshift; j--)
                         Tix[j]= sum(Lix[j], prod(factor, Tix[j-Tshift]));
+                }else{
+                    for(j=Tix[-1]; j>=Tshift; j--)
+                        Tix[j]=             prod(factor, Tix[j-Tshift]);
+                    for(; j>Lix[-1]; j--)
+                        Tix[j]= 0;
                 }
                 for(; j>0; j--)
                     Tix[j]= Lix[j];
@@ -404,53 +384,49 @@ static int rsBerlekampMassey(GFF4Element
 
 int EXT(rsDecode)(GFF4Element *data, int *erasure, int erasureCount, int parityCount){
     int i;
-    int errorCount, gfftEval, phantomErrorCount;
+    int errorCount, gfftEval;
     const int dataCount= SIZE - 1 - parityCount;
     GFF4Element erasureLocator[erasureCount + 1];
-    GFF4Element erasureSynthSrc[erasureCount*3 + 1];
     GFF4Element locator0[SIZE]; //[parityCount + 1 + 1];
     GFF4Element locator1[SIZE]; //[parityCount + 1 + 1];
     GFF4Element *locators[2]={locator0, locator1};
     GFF4Element omega0[SIZE - 1 + 1];
     GFF4Element omega1[SIZE - 1 + 1];
     GFF4Element *omegas[2]={omega0, omega1}; //FIXME clean this shit up
-    GFF4Element *fErrorLocator, *errorLocator, *omega, *syn, *psi;
+    GFF4Element *errorLocator, *omega, *syn, *psi;
 
     syn= omegas[1];
     
     /* kill erased symbols */
     for(i=0; i<erasureCount; i++){
-//printf("%d\n", erasure[i]);
         data[ erasure[i] ]=0;
     }
 
     EXT(getSyndrom)(syn, data, parityCount);
+
     for(i=1; i<=parityCount; i++){
         if(syn[i+1]) break;
     }
     if(i>parityCount)
         return 0;
 
-    phantomErrorCount=0;
     //FIXME check truncated symbols syndrom
 
-    for(i=0; i<erasureCount; i++){
-        SET_POLY1(erasureSynthSrc + 3*i, 1, neg(EXT(exp)[erasure[i]]));
+    if(erasureCount){
+        EXT(synthPoly)(erasureLocator, erasure, erasureCount);
+        EXT(partialProdPoly)(syn, erasureLocator, syn, syn[0]);
     }
-    EXT(synthPoly)(erasureLocator, erasureSynthSrc, erasureCount);
-//  EXT(printPoly)(erasureLocator, erasureCount);
-#if 1
+#if 0
 for(i=0; i<erasureCount; i++){
     int eval= EXT(evalPoly)(erasureLocator, inv(EXT(exp)[erasure[i]]));
     assert( eval==0);
 }
 #endif
-    EXT(partialProdPoly)(syn, erasureLocator, syn, syn[0]);
-START_TIMER
+//START_TIMER
     i= rsBerlekampMassey(locators, omegas, parityCount, erasureCount);
 //    i= rsEuclid(locators, omegas, parityCount, erasureCount);
 //    i= rsSchoenhage(locators, omegas, parityCount, erasureCount);
-STOP_TIMER
+//STOP_TIMER
 //printf("errorCount %d \n", errorCount);
 
     if(i<0) return -1;
@@ -458,20 +434,23 @@ STOP_TIMER
 
     omega= omegas[i];
     errorLocator= locators[i];
-    fErrorLocator= omegas[1-i]; //reuse some unused space
     psi= locators[1-i]; //reuse some unused space
 
-    gfftEval= (errorCount>20) + (errorCount>1024); //FIXME finetune thresholds
-    if(M) gfftEval=0;
-    if(gfftEval){
-        memcpy(fErrorLocator, errorLocator+1, (errorCount+1)*sizeof(GFF4Element));
-        memset(fErrorLocator + errorCount+1, 0, (SIZE - errorCount-2)*sizeof(GFF4Element));
-        EXT(gfft)(fErrorLocator, fErrorLocator, SHIFT);
-    }
+    if(M)                                   gfftEval= 0;
+    else if(errorCount + erasureCount>1024) gfftEval= 2;
+    else                                    gfftEval= errorCount>20;
 
-    EXT(prodPoly)(psi, errorLocator, erasureLocator);
+    if(erasureCount)
+        EXT(prodPoly)(psi, errorLocator, erasureLocator);
+    else
+        memcpy(psi, errorLocator, (errorLocator[0]+2)*sizeof(GFF4Element));
     EXT(getDerivative)(psi, psi);
 
+    if(gfftEval){
+        memset(errorLocator + errorCount+2, 0, (SIZE - errorCount-2)*sizeof(GFF4Element));
+        errorLocator++;
+        EXT(gfft)(errorLocator, errorLocator, SHIFT);
+    }
     if(gfftEval>1){
         memset(omega + omega[0] + 2, 0, (SIZE - omega[0] - 2)*sizeof(GFF4Element));
         memset(psi   + psi  [0] + 2, 0, (SIZE - psi  [0] - 2)*sizeof(GFF4Element));
@@ -487,7 +466,7 @@ STOP_TIMER
             GFF4Element X, nom, denom;
 
             if(gfftEval){ //FIXME optimize
-                if(fErrorLocator[i])
+                if(errorLocator[i])
                     continue;
                     
                 r= i ? MINUS1 -  bitReverse(i) : 0; //FIXME i=0 <->MINUS1 ?
@@ -511,12 +490,8 @@ STOP_TIMER
 
             assert(r>=0 && r<=SIZE-2);
             e= prod(nom, inv(denom));
-            if(data[r]==0 && e==MINUS1){
-                if(r < dataCount)
-                    return -1;
-                else
-                    phantomErrorCount++;
-            }
+            if(data[r]==0 && e==MINUS1 && r < dataCount)
+                return -1;
             data[r]= sum(data[r], e);
         }
     }
@@ -542,6 +517,6 @@ STOP_TIMER
         data[r]= prod(nom, inv(denom));
     }
 
-    return erasureCount + errorCount - phantomErrorCount;
+    return erasureCount + errorCount;
 }
 

Modified: trunk/noe/rs.h
==============================================================================
--- trunk/noe/rs.h	(original)
+++ trunk/noe/rs.h	Tue Jul 10 00:53:52 2007
@@ -18,19 +18,17 @@
 
 typedef struct noe_RsEncoder{
     int parityCount;
-    GFF4Element *parityLocator;
-    GFF4Element *parityFactor;
+    int tLocation;
+    GFF4Element *tPoly;
 }noe_RsEncoder;
 
-void EXT(getRsGenerator)(GFF4Element *dst, int first, int order);
-
 /**
  * gets the syndroms.
  * @param src a 2^16 entry long polynom
  * @param syn the syndrom polynom will be stored here
  */
 void EXT(getSyndrom)(GFF4Element *syn, GFF4Element *src, int order);
-noe_RsEncoder *EXT(getRsEncoder)(int parityOrder);
+noe_RsEncoder *EXT(getRsEncoder)(int parityOrder, int tLocation);
 void EXT(freeRsEncoder)(noe_RsEncoder *encoder);
 
 /**
@@ -39,3 +37,6 @@ void EXT(freeRsEncoder)(noe_RsEncoder *e
 void EXT(rsEncode)(noe_RsEncoder *encoder, GFF4Element *data);
 
 int EXT(rsDecode)(GFF4Element *data, int *erased, int erasedCount, int parityCount);
+#if M==0
+int EXT(rsTransform)(noe_RsEncoder *e, GFF4Element *data, int encode);
+#endif

Modified: trunk/noe/test.c
==============================================================================
--- trunk/noe/test.c	(original)
+++ trunk/noe/test.c	Tue Jul 10 00:53:52 2007
@@ -124,6 +124,7 @@ int main(){
 #if 1
     srand (31415);
 
+#if 0 // only func using getRsGenerator so droped with it
     printf("\nTesting generator polynoms");
     for(i=1; i<100; i+=5){
         int n= SIZE - 1;
@@ -151,9 +152,9 @@ int main(){
         }
         printf("."); fflush(stdout);
     }
-
+#endif
     printf("\nTesting encoder");
-    for(i=1; i<1000 && i<SIZE/2; i+=i+5){
+    for(i=1; i<SIZE/2; i+=i+5){
         int n= SIZE - 1;
         int k= SIZE - 1 - i;
         GFF4Element syn[n+1];
@@ -162,17 +163,19 @@ int main(){
         GFF4Element code[n];
         noe_RsEncoder *encoder;
         int pass=5;
+        int tLoc= rand() % k;
 
-        encoder= EXT(getRsEncoder)(i);
+        encoder= EXT(getRsEncoder)(i, tLoc);
 
         for(pass=5; pass; pass--){
             for(j=0; j<k; j++)
                 data[j]= rand() & MASK;
+            data[tLoc] %= SIZE - 1 - n + k;
 
             memcpy(code, data, n*sizeof(GFF4Element));
-
+{START_TIMER
             EXT(rsEncode)(encoder, code);
-
+STOP_TIMER}
             EXT(getSyndrom)(syn, code, i);
             //FIXME check that code contains data, but its not touched so ...
 
@@ -182,6 +185,27 @@ int main(){
             for(j=0; j<i; j++){
                 if(syn[j+1+1]) printf("FAIL generator:%d coefficient:%d is %X\n", i, j, syn[j+1]);
             }
+#if M==0
+            if(EXT(rsTransform)(encoder, code, 1) < 0)
+                printf("transform failure\n");
+            EXT(getSyndrom)(syn, code, i);
+            for(j=0; j<i; j++)
+                if(syn[j+1+1]) printf("FAILT generator:%d coefficient:%d is %X\n", i, j, syn[j+1]);
+            for(j=0; j<n; j++){
+                if(code[j] == MINUS1)
+                    printf("illegal symbol detected\n");
+            }
+
+            if(EXT(rsTransform)(encoder, code, 0) < 0)
+                printf("inv-transform failure\n");
+            EXT(getSyndrom)(syn, code, i);
+            for(j=0; j<i; j++)
+                if(syn[j+1+1]) printf("FAILTT generator:%d coefficient:%d is %X\n", i, j, syn[j+1]);
+            for(j=0; j<k; j++){
+                if(data[j] != code[j])
+                    printf("data changed\n");
+            }
+#endif
             printf("."); fflush(stdout);
         }
         EXT(freeRsEncoder)(encoder);
@@ -192,13 +216,13 @@ int main(){
         int n= SIZE - 1;
         int k= SIZE - 1 - i;
 //      GFF4Element syn[n];
-        GFF4Element erased[i];
+        int erased[i];
         GFF4Element data[k];
         GFF4Element code[n];
         noe_RsEncoder *encoder;
         int pass=5;
 
-        encoder= EXT(getRsEncoder)(i);
+        encoder= EXT(getRsEncoder)(i, 0);
 
         for(pass=5; pass; pass--){
             int erasedCount, errorCount;
@@ -281,53 +305,6 @@ int main(){
     }
 #endif
     printf("\n");
-return 0;
-    printf("Testing erasure decoding\n");
-    for(i=1; i<20; i++){
-        int n= SIZE - 1;
-        int k= SIZE - 1 - i;
-        GFF4Element syn[n+1];
-        GFF4Element gen[i+1+1];
-        GFF4Element data[k+1];
-        GFF4Element code[n+1];
-        GFF4Element locator[i+1+1];
-        GFF4Element omega[2*i+1];
-
-        data[0]= k-1;
-        for(j=0; j<k; j++)
-            data[j+1]= rand() & MASK;
-        
-        EXT(getRsGenerator)(gen, 1, i);
-        EXT(prodPoly)(code, gen, data);
-        assert(code[0]==n);
-        
-        locator[0] = 0;
-        locator[1] = 1;
-        for(j=0; j<i; j++){
-            int index;
-            GFF4Element locationFactor[3];
-            
-            if((rand()&7)==0 && j) continue;
-            
-            index= rand()%MASK;
-            SET_POLY1(locationFactor, 1, inv(EXT(exp)[index]));
-            EXT(prodPoly)(locator, locator, locationFactor);
-            
-            code[index] = 0;
-        }
-        
-        EXT(getSyndrom)(syn, code+1, i);
-        EXT(prodPoly)(omega, syn, locator);
-//      EXT(getDeriative)(locator, locator, locatorOrder);
-        
-
-//      EXT(printPoly)(gen, i);
-//      EXT(printPoly)(syn, i-1);
-        
-        for(j=0; j<i; j++){
-            if(syn[j+1]) printf("FAIL generator:%d coefficient:%d\n", i, j);
-        }
-    }
 
     return 0;
 }



More information about the Mndiff-dev mailing list