code:python_reedsolomon
This is an old revision of the document!
#!/usr/bin/python mm=7 # RS code over GF(2**4) - change to suit tt=16 # number of errors that can be corrected nn=2**mm -1 #length of codeword kk = nn-2*tt primative_polynomials=[[1,1,0,1],[1,1,0,0,1],[1,0,1,0,0,1],[1,1,0,0,0,0,1],[1,0,0,1,0,0,0,1],[1,0,1,1,1,0,0,0,1],[1,0,0,0,1,0,0,0,0,1]] pp=primative_polynomials[mm-3] # specify irreducible polynomial coeffts m=5 alpha_to=[0]*(nn+1) index_of=[0]*(nn+1) gg=[0]*(nn-kk+1) recd=[0]*(nn) data=[0]*(kk) bb=[0]*(nn-kk) def generate_gf(): # generate GF(2**mm) from the irreducible polynomial p(X) in pp[0]..pp[mm] # lookup tables: index->polynomial form alpha_to[] contains j=alpha**i; # polynomial form -> index form index_of[j=alpha**i] = i # alpha=2 is the primitive element of GF(2**mm) mask = 1 alpha_to[mm] = 0 for i in range(mm): alpha_to[i] = mask index_of[alpha_to[i]] = i if pp[i]!=0: alpha_to[mm] ^= mask mask <<= 1 index_of[alpha_to[mm]] = mm mask >>= 1 for i in range(mm+1,nn): if (alpha_to[i-1] >= mask): alpha_to[i] = alpha_to[mm] ^ ((alpha_to[i-1]^mask)<<1) else: alpha_to[i] = alpha_to[i-1]<<1 index_of[alpha_to[i]] = i index_of[0] = -1 ; def gen_poly(): # Obtain the generator polynomial of the tt-error correcting, length # nn=(2**mm -1) Reed Solomon code from the product of (X+alpha**i), i=1..2*tt gg[0] = 2 # primitive element alpha = 2 for GF(2**mm) gg[1] = 1 # g(x) = (X+alpha) initially for i in range(2,nn-kk+1): gg[i] = 1 for j in range(i-1,0,-1): if (gg[j] != 0): gg[j] = gg[j-1]^ alpha_to[(index_of[gg[j]]+i)%nn] else: gg[j] = gg[j-1] gg[0] = alpha_to[(index_of[gg[0]]+i)%nn] ; # gg[0] can never be zero # convert gg[] to index form for quicker encoding */ for i in range(0,nn-kk): gg[i] = index_of[gg[i]] ; def encode_rs(): # take the string of symbols in data[i], i=0..(k-1) and encode systematically # to produce 2*tt parity symbols in bb[0]..bb[2*tt-1] # data[] is input and bb[] is output in polynomial form. # Encoding is done by using a feedback shift register with appropriate # connections specified by the elements of gg[], which was generated above. # Codeword is c(X) = data(X)*X**(nn-kk)+ b(X) for i in range(0,nn-kk): bb[i] = 0 for i in range(kk-1,-1,-1): feedback = index_of[data[i]^bb[nn-kk-1]] if feedback != -1: for j in range(nn-kk-1,0,-1): if gg[j] !=-1: bb[j] = bb[j-1]^alpha_to[(gg[j]+feedback)%nn] else: bb[j] = bb[j-1] bb[0] = alpha_to[(gg[0]+feedback)%nn] else: for j in range(nn-kk-1,0,-1): bb[j] = bb[j-1] bb[0] = 0 # generate the Galois Field GF(2**mm) generate_gf() print "Look-up tables for GF(2**%2d)\n",mm print " i alpha_to[i] index_of[i]\n" for i in range(0,nn): print i,alpha_to[i],index_of[i] print "\n\n" #compute the generator polynomial for this RS code gen_poly() #for known data, stick a few numbers into a zero codeword. Data is in # polynomial form. # for example, say we transmit the following message (nothing special!) s=raw_input("enter string ") d=0 for i in s: data[d]=ord(i) d=d+1 for i in range(len(s),kk): data[i]=0 # encode data[] to produce parity in bb[]. Data input and parity output # is in polynomial form encode_rs() # put the transmitted codeword, made up of data plus parity, in recd[] */ for i in range(0,nn-kk): recd[i] = bb[i] for i in range(0,kk): recd[i+nn-kk] = data[i] print recd print len(recd) s='' for i in range(len(recd)-1,-1,-1): if i>40 or i<17: s+=chr(recd[i]) else: s+=" " print s from reedsolomon import Codec c=Codec(nn,kk,mm) print c.decode(s,[86, 87, 88, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101])
code/python_reedsolomon.1202313967.txt.gz · Last modified: 2008/07/19 23:31 (external edit)