## Python problem: class returns by reference, not by value

I made some changes to my python-class Gold.

The problem is that in python an attribute of a class is returned by reference. This might be nasty in some situations.

if I modify prn I modify the internal list prnRegAnalogue because they are identical: they have the same address!

#import old class Gold goldutil myGold1 = Gold(13) prn = myGold1. prnRegAnalogue prn[0] = 100 prn2 = myGold1. prnRegAnalogue

if I modify prn I modify the internal list prnRegAnalogue because they are identical: they have the same address!

So now prn[0] and prn2[0] both have the same value 100.

Therefore I changed my class. Now I hide prnRegAnalogue by prefixing the name with a double-underscore: __ prnRegAnalogue and I defined a method prnReg() to return a copy of prnRegAnalogue. The effect is that now the list prnRegAnalogue is called by value.

# import modified class Gold goldutil myGold1 = Gold(13) prn = myGold1. prnRegAnalogue prn[0] = 100 prn2 = myGold1. prnRegAnalogue

Now prn2[0] has not been changed.

The new class now reads:

class Gold: # Gold-sequence goldutil.py # Nov-2013 Kees de Groot # call: # myprn1 = Gold(3, 7) for code phase selection = (3,7) # or # myprn2 = Gold(4) for the 4th from the list # = satellite ID number = GPS PRN signal number # start with (0,0) to start numbering with 1 # note that PRN 34 and 37 are identical prnlist = [ (0,0), (2,6), (3,7), (4,8), (5,9), (1,9), (2,10), (1,8), (2,9), (3,10), (2,3), (3,4), (5,6), (6,7), (7,8), (8,9), (9,10), (1,4), (2,5), (3,6), (4,7), (5,8), (6,9), (1,3), (4,6), (5,7), (6,8), (7,9), (8,10), (1,6), (2,7), (3,8), (4,9), (5,10), (4,10), (1,7), (2,8), (4,10)] def __init__(self, a, b = 0): self.reg1 = [1]*11 # don't use first bit-position [0] self.reg2 = [1]*11 # so count bit-positions from 1 to 10 self.__prnRegBinary = [] # contains 0 and 1 self.__prnRegAnalogue = [] # contains -1 and +1 if b == 0: # a = GPS PRN signal number or satellite ID (self.prn1, self.prn2) = self.prnlist[a] else: # direct code phase selection (self.prn1, self.prn2) = (a, b) for i in range(1024): G1 = self.reg1[10] G2 = self.reg2[self.prn1] ^ self.reg2[self.prn2] out = G1 ^ G2 self.__prnRegBinary.append(out) # values 0 and 1 self.__prnRegAnalogue.append(2*out-1) # values -1 and +1 val1 = self.reg1[3] ^ self.reg1[10] val2 = (self.reg2[2] ^ self.reg2[3] ^ self.reg2[6] ^ self.reg2[8] ^ self.reg2[9] ^ self.reg2[10]) # shift one position to the right for j in range(9): # j = 0,1,2,3,4,5,6,7,8 k = 10 - j # k = 10,9,8,7,6,5,4,3,2 self.reg1[k] = self.reg1[k-1] self.reg2[k] = self.reg2[k-1] self.reg1[1] = val1 self.reg2[1] = val2 # return list by value, not by reference def prnReg(self): alist = [] for val in self.__prnRegAnalogue: alist.append(val) return alist

## Correlation with fft

A very interesting property of the Fourier transform is that convolution in the time-domain can be done by multiplication in the frequency-domain.

That means that to calculate a @ b, in which @ is meant to stand for convolution, you can do that calculation as follows:

ifft( fft(a) * fft(b) )which means do an fft of a, do an fft of b, multiply and inverse-fft the result.

But the result of fft(a) has to be complex conjugated first, which means that all imaginary values has to change sign.

So, the procedure is: ifft(conjugate(fft(a)) * fft(b))

It sound strange perhaps, but this complex (sic) procedure might be very efficient to calculate because in practice the part conjugate(fft(a)) has to be calculated only once during initialization. The remaining calculation boils down to

c = conjugate(fft(a)

…

ifft(c * fft(b))

To calculate the correlation between a and b.

I made a small python program to experiment with this procedure:

This program might look daunting, but there is really a lot of comment and not so much code.import goldutil as gl # home-made Gold-sequence generator import pylab as pl # module for plotting import numpy as np # module for fft # CorrelFFTdemo2s.py # Simulation of message constructed using Gold-sequences # - Define an 8-bit message: eg 0 0 1 1 0 1 0 0 1 = -1, 1, 1, -1, 1, -1, -1, 1 # - Construct an array with 9*1024 samples of the same Gold-sequence, # eg. satellite id 12 # - Multiply the first 1024 samples with 1, the 2n set with -1, # the 3rd set with -1, etc to code our message # - Now try to decode our message using correlation implemented with FFT # nov 2013 Kees de Groot messageBin = (0, 0, 1, 1, 0, 1, 0, 0, 1) print "messageBin = ", messageBin messageAna = [] for val in messageBin: messageAna.append(2*val-1) print "messageAna = ", messageAna sat1 = 12 # this is the sat we are looking for myGold1 = gl.Gold(sat1) pprn1 = myGold1.prnReg() # get list with 1024 bits of Goldsequence # code message into codedmsg by multiplying codedmsg = [] for mesbit in messageAna: for prnbit in pprn1: codedmsg.append(mesbit * prnbit) codedmsg = [0]*200 + codedmsg # add some zeroes to improve readability of plot # decode codedmsg by correlation with origional pprn1 Gold-sequence # cross-correlation with the right Gold-sequence # correlate codedmsg with pprn1 and put result in ac # once-only calculate complex conjugate of prn-sequence pprn1 # problem is that I want the fft to have the same length as codedmsg # so, just add zeros at the end of pprn1 numzeros = len(codedmsg) - len(pprn1) pprn1 = pprn1 + [0]*numzeros fftprn = np.fft.rfft(pprn1) # complex conjugate (there must be a better method) for i in range(len(fftprn)): fftprn[i] = fftprn[i].conjugate() fftcodedmsg = np.fft.rfft(codedmsg) # convolution (correlation) == multiplying with complex conjugate # after fft and inverse-fft the result # multiply fftcodedmsg with fftprn for i in range(len(fftcodedmsg)): fftcodedmsg[i] *= fftprn[i] # inverse-fft ac = np.fft.irfft(fftcodedmsg) pl.figure(1) pl.title("correlation, calculated with fft-method") pl.ylim(-1500, 1500) pl.stem(ac) pl.show()

The output is:

GraphicsPython 2.7.3 (default, Apr 10 2012, 23:31:26) [MSC v.1500 32 bit (Intel)] on win32 Type "copyright", "credits" or "license()" for more information. >>> ================================ RESTART ================================ >>> messageBin = (0, 0, 1, 1, 0, 1, 0, 0, 1) messageAna = [-1, -1, 1, 1, -1, 1, -1, -1, 1]

It is obvious that this is exactly the same result as in the previous blog-message where I programmed plain convolution.

## Correlation with fft, block-mode

I expect to handle my real-time samples in blocks of 1024 chips. A chip is one bit of a Gold-sequence. In this example I simulated this procedure. In a loop I correlate 1024 chips from my input and decode the message bit by bit.

with output:import goldutil as gl # home-made Gold-sequence generator import pylab as pl # module for plotting import numpy as np # module for fft # CorrelFFTblock.py # Simulation of message constructed using Gold-sequences # - Define an 8-bit message: eg 0 0 1 1 0 1 0 0 1 = -1, 1, 1, -1, 1, -1, -1, 1 # - Construct an array with 9*1024 samples of the same Gold-sequence, # eg. satellite id 12 # - Multiply the first 1024 samples with 1, the 2n set with -1, # the 3rd set with -1, etc to code our message # - Now try to decode our message using correlation implemented with FFT # do this piece-wise in blocks of 1024 bits # nov 2013 Kees de Groot def minmax(ac): minval = min(ac) maxval = max(ac) if abs(minval) > abs(maxval): result = 0 else: result = 1 return result messageBin = (0, 0, 1, 1, 0, 1, 0, 0, 1) print "messageBin = ", messageBin messageAna = [] for val in messageBin: messageAna.append(2*val-1) print "messageAna = ", messageAna sat1 = 12 # this is the sat we are looking for myGold1 = gl.Gold(sat1) pprn1 = myGold1.prnReg() # get list with 1024 bits of Goldsequence # code message into codedmsg by multiplying codedmsg = [] for mesbit in messageAna: for prnbit in pprn1: codedmsg.append(mesbit * prnbit) codedmsg = [0]*200 + codedmsg + [0]*824 # add some zeroes # decode codedmsg by correlation with origional pprn1 Gold-sequence # cross-correlation with the right Gold-sequence # correlate codedmsg with pprn1 and put result in ac # once-only calculate complex conjugate of prn-sequence pprn1 fftprn = np.fft.rfft(pprn1) # complex conjugate (there must be a better method) for i in range(len(fftprn)): fftprn[i] = fftprn[i].conjugate() # the above calculation has to be done once during initialization # now take 1024 bit blocks and correlate using fft N = len(codedmsg) / len(pprn1) print "N = ", N for n in range(N): fftcodedmsg = np.fft.rfft(codedmsg[n*1024:(n+1)*1024]) for i in range(len(fftcodedmsg)): fftcodedmsg[i] *= fftprn[i] ac = np.fft.irfft(fftcodedmsg) print "block ", n, ", message-bit = ", minmax(ac) pl.figure(n+1) pl.title("correlation, calculated with fft-method " + "block " + str(n)) pl.ylim(-1500, 1500) pl.stem(ac) pl.show()

and a lot of nice graphics:Python 2.7.3 (default, Apr 10 2012, 23:31:26) [MSC v.1500 32 bit (Intel)] on win32 Type "copyright", "credits" or "license()" for more information. >>> ================================ RESTART ================================ >>> messageBin = (0, 0, 1, 1, 0, 1, 0, 0, 1) messageAna = [-1, -1, 1, 1, -1, 1, -1, -1, 1] N = 10 block 0 , message-bit = 0 block 1 , message-bit = 0 block 2 , message-bit = 1 block 3 , message-bit = 1 block 4 , message-bit = 0 block 5 , message-bit = 1 block 6 , message-bit = 0 block 7 , message-bit = 0 block 8 , message-bit = 1 block 9 , message-bit = 1

## Conclusion

To find a message coded with a Gold-sequence one uses correlation. One method is to use convolution. It is a simple process but might cost a lot of time.

Another method is to make use of a property of the Fourier-transform:

Convolution in the time-domain is the same as multiplication in the frequency-domainI experimented with this method and showed it to be possible.

Now I am confident to use this method with real samples from the bladeRF.

Keep tuned!