ineff13.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/afs/cern.ch/work/i/innocent/public/code/s7/CMSSW_9_3_0_pre2\n",
      "/afs/cern.ch/user/i/innocent/pyTools\n"
     ]
    }
   ],
   "source": [
    "%%bash\n",
    "echo $CMSSW_BASE\n",
    "echo $PWD"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Linux vinavx3.cern.ch 3.10.0-514.10.2.el7.x86_64 #1 SMP Fri Mar 3 00:04:05 UTC 2017 x86_64 x86_64 x86_64 GNU/Linux\n",
      "innocent\n",
      "root://eoscms.cern.ch///store/relval/CMSSW_9_3_0_pre3/RelValSingleMuPt100/GEN-SIM-RECO/92X_upgrade2017_realistic_v10_resub-v1/00000/3C1A3236-AE7C-E711-ACA8-0025905A60FE.root\n",
      "root://eoscms.cern.ch///store/relval/CMSSW_9_3_0_pre3/RelValSingleMuPt100/GEN-SIM-RECO/92X_upgrade2017_realistic_v10_resub-v1/00000/3C1A3236-AE7C-E711-ACA8-0025905A60FE.root (1 runs, 60 lumis, 6000 events, 1208929930 bytes)\n"
     ]
    }
   ],
   "source": [
    "%%bash\n",
    "uname -a\n",
    "whoami\n",
    "\n",
    "edmFileUtil root://eoscms.cern.ch///store/relval/CMSSW_9_3_0_pre3/RelValSingleMuPt100/GEN-SIM-RECO/92X_upgrade2017_realistic_v10_resub-v1/00000/3C1A3236-AE7C-E711-ACA8-0025905A60FE.root\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "import math\n",
    "import time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "application/javascript": [
       "IPython.OutputArea.prototype._should_scroll = function(lines) {\n",
       "    return false;\n",
       "}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%javascript\n",
    "IPython.OutputArea.prototype._should_scroll = function(lines) {\n",
    "    return false;\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Welcome to JupyROOT 6.10/01\n"
     ]
    }
   ],
   "source": [
    "from DataFormats.FWLite import Handle, Events"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "singleMU100 = \"root://eoscms.cern.ch///store/relval/CMSSW_9_3_0_pre3/RelValSingleMuPt100/GEN-SIM-RECO/92X_upgrade2017_realistic_v10_resub-v1/00000/3C1A3236-AE7C-E711-ACA8-0025905A60FE.root\"\n",
    "\n",
    "label = \"generalTracks\"\n",
    "quality = \"highPurity\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def nt(fevents) :\n",
    "    events = Events(fevents)\n",
    "    tracks = Handle(\"std::vector<reco::Track>\")\n",
    "    nt = []\n",
    "    for i in range(0,1000) : # events.size()):\n",
    "      a= events.to(i)\n",
    "      if (i%500==0) : print \"Event\", i \n",
    "      a=events.getByLabel(label, tracks)\n",
    "      nt.append(tracks.product().size())\n",
    "    print len(nt)\n",
    "    return nt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Event 0\n",
      "Event 500\n",
      "1000\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEZCAYAAABmTgnDAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAF1JJREFUeJzt3Xm0bGV95vHvw6gIiKgRAbmiLCTSjWi7CA5pammciJFE\n0QBRAQdo7WiykjYkxjQH7Da26c5yQE2MSqIJgh0HUDDBFg+KgHFCFEGuLgREwIF5UBF//cd5r1aK\nc++uM9SpXdzvZ61at+rdVXv/ap9966m933fvSlUhSdKmbDHtAiRJ/WdYSJI6GRaSpE6GhSSpk2Eh\nSepkWEiSOhkW6o0kJyc5ccrLvyHJhVNa9oree5IjkvzL6lUl/ZJhoY1K8p0k1ye571DbS5N8erqV\nrb4kTwaeCuxaVQcuMv3IJJ+dTnXjqapTquqZ065jKZL8PMkjpl2HuhkW2pRq28gfLtLea0mWum0/\nHPhOVf14Y7Pset/LWKZmYFvSAjdudfkr4I+T7Dg6Icm69s1wi6G2Tyd5Cb/8Nn5ekr9OcmOS9Ume\nmOSoJFcluS7Ji0dm++AkZye5pc1rj6F579Om/SjJpUmePzTt5CTvSHJmkluBwSL1PjTJ6e31lyd5\nWWt/CfB3wBPaco8fed0+wDvb9FuT3LCxZSY5OMmXk9yc5MpF5vXkJJ9r6+PKRd4/SXZIck6SN7fH\nBye5pNV2dZI/WuwPNbr30/42x7b3ekOSkxb9Cy88N0n+NMm3kvwgyalJdmrTzkryypHnX5Tkt8f8\nu5yU5OOt/guS7NmmndtC+OI27fn3rEy9UVXevC16A64AngL8M/D61vZS4Jx2fx1wN7DF0Gs+Dbyk\n3T8S+Cnw4vah8HrgauBtwNbA04BbgO3a808Gbgae1Ka/Gfhsm7YdcNXQvB4D/ADYZ+i1NwIHtsfb\nLPJ+PjO07McA3wcGQ7V+ZhPr4h7TF1sm8J+Bfdvj/wBcCzxnaH3dArwA2BJ4ALDf0LxOBHYGPg+c\nMLSc7wFPbPfvD+w/To3Az4EzgB2Ah7X3+/SNvPYPgPOBh7b1807glDbtRcB5Q899NHADsNWYf5cf\nAP+pfTn9xw3zHapxz2lv6966b+5ZaBzHA7+f5IHLeO0VVfW+WvhkOA3YtX0Q3lVVn2xhstfQ88+s\nqs9V1V3AnwMHJtkNePbwvKrqq8CHgOFvo6dX1YUsfAn66XARSXYHngAc15b9VeDd7UNuJf7dMqvq\nM1V1SXv8deBU4KD23MOBT1bVB6vq7qq6saouHprXbsC5wGlVNbxH8lNg3yQ7VNXNVXXREur7y6q6\ntaqubkG+/0aedyzw51V1bVv3JwKHtr3GjwCPSfKw9twjgA9X1c/G/Lt8pKq+VFU/B/5pkRqyhPej\nKTEs1Kl9+H0c+LNlvPz6oft3tvn9cKRt+6HHVw8t9/b2zX3X9q38wHY45YYkN7YPrYcs9tpF7Arc\nUFV3DLVd2T6gV+LfLTPJAe0Q0veT3NQ+hB/UJj8M+PYm5vWbwH2Avx1pf16bdmU7NHePDvhNGF7/\nd4ys62HrgI9sWL/AN4C7gIdU1W3AWcBh7bmHtz0Exvy7XDdmDeoxw0LjmgNePvLhenv7d7uhtl1W\nuJwN315Jsn07VPO99qE8X1U7t9sDqmrHqvr9odduqrP0e8DOSe431LYHcM2YdW1s3qPtpwAfBXar\nqp3aB/+Gb85Xj+xFjXoX8C/AJ4ZHoLVv5b8NPBg4HfjgmDUvxVXAs0bW7/2q6to2/QPAES2otq2q\n+aH31PV30b2AYaGxVNW322GkVw+1/bB92L4wyRato/iRHbPqOuRwcOsE36b1cVxYVde0PZu9k7ww\nyVZJtk7y+CSPGrP+77Zj8n+ZZNsk+7X+l/ePtQIWvqHvnmTrjudtD9xYVXclOaB9y97gn4CnJjk0\nyZZJdk7ymJE6XwV8E/h4kvu093lEkh2r6m7g1tZPtNr+FnjDhgEFSR6c5DlD089qexEntu1ggxX9\nXdpeh0NnZ0DvwqKNrHhnktOSvHTa9WzmRr81n9j2IobbXw78CfBD4FeBzy1xnjVy/5S2F/Mj4LHA\nC1n4EL0NeHo7FPK9dnsjsO0S3s/hwJ7ttR8C/qKqxj1n5BzgEuC6JN/fxPNeCbw+yc3A64Y/WFu/\nwcHAf2sdxF8B9ltkHscA3217KNu2DuYr2mGtY0YCaFM2ta5HvaXttZzdaj8fOGCo9p8CH27nopwy\n1L7Sv8sc8L52COvQMV+jKUj19MePkgQ4tap+d9q1SNLmbuJ7Fkne084Cvnik/ZlJLmtjwI8bmfZb\nwJltJIkkacomvmfRLqNwG/C+qtqPX57pennbpf0e8AXgsKq6bOS1p1fVIRMtUJLUaatJL6Cqzkuy\nbqT5AGB9VV3JQiicChwCXJbkIOC5bQjhve4aRJI0iyYeFhux28j49O9u6EyrqnPbiUmSpJ6YVlis\nSJJ+9spLUs9V1bLOmJ/W0Nlr2glRG+y+hJOjwGtarert+OOPn3oN95bbkUcWhxwy3fV5+eXFXntN\nf12s1s3tc/VuK7FWYZGRk7G+AOzVrlq6TRujfcYa1SJJWqK1GDp7SjvBZ+92Weqj25morwLObic6\nnVpVly5lvnNzc8zPz4/xTEnavM3PzzM3N7eiefT2pLxNSVKzWHdfzc/PMxjc4+cftAxHHQW77DLP\nG984vfW5fj0cfPDCv/cGbp+rJwmz1mehHvE/4uraZx/X52py++wHw0KS1Glmw8I+C0kaj30WUs8c\ndRQMBgv/Tsu9rc9Cq8c+C0nSRBkWkqROMxsW9llI0njss5B6xj4L9Zl9FpKkiTIsJEmdDAtJUqeZ\nDQs7uCVpPHZwSz1jB7f6zA5uSdJEGRaSpE6GhSSpk2EhSeo0s2HhaChJGo+joaSecTSU+szRUJKk\niTIsJEmdDAtJUifDQpLUybCQJHWa2bBw6Kwkjcehs1LPOHRWfebQWUnSRBkWkqROhoUkqZNhIUnq\nZFhIkjoZFpKkToaFJKmTYSFJ6jSzYeEZ3JI0Hs/glnrGM7jVZ57BLUmaKMNCktTJsJAkdTIsJEmd\nDAtJUifDQpLUybCQJHUyLCRJnQwLSVInw0KS1MmwkCR1mtmw8EKCkjQeLyQo9YwXElSfeSFBSdJE\nGRaSpE6GhSSpk2EhSepkWEiSOhkWkqROhoUkqZNhIUnqZFhIkjoZFpKkToaFJKmTYSFJ6mRYSJI6\nGRaSpE6GhSSp01bTLmAxSQ4BfhPYAXhvVX1y2jVJ0uasl2FRVacDpyfZCfgrwLCQpClak8NQSd6T\n5PokF4+0PzPJZUkuT3LcIi99HfD2tahRkrRxa9VncTLwjOGGJFsAJ7X2fYHDk+wzNP2NwFlVddEa\n1ShJ2og1CYuqOg+4caT5AGB9VV1ZVXcBpwKHsBAUrwKeChya5Ji1qFGStHHT7LPYDbh66PF3W4BQ\nVW8D3rapF8/Nzf3i/mAwYDAYTLBUSZo98/PzzM/Pr8q8etnBPY7hsJAk3dPoF+kTTjhh2fOa5nkW\n1wB7DD3evbVJknpmLcMi7bbBF4C9kqxLsg1wGHDGuDObm5tbtd0rSbo3m5+fX/HRmFTVqhW00YUk\npwAD4IHA9cDxVXVykmcBb26h9Z6qeuOY86u1qFtaqqOOgsFg4d9pWb8eDj544V9pWBKqKmM89R7W\npM+iqo7YSPsngE+sRQ2SpOXz2lCSpE4zGxb2WUjSeGamz2K12WehvrLPQn22kj6Lmd2zkCStHcNC\nktTJsJAkdZrZsLCDW5LGYwe31DN2cKvP7OCWJE2UYSFJ6mRYSJI6zWxY2MEtSeOxg1vqGTu41WcT\n7+BO8sgk27b7gySvTrLTchYoSZo94x6G+hBwd5K9gHcBDwNOmXBtkqSeGDcsfl5VPwN+B3hbVb0G\neOiEa5Mk9cS4YXFXksOBI4GPt7atJ1hXJzu4JWk8a9bBneTRwH8BLqiqDyTZE3hBVf2vFS19mezg\nVl/Zwa0+W4ufVX1aVb16w4OquiLJj5ezQEnS7Bn3MNSRi7RN8buTJGktbXLPovVTHAHsmeSMoUk7\nADdMvjxJUh90HYY6H7gWeBDwf4babwUunnBtkqSe2GRYVNWVwJXAE9auJElS34x7Bvdzk6xPcnOS\nW5LcmuSWyZcnSeqDcUdDvQn4raq6dML1jG1ubo7BYMBgMJh2KZLUa/Pz8ys+L23c8yw+V1VPWtGS\nVpHnWaivPM9CfbYW51l8MclpwEeBn2xorKoPL2ehkqTZMm5Y7AjcATx9qK0Aw0KSNgNjhUVVHT35\nUiRJfTXuaKi9k3wqydfb4/2SvG7i1UmSemHcy338HfBnwF0s7GlcDBw22dIkSX0xblhsV1X/NtL2\nswnUI0nqoXHD4odJHtk6tUlyaLsMiCRpMzDuaKj/2n5OdZ8k1wBXAL834do2yZPyJGk8a3lS3pZV\ndXeS+wFbVNWtK1rqCnlSnvrKk/LUZys5KW/cw1BXJHkXcCBw23IWJEmaXeOGxT7A/2uHo65IclKS\nJ0+4NklST4wVFlV1R1V9sKqeCzy2ndF97uTLkyT1wbh7FiQ5KMk7gC8B9wFeMNnSJEl9MdZoqCTf\nAb4CfBB4TVXdPvnSJEl9Me7Q2f2qyh87kqTN1LiHoXbx2lCStPny2lCSpE5eG0qS1MlrQ0mSOs3s\ntaEkSWtn3LA4oqp+o10bass+jIzyQoKSNJ6JX0gwyXHAZ4B3VtX+re3LVfW4FS11hbyQoPrKCwmq\nz1ZyIcGuPYvLgOcDj0jy2fb4gUkeVVXfXF65kqRZ09XBfRPwWuBbwAB4S2v/0yTnr0F9kqQe6Nqz\neAbw34FHAn8NXAzcXlVHr1F9kqQe2OSeRVW9tqqeCnwHeD+wJfDgJOcl+djalSlJmqZxR0P9a1V9\nEfhikldU1ZOTPGjCtUmSemLc37P4k6GHR7W2H06sKklSr4z9exYbVNVXJ1OKJKmvlhwWkqTNj2Eh\nSepkWEiSOhkWkqROhoUkqZNhIUnqZFhIkjoZFpKkTr0LiyR7Jnl3kg9OuxZJ0oLehUVVXVFVL5t2\nHZKkX5p4WCR5T5Lrk1w80v7MJJclubz9Ip8kqafWYs/i5Pa7GL+QZAvgpNa+L3B4kn1GXresn/6T\nJK2+iYdFVZ0H3DjSfACwvqqurKq7gFOBQ1gIkp2TvBPY3z0OSeqHcX/PYrXtBlw99Pi7LUCoqhuA\nV3TNYG5u7hf3B4MBg8FgQqVK0myan59nfn5+VeY1rbBYseGwkCTd0+gX6RNOOGHZ85rWaKhrgD2G\nHu/e2iRJPbRWYZGRDusvAHslWZdkG+Aw4Iw1qkWStERrMXT2FOB8YO8kVyU5uqruBl4FnA1cApxa\nVZcuZb5zc3OrdixOku7N5ufnV3zoPlW1agWtlSQ1i3Xr3u+oo2AwWPh3Wtavh4MPXvhXGpaEqlrW\naQm9O4NbktQ/hoUkqdPMhoV9FpI0HvsspJ6xz0J9Zp+FJGmiDAtJUqeZDQv7LCRpPPZZSD1jn4X6\nzD4LSdJEGRaSpE6GhSSp08yGhR3ckjQeO7ilnrGDW31mB7ckaaIMC0lSJ8NCktTJsJAkdZrZsHA0\nlCSNx9FQUs84Gkp95mgoSdJEGRaSpE6GhSSpk2EhSeq01bQLWK4TT5x2BfCyl8Guu067CkmavJkN\ni099ao516wY8/OGDqSz/k5+E+94XXvOaqSxeksY2Pz+/4lMNZjYszj13ZWOGV+rOO6e6eEka22Aw\nYDAYcMIJJyx7HvZZSJI6GRaSpE6GhSSpk2EhSepkWEiSOhkWkqROMzt0VpJmxU03wbe/Pe0qVmZm\nw2Jubu4XY4clqc+OPRa++EXYaafpLP/WW+f50Y8205PyVvpDHpK0Vm6/Hd7yFnj2s6dVwYALLhjw\nxCd6Up4kaYIMC0lSJ8NCktTJsJAkdTIsJEmdDAtJUifDQpLUybCQJHUyLCRJnQwLSVInw0KS1Glm\nw2Jubo75+ZVdGEuSNgfz8/O8+90ru56eFxKUpHu5wWDAttsOeO97vZCgJGmCDAtJUifDQpLUybCQ\nJHUyLCRJnQwLSVInw0KS1MmwkCR1MiwkSZ0MC0lSJ8NCktTJsJAkdTIsJEmdDAtJUifDQpLUybCQ\nJHXq3Y8fJdkOeAfwE+Dcqjpl2jVJ0uauj3sWzwX+b1UdCzxn2sVsDvx52tV12WWuz9Xk9tkPEw+L\nJO9Jcn2Si0fan5nksiSXJzluaNLuwNXt/t2Trk/+Z1xthsXqcvvsh7XYszgZeMZwQ5ItgJNa+77A\n4Un2aZOvboEBkDWoT5LUYeJhUVXnATeONB8ArK+qK6vqLuBU4JA27SPAoUneDnxs0vVJkrqlqia/\nkGQd8LGq2q89fh7wjKo6pj1+IXBAVb16zPlNvmhJuheqqmUdsendaKhxLPfNSpKWZ1qjoa4B9hh6\nvHtrkyT10FqFRUY6q78A7JVkXZJtgMOAM9aoFknSEq3F0NlTgPOBvZNcleToqrobeBVwNnAJcGpV\nXbrIaxcddjvynLcmWZ/koiT7T/r9zKqudZnkoCQ3Jflyu71u7aucHUl2T3JOkkuSfC3Jov1tbp/d\nxlmXbp/jS7Jtks8n+Upbp2/YyPOWtm1WVW9vwJOB/YGLNzL9WcCZ7f6vARdOu+a+3sZYlwcBZ0y7\nzlm5AbsA+7f72wPfBPYZeY7b5+qtS7fPpa3T7dq/WwIXAk8amb7kbbOPZ3D/wkaG3Q47BHhfe+7n\ngfsnecjaVTg7xliXeF7L+Krquqq6qN2/DbgU2G3kaW6fYxhzXeL2Ob6quqPd3bYdQRr9v7/kbbPX\nYTGG3YbO9qZ1ki+2kWk8T2i7pGcmefS0i5kVSR7e9to+PzLJ7XOJNrEucfscX5ItknwFuA6Yr6pv\njDxlydvmTA6d1UR8Cdijqu5I8izgo8De0y6q75JsD/wz8AftW7GWqWNdun0uQVX9HHhskh2Bs5Mc\nVFXnrmSes75ncQ3wsKHHDsFdpqq6bcOua1V9Atg6yc7TrqvPkmzVPtzeX1WnL/IUt88xda1Lt8/l\nqapbgDOBx49MWvK2OQthMTrsdtgZwItZ2NgOBG6qquvXtryZstF1OXy8MskB7ez+G9a0utnzXuAb\nVfWWjUx3+xzfJtel2+f4kjwoyf3b/fsCTwMuGnnakrfNXh+GasNuB8ADk1wFHA9ssxCY9a6qOivJ\nwUm+BdwOHD3tmvuqa12263G9ArgLuBP43WnX3GdJngT8HvC1dmy4gNcC69w+l2acden2uSQPBf4h\nSdoOwfur6lNJjl3Jtrkm14aSJM22WTgMJUmaMsNCktTJsJAkdTIsJEmdDAtJUifDQpLUybCQRiR5\nQ7sk9iFJjmttRybZZZXmvy7J11ZjXtJaMSyke/q1diG7g4DPtLajNnahtSTL+X/kCU6aKb0+g1ta\nS0neBDwDeHj7wa69gKck+VC7ts4/JrkTeGK7jPZpwG8Ab2oXbDsG2Br4FvCiqvpxkl8B/gZ4RAuI\nVwDXDi3zEe2aSC9vZyaf3OaxBfC8qvr2dNeKtMAzuKUhSR4PvAj4o3Zp519v7ecAf1xVX2mPrwDe\nXlX/uz1+QFXd2O6/Hriuqt6e5FTg/Kp6a7v8wvbAzsDHgEOBU4EXV9XXk7wVuKCqPtAurLdlVf1k\nyqtEAvcspHt4HHAx8KvAZUPti12E8bSh+/8xyf8AdgLuB/xra39KCx9q4ZvZre1qqb/SLrP93Kra\nsJwLgD9Psjvwkar61mTfqjQ+w0Ja2Bt4DPD37VLNP2gf+CT5cjvstJjbh+7/PfCctodwZOvvYBN9\nEzcDVwG/viGU2h7FhcCzgbOSHFNV86v6RqVlsoNbWvig/mpVPRb4ZlU9GjgHeHpVPa6qfgzcAuy4\niVlsD1yXZOt2BdUNPgW8kl/+etmGefwE+B3gxUkOb9P3rKorquptwOnAfhN909ISuGchNUkeNPRb\nxY+qqm8OTf4H4G+S3NH2NEb3GP4C+Dfg+20k1Q6t/Q+BdyV5KfCz1sF9HQsBdWeSZ7dfMrsV2DfJ\ni9pluK8F/ufavHOpmx3ckqROHoaSJHUyLCRJnQwLSVInw0KS1MmwkCR1MiwkSZ0MC0lSJ8NCktTp\n/wPLdbM8r8lF0gAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f80ac1fcf50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning in <TClass::Init>: no dictionary for class edm::Wrapper<edm::DetSetVector<CTPPSPixelCluster> > is available\n",
      "Warning in <TClass::Init>: no dictionary for class edm::DetSetVector<CTPPSPixelCluster> is available\n"
     ]
    }
   ],
   "source": [
    "nt = nt(singleMU100)\n",
    "fig, ax = plt.subplots()\n",
    "ax.hist(nt,histtype='step', fill=False)\n",
    "ax.set_title('Number of tracks in event')\n",
    "ax.set_xlabel('#tracks')\n",
    "ax.set_ylabel('#events')\n",
    "plt.yscale('log', nonposy='clip')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAEACAYAAABYq7oeAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFsVJREFUeJzt3X+QXXWZ5/H304ZfSqDDjw5WEtMMOBBYdqPjBkdE2g3D\nDxeBWqoC2WUx4vIPKpTsuiRQkvDHyiBFTWQo1kJYEkYtJrJTCluOBopprVnHwCAIQgYiVBLSS5IB\nCWEKS4E8+0dO4iW325Duc/ve8+X9qrqVe799zu3z6dPpJ+d5bt9EZiJJUqu+bh+AJKn3WBwkSW0s\nDpKkNhYHSVIbi4MkqY3FQZLUZq/FISLujIgtEfFEy9q0iFgdEc9ExI8i4tCWjy2JiHURsTYizmhZ\n/3BEPBERz0bE8g7lkSTV4J1cOdwFnLnH2mLgwcw8DngIWMLOAnACsACYA5wN3BYRUe3zP4HPZeYf\nA38cEXs+pySpR+y1OGTm3wOv7LF8HrCyur8SOL+6fy5wT2a+mZnrgXXAvIg4CpiamY9U293dso8k\nqceMd+YwkJlb2Fk8NgMD1foM4IWW7UaqtRnAppb1TdWaJKkH1TWQ9j04JKkgU8a535aImJ6ZW6qW\n0dZqfQSY1bLdzGptrPVRRYTFRpLGITPjHWy2V+/0yiGq2y73AYuq+58Bvt+yflFE7B8RRwPHAg9X\nradXI2JeNaC+pGWfMSSZZd6WLl3a9WMwm/nMV96tTnu9coiI7wBDwOERsRFYCvw58N2IuBTYUL1C\nicx8OiJWAU8DbwCX5++P+PPACuBA4AeZ+cNakzTI+vXru30IHVNyNszXeKXnq9Nei0Nm/scxPnT6\nGNvfANwwyvqjwEnjOUhJ0uTyN6S7YNGiRe9gq2YqORvma7zS89Up6u5T1WHnQDrpwUOTpJ4VEUz2\nQFo1Gh4e7vYhdEzJ2TBfm8HBQSLC2yTfBgcHO/Y9sMt4X8oqSWzYsKH2V8lo737/rkQd/By9eGJt\nK0nNULUxun0Y7zpjfd1tK0mSOsri0AUl961Lzob59C5icZCkSbJy5UpOPfXUbh/GO2Jx6IKhoaFu\nH0LHlJwN8zXO4OAg06dP5ze/+c3utTvvvJNPfvKTo26/YcMG+vr62LFjR8eOaTKGyXWwOEgqVkSw\nY8cOli9f3rY+mszc65D9rbfeqv04e5HFoQtK7uuWnA3zNdKXv/xlbr75ZrZv377XbU877TQA+vv7\nOeSQQ1izZg0rV67k4x//OFdddRVHHHEE119/Pc8//zzz58/niCOOYGBggIsvvvhtz79p0yYuuOAC\nBgYGOPLII7niiivGPLZPfOITvPbaazz33HMMDQ3R39/PwMAACxcurPGrsO8sDpKK9pGPfIShoSFu\nuummvW77k5/8BIDt27ezfft2Tj75ZADWrFnDsccey9atW7n22mvJTK655ho2b97M2rVr2bRpE8uW\nLQNgx44dnHPOORx99NFs3LiRkZERLrroord9nszksssu45e//CUPPPAAU6dO5Stf+Qpnnnkm27Zt\nY9OmTXzxi1/syNfjnbI4dEFpfd1WJWfDfOMSUc9tIq6//npuvfVWXn755Xe0/Z5tpRkzZnD55ZfT\n19fHAQccwDHHHMP8+fOZMmUKhx9+OF/60pf48Y9/DFUhefHFF/na177GgQceyP7778/HPvax3c/1\nu9/9joULF7Jt2zbuv/9+DjjgAAD2228/NmzYwMjISNs+3WBxkNRRmfXcJuLEE0/knHPO4YYb2t4w\n+h2ZNWvW2x5v3bqVhQsXMnPmTPr7+7n44ot56aWXoGopzZ49m76+0X+8/upXv+K+++5j6dKlTJny\n+zepuOmmm9ixYwfz5s3jpJNO4q677hrXsdbF4tAFJfZ1dyk5G+ZrtGXLlvHNb36TkZEx/xPKMQfV\ne65fc8019PX18dRTT7Ft2za+9a1v7b7amDVrFhs3bhzzFU8nnHACd911F2eddRbPPvvs7vWBgQFu\nv/12RkZG+MY3vsHll1/O888/P860E2dxkPSucMwxx3DhhRdyyy23jLnNkUceSV9fH88999wffK7X\nXnuNgw8+mKlTpzIyMvK2eca8efN4//vfz+LFi3n99df57W9/y09/+tO37X/hhRfy1a9+ldNPP313\nAbj33nt3F67+/n76+vrGvPqYDBaHLii5b11yNszXOHv+i/+6667j9ddfH/MK4aCDDuLaa6/llFNO\n4bDDDuPhhx8edbulS5fy6KOP0t/fz6c//WkuuOCC3R/r6+vj/vvvZ926dXzgAx9g1qxZrFq1qu05\nLrnkEq677jrmz5/Pxo0beeSRRzj55JM55JBDOP/887nlllsm5d1Xx+Ib70kaN994rzt8471CldzX\nLTkb5tO7iMVBktTGtpKkcbOt1B22lSRJXWFx6IKS+7olZ8N8ehexOEiS2jhzkDRug4ODbNiwoduH\n8a4ze/Zs1q9f37Ze58zB4iBJhXAg3XAl93VLzob5Gq/0fHWyOEiS2thWkqRC2FaSJHWUxaELSu57\nlpwN8zVe6fnqZHGQJLVx5iBJhXDmIEnqKItDF5Tc9yw5G+ZrvNLz1cniIElqM6GZQ0QsAS4G3gKe\nBD4LvA/4a2A2sB5YkJmvtmx/KfAmcGVmrh7jeZ05SNI+6omZQ0TMBi4DPpSZ/xqYAiwEFgMPZuZx\nwEPAkmr7E4AFwBzgbOC2GOt/+ZYkddVE2krbgd8B74uIKcBBwAhwHrCy2mYlcH51/1zgnsx8MzPX\nA+uAeTVkaJyS+54lZ8N8jVd6vjqNuzhk5ivAzcDGqii8mpkPAtMzc0u1zWZgoNplBvBCy1OMVGuS\npB4zZbw7RsQfAV+qZguvAt+NiP8E7DkpGOfkYBHLlg0C0N/fz9y5cxkaGoKW6t/Ux7vWeuV46nw8\nNDTUU8djPvOVnG94eJgVK1ZA9X9r1GncA+mIWAD8WWZeVj3+z8BHgX8HDGXmlog4Cvi7zJwTEYt3\nXkzkjdX2PwSWZuaaUZ7bgbQk7aOeGEgDzwAfjYgDq8HyfOBp4D5gUbXNZ4DvV/fvAy6KiP0j4mjg\nWODhGjI0Tsl9z5KzYb7GKz1fncbdVsrMX0TE3cCj1UtZHwNuB6YCqyLiUmBD9QolMvPpiFhVFZA3\ngMuzF9+7Q5LkeytJUil6pa0kSSqUxaELSu57lpwN8zVe6fnqZHGQJLVx5iBJhXDmIEnqKItDF5Tc\n9yw5G+ZrvNLz1cniIElq48xBkgrhzEGS1FEWhy4oue9ZcjbM13il56uTxUGS1MaZgyQVwpmDJKmj\nLA5dUHLfs+RsmK/xSs9XJ4uDJKmNMwdJKoQzB0lSR1kcuqDkvmfJ2TBf45Wer04WB0lSG2cOklQI\nZw6SpI6yOHRByX3PkrNhvsYrPV+dLA6SpDbOHCSpEM4cJEkdZXHogpL7niVnw3yNV3q+OlkcJElt\nnDlIUiGcOUiSOsri0AUl9z1Lzob5Gq/0fHWyOEiS2jhzkKRCOHOQJHWUxaELSu57lpwN8zVe6fnq\nZHGQJLWZ0MwhIg4F7gD+FbADuBR4FvhrYDawHliQma9W2y+ptnkTuDIzV4/xvM4cJGkf9dLM4evA\nDzJzDvBvgH8CFgMPZuZxwEPAkuqgTwAWAHOAs4HbIqKWEJKkeo27OETEIcCpmXkXQGa+WV0hnAes\nrDZbCZxf3T8XuKfabj2wDphXS4qGKbnvWXI2zNd4peer00SuHI4GXoqIuyLi5xFxe0S8F5iemVvY\nWTA2AwPV9jOAF1r2H6nWJEk9Ztwzh4j4E+BnwJ9m5j9GxF8ArwFfyMzDWrZ7OTMPj4i/BP4hM79T\nrd9RtaT+ZpTnduYgSfuozpnDlAnsuwl4ITP/sXr8v6t5w5aImJ6ZWyLiKGBr9fERYFbL/jOrtTEs\nYtmyQQD6+/uZO3cuQ0ND0HJp6GMf+9jH7+bHw8PDrFixAoDBwZ0/L+sy0Vcr/Ri4LDOfjYilwHur\nD/06M2+MiKuBaZm5uBpIfxs4uWonPQB8MEc5gNKvHIaHh3ef6NKUnA3zNV7p+XrlygHgCuDbEbEf\n8DzwWeA9wKqIuBTYUL1Cicx8OiJWAU8DbwCXj1YYJEnd53srSVIheun3HCRJBbI4dEHJr7UuORvm\na7zS89XJ4iBJauPMQZIK4cxBktRRFocuKLnvWXI2zNd4peerk8VBktTGmYMkFcKZgySpoywOXVBy\n37PkbJiv8UrPVyeLgySpjTMHSSqEMwdJUkdZHLqg5L5nydkwX+OVnq9OFgdJUhtnDpJUCGcOkqSO\nsjh0Qcl9z5KzYb7GKz1fnSwOkqQ2zhwkqRDOHCRJHWVx6IKS+54lZ8N8jVd6vjpZHCRJbZw5SFIh\nnDlIkjrK4tAFJfc9S86G+Rqv9Hx1sjhIkto4c5CkQjhzkCR1lMWhC0rue5acDfM1Xun56mRxkCS1\nceYgSYVw5iBJ6iiLQxeU3PcsORvma7zS89XJ4iBJauPMQZIK0VMzh4joi4ifR8R91eNpEbE6Ip6J\niB9FxKEt2y6JiHURsTYizpjo55YkdUYdbaUrgadbHi8GHszM44CHgCXsLAwnAAuAOcDZwG0RUUuF\na5qS+54lZ8N8jVd6vjpNqDhExEzgU8AdLcvnASur+yuB86v75wL3ZOabmbkeWAfMm8jnlyR1xoRm\nDhHxXeB/AIcC/zUzz42IVzJzWss2v87MwyLiL4F/yMzvVOt3AD/IzL8Z5XmdOUjSPuqJmUNE/Htg\nS2Y+Dvyhg/FHvCQ1zJQJ7HsKcG5EfAo4CJgaEX8FbI6I6Zm5JSKOArZW248As1r2n1mtjWERy5YN\nAtDf38/cuXMZGhqClr5hUx8vX768qDytj1t7ur1wPOYzX8n5hoeHWbFiBQCDg4PUqZaXskbEaS1t\npa8BL2fmjRFxNTAtMxdXA+lvAycDM4AHgA/mKAdQeltpeHh494kuTcnZMF/jlZ6vzrZSJ4rDYcCq\n6iphA7AgM7dV2y0BPge8AVyZmavHeL6ii4MkdULPFYe6WRwkad/1xEBa41fya61Lzob5Gq/0fHWy\nOEiS2thWkqRC2FaSJHWUxaELSu57lpwN8zVe6fnqZHGQJLVx5iBJhXDmIEnqKItDF5Tc9yw5G+Zr\nvNLz1cniIElq48xBkgrhzEGS1FEWhy4oue9ZcjbM13il56uTxUGS1MaZgyQVwpmDJKmjLA5dUHLf\ns+RsmK/xSs9XJ4uDJKmNMwdJKoQzB0lSR1kcuqDkvmfJ2TBf45Wer04WB0lSG2cOklQIZw6SpI6y\nOHRByX3PkrNhvsYrPV+dLA6SpDbOHCSpEM4cJEkdZXHogpL7niVnw3yNV3q+OlkcJEltnDlIUiGc\nOUiSOsri0AUl9z1Lzob5Gq/0fHWyOEiS2jhzkKRCOHOQJHXUuItDRMyMiIci4qmIeDIirqjWp0XE\n6oh4JiJ+FBGHtuyzJCLWRcTaiDijrhBNU3Lfs+RsmK/xSs9Xp4lcObwJXJWZJwJ/Cnw+Io4HFgMP\nZuZxwEPAEnYWhhOABcAc4Gzgtoio5fJHklSv2mYOEfE94NbqdlpmbomIo4DhzDw+IhYDmZk3Vtv/\nLbAsM9eM8lzOHCRpH/XczCEiBoG5wM+A6Zm5hZ2VYDMwUG02A3ihZbeRak2S1GOmTPQJIuJg4F7g\nysz8l53/6n+bcf77fxHLlg0C0N/fz9y5cxkaGoKWvmFTHy9fvryoPK2PW3u6vXA85jNfyfmGh4dZ\nsWIFAIODg9RpQm2liJgC/B/gbzPz69XaWmCopa30d5k5Z5S20g+Bpe/GttLw8PDuE12akrNhvsYr\nPV+dbaWJFoe7gZcy86qWtRuBX2fmjRFxNTAtMxdXA+lvAydX7aQHgA/mKAdQenGQpE7oieIQEacA\nPwGerFpHCVwDPAysAmYBG4AFmbmt2mcJ8DngjaoNtXqM57Y4SNI+6omBdGb+38x8T2bOzcwPZeaH\nM/OHmfnrzDw9M4/LzDN2FYZqnxsy89jMnDNWYXg3KPm11iVnw3yNV3q+Ovkb0pKkNr63kiQVoifa\nSpKkclkcuqDkvmfJ2TBf45Wer04WB0lSG2cOklQIZw6SpI6yOHRByX3PkrNhvsYrPV+dLA6SpDbO\nHCSpEM4cJEkdZXHogpL7niVnw3yNV3q+OlkcJEltnDlIUiGcOUiSOsri0AUl9z1Lzob5Gq/0fHWy\nOEiS2jhzkKRCOHOQJHWUxaELSu57lpwN8zVe6fnqZHGQJLVx5iBJhXDmIEnqKItDF5Tc9yw5G+Zr\nvNLz1cniIElq48xBkgrhzEGS1FEWhy4oue9ZcjbM13il56uTxUGS1MaZgyQVwpmDJKmjLA5dUHLf\ns+RsmK/xSs9XJ4uDJKmNMwdJKoQzB0lSR016cYiIsyLinyLi2Yi4erI/fy8oue9ZcjbM13il56vT\npBaHiOgDbgXOBE4EFkbE8ZN5DL3g8ccf7/YhdEzJ2TBf45Wer06TfeUwD1iXmRsy8w3gHuC8ST6G\nrtu2bVu3D6FjSs6G+Rqv9Hx1muziMAN4oeXxpmpNktRDHEh3wfr167t9CB1TcjbM13il56vTpL6U\nNSI+CizLzLOqx4uBzMwb99jOF7FK0jjU9VLWyS4O7wGeAeYDLwIPAwszc+2kHYQkaa+mTOYny8y3\nIuILwOqqpXWnhUGSek9P/oa0JKm7emogXcovyEXE+oj4RUQ8FhEPV2vTImJ1RDwTET+KiENbtl8S\nEesiYm1EnNHVgx9FRNwZEVsi4omWtX3OExEfjognqvO7vBtZ9jRGtqURsSkifl7dzmr5WGOysfO4\nZkbEQxHxVEQ8GRFXUNb52zPfFynkHEbEARGxpvo58lREfJXJPHeZ2RO3qlD9CpgN7Ac8Dhzf7eMa\nZ5bngWl7rN0I/Pfq/tXAn1f3TwAeq1p8g9XXILqdYY9j/zgwF3hiInmANcC/re7/ADizR7MtBa4a\nZds5TcpWHctRwNzq/sHVzO/4gs7fWPmKOIfAe6s/3wP8DDhlss5dL105lPQLcjHKVdl5wMrq/krg\n/Or+ucA9mflmZq4H1lVfi56RmX8PvLLH8j7liYijgKmZ+Ui13d0t+3TNGNmozuGezmtSNnbm25yZ\nj1f3/wVYC8ws6PyNlm/X7041/hxm5uvV3QOqnymvTNa566XiUNIvyCXwQEQ8EhH/pVqbnplbqL6h\ngYFqfc/cIw3JPbCPeWZU53SXXj+/X4iIxyPijpbL9kZni4jB6irpZ+P4fuz5jC351lRLjT+HEdEX\nEY8Bm4HhzHx6ss5dLxWHkpySmR8GPgV8PiJOrQpGq9JeCVBSntuAP8rMudVfypu7fUATFREHA/cC\nV1b/wi7q+3GUfEWcw8zckZkfqq72To2Iock6d71UHEaAD7Q8nlmtNU5mvlj9+c/A96o20ZaImM7O\nb+SjgK3V5iPArJbdm5J7X/M0Jmdm/nPm7pfxfbOlzdfIbBExpfrB+VeZ+f1quZjzN1q+0s5hZm6v\nZgUfmaxz10vF4RHg2IiYHRH7AxcB93X7oPZVRLy3+lcMEfE+4AzgySrLomqzzwC7/pLeB1wUEftH\nxNHAsdUvB/aa2KOHu095qsvfVyNiXkQEcEnLPt32tmzVX7hd/gPwy+p+E7MB/C/g6cz8estaSeev\nLV8J5zAijtjVDouIg4A/qwbOk3Puuj2N32Myf1b1aoN1wOJuH884MxxdvdLqsaooLK7WDwMerPKt\nBvpb9llSvbJgLXBGtzOMkuk7wP8DfgtsBD4LTNvXPMCfVF+TdcDXu53rD2S7G3iiOo/fq3q8jctW\nHdcpwFst35M/r/6e7fP3Yy9m/AP5Gn8OgZOqPI8BvwD+W47zZ8l4svlLcJKkNr3UVpIk9QiLgySp\njcVBktTG4iBJamNxkCS1sThIktpYHCRJbSwOkqQ2/x/Htrbxmw5PrQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f8081ddb490>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "y,x = np.histogram(nt,np.linspace(0.,3000.,300))\n",
    "x.resize(len(y))\n",
    "plt.step(x,y,where='post',label='N tracks')\n",
    "plt.grid(True)\n",
    "plt.legend(loc='upper right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# testBit() returns a nonzero result, 2**offset, if the bit at 'offset' is one.\n",
    "def testBit(int_type, offset):\n",
    "    mask = 1 << offset\n",
    "    return(int_type & mask)\n",
    "\n",
    "def onlyBit(int_type, offset):\n",
    "    mask = 1 << offset\n",
    "    return(int_type ^ mask)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def hpa(pattern) :\n",
    "    pxp = 0\n",
    "    mul = 1\n",
    "    for i in range(1,5) :\n",
    "        if pattern.hasValidHitInPixelLayer(1,i) : pxp+=mul\n",
    "        mul *=10\n",
    "    for i in range(1,4) :\n",
    "        if pattern.hasValidHitInPixelLayer(2,i) : pxp+=mul\n",
    "        mul *=10\n",
    "    return pxp, 10*pattern.numberOfValidTIBLayersWithMonoAndStereo(0)+pattern.numberOfValidTIDLayersWithMonoAndStereo(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def loadTk(fevents) :\n",
    "    print fevents\n",
    "    events = Events(fevents)\n",
    "    print events.size()\n",
    "    tracks = Handle(\"std::vector<reco::Track>\")\n",
    "    evidsA = []\n",
    "    tkParsA = []\n",
    "    tkHitsA = []\n",
    "    for i in range(0,events.size()) :\n",
    "      if (len(tkParsA) > 1000000) : break\n",
    "      a= events.to(i)\n",
    "      id = events.object().id()\n",
    "      evid = '{:d}:{:d}:{:d}'.format(int(id.run()),int(id.luminosityBlock()), int(id.event()))\n",
    "      if (i%500==0) : \n",
    "        print \"Event\", i , evid, len(tkParsA)\n",
    "      a=events.getByLabel(label, tracks)\n",
    "      for tk in tracks.product() :\n",
    "        # if (not tk.quality(tk.qualityByName(quality))) : continue\n",
    "        pattern = tk.hitPattern()\n",
    "        p,t = hpa(pattern)\n",
    "        evidsA.append(evid)\n",
    "        tkParsA.append([tk.eta(),tk.phi(),tk.pt(),tk.dz(),tk.ndof(),tk.chi2()])\n",
    "        tkHitsA.append([pattern.numberOfValidHits(),pattern.numberOfValidPixelHits(),tk.originalAlgo()-4,tk.algoMaskUL(),\n",
    "                       p,t])\n",
    "\n",
    "    print len(tkParsA)\n",
    "    tkPars = np.array(zip(*tkParsA), dtype=np.float)\n",
    "    print len(tkPars)\n",
    "    tkHits = np.array(zip(*tkHitsA), dtype=np.int)\n",
    "    print len(tkHits)\n",
    "    evids = np.array(evidsA,dtype=np.str)\n",
    "    return (evids,tkPars,tkHits)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def loadFile() :\n",
    "    loaded = np.load('/tmp/singleMuon100.npz')\n",
    "    return loaded['evids'],loaded['tkPars'], loaded['tkHits']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "root://eoscms.cern.ch///store/relval/CMSSW_9_3_0_pre3/RelValSingleMuPt100/GEN-SIM-RECO/92X_upgrade2017_realistic_v10_resub-v1/00000/3C1A3236-AE7C-E711-ACA8-0025905A60FE.root\n",
      "6000\n",
      "Event 0 1:1:2 0\n",
      "Event 500 1:3:205 998\n",
      "Event 1000 1:11:1004 1999\n",
      "Event 1500 1:18:1707 2995\n",
      "Event 2000 1:20:1902 3993\n",
      "Event 2500 1:49:4803 4992\n",
      "Event 3000 1:61:6006 5991\n",
      "Event 3500 1:68:6703 6991\n",
      "Event 4000 1:71:7007 7989\n",
      "Event 4500 1:77:7603 8990\n",
      "Event 5000 1:84:8303 9993\n",
      "Event 5500 1:83:8202 10993\n",
      "11995\n",
      "6\n",
      "6\n"
     ]
    }
   ],
   "source": [
    "fromAOD = True\n",
    "evids, tkPars, tkHits = loadTk(singleMU100) if fromAOD else loadFile()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "if fromAOD :\n",
    "    np.savez_compressed('/tmp/singleMU100',evids=evids, tkPars=tkPars, tkHits=tkHits)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-rw-r--r--. 1 innocent zh 480603 Sep 10 12:04 /tmp/singleMU100.npz\n"
     ]
    }
   ],
   "source": [
    "%%bash\n",
    "ls -l /tmp/singleMU100*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def nhit(var,tkPars, tkHits) :\n",
    "    hp = np.greater(tkPars[2],4.)\n",
    "    notQ = np.equal(tkHits[2],18)|np.equal(tkHits[2],1)|np.equal(tkHits[2],2)|np.equal(tkHits[2],3)\n",
    "    quad = np.equal(tkHits[2],0)|np.equal(tkHits[2],19)|np.equal(tkHits[2],20)\n",
    "    mu = np.not_equal(testBit(tkHits[3],4+9),0)\n",
    "    nbins=62\n",
    "    yn,x = np.histogram(tkPars[var],np.linspace(-3.1,3.1,nbins),weights=1.*tkHits[0]*hp)\n",
    "    ya,x = np.histogram(tkPars[var],np.linspace(-3.1,3.1,nbins),weights=1.*hp)\n",
    "    x.resize(len(ya))\n",
    "    return x,yn/ya,ya\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def printMu(evids, tkPars, tkHits) :\n",
    "    mu = np.equal(tkHits[2],10)\n",
    "#    mu = np.not_equal(testBit(tkHits[3],4+10),0)\n",
    "    for i in range(0,len(mu)) :\n",
    "        if (mu[i]) : \n",
    "            print \"'{:s}-{:s}',\".format(evids[i],evids[i])\n",
    "            #print evids[i],tkPars[0][i],tkPars[1][i],tkPars[2][i],tkPars[3][i],\\\n",
    "            #tkHits[0][i],tkHits[1][i],tkHits[4][i],tkHits[5][i]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "'1:8:735-1:8:735',\n",
      "'1:10:912-1:10:912',\n",
      "'1:2:129-1:2:129',\n",
      "'1:2:147-1:2:147',\n",
      "'1:5:485-1:5:485',\n",
      "'1:11:1007-1:11:1007',\n",
      "'1:12:1146-1:12:1146',\n",
      "'1:22:2163-1:22:2163',\n",
      "'1:21:2012-1:21:2012',\n",
      "'1:21:2090-1:21:2090',\n",
      "'1:29:2859-1:29:2859',\n",
      "'1:29:2862-1:29:2862',\n",
      "'1:29:2879-1:29:2879',\n",
      "'1:71:7015-1:71:7015',\n",
      "'1:72:7192-1:72:7192',\n",
      "'1:73:7210-1:73:7210',\n",
      "'1:74:7340-1:74:7340',\n",
      "'1:77:7674-1:77:7674',\n",
      "'1:78:7797-1:78:7797',\n",
      "'1:86:8536-1:86:8536',\n",
      "'1:86:8543-1:86:8543',\n",
      "'1:82:8174-1:82:8174',\n",
      "18.0 12 2"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/cvmfs/cms.cern.ch/slc7_amd64_gcc630/external/py2-pippkgs/5.0-ghjeda2/lib/python2.7/site-packages/ipykernel_launcher.py:10: RuntimeWarning: invalid value encountered in divide\n",
      "  # Remove the CWD from sys.path while we load stuff.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhQAAAEPCAYAAAAag10gAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X2UXGWV7/Hf7hAM0ISuBi7RBGm4Js7lLcAAsi4amkQc\nwBEM4AtCpKMXHUUIjDgiCi69g8bMXQws7yDOgphADIgQeXNAZEjhiIOjkIYEInh1eBM7UbojxAAD\nyb5/9EmodFd1vT51znP6+1mrFl3ve/eprmyevc855u4CAABoRkfaAQAAgPhRUAAAgKZRUAAAgKZR\nUAAAgKZRUAAAgKZRUAAAgKYFLSjMbJqZ3Wdmj5nZajM7r+S+c81sbXL7wpBxAACAsHYI/PqvS/pb\nd+83s05JD5nZPZKmSHqfpIPc/XUz2yNwHAAAIKCgBYW7D0gaSH7eaGZrJU2V9AlJC9399eS+P4aM\nAwAAhNW2GQoz65F0iKSfS5ohaZaZPWhmK83s8HbFAQAAWi90y0MaLiY6Jd0saUGyUrGDpIK7H2Vm\nR0i6SdJ+7YgFAAC0XvCCIikebpZ0vbvfltz8rKQVGm53/MLMtpjZ7u7+wojncqIRAGiAu1ujz91p\np50GXnnllb1aGxHyYNKkSetefvnlKeXua0fLY7Gkx939ypLbbpU0W8NFwwxJE0cWE1u5+3aXL3/5\ny6Nuy8OFvOK6kFdcl/GWV7NeeeWVvdLOjUs2L2MVmkFXKMzsaElnSFptZqskuaSLJX1H0mIzWy3p\nVUkfrfU1n3rqqZAhp4a84kJecSEvILzQe3k8IGlChbvnhXxvAADQPtEdKbOvry/tEIIgr7iQV1zI\nCwjPvAX9tlDMzLMcHwBkkZmpmaFMvnvba/78+dp777311a9+teWvfeCBB+qqq67SrFmzWvJ6Y322\noluhKBaLaYcQBHnFhbziQl4Yr9asWbOtmLj66qu1//77q6urS7NmzdJvfvMbSdKjjz6q97znPdp9\n9921zz77aOnSpQ29V3QFBQAAW23evDntEKKxYcMG3XrrrRocHNTMmTO1YMECSdKzzz6rz3zmMxoY\nGNCNN96oT37ykxoYGKj79aMrKHp7e9MOIQjyigt5xYW88mXffffVokWLNHPmTHV2dmrz5s3q6OjQ\nb3/7222PmT9/vi699FJJ0v3336+9995bl19+ufbaay9NnTpVS5YsKfvaQ0ND2nvvvfXDH/5QkvTn\nP/9Z06dP17Jly/TLX/5SU6ZMUWk7aMWKFTr00EMlSYVCQZMnT9bkyZPV2dmpjo4OPfPMM3XltnTp\nUr3rXe/a7rbS3ObPn69zzjlHJ554onbddVfNmjVLAwMDOv/881UoFLT//vvrkUce2e53dd9990mS\nLrroIs2YMUMdHR06+uijtW7dOknSe9/7Xp100kmaOHGijjzySE2cOFEvvFD2SA5jiq6gAADgxhtv\n1F133aUNGzZowoQJMht7ZGRgYEAvvfSSnn/+eV1zzTU655xz9Kc//WnU4wqFghYvXqyzzz5bf/jD\nH3T++efrsMMO05lnnqnDDz9ce+yxh+65555tj1+2bJnOOussKSlGXnzxRb344otasGCBjjnmGE2d\nOrXu3EbmMvL697//fX3ta1/TCy+8oIkTJ+qoo47SEUccocHBQZ166qm64IILxnz99evX65JLLtH8\n+fNH3XfBBRfogAMO0AEHHFB33NEVFHntGZJXXMgrLuTVematuTRqwYIFestb3qI3velNUnIQxLHs\nuOOOuuSSSzRhwgSdcMIJ6uzs1BNPPFH2sccdd5w+8IEPaM6cObr77rt19dVXb7tv3rx5uv766yVJ\ng4OD+tGPfqTTTz99u+d/73vf0w033KBbbrlFEyZUOnJC7UbmNnfuXB1yyCHacccdNXfuXO2yyy46\n44wzZGb60Ic+pP7+/oqv9dprr+n444/X+973Pn3605/e7r5Fixbpxz/+se64446G4oyuoAAApM+9\nNZdGTZs2ra7H77777uroeOOfvJ133lkbN26s+Pizzz5ba9asUV9fnwqFwrbbzzzzTN155516+eWX\nddNNN2nWrFnaa683Dh65atUqnXvuubr11lvV3d1dd161KH2/nXbaadT1sfIqFovauHGjLr/88lH3\nXXnllbrmmmu05557NhRXdAVFXnuG5BUX8ooLeeXPyDbAzjvvrE2bNm273shQ4VZbtmzRJz7xCZ11\n1lm66qqrtpvNmDp1qo466ijdcsstWrZsmebNe+MYjevXr9fcuXP1rW99SwcffHBD773LLru0LI9y\nfv/732vKlLKn4tDAwIDe/OY3N/za0RUUAACMdOihh2r58uXasmWL7r77bt1///0Nv9Zll12mjo4O\nLV68WBdeeKHmzZu3Xdth3rx5WrRokdasWaNTTjlFSvY2Oe200zRv3jydeuqpo16zo6NDP/nJT6q+\n98yZM/XYY4/p0Ucf1auvvqqvfOUrVedDRhqr/fPBD35Qt912W9n7nn/+ee277751vVep6AoKeqFx\nIa+4kFdc8ppXNeX+gb3iiit0++23q1Ao6IYbbtDcuXPrfg1Jevjhh3XFFVfo+uuvl5np85//vDo6\nOrRw4cJtjznllFP09NNP65RTTtGkSZMkSc8995weeOABXXHFFZo8ebJ23XVXTZ48Wc8995yeffZZ\nTZ48WQcddFDVWKZPn65LL71Uc+bM0YwZM0bt8VGL0tcbmeeKFStGzXxs9ba3vU1PP/103e+37b2y\nfDS0ckdrKxaLuVzmI6+4kFdcxlteHCkzvOnTp+vb3/62Zs+eXfWx3/3ud/X444/rsssua0tsIY31\n2YquoMiz7m5paGj07YWCNDiYRkQAYkRBEdaKFSt00UUX6cknn0w7lLYb67MV9GyjqM/QUPmp52Z2\nrQIAtM6xxx6rtWvXatmyZWmHkjnMUGQEecWFvOJCXmiVlStXamBgQO9+97vTDiVzoisoAABA9jBD\nkSFmlVse4+jXAGyn0mxROcwbDWOGAqEwQwEgWpVmi8ph3ghIT3Qtj7z2DMkrLuQVF/Kqz6RJk9aZ\nmbhwGXmZNGnSukqfG1YoAADbefnll8sfmxkYAzMUGcIMBTBaPZ9//laGWZMzFEAjomt5jEeFwuhT\n/gY6iR0AAA2JrqAYj73QwcHRp/ytdeo9beNxe8Us9rwqFd+x51VJXvNCnIIWFGY2zczuM7PHzGy1\nmZ034v7PmtkWM+P/t4Fxrrt7dDFgNlwk1Crm4huIXdAZCjObImmKu/ebWaekhySd7O6/MrNpkq6R\n9HZJf+nuo/YeZ4aiNY8FYhDqMz0e/1aYoUAagq5QuPuAu/cnP2+UtFbS1OTuf5T0uZDvDwAA2qNt\nMxRm1iPpEEk/N7OTJD3r7qvrfZ289gzJKy7kFRfyAsJry3EoknbHzZIWSNos6WJJx5U+pNJz+/r6\n1NPTI0nq6uqSJPX29kolf0x5uS4VVSxmJ55WXO/v789UPFzP7vYK9fmXsvP7DbW9isWilixZIknb\nvi+Bdgt+HAoz20HSnZLucvcrzexASfdK2pQUEtMk/U7Ske6+fsRzmaFowWOBGDBD0TrMUCAN7Sgo\nrpP0R3f/2wr3/6ekw9x91Cw2BUVrHgvEgIKidSgokIbQu40eLekMSbPNbJWZPWxmx494mI/V8hgp\nLz3D0bvIFevaPS4WedleI5FXXMgLCC/oDIW7PyBpQpXH7BcyhqwaeQbF4d5xmhEBANA4zuWRkmaX\nYcfjMi7yLdRnuru7/MGtCoXhA2HlES0PpIGCIiUUFMD22v2ZzvPfEAUF0sC5PDKCvOJCXnEhLyC8\n6AoKDCt3EiTOQgoASAstj5SwixywPVoerUPLA2lghQIAADQtuoIirz3DVuVVrhWSZhuE7RUX8opL\nXvNCnNpyLg+0T7nd4IyFTwBAYMxQpKSd/ds894qRH8xQtA4zFEhDdC0PAACQPdEVFHntGZJXXMgr\nLuQFhBddQQEAALKHGYqUMEMBbI8ZitZhhgJpYIUCAAA0LbqCIq89Q/KKC3nFhbyA8KIrKAAAQPYw\nQ5ESZiiA7TFD0TrMUCANrFAE1t1d/qyghULakQEA0DrRFRSx9QyHhob/L2jkZeQhsmPLq1bkFRfy\nikte80KcoisoAABA9jBDEVgW+rRZiAGohhmK1mGGAmlghQIAADQtuoIirz1D8ooLecWFvIDwghYU\nZjbNzO4zs8fMbLWZnZvcvsjM1ppZv5ndYmaTQ8YBAADCCjpDYWZTJE1x934z65T0kKSTJU2TdJ+7\nbzGzhZLc3b9Q5vnMUOQkBqCadn9Ou7uH98IqVSiM3gMrRsxQIA07hHxxdx+QNJD8vNHM1kqa6u73\nljzsQUmnhowDAEYqVzgY/wQDDWvbDIWZ9Ug6RNLPR9z1MUl31fo6ee0ZkldcyCsu5AWEF3SFYquk\n3XGzpAXuvrHk9i9Kes3dl1d6bl9fn3p6eiRJXV1dkqTe3l6p5I8p69el8fv+/f39qf/+uR7H9pKK\nKhbT/32k/ffayPYqFotasmSJJG37vgTaLfhxKMxsB0l3SrrL3a8sub1P0tmSZrv7qxWeywxFTmIA\nqsnC5zQLMbQCMxRIQztWKBZLenxEMXG8pM9JmlWpmAAAAPEIvdvo0ZLOkDTbzFaZ2cNmdoKkb0rq\nlPTj5Laran3NvPYMySsu5BUX8gLCC72XxwOSJpS5a3rI9wUAAO3FuTwCy0JPNgsxANVk4XOahRha\ngRkKpCG6Q28DAIDsia6gyGvPkLziQl5xIS8gvOgKCgAAkD3MUASWhZ5sFmIAqsnC5zQLMbQCMxRI\nAysUAACgadEVFHntGZJXXMgrLuQFhBddQQEAALKHGYrAstCTzUIMQDVZ+JxmIYZWYIYCaWCFAgAA\nNC26giKvPUPyigt5xYW8gPCiKygAAED2MEMRWBZ6st3d0tDQ6NsLBWlwMI2IgNGy8LeShRhagRkK\npCHo2UaRDZWKBuPrBikpV+QWCmlFA6AVomt55LVnSF5xIa/mDA0NrwSUXkKulrG9gPCiKyiyrLt7\n+P/6Sy/8XxcAYDxghqKFYuu/xhYv8iOrn72sxlUvZiiQBlYoAABA06IrKPLaMySvuJBXXGrNq1AY\n3bY0G25nZlFetxfixF4eAJBgjyigccxQtFBs/dfY4kV+xPbZiy9eZijQftG1PAAAQPZEV1DktWdI\nXnEhr7iQFxBe0ILCzKaZ2X1m9piZrTaz85LbC2Z2j5k9YWY/MrPdQsYBAADCCjpDYWZTJE1x934z\n65T0kKSTJc2X9IK7LzKzz0squPtFZZ7PDEVAscWL/IjtsxdfvMxQoP3qWqEws7lJYVATdx9w9/7k\n542S1kqalhQVS5OHLZX0/noDBwAA2VFzQWFm/13STZLObOSNzKxH0iGSHpS0l7uvU1J0SPpvtb5O\nFnqG5Q6x3exhtrOQVwjkFRfyikte80Kc6jkOxXxJ35D0MUlX1/MmyarGzZIWuPtGMxu5eFhxMbGv\nr089PT2SpK6uLklSb2+vVPLH1O7rQ0O9ci9/f7GYfny1XpfCx9vf35+ZfLnO9srb9a3bq1gsasmS\nJZK07fsSaLeaZijMbIKkxyUdLmmFpAvd/ZGa3sBsB0l3SrrL3a9Mblsrqdfd1yVzFivd/X+UeW4m\nZyhi66dWkpc8EJ/YPnvxxcsMBdqv1pbHiZIedPeXJC2W9PE63mOxpMe3FhOJ2yX1JT+fJem2Ol4P\nAABkTK0FxcclXZv8/ANJ7zWzHas9ycyOlnSGpNlmtsrMHjaz45PWyXFm9oSkOZIW1hpwXnuG5BUX\n8ooLeQHhVZ2hMLMuSV3u/hMND1G+YmY3S5ot6e6xnuvuD0iaUOHudzcaNAAAyBbO5dGA2PqpleQl\nD8Qnts9efPEyQ4H2q6nlYWZHm9kuyc9nmtnlZrZP8OgAAEAUap2h+JakTWY2U9JnJf1G0nWBYysr\nrz1D8ooLecWFvIDwai0oXk96DydL+r/u/k+Sdg0cGwAAiEStx6G4PxnAnC9plqT1kh5x94OCBscM\nRVB5yQPxie2zF1+8zFCg/WpdofiQpFclfTw5VPY0Sf8QODYAABCJWguKC9z9cnf/Nw3vDvqMpAPC\nhlZeXnuG5BUX8opLs3kVCqPP3dPd3bLwGpbX7YU41VpQHFfmthNaHAsAZNLg4HDLo/QyNJR2VEC2\njDlDYWafkvRpSfsle3Zstaukn7n7GUGDY4YiqLzkgfjk4bOX5RyYoUAaqhUUu0kqSPq6pItK7nrJ\n3QeDB0dBEVRe8kB88vDZy3IOFBRIQ7WWh7v7U5LOkfRSyUVmlkoHMa89Q/KKC3nFhbyA8Kqdy2O5\npL+W9JAkl1Ra8XrSCgEAAOMc5/JoQJaXOuuRlzwQnzx89rKcAy0PpKHq2Ua3MrOpkvYpfc7WM5AC\nAIDxrdaTg31D0gOSviTpc8nlwvDhjZbXniF5xYW84kJeQHi1rlC8X9Lb3f3VwPEAAIAI1Xouj7sk\nfcDdN7YlqjfelxmKgPKSB+KTh89elnNghgJpGHOFwsy+mezNsUlSv5n9a3JOD2l4huK8tkQJAAAy\nrdoMxS+TXUZvl/S/Jf0sub710nZ57RmSV1zIKy7kBYQ35gqFuy9tXygAACBW1Q69vTJpeZTj7j4n\nWGTMUASXlzwQnzx89rKcAzMUSEO1vTzK7Rp6lKS/k7Q+UEwAACAyY85QuPtDWy+SOiV9Q9Lpkv7G\n3Y9oX5hvyGvPkLziQl5xIS8gvKoHtjKzvzKzf5N0iaTL3P2d7n5XLS9uZtea2Toze7TktiPM7D/M\nbFXy38ObTQLx6u4eXjqu5dKdyunoAAC1qDZD8QtJe0r6B0n/PvJ+d394zBc3e6ekjZKuc/eD9cZc\nxtfd/R4zO0HS37n7sRWezwxFQFnIo54YshBvs7q7paGh2h5bKEiDg6EjSkcetmWWc2CGAmmoNkPx\n56QgOE3SqWXONjp7rCe7+0/NbJ8RN/9e0m7Jz12SftdY6EDlf6Dr+ce43GuE+sd8aKi+AgoAYlFt\nhqLX3Y9NLrNLfj7W3ccsJsZwkaTLzewZSYskfaGeJ+e1Z0hejdn6D/TIi1R7y6Tca1RbRWB7xYW8\ngPCqnsvDzHaWNN3dHym57a2SNrt7I6sL10o6191vNbPTJC2WdFylB/f19amnp0eS1NXVJUnq7e2V\nSv6Y2n1dSvf9W3VdKqpYDPt+/f39Lft9dnYWk/9rfyP+zs7yzx8cHP18s/L5tnt71vp+7dg+jWwv\nrqe3fSptr2KxqCVLlkjStu9LoN2qnsvDzCZK+pWkg939z8lt90i62N1/WfUNhlsed5TMULzo7pNL\n7v+Tu+9W4bnMUASUhTzaGUOl9yp3e6i4xtvMSCV5yC3LOTBDgTRU3cvD3V+T9ANJH9QbqxN71lJM\nJGzE7MWvzeyY5LXmSHqy0eCBehQK5fceKRTCvF+5PVjqea9K8dazB0y5GNhbBkAIVQuKxDWS5ic/\nf1TSd2p5kpktT87/McPMnjGz+ZI+IWmRma2S9PfJ9ZrltWdIXuENDpaft2hk+LKWvMrNZtTzXpXi\nLXepNPNR73xIlrZXK5EXEF7VGQoNr1L8yobNkPRhSe+q8XkfqXDXO+qKEsCYtq5mlLsdANqh6gzF\ntgea9Un6mKTfufvpwSNjhiK4LOSRhRjKaUVcec4tDzE0K8s5MEOBNNTa8pCkmyTNTPbSQA6U69HT\nXwcANKLmgsLdN7n7bu5+b9iQxpbXnmEaeZXr0dd6FMdasb3iQl5xyWteiFNNMxTAeFNpJqHSY/N6\niGwAqFXNMxRpYIai/dqdWx5+l604/Hc7ZeF3noUYmpXlHJihQBpYoQCalMWiAQDarZ6hzEzIa8+Q\nvOJCXnEhLyC86AoKAACQPcxQNCDLvdNmMUORf1n4nWchhmZlOQdmKJAGViiAcYbjjwAIIbqCIq89\nQ/KKS8x5jXX8kZjzGgt5AeFFV1AgTuXOehnyTJ8AgPZihqIBWe6dNitUbnn+neUBszP1y3IOzFAg\nDaxQAACApkVXUOS1Z0hecSGvuJAXEF50BQUAAMgeZigakOXeabOYoRifmKGoX5ZzYIYCaWCFAkAw\n7N0DjB/RFRR57RmSV1zIqzZDQ6OPeeHe/hOqsb2A8KIrKAAAQPYwQ9GALPdOm8UMxfjEdq9flnNj\nhgJpYIUCQNnze3CODwD1iK6gyGvPkLzikre8tp7fY+XKYtlzfMQub9trq7zmhTgFLSjM7FozW2dm\nj464/VwzW2tmq81sYcgYAABAeEFnKMzsnZI2SrrO3Q9ObuuVdLGkE939dTPbw93/WOH5zFC0Gb10\nlGp2u+V5u2c5N2YokIagKxTu/lNJIxdNPyVpobu/njymbDEBAADikcYMxQxJs8zsQTNbaWaH1/Pk\nvPYMY8yr0kGLtr8Uc3kgoxi3Vy3IKy55zQtx2iGl9yy4+1FmdoSkmyTtV+nBfX196unpkSR1dXVJ\nknp7e6WSP6Z2X5fSff+s5Dc0VNTKlWO/Xn9/v84/P1v5cX3s7dXqz3te/16koorFbGyvYrGoJUuW\nSNK270ug3YIfh8LM9pF0R8kMxb9I+oa7359c/3+S3uHuL5R5LjMUbVZPbnn+PWAYMxSVZTk3ZiiQ\nhna0PCy5bHWrpNka/tDPkDSxXDEBAADiEXq30eWSfiZphpk9Y2bzJS2WtJ+ZrZa0XNJH63nNvPYM\ns5JXpQMcNXqCp6zk1WrkFRfyAsILOkPh7h+pcNe8kO+LxrX7pE0AgHzgXB4NyHLvFGglZigqy3Ju\nzFAgDdEdehsAAGRPdAVFXnuG5BUX8hqt3HFJsnL8EbYXEF4ax6EAEImtQ7ojbys3azM0lN0WAIDw\nmKFoQJZ7p0Bo3d3lz0JaqdDIqyx/DzBDgTRQUDQgy18kANojy98DFBRIAzMUGUFecSGvuJAXEF50\nBQUAAMgeWh4NyPJSJ4D2yPL3AC0PpIEVCgAA0LToCoq89gzJKy7kFRfyAsKLrqAAAADZwwxFA7Lc\nOwXQHln+HmCGAmlghQIAADQtuoKi3T3Ddp2fIK+9UPKKC3nFJa95IU6cy6MKzk8AAEB1zFBUjYGC\nAsBoWf5uYIYCaYiu5QEAALInuoIirz1D8ooLecWFvIDwoisoAABA9jBDUTWG7PZJAaQny98NzFAg\nDaxQAACApkVXUOS1Z0hecSGvuJAXEF7QgsLMrjWzdWb2aJn7PmtmW8ysO2QMABBCoTD6oHdmwwfD\nA8ajoDMUZvZOSRslXefuB5fcPk3SNZLeLukv3X2wwvOZoQAQlSx8ZzBDgTQEXaFw959KGipz1z9K\n+lzI9wYAAO3T9hkKMztJ0rPuvrqR5+e1Z0hecSGvuJAXEF5bz+VhZjtJuljScaU3j/Wcvr4+9fT0\nSJK6urokSb29vVLJH1Po61J73y9P1/v7+zMVD9fZXnn6vti6vYrFopYsWSJJ274vgXYLfhwKM9tH\n0h3ufrCZHSjpXkmbkkJimqTfSTrS3deXeS4zFACikoXvDGYokIZ2rFDY1lUId18jacq2O8z+U9Jh\n7l5uzgIAAEQi9G6jyyX9TNIMM3vGzOaPeIhXa3mMlNeeIXnFhbziQl5AeEFXKNz9I1Xu3y/k+wMA\ngPbgXB5VY0i/HwogHln4zmCGAmmI7tDbAAAge6IrKPLaMySvuJBXXMgLCC+6ggIAAGQPMxRVY0i/\nHwogHln4zmCGAmlghQIAADQtuoIirz1D8ooLecWFvIDwoisoAABA9jBDUTWG9PuhAOKRhe8MZiiQ\nBlYoAABA06IrKPLaMySvuJBXXMgLCC+6ggIAAGQPMxRVY0i/HwogHln4zmCGAmlghQIAADQtuoIi\nrz1D8ooLecWFvIDwoisoAABA9jBDUTWG9PuhAOKRhe8MZiiQBlYoAABA06IrKEL1DM3KXwqFIG83\nSl57oeQVF/KKS17zQpx2SDuArEh7iRIAgJgxQwEALcQMBcar6FoeAAAge4IWFGZ2rZmtM7NHS25b\nZGZrzazfzG4xs8n1vGZee4bkFRfyiks78yoUys9jhZDX7YU4hV6h+I6kvxpx2z2SDnD3QyT9WtIX\n6nnB/v7+1kaYEeQVF/KKSzvzGhwcbnmMvISQ1+2FOAUtKNz9p5KGRtx2r7tvSa4+KGlaPa+5YcOG\nlsaYFeQVF/KKC3kB4aU9Q/ExSXelHAMAAGhSagWFmX1R0mvuvrye5z311FPhgkoRecWFvOJCXkB4\nwXcbNbN9JN3h7geX3NYn6WxJs9391TGeyz6jANAAdhtFu7XjwFaWXIavmB0v6XOSZo1VTIg/CAAA\nohF0hcLMlkvqlbS7pHWSvizpYkk7SnohediD7v7pYEEAAIDgMn2kTAAAEIe09/JoiJl91cweSQ6O\nda+Z1bXraVY1e9CvrDKz08xsjZltNrPD0o6nWWZ2vJn9ysyeNLPPpx1PK5Q7CF0emNk0M7vPzB4z\ns9Vmdl7aMbWCmb3JzH5uZquS3L6WdkxAlCsUZtbp7huTn8+VNNPd/1facTXLzN4t6T5332JmC4fH\nSLyuA39lkZm9XdIWSd+WdKG7P5x2TI0ysw5JT0qaI+l5Sb+Q9GF3/1XasTXDzN4paaOk60oHqGNn\nZlMkTXH3fjPrlPSQpJNj314azm1nd99kZhMkPSDps+7+QNpxYfyKcoViazGR2EXSH1MMp2WaPehX\nVrn7E+7+69Lh3IgdKenX7v60u78m6UZJJ6cdVLPKHYQuD9x9wN379cb3xlpJU9OOqxXcfVPy45uS\n7/LcbT/EJcqCQsPV+d+b2TOS+iR9Pe14AuCgX9k0VdKzJdefy8s/UHlnZj2SDpH087RjaQUz6zCz\nVZIGJBXd/fG0Y8L41o7dRhtiZj+WtFfpTZJc0hfd/Q53/5KkLyU97CskzU8x3JpVy0tNHPQrTbXk\nBaQlaXfcLGnBiBXOaCWrmYcms1b3mNkx7n5/2nFh/MpsQeHux9X40OWS/iVwOC1TLa/koF8nSprd\nvqiaV8f2it3vJL215Pq05DZklJntkBQT17v7bWnH02ru/qKZ/VDS4ZIoKJCaKFseZva2kqvvl5SL\nU+6VHPTAtGTuAAACIElEQVTrpGoH/YpY7HMUv5D0NjPbx8x2lPRhSbenHVSLWA62TzmLJT3u7lem\nHUirmNkeZrZb8vNOko7Ly/cg4hXrXh43S5ohabOk30r6lLuvTzuuZpnZr/N40C8ze7+kb0raQ9IG\nSf3ufkLacTUqKfyuTArya919YdoxNavcQejc/Ttpx9UsMzta0k8krU5acC7pYne/O+3YmmFmB0la\nmhSAHcnqy/9JOy6Mb1EWFAAAIFuibHkAAIBsoaAAAABNo6AAAABNo6AAAABNo6AAAABNo6AAAABN\no6AA6mRm0Z8BFgBajeNQAHUys5fcfde04wCALMnsuTyALDCzMySdJ2mipP+Q9KKknczsYUmPufs8\nM/tBck6PSZKudPdr0o4bANqNggKowMz+QtKHJP1Pd99sZv+UHMJ5k7sfVvLQ+e6+wcwmSfqFmd3i\n7kMphg4AbUdBAVQ2R9JhSZFgyQpEuXPGnJ+cr0TJSsX0ZDUDAMYNCgqgMpO01N2/uN2NZheW/HxM\ncqr5d7j7q2a2Mik8AGBcYS8PoLJ/lXSame2p4eKhYGZvlfRfZjYhecxukoaSYuIvJB2VbsgAkA5W\nKIAK3H2tmX1J0j1m1iHpvySdI+mfJa02s4ckfUzS35jZY5KekPTvaccNAGlgt1EAANA0Wh4AAKBp\nFBQAAKBpFBQAAKBpFBQAAKBpFBQAAKBpFBQAAKBpFBQAAKBpFBQAAKBp/x85VtMrJrjPkQAAAABJ\nRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f8081c48350>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    }
   ],
   "source": [
    "printMu(evids,tkPars, tkHits)\n",
    "\n",
    "print tkPars[4][0],tkHits[0][0],tkHits[1][0],\n",
    "eta=0\n",
    "phi=1\n",
    "pt=2\n",
    "\n",
    "var = eta\n",
    "x,y,ya = nhit(eta,tkPars, tkHits)\n",
    "\n",
    "plt.step(x,y,where='post',label='run xyz, lumi??')\n",
    "plt.grid(True)\n",
    "plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)\n",
    "plt.xlim(-3.1, 3.1)\n",
    "plt.xlabel('eta')\n",
    "plt.ylabel('<Nhits>')\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}