In [1]:
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
import matplotlib
# showing figures inline
%matplotlib inline
In [2]:
# plotting options
font = {'size' : 20}
plt.rc('font', **font)
plt.rc('text', usetex=True)
matplotlib.rc('figure', figsize=(18, 6) )
In [3]:
# max. numbers of errors and/or symbols
max_errors = 1e2
max_syms = 1e5
# Eb/N0
EbN0_db_min = 0
EbN0_db_max = 30
EbN0_db_step = 3
# initialize Eb/N0 array
EbN0_db_range = np.arange( EbN0_db_min, EbN0_db_max, EbN0_db_step )
EbN0_range = 10**( EbN0_db_range / 10 )
# constellation points
mod_points_bpsk = [-1, 1]
In [4]:
# initialize BER array
# theoretical ber for bpsk as on slides
ber_bpsk = np.zeros_like( EbN0_db_range, dtype=float )
ber_bpsk_theo = 1 - stats.norm.cdf( np.sqrt( 2 * EbN0_range ) )
# ber in fading channel
ber_fading = np.zeros_like( EbN0_db_range, dtype=float )
ber_fading_theo = 1 / ( 4 * 10**(EbN0_db_range / 10 ) )
# ber when applying channel inversion
ber_inverted = np.zeros_like( EbN0_db_range, dtype=float )
# loop for snr
for ind_snr, val_snr in enumerate( EbN0_range ):
# initialize error counter
num_errors_bpsk = 0
num_errors_fading = 0
num_errors_inverted = 0
num_syms = 0
# get noise variance
sigma2 = 1. / ( val_snr )
# loop for errors
while ( num_errors_bpsk < max_errors and num_syms < max_syms ):
# generate data and modulate by look-up
d = np.random.randint( 0, 2)
s = mod_points_bpsk[ d ]
###
# bpsk without fading
###
# add noise
noise = np.sqrt( sigma2 / 2 ) * ( np.random.randn() + 1j * np.random.randn() )
r_bpsk = s + noise
# demod
d_est_bpsk = int( np.real( r_bpsk ) > 0 )
###
# bpsk with slow flat fading
###
# flat fading
h = 1/np.sqrt(2) * ( np.random.randn() + 1j * np.random.randn() )
r_flat = h * s + noise
# matched filter and inverting channel
y_mf = np.conjugate( h / np.abs(h) )* r_flat
y_inv = r_flat / h
# demodulate symbols
d_est_flat = int( np.real( y_mf ) > 0 )
d_est_inv = int( np.real( y_inv ) > 0 )
###
# count errors
num_errors_bpsk += int( d_est_bpsk != d )
num_errors_fading += int( d_est_flat != d )
num_errors_inverted += int( d_est_inv != d )
# increase counter
num_syms += 1
# ber by relative amount of errors
ber_bpsk[ ind_snr ] = num_errors_bpsk / num_syms
ber_fading[ ind_snr ] = num_errors_fading / ( num_syms * 1.0 )
ber_inverted[ ind_snr ] = num_errors_inverted / ( num_syms * 1.0 )
# show progress
print('Eb/N0 planned (dB) = {:2.1f}\n'.format( 10*np.log10(val_snr) ) )
In [5]:
# plot bpsk results using identical colors for theory and simulation
ax_sim = plt.plot( EbN0_db_range, ber_bpsk, marker = 'o', mew=4, ms=18, markeredgecolor = 'none', linestyle='None', label='AWGN, sim.' )
color_sim = ax_sim[0].get_color()
plt.plot(EbN0_db_range, ber_bpsk_theo, linewidth = 2.0, color = color_sim, label='AWGN, theo.')
# plot slow flat results using identical colors for theory and simulation
ax_sim = plt.plot( EbN0_db_range, ber_fading , marker = 'D', mew=4, ms=18, markeredgecolor = 'none', linestyle='None', label = 'Fading, sim.' )
color_sim = ax_sim[0].get_color()
plt.plot(EbN0_db_range, ber_fading_theo, linewidth = 2.0, color = color_sim, label='Fading, theo.')
# plot ber when using channel inversion
ax_sim = plt.plot( EbN0_db_range, ber_inverted , marker = 'v', mew=4, ms=18, markeredgecolor = 'none', linestyle='None', label = 'Fading, inv., sim.' )
plt.yscale('log')
plt.grid(True)
plt.legend(loc='lower left')
plt.ylim( (1e-7, 1) )
plt.xlabel('$E_b/N_0$ (dB)')
plt.ylabel('BER')
Out[5]:
In [ ]: