From 96d3003311ff72e0f07b45eb36eca73a317f66ab Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Mon, 16 Jan 2023 17:00:45 +0100
Subject: [PATCH 01/11] Update notebook

---
 notebooks/jec.ipynb | 182 ++++++++++++++++++++++++++++++--------------
 1 file changed, 124 insertions(+), 58 deletions(-)

diff --git a/notebooks/jec.ipynb b/notebooks/jec.ipynb
index 28fb1e17..2993d7a9 100644
--- a/notebooks/jec.ipynb
+++ b/notebooks/jec.ipynb
@@ -10,7 +10,7 @@
      "name": "stderr",
      "output_type": "stream",
      "text": [
-      "/work/dvalsecc/miniconda3/envs/pocket-coffea/lib/python3.9/site-packages/coffea/util.py:154: FutureWarning: In coffea version v0.8.0 (target date: 31 Dec 2022), this will be an error.\n",
+      "/work/mmarcheg/miniconda3/envs/pocket-coffea/lib/python3.8/site-packages/coffea/util.py:154: FutureWarning: In coffea version v0.8.0 (target date: 31 Dec 2022), this will be an error.\n",
       "(Set coffea.deprecations_as_errors = True to get a stack trace now.)\n",
       "ImportError: coffea.hist is deprecated\n",
       "  warnings.warn(message, FutureWarning)\n"
@@ -27,7 +27,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 26,
+   "execution_count": 2,
    "id": "b0fbce2b",
    "metadata": {},
    "outputs": [],
@@ -39,17 +39,17 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 4,
+   "execution_count": 3,
    "id": "fc676d4c",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "text/plain": [
-       "<Array [[0.0635, 0.04, ... 0.336, 0.158]] type='100 * var * float32[parameters={...'>"
+       "<Array [[0.0635, 0.04, ... -0.00543, -0.664]] type='400000 * var * float32[param...'>"
       ]
      },
-     "execution_count": 4,
+     "execution_count": 3,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -60,17 +60,17 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 5,
+   "execution_count": 4,
    "id": "7f52a8fb",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "text/plain": [
-       "<Array [0.937, 0.96, 0.958, ... 0.664, 0.842] type='1098 * float32'>"
+       "<Array [0.937, 0.96, 0.958, ... 1.01, 1.66] type='4343693 * float32'>"
       ]
      },
-     "execution_count": 5,
+     "execution_count": 4,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -82,20 +82,18 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 6,
+   "execution_count": 5,
    "id": "3628b6ad",
    "metadata": {},
    "outputs": [
     {
      "data": {
-      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAARK0lEQVR4nO3dbYxcV33H8e+POAEEFOdh61q2i6mwQGlVglmlRlSIYlHlocKRClFQRUzkylWbtiAqFZcXRVR9Ed5ASVsFWYTWqSiQBmjcEGgtJwj1RQIbCAESaJaIyLaSeAnEPKSAQv99sccwWdae2d3ZXe/h+5FGc+655+78T67z890zd8apKiRJfXnGahcgSRo/w12SOmS4S1KHDHdJ6pDhLkkdWrfaBQBccMEFtXXr1tUuQ5LWlHvuuedbVTUx374zIty3bt3K1NTUapchSWtKkodPtc9lGUnqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOjQ03JO8OMm9A4/vJnlrkvOSHEryYHs+t41PkuuTTCe5L8n25Z+GJGnQ0HCvqq9X1UVVdRHwcuBJ4BPAPuBwVW0DDrdtgEuBbe2xF7hhGeqWJJ3GQj+huhP4RlU9nGQX8OrWfwD4DPB2YBdwU83+KyB3JVmfZGNVPTKmmrWKtu775KKP/eZ1l4+xEkmns9A196uAD7f2hoHAfhTY0NqbgCMDxxxtfU+TZG+SqSRTMzMzCyxDknQ6I4d7knOA1wH/Nndfu0pf0L/XV1X7q2qyqiYnJub93htJ0iIt5Mr9UuALVfVY234syUaA9ny89R8Dtgwct7n1SZJWyELC/Y38bEkG4CCwu7V3A7cO9F/d7prZAZxwvV2SVtZIb6gmeQ7wWuCPBrqvA25Osgd4GLiy9d8OXAZMM3tnzTVjq1aSNJKRwr2qfgCcP6fvcWbvnpk7toBrx1KdJGlR/ISqJHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6NFK4J1mf5JYkX0vyQJJXJDkvyaEkD7bnc9vYJLk+yXSS+5JsX94pSJLmGvXK/X3Ap6vqJcBLgQeAfcDhqtoGHG7bAJcC29pjL3DDWCuWJA01NNyTPB94FXAjQFX9uKqeAHYBB9qwA8AVrb0LuKlm3QWsT7JxzHVLkk5jlCv3FwIzwD8l+WKSDyR5DrChqh5pYx4FNrT2JuDIwPFHW9/TJNmbZCrJ1MzMzOJnIEn6OaOE+zpgO3BDVb0M+AE/W4IBoKoKqIW8cFXtr6rJqpqcmJhYyKGSpCFGCfejwNGqurtt38Js2D92crmlPR9v+48BWwaO39z6JEkrZGi4V9WjwJEkL25dO4H7gYPA7ta3G7i1tQ8CV7e7ZnYAJwaWbyRJK2DdiOP+DPhQknOAh4BrmP2L4eYke4CHgSvb2NuBy4Bp4Mk2VpK0gkYK96q6F5icZ9fOecYWcO3SypIkLYWfUJWkDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nq0EjhnuSbSb6c5N4kU63vvCSHkjzYns9t/UlyfZLpJPcl2b6cE5Ak/byFXLn/TlVdVFWTbXsfcLiqtgGH2zbApcC29tgL3DCuYiVJo1nKsswu4EBrHwCuGOi/qWbdBaxPsnEJryNJWqBRw72A/0pyT5K9rW9DVT3S2o8CG1p7E3Bk4Nijre9pkuxNMpVkamZmZhGlS5JOZd2I4367qo4l+WXgUJKvDe6sqkpSC3nhqtoP7AeYnJxc0LGSpNMb6cq9qo615+PAJ4CLgcdOLre05+Nt+DFgy8Dhm1ufJGmFDA33JM9J8ryTbeB3ga8AB4Hdbdhu4NbWPghc3e6a2QGcGFi+kSStgFGWZTYAn0hycvy/VtWnk3weuDnJHuBh4Mo2/nbgMmAaeBK4ZuxVS5JOa2i4V9VDwEvn6X8c2DlPfwHXjqU6SdKi+AlVSeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1aORwT3JWki8mua1tvzDJ3Ummk3w0yTmt/5lte7rt37pMtUuSTmEhV+5vAR4Y2H438N6qehHwHWBP698DfKf1v7eNkyStoJHCPclm4HLgA207wGuAW9qQA8AVrb2rbdP272zjJUkrZNQr978D/hL4v7Z9PvBEVT3Vto8Cm1p7E3AEoO0/0cY/TZK9SaaSTM3MzCyueknSvIaGe5LfA45X1T3jfOGq2l9Vk1U1OTExMc4fLUm/8NaNMOaVwOuSXAY8C/gl4H3A+iTr2tX5ZuBYG38M2AIcTbIOeD7w+NgrlySd0tAr96r6q6raXFVbgauAO6rqD4A7gde3YbuBW1v7YNum7b+jqmqsVUuSTmsp97m/HXhbkmlm19RvbP03Aue3/rcB+5ZWoiRpoUZZlvmpqvoM8JnWfgi4eJ4xPwTeMIbaJEmL5CdUJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUoaHhnuRZST6X5EtJvprkXa3/hUnuTjKd5KNJzmn9z2zb023/1mWegyRpjlGu3H8EvKaqXgpcBFySZAfwbuC9VfUi4DvAnjZ+D/Cd1v/eNk6StIKGhnvN+n7bPLs9CngNcEvrPwBc0dq72jZt/84kGVfBkqThRlpzT3JWknuB48Ah4BvAE1X1VBtyFNjU2puAIwBt/wng/DHWLEkaYqRwr6qfVNVFwGbgYuAlS33hJHuTTCWZmpmZWeqPkyQNWNDdMlX1BHAn8ApgfZJ1bddm4FhrHwO2ALT9zwcen+dn7a+qyaqanJiYWFz1kqR5jXK3zESS9a39bOC1wAPMhvzr27DdwK2tfbBt0/bfUVU1xpolSUOsGz6EjcCBJGcx+5fBzVV1W5L7gY8k+Vvgi8CNbfyNwL8kmQa+DVy1DHVLkk5jaLhX1X3Ay+bpf4jZ9fe5/T8E3jCW6iRJi+InVCWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6NDTck2xJcmeS+5N8NclbWv95SQ4lebA9n9v6k+T6JNNJ7kuyfbknIUl6ulGu3J8C/qKqLgR2ANcmuRDYBxyuqm3A4bYNcCmwrT32AjeMvWpJ0mkNDfeqeqSqvtDa3wMeADYBu4ADbdgB4IrW3gXcVLPuAtYn2TjuwiVJp7agNfckW4GXAXcDG6rqkbbrUWBDa28CjgwcdrT1zf1Ze5NMJZmamZlZaN2SpNMYOdyTPBf4GPDWqvru4L6qKqAW8sJVtb+qJqtqcmJiYiGHSpKGGCnck5zNbLB/qKo+3rofO7nc0p6Pt/5jwJaBwze3PknSChnlbpkANwIPVNV7BnYdBHa39m7g1oH+q9tdMzuAEwPLN5KkFbBuhDGvBN4EfDnJva3vHcB1wM1J9gAPA1e2fbcDlwHTwJPANeMsWJI03NBwr6r/BnKK3TvnGV/AtUusS5K0BH5CVZI6ZLhLUocMd0nq0ChvqKojW/d9crVLkLQCvHKXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR3yK3+1YpbydcPfvO7yMVYi9c8rd0nq0NBwT/LBJMeTfGWg77wkh5I82J7Pbf1Jcn2S6ST3Jdm+nMVLkuY3ypX7PwOXzOnbBxyuqm3A4bYNcCmwrT32AjeMp0xJ0kIMDfeq+izw7Tndu4ADrX0AuGKg/6aadRewPsnGMdUqSRrRYtfcN1TVI639KLChtTcBRwbGHW19PyfJ3iRTSaZmZmYWWYYkaT5LfkO1qgqoRRy3v6omq2pyYmJiqWVIkgYsNtwfO7nc0p6Pt/5jwJaBcZtbnyRpBS023A8Cu1t7N3DrQP/V7a6ZHcCJgeUbSdIKGfohpiQfBl4NXJDkKPBO4Drg5iR7gIeBK9vw24HLgGngSeCaZahZkjTE0HCvqjeeYtfOecYWcO1Si5IkLY2fUJWkDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHRr6lb/SmWDrvk8u+thvXnf5GCuR1gav3CWpQ165q3te9esXkVfuktQhr9yl01jKVT945a/V45W7JHXIcJekDi3LskySS4D3AWcBH6iq65bjdX5RLXWpQCvHN3O1WsZ+5Z7kLOAfgUuBC4E3Jrlw3K8jSTq15bhyvxiYrqqHAJJ8BNgF3L8MryV1a63+huZvHGeG5Qj3TcCRge2jwG/NHZRkL7C3bX4/ydcX+XoXAN9a5LFnGudy5ullHrBCc8m7l/sVAM/LSS841Y5VuxWyqvYD+5f6c5JMVdXkGEpadc7lzNPLPMC5nKmWay7LcbfMMWDLwPbm1idJWiHLEe6fB7YleWGSc4CrgIPL8DqSpFMY+7JMVT2V5E+B/2T2VsgPVtVXx/06A5a8tHMGcS5nnl7mAc7lTLUsc0lVLcfPlSStIj+hKkkdMtwlqUNrJtyTXJLk60mmk+ybZ/8zk3y07b87ydZVKHMkI8zlzUlmktzbHn+4GnUOk+SDSY4n+cop9ifJ9W2e9yXZvtI1jmqEubw6yYmBc/LXK13jKJJsSXJnkvuTfDXJW+YZsybOy4hzWSvn5VlJPpfkS20u75pnzHgzrKrO+Aezb8x+A/g14BzgS8CFc8b8CfD+1r4K+Ohq172EubwZ+IfVrnWEubwK2A585RT7LwM+BQTYAdy92jUvYS6vBm5b7TpHmMdGYHtrPw/4n3n+fK2J8zLiXNbKeQnw3NY+G7gb2DFnzFgzbK1cuf/0Kw2q6sfAya80GLQLONDatwA7k2QFaxzVKHNZE6rqs8C3TzNkF3BTzboLWJ9k48pUtzAjzGVNqKpHquoLrf094AFmPzU+aE2clxHnsia0/9bfb5tnt8fcu1nGmmFrJdzn+0qDuSf5p2Oq6ingBHD+ilS3MKPMBeD326/MtyTZMs/+tWDUua4Vr2i/Vn8qya+vdjHDtF/rX8bsVeKgNXdeTjMXWCPnJclZSe4FjgOHquqU52UcGbZWwv0XzX8AW6vqN4FD/Oxvc62eLwAvqKqXAn8P/PvqlnN6SZ4LfAx4a1V9d7XrWYohc1kz56WqflJVFzH7qf2Lk/zGcr7eWgn3Ub7S4KdjkqwDng88viLVLczQuVTV41X1o7b5AeDlK1TbuHXzVRRV9d2Tv1ZX1e3A2UkuWOWy5pXkbGbD8ENV9fF5hqyZ8zJsLmvpvJxUVU8AdwKXzNk11gxbK+E+ylcaHAR2t/brgTuqvTNxhhk6lznrn69jdq1xLToIXN3uztgBnKiqR1a7qMVI8isn1z+TXMzs/ztn3MVDq/FG4IGqes8phq2J8zLKXNbQeZlIsr61nw28FvjanGFjzbA18Q9k1ym+0iDJ3wBTVXWQ2T8E/5Jkmtk3xq5avYpPbcS5/HmS1wFPMTuXN69awaeR5MPM3q1wQZKjwDuZfaOIqno/cDuzd2ZMA08C16xOpcONMJfXA3+c5Cngf4GrztCLh1cCbwK+3NZ3Ad4B/CqsufMyylzWynnZCBzI7D9m9Azg5qq6bTkzzK8fkKQOrZVlGUnSAhjuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUP/D2YsydU15EV0AAAAAElFTkSuQmCC\n",
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGsCAYAAAAPJKchAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAfIklEQVR4nO3de2yV533A8d8BgkkX7ISm2AYcoKGhuRpKuZisBTYaj7Io1rSMRVphbcnWCqpQtlWwS6K0m5woIWHqWGiUBSuNGCnNAA3aBAoFlOCo4qZClrKSC5AEO4mS2MHtHGa/+yOKOzeY+DiYx5fPR3r/8OvnPec5j46Ov3rPe45zWZZlAQCQyIDUEwAA+jcxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJNWrYmT37t1x4403xogRIyKXy8XGjRvzvo0sy+Lee++NK664IgoKCmLkyJHxT//0T+d+sgBApwxKPYF8NDU1RXl5eXzlK1+JP/qjP+rSbdx2222xdevWuPfee+Paa6+NN998M958881zPFMAoLNyvfUf5eVyudiwYUNUVVW17Wtubo6/+7u/i3//93+Pt99+O6655pq4++67Y+bMmRER8dxzz8V1110Xhw8fjvHjx6eZOADQTq96m+bDLF68OGpra2PdunXx85//PG6++eb4gz/4g/jlL38ZERH/+Z//GZ/85Cdj8+bNMXbs2BgzZkwsXLjQmREASKjPxMjx48djzZo1sX79+vjc5z4Xl19+efz1X/91/O7v/m6sWbMmIiJeeOGFOHbsWKxfvz4eeeSRqKmpiX379sUf//EfJ549APRfveqakbM5dOhQtLS0xBVXXNFuf3Nzc3z84x+PiIjW1tZobm6ORx55pG3cv/3bv8WkSZPiyJEj3roBgAT6TIycOnUqBg4cGPv27YuBAwe2+91FF10UERGlpaUxaNCgdsFy5ZVXRsR7Z1bECACcf30mRiZOnBgtLS3x2muvxec+97kzjrn++uvjf//3f+P555+Pyy+/PCIi/vu//zsiIkaPHn3e5goA/Eav+jTNqVOn4ujRoxHxXnzcd999MWvWrBg2bFhcdtll8Wd/9mfx9NNPx4oVK2LixInx+uuvx/bt2+O6666LuXPnRmtra0yePDkuuuiiWLlyZbS2tsaiRYuisLAwtm7dmvjRAUD/1KtiZOfOnTFr1qwP7F+wYEHU1NTE6dOn4x//8R/jkUceiVdeeSUuvfTSmDZtWtx5551x7bXXRkTEq6++Gt/4xjdi69at8Tu/8zsxZ86cWLFiRQwbNux8PxwAIHpZjAAAfU+f+WgvANA7iREAIKle8Wma1tbWePXVV2Po0KGRy+VSTwcA6IQsy+Kdd96JESNGxIABHZ//6BUx8uqrr0ZZWVnqaQAAXXDixIkYNWpUh7/vFTEydOjQiHjvwRQWFiaeDQDQGY2NjVFWVtb2d7wjvSJG3n9rprCwUIwAQC/zYZdYuIAVAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJDUoNQTgJ5szLIt3XK7L901t1tuF6A3cmYEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApPKKkerq6pg8eXIMHTo0hg8fHlVVVXHkyJGzHlNTUxO5XK7dNmTIkI80aQCg78grRnbt2hWLFi2KZ555JrZt2xanT5+OG264IZqams56XGFhYZw8ebJtO3bs2EeaNADQdwzKZ/ATTzzR7ueampoYPnx47Nu3Lz7/+c93eFwul4uSkpKuzRAA6NM+0jUjDQ0NERExbNiws447depUjB49OsrKyuKmm26KZ5999qzjm5ubo7Gxsd0GAPRNXY6R1tbWWLJkSVx//fVxzTXXdDhu/Pjx8fDDD8emTZvi0UcfjdbW1pg+fXq8/PLLHR5TXV0dRUVFbVtZWVlXpwkA9HC5LMuyrhz49a9/PX784x/HU089FaNGjer0cadPn44rr7wybrnllvjOd75zxjHNzc3R3Nzc9nNjY2OUlZVFQ0NDFBYWdmW60CVjlm3pltt96a653XK7AD1JY2NjFBUVfejf77yuGXnf4sWLY/PmzbF79+68QiQi4oILLoiJEyfG0aNHOxxTUFAQBQUFXZkaANDL5PU2TZZlsXjx4tiwYUPs2LEjxo4dm/cdtrS0xKFDh6K0tDTvYwGAvievMyOLFi2KtWvXxqZNm2Lo0KFRV1cXERFFRUVx4YUXRkTE/PnzY+TIkVFdXR0REd/+9rdj2rRpMW7cuHj77bfjnnvuiWPHjsXChQvP8UMBAHqjvGLkgQceiIiImTNnttu/Zs2a+PM///OIiDh+/HgMGPCbEy5vvfVW3HrrrVFXVxeXXHJJTJo0Kfbs2RNXXXXVR5s5ANAndPkC1vOpsxfAwLnmAlaAruvs32//mwYASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUnnFSHV1dUyePDmGDh0aw4cPj6qqqjhy5MiHHrd+/fr49Kc/HUOGDIlrr702fvSjH3V5wgBA35JXjOzatSsWLVoUzzzzTGzbti1Onz4dN9xwQzQ1NXV4zJ49e+KWW26Jr371q3HgwIGoqqqKqqqqOHz48EeePADQ++WyLMu6evDrr78ew4cPj127dsXnP//5M46ZN29eNDU1xebNm9v2TZs2LSZMmBCrV6/u1P00NjZGUVFRNDQ0RGFhYVenC3kbs2xLt9zuS3fN7ZbbBehJOvv3+yNdM9LQ0BAREcOGDetwTG1tbcyePbvdvsrKyqitre3wmObm5mhsbGy3AQB9U5djpLW1NZYsWRLXX399XHPNNR2Oq6uri+Li4nb7iouLo66ursNjqquro6ioqG0rKyvr6jQBgB6uyzGyaNGiOHz4cKxbt+5cziciIpYvXx4NDQ1t24kTJ875fQAAPcOgrhy0ePHi2Lx5c+zevTtGjRp11rElJSVRX1/fbl99fX2UlJR0eExBQUEUFBR0ZWoAQC+T15mRLMti8eLFsWHDhtixY0eMHTv2Q4+pqKiI7du3t9u3bdu2qKioyG+mAECflNeZkUWLFsXatWtj06ZNMXTo0LbrPoqKiuLCCy+MiIj58+fHyJEjo7q6OiIibrvttpgxY0asWLEi5s6dG+vWrYu9e/fGgw8+eI4fCgDQG+V1ZuSBBx6IhoaGmDlzZpSWlrZtjz32WNuY48ePx8mTJ9t+nj59eqxduzYefPDBKC8vjx/+8IexcePGs170CgD0H3mdGenMV5Ls3LnzA/tuvvnmuPnmm/O5KwCgn/C/aQCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAklXeM7N69O2688cYYMWJE5HK52Lhx41nH79y5M3K53Ae2urq6rs4ZAOhD8o6RpqamKC8vj1WrVuV13JEjR+LkyZNt2/Dhw/O9awCgDxqU7wFz5syJOXPm5H1Hw4cPj4svvjjv4wCAvu28XTMyYcKEKC0tjS984Qvx9NNPn3Vsc3NzNDY2ttsAgL6p22OktLQ0Vq9eHY8//ng8/vjjUVZWFjNnzoz9+/d3eEx1dXUUFRW1bWVlZd09TQAgkVyWZVmXD87lYsOGDVFVVZXXcTNmzIjLLrssvv/975/x983NzdHc3Nz2c2NjY5SVlUVDQ0MUFhZ2dbqQtzHLtnTL7b5019xuuV2AnqSxsTGKioo+9O933teMnAtTpkyJp556qsPfFxQUREFBwXmcEQCQSpLvGTl48GCUlpamuGsAoIfJ+8zIqVOn4ujRo20/v/jii3Hw4MEYNmxYXHbZZbF8+fJ45ZVX4pFHHomIiJUrV8bYsWPj6quvjv/5n/+Jhx56KHbs2BFbt249d48CAOi18o6RvXv3xqxZs9p+Xrp0aURELFiwIGpqauLkyZNx/Pjxtt+/++678Vd/9VfxyiuvxMc+9rG47rrr4ic/+Um72wAA+q+PdAHr+dLZC2DgXHMBK0DXdfbvt/9NAwAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASeUdI7t3744bb7wxRowYEblcLjZu3Pihx+zcuTM+85nPREFBQYwbNy5qamq6MFUAoC/KO0aampqivLw8Vq1a1anxL774YsydOzdmzZoVBw8ejCVLlsTChQvjySefzHuyAEDfMyjfA+bMmRNz5szp9PjVq1fH2LFjY8WKFRERceWVV8ZTTz0V999/f1RWVuZ79wBAH9Pt14zU1tbG7Nmz2+2rrKyM2traDo9pbm6OxsbGdhsA0Dd1e4zU1dVFcXFxu33FxcXR2NgYv/71r894THV1dRQVFbVtZWVl3T1NACCRHvlpmuXLl0dDQ0PbduLEidRTAgC6Sd7XjOSrpKQk6uvr2+2rr6+PwsLCuPDCC894TEFBQRQUFHT31ACAHqDbz4xUVFTE9u3b2+3btm1bVFRUdPddAwC9QN4xcurUqTh48GAcPHgwIt776O7Bgwfj+PHjEfHeWyzz589vG/+1r30tXnjhhfjWt74Vv/jFL+Jf//Vf4wc/+EF885vfPDePAADo1fKOkb1798bEiRNj4sSJERGxdOnSmDhxYtx+++0REXHy5Mm2MImIGDt2bGzZsiW2bdsW5eXlsWLFinjooYd8rBcAiIiIXJZlWepJfJjGxsYoKiqKhoaGKCwsTD0d+pExy7Z0y+2+dNfcbrldgJ6ks3+/e+SnaQCA/kOMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUoNSTwA+qjHLtqSeAgAfgTMjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJBUl2Jk1apVMWbMmBgyZEhMnTo1fvazn3U4tqamJnK5XLttyJAhXZ4wANC35B0jjz32WCxdujTuuOOO2L9/f5SXl0dlZWW89tprHR5TWFgYJ0+ebNuOHTv2kSYNAPQdg/I94L777otbb701vvzlL0dExOrVq2PLli3x8MMPx7Jly854TC6Xi5KSko82U+hDxizb0m23/dJdc7vttgG6Q15nRt59993Yt29fzJ49+zc3MGBAzJ49O2prazs87tSpUzF69OgoKyuLm266KZ599tmz3k9zc3M0Nja22wCAvimvGHnjjTeipaUliouL2+0vLi6Ourq6Mx4zfvz4ePjhh2PTpk3x6KOPRmtra0yfPj1efvnlDu+nuro6ioqK2raysrJ8pgkA9CLd/mmaioqKmD9/fkyYMCFmzJgR//Ef/xGf+MQn4nvf+16HxyxfvjwaGhrathMnTnT3NAGARPK6ZuTSSy+NgQMHRn19fbv99fX1nb4m5IILLoiJEyfG0aNHOxxTUFAQBQUF+UwNAOil8jozMnjw4Jg0aVJs3769bV9ra2ts3749KioqOnUbLS0tcejQoSgtLc1vpgBAn5T3p2mWLl0aCxYsiM9+9rMxZcqUWLlyZTQ1NbV9umb+/PkxcuTIqK6ujoiIb3/72zFt2rQYN25cvP3223HPPffEsWPHYuHChef2kQAAvVLeMTJv3rx4/fXX4/bbb4+6urqYMGFCPPHEE20XtR4/fjwGDPjNCZe33norbr311qirq4tLLrkkJk2aFHv27Imrrrrq3D0KAKDXymVZlqWexIdpbGyMoqKiaGhoiMLCwtTToYfpzu/s6I18zwjQU3T277f/TQMAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSg1JPADi3xizb0i23+9Jdc7vldgGcGQEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUr4OHuiU7vqa+QhfNQ/9nTMjAEBSYgQASEqMAABJuWYESK67rkdxLQr0Ds6MAABJOTMC9Fk+AQS9gzMjAEBSXTozsmrVqrjnnnuirq4uysvL47vf/W5MmTKlw/Hr16+Pf/iHf4iXXnopPvWpT8Xdd98dX/ziF7s8aYDUXOcC507eMfLYY4/F0qVLY/Xq1TF16tRYuXJlVFZWxpEjR2L48OEfGL9nz5645ZZborq6Ov7wD/8w1q5dG1VVVbF///645pprzsmDoOfrztPl0Jd4a4n+KJdlWZbPAVOnTo3JkyfHv/zLv0RERGtra5SVlcU3vvGNWLZs2QfGz5s3L5qammLz5s1t+6ZNmxYTJkyI1atXd+o+Gxsbo6ioKBoaGqKwsDCf6dJDiBFIT4xwvnX273deZ0befffd2LdvXyxfvrxt34ABA2L27NlRW1t7xmNqa2tj6dKl7fZVVlbGxo0bO7yf5ubmaG5ubvu5oaEhIt57UPROrc2/Sj0F6Pcu++b61FPI2+E7K1NPgY/g/b/bH3beI68YeeONN6KlpSWKi4vb7S8uLo5f/OIXZzymrq7ujOPr6uo6vJ/q6uq48847P7C/rKwsn+kC0MsVrUw9A86Fd955J4qKijr8fY/8aO/y5cvbnU1pbW2NN998Mz7+8Y9HLpc7Z/fT2NgYZWVlceLECW//fAhrlR/r1XnWqvOsVedZq87rzrXKsizeeeedGDFixFnH5RUjl156aQwcODDq6+vb7a+vr4+SkpIzHlNSUpLX+IiIgoKCKCgoaLfv4osvzmeqeSksLPRk7SRrlR/r1XnWqvOsVedZq87rrrU62xmR9+X1PSODBw+OSZMmxfbt29v2tba2xvbt26OiouKMx1RUVLQbHxGxbdu2DscDAP1L3m/TLF26NBYsWBCf/exnY8qUKbFy5cpoamqKL3/5yxERMX/+/Bg5cmRUV1dHRMRtt90WM2bMiBUrVsTcuXNj3bp1sXfv3njwwQfP7SMBAHqlvGNk3rx58frrr8ftt98edXV1MWHChHjiiSfaLlI9fvx4DBjwmxMu06dPj7Vr18bf//3fx9/+7d/Gpz71qdi4cWOP+I6RgoKCuOOOOz7wlhAfZK3yY706z1p1nrXqPGvVeT1hrfL+nhEAgHPJ/6YBAJISIwBAUmIEAEhKjAAASfX5GFm1alWMGTMmhgwZElOnTo2f/exnZx2/fv36+PSnPx1DhgyJa6+9Nn70ox+dp5mml89a1dTURC6Xa7cNGTLkPM42nd27d8eNN94YI0aMiFwud9b/s/S+nTt3xmc+85koKCiIcePGRU1NTbfPsyfId6127tz5gedVLpc767+P6Cuqq6tj8uTJMXTo0Bg+fHhUVVXFkSNHPvS4/via1ZW16q+vWQ888EBcd911bV9oVlFRET/+8Y/PekyK51SfjpHHHnssli5dGnfccUfs378/ysvLo7KyMl577bUzjt+zZ0/ccsst8dWvfjUOHDgQVVVVUVVVFYcPHz7PMz//8l2riPe+re/kyZNt27Fjx87jjNNpamqK8vLyWLVqVafGv/jiizF37tyYNWtWHDx4MJYsWRILFy6MJ598sptnml6+a/W+I0eOtHtuDR8+vJtm2HPs2rUrFi1aFM8880xs27YtTp8+HTfccEM0NTV1eEx/fc3qylpF9M/XrFGjRsVdd90V+/bti71798bv/d7vxU033RTPPvvsGccne05lfdiUKVOyRYsWtf3c0tKSjRgxIquurj7j+D/5kz/J5s6d227f1KlTs7/8y7/s1nn2BPmu1Zo1a7KioqLzNLueKyKyDRs2nHXMt771rezqq69ut2/evHlZZWVlN86s5+nMWv30pz/NIiJ76623zsucerLXXnsti4hs165dHY7pz69Z/19n1spr1m9ccskl2UMPPXTG36V6TvXZMyPvvvtu7Nu3L2bPnt22b8CAATF79uyora094zG1tbXtxkdEVFZWdji+r+jKWkVEnDp1KkaPHh1lZWVnLe3+rr8+rz6KCRMmRGlpaXzhC1+Ip59+OvV0kmhoaIiIiGHDhnU4xnPrPZ1ZqwivWS0tLbFu3bpoamrq8F+ypHpO9dkYeeONN6KlpaXtm2HfV1xc3OH7z3V1dXmN7yu6slbjx4+Phx9+ODZt2hSPPvpotLa2xvTp0+Pll18+H1PuVTp6XjU2Nsavf/3rRLPqmUpLS2P16tXx+OOPx+OPPx5lZWUxc+bM2L9/f+qpnVetra2xZMmSuP7668/6bdX99TXr/+vsWvXn16xDhw7FRRddFAUFBfG1r30tNmzYEFddddUZx6Z6TuX9dfAQ8d4/QPz/ZT19+vS48sor43vf+1585zvfSTgzerPx48fH+PHj236ePn16PP/883H//ffH97///YQzO78WLVoUhw8fjqeeeir1VHq8zq5Vf37NGj9+fBw8eDAaGhrihz/8YSxYsCB27drVYZCk0GfPjFx66aUxcODAqK+vb7e/vr4+SkpKznhMSUlJXuP7iq6s1W+74IILYuLEiXH06NHumGKv1tHzqrCwMC688MJEs+o9pkyZ0q+eV4sXL47NmzfHT3/60xg1atRZx/bX16z35bNWv60/vWYNHjw4xo0bF5MmTYrq6uooLy+Pf/7nfz7j2FTPqT4bI4MHD45JkybF9u3b2/a1trbG9u3bO3yvrKKiot34iIht27Z1OL6v6Mpa/baWlpY4dOhQlJaWdtc0e63++rw6Vw4ePNgvnldZlsXixYtjw4YNsWPHjhg7duyHHtNfn1tdWavf1p9fs1pbW6O5ufmMv0v2nOrWy2MTW7duXVZQUJDV1NRk//Vf/5X9xV/8RXbxxRdndXV1WZZl2Ze+9KVs2bJlbeOffvrpbNCgQdm9996bPffcc9kdd9yRXXDBBdmhQ4dSPYTzJt+1uvPOO7Mnn3wye/7557N9+/Zlf/qnf5oNGTIke/bZZ1M9hPPmnXfeyQ4cOJAdOHAgi4jsvvvuyw4cOJAdO3Ysy7IsW7ZsWfalL32pbfwLL7yQfexjH8v+5m/+JnvuueeyVatWZQMHDsyeeOKJVA/hvMl3re6///5s48aN2S9/+cvs0KFD2W233ZYNGDAg+8lPfpLqIZw3X//617OioqJs586d2cmTJ9u2X/3qV21jvGa9pytr1V9fs5YtW5bt2rUre/HFF7Of//zn2bJly7JcLpdt3bo1y7Ke85zq0zGSZVn23e9+N7vsssuywYMHZ1OmTMmeeeaZtt/NmDEjW7BgQbvxP/jBD7IrrrgiGzx4cHb11VdnW7ZsOc8zTieftVqyZEnb2OLi4uyLX/xitn///gSzPv/e//jpb2/vr8+CBQuyGTNmfOCYCRMmZIMHD84++clPZmvWrDnv804h37W6++67s8svvzwbMmRINmzYsGzmzJnZjh070kz+PDvTOkVEu+eK16z3dGWt+utr1le+8pVs9OjR2eDBg7NPfOIT2e///u+3hUiW9ZznVC7Lsqx7z70AAHSsz14zAgD0DmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgqf8DgyUK1RABovAAAAAASUVORK5CYII=\n",
       "text/plain": [
-       "<Figure size 432x288 with 1 Axes>"
+       "<Figure size 640x480 with 1 Axes>"
       ]
      },
-     "metadata": {
-      "needs_background": "light"
-     },
+     "metadata": {},
      "output_type": "display_data"
     }
    ],
@@ -105,7 +103,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 7,
+   "execution_count": 6,
    "id": "6c2a4a1a",
    "metadata": {},
    "outputs": [],
@@ -115,17 +113,17 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 8,
+   "execution_count": 7,
    "id": "a7231e19",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "text/plain": [
-       "<Array [1.07, 1.04, 1.04, ... 1.21, 1.51, 1.19] type='1098 * float32'>"
+       "<Array [1.07, 1.04, 1.04, ... 0.995, 0.601] type='4343693 * float32'>"
       ]
      },
-     "execution_count": 8,
+     "execution_count": 7,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -137,20 +135,18 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 9,
+   "execution_count": 8,
    "id": "abfe684d",
    "metadata": {},
    "outputs": [
     {
      "data": {
-      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAOzUlEQVR4nO3db4xcV33G8e9DnAAqFEO8dSPbYVMRqUqrAqkVjKiqlIgqfyocqQEFVcREriy1qQqiVTG8KKLqC/OGlLQVyCKoDqKQNEDjhtDWCkGoUhPYQAj5A2WJEsVWwCYkBpRCZfrriz1Oh2XXM7s7u+s5+X6k0Zx7zpmZ38l1Ht+9e+c6VYUkqS/PW+8CJEnjZ7hLUocMd0nqkOEuSR0y3CWpQxvWuwCATZs21fT09HqXIUkT5d577/1eVU0tNHZahPv09DQzMzPrXYYkTZQkjy025mkZSeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUodGCvckjyb5epL7ksy0vpclOZTkW+35pa0/SW5IMpvk/iQXruYCJEk/bynfUP2dqvrewPZe4M6q2pdkb9t+F3AZcH57vAb4UHtWB6b3fnbZr3103xVjrETSqazktMxO4EBrHwCuHOi/qebcDWxMcs4KPkeStESjhnsB/57k3iR7Wt/mqnqitb8DbG7tLcDjA6893PokSWtk1NMyv1VVR5L8EnAoyTcGB6uqkizpH2Ntf0nsATj33HOX8lJJ0hAjHblX1ZH2fBT4DHAR8N2Tp1va89E2/QiwbeDlW1vf/PfcX1Xbq2r71NSCd6yUJC3T0HBP8gtJXnyyDfwu8ABwENjVpu0Cbmvtg8A17aqZHcDxgdM3kqQ1MMppmc3AZ5KcnP+PVfWvSb4M3JJkN/AY8OY2/w7gcmAWeAa4duxVS5JOaWi4V9UjwCsX6H8SuGSB/gKuG0t1kqRl8RuqktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA6NHO5Jzkjy1SS3t+3zktyTZDbJzUnOav3Pb9uzbXx6lWqXJC1iKUfubwceHth+P3B9Vb0CeArY3fp3A0+1/uvbPEnSGhop3JNsBa4APtK2A7weuLVNOQBc2do72zZt/JI2X5K0RkY9cv8b4C+A/23bZwNPV9WJtn0Y2NLaW4DHAdr48Tb/ZyTZk2QmycyxY8eWV70kaUFDwz3J7wFHq+recX5wVe2vqu1VtX1qamqcby1Jz3kbRpjzOuCNSS4HXgD8IvBBYGOSDe3ofCtwpM0/AmwDDifZALwEeHLslUuSFjU03Kvq3cC7AZJcDPx5Vf1Bkn8CrgI+CewCbmsvOdi2/7ONf76qauyVa+JM7/3ssl/76L4rxliJ1L+VXOf+LuCdSWaZO6d+Y+u/ETi79b8T2LuyEiVJSzXKaZlnVdUXgC+09iPARQvM+THwpjHUJklaJr+hKkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SerQ0HBP8oIkX0rytSQPJnlf6z8vyT1JZpPcnOSs1v/8tj3bxqdXeQ2SpHlGOXL/CfD6qnol8Crg0iQ7gPcD11fVK4CngN1t/m7gqdZ/fZsnSVpDQ8O95vyobZ7ZHgW8Hri19R8ArmztnW2bNn5JkoyrYEnScCOdc09yRpL7gKPAIeDbwNNVdaJNOQxsae0twOMAbfw4cPYC77knyUySmWPHjq1oEZKknzVSuFfVT6vqVcBW4CLgV1f6wVW1v6q2V9X2qamplb6dJGnAkq6WqaqngbuA1wIbk2xoQ1uBI619BNgG0MZfAjw5jmIlSaMZ5WqZqSQbW/uFwBuAh5kL+avatF3Aba19sG3Txj9fVTXGmiVJQ2wYPoVzgANJzmDuL4Nbqur2JA8Bn0zy18BXgRvb/BuBjyWZBb4PXL0KdUuSTmFouFfV/cCrF+h/hLnz7/P7fwy8aSzVSZKWxW+oSlKHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjo0NNyTbEtyV5KHkjyY5O2t/2VJDiX5Vnt+aetPkhuSzCa5P8mFq70ISdLPGuXI/QTwZ1V1AbADuC7JBcBe4M6qOh+4s20DXAac3x57gA+NvWpJ0ikNDfeqeqKqvtLaPwQeBrYAO4EDbdoB4MrW3gncVHPuBjYmOWfchUuSFrekc+5JpoFXA/cAm6vqiTb0HWBza28BHh942eHWN/+99iSZSTJz7NixpdYtSTqFkcM9yYuATwHvqKofDI5VVQG1lA+uqv1Vtb2qtk9NTS3lpZKkIUYK9yRnMhfsH6+qT7fu75483dKej7b+I8C2gZdvbX2SpDUyytUyAW4EHq6qDwwMHQR2tfYu4LaB/mvaVTM7gOMDp28kSWtgwwhzXge8Ffh6kvta33uAfcAtSXYDjwFvbmN3AJcDs8AzwLXjLFiSNNzQcK+q/wCyyPAlC8wv4LoV1iVJWgG/oSpJHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nq0NBwT/LRJEeTPDDQ97Ikh5J8qz2/tPUnyQ1JZpPcn+TC1SxekrSwUY7c/wG4dF7fXuDOqjofuLNtA1wGnN8ee4APjadMSdJSDA33qvoi8P153TuBA619ALhyoP+mmnM3sDHJOWOqVZI0ouWec99cVU+09neAza29BXh8YN7h1vdzkuxJMpNk5tixY8ssQ5K0kBX/QrWqCqhlvG5/VW2vqu1TU1MrLUOSNGC54f7dk6db2vPR1n8E2DYwb2vrkyStoeWG+0FgV2vvAm4b6L+mXTWzAzg+cPpGkrRGNgybkOQTwMXApiSHgfcC+4BbkuwGHgPe3KbfAVwOzALPANeuQs2SpCGGhntVvWWRoUsWmFvAdSstSpK0Mn5DVZI6ZLhLUoeGnpZRX6b3fna9S5C0Bjxyl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pB3hZxAz8U7O65kzY/uu2KMlUiTwSN3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR3y9gPr5Ll4CwFJa8cjd0nqkOEuSR0y3CWpQ4a7JHVoVX6hmuRS4IPAGcBHqmrfanyONArvBa/norEfuSc5A/h74DLgAuAtSS4Y9+dIkha3GkfuFwGzVfUIQJJPAjuBh1bhs9aVlzP2b6X72CN/rZfVCPctwOMD24eB18yflGQPsKdt/ijJN5f5eZuA7y3ztacb13L6WdE68v4xVrJyvewTcC0nvXyxgXX7ElNV7Qf2r/R9ksxU1fYxlLTuXMvpp5d1gGs5Xa3WWlbjapkjwLaB7a2tT5K0RlYj3L8MnJ/kvCRnAVcDB1fhcyRJixj7aZmqOpHkT4B/Y+5SyI9W1YPj/pwBKz61cxpxLaefXtYBruV0tSprSVWtxvtKktaR31CVpA4Z7pLUoYkJ9ySXJvlmktkkexcYf36Sm9v4PUmm16HMkYywlrclOZbkvvb4w/Woc5gkH01yNMkDi4wnyQ1tnfcnuXCtaxzVCGu5OMnxgX3yl2td4yiSbEtyV5KHkjyY5O0LzJmI/TLiWiZlv7wgyZeSfK2t5X0LzBlvhlXVaf9g7hez3wZ+BTgL+Bpwwbw5fwx8uLWvBm5e77pXsJa3AX+33rWOsJbfBi4EHlhk/HLgc0CAHcA9613zCtZyMXD7etc5wjrOAS5s7RcD/7XAn6+J2C8jrmVS9kuAF7X2mcA9wI55c8aaYZNy5P7sLQ2q6n+Ak7c0GLQTONDatwKXJMka1jiqUdYyEarqi8D3TzFlJ3BTzbkb2JjknLWpbmlGWMtEqKonquorrf1D4GHmvjU+aCL2y4hrmQjtv/WP2uaZ7TH/apaxZtikhPtCtzSYv5OfnVNVJ4DjwNlrUt3SjLIWgN9vPzLfmmTbAuOTYNS1TorXth+rP5fk19a7mGHaj/WvZu4ocdDE7ZdTrAUmZL8kOSPJfcBR4FBVLbpfxpFhkxLuzzX/AkxX1W8Ah/j/v821fr4CvLyqXgn8LfDP61vOqSV5EfAp4B1V9YP1rmclhqxlYvZLVf20ql7F3Lf2L0ry66v5eZMS7qPc0uDZOUk2AC8BnlyT6pZm6Fqq6smq+knb/Ajwm2tU27h1cyuKqvrByR+rq+oO4Mwkm9a5rAUlOZO5MPx4VX16gSkTs1+GrWWS9stJVfU0cBdw6byhsWbYpIT7KLc0OAjsau2rgM9X+83EaWboWuad/3wjc+caJ9FB4Jp2dcYO4HhVPbHeRS1Hkl8+ef4zyUXM/b9z2h08tBpvBB6uqg8sMm0i9ssoa5mg/TKVZGNrvxB4A/CNedPGmmHrdlfIpahFbmmQ5K+Amao6yNwfgo8lmWXuF2NXr1/FixtxLX+a5I3ACebW8rZ1K/gUknyCuasVNiU5DLyXuV8UUVUfBu5g7sqMWeAZ4Nr1qXS4EdZyFfBHSU4A/w1cfZoePLwOeCvw9XZ+F+A9wLkwcftllLVMyn45BziQuX/M6HnALVV1+2pmmLcfkKQOTcppGUnSEhjuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUP/B+tC2+0i88UDAAAAAElFTkSuQmCC\n",
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiwAAAGsCAYAAAD+L/ysAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAArDUlEQVR4nO3de3BUZZ7/8U8nmA6oaUDITSMJF8OIkDAoMag/YG0N2SxldmoUqBmJGdEdFi2YeJnEVdDV3aAiojNZMioQmFkFWQVrhIkw0UAhQYpASnCQBQw3SYeL0k16xuCkz+8Py3ZbEsgJJv3Qeb+qTmk/53uefM+ppvtTp0+fdliWZQkAAMBgUeFuAAAA4HwILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeBEXWDZu3KiJEycqOTlZDodDq1evtj2HZVmaN2+errnmGjmdTl155ZX6j//4jx++WQAA0C49wt3AD83v9ysjI0O/+MUv9JOf/KRDc8ycOVPr1q3TvHnzNHz4cH3xxRf64osvfuBOAQBAezki+ccPHQ6HVq1apfz8/OBYc3Oz/u3f/k1vvPGGTp06peuuu07PPvusxo0bJ0navXu3RowYoV27dik9PT08jQMAgBAR95HQ+TzwwAOqqanR8uXL9fHHH+vOO+/UhAkTtHfvXknSH//4Rw0cOFDvvvuu0tLSlJqaqmnTpnGGBQCAMOpWgeXQoUNasmSJVq5cqVtuuUWDBg3Sww8/rJtvvllLliyRJH322Wc6ePCgVq5cqWXLlqmiokK1tbX66U9/GubuAQDoviLuGpZz2blzp1paWnTNNdeEjDc3N+uKK66QJAUCATU3N2vZsmXBukWLFmnUqFHas2cPHxMBABAG3SqwNDU1KTo6WrW1tYqOjg5Zd9lll0mSkpKS1KNHj5BQ86Mf/UjSN2doCCwAAHS9bhVYRo4cqZaWFh07dky33HJLqzU33XST/v73v2v//v0aNGiQJOl///d/JUkDBgzosl4BAMB3Iu5bQk1NTdq3b5+kbwLK/PnzNX78ePXt21dXX321fv7zn+vDDz/UCy+8oJEjR+r48eOqqqrSiBEjlJeXp0AgoBtuuEGXXXaZFixYoEAgoBkzZiguLk7r1q0L894BANA9RVxgqa6u1vjx488aLygoUEVFhb7++ms988wzWrZsmT7//HP169dPN954o5566ikNHz5cknT06FE9+OCDWrdunS699FLl5ubqhRdeUN++fbt6dwAAgCIwsAAAgMjTrb7WDAAALk4EFgAAYLyI+JZQIBDQ0aNHdfnll8vhcIS7HQAA0A6WZen06dNKTk5WVNS5z6FERGA5evSoUlJSwt0GAADogMOHD+uqq646Z01EBJbLL79c0jc7HBcXF+ZuAABAe/h8PqWkpATfx88lIgLLtx8DxcXFEVgAALjItOdyDi66BQAAxiOwAAAA4xFYAACA8WwFltLSUt1www26/PLLFR8fr/z8fO3Zs+e8261cuVJDhw5VbGyshg8frrVr14astyxLs2fPVlJSknr27Cm32629e/fa2xMAABCxbAWWDRs2aMaMGdqyZYvWr1+vr7/+Wrfffrv8fn+b22zevFlTpkzRvffeqx07dig/P1/5+fnatWtXsOa5557Tyy+/rPLycn300Ue69NJLlZOTo6+++qrjewYAACLGBf2W0PHjxxUfH68NGzbo//2//9dqzaRJk+T3+/Xuu+8Gx2688UZlZmaqvLxclmUpOTlZDz30kB5++GFJktfrVUJCgioqKjR58uTz9uHz+eRyueT1evmWEAAAFwk7798XdA2L1+uVpHP+inFNTY3cbnfIWE5OjmpqaiRJ9fX18ng8ITUul0tZWVnBmu9rbm6Wz+cLWQAAQOTqcGAJBAKaNWuWbrrpJl133XVt1nk8HiUkJISMJSQkyOPxBNd/O9ZWzfeVlpbK5XIFF+5yCwBAZOtwYJkxY4Z27dql5cuX/5D9tEtJSYm8Xm9wOXz4cJf3AAAAuk6H7nT7wAMP6N1339XGjRvPe+//xMRENTY2how1NjYqMTExuP7bsaSkpJCazMzMVud0Op1yOp0daR0AAFyEbJ1hsSxLDzzwgFatWqX3339faWlp590mOztbVVVVIWPr169Xdna2JCktLU2JiYkhNT6fTx999FGwBgAAdG+2zrDMmDFDr7/+ut555x1dfvnlwWtMXC6XevbsKUmaOnWqrrzySpWWlkqSZs6cqbFjx+qFF15QXl6eli9frm3btumVV16R9M3vB8yaNUvPPPOMhgwZorS0ND3xxBNKTk5Wfn7+D7irAADgYmUrsCxcuFCSNG7cuJDxJUuW6J577pEkHTp0SFFR3524GTNmjF5//XU9/vjjeuyxxzRkyBCtXr065ELdRx99VH6/X/fff79OnTqlm2++WZWVlYqNje3gbgEAgEhyQfdhMQX3YQEA4OJj5/27QxfdAvhGavGaTpv7wNy8TpsbAC42/PghAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADCe7cCyceNGTZw4UcnJyXI4HFq9evU56++55x45HI6zlmHDhgVrnnzyybPWDx061PbOAACAyGQ7sPj9fmVkZKisrKxd9S+99JIaGhqCy+HDh9W3b1/deeedIXXDhg0Lqdu0aZPd1gAAQITqYXeD3Nxc5ebmtrve5XLJ5XIFH69evVpffvmlCgsLQxvp0UOJiYl22wEAAN1Al1/DsmjRIrndbg0YMCBkfO/evUpOTtbAgQP1s5/9TIcOHWpzjubmZvl8vpAFAABEri4NLEePHtWf/vQnTZs2LWQ8KytLFRUVqqys1MKFC1VfX69bbrlFp0+fbnWe0tLS4Jkbl8ullJSUrmgfAACESZcGlqVLl6p3797Kz88PGc/NzdWdd96pESNGKCcnR2vXrtWpU6f05ptvtjpPSUmJvF5vcDl8+HAXdA8AAMLF9jUsHWVZlhYvXqy7775bMTEx56zt3bu3rrnmGu3bt6/V9U6nU06nszPaBAAABuqyMywbNmzQvn37dO+99563tqmpSfv371dSUlIXdAYAAExnO7A0NTWprq5OdXV1kqT6+nrV1dUFL5ItKSnR1KlTz9pu0aJFysrK0nXXXXfWuocfflgbNmzQgQMHtHnzZv3zP/+zoqOjNWXKFLvtAQCACGT7I6Ft27Zp/PjxwcdFRUWSpIKCAlVUVKihoeGsb/h4vV699dZbeumll1qd88iRI5oyZYpOnjyp/v376+abb9aWLVvUv39/u+0BAIAI5LAsywp3ExfK5/PJ5XLJ6/UqLi4u3O2gG0ktXtNpcx+Ym9dpcwOACey8f/NbQgAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgvB52N9i4caOef/551dbWqqGhQatWrVJ+fn6b9dXV1Ro/fvxZ4w0NDUpMTAw+Lisr0/PPPy+Px6OMjAz95je/0ejRo+22B0SM1OI1nTLvgbl5nTIvAHQm22dY/H6/MjIyVFZWZmu7PXv2qKGhIbjEx8cH161YsUJFRUWaM2eOtm/froyMDOXk5OjYsWN22wMAABHI9hmW3Nxc5ebm2v5D8fHx6t27d6vr5s+fr/vuu0+FhYWSpPLycq1Zs0aLFy9WcXGx7b8FAAAiS5ddw5KZmamkpCTddttt+vDDD4PjZ86cUW1trdxu93dNRUXJ7Xarpqam1bmam5vl8/lCFgAAELk6PbAkJSWpvLxcb731lt566y2lpKRo3Lhx2r59uyTpxIkTamlpUUJCQsh2CQkJ8ng8rc5ZWloql8sVXFJSUjp7NwAAQBjZ/kjIrvT0dKWnpwcfjxkzRvv379eLL76o3//+9x2as6SkREVFRcHHPp+P0AIAQATr9MDSmtGjR2vTpk2SpH79+ik6OlqNjY0hNY2NjSHfIvq/nE6nnE5np/cJAADMEJb7sNTV1SkpKUmSFBMTo1GjRqmqqiq4PhAIqKqqStnZ2eFoDwAAGMb2GZampibt27cv+Li+vl51dXXq27evrr76apWUlOjzzz/XsmXLJEkLFixQWlqahg0bpq+++kqvvfaa3n//fa1bty44R1FRkQoKCnT99ddr9OjRWrBggfx+f/BbQwAAoHuzHVi2bdsWciO4b68lKSgoUEVFhRoaGnTo0KHg+jNnzuihhx7S559/rl69emnEiBH685//HDLHpEmTdPz4cc2ePVsej0eZmZmqrKw860JcAADQPTksy7LC3cSF8vl8crlc8nq9iouLC3c76EY66260nYk73QIwhZ33b35LCAAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMZzuwbNy4URMnTlRycrIcDodWr159zvq3335bt912m/r376+4uDhlZ2frvffeC6l58skn5XA4QpahQ4fabQ0AAEQo24HF7/crIyNDZWVl7arfuHGjbrvtNq1du1a1tbUaP368Jk6cqB07doTUDRs2TA0NDcFl06ZNdlsDAAARqofdDXJzc5Wbm9vu+gULFoQ8/s///E+98847+uMf/6iRI0d+10iPHkpMTLTbDgAA6Aa6/BqWQCCg06dPq2/fviHje/fuVXJysgYOHKif/exnOnToUJtzNDc3y+fzhSwAACBydXlgmTdvnpqamnTXXXcFx7KyslRRUaHKykotXLhQ9fX1uuWWW3T69OlW5ygtLZXL5QouKSkpXdU+AAAIgy4NLK+//rqeeuopvfnmm4qPjw+O5+bm6s4779SIESOUk5OjtWvX6tSpU3rzzTdbnaekpERerze4HD58uKt2AQAAhIHta1g6avny5Zo2bZpWrlwpt9t9ztrevXvrmmuu0b59+1pd73Q65XQ6O6NNAABgoC45w/LGG2+osLBQb7zxhvLy8s5b39TUpP379yspKakLugMAAKazfYalqakp5MxHfX296urq1LdvX1199dUqKSnR559/rmXLlkn65mOggoICvfTSS8rKypLH45Ek9ezZUy6XS5L08MMPa+LEiRowYICOHj2qOXPmKDo6WlOmTPkh9hEAAFzkbJ9h2bZtm0aOHBn8SnJRUZFGjhyp2bNnS5IaGhpCvuHzyiuv6O9//7tmzJihpKSk4DJz5sxgzZEjRzRlyhSlp6frrrvu0hVXXKEtW7aof//+F7p/AAAgAjgsy7LC3cSF8vl8crlc8nq9iouLC3c76EZSi9eEuwXbDsw9/8eyANAV7Lx/81tCAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA824Fl48aNmjhxopKTk+VwOLR69erzblNdXa0f//jHcjqdGjx4sCoqKs6qKSsrU2pqqmJjY5WVlaWtW7fabQ0AAEQo24HF7/crIyNDZWVl7aqvr69XXl6exo8fr7q6Os2aNUvTpk3Te++9F6xZsWKFioqKNGfOHG3fvl0ZGRnKycnRsWPH7LYHAAAikMOyLKvDGzscWrVqlfLz89us+fWvf601a9Zo165dwbHJkyfr1KlTqqyslCRlZWXphhtu0G9/+1tJUiAQUEpKih588EEVFxeftw+fzyeXyyWv16u4uLiO7g5gW2rxmnC3YNuBuXnhbgEAJNl7/+70a1hqamrkdrtDxnJyclRTUyNJOnPmjGpra0NqoqKi5Ha7gzXf19zcLJ/PF7IAAIDI1emBxePxKCEhIWQsISFBPp9Pf/vb33TixAm1tLS0WuPxeFqds7S0VC6XK7ikpKR0Wv8AACD8LspvCZWUlMjr9QaXw4cPh7slAADQiXp09h9ITExUY2NjyFhjY6Pi4uLUs2dPRUdHKzo6utWaxMTEVud0Op1yOp2d1jMAADBLp59hyc7OVlVVVcjY+vXrlZ2dLUmKiYnRqFGjQmoCgYCqqqqCNQAAoHuzHViamppUV1enuro6Sd98bbmurk6HDh2S9M3HNVOnTg3W//KXv9Rnn32mRx99VJ9++qn+67/+S2+++aZ+9atfBWuKior06quvaunSpdq9e7emT58uv9+vwsLCC9w9AAAQCWx/JLRt2zaNHz8++LioqEiSVFBQoIqKCjU0NATDiySlpaVpzZo1+tWvfqWXXnpJV111lV577TXl5OQEayZNmqTjx49r9uzZ8ng8yszMVGVl5VkX4gIAgO7pgu7DYgruw4Jw4T4sANBxRt2HBQAA4EIRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA43UosJSVlSk1NVWxsbHKysrS1q1b26wdN26cHA7HWUteXl6w5p577jlr/YQJEzrSGgAAiEA97G6wYsUKFRUVqby8XFlZWVqwYIFycnK0Z88excfHn1X/9ttv68yZM8HHJ0+eVEZGhu68886QugkTJmjJkiXBx06n025rAAAgQtk+wzJ//nzdd999Kiws1LXXXqvy8nL16tVLixcvbrW+b9++SkxMDC7r169Xr169zgosTqczpK5Pnz4d2yMAABBxbAWWM2fOqLa2Vm63+7sJoqLkdrtVU1PTrjkWLVqkyZMn69JLLw0Zr66uVnx8vNLT0zV9+nSdPHnSTmsAACCC2fpI6MSJE2ppaVFCQkLIeEJCgj799NPzbr9161bt2rVLixYtChmfMGGCfvKTnygtLU379+/XY489ptzcXNXU1Cg6OvqseZqbm9Xc3Bx87PP57OwGAAC4yNi+huVCLFq0SMOHD9fo0aNDxidPnhz8/+HDh2vEiBEaNGiQqqurdeutt541T2lpqZ566qlO7xcAAJjB1kdC/fr1U3R0tBobG0PGGxsblZiYeM5t/X6/li9frnvvvfe8f2fgwIHq16+f9u3b1+r6kpISeb3e4HL48OH27wQAALjo2AosMTExGjVqlKqqqoJjgUBAVVVVys7OPue2K1euVHNzs37+85+f9+8cOXJEJ0+eVFJSUqvrnU6n4uLiQhYAABC5bH9LqKioSK+++qqWLl2q3bt3a/r06fL7/SosLJQkTZ06VSUlJWdtt2jRIuXn5+uKK64IGW9qatIjjzyiLVu26MCBA6qqqtIdd9yhwYMHKycnp4O7BQAAIonta1gmTZqk48ePa/bs2fJ4PMrMzFRlZWXwQtxDhw4pKio0B+3Zs0ebNm3SunXrzpovOjpaH3/8sZYuXapTp04pOTlZt99+u55++mnuxQIAACRJDsuyrHA3caF8Pp9cLpe8Xi8fD6FLpRavCXcLth2Ym3f+IgDoAnbev/ktIQAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwXocCS1lZmVJTUxUbG6usrCxt3bq1zdqKigo5HI6QJTY2NqTGsizNnj1bSUlJ6tmzp9xut/bu3duR1gAAQASyHVhWrFihoqIizZkzR9u3b1dGRoZycnJ07NixNreJi4tTQ0NDcDl48GDI+ueee04vv/yyysvL9dFHH+nSSy9VTk6OvvrqK/t7BAAAIo7twDJ//nzdd999Kiws1LXXXqvy8nL16tVLixcvbnMbh8OhxMTE4JKQkBBcZ1mWFixYoMcff1x33HGHRowYoWXLluno0aNavXp1h3YKAABEFluB5cyZM6qtrZXb7f5ugqgoud1u1dTUtLldU1OTBgwYoJSUFN1xxx365JNPguvq6+vl8XhC5nS5XMrKympzzubmZvl8vpAFAABELluB5cSJE2ppaQk5QyJJCQkJ8ng8rW6Tnp6uxYsX65133tEf/vAHBQIBjRkzRkeOHJGk4HZ25iwtLZXL5QouKSkpdnYDAABcZDr9W0LZ2dmaOnWqMjMzNXbsWL399tvq37+/fve733V4zpKSEnm93uBy+PDhH7BjAABgGluBpV+/foqOjlZjY2PIeGNjoxITE9s1xyWXXKKRI0dq3759khTczs6cTqdTcXFxIQsAAIhctgJLTEyMRo0apaqqquBYIBBQVVWVsrOz2zVHS0uLdu7cqaSkJElSWlqaEhMTQ+b0+Xz66KOP2j0nAACIbD3sblBUVKSCggJdf/31Gj16tBYsWCC/36/CwkJJ0tSpU3XllVeqtLRUkvTv//7vuvHGGzV48GCdOnVKzz//vA4ePKhp06ZJ+uYbRLNmzdIzzzyjIUOGKC0tTU888YSSk5OVn5//w+0pAAC4aNkOLJMmTdLx48c1e/ZseTweZWZmqrKyMnjR7KFDhxQV9d2Jmy+//FL33XefPB6P+vTpo1GjRmnz5s269tprgzWPPvqo/H6/7r//fp06dUo333yzKisrz7rBHAAA6J4clmVZ4W7iQvl8PrlcLnm9Xq5nQZdKLV4T7hZsOzA3L9wtAIAke+/f/JYQAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHg9wt0A0BVSi9eEuwUAwAXgDAsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeAQWAABgPAILAAAwHoEFAAAYr0OBpaysTKmpqYqNjVVWVpa2bt3aZu2rr76qW265RX369FGfPn3kdrvPqr/nnnvkcDhClgkTJnSkNQAAEIFsB5YVK1aoqKhIc+bM0fbt25WRkaGcnBwdO3as1frq6mpNmTJFH3zwgWpqapSSkqLbb79dn3/+eUjdhAkT1NDQEFzeeOONju0RAACIOLYDy/z583XfffepsLBQ1157rcrLy9WrVy8tXry41fr//u//1r/+678qMzNTQ4cO1WuvvaZAIKCqqqqQOqfTqcTExODSp0+fju0RAACIOLYCy5kzZ1RbWyu32/3dBFFRcrvdqqmpadccf/3rX/X111+rb9++IePV1dWKj49Xenq6pk+frpMnT7Y5R3Nzs3w+X8gCAAAil63AcuLECbW0tCghISFkPCEhQR6Pp11z/PrXv1ZycnJI6JkwYYKWLVumqqoqPfvss9qwYYNyc3PV0tLS6hylpaVyuVzBJSUlxc5uAACAi0yPrvxjc+fO1fLly1VdXa3Y2Njg+OTJk4P/P3z4cI0YMUKDBg1SdXW1br311rPmKSkpUVFRUfCxz+cjtAAAEMFsnWHp16+foqOj1djYGDLe2NioxMTEc247b948zZ07V+vWrdOIESPOWTtw4ED169dP+/bta3W90+lUXFxcyAIAACKXrcASExOjUaNGhVww++0FtNnZ2W1u99xzz+npp59WZWWlrr/++vP+nSNHjujkyZNKSkqy0x4AAIhQtr8lVFRUpFdffVVLly7V7t27NX36dPn9fhUWFkqSpk6dqpKSkmD9s88+qyeeeEKLFy9WamqqPB6PPB6PmpqaJElNTU165JFHtGXLFh04cEBVVVW64447NHjwYOXk5PxAuwkAAC5mtq9hmTRpko4fP67Zs2fL4/EoMzNTlZWVwQtxDx06pKio73LQwoULdebMGf30pz8NmWfOnDl68sknFR0drY8//lhLly7VqVOnlJycrNtvv11PP/20nE7nBe4eAACIBA7LsqxwN3GhfD6fXC6XvF4v17OgVanFa8LdgjEOzM0LdwsAIMne+ze/JQQAAIxHYAEAAMYjsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMF6X/vghcC7cKwUA0BYCC9DNdGYw5KZ0ADoLHwkBAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMYjsAAAAOMRWAAAgPG4NT9s4zd/AABdjTMsAADAeAQWAABgPAILAAAwHoEFAAAYj8ACAACMR2ABAADGI7AAAADjEVgAAIDxCCwAAMB4BBYAAGA8AgsAADAegQUAABiPwAIAAIxHYAEAAMbrUGApKytTamqqYmNjlZWVpa1bt56zfuXKlRo6dKhiY2M1fPhwrV27NmS9ZVmaPXu2kpKS1LNnT7ndbu3du7cjrQEAgAjUw+4GK1asUFFRkcrLy5WVlaUFCxYoJydHe/bsUXx8/Fn1mzdv1pQpU1RaWqp/+qd/0uuvv678/Hxt375d1113nSTpueee08svv6ylS5cqLS1NTzzxhHJycvSXv/xFsbGxF76X3VBq8Zpwt4BuqLOedwfm5nXKvAAuHg7Lsiw7G2RlZemGG27Qb3/7W0lSIBBQSkqKHnzwQRUXF59VP2nSJPn9fr377rvBsRtvvFGZmZkqLy+XZVlKTk7WQw89pIcffliS5PV6lZCQoIqKCk2ePPm8Pfl8PrlcLnm9XsXFxdnZnYhFYEEkIbAAkcnO+7etMyxnzpxRbW2tSkpKgmNRUVFyu92qqalpdZuamhoVFRWFjOXk5Gj16tWSpPr6enk8Hrnd7uB6l8ulrKws1dTUtBpYmpub1dzcHHzs9XolfbPj+Eag+a/hbgH4wfBvG4hM3/7bbs+5E1uB5cSJE2ppaVFCQkLIeEJCgj799NNWt/F4PK3Wezye4Ppvx9qq+b7S0lI99dRTZ42npKS0b0cAXFRcC8LdAYDOdPr0ablcrnPW2L6GxQQlJSUhZ20CgYC++OILXXHFFXI4HD/o3/L5fEpJSdHhw4f5uOk8OFbtx7FqP46VPRyv9uNYtV9nHSvLsnT69GklJyeft9ZWYOnXr5+io6PV2NgYMt7Y2KjExMRWt0lMTDxn/bf/bWxsVFJSUkhNZmZmq3M6nU45nc6Qsd69e9vZFdvi4uJ4QrcTx6r9OFbtx7Gyh+PVfhyr9uuMY3W+MyvfsvW15piYGI0aNUpVVVXBsUAgoKqqKmVnZ7e6TXZ2dki9JK1fvz5Yn5aWpsTExJAan8+njz76qM05AQBA92L7I6GioiIVFBTo+uuv1+jRo7VgwQL5/X4VFhZKkqZOnaorr7xSpaWlkqSZM2dq7NixeuGFF5SXl6fly5dr27ZteuWVVyRJDodDs2bN0jPPPKMhQ4YEv9acnJys/Pz8H25PAQDARct2YJk0aZKOHz+u2bNny+PxKDMzU5WVlcGLZg8dOqSoqO9O3IwZM0avv/66Hn/8cT322GMaMmSIVq9eHbwHiyQ9+uij8vv9uv/++3Xq1CndfPPNqqysNOIeLE6nU3PmzDnrIyicjWPVfhyr9uNY2cPxaj+OVfuZcKxs34cFAACgq/FbQgAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AIqmsrEypqamKjY1VVlaWtm7des76lStXaujQoYqNjdXw4cO1du3aLuo0/Owcq4qKCjkcjpDFhG9+dYWNGzdq4sSJSk5OlsPhCP521rlUV1frxz/+sZxOpwYPHqyKiopO79MEdo9VdXX1Wc8rh8PR5k95RJLS0lLdcMMNuvzyyxUfH6/8/Hzt2bPnvNt1x9esjhyr7vqatXDhQo0YMSJ4U7js7Gz96U9/Ouc24XhOdfvAsmLFChUVFWnOnDnavn27MjIylJOTo2PHjrVav3nzZk2ZMkX33nuvduzYofz8fOXn52vXrl1d3HnXs3uspG/uitjQ0BBcDh482IUdh4/f71dGRobKysraVV9fX6+8vDyNHz9edXV1mjVrlqZNm6b33nuvkzsNP7vH6lt79uwJeW7Fx8d3Uofm2LBhg2bMmKEtW7Zo/fr1+vrrr3X77bfL7/e3uU13fc3qyLGSuudr1lVXXaW5c+eqtrZW27Zt0z/8wz/ojjvu0CeffNJqfdieU1Y3N3r0aGvGjBnBxy0tLVZycrJVWlraav1dd91l5eXlhYxlZWVZ//Iv/9KpfZrA7rFasmSJ5XK5uqg7c0myVq1adc6aRx991Bo2bFjI2KRJk6ycnJxO7Mw87TlWH3zwgSXJ+vLLL7ukJ5MdO3bMkmRt2LChzZru/Jr1f7XnWPGa9Z0+ffpYr732WqvrwvWc6tZnWM6cOaPa2lq53e7gWFRUlNxut2pqalrdpqamJqReknJyctqsjxQdOVaS1NTUpAEDBiglJeWcib27667PqwuRmZmppKQk3Xbbbfrwww/D3U5YeL1eSVLfvn3brOG59Y32HCuJ16yWlhYtX75cfr+/zZ/HCddzqlsHlhMnTqilpSV4l95vJSQktPl5uMfjsVUfKTpyrNLT07V48WK98847+sMf/qBAIKAxY8boyJEjXdHyRaWt55XP59Pf/va3MHVlpqSkJJWXl+utt97SW2+9pZSUFI0bN07bt28Pd2tdKhAIaNasWbrppptC7hz+fd31Nev/au+x6s6vWTt37tRll10mp9OpX/7yl1q1apWuvfbaVmvD9ZyyfWt+oL2ys7NDEvqYMWP0ox/9SL/73e/09NNPh7EzXMzS09OVnp4efDxmzBjt379fL774on7/+9+HsbOuNWPGDO3atUubNm0KdyvGa++x6s6vWenp6aqrq5PX69X//M//qKCgQBs2bGgztIRDtz7D0q9fP0VHR6uxsTFkvLGxUYmJia1uk5iYaKs+UnTkWH3fJZdcopEjR2rfvn2d0eJFra3nVVxcnHr27Bmmri4eo0eP7lbPqwceeEDvvvuuPvjgA1111VXnrO2ur1nfsnOsvq87vWbFxMRo8ODBGjVqlEpLS5WRkaGXXnqp1dpwPae6dWCJiYnRqFGjVFVVFRwLBAKqqqpq87O77OzskHpJWr9+fZv1kaIjx+r7WlpatHPnTiUlJXVWmxet7vq8+qHU1dV1i+eVZVl64IEHtGrVKr3//vtKS0s77zbd9bnVkWP1fd35NSsQCKi5ubnVdWF7TnXqJb0XgeXLl1tOp9OqqKiw/vKXv1j333+/1bt3b8vj8ViWZVl33323VVxcHKz/8MMPrR49eljz5s2zdu/ebc2ZM8e65JJLrJ07d4ZrF7qM3WP11FNPWe+99561f/9+q7a21po8ebIVGxtrffLJJ+HahS5z+vRpa8eOHdaOHTssSdb8+fOtHTt2WAcPHrQsy7KKi4utu+++O1j/2WefWb169bIeeeQRa/fu3VZZWZkVHR1tVVZWhmsXuozdY/Xiiy9aq1evtvbu3Wvt3LnTmjlzphUVFWX9+c9/DtcudJnp06dbLpfLqq6uthoaGoLLX//612ANr1nf6Mix6q6vWcXFxdaGDRus+vp66+OPP7aKi4sth8NhrVu3zrIsc55T3T6wWJZl/eY3v7GuvvpqKyYmxho9erS1ZcuW4LqxY8daBQUFIfVvvvmmdc0111gxMTHWsGHDrDVr1nRxx+Fj51jNmjUrWJuQkGD94z/+o7V9+/YwdN31vv3q7feXb49PQUGBNXbs2LO2yczMtGJiYqyBAwdaS5Ys6fK+w8HusXr22WetQYMGWbGxsVbfvn2tcePGWe+//354mu9irR0nSSHPFV6zvtGRY9VdX7N+8YtfWAMGDLBiYmKs/v37W7feemswrFiWOc8ph2VZVueewwEAALgw3foaFgAAcHEgsAAAAOMRWAAAgPEILAAAwHgEFgAAYDwCCwAAMB6BBQAAGI/AAgAAjEdgAQAAxiOwAAAA4xFYAACA8QgsAADAeP8f72oB3J5OJH8AAAAASUVORK5CYII=\n",
       "text/plain": [
-       "<Figure size 432x288 with 1 Axes>"
+       "<Figure size 640x480 with 1 Axes>"
       ]
      },
-     "metadata": {
-      "needs_background": "light"
-     },
+     "metadata": {},
      "output_type": "display_data"
     }
    ],
@@ -160,35 +156,38 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 10,
+   "execution_count": 9,
    "id": "4a38abae",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "text/plain": [
-       "(array([  0.,   0.,   0.,   0.,   9., 169., 704., 118.,  34.,  25.,  19.,\n",
-       "          9.,   4.,   4.,   2.,   0.,   1.,   0.,   0.,   0.]),\n",
-       " array([0.  , 0.15, 0.3 , 0.45, 0.6 , 0.75, 0.9 , 1.05, 1.2 , 1.35, 1.5 ,\n",
-       "        1.65, 1.8 , 1.95, 2.1 , 2.25, 2.4 , 2.55, 2.7 , 2.85, 3.  ],\n",
-       "       dtype=float32),\n",
+       "(array([0.000000e+00, 0.000000e+00, 0.000000e+00, 2.050000e+02,\n",
+       "        2.136400e+04, 7.107010e+05, 2.630348e+06, 5.204000e+05,\n",
+       "        1.720750e+05, 9.079700e+04, 6.480300e+04, 4.666500e+04,\n",
+       "        3.338600e+04, 2.284300e+04, 1.320600e+04, 8.632000e+03,\n",
+       "        4.469000e+03, 1.948000e+03, 1.068000e+03, 5.150000e+02]),\n",
+       " array([0.        , 0.15000001, 0.30000001, 0.44999999, 0.60000002,\n",
+       "        0.75      , 0.89999998, 1.04999995, 1.20000005, 1.35000002,\n",
+       "        1.5       , 1.64999998, 1.79999995, 1.95000005, 2.0999999 ,\n",
+       "        2.25      , 2.4000001 , 2.54999995, 2.70000005, 2.8499999 ,\n",
+       "        3.        ]),\n",
        " <BarContainer object of 20 artists>)"
       ]
      },
-     "execution_count": 10,
+     "execution_count": 9,
      "metadata": {},
      "output_type": "execute_result"
     },
     {
      "data": {
-      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAARK0lEQVR4nO3dbYxcV33H8e+POAEEFOdh61q2i6mwQGlVglmlRlSIYlHlocKRClFQRUzkylWbtiAqFZcXRVR9Ed5ASVsFWYTWqSiQBmjcEGgtJwj1RQIbCAESaJaIyLaSeAnEPKSAQv99sccwWdae2d3ZXe/h+5FGc+655+78T67z890zd8apKiRJfXnGahcgSRo/w12SOmS4S1KHDHdJ6pDhLkkdWrfaBQBccMEFtXXr1tUuQ5LWlHvuuedbVTUx374zIty3bt3K1NTUapchSWtKkodPtc9lGUnqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOjQ03JO8OMm9A4/vJnlrkvOSHEryYHs+t41PkuuTTCe5L8n25Z+GJGnQ0HCvqq9X1UVVdRHwcuBJ4BPAPuBwVW0DDrdtgEuBbe2xF7hhGeqWJJ3GQj+huhP4RlU9nGQX8OrWfwD4DPB2YBdwU83+KyB3JVmfZGNVPTKmmrWKtu775KKP/eZ1l4+xEkmns9A196uAD7f2hoHAfhTY0NqbgCMDxxxtfU+TZG+SqSRTMzMzCyxDknQ6I4d7knOA1wH/Nndfu0pf0L/XV1X7q2qyqiYnJub93htJ0iIt5Mr9UuALVfVY234syUaA9ny89R8Dtgwct7n1SZJWyELC/Y38bEkG4CCwu7V3A7cO9F/d7prZAZxwvV2SVtZIb6gmeQ7wWuCPBrqvA25Osgd4GLiy9d8OXAZMM3tnzTVjq1aSNJKRwr2qfgCcP6fvcWbvnpk7toBrx1KdJGlR/ISqJHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6NFK4J1mf5JYkX0vyQJJXJDkvyaEkD7bnc9vYJLk+yXSS+5JsX94pSJLmGvXK/X3Ap6vqJcBLgQeAfcDhqtoGHG7bAJcC29pjL3DDWCuWJA01NNyTPB94FXAjQFX9uKqeAHYBB9qwA8AVrb0LuKlm3QWsT7JxzHVLkk5jlCv3FwIzwD8l+WKSDyR5DrChqh5pYx4FNrT2JuDIwPFHW9/TJNmbZCrJ1MzMzOJnIEn6OaOE+zpgO3BDVb0M+AE/W4IBoKoKqIW8cFXtr6rJqpqcmJhYyKGSpCFGCfejwNGqurtt38Js2D92crmlPR9v+48BWwaO39z6JEkrZGi4V9WjwJEkL25dO4H7gYPA7ta3G7i1tQ8CV7e7ZnYAJwaWbyRJK2DdiOP+DPhQknOAh4BrmP2L4eYke4CHgSvb2NuBy4Bp4Mk2VpK0gkYK96q6F5icZ9fOecYWcO3SypIkLYWfUJWkDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nq0EjhnuSbSb6c5N4kU63vvCSHkjzYns9t/UlyfZLpJPcl2b6cE5Ak/byFXLn/TlVdVFWTbXsfcLiqtgGH2zbApcC29tgL3DCuYiVJo1nKsswu4EBrHwCuGOi/qWbdBaxPsnEJryNJWqBRw72A/0pyT5K9rW9DVT3S2o8CG1p7E3Bk4Nijre9pkuxNMpVkamZmZhGlS5JOZd2I4367qo4l+WXgUJKvDe6sqkpSC3nhqtoP7AeYnJxc0LGSpNMb6cq9qo615+PAJ4CLgcdOLre05+Nt+DFgy8Dhm1ufJGmFDA33JM9J8ryTbeB3ga8AB4Hdbdhu4NbWPghc3e6a2QGcGFi+kSStgFGWZTYAn0hycvy/VtWnk3weuDnJHuBh4Mo2/nbgMmAaeBK4ZuxVS5JOa2i4V9VDwEvn6X8c2DlPfwHXjqU6SdKi+AlVSeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1aORwT3JWki8mua1tvzDJ3Ummk3w0yTmt/5lte7rt37pMtUuSTmEhV+5vAR4Y2H438N6qehHwHWBP698DfKf1v7eNkyStoJHCPclm4HLgA207wGuAW9qQA8AVrb2rbdP272zjJUkrZNQr978D/hL4v7Z9PvBEVT3Vto8Cm1p7E3AEoO0/0cY/TZK9SaaSTM3MzCyueknSvIaGe5LfA45X1T3jfOGq2l9Vk1U1OTExMc4fLUm/8NaNMOaVwOuSXAY8C/gl4H3A+iTr2tX5ZuBYG38M2AIcTbIOeD7w+NgrlySd0tAr96r6q6raXFVbgauAO6rqD4A7gde3YbuBW1v7YNum7b+jqmqsVUuSTmsp97m/HXhbkmlm19RvbP03Aue3/rcB+5ZWoiRpoUZZlvmpqvoM8JnWfgi4eJ4xPwTeMIbaJEmL5CdUJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR0y3CWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUoaHhnuRZST6X5EtJvprkXa3/hUnuTjKd5KNJzmn9z2zb023/1mWegyRpjlGu3H8EvKaqXgpcBFySZAfwbuC9VfUi4DvAnjZ+D/Cd1v/eNk6StIKGhnvN+n7bPLs9CngNcEvrPwBc0dq72jZt/84kGVfBkqThRlpzT3JWknuB48Ah4BvAE1X1VBtyFNjU2puAIwBt/wng/DHWLEkaYqRwr6qfVNVFwGbgYuAlS33hJHuTTCWZmpmZWeqPkyQNWNDdMlX1BHAn8ApgfZJ1bddm4FhrHwO2ALT9zwcen+dn7a+qyaqanJiYWFz1kqR5jXK3zESS9a39bOC1wAPMhvzr27DdwK2tfbBt0/bfUVU1xpolSUOsGz6EjcCBJGcx+5fBzVV1W5L7gY8k+Vvgi8CNbfyNwL8kmQa+DVy1DHVLkk5jaLhX1X3Ay+bpf4jZ9fe5/T8E3jCW6iRJi+InVCWpQ4a7JHXIcJekDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6NDTck2xJcmeS+5N8NclbWv95SQ4lebA9n9v6k+T6JNNJ7kuyfbknIUl6ulGu3J8C/qKqLgR2ANcmuRDYBxyuqm3A4bYNcCmwrT32AjeMvWpJ0mkNDfeqeqSqvtDa3wMeADYBu4ADbdgB4IrW3gXcVLPuAtYn2TjuwiVJp7agNfckW4GXAXcDG6rqkbbrUWBDa28CjgwcdrT1zf1Ze5NMJZmamZlZaN2SpNMYOdyTPBf4GPDWqvru4L6qKqAW8sJVtb+qJqtqcmJiYiGHSpKGGCnck5zNbLB/qKo+3rofO7nc0p6Pt/5jwJaBwze3PknSChnlbpkANwIPVNV7BnYdBHa39m7g1oH+q9tdMzuAEwPLN5KkFbBuhDGvBN4EfDnJva3vHcB1wM1J9gAPA1e2fbcDlwHTwJPANeMsWJI03NBwr6r/BnKK3TvnGV/AtUusS5K0BH5CVZI6ZLhLUocMd0nq0ChvqKojW/d9crVLkLQCvHKXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHTLcJalDhrskdchwl6QOGe6S1CHDXZI6ZLhLUocMd0nqkOEuSR3yK3+1YpbydcPfvO7yMVYi9c8rd0nq0NBwT/LBJMeTfGWg77wkh5I82J7Pbf1Jcn2S6ST3Jdm+nMVLkuY3ypX7PwOXzOnbBxyuqm3A4bYNcCmwrT32AjeMp0xJ0kIMDfeq+izw7Tndu4ADrX0AuGKg/6aadRewPsnGMdUqSRrRYtfcN1TVI639KLChtTcBRwbGHW19PyfJ3iRTSaZmZmYWWYYkaT5LfkO1qgqoRRy3v6omq2pyYmJiqWVIkgYsNtwfO7nc0p6Pt/5jwJaBcZtbnyRpBS023A8Cu1t7N3DrQP/V7a6ZHcCJgeUbSdIKGfohpiQfBl4NXJDkKPBO4Drg5iR7gIeBK9vw24HLgGngSeCaZahZkjTE0HCvqjeeYtfOecYWcO1Si5IkLY2fUJWkDhnuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUOGuyR1yHCXpA4Z7pLUIcNdkjpkuEtShwx3SeqQ4S5JHRr6lb/SmWDrvk8u+thvXnf5GCuR1gav3CWpQ165q3te9esXkVfuktQhr9yl01jKVT945a/V45W7JHXIcJekDi3LskySS4D3AWcBH6iq65bjdX5RLXWpQCvHN3O1WsZ+5Z7kLOAfgUuBC4E3Jrlw3K8jSTq15bhyvxiYrqqHAJJ8BNgF3L8MryV1a63+huZvHGeG5Qj3TcCRge2jwG/NHZRkL7C3bX4/ydcX+XoXAN9a5LFnGudy5ullHrBCc8m7l/sVAM/LSS841Y5VuxWyqvYD+5f6c5JMVdXkGEpadc7lzNPLPMC5nKmWay7LcbfMMWDLwPbm1idJWiHLEe6fB7YleWGSc4CrgIPL8DqSpFMY+7JMVT2V5E+B/2T2VsgPVtVXx/06A5a8tHMGcS5nnl7mAc7lTLUsc0lVLcfPlSStIj+hKkkdMtwlqUNrJtyTXJLk60mmk+ybZ/8zk3y07b87ydZVKHMkI8zlzUlmktzbHn+4GnUOk+SDSY4n+cop9ifJ9W2e9yXZvtI1jmqEubw6yYmBc/LXK13jKJJsSXJnkvuTfDXJW+YZsybOy4hzWSvn5VlJPpfkS20u75pnzHgzrKrO+Aezb8x+A/g14BzgS8CFc8b8CfD+1r4K+Ohq172EubwZ+IfVrnWEubwK2A585RT7LwM+BQTYAdy92jUvYS6vBm5b7TpHmMdGYHtrPw/4n3n+fK2J8zLiXNbKeQnw3NY+G7gb2DFnzFgzbK1cuf/0Kw2q6sfAya80GLQLONDatwA7k2QFaxzVKHNZE6rqs8C3TzNkF3BTzboLWJ9k48pUtzAjzGVNqKpHquoLrf094AFmPzU+aE2clxHnsia0/9bfb5tnt8fcu1nGmmFrJdzn+0qDuSf5p2Oq6ingBHD+ilS3MKPMBeD326/MtyTZMs/+tWDUua4Vr2i/Vn8qya+vdjHDtF/rX8bsVeKgNXdeTjMXWCPnJclZSe4FjgOHquqU52UcGbZWwv0XzX8AW6vqN4FD/Oxvc62eLwAvqKqXAn8P/PvqlnN6SZ4LfAx4a1V9d7XrWYohc1kz56WqflJVFzH7qf2Lk/zGcr7eWgn3Ub7S4KdjkqwDng88viLVLczQuVTV41X1o7b5AeDlK1TbuHXzVRRV9d2Tv1ZX1e3A2UkuWOWy5pXkbGbD8ENV9fF5hqyZ8zJsLmvpvJxUVU8AdwKXzNk11gxbK+E+ylcaHAR2t/brgTuqvTNxhhk6lznrn69jdq1xLToIXN3uztgBnKiqR1a7qMVI8isn1z+TXMzs/ztn3MVDq/FG4IGqes8phq2J8zLKXNbQeZlIsr61nw28FvjanGFjzbA18Q9k1ym+0iDJ3wBTVXWQ2T8E/5Jkmtk3xq5avYpPbcS5/HmS1wFPMTuXN69awaeR5MPM3q1wQZKjwDuZfaOIqno/cDuzd2ZMA08C16xOpcONMJfXA3+c5Cngf4GrztCLh1cCbwK+3NZ3Ad4B/CqsufMyylzWynnZCBzI7D9m9Azg5qq6bTkzzK8fkKQOrZVlGUnSAhjuktQhw12SOmS4S1KHDHdJ6pDhLkkdMtwlqUP/D2YsydU15EV0AAAAAElFTkSuQmCC\n",
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGsCAYAAAAPJKchAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAfIklEQVR4nO3de2yV533A8d8BgkkX7ISm2AYcoKGhuRpKuZisBTYaj7Io1rSMRVphbcnWCqpQtlWwS6K0m5woIWHqWGiUBSuNGCnNAA3aBAoFlOCo4qZClrKSC5AEO4mS2MHtHGa/+yOKOzeY+DiYx5fPR3r/8OvnPec5j46Ov3rPe45zWZZlAQCQyIDUEwAA+jcxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJNWrYmT37t1x4403xogRIyKXy8XGjRvzvo0sy+Lee++NK664IgoKCmLkyJHxT//0T+d+sgBApwxKPYF8NDU1RXl5eXzlK1+JP/qjP+rSbdx2222xdevWuPfee+Paa6+NN998M958881zPFMAoLNyvfUf5eVyudiwYUNUVVW17Wtubo6/+7u/i3//93+Pt99+O6655pq4++67Y+bMmRER8dxzz8V1110Xhw8fjvHjx6eZOADQTq96m+bDLF68OGpra2PdunXx85//PG6++eb4gz/4g/jlL38ZERH/+Z//GZ/85Cdj8+bNMXbs2BgzZkwsXLjQmREASKjPxMjx48djzZo1sX79+vjc5z4Xl19+efz1X/91/O7v/m6sWbMmIiJeeOGFOHbsWKxfvz4eeeSRqKmpiX379sUf//EfJ549APRfveqakbM5dOhQtLS0xBVXXNFuf3Nzc3z84x+PiIjW1tZobm6ORx55pG3cv/3bv8WkSZPiyJEj3roBgAT6TIycOnUqBg4cGPv27YuBAwe2+91FF10UERGlpaUxaNCgdsFy5ZVXRsR7Z1bECACcf30mRiZOnBgtLS3x2muvxec+97kzjrn++uvjf//3f+P555+Pyy+/PCIi/vu//zsiIkaPHn3e5goA/Eav+jTNqVOn4ujRoxHxXnzcd999MWvWrBg2bFhcdtll8Wd/9mfx9NNPx4oVK2LixInx+uuvx/bt2+O6666LuXPnRmtra0yePDkuuuiiWLlyZbS2tsaiRYuisLAwtm7dmvjRAUD/1KtiZOfOnTFr1qwP7F+wYEHU1NTE6dOn4x//8R/jkUceiVdeeSUuvfTSmDZtWtx5551x7bXXRkTEq6++Gt/4xjdi69at8Tu/8zsxZ86cWLFiRQwbNux8PxwAIHpZjAAAfU+f+WgvANA7iREAIKle8Wma1tbWePXVV2Po0KGRy+VSTwcA6IQsy+Kdd96JESNGxIABHZ//6BUx8uqrr0ZZWVnqaQAAXXDixIkYNWpUh7/vFTEydOjQiHjvwRQWFiaeDQDQGY2NjVFWVtb2d7wjvSJG3n9rprCwUIwAQC/zYZdYuIAVAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJDUoNQTgJ5szLIt3XK7L901t1tuF6A3cmYEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApPKKkerq6pg8eXIMHTo0hg8fHlVVVXHkyJGzHlNTUxO5XK7dNmTIkI80aQCg78grRnbt2hWLFi2KZ555JrZt2xanT5+OG264IZqams56XGFhYZw8ebJtO3bs2EeaNADQdwzKZ/ATTzzR7ueampoYPnx47Nu3Lz7/+c93eFwul4uSkpKuzRAA6NM+0jUjDQ0NERExbNiws447depUjB49OsrKyuKmm26KZ5999qzjm5ubo7Gxsd0GAPRNXY6R1tbWWLJkSVx//fVxzTXXdDhu/Pjx8fDDD8emTZvi0UcfjdbW1pg+fXq8/PLLHR5TXV0dRUVFbVtZWVlXpwkA9HC5LMuyrhz49a9/PX784x/HU089FaNGjer0cadPn44rr7wybrnllvjOd75zxjHNzc3R3Nzc9nNjY2OUlZVFQ0NDFBYWdmW60CVjlm3pltt96a653XK7AD1JY2NjFBUVfejf77yuGXnf4sWLY/PmzbF79+68QiQi4oILLoiJEyfG0aNHOxxTUFAQBQUFXZkaANDL5PU2TZZlsXjx4tiwYUPs2LEjxo4dm/cdtrS0xKFDh6K0tDTvYwGAvievMyOLFi2KtWvXxqZNm2Lo0KFRV1cXERFFRUVx4YUXRkTE/PnzY+TIkVFdXR0REd/+9rdj2rRpMW7cuHj77bfjnnvuiWPHjsXChQvP8UMBAHqjvGLkgQceiIiImTNnttu/Zs2a+PM///OIiDh+/HgMGPCbEy5vvfVW3HrrrVFXVxeXXHJJTJo0Kfbs2RNXXXXVR5s5ANAndPkC1vOpsxfAwLnmAlaAruvs32//mwYASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUnnFSHV1dUyePDmGDh0aw4cPj6qqqjhy5MiHHrd+/fr49Kc/HUOGDIlrr702fvSjH3V5wgBA35JXjOzatSsWLVoUzzzzTGzbti1Onz4dN9xwQzQ1NXV4zJ49e+KWW26Jr371q3HgwIGoqqqKqqqqOHz48EeePADQ++WyLMu6evDrr78ew4cPj127dsXnP//5M46ZN29eNDU1xebNm9v2TZs2LSZMmBCrV6/u1P00NjZGUVFRNDQ0RGFhYVenC3kbs2xLt9zuS3fN7ZbbBehJOvv3+yNdM9LQ0BAREcOGDetwTG1tbcyePbvdvsrKyqitre3wmObm5mhsbGy3AQB9U5djpLW1NZYsWRLXX399XHPNNR2Oq6uri+Li4nb7iouLo66ursNjqquro6ioqG0rKyvr6jQBgB6uyzGyaNGiOHz4cKxbt+5cziciIpYvXx4NDQ1t24kTJ875fQAAPcOgrhy0ePHi2Lx5c+zevTtGjRp11rElJSVRX1/fbl99fX2UlJR0eExBQUEUFBR0ZWoAQC+T15mRLMti8eLFsWHDhtixY0eMHTv2Q4+pqKiI7du3t9u3bdu2qKioyG+mAECflNeZkUWLFsXatWtj06ZNMXTo0LbrPoqKiuLCCy+MiIj58+fHyJEjo7q6OiIibrvttpgxY0asWLEi5s6dG+vWrYu9e/fGgw8+eI4fCgDQG+V1ZuSBBx6IhoaGmDlzZpSWlrZtjz32WNuY48ePx8mTJ9t+nj59eqxduzYefPDBKC8vjx/+8IexcePGs170CgD0H3mdGenMV5Ls3LnzA/tuvvnmuPnmm/O5KwCgn/C/aQCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAklXeM7N69O2688cYYMWJE5HK52Lhx41nH79y5M3K53Ae2urq6rs4ZAOhD8o6RpqamKC8vj1WrVuV13JEjR+LkyZNt2/Dhw/O9awCgDxqU7wFz5syJOXPm5H1Hw4cPj4svvjjv4wCAvu28XTMyYcKEKC0tjS984Qvx9NNPn3Vsc3NzNDY2ttsAgL6p22OktLQ0Vq9eHY8//ng8/vjjUVZWFjNnzoz9+/d3eEx1dXUUFRW1bWVlZd09TQAgkVyWZVmXD87lYsOGDVFVVZXXcTNmzIjLLrssvv/975/x983NzdHc3Nz2c2NjY5SVlUVDQ0MUFhZ2dbqQtzHLtnTL7b5019xuuV2AnqSxsTGKioo+9O933teMnAtTpkyJp556qsPfFxQUREFBwXmcEQCQSpLvGTl48GCUlpamuGsAoIfJ+8zIqVOn4ujRo20/v/jii3Hw4MEYNmxYXHbZZbF8+fJ45ZVX4pFHHomIiJUrV8bYsWPj6quvjv/5n/+Jhx56KHbs2BFbt249d48CAOi18o6RvXv3xqxZs9p+Xrp0aURELFiwIGpqauLkyZNx/Pjxtt+/++678Vd/9VfxyiuvxMc+9rG47rrr4ic/+Um72wAA+q+PdAHr+dLZC2DgXHMBK0DXdfbvt/9NAwAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASeUdI7t3744bb7wxRowYEblcLjZu3Pihx+zcuTM+85nPREFBQYwbNy5qamq6MFUAoC/KO0aampqivLw8Vq1a1anxL774YsydOzdmzZoVBw8ejCVLlsTChQvjySefzHuyAEDfMyjfA+bMmRNz5szp9PjVq1fH2LFjY8WKFRERceWVV8ZTTz0V999/f1RWVuZ79wBAH9Pt14zU1tbG7Nmz2+2rrKyM2traDo9pbm6OxsbGdhsA0Dd1e4zU1dVFcXFxu33FxcXR2NgYv/71r894THV1dRQVFbVtZWVl3T1NACCRHvlpmuXLl0dDQ0PbduLEidRTAgC6Sd7XjOSrpKQk6uvr2+2rr6+PwsLCuPDCC894TEFBQRQUFHT31ACAHqDbz4xUVFTE9u3b2+3btm1bVFRUdPddAwC9QN4xcurUqTh48GAcPHgwIt776O7Bgwfj+PHjEfHeWyzz589vG/+1r30tXnjhhfjWt74Vv/jFL+Jf//Vf4wc/+EF885vfPDePAADo1fKOkb1798bEiRNj4sSJERGxdOnSmDhxYtx+++0REXHy5Mm2MImIGDt2bGzZsiW2bdsW5eXlsWLFinjooYd8rBcAiIiIXJZlWepJfJjGxsYoKiqKhoaGKCwsTD0d+pExy7Z0y+2+dNfcbrldgJ6ks3+/e+SnaQCA/kOMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUoNSTwA+qjHLtqSeAgAfgTMjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJBUl2Jk1apVMWbMmBgyZEhMnTo1fvazn3U4tqamJnK5XLttyJAhXZ4wANC35B0jjz32WCxdujTuuOOO2L9/f5SXl0dlZWW89tprHR5TWFgYJ0+ebNuOHTv2kSYNAPQdg/I94L777otbb701vvzlL0dExOrVq2PLli3x8MMPx7Jly854TC6Xi5KSko82U+hDxizb0m23/dJdc7vttgG6Q15nRt59993Yt29fzJ49+zc3MGBAzJ49O2prazs87tSpUzF69OgoKyuLm266KZ599tmz3k9zc3M0Nja22wCAvimvGHnjjTeipaUliouL2+0vLi6Ourq6Mx4zfvz4ePjhh2PTpk3x6KOPRmtra0yfPj1efvnlDu+nuro6ioqK2raysrJ8pgkA9CLd/mmaioqKmD9/fkyYMCFmzJgR//Ef/xGf+MQn4nvf+16HxyxfvjwaGhrathMnTnT3NAGARPK6ZuTSSy+NgQMHRn19fbv99fX1nb4m5IILLoiJEyfG0aNHOxxTUFAQBQUF+UwNAOil8jozMnjw4Jg0aVJs3769bV9ra2ts3749KioqOnUbLS0tcejQoSgtLc1vpgBAn5T3p2mWLl0aCxYsiM9+9rMxZcqUWLlyZTQ1NbV9umb+/PkxcuTIqK6ujoiIb3/72zFt2rQYN25cvP3223HPPffEsWPHYuHChef2kQAAvVLeMTJv3rx4/fXX4/bbb4+6urqYMGFCPPHEE20XtR4/fjwGDPjNCZe33norbr311qirq4tLLrkkJk2aFHv27Imrrrrq3D0KAKDXymVZlqWexIdpbGyMoqKiaGhoiMLCwtTToYfpzu/s6I18zwjQU3T277f/TQMAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgKTECACQlRgCApMQIAJCUGAEAkhIjAEBSg1JPADi3xizb0i23+9Jdc7vldgGcGQEAkhIjAEBSYgQASEqMAABJiREAICkxAgAkJUYAgKTECACQlBgBAJISIwBAUr4OHuiU7vqa+QhfNQ/9nTMjAEBSYgQASEqMAABJuWYESK67rkdxLQr0Ds6MAABJOTMC9Fk+AQS9gzMjAEBSXTozsmrVqrjnnnuirq4uysvL47vf/W5MmTKlw/Hr16+Pf/iHf4iXXnopPvWpT8Xdd98dX/ziF7s8aYDUXOcC507eMfLYY4/F0qVLY/Xq1TF16tRYuXJlVFZWxpEjR2L48OEfGL9nz5645ZZborq6Ov7wD/8w1q5dG1VVVbF///645pprzsmDoOfrztPl0Jd4a4n+KJdlWZbPAVOnTo3JkyfHv/zLv0RERGtra5SVlcU3vvGNWLZs2QfGz5s3L5qammLz5s1t+6ZNmxYTJkyI1atXd+o+Gxsbo6ioKBoaGqKwsDCf6dJDiBFIT4xwvnX273deZ0befffd2LdvXyxfvrxt34ABA2L27NlRW1t7xmNqa2tj6dKl7fZVVlbGxo0bO7yf5ubmaG5ubvu5oaEhIt57UPROrc2/Sj0F6Pcu++b61FPI2+E7K1NPgY/g/b/bH3beI68YeeONN6KlpSWKi4vb7S8uLo5f/OIXZzymrq7ujOPr6uo6vJ/q6uq48847P7C/rKwsn+kC0MsVrUw9A86Fd955J4qKijr8fY/8aO/y5cvbnU1pbW2NN998Mz7+8Y9HLpc7Z/fT2NgYZWVlceLECW//fAhrlR/r1XnWqvOsVedZq87rzrXKsizeeeedGDFixFnH5RUjl156aQwcODDq6+vb7a+vr4+SkpIzHlNSUpLX+IiIgoKCKCgoaLfv4osvzmeqeSksLPRk7SRrlR/r1XnWqvOsVedZq87rrrU62xmR9+X1PSODBw+OSZMmxfbt29v2tba2xvbt26OiouKMx1RUVLQbHxGxbdu2DscDAP1L3m/TLF26NBYsWBCf/exnY8qUKbFy5cpoamqKL3/5yxERMX/+/Bg5cmRUV1dHRMRtt90WM2bMiBUrVsTcuXNj3bp1sXfv3njwwQfP7SMBAHqlvGNk3rx58frrr8ftt98edXV1MWHChHjiiSfaLlI9fvx4DBjwmxMu06dPj7Vr18bf//3fx9/+7d/Gpz71qdi4cWOP+I6RgoKCuOOOOz7wlhAfZK3yY706z1p1nrXqPGvVeT1hrfL+nhEAgHPJ/6YBAJISIwBAUmIEAEhKjAAASfX5GFm1alWMGTMmhgwZElOnTo2f/exnZx2/fv36+PSnPx1DhgyJa6+9Nn70ox+dp5mml89a1dTURC6Xa7cNGTLkPM42nd27d8eNN94YI0aMiFwud9b/s/S+nTt3xmc+85koKCiIcePGRU1NTbfPsyfId6127tz5gedVLpc767+P6Cuqq6tj8uTJMXTo0Bg+fHhUVVXFkSNHPvS4/via1ZW16q+vWQ888EBcd911bV9oVlFRET/+8Y/PekyK51SfjpHHHnssli5dGnfccUfs378/ysvLo7KyMl577bUzjt+zZ0/ccsst8dWvfjUOHDgQVVVVUVVVFYcPHz7PMz//8l2riPe+re/kyZNt27Fjx87jjNNpamqK8vLyWLVqVafGv/jiizF37tyYNWtWHDx4MJYsWRILFy6MJ598sptnml6+a/W+I0eOtHtuDR8+vJtm2HPs2rUrFi1aFM8880xs27YtTp8+HTfccEM0NTV1eEx/fc3qylpF9M/XrFGjRsVdd90V+/bti71798bv/d7vxU033RTPPvvsGccne05lfdiUKVOyRYsWtf3c0tKSjRgxIquurj7j+D/5kz/J5s6d227f1KlTs7/8y7/s1nn2BPmu1Zo1a7KioqLzNLueKyKyDRs2nHXMt771rezqq69ut2/evHlZZWVlN86s5+nMWv30pz/NIiJ76623zsucerLXXnsti4hs165dHY7pz69Z/19n1spr1m9ccskl2UMPPXTG36V6TvXZMyPvvvtu7Nu3L2bPnt22b8CAATF79uyora094zG1tbXtxkdEVFZWdji+r+jKWkVEnDp1KkaPHh1lZWVnLe3+rr8+rz6KCRMmRGlpaXzhC1+Ip59+OvV0kmhoaIiIiGHDhnU4xnPrPZ1ZqwivWS0tLbFu3bpoamrq8F+ypHpO9dkYeeONN6KlpaXtm2HfV1xc3OH7z3V1dXmN7yu6slbjx4+Phx9+ODZt2hSPPvpotLa2xvTp0+Pll18+H1PuVTp6XjU2Nsavf/3rRLPqmUpLS2P16tXx+OOPx+OPPx5lZWUxc+bM2L9/f+qpnVetra2xZMmSuP7668/6bdX99TXr/+vsWvXn16xDhw7FRRddFAUFBfG1r30tNmzYEFddddUZx6Z6TuX9dfAQ8d4/QPz/ZT19+vS48sor43vf+1585zvfSTgzerPx48fH+PHj236ePn16PP/883H//ffH97///YQzO78WLVoUhw8fjqeeeir1VHq8zq5Vf37NGj9+fBw8eDAaGhrihz/8YSxYsCB27drVYZCk0GfPjFx66aUxcODAqK+vb7e/vr4+SkpKznhMSUlJXuP7iq6s1W+74IILYuLEiXH06NHumGKv1tHzqrCwMC688MJEs+o9pkyZ0q+eV4sXL47NmzfHT3/60xg1atRZx/bX16z35bNWv60/vWYNHjw4xo0bF5MmTYrq6uooLy+Pf/7nfz7j2FTPqT4bI4MHD45JkybF9u3b2/a1trbG9u3bO3yvrKKiot34iIht27Z1OL6v6Mpa/baWlpY4dOhQlJaWdtc0e63++rw6Vw4ePNgvnldZlsXixYtjw4YNsWPHjhg7duyHHtNfn1tdWavf1p9fs1pbW6O5ufmMv0v2nOrWy2MTW7duXVZQUJDV1NRk//Vf/5X9xV/8RXbxxRdndXV1WZZl2Ze+9KVs2bJlbeOffvrpbNCgQdm9996bPffcc9kdd9yRXXDBBdmhQ4dSPYTzJt+1uvPOO7Mnn3wye/7557N9+/Zlf/qnf5oNGTIke/bZZ1M9hPPmnXfeyQ4cOJAdOHAgi4jsvvvuyw4cOJAdO3Ysy7IsW7ZsWfalL32pbfwLL7yQfexjH8v+5m/+JnvuueeyVatWZQMHDsyeeOKJVA/hvMl3re6///5s48aN2S9/+cvs0KFD2W233ZYNGDAg+8lPfpLqIZw3X//617OioqJs586d2cmTJ9u2X/3qV21jvGa9pytr1V9fs5YtW5bt2rUre/HFF7Of//zn2bJly7JcLpdt3bo1y7Ke85zq0zGSZVn23e9+N7vsssuywYMHZ1OmTMmeeeaZtt/NmDEjW7BgQbvxP/jBD7IrrrgiGzx4cHb11VdnW7ZsOc8zTieftVqyZEnb2OLi4uyLX/xitn///gSzPv/e//jpb2/vr8+CBQuyGTNmfOCYCRMmZIMHD84++clPZmvWrDnv804h37W6++67s8svvzwbMmRINmzYsGzmzJnZjh070kz+PDvTOkVEu+eK16z3dGWt+utr1le+8pVs9OjR2eDBg7NPfOIT2e///u+3hUiW9ZznVC7Lsqx7z70AAHSsz14zAgD0DmIEAEhKjAAASYkRACApMQIAJCVGAICkxAgAkJQYAQCSEiMAQFJiBABISowAAEmJEQAgqf8DgyUK1RABovAAAAAASUVORK5CYII=\n",
       "text/plain": [
-       "<Figure size 432x288 with 1 Axes>"
+       "<Figure size 640x480 with 1 Axes>"
       ]
      },
-     "metadata": {
-      "needs_background": "light"
-     },
+     "metadata": {},
      "output_type": "display_data"
     }
    ],
@@ -206,7 +205,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 16,
+   "execution_count": 10,
    "id": "8e051b62-1e63-467d-9f43-6a48813a0e2b",
    "metadata": {},
    "outputs": [],
@@ -218,7 +217,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 20,
+   "execution_count": 11,
    "id": "4ea29864-913c-4dbe-8751-4be371b6d396",
    "metadata": {},
    "outputs": [],
@@ -228,7 +227,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 21,
+   "execution_count": 12,
    "id": "4b09a818-f03a-45ef-996e-bcd431e9a8f4",
    "metadata": {},
    "outputs": [],
@@ -239,7 +238,28 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 22,
+   "execution_count": 16,
+   "id": "5e951704",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "<coffea.jetmet_tools.CorrectedJetsFactory.CorrectedJetsFactory at 0x7f16a65fc190>"
+      ]
+     },
+     "execution_count": 16,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "jmestuff['fatjet_factory']['2018']"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 17,
    "id": "6b9387c1-28d0-4e28-bf35-bebce40cfea7",
    "metadata": {},
    "outputs": [],
@@ -254,7 +274,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 23,
+   "execution_count": 18,
    "id": "10b431d3-6bab-4a0e-a474-2f5257e8623e",
    "metadata": {},
    "outputs": [],
@@ -265,7 +285,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 27,
+   "execution_count": 19,
    "id": "0f0255f9-238c-461c-a003-1113e353996c",
    "metadata": {},
    "outputs": [
@@ -273,41 +293,87 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "CPU times: user 3.87 s, sys: 616 ms, total: 4.49 s\n",
-      "Wall time: 4.48 s\n"
+      "CPU times: user 2.3 s, sys: 165 ms, total: 2.46 s\n",
+      "Wall time: 7.43 s\n"
      ]
     }
    ],
    "source": [
     "%%time\n",
-    "j = jmestuff[\"jet_factory\"][\"2018\"].build(add_jec_variables(events.Jet,events.fixedGridRhoFastjetAll), jec_cache)"
+    "#j = jmestuff[\"jet_factory\"][\"2018\"].build(add_jec_variables(events.Jet,events.fixedGridRhoFastjetAll), jec_cache)\n",
+    "j = jmestuff[\"fatjet_factory\"][\"2018\"].build(add_jec_variables(events.FatJet,events.fixedGridRhoFastjetAll), jec_cache)"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 28,
+   "execution_count": 20,
    "id": "0f86f278-ab3e-4e3f-9c18-f2820cb4f840",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "text/plain": [
-       "['JES_FlavorQCD',\n",
+       "['JES_AbsoluteStat',\n",
+       " 'JES_AbsoluteScale',\n",
+       " 'JES_AbsoluteSample',\n",
+       " 'JES_AbsoluteFlavMap',\n",
+       " 'JES_AbsoluteMPFBias',\n",
+       " 'JES_Fragmentation',\n",
+       " 'JES_SinglePionECAL',\n",
+       " 'JES_SinglePionHCAL',\n",
+       " 'JES_FlavorQCD',\n",
+       " 'JES_TimePtEta',\n",
+       " 'JES_RelativeJEREC1',\n",
+       " 'JES_RelativeJEREC2',\n",
+       " 'JES_RelativeJERHF',\n",
+       " 'JES_RelativePtBB',\n",
+       " 'JES_RelativePtEC1',\n",
+       " 'JES_RelativePtEC2',\n",
+       " 'JES_RelativePtHF',\n",
        " 'JES_RelativeBal',\n",
-       " 'JES_HF',\n",
-       " 'JES_BBEC1',\n",
-       " 'JES_EC2',\n",
-       " 'JES_Absolute',\n",
-       " 'JES_Absolute_2018',\n",
-       " 'JES_HF_2018',\n",
-       " 'JES_EC2_2018',\n",
-       " 'JES_RelativeSample_2018',\n",
-       " 'JES_BBEC1_2018',\n",
+       " 'JES_RelativeSample',\n",
+       " 'JES_RelativeFSR',\n",
+       " 'JES_RelativeStatFSR',\n",
+       " 'JES_RelativeStatEC',\n",
+       " 'JES_RelativeStatHF',\n",
+       " 'JES_PileUpDataMC',\n",
+       " 'JES_PileUpPtRef',\n",
+       " 'JES_PileUpPtBB',\n",
+       " 'JES_PileUpPtEC1',\n",
+       " 'JES_PileUpPtEC2',\n",
+       " 'JES_PileUpPtHF',\n",
+       " 'JES_PileUpMuZero',\n",
+       " 'JES_PileUpEnvelope',\n",
+       " 'JES_SubTotalPileUp',\n",
+       " 'JES_SubTotalRelative',\n",
+       " 'JES_SubTotalPt',\n",
+       " 'JES_SubTotalScale',\n",
+       " 'JES_SubTotalAbsolute',\n",
+       " 'JES_SubTotalMC',\n",
        " 'JES_Total',\n",
+       " 'JES_TotalNoFlavor',\n",
+       " 'JES_TotalNoTime',\n",
+       " 'JES_TotalNoFlavorNoTime',\n",
+       " 'JES_FlavorZJet',\n",
+       " 'JES_FlavorPhotonJet',\n",
+       " 'JES_FlavorPureGluon',\n",
+       " 'JES_FlavorPureQuark',\n",
+       " 'JES_FlavorPureCharm',\n",
+       " 'JES_FlavorPureBottom',\n",
+       " 'JES_TimeRunA',\n",
+       " 'JES_TimeRunB',\n",
+       " 'JES_TimeRunC',\n",
+       " 'JES_TimeRunD',\n",
+       " 'JES_CorrelationGroupMPFInSitu',\n",
+       " 'JES_CorrelationGroupIntercalibration',\n",
+       " 'jet_energy_uncertainty_CorrelationGroupbJES',\n",
+       " 'JES_CorrelationGroupbJES',\n",
+       " 'JES_CorrelationGroupFlavor',\n",
+       " 'JES_CorrelationGroupUncorrelated',\n",
        " 'JES_jes']"
       ]
      },
-     "execution_count": 28,
+     "execution_count": 20,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -662,7 +728,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.9.15"
+   "version": "3.8.15"
   }
  },
  "nbformat": 4,
-- 
GitLab


From a6356bf358e2c472011b5723ae068a82481f8501 Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Fri, 10 Feb 2023 16:32:17 +0100
Subject: [PATCH 02/11] Temporary solution for plotting ggH

---
 pocket_coffea/utils/plot_utils.py | 78 ++++++++++++++++++++++++++++++-
 scripts/plot/make_plots.py        |  4 +-
 2 files changed, 79 insertions(+), 3 deletions(-)

diff --git a/pocket_coffea/utils/plot_utils.py b/pocket_coffea/utils/plot_utils.py
index ae7dc0fa..df6004f5 100644
--- a/pocket_coffea/utils/plot_utils.py
+++ b/pocket_coffea/utils/plot_utils.py
@@ -277,6 +277,9 @@ def plot_data_mc_hist1D(
     samples = h.keys()
     samples_data = list(filter(lambda d: 'DATA' in d, samples))
     samples_mc = list(filter(lambda d: 'DATA' not in d, samples))
+    
+    print("********************")
+    print(samples_mc)
 
     h_mc = h[samples_mc[0]]
 
@@ -292,6 +295,8 @@ def plot_data_mc_hist1D(
     for year in years:
 
         for cat in categories:
+            if 'noreweight' in cat:
+                continue
             if only_cat:
                 if isinstance(only_cat, list):
                     if not cat in only_cat:
@@ -370,7 +375,7 @@ def plot_data_mc_hist1D(
                             [
                                 h[d][slicing_mc]
                                 for d in samples_mc
-                                if d.endswith(f'_{f}')
+                                if (d.endswith(f'_{f}') & ('GluGluH' not in d))
                             ]
                         )
                     )
@@ -382,7 +387,7 @@ def plot_data_mc_hist1D(
                             [
                                 h[d][slicing_mc_nominal]
                                 for d in samples_mc
-                                if d.endswith(f'_{f}')
+                                if (d.endswith(f'_{f}') & ('GluGluH' not in d))
                             ]
                         )
                     )
@@ -392,6 +397,57 @@ def plot_data_mc_hist1D(
                 nevents = {
                     f: round(sum(dict_mc_nominal[f].values()), 1) for f in flavors
                 }
+                print("************")
+                print(histname)
+                print(samples_mc)
+                dict_gghbb = {
+                    f: stack_sum(
+                        hist.Stack.from_iter(
+                            [
+                                h[d][slicing_mc]
+                                for d in samples_mc
+                                if (d.endswith(f'_{f}') & ('GluGluHToBB' in d))
+                            ]
+                        )
+                    )
+                    for f in flavors
+                }
+                dict_gghbb_nominal = {
+                    f: stack_sum(
+                        hist.Stack.from_iter(
+                            [
+                                h[d][slicing_mc_nominal]
+                                for d in samples_mc
+                                if (d.endswith(f'_{f}') & ('GluGluHToBB' in d))
+                            ]
+                        )
+                    )
+                    for f in flavors
+                }
+                dict_gghcc = {
+                    f: stack_sum(
+                        hist.Stack.from_iter(
+                            [
+                                h[d][slicing_mc]
+                                for d in samples_mc
+                                if (d.endswith(f'_{f}') & ('GluGluHToCC' in d))
+                            ]
+                        )
+                    )
+                    for f in flavors
+                }
+                dict_gghcc_nominal = {
+                    f: stack_sum(
+                        hist.Stack.from_iter(
+                            [
+                                h[d][slicing_mc_nominal]
+                                for d in samples_mc
+                                if (d.endswith(f'_{f}') & ('GluGluHToCC' in d))
+                            ]
+                        )
+                    )
+                    for f in flavors
+                }
             else:
                 dict_mc = {d: h[d][slicing_mc] for d in samples_mc}
                 dict_mc_nominal = {d: h[d][slicing_mc_nominal] for d in samples_mc}
@@ -406,6 +462,15 @@ def plot_data_mc_hist1D(
                     dict_mc_nominal[sample] = histo_reweighted
             stack_mc = hist.Stack.from_dict(dict_mc)
             stack_mc_nominal = hist.Stack.from_dict(dict_mc_nominal)
+            if flavorsplit == '5f':
+                stack_gghbb = hist.Stack.from_dict(dict_gghbb)
+                stack_gghbb_nominal = hist.Stack.from_dict(dict_gghbb_nominal)
+                stack_gghcc = hist.Stack.from_dict(dict_gghcc)
+                stack_gghcc_nominal = hist.Stack.from_dict(dict_gghcc_nominal)
+                hist_gghbb = stack_sum(stack_gghbb)
+                hist_gghbb_nominal = stack_sum(stack_gghbb_nominal)
+                hist_gghcc = stack_sum(stack_gghcc)
+                hist_gghcc_nominal = stack_sum(stack_gghcc_nominal)
 
             if not is_mc_only:
                 # Sum over eras if era axis exists in data histogram
@@ -463,6 +528,12 @@ def plot_data_mc_hist1D(
                     x, h_data.values(), yerr=np.sqrt(h_data.values()), **opts_data
                 )
                 # stack_data.plot(stack=True, color='black', ax=ax)
+            if flavorsplit == '5f':
+                sf = 100
+                hist_gghbb_nominal = sf * hist_gghbb_nominal
+                hist_gghcc_nominal = sf * hist_gghcc_nominal
+                hist_gghbb_nominal.plot(stack=False, histtype='step', ax=ax, color='red', label=f'ggHbb x {sf}')
+                hist_gghcc_nominal.plot(stack=False, histtype='step', ax=ax, color='green', label=f'ggHcc x {sf}')
             if rebinning:
                 syst_err_up, syst_err_down = get_systematic_uncertainty(
                     stack_mc,
@@ -515,6 +586,7 @@ def plot_data_mc_hist1D(
             rax.tick_params(axis='y', labelsize=fontsize)
 
             if log:
+                ax.set_yscale("log")
                 ax.set_xlim(0, 1)
                 exp = math.floor(
                     math.log(max(stack_sum(stack_mc_nominal).values()), 10)
@@ -529,6 +601,7 @@ def plot_data_mc_hist1D(
                     if histname in config.plot_options.keys():
                         if 'xlim' in config.plot_options[histname].keys():
                             ax.set_xlim(*config.plot_options[histname]['xlim'])
+                        """
                         if 'scale' in config.plot_options[histname].keys():
                             ax.set_yscale(config.plot_options[histname]['scale'])
                             if config.plot_options[histname]['scale'] == 'log':
@@ -539,6 +612,7 @@ def plot_data_mc_hist1D(
                                 )
                                 ax.set_ylim((0.01, 10 ** (exp + 2)))
                                 # ax.legend(handles, labels, loc="upper right", fontsize=fontsize, ncols=2)
+                        """
                         if 'ylim' in config.plot_options[histname].keys():
                             if isinstance(config.plot_options[histname]['ylim'], tuple):
                                 ax.set_ylim(*config.plot_options[histname]['ylim'])
diff --git a/scripts/plot/make_plots.py b/scripts/plot/make_plots.py
index e998b553..76d1edef 100644
--- a/scripts/plot/make_plots.py
+++ b/scripts/plot/make_plots.py
@@ -112,7 +112,9 @@ if not os.path.exists(config.plots):
 def make_plots(entrystart, entrystop):
     _accumulator = slice_accumulator(accumulator, entrystart, entrystop)
     for (histname, h) in _accumulator['variables'].items():
-        plot_data_mc_hist1D(h, histname, config, flavorsplit=None, mcstat=True, stat_only=False)
+        #if ('FatJetGood_DDX' in histname) | ('FatJetGood_particleNetMD_Xqq' in histname) | ('FatJetGood_particleNet_H4qvsQCD' in histname):
+        #    continue
+        plot_data_mc_hist1D(h, histname, config, flavorsplit='5f', mcstat=True, stat_only=False, log=True)
 
 # Filter dictionary of histograms with `args.only`
 accumulator['variables'] = { k : v for k,v in accumulator['variables'].items() if args.only in k }
-- 
GitLab


From b63fe558a07113ab6c698ea379f11b4e0706b546 Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Thu, 16 Mar 2023 11:04:59 +0100
Subject: [PATCH 03/11] Nicer log plots

---
 pocket_coffea/utils/plot_utils.py | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/pocket_coffea/utils/plot_utils.py b/pocket_coffea/utils/plot_utils.py
index 00c0d4ab..09d688ae 100644
--- a/pocket_coffea/utils/plot_utils.py
+++ b/pocket_coffea/utils/plot_utils.py
@@ -728,7 +728,7 @@ def plot_data_mc_hist1D(
                 exp = math.floor(
                     math.log(max(stack_sum(stack_mc_nominal).values()), 10)
                 )
-                ax.set_ylim((0.01, 10 ** (exp + 2)))
+                ax.set_ylim((0.01, 10 ** (exp + 3)))
                 # if flavorsplit == '5f':
                 #    ax.legend(handles, labels, loc="upper right", fontsize=fontsize, ncols=2)
                 # else:
@@ -747,7 +747,7 @@ def plot_data_mc_hist1D(
                             handles_new.append(handles[i])
                         labels = labels_new
                         handles = handles_new
-                        ax.legend(handles, labels, fontsize=fontsize, ncols=2)
+                        ax.legend(handles, labels, fontsize=fontsize, ncols=2, loc="upper right")
                     if ("variables" in config.plot_options) & (histname in config.plot_options["variables"].keys()):
                         cfg_plot = config.plot_options["variables"][histname]
                         if 'xlim' in cfg_plot.keys():
@@ -760,7 +760,7 @@ def plot_data_mc_hist1D(
                                         max(stack_sum(stack_mc_nominal).values()), 10
                                     )
                                 )
-                                ax.set_ylim((0.01, 10 ** (exp + 2)))
+                                ax.set_ylim((0.01, 10 ** (exp + 3)))
                                 # ax.legend(handles, labels, loc="upper right", fontsize=fontsize, ncols=2)
                         if 'ylim' in cfg_plot.keys():
                             if isinstance(cfg_plot['ylim'], tuple):
-- 
GitLab


From 05a3c61438e2cc6a828fe0026456ae3d20005ca6 Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Mon, 20 Mar 2023 23:31:20 +0100
Subject: [PATCH 04/11] Update genWeightsProcessor for skim

---
 pocket_coffea/workflows/genweights.py | 9 +++++++++
 1 file changed, 9 insertions(+)

diff --git a/pocket_coffea/workflows/genweights.py b/pocket_coffea/workflows/genweights.py
index 1a404377..00fe5d08 100644
--- a/pocket_coffea/workflows/genweights.py
+++ b/pocket_coffea/workflows/genweights.py
@@ -48,4 +48,13 @@ class genWeightsProcessor(BaseProcessorABC):
             if self._hasSubsamples:
                 raise Exception("This processor cannot compute the sum of genweights of subsamples.")
 
+        self.process_extra_before_skim()
+        self.skim_events()
+        if not self.has_events:
+            return self.output
+
+        if self.cfg.save_skimmed_files:
+            self.export_skimmed_chunk()
+            return self.output
+
         return self.output
-- 
GitLab


From 35de55d76ed059c8e85c13e07699bc554acfe966 Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Tue, 21 Mar 2023 01:13:54 +0100
Subject: [PATCH 05/11] Updates for skim

---
 pocket_coffea/workflows/genweights.py | 48 ++++++++++++++++++++++++++-
 1 file changed, 47 insertions(+), 1 deletion(-)

diff --git a/pocket_coffea/workflows/genweights.py b/pocket_coffea/workflows/genweights.py
index 00fe5d08..dc0ffab7 100644
--- a/pocket_coffea/workflows/genweights.py
+++ b/pocket_coffea/workflows/genweights.py
@@ -1,9 +1,16 @@
+import numpy as np
 import awkward as ak
 
 import copy
 
+from coffea.lumi_tools import LumiMask
+from coffea.analysis_tools import PackedSelection
+
 from .base import BaseProcessorABC
 from ..utils.configurator import Configurator
+from ..parameters.event_flags import event_flags, event_flags_data
+from ..parameters.lumi import goldenJSON
+from ..parameters.btag import btag
 
 class genWeightsProcessor(BaseProcessorABC):
     def __init__(self, cfg: Configurator) -> None:
@@ -11,7 +18,8 @@ class genWeightsProcessor(BaseProcessorABC):
         self.output_format = {
             "sum_genweights": {},
             "cutflow": {
-                "initial": {s: 0 for s in self.cfg.datasets}
+                "initial": {s: 0 for s in self.cfg.datasets},
+                "skim": {s: 0 for s in self.cfg.datasets}
             }
         }
 
@@ -19,9 +27,14 @@ class genWeightsProcessor(BaseProcessorABC):
         self._dataset = self.events.metadata["dataset"]
         self._sample = self.events.metadata["sample"]
         self._year = self.events.metadata["year"]
+        self._btag = btag[self._year]
         self._isMC = self.events.metadata["isMC"] == "True"
         if self._isMC:
+            self._era = "MC"
             self._xsec = self.events.metadata["xsec"]
+        else:
+            self._era = self.events.metadata["era"]
+            self._goldenJSON = goldenJSON[self._year]
 
         # Check if the user specified any subsamples without performing any operation
         if self._sample in self._subsamplesCfg:
@@ -29,6 +42,39 @@ class genWeightsProcessor(BaseProcessorABC):
         else:
             self._hasSubsamples = False
 
+    def skim_events(self):
+        self._skim_masks = PackedSelection()
+        mask_flags = np.ones(self.nEvents_initial, dtype=np.bool)
+        flags = event_flags[self._year]
+        if not self._isMC:
+            flags += event_flags_data[self._year]
+        for flag in flags:
+            mask_flags = getattr(self.events.Flag, flag)
+        self._skim_masks.add("event_flags", mask_flags)
+
+        self._skim_masks.add("PVgood", self.events.PV.npvsGood > 0)
+
+        if not self._isMC:
+            mask_lumi = LumiMask(self._goldenJSON)(
+                self.events.run, self.events.luminosityBlock
+            )
+            self._skim_masks.add("lumi_golden", mask_lumi)
+
+        for skim_func in self._skim:
+            mask = skim_func.get_mask(
+                self.events,
+                year=self._year,
+                sample=self._sample,
+                btag=self._btag,
+                isMC=self._isMC,
+            )
+            self._skim_masks.add(skim_func.id, mask)
+
+        self.events = self.events[self._skim_masks.all(*self._skim_masks.names)]
+        self.nEvents_after_skim = self.nevents
+        self.output['cutflow']['skim'][self._dataset] += self.nEvents_after_skim
+        self.has_events = self.nEvents_after_skim > 0
+
     def apply_object_preselection(self, variation):
         pass
 
-- 
GitLab


From 22cc6e329ec02c2e6c5b81b82f35be0211dc9901 Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Tue, 21 Mar 2023 13:28:35 +0100
Subject: [PATCH 06/11] Documenting hadd script

---
 scripts/README.md             | 85 +++++++++++++++++++++++++++++++++++
 scripts/hadd_skimmed_files.py |  9 ++--
 2 files changed, 91 insertions(+), 3 deletions(-)
 create mode 100644 scripts/README.md

diff --git a/scripts/README.md b/scripts/README.md
new file mode 100644
index 00000000..940cc5a1
--- /dev/null
+++ b/scripts/README.md
@@ -0,0 +1,85 @@
+# PocketCoffea scripts
+
+The current guide is a documentation of how to run the scripts of the PocketCoffea package.
+
+## Merge skimmed files with hadd
+
+In order to merge the skimmed n-tuples, the `hadd` command from ROOT needs to be used.
+To source ROOT inside the `pocket_coffea` environment, run the following command:
+
+```
+source /cvmfs/sft.cern.ch/lcg/views/LCG_102/x86_64-centos7-gcc11-opt/setup.sh
+```
+
+To run the `hadd` script run the following command:
+
+```
+cd /path/to/PocketCoffea
+python scripts/hadd_skimmed_files.py -fl /path/to/output/output_all.coffea -o /path/to/outputdir
+```
+
+## Build dataset
+
+The first step is to build the json files containing the list of files of the UL datasets to process. In our case, we need to include the `TTToSemiLeptonic` and `TTTo2L2Nu` MC datasets and `SingleMuon` as data dataset.
+
+A step-by-step detailed guide on the creation of the dataset file can be found at this [link](https://pocketcoffea.readthedocs.io/en/latest/examples.html).
+
+### Append sum of genweights to datasets definitions
+
+In order to update the `datasets_definitions.json` and include the `sum_genweights` metadata in the datasets definitions, a dedicated script needs to be run with the following command:
+
+```
+cd /path/to/PocketCoffea
+scripts/dataset/append_genweights.py --cfg dataset/dataset_definitions.json -i output/genweights/genweights_2016_PreVFP/output_all.coffea --overwrite
+scripts/dataset/append_genweights.py --cfg dataset/dataset_definitions.json -i output/genweights/genweights_2016_PostVFP/output_all.coffea --overwrite
+```
+
+N.B.: the argument `-i output/genweights/genweights_201*/output_all.coffea` should be the path to the Coffea output file where the genweights are stored. The option `--overwrite` automatically updates the config file passed as `--cfg` argument with the genweights, in this case `dataset/dataset_definitions.json`.
+
+### Build dataset with sum of genweights
+
+Now the steps of the [Build dataset](#build-dataset) section need to be repeated to generate the json datasets with the additional metadata `sum_genweights`.
+
+## Run the analysis
+
+In order to run the analysis workflow and produce the output histograms, run the following command:
+
+```
+cd /path/to/PocketCoffea
+runner.py --cfg /path/to/config/config.py --full
+```
+
+N.B.: the argument `--full` will process all the datasets together at once and save the output in a single output file, `output_all.coffea`. Otherwise, the datasets are processed separately and an output file is saved for each dataset.
+
+## Accumulate output files
+
+If the output of the [previous step](#run-the-analysis) has been produced without the argument `--full`, the output files need to be merged in a single output file `output_all.coffea`. If the output has been produced with the argument `--full` and the output file `output_all.coffea` is already existing, skip this step and continue with the [next one](#produce-datamc-plots).
+
+Once the Coffea output files are produced, one needs to merge the files into a single file by using the script `accumulate_files.py` by running this command:
+
+```
+cd /path/to/output
+python accumulate_files.py -i output1.coffea output2.coffea output3.coffea -o output_all.coffea
+```
+
+## Produce data/MC plots
+
+To produce plots from a Coffea output file execute the plotting script:
+
+```
+cd /path/to/PocketCoffea
+make_plots.py --cfg /path/to/config/config.py -i /path/to/output/output_all.coffea
+```
+
+## Run trigger SF script
+
+To run the script that computes the single electron trigger SF and produces the validation plots, run:
+
+```
+cd /path/to/PocketCoffea
+python scripts/plot/trigger_efficiency.py --cfg config/semileptonic_triggerSF/semileptonic_triggerSF_2017.py --save_plots
+```
+
+If the additional argument `--save_plots` is not specified, only the scale factor maps are saved without saving the plots.
+
+The output plots are saved in `/path/to/output/plots/trigger_efficiency` and `/path/to/output/plots/trigger_scalefactor`, while the 1D and 2D scale factor maps are saved in the folder specified in `workflow_options['output_triggerSF']` in the config file.
diff --git a/scripts/hadd_skimmed_files.py b/scripts/hadd_skimmed_files.py
index 2280fcca..6d870ba2 100644
--- a/scripts/hadd_skimmed_files.py
+++ b/scripts/hadd_skimmed_files.py
@@ -18,7 +18,7 @@ parser.add_argument(
 )
 parser.add_argument("-o", "--outputdir", required=False, type=str, help="Output folder")
 parser.add_argument(
-    "--only-datasets", nargs="+", type=str, help="Restring list of datasets"
+    "--only-datasets", nargs="+", type=str, help="Restricting list of datasets"
 )
 parser.add_argument("-f", "--files", type=int, help="Limit number of files")
 parser.add_argument("-e", "--events", type=int, help="Limit number of files")
@@ -43,6 +43,9 @@ groups_metadata = {}
 for dataset in df["skimmed_files"].keys():
     if args.only_datasets and dataset not in args.only_datasets:
         continue
+    outputdir_dataset = os.path.join(args.outputdir, dataset)
+    if not os.path.exists(outputdir_dataset):
+        os.makedirs(outputdir_dataset)
     groups_metadata[dataset] = defaultdict(dict)
     nevents_tot = 0
     nfiles = 0
@@ -54,7 +57,7 @@ for dataset in df["skimmed_files"].keys():
         if (args.files and (nfiles + 1) > args.files) or (
             args.events and (nevents_tot + nevents) > args.events
         ):
-            outfile = f"{args.outputdir}/{dataset}/{dataset}_{ngroup}.root"
+            outfile = os.path.join(outputdir_dataset, f"{dataset}_{ngroup}.root")
             workload.append((outfile, group[:]))
             groups_metadata[dataset]["files"][outfile] = group[:]
             group.clear()
@@ -68,7 +71,7 @@ for dataset in df["skimmed_files"].keys():
 
     # add last group
     if len(group):
-        outfile = f"{args.outputdir}/{dataset}/{dataset}_{ngroup}.root"
+        outfile = os.path.join(outputdir_dataset, f"{dataset}_{ngroup}.root")
         workload.append((outfile, group[:]))
         groups_metadata[dataset]["files"][outfile] = group[:]
 
-- 
GitLab


From 2cdc97f3411628950176da20d345a999ab807b1f Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Wed, 22 Mar 2023 03:06:57 +0100
Subject: [PATCH 07/11] Update object preselection for mutag

---
 pocket_coffea/parameters/object_preselection.py | 3 ---
 1 file changed, 3 deletions(-)

diff --git a/pocket_coffea/parameters/object_preselection.py b/pocket_coffea/parameters/object_preselection.py
index e93c0453..291f816d 100644
--- a/pocket_coffea/parameters/object_preselection.py
+++ b/pocket_coffea/parameters/object_preselection.py
@@ -89,9 +89,6 @@ object_preselection = {
             "eta": 2.4,
             "msd": 40,
             "jetId": 2,
-            "nsubjet": 2,
-            "nmusj": 1,
-            "dimuon_pt_ratio": 0.6,
         },
     },
 }
-- 
GitLab


From ff7cb05bfd6a8472aa4849c1965e29b32eaa695e Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Fri, 24 Mar 2023 02:02:36 +0100
Subject: [PATCH 08/11] Define fatjet selection

---
 pocket_coffea/lib/jets.py | 14 ++++++++------
 1 file changed, 8 insertions(+), 6 deletions(-)

diff --git a/pocket_coffea/lib/jets.py b/pocket_coffea/lib/jets.py
index 79051fde..07e37cee 100644
--- a/pocket_coffea/lib/jets.py
+++ b/pocket_coffea/lib/jets.py
@@ -226,7 +226,7 @@ def jet_selection(events, Jet, finalstate):
     # Only jets that are more distant than dr to ALL leptons are tagged as good jets
     leptons = events["LeptonGood"]
     # Mask for  jets not passing the preselection
-    presel_mask = (
+    mask_presel = (
         (jets.pt > cuts["pt"])
         & (np.abs(jets.eta) < cuts["eta"])
         & (jets.jetId >= cuts["jetId"])
@@ -234,18 +234,20 @@ def jet_selection(events, Jet, finalstate):
     # Lepton cleaning
     if "dr" in cuts.keys():
         dR_jets_lep = jets.metric_table(leptons)
-        lepton_cleaning_mask = ak.prod(dR_jets_lep > cuts["dr"], axis=2) == 1
+        mask_lepton_cleaning = ak.prod(dR_jets_lep > cuts["dr"], axis=2) == 1
 
     if Jet == "Jet":
-        jetpuid_mask = (jets.puId >= cuts["puId"]["value"]) | (
+        mask_jetpuid = (jets.puId >= cuts["puId"]["value"]) | (
             jets.pt >= cuts["puId"]["maxpt"]
         )
-        good_jets_mask = presel_mask & lepton_cleaning_mask & jetpuid_mask
+        mask_good_jets = mask_presel & mask_lepton_cleaning & mask_jetpuid
 
     elif Jet == "FatJet":
-        raise NotImplementedError()
+        # Apply the msd and preselection cuts
+        mask_msd = (events.FatJet.msoftdrop > cuts["msd"])
+        mask_good_jets = mask_presel & mask_msd
 
-    return jets[good_jets_mask], good_jets_mask
+    return jets[mask_good_jets], mask_good_jets
 
 
 def btagging(Jet, btag):
-- 
GitLab


From 6fb0e4bd7757fc5ffd47ed7a94516c4f75ecb5ca Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Fri, 24 Mar 2023 02:43:50 +0100
Subject: [PATCH 09/11] Define FatJetGood

---
 pocket_coffea/utils/configurator.py             | 12 ++++++++++++
 pocket_coffea/workflows/tthbb_base_processor.py |  4 ++++
 2 files changed, 16 insertions(+)

diff --git a/pocket_coffea/utils/configurator.py b/pocket_coffea/utils/configurator.py
index a3c626a1..8b5bac46 100644
--- a/pocket_coffea/utils/configurator.py
+++ b/pocket_coffea/utils/configurator.py
@@ -209,6 +209,18 @@ class Configurator:
             self.subsamples_list += sam
         self.total_samples_list = list(set(self.samples + self.subsamples_list))
 
+        for key, subscfg in self.subsamples_cuts.items():
+            if isinstance(subscfg, dict):
+                # Convert it to StandardSelection
+                self.subsamples[key] = StandardSelection(subscfg)
+            elif isinstance(subscfg, StandardSelection):
+                self.subsamples[key] = subscfg
+            elif isinstance(subscfg, CartesianSelection):
+                self.subsamples[key] = subscfg
+        # Unique set of cuts
+        logging.info("Subsamples:")
+        logging.info(self.subsamples)
+
     def filter_dataset(self, nfiles):
         filtered_dataset = {}
         for sample, ds in self.fileset.items():
diff --git a/pocket_coffea/workflows/tthbb_base_processor.py b/pocket_coffea/workflows/tthbb_base_processor.py
index 8602188b..5c035380 100644
--- a/pocket_coffea/workflows/tthbb_base_processor.py
+++ b/pocket_coffea/workflows/tthbb_base_processor.py
@@ -66,6 +66,10 @@ class ttHbbBaseProcessor(BaseProcessorABC):
                 self.events.ElectronGood, self.events.MuonGood
             )
 
+        self.events["FatJetGood"], self.jetGoodMask = jet_selection(
+            self.events, "FatJet", self.cfg.finalstate
+        )
+
     def count_objects(self, variation):
         self.events["nMuonGood"] = ak.num(self.events.MuonGood)
         self.events["nElectronGood"] = ak.num(self.events.ElectronGood)
-- 
GitLab


From d90c609c837ff6465cf19906ac133b25e27a4055 Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Fri, 24 Mar 2023 02:52:28 +0100
Subject: [PATCH 10/11] Treat subsamples as StandardSelection or
 CartesianSelection

---
 config/test_subsamples.py           | 107 +++++++++++++++++-----------
 pocket_coffea/lib/hist_manager.py   |   4 +-
 pocket_coffea/utils/configurator.py |  39 +++++-----
 pocket_coffea/workflows/base.py     |  66 +++++------------
 4 files changed, 105 insertions(+), 111 deletions(-)

diff --git a/config/test_subsamples.py b/config/test_subsamples.py
index a75afa41..cdc64b6f 100644
--- a/config/test_subsamples.py
+++ b/config/test_subsamples.py
@@ -1,12 +1,65 @@
 from pocket_coffea.parameters.cuts.preselection_cuts import *
 from pocket_coffea.workflows.tthbb_base_processor import ttHbbBaseProcessor
-from pocket_coffea.lib.cut_functions import get_nObj_min, get_nObj_eq, get_nBtag, get_HLTsel, get_nObj_less
+from pocket_coffea.lib.cut_functions import get_nObj_min, get_nObj_eq, get_HLTsel, get_nObj_less
 from pocket_coffea.parameters.histograms import *
 from pocket_coffea.lib.columns_manager import ColOut
 from pocket_coffea.parameters.btag import btag_variations
-from pocket_coffea.lib.weights_manager import WeightCustom
-import numpy as np
 
+def ptmsd(events, params, **kwargs):
+    # Mask to select events with a fatjet with minimum softdrop mass and maximum tau21
+    #return (events.FatJetGood[:,0].pt > params["pt"]) & (events.FatJetGood[:,0].msoftdrop > params["msd"])
+    mask = (events.FatJetGood.pt > params["pt"]) & (events.FatJetGood.msoftdrop > params["msd"])
+
+    assert not ak.any(ak.is_none(mask, axis=1)), f"None in ptmsd\n{events.FatJetGood.pt[ak.is_none(mask, axis=1)]}"
+
+    return ak.where(~ak.is_none(mask, axis=1), mask, False)
+
+def get_ptmsd(pt, msd, name=None):
+    if name == None:
+        name = f"pt{pt}msd{msd}"
+    return Cut(
+        name=name,
+        params= {"pt" : pt, "msd" : msd},
+        function=ptmsd,
+        collection="FatJetGood"
+    )
+
+def flavor_mask(events, params, **kwargs):
+    mask = {
+        "l"  : events.FatJetGood.hadronFlavour < 4,
+        "c"  : events.FatJetGood.hadronFlavour == 4,
+        "b"  : events.FatJetGood.hadronFlavour == 5,
+        "cc" : abs(events.FatJetGood.hadronFlavour == 4) & (events.FatJetGood.nBHadrons == 0) & (events.FatJetGood.nCHadrons >= 2),
+        "bb" : abs(events.FatJetGood.hadronFlavour == 5) & (events.FatJetGood.nBHadrons >= 2)
+    }
+
+    if params["flavor"] in ["bb", "cc"]:
+        return mask[params["flavor"]]
+    elif params["flavor"] == "b":
+        return mask[params["flavor"]] & ~mask["bb"]
+    elif params["flavor"] == "c":
+        return mask[params["flavor"]] & ~mask["cc"]
+    elif params["flavor"] == "l":
+        return mask[params["flavor"]] & ~mask["bb"] & ~mask["cc"] & ~mask["b"] & ~mask["c"]
+    else:
+        raise NotImplementedError
+
+def get_flavor(flavor):
+    return Cut(
+        name=flavor,
+        params={"flavor": flavor},
+        function=flavor_mask,
+        collection="FatJetGood"
+    )
+
+samples = ["TTToSemiLeptonic",
+           "ttHTobb",
+           "DATA_SingleMuon"
+          ]
+
+subsamples = {}
+for s in filter(lambda x: 'DATA' not in x, samples):
+    subsamples[s] = {f"{s}_{f}" : [get_flavor(f)] for f in ['l', 'c', 'b', 'cc', 'bb']}
 
 cfg =  {
     "dataset" : {
@@ -16,35 +69,17 @@ cfg =  {
                   "datasets/DATA_SingleMuon_local.json",
                   "datasets/DATA_SingleEle.json"],
         "filter" : {
-            "samples": [
-                # "TTToSemiLeptonic",
-                "TTbbSemiLeptonic",
-                "ttHTobb",
-                "DATA_SingleMu",
-                # "DATA_SingleEle"
-            ],
+            "samples": samples,
             "samples_exclude" : [],
             "year": ['2018']
         },
-        "subsamples": {
-            "TTbbSemiLeptonic":
-            {
-                "tt<3b": [get_nObj_less(3, coll="BJetGood")],
-                "tt+3b": [get_nObj_eq(3, coll="BJetGood")],
-                "tt+4b": [get_nObj_eq(4, coll="BJetGood")],
-                "tt>4b": [get_nObj_min(5, coll="BJetGood")],
-            },
-            "ttHTobb":{
-                "ttH<4b": [get_nObj_less(4, coll="BJetGood")],
-                "ttH+4b": [get_nObj_min(4, coll="BJetGood")],
-            }
-        }
+        "subsamples": subsamples
     },
     
 
     # Input and output files
     "workflow" : ttHbbBaseProcessor,
-    "output"   : "output/test_subsamples",
+    "output"   : "output/test/test_subsamples",
     "worflow_options" : {},
 
     "run_options" : {
@@ -67,16 +102,15 @@ cfg =  {
     },
 
     # Cuts and plots settings
-    "finalstate" : "semileptonic",
-    "skim": [get_nObj_min(4, 15., "Jet"),
+    "finalstate" : "mutag",
+    "skim": [get_nObj_min(1, 400., "FatJet"),
              get_HLTsel("semileptonic")],
-    "preselections" : [semileptonic_presel],
+    "preselections" : [semileptonic_presel, get_nObj_min(1, 450., "FatJetGood")],
     "categories": {
         "baseline": [passthrough],
+        "pt450msd40" : [get_ptmsd(450., 40.)],
     },
 
-    
-
     "weights": {
         "common": {
             "inclusive": ["genWeight","lumi","XS",
@@ -111,23 +145,16 @@ cfg =  {
 
    "variables":
     {
-
+        **fatjet_hists(coll="FatJetGood", fields=["pt", "eta", "phi",]),
+        **fatjet_hists(coll="FatJetGood", fields=["pt", "eta", "phi",], pos=0),
         **ele_hists(coll="ElectronGood", pos=0),
         **muon_hists(coll="MuonGood", pos=0),
     },
 
     "columns": {
         "common": {
-            "inclusive": [ ColOut("JetGood", ["pt", "eta","phi"]),
-                           ColOut("ElectronGood", ["pt","eta","phi"])]
+            "inclusive": [ ColOut("FatJetGood", ["pt", "eta", "phi"])]
         },
-        "bysample":{
-            "ttHTobb": {
-                "inclusive": [ColOut("MuonGood", ["pt", "eta", "phi"])]
-            },
-            "ttH+4b":{
-                "inclusive": [ColOut("BJetGood", ["pt", "eta", "phi"])]
-            }
-        }
+        "bysample" : { subs : { "inclusive" : [ColOut("FatJetGood", ["hadronFlavour", "nBHadrons", "nCHadrons"])] } for subs in subsamples }
     }
 }
diff --git a/pocket_coffea/lib/hist_manager.py b/pocket_coffea/lib/hist_manager.py
index 5e1dd11d..4f65ad18 100644
--- a/pocket_coffea/lib/hist_manager.py
+++ b/pocket_coffea/lib/hist_manager.py
@@ -258,7 +258,7 @@ class HistManager:
         weights_manager,
         categories,
         shape_variation="nominal",
-        subsamples_masks=None,  # This is a dictionary with name:ak.Array(bool)
+        subsamples=None,  # This is a dictionary with name:ak.Array(bool)
         custom_fields=None,
     ):
         '''
@@ -367,7 +367,7 @@ class HistManager:
             # Mask the events, the weights and then flatten and remove the None correctly
             for category, cat_mask in categories.get_masks():
                 # loop directly on subsamples
-                for subsample, subs_mask in subsamples_masks.items():
+                for subsample, subs_mask in subsamples.get_masks():
                     # logging.info(f"\t\tcategory {category}, subsample {subsample}")
                     mask = cat_mask & subs_mask
                     # Skip empty categories and subsamples
diff --git a/pocket_coffea/utils/configurator.py b/pocket_coffea/utils/configurator.py
index 8b5bac46..454b0741 100644
--- a/pocket_coffea/utils/configurator.py
+++ b/pocket_coffea/utils/configurator.py
@@ -186,37 +186,38 @@ class Configurator:
 
     def load_subsamples(self):
         # subsamples configuration
-        # Dictionary with subsample_name:[list of Cut ids]
-        self.subsamples_cuts = {
-            key: subscfg
-            for key, subscfg in self.dataset.get("subsamples", {}).items()
-            if key in self.samples
-        }
-        # Map of subsamples names
+        subsamples_dict = self.dataset.get("subsamples", {})
         self.subsamples = {}
+        self.subsamples_names = {}
         self.has_subsamples = {}
+        # Save list of subsamples for each sample
         for sample in self.samples:
-            if sample in self.subsamples_cuts:
-                self.subsamples[sample] = list(self.subsamples_cuts[sample].keys())
+            if sample in subsamples_dict.keys():
+                self.subsamples_names[sample] = list(subsamples_dict[sample].keys())
                 self.has_subsamples[sample] = True
             else:
-                self.subsamples[sample] = []
                 self.has_subsamples[sample] = False
 
         # Complete list of samples and subsamples
         self.subsamples_list = []
-        for sam in self.subsamples.values():
+        for sam in self.subsamples_names.values():
             self.subsamples_list += sam
         self.total_samples_list = list(set(self.samples + self.subsamples_list))
 
-        for key, subscfg in self.subsamples_cuts.items():
-            if isinstance(subscfg, dict):
-                # Convert it to StandardSelection
-                self.subsamples[key] = StandardSelection(subscfg)
-            elif isinstance(subscfg, StandardSelection):
-                self.subsamples[key] = subscfg
-            elif isinstance(subscfg, CartesianSelection):
-                self.subsamples[key] = subscfg
+        # Now saving the subsamples
+        for sample in self.samples:
+            if sample in subsamples_dict:
+                subscfg = subsamples_dict[sample]
+                if isinstance(subscfg, dict):
+                    # Convert it to StandardSelection
+                    self.subsamples[sample] = StandardSelection(subscfg)
+                elif isinstance(subscfg, StandardSelection):
+                    self.subsamples[sample] = subscfg
+                elif isinstance(subscfg, CartesianSelection):
+                    self.subsamples[sample] = subscfg
+            else:
+                # if there is no configured subsample, the full sample becomes its subsample
+                self.subsamples[sample] = StandardSelection({sample: [passthrough]})
         # Unique set of cuts
         logging.info("Subsamples:")
         logging.info(self.subsamples)
diff --git a/pocket_coffea/workflows/base.py b/pocket_coffea/workflows/base.py
index bbcbc9fd..3d6f4d36 100644
--- a/pocket_coffea/workflows/base.py
+++ b/pocket_coffea/workflows/base.py
@@ -55,7 +55,7 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
         self._categories = self.cfg.categories
 
         # Subsamples configurations: special cuts to split a sample in subsamples
-        self._subsamplesCfg = self.cfg.subsamples_cuts
+        self._subsamples = self.cfg.subsamples
 
         # Weights configuration
         self.weights_config_allsamples = self.cfg.weights_config
@@ -116,15 +116,7 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
             self._era = self.events.metadata["era"]
             self._goldenJSON = goldenJSON[self._year]
         # Loading metadata for subsamples
-        if self._sample in self._subsamplesCfg:
-            self._subsamples = self._subsamplesCfg[self._sample]
-            self._subsamples_names = list(self._subsamples.keys())
-            self._hasSubsamples = True
-        else:
-            # if there is no configured subsample, the full sample becomes its subsample
-            self._subsamples = {self._sample: [passthrough]}
-            self._subsamples_names = [self._sample]
-            self._hasSubsamples = False
+        self._hasSubsamples = self.cfg.has_subsamples[self._sample]
 
     def load_metadata_extra(self):
         '''
@@ -276,28 +268,12 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
             isMC=self._isMC,
         )
 
-        # Defining the subsamples cut
-        # saving all the cuts in a single selector
-        self._subsamples_masks = PackedSelection()
-        self._subsamples_cuts_ids = []
-        # saving the map of cut ids for each subsample
-        self._subsamples_map = defaultdict(list)
-        self._subsamples_masks = PackedSelection()
-
-        for subs, cuts in self._subsamples.items():
-            for cut in cuts:
-                if cut.id not in self._subsamples_cuts_ids:
-                    self._subsamples_masks.add(
-                        cut.id,
-                        cut.get_mask(
-                            self.events,
-                            year=self._year,
-                            sample=self._sample,
-                            isMC=self._isMC,
-                        ),
-                    )
-                    self._subsamples_cuts_ids.append(cut.id)
-                self._subsamples_map[subs].append(cut.id)
+        self._subsamples[self._sample].prepare(
+            events=self.events,
+            year=self._year,
+            sample=self._sample,
+            isMC=self._isMC,
+        )
 
     def define_common_variables_before_presel(self, variation):
         '''
@@ -377,11 +353,7 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
 
             # If subsamples are defined we also save their metadata
             if self._hasSubsamples:
-                for subs in self._subsamples_names:
-                    # Get the mask
-                    subsam_mask = self._subsamples_masks.all(
-                        *self._subsamples_map[subs]
-                    )
+                for subs, subsam_mask in self._categories.get_masks():
                     mask_withsub = mask_on_events & subsam_mask
                     self.output["cutflow"][category][subs] = ak.sum(mask_withsub)
                     if self._isMC:
@@ -410,7 +382,7 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
         self.hists_managers = HistManager(
             self.cfg.variables,
             self._sample,
-            self._subsamples_names,
+            self._subsamples[self._sample].keys(),
             self._categories,
             variations_config=self.cfg.variations_config[self._sample],
             custom_axes=self.custom_axes,
@@ -434,21 +406,17 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
         throught the HistManager.
         '''
         # Filling the autofill=True histogram automatically
-        subs_masks = {}
-        for subs in self._subsamples_names:
-            # Get the mask
-            subs_masks[subs] = self._subsamples_masks.all(*self._subsamples_map[subs])
         # Calling hist manager with the subsample masks
         self.hists_managers.fill_histograms(
             self.events,
             self.weights_manager,
             self._categories,
-            subsamples_masks=subs_masks,
+            subsamples=self._subsamples[self._sample],
             shape_variation=variation,
             custom_fields=self.custom_histogram_fields,
         )
         # Saving the output
-        for subs in self._subsamples_names:
+        for subs in self._subsamples[self._sample].keys():
             # Saving the output for each subsample
             for var, H in self.hists_managers.get_histograms(subs).items():
                 self.output["variables"][var][subs] = H
@@ -466,7 +434,7 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
         If Subsamples are defined a columnsmager is created for each of them.
         '''
         self.column_managers = {}
-        for subs in self._subsamples_names:
+        for subs in self._subsamples[self._sample].keys():
             self.column_managers[subs] = ColumnsManager(
                 self.cfg.columns[subs], subs, self._categories
             )
@@ -483,14 +451,12 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
         # TODO Fill column accumulator for different variations
         if self._hasSubsamples:
             # call the filling for each
-            for subs in self._subsamples_names:
-                # Get the mask
-                subsam_mask = self._subsamples_masks.all(*self._subsamples_map[subs])
+            for subs in self._subsamples[self._sample].keys():
                 # Calling hist manager with a subsample mask
                 self.output["columns"][subs] = self.column_managers[subs].fill_columns(
                     self.events,
                     self._categories,
-                    subsample_mask=subsam_mask,
+                    subsample_mask=self._subsamples[self._sample].get_mask(subs),
                 )
         else:
             self.output["columns"][self._sample] = self.column_managers[
@@ -658,7 +624,7 @@ class BaseProcessorABC(processor.ProcessorABC, ABC):
             # This is computed before any preselection
             self.output['sum_genweights'][self._sample] = ak.sum(self.events.genWeight)
             if self._hasSubsamples:
-                for subs in self._subsamples_names:
+                for subs in self._subsamples[self._sample].keys():
                     self.output['sum_genweights'][subs] = self.output['sum_genweights'][
                         self._sample
                     ]
-- 
GitLab


From 4fe74a09ba85ad2a14bec6913e203e2f4c07c3a5 Mon Sep 17 00:00:00 2001
From: mmarchegiani <matmarcheg@gmail.com>
Date: Fri, 24 Mar 2023 03:24:02 +0100
Subject: [PATCH 11/11] Undo changes inherited from other branches

---
 .../parameters/object_preselection.py         |  3 +
 pocket_coffea/utils/plot_utils.py             | 81 ++----------------
 pocket_coffea/workflows/genweights.py         | 57 +------------
 scripts/README.md                             | 85 -------------------
 scripts/hadd_skimmed_files.py                 |  9 +-
 5 files changed, 12 insertions(+), 223 deletions(-)
 delete mode 100644 scripts/README.md

diff --git a/pocket_coffea/parameters/object_preselection.py b/pocket_coffea/parameters/object_preselection.py
index 291f816d..e93c0453 100644
--- a/pocket_coffea/parameters/object_preselection.py
+++ b/pocket_coffea/parameters/object_preselection.py
@@ -89,6 +89,9 @@ object_preselection = {
             "eta": 2.4,
             "msd": 40,
             "jetId": 2,
+            "nsubjet": 2,
+            "nmusj": 1,
+            "dimuon_pt_ratio": 0.6,
         },
     },
 }
diff --git a/pocket_coffea/utils/plot_utils.py b/pocket_coffea/utils/plot_utils.py
index 46d3081d..00c0d4ab 100644
--- a/pocket_coffea/utils/plot_utils.py
+++ b/pocket_coffea/utils/plot_utils.py
@@ -468,9 +468,6 @@ def plot_data_mc_hist1D(
     samples = h.keys()
     samples_data = list(filter(lambda d: 'DATA' in d, samples))
     samples_mc = list(filter(lambda d: 'DATA' not in d, samples))
-    
-    print("********************")
-    print(samples_mc)
 
     h_mc = h[samples_mc[0]]
 
@@ -490,8 +487,6 @@ def plot_data_mc_hist1D(
     for year in years:
 
         for cat in categories:
-            if 'noreweight' in cat:
-                continue
             if only_cat:
                 if isinstance(only_cat, list):
                     if not cat in only_cat:
@@ -571,7 +566,7 @@ def plot_data_mc_hist1D(
                             [
                                 h[d][slicing_mc]
                                 for d in samples_mc
-                                if (d.endswith(f'_{f}') & ('GluGluH' not in d))
+                                if d.endswith(f'_{f}')
                             ]
                         )
                     )
@@ -583,7 +578,7 @@ def plot_data_mc_hist1D(
                             [
                                 h[d][slicing_mc_nominal]
                                 for d in samples_mc
-                                if (d.endswith(f'_{f}') & ('GluGluH' not in d))
+                                if d.endswith(f'_{f}')
                             ]
                         )
                     )
@@ -593,57 +588,6 @@ def plot_data_mc_hist1D(
                 nevents = {
                     f: round(sum(dict_mc_nominal[f].values()), 1) for f in flavors
                 }
-                print("************")
-                print(histname)
-                print(samples_mc)
-                dict_gghbb = {
-                    f: stack_sum(
-                        hist.Stack.from_iter(
-                            [
-                                h[d][slicing_mc]
-                                for d in samples_mc
-                                if (d.endswith(f'_{f}') & ('GluGluHToBB' in d))
-                            ]
-                        )
-                    )
-                    for f in flavors
-                }
-                dict_gghbb_nominal = {
-                    f: stack_sum(
-                        hist.Stack.from_iter(
-                            [
-                                h[d][slicing_mc_nominal]
-                                for d in samples_mc
-                                if (d.endswith(f'_{f}') & ('GluGluHToBB' in d))
-                            ]
-                        )
-                    )
-                    for f in flavors
-                }
-                dict_gghcc = {
-                    f: stack_sum(
-                        hist.Stack.from_iter(
-                            [
-                                h[d][slicing_mc]
-                                for d in samples_mc
-                                if (d.endswith(f'_{f}') & ('GluGluHToCC' in d))
-                            ]
-                        )
-                    )
-                    for f in flavors
-                }
-                dict_gghcc_nominal = {
-                    f: stack_sum(
-                        hist.Stack.from_iter(
-                            [
-                                h[d][slicing_mc_nominal]
-                                for d in samples_mc
-                                if (d.endswith(f'_{f}') & ('GluGluHToCC' in d))
-                            ]
-                        )
-                    )
-                    for f in flavors
-                }
             else:
                 dict_mc = {d: h[d][slicing_mc] for d in samples_mc}
                 dict_mc_nominal = {d: h[d][slicing_mc_nominal] for d in samples_mc}
@@ -672,15 +616,6 @@ def plot_data_mc_hist1D(
                     dict_mc_nominal[sample] = histo_reweighted
             stack_mc = hist.Stack.from_dict(dict_mc)
             stack_mc_nominal = hist.Stack.from_dict(dict_mc_nominal)
-            if flavorsplit == '5f':
-                stack_gghbb = hist.Stack.from_dict(dict_gghbb)
-                stack_gghbb_nominal = hist.Stack.from_dict(dict_gghbb_nominal)
-                stack_gghcc = hist.Stack.from_dict(dict_gghcc)
-                stack_gghcc_nominal = hist.Stack.from_dict(dict_gghcc_nominal)
-                hist_gghbb = stack_sum(stack_gghbb)
-                hist_gghbb_nominal = stack_sum(stack_gghbb_nominal)
-                hist_gghcc = stack_sum(stack_gghcc)
-                hist_gghcc_nominal = stack_sum(stack_gghcc_nominal)
 
             if not is_mc_only:
                 # Sum over eras if era axis exists in data histogram
@@ -731,12 +666,6 @@ def plot_data_mc_hist1D(
                     x, h_data.values(), yerr=np.sqrt(h_data.values()), **opts_data
                 )
                 # stack_data.plot(stack=True, color='black', ax=ax)
-            if flavorsplit == '5f':
-                sf = 100
-                hist_gghbb_nominal = sf * hist_gghbb_nominal
-                hist_gghcc_nominal = sf * hist_gghcc_nominal
-                hist_gghbb_nominal.plot(stack=False, histtype='step', ax=ax, color='red', label=f'ggHbb x {sf}')
-                hist_gghcc_nominal.plot(stack=False, histtype='step', ax=ax, color='green', label=f'ggHcc x {sf}')
             if rebinning:
                 syst_err2_up, syst_err2_down = get_systematic_uncertainty(
                     stack_mc,
@@ -799,7 +728,7 @@ def plot_data_mc_hist1D(
                 exp = math.floor(
                     math.log(max(stack_sum(stack_mc_nominal).values()), 10)
                 )
-                ax.set_ylim((0.01, 10 ** (exp + 3)))
+                ax.set_ylim((0.01, 10 ** (exp + 2)))
                 # if flavorsplit == '5f':
                 #    ax.legend(handles, labels, loc="upper right", fontsize=fontsize, ncols=2)
                 # else:
@@ -818,7 +747,7 @@ def plot_data_mc_hist1D(
                             handles_new.append(handles[i])
                         labels = labels_new
                         handles = handles_new
-                        ax.legend(handles, labels, fontsize=fontsize, ncols=2, loc="upper right")
+                        ax.legend(handles, labels, fontsize=fontsize, ncols=2)
                     if ("variables" in config.plot_options) & (histname in config.plot_options["variables"].keys()):
                         cfg_plot = config.plot_options["variables"][histname]
                         if 'xlim' in cfg_plot.keys():
@@ -831,7 +760,7 @@ def plot_data_mc_hist1D(
                                         max(stack_sum(stack_mc_nominal).values()), 10
                                     )
                                 )
-                                ax.set_ylim((0.01, 10 ** (exp + 3)))
+                                ax.set_ylim((0.01, 10 ** (exp + 2)))
                                 # ax.legend(handles, labels, loc="upper right", fontsize=fontsize, ncols=2)
                         if 'ylim' in cfg_plot.keys():
                             if isinstance(cfg_plot['ylim'], tuple):
diff --git a/pocket_coffea/workflows/genweights.py b/pocket_coffea/workflows/genweights.py
index dc0ffab7..1a404377 100644
--- a/pocket_coffea/workflows/genweights.py
+++ b/pocket_coffea/workflows/genweights.py
@@ -1,16 +1,9 @@
-import numpy as np
 import awkward as ak
 
 import copy
 
-from coffea.lumi_tools import LumiMask
-from coffea.analysis_tools import PackedSelection
-
 from .base import BaseProcessorABC
 from ..utils.configurator import Configurator
-from ..parameters.event_flags import event_flags, event_flags_data
-from ..parameters.lumi import goldenJSON
-from ..parameters.btag import btag
 
 class genWeightsProcessor(BaseProcessorABC):
     def __init__(self, cfg: Configurator) -> None:
@@ -18,8 +11,7 @@ class genWeightsProcessor(BaseProcessorABC):
         self.output_format = {
             "sum_genweights": {},
             "cutflow": {
-                "initial": {s: 0 for s in self.cfg.datasets},
-                "skim": {s: 0 for s in self.cfg.datasets}
+                "initial": {s: 0 for s in self.cfg.datasets}
             }
         }
 
@@ -27,14 +19,9 @@ class genWeightsProcessor(BaseProcessorABC):
         self._dataset = self.events.metadata["dataset"]
         self._sample = self.events.metadata["sample"]
         self._year = self.events.metadata["year"]
-        self._btag = btag[self._year]
         self._isMC = self.events.metadata["isMC"] == "True"
         if self._isMC:
-            self._era = "MC"
             self._xsec = self.events.metadata["xsec"]
-        else:
-            self._era = self.events.metadata["era"]
-            self._goldenJSON = goldenJSON[self._year]
 
         # Check if the user specified any subsamples without performing any operation
         if self._sample in self._subsamplesCfg:
@@ -42,39 +29,6 @@ class genWeightsProcessor(BaseProcessorABC):
         else:
             self._hasSubsamples = False
 
-    def skim_events(self):
-        self._skim_masks = PackedSelection()
-        mask_flags = np.ones(self.nEvents_initial, dtype=np.bool)
-        flags = event_flags[self._year]
-        if not self._isMC:
-            flags += event_flags_data[self._year]
-        for flag in flags:
-            mask_flags = getattr(self.events.Flag, flag)
-        self._skim_masks.add("event_flags", mask_flags)
-
-        self._skim_masks.add("PVgood", self.events.PV.npvsGood > 0)
-
-        if not self._isMC:
-            mask_lumi = LumiMask(self._goldenJSON)(
-                self.events.run, self.events.luminosityBlock
-            )
-            self._skim_masks.add("lumi_golden", mask_lumi)
-
-        for skim_func in self._skim:
-            mask = skim_func.get_mask(
-                self.events,
-                year=self._year,
-                sample=self._sample,
-                btag=self._btag,
-                isMC=self._isMC,
-            )
-            self._skim_masks.add(skim_func.id, mask)
-
-        self.events = self.events[self._skim_masks.all(*self._skim_masks.names)]
-        self.nEvents_after_skim = self.nevents
-        self.output['cutflow']['skim'][self._dataset] += self.nEvents_after_skim
-        self.has_events = self.nEvents_after_skim > 0
-
     def apply_object_preselection(self, variation):
         pass
 
@@ -94,13 +48,4 @@ class genWeightsProcessor(BaseProcessorABC):
             if self._hasSubsamples:
                 raise Exception("This processor cannot compute the sum of genweights of subsamples.")
 
-        self.process_extra_before_skim()
-        self.skim_events()
-        if not self.has_events:
-            return self.output
-
-        if self.cfg.save_skimmed_files:
-            self.export_skimmed_chunk()
-            return self.output
-
         return self.output
diff --git a/scripts/README.md b/scripts/README.md
deleted file mode 100644
index 940cc5a1..00000000
--- a/scripts/README.md
+++ /dev/null
@@ -1,85 +0,0 @@
-# PocketCoffea scripts
-
-The current guide is a documentation of how to run the scripts of the PocketCoffea package.
-
-## Merge skimmed files with hadd
-
-In order to merge the skimmed n-tuples, the `hadd` command from ROOT needs to be used.
-To source ROOT inside the `pocket_coffea` environment, run the following command:
-
-```
-source /cvmfs/sft.cern.ch/lcg/views/LCG_102/x86_64-centos7-gcc11-opt/setup.sh
-```
-
-To run the `hadd` script run the following command:
-
-```
-cd /path/to/PocketCoffea
-python scripts/hadd_skimmed_files.py -fl /path/to/output/output_all.coffea -o /path/to/outputdir
-```
-
-## Build dataset
-
-The first step is to build the json files containing the list of files of the UL datasets to process. In our case, we need to include the `TTToSemiLeptonic` and `TTTo2L2Nu` MC datasets and `SingleMuon` as data dataset.
-
-A step-by-step detailed guide on the creation of the dataset file can be found at this [link](https://pocketcoffea.readthedocs.io/en/latest/examples.html).
-
-### Append sum of genweights to datasets definitions
-
-In order to update the `datasets_definitions.json` and include the `sum_genweights` metadata in the datasets definitions, a dedicated script needs to be run with the following command:
-
-```
-cd /path/to/PocketCoffea
-scripts/dataset/append_genweights.py --cfg dataset/dataset_definitions.json -i output/genweights/genweights_2016_PreVFP/output_all.coffea --overwrite
-scripts/dataset/append_genweights.py --cfg dataset/dataset_definitions.json -i output/genweights/genweights_2016_PostVFP/output_all.coffea --overwrite
-```
-
-N.B.: the argument `-i output/genweights/genweights_201*/output_all.coffea` should be the path to the Coffea output file where the genweights are stored. The option `--overwrite` automatically updates the config file passed as `--cfg` argument with the genweights, in this case `dataset/dataset_definitions.json`.
-
-### Build dataset with sum of genweights
-
-Now the steps of the [Build dataset](#build-dataset) section need to be repeated to generate the json datasets with the additional metadata `sum_genweights`.
-
-## Run the analysis
-
-In order to run the analysis workflow and produce the output histograms, run the following command:
-
-```
-cd /path/to/PocketCoffea
-runner.py --cfg /path/to/config/config.py --full
-```
-
-N.B.: the argument `--full` will process all the datasets together at once and save the output in a single output file, `output_all.coffea`. Otherwise, the datasets are processed separately and an output file is saved for each dataset.
-
-## Accumulate output files
-
-If the output of the [previous step](#run-the-analysis) has been produced without the argument `--full`, the output files need to be merged in a single output file `output_all.coffea`. If the output has been produced with the argument `--full` and the output file `output_all.coffea` is already existing, skip this step and continue with the [next one](#produce-datamc-plots).
-
-Once the Coffea output files are produced, one needs to merge the files into a single file by using the script `accumulate_files.py` by running this command:
-
-```
-cd /path/to/output
-python accumulate_files.py -i output1.coffea output2.coffea output3.coffea -o output_all.coffea
-```
-
-## Produce data/MC plots
-
-To produce plots from a Coffea output file execute the plotting script:
-
-```
-cd /path/to/PocketCoffea
-make_plots.py --cfg /path/to/config/config.py -i /path/to/output/output_all.coffea
-```
-
-## Run trigger SF script
-
-To run the script that computes the single electron trigger SF and produces the validation plots, run:
-
-```
-cd /path/to/PocketCoffea
-python scripts/plot/trigger_efficiency.py --cfg config/semileptonic_triggerSF/semileptonic_triggerSF_2017.py --save_plots
-```
-
-If the additional argument `--save_plots` is not specified, only the scale factor maps are saved without saving the plots.
-
-The output plots are saved in `/path/to/output/plots/trigger_efficiency` and `/path/to/output/plots/trigger_scalefactor`, while the 1D and 2D scale factor maps are saved in the folder specified in `workflow_options['output_triggerSF']` in the config file.
diff --git a/scripts/hadd_skimmed_files.py b/scripts/hadd_skimmed_files.py
index 6d870ba2..2280fcca 100644
--- a/scripts/hadd_skimmed_files.py
+++ b/scripts/hadd_skimmed_files.py
@@ -18,7 +18,7 @@ parser.add_argument(
 )
 parser.add_argument("-o", "--outputdir", required=False, type=str, help="Output folder")
 parser.add_argument(
-    "--only-datasets", nargs="+", type=str, help="Restricting list of datasets"
+    "--only-datasets", nargs="+", type=str, help="Restring list of datasets"
 )
 parser.add_argument("-f", "--files", type=int, help="Limit number of files")
 parser.add_argument("-e", "--events", type=int, help="Limit number of files")
@@ -43,9 +43,6 @@ groups_metadata = {}
 for dataset in df["skimmed_files"].keys():
     if args.only_datasets and dataset not in args.only_datasets:
         continue
-    outputdir_dataset = os.path.join(args.outputdir, dataset)
-    if not os.path.exists(outputdir_dataset):
-        os.makedirs(outputdir_dataset)
     groups_metadata[dataset] = defaultdict(dict)
     nevents_tot = 0
     nfiles = 0
@@ -57,7 +54,7 @@ for dataset in df["skimmed_files"].keys():
         if (args.files and (nfiles + 1) > args.files) or (
             args.events and (nevents_tot + nevents) > args.events
         ):
-            outfile = os.path.join(outputdir_dataset, f"{dataset}_{ngroup}.root")
+            outfile = f"{args.outputdir}/{dataset}/{dataset}_{ngroup}.root"
             workload.append((outfile, group[:]))
             groups_metadata[dataset]["files"][outfile] = group[:]
             group.clear()
@@ -71,7 +68,7 @@ for dataset in df["skimmed_files"].keys():
 
     # add last group
     if len(group):
-        outfile = os.path.join(outputdir_dataset, f"{dataset}_{ngroup}.root")
+        outfile = f"{args.outputdir}/{dataset}/{dataset}_{ngroup}.root"
         workload.append((outfile, group[:]))
         groups_metadata[dataset]["files"][outfile] = group[:]
 
-- 
GitLab