|
|
{ "cells": [ { "cell_type": "code", "execution_count": 43, "metadata": {}, "outputs": [], "source": [ "import uproot\t\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from mpl_toolkits import mplot3d\n", "import awkward as ak\n", "from scipy.optimize import curve_fit\n", "from methods.fit_linear_regression_model import fit_linear_regression_model\n", "import sklearn\n", "%matplotlib inline" ] }, { "cell_type": "code", "execution_count": 44, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "40402 10099\n", "50501\n" ] } ], "source": [ "file = uproot.open(\n", " \"tracking_losses_ntuple_B_default_radlength_endVelo.root:PrDebugTrackingLosses.PrDebugTrackingTool/Tuple;1\"\n", ")\n", "\n", "# selektiere nur elektronen von B->K*ee\n", "allcolumns = file.arrays()\n", "found = allcolumns[(allcolumns.isElectron) & (~allcolumns.lost) &\n", " (allcolumns.fromB)] # B: 9056\n", "lost = allcolumns[(allcolumns.isElectron) & (allcolumns.lost) &\n", " (allcolumns.fromB)] # B: 1466\n", "\n", "electrons = allcolumns[(allcolumns.isElectron) & (allcolumns.fromB)]\n", "\n", "print(ak.num(found, axis=0), ak.num(lost, axis=0))\n", "print(ak.num(electrons, axis=0))\n", "# ak.count(found, axis=None)" ] }, { "cell_type": "code", "execution_count": 45, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "stretch factor: 0.24996287312509283\n" ] } ], "source": [ "rad_length_found = ak.to_numpy(found[\"rad_length_frac\"])\n", "eta_found = ak.to_numpy(found[\"eta\"])\n", "rad_length_lost = ak.to_numpy(lost[\"rad_length_frac\"])\n", "eta_lost = ak.to_numpy(lost[\"eta\"])\n", "\n", "stretch_factor = ak.num(eta_lost, axis=0) / ak.num(eta_found, axis=0)\n", "print(\"stretch factor: \", stretch_factor)" ] }, { "cell_type": "code", "execution_count": 46, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAj8AAAHLCAYAAAAnR/mlAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAyUUlEQVR4nO3dUYgb+4Hv+Z8nuafjTONT3c4+ZOFCXIJ5GC5cU+pc2HCJD1gisHnZpSUb9iGZ5dAqQr95GQnPi+2XEdXMa0NKZpizeVg4rTK7CzsPM6ozHIdlYNdWjQPLstzdLgfykAub010xfeP0mUx0H5wqt1rqbnVJJam6vh9oznGVqv5/ValVv/7//1X/a4PBYCAAAICC+KNFVwAAAGCeCD8AAKBQCD8AAKBQCD8AAKBQCD8AAKBQCD8AAKBQCD8AAKBQCD8AAKBQCD8otCAI1Ol0Fl2NwouiSJ1OR2EYLroquCR+h5BHhB8UUhiGqtfrKpfLcl13aPna2tpMv8yz2OekgiBQq9VSuVxWuVyee/mT6HQ6unXrlmzbVhRFU+0rfr/1el2lUkk7OzuzqWRKizz3WRv3O+T7vkqlkq5du6Zr166pXC7L87yRbT3PU7lcTl4TBMG5ZeXhc4x8IfygkEzTVLfbHVkeRZGiKNL+/n6q/YZhOHIBn3af07AsS/fv37/w4jJPp49Ro9FQo9GYer9BEOju3btyHEfdble1Wm2ux3zZzn3Wxv0OVSoV9fv95N+2batWq41sW6vVZNu2DMNQv9+XZVnnlrWMn2PkG+EHOMGyLA0GAzmOk2r7er2ug4ODme5zWhddWOZt3DG6efPm1Pttt9taX19P/u04zlCrXtaW8dwvgmEYajabkjT2D4xYv9/Xw4cPJ97vsn2OkW+EH2BG6vU6f5leIMtjtMhjz7kfFoca3/fPPC57e3szafED0iD8IBc8z1O1WpXv++p0OlpbW5Nt28n6nZ0d2badjAsYN8YiiiLZtp38jBsP4nme6vW66vX6yLrzyvA8L/mSt2176GJ43j7jOrVaLVWr1eQ9nlzf6XSSsRO+7ydjJcbt77KCIFC9Xle1WlWpVFKr1Zqq7PjYlkqlkUBw3jGKHRwcJGVO+h47nY7q9brCMEzGodTr9eQ4zvqzE5+nk+d3mc/9LM/xJL9D0rvWn7i7q91uj6z3PE+VSkWGYUxUz4tcdCyBEQNgyXW73YFpmgNJg0ajMWg2mwPLsgaWZQ0Gg8Gg2WwOTn6Ue73eQNKg1+sly/b39weGYQwtcxxnICnZz/7+frKsUqkM1WGSMuLX7O/vD5V71j77/f7AMIxBv99PlrmuO5A0cBwn2b5WqyXbN5vNQb/fHzQajaHXXeTk+zxZ/sk6dbvd5BinKduyrEGz2Uy2lTSQNDBNM9nnuGM0GAyGjpHjOKneo2maA9M0h5bN6rNjmuZQnQ3DGBiGkfx7Wc/9LM/xJL9DJ538DBweHg6tsyxraD8X1TN21uf4omMJnEb4QS7EX7LxxfWkSqUydCE6PDwceW2tVhu5AA0G479Mx12sJinjrAv7Wfu0LGtsnSzLGtpPfCE4/UU+bp9nGfc+LcsaumAMBu8u6icvVpOWHb/u5P7iC+jJ43FR+HFdN/V7HBd+Tu477WfHsqyR9x/vMz5Oy3ruZ3mOL/M7FKtUKiPHMw5Rl63nWWVNeiyBk+j2Qi7EzePf/va3R9Z1u92hO0xevnwpScmdN2EYJl0faV1UxmWFYaggCMYO4oy7ZE4P1j3ZRRA7PcD2suW32+2kW+ZkN0f8/iYt+8WLFyPr4/1dZizMxsbGSLnTPvtn2s9OEASqVCpD2zWbTQ0Gg7HH5SLzOvezPMdpf4firquTXYmO4wwNdL5sPce9x8scS0CSvrroCgCXMe4L2jAMGYYhz/P06aefjlzk4ounaZpTlXteGZd1XiCIA0CWD/yLyz/vbpzLiI+H7/sjF6JluUsnzWcnPk5pQs5Z5nXuZ3mO0/4OVSoVmaapMAzV6XTUaDS0t7en169fz6Sei/49Qn7R8oPcC8NQ5XJZYRiq2+0mt9meXC+lbyWZpIy0xrUcxRfak7dtz1p8TGZ1YajVaqpUKmq32/J9X1EUyXEcNZvNqUJn1ib97GRxAc363M+y7tP8DsWtP47jyPM8bWxsDIXJWdRzUb9HyC/CD3KvWq1qfX39zEASX3xPdm/MuozLiltDxt2REn+Rl0qlmZQ1TnxMxj1996x6XaTb7apSqSTTHTiOs/TPt7novMbn6axWiTQX7Hmd+1me42l+hxqNRtJ9ubW1NXSn3bT1XPTvEfKL8INci29vPvmXZPylF/+VGjd/dzqdsX8hXjRuZ5IyLrM/6d0XvmVZyb5PevnypQzDyPQZKPEYllarNdJ1kHYqhnq9nrSeNJvNc7u7Th+jL774QtJ0rXMHBweX2v6yn53TF9hWqzXSqrBM536W53ja36GTY3xOP/F5mnou+vcI+UX4QS7EX66nv2Tji4/neep0Oup0OkkzexAEyV+T8V/25XJZvu8rDMPkdWEYJs8rGRdqJikjiqLkL0zXdZMBomftU3rXmmAYxtBfwnF30dOnT5OL8jSB4Cwnn8JbLpdVr9e1s7OjarWq/f395II0adlxOIj30+l0hp5/EzvrGJ0liqKp3/+0n5249aparaperyfPAyqVSsk5WsZzP8tzfHpf5/0OjRMHkHFBZNJ6nmXSYwkMWfTtZsBFTj6rxTTNkduhXdcdGIYxME0zuaW20WgMDMMYusXWdd1kP5ZlJc9vaTabg/39/UG/30+eeaI/3Pob32Y7aRmWZQ0Mw0ieT3LePgeDd7dWx7cQNxqNQaPRGLrlt9/vJ7fsmqY56PV6g8PDw+Q2co25Rfmkk89sGVe+4zhnHtvLlB0fy3j5yR/TNIfKPH2Mut1ucltzfMvz6XLG3aZ+1ntsNBrJM2Rm9dnpdrvJsTj9jJqz3teiz31sVuc4Pl7n/Q6dp9FonPuai+p53uf4omMJnHZtMBgMsotWAIogCAJ9+umnevjwoQ4ODoZabLrdrkql0szGSwHAtLjVHcBU4jumDg8Pk1vHTzJNk6kGACwVxvwAmEo80HRra2tojE/8bBfXdRl0CmCp0O0FYGo7Oztqt9tDg4oty5LjOBcOWAWAeSP8AJiZeKzPMj/YEAAIPwAAoFAY8wMAAAqlEHd7/epXv9Lf/d3f6Vvf+pauX7++6OoAAIAJvH37Vj//+c/1ve99T9/4xjdmtt9ChJ+//du/1Z/92Z8tuhoAACCFTz75RD/84Q9ntr9ChJ9vfvObkqSnT5+eO9/QNDY3N/Xs2bNM9n2Vyjg6OtKdO3f0/Plzra6uZlZO1u+DczG5q3CsOBfFKYNzsVxlBEGgra2t5Do+K5mFnx/96Efa2NjQxx9/nFURE/va174mSfqTP/mTzMLP9evXM9v3VSrjzZs3kqTbt2/rxo0bmZWT9fvgXEzuKhwrzkVxyuBcLFcZR0dHkt5fx2clkwHPr1+/luu6PM4eAAAsnUxafm7duiXXdXnWBwAAWDqZdXttbW1ltWsAAIDUMgs/3/ve9+Q4jm7fvp1VEZf2gx/8QF//+tfHrtve3tb29vacawQAQLHt7u5qd3d37Lrf/OY3mZSZOvzcv3//zHVRFMn3fe3t7S1V+PnJT36i7373u4uuBgAA+IPzGh9++tOf6s6dOzMvM3X46Xa7E73mL//yL9MWkSvzaDW6KmXMQ9bvg3MxuatwrDgXxSsja1flOOX1XKSe2+vevXtyHEfr6+sj6/b399XpdPTjH/946grOQpwcnz9/TsvPgr1580Yffvihfv3rX2d6GykuxrlYHpyL5cG5WC5ZXb9T3+pu27Zu3bqlDz/8cOTHsiyVy2X9xV/8xcwqCgAAMAupw8/du3fPXW+aplzXTbt7AACATKQe8/Pq1asz14VhqFarlXbXAAAAmUkdfizL0rVr185cPxgMtLOzk3b3AAAAmUgdfgzD0L1792QYxsi6mzdvyrKsC7vGAAAA5i11+Hn69Kk2NzdnWRcAAIDMpR7wTPABAAB5lNn0Fj/60Y+0sbGhjz/ "text/plain": [ "<Figure size 640x480 with 1 Axes>" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.hist(\n", " rad_length_lost,\n", " bins=100,\n", " density=True,\n", " alpha=0.5,\n", " color=\"darkorange\",\n", " histtype=\"bar\",\n", " label=\"lost\",\n", " range=[0, 1],\n", ")\n", "plt.hist(\n", " rad_length_found,\n", " bins=100,\n", " density=True,\n", " alpha=0.5,\n", " color=\"blue\",\n", " histtype=\"bar\",\n", " label=\"found\",\n", " range=[0, 1],\n", ")\n", "plt.xlim(0, 1)\n", "# plt.yscale(\"log\")\n", "plt.title(\"radiation length fraction endVelo\")\n", "plt.xlabel(f\"$x/X_0$\")\n", "plt.ylabel(\"a.u.\")\n", "\n", "plt.legend()\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 47, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABiMAAAL7CAYAAACfnvZxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAAC5AUlEQVR4nOz9f5hb1X3o+3+2h+BOsI08pN9CnPZgTXsvpvnBlUwI+BxcQLopTdIfsWRyn5O0+Z6ClHbSQk0zKt+ee4bp87QTTYoLaaeNZN+Wc5KeE1tye9vclOZIQMx5jEsyUmlIMKftyGmLA+2lM/KPMDFhvL9/TCRmZq+1RtqjvbUlvV/Powe8PuvX/qVZey/tvS3btm0BAAAAAAAAAADwyKZudwAAAAAAAAAAAPQ3JiMAAAAAAAAAAICnmIwAAAAAAAAAAACeYjICAAAAAAAAAAB4iskIAAAAAAAAAADgKSYjAAAAAAAAAACAp5iMAAAAAAAAAAAAnmIyAgAAAAAAAAAAeIrJCAAAAAAAAAAA4CkmIwAAAAKmWq3K9PS0RKNRicfj3e5Oz6jVapLP5/tynZXLZUmn07J9+3YpFovd7g56FN8t3gjieu3n70MAANC7mIwAAAAIkGq1KkeOHJFMJiPValWZp1gsyvbt22V6etqXPvndXruq1arE43EZHR2VdDot5XK5213qqGKxKLlcTvL5vNTrdU/amJ6elmQyKdu3b5ft27dLOp32rK1OC/r+GRStfLegfab1WiwWJRqNimVZqz7xeFy7v1arVUkmk6vyJ5NJqdVqLfenn78PAQBAb7Ns27a73QkAAACslkwmpVgsSiwWk1KptCoWjUalWq1KOByWubm5jrTXuPAcCoUcMS/a80I+n5d0Oi0iIr06xDVtB8uyRESkUChIIpHoWHt33HGHxGIxyWazUq/XJRqNSq1Wk2w2K+Pj4x1pZ6P6Yf8MCtN3S6eZtlu/Ma3XeDzenBSIRCJSqVTWrW96eloymYzr47Afvg8BAED/4c4IAACAABoZGdHGstls8+Jxp5h+eetFe14Ih8Pd7sKGtfML6E5o/Jq7cdEyFArJ3NyclEqlwExEiPTH/hkUpu+WTvN7f+4m03otFArN/2/1rpR//dd/lVAo5Po47IfvQwAA0H8u63YHAAAA0J5YLCaxWKxj9WUyGeOjPDrdHtTW2w5eyOfzIuK8kBqk7c3+2Zu6sT8HVSgUklQq1Tzepqen151kyOfzkkql/OgeAACAb7gzAgAAYIAVi0WetR8A3dgOvfCLdfbP3sR2c1p5587U1JQxb7lclnq9Lg888IDX3QIAAPAVkxEAAAAtKJfLkkwmm79sTafTzReRrlQsFiUej0symZR4PC7RaLRZRieTycjo6KhEo9GW8tdqNclkMhKNRrV5WulHsViUe+65p/nvO+64o/kC43bbE1n+tW+jrdHRUYnH49plqdfrks/nV+XJ5/MyOjoqlmU13wPQafV6XdLp9Kp2VBdNO9G/xvpobNt0Oi3FYtHxYuhWt8NamUxGtm/f3vb6yufzsn379lXbc+fOnc02V64PL/f7letn+/btEo/HVy1Dv+6fra6rTvSh3e+W9bRy/LSzP7vZbzbSt0a+bqzXUCjUfNdLvV6XYrGozZvNZiWRSDjetdHqMrajneMCAABgw2wAAABoZbNZOxwO2yJii4idy+XsWCxmh0KhZtrc3Jxt27Ydi8VsEbEXFhaa5cfHx20RsbPZrKPuhYUFOxKJ2LFYbFWZRCLRrDsWizXTK5VKsw0RsUOhkLLP7fRjYWGhWV+lUlkVa7W9SqVih8NhO5FIrKq3sRxrl69UKtmRSKRZbzabtROJhB2JROxUKrVqfa8st55SqdQsp4tHIhG7VCrZtm3bc3NzzeWLRCId7V8kEllVZ6FQaJZZ+Vm5vnTboWHlPtjYb8bHx12vr7m5Occ+3ODlfj83N2eHw2E7lUo1y6zsS6FQaGm99Nr+2c662mgf2v1uaUWrx0+j/fX253b3m070rdvrdeV+Hg6HjXkay9LuMq7Mb/o+bPe4AAAA6AQmIwAAANax8sJaOBxuXgzKZrN2KpWybfuNCz9rL4g2LiypLpRGIhHtBdTGRV/Vha1UKqWts91+tHLR0NReo6/rTYyolqOxjKFQaNWFt4WFhWasnYuSpotvjeVfu4wrl398fLwj/WtcUF15Ud2237hgGQ6H7YWFhVUX+tqZjAiHw448jf6sXQYT02TE2j51er9feQG0obGvrLxI24/7Zzvrym0f3H636LR7/Ky33dysi071zba7u15XTqCtnXCwbbs5wbjRZVxvMsLtcQEAALARPKYJAABgHaFQqPm4jJUvyx0fH5dcLrcqbzgcXvXvxouB6/X6qsfz5PN5qVar2meCm17Iu97jaNrpRytM7WUyGanX69oXrTaek14ulx2PJWn0MZVKrVreUCgku3fvFhGRubm5tvpq6mcoFJJIJLIqfWXa2keTuO1fYznXboN0Oi0ib7yrYe0jWFqVzWYdy9HoXyffA+HFfp/JZLT7fTablVAo1PbLqHtt/2zn2HTTh418t+i4OX5a0YnvKb+O7U6t10wm0/z/le+RWNlO47tiZZlOrv+NHBcAAAAbcVm3OwAAANBLdBc+Y7GYzM3NOS6urbw4PD8/37y427jYs/biUkPjoly7sXb70QpTe40LVTfeeKMyHolEJBwOS61WkyNHjjSfmb6eRv87dXG90c/R0VFHrNFG4wJoK+vG1L9G2vz8/Kr0lRcr290GG+lPJ3Rqv29cMFXt95FIRBYWFtruW6/sn50+NnV92Mh3i06nj59Orgu/ju1OrddYLCaRSESq1aqUy2WpVqvNOhvvlFk7SdDpZfTquAAAAFgPkxEAAAAdsvLCWrFYlCNHjjgutjU0LiDp4n71Y6NauRjbuKjlxQupW9FoNxaLSalU8rVd3a+lvdoe3dDq/lar1aRer/u67EHbP/04Njv93eLV8dOJdeHnsd3J9frAAw9IMpkUEZGpqSkpFArN/1/74movljFoxwUAABgcPKYJAACgg/L5vIyOjsr8/LwUCgXlIz28+vV6u/3YqLW/YtZpXLwz5fFSo10/1ruINH/VvPZRRo3HzvTjr4zb2e/92g+Cun96eWx6sY97efxsdF34dWx3uv5EItHc74rFYvOiv+oxUJ1exqAeFwAAYDBwZwQAAECHxONxKZfLysePrLTyWei1Wq3jv4xutR8btbJu07PzG7/y7dbdAI3Hp/g1GZHL5aRWq0m5XJZkMinZbFZGRkYkmUxKJBJp/gq6X7S6v619H0AnH1OlEsT90+tj04vvFq+On06sC7+ObS/WayaTab4bovEIqHA47HgMVKeXMYjHBQAAGBzcGQEAANABmUxGyuWypFKpdS/erLzY1OlHYLTTj05otFEul7V5GhfRGi+G9dvax8Ho1Ot1Vy/hVSmVShKJRGRkZEQymUzzwmOlUulI/UHhdr83redardbxC69B2D/9ODa9+G7x4vjp1Lrw69j2Yr2mUqlV71LJ5/OrXm7d4MUyBum4AAAAg4XJCAAAgA4wvWBUpXExaGpqqqv9WMnN4zgav+ytVqvaC8iNi3eqC21+CIVCzfV9zz33rPqV80rJZFL7jod2JZNJKRQKksvlmv9t9fFMvfRYlHb3t8ZF3ampKe3+kk6nlReoe33/3Mix2Y5Of7ds9PhRbbdOrQs/j20vvrPXPpJp7YurRbxZxiAdFwAAYLAwGQEAANABjQtuR44cWZW+8sLVyotIjcdy1Ov15otMV+Y7evSoiKgfzdFoS3VRqt1+qF6UurZdU3vj4+PNC2WNZVqpcbFrZb5WdPqC/Mr1HY1GVy1rtVqVaDQqkUik5T6a+pdOp9t+lEsr26HTdBc127GR/T4ajUo+n2/Gy+WyjI6Orrr42U/7Z7vrym0fNvLdotPu8bPeduvkuvDr2PZivY6Pjzf/XzURoWq7E8vo1XEBAACwLhsAAADrEhFbROxUKqWMp1KpZp5wOGwnEgk7HA7bhULBkT43N6csMz4+bo+Pj9vhcNiORCLNWCgUskulUrO "text/plain": [ "<Figure size 2000x800 with 3 Axes>" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "nbins = 100\n", "vmax = 80\n", "\n", "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(20, 8))\n", "\n", "a0 = ax0.hist2d(\n", " rad_length_found,\n", " eta_found,\n", " density=False,\n", " bins=nbins,\n", " cmap=plt.cm.jet,\n", " cmin=1,\n", " vmax=vmax,\n", " range=[[0, 0.6], [2, 5]],\n", ")\n", "ax0.set_xlabel(f\"$x/X_0$\")\n", "ax0.set_ylabel(f\"$\\eta$\")\n", "ax0.set_title(f\"found $\\eta$ rad_length_frac\")\n", "\n", "a1 = ax1.hist2d(\n", " rad_length_lost,\n", " eta_lost,\n", " density=False,\n", " bins=nbins,\n", " cmap=plt.cm.jet,\n", " cmin=1,\n", " vmax=vmax * stretch_factor,\n", " range=[[0, 0.6], [2, 5]],\n", ")\n", "ax1.set_xlabel(f\"$x/X_0$\")\n", "ax1.set_ylabel(f\"$\\eta$\")\n", "ax1.set_title(f\"lost $\\eta$ rad_length_frac\")\n", "# ax1.set(xlim=(0,4000), ylim=(-1000,1000))\n", "\n", "plt.suptitle(\"radiation length fraction and eta endVelo\")\n", "plt.colorbar(a0[3], ax=ax1)\n", "\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 50, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Parameterisation for rad_length_frac:\n", "intercept= 0.0\n", "coef= {'1': 0.19436659767031483, 'x': -0.0002934072350452474, 'y': -0.00011050606423805956, 'tx': 0.22076931190684765, 'ty': 0.07813480395105185, 'qop': -11.837351954632934, 'x^2': -2.0111718729828063e-06, 'x y': 9.09756297060689e-06, 'x tx': -0.0071441555126283435, 'x ty': 0.01601317688273909, 'x qop': 0.6215416537184049, 'y^2': -6.877429438893236e-06, 'y tx': -0.03246512737106646, 'y ty': 0.009236443466565325, 'y qop': -0.7656982001003853, 'tx^2': 3.2721253340872036, 'tx ty': 6.990115604991107, 'tx qop': -728.51712427984, 'ty^2': -4.166569974702085, 'ty qop': 587.2224646668662, 'qop^2': 6.145749766573158, 'x^3': -1.6562683153707752e-07, 'x^2 y': -1.2001632657575304e-06, 'x^2 tx': 0.00040230376049067056, 'x^2 ty': -0.0003226579491818191, 'x^2 qop': 0.002985359658961018, 'x y^2': -5.185010866064232e-07, 'x y tx': 0.0029804615695217043, 'x y ty': 0.0012164498083769432, 'x y qop': 0.0029776397974184527, 'x tx^2': -0.3154073863149152, 'x tx ty': -0.4010820725541255, 'x tx qop': -2.1214929005307828, 'x ty^2': -0.6824284368708862, 'x ty qop': -14.944686049485084, 'x qop^2': -110.18157556408426, 'y^3': 2.4713398261155817e-07, 'y^2 tx': -8.884404196152218e-05, 'y^2 ty': -0.0005372353581591938, 'y^2 qop': -0.0036932151114263956, 'y tx^2': -1.5419692025080922, 'y tx ty': -0.11860718642094348, 'y tx qop': 11.982547200064268, 'y ty^2': 0.38784142586108383, 'y ty qop': 3.211190671760084, 'y qop^2': 46.063719746472096, 'tx^3': 79.89361864801454, 'tx^2 ty': 468.84573407547805, 'tx^2 qop': -100.98490043430814, 'tx ty^2': 186.00390862227553, 'tx ty qop': 7.2472113332685115, 'tx qop^2': 0.3801287848125358, 'ty^3': -92.81882058287101, 'ty^2 qop': 31.69467333269877, 'ty qop^2': 0.2896764307992765, 'qop^3': -0.004367926967413577, 'x^4': 4.0132150047611503e-08, 'x^3 y': -2.8422402209571374e-09, 'x^3 tx': -0.00012178298386800179, 'x^3 ty': -3.6970947007830546e-05, 'x^3 qop': -0.00037942687285408283, 'x^2 y^2': 7.386681488696922e-09, 'x^2 y tx': 4.123194972294186e-05, 'x^2 y ty': -5.2067299577629456e-05, 'x^2 y qop': 8.539712169408631e-06, 'x^2 tx^2': 0.1375456025515689, 'x^2 tx ty': 0.07276853278925399, 'x^2 tx qop': 0.5913446466084463, 'x^2 ty^2': 0.05481754041796982, 'x^2 ty qop': 0.12778897954363583, 'x^2 qop^2': 7.599038189537376, 'x y^3': 4.017000065914544e-09, 'x y^2 tx': 3.555973574975724e-05, 'x y^2 ty': 2.5015913422521407e-05, 'x y^2 qop': 0.0006923060020653793, 'x y tx^2': -0.0728178744632082, 'x y tx ty': -0.030849602066386872, 'x y tx qop': -0.053150965744121384, 'x y ty^2': -0.054271488256181, 'x y ty qop': -0.3845284787500334, 'x y qop^2': 3.5382762706451607, 'x tx^3': -68.2763499577661, 'x tx^2 ty': -35.09421296048657, 'x tx^2 qop': -209.5850008849592, 'x tx ty^2': -25.837301721685865, 'x tx ty qop': -72.99376890218326, 'x tx qop^2': -20.075729532656002, 'x ty^3': 24.168201409024896, 'x ty^2 qop': 77.70073016158322, 'x ty qop^2': -33.972956067372266, 'x qop^3': 4.163350823159558, 'y^4': 2.935763632194721e-09, 'y^3 tx': -3.516699649708244e-05, 'y^3 ty': -8.2068742131014e-06, 'y^3 qop': -1.523748721155016e-05, 'y^2 tx^2': -0.009411124024499173, 'y^2 tx ty': 0.06338416840865246, 'y^2 tx qop': -0.6891109275638628, 'y^2 ty^2': 0.008253872692153719, 'y^2 ty qop': 0.02346960416570875, 'y^2 qop^2': -0.4386626212085398, 'y tx^3': 32.94713780622537, 'y tx^2 ty': 19.368744361202825, 'y tx^2 qop': 17.426884247940652, 'y tx ty^2': -27.528264177777427, 'y tx ty qop': 341.8447125490805, 'y tx qop^2': -32.944204297917075, 'y ty^3': -3.4360959680136522, 'y ty^2 qop': -10.352300738708175, 'y ty qop^2': 79.19546293155958, 'y qop^3': 0.22012908047276727, 'tx^4': 12563.778412336536, 'tx^3 ty': 749.5755851917718, 'tx^3 qop': -9.94600150329609, 'tx^2 ty^2': 1299.398342139433, 'tx^2 ty qop': 0.1819109935589588, 'tx^2 qop^2': -0.05652058959156786, 'tx ty^3': 402.1468237499123, 'tx ty^2 qop': 0.5288025377553879, 'tx ty qop^2': -0.09503849406227534, 'tx qop^3': 0.005526136730138975, 'ty^4': 490.8007853758815, 'ty^3 qop': -0.032786661645227415, 'ty^2 qop^2': 0.21306953809390777, 'ty qop^3': 0.00041305791547802015, 'qop^4': 5.0463795291232325e-05}\ "r2 score= 0.052596793429776634\n", "RMSE = 0.12333913500444761\n", "\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAikAAAHLCAYAAAAN9wgKAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABqTElEQVR4nO3de3xTdZ4//le5lUuhoaCIItBE7hchKY5zE50mjndG2sDud247O9rIdL+zu442sisCOhJSndmd2emMaZ0dd2b2+xtJcAQv45ig4owzjjQRtVwEegoKIgJpKIVeaHN+f+CJTZq0Sc5Jzkn6ej4eeejJub1PE3rePZ/P5/0pEEVRBBEREZHGDFM7ACIiIqJ4mKQQERGRJjFJISIiIk1ikkJERESaxCSFiIiINIlJChEREWkSkxQiIiLSJCYpREREpElMUog0LBAIoL6+fsiePxNCoRDq6+shCILaoRDRIJik0JAVCARgt9thtVphMBhQW1urdkgRgiDAarXCZDLB5XLlxPmln6fJZILJZMpwhOmpr69HaWkpbDYbQqGQYsfV8neJKJeNUDsAIjUEAgGUl5ejtbUVAGC329Hc3KxyVJ/R6/Vwu90oKCjImfMbjUYAQG1tbeT/1SYIAkpKSqDT6QAAVVVVaG5uVjSJ0Pp3iSiXMUmhIcnhcKCkpCSy7HQ6VYwmf2glOZFYrVa43e5IkgIAkyZNUvQc/C4RZQ6be2hICgQCaodAGWa1WrPyOfO7RJQ5TFJoSKmvr4fVaoUgCJF+F1arFT6fL7JNKBSCzWaD3W6HxWKBxWKJWu/xeDBx4kQUFBREblA+nw9WqxUFBQWwWq2R49TX18NkMsHj8cDn88FkMkVt05d0XumVSpOEx+OJxFlfX4+JEyfCZrNF1tfW1kauyWQyxe0MK+f8yQgEArBarbBYLDAYDLDb7VHnTuVnBSASp8Fg6JeQeDyeyLLNZoubsASDwcg5BzpPIoN9l5T+TKTvIpMiGlJEoiFIr9eLer2+3/t+v1/U6XSi3++PvOdyuUQAotPpjLxXVVUlAojarrm5WQQgVlZWRpYrKytFAKLZbBZrampEv98f2bfv8Zqbm0WdTid6vd7Ie06nUwQgGo3GAa/F7XaLer1eBCBWVVWJNTU1otFojOxXU1Mj9v2n7vV6RQBR55Jz/ljx9vH7/aLZbI6KWYpXOn+yPytRFEWj0SjW1NRE9gUgAhD1en3kmNJ1Nzc3R+0rXZfZbBadTueA50lGvO+SUp+JXq+Pil+n04k6nS7lGIlyFZMUGpISJSlGozHqZtr3/b43POkm0zdJaW1tjUpSRPGzm3HszU+6SUoqKyvjnjfZJEG68Uo37r7MZnPUjU2Ks++2cs8/2D5GozHqZyWKF2+4AMTW1lZRFJP/WUnb9T2elGT0vaEPlqS4XK4Bz5OsRN8luZ+J0Wjs97OQjin9zIjyHZt7iD4lCAICgUDczp/SY/p0hwP37bgpCQaDkfNKTQPpko6/bNmyfuvcbjf8fn9kubGxEQAiQ3CVOP9ApJ+rw+GINIn0bVqR4pEM9LMCgF27dvVbLx0vlaaQsrKyfudVsnaK3M8kEAjAbDZH7VdTUwNRFOP+jIjyEUf3EH1qoBucdEPLRAEw6Zh6vV72seLdvHQ6HXQ6HTweD55++ul+N00lzx+P9HN1u92KHE+K3+fz9UsotTa6CEjvM5F+ZkxGaKjjkxSiGPGKfEk3i75DTZUiJQl9nxYofXyTyQRBEOB2u1FTU5P18/f9r1yVlZUwm81wOBzw+XwIhUJwOp2oqanJWKKltGQ/E1bFpaGOSQrRp6S/wvuO5JFIiYvBYFD8vNKNte/jfyVZLBaUlJT0uxFm6/zS8T0eT9z18X7eg3G73TCbzZGy/U6nM6fqkwz2mUjfxURPn5i80FDBJIWGpGAw2O/JgV6vh9FojAwp7auxsRE6nQ5VVVUAPisI1nc76f9TLbcuNSXV19fH3VdO+XbpWvo2G0jHk64/k+cHEOlXYbfb+zWppTsvkFSkraamBjU1NQM288TGf/r0aQDKPTmK910aSKqfSWwSZ7fbM/JEj0iLmKQQ9SFVJ+1bz0JqTmhoaIjcWKSbot1uj9TBkDrV+ny+SCfUZG5eOp0u8he1yWSCz+eDIAiROiKCIAxas0S6ycXekKWbmcfjQX19Perr6yPHDQQCkacbcs+fyvVZrVbU1tbCYrGgubk5ksQke6OXbtzScerr66Pqokikp14ulyvSOXggoVBI0SYvuZ+J9GTIYrHAarVG6qkYDAb2VaGhQ+3hRUTZ1LcmBj6tYdG3NoUoXhwOKg3JraqqEquqqvoNnxXFi8NBpboV0tBRvV4fqfHh9/sjQ5f1er3o9XrF1tbWqPP3HWLqcrkitTWMRmOkTkZNTU2/YbR99a3Jodfr+w2tdblcok6nixxLFC8O2e0bdyrnl2KP94q9tr5DZZ1OZ8I4U/lZSXHFO79er486p9FoFHU6XaR2itvtjgx9loZFx56n789koGsd6Luk1GfidrsjPxej0djvu0qU7wpEURQzlwIRESkrEAjg6aefxtq1axEMBqOegLjdbhgMhoR9PYgot3AIMhHlDGlUTGtra2QYb196vT6tjrhEJI8gCBkZXZeXSYpUOEqv1yMUCsFisaCysjKlY/Rtw66qqmIbMJEGSJ2T7777bqxduzbSN0gQBPh8PjQ3N+fUKB+iXFVQUBC1bDQao0YIKnEfBoC8a+4RBAEGgwF+vz/yC0yazEwamTHY/na7HaFQCC6XK2fqLhANFbW1tXA4HFEdUo1GI5xOZ78KrUSkvPr6evj9/qiSDGazOeqPBjn34b7yLkmRRlV4vd7Ie/X19bDZbBjsUgOBAMrLy7Fq1aq0y58TUXZIfVH4hwRRdlkslqh7bLz1QHr34Vh5NQQ5FApFDf+U9K05MNC+5eXl0Ov1TFCIcoBOp2OCQpRlHo8HjY2NsFqtce+pcu7D8eRVkiJN0hX7i0t63DRQ5ic18bA9m4iIKD6v14tQKASPxwObzYaJEydGdVaXcx+OJ686zkqd6hJ1ch2olLSU3Xm9XtjtdgiCgLKysgH7pZw6dQovvPACpk6ditGjR6cd96hRozBq1Ki09ycioszq7u5Gd3d32vt3dnbi+PHjuPXWWzF58mQFI4v2wQcf4NSpU7LjBYCpU6di1qxZUe+5XC64XC4EAgG4XC7U19dHCjPq9XpZ9+G4MlWA5cyZM2JLS4u4Y8cOsaWlJVOniVJTUxMpshQLnxZVisfv90cVsBLF6IJRfYtD9fXUU08NWOyJL7744osvvvq+nnrqKcXuebGOHDkijlQw1pEjR4pHjhwZ8Jxut1sEIFZWVoqimP59OBHFn6S88sor2Lx5M3bs2AEAEEUxMlTJaDSitrYWN9xwg9KnBfBZGexEpa0TPRGRMjubzRbZRuqbYrFY4HA44jYDTZ06FQDw05/+FFdfffWg8X3zm9/Eb37zm37vx3uSUlFRga1btw56zHS2T3bb9vZ2LF++HDt37kRRUZFqcaS6bSrbp3qNWvhcUt02k9eohZ8Hv6fZiyOT22r9e5royUSi3+sA0HkhjF+/dQzHQh24piiE2ofuj9w3MuHUqVO4AKASQAmAXjnHArDtwgWcOnUK06dPT7hdZWUlKisrI9NSpHsfTkTRJGXr1q1wOBxYvXp1v0mwgsEg/H4/7r//ftxzzz246667lDw1gM8uPtGEaIl+OIkeS0nDGRM9npKaeK6++mpcd911g8Y3bty4pLYDgDFjxgw4aZqc7ZPdtq2tDQCwZMkSTJgwQbU4Ut02le1TvUYtfC6pbpvJa9TCz4Pf0+zFkcltc/V7muj3+pnzF7B+exN6JhXi5/+wACcP7kYtIKtrQLIuB3CFzGOMS2Fbi8US6ZeS7n04EUWTlMbGxkinmXjKy8tRU1ODBx54QMnTRki9h2OTCmnZZDINuF9zc3Pc9ZxxlIiIkvXJ2U489OwenOvugWPlYpROHoeTB9WOKrOk+2i69+FEFB3dI01fP5hMdRrS6XQwGo39eg9LGd6qVasS7mc2m/uV05YywVR "text/plain": [ "<Figure size 640x480 with 2 Axes>" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "rad_length_frac = found[\"rad_length_frac\"]\n", "# @ z = 9400.mm or 770.mm\n", "state = 1\n", "\n", "if state == 1:\n", " slopex = found[\"ideal_state_770_tx\"]\n", " slopey = found[\"ideal_state_770_ty\"]\n", " x = found[\"ideal_state_770_x\"]\n", " y = found[\"ideal_state_770_y\"]\n", " qop = found[\"ideal_state_770_qop\"]\n", "elif state == 2:\n", " slopex = found[\"ideal_state_9410_tx\"]\n", " slopey = found[\"ideal_state_9410_ty\"]\n", " x = found[\"ideal_state_9410_x\"]\n", " y = found[\"ideal_state_9410_y\"]\n", " qop = found[\"ideal_state_9410_qop\"]\n", "\n", "data = ak.zip(\n", " {\n", " \"rad_length_frac\": rad_length_frac,\n", " \"x\": x,\n", " \"y\": y,\n", " \"tx\": slopex,\n", " \"ty\": slopey,\n", " \"qop\": qop,\n", " }\n", ")\n", "lin_reg, features, xx0_test, xx0_predicted = fit_linear_regression_model(\n", " data,\n", " \"rad_length_frac\",\n", " [\"x\", \"y\", \"tx\", \"ty\", \"qop\"],\n", " 4,\n", " include_bias=True,\n", ")\n", "\n", "nbins = 100\n", "vmax = 50\n", "\n", "a0 = plt.hist2d(\n", " xx0_test,\n", " xx0_predicted,\n", " density=False,\n", " bins=nbins,\n", " cmap=plt.cm.jet,\n", " cmin=1,\n", " vmax=vmax,\n", " range=[[-0.1, 0.6], [-0.1, 0.6]],\n", ")\n", "plt.plot([-0.1, 0.6], [-0.1, 0.6], marker=\"\", alpha=0.8)\n", "plt.xlabel(f\"True $x/X_0$\")\n", "plt.ylabel(f\"Predicted $x/X_0$\")\n", "plt.title(f\"found rad_length_frac\")\n", "# ax1.set(xlim=(0,4000), ylim=(-1000,1000))\n", "\n", "plt.colorbar(a0[3])\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 51, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAj8AAAHLCAYAAAAnR/mlAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAx5klEQVR4nO3dTYgjeWLn/V+u/XROeZLqqKy5zMKCOwQ+F6Hsi/0wNVARDMxcFkqqunkMQ0kMeVjoxRlOXyrrYhGJrwkTSgyN2Uulonaf52CDLbVxDXjg2UqFy5eHZ5dVlMGHWXimldFF4prstUd7qImoVOollXpJhTK+H0i6S1JE/BWhVPzy/7rW6/V6AgAAyIl/s+wCAAAAXCfCDwAAyBXCDwAAyBXCDwAAyBXCDwAAyBXCDwAAyBXCDwAAyBXCDwAAyBXCD3ItDEPV6/VlFyP34jhWvV5XFEXLLgquiN8hrCLCD3IpiiKVy2UVi0X5vt/3+J07d+b6Zb6IfU4qDEO5rqtisahisXjtx59EvV7XJ598omq1qjiOZ9pX8n7L5bIKhYL29/fnU8gpLfPaL9qw36FWq6VCoaC1tTWtra2pWCwqCIKBbYMgULFYTF8ThuHYY63C5xirhfCDXDJNU41GY+DxOI4Vx7E6nc5U+42iaOAGPus+Z2FZlh4/fnzpzeU6XTxHlUpFlUpl5v2GYagHDx7I8zw1Gg2VSqVrPedZu/aLNux3yLZttdvt9N/ValWlUmlg21KppGq1KsMw1G63ZVnW2GNl8XOM1Ub4Ac6xLEu9Xk+e5021fblcVrfbnes+Z3XZjeW6DTtHd+/enXm/tVpNm5ub6b89z+ur1Vu0LF77ZTAMQzs7O5I09A+MRLvd1u7u7sT7zdrnGKuN8APMSblc5i/TSyzyHC3z3HPt+yWhptVqjTwvR0dHc6nxA6ZB+MFKCIJAjuOo1WqpXq/rzp07qlar6fP7+/uqVqtpv4BhfSziOFa1Wk1/hvUHCYJA5XJZ5XJ54LlxxwiCIP2Sr1arfTfDcftMyuS6rhzHSd/j+efr9Xrad6LVaqV9JYbt76rCMFS5XJbjOCoUCnJdd6ZjJ+e2UCgMBIJx5yjR7XbTY076Huv1usrlsqIoSvuhlMvl9DzO+7OTXKfz1zfL136e13iS3yHpfe1P0txVq9UGng+CQLZtyzCMicp5mcvOJTCgB2Rco9HomabZk9SrVCq9nZ2dnmVZPcuyer1er7ezs9M7/1FuNps9Sb1ms5k+1ul0eoZh9D3meV5PUrqfTqeTPmbbdl8ZJjlG8ppOp9N33FH7bLfbPcMweu12O33M9/2epJ7neen2pVIp3X5nZ6fXbrd7lUql73WXOf8+zx//fJkajUZ6jqc5tmVZvZ2dnXRbST1JPdM0030OO0e9Xq/vHHmeN9V7NE2zZ5pm32Pz+uyYptlXZsMweoZhpP/O6rWf5zWe5HfovPOfgZOTk77nLMvq289l5UyM+hxfdi6Biwg/WAnJl2xycz3Ptu2+G9HJycnAa0ul0sANqNcb/mU67GY1yTFG3dhH7dOyrKFlsiyrbz/JjeDiF/mwfY4y7H1altV3w+j13t/Uz9+sJj128rrz+0tuoOfPx2Xhx/f9qd/jsPBzft/TfnYsyxp4/8k+k/OU1Ws/z2t8ld+hhG3bA+czCVFXLeeoY016LoHzaPbCSkiqxz/99NOB5xqNRt8Ik+PjY0lKR95EUZQ2fUzrsmNcVRRFCsNwaCfOpEnmYmfd800EiYsdbK96/FqtljbLnG/mSN7fpMd+9erVwPPJ/q7SF2Zra2vguLPO/TPrZycMQ9m23bfdzs6Oer3e0PNymeu69vO8xtP+DiVNV+ebEj3P6+vofNVyDnuPVzmXgCT95rILAFzFsC9owzBkGIaCINDz588HbnLJzdM0zZmOO+4YVzUuECQBYJET/iXHHzca5yqS89FqtQZuRFkZpTPNZyc5T9OEnFGu69rP8xpP+ztk27ZM01QURarX66pUKjo6OtKbN2/mUs5l/x5hdVHzg5UXRZGKxaKiKFKj0UiH2Z5/Xpq+lmSSY0xrWM1RcqM9P2x73pJzMq8bQ6lUkm3bqtVqarVaiuNYnudpZ2dnptC5aJN+dhZxA130tZ9n2Wf5HUpqfzzPUxAE2tra6guT8yjnsn6PsLoIP1h5juNoc3NzZCBJbr7nmzfmfYyrSmpDho1ISb7IC4XCXI41THJOhs2+O6pcl2k0GrJtO13uwPO8zM9vc9l1Ta7TqFqJaW7Y13Xt53mNZ/kdqlQqafPlkydP+kbazVrOZf8eYXURfrDSkuHN5/+STL70kr9Sk+rver0+9C/Ey/rtTHKMq+xPev+Fb1lWuu/zjo+PZRjGQudASfqwuK470HQw7VIM5XI5rT3Z2dkZ29x18Rx9+eWXkmarnet2u1fa/qqfnYs3WNd1B2oVsnTt53mNZ/0dOt/H5+KMz7OUc9m/R1hdhB+shOTL9eKXbHLzCYJA9Xpd9Xo9rWYPwzD9azL5y75YLKrVaimKovR1URSl85UMCzWTHCOO4/QvTN/30w6io/Ypva9NMAyj7y/hpLno8PAwvSnPEghGOT8Lb7FYVLlc1v7+vhzHUafTSW9Ikx47CQfJfur1et/8N4lR52iUOI5nfv+zfnaS2ivHcVQul9P5gAqFQnqNsnjt53mNL+5r3O/QMEkAGRZEJi3nKJOeS6DPsoebAZc5P1eLaZoDw6F93+8ZhtEzTTMdUlupVHqGYfQNsfV9P92PZVnp/C07Ozu9TqfTa7fb6Zwn+vXQ32SY7aTHsCyrZxhGOj/JuH32eu+HVidDiCuVSq9SqfQN+W232+mQXdM0e81ms3dycpIOI9eQIcrnnZ+zZdjxPc8beW6vcuzkXCaPn/8xTbPvmBfPUaPRSIc1J0OeLx5n2DD1Ue+xUqmkc8jM67PTaDTSc3FxjppR72vZ1z4xr2ucnK9xv0PjVCqVsa+5rJzjPseXnUvgorVer9dbXLQCkAdhGOr58+fa3d1Vt9vtq7FpNBoqFApz6y8FALNiqDuAmSQjpk5OTtKh4+eZpslSAwAyhT4/AGaSdDR98uRJXx+fZG4X3/fpdAogU2j2AjCz/f191Wq1vk7FlmXJ87xLO6wCwHUj/ACYm6SvT5YnNgQAwg8AAMgV+vwAAIBcycVor1/84hf6q7/6K/32b/+2bt26teziAACACbx7907/+I//qO9973v61re+Nbf95iL8/MVf/IX+4A/+YNnFAAAAU/j888/1wx/+cG77y0X4+fa3vy1JOjw8HLve0CwePnyoFy9eLGTfN+kYp6enun//vl6+fKmNjY2FHWfR74NrMbmbcK64Fvk5BtciW8cIw1BPnjxJ7+PzsrDw8+Mf/1hbW1v60Y9+tKhDTOwb3/iGJOl3fud3FhZ+bt26tbB936RjvH37VpJ079493b59e2HHWfT74FpM7iacK65Ffo7BtcjWMU5PTyV9uI/Py0I6PL9580a+7zOdPQAAyJyFhJ9PPvlEvu/r6OjoytsGQaBisai1tTUVCoWh0+KHYZiurlytVi9dGRoAACCxsGavJ0+eXHmber2udrstz/MkSa7rynEcdTqddNK0ZB2hdrudVrUVCgV1u12m0AcAAJda2Dw/3/ve9/T69esrbRPHsXzfl23bsm1bh4eHktS3XlC1WpVt231tjEkNEAAAwGWmrvl5/PjxyOfiOFar1dLR0ZHu3bs38T4v9hFKVodOgk6y36RmKLG1tSXpfc0RtT8AAGCcqcNPo9GY6DV/8id/Mu0hFASBPM9Lm7yOj48laWDdoCQcNZvNpYWf7e1tjpEhi34fXIvJ3YRzxbXI3zEW7aacp1W9FlOv7fXo0SN5nqfNzc2B5zqdjur1un7yk59MXTDXdVWv13V4eKhSqSTpfc1OtVpVs9kcWCl6bW1NlmWp3W4P7OunP/2p7t+/r7/8y7/U7/3e701dpvX1da2vr0+9Pd4PI/3444/11VdfLXQYKS7HtcgOrkV2cC3m6+zsTGdnZ1Nv/3d/93f6/ve/r5cvX+o73/nO3Mo1dc1PtVrVJ598MvQ "text/plain": [ "<Figure size 640x480 with 1 Axes>" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.hist(\n", " xx0_test,\n", " bins=100,\n", " density=True,\n", " alpha=0.5,\n", " color=\"darkorange\",\n", " histtype=\"bar\",\n", " label=\"test\",\n", " range=[0, 1],\n", ")\n", "plt.hist(\n", " xx0_predicted,\n", " bins=100,\n", " density=True,\n", " alpha=0.5,\n", " color=\"blue\",\n", " histtype=\"bar\",\n", " label=\"predicted\",\n", " range=[0, 1],\n", ")\n", "plt.xlim(0, 1)\n", "# plt.yscale(\"log\")\n", "plt.title(\"radiation length fraction endVelo\")\n", "plt.xlabel(f\"$x/X_0$\")\n", "plt.ylabel(\"a.u.\")\n", "\n", "plt.legend()\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 49, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Parameterisation for rad_length_frac:\n", "intercept= 0.0\n", "coef= {'1': 0.24562140630345797, 'x': -0.0006147519300623754, 'y': 0.00013149294035773212, 'tx': 0.40993139584024857, 'ty': -0.17613130993997955, 'qop': -2.431621847221756, 'x^2': -2.537929422256502e-05, 'x y': 1.6503609236837168e-05, 'x tx': 0.014114349319941343, 'x ty': -0.015914217992780068, 'x qop': -0.3138750804691517, 'y^2': 8.004082295422125e-06, 'y tx': -0.0007103518977325451, 'y ty': -0.012102735170731998, 'y qop': -0.01613366131755978, 'tx^2': -1.1805199392607937, 'tx ty': 2.630285174738731, 'tx qop': 0.2259863956524285, 'ty^2': -0.9000770143303645, 'ty qop': 0.035317808178541486, 'qop^2': 0.506257028947111, 'x^3': -1.054613771280784e-06, 'x^2 y': -7.542919321773425e-07, 'x^2 tx': 0.0016046817739745617, 'x^2 ty': 0.0041420285105903745, 'x^2 qop': 0.021289920901612722, 'x y^2': 1.6630958165414549e-06, 'x y tx': -0.0029322410866066016, 'x y ty': 0.00026646565396820643, 'x y qop': 0.008117926325135527, 'x tx^2': -0.6007987905052633, 'x tx ty': -2.98887506065023, 'x tx qop': -13.962516223136573, 'x ty^2': -0.8388571705095013, 'x ty qop': -3.036359977877559, 'x qop^2': 0.055324587199241985, 'y^3': -8.030979246897996e-08, 'y^2 tx': -0.002767856037518542, 'y^2 ty': 0.00010106006074174414, 'y^2 qop': 0.0087584836303898, 'y tx^2': 2.52086358169213, 'y tx ty': 1.771954217213756, 'y tx qop': -3.4960101511969945, 'y ty^2': -0.02778479572722497, 'y ty qop': -8.306825049361231, 'y qop^2': 0.011633837171348742, 'tx^3': -0.1376005395953692, 'tx^2 ty': 0.08788148993360365, 'tx^2 qop': -0.03627295634145932, 'tx ty^2': 0.03876102551307681, 'tx ty qop': -0.008986008881728739, 'tx qop^2': 9.888031370307247e-05, 'ty^3': -0.3359335900376418, 'ty^2 qop': -0.02333508540052829, 'ty qop^2': 9.573184641700597e-06, 'qop^3': -4.597612275351735e-06, 'x^4': -6.108581640518196e-09, 'x^3 y': -6.564809273967853e-08, 'x^3 tx': 1.5150502404637223e-05, 'x^3 ty': 3.104184179880676e-05, 'x^3 qop': 8.671096332828199e-05, 'x^2 y^2': -7.650740829062697e-09, 'x^2 y tx': 0.00011735105575993998, 'x^2 y ty': -3.2175113096499786e-05, 'x^2 y qop': 6.774919445144169e-05, 'x^2 tx^2': -0.011724689654149013, 'x^2 tx ty': -0.043137944250594654, 'x^2 tx qop': -0.0234816944437721, 'x^2 ty^2': -0.015910369617829878, 'x^2 ty qop': 0.47749147066068204, 'x^2 qop^2': 10.59413029512916, 'x y^3': 3.234245482541809e-08, 'x y^2 tx': 5.3043264426122505e-05, 'x y^2 ty': -5.8168920071072217e-05, 'x y^2 qop': -0.00018640913107006996, 'x y tx^2': -0.06797789037420698, 'x y tx ty': 0.044976962291121, 'x y tx qop': -0.5074481714386293, 'x y ty^2': 0.03709689469083238, 'x y ty qop': 0.42061549746065124, 'x y qop^2': 5.878761197856953, 'x tx^3': 2.928172660629093, 'x tx^2 ty': 15.518053547545511, 'x tx^2 qop': -0.4179045108625407, 'x tx ty^2': 2.863496001226073, 'x tx ty qop': 0.03316794782068423, 'x tx qop^2': 0.011365298564370092, 'x ty^3': -8.111988620415804, 'x ty^2 qop': -0.21644263264392571, 'x ty qop^2': 0.00843273298050687, 'x qop^3': -0.00016877332507871679, 'y^4': 6.023639254060242e-09, 'y^3 tx': -1.594750695611824e-05, 'y^3 ty': -1.437326358022517e-05, 'y^3 qop': -2.1442523502287515e-05, 'y^2 tx^2': -0.04753636763140216, 'y^2 tx ty': 0.018902432682699272, 'y^2 tx qop': -0.31014755455483073, 'y^2 ty^2': 0.011363729381750311, 'y^2 ty qop': 0.01681916414935703, 'y^2 qop^2': -5.242995516210874, 'y tx^3': 12.07113658017565, 'y tx^2 ty': 2.5063182953688687, 'y tx^2 qop': 0.05008258889078114, 'y tx ty^2': -5.840846904263012, 'y tx ty qop': -0.2222888950764648, 'y tx qop^2': 0.008403354774191822, 'y ty^3': -2.8576703684113554, 'y ty^2 qop': 0.156828601288708, 'y ty qop^2': -0.010489395507591974, 'y qop^3': 1.6753570524717795e-05, 'tx^4': -0.01910867380555672, 'tx^3 ty': 0.0837123623190646, 'tx^3 qop': -0.0014826028602247271, 'tx^2 ty^2': 0.0015368291104171875, 'tx^2 ty qop': 8.478961198550636e-05, 'tx^2 qop^2': 1.2316906385746789e-05, 'tx ty^3': -0.0271004544137963, 'tx ty^2 qop': -0.0008124467158542836, 'tx ty qop^2': 1.1953538675306295e-05, 'tx qop^3': -1.6561432689863052e-07, 'ty^4': -0.051346978654440316, 'ty^3 qop': 0.0005598437003601672, 'ty^2 qop^2': -1.806470433082 "r2 score= 0.013501590753704495\n", "RMSE = 0.2590497099709809\n", "\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhoAAAHLCAYAAACOMzExAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABmd0lEQVR4nO3dfXxTZZ43/k95KsUCoYCIINITQHkQIYEZddZhnCao4wMjbfCe+575zeyONMN279ldZ2xkV+TBkZCqs7sz2xnT6s7sPbs7a5O6gs806KAzjkoTUVtAoKc8qiikoRRKS5vz+6Oe2KRJmudzknzer1dekJycc745bZNvrut7XVeBJEkSiIiIiNJghNIBEBERUe5iokFERERpw0SDiIiI0oaJBhEREaUNEw0iIiJKGyYaRERElDZMNIiIiChtmGgQERFR2jDRIMpyHo8HdXV1SoeRUj6fD3V1dRBFUelQiChJTDQoZ3k8HlgsFuj1euj1eqXDSTlRFGEymaDX62G322PaJxuuSV1dHUpLS2E2m+Hz+VJyTPl1m0wmaLVa1NTUpOS4RDS8UUoHQJQuOp0OAFBTUxP4fzqIooiSkhJoNJq0nSMcQRDgcDhQUFAQ8z6ZuibxCL1+lZWVaGtrS1ky4PF4UFZWho6ODgCAxWJBW1tbSo5NRMNjiwbltEx8mJpMJni93rSfJ1XUkmDIwl2/yZMnp+z4VqsVJSUlgfs2my3mFiAiSh4TDaIkmEwmeDwepcPIWpm4fvz5ECmLiQblLZ/PB7PZDIvFAqPRCKPRCJfLFfE5ZrMZWq02UHjpdDoDH2Jms3nYD02n0xk4R11dHSZNmgSz2RzYXlNTEziXXq8PW+ApxyPfUl1r4PF4YDKZYDQaodVqYbFYgs5dV1cHvV4Pp9MJl8sFvV6PgoICmEymsMeT49RqtUOuTyzXz+v1Bs4Z7Tzh1NXVwWQyQRTFQD2LyWQK/IxT/fOQf4eY2BCFkIhyHABJp9MFPeZ2uyWNRiO53e7AY3a7XQIg2Wy2wGMVFRVSdXV10HMGb6+urpYASG1tbVFjcDgckiAIEgCpsrJSqq6ulnQ6XSAu+TiypqYmCYDU1NQUeKytrU3SaDRBj9lstrCvbziRronBYAiKWY5XPn9FRYUEQDIYDFJ1dbXkdrulysrKIddNkiRJp9MFrl1bW5sEQAIgCYIQOGak6ye/LoPBINlstqjnGY4gCJIgCEGPpernIQhCUOwajUbSaDRxxUeU65hoUM4L96Gq0+mCPlQHPz74g0+j0QQlGpIkJZRoyPsBGHI8SZIkg8EQ9AHV0dEx5LkVFRVhY05VoqHT6YISL0kaeP0ApI6ODkmSvkw+Qj/s5aRAJj9v8PHkRGHwtRou0bDb7VHPE4twicbgcyT689DpdEOug3xM+XoRkSSx64TyjiiK8Hg8YYsi5aZzuVhQEATU1NQEdVFUV1cndF55VMXy5cuHbHM4HHC73YH7zc3NABAY3imKYqCpPx3ka2K1WgNdDIO7KeR4ZOFG2Awu6NyzZ8+Q7fLx4ulaWLZs2ZDzpmpujWR/Hh6PBwaDIWi/6upqSJKU8RFIRGrG4a2Ud6J90MkfbPKHmcPhgF6vh8Vigd1uh8PhSHrURrgPIY1GA41GA6fTiWeeeWbIh58cjyAISZ07EvmaOByOlBxPjt/lcg25Xmob9ZLIz0O+XkwoiIbHFg3KW+Emg5I/OOThkIIgoL29HQaDAaIoRiwKTJZ8bFEU4XA4hrSayIlGuobRysdPVWtBRUUFDAYDrFYrXC4XfD4fbDYbqqur05YspVKsPw/OXEo0PCYalHfkb9ShI0yAL5MPrVYLYOCDRKPRoKmpKfBtf/DIhFQxGo0oKSmJ2C0jfzgPbs5PJfn4Tqcz7PZw12o4DocDBoMhMEW6zWaDzWZLKs5MGe7nIf8ORWoBYgJC9CUmGpR3BEGATqcLDHscrLm5GRqNBpWVlQAQ9MFYUVERqN0I3S+ZqbLlOAY3w8vHk1sw5C6durq6sOdKdqpuudbAYrEM6VpKtAXHZDIFWgOqq6ujdpmExn/mzBkAqWnB8Xq9cR0n3p9HaBJmsViCJggjyndMNCgvORwOaDSaoNYJuXm/vr4+8CHT0NAQlFT4fD4IghBoAZBbPux2e6BgMxL5wyr0Q1X+UHI6nairq0NdXV1g/gqPxxM4pvztWq/Xw+VyQRTFwPNEUUxqTg2NRhN0fJPJhJqaGhiNRrS1tQUSkVg/sOUPYPk4dXV1QfNmyOK5fsDAtUtV91GyPw85CTUajTCZTIH5NrRaLWs3iAZTetgLUboMnnsBXwzJHDzssKOjIzBktLKyUqqsrBwyvNNgMEiCIEjV1dVSdXW1VFFRMWTook6nkzQaTWBuiHAGz9sgCMKQYZt2u13SaDSBc0nSwHDQ0OG1drs9cBydTheYy6G6ujpoiKj8msPdol0Tm80WMU632x0Y/isIgtTU1CR1dHQMOZ4kfTnHRLjzC4IQdM7Q6+dwOALDauUht6HnGXxNYn2tlZWVgXkwUvXzcDgcgWui0+mC5tkgogEFkiRJ6U1liCjfeDwePPPMM1i/fj28Xm9QS4TD4YBWq014mDARZRcObyWilJJHbHR0dASGiQ4mCEJCxaVElJ1UmWj4fD5YrVYAiLlKXZ5sSBAE+Hw+GI1GVFRUpDNMIgpDrmlZu3Yt1q9fHygCFUURLpcLbW1tWTP6hIiSp7pEw+VywW63w+l0Bir/hyN/g3K73YE3Na1WC6/XG/MxiCg1DAYDbDYbrFZrUHGnTqeDzWbj3yRRnlFtjUZBQQEqKysDwwmjkadlbmpqCjxWV1cHs9kMlb48orwg12ZkwyRdRJQeWT+81efzweVyDVkDYvA4dyJShkajYZJBlOeyPtGQFzsKfTOTu1AGt3IQERFRZqmuRiNecuFZpAlyhpsK+PTp03jxxRcxffp0jB07NuE4xowZgzFjxiS8PxERpV9vby96e3sT3v/ixYv45JNPcMcdd2DKlCkpjCzYsWPHcPr06aSPM2XKFMyaNSsFESUubYlGZ2cnvF4vRFGEIAiYPXt2Ws7T1tYGABGn/B1uauZXX30VP/jBD1IcFRER5bLf/va3+P73v5+WYx87dgxzrr4al1JwrHHjxmH//v2KJhspTzRee+01bNu2Dbt27QIASJKEgoICAAPdGTU1NbjllltSdj55CuNI0xIP1z8sJ0D19fWYN29ewnGMGTMG3/nOd9DY2JjwMQCgq6sLK1aswO7du1FcXJzUscrLy5OOR23H4fWJjtdneLxG0eX69fn9738fd4vG+yfPocHzCa6ddhkWjv4Mf/93f4fp06cnHU8kp0+fxiUAFQAuT+I4nwFwXriA06dP506i0djYCKvVinvvvXfIwkJerxdutxsPPPAAfvSjH+G+++5LyTnlRCJSy8VwiUZRURGAgSQo2qJPsSgqKkr6GJ2dnQCAJUuWYMKECYrHo7bj8PpEx+szPF6j6HL9+txwww1x7fPSh5+g6YM23LNyLv62bC7+9Mc3ASCprvZYXQlgRhL7q6U2IqVxNDc3B4ozwykrK0N1dTUefPDBlJ1THl0SWosh39fr9Sk7FxER5QdJkuBoPoHfvX0Ud10/Hff9hYARIwqUDisrpXTUyeTJk2N6XioLaDQaDXQ63ZDRJfIUx2vWrEnZuYiIKPdJkoTf/OkIfvf2Ufzvr87C2puVSTJGYqA1INHbyIxHHF5KE43Dhw9j/fr1OHr0aNjtR44cweOPPx4o4IwkWgGnKIrQarVBayXU19cHls2W2Ww22Gw2LtdMREQx8/sl/PK1w/if905i7dcFfOcrswJ1hpk2CsDoJG5Z23Xy2muv4Zvf/GbYbU8++SRMJhNKS0tRUFAAjUaDkpKSwOqNAFBRUYFnnnkm4vE9Hk9gNtCGhgYYjUYYDIZAwiDPNDg4GdHpdHC73bBYLBAEAaIowmKxcKpjIiKKWW+fH0/s/Ahvi2dwv3Eebrk2mVJMksWdaJhMJpw5cybidofDgfb2drjdbjQ3N8Pn80Gr1UIQBOh0OpSWlkY9vk6ng91ujzj1uE6nQ0d "text/plain": [ "<Figure size 640x480 with 2 Axes>" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "rad_length_frac = lost[\"rad_length_frac\"]\n", "# @ z = 9400.mm or 770.mm\n", "state = 1\n", "\n", "if state == 1:\n", " slopex = lost[\"ideal_state_770_tx\"]\n", " slopey = lost[\"ideal_state_770_ty\"]\n", " x = lost[\"ideal_state_770_x\"]\n", " y = lost[\"ideal_state_770_y\"]\n", " qop = lost[\"ideal_state_770_qop\"]\n", "elif state == 2:\n", " slopex = lost[\"ideal_state_9410_tx\"]\n", " slopey = lost[\"ideal_state_9410_ty\"]\n", " x = lost[\"ideal_state_9410_x\"]\n", " y = lost[\"ideal_state_9410_y\"]\n", " qop = lost[\"ideal_state_9410_qop\"]\n", "\n", "data = ak.zip({\n", " \"rad_length_frac\": rad_length_frac,\n", " \"x\": x,\n", " \"y\": y,\n", " \"tx\": slopex,\n", " \"ty\": slopey,\n", " \"qop\": qop,\n", "})\n", "lin_reg, features, xx0_test, xx0_predicted = fit_linear_regression_model(\n", " data,\n", " \"rad_length_frac\",\n", " [\"x\", \"y\", \"tx\", \"ty\", \"qop\"],\n", " 5,\n", " include_bias=True,\n", ")\n", "\n", "nbins = 100\n", "vmax = 50\n", "\n", "a0 = plt.hist2d(\n", " xx0_test,\n", " xx0_predicted,\n", " density=False,\n", " bins=nbins,\n", " cmap=plt.cm.jet,\n", " cmin=1,\n", " vmax=vmax * stretch_factor,\n", " range=[[-0.1, 1.0], [-0.1, 1.0]],\n", ")\n", "plt.plot([-0.1, 1.0], [-0.1, 1.0], marker=\"\", alpha=0.8)\n", "plt.xlabel(f\"True $x/X_0$\")\n", "plt.ylabel(f\"Predicted $x/X_0$\")\n", "plt.title(f\"lost rad_length_frac\")\n", "# ax1.set(xlim=(0,4000), ylim=(-1000,1000))\n", "\n", "plt.colorbar(a0[3])\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 58, "metadata": {}, "outputs": [ { "data": { "text/html": [ "<pre>[nan,\n", " nan,\n", " nan,\n", " nan,\n", " nan,\n", " nan,\n", " nan,\n", " nan,\n", " nan]\n", "-----------------\n", "type: 9 * float64</pre>" ], "text/plain": [ "<Array [nan, nan, nan, nan, nan, nan, nan, nan, nan] type='9 * float64'>" ] }, "execution_count": 58, "metadata": {}, "output_type": "execute_result" } ], "source": [ "endT_x = lost[np.isnan(lost.ideal_state_9410_x)][\"ideal_state_9410_x\"]\n", "endT_x" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "tuner", "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.12" } }, "nbformat": 4, "nbformat_minor": 2}
|
