{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Recover time information from multiple antennae"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "rng = np.random.default_rng(12345)\n",
    "from scipy import signal as sgl\n",
    "\n",
    "import os\n",
    "import sys\n",
    "# Append parent directory to import path so lib can be found\n",
    "sys.path.append(os.path.dirname(os.path.abspath(os.getcwd())))\n",
    "\n",
    "from lib import TravelSignal\n",
    "from lib.sampling import Digitizer\n",
    "from lib.location import Receiver, Emitter\n",
    "\n",
    "import lib.sampling as smp\n",
    "import lib.location as loc\n",
    "import lib.util as util\n",
    "\n",
    "km = 1e3\n",
    "ns = 1e9"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Setup emitters and Antenna\n",
    "############################\n",
    "if True:\n",
    "    emitter = Emitter([1,1])*km\n",
    "else:\n",
    "    emitter = Emitter([0,0])*km\n",
    "\n",
    "# reference antenna\n",
    "ref_ant = Receiver(np.array([2,2])*km)\n",
    "\n",
    "# antenna circle\n",
    "antennae = [\n",
    "    Receiver( 5*km * util.rot_vector(np.pi*2/6* n)) for n in range(6)\n",
    "]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 14880 characters]\n",
      "text/plain": [
       "<Figure size 288x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "## Show Geometry\n",
    "#################\n",
    "if True:\n",
    "    fig, ax = plt.subplots(1,1, figsize=(4,4))\n",
    "    loc.plot_geometry(ax, [emitter], antennae)\n",
    "\n",
    "    if True:\n",
    "        ax.plot(*ref_ant.x, \"x\")\n",
    "        ax.annotate(\"$\\mathrm{A}_\\mathrm{ref}$\", ref_ant.x);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 17812 characters]\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Setup signal\n",
    "########################\n",
    "\n",
    "fc = 50e6# MHz\n",
    "signal_t_length = 200e-9\n",
    "signal_pulse_center = 125e-9 # 100 ns\n",
    "signal_sample_rate = 5e9 # Hz\n",
    "signal_time = util.sampled_time(\n",
    "    sample_rate=signal_sample_rate,\n",
    "    start=0,\n",
    "    end=signal_t_length\n",
    ")\n",
    "\n",
    "analog_signal = TravelSignal(\n",
    "    sgl.gausspulse(signal_time - signal_pulse_center, fc=fc),\n",
    "    signal_sample_rate,\n",
    "    t_0 = 0,\n",
    "    periodic=False\n",
    ")\n",
    "\n",
    "del fc, signal_t_length, signal_pulse_center, signal_sample_rate\n",
    "\n",
    "if True:\n",
    "    fig, ax = plt.subplots()\n",
    "    ax.plot(signal_time, analog_signal(signal_time))\n",
    "    plt.show()\n",
    "    \n",
    "    \n",
    "## the emitted signal\n",
    "emitted = emitter.emit(analog_signal)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Spatial time offset emitter-antenna: 4717.31 ns\n"
     ]
    },
    {
     "data": {
      "image/png": "[encoded data, 22752 characters]\n",
      "text/plain": [
       "<Figure size 1152x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Show how the signal is emitted\n",
    "# and received by a reference antenna\n",
    "########################\n",
    "\n",
    "# analog time corresponds to a very high sampled time range\n",
    "#  (to fully show the signal)\n",
    "# that is sufficiently long to show both the emitted and received signal\n",
    "if not True:\n",
    "    t_start = 0\n",
    "else:\n",
    "    t_start = -50e-8\n",
    "\n",
    "if not True:\n",
    "    t_end = 5e-7 # ms\n",
    "else:\n",
    "    t_end = 50e-7 # ms\n",
    "t_offset = 0 # 100 ns\n",
    "sample_rate = 1e10 # 1/s\n",
    "\n",
    "analog_time = util.sampled_time(sample_rate=sample_rate, start=t_start, end=t_end, offset=t_offset)\n",
    "\n",
    "\n",
    "if True:\n",
    "    # Figure\n",
    "    fig, ax = plt.subplots(1, 1, figsize=(16,4))\n",
    "    ax.set_title(\"Signal as emitted from the emitter, and received by the reference antenna\")\n",
    "    ax.set_ylabel(\"Amplitude\")\n",
    "    ax.set_xlabel(\"Time (ns)\")\n",
    "\n",
    "    ## Emitter\n",
    "    received = emitter.recv(emitted)\n",
    "\n",
    "    ax.plot(analog_time * ns,   received(analog_time), label=\"Emitter\")\n",
    "\n",
    "\n",
    "    ## Reference antenna\n",
    "    received = ref_ant.recv(emitted)\n",
    "\n",
    "    ax.plot(analog_time * ns,   received(analog_time), label=\"Reference Antenna\")\n",
    "    print(\"Spatial time offset emitter-antenna: {:.6g} ns\".format(analog_signal.spatial_time_offset(emitter.x, ref_ant.x)*ns))\n",
    "\n",
    "    ax.legend(loc='upper center')\n",
    "    plt.show();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 21092 characters]\n",
      "text/plain": [
       "<Figure size 864x432 with 8 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Show absolute timing over all antennae\n",
    "########################################\n",
    "\n",
    "if True:\n",
    "    longtime = util.sampled_time(\n",
    "        sample_rate=sample_rate,\n",
    "        start=t_start,\n",
    "        end=5*t_end,\n",
    "        offset=t_offset\n",
    "    )\n",
    "\n",
    "    ylabel_kw = {\"rotation\": \"horizontal\", \"va\":\"center\", \"ha\":\"center\", \"labelpad\": 30}\n",
    "    \n",
    "    fig, axs = plt.subplots(2+len(antennae), 1, sharex=True, figsize=(12,6), gridspec_kw={\"hspace\":0})\n",
    "    axs[0].set_title(\"Traces of Emitter and Antennae\")\n",
    "\n",
    "    # Emitter\n",
    "    i = 0\n",
    "    axs[i].set_ylabel(\"Emitter\".format(), **ylabel_kw)\n",
    "    axs[i].set_xlabel(\"Time (ns)\")\n",
    "    axs[i].plot(longtime * ns, emitted(longtime), label=\"Emitter\")\n",
    "    axs[i].plot(signal_time * ns, emitted(signal_time), label=\"Signal time\", zorder=2.5)\n",
    "\n",
    "    # Reference Antenna\n",
    "    i +=1\n",
    "    axs[i].set_ylabel(\"Ref Antenna\", **ylabel_kw)\n",
    "    axs[i].plot(longtime * ns, ref_ant.recv(emitted)(longtime), label=\"Ref Antenna\")\n",
    "    axs[i].set_xlabel(\"Time (ns)\")\n",
    "\n",
    "    # Antenna\n",
    "    for j, ant in enumerate(antennae):\n",
    "        i +=1\n",
    "        axs[i].set_ylabel(\"Antenna {}\".format(j), **ylabel_kw)\n",
    "        axs[i].plot(longtime * ns, ant.recv(emitted)(longtime), label=\"Antenna {}\".format(j))\n",
    "        axs[i].set_xlabel(\"Time (ns)\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Time delay ref ant -- emitter (ns): 9434.617346998737\n",
      "T_threshold Emitter (ns): 123.0\n",
      "T_threshold Antenna (ns): 9680.617346998737\n",
      "Diff thresholds (ns): 9557.617346998737\n"
     ]
    },
    {
     "data": {
      "image/png": "[encoded data, 15472 characters]\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Determine time from signal thresholds for emitter and reference receiver\n",
    "######################################\n",
    "\n",
    "threshold = 0.8\n",
    "\n",
    "if True:\n",
    "    mysignal = analog_signal\n",
    "else:\n",
    "    mysignal = emitted\n",
    "    \n",
    "received = ref_ant.recv(mysignal)\n",
    "#received.periodic=True\n",
    "\n",
    "## Recover signal time from signal\n",
    "signal_time = util.sampled_time(sample_rate=mysignal.sample_rate, start=0, end=mysignal.time_length)\n",
    "\n",
    "## Find time offset and recover approximate signal time at antenna\n",
    "# Note: no idea where 2 comes from right now\n",
    "time_offset = -2 * analog_signal.total_time_offset(x_f = ref_ant.x, x_0 = emitter.x, velocity = None)\n",
    "ref_ant_signal_time = 2*signal_time + time_offset\n",
    "\n",
    "\n",
    "# Find t_threshold for emitter\n",
    "raw_emit = mysignal(signal_time)\n",
    "t_thres_emit_idx = util.detect_edges(threshold, raw_emit, rising=True, falling=False)\n",
    "\n",
    "if t_thres_emit_idx:\n",
    "    t_thres_emit = signal_time[t_thres_emit_idx[0]]\n",
    "else:\n",
    "    t_thres_emit = None\n",
    "del t_thres_emit_idx\n",
    "    \n",
    "    \n",
    "# Find t_threshold for receiver\n",
    "raw_ant = received(ref_ant_signal_time)\n",
    "t_thres_ant_idx = util.detect_edges(threshold, raw_emit, rising=True, falling=False)\n",
    "\n",
    "if t_thres_ant_idx:\n",
    "    t_thres_ant = ref_ant_signal_time[t_thres_ant_idx[0]]\n",
    "else:\n",
    "    t_thres_ant = None\n",
    "del t_thres_ant_idx\n",
    "\n",
    "print(\"Time delay ref ant -- emitter (ns):\", time_offset*ns)\n",
    "print(\"T_threshold Emitter (ns):\", t_thres_emit*ns)\n",
    "print(\"T_threshold Antenna (ns):\", t_thres_ant*ns)\n",
    "print(\"Diff thresholds (ns):\", (t_thres_ant - t_thres_emit)*ns)\n",
    "\n",
    "## Make another plot\n",
    "fig, ax = plt.subplots()\n",
    "axs[i].set_xlabel(\"Time (ns)\")\n",
    "# Emitter\n",
    "if not True:\n",
    "    ax.plot(signal_time *ns,  raw_emit, label='Emitter')\n",
    "    ax.plot(t_thres_emit*ns, mysignal(t_thres_emit), '*')\n",
    "    ax.annotate(\"t_1\", xy=(t_thres_emit*ns, mysignal(t_thres_emit)))\n",
    "\n",
    "# Antenna\n",
    "if True:\n",
    "    ax.plot(ref_ant_signal_time*ns,   raw_ant, label=\"Reference Antenna\")\n",
    "    ax.plot(t_thres_ant*ns, received(t_thres_ant), '*')\n",
    "    ax.annotate(\"t_2\", xy=(t_thres_ant*ns, received(t_thres_ant)))\n",
    "\n",
    "plt.show();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.10.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
