UKHAS Wiki

UK High Altitude Society

User Tools

Site Tools


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)

Donate Powered by PHP Valid HTML5 Valid CSS Driven by DokuWiki