2022-07-11 11:29:38 +02:00
|
|
|
{
|
|
|
|
"cells": [
|
|
|
|
{
|
|
|
|
"cell_type": "markdown",
|
|
|
|
"metadata": {},
|
|
|
|
"source": [
|
|
|
|
"# Beacon Sync\n",
|
|
|
|
"\n",
|
|
|
|
"Synchronise two delta peaks, by using an intermediate beacon that was sent out together with it."
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
|
|
|
"execution_count": 1,
|
|
|
|
"metadata": {},
|
|
|
|
"outputs": [],
|
|
|
|
"source": [
|
|
|
|
"import numpy as np\n",
|
|
|
|
"import matplotlib.pyplot as plt\n",
|
|
|
|
"import scipy.signal as signal\n",
|
|
|
|
"import scipy.fft as ft"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
|
|
|
"execution_count": 2,
|
|
|
|
"metadata": {},
|
|
|
|
"outputs": [],
|
|
|
|
"source": [
|
|
|
|
"# Monkey patch correlation_lags\n",
|
|
|
|
"if not hasattr(signal, 'correlation_lags'):\n",
|
|
|
|
" def correlation_lags(in1_len, in2_len, mode='full'):\n",
|
|
|
|
" r\"\"\"\n",
|
|
|
|
" Calculates the lag / displacement indices array for 1D cross-correlation.\n",
|
|
|
|
" Parameters\n",
|
|
|
|
" ----------\n",
|
|
|
|
" in1_size : int\n",
|
|
|
|
" First input size.\n",
|
|
|
|
" in2_size : int\n",
|
|
|
|
" Second input size.\n",
|
|
|
|
" mode : str {'full', 'valid', 'same'}, optional\n",
|
|
|
|
" A string indicating the size of the output.\n",
|
|
|
|
" See the documentation `correlate` for more information.\n",
|
|
|
|
" See Also\n",
|
|
|
|
" --------\n",
|
|
|
|
" correlate : Compute the N-dimensional cross-correlation.\n",
|
|
|
|
" Returns\n",
|
|
|
|
" -------\n",
|
|
|
|
" lags : array\n",
|
|
|
|
" Returns an array containing cross-correlation lag/displacement indices.\n",
|
|
|
|
" Indices can be indexed with the np.argmax of the correlation to return\n",
|
|
|
|
" the lag/displacement.\n",
|
|
|
|
" Notes\n",
|
|
|
|
" -----\n",
|
|
|
|
" Cross-correlation for continuous functions :math:`f` and :math:`g` is\n",
|
|
|
|
" defined as:\n",
|
|
|
|
" .. math::\n",
|
|
|
|
" \\left ( f\\star g \\right )\\left ( \\tau \\right )\n",
|
|
|
|
" \\triangleq \\int_{t_0}^{t_0 +T}\n",
|
|
|
|
" \\overline{f\\left ( t \\right )}g\\left ( t+\\tau \\right )dt\n",
|
|
|
|
" Where :math:`\\tau` is defined as the displacement, also known as the lag.\n",
|
|
|
|
" Cross correlation for discrete functions :math:`f` and :math:`g` is\n",
|
|
|
|
" defined as:\n",
|
|
|
|
" .. math::\n",
|
|
|
|
" \\left ( f\\star g \\right )\\left [ n \\right ]\n",
|
|
|
|
" \\triangleq \\sum_{-\\infty}^{\\infty}\n",
|
|
|
|
" \\overline{f\\left [ m \\right ]}g\\left [ m+n \\right ]\n",
|
|
|
|
" Where :math:`n` is the lag.\n",
|
|
|
|
" Examples\n",
|
|
|
|
" --------\n",
|
|
|
|
" Cross-correlation of a signal with its time-delayed self.\n",
|
|
|
|
" >>> from scipy import signal\n",
|
|
|
|
" >>> from numpy.random import default_rng\n",
|
|
|
|
" >>> rng = default_rng()\n",
|
|
|
|
" >>> x = rng.standard_normal(1000)\n",
|
|
|
|
" >>> y = np.concatenate([rng.standard_normal(100), x])\n",
|
|
|
|
" >>> correlation = signal.correlate(x, y, mode=\"full\")\n",
|
|
|
|
" >>> lags = signal.correlation_lags(x.size, y.size, mode=\"full\")\n",
|
|
|
|
" >>> lag = lags[np.argmax(correlation)]\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
"\n",
|
|
|
|
" # calculate lag ranges in different modes of operation\n",
|
|
|
|
" if mode == \"full\":\n",
|
|
|
|
" # the output is the full discrete linear convolution\n",
|
|
|
|
" # of the inputs. (Default)\n",
|
|
|
|
" lags = np.arange(-in2_len + 1, in1_len)\n",
|
|
|
|
" elif mode == \"same\":\n",
|
|
|
|
" # the output is the same size as `in1`, centered\n",
|
|
|
|
" # with respect to the 'full' output.\n",
|
|
|
|
" # calculate the full output\n",
|
|
|
|
" lags = np.arange(-in2_len + 1, in1_len)\n",
|
|
|
|
" # determine the midpoint in the full output\n",
|
|
|
|
" mid = lags.size // 2\n",
|
|
|
|
" # determine lag_bound to be used with respect\n",
|
|
|
|
" # to the midpoint\n",
|
|
|
|
" lag_bound = in1_len // 2\n",
|
|
|
|
" # calculate lag ranges for even and odd scenarios\n",
|
|
|
|
" if in1_len % 2 == 0:\n",
|
|
|
|
" lags = lags[(mid-lag_bound):(mid+lag_bound)]\n",
|
|
|
|
" else:\n",
|
|
|
|
" lags = lags[(mid-lag_bound):(mid+lag_bound)+1]\n",
|
|
|
|
" elif mode == \"valid\":\n",
|
|
|
|
" # the output consists only of those elements that do not\n",
|
|
|
|
" # rely on the zero-padding. In 'valid' mode, either `in1` or `in2`\n",
|
|
|
|
" # must be at least as large as the other in every dimension.\n",
|
|
|
|
"\n",
|
|
|
|
" # the lag_bound will be either negative or positive\n",
|
|
|
|
" # this let's us infer how to present the lag range\n",
|
|
|
|
" lag_bound = in1_len - in2_len\n",
|
|
|
|
" if lag_bound >= 0:\n",
|
|
|
|
" lags = np.arange(lag_bound + 1)\n",
|
|
|
|
" else:\n",
|
|
|
|
" lags = np.arange(lag_bound, 1)\n",
|
|
|
|
" return lags\n",
|
|
|
|
"\n",
|
|
|
|
" signal.correlation_lags = correlation_lags"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
|
|
|
"execution_count": 3,
|
|
|
|
"metadata": {},
|
|
|
|
"outputs": [],
|
|
|
|
"source": [
|
|
|
|
"### signal generation\n",
|
|
|
|
"def fft_bandpass(signal, band, samplerate):\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" Simple bandpassing function employing a FFT.\n",
|
|
|
|
"\n",
|
|
|
|
" Parameters\n",
|
|
|
|
" ----------\n",
|
|
|
|
" signal : arraylike\n",
|
|
|
|
" band : tuple(low, high)\n",
|
|
|
|
" Frequencies for bandpassing\n",
|
|
|
|
" samplerate : float\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" signal = np.asarray(signal)\n",
|
|
|
|
"\n",
|
|
|
|
" fft = ft.rfft(signal)\n",
|
|
|
|
" freqs = ft.rfftfreq(signal.size, 1/samplerate)\n",
|
|
|
|
" fft[(freqs < band[0]) | (freqs > band[1])] = 0\n",
|
|
|
|
" \n",
|
|
|
|
" return ft.irfft(fft, signal.size), (fft, freqs)\n",
|
|
|
|
"\n",
|
|
|
|
"def deltapeak(timelength=1e3, samplerate=1, offset=None, peaklength=1):\n",
|
|
|
|
" N_samples = int(timelength * samplerate)\n",
|
|
|
|
" if offset is None:\n",
|
|
|
|
" offset = (np.random.random(1)*N_samples).astype(int) % N_samples\n",
|
|
|
|
" elif isinstance(offset, (tuple, list)):\n",
|
|
|
|
" offset_min = offset[0]\n",
|
|
|
|
" offset_max = offset[-1]\n",
|
|
|
|
" \n",
|
|
|
|
" offset = (np.random.random(1)*(offset_max - offset_min)+offset_min).astype(int) % N_samples\n",
|
|
|
|
" \n",
|
|
|
|
" position = (offset + np.arange(0, peaklength)).astype(int) % N_samples\n",
|
|
|
|
" \n",
|
|
|
|
" signal = np.zeros(N_samples)\n",
|
|
|
|
" signal[position] = 1\n",
|
|
|
|
" \n",
|
|
|
|
" return signal, position\n",
|
|
|
|
"\n",
|
|
|
|
"def sin_delay(f, t, t_delay=0, phase=0):\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" return np.sin(2*np.pi*f * (t + t_delay) + phase)\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"\n",
|
|
|
|
"def annotate_width(ax, name, x1, x2, y, text_kw={}, arrow_kw={}):\n",
|
|
|
|
" default_arrow_kw = dict(\n",
|
|
|
|
" xy = (x1, y),\n",
|
|
|
|
" xytext = (x2,y),\n",
|
|
|
|
" arrowprops = dict(\n",
|
|
|
|
" arrowstyle=\"<->\",\n",
|
|
|
|
" shrinkA=False,\n",
|
|
|
|
" shrinkB=False\n",
|
|
|
|
" ),\n",
|
|
|
|
" )\n",
|
|
|
|
"\n",
|
|
|
|
" default_text_kw = dict(\n",
|
|
|
|
" va='bottom',\n",
|
|
|
|
" ha='center',\n",
|
|
|
|
" xy=((x1+x2)/2, y)\n",
|
|
|
|
" )\n",
|
|
|
|
"\n",
|
|
|
|
" an1 = ax.annotate(\"\", **{**default_arrow_kw, **arrow_kw})\n",
|
|
|
|
" an2 = ax.annotate(name, **{**default_text_kw, **text_kw})\n",
|
|
|
|
"\n",
|
|
|
|
" return [an1, an2]\n",
|
|
|
|
"\n",
|
|
|
|
"def time2phase(time, frequency=1):\n",
|
|
|
|
" return 2*np.pi*frequency*time\n",
|
|
|
|
"\n",
|
|
|
|
"def phase2time(phase, frequency=1):\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" return phase/(2*np.pi*frequency)\n",
|
|
|
|
"\n",
|
|
|
|
"def time_roll(a, samplerate, time_shift, *roll_args, **roll_kwargs):\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" Like np.roll, but use samplerate and time_shift to approximate\n",
|
|
|
|
" the offset to roll.\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" shift = np.rint(time_shift*samplerate).astype(int)\n",
|
|
|
|
" return np.roll(a, shift, *roll_args, **roll_kwargs)"
|
2022-07-11 11:29:38 +02:00
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 4,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"name": "stdout",
|
|
|
|
"output_type": "stream",
|
|
|
|
"text": [
|
2022-08-01 15:17:24 +02:00
|
|
|
"Beacon period [ns]: 14.285714285714285\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"Beacon initial [ns]: 4.4\n",
|
|
|
|
"Beacon initial [phase]: 1.9352210746113125\n",
|
|
|
|
"Beacon initial [idx]: 22.0\n",
|
|
|
|
"Beacon difference [ns]: 8.571428571428571\n",
|
|
|
|
"Beacon difference [phase]: 3.7699111843077517\n",
|
|
|
|
"Impulse offsets [ns]: [ 75.4 149.6]\n",
|
|
|
|
"Time difference Impulses [ns]: 74.20000000000002\n"
|
2022-07-11 11:29:38 +02:00
|
|
|
]
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"us = 1e3 # ns\n",
|
|
|
|
"ns = 1/us # us\n",
|
|
|
|
"\n",
|
|
|
|
"\n",
|
|
|
|
"band = (30, 80) # MHz\n",
|
|
|
|
"samplerate = 5000 # MHz\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"timelength = 0.2 # us\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"\n",
|
|
|
|
"time = np.arange(0, timelength, 1/samplerate)\n",
|
|
|
|
"\n",
|
|
|
|
"# generate beacons\n",
|
|
|
|
"if True: # in-band\n",
|
|
|
|
" f_beacon = 70 # MHz\n",
|
|
|
|
"else: # under band\n",
|
|
|
|
" f_beacon = 20 # MHz\n",
|
|
|
|
"\n",
|
|
|
|
"beacon_amplitude = 0.1\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"beacon_init_phase = time2phase(4.4*ns, f_beacon)\n",
|
|
|
|
"beacon_phase_offset = 1.2*np.pi\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"\n",
|
|
|
|
"beacons = np.array([\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" beacon_amplitude * sin_delay(f_beacon, time, t_delay=0, phase=-beacon_init_phase),\n",
|
|
|
|
" beacon_amplitude * sin_delay(f_beacon, time, t_delay=0, phase=-beacon_init_phase-beacon_phase_offset)\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"])\n",
|
|
|
|
"\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"# generate impulses\n",
|
|
|
|
"impulses = []\n",
|
|
|
|
"impulses_offsets = []\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"impulses_def_offsets = [\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" (0.3*len(time),0.4*len(time)),\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
" (0.7*len(time),0.8*len(time)),\n",
|
|
|
|
" ]# random offsets in interval\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"for i in range(2):\n",
|
|
|
|
" offset = None\n",
|
|
|
|
" if impulses_def_offsets:\n",
|
|
|
|
" if len(impulses_def_offsets) == 1:\n",
|
|
|
|
" offset = impulses_def_offsets[0]\n",
|
|
|
|
" else:\n",
|
|
|
|
" offset = impulses_def_offsets[i]\n",
|
|
|
|
" orig_imp, imp_offset = deltapeak(timelength, samplerate, offset=offset, peaklength=1)\n",
|
|
|
|
"\n",
|
|
|
|
" ## Bandpass it\n",
|
|
|
|
" imp, _ = fft_bandpass(orig_imp, band, samplerate)\n",
|
|
|
|
" imp /= np.max(imp)\n",
|
|
|
|
" \n",
|
|
|
|
" impulses.append(imp)\n",
|
|
|
|
" impulses_offsets.append(imp_offset/samplerate)\n",
|
|
|
|
"\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"impulses = np.array(impulses)\n",
|
|
|
|
"impulses_offsets = np.array(impulses_offsets)\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"print(\"Beacon period [ns]:\", 1/f_beacon/ns)\n",
|
|
|
|
"print(\"Beacon initial [ns]:\", phase2time(beacon_init_phase, f_beacon) /ns)\n",
|
|
|
|
"print(\"Beacon initial [phase]:\", beacon_init_phase)\n",
|
|
|
|
"print(\"Beacon initial [idx]:\", phase2time(beacon_init_phase, f_beacon)*samplerate)\n",
|
|
|
|
"print(\"Beacon difference [ns]:\", phase2time(beacon_phase_offset, f_beacon)/ns)\n",
|
|
|
|
"print(\"Beacon difference [phase]:\", beacon_phase_offset)\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"print(\"Impulse offsets [ns]:\", impulses_offsets[:,0]/ns)\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"print(\"Time difference Impulses [ns]: {}\".format( (impulses_offsets[1,0]-impulses_offsets[0,0])/ns ))"
|
2022-07-11 11:29:38 +02:00
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 15:17:24 +02:00
|
|
|
"execution_count": 5,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [],
|
|
|
|
"source": [
|
|
|
|
"full_signals = impulses + beacons"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 6,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"data": {
|
2022-08-01 17:03:26 +02:00
|
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAroAAAEHCAYAAABBQeUvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOydd3xb1dmAnyt57z0SJ3Ycx3H2MtmLDEYYoS0U2kIKZbQQCoUOKFCgDVAos6wCAcoofFBGGAkNBMgki4QM24lXllcc772l8/2hSJFtyb6yJUvXOc/vZ3Cke68evz4+573nnqEIIZBIJBKJRCKRSAYbOncLSCQSiUQikUgkrkAmuhKJRCKRSCSSQYlMdCUSiUQikUgkgxKZ6EokEolEIpFIBiUy0ZVIJBKJRCKRDEpkoiuRSCQSiUQiGZR4ueKiUVFRIikpyRWXdjoNDYcBCAoa42YTdbS0tADg5+fnZpPe0ZIraMtXK65NTTkA6HSJgOf7gnZia0ZLvlpyBW35askVtOWrJVfQnq8z2Lt3b4UQItrWey5JdJOSktizZ48rLu10tm+fAcDs2bvcbKKO3NxcAFJTU91s0jtacgVt+WrFdd++hQAEBr4CeL4vaCe2ZrTkqyVX0JavllxBW75acgXt+ToDRVFO2HtPDl2QSCTceuutJCYmOvWa7QYjVY1tnKxtJrO41qnXlkgkEolEDTLRlUjOco4dO8amTZtoa2ujvr7eKdcsq29h+fPfkVdWT0FVE7/9v308vj6bDoPRKdeXSCQSiUQNiiu2AE5PTxdaGbogkZztrFixgmXLlrF69WoefvhhZs6c2a/rNbcZuPKVHeSdauD58/9GiL83G8tf5sVNR7huThIPXDLOSeYSiUQikYCiKHuFEOm23pM9uhLJWUxWVhaZmZlceeWVjBkzhqysrH5f87H12WQU1/Lsz6YQEeiDl07hTxekce3sJP793XG25pU7wVwikUgkkt6Ria5EchZz7733smrVKhRFYcyYMWRmZgJw9OhRrr/+ei6//HKHrpdTWs/bO0/wixnDWTo2ttN7d1+YRkpMEPeuyaRdDmGQSCQSyQBw1ie627fPsKy8oAVyc3MtMyo9HS25grZ8neG6a9cuvvzyS1auXElSUhKrVq2y9OgmJyfz2muvOXzNx9ZnE+Trxe+Xju7mW3DsCH++MI2CqiY+2lvUL3dXoqVyANry1ZIraMtXS66gLV8tuYL2fF1NnxJdRVFecbaIRCIZWO655x7Wrl3L8ePHOX78OAcOHLD06PaF3FP1fJtdxq/mjCA80MfmMYvSYpg8LIznvs2XE9MkEolE4nLsJrqKokTY+YoElg2go0QicTIbNmygtbWVxYsXW16LjY2lsbGRqqqqPl3z1a1H8fPWcc0s+8uUKYrCLQtHUlzTzNeHy/r0ORKJRCKRqKWnHt1yYA+w1+prz+mvGNerSSQSV7F06VK2bdvW7fXa2loiIiKorKzkN7/5Dfv27ePvf/97r9erbW7n0/0l/HhqAhF2enPNLEqLIT7Uj3d22V3fWyKRSCQSp9DTzmhHgcVCiIKubyiKUug6pYElJuYadys4RGRkpLsVVKMlV9CWb1fXxsZGXnrpJYxGI3/84x+dcv2XXnpJ9fGfHyihseoUB//zBnuH/4lp06bZ9fXS6/j59OE8uSGX4xWNJEUF9tvXmWipHIC2fLXkCtry1ZIraMtXS66gPV9X01OP7jNAuJ33/uECF7eQknIrKSm3ultDNZGRkS4vxC+//DLx8fFMnjyZlJQULrvsMtra2hy+Tl9dV69ezeTJk5k8eTI6nc7y/Z133unwtRxhIGLrLMyujY2NPPnkk4wcOZKdO3dyySWXuMXng71FjB2ZyPkL57J8+XIuvfRS9u7d283XzOXpCSgKfLq/xB26PaKlcgDa8tWSK2jLV0uuoC1fLbmC9nxdjdwwQtKNlStXkp6eznXXXYfRaCQ1NZX33nuP9HSbazG7jOLiYmbPns2JE/IRty3279/P9OnT8fHxISUlhcBA9/SMtrQbyCqpIyHcn9gQP4xGIyUlJRQUFDB9ehgvvTSJKVM2dTvvqld2UFbfyjd3LkBRlIEXl0gkEsmgoKcNI3oaunBWkJ//PIBmenUrKysB1z6ayMjI4MYbbwQgPz8fIQSpqakOX8eW6xVXXEFsbCz79++nsLCQd955h1deeYWdO3cyb968TktaZWZmMmHChH7+NP3z9VQqKyuJj4/n9ttv56233iIgIIBrr72WsWPHDrjLRz8UUbyrgCevnkaA3sCaNWvIyclh0aJF/OIXlRZf6Bzb5ZOH8uePM8gqqWP80NAB97aHlsoBaMtXS66gLV8tuYK2fLXkCtrzdTVnfaJbVvY2IBNda7KyslixYgXt7e0UFRWxbt06QkJCejxn3rx51NfXd3qttbWVu+66i2uvvdbyWkZGBrNmzeL555/n/vvv5/rrr2fTpk1ERkYSGxvLiy++iK+vr+XY8ePHO/3ns4eWKgez6+OPP85DDz3Ea6+9xt/+9jdWrFjBww8/PKAuj+3bxoxZCUwbn8i0adNYuHAhW7ZsYfz48ezbt7CTr3VsLxwfx18+yWRdxkmZ6PYDLflqyRW05aslV9CWr5ZcQXu+ruasT3QlnSksLCQmJoaDBw8C8NZbb7Fq1SpefvllHn74YWpra/nwww+7nbd169Zur3VdsLqlpYWamhp+97vfAeDv78/1119PfHw8AAEBAfj4nJmxn5mZydKlS532sw1WfH19ueWWW7j++uuprq4e0M8uqm7iYFEtd1+YRnx8PHv37iUhIUHVuWEBPsxIjmDDoVPcdUGai00lEolEcjaiasMIRVHCFUWZrijKfPOXq8Uk7uHgwYOdHn9PmjSJsrKyXnfKmjdvnmXSmPlr+fLlbN++3XJMVlYWU6dORaczFbsDBw4wY4ZpV7qioiKGDBnSaazmQPfoah1fX1/i4uIG9DPXZ5YCcMG4OPR6veok18zSMbHklzVwrKLRFXoSiUQiOcvptUdXUZQbgNuBBGA/MBPYASxyrZrEHWRkZDBmzBgAhBC8+eabLFmypNfz1PToZmRkMGnSJMu/Dx48yMSJEwFT0mv+HsBoNJKXl0damuzp82TWZ5aSFhfc5yXCloyN5cHPD7HhUCk3zR/pZDuJRCKRnO2oGbpwO3AOsFMIca6iKGnAX12rNXDMnr3L3QoO0ZdJYY6QkZHB5s2bWbduHYqiMHPmTJ544ok+Xaura0ZGBtOnTwdMwxiam5sJDzetYGed9IJpElxCQoJlvO5A4OrYOhNPcK1qbGNvQTW3LRrV67H2fBPCAxgbH8KGQ6c8JtH1hNg6gpZ8teQK2vLVkitoy1dLrqA9X1ejJtFtEUK0KIqCoii+QohsRVFGu9xM4hbeeecdm69XVlZy7733WnbK+vOf/+zwtZ988knL935+fhw7dszy767XS01N5dChQw5/hmTg2JpXjhBwblr/NkpcOjaW577No6KhlaiggbuxkUgkEsngR02iW6QoShjwCbBBUZRqwPNWeXciFRUVNDc3M2zYMHereAyO7pQlcR6VlZXU19eTlJTkbpVObMopJyLQh4n9XDFh6dhY/vlNHltyy/nxVMfG+EokEolE0hO9TkYTQvxICFEjhHgQ+AvwGnCZq8UGiu3bZ7B9u2lCVHl5OXfffTejR49m/fr1bjazTW5ubrexr56KllzBc30zMjIsG3jk5+cD7nc1GgVbcsuZPyoKna73zR568h0bH0JUkA+bc8udrdkn3B1bR9GSr5ZcQVu+WnIFbflqyRW05+tqVC0vpiiKHogFzM+a44ACV0kNNDU17Vx++eV89dVXzJ07l8cff5zo6Gg+//xzd6t1o7i4GICcnBw3m/SOllzBs32fe+45Pv/8c6ZOncr06dO5/fbbGT3afSOIMoprqWxsY+Ho/g1bANDpFOaNimZzbjlGo1CVOEskEolEogY1qy78FngAOAUYT78sgIl2T9IYOTlNrF27loCAABoaGlizZo27lezS2Ghahsld2706gpZcwfN9m5q
|
2022-07-11 11:29:38 +02:00
|
|
|
"text/plain": [
|
|
|
|
"<Figure size 864x288 with 3 Axes>"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
"metadata": {
|
|
|
|
"needs_background": "light"
|
|
|
|
},
|
|
|
|
"output_type": "display_data"
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"# Make a figure showing two signals with a beacon per signal\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"colors = ['y','g']\n",
|
|
|
|
"multiplier_name = ['m','n']\n",
|
|
|
|
"\n",
|
|
|
|
"\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"fig, axes = plt.subplots(3,1, sharex=True, figsize=(12,4))\n",
|
|
|
|
"axes[-1].set_xlabel(\"Time [ns]\")\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"axes[-1].set_yticks([],[])\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"for i in range(0, 2):\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
" axes[i].set_yticks([],[])\n",
|
|
|
|
" axes[i].set_ylabel(\"Antenna {:d}\".format(i+1))\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" axes[i].plot(time/ns, impulses[i])\n",
|
|
|
|
" axes[i].plot(time/ns, beacons[i], marker='.')\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
" if not True:\n",
|
|
|
|
" axes[i].plot(time/ns, full_signals[i])\n",
|
|
|
|
"\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"\n",
|
|
|
|
"# indicate timing of pulses\n",
|
|
|
|
"for i, impulse_offset in enumerate(impulses_offsets):\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" kwargs = dict(color=colors[i])\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" [ax.axvline(impulse_offset/ns, **kwargs) for ax in (axes[i], axes[-1])]\n",
|
|
|
|
"\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"# indicate timing of the beacons\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"# and annotate ticks and impulse widths\n",
|
|
|
|
"tmp_beacon_phases = beacon_init_phase + np.arange(0,2)*beacon_phase_offset\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"if not True: # mod phases\n",
|
|
|
|
" tmp_beacon_phases %= 2*np.pi\n",
|
|
|
|
"tmp_beacon_offsets = phase2time(tmp_beacon_phases, f_beacon)\n",
|
|
|
|
"\n",
|
|
|
|
"\n",
|
|
|
|
"A = np.empty(2)\n",
|
|
|
|
"B = np.empty(2)\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"for i in range(0,2):\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" kwargs = dict(color=colors[i], ls=(0, (3,2)))\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" tick_kwargs = dict(color='k', alpha=0.2)\n",
|
|
|
|
"\n",
|
|
|
|
" # indicate every period of the beacon\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" beacon_ticks = tmp_beacon_offsets[i] + [(n)*1/f_beacon for n in range(1+int((time[-1] - time[0]) * f_beacon))]\n",
|
|
|
|
"\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" [axes[i].axvline(tick/ns, **{**kwargs, **tick_kwargs}) for tick in beacon_ticks]\n",
|
|
|
|
"\n",
|
|
|
|
" # reference period in beacon\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
" [ax.axvline(tmp_beacon_offsets[i]/ns, **kwargs) for ax in (axes[i], axes[-1])]\n",
|
|
|
|
"\n",
|
|
|
|
" # annotate width between impulse and closest beacon tick\n",
|
|
|
|
" # and closest beacon tick and reference tick\n",
|
|
|
|
" closest_beacon_tick_id = np.argmin(np.abs(beacon_ticks-impulses_offsets[i]))\n",
|
|
|
|
" if closest_beacon_tick_id != 0 and beacon_ticks[closest_beacon_tick_id] > impulses_offsets[i]:\n",
|
|
|
|
" closest_beacon_tick_id -= 1\n",
|
|
|
|
" closest_beacon_tick = beacon_ticks[closest_beacon_tick_id]\n",
|
|
|
|
"\n",
|
|
|
|
" annotate_width(axes[i], f\"$A_{i+1}$\", closest_beacon_tick/ns, impulses_offsets[i]/ns, 0.7)\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" annotate_width(axes[i], f\"$B_{i+1}={multiplier_name[i]}T$\", closest_beacon_tick/ns, tmp_beacon_offsets[i]/ns, 0.4)\n",
|
|
|
|
"\n",
|
|
|
|
" A[i] = closest_beacon_tick - impulses_offsets[i]\n",
|
|
|
|
" B[i] = closest_beacon_tick - tmp_beacon_offsets[i]\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"\n",
|
|
|
|
"# annotate width between beacon reference periods\n",
|
|
|
|
"annotate_width(axes[-1], \"$t_\\phi$\", tmp_beacon_offsets[0]/ns, tmp_beacon_offsets[-1]/ns, 0.4)\n",
|
|
|
|
"\n",
|
|
|
|
"# annotate width between pulses\n",
|
|
|
|
"annotate_width(axes[-1], \"$\\Delta t$\", impulses_offsets[0]/ns, impulses_offsets[-1]/ns, 0.4)\n",
|
|
|
|
"\n",
|
|
|
|
"\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"fig.show()\n",
|
|
|
|
"if False:\n",
|
|
|
|
" fig.savefig('figures/08_beacon_sync_timing_outline.pdf')"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
|
|
|
"execution_count": 7,
|
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"name": "stdout",
|
|
|
|
"output_type": "stream",
|
|
|
|
"text": [
|
|
|
|
"\\Delta t = 74.20000000000002ns\n",
|
|
|
|
"\\Delta A = 5.7999999999999785ns\n",
|
|
|
|
"t_phi = 8.571428571428573ns\n",
|
|
|
|
"B_1 = 57.14285714285714ns = 4.0T\n",
|
|
|
|
"B_2 = 128.57142857142856ns = 9.0T\n",
|
|
|
|
"kT = 71.42857142857142ns = 5.0T\n"
|
|
|
|
]
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"t_phi = (tmp_beacon_offsets[1]-tmp_beacon_offsets[0])\n",
|
|
|
|
"Delta_A = (A[1] - A[0])\n",
|
|
|
|
"\n",
|
|
|
|
"print(\"\\Delta t = {}ns\".format( (impulses_offsets[1][0] - impulses_offsets[0][0])/ns ))\n",
|
|
|
|
"print(\"\\Delta A = {}ns\".format( Delta_A/ns ))\n",
|
|
|
|
"print(\"t_phi = {}ns\".format( t_phi/ns ))\n",
|
|
|
|
"print(\"B_1 = {}ns = {}T\".format( (B[0])/ns, (B[0]*f_beacon) ))\n",
|
|
|
|
"print(\"B_2 = {}ns = {}T\".format( (B[1])/ns, (B[1]*f_beacon) ))\n",
|
|
|
|
"print(\"kT = {}ns = {}T\".format( (B[1]-B[0])/ns, (B[1]-B[0])*f_beacon ))\n",
|
|
|
|
"\n"
|
2022-08-01 15:17:24 +02:00
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "markdown",
|
|
|
|
"metadata": {},
|
|
|
|
"source": [
|
2022-08-01 17:03:26 +02:00
|
|
|
"$\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"\\Delta t = (A_2 + B_2) - (A_1 + B_1) + t_\\phi\\\\\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"\\quad = (A_2 - A_1) + (B_2 - B_1) + t_\\phi\\\\\n",
|
|
|
|
"\\quad = (A_2 - A_1) + (nT - mT) + t_\\phi\\\\\n",
|
|
|
|
"\\quad = \\Delta A + (kT) + t_\\phi\n",
|
|
|
|
"$\n",
|
|
|
|
"\n",
|
|
|
|
", where $\\Delta A < T$ and $k \\in \\mathbb{Z}$ and $t_\\phi$ is minimisable.\n",
|
|
|
|
"\n",
|
|
|
|
"Then $\\Delta t$ can be determined by iteratively summing the signals, changing $k$, and finding the $k$ belonging to the maximum of the sums."
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
|
|
|
"execution_count": 8,
|
|
|
|
"metadata": {},
|
|
|
|
"outputs": [],
|
|
|
|
"source": [
|
|
|
|
"def find_best_integer_periods_sum(samplerate, f_beacon, ref_impulse, impulse, k_step=1):\n",
|
|
|
|
" max_k = int( len(ref_impulse)*f_beacon/samplerate )\n",
|
|
|
|
" ks = np.arange(-max_k/2, max_k/2, step=k_step)\n",
|
|
|
|
" \n",
|
|
|
|
" maxima = np.empty(len(ks))\n",
|
|
|
|
" \n",
|
|
|
|
" best_i = 0\n",
|
|
|
|
" \n",
|
|
|
|
" for i,k in enumerate(ks, 0):\n",
|
|
|
|
" augmented_impulse = time_roll(impulse, samplerate, k/f_beacon)\n",
|
|
|
|
" \n",
|
|
|
|
" maxima[i] = max(ref_impulse + augmented_impulse)\n",
|
|
|
|
" \n",
|
|
|
|
" if maxima[i] > maxima[best_i]:\n",
|
|
|
|
" best_i = i\n",
|
|
|
|
" \n",
|
|
|
|
" return ks[best_i], (ks, maxima)"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
|
|
|
"execution_count": 9,
|
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"name": "stdout",
|
|
|
|
"output_type": "stream",
|
|
|
|
"text": [
|
|
|
|
"Best k: -5\n",
|
|
|
|
"Maximum: 2.0\n"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"data": {
|
|
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAtkAAAEYCAYAAACJEeoSAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzde5Qc93XY+e+tV7/njTdAggRpUQIkAhFkWcfWw5F1aEkrQ1jtOaRoKY59Eq031m6ieFfyOjqWuHHWkTcnyjpaR9EyoSVlIzkRbTDiOhspR4QoMREpSoZJUNSQIAiQAw6AwTz73V1Vv/2juoaDwTy6e6a7BuD9nDMH6Omq7l9X93TfvnV/9yfGGJRSSimllFKbx0p6AEoppZRSSt1oNMhWSimllFJqk2mQrZRSSiml1CbTIFsppZRSSqlNpkG2UkoppZRSm0yDbKWUUkoppTaZBtlKKaWUUkptMg2ylVJKKaWU2mQaZCuluiYirxORvxSRooj8Tytc/4yIvCuBcZ0TkV/q9/2uZb1j0emYReRlETmyKYNb/T6uen7Xe75VRI+bUgo0yFbquiAivyAi/0VE5kVkRkQeE5G3JD0u4JPASWNMwRjzR8uvNMYcNMac7P+wtp6lx2KjXwJEZBjYBTy77Pe/JiKl1k9NRIIll+dEJNXhXS1/ftd8vjsY/6Z8CWrdTrX1+GZF5P8VkX2bdLsbGd+WPm6r3PbJ1msmfr2Mr7Ht31zpC2MPXn9KXdc0yFZqixORAeBh4J8DI8Ae4D6gnuS4Wm4Gnkl6EK9BbwTOGGNqS39pjPmyMSZvjMkD/zvwcHzZGDNkjOn0NbP8+d2Kz/cHWo93F3CJ6O8kadfDcVvJx5e8Xl63/EoR+e9F5PirF+VjSy734vWn1HVNg2yltr6fATDGfM0YExhjqsaYbxljngIQESMit8Ubi8ifiMjvL7l8TkT+FxF5SkTKIvKvRGSHiPzH1unr/9zKjF5DRF7fynDNtcodfmXJdd8BfhH4QitT9TMr7L+Yeet0HK3t/1cR+UkrS/mAiKTbeczLxvApEbnQuo9xEXl36/e7ReRBEZkSkReXn8Zfbb9l2/y6iHxzyeUzIvLvllx+WUQOLz0WIvJV4Cbgm63j9sklN3m4dXzmReRP48e7gjcBp1u3mxWRfysifyYi+aW3BfzVKvvH4+vk+b3m+e7m2K70+Ns51utpfeH4BvCGZY9xrbFcc7/rPD/rHruNHLe1xtvuuHrsXwMHgL9HFESHwEOrbLvu60+pG54xRn/0R3+28A8wAEwDXwbeCwwvu94Aty25/CfA7y+5fA74AbCDKAt+GfgxcARIAd8BPrPC/brAGeB3AQ/460AReN2SbU4Cf2uNsZ8DfqmbcbS2Pw3sI8rgPxY/rjYf8y8BrwNeBna3fr+fKEiwgB8Bv9d6bLcCZ4G7WtutuN8Kj+9WYK51e7uA88CFJdfNAtYqx+KXVjhWTwC7W4/3WeA3VzmuXwI+A9zSOoafAWTZNmeB/3aN56bj53fp5W6P7QrHoq1j3cbrK0v0N/KVJdevOpa17nel56eTY9fNcVtvvGu8bv4Y+ONNeJ85CUwBV4j+1t61yuP+beAU8F+AX1/+umv39ac/+vNa+NFMtlJbnDFmAfgFosDy/wamROQ/iMiODm7mnxtjLhljLgDfAx43xvyliU7f/jlRoLvczwF54B8bYxrGmO8Qla18eAMPp9NxfMEY87IxZgb4R13cd0AUwL9BRFxjzDljzAvAW4Btxpj/rfXYzhId23vW2e8qrf2KRFm7dwL/CbggIne0Ln/PGBN2MN4/Msa80nq832zd7kreSBTUfwe4zxhznzHGxFdKVGK0nygYWs1Gn99uj227t9OuEyIyBywA7wH+jyXXrTWWjdzvRo7dWvfb6bHDGPN3jDF/p81xr+VTREH9HqIvcd8UkQPLtvkN4EXgnwH/gOiLwLHlN9Tm60+pG54G2UpdB4wxzxpj/qYxZi9wiCjb+c86uIlLS/5fXeFynmvtBl5eFiSeJ/oQ7lan43h52X3v7uTOjDFniE5tfxa4LCJfF5HdRDWyu1un+udaQdrvEmXZ19pvJd8F3gW8o/X/k0QB9jtblztxccn/K6zwvIiIEL0GjgNfNMasdLr+TqLg/8U17mtDz2+3x7aD22nXB40xQ0SB68eB74rIztZ1q45lg/fb9bFb5347OnbtaJW0mFV+vr9kXI8bY4rGmLox5stE2ez3LRv7vzTG/NmrF82/NMacWOFu23n9KXXD0yBbqeuMMeanROURh1q/qhCdKo/tXL5Pl14B9onI0veJm4ALm3T77VjaKeKm1pigg8dsjPm3xphfIApgDPA5ouD9RRNNxop/CsaY962z30riIPvtrf9/l/WDbLPK79txS+vfXwJ+W0SOrrDNYeCppdntFWz4+e322LLs8XdwrNcaS9AKAAOiMz+sN5Y17ne952dDx26N+13v2HX8ujHGvMsYI6v8/MJauwKyym3+iVm7a1A7rz+lbngaZCu1xYnIHSLy2yKyt3V5H9Fp6R+0NjkF3Csitoj8MlFwtxkeB8rAJ0XElahl1weAr2/S7bfjt0Rkr4iMEGX0/rT1+7Yes0T9if+6RK3DakTZ8oCo9nmhNQEt07qdQ9Jqi7jGfiv5LtHEtowxZoKoDOaXgVHgL1fZ5xLRqfluvIkogHka+Bjw5yKya9k2h1n/VP2Gnt9uj23L4uPv8FivNR4RkWPAMK+2Nlx1LOvc73rPT9fHbp37Xe/YbeR1s9aYhkTkLhFJi4gjIr9KdGbmP3V5k+28/pS64WmQrdTWVwTeCjwuImWi4Po00QQkgL9L9AE/B/wqsNLp244ZYxrArxBNtrxCNMHqb7Qy6f3yb4FvEU2iOgvEHUTafcwp4B8Tjf8isB34XWNM0Nr/MNEp7SvA/cDgWvutdAfGmOeAElFwHdfQnwUea93PSv4A+HSrJOB/XvMIXOuNwFOt+zpBVD97Qq7uRHIn6wQ5m/D8dntsYcnjB+5e6XbiDSXqPrPisW/5poiUiGqy/xHwa8aYZ1qPca2xrPUcr/n8bPDYrXq/bRy7a8YlIl8UkS+2cb9rcYn+tuKJj/8jURnOqr2y17Hu60+p1wLRszlKqa1IRM4RdWT4z0mPRSmllOqUZrKVUkoppZTaZBpkK6WUUkoptcm0XEQppZRSSqlNpplspZRSSimlNpmT9AA229jYmNm/f3/Sw1BKKaWUUje4H/3oR1eMMdtWuu6GC7L379/Pk08+mfQwlFJKKaXUDU5Ezq92nZaLKKWUUkoptck0yFZKKaWUUmqTaZCtlFJKKaXUJrvharKVUkoppVTvNZtNJiYmqNVqSQ+l59LpNHv37sV13bb30SBbKaWUUkp1bGJigkKhwP79+xGRpIfTM8YYpqenmZiY4JZbbml7v8TKRURkn4g8IiLPisgzIvJ3V9hGROSPROSMiDwlIn8tibEqpZRSSqmr1Wo1RkdHb+gAG0BEGB0d7Thjn2Qm2wd+2xjzYxEpAD8SkW8bY36yZJv3Are3ft4K/IvWv0oppZRSKmE3eoAd6+ZxJhZkG2MmgcnW/4si8iywB1gaZB8DvmKitd9/ICJDIrKrte+KLi7U+Nz/99NeDl0pdR2ohfO8VPsBO7xDDDp7kh6OUkrdcN65vcnkfDXpYfTNfLXZUYy5JbqLiMh+4Ajw+LKr9gAvL7k80frd8v0/JiJPisiTlUqlV8NUSl0nfFPn2zOf4Yni/fzHmd+hFFxKekhKKaV64OXz5/nFtx3d0G382b/7Or/0C29d/NkznOP0U3+14bElPvFRRPLAg8DfM8YsLL96hV3MNb8w5kvAlwCOHj1qPvXLd2z6OJVS149/85N/w/zlC/zOz/4Of/TjP6JR+As+9Y7PJT0spZS6oTz77LPsGswkOob6QBrHkg2N47f+9q/zW3/71wF4+umnOXbsGO95+89ds91cxmV5jPk7a9xuoplsEXGJAuz/xxjzZytsMgHsW3J5L/BKP8amlLo+GWP4xnPf4E1jb+JXX/+rHLvtGN8+/23
|
|
|
|
"text/plain": [
|
|
|
|
"<Figure size 864x288 with 1 Axes>"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
"metadata": {
|
|
|
|
"needs_background": "light"
|
|
|
|
},
|
|
|
|
"output_type": "display_data"
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"# Make figure showing the sums depending on k\n",
|
|
|
|
"ref_impulse = impulses[0]\n",
|
|
|
|
"my_impulse = impulses[1]\n",
|
|
|
|
"\n",
|
|
|
|
"# remove 'already determined' offsets\n",
|
|
|
|
"if True:\n",
|
|
|
|
" # $t_\\phi$ offset\n",
|
|
|
|
" my_impulse = time_roll(my_impulse, samplerate, -t_phi)\n",
|
|
|
|
"\n",
|
|
|
|
" # $\\Delta A$ offset\n",
|
|
|
|
" my_impulse = time_roll(my_impulse, samplerate, Delta_A)\n",
|
|
|
|
"\n",
|
|
|
|
"best_k, (ks, maxima) = find_best_integer_periods_sum(samplerate, f_beacon, ref_impulse, my_impulse)\n",
|
|
|
|
"print(\"Best k: {:0g}\".format(best_k))\n",
|
|
|
|
"print(\"Maximum: {}\".format(maxima[np.where(ks == best_k)][0]))\n",
|
|
|
|
"\n",
|
|
|
|
"\n",
|
|
|
|
"# Make figure\n",
|
|
|
|
"fig, axes = plt.subplots(1, 1, sharex=True,figsize=(12,4))\n",
|
|
|
|
"if not hasattr(axes, 'ndim'):\n",
|
|
|
|
" axes = [axes]\n",
|
|
|
|
"\n",
|
|
|
|
"axes[0].set_title(\"Sum of impulses with $kT$ offsets. Best offset: ${:.0f}*T$\".format(best_k))\n",
|
|
|
|
"axes[-1].set_xlabel(\"Time [ns]\")\n",
|
|
|
|
"\n",
|
|
|
|
"if not True:\n",
|
|
|
|
" i=0\n",
|
|
|
|
" axes[i].set_ylabel(\"Reference\")\n",
|
|
|
|
" axes[i].plot(time/ns, ref_impulse, label=\"reference\")\n",
|
|
|
|
" axes[i].plot(time/ns, my_impulse, label='impulse')\n",
|
|
|
|
" axes[i].legend()\n",
|
|
|
|
"\n",
|
|
|
|
"axes[-1].set_ylabel(\"Sum\")\n",
|
|
|
|
"\n",
|
|
|
|
"best_maximum = np.max(maxima)\n",
|
|
|
|
"axes[-1].axhline(best_maximum, alpha=0.7)\n",
|
|
|
|
"\n",
|
|
|
|
"for i, k in enumerate(ks, 0):\n",
|
|
|
|
" sample_offset = int(k*1/f_beacon*samplerate)\n",
|
|
|
|
" augmented_impulses = np.roll(my_impulse, sample_offset)\n",
|
|
|
|
" \n",
|
|
|
|
" summed_impulse = ref_impulse + augmented_impulses\n",
|
|
|
|
" if True or k%2 == 1:\n",
|
|
|
|
" axes[-1].plot(time/ns, summed_impulse, label='k={:.0f}'.format(k),\n",
|
|
|
|
" alpha=0.1 + 0.9*1/(1+4*abs(best_maximum-maxima[i]))\n",
|
|
|
|
" )\n",
|
|
|
|
" \n",
|
|
|
|
"axes[-1].legend()\n",
|
|
|
|
"fig.show()\n",
|
|
|
|
"\n",
|
|
|
|
"del ref_impulse\n",
|
|
|
|
"del my_impulse"
|
2022-07-11 11:29:38 +02:00
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "markdown",
|
|
|
|
"metadata": {},
|
|
|
|
"source": [
|
|
|
|
"## 1. Solve it"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 10,
|
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"data": {
|
|
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAEWCAYAAABIVsEJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dd3iUZdbH8e9JT+iQECAQQiD0TgTBAggoqIBtVWxgw+7q6u7q6qrru7rqKrZ1VexrR9eCKyJFsNKC9JaE0EJLaKEGUs77xzzZHWMaYSbPTHI+1zVXZp76mwTmzP2U+xZVxRhjjDleIW4HMMYYE5ysgBhjjKkWKyDGGGOqxQqIMcaYarECYowxplqsgBhjjKkWKyAmqInI5SIyw0fbmisi11VjvSQRUREJ80WOYHcifxMReUhE3vF1JuMfVkCMq5wP7b0iElmd9VX1XVU909e5TPXZ36TusAJiXCMiScBpgAJjXA1jfMJaYXWLFRDjpquA+cCbwPiKFhSRCSKSJSIHRGSDiFzuNf0Hr+VURG4UkQynZfOCiIgzL1REnhKRXc42bq3o0JOIXCMia5ztfC0ibSt5P9eIyDYR2S4id3ltJ0RE7hGR9SKyW0SmiEhTr/kficgOEckTke9EpJvXvGgn8yZn/g8iEu3MGyMiq0Rkn9OS6+K13kYRuVtEljvrfSgiURX8bn8UkeedZdeKyDCv+Y1E5DXnfW0Vkb+KSGipdZ8WkT3AQ2X8TQaJyCJn24tEZJDXvHYi8q3zd50JxFbyOzYBxAqIcdNVwLvO4ywRiS9rIRGpBzwHjFLVBsAgYGkF2z0XOAnoBVwMnOVMvx4YBfQG+gLnlbcBETkP+BNwARAHfA+8X8n7GQqkAGcC94jIcGf67c6+BgOtgL3AC17rfeWs1xz4Gc/vo8STQD8877kp8AegWEQ6OnnucPJNA74QkQivdS8GRgLtgJ7AhAqyDwCy8HyAPwh84lXk3gIKgQ5AH+f9XVfGus2BR7w36mzjSzx/v2bAJOBLEWnmLPIesNjZ7/9RyRcJE2BU1R72qPEHcCpQAMQ6r9cCd5azbD1gH3AhEF1q3gTgB6/XCpzq9XoKcI/z/BvgBq95w53lw5zXc4HrnOdfAdd6LRsCHAbalpEvydlOZ69pTwCvOc/XAMO85rV03ntYGdtq7GyrkbPPI0CvMpb7MzClVL6twBDn9UbgilJ5Xirn9zsB2AaI17SFwJVAPHDU+/cOjAPmeK27uby/ibONhaXmz3OWScRTmOp5zXsPeMftf5/2qNrDWiDGLeOBGaq6y3n9HuV8+1TVQ8AlwI3AdhH5UkQ6V7DtHV7PDwP1neetgC1e87yfl9YWeNY5PLQP2AMIkFDBOt7b2+Tsr2Rbn3ptaw1QBMQ7h9Uecw5v7cfzwQ+eb+SxQBSwvox9tXL2AYCqFjv7985X3u+hLFvV+QQvlb8tEI7n916S/2U8rY2y3neFOb22neDM2+v8fb3nmSBhBcTUOOcY/sXAYOfY/w7gTqCXiPQqax1V/VpVR+D59r4WeKUau94OtPZ63aaCZbfgaa009npEq+pPFazjvb1EPN/qS7Y1qtS2olR1K3AZMBZPa6gRntYMeIrVLiAfaF/Gvrbh+XD3LOw5z9MGTyukOhJKzhWVyr8FTwsk1it7Q1Xt5rVsRV16/yKn17a34vl7NHEOUXrPM0HCCohxw3l4voF3xXM+ojfQBc95hqtKLywi8c4J43p4PswOOusfrynAb0UkQUQaA3+sYNmXgHtLTmg7J5J/U8n2/ywiMc46VwMfem3rkZKT8CISJyJjnXkNnPe0G4gBHi3ZmNOqeB2YJCKtnNbKQPFc8jwFOEdEholIOHCXs52KClxFmgO3i0i48z67ANNUdTswA3hKRBo6FwS0F5HBVdzuNKCjiFwmImEicgmev/t/VHUTkAb8RUQiRORUYHQ18xsXWAExbhgPvKGqm1V1R8kD+AdweRlXRYXg+YDchudQ0mDg5mrs9xU8H4bLgSV4PtwKKaMYqeqnwOPAB86hpZV4TsBX5FsgE5gNPKmqJTfTPQtMBWaIyAE8V54NcOb9C89hm63Aameet7uBFcAiPO/9cSBEVdcBVwDP42mpjAZGq+qxSn8LZVuA50T+Ljwnwi9S1d3OvKuACCffXuBjPC3BSjnbOBfP3283nosAzvU6dHkZnt/FHjwn7/9VzfzGBfLLw57G1B0iMgrPieXKLs+t1URkAp6LB051O4sJLtYCMXWGc0/F2c6hlAQ833g/dTuXMcHKCoipSwT4C57DMEvwXA31gKuJjAlidgjLGGNMtVgLxBhjTLXUqY7PYmNjNSkpye0YxhgTVBYvXrxLVeNKT69TBSQpKYm0tDS3YxhjTFARkTJ7CLBDWMYYY6rFCogxxphqsQJijDGmWqyAGGOMqRYrIMYYY6rF1QIiIq+LSI6IrCxnvojIcyKS6QzN2ddr3njxDFuaISI2ipkxxtQwt1sgb+IZcrM8o/D0EJoCTARehP8Ok/kgnl48+wMPikgTvyY1xhjzC67eB6Kq34lIUgWLjAX+5YyUNl9EGotIS2AIMFNV9wCIyEw8haiyMauNH6gqWbsOsXXvEfbnF5B3pID9RwopVqVFwyhaNo4ioXE0LRpFERkW6nZcY4yPBPqNhAn8crjMbGdaedN/RUQm4mm9kJhog535ytod+5m3fjcLN+xh4YY97D5U+TAUoSFCv7ZNGNqpOWd0bk7H+Pr8chA8Y0wwCfQCUtani1Yw/dcTVScDkwFSU1Ot58gTkF9QxBfLtvH2/E0sz84DoHWTaAZ3imNAu6Ykx9WnUXQ4jaLDaRgVjghsz8tn274jbNt3hPW5h/guPZfHp6/l8elrSWgczYX9WjNhUBJN60W4/O6MMccr0AtINr8cZ7o1nlHpsvEcxvKePrfGUtUxOQfyefX7DUxJ28K+wwV0aF6fv4zpxoiu8bRqHF3huu1i69Eu9n9DXt8zqjM78vKZsy6HGat28NzsDCZ/t55LT0rkutPa0bpJjL/fjjHGR1zvzt05B/IfVe1exrxzgFuBs/GcMH9OVfs7J9EXAyVXZf0M9Cs5J1Ke1NRUtb6wqq6gqJi3ftrIM7MyOFJQxJld47lyYFsGJjfz2aGnzJwDvPxtFp8u2YoCv+nXmj+O7EwTa5EYEzBEZLGqpv5qupsFRETex9OSiAV24rmyKhxAVV8Sz6fUP/CcID8MXK2qac661wB/cjb1iKq+Udn+rIBU3Q8Zu3joi1Vk5hxkcMc4HhjdlfZx9f22v237jvDK91m8PW8TjaLDeXBMN0b3bGnnSIwJAAFZQGqaFZDKHTxayJ8/W8mnS7aS2DSGB87tyrAuzWvsg3zN9v3c8+/lLMvOY2inOP56fg8SKjlMZozxLysgWAGpzNod+7n53Z/ZuOsQt52Rwk1D2hMVXvOX3RYVK2/+tJEnv15HiMBTF/dmZPcWNZ7DGONRXgFx+0ZCEyA+StvCeS/8yIH8Qt697mTuHNHRleIBnst9rz21HTPuPJ0O8Q248Z3FTJqZTnFx3fmyY0wwsAJSxxUUFXPPv5fz+4+X06dNE768/VQGtm/mdiwA2jSN4cOJJ3NRv9Y8NzuDiW+ncSC/wO1YxhiHFZA67MixIib+K40PFm3hlqHteee6ATRvEOV2rF+ICg/l7xf15C9jujFnXS7nvfAjm3YfcjuWMQYrIHVW3uECrnxtAd+m5/K3C3rw+7M6ExoSmFc8iQjjByXxzrUD2HPoGBe/PI/MnINuxzKmzrMCUgft3J/PJZPnsTw7jxcu68u4/sHRxcvA9s34YOJAiorh0snzWLtjv9uRjKnTrIDUMdl7D3PRSz+xZc9h3rj6JEb1aOl2pOPSqUUDPrzhZMJCQrh08nxWOF2qGGNqnhWQOiT3wFGufG0heYcLeO/6kzmlQ6zbkaqlfVx9ptwwkHoRYVz26nx+3rzX7UjG1ElWQOqIvCMFjH99ITvy8nnj6pPo1aax25FOSGKzGKbcOJCm9SKY8PpCMnYecDuSMXWOFZA64MixIq57axEZOQd
|
|
|
|
"text/plain": [
|
|
|
|
"<Figure size 432x288 with 1 Axes>"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
"metadata": {
|
|
|
|
"needs_background": "light"
|
|
|
|
},
|
|
|
|
"output_type": "display_data"
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"# single out one period of the beacon\n",
|
|
|
|
"beacon_samplerate = samplerate # MHz\n",
|
|
|
|
"beacon_time = np.arange(0, 1/f_beacon, 1/beacon_samplerate)\n",
|
|
|
|
"ref_beacon = sin_delay(f_beacon, beacon_time, phase=0, t_delay=0)\n",
|
|
|
|
"\n",
|
|
|
|
"# .. and show beacon period\n",
|
|
|
|
"fig, ax = plt.subplots()\n",
|
|
|
|
"ax.set_title(\"A single beacon period\")\n",
|
|
|
|
"ax.set_xlabel(\"Time [ns]\")\n",
|
|
|
|
"ax.set_ylabel(\"Amplitude\")\n",
|
|
|
|
"ax.plot(beacon_time/ns,ref_beacon)\n",
|
|
|
|
"fig.show()"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
|
|
|
"execution_count": 11,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [],
|
|
|
|
"source": [
|
|
|
|
"### correlation\n",
|
|
|
|
"def correlation_and_lag(sig1, sig2, mode=\"full\", normalise=False):\n",
|
|
|
|
" corr = signal.correlate(sig1, sig2, mode=mode)\n",
|
|
|
|
" if normalise:\n",
|
|
|
|
" corr /= np.max(corr)\n",
|
|
|
|
"\n",
|
|
|
|
" lags = signal.correlation_lags(sig1.size, sig2.size, mode=mode)\n",
|
|
|
|
" \n",
|
|
|
|
" return corr, lags\n",
|
|
|
|
"\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"def find_best_lag(sig1, sig2, fix_one_short=False, fix_positive=False, subtract_means=True, **corr_kwargs):\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" if subtract_means:\n",
|
|
|
|
" sig1 -= np.mean(sig1)\n",
|
|
|
|
" sig2 -= np.mean(sig2)\n",
|
|
|
|
" \n",
|
|
|
|
" corr, lags = correlation_and_lag(sig1, sig2, **corr_kwargs)\n",
|
|
|
|
" lag_id = corr.argmax()\n",
|
|
|
|
"\n",
|
|
|
|
" lag = lags[lag_id]\n",
|
|
|
|
" \n",
|
|
|
|
" if fix_one_short:\n",
|
|
|
|
" # for some reason it is always one short\n",
|
|
|
|
" if lag > 0:\n",
|
|
|
|
" lag += 1\n",
|
|
|
|
" elif lag < 0:\n",
|
|
|
|
" lag -= 1\n",
|
|
|
|
"\n",
|
|
|
|
" # turn negative lag into positive\n",
|
|
|
|
" if fix_positive and lag < 0:\n",
|
|
|
|
" lag += len(sig2)\n",
|
|
|
|
" \n",
|
|
|
|
" return lag, (corr, lags)\n",
|
|
|
|
"\n",
|
|
|
|
"def find_beacon_phase_delay(samplerate, f_beacon, reference_beacon, delayed_beacon, **lag_kwargs):\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" Return phase delay of `beacon` with respect to `reference_beacon`.\n",
|
|
|
|
" Note that the returned value can be off by a multiple of $2\\pi$.\n",
|
|
|
|
" \n",
|
|
|
|
" Parameters\n",
|
|
|
|
" ==========\n",
|
|
|
|
" samplerate : float\n",
|
|
|
|
" Samplerate of both reference_beacon and delayed_beacon\n",
|
|
|
|
" f_beacon : float\n",
|
|
|
|
" Frequency of the beacons\n",
|
|
|
|
" reference_beacon : ndarray\n",
|
|
|
|
" The beacon to use as a reference\n",
|
|
|
|
" beacon : ndarray\n",
|
|
|
|
" The beacon to find the delay for\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" \n",
|
|
|
|
" calc_lag, _ = find_best_lag(reference_beacon, delayed_beacon, **lag_kwargs)\n",
|
|
|
|
" \n",
|
|
|
|
" time = calc_lag / samplerate\n",
|
|
|
|
" \n",
|
|
|
|
" return 2*np.pi* f_beacon * time"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "markdown",
|
|
|
|
"metadata": {},
|
|
|
|
"source": [
|
|
|
|
"##### 1.1 Beacon Phase Delays"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 12,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"name": "stdout",
|
|
|
|
"output_type": "stream",
|
|
|
|
"text": [
|
2022-08-01 17:03:26 +02:00
|
|
|
"Beacon delays [ns] \\pm k*14.285714285714285ns: [0. 5.6]\n"
|
2022-07-11 11:29:38 +02:00
|
|
|
]
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"beacon_phase_delays = np.array([\n",
|
|
|
|
" find_beacon_phase_delay(beacon_samplerate, f_beacon, beacons[0], beacon)\n",
|
|
|
|
" for beacon in beacons\n",
|
|
|
|
"])\n",
|
|
|
|
"\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"beacon_time_delays = phase2time(beacon_phase_delays, f_beacon)\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"\n",
|
|
|
|
"print(\"Beacon delays [ns] \\pm k*{}ns: {}\".format(1/f_beacon/ns, beacon_time_delays/ns))"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 13,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
2022-08-01 15:17:24 +02:00
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"name": "stdout",
|
|
|
|
"output_type": "stream",
|
|
|
|
"text": [
|
|
|
|
"0 0.0\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"1 2.463008640414398\n"
|
2022-08-01 15:17:24 +02:00
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"data": {
|
2022-08-01 17:03:26 +02:00
|
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEpCAYAAACN9mVQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOydd3gU17m432/Ve0ESSEiIXoQKCNFcwMa42xhcwI5LHCdx4vSb+0vim1zn+jrOTbeT2Cl2HNtxBzeMbdzihk0zXSB6EaiAUAH1unt+f5wVLGIlrbRV4rzPM4+0M2fO+WZndr5TviJKKQwGg8Fg6A6LvwUwGAwGQ2BjFIXBYDAYesQoCoPBYDD0iFEUBoPBYOgRoygMBoPB0CNGURgMBoOhR4yiMBgMBkOPGEVhMBgMhh4xisJwBiLytIg86GLZYhGZ722ZHNpTItIoIr90s56PRKRFRD73lGxd6i8RkaneqNtg8AdGUbiJ/WXZLCINInJCRN4WkQx/yzWIyVNK/cydCpRS84BvekieMxCRBCAV2NVLuTtF5KJujn1HRDaKSKuIPN1DHePsCu+5bo6Hicg/ReSwiNSLyBYRubJLmZEistL+7B4TkUdFJNjh+Cf2Nhrs2x6HYw1dNquIPNKX63B2De7K3ZPMA1luf2IUhWe4VikVjX5BVACP9FLeMIAQkftF5H4Xi+cA+5VSLd3U9Q0RWXT6o9zt8LmTcuBB4Mle2voLsKGH48FACTAXiAPuA5aJyEiHMn8FjqOf3Sn2st/qUs93lFLR9m1C506HfdHAUKAZeLmP1+HsGjwht1OZB7Lc/sQoCg9ifzm8AmR17hORNBF5VUQqReSQiHzP4di9InLA3vvY6fjCEJEMEXnNfl61iDzqcGySvfdxUkSKRGSBw7FiEfl/IlIoIrUislREwruTWUSmishmuwxLgXCHY93K7qSenq7lRyLyapfyj4jIH+3//0REyuzn7hGRS3r8ol281v7W6ya5wA57+5Ei8oL9Pkbbjz8JjAF+APwfYAPecKxAKfWaUmo5UN1dIyJyM3AS+LC7MkqpRqXU/UqpYqWUTSn1FnAImOZQbBSwTCnVopQ6BrwLTO7bJQNwI/oF+Jmr19HdNRi5Aw+jKDyIiEQCS4B19s8W4E1gGzAcuAT4gYhcbj/lAHAhuvfxv8BzIpIqIkHAW8BhYKT93JfsdYbY63wfSAG+CzwvIo69j8XAFeiHMhe4sxt5Q4HlwLNAIrpXdYOLsnfF6bXYjz0HXCEi8fa6g+3f07N2ub8DTFdKxQCXA8XdtOEMp9fqgXr7Sy6wXURGAZ8De4AblFINDmWUw1+rw2eXEJFY4AHgP/t43lBgPFDksPtPwM12pTYcuBL98nLkVyJSJSKrpZvpMuDLwDPKxSijfbmGfsrtiswDWW7fopQymxsb+uXTgO5hdKCHrTn2YzOBI13K/xfwVDd1bQWuA2YDlUCwkzIXAscAi8O+F4H7HeS5zeHYb4G/d9PeHLu84rBvDXrY3avs9rbm93QtDp/fAb5u//8aYKf9/7HoHt18IKSX71oBY7u07/Rae6sXrVA+76Gtt+z39CTQYt86P7/Vw3lrgb+je5LXOTn+DeB6e/sX2z8v7KauB4Gnnez/E/AT+//3A8+58JyGAP8GHuuyfxKwyf7sKuDpLs/DTCAGCEO/VOuBMV3qGIFWeKNcvQ5Xr6E/crsi80CW2x+bGVF4hoVKqXj0Df4O8KmIDAMygTT7FNFJETkJ/BQ9L4qI3CEiWx2OZQNJQAZwWCnV4aStNKBEKWVz2HcY3evv5JjD/01ANM5JA8qU/Sl1qIveZO9KD9fSyb+A2+z/34YexaCU2o+ehrkfOC4iL4lIWjfyOsPptbpbr1LqGqVUvP2+/hr4dednpdQ1zs4REUFf9yK0wnqjaxml1GNKqddOf1SPKT3N4RIiMgWt/B7uwzkW9Pfdhn4+Hfe/B7wGRKHvVwLwGwd51yul6pVSrUqpfwGrgau6NHEHWuke8uQ19FduF2UeyHL7HKMoPIhSymp/CViBC9ALW4ccXjDxSqkYpdRVIpIJ/AP9IA2xv5B2AGI/b4Q4WJ84UA5k2B+6TkYAZf0Q+Sgw3P6Cc6yLnmTvWkkv19LJciBXRLLRI4rnOw8opV5QSl2AVk4KhxeVO3ir3h4YZf87H/hPESnoQbanlVKf9KONi9DTkUdE5Bjw/4AbRGSzs8L2e/tPtIK/QSnV7nA4Ed0pedT+cqoGnqLnl5PizPsK+oX7L09eg4fldibzQJbb5xhF4UFEcx26l7AL+AKosy+qRohIkIhki8h0dI9CoaeYEJGvoHuj2M87CvxaRKJEJFxEzrcfWw80Aj8WkRD7POa12Ncw+sha9BD4eyISLCLXAzMcZOhO9q70dC3AGQv9LwBfKKWO2MtOEJF5IhKGnt5pRitat/BkvUovUN7vQtFcoFAptR24G3jdYZ3GZez3IhwIAoLs97+z0/A4ejF8in37O/A2eg3GGX9DT3dcq5RqdjyglKpCT5HdY28zHj3lsc0uR7yIXN7Zvojcip6ufM9B1vPQo1lHq6HersOVa+iX3K7IPJDl9hv+nvsa6Bt6nrwZvU5Rj+5J3+pwPA29hnAMOIFe6J5vP/ZLoAaoAh4CPgW+Zj82At0Lr7Yf/7NDnZPtZWuBncCiLvLMd/h8Pz3MYQMFwBa77Evt24O9yd61rZ6uxaH8BWiF8hWHfblopVRvP/8tIK0bWZ2tUTi91t7qpfc1infs99TZ9k4359wH/K3L5/VAeB+fqfvt1+q43d9D2ee6yP1T+/+dI6mWLvI7Pp9TgE/s97cK/eJMsR9LRpuA1qPXZtYBl3Zp/zHgWXeuw8k19FtuV2QeyHL7a+tcRDEYvI6IjAB2A8OUUnX9OL8FaEUrzfvckOMDYBZ6ZOMLk1mDYUBjFIXBJ9jXVB4CYpVSd/lbHoPB4DrOFksNBo8iIlFoj/XDaJ8Hg8EwgDAjCoPBYDD0iLF6MhgMBkOPGEVhMBgMhh4xisJgMBgMPWIUhcFgMBh6xCgKQ0AjHspq521EZLycToLzNX/LYzB4EqMoDAOBM7LaicgUEdkkIk32v1P6W7GDIurMKvZEL+VvFpFd9nMOiMiFAEqpvUonwvmsp/O9iZzOkJbhsG++iBT7SybD4MAoCsOAQnQOjTfQOS4S0EHd3rDv7y956nRWsW5HAyJyKTqw4FfQ4aDnAAfdaNcbNKJDhxgMHsMoCsNA4yK0o+gflY7A+Wd0hM15Pmj7f4EHlFLrlM5gVqaU6k/UXm/yZ+AWERnr7KD0nB3xaRH5i+i87/Uisl5ExjgcnygiH4hIjeiMgYt9cD2GAMAoCsNAYzI6Qqujp2gh7qWTXCU60f1rcmZ+41OIzjpYACSLyH4RKRWRR0UkwtVGRKeLPdnd1tdy3VCGDvl+v5P2XcmOeAtaISYA+9HBHju96z9AR/9NsZf7q4gEfBpPg/sYRWEYaESjo+Y6UoueCuoPc9E5Biaic328Jc7zgAxFZy27EZ1lcAowFfhvVxtSSjkmPzpr62u5HvgVcK2Tl/gs9Pf3a6VUm1LqI3RU3VscyrymlPpC6aRZz9uvE3QOkWKl1FNKqQ6l1GbgVfv3YRjkGEVhGGg0ALFd9sWiwzP3GaXUKvtL8yTwfXTyoUlOinbmFnhEKXVU6dwCDxEgGcgcUUpVAo+iczs74k52xExgZpeRza3AMI8KbwhIjKIwDDSK0JnyHDN/5XJmAnt3cJpVTCl1Aii1H+8XIvJTB+uqs7a+luuF36Fzck9z2OdOdsQS4NMuo5topdQ9LspjGMAYRWEYaHyCzlT3PREJE5HOnMQf9bUiEZlsN7UNEpFo4A/ol+aubk55CviuiKSISAI6J/dbrranlPo/B+uqs7a+luulrZP26/mxw253siO+BYwXkdvt54aIyHQRcTb6MgwyjKIwDCiUUm3AQnS+45PAXcBC+/7O3vg7LlY3FJ3Rrw5t5joSuEbZ8xw7qesX6Cxke9HKZAv2xd4A5U84pH+1f0cLgCvR2dX+CtyhlNrdW0V
|
2022-08-01 15:17:24 +02:00
|
|
|
"text/plain": [
|
|
|
|
"<Figure size 432x288 with 1 Axes>"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
"metadata": {
|
|
|
|
"needs_background": "light"
|
|
|
|
},
|
|
|
|
"output_type": "display_data"
|
|
|
|
}
|
|
|
|
],
|
2022-07-11 11:29:38 +02:00
|
|
|
"source": [
|
|
|
|
"# Make a figure showing the corrected beacons\n",
|
|
|
|
"fig, ax = plt.subplots(1,1, sharex=True)\n",
|
|
|
|
"ax.set_xlabel(\"Time [ns]\")\n",
|
|
|
|
"ax.set_ylabel(\"Amplitude [au]\")\n",
|
|
|
|
"ax.set_title(\"Beacon delays [ns] $\\pm$ $k*{}$\\n{} == {}\".format(1/f_beacon/ns, beacon_time_delays/ns, None))\n",
|
|
|
|
"\n",
|
|
|
|
"for i, _ in enumerate(beacons):\n",
|
|
|
|
" l = ax.plot(time/ns, beacons[i], label=\"ch {}\".format(i), ls ='--', alpha=0.5)\n",
|
|
|
|
" \n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" print(i, beacon_phase_delays[i], )\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
" if True:\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" corrected_beacon = sin_delay(f_beacon, beacon_time, phase=beacon_phase_delays[i], t_delay=0)\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" else:\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
" corrected_beacon = sin_delay(f_beacon, beacon_time, t_delay=beacon_time_delays[i], phase=0)\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" \n",
|
2022-08-01 15:17:24 +02:00
|
|
|
" ax.plot(beacon_time/ns, 2*beacon_amplitude*corrected_beacon, label='ch {} corrected'.format(i), color=l[0].get_color())\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
" \n",
|
|
|
|
" # indicate start of uncorrected beacons\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
" ax.axvline(phase2time((2*np.pi+beacon_time_delays[i])%2*np.pi, f_beacon)/ns, color=l[0].get_color())\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"\n",
|
|
|
|
"ax.legend(ncol=2)\n",
|
|
|
|
"ax.margins(y=0.3)\n",
|
|
|
|
"if True:\n",
|
|
|
|
" ax.set_xlim(time[0]/ns - 1, time[2*samplerate//f_beacon]/ns)\n",
|
|
|
|
"\n",
|
|
|
|
"fig.show()\n",
|
|
|
|
"\n",
|
|
|
|
"\n",
|
2022-08-01 17:03:26 +02:00
|
|
|
"#ax.plot((double_signal_time) * ns, signal_2(double_signal_time + calc_shift), 'r--', label='Recovered')"
|
2022-07-11 11:29:38 +02:00
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "markdown",
|
|
|
|
"metadata": {},
|
|
|
|
"source": [
|
|
|
|
"##### 1.2 Impulse vs beacon delays\n",
|
|
|
|
"\n",
|
|
|
|
"Find the delay within a single beacon period"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 14,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [],
|
|
|
|
"source": [
|
|
|
|
"#def find_beacon_impulse_delay(samplerate, f_beacon, impulse, init_phase=0):\n",
|
|
|
|
"def find_beacon_impulse_phase_delay(samplerate, f_beacon, reference_beacon, impulse, **lag_kwargs):\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" Return phase delay of `beacon` with respect to `reference_beacon`.\n",
|
|
|
|
" Note that the returned value can be off by a multiple of $2\\pi$.\n",
|
|
|
|
" \n",
|
|
|
|
" Parameters\n",
|
|
|
|
" ==========\n",
|
|
|
|
" samplerate : float\n",
|
|
|
|
" Samplerate of both reference_beacon and delayed_beacon\n",
|
|
|
|
" f_beacon : float\n",
|
|
|
|
" Frequency of the beacons\n",
|
|
|
|
" reference_beacon : ndarray\n",
|
|
|
|
" The beacon to use as a reference\n",
|
|
|
|
" beacon : ndarray\n",
|
|
|
|
" The beacon to find the delay for\n",
|
|
|
|
" \"\"\"\n",
|
|
|
|
" \n",
|
|
|
|
" calc_lag, _ = find_best_lag(reference_beacon, impulse, **lag_kwargs)\n",
|
|
|
|
" \n",
|
|
|
|
" return 2*np.pi* f_beacon * calc_lag / samplerate"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 15,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"name": "stdout",
|
|
|
|
"output_type": "stream",
|
|
|
|
"text": [
|
2022-08-01 17:03:26 +02:00
|
|
|
"Beacon Impuls delays: $[-454.90261624 -921.11496603]$ns\n"
|
2022-07-11 11:29:38 +02:00
|
|
|
]
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"impulse_beacon_phase_delays = np.empty( len(impulses) )\n",
|
|
|
|
"\n",
|
|
|
|
"for i, _ in enumerate(impulses):\n",
|
|
|
|
" impulse_beacon_phase_delays[i] = find_beacon_impulse_phase_delay(\n",
|
|
|
|
" beacon_samplerate, f_beacon, \n",
|
|
|
|
" ref_beacon, impulses[i]\n",
|
|
|
|
" )\n",
|
|
|
|
"\n",
|
|
|
|
"print(\"Beacon Impuls delays: ${}$ns\".format(impulse_beacon_phase_delays/f_beacon/ns))"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"cell_type": "code",
|
2022-08-01 17:03:26 +02:00
|
|
|
"execution_count": 16,
|
2022-07-11 11:29:38 +02:00
|
|
|
"metadata": {},
|
|
|
|
"outputs": [
|
|
|
|
{
|
|
|
|
"data": {
|
|
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEKCAYAAAAfGVI8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAVt0lEQVR4nO3df6zdd33f8ecLQxItUDDYk5DjEENdkhS1CRwZKqTBBEkc/rCRilqnQyQsnTVGYIOtUhCagsw/FDbRsaUlpvWASsOBVNpuO1CaEhhbh8HHIwuxK8PFUOKZKQaHaF2yZE7e++N8s3t8fa/v19fH9zj383xIRz7fz/fz+Z7P/eje8/L35ydVhSSpXc+bdgckSdNlEEhS4wwCSWqcQSBJjTMIJKlxBoEkNW7JIEiyJ8kjSR5aZH2SfCrJbJIHk7x2bN3NSb7fvW6eZMclSZPRZ4/gs8DWM6y/EdjcvXYCfwCQ5KXAHcDrgS3AHUnWnktnJUmTt2QQVNU3gBNnqLId+HyN7ANekuTlwA3AfVV1oqoeBe7jzIEiSZqC509gGxuAh8eWj3Zli5WfJslORnsTXHrppa+78sorJ9AtSWrHgQMHflpV65fTdhJBkAXK6gzlpxdW7QZ2AwwGgxoOhxPoliS1I8lfL7ftJK4aOgpsHFu+DDh2hnJJ0gVkEkEwA7yru3roDcBjVfUT4F7g+iRru5PE13dlkqQLyJKHhpJ8AXgzsC7JUUZXAr0AoKo+DXwZeBswCzwOvLtbdyLJR4H93aZ2VdWZTjpLkqZgySCoqpuWWF/AexdZtwfYs7yuSZJWgncWS1LjDAJJapxBIEmNMwgkqXEGgSQ1ziCQpMYZBJLUOINAkhpnEEhS4wwCSWqcQSBJjTMIJKlxBoEkNc4gkKTGGQSS1DiDQJIa1ysIkmxNcjjJbJLbF1j/ySQPdK/vJfn52Lqnx9bNTLLzkqRz12eqyjXAncB1jCak359kpqoOPVunqj4wVv99wLVjm3iiqq6ZXJclSZPUZ49gCzBbVUeq6ilgL7D9DPVvAr4wic5Jks6/PkGwAXh4bPloV3aaJK8ANgH3jxVfkmSYZF+Sty/SbmdXZ3j8+PGeXZckTUKfIMgCZbVI3R3APVX19FjZ5VU1AH4L+L0krzptY1W7q2pQVYP169f36JIkaVL6BMFRYOPY8mXAsUXq7mDeYaGqOtb9ewT4OqeeP5AkTVmfINgPbE6yKclFjL7sT7v6J8mrgbXAN8fK1ia5uHu/DngjcGh+W0nS9Cx51VBVnUxyG3AvsAbYU1UHk+wChlX1bCjcBOytqvHDRlcBdyV5hlHofGz8aiNJ0vTl1O/t6RsMBjUcDqfdDUl6TklyoDsfe9a8s1iSGmcQSFLjDAJJapxBIEmNMwgkqXEGgSQ1ziCQpMYZBJLUOINAkhpnEEhS4wwCSWqcQSBJjTMIJKlxBoEkNc4gkKTGGQSS1LheQZBka5LDSWaT3L7A+luSHE/yQPf67bF1Nyf5fve6eZKdlySduyWnqkyyBrgTuI7RRPb7k8wsMOXk3VV127y2LwXuAAZAAQe6to9OpPeSpHPWZ49gCzBbVUeq6ilgL7C95/ZvAO6rqhPdl/99wNbldVWSdD70CYINwMNjy0e7svl+PcmDSe5JsvFs2ibZmWSYZHj8+PGeXZckTUKfIMgCZfNnvP9T4Iqq+hXgL4DPnUVbqmp3VQ2qarB+/foeXZIkTUqfIDgKbBxbvgw4Nl6hqn5WVU92i58BXte3rSRpuvoEwX5gc5JNSS4CdgAz4xWSvHxscRvwV937e4Hrk6xNsha4viuTJF0glrxqqKpOJrmN0Rf4GmBPVR1MsgsYVtUM8P4k24CTwAnglq7tiSQfZRQmALuq6sR5+DkkScuUqtMO2U/VYDCo4XA47W5I0nNKkgNVNVhOW+8slqTGGQSS1DiDQJIaZxBIUuMMAklqnEEgSY0zCCSpcQaBJDXOIJCkxhkEktQ4g0CSGmcQSFLjDAJJapxBIEmNMwgkqXG9giDJ1iSHk8wmuX2B9R9McqibvP6rSV4xtu7pJA90r5n5bSVJ07XkDGVJ1gB3AtcxmoN4f5KZqjo0Vu07wKCqHk/yHuDjwG92656oqmsm3G9J0oT02SPYAsxW1ZGqegrYC2wfr1BVX6uqx7vFfYwmqZckPQf0CYINwMNjy0e7ssXcCnxlbPmSJMMk+5K8faEGSXZ2dYbHjx/v0SVJ0qQseWgIyAJlC050nOSdwAB401jx5VV1LMkrgfuTfLeqfnDKxqp2A7thNGdxr55Lkiaizx7BUWDj2PJlwLH5lZK8FfgwsK2qnny2vKqOdf8eAb4OXHsO/ZUkTVifINgPbE6yKclFwA7glKt/klwL3MUoBB4ZK1+b5OLu/TrgjcD4SWZJ0pQteWioqk4muQ24F1gD7Kmqg0l2AcOqmgE+AbwQ+FISgB9X1TbgKuCuJM8wCp2PzbvaSJI0Zam6sA7JDwaDGg6H0+6GJD2nJDlQVYPltPXOYklqnEEgSY0zCCSpcQaBJDXOIJCkxhkEktQ4g0CSGmcQSFLjDAJJapxBIEmNMwgkqXEGgSQ1ziCQpMYZBJLUOINAkhpnEEhS43oFQZKtSQ4nmU1y+wLrL05yd7f+W0muGFv3oa78cJIbJtd1SdIkLBkESdYAdwI3AlcDNyW5el61W4FHq+oXgU8Cv9u1vZrRHMe/DGwFfr/bniTpAtFnj2ALMFtVR6rqKWAvsH1ene3A57r39wBvyWjy4u3A3qp6sqp+CMx225MkXSCWnLwe2AA8PLZ8FHj9YnW6ye4fA17Wle+b13bD/A9IshPY2S0+meShXr1f/dYBP512Jy4QjsUcx2KOYzHn1ctt2CcIskDZ/BnvF6vTpy1VtRvYDZBkuNwJmFcbx2KOYzHHsZjjWMxJMlxu2z6Hho4CG8eWLwOOLVYnyfOBFwMneraVJE1RnyDYD2xOsinJRYxO/s7MqzMD3Ny9fwdwf1VVV76ju6poE7AZ+PZkui5JmoQlDw11x/xvA+4F1gB7qupgkl3AsKpmgD8C/jjJLKM9gR1d24NJvggcAk4C762qp5f4yN3L/3FWHcdijmMxx7GY41jMWfZYZPQfd0lSq7yzWJIaZxBIUuMMAklqnEEgSY0zCCSpcQaBJDXOIJCkxhkEktQ4g0CSGmcQSFLjDAJJalyfqSr3JHlkscliMvKpbl7iB5O8dmzdzUm+371uXqi9JGm6+uwRfJbRfMOLuZHR46U3M5pl7A8AkrwUuIPRbGZbgDuSrD2XzkqSJm/JIKiqbzB6tPRitgOfr5F9wEuSvBy4Abivqk5U1aPAfZw5UCRJU9BnqsqlLDSn8YYzlJ9mfM7iSy+99HVXXnnlBLolSe04cODAT6tq/XLaTiIIzmm+Yjh1zuLBYFDD4bKn3pSkJiX56+W2ncRVQ4vNS+x8xZL0HDCJIJgB3tVdPfQG4LGq+gmjqS2vT7K2O0l8fVcmSbqALHloKMkXgDcD65IcZXQl0AsAqurTwJeBtwGzwOPAu7t1J5J8FNjfbWpXVZ3ppLMkaQr6TF5/0xLrC3jvIuv2AHuW1zVJ0krwzmJJapxBIEmNMwgkqXEGgSQ1ziCQpMYZBJLUOINAkhpnEEhS4wwCSWqcQSBJjTMIJKlxBoEkNc4gkKTGGQSS1DiDQJIaZxBIUuN6BUGSrUkOJ5lNcvsC6z+Z5IHu9b0kPx9b9/TYuplJdl6SdO76TFW5BrgTuI7RhPT7k8xU1aFn61TVB8bqvw+4dmwTT1TVNZPrsiRpkvrsEWwBZqvqSFU9BewFtp+h/k3AFybROUnS+dcnCDYAD48tH+3KTpPkFcAm4P6x4kuSDJPsS/L2Rdrt7OoMjx8/3rPrkqRJ6BMEWaCsFqm7A7inqp4eK7u8qgbAbwG/l+RVp22sandVDapqsH79+h5dkiRNSp8gOApsHFu+DDi2SN0dzDssVFXHun+PAF/n1PMHkqQp6xME+4HNSTYluYjRl/1pV/8keTWwFvjmWNnaJBd379cBbwQOzW8rSZqeJa8aqqq
|
|
|
|
"text/plain": [
|
|
|
|
"<Figure size 432x288 with 2 Axes>"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
"metadata": {
|
|
|
|
"needs_background": "light"
|
|
|
|
},
|
|
|
|
"output_type": "display_data"
|
|
|
|
}
|
|
|
|
],
|
|
|
|
"source": [
|
|
|
|
"# Make a figure showing the corrected beacons\n",
|
|
|
|
"fig, axes = plt.subplots(len(impulses),1, sharex=True)\n",
|
|
|
|
"axes[-1].set_xlabel(\"Time [us]\")\n",
|
|
|
|
"ax.set_title(\"Beacon Impuls delays: ${}$ns\".format(impulse_beacon_phase_delays/f_beacon/ns))\n",
|
|
|
|
"\n",
|
|
|
|
"for i, beacon in enumerate(beacons):\n",
|
|
|
|
" ax.set_ylabel(\"Amplitude [au]\")\n",
|
|
|
|
" ax.plot(time, beacon, label=\"ch {}\".format(i))\n",
|
|
|
|
"\n",
|
2022-08-01 15:17:24 +02:00
|
|
|
"ax.plot(beacon_time, corrected_beacon, label='phase corrected (overlaps ch 0)')\n",
|
2022-07-11 11:29:38 +02:00
|
|
|
"\n",
|
|
|
|
"ax.legend()\n",
|
|
|
|
"\n",
|
|
|
|
"fig.show()"
|
|
|
|
]
|
|
|
|
},
|
|
|
|
{
|
|
|
|
"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.7.6"
|
|
|
|
}
|
|
|
|
},
|
|
|
|
"nbformat": 4,
|
|
|
|
"nbformat_minor": 4
|
|
|
|
}
|