diff --git a/coherency_tests.ipynb b/coherency_tests.ipynb
index b2251501f7f3b5d7057b8657927987de41af761b..df9c06d9dfbcc11f29dc45833eb4a222f3f18a13 100644
--- a/coherency_tests.ipynb
+++ b/coherency_tests.ipynb
@@ -22,7 +22,7 @@
     {
      "data": {
       "text/plain": [
-       "(1842842, 5)"
+       "(921421, 5)"
       ]
      },
      "execution_count": 2,
@@ -65,17 +65,17 @@
      "output_type": "stream",
      "text": [
       "Building the KDtree v2:\n",
-      "\tTotal time: 1.741s\n",
+      "\tTotal time: 1.607s\n",
       "\n",
       "knn search:\n",
-      "\tTotal time: 18.676s\n",
+      "\tTotal time: 8.044s\n",
       "\n",
       "ID estimation:\n",
-      "\tID value: 3.920865\n",
-      "\tTotal time: 0.391s\n",
+      "\tID value: 3.862658\n",
+      "\tTotal time: 0.180s\n",
       "\n",
       "Density and k* estimation:\n",
-      "\tTotal time: 5.015s\n",
+      "\tTotal time: 2.085s\n",
       "\n"
      ]
     }
@@ -83,7 +83,7 @@
    "source": [
     "data.compute_distances(299)\n",
     "data.compute_id_2NN()\n",
-    "#data.id = 4\n",
+    "data.id = 4\n",
     "data.compute_density_kstarNN()"
    ]
   },
@@ -97,7 +97,7 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "3.920865231328582\n"
+      "4\n"
      ]
     }
    ],
@@ -116,7 +116,7 @@
     {
      "data": {
       "text/plain": [
-       "0.0001499206205206319"
+       "8.791685787264807e-08"
       ]
      },
      "execution_count": 6,
@@ -137,8 +137,8 @@
     {
      "data": {
       "text/plain": [
-       "array([-3.64725383, -4.53602916, -4.30573383, ...,  7.18148125,\n",
-       "       -6.25416517, -3.54944958])"
+       "array([-3.06640128, -3.86980622, -3.66751667, ..., -0.64454652,\n",
+       "       -3.35926525, -7.39476774])"
       ]
      },
      "execution_count": 7,
@@ -159,8 +159,8 @@
     {
      "data": {
       "text/plain": [
-       "array([-3.64714384, -4.53592253, -4.30561972, ...,  7.18179655,\n",
-       "       -6.25407934, -3.54934192])"
+       "array([-3.06640124, -3.86980629, -3.66751671, ..., -0.64454651,\n",
+       "       -3.35926533, -7.39476776])"
       ]
      },
      "execution_count": 8,
@@ -182,13 +182,13 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "-0.5825421810150146 -1.2760780561214355\n"
+      "-8.463785171508789 -8.463785648321064 max diff 4.768122749965187e-07\n"
      ]
     },
     {
      "data": {
       "text/plain": [
-       "1789174"
+       "674003"
       ]
      },
      "execution_count": 9,
@@ -199,20 +199,20 @@
    "source": [
     "\n",
     "a = np.argmax(np.abs(den_comp - den_gt))\n",
-    "print(den_comp[a], den_gt[a])\n",
+    "print(den_comp[a], den_gt[a], \"max diff\", np.abs(den_comp[a] - den_gt[a]))\n",
     "a"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": 10,
-   "id": "10dcc6ee-af0f-412d-8be9-d7500be6cd08",
+   "id": "cd441b5f-f7d8-476c-9b8f-2c1a7f92380f",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "text/plain": [
-       "0.6935358751064209"
+       "14"
       ]
      },
      "execution_count": 10,
@@ -221,34 +221,44 @@
     }
    ],
    "source": [
-    "np.abs(den_comp - den_gt)[a]"
+    "data.kstar[a]"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": 11,
-   "id": "19b1c9e6-bf00-44f9-822f-f13412263ce5",
+   "id": "1338cf31-e49d-47da-bbb4-1225db9aefd7",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "text/plain": [
-       "130"
+       "[<matplotlib.lines.Line2D at 0x7485f7bad5d0>]"
       ]
      },
      "execution_count": 11,
      "metadata": {},
      "output_type": "execute_result"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhYAAAGsCAYAAACB/u5dAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABGLklEQVR4nO3df3xU1Z0//tckIckkkEkg/EggIQGEqPyQH4II/qDQIvJVEWtdpS1ru65a7OraWonb1rWrJrptH239orXa6nargi6gVgWLICCiEIKB8FOBhAQIxEAyIcmQhMz5/BFmyEzu/L73nnvnvp6PRx4PMnNz72HuzL3vOed93scmhBAgIiIiUkGC7AYQERFR/GBgQURERKphYEFERESqYWBBREREqmFgQURERKphYEFERESqYWBBREREqmFgQURERKphYEFERESqYWBBREREqpEWWGzevBk33XQTcnNzYbPZ8Pbbb2t6vIKCAthstl4/S5Ys0fS4REREViItsGhtbcWECROwbNkyXY5XVlaGuro678+6desAALfffrsuxyciIrICaYHFvHnz8OSTT+LWW29VfL69vR0//elPMXToUKSnp2PatGnYuHFj1McbOHAghgwZ4v157733MHLkSFx33XVR75OIiIh8GTbH4oEHHsBnn32G5cuXY/fu3bj99ttxww034Kuvvop53x0dHfjb3/6GH/zgB7DZbCq0loiIiADAZoRl0202G1avXo0FCxYAAGpqajBixAjU1NQgNzfXu92cOXMwdepUPP300zEd780338Rdd93Va/9EREQUG0P2WFRWVqKrqwujR49G3759vT+bNm3C4cOHAQAHDhxQTMbs+bN06VLF/f/5z3/GvHnzGFQQERGpLEl2A5S0tLQgMTER5eXlSExM9Hmub9++AIARI0Zg//79QfczYMCAXo8dPXoUH330EVatWqVeg4mIiAiAQQOLiRMnoqurC/X19bjmmmsUt0lOTkZRUVHE+37llVcwaNAgzJ8/P9ZmEhERkR9pgUVLSwsOHTrk/b2qqgoVFRXo378/Ro8ejUWLFuH73/8+fvOb32DixIn4+uuvsX79eowfPz7qoMDtduOVV17B4sWLkZRkyJiKiIjI1KQlb27cuBGzZs3q9fjixYvx6quvorOzE08++ST++te/4vjx48jOzsZVV12FJ554AuPGjYvqmP/4xz8wd+5cHDx4EKNHj471v0BERER+Igos/vM//xNPPPGEz2NjxozBgQMHVG8YERERmU/E4wGXX345Pvroo4s74JACERERXRBxVJCUlIQhQ4Zo0RYiIiIyuYgDi6+++gq5ublITU3F9OnTUVJSgvz8/IDbt7e3o7293fu72+3GmTNnMGDAAFa9JCIiMgkhBM6ePYvc3FwkJAQugxVRjsWaNWvQ0tKCMWPGoK6uDk888QSOHz+OPXv2oF+/fop/o5SXQUREROZUW1uLYcOGBXw+plkhTU1NGD58OH7729/ihz/8oeI2/j0WTqcT+fn5qK2tRUZGRrSHJiIiIh01NzcjLy8PTU1NcDgcAbeLKfMyMzMTo0eP9qlH4S8lJQUpKSm9Hs/IyGBgQUREZDKh0hhiWiukpaUFhw8fRk5OTiy7ISIiojgRUWDx05/+FJs2bUJ1dTW2bt2KW2+9FYmJibjzzju1ah8RERGZSERDIceOHcOdd96J06dPY+DAgZg5cyY+//xzDBw4UKv2ERERkYlEFFgsX75cq3YQERFRHIgpx4KIiIioJwYWREREpBoGFkRERKQaBhZERESkGgYWREREpBoGFkRERKQaBhZxpM7pwtbDDahzumQ3hYgMiNcI0kNMa4XEkzqnC1UNrSjMTkeOwx713wJQ/HfPfdY5XdhRfQY2mw15WXa0dnQhPTkRrR1d3r/zPD95eJbP7/Y+CTjS0IqpBf0xIS/Lu88XNx9GyQcHvL9/e+JQjM93oKqhDSOy0zF+mKPXcTxt6tmeycOzev3/65wu/F95LfYeb0bBgDSMHZapuJ3S67Kj+gyaXJ3ISkv2/l8DvcaxnAOr8bxWro7zqKhtwqB+qZhz2WBDvG67ahuxvfoMphb0x6CM1ICfgZ6Pe373/xz4/63SZy09ORE1Z9p8Pk+h3kNKx+u5fTj/h2D77vn59rTN81lW+j+t23cSDS0dmF00yOdzHchvPjyAt8qPIS05EbmOVPRN7YPZlw7C0Kw0n9d0R/UZ7DnhxMG6s2ht78KOmkYIAdgA3Dk1D1ePyu51DQp0jTDiZzPQtddIbbSimFY3jUZzczMcDgecTqdhFiFbUVaD4lWVcAsgwQaULByHO67M9z4f7CLU8289y7II+P675z5XlNVg6cpKqPGiX3NJNp799ni8W3ECJWsOhP6DHmwA7rmmEAP6paD0gwPe9tgAlN528f+/oqwGj66sVNxH8bwi3HxFruKH+Y3t3a+LEqXXONQ5oIt6vlb+nrkt8OsWKoBUw0/erMDKnccVn1t81XA8sWBsr3N9w+VDsGbvSfS8Eil9fgD4/L9tF55XYgOwdF4R7r1uJADfQOHAybOKr99dU/Pw49mX4NcfHlT8PwR6X/a8Pry3uw4vfVIV8PXp2eZvXToY/fv2wfKyYz7beD7Xgc7Ppb9YA1enO+gxFk4ailU7j8d0nUmwAY/eUISGlna8vKUKwmCfzUDXXr3aGCo4jWZf/oF2sIA3nABUbeHevy0fWNQ5XZhRusHnIpNos2HL0lnIcdi7A4FVlT4XPc8bt2hIPyx4fivCeQUTADw8dzR+/eGXqv8f1GYDsLX4GwCA6SUbwv4bT0BS53SF/Luer/Gu2sZer2MCgE+Lv8FvHn6U3q89JdiAT5f2ft2e//gQnv3woPd3/wDS/xjRfPPbVduIW5ZtDbpNXlYqjjedC9j+QHreOCKx4IocfFXfir0nmiP8S2U937dA8CAvFjYbUKpwc7z7L9vx8Zdfq3uwCPm/BjKE+hxo3cblF7449Tx8uAGN/+fLP9C+deJQrP7ieK8vWf5B+22ThuI337nCJ8Dp2Tumxf893Pu35YdCqhpae705u4RAdUMbgO5vSP6Bg1sg4l4HN2CKoALovni/sqUa1xeFvwaMQPdrcu3ogSgNo/fE8xpv/vJrxdfSfaENj82/NJKmx71XtlQFvYm5BbzvXc/F692KEz5BBeB7vnpegGLpOdpefSbkNrWN58Lal79o79tvV9RF+ZfKPO9bz7dKLYIKABACKPY7P3VOl/SgAvB9DWRRum731CUEPtp3CufOd6n+7b7O6ULxaoVrlgAeW7Wn12eq59/9ZUsV/nzhM+zpESpde8B7j3EL+AQPnn1m903u1Yu2cudx5DhS8fzGw71ei2BfHPRg+cCiMDu9V5eqDUBBdlrQN6+u3TwSvLzlCOaPHxLR3wgAz6w5gHcqToS1fVtHZ9AA7eUtR3D3zAL2WlxQ53QF7Wb3+HBvHf7ns6PesfRAr68ngPQEb/43Sv8LZaiejKkF/aP7j5lIos2Gguw0AKFvbrHyD67LjzZqd7AI9HwNYhVt71hhdjoSbAj6+v/inb3ef3u+3atx7KqG1oC91IGCLqUhcLdAWEPYXULg4wPKAeWyjw8rfr4Fur8UBwpytMZZIUou9Lt+eqhBbjskcgvgt/+IvIfl7TCDCgDYVesMGqC5BfD+7jpmsF9Q1dAaVkD76taj3gtfqO1f3nLE+/oG671bUVaDGaUbcNdL23B1yQY8/f6+XudlQl4WJgxzhPefMSEbgKcXjvVeqD03Ny15zs+Ksho88PoX2h4sDP6vQSxe3VqFq0u631MzSjdgRVlN2H+b47B7827CsXLnceyqvRiY9Xw/R3psT6KokgQbvEGXZwbOWztiy6tLADArQO9xqOunp/dSb5YPLJQu1kIAz60/hGUfH5bSJqPY9JV2gVWizYaBGSkht3vy/f0Rf/DjUZ3ThdMt7arfyHpefJRulIk2G9KSE3x6MgSAP31SpXhephbGb6+FAHDt6IsX+EhvbtFwC+CNbd03JiO4NKefKt3rdU4X/vPdfd5rr6d3LJIvET3bMTY3dL7eS5urvMf275krXlkZ9rE3BxmO+peZIwAAT72/zxu4PPJ/sSXrCwANLR2YPsL3s3X9mOBD1T2DHL1ZPrAI9K3j9e3WvpFp7aYJORg/NLxvt9FcdOKJ59vVj9+oCCtR2F+oWCQtufsy4H+j9Hw7be3oUuxy9j8vP3mzIqyhGiO7akTwwOiVLdU+v+sxhv2HDYcMM/S6r+6sKp/DqobWXo/1zG2LVGZan5DbfLCnzjv84f9+9gw7heIJSpTYAAzol4yrSzbgpU+C50JFQqD7c7Zg4jDvYzNHDcCmg8HzbUoWjpM2jGz5wCLHYcetE4fKboblvF1xArc+vxUF/cOLqGO56JiZ/7eraK5VS28sCvp8W8fFqYs9b5S3ThyKO67MD9rl7zkvu2obA04zNZPrRg8K+nzPoSMAluxJU+NzqDSc4Mlti0Zqn8SQ23h65wINZby85Qh21TYGLSAWLK/mrmn5eGbNAU2CwC4h8PXZi4nPWw6dDnqc+68fIXVKsOUDizqnC6u/MP8F0YzcAqg+E95FSs2EMTOJNUHwhrFDcPOE3IDPB7uYp6d053YH6/L3nJdwZoSYQfO5zqDP9xw6CvbtNV7FcvMPa+ca8gwN5DjsmF3UO4B0C2DB81uD5l14kv2VjBiYrlkyb4IN6JMY/u1aqUdIT5YPLLTO7KbYqZkwBpirrHGkCYID+yb7/H7d6IF4ZUvswxNK334SbBfPS7zMCGlq7Qi5zYd7u6ewWvHacWVhf1U+h0o3PqFhsqHN5js0EGgqfc9pn0rDrzkOO2aOylb825yMVPUa7Oeblw7Goa9bwt7+wz2npF7fLB9Y6JHZTbHxT5qLRSzZ4DJEmiCY3c/34uZs6wya9yAQ/cX8459e7w04JuRlYVocJG42tIQOLP7ns6Ooc7qCzg6IVyMH9lVlP0qvnZa9kluXfsMnOM5KSw6ydbdAw69jhvRT3H5A3xT0TQk9JBONf+w7hbd2HAu94QWxfK7VYPnAIsdhx6M3BB+DpvAMy9QuYn9u/aGY9xGoToPRey4iGStNTvL9SDe0tgcdiw0nc9zTw+NviMP3fH9/ekG4zTSsdftPhdzG88163wmnDi2KP54Eyp4SbTZVeyX9+e+3sS10ABlpoHO6pR0t7V0Rty0ckXaMaTpkFQbLF8gCgHFxPPdeT8eaoquqGI7lZTX48exRMV14gtVpiJciXB3nfdeQyE5PCVpIKFTmeLCS1W/tOIbvXjXc+3s4F+t48dLmw9gQIis/VouvHo7/2XpU02PoLdD7Se8S4V+fbQ/6fDSBTl2zdte/SM0dK3dBQsv3WADBC56QMahR7CVQnYZ4Sgp1+t3cHWl9oq610Np+PmjJ6l++49vbE+piHU+0DioAoDbMxGY9tXWcj/pvg5VA1/smOLBf8Bo6W5bOinhWxZGvtUmYjGaofvywTNXbEQkGFgBWx8E0uXinRrGXXnUabOomhRrBCafvt6bDX7egKMCYMNBd9jfQUNDZc+fDXpcECH2xpshsCFDGWaaz56IPLMJNdo0mufp4U3jbBhrWA4BhWRevA9FcE17fpk2+lhmH6i0/FFLndPVaoImM59F5RaoEAHdcme9dBv7Hs0YZYvlnLb38SRX+HCR50xMcKL22/VKDXx78g71wEuLI3EK9J4IJp2e452rSkSyCt7/ubNDn/RcAU+oEGJCejGONgQMUT25Ia7tycKXVBKHWGHqJZLF8YCF7vm88uWlCDv6+S93VJAFgyayRuPfakarvN8MeulpfPIg2eTM9JQmX5fTDvgAX7V/dEl+9PRTa7mO+CauRLOSV47AjLTkRbR3KCY51Tpc3qABCrxYaieklG3x+DxUEbD3c4PN/6pkbovckwj+okLiuN8sHFsyvUM+TC8ZpElj0TBAkdQXrCWptPw+bLfBl9PYpw3x+t1LyplVVNbRiV20jJuRl+dxsw+1dSElKCBhYKK0aKiu5+q6Xtnn/T9eOHhhz9VursXyORY7DjkyLfHPV2qqd4c+zJvlGZKcF7Qla9cVx7D3RHNa+6pwurOT5t4Qd1Y0Bp26HKokdjFJVSz2Tq/1nVHn+T+VHGy1XCC1Wlu+xALq7fJtcwUv5Umi/+vs+2U0gBYk2G7oUVi8bPkCd3roVZTXevBWKf1MKsgJO3V7w/Nag+RHuAHfoOqcLOQ47bps8DP9X3h2gal3bwp9/YAF0/58guhO9o1kA0Kos32MBAGdDrA9A4dHqc1dvoPnhZrRl6SxN9isEsKu20TBLepP2BvZNxoS8rIAVi0OVxD7f1fvmDQA7jzYCAKb2qN4azZTPcCn1UvsXlwMu3CBtwLcnXxz2Y6Hm0CwfWLy46TCaY5hCRdq79fmthi+9bWSBvvGd6+yKac2UZ9YewIJlWznmbCFft3TgxU2HkeOw4/GbLvM+rnSzVSqJLQLclZV6A0L1VMRSMVeph9o/sLChezn1B17/wqec9j9dmRf1cfUi+8uypQOLOqcLpWsOyG4GhWCW0ttm89mRMzGtmfLKp9UMKiyoZM0BvLj5MK4eeXExrpcXTw6r+FyfBOVbzuSCrIja8Mb2GlxduiH0hhpwmmDYvLGVgYU0VQ2tvDCaRKAFgWJhhjHTYAV91MLAjSJV+sEBzP3dZu/vh75uRfGNl3p/TwhQfC4hQBnJiEpnO114rMe0VK0E2v0He05qe2AVZKXLnZBg6eRNzzghM37VoeVrGW+lt8MRbJ0OtcXbmimkLQHfwPyZNQfwzgMzvL+vefAajBmSocmx+YXQ+CzdYxHpktQU3NqHrtVkv4G+/cSzXbWNWKpTUAGoUzKdrMst4FO1cnCGdisds/ZQaBwKkeza0QNlNyFufHpImy77+68fGfelt3tavr0GtyzbqutQzb/MHIEch53DIRSVBJvvWhtaynHY0T+d5eODkT0UYvnAovzCNCeK3X+9p00dixc2HrbMDa/O6ULxanWnb4bz2t09swAAS9xTePwzJR6dV6RpL4W/vimWHsUPqV8qAwuphBky+ExCq257NZZMNwulssaxemVLddDne94k2M1M4fB/i948Idfn97KqM+b7MhCkfD1FxvKBxZSC/qE3orAESPhWRVqyNd6qWtzYX95yJOjzAhcDNyvlsVD0Qn0a7/nfcp9pzJ7ZTYEqb8ZST0Utp1vapR4/nli+P4kXUvXcNmko3io/rsm+2zqUK/aZSTgrQeY47JhamIXtVeoN0YXqSbLijBuKjdKn8ZRfhVzPNOamtg6UrjkYdCaHZ9GvhZOGBdlKW8GWTKfIWONrIOniysIBmuw3HmYsrCirwYzSDWEVpBo5sK+qx04M0cVrtRk3pL765nOoPdN7uLJLiJBBhYdbgAvZqcQmufA4AwtSzdBMbW5OP7p+pCY3PqHTbPg6p8tn6qjeBameXjg26PNWmnFD2jjW6EJe/97BfwIiW0OIKW/xgYEFqeb1bdqs57HoquGa7FcvSgmZWlQSDSRU4MB1WChWw7LsirNCHp1XFNF3Z5n5k3pNl7UCBhakmvcr6zTZ77u7TmiyX70oJWQmQL/hnVA9IyznTbEapBBU9Em04d7rRmLmqGyFv+gtwQbcJjHHglNY1cPAggzvmTUHNLnx6dXtmuOwY1qh7+wjAWDzl1/rcvwZIRZr0rP3hKzDdqH7YcyQft7HHKmBb973Xz/SZ9l0vXWcN3+CuMfR061SvywwsCDDi4c6FjkO3290AkDxykrsqtW+QBtnhZBRBFqEDAD+uPEImto6dGyNr7bOLmnHVtvystqoVy1WAwMLMrx4mBVy9tz5Xo+5ASx4fqv0HIcFE3M5K4RiosYshC4hcLolvMBCi1yMk85zoTcyEZmrFjOwIMN7dF6R6W98/QJ0AQsDLFn+9hcnmGNBhtAnkdUv1SRrmNPygcVP3qyQ3QQKwb9csBmlBUkMk53jIPv4RB6dXb3H7TxVOxn8Rk7WMKel02B31TZi5U5tKkUShav3h1/fb23MsSAjSLTZMMBv1dIVZTUovlADJsEGlCwcp1ndlWFZ9riqvplos0krfmfpHovt1WdkN4HCUN8cX2OfPSXY5Fe+lH18IqD7fZiZdjGwqHO6vEEFoH3OQLwtxb5l6Sxpxe8sHVhM5QJkpvDmjvgt8/vCoknSK1/KPj5Zhy1I1qX/+7CqobXXjKYuIVBerf1Mqngg88uCpQOLCXlZuG3SUNnNoBDe2F5j+vHVQDUzBvRN0fzYodYK6cnsrzPFj0Ar/f7b8i8UZ1mRcVg6sACA33znCtlNoBDioY6FTFuWzgprO89CaURGkOOwIyu9T6/H3QI40yqv3gWFZvnAgswhLTm6tyozysPrEvUfzyYyAj3LbHfGUeVN2Sw9KwRg169ZtHVE/qHXM6M8NGPfsZXGs4lk2lXbiHM6VsPcf/KsbseKd5bvsahqaJXdBAohmsqb4WSU8z56UWF2OoJUWybSxfaqizP1blm2FV+fVR7yiLcZHPHG8oFFoAQhMo5oKm8Gyijvmashe6EeI8lx2FGycJzsZpCF/e/n1Vi5M/QMsKkFWQEr2ZIxWD6w4Px944um8qZSwOhfCOqN7XIX6jEaTjslmX7x9t6wehFT+iTqXELOnGRe1ywfWFB8ynHYMWrQxeDCU4XOn8yFeoiItCLzusbAggzvb58fjervBmdcXKrcU4VOKaeGa2UQmUe8rUKqFZnXtZgCi9LSUthsNjz00EMqNUd//KZqfMs+PowXNx+OaR+eIa9whkhk02JJ6GD4GSA1nPIrvS+EwNbDDWhtV7eY1Vf1LTjfxamhoci8rkWdAVNWVoYXX3wR48ePV7M9uuOsEHN4Zs0B3DwhN+acGP+/12utjkCVN43g/d11mD8+B5u//Fp2U8ikHltViY1f1vs81tklcNdL2zTJh+jk3OiQZK4BFFWPRUtLCxYtWoSXXnoJWVlZardJV5wVYg5aVd98/9+usXzS4pPv78fVJRuwdGWl7KaQSW04WB+wDkokIcClQ/qFtV0fzo0OSeZ1LarAYsmSJZg/fz7mzJkTctv29nY0Nzf7/BBFKppaFuHomYdhZQKs60HyJSeFviVdMqgvkhKZHmhkEZ+d5cuXY+fOnSgpKQlr+5KSEjgcDu9PXl5exI3UEodCzOFfZo4w1NRglgonipxbhSGMIQ5+GTC6iAKL2tpaPPjgg3jttdeQmhreyS0uLobT6fT+1NbWRtVQrXAoxBzunlmgyn7853a//cXxqPYxo3QD7nppW8x1MPRO1AzGQE2hONXk6ox5H0bOV6JuEQUW5eXlqK+vx6RJk5CUlISkpCRs2rQJf/jDH5CUlISurt513VNSUpCRkeHzQxSpd3ediHkfnjLfPT35/r6Ieh3CKRVuVktvLJLdBCKKAxEFFrNnz0ZlZSUqKiq8P1OmTMGiRYtQUVGBxMRErdqpGQ6FmMMzaw7EfPNWKvMdaVJoOKXCzah4XhHuvXak7GaQxZ3mcuhxIaLppv369cPYsb7VC9PT0zFgwIBej5sFh0LMwRMAxJJn4Vloq2dgEGlSqNI+wpkvbvTu23uvY1BB0bNBneTfY43m7/kjVt40VEIgBRZNsRf/m7nSQls/n39ZRO8B/33oVQdDa1wvhWLxzcsGB3yOX97kkDk8G3NgsXHjRvzud79ToSlEytS8efvP7V4wcWhM+3hqwTjV54vLSKKMlzwRMh4thpttRsp6NiiZCyxavseCjG/1j642bBGrrPTksLYTEXQUy7hmxkOeCMmzbt8pVfaTae8Tcpumtg6W9A6DzMRyLmpPRIZbL4XMRa0UonCmo+45wSKL4fJ8YdB7qJY9FmR4C57fqlqXnhp1LOJRPOSJEJEvWV8YGFiQ4QkBLF1VGXOXnhp1LOKVUYeaiCg6iTabtC8MDCzIFIQAdh5tjGkf4daxYLluIjK7n90wRtoXBuZYkGnEWgsinDoWK8pqsHRVJYTofq5kofqzPoiItPbs2oO4+Ypc9ljIwG+l5mADMLkgK6Z9KNWxeOzGS70fPM9QiYjDct1EZC0yZ3pZPrB4bv1XsptAYVg4aajmdSzitVx3uF7cfFh2E4hIJTJnelk6sKhzuvD6dmOttkrKVu08rnnPgWeopCerTMN8cdNhlHxwQHYziEglMmd6WTqw4AJk5iEAlFfHlrwZSo7DjgfnXOL9XWZWtd5K1zCoIIonMnPDLB1YKH1DJeNqcmm/8uG8sTnef29ZOku1D6fMRcjC6ekx+BppRGQilg4slJL5yLgy7eGVz1aLmj0VgW/c2ke2O6rPaH4MIiIPSwcWZB7RzAox+lLlgdhUDja4YBMR6cnSgUWd04VHV1aG3pCku+eaEZbIddDC5OGxTdMlIoqEpQMLJm+agw3A3TMLZDfDtBiQEZGeLB1YFGany24ChWHJrJG8ORIRmYSlA4schx3F84pkN4NCuDQnQ3YTTG3r4QbZTSAiC7F0YAEA9143UnYTKIT9dc2ym2Bqd720TXYTiMhCLB9YkPE9v/Ew1+sgIjIJBhZkeEpLmxMRkTExsCDD81/aXE1mrXVBRGRUDCzI8B6dV2TgWSHhRSYMYIjIKhhYkOHdMmFo6I2IiMgQGFiQpQmdlt9q6zivy3GIiABgRVmNtGMzsCDS2IqyGqzZczLs7bm0BxHF6rFVe6TNpmNgQYZX33xOdhOiVud0oXgV16MhIn11CSFtNh0DCzK8Bc9/GlW3nl7DHMFUNbTCLb8ZRGRBu483STkuAwsyPLeIrlsv1pkYQoWpHIXZ6Ujg0AYRSfDsmoNShkMYWJAp6NWtp3YMkOOwo2ThuMDHY9BBRBqRNRzCwIJMQa0iWf5DKu9UnIh5n6HccWW+5scgIvKXaLNpVlwwGAYWZApqFMlSSqQs+WA/1yEhorg0d+xgKcUFGViQ4T0waxTuvTa2VWjrnC7FREquQ0JE8WrtnpPMsZCB31aN71uXDY7q7+rPtnv/PaN0AyqPOXslUmq5DgkRkUyyvjhZPrAoP9oouwkUwu7jzoj/xtND4eEWwLNrD+LReUU+2xXfeKmB1yEhIooecywkUWNKIWnr52/vibiORc+gwqNLCIwfmunz2C1X5MbSNE1woggRxSrBBjy9cCxzLGTI789ucDNYurIyomGrwuz0Xo/Jit4jZeMcVCKKkczvzJYPLFo7umQ3gcIgAJRXhz9s5R+lJ9psytE7O6yIKA4JyFsvxPKBhdI3WzKmSL7I+w+d/OyGMawnQUSWwgJZkjBxzxxsNmDS8KywtlWqV/Hs2tClbeucLlTUNoW1fyIio5M1/Juk+xGJIpRgA0oWjgs7CFSqV+GJ3APtY0VZDYpXVYZcMOyN7TV4bPXFoOXTQw24YWxOWO0KR53TxcCFiGIWcPhXB5bvsSDje/tHMyIaxlBa+CtY5O7p4fAPKvxv8HVOFx5bVemTFPXathrVAoEXNx/G1aUb8OHeU6rsj4isS+bwr+UDi2iW4yZjU1r4K1DkLhB4afOjp33HJqsaWnvleqpVgObFTYdR8sEBqZncRBQ/whn+1YqlAwulsXgyntrGyD8c/pF6sMg90NLmwwek9drOnxqVO+ucLpSuORDTPoiIepKVuAlYPLAI9E2VjCUvS9sxwkBLm/v3cOQ47Lh0SD+fxxZNy495DFOpJ4SIKBYy6/ZYOrAI9E2VjGVQRqrmxwh3LDI30zeImDEqO+Zj831IRGpbMDFX2qxHSwcWgb6pkrHUN5+T3QTN2MD3IRGp7+0vTjDHQhYWTTK+aHIs9KJWsiXfh0SkJuZYEAWhdY4FEVG8USOxPOpjSzkqUQT+uPmIlONy6icRmdWVBf2ZY0EUyAeVddhVG/4CZMH41y15t+KEKvslIjKSsuozzLEgCmZHBCubBqJUt6RkzX6W0CaiuKNW8b5oWD6w4E3FHKYUhLcAWTBKdUtkfviIiLTCOhYSVTW0ym4ChXDjuBxMyIs9sFCqFxEswWnFjtqYj0lEJIOsBcgABhaKZZrJWK7Iy1RlP0r1IpbOuxQ5Drtiz9XPV1eyR4uITEnmFHbLBxayIjoKX6mKeRD+H7abJ+QCUO656uIwCRFRxCwfWHB1U+PTIw9CqecqUeI8cCIis7J0YMHVTc0jLVnbt6pSz9WTt45jjxYRUYQsHVhwdVPzOKZRWe9g65DcMSVPk2MSEcUzSwcWTNw0j08PnVZlP/5DX7cs+5TDYUREKooosHjhhRcwfvx4ZGRkICMjA9OnT8eaNWu0apvmchx2zC4aKLsZFIblZTUxJ3DWOV1YutJ36EsAKF7F2R9ERGqJKLAYNmwYSktLUV5ejh07duAb3/gGbrnlFuzdu1er9mnuujEMLMwg0gRO/16IFWU1qGpohdLIF4tkERGpJ6LA4qabbsKNN96ISy65BKNHj8ZTTz2Fvn374vPPP9eqfZqbMCxTdhMoDJFUkVNKyn1s1R6kJyfCprC9zFUAiYjiTdQ5Fl1dXVi+fDlaW1sxffr0gNu1t7ejubnZ58dIajVKCiT1JNgiqyKnlJTbJQTaOtwovc23QJYNQMlCzv4gIlJLxIFFZWUl+vbti5SUFNx3331YvXo1LrvssoDbl5SUwOFweH/y8oyVaS+4Nrbhvb1kRkRV5Aqz03v1TNgu9Er47+flxZOlVqiz2ZT6UIiIzCviwGLMmDGoqKjAtm3bcP/992Px4sXYt29fwO2Li4vhdDq9P7W1xlp/YUpBf9lNoFDUiP0C7OOev5YHnBXChE4ioshFHFgkJydj1KhRmDx5MkpKSjBhwgT8/ve/D7h9SkqKdxaJ58dI2AVufAuej2xKqFKSpkB3gqZ/sOAW3fkXSkHENc9+zKmoRGRKMq9dMdexcLvdaG9vV6MtUvBbqfEFu/krKcxOh/8Igyf5U3lNEKE4KySc43Ikg4iMqHilvGn0EQUWxcXF2Lx5M6qrq1FZWYni4mJs3LgRixYt0qp9muOy6eYQ6OavJMdhxz0zR/g85kn+VF4TJPCMk0iOS0RkFG4Ar2yplnLsiAKL+vp6fP/738eYMWMwe/ZslJWV4cMPP8Q3v/lNrdqnufTkRNlNoDBEMt0UAGYVDfL53ZOg6T/0FWrGSaTHJSIyipe3HJHSa5EUycZ//vOftWqHNLuONcluAoVgQ2TTTSOx6v6rcUV+VtjH5dAHEZmFp/if3rmEll4rBACq2M1tfBrezAdlpAZ87kfXj5Q6FZWIKBY2yCn+Z/nAYkB6suwmUAgiwuRNtWTY+4TchmVQiIh8WT6wKOAKp6bAJEoiosh4ptnrzfKBxeThyuPrZDxpyeG/XUWAilj+c7v/vutETG0iIjIqWcnnlg8sWCDLPNo63DH9vdKy6aVrDrCWCRHFnUSbTbOk91AsH1iQOaixAmn50UbFipw7jzYqbq/HDBCuVUNEWtiydJa05HMGFmQKj84rijnyDnQT572diOKNzN54BhZkeA/MGoV7rx0Z836mFPTvveopgMkFzLMhIlILAwsyvO9eNVyV/eQ47Ci9bZzPY2r0hBARGY3M3DEGFmQp/mOON03IldQSX0wgJSI1yZyez8CCSCLbhQxRLoZHRGqSucYRAwuyFP86Fq99flRSS3wprbpKRBQtJm9KxC5o43s3miJWCjM96pwuFK/yrWPx/MbDeHHzYcVd2MJYpKSxrSPytilgngcRqcn/S5SeLB9YsAva+ErX7I84AFSaQVrV0Aq3whPPxFAk6z9W75H6ASYiUiJjfSUPywcW6cmJsptAIXiW/o1VYXa6Yh9EJPs/19nl87uA3A8wEZESmesrWT6waO3oCr0RSRfJOiGB5DjsWDqvqNfjkdTTb1N4v3CBNCIyGlnrhAAMLFB5zCm7CRSGWNcJ8bj3Ot9CWwk2RFRPP02hh0vmB5iISImsdUIAiwcWdU4Xnll7QHYzKAy7jzepsh//fIj7rhsRUT39lKTegYXMDzARkRJZ64QAFg8sAiXzkfE8u+ZgzHkMSrNC/rjpSMz7lfkBJiIyGksHFoGS+ch41MhjUAokgyVu6rG6KRFRvLF0YJHjsOOGsUNkN4PCEGkeg9KKpYXZ6UjwCxbUWI49mDqnC1sPN2i2fyIio7F0YFHndGHNnpOym0EhRJpgGUiOw46Shb6LkP1s7hjN8iNWlNVgRukG3PXSNk32T0RkRJYOLFgcyxxe+v5k1fIY/Pfz/2m0CJknn4M5PERkNZYOLLg+gznc89dy01W3ZGIwEVmVpQOLHIcdt00aKrsZFIJbmK+6JRODiciqLB1Y1DldWP3FcdnNoDBoVd1SKclTDTkOO24cl6PNzomIDMzSgQW7q83DjNUtJ+Znym4CEZHuLB1YsLvaHNSaFQL0rrz53u7IlmRXu7YF339EFG8sHViwjoU5rP7R1RHPChEKC6crVd58du1B7KptjKl9RER0kaUDizqnCx/uZR0LoxuUkarKfpSGvgSABc9vNd2sEyIio7J0YMEci/gVbuVNz7aPrdqjfaOIiCzA0oFFoBsNGYtNpUwEpcqbHl0aTA+xcbERIrIgSwcWOQ475l7OHAujU8qXiJSnBkagXI1EBgFERKqwdGDBtUKsY0bphoB5FLYLs06IiCh2SbIboJY6pws7qs+gydWJrLRkHKxrxod7TyIhwYbzXQJdwg1XRxda288jMSEBhQPScersOdnNpjB88zcbkZtpx8T8LNw1LR8T8rJC/s3/flbt87tbAMUrK3Ht6IG9tv3GmIHYfayp1+PPrf8SGw/U45EbxmB71Rm8taMW9Wfbe2330ubDeKusBiec55Bh74PBGakYkJ6MhpaOkO38579sR25mKs60ht6WiChcdU6XZgsshmITQqvag8qam5vhcDjgdDqRkZGhyj5XlNVg6cpKFTrMyQxumzQUv/nOFahzulDV0IrC7HSfD9D9fysP2BP1r9eMwJ8+OaJXU4mIpHnmtnGqLeAIhH//Nn2PRZ3TxaDCYlbuPI4cRyqe33gYbtFdQKtkYfcHaFdtY9DhLQYVRGQVSy/00urdc2H6HIuqhlYGFRb0/3982DtVuOciZR/tPyW3YUREBiEAlFfrXwDQ9IEFy3ITcHGRstSkRNlNISIyjCaX/vlbpg8schx2LOTS5wSgIDsN+QPMtVAZEZGWtnzZoPsxTR9YcOlz8qhvPof8/gwsiIg8Dp46q/sxTR9YsCw3eeyobkRrR5fsZhARGcb1ClPstWb6wKIwO112E0hn/VKU37ZTCrL4fiAi6mFKYX/dj2n6wIKsx56chNv88mpumzQ0rMJZRERWom+lqm6mDyyqGlplN4F0dr4r8CeF7wcioovy+utffdP0gQW7vq2n9VwnVu70TdhdufM4dtU24tOv9M+AJiIyqmONLt2PafrAAgDrWFhMu1v58T+sP4RlGw/r2xgiIgOrPN6k+zFNH1iw8iZ5bDhQL7sJRESGUt3QpvsxTR9YpCez0iJ1Y4BJROSrIFv/2j6mDyxYt4A8EmzA4unDZTeDiMgwWtrP635M0wcWrg79XzQyph9dPxLXjs6W3QwiIsM4yqGQyB3h9EK6YMaogcyzICLq4VpW3ozc1AL9q4qR8dhs3WOJV+Rlym4KEZFh5DhYxyJiE/KykJuRKrsZJNuFzE17cpLcdhARGcjRM/r36ps+sACATneAwgZkGQLd06qEjPq1REQGdbqlQ/djmj6wWL//JL6W8MKRsSRcGAo50XROdlOIiAxDRnVq0wcWTNYjAPjVLWMBACVrDkhuCRGRceRm6p8qYPrAgsl61pRo8y3kfvuUYVyAjIjIDytvRmFolv5VxUi+703P9/ldCC5IR0Tkb0pBlu7HNH1gwZuJNW086DsEdtJ5DjkOO749caikFhERGc8gCbMmIwosSkpKcOWVV6Jfv34YNGgQFixYgIMHD2rVtrDkOOwonlcktQ2kv+rTvksBz/r1Rqwoq8GY3H6SWkREZDzl1Y26HzOiwGLTpk1YsmQJPv/8c6xbtw6dnZ341re+hdZWuWPb57s43dTqBIDiVZXItPeR3RQiIsPwS0fTRUTVhNauXevz+6uvvopBgwahvLwc1157raoNi0TlsSZpxybjcAvA6eLaMUREHsOy9K+8GVOZQqfTCQDo3z9wWe329na0t7d7f29ubo7lkIoKmGdB6K5lUShhiWAiIqM61ujChDx9EzijTt50u9146KGHMGPGDIwdOzbgdiUlJXA4HN6fvLy8aA8ZEGeGENC9uilLehMRXXSm1USVN5csWYI9e/Zg+fLlQbcrLi6G0+n0/tTW1kZ7yICaXKy8Sd2rm3KWEBHRRQdPqj9KEEpUX+8eeOABvPfee9i8eTOGDRsWdNuUlBSkpKRE1bhwnWnt1HT/ZHyekt4yVvIjIjKqo6f1L5AVUWAhhMCPf/xjrF69Ghs3bkRhYaFW7YrICH5Ltbxf3TIWOQ476pyu0BsTEVnEtaMH6n7MiAKLJUuW4PXXX8c777yDfv364eTJkwAAh8MBu13eN8W2ji5pxyZj+Pbk7p6zj/afktwSIiLjGDs0U/djRpRj8cILL8DpdOL6669HTk6O92fFihVatS+kOqeLC08RXv6kCgBQecwpuSVERMaRlqx/ge2Ih0KMhgtPEQD8+h8H0SfJhvTkRNlNISIyjLYO/QtImn6tEFcHCyJRt2fWHMDMS7JlN4OIyDDaOvSf3GD6wGLLV6dlN4EMwi0Al4TonIjIqHZLGB42fWDRyh4LuiDRZkMja5oQEXll99W23IMS0wcWY4dmyG4CGYDNBjy9cCyy0pJlN4WIyDDGD3PofkzTBxbfvGyI7CaQAdw0Pgd3XJmPycP1rYlPRGRktWf0r+1j+sAix2HHlOGZsptBkv19dx3qnC5W3iQi6kHGkhdxsWJT35S4+G9QDIQAqhvasO8E61gQEXlk2vUfHjZ9jwUATCkMvGw7WYMN3WuFbDhQL7spRESGYZdQICsuAgtXO0t6W92cSwchx2HHN4oGyW4KEZFhfHpI/5IMpg8s6pwuPL/xsOxmkGQf7a/HirIazL6UybxERB6t7fqXZDB9YFHV0ArjFRonvQkAj63aw9VNiYh6yB+QpvsxTR9YFHLJdLqgSwisLD8muxlERIYxvL/+90jTBxZEHok2G97ddUJ2M4iIDGNygf61fUwfWHB1U/K47/oR+PJUi+xmEBEZwnWjs6XU9jF9YFGYnQ6b7EaQIbCcNxHRRbK+eJs+sMhx2HHn1DzZzSADuLIgi0EmEdEFNWdcUhLaTR9YAMD0kQNkN4EkswEYlJGKmZdky24KEZFhyEhoj4vA4lwnC2RZnQBQXt2IT75qkN0UIiLDqGrQP+8sLgKLrYf1ryxGxmKzARwHISLyddUI/Xv04yKwOH1W/9XbyFi+M2UYl0wnIvIzLIt1LCJW53Rh8yF2f1vdVSMGIMdhx11M5CUiAnBxcUa9mT6w+Gj/KdlNIAMQF+q6/3j2JXIbQkRkEJlpfVjHIhpHvmaBLLpIxoeIiMiIGts6Od00GoUSunnIuHbVNspuAhGRYew8qv810fSBxTcv4zLZdHEohENjREQXCQnLf5s+sMhx2FE8r0h2M8ggBvVLld0EIiLD4CJkUbr5ilzZTSDJGtu6pxzPuWyw5JYQERlHffM53Y8ZF4HFR/vY/W11T76/HyvKapi8SUTUw45q/XMsknQ/ogbqz+ofkZHxFK+qRNGQfrKbQURkGBl2/W/zcdFj0dx2XnYTyADcAli/v152M4iIDONYI6ebRqzO6cL/fH5UdjPIABJsQMd5LkhHROSRnKT/bd70gUW5hDm6ZEwzR2Wjj4QPERGRUWWk9tH9mKa/CgsZk3TJkLYcakB+fxZMIyLy6J+erPsxTR9YTCnoL7sJZBBuAThdzLchIvIYlsW1QoiilmizscQ7EVEPbR1u3Y9p+sCiqoGLkFH38sBPLxwLV6f+HyIiIqNKS2byZsTSkxNlN4EMwJNpw5wbIqKLON00CrUSXjQypuJVlUzeJCLqgYuQReFMa7vsJpBBuAVQe4aBJhGRh51DIZHrn54iuwlkEAk2dCdbEBERAKC6oU33Y5o+sMiTMJWGjOnReUV8PxAR9cC1QqLAHAvyGJZp5/uBiKiHZgm1fUy/umk1p5vSBY1tHWhq65TdDCIiwyiQUNvH9D0Wp1uYvEndMu3JfD8QEfWQlsy1QiI2YlBf2U0gg5hckIUB/ZjMS0Tk0dahfy+u6QOL8UMdsptABnDbpKHIcdhhY30sIiKvTw+d1v2Ypg8sWju6ZDeBDGBAejLqnC40cCiEiMhLRoEs0ydvFmanw4aLJZ3Jmv70SRX+9EkVpgzPlN0UIiLDmHnJAN2PafoeixwH6xbQRTuONsluAhGRYdQ5z+l+TNMHFj97s4K9FURERApe3nxY92OaPrB4d/cJ2U0gIiIypOoz7LGImNvN/goiIiKjMH9gwbiCiIhIkYybvOkDi/MMLIiIiBQlSljx2fSBBVfJJiIiCoCBReSSTV+Jg4iISBtut/7HNH1gkWQz/X+BiIhIEzJqU5v+rtzaKSEcIyIiMgEZ6QKmDyyIiIhIWaKEu3zEh9y8eTNuuukm5Obmwmaz4e2339agWeGzJzF9k4iISEmSGQKL1tZWTJgwAcuWLdOiPRHL7pciuwlERESGlJig/5fviOdUzJs3D/PmzdOiLVFpOXdedhOIiIgMqUNCsSfNJ2u2t7ejvb3d+3tzc7Oq+292MbAgIiJSImN+g+ajLyUlJXA4HN6fvLw8VfcvYyoNERERKdM8sCguLobT6fT+1NbWqrr/ZOZuEhERGYbmQyEpKSlISdEuwbJPcgI62lnLgoiIyF+6hG/fpq9j0XmeQQUREZGSdAnrXkR8xJaWFhw6dMj7e1VVFSoqKtC/f3/k5+er2rhwSJhJQ0REZArN7Z26HzPiwGLHjh2YNWuW9/eHH34YALB48WK8+uqrqjUsXIm2BADstSAiIvJ3Tv+4IvLA4vrrr4cQ+s+LDSQtJZHrhRARERmE6XMsEqQssUJERERKTB9YnGnrkN0EIiIiusD0gQVHQYiIiIzD9IEFERERKZORLGD6wEL/GbpERETm0McMy6YbjT2FyZtERERKEhhYRK6zyzhTX4mIiIzknIQFwE0fWMh40YiIiEiZ6QMLIiIiUpYo4ZimDyyYYUFERKSsS8IxTR9YSFi4jYiIiAIwfWDRR0bKKxERESky/V25S7D0JhERkVGYPrCQsSQsERGRGcjIFjB9YMHkTSIiImVJSfrfJU0fWHAghIiISJmMIpKmDyxSOSuEiIhIkVtCcWrTBxaJCRwMISIiMgrTBxYd57lWCBERkZI0Ccubmj6wcDPJgoiISFG6hHwB0wcWiab/HxAREWmjVcJKnaa/LQuOhBARESlq7dS/W9/0gUUnAwsiIiLDMH1gwTkhRERExmH6wIIdFkRERMZh+sCCiIiIjIOBBREREamGgQURERGphoEFERERqYaBBREREamGgQURERGphoEFERERqYaBBREREamGgQURERGphoEFERERqYaBBREREamGgQURERGpxvSBBVc3JSIiMg7TBxZc3ZSIiMg4TB9YEBERkXEwsCAiIiLVmD6wYI4FERGRcZg+sGCOBRERkXGYPrAgIiIi42BgQURERKphYEFERESqMX1gkSi7AURERORl+sCiS3YDiIiIyMv0gQUREREZBwMLIiIiUg0DCyIiIlINAwsiIiJSDQMLIiIiUg0DCyIiIlINAwsiIiJSDQMLIiIiUg0DCyIiIlINAwsiIiJSDQMLIiIiUk1UgcWyZctQUFCA1NRUTJs2Ddu3b1e7XURERGRCEQcWK1aswMMPP4zHH38cO3fuxIQJEzB37lzU19dr0T4iIiIykYgDi9/+9re45557cPfdd+Oyyy7DH//4R6SlpeEvf/mLFu0jIiIiE4kosOjo6EB5eTnmzJlzcQcJCZgzZw4+++wzxb9pb29Hc3Ozzw8RERHFp4gCi4aGBnR1dWHw4ME+jw8ePBgnT55U/JuSkhI4HA7vT15eXvStJSIiIkPTfFZIcXExnE6n96e2tlbV/f941khV90dERETRiyiwyM7ORmJiIk6dOuXz+KlTpzBkyBDFv0lJSUFGRobPj5p+MrcI6cmJqu6TiIgoHlSXztf9mBEFFsnJyZg8eTLWr1/vfcztdmP9+vWYPn266o0L195f3aBLz4VN8yMQERGpQ0ZQAQBJkf7Bww8/jMWLF2PKlCmYOnUqfve736G1tRV33323Fu0L20/mFuEnc4uktoGIiMjqIg4s7rjjDnz99df45S9/iZMnT+KKK67A2rVreyV0EhERkfXYhBBCzwM2NzfD4XDA6XSqnm9BRERE2gj3/s21QoiIiEg1DCyIiIhINQwsiIiISDUMLIiIiEg1DCyIiIhINQwsiIiISDUMLIiIiEg1DCyIiIhINQwsiIiISDURl/SOlafQZ3Nzs96HJiIioih57tuhCnbrHlicPXsWAJCXl6f3oYmIiChGZ8+ehcPhCPi87muFuN1unDhxAv369YPNpt5C5M3NzcjLy0NtbS3XIJGI58EYeB6MgefBGHge1CGEwNmzZ5Gbm4uEhMCZFLr3WCQkJGDYsGGa7T8jI4NvHAPgeTAGngdj4HkwBp6H2AXrqfBg8iYRERGphoEFERERqSZuAouUlBQ8/vjjSElJkd0US+N5MAaeB2PgeTAGngd96Z68SURERPErbnosiIiISD4GFkRERKQaBhZERESkGgYWREREpJq4CSyWLVuGgoICpKamYtq0adi+fbvsJhlSSUkJrrzySvTr1w+DBg3CggULcPDgQZ9tzp07hyVLlmDAgAHo27cvbrvtNpw6dcpnm5qaGsyfPx9paWkYNGgQHnnkEZw/f95nm40bN2LSpElISUnBqFGj8Oqrr/ZqT6jzFk5bzK60tBQ2mw0PPfSQ9zGeA30cP34c3/3udzFgwADY7XaMGzcOO3bs8D4vhMAvf/lL5OTkwG63Y86cOfjqq6989nHmzBksWrQIGRkZyMzMxA9/+EO0tLT4bLN7925cc801SE1NRV5eHp599tlebXnrrbdQVFSE1NRUjBs3Dh988IHP8+G0xYy6urrwi1/8AoWFhbDb7Rg5ciT+67/+y2c9Cp4HkxFxYPny5SI5OVn85S9/EXv37hX33HOPyMzMFKdOnZLdNMOZO3eueOWVV8SePXtERUWFuPHGG0V+fr5oaWnxbnPfffeJvLw8sX79erFjxw5x1VVXiauvvtr7/Pnz58XYsWPFnDlzxBdffCE++OADkZ2dLYqLi73bHDlyRKSlpYmHH35Y7Nu3Tzz33HMiMTFRrF271rtNOOctVFvMbvv27aKgoECMHz9ePPjgg97HeQ60d+bMGTF8+HDxz//8z2Lbtm3iyJEj4sMPPxSHDh3yblNaWiocDod4++23xa5du8TNN98sCgsLhcvl8m5zww03iAkTJojPP/9cfPLJJ2LUqFHizjvv9D7vdDrF4MGDxaJFi8SePXvEG2+8Iex2u3jxxRe923z66aciMTFRPPvss2Lfvn3i5z//uejTp4+orKyMqC1m9NRTT4kBAwaI9957T1RVVYm33npL9O3bV/z+97/3bsPzYC5xEVhMnTpVLFmyxPt7V1eXyM3NFSUlJRJbZQ719fUCgNi0aZMQQoimpibRp08f8dZbb3m32b9/vwAgPvvsMyGEEB988IFISEgQJ0+e9G7zwgsviIyMDNHe3i6EEOJnP/uZuPzyy32Odccdd4i5c+d6fw913sJpi5mdPXtWXHLJJWLdunXiuuuu8wYWPAf6ePTRR8XMmTMDPu92u8WQIUPEf//3f3sfa2pqEikpKeKNN94QQgixb98+AUCUlZV5t1mzZo2w2Wzi+PHjQgghnn/+eZGVleU9L55jjxkzxvv7d77zHTF//nyf40+bNk3ce++9YbfFrObPny9+8IMf+Dy2cOFCsWjRIiEEz4MZmX4opKOjA+Xl5ZgzZ473sYSEBMyZMwefffaZxJaZg9PpBAD0798fAFBeXo7Ozk6f17OoqAj5+fne1/Ozzz7DuHHjMHjwYO82c+fORXNzM/bu3evdpuc+PNt49hHOeQunLWa2ZMkSzJ8/v9frxHOgj3fffRdTpkzB7bffjkGDBmHixIl46aWXvM9XVVXh5MmTPv93h8OBadOm+ZyHzMxMTJkyxbvNnDlzkJCQgG3btnm3ufbaa5GcnOzdZu7cuTh48CAaGxu92wQ7V+G0xayuvvpqrF+/Hl9++SUAYNeuXdiyZQvmzZsHgOfBjHRfhExtDQ0N6Orq8rnAAsDgwYNx4MABSa0yB7fbjYceeggzZszA2LFjAQAnT55EcnIyMjMzfbYdPHgwTp486d1G6fX2PBdsm+bmZrhcLjQ2NoY8b+G0xayWL1+OnTt3oqysrNdzPAf6OHLkCF544QU8/PDDeOyxx1BWVoZ/+7d/Q3JyMhYvXuz9/ym9Pj1f40GDBvk8n5SUhP79+/tsU1hY2GsfnueysrICnque+wjVFrNaunQpmpubUVRUhMTERHR1deGpp57CokWLAIT3f+d5MBbTBxYUvSVLlmDPnj3YsmWL7KZYSm1tLR588EGsW7cOqampsptjWW63G1OmTMHTTz8NAJg4cSL27NmDP/7xj1i8eLHk1lnHm2++iddeew2vv/46Lr/8clRUVOChhx5Cbm4uz4NJmX4oJDs7G4mJib2y1E+dOoUhQ4ZIapXxPfDAA3jvvffw8ccf+yxjP2TIEHR0dKCpqcln+56v55AhQxRfb89zwbbJyMiA3W4P67yF0xYzKi8vR319PSZNmoSkpCQkJSVh06ZN+MMf/oCkpCQMHjyY50AHOTk5uOyyy3weu/TSS1FTUwPg4usY6vWpr6/3ef78+fM4c+aMKueq5/Oh2mJWjzzyCJYuXYp/+qd/wrhx4/C9730P//7v/46SkhIAPA9mZPrAIjk5GZMnT8b69eu9j7ndbqxfvx7Tp0+X2DJjEkLggQcewOrVq7Fhw4ZeXYOTJ09Gnz59fF7PgwcPoqamxvt6Tp8+HZWVlT4f5HXr1iEjI8N7oZ4+fbrPPjzbePYRznkLpy1mNHv2bFRWVqKiosL7M2XKFCxatMj7b54D7c2YMaPXVOsvv/wSw4cPBwAUFhZiyJAhPv/35uZmbNu2zec8NDU1oby83LvNhg0b4Ha7MW3aNO82mzdvRmdnp3ebdevWYcyYMcjKyvJuE+xchdMWs2pra0NCgu+tKDExEW63GwDPgynJzh5Vw/Lly0VKSop49dVXxb59+8S//uu/iszMTJ+Meep2//33C4fDITZu3Cjq6uq8P21tbd5t7rvvPpGfny82bNggduzYIaZPny6mT5/ufd4z1fFb3/qWqKioEGvXrhUDBw5UnOr4yCOPiP3794tly5YpTnUMdd5CtSVe9JwVIgTPgR62b98ukpKSxFNPPSW++uor8dprr4m0tDTxt7/9zbtNaWmpyMzMFO+8847YvXu3uOWWWxSnOU6cOFFs27ZNbNmyRVxyySU+0xybmprE4MGDxfe+9z2xZ88esXz5cpGWltZrmmNSUpL49a9/Lfbv3y8ef/xxxWmOodpiRosXLxZDhw71TjddtWqVyM7OFj/72c+82/A8mEtcBBZCCPHcc8+J/Px8kZycLKZOnSo+//xz2U0yJACKP6+88op3G5fLJX70ox+JrKwskZaWJm699VZRV1fns5/q6moxb948YbfbRXZ2tvjJT34iOjs7fbb5+OOPxRVXXCGSk5PFiBEjfI7hEeq8hdOWeOAfWPAc6OPvf/+7GDt2rEhJSRFFRUXiT3/6k8/zbrdb/OIXvxCDBw8WKSkpYvbs2eLgwYM+25w+fVrceeedom/fviIjI0Pcfffd4uzZsz7b7Nq1S8ycOVOkpKSIoUOHitLS0l5tefPNN8Xo0aNFcnKyuPzyy8X7778fcVvMqLm5WTz44IMiPz9fpKamihEjRoj/+I//8JkWyvNgLlw2nYiIiFRj+hwLIiIiMg4GFkRERKQaBhZERESkGgYWREREpBoGFkRERKQaBhZERESkGgYWREREpBoGFkRERKQaBhZERESkGgYWREREpBoGFkRERKQaBhZERESkmv8H3LglSiH2760AAAAASUVORK5CYII=",
+      "text/plain": [
+       "<Figure size 640x480 with 1 Axes>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
     }
    ],
    "source": [
-    "data.kstar[a]"
+    "plt.plot(np.abs(den_comp - den_gt), \".-\")"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "cd441b5f-f7d8-476c-9b8f-2c1a7f92380f",
+   "id": "a309b4e5-0ae4-4b55-a971-6fb05320c49b",
    "metadata": {},
    "outputs": [],
    "source": []
diff --git a/src/common/common.c b/src/common/common.c
index 1df7043f4ab5d953a2e3ed70d84f1b634a486b8f..fbf60c5b64a1d93755ef2947ba9add1fc9ecdbf4 100644
--- a/src/common/common.c
+++ b/src/common/common.c
@@ -20,6 +20,7 @@ void get_context(global_context_t* ctx)
     ctx -> n_halo_points_send = NULL;
     ctx -> halo_datapoints  = NULL;
     ctx -> local_datapoints = NULL;
+    ctx -> __recv_heap_buffers = NULL;
 }
 
 void free_context(global_context_t* ctx)
@@ -39,14 +40,17 @@ void free_context(global_context_t* ctx)
     {
         for(int i = 0; i < ctx -> world_size; ++i) 
         {
+            /*
             for(int j = 0; j < ctx -> n_halo_points_recv[i]; ++j)
             {
                 FREE_NOT_NULL(ctx -> halo_datapoints[i][j].ngbh.data);
             }
+            */
             FREE_NOT_NULL(ctx -> halo_datapoints[i]);
         }
     }
     FREE_NOT_NULL(ctx -> halo_datapoints);
+    FREE_NOT_NULL(ctx -> __recv_heap_buffers);
 
     if(ctx -> idx_halo_points_recv)
     {
diff --git a/src/common/common.h b/src/common/common.h
index 1ca7b7f634b5891263728381d7a94e5ed4d8e5df..59c4eca3f456330eb6f3da801c94476490ffd986 100644
--- a/src/common/common.h
+++ b/src/common/common.h
@@ -29,9 +29,12 @@ typedef struct datapoint_info_t {
 	#define float_t double
 #endif
 
+
 #define MY_TRUE  1
 #define MY_FALSE 0
 
+#define CHECK_ALLOCATION(x) if(!x){printf("[!!!] %d rank encountered failed allocation at line %s ", ctx -> mpi_rank, __LINE__ ); exit(1);};
+
 #define DB_PRINT(...) printf(__VA_ARGS__)
 #ifdef NDEBUG
 	#undef DB_PRINT(...)
@@ -126,6 +129,7 @@ struct global_context_t
     int* rank_n_points;
 	char processor_mame[MPI_MAX_PROCESSOR_NAME];
 	MPI_Comm mpi_communicator;
+    heap_node* __recv_heap_buffers;
 };
 
 struct pointset_t
diff --git a/src/tree/tree.c b/src/tree/tree.c
index adfd107b54985e03503c49d1bb5936ab6e29ebea..062a199a988a5db79fbfb428b7b727b0a1e59ff0 100644
--- a/src/tree/tree.c
+++ b/src/tree/tree.c
@@ -2366,12 +2366,20 @@ void find_foreign_nodes(global_context_t* ctx, datapoint_info_t* dp, datapoint_i
     {
         for(int j = 0; j < n_heap_to_recv[i]; ++j)
         {
+            /*
             foreign_dp[i][j].array_idx = array_indexes_to_request[i][j];
             init_heap(&(foreign_dp[i][j].ngbh));
             allocate_heap(&(foreign_dp[i][j].ngbh), k);
             foreign_dp[i][j].ngbh.N     = k;
             foreign_dp[i][j].ngbh.count = k;
             memcpy(foreign_dp[i][j].ngbh.data, heap_buffer_to_recv + k * (j + rdispls[i]), k * sizeof(heap_node));
+            */
+
+            foreign_dp[i][j].array_idx = array_indexes_to_request[i][j];
+            //init_heap(&(foreign_dp[i][j].ngbh));
+            foreign_dp[i][j].ngbh.N     = k;
+            foreign_dp[i][j].ngbh.count = k;
+            foreign_dp[i][j].ngbh.data =  heap_buffer_to_recv + k * (j + rdispls[i]);
 
             if(foreign_dp[i][j].ngbh.data[0].array_idx != array_indexes_to_request[i][j])
             {
@@ -2402,7 +2410,8 @@ void find_foreign_nodes(global_context_t* ctx, datapoint_info_t* dp, datapoint_i
     free(count_to_request);
     free(capacities);
 
-    free(heap_buffer_to_recv);
+    /* free(heap_buffer_to_recv); this needs to be preserved*/
+    ctx -> __recv_heap_buffers = heap_buffer_to_recv;
     free(heap_buffer_to_send);
     free(idx_buffer_to_send);
     free(idx_buffer_to_recv);
@@ -2623,7 +2632,7 @@ void compute_density_kstarnn(global_context_t* ctx, const float_t d, int verbose
             /* TO REMOVE
             if(local_datapoints[i].array_idx == 17734) printf("%lu ksel i %lu j %lu tmp_dp %lu di %lf fj %lf vvi %lf vvj %lf\n", ksel, i, jj, tmp_dp.array_idx, 
                         sqrt(local_datapoints[i].ngbh.data[ksel].value), sqrt(tmp_dp.ngbh.data[ksel].value), vvi, vvj);
-                        */
+            */
 
             vp = (vvi + vvj)*(vvi + vvj);
             dL = -2.0 * ksel * log(4.*vvi*vvj/vp);
@@ -2689,7 +2698,7 @@ void simulate_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)
 
         /* 8M points */
         
-        //data = read_data_file(ctx,"../norm_data/std_g0144846_Me14_091_0001",MY_TRUE);
+        // data = read_data_file(ctx,"../norm_data/std_g0144846_Me14_091_0001",MY_TRUE);
 
         //88M 
         //data = read_data_file(ctx,"../norm_data/std_g5503149_091_0000",MY_TRUE);
@@ -2701,7 +2710,7 @@ void simulate_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)
 
         //ctx -> n_points = 48*5*2000;
         ctx->n_points = ctx->n_points / ctx->dims;
-        // ctx->n_points = (ctx->n_points * 0.5) / 10;
+        ctx->n_points = (ctx->n_points * 0.1) / 10;
         // ctx -> n_points = ctx -> world_size * 1000;
 
         //ctx -> n_points = 10000000 * ctx -> world_size;
@@ -2817,7 +2826,7 @@ void simulate_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)
     elapsed_time = TIME_STOP;
     //id = 3.920865231328582;
     //id = 4.008350298212649;
-    //id = 4.;
+    id = 4.;
     LOG_WRITE("ID estimate", elapsed_time)
 
     MPI_DB_PRINT("ID %lf \n",id);