Commit 60b8c173 authored by Kelly Chang's avatar Kelly Chang
parents 18a92d19 ac651fcd
......@@ -2,7 +2,7 @@
"cells": [
{
"cell_type": "code",
"execution_count": 130,
"execution_count": 33,
"id": "dbef8759",
"metadata": {
"id": "dbef8759"
......@@ -21,12 +21,13 @@
"from numpy import linalg as la\n",
"from scipy.stats import gaussian_kde\n",
"import seaborn as sns\n",
"import pywt"
"import pywt\n",
"from collections import Counter"
]
},
{
"cell_type": "code",
"execution_count": 364,
"execution_count": 2,
"id": "9ed20f84",
"metadata": {
"id": "9ed20f84"
......@@ -43,7 +44,7 @@
" image = tiff_list[i]\n",
" image = Image.open(image) #Open the image and read it as an Image object\n",
" image = np.array(image)[1:,:] #Convert to an array, leaving out the first row because the first row is just housekeeping data\n",
" image_int = image.astype(np.int_)\n",
" image_int = image.astype(int)\n",
" \n",
" A = np.array([[3,0,-1],[0,3,3],[1,-3,-4]]) # the matrix for system of equation\n",
" \n",
......@@ -59,7 +60,8 @@
" y = np.vstack((y0,y1,y2))\n",
" \n",
" # use numpy solver to solve the system of equations all at once\n",
" predict = np.linalg.solve(A,y)[-1]\n",
" #predict = np.linalg.solve(A,y)[-1]\n",
" predict = np.floor(np.linalg.solve(A,y)[-1]).astype(int)\n",
" #predict = []\n",
" \n",
" # flatten the neighbor pixels and stack them together\n",
......@@ -71,7 +73,7 @@
" \n",
" # calculate the difference\n",
" diff = np.max(neighbor,axis = 1) - np.min(neighbor, axis=1)\n",
" \n",
" diff = np.pad(diff.reshape(510,638), pad_width=1)\n",
" \n",
" # flatten the image to a vector\n",
" small_image = image_int[1:-1,1:-1]\n",
......@@ -98,19 +100,19 @@
},
{
"cell_type": "code",
"execution_count": 202,
"execution_count": 3,
"id": "ba2881d9",
"metadata": {},
"outputs": [],
"source": [
"scenes = file_extractor()\n",
"images = image_extractor(scenes)\n",
"num_images = im_distribution(images, \"_9\")"
"num_images = im_distribution(images, \"11\")"
]
},
{
"cell_type": "code",
"execution_count": 365,
"execution_count": 4,
"id": "11e95c34",
"metadata": {},
"outputs": [],
......@@ -120,7 +122,7 @@
},
{
"cell_type": "code",
"execution_count": 402,
"execution_count": 5,
"id": "434e4d2f",
"metadata": {},
"outputs": [],
......@@ -135,7 +137,7 @@
" error (array): matrix of errors computed in encoding. Same \n",
" shape as the original image (512, 640) in this case\n",
" A (array): Matrix used for the system of equations to create predictions\n",
" Returns: \n",
" Returns: cd cdcd\n",
" image (array): The reconstructed image\n",
" \"\"\"\n",
" new_e = error.copy()\n",
......@@ -146,62 +148,316 @@
" z0, z1, z2, z3 = new_e[r-1][c-1], new_e[r-1][c], new_e[r-1][c+1], new_e[r][c-1]\n",
" y = np.vstack((-z0+z2-z3, z0+z1+z2, -z0-z1-z2-z3))\n",
"\n",
" if r == 345 and c == 421:\n",
" print(new_e[r][c])\n",
" print(np.linalg.solve(A,y)[-1])\n",
" print(new_e[r][c] + np.linalg.solve(A,y)[-1])\n",
" print(np.ceil(new_e[r][c]) + np.floor(np.linalg.solve(A,y)[-1]))\n",
" \n",
" \n",
" #Real solution that works, DO NOT DELETE\n",
" new_e[r][c] = int(np.ceil(new_e[r][c] + np.linalg.solve(A,y)[-1])) \n",
" #new_e[r][c] = int(np.ceil(new_e[r][c] + np.linalg.solve(A,y)[-1]))\n",
" \n",
" #new_e[r][c] = np.ceil(new_e[r][c]) + np.floor(np.linalg.solve(A,y)[-1])\n",
" new_e[r][c] = np.round(new_e[r][c] + np.linalg.solve(A,y)[-1], 1)\n",
" \n",
" return new_e.astype(int)\n",
" "
" "
]
},
{
"cell_type": "code",
"execution_count": 403,
"id": "7f395ab2",
"execution_count": 6,
"id": "4f4a5a35",
"metadata": {},
"outputs": [],
"source": [
"new_error = reconstruct(err, A)"
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "6d95ffce",
"metadata": {},
"outputs": [],
"source": [
"first = []\n",
"second = []\n",
"third = []\n",
"fourth = []\n",
"\n",
"for i in range(1,diff.shape[0]-1):\n",
" for j in range(1,diff.shape[1]-1):\n",
" if diff[i][j] <= 50:\n",
" first.append(np.abs(err[i][j]))\n",
" elif diff[i][j] > 50 and diff[i][j] <= 100:\n",
" second.append(np.abs(err[i][j]))\n",
" elif diff[i][j] > 100 and diff[i][j] <= 200:\n",
" third.append(np.abs(err[i][j]))\n",
" else:\n",
" fourth.append(np.abs(err[i][j]))\n"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "1c848109",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAD4CAYAAADsKpHdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAQeklEQVR4nO3dcayddX3H8fdnrSJgigUKqW2zW0OjQjOnNKzqYpbVhCrG8gckdxmj2Zo0IUzRmLh2/mH2R5OSGRGSwdKAUpAITWWjkeAkRbMsYWUXMUKpHXeW0Uql14HItoAWv/vj/G52erm999zb0nMPfb+Sk/Oc7/P8nvP70nv43Od5zjk3VYUkSb/T7wlIkuYGA0GSBBgIkqTGQJAkAQaCJKmZ3+8JzNb5559fQ0ND/Z6GJA2Uxx9//BdVtWiydQMbCENDQ4yMjPR7GpI0UJL85/HWecpIkgQYCJKkxkCQJAEGgiSpMRAkSYCBIElqDARJEmAgSJIaA0GSBAzwJ5VPxNCmB/v23M9uvaJvzy1JU/EIQZIEGAiSpMZAkCQBBoIkqTEQJEmAgSBJagwESRJgIEiSGgNBkgQYCJKkxkCQJAEGgiSpMRAkSYCBIElqDARJEmAgSJKangIhyeeT7E3yVJJvJXlHknOTPJzkmXa/sGv7zUlGk+xPcnlX/dIkT7Z1tyRJq5+R5L5W35Nk6KR3Kkma0rSBkGQJ8FlgVVWtBOYBw8AmYHdVrQB2t8ckubitvwRYC9yaZF7b3W3ARmBFu61t9Q3AS1V1EXATcONJ6U6S1LNeTxnNB85MMh84C3geWAdsb+u3A1e25XXAvVX1WlUdAEaBy5IsBhZU1aNVVcBdE8aM72snsGb86EGSdGpMGwhV9TPgK8BzwGHg5ar6HnBhVR1u2xwGLmhDlgAHu3ZxqNWWtOWJ9WPGVNVR4GXgvIlzSbIxyUiSkbGxsV57lCT1oJdTRgvp/Aa/HHg3cHaSa6YaMkmtpqhPNebYQtW2qlpVVasWLVo09cQlSTPSyymjjwMHqmqsqn4D3A98BHihnQai3R9p2x8ClnWNX0rnFNOhtjyxfsyYdlrqHODF2TQkSZqdXgLhOWB1krPaef01wD5gF7C+bbMeeKAt7wKG2zuHltO5ePxYO630SpLVbT/XThgzvq+rgEfadQZJ0ikyf7oNqmpPkp3AD4GjwBPANuCdwI4kG+iExtVt+71JdgBPt+2vr6rX2+6uA+4EzgQeajeAO4C7k4zSOTIYPindSZJ6Nm0gAFTVl4EvTyi/RudoYbLttwBbJqmPACsnqb9KCxRJUn/4SWVJEmAgSJIaA0GSBBgIkqTGQJAkAQaCJKkxECRJgIEgSWoMBEkSYCBIkhoDQZIEGAiSpMZAkCQBBoIkqTEQJEmAgSBJagwESRJgIEiSGgNBkgQYCJKkxkCQJAEGgiSpMRAkSYCBIElqDARJEmAgSJIaA0GSBBgIkqTGQJAkAQaCJKkxECRJgIEgSWoMBEkSYCBIkhoDQZIEGAiSpMZAkCQBBoIkqekpEJK8K8nOJD9Jsi/Jh5Ocm+ThJM+0+4Vd229OMppkf5LLu+qXJnmyrbslSVr9jCT3tfqeJEMnvVNJ0pR6PUK4GfhuVb0P+ACwD9gE7K6qFcDu9pgkFwPDwCXAWuDWJPPafm4DNgIr2m1tq28AXqqqi4CbgBtPsC9J0gxNGwhJFgAfA+4AqKpfV9UvgXXA9rbZduDKtrwOuLeqXquqA8AocFmSxcCCqnq0qgq4a8KY8X3tBNaMHz1Ikk6NXo4Q3gOMAd9I8kSS25OcDVxYVYcB2v0FbfslwMGu8YdabUlbnlg/ZkxVHQVeBs6bOJEkG5OMJBkZGxvrsUVJUi/m97jNh4DPVNWeJDfTTg8dx2S/2dcU9anGHFuo2gZsA1i1atUb1g+CoU0P9uV5n916RV+eV9Lg6OUI4RBwqKr2tMc76QTEC+00EO3+SNf2y7rGLwWeb/Wlk9SPGZNkPnAO8OJMm5Ekzd60gVBVPwcOJnlvK60BngZ2AetbbT3wQFveBQy3dw4tp3Px+LF2WumVJKvb9YFrJ4wZ39dVwCPtOoMk6RTp5ZQRwGeAe5K8Hfgp8Od0wmRHkg3Ac8DVAFW1N8kOOqFxFLi+ql5v+7kOuBM4E3io3aBzwfruJKN0jgyGT7AvSdIM9RQIVfUjYNUkq9YcZ/stwJZJ6iPAyknqr9ICRZLUH35SWZIEGAiSpMZAkCQBBoIkqTEQJEmAgSBJagwESRJgIEiSGgNBkgQYCJKkxkCQJAEGgiSpMRAkSYCBIElqDARJEmAgSJIaA0GSBBgIkqTGQJAkAQaCJKkxECRJgIEgSWoMBEkSYCBIkhoDQZIEGAiSpMZAkCQBBoIkqTEQJEmAgSBJagwESRJgIEiSGgNBkgQYCJKkxkCQJAEGgiSpMRAkScAMAiHJvCRPJPlOe3xukoeTPNPuF3ZtuznJaJL9SS7vql+a5Mm27pYkafUzktzX6nuSDJ3EHiVJPZjJEcINwL6ux5uA3VW1AtjdHpPkYmAYuARYC9yaZF4bcxuwEVjRbmtbfQPwUlVdBNwE3DirbiRJs9ZTICRZClwB3N5VXgdsb8vbgSu76vdW1WtVdQAYBS5LshhYUFWPVlUBd00YM76vncCa8aMHSdKp0esRwteALwK/7apdWFWHAdr9Ba2+BDjYtd2hVlvSlifWjxlTVUeBl4HzJk4iycYkI0lGxsbGepy6JKkX0wZCkk8BR6rq8R73Odlv9jVFfaoxxxaqtlXVqqpatWjRoh6nI0nqxfwetvko8OkknwTeASxI8k3ghSSLq+pwOx10pG1/CFjWNX4p8HyrL52k3j3mUJL5wDnAi7PsSZI0C9MeIVTV5qpaWlVDdC4WP1JV1wC7gPVts/XAA215FzDc3jm0nM7F48faaaVXkqxu1weunTBmfF9Xted4wxGCJOnN08sRwvFsBXYk2QA8B1wNUFV7k+wAngaOAtdX1ettzHXAncCZwEPtBnAHcHeSUTpHBsMnMC9J0izMKBCq6gfAD9ryfwFrjrPdFmDLJPURYOUk9VdpgSJJ6g8/qSxJAgwESVJjIEiSAANBktQYCJIkwECQJDUGgiQJMBAkSY2BIEkCDARJUmMgSJIAA0GS1JzIt51qgAxterBvz/3s1iv69tySeucRgiQJMBAkSY2BIEkCDARJUmMgSJIAA0GS1BgIkiTAQJAkNQaCJAkwECRJjYEgSQIMBElSYyBIkgADQZLUGAiSJMBAkCQ1BoIkCTAQJEmNgSBJAgwESVJjIEiSAANBktQYCJIkwECQJDUGgiQJ6CEQkixL8v0k+5LsTXJDq5+b5OEkz7T7hV1jNicZTbI/yeVd9UuTPNnW3ZIkrX5GkvtafU+SoTehV0nSFHo5QjgKfKGq3g+sBq5PcjGwCdhdVSuA3e0xbd0wcAmwFrg1yby2r9uAjcCKdlvb6huAl6rqIuAm4MaT0JskaQamDYSqOlxVP2zLrwD7gCXAOmB722w7cGVbXgfcW1WvVdUBYBS4LMliYEFVPVpVBdw1Ycz4vnYCa8aPHiRJp8aMriG0UzkfBPYAF1bVYeiEBnBB22wJcLBr2KFWW9KWJ9aPGVNVR4GXgfMmef6NSUaSjIyNjc1k6pKkafQcCEneCXwb+FxV/WqqTSep1RT1qcYcW6jaVlWrqmrVokWLppuyJGkGegqEJG+jEwb3VNX9rfxCOw1Euz/S6oeAZV3DlwLPt/rSSerHjEkyHzgHeHGmzUiSZq+XdxkFuAPYV1Vf7Vq1C1jfltcDD3TVh9s7h5bTuXj8WDut9EqS1W2f104YM76vq4BH2nUGSdIpMr+HbT4K/BnwZJIftdpfA1uBHUk2AM8BVwNU1d4kO4Cn6bxD6fqqer2Nuw64EzgTeKjdoBM4dycZpXNkMHxibUmSZmraQKiqf2Hyc/wAa44zZguwZZL6CLBykvqrtECRJPWHn1SWJAEGgiSpMRAkSYCBIElqDARJEmAgSJIaA0GSBBgIkqTGQJAkAQaCJKkxECRJgIEgSWoMBEkSYCBIkppe/h6CdEKGNj3Yl+d9dusVfXleaVB5hCBJAgwESVJjIEiSAANBktQYCJIkwECQJDUGgiQJMBAkSY2BIEkCDARJUmMgSJIAA0GS1BgIkiTAQJAkNQaCJAkwECRJjYEgSQL8i2l6C/MvtUkz4xGCJAkwECRJjYEgSQIMBElSYyBIkoA59C6jJGuBm4F5wO1VtbXPU5JmpV/vbgLf4aQTMyeOEJLMA/4O+ARwMfAnSS7u76wk6fQyV44QLgNGq+qnAEnuBdYBT/d1VtKA8bMXOhFzJRCWAAe7Hh8C/mDiRkk2Ahvbw/9Osn+Wz3c+8ItZjp0r7GHueCv0cUI95MaTOJPZO+3/HXr0u8dbMVcCIZPU6g2Fqm3AthN+smSkqlad6H76yR7mjrdCH/YwN/S7hzlxDYHOEcGyrsdLgef7NBdJOi3NlUD4N2BFkuVJ3g4MA7v6PCdJOq3MiVNGVXU0yV8C/0Tnbadfr6q9b+JTnvBppznAHuaOt0If9jA39LWHVL3hVL0k6TQ0V04ZSZL6zECQJAGnYSAkWZtkf5LRJJv6PZ9eJFmW5PtJ9iXZm+SGVj83ycNJnmn3C/s91+kkmZfkiSTfaY8Hqock70qyM8lP2r/Hhwewh8+3n6OnknwryTvmeg9Jvp7kSJKnumrHnXOSze01vj/J5f2Z9bGO08Pftp+lHyf5hyTv6lp3yns4rQJhgL8i4yjwhap6P7AauL7NexOwu6pWALvb47nuBmBf1+NB6+Fm4LtV9T7gA3R6GZgekiwBPgusqqqVdN7EMczc7+FOYO2E2qRzbq+NYeCSNubW9trvtzt5Yw8PAyur6veAfwc2Q/96OK0Cga6vyKiqXwPjX5Exp1XV4ar6YVt+hc7/hJbQmfv2ttl24Mq+TLBHSZYCVwC3d5UHpockC4CPAXcAVNWvq+qXDFAPzXzgzCTzgbPofOZnTvdQVf8MvDihfLw5rwPurarXquoAMErntd9Xk/VQVd+rqqPt4b/S+QwW9KmH0y0QJvuKjCV9msusJBkCPgjsAS6sqsPQCQ3ggj5OrRdfA74I/LarNkg9vAcYA77RTnvdnuRsBqiHqvoZ8BXgOeAw8HJVfY8B6qHL8eY8qK/zvwAeast96eF0C4SeviJjrkryTuDbwOeq6lf9ns9MJPkUcKSqHu/3XE7AfOBDwG1V9UHgf5h7p1am1M6zrwOWA+8Gzk5yTX9nddIN3Os8yZfonBq+Z7w0yWZveg+nWyAM7FdkJHkbnTC4p6rub+UXkixu6xcDR/o1vx58FPh0kmfpnKr74yTfZLB6OAQcqqo97fFOOgExSD18HDhQVWNV9RvgfuAjDFYP444354F6nSdZD3wK+NP6/w+G9aWH0y0QBvIrMpKEznnrfVX11a5Vu4D1bXk98MCpnluvqmpzVS2tqiE6/90fqaprGKwefg4cTPLeVlpD5yvaB6YHOqeKVic5q/1craFzTWqQehh3vDnvAoaTnJFkObACeKwP85tWOn8Y7K+AT1fV/3at6k8PVXVa3YBP0rma/x/Al/o9nx7n/Id0Dhd/DPyo3T4JnEfn3RXPtPtz+z3XHvv5I+A7bXmgegB+Hxhp/xb/CCwcwB7+BvgJ8BRwN3DGXO8B+Badax6/ofPb84ap5gx8qb3G9wOf6Pf8p+hhlM61gvHX9d/3swe/ukKSBJx+p4wkScdhIEiSAANBktQYCJIkwECQJDUGgiQJMBAkSc3/AZOjkH9qMpzuAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"3.499999999992724\n",
"[13644.5]\n",
"[13648.]\n",
"[13648.]\n"
"124\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAD4CAYAAADsKpHdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAQfElEQVR4nO3df6zd9V3H8edr7WToBuNHIU2LXibVCCSO0SDJ3LLIIh3MgQqmi0oTmzQSlmxR44pLdP5BAhqHIQoLykLBOahsC80IcaRsLiYIXhgMCqt0g0Gl0g6QsSho2ds/zvua09t7b29vS88pfT6Sk/M97/P9fM/7fM6B1/1+v+ecpqqQJOkto25AkjQeDARJEmAgSJKagSBJAgwESVJbPOoGFurEE0+siYmJUbchSYeVBx988PtVtWSm+w7bQJiYmGBycnLUbUjSYSXJ92a7z0NGkiTAQJAkNQNBkgQYCJKkZiBIkgADQZLUDARJEmAgSJKagSBJAg7jbyofiIn1d43ssZ+++sKRPbYkzcU9BEkSYCBIkpqBIEkCDARJUjMQJEmAgSBJagaCJAkwECRJzUCQJAEGgiSpGQiSJMBAkCQ1A0GSBBgIkqRmIEiSAANBktQMBEkSYCBIkpqBIEkCDARJUjMQJEmAgSBJagaCJAkwECRJbd6BkGRRkm8m+UrfPj7JPUme7Ovjhta9Msm2JFuTnD9UPzvJo33fdUnS9aOS3N71+5NMHMTnKEmah/3ZQ/g48MTQ7fXA5qpaAWzu2yQ5HVgNnAGsAq5PsqjH3ACsA1b0ZVXX1wIvVdVpwLXANQt6NpKkBZtXICRZDlwI/O1Q+SJgQy9vAC4eqt9WVa9V1VPANuCcJEuBY6rqvqoq4JZpY6a2dQdw3tTegyTp0JjvHsJfAn8I/GiodnJV7QDo65O6vgx4dmi97V1b1svT63uMqardwMvACdObSLIuyWSSyV27ds2zdUnSfOwzEJJ8GNhZVQ/Oc5sz/WVfc9TnGrNnoerGqlpZVSuXLFkyz3YkSfOxeB7rvBf4SJILgLcBxyT5O+D5JEurakcfDtrZ628HThkavxx4ruvLZ6gPj9meZDFwLPDiAp+TJGkB9rmHUFVXVtXyqppgcLL43qr6LWATsKZXWwPc2cubgNX9yaFTGZw8fqAPK72S5Nw+P3DZtDFT27qkH2OvPQRJ0htnPnsIs7ka2JhkLfAMcClAVW1JshF4HNgNXFFVr/eYy4GbgaOBu/sCcBNwa5JtDPYMVh9AX5KkBdivQKiqrwNf7+UXgPNmWe8q4KoZ6pPAmTPUX6UDRZI0Gn5TWZIEGAiSpGYgSJIAA0GS1AwESRJgIEiSmoEgSQIMBElSMxAkSYCBIElqBoIkCTAQJEnNQJAkAQaCJKkZCJIkwECQJDUDQZIEGAiSpGYgSJIAA0GS1AwESRJgIEiSmoEgSQIMBElSMxAkSYCBIElqBoIkCTAQJEnNQJAkAQaCJKkZCJIkwECQJDUDQZIEGAiSpGYgSJIAA0GS1PYZCEneluSBJI8k2ZLkT7t+fJJ7kjzZ18cNjbkyybYkW5OcP1Q/O8mjfd91SdL1o5Lc3vX7k0y8Ac9VkjSH+ewhvAb8UlX9PPBuYFWSc4H1wOaqWgFs7tskOR1YDZwBrAKuT7Kot3UDsA5Y0ZdVXV8LvFRVpwHXAtcc+FOTJO2PxftaoaoK+GHffGtfCrgI+EDXNwBfBz7Z9duq6jXgqSTbgHOSPA0cU1X3ASS5BbgYuLvHfLq3dQfwV0nSj/2mMrH+rpE87tNXXziSx5V0+JjXOYQki5I8DOwE7qmq+4GTq2oHQF+f1KsvA54dGr69a8t6eXp9jzFVtRt4GThhhj7WJZlMMrlr1655PUFJ0vzMKxCq6vWqejewnMFf+2fOsXpm2sQc9bnGTO/jxqpaWVUrlyxZso+uJUn7Y78+ZVRV/8ng0NAq4PkkSwH6emevth04ZWjYcuC5ri+fob7HmCSLgWOBF/enN0nSgZnPp4yWJHlnLx8NfBD4NrAJWNOrrQHu7OVNwOr+5NCpDE4eP9CHlV5Jcm5/uuiyaWOmtnUJcO+b8fyBJI2zfZ5UBpYCG/qTQm8BNlbVV5LcB2xMshZ4BrgUoKq2JNkIPA7sBq6oqtd7W5cDNwNHMziZfHfXbwJu7RPQLzL4lJIk6RCaz6eMvgWcNUP9BeC8WcZcBVw1Q30S2Ov8Q1W9SgeKJGk0/KayJAkwECRJzUCQJAEGgiSpGQiSJMBAkCQ1A0GSBBgIkqRmIEiSAANBktQMBEkSYCBIkpqBIEkCDARJUjMQJEmAgSBJagaCJAkwECRJzUCQJAEGgiSpGQiSJMBAkCQ1A0GSBBgIkqRmIEiSAANBktQMBEkSYCBIkpqBIEkCDARJUjMQJEmAgSBJagaCJAkwECRJzUCQJAHzCIQkpyT5WpInkmxJ8vGuH5/kniRP9vVxQ2OuTLItydYk5w/Vz07yaN93XZJ0/agkt3f9/iQTb8BzlSTNYT57CLuB36+qnwPOBa5IcjqwHthcVSuAzX2bvm81cAawCrg+yaLe1g3AOmBFX1Z1fS3wUlWdBlwLXHMQnpskaT/sMxCqakdVPdTLrwBPAMuAi4ANvdoG4OJevgi4rapeq6qngG3AOUmWAsdU1X1VVcAt08ZMbesO4LypvQdJ0qGxX+cQ+lDOWcD9wMlVtQMGoQGc1KstA54dGra9a8t6eXp9jzFVtRt4GThhhsdfl2QyyeSuXbv2p3VJ0j7MOxCSvB34IvCJqvrBXKvOUKs56nON2bNQdWNVrayqlUuWLNlXy5Kk/TCvQEjyVgZh8Pmq+lKXn+/DQPT1zq5vB04ZGr4ceK7ry2eo7zEmyWLgWODF/X0ykqSFm8+njALcBDxRVZ8ZumsTsKaX1wB3DtVX9yeHTmVw8viBPqz0SpJze5uXTRszta1LgHv7PIMk6RBZPI913gv8NvBokoe79kfA1cDGJGuBZ4BLAapqS5KNwOMMPqF0RVW93uMuB24Gjgbu7gsMAufWJNsY7BmsPrCnJUnaX/sMhKr6Z2Y+xg9w3ixjrgKumqE+CZw5Q/1VOlAkSaPhN5UlSYCBIElqBoIkCTAQJEnNQJAkAQaCJKkZCJIkwECQJDUDQZIEGAiSpGYgSJIAA0GS1AwESRJgIEiSmoEgSQIMBElSMxAkSYCBIElqBoIkCTAQJEnNQJAkAQaCJKkZCJIkwECQJDUDQZIEGAiSpLZ41A3o0JhYf9fIHvvpqy8c2WNLmj/3ECRJgIEgSWoGgiQJMBAkSc1AkCQBBoIkqRkIkiTAQJAkNQNBkgTMIxCSfC7JziSPDdWOT3JPkif7+rih+65Msi3J1iTnD9XPTvJo33ddknT9qCS3d/3+JBMH+TlKkuZhPnsINwOrptXWA5uragWwuW+T5HRgNXBGj7k+yaIecwOwDljRl6ltrgVeqqrTgGuBaxb6ZCRJC7fPQKiqbwAvTitfBGzo5Q3AxUP126rqtap6CtgGnJNkKXBMVd1XVQXcMm3M1LbuAM6b2nuQJB06Cz2HcHJV7QDo65O6vgx4dmi97V1b1svT63uMqardwMvACTM9aJJ1SSaTTO7atWuBrUuSZnKwTyrP9Jd9zVGfa8zexaobq2plVa1csmTJAluUJM1koYHwfB8Goq93dn07cMrQesuB57q+fIb6HmOSLAaOZe9DVJKkN9hCA2ETsKaX1wB3DtVX9yeHTmVw8viBPqz0SpJz+/zAZdPGTG3rEuDePs8gSTqE9vkP5CT5AvAB4MQk24E/Aa4GNiZZCzwDXApQVVuSbAQeB3YDV1TV672pyxl8Yulo4O6+ANwE3JpkG4M9g9UH5ZlJkvbLPgOhqj46y13nzbL+VcBVM9QngTNnqL9KB4okaXT8prIkCTAQJEnNQJAkAQaCJKkZCJIkwECQJDUDQZIEGAiSpGYgSJIAA0GS1AwESRJgIEiSmoEgSQIMBElSMxAkSYCBIElqBoIkCTAQJEnNQJAkAQaCJKkZCJIkwECQJLXFo25Ab34T6+8ayeM+ffWFI3lc6XDlHoIkCTAQJEnNQJAkAQaCJKkZCJIkwECQJDUDQZIEGAiSpGYgSJIAA0GS1AwESRJgIEiSmj9upzctf1RP2j9js4eQZFWSrUm2JVk/6n4k6UgzFoGQZBHw18CHgNOBjyY5fbRdSdKRZVwOGZ0DbKuq7wIkuQ24CHh8pF1JCzCqQ1Xg4SodmHEJhGXAs0O3twO/MH2lJOuAdX3zh0m2LvDxTgS+v8Cxh5J9Hlxv+j5zzUHuZG6Hy3zC4dProejzp2a7Y1wCITPUaq9C1Y3AjQf8YMlkVa080O280ezz4LLPg+tw6RMOn15H3edYnENgsEdwytDt5cBzI+pFko5I4xII/wqsSHJqkh8DVgObRtyTJB1RxuKQUVXtTvIx4B+BRcDnqmrLG/iQB3zY6RCxz4PLPg+uw6VPOHx6HWmfqdrrUL0k6Qg0LoeMJEkjZiBIkoAjMBDG9ScykpyS5GtJnkiyJcnHu/7pJP+e5OG+XDAGvT6d5NHuZ7Jrxye5J8mTfX3ciHv82aE5ezjJD5J8YhzmM8nnkuxM8thQbdb5S3Jlv1+3Jjl/xH3+eZJvJ/lWki8neWfXJ5L899C8fnbEfc76Oo/ZfN4+1OPTSR7u+mjms6qOmAuDE9bfAd4F/BjwCHD6qPvq3pYC7+nldwD/xuBnPD4N/MGo+5vW69PAidNqfwas7+X1wDWj7nPa6/4fDL6QM/L5BN4PvAd4bF/z1++BR4CjgFP7/btohH3+MrC4l68Z6nNieL0xmM8ZX+dxm89p9/8F8MejnM8jbQ/h/38io6r+B5j6iYyRq6odVfVQL78CPMHgG9yHi4uADb28Abh4dK3s5TzgO1X1vVE3AlBV3wBenFaebf4uAm6rqteq6ilgG4P38Uj6rKqvVtXuvvkvDL4zNFKzzOdsxmo+pyQJ8BvAFw5FL7M50gJhpp/IGLv/6SaZAM4C7u/Sx3oX/XOjPhTTCvhqkgf750QATq6qHTAIN+CkkXW3t9Xs+R/auM0nzD5/4/ye/R3g7qHbpyb5ZpJ/SvK+UTU1ZKbXeVzn833A81X15FDtkM/nkRYI8/qJjFFK8nbgi8AnquoHwA3ATwPvBnYw2K0ctfdW1XsY/DrtFUneP+qGZtNfdPwI8A9dGsf5nMtYvmeTfArYDXy+SzuAn6yqs4DfA/4+yTGj6o/ZX+exnE/go+z5R8tI5vNIC4Sx/omMJG9lEAafr6ovAVTV81X1elX9CPgbDtHu7Vyq6rm+3gl8mUFPzydZCtDXO0fX4R4+BDxUVc/DeM5nm23+xu49m2QN8GHgN6sPePchmBd6+UEGx+Z/ZlQ9zvE6j+N8LgZ+Dbh9qjaq+TzSAmFsfyKjjyHeBDxRVZ8Zqi8dWu1Xgcemjz2UkvxEkndMLTM4yfgYg3lc06utAe4cTYd72eMvr3GbzyGzzd8mYHWSo5KcCqwAHhhBf8DgU3rAJ4GPVNV/DdWXZPDvmpDkXQz6/O5oupzzdR6r+WwfBL5dVdunCiObz0N9FnvUF+ACBp/g+Q7wqVH3M9TXLzLYdf0W8HBfLgBuBR7t+iZg6Yj7fBeDT2k8AmyZmkPgBGAz8GRfHz8Gc/rjwAvAsUO1kc8ng4DaAfwvg79Y1841f8Cn+v26FfjQiPvcxuAY/NR79LO97q/3++ER4CHgV0bc56yv8zjNZ9dvBn532rojmU9/ukKSBBx5h4wkSbMwECRJgIEgSWoGgiQJMBAkSc1AkCQBBoIkqf0fJxUqGvVsO1sAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"181\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAD4CAYAAADsKpHdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAATH0lEQVR4nO3df4xc13mf8ecbMmFkO3Qka6Wwu3SXiYk0FNHa1kJhYyBwwCRio8BUAQtYo4mIlgBbgWmdNkVKJkDUfwjIyA+1AiICjKWQSl3RhGJDRFy5FugERgFZyspWQlEMq22oSmsy4qZ2HTaFmZB9+8ccwuPl7HI5Q+1wtc8HGMyd955z75mLC373njszTFUhSdJ3DXsAkqQbg4EgSQIMBElSYyBIkgADQZLUrB72APp166231vj4+LCHIUnLyosvvviXVTXSa92yDYTx8XGmpqaGPQxJWlaS/M/51jllJEkCDARJUmMgSJIAA0GS1BgIkiTAQJAkNQaCJAkwECRJjYEgSQKW8TeVBzG+53ND2/drD90ztH1L0kK8QpAkAQaCJKkxECRJgIEgSWoMBEkSYCBIkhoDQZIELCIQkjye5FySl3us+7dJKsmtXbW9SaaTnEpyd1f9ziTH27pHkqTV1yT5dKs/n2T8Or03SdI1WMwVwkFg29xikvXATwGvd9U2AZPAHa3Po0lWtdX7gV3Axva4vM2dwDeq6n3Aw8An+nkjkqTBXDUQqupLwNd7rHoY+GWgumrbgcNVdaGqTgPTwF1J1gFrq+q5qirgCeDerj6H2vJTwNbLVw+SpKXT1z2EJB8BvlZVfzJn1SjwRtfrmVYbbctz69/Rp6ouAt8E3jPPfnclmUoyNTs728/QJUnzuOZASPIO4FeBX+u1uketFqgv1OfKYtWBqpqoqomRkZHFDFeStEj9XCH8ELAB+JMkrwFjwFeS/ACdv/zXd7UdA860+liPOt19kqwG3k3vKSpJ0lvomgOhqo5X1W1VNV5V43T+Qf9gVf0FcBSYbJ8c2kDn5vELVXUWOJ9kS7s/cD/wdNvkUWBHW/4o8MV2n0GStIQW87HTJ4HngB9OMpNk53xtq+oEcAR4Bfg8sLuqLrXVDwCfpHOj+X8Az7T6Y8B7kkwD/wbY0+d7kSQN4Kr/H0JVfewq68fnvN4H7OvRbgrY3KP+LeC+q41DkvTW8pvKkiTAQJAkNQaCJAkwECRJjYEgSQIMBElSYyBIkgADQZLUGAiSJMBAkCQ1BoIkCTAQJEmNgSBJAgwESVJjIEiSAANBktQYCJIkwECQJDUGgiQJWEQgJHk8ybkkL3fVfj3JnyX50ySfTfL9Xev2JplOcirJ3V31O5Mcb+seSZJWX5Pk063+fJLx6/sWJUmLsZgrhIPAtjm1Z4HNVfX3gf8O7AVIsgmYBO5ofR5Nsqr12Q/sAja2x+Vt7gS+UVXvAx4GPtHvm5Ek9e+qgVBVXwK+Pqf2haq62F5+GRhry9uBw1V1oapOA9PAXUnWAWur6rmqKuAJ4N6uPofa8lPA1stXD5KkpXM97iH8M+CZtjwKvNG1bqbVRtvy3Pp39Gkh803gPb12lGRXkqkkU7Ozs9dh6JKkywYKhCS/ClwEPnW51KNZLVBfqM+VxaoDVTVRVRMjIyPXOlxJ0gL6DoQkO4CfBf5JmwaCzl/+67uajQFnWn2sR/07+iRZDbybOVNUkqS3Xl+BkGQb8O+Aj1TV/+1adRSYbJ8c2kDn5vELVXUWOJ9kS7s/cD/wdFefHW35o8AXuwJGkrREVl+tQZIngQ8DtyaZAR6k86miNcCz7f7vl6vqX1TViSRHgFfoTCXtrqpLbVMP0PnE0k107jlcvu/wGPB7SabpXBlMXp+3Jkm6FlcNhKr6WI/yYwu03wfs61GfAjb3qH8LuO9q45AkvbX8prIkCTAQJEmNgSBJAgwESVJjIEiSAANBktQYCJIkwECQJDUGgiQJMBAkSY2BIEkCDARJUmMgSJIAA0GS1BgIkiTAQJAkNQaCJAkwECRJjYEgSQIWEQhJHk9yLsnLXbVbkjyb5NX2fHPXur1JppOcSnJ3V/3OJMfbukeSpNXXJPl0qz+fZPw6v0dJ0iIs5grhILBtTm0PcKyqNgLH2muSbAImgTtan0eTrGp99gO7gI3tcXmbO4FvVNX7gIeBT/T7ZiRJ/btqIFTVl4CvzylvBw615UPAvV31w1V1oapOA9PAXUnWAWur6rmqKuCJOX0ub+spYOvlqwdJ0tLp9x7C7VV1FqA939bqo8AbXe1mWm20Lc+tf0efqroIfBN4T6+dJtmVZCrJ1OzsbJ9DlyT1cr1vKvf6y74WqC/U58pi1YGqmqiqiZGRkT6HKEnqZXWf/d5Msq6qzrbpoHOtPgOs72o3Bpxp9bEe9e4+M0lWA+/myimqt43xPZ8byn5fe+ieoexX0vLR7xXCUWBHW94BPN1Vn2yfHNpA5+bxC21a6XySLe3+wP1z+lze1keBL7b7DJKkJXTVK4QkTwIfBm5NMgM8CDwEHEmyE3gduA+gqk4kOQK8AlwEdlfVpbapB+h8Yukm4Jn2AHgM+L0k03SuDCavyzuTJF2TqwZCVX1snlVb52m/D9jXoz4FbO5R/xYtUCRJw+M3lSVJgIEgSWoMBEkSYCBIkhoDQZIEGAiSpMZAkCQBBoIkqTEQJEmAgSBJagwESRJgIEiSGgNBkgQYCJKkxkCQJAEGgiSpMRAkSYCBIElqDARJEjBgICT510lOJHk5yZNJvjfJLUmeTfJqe765q/3eJNNJTiW5u6t+Z5Ljbd0jSTLIuCRJ167vQEgyCvwrYKKqNgOrgElgD3CsqjYCx9prkmxq6+8AtgGPJlnVNrcf2AVsbI9t/Y5LktSfQaeMVgM3JVkNvAM4A2wHDrX1h4B72/J24HBVXaiq08A0cFeSdcDaqnquqgp4oquPJGmJ9B0IVfU14DeA14GzwDer6gvA7VV1trU5C9zWuowCb3RtYqbVRtvy3PoVkuxKMpVkanZ2tt+hS5J6GGTK6GY6f/VvAP4O8M4kP7dQlx61WqB+ZbHqQFVNVNXEyMjItQ5ZkrSAQaaMfhI4XVWzVfW3wGeAHwPebNNAtOdzrf0MsL6r/xidKaaZtjy3LklaQoMEwuvAliTvaJ8K2gqcBI4CO1qbHcDTbfkoMJlkTZINdG4ev9Cmlc4n2dK2c39XH0nSElndb8eqej7JU8BXgIvAV4EDwLuAI0l20gmN+1r7E0mOAK+09rur6lLb3APAQeAm4Jn2kCQtob4DAaCqHgQenFO+QOdqoVf7fcC+HvUpYPMgY5EkDcZvKkuSAANBktQYCJIkwECQJDUGgiQJMBAkSY2BIEkCDARJUmMgSJIAA0GS1BgIkiTAQJAkNQaCJAkwECRJjYEgSQIMBElSYyBIkgADQZLUGAiSJGDAQEjy/UmeSvJnSU4m+YdJbknybJJX2/PNXe33JplOcirJ3V31O5Mcb+seSZJBxiVJunaDXiH8R+DzVfX3gH8AnAT2AMeqaiNwrL0mySZgErgD2AY8mmRV285+YBewsT22DTguSdI16jsQkqwFfhx4DKCq/qaq/jewHTjUmh0C7m3L24HDVXWhqk4D08BdSdYBa6vquaoq4ImuPpKkJTLIFcIPArPA7yb5apJPJnkncHtVnQVoz7e19qPAG139Z1pttC3PrV8hya4kU0mmZmdnBxi6JGmuQQJhNfBBYH9VfQD4a9r00Dx63ReoBepXFqsOVNVEVU2MjIxc63glSQsYJBBmgJmqer69fopOQLzZpoFoz+e62q/v6j8GnGn1sR51SdIS6jsQquovgDeS/HArbQVeAY4CO1ptB/B0Wz4KTCZZk2QDnZvHL7RppfNJtrRPF93f1UeStERWD9j/XwKfSvI9wJ8D/5ROyBxJshN4HbgPoKpOJDlCJzQuArur6lLbzgPAQeAm4Jn2kCQtoYECoapeAiZ6rNo6T/t9wL4e9Slg8yBjkSQNxm8qS5IAA0GS1BgIkiTAQJAkNQaCJAkwECRJjYEgSQIMBElSYyBIkgADQZLUGAiSJMBAkCQ1BoIkCTAQJEmNgSBJAgwESVJjIEiSAANBktQYCJIk4DoEQpJVSb6a5A/a61uSPJvk1fZ8c1fbvUmmk5xKcndX/c4kx9u6R5Jk0HFJkq7N9bhC+Dhwsuv1HuBYVW0EjrXXJNkETAJ3ANuAR5Osan32A7uAje2x7TqMS5J0DQYKhCRjwD3AJ7vK24FDbfkQcG9X/XBVXaiq08A0cFeSdcDaqnquqgp4oquPJGmJrB6w/38Afhn4vq7a7VV1FqCqzia5rdVHgS93tZtptb9ty3PrV0iyi86VBO9973sHHPrKMr7nc0Pb92sP3TO0fUtavL6vEJL8LHCuql5cbJcetVqgfmWx6kBVTVTVxMjIyCJ3K0lajEGuED4EfCTJzwDfC6xN8p+AN5Osa1cH64Bzrf0MsL6r/xhwptXHetQlSUuo7yuEqtpbVWNVNU7nZvEXq+rngKPAjtZsB/B0Wz4KTCZZk2QDnZvHL7TppfNJtrRPF93f1UeStEQGvYfQy0PAkSQ7gdeB+wCq6kSSI8ArwEVgd1Vdan0eAA4CNwHPtIckaQldl0Coqj8C/qgt/y9g6zzt9gH7etSngM3XYyySpP74TWVJEmAgSJIaA0GSBBgIkqTGQJAkAQaCJKkxECRJgIEgSWoMBEkSYCBIkhoDQZIEGAiSpMZAkCQBBoIkqTEQJEmAgSBJagwESRJgIEiSGgNBkgQMEAhJ1if5wyQnk5xI8vFWvyXJs0lebc83d/XZm2Q6yakkd3fV70xyvK17JEkGe1uSpGs1yBXCReCXqupHgC3A7iSbgD3AsaraCBxrr2nrJoE7gG3Ao0lWtW3tB3YBG9tj2wDjkiT1oe9AqKqzVfWVtnweOAmMAtuBQ63ZIeDetrwdOFxVF6rqNDAN3JVkHbC2qp6rqgKe6OojSVoi1+UeQpJx4APA88DtVXUWOqEB3NaajQJvdHWbabXRtjy33ms/u5JMJZmanZ29HkOXJDUDB0KSdwG/D/xiVf3VQk171GqB+pXFqgNVNVFVEyMjI9c+WEnSvAYKhCTfTScMPlVVn2nlN9s0EO35XKvPAOu7uo8BZ1p9rEddkrSEBvmUUYDHgJNV9Vtdq44CO9ryDuDprvpkkjVJNtC5efxCm1Y6n2RL2+b9XX0kSUtk9QB9PwT8PHA8yUut9ivAQ8CRJDuB14H7AKrqRJIjwCt0PqG0u6outX4PAAeBm4Bn2kOStIT6DoSq+m/0nv8H2DpPn33Avh71KWBzv2ORJA3ObypLkgADQZLUGAiSJMBAkCQ1g3zKSFqU8T2fG8p+X3vonqHsV1quvEKQJAEGgiSpMRAkSYCBIElqDARJEmAgSJIaA0GSBBgIkqTGQJAkAQaCJKkxECRJgIEgSWoMBEkS4K+d6m3MX1mVrs0Nc4WQZFuSU0mmk+wZ9ngkaaW5IQIhySrgt4F/BGwCPpZk03BHJUkry40yZXQXMF1Vfw6Q5DCwHXhlqKOS+jCsqSpwukqDuVECYRR4o+v1DPCjcxsl2QXsai//T5JTfe7vVuAv++z7duex6W1ZHJd8Yii7XRbHZghu1OPyd+dbcaMEQnrU6opC1QHgwMA7S6aqamLQ7bwdeWx687jMz2PT23I8LjfEPQQ6VwTru16PAWeGNBZJWpFulED4Y2Bjkg1JvgeYBI4OeUyStKLcEFNGVXUxyS8A/xVYBTxeVSfewl0OPO30Nuax6c3jMj+PTW/L7rik6oqpeknSCnSjTBlJkobMQJAkASswEPyJjG9L8lqS40leSjLVarckeTbJq+355mGPcykkeTzJuSQvd9XmPRZJ9rZz6FSSu4cz6rfePMfl3yf5WjtvXkryM13rVspxWZ/kD5OcTHIiycdbfVmfMysqEPyJjJ5+oqre3/V56T3AsaraCBxrr1eCg8C2ObWex6KdM5PAHa3Po+3cejs6yJXHBeDhdt68v6r+C6y443IR+KWq+hFgC7C7vf9lfc6sqECg6ycyqupvgMs/kaFv2w4casuHgHuHN5SlU1VfAr4+pzzfsdgOHK6qC1V1Gpimc2697cxzXOazko7L2ar6Sls+D5yk84sLy/qcWWmB0OsnMkaHNJYbQQFfSPJi+1kQgNur6ix0TnrgtqGNbvjmOxaeR/ALSf60TSldnhZZkcclyTjwAeB5lvk5s9ICYVE/kbGCfKiqPkhnCm13kh8f9oCWiZV+Hu0Hfgh4P3AW+M1WX3HHJcm7gN8HfrGq/mqhpj1qN9yxWWmB4E9kdKmqM+35HPBZOpewbyZZB9Cezw1vhEM337FY0edRVb1ZVZeq6v8Bv8O3pz5W1HFJ8t10wuBTVfWZVl7W58xKCwR/IqNJ8s4k33d5Gfhp4GU6x2NHa7YDeHo4I7whzHcsjgKTSdYk2QBsBF4YwviG4vI/eM0/pnPewAo6LkkCPAacrKrf6lq1rM+ZG+KnK5bKEH4i40Z2O/DZznnNauA/V9Xnk/wxcCTJTuB14L4hjnHJJHkS+DBwa5IZ4EHgIXoci6o6keQInf+v4yKwu6ouDWXgb7F5jsuHk7yfzpTHa8A/h5V1XIAPAT8PHE/yUqv9Csv8nPGnKyRJwMqbMpIkzcNAkCQBBoIkqTEQJEmAgSBJagwESRJgIEiSmv8PO4PnJtPRJ2IAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"216\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAD4CAYAAAAAczaOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAATFUlEQVR4nO3df6ye533X8fcHe/X6Y9ESfBI828Xu5I050aDpIQQKVSGweM1UZ39EcqRSCyJZRN7oEGPY9I/sH0veGIVVIpFMG+pCiWV1HbGoOhqZbdGkLOakS5o4mRd3DvFpvPiUCBZAcpf0yx/PFfZw8hwfn+c5Pic+1/slHd33872v+7muy7f1Ofe5nx93qgpJUh/+zGoPQJK0cgx9SeqIoS9JHTH0Jakjhr4kdWT9ag9gMRs3bqxt27at9jAk6Zry1FNPfaeqpubX3/Ghv23bNmZmZlZ7GJJ0TUny30bVvbwjSR0x9CWpI4a+JHXE0Jekjhj6ktQRQ1+SOmLoS1JHDH1J6oihL0kdecd/IncS2w58dVX6fenwXavSryQtxjN9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1BFDX5I6YuhLUkcMfUnqyKKhn+ThJBeTPDev/rNJziQ5neSXh+oHk5xt2+4cqn8oybNt22eTZHmnIklazJWc6X8B2DVcSPI3gd3Aj1fVzcCvtPpOYA9wc9vnwSTr2m4PAfuAHe3n/3tOSdLVt2joV9XjwGvzyvcDh6vqUmtzsdV3A8eq6lJVnQPOArcl2QRcV1VPVFUBXwTuXqY5SJKu0LjX9H8E+BtJnkzy20n+cqtvBs4PtZtttc1tfX59pCT7kswkmZmbmxtziJKk+cYN/fXA9cDtwD8Bjrdr9KOu09dl6iNV1ZGqmq6q6ampqTGHKEmab9zQnwW+UgOngO8BG1t961C7LcArrb5lRF2StILGDf3/CPwtgCQ/ArwL+A5wAtiTZEOS7QxesD1VVReA15Pc3v4i+CTw6KSDlyQtzaI3UUnyCPBRYGOSWeAB4GHg4fY2zu8Ce9sLtKeTHAeeB94A9lfVm+2p7mfwTqB3A19rP5KkFbRo6FfVvQts+sQC7Q8Bh0bUZ4BbljQ6SdKy8hO5ktQRQ1+SOmLoS1JHDH1J6oihL0kdMfQlqSOGviR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdWTR0E/ycJKL7YYp87f9fJJKsnGodjDJ2SRnktw5VP9Qkmfbts+2O2hJklbQlZzpfwHYNb+YZCvwd4CXh2o7gT3AzW2fB5Osa5sfAvYxuIXijlHPKUm6uhYN/ap6HHhtxKZ/CfwCUEO13cCxqrpUVeeAs8BtSTYB11XVE+22il8E7p508JKkpRnrmn6SjwPfrqpn5m3aDJwfejzbapvb+vy6JGkFLXqP3PmSvAf4NPATozaPqNVl6gv1sY/BpSDe//73L3WIkqQFjHOm/8PAduCZJC8BW4BvJPlzDM7gtw613QK80upbRtRHqqojVTVdVdNTU1NjDFGSNMqSQ7+qnq2qG6tqW1VtYxDot1bVHwEngD1JNiTZzuAF21NVdQF4Pcnt7V07nwQeXb5pSJKuxJW8ZfMR4AngR5PMJrlvobZVdRo4DjwP/Aawv6rebJvvBz7H4MXdbwFfm3DskqQlWvSaflXdu8j2bfMeHwIOjWg3A9yyxPFJkpaRn8iVpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1BFDX5I6YuhLUkcMfUnqiKEvSR0x9CWpI4a+JHXkSu6c9XCSi0meG6r98yS/n+SbSX49yQ8ObTuY5GySM0nuHKp/KMmzbdtn220TJUkr6ErO9L8A7JpXewy4pap+HPgD4CBAkp3AHuDmts+DSda1fR4C9jG4b+6OEc8pSbrKFg39qnoceG1e7etV9UZ7+LvAlra+GzhWVZeq6hyD++HelmQTcF1VPVFVBXwRuHuZ5iBJukLLcU3/7/OnNznfDJwf2jbbapvb+vz6SEn2JZlJMjM3N7cMQ5QkwYShn+TTwBvAl94qjWhWl6mPVFVHqmq6qqanpqYmGaIkacj6cXdMshf4KeCOdskGBmfwW4eabQFeafUtI+qSpBU01pl+kl3APwU+XlX/Z2jTCWBPkg1JtjN4wfZUVV0AXk9ye3vXzieBRyccuyRpiRY900/yCPBRYGOSWeABBu/W2QA81t55+btV9Q+q6nSS48DzDC777K+qN9tT3c/gnUDvZvAawNeQJK2oRUO/qu4dUf78ZdofAg6NqM8AtyxpdJKkZTX2NX0tbNuBr65a3y8dvmvV+pb0zufXMEhSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1BFDX5I6YuhLUkcMfUnqiKEvSR0x9CWpI4a+JHXE0Jekjiwa+kkeTnIxyXNDtRuSPJbkxba8fmjbwSRnk5xJcudQ/UNJnm3bPtvuoCVJWkFXcqb/BWDXvNoB4GRV7QBOtsck2QnsAW5u+zyYZF3b5yFgH4NbKO4Y8ZySpKts0dCvqseB1+aVdwNH2/pR4O6h+rGqulRV54CzwG1JNgHXVdUT7SbqXxzaR5K0Qsa9pn9Tu9k5bXljq28Gzg+1m221zW19fn2kJPuSzCSZmZubG3OIkqT5lvuF3FHX6esy9ZGq6khVTVfV9NTU1LINTpJ6N27ov9ou2dCWF1t9Ftg61G4L8EqrbxlRlyStoHFD/wSwt63vBR4dqu9JsiHJdgYv2J5ql4BeT3J7e9fOJ4f2kSStkPWLNUjyCPBRYGOSWeAB4DBwPMl9wMvAPQBVdTrJceB54A1gf1W92Z7qfgbvBHo38LX2I0laQYuGflXdu8CmOxZofwg4NKI+A9yypNFJkpaVn8iVpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1BFDX5I6YuhLUkcMfUnqiKEvSR0x9CWpIxOFfpJ/lOR0kueSPJLk+5PckOSxJC+25fVD7Q8mOZvkTJI7Jx++JGkpxg79JJuBfwhMV9UtwDpgD3AAOFlVO4CT7TFJdrbtNwO7gAeTrJts+JKkpZj08s564N1J1gPvYXCz893A0bb9KHB3W98NHKuqS1V1DjgL3DZh/5KkJRg79Kvq28CvMLhH7gXgf1bV14Gb2o3Qacsb2y6bgfNDTzHbam+TZF+SmSQzc3Nz4w5RkjTPJJd3rmdw9r4d+CHgvUk+cbldRtRqVMOqOlJV01U1PTU1Ne4QJUnzTHJ5528D56pqrqr+BPgK8NeAV5NsAmjLi639LLB1aP8tDC4HSZJWyCSh/zJwe5L3JAlwB/ACcALY29rsBR5t6yeAPUk2JNkO7ABOTdC/JGmJ1o+7Y1U9meTLwDeAN4DfA44A7wOOJ7mPwS+Ge1r700mOA8+39vur6s0Jxy9JWoKxQx+gqh4AHphXvsTgrH9U+0PAoUn6lCSNz0/kSlJHDH1J6oihL0kdMfQlqSOGviR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1JGJQj/JDyb5cpLfT/JCkr+a5IYkjyV5sS2vH2p/MMnZJGeS3Dn58CVJSzHpmf6vAr9RVX8B+IsMbpd4ADhZVTuAk+0xSXYCe4CbgV3Ag0nWTdi/JGkJxg79JNcBHwE+D1BV362q/wHsBo62ZkeBu9v6buBYVV2qqnPAWeC2cfuXJC3dJGf6HwDmgH+b5PeSfC7Je4GbquoCQFve2NpvBs4P7T/bam+TZF+SmSQzc3NzEwxRkjRsktBfD9wKPFRVHwT+N+1SzgIyolajGlbVkaqarqrpqampCYYoSRo2SejPArNV9WR7/GUGvwReTbIJoC0vDrXfOrT/FuCVCfqXJC3R2KFfVX8EnE/yo610B/A8cALY22p7gUfb+glgT5INSbYDO4BT4/YvSVq69RPu/7PAl5K8C/hD4O8x+EVyPMl9wMvAPQBVdTrJcQa/GN4A9lfVmxP2L0lagolCv6qeBqZHbLpjgfaHgEOT9ClJGp+fyJWkjhj6ktQRQ1+SOmLoS1JHJn33jt5hth346qr0+9Lhu1alX0lL45m+JHXE0Jekjhj6ktQRQ1+SOmLoS1JHDH1J6oihL0kdMfQlqSOGviR1xNCXpI5MHPpJ1rUbo/+n9viGJI8lebEtrx9qezDJ2SRnktw5ad+SpKVZjjP9TwEvDD0+AJysqh3AyfaYJDuBPcDNwC7gwSTrlqF/SdIVmij0k2wB7gI+N1TeDRxt60eBu4fqx6rqUlWdA84Ct03SvyRpaSY90/9XwC8A3xuq3VRVFwDa8sZW3wycH2o322pvk2RfkpkkM3NzcxMOUZL0lrFDP8lPARer6qkr3WVErUY1rKojVTVdVdNTU1PjDlGSNM8k36f/YeDjST4GfD9wXZJ/D7yaZFNVXUiyCbjY2s8CW4f23wK8MkH/kqQlGvtMv6oOVtWWqtrG4AXa/1JVnwBOAHtbs73Ao239BLAnyYYk24EdwKmxRy5JWrKrceesw8DxJPcBLwP3AFTV6STHgeeBN4D9VfXmVehfkrSAZQn9qvot4Lfa+n8H7lig3SHg0HL0KUlaOj+RK0kdMfQlqSOGviR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1BFDX5I6YuhLUkcmuUfu1iS/meSFJKeTfKrVb0jyWJIX2/L6oX0OJjmb5EySO5djApKkKzfJmf4bwD+uqh8Dbgf2J9kJHABOVtUO4GR7TNu2B7gZ2AU8mGTdJIOXJC3NJPfIvVBV32jrrwMvAJuB3cDR1uwocHdb3w0cq6pLVXUOOAvcNm7/kqSlW5Zr+km2AR8EngRuqqoLMPjFANzYmm0Gzg/tNttqkqQVMnHoJ3kf8GvAz1XVH1+u6YhaLfCc+5LMJJmZm5ubdIiSpGai0E/yfQwC/0tV9ZVWfjXJprZ9E3Cx1WeBrUO7bwFeGfW8VXWkqqaranpqamqSIUqShqwfd8ckAT4PvFBVnxnadALYCxxuy0eH6v8hyWeAHwJ2AKfG7V/vLNsOfHXV+n7p8F2r1rd0rRk79IEPA38XeDbJ0632zxiE/fEk9wEvA/cAVNXpJMeB5xm882d/Vb05Qf+SpCUaO/Sr6ncYfZ0e4I4F9jkEHBq3T0nSZPxEriR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOTPLdO9I7wmp92Ztf9KZrkWf6ktQRQ1+SOmLoS1JHDH1J6oihL0kdWfF37yTZBfwqsA74XFUdXukxSMvBW0TqWrSiZ/pJ1gH/GvhJYCdwb5KdKzkGSerZSp/p3wacrao/BEhyDNjN4L65kq7Qav6VsRr8y2b5rHTobwbODz2eBf7K/EZJ9gH72sP/leTMmP1tBL4z5r7Xqt7m3Nt8ob85b8wvdTVfWJ5j/OdHFVc69EfdSL3eVqg6AhyZuLNkpqqmJ32ea0lvc+5tvtDfnHubL1zdOa/0u3dmga1Dj7cAr6zwGCSpWysd+v8V2JFke5J3AXuAEys8Bknq1ope3qmqN5L8DPCfGbxl8+GqOn0Vu5z4EtE1qLc59zZf6G/Ovc0XruKcU/W2S+qSpDXKT+RKUkcMfUnqyJoM/SS7kpxJcjbJgdUez9WS5KUkzyZ5OslMq92Q5LEkL7bl9as9zkkkeTjJxSTPDdUWnGOSg+24n0ly5+qMenwLzPcXk3y7Heenk3xsaNs1PV+AJFuT/GaSF5KcTvKpVl+Tx/ky812Z41xVa+qHwQvE3wI+ALwLeAbYudrjukpzfQnYOK/2y8CBtn4A+KXVHueEc/wIcCvw3GJzZPDVHs8AG4Dt7f/ButWewzLM9xeBnx/R9pqfb5vHJuDWtv4DwB+0ua3J43yZ+a7IcV6LZ/r/76sequq7wFtf9dCL3cDRtn4UuHv1hjK5qnoceG1eeaE57gaOVdWlqjoHnGXw/+GascB8F3LNzxegqi5U1Tfa+uvACww+vb8mj/Nl5ruQZZ3vWgz9UV/1cLl/0GtZAV9P8lT76gqAm6rqAgz+cwE3rtrorp6F5riWj/3PJPlmu/zz1mWONTffJNuADwJP0sFxnjdfWIHjvBZD/4q+6mGN+HBV3crgW0v3J/nIag9ola3VY/8Q8MPAXwIuAP+i1dfUfJO8D/g14Oeq6o8v13RE7Zqb94j5rshxXouh381XPVTVK215Efh1Bn/yvZpkE0BbXly9EV41C81xTR77qnq1qt6squ8B/4Y//dN+zcw3yfcxCMAvVdVXWnnNHudR812p47wWQ7+Lr3pI8t4kP/DWOvATwHMM5rq3NdsLPLo6I7yqFprjCWBPkg1JtgM7gFOrML5l9VbwNT/N4DjDGplvkgCfB16oqs8MbVqTx3mh+a7YcV7tV7Kv0qvjH2Pwivi3gE+v9niu0hw/wOAV/WeA02/NE/izwEngxba8YbXHOuE8H2Hwp+6fMDjjue9ycwQ+3Y77GeAnV3v8yzTffwc8C3yzBcCmtTLfNoe/zuByxTeBp9vPx9bqcb7MfFfkOPs1DJLUkbV4eUeStABDX5I6YuhLUkcMfUnqiKEvSR0x9CWpI4a+JHXk/wKcuwLFfq6utQAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"251\n"
]
}
],
"source": [
"new_error = reconstruct(err, A)"
"\n",
"plt.hist(first)\n",
"plt.show()\n",
"print(np.max(first))\n",
"plt.hist(second)\n",
"plt.show()\n",
"print(np.max(second))\n",
"plt.hist(third)\n",
"plt.show()\n",
"print(np.max(third))\n",
"plt.hist(fourth)\n",
"plt.show()\n",
"print(np.max(fourth))"
]
},
{
"cell_type": "code",
"execution_count": 420,
"id": "06ccaf8e",
"execution_count": 9,
"id": "139938e3",
"metadata": {},
"outputs": [],
"source": [
"class NodeTree(object):\n",
" def __init__(self, left=None, right=None):\n",
" self.left = left\n",
" self.right = right\n",
"\n",
" def children(self):\n",
" return self.left, self.right\n",
"\n",
" def __str__(self):\n",
" return self.left, self.right\n",
"\n",
"\n",
"def huffman_code_tree(node, binString=''):\n",
" '''\n",
" Function to find Huffman Code\n",
" '''\n",
" if type(node) is str:\n",
" return {node: binString}\n",
" (l, r) = node.children()\n",
" d = dict()\n",
" d.update(huffman_code_tree(l, binString + '0'))\n",
" d.update(huffman_code_tree(r, binString + '1'))\n",
" return d\n",
"\n",
"\n",
"def make_tree(nodes):\n",
" '''\n",
" Function to make tree\n",
" :param nodes: Nodes\n",
" :return: Root of the tree\n",
" '''\n",
" while len(nodes) > 1:\n",
" (key1, c1) = nodes[-1]\n",
" (key2, c2) = nodes[-2]\n",
" nodes = nodes[:-2]\n",
" node = NodeTree(key1, key2)\n",
" nodes.append((node, c1 + c2))\n",
" nodes = sorted(nodes, key=lambda x: x[1], reverse=True)\n",
" return nodes[0][0]"
]
},
{
"cell_type": "code",
"execution_count": 263,
"id": "e477d0c8",
"metadata": {},
"outputs": [],
"source": [
"def enc_experiment(images, plot=True):\n",
" origin, predict, diff, error, A = plot_hist(images, 2)\n",
" image = Image.open(images[0]) #Open the image and read it as an Image object\n",
" image = np.array(image)[1:,:] #Convert to an array, leaving out the first row because the first row is just housekeeping data\n",
" image = image.astype(int)\n",
" new_error = np.copy(image)\n",
" #new_error[1:-1,1:-1] = np.reshape(error[1:-1,1:-1],(510, 638))\n",
" new_error[1:-1, 1:-1] = error[1:-1, 1:-1]\n",
" keep = new_error[0,0]\n",
" new_error[0,:] = new_error[0,:] - keep\n",
" new_error[-1,:] = new_error[-1,:] - keep\n",
" new_error[1:-1,0] = new_error[1:-1,0] - keep\n",
" new_error[1:-1,-1] = new_error[1:-1,-1] - keep\n",
" new_error[0,0] = keep\n",
" new_error = np.ravel(new_error)\n",
" if plot:\n",
" plt.hist(new_error[1:],bins=100)\n",
" plt.show()\n",
" \n",
" #ab_error = np.abs(new_error)\n",
" #string = [str(i) for i in ab_error]\n",
" string = [str(i) for i in new_error]\n",
" #string = [str(i) for i in np.arange(0,5)] + [str(i) for i in np.arange(0,5)] + [str(i) for i in np.arange(0,2)]*2\n",
" freq = dict(Counter(string))\n",
" freq = sorted(freq.items(), key=lambda x: x[1], reverse=True)\n",
" \n",
" node = make_tree(freq)\n",
" encoding_dict = huffman_code_tree(node)\n",
" #encoded = [\"1\"+encoding[str(-i)] if i < 0 else \"0\"+encoding[str(i)] for i in error]\n",
" #print(time.time()-start)\n",
" encoded = new_error.reshape((512,640)).copy().astype(str)\n",
" \n",
" #encoded = np.zeros_like(new_error.reshape((512,640))).astype(str)\n",
" print(encoded)\n",
" for i in range(encoded.shape[0]):\n",
" for j in range(encoded.shape[1]):\n",
" if i == 0 and j == 0:\n",
" encoded[i][j] = encoded[i][j]\n",
" else:\n",
" #print(encoding_dict[encoded[i][j]])\n",
" encoded[i][j] = encoding_dict[encoded[i][j]]\n",
" #print(encoded[i][j])\n",
" \n",
" return encoding_dict, encoded, new_error.reshape((512,640))\n",
" #print(encoding)"
]
},
{
"cell_type": "code",
"execution_count": 264,
"id": "fd8e96a5",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"518"
]
},
"execution_count": 420,
"metadata": {},
"output_type": "execute_result"
"name": "stdout",
"output_type": "stream",
"text": [
"[['22541' '-10' '14' ... '32' '48' '33']\n",
" ['7' '67' '-21' ... '-1' '-1' '77']\n",
" ['7' '-15' '-3' ... '10' '-45' '58']\n",
" ...\n",
" ['49' '82' '-2' ... '-64' '5' '151']\n",
" ['27' '-33' '18' ... '47' '-16' '208']\n",
" ['17' '0' '-5' ... '138' '207' '226']]\n"
]
}
],
"source": [
"encode_dict, encoding, error = enc_experiment(images, plot=False)"
]
},
{
"cell_type": "code",
"execution_count": 265,
"id": "7ebd8dd8",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"11110111110\n",
"92\n",
"11110111110001\n"
]
}
],
"source": [
"e = np.round(err, 1)\n",
"len(np.unique(e[1:-1, 1:-1]))"
"print(encoding[0][100])\n",
"print(error[0][100])\n",
"print(encode_dict['92'])"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "a6a579b7",
"metadata": {},
"outputs": [],
"source": [
"def decoder(A, encoded_matrix, encoding_dict):\n",
" \"\"\"\n",
" Function that accecpts the prediction matrix A for the linear system,\n",
" the encoded matrix of error values, and the encoding dicitonary.\n",
" \"\"\"\n",
" error_matrix = encoded_matrix.copy()\n",
" for i in range(error_matrix.shape[0]):\n",
" for j in range(error_matrix.shape[1]):\n",
" if i == 0 and j == 0:\n",
" error_matrix[i][j] = encoded_matrix[i][j]\n",
" else:\n",
" error_matrix[i][j] = encoding_dict."
]
}
],
......
......@@ -932,7 +932,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.1"
"version": "3.8.11"
}
},
"nbformat": 4,
......
......@@ -2,7 +2,7 @@
"cells": [
{
"cell_type": "code",
"execution_count": 626,
"execution_count": 33,
"id": "dbef8759",
"metadata": {
"id": "dbef8759"
......@@ -21,12 +21,13 @@
"from numpy import linalg as la\n",
"from scipy.stats import gaussian_kde\n",
"import seaborn as sns\n",
"import pywt"
"import pywt\n",
"from collections import Counter"
]
},
{
"cell_type": "code",
"execution_count": 654,
"execution_count": 2,
"id": "9ed20f84",
"metadata": {
"id": "9ed20f84"
......@@ -99,7 +100,7 @@
},
{
"cell_type": "code",
"execution_count": 647,
"execution_count": 3,
"id": "ba2881d9",
"metadata": {},
"outputs": [],
......@@ -111,7 +112,7 @@
},
{
"cell_type": "code",
"execution_count": 655,
"execution_count": 4,
"id": "11e95c34",
"metadata": {},
"outputs": [],
......@@ -121,7 +122,7 @@
},
{
"cell_type": "code",
"execution_count": 537,
"execution_count": 5,
"id": "434e4d2f",
"metadata": {},
"outputs": [],
......@@ -159,30 +160,18 @@
},
{
"cell_type": "code",
"execution_count": 606,
"id": "bf427edd",
"execution_count": 6,
"id": "4f4a5a35",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"-19\n",
"[22286.]\n",
"22285.0\n",
"[22266.]\n",
"22266.999999999996\n"
]
}
],
"outputs": [],
"source": [
"new_error = reconstruct(err, A)"
]
},
{
"cell_type": "code",
"execution_count": 669,
"id": "5edcf208",
"execution_count": 7,
"id": "6d95ffce",
"metadata": {},
"outputs": [],
"source": [
......@@ -205,8 +194,8 @@
},
{
"cell_type": "code",
"execution_count": 673,
"id": "953375c5",
"execution_count": 8,
"id": "1c848109",
"metadata": {},
"outputs": [
{
......@@ -304,22 +293,172 @@
},
{
"cell_type": "code",
"execution_count": 679,
"id": "90387cb9",
"execution_count": 9,
"id": "139938e3",
"metadata": {},
"outputs": [],
"source": [
"class NodeTree(object):\n",
" def __init__(self, left=None, right=None):\n",
" self.left = left\n",
" self.right = right\n",
"\n",
" def children(self):\n",
" return self.left, self.right\n",
"\n",
" def __str__(self):\n",
" return self.left, self.right\n",
"\n",
"\n",
"def huffman_code_tree(node, binString=''):\n",
" '''\n",
" Function to find Huffman Code\n",
" '''\n",
" if type(node) is str:\n",
" return {node: binString}\n",
" (l, r) = node.children()\n",
" d = dict()\n",
" d.update(huffman_code_tree(l, binString + '0'))\n",
" d.update(huffman_code_tree(r, binString + '1'))\n",
" return d\n",
"\n",
"\n",
"def make_tree(nodes):\n",
" '''\n",
" Function to make tree\n",
" :param nodes: Nodes\n",
" :return: Root of the tree\n",
" '''\n",
" while len(nodes) > 1:\n",
" (key1, c1) = nodes[-1]\n",
" (key2, c2) = nodes[-2]\n",
" nodes = nodes[:-2]\n",
" node = NodeTree(key1, key2)\n",
" nodes.append((node, c1 + c2))\n",
" nodes = sorted(nodes, key=lambda x: x[1], reverse=True)\n",
" return nodes[0][0]"
]
},
{
"cell_type": "code",
"execution_count": 263,
"id": "e477d0c8",
"metadata": {},
"outputs": [],
"source": [
"def enc_experiment(images, plot=True):\n",
" origin, predict, diff, error, A = plot_hist(images, 2)\n",
" image = Image.open(images[0]) #Open the image and read it as an Image object\n",
" image = np.array(image)[1:,:] #Convert to an array, leaving out the first row because the first row is just housekeeping data\n",
" image = image.astype(int)\n",
" new_error = np.copy(image)\n",
" #new_error[1:-1,1:-1] = np.reshape(error[1:-1,1:-1],(510, 638))\n",
" new_error[1:-1, 1:-1] = error[1:-1, 1:-1]\n",
" keep = new_error[0,0]\n",
" new_error[0,:] = new_error[0,:] - keep\n",
" new_error[-1,:] = new_error[-1,:] - keep\n",
" new_error[1:-1,0] = new_error[1:-1,0] - keep\n",
" new_error[1:-1,-1] = new_error[1:-1,-1] - keep\n",
" new_error[0,0] = keep\n",
" new_error = np.ravel(new_error)\n",
" if plot:\n",
" plt.hist(new_error[1:],bins=100)\n",
" plt.show()\n",
" \n",
" #ab_error = np.abs(new_error)\n",
" #string = [str(i) for i in ab_error]\n",
" string = [str(i) for i in new_error]\n",
" #string = [str(i) for i in np.arange(0,5)] + [str(i) for i in np.arange(0,5)] + [str(i) for i in np.arange(0,2)]*2\n",
" freq = dict(Counter(string))\n",
" freq = sorted(freq.items(), key=lambda x: x[1], reverse=True)\n",
" \n",
" node = make_tree(freq)\n",
" encoding_dict = huffman_code_tree(node)\n",
" #encoded = [\"1\"+encoding[str(-i)] if i < 0 else \"0\"+encoding[str(i)] for i in error]\n",
" #print(time.time()-start)\n",
" encoded = new_error.reshape((512,640)).copy().astype(str)\n",
" \n",
" #encoded = np.zeros_like(new_error.reshape((512,640))).astype(str)\n",
" print(encoded)\n",
" for i in range(encoded.shape[0]):\n",
" for j in range(encoded.shape[1]):\n",
" if i == 0 and j == 0:\n",
" encoded[i][j] = encoded[i][j]\n",
" else:\n",
" #print(encoding_dict[encoded[i][j]])\n",
" encoded[i][j] = encoding_dict[encoded[i][j]]\n",
" #print(encoded[i][j])\n",
" \n",
" return encoding_dict, encoded, new_error.reshape((512,640))\n",
" #print(encoding)"
]
},
{
"cell_type": "code",
"execution_count": 264,
"id": "fd8e96a5",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"15"
]
},
"execution_count": 679,
"metadata": {},
"output_type": "execute_result"
"name": "stdout",
"output_type": "stream",
"text": [
"[['22541' '-10' '14' ... '32' '48' '33']\n",
" ['7' '67' '-21' ... '-1' '-1' '77']\n",
" ['7' '-15' '-3' ... '10' '-45' '58']\n",
" ...\n",
" ['49' '82' '-2' ... '-64' '5' '151']\n",
" ['27' '-33' '18' ... '47' '-16' '208']\n",
" ['17' '0' '-5' ... '138' '207' '226']]\n"
]
}
],
"source": []
"source": [
"encode_dict, encoding, error = enc_experiment(images, plot=False)"
]
},
{
"cell_type": "code",
"execution_count": 265,
"id": "7ebd8dd8",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"11110111110\n",
"92\n",
"11110111110001\n"
]
}
],
"source": [
"print(encoding[0][100])\n",
"print(error[0][100])\n",
"print(encode_dict['92'])"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "a6a579b7",
"metadata": {},
"outputs": [],
"source": [
"def decoder(A, encoded_matrix, encoding_dict):\n",
" \"\"\"\n",
" Function that accecpts the prediction matrix A for the linear system,\n",
" the encoded matrix of error values, and the encoding dicitonary.\n",
" \"\"\"\n",
" error_matrix = encoded_matrix.copy()\n",
" for i in range(error_matrix.shape[0]):\n",
" for j in range(error_matrix.shape[1]):\n",
" if i == 0 and j == 0:\n",
" error_matrix[i][j] = encoded_matrix[i][j]\n",
" else:\n",
" error_matrix[i][j] = encoding_dict."
]
}
],
"metadata": {
......@@ -343,7 +482,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.1"
"version": "3.8.11"
}
},
"nbformat": 4,
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment