{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "e8ab5a15",
   "metadata": {},
   "source": [
    "### Fits two lines to par0(x-par1) and par2(x-par3) with MInuit\n",
    "### Adds a soft constraint on par1-par3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "a12892ea",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import math\n",
    "import matplotlib.pyplot as plt\n",
    "import scipy.stats as stats\n",
    "from iminuit import Minuit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "afea4692",
   "metadata": {},
   "outputs": [],
   "source": [
    "# This is a nice utility to put a ROOT style stat box on a histogram plot\n",
    "def statBox(ax, entries, binEdges, x=0.96, y=0.98, fontsize='medium'):\n",
    "    \"\"\"\n",
    "    Put a stat box on the histogram at coord x and y\n",
    "    font = medium is appropriate for 1x1.  Other choices are\n",
    "    size in points, 'xx-small', 'x-small', 'small', 'medium', 'large', 'x-large', 'xx-large'\n",
    "    \"\"\"\n",
    "    en = len(entries)                   # number of entries\n",
    "    ov = (entries>binEdges[-1]).sum()   # overflows\n",
    "    uf = (entries<binEdges[0]).sum()    # underflows\n",
    "    mn = entries.mean()                 # mean\n",
    "    sd = entries.std()                  # standard deviation\n",
    "    textstr = 'N=%i \\nOverflow=%i \\nUnderflow=%i \\n$\\mu=%.2f$ \\n$\\sigma=%.2f$' % (en, ov, uf, mn, sd)\n",
    "    props = dict(boxstyle='round', facecolor='white', alpha=0.5)\n",
    "    ax.text(x, y, textstr,\n",
    "            transform=ax.transAxes,\n",
    "            bbox=props,\n",
    "            fontsize=fontsize,\n",
    "            horizontalalignment='right',\n",
    "            verticalalignment='top'  ) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "9f942356",
   "metadata": {},
   "outputs": [],
   "source": [
    "class myFit:\n",
    "    \"\"\"\n",
    "    Wrapper around a minuit chisq with a Gaussian constraint.\n",
    "    y1fit = par[0]* (x1-par[1])\n",
    "    y2fit = par[2]* (x2-par[3])\n",
    "    chisq = (y1fit-y1meas)**2/s1**2 + (y2fit-y2meas)**2/s2**2 + (par[1]-par[2])**2/sigma**2\n",
    "    Where sigma is something small enough to be an effective constraint\n",
    "    \n",
    "    Instantiate as \n",
    "    f = myFit(x1,x2,y1meas,y2meas,s1,s2,sigma=1e-4)\n",
    "    (default sigma is 1 micron)\n",
    "    \n",
    "    Fit as (MINOS not needed for this)\n",
    "    m = Minuit(f.CHI2, pinit, name=(\"p0\",\"p1\",\"p2\",\"p3\"))\n",
    "    m.errordef = Minuit.CHISQUARED\n",
    "    m.migrad()\n",
    "    \"\"\"\n",
    "    \n",
    "    def __init__(self,x1,x2,y1m,y2m,s1,s2,sigma=1e-4):\n",
    "        self.x1  = x1\n",
    "        self.x2  = x2\n",
    "        self.y1m = y1m\n",
    "        self.y2m = y2m\n",
    "        self.s1  = s1\n",
    "        self.s2  = s2\n",
    "        self.s   = sigma\n",
    "        \n",
    "    def CHI2(self, par):\n",
    "        y1fit = par[0]* (self.x1-par[1])\n",
    "        y2fit = par[2]* (self.x2-par[3])\n",
    "        t1    = ( (y1fit-self.y1m)**2/(self.s1*self.s1) ).sum()\n",
    "        t2    = ( (y2fit-self.y2m)**2/(self.s2*self.s2) ).sum()\n",
    "        t3    = (par[1]-par[3])**2/(self.s*self.s)\n",
    "        return t1+t2+t3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "5f7a60ff",
   "metadata": {},
   "outputs": [],
   "source": [
    "# --------------------------------\n",
    "# Here are the needed constants\n",
    "#---------------------------------\n",
    "nDetectors = 4\n",
    "w          = 0.005  # 50 micron is strip width\n",
    "s          = np.full(nDetectors, w/np.math.sqrt(12))  # resolution in each hit\n",
    "xdet       = np.array([2., 3., 5., 7.]) # x coordinates of detectors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "b5c9e119",
   "metadata": {},
   "outputs": [],
   "source": [
    "# read all data into numpy arrays\n",
    "data = np.loadtxt(\"straightTracks.txt\")\n",
    "xv   = data[:,0]   # true xverteces (called X0 in the exercise pdf)\n",
    "yv   = data[:,1]   # true yverteces\n",
    "nev  = len(xv)     # number of pairs of tracks\n",
    "trk1 = data[:,[2,3,4,5]]   # hits for track number 1\n",
    "trk2 = data[:,[6,7,8,9]]   # hits for track number 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "d3764c97",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now do the constrained fit (store results in arrays)\n",
    "fitVertex1 = np.empty(nev)  \n",
    "fitVertex2 = np.empty(nev)  \n",
    "pull1      = np.empty(nev)\n",
    "pull2      = np.empty(nev)\n",
    "chiProb    = np.empty(nev)\n",
    "x1d   = xdet.copy()  # not sure if the copy is necessary, but it can't hurt\n",
    "x2d   = xdet.copy()\n",
    "s1d   = s.copy()\n",
    "s2d   = s.copy()\n",
    "igood = 0\n",
    "for i in range(nev):  # loop over pairs\n",
    "    y1 = trk1[i][:]   # coordinates of trk 1 for this pair\n",
    "    y2 = trk2[i][:]   # coordinates of trk 2 for this pair\n",
    "\n",
    "    # guess parameters\n",
    "    slope1 = ( y1[3] - y1[0] ) / ( x1d[3] - x1d[0] )\n",
    "    slope2 = ( y2[3] - y2[0] ) / ( x2d[3] - x2d[0] )\n",
    "    inter1 = y1[3] - slope1 * x1d[3]\n",
    "    inter2 = y2[3] - slope2 * x2d[3]\n",
    "    guess = np.array( [slope1, inter1, slope2, inter2] )\n",
    "\n",
    "    # Now minuit\n",
    "    f = myFit(x1d,x2d,y1,y2,s1d,s2d,sigma=1e-4)\n",
    "    m = Minuit(f.CHI2, guess, name=(\"p0\",\"p1\",\"p2\",\"p3\"))\n",
    "    m.errordef = 1.\n",
    "    m.migrad()\n",
    "    \n",
    "    if m.valid:    # fit worked\n",
    "    \n",
    "        # The fitted vertex in microns and the pull \n",
    "        thisFitVertex1     = 10000*m.values['p1']\n",
    "        thisPull1          = (m.values['p1']-xv[i])/m.errors['p1']\n",
    "        fitVertex1[igood]  = thisFitVertex1\n",
    "        pull1[igood]       = thisPull1\n",
    "        thisFitVertex2     = 10000*m.values['p3']\n",
    "        thisPull2          = (m.values['p3']-xv[i])/m.errors['p3']\n",
    "        fitVertex2[igood]  = thisFitVertex2\n",
    "        pull2[igood]       = thisPull2\n",
    "    \n",
    "        # chiprob is like TMath::Prob in ROOT.\n",
    "        # Probability that an observed Chi-squared exceeds\n",
    "        # the value chisq by chance, even for a correct model\n",
    "        ndof  = len(y1) + len(y2) - (len(guess)-1)\n",
    "        chisq = m.fval\n",
    "        thisChiProb     = 1. - stats.chi2.cdf(chisq, ndof)\n",
    "        chiProb[igood]  = thisChiProb\n",
    "    \n",
    "        igood = igood+1\n",
    "        \n",
    "    else:\n",
    "        \n",
    "        print(\"Fit on event \",i,\" failed\")\n",
    "    \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "2135b549",
   "metadata": {},
   "outputs": [],
   "source": [
    "# If some fits failed, take them out from the list\n",
    "if igood != nev:\n",
    "    fitVertex1 = fitVertex1[0:igood]\n",
    "    fitVertex2 = fitVertex2[0:igood]\n",
    "    pull1      = pull1[0:igood]\n",
    "    pull2      = pull2[0:igood]\n",
    "    chiprob    = chiprob[0:igood]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "826de1fb",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Sanity check: They need to be about the same\n",
    "np.testing.assert_allclose(fitVertex1, fitVertex2, rtol=0.01, atol=1)\n",
    "np.testing.assert_allclose(pull1,          pull2, rtol=1, atol=0.01)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "ab1abcd8",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEGCAYAAABrQF4qAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAmJ0lEQVR4nO3deXSV1f3v8feXSZmnBophCBpEAoQpAraSMijyQxaDIsOiFZXfpaulrSPCT6+1VK2g/Jy9ernFKygySFW4LSIyWBAHipikzCAkClJBJgUEJdn3j+fJ6SFkOEnOSU4ePq+1zsp59jPtDck3O/vs/X3MOYeIiARXtcqugIiIxJYCvYhIwCnQi4gEnAK9iEjAKdCLiARcjYq82Y9+9COXlJRUkbcUEanyPvnkk6+dcwllPb9CA31SUhIbN26syFuKSAX78ssv+eCDDzh27Bh5eXmVXZ24UrNmTRISEhgwYAC1a9eO+DwzyynPfSs00ItIsO3fv5/XXnuN9PR0evXqRbVqGh0O9/3337N582ZeeeUVfvGLX5Qq2JeHAr2IRM1bb73F4MGD6dixY2VXJW4lJSWxdOlS1q1bx8CBAyvknvp1KyJR4ZzjyJEjtG/fvrKrEtfMjPbt2/P1119X2D0V6EUkKpxzOOeoUUMDBSWpVasWZ8+erbD7KdCLSJVmZtx9992h7ZkzZ/KHP/wh4vMHDRpEo0aNGDJkyDnle/fupVevXiQnJzN69Gi+//57AM6cOcPo0aNJTk6mV69eZGdnh8559NFHSU5Opn379rzzzjvlalc0KdCLSJV20UUX8cYbb5R5KGTy5Mm88sor55VPmTKFO++8k927d9O4cWNmz54NwOzZs2ncuDG7d+/mzjvvZMqUKQBs3bqVBQsWsGXLFpYvX86vf/1rcnNzy96wKFKgF5EqrUaNGkycOJEnn3yyTOcPGDCA+vXrn1PmnGP16tWMHDkSgPHjx/PWW28BsGTJEsaPHw/AyJEjWbVqFc45lixZwpgxY7joooto27YtycnJbNiwoewNiyINpolIlTdp0iRSU1O59957zymfN28ejz/++HnHJycns3jx4iKvd/jwYRo1ahT6vKFly5bs378f8KaQtmrVCvB+yTRs2JDDhw+zf/9+evfuHbpG+DmVTYFeRKq8Bg0acPPNN/PMM8+cMzd93LhxjBs3rhJrFh80dCOCN7fZzAp9KW1H1XDHHXcwe/ZsTp48GSqbN28eXbt2Pe+VPyRTlKZNm3Ls2LHQzJh9+/aRmJgIQGJiIl988QUAZ8+e5fjx4zRt2vSc8oLnVDYFehEgJycnND2w4Csnp1yrz6WCNGnShFGjRoU+NAWvR5+RkXHeq7hhG/Bm8vTr1y903Jw5cxg2bBgAQ4cOZc6cOQAsXryY/v37Y2YMHTqUBQsWcObMGfbu3cuuXbvo2bNnjFpbOgr0IhIYd999d6ln3/Tp04ebbrqJVatW0bJly9C0yBkzZvDEE0+QnJzM4cOHmTBhAgATJkzg8OHDJCcn88QTTzB9+nQAOnbsyKhRo0hJSWHQoEE8//zzVK9ePboNLKOIxujNrBHwZ6AT4IDbgB3AQiAJyAZGOeeOxqKSIiJFOXHiROh98+bNOXXqVKnOX7duXaHll156aaGzZi6++GJef/31Qs+5//77uf/++0t1/4oQaY/+aWC5c+4KoAuwDZgKrHLOtQNW+dsiIhJnSgz0ZtYQSAdmAzjnvnfOHQOGAXP8w+YAw2NTRRERKY9Ihm7aAoeA/2tmXYBPgNuB5s65A/4x/wKaF3aymU0EJgK0bt263BUWkWDYt28fkyZNYuvWreTl5TFkyBAef/xxatWqFZXrT548mWXLljF48GDq1q1LvXr1uOeee6Jy7aJUr16dzp07A168W7p0aUzvF6lIhm5qAN2BF5xz3YCTFBimcc45vLH78zjnZjnn0pxzaQkJZX5AiogEiHOOG264geHDh7Nr1y527tzJiRMnojK+nT8lctasWWRlZRW6YCpWateuHZrZEy9BHiIL9PuAfc65j/3txXiB/yszawHgfz0YmyqKSNCsXr2aiy++mFtvvRXwesJPPvkkL730EqdOnaJ3795s2bIldHzfvn3ZuHEjJ0+e5LbbbqNnz55069aNJUuWAPDyyy8zdOhQ+vfvz4ABAxg6dCgnTpygR48eLFy48Jx7Z2Rk0Lt3b1JTUxkxYgRHjx7l4MGD9OjRA4DMzEzMjM8//xyAyy67rNQf8MabEgO9c+5fwBdmlp9kegCwFVgKjPfLxgNLYlJDEQmcLVu2hAJrvgYNGtC6dWt2797N6NGjWbRoEQAHDhzgwIEDpKWl8cgjj9C/f382bNjAmjVrmDx5cmiB1KZNm1i8eDF///vfWbp0aah3PXr06HPuc/PNNzNjxgyysrLo3Lkz06ZNo1mzZpw+fZpvvvmGdevWkZaWxrp168jJyaFZs2bUqVMnosVXp0+fJi0tjd69e4dy48SDSFMg/BaYZ2a1gD3ArXi/JBaZ2QQgBxgVmyqKyIVm1KhRDBw4kGnTprFo0aJQMF2xYgVLly5l5syZgBdY83ve1157LU2aNCn2usePH+fYsWP87Gc/A7xkZTfddBMAP/nJT1i/fj1r167lvvvuY/ny5Tjn6NOnDxBZOoWcnBwSExPZs2cP/fv3p3Pnzlx22WVl/4eIkogCvXMuA0grZNeAqNZGRC4IKSkp561O/eabb/j8889JTk6mTp06NG3alKysLBYuXMiLL74IeGP7f/nLX857itXHH39M3bp1y1Wn9PT0UC9+2LBhzJgxAzPj+uuvByJLkJaf8uDSSy+lb9++fPrpp3ER6LUyVkQq3IABAzh16hRz584FIDc3l7vvvptbbrmFOnXqADB69Ggee+wxjh8/TmpqKgDXXXcdzz77LN78D/j0009Ldd+GDRvSuHHj0CKpV155JdS779OnD6+++irt2rWjWrVqNGnShGXLlnH11VcDJadTOHr0KGfOnAHg66+/Zv369aSkpJTnnylqFOhFpMKZGW+++Savv/467dq14/LLL+fiiy/mT3/6U+iYkSNHsmDBAkaN+veo8AMPPMAPP/xAamoqHTt25IEHHij1vefMmcPkyZNJTU0lIyOD3//+94CX2M45R3p6OgBXX301jRo1onHjxhFdd9u2baSlpdGlSxf69evH1KlT4ybQW/5vxoqQlpbmNm7cWGH3E4mUmVHUz0Jx++Tf8vLyePjhh0OBU4q2Z88e3n//fW6++eaIjjezT5xzhQ2fR0Q9ehGRgFOgFxEJOAV6EZGAU6AXkZjIzs6mU6dO55T94Q9/CM2Bj1T+qtjS2L59O127dqVbt2589tln1KtXr1Tnl8UDDzxAamoqXbt2ZeDAgXz55Zcxv2ekFOhFJFByc3N56623GDlyZIXOY588eTJZWVlkZGQwZMgQ/vjHP1bIfSOhQC8ilaJv375MmTKFnj17cvnll4fmtn/33XeMGTOGDh06MGLECL777rvQOStWrOCqq66ie/fu3HTTTaGHjiQlJTFlyhS6d+/OwoULeeqpp3jhhRfo16/fOfd0zjF58mQ6depE586dQ3lwJk2aFEpCNmLECG677TYAXnrppYgTrTVo0CD0/uTJk5hZGf9loi/SFAgiIlF39uxZNmzYwLJly5g2bRorV67khRdeoE6dOmzbto2srCy6d+8OeIuQHn74YVauXEndunVDj/rLn87ZtGlTNm3aBMDOnTsLTUv8xhtvkJGRQWZmJl9//TVXXnkl6enp9OnTh3Xr1jF06FD279/PgQNeBvZ169YxZswYwFtQ9e23357XhpkzZ3LNNdcA3hOm5s6dS8OGDVmzZk1s/tHKQIFeRGKiqB5tePkNN9wAQI8ePcjOzgZg7dq1/O53vwMgNTU1tCr2o48+YuvWrfz0pz8F4Pvvv+eqq64KXatg8rLCvP/++4wdO5bq1avTvHlzfvazn/GPf/yDPn368NRTT7F161ZSUlI4evQoBw4c4MMPP+SZZ54Bin7kYLhHHnmERx55hEcffZTnnnuOadOmlXhORVCgF5GYaNq0KUePnvsY6SNHjtC2bdvQ9kUXXQR4aYrz88gXxTnHtddey/z58wvdX55cN4mJiRw7dozly5eTnp7OkSNHWLRoEfXq1aN+/fpAZD36fOPGjWPw4MFxE+g1Ri8iMVGvXj1atGjB6tWrAS/IL1++PJQ7pijp6em89tprAGzevJmsrCwAevfuzfr169m9ezfgjYPv3LmzVHXq06cPCxcuJDc3l0OHDrF27Vp69uwZuv5TTz0VGsqZOXNmKHMleD36wnLd5Af5Xbt2hY5dsmQJV1xxRanqFkvq0YtIzMydO5dJkyZx1113AfDggw+WOAvmV7/6FbfeeisdOnSgQ4cOobz1CQkJvPzyy4wdOzaUPOzhhx/m8ssvj7g+I0aM4MMPP6RLly6YGY899hg//vGPAe+XwIoVK0hOTqZNmzYcOXLknEBfkqlTp7Jjxw6qVatGmzZtQhk344Fy3cgFJSkpiZycnPPK27RpExojLki5biKjXDeRq+hcN+rRywUlJydHQVsuOBqjFxEJOAV6EZGAU6AXEQk4BXoRkYBToBcRCXPbbbfRrFmz8zJvFvTkk0/SsWNHOnXqxNixYzl9+nRoX25uLt26dWPIkCGxrm5EFOhFRMLccsstLF++vNhj9u/fzzPPPMPGjRvZvHkzubm5LFiwILT/6aefpkOHDrGuasQU6EUkEPr27cv27dsBOHz4cIk98qKkp6fTpEmTEo87e/Ys3333HWfPnuXUqVNccsklAOzbt4+//e1v/Od//meZ7h8LmkcvIoGwe/fu0CrZrKwsOnfufM7+0uSqKUliYiL33HMPrVu3pnbt2gwcOJCBAwcCcMcdd/DYY48Veq/KokAvIlVeTk4OiYmJVKvmDVJkZWWFsl7miyT7ZKSOHj3KkiVL2Lt3L40aNeKmm27i1VdfpVGjRjRr1owePXrw3nvvRe1+5RVRoDezbOBbIBc465xLM7MmwEIgCcgGRjnnjhZ1DRGRWMnMzDwnsH/yySfnpS2OZo9+5cqVtG3bloSEBMBLt/zBBx/QsGFDli5dyrJlyzh9+jTffPMNP//5z3n11VfL0KroKU2Pvp9z7uuw7anAKufcdDOb6m9PiWrtREQikJGREZr1smvXLpYsWcLDDz98zjHR7NG3bt2ajz76iFOnTlG7dm1WrVpFWloav/3tb3n00UcBeO+995g5c2alB3ko34exw4A5/vs5wPBy10ZEpAwyMzPJy8ujS5cu/PGPfyQlJYU5c+aUfGIhxo4dy1VXXcWOHTto2bIls2fPDu0bPHgwX375Jb169WLkyJF0796dzp07k5eXx8SJE6PVnKiLKHulme0FjgIO+N/OuVlmdsw518jfb8DR/O0C504EJgK0bt26R2GZA0UqSlkyUUY7e2VRGTSh+Cya8a4ys1e2a9eOTZs2hR4SEu/iNXvl1c65/WbWDHjXzLaH73TOOTMr9CfBOTcLmAVemuKyVlQkKIrLoBlPD5SuKr799lvMrMoE+coQ0dCNc26///Ug8CbQE/jKzFoA+F8PxqqSIiJFqV+/fqmfNHWhKTHQm1ldM6uf/x4YCGwGlgLj/cPGA0tiVUkJtqSkJMzsvFdSUlJlV00kECIZumkOvOn/SVkDeM05t9zM/gEsMrMJQA4wKnbVlCAraihDwxgi0VFioHfO7QG6FFJ+GBgQi0qJiEj0KNeNiEjAKdCLiAScAr2IBF5ROeaXL19O+/btSU5OZvr06efsS0pKonPnznTt2pW0tPOnsH/xxRf069ePlJQUOnbsyNNPPx3ad+zYMUaOHMkVV1xBhw4d+PDDD2PTsAgp0ItI4BWWYz43N5dJkybx9ttvs3XrVubPn8/WrVvPOWbNmjVkZGSwcePG865Zo0YN/vu//5utW7fy0Ucf8fzzz4fOv/322xk0aBDbt28nMzOz0nPTK9CLSNzKzMwkPT2dlJQUqlWrhpmVaeVtYTnmN2zYQHJyMpdeeim1atVizJgxLFkS+SzxFi1a0L17d8Cby9+hQwf279/P8ePHWbt2LRMmTACgVq1aNGrUqNR1jialKRaRuHT69GlGjx7N3Llz6dmzJw888ACnT59m2rRpoWPKk5Fy//79tGrVKrTdsmVLPv7449C2mTFw4EDMjF/+8pfF5rLJzs7m008/pVevXuzZs4eEhARuvfVWMjMz6dGjB08//TR169YtTfOjSoFeROLSypUr6d69Oz179gQgNTWV5cuXn7O+IpoZKQt6//33SUxM5ODBg1x77bVcccUVpKenn3fciRMnuPHGG3nqqado0KABZ8+eZdOmTTz77LP06tWL22+/nenTp/PQQw/FrK4l0dCNiMSlzZs3n/OUqE2bNoWGSvL16dOHrl27nvdauXJliddPTEzkiy++CG3v27ePxMTEc/YDNGvWjBEjRrBhw4bzrvHDDz9w4403Mm7cOG644QbA+8ugZcuW9OrVC4CRI0eyadOmUrQ8+tSjF5G41LRpU1avXg3Azp07eeONN/jggw/OOaY8Pforr7ySXbt2sXfvXhITE1mwYAGvvfYaACdPniQvL4/69etz8uRJVqxYcd5nA845JkyYQIcOHbjrrrtC5T/+8Y9p1aoVO3bsoH379qxatYqUlJQy1zMa1KMXkbg0duxYTpw4QadOnZg4cSLz58+nadOmZb5WwRzzNWrU4LnnnuO6666jQ4cOjBo1io4dOwLw1VdfcfXVV9OlSxd69uzJ9ddfz6BBg4B/56Rfv349r7zyCqtXrw79JbFs2TIAnn32WcaNG0dqaioZGRncd9990flHKaOI8tFHS1pamitsmpJc2IrK9x7tPPBlvWa061Hc9WLR5oqSl5fHQw89xIMPPljZVYl7n332GevXr6+wfPTq0YtIVJgZtWrV4tSpU5Vdlbh38uRJ6tSpU2H3U6AXkagwM5KTk3n77bfJy8ur7OrErW+//Za1a9eSnJxcYffUh7EiEjXDhw9nwYIF/PnPf6Z58+ZUr169sqsUV77//nuys7O58sor6dq1a4XdV4FepByKev5rVX72a3nUrFmTsWPHkpOTw7Fjx9SzL6BmzZr07t2bSy65pELvq0AvUg56aMr5atSowWWXXVbZ1ZAwGqMXEQk4BXoRkYBToBcRCTgFehGRgFOgFxEJOAV6EZGAU6AXEQk4BXoRkYCLONCbWXUz+9TM/upvtzWzj81st5ktNLNasaumyIWhTZs2mNl5r6SkpMqumlRhpenR3w5sC9ueATzpnEsGjgITolkxkQtRdnY2zrnzXoWlWRCJVESB3sxaAtcDf/a3DegPLPYPmQMMj0H9RESknCLt0T8F3AvkZyhqChxzzp31t/cBiYWch5lNNLONZrbx0KFD5amriIiUQYmB3syGAAedc5+U5QbOuVnOuTTnXFpCQkJZLiEiIuUQSfbKnwJDzWwwcDHQAHgaaGRmNfxefUtgf+yqKSIiZVVij94591/OuZbOuSRgDLDaOTcOWAOM9A8bDyyJWS1FRKTMyjOPfgpwl5ntxhuznx2dKkkQJSUlFTpt0Mxo06ZNZVevWEVNeazIuhdXB029lJJYRT5xPi0tzW3cuLHC7ifxw8wKfUBHtM+pjGuW9j7RrkNFtUkqj5l94pxLK+v5WhkrIhJwCvQSOFV5mEgkFvTMWAmcop7jKnKhUo9eRCTgFOilStLwjEjkNHQjVZKGZ0Qipx69iEjAKdCLiAScAr2ISMAp0IuIBJwCvYhIwGnWjcSt/EReRe0Tkcgo0Evcys7OruwqiASChm5ERAJOgV5EJOAU6EVEAk6BXkQk4BToRUQCTrNuRGJAU0MlnijQi8SApoZKPNHQjYhIwCnQS6kV99APPQxEJP5o6EZKTQ/9EKla1KMXEQm4EgO9mV1sZhvMLNPMtpjZNL+8rZl9bGa7zWyhmdWKfXVFpDSKGmZLSkqq7KpJBYqkR38G6O+c6wJ0BQaZWW9gBvCkcy4ZOApMiFktRaRM8ofZCr5ycnIqu2pSgUoM9M5zwt+s6b8c0B9Y7JfPAYbHooIiIlI+EY3Rm1l1M8sADgLvAp8Bx5xzZ/1D9gGJMamhiIiUS0SB3jmX65zrCrQEegJXRHoDM5toZhvNbOOhQ4fKVksRESmzUs26cc4dA9YAVwGNzCx/emZLYH8R58xyzqU559ISEhLKU1cRESmDSGbdJJhZI/99beBaYBtewB/pHzYeWBKjOoqISDlEsmCqBTDHzKrj/WJY5Jz7q5ltBRaY2cPAp8DsGNZTRETKqMRA75zLAroVUr4Hb7xeRETimFbGiogEnAK9iJyjuKR1WlFbNSmpmYico7ikdUU9TEXim3r0IiIBp0AvIhJwCvRSqOLGafUQEZGqRWP0Uig9XEQkONSjFxEJOAV6EZGAU6AXEQk4BXoRkYBToL+AaWaNyIVBs24uYJpZI3JhUI9eRCTgFOhFRAJOQzciVVybNm2KTDZW1GctZTlHqi4FepEqLjs7u0LOkapLQzciIgGnQC8iEnAK9CIiAadALyIScAr0IiIBp0AvIhJwCvQiIgGnQC8iEnAlBnoza2Vma8xsq5ltMbPb/fImZvaume3yvzaOfXWlKMVlokxKSqrs6klA5K+o1fdZ1RLJytizwN3OuU1mVh/4xMzeBW4BVjnnppvZVGAqMCV2VZXiFJeJsqil7iKlVdyKWn2fxa8Se/TOuQPOuU3++2+BbUAiMAyY4x82BxgeozqKiEg5lGqM3sySgG7Ax0Bz59wBf9e/gOZFnDPRzDaa2cZDhw6Vp64iIlIGEQd6M6sH/AW4wzn3Tfg+540ZFDpu4Jyb5ZxLc86lJSQklKuyIiJSehEFejOriRfk5znn3vCLvzKzFv7+FsDB2FRRRETKI5JZNwbMBrY5554I27UUGO+/Hw8siX71RESkvCKZdfNT4BfAP80swy+7D5gOLDKzCUAOMComNRQRkXIpMdA7594Hipo3NSC61RERkWjTylgRkYBToBcRCTgFehGRgFOgFxEJOAV6EZGAU6AXEQk4BXoRkYBToBeRqCgqV73y1Fe+SFbGioiUqKhc9cpTX/nUoxcRCTgFehGRgFOgFxEJOAV6EZGAU6AXEQk4BXoRkYBToBcRCTgFehGRgFOgjzNJSUmFri4s6dWmTZsir1nUisXizhGpbEX9LGilbelpZWycycnJwTkX1WsWtWJRJJ4V9bOglbalpx69iEjAKdCLiAScAr2ISMAp0IuIBJwCvYhIwJUY6M3sJTM7aGabw8qamNm7ZrbL/9o4ttUUkaqqqOm9mipZcSLp0b8MDCpQNhVY5ZxrB6zyt0VEzpOdnY1zrtBXTk5OZVfvglBioHfOrQWOFCgeBszx388Bhke3WiIiEi1lHaNv7pw74L//F9C8qAPNbKKZbTSzjYcOHSrj7YKnqFV/Wq0qItFW7pWxzjlnZkUu5XTOzQJmAaSlpUV3yWcVFosVsCIihSlrj/4rM2sB4H89GL0qiYhINJU10C8FxvvvxwNLolMdERGJtkimV84HPgTam9k+M5sATAeuNbNdwDX+tohIqRQ39VKfV0VPiWP0zrmxRewaEOW6iMgFRplVK4ZWxoqIBJwCfTnF4kEhIiLRpAePlJOmSYpIvFOPXkQk4BToRUQCToFeRCTgFOhFRAJOgV5EJOACG+iLmvZY3IMOipsqqQckiMQHPcik9AI7vbKoaY9mVupzSjpPRCpOcatp9XNauMD26EVExHPBBfqyJlEq6jytcBWJH2UZ1rkQhmwDO3RTlLImUVLyJZH4V5ZhnQthyPaC69GLiFxoFOhFRAJOgV5EJOAU6EVEAk6BXkQk4OIm0JdlJauISKTKMkU6KKtw42Z6ZVlWsoqIRKosU6SDsgo3bnr0IiISG1U60Be3ok0rVkWkKikunpVX3AzdlIWe1yoiQRHLFbpVukcvIiIlU6AXESmDaCdQi+Vwc7mGbsxsEPA0UB34s3NuelRqJSIS56KdQC2WytyjN7PqwPPAfwApwFgzS4lWxUREJDrKM3TTE9jtnNvjnPseWAAMi061REQkWsozdJMIfBG2vQ/oVfAgM5sITPQ3z5jZ5qIuWNSfO8V94hxnixZ+BHxd2ZWIkSC3DdS+qi7u2leWeFaM9uWpS8ynVzrnZgGzAMxso3MuLdb3rCxBbl+Q2wZqX1V3IbSvPOeXZ+hmP9AqbLulXyYiInGkPIH+H0A7M2trZrWAMcDS6FRLRESipcxDN865s2b2G+AdvOmVLznntpRw2qyy3q+KCHL7gtw2UPuqOrWvGKYUAiIiwaaVsSIiAadALyIScDEJ9Gb2kJllmVmGma0ws0v8cjOzZ8xst7+/e9g5481sl/8aH4t6RYuZPW5m2/02vGlmjcL2/Zffvh1mdl1Y+SC/bLeZTa2UikfIzG4ysy1mlmdmaQX2Vfn2FVSV657PzF4ys4Ph61TMrImZvev/TL1rZo398iJ/DuORmbUyszVmttX/vrzdLw9K+y42sw1mlum3b5pf3tbMPvbbsdCf9IKZXeRv7/b3J5V4E+dc1F9Ag7D3vwNe9N8PBt4GDOgNfOyXNwH2+F8b++8bx6JuUWrfQKCG/34GMMN/nwJkAhcBbYHP8D6oru6/vxSo5R+TUtntKKZ9HfAWaLwHpIWVB6J9BdpaZeteoB3pQHdgc1jZY8BU//3UsO/TQn8O4/UFtAC6++/rAzv978WgtM+Aev77msDHfr0XAWP88heBX/nvfx0WU8cAC0u6R0x69M65b8I26wL5n/gOA+Y6z0dAIzNrAVwHvOucO+KcOwq8CwyKRd2iwTm3wjl31t/8CG8NAXjtW+CcO+Oc2wvsxksVUaXSRTjntjnndhSyKxDtK6Aq1z3EObcWOFKgeBgwx38/BxgeVl7Yz2Fccs4dcM5t8t9/C2zDW5kflPY559wJf7Om/3JAf2CxX16wffntXgwMsBKW28ZsjN7MHjGzL4BxwO/94sLSJiQWU14V3IbXe4Bgti9cENtXletekubOuQP++38Bzf33VbbN/jBFN7xeb2DaZ2bVzSwDOIjX0f0MOBbWoQxvQ6h9/v7jQNPirl+e7JUrzWxzIa9hfgXud861AuYBvynrfSpLSe3zj7kfOIvXxiolkvZJcDjv7/wqPZfazOoBfwHuKDBqUOXb55zLdc51xRsd6AlcEc3rl2fB1DURHjoPWAY8SNFpE/YDfQuUv1fWukVDSe0zs1uAIcAA/5sMik8LEVfpIkrx/xeuyrSvFIKcyuMrM2vhnDvgD10c9MurXJvNrCZekJ/nnHvDLw5M+/I5546Z2RrgKrwhpxp+rz28Dfnt22dmNYCGwOHirhurWTftwjaHAdv990uBm/1PxXsDx/0/vd4BBppZY/+T84F+WVwy74Er9wJDnXOnwnYtBcb4n4q3BdoBGwhOuoggtq8q170kS4H8GWzjgSVh5YX9HMYlf/x5NrDNOfdE2K6gtC/B/Jl7ZlYbuBbvc4g1wEj/sILty2/3SGB1WGezcDH6FPkvwGYgC/h/QGLYp8vP440//ZNzZ3Tchvfh3m7g1ljUK4rt2403Rpbhv14M23e/374dwH+ElQ/Gmy3wGXB/ZbehhPaNwBsTPAN8BbwTpPYV0t4qW/ewNswHDgA/+P93E/DGbVcBu4CVQBP/2CJ/DuPxBVyNNyyTFfYzNzhA7UsFPvXbtxn4vV9+KV5HajfwOnCRX36xv73b339pSfdQCgQRkYDTylgRkYBToBcRCTgFehGRgFOgFxEJOAV6EZGAU6CXqDOz35nZNjObZ2ZDzc8IaWbDzSwl7LhbzM9sWoprJ1lYhsaw8jQ/819+hr/LzGyPmTXwtwvNulngGmZmq/PPibA+ofZVJDP7jZndVtH3lapJ0ysl6sxsO3CNc25fgfKXgb865xb72+8B9zjnIn7CvZ/r5K/OuU6F7PtfwD7n3J/MbDkwxzk33//lMh9vafkleHOuL3fO5RY4/3q/3ndG3NhS8Bf+mHMuLwrXqgOsd851K3/NJOjUo5eoMrMX8RZ6vG1md/q99ufM7CfAUOBx855TMAVIA+b527XNrIeZ/d3MPjGzd/xl7fjlmWaWCUwq5vb3Af/DzO7FSyM93y8vKutmQePwVx/6fzlsN7OXzWyn/9fJNWa23rz85z39424xs+f8983Nez5Bpv/6iX+dHWY2F28xTCvznmew2cz+aWaj/XP7mtl7ZrbYv+88/xcDZjbdvFzsWWY2E8B5K7Kz8+shUqzKXhWmV/BeQDbwI//9LcBz/vuXgZFhx72Hv2oRLzXrB0CCvz0a74Hz4K0YTPffP05YzvVC7v1LIBdoH1b2HPDzsO3Z4fUIK88B6vvvk/AS1nXG6xB9AryEt+pyGPBWIe1biJdwC7w89w396+QBvf3yG/GyE1bHy7b4OV6+9b54WQhb+vf7EG9FaFO8Vcj5f303Cqvv/cDdlf3/rVf8v9Sjl3jRHugEvGteutb/CbT0c4A0cl6+dYBXSrjOf+ClbUgp4bjCNHFevvN8e51z/3TeUMsWYJVzzuEtq08q5Pz+wAsQykZ43C/PcV5edPCC93x//1fA34Er/X0bnHP7/Ptl+Pc4DpwGZpvZDUB4bqWDeENRIsUqc/ZKkSgzYItz7qpzCsMe01jiBcyG4PWirwPeNLN3nDfEEWk2w7NmVs39ewz9TNi+vLDtPEr3s3MywuPC75eLN/x01h+eGYCXwOo3eL9QwMt58l0p6iEXKPXopSJ9i/couMK2dwAJZnYVeGlpzayjc+4YcMzMrvaPG1fYhf2sf08Ak5xz/8Qba7/f311U1s2CduB9vlBWq4Bf+fWpbmYNCzlmHTDa35+A9wjAwuqCf516QEPn3DLgTqBL2O7L8cb9RYqlQC8VaQEw2cw+NbPL8MbsX/SHaqrj9Vhn+B+6ZgA/8c+7FXjeP66oR6Y9ALzpnNvqb/8BGGtm7ZxzW/Cev7kVWI73yyC3kGv8jXOfi1BatwP9zOyfeGP6hQ0fvYn3mUMmsBq41zn3r2KuWR/4q5llAe8Dd4Xt+yneeL9IsTS9UsTnz/KZ65y7trLrUhIz6wbc5Zz7RWXXReKfevQiPuc9nOL/lGbBVCX6Ed5fMSIlUo9eRCTg1KMXEQk4BXoRkYBToBcRCTgFehGRgFOgFxEJuP8PQmKYHaGYCpAAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the fitted value\n",
    "ax33 = plt.subplot(111)\n",
    "con33, bins33, _ = ax33.hist(fitVertex1, np.linspace(-300.,300.,51), histtype='step', color='black')\n",
    "statBox(ax33, fitVertex1, bins33)\n",
    "ax33.set_xlim(bins33[0], bins33[-1])\n",
    "_ = ax33.set_xlabel(\"fitted X0 (microns)\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "2e2fda18",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEJCAYAAACdePCvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAnW0lEQVR4nO3deXRV9bn/8ffDIILIKFJMIEGClsGIMURohTKoKLpA/CHD9VYU/bFuS29br0VsXaJ0VOtC622rpeIVWwsoReH2p5FJL1oVihhSwIG0JAqXVkaVSSE8vz/2zvEkZD4ZTrI/r7XOyt7fPX2/Sc5+9ncPzzZ3R0REoqdFY1dAREQahwKAiEhEKQCIiESUAoCISEQpAIiIRJQCgIhIRFUZAMzsCTP7yMy2xJX93MzeNbN8M3vOzDrFTfu+mRWY2XtmNqae6i0iIgmqTg/gSeDKMmWrgIHungm8D3wfwMz6A1OAAeEyvzazlnVWWxERqTOtqprB3deZWXqZspVxo28CE8Ph8cBid/8M2GFmBUAO8EZl2zjrrLM8PT29sllERKSMt956a6+7d6vt8lUGgGqYDiwJh1MIAkKJnWFZpdLT09m4cWMdVEVEktWRI0dYs2YNe/fu5fjx441dnaRiZnTs2JGcnBxqcjBsZkWJbDehAGBmdwEngKdrsewMYAZAr169EqmGiCS5I0eOsHDhQtLS0hg5ciStW7du7CollZMnT7J3716effZZrrvuOvr06dMg2611ADCzm4BrgNH+RUKhXUDPuNlSw7JTuPt8YD5Adna2EhKJNGMvv/wyvXr14qqrrsLMGrs6Salnz5506dKFJUuWMGvWrAb5PdXqNlAzuxK4Axjn7kfiJq0ApphZGzPrDfQFNiReTRFpyvbu3cuXv/xl7fyrkJaWxsmTJzl69GiDbK86t4EuIriIe76Z7TSzW4BfAmcCq8wsz8weA3D3rcAzwDYgF5jp7sX1VnsRaRKKi4sb/LSPmXH77bfHxh988EHuvffeai9/5ZVX0qlTJ6655ppS5Tt27OCSSy4hIyODyZMn8/nnnwPw2WefMXnyZDIyMrjkkksoLCyMLfOzn/2MjIwMzj//fF566aVKt3vaaadx4sSJatczEVUGAHef6u493L21u6e6+wJ3z3D3nu4+KPz8W9z8P3H3Pu5+vru/WL/VFxEpX5s2bVi2bBl79+6t1fKzZs3id7/73Snls2fP5rbbbqOgoIDOnTuzYMECABYsWEDnzp0pKCjgtttuY/bs2QBs27aNxYsXs3XrVnJzc/nmN79JcXFyHBfrSWARaZZatWrFjBkzeOihh2q1/OjRoznzzDNLlbk7a9euZeLE4M73adOm8fzzzwOwfPlypk2bBsDEiRNZs2YN7s7y5cuZMmUKbdq0oXfv3mRkZLBhQ3KcGa+L20BFRJLSzJkzyczM5I477ihV/vTTT/Pzn//8lPkzMjJYunRphevbt28fnTp1olWrYNeZmprKrl3BfS67du2iZ8/gHphWrVrRsWNH9u3bx65duxgyZEhsHfHLNDYFABFptjp06MCNN97II488Qtu2bWPlN9xwAzfccEMj1iw5KACISLP23e9+l6ysLG6++eZYWW17AF27duXgwYOcOHGCVq1asXPnTlJSgmddU1JS+PDDD0lNTeXEiRN8/PHHdO3aNVZeIn6ZxqZrACJlpKenY2axj9KUNG1dunRh0qRJsYu1EPQA8vLyTvlUtvOH4M6ikSNHxuZbuHAh48ePB2DcuHEsXLgQgKVLlzJq1CjMjHHjxrF48WI+++wzduzYwfbt28nJyamn1taMAoBIGUVFRbh77FNUlNDT9pIEbr/99hrfDTRs2DCuv/561qxZQ2pqauz2zfvvv5958+aRkZHBvn37uOWWWwC45ZZb2LdvHxkZGcybN4/77rsPgAEDBjBp0iT69+/PlVdeya9+9StatkyOHJk6BSQizdKhQ4diw927d+fIkSOVzH2qV199tdzyc889t9y7eE4//XSeffbZcpe56667uOuuu2q0/YagHoCISEQpAIiIRJQCgIhIRCkAiEhS2blzJ+PHj6dv37706dOH73znO7F8O3Vh1qxZDBgwgFmzZnHvvffy4IMP1tm6K7Jw4UL69u1L3759Y3cKJQMFABFJGu7Oddddx7XXXsv27dt5//33OXToUJ1cQC1JsDZ//nzy8/PLfQ6gPuzfv5+5c+eyfv16NmzYwNy5czlw4ECDbLsqCgAikjTWrl3L6aefHntoq2XLljz00EM88cQTHDlyhCFDhrB169bY/CNGjGDjxo0cPnyY6dOnk5OTw0UXXcTy5csBePLJJxk3bhyjRo1i9OjRjBs3jkOHDnHxxRezZMmSUtvOy8tjyJAhZGZmMmHCBA4cOMBHH33ExRdfDMDmzZsxMz744AMA+vTpU607i1566SUuv/xyunTpQufOnbn88svJzc2tk99XohQARCRpbN26NbbDLdGhQwd69epFQUEBkydP5plnngFg9+7d7N69m+zsbH7yk58watQoNmzYwMsvv8ysWbM4fPgwAJs2bWLp0qX8z//8DytWrKBt27bk5eUxefLkUtu58cYbuf/++8nPz+eCCy5g7ty5nH322Rw7doxPPvmEV199lezsbF599VWKioo4++yzadeuHU8//TSDBg065VOSMC4+RxAoF5CISK1MmjSJK664grlz5/LMM8/EdrIrV65kxYoVsfP5x44dix2plxx9V+bjjz/m4MGDfO1rXwOCLJ/XX389AF/5ylf485//zLp16/jBD35Abm4u7s6wYcOApp1XSAFARJJG//79T0nH8Mknn/DBBx+QkZFBu3bt6Nq1K/n5+SxZsoTHHnsMCK4d/PGPf+T8888vtez69es544wzEqrT8OHDY0f948eP5/7778fMuPrqq4Gq8wqlpKTwyiuvxMp37tzJiBEjEqpTXdEpIBFJGqNHj+bIkSM89dRTQPAmsdtvv52bbrqJdu3aATB58mQeeOABPv74YzIzMwEYM2YM//mf/0nJ68nffvvtGm23Y8eOdO7cOfb07+9+97tYb2DYsGH8/ve/p2/fvrRo0YIuXbrwwgsvcOmllwJV5xUaM2YMK1eu5MCBAxw4cICVK1cyZsyYBH9TdUMBQESShpnx3HPP8eyzz9K3b1/OO+88Tj/9dH7605/G5pk4cSKLFy9m0qRJsbK7776b48ePk5mZyYABA7j77rtrvO2FCxcya9YsMjMzycvLY86cOUCQHNDdGT58OACXXnopnTp1onPnztVab5cuXbj77rsZPHgwgwcPZs6cOVWekmooVhIxG1N2drZv3LixsashAgQ7ofjvRdlxqbknnniCyy67jF69ejV2VZLevHnzuPXWW+nQoUOV85rZW+6eXdttqQcgIhJRCgAiIhGlACAiDa6wsJCBAweWKqtNWoaSB8Fq4t1332XQoEFcdNFF/O1vf6N9+/Y1Wr42nn32WQYMGECLFi1qXN/6pAAgIpFRXFzM888/z8SJE3n77bfp06dPg2x34MCBLFu2LHYhOVkoAIhI0hkxYgSzZ88mJyeH8847L3Z75tGjR5kyZQr9+vVjwoQJHD16NLbMypUrGTp0KFlZWVx//fWxF8Kkp6cze/ZssrKyWLJkCQ8//DCPPvooI0eOLLVNd2fWrFkMHDiQCy64IJYqYubMmaxYsQKACRMmMH36dCC4sF3dHEX9+vU75RmFZKAHwUQkKZ04cYINGzbwwgsvMHfuXFavXs2jjz5Ku3bteOedd8jPzycrKwuAvXv38uMf/5jVq1dzxhlnxF7bWHIrZ9euXdm0aRMA77//Pu3bt+d73/teqe0tW7aMvLw8Nm/ezN69exk8eDDDhw9n2LBhvPrqq4wbN45du3axe/duIHhj2JQpU4DgWYFPP/30lDY8+OCDXHbZZfX2O0qUAoCINDgzq7L8uuuuA+Diiy+msLAQgHXr1vHtb38bgMzMzNiDYG+++Sbbtm3jq1/9KgCff/45Q4cOja2rbN6f8rz22mtMnTqVli1b0r17d772ta/xl7/8hWHDhvHwww+zbds2+vfvz4EDB9i9ezdvvPEGjzzyCFDx6yOTnQKAiDS4rl27npISef/+/fTu3Ts23qZNGyDICFqSyrki7s7ll1/OokWLyp2eSDqIlJQUDh48SG5uLsOHD2f//v0888wztG/fnjPPPBNouj2AKq8BmNkTZvaRmW2JK+tiZqvMbHv4s3NYbmb2iJkVmFm+mWXVZ+VFpGlq3749PXr0YO3atUCw88/NzY2lV6jI8OHD+cMf/gDAli1byM/PB2DIkCH8+c9/pqCgAIDDhw/z/vvv16hOw4YNY8mSJRQXF7Nnzx7WrVtHTk5ObP0PP/xw7JTQgw8+GEsGB0EPoLx0EMm884fqXQR+EriyTNmdwBp37wusCccBrgL6hp8ZwKN1U00RaW6eeuopfvSjHzFo0CBGjRrFPffcU+VdOd/4xjc4dOgQ/fr1Y86cObHU0d26dePJJ59k6tSpZGZmMnToUN59990a1WfChAlkZmZy4YUXMmrUKB544AG+9KUvAUFwOHHiBBkZGWRlZbF///5SAaAqzz33HKmpqbzxxhtcffXVSZMLqFqpIMwsHfiTuw8Mx98DRrj7bjPrAbzi7ueb2W/C4UVl56ts/UoFIclEqSDqnlJBVF9TSAXRPW6n/g+gezicAnwYN9/OsExERJJMws8BeHBoVOPDIzObYWYbzWzjnj17Eq2GiIjUUG0DwD/DUz+EPz8Ky3cBPePmSw3LTuHu8909292zu3XrVstqiIhIbdU2AKwApoXD04DlceU3hncDDQE+rur8v4iINI4qnwMws0XACOAsM9sJ3APcBzxjZrcARUDJmxleAMYCBcAR4OZ6qLOIiNSBKnsA7j7V3Xu4e2t3T3X3Be6+z91Hu3tfd7/M3feH87q7z3T3Pu5+gbvr1h4RaTZyc3M5//zzycjI4L777qtwvunTp3P22WeXynh67NgxcnJyuPDCCxkwYAD33HNPQ1S5UkoGJyJSDcXFxcycOZMXX3yRbdu2sWjRIrZt21buvDfddBO5ubmlytq0acPatWvZvHkzeXl55Obm8uabbzZE1SukACAizd6IESNiD4bt27fvlHcRVMeGDRvIyMjg3HPP5bTTTmPKlCksX7683HmHDx9+ynt/zSz27oHjx49z/PjxCnMiNRTlAhKRZq+goIDzzjsPgPz8fC644IJS06uTy2fXrl307PnFTY6pqamsX7++RvUoLi7m4osvpqCggJkzZ3LJJZfUtCl1SgFARJq1oqIiUlJSaNEiOOGRn58fyyJaoqGyebZs2ZK8vDwOHjzIhAkT2LJlS616I3VFp4BEaiA9PR0zw8xIT09v7OpINWzevLnUDv+tt946JQAMGzaMQYMGnfJZvXp1bJ6UlBQ+/PCLRAc7d+4kJaV2iQ46derEyJEjT7lO0NDUAxCpgaKiolheoMY+fyvVk5eXx7FjxwDYvn07y5cv58c//nGpearTAxg8eDDbt29nx44dpKSksHjx4lhm0urYs2cPrVu3plOnThw9epRVq1Yxe/bsmjWmjqkHICLN2ubNmzl58iQXXnghP/zhD+nfvz8LFy6s8XpatWrFL3/5S8aMGUO/fv2YNGkSAwYMiE0fO3Ys//u//wvA1KlTGTp0KO+99x6pqaksWLCA3bt3M3LkSDIzMxk8eDCXX34511xzTZ21szbUAxCRZi0/P59NmzbFXt6SiLFjxzJ27Nhyp73wwgux4YpeTPP2228nXIe6pB6AiDRbn376KWZWJzv/5kgBQESarTPPPLPGbwaLEgUAEZGIUgAQEYkoBQARkYhSABARiSgFABGRiFIAEKmltLS0WFoIpYZoPirL2//QQw8xYMAABg4cyNSpU2NPGJdV3nsD3nvvvVJpJjp06MDDDz/cEE2qkB4EE6mlwsLCUuNKDdE8lOTtb9++PcePH+fSSy/lqquuomfPnjzyyCNs27aNtm3bMmnSJBYvXsxNN91UavmS9wasWrWK1NRUBg8ezLhx4+jfvz95eXmxeVJSUpgwYULDNzCOegAi0mxs3ryZ4cOH079/f1q0aIGZMWfOnBqto7K8/SdOnODo0aOcOHGCI0eOcM4555yyfHXeG7BmzRr69OlDWlpaLVtaN9QDEJFm4dixY0yePJmnnnqKnJwc7r77bo4dO8bcuXNj81Qn7z9UnLf/e9/7Hr169aJt27ZcccUVXHHFFaesqzrvDVi8eDFTp05NuM2JUgAQkWZh9erVZGVlkZOTA0BmZia5ubmlTs1VN+9/eXn7U1JSWL58OTt27KBTp05cf/31/P73v+df//Vfa1TPzz//nBUrVvCzn/2sRsvVBwUAEWkWtmzZUupNX5s2bSIrK6vUPNXtAZSIz9uflpZG79696datGwDXXXcdr7/++ikBoKr3Brz44otkZWXRvXv32jW0DikAiEiz0LVrV9auXQvA+++/z7Jly3j99ddLzVOdHkBFefu7devGm2++yZEjR2jbti1r1qwhOzv7lOWrem/AokWLkuL0DygAiEgzMXXqVFasWMHAgQM566yzWLRoEV27dq3xenbv3s20adMoLi7m5MmTTJo0KZa3f+LEiWRlZdGqVSsuuugiZsyYEVtu7NixPP7445xzzjmx9wYUFxczffr02HsDDh8+zKpVq/jNb35TN41OkJW83agxZWdn+8aNGxu7GiJAcBdI/PcifrzstMqWky888cQTXHbZZfTq1auxq5L05s2bx6233kqHDh2qnNfM3nL3U7sh1aTbQEWkQSg4Vk9D/p4UAESk3rVt25bDhw83djWS3smTJzl69Cht27ZtkO0pAIhIvcvIyOCVV15REKiEu/PSSy+RmppK69atG2SbCV0ENrPbgFsBB/4K3Az0ABYDXYG3gK+7++cJ1lOkXqWnp1NUVATQ6E9nNkfZ2dkcOnSI3/72t/Tq1YvTTjutsauUVE6ePMnevXspLi7m61//eoNtt9YXgc0sBXgN6O/uR83sGeAFYCywzN0Xm9ljwGZ3f7SydekisDS26l7c1UXg2nN3du3axd69ezl+/HhjVyepmBkdO3YkLS2tRsEx0YvAid4G2gpoa2bHgXbAbmAU8C/h9IXAvUClAUBEmj8zIzU1ldTU1MauioRqfQ3A3XcBDwIfEOz4PyY45XPQ3U+Es+0EUspfg4iINKZaBwAz6wyMB3oD5wBnAFfWYPkZZrbRzDbu2bOnttUQqXfxef91fUCak0ROAV0G7HD3PQBmtgz4KtDJzFqFvYBUYFd5C7v7fGA+BNcAEqiHSL0qm/dfpLlI5DbQD4AhZtbOgnR7o4FtwMvAxHCeacDyCpYXEZFGlMg1gPXAUmATwS2gLQiO6GcD/2FmBQS3gi6og3qKiEgdS+guIHe/B7inTPHfgZxE1isiIvVPTwKLiESUAoCISEQpAIiIRJQCgIhIRCkAiIhElAKAiEhEKQCIiESUAoCISEQpAIiIRJQCgERWenq6snxKpCX6QhiRJquoqEhv8JJIUw9ARCSiFABERCJKAUBEJKIUAEREIkoBQEQkohQAREQiSgFARCSiFABERCJKAUBEJKIUAEREIkoBQEQkohQAJDLik78pAZyIksFJhCj5m0hp6gGIiESUAoCISEQpAIiIRFRCAcDMOpnZUjN718zeMbOhZtbFzFaZ2fbwZ+e6qqyIiNSdRHsAvwBy3f3LwIXAO8CdwBp37wusCcdFRCTJ1DoAmFlHYDiwAMDdP3f3g8B4YGE420Lg2sSqKCIi9SGRHkBvYA/wX2b2tpk9bmZnAN3dfXc4zz+A7olWUkRE6l4iAaAVkAU86u4XAYcpc7rHg5uuy73x2sxmmNlGM9u4Z8+eBKohIiK1kUgA2AnsdPf14fhSgoDwTzPrARD+/Ki8hd19vrtnu3t2t27dEqiGiIjURq0DgLv/A/jQzM4Pi0YD24AVwLSwbBqwPKEaiohIvUg0FcS/A0+b2WnA34GbCYLKM2Z2C1AETEpwGyIiUg8SCgDungdklzNpdCLrFRGR+qcngUVEIkoBQEQkohQAREQiSgFARCSiFABERCJKAUBEJKIUAETqSFpaWux9w+np6Y1dHZEq6Z3AInWksLAwNmxmjVcRkWpSD0BEJKIUAEREIkoBQEQkohQAREQiSgFARCSiFABERCJKAUBEJKIUAEREIkoBQEQkohQAREQiSgFARCSiFABERCJKAUCatfT09FiGzrS0tMaujkhSUTZQadaKiopw98auhkhSUg9ARCSiFABERCJKAUBEJKIUAEREIkoBQEQkohIOAGbW0szeNrM/heO9zWy9mRWY2RIzOy3xaoqISF2rix7Ad4B34sbvBx5y9wzgAHBLHWxDRETqWEIBwMxSgauBx8NxA0YBS8NZFgLXJrINkaYoLS0t9gBaenp6Y1dHpFyJPgj2MHAHcGY43hU46O4nwvGdQEqC2xBpcgoLC2PDwXGRSPKpdQ/AzK4BPnL3t2q5/Awz22hmG/fs2VPbaoiISC0lcgroq8A4MysEFhOc+vkF0MnMSnoWqcCu8hZ29/nunu3u2d26dUugGiIiUhu1DgDu/n13T3X3dGAKsNbdbwBeBiaGs00DlidcSxERqXP18RzAbOA/zKyA4JrAgnrYhoiIJKhOsoG6+yvAK+Hw34GculiviIjUHz0JLFLP4m8J1W2hkkz0PgCRehZ/SyjotlBJHuoBiIhElAKAiEhEKQCIiESUAoCISEQpAEizk56eHrvjJi0trbGrI5K0dBeQNDtFRUW4e2NXQyTpqQcgIhJRCgAiIhGlACAiElEKACIiEaUAICISUQoAIiIRpQAgIhJRCgAiIhGlACAiElEKACIiEaUAICISUQoAIiIRpQAgIhJRCgAiIhGlACAiElEKACIiEaUAINLA0tLSYm8sS09Pb+zqSITpjWAiDaywsDA2bGaNVxGJPPUAREQiqtYBwMx6mtnLZrbNzLaa2XfC8i5mtsrMtoc/O9dddUVEpK4k0gM4Adzu7v2BIcBMM+sP3Amscfe+wJpwXEREkkytA4C773b3TeHwp8A7QAowHlgYzrYQuDbBOoqISD2ok2sAZpYOXASsB7q7++5w0j+A7nWxDRERqVsJBwAzaw/8Efiuu38SP83dHfAKlpthZhvNbOOePXsSrYaIiNRQQgHAzFoT7PyfdvdlYfE/zaxHOL0H8FF5y7r7fHfPdvfsbt26JVINERGphUTuAjJgAfCOu8+Lm7QCmBYOTwOW1756Is2bHgqTxpRID+CrwNeBUWaWF37GAvcBl5vZduCycFykTqWnpzeLHWdhYSHujrtTVFTU2NWRiKn1k8Du/hpQ0WOMo2u7XpHqKCoqIrjEpKdpRWpLTwKLiESUAoCISEQpGZw0eSUXUuPHRaRqCgDS5MVn1xSR6tMpIBGRiFIAEBGJKAUAEZGIUgAQEYkoBQARkYhSABARiSgFABGRiFIAEBGJKAUAEZGIUgCQJiE+/bOZKd2DSB1QKghpEuLTP4tI3VAPQCRJxL8drKm/6EaaBgUASVrxp32icMon/u1gZd8Q1lzegCbJRaeAJGnptM8X9AY0qQ/qAUijKntxV0e3Ig1HPQBpVGWP8nV0K9Jw1AOQpBJ/ITQK5/1ro+zFYvWgpLbUA5Ckord7Va2y35F6UFIT6gGIJKna9Ibil1FvQKqiHoBIkqpNbyh+GfUGpCrqAUiDi9r9/clId18JKABIqL53CPHrB2IPO+mcf/2p7MnikruvSj6AgkEE6RSQAPV/O6Ye6mp4ZYNrZX9TnTqKpnrrAZjZlWb2npkVmNmd9bWdZBa1x/fL9iKUvTO6ova/31TVSw/AzFoCvwIuB3YCfzGzFe6+rT62l6yi9vi+jvKlRNT+95uq+uoB5AAF7v53d/8cWAyMr6dt1UhzPTKp6Oi7sjZWdjG2soeNdJTfNFX3ttLqPmiWDN+lKF7MLns9LRH1dQ0gBfgwbnwncEk9batGmuuRSUVH35W1sbIjdl2cbX6q+zet7oNmyfBdimIqkbr8vTfaRWAzmwHMCEc/M7MtDbjtcofrcVtnAXsbaFvVLq9qWg2cBeytixUlqebcvhq3raLvT9n/pYb6nlW2XZrp3y6ujecnsp76CgC7gJ5x46lhWYy7zwfmA5jZRnfPrqe6NDq1r2lrzu1rzm2DaLQvkeXr6xrAX4C+ZtbbzE4DpgAr6mlbIiJSC/XSA3D3E2b2LeAloCXwhLtvrY9tiYhI7dTbNQB3fwF4oZqzz6+veiQJta9pa87ta85tA7WvUqb7tkVEokm5gEREIqrBA4CZ/cjM8s0sz8xWmtk5YbmZ2SNh6oh8M8uKW2aamW0PP9Maus7VZWY/N7N3w/o/Z2ad4qZ9P2zbe2Y2Jq68yaTMMLPrzWyrmZ00s+wy05p8+8pqynUvYWZPmNlH8bdZm1kXM1sVfp9WmVnnsLzC72CyMrOeZvaymW0L/ze/E5Y3+Taa2elmtsHMNodtmxuW9zaz9WEbloQ32mBmbcLxgnB6epUbic8I2BAfoEPc8LeBx8LhscCLgAFDgPVheRfg7+HPzuFw54audzXbdgXQKhy+H7g/HO4PbAbaAL2BvxFcHG8ZDp8LnBbO07+x21FJ+/oR3Hf8CpAdV94s2lemrU227mXaMRzIArbElT0A3BkO3xn3f1rudzCZP0APICscPhN4P/x/bPJtDOvYPhxuDawP6/wMMCUsfwz4Rjj8zbj96RRgSVXbaPAegLt/Ejd6BlByEWI88JQH3gQ6mVkPYAywyt33u/sBYBVwZYNWuprcfaW7nwhH3yR4/gGCti1298/cfQdQQJAuI2lTZpTH3d9x9/fKmdQs2ldGU657jLuvA/aXKR4PLAyHFwLXxpWX9x1MWu6+2903hcOfAu8QZCJo8m0M63goHG0dfhwYBSwNy8u2raTNS4HRVsUTeI1yDcDMfmJmHwI3AHPC4vLSR6RUUp7sphMcaUDza1tZzbF9TbnuVenu7rvD4X8A3cPhJt3m8JTHRQRHys2ijWbW0szygI8IDn7/BhyMO9CMr3+sbeH0j4Gula2/XgKAma02sy3lfMaHlbvL3XsCTwPfqo861Jeq2hbOcxdwgqB9TUp12ifNhwfnC5r8rYBm1h74I/DdMmcZmnQb3b3Y3QcRnE3IAb5cl+uvrwfBLqvmrE8TPCtwDxWnj9gFjChT/krClaylqtpmZjcB1wCjw388qDw1RqUpMxpaDf528ZpM+2qgynQmTdg/zayHu+8OT398FJY3yTabWWuCnf/T7r4sLG5WbXT3g2b2MjCU4LRVq/AoP77+JW3baWatgI7AvsrW2xh3AfWNGx0PvBsOrwBuDK/SDwE+DrtwLwFXmFnn8Er+FWFZ0jGzK4E7gHHufiRu0gpgSniVvjfQF9hA80mZ0Rzb15TrXpUVQMnddNOA5XHl5X0Hk1Z4jnsB8I67z4ub1OTbaGbdLLyT0MzaErxf5R3gZWBiOFvZtpW0eSKwNu4gtHyNcGX7j8AWIB/4byAl7or3rwjOcf2V0neZTCe4sFgA3NzQda5B2woIzsHlhZ/H4qbdFbbtPeCquPKxBHcu/A24q7HbUEX7JhCcc/wM+CfwUnNqXzntbbJ1j2vDImA3cDz8291CcF54DbAdWA10Ceet8DuYrB/gUoLTO/lx37uxzaGNQCbwdti2LcCcsPxcggOsAuBZoE1Yfno4XhBOP7eqbehJYBGRiNKTwCIiEaUAICISUQoAIiIRpQAgIhJRCgAiIhGlACAiElEKABITptbdYWZdwvHO4Xh6OF5lWm4zGxamrs0zsxQzWxqWDzKzsXHzjTCzr9SijoVmdlYN5s8O61OSMrePmf3dzDqE4+WmsS6zjh/UtJ41YWYXmdmCGi7zuJn1r686VbLd1eEDmdIM6DkAKcXM7gAy3H2Gmf0GKHT3n4VBYSOQTfDgzVvAxR5kaI1f/jHgNXf/fZnymwgeuvlWOH4vcMjdH6xh/QrD9eytwTK/Bna6+0/NLBdY6O6Lwh3oIoIcK+cQPDB0nrsXl1n+kLu3L2e9RvAdOlmTNpSznmeBH7v75kTWU8n6S9IG1MW6pgGp7v6TulifNLLGftpNn+T6EKSczQe+C2wFWoflU4HfxM33G2BqmWVvJUg9vIMgz1M6wROMpwEfAHsIntScTZChcVc4PgzoRvCU+F/Cz1fDdXYFVoZ1eRwoAs6qYZs6hXW6A1gdV/594Ptx4y8BQ8ssex9QHNazpE3vAU+FdUojCGQl808EngyHy21TmfWfCbwXN34vQUrfV8O2XkeQ2/6vQG7c3+MVwqdYCdKjbyJ4Z8GauPX8DvgzQZBLB9aGf9s1QK9wvieBR4DXCd61MTEs7wGsC9u9BRgWlncm7t0C+jTtT729FF6aJnc/bmazCHY2V7j78XBSlWl03f1xM7sU+JO7Ly05deTun5vZHEr3ANoS1wMwsz8AD7n7a2bWi2Bn3I8gUeBr7v5DM7uaIJVBTdt00MzuA35N8LKQEikE722orE13mtm3PMjIWJJyuC8wzYN88ljFKdd/UUGb4mUT7GDj9QFGhnV9A/g/7n6HmT0HXA08XzKjmXUDfgsMd/fY6btQf+BSdz9qZv9N0PNZaGbTCXb614bz9SBIqfBlgnwyS4F/IUj18RMzawm0C38fB8KcT13dvdJEY5L8FACkPFcR5I8ZSJCDvCFcBvSP25l2CFP8Dic4Csbd/5+ZHahg+apcRZC/qD/BEXwiikp2/lUot03+xUs+INj57imz3IthIP4rwZvJcsPyvxIcyccbAqzz4EU8uHv8y19WuPvRcHgo4e+RoGfwQNx8z3twGmubmZXkzf8L8IQFmTafd/e8uPk/IjhlpgDQxOkisJRiZoMIsg4OAW6zL96WVN9pdFsAQ9x9UPhJKbOjrKzOE8KLznlW5l3F4fRrCFLjjgF+bmbtwkm1bdPhMuPxF9JOjxuuTpuOllkGgmR7hDvl4+5esv6T1OygrWw9K/JZ3LCF215HEHx3AU+a2Y1x85we1luaOAUAiQkvaj5K8FKND4CfAyUXaRNNy/0pwfnuisZXAv8eV5dB4eA6gtMRmNlVBOegS3H35+J2shvLtKktMA+Y6e5/JUide1c4uaI01mUdD4+EK/JPM+tnZi0IMqZW1aZ47wAZlay7Km8Cw8P6U+YUULzXCVJaQ/AmvlcrW6mZpQH/dPffElx7yQrLDfgSUJhAnSVJKABIvP8LfODuJad9fg30M7OvhacWfsQXFzR/WOZ0Q1VeJjgdkmdmkwlSgZccuQ8Dvg1km1m+mW0D/i1cbi7BDm4rwSmMD2rYpruB59x9Wzh+LzDVzPq6+1aCF2xvIzjNMtPL3AEUmg/km1lFb3i7E/gTwU42Prd8RW2Kcfd3gY5mdmbZadXh7nuAGcAyM9sMLKlg1n8HbjazfODrwHeqWPUIYLOZvQ1MJrieAXAx8KbX0V1F0rh0G6hIIzOz24BP3f3xxq5LVczsFwTXFtY0dl0kceoBiDS+Ryl9Hj6ZbdHOv/lQD0BEJKLUAxARiSgFABGRiFIAEBGJKAUAEZGIUgAQEYmo/w9G9B3JzgdGLAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the difference between the fitted and true vertex\n",
    "ax3 = plt.subplot(111)\n",
    "dxfit = fitVertex1 - 10000*xv\n",
    "con3, bins3, _ = ax3.hist(dxfit, np.linspace(-300.,300.,101), histtype='step', color='black')\n",
    "statBox(ax3, dxfit, bins3)\n",
    "ax3.set_xlim(bins3[0], bins3[-1])\n",
    "_ = ax3.set_xlabel(\"X0 fitted - X0 true (microns)\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "f03cb9af",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXMAAAEGCAYAAACXVXXgAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAiWklEQVR4nO3de3RV5bnv8e8jiIigXBqQBmVREjlcjIjhViUiiFDaAWKRy/BUrIzN2N3sU2stheoApdq9pWUrdddtB6O4jVa5aFE4PRoR0IK2guGWQlRATQQaJXJRLiICz/ljzawdIMlaSVayksnvM8YarDnXO+d8VpSHN+985/OauyMiIo3beakOQEREak/JXEQkBJTMRURCQMlcRCQElMxFREKgaX1e7Bvf+IZHIpH6vKSISKO3YcOGz9w9rao29ZrMI5EI+fn59XlJEaln//jHP/jrX//KwYMHOXXqVKrDaVDOP/980tLSGDp0KBdeeGHCx5lZcbw29ZrMRSTc9uzZw3PPPUdOTg79+/fnvPM0klve8ePH2bp1K8888ww/+MEPqpXQ41EyF5Gkeemllxg5ciQ9e/ZMdSgNViQSYfny5axdu5abbropaefVP5sikhTuzv79++nWrVuqQ2nQzIxu3brx2WefJfW8SuYikhTujrvTtKl+4Y+nWbNmnDhxIqnnVDIXkUbNzLjnnnti23PnzuWBBx5I+PgRI0bQunVrvve97522/6OPPqJ///5kZGQwfvx4jh8/DsBXX33F+PHjycjIoH///hQVFcWO+fd//3cyMjLo1q0br776aq2+V3UpmYtIo3bBBRewdOnSGg9bTJs2jWeeeeas/dOnT+fuu+9m586dtGnThgULFgCwYMEC2rRpw86dO7n77ruZPn06AIWFhSxatIht27aRl5fHv/zLv3Dy5Mmaf7FqUjIXkUatadOmTJkyhUcffbRGxw8dOpRWrVqdts/dWb16NWPHjgVg0qRJvPTSSwAsW7aMSZMmATB27FhWrVqFu7Ns2TImTJjABRdcQJcuXcjIyGD9+vU1/2LVpMEtEWn0pk6dSlZWFj//+c9P2//ss8/ym9/85qz2GRkZvPDCC5Web9++fbRu3To2/t+pUyf27NkDRKdfXnbZZUD0H5JLLrmEffv2sWfPHgYMGBA7R/lj6oOSuYg0ehdffDG33347jz322Glzt2+77TZuu+22FEZWfzTMIqETiUQws1q/VHqicfnJT37CggULOHLkSGzfs88+S+/evc96lQ2fVKZdu3YcPHgwNuNk9+7dpKenA5Cens6uXbsAOHHiBJ9//jnt2rU7bf+Zx9QHJXMJneLi4tg0udq8iovjPkEtDUjbtm0ZN25c7EYlRHvmmzdvPutV1RALRGfI3HDDDbF2ubm5jB49GoBRo0aRm5sLwAsvvMCQIUMwM0aNGsWiRYv46quv+Oijj9ixYwf9+vWro297toSSuZndbWbbzGyrmS00s+Zm1sXM1pnZTjNbbGbN6jpYEZGq3HPPPdWe1TJo0CBuvfVWVq1aRadOnWJTCufMmcMjjzxCRkYG+/btY/LkyQBMnjyZffv2kZGRwSOPPMLDDz8MQM+ePRk3bhw9evRgxIgRPP744zRp0iS5X7AKccfMzSwd+DHQw92/NLMlwARgJPCouy8ys98Dk4En6jRaEZEzHD58OPa+Q4cOHD16tFrHr127tsL93/rWtyqcjdK8eXOef/75Co+57777uO+++6p1/WRJdJilKXChmTUFWgAlwBCg7HeVXODmpEcnIiIJiZvM3X0PMBf4mGgS/xzYABx097LnUXcDFY70m9kUM8s3s/zS0tLkRC0iIqeJm8zNrA0wGugCfBO4CBiR6AXcfb67Z7t7dlpalbXVReQcsnv3bkaPHk1mZiZdu3blrrvuij0ynwzTpk2jZ8+eTJs2jQceeIC5c+cm7dyVyc3NJTMzk8zMzNhN0vqSyDDLjcBH7l7q7l8DS4FrgdbBsAtAJ6D+ZseLSKPm7txyyy3cfPPN7Nixg+3bt3P48OGkjDeXTSecP38+BQUFFT40VBf279/P7NmzWbduHevXr2f27NkcOHCgXq4NiSXzj4EBZtbCzAwYChQCrwNlkzUnAcvqJkQRCZvVq1fTvHlzfvjDHwLQpEkTHn30UZ588kmOHj3KgAED2LZtW6z94MGDyc/P58iRI9x5553069ePq6++mmXLomnnqaeeYtSoUQwZMoShQ4cyatQoDh8+zDXXXMPixYtPu/bmzZsZMGAAWVlZjBkzhgMHDrB3716uueYaALZs2YKZ8fHHHwPQtWvXhG6qvvrqqwwbNoy2bdvSpk0bhg0bRl5eXlJ+XolIZMx8HdEbnRuBvwfHzAemAz81s51AO2BBpScRESln27ZtseRZ5uKLL+byyy9n586djB8/niVLlgBQUlJCSUkJ2dnZ/OpXv2LIkCGsX7+e119/nWnTpsUeEtq4cSMvvPACf/nLX1i+fDkXXnghmzdvZvz48add5/bbb2fOnDkUFBRw5ZVXMnv2bNq3b8+xY8f44osvWLt2LdnZ2axdu5bi4mLat29PixYt4j6AVP4xf2igj/O7+/3A/Wfs/hCovxnxInLOGDduHDfddBOzZ89myZIlsYS5YsUKli9fHhv/PnbsWKwHXdYrrsrnn3/OwYMHuf7664FoAa1bb70VgG9/+9u89dZbrFmzhnvvvZe8vDzcnUGDBgENvzSAarOISL3r0aPHWU9hfvHFF3z88cdkZGTQokUL2rVrR0FBAYsXL+b3v/89EB1r/9Of/nTWakbr1q3joosuqlVMOTk5sd746NGjmTNnDmbGd7/7XSB+0a709HTeeOON2P7du3czePDgWsVUHXqcX0Tq3dChQzl69ChPP/00ACdPnuSee+7hjjvuoEWLFgCMHz+eX//613z++edkZWUBMHz4cP7zP/8Tdwdg06ZN1bruJZdcQps2bWIPCj3zzDOxXvqgQYP44x//SGZmJueddx5t27bl5Zdf5rrrrgPilwYYPnw4K1as4MCBAxw4cIAVK1YwfPjwWv6kEqdkLiL1zsx48cUXef7558nMzOSKK66gefPm/Nu//VuszdixY1m0aBHjxo2L7Zs5cyZff/01WVlZ9OzZk5kzZ1b72rm5uUybNo2srCw2b97MrFmzgGiBNncnJycHgOuuu47WrVvTpk2bhM7btm1bZs6cSd++fenbty+zZs2KO+yTTFb2L1x9yM7O9vz8/Hq7npybzIxk/H+drPOcK06dOsVDDz0US45SuQ8//JA333yT22+/PaH2ZrbB3bOraqOeuUgdileOV2V2JVl0A1SkDpWV461M9NENkdpTz1xEJASUzEWkThQVFdGrV6/T9tWkRkrZ05/V8d5779G7d2+uvvpqPvjgA1q2bFmt42ti//79DBs2jMzMTIYNG1avj/KDkrmIhMzJkyd56aWXGDt2LJs2baJr1671ct2HH36YoUOHsmPHDoYOHRpbtKK+KJmLSEoMHjyY6dOn069fP6644orY3O8vv/ySCRMm0L17d8aMGcOXX34ZO2bFihUMHDiQPn36cOutt8YWpohEIkyfPp0+ffqwePFi5s2bxxNPPMENN9xw2jXdnWnTptGrVy+uvPLKWN2WqVOnsnz5cgDGjBnDnXfeCcCTTz6ZcPGvZcuWMWnSJCD6ZOlLL71U8x9ODegGqIikzIkTJ1i/fj0vv/wys2fPZuXKlTzxxBO0aNGCd999l4KCAvr06QPAZ599xkMPPcTKlSu56KKLYsu6lU2FbNeuHRs3bgRg+/bttGzZkp/97GenXW/p0qVs3ryZLVu28Nlnn9G3b19ycnIYNGgQa9euZdSoUezZs4eSkhIgugrRhAkTgOhDRYcOHTrrO8ydO5cbb7yRTz/9lI4dOwJw6aWX8umnn9bND60SSuYiNRSJROIu+ty5c+d6iqbhqWymTvn9t9xyCwDXXHMNRUVFAKxZs4Yf//jHAGRlZcWe/nz77bcpLCzk2muvBeD48eMMHDgwdq4zC2pV5M0332TixIk0adKEDh06cP311/POO+8waNAg5s2bR2FhIT169ODAgQOUlJTwt7/9jcceewyofHm5yr5jfc9UUjIXqaF40w7Pde3atTvrJuD+/fvp0qVLbPuCCy4AoiVwy+qQV8bdGTZsGAsXLqzw89rUZklPT+fgwYPk5eWRk5PD/v37WbJkCS1btqRVq1ZA/J55hw4dKCkpoWPHjpSUlNC+ffsax1MTGjMXkTrRsmVLOnbsyOrVq4FoIs/Ly4vVOqlMTk4Ozz33HABbt26loKAAgAEDBvDWW2+xc+dOAI4cOcL27durFdOgQYNYvHgxJ0+epLS0lDVr1tCvX7/Y+efNmxcbdpk7d26sYiJEe+YV1Wa58cYbARg1alRsdaHc3FxGjx5drdhqS8lcROrM008/zYMPPkjv3r0ZMmQI999/f9zZJT/60Y84fPgw3bt3Z9asWbG652lpaTz11FNMnDiRrKwsBg4cyHvvvVeteMaMGUNWVhZXXXUVQ4YM4de//jWXXnopEE30J06cICMjgz59+rB///7Tknk8M2bM4LXXXiMzM5OVK1cyY8aMasVWW6rNIqFTX7VZknGdMNV/UW2WxKWkNouZdTOzzeVeX5jZT8ysrZm9ZmY7gj8TKy0mIiJJl8iyce+7e2937w1cAxwFXgRmAKvcPRNYFWyLiEgKVHfMfCjwgbsXA6OB3GB/LnBzEuMSEZFqqG4ynwCUzQvq4O4lwftPgA4VHWBmU8ws38zyS0tLaximyP+IV1b2XJ7bLeeuhJO5mTUDRgHPn/mZR+/gVHgXx93nu3u2u2enpaXVOFCRMmXzuyt7lT18InIuqU7P/DvARncve0b1UzPrCBD8uTfZwYmI1Le8vDy6detGRkZGlcWyKmq3a9cubrjhBnr06EHPnj357W9/W19hVyuZT+R/hlgAlgOTgveTgGXJCkpEJBVOnjzJ1KlTeeWVVygsLGThwoUUFhYm3K5p06b8x3/8B4WFhbz99ts8/vjjFR5fFxJK5mZ2ETAMWFpu98PAMDPbAdwYbIuIpMTgwYNjDxHt27fvrFrqiVi/fj0ZGRl861vfolmzZkyYMIFly87up1bWrmPHjrHCYK1ataJ79+7s2bOndl8sQQnVZnH3I0C7M/btIzq7RUQk5Xbu3MkVV1wBQEFBAVdeeeVpn8errQKwZ88eLrvssthnnTp1Yt26dWcdk0i7oqIiNm3aRP/+/Wv+papBhbZEpNErLi4mPT2d886LDjYUFBTEqi2WqU7Vw9o6fPgw3//+95k3bx4XX3xxvVxTyVxEGr0tW7aclrw3bNhwVkncRHrm6enp7Nq1K/bZ7t27SU9PP+uYqtp9/fXXfP/73+e2226LlfitD0rmItLobd68mWPHjgGwY8cOli1bxkMPPXRam0R65n379mXHjh189NFHpKens2jRolgFx0TauTuTJ0+me/fu/PSnP03Ol0uQqiaKSKO3ZcsWTp06xVVXXcUvf/lLevToEStHWx1Nmzbld7/7HcOHD6d79+6MGzeOnj17xj4fOXIk//jHPypt99Zbb/HMM8+wevVqevfuTe/evXn55ZeT+VUrj71eriIiUocKCgrYuHFjbCGJ2hg5ciQjR46s8LPyibmidtddd13KqmCqZy4ijdqhQ4cws6Qk8sZMyVxEGrVWrVpVe8WhMFIyFxEJASVzEZEQ0A1QkUp07twZM6vyc5GGQslcpBIqpSuNiYZZRERCQMlcROQMd955J+3bt49bebGy2ueJ1kRPJiVzEZEz3HHHHeTl5VXZprKa5onWRE82JXMRCY0tW7aQk5NDjx49OO+88zAzZs2aVe3z5OTk0LZt2yrbVFbTPNGa6MmmG6AiEgrHjh1j/PjxPP300/Tr14+ZM2dy7NgxZs+eHWuTSOXERFVW0zzRmujJllAyN7PWwB+AXkQXbr4TeB9YDESAImCcux+oiyBFROJZuXIlffr0oV+/fgBkZWWRl5d32vTS+qxpXt8S7Zn/Fshz97Fm1gxoAdwLrHL3h81sBjADmF5HcYqEUry57GVtNE0yvq1bt562utDGjRtjS7iVSWbPvLKa5onWRE+2uMnczC4BcoA7ANz9OHDczEYDg4NmucAbKJmLVEsiSTpespeodu3asXr1agC2b9/O0qVL+etf/3pam2T2zCurad6tW7eEaqInWyI3QLsApcB/m9kmM/tDsMBzB3cvCdp8AnSo6GAzm2Jm+WaWX1pampyoJbQikQhmVuVLT15KRSZOnMjhw4fp1asXU6ZMYeHChbRr1y7+gZWca+DAgbz//vt06tSJBQsWAP9Tzxwqr30eryZ6XbF4tXfNLBt4G7jW3deZ2W+BL4D/4+6ty7U74O5tqjpXdna25+fn1z5qCS0zS1k96IaqsfxMTp06xYMPPsj999+f6lAavA8++IC33nqL22+/PaH2ZrbB3bOrapNIz3w3sNvdy27HvgD0AT41s47BhToCexOKSkRCycxo1qwZR48eTXUoDd6RI0do0aJFUs8ZN5m7+yfALjPrFuwaChQCy4FJwb5JQN1PpBSRBsvMyMjI4JVXXuHUqVOpDqfBOnToEGvWrCEjIyOp5407zAJgZr2JTk1sBnwI/JDoPwRLgMuBYqJTE/dXdR4Ns0g8jWVIoT41pp/J119/zaJFi/jyyy/p0KEDTZo0SXVIDcrx48cpKiqib9++DBo0KOHjEhlmSSiZJ4uSucTTmBJXfWlsP5MTJ05QXFzMwYMH1UM/w/nnn0/79u355je/Wa3jEknmegJURJKqadOmdO3aNdVhnHNUm0VEJASUzEVEQkDJXEQkBJTMpV7Fe8JTT3eK1IxugEq9Ki4ublQzM0QaC/XMRURCQMlcRCQElMxFREJAyVxEJASUzEVEQkDJXEQkBJTMJWm0SpBI6mieuSSN5pCLpI565iIiIaBkLiISAgkNs5hZEXAIOAmccPdsM2sLLAYiQBHRlYYO1E2YIiJSler0zG9w997lVruYAaxy90xgVbAtIiIpUJthltFAbvA+F7i51tGIiEiNJJrMHVhhZhvMbEqwr4O7lwTvPwE6JD06ERFJSKJTE69z9z1m1h54zczeK/+hu7uZVTgnLUj+UwAuv/zyWgUrIiIVS6hn7u57gj/3Ai8C/YBPzawjQPDn3kqOne/u2e6enZaWlpyoRUTkNHGTuZldZGatyt4DNwFbgeXApKDZJGBZXQUpci7r3LlzlU/VRiKRVIcoDUAiwywdgBfNrKz9c+6eZ2bvAEvMbDJQDIyruzBFzl1FRUVVfh783ZRzXNxk7u4fAldVsH8fMLQughIRkerRE6AiIiGgZC4iEgJK5iKNXLwbpLpJem5QCVyRRi7eDVLQTdJzgXrmIiIhoGQuIhICSuYiIiGgZC4iEgJK5iIiIaBkLiISAkrmIiIhoGQuIhICSuYiIiGgZC4iEgJK5iIiIaBkLiISAkrmIiIhkHAyN7MmZrbJzP4cbHcxs3VmttPMFptZs7oLU0REqlKdnvldwLvltucAj7p7BnAAmJzMwEREJHEJJXMz6wR8F/hDsG3AEOCFoEkucHMdxCciIglItGc+D/g5cCrYbgccdPcTwfZuIL2iA81sipnlm1l+aWlpbWIVEZFKxE3mZvY9YK+7b6jJBdx9vrtnu3t2WlpaTU4hIiJxJLJs3LXAKDMbCTQHLgZ+C7Q2s6ZB77wTsKfuwpRUi0QiFBcXV9mmc+fO9RSNiJwpbs/c3X/h7p3cPQJMAFa7+23A68DYoNkkYFmdRSkpV1xcjLtX+UpkLUoRqRu1mWc+Hfipme0kOoa+IDkhiYhIdSUyzBLj7m8AbwTvPwT6JT8kERGpLj0BKiISAkrmIueAzp07Y2ZVviKRSKrDlFqo1jCLiDROidycjj4LKI2VeuYiIiGgZC4iEgJK5iIiIaBkLiISAkrmIiIhoGQuIhICSuYiIiGgZC4iEgJK5iIiIaBkLiISAkrmIiIhoGQuIhICSuYiIiGQyILOzc1svZltMbNtZjY72N/FzNaZ2U4zW2xmzeo+XBERqUgiPfOvgCHufhXQGxhhZgOAOcCj7p4BHAAm11mUIiJSpUQWdHZ3Pxxsnh+8HBgCvBDszwVurosARUQkvoTGzM2siZltBvYCrwEfAAfd/UTQZDeQXsmxU8ws38zyS0tLkxCyiIicKaFk7u4n3b030InoIs7/K9ELuPt8d8929+y0tLSaRSkiIlWq1mwWdz8IvA4MBFqbWdmyc52APckNTUREEpXIbJY0M2sdvL8QGAa8SzSpjw2aTQKW1VGMUg8ikUiVi/127tw51SGKSBUSWdC5I5BrZk2IJv8l7v5nMysEFpnZQ8AmYEEdxil1rLi4GHdPdRgiUkNxk7m7FwBXV7D/Q6Lj59LARSIRiouLq2yjnrdI45ZIz1waOfW6RcJPj/OLiISAkrmISAgomYuIhICSuYhICCiZi4iEgJK5iEgIKJmLiISAkrmISAgomYuIhICSuYhICCiZi4iEgJK5iEgIKJmLiISAkrmISAgomYuIhEAiy8ZdZmavm1mhmW0zs7uC/W3N7DUz2xH82abuwxURkYok0jM/Adzj7j2AAcBUM+sBzABWuXsmsCrYFhGRFIibzN29xN03Bu8PEV3MOR0YDeQGzXKBm+soRhERiaNaY+ZmFiG6Hug6oIO7lwQffQJ0qOSYKWaWb2b5paWltYlVREQqkXAyN7OWwJ+An7j7F+U/8+gCkxUuMunu8909292z09LSahWsiIhULKFkbmbnE03kz7r70mD3p2bWMfi8I7C3bkIUEZF4EpnNYsAC4F13f6TcR8uBScH7ScCy5IcnIiKJaJpAm2uBHwB/N7PNwb57gYeBJWY2GSgGxtVJhCIiElfcZO7ubwJWycdDkxuOiDRUkUiE4uLiKtt07tyZoqKi+glITpNIz1xEhOLiYqJzHSoXHZWVVNDj/CIiIaBkLiISAkrmIiIhoDFzEQGiNy+rGvPu3LlzPUYj1aVkLiIAmoXSyGmYRUQkBJTMRURCQMlcRCQElMxFREJAyVxEJASUzEVEQkDJPAQikQhmVulL84NFwk/zzEMgkQJIIhJu6pmLiISAkrmISAgksmzck2a218y2ltvX1sxeM7MdwZ9t6jZMERGpSiI986eAEWfsmwGscvdMYFWwLSLnuLJiXZW9IpFIqkMMrUSWjVtjZpEzdo8GBgfvc4E3gOnJDExEGp94xbq0ElHdqemYeQd3LwnefwJ0SFI8IiJSA7W+AerROXGVzoszsylmlm9m+aWlpbW93Dkn3hxyzSMXEaj5PPNPzayju5eYWUdgb2UN3X0+MB8gOztbk6GrSXPIRSQRNe2ZLwcmBe8nAcuSE46IiNREIlMTFwJ/A7qZ2W4zmww8DAwzsx3AjcG2iIikSCKzWSZW8tHQJMciIiI1pCdARURCQMlcRCQElMxFREJAyVxEJASUzEVEQkDJPMW0SpCIJINWGkoxPeEpIsmgnrmINCiJ1CNSKd2zqWcuIg1KIr+tqpTu2dQzFxEJAfXMRaTelK1EFK+NVJ+SuYjUm3grEUnNaZilhhK5SZPIS70QkerTWqNnU8+8hjSlUCR1tNbo2dQzFxEJASVzEZEQUDIXEQmBWiVzMxthZu+b2U4zm5GsoBoC1UwRabzi3SAN403SGt8ANbMmwOPAMGA38I6ZLXf3wmQFl0q6wSnSeCUyBTJsN0lr0zPvB+x09w/d/TiwCBidnLBERKQ6ajM1MR3YVW57N9D/zEZmNgWYEmx+ZWZba3HN+vIN4LNG8C/3N4DPUh1EHI0hRlCcydYo4jSzRhEn0C1egzqfZ+7u84H5AGaW7+7ZdX3N2lKcydMYYgTFmWyKM7nMLD9em9oMs+wBLiu33SnYJyIi9aw2yfwdINPMuphZM2ACsDw5YYmISHXUeJjF3U+Y2b8CrwJNgCfdfVucw+bX9Hr1THEmT2OIERRnsinO5Iobp2n6nYhI46cnQEVEQkDJXEQkBFKWzM3sHjPzYJ5ng2JmD5pZgZltNrMVZvbNVMdUETP7jZm9F8T6opm1TnVMFTGzW81sm5mdMrMGNw2sMZSlMLMnzWxvQ35Ow8wuM7PXzaww+O99V6pjqoiZNTez9Wa2JYhzdqpjqoqZNTGzTWb256rapSSZm9llwE3Ax6m4fgJ+4+5Z7t4b+DMwK8XxVOY1oJe7ZwHbgV+kOJ7KbAVuAdakOpAzlStL8R2gBzDRzHqkNqoKPQWMSHUQcZwA7nH3HsAAYGoD/Vl+BQxx96uA3sAIMxuQ2pCqdBfwbrxGqeqZPwr8HGiQd1/d/YtymxfRcONc4e4ngs23ic71b3Dc/V13fz/VcVSiUZSlcPc1wP5Ux1EVdy9x943B+0NEE1B6aqM6m0cdDjbPD14N8u+4mXUCvgv8IV7bek/mZjYa2OPuW+r72tVhZr8ys13AbTTcnnl5dwKvpDqIRqiishQNLgE1NmYWAa4G1qU4lAoFQxebgb3Aa+7eIOME5hHt+J6K17BOHuc3s5XApRV8dB9wL9EhlpSqKkZ3X+bu9wH3mdkvgH8F7q/XAAPx4gza3Ef0V9xn6zO28hKJU84NZtYS+BPwkzN+y20w3P0k0Du4z/SimfVy9wZ1P8LMvgfsdfcNZjY4Xvs6SebufmNF+83sSqALsCUoYtUJ2Ghm/dz9k7qIpTKVxViBZ4GXSVEyjxenmd0BfA8Y6il8aKAaP8+GRmUpksjMzieayJ9196Wpjicedz9oZq8TvR/RoJI5cC0wysxGAs2Bi83sj+7+vytqXK/DLO7+d3dv7+4Rd48Q/ZW2T30n8njMLLPc5mjgvVTFUhUzG0H0V7BR7n401fE0UipLkSQW7aEtAN5190dSHU9lzCytbOaXmV1IdE2GBvd33N1/4e6dglw5AVhdWSIHzTOvzMNmttXMCogOCTXIKVbA74BWwGvBNMrfpzqgipjZGDPbDQwE/p+ZvZrqmMoEN5DLylK8CyxJoCxFvTOzhcDfgG5mttvMJqc6pgpcC/wAGBL8/7g56FU2NB2B14O/3+8QHTOvctpfY6DH+UVEQkA9cxGREFAyFxEJASVzEZEQUDIXEQkBJXMRkRBQMpcGx8xOBtPatprZ82bWIk77N8qqMZpZUXUqcQYVHd8NHhwpvz87qKjXLNjuamYfmtnFwfYvgiqL75vZ8HjnDs73WLB/sJl9O9EYRRKhZC4N0Zfu3tvdewHHgX+uw2tNBv7J3W8ov9Pd84G/AD8Ldj1OtDTBF0ElwAlAT6JPDv5XUH2x0nO7e767/zjYPxhQMpekUjKXhm4tkBH0ZmMPdpjZ74JSBgkxs4lm9vegtz8n2DcLuA5YYGa/qeCwe4F/MrOfA03dfWGwfzSwyN2/cvePgJ1Eqy+Wv95p5y6LPyhA9c/A3cFvH4MS/Q4iVamT2iwiyWBmTYnWGc+r5Xm+CcwBrgEOACvM7GZ3/6WZDQF+FvTETxPU7XgY+C+itc7LpBMtOVzmrEqLZ567rFCSuxcFT+oedve5tfleIuWpZy4N0YVBedJ8oguYLKjl+foCb7h7afD4/rNAToLHfgf4lNOTuUiDo565NERfBqs8xZjZCU7vfDSv6yCCEqSXAMOJlkl9NShopkqL0uCoZy6NRTHQw8wuCCreDa3GseuB683sG8GNyolEb25WKqim9wgw1d3/DiwjWo8folUVJwSxdAEyg2sk6hDRAmkiSaNkLo2Cu+8ClhCtOb0E2FSNY0uAGcDrwBZgQwILZswEXnT3wmD7AaLrg2YGVRWXAIVEx/OnBosdJOr/AmN0A1SSSVUTRURCQD1zEZEQUDIXEQkBJXMRkRBQMhcRCQElcxGREFAyFxEJASVzEZEQ+P9nOhyAT9AivAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the pull\n",
    "ax4 = plt.subplot(111)\n",
    "con4, bins4, _ = ax4.hist(pull1, np.linspace(-4.,4.,40), histtype='step', color='black')\n",
    "statBox(ax4, pull1, bins3)\n",
    "ax4.set_xlim(bins4[0], bins4[-1])\n",
    "_ = ax4.set_xlabel(\"Pull of X0 fit\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "57e3ab4b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAEGCAYAAABvtY4XAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAjf0lEQVR4nO3de3hU9bX/8fcCRERELgKlQRMlkRIwYggXW0EEUQ/1cFHkcmjFyjk8p/V3Wq1S/NWfVFpaRamibY89/IrHiMrFSyGnRQRBCtIKjRIioBWqoCBV5CI3EYF1/pidIUBCZiYzmWTzeT3PPNn3veZLWPnOd/Ze29wdEREJn3rpDkBERFJDCV5EJKSU4EVEQkoJXkQkpJTgRURCqkFNnuy8887zrKysmjyliEid98Ybb3zq7q3i3a9GE3xWVhbFxcU1eUoRqWEfffQRf/7zn9m9ezdHjx5Ndzi1yhlnnEGrVq3o168fZ511Vsz7mdnmRM5XowleRMJt69atPPvss/Tu3ZsePXpQr55Ggcs7dOgQa9euZcaMGXz729+OK8knQgleRJJm7ty5DBgwgE6dOqU7lForKyuLoqIili9fzjXXXJPSc+nPq4gkhbuzc+dOOnTokO5QajUzo0OHDnz66acpP5cSvIgkhbvj7jRooIGBqjRs2JDDhw+n/DxK8CJSp5kZd955Z3R+ypQp3HfffTHvf91119GsWTOuv/7645a///779OjRg+zsbIYPH86hQ4cA+OKLLxg+fDjZ2dn06NGDTZs2Rfe5//77yc7OpkOHDrz88svVel/JoAQvInXamWeeyYsvvpjwkMe4ceOYMWPGScvHjx/PHXfcwcaNG2nevDnTp08HYPr06TRv3pyNGzdyxx13MH78eADWr1/PrFmzWLduHQsWLOB73/seR44cSfyNJYESvIjUaQ0aNGDs2LE88sgjCe3fr18/zjnnnOOWuTtLlixh6NChAIwePZq5c+cCMG/ePEaPHg3A0KFDWbx4Me7OvHnzGDFiBGeeeSYXXngh2dnZrFq1KvE3lgQaLBOROu+2224jLy+PH/3oR8ctf+aZZ3jooYdO2j47O5vnn3++0uPt2LGDZs2aRb9PaNeuHVu3bgUil4Kef/75QOSPy7nnnsuOHTvYunUrPXv2jB6j/D7pogQvInVe06ZNufnmm3nssceOu7Z81KhRjBo1Ko2RpVeoh2iysrIws7hfKqcgUvfcfvvtTJ8+nf3790eXPfPMM3Tp0uWkV9nQS2VatmzJ7t27o1e6bNmyhYyMDAAyMjL48MMPATh8+DCfffYZLVu2PG75ifukS6gT/ObNm6OXbsXz2rw5obuCRSSNWrRowbBhw6JfhkKkB19SUnLS61TDMxC5Mueqq66KbldYWMigQYMAGDhwIIWFhQA8//zz9O3bFzNj4MCBzJo1iy+++IL333+fDRs20L179xS929iEOsGLyOnlzjvvjPtqml69enHTTTexePFi2rVrF728cfLkyTz88MNkZ2ezY8cOxowZA8CYMWPYsWMH2dnZPPzwwzzwwAMAdOrUiWHDhpGbm8t1113Hb37zG+rXr5/cNxinmMbgzawZ8DugM+DArcDfgNlAFrAJGObuu1IRpIhIZfbt2xedbtOmDQcOHIhr/+XLl1e4/KKLLqrwKphGjRrx3HPPVbjPPffcwz333BPX+VMp1h78o8ACd/8acCnwNnA3sNjdc4DFwbyIiNQSVSZ4MzsX6A1MB3D3Q+6+GxgEFAabFQKDUxOiiIgkIpYe/IXAduC/zWy1mf3OzM4G2rj7tmCbfwBtKtrZzMaaWbGZFW/fvj05UadYZmamrr4RSbEtW7YwaNAgcnJyaN++PT/4wQ+i5QCSYdy4cXTq1Ilx48Zx3333MWXKlKQduzKFhYXk5OSQk5MT/SI2nWJJ8A2AfOBxd78M2M8JwzHu7kTG5k/i7tPcvcDdC1q1ivuBJGmxadMmXX0jkkLuzg033MDgwYPZsGED7777Lvv27UvK+HXZpY3Tpk2jtLS0whudUmHnzp1MnDiRlStXsmrVKiZOnMiuXen9WjKWBL8F2OLuK4P554kk/I/NrC1A8POT1IQoImGzZMkSGjVqxHe+8x0A6tevzyOPPMITTzzBgQMH6NmzJ+vWrYtu36dPH4qLi9m/fz+33nor3bt357LLLmPevHkAPPnkkwwcOJC+ffvSr18/Bg4cyL59++jatSuzZ88+7twlJSX07NmTvLw8hgwZwq5du/jkk0/o2rUrAGvWrMHM+OCDDwBo3759TF/cvvzyy/Tv358WLVrQvHlz+vfvz4IFC5LSXomqMsG7+z+AD82srMhzP2A9UASMDpaNBualJEIRCZ1169ZFE2qZpk2bcsEFF7Bx40aGDx/OnDlzANi2bRvbtm2joKCAn//85/Tt25dVq1bx6quvMm7cuOiNTW+++SbPP/88f/rTnygqKuKss86ipKSE4cOHH3eem2++mcmTJ1NaWsoll1zCxIkTad26NQcPHmTPnj0sX76cgoICli9fzubNm2ndujWNGzeu8qap8iUMoG6VKvgP4Bkzawi8B3yHyB+HOWY2BtgMDEtNiFKZrKysuIeFMjMzjytvKnIq8f6OZWVlcfToUZo1a8btt9+e8HmHDRvGNddcw8SJE5kzZ040iS5cuJCioqLoePrBgwejPe2y3vOpfPbZZ+zevZsrr7wSiBQRu+mmmwD4+te/zooVK1i2bBk//vGPWbBgAe5Or169gLpZ9iCmBO/uJUBBBav6JTUaiUvZnbrxMLMURSNhFM/v2NGjR5k0aRITJkyosh57bm7uSXeT7tmzhw8++IDs7GwaN25My5YtKS0tZfbs2fz2t78FImP3L7zwwklPjVq5ciVnn3127G+sAr1794722gcNGsTkyZMxM775zW8CVRcuy8jIYOnSpdHlW7ZsoU+fPtWKqbp0J6uI1Lh+/fpx4MABnnrqKQCOHDnCnXfeyS233ELjxo0BGD58OA8++CCfffYZeXl5AFx77bX86le/iv7RWb16dVznPffcc2nevHn05qYZM2ZEe/O9evXi6aefJicnh3r16tGiRQvmz5/PFVdcAVRd9uDaa69l4cKF7Nq1i127drFw4UKuvfbaarZU9SjBi0iNMzN+//vf89xzz5GTk8PFF19Mo0aN+MUvfhHdZujQocyaNYthw46N/t577718+eWX5OXl0alTJ+699964z11YWMi4cePIy8ujpKSECRMmAJHhJXend+/eAFxxxRU0a9aM5s2bx3TcFi1acO+999KtWze6devGhAkTqhwySjWL9yN+dRQUFHhxcXGNnc/M4h7CCPv5ajpGqdvi+X05cYgmnsfmnW7ee+89XnvtNW6++eaYtjezN9y9omHyU1IPXkQkpJTgRURCSgleRCSklOBFJCU2bdpE586dj1uWSE2YsrtY4/HOO+/QpUsXLrvsMv7+97/TpEmTuPZPxM6dO+nfvz85OTn0798/7WUKQAleRELmyJEjzJ07l6FDh7J69Wrat29fI+d94IEH6NevHxs2bKBfv37RB4GkkxK8iKRFnz59GD9+PN27d+fiiy+OXpv++eefM2LECDp27MiQIUP4/PPPo/ssXLiQyy+/nPz8fG666abowz6ysrIYP348+fn5zJ49m6lTp/L4449z1VVXHXdOd2fcuHF07tyZSy65JFqn5rbbbqOoqAiAIUOGcOuttwLwxBNPxFwAbd68eYweHaneMnr0aObOnZt44yRJrKUKRESS7vDhw6xatYr58+czceJEXnnlFR5//HEaN27M22+/TWlpKfn5+QB8+umnTJo0iVdeeYWzzz47+ki9suvYW7ZsyZtvvgnAu+++S5MmTbjrrruOO9+LL75ISUkJa9as4dNPP6Vbt2707t2bXr16sXz5cgYOHMjWrVvZti1SCX358uWMGDECiNwItXfv3pPew5QpU7j66qv5+OOPadu2LQBf+cpX+Pjjj1PTaHFQgheRlKisLEb55TfccAMAXbt2jdZIWrZsGd///vcByMvLi97F+vrrr7N+/Xq+8Y1vAHDo0CEuv/zy6LFOLCpWkddee42RI0dSv3592rRpw5VXXslf//pXevXqxdSpU1m/fj25ubns2rWLbdu28Ze//IXHHnsMqPzRfpW9x9pQFkQJXkRSomXLlid90bhz504uvPDC6PyZZ54JRMoFl9Vxr4y7079/f2bOnFnh+urUosnIyGD37t0sWLCA3r17s3PnTubMmUOTJk0455xzgKp78G3atGHbtm20bduWbdu20bp164TjSRaNwYtISjRp0oS2bduyZMkSIJLcFyxYEK3tUpnevXvz7LPPArB27VpKS0sB6NmzJytWrGDjxo0A7N+/n3fffTeumHr16sXs2bM5cuQI27dvZ9myZXTv3j16/KlTp0aHbKZMmRKtJAmRHnxFtWiuvvpqAAYOHBh9ilNhYSGDBg2KK7ZUUIIXkZR56qmn+NnPfkaXLl3o27cvP/nJT6q8quW73/0u+/bto2PHjkyYMCFaN75Vq1Y8+eSTjBw5kry8PC6//HLeeeeduOIZMmQIeXl5XHrppfTt25cHH3yQr3zlK0Ak+R8+fJjs7Gzy8/PZuXPncQm+KnfffTeLFi0iJyeHV155hbvvvrvqnVJMtWhOs/OpFo3EQ7VoUkO1aEREpFqU4EVEQkoJXkQkpJTgRURCSgleRCSklOBFRMpZsGABHTp0IDs7u8qCYUeOHOGyyy7j+uuvB+Bvf/sbXbp0ib6aNm3K1KlTayDqiulO1iTKzMys0duTMzMzE9onkRgzMzOjt5KLhNWRI0e47bbbWLRoEe3ataNbt24MHDiQ3NzcCrd/9NFH6dixI3v27AGgQ4cOlJSURI+VkZHBkCFDair8k6gHn0SbNm3C3WvslUjCTTTGzZs3J7/BRJKoT58+0RufduzYcVIt+lisWrWK7OxsLrroIho2bMiIESOYN29ehdtu2bKFP/7xj/zrv/5rhesXL15M+/btE+qIJYt68CISChs3buTiiy8GoLS0lEsuueS49VXVkgHYunUr559/fnRdu3btWLlyZYXnu/3223nwwQcrPCbArFmzGDlyZELvJVmU4EWkztu8eTMZGRnUqxcZlCgtLY1WoSwTTzXIqvzhD3+gdevWdO3alaVLl560/tChQxQVFXH//fcn7ZyJiCnBm9kmYC9wBDjs7gVm1gKYDWQBm4Bh7p7+Z1SJyGlnzZo1xyX0N95446TywbH04DMyMvjwww+j67Zs2UJGRsZJ+6xYsYKioiLmz5/PwYMH2bNnD9/61rd4+umnAXjppZfIz8+nTZs2SXl/iYqnB3+Vu39abv5uYLG7P2Bmdwfz45ManYhIDEpKSjh48CAAGzZsYN68eUyaNOm4bWLpwXfr1o0NGzbw/vvvk5GRwaxZs6KVLcu7//77o73zpUuXMmXKlGhyB5g5c2bah2egel+yDgIKg+lCYHC1oxERScCaNWs4evQol156KT/96U/Jzc2Nlu6NR4MGDfj1r3/NtddeS8eOHRk2bBidOnWKrh8wYAAfffTRKY+xf/9+Fi1aFH2YSTrF2oN3YKGZOfBf7j4NaOPu24L1/wAq/CxiZmOBsQAXXHBBNcOV00VWVlZCV+7ocs7TU2lpKW+++Wb04RzVMWDAAAYMGFDhuvnz55+0rE+fPvTp0yc6f/bZZ7Njx45qx5EMsSb4K9x9q5m1BhaZ2XFFmN3dg+R/kuCPwTSIlAuuVrRy2ti8eXNCZY1rw2PSpGbt3bsXM0tKcg+bmIZo3H1r8PMT4PdAd+BjM2sLEPz8JFVBiohU5pxzzon7yU6niyoTvJmdbWbnlE0D1wBrgSJgdLDZaKDiuwFERCQtYhmiaQP8Pvjo2wB41t0XmNlfgTlmNgbYDAxLXZgiIhKvKhO8u78HXFrB8h1Av1QEJSIi1adaNCIiIaUELyISUkrwIiIVqKou/IcffshVV11Fbm4unTp14tFHH42uqy114VVsTETkBLHUhW/QoAG//OUvyc/PZ+/evXTt2pX+/fuTm5tba+rCqwcvIqGyZs0aevfuTW5uLvXq1cPMmDBhQlzHiKUufNu2bcnPzwci1+J37NiRrVu3nnSsdNaFVw9eRELj4MGDDB8+nKeeeoru3btz7733cvDgQSZOnBjdJtl14SHyIJ3Vq1fTo0ePk9alsy68EryIhMYrr7xCfn4+3bt3ByAvL48FCxYcV8IimXXhAfbt28eNN97I1KlTadq06XHr0l0XXgleUqo6RcPkZCrCdmpr16497klOb775ZnQYpUwy68J/+eWX3HjjjYwaNarC6pHprguvBC8plWjRMKmYirCdWsuWLVmyZAkA7777Li+++CJ//vOfj9smWXXh3Z0xY8bQsWNHfvjDH1Z4nHTXhdeXrCISGiNHjmTfvn107tyZsWPHMnPmTFq2bBn3cU5VF76sJvyKFSuYMWMGS5YsiV4OWb6ccG2oC68evMQkMzMzoV5gTQ+1VCfO02EII9XS/WmtSZMm/M///E9SjlVZXfiyJP7Vr371lO/3VHXha6qdlOAlJnUl+SUa5+kyhJFKZkbDhg05cOBAukOp9fbv30/jxo1Tfh4N0YhIUpgZ2dnZvPTSS2nvyddme/fuZdmyZWRnZ6f8XOrBi0jSDB48mFmzZrFo0SLy8/OpX79+ukOqVQ4dOsSmTZvo1q0bXbp0Sfn5rCb/0hYUFHhxcXGNnc/M1JOQmNSV35VE46zJ/Q4fPswZZ5xBcXExR48ejfucYXbGGWfQunVrvvrVr8a1n5m94e4F8Z5PPXgRSaoGDSJppWvXrmmORDQGLyISUurBi5wG6splrpJcSvAip4G6cpmrJJeGaEREQkoJXkQkpJTgRURCSgleRCSklOBFREIq5gRvZvXNbLWZ/SGYv9DMVprZRjObbWYNUxemiIjEK54e/A+At8vNTwYecfdsYBcwJpmBiYhI9cSU4M2sHfBN4HfBvAF9geeDTQqBwSmIT0REEhRrD34q8COgrHJQS2C3ux8O5rcAJz+wEDCzsWZWbGbF27dvr06sIiIShyoTvJldD3zi7m8kcgJ3n+buBe5e0KpVq0QOISIiCYilVME3gIFmNgBoBDQFHgWamVmDoBffDtiaujBFRCReVfbg3f3/uns7d88CRgBL3H0U8CowNNhsNDAvZVGKiEjcqnMd/Hjgh2a2kciY/PTkhCQiIskQVzVJd18KLA2m3wO6Jz8kERFJBt3JKiISUkrwIhx7IEa8r6ysrHSHLlIpPfBDhMQfiJHIU5JEaop68CIiIaUELyISUkrwIiIhpQQvIhJSSvAiIiGlBC8itUZWVladuFy1rsSpyyRFpNbYvHkz7h73fjV9uWpdiVM9eBGRkFKCF6kG3QErtZmGaESqQXfASm2mHryISEgpwYtI0iU6dJWZmVmjcSZ6NUxNx5koDdGISNIlOnRV0xK9GqauUA9eRCSklOBFREKqTiT4sI+TiYikQp0Ygw/7OJmISCrUiR68iIjEr0704EXCpuwywkT2k5OpPSumBC+SBnXlMsK6Qu1ZMQ3RiIiEVJUJ3swamdkqM1tjZuvMbGKw/EIzW2lmG81stpk1TH24IiISq1h68F8Afd39UqALcJ2Z9QQmA4+4ezawCxiTsihFRCRuVSZ4j9gXzJ4RvBzoCzwfLC8EBqciQBERSUxMY/BmVt/MSoBPgEXA34Hd7n442GQLkJGSCEVEJCExJXh3P+LuXYB2QHfga7GewMzGmlmxmRVv3749sShFRCRucV1F4+67gVeBy4FmZlZ2mWU7YGsl+0xz9wJ3L2jVqlV1YhURkTjEchVNKzNrFkyfBfQH3iaS6IcGm40G5qUoRhERSUAsNzq1BQrNrD6RPwhz3P0PZrYemGVmk4DVwPQUxikiInGqMsG7eylwWQXL3yMyHi8iIrWQ7mQVEQkpJXgRkZBSghcRCSkleBGRkFKCFxEJKSV4EZGQUoIXEQkpJXgRkZBSghcRCSkleBGRkFKCFxEJKSV4EZGQUoIXEQkpJXgRkZBSghcRCSkleBGRkFKCFxEJKSV4EZGQUoIXEQkpJXgRkZBSghcRCSkleBGRkFKCFxEJKSV4EZGQUoIXEQmpKhO8mZ1vZq+a2XozW2dmPwiWtzCzRWa2IfjZPPXhiohIrGLpwR8G7nT3XKAncJuZ5QJ3A4vdPQdYHMyLiEgtUWWCd/dt7v5mML0XeBvIAAYBhcFmhcDgFMUoIiIJiGsM3syygMuAlUAbd98WrPoH0KaSfcaaWbGZFW/fvr06sYqISBxiTvBm1gR4Abjd3feUX+fuDnhF+7n7NHcvcPeCVq1aVStYERGJXUwJ3szOIJLcn3H3F4PFH5tZ22B9W+CT1IQoIiKJiOUqGgOmA2+7+8PlVhUBo4Pp0cC85IcnIiKJahDDNt8Avg28ZWYlwbIfAw8Ac8xsDLAZGJaSCEVEJCFVJnh3fw2wSlb3S244IiKSLLqTVUQkpJTgRURCSgleRCSklOBFREJKCV5EJKSU4EVEQkoJXkQkpJTgRURCSgleRCSklOBFREJKCV5EJKSU4EVEQkoJXkQkpJTgRURCSgleRCSklOBFREJKCV5EJKSU4EVEQkoJXkQkpJTgRURCSgleRCSklOBFREJKCV5EJKSU4EVEQqrKBG9mT5jZJ2a2ttyyFma2yMw2BD+bpzZMEZG6LzMzEzOL+5WoWHrwTwLXnbDsbmCxu+cAi4N5ERE5hU2bNuHucb8SVWWCd/dlwM4TFg8CCoPpQmBwwhGIiEhKJDoG38bdtwXT/wDaVLahmY01s2IzK96+fXuCpxMRkXhV+0tWj3x+qPQzhLtPc/cCdy9o1apVdU8nIiIxSjTBf2xmbQGCn58kLyQREUmGRBN8ETA6mB4NzEtOOCIikiyxXCY5E/gL0MHMtpjZGOABoL+ZbQCuDuZFRKQWaVDVBu4+spJV/ZIci4iIJJHuZBURCakqe/DJ9NZbbyV0V1ZmZmYKohERCbcaTfCHDh2q1l1ZIiISOw3RiIiElBK8iEhIKcGLiISUEryISEgpwYuIhJQSvIhISCnBi4iElBK8iEhIKcGLiISUEryISEgpwYuIhJQSvIhISCnBi4iElBK8iEhIKcGLiISUEryISEgpwYuIhJQSvIhISCnBi4iElBK8iEhIKcGLiISUEryISEhVK8Gb2XVm9jcz22hmdycrKBERqb6EE7yZ1Qd+A/wTkAuMNLPcZAUmIiLVU50efHdgo7u/5+6HgFnAoOSEJSIi1dWgGvtmAB+Wm98C9DhxIzMbC4wNZr8ws7XVOGeYnAd8mu4gagm1xTFqi2PUFsd0SGSn6iT4mLj7NGAagJkVu3tBqs9ZF6gtjlFbHKO2OEZtcYyZFSeyX3WGaLYC55ebbxcsExGRWqA6Cf6vQI6ZXWhmDYERQFFywhIRkepKeIjG3Q+b2f8BXgbqA0+4+7oqdpuW6PlCSG1xjNriGLXFMWqLYxJqC3P3ZAciIiK1gO5kFREJKSV4EZGQSkmCr6qEgZmdaWazg/UrzSwrFXGkWwzt8EMzW29mpWa22Mwy0xFnTYi1rIWZ3WhmbmahvTwulrYws2HB78Y6M3u2pmOsKTH8H7nAzF41s9XB/5MB6YizJpjZE2b2SWX3ClnEY0FblZpZfpUHdfekvoh84fp34CKgIbAGyD1hm+8Bvw2mRwCzkx1Hul8xtsNVQONg+rthbIdY2yLY7hxgGfA6UJDuuNP4e5EDrAaaB/Ot0x13GttiGvDdYDoX2JTuuFPYHr2BfGBtJesHAC8BBvQEVlZ1zFT04GMpYTAIKAymnwf6mZmlIJZ0qrId3P1Vdz8QzL5O5F6CMIq1rMXPgMnAwZoMrobF0hb/BvzG3XcBuPsnNRxjTYmlLRxoGkyfC3xUg/HVKHdfBuw8xSaDgKc84nWgmZm1PdUxU5HgKyphkFHZNu5+GPgMaJmCWNIplnYobwyRv85hVGVbBB83z3f3P9ZkYGkQy+/FxcDFZrbCzF43s+tqLLqaFUtb3Ad8y8y2APOB/6iZ0GqleHNK6ksVSNXM7FtAAXBlumNJBzOrBzwM3JLmUGqLBkSGafoQ+VS3zMwucffd6QwqTUYCT7r7L83scmCGmXV296PpDqwuSEUPPpYSBtFtzKwBkY9eO1IQSzrFVMrBzK4G7gEGuvsXNRRbTauqLc4BOgNLzWwTkfHFopB+0RrL78UWoMjdv3T394F3iST8sImlLcYAcwDc/S9AIyJFyE5HcZeHSUWCj6WEQREwOpgeCizx4FuEEKmyHczsMuC/iCT3sI6zQhVt4e6fuft57p7l7llEvo8Y6O4JFViq5WL5/zGXSO8dMzuPyJDNezUYY02JpS0+APoBmFlHIgl+e41GWXsUATcHV9P0BD5z922n2iHpQzReSQkDM/spUOzuRcB0Ih+1NhL5UmFEsuNItxjb4SGgCfBc8B3zB+4+MG1Bp0iMbXFaiLEtXgauMbP1wBFgnLuH7RNurG1xJ/D/zewOIl+43hLCziAAZjaTyB/284LvHH4CnAHg7r8l8h3EAGAjcAD4TpXHDGlbiYic9nQnq4hISCnBi4iElBK8iEhIKcGLiISUEryISEgpwctJzOyImZWY2Voze87MGsex7y1m9us4z7evkuU/DW4Ew8yWlt34ZGbzzaxZ8PpePOeqIo6HguqND1Uz7n83s5uTFVcqVBa7hItKFUhFPnf3LgBm9gzw70RKCRAsaxDUEEopd59QyfIBQRxZRCqT/meSTjkWaOHuR6pzkOCa5RpXU/8uUneoBy9VWQ5km1kfM1tuZkXAejNrZGb/bWZvBbW6ryq3z/lBj3uDmf2kbKGZzTWzN4Je8tjyJzGzR4Lli82sVbDsSTMbemJAZrYpuMPzAaB98GnjITN7yswGl9vuGTMbdMK+Fmy7Noh9eLC8iMhNZ2+ULSu3T5Ny77XUzG4st+7nZrYmKArWJlh2n5ndFUx/347V/J8VLGtpZguD9/s7M9tsZueZWZaVqwVuZneZ2X3B9L+Z2V+Dc71Q9qkqaKPfmtlK4EEza29mC4J2Xm5mXwu2u9DM/hK8h0mn+geXEEl3DWS9at8L2Bf8bADMI1Krvg+wH7gwWHcnkTsPAb5G5JbyRkQKhm0jUh30LGAtQW13Ir1jyi1vGcw7MCqYngD8Oph+EhgaTC8td5xNROqRZFGudjaRYm1zg+lzgfeBBie8txuBRUTunGwTxN22/PuuoD0mA1PLzTcvF/c/B9MPAv8vmL4PuCuY/gg4M5huFvx8DJgQTH8zOE5F7+cu4L5gumW55ZOA/yjXRn8A6gfzi4GcYLoHkTIgENzmHkzfVtl71StcL/XgpSJnmVkJUEwkAU4Plq/ySPErgCuApwHc/R1gM5GaKQCL3H2Hu38OvBhsC/B9M1tDpNbM+RwroHUUmB1MP11u+7i4+5+I1DZpRaQK4Qt+8pDFFcBMdz/i7h8DfwK6VXHoq4HflDvPrmDyEJHkCvAGkQR9olLgGYtUDC2LpTfH2u6PwK4K9jtR56BH/hYwCuhUbt1z7n7EzJoAXydS+qKESJ2jsnrh3wBmBtMzYjifhIDG4KUi0TH4MhaplbM/xv1PrH/hZtaHSKK83N0PmNlSIj3+WPaPx1PAt4jUN6qyVkc1fenuZbEeoeL/T98kktD/GbjHzC45xfEOc/ywafn2eRIY7O5rzOwWgmJkgbJ/l3rA7hP/7cpRXZLTjHrwkqjlRHqSmNnFwAXA34J1/c2shZmdBQwGVhAZMtkVJPevESkJXKYekaqiAP8CvBZjDHuJlBou70ngdgB3X19J3MPNrH7Q0+8NrKriPIuIDGsAYGbNYwnOInXuz3f3V4HxRNqgCZHHEv5LsM0/AWXH+xhoHYzRnwlcX+5w5wDbzOwMgnY/kbvvAd43s5uCY5uZXRqsXsGxon4V7i/howQvifpPoF4wZDCbSJW/snr2q4AXiAxPvOCRsr8LgAZm9jaRL0dfL3es/UD34AvGvsBPYwnAIxUWVwRfmD4ULPsYeBv470p2+30Q1xpgCfAjd/9HFaeaBDQPzrOGyLN0Y1EfeDpoo9XAYx55aMdEoLeZrQNuIDIMhrt/SeS9ryLyR+Wdcse6F1hJJFGXX36iUcCYIM51HHsE3g+A24JYTvkUIAkPVZOUUAmuLnkLyHf3z9IdTyws8pCTAnf/NN2xSLioBy+hYZGbot4GflVXkrtIKqkHLyISUurBi4iElBK8iEhIKcGLiISUEryISEgpwYuIhNT/ArEVEwhEG9ZnAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the chisquared probability \n",
    "probAx = plt.subplot(111)\n",
    "probContents, probBins, _ = probAx.hist(chiProb, np.linspace(0.,1.,25), histtype='step', color='black')\n",
    "statBox(probAx, chiProb, probBins)\n",
    "probAx.set_xlim(probBins[0], probBins[-1])\n",
    "_ = probAx.set_xlabel(\"Probability of chisquared\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c604939d",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.9.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}