gnuradio / grdigital / examples / berawgn.py @ master
History  View  Annotate  Download (4.78 KB)
1 
#!/usr/bin/env python


2 
#

3 
# Copyright 2012,2013 Free Software Foundation, Inc.

4 
#

5 
# This file is part of GNU Radio

6 
#

7 
# GNU Radio is free software; you can redistribute it and/or modify

8 
# it under the terms of the GNU General Public License as published by

9 
# the Free Software Foundation; either version 3, or (at your option)

10 
# any later version.

11 
#

12 
# GNU Radio is distributed in the hope that it will be useful,

13 
# but WITHOUT ANY WARRANTY; without even the implied warranty of

14 
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the

15 
# GNU General Public License for more details.

16 
#

17 
# You should have received a copy of the GNU General Public License

18 
# along with GNU Radio; see the file COPYING. If not, write to

19 
# the Free Software Foundation, Inc., 51 Franklin Street,

20 
# Boston, MA 021101301, USA.

21 
#

22  
23 
"""

24 
BER simulation for QPSK signals, compare to theoretical values.

25 
Change the N_BITS value to simulate more bits per Eb/N0 value,

26 
thus allowing to check for lower BER values.

27 

28 
Lower values will work faster, higher values will use a lot of RAM.

29 
Also, this app isn't highly optimizedthe flow graph is completely

30 
reinstantiated for every Eb/N0 value.

31 
Of course, expect the maximum value for BER to be one order of

32 
magnitude below what you chose for N_BITS.

33 
"""

34  
35  
36 
import math 
37 
import numpy 
38 
from gnuradio import gr, digital 
39 
from gnuradio import analog 
40 
from gnuradio import blocks 
41 
import sys 
42  
43 
try:

44 
from scipy.special import erfc 
45 
except ImportError: 
46 
print "Error: could not import scipy (http://www.scipy.org/)" 
47 
sys.exit(1)

48  
49 
try:

50 
import pylab 
51 
except ImportError: 
52 
print "Error: could not import pylab (http://matplotlib.sourceforge.net/)" 
53 
sys.exit(1)

54  
55 
# Best to choose powers of 10

56 
N_BITS = 1e7

57 
RAND_SEED = 42

58  
59 
def berawgn(EbN0): 
60 
""" Calculates theoretical bit error rate in AWGN (for BPSK and given Eb/N0) """

61 
return 0.5 * erfc(math.sqrt(10**(float(EbN0)/10))) 
62  
63 
class BitErrors(gr.hier_block2): 
64 
""" Two inputs: true and received bits. We compare them and

65 
add up the number of incorrect bits. Because integrate_ff()

66 
can only add up a certain number of values, the output is

67 
not a scalar, but a sequence of values, the sum of which is

68 
the BER. """

69 
def __init__(self, bits_per_byte): 
70 
gr.hier_block2.__init__(self, "BitErrors", 
71 
gr.io_signature(2, 2, gr.sizeof_char), 
72 
gr.io_signature(1, 1, gr.sizeof_int)) 
73  
74 
# Bit comparison

75 
comp = blocks.xor_bb() 
76 
intdump_decim = 100000

77 
if N_BITS < intdump_decim:

78 
intdump_decim = int(N_BITS)

79 
self.connect(self, 
80 
comp, 
81 
blocks.unpack_k_bits_bb(bits_per_byte), 
82 
blocks.uchar_to_float(), 
83 
blocks.integrate_ff(intdump_decim), 
84 
blocks.multiply_const_ff(1.0/N_BITS),

85 
self)

86 
self.connect((self, 1), (comp, 1)) 
87  
88 
class BERAWGNSimu(gr.top_block): 
89 
" This contains the simulation flow graph "

90 
def __init__(self, EbN0): 
91 
gr.top_block.__init__(self)

92 
self.const = digital.qpsk_constellation()

93 
# Source is N_BITS bits, nonrepeated

94 
data = map(int, numpy.random.randint(0, self.const.arity(), N_BITS/self.const.bits_per_symbol())) 
95 
src = blocks.vector_source_b(data, False)

96 
mod = digital.chunks_to_symbols_bc((self.const.points()), 1) 
97 
add = blocks.add_vcc() 
98 
noise = analog.noise_source_c(analog.GR_GAUSSIAN, 
99 
self.EbN0_to_noise_voltage(EbN0),

100 
RAND_SEED) 
101 
demod = digital.constellation_decoder_cb(self.const.base())

102 
ber = BitErrors(self.const.bits_per_symbol())

103 
self.sink = blocks.vector_sink_f()

104 
self.connect(src, mod, add, demod, ber, self.sink) 
105 
self.connect(noise, (add, 1)) 
106 
self.connect(src, (ber, 1)) 
107  
108 
def EbN0_to_noise_voltage(self, EbN0): 
109 
""" Converts Eb/N0 to a singlesided noise voltage (assuming unit symbol power) """

110 
return 1.0 / math.sqrt(2.0 * self.const.bits_per_symbol() * 10**(float(EbN0)/10)) 
111  
112  
113 
def simulate_ber(EbN0): 
114 
""" All the work's done here: create flow graph, run, read out BER """

115 
print "Eb/N0 = %d dB" % EbN0 
116 
fg = BERAWGNSimu(EbN0) 
117 
fg.run() 
118 
return numpy.sum(fg.sink.data())

119  
120 
if __name__ == "__main__": 
121 
EbN0_min = 0

122 
EbN0_max = 15

123 
EbN0_range = range(EbN0_min, EbN0_max+1) 
124 
ber_theory = [berawgn(x) for x in EbN0_range] 
125 
print "Simulating..." 
126 
ber_simu = [simulate_ber(x) for x in EbN0_range] 
127  
128 
f = pylab.figure() 
129 
s = f.add_subplot(1,1,1) 
130 
s.semilogy(EbN0_range, ber_theory, 'g.', label="Theoretical") 
131 
s.semilogy(EbN0_range, ber_simu, 'bo', label="Simulated") 
132 
s.set_title('BER Simulation')

133 
s.set_xlabel('Eb/N0 (dB)')

134 
s.set_ylabel('BER')

135 
s.legend() 
136 
s.grid() 
137 
pylab.show() 
138 