ifft去噪还原.ipynb · 杨皓宇/逆快速傅里叶变换ifft去噪还原

加入 Gitee

与超过 1200万开发者一起发现、参与优秀开源项目，私有仓库也完全免费：）

文件

该仓库未声明开源许可证文件（LICENSE），使用请关注具体项目描述及其代码上游依赖。

克隆/下载

ifft去噪还原.ipynb 27.83 KB

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "fb8b43f4-cb4c-4554-a8e6-1fcaec7d491c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABHcAAAD4CAYAAABmBgoFAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAUBklEQVR4nO3df4xlZ3kf8O/D7kI3BLqJPEX22s5CRSwhaFk6AiErKCUJNjQCl/6InTaiSaQlEqlArdzg9o/SSBU0bqLkrzQu0FIVDLTYBKEEQ0VolIofnrVdjG2cENfArl08FG2AdkXw5ukfe9dar2d27njn7r3vnc9HWu3cd889fp55/d475zvnnlPdHQAAAADG9Ix5FwAAAADA0yfcAQAAABiYcAcAAABgYMIdAAAAgIEJdwAAAAAGtncWO73kkkv60KFDs9g1AAAAwK509OjRb3b3yrnjMwl3Dh06lLW1tVnsGgAAAGBXqqqvbjTuY1kAAAAAAxPuAAAAAAxMuAMAAAAwMOEOAAAAwMCEOwAAAAADm+puWVV1IMm7k7w4SSf5he7+7AzrYkF89O7jufmOB/PIiZO57MD+3HjNVbnu8MF5lwUAAABMTHsr9N9K8onu/rtV9cwkPzDDmlgQH737eG667d6c/P6pJMnxEydz0233JomABwAAABbElh/LqqrnJnlVkvckSXf/eXefmHFdLICb73jwiWDnjJPfP5Wb73hwThUBAAAA55rmmjsvSLKe5D9U1d1V9e6qeva5G1XVkapaq6q19fX1HS+Ui++REye3NQ4AAABcfNOEO3uTvCzJb3f34ST/N8nbz92ou2/p7tXuXl1ZWdnhMpmHyw7s39Y4AAAAcPFNE+4cS3Ksuz8/efxfczrsYcndeM1V2b9vz5PG9u/bkxuvuWpOFQEAAADn2jLc6e7/neTrVXXmiP4nktw/06pYCNcdPph3vvEleeae0/+bHDywP+9840tcTBkAAAAWyLR3y/rHSd4/uVPWQ0l+fnYlsUiuO3wwt37ha0mSD735lXOuBgAAADjXVOFOd9+TZHW2pQAAAACwXdNccwcAAACABSXcAQAAABiYcAcAAABgYMIdAAAAgIEJdwAAAAAGJtwBAAAAGJhwBwAAAGBgwh0AAACAgQl3AAAAAAYm3AEAAAAYmHAHAAAAYGDCHQAAAICBCXcAAAAABibcAQAAABiYcAcAAABgYMIdAAAAgIEJdwAAAAAGJtwBAAAAGJhwBwAAAGBgwh0AAACAgQl3AAAAAAYm3AEAAAAYmHAHAAAAYGDCHQAAAICB7Z1mo6p6OMl3kpxK8nh3r86yKAAAAACmM1W4M/E3u/ubM6sEAAAAgG3zsSwAAACAgU0b7nSST1bV0ao6stEGVXWkqtaqam19fX3nKgQAAABgU9OGO1d398uSvDbJW6rqVedu0N23dPdqd6+urKzsaJEAAAAAbGyqcKe7H5n8/ViS25O8fJZFAQAAADCdLcOdqnp2VT3nzNdJXpPkS7MuDAAAAICtTXO3rOclub2qzmz/ge7+xEyrAgAAAGAqW4Y73f1Qkr9+EWoBAAAAYJvcCh0AAABgYMIdAAAAgIEJdwAAAAAGJtwBAAAAGJhwBwAAAGBgwh0AAACAgQl3AAAAAAYm3AEAAAAYmHAHAAAAYGDCHQAAAICBCXcAAAAABibcAQAAABiYcAcAAABgYMIdAAAAgIEJdwAAAAAGJtwBAAAAGJhwBwAAAGBgwh0AAACAgQl3AAAAAAYm3AEAAAAYmHAHAAAAYGDCHQAAAICBCXcAAAAABibcAQAAABjY1OFOVe2pqrur6uOzLAgAAACA6W3nzJ23JnlgVoUAAAAAsH1ThTtVdXmSv5Xk3bMtBwAAAIDtmPbMnd9M8s+S/MVmG1TVkapaq6q19fX1nagNAAAAgC1sGe5U1U8neay7j55vu+6+pbtXu3t1ZWVlxwoEAAAAYHPTnLlzdZLXV9XDST6Y5NVV9Z9nWhUAAAAAU9ky3Onum7r78u4+lOT6JJ/u7n8488oAAAAA2NJ27pYFAAAAwILZu52Nu/szST4zk0oAAAAA2DZn7gAAAAAMTLgDAAAAMDDhDgAAAMDAhDsAAAAAAxPuAAAAAAxMuAMAAAAwMOEOAAAAwMCEOwAAAAADE+4AAAAADEy4AwAAADAw4Q4AAADAwIQ7AAAAAAMT7gAAAAAMTLgDAAAAMDDhDgAAAMDAhDsAAAAAAxPuAAAAAAxMuAMAAAAwMOEOAAAAwMCEOwAAAAADE+4AAAAADEy4AwAAADAw4Q4AAADAwIQ7AAAAAAPbMtypqr9UVV+oqv9ZVfdV1b+6GIUBAAAAsLW9U2zzvSSv7u7vVtW+JH9UVb/f3Z+bcW0AAAAAbGHLcKe7O8l3Jw/3Tf70LIsCAAAAYDpTXXOnqvZU1T1JHkvyqe7+/AbbHKmqtapaW19f3+EyAQAAANjIVOFOd5/q7pcmuTzJy6vqxRtsc0t3r3b36srKyg6XCQAAAMBGtnW3rO4+keQzSa6dRTEAAAAAbM80d8taqaoDk6/3J/nJJF+ecV0AAAAATGGau2VdmuR9VbUnp8OgD3f3x2dbFgAAAADTmOZuWV9Mcvgi1AIAAADANm3rmjsAAAAALBbhDgAAAMDAhDsAAAAAAxPuAAAAAAxMuAMAAAAwMOEOAAAAwMCEOwAAAAADE+4AAAAADEy4AwAAADAw4Q4AAADAwIQ7AAAAAAMT7gAAAAAMTLgDAAAAMDDhDgAAAMDAhDsAAAAAAxPuAAAAAAxMuAMAAAAwMOEOAAAAwMCEOwAAAAADE+4AAAAADEy4AwAAADAw4Q4AAADAwIQ7AAAAAAMT7gAAAAAMbMtwp6quqKo/qKoHquq+qnrrxSgMAAAAgK3tnWKbx5P80+6+q6qek+RoVX2qu++fcW0AAAAAbGHLM3e6+9Huvmvy9XeSPJDk4KwLAwAAAGBr27rmTlUdSnI4yec3+LcjVbVWVWvr6+s7VB4AAAAA5zN1uFNVP5jkI0ne1t3fPvffu/uW7l7t7tWVlZWdrBEAAACATUwV7lTVvpwOdt7f3bfNtiQAAAAApjXN3bIqyXuSPNDdvzH7kgAAAACY1jRn7lyd5OeSvLqq7pn8ed2M6wIAAABgClveCr27/yhJXYRaAAAAANimbd0tCwAAAIDFItwBAAAAGJhwBwAAAGBgwh0AAACAgQl3AAAAAAYm3AEAAAAYmHAHAAAAYGDCHQAAAICBCXcAAAAABibcAQAAABiYcAcAAABgYMIdAAAAgIEJdwAAAAAGJtwBAAAAGJhwBwAAAGBgwh0AAACAgQl3AAAAAAYm3AEAAAAYmHAHAAAAYGDCHQAAAICBCXcAAAAABibcAQAAABiYcAcAAABgYMIdAAAAgIFtGe5U1Xur6rGq+tLFKAgAAACA6U1z5s5/THLtjOsAAAAA4GnYMtzp7j9M8q2LUAsAAAAA27Rj19ypqiNVtVZVa+vr6zu1WwAAAADOY8fCne6+pbtXu3t1ZWVlp3YLAAAAwHm4WxYAAADAwIQ7AAAAAAOb5lbotyb5bJKrqupYVf3i7MsCAAAAYBp7t9qgu2+4GIUAAAAAsH0+lgUAAAAwMOEOAAAAwMCEOwAAAAADE+4AAAAADEy4AwAAADAw4Q4AAADAwIQ7AAAAAAMT7gAAAAAMTLgDAAAAMDDhDgAAAMDAhDsAAAAAAxPuAAAAAAxMuAMAAAAwMOEOAAAAwMCEOwAAAAADE+4AAAAADEy4AwAAADAw4Q4AAADAwIQ7AAAAAAMT7gAAAAAMTLgDAAAAMDDhDgAAAMDAhDsAAAAAAxPuAAAAAAxs7zQbVdW1SX4ryZ4k7+7ud820qjn76N3Hc/MdD+aREydz2YH9ufGaq5JkqrHrDh+8oOcv0th1hw/u+Pdnnr2MWLdexhlbpl6WvT+9LOaYXsYZW6b+9LKYY3pZzLFl6mXZ+9PL1vtcRtXd59+gak+SP07yU0mOJbkzyQ3dff9mz1ldXe21tbWdrPOi+ejdx3PTbffm5PdPPTG27xmVVPL9U33esf379uTv/I2D+cjR40/r+Ys0tn/fnrzzjS/JdYcP5md+57NJkg+9+ZUX9P2ZZy/LNC96WbyxZepl2fvTy2KO6WWcsWXqTy+LOaaXxRxbpl6WvT+9bL3PM8e5o6qqo929+pTxKcKdVyZ5R3dfM3l8U5J09zs3e87I4c7V7/p0jp84mTd/8Xfzgj87vu3nV1W2+p6O4ll79+TwlQdy/6PfTpK86NLn5u6vncj3Hj+1xTMXzzLNi14W0zL1spFl6k8vi0kv41im/vSymPSymJapl40sU396eaqH/vLB/M5fe0OS5OCB/fkfb3/1Be9zXjYLd6b5WNbBJF8/6/GxJK/Y4D9wJMmRJLnyyiufZpnz98iJkxf0/GVZREmeCHFedOlznzI2mmWaF70spmXqZSPL1J9eFpNexrFM/ellMellMS1TLxtZpv70cn4Xesy/qKYJd2qDsad8h7v7liS3JKfP3LnAuubmsgP7c/zEySdSve3aU5VTS7KYDh7Yn585J9H82cmZTaNZpnnRy2Japl42skz96WUx6WUcy9SfXhaTXhbTMvWykWXqTy/nd9mB/Tu6v0Uxzd2yjiW54qzHlyd5ZDblzN+N11yV/fv2PGls3zMq+/bUlmP79+3JDa+44mk/f5HG9u/b88QFqM52Id+fefayTPOil8UbW6Zelr0/vSzmmF7GGVum/vSymGN6WcyxZepl2fvTy9b73Og4dxlMc+bOnUleWFXPT3I8yfVJfnamVc3RmQsrXciVt1d/5IcX5uriF9rLTn9/5tnLMs2LXhZvbJl6Wfb+9LKYY3oZZ2yZ+tPLYo7pZTHHlqmXZe9PL1vvcxlteUHlJKmq1yX5zZy+Ffp7u/tfn2/7kS+oDAAAALCILuSCyunu30vyezteFQAAAAAXZJpr7gAAAACwoIQ7AAAAAAMT7gAAAAAMTLgDAAAAMLCp7pa17Z1WrSf56o7v+OK7JMk3510Ec2P+dy9zv7uZ/93L3O9u5n/3Mve7l7nf3Uad/x/p7pVzB2cS7iyLqlrb6BZj7A7mf/cy97ub+d+9zP3uZv53L3O/e5n73W3Z5t/HsgAAAAAGJtwBAAAAGJhw5/xumXcBzJX5373M/e5m/ncvc7+7mf/dy9zvXuZ+d1uq+XfNHQAAAICBOXMHAAAAYGDCHQAAAICBCXc2UVXXVtWDVfWVqnr7vOthdqrqiqr6g6p6oKruq6q3TsbfUVXHq+qeyZ/XzbtWZqOqHq6qeyfzvDYZ++Gq+lRV/cnk7x+ad53srKq66qz1fU9Vfbuq3mbtL6+qem9VPVZVXzprbNO1XlU3TX4OeLCqrplP1eyETeb+5qr6clV9sapur6oDk/FDVXXyrNeAfze3wtkRm8z/pq/11v7y2GTuP3TWvD9cVfdMxq39JXKeY7ylfd93zZ0NVNWeJH+c5KeSHEtyZ5Ibuvv+uRbGTFTVpUku7e67quo5SY4muS7J30/y3e7+t/Osj9mrqoeTrHb3N88a+7Uk3+rud00C3h/q7l+ZV43M1uR1/3iSVyT5+Vj7S6mqXpXku0n+U3e/eDK24VqvqhcluTXJy5NcluS/JfnR7j41p/K5AJvM/WuSfLq7H6+qf5Mkk7k/lOTjZ7ZjfJvM/zuywWu9tb9cNpr7c/7915P8WXf/qrW/XM5zjPePsqTv+87c2djLk3ylux/q7j9P8sEkb5hzTcxIdz/a3XdNvv5OkgeSHJxvVSyANyR53+Tr9+X0mwHL6yeS/Gl3f3XehTA73f2HSb51zvBma/0NST7Y3d/r7v+V5Cs5/fMBA9po7rv7k939+OTh55JcftEL46LYZO1vxtpfIueb+6qqnP5l7q0XtSguivMc4y3t+75wZ2MHk3z9rMfH4mB/V5gk9oeTfH4y9MuT07Xf62M5S62TfLKqjlbVkcnY87r70eT0m0OSvzK36rgYrs+Tf7iz9nePzda6nwV2l19I8vtnPX5+Vd1dVf+9qn5sXkUxcxu91lv7u8ePJflGd//JWWPW/hI65xhvad/3hTsbqw3GfH5tyVXVDyb5SJK3dfe3k/x2kr+a5KVJHk3y6/Orjhm7urtfluS1Sd4yOYWXXaKqnpnk9Un+y2TI2ifxs8CuUVX/IsnjSd4/GXo0yZXdfTjJP0nygap67rzqY2Y2e6239nePG/LkX+xY+0tog2O8TTfdYGyotS/c2dixJFec9fjyJI/MqRYugqral9OL/v3dfVuSdPc3uvtUd/9Fkn+fwU7LY3rd/cjk78eS3J7Tc/2NyWd1z3xm97H5VciMvTbJXd39jcTa34U2W+t+FtgFqupNSX46yT/oyYUoJ6fk/5/J10eT/GmSH51flczCeV7rrf1doKr2Jnljkg+dGbP2l89Gx3hZ4vd94c7G7kzywqp6/uQ3utcn+dica2JGJp+3fU+SB7r7N84av/Sszf52ki+d+1zGV1XPnlxkLVX17CSvyem5/liSN002e1OS351PhVwET/rNnbW/62y21j+W5PqqelZVPT/JC5N8YQ71MSNVdW2SX0ny+u7+f2eNr0wusp6qekFOz/1D86mSWTnPa721vzv8ZJIvd/exMwPW/nLZ7BgvS/y+v3feBSyiyV0TfjnJHUn2JHlvd98357KYnauT/FySe8/cCjHJP09yQ1W9NKdPx3s4yZvnURwz97wkt59+/c/eJB/o7k9U1Z1JPlxVv5jka0n+3hxrZEaq6gdy+s6IZ6/vX7P2l1NV3Zrkx5NcUlXHkvzLJO/KBmu9u++rqg8nuT+nP7LzlpHumMGTbTL3NyV5VpJPTd4DPtfdv5TkVUl+taoeT3IqyS9197QX42UBbTL/P77Ra721v1w2mvvufk+eeq29xNpfNpsd4y3t+75boQMAAAAMzMeyAAAAAAYm3AEAAAAYmHAHAAAAYGDCHQAAAICBCXcAAAAABibcAQAAABiYcAcAAABgYP8fojjP6mvtnQQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 1440x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "def denoise(arr, thresh):\n",
    "    mask=(arr > thresh) \n",
    "    arr = mask*fd\n",
    "    return arr\n",
    "def discard(arr, thresh):\n",
    "    arr[thresh:] = 0\n",
    "    return arr  \n",
    "srate=200\n",
    "fd=np.random.rand(srate) #产生一个列表，其中更包含200个0到1的随机数,模拟200个噪音信号\n",
    "\n",
    "fd[13]=6  # 模拟在13hz有一个振幅为6的信号\n",
    "fd[100]=9 # 模拟在100hz有一个振幅为9的信号\n",
    "fd[188]=3 # 模拟在188hz有一个振幅为3的信号\n",
    "\n",
    "\n",
    "fd=denoise(fd,1)\n",
    "fd=discard(fd,50)\n",
    "\n",
    "plt.figure(figsize=(20,4))\n",
    "plt.stem(fd)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "349cf0ca-b4d1-4663-847a-cc99c7b566c7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "5bcb0cb1b7434ec79da1ad57f2354b8c",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(FloatSlider(value=0.5, continuous_update=False, description='频率（Hz）', max=2.0, min=0.1, …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "<function __main__.draw_circle(f, t)>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import ipywidgets\n",
    "import math\n",
    "\n",
    "# 定义绘制圆形轨迹的函数\n",
    "def draw_circle(f, t):\n",
    "    srate = 10000  # 设置采样率\n",
    "    t = np.linspace(0, t, srate)  # 根据采样率在时间范围内生成等间隔的点\n",
    "    \n",
    "    # 根据频率 f 和时间 t，计算复平面上的圆形轨迹点\n",
    "    p = np.exp(-1j * 2 * math.pi * f * t)\n",
    "    # 绘制圆形轨迹点，横坐标为实数部分，纵坐标为虚数部分\n",
    "    plt.plot(np.real(p), np.imag(p), \"-\",color = 'blue')\n",
    "    rp = np.real(p)\n",
    "    ip = np.imag(p)\n",
    "    plt.plot(rp[-1],ip[-1],\"o\",color = 'yellow')\n",
    "    \n",
    "    plt.axis('square')  # 设置坐标轴比例相等\n",
    "\n",
    "    plt.grid()  # 显示网格线\n",
    "    plt.show()  # 显示图像\n",
    "\n",
    "# 创建一个滑动条，用于调节画圆的频率，范围从 0.1 Hz 到 2 Hz，步长为 0.01 Hz\n",
    "f = ipywidgets.FloatSlider(\n",
    "    min=0.1, max=2, step=0.01, value=0.5, \n",
    "    description='频率（Hz）',\n",
    "    continuous_update=False)\n",
    "\n",
    "# 创建一个滑动条，用于调节画圆的总时间，范围从 0 秒到 10 秒，步长为 0.1 秒\n",
    "t = ipywidgets.FloatSlider(\n",
    "    min=0, max=10, step=0.1, value=0, \n",
    "    description='时间（秒）',\n",
    "    continuous_update=False)\n",
    "\n",
    "# 使用 ipywidgets 的交互功能，将滑动条与绘制圆形轨迹的函数关联\n",
    "ipywidgets.interact(draw_circle, t=t, f=f)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "d483a993-0f26-4954-b7cf-8d0c1c919750",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "78078e0eaedd45d1a154f9ef0c586d5b",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(FloatSlider(value=4.0, continuous_update=False, description='波形1 频率', min=1.0, step=1.0)…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "<function __main__.draw_cos(f1, f2, amp1, amp2)>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import math\n",
    "import matplotlib.pyplot as plt\n",
    "import ipywidgets\n",
    "def draw_cos(f1,f2,amp1,amp2):\n",
    "    ####### 把这里函数完整 ##########\n",
    "    # 1. 设置显示时间为1秒，采样率为2000\n",
    "    srate = 2000\n",
    "    # 生成时间序列\n",
    "    t = np.linspace(0, 1, srate)  #1s内，采样2000次\n",
    "    # 2. 定义两个cos函数，f1, amp1为第一个cos函数的频率和振幅；f2, amp2为第二个cos函数的频率和振幅\n",
    "    s1 = amp1 * np.cos(2 * math.pi * f1 * t)         # 振幅为f1\n",
    "    s2 = amp2 * np.cos(2 * math.pi * f2 * t)         # 振幅为f2\n",
    "    # 3. 以一个合适的画板大小 figsize 显示图像\n",
    "    plt.figure(figsize=(20, 4))  # 创建自定义图像,figsize:指定figure的宽和高\n",
    "    # 4. plot两个图像，横轴为时间，纵轴为对应的函数值\n",
    "   # plt.plot(t, s1,color='blue')\n",
    "   # plt.plot(t, s2,color='red', linewidth=2) \n",
    "    plt.plot(t, s1)\n",
    "    plt.plot(t, s2,linewidth=2) \n",
    "    # 5. 用xlim方法设置显示的范围为：频率最低的那个函数的2个周期\n",
    "    T = 1/min(f1,f2)\n",
    "    plt.xlim([0, 2*T])\n",
    "\n",
    "f1 = ipywidgets.FloatSlider(min=1, max=100, step=1, value=4, continuous_update=False,description=\"波形1 频率\")\n",
    "f2 = ipywidgets.FloatSlider(min=1, max=200, step=1, value=2, continuous_update=False,description=\"波形2 频率\")\n",
    "amp1 = ipywidgets.FloatSlider(min=0, max=10, step=0.1, value=1, continuous_update=False,description=\"波形1 振幅\")\n",
    "amp2 = ipywidgets.FloatSlider(min=0, max=10, step=0.1, value=1, continuous_update=False,description=\"波形2 振幅\")\n",
    "ipywidgets.interact(draw_cos, f1=f1,f2=f2,amp1=amp1,amp2=amp2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "e332b14e-99bd-4ea6-9458-3b2061dc5882",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0IAAAEvCAYAAABoq3HZAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAhgklEQVR4nO3df5BlZZnY8e9Dz2BGYDO4Njg0rGA521mUyKwUSihTF10FJu4ykk0CSVbctWp0I1VLlSEOscrVtRLYUO4mWyhkdqXArOFHRRhRRgci3qBVIAIzMLDQy8iyMj1TTKE20tCBmebNH32aXJp7+8c99+d5v5+qW33v+XHP+5x73vPcp++5742UEpIkSZKUk8P63QBJkiRJ6jULIUmSJEnZsRCSJEmSlB0LIUmSJEnZsRCSJEmSlB0LIUmSJEnZWdXvBjSzdu3a9Pa3v73fzeiJF154gSOOOKLfzegJY60mY62uBx544NmU0mi/2zGIzFPVZKzVlFOskFe8ZfPUQBZCxx57LPfff3+/m9ET9XqdWq3W72b0hLFWk7FWV0T8fb/bMKjMU9VkrNWUU6yQV7xl85SXxkmSJEnKjoWQJEmSpOxYCEmSJEnKjoWQJEmSpOxYCEmSJEnKjoWQJEmSpOxYCEmSJEnKzpKFUEScEBHfj4jHIuLRiPijYvqbIuLOiHii+Ht0i/XPiYiJiNgTEVs6HYAkKW/mKUlSO5bzidAh4NMppd8A3gt8KiJOBrYA30sprQe+Vzx+jYgYAb4MnAucDFxYrDvwtu2c5Mwr7uKkLbdz5hV3sW3nZL+bJKmiPN+UlmWekiSVs2QhlFLan1J6sLj/PPAYMAacB1xfLHY9sKnJ6qcDe1JKT6aUXgZuLNYbaNt2TnLZLbuZnJohAZNTM1x2y27fnEjqOM835eWYpyRJ5a3oO0IRcSKwAfgRcGxKaT/MJSHgmCarjAFPNzzeW0wbaFfumGDm4Oxrps0cnOXKHRN9apGkqvJ801m55ClJUnmrlrtgRBwJfAO4JKX0y4hY1mpNpqUWz78Z2AwwOjpKvV5fbtM6bnJqpuX0Trdrenq6r7H2krFWk7GW08vzTdXllKd6yT5eTcZaXbnFW8ayCqGIWM1ccvl6SumWYvIzEbEupbQ/ItYBB5qsuhc4oeHx8cC+ZttIKW0FtgKMj4+nWq22vAi6YOzeu5q+ORlbu4ZOt6ter3f8OQeVsVaTsZbTy/NNleWWp3rJPl5NxlpducVbxnJGjQvgq8BjKaU/a5h1G3BRcf8i4JtNVv8xsD4iToqIw4ELivUG2qVnj7Nm9chrpq1ZPcKlZ4/3qUWSqsrzTXk55ilJUnnL+Y7QmcDvAe+PiF3FbSNwBfDBiHgC+GDxmIg4LiK2A6SUDgEXAzuY+/LqzSmlR7sQR0dt2jDG5eefwuEjc7tnbO0aLj//FDZt8LJxSZ3l+aYjsstTkqTylrw0LqX0Q5pfQw3wgSbL7wM2NjzeDmxvt4H9smnDGDfc91MAbvrEGX1ujaQq83xTTq55SpJUzopGjZMkSZKkKrAQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpQdCyFJkiRJ2bEQkiRJkpSdVUstEBHXAh8GDqSU3llMuwkYLxZZC0yllE5tsu5TwPPALHAopXRaR1otSVLBPCVJaseShRBwHXAV8LX5CSmlfzV/PyK+BDy3yPpnpZSebbeBkiQt4TrMU5KkFVqyEEop3R0RJzabFxEB/Evg/R1ulyRJy2KekiS1o+x3hN4HPJNSeqLF/ATcEREPRMTmktuSJGmlzFOSpKaWc2ncYi4Eblhk/pkppX0RcQxwZ0Q8nlK6u9mCRQLaDDA6Okq9Xi/ZtPKmpmYAutqW6enpgYi1F4y1moy1M3pxvslUpfNUL9jHq8lYqyu3eMtouxCKiFXA+cC7Wy2TUtpX/D0QEbcCpwNNE0xKaSuwFWB8fDzVarV2m9YxV0/cA0CtdkbXtlGv1xmEWHvBWKvJWDujF+eb3OSQp3rBPl5NxlpducVbRplL434LeDyltLfZzIg4IiKOmr8PfAh4pMT2JElaCfOUJKmlJQuhiLgBuAcYj4i9EfHxYtYFLLjcICKOi4jtxcNjgR9GxEPAfcDtKaXvdq7pkiSZpyRJ7VnOqHEXtpj+sSbT9gEbi/tPAu8q2T5JkhZlnpIktaPsqHGSJEmSNHQshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnYshCRJkiRlx0JIkiRJUnaWLIQi4tqIOBARjzRM+3xETEbEruK2scW650TERETsiYgtnWy4JElgnpIktWc5nwhdB5zTZPqfp5ROLW7bF86MiBHgy8C5wMnAhRFxcpnGSpLUxHWYpyRJK7RkIZRSuhv4eRvPfTqwJ6X0ZErpZeBG4Lw2nkeSpJbMU5KkdpT5jtDFEfFwcUnC0U3mjwFPNzzeW0yTJKkXzFOSpJZWtbne1cAXgVT8/RLwBwuWiSbrpVZPGBGbgc0Ao6Oj1Ov1NpvWOVNTMwBdbcv09PRAxNoLxlpNxtoZvTjfZCaLPNUL9vFqMtbqyi3eMtoqhFJKz8zfj4i/BL7dZLG9wAkNj48H9i3ynFuBrQDj4+OpVqu107SOunriHgBqtTO6to16vc4gxNoLxlpNxtoZvTjf5CSXPNUL9vFqMtbqyi3eMtq6NC4i1jU8/AjwSJPFfgysj4iTIuJw4ALgtna2J0nSSpinJElLWfIToYi4AagBb46IvcAfA7WIOJW5SwieAj5RLHsc8FcppY0ppUMRcTGwAxgBrk0pPdqNICRJ+TJPSZLasWQhlFK6sMnkr7ZYdh+wseHxduB1Q5ZKktQp5ilJUjvKjBonSZIkSUPJQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGXHQkiSJElSdiyEJEmSJGVnyUIoIq6NiAMR8UjDtCsj4vGIeDgibo2ItS3WfSoidkfEroi4v4PtliQJME9JktqznE+ErgPOWTDtTuCdKaV/DPwtcNki65+VUjo1pXRae02UJGlR12GekiSt0JKFUErpbuDnC6bdkVI6VDy8Fzi+C22TJGlJ5ilJUjsipbT0QhEnAt9OKb2zybxvATellP66yby/A34BJOC/p5S2LrKNzcBmgNHR0XfffPPNy42hay7/0QwAl71nTde2MT09zZFHHtm15x8kxlpNxtoZvTjfrNRZZ531wLB8SpJrnuoF+3g1GWt15RRv2Ty1qszGI+KzwCHg6y0WOTOltC8ijgHujIjHi//cvU6RfLYCjI+Pp1qtVqZpHXH1xD0A1GpndG0b9XqdQYi1F4y1moy1M3pxvslR1fNUL9jHq8lYqyu3eMtoe9S4iLgI+DDwb1KLj5VSSvuKvweAW4HT292eJEkrYZ6SJC2mrUIoIs4BPgP8TkrpxRbLHBERR83fBz4EPNJsWUmSOsk8JUlaynKGz74BuAcYj4i9EfFx4CrgKOYuI9gVEdcUyx4XEduLVY8FfhgRDwH3AbenlL7blSgkSdkyT0mS2rHkd4RSShc2mfzVFsvuAzYW958E3lWqdZIkLcE8JUlqR9vfEZIkSZKkYWUhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7FkKSJEmSsmMhJEmSJCk7SxZCEXFtRByIiEcapr0pIu6MiCeKv0e3WPeciJiIiD0RsaWTDZckCcxTkqT2LOcToeuAcxZM2wJ8L6W0Hvhe8fg1ImIE+DJwLnAycGFEnFyqtRWybeckZ15xFx/77gucecVdbNs52e8mSWrB/jrwrsM81XEe99Jwsc+u3JKFUErpbuDnCyafB1xf3L8e2NRk1dOBPSmlJ1NKLwM3Futlb9vOSS67ZTeTUzMATE7NcNktuz1gpQFkfx185qnO87iXhot9tj3tfkfo2JTSfoDi7zFNlhkDnm54vLeYlr0rd0wwc3D2NdNmDs5y5Y6JPrVIUiv216FlnirB414aLvbZ9qzq4nNHk2mp5cIRm4HNAKOjo9Tr9S41a/mmiqq6022Zr9abTR+EuLtlenq60vE1Mtbq6FV/7db5Rosa+jzVLeap6jPWasm1z5bVbiH0TESsSyntj4h1wIEmy+wFTmh4fDywr9UTppS2AlsBxsfHU61Wa7NpnXP1xD0A1GpndPR5x+69q+kBO7Z2DYMQd7fU6/VKx9fIWKujV/21W+ebjGWRp7rFPFV9xlotufbZstq9NO424KLi/kXAN5ss82NgfUScFBGHAxcU62Xv0rPHWbN65DXT1qwe4dKzx/vUIkmt2F+HlnmqBI97abjYZ9uznOGzbwDuAcYjYm9EfBy4AvhgRDwBfLB4TEQcFxHbAVJKh4CLgR3AY8DNKaVHuxPGcNm0YYzLzz+Fw0fmdv/Y2jVcfv4pbNrgpenSoLG/Dj7zVOd53EvDxT7bniUvjUspXdhi1geaLLsP2NjweDuwve3WVdimDWPccN9PmZqaYsdn3t/v5khahP11sJmnusPjXhou9tmVa/fSOEmSJEkaWhZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpO20XQhExHhG7Gm6/jIhLFixTi4jnGpb5XOkWS5K0DOYpSdJiVrW7YkppAjgVICJGgEng1iaL/iCl9OF2tyNJUjvMU5KkxXTq0rgPAD9JKf19h55PkqROMk9Jkl6jU4XQBcANLeadEREPRcR3IuIdHdqeJEkrYZ6SJL1G25fGzYuIw4HfAS5rMvtB4K0ppemI2AhsA9a3eJ7NwGaA0dFR6vV62aaVNjU1A9C1tkxNzTA7OzsQsfbC9PS0sVZQLrF2u792+3yTsyrnqW4zT1WXsVZTbn22rNKFEHAu8GBK6ZmFM1JKv2y4vz0ivhIRb04pPdtk2a3AVoDx8fFUq9U60LRyrp64B4Ba7YyuPf/U1BSDEGsv1Ot1Y62gXGLtdn/t9vkmc5XNU91mnqouY62m3PpsWZ24NO5CWlxuEBFviYgo7p9ebO9nHdimJEnLZZ6SJL1OqU+EIuKNwAeBTzRM+yRASuka4HeBP4yIQ8AMcEFKKZXZpiRJy2WekiS1UqoQSim9CPzqgmnXNNy/CriqzDYkSWqXeUqS1EqnRo2TJEmSpKFhISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpOxZCkiRJkrJjISRJkiQpO6UKoYh4KiJ2R8SuiLi/yfyIiL+IiD0R8XBE/GaZ7UmStBLmKUlSK6s68BxnpZSebTHvXGB9cXsPcHXxV5KkXjFPSZJep9uXxp0HfC3NuRdYGxHrurxNSZKWyzwlSZkqWwgl4I6IeCAiNjeZPwY83fB4bzFNkqReME9Jkpoqe2ncmSmlfRFxDHBnRDyeUrq7YX40WSc1e6IiQW0GGB0dpV6vl2xaeVNTMwBda8vU1Ayzs7MDEWsvTE9PG2sF5RJrt/trt883Gat0nuo281R1GWs15dZnyypVCKWU9hV/D0TErcDpQGOC2Quc0PD4eGBfi+faCmwFGB8fT7VarUzTOuLqiXsAqNXO6NrzT01NMQix9kK9XjfWCsol1m73126fb3JV9TzVbeap6jLWasqtz5bV9qVxEXFERBw1fx/4EPDIgsVuAz5ajMrzXuC5lNL+tlsrSdIymackSYsp84nQscCtETH/PP8zpfTdiPgkQErpGmA7sBHYA7wI/H655kqStGzmKUlSS20XQimlJ4F3NZl+TcP9BHyq3W1IktQu85QkaTHdHj5bkiRJkgaOhZAkSZKk7FgISZIkScqOhZAkSZKk7FgISZIkScqOhZAkSZKk7FgISZIkScqOhZAkSZKk7FgISZIkScqOhZAkSZKk7FgISZIkScqOhZAkSZKk7FgISZIkScqOhZAkSZKk7FgISZIkScqOhZAkSZKk7FgISZIkScrOqn43QBpG23ZOcuWOCfZNzXDc2jVcevY4mzaMdWx5SZLKME9JS7MQklZo285JLrtlNzMHZwGYnJrhslt2AzRNGitdXpKkMsxT0vJ4aZy0QlfumHg1WcybOTjLlTsmOrK8JEllmKek5bEQklZo39RMV6dLklSGeUpaHgshaYWOW7umq9MlSSrDPCUtj4WQtEKXnj3OmtUjr5m2ZvUIl5493pHlJUkqwzwlLU/bhVBEnBAR34+IxyLi0Yj4oybL1CLiuYjYVdw+V665Uv9t2jDG5eefwuEjc91nbO0aLj//lJZfKF3p8pI6wzylXJmnpOUpM2rcIeDTKaUHI+Io4IGIuDOl9DcLlvtBSunDJbYjDZxNG8a44b6fAnDTJ87o+PKSOsI8pWyZp6Sltf2JUEppf0rpweL+88BjgP86kCQNBPOUJGkxHfmOUEScCGwAftRk9hkR8VBEfCci3tGJ7UmStBLmKUnSQqV/UDUijgS+AVySUvrlgtkPAm9NKU1HxEZgG7C+xfNsBjYDjI6OUq/XyzattKli2MhutWVqaobZ2dmBiLUXpqenKxXrYsdHs1i7fTz1S9Ve11a63V+renwMgirnqW4zTw0389Scqr2ui8mtz5ZVqhCKiNXMJZevp5RuWTi/MeGklLZHxFci4s0ppWebLLsV2AowPj6earVamaZ1xNUT9wBQq3XnWtmrJ+5hamqKQYi1F+r1eqViXez4aBZrt4+nfqna69pKt/trVY+Pfqt6nuo289RwM0/Nqdrrupjc+mxZZUaNC+CrwGMppT9rscxbiuWIiNOL7f2s3W1KkrRc5ilJ0mLKfCJ0JvB7wO6I2FVM+4/ArwGklK4Bfhf4w4g4BMwAF6SUUoltSpK0XOYpSVJLbRdCKaUfArHEMlcBV7W7DUmS2mWekiQtpiOjxkmSJEnSMLEQkiRJkpQdCyFJkiRJ2bEQkiRJkpSd0j+oKvXbtp2TXLljgn1TMxy3dg2Xnj3Opg1j/W7WQHDfSFL/eS5uzX2jfrIQ0lDbtnOSy27ZzczBWQAmp2a47JbdANmfSN03ktR/notbc9+o37w0TkPtyh0Tr55A580cnOXKHRN9atHgcN9IUv95Lm7NfaN+sxDSUNs3NbOi6Tlx30hS/3kubs19o36zENJQO27tmhVNz4n7RpL6z3Nxa+4b9ZuFkIbapWePs2b1yGumrVk9wqVnj/epRYPDfSNJ/ee5uDX3jfrNwRI01Oa/TPkf/tfDvDz7CmOOOPMq940k9Z/n4tbcN+o3CyENvU0bxrjhvp8CcNMnzuhzawaL+0aS+s9zcWvuG/WTl8ZJkiRJyo6fCDUx/+Nek1MzHD5yGNt2TvoxraSuevb5l3j6FzOctOV2f1RQkqQe8BOhBeZ/3GuyGLrx5dlXuOSmXWz4kzvYtnOyz62TVDXbdk5y6hfu4CfPvsDLs6+Q+P8/Kug5R5Kk7rEQWqDZj3sB/OLFg74xkdRR8/94mZo5+Lp5/qigJEndZSG0wGI/4uUbE0md1OofL/P8UUFJkrrHQmiBpX7EyzcmkjplqfOJPyooSVL3WAgtcNY/GiUWmf8P16zuWVskVdvaNy5+Pnnx5UNejitJUpc4alyDbTsn+cYDk6RFlnmheGNSZjSnxlHpAE7ccrs/IiYNqMb+GsCGP7mDqRcPlh7ZbdvOSab/76FFl5n/biLguUGSpA7zE6EGX/jWo4terw9wcDbx6Zsfavu/tAtHpZvnKFHS4FnYXxNzxUknRnb7wrce5eAri/3bZc7MwVk+f9ujbW1DkiS1ln0hND907YlbbucXL75+5KZmZlPikpt28Y7PfXdFb4K27Zzk0zc/1LLY8g2PNFiW+ufISvtsO+cbgKmZgys+30iSpMVlc2nctp2TfP62R5sOU9uuF16e5ZKbdnHJTbs69pzzb3j+00dO8VIYqU/mPgl6mJmDryy57NTMQU7ccnvX29TqfHP0G1fzx7/9Ds8XkiStUKS09KUZLVeOOAf4b8AI8FcppSsWzI9i/kbgReBjKaUHl3reN6xbn9Zd9F/bbpfyMz/ARftHc/nt92vbS+n3vtHw2n/9Jby0/4nFxo8ZeOap9h0WsIyrN7UMg7AvV5qnepXXBmHfVMkgvx/phrJ5qu1L4yJiBPgycC5wMnBhRJy8YLFzgfXFbTNwdbvbkxaT6G/HH+STTr/3jVZm1WFDXXcMFPNUOb457ZxB2JcrbUKvmjwI+6ZK3J0rU+bSuNOBPSmlJwEi4kbgPOBvGpY5D/hamvvY6d6IWBsR61JK+xd74uOfP8B//sFXSjStcw6LuTclr5T45EySluOwCN6w6rAlB23plY/2uwHlZZGnJClXZfNUmcESxoCnGx7vLaatdBkAImJzRNwfEfeXaFPHjBwG6448jF8/eu629h/4X1pJ3dF4vnnrUXi+6ZxK5ylJUjllPhFqlqkXfmyynGXmJqa0FdgKc9def+Z9/65E09p3xOEjTQcqeBcr+wL1cvzb9/4ak5OTfP/pwfjvr6TFnXXCCGNjY/z1vT/tyPP18nzTlicv6d+2O6OSeUqSVCiZp8oUQnuBExoeHw/sa2OZgbCckZc2bRhj04ax0iPQNW6rXv8Z5/2T9R0f0U5S58wXLGufe4Ja7RROe+ubOnYOaKVT55vMVSpPSZI6q0wh9GNgfUScBEwCFwD/esEytwEXF9dlvwd4bqnrrrutE0PNzr9B6ZTG5/NNT3veuHruKs8Xe/zf8/nRbgZ5lJZ+7ZthNv+6jq1dw6Vnj7/aP+v1J4DOnwMW02pbniuWZSjz1KBwNK/O6WeuWOm2e93WYcijw+SNqw9jdvYVXjLlL0vbhVBK6VBEXAzsYG5Y0mtTSo9GxCeL+dcA25kbknQPc8OS/v5ynvvEXzmMiSv+WbtNG2q9fIPVa/V6nVqt1u9m9ISxqtt6ca6IP/3wA13dQJeZpzojpz5urNWUU6yQV7xl81SpH1RNKW1nLok0Trum4X4CPlVmG5Iktcs8JUlqpcyocZIkSZI0lCIN4O/jRMTzwES/29Ejbwae7XcjesRYq8lYq2s8pXRUvxsxiMxTlWWs1ZRTrJBXvKXyVKlL47poIqV0Wr8b0QsRcb+xVo+xVlNOscJcvP1uwwAzT1WQsVZTTrFCXvGWzVNeGidJkiQpOxZCkiRJkrIzqIXQ1n43oIeMtZqMtZpyihXyi3clcto3xlpNxlpdOcVbKtaBHCxBkiRJkrppUD8RkiRJkqSu6WshFBH/IiIejYhXIuK0huknRsRMROwqbtc0zHt3ROyOiD0R8RcREf1p/cq0irWYd1kRz0REnN0wfShjbRQRn4+IyYbXcmPDvKZxD7uIOKeIaU9EbOl3ezotIp4qjstd86O1RMSbIuLOiHii+Ht0v9vZjoi4NiIORMQjDdNaxjbMx3CLWLPrr0sxT706r7J5CvI79s1T5qlh0JM8lVLq2w34DWAcqAOnNUw/EXikxTr3AWcAAXwHOLefMXQg1pOBh4A3ACcBPwFGhjnWBXF/Hvj3Taa3jHuYb8BIEcvbgMOLGE/ud7s6HONTwJsXTPsvwJbi/hbgT/vdzjZj+6fAbzaef1rFNuzHcItYs+qvy9xP5qmK56kijmyOffOUearfMZSMtaN9ta+fCKWUHkspLfsH6SJiHfArKaV70lzUXwM2dat9nbRIrOcBN6aUXkop/R2wBzh9mGNdpqZx97lNnXA6sCel9GRK6WXgRuZirbrzgOuL+9czpMdqSulu4OcLJreKbaiP4RaxtjLUsZZhngLyzVNQzWPfPGWeGgq9yFOD/B2hkyJiZ0T8n4h4XzFtDNjbsMzeYtowGwOebng8H1OVYr04Ih4uPuKc/7i2VdzDrqpxNUrAHRHxQERsLqYdm1LaD1D8PaZvreu8VrFV9bXOqb+WZZ6qVqy5HPtVjGkh89Scqr7WHeurq7rRukYR8b+BtzSZ9dmU0jdbrLYf+LWU0s8i4t3Atoh4B3MfvS80MMPetRlrq5gGOtZGi8UNXA18kbm2fxH4EvAHDFF8K1TVuBqdmVLaFxHHAHdGxOP9blCfVPG1zq2/AuapQqXzFJirGlQxpoXMU3Oq+Fp3tK92vRBKKf1WG+u8BLxU3H8gIn4C/Dpz1d3xDYseD+zrRDs7oZ1YmYvphIbH8zENdKyNlht3RPwl8O3iYau4h11V43pVSmlf8fdARNzK3EfPz0TEupTS/uJymQN9bWRntYqtcq91SumZ+fuZ9FfAPLUMQ5+nwFzVoIoxvYZ5yjzFMmMdyEvjImI0IkaK+28D1gNPFh/3PR8R7y1Gpvko0Oo/WMPiNuCCiHhDRJzEXKz3VSXWokPO+wgwP/JH07h73b4u+DGwPiJOiojDgQuYi7USIuKIiDhq/j7wIeZe09uAi4rFLmIIj9VFtIqtcsdwhv21beap6uQpyO7YN0+Zp4ZWx/tqv0aCKEZ4+AhzFdxLwDPAjmL6PwceZW70hweB325Y57Qi6J8AV1H8KOyg31rFWsz7bBHPBA0j7gxrrAvi/h/AbuDh4iBdt1Tcw34DNgJ/W8T22X63p8Oxva3olw8VffSzxfRfBb4HPFH8fVO/29pmfDcwd8nTwaK/fnyx2Ib5GG4Ra3b9dRn7yTy1yOs/rLE2iT2rY988ZZ4ahlsv8lQUK0qSJElSNgby0jhJkiRJ6iYLIUmSJEnZsRCSJEmSlB0LIUmSJEnZsRCSJEmSlB0LIUmSJEnZsRCSJEmSlB0LIUmSJEnZ+X9DAdT5a23JRAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1008x360 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#导入库\n",
    "import numpy as np\n",
    "#用于进行数学计算\n",
    "import matplotlib.pyplot as plt\n",
    "#用于绘制图形\n",
    "import math\n",
    "#用于数学计算\n",
    "#定义一个函数，用于生成正弦波\n",
    "def getSin(amp, freq, phase, sampleList):\n",
    "    return amp * np.sin(-2 * math.pi * freq * sampleList + phase)\n",
    "# 定义一个函数，用于生成余弦波\n",
    "def getCos(amp, freq, phase, sampleList):\n",
    "    return amp * np.cos(-2 * math.pi * freq * sampleList + phase)\n",
    "# 定义去噪函数，将振幅小于阈值的噪声去除\n",
    "def denoise(arr, thresh):\n",
    "    mask = arr > thresh\n",
    "    return mask * arr\n",
    "# 设置采样率为 3000\n",
    "srate = 3000\n",
    "# 在0到 1的范围内生成等间距的采样点\n",
    "t = np.linspace(0, 1, srate)\n",
    "#生成四个正弦和余弦波\n",
    "s1 = getSin(amp=1.5, freq=30, phase=0, sampleList=t)\n",
    "s2 = getCos(amp=3, freq=5, phase=0, sampleList=t)\n",
    "s3 = getSin(amp=10, freq=100, phase=0, sampleList=t)\n",
    "s4 = getCos(amp=20, freq=120, phase=0, sampleList=t)\n",
    "# 将四个正弦和余弦波相加，得到混合信号\n",
    "mixed_signal = s1 + s2 + s3 + s4\n",
    "# 使用 numpy 的 fft 函数计算傅里叶系数\n",
    "fCoefs = np.fft.fft(mixed_signal, srate)\n",
    "# 计算振幅列表\n",
    "amp_list = 2 * np.abs(fCoefs / srate)\n",
    "# 将振幅数组和频率数组进行平移，使 0 Hz 位于中心\n",
    "freqs = np.fft.fftfreq(len(amp_list), 1/srate)\n",
    "amp_shifted = np.fft.fftshift(amp_list)\n",
    "freq_shift = np.fft.fftshift(freqs)\n",
    "# 创建两个子图，分别用于显示原始振幅频谱和处理后的振幅频谱\n",
    "fig, ax = plt.subplots(1, 2, figsize=(14, 5))\n",
    "# 绘制原始振幅频谱\n",
    "ax[0].stem(freq_shift, amp_shifted)\n",
    "ax[0].set_xlim([-150, 150])\n",
    "ax[0].set_ylim([-0.5, 21])\n",
    "ax[0].grid()\n",
    "# 对振幅频谱进行去噪处理\n",
    "amp_shifted[(freq_shift > 110) | (freq_shift < -110)] = 0\n",
    "amp_shifted = denoise(amp_shifted, 1)\n",
    "# 绘制处理后的振幅频谱\n",
    "ax[1].stem(freq_shift, amp_shifted)\n",
    "# 设置第二个子图的坐标轴范围\n",
    "ax[1].set_xlim([-150, 150])\n",
    "ax[1].set_ylim([-0.5, 21])\n",
    "# 给第二个子图添加网格线\n",
    "ax[1].grid()\n",
    "# 显示图形\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ef77bb12-d71b-4fe6-b361-5cf388d3d006",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.9.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}

一键复制编辑原始数据按行查看历史

提交于 2023-04-16 05:15 . ifft去噪还原

马建仓 AI 助手

尝试更多

代码解读

代码找茬

代码优化

https://gitee.com/yang-haoy/ifft_denoise.git

[email protected]:yang-haoy/ifft_denoise.git

yang-haoy

ifft_denoise

逆快速傅里叶变换ifft去噪还原

master

杨皓宇/逆快速傅里叶变换ifft去噪还原

简介

发行版

贡献者

近期动态

杨皓宇/逆快速傅里叶变换ifft去噪还原 .gitee-modal { width: 500px !important; }

简介

发行版

贡献者

近期动态

搜索帮助

杨皓宇/逆快速傅里叶变换ifft去噪还原