exploratory_data_analysis.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:22:43.828531Z",
     "iopub.status.busy": "2023-05-29T17:22:43.828145Z",
     "iopub.status.idle": "2023-05-29T17:22:45.149390Z",
     "shell.execute_reply": "2023-05-29T17:22:45.148594Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np \n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt \n",
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:22:45.153899Z",
     "iopub.status.busy": "2023-05-29T17:22:45.153304Z",
     "iopub.status.idle": "2023-05-29T17:38:57.951193Z",
     "shell.execute_reply": "2023-05-29T17:38:57.950543Z"
    }
   },
   "outputs": [],
   "source": [
    "import uproot\n",
    "from pathlib import Path\n",
    "from tqdm.auto import tqdm\n",
    "\n",
    "# load utils\n",
    "df_pmt_id_conversion = pd.read_csv(\"/home/ferracci/new_dataset/utils/PMT_ID_conversion.csv\")\n",
    "df_pmt_position = pd.read_csv(\"/home/ferracci/new_dataset/utils/PMTPos_CD_LPMT.csv\")\n",
    "df_spmt_id_conversion = pd.read_csv(\"/home/ferracci/new_dataset/utils/SPMT_ID_conversion.csv\")\n",
    "df_spmt_position = pd.read_csv(\"/home/ferracci/new_dataset/utils/PMTPos_CD_SPMT.csv\")\n",
    "\n",
    "# computes PMTs positions \n",
    "x = np.array(df_pmt_position['x']/1000).reshape((-1, 1))\n",
    "y = np.array(df_pmt_position['y']/1000).reshape((-1, 1))\n",
    "z = np.array(df_pmt_position['z']/1000).reshape((-1, 1))\n",
    "pos = np.hstack((x, y, z))\n",
    "\n",
    "x_s = np.array(df_spmt_position['x']/1000).reshape((-1, 1))\n",
    "y_s = np.array(df_spmt_position['y']/1000).reshape((-1, 1))\n",
    "z_s = np.array(df_spmt_position['z']/1000).reshape((-1, 1))\n",
    "pos_s = np.hstack((x_s, y_s, z_s))\n",
    "\n",
    "pmt_id_raw_to_id_map = dict(zip(df_pmt_id_conversion['CdID'], df_pmt_id_conversion['PMTID']))\n",
    "pmt_id_to_pos_map = dict(zip(df_pmt_position['PMTID'], pos))\n",
    "\n",
    "spmt_id_raw_to_id_map = dict(zip(df_spmt_id_conversion['CdID'], df_spmt_id_conversion['SPMTID']))\n",
    "spmt_id_to_pos_map = dict(zip(df_spmt_position['SPMTID'], pos_s))\n",
    "\n",
    "files = list(Path(\"/mnt/data/train_flat_dataset_processed/\").glob(\"*\"))\n",
    "\n",
    "pmt_pos, charge, fht = np.array([]), np.array([]), np.array([])\n",
    "spmt_pos, charge_s, fht_s = np.array([]), np.array([]), np.array([])\n",
    "targets_dataframe = []\n",
    "\n",
    "# compute PMTs positions, charges and first-hit-times for roughly 500k events\n",
    "for filename in tqdm(files[:125]):\n",
    "    hits = np.load(filename, allow_pickle=True)[\"hits\"]\n",
    "    targets = np.load(filename, allow_pickle=True)[\"primaries\"]\n",
    "\n",
    "    pmt_id_raw, b, c = hits[0, :], hits[1, :], hits[2, :]\n",
    "    spmt_id_raw, e, f = hits[3, :], hits[4, :], hits[5, :]\n",
    "\n",
    "    # maps raw ids (CdID) to standard ids (PMTID)\n",
    "    pmt_id = np.array([np.array([pmt_id_raw_to_id_map[n] for n in event]) for event in pmt_id_raw], dtype=object)\n",
    "    spmt_id = np.array([np.array([spmt_id_raw_to_id_map[n] for n in event]) for event in spmt_id_raw], dtype=object)\n",
    "\n",
    "    # maps standard ids (PMTID) to PMTs positions\n",
    "    a = np.array([np.array([pmt_id_to_pos_map[n] for n in event]) for event in pmt_id], dtype=object)\n",
    "    d = np.array([np.array([spmt_id_to_pos_map[n] for n in event]) for event in spmt_id], dtype=object)\n",
    "\n",
    "    # append to numpy arrays\n",
    "    pmt_pos, charge, fht = np.append(pmt_pos, a), np.append(charge, b), np.append(fht, c)\n",
    "    spmt_pos, charge_s, fht_s = np.append(spmt_pos, d), np.append(charge_s, e), np.append(fht_s, f)\n",
    "\n",
    "    # retrieve targets_dataframes\n",
    "    t_df = pd.DataFrame({\"Qedep\": targets[:, 1], \"Redep\": targets[:, 5]})\n",
    "    targets_dataframe.append(t_df)\n",
    "\n",
    "targets_dataframe = pd.concat(targets_dataframe, axis=0, ignore_index=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Event 3D Plot - Large PMTs and Small PMTs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:38:57.954864Z",
     "iopub.status.busy": "2023-05-29T17:38:57.954237Z",
     "iopub.status.idle": "2023-05-29T17:39:00.374378Z",
     "shell.execute_reply": "2023-05-29T17:39:00.373804Z"
    }
   },
   "outputs": [],
   "source": [
    "plt.style.use(\"dark_background\")\n",
    "plt.rc(\"text\", usetex=True)\n",
    "plt.rc(\"font\", family=\"cm\")\n",
    "plt.rcParams[\"grid.color\"] = (0.5, 0.5, 0.5, 0.2)\n",
    "\n",
    "fig, axes = plt.subplots(nrows=2, ncols=2, figsize=(11,9), subplot_kw={\"projection\": \"3d\"}, dpi=200)\n",
    "\n",
    "# charge 3d\n",
    "ax = axes[0, 0]\n",
    "\n",
    "scatter = ax.scatter(pmt_pos[5][:, 0], pmt_pos[5][:, 1], pmt_pos[5][:, 2], s=0.1, alpha=1, c=charge[5], cmap=\"YlOrRd\", vmin=0, vmax=30)\n",
    "\n",
    "ax.set_xlabel(\"$x$\", fontsize=12)\n",
    "ax.set_ylabel(\"$y$\", fontsize=12)\n",
    "ax.set_zlabel(\"$z$\", fontsize=12)\n",
    "\n",
    "cbar = fig.colorbar(scatter, ax=ax, shrink=0.5, aspect=10, ticks=[0, 10, 20, 30])\n",
    "cbar.ax.set_title(\"Charge\", fontsize=12)\n",
    "cbar.outline.set_edgecolor('w')\n",
    "cbar.outline.set_linewidth(0.5) \n",
    "\n",
    "ax.view_init(elev=10, azim=50)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=9)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=9)\n",
    "ax.set_xticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_yticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_zticks(ticks=[-20, -10, 0, 10, 20])\n",
    "\n",
    "ax.xaxis.pane.fill = False\n",
    "ax.yaxis.pane.fill = False\n",
    "ax.zaxis.pane.fill = False\n",
    "\n",
    "ax.xaxis.pane.set_edgecolor('k')\n",
    "ax.yaxis.pane.set_edgecolor('k')\n",
    "ax.zaxis.pane.set_edgecolor('k')\n",
    "\n",
    "# fht 3d\n",
    "ax = axes[0, 1]\n",
    "scatter = ax.scatter(pmt_pos[5][:, 0], pmt_pos[5][:, 1], pmt_pos[5][:, 2], s=0.1, alpha=1, c=fht[5], cmap=\"PuBuGn\", vmin=0, vmax=1000)\n",
    "\n",
    "ax.set_xlabel(\"$x$\", fontsize=12)\n",
    "ax.set_ylabel(\"$y$\", fontsize=12)\n",
    "ax.set_zlabel(\"$z$\", fontsize=12)\n",
    "\n",
    "cbar = fig.colorbar(scatter, ax=ax, shrink=0.5, aspect=10, ticks=[0, 200, 400, 600, 800, 1000])\n",
    "cbar.ax.set_title(\"FHT\", fontsize=12)\n",
    "cbar.outline.set_edgecolor('w')\n",
    "cbar.outline.set_linewidth(0.5)\n",
    "\n",
    "ax.view_init(elev=10, azim=50)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=9)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=9)\n",
    "ax.set_xticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_yticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_zticks(ticks=[-20, -10, 0, 10, 20])\n",
    "\n",
    "ax.xaxis.pane.fill = False\n",
    "ax.yaxis.pane.fill = False\n",
    "ax.zaxis.pane.fill = False\n",
    "\n",
    "ax.xaxis.pane.set_edgecolor('k')\n",
    "ax.yaxis.pane.set_edgecolor('k')\n",
    "ax.zaxis.pane.set_edgecolor('k')\n",
    "\n",
    "plt.subplots_adjust(wspace=0.1, hspace=0)\n",
    "\n",
    "# charge 3d\n",
    "ax = axes[1, 0]\n",
    "\n",
    "scatter = ax.scatter(spmt_pos[5][:, 0], spmt_pos[5][:, 1], spmt_pos[5][:, 2], s=0.5, alpha=1, c=charge_s[5], cmap=\"YlOrRd\", vmin=0, vmax=3)\n",
    "\n",
    "ax.set_xlabel(\"$x$\", fontsize=12)\n",
    "ax.set_ylabel(\"$y$\", fontsize=12)\n",
    "ax.set_zlabel(\"$z$\", fontsize=12)\n",
    "\n",
    "cbar = fig.colorbar(scatter, ax=ax, shrink=0.5, aspect=10, ticks=[0., 1., 2., 3.0])\n",
    "cbar.ax.set_title(\"Charge\", fontsize=12)\n",
    "cbar.outline.set_edgecolor('w')\n",
    "cbar.outline.set_linewidth(0.5)\n",
    "\n",
    "ax.view_init(elev=10, azim=50)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=9)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=9)\n",
    "ax.set_xticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_yticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_zticks(ticks=[-20, -10, 0, 10, 20])\n",
    "\n",
    "ax.xaxis.pane.fill = False\n",
    "ax.yaxis.pane.fill = False\n",
    "ax.zaxis.pane.fill = False\n",
    "\n",
    "ax.xaxis.pane.set_edgecolor('k')\n",
    "ax.yaxis.pane.set_edgecolor('k')\n",
    "ax.zaxis.pane.set_edgecolor('k')\n",
    "\n",
    "# fht 3d\n",
    "ax = axes[1, 1]\n",
    "scatter = ax.scatter(spmt_pos[5][:, 0], spmt_pos[5][:, 1], spmt_pos[5][:, 2], s=0.5, alpha=1, c=fht_s[5], cmap=\"PuBuGn\", vmin=0, vmax=1000)\n",
    "\n",
    "ax.set_xlabel(\"$x$\", fontsize=12)\n",
    "ax.set_ylabel(\"$y$\", fontsize=12)\n",
    "ax.set_zlabel(\"$z$\", fontsize=12)\n",
    "\n",
    "cbar = fig.colorbar(scatter, ax=ax, shrink=0.5, aspect=10, ticks=[0, 200, 400, 600, 800, 1000])\n",
    "cbar.ax.set_title(\"FHT\", fontsize=12)\n",
    "cbar.outline.set_edgecolor('w')\n",
    "cbar.outline.set_linewidth(0.5)\n",
    "\n",
    "ax.view_init(elev=10, azim=50)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=9)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=9)\n",
    "ax.set_xticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_yticks(ticks=[-20, -10, 0, 10, 20])\n",
    "ax.set_zticks(ticks=[-20, -10, 0, 10, 20])\n",
    "\n",
    "ax.xaxis.pane.fill = False\n",
    "ax.yaxis.pane.fill = False\n",
    "ax.zaxis.pane.fill = False\n",
    "\n",
    "ax.xaxis.pane.set_edgecolor('k')\n",
    "ax.yaxis.pane.set_edgecolor('k')\n",
    "ax.zaxis.pane.set_edgecolor('k')\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/scatter_event.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.4)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### $R_{cht}$ vs. $R_{cc}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:00.384409Z",
     "iopub.status.busy": "2023-05-29T17:39:00.384095Z",
     "iopub.status.idle": "2023-05-29T17:39:36.000469Z",
     "shell.execute_reply": "2023-05-29T17:39:35.999838Z"
    }
   },
   "outputs": [],
   "source": [
    "# compute center of charge and center of first hit time for all events\n",
    "x_cc, y_cc, z_cc = np.zeros(len(pmt_pos)), np.zeros(len(pmt_pos)), np.zeros(len(pmt_pos))\n",
    "for i in range(len(pmt_pos)):\n",
    "    x_cc[i] = np.sum(pmt_pos[i][:, 0] * charge[i]) / np.sum(charge[i])\n",
    "    y_cc[i] = np.sum(pmt_pos[i][:, 1] * charge[i]) / np.sum(charge[i])\n",
    "    z_cc[i] = np.sum(pmt_pos[i][:, 2] * charge[i]) / np.sum(charge[i])\n",
    "R_cc = np.sqrt(x_cc**2 + y_cc**2 + z_cc**2)\n",
    "\n",
    "x_cht, y_cht, z_cht = np.zeros(len(pmt_pos)), np.zeros(len(pmt_pos)), np.zeros(len(pmt_pos))\n",
    "for i in range(len(pmt_pos)):\n",
    "    x_cht[i] = np.sum(pmt_pos[i][:, 0] / (fht[i] + 50)) / np.sum(1 / (fht[i] + 50))\n",
    "    y_cht[i] = np.sum(pmt_pos[i][:, 1] / (fht[i] + 50)) / np.sum(1 / (fht[i] + 50))\n",
    "    z_cht[i] = np.sum(pmt_pos[i][:, 2] / (fht[i] + 50)) / np.sum(1 / (fht[i] + 50))\n",
    "R_cht = np.sqrt(x_cht**2 + y_cht**2 + z_cht**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:36.003946Z",
     "iopub.status.busy": "2023-05-29T17:39:36.003448Z",
     "iopub.status.idle": "2023-05-29T17:39:36.006937Z",
     "shell.execute_reply": "2023-05-29T17:39:36.006450Z"
    }
   },
   "outputs": [],
   "source": [
    "Redep = targets_dataframe[\"Redep\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:36.009301Z",
     "iopub.status.busy": "2023-05-29T17:39:36.008847Z",
     "iopub.status.idle": "2023-05-29T17:39:41.692442Z",
     "shell.execute_reply": "2023-05-29T17:39:41.691688Z"
    }
   },
   "outputs": [],
   "source": [
    "# scatter of R_cht against R_cc --> non linearity \n",
    "plt.style.use(\"default\")\n",
    "plt.rc(\"text\", usetex=True)\n",
    "plt.rc(\"font\", family=\"cm\")\n",
    "plt.rcParams[\"grid.color\"] = (0.5, 0.5, 0.5, 0.2)\n",
    "\n",
    "fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(7,5), dpi=200)\n",
    "\n",
    "scatter = ax.scatter(R_cc, R_cht, s=0.3, alpha=0.9, c=Redep, cmap=\"inferno\")\n",
    "\n",
    "ax.set_xlabel(\"$R_{cc}$ [m]\", fontsize=15)\n",
    "ax.set_ylabel(\"$R_{cht}$ [m]\", fontsize=15)\n",
    "\n",
    "ax.set_xlim((0, 12))\n",
    "ax.set_ylim((0, 10))\n",
    "\n",
    "cbar = fig.colorbar(scatter, ax=ax, shrink=0.8, aspect=10)\n",
    "cbar.ax.set_title(\"$R$ [m]\", fontsize=15)\n",
    "cbar.outline.set_edgecolor('k')\n",
    "cbar.outline.set_linewidth(0.5)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "\n",
    "ax.grid()\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/scatter_Rcht_vs_R_cc.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.2);"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### AccumCharge vs. $R^3$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:41.697607Z",
     "iopub.status.busy": "2023-05-29T17:39:41.696957Z",
     "iopub.status.idle": "2023-05-29T17:39:44.343083Z",
     "shell.execute_reply": "2023-05-29T17:39:44.342267Z"
    }
   },
   "outputs": [],
   "source": [
    "accum_charge = np.zeros(len(pmt_pos))\n",
    "for i in range(len(pmt_pos)):\n",
    "    accum_charge[i] = np.sum(charge[i]) + np.sum(charge_s[i])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:44.346230Z",
     "iopub.status.busy": "2023-05-29T17:39:44.345975Z",
     "iopub.status.idle": "2023-05-29T17:39:44.356573Z",
     "shell.execute_reply": "2023-05-29T17:39:44.355929Z"
    }
   },
   "outputs": [],
   "source": [
    "targets_dataframe[\"average_charge\"] = accum_charge / targets_dataframe[\"Qedep\"]\n",
    "targets_cut = targets_dataframe[(targets_dataframe[\"average_charge\"] > 1600) & (targets_dataframe[\"average_charge\"] < 2800)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:44.359487Z",
     "iopub.status.busy": "2023-05-29T17:39:44.358861Z",
     "iopub.status.idle": "2023-05-29T17:39:45.049255Z",
     "shell.execute_reply": "2023-05-29T17:39:45.048555Z"
    }
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(7,5), dpi=200)\n",
    "\n",
    "hist2d = ax.hist2d(np.array(targets_cut[\"Redep\"])**3, np.array(targets_cut[\"average_charge\"]), bins=100, cmap=\"inferno\")\n",
    "\n",
    "ax.set_xlabel(\"$R^3$ [$\\mathrm{m}^3$]\", fontsize=15)\n",
    "ax.set_ylabel(\"Accumulated charge per MeV [PE/MeV]\", fontsize=15)\n",
    "\n",
    "cbar = fig.colorbar(hist2d[-1], ax=ax, shrink=0.8, aspect=10)\n",
    "cbar.ax.set_title(\"Count\", fontsize=15)\n",
    "cbar.outline.set_edgecolor('k')\n",
    "cbar.outline.set_linewidth(0.5)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/scatter_charge_vs_R3.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.2);"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Charge and FHT distributions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:45.052204Z",
     "iopub.status.busy": "2023-05-29T17:39:45.051519Z",
     "iopub.status.idle": "2023-05-29T17:39:45.696873Z",
     "shell.execute_reply": "2023-05-29T17:39:45.696181Z"
    }
   },
   "outputs": [],
   "source": [
    "# visualize average pdf for fht distribution with errorbars\n",
    "test = targets_dataframe[(targets_dataframe[\"Qedep\"] > 1.9) & (targets_dataframe[\"Qedep\"] < 2.1) & \n",
    "                         (targets_dataframe[\"Redep\"] > 15.9) & (targets_dataframe[\"Redep\"] < 16.1)]\n",
    "fht_E2_R16 = fht[test.index.values]\n",
    "\n",
    "# remove artificial bar at fht=0 due to how the trigger algorithm works\n",
    "for i in range(len(fht_E2_R16)):\n",
    "    fht_E2_R16[i] = fht_E2_R16[i][fht_E2_R16[i] != 0]\n",
    "\n",
    "n_bins = 80\n",
    "hist = np.zeros((len(fht_E2_R16), n_bins))\n",
    "hist[0, :], bins = np.histogram(fht_E2_R16[0], bins=n_bins, range=(0, 1200))\n",
    "for i in range(len(fht_E2_R16)-1):\n",
    "    hist[i+1, :], _ = np.histogram(fht_E2_R16[i+1], bins=bins)\n",
    "mean_hist = np.mean(hist, axis=0)\n",
    "std_hist = np.std(hist, axis=0)\n",
    "\n",
    "fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(7,5), dpi=200)\n",
    "\n",
    "ax.bar(bins[:-1], mean_hist, width=(bins[1]-bins[0])/1.25, color=\"purple\")\n",
    "ax.errorbar(bins[:-1], mean_hist, yerr=std_hist, fmt=\".\", markersize=0, capsize=1.5, ecolor=\"k\", elinewidth=0.75)\n",
    "\n",
    "ax.set_xlabel(\"t [ns]\", fontsize=15)\n",
    "ax.set_ylabel(\"Number of fired channels per 15 ns\", fontsize=15)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "\n",
    "ax.grid()\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/cdf_fht_error.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.2);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:45.699763Z",
     "iopub.status.busy": "2023-05-29T17:39:45.699066Z",
     "iopub.status.idle": "2023-05-29T17:39:45.705885Z",
     "shell.execute_reply": "2023-05-29T17:39:45.705279Z"
    }
   },
   "outputs": [],
   "source": [
    "# let's now compare the distributions for different energies at the center of the detector\n",
    "df = targets_dataframe[targets_dataframe[\"Redep\"] < 3]\n",
    "E_L = df[(df[\"Qedep\"] > 0.9) & (df[\"Qedep\"] < 1.1)]\n",
    "E_C = df[(df[\"Qedep\"] > 4.9) & (df[\"Qedep\"] < 5.1)]\n",
    "E_R = df[(df[\"Qedep\"] > 8.9) & (df[\"Qedep\"] < 9.1)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:45.708573Z",
     "iopub.status.busy": "2023-05-29T17:39:45.707997Z",
     "iopub.status.idle": "2023-05-29T17:39:45.713960Z",
     "shell.execute_reply": "2023-05-29T17:39:45.713406Z"
    }
   },
   "outputs": [],
   "source": [
    "charge_E_L, fht_E_L = charge[E_L.index.values], fht[E_L.index.values]\n",
    "charge_E_C, fht_E_C = charge[E_C.index.values], fht[E_C.index.values]\n",
    "charge_E_R, fht_E_R = charge[E_R.index.values], fht[E_R.index.values]\n",
    "\n",
    "for i in range(len(fht_E_L)):\n",
    "    fht_E_L[i] = fht_E_L[i][fht_E_L[i] != 0]\n",
    "for i in range(len(fht_E_C)):\n",
    "    fht_E_C[i] = fht_E_C[i][fht_E_C[i] != 0]\n",
    "for i in range(len(fht_E_R)):\n",
    "    fht_E_R[i] = fht_E_R[i][fht_E_R[i] != 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:45.716700Z",
     "iopub.status.busy": "2023-05-29T17:39:45.716103Z",
     "iopub.status.idle": "2023-05-29T17:39:45.720332Z",
     "shell.execute_reply": "2023-05-29T17:39:45.719578Z"
    }
   },
   "outputs": [],
   "source": [
    "def get_mean_hist(array, xlim, n_bins):\n",
    "    hist = np.zeros((len(array), n_bins))\n",
    "    hist[0, :], bins = np.histogram(array[0], bins=n_bins, range=xlim)\n",
    "    for i in range(len(array)-1):\n",
    "        hist[i+1, :], _ = np.histogram(array[i+1], bins=bins)\n",
    "    mean_hist = np.mean(hist, axis=0)\n",
    "    return mean_hist, bins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:45.722742Z",
     "iopub.status.busy": "2023-05-29T17:39:45.722551Z",
     "iopub.status.idle": "2023-05-29T17:39:45.725945Z",
     "shell.execute_reply": "2023-05-29T17:39:45.725385Z"
    }
   },
   "outputs": [],
   "source": [
    "def plot_hist(ax, array, xlim, n_bins, color, alpha):\n",
    "    hist, bins = get_mean_hist(array, xlim, n_bins)\n",
    "    ax.bar(bins[:-1], hist, width=bins[1]-bins[0], color=color, alpha=alpha)\n",
    "    return ax"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:45.728635Z",
     "iopub.status.busy": "2023-05-29T17:39:45.728060Z",
     "iopub.status.idle": "2023-05-29T17:39:45.733313Z",
     "shell.execute_reply": "2023-05-29T17:39:45.732758Z"
    }
   },
   "outputs": [],
   "source": [
    "def plot_cdf(ax, array, xlim, n_bins, color, plot_type, legend):\n",
    "    \"\"\" \n",
    "    Displays the average cumulative distribution function plot, with a vertical line at the average mean value.\n",
    "\n",
    "    Parameters:\n",
    "        ax (matplotlib.axes): matplotlib ax on which to work\n",
    "        array (ndarray): array of objects containing samples from the distribution\n",
    "        xlim (tuple): range to use for the histogram\n",
    "        n_bins (int): number of bins to use for the histogram \n",
    "        color (string): specify color for the plot\n",
    "        plot_type (string): \"R\" for plot at different radii, \"E\" for plot at different energies\n",
    "        legend (float): average radius/energy value (depending on plot_type) to display in the legend\n",
    "\n",
    "    Returns: \n",
    "        matplotlib.axes\n",
    "    \"\"\"\n",
    "    hist, bins = get_mean_hist(array, xlim, n_bins)\n",
    "    cumulative_hist = np.cumsum(hist)\n",
    "    cumulative_hist_norm = cumulative_hist / np.sum(hist)\n",
    "\n",
    "    if plot_type == \"R\":\n",
    "        ax.plot(bins[:-1], cumulative_hist_norm, linewidth=1, color=color, label=f\"R = {legend} m\")\n",
    "    elif plot_type == \"E\":\n",
    "        ax.plot(bins[:-1], cumulative_hist_norm, linewidth=1, color=color, label=f\"E = {legend} MeV\")\n",
    "    else:\n",
    "        raise ValueError(\"Only 'R' and 'E' can be passed to plot_type\")\n",
    "\n",
    "    mean_value = 0\n",
    "    for i in range(len(array)):\n",
    "        mean_value += np.mean(array[i])\n",
    "    mean_value /= len(array)\n",
    "    ax.axvline(mean_value, color=color, linestyle=\"dashed\", linewidth=1)\n",
    "\n",
    "    return ax"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:45.735377Z",
     "iopub.status.busy": "2023-05-29T17:39:45.735200Z",
     "iopub.status.idle": "2023-05-29T17:39:47.250138Z",
     "shell.execute_reply": "2023-05-29T17:39:47.249578Z"
    }
   },
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(8,5), dpi=150)\n",
    "\n",
    "ax1 = fig.add_axes([0, 0, 1, 1])\n",
    "plot_cdf(ax1, charge_E_L, xlim=(0, 6), n_bins=500, color=\"darkmagenta\", plot_type=\"E\", legend=round(np.mean(E_L[\"Qedep\"]), 1))\n",
    "plot_cdf(ax1, charge_E_C, xlim=(0, 6), n_bins=500, color=\"cornflowerblue\", plot_type=\"E\", legend=round(np.mean(E_C[\"Qedep\"]), 1))\n",
    "plot_cdf(ax1, charge_E_R, xlim=(0, 6), n_bins=500, color=\"seagreen\", plot_type=\"E\", legend=round(np.mean(E_R[\"Qedep\"]), 1))\n",
    "ax1.set_xlabel(\"nPE\", fontsize=15)\n",
    "ax1.set_ylabel(\"CDF(nPE)\", fontsize=15)\n",
    "ax1.set_xlim((0, 6))\n",
    "ax1.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax1.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax1.grid()\n",
    "ax1.legend(loc=(0.37, 1.01), ncol=3, frameon=False, prop={'size': 12})\n",
    "\n",
    "ax2 = fig.add_axes([0.425, 0.1, 0.55, 0.55])\n",
    "plot_hist(ax2, charge_E_R, xlim=(0, 6), n_bins=100, color=\"seagreen\", alpha=0.9)\n",
    "plot_hist(ax2, charge_E_C, xlim=(0, 6), n_bins=100, color=\"cornflowerblue\", alpha=0.9)\n",
    "plot_hist(ax2, charge_E_L, xlim=(0, 6), n_bins=100, color=\"darkmagenta\", alpha=0.9)\n",
    "ax2.set_xlim((0, 6))\n",
    "ax2.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax2.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax2.grid()\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/cdf_charge_center.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.2);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:47.253010Z",
     "iopub.status.busy": "2023-05-29T17:39:47.252496Z",
     "iopub.status.idle": "2023-05-29T17:39:48.747288Z",
     "shell.execute_reply": "2023-05-29T17:39:48.746722Z"
    }
   },
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(8,5), dpi=150)\n",
    "\n",
    "ax1 = fig.add_axes([0, 0, 1, 1])\n",
    "plot_cdf(ax1, fht_E_L, xlim=(0, 1250), n_bins=500, color=\"darkmagenta\", plot_type=\"E\", legend=round(np.mean(E_L[\"Qedep\"]), 1))\n",
    "plot_cdf(ax1, fht_E_C, xlim=(0, 1250), n_bins=500, color=\"cornflowerblue\", plot_type=\"E\", legend=round(np.mean(E_C[\"Qedep\"]), 1))\n",
    "plot_cdf(ax1, fht_E_R, xlim=(0, 1250), n_bins=500, color=\"seagreen\", plot_type=\"E\", legend=round(np.mean(E_R[\"Qedep\"]), 1))\n",
    "ax1.set_xlabel(\"t [ns]\", fontsize=15)\n",
    "ax1.set_ylabel(\"CDF(t)\", fontsize=15)\n",
    "ax1.set_xlim((0, 1250))\n",
    "ax1.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax1.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax1.grid()\n",
    "ax1.legend(loc=(0.37, 1.01), ncol=3, frameon=False, prop={'size': 12})\n",
    "\n",
    "ax2 = fig.add_axes([0.425, 0.1, 0.55, 0.55])\n",
    "plot_hist(ax2, fht_E_R, xlim=(0, 1200), n_bins=100, color=\"seagreen\", alpha=0.9)\n",
    "plot_hist(ax2, fht_E_C, xlim=(0, 1200), n_bins=100, color=\"cornflowerblue\", alpha=0.9)\n",
    "plot_hist(ax2, fht_E_L, xlim=(0, 1200), n_bins=100, color=\"darkmagenta\", alpha=0.9)\n",
    "ax2.set_xlim((0, 1200))\n",
    "ax2.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax2.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax2.grid()\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/cdf_fht_center.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.2);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:48.750083Z",
     "iopub.status.busy": "2023-05-29T17:39:48.749557Z",
     "iopub.status.idle": "2023-05-29T17:39:48.757041Z",
     "shell.execute_reply": "2023-05-29T17:39:48.756559Z"
    }
   },
   "outputs": [],
   "source": [
    "# let's now compare the distributions for different radii at E = 1MeV\n",
    "df = targets_dataframe[(targets_dataframe[\"Qedep\"] > 0.75) & (targets_dataframe[\"Qedep\"] < 1.25)]\n",
    "R2 = df[(df[\"Redep\"] > 1.55) & (df[\"Redep\"] < 2.3)]\n",
    "R9 = df[(df[\"Redep\"] > 8.8) & (df[\"Redep\"] < 9.1)]\n",
    "R17 = df[(df[\"Redep\"] > 16.9) & (df[\"Redep\"] < 17.05)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:48.759650Z",
     "iopub.status.busy": "2023-05-29T17:39:48.759109Z",
     "iopub.status.idle": "2023-05-29T17:39:48.766326Z",
     "shell.execute_reply": "2023-05-29T17:39:48.765838Z"
    }
   },
   "outputs": [],
   "source": [
    "charge_R2, fht_R2 = charge[R2.index.values], fht[R2.index.values]\n",
    "charge_R9, fht_R9 = charge[R9.index.values], fht[R9.index.values]\n",
    "charge_R17, fht_R17 = charge[R17.index.values], fht[R17.index.values]\n",
    "\n",
    "for i in range(len(fht_R2)):\n",
    "    fht_R2[i] = fht_R2[i][fht_R2[i] != 0]\n",
    "for i in range(len(fht_R9)):\n",
    "    fht_R9[i] = fht_R9[i][fht_R9[i] != 0]\n",
    "for i in range(len(fht_R17)):\n",
    "    fht_R17[i] = fht_R17[i][fht_R17[i] != 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:48.768801Z",
     "iopub.status.busy": "2023-05-29T17:39:48.768350Z",
     "iopub.status.idle": "2023-05-29T17:39:50.541513Z",
     "shell.execute_reply": "2023-05-29T17:39:50.540964Z"
    }
   },
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(8,5), dpi=150)\n",
    "\n",
    "ax1 = fig.add_axes([0, 0, 1, 1])\n",
    "plot_cdf(ax1, charge_R2, xlim=(0, 6), n_bins=500, color=\"darkmagenta\", plot_type=\"R\", legend=round(np.mean(R2[\"Redep\"]), 1))\n",
    "plot_cdf(ax1, charge_R9, xlim=(0, 6), n_bins=500, color=\"cornflowerblue\", plot_type=\"R\", legend=round(np.mean(R9[\"Redep\"]), 1))\n",
    "plot_cdf(ax1, charge_R17, xlim=(0, 6), n_bins=500, color=\"seagreen\", plot_type=\"R\", legend=round(np.mean(R17[\"Redep\"]), 1))\n",
    "ax1.set_xlabel(\"nPE\", fontsize=15)\n",
    "ax1.set_ylabel(\"CDF(nPE)\", fontsize=15)\n",
    "ax1.set_xlim((0, 6))\n",
    "ax1.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax1.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax1.grid()\n",
    "ax1.legend(loc=(0.43, 1.01), ncol=3, frameon=False, prop={'size': 12})\n",
    "\n",
    "ax2 = fig.add_axes([0.425, 0.1, 0.55, 0.55])\n",
    "plot_hist(ax2, charge_R2, xlim=(0, 6), n_bins=100, color=\"darkmagenta\", alpha=0.9)\n",
    "plot_hist(ax2, charge_R9, xlim=(0, 6), n_bins=100, color=\"cornflowerblue\", alpha=0.9)\n",
    "plot_hist(ax2, charge_R17, xlim=(0, 6), n_bins=100, color=\"seagreen\", alpha=0.9)\n",
    "ax2.set_xlim((0, 6))\n",
    "ax2.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax2.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax2.grid()\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/cdf_charge_energy.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.2);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:50.544554Z",
     "iopub.status.busy": "2023-05-29T17:39:50.543968Z",
     "iopub.status.idle": "2023-05-29T17:39:52.066568Z",
     "shell.execute_reply": "2023-05-29T17:39:52.066050Z"
    }
   },
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(8,5), dpi=150)\n",
    "\n",
    "ax1 = fig.add_axes([0, 0, 1, 1])\n",
    "plot_cdf(ax1, fht_R2, xlim=(0, 1250), n_bins=500, color=\"darkmagenta\", plot_type=\"R\", legend=round(np.mean(R2[\"Redep\"]), 1))\n",
    "plot_cdf(ax1, fht_R9, xlim=(0, 1250), n_bins=500, color=\"cornflowerblue\", plot_type=\"R\", legend=round(np.mean(R9[\"Redep\"]), 1))\n",
    "plot_cdf(ax1, fht_R17, xlim=(0, 1250), n_bins=500, color=\"seagreen\", plot_type=\"R\", legend=round(np.mean(R17[\"Redep\"]), 1))\n",
    "ax1.set_xlabel(\"t [ns]\", fontsize=15)\n",
    "ax1.set_ylabel(\"CDF(t)\", fontsize=15)\n",
    "ax1.set_xlim((0, 1250))\n",
    "ax1.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax1.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax1.grid()\n",
    "ax1.legend(loc=(0.43, 1.01), ncol=3, frameon=False, prop={'size': 12})\n",
    "\n",
    "ax2 = fig.add_axes([0.425, 0.1, 0.55, 0.55])\n",
    "plot_hist(ax2, fht_R2, xlim=(0, 1200), n_bins=100, color=\"darkmagenta\", alpha=0.9)\n",
    "plot_hist(ax2, fht_R9, xlim=(0, 1200), n_bins=100, color=\"cornflowerblue\", alpha=0.9)\n",
    "plot_hist(ax2, fht_R17, xlim=(0, 1200), n_bins=100, color=\"seagreen\", alpha=0.9)\n",
    "ax2.set_xlim((0, 1200))\n",
    "ax2.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax2.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "ax2.grid()\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/cdf_fht_energy.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.2);"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Features Exploration"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:52.069362Z",
     "iopub.status.busy": "2023-05-29T17:39:52.068851Z",
     "iopub.status.idle": "2023-05-29T17:39:52.498682Z",
     "shell.execute_reply": "2023-05-29T17:39:52.497912Z"
    }
   },
   "outputs": [],
   "source": [
    "# 1% of the dataset is enough for basic exploration \n",
    "features_dataframe = pd.read_csv(\"/mnt/ferracci/features_dataframe_new.csv.gz\", nrows=10000)\n",
    "targets_dataframe = pd.read_csv(\"/mnt/ferracci/targets_dataframe_new.csv.gz\", nrows=10000)\n",
    "Qedep = np.array(targets_dataframe[\"Qedep\"])\n",
    "Redep = np.array(targets_dataframe[\"Redep\"])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:52.501840Z",
     "iopub.status.busy": "2023-05-29T17:39:52.501440Z",
     "iopub.status.idle": "2023-05-29T17:39:52.565547Z",
     "shell.execute_reply": "2023-05-29T17:39:52.565036Z"
    }
   },
   "outputs": [],
   "source": [
    "corr_targets_dataframe = features_dataframe.corrwith(targets_dataframe[\"Qedep\"])\n",
    "print(corr_targets_dataframe.sort_values(ascending=False)[:12])\n",
    "print(corr_targets_dataframe.sort_values(ascending=True)[:12])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:52.568388Z",
     "iopub.status.busy": "2023-05-29T17:39:52.567820Z",
     "iopub.status.idle": "2023-05-29T17:39:53.507783Z",
     "shell.execute_reply": "2023-05-29T17:39:53.507238Z"
    }
   },
   "outputs": [],
   "source": [
    "plt.style.use(\"default\")\n",
    "plt.rc(\"text\", usetex=True)\n",
    "plt.rc(\"font\", family=\"cm\")\n",
    "plt.rcParams[\"grid.color\"] = (0.5, 0.5, 0.5, 0.2)\n",
    "\n",
    "fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(7,5), dpi=150)\n",
    "\n",
    "ax = sns.heatmap(features_dataframe.corr(), vmin=-1, vmax=1, xticklabels=False, yticklabels=False)\n",
    "ax.figure.axes[-1].yaxis.label.set_size(12)\n",
    "\n",
    "fig.savefig(\"/home/ferracci/new_dataset/images/features_corr.png\", dpi=300, bbox_inches=\"tight\", pad_inches=0.1);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:53.510602Z",
     "iopub.status.busy": "2023-05-29T17:39:53.510071Z",
     "iopub.status.idle": "2023-05-29T17:39:54.451062Z",
     "shell.execute_reply": "2023-05-29T17:39:54.450525Z"
    }
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(nrows=1, ncols=2, figsize=(11,5), dpi=150)\n",
    "\n",
    "ax = axes[0]\n",
    "scatter = ax.scatter(features_dataframe[\"accum_charge\"], Qedep, s=0.3, alpha=0.9, c=Redep, cmap=\"inferno\")\n",
    "\n",
    "ax.set_xlabel(\"AccumCharge\", fontsize=15)\n",
    "ax.set_ylabel(\"$E_{kin}$ [MeV]\", fontsize=15)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "\n",
    "cbar = fig.colorbar(scatter, ax=ax, shrink=0.8, aspect=10)\n",
    "cbar.ax.set_title(\"$R$ [m]\", fontsize=15)\n",
    "cbar.outline.set_edgecolor('k')\n",
    "cbar.outline.set_linewidth(0.5)\n",
    "\n",
    "ax.grid()\n",
    "\n",
    "ax = axes[1]\n",
    "scatter = ax.scatter(features_dataframe[\"nPMTs\"], Qedep, s=0.3, alpha=0.9, c=Redep, cmap=\"inferno\")\n",
    "\n",
    "ax.set_xlabel(\"nPMTs\", fontsize=15)\n",
    "ax.set_ylabel(\"$E_{kin}$ [MeV]\", fontsize=15)\n",
    "\n",
    "ax.tick_params(axis=\"both\", which=\"major\", labelsize=12)\n",
    "ax.tick_params(axis=\"both\", which='minor', labelsize=12)\n",
    "\n",
    "cbar = fig.colorbar(scatter, ax=ax, shrink=0.8, aspect=10)\n",
    "cbar.ax.set_title(\"$R$ [m]\", fontsize=15)\n",
    "cbar.outline.set_edgecolor('k')\n",
    "cbar.outline.set_linewidth(0.5)\n",
    "\n",
    "ax.grid();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:54.455300Z",
     "iopub.status.busy": "2023-05-29T17:39:54.454857Z",
     "iopub.status.idle": "2023-05-29T17:39:54.458221Z",
     "shell.execute_reply": "2023-05-29T17:39:54.457751Z"
    }
   },
   "outputs": [],
   "source": [
    "def plot_feature_distribution(ax, feature, feature_name, bins):\n",
    "    ax.hist(feature, bins=bins, color=\"steelblue\", density=True)\n",
    "    ax.set_xlabel(feature_name, fontsize=15)\n",
    "    ax.grid()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:54.460647Z",
     "iopub.status.busy": "2023-05-29T17:39:54.460176Z",
     "iopub.status.idle": "2023-05-29T17:39:56.358465Z",
     "shell.execute_reply": "2023-05-29T17:39:56.357929Z"
    }
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(nrows=2, ncols=4, figsize=(15,8), dpi=100)\n",
    "\n",
    "ax = axes[0, 0]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_std\"], \"$pe_{std}$\", 50)\n",
    "\n",
    "ax = axes[0, 1]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_skew\"], \"$pe_{skewness}$\", 50)\n",
    "\n",
    "ax = axes[0, 2]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_kurtosis\"], \"$pe_{kurtosis}$\", 50)\n",
    "\n",
    "ax = axes[0, 3]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_entropy\"], \"$pe_{entropy}$\", 50)\n",
    "\n",
    "ax = axes[1, 0]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_cht_std\"], \"$ht_{std}$\", 50)\n",
    "\n",
    "ax = axes[1, 1]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_cht_skew\"], \"$ht_{skewness}$\", 50)\n",
    "\n",
    "ax = axes[1, 2]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_cht_kurtosis\"], \"$ht_{kurtosis}$\", 50)\n",
    "\n",
    "ax = axes[1, 3]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe_cht_entropy\"], \"$ht_{entropy}$\", 50);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:56.361251Z",
     "iopub.status.busy": "2023-05-29T17:39:56.360749Z",
     "iopub.status.idle": "2023-05-29T17:39:59.355059Z",
     "shell.execute_reply": "2023-05-29T17:39:59.354526Z"
    }
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(nrows=3, ncols=4, figsize=(15,12), dpi=100)\n",
    "\n",
    "ax = axes[0, 0]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe10_cc\"], \"$pe_{10\\%, cc}$\", 50)\n",
    "\n",
    "ax = axes[0, 1]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe40_cc\"], \"$pe_{40\\%, cc}$\", 50)\n",
    "\n",
    "ax = axes[0, 2]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe70_cc\"], \"$pe_{70\\%, cc}$\", 50)\n",
    "\n",
    "ax = axes[0, 3]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe95_cc\"], \"$pe_{95\\%, cc}$\", 50)\n",
    "\n",
    "ax = axes[1, 0]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe10_cht\"], \"$pe_{10\\%, cht}$\", 50)\n",
    "\n",
    "ax = axes[1, 1]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe40_cht\"], \"$pe_{40\\%, cht}$\", 50)\n",
    "\n",
    "ax = axes[1, 2]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe70_cht\"], \"$pe_{70\\%, cht}$\", 50)\n",
    "\n",
    "ax = axes[1, 3]\n",
    "plot_feature_distribution(ax, features_dataframe[\"pe95_cht\"], \"$pe_{95\\%, cht}$\", 50)\n",
    "\n",
    "ax = axes[2, 0]\n",
    "plot_feature_distribution(ax, features_dataframe[\"ht10_5\"], \"$ht_{10\\%-5\\%}$\", 50)\n",
    "\n",
    "ax = axes[2, 1]\n",
    "plot_feature_distribution(ax, features_dataframe[\"ht40_35\"], \"$ht_{40\\%-35\\%}$\", 50)\n",
    "\n",
    "ax = axes[2, 2]\n",
    "plot_feature_distribution(ax, features_dataframe[\"ht70_65\"], \"$ht_{70\\%-65\\%}$\", 50)\n",
    "\n",
    "ax = axes[2, 3]\n",
    "plot_feature_distribution(ax, features_dataframe[\"ht95_90\"], \"$ht_{95\\%-90\\%}$\", 50);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-05-29T17:39:59.357791Z",
     "iopub.status.busy": "2023-05-29T17:39:59.357277Z",
     "iopub.status.idle": "2023-05-29T17:40:01.044207Z",
     "shell.execute_reply": "2023-05-29T17:40:01.043666Z"
    }
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(nrows=2, ncols=3, figsize=(12,7), dpi=150)\n",
    "\n",
    "ax = axes[0, 0]\n",
    "plot_feature_distribution(ax, features_dataframe[\"gamma_x_cc\"], \"$\\gamma_{x,cc}$\", 100)\n",
    "\n",
    "ax = axes[0, 1]\n",
    "plot_feature_distribution(ax, features_dataframe[\"gamma_y_cc\"], \"$\\gamma_{y,cc}$\", 100)\n",
    "\n",
    "ax = axes[0, 2]\n",
    "plot_feature_distribution(ax, features_dataframe[\"gamma_z_cc\"], \"$\\gamma_{z,cc}$\", 100)\n",
    "\n",
    "ax = axes[1, 0]\n",
    "plot_feature_distribution(ax, features_dataframe[\"gamma_x_cht\"], \"$\\gamma_{x,cht}$\", 100)\n",
    "\n",
    "ax = axes[1, 1]\n",
    "plot_feature_distribution(ax, features_dataframe[\"gamma_y_cht\"], \"$\\gamma_{y,cht}$\", 100)\n",
    "\n",
    "ax = axes[1, 2]\n",
    "plot_feature_distribution(ax, features_dataframe[\"gamma_z_cht\"], \"$\\gamma_{z,cht}$\", 100);"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "deep_learning",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.6"
  },
  "widgets": {
   "application/vnd.jupyter.widget-state+json": {
    "state": {
     "1f0781ec49504a81adbec182c3952aec": {
      "model_module": "@jupyter-widgets/controls",
      "model_module_version": "2.0.0",
      "model_name": "HTMLModel",
      "state": {
       "_dom_classes": [],
       "_model_module": "@jupyter-widgets/controls",
       "_model_module_version": "2.0.0",
       "_model_name": "HTMLModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/controls",
       "_view_module_version": "2.0.0",
       "_view_name": "HTMLView",
       "description": "",
       "description_allow_html": false,
       "layout": "IPY_MODEL_4d2013fb39f64c77a1d97c1525ecc678",
       "placeholder": "​",
       "style": "IPY_MODEL_dbc44638940846c0b6ed6ef26b4e04e2",
       "tabbable": null,
       "tooltip": null,
       "value": "100%"
      }
     },
     "41de841ebfb248d38d8f50e41fbffc6b": {
      "model_module": "@jupyter-widgets/controls",
      "model_module_version": "2.0.0",
      "model_name": "ProgressStyleModel",
      "state": {
       "_model_module": "@jupyter-widgets/controls",
       "_model_module_version": "2.0.0",
       "_model_name": "ProgressStyleModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/base",
       "_view_module_version": "2.0.0",
       "_view_name": "StyleView",
       "bar_color": null,
       "description_width": ""
      }
     },
     "4d2013fb39f64c77a1d97c1525ecc678": {
      "model_module": "@jupyter-widgets/base",
      "model_module_version": "2.0.0",
      "model_name": "LayoutModel",
      "state": {
       "_model_module": "@jupyter-widgets/base",
       "_model_module_version": "2.0.0",
       "_model_name": "LayoutModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/base",
       "_view_module_version": "2.0.0",
       "_view_name": "LayoutView",
       "align_content": null,
       "align_items": null,
       "align_self": null,
       "border_bottom": null,
       "border_left": null,
       "border_right": null,
       "border_top": null,
       "bottom": null,
       "display": null,
       "flex": null,
       "flex_flow": null,
       "grid_area": null,
       "grid_auto_columns": null,
       "grid_auto_flow": null,
       "grid_auto_rows": null,
       "grid_column": null,
       "grid_gap": null,
       "grid_row": null,
       "grid_template_areas": null,
       "grid_template_columns": null,
       "grid_template_rows": null,
       "height": null,
       "justify_content": null,
       "justify_items": null,
       "left": null,
       "margin": null,
       "max_height": null,
       "max_width": null,
       "min_height": null,
       "min_width": null,
       "object_fit": null,
       "object_position": null,
       "order": null,
       "overflow": null,
       "padding": null,
       "right": null,
       "top": null,
       "visibility": null,
       "width": null
      }
     },
     "6e90325ca1ef4289a5bc545304dad3ad": {
      "model_module": "@jupyter-widgets/base",
      "model_module_version": "2.0.0",
      "model_name": "LayoutModel",
      "state": {
       "_model_module": "@jupyter-widgets/base",
       "_model_module_version": "2.0.0",
       "_model_name": "LayoutModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/base",
       "_view_module_version": "2.0.0",
       "_view_name": "LayoutView",
       "align_content": null,
       "align_items": null,
       "align_self": null,
       "border_bottom": null,
       "border_left": null,
       "border_right": null,
       "border_top": null,
       "bottom": null,
       "display": null,
       "flex": null,
       "flex_flow": null,
       "grid_area": null,
       "grid_auto_columns": null,
       "grid_auto_flow": null,
       "grid_auto_rows": null,
       "grid_column": null,
       "grid_gap": null,
       "grid_row": null,
       "grid_template_areas": null,
       "grid_template_columns": null,
       "grid_template_rows": null,
       "height": null,
       "justify_content": null,
       "justify_items": null,
       "left": null,
       "margin": null,
       "max_height": null,
       "max_width": null,
       "min_height": null,
       "min_width": null,
       "object_fit": null,
       "object_position": null,
       "order": null,
       "overflow": null,
       "padding": null,
       "right": null,
       "top": null,
       "visibility": null,
       "width": null
      }
     },
     "72e18cc4c0124f3c8c02806a2cca2cb5": {
      "model_module": "@jupyter-widgets/controls",
      "model_module_version": "2.0.0",
      "model_name": "HBoxModel",
      "state": {
       "_dom_classes": [],
       "_model_module": "@jupyter-widgets/controls",
       "_model_module_version": "2.0.0",
       "_model_name": "HBoxModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/controls",
       "_view_module_version": "2.0.0",
       "_view_name": "HBoxView",
       "box_style": "",
       "children": [
        "IPY_MODEL_1f0781ec49504a81adbec182c3952aec",
        "IPY_MODEL_88769189c1ff4d76ad43014c1644263f",
        "IPY_MODEL_d0953608d0f141239d4c65456d47dbd5"
       ],
       "layout": "IPY_MODEL_e784cd0ff4854c4d86c2cce64ce7b9b3",
       "tabbable": null,
       "tooltip": null
      }
     },
     "88769189c1ff4d76ad43014c1644263f": {
      "model_module": "@jupyter-widgets/controls",
      "model_module_version": "2.0.0",
      "model_name": "FloatProgressModel",
      "state": {
       "_dom_classes": [],
       "_model_module": "@jupyter-widgets/controls",
       "_model_module_version": "2.0.0",
       "_model_name": "FloatProgressModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/controls",
       "_view_module_version": "2.0.0",
       "_view_name": "ProgressView",
       "bar_style": "success",
       "description": "",
       "description_allow_html": false,
       "layout": "IPY_MODEL_6e90325ca1ef4289a5bc545304dad3ad",
       "max": 125,
       "min": 0,
       "orientation": "horizontal",
       "style": "IPY_MODEL_41de841ebfb248d38d8f50e41fbffc6b",
       "tabbable": null,
       "tooltip": null,
       "value": 125
      }
     },
     "bad6749bc685412190bae30bf886a037": {
      "model_module": "@jupyter-widgets/base",
      "model_module_version": "2.0.0",
      "model_name": "LayoutModel",
      "state": {
       "_model_module": "@jupyter-widgets/base",
       "_model_module_version": "2.0.0",
       "_model_name": "LayoutModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/base",
       "_view_module_version": "2.0.0",
       "_view_name": "LayoutView",
       "align_content": null,
       "align_items": null,
       "align_self": null,
       "border_bottom": null,
       "border_left": null,
       "border_right": null,
       "border_top": null,
       "bottom": null,
       "display": null,
       "flex": null,
       "flex_flow": null,
       "grid_area": null,
       "grid_auto_columns": null,
       "grid_auto_flow": null,
       "grid_auto_rows": null,
       "grid_column": null,
       "grid_gap": null,
       "grid_row": null,
       "grid_template_areas": null,
       "grid_template_columns": null,
       "grid_template_rows": null,
       "height": null,
       "justify_content": null,
       "justify_items": null,
       "left": null,
       "margin": null,
       "max_height": null,
       "max_width": null,
       "min_height": null,
       "min_width": null,
       "object_fit": null,
       "object_position": null,
       "order": null,
       "overflow": null,
       "padding": null,
       "right": null,
       "top": null,
       "visibility": null,
       "width": null
      }
     },
     "d0953608d0f141239d4c65456d47dbd5": {
      "model_module": "@jupyter-widgets/controls",
      "model_module_version": "2.0.0",
      "model_name": "HTMLModel",
      "state": {
       "_dom_classes": [],
       "_model_module": "@jupyter-widgets/controls",
       "_model_module_version": "2.0.0",
       "_model_name": "HTMLModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/controls",
       "_view_module_version": "2.0.0",
       "_view_name": "HTMLView",
       "description": "",
       "description_allow_html": false,
       "layout": "IPY_MODEL_bad6749bc685412190bae30bf886a037",
       "placeholder": "​",
       "style": "IPY_MODEL_d169b66a2e3444e5ad4dc3c938131543",
       "tabbable": null,
       "tooltip": null,
       "value": " 125/125 [16:12&lt;00:00,  7.66s/it]"
      }
     },
     "d169b66a2e3444e5ad4dc3c938131543": {
      "model_module": "@jupyter-widgets/controls",
      "model_module_version": "2.0.0",
      "model_name": "HTMLStyleModel",
      "state": {
       "_model_module": "@jupyter-widgets/controls",
       "_model_module_version": "2.0.0",
       "_model_name": "HTMLStyleModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/base",
       "_view_module_version": "2.0.0",
       "_view_name": "StyleView",
       "background": null,
       "description_width": "",
       "font_size": null,
       "text_color": null
      }
     },
     "dbc44638940846c0b6ed6ef26b4e04e2": {
      "model_module": "@jupyter-widgets/controls",
      "model_module_version": "2.0.0",
      "model_name": "HTMLStyleModel",
      "state": {
       "_model_module": "@jupyter-widgets/controls",
       "_model_module_version": "2.0.0",
       "_model_name": "HTMLStyleModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/base",
       "_view_module_version": "2.0.0",
       "_view_name": "StyleView",
       "background": null,
       "description_width": "",
       "font_size": null,
       "text_color": null
      }
     },
     "e784cd0ff4854c4d86c2cce64ce7b9b3": {
      "model_module": "@jupyter-widgets/base",
      "model_module_version": "2.0.0",
      "model_name": "LayoutModel",
      "state": {
       "_model_module": "@jupyter-widgets/base",
       "_model_module_version": "2.0.0",
       "_model_name": "LayoutModel",
       "_view_count": null,
       "_view_module": "@jupyter-widgets/base",
       "_view_module_version": "2.0.0",
       "_view_name": "LayoutView",
       "align_content": null,
       "align_items": null,
       "align_self": null,
       "border_bottom": null,
       "border_left": null,
       "border_right": null,
       "border_top": null,
       "bottom": null,
       "display": null,
       "flex": null,
       "flex_flow": null,
       "grid_area": null,
       "grid_auto_columns": null,
       "grid_auto_flow": null,
       "grid_auto_rows": null,
       "grid_column": null,
       "grid_gap": null,
       "grid_row": null,
       "grid_template_areas": null,
       "grid_template_columns": null,
       "grid_template_rows": null,
       "height": null,
       "justify_content": null,
       "justify_items": null,
       "left": null,
       "margin": null,
       "max_height": null,
       "max_width": null,
       "min_height": null,
       "min_width": null,
       "object_fit": null,
       "object_position": null,
       "order": null,
       "overflow": null,
       "padding": null,
       "right": null,
       "top": null,
       "visibility": null,
       "width": null
      }
     }
    },
    "version_major": 2,
    "version_minor": 0
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}