23 February 2016 GNU Radio / Simulation

This IPython notebook generates channel coefficients with GNU Radio’s flat fader block and plots their autocorrelation at different points in time. The goal is to show how autocorrelation properties change, highlighting a potential bug in the model. See also the corresponding thread on the GNU Radio mailing list.

In [1]:

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
plt.style.use('bmh')
plt.rcParams['figure.figsize'] = (10, 10.0*2/3)
plt.rcParams['font.size'] = 14
plt.rcParams['lines.linewidth'] = '3'


N = 12              # N sinusoids
v = 400 / 3.6       # speed (km/h to m/s)
f = 6e9             # freq
c = 3e8             # speed of light
f_d = f * v / c     # max. Doppler
samp_rate = 1e6     # 10 MHz
N_samp = int(1000 * samp_rate)  # samples to generate
fDTs = f * v / c / samp_rate
seed = np.random.randint(10000)

Generate Channel CoefficientsÂ¶

In [2]:

from gnuradio import gr, blocks, channels, analog
sig = analog.sig_source_c(0, analog.GR_CONST_WAVE, 0, 0, 1)
head = blocks.head(gr.sizeof_gr_complex*1, N_samp)
chan = channels.fading_model(N, fDTs, False, 4.0, seed)
snk = blocks.file_sink(gr.sizeof_gr_complex*1, "/home/basti/channel.bin", False)
snk.set_unbuffered(True)

## do simulations
tb = gr.top_block()
tb.connect(sig, head)
tb.connect(head, chan)
tb.connect(chan, snk)
tb.run()

Autocorrelation over TimeÂ¶

In [3]:

import os

n_w = 4
w = int(n_w * c / f / v * samp_rate)

def autocorr(x, t=1):
    return np.corrcoef(np.array([x[0:len(x)-t], x[t:len(x)]]))

# resolution of autocorrelation
gaps = np.round(np.linspace(0, w, num=80))
gaps = gaps.astype(int)

## analytical
import scipy.special
bessel = map(lambda x: scipy.special.jv(0, 2 * np.pi * f_d / samp_rate * x), gaps)

N_items = 1000000

i = 0
while i < (N_samp - N_items):
    file = open("/home/basti/channel.bin", "rb")
    file.seek(i * 8, os.SEEK_SET)
    chan = np.fromfile(file, dtype=np.complex64, count=N_items)
    file.close()
    acorr = map(lambda x: np.abs(autocorr(chan, x)[0, 1]), gaps)
    
    fig, ax = plt.subplots()

    ax.plot(gaps/c*f*v/samp_rate, acorr, label='Simulations')
    ax.plot(gaps/c*f*v/samp_rate, np.abs(bessel), linestyle=':', label='Analytical')
    ax.axvline(0.4, linewidth=1, color='gray');
    ax.legend(loc=1)
    plt.show()
    
    i += 100000000

Values over timeÂ¶

To get a better understanding we plot the variables over time.

In [1]:

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
plt.style.use('bmh')
plt.rcParams['figure.figsize'] = (18, 15)
plt.rcParams['font.size'] = 14
plt.rcParams['lines.linewidth'] = '3'

d_step = .1
fDTs = 1
N_samp = 1600
n = 0
d_theta = 0
d_N = 8

thetas = []
alphas = []
dopplers = []
coeffs = []
diffs = []

def update_theta():
    global d_theta
    global d_step
    d_theta += np.random.random() * d_step

    if (d_theta > np.pi):
        d_theta = np.pi
        d_step = -d_step
    elif (d_theta < -np.pi):
        d_theta = -np.pi
        d_step = -d_step


tmp = 0

for d_m in range(N_samp):
    alpha_n = (2*np.pi*n - np.pi + d_theta)/(4*d_N)
    doppler =      2*np.pi*fDTs*d_m*np.cos(alpha_n)
    coeff = np.cos(2*np.pi*fDTs*d_m*np.cos(alpha_n))

    diffs.append(doppler - tmp)
    tmp = doppler

    thetas.append(d_theta)
    alphas.append(alpha_n)
    dopplers.append(doppler)
    coeffs.append(coeff)
    update_theta()

fig, ax = plt.subplots(4, sharex=True)
fig.set_size_inches(18.5, 15.5)

ax[0].plot(thetas)
ax[0].set_title("d_theta")
ax[1].plot(alphas)
ax[1].set_title("alpha_n")
ax[2].plot(diffs)
ax[2].set_title("phase diff")
ax[3].plot(coeffs)
ax[3].set_title("coeff")
ax[3].set_xlabel("d_m (time)");