Statistics
| Branch: | Tag: | Revision:

gnuradio / gr-digital / 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 02110-1301, 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 optimized--the 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, non-repeated
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 single-sided 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, 'b-o', 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