{ "cells": [ { "cell_type": "markdown", "id": "8136f6c4", "metadata": {}, "source": [ "# Beam Sweeping and Probing Along a UE Trajectory\n", "\n", "This notebook demonstrates how **NeoRadium** can be used to model beam sweeping and beam probing in a DeepMIMO scenario with a UE moving along a predefined trajectory. It creates a trajectory-based channel, configures CSI-RS resources and CSI reports for sweeping and probing, and visualizes how the selected beam direction evolves as the UE moves through the environment.\n", "\n", "The example is intended as an illustrative API demonstration rather than a fully standards-faithful implementation of practical 5G NR beam management. In a real system, beam management procedures, CSI-RS resource design, CSI reporting, triggering, timing relationships, and scheduling are subject to additional 3GPP constraints and implementation-specific details. Here, the workflow is intentionally simplified so the beam sweeping and probing APIs in **NeoRadium** can be demonstrated clearly and compactly." ] }, { "cell_type": "code", "execution_count": 1, "id": "0d0aaa42", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import time\n", "import matplotlib\n", "from IPython.display import HTML, Markdown, display\n", "\n", "from neoradium import DeepMimoData, TrjChannel, BandwidthPart, AntennaPanel, PDSCH, random\n", "from neoradium import CsiRs, CsiRsSet, CsiRsConfig, CsiReport, CsiReportMan\n", "from neoradium.utils import toDb, toLinear" ] }, { "cell_type": "code", "execution_count": 2, "id": "05aeedf9", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "DeepMimoData Properties:\n", " Scenario: asu_campus_3p5\n", " Version: 4.0.0a3\n", " UE Grid: rx_grid\n", " Grid Size: 411 x 321\n", " Base Station: BS (at [166. 104. 22.])\n", " Total Grid Points: 131,931\n", " UE Spacing: [1. 1.]\n", " UE bounds (xyMin, xyMax) [-225.55 -160.17], [184.45 159.83]\n", " UE Height: 1.50\n", " Carrier Frequency: 3.5 GHz\n", " Num. paths (Min, Avg, Max): 0, 6.21, 10\n", " Num. total blockage: 46,774\n", " LOS percentage: 19.71%\n", "\n" ] } ], "source": [ "# Replace this with the folder on your computer where you store DeepMIMO scenarios\n", "dataFolder = \"/data/RayTracing/DeepMIMO/Scenarios/V4/\"\n", "DeepMimoData.setScenariosPath(dataFolder)\n", "\n", "# Create a DeepMimoData object\n", "dmData = DeepMimoData(\"asu_campus_3p5\")\n", "dmData.print()" ] }, { "cell_type": "code", "execution_count": 3, "id": "b1fd2c2a-9c42-46c2-ac0a-55af9d0c053e", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Trajectory Properties:\n", " start (x,y,z): (35.45, 89.83, 1.50)\n", " No. of points: 12,782\n", " curIdx: 0 (0.00%)\n", " curSpeed: [14.93 0. 0. ]\n", " Total distance: 191.37 meters\n", " Total time: 12.781 seconds\n", " Average Speed: 14.973 m/s\n", " Carrier Frequency: 3.5 GHz\n", " Paths (Min, Avg, Max): 5, 9.09, 10\n", " Totally blocked: 0\n", " LOS percentage: 100.00%\n", "\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAqQAAAG/CAYAAACOv8VsAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjExLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvctoD+AAAAAlwSFlzAAAPYQAAD2EBqD+naQAAXmVJREFUeJzt3Qd4U+Xbx/G77L1k7w2yRAQBeVkCLgQUkI2A4N7gwgmKgAsHijhQhiCCIktBmSoi/AEHoGyQJXuvMkre637qiUlI27RNcnKS7+e6QtvkNDknJzS/3M+Kc7lcLgEAAABsksGuBwYAAAAUgRQAAAC2IpACAADAVgRSAAAA2IpACgAAAFsRSAEAAGArAikAAABsRSAFAACArTLZ+/AAItW5c+fk008/ldWrV8upU6fkueeekwoVKti9WwiD1157TU6fPi0vvPACzzeAsIhjpSYgOBYsWCATJkww37/44otSunTpS7Z5++235bfffpO8efOa7yPZjTfeKJs3b5ZHHnlEcubMKa1bt5ZChQoledxPPfWUVK1aNdn73Lp1q8ycOdPcry4SV7FiRWnXrp2UL1/e7/YahufMmSM7duwwz1n16tXN9rly5fK7/TfffCOzZ8+W999/Xy5evCh33HGHub53797SrFkzr23/+usvefXVV+Xuu++Whg0bmuus36lfv77ce++9yR7L2bNnzbGsWLFCjh49ap6bJk2ayHXXXSdxcXGXbH/y5EmZPHmybNy4Uc6cOWPC/Q033JDicxYqb731lhw8eFCGDBlyyW0tW7Y0x7Ry5Upb9g1A7KHJHgiSdevWybhx48xFK4v+Askzzzxjbv/iiy8i+nnfuXOnzJ0714TM+++/3wQ6f2HU87j37t2b5P1p0NP7qlKliixatEgqV64sl19+ufz444/mOr1Nt7FoWL3rrrvkqquukt9//91so4F02LBhUqJECRNS/XnnnXdk9+7d7se0zoceQ0JCgte2//zzj7lty5YtXvup1/3www/JPj9LliyRSpUqycCBAyVHjhxSt25dOX/+vPTo0cN8v23bNq/t582bJ2XKlJE333xT8uTJY479f//7n9SoUcP8jh3mz58v06dPt+WxAcAXTfZAkGm1TUPN888/71UpmzJligktGlg08EUyrUiqIkWKBOX+nnzySXnjjTdk4sSJ0rVrV/f1DzzwgHz++efSrVs3Exi1qVhNmzZNPvroI1PBfPzxx93bDxgwQJ5++mn3/nk6fvy4LF68WEaNGuV1vVaqtRo6ZswYE3LTa82aNXL99ddL48aNTaDLli2b1/41bdpUmjdvLr/++qsUKFDAHNftt98uRYsWlVWrVrm3v++++0xQHjx4cLr3CQCcjkAKBJlWE7UZWKt/Gk4sWjVt06aNqZT6BlLd9pNPPjHfa4jVJulatWpJ586dTUXNopXV77//3oQrrTTOmDHDXK/N6a1atQq4b+jUqVNl+fLlcuHCBalWrZqp0uXLl8/dlKvVUev7L7/80oQ67YaQFtpEPWLECBM6PcOoRa/TZnbd5s477zTVU20GVzfddJPXtpkyZTIh9fDhw5fcj1ZN9Xj0Ofakz8umTZtk0KBB0r17d9P9ID0effRRU8EdP368Vxi1ArwGaW26Hzp0qLz++uvmXGv1uGPHjpds36hRI9OMH2ifTv2Qo6+Bn376SXLnzi1dunSR2rVre2372WefmeqnypAhgzmvDRo0kFtvvVUyZ85srn/22Wfljz/+MK9Ffb2qjBkzmteVr6+//tp0y9BKsL4etWrtSZ9zDeZ6zrSvsZ4/PadJVdQBwB+a7IEgq1OnjgmTY8eOdV+nfSa1mbdPnz5+f6dYsWKmj6NetPKmwUbDzBVXXCFHjhxxb6dv+lp9fe+990zztPa91BCgFbvHHnssxX3TPoP16tUzoSp//vxStmxZGT16tGkS16ZxVbNmTXfI0e91n7Sqm1YauLQp3Ao+/uhtuo3VlUGbt5U2a/ujlUdf2p9Tg1fhwoX9Bro9e/aYKm167N+/34QzDcr+Hkfp+dPzopVf69xqENTqqGe3hOSOxdd3330ns2bNMv1aFy5caM7Xn3/+acKhBmNP2pXAei1dc801kjVrVnnwwQfNa8R6fH0NXHbZZSacW9t6fniy6GtKg771+tDn17M7g1Z/9X719aT3pV0Q/v77bxO0165dG8AzCgD/0kFNANJv5MiRLv0vtWLFCteIESNcuXLlcp08edLc9vTTT7uKFSvmunDhguv66693FSlSJMX7O3bsmKtQoUKuJ5980n3dgAEDXHFxca7777/fa9tXX33VPPbixYuTvc8ePXq4smTJ4tq0aZP7uhMnTrgqVqzouvzyy10JCQnmukWLFpn7mzVrVsDHrb/jT8eOHc3t+/btS/I+9uzZY7a57bbb3Puk+5MxY0ZX165dXePGjXP9/vvvrosXL/r9/fPnz7vy5ctnngfP6/Q++/bta37u1KmTOSd79+41P8+bN8/cPmHChEt+p3Pnzn4fZ/78+eb2wYMHJ/uctGvXzmx36NAh8/OQIUPMz7Vq1XINHz7c3I+e30C1aNHClTNnTtc777zjdX23bt1c2bNnN89fcv766y/z+J9//rn7utatW7uqV6+e7OO9//777uv0tVG5cmVX8+bN3dfpa13v97vvvvP6fT1/Bw8eDPj4AIAKKRAC2gSuo7C1uVurUlrF0uu0WTQpWkXVSp5WwbSS+tBDD5nme6ty6fEhUu655x6v67QvplbCtI9mck31uj/t27c3o9st2j1A+zLq4CSt4gWbNgur5JrKrVHz1rb6s1aDtcuANs9rP1Kt2mpVUvuQavO1J63a6ahwHYGfFB0Qpc+BNt2H8lj8HY8OZtMuFlpp1sq5jq7XCqVWF3UmgUBYA708aZ9VHbGvfW496XOnfVN1e60+v/LKK6a7g+9rKTn62uvXr5/7Z23+1xkEtG+sxXo9a99dbbr3PH49PgAIFIEUCAHtP6fTJmn40P58u3btSrbJWqdN0pHX+sau0wFps682o2r/Pw1avjwDpcqePbuULFnSTKuUFB19Hh8fb5pffVlTD3mOOA8WK5gcOnQoyW2s2zxDjIY+Ddran3Xfvn1m37QfqgZL364P2pdWj0H7LyZFm9E1yH/88ceyYcOGkB2L5+2ezfF6PrWPpgZ/Paf6IUXDqJ7r7du3p/jY2o9XP3R4so7X87zrjAXaVK/X6TRZ2hSvj62B1N9rKSnanUN/x/d1fezYMROC1ZVXXmk+QOk50W4m2k9Wu5Po+QKA1CCQAiGiAVQrdzoYSOe11MFDSVW+tNKlAVbn0dR+ezoXptWv0h+t9PnSiqxvgPBkDWjx97saVD23CSY9drVs2bIkt7Fus7ZNKlDq3K0asHRQ1okTJ7z6jyZXHbXooCANujrqPy00gOlzlNyxaKVQ5+/UMJjUfKm6DzrwR6eB0tkBtHKdEj2/SV1nnXcdQKWDvnQ6Ku1r/PDDD5vwftttt7nPcaB8w69VJVWeU2jprAb6YeHll182g7Z0AQX9wKQfrgAgUARSIERuvvlmU1H7+eefk62OavPzgQMHLhk4pOFCB4j4o5Pre9JR3Dqvpg6CSorO36n742+yc+s6HYwVbFrV1GCmYVLDty+9Tgdo6TbWKHzdH88mYE/abK+/o5U6pc3QWmEMJJDq8WsFUSuqOsgstTRI6kh9DVtJNX9rtwkdPGY1r2tVUquiSR2L8hy4lhR9Pej9erK6WFjnXafD0ufG97Wkszj40hCb1Aee1LKqzzrCX2c00Kb8kSNHBuW+AcQGAikQIlpJ07lHdbonDWXJhZxy5cqZydOtypNWvrRq6m8eUK1SffDBB2Z0vdJQoRU/fTzPPn/++gRqX1F9HB2xbdGwpFUurdDqCO1g0xCoKyf98ssvZu5Nz0qfVmt1n/Q2He1vNYnrNEItWrS4pGldg7g24WsA0y4KSsOlPk/JVVc96cpTpUqVMrMYpIX+njZnd+jQQdavX+91m47A176/1157rTlWqx+pNqFrxdIzAGrFUiukel502q6U6HOj/V+t+9AKsS7tqR80brnlFnOddsfQ14E1bZfSDzvaF1e7dXjS39MPMkkF/0Do+fAN9lYV1XO6MgBICfOQAiGkE6QHQsNYp06dTDOvTp2jc0Q+8cQTZpoh36ZWDTA6z+nVV19t+k3qpO9aPdM+iSkFSp1/Uqc/0oFNOjWPhmHtVqDT+fhOH5Raw4cP95rqyqLBWgd06UpLGgZ1Sqf/+7//M7dp9ViDkoZKz/lDdfCMBh19PjRoav9J7QOr2+vPWonzbK7X37WCUEq0Wfmll15Ktmqtc7T6u12n9NLAaYVrrSjr86hTO2mztZ43rYxqs7nVjK7Ta2mTuR67Nmfr/us51WZ/DZc6b6m1dGlyNIBrENfBXXqedR+1GqrPnRU2CxYsaOZz7d+/v6me6n7p9FD6GLq9p759+5q5b/Xca3cSDbL+5iFN6bnU49LqvJ4rrYzqHKnaH9rfkqQAkBTWsgeCRKtlGjI0HCU3wlgrlDpyXCcZ9x0Io7+vFUQNKBomtNKlVVOrgqb9S7XapVUtvQ9rTkidiD01o5q1aVdHYlsT4+t8o560cqaPrWuaW5XIlI47KToxvVbjlAYorappeLPWstd+mf7Wfrf2QwO3hmht0tewroO+LBrENazqxPq+VUa9f61K6sAfrVB60iCoTev63Opzp03Onr+TFH0srX5adPCOPo/aD1TDoIa7pCqDel51lSfthqGPq1VabVrPkiWLpMRzbXk9dxou9fnQ/rQ6Yb0v3UYDqYZM3UYn0Z80aZKpxHuGX91/Ddfa/UHPga4opXTxBa1ea7cTTzoIS0fZ6wcMz/7KOoBKj02fPw3LGk4BIDUIpICDeAZSiBnRrd0VtG+l7ypI0cQzkAJANKIPKQDH0iZs7acbzWEUAGIBfUgBOJbVFxUA4GwEUsBBunTpEpKpmRDZdICbv/ljASBa0IcUAAAAtqIPKQAAAGxFIAUAAICt6EOaBjqHoU4ErXP7JTV/IgAASKRz1OrqYsWLF79kEQu9Taeys1aqQ3TQhTJ0vuJAcxKBNA00jOqk1gAAIHC6mIXnYhs6WE8Xvjh9+jRPYxTShTt0kZdAFgAhkKaBVkat/1is1wwAQPJ0NTMt5Fjvn1Zr47Zt20wlTSunGlpodYwOWvXWDxsHDhww51hXcEtpeWcCaRpY/2E0jBJIAQBI3fun0sCioVSDqr8lcOFs2bNnN8sXb9++3ZzrlBYwYVATAACwTUqVM8TGueVVAAAAAFvRZB9sETLqPhx74QrDY0SKyDirkcMVpufKFYXnLpT/byLtWGPWoPA+nGtQ5L1WtA8hkBoEUiAFvMnb91zp/aX3bU1/n3MIOMuOHTvk4MGDYXmsggULSunSpQPe/siRI3Lbbbe5f9Z+kjpzQOfOnaVly5aXbDd8+HCpW7duksF9+vTpMmvWLNm7d6/ZlxtvvFE6depkBnv5+v777832+vxcdtllcsUVV0jfvn0lb9685va3335btmzZIu+88477d/bv3y/33XefVKhQQYYNGyYjR440j+erRIkSMm7cOPf9WNvo8eXLl0+qV68u7du3l2rVqkkoEEiBFBBoAkPoAxAMGraqVKki8fHxYXlCdbDNhg0bAg6lZ8+elQULFphwp0FTB+wsW7ZMrr/+epk4caJ06dLFa7ukgrX+nga8lStXyuOPPy4dO3Y0I9Kffvpp+eCDD+Tbb791D/bSeVq7du0qixYtkgcffFBuuukmOXr0qKxevVquvPJK+euvv8xxrFu3Tn7//Xf3Y/z9999y3XXXmcD86aefmj6duo1OX+kZWpXnwDLPbXTgmY6W18fWx3r44Yfl1VdflWAjkALpCKWuGAxk4W6IC0aVFIBzaIALVxhV+lj6mKmpkqratWu7K6IaEJcuXSrTpk1zB9KUvPTSSyawrl271lQvLVpVvfzyy+WJJ56Qd999172tViw1gFauXNnrfp599llTxfT1559/mjBav359+fzzzyVr1qzu23SGIM9qrj++23Tv3t1Ubm+44QapVauW9OjRQ4KJQU0AAADpcOrUKdm6dauULVs2oO216vj+++9Lnz59vMKoKly4sKlCfvzxx3LmzBmz7XvvvWea5n3DqBUcfZv3tWLbpEkT0/w/depUrzCaHhpwW7RoIR999JEEGxVSII1V0lip2rlSqFTGSmUYADwNHDhQXn/9dTl//rypRrZq1Uqef/75gJ4knZvz0KFDUq9ePb+3161b1zT56/1qv1Ld9qqrrgrovjUYa2VT+40m1bSuTfK+FVL9+amnnkrx/nWftUtBTAdSXedWS9ZffPGF6XyrLwTfsru/UvlDDz0k1157rftn/bQxduxYUyrXPhdagta+H9EkHP0eoz2I+AtfsRJC/YnlYwcAXzqISYOj9u/U4Dh06FDTNH7nnXem+GRp5VNZg5F86SAiaztrWVXrupRotVQvu3btMrnJ3+AozVC+4VOX+Ax0wvtQdKlwTCDVTyC69JT2W9DS+KZNmy7ZRl8UM2bMkFdeecWrrF2xYkWv7fr16yffffed+XSjo+Datm1rRp3dddddYTkWAADgbJ59SLVfpWaQRx99VHr16pXi2u0a/nTVKq2U+vP333+7g2P+/PnNttZ1KSlTpowZJa+Ftm7dupmBVpkyece9QPqQJmX37t1StGhRidlAqgl/yZIlZqTYI488Yr5PivabaNCggd/b/vjjDzPSbPHixdK0aVNznY46008K+iIKVj8LAAAQOzQ8asFMC11FihRJdlsNmQ0bNpRJkyaZEOtLQ2TVqlWlfPny5udrrrnGVF91W8/lV5OiA5nmz59v+nzq6Hz9Xd9QmhbajeCbb74xfVNjNpBqaNQwGgjtM6HBUjsK9+7d26tCOmfOHNMfQ0Or54g2HaWmnYCtkBoNmK4IiC10q0heXJgnrEfs0Orol19+aUbqpxRGPbNK8+bN5YUXXjB9T62mdR1ZP3fuXBP8LK+99po0a9bMtOwOHjzYXTzT5vxRo0bJAw88cMla8dqdQLsmat9WzTlTpkzxOxo/UPv27ZN7773XPKbuR8wG0kBpCL366qtNOVzDZ82aNc0kslYfUZ3jS4Ot5ycMa6oHvc1fINVPBHqxHD9+XJyCUIpowNRPCNaKRoRSBHtQk4bRjRs3Sq5cuUzoS2o7T2+88YY0atTI5BQdfPThhx+arobaLK/55OuvvzbdACxaTdWKp85BqgOKdJL6Y8eOmZCoI/WTCpo6b+jChQtN83yHDh1MaE5qUJMGYu3OaLG20XE3OqhKj1HD7c8//yzlypWTYIuqQKodbXVCWH1RKG2C18v9998vmzdvNtdpsPSc/FXppwo9EUl10tXJb/UTCZAcwj8Q+QilkU9bMfV9OZwT4+tjBqpAgQIyb94898+aH4oXL25aYz0HEPlu58kqhOkUSjop/4YNG0y41NWXdCUkf83yjRs3NhlH+51aKzVpEc6zq6F2aTxx4oTX7+nYG518XwPl4cOHzTY6Cb+/lmjP+7G20aZ+HVClixVozgqVqAqk+kKwwqilXbt2Mn78eNOnQ/ts6Ig2/d6TrnagI9H0dn/0003//v29KqSlSpUK0VEA8EVTNEKx9jvV0sikYU0DWqQuHaoDlgIZEBTodqpKlSrmEuigJb34o/1O/dHjs45RByQltZ3n/aS0TbBFVSD1R8OmftKwPrXouq9a7tZPELlz5zbXWcts6acIf/TTB4OdACC6UC2NXJ4BCrEhqlZq0s67OiGsZf/+/abfhnYa1ikOrIqpfmrRaZ6U9o3Qvhw64awu1RWt1SXPCwDg37+PDHQCIoKjKqS6lJb2nVizZo3pB3HLLbeY6ydPnmz6gGhzffv27U0/CO3vsGrVKjP90yeffOK+D+1zMWHCBOnZs6fpNKwV1HPnzpmOxQAAAAg/RwVSrW76G+FujS7Tebc0hOqKCQcOHDDzd/kbCab3ox2CV6xYYZriNbSmNIktEIuoqCPWqqT0KwXs4ahA6rn8Z1K0r2hSfUE9aQVVpy8A4B9hFLGIfqWAPaKqDymA4CCMIpbRrxQIP0dVSAHEpvROjO/vd1NefA+xjKmhgPCiQgog5hBGESiqpUB4EEhjDG/EcNo0ZcGcrkxf//wfQKpfk0wNBYQcTfYxhDdihILLQa8xJ+0rIgtN+GG2Y4dImFZqEl02NBWT8J8+fVpGjRpl5je/6667vG7T6SR1eU1dh97aTpfgLFu2bJL3pytSLVy4UPbu3WtWjdIB3Lqaki9du17nW9ftdAahpk2bXrIUupMRSGNAXBQNiomL8EE5DAaKXARQnqNopVNVRVUVV8OoLqMZprXsJVs2kQ0bAg6lOv3k448/br7XJTyvv/56921jxowxoVIDqbVdjRo1kgyk77zzjjz99NPSsGFDszjPvHnzpF+/fvLMM8+Yi0Wv79Spk9SsWdNc5s6dKw888ICMGDHCTGUZDQikSBcCWPSKc+gApqTuE4BDaGU0XGFU6WPpY6ZyqVINhk8++aSZQlIX5EmtL774wiz4M3XqVFNFtSxevNjcpwbbu+++21ynX7t16ybvvfeee7tDhw7Jpk2bJFrQhxRAxLH6ehIkAUQqrWBu2bJFJk6cmKbfHzRokLRu3dorjKpmzZrJ7bffbm7X5c0TEhLMYj4agD3pypO6sE+0IJACAACkUpEiRaR///7y3HPPydmzZ1P1u/v27ZP169fLjTfe6Pf2G2+80fQV1W10wR/tV/riiy/K6NGjZePGjeJyRV/7JIEUQESjSgogUmkf0TNnzng1pQdCBzKpEiVK+L29ZMmSXttp8772IR0yZIhUqVJFChcubAZU6TLp0YJAGgNCNYUOAMCZeB8Ijly5cpkK6csvvyxHjx4N+Pd0+XKlVVB/9uzZ47Vd/vz55a233pJdu3bJ9u3bZejQoTJz5kxp06aNRAsCKQAg1Qg0gLgHHGlgHD58eMBPiVZGdeT9okWL/N6+ePFic5/VqlW75LbSpUvLnXfeaULp8uXLo6ZKSiCNIQwUARDMvycARDJnzmwqpG+//bbs3r074KfkqaeeMiPsNXx6Wr16tXz44YdmBH+mTImTIf3yyy/iSwc6ZcuWzV1FdTqmfYoBvHHA6UIxFRQABIv273z99ddl5cqVcsUVV3jd9tVXX8natWu9ruvcubOprGrzuw5guu2228w8pNu2bZNJkyaZuUituU6Vfn/hwgWpX7++FCpUyITWGTNmmMfUQBwNCKQAAAABypkzpwwYMEBKlSrlvi4uLs6syqSDj+rVq+e1nb++oufOnTNftdm9T58+8u2335qR9zpgSZvhq1ev7rX9kiVLZNWqVfLzzz+bJnqdGuqVV14xKzZFizhXNM4dEGK6+kLevHnNMl66dJiXuMioR9q1Fy4H7TMv/ORFxivZma8tp77uIul5ilgRuiJSWlZqCuX59o0W/t434+PjTUVQQ5U2PTtlpSYELslz7AcVUoQVb3iIpgAHIAQ0GGpAjNC17BEaBFIEFYET6REL4ZP+sEAANCASEmMKgRQxKy7GgxEAAJGCaZ8AAABgKwIpAAAAbEUgBQAAgK0IpAAAALAVg5oARMwsCf4enwFmABD9CKRAjLE7dAIA4ItACgAAIsqOYzvk4OnwTIxfMEdBKZ038InxdU359evXS9myZSVXrlzJbnv06FHZv3+/XHbZZeaSlNOnT5vlRUuWLClZsmSRWEQgBQDAwW7OKSI3ici3EjVhtMq7VST+QniWDs2WKZtseGBDwKH04MGDUrNmTZkzZ47ccMMNfrf57bff5JFHHpGVK1dKiRIlZM+ePVKtWjV588035ZprrnFvt3PnTrn77rtl4cKFZjtdp75Vq1by8ssvS9WqVSWWMKgJAAAHh9GviolIdxFpLVFBK6PhCqNKHyuY1djVq1dL48aN5fLLLzfV0Y0bN5oQ27JlS7n22mtl6dKl7m179OghZ86cMYF1y5YtcuTIEbn99tvlf//7n8QaKqSIWQyWAeB01bOIZLE6hnf79+s3Nu4QZMCAAVK6dGkZNWqUZMiQWPfLmjWrqXouWbJEHnzwQVm1apUkJCSYcPr+++9L/vz5zXYZM2aUdu3axeSzSIUUAACHeuWIyDOexb1u0VMpdaITJ06Y5vfevXu7w6invn37yq+//mqa6jV8litXTj777DPZtm2bxDoCKQAADjb0iIhM8biCUGqbXbt2ycWLF6V8+fJ+b7eu3759u/mqYXTfvn3mer10795dpk2bJrGIQAoAgNPNIJRGAmuEvI6a9+fUqVNe21199dWybt06WbNmjTz11FMSHx8vHTt2lPvuu09iDYEUEdWnM9QXAIhahFLblSlTRnLnzi1//PGH39v1+kyZMpkBT55q1Kghd911l3z11VcyaNAgGT16tJw8eVJiCYEUMYVwCiCqEUptpWGzT58+8tFHH5l5RX37l7799tvStWtXE1qVVkR9lShRQuLi4sTliq0yCqPsAUQ0HUAcW3+WgSCEUtXp36+Mvg8J7Qe6du1ar+t0kNLQoUPNtE1NmjQx1U6thuqgpZdeekkKFiwob731lnt7nWu0V69e0qBBAylUqJCZMurZZ5+VTp06uUNrrCCQRil/b+AsGQkAMYJQGjKZM2eW6tWry8iRI83F06effir16tWTxYsXy5gxY2TChAlm0JKu0tSzZ0+55557JEeOHO7tV6xYIR9++KG5H50Uv1ixYibE6ij9WEMgjVKRHD6pdgFAGDg0lOpSnrp6UjhXatLHDJSGS9/KqC+dd1QHJqU0OEmros8880zAjx3NCKQAEGZxg9L3+650/j5iiANDqS7hqUt5Rupa9ggNAiliVnJVZKq4iOZAixiTQijNllkkT3aR42dE4s9LRNCASEiMLQRSBJXLId0GACDWQ+m7hUSK5hW5pa5IxgwiCRdFpq8UGTFHZOlGG/cVMYlpnxAyzAUKAJE7JdT9rf4Lo0q/tr1K5KfnRe5uYdteIkYRSBH14pK4AEBMhtIT//1ohVFL5owiGeJERvURuaZy2PcOMYxACgBArCgiIgFMb6nN94/eGI4dAhzYh1RXLfjuu+/kiy++kOLFi8vLL7/sd5vPP/9cFixYINmyZTOTyzZt2jTV2zhdcoNyqA4CQIy6JrDNtFJ6a93EAU+RMtAJ0c0xFdLz589LpUqVZMSIEWbFAw2m/tx7773Sv39/szJCvnz5pFWrVjJ27NhUbwMgejBrAvCv+oE/E9qcr6PvgXBwTIU0Y8aM8v3330v58uXlkUcekSVLllyyjU5U+8EHH8j8+fOlRYvEHtlZsmSRxx57TLp162a+D2QbpIw3eABwYHN9qcA312Z7nQoKCAfHBNIMGTKYMJqcb7/91qyg0Lx5c/d1nTt3NstwLVu2zKwrG8g20Y4wCQAxqEEqtr0oknGlyBm7mutP7RA5G56J8SVrQZGcTIxvN8cE0kBs2bJFSpYsacKrpWzZsubr1q1bTdgMZBtfZ8+eNRfL8ePHQ3wkCEcApy8tAnmdpAWvLTi9ud506Dsm9oXRWVVELoZn6VDJkE2kzYaAQ6muOX/llVd6rW2vuUKLW/fff7/Exf33F2DWrFny/vvvm3yhxbCWLVvKQw89ZL6HQ/uQBiI+Pl5y5crldZ0OWtLmfr0t0G18DRs2TPLmzeu+lCqVijYPIMrmkbXjAiAIzfVlAtjuosf3rUSktQ3PvFZGwxVGlT5WKqqxCQkJsnv3bhk+fLhpWdUB0jo25fHHH5d33nnHvd3kyZOlQ4cOcsMNN8j06dPlzTffNHmjb9++IToQZ4uqCqmGxcOHD3tdd/ToUfPi0cFLgW7ja+DAgWYQlGeFlFDqfLEadGL1uNOL5ToRVc31+ja46d8poKr+W57SMLpSRE6KyLX+lxnFfwoWLGgqo0q7FI4fP15++uknefjhh81148aNk1tuucVURC1XX321yRuI8kBaq1Yt+eijj+TUqVOSM2dOc92aNWvM15o1awa8ja+sWbOaC4DoDYE0syPqm+s1hK4QkeUissHjtswioqPpdQCT1Wf0oPcyowahNEnaJP/bb7+ZQdcWzQ3aTVBbX7UyatEWWUR5IG3Xrp08+uijpr+GjprX+Ua1RF67dm2pXr16wNsAsYoqIBCFNP/ozIZJrU+vIfS8nxWdFKE0ST179jShU6el1H6lOlPPk08+6dW6etNNN0mZMmXM9JL169c3zfc6hSUcHki1f8bOnTvNpxA9+V26dHGXxfVFUahQIfnkk0+kT58+MmPGDNMUr83zc+bMcd9HINsA4UIABBByCcmE0eQQSpOlxaxmzZrJhQsX5M8//zT9SJ999lkZOnSouV0D6KZNm+Trr7+WH3/80Wyvzfl60e/hLc6lJUKHmDt3rgmQvjp27CiZMv2XrTWsakdjDamNGzeW7Nkvndk3kG2Son1ItS/qsWPHJE+ePN43eoyuA1JCIAUQDK5Qdm9p51EpVZMCaL73iRb+3je1KVsXuilXrpxXk7Yc/lVk7lUSVjesEilQJ6BN9+7dK8WKFTOFLK14Wt577z3TZK/HmCNHDr+/q4v7DBgwQJYuXSoNGzaUaBef1Dl2eoXU88QnR6ugbdq0Sfc2iE6EQESTOImTMlJGckkuOSknZbtsFxdD1xBMVEoDojP4aLX09OnTSQbSrl27mkC6a9cuXqNODqRwLkIgEHyXy+Vyg9wgeSWv+7pjckzmylxZJ+t4yhE8hNJkaavr6NGjzfykOvpePfHEE9KgQQO57rrrTFjVwdQ6jaSGVb0e3gikUYwQCER3GO3k1Y6aKI/kMddPkSmEUgQXodTvoCatih45csQsR/7uu++6b7/tttvktddeM2NWdFly7bagA6hnzpzJ1JF+EEiDjBAIIBzN9FoZtb73vU2b7PX29bKe5ns4L5TqUp66elI4V2rSxwxQ4cKFzQBrz2mc9Drf6Zzq1asnU6ZMMd8fOnTIzHXOlE9JI5ACgMNon1GrmT5btjPSvftn4nLFXXJp78oqpy7Gy0WXuC8J/37dckTkxR9ETp6z+2jgOKEOpbqEpy7lGaFr2evS49aE+IFiqdCUEUgBwGF0AJMlY8YEKVlyt9/tyqVwP20ri6zZ/19Y1YsOjvb8WRfv8b399HmRj38V2XAoyAcG54bSJiLyvZ/5TNNKA2IqQiKcj0AKAA6jo+kt8fFZ5fPPu0pcnMvjctF8/SnuBzkUd1AyxIm5ZMyQ+DVbJpFnG4tUKZh4SYs764iM+0PkXILI+Ysi5//9uvmwyBdrWaI2VrqLPX1QpHtukebnRPY/89/1jplPEhGDQAoADqNTO+loeh3AlJCQWTZurOJ1u/YhPS7HZaQcSjIYfPmXSOtKIlkyJoZUveg0yhmSuWhvVf3aqrxIo9IiD+lylH50qS6yZOd/IVW/nk0Qmb9V5J8TwX8+YJ+hR0TeOipymgSKdCKQAoDDaODUqZ10NL1+7zmwyZqDVG9Pbj7SvSdFxvyWtscftkTkjitFSuYRyZxBJHPGxK85M4t0ryXSrmrixdeRMyLPLBQ5Eu8dVg+cFln5T9r2BfYjjCIYCKQA4EA6z6hO7eQ7D6lWRkM9D6k2049e6f+2D38V6XWFSPZM/wVV/Vr5MpFqhURGtfb/exNXi3z8W+J96+XshcQ+qvo9gOhHIAUAh9LQqVM7RdJKTct2JV58ab/VQc1ErijiXVXVr3qdVlb14mn7UZFe00V2Hf8vqB4/K3LmQtgOB0CYEEgBwME0fP4tf0uki78g8tR8/7ddW05kUFOR/NkT+7RmzZj4fZl8Iot7e2+rldOnF4p8sDLxPnUaKwDORyAFANhq4bbEi6f82UQ+bivSpExiJdUE1UyJlzeuS7worZpO/VPkrtmJ01EBcCYCKQAg4ujApw6Ji9x4ebi+yMvXiuTMkvizBlVt6r/1cpFj8SK/7xXpN4vR/HCe7du3y0svvSQjR46U7NmzS6whkAIAHOPt5SLvLE/sk6qX2kVFJnUQKZpLJEdmkWK5Rf68T+TvoyJ7Tog8u0jk1z127zVSb4eIhGmlJtHJeAObhH/69Okye/bsZLd55ZVX0rQy04EDB2TMmDHy+uuvBy2Qbt26VYYOHSrvvfeeZM2aVSIZgRQA4CjabVQHNull0d8iZd5KnIKqSM7EZn4dza9BVS/a5L98d2J/0zPnRb5cJzJ5rd1HgJTDqM6tG6a17CWbiGwIKJSWKlVKGjRo4P75iSeekP/7v/+Ttm3buq9La/ArW7asfPTRR5IjRw4Jlv3795uQ+9ZbbxFIAQAIJe1HuvVI4qX2aJE6xUTyZBXp31DkhoqJg6YsHaqJdK8psvuEyB97RT5YlbgkKiLJwTCGUfn3sQ4GFEivuuoqc7EMGjRI6tSpI/369TM/b9iwQR555BEZNmyYfPzxx6YZ/rHHHjMV08cff9xskzlzZilfvrx0795dihcv7r6vU6dOybJly8z1FpfLJTNmzJAffvhBsmXLJs2aNZPrr7/ea590m1mzZrm36dy5s9SqVctUXHU/1AMPPCCZMmWShg0bSt++feXs2bPy2WefyW+//SZ58+aVjh07ypVXXum+Tz2O1157zes4+vfvL2+//bb06tVLrr76ave2x44dM8eo4bxSpUppPAciGdL8mwAARBidbF8rovO2itw8SaTleJEuX4r0ni7y5rLEbW6uLHL3VYlzoi7qlbgMqvZFBdJrz549piJ5zTXXyOHDh014zZ07t6lOamVVLzVq1DBBsFq1arJly5ZLmuzPnj3rDpoaFJ977jkpWrSo5MmTR+68804T/CwXL16U9u3bm+u1spozZ04TOH/55RfzmNWrVzfb1atXzzx2hQoVJCEhQa699loTOEuXLm32UwPmtGnTkj0ODa4HDx6UESNGeB3zpEmTZObMmabCmx402QMAopJOCbXAc/T+HyJf/SVSr4RI4ZwijzZIbNLXi648pX1Tv16f2LwPpMfzzz8vPXv29LrOqqKq+++/X7p06SJvvPGGjBo1yu99aNBbuXKl/PXXXyZoKg2fGmTvvfdeKVeunIwfP17mzp1rttGflVZijxw5YgKsdiXQKqdWNXPlymVuHzt2rKxevVq2bdsmBQtq/1kxFdyHH35Y2rRpYyq4SR2HHoNuo/efP39+c90nn3xitvH8vbQgkAIAYsbPOxMvasqfIn1qJ4bRBiUTL0t3itw0UeRYYpEKSJMbb7zxkut27twpkydPlh07dsiZM2dMdXTfvn1J3sfs2bNNyBswYICpllqXjBkzmkCpAfTbb7+Vli1busOo0t8pXLhwkve7ePFi8ztWGFU9evSQl19+2QyCqlKlSpLHob+n3QwmTpxougGsXbvWhOZx48ZJetFkH2yD/r0AACKaThH18FyR+h+LfPJb4jym15QSWdhL5LLYm3UHQZQvXz6vn1etWmWC3u+//27Cozaf6wCpEydOJHkfBw8eNE31devWNU3u2qxev359U1HVZn+llcoiRYqkeqCTZxhV1s++Adn3OOLi4uSOO+4wVVGlzfq6T1q1TS8qpABiQ7g/KPLB1DH+OiDSd6bIW8tE5vVMHBT1Q2+R7tNEVu9LHNUPpMenn35qBiNpZdGilUWtlialRIkSsn79eq+mfl8lS5b06ofqSwOkLw3C2lzvSQctKe1TmpI+ffqYwVwrVqwwxzNkyBAJBiqkAACIyJr9Ik3Hiuw6LlK9sMjv94isuDNxjlMgPbSp3bMa+vfff5vm++T06NFDli9fLlOnTr2kKV+b/FXXrl3N6PoFCxa4b9+1a5cJsqpAgQLmqw5MsuhAqYULF8off/zh3jfty6oDl8qUKZPisWhQ1qb822+/3cwMoH1hg4EKKRCJBkVu5c2VzONf+lk8htF9x5E2HBJp/KnI5x0S+5ReVVzkx94iLcaL7Dxu997Bqe677z4zX2mjRo2kWLFiJkTq1E/JadmypRkJrwOG3n33XfN72ne0YsWKct11iWvn6tfBgwdL69atpUmTJmZC/c2bN7tHzOuo+po1a8ott9xipnXSUfM6Cv+ee+6Rxo0bS6tWrUx1VPu3zpkzx29F1R+t2up9aijVwVPBEOfSaIxUOX78uJn+QOfe8j0R7pNJcx0CCHDu101qnq1BkbHP/hBI/eBvgWOVyyey4HaRcvlFth9NDKVbjti9V87gesGV4vtmfHy8aTrWPpU6f6YTJsb39fnnn5u+oToXqTVd0jfffGNCn2+402mdNIheuHDBBMOjR4+a42/Xrp25XZvAtZ+oPlc6VZRF+3UuWbLETNdUu3ZtqVy58iX7oU3/S5cuNc+tBlNrRL3SaqpWRPV+NATrXKZKByRZ85A2b97c6zGTOw7r8bSaqgOkmjZtmuTzk/Q5vhSBNFSB1B/emKJKWkNbJATSYO67JwJpEvi/71glcieG0ioFRc5eEJmzWWTA94mT8CNUgTRylw4NpZkzZ5omeG0Gj3Q6HdRXX31lQm1yuSc1gZQmeyCMgS0cnLzvUYvme8fSFZ2ajBX5vofIFUVFbqkqcnWJxAn314UrL8Wk0hEREsNl6NChMnr0aDNgKJJphfeDDz4wXQKmTJkScBN/IAikiCnRENai4RhiEqHUsfafEqn3kUiL8iKvt0oc8KSj8FtNEPkj6WkkgYDVrl3bzOVpNadHKp1iSvud6mpRus/BRCBFVCO8IaIQSh29JOnczSL/2y3yXQ+RujrYqY/IyP8lzmFKEz7S46abbnLEE1i1alVzCQWmfQKAcKLC7WiHzyQOblqyQyRPVpFnGidODaUBFUDaEUiBMIvzcwHgHMfPilz/mchj34usPyhSILvI/J6JqzwBSBsCKRAmhE8geugyo2/8kti3dPHfInmzJTblN015XnEAfhBIAQBIo5PnRG6aKPLdZpFcWUTmdBe5rZpIlow8pUBqEEgBAEiHMxdE2k0WmblBJHtmkSm3iezuL9KEaikQMAIpEAFcqbgAiDxnE0Q6ThH56NfEnwvmEPmmG/1KEbgjR47Il19+KefPn4/Jp41pnxBVnDDNE6ESiN6poe6aJfLgtyIzu4pcVyGxCV8n0V/xj9175yy6otPp06fD8lg5cuQwq0gF4s8//5R169alOIWT3mdqbdmyRW677TYTTPPlyxfw7+3evVvWr18vWbJkkSuuuMJrBcnDhw+bZUNvvfVWyZgx/f1Ign1/ngikcBwnhE5/CKJA7FRLb5mcWCFtXi5xsJNOFfXbXrv3zDlh9N133zVrvodDpkyZ5IEHHggolK5Zs8ZUMS3ffvutlC1bVqpVq+a+Ttd2T0sgLVCggHTo0MEEy0C4XC655557ZMKECVK/fn3JkCGDCaadOnWSN99802yzceNGE3JPnDjhtb59WgX7/jwRSBFRnBo2gVj5/xkT05QNCk6/0jafi8ztIfJ/pUXm9RRpPk5kzf5g7GB008pouMKo0sfSxwwkkHbp0sVcLCVLljQBcNCgxBfNgQMHzPKaGix///132b59uzRu3NgE1G+++cZskzlzZilfvrzUqFHDa+nN/Pnzm/vW2wMxfvx4c9HHqVKlirnu7NmzMmbMGPP9qVOnZNGiReb76dOnm7Xky5QpYx43pX3xdxx16tTxe3/16tWTYCCQIqwInICz8H827U6dTxyB/31PkQYlRebfLtL4U5GNh4J4ghBRtElfK4jt2rWTrVu3SqVKlczKRoULF5bJkye7Q+OqVaukYsWKpsJqVRq3pLLJXqu1GgitMKqyZs0q9913n/leq5jz5s0z33/11Vemib1JkyZSoUKFFPfF33Fcdtllfu+PQIqIxJsX4Hz8Pw6eE+dEbvhM5Kc+IjWLiDzXRKTn1xL91eQXJKaVK1fOVBE9eTb1x8fHm6Z9bVp/7rnn0vQYV199tbzxxhvy0ksvSY8ePcxjeipatKgMHTpUGjZsaJr1PZvYA90X3+NI6v6CgQopUoU3KiB68f87NI6dFRn8g8iXnUQ6VhN55WeRtZHUdE9XqaB75JFHkqxq7tixQ86cOWOqm8uWLUvzY3Tq1MkMaBoxYoQ8//zzpgrbokULeeKJJ6R27dop/n4g+5LUcYQCgRReeEMCYgv/58Nj9kaR5btE6pcUWXC7SNOxicuOIjoVL17c6+dDhw7J9ddfL7t27ZJatWqZkfAbNmww/TDT49FHHzUXHWy0dOlSGT16tKmcrlixwoy49yc1++J7HKFEII0xvPkA4O+APSPvb5iYGEbrFPsvlG4+zOsxGnkOEFJaxdTrdu7c6R609Pjjj7sHCaVX5cqVzaVz584mRE6ZMiXJQJqaffE9jlAikEYZ3mgA8LchMh2NF2k1QWRRL5FaRUQW3i7SZKzI30ft3jOE2p49e0xgtAKgTn4/a9asZPthulwuM3joyiuvNAORfGmFs0SJEl6hUX8nISHBPRepdf86eMn6Pi37YvF3f8ESVYFUn6C77777kuv79esn//d//+d13bRp02TBggWmRK0jyRo0aBDGPQUAxKLDZxInyl/cW6RaocRw2uRTkZ3H7d4zhFKbNm3MlE464r5YsWIybtw42bt3rxndnpSEhASTT0aOHGnmSfWlYfXDDz+U9u3bm1Hwx48fl7Fjx5rpq26//XazjU7ppCP2X3jhBTNwSedMTcu+WPzdX7BG2UfV0qGa8vWJLVWqlDRr1sx90ZFmnh566CETXPVEKJ0jbOLEieKUCmhyFwDw9zcDkePA6cSJ8nX6p7L5EqeDyhxV78axo3Xr1l6T4hcqVMjM3amT1HvSlY00QOp8nkuWLJE+ffrIxx9/bAYh+Yr7t+Kp96H3lVRQfPjhh80I+OzZs5sCm07VpPe7du1ad77R+U91qia9z6lTp8ovv/wS0L4kdRz+7i9Y4lxa340SJ0+elNy5c5snKKmK519//WUmgJ0zZ47p1Kt0moMPPvjAjFYLZEJa/RSin0B0NQnPJbpS7G+RwpsCbxpwupiYND0YQhwQQ/m3JCbOcZgCfIncIlsfFsmSMXFqqO+2SNQcr2+08Pe+qdMNbdu2zUwt5DmgJpJXagqlRYsWyXXXXWdae32DoFMldY6jvsneMmrUKFMp1T4X3bp18xolphO/6moIrVq1cl/XtWtXGTJkiCxfvvySpv1gInAC4O8LLLtPiMzdLNK2isj7rRMnzdfrYp0GQw2IkbiWfagsWLDAzCmqLbbREkZTK+oCqZapNYDqV62Cvvjii2aJLD3JavPmzaZJ3/OEW5PJ6ioJ/gKpflrRi+cnPQCIBHzQdbZ7Zif2Ja1YQGRhr8SR93tPhncfXBFY+daAaHdIDKcffvhBatasKf3795dYFVWBVMvB2neiQIEC7v4VOnGs9hfVpnqlE8DmzJnT6/e0/4UugaW3+TNs2DAZPHhwGI4AAFJGCI0ee06KXDtO5Mc+IpUvE5nfU6TZOJGD4SkOIkK8+OKLEuuiqi6s/UCsMGrREWrr1q0zfVKU9l05etR7jg29TUezJfVpbODAgWYb66JzdwFAuDFAKTrpCHsNpbuOi1QvLDKvp0i+9M2XDjhOVFVI/bH6oFgdrLUkPmbMGHO99htRWlVVOtjJn6xZs5oLAERDJTTSmmcjTlqf93Scr21HE0feL+kjUruoSP+GIs8HZ850wBGiqkL6888/m7m0PPt6vv3229KoUSMzb5Zq166dGQn/0Ucfubd75513TBjVsAoAdqMSGpt0GqhXlyZ+3+sKkfwxUiWNosl+kI5zG1UVUm1213m0ihQpYgKoBlQdsDR+/Hj3NnqbTvGk/Upnz54tR44cMU3wOgAKAOxCv1Coz1aL9G8gUjqvyNweiZPonzgXhtefDZVza5pFbbHUsRyI3lbqQKbUjKp5SJWOhtfpmw4cOGBWFKhdu7bfuUH/+ecfE1i1Kb558+Zm/tJApXUe0qh6ogE/aApOfdOuHUGU8xQiQTqXOur+h94iBXOI/LRd5IaJIqfPS8j4ew3GhWEeUmsZSx3XUbhwYdONLpxrpyN09PxrGN2/f78pEFoT9cdUIA0HAingH28lgbH7jy7nKUSC+OFC+5HqsqI6uEmrpj2/lrAE0mC9NgINpLqddrXzHWyM6KBhVFfLDOSDRlQ12QMAEA1+3yvSZ4bI151FOlcXuXt2aKukdtGgotUzrZDq8t+IHtpMr1NqBopACgBABNJVnA6cEimUU2R6Z5E2n4ucTZCopMElNeEF0SeqRtkDABAt4i+ItJ0scvKcSKsKIlNuE8nEuzaiFC9tAAAC5ErmEgrLdoncPEnkzPnENe/H38KpQnQikAJAjEkuVNk94Mqp4kI4W8IP20XaT0n8vmtNkRqFQ/dYgF0IpAAAOKA/6Yrdid8/1cjuvQGCj0AKAECEVkc9PbdIJOGiSPdaInWLh+cxgXAhkAIAEOFhVH23RWTFP4nfF8sVvscFwoFACgCAg0beq2xM2ogoQyAFACDCq6MWHW2vKhYIzf0zsA12IZACAOCAMKrmbU38OriZyOUF7dkHIBQIpAAAOMRbyxLnJs2cUaRBSbv3BggeAikAAA6hTer7TiZ+nzGE7+A03SPc6BYNAECEN9V7Svh39YJwLCOqDxUX+ocBqJACAOAkFy4mfs1IUkQUoUIKAEHE0psINZ0cP9RN9kC48XIGgCCiaIVwNdlTIUU0oULqYLzxhRaVLgCRiAopohEVUgAAHIQKKaIRgRQAAAc5/28fUibGRzQhkAIA4CAz1id+7VVb5KZKdu8NEBwEUgAAHGTOZpHP1yR+X6eY3XsDBAeBFAAAhzkSn/iVkfaIFoyyBxDxsxswo4T95yDUOMdpmxw/HKs1AeFAIAUQNEyVhagOzRGwbOglI+0JpIgSBFIAgO2ojKdxLlJKy4gSfLYCAMBB1VFFkz2iDRVSAEBEogtI0miyR7QhkAKAB0IQIr06qk6fT/zauLRI9kwiZy7YvUdA+tBkD8ARITFcF8AJJq0ROXQ6cR7ShxvYvTdA+hFIAQBwUHVU/X1U5L0Vid+XymP33gDpRyAFAMDBa9pnYKQ9ogCBFAAAB7r4bx8TAimiAYEUAAAHB1Imx0c0IJACAODgyfF5I0c04HUMAIADUSFFNCGQAgDg4Mnx6UOKaMDE+ACAgMTkYO7UTPs0yKYKaUyeGEQbAikAIDbDZrT0IeXkIQrQZA8AgAPRhxTRhAppiIVyKUKWOQSA2HUuIfFrtUIiOTL/t759KN5rKMIi1KiQAgDgQHM2ixw8nRhIn29q994A6UMgBQDAgf45ITL4h8Tvq1xm994A6ROzTfbff/+9LFiwQLJlyyYdOnSQWrVq2b1LABBSNLuG2KDwj8I/eS7xa2bKS3C4gF/C8+fPly1btkg0GDhwoHTu3Fni4uJk7969UrduXZk2bVrQ+tp4XgAACJUL/460z0QgRaxUSLdu3Sq33HKLjBgxQu666y5xqo0bN8qrr74qX3/9tbRt29ZclzdvXnnggQfMz5kyxWzRGADgMARSRIuAP1NpCH377bflsccek9atW5vKohPNnj1b8uTJY47B0rNnT9mzZ4+sWLHC1n0DACA1CKSIFqkq8vft21dWr14tp0+flho1asjnn38umzdv9rocOnRIItmmTZukdOnSkjFjRvd15cuXd9/mz9mzZ+X48eNeFwAA7EYgRbRIda+TsmXLysKFC6VBgwbSrVs3qVSpktdl2LBhEsnOnDkjuXPn9rouZ86cJqBq0PZHj0mb9a1LqVKlwrS3AAAk7dC/b1tXFRepX4JnCs6V6g6Thw8flvvuu8+MUNe+mPXq1fO6PdLDWq5cueTo0aNe12nFMyEhwTTlJzUIqn///l7bR/pxAgCi3887RWZuEGlbRWRwM5EbJtq9R0AYAuncuXNNs33RokVl1apVUq1aNXGa6tWry9ixYyU+Pt5M+aTWrVtnviZ1PFmzZjUXALE35RGzZSDSlw99f2ViIC2Yw+69AcLQZD9p0iRp06aN3HHHHbJs2TJHhlHVrl07uXDhgowbN8593XvvvSdVqlSRK664wtZ9A2IxJIbrAkSr8/8uIZr5v6ERQPRWSLNnzy5LliyR+vXri5MVL15cRo4cKQ8//LCp+GoXBB2o9e2335p5SQEAcJJQD2zinRERFUhvvfVWiRZ33nmnXHvttfLjjz+apvjrr79eLruMddcAAM5z/mLwV2sihCLcYnYW+AoVKpgLAADR0GSfJ2tiKLUCKuAkLDYGAICDbTwkcuKsSJFcIq9fF5z7ZDAfwo1ACgCAgx07K9JvVuL3t1a1e2+AtCGQAgDgcH/uT/yaNWY74sHpeOkCAOBw5/7tR5olRFM/0YSPUCOQAgDgcGdDEEgJoQgnmuwdLJyTisfiBQCcQgc1qRyZRdpUtntvgNQjkAIA4HBH4hOXEFWjb7Z7b4DUI5ACABAFhvyY+LVwTrv3BEg9AikAAFHARadPOBiDmgDAB+/rcCJdoemXnSIJvIDhQARSAPChg9p4T4fTHDwtcs0ndu8FkDY02QOAD8IoAIQXgRQAAAC2oskeAABcIj3zMdPKgNSiQgoAAABbEUgB2IYqCgBA0WQPAKloqiREIxawfDLCjQopAAAAbEUgBQAAgK0IpIAfNMsCABA+BFIAAOCFD+UINwIpAAAAbMUoewAI8+hjqk8A4I0KKQAAAGxFIAUAAJegko9wIpACMdrsbF0AICmEUoQLfUiBIIqL4DeQSNg3AM6jf1P4+4FQo0IKxAgqHQCASEUgBaJMapviqXwAAOxGIHV4xSu1F8QOK2jSXxRAsn8rBvH8wH70IY1iVL7AawAA4AQEUgAmuNpVQQ/kcQnWABDdaLIHAACArQikAAAgWYxBQKgRSAEYNIsDAOxCIAUAAICtGNQEADEiuWZXKuShf655joGkUSEF4Mb69kBwED6B1KFCCgBAkDD4B0gbKqQAAIQYqyEBySOQAgAAwFYEUgAAANiKPqQAHNUvj8EioX+OASDcoiqQnjt3Tp5++ulLru/SpYvUrVvX67off/xRFixYINmyZZNbb71VqlatGsY9BSIToQQILvqOAjHYZK+B9I033jBfixYt6r5kz57da7tBgwbJzTffLMeOHZP169fLFVdcIbNmzbJtv4FAw2JyFwAAnCqqKqSWbt26SYMGDfzetmXLFhkyZIh88cUX0qFDB3NdgQIF5N5775WbbrpJMmbMGOa9RSTSZmFCHgAA4RFVFVLLpEmT5LnnnpOxY8fK0aNHvW7TSmiuXLmkXbt27ut69+4tu3fvlpUrV9qwtwAAALEt6gJp3rx55fTp0+b7UaNGmb6hv/76q/v2DRs2SOnSpSVTpv+KwxUqVDBfN27c6Pc+z549K8ePH/e6AJGGii4AwKkiusleg+Xzzz+f7DY1atQwFU6VNWtWWbt2rZQsWdL8/OKLL0rr1q2lX79+7lCq95knTx6v+9CKqTbVnzp1yu9jDBs2TAYPHhyko0I0s3P9akafA5GFAU1AlATSuLg4MygpOfnz53d/nzlzZncYtX6/V69e0rVrVzl58qQJnnrxbcY/ceKEJCQkSO7cuf0+xsCBA6V///7un7VCWqpUqXQcGRBc0R5G46L0eKhqA4ADAqmOjn/sscfSdR8XLlwQl8tlRt6ratWqyfjx400zvFZUlY60V5dffrnf+9DtrG2dxBWDb/QAAMB5oqoP6W+//ebVv1ND5+jRo6VOnTpmJL3SwUwaTidOnOjeTrepWLGiXHnlleIkcUG4AAAA2C2iK6SptW/fPjPlU82aNSVfvnxm4nttxp86dap7G23S17lKH3jgAZk3b54cPnxYfvnlF5k9e7Zp4gfSg1cQAPO3YBDPA5AacS5tz44iOtm9BtEDBw5I+fLlpXnz5l4j6i3r1q2TxYsXm6Z4nX80pb6qnrQKq6P59bF8B0hJGEMt4Se0XBF2TlwhfuxI/kMQra/1SH7OEVuB1DUouP8PfaNFsu+bQLRVSJW+4Nu3b5/idtpfNKk+owAAAAifqOpDCgAAAOeJugopACC6xDms+4PTmuuBSECFFAAAALYikAJRIpIqREAwMD0dEDsIpACieiQ7ACDy0YcUAOB4ofhARasDED4EUgBAulFh//d5YEATkCY02QMAAMBWBFIgiqS3iZFBJAAAOxBIAQAAYCv6kAKwFX0PES3oPwqkHRVSAAAA2IpACkQZpqoBADgNgRTAJWhGBwCEE4EUAAAAtiKQAgAQ4wOanL7/cD5G2cdQX8G4KJsvEwAARAcCKQAgzeJi5NgYLAiEFk32MYQ/qLEjveea1woA/i4gnAikAAAAsBVN9gDcqIzGnmhucg/r85TUoCAGCwEBoUIKAAAAW1EhBaK42pma6hfVUQAWKucINyqkAAAAsBWBFLbg03dkoToKALATgRQAAAC2IpACAADAVgRSIMbRXA8AsBuBFIhihE3wGgHgBARSIIYRWAEAkYB5SGOMvwDCiHcAQJznqlKsMIUwI5AiqFWyuAid9B3+n0MAACIBgRRBRciJ3HNCgAcARCr6kAIxiA8OAIBIQoUUAEKI8A8AKaNCCsQYAhIAINIQSAEAAGArAikAAABsRR9SIEbQVM/zhvD/n2F2CyAwBNIoxh9CAKEMXvyNSeH5YXJ5IGAEUgBIIXRRXfaPBSqiF695hBt9SAEAAGArAikAAABs5bgm+02bNsmUKVOkQIECcu+99/rdZvny5bJw4ULJli2btGvXTsqXL5+mbQCEHs2+QGT2f3XRBxZh5JgK6YULF6RVq1bSunVrmTp1qowZM8bvdsOGDZMWLVrI9u3bZdmyZVK9enWZO3duqrcBACCtGNAEpE6cy+VyRN/lhIQEU9Fs2bKlPProo7JkyRJZuXKl1zbbtm2TypUry2effSadO3c21z344IMyc+ZMc1uGDBkC2iYlx48fl7x588qxY8ckT5483jfGRc6408jZE2dyRdhz7oj/qFH6Wo3m5z7az51d5z0aAmm6KqQ+0SLZ903ASRXSjBkzmgppXDKBT0Nljhw5pH379u7r7rjjDtmxY4c7vAayDQAAAMLHMYE0EBs2bJAyZcpI5syZ3ddVqlTJfVug2/g6e/as+XTneQEAAIDDBzWdPn1ahg4dmuw21apVk27dugV8n6dOnbqkKSBXrlymuqq3BbqNL+1zOnjw4ID3AwAAAA6okGrTu45wT+7iWcUMhAZL7Z/i6cSJE6b/qd4W6Da+Bg4caH7HuuzcuTPVxwsAAIAIq5Bmz55dnn322aDeZ9WqVWXChAly7tw5yZIli7lu48aN7tsC3cZX1qxZzQVAZGGwESJRNAxoAsItqvqQtm3bVuLj4808pZaPPvpIypUrJ3Xq1Al4GwDhD5ZpuQChfE0CCB9HTYw/cuRI2bdvn/zyyy/yzz//uCuszz//vKl26mAl7e959913y4IFC+Tw4cMyf/58M7Lems4pkG0AAPYGt9ROJ0WABJzNUYFUm821b2mbNm28rvecCmrAgAFm0vtFixaZ7d99910pVaqU1/aBbAMAsA8reAGxxTET40cSJsaPDUyMD9grNVXSSHoji5Y+pEyMj3ByVIUUABA7IilkxloYBcKNTpMAAACwFYEUAAAAtiKQAgAAwFYEUgAAANiKQU1AiAdQMH0NAADJI5DGSMBK7STTSP9zDiC2MMIeSDsCaYwgKAEAgEhFH1IAAADYikAKAAAAWxFIAQAAYCsCKQAAAGxFIAUAIJ0YYQ+kD6PsgTBglgMAAJJGhRQAAAC2IpACAADAVgRSAAAA2IpACgAAAFsRSAEASAdG2APpRyAFAACArQikAAAAsBWBFAAAALYikAIAAMBWBFIAANKIAU1AcBBIAQAAYCsCKQAAAGxFIAUAAICtCKQAAACwFYEUAAAAtspk78MDQMrikrjexZMHGzHCHggeAikAALFuUOAfBAPBh0WkFk32AAAAsBWBFAAAALYikAIAAMBWBFIAAFKJAU1AcBFIAQAAYCsCKQAAAGxFIAUAAICtCKQAAACwFYEUAAAAtiKQAgAAwFYsHQoAQICY7gkIDQJpkPmu/ct6vkBopWe97VDj/z8ARGkg3bVrl3z11VeSL18+6dWrl9dt58+flzfeeOOS32ndurXUrFnT67q1a9fKokWLJFu2bOb24sWLh3zfAQAA4OA+pAkJCXLrrbdKo0aN5MMPP5SRI0dess3Zs2dl4MCB8tdff8nRo0fdl3Pnznlt9/bbb8vVV18ty5cvl+nTp0vlypVl8eLFYTwaANFeHUV4XwdJXQA4g2MqpC6XS26//XaZOnWqPPbYY7JkyZIkt73vvvukQYMGfm/buXOnPPHEEybUWhXWO++8U/r16yebNm2SuLjg/glL6t5oygOA4P9tTW47/u4CkcsxFdJMmTKZCql+Tck333xjqqCzZ8+W+Ph4r9tmzJghWbJkka5du7qvu/vuu2XLli3y66+/hmTfgWirOoX6AqT0mkzrazo9GNAEhI5jAmmgsmbNKmvWrJH169dL//79pUaNGuZ7y7p166Rs2bImlFqqVKnivs0f7Qpw/PhxrwvsR5iJ3ucxkvYFkSFYrwNeU0Bksq3J/syZM6aKmRzt29m+ffuA71ND5u+//y5Vq1Y1P2vf0VatWknfvn3l559/NtedPHlS8ubN6/V7uXPnlowZM5rb/Bk2bJgMHjxYoukPtCsC9oPmMwCh/BuT0n3yNwiIHBns7BPqOfDI3+XUqVOpuk8NpFYYtX6+6667ZNmyZXL69GlzXc6cOS+pcGoQ1UFTeps/OlDq2LFj7ov2QwUAAIDDK6Q5cuSQ4cOHh/xxtPJ58eJFE0j1MbXqOnHiRLlw4YK7P+rmzZvNV70tqW4AegEAhBddNoDYEFV9SLWvqPb3tGgQHTt2rFSrVk0KFixormvbtq2pvH799dfu7T799FMpWbKk1KtXT5yE5qbQYgAOEN1hlLALRA7HTPukxowZIwcOHJCVK1fK3r173RXWAQMGSObMmc38ox07dpQmTZqYifO/++472bNnj5lI31K+fHl54YUX5I477pAff/xRDh8+LF9++aVMmzZNMmSIqnwe04MWACCYGGEPhFacSztzOsRbb71lgqivF1980T1q/p9//jHTPWlw1fCpFVF/fUN1kNPChQtNU7xOJ1WpUqWA90P7oOrAKO1PmidPHq/bAp3HNFhPupPDlyuKjiVQkfifLRaed7tE4vl2irgIO09JBlK9PqnbnCTIx+AbLZJ73wQcF0gjBYE0OAikkYFAGjqx/sc1Lh3Pgx2vS1dqA2mg1zkBgRQ2c1STPQAgctHtxifYOTWcAjYgkAIAQhpCY3bOT8IpEDBG8QAAUi2WVjwKuLk+OVRLgWQRSAEAqZKWIBor4RVA2tBkH2ThHiMWTc1g0XQsTsLzjkh7zUTaa9Lv/rwg0SXajgeOQ4UUAAAAtiKQAgAAwFYEUgAAANiKQAoAAABbMagpHQOXdMUmAACQPOv9ksUhkRQCaRqcOHHCfC1VqlRafh0AgJh9/9Q17QFfrGWfBhcvXpR//vlHcufOLXFxcbZ+4tRQvHPnTsmTJ49EE47NmThvzsR5cyYnnTetjGoYLV68uGTIQG9BXIoKaRrof6aSJUtKpNA/RJH+xyitODZn4rw5E+fNmZxy3qiMIjl8TAEAAICtCKQAAACwFYHUwbJmzSovvPCC+RptODZn4rw5E+fNmaL5vCH2MKgJAAAAtqJCCgAAAFsRSAEAAGArAikAAABsxTykDqATCk+fPl2+/vprMyF/xYoV5cEHH5Tq1at7bXfkyBEZNmyYrFixQvLnzy/9+vWTm266KdXbhNuhQ4dk7Nix5hivueYaeeWVV7xu3759u3Tu3PmS3xs5cqTUq1fP/fOxY8dk+PDhsmzZMjPfXZ8+faRdu3ZiJ32+x40bJ9OmTZOrrrpK3nzzzUu20cmiX331Vfn555/NYgu33367dOjQIdXb2O2NN96QqVOnel1XtmxZmTx5std1S5YskXfffVf27t0rNWrUkIEDB0qJEiUkUsXHx8uIESNkwYIFki1bNvNa1OffaT744AP59NNPva677LLL5JtvvvG6buXKlfLWW2+ZydarVq0qTz31lJQrV04iyZkzZ+SLL76QiRMnSq5cuczfRl8XLlwwfyO+/fZbyZgxo/lbcPfdd3tNyh7INuF29uxZ8/9owoQJkilTpkvOj2rVqpV7xUDLPffcI71793b/nJCQIKNGjZLZs2ebn2+++Wa57777zHECkYhA6gCPPfaY/P3333LLLbeYVS6++uorqVOnjvz4449Sv359s8358+fl2muvlRw5cpg3kL/++sv8cZ00aZLcdtttAW8Tbvqm16BBA/Mmr8F7w4YNft98li9fLjNnzpRChQq5r69SpYrXH9/rrrvOfP/MM8/Ipk2bpGPHjuYNuEePHmKH/fv3y5VXXmn2Q99Y1q1bd8k2esz6geDUqVNmtKyG727dusl7771nPiwEuk0k2LZtmxnt6/mBInv27F7b/PTTT9KiRQvzmtZQp8G0UaNG8scff0TspNldunQx5+7ll182HzDuv/9+88FQ/w85if5fO3funAkplixZsnht89tvv0njxo1NKOvevbuMGTNGGjZsaM5PkSJFJFLUrFnTfHjVvwdLly71u40ew/fffy+vvfaaOe7+/fvLli1b5PXXX0/VNuGmf9tr165tnm/dN3+0oPDkk09K8+bN3df5LmX90EMPmQ/Ceiy6ouCjjz5q/r7q/zkgIrkQ8U6cOHHJdQ0bNnT17NnT/fP48eNdmTJlcu3fv9993f333++qWLFiqrYJt7Nnz7rOnDljvu/QoYOrXbt2l2yzbt06l75Ud+7cmeT9TJkyxZUhQwbXrl273NcNGDDAVapUKdfFixdddjh37pzr9OnT5vvu3bu7rr/++ku2mTlzpjm2rVu3uq975plnXEWKFHFduHAh4G0igb6W/J0/T02aNDHn2aLnPl++fK7hw4e7ItHSpUvNc79ixQr3dW+++aYrZ86crpMnT7qcRF8zTZs2TXabtm3bulq0aOH++fz5866SJUu6nnrqKVckOXr0qPn60ksvucqUKXPJ7Rs3bjTnbc6cOe7rJkyYYP7+7d27N+Bt7Dy21157zfwf9ydv3ryuqVOnJnkff//9t/l7+PXXX7uv0+31uh07doRgr4H0ow+pA2iTlL/r9BO9RZsTtdLoWUHU6ufmzZtNRS3QbcJNKzTaDBqIO+64w1TXtNnJt9qox6aVBc+mXz02rQpt3LhR7JA5c+ZLKoS+dL+12dqzSVT3e9++fbJ27dqAt4kUq1atMlX4W2+91TTha/OjZ9O3djlo06aN+zo991rZnj9/vkQife6LFi0qdevW9XrutVqtXUOcZv369dKyZUtp27atDB06VE6fPu2+TSvxCxcu9Do/WtnX6nyknZ+UqulW9wo9Vs/zpk30ixcvDngbOwTaUqCVz2bNmkmvXr0uqaQuWrTIdDu44YYb3Ndpk71WSvUcA5GIQOpA+kaof0y1Cd+igVKb8z1ZP1thM5BtIpWGHO0Tqk292oR/xRVXmC4LFqceWzSdNw3f2j1Cm7L1jX306NHStGlT8wavdu/ebbpW+DuWSDoOT/6ee+tDT6Tuc3If/rQZXv8PaRcZ7apz9dVXmw8K6vDhw3Ly5ElHnZ+k6P4WLlzYBGqL9r3WD/Ke/69S2iZSad9e/Xuo3ZNKly5tPmB49k/X/S9QoIDXh339XscNRPqxIXbRh9QGf/75p/Tt2zfZbfSNQwcu+dK+pDqYRfu16cWi/UN9V+uwqnN6W6DbpJf23ezZs2ey22h/1QEDBgR8n+XLlzcVGv10r2688UbTl++JJ55wV6l0//XNJJTHltTgKk8axHSQTqAi5bz5o4PfZsyYkew2OrCkTJky5vshQ4Z47aeG0cqVK5tBTRpUrX31dyyhPI708Pfca+Vbq0+Rus9J0T6HnseirQ06QPLjjz+WBx54wJHnJzXnzfdYAtkmUv3www/ufdcBTjlz5pSnn35a7r33XhM8nXxsiF0EUhvoG7iOYk2Ob5VC7dixw1QKdSCTjtz2pJ+GtcLhO3rdGkkb6DbppdWjlI5Nm0BTw3fghdJmtscff9z9cziOTQcZpHRsWnFJDd3vrVu3pnjeUtomFDR8ew6a8MdzoIvvG6B2MahQoYIZEKOBVI9D+TtPoTyO9PD3ujp69KhcvHgxYvc5Kb7nR/8falcQPT8qX758Jmg76fyk5rxplwS9Lrm/h77bOOVc6t9D/SCsg5a09cjfsTn1XCJ2EEhtoE1C2pczNbQvpIaDWrVqmaqUZzOT0v6TOq2L/kG1Kok6Ml0/EVuj0QPZJr10BH9qjy0tdMogrQpY9Ni0T5U2D1vPjR6bVrN8p8dKK608BPvYdL91uivta2m9yeh+69QsOpI40G1CQSvTekkrbZ4/cOCA+zxpWNcPLDpC2LO7iR6L5/RdkUSfe50WSCvy2txp7a/SGRScTv8fWa9pfX1ffvnl5vx4Th+kx+u0Y9XzpuFLZ36w+l5r/2Z9TVrHEsg2TjqPyvq/psem/Zy1r72eU6X9zbXPsNOODTEkCAOjEGI6crxChQquW265xYzc9mfz5s2uzJkzu95//33z8+HDh12VK1d29e3bN1Xb2CmpUfafffaZa/369e6fV61aZUZm66hui44czZYtmxkBrY4dO+aqXr26Gd0eCZIaZb9nzx4zYnvYsGHuGRVq167tNRI9kG3spqP9hw4d6oqPj3f//MQTT5gRy2vWrHFv9/zzz7uKFi3qHuk7ffp0V1xcnBnNHomOHz/uKliwoJmxwZoVQkeqpzRaPRLp+Tl16pT5XmeeGDJkyCXP/YgRI1z58+d3bdiwwfy8YMECMzLbcyR6JElqlL2ep7Jly7p69+5tfk5ISDB/W2rWrOmedSOQbeyU1Cj7RYsWuebPn+/+effu3ebvQd26dd3X6f8//dvetWtXcyx66dSpk7kukmbmADwRSB1A/6joZ4c6deq46tev775Yf0gtkydPduXJk8dVvnx5V44cOVytWrUywSy124Rby5YtzfEUKFDAvBnq93qdRd8w9dhLly5t/qBmzZrV9dBDD7mni7J89dVXJqiWK1fOlStXLlezZs1M6LbTTTfdZI5HQ41O1aLf69RHnmbNmmWOXd8cc+fO7WrUqJHrwIEDqd7Gbi+88II5zho1arguu+wy8xqbPXu21zYaArp06WI+PFSqVMmVPXt21zvvvOOKZD/88IOrWLFiZvojfX3qm//27dtdTvPqq6+6Chcu7KpWrZr5qsej06V50rCiH1D1/5j1f+3ll192RZp7773X/F/SY8iSJYv7b+I///zj9cFV/7/oB6BChQq5qlSpYqaQ8xTINuH2yCOPmGPRv3daQLCOTadyUvrau/nmm83+6oduPUdt2rS5ZFq81atXmyn9NNTq+dbv9TogUsXpP3ZXaZHyQCGrz6Bv07/2AfOkI9C1H5E2L1qDTXwFsk046cow1khsiza7e061Y43SPn78uGleS2qqKB0xrFPb6NQpkbC6zK+//uo1PZfSfno6utmTNsfrfuvArKSayQPZxm46YEKn2dL+iNoP2uoa4ksnltdpq3RQje9gtEikr09t/tQuEzpQy6n0OPTvif7tKFmyZJLnR8+NniP9P6TnMtLooh76t8CXNkd79q/Uvr563vT/nI5M93e8gWwTTvp/XPsp+9LJ8j3/7mk3kl27dpm/4Xny5PF7X3psen9Kj83OFaiAlBBIAQAAYCs+LgEAAMBWBFIAAADYikAKAAAAWxFIAQAAYCsCKQAAAGxFIAUAAICtCKQAAACwFYEUgKP99NNPMm/evEuu3759u0yePNms6Q0AiGwEUgCOpqvR3HjjjbJ48WKvFaNuvfVWmTFjhuTMmdPW/QMApIyVmgA43qOPPmrC5+rVq82ymE899ZRMmjTJ/ByJS18CALwRSAE4Xnx8vFnHvEmTJtK9e3dp0aKFfP/999K8eXO7dw0AEAACKYCo8L///U8aNWpkKqK9e/eW1157ze5dAgAEiEAKIGpohXTp0qWya9cuKVq0qN27AwAIEIOaAEQFHVGvVdKKFSvKgAED7N4dAEAqUCEF4Hg7d+6UWrVqyQsvvCDXXnut1KtXzwxq6tChg927BgAIAIEUgKO5XC4ziCljxoxmIFNcXJwMHTpU3nrrLVm7dq0ULlzY7l0EAKSAQArA0XTw0vDhw2XNmjVSvHhxc11CQoJcc801UqJECZk2bZrduwgASAF9SAE41smTJ00VdPz48e4wqrRaOm7cOMmSJYusW7fO1n0EAKSMCikAAABsRYUUAAAAtiKQAgAAwFYEUgAAANiKQAoAAABbEUgBAABgKwIpAAAAbEUgBQAAgK0IpAAAALAVgRQAAAC2IpACAADAVgRSAAAA2IpACgAAALHT/wOa/PFUeTszCQAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "random.setSeed(123) # Make results reproducible\n", "txBearingAngle = -135 # Pointing TX panel to the center of map\n", "\n", "bwp = BandwidthPart(numRbs=24, spacing=15) # Create a bandwidth part\n", "\n", "# Create a trajectory passing through the following points\n", "trjPoints = [[35, 90], [57, 89], [76, 86], [91, 83], [107, 78], [118, 70], \n", " [128, 62], [139, 49], [147, 38], [149, 22], [151, 2], [154, -14]]\n", "trajectory = dmData.trajectoryFromPoints(trjPoints, bwp, speedMps=15)\n", "\n", "trajectory.print() # Print the trajectory information\n", "ax = dmData.drawMap(\"LOS-NLOS\", trajectory) # Draw the map with the trajectory\n", "dmData.drawBsPanel(ax, txBearingAngle) # Draw the base station antenna panel" ] }, { "cell_type": "code", "execution_count": 4, "id": "490ba356-c416-4413-82f9-cfdeb413c384", "metadata": {}, "outputs": [], "source": [ "def createCsiRS():\n", " # Create a CSI-RS configuration for beam sweeping and beam probing\n", " # Beam-sweeping:\n", " # - 8 beams and 8 CSI-RS objects with resource IDs 1 to 8\n", " # - CSI-RS resources at symbol 4 and resource elements 2 to 9\n", " # - One periodic NZP CSI-RS resource set (set ID=1) containing all 8 sweeping resources\n", " # with period 20\n", " # - One periodic CSI report (ID=11) with period 20 and offset=1\n", " numSweep = 8\n", " sweepResources = []\n", " for i in range(numSweep):\n", " sweepResources += [ CsiRs(resourceId=i+1, numPorts=1, symbols=[4],\n", " freqMap=\"\".join([str(int(x)) for x in np.eye(12)[i+2]])[::-1]) ] # REs 2 to 9\n", " sweepSet = CsiRsSet(\"NZP\", bwp, resourceType=\"periodic\", rsId=1, period=20, csiRsList=sweepResources)\n", " sweepReport = CsiReport(sweepSet, reportId=sweepSet.rsId+10, reportType='periodic', \n", " period=20, offset=1, quantity=\"Cri\") # CSI report for beam sweeping\n", " \n", " # Beam-probing:\n", " # - 4 beams and 4 CSI-RS objects with resource IDs 9 to 12\n", " # - CSI-RS resources at symbol 4 and resource elements 4 to 7\n", " # - One aperiodic NZP CSI-RS resource set (set ID=2) containing all 4 probing resources\n", " # triggered when sweeping CRI is received.\n", " # - One aperiodic CSI report (ID=12) triggered when sweeping CRI is received.\n", " numProbe = 4\n", " probeResources = []\n", " for i in range(numProbe):\n", " probeResources += [ CsiRs(resourceId=i+9, numPorts=1, symbols=[4],\n", " freqMap=\"\".join([str(int(x)) for x in np.eye(12)[i+4]])[::-1]) ] # REs 4 to 7\n", " probeSet = CsiRsSet(\"NZP\", bwp, resourceType=\"aperiodic\", rsId=2, csiRsList=probeResources)\n", " probeReport = CsiReport(probeSet, reportId=probeSet.rsId+10, \n", " reportType='aperiodic', quantity=\"Cri\") # CSI report for beam probing\n", " \n", " csiRsConfig = CsiRsConfig([sweepSet, probeSet]) # CSI-RS config\n", " csiReportMan = CsiReportMan([sweepReport, probeReport]) # CSI report manager\n", " return csiRsConfig, csiReportMan\n", "\n", "def processFeedback(csiReportMan):\n", " # Process CSI feedback:\n", " sweepReport, probeReport = csiReportMan.csiReports\n", " criSweep, criProbe, rsrpProbe = None, None, None\n", "\n", " csiReportInfo = csiReportMan.getFeedback() # Get all available CSI reports from CsiReport objects\n", " for reportId, csiFeedback in csiReportInfo.items(): # Get the CSI feedback for each report\n", " if reportId == sweepReport.reportId: # Sweeping report\n", " criSweep = csiFeedback.cri.cri # CSI-RS resource ID of the best beam\n", " probeReport.csiRsSets[0].trigger() # Trigger Probing CSI-RS resource set\n", " probeReport.trigger() # Trigger Probing report\n", " \n", " elif reportId == probeReport.reportId: # Probing report\n", " criProbe = csiFeedback.cri.cri # CSI-RS resource ID of the best beam\n", " rsrpProbe = csiFeedback.cri.rsrp # RSRP of the best beam\n", " else:\n", " print(f\"Unhandled report: {reportId}\")\n", "\n", " return criSweep, criProbe, rsrpProbe" ] }, { "cell_type": "code", "execution_count": 5, "id": "66a4df45-1640-4486-a561-4807173247b4", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Simulating communication over 12782 slots along the trajectory ...\n", " Slot Blk Errors BLER Beam Ang RSRP (dB) Trj. Time Exe. Time\n", "------ ---------- ------ -------- --------- --------- ---------\n", "12,782 0 0.00 38.21 -55.59 12.78 1560.99 " ] } ], "source": [ "random.setSeed(123) # Make results reproducible\n", "numSlots = trajectory.numPoints # Total number of slots = number of points on the trajectory\n", "\n", "# Note:\n", "# Since our goal is to show that different UE locations along a trajectory experience different \n", "# beam RSRPs, we should not normalize the channel gain independently at each location. So, we set \n", "# the 'normalizeGains' to 'False' in our channel model below and use a fixed noise power.\n", "\n", "# Calculating the noise power:\n", "k = 1.380649e-23 # Boltzmann constant (joules per kelvin)\n", "tempK = 290.0 # Temperature in kelvin\n", "nf = 9 # Receiver noise figure (dB)\n", "noiseVarFreq = k*tempK*bwp.spacing*1000*toLinear(nf)\n", "\n", "# Creating a trajectory-based channel model:\n", "channel = TrjChannel(bwp, trajectory, \n", " normalizeGains = False,\n", " txAntenna = AntennaPanel([1,4], polarization='x'), # 8 TX antennas\n", " txOrientation = [txBearingAngle,0,0], # BS antenna orientation\n", " rxAntenna = AntennaPanel([1,1], polarization='x', # 2 RX antennas\n", " beamWidth=[65,360])) # Omnidirectional\n", "\n", "csiRsConfig, csiReportMan = createCsiRS()\n", "# csiRsConfig.print() # Uncomment to print the CSI-RS configuration\n", "# csiReportMan.print() # Uncomment to print the CSI reporting information\n", "\n", "# Create a simple PDSCH object for end-to-end transmission\n", "pdsch = PDSCH(bwp, numLayers=1, modulation=\"16QAM\", csiRsConfig=csiRsConfig)\n", "pdsch.setDMRS(additionalPos=1) # DMRS configuration\n", "\n", "ldpc = pdsch.getLdpcCodec(coderates=490/1024)\n", "\n", "minMse, maxMse = 100, 0\n", "t0 = time.monotonic()\n", "print(f\"\\nSimulating communication over {numSlots} slots along the trajectory ...\")\n", "print(\" Slot Blk Errors BLER Beam Ang RSRP (dB) Trj. Time Exe. Time\")\n", "print(\"------ ---------- ------ -------- --------- --------- ---------\")\n", "\n", "txBlockErrors = 0\n", "sweepWs, sweepBeams = channel.txAntenna.getSweepingBeams(numTheta=1, numPhi=8)\n", "curSweepCri, curProbeCri, curProbeRsrp = None, None, None\n", "probePhiLocal, probePhiGlobal = None, None\n", "channel.restart()\n", "blefs = np.zeros(numSlots, dtype=np.int8) # Block error flags for each slot\n", "beamAngles = np.zeros((numSlots,2), dtype=np.float32) # Beam angles for each slot (theta, phi)\n", "beamRsrps = np.zeros(numSlots, dtype=np.float32) # Beam RSRPs\n", "for slotNo in range(numSlots):\n", " \n", " criSweep, criProbe, rsrpProbe = processFeedback(csiReportMan) # Process CSI feedback\n", " if criSweep is not None: curSweepCri = criSweep # Update current Sweeping CRI \n", " if criProbe is not None:\n", " curProbeCri = criProbe # Update current Probing CRI\n", " curProbeRsrp = rsrpProbe # Update current Probing RSRP\n", " probePhiLocal = probeBeams[1][criProbe-9]\n", " _, probePhiGlobal = AntennaPanel.local2Global( probeBeams[0][criProbe-9],\n", " probePhiLocal,\n", " channel.txOrientation)\n", " \n", " pdsch.initGrid() # Initialize PDSCH's internal resource grid \n", " numBits = pdsch.getBitCapacity()[0] # Capacity of the resource grid for PDSCH data\n", " txBlock = random.bits(ldpc.txBlockSizes[0]) # Create a random transport block\n", " \n", " rateMatchedCodeBlocks = ldpc.encode(txBlock, numBits) # LDPC rate-matching and encoding\n", " pdsch.setPdschData(rateMatchedCodeBlocks) # Populates the PDSCH's internal grid. \n", " \n", " channelMatrix = channel.getChannelMatrix() # Get the channel matrix\n", " precoder = pdsch.getPrecodingMatrix(channelMatrix) # Get the precoder matrix from PDSCH object\n", "\n", " txGrid = bwp.createGrid(len(channel.txAntenna)) # Create the transmitted resource grid\n", " pdsch.precodeTo(txGrid, precoder) # Perform the precoding and put data in `txGrid`\n", "\n", " # Process CSI-RS:\n", " csiRsResources = csiRsConfig.getResources() # Get all CSI-RS resources for the current slot, if any\n", " for csiSetId, setResources in csiRsResources.items():\n", " if csiSetId == 1: # Beam Sweeping:\n", " for resourceId, (lIdx, kIdx, sweepReValues) in setResources.items():\n", " b = resourceId-1 # Sweeping resource ID to beam index\n", " w = sweepWs[:,b:b+1] # Sweeping vector (nt x 1)\n", " txGrid[:,lIdx, kIdx] = (w @ sweepReValues, # Update txGrid with precoded CSI-RS values\n", " \"CSIRS_NZP\", # Set RE's type to NZP CSI-RS\n", " resourceId) # Set RE's object ID to CSI-RS resource ID\n", " \n", " elif csiSetId == 2: # Beam probing:\n", " b = curSweepCri-1 # Beam index of current sweeping CRI \n", " theta0, phi0 = sweepBeams[0][b], sweepBeams[1][b] # Best beam angles derived from Sweeping CRI\n", " probeWs, probeBeams = channel.txAntenna.getProbingBeams(theta0, phi0, 4, polStrategy='equal')\n", " for resourceId, (lIdx, kIdx, probeReValues) in setResources.items():\n", " b = resourceId - 9 # Probing resource ID to beam index\n", " w = probeWs[:,b:b+1] # Precoding vector for this beam. nt x 1\n", " txGrid[:,lIdx, kIdx] = (w @ probeReValues, # Update txGrid with precoded CSI-RS values\n", " \"CSIRS_NZP\", # Set RE's type to NZP CSI-RS\n", " resourceId) # Set RE's object ID to CSI-RS resource ID\n", "\n", " rxGrid = txGrid.applyChannel(channelMatrix) # Apply the channel in frequency domain\n", " noisyRxGrid = rxGrid.addNoise(noiseVar=noiseVarFreq) # Add noise\n", "\n", " csiReportMan.processRxGrid(noisyRxGrid, csiRsResources) # UE-side processing of the received CSI-RS\n", "\n", " effChanMat = channel.getEffChannel(channelMatrix, precoder) # Get effective channel matrix\n", " eqGrid, llrScales = pdsch.equalize(noisyRxGrid, effChanMat) # Equalize the received grid\n", " llrs = pdsch.getLLRs(eqGrid, llrScales) # Demodulate and get LLRs\n", " decodedTxBlocks, crcMatch = ldpc.decode(llrs) # LDCP rate-recovery and decoding\n", " txBlockErrors += 0 if crcMatch[0][0] else 1 # crcMatch[0][0] -> CRC match for txBlock\n", "\n", " blefs[slotNo] = 0 if crcMatch[0][0] else 1 # Block error flag for this slot\n", " beamAngles[slotNo] = np.nan if probePhiGlobal is None else probePhiGlobal # Current beam angles based on probing CRI\n", " beamRsrps[slotNo] = np.nan if curProbeRsrp is None else curProbeRsrp\n", " dt = time.monotonic()-t0 # Total time spent so far\n", " \n", " print(f\"\\r{slotNo+1:^6,d} {txBlockErrors:^10,d} {100*txBlockErrors/(slotNo+1):^6.2f} \" \n", " f\"{\"N/A\" if probePhiLocal is None else str(np.round(probePhiLocal,2)):^8s} \"\n", " f\"{\"N/A\" if curProbeRsrp is None else str(np.round(curProbeRsrp,2)):^9s} \"\n", " f\"{channel.trajectory.cur.time:^9.2f} {dt:^9.2f}\", end='')\n", " channel.goNext() # Move to the next slot and trajectory point\n" ] }, { "cell_type": "code", "execution_count": 6, "id": "63567763-0aae-4a08-8f90-61be9efe3234", "metadata": {}, "outputs": [ { "data": { "text/markdown": [ "![demo](AnimateCRI.gif)" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "arrow, prevRsrp = None, None\n", "\n", "# Since we are averaging the angles, we need to make sure they are continuous:\n", "# Local angle range: -90 to 90\n", "# Global angle range: 135 .. -45\n", "# Adjusted global angle range: 135 .. 315 (Add 360 to the negative values)\n", "beamAngles = np.array(beamAngles)\n", "beamAngles[beamAngles<0] += 360 # Range: 135 .. 315\n", "\n", "# Callback used to initialize and update the scenario map and the graphs below it\n", "def handleGraph(request, ax, trajectory, points=None):\n", " global arrow, prevRsrp\n", " if request==\"Config\":\n", " ax[0].set_xlim(0,trajectory.numPoints)\n", " ax[0].set_ylim(-65,-50)\n", " ax[0].set_title(\"RSRP of the Best Beam (dB)\")\n", " ax[0].set_xlabel(\"Trajectory Points\")\n", " ax[0].grid()\n", "\n", " elif request==\"ConfigMap\":\n", " ax.set_title(\"Beam Tracking Along a UE Trajectory Using CSI Feedback\")\n", " dmData.drawBsPanel(ax, txBearingAngle) # Draw TX antenna panel\n", "\n", " # Create the arrow patch\n", " arrow = dmData.drawBeamArrow(ax, txBearingAngle, color=\"cyan\", length=50)\n", " arrow.set_animated(True)\n", "\n", " elif request==\"Draw\":\n", " # Draw the RSRP plot below the map\n", " p0, p1 = points\n", " segmentRsrps = beamRsrps[p0:p1]\n", " meanRsrp = segmentRsrps[~np.isnan(segmentRsrps)].mean()\n", " if prevRsrp is None: prevRsrp = meanRsrp\n", " ax[0].plot([p0,p1], [prevRsrp, meanRsrp], 'blue', markersize=1)\n", " prevRsrp = meanRsrp\n", " \n", " elif request==\"DrawOnMap\":\n", " p0, p1 = points\n", " segmentAngles = beamAngles[p0:p1]\n", " meanPhi = segmentAngles[~np.isnan(segmentAngles)].mean()\n", " # Now that we have computed the mean angle, we can remove the 360 to get \n", " # the angle back in the correct range\n", " if meanPhi>180: meanPhi -= 360\n", " # Update the arrow direction\n", " dmData.drawBeamArrow(ax, meanPhi, color=\"cyan\", length=50, arrow=arrow)\n", " return (arrow,)\n", "\n", "# Create the animation and display it below\n", "anim = dmData.animateTrajectory(trajectory, numGraphs=1, pointsPerFrame=200, \n", " graphCallback=handleGraph, fileName='AnimateCRI.gif',\n", " lastFrameDur=3000) # Freeze on last frame for 3 seconds\n", "display(Markdown(\"![demo](AnimateCRI.gif)\"))\n", "\n", "# Alternatively, the following code provides better control over running \n", "# the animation. Note that for this method to work, you should not pass a \n", "# 'fileName' to the 'animateTrajectory' function.\n", "# # Increase the animation memory limit for HTML-based animation display\n", "# matplotlib.rcParams['animation.embed_limit'] = 100000000\n", "# anim = dmData.animateTrajectory(trajectory, numGraphs=1, pointsPerFrame=200, \n", "# graphCallback=handleGraph)\n", "# HTML(anim.to_jshtml())" ] }, { "cell_type": "code", "execution_count": null, "id": "c92bdc1e-eb87-42be-8f6a-bd47b9f0f752", "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.12.10" } }, "nbformat": 4, "nbformat_minor": 5 }