JPEG_Detection.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "end_time": "2024-11-12T09:47:54.883955900Z",
     "start_time": "2024-11-12T09:47:54.859142800Z"
    }
   },
   "outputs": [],
   "source": [
    "import math\n",
    "import os\n",
    "import cv2\n",
    "import numpy as np\n",
    "from PIL import Image\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "from scipy.signal import find_peaks\n",
    "from scipy.signal import savgol_filter\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "end_time": "2024-11-12T09:47:54.884953200Z",
     "start_time": "2024-11-12T09:47:54.867148Z"
    }
   },
   "outputs": [],
   "source": [
    "def zero_pad(Y):\n",
    "    mx, my = Y.shape\n",
    "    padded_x = ((mx + 7) // 8) * 8\n",
    "    padded_y = ((my + 7) // 8) * 8\n",
    "    img = np.zeros((padded_x, padded_y), dtype=np.uint8)\n",
    "    img[:mx, :my] = Y\n",
    "\n",
    "    return img\n",
    "\n",
    "def DCT(img):\n",
    "    w, h = img.shape\n",
    "    num_blocks = (w // 8) * (h // 8)\n",
    "    dct_output = np.zeros((8, 8, num_blocks), dtype=float)\n",
    "\n",
    "    for ci, (i, j) in enumerate((x, y) for x in range(0, w, 8) for y in range(0, h, 8)):\n",
    "        dct_output[:, :, ci] = cv2.dct(img[i:i+8, j:j+8].astype(np.float32))  # Ensure input is float for DCT\n",
    "\n",
    "    return dct_output\n",
    "\n",
    "\n",
    "def hist(DC, N=2000, show=False):\n",
    "    hcount, bin_edges = np.histogram(DC, bins=N, density=True)\n",
    "    if show:\n",
    "        plt.figure(figsize=(10, 6))\n",
    "        sns.histplot(x=bin_edges[:-1], weights=hcount, bins=N, kde=True, stat=\"density\")\n",
    "        plt.title(\"DC Coefficient Histogram\")\n",
    "        plt.xlabel(\"DC Coefficient Value\")\n",
    "        plt.ylabel(\"Density\")\n",
    "        plt.show()\n",
    "    return hcount, bin_edges\n",
    "\n",
    "def FFT(hcount):\n",
    "  hcount -= np.mean(hcount)\n",
    "  x = np.fft.fft(hcount.astype(float))\n",
    "  x = abs(np.fft.fftshift(x))\n",
    "  x = savgol_filter(x, 11, 2)\n",
    "  return x\n",
    "\n",
    "def quality(x, N=2000):\n",
    "    ref = np.max(x)*0.16+np.mean(x)\n",
    "    peaks, _ = find_peaks(x, height=ref, distance=15)\n",
    "    uncompressed = False\n",
    "    avg_dist = 0\n",
    "    if len(peaks) <= 1:\n",
    "        uncompressed = True\n",
    "    else:\n",
    "        distances = np.diff(peaks)\n",
    "        avg_dist = np.bincount(distances).argmax()\n",
    "\n",
    "    return uncompressed, avg_dist, peaks\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-11-12T09:50:45.387423700Z",
     "start_time": "2024-11-12T09:47:54.886954100Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "10\n",
      "20\n",
      "30\n",
      "40\n",
      "50\n",
      "60\n",
      "70\n",
      "80\n",
      "90\n",
      "100\n"
     ]
    }
   ],
   "source": [
    "def compress(path, quality=100):\n",
    "    img = Image.open(path)\n",
    "    img = img.convert('RGB')\n",
    "    img.save('temp_jpeg.jpg',quality = quality)\n",
    "    return Image.open(\"temp_jpeg.jpg\")\n",
    "\n",
    "# Generating training data\n",
    "image_number = 700 # number of images\n",
    "quality_peak_list=[]\n",
    "\n",
    "for q in range(10,101,10):\n",
    "    l=[]\n",
    "    print(q)\n",
    "    for i in range(0, image_number):              # each image has 4 versions\n",
    "        file_name = f\"{i + 1:04}x4w1.png\"\n",
    "        path = os.path.join(\"./DIV2K_train_LR_wild\", file_name)\n",
    "        \n",
    "        img = compress(path, q)\n",
    "\n",
    "        img = img.convert('YCbCr')\n",
    "        img = zero_pad(np.array(img)[:,:,0])\n",
    "        img = img.astype(np.float64)\n",
    "        dct = DCT(img)[0,0,:]\n",
    "        hcount, bin_edges = hist(dct)\n",
    "        fft = FFT(hcount)\n",
    "        _, avg_peak_dist,_ = quality(fft)\n",
    "        l.append(avg_peak_dist)\n",
    "    quality_peak_list.append(l)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "end_time": "2024-11-12T09:50:45.903663300Z",
     "start_time": "2024-11-12T09:50:45.389933600Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 25. 108. 267. 254.  32.   7.   0.   1.   0.   0.   0.   0.   0.   1.\n",
      "    0.   0.   0.   0.   0.   0.   0.   0.   1.   0.   0.   0.   0.   0.\n",
      "    0.   0.   0.   0.   0.   0.   0.   0.   0.   0.   0.   4.]\n",
      " [ 25. 108. 267. 254.  33.  12.  13.  50. 160. 319. 147.   0.   0.   1.\n",
      "    1.   0.   1.   0.   0.   0.   0.   0.   1.   0.   0.   0.   0.   0.\n",
      "    0.   0.   0.   0.   0.   0.   0.   0.   0.   0.   0.   8.]\n",
      " [ 26. 109. 270. 255.  38.  12.  14.  51. 160. 325. 160.  67. 211. 338.\n",
      "   50.   0.   1.   0.   0.   0.   0.   0.   1.   0.   0.   0.   0.   0.\n",
      "    0.   0.   0.   0.   0.   0.   0.   0.   0.   0.   0.  12.]\n",
      " [ 33. 110. 271. 256.  38.  14.  14.  53. 160. 328. 164.  68. 218. 355.\n",
      "  130. 247. 288.  35.   1.   0.   0.   1.   1.   0.   0.   0.   0.   0.\n",
      "    0.   0.   0.   0.   0.   0.   0.   0.   0.   0.   0.  15.]\n",
      " [ 38. 114. 273. 257.  38.  14.  14.  53. 160. 328. 164.  68. 219. 355.\n",
      "  136. 261. 361. 259. 309.  57.   0.   1.   1.   1.   0.   0.   0.   0.\n",
      "    0.   0.   0.   0.   0.   0.   0.   0.   0.   0.   0.  19.]\n",
      " [ 45. 118. 275. 258.  39.  15.  14.  53. 160. 328. 165.  68. 219. 355.\n",
      "  137. 261. 368. 275. 381. 277. 307.  53.   1.   1.   2.   0.   0.   1.\n",
      "    0.   1.   0.   0.   0.   0.   0.   0.   0.   0.   0.  23.]\n",
      " [ 58. 127. 281. 260.  40.  16.  14.  53. 160. 329. 165.  68. 219. 355.\n",
      "  137. 261. 368. 276. 385. 287. 347. 178. 297. 183.   5.   0.   1.   1.\n",
      "    1.   1.   0.   0.   0.   0.   0.   0.   0.   0.   0.  27.]\n",
      " [ 73. 141. 290. 264.  44.  17.  14.  53. 160. 329. 165.  68. 219. 355.\n",
      "  137. 261. 368. 276. 385. 287. 347. 178. 298. 188.  16.  42. 146. 305.\n",
      "  135.   1.   0.   1.   0.   0.   2.   0.   0.   0.   0.  35.]\n",
      " [ 92. 157. 297. 271.  52.  19.  14.  54. 162. 330. 165.  68. 219. 355.\n",
      "  137. 261. 368. 276. 385. 287. 347. 178. 299. 188.  17.  42. 146. 307.\n",
      "  136.   2.   6.  13.  62. 204. 233.  32.   0.   0.   0. 119.]\n",
      " [121. 177. 304. 279.  60.  23.  14.  55. 164. 332. 165.  68. 219. 357.\n",
      "  137. 261. 368. 276. 385. 287. 347. 178. 299. 188.  17.  42. 146. 308.\n",
      "  136.   2.   7.  13.  62. 204. 233.  32.   0.   0.   0. 734.]]\n",
      "[1.7080502  2.81304283 2.91803545 3.02302808 3.12802071 3.23301334\n",
      " 3.33800596 3.44299859 3.54799122 3.65298384 3.75797647 3.8629691\n",
      " 3.96796172 4.07295435 4.17794698 4.28293961 4.38793223 4.49292486\n",
      " 4.59791749 4.70291011 4.80790274 4.91289537 5.01788799 5.12288062\n",
      " 5.22787325 5.33286587 5.4378585  5.54285113 5.64784376 5.75283638\n",
      " 5.85782901 5.96282164 6.06781426 6.17280689 6.27779952 6.38279214\n",
      " 6.48778477 6.5927774  6.69777003 6.80276265 7.90775528]\n"
     ]
    },
    {
     "data": {
      "text/plain": "<Figure size 1000x600 with 1 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA08AAAIhCAYAAACWt4GEAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAACCo0lEQVR4nO3dd3RU1f7+8efMpCeEQBICBKRLTwgdgUv10kSaBUQQUQEB8YoKckGRrnTpRSxgARFFVGygIipFIwQRUJrcUARCCDUhZeb3B7/Ml5AEJjMTZhLer7VYi5wze89nzpwpz5x99jGsVqtVAAAAAIAbMrm7AAAAAAAoCAhPAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB0ITwAAAABgB8ITAAAAANiB8AQAAAAAdiA8AU7gGtP/x93bwt33D9fgecxfbF/k5FbuF+yDKOgIT3BYnz59VLVq1Sz/atWqpZYtW2rcuHE6d+6cy+9z27Ztqlq1qrZt2+Z0X3369FGfPn1sf7du3VovvPCC3e03btyokSNHOl2Hux09elRVq1bVRx99JCnv2zg1NVWTJ0/Wp59+etPbVq1aVXPnznXofm4kNjZWAwYMsP19/WPKb59//rlatWqlWrVq6aWXXsrX+/r2229VtWrVHNf9+OOP6tGjh6Kjo9W6dWstW7bMo76oXP+cz507N8tjuf55dMV93ErXv6c4Y9CgQVq9erVL+sq0YMECLVu2zKV9OutWvlYdeZ1ev4/m9XOiIFi9erVeffXVPLVx5HWW02fF9a+Z/PqMcMS1tTjjo48+UtWqVXX06FEXVPV/zp49q5YtWyo+Pt6l/eLmvNxdAAq2GjVqaOzYsba/09LS9Mcff2jmzJnau3ev3n//fRmG4cYK7Tdv3jwFBQXZffu33nor/4pxo5o1a2rVqlWqXLmyXbc/deqU3n77bU2ZMuWmt121apVKlizpbInZrF69WgcPHrT9XaJECa1atUp33HGHy+8rJ+PHj1f58uX1yiuvKCIiIt/uZ9u2bXr22WdzXLdz504NGjRIHTp00NNPP63Y2FhNmzZNGRkZTgeS/HL//ferefPmtr+vfx5vVx999JFOnjypHj16uLTf1157TUOHDnVpn866la/VW/U6LWgWLlyohg0b5vv92PNZkV+fEe7UsmVLrVq1SiVKlHBpv8WKFVO/fv303//+V8uXLy8w37UKA8ITnBIUFKQ6depkWdagQQNdunRJc+bMUVxcXLb1nqpGjRruLsEj5PScusqt2hd8fHxu6X6XlJSkpk2bqlGjRvnS/8WLF7V06VItXbpURYoU0eXLl7PdZu7cuapevbqmTZsmSfrXv/6l9PR0LVq0SH379pWfn1++1OaMkiVLFrovSs5KSUnR9OnTNXbsWJlMhX9wyK18reb36xTOKyjfF/KiePHiKl68eL70/dBDD2nhwoX65ptv9O9//ztf7gPZFf53ZrhFrVq1JEnHjx+3LduwYYO6d++u2rVrq2nTppo4cWK2L4EbNmzQQw89pJiYGNWqVUvt27fXu+++m+v9pKamqn///mrUqJH27t2b6+2OHz+uoUOHql69emratKnefPPNbLe5fjjGZ599pnvvvVdRUVFq3LixnnvuOZ08eVLS1aEG27dv1/bt27MMK9i3b5+GDh2qxo0bq2bNmmrevLkmTpyolJQUW79Vq1bVu+++q9GjR6thw4aKiYnR008/rYSEhCz1rF27Vt26dVN0dLRatmypGTNmKDU11bb+r7/+0sCBA1W3bl3VrVtXQ4YMsevw/ddff217XN26ddO+ffuyrL9+qERKSopefvll/etf/7I9J5lDf44ePao2bdpIkkaNGqXWrVtLkl544QU98sgjGjt2rOrWrauOHTsqIyMjx2EQBw4c0EMPPaTatWvr7rvv1ooVK2zrchvS88ILL2S5r48//ljHjh2z3Tandn///beGDRumpk2bqk6dOurTp49iY2Oz3dcXX3yhYcOGKSYmRg0bNtSYMWNyDCvXbitJmj9/fpahGT/99JMeeugh1atXT40aNdKzzz6rEydO2Np+9NFHqlGjhlavXq2mTZuqYcOGOnDgQI738+GHH+qDDz7QSy+9pIcffjjb+tTUVG3btk133313luXt2rXTpUuXsjzO6125ckWTJk3SXXfdpZiYGI0YMULvvfdelqFKOQ1Hy2lIzS+//KLHHntMDRo0UK1atdS6dWvNnTtXFoslx/u+dkhUTs9jjx491LNnz2zt+vXrp0cffTTXx3S9G71Wrly5onr16mUbtpSenq7GjRtr4sSJtmWrV69Wp06dbMOT586dq4yMjFzv96efftIDDzygmJgYNWjQQE8++eRNj6ytWbNGV65cUatWrbIs37Rpk3r27Kk6deqoWbNmeumll3T+/Hnb+ptt+8ztPG/evCzPrT3vIwcPHtQTTzyhunXr6q677tKsWbM0atSoLPvElStXNH/+fLVv3161a9fWv//9by1ZsiTLc9+nTx8999xzGjZsmOrUqaNHH300x9fq8ePHNXz4cDVs2FDR0dF65JFHtGfPniw13ej9+XrOvE5v5sKFC5oyZYratm2r2rVr65577tGHH35oW9+tWzc9+eSTWdq0bdtWLVu2zLJs8ODBeuyxx25Y/5YtW9SnTx9FRUWpZcuWWr16tU6dOqWhQ4cqJiZGLVq0yDYi4mafSa1bt9axY8f08ccf27ZL5jCzuLg4devWTVFRUercubO+/PLLG26L33//XY899pgaNWqkunXratCgQdq/f7+k3D8rrpfXz4jcHDp0SEOHDlXDhg3VoEEDDRw4MMtr72bPW05OnTqlUaNGqUWLFoqKitJ9992njRs3Zqt/3rx56t69u6KiojRv3rxsw/ZeeOEF9evXT2vWrFG7du1Uq1YtdenSRT/88EOWvnbs2KHevXurTp06atmypd5++23169cvy/cUHx8ftWvXTosXL77pNoHrEJ6QLw4fPixJKlu2rCTp008/1ZAhQ1SxYkXNnz9fQ4cO1bp16zR48GDbORnff/+9hgwZopo1a2rBggWaO3euypYtq/HjxysuLi7bfaSnp+uZZ57R7t279cYbb6h69eo51nL58mU9/PDD+uuvvzRhwgS9+OKLWr16tXbs2JFr/bGxsRoxYoT+/e9/a+nSpRo1apS2bt1qGzI1duxY1ahRQzVq1NCqVatUs2ZNnTp1Sr1791ZycrJeeeUVLV26VJ06ddKKFSu0fPnyLP3PmjVLFotFM2fO1IgRI/Tdd99p8uTJtvXvvvuuRo4cqZo1a2revHkaMGCAVqxYYfsid/jwYfXs2VNnzpzRq6++qkmTJik+Pl69evXSmTNncn1c3377rYYNG6aqVatq/vz56tChg55//vlcby9JkydP1g8//KCRI0dq2bJlatOmjaZOnao1a9aoRIkSmjdvniTpySeftP1fkn799VedOHFC8+fP17PPPiuz2Zxj/1OmTFGdOnW0cOFC2wf722+/fcOarjV48GC1aNFC4eHhWrVqVbYvJdLVD9/u3bvr6NGjGjNmjKZPny7DMPTII49o+/btWW47duxYRUZGasGCBXrsscf04YcfauHChTned+YQR0m67777bEMz1q5dq/79+6tUqVKaOXOmRo0apR07dujBBx/M8vxkZGTojTfe0KRJkzRq1ChVqlQpx/tp3bq1vv322xyDhCTFx8crLS1N5cuXz7K8XLlykv7v9ZiT5557TqtXr9aAAQM0e/ZsJSUlafbs2bnePjf79u1Tv379FBISolmzZmnhwoWqX7++5s2bpy+++OKm7XN6Hu+77z7t2LFDR44csd3uxIkT2rZtm7p3725XXTd7rfj6+qpdu3b64osvspwf9tNPP+ns2bPq0qWLJGnx4sV68cUX1aRJEy1atEi9e/fW0qVL9eKLL+Z4v/Hx8Ro8eLBq1aqlhQsXatKkSTp8+LAGDBiQa5iUpHXr1qlly5by8fGxLfvuu+80cOBAhYaGavbs2Xruuee0YcMGPfPMM5Ls2/bX76f2bBtJSkxM1MMPP6wTJ05oypQpGjNmjL788kt99tlntvqsVqsGDRqk119/Xffff78WLVqk9u3ba/bs2VmGdUvSF198ocDAQC1cuFCPP/54tsefmJionj176o8//tCLL76oGTNmyGKxqHfv3rYvvzd7f76es6/T3KSkpOihhx7Sp59+qscff1wLFixQvXr1NHr0aC1atEiS1KJFC23fvt0Wso8ePar4+HidOHHCFlLT0tK0ZcuWHN+7rjV8+HC1bt1aixcvVoUKFTR27Fj17dtXVapU0YIFCxQVFaUpU6Zo165dkmTXZ9K8efMUHh6uFi1aZBtaNnDgQLVp00bz5s1ThQoV9J///EebNm3KsbatW7eqV69ekq5+ZkycOFEnTpxQz549dfDgwRt+VtxMXj8jTp48qQcffFB///23Xn75ZU2bNk0JCQl65JFHlJSUZNfzdr2EhATdd999+vXXX/XMM89o7ty5ioyM1JAhQ7Ru3bost120aJE6d+6sOXPmqF27djn2t3v3bi1btkzDhg3T/PnzZTab9dRTT9nOFT948KD69esnSZo5c6aeeuopLVmyJMcfwtq3b6/du3ff8H0ersWwPTjFarUqPT3d9ve5c+e0fft2LVy40Hb0yGq1avr06WrevLmmT59uu2358uXVr18/bdq0SS1bttSBAwfUrVs3jR492nabmJgYNWrUSNu2bVN0dLRtucVi0QsvvKBt27bpzTffVM2aNXOt8eOPP9bx48f12Wef2c7jiY6OzvYr/bViY2Pl5+enAQMG2L7EhISE6Pfff5fValXlypVt50dlDjPYuXOnqlevrtdee8227q677tJPP/2kbdu2ZTnv5M4778wy7nvXrl22X/UsFovmz5+vtm3bZvnVOzk5WZ9//rnS0tI0b948+fv766233rLdV5MmTdS2bVu9/vrruU5kMX/+fEVFRdmGdmWebzJjxoxct8X27dvVtGlTderUSZLUqFEjBQQEKDQ0VD4+PrbQescdd2QZ+pienq7x48ffdFjWAw88oBEjRkiSmjVrppMnT2rx4sV2n3h/xx13qHjx4lmG/1x/pGjevHny8fHR8uXLbdurZcuWuueeezR16tQsvzi2aNHCtv2aNGmin376Sd9//32OX8yuHeJYsmRJ1alTRxaLRdOnT1ezZs2ybNfMI3DLli2zPV7p6sQAN/vSdLPzQS5cuGCr51qBgYGSrg77y8lff/2lr7/+Wi+99JJ69+4t6epz0KlTpzxP+LJv3z7dddddmjZtmm24WdOmTfXtt99q27Zttv0nNzk9j/fcc49eeeUVffLJJxo2bJgk6ZNPPlFgYOANX7/Xsue10qVLF61Zs0axsbGqX7++pKuTC1SsWFG1a9fWhQsXtGDBAj344IMaM2aMbTuFhIRozJgxevTRR1WlSpUs97tr1y6lpKRo4MCBtvNrSpYsqY0bN+ry5cs5nl958eJF/f777+rQoUOW5ZlDMufNm2c7r8HHx0evvfaaEhIS7Nr21++n9m6bFStW6NKlS1q7dq3tcURHR2f5UvjDDz/o559/1syZM23Pc9OmTeXn56fXXnvN9gVfkry9vTVu3Djb++r1J9G//fbbSkpK0vvvv6/IyEhJV4egduzYUa+99prmzJlz0/fn68/9cMXrNCcfffSR/vrrL61cuVIxMTGSrr6npqena8GCBerZs6datmyphQsXateuXYqJidGWLVtUvnx5JSQk6JdfflHZsmUVGxury5cvZzvaeL0ePXrYjrgGBATogQceUFRUlJ5++mlJUrVq1fT111/rt99+U1RUlP7666+bfibVqFFDPj4+Kl68eLYhc3369NGQIUNsj6tbt26aP3++WrRoka22GTNmqFy5clqyZInth7JmzZrp7rvv1pw5c/Taa6/l+llxMzf6jMhpaOtbb72l1NRUvfnmmwoPD7dtm169eikuLk7Hjh276fMWEhKSpc8333xTiYmJ+uqrr2z7ZYsWLdSvXz9NnTpV99xzj62W+vXrZzky/vvvv2er8cKFC/roo49s7+0BAQF6+OGHtXXrVtuRpCJFiuj111+Xv7+/JKlixYo5/oBWu3ZtSdKWLVtUoUIFu7crHMeRJzjll19+Uc2aNW3/7rrrLg0fPly1atXSjBkzZBiGDh06pH/++UetW7dWenq67V+DBg0UFBSkn376SZL0+OOP65VXXtGlS5e0e/durV+/3nYo+trhapI0ffp0ffrpp+rbt6/tjSM3v/76q+64444sEyCUKlXqhmOrGzRooOTkZN1zzz2aMWOGfv31VzVr1kxDhw7N9aTMZs2a6Z133pGvr68OHDigjRs3auHChUpMTMxW//X3XbJkSSUnJ0u6+mvwmTNnsn05fOyxx/TRRx/J29tbW7duVcOGDeXn52fbnkFBQapfv75+/vnnHOtLSUnRH3/8ke0D+vovatdr1KiRPvjgAz3xxBN65513FB8fryFDhtz0C39ISIhd57N07Ngxy9933323zpw5o0OHDt20rb22b9+uVq1aZfnC6uXlpU6dOmn37t26dOmSbXlOz01uw/ZycvjwYZ0+fVr33HNPluV33HGHYmJish3pyu2IaV7c6EiGpFzPnfn1118lyTacRpLMZnO22u3RtWtXLV26VGlpadq3b5+++uorzZkzRxkZGUpLS8tzf5JUpEgR/fvf/87yy+7HH3+sjh072n0Olz2vlYYNG6p06dL6/PPPJV0dgrZhwwbbUacdO3YoJSUl23tY5tCjzPewa0VHR8vX11f33XefJk2apM2bN6tatWp65plncp2Y5sSJE8rIyFCZMmVsy1JSUrRnzx61bds2y3tPx44d9dVXXyksLMzhbW/Pttm6datiYmKyTLAQGRlp+9IpXX19eXl5qX379ln6v/fee23rM1WsWDHLUbXrbdmyRdWrV1dERIStJpPJpH/961+2mhx5f75eXl+nOdm+fXu2bZH5uK9cuaK4uDhFRUWpWLFiWbZno0aNFB0drV9++UXS1fBZpUqVLM97Tq69n9DQUEnK8qNisWLFJP3fjyl5+UzKSbdu3Wz/NwxDd999t+1HgWtdvnzZFvqvHWEQHBysVq1a2bUtbySvnxGxsbGqU6eOLThJV9/Hv/vuO9uRwJs9b9fbvn27YmJibMHp2janT5/OUos97+nFixfP8qNY5mdl5veArVu36l//+pctOEnK8f6lq++TwcHBLp/ND7njyBOcUrNmTY0bN07S1TdXX19flSpVKsuXg6SkJEnSuHHjbLe91qlTpyRdHa4xduxYbdiwQYZhqFy5crZfga+fbvnw4cNq0KCB3n77bT344IM3nDnp3Llztg+Va4WHh2c7zyhTTEyMlixZorfeektvvvmmlixZorCwMA0aNCjXIyKZw/DeffddXb58WaVKlVJUVJR8fX2z3fbaN0Tp6pfbzMeYub0yPxxzkpSUpPXr12v9+vXZ1uV2Yuq5c+dktVqzbYubzQA0evRolSxZUuvWrdOECRM0YcIExcTE6OWXX1a1atVybZd51ONmwsLCsvyd+bjPnTvnskkOzp07l+1+Mu/barVmOTJzo+fGHpnPX273d/25GwEBAXb3nZsiRYpIUpYQKP3fEafcvqxnHl26fp9xZCaylJQUTZgwQZ988onS09NVpkwZxcTEyMvLy6np0u+77z6tW7dOv/76q8xms/7+++88Tatsz2vFMAx17txZq1ev1pgxY/Tdd9/p8uXL6ty5s60PSbnOWpj5HnatMmXK6J133tGSJUv04Ycfavny5QoODtZDDz2k//znPzl+yc/80nvtPpH5ur3R+4Gj296ebZOYmJjjkf2wsDDb+2fme+z1Q3Mzv7xmPi7p5u8LSUlJOnLkSK6jCZKTkx16f87pfjIfR06P7frXaU7OnTuX5Qv6te0l6fz587bgt2XLFg0ZMkRbt27Vf//7X5UuXdo2Ff3mzZtvetRJyvl1fP371bXy8pmUk+s/G0JDQ2W1WrOcayddfX6tVmuu2/La598RN/qMyElSUtINg6g9z1tObTJPQ7hZG3ve069/3jLfDzJ/CEtMTMzxNZ/TNs7sL7cRBnA9whOcEhgYeNMjP8HBwZKkESNG5DgdatGiRSVdPffi0KFDeuuttxQTEyMfHx8lJyfrgw8+yNZmwoQJatiwoTp06KBx48ZpwYIFud5/sWLFspwzkSnzwzM3zZs3V/PmzZWcnKytW7dq+fLlmjhxoqKjoxUVFZXt9pkf5uPGjdO///1v2xfa++6774b3c73M7ZWYmJhl+dmzZ7Vnzx7FxMSoSJEiuuuuu3I8ad7LK+eXdUhIiEwmU7bAeLPt4OPjoyeffFJPPvmkjh8/ru+++04LFizQs88+a/ul3hnXfwBm1hcaGmr7QLn+pPy8HAmSru5jOQXl06dPS7q6j+T0BdgRmcM9cru/nIK8s+644w6ZzeZs+/n//vc/Scr1XKrMWhISElS6dGnb8rNnz2a77c2eg0mTJumrr77S7Nmzddddd9m+QDRp0iSPjyarhg0b6o477tCXX34pk8mkihUr5mlGLntfK126dNHixYu1bds2rV+/Xg0aNLD9ypv5mpw+fXq288qk3L/QZJ4wnpqaqtjYWK1atUqLFi1StWrVcjzim/l8XPtFLCgoSIZhZHs/uHLlirZu3aro6GjNmDHDoW1vz7YpWbJkjvvytecEFS1aVGfPnlVGRkaWAJX5msrLPl+kSBE1bNgw1yFzmUet8vr+fD1XvE6LFi2a42fLte8r0tUhwiNGjNCuXbuUkJBgO9I5a9Ys7dixQ3/99Zdefvnlm95fXjn7mZSUlJRl305ISJDZbM42pK1IkSIyDCPXbXn97fPqRp8ROSlSpEi214t09ahmmTJl7H7erlW0aFHbenvbOONGr7uKFStmW37+/Pl8+WxBzhi2h3xXsWJFhYaG6ujRo6pdu7btX0REhGbMmGH7hS82Nlb//ve/1ahRI9sHZObsM9cPSwoLC1N4eLiGDx+ujRs33vCE9MaNG+vo0aNZxh0nJiZq586dubZ59dVX1aNHD1mtVvn7+6tVq1a282AyZxC8fihUbGysKleurB49etg+pE6ePKm//vrrpsOqrlWxYkUVK1ZM3333XZbln3zyiQYMGKC0tDTbzGzVq1e3bc9atWrprbfe0jfffJNjv76+voqJidHXX3+d5dfob7/9NtdaUlJS1K5dO73xxhuSpNKlS6t3797q1KmTbTvkNhGEvb7//vssf3/++ecqVaqUypUrZ/ul9dpZtNLS0mwnRGe62ZTODRo00HfffZfll7mMjAx9/vnnql279g2HEeVVhQoVFB4enuWEeunqBAI7d+5U3bp1XXZfmXx9fVW/fn198803WZ7br776SkWKFMn1y2Tjxo1lGEa2WbSu3/eCgoL0zz//ZFl2/YnLsbGxatSokdq2bWv78r57924lJibavf/n9DwahqHu3btrw4YN+vbbb7MMJbKHva+VSpUqqWbNmvr888+1adMm25Az6erQKG9vb508eTLLe5iXl5dmzpyZ43CZt956S61atVJqaqp8fHzUpEkTTZgwQVLWWUivFRERIbPZnGVbBwYGqnr16tmekx9++EEDBgzQqVOn7N72129fe7ZNgwYNtHPnzixfHE+dOpXl/bNhw4ZKT0/Pth9lDresV69ejo83Jw0bNtThw4dVoUKFLNv6k08+0Ycffiiz2WzX+/PNuOJ12qBBAx07dizb5EPr1q2Tt7e37XXXrFkzWa1W20QP4eHhql27tgICAjRt2jQVK1Ys2xAyV7D3Mym3988NGzbY/m+1WvX111+rXr162d4vAwICVKtWLX3xxRdZfmS5cOGCvv/+e9vz7+hnxY0+I3JSv359xcXFZQlQZ86c0eOPP65NmzbZ/bxdq0GDBtqxY4eOHTuWrU14eHiutTiqQYMG2rx5s65cuWJbtmfPnhzfa86dO6fk5OQsP4Ahf3HkCfnObDbrmWee0UsvvSSz2axWrVrp/PnzWrBggU6ePGkbnhEVFaVPP/1UNWvWVMmSJfXbb79pyZIlMgzDNg74ej179tTatWttUy1nHsW6VpcuXbR8+XINHTrUdr7BwoULb/iFrnHjxnrzzTf1wgsv6N5771VaWppef/11hYSEqHHjxpKu/hq9Y8cObdmyRTVq1FBUVJQWLFigJUuWqE6dOjpy5IgWL16s1NTUXOvPbXs99dRTGj9+vEJDQ9W6dWsdPnxYc+bMUe/evVW0aFENHjxYPXv21MCBA9WrVy/5+vpq1apV2rBhg+bMmZNr38OHD9cjjzyioUOH6sEHH9Thw4dznV1Ikvz8/Gwz/nl7e6tq1ao6fPiwPv74Y9sJ45kfylu2bFGlSpWyjMG3x4oVKxQYGKgaNWro888/1+bNmzV16lQZhqGiRYsqJiZGK1asULly5VS0aFEtX75cKSkpWYZGBAcHKyEhQZs2bcpxvPnQoUP1ww8/qG/fvhowYIC8vb1t52+9/vrrear3Zkwmk4YPH65Ro0bp2Wef1b333quzZ89q3rx5Klq0aJ6m2M6LJ598Uo8++qiefvpp9ejRQzt27NCyZcv07LPP5jq0p3z58urZs6dmzZql9PR01ahRQ+vWrcv2paJVq1b69ttvNWXKFLVu3Vq//vqr1q5dm+U2UVFR+uKLL/T++++rUqVK2rdvnxYuXHjD1+/1rn8eM4cNde/e3TZ9ceZ5SPbKy2ulS5cuevXVV7Odv1OsWDE9/vjjeu2113Tx4kU1atRIJ0+e1GuvvSbDMHIcvtq4cWNNnz5dQ4YM0cMPPyyz2ayVK1fKx8cn1yFaAQEBqlu3rmJjY20zbUnSsGHD9OSTT2r48OHq2rWrEhISNHPmTLVt21Z33nmn3ds+ODhYv/32m3755RfVr1/frm3Tt29fvfvuu3rsscdskwcsWLBAaWlptiPD//rXv9SoUSONGTNGJ0+eVLVq1bR9+3YtXbpU3bp1s/uC29LVaeg/+eQT9evXT/3791exYsW0fv16ffDBBxo1apRt297s/flmXPE67d69u9577z0NGTJEw4YNU5kyZfTtt99qzZo1Gjp0qO2IZXBwsGJiYrRhwwY9+OCDkq4e2atfv75++OEHdenSJV+u6WXvZ1JwcLD27Nmj7du3ZwkOU6dO1ZUrV1ShQgXbBaxzm+Xu2Wef1WOPPaYBAwbooYceUlpampYsWaLU1FTbfuPoZ8WNPiNy0q9fP61du1aPP/64Bg4cKG9vby1cuFAlS5ZU586d5ePjY9fzdq1HH31U69atU79+/TR06FCFhIRo7dq12rp1qyZPnuzy52/QoEFav369Hn/8cfXv31/nz5/Xa6+9JpPJlO1xZ/6Q1axZM5fWgNwRnnBL3H///QoMDNTrr7+uVatW2b4kTJ8+3TaO+JVXXrGdUyNd/WI3btw42/kOOTGZTBo/frx69OihV199Nct035l8fHz09ttva/LkyZo0aZIMw9ADDzygsmXL5jodbYsWLTR9+nS98cYbtpOQ69Wrp+XLl9uGIPTu3Vu7d+/WE088oSlTpmjgwIE6e/asli9frvnz56tUqVLq0qWLDMPQ4sWLdf78+RzflHPSu3dvBQQEaNmyZbYrrj/xxBN64oknJF2dOejdd9/VrFmzNGLECFmtVt15552aP39+lpP/r1e/fn0tXbpUM2fO1NChQ1WmTBlNnjxZgwYNyrXN+PHjNXv2bL3xxhs6ffq0QkNDdd9999lmeAoKCtKjjz6qVatWadOmTTmePH8jEydO1Ouvv67Zs2erbNmyWWbskv5vvxgzZoyCgoJ03333qV69erbzBaSrX2I2bdpk+zC8/gTjKlWq6L333rNNR2wYhqKiorR8+XLbeXWu1L17dwUGBmrx4sUaMmSIgoKC1Lx5cw0fPjzHsfau0KRJE82dO1dz5szRkCFDFBERoREjRqh///43bPfSSy8pLCxMK1as0Pnz59WiRQs99NBDWa6v1qNHD/3vf//Txx9/rJUrV6pBgwaaM2eObWpi6eq1S9LS0jR79mylpqaqTJkyevLJJ3XgwAF9++23N7weUqbrn8fMc4wiIiJUrVo1hYWF5fl8rLy8VjJnX2zVqpXti16m//znPwoPD9d7772n119/XUWLFlWTJk00fPjwbLfNvN9FixZp/vz5Gj58uDIyMlSrVi298cYbOQ67ydSuXTvNnTtXV65csZ2b0qpVKy1atEjz5s3TkCFDVLx4cXXu3FlPPfWUJPu2vdls1qBBg7RgwQI98cQTWr9+vV3bJjg4WMuXL9ekSZM0YsQIBQYG6qGHHpK/v7/tB4zM97g5c+borbfeUmJiosqUKaPhw4fn+ceCiIgIrVy5UjNmzNDLL7+sK1euqHz58po0aZJtuJk978/2cPZ16u/vrxUrVmjGjBm2YF2xYsUstWZq0aKFfvnllywX6G3UqJF++OGHm06+4yh7P5P69++vyZMn67HHHstyDcSXX35ZixcvVnx8vGrUqKE33ngj1/fLJk2a6M0339ScOXM0fPhw+fj4qH79+nr11VdtMy06+llxs8+I65UqVUrvvfeepk2bphdeeEE+Pj5q1KiRZs2aZfuB1d7nLVN4eLjef/99zZgxQxMnTlRaWpqqVaumBQsW3PAz11HlypXTsmXLNHXqVA0bNkyhoaEaOHCgFi5cmO28wR9++EFRUVE5TiaB/GFYnTmTFwBQ6MydO1fz5s3Tn3/+6e5SJF0datSqVSvNmTNHbdu2dXc5+So5OVlt27bV888/r65du7q7HMXFxSkpKSnL9NTp6elq2bKlOnXqZDsahMLjo48+0qhRo7Rx48abzgCI/LFlyxZ5e3tnCavnz5/XXXfdpREjRqhv376Srp572rx5c7366quF/r3Rk3DkCQDgkfbu3auNGzfqq6++Uvny5W1Tgxdm/v7+euqpp7Rs2TJ17tzZ6XMKnXX8+HE988wzGjJkiBo2bKjk5GStWrVKFy5c0AMPPODW2oDC6o8//rAdxatZs6aSkpL05ptvqkiRIlmm11+5cqWqVKmSL0e/kDvCEwDAI125ckVvvvmmIiIiNHPmzHw5L8QT9ezZUxs3btTq1atzvCjmrdShQwclJSXpvffe07Jly+Tt7a3o6Gi98847uc7iCMA5/fv3V2pqqt5//32dOHFCAQEBatiwoaZMmZLlMgJvvfWWVqxYYff1zeAaDNsDAAAAADvcHj/jAQAAAICTCE8AAAAAYAfCEwAAAADYgfAEAAAAAHYgPAEAAACAHW77qcrPnLkg5ht0nmFIoaFF2J7wOOyb8GTsn/Bk7J/wZK7ePzP7u5nbPjxZreINwYXYnvBU7JvwZOyf8GTsn/Bkt3r/ZNgeAAAAANiB8AQAAAAAdiA8AQAAAIAdCE8AAAAAYAfCEwAAAADYgfAEAAAAAHYgPAEAAACAHQhPAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB0ITwAAAABgB8ITAAAAANiB8AQAAAAAdiA8AQAAAIAdCE8AAAAAYAfCEwAAAADYgfAEAAAAAHbwcncBAAAAAAoWk8mQyWQ41YfFYpXFYnVRRbcG4QkAAACA3UwmQ0WLBcrLyfCUbrHq3NlLBSpAEZ4AAAAA2M1kMuRlMjR4zxHtv5TiUB9VAv20oEY5mUwG4QkAAABA4bb/Uop+v5js7jJuKSaMAAAAAAA7EJ4AAAAAwA6EJwAAAACwA+EJAAAAAOxAeAIAAAAAOxCeAAAAAMAOhCcAAAAAsAPhCQAAAADsQHgCAAAAADsQngAAAADADoQnAAAAALAD4QkAAAAA7OC28PTRRx+patWq2f5Vq1ZNkrRnzx7df//9io6OVo8ePbR79+4s7T/77DO1bdtW0dHRGjJkiBITE93xMAAAAADcJtwWnjp27Kgff/zR9u/7779XuXLl1LdvX12+fFkDBgxQ/fr19dFHHykmJkYDBw7U5cuXJUm7du3S6NGjNXToUK1atUrnz5/XqFGj3PVQAAAAANwG3Bae/Pz8FB4ebvu3bt06Wa1WPffcc1q/fr18fX01YsQIVapUSaNHj1ZgYKC+/PJLSdI777yjDh06qGvXrqpWrZqmTp2qTZs2KT4+3l0PBwAAAEAh5xHnPCUlJWnp0qV69tln5ePjo7i4ONWrV0+GYUiSDMNQ3bp1tXPnTklSXFyc6tevb2tfqlQplS5dWnFxce4oHwAAAMBtwMvdBUjS+++/rxIlSqh9+/aSpNOnT6ty5cpZbhMaGqr9+/dLkk6dOqUSJUpkW//PP//k+b7/fz6DkzK3I9sTnoZ9E56M/ROejP0Tt4oj+5ir9097+3F7eLJarVq9erUef/xx27Lk5GT5+PhkuZ2Pj49SU1MlSSkpKTdcnxehoUUcqBq5YXvCU7FvwpOxf8KTsX8iPxUrFuhU+1u9f7o9PP3+++86efKkOnXqZFvm6+ubLQilpqbKz8/vhuv9/f3zfP9nzlyQ1epA4cjCMK7uvGxPeBr2TXgy9k94MvZP5MZsNjkdejKdPXtJGRmWPLdz9f6Z2d/NuD08bd68WfXr11fRokVtyyIiIpSQkJDldgkJCbahermtDw8Pz/P9W63iDcGF2J7wVOyb8GTsn/Bk7J/Ib87sX7d6/3T7hBG7du1S3bp1syyLjo7Wjh07ZP3/W8Jqteq3335TdHS0bX1sbKzt9idOnNCJEyds6wEAAADA1dwenvbv359tcoj27dvr/PnzmjRpkg4cOKBJkyYpOTlZHTp0kCT16tVLn3zyiVavXq19+/ZpxIgRatmypcqWLeuOhwAAAADgNuD28JSQkKDg4OAsy4KCgrR48WLFxsaqe/fuiouL05IlSxQQECBJiomJ0fjx4zV//nz16tVLRYsW1ZQpU9xRPgAAAIDbhNvPedq1a1eOy6OiovTxxx/n2q579+7q3r17fpUFAAAAAFm4/cgTAAAAABQEhCcAAAAAsAPhCQAAAADsQHgCAAAAADsQngAAAADADoQnAAAAALAD4QkAAAAA7EB4AgAAAAA7EJ4AAAAAwA6EJwAAAACwA+EJAAAAAOxAeAIAAAAAOxCeAAAAAMAOhCcAAAAAsAPhCQAAAADsQHgCAAAAADsQngAAAADADoQnAAAAALAD4QkAAAAA7EB4AgAAAAA7EJ4AAAAAwA6EJwAAAACwA+EJAAAAAOxAeAIAAAAAOxCeAAAAAMAOhCcAAAAAsAPhCQAAAADsQHgCAAAAADsQngAAAADADoQnAAAAALAD4QkAAAAA7EB4AgAAAAA7EJ4AAAAAwA6EJwAAAACwA+EJAAAAAOxAeAIAAAAAOxCeAAAAAMAOhCcAAAAAsAPhCQAAAADsQHgCAAAAADsQngAAAADADoQnAAAAALAD4QkAAAAA7EB4AgAAAAA7EJ4AAAAAwA6EJwAAAACwA+EJAAAAAOxAeAIAAAAAO7g1PKWmpmrcuHFq0KCB7rrrLs2cOVNWq1WStGfPHt1///2Kjo5Wjx49tHv37ixtP/vsM7Vt21bR0dEaMmSIEhMT3fEQAAAAANwm3BqeJk6cqJ9//lnLli3TjBkz9MEHH2jVqlW6fPmyBgwYoPr16+ujjz5STEyMBg4cqMuXL0uSdu3apdGjR2vo0KFatWqVzp8/r1GjRrnzoQAAAAAo5LzcdcdJSUlas2aN3nzzTUVFRUmS+vfvr7i4OHl5ecnX11cjRoyQYRgaPXq0fvjhB3355Zfq3r273nnnHXXo0EFdu3aVJE2dOlWtWrVSfHy8ypYt666HBAAAAKAQc9uRp9jYWAUFBalhw4a2ZQMGDNCUKVMUFxenevXqyTAMSZJhGKpbt6527twpSYqLi1P9+vVt7UqVKqXSpUsrLi7ulj4GAAAAALcPtx15io+PV2RkpNauXatFixYpLS1N3bt315NPPqnTp0+rcuXKWW4fGhqq/fv3S5JOnTqlEiVKZFv/zz//5LmO/5/P4KTM7cj2hKdh34QnY/+EJ2P/xK3iyD7m6v3T3n7cFp4uX76sI0eOaOXKlZoyZYpOnz6tl156Sf7+/kpOTpaPj0+W2/v4+Cg1NVWSlJKScsP1eREaWsTxB4Fs2J7wVOyb8GTsn/Bk7J/IT8WKBTrV/lbvn24LT15eXrp48aJmzJihyMhISdLx48f1/vvvq1y5ctmCUGpqqvz8/CRJvr6+Oa739/fPcx1nzlzQ/5/gD04wjKs7L9sTnoZ9E56M/ROejP0TuTGbTU6Hnkxnz15SRoYlz+1cvX9m9nczbgtP4eHh8vX1tQUnSapQoYJOnDihhg0bKiEhIcvtExISbEP1IiIiclwfHh6e5zqsVvGG4EJsT3gq9k14MvZPeDL2T+Q3Z/avW71/um3CiOjoaF25ckWHDx+2LTt06JAiIyMVHR2tHTt22K75ZLVa9dtvvyk6OtrWNjY21tbuxIkTOnHihG09AAAAALia28JTxYoV1bJlS40aNUr79u3T5s2btWTJEvXq1Uvt27fX+fPnNWnSJB04cECTJk1ScnKyOnToIEnq1auXPvnkE61evVr79u3TiBEj1LJlS6YpBwAAAJBv3HqR3OnTp+uOO+5Qr169NHLkSPXu3Vt9+vRRUFCQFi9erNjYWHXv3l1xcXFasmSJAgICJEkxMTEaP3685s+fr169eqlo0aKaMmWKOx8KAAAAgELOsFpv71GsCQmcBOkKhiGFhRVhe8LjsG/Ck7F/wpOxfyI3Xl5XJ4y4+5c/9fvFZIf6qB3kr28aVNXZs5eUnu7YhBGu3D8z+7sZtx55AgAAAICCgvAEAAAAAHYgPAEAAACAHQhPAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB0ITwAAAABgB8ITAAAAANiB8AQAAAAAdiA8AQAAAIAdCE8AAAAAYAfCEwAAAADYgfAEAAAAAHYgPAEAAACAHQhPAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB0ITwAAAABgB8ITAAAAANiB8AQAAAAAdiA8AQAAAIAdCE8AAAAAYAfCEwAAAADYgfAEAAAAAHYgPAEAAACAHQhPAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB0ITwAAAABgB8ITAAAAANiB8AQAAAAAdiA8AQAAAIAdCE8AAAAAYAfCEwAAAADYgfAEAAAAAHYgPAEAAACAHQhPAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB0ITwAAAABgB8ITAAAAANiB8AQAAAAAdiA8AQAAAIAdCE8AAAAAYAe3hqdvvvlGVatWzfJv2LBhkqQ9e/bo/vvvV3R0tHr06KHdu3dnafvZZ5+pbdu2io6O1pAhQ5SYmOiOhwAAAADgNuHW8HTgwAG1atVKP/74o+3fxIkTdfnyZQ0YMED169fXRx99pJiYGA0cOFCXL1+WJO3atUujR4/W0KFDtWrVKp0/f16jRo1y50MBAAAAUMi5NTwdPHhQd955p8LDw23/goODtX79evn6+mrEiBGqVKmSRo8ercDAQH355ZeSpHfeeUcdOnRQ165dVa1aNU2dOlWbNm1SfHy8Ox8OAAAAgELM7eGpfPny2ZbHxcWpXr16MgxDkmQYhurWraudO3fa1tevX992+1KlSql06dKKi4u7FWUDAAAAuA15ueuOrVarDh8+rB9//FGLFy9WRkaG2rdvr2HDhun06dOqXLlyltuHhoZq//79kqRTp06pRIkS2db/888/ea7j/+czOClzO7I94WnYN+HJ2D/hydg/cas4so+5ev+0tx+3hafjx48rOTlZPj4+mj17to4ePaqJEycqJSXFtvxaPj4+Sk1NlSSlpKTccH1ehIYWcfxBIBu2JzwV+yY8GfsnPBn7J/JTsWKBTrW/1fun28JTZGSktm3bpqJFi8owDFWvXl0Wi0XPP/+8GjZsmC0Ipaamys/PT5Lk6+ub43p/f/8813HmzAVZrY4/DlxlGFd3XrYnPA37JjwZ+yc8GfsncmM2m5wOPZnOnr2kjAxLntu5ev/M7O9m3BaeJCkkJCTL35UqVdKVK1cUHh6uhISELOsSEhJsQ/UiIiJyXB8eHp7nGqxW8YbgQmxPeCr2TXgy9k94MvZP5Ddn9q9bvX+6bcKIzZs3q1GjRkpOTrYt27t3r0JCQlSvXj3t2LFD1v+/JaxWq3777TdFR0dLkqKjoxUbG2trd+LECZ04ccK2HgAAAABczW3hKSYmRr6+vhozZowOHTqkTZs2aerUqXr88cfVvn17nT9/XpMmTdKBAwc0adIkJScnq0OHDpKkXr166ZNPPtHq1au1b98+jRgxQi1btlTZsmXd9XAAAAAAFHJuC09BQUFatmyZEhMT1aNHD40ePVoPPvigHn/8cQUFBWnx4sWKjY1V9+7dFRcXpyVLliggIEDS1eA1fvx4zZ8/X7169VLRokU1ZcoUdz0UAAAAALcBw2q9vUexJiRwEqQrGIYUFlaE7QmPw74JT8b+CU/G/onceHldnTDi7l/+1O8Xk2/eIAe1g/z1TYOqOnv2ktLTHZswwpX7Z2Z/N+PWi+QCAAAAQEFBeAIAAAAAOxCeAAAAAMAOhCcAAAAAsAPhCQAAAADsQHgCAAAAADsQngAAAADADoQnAAAAALBDnsNTenq63n//fR0/flyS9Nprr6lTp056/vnnlZSU5Or6AAAAAMAj5Dk8vfLKK1qwYIHOnz+vDRs2aOnSperSpYtOnDihCRMm5EeNAAAAAOB2XnltsH79ei1YsEDVqlXT0qVL1axZMw0YMECtWrVSz54986NGAAAAAHC7PB95Sk5OVmhoqNLT0/XDDz+oVatWkiSLxSIvrzxnMQAAAAAoEPKcdurWratp06YpKChIycnJatu2rfbt26cJEyaocePG+VEjAAAAALhdno88TZw4UWlpafrjjz80ZcoUhYaG6osvvlBoaKjGjh2bHzUCAAAAgNvl+chTqVKltHDhwizLnnnmGZcVBAAAAACeyKHrPMXGxmrYsGG2WfaWLFmizz//3NW1AQAAAIDHyHN4+vrrrzVgwABFRkbq8OHDSk9Pl5eXl1544QW99957+VEjAAAAALhdnsPTvHnz9PLLL2vkyJEym82SpP79+2vy5Ml68803XV4gAAAAAHiCPIenI0eOqE6dOtmWR0VF6eTJk66oCQAAAAA8Tp7DU+XKlbV58+Zsyz/++GNVrlzZJUUBAAAAgKfJ82x7o0aN0qBBg7R161alpaVp0aJFOnLkiHbv3p1tFj4AAAAAKCzyfOSpfv36+uKLL1SpUiW1bt1aSUlJqlOnjtavX68mTZrkR40AAAAA4HZ5PvIkSeHh4Xr66addXQsAAAAAeKw8h6fWrVvLMIxsyw3DkLe3t8LDw9WhQwf16tXLJQUCAAAAgCfIc3h6+OGHNW/ePD388MOqU6eOrFardu/erRUrVqhHjx4qUaKEFi5cqIsXL+qJJ57Ij5oBAAAA4JbLc3hau3atJkyYoE6dOtmWtWnTRlWrVtWiRYu0du1aVa9eXWPGjCE8AQAAACg08jxhxP/+9z9Vq1Yt2/IqVaro0KFDkqTy5cvrzJkzzlcHAAAAAB4iz+GpTp06mjt3ri5fvmxbdvnyZc2fP19RUVGSpE2bNqlcuXKuqxIAAAAA3CzPw/YmTJigQYMGqXnz5ipfvrysVquOHDmiUqVKac6cOfrxxx81efJkvfbaa/lRLwAAAAC4RZ7DU9myZbVu3Tpt2bJFf/31l8xms6pUqaImTZrIMAyFhIRo06ZNKl68eH7UCwAAAABu4dB1nsxms5o1a6ZmzZrZlqWmpmrv3r2Kjo52WXEAAAAA4CnyHJ5+++03jRs3TgcOHJDFYsmyzmw2a/fu3S4rDgAAAAA8RZ4njJg4caIiIyO1aNEi+fv7a+7cuRozZoxCQkI0derU/KgRAAAAANwuz0ee9u/fr2nTpqlSpUqqWbOmvL291bt3b4WGhmrp0qXq2LFjftQJAAAAAG6V5yNP/v7+MpvNkqSKFSvqzz//lCRFRUXp8OHDrq0OAAAAADxEnsNT48aNNWPGDJ08eVIxMTFav369kpKS9O233yo4ODg/agQAAAAAt8tzeBo9erTOnTunr7/+Wp06dVJQUJAaN26sKVOmaMiQIflRIwAAAAC4XZ7PeYqIiNDy5cttf69YsUIHDhxQcHCwIiIiXFocAAAAAHgKh67ztG/fPh06dEipqanZ1nXt2tXZmgAAAADA4+Q5PE2fPl2vv/66QkND5evrm2WdYRiEJwAAAACFUp7D06pVqzRp0iT16NEjP+oBAAAAAI+U5wkjihQpotq1a+dHLQAAAADgsfJ85GnkyJEaP368hg0bptKlS8tkypq/Spcu7bLiAAAAAMBT5Dk8paSk6I8//lDfvn1lGIZtudVqlWEY2rt3r0sLBAAAAABPkOfwNG3aND3wwAN64IEH5Ofnlx81AQAAAIDHyXN4Sk1N1cMPP6yyZcvmRz0AAAAA4JHyPGFE//79tXjxYl25ciU/6gEAAAAAj5TnI08//fSTdu7cqbVr1yosLExmsznL+o0bN7qsOAAAAADwFHkOT927d1f37t1dXsiAAQNUvHhxvfLKK5KkPXv2aOzYsfrrr79UuXJljRs3TrVq1bLd/rPPPtPs2bN1+vRpNWvWTBMmTFDx4sVdXhcAAAAASA6Ep27durm8iM8//1ybNm2y9X358mUNGDBAnTt31iuvvKL3339fAwcO1DfffKOAgADt2rVLo0eP1rhx41StWjVNmjRJo0aN0uLFi11eGwAAAABIdoanvn37at68eQoODlafPn2yTFF+veXLl+epgKSkJE2dOjXLhXfXr18vX19fjRgxQoZhaPTo0frhhx/05Zdfqnv37nrnnXfUoUMHde3aVZI0depUtWrVSvHx8UxkAQAAACBf2BWeGjZsKG9vb0lSo0aNXFrAq6++qi5duujUqVO2ZXFxcapXr54tpBmGobp162rnzp3q3r274uLi9MQTT9huX6pUKZUuXVpxcXGEJwAAAAD5wq7wNHTo0Bz/76wtW7bo119/1aeffqqXX37Ztvz06dOqXLlyltuGhoZq//79kqRTp06pRIkS2db/888/ea7hBgfRkAeZ25HtCU/DvglPxv4JT8b+iVvFkX3M1funvf3k+ZwnV7ly5YrGjh2rl156KdvFdpOTk+Xj45NlmY+Pj1JTUyVJKSkpN1yfF6GhRfLcBrlje8JTsW/Ck7F/wpOxfyI/FSsW6FT7W71/ui08zZs3T7Vq1VLz5s2zrfP19c0WhFJTU20hK7f1/v7+ea7jzJkLslrz3AzXMYyrOy/bE56GfROejP0Tnoz9E7kxm01Oh55MZ89eUkaGJc/tXL1/ZvZ3M3aFp0uXLikw0DUbKNPnn3+uhIQExcTESJItDH311Ve65557lJCQkOX2CQkJtqF6EREROa4PDw/Pcx1Wq3hDcCG2JzwV+yY8GfsnPBn7J/KbM/vXrd4/TfbcqFWrVjpx4oQkadSoUbp48aLTd7xixQp9+umnWrt2rdauXavWrVurdevWWrt2raKjo7Vjxw5Z//+WsFqt+u233xQdHS1Jio6OVmxsrK2vEydO6MSJE7b1AAAAAOBqdh15slgs+umnn9SkSROtXbtWDz/8sIoVK5bjbUuXLm3XHUdGRmb5O/PIVrly5RQaGqoZM2Zo0qRJ6tmzp1auXKnk5GR16NBBktSrVy/16dNHderUUe3atTVp0iS1bNmSmfYAAAAA5Bu7wtMjjzyiMWPG2KYOv++++yTJdmTIMAxZrVYZhqG9e/c6XVRQUJAWL16ssWPH6oMPPlDVqlW1ZMkSBQQESJJiYmI0fvx4zZkzR+fOnVPTpk01YcIEp+8XAAAAAHJjWK32jRI8f/68Lly4oDZt2mj16tUqXrx4jre7/oiSp0tI4CRIVzAMKSysCNsTHod9E56M/ROejP0TufHyujphxN2//KnfLyY71EftIH9906Cqzp69pPR0xyaMcOX+mdnfzdg9215wcLCCg4O1ceNGlS5dWikpKTpy5IgsFovuuOMOBQUFOVUwAAAAAHiyPE9VXqJECU2ZMkXvvfee0tPTr3bi5aXOnTtr3Lhx2a6/BADA9by9zfLxMTvVR2pqhtLSMlxUEQAAN5fn8PTqq69q06ZNWrhwoWJiYmSxWLRjxw5NnDhRs2bN0siRI/OjTgBAIWEyGSpa1N92Hq2j/P2tSky8JIuF8UQAgFsjz+Hps88+02uvvaZGjRrZlrVo0UK+vr567rnnCE8AgBsymQwZhqHkA2dlSUl3rA8/L/lXLiaTySA8AQBumTyHJ6vVqtDQ0GzLixcvrkuXLrmkKABA4Wa1WOVfOedLXuSlDwAAbqU8h6fGjRtr+vTpmj59um2SiPPnz2vmzJlZjkYBAJAbw2Ro88rlOn/qpEPtg0tEqHnPvi6uCgCAG8tzePrvf/+rvn37qnnz5qpQoYIk6fDhwypbtqwWLlzo8gIBAIXT3ztjderwQYfalqhQifAEALjl8hyeIiIi9Nlnn+mHH37QoUOH5OvrqwoVKqhp06YymUz5USMAAAAAuF2ew5MkeXt7q02bNmrTpo2r6wEAAAAAj8ShIgAAAACwA+EJAAAAAOyQ52F7//vf/3THHXfkRy0AgHzk42OWj49Do7VtUlPTlZqa4aKKAAAoWPL8KdqrVy8tXrxYtWrVyo96AAD5wGQyFBzsL8MwnOrHz89biYmXuDAtAOC2lOfwFBYWpjNnzuRHLQCAfGIyGTIMQ4fjflPKpYsO9eEXGKQK0XVlMhmEJwDAbSnP4alGjRoaPHiwateurcjISPn4+GRZP2XKFJcVBwBwrR/ff9upaytViK7r4ooAACg4HBr8fu+997q6DgDALVA8sqxb2gIAUBjkOTxxZAkACiaLJUOdnnrO6T4AALhdOXTkKTY2Vm+//baOHDmiRYsW6dNPP1VkZKQ6derk6voAAC5iMpm19ZODOp+Q7FD74DB/Ne5SycVVAQBQcOQ5PH399dcaNWqUHnjgAX3//fdKT0+Xl5eXXnjhBZ07d04PPfRQftQJALeUyWTIZHJuZjqLxepxEysc2X1GCfGOTRgRVjaI8AQAuK3lOTzNmzdPL7/8sjp37qyVK1dKkvr376/w8HDNmTOH8ASgwDOZDIUUC5TZyfCUYbEq6SzTegMAUFiY8trgyJEjqlOnTrblUVFROnnypCtqAgC3MpkMybncdJUhp49eAQAAz5HnI0+VK1fW5s2bsx1h+vjjj1W5cmWXFQYA7mQ2DCXv2ydLsmPnB5n8/eVfrZqLqwIAAO6U5/A0atQoDRo0SFu3blVaWpoWLVqkI0eOaPfu3Vq4cGF+1AgAbvHPf0crZc8eh9r61aihCh+tcXFFAADAnfI8bK9+/fr64osvVKlSJbVu3VpJSUmqU6eO1q9fryZNmuRHjQAAAADgdg5NVR4eHq6nn35aFy5ckLe3t/z8/FxdFwAAt4S3t1k+Pman+0lNzVBaGtfBAoDCLM/hKS0tTYsXL9bKlSt15swZSVLJkiXVr18/PfLIIy4vEACA/GIyGSpa1F+G4fzEHv7+ViUmMrsiABRmeQ5PEyZM0ObNm/Xcc8+pRo0aslgs2rVrl+bMmaMzZ85o+PDh+VEnAAAuZzIZMgxDl+JOyXIpzfF+Ar0VGF1CJpNBeAKAQizP4enzzz/X4sWLVb9+fduyatWqKTIyUsOHDyc8AQAKnIubjirt+CWH23uXDlRgdAkXVgQA8ER5njAiKChIXl7ZM1eRIkVyXA4AAAAAhYFdaef48eO2//ft21cjR47U6NGjVbt2bZnNZv31118aP368nnrqqXwrFABuNZ+KFd3SFgAAeCa7wlPr1q1tJ9NarVfHcg8YMCDbsnHjxqlnz575UScA3FLWjAxFTp/mdB8AAKDwsCs8bdy4Mb/rAACPYpjNmvPbHB27eMyh9pFBkRpWd5iLqwIAAO5kV3iKjIzMtiwhIUGpqakuLwgAPMWPx37U3sS9DrWtXrw64QkAgEImzzM8fPnllxo7dqzOnz+fZbnVapVhGNq717EvGgBuPZPJkMnk/PVtLBYr0zMDAIBCL8/hacqUKerYsaMefvhh+fn55UdNAG4Bk8lQ0WKB8nJBeEq3WHXuLBcHBQAAhVuew9Ply5fVt29fVahQIT/qAXCLXL04qGv6MgxxcVAAAFDo5fk6Tw899JDefPNNzncCCgGTLB7VDwAAgCfL85Gn9u3b65FHHtHatWsVFhZmm648EzPzAQWHYZiVvH+/rMnJjvfh7y//KlVcWBUAAIBnynN4ev7551WlShXdc889nPMEFHDWjAyXBB+uZwQAAG4HeQ5PR48e1cKFC1W2bNn8qAfALeTstYwkrmcEAABuH3kOT61atdLPP/+sBx98MD/qAXCLHb94XIfPHXa4vSEXzToBAADg4fIcniIjIzVp0iStXbtWZcuWldlszrJ+ypQpLisOQP7KsGTolX+94pJ+AAAACrs8h6czZ86oU6dO+VELgFvMbDJLG8dLSf9zvJOQO2Ru85LrigIAAPBQDl0kF0AhcmCDdCLO8faloiXCEwAAuA3kOTzNmzfvhuuHDh3qcDEAAAAA4KnyHJ62bduW5e+MjAwdPXpU58+fV7t27VxWGAAAAAB4kjyHpxUrVuS4fPLkydkumAsAQG6KRzp+yQtn2gIA4Kg8h6fc9OnTR927d9eoUaNc1SUAoJCyWDLU6annnO4DyAuTyZDJ5PwPvRaLVRaL1QUVAShoXBaeNm3aJF9f3zy1OXLkiMaPH6/ffvtNRYsW1cMPP6zHH39ckhQfH68XX3xRO3fuVOnSpfXf//5XzZo1s7X9+eefNXnyZMXHxys6OlqTJk3iwr0AUECYTGZt/eSgzickO9Q+OMxfjbtUclk9XiUC3Noe+c9kMlS0WKC8XBCe0i1WnTt7iQAF3IbyHJ5at26dbXjepUuXdO7cOY0cOdLufiwWiwYMGKDatWvr448/1pEjRzR8+HBFRETonnvu0ZAhQ3TnnXdqzZo12rBhg4YOHar169erdOnSOn78uIYMGaKnnnpKzZs31/z58zV48GCtW7eOoYMu4Mwvc2azyfZ/fpkDcCNHdp9RQvxFh9qGlQ1yWXiyWqwK7VnNJf3Ac5lMhlz1FcEwrvbHZxxw+8lzeHrqqaey/G0Yhry9vVWrVi2VK1fO7n4SEhJUvXp1vfzyywoKClL58uXVpEkTxcbGKiwsTPHx8Vq5cqUCAgJUqVIlbdmyRWvWrNFTTz2l1atXq1atWurfv7+kq9OnN23aVNu3b1ejRo3y+pBwDWd/mStWLND2f36ZA1AQGCZDm1cu1/lTJx3uI7hEhJr37OvCqpAfTLJIMruoHwC3ozyHp27durnkjkuUKKHZs2dLkqxWq3777Tf98ssvGjt2rOLi4lSjRg0FBPzfMIh69epp586dkqS4uDjVr1/fts7f3181a9bUzp07CU9OMpkMeZkMDd5zRPsvpTjcT5VAPy2oUY5f5gAUCH/vjNWpwwcdbl+iQiXCUwFgGGZt3LhRSUlJDvcREhKiNm3auK4oAAVKnsPT+fPn9cYbb+j3339Xenq6rNasX4yXL1+e5yJat26t48ePq1WrVmrXrp0mT56sEiVKZLlNaGio/vnnH0nS6dOnb7g+Lxjll7P9l1L0+0XHzkW4Htv49uGJz3VmTe6szRO3i7M86TF5Ui1S3uq5fv80m01Zhj87KiPDoowMjo7k5MCBAzpx4oTD7UuVKmULT56277maJ7x/4vbgyD7m6v3T3n7yHJ5GjBih33//XZ07d1ZQUFBem+dozpw5SkhI0Msvv6wpU6YoOTlZPj4+WW7j4+Oj1NRUSbrp+rwIDS3ieOG4qWuH8aFw8/Tn2l2vdU/fLo7wpMfkSbVIjteTuX9arRYZhvPhyVX9IHeetu/lJ74rIT85+1q61ftnnsPTzz//rHfeeUdRUVEuK6J27dqSpCtXrui5555Tjx49lJyc9ahHamqq/Pz8JEm+vr7ZglJqaqqCg4PzfN9nzlyQlVFlNmazScWKBapKoJ9T/WS2P3v2Er9+eqjM59pVPPW5Noyrb6x5ea27ctt4ynbxpMdUWGtxpJ5r90+z2aSiRZ177/0/ViUlXVJ6uvv3PU/h7ue6IHLk/RO3B094H3f1/pnZ383kOTxFRETIZHL+16yEhATt3LlTbdu2tS2rXLmy0tLSFB4erkOHDmW7feZQvYiICCUkJGRbX7169TzXYbWKN4TrZFitWlDD/sk/btSPxPa9nXjyc+3O17onbxdHedJj8qRaJMfqydw/rVbXzAhntRp8vt0Ct8v2ZV9CfnNm/7rV+6dDw/ZefvllDRs2TOXKlZO3t3eW9aVLl7arn6NHj2ro0KHatGmTIiIiJEm7d+9W8eLFVa9ePb3xxhtKSUmxHW2KjY1VvXr1JEnR0dGKjY219ZWcnKw9e/Zo6NCheX04yIHZMHTg4AylJMc73Ieff1lVrvSsC6sCgMLPZDIp6avDykh0fMIec3E/hbSr4MKqAACZHJ6qfMCAAVmuqWS1WmUYhvbu3WtXP7Vr11bNmjX13//+V6NGjdKxY8c0bdo0DRo0SA0bNlSpUqU0atQoDR48WN9995127dqlKVOmSJJ69OihZcuWacmSJWrVqpXmz5+vMmXKMNOeCyWe2aQLF/9wuH2RoJoS4QkA8uzKn2eVdvySw+29SwdKhCcAyBd5Dk8bN250yR2bzWYtWLBAEyZM0IMPPih/f3/16dNHffv2lWEYWrBggUaPHq3u3burXLlymj9/vu2oVpkyZTR37lxNnjxZ8+fPV0xMjObPn88Fcl0oINC5i0862x4AAADwNHkOT5GRkS6784iICM2bNy/HdeXKldM777yTa9sWLVqoRYsWLqsF/8dqzVCtmrNc0g8Az1K8lOMn+DrTFgCAwiDP4QmFHxcRBAoni8Wqu/vXdLoPT+Ln553teoN5wYgFAEBeEJ6QI1deRBCAZzCZDJ2aNVtpx4461N47soxKPPMf1xbloIBgH1mtVvn7+9z8xjdhtVoVULSYC6oCABR2hCcAuI1c2rxZKXv2ONTWr0YNyUPCk2+AlwzD0KW4U7JcSnO4H1OgtwKjS8gvkCGJAICbIzwhR2FhYW5tDwD2uLjpqNMz0wVGl3BhRQCAwozwhGwsFot69Ojhkn7MZucuqGyxWD3uHAsAAADcnghPyMYVF2n0Kh2ooH+VVXCwv1O1pFusOnf2EgEKAAAAbkd4Qo6cvUijvySzYWjwniPaf8mxEFYl0E8LapSTyWQQngAP4kkz3HmVCHBrewDA7YXwhBw5+4XCXNxPkrT/Uop+v5jsipIAuJk5LMy1M9wFO96PyWTIarEqtGc152uxWJVyyfEfiwAAtw/CE7Jx1RcSAIWLuXgxySrJFQeOrJJ/kLfDzS0WqwyToc0rl+v8qZMO91OqanXVbXePLp8763AfAIDbB+EJ2VhlkSGzS/qqEujnlrYA8kF6uksCS3CJCDXv2dclw3H/3hmrU4cPOtdJu3ucrgMAcHsgPCEHzs2QlynDatWCGuWc7gOAZ3E2sJSoUEnNe/Z1YUUAANwahCdkY1gtkguOPJkNQwcOzlBKcrxD7f38y6pypWedrgMAAABwBcITsjHMZq35a41OJ592uI8qIVXUplwbJZ7ZpAsX/3CojyJBNSXCEwAAADwE4QnZZFgy1ONO5y+SK0kBgZXc0hYAAABwNcITsjGbzHr758M6ed7xi+TWLB2sjrVLqlbNWU7VYrVmONUeAAAAcBXCE3L0wa9H9cfx8w63vze6lDpFRWrjxo1KSkpyqI+QkBC1adPG4RqAwsLX1ywfH8en9ZYkF16XFgCA2xbhCfkqKSlJCQkJ7i4DKLBMJkNFivjLcEH6sVqtMoeFOd1P8ciybm0PAIC7EJ6QL7xMJlksFvXo4dy5UxaLxUUVAQWTyWTIMAydOHRWacmOD2P19jerVMViMgcXcbwYLy9ZLBZ1euo5x/v4/ywWi65cTne6HwAAbiXCE/KNIed/KXdFH0Bh8MP7+5UQf9Hh9mFlg/Tg6IbOFZGeLpPJpK2fHNT5hGSHuylZqaiiWpbV5fOpztUDAMAtRnhCvrHKIsPJ60VZxZEnwNMc2X3GqSAnSVEtGboHACh4TO4uAIWXYfWMPgAAAABX4MgT8kW6xSLDbNancceUeMmxoTnFA33UOTrSxZUBAAAAjiE8IV8t2nTI4SnPa5YOJjwBAADAYzBsDwAAAADswJEnACgAipcKdGt7AABAeAIAj2exWHV3/5ou6Sfj/AUXVAQAwO2J8AQAHs5kMnRq1mylHTvqcB/+MTEq3ru3MhISXFgZAAC3F8ITABQAlzZvVsqePc510ru3a4oBAOA2xYQRAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB2YMAIAAACFjtlskpeXc8cJ0tMtysiwuKgiFAaEJwAAABQqJpOhYsX8ZRjOhSer1aLExMuyWKwuqgwFHeEJAAAAhYrJZMgwTDpwcIZSkuMd6sPPv6wqV3pWJpNBeIIN4QkAAACFUuKZTbpw8Q+H2hYJqilVetbFFaGgY8IIAAAAALAD4QkAAAAA7EB4AgAAAAA7EJ4AAAAAwA5MGAGgUPH1NcvHxzvb8qAgP7v7MAxXVgQAAAoLwhOAQsNkMlSkiL+MHNKPn1/2QHUjVqtVof6hrioNAAAUAgzbA1BomEyGrC66FIfVKhXzLeaazgAAQKHAkScAhYrJZOjUrNlKO3bU4T68I8uoxDP/UYY1w4WVAQCAgo7wBKDQubR5s1L27HG4vV+NGtIz/3FdQQAAoFAgPAEodHwqVnRrewAAUDgRngAUKtaMDEVOn+aSfs6nnndBRQAAoLAgPAEoVAyzWXN+m6NjF4853Eed8DrqVb2XziSfcWFlAACgoHNreDp58qQmTZqkrVu3ytfXVx07dtTw4cPl6+ur+Ph4vfjii9q5c6dKly6t//73v2rWrJmt7c8//6zJkycrPj5e0dHRmjRpksqWLevGRwPAU/x47EftTdzrVB+9qvdyUTUAAKCwcNtU5VarVcOGDVNycrLeffddzZo1S999951mz54tq9WqIUOGKCwsTGvWrFGXLl00dOhQHT9+XJJ0/PhxDRkyRN27d9eHH36o4sWLa/DgwbK6ao5iAAAAuIXZbJKvr5dT/7y9ze5+GCik3Hbk6dChQ9q5c6d++uknhYWFSZKGDRumV199Vf/6178UHx+vlStXKiAgQJUqVdKWLVu0Zs0aPfXUU1q9erVq1aql/v37S5KmTJmipk2bavv27WrUqJG7HhIAAACcYDIZKlbMX4bBpUjhmdwWnsLDw/X666/bglOmixcvKi4uTjVq1FBAQIBteb169bRz505JUlxcnOrXr29b5+/vr5o1a2rnzp2EJwAAgALKZDJkGCYdODhDKcnxDvdTtGg9lS3bx4WVAVe5LTwFBwerefPmtr8tFoveeecdNW7cWKdPn1aJEiWy3D40NFT//POPJN10fV4YhgPF45bjeSoYeJ5yxna5veTl+c68bX7tI+x7+auwb9/83j9vJPHMJl24+IdTfbgyPBX259rdHNm+rt4/7e3HY2bbmzZtmvbs2aMPP/xQb731lnx8fLKs9/HxUWpqqiQpOTn5huvzIjS0iONF45YoVizQ3SXADjxPOWO73F4cfb7z47OIfS9/3U7b93b/rnQ7Pdfu4Oz2vdX7p0eEp2nTpuntt9/WrFmzdOedd8rX11dJSUlZbpOamio/Pz9Jkq+vb7aglJqaquDg4Dzf95kzF8Q8E//HbDZ53JvE2bOXlJFhcXcZhY6rn2tPeJ4K4/7riY8Jucvr820YVz/4z5y5IJOp8L0mPUlhfM/Lb9fun7fqu5InvufdDs91XrnyeXJ0+7p6/8zs72bcHp4mTJig999/X9OmTVO7du0kSRERETpw4ECW2yUkJNiG6kVERCghISHb+urVq+f5/q1WEZ4KAJ6jgoHnKWdsl9uLI893fu0j7Hv563bZvnxX4vHnN2e2763eP906lcm8efO0cuVKzZw5U506dbItj46O1h9//KGUlBTbstjYWEVHR9vWx8bG2tYlJydrz549tvUAAAAA4GpuC08HDx7UggUL9MQTT6hevXo6ffq07V/Dhg1VqlQpjRo1Svv379eSJUu0a9cu3XfffZKkHj166LffftOSJUu0f/9+jRo1SmXKlGGmPQAAAAD5xm3haePGjcrIyNDChQvVrFmzLP/MZrMWLFig06dPq3v37lq3bp3mz5+v0qVLS5LKlCmjuXPnas2aNbrvvvuUlJSk+fPny2AqFAAAAAD5xG3nPA0YMEADBgzIdX25cuX0zjvv5Lq+RYsWatGiRX6UBgBAgefr6yVvb8dPBEhLsyg9PcOFFQFAwef2CSMAAIDrmIK8ZbVaFRDg61Q/VqtViYmXZLFwpjwAZCI8AQBQiJj8vWQYhrZv364LFy441EeRIkXUsGFDmUwG4QkArkF4AgCgEIqPj892WQ97hYWFqWHDhi6uCAAKPsITAACFicmQxWJRjx49nOrGYuGioABwPcITAACFicUqk8mkpK8OKyMx5ea3z4G5uJ9C2lVwcWEAUPARngAAcBE/v6uTNeRVQICPyy+3ceXPs0o7fsmhtt6lAyXCEwBkQ3gqREwmQyaTcx++ZrPbLv0FAAVWQNFislqt8vf3cay9kzPjAQBuDcJTIWEyGSpaLFBeToYnScqwWhUexAc5ANjLLzBQhmHoUtwpWS6lOdSHV5i//O8s7uLKAACuRHgqJEwmQ14mQ4P3HNH+S46NcZekKoF+WlCjnIL92TUAIK8ubjrq8FA5/+gwwhMAeDi+IRcy+y+l6PeLye4uAwAAACh0OMEFAAAAAOxAeAIAAAAAOxCeAAAAAMAOnPMEFFDOTk3v6mnpAwN95MDlbbJITU1TSkqGawoCAABwMcITUAC5amr6DKtV5qASzhUTVEJWq1U+Pt7O9SPJx8dLqamXZLE4mcIAAADyAeEJKIAyp6af+fc/OpqS6lAfZfx8NLx8SckvxLli/EJkGIaS9+2TJdnxmR5N/v7yr1ZNJpNBeAIAAB6J8AQUUBlW69Xw42QfZpNr3gYSX1+m1EOHHG7vU7GiIqdPc0ktAAAA+YHwBBRQJlkkmZ3v43KCk4V4yZqR4ZLgY83gfCcAAOC5CE+FTJVAP7e2x61jGGYp4YCUnuJYB75FZRQrK1085VwhlnQZZrPm/DZHxy4ec7ibyKBIDas7zLlaAAAA8hHhqRDJsFq1oEY5l/TjZWIW+wJhTX/pRJxjbWvfL/V43WWl/HjsR+1N3Otw++rFqxOeAACARyM8FSJmw9CBgzOUkhzvcB9+/mVVudKzSrdYXFgZAAAAUPARngqZxDObdOHiHw63LxJUU6r0rAsrAgAAAAoHxmYBAAAAgB0ITwAAAABgB4btAXCJikUrurU9gMLNy8ssb2/Hf/M1mQwXVgPgdkV4AuAck5cyLBl65V+vON1VhoXrPAHIzmQyFBLiL8MgAAFwL8ITAOdY0mU2maWN46Wk/zneT8gdMrd5yXV1ASg0TCZDhmFo+/btunDhgkN9REREqFatWi6uDMDthvBUyAQEVnJre9zGDmxw/JpTklQqWiI8AbiBHTt26MSJEw61rV27NuEJgNMIT4WI1ZqhWjVnuaQfLpILAAAAZEV4KkQMw6yNGzcqKSnJ4T5CQkLUpk0bLpILAAAAXIfwVMgcOHDA4SENklSqVCm1adPGhRUBt7eAAB9ZrY635/x4wDN5e5tlNjs3SiM93aKMDH6sBAoSwhMA5INQ/1BZrVb5+no73ZfVapU5LMwFVQFwVlBQkKxWi4KC/Jzuy2q1KDHxsiwWJ35hAXBLEZ4AIB8E+wTLMAydOHRWacmOT8Hu7W9WqYrFZA4u4sLqADgqICBAkqvCjlUmk0F4AgoQwhMA5KMf3t+vhPiLDrcPKxukB0c3dGFFAJxhsVhceo4xgIKF8AQA+ah4qUC3tgeQPzjHGLg9EZ4AIB+YDbMsFqvu7l/T6b4sFqsyzjt2YVAAAOA6hKdCJszJk8qdbQ/gqgxrhkwmQ6dmzVbasaMO9+MfE6PivXsrIyHBhdUBAABHEJ4KEYvFoh49erikHy6SC7jGpc2blbJnj3Od9O7tmmIAAIBTCE+FiMlkUtJXh5WRmOJwH+bifgppV4GL5AIAAADXITx5CJPJkMnk+NUwMy/Ud+XPs0o7fsnhfrxLB0rtKjjcHgAAACisCE8ewGQyVLRYoLycCE+SlGG1yhTk/AU5AQAAgJupEuj4xaKdaetOhCcPYDIZ8jIZGrzniPZfcmzIXZVAPy2oUU7ekUGyXExzuBavEgEOtwUAAMDtIcNq1YIa5Zzuo6AhPHmQ/ZdS9PvFZKf6CGlXwelhd1aLlQkjAAAAkCuzYejAwRlKSY53qL2ff1lVrvSsi6vKf4SnQmb/9p+VEH/E4faBIcUV1aYdE0YAAADghhLPbNKFi3841LZIUE2J8AR3q9LwLlVpeJdTfTBVOQAAAJAd4amQ2fV9vP45eM7h9sFh/mrcpRJHngAAAIDrEJ4KmcvnUnX2n8sOtzcM52b8AwAAAAorwlMh07hLJTXuUsmpPixMGAEAAABk4xHfkFNTU3XPPfdo27ZttmXx8fHq16+f6tSpo44dO+rHH3/M0ubnn3/WPffco+joaPXt21fx8Y7N9IHsrLIo8VKqu8sAAAAAPIrbjzxduXJFzz77rPbv329bZrVaNWTIEN15551as2aNNmzYoKFDh2r9+vUqXbq0jh8/riFDhuipp55S8+bNNX/+fA0ePFjr1q277YedfRp3THtPnHe4ffVSweocHanTF6+4sCoAAACg4HNreDpw4ICeffZZWa+7QNbWrVsVHx+vlStXKiAgQJUqVdKWLVu0Zs0aPfXUU1q9erVq1aql/v37S5KmTJmipk2bavv27WrUqJE7HorH+GbPSa2LO+Fw+3ujS6lzdKQLKwIAAAAKB7cO28sMO6tWrcqyPC4uTjVq1FBAQIBtWb169bRz507b+vr169vW+fv7q2bNmrb1AAAAAOBqbj3y9NBDD+W4/PTp0ypRokSWZaGhofrnn3/sWp8Xt/kovwKD5+n2wvMNeA5ej/nPU7dxZl2eWt+txDbIX45sX1fvn/b24/ZznnKSnJwsHx+fLMt8fHyUmppq1/q8CA0t4nihuCWKFQt0dwm4hXi+Ac/B6zH/FYRtfLt/VyoIz1FB5uz2vdX7p0eGJ19fXyUlJWVZlpqaKj8/P9v664NSamqqgoOD83xfZ85c0HWnXN1yZrOJF+YNnD17SRkZXLT3WoV5n3Hm+S7M2wVwB095/y3Mr21P2cY5MYyrX0xv5XclT3yuPfk5chdXPk+Obl9X75+Z/d2MR4aniIgIHThwIMuyhIQE21C9iIgIJSQkZFtfvXr1PN+X1Sq3hyfcHM/R7YXnG/AcvB7zn6dvY74r8fjzmzPb91bvnx5xnafrRUdH648//lBKSoptWWxsrKKjo23rY2NjbeuSk5O1Z88e23oAAAAAcDWPPPLUsGFDlSpVSqNGjdLgwYP13XffadeuXZoyZYokqUePHlq2bJmWLFmiVq1aaf78+SpTpkyBn6a8SqCfW9oCAAAAuDmPDE9ms1kLFizQ6NGj1b17d5UrV07z589X6dKlJUllypTR3LlzNXnyZM2fP18xMTGaP39+gb5AbobVqgU1yjndh5fJIw8mAgAAAAWex4SnP//8M8vf5cqV0zvvvJPr7Vu0aKEWLVrkd1m3jNkwdODgDKUkxzvU3s+/rCpXelbpFk5oBAAAAPKDx4QnSIlnNunCxT8calskqKZU6VkXVwQAAAAgE2O8AAAAAMAOhCcAAAAAsAPhCQAAAADsQHgCAAAAADsQngAAAADADoQnAAAAALAD4QkAAAAA7MB1njxIQGAlt7QFAAAAcHOEJw9htWaoVs1ZTvfhZeJgIgAAAJAfCE8ewjDM2rhxo5KSkhxqHxISojZt2ijdYnFtYQAAAAAkEZ48yoEDB3TixAmH2pYqVUpt2rRxcUXArRcY6COr1bG2huHaWgAAAK5FeALgGYJKyGq1ysfH26lurFarQv1DXVQUAADA/+EEGQCeISBMcsWwU4auAgCAfMKRJwCewZIuw2zWnN/m6NjFYw51USe8jnpV76UzyWdcXBwAAADhCYCH+fHYj9qbuNfh9r2q93JhNQAAe5nNJnl5OTeoyWTi5FV4NsITPJ7Z7PzoUovFKovFwVkIAADADZlMhooV85dhcEYICjfCkwcJCwtzS1tPFR7kqwyrVcHB/k73lW6x6tzZSwQoAADygclkyDBMOnBwhlKS4x3up2jReipbto8LKwNci/DkISwWi3r06OF0H4XpIrnB/l4yG4YG7zmi/ZdSHO6nSqCfFtQoJ5PJIDwBhUzxyLJOtQ8uEeGiSgBIUuKZTbpw8Q+n+iA8wZMRnjyEyWRS0leHlZHoWEgwF/dTSLsKhfIiufsvpej3i8nuLsNlGBMOOO/qjyEWdXrqOXeXAgC4jRCePMiVP88q7fglh9p6lw6U2lVwcUVwNVeOCbdaLTKCSrigKqDgsVisMplM2vrJQZ1PcPzHlZKViiqqpXNHr67lVSLA4bbm4n4uqwMAkD8IT/B4VQKd+0LhbHtXunrEyFVDB61Xr40E3MaO7D6jhPiLTvXhivBkmL1ktVgV2rOa030BADwX4Qkey8tkUobVqgU1yjndV4bVg851skgyu6gfeDyfihWdau8dWcZFlSA/WTPSZZgMbV65XOdPnXSoj1JVq6tuu3tcXBkAwJUIT/BY6RaLzIbh9Mw9fv5lVbnSsy6szDkWw0XZyZDMlxNc0BPyg9kwy5qRocjp09xdCm6hv3fG6tThg453QHgCAI9GeILHc3bmniJBNSUPCk9mk1naOF5K+p/jnZRtKHPDAdLFU64rDC6VYc2QYTZrzm9zdOziMYf7qRNehwv/AgDgIQhPHsSZE42daevpAgIrubV9vjiwQToR51wfDQe4phbkqx+P/ai9iXud6oPwBACAZyA8eQhXnGhstVgL1XWevEwmWa0ZqlVzltN9Wa0ZLqgIAAAAtzPCk4dw9kTj4BIRat6zb6G6zlO6xSLDMGvjxo1KSkpyuJ+QkBC1adPGdYUBAOzm5WWWtzfXtgNQOBCePIgzJxqXqFBJzXv2dXFFniEpKUkJCUyMAAAFjclkKCTEX4ZB+AFQOBCePEjxSMevNeJMW0/lZTLJYrGoR48eTvdlKURH5ACgoDCZDBmGoe3bt+vChQsO9xMREaFatWq5sDIAcAzhyUNYLBnq9NRzTvZhKVTnPAEACocdO3boxIkTDrevXbs24QmARyA8eQiTyaytnxzU+YRkh9oHh/mrcZdKhe6cJ5PJpC8Pf6nElESH+ynuV1ztK7R3YWXA7a14qUCn2geH+buoEgCAOzkzo7FHzoZsB8KTB7lwJkVn/7nsUNvCPJ78+/jvdejcIYfbVyxakfAEuIKXlywWq+7uX9PdlQAA3MwVMyIXxNmQCU8ewhVfSCyFcKryDEuGXvnXK073lWEpeC9OwOOkp8tkMnRq1mylHTvqcDf+MTEq3ru3CwsDgPzj7W2W2ezc96v0dIsyMgrP6CBJTs+IXFBnQyY8eQhnv5B4R5ZRiWf+U6iG7UmSq6Jg4YmUgPtd2rxZKXv2ONcJ4QmAh/PxCZPValFQkJ/TfVmtFiUmXpbFYnVBZZ7jwIEDDp/PWKpUKcITnOPMFxK/GjWkZ/7j2oLcLN1ikWEya9pX+xSf6NhwRkkqWzxAz7dz7gLEAADg9uLlFSzDMOnAwRlKSY53uB8//7KqXOlZmUxGoQtPtyPCkwfxqVjRLW093fd/ntYfx8873L5m6WDCEwAAcEjimU26cPEPh9sXCaopVXrWhRV5jrCwMLe0dSfCk4ewZmQocvo0p/soTOc8AQAAwDO54lqcBfE6nIQnD2ExJLML+ki8lOqSegojZ0/2tFisHG4HAACQZDKZdHH7CWVccOy7p7mIj4IalnJxVfmP8OQhzCaz3v75sE6eT3GoffVSweocHanTF6+4uLKCLzzIVxlWq4KDnbu2TLrFqnNnLxGgAADAbc9qsTodfqwF8DsV4cmDfPDrUYfP7bk3upQ6R0e6uKLCIdjfS2bD0My//9HRFMd+HSnj56Ph5UtysicAAIAkw2Ro18avdCkp0aH2gSHFFdWmnYuryn+EJxR6XiaTMqxWDS9f0ql+MqyEJgAAAEmyWDKcDj+WAngdTsITCr10i0Vmw9CxY6uUmnrKoT58fEooMvJBF1cGAJ7N19dL3t6O/3BkMhkurAaAJzGZzNr6yUGdT0h2qH1wmL8ad6nk4qryH+EJhZ6XySSrNcPp8GO1ZigoyEfOHIAy+B4BoAAwBXnLarUqIMDX3aUA8GAXzqTo7D+OXYvTKKBfighPuD1YTZLTr1GTvL2dnRNRslqtMoJKON0PAOQXk7+XDMPQ9u3bdeHCBYf7iYiIUK1atVxYGfKL2WySl1fOs9L6+Nz86yJHGW8/FotVd/ev6XQfBQ3hCYVeusUiw2To4patslxwbEIOU5FgBTVp7JJ6LFaL09PSA8CtsGPHDp04ccLh9rVr1yY8FQAmk6FixfxlGDmHJ2dnq0XhZDIZOvvBaqWfPu1Qe6/wcBV74H4XV5X/CE8o9LxMJlkzMpwOP9aMDBlb5kondzveSdmGMjccIF107NwrAABczWQyZBgmHTg4QynJ8Q71UbRoPZUt28fFlcGTWTMynA4/1gwmjAA8TrrFIsNslrYvkS7841gnETVl1Ooh/fGRdCLOuYIaDnCuPQAA+SDxzCZduPiHw+0JT7cXw2zWmr/W6HSyY0eewv3D1ePOHi6uKv8V6PB05coVjRs3Tl9//bX8/PzUv39/9e/f391lwVPFb5cS/nK8fa2C9wIHABRuNzpXyV6crwRHZFgynA4/GUxVfmtNnTpVu3fv1ttvv63jx49r5MiRKl26tNq3b+/u0uBBvEwmWS0ZMnq87u5SAKBACQsLc6p9SEiIawpBjm52rhKQn8wmsxT7tnTRwVE9QSVlrveIa4u6BQpseLp8+bJWr16tpUuXqmbNmqpZs6b279+vd999l/CELK5OGGHWtK/2KT7Rsek065UrpkfuquDiygDAQ5kMWSwW9ejBEXdP5opzlSTOV4JjrJYMGU6GHytHnm6dffv2KT09XTExMbZl9erV06JFi2SxWGQy8SsMsvr+z9P647hjs+1Jcm14CrvTufYhd7imDtxQxaIVHW4bGRTpdB/X9uMqPhUdr8c7sowLK5GKlwp0uG1wmGfO/lU8sqzDbYNLREiSvEoEOFWDubif0/2YQ3xlMpl0cfsJZVxIdbgf74gABdQKd7j99Zw5EpZ5FMxVR9Oc7Sezvb+/j6wOXkAw8zo5zp6rJF09Xykg0PELlvr5O77v56fC+Jg8hbM/TJctHqDn21VzcVX5r8CGp9OnT6tYsWLy8fGxLQsLC9OVK1eUlJSk4sWL29WPySSnLnrqSjVLB8vfx7FJrCuFBzndh6f144m1qGSU5O3EF5sipSRLhuSqIYTO1BNa2fk+8qGf6sWry9/LsS/GFYIrON2HJJUIKKEMS4Ze+dcrDvchySV9ZPKtXl2Gv+OPyatEhKwZGYqcPs3ttfhUruyS64NIUvgdReTl4/iPZcVKXt1nS5SvJG9fxy8IG1g8VBaLRZ2ees7hPiTJarEqtKfzXyZc0Y/VYlVQw1JO1yJJJUuWlLe3t8PtixQp4pIjYa46muaqfqxWq/z8HN8umYKK1JDJ5Phr0tcvQlZrhmrVnOX2WvwDKrmkH096TAGBVz+bgoJ8bnLLG0tLS3PqdeTqfiTp9IUrOn4uxaG2ftdcO9ORYx6Z19h11Xd5e6/Za1gd/bnDzdauXavXXntN3333nW1ZfHy82rZtq02bNqlkyZJurA4AAABAYVNgx7b5+voqNTXrUILMv/38/NxREgAAAIBCrMCGp4iICJ09e1bp6em2ZadPn5afn5+Cg4PdWBkAAACAwqjAhqfq1avLy8tLO3futC2LjY1V7dq1mSwCAAAAgMsV2JTh7++vrl276uWXX9auXbu0YcMGvfHGG+rbt6+7SwMAAABQCBXYCSMkKTk5WS+//LK+/vprBQUF6bHHHlO/fv3cXRYAAACAQqhAhycAAAAAuFUK7LA9AAAAALiVCE8AAAAAYAfCEwAAAADYgfAEpxw5ckSPPfaYYmJi1LJlS73++uvuLgnI0YABA/TCCy+4uwzA5ptvvlHVqlWz/Bs2bJi7ywIkSampqRo3bpwaNGigu+66SzNnzhSnycMTfPTRR9neO6tWrapq1ardkvv3uiX3gkLJYrFowIABql27tj7++GMdOXJEw4cPV0REhDp37uzu8gCbzz//XJs2bVK3bt3cXQpgc+DAAbVq1UoTJkywLfP19XVjRcD/mThxorZt26Zly5bp0qVLeuaZZ1S6dGn17NnT3aXhNtexY0c1b97c9nd6eroeeeQRtWzZ8pbcP+EJDktISFD16tX18ssvKygoSOXLl1eTJk0UGxtLeILHSEpK0tSpU1W7dm13lwJkcfDgQd15550KDw93dylAFklJSVqzZo3efPNNRUVFSZL69++vuLg4whPczs/PT35+fra/Fy9eLKvVqueee+6W3D/hCQ4rUaKEZs+eLUmyWq367bff9Msvv2js2LHuLQy4xquvvqouXbro1KlT7i4FyOLgwYO666673F0GkE1sbKyCgoLUsGFD27IBAwa4sSIgZ0lJSVq6dKkmTpwoHx+fW3KfnPMEl2jdurUeeughxcTEqF27du4uB5AkbdmyRb/++qsGDx7s7lKALKxWqw4fPqwff/xR7dq1U9u2bTV9+nSlpqa6uzRA8fHxioyM1Nq1a9W+fXu1adNG8+fPl8VicXdpQBbvv/++SpQoofbt29+y+yQ8wSXmzJmjRYsWae/evZoyZYq7ywF05coVjR07Vi+99FKWw/uAJzh+/LiSk5Pl4+Oj2bNna+TIkfr00081depUd5cG6PLlyzpy5IhWrlypKVOmaOTIkVqxYoXeeustd5cG2FitVq1evVoPP/zwLb1fhu3BJTLPJ7ly5Yqee+45jRgx4pYdPgVyMm/ePNWqVSvLSaWAp4iMjNS2bdtUtGhRGYah6tWry2Kx6Pnnn9eoUaNkNpvdXSJuY15eXrp48aJmzJihyMhISVcD//vvv6/+/fu7uTrgqt9//10nT55Up06dbun9Ep7gsISEBO3cuVNt27a1LatcubLS0tJ08eJFFS9e3I3V4Xb3+eefKyEhQTExMZJkGw711VdfaceOHe4sDZAkhYSEZPm7UqVKunLlis6dO8f7J9wqPDxcvr6+tuAkSRUqVNCJEyfcWBWQ1ebNm1W/fn0VLVr0lt4vw/bgsKNHj2ro0KE6efKkbdnu3btVvHhxPvjhditWrNCnn36qtWvXau3atWrdurVat26ttWvXurs0QJs3b1ajRo2UnJxsW7Z3716FhITw/gm3i46O1pUrV3T48GHbskOHDmUJU4C77dq1S3Xr1r3l90t4gsNq166tmjVr6r///a8OHDigTZs2adq0aRo0aJC7SwMUGRmpcuXK2f4FBgYqMDBQ5cqVc3dpgGJiYuTr66sxY8bo0KFD2rRpk6ZOnarHH3/c3aUBqlixolq2bKlRo0Zp37592rx5s5YsWaJevXq5uzTAZv/+/apcufItv1+G7cFhZrNZCxYs0IQJE/Tggw/K399fffr0Ud++fd1dGgB4tKCgIC1btkyTJ09Wjx49FBgYqJ49exKe4DGmT5+uCRMmqFevXvL391fv3r3Vp08fd5cF2CQkJCg4OPiW369htVqtt/xeAQAAAKCAYdgeAAAAANiB8AQAAAAAdiA8AQAAAIAdCE8AAAAAYAfCEwAAAADYgfAEAAAAAHYgPAEAAACAHQhPAAAAAGAHwhMAIJuqVatq27ZtDrV94YUX9MILL0iS5s6dqz59+ty0jdVq1bvvvuvQ/d0KM2fO1OrVqyU5t22cde19nzlzRl988UWO65zVp08fzZ0796a3mzVrlj744AOX3CcAFASEJwBAvunfv79dX8J/+eUXjR8//hZUlHeHDh3SN998o27durm7FP3444+KiYmRJE2fPl2bNm3Kl/uZO3eu+vfvf9PbPfbYY1q8eLHOnj2bL3UAgKchPAEA8k1gYKBCQkJuejur1Zr/xTho6dKl6tatm7y8vNxdisLDw+Xj4yMpf7dZSEiIAgMDb3q74OBgNWvWTO+9916+1QIAnoTwBAAF0NGjR1W1alV9+umnat68uerXr6+JEycqPT3ddptvvvlGHTt2VHR0tO677z5t377dtu7ixYsaNWqUmjRpolq1aql9+/basGFDjve1YcMGRUVFafPmzTmu//XXX9W1a1dFRUXp6aefVnJysm3dtcP20tLSNGbMGDVq1EgxMTEaNGiQTp48qaNHj6pv376S/m/oWWpqqqZMmaLmzZurZs2aat26tVatWmXrt3Xr1nr33Xf1wAMPqHbt2urSpYt2795tW3/kyBE99thjiomJUcuWLbV8+XLbur/++kt9+vRRVFSU2rVrd8PhgufPn9f69evVpk2bHNdfuXJF06ZNU4sWLVSnTh0NGjRIJ06csK2Pj49Xv379FB0drc6dO2vZsmVq3bq1bf3q1avVvn171apVS40aNdK4ceOUkZEh6f+GP957771q0qSJ/v77b9v2mTt3rj7++GN9/PHHWfr79ddf1blzZ9WuXVsPP/ywjh07Jknatm2bWrdurQ8//FBNmzZVgwYNtHTpUv3yyy9q3769YmJiNGLECFksFknZh+29+eabat26tWJiYvTYY48pPj4+y3OxatUqW1sAKMwITwBQgM2bN0+zZs3SvHnz9PXXX9u+8O7bt08jR47Uk08+qXXr1unee+/VE088oSNHjkiSJk2apMOHD+uNN97QZ599pvr162v06NFKTU3N0v9vv/2m559/Xq+88oqaN2+e7f4TExM1cOBA3XXXXVq7dq0qV66sL7/8Msda3333Xf3yyy9644039OGHH+rSpUuaPHmySpUqZas7c1jakiVL9P3332vu3Ln68ssv1bVrV02YMEEJCQm2/ubOnasBAwZo3bp1KlKkiCZOnCjpaqDp37+/AgMD9cEHH+ill17SrFmz9N133yklJUVPPPGE6tWrp3Xr1mnkyJFasGCB1q5dm2PN27dvV0hIiCpVqpTj+rFjx+qbb77Rq6++qpUrVyo9PV2DBw+WxWJRenq6Bg4cqODgYK1Zs0YDBgzQvHnzsvQ9ceJEDR8+XF9++aXGjRunDz/8UBs3brTd5pNPPtF//vMfLV68WOXLl7ct79+/vzp06KAOHTroww8/tC1fvXq1xowZow8//FDnzp3T9OnTbetOnTqlDRs2aMWKFRo0aJBmzpypyZMn65VXXtHMmTO1fv36LPedaeXKlZo3b56ee+45ffzxxwoMDNTTTz9tW9+4cWMlJCTor7/+ynEbAUBh4v4xCAAAhz3//POqX7++JOnpp5/W9OnT9Z///EfLli3TAw88oM6dO0uS+vbtq19++UXvv/++XnjhBTVo0ECPPvqo7rzzTklXv4yvXr1aZ86cUalSpSRdPddn9uzZGjlypDp27Jjj/X/xxRcqXry4nn/+eRmGoaeeeirX83COHj0qX19fRUZGKiQkRK+88oqSkpJkNptVtGhRSVeHpUlStWrV1LhxY9WpU0eSNGjQIM2fP19///23wsLCJEndunVT27ZtJUmPPvqo7Qv9jz/+qMTERE2ePFlBQUGqUqWKxowZI5PJpE8//VShoaH6z3/+I0kqX768jh07puXLl6tr167Zat6zZ0+uwencuXP65JNPtHTpUjVu3FjS1fOQWrZsqZ9++kmGYejEiRP64IMPFBQUpMqVK+uvv/7S559/LkkKCAjQpEmT9O9//1uSVKZMGb355pvav3+/bVnt2rWzHFnKFBgYKD8/P0lS8eLFbcuffPJJNWrUSJJ03333aeXKlbZ1aWlpGjlypCpUqKDSpUtr6tSp6t27t20bV69eXYcOHcp2X6tWrVK/fv1s+8BLL72kZcuWKSUlRX5+fvL19VXZsmW1Z88eVatWLcdtBQCFBeEJAAqwunXr2v5fq1YtJSYm6uzZszp48KC++OKLLEPd0tLS1KxZM0lS165dtWHDBn3wwQc6dOiQ/vjjD0myDRmTrh6dSk9Pt4WpnBw4cEDVqlWTYRi2ZbVr184ydC/Tgw8+qM8//1zNmjVTw4YN1bZtW3Xv3j3Hftu2bauffvpJr7zyig4dOqQ9e/Zkq+/aIzFBQUFKS0uTJB0+fFgVKlRQUFCQbX2PHj0kSa+++qr27dtnm3Qhs0+z2ZxjHYmJiSpWrFiO6/7++29ZLBZFR0fbloWEhKhChQo6ePCgMjIystVRp04dW3iqVauW/Pz8NGfOHB04cEB//vmnjhw5YnuOJCkyMjLH+87NHXfcYft/kSJFdOXKlSzry5YtK0m24HVt/35+ftmOPEpXt2fNmjVtf4eFhWnkyJFZbhMSEqIzZ87kqVYAKIgITwBQgHl7e9v+n3nOiWEYysjI0BNPPJHtaErml+YRI0Zox44d6tKli3r16qXw8HA9+OCDWW7bs2dPeXt7a+LEiWrSpIltooLrXT9xgbe3d47hqUqVKvr222/1/fff6/vvv9fMmTP12Wef5XjO0axZs7R69Wp1795dXbt21dixY7Mdgbn2sV/rRhM7pKenq0mTJnrppZdyvc21MrdlTnx9fXNcnpGRIYvFIrPZnG3bXPv35s2bNWTIEHXt2lXNmzfXkCFDNG7cOLvuIzcm041H41+/bW52+5za5MRisdjVFwAUdLzTAUABtnfvXtv/d+/erRIlSqhYsWKqUKGCjh49qnLlytn+rVq1Sj/88IMuXryozz77TLNmzdKwYcN0991369y5c5Kyfrm/++67NWTIECUnJ2vJkiU53n+VKlW0Z8+eLAHj2pqutXbtWn333Xfq0KGDXn31Vb3++uuKjY3VmTNnshy5kq6eZ/Piiy/queeeU8eOHW1hzJ4Z5sqXL68jR45kCXCvvvqqJk6cqAoVKujw4cMqU6aMbbvs3LlTK1asyLGv0NBQJSUl5biubNmy8vLy0s6dO23Lzp49qyNHjqhChQqqUqWK/v77b128eNG2PvMIn3T1/KQePXpo/Pjxuv/++1WpUiX973//s3sWveu3WX4pV66c9u3bZ/v77Nmzaty4sY4ePZplWeZwSgAozAhPAFCATZo0Sb///rt+/vlnvfbaa+rdu7ckqV+/flq/fr2WL1+u//3vf3rrrbf01ltvqXz58vLx8ZG/v7++/vprHT16VJs3b7ZdY+n6YVtBQUEaPny4li5dmuXLcqZOnTopOTlZkyZN0qFDh2yBKCcXLlzQpEmTtGXLFsXHx+vTTz9VyZIlVaxYMfn7+0u6GgCvXLmikJAQfffdd4qPj9evv/6qESNG5FhfTpo1a6awsDC99NJLOnjwoDZu3KiVK1eqWbNmuvfee5WSkmJbt2nTJk2aNEmhoaE59lWjRo1cJ0IIDAzU/fffrwkTJmjbtm3at2+fnn/+eZUsWVJNmzZVkyZNVKpUKb344os6ePCgvvzyyyyz/oWEhGjHjh36888/tX//fr3wwgs6ffq0XY9Rkvz9/XXs2DGdPHnSrts7qk+fPnr77be1YcMGHT58WGPHjlWZMmVUpkwZSVdnbjx27FiWoX0AUFgRngCgAOvYsaMGDhyo4cOH6/7779eAAQMkXT23ZurUqXrvvffUsWNHffDBB5oxY4YaNGggHx8fTZs2TV999ZU6deqkV155RU8++aTCw8NzPGrUrVs33XnnnbbZ7K5VtGhRvf766/r999/VpUsX/fzzz+rSpUuOtfbu3Vtdu3bV888/r44dO2rPnj1auHChzGazqlatqqZNm6pnz57atGmTJk+erL1796pTp04aNWqU2rdvr6ioqFyPal3Ly8tLCxYs0KlTp9StWzdNmjRJI0aMUMuWLRUUFKSlS5fq77//VteuXTVmzBj17t1bAwcOzLGvRo0a6fz58zp8+HCO60eOHKm77rpLw4YNU69eveTr66u33npLPj4+MplMmjt3rk6ePKkuXbpowYIF6t69u2244dChQxUaGqoHH3xQjz76qHx9fdWrVy+7HqMkdenSRYcPH9a9996br9d86tKli/r3769x48ape/fuunLliubMmWNbv2PHDpUsWVKVK1fOtxoAwFMYVk++MiEAIEdHjx5VmzZttHHjRtsRAOSPF154QWXKlNHQoUPz1O7MmTPas2dPlineX3/9dW3atCnXYYIF0ahRo1S2bFkNHjzY3aUAQL7jyBMAADfw+OOPa926dbbZ/PLiySef1Hvvvadjx47p559/1ttvv6327dvnQ5XucfbsWf3000/q1auXu0sBgFuC8AQAwA1UrlxZd999tz766KM8tQsNDdXs2bP1/vvvq3379ho9erQefvhhPfTQQ/lU6a33xhtv6Mknn8x1OncAKGwYtgcAAAAAduDIEwAAAADYgfAEAAAAAHYgPAEAAACAHQhPAAAAAGAHwhMAAAAA2IHwBAAAAAB2IDwBAAAAgB0ITwAAAABgh/8HfNJ+ESyQbY0AAAAASUVORK5CYII="
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": "9"
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def find_bin(avg_dist):\n",
    "    return np.digitize(logaritmare(avg_dist),bins)\n",
    "\n",
    "def logaritmare(x):\n",
    "    if x==0:\n",
    "        x=1000\n",
    "    return math.log(x)\n",
    "\n",
    "vfunc = np.vectorize(logaritmare)\n",
    "\n",
    "X=[np.array(vfunc(x)) for x in quality_peak_list]\n",
    "\n",
    "plt.figure(figsize=(10, 6))\n",
    "\n",
    "nb,bins,_ = plt.hist(X,bins=40, stacked=True ,histtype='bar',cumulative=False)\n",
    "\n",
    "#floating point error: when an image gets a value exactly at the boundary, floating point precision might cause it to overshoot the bin\n",
    "bins[-1]+=1\n",
    "bins[0]-=1\n",
    "\n",
    "\n",
    "plt.title(\"Peak distance distribution for 10 quality levels (categories follow matplotlib coloring)\")\n",
    "plt.xlabel(\"peak distance (logarithmic)\")\n",
    "plt.ylabel(\"number of images\")\n",
    "plt.show()\n",
    "\n",
    "\n",
    "def most_likely_category(avg_dist):\n",
    "    # when I get a new image see which bin it is in\n",
    "    # P(cat1| bin x) = P(binx|cat1)*P(cat1) / P(binx)\n",
    "    # greatest prob wins\n",
    "    # p1 = nr_binx_cat1/nr_cat1 * nr_cat1/ nr_all_cat / nr_binx_total * nr_all_cat\n",
    "    # p1 = nr_binx_cat1/nr_binx_total\n",
    "    # since nr_binx_total is the same for all probabilities, we can just compare the number of images in each bin\n",
    "    bin=find_bin(avg_dist)\n",
    "    candidate_bin = nb[0][bin-1]\n",
    "    ind=0\n",
    "    for i in range(1,10):\n",
    "        if nb[i][bin-1]-nb[i-1][bin-1] > candidate_bin:\n",
    "            candidate_bin = nb[i][bin-1]-nb[i-1][bin-1]\n",
    "            ind=i\n",
    "    return ind+1\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "end_time": "2024-11-12T09:50:56.045842800Z",
     "start_time": "2024-11-12T09:50:45.896667100Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "85\n",
      "40 40\n",
      "0\n",
      "100 100\n",
      "63\n",
      "30 30\n",
      "131\n",
      "60 60\n",
      "147\n",
      "70 70\n",
      "244\n",
      "80 80\n",
      "55\n",
      "30 30\n",
      "147\n",
      "70 70\n",
      "603\n",
      "90 90\n",
      "113\n",
      "50 60\n",
      "91\n",
      "40 50\n",
      "67\n",
      "30 40\n",
      "129\n",
      "70 60\n",
      "129\n",
      "70 60\n",
      "129\n",
      "70 60\n",
      "99\n",
      "60 50\n",
      "506\n",
      "90 90\n",
      "0\n",
      "100 100\n",
      "151\n",
      "70 70\n",
      "94\n",
      "50 50\n",
      "558\n",
      "90 90\n",
      "84\n",
      "40 40\n",
      "168\n",
      "70 70\n",
      "129\n",
      "60 60\n",
      "269\n",
      "80 80\n",
      "20\n",
      "10 10\n",
      "60\n",
      "30 30\n",
      "41\n",
      "20 20\n",
      "0\n",
      "100 100\n",
      "22\n",
      "10 10\n",
      "67\n",
      "30 40\n",
      "90\n",
      "40 50\n",
      "140\n",
      "70 70\n",
      "18\n",
      "10 10\n",
      "0\n",
      "100 100\n",
      "70\n",
      "40 40\n",
      "102\n",
      "50 50\n",
      "102\n",
      "50 50\n",
      "546\n",
      "90 90\n",
      "273\n",
      "80 80\n",
      "29\n",
      "90 20\n",
      "21\n",
      "10 10\n",
      "42\n",
      "20 20\n",
      "167\n",
      "70 70\n",
      "0\n",
      "100 100\n",
      "114\n",
      "60 60\n",
      "114\n",
      "60 60\n",
      "54\n",
      "30 30\n",
      "0\n",
      "100 100\n",
      "33\n",
      "20 20\n",
      "0\n",
      "100 100\n",
      "222\n",
      "80 80\n",
      "232\n",
      "80 80\n",
      "87\n",
      "50 40\n",
      "69\n",
      "40 40\n",
      "87\n",
      "50 40\n",
      "268\n",
      "80 80\n",
      "538\n",
      "90 90\n",
      "20\n",
      "10 10\n",
      "80\n",
      "40 40\n",
      "60\n",
      "30 30\n",
      "126\n",
      "60 60\n",
      "0\n",
      "90 100\n",
      "102\n",
      "50 50\n",
      "18\n",
      "40 10\n",
      "22\n",
      "10 10\n",
      "21\n",
      "70 10\n",
      "18\n",
      "60 10\n",
      "17\n",
      "10 10\n",
      "19\n",
      "30 10\n",
      "19\n",
      "60 10\n",
      "19\n",
      "30 10\n",
      "39\n",
      "20 20\n",
      "59\n",
      "30 30\n",
      "23\n",
      "90 10\n",
      "23\n",
      "100 10\n",
      "121\n",
      "60 60\n",
      "20\n",
      "10 10\n",
      "37\n",
      "20 20\n",
      "0\n",
      "100 100\n",
      "63\n",
      "30 30\n",
      "106\n",
      "50 50\n",
      "63\n",
      "30 30\n",
      "85\n",
      "40 40\n",
      "38\n",
      "20 20\n",
      "38\n",
      "20 20\n",
      "75\n",
      "40 40\n",
      "38\n",
      "20 20\n",
      "138\n",
      "60 70\n",
      "45\n",
      "20 20\n",
      "299\n",
      "80 80\n",
      "112\n",
      "50 60\n",
      "122\n",
      "60 60\n",
      "59\n",
      "30 30\n",
      "122\n",
      "60 60\n",
      "263\n",
      "80 80\n",
      "96\n",
      "50 50\n",
      "156\n",
      "70 70\n",
      "78\n",
      "40 40\n",
      "522\n",
      "90 90\n",
      "52\n",
      "30 30\n",
      "468\n",
      "90 90\n",
      "108\n",
      "60 50\n",
      "52\n",
      "30 30\n",
      "0\n",
      "100 100\n",
      "88\n",
      "40 40\n",
      "23\n",
      "30 10\n",
      "109\n",
      "50 50\n",
      "125\n",
      "60 60\n",
      "271\n",
      "80 80\n",
      "20\n",
      "10 10\n",
      "82\n",
      "40 40\n",
      "0\n",
      "100 100\n",
      "127\n",
      "60 60\n",
      "20\n",
      "10 10\n",
      "549\n",
      "90 90\n",
      "534\n",
      "90 90\n",
      "160\n",
      "70 70\n",
      "534\n",
      "90 90\n",
      "101\n",
      "50 50\n",
      "0\n",
      "100 100\n",
      "0\n",
      "90 100\n",
      "0\n",
      "100 100\n",
      "0\n",
      "90 100\n",
      "0\n",
      "100 100\n",
      "83\n",
      "40 40\n",
      "0\n",
      "90 100\n",
      "21\n",
      "10 10\n",
      "90\n",
      "40 50\n",
      "299\n",
      "80 80\n",
      "23\n",
      "10 10\n",
      "90\n",
      "40 50\n",
      "545\n",
      "90 90\n",
      "126\n",
      "60 60\n",
      "20\n",
      "10 10\n",
      "0\n",
      "100 100\n",
      "461\n",
      "90 90\n",
      "230\n",
      "80 80\n",
      "230\n",
      "80 80\n",
      "230\n",
      "80 80\n",
      "90\n",
      "40 50\n",
      "180\n",
      "70 70\n",
      "65\n",
      "30 30\n",
      "90\n",
      "40 50\n",
      "0\n",
      "100 100\n",
      "21\n",
      "10 10\n",
      "104\n",
      "50 50\n",
      "128\n",
      "60 60\n",
      "122\n",
      "60 60\n",
      "0\n",
      "100 100\n",
      "0\n",
      "100 100\n",
      "122\n",
      "60 60\n",
      "561\n",
      "90 90\n",
      "129\n",
      "60 60\n",
      "0\n",
      "100 100\n",
      "106\n",
      "50 50\n",
      "22\n",
      "10 10\n",
      "85\n",
      "40 40\n",
      "85\n",
      "40 40\n",
      "289\n",
      "80 80\n",
      "0\n",
      "100 100\n",
      "73\n",
      "40 40\n",
      "0\n",
      "100 100\n",
      "55\n",
      "30 30\n",
      "467\n",
      "90 90\n",
      "234\n",
      "80 80\n",
      "0\n",
      "100 100\n",
      "70\n",
      "40 40\n",
      "120\n",
      "60 60\n",
      "0\n",
      "100 100\n",
      "19\n",
      "10 10\n",
      "78\n",
      "40 40\n",
      "20\n",
      "10 10\n",
      "119\n",
      "60 60\n",
      "257\n",
      "80 80\n",
      "154\n",
      "70 70\n",
      "138\n",
      "70 70\n",
      "51\n",
      "30 30\n",
      "106\n",
      "60 50\n",
      "51\n",
      "30 30\n",
      "94\n",
      "60 50\n",
      "204\n",
      "80 80\n",
      "31\n",
      "20 20\n",
      "204\n",
      "80 80\n",
      "138\n",
      "70 70\n",
      "138\n",
      "70 70\n",
      "33\n",
      "20 20\n",
      "0\n",
      "100 100\n",
      "56\n",
      "30 30\n",
      "75\n",
      "40 40\n",
      "0\n",
      "100 100\n",
      "150\n",
      "70 70\n",
      "418\n",
      "90 90\n",
      "46\n",
      "30 20\n",
      "97\n",
      "60 50\n",
      "16\n",
      "10 100\n",
      "79\n",
      "40 40\n",
      "529\n",
      "90 90\n",
      "20\n",
      "10 10\n",
      "40\n",
      "20 20\n",
      "145\n",
      "70 70\n",
      "54\n",
      "30 30\n",
      "35\n",
      "20 20\n",
      "71\n",
      "40 40\n",
      "0\n",
      "100 100\n",
      "113\n",
      "60 60\n",
      "245\n",
      "80 80\n",
      "19\n",
      "10 10\n",
      "0\n",
      "100 100\n",
      "0\n",
      "100 100\n",
      "0\n",
      "100 100\n",
      "118\n",
      "60 60\n",
      "472\n",
      "90 90\n",
      "472\n",
      "90 90\n",
      "71\n",
      "40 40\n",
      "35\n",
      "20 20\n",
      "42\n",
      "20 20\n",
      "105\n",
      "50 50\n",
      "563\n",
      "90 90\n",
      "84\n",
      "40 40\n",
      "132\n",
      "60 60\n",
      "132\n",
      "60 60\n",
      "132\n",
      "60 60\n",
      "170\n",
      "70 70\n",
      "104\n",
      "50 50\n",
      "166\n",
      "70 70\n",
      "62\n",
      "30 30\n",
      "21\n",
      "10 10\n",
      "80\n",
      "40 40\n",
      "123\n",
      "60 60\n",
      "60\n",
      "30 30\n",
      "101\n",
      "50 50\n",
      "268\n",
      "80 80\n",
      "80\n",
      "40 40\n",
      "101\n",
      "50 50\n",
      "40\n",
      "20 20\n",
      "160\n",
      "70 70\n",
      "40\n",
      "20 20\n",
      "160\n",
      "70 70\n",
      "160\n",
      "70 70\n",
      "479\n",
      "90 90\n",
      "72\n",
      "40 40\n",
      "239\n",
      "80 80\n",
      "36\n",
      "20 20\n",
      "0\n",
      "100 100\n",
      "0\n",
      "100 100\n",
      "0\n",
      "100 100\n",
      "0\n",
      "100 100\n",
      "17\n",
      "10 10\n",
      "17\n",
      "10 10\n",
      "34\n",
      "20 20\n",
      "68\n",
      "40 40\n",
      "34\n",
      "60 20\n",
      "86\n",
      "40 40\n",
      "46\n",
      "70 20\n",
      "569\n",
      "90 90\n",
      "120\n",
      "60 60\n",
      "260\n",
      "80 80\n",
      "78\n",
      "40 40\n",
      "19\n",
      "10 10\n",
      "61\n",
      "30 30\n",
      "103\n",
      "50 50\n",
      "127\n",
      "60 60\n",
      "103\n",
      "50 50\n",
      "0\n",
      "100 100\n",
      "85\n",
      "40 40\n",
      "21\n",
      "10 10\n",
      "107\n",
      "50 50\n",
      "0\n",
      "90 100\n",
      "143\n",
      "70 70\n",
      "0\n",
      "100 100\n",
      "237\n",
      "80 80\n",
      "241\n",
      "80 80\n",
      "36\n",
      "20 20\n",
      "484\n",
      "90 90\n",
      "19\n",
      "10 10\n",
      "84\n",
      "40 40\n",
      "169\n",
      "70 70\n",
      "84\n",
      "40 40\n",
      "84\n",
      "40 40\n",
      "59\n",
      "30 30\n",
      "526\n",
      "90 90\n",
      "0\n",
      "100 100\n",
      "79\n",
      "40 40\n",
      "0\n",
      "100 100\n",
      "45\n",
      "20 20\n",
      "600\n",
      "90 90\n",
      "139\n",
      "60 70\n",
      "128\n",
      "60 60\n",
      "62\n",
      "30 30\n",
      "168\n",
      "70 70\n",
      "105\n",
      "50 50\n",
      "124\n",
      "60 60\n",
      "0\n",
      "100 100\n",
      "41\n",
      "20 20\n",
      "80\n",
      "40 40\n",
      "119\n",
      "60 60\n",
      "119\n",
      "60 60\n",
      "119\n",
      "60 60\n",
      "260\n",
      "80 80\n",
      "98\n",
      "50 50\n",
      "39\n",
      "20 20\n",
      "20\n",
      "10 10\n",
      "258\n",
      "80 80\n",
      "44\n",
      "30 20\n",
      "19\n",
      "10 10\n",
      "64\n",
      "40 30\n",
      "102\n",
      "60 50\n",
      "38\n",
      "20 20\n",
      "38\n",
      "20 20\n",
      "118\n",
      "60 60\n",
      "96\n",
      "50 50\n",
      "18\n",
      "100 10\n",
      "18\n",
      "90 10\n",
      "121\n",
      "60 60\n",
      "101\n",
      "50 50\n",
      "60\n",
      "30 30\n",
      "40\n",
      "20 20\n",
      "20\n",
      "10 10\n",
      "40\n",
      "20 20\n",
      "102\n",
      "50 50\n",
      "20\n",
      "10 10\n",
      "544\n",
      "90 90\n",
      "273\n",
      "80 80\n",
      "0\n",
      "100 100\n",
      "61\n",
      "40 30\n",
      "0\n",
      "100 100\n",
      "75\n",
      "50 40\n",
      "20\n",
      "10 10\n",
      "94\n",
      "50 50\n",
      "250\n",
      "80 80\n",
      "20\n",
      "10 10\n",
      "124\n",
      "60 60\n",
      "124\n",
      "60 60\n",
      "59\n",
      "30 30\n",
      "59\n",
      "30 30\n",
      "66\n",
      "40 40\n",
      "34\n",
      "20 20\n",
      "17\n",
      "10 10\n",
      "66\n",
      "40 40\n",
      "20\n",
      "10 10\n",
      "548\n",
      "90 90\n",
      "164\n",
      "70 70\n",
      "126\n",
      "60 60\n",
      "518\n",
      "90 90\n",
      "259\n",
      "80 80\n",
      "19\n",
      "10 10\n",
      "78\n",
      "40 40\n",
      "116\n",
      "60 60\n",
      "94\n",
      "50 50\n",
      "56\n",
      "30 30\n",
      "502\n",
      "90 90\n",
      "494\n",
      "90 90\n",
      "0\n",
      "100 100\n",
      "0\n",
      "100 100\n",
      "37\n",
      "20 20\n",
      "0\n",
      "100 100\n",
      "287\n",
      "80 80\n",
      "132\n",
      "60 60\n",
      "287\n",
      "80 80\n",
      "88\n",
      "40 40\n",
      "108\n",
      "50 50\n",
      "0\n",
      "90 100\n",
      "178\n",
      "70 70\n",
      "79\n",
      "40 40\n",
      "39\n",
      "20 20\n",
      "39\n",
      "20 20\n",
      "263\n",
      "80 80\n",
      "79\n",
      "40 40\n",
      "530\n",
      "90 90\n",
      "40\n",
      "20 20\n",
      "79\n",
      "40 40\n",
      "17\n",
      "10 10\n",
      "217\n",
      "80 80\n",
      "64\n",
      "40 30\n",
      "0\n",
      "100 100\n",
      "23\n",
      "10 10\n",
      "181\n",
      "70 70\n",
      "138\n",
      "60 70\n",
      "181\n",
      "70 70\n",
      "64\n",
      "30 30\n",
      "22\n",
      "10 10\n",
      "0\n",
      "90 100\n",
      "22\n",
      "10 10\n",
      "467\n",
      "90 90\n",
      "0\n",
      "100 100\n",
      "467\n",
      "90 90\n",
      "107\n",
      "60 50\n",
      "93\n",
      "50 50\n",
      "149\n",
      "70 70\n",
      "494\n",
      "90 90\n",
      "113\n",
      "60 60\n",
      "23\n",
      "10 10\n",
      "42\n",
      "20 20\n",
      "291\n",
      "80 80\n",
      "66\n",
      "30 40\n",
      "59\n",
      "30 30\n",
      "157\n",
      "70 70\n",
      "121\n",
      "60 60\n",
      "39\n",
      "20 20\n",
      "average error: 3.075 \n"
     ]
    }
   ],
   "source": [
    "number_of_tests = 100\n",
    "testing_cutoff = 700 # images from 700 onwards are for testing\n",
    "\n",
    "error=0\n",
    "\n",
    "for i in range(0, number_of_tests):\n",
    "        for j in range(1, 5):\n",
    "            \n",
    "            file_name = f\"{testing_cutoff + i + 1:04}x4w1.png\"\n",
    "            path = os.path.join(\"DIV2K_train_LR_wild\", file_name)\n",
    "            \n",
    "            qual_cat = np.random.randint(1,11)\n",
    "        \n",
    "            img = compress(path, qual_cat*10)\n",
    "            img = img.convert('YCbCr')\n",
    "            img = zero_pad(np.array(img)[:,:,0])\n",
    "            img = img.astype(np.float64)\n",
    "            dct = DCT(img)[0,0,:]\n",
    "            hcount, bin_edges = hist(dct)\n",
    "            fft = FFT(hcount)\n",
    "            _, avg_peak_dist,_ = quality(fft)\n",
    "            print(avg_peak_dist)\n",
    "            cat=most_likely_category(avg_peak_dist)\n",
    "            print(qual_cat*10, cat*10)\n",
    "            error += np.abs(cat*10-qual_cat*10)\n",
    "            \n",
    "            \n",
    "print(f\"average error: {error / number_of_tests / 4} \") # we do 4 samples for each test"
   ]
  }
 ],
 "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.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}