IUE Pulsamplitudenmodulation.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Copyright 2019 Institut für Nachrichtentechnik, RWTH Aachen University\n",
    "%matplotlib notebook\n",
    "\n",
    "from ipywidgets import interact, interactive, fixed, interact_manual\n",
    "\n",
    "import ipywidgets as widgets\n",
    "from IPython.display import Markdown, Latex\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import signal # convolution, \n",
    "from scipy.io import wavfile # wavfile\n",
    "\n",
    "from ient_nb.ient_plots import *\n",
    "from ient_nb.ient_signals import *\n",
    "from ient_nb.ient_audio import *\n",
    "from ient_nb.ient_transforms import *\n",
    "\n",
    "def convolution(s, h):\n",
    "    # Convolve s and h numerically\n",
    "    g = signal.convolve(s, h, mode='same')*deltat; #g = g[0:len(s)];\n",
    "    return g\n",
    "\n",
    "def ient_ideal_lowpass(s, fg):\n",
    "    return convolution(s, 2*fg*si(2*np.pi*fg*t))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div>\n",
    "    <img src=\"figures/rwth_ient_logo@2x.png\" style=\"float: right;height: 5em;\">\n",
    "</div>\n",
    "\n",
    "# Pulsamplitudenmodulation\n",
    "\n",
    "## Übersicht\n",
    "\n",
    "![Blockdiagramm](figures/pam_block_diagram.png)\n",
    "\n",
    "Hier Rauschprozess $n(t)$ (Kanal) vernachlässigt (ansonsten $g(t)\\rightarrow y(t)=g(t)+n_e(t)$ mit $n_e(t)=n(t)\\ast s(-t)$)\n",
    "\n",
    "## Nutzsignal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fs = 400000 # interne Abtastrate\n",
    "fs0, data = wavfile.read('data/krawehl.wav')\n",
    "f = 0.99*data/np.max(np.abs(data)) # Normalisieren\n",
    "f = signal.resample(f, int(len(f)/fs0*fs)) # Etwas hochtasten für schönere Plots\n",
    "(t, deltat) = np.linspace(-len(f)/fs/2, len(f)/fs/2, len(f), retstep=True) # Zeitachse in Sekunden\n",
    "\n",
    "F, f_ax = ient_dft(f, fs) # Fourier-Transformation\n",
    "\n",
    "# Plot\n",
    "fig,axs=plt.subplots(2,1, figsize=(8,4)); \n",
    "ax = axs[0]\n",
    "ax.plot(1000*t, f, 'rwth');\n",
    "ax.set_xlabel(r'$\\rightarrow t$ [ms]', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow f(t)$'); ient_axis(ax);\n",
    "\n",
    "ax = axs[1]\n",
    "ax.plot(f_ax, np.abs(F), 'rwth');\n",
    "ax.set_xlabel(r'$\\rightarrow f$ [Hz]', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow |F(f)|$'); \n",
    "ax.set_xlim([-10E3, 10E3]); ient_axis(ax);\n",
    "ient_audio_play(f, fs, r'$f(t)$')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Sender\n",
    "\n",
    "Trägersignal $s(t)=\\frac{1}{\\sqrt{t_0}}\\mathrm{rect}\\left(\\frac{t}{t_0}\\right)$ mit Breite $t_0 = 0{,}5 T$ (oder $t_0 = 1{,}25 T$ für Aufgabe 11) und $T=125\\mu\\mathrm{s}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "T = 125E-6\n",
    "t0 = 0.5*T\n",
    "#t0 = 1.25*T # Aufgabe 11: '#' entfernen und alle Cells erneut laufen lassen\n",
    "\n",
    "# Trägersignal\n",
    "s = lambda t,t0: 1/np.sqrt(t0)*rect(t/t0)\n",
    "phi_ss = convolution(s(-t,t0),s(t,t0))\n",
    "\n",
    "# Plot\n",
    "fig, axs = plt.subplots(1,2,figsize=(8,4)); \n",
    "ax = axs[0]; ax.plot(t*1E6, s(t,t0), 'rwth'); \n",
    "ax.set_xlabel(r'$\\rightarrow t$ [$\\mu$s]'); ax.set_ylabel(r'$\\uparrow s(t)$'); \n",
    "ax.set_xlim([-2E6*t0, 2E6*t0]); ient_axis(ax);\n",
    "ient_annotate_distance(ax, r'$t_0$', (-t0/2*1E6,10), (t0/2*1E6,10));\n",
    "\n",
    "ax = axs[1]; ax.plot(t*1E6, phi_ss, 'rwth'); \n",
    "ax.set_xlabel(r'$\\rightarrow t$ [$\\mu$s]', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow s(t)\\ast s(-t)$'); \n",
    "ax.set_xlim([-2E6*t0, 2.5E6*t0]); ient_axis(ax);\n",
    "ient_annotate_distance(ax, r'$2t_0$', (-t0*1E6,.1), (t0*1E6,.1)); fig.tight_layout();\n",
    "\n",
    "Markdown(\"$t_0 = {:.3f}\\ \\mu\\\\mathrm{{s}}$\".format(t0*1E6))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Tiefpassfilterung (ideal, Grenzfrequenz $f_\\mathrm{g} = \\frac{1}{2T}$) und ideale Abtastung"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Tiefpassfilterung\n",
    "fg = 1/(2*T)\n",
    "fTP = ient_ideal_lowpass(f, fg)\n",
    "\n",
    "# Abtastung\n",
    "nT_idx = np.arange(0, len(t), int(fs*T)); nT = t[nT_idx]; \n",
    "fa_nT = fTP[nT_idx];\n",
    "fa = np.zeros_like(t); fa[nT_idx] = fa_nT\n",
    "\n",
    "# Plot\n",
    "ient_audio_play(f, fs, r'$f(t)$')\n",
    "ient_audio_play(fTP, fs, r'$f(t)\\ast h_\\mathrm{LP}(t)$')\n",
    "\n",
    "display(Markdown(\"$f_g = {:.1f}\\ \\\\mathrm{{Hz}}$\".format(fg)))\n",
    "\n",
    "fig,ax=plt.subplots(); \n",
    "ax.plot(1000*t, fTP, 'rwth', label=r'$f(t)$');\n",
    "ient_plot_dirac(ax, 1000*nT, fa_nT, 'rot', label=r'$f_\\mathrm{a}(t)$')\n",
    "ax.set_xlabel(r'$\\rightarrow t$ [ms]', bbox=ient_wbbox); ax.legend();\n",
    "ax.set_xlim([0, 10]); ient_axis(ax);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sendesignal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sender_carrier(s, fa, t0):\n",
    "    m = convolution(s(t, t0), fa/deltat)\n",
    "    return m\n",
    "\n",
    "m = sender_carrier(s, fa, t0)\n",
    "\n",
    "# Plot\n",
    "fig,ax=plt.subplots(); \n",
    "ax.plot(1000*t, m, 'rwth');\n",
    "ax.set_xlabel(r'$\\rightarrow t$ [ms]', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow m(t)$');\n",
    "ax.set_xlim([0, 10]); ient_axis(ax);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Kanal\n",
    "\n",
    "Hier wird $n(t)=0$ angenommen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = np.zeros_like(t)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Empfänger\n",
    "\n",
    "### Korrelationsfilter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Schicke m(t) und n(t) getrennt durch Korrelationsfilter\n",
    "def receiver_filter(m, n, s, t0):\n",
    "    h = s(-t, t0)\n",
    "    g = convolution(h, m)\n",
    "    ne = convolution(h, n)\n",
    "    y = g+ne\n",
    "    return (y,g,ne,h)\n",
    "\n",
    "y,g,ne,h = receiver_filter(m, n, s, t0)\n",
    "\n",
    "# Plot\n",
    "fig,ax = plt.subplots(); \n",
    "ax.plot(1000*t, y, 'rwth');\n",
    "ax.set_xlabel(r'$\\rightarrow t$ [ms]', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow y(t)$');\n",
    "ax.set_xlim([0, 10]);  ient_axis(ax);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Abtastung und Rekonstruktion"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def receiver_sampling_lp(y):\n",
    "    # Abastung\n",
    "    ya_nT = y[nT_idx]\n",
    "    ya = np.zeros_like(t)\n",
    "    ya[nT_idx] = ya_nT\n",
    "\n",
    "    # Rekonstruktion: Tiefpass und Multiplikation mit T\n",
    "    fe = ient_ideal_lowpass(ya/deltat, fg)*T\n",
    "    return (fe, ya, ya_nT)\n",
    "\n",
    "fe, ya, ya_nT = receiver_sampling_lp(y)\n",
    "\n",
    "# Plot\n",
    "fig,ax=plt.subplots(); \n",
    "ient_plot_dirac(ax, 1000*nT, fe[nT_idx], 'black-50')\n",
    "ax.plot(1000*t, fe, 'rwth');\n",
    "ax.set_xlabel(r'$\\rightarrow t$ [ms]', bbox=ient_wbbox);  ax.set_ylabel(r'$\\uparrow f_\\mathrm{e}(t)$');\n",
    "ax.set_xlim([0, 10]); ient_axis(ax);\n",
    "\n",
    "ient_audio_play(fTP, fs, r'Sender: $f(t)\\ast h_\\mathrm{LP}(t)$')\n",
    "ient_audio_play(fe, fs, r'Empfänger: $f_\\mathrm{e}(t)$')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Kompensation am Empfänger\n",
    "\n",
    "In diesem Abschnitt geht es um den in Aufgabe 11 beschriebenen Fall, in dem $t_0 > T$ gilt."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def receiver_compensation_filter(t0, T):\n",
    "    if t0 > T:\n",
    "        # Kompensationsfilter\n",
    "        H2 = 1/(1+0.4*np.cos(2*np.pi*f_ax*T))*rect(f_ax/(2*fg))\n",
    "        print(\"Kompensiere\")\n",
    "    else:\n",
    "        H2 = np.ones_like(F)\n",
    "        print(\"Kompensation nicht nötig\")\n",
    "    return H2\n",
    "\n",
    "H2 = receiver_compensation_filter(t0, T)\n",
    "\n",
    "# Plot\n",
    "fig,ax=plt.subplots(); ax.plot(f_ax, np.abs(H2), 'rwth'); \n",
    "ax.set_xlabel(r'$\\rightarrow f$ [Hz]', bbox=ient_wbbox);  ax.set_ylabel(r'$\\uparrow |H_2(f)|$');\n",
    "ax.set_xlim([-2*fg, 2*fg]); ient_axis(ax);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Filterung und Vergleich von $F_\\mathrm{e}(f)$ mit $F(f)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def receiver_compensate(fe, H2):\n",
    "    # Fourier-Transformation\n",
    "    Fe,_ = ient_dft(fe, fs)\n",
    "\n",
    "    # Kompensation (im Frequenzbereich)\n",
    "    Fe2 = Fe*H2\n",
    "    \n",
    "    # Inverse Fourier-Transformation\n",
    "    fe2 = np.real(ient_idft(Fe2, len(t)))\n",
    "    \n",
    "    return (fe2, Fe2, Fe)\n",
    "\n",
    "fe2, Fe2, Fe = receiver_compensate(fe, H2)\n",
    "\n",
    "# Plot\n",
    "fig, axs = plt.subplots(2,2, figsize=(8,4)); \n",
    "ax = axs[0,0]; \n",
    "ax.plot(f_ax, np.abs(F),  'rwth', label=r'$|F(f)|$'); \n",
    "ax.plot(f_ax, np.abs(Fe), 'grun', label=r'$|F_\\mathrm{e}(f)|$'); \n",
    "ax.set_xlabel(r'$\\rightarrow f$ [Hz]', bbox=ient_wbbox); \n",
    "ax.set_xlim([-fg, fg]); ax.legend(); ient_axis(ax);\n",
    "\n",
    "ax = axs[0,1]; \n",
    "ax.plot(f_ax, np.abs(F),  'rwth'); \n",
    "ax.plot(f_ax, np.abs(Fe), 'grun'); \n",
    "ax.set_xlabel(r'$\\rightarrow f$ [Hz]', bbox=ient_wbbox); \n",
    "ax.set_xlim([0, fg/10]); ient_axis(ax);\n",
    "\n",
    "ax = axs[1,0]; \n",
    "ax.plot(f_ax, np.abs(F), 'rwth', label=r'$|F(f)|$'); \n",
    "ax.plot(f_ax, np.abs(Fe2), 'grun', label=r'$|F_\\mathrm{e}(f) \\cdot H_2(f)|$'); \n",
    "ax.set_xlabel(r'$\\rightarrow f$ [Hz]', bbox=ient_wbbox); \n",
    "ax.set_xlim(axs[0,0].get_xlim()); ax.set_ylim(axs[0,0].get_ylim()); ax.legend(); ient_axis(ax);\n",
    "\n",
    "ax = axs[1,1]; \n",
    "ax.plot(f_ax, np.abs(F), 'rwth'); \n",
    "ax.plot(f_ax, np.abs(Fe2), 'grun'); \n",
    "ax.set_xlabel(r'$\\rightarrow f$ [Hz]', bbox=ient_wbbox); \n",
    "ax.set_xlim(axs[0,1].get_xlim()); ax.set_ylim(axs[0,1].get_ylim()); ient_axis(ax);\n",
    "\n",
    "ient_audio_play(fTP, fs, r'$f(t)\\ast h_\\mathrm{LP}(t)$')\n",
    "ient_audio_play(fe2, fs, r'$f_\\mathrm{e}(t)\\ast h_2(t)$')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Interaktive Demo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fix,axs = plt.subplots(2,2, figsize=(8,8)); \n",
    "@interact(t0byT = widgets.FloatSlider(min=0.5,max=1.25,step=0.25,value=0.5,description=r'$t_0/T$:', continuous_update=False))\n",
    "def update_plot(t0byT):\n",
    "    t0 = t0byT*T\n",
    "    # Sender: Träger mit Breite t0\n",
    "    m = sender_carrier(s, fa, t0)\n",
    "    # Empfänger: Korrelationsfilter, Abtastung und Tiefpass\n",
    "    y,g,ne,h = receiver_filter(m, n, s, t0)\n",
    "    fe, ya, ya_nT = receiver_sampling_lp(y)\n",
    "    # Empfänger: Kompensation\n",
    "    H2 = receiver_compensation_filter(t0, T)\n",
    "    fe2, Fe2, Fe = receiver_compensate(fe, H2)\n",
    "    \n",
    "    # Plot\n",
    "    if not axs[0,0].lines:\n",
    "        ax = axs[0,0]; ax.plot(t/T, m, 'rwth');\n",
    "        ax.set_xlabel(r'$\\rightarrow t/T$', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow m(t)$', bbox=ient_wbbox);\n",
    "        ax.set_xlim([0,50]); ax.grid(); ient_axis(ax)\n",
    "        ax = axs[0,1]; ax.plot(t/T, y, 'rwth');\n",
    "        ax.set_xlabel(r'$\\rightarrow t/T$', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow y(t)$', bbox=ient_wbbox);\n",
    "        ax.set_xlim(axs[0,0].get_xlim()); ax.grid(); ient_axis(ax)\n",
    "        ax = axs[1,0]; ax.plot(t/T, fe, 'rwth');\n",
    "        ax.set_xlabel(r'$\\rightarrow t/T$', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow f_\\mathrm{e}(t)$', bbox=ient_wbbox);\n",
    "        ax.set_xlim(axs[0,0].get_xlim()); ax.grid(); ient_axis(ax)\n",
    "        ax = axs[1,1]; ax.plot(t/T, fe, 'rwth');\n",
    "        ax.set_xlabel(r'$\\rightarrow t/T$', bbox=ient_wbbox); ax.set_ylabel(r'$\\uparrow f_{\\mathrm{e}}(t)\\ast h_2(t)$', bbox=ient_wbbox);\n",
    "        ax.set_xlim(axs[0,0].get_xlim()); ax.grid(); ient_axis(ax)\n",
    "    else:\n",
    "        axs[0,0].lines[0].set_ydata(m);  axs[0,1].lines[0].set_ydata(y);\n",
    "        axs[1,0].lines[0].set_ydata(fe); axs[1,1].lines[0].set_ydata(fe2);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This notebook is provided as [Open Educational Resource](https://en.wikipedia.org/wiki/Open_educational_resources) (OER). Feel free to use the notebook for your own purposes. The code is licensed under the [MIT license](https://opensource.org/licenses/MIT). \n",
    "\n",
    "Please attribute the work as follows: \n",
    "*Christian Rohlfing, Übungsbeispiele zur Vorlesung \"Informationsübertragung\"*, gehalten von Jens-Rainer Ohm, 2019, Institut für Nachrichtentechnik, RWTH Aachen University."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}