Playing arround with plotting OFDM carriers with matplotlib

nbiot_ofdm_visualisation.py

from matplotlib import pyplot as plt
import numpy as np
from scipy.optimize import fsolve
from matplotlib.ticker import EngFormatter

def get_lobe(carrier, df, x_init, ncarriers=12,):
  "Calculating the side lobe signal strength."
  # Sinc function is by default scaled with pu cos is not
  carrier = carrier - ncarriers/2
  func = lambda z : np.sinc(z) - np.cos(z*np.pi)
  x_init = x_init / (carrier * df)
  x_solution = fsolve(func, x_init)
  x_scaled = (carrier + x_solution) * df
  y = np.sinc(x_solution)
  y_db = 10*np.log10(y**2)
  return [x_scaled, y, y_db]

w = float(1)/15e3
carriers = 12
x = np.linspace(-15e3*(carriers+2),15e3*(carriers+2),1000)
y = np.zeros((carriers, len(x)))

tt = np.zeros(len(x))

# Plot and subplots
f, (ax1, ax2, ax3) = plt.subplots(3, 1)

# Sexy formatter engineering style
formatter0 = EngFormatter(unit='Hz')

# Calculate and plot all different carriers
for i in xrange(-carriers/2, carriers/2+1):
  if i != 0:
    y[i] = np.sinc(x*w - i)
    tt = tt + y[i]
    ax1.plot(x, y[i], '-', linewidth=1, label="Carrier " + str(i + carriers/2))

# Plot the sum signal in an extra subplot
tt_max = np.max(tt)
tt_db = 10*np.log10((tt/tt_max)**2)
ax2.plot(x, tt_db, '-', linewidth=1, label="Sum signal")
ax2.set(ylim=(-60, 3))

# Reference plot for finding maximas
b = np.sinc(x*w)
c = np.cos(x*np.pi*w)
ax3.plot(x,b,label="Sinc")
ax3.plot(x, c, label="Cosine")

# Sidelobe
lx, ly, ly_db = get_lobe(12, 15e3, 112e3)
# Sum signal has another value for the side lobe
ax = min(x.tolist(), key=lambda x:abs(x-lx))
ax = x.tolist().index(ax)
ay = 10*np.log10( (tt[ax]/tt_max)**2 )

ax2.annotate("First side lobe, " + str(round(ay,1)) + "dB",
            xy=(lx, ay), xycoords='data',
            xytext=(-15, 25), textcoords='offset points',
            arrowprops=dict(facecolor='black', shrink=0.05),
            horizontalalignment='right', verticalalignment='bottom')

# Plot formatting
for ax in (ax1, ax2, ax3):
  ax.xaxis.set_major_formatter(formatter0)
  ax.set_ylabel("H(f)")
  ax.set_xlabel("f in kHz")
  ax.legend()
  ax.grid(True)

ax3.set_title("Maximas of the sinc fucntion")
ax2.set_title("Sum signal")
ax1.set_title("Different subcarriers")

plt.show()