In [ ]:
from bokeh.plotting import output_notebook
output_notebook()

In [ ]:
from simple_precoded_srs import *
from pyphysim.reference_signals.zadoffchu import get_shifted_root_seq

In [ ]:
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Scenario Description xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# 3 Base Stations, each sending data to its own user while interfering
# with the other users.

# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Configuration xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
num_prbs = 25  # Number of PRBs to simulate
Nsc = 12 * num_prbs  # Number of subcarriers
Nzc = 149  # Size of the sequence
u1 = 1  # Root sequence index of the first user
u2 = 2  # Root sequence index of the first user
u3 = 3  # Root sequence index of the first user
numAnAnt = 4  # Number of Base station antennas
numUeAnt = 2  # Number of UE antennas

num_samples = 1  # Number of simulated channel samples (from
# Jakes process)

# Channel configuration
speedTerminal = 0 / 3.6  # Speed in m/s
fcDbl = 2.6e9  # Central carrier frequency (in Hz)
timeTTIDbl = 1e-3  # Time of a single TTI
subcarrierBandDbl = 15e3  # Subcarrier bandwidth (in Hz)
numOfSubcarriersPRBInt = 12  # Number of subcarriers in each PRB
L = 16  # The number of rays for the Jakes model.

# Dependent parameters
lambdaDbl = 3e8 / fcDbl  # Carrier wave length
Fd = speedTerminal / lambdaDbl  # Doppler Frequency
Ts = 1. / (Nsc * subcarrierBandDbl)  # Sampling time

# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Generate the root sequence xxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
a_u1 = get_extended_ZF(calcBaseZC(Nzc, u1), Nsc / 2)
a_u2 = get_extended_ZF(calcBaseZC(Nzc, u2), Nsc / 2)
a_u3 = get_extended_ZF(calcBaseZC(Nzc, u3), Nsc / 2)

print("Nsc: {0}".format(Nsc))
print("a_u.shape: {0}".format(a_u1.shape))

# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Create shifted sequences for 3 users xxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# We arbitrarely choose some cyclic shift index and then we call
# zadoffchu.getShiftedZF to get the shifted sequence.
shift_index = 4
r1 = get_shifted_root_seq(a_u1, shift_index, 8)
r2 = get_shifted_root_seq(a_u2, shift_index, 8)
r3 = get_shifted_root_seq(a_u3, shift_index, 8)

# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Generate channels from users to the BS xxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
jakes_all_links = np.empty([3, 3], dtype=object)
tdlchannels_all_links = np.empty([3, 3], dtype=object)
impulse_responses = np.empty([3, 3], dtype=object)
# Dimension: `UEs x ANs x num_subcarriers x numUeAnt x numAnAnt`
freq_responses = np.empty([3, 3, Nsc, numUeAnt, numAnAnt], dtype=complex)

for ueIdx in range(3):
    for anIdx in range(3):
        jakes_all_links[ueIdx, anIdx] = JakesSampleGenerator(
            Fd, Ts, L, shape=(numUeAnt, numAnAnt))

        tdlchannels_all_links[ueIdx, anIdx] = TdlChannel(
            jakes_all_links[ueIdx, anIdx],
            tap_powers_dB=COST259_TUx.tap_powers_dB,
            tap_delays=COST259_TUx.tap_delays)

        tdlchannels_all_links[ueIdx, anIdx]._generate_impulse_response(
            num_samples)

        impulse_responses[ueIdx, anIdx] \
            = tdlchannels_all_links[ueIdx, anIdx].get_last_impulse_response()

        freq_responses[ueIdx, anIdx] = \
            impulse_responses[ueIdx, anIdx].get_freq_response(Nsc)[:, :, :, 0]

# xxxxxxxxxx Channels in downlink direction xxxxxxxxxxxxxxxxxxxxxxxxxxx
# Dimension: `Nsc x numUeAnt x numAnAnt`
dH11 = freq_responses[0, 0]
dH12 = freq_responses[0, 1]
dH13 = freq_responses[0, 2]
dH21 = freq_responses[1, 0]
dH22 = freq_responses[1, 1]
dH23 = freq_responses[1, 2]
dH31 = freq_responses[2, 0]
dH32 = freq_responses[2, 1]
dH33 = freq_responses[2, 2]

# xxxxxxxxxx Principal dimension in downlink direction xxxxxxxxxxxxxxxx
sc_idx = 124  # Index of the subcarrier we are interested in
[dU11, dS11, dV11_H] = np.linalg.svd(dH11[sc_idx])
[dU22, dS22, dV22_H] = np.linalg.svd(dH22[sc_idx])
[dU33, dS33, dV33_H] = np.linalg.svd(dH33[sc_idx])

# xxxxxxxxxx Users precoders xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# Users' precoders are the main column of the U matrix
F11 = dU11[:, 0].conj()
F22 = dU22[:, 0].conj()
F33 = dU33[:, 0].conj()

# xxxxxxxxxx Channels in uplink direction xxxxxxxxxxxxxxxxxxxxxxxxxxxxx
int_path_loss = 0.05  # Path loss of the interfering links
int_g = math.sqrt(int_path_loss)  # Gain of interfering links
dir_d = 1.0                   #  Gain of direct links

pl = np.array([[  2.21e-08,   2.14e-09,   1.88e-08],
               [  3.45e-10,   2.17e-08,   4.53e-10],
               [  4.38e-10,   8.04e-10,   4.75e-08]])
# pl = np.array([[  1,   0.1,   0.1],
#                [  0.1,   1,   0.1],
#                [  0.1,   0.1,   1]])


# Dimension: `Nsc x numAnAnt x numUeAnt`
uH11 = math.sqrt(pl[0, 0]) * np.transpose(dH11, axes=[0, 2, 1])
uH12 = math.sqrt(pl[0, 1]) * np.transpose(dH12, axes=[0, 2, 1])
uH13 = math.sqrt(pl[0, 2]) * np.transpose(dH13, axes=[0, 2, 1])
uH21 = math.sqrt(pl[1, 0]) * np.transpose(dH21, axes=[0, 2, 1])
uH22 = math.sqrt(pl[1, 1]) * np.transpose(dH22, axes=[0, 2, 1])
uH23 = math.sqrt(pl[1, 2]) * np.transpose(dH23, axes=[0, 2, 1])
uH31 = math.sqrt(pl[2, 0]) * np.transpose(dH31, axes=[0, 2, 1])
uH32 = math.sqrt(pl[2, 1]) * np.transpose(dH32, axes=[0, 2, 1])
uH33 = math.sqrt(pl[2, 2]) * np.transpose(dH33, axes=[0, 2, 1])

# Compute the equivalent uplink channels
uH11_eq = uH11.dot(F11)
uH12_eq = uH12.dot(F22)
uH13_eq = uH13.dot(F33)
uH21_eq = uH21.dot(F11)
uH22_eq = uH22.dot(F22)
uH23_eq = uH23.dot(F33)
uH31_eq = uH31.dot(F11)
uH32_eq = uH32.dot(F22)
uH33_eq = uH33.dot(F33)

# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Compute Received Signals xxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# Calculate the received signals
comb_indexes = np.r_[0:Nsc:2]
Y1_term11 = uH11_eq[comb_indexes] * r1[:, np.newaxis]
Y1_term12 = uH12_eq[comb_indexes] * r2[:, np.newaxis]
Y1_term13 = uH13_eq[comb_indexes] * r3[:, np.newaxis]
Y1 = Y1_term11 + Y1_term12 + Y1_term13

Y2_term21 = uH21_eq[comb_indexes] * r1[:, np.newaxis]
Y2_term22 = uH22_eq[comb_indexes] * r2[:, np.newaxis]
Y2_term23 = uH23_eq[comb_indexes] * r3[:, np.newaxis]
Y2 = Y2_term21 + Y2_term22 + Y2_term23

Y3_term31 = uH31_eq[comb_indexes] * r1[:, np.newaxis]
Y3_term32 = uH32_eq[comb_indexes] * r2[:, np.newaxis]
Y3_term33 = uH33_eq[comb_indexes] * r3[:, np.newaxis]
Y3 = Y3_term31 + Y3_term32 + Y3_term33

In [ ]:
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Estimate the equivalent channel xxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
(uH11_eq_est, uH12_eq_est, uH13_eq_est, uH21_eq_est, uH22_eq_est,
     uH23_eq_est, uH31_eq_est, uH32_eq_est, uH33_eq_est
    ) = estimate_channels_remove_only_direct(Y1, Y2, Y3, r1, r2, r3, Nsc, comb_indexes)

In [ ]:
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxx Plot the true and estimated channels xxxxxxxxxxxxxxxx
# xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
#session = bp.Session(load_from_config=False)
#bp.output_server(docname='simple_precoded_srs_AN1', session=session)

p1 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH11_eq, uH11_eq_est,
    title='Direct Channel from UE1 to AN1')
#bp.show(p1)

p2 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH12_eq, uH12_eq_est,
    title='Interfering Channel from UE2 to AN1')
#bp.show(p2)

p3 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH13_eq, uH13_eq_est,
    title='Interfering Channel from UE3 to AN1')
#bp.show(p3)

tab1 = bw.Panel(child=p1, title="UE1 to AN1")
tab2 = bw.Panel(child=p2, title="UE2 to AN1")
tab3 = bw.Panel(child=p3, title="UE3 to AN1")
tabs_an1 = bw.Tabs(tabs=[tab1, tab2, tab3])
#bp.show(tabs_an1)

In [ ]:
#bp.output_server(docname='simple_precoded_srs_AN2', session=session)
p1 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH21_eq, uH21_eq_est,
    title='Interfering Channel from UE2 to AN2')
#bp.show(p1)

p2 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH22_eq, uH22_eq_est,
    title='Direct Channel from UE2 to AN2')
#bp.show(p2)

p3 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH23_eq, uH23_eq_est,
    title='Interfering Channel from UE3 to AN2')
#bp.show(p3)

tab1 = bw.Panel(child=p1, title="UE1 to AN2")
tab2 = bw.Panel(child=p2, title="UE2 to AN2")
tab3 = bw.Panel(child=p3, title="UE3 to AN2")
tabs_an2 = bw.Tabs(tabs=[tab1, tab2, tab3])
#bp.show(tabs_an2)

In [ ]:
#bp.output_server(docname='simple_precoded_srs_AN3', session=session)
p1 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH31_eq, uH31_eq_est,
    title='Interfering Channel from UE1 to AN3')
p2 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH32_eq, uH32_eq_est,
    title='Interfering Channel from UE2 to AN3')
p3 = plot_true_and_estimated_channel_with_bokeh_all_antennas(
    uH33_eq, uH33_eq_est,
    title='Direct Channel from UE3 to AN3')
#bp.show(p3)

tab1 = bw.Panel(child=p1, title="UE1 to AN3")
tab2 = bw.Panel(child=p2, title="UE2 to AN3")
tab3 = bw.Panel(child=p3, title="UE3 to AN3")
tabs_an3 = bw.Tabs(tabs=[tab1, tab2, tab3])
#bp.show(tabs_an3)

In [ ]:
tabs1 = bw.Panel(child=tabs_an1, title="AN1")
tabs2 = bw.Panel(child=tabs_an2, title="AN2")
tabs3 = bw.Panel(child=tabs_an3, title="AN3")
tabs_all = bw.Tabs(tabs=[tabs1, tabs2, tabs3])
bp.show(tabs_all)

In [ ]: