mirror of
https://gitlab.science.ru.nl/mthesis-edeboone/m-thesis-introduction.git
synced 2024-11-13 10:03:32 +01:00
411 lines
201 KiB
Text
411 lines
201 KiB
Text
|
{
|
||
|
"cells": [
|
||
|
{
|
||
|
"cell_type": "markdown",
|
||
|
"metadata": {},
|
||
|
"source": [
|
||
|
"# Translation of a discrete 1D signal\n",
|
||
|
"\n",
|
||
|
"\n",
|
||
|
"#### wave equation\n",
|
||
|
"\n",
|
||
|
"$$\n",
|
||
|
"\\frac{\\partial^2 u}{\\partial t^2}\n",
|
||
|
" = \n",
|
||
|
" v^2 \\frac{\\partial^2 u}{\\partial x^2}\n",
|
||
|
"$$\n",
|
||
|
"\n",
|
||
|
"\n",
|
||
|
"##### simple sine solution\n",
|
||
|
"\n",
|
||
|
"$$\n",
|
||
|
"u(x, t) = \\sin\\left( \\frac{2 \\pi x}{\\lambda} -2\\pi f t\\right)\n",
|
||
|
"= \\sin\\left( 2 \\pi f \\, ( \\frac{x}{v(f)} -t) \\right)\n",
|
||
|
"$$\n",
|
||
|
"\n",
|
||
|
"$$\\lambda = \\dfrac{v}{f}$$\n",
|
||
|
"##### plane wave per mode\n",
|
||
|
"\n",
|
||
|
"\n",
|
||
|
"\n",
|
||
|
"$$ \n",
|
||
|
"u_f(x,t) \n",
|
||
|
"= \n",
|
||
|
" A \\exp(- 2\\pi i f \\left[ \\frac{x}{v(f)} - t \\right] )\n",
|
||
|
" +\n",
|
||
|
" B \\exp(- 2\\pi i f \\left[ \\frac{x}{v(f)} + t \\right] )\n",
|
||
|
"$$\n"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "code",
|
||
|
"execution_count": 1,
|
||
|
"metadata": {},
|
||
|
"outputs": [],
|
||
|
"source": [
|
||
|
"%matplotlib inline\n",
|
||
|
"\n",
|
||
|
"import numpy as np\n",
|
||
|
"import scipy.fft as ft\n",
|
||
|
"import matplotlib.pyplot as plt\n",
|
||
|
"import matplotlib.gridspec as gridspec\n",
|
||
|
"import matplotlib.ticker as tck\n",
|
||
|
"import scipy.interpolate as interp\n",
|
||
|
"\n",
|
||
|
"rng = np.random.default_rng()"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "code",
|
||
|
"execution_count": 2,
|
||
|
"metadata": {},
|
||
|
"outputs": [],
|
||
|
"source": [
|
||
|
"class TravelSignal:\n",
|
||
|
" \"\"\"\n",
|
||
|
" Model an arbitrary signal that can be translated to another position and time.\n",
|
||
|
" \"\"\"\n",
|
||
|
"\n",
|
||
|
" def __init__(self, signal, sample_rate, t_0 = 0, x_0 = 0, periodic=True):\n",
|
||
|
" assert periodic, \"Non-periodic signal is not yet supported\"\n",
|
||
|
"\n",
|
||
|
" self.raw = signal\n",
|
||
|
"\n",
|
||
|
" self.sample_rate = sample_rate # Hz\n",
|
||
|
" self.sample_length = len(self.raw)\n",
|
||
|
" self.time_length = self.sample_length*sample_rate # s\n",
|
||
|
" \n",
|
||
|
" self.x_0 = x_0\n",
|
||
|
" self.t_0 = t_0\n",
|
||
|
"\n",
|
||
|
" # choose interpolation method\n",
|
||
|
" if not interp1d_kw:\n",
|
||
|
" self.interp_f = None\n",
|
||
|
"\n",
|
||
|
" else:\n",
|
||
|
" # offload interpolation to scipy.interpolate\n",
|
||
|
" import scipy.interpolate as interp\n",
|
||
|
"\n",
|
||
|
" interp1d_kw_defaults = {\n",
|
||
|
" \"copy\": False,\n",
|
||
|
" \"kind\": 'linear',\n",
|
||
|
" \"assume_sorted\": True,\n",
|
||
|
" \"bounds_error\": True\n",
|
||
|
" }\n",
|
||
|
"\n",
|
||
|
" if self.periodic:\n",
|
||
|
" interp1d_kw_defaults['bounds_error'] = False\n",
|
||
|
" interp1d_kw_defaults['fill_value'] = (self.raw[-1], self.raw[0])\n",
|
||
|
" \n",
|
||
|
" # merge kwargs\n",
|
||
|
" if interp1d_kw is not True:\n",
|
||
|
" interp1d_kw = { **interp1d_kw_defaults, **interp1d_kw }\n",
|
||
|
"\n",
|
||
|
" self.interp_f = interp.interp1d(\n",
|
||
|
" np.arange(0, self.sample_length),\n",
|
||
|
" self.raw,\n",
|
||
|
" **interp1d_kw\n",
|
||
|
" )\n",
|
||
|
" else:\n",
|
||
|
" self.interp_f = None\n",
|
||
|
" \n",
|
||
|
" def __len__(self):\n",
|
||
|
" return self.sample_length\n",
|
||
|
" \n",
|
||
|
" def __call__(self, t_f = None, x_f = None, **kwargs):\n",
|
||
|
" \"\"\"\n",
|
||
|
" Allow this class to be used as a function.\n",
|
||
|
" \"\"\"\n",
|
||
|
" return self._translate(t_f, x_f, **kwargs)[0]\n",
|
||
|
" \n",
|
||
|
" def _translate(self, t_f = None, x_f = None, t_0 = None, x_0 = None, velocity = None):\n",
|
||
|
" \"\"\"\n",
|
||
|
" Translate the signal from (t_0, x_0) to (t_f, x_f) with optional velocity.\n",
|
||
|
" \n",
|
||
|
" Returns the signal at (t_f, x_f)\n",
|
||
|
" \"\"\"\n",
|
||
|
" \n",
|
||
|
" if t_0 is None:\n",
|
||
|
" t_0 = self.t_0\n",
|
||
|
" \n",
|
||
|
" if velocity is None:\n",
|
||
|
" velocity = 1\n",
|
||
|
"\n",
|
||
|
"\n",
|
||
|
" ## spatial offset\n",
|
||
|
" if x_f is None:\n",
|
||
|
" spatial_time_offset = 0\n",
|
||
|
" else:\n",
|
||
|
" x_f = np.asarray(x_f)\n",
|
||
|
" if x_0 is None:\n",
|
||
|
" x_0 = self.x_0\n",
|
||
|
"\n",
|
||
|
" spatial_time_offset = np.sum(np.sqrt( (x_f - x_0)**2 )/velocity)\n",
|
||
|
"\n",
|
||
|
" ## temporal offset\n",
|
||
|
" if t_f is None:\n",
|
||
|
" temporal_time_offset = 0\n",
|
||
|
" else:\n",
|
||
|
" t_f = np.asarray(t_f)\n",
|
||
|
" \n",
|
||
|
" if t_0 is None:\n",
|
||
|
" t_0 = self.t_0\n",
|
||
|
" \n",
|
||
|
" temporal_time_offset = t_f - t_0\n",
|
||
|
"\n",
|
||
|
" # total offset\n",
|
||
|
" total_time_offset = spatial_time_offset + temporal_time_offset\n",
|
||
|
"\n",
|
||
|
" # periodic signal\n",
|
||
|
" n_offset = (total_time_offset * sample_rate ) % self.sample_length\n",
|
||
|
"\n",
|
||
|
" # offload to scipy interpolation\n",
|
||
|
" if self.interp_f:\n",
|
||
|
" amplitude = self.interp_f(n_offset)\n",
|
||
|
" else:\n",
|
||
|
" # self written interpolation\n",
|
||
|
" n_offset_eps, n_offset_int = np.modf(n_offset)\n",
|
||
|
" n_offset_int = n_offset.astype(int)\n",
|
||
|
"\n",
|
||
|
" if True:\n",
|
||
|
" amplitude = (1-n_offset_eps) * self.raw[n_offset_int] \\\n",
|
||
|
" + n_offset_eps * self.raw[(n_offset_int + 1) % self.sample_length]\n",
|
||
|
"\n",
|
||
|
" else:\n",
|
||
|
" amplitude = self.raw[n_offset_int]\n",
|
||
|
"\n",
|
||
|
" \n",
|
||
|
" return amplitude, total_time_offset"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "markdown",
|
||
|
"metadata": {},
|
||
|
"source": [
|
||
|
"# Testing:"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "markdown",
|
||
|
"metadata": {},
|
||
|
"source": [
|
||
|
"### Interpolation: Sine wave"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "code",
|
||
|
"execution_count": 3,
|
||
|
"metadata": {},
|
||
|
"outputs": [
|
||
|
{
|
||
|
"data": {
|
||
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA8MAAAEGCAYAAABBz0cYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeXxU5dXA8d+9syQhIQQSAgESSCAgSyBsYQcBxQW1rbRWigtVUKu41L3VVuurtmqpG1UrUaQuUPcFFdn3JSRh2CFAFrKHQBYSSGa59/3jJpMFgihhZhLO108+MHfuTB5MZuae5znPOYqu6wghhBBCCCGEEBcT1dsDEEIIIYQQQgghPE2CYSGEEEIIIYQQFx0JhoUQQgghhBBCXHQkGBZCCCGEEEIIcdGRYFgIIYQQQgghxEXH7O0BNLewsDC9R48e3h6GEEIIIYQQQogLICUlpVjX9Y7n+zytLhju0aMHycnJ3h6GEEIIIYQQQogLQFGUrOZ4HkmTFkIIIYQQQghx0ZFgWAghhBBCCCHERUeCYSGEEEIIIYQQFx0JhoUQQgghhBBCXHQkGBZCCCGEEEIIcdHxajCsKMq7iqIUKYqyu4n7L1UUpUxRFFvN1189PUYhhBBCCCGEEK2Pt1eG3wOu/JFz1uu6Hl/z9cyPPmNFIWQnNX1/dhKsn3v2c4QQQgghhBBCeNa5xGrZSXRtq3Rujm/n1T7Duq6vUxSlR7M+aXk+LLgKBt8MHaJBNYNqAdUEZTmw6XXQnGCywBXPQ2QCWIPALxj8gsDsDznbIHM99Bhn3C+EEEIIIYQQ4ufLTqqLsboOA0clVJ+A6gqwnzBisGV/AZfDiOFG3QMhkeByGvGb5oCSTNj+AZ2DlK7NMSSvBsPnaJSiKDuAPOBhXdf3ND5BUZQ7gDsAhkaoxv+slAVnf1aXHb57+Ax3qIBW+8TQLQHCYiEwHAI7QlA4BIbBiQIoPgixUyBqxHn9A4UQQgghhBCixcraDGlLIbQXBHWCyqNQWQSVxcbfjx2C3FRAP7fn0xyw8ZULOmTw/WA4Feiu63qFoihXA18CsY1P0nX9beBtgGFdzDrmALjpU+gy2JhZ0FzG/9CcZPjsNmN2QTXDlGcguJsxI2GvgOpyOLQKsjYCOui6MftQkmX8EHXX6SNc/09oFwlhvY2Zi3aREBJl/HnyGBzdD9HjZYVZCCGEEEII0TJlrIcD30NIN7AEQlk2lGYbfxYfgsrCMz/OHABBHY3FSncgrBirw7GXgV9bsLY1/izLhh+eqMniNcO0BdBtmJHRq5qM+C3PBh9MQ6fiHKPqs1N0vVme5+cPwEiTXqLr+oBzODcTGKbrenFT5wy7JFJPXv5Z08Fn/eX5M52TnQQLrzNWjk1WuPVr4zxNg6pSIyje/G9I/S/GD1SBjn2M9OqybCMAPn3k0H0MRI00guawWOPLr+2Pj0cIIYQQQgghLpTaeKT7WGjf3ch+LU6r+7Ng9+nBrqJC2y7GYqC9wjgH3Tg+dCaMvs/IqvULqvseZ4qxmhrL2WKj7CS69R+Zm1OudTvff7pPB8OKonQGCnVd1xVFSQA+xVgpbnLQw4YN05OTk89vUD83YAawVxp7kze8Cjs+wj0DEtABqsoari63CYVTJcYKtGqGa1+BuN+A2e/8xi+EEEIIIYQQTTl5HAp3w/5vIWn+mTNgLW2MtGddg8I9uIPdUXNg8l+NFVto3kD3HCmKkqLr+rDzehK8HAwrirIIuBQIAwqBpwALgK7rbymKMgf4A+AETgEP6rq+6WzP2SzB8Ln4OQFzRDyUZNTNsuz5Egp2NHycaoawPtA5DiIGGn867cZ5snoshBBCCCGEOFdHtsL+JcbCnLMKCnYZX2VHznCyAr2vhBF3GlmsbbuAqp5bsOvhbNdWEQxfCB4Lhs/FTwqYzTD+UXCcgoKdxi/pifyG5ysmGHEXxP3aCJJrZ2OEEEIIIYQQoqoMclMgexscXA652xreH9YbOtcsuHWOM/bnfnyrR1d1m4MEw03wqWD4XJztF6viKKx4Cmz10q1rmf2NAmHdhhkVr1UzHN3nM7+gQgghhBBCiAsoawvs/RIsAUZdo5xko3hvbV2jwI5GRWcw0psnPAaXPn768/hYoHsuJBhuQosLhn9M47SEaYnG33OSIScJ8ncYt2spJhj3kLFxvV2ztN8SQgghhBBCeJumGft8MzfAvq/hyOa6+6xBEDXKCGa7DYeuQ+DogXPby9sCSTDchFYXDMPZZ2uc1bD0cUhewGmrxx1ijMf0GAfR46D0SIub9RFCCCGEEOKik50EGeugXTc4VWpcw2dtNIrvAvi3h6qavysmmPgEjH/ozM/TCq//mysY9vU+wwKMX9ymfnnNfjBoOtgW1c36TP2X8eLIWG8U6UpdWHOyYvxhssCMTyFmgkeGL4QQQgghhDgH5fmQ9DZsfMWo4lwrJAr6TDUWuHqMhfK8hqu+0ePO/HxniyOErAy3Gk3N+mguoyDXqmfh0Iq644rJeNH0uhx6XWb0SlYUz49bCCGEEEKIi5XTDtlbjev0QyuMNOgGVBg9B6b83+mPbaWrvudC0qSbcNEGwz+m/t5j1QR9r4XCvUbRLYB2kdBrMoR0N1Kve02+6F5UQgghhBBCXHD7lsCOj4z05/ydYD8BqgWiRhqLVEGdYckDrXKvb3ORNGnx00QmGC+kxrNHpdlweKVRen3Hx+A8aRxf9yKMvh9G3Q1B4d4btxBCCCGEEC2ZrkPRPqPf787/wbFDNXco0OdqGDwDoseDX9u6x4TGXLSrvp4kK8OiztqXYPXzQL39CSgQOQIumWp8hfa8qFMyhBBCCCGEaFLtdXLUaOP2/iWw/1soyQAUCO4K5bmAbmxbnPSE0QlG/CSyMiyaX8wEWD+3LiXjmlegNMt4ES//i/EV0t14AesamPwkbUMIIYQQQgiAzI3w/i8btj01WSF6Aoy5H/pcZXR3qV/4qkcTha+ER8jKsGioqVXfkiw48B1seQtKM+uOR4+Hq16C8Es8PlQhhBBCCCG8yuU0rp13fwa7PgFnVd19A6YZi0v+wQ0fI1mW500KaDVBguELLDsJFl5rVL5TqCv5Ht4PBlwP/a+XVGohhBBCCNH6uFOgx4DuMgLgvV/ByWKwBhlbCzPXG91cpPDVBSXBcBMkGPaA+oFuSJTxJrD7c8jeYtwf2gtKMiWVWgghhBBCtA5HthoLQi47UBM/mQOgz5XGYlDs5WAJkAUhD5FguAk9+vfQv1z9JfHh8QDYimwkFyYzrNMw4sPjT7t9pnPEz1SWA3u+hM3/hhN5dcf7/RJ+8W/wC/Le2IQQQgghhPipSrONCtBb3oCTx+qO978erntdrm/P08+J1WxFNoZfMjzXcdzR7Xy/f6sLhttEt9Fjn44lLiwOk2rCdtSGpmv4m/x5dPij/D3p7zg0BybFxPWx1+Nv8mfxgcW4NBdWk5Xnxz6PQ3OQVpLG2G5jGdbJmHA4lx+UqJGdBAuvMVKpAdDBEgj9roNBNxozZarJq0MUQgghhBDijKpPwL5vwPYRZG4AdOgcB0X7azIfJQW6sXMJYrfmbWVj3kZi28fSrW03CioL+OvGv2J32VEVlTFdx7AhdwMu3YVJMXH3oLtJiEhg1rJZVLuqURWV6OBoMsoySHs6jVMZp5TzHXerC4YDogP0Xk/3IrxNOJqmUVxVDIBJMTEiYgSb8jY1+ViTYiIiMIKcihz3sTbmNnRr240j5UeMH5SqMiZiDBvzN+LSXFhUC69Neo0xXcdIwFxfbYpI97HGm8aORbDnC6guN0rKD7wBOvU3KupJGokQQgghhPCmrM2w/QM4kQ9HNoPjJHSIgUHTjevW9j0u6hTos8U5ALf/cLt7wfGq6Ksot5ezMXcjmq5hNVkBqHJVNXjOAWED2Fu8F61BW9c6oyJGkRCRwOupr7vPCbIEUeGo4NDTh5olGG51rZUUFPxN/sydMBeA2ctm49AcWFQLl0VdRmphqvv2G5e9wSnnKR5c8yBOzYlFtdCnQx9yK3LR0VFQ6BXSC1VVsbvsaGhomsb
|
||
|
"text/plain": [
|
||
|
"<Figure size 1152x288 with 1 Axes>"
|
||
|
]
|
||
|
},
|
||
|
"metadata": {
|
||
|
"needs_background": "light"
|
||
|
},
|
||
|
"output_type": "display_data"
|
||
|
}
|
||
|
],
|
||
|
"source": [
|
||
|
"# Sine wave\n",
|
||
|
"sample_rate = 1e3 # Hz\n",
|
||
|
"f_sine = 10\n",
|
||
|
"time = np.arange(0, 1, 1/(sample_rate)) #s\n",
|
||
|
"time2 = np.arange(0, 1, 1/(sample_rate*2)) #s\n",
|
||
|
"phase_offset = 0\n",
|
||
|
"\n",
|
||
|
"signal_func = lambda t: np.cos(2*np.pi*f_sine*t + phase_offset)\n",
|
||
|
"\n",
|
||
|
"mysignal = TravelSignal(signal_func(time), sample_rate)\n",
|
||
|
"\n",
|
||
|
"fig, ax = plt.subplots(figsize=(16,4))\n",
|
||
|
"ax.plot(time, signal_func(time), '.-', label=\"Raw signal\")\n",
|
||
|
"ax.plot(time, mysignal(time)+0.5, '.-', label=\"MySignal +0.5\")\n",
|
||
|
"ax.plot(time2, mysignal(time2)+0.25, '.--', label=\"MySignal(interp) +0.25\")\n",
|
||
|
"ax.set_ylabel(\"Amplitude\")\n",
|
||
|
"ax.set_xlabel(\"Time\")\n",
|
||
|
"ax.set_xlim(0, 1/f_sine)\n",
|
||
|
"ax.legend()\n",
|
||
|
"plt.show();"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "markdown",
|
||
|
"metadata": {},
|
||
|
"source": [
|
||
|
"### Extrapolation: Sine wave"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "code",
|
||
|
"execution_count": 4,
|
||
|
"metadata": {},
|
||
|
"outputs": [
|
||
|
{
|
||
|
"data": {
|
||
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA7gAAAEGCAYAAABPZjl3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOx9d3xU15n2c2ZGXUIIFapQAUx1EEgWuIMbYDvYad4Y23G8STZlv+xmd3+/xHG+2BuvkzjJ2kn227SNU1xiO8YtcsNUU0wREgxdoIKkEQJGEpJQn3Lv98e5fdqd0dwi9jz++TfnltF9ec88M/fc932fl/A8DwYGBgYGBgYGBgYGBgaGiQ6H1QYwMDAwMDAwMDAwMDAwMCQDbIHLwMDAwMDAwMDAwMDAcEWALXAZGBgYGBgYGBgYGBgYrgiwBS4DAwMDAwMDAwMDAwPDFQG2wGVgYGBgYGBgYGBgYGC4IuCy2oBko6CggC8tLbXaDAYGBgYGBgYGBgYGBgYDUF9f383zfGG4Y1fcAre0tBR1dXVWm8HAwMDAwMDAwMDAwMBgAAghbZGOsRRlBgYGBgYGBgYGBgYGhisCbIHLwMDAwMDAwMDAwMDAcEWALXAZGBgYGBgYGBgYGBgYrgiwBS4DAwMDAwMDAwMDAwPDFQG2wGVgYGBgYGBgYGBgYGC4ImDpApcQ8kdCiJcQcjzC8VWEkH5CiFv4/3GzbWRgYGBgYGBgYGBgYGCYGLA6gvtnAGtjnLOb5/kK4f8nTbDJWNT9GXjxU/TVjvDUAq9uAH5/i31trPszte/VB6i9dgPzYXKw5Qngv5bRVzuCzfP4wXyYHHhqgd3P2Nc+BgYGBgYGE2FpH1ye53cRQkqttMFU1P0ZePef6bh5O9B7Frj9B5aapIKnFvjjWoAP0u1z9fS16ouWmRSCLU8AH/9C3j79AfD3m4DiautsUiKcD+02z8rPIWA/HwLAG18Bjr1Gxx//ArjcCXzm99bapATjyvjBuJIcHPwj8P6/ATwAVxrwcI297GNgYGBgYDAZVkdw9eBaQsgRQsgHhJDF4U4ghPwDIaSOEFLX1dVltn36cepv6u2Pf2mvJ+57fiHfbIrY/2trbAkHTy31mRJ8EDjysjX2hMORV0J9aLd53vtf6m0+qF4IWQ1Prby4FXHsNXv5cM/PGVfGi4nAlXA+tBtX3vtXgOcAcEBwDGjdbbVVDAwMDAwMlsLuC9xDAEp4nl8K4P8BeDvcSTzP/w/P81U8z1cVFhaaamBcWHiPZgdvn5slTy1w+v3Q/d1n7HPDeeQV0DCFBl1nTDclIjqPhNnJ22dh4akFLjWH7m/4wD7zvPvn4ffbiisfhO5nXIkP5w6H2WkzrvS2hO63FVeehXqeCVB6o1XWMDAwMDAw2AK2XuDyPH+Z5/lBYfw+gBRCSIHFZiWOqi8CJdep99nlZinSDbGdFuHnj4bf37bPHj701AKd9eGP2WVhsSfC4hGcPRYWnlrgTJjFI8C4Eg/CLh5hL66cj2CjXbiy+9kIB+zElU3qffPXsvRkBgYGBob/9bD1ApcQMo0QQoRxNai9PdZaNU7c9gOo3W6Tm6XL5yMfO73J+ptiT61c5xgCm/jwyCuRj7Xvt4cPT2+KfHzQBun9ERePgG3mub8z8jG7cOW8O8JBm/hwInDlzIeRj9uSKw7g+m9ZZQ0DAwMDA4NtYHWboFcA7AMwnxDSQQj5EiHka4SQrwmnfBbAcULIEQD/BeDzPM9HuvudGCiuBhbcqdlJLDFFhZFezQ6FTTxnfV1X625EXvjAHjecfR2aHXb0IafYYYPPnRZaemdP1ZxgA5uHL2l22HGebc6V3jbNDjv6cIJxZcE6Fr1lYGBgYGCAxQtcnufv53l+Os/zKTzPz+J5/g88z/+W5/nfCsf/m+f5xTzPL+V5fiXP83uttDdpmHu7envaUmvsEOGpBToOytuOFODqzylO4IGMfNPNUiFjinr76vsA4pS3G7dYG/Xx1AIt2+VtW/pQc/2rPwc4FELqVvsQACbNUG+X3azetgNXlGnodpzn9MnqbTty5exH8rYdfTgRuJKl0ZuYe4c1djAwMDAwMNgMtk5RvmIx0gNVROBCpHRCk9C6W6FmSoDlDwJFCxQnEMFmC9H2sWLDQe1b9qC8K+i3NurTuhvgAsKG0odE3mf1PA8r55BQ+5Z/Qd5ltQ8BoHWPYsMBDHeDcSVOaH1YtABYukHeZfU8t+4GuHA+tBFXhpRRbptypU35vNdh/eeOgYGBgYHBJmALXCtQeiPgTJG3D79sbTRAFa3ggWkVofusjKh4aoETb8nbzhTqwxnLFSdx1toYzoelNyqiPrz186xKWxXmVBURtdiHnlr14syZQpXH7cSV9DzFhk25omxHJnJlVpXiJMaVmAj6FRs25YpygSvOMwMDAwMDAwNb4FqC4mp7RR/b9yk2hEiAnaLM2ujosg3UhyOaiKSVEQxtdHSkR5jnB+TdVs9zqyYKHjLPFkfOtNHRZRuo8riduKLNJLAlVzQ+tB1XuhUbCq4s/by82+p5jvmdaAOuSDXCinlmYGBgYGBgYAtcy2CXaICnFjj+urwtRgLsFGUOF/EJt9/KiIpDSSWFLdOXKfZbPM9nd8rbynm2S+RMVWetmGc7ceWkohU340pi4JXiTQpbZtokyuypVT/IUHFFrGW2mCvhMgkYGBgYGBgYALAFrnVQRgOIhfVTkSI+2igzF7AuoqLyjcJX2oiKlVGp7ibFhtZGAVbPszY6Ks7zJ+6Tz7Nynvs9io0I82y1DxlXxo+eCcyVJZ+Rz7NynnvPKjZY/S0DAwMDA4MSbIFrFVRRH2JdtKL0RoAIN77OVLUYzdL7Id0UO5zW1XgVr5THylqz0hsBVxodEwt9CACBEcEOB7VJaaOoYOtwWefD0hshzaV2niu/KI+tnGep7pGE+pBxRR+Udel25Yp/WLDDzlwRoJ3nqi/JY0u54hMGGq4wMDAwMDAwsAWuZSiuBm55nI55Dtj0qHXpblI/xTC9M4kNPiIDFxQbChuLq4G1PxF2W+hDTy1wUhD2IQ5g7dPqejhxUWRlC2eOg+y7cPOsiO5ZAU8tsP83dOzQ+LC4Grj5O3RsOVfE9NpoXLGwZ+rlcwobtFx5mo6t5krD+3QcjSvR+vgaDa3AlBKEQBUJtwKeWuDgH+jY4Qz1IQMDAwMDw/9y2GD18r8Y0lN4HgiMWZPu1rob0k2cNuWudbd8Qx/0A0deNt08AOrrckG1jcrUPCt9KIpg8ZzaJqUPuYB1Pqz7gzzW+rB1t7z4tkrcp3U3wAkLC56PknJpNVcEROWKz7p5PrZRGPBhuHJJHlvpQzH9NxpXghZy5eBz8jgcV2AnrnAsPZmBgYGBgUEDtsC1EpnKNEGLRFVSs+Uxr7HBDgJEnlqgeYe8rU1dzLCBD5XXtKsPT7wpb080HwJAZoFiwyIbnenyOOw8WyxA5KkFzioWXBNtnu3ClVM18vZE8yEDAwMDAwMDW+BaCjuIqvQ0KzY0NhRXA1d/Vt62QlRFGdUJ1w7DDj4c6lJshPGh1QJEkcSRRKgEiCzyYSRxpHDHrZrnSxHEkQDqz8Wfkrct40qseRYPW+TDQWW5QRgfVijaajGuhEcsrjAwMDAwMPwvB1vgWgk7iKpIadKO8GIllY/IYytEVaIJvojHRQEiy8RzxDq9CIIvVgsQld4o14dG9GEqHVvlw5Ib5LFSHElE6Y2Qvq6s4kpgTBhE4IrVAkTRhMTE4w6RKw5r5lm5eAznwwqFzZZxJYYPreZK+mRhwASmGBgYGBgYwoEtcK0GsVCQxlMLHH6RjrXCPiLsIKoSTRypuBpY+Q065iwQz/HUArt+RsfhRHNEWC3iFE0cqbgauP1J+TwrBIiGuxUbEQSGrObK0b/ScSSuSCnKgG25co2wCOeC1nDl41/QcVSuWCzWFU10r7gauOX7wmGLvm8+/C4dE8IEphgYGBgYGMKALXCtROtuQd0WFqbjieJIEYR9oolQmYGWj+SxVvBFBFHUPgZ95trYulutuhrJh7yFPlTNYQQ
|
||
|
"text/plain": [
|
||
|
"<Figure size 1152x288 with 1 Axes>"
|
||
|
]
|
||
|
},
|
||
|
"metadata": {
|
||
|
"needs_background": "light"
|
||
|
},
|
||
|
"output_type": "display_data"
|
||
|
}
|
||
|
],
|
||
|
"source": [
|
||
|
"# Sine wave\n",
|
||
|
"sample_rate = 1e3 # Hz\n",
|
||
|
"f_sine = 10\n",
|
||
|
"time = np.arange(0, 1, 1/(sample_rate)) #s\n",
|
||
|
"time2 = np.arange(3/4, 2, 1/(sample_rate*2)) #s\n",
|
||
|
"phase_offset = 0\n",
|
||
|
"\n",
|
||
|
"signal_func = lambda t: np.cos(2*np.pi*f_sine*t + phase_offset)\n",
|
||
|
"\n",
|
||
|
"mysignal = TravelSignal(signal_func(time), sample_rate)\n",
|
||
|
"\n",
|
||
|
"fig, ax = plt.subplots(figsize=(16,4))\n",
|
||
|
"ax.plot(time, signal_func(time), '.-', label=\"Raw signal\")\n",
|
||
|
"ax.plot(time, mysignal(time)+0.5, '.-', label=\"MySignal +0.5\")\n",
|
||
|
"ax.plot(time2, mysignal(time2)+0.25, '.--', label=\"MySignal(interp) +0.25\")\n",
|
||
|
"ax.set_ylabel(\"Amplitude\")\n",
|
||
|
"ax.set_xlabel(\"Time\")\n",
|
||
|
"ax.legend()\n",
|
||
|
"plt.show();"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "markdown",
|
||
|
"metadata": {},
|
||
|
"source": [
|
||
|
"### Interpolation: Delta function"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "code",
|
||
|
"execution_count": 6,
|
||
|
"metadata": {},
|
||
|
"outputs": [
|
||
|
{
|
||
|
"name": "stdout",
|
||
|
"output_type": "stream",
|
||
|
"text": [
|
||
|
"P Raw S1 S2\n",
|
||
|
"a 1.0 1.0 0.9130434782608696\n",
|
||
|
"T 0.4666666666666667 0.4666666666666667 0.4782608695652174\n",
|
||
|
"[0. 0.03333333 0.06666667 0.1 0.13333333 0.16666667\n",
|
||
|
" 0.2 0.23333333 0.26666667 0.3 0.33333333 0.36666667\n",
|
||
|
" 0.4 0.43333333 0.46666667 0.5 0.53333333 0.56666667\n",
|
||
|
" 0.6 0.63333333 0.66666667 0.7 0.73333333 0.76666667\n",
|
||
|
" 0.8 0.83333333 0.86666667 0.9 0.93333333 0.96666667] [0. 0.01449275 0.02898551 0.04347826 0.05797101 0.07246377\n",
|
||
|
" 0.08695652 0.10144928 0.11594203 0.13043478 0.14492754 0.15942029\n",
|
||
|
" 0.17391304 0.1884058 0.20289855 0.2173913 0.23188406 0.24637681\n",
|
||
|
" 0.26086957 0.27536232 0.28985507 0.30434783 0.31884058 0.33333333\n",
|
||
|
" 0.34782609 0.36231884 0.37681159 0.39130435 0.4057971 0.42028986\n",
|
||
|
" 0.43478261 0.44927536 0.46376812 0.47826087 0.49275362 0.50724638\n",
|
||
|
" 0.52173913 0.53623188 0.55072464 0.56521739 0.57971014 0.5942029\n",
|
||
|
" 0.60869565 0.62318841 0.63768116 0.65217391 0.66666667 0.68115942\n",
|
||
|
" 0.69565217 0.71014493 0.72463768 0.73913043 0.75362319 0.76811594\n",
|
||
|
" 0.7826087 0.79710145 0.8115942 0.82608696 0.84057971 0.85507246\n",
|
||
|
" 0.86956522 0.88405797 0.89855072 0.91304348 0.92753623 0.94202899\n",
|
||
|
" 0.95652174 0.97101449 0.98550725]\n"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"data": {
|
||
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA7gAAAHgCAYAAACPe8HoAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeXhU5f3//+eZyb4vhAAJEMCAgawk7Ms3gCi4gKC41A2o4gf7kYq1Sm1VrLVuqK3a1rog+ik/l4q4axGBCghCkH3fAoQlZN/3Ob8/IiMhCyTMZAh5Pa6Ly8w597nv18RchPfc59y3YZomIiIiIiIiIm2dxdUBRERERERERBxBBa6IiIiIiIhcFFTgioiIiIiIyEVBBa6IiIiIiIhcFFTgioiIiIiIyEVBBa6IiIiIiIhcFNxcHaC5OnToYEZFRbk6hoiIiIiIiDjBhg0bsk3TDGvJtW2uwI2KiiItLc3VMURERERERMQJDMM41NJrdYuyiIiIiIiIXBRU4IqIiIiIiMhFQQWuiIiIiIiIXBTa3DO4IiIiIiLSOqqqqsjIyKC8vNzVUeQi5OXlRWRkJO7u7g7rUwWuiIiIiIg0KCMjA39/f6KiojAMw9Vx5CJimiY5OTlkZGTQo0cPh/WrW5RFRERERKRB5eXlhIaGqrgVhzMMg9DQUIffHaACV0REREREGqXiVpzFGT9bKnBFREREROSCZbVaSUxMJDY2lmuuuYb8/Hynj/nqq6/yzjvvOLzf1NRU0tLSHN6v/EwFroiIiIiIXLC8vb3ZtGkT27ZtIyQkhL/97W9OH/N//ud/uP32250+jjieClwREREREXGYDYfy+NvyfWw4lOfwvocMGcLRo0cBKC4uZsyYMfTv35+4uDg++eQTAJ599lleeuklAGbPns3o0aMB+Pbbb7n11lvr9Tlnzhz69u1LfHw8DzzwAABz585l3rx5AKxfv574+HiGDBnCb3/7W2JjYwFYsGABkydPZty4cURHR/Pggw/a+5w5cyYpKSn069ePxx57zOHfB2mcVlEWEREREZGzevyz7ew4Vthkm6LyKnadKMJmgsWASzv54+/V+BYwfbsE8Ng1/c5p/JqaGr799lt++ctfArVbzCxevJiAgACys7MZPHgwEyZMYOTIkTz//PPMmjWLtLQ0KioqqKqqYtWqVYwYMaJOn7m5uSxevJhdu3ZhGEaDtz9PmzaN1157jaFDhzJnzpw65zZt2sTGjRvx9PSkT58+3HvvvXTt2pUnn3ySkJAQampqGDNmDFu2bCE+Pv6c3qecH83gioiIiIiIQxSWV2Mza7+2mbWvz1dZWRmJiYmEhoaSm5vL2LFjgdptZh5++GHi4+O57LLLOHr0KJmZmSQnJ7NhwwaKiorw9PRkyJAhpKWlsXLlynoFbkBAAF5eXtx555189NFH+Pj41Dmfn59PUVERQ4cOBeAXv/hFnfNjxowhMDAQLy8v+vbty6FDhwD44IMP6N+/P0lJSWzfvp0dO3ac9/dBzo1mcEVERERE5KzOZaZ1w6E8bnljLVXVNtzdLPz1piSSuwef17innsEtKCjg6quv5m9/+xuzZs1i4cKFZGVlsWHDBtzd3YmKiqK8vNz+9VtvvcXQoUOJj49n+fLl7N+/n5iYmDp9u7m5sW7dOr799lvee+89XnnlFZYtW2Y/b5pmk9k8PT3tX1utVqqrqzl48CDz5s1j/fr1BAcHM3XqVIdvhSON0wyuiIiIiIg4RHL3YBbeOZj7L+/DwjsHn3dxe7rAwEBeeukl5s2bR1VVFQUFBXTs2BF3d3eWL19unz0FGDlyJPPmzWPkyJGMGDGCV199lcTExHrb0hQXF1NQUMCVV17JX/7yFzZt2lTnfHBwMP7+/qxduxaA995776w5CwsL8fX1JTAwkMzMTL766isHvHs5V5rBFRERERERh0nuHuzQwvZ0SUlJJCQk8N5773HLLbdwzTXXkJKSQmJiIpdeeqm93YgRI3jyyScZMmQIvr6+eHl51bs9GaCoqIiJEydSXl6OaZq8+OKL9dq8+eab3HXXXfj6+pKamkpgYGCTGRMSEkhKSqJfv3707NmTYcOGnf8bl3NmnG3a/UKTkpJiau8oERERERHn27lzZ73betub4uJi/Pz8AHj66ac5fvw4f/3rX12c6uLR0M+YYRgbTNNMaUl/msEVERERERFpxBdffMFTTz1FdXU13bt3Z8GCBa6OJE1QgSsiInIOTvz5z1Ts3NXq43rGXEqnhx9u9XFFRKTWjTfeyI033ujqGHKOtMiUiIiIiIiIXBQ0gysiInIONIsqIiJy4dMMroiIiIiIiFwUNIMrIiJyDlZ+sIfsI8XNvq5DVz9G3NDbCYlERETkTE6bwTUMY75hGCcNw9jWyHnDMIyXDMPYZxjGFsMw+jsri4iIiIiItE2GYXDbbbfZX1dXVxMWFsbVV1/d5HW7d+8mNTWVxMREYmJimDFjBgBpaWnMmjXL4Tnnzp3LvHnzWnTt22+/TXR0NNHR0bz99tuN9h8REUFiYiKJiYl8+eWX5xP3ouXMGdwFwCvAO42cHw9E//RnEPCPn/4rIiJywfEZFMaxMAuDe4aS3D3Y1XGcasOhPNYeyGkX71VELny+vr5s27aNsrIyvL29+eabb4iIiDjrdbNmzWL27NlMnDgRgK1btwKQkpJCSkqLtlg9b6mpqSxYsICoqCj7sdzcXB5//HHS0tIwDIPk5GQmTJhAcHD9v39nz57NAw880IqJ2x6nzeCapvkdkNtEk4nAO2attUCQYRidnZVHRESkpTYcyuOx1//Ny0u2ccsba9lwKM/VkZxmw6E87nx9BR9/s+yif68i0naMHz+eL774AoB3332Xm2++GQCbzUZ0dDRZWVn215dccgnZ2dkcP36cyMhIex9xcXEArFixwj77m5WVxdixY+nfvz9333033bt3Jzs7m/T0dGJiYrjrrrvo168fl19+OWVlZQC8/vrrDBgwgISEBK677jpKS0vP67395z//YezYsYSEhBAcHMzYsWP5+uuvz6vP9syVi0xFAEdOe53x07F6DMOYYRhGmmEYaad+eEVERFrL2l1HWGR9mJnWT6iqtrH2QI6rIznN2gM5/Ir3Wez+CNXV1Rf1exWRZvpqDrx1lWP/fDXnnIa+6aabeO+99ygvL2fLli0MGlR746fFYuHWW29l4cKFACxdupSEhAQ6dOjA7NmzGT16NOPHj+fFF18kPz+/Xr+PP/44o0eP5scff2TSpEkcPnzYfm7v3r386le/Yvv27QQFBbFo0SIAJk+ezPr169m8eTMxMTG8+eab5/VtPXr0KF27drW/joyM5OjRow22feWVV4iPj2f69Onk5ekDyIa4ssA1GjhmNtTQNM3XTNNMMU0zJSwszMmxRERE6uoXVImnUcU4yzrc3WpvU75YDY4KZrx1PX5GOSHW8ov6vYpI2xEfH096ejrvvvsuV155ZZ1z06dP5513ap+KnD9/PtOmTQNg2rRp7Ny5kylTprBixQoGDx5MRUVFnWtXrVrFTTfdBMC4cePq3Bbco0cPEhMTAUhOTiY9PR2Abdu2MWLECOLi4li4cCHbt29vMvtbb71lf242LS2NK6+8ksTERCZNmgSAadYvgQyjfqk0c+ZM9u/fz6ZNm+jcuTO/+c1vmhy3vXLlKsoZQNfTXkcCx1yURUREpFGRXuUA9LFk8I/Lgy/q51KDivcSYWQD8MyVERf1exWRZhr/tEuHnzBhAg888AArVqwgJ+fnu0u6du1KeHg4y5Yt44cffrDP5gJ06dKF6dOnM336dGJjY9m2re76tw0Vl6d4enrav7ZarfZblKdOncrHH39MQkICCxYsYMWKFU3mnjZtmr3obugZ3MjIyDp9ZGRkkJqaWq+f8PBw+9d33XXXWRfZaq9cOYP7KXD7T6spDwYKTNM87sI8IiIiDSrLP2n/2vvgEhcmcb6T6z+yf93Z/fyeKxMRcaTp06fz6KOP2p+lPd2dd97Jrbfeyg033IDVagXg66+/pqqqCoATJ06Qk5NTb3Gq4cO
|
||
|
"text/plain": [
|
||
|
"<Figure size 1152x576 with 1 Axes>"
|
||
|
]
|
||
|
},
|
||
|
"metadata": {
|
||
|
"needs_background": "light"
|
||
|
},
|
||
|
"output_type": "display_data"
|
||
|
}
|
||
|
],
|
||
|
"source": [
|
||
|
"## 1d signal: Delta peak with time interval\n",
|
||
|
"sample_rate = 3e1 # Hz\n",
|
||
|
"interp_sample_rate = 2.3*sample_rate\n",
|
||
|
"time = np.arange(0, 1, 1/sample_rate) #s\n",
|
||
|
"time2 = np.arange(0, 1, 1/(interp_sample_rate)) #s\n",
|
||
|
"\n",
|
||
|
"signal = np.zeros(len(time))\n",
|
||
|
"signal[len(signal)//2-1] = 1\n",
|
||
|
"signal[len(signal)//2] = 0.75\n",
|
||
|
"signal[len(signal)//2+1] = 1\n",
|
||
|
"\n",
|
||
|
"mysignal = TravelSignal(signal, sample_rate+20)\n",
|
||
|
"time_signal, time_offsets = mysignal._translate(time)\n",
|
||
|
"time2_signal, time2_offsets = mysignal._translate(time2)\n",
|
||
|
"\n",
|
||
|
"raw_max_idx = np.argmax(signal)\n",
|
||
|
"my_max_idx = np.argmax(mysignal(time))\n",
|
||
|
"my2_max_idx = np.argmax(mysignal(time2))\n",
|
||
|
"\n",
|
||
|
"\n",
|
||
|
"fig, ax = plt.subplots(1,1,sharex=True, figsize=(16,8))\n",
|
||
|
"ax.plot(time, mysignal.raw, '.-', label=\"Raw signal\")\n",
|
||
|
"ax.plot(time, time_signal, '-', label=\"MySignal +0.5\")\n",
|
||
|
"ax.plot(time2, time2_signal, '--', label=\"MySignal(interp) +0.25\")\n",
|
||
|
"ax.plot(np.array([0, 1/sample_rate]) + time[raw_max_idx], [1.02, 1.02], label=\"Orig Samplerate\")\n",
|
||
|
"ax.plot(np.array([0, 1/interp_sample_rate]) + time[raw_max_idx], [1.01, 1.01], label=\"Interp Samplerate\")\n",
|
||
|
"ax.set_ylabel(\"Amplitude\")\n",
|
||
|
"ax.set_xlabel(\"Time\")\n",
|
||
|
"ax.set_xlim(0, 1)\n",
|
||
|
"ax.legend()\n",
|
||
|
"\n",
|
||
|
"\n",
|
||
|
"print(\"P\", \"Raw\", \"S1\", \"S2\")\n",
|
||
|
"print(\"a\", signal[raw_max_idx], mysignal(time)[my_max_idx], mysignal(time2)[my2_max_idx])\n",
|
||
|
"print(\"T\", time[raw_max_idx], time[my_max_idx], time2[my2_max_idx])\n",
|
||
|
"\n",
|
||
|
"print(time_offsets, time2_offsets)\n",
|
||
|
"plt.show();"
|
||
|
]
|
||
|
},
|
||
|
{
|
||
|
"cell_type": "code",
|
||
|
"execution_count": null,
|
||
|
"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.7.6"
|
||
|
}
|
||
|
},
|
||
|
"nbformat": 4,
|
||
|
"nbformat_minor": 4
|
||
|
}
|