Skip to content

Instantly share code, notes, and snippets.

@campagnola
Created June 21, 2019 21:40
Show Gist options
  • Save campagnola/28c5fafd6a7adf7da79b0e5df1f19f09 to your computer and use it in GitHub Desktop.
Save campagnola/28c5fafd6a7adf7da79b0e5df1f19f09 to your computer and use it in GitHub Desktop.
Very simple, artifact-free 60Hz hum removal (using scipy minimize)
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"from scipy.optimize import minimize\n",
"from matplotlib import pyplot as plt"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"def line_noise(x, t):\n",
" return x[0] * np.sin(x[1] + (60.0 * 2 * np.pi * t))\n",
"\n",
"def err(x, t, v):\n",
" return (v - line_noise(x, t)).std()"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAD8CAYAAABdCyJkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztnXmcFMX5/z8PLCw37HKJXMupnAKu4A0IgiICEaOAGqLibfwa44FigopR4pHwi4mJBIlo4gGKAUXFFUFALhfkkvuGlXMXOQSWY5/fH9OjvTPdM91d1dU9s/V+vea1Mz3VVbXT3fVUPfUcxMzQaDQajcYr5YLugEaj0WhSGy1INBqNRiOEFiQajUajEUILEo1Go9EIoQWJRqPRaITQgkSj0Wg0QmhBotFoNBohtCDRaDQajRBakGg0Go1GiIygO6CCOnXqcE5OTtDd0Gg0mpRi6dKlB5i5brJyZUKQ5OTkID8/P+huaDQaTUpBRNudlNOqLY1Go9EIoQWJRqPRaITQgkSj0Wg0QmhBotFoNBohtCDRaDQajRBakGg0Go1GCC1INBqNRiOEFiQpRkkJY3L+Tpw6UxJ0VzQajQaAFiQpx9RvC/Do+ysxfu6WoLui0Wg0ALQgSTl+OHYSAFB49GTAPdFoNJoIWpCkKERB90Cj0WgiaEGSYjAH3QONRqMpjRYkGo1GoxFCCxKNRqPRCKEFSYrBiOi29BaJRqMJC1qQpCh6s12T6hwtPo31e44E3Y2U4MWZ65AzckbQ3bBFSJAQUTYR5RHRRuNvlk254UaZjUQ03HS8IhGNJ6INRLSOiAYbxx8iojVEtJKIZhFRU9M5Z4houfGaLtJ/Tfozf+MBfLZ6Dw4dP4Ut+48G3R2NiV+9vhh9x80Nuhspwd9nbw66CwkRzZA4EsAsZh5LRCONz4+ZCxBRNoDRAHIBMIClRDSdmQ8CGAVgHzO3JqJyALKN074FkMvMx4joHgAvALjR+O44M3cS7LemjHDz64sBADm1q2Bb4TFsG3tNwD3SRFm244egu6CRhKhqayCAScb7SQAGWZTpCyCPmYsM4ZEH4Crju9sAPA8AzFzCzAeM97OZ+ZhRZhGARoL9TBu0+a83thUeS15Ig1fnbEKvl+cE3Q2NieMnzwTdhaSICpL6zLzbeL8HQH2LMg0B7DR93gWgIRHVMj6PIaJlRDSFiKzOvx3Ap6bPlYgon4gWEZGV4AIAENGdRrn8/fv3O/+PBNlZdAyHjp/yrf6oHCG9SVKmWLD5AFbu8n8G/8Jn67F5/4++t6NxzpkUmD0mFSRE9AURrbZ4DTSXY2bGz+OcEzIQWWksYOYuABYCeCmm7ZsRUYm9aDrclJlzAQwDMI6IWlhVzszjmTmXmXPr1q3roltiXPbCbFzz13m+t+O3GDl28jT+NXcLSkrCfxOXBYb9azEG/O1rZe2dTpOgoAs2HUDOyBlYXXAo6K544oOluzDqw1VBdyMpSfdImLm33XdEtJeIGjDzbiJqAGCfRbECAD1MnxsBmAOgEMAxAFON41MQWX1E6+6NyB5Kd2YuNvWnwPi7hYjmAOgMIBQ7UYVHI93cdfB4wD0R58n/rcbUZQVomFUZ/To0CLo7GsVsK/wRLetVD7obwnyxNjIkLdpSiPYNawbcG/f8bsqKoLvgCFHV1nQAUSus4QCmWZSZCaAPEWUZVl19AMw0VjAf4Wch0wvAGgAgos4AXgMwgJl/Ek5GHZnG+zoALomeEwb+u3hH0F2Qwvo9RzB1WQEA4MSp8Otnw8Dsdfv0b+WCo8Wng+6CRiKiVltjAUwmotsBbAdwAwAQUS6Au5l5BDMXEdEYAN8Y5zzDzEXG+8cAvEVE4wDsB3CrcfxFANUATDH2AnYw8wAAbQC8RkQliAjBscwcGkFSTsG2hZ/7L1HW7E5NNUBQrC44hFvf+AZDuzbB89d1CLo7ntlz6ISSdhZsOoBhExb72sbm/Ucx6G9f44gWWEoQEiTMXIjISiL2eD6AEabPEwFMtCi3HcDlFsct1WnMvABAKJ/UpdsP4qXPN/jezj/mGFo8H4XWqTP+7oss3V6E02cY3ZrXBgCM+XgNVhUcwuS7LvK1Xb84aIT231GU2pvUg/6uZg/mm20HfW/jkSkrSgmRFNivTmlEVyQagxkrdycvlCKc8XmDffA/FgIAbsxtjD9d3xGvz9/qa3t+Ex2kKMUD1+w5rGZFomLlHuujwq7sgJxzpoRx4Ggx6teo5Ev9qYIOkeITxadTV19unr35OZN7L39nqc/j527G3A3qTLVloy2yk/O/bwvwcp6/K/cfFaqzXpi5Dt2em4WvNuxHzsgZ2OGzv1LOyBnYdTB8PlFakPjEWwu3S6urpITjzBdPnDyDbQf8UaV89/3PbW054G9YkT9/vv6n9899sg6/mrjE1/b8ICpr522MmJq+uXBbgL0JNw++t9z3Nr5cZ2U86g+vfRVJef1PQ+Wct3av721+s60oeSHFaEHiEzL3GcbP24L+r8xHvukGmrRwO3q8NEdaG2bM1md+x/j565ebfK1fBSUxy7bo4OIHBT8cxx+mrcZNExb51oZqFmw6EHQXhFG5Gg3jfo8WJBJYur0IE78uredfs/uwtPq/+z5SV8EPqe+fopJ9inT+Pqnff8J83Xf/cBxvLtyOrzcV+tuoQoZNWOy7OiqMg69Xwvi/aEEigdsn5ccd+2jF90ra9ntjPJVRoSY7dOwU/re8oNQx2QL/krFf/vT+xvFqViJT8ncpvbdkOvGqWh0cOvazKb7KFcnuQ+GbUGpBIoEfjvnr25HoHj1dkh6hLPxgtwK/iAfe/RbTlquZNADqJg6vzd2Cf3+tzprujQXbpNVlNWP341frNObzn95HV4isYLmgws3ALdr8N8WRed++vXgHmtetKq/CgFHxUIdxdiiLA0dPKmvr2x3++5bIJkgV055DJ5CZUQ5ZVSsG1wkTekWi+YknPlyFIYpUJyo4fCK9vZr9FpSpas5s1e93luzAydP+r95j1Zx+ceHzs9B5TJ6StpygBUkKUXxKq7FEmS3RNPS/i7djq08m2E7wW8ulwnEwit9pEbYXHsOrc/y3EFxdIM/IJpXQgsRHZKV2jT5jj36wUkp9QaJC3ZQImT4eoz5c7Xs4GVVYpQsoPlWCvDX++0UAciP+nLCZcPm9l1mW0YLER2I9t72icjPXb8bP9c/Hoqwhc/D9x1fx/kIT5m/FHW/mY833/s+yZS5IHk6R0OvphBYkKU4YbcoT8cmqYGOSFSvQk6tC5qWflkC3f1xBePyiH/3f2JclrI4Wn8Ytr/sbvTjV0FZbPlKgIMHVnsMn0KyOv5ZWS7cfRPVKGWhdXyzR0aZ9R7Bil/8h6hdsOoCVNhnxFmxOH0c+mWzYa6+GrVje//mmClNtWUE1P/9uD+ZtTH1vfJnoFYkgexN4T6vIHRKN8eMng/+xAH3+Mle4nkSDlUyGTViMsZ+uU9KW36gMQGhHRvkUNd+KQZUVWiqkxpWNFiSCjPpwte13KmYtK3b9kLyQJmUZ/I8Ftt+pshjLkGS+NXu9umCKVqgSh+mSKdUNWpAIErQV0ro9R6TUc0hbtISSRNf3tje+sf1OJrJm8rf+2//+WlmfRVkSwqi5iUilVBRakAgSG/k1Vbn8xdm+t5EmP1Vo0HHW4vnvYvv0DSsl7M8dOn4KD0323yqMmbF0e+p4+2tBIgAzY32SFYHfiW5koWI/RxNOnv90bdBdkMb3Pm/am3P1+MnEr7dh2L9SxzJMCxIBXp+/NemN2/+VeYp6o3GKCpVQIhWLLGSpVZPlT0knk+lUIW/NnqC74AohQUJE2USUR0Qbjb9ZNuWGG2U2EtFw0/GKRDSeiDYQ0ToiGmwc/zUR7Sei5cZrRLK6guBrBwl50j3eUyqiIoNe8yc+8b0NVYz5eE3QXXCMCudJFaSaGljUj2QkgFnMPJaIRhqfHzMXIKJsAKMB5CLiQ7WUiKYz80EAowDsY+bWRFQOQLbp1PeY+X4XdWlCDvudAUohQRtZqGTRltTZpP5qw35/Gyg7l90VoqqtgQAmGe8nARhkUaYvgDxmLjIG/DwAVxnf3QbgeQBg5hJmTjbFT1SXUg6fOIXZ632+adOMPQqczjQajXpEBUl9Zo7GvNgDoL5FmYYAzEGndgFoSES1jM9jiGgZEU0hIvP5g4loJRG9T0SNE9Ul+D944rfvLndcdtEW797UPxxTlxPCb56dkT6bumHA741lMycUhEnRpC5JBQkRfUFEqy1eA83lOLLWd7PwywDQCMACZu4CYCGAl4zvPgKQw8wdEVl1TLKuImG/7ySifCLK379f/sphR5Fza6wlW72rBlR5g2vko8rCRwXLBBNPhcUnQoURRFkkqSBh5t7M3N7iNQ3AXiJqAADGX6tdzAIAjU2fGxnHCgEcAzDVOD4FQBejzUJmLjaOTwBwfpK6rPo9nplzmTm3bt26yf5N17hx0hJRpzvRxR85oU133SK6x+FEwC8UjOu1dHvyCchRVSFUBMffnS4mXn7y+2n2kShksmCzWFSLVBN3oqqt6QCillPDAUyzKDMTQB8iyjKsuvoAmGmsYD4C0MMo1wvAGuAnoRRlAIC1ieoS/B88ISsAnAxSyXEpXVizO/lqQ9RhsNBBqtv/e+dboTZUMWNlOMxZVYUvEZ1EpBqiVltjAUwmotsBbAdwAwAQUS6Au5l5BDMXEdEYAFHj/WeYOTrVegzAW0Q0DsB+ALcaxx8gogEATgMoAvBrAEhSV2gJu7WS09n57kPH0aBmZZ97I4ZTB1Bm/4P4iWpRnIRWX79XTogcv/lacIbuhNNn/Pd3UfUkhyFYpxuEBAkzFyKykog9ng9ghOnzRAATLcptB3C5xfHHATxu06ZlXaoJUz7rTfuOonvrup7TlU7J3+Wo3LWvzEf+k1d6akMVMxTlO1GxIn3cQRRZ0ftQVZiVYyf9Hxhfn7/V9zZUoSK/vEy0Z3sa8OyMtUIP0do9zpy4DjhQtQTJmRLGa3OdhdUXHT6dPOiiK1EnC0VRgTbZYRZP0d/LSQgeUWFTqCA5llM1sqiADtNE1QlakChAhe+aSMKodPGt+/DbAsd5uUU32//0WfJ8JxPnbxNqQwVHFUVeKHEwwZ69Lvx+WX/O2+Co3KsK8gSFCS1IPOJVjeSWNBnjlaAiJSwQsZRyMvs9cLQ4aRlRVOnSl+8Uy3vjZIb+2AcrhdrQBIcWJB5Zu9t5TB8VwuAbAV8VFYQtpIhIb54MUQY8FeocAHhx5nqh852o+VSZMoftXrTCqf/YRyu+97knztCCJOQMGb/IUbk9CVL+hoG3FtnniQgCkbFk2Q6dldItCgyqHO8Wvb0k3BkMv3Xh/PmbkJh/a0GiUTJDm7vBf/NPVfuTqWaaGQb8DvVTeLQYr81NHA4/yuff7fW1L6J8/0O4J4VWaEGiUaJ6Kx+yO22DgP9F+BUj4eO0z2bG012oeHyPEFwGCdnjnaZ4nPF/sNSZf0cq4CY50sXPz/KxJxGG/cuZytCKVNCxlzU+WObuWdHXUC5akISY9xUJEhXP1BwXIfdVRLUVsfNPpyEo7FEXnLK6wF1Cq3kb/Ve1eiUVr4kWJCEm30HQPo03fjzp3VQ4xXzFEvLcJ8n9YdKRVPMcDztakIQYVavvsFlUeUWVN7CKyzJjpfNQL1sP/OhjT8Q4WnwaD09ZEXQ3pKDVYfZoQaIAr7efioExb024LVjCiJucFvnb3K8qT54uwX1vL3Nc3o1PUxRmduSdH3uOW974eqsrFa2KmFyAt2crX0fZtkULkhCjIjDgthDPZsOKm+H0iIcQJG515G97CI2+70gx/uEyjMdLn4s5JTrhyj/P9b0Nr2h1mD1akCjAa3gJFZtubmdmKkJ1q8KzqsLNaQpWlfM3qdk4fmuh/yrQgh+O+95G2EnFJI5akCjAq4VImJJnRRENlREmvMiRLmPycMSFQ+IKD5MIFdfdSwspOL7ZEuboug+ExFvdDVqQhJggzQD7tK1veTzMemK3aiQvv66TZFNmxn2x0UMrClA0kKqIaKAJHi1IQsrR4tM4dSY4QVI10zrnWYgnchj7qbvN4837nQXGc8KgTmdLq0sFXlY9Xq79Eg/GBumCkxws6YIWJCGlSFESKbtw+JUqWN8aYVYJuGWPRMfHcUM6234Xxk1aL9fxsKLcJWHF6icbM7Ad6tfItCx/3tOfu26j0Cb1wL09WriuSyVakJRx7MaTJ/q1sSmfPpJE1XrvxGk1eVI0/nLMwom1b7uzcFW7s6S1YRdx4dGrzpXWhh9oQRJS7PZHxt3YSUn71StVsP5CkRz5XoH1TlgdzPZ7SIi11GUUhKCnAzdf2CTQ9r1MiEa8me9DT9IDIUFCRNlElEdEG42/WTblhhtlNhLRcNPxikQ0nog2ENE6IhpsHP8LES03XhuI6AfTOWdM300X6b/Gmu6t69p+53c48Cj7j6jILhjOlcL9LpwRo6zb4y6asaoMn3Y0r1Mt0PY1chFdkYwEMIuZWwGYZXwuBRFlAxgNoBuArgBGmwTOKAD7mLk1gLYAvgIAZv4tM3di5k4AXgEw1VTl8eh3zDxAsP+hxW6y3LR2Fdtz7vuv+wHIajx5pO85tuU37D0qTecvc4PaanXxWBJ1gBvv8UQsfPwKAMD8x3pilI1K0A2HFWzSlgt6SZIAWSvF+jUycVf35lLq0iRGVJAMBDDJeD8JwCCLMn0B5DFzETMfBJAH4Crju9sAPA8AzFzCzFa2gkMBvCPYz7ShcxPLRR8AYMYq5/GZEtG+Yc2E358ukSNI/nR9R6x6qg9mP9xDuC4rR7Z7erTAi9d3xMe/uRQP92kt3IYdDWpWBgA0yqqCvhb6cre+JF6GUbdjb9B7XYkWRK/P3yqljXt7tMSjfc/FgpFXWHRAShOlyKpiow4uA4gKkvrMHB299gCwcj5oCGCn6fMuAA2JqJbxeQwRLSOiKURU6nwiagqgGYAvTYcrEVE+ES0iIivBFRjZVSsG3QXXJHqeep5jr+JyQ87IGZbHMzPKo3qlCmhWp6pwG3YD6S9zG6N9w5q4sq31hqjsPOFNLFaMt7y+xF0l4dy6kUoiwTdtubs85MsSpKYtX45wdq3KFu3L+ZGzTM98oyx7bUG6k1SQENEXRLTa4jXQXI4jV8bN1ckA0AjAAmbuAmAhgJdiygwB8D4zm5XZTZk5F8AwAOOIyNIujojuNARO/v79/mdEq1CepMyso1j9kCv+0Eda/U7o0KiW5XEVs1m3j/mbC7clqc+6xlMhNM1VQohVW263b657dYHrem57I99V8E07KphSf1q1d1wgXUEsD13p36palKSChJl7M3N7i9c0AHuJqAEAGH/3WVRRAKCx6XMj41ghgGP4ef9jCoAuMecOQYxai5kLjL9bAMwBYGnAz8zjmTmXmXPr1pUzs07ERS3qoGZlf5e2NVUvnW1mbX7s0z56Vel9meJT7h7Af82Tow5JVVRFzVWBKkOAwyfk7kVZPS5unRKLEhizPNCrldsuKUNUtTUdQNQKaziAaRZlZgLoQ0RZxiZ7HwAzjRXMRwB6GOV6AVgTPYmIzgWQhchKJXosi4gyjfd1AFxiPkcVVsti2aakqkxTg7beiXJpyzqlPv/dZWTaZPi1iqpcobwv9bpFdoKqD+65yPK4SGZJM4k2+/0wBKhS0f/rZLXqdft4PT51laTeqEVUkIwFcCURbQTQ2/gMIsologkAwMxFAMYA+MZ4PWMcA4DHADxFRCsB3ALgd6a6hwB4l0uPqG0A5BPRCgCzAYxlZuWC5II/fmH7nQx9v0pk7xF4pWJG6Vtxi8TwJX7S7uwa0utMNFR/+bvuvjfy3xHdcH7TbMvvnpr+nZTmh3Rtglb1qlnOsr/d4S1adiL+PixW2SEfGZMVFRZ7fmAdUMkhzFyIyEoi9ng+gBGmzxMBTLQotx3A5TZ1P2VxbAGADt57LIcDCcKXTL7rImzcewR/+mwdVuw69NPxXQePudqM++ssb8H+jhafRjWbOFlWJIzmazOd2ne42HJTWYRz6lcv3YYCPxIZvD78AmVtdWpcC83r+u9/kWgx/Nai7RgzqL1Q/U8PaIdKFcoj76GIUPR6r6cCCzcXYlDnho7Lp6qdhfZsl0T04atbPRMXt6wTNwj/++ttrur7XwLLlUT+FyoSVT34nvww17EqtjDGp7LCj30rO7VmmH0/3HDzhU1LfZZppBIkFcrHX6B3lrhLOublEochR5AWJJJ4sHfMEj1mMJA5CIwb0hlT7rbWYXtJu+qWkyG4cWXhRoet6oHdVnjM8vhFLWpLa2NVwaHkhXyifMzDkGrqYDv6dWggXMfm/e4ngoUuUxv4gRYkEhhyQWPk5pTWKcd+Lid5U9tO/fDWInlZ7FItNLrfTJDkKOeFuY/0xENXRizb2jcU35cZ8eY3wnWkEqK5fZzkocmwmC2qsGXp9twsnAp4cqcFiU/ErVAUqSVkNtO8bjVsfq5f3PGQxjpMiIwHeu9heWHn3VCOIo6O0Zn8B/dcjLHXiW0VJrqG0UH39kubCbURdtzcx7Ks1fwi6P5pQeITsdFz/bBEsUKFOe933ztXn20v9H/Pxgl2g0bQoUKcEHtNMzPK4+IWdWxKOyPRsBN1svt9/7a4MbdxgpLhRrYWwI47Lm+O8xrVxIBOzjfVrdhho9JMBYSstjQRnNyvS7aqyRQXEreQn+jzl7mWx889q7rl8TATppWYbKu5KKP6tUG3Zj+rZUtsfKb8nrAcOn5KyMH3tkua4Ze5jST2KJ561SMJrRplVcG0+y8Vru/bndahXmItGsOIXpFIwdlDNXeDvFAtdpY9Pwr6hfzvvkuEzo+l2Mb6StVsMSyIhONwGvPMlae2TXfuuLx5KSFRt3p89r9vttnHtpKFl+yCZv5wbVtkZiR2QnRzRSZ+Hb8/9quLmsYdi1Npu8DOUvG9uy70XKcqtCCRQMt6zmz7veSZsMPuIdiwV8yRr1Pj0vG1ypcj9O8obo0SS1sfHPkSEbTcOiIg4P/m0Jlu3W7nOUmcWt79n8XAeLQ49Zzm2jSIv9/cRI+YvzE+MLnVqqxjo8SRsxNh59Veq8rPgSEf6XsOrj0v3ggm6PtbCxIJ3HpxjqNyIpqRs2pUEjjbGqc62WHd4rPZia6uwrIeeeyDlUraiTV5dUMlh2FY/AirYzWrn5K/S3o7ftOgZmU0r1vazNjNIvGEy9hvUbYe+NHxavS0g3L39WyJV4ZahhcMFC1IJFBOgadYq/ryPZqveHlOqc9u4hE5MYdMhJNHq1gg1/mzMd7XtWz07Z99t8dzG25Q4Uwoargz6baujso5vfZhS2XcP8bPY/4m55Ohjfu8rfT3Hi7GK19u8nRuKqEFScg4U8KYtrwg7vhzvyht7nm2kUzpinPreW4rdgZkuz8S0Hjw5883SKurXo1K+PJ33XFBjn1isGSIDIz3eshe6RZRX4lEKZZLt+OMU2fCJUhi9T+/fW+FUHV1qjnLP/TNNvmGNtd1KW0hds6Tn2HfkWDM0wEtSJTixNb7P4u24//eXR53vHF2aSudJrWrYNHjvXBfT8t0LFKx6vXHK8WyMToZk3cf8v5gWFXfvG41VK4oz1Dxxes7Oi47Z738nDhxG+Gqxm2H7VhZfAWK5P50beYs0oAfv0Mfi0RtQZoPa0GiECf304GjzoMVnlWzkuWG3z6PjnN22herfn+xdq+nNoDIxufj/RLnUwf8GRebZns3m90fc21iIxZHqV7Ju7B6/tO1jstOuat0mBwR1ZYfBhUzFagN3exdqBBrVs+KiByZ8YBzs+LjHvdxZKAFiY/0ilE7iaoerLC6Se94a6mnuuwsP2T3e+o9F6NOtXiz0rh2BZ5AN1sSTsPKfLKq9MBoZ8K86qm+6NfBOrVvMl77aovjsjkxMaqsVKJO8cMc22pl7YY9DlakW1zEpqrqIiq2V6xuWa8/bdsGNdDubOdWYEEuALUgcUnsDKh1gk3wV28ubbbp5ELLuBlW7PTqRW99x8u8QS9rVQeVHW7qizTr5twXPvOWFCpok8tYpiz1bk3lxhjAjwmRFbe9kTweWGYF50NYg5pyLR+d/mZe75PEIWrCpTbUgsQlsaqnRAlzYk0n/ZgxOPVhESGoW9bpiiSoMNpWeupUJSyZMs04SVPb6+WvhNpY4yLcj5mRV5+LJg7VpMc85m1PdEmsVpBBXkItSFzywDulc3HExtRKhJNNN7c3g8w88baqLY8SUNT8s/iUMwHxxIdi6Um9Pn92eySROoN5qr160Lu578K2hy7CsAmLPJ13d/cWjoXvroPHk5b5wGI12bttfdvysdHFg0YLEpcsiwm+aBVCwg4nz5+qh3TBJgtPXcltzHHptBjrtDdr3T5H532wzPvegF94Uf842RNIRsEPyQet1Ra5SGQLvh+OBZ8jwwle9obcGlM4acEqFl+NBJPU7KrOTI9VoYM2ClChPLnyWA6Tg9awCYvjjtWz8Z732uuDLp0WrX7JLfuPekovW9fBZn6Uwyfk5633cql/+57Y5rRT+r8yP+6Ym/E0f3vyWFtPf7Qm7ljeby8P1LfEahYve/LUwSJEyr4jxUkDXcowEQ4ykrVekQjg9sKFR4xYU9kmFIedAPQiGHNq22fDG9o1PhTLrLXJVyVW/jl921mrBVTpkZ0695mRkZzIa5Rp2Z73x07GC+dW9au7irHm5VqNu7ETPrjnYsvvGtaqbNFG8kbchAOqX6MSto29Ju54Mr+rPRJy3agygrBCCxIBzrgcSP1akFhlZvOC3erKrt/JAkTuOxLvE3Pn5c1tyz89oB2qxZho/vET534VZoLePLYKEpgMGV3+YJk3yy2Ze21BMqhzQ5zf1Hn0Aie/uYjPVJStBxKbKW/Y6zzgZhgRFiRElE1EeUS00fhreRWJaLhRZiMRDTeOVSei5abXASIaZ3yXSUTvEdEmIlpMRDmmuh43jq8nor6i/4NXgs5KFmVI1/jkQzLVaHZVJfv/x34ab1abmWCDulw5wuAuYsmBkuF1rE42EIQFryoSq9m6CNtdelnLmgy5xUmzh2Osx7z0NNllOeHQsESkDT+RsSIZCWAWM7cCMMv4XAqcoAsPAAAgAElEQVQiygYwGkA3AF0BjCaiLGY+wsydoi8A2wFMNU67HcBBZm4J4C8A/mTU1RbAEADtAFwF4FUich5tsIwQaxQgQjsJOcKj2O3DRJGxkqjhwbP80LHEpqaxg4kfeNFx147ZdPU6t6kheUVyxOW+k5UATHYr7D6U3LBABmtjwvNHM0jKRMbEL8hprYxfZCCAScb7SQAGWZTpCyCPmYuY+SCAPESEwE8QUWsA9QDMs6j3fQC9KDLKDATwLjMXM/NWAJsQEU5lFqsByGvYaysa1KyM+Y/1jDu+N8AgcYlonSCjnN2Aed4zYomUpOBBhsoYPM49qzoGCaaJjcWJ9ZgZq+uSTLBe9PyXrtqwYu/h5CGJ1seonZw61JpZuOWA1MR2AHBpy9LploM05pEhSOozc3QnaQ8Aq13OhgB2mj7vMo6ZGQLgPf751/jpHGY+DeAQgNoO6wIR3UlE+USUv3+//IB5XrHLghbFy4bZvQoCN1rlpfhslbtYSv/Xy1v2uC/Xieuoo9zXs6Wn81RsuVg18fFv3KVwdWspB0SMHFSkQkjEVJsN8jDy2FXJ48TFsmhLEX41cYnUfsTek6/Pj8/iqApHgoSIviCi1RavgeZyhhDwKhaHAHjH47lxMPN4Zs5l5ty6dd1b0FghmsYW8GdAalDTwhpFchtWs50jLjPltfDohX/bG/lYt8e5B3Ki37hShfJolOV+P8Ctd7LVQzA1yUZ4bNbCFwZ3RPuGiWMtxVraecmb0bddYg/99hJVm3Y0r1sN13SQHzjSD2o7DB+vGu+hkcRxJEiYuTczt7d4TQOwl4gaAIDx18peswCAeUe4kXEMxnnnAchg5qVW5xBRBoCaAAqT1eUn7UbPLPW5dxv3uUCSDfBWy/l7evi/4kiG1cAYG8QwGSLCzc3mbbIVflcPXsFDxnvzgDbz0OTE+S9iPaBrVFbj5nVWkhhUH//mMuE2Fj3eK3mhmBtkR9Ex6Y6NH93vboWnChlaKT/8oZwiQ7U1HcBw4/1wANMsyswE0IeIsgyrrj7GsShDEb8aMdd7PYAvjRXPdABDDKuuZgBaAZC7ZnTIAA96ZS/3i1eVkBuSrZRk3OhOvIjtw7SItx/l+cEdkhcSRI6+OvnvFSYn10QkE1aA9X87zyJXughWDoNuaVXPfg9OJWG69DKmPGMBTCai2xGxuroBAIgoF8DdzDyCmYuIaAyAaDjPZ5jZ7Dl1A4B+MfW+DuAtItoEoAgR1ReY+TsimgxgDYDTAO5j5kAC8fvxEAflVJRsyEoUV8pxGz6o9bzElsrMKA8isQfRaXY8jXOsLPZCNFb+hJuwSGUFYUHCzIUA4tatzJwPYITp80QAE23qiPNSY+YTAH5pU/6PAP7oscvSUJUBLgyBWVXF9rG31LHxrvfYTjki1w6lZj681yYtsYEXM+b9Fg6cyXD7H/w5T176Yo08QpdN0iXasz0Jp86U4PGpKy3NGbOquB9ck90v32xLHsfIES7HsViPcj9wYhgkQ2g6eSRFjZRiUx/H9UHRwOC2mb/O2uhPRyQgeul/dVFTKf0Igh89hpoPC1qQJGHB5kK8s2QnRn6wMu67WDtuJ9zy+mIcPmFv7eQ1VlIsL3++wdVg9tehnZOWmf1wD4EeASJDhd2/8tLn6731JAzLvCQ46WKQ8ZVkY/X/urmHnYZGUWGFJoqTYLD9QmTlpgWJABkePFwXby3CRyu+d9dOOfftLN1+EMU2PitWwQGd6H2b1YkPuOgmu6CTgdEuHardcPKPOZvjjiVySIxS3mdBkj7Duzqsrohoul4rbrkw/CuXrx+7ImmZoRahkYJCC5IAOO5iGbtt7DWuQtWbsZvM/cdhjnInvGoxkNvh5L+w+1et/pd5G60dTZ8a0DZpO0HFdkoXEu21LNpS6KlON/nJrXC6eHETikZmhAg3OLJyC9GqWguSAIiN3eMXdht4XlN/iiJy48+y8G7f/UN8iJZ5j/a09MKPJWhPbic4WSUGtUebaK9lcv5O2+8SMfj8Rl67AwA45yyHZrkuLv2Yj+Pzqmji0YIkAFRNJNxYJakYkBxttts85VMtsiBa/Y7JNsGjeF3lqaRFneSRAKwum+yYTqmClerVCjdX/r+Ld5T6/OXvurs4u+ygBUkA7Ciy9tIWyUkw3MJihcUjU0slS7IJscgKp0sT5zkrZPL8J2stvbWnu9w3i9LBIoSKlzAp6UCiFAVmnDjGHjp+Cq9YrLqcZut05MmfRmhBEgBWllklJYw+f5nruc7hF+fEHXthpvON8BYe0tm6RfbgvW638/hbsVzasrbEnjjntblb8OyM+GRdD7zzraf6rKzt3Fg6ecnkGFacTixibVesjE+e/ug7vCzgc1O/RtlyWtSCxGcSZQQ0s3yX/IBrscvyKLHP2x/6t/UUGtsPKlVwfktOEIh26nSoPXWmBH09CPhEY7mMlLpRRP1//q+3/PA7O21W3MlQ5XtTr3rpjWyrGGjHisX2EcO0Ea4CLUh85ol+bRyV8xLqw4zTG/e+t5fhhc9K+16E6Z63Wln5gdMxa/+R4rh8FKJEf25mxjYfMi9arXjs8EPF59WpVpXdQOymvFtzfL9JxdWMFiRpgtPZ3IyVu5MXCpBKFdSsjPwOSeFEX/+fxTvQ46U5WLrd2gnVg/uQUnJGzojb1xPxnHd6SU6cOoM3vi69Gn24T2vP7QLApAXbSvdFoSeQ3fVPJUJ+qwZPdIDeuFfuBua+w6VNV/ccDme2wVjmPRqfKdEKN74yQeB00PI6nLQ7uwbGDGqPyXddFPdddPX47Y7IzH3rgXhVUPVKGaheSW76Wz+I3dcTieVlNxmKjQTx0sz1eOqj0ma5F7VwH2XCzOjp3wmdL8LgfywMrG1ZaEHikNiB3m3muli6PjcL+0ypaq3MW91QR1FE0oa1nCWFOnDUfQBClVQoH6/Pmzh/K3YdLD2oe9XbExFuubApWnlM5nV1+8TJptIRO3P12P2KvRbBLR37kEB+0rcg6XlOaWOJPYeCmZBqQeIRN5vCK/7Qx/J4kYe0qHbUEJi9unmwwrSfIkL9GvGew898vCYuHaqoBszK5Hl1wSF8tvrnpGBWwsqN97UVZwT33ILArstO1EwqNurvU5DS2i2xquC1ApaMImhBooCaVawHefO9H2SCIjde3kFao5w6UyJslJCMw8dLZ5nz47Js3HcUd/9naUJhYRdzzCkLN3sLU2JH/pO9PZ+74dmrHZWzu7ayr0HqiVjnLNkWzH6LFiSeSZOpOYAbcsWCv8UOAKM+XIW+40rrzv8+rIvj+q441zqFcatRn+I379r7W6x+uq/jNuzhmE/BDDs1K4vtj5wqkeuNWqead9Wp0ygCKZ6SIxBi53WbAnJG1YJEI2wpFavb/u/iHXHxvK7p6DzkdaJhJ5HVmRufCruF1YGjJ7HRZIl04pS3SMlBczIm8vMjU1bgwQRC2E+cTrlqVXUmPFWY66aKULs4xsggb018TDoVaEESILJv1jED21ke/3rTAXz3/SG5jfnI9YLB+0T5etPPecJve+ObuO9V2PmLahBj90imLN2F/y0Pxl/C6f8iss/nhmTdyd9WhM8tBmQZz2vus19IdUh1moPFb7QgKQPcNGExrvnr/KC74ZizHVqGiWE/nJjHC6vMmKlAqqZuHdq1SWBtby+MOIc+9ZG1KXDs7N8LB44WS7VobJgV/6zkjJyhfM9VC5I0IlORM18ssu9ZJ6bM7y6xDv/ilB7nhCfGlOjP5zTqbVC4MdD446D2cce++z65JZKMW/DlzyM+MHb3cz1JK9FiC3WpV+xWcaqN9oQECRFlE1EeEW00/lqus4houFFmIxENN45VJ6LlptcBIhpnfPcQEa0hopVENIuImprqOmM6Z7pI/0UIoxns5a3cD47nurC/t+OyF77E6gJ5qjMnviojp64SaiPRvpCqyVyie8jN7ZUCEfEdY2VBeMeb+VLbsLu8qpJYmTOXBmWuKxvRFclIALOYuRWAWcbnUhBRNoDRALoB6ApgNBFlMfMRZu4UfQHYDmCqcdq3AHKZuSOA9wG8YKryuOm8AYL9Dw0yxq5sC58Fs77f6Tlu2Xu4GP1fmY/TEnW/YeasGnJVb1OX7ZJaHwDc/3YwG+t+UlLCvpp/q5rEmy0Br/5/83xpQ/XcQlSQDAQwyXg/CcAgizJ9AeQxcxEzHwSQB+AqcwEiag2gHoB5AMDMs5k56mK8CECwu68WyLhQ0RtqZ9ExzFnvTzKimyYsTvj96TPyHp/TitbTVvk8/MBOz/zyDedJbWfRlnjb/0wXDq9e+eh+segMqmn+xCe45hXxvT67Z1f2SvTiFtapClJ0+yohondrfWaO2mPuAVDfokxDAObcm7uMY2aGAHiPrZ/c2wF8avpciYjyiWgREVkJLiXIuBc+WBoJizLgb3I2wr2o237fP3l+87ChKvzKil3W6jpRHw8nuImCbLcH0f3F2Rg9bTW27Lf2LWjfsIaXrtni5wZv8emI2slOFSSj6WT9F3USjfLqnM1S6kmEalmVVJAQ0RdEtNriNdBczhACXvs/BMA7Fm3fDCAXwIumw02ZORfAMADjiMgybgER3WkInPz9+8OZenSiEcH04LFTSUr6R5sG4nskUVTNtGLD4Msm+m+csXDqu+OyZr62HcVJ3vlkbC88hkkLt+OKl7+S0KPkTPPRvPhtm9w6MokNDhmL0zhzUeyeh7CFrZdBUhHLzLaxEYhoLxE1YObdRNQAwD6LYgUAepg+NwIwx1THeQAymHlpTN29AYwC0J2Zf5qCMnOB8XcLEc0B0BlAnIhn5vEAxgNAbm6u9CEuXfY3M8r7r0JpXV9u9kUrG38v1K5aEYU28c5yRs5AS4uAixfkZEtpWyZe70VZ4W7u/e9S/G1oF19NpZ+OifbrB17zqNiRqibYXhAdRaYDGG68Hw5gmkWZmQD6EFGWYdXVxzgWZShiViNE1BnAawAGMPM+0/EsIso03tcBcAkA3+6wE6fO4Nf/jndIk8m8jfarJbfPedDCjcFSrbf8pmsza6EQVXFYhZvoaRO+JRF2jqLpwier9uDgsZPYnCK54mtWroDG2farC1njfxjEyL7DJ5SogkUFyVgAVxLRRgC9jc8golwimgAAzFwEYAyAb4zXM8axKDcgXq31IoBqAKbEmPm2AZBPRCsAzAYwlpl9EyQ/JFA5yZrNPfb+Svs2pLTwM4eOy1GhJfLs7m+xGSoaydYvHr3qXNfnVPCwgqtnEWkYkHd9f9EldstRPQxg6rdiqRCiXNbKg+Ofi1E7o3w5zHv0CsvvVuyUmPI6QZ/8dhhcu/swnv9kLbo+Nwu5z37ha1uAA9VWIpi5EEAvi+P5AEaYPk8EMNGmjrik5nbqNGZeAKCD1/66JVHAvooOMuA54fsE+QNkR9qV9ZBMv/9SdHtuVtzxuRusTY3D6HMD2F/DoFUSl7Z0N5De070FerSuh35/9ceUVDWvD78ArZ/8NHlBHxj496/RtoEcI4TqleyH17eX7ED/Dmdbfte6vvi+pdWEzk+0Z7sHHujVyvXG27Bu7kM/1JJsHRSba8MrVrk8AODu/yy1PO5F6Hr5vdxSyaZfKqxqElGnmjvfHiJSYi6cCJmyN8ODh2XVTHlRHdZIchJ88Zf2ZuKTv9mJ85753PK7v9/kPFJ2WJBjz1bG8JK9zu2z8eyg9q5npkHmCkmEl8BysoWoFbVtQqMnUml6we0gO+KyuEW6JgH/ub2bEqMRtyRy9k1kWq4qeKVMtCDxgJfZV7uza7oqf/OFTZMXShGe6NfG9Tlh2KiURbuz3alKUiFMfSz7LdLfAt6CMLo2Mgnn/AlAJLXB0eLTyQsahPl/SUT4xHia8ovOwW+IBoWXDeo70mhW3ji7iqvyqTiW2O3PqBgYnSbOEsGTAYAHUvHaA1qQJESm3rdcqk41AkJGDLCwM2WpdYwtOyuvoFkw0trSKREq7noVz1Z/F4nZzLjtmRf19B9/ER8xWTVakHjAS/rVdIrQqimb1FCwbwW4H0xVzNGURYT2cE6z2sGnEdCCRBEqZk1aVmnCRroM8mHes2sagnw0WpB4wMuNq+KBssrloJGHKj25F8oruMG8mOVq/Kd2CNTAWpAogohc+55owkWY97ma1na3oe+FShXK4+07urk6J6xRDZTh2gLN/e8VhpBeWpAkQPb1+eCeiyXXKMbIq92HCCnLhOB5tUWVD5GMvOWpiKrBukL51BS8WpCEkCoV1eRev7u7ZQT+lOSskFo6RRlxqbPw816cN8NMiBdxrvCq1kuTfz8pWpAkwC5Wjt+zkw4N3TkvaoC5j/YMugtSuKqd+6gJGv+5yCbboSaCFiQJqGjjSJco+m0iqiUI4qYRQ1YQzaDxYlquCS//vvWCoLughPR4+nzCalm+4g99PDuMVXOYqlNEHVDHJn5UUIgktfrkgcsk9kQcv0N/R9rwvQmlqFDtnNfY/xW814jQ5zfNxoO9WzkuX1mRWls2WpC4pGYV/52yRCxdZj6oZvD93ZWtHZX78N5LPLfR1mGMqiwF1wRwH97djNPJQZrJESVGAF7TEjvdtwKAEoEL0yjLmUXdyKvPxaRbu7quPwyrWK1rSUBQposiz55dRFvZOG2nqsNVmAg1FXlc33m59/hflSo4G+zSbUUSZtq4yDsikqPG6Uo2lY1f9IokhKSLpYsqVI29IrNrp6Hhz2sUfkMLN46ZdumMkxG2PS+namkrysLkIFxXSwMgNZy4wrCcTiWcrpouFlCfqeIehzPnpU/2Rr8O3oIdOvXUv6u7mijRdsncwkAlj6o9mWhBkoCgVgapsCIJk7ALT0+CZcXoPpj/mP9m0E6FnQo1a5uz5KTFTWXKlSPXEQek9yHQ1jWWqLC8Eg2blArCzglXtq0fdBekUbNyBccbuyKqGhVc2NyZSiwVVsYi+ytOCXpipwVJCBkzyP/8Ai9cb59P2glh0vuKdOWfN5+PC3LSy5s8GW+P6Ia8hy4PuhsJ+ftNXdC5Sa2k5a5u7011BqSXhZyK5F6JEBIkRJRNRHlEtNH4a/lEEtFwo8xGIhpuHKtORMtNrwNENM747tdEtN/03YhEdaUbKmaLqbCpq4Ly5Qi/uign6G4o5cLmtdGgZrgDiFapmIFW9ZL7IDm1hguSSxTse1UO+HcQXZGMBDCLmVsBmGV8LgURZQMYDaAbgK4ARhNRFjMfYeZO0ReA7QCmmk59z/T9hER1Cf4PtqSq9saJmWqr+tWF2siuqsbk1gkXNtPhK9yQKmpJv1e9Tn+Gz38rtnpTMTGsXzNYR2RRQTIQwCTj/SQAgyzK9AWQx8xFzHwQQB6Aq8wFiKg1gHoArBM/u6hLo4bqlcIjSLyamEZRNbDWqRZ83ghAXaRgUfxWPTmtv7XgpEsF9aoHa1UmKkjqM/Nu4/0eAFY7lw0B7DR93mUcMzMEkRWI+doOJqKVRPQ+ETV2URcAgIjuJKJ8Isrfv3+/w38neGY+GG7ddRgRHXCcmJoO7tJIsBUg/8krhevQaMJIUkFCRF8Q0WqL10BzOUMIeH2mhwB4x/T5IwA5zNwRkVXHJMuzEsDM45k5l5lz69at67Fb6jnnrPDPfpxSPeSWQVGcZJa8/4qWCnqiMRMmgw4RalWpgKFdGycsk+rWg0kFCTP3Zub2Fq9pAPYSUQMAMP7us6iiAID5V2xkHINx3nkAMph5qanNQmYuNj5OAHC+k7pkkyoqgLCSE4Jc0k5wsiJpliL/SzpxzlneA346oXldNdeUiPD8dR0Tlnk8xZPMiaq2pgOIWk4NBzDNosxMAH2IKMvYGO9jHIsyFKVXI1GhFGUAgLUO6/IVJ+aIYSBVVgKyEBX3mRW0FXwY8TtOW5cm4TH7TvXFl+gTNBbAlUS0EUBv4zOIKJeIJgAAMxcBGAPgG+P1jHEsyg2IESQAHiCi74hoBYAHAPzaYV2+IhLJViV3Jgkb8awCPxVAnbOYaCuXtKjjKoCfRqMpjZAgMVRQvZi5laECKzKO5zPzCFO5iczc0nj9O6aO5sy8LubY48zcjpnPY+ae5u8T1ZUKqDAFTBZW++YLm/reB1mseqqP722UK0f4TQj2QEStzzTu+eqRHkocgJOR6vtBek2fAD92SNJl10XV/xEmM2O/mXzXRVLqCdrLWRWv3XJ+8kJJaFq7Km5JoYlVWNGCROOJrs2yQ/MA2mRELrM8PaBdoO13b10Xtav66zPz+/5t0Tet8tun9pJEP4KqSZPJYkb5cqFQCQBA/45nC9eRSLXQoGZ4Q4hbUSVButYJv8r1vf0n+rXB0t/76zNTNUVT0qYrWpCEiFSxCgsbFXxekky7PzWMLKI0rGUdR6t21YrorcBfQYWRhbbMDxdakCgm0aD3zh0X+t5+6/r+2uabqaso7a/fBB1+wi3dmlvHHkunwbdcOv0zaYAWJAnw415t08Dec11FJNPaVdUM7i9e3xHjbuzsezsP9m4lpZ5UyGshSjo52KaiQcFXj/Sw/U5bbWlc8epN52PSbV0Da1/VgPnL3MaoWcV/i6sODdMjHL6KMOAyh9737rRfPZ+lIC1tKgqSprXTNzqCFiQO6d/RewIdMzUrV0D31sHF/pKdSU1GMEMRUn0mF8VJ6H9RZC5I7NRnHRrWRK0qciy2ftG5oe1zl04x6QB/bbZaKAgFowWJQ/7Qv23QXZDCSzeIZUaMZUiSYHSa8BB0Ola3VKmYgb8N62L53bk6V7tj2p3t/6pdC5IEmHXK9RQs11VgZ9HjFasVQQpqHQJf2ahoXsV1ScVrHwb8jCtWs7L/KmYtSAIinfOEq9zUvbSV/2lM0wUV18WPNtLIRsCS8becL2WC98AVLdElIBcCLUgCokaahP6IbhK3rFcN9WtELMJUPfeP9D1HmqWb6qCNfVI8/4QdqboiaXd2/PW/rnNDzH64h+9t95Hkof9Qn3MwNaDAslqQBES6zLI6NKqJFwZ3xNR7L0beQ92Vtn1fT3mBFlvWU+dfAwB/uLb0nluvc+tJrT+3afyKt5zkp33JqF64pkPpzXA/rKli1Y7TfXAQfXpAO+TUrlLq2IXNa+s8NA7RgiQggtbJy+SGCxqjRqUKyMyI3E6pKiRfH+5/+JAojbJKD1rnNZarkvjT9fGJlGRvtterXikul4sK9VnHRvLVN7k52ZjzSM9Sx85PE/XzwE7iIYSSUbYyIIUIv5P2BEFUOKaadVCUoOz8/VAHWU1UGmXJNbQA4q91al75eLaNvSboLkjhus4NkZvjf3oCvSIJiE6SZ6ChQvJocvulzTD8IrWRhlVYurwyNOL5PzZJGlZvxEsSFftyOnRJ2ST9psU+cHV7+eGq00iz9RMZxtT6WgnReM383vDhmbRwu9R6Y4mOgY2yKmP2wz3QatSn0ttY+PgVOFMSufrXnnc2rj3PH7WD1YpExRgvex8GiEQsfmPBNszfdEB+5WmOqnFGC5IkLH2yd5lKriRCRvlyWPb7K1G9kj+3VXbViij68aQvdQM/L6QyypFvEYUb1JSvXrKiRNEIEiuc/FiR9G5bH73b1kfOyBnS69bIQau2klC7WiYqZuifySnZVSv6NgjPfPByfPybS9E1Jxuv3mTt8SxCtN+10yBqcVCxqFSoBFOZTAVjyZJRvTCqXxvf2zGjVyQBwelktqWIutUzUbd6JibfLSclbSyNs6tg7HUd0KtN6vt4tKhbFZe1qoN5Gw+gYkY5nDxdIj2qAVB6O2xUvza44QIdMicR1TIzUHzav1U1ELGmu65LQ7w2dwvuuMz/GG6AoCAhomwA7wHIAbANwA3MfNCi3HAATxofn2XmSURUHcA8U7FGAP7DzA8S0V8ARG3xqgCox8y1jLrOAFhlfLeDmQeI/A9Bc+slObinRwtUqahlehgY0rVJ0F2QAhHhtkuaYd7GA7i4RW3cmNsYV7SR66sSyx0KAk+mOu/fczFmrd2LZ2es9bWd2tUykf9kb1/bMCO6zhoJYBYztwIwy/hcCkPYjAbQDUBXAKOJKIuZjzBzp+gLwHYAUwGAmX9rOv5K9LjBcdN5KStESsXxql4J1SSbA3/8m0txXeeGUuvUpC5Xd2iAzIzUT0+b7XMueL9pVqcqRphWCff2aBFgb+QhOnoNBNDDeD8JwBwAj8WU6Qsgj5mLAICI8gBcBeCdaAEiag2gHkqvUKIMRUQQpRV+q7baN6yJVvXTK9S2JnzUULgnMu/Rnr4ZckQZ2rVJnIe7nzx61bnK2vIT0RVJfWbebbzfA8BKudwQwE7T513GMTNDALzHMaMrETUF0AzAl6bDlYgon4gWEdEgod6HgFR13tNoAODhPucoa6txdhVpuU7seP66Drire3qsElSSVLwT0RcArBwpRpk/MDMTkddp9hAAt9gcf5+Zz5iONWXmAiJqDuBLIlrFzJst+n0ngDsBoEmT9NB7azRe8HPxW7li6qvLNOIkFSTMbLtjQ0R7iagBM+8mogYA9lkUK8DP6i8gsqk+x1THeQAymHmpxblDANwX058C4+8WIpoDoDOAOEHCzOMBjAeA3NzcMmkiVRbykGsSoGixe+slOekdqUGTFFHV1nQAw433wwFMsygzE0AfIsoioiwAfYxjUYbCtF8ShYjOBZAFYKHpWBYRZRrv6wC4BMAawf8hUFR4G9/VXVvTeMVvnXw6MPradhjYSRt2uGHJE72w5IleQXdDGqJPyVgAk4nodkSsrm4AACLKBXA3M49g5iIiGgPgG+OcZ6Ib7wY3AOhnUfcQAO/G7Ju0AfAaEZUgIgTHMnNKCxI/1Q56/0WMJaN6pbSl0wU52TinfnU80lfdPobGGemScTWKkCBh5kIAcWKVmfMBjDB9nghgok0dltNlZn7K4tgCAB08drfMcVGL2gCA7q3rokPDmli161DAPUot6lVP7Ye9WmYGZv728qC7oSkD6HV7wPip2urUuFapcI0S5aMAAAZ8SURBVNj9JQdT1Gg0GkDH2goMHSFFo9GkC1qQaDQajUYILUgCRm+HazSaVEcLkoDIKB8RIRV0iHqNRpPi6M32gBjWrQn2HDqB+3u2DLorGo1GI4QWJAGRmVEejytOPqPRaDR+oPUqGo1GoxFCCxKNRqPRCKEFiUaj0WiE0IJEo9FoNEJoQaLRaDQaIbQg0Wg0Go0QWpBoNBqNRggtSDQajUYjBHEZCENLRPsRSbzllToADkjqjkx0v9yh++UO3S93pGO/mjJz3WSFyoQgEYWI8pk5N+h+xKL75Q7dL3fofrmjLPdLq7Y0Go1GI4QWJBqNRqMRQgsSZ4wPugM26H65Q/fLHbpf7iiz/dJ7JBqNRqMRQq9INBqNRiNEmRMkRHQVEa0nok1ENNLi+0wies/4fjER5Zi+e9w4vp6I+jqt089+EdGVRLSUiFYZf68wnTPHqHO58aqnsF85RHTc1PY/Teecb/R3ExH9lYhcZxwW6NdNpj4tJ6ISIupkfKfi97qciJYR0Wkiuj7mu+FEtNF4DTcdV/F7WfaLiDoR0UIi+o6IVhLRjabv3iCirabfq5OqfhnfnTG1Pd10vJlxzTcZ90BFVf0iop4x99cJIhpkfKfi93qIiNYY12oWETU1fefb/QVmLjMvAOUBbAbQHEBFACsAtI0pcy+AfxrvhwB4z3jf1iifCaCZUU95J3X63K/OAM423rcHUGA6Zw6A3IB+rxwAq23qXQLgQkRS1n8K4GpV/Yop0wHAZsW/Vw6AjgDeBHC96Xg2gC3G3yzjfZbC38uuX60BtDLenw1gN4Baxuc3zGVV/l7Gd0dt6p0MYIjx/p8A7lHZr5hrWgSgisLfq6epvXvw8/Po2/3FzGVuRdIVwCZm3sLMJwG8C2BgTJmBACYZ798H0MuQ0AMBvMvMxcy8FcAmoz4ndfrWL2b+lpm/N45/B6AyEWW6bF96v+wqJKIGAGow8yKO3MVvAhgUUL+GGufKImm/mHkbM68EUBJzbl8AecxcxMwHAeQBuErV72XXL2bewMwbjfffA9gHIKmDmt/9ssO4xlcgcs2ByD2g7PeK4XoAnzLzMZfti/Rrtqm9RQAaGe/9vL/KnCBpCGCn6fMu45hlGWY+DeAQgNoJznVSp5/9MjMYwDJmLjYd+7exjP69hyWraL+aEdG3RPQVEV1mKr8rSZ1+9yvKjQDeiTnm9+/l9lxVv1dSiKgrIjPhzabDfzTUKH/xMIER7VclIsonokVR9REi1/gH45p7qVNGv6IMQfz9pfL3uh2RFUaic2XcX2VOkKQtRNQOwJ8A3GU6fBMzdwBwmfG6RWGXdgNowsydATwE4G0iqqGw/YQQUTcAx5h5telwkL9XqDFmrm8BuJWZo7PwxwGcC+ACRFQmjynuVlOOeGwPAzCOiFoobt8W4/fqAGCm6bCy34uIbgaQC+BFv9owU9YESQGAxqbPjYxjlmWIKANATQCFCc51Uqef/QIRNQLwIYBfMfNPs0VmLjD+HgHwNiJLYyX9MlSAhUb7SxGZxbY2yjcyna/89zKImy0q+r3cnqvq97LFmADMADCKmRdFjzPzbo5QDODfUPt7ma/XFkT2tzojco1rGdfcdZ0y+mVwA4APmfmUqb9Kfi8i6g1gFIABJu2En/dXmdtsz0Bkk6kZft6sahdT5j6U3qSdbLxvh9Kb7VsQ2fxKWqfP/apllL/Oos46xvsKiOiM71bYr7oAyhvvmxs3ZzZbb+71U9Uv43M5oz/NVf9eprJvIH6zfSsiG6FZxntlv1eCflUEMAvAgxZlGxh/CcA4AGMV9isLQKbxvg6AjTA2ngFMQenN9ntV9ct0fBGAnqp/L0SE6WYYBhIq7i9mLluCxPjR+gHYYPzYo4xjzyAivQGgknEjbjJ+YPNgM8o4bz1Mlg1WdarqF4AnAfwIYLnpVQ9AVQBLAaxEZBP+/8EY2BX1a7DR7nIAywBca6ozF8Bqo86/wXCMVXgdewBYFFOfqt/rAkT00D8iMnv+znTubUZ/NyGiQlL5e1n2C8DNAE7F3F+djO++BLDK6Nt/AFRT2K+LjbZXGH9vN9XZ3Ljmm4x7IFPxdcxBZKJSLqZOFb/XFwD2mq7VdBX3l/Zs12g0Go0QZW2PRKPRaDSS0YJEo9FoNEJoQaLRaDQaIbQg0Wg0Go0QWpBoNBqNRggtSDQajUYjhBYkGo1GoxFCCxKNRqPRCPH/ATPwhC1cHnOtAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"sample_rate = 50e3\n",
"t = np.arange(10000) / sample_rate\n",
"\n",
"x = [5e-3, -1.23] # noise amplitude, phase\n",
"v = -70e-3 + np.random.normal(size=len(t), scale=1e-3) + line_noise(x, t)\n",
"\n",
"plt.plot(t, v);"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"('guessed:', array([ 0.00501016, -1.23262113]))\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAD8CAYAAABkbJM/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztnWmYFNXVgN8zDLsoq4hsowgqLiAO4C4KAkoMRlww0YwLGo35YmJMxC0YSSIuiWYxGqIY1LjHiDsOCMYoIAMiIIKgIovsICDIOuf70dVQ01PdXdVVvcx43ufpp7tu3bp1upZ77nLOuaKqGIZhGEY6ivItgGEYhlEzMIVhGIZh+MIUhmEYhuELUxiGYRiGL0xhGIZhGL4whWEYhmH4whSGYRiG4QtTGIZhGIYvTGEYhmEYvijOtwBR0rJlSy0pKcm3GIZhGDWKGTNmrFXVVuny1SqFUVJSQkVFRb7FMAzDqFGIyBd+8tmQlGEYhuGLUApDRJqLSLmILHS+myXJV+bkWSgiZa70eiIyWkQ+EZH5IjLEte8CEZknIh+JyJNh5DQMwzDCE7aHMRyYqKqdgYnOdhVEpDkwAugN9AJGuBTLLcBqVe0CdAXedo7pDNwEnKiqRwA/CymnYRiGEZKwCmMwMNb5PRY4xyPPAKBcVder6gagHBjo7LscuBNAVStVda2TfiXwgJMfVV0dUk7DMAwjJGEVRmtVXeH8Xgm09sjTFljq2l4GtBWRps72SBGZKSLPiUj8+C5AFxF5V0SmishADMMwjLySVmGIyAQRmevxGezOp7GVmIKsxlQMtAPeU9UewBTgXte+zkAf4CLgHy4FkyjfVSJSISIVa9asCXB6wzAMIwhpzWpVtV+yfSKySkTaqOoKEWkDeA0dLSdW8cdpB0wG1gFbgRec9OeAK5zfy4BpqroT+FxEPiGmQKZ7yDcaGA1QWlpqywcahmFkibBDUi8BcaunMmCcR57xQH8RaeZMdvcHxjs9kpfZq0z6AvOc3y/G00WkJbEhqs9CymoUALt2V/Ls9KXsrjTdbhg1jbAKYxRwhogsBPo524hIqYg8DKCq64GRxHoH04E7nDSAG4HbRWQ2cAnwCyd9PLBOROYBk4Bfquq6kLIaBcA/31vMr/49m6enL8m3KIZhBCSUp7dTiff1SK8Ahrm2xwBjPPJ9AZzika7A9c7HqEWs27IDgK+27syzJIZhBMU8vQ3DMAxfmMIwDMMwfGEKwzAMw/CFKQyHN+au5JJHpuVbDMMwjIKlVoU3D8PVT8zItwiGYRgFjfUwDMMwCpDdlcr1z85i/spN+RZlD6YwDMMwCpDP137NCzOXc+2/ZuZblD2YwjCMLHPR6Kn0+t2EfIthGKGxOQwjp+i3MCLIlM8sSMG3hTWbt9OicT2KiiTfomQF62EUKGPfW8xlj76fbzGyhtTO98n4FrNi4zf0/N0E/vzWwnyLkjVMYRQoI176iEkLam+49m9jT8Oo3azcuA2gVr+3pjCMnGI9C8OouZjCMAzDKGAKqTNuCsMoCG59cQ7XPf1BvsUwDCMFpjCMguCJqUsYN+vLfIthGAVHIY3imsIwDMMwfGEKwzAMw/BFKIUhIs1FpFxEFjrfzZLkK3PyLBSRMld6PREZLSKfiMh8ERnipN8nIrOczyci8lUYOY2aw/Kvvgl1/I5dlVTaeuGGkRXC9jCGAxNVtTMw0dmugog0B0YAvYFewAiXYrkFWK2qXYCuwNsAqvpzVe2uqt2BvwAvhJTTqCH0uWcSUz7N3DO6y62v8382eW7UIgqp+RNWYQwGxjq/xwLneOQZAJSr6npV3QCUAwOdfZcDdwKoaqWqrvU4/iLgqZByGgVCOoe9nbuVi/4xlXcXeT0K/nh19oqMj40KVeWmF2YzzcKCFDRzl29kxhfroy00Mq/UQprujhFWYbRW1fjbuRJo7ZGnLbDUtb0MaCsiTZ3tkSIyU0SeE5Eqx4tIR+Ag4K2QcmbMotWb+WLdlnydPhBfhhzOySXpHPhWbdqWG0GyyFPvL+XC0VPzLUbemfLpOhavLcx36Dt/+R9DHpySbzFqDGkVhohMEJG5Hp/B7nyqqgTrPRUD7YD3VLUHMAW4NyHPUOB5Vd2dQr6rRKRCRCrWrIneJb/fH//LqfdMjrzcZRu2MmfZxipp4z9ayXMVS5MckZppn63jhFFv8Z8PlkUhnmFExkX/mEqfeyfnW4zcEVk4g0IajIqRVmGoaj9VPdLjMw5YJSJtAJzv1R5FLAfau7bbOWnrgK3snZ94DuiRcOxQ0gxHqepoVS1V1dJWrVql+zsFw0l3TeLsv/6vStqPHp/BL5+fnVF5C1ZtBmDmFzXDPsBiSeWG3ZXK+i078i2GEYJCGpgKOyT1EhC3eioDxnnkGQ/0F5FmzmR3f2C80yN5Gejj5OsLzIsfJCKHAc2I9TyMWoLFksotd78xnx4jy9kQQmms3rSN7buSdvKNCKgp1ziswhgFnCEiC4F+zjYiUioiDwOo6npgJDDd+dzhpAHcCNwuIrOBS4BfuMoeCjztKJackePTGbWUQnmMxn+0EoCvvtmZcRm9fj+xoFZ9K3gyuPm9fj+RHz9R+Nc41AJKqrqOWM8gMb0CGObaHgOM8cj3BXBKkrJvDyNbLlmybitrt2ynRwdPNxQjAwqlwq0thG0ITfjYa7TZiJKJ81ezfssOGtQtolG9vVVzIb0K5ukdAafcM4lz//ZevsUwjGqIjQFmxLINWxk3a3lmB4e45j1GlvOdv8TnNgvv3pnCMIwEKiuVse8tZtvOwh9T9kshtVILhdWbk5tuD3nwPa57elYOpdnLZ2sK0wQZTGEYOWT+yk08OPlTX3lTNdLmLNvIM9OXVEl7ZfaXfLXV38TuotVfp3TWenn2l4x46SPun1Dzl9osvDZq4dDrdxOT7lu1aXsOJUlG4an5UHMYtRHV6Cx5Xp+zggb16nDaoftHU2ANZ2JE4+Bxc+QLe3YAYsMHP3nyA046pKWv4/v98W0AFo8a5Ln/6+27ANiY4UTxXW/M960YE9m+azf1i+tkdGwqbE4oh0R8sQtJ6VsPI4tc86+ZXPbo9Jycq5AqhJ27K9n4zU4enPzpnm6/qnL/hE+ycr7tuyqBwvF0z1RZLFi5mUNvfYPX5kQY2iQHtc2mbTu58rEK1n4drFW+fsuOzOcJjLxgCqNAWLp+a75FiIxr/zWTbr95k7vemM/Pn4mNA2/evoudu/dqNXcvzmuuoJAUYK6Yszzm+T/h41V5liQYz7y/lPJ5q3gooKK8+vEZXPf0LP42eRGfrfkaiM0rdL/jTeav3OR5zIKVm9ntEY14xhcbeD1KRRuGiA0NCulVMIWRY8rnreKCv1f3RVy6IZzCKCRjmDfn7a3wvt6efuL4J0/mLrrs5m3ph5lqp7Kq+qde/GA5i1ZvzpMsMb7cGOsR3v3GAgb+6R0gNmz51dad/PPdxdXyf7JqMwPu/y/3lVfvqQ558D2uqXW+IrGX+rM1WwrGAMMURgaoasYRLq95Ygbvfx5xdMwA7NxdySP/+5yduyt95V+wcjMlw19l1tKvqFi8fk9L0Dc+at+oWtR+6vmfBbB8EYGNW3fy8Duf1ViHzmTtiJ89M4t+f/xv2uOTeYiv37Ij8DVZsLKqgnIfvmOX9/P48YpNlAx/ldnLvtoTkHLW0poR/iZKfvvqvPSZcoApjAT8vAJPTFvCkAen8KbjRVtozF8Ze8kqFldXTI+++zkjX5nHY1O+8FXWW/NjE9Wvz13BeQ9N4fQ/vB2JjFHWv0E6V58HiJqqCsNfmM1vX/2Y6Ys3BBcsDZMWrKbTza+xyUevJx1vfrSSb3ZE3wo9ZmR5tbRFq7+mx8jyas/QrjQLVw24P72CSmSi05gYX6DvmidZaFwsXlsYQ9amMDLg09WxVvayDYUxyQpVn9F3PomtJfHG3Oov2eZtMQugLY4lUNblyvC4IENsfs6xdccunpy2JH3GhHPHLaX89siC8OeJC9ldqSxcFXueMu3FfPTlRq56fAa3vjiXHbsq05bzdch7Hw9V/t9PqkaH/ud7iwOVk2p1Ra+/UEM7ebUKUxghSXxpMiZHL0O8corz2JTF/OHNBaHL3bRtJ7f8Z07SVu7UEKvoJWPZhq1s3bG38kulY0a+8jE3/2cOn/noYeS7YpKApk3xRsDMJRvocuvrnr1H93/amWT4J5ZPOeaON3nqfX/KFTKbP0umgN1FPf3+Eh56+7NIzpdLPly2kbfmRzHMWvVB7PfHtznTmevJF6YwQvKPd6o/0Pngb5MXpc2za3clbzhd+/g79+txH/GXt9Ifm46/vrWIf01bwr+mVa2s4hXVlh3BWrXXP/th0nHtOCfdNYmyMe/vPVeKvOsCmnxC4VdMsHdOCtiz0NdLH365Z3/Q0CCVChu2xpR/Nqn0oZWHvzBnT28oqBKvWLye9z7NfNXGZKzZvN3XgmqX/7Mi1HkWrf6aW1+cWy3t4xXe1mO5whz3MiDMBGi6cd4glAx/FYBpN/f19ExNPNPuLDadK53/lXiKuKlo0FYzxFbda9+8UZW0DxMmPKcv3lBQjk1+uPKxCurWiUbqR9/9nPJ56Vuz7tuSLX+YIGS7F3feQzFLxGTOmZnS83cTkpb7vQziyb2XZCniHz1ewacFGCLEehgJ+FEGY31OGPvhldlfcvhtb7A9xBj5B0uqVqLJGpWp/tpXW3fwnb+8k7T1tC3NhOrDTis3SrzkHfzAu5GVv7tSeX3OisANgCenLaFk+KtVhsPiPOLjOpTPW8Vrc9zzS5nXnnGjhGSleD0Kfp7fVBIl7vtqa/hJ+0zOu3nbTj4NarUXkAcmhe99p+L7D0/j+w9Py+o5osQURp6587X5fLNzN2s3ew+ZTF6w2tNefpcPBaMaGwf2qtgSeX3uSuYurx7rSZ3XNFkls33X7j09HSjcYZw3PVrho//7Gdf8ayavzA7m8BUf/lv3dXWT05Gv+Dd/jPpSpdJ7fnViEJni9/qvHpXqKXdP4qYXMh/Wiq8gmYoL/j6VvhFZ7SXjnvH+5vf+/nZm3v01DVMYPvFaDSto5egVHTPeuk22NOulj073tJePd429iFcOUz5bx/AX5nDHy9UrMRGY92X48dCvt1VVRskqpsRr5efaxZXV6s3b+GBJerPWoD2FlY7jWCbzG5mgqlXWXJ/p6hlu3bGryjOwa3clZ/3pHSYtSB1/a3bCuvCJ+LnON70wh5v/M4fBD7zLw/9LPicX94Pww5L1W1NOnCfrGcTFXezR0/3NSx9V2Y6P52caiiVK7nx9fr5FyAmmMHxy2G1vVEsLWuGmio4ZhLnLN7LBNQxw9RMzPPPFTWfXfr2DFRu/qTJ0JSKc9efqFhdhx5b9eFJDdUXjxX8XrmXo6Cn0+t3EjMaHs0XiRP6azdu59smZaU2V3120jp8/86HnvjUJPcy1X+9g3opNDP936jXet3oMFa7fsoOuv36jipLVFIM9T72/hCenLeHDpV/x+9diFZ/Xc9D799E8vxAzkghKMgu3u97IfWVdMvxVlmUQnWH2sq98WVau9ei9JiMTJ8pMMYXhE6/78dyMZdUTPYh6rd4rH0ttgbHDGa5a4opPddJdk7joH1OTHpPJ8MjWHbuqTaT/2aMi8JqU/ZuPVuFtL85l6mf+veLdFkGDH3iXOWla30H57SvzKBn+6h7/gd++8jEAfyz/hFdnr+DFNIH0/CpTSF3Bp+PztVvYumM3/3jns4yMDdIRpHLq/fvkPeFUpJI6k0WhKivVl3UTxOJS+WH+iuChVb7713f5ocuyLxmJkZKT1SFL1m2lx8hyHn4n+jlEL0IpDBFpLiLlIrLQ+fZco1REypw8C0WkzJVeT0RGi8gnIjJfRIY46R1EZJKIfCAis0XkrDByBiEbevrQW9/IaWgJL8csr4Bt6diyfdce66dEdu6upOuvx3N7wjCBF7kyBXR7cX+49CvueCW9bEFInNj/6puqrcDNaXpNmczvRFXhZ+fxSy9bputKbEiYSHeLn4nv05/fWsip90xO24iYtGA1Qx6Mpjebaqgz6DBoskgD8Rh06YYuoyJsD2M4MFFVOwMTne0qiEhzYATQG+gFjHAplluA1araBegKxGewbgWeVdVjgKHA30LKGRmPT1lcZfudhf5svUcFGOP0M6GdjKgmnVWVI0aM5943vU0w445XVa19CpN0oTf+NLHqQkledaufCjfIPc41yWJC5YtkXt5RG02s2rSNh97+lGlOTzW+lkoyMonekMxz3mu4ME6yYWS/rNy4LZKQMkEJqzAGA2Od32OBczzyDADKVXW9qm4AyoGBzr7LgTsBVLVSVeO1rwL7Or/3A76kANj4zU5uG5dZq/Xv//Xv4PfJquoTgpWV6tsyKgqemb406b5CiZwZJ9kQxTqnkrw+ybxBnMTWbNWy/Zw/fR4nZ9I9U1ye8CLR3Ed3GYVmuplusj4qysa8z6jX5/PZ2r3v1ONTFkeynED8vpclGWK65JHk1zzsin7H3TmR0++dHKqMTAirMFqratwmcSXQ2iNPW8Bd+ywD2opIU2d7pIjMFJHnRCR+/O3AxSKyDHgN+L9kAojIVSJSISIVa9ZkFqbD7xxDt9+8mXL/V1t3VLGCyRQvL9jrn53FIbe87uv4KBppqUJo1JRFb+JrI6eKWeRF/Pq9+MFy3l2UPqTJxgj8EIYnMUEN0+KetGA1W3dmM2aYJh2yjPxMPk/j9nIHmL+y+jzDbeM+4uS7J0UhFpB8zmPxuuwGDAwyMR4VaRWGiEwQkbken8HufBobpA/y9BQD7YD3VLUHMAW419l3EfBPVW0HnAU8LiKesqrqaFUtVdXSVq1aBTj9Xm5+YW76TD742TOz+PkzHwZyJvJ7wV6cFW0nK1lF5Cea67QAE9FRkknFPGfZRl+TnW9/soZT7p5UpfHgHlJIVXG/6nPhniCV/67d4SvibTsrWbo+uwEyT/vD5KyWH+chn34OP33Ke22VOllwEErVC0/Hyo3+TZT98N6n6yKKX5WatApDVfup6pEen3HAKhFpA+B8e828LAfau7bbOWnrgK3AC076c0AP5/cVwLPO+acADQB/CzZngHtti5dCVMzxh2D7zugjm2ab+Pv0/uL1fPRl6uGCFz7ITw/jD+XBgySe/df/pRxLjvPrcXNZsn4rK77yPzaczQnpe52AkCsirFgufnha5CaoX0TYivbyF4qKoiJ/9ypZrt2VWs1wxMsZ1C87shD9+MlpmSswv4QdknoJiFs9lQHjPPKMB/qLSDNnsrs/MN7pkbwM9HHy9QXiT8wSZxsROZyYwogoLGxqHpsaPuzH/z01k35/zK4HajIybUi5u8/LNnyTE6uu6R7rdaQiE0uvTLj7jeqKyWvMe87yjawPMJkc5NYkLjYUBf9btDZSJ7cJH2dmmZOsJTzm3eyZhtbxqTC82Lm7kk43v8YpHsNYYUzm73zt44yPzRdhFcYo4AwRWQj0c7YRkVIReRhAVdcDI4HpzucOJw3gRuB2EZkNXAL8wkn/BXCliHwIPAVcqjle8mzX7kr+kmA945dP12xh0erMY9yk+qc7d1emrHgyvUruikQ1tlRmKoKexyvybJQxiMKS7v8km2M47yH/JpiZ+A8UGtt27k677sXDKSI4q2rKSK6zlmZnMrwoxLWP1wNec2HphqWG/3s2S5L0woIYwiQj16H4Q0WrVdV1OD2BhPQKYJhrewwwxiPfF8ApHunzgBPDyBaEKtdclcenLGZ3pfIHj7WDUxFVhZDKaSuqRZvSeVqnWzfi9IBj116hTwpx2dOgt/CzABFFH0kSdiPxnM/PWEaX1vsEEyRH/LH8E95LsbbJqk3b+O2ryVvO6VbO+2prdiZyQ3QwUvak0s01PT19KU9PXxpp1Nz/uUz5d+zOrcWihTdPYO6Xm/hwWWams5k4qHmNUaeqR2cv+4ov04xr+6n05qWQdaGPwG9BzAKTrcH8ZUDrpWRck6FNuzugYhjPar8E8VovVNJV6OmGDdflwB/k2Yrqrf4wQ1Kp3pV8NHnWbdn77oVddyMopjCoOracq3HyTLnu6Vkp9785byUt96lfJc1LgaR6fYL2rNLhDsHt5vaX57ElwDrUydZ98DKdzJRshNJIh9cj59VoeGPuSkpLmlW7vzWJenWiiUaUqnf6hMc8ZD7uq5tJC1bzzPvRTEp/mKVhOz+YwiA/rYRssWrT9mqtf693q1DG0/2GjwZYnSQEfBRE4UuRKec7i/2k4uvtu7j6iRmUtGjEE8N68+rsFfQ93MvtKTq8wk1kexTRb/k/SOGI6OUUmO/H/bJHp0dW1tgpiz3TN36T/d6bKYwCpDYpsJrCJmdOJ98VSzJ2O2Pli9dt5aS7YtY6D0awBkOqHrVXJec34GYyomqopJpHySX5mIdLds4gFnuZYtFqv6UkGyYyChOvOZYorMz63BudxzOEV7jZmkvKVo861QR/tkh2hXKhu0xhFCCFaD30bSHKkBFhcHsCqyqPvrs4K+fJtid4Il/nIWAepPeB+elTH7B607aC7WH6wSukUNTYkFQBsiSCwGhGzWazKwJqj5HlKQMk1iRuT+PNna+20ksffkm94iKO6dA0feY8I3j3MhLXpskG1sOg8EJ5eJkFhqEmt5qM1NF0DX8U+ajpanrHvjIH1Zj1MIBtEa+IF5aa/uAa315SWS/5IZ+Nm3/PXMa/Z4aPNp1tklUPu3KgMayHUYBErTAS14s2jGwRxPM9l0QR/bfQyUI8w2qYwihAorYUSeZpbRiFRrZ611E6d+abZNcoF8YypjAKkNoQQsIwMiFMyPBvO7mwkjKFQTQr1BmGYeSTXAy6mcIwDMMwfGEKwzAMoxaQi5ESUxiGYRi1ABuSMgzDMAqGUApDRJqLSLmILHS+myXJV+bkWSgiZa70eiIyWkQ+EZH5IjLESe8oIhNFZLaITBaRdmHk9PE/slm8YRhG1qkJwQeHAxNVtTMw0dmugog0B0YAvYFewAiXYrkFWK2qXYCuwNtO+r3AY6p6NHAHcGdIOQ3DMGo1NcEPYzAw1vk9FjjHI88AoFxV16vqBqAcGOjsuxxHGahqparGF6vtCrzl/J7knCdr5CKOvGEYRjapCXMYrVV1hfN7JeC1BFhbwB1NbxnQVkTiYSFHishMEXlOROLHfwic6/z+HtBERFqElNUwDKP2UghDUiIyQUTmenyqtPo11h8KInIx0A54T1V7AFOIDUUB3ACcKiIfAKcCywHPCIEicpWIVIhIxZo1awKc3jAMo/aQix5G2mi1qtov2T4RWSUibVR1hYi0AbyWcVsO9HFttwMmA+uArcALTvpzwBXOOb/E6WGIyD7AEFX1DIikqqOB0QClpaW1P8KYYRhGngg7JPUSELd6KgPGeeQZD/QXkWbOZHd/YLzTI3mZvcqkLzAPQERaikhctpuAMSHlNAzDMEISVmGMAs4QkYVAP2cbESkVkYcBVHU9MBKY7nzucNIAbgRuF5HZwCXAL5z0PsACEfmE2LzI70LKaRiGYYQk1AJKqrqOWM8gMb0CGObaHoNHL0FVvwBO8Uh/Hng+jGyGYRhGtJint2EYRi2gJvhhGIZhGAVATfDDMAzDML4lmMIwDMOoBdSEWFKGYRhGAWBLtBqGYRgFgykMwzCMWoBNehuGYRgFgykMwzAMwxemMAzDMAxfmMIwDMOoDZhZrWEYhlEomMIwDMOoBWgOuhimMAzDMGoB5ultGIZhFAymMAzDMGoB5rhnGIZh+MLWwzAMwzAKhlAKQ0Sai0i5iCx0vpslyVfm5FkoImVOWhMRmeX6rBWR+5199UXkGRFZJCLTRKQkjJyGYRi1nZowJDUcmKiqnYGJznYVRKQ5MALoDfQCRohIM1XdrKrd4x/gC+AF57ArgA2qeghwH3BXSDkNwzCMkIRVGIOBsc7vscA5HnkGAOWqul5VNwDlwEB3BhHpAuwPvONR7vNAXxGRkLIahmEYIQirMFqr6grn90qgtUeetsBS1/YyJ83NUOAZ3Ttrs+cYVd0FbARaeAkgIleJSIWIVKxZsyazf2EYhlHDyYUfRnG6DCIyATjAY9ct7g1VVRHJVOShwCWZHKiqo4HRAKWlpbkYxjMMw/hWklZhqGq/ZPtEZJWItFHVFSLSBljtkW050Me13Q6Y7CqjG1CsqjMSjmkPLBORYmA/YF06WQ3DMIzsEXZI6iWgzPldBozzyDMe6C8izRwrqv5OWpyLgKdSlHse8JbmwsjYMAzDSEraHkYaRgHPisgVxKycLgAQkVLgalUdpqrrRWQkMN055g5VXe8q4wLgrIRyHwEeF5FFwHpiQ1aGYRhGHgmlMFR1HdDXI70CGObaHgOMSVLGwR5p24Dzw8hmGIZhRIt5ehuGYRi+MIVhGIZh+MIUhmEYhuELUxiGYRiGL0xhGIZhGL4whWEYhmH4whSGYRiG4QtTGIZhGIYvTGEYhmEYvjCFYRiGYfjCFIZhGIbhC1MYhmEYhi9MYRiGYRi+MIVhGIZh+MIUBiCSbwkMwzAKH1MYhmEYhi9MYQC2+KthGEZ6QikMEWkuIuUistD5bpYkX5mTZ6GIlDlpTURkluuzVkTud/adIiIzRWSXiJwXRkZ//yPbZzAMw6j5hO1hDAcmqmpnYKKzXQURaQ6MAHoDvYARItJMVTeravf4h9ia4C84hy0BLgWeDCmfYRjGt4IurffJ+jnCKozBwFjn91jgHI88A4ByVV2vqhuAcmCgO4OIdAH2B94BUNXFqjobqAwpn2EYxrcCIftDJWEVRmtVXeH8Xgm09sjTFljq2l7mpLkZCjyjarMJhmEYhUpxugwiMgE4wGPXLe4NVVURybTCHwpcksmBInIVcBVAhw4dMjy9YRiGkY60CkNV+yXbJyKrRKSNqq4QkTbAao9sy4E+ru12wGRXGd2AYlWd4VfoBPlGA6MBSktLrYdiGMa3kmv6dMr6OcIOSb0ElDm/y4BxHnnGA/1FpJljRdXfSYtzEfBUSDkMwzC+1ZzcuWXWzxFWYYwCzhCRhUA/ZxsRKRWRhwFUdT0wEpjufO5w0uJcQILCEJGeIrIMOB/4u4h8FFLOlJhVrWEYRnrSDknbXTSVAAAen0lEQVSlQlXXAX090iuAYa7tMcCYJGUc7JE2ndjQlWEYhuGD4jrZ98M2T2/DMIxawH4N62b9HKYwDMMwDF+YwgDEYoMYhmGkxRSGYRiGw++/d1S+RShoTGEYhmE4fL+3Of+mwhRGjmnaKPsTU4ZRm9mnfijjTiMEpjByTFENmi+ZelM1i2nDCMQRB+4beZkXH9cx8jINf5jCILeOe7k6V9umDUMdf1Tb/ThgvwYRSfPt4MqTD8q3CAVHm/3CPYeJTLnpdOrWqTmNrtqGKYyQ3D3k6ED5c2WR9Z8fnxDq+D9c0C0iSbJLLtYA8MuZR7XJtwi1njb7NbTIDHnEFAbhVty7oGf7QPmLcvS0779vuN5BMjHfHX56qHJrMz06NKNrm6pDML8aeGiepKm91NYIow9dfGy+RUiLKQygaaN6oY4/7dBWvvMW2hxGUO/Q/ZvUp7Sj50q8eUEQ2hTQ0Nk5xxyYs3MdW0D3wQhPgVUNnpjCIPywRqdW/o/P5UPhR5F9OKK/Z3oqOWtrCy8ZxbnqFgakZ0nzfItghKROgT5byTCFkYYWjdP3PnZV+q9Cc9nDuOu8YPMrfggr/dndctcCj4pnrz4+42OzuYZkTWiRhqFJg9pvPtuoXp2U+y8sDTbknW1MYURAy338D2kV5fCKF4c6WfLaqJBW0s12pfmnod1pUJz6pXaTi3WVg3JGV6+Vkwuf7u2b5luEnOL15Pz8jC45lyMVpjAi4Een+l/pKkyFMvO2MzI+NihBK+KLj6udHrL7N2lAxxaNgNQV2D9+WJo1Gf59jbfFW+GppuqEUei7A/Tcc0mzNM63V51SbcWGpLgvT5h/++9rMu8FB8EUBuGHDeoGiEOfbsjy8hOT2/I3b1yPxaMG0SvPY9depsF+5nF+/72jQldyJU7lnStEoHH9YhaPGsSAI7yWto9x4iEtgJjVXLtme30Pfnh8dSezoCujFeLk9qGtm2T9HFec5P0u5FtRjv/ZKSn3Zypfsnqo5T71aJLGu/3YjrmpE0xhpCHqNo4CPzntkIhLzYyXf3JS4GNExHNI6pQu6SfYjzs4/EMdZcz/+y7sxo0DDwP8tYT9zGft17Auj5T13LPdpEHdav87WcWwbw0as3/8il5ZLf+XAw7ltEP3z+o5MiWdyXqmdUayZ3Dazf344Ne5G11IRSiFISLNRaRcRBY6355NIREpc/IsFJEyJ62JiMxyfdaKyP3OvutFZJ6IzBaRiSJSq2IB3DDgUBaPGuS5L0pHtFR1YJHAUe32y+jYRD76zQAObtnYV950L1Pa1rdIZJYlJx3SinrF/l+B845tx0/7dk6bTxP+ZaKCaOZD8fghn5PeTRr4U9wt96m/5/ekG/oEOkey/+e3Qu68fzTv0qGtm3D+sdlfADRZQ6JOkeRkNT0/hJViODBRVTsDE53tKohIc2AE0BvoBYwQkWaqullVu8c/wBfAC85hHwClqno08Dxwd0g5U5LuxZtw/am+yzowhE/Ad7sdyIU929MnjTlsYoWUDnfAw6NTKAk3qTzSE88eZS8s3fBgk/rFeZt0LyoSBnf3tvJKNTeVKG3Dut6vXbb+VTb0il9lVc8J43H72V05yGejotCoUyTcc77/yAdBns+athZPWIUxGBjr/B4LnOORZwBQrqrrVXUDUA4MdGcQkS7A/sA7AKo6SVW3Orunkuf1vQ8J0FK5YcChzLnd27chHUe13Q8R4eGIJ1C9nt9MH9RkL4Pf8hJznXtM20Dnv+/C7lW2i0PEFVJ0z/9JZu4c5N4nPY/fCqQw53gLjlR3PD7W//p1J3vu/+g3A7Lufd/C1asKQk3wyQirMFqr6grn90rAy36vLbDUtb3MSXMzFHhGvd+sK4DXQ8qZMQ0CDFlArOWVqrvup+5I1/30a2nlVQf6rbu8zrBvg2KK6xRF6ltQp0ioF6C73apJ1ZfxwR9EE07B6//+/ntHVRlSAajjQzEmXp8CskKOjCgaxi9ee2LkJ3/lpydx/4XdObyNd5TcxvWLaVQ3Zibd7/D9OSdJjzFTZtzaz9dcVxz3XxHg7gTfqfoB659sk1YaEZkgInM9PoPd+ZzKPtNXYyjwlMe5LwZKgXtSyHeViFSISMWaNWsyOnmqF/rJK4/LqMyw/OWiYyIry62H/b7oDT0ciq7Pgk24CL7HTOJxmty9mfbNs2c15TX017FFI24ceBhDevjv9A7tVd3kuPznp/CvYb0TzrcXP97lfhsOmQ4FpZpP8zp33TrC1QFMzLOxrkXHFo05x2evtW3Thp7Wfb/57hE8c1Vm73263oXbgg7Sm9VGNd8VFWkVhqr2U9UjPT7jgFUi0gbA+V7tUcRywO2u2M5JwzmuG1CsqjPcB4lIP+AW4Luquj2FfKNVtVRVS1u18h/TyS8lEY+7uivs3gc1T7ovlUe03zkMr5f6qLbp5zAevawnrRMsQfp3bc0QZ+Ivm3MWiib1bn0tyTBDFJx2WMwi59ITSlLmExGu6dOJA5v6n6s6z2PCtHPrJpx4SPIJ/iiv8Q0DMhuC6ZPCSsmr4dGwbp2kczxREtXATeI1js/1lZ1QQu+DW/gq43mPKABB7t0+CZZx7v9WiAtFhe3vvASUOb/LgHEeecYD/UWkmWNF1d9Ji3MRCb0LETkG+DsxZeGlhLLCQxf3CF1GEMe8rI9UOKLUr7u3t3BD//SVx6mdqyve0T8sTTrUlviPK27t51tEL+4672hO6OTvhQ1LvLfSrllDFo8axO3fPSIn502kQd2iwJP5fse8g/gJuUklj9eZRSS0ZVOUJDunu4fqR65U81ilAX2i3NfnshNLGNpzb++z8GcwwiuMUcAZIrIQ6OdsIyKlIvIwgKquB0YC053PHU5anAuoPhx1D7AP8JxjcvtSSDl9UhNumX/2a1iXXw44lKdd3Ws/8XnCjk+7x/3dQ2vJKpQq3fJc1iyaulJMJUvi33D/L6/j3D4WR7er7jEuSJXKy88tCOJRnAmpHK29jAQ8722KhylKA6GSFo18WwAGJeg8QjeP++vFD48vqXId0z36iaHz4+RyHfJQfR5VXQdUW8dTVSuAYa7tMcCYJGVUe+pVNVwTNY+ka83legL02gQnwShOf+ugw7nx+dl8vm4LqqnLbNWkPiUtGrF43VbP/SJQ0qIxC1ZtrpJ+6AFNeO/TdRFI64+gIVv8XEd3D+DtX57G5m272FlZWcVnZd8GxWzatouDWjZm8botgWRoXL+Y5o3rsX7LjirpT115HPdN+IT3P1+f5Mi9XHZiCY++uzjQeSGJcghaRsL2/JEDua/8E/7+388889815CjP9Hd+dVrGc1nu9zEqR8FDD/DnBa+qgU3kEzm1Syt+/z3v65INCmsKPk8c54xX7r9v8gmr4Wce5qusgUcmDx+RDr+mqUccGL4llepM6eToWdKct27os8faJCzHdNjbIou/wDedeTiPXZ5db+J0BAlr7tUQcDuONWtcjw4tGtGp1T5Vrm987mRYhsu7evWQju/Ugo4BK9AgDozg/Yx0aBFsvi+xjAZ16+xRson/a/GoQVzY07slnUxZpLt77lOc26MtdyVZPTNXrhKZnCbXQ32mMIiF6njnV6fRvlnyl+zqUzvxGx/j2+nGlb0evk6tgr1oN591eKD8XuTyQfu/02Pe0a1TKORE6hUX+Qo3AnCd432dKg6XF1rld/UrcuaRuVtyNdOw91Hcx8WjBvl6tlNRr04RY8qC+Q8J8N9fnlY1zeMyHHGg91BMpnido23ThoGVZng5at4QuCkMYh687Zs3SmuVUHZCCYf57G76wnnb6wcInw3BW4MQC1x4cEDFFJSDWzbmp6dXj5M15Nh2LB41iEb1iquZFZYmBE0L0kWPW+QM7dWBxaMG8euzu6aNJJqI10vb2DEprpPCKTCqVz1eTsyJMKJCiUaRBJHnjK6tAzusiUCHFMEkRYTnrz6+mvlxlMSft1T381wPE+orM+4RVt1OHAotdCViCsNFw3p1WPi7M7NSdrx1nSq+TZBHpeU+9QJFrZ152xm89Ys+GZ3LL2/d0Ifr01hh/axfF0ZfstfZ7vzSdhl3+f8QIFyDF5pm0jsX3DKoKxeWtueso9pkNJ496txwi2SlDmniX56orMviBhPNG8e+S0uaV1tCOSt1aopCj/Mwsb1lUNek8eD8EvTZ83bEze3zawojgUxNENMRd3ob6DHMcdZRsXmPVHb5iVTcegZX98mulUwUJD7jdesU0f+IA/aEcEhsUQV5/v1kffJKf63TsAsfZVqJtWpSn7vOOzpwLzNOmDkzCB6XLBmJHvjJcId797rmPzy+hPsv7M7Qntlbaa7NfrFebscWjbJmhNK9fVO6hbTaKsS+himMgNx7fjf6HhbcmuKC0vY8dnkvhvTY64Uaf1l7ljRn8ahBvq0ragNv/+o0X9FLvQI/BnmRTujUkvdvqWbIB4SLQ5UNoqy8gv4zr/AscXnaJATUvP3srhlKBXcMPpK2TWMVtpeSrVMknHNMW4qyGFfpjK6t+dew3lXmvKI422muoKEvXnsit32n+nWqEgoksbEEnOKK1pz4OBTCaJUpDA/6HZ5cIRzZdj8eubRn0v3JEBFO6dLKc4zSz7hlVOsPFBcJ5x/bLpKQJxc4rcBUsaCS/bfmjeu5QlYk//9RBP/bv8neCu+v3z+GXiXNeejiHlX8RaJqaUdV1tBe4VrYQSUY3P3Aaj469Z2ouonhMy4NaFyQjgZJovcmI1PFeq6rsXbiIS0pKpJIjT8evSz8O7r/vg18By89pkPTSAxgglB4vucFwMNlPSkZ/mrWzxPkwe99UPVx1ExeHJFgoZpTcdugrtw48LCUk/BBK4OwpFO+3zn6QL5z9N7wFdmcZAxatvt2dogoRpbf3kBxnSKu69uZ37768Z605o3q8djlvejWvindfvNm4HM/Oaw3b85bBVR9VuPj7vHL89Yv+rBswzeByw9KqjmvMI/BBaXtMrJyS5x/2GsA4Y///DjD4I0hMIVRAPh51gqhO5pIUZHQoCj12LuXV3Mqcj0FnfGkYTZuiE9RLj2hhO8c7c/kt1Eay79Uf//wNvv6Nm1O5NDWTTjhkJZ7FIYXcYV6YNOGHNi0YdJ8UeGpwCMYB7z7vMwbYKki2xbgK29DUsm4/8LuPBLQrtzIjLYBAvkB9OgQW9gxyhfKPQHb9/BYlP4gjnvZxl1x3/7dIwLHMHr1pyfRrX0w5d08S5FSCzHae1ijBy/S/U8R8QxM6ff4fGAKIwnnHNN2T8WRjANcEV1vHHgYz/6oeuTKVIR9IGrLOgtXn9pp7wSkx3867dBWlLjs9R+5tJSXf3JS1patvOf8o3l3+Ok0yNCTPdP74p77UIUzjzyAR535sqBe74kyHHHgfpx0SLCAju7YRfH1JZqE6LHE+b4T7r1phOuze9HRh+d5vl+hoiKpZkSTTHVlQ6kFxYakMuTd4adXcfS7pk/VdQB+eHxH2jZtyJ2vz09blp/HIF+PypPDevPvmcvTZwxBcZ0iBhzRmjHvfu65P3EysUmDuknXI890iMldWdcvrrPHkicZ5x/bjmemL2HVpqSR9wPfs0TRH7w42OJQD13cg5c/XFElzS3Dvj7X4QY4oVOLKpZKZ3drw8crNvHn7/tbpyXViN1PTj+EH592SMYrzPkdDbzvwm68//l6rhhbkaIsSVrm6EuOrRZ+PAiJ9zOdst1zXMZnzD7Ww8iQtk0bsl+KFtIdg4/kRwEWk4mCTB3ZRl9yLDf0914c6YRDWvKHC6KZJPdDlNZK2aR980ZMu7lfzsNJpGLgkW144AfJQ/RfftJBjEiYBL/4OH+RTn90Sieeueo43wH64pVlfDGjU10mpyKSk+VImzSom3aU4MqTD+KiXu25/KTqll/9jziAEzr5941KRWnHZky7pa9nLyFuCdh8nwQHxYR8hTCPWThP+7eQEWd3pVv7phzpY1EjPwxJMR6aiv5HHMBPTu/Mz/p1jkSOTIjaWumG/l0ChW8ohO5+74ODzUsEpW6dIi5zmcSOvbwXh+zvz/enTpH4WlQo8TZ2b9+UxaMGea5sVwg0aVCXO889OuuLFQ07+WAa1fM+xw0DDuXJYb33zM0VMqYw8sjR7Zoy7toTfY2Ve1WoUbfFf9avS+hwB2GJal7mol4duGVQ5g5mfvFSM4cd0ITvHdOWPw3tHqis0ZeUMnJwLMRGEK//TDk1wQKqtsyJFQrxJW4fuvjYlB75desUcYKP+31d385572XYHIZREET1IuQ6eJtXHVtcp4j7LgymLCC2vsUlx5dwyfElEcgVrvY35RGepo3qRdoA63t4az6/c1BOfMSSYT0Mo6DIVz1VU+ZOvPjT0O5V1hRxk00FGlX0AaPmYAqjhuD12uc70mqU5H8GoeYyuHvbSLx++7pC4vjRMyd7rP2eCzIN1GiEJ5TCEJHmIlIuIgudb89ZGxEpc/IsFJEyJ62Js153/LNWRO539l0tInOc9P+JSPYHowuQB76f3OKlthFfS8HveshRUwiT3vnm4Fb7ZHXtiai47MQSftq3M2//sg8Tf1E9OGUiw888jI4p1t0oNOKx2U7LIMhptgnbwxgOTFTVzsBEZ7sKItIcGAH0BnoBI0SkmapuVtXu8Q/wBfCCc9iTqnqUk3438MeQctZIBvkM/1AbOKhlY1776cncdJa/pXCDcuugw3nx2tzH3skX8Qntw9vUvgjIDerW4fozutCxRWNf1ldXn9qJtxNW9ssXqTy74zSoW4d3h5/OPeeHW+skG4Sd9B4M9HF+jwUmAzcm5BkAlKvqegARKQcGAk/FM4hIF2B/4B0AVd3kOr4xhe3LkhPybR2RC7pGvBSnm2EnF/7aIVEyuHtbTj9sf5oEcNYzss//nX4IPUua88CkRdVWn3STznE0X4RVGK1VNe5auhLw8pJpCyx1bS9z0twMBZ5R16C8iFwLXA/UA05PJoCIXAVcBdChgz8nJCP71CsuYseuypyft2OLRqzfsoPiooAhszNtkxRwU8aUReEhIhzfqQXHdwoWpqVQSKswRGQC4GVEfIt7Q1VVRDJ9fYYClySU9wDwgIh8H7gVKPM6UFVHA6MBSktLC/j1zYwD92vAlxu3oeq/l+E3BEE2eesXp/L52i05P+8jZT2Z8cUG9gu4tndYctUDjK83bhj5IG3Noqr9ku0TkVUi0kZVV4hIG2C1R7bl7B22AmhHbOgqXkY3oFhVZyQ5zdPAg+nkrK08f80JzPhig+8VyO4achS9PNbOyDXtmjWiXbPcTzQ2b1yPM7qmDgfhRU2Y9P7oNwMyWnchm/xq4KGMn7sy32IYOSLspPdL7G35lwHjPPKMB/qLSDPHiqq/kxbnIlzzGQAi4o5RMQhYGFLOGsuBTRtydrcDPfd5dacu7NnBtZKdUZtoXL+YhlnuYXRpHZsk/0Hvjmlyxvhxn0MY95OTsimSUUCEHbsYBTwrIlcQs3K6AEBESoGrVXWYqq4XkZHAdOeYO+IT4A4XAGcllPsTEekH7AQ2kGQ4yjCM4Aw84gB27PaeX2rVpH7ew8MYhUsohaGq64C+HukVwDDX9hhgTJIyqpmvqOp1YeQyDCM5D10SLGy6EWNIj3bMXb4x32JwwL4NuC5PgULzPztqZEwtcvTOG+f2aMfLs1dwRYDItm6ObrcfFV9sKLi5hXwRX1Tsgp7t8yxJ9OQyzH8qpt5crY2eM0xh1ALaNm3I8DOz4/BW22neuB7jQjj0jbmsJ5+u/pq6WVr9r6YRdcA9o7Cwp7wWcMSB+yadGDeyy74N6nJMDVjHwDCiwBSGYRiG4QtTGDUam8QwDCN3mMKoBdh8q2EYucAUhmEYhuELUxg1GDOrNQwjl5jCqAXUhDhIhmHUfExhZJn4QjaGYRg1HXPcyzKPXtqTyiyNHdVxItjWKza9bxhG9jGFkWWKioSiLA0Z9T28Ndf06cRV37LV5AzDyA+mMGowdYqEGwdaSBDDMHKDjWUYhmEYvjCFYRiGYfjCFIZhGIbhC1MYhmEYhi9CKQwRaS4i5SKy0Pn2jPMsImVOnoUiUuakNRGRWa7PWhG5P+G4ISKizpKvhmEYRh4J28MYDkxU1c7ARGe7CiLSHBgB9AZ6ASNEpJmqblbV7vEPsTXBX3Ad1wS4DpgWUkbDMAwjAsIqjMHAWOf3WOAcjzwDgHJVXa+qG4ByYKA7g4h0AfYH3nEljwTuAraFlNEwDMOIgLAKo7WqrnB+rwRae+RpCyx1bS9z0twMBZ5RjblEi0gPoL2qvhpSPsMwDCMi0jruicgE4ACPXbe4N1RVRSTTGBhDgUuc8xUBfwQu9XOgiFwFXOVsfi0iCzKUoSWwNsNjs4nJFQyTKziFKpvJFYwwcnX0kymtwlDVfsn2icgqEWmjqitEpA2w2iPbcqCPa7sdMNlVRjegWFVnOElNgCOByRJbGegA4CUR+a6qVnjINxoYne5/pENEKlS14CbXTa5gmFzBKVTZTK5g5EKusENSLwFlzu8yYJxHnvFAfxFp5lhR9XfS4lwEPBXfUNWNqtpSVUtUtQSYCngqC8MwDCN3hFUYo4AzRGQh0M/ZRkRKReRhAFVdT2wCe7rzucNJi3MBLoVhGIZhFCahgg+q6jqgr0d6BTDMtT0GGJOkjJShVlW1TxgZAxB6WCtLmFzBMLmCU6iymVzByLpcorbOp2EYhuEDCw1iGIZh+KLWKgwRGSgiC0RkkYh4eaDXF5FnnP3TRKTEte8mJ32BiAzwW2Y25RKRM0RkhojMcb5Pdx0z2SkzHmZl/xzKVSIi37jO/ZDrmGMdeReJyJ/FMXvLkVw/SAg9Uyki3Z19ubhep4jITBHZJSLnJeyrFirHSc/F9fKUS0S6i8gUEflIRGaLyIWuff8Ukc9d16t7ruRy9u12nfslV/pBzj1f5DwD9XIll4iclvB8bRORc5x9ubhe14vIPOdeTRSRjq59WXu+UNVa9wHqAJ8CBwP1gA+Brgl5fgw85PyOOw4CdHXy1wcOcsqp46fMLMt1DHCg8/tIYLnrmMlAaZ6uVwkwN0m57wPHAQK8DpyZK7kS8hwFfJrj61UCHA08BpznSm8OfOZ8N3N+N8vh9UomVxegs/P7QGAF0NTZ/qc7by6vl7Pv6yTlPgsMdX4/BFyTS7kS7ul6oFEOr9dprvNdw973MWvPl6rW2h5GL2CRqn6mqjuAp4mFMXHjDmvyPNDX0biDgadVdbuqfg4scsrzU2bW5FLVD1T1Syf9I6ChiNQPeP7I5UpWoMT8cvZV1akae1ofwzt0TC7kusg5NirSyqWqi1V1NlCZcKxnqJxcXa9kcqnqJ6q60Pn9JTGfqlYBzx+5XMlw7vHpxO45JA9NlAu5zgNeV9WtAc8fRq5JrvNNJebfBtl9vmqtwvATjmRPHlXdBWwEWqQ41k+Z2ZTLzRBgpqpud6U96nR/b8ugqxlWroNE5AMReVtETnblX5amzGzLFedCqptuZ/t6BT02V9crLSLSi1jL9lNX8u+c4Y/7MmiohJWrgYhUiMjU+LAPsXv8lXPPMykzCrniDKX685XL63UFsR5DqmOjeL5qrcKotYjIEcSCMv7IlfwDVT0KONn5XJJDkVYAHVT1GOB64EkR2TeH50+JiPQGtqrqXFdyPq9XQeO0RB8HLlPVeKv6JuAwoCexoY4bcyxWR415MH8fuF9EOuX4/ElxrtdRVHVGztn1EpGLgVLgnmydw01tVRjLgfau7XZOmmceESkG9gPWpTjWT5nZlAsRaQf8B/ihqu5p/anqcud7M/AksS5tTuRyhu7WOeefQaxV2sXJ3851fM6vl0O11l+OrlfQY3N1vZLiKPpXgVtUdWo8XVVXaIztwKPk9nq579dnxOafjiF2j5s69zxwmVHI5XAB8B9V3emSNyfXS0T6EYvp913XaEM2n69aO+ldTGyy5yD2ThodkZDnWqpOlj7r/D6CqpPenxGbhEpbZpblaurkP9ejzJbO77rExnSvzqFcrYA6zu+DnYewuXpPsp2VK7mc7SJHnoNzfb1cef9J9Unvz4lNSDZzfufseqWQqx6xNW1+5pG3jfMtwP3AqBzK1Qyo7/xuCSzEmQAGnqPqpPePcyWXK30qcFqurxcxpfkpjqFCLp4vVa2dCsO5OGcBnzgX9RYn7Q5i2higgfPALXIupLtSucU5bgEuSwKvMnMlF3ArsAWY5frsDzQGZgCziU2G/wmnAs+RXEOc884CZgJnu8osBeY6Zf4Vx1E0h/exDzA1obxcXa+exMaJtxBrDX/kOvZyR95FxIZ+cnm9POUCLgZ2Jjxf3Z19bwFzHNmeAPbJoVwnOOf+0Pm+wlXmwc49X+Q8A/VzfB9LiDVIihLKzMX1mgCsct2rl3LxfJmnt2EYhuGL2jqHYRiGYUSMKQzDMAzDF6YwDMMwDF+YwjAMwzB8YQrDMAzD8IUpDMMwDMMXpjAMwzAMX5jCMAzDMHzx/0Fs79DaozsRAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# now try to guess the parameters in x from the mixed signal:\n",
"x0 = [10e-3, 0.0] # initial guess\n",
"x1 = minimize(err, x0, (t, v)).x\n",
"print(\"guessed:\", x1)\n",
"\n",
"v1 = v - line_noise(x1, t)\n",
"plt.plot(t, v1);"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"language": "python",
"name": "python2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.15"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment