{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Effective Area ##\n",
    "\n",
    "Here, we describe the meaning and howto of the neutrino effective area"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ngen = 200000000.0\n"
     ]
    }
   ],
   "source": [
    "ROOT.gStyle.SetOptStat(0)\n",
    "\n",
    "f = EventFile(\"../evtfiles/numu_jgandalf.root\")\n",
    "ngen = float( f.header.get_field(\"genvol\",\"numberOfEvents\") )\n",
    "print \"ngen =\",ngen\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " /pi1/data/shanidze/km3net/TDR/DETECTORS/km3net_wpd_V2.det\n"
     ]
    }
   ],
   "source": [
    "print f.header.get_line(\"detector\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    " The $w^2$ weights in the KM3NeT Monte Carlo event fiels are defined in such a way, that the number\n",
    " of detected neutrino events can be computed like\n",
    "$$\n",
    "N = \\frac{1}{N_{\\rm gen}} \\sum_{i \\in \\rm bin} w^2_i \\Phi(E_i,\\theta_i,\\phi_i),\n",
    "$$\n",
    " where the flux \n",
    " $\\Phi \\equiv \\frac{d N_\\nu}{d E ~ d \\Omega}$ is in units \n",
    " $GeV^{-1} sr^{-1} s^{-1} m^{-2}$. $N_{gen}$ is the total number of neutrino interactions generated.\n",
    " \n",
    " The sum is over events in some bin of energy and/or direction, $N$ is the \n",
    " the number of events in that bin. Usually, we draw all bins using a single command:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAF0CAIAAAC/mplpAAAABmJLR0QAAAAAAAD5Q7t/AAAUF0lE\nQVR4nO3dXZajOBKGYWnOrCsluk9tZXoVnUDvomsrVZmI3hhzEVMaGQSJMYIwfp+LOi6nbWQbPkvi\nJ+wwDAYAtPrX2Q0AgCWEFADVCCkAqhFSAFQjpACoRkgBUI2QAqAaIQVANUIKF/fz509rrbX258+f\nhy20aRrv/WGLuzZCChf37du3YRh+/Pjx7du3Y5bovW/bdvPTN6fqVZORkNpTCMF7H0K494lx9dr2\n9GXv7+/OuY+Pj+mfQgibF9c0zcKf5K8+J31YXHrTNLJleu8XXllMPyX5AOXzHz34x48fxph//vnn\n7e0t+2ryLFn6mkWPFjH90kMIdV0vv86Cban6YDKqNuArzrm6rtc80hjjnOu67q7X77rOGCOLiDd2\nZIx5e3vL/sk555zb8JqyEWbfqbyduGj5TKL03U0fVte1c05uL7+j0aobnxJbJa8ji3t7e5v7BOSN\nyKJjsnz5xrP3pEuUN7LwOgt+/PgxDMOff/451+aFtm1eqGbXDKn1MTH3yPT+0TbTdd3cs6bb7cfH\nx5cLHW3Vj4dU+voyZPj58+foMbJFZePgy09PGjwXUmnwLcRNfNdzm/3n52f64HRZo6fIS8kDpq9m\njJHNfq4Now982ubR2xw9Jfsel/MifcF4e2WqDrefxl7JqNnVQipuP3HVkW8xPiD9RuMj5R55btx6\n05/H9Id6rhM6erX0dSSqRq88bfCoS5WuiHHpcfPIxlm6CHl6/O8op7qui52IbGMWOhTxU8qGVHr/\nQkhJz2XujdR1HUMq26T0WaNgGjVs4R2N1g2RRkC6koxyYW5x8YnTN55+O9L+6YsPK1J17u0QUs/B\n3CaOrHBxJcjeGVfxGFLDZAg2zbvs5mduxxppNg3Jj96anlQ2pGKzs8O0NDjiQud6Utl3kf53Ll/i\nK2c3zvS9DIvDPXn65+enSTpNn7eG+c984bvIpt6at98lhvkO2tztL6VfUJpN0/bPxVB2/YwIqScQ\nJ0rEtAOS5oVzLj4yZoe57QWMQkpef2EDSFe79/d3uVNmrD8+PrK/28M9ITUkv8bZpce2xfciIbXQ\n4PjKaYLHftaQbLrD7XY714NYH1LDMIxCKu1rpN/dqEnD5KclbUk2wVe+/TQdsitJfOJo3VgjXQFG\nt1e+yEIH9sIuGFKpNFyG23VrJIZUfLXsihgX4XIT5OkGnM5G7RhSaXdvuvTRXMmGkBpJP6v4gHR6\ne/pjPvoMs02NY71ps9MnzjUpfeK0L3PXNj99ZLq2TFeS+EQz351cs7i522teJF23X8EFD0GI763r\nurgj3PzaZS7/jlaL+MhlsrNZHt/3fVVVCw/u+15uSGfhEfGljDGy0L7vlw8d2HxgQfqZjFI77npv\n21aa9GUz5rRtm+7In9t3Lo9JUyCNLedc27Zt27rbuZ71vPfTtxD/u7CSpOvVkYcmySEjcQ285FFR\nGQWC70zmdk4q/bk2kw5Rl5uTmr6Uyc1JZbtFZjIlNEzmpKZtnutJuWT2YdQlzL5UdqHre1Kj/5qv\nfqvNujmpdMQUh43pYyTEnXMy6IuZHod7c02KgTVqhlk9JzXczg0Nt6Pp7EoyeuJd3ZnHe1J3fUGX\ncbWQSkd86VeYrnAi/flNcy0+wNzOZE1HHwvDveF2sDCaQR/JhtRo6BqXnsbudB1N35S05K6Qms7L\nLMh+AgsfQvppj1o+6mzKIC7dGueaNHfnXbM2o45Y+nVMV5LRs+5d0IMhddcXdBnXfJ/rV52VjxyN\nOFY+6+PjI3uc1HobpkjXN6/QK9w1wzKSHhu1rUl37W4bLWLDn07x+Ff8XOww83MHbBNCqKrqrPXK\nWtt13atM1ryGC06c41ze+7qu1+yI2F0IwTlHQl0MPSkAqtGTAqAaIQVANUIKgGqEFADVCCkAqhFS\nAFQjpACoRkgBUI2QAqAaIQVANUIKgGqEFADV/n3KUq21pywXwCNOuR7BOSFlTnq3ADY7q2/BcA+A\naoQUANVOG+6Nuo6M/gBkMScFQDWGewBU2zmkvPfe+/Qi/J+fn3LntAY6AHxpz5CSaSbvfdu2cvvz\n8/P333+XMop//fUXOQXgbnsV8EuLMkqd1a7r0jqR7+/vJqnTu9dyARzjrM12t4nzabEz772ElPy3\n7/tRPWsA+NrusSddqrqu0zslnmLN8eNbBeBBZ22YOx+CIFNRaQFbmZaShPrtt9/SGNp30QAuac+Q\nstY650II8R5JqPf397Ztd1wQgNexW5n1pmnatpUpc+G9l9iSKXMhnSlrKe/+Qqz9Pgz/ObsVeNRZ\nm+1uPSnpQFVVFe+RwOr7XoZ7gmwCcJedQ2qESALwIE4wBqDaaefujfYyntUM6GHt97ObAI04wRiA\naoQUANUIKQCqEVIAVCOkAKjGIQgAVOMa5wBUY7gHQLXTelK4qtwxmX9ka99y1jHWIKSws2n0WEse\nYTuGewBUY+8eANXYu4dzMHWFlZiTwjmYusJKxSsYhxDkzuxV8QBgWdkKxiGEeEHh9MrCALDSbiEl\nvacQQtM0cnVzuS31Y0IIzrm0hwUAa+wWUt57KQua3tP3faxsLD2svRYH4EXsWWZd8khqW40Ca8pm\nd+T8wr4/AKJ4BeM5xBCANcpWMHbOxdsyLbXj4gC8gp0nzpumCb+Y23movu+ZOAdwr7IVjCWzZAzo\nnIuT6HgpHKKJRxxR3F2O57xZ6klF5aGZtd+JM83O2myPOC0m24HiBGMAa3CCMQDVuJ4Ung8F2V8K\nIQUtmJBCFiEFQDVCCoBqhBQA1bjGOQDVOAQBgGoM9wCoRkgBUI2QwnYcVIkDEFIAVGPvHgDV2LsH\nQDUqGEO1mWmvTEF2Tv27qv1DSioYx2tISfU9Y0x6J7BSNnooyP5Sdp44b5qm7/v431EFYyqtA7jX\nztViRvfECsbyV7nk+Y5LBHB5e/akhmEYTYdLEeO0eMyOiwPwCva/srq1tuu6mEdp9youa7l8sWHf\n35M4q3QCJRtOcVYhhrIHc3rvnXPSwxqVtBoWFW0VgCdSNqT6vo/BJEO/oosDcD3Fj5OKM+WUWX9q\n649XMhyyhF2VDamu66qqijNQjOOeF8cr4Sz7h1SaRN77YRjYtQdgMyoYA1CNE4wBqMb1pACoRkjh\n+TBb/1IIKQCqEVIAVCOkAKjGNc4BqHZaT4ozilEaFbeugeEeANUIKWzHoQA4ACEFQDVCCoBq7N0D\noBonGANQjeEeANX2Dynv/ai4nvdeyhrvviwAl1e2grH5NffkvW/blsqgAO51XAVjY0wIgYsIA7jL\nniElc+FpVEkqyUBPBn07Lg7AKyhbwTgGlnOu7/vp/XPY94e7zFXcMubv0V0cJb/ZWRWMix+CEId7\n3vuqquKbJIawIypuXVjxQxDSCsallwVswMUSlCsbUnVdt20rt9u2pYIxgHuVDSnZu2etlUkoDkEA\ncK+yFYzNr2Di4AMA2xx07t40oTjBGMAanGAMQDVOMAagGiEFQDVCCoBqhBQA1QgpXBbnxFwD1zgH\noBqHIOCFzF0sIXtVDjpiSpwWUtDM2u+X3ES5WMIzYk4KgGqEFADVCCkAqrF3D4Bq7N0DoBrDPQCq\n7RxS0/LFxpimabjiHYBt9hzuTcsXi3iZcwC4124hNVdK78sSewCwYLeQmpYvNsZI7eK0Zgy0WX+m\nCIdl4xQF9+6FENq2HYZBompkuYfFvr/DcKbI67zTJ1Vw715VVVK+OBaMSf86LCrXKgDPpexxUn3f\nV1Ult9Ma6wCwUsGeVOwW1XVtGMEB2ISDOQGoZk/p4Fh7znKx0lWvJ4VHnLXZcoIxANU4wRiAasxJ\nAVCNkAKgGhPnAFY5a7OlJwVANUIKgGocggBAtdN6UpxRjOcyc00bFMdwD4BqhBQA1QgpAKoRUgBU\nY+8eANU4wRiAagz3AKi2f0iNihiHELz33vtszRgAWLbzcG9UxDiEUFVVXddS3sr8qsQHKJc7dDNT\ni9BQEau8PU9rjnPhXdd5740x8q90rJqmkTJ8hqsg4AlxSeUrXD54WsRYBnpyO4TgnNtxcQBeQdm9\ne3Fw573v+77ruvgnKhgDWKP43r0QgrVWEir2qgwVjE/CWbJ4OmVDKk6cD8OQJhQArFR2uFdVlXMu\nPSiBqAJwl+JHnPd9X1VV/C9DOQB32T+k0hgikgA8iNNiAKjGVRAAqMZVEIBVXvxw8xOdFlIobeaQ\nqMwJaGx+0IyQuqxs9FhLJOHJMHEOQDVCCoBq7N0DoBp79wCoxnAPgGqEFADVCCkAqp1z0WKucQ48\nnbM2W3pSAFTjEAQAqnEIAgDVjqtgnN4JACvtHFLZCsZyO72IMHBV1OPZ3Z7DvWkpvaZpnHPSh/Le\nN01DmXUAd9mzJzUtmdf3fSwP471v23bHxQF4BVr27o0wrQ5AsHcPgGplD+Z0zsXbIYT0vwCwRtmQ\nSueh+r5n1hzAvcqGlOzds9Zaa6XeetHFAbieI84YlOM5b5bKCca4hNxRUX8Y8/f0kReof3HWZstV\nEIA9Wfv9AnmUddZmq+UQBDILQBaHIABQjetJAVCNkAKgGiEFQDVCCoBq7N0DoBrHSQFYhWoxAJBB\nSAFQjZACoBohBUA1QgqAahyCAEA1TjAGoFrx4Z6UL+bCwQC2KRtSMqaTK51z7WAAGxQc7knhYvnX\nGENlUAAblD3O3Vpb17X3vqqqWG/dcFoM8ISueY1z733f93K767o44lsuX2yYVgf0ueC5e03T9H3f\ndd0wDHVdV1WV/nVYVK5VAJ5LwWiUflM6xIudKYZ7wNPVlblgT0rGehJScggCO/gA3Kvg3r2maUII\ncZTXdV25ZQG4qrJHnMejEOhDAdjmiBOMSSgAm3GCMQDVOMEYgGqnhRTwOqz9nrv7j+lBzc91UMIx\nCKkreLojbl5N9tuxlkhahStzAlCNkAKgGnv3AKjG3j0AqjHcA6Aae/eeT25/dmZntmF/Ni6BkHo+\n0+hhZzYujOEecA5+V1YipACoxiEIAFTjEAQAqhUf7jVNQwVjAJuVDalYu5gKxkUxBYsLK1vBuO97\nGdbRmQKwTcGeVLzAucRTrG0FAOsVD6mqqqRmzGi4ZxeVaxWA51J8754M99LaVun9ALCsbHHQ0W1G\nfADuVTCk0qkoKhgD28xcH/2FlB3u1XVNBWMAj7AHzA1NKxhbe8RygQvQU2XjrM2WCsYAVOMEYwCq\ncYIxANW4MiegBYWOswgpQAsKHWdxZU4AqhFSAFRj7x4A1di7B0A1hnsAVDvnOHdOiwGezpVPiwFQ\n2oUvlkBIAVCNkFLtwj+PwEocggBANQ5BAKDaQcM9CsAA2OaInhQJtVJuBipzBrx5sZPgMbL+Ygnm\nEqtK8ZCicPF60/WJM+Ax9WoXSyg73AshtG1b13XRpQC4sLIhVVVVXdfZa5xTwRjAGsWLg3rvpfTe\nqDLosKhcqwA8l4In40w7RF3XSXJx7h6wrwMqX5212R5Ud6+qqnRBhBTwdDjBGAAyuFQLgFXoSQFA\nBicYA1CNE4wBqMZwD4BqhBQA1QgpAKoRUgBUY+8eANXYuwdANYZ7AMZUlSkipIAXoip9ViKkAKhG\nSAFQ7bSJcwBKrC5T9PcBjZniEATgsjYXScvWnrHWGnNCQRoOQQAu6xpF0orPSXnvvfdU3wOwTdmQ\nkjGd975tWwpVAdigYEhJ7ymE0DRN13VmUtUKgE6qhoQF56SmNUGzVUJf0wEFiICpZ1zryoaUpFLT\nNNNi68ujP6bVAYjie/ckjOq6Hs2dE0MA1igbUtZa59yLT0XNnC216nAVAAVDSrpOTdPEkHrNOals\n7jzj4SrAKQqGlGRTVVXxnq7rXjOnAGxGBWMAq1DBGAAyOMEYgGqcYAxANYZ7AFQjpACoRkgBUI2Q\nAqAae/cAqMbePQCqMdzb2TMWXwQ0o6TVdusvb8C5xMBmhFTGynOUzrq8gfIzHzU3T3PbjPrmnYXh\nHgDVuApCBs17hObmaW6boXkzOAQBgGocggBANaUVjJdryXxZZ/TBpz/y4o8/4NzmPf7uHny65uaV\n/uq/pHm7KOe4CsZUWgewQfFrnMcqDM+eU5TzBE5RPKSEdKbKLWt3uQM1/3+UZvwrsQWUVnCfohSG\niVFlrY3VYtSOfgEsuNohCGnvSaIq1rNi1x6AlTjiHIBqZQ8hjUO80dAPAFYq25Oq67qqKmtt3/ej\neXSdgbXtqK5jhBA0Ny9SOOHovbe/6Kyh3TSNzm+2aRp76/gtt+wR503TNE0jW1d6Z9/3RZe7jbXW\nOSdTaW3baps4q6oqNi+EoDPlFSaUMabve/nozm5Inve+7/u6rhV+s+mHFkI4Z8sdjhWX23XdwYte\nVtd1/DS6rtPWwrR56W1VpGEK26awSZGsbPG2c+7U5ixxztV1ffxyj544l6UevNA1vPdxG4v3nNSW\njKZpZG2WX1rn3NktGgshtG07+gw1iHuWdY714jHPMtZT1Y1KSTfqnAHp8bkoIaWqn5KSzeyUX4wv\nxQhQ+OnJh5b2C5SIX6jOjl78vZEbantSJ24UhNQNWV10JlQka/PZrbghTeq6ToJA55c7qBzIp9+m\nwogX584wEFL/p/l3TDoCclvhqjztoev5fruuSxujqm3DZPvX1jxx1myU4GDO/5HBtuyL1LaHRYwO\n31clrk8xQPXM/oQQqqqSD02+ZT1tM7+apLZ5ou/7M1t1SjQafT8X06lozS3U1rZIYS9vUP/RpXsb\nFDbv9O9U9TWVMTI6BRJ3GR2vp43y5p2IkAKgGnNSAFQjpACoRkgBUI2QAqAaIQVANUIKgGqEFADV\nCCkAqhFSAFQjpACoRkgBUI2QAqAaIQVANUIKgGqEFADVCCkAqhFSAFQjpACoRkgBUI2QAqAaIQVA\nNUIKgGqEFADV/gt0ZenjXC9eLwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "c1 = ROOT.TCanvas(\"c1\",\"c1\",400,400)\n",
    "\n",
    "hevts  = TH1D(\"hevts\",\"events for flux of 1e-4 (E/GeV)^{-2} GeV^{-1} m^{-2} s^{-1}\", 24, 1,7 )    \n",
    "\n",
    "f.roottree().Draw( \"log10(mc_trks[0].E)>>hevts\",\"w[1]*1e-4/\"+str(ngen)+\"*mc_trks[0].E**-2\" ,\"goff\")\n",
    "\n",
    "hevts.Draw()\n",
    "c1.Draw()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    " The plot above simply shows the number of detected events in each log(E) bin.\n",
    " \n",
    " \n",
    " The effective Area is defined as the ratio between\n",
    " the rate of detected events and the total flux of neutrinos. \n",
    " For the (average) effective area in a bin this gives:\n",
    "$$ \n",
    "A^{\\rm eff} \\equiv \n",
    "\\frac{N}{m \\int_{\\rm bin} \\Phi dE ~ d \\Omega} =\n",
    "\\frac{1}{N_{\\rm gen}}  \\frac{\\sum_{i \\in \\rm{bin}} w^2_i \\Phi(E_i,\\theta_i,\\phi_i)}{ m \\int_{bin} \\Phi dE ~ d \\Omega},\n",
    "$$\n",
    " where $m = 365\\times24\\times3600$ accounts for the fact that the \n",
    " numerator is the number of events per year while the integral in the numerator is \n",
    " the number of neutrinos per second.\n",
    " This relation is true for any $\\Phi$. It is often convenient to use \n",
    " logarithmic bins for the energy and bins of equal solid angle for the direction.\n",
    " In that case an isotropic flux of $\\Phi = E^{-1}$ \n",
    " is a convenient choice as the integral is the same in each bin, namely:\n",
    "$$\n",
    "\\int_{E_{\\rm bin}} \\int_{\\Omega_{\\rm bin}} E^{-1} dE ~ d \\Omega = \\Omega_{\\rm bin} \\log(10) W,\n",
    "$$\n",
    " where $\\log$ is the natural logarithm and $W$ is the width of the 10-logarithmic energy bin \n",
    " [$\\log_{10}(E_{\\rm bin}^{\\rm max})- \\log_{10}(E_{\\rm bin}^{\\rm min}) $].\n",
    " \n",
    " Hence, the effictive area in $m^2$ can be computed as:\n",
    "$$\n",
    "A^{\\rm eff} = \\frac{1}{ m N_{\\rm gen} \\log(10) W \\Omega_{\\rm bin}}\n",
    " \\sum_{i \\in \\rm bin} w^2_i E_i^{-1} .\n",
    "$$\n",
    "\n",
    " It is intesting to note that the effecitive area is proportional to\n",
    " the number of events detected from a physical $E^{-1}$ neutrino \n",
    " spectrum.\n",
    "\n",
    " Hence, the effective area can be plotted with a few lines in aanet:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAF0CAIAAAC/mplpAAAABmJLR0QAAAAAAAD5Q7t/AAARGElE\nQVR4nO3dW5qqyLYGUDjf6RdYjTndELtRjdkELeM8RBWbhabphcsExnhYn+ZSIzJN/pwxCbXs+74A\niOp/tp4AwDNCCghNSAGhCSkgNCEFhCakjq8sy7Is8+WmacqyrOt6loed5XHgOSF1cDlHrtdrURQp\npdvtVlVV0zSfPZpgYn3/u/UEWFbXdUVRjFOpaZpZgsYOO9ahkjqIlFJd17nSSSnlLw5hlL94uVyK\nomiaZsish/fKhoXhcOP8aF3X5QvD44xvM9xy+ErOxF9LsPHNxjPJV/Nknszt/tv5uFoknJ79a9s2\nP5tVVeULbdv2fT9crapquE1RFNfrte/7vAa8v1f/b4lUVVX+r6qqJo+WbzNcGP8i5Ye9n8CT37f7\nm+UZjuc83HdyyzyHh99OfhD2TkgdwfhY7f9MjfHlfMyPk2hyr3x1nDL9vwlyf5efbj+MmIcbJ8Vk\nxCfzz1eHkBr+azyZydCTW/obfBiWe7uXF0fjxdTQJn9yr7waGq+JqqrKDazb7TZ+wJTSUKQ8eajx\ncOMJPBxiou/74e730x4qo6Iouq4bXx0PnX+hf3oQdmzTiGQGPyXIfX0xrqTGh/rk96H4od7pf6ik\nxqPkyeQvPhliom3byY3HldSkFrs3TON+xA9/pkTi7N5BXK/XSWf6eaO6ruuu6yZNn29Gv91uwxaH\n8X+9MkTu6OdvIXe+n9z4py0U+V4vPgh7snVK8q37cmNczoyf5XElNWkk9X2f2+T9Xd9nfK/ih0pq\nuDqeyZMhxu5vVvxWSd3P7ddbsl+exSMY0qFt23zM3y/E+keN8+EIz/can1PLJwQnresnIfXw/N14\niEmEDYbhxo/zU0gN391kbg8fREgdg2fxIH4qkJ+E1GQhNk6cSZ9rkmvFo9bVJBx/HWJsvEIc73u4\nD6knc7vvasmpYyh7+4aP4n4P5+v3eniXh//1n//8pyiKv/76a5YhJrcZbvbrXR7e4N0HYReEFBCa\nfVJAaEIKCE1IAaEJKSA0IQWEJqSA0IQUEJqQAkITUkBoQgoITUgBoQkpIDQhBYQmpIDQhBQQmpAC\nQhNSQGjbfKSVjxuCPdrkjXw3+9w9b1sM+7JVbbHsci+lVJZlWZYPP80R4FfLfhBD/izZpmnK8o+B\nJleB+LY6bBevpHJOTT4TDeBFi5/dyyHVdd3SAwGH9G1I1XU9fCJjMSqd8hfLskwpaUgBH/sqpJqm\nGZdIKaXL5ZIv5wtt214ul7IsLfeAz3zeCRvOR7Ztmz/JevLZ1rll/tN9Nc5hX/bXOO/7fjLjruty\nThVFUdf17XZ7cvfyqY9nBRyMzZxAaJuF1KRcklnAQ3OG1Lg7nlJ63iyXSsAr5twnNe5DdV33fOeB\nJhTEVJZ/bz2FP8xZSTVNk1+sVxRFVVVDE/0hlRQEVJZ/9/3/bT2LP8x/TjHv5/xlVFsQIJ7nCbXV\nYbvRqHfrO5kF2/q1htoqpGxBACKu8gabhRQQQW6Th02owj4pOLPIBdTAcg9OahcJVfi0GDinvSRU\nYbkHJ7SjhCos9+Bs9pVQhbN7cAyvv5ZlXwlVCCnYu/h7CL6kJwV7dfh4yvSkYJd211r6mOUe7MxJ\nCqiBkII9OU8BNRBSsA9nK6AGQgqiO208Zc7uQVwnj6fM2T0IZ9iZefJ4yiz3IITxlnHZNCakYEuK\npl8JKdiAbHqdkIJV6YW/S0jBSsTTZ4QULMvK7kv2ScFSlE6zsE8K5ieeZmS5B3MST7MTUjAP8bQQ\nIQXfEk+LElLwOfG0AiEFnxBPq1n8E4zLsizLsmmapQeCdZTl3/ntMSXUOspFtwLUdV3XddM0ZVm2\nbVvX9T+jlsuOC0s4efW01WG77KgppbquU0qXy2U8kJBid0745uITWx223y73cgYNV3MqDV/MZdTl\ncmnb9suBYCvD+m7riZzUV43zpmm6rhuu5oqpqqqiKHLp1DRNSknRxE6dfH0XxOf12/Diu6HZlP8d\naqi6rm+323B7PSl2RDzd2+qw/bySytMdv06467rr9Zov54R68i1NXmD88MFhfeIpGi8whn+Ip5i8\nVQuIp9DmDKncMs9SSuOr96QSQThzF9ycO87HnfKu6+wyJz4JFd+clVTecJDXcVVVDefyHrLcY1uW\neHsx/znFvJ/zl1FtQWBTCqgPHPNlMT+Oerf/QGaxGgn1mf3tk/qSVGJ9lnh75P2kOAsF1E7ZJ8Up\nSKj9stzj4Czx9s5yjyNTQB3A4m8fDFuRUMegJ8UBWeIdiZ4UhyKejkdPiuOwvjskIcURKKAOTE+K\n3VNAHZueFDumgDoDyz12STydh5Bif6zvTkVIsScKqBMSUuyGAuqcnN1jBxRQZ+bsHqGJJyz3iMv6\njkJIEZMCioGQIhwFFGNCikAUUNwTUkShgOIh78xJCBKKn9gnxfYkFE/YJ8WWNKH4lZ4Um1FA8Qo9\nKbYhoXiRkGIDEorXCSnWJqF4i54U69Em5wNCipUooPiM5R5rkFB8bOaQqus6pTRcTSnVdT35Imcj\nofjGnCHVNE3XdcPVlNLlcsmXhwucjYTiS+VcO7+Hl7m0bVvXdVEU+d9cQ+V6qmma4cZ2nJ+BhDqS\nrQ7b2Sqpvu8n30DXdTmniqKo6/p2u801FrsgoZhFlBcYT6iz9k5CMRcvMGZ+EooZLRhSVVUNl1NK\n46uFt2o5LgnFvBbcJzXuQ3VdN3TNs/5Py02DNUkoZrdgSDVNU1VVWZZlWVZVNTTRAV63+DnFvJ9z\nOupd11wxdQDKqGPbagvCRqPaJ3U4Eurwtjpso2xBkFm7JqFYji0IfEtCsSjvgsBXJBRLs9zjcxKK\nFVju8SEJxTos9/iEhGI1QgoITU+KtymjWJOeFO+RUKzMco83SCjWJ6R4lYRiE3pSvERCsRU9KSA0\nyz1+p4xiQ0KKX0gotiWkeEZCsTkhxY8kFBE4uweE5uwejymjCMJyjwckFHEIKaYkFKEIKf4goYhG\nSPFfEoqAhBQQmpDiH8ooYrJPiqKQUARmnxQSitAs985OQhGckDo1CUV8Quq8JBS7IKROSkKxF0IK\nCE1InZEyih0RUqcjodiXlUKqruuU0jpj8YSEYnfWCKmmabquW2EgnpNQ7NHiO84nL38BeMvilVTf\n914BE4Eyip2K8gLjCbk2LwnFfnmBMRBalEpKZi1HGcWuqaQOTkKxdyopILSVQuo+g6TSCpRRHICX\nxRyWhOIYLPeA0DTOj0kZxWFY7h2QhOJILPeA0Cz3jkYZxcFY7h2KhOJ4hBQQmp7UcSijOCQ9qYOQ\nUByV5R4QmpA6AmUUB6YntXsSimPTkwJCs9zbN2UUhyekgNCE1I4pozgDjfO9klCchMY5EJrl3i4p\nozgPIbU/EopTEVJAaEJqZ5RRnI2QAkITUnuijOKE7JPaDQnFOdknBYRmubcPyihOS0jtgITizIQU\nEJqQik4ZxckJKSA0IRWaMgrWC6nJxih+JaGgWCekUkoSCvhMudqmyrL871jjyzykjCKarQ7bDyup\nuq5TSsPVlFJd18MX08gMcwRO7JOXxTRN03XdcDWldLlcqqoqiuJyufR9X9f1XPM7J2UUDN4Oqfvu\nUtM0VVXloqmu66ZpmqaZY24A7y/3+r6frEu7rhtKp7qub7fbT3ccXy2fendWR6KMgjHvggCE5v2k\nYlFGwcQMIZVb5llKaXz1CakEvGKGzZzjPlTXdS92zTWh7imj4N4MlVTTNMOe8qqqXtx/oJICXjHb\nFtK8n/PVUe04v6OMIritDtvZGufvbuDUOAdeYQtCCMoo+In3kwJCs09qe8ooeMJyDwjNcm9jyih4\nznIPCM1yb0vKKPiV5R4QmpDajDIKXiGkgNA0zrehjIIXaZwDoVnuAaEJqQ1Y68HrhBQQmsb52pRR\n8BaNcyA0y71VKaPgXUIKCE1IrUcZBR8QUkBoQgoITUitxFoPPmOfFBCafVJrUEbBxyz3gNCE1OKU\nUfANIQWEJqSWpYyCLwkpIDQhtSBlFHxv8ZAqy7Isy6Zplh4IOKRlQ6qu6+v12vf97XZLKS06FnBI\n5aKbKlNKdV2nlC6Xy3igslx23Ais9TiYrQ7btyupHDrD1RxDwxfTSL5x0zSXy6Vt27lmDJzKey+L\naZqm67rhai6RqqoqiiLXSnVdT26fUjp80XRPGQVzeaN+G14S3LZtDqP871A05brp4V3G9ypOsNwT\nUhzPDpZ7fd9Ppth13ZA7dV3fbreHd8kmRVb51AffSRwSCmbkXRCA0Lyf1MyUUTCvr0Iqt8yzlNL4\n6q8Ok0rAor7azDnuQ3Vd99a28iM1oYDlfFVJ5R0GOWKqqpq0xp87ZCVlrQezezukJuGS922+FU8A\nr5uhcf5ZQh21cQ7MyxYEIDTvJzUbDSlYgn1SQGiWe0BolntAaJZ789CQgoVY7gGhWe4BoQmpGVjr\nwXKEFBCaxjkQmsY5EJrl3rc0pGBRQgoITUgBoWmcA6FpnH9FQwqWZrkHhCakgNCE1Oes9WAFQgoI\nTUgBoQkpIDT7pD6kIQXrsE8KCM1yDwhNSAGhCalPaEjBaoQUENqyIZVSKsuyLMumaRYdCDiqctGz\nbHVd13XdNE1Z/jHQ5OruWO5xQlsdtotXUjmnqqpadKA1SShY08z7pFJKw+W6rot/i6nb7TbvQMBJ\nvF1J1XU9TqKhVspfrEeKoijLMqWkIQV87L2Qapqm67rhakrpcrnky8OFsbZtL5dLWZaHWe5Z68HK\n3uiEDa+2a9t2WMoV/y7xhh75iw+108a5kOK0dtA47/t+MsWu63JOFUXxbuOpfOr1xwGOzQuMgdDs\nOH+DtR6s76tKatwOTym91R3f+/tJAev4qpIa96G6rntrq0H/p2+mARzYV5VU0zT51XlFUVRVNTTR\nX6GSAl4xwznFvJ/zvVF3uAVBQ4qT2+qw3WjUu00Gu8ssOJsd7JOa1756UmX599ZTgJOyBeF3Fnqw\nIR9p9QsJBduy3HtgCNCYCRX8ZUORpxd5bkX46W3Fcu9HMRMKzkZIPSahIAg9qal8Ik9CQRB6Un9Q\nQEE0QZd7zzuIv/YXP7v7iwn15ei/3uDL7unSo595eks/9b/a5LjYXNCQWp8aCmLSkyoKCQWBeWdO\nCQWhRXmBMRDfid4FAeBFGudAaEIKCG2bkJp8Vnscb33E6cqGT7SPOb1BwIZjXdfDRzq++y6y62ia\nJuYz2zTN5DMx1z9yNzi7N/ms9jjyx8HnT5e43W7RunWXy2WYXkopZsoHTKiiKLque/c9+NdU13XX\nddfrNeAzO/6hpZS2OXL7dQ3jtm278tDPXa/X4afRtm20GY6nN74cSp5YwLkFnNIg/7INl6uq2nQ6\nz1RVdb1e1x937eVeHnXlQV9R1/VwjA1f2WguDzRNk3+b81/atz7icB0ppdvtNvkZRpALk7ziC/Wc\nZnl6KaW81gtVRo3lMmqbBen6uZhDKlSdMpYPs03+YvxqiICAP738QxvXBUEMT2jMQm/4e5MvhK2k\nNjwohNQf8q9LzIQa5N/mrWfxhzyltm1zEMR8cvuQC/nxsxkw4rNtOwxC6r8i/x3LhUC+HPBX+b5C\nj/P8tm07nkyoufV3x3+06WVbdaMy+6T+kRfb+TOZo51hyYZPtA84t+H3aQjQON2flNLlcsk/tPws\nx5lb8e+Uwk4v67puy1ltEo1FvD8X963oyDOMNrdBwCqvD/+jG59tCDi9zZ9Tr93bk+FE1cbz2Ke8\nG3brWfwo+PQ2JKSA0PSkgNCEFBCakAJCE1JAaEIKCE1IAaEJKSA0IQWEJqSA0IQUEJqQAkITUkBo\nQgoITUgBoQkpIDQhBYQmpIDQhBQQmpACQhNSQGhCCghNSAGhCSkgtP8HOAeOlob3K2QAAAAASUVO\nRK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "c2 = ROOT.TCanvas(\"c2\",\"c2\",400,400)\n",
    "\n",
    "haeff  = TH1D(\"haeff\",\"effective area\", 24, 1,7 )\n",
    "f.roottree().Draw( \"log10(mc_trks[0].E)>>+haeff\",\"w[1]/mc_trks[0].E\" ,\"goff\")\n",
    "haeff.Scale ( 1 / ( ngen * log(10) * haeff.GetBinWidth(1) * 3600.0 * 24 * 365 * 4*pi ))\n",
    "\n",
    "haeff.Draw(\"Lhist\")\n",
    "c2.SetLogy()\n",
    "c2.Draw()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note the $\\Omega_{bin} = 4\\pi$ is the solid angle of the region of the sky we *average*\n",
    " over: it is not by definition equal the generation phase space of \n",
    " the production. We may for example compute the average effective area for\n",
    " up- and downgoing neutrinos with (note the 2pi, regardless of whether the files contains a production of upgoing or \n",
    " up- and downgoing neutrinos.) :\n",
    " \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning in <TFile::Append>: Replacing existing TH1: haeff_up (Potential memory leak).\n",
      "Warning in <TFile::Append>: Replacing existing TH1: haeff_dn (Potential memory leak).\n",
      "Warning in <TCanvas::Constructor>: Deleting canvas with same name: c3\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAF0CAIAAAC/mplpAAAABmJLR0QAAAAAAAD5Q7t/AAAUFklE\nQVR4nO3dbXajOBOGYWnOuy/ELGZ6GRbLyCxmkFbG+6MSRQYbfwEqifs6ffrYadsoif10lRBgp2ky\nAKDVX6UHAABrCCkAqhFSAFQjpACoRkgBUI2Qap+11lort7331lrn3CYvu8nrAOsIqcZJjlwuF2NM\nCGEYhq7rvPfvvRrBhOP9r/QAsK8YozEmTyXv/SZBwwo7HINKqhEhBOecVDohBPliCiP5Yt/3xhjv\nfcqsm88SqTFMD5ZXizHKjfQ6+WPSI9NXJBMflmD5w/KRyF0ZzMrYlt/O29Ui1JlQv3Ec5bfZdZ3c\nGMdxmqZ0t+u69BhjzOVymaZJesDls6afEqnrOvmnrutmryaPSTfyN5K87HIAK++35cNkhPmY03Nn\nj5Qx3Px25EVQO0KqBflndbpOjfy2fObzJJo9S+7mKTP9JMjyKfcen7Yom8uTYrbFlfHL3RRS6Z/y\nwcw2PXsk/wc3g3avetIc5c1UmiZfeZZ0Q3lP1HWdTGANw5C/YAghFSkrL5VvLh/AzU3MTNOUnr4c\ndqqMjDExxvxuvml5Q997EVSsaERiA/cSZFlf5JVU/lGfvR/MnXpnulNJ5VuRwcgXVzYxM47j7MF5\nJTWrxZbSMJZbfPNnCk3Yu9eIy+Uym5len6h2zsUYZ5M+n2x9GIa0xCH/p2c2ITP68i3IzPfKg+8t\noZBnPfkiqEnplMSnluVGXs7kv+W8kppNJE3TJNPk02LeJ3+WuVNJpbv5SFY2kVs+zDyqpJZje/hI\n1IvfYgtSOozjKJ/5ZSM23Zo4T59weVa+T012CM6mrldC6ub+u3wTswhL0uby17kXUum7m43t5osQ\nUm3gt9iIewXySkjNGrE8cWbzXLNcM7emrmbh+HATubxDzNc9LENqZWzLWS1yqg12Yt1wK5ZrOJ9/\n1s2n3Pyn//77zxjz999/b7KJ2WPSwx4+5eYDXn0RVIGQAqAa66QAqEZIAVCNkAKgGiEFQDVCCoBq\nhBQA1QgpAKoRUgBUI6QAqEZIAVCNkAKgGiEFQDVCCoBqhBQA1QgpAKoRUgBUI6QAqFbmklZcbgio\nUZET+Ra77h6nLQbqUqq22LfdCyFYa621N6/mCAAP7XshBrmWrPfe2qsNze4C0K/Ux3b3SkpyanZN\nNACaWaNo1nj3vXsSUjHGvTcEYBPW2MkoanQ+DSnnXLoio8lKJ/mitTaEwIQUUAttCWU+DCnvfV4i\nhRD6vpfbcmMcx77vrbW0e4B+ChPKfDJxnvZHjuMoV7KeXdtapszvPZeJc0CVhwlV38T5NE2zEccY\nJaeMMc65YRhWnm5XvT0qAG/QWUMJFnMCZ6c5oUzBkJqVS2QWUITyhDLbhlQ+Ox5CWJ8sJ5UAPGPL\ndVL5PFSMcX3lAZNQQHE3yyhr/y0ymHu2rKS893KwnjGm67o0iX4TlRRQ1r2EmqZ/ioznnu33Kcp6\nzgdbZQkCUNQbCVXqY1toq4v+jswCDvNeDVUqpFiCAJxLLV1eUiykAGzo+fMWzBJKpsnVJpRhnRRQ\nO4mn99Y6aS6gEto9oGKfLMWsIqEM7R5QqU8KKFNPQhnaPaBGHx7LUlFCGdo9oDqnSihDuwdUZKXF\ne/5YlroSyhBSQC3uFVD61xB8iDkpQLt7BVTz8SSYkwL0Wu/vmo8nQbsHqJOWj5fp76w1mmoIQgpQ\nIT+uZWXn3b4FlEzCaEooQ0gBZa0XTVeP3LWAUhlPgpACynh+yfhp40mwdw84GvH0EvbuAcd5Mp7S\nysyTx5Og3QOO8DCe8iXjr2XTG9cxqSSeBCEF7G7laLs3i6Y8mKpKnDcQUsCO1heLm5ey6UzBlCOk\ngF1sdixLyqYzBVOOkAI2tnE8nTWbEkIK2NJy+un9zu708SRYJwVsZpZQ76xyIp4WWCcFbGDW4hFP\nG6LdAz6VF1DE0+YIKeAjKaGIp50QUsCbUov3WjyddbnT22yRuSFry2wX2IoUUA/i6eYBK9W+80t9\nbAkp4GXWWGO/TB5PbeXRTaU+tn/tvQFrrbXWe7/3hoADWGOtsZM1k/kzmT/G2u8/03TjD7awb0g5\n5y6XyzRNwzCEEHbdFrAv+xNP1pBHR9q3fgshOOdCCH3f5xui3UNdnj9NXcNqbfckg9JdSaX0Reec\n977v+3EcP9wQUIS0c8Z+nTyhCvpoCYL3PsaY7krF1HWdMUZKJ+99CIGiCdX5vT7Cd3NXdjin9n79\nlg6+G8fROWeMkb9TDeWcG4YhPT49zNDuQavf60rZr8n8YbIpV+pj+34lJcPNjxOOMV4uF7ktCbXy\nLdnVc54SYTjYfGUmCaUGBxjj7PJ4subfafpH2yV8T45TteC8bsSTUXeRcWwZUjJlLkII+d0lUgll\n5QcG/8aTaW2ZeAO2XMyZz5THGFllDp2+V40vE4plmSptWUnJggPp47quS/vybqLdQxG3z6xCi6fY\n9vsUZT3ng62yBAEl3CigDAn1rNOdBWH2FTILu8qPayGh3lPfOqkPkUo4zN0WzzBNXgHOzImWUUA1\ngHVSaNbs+ggkVKVo99AmWrxm0O6hQezFawkhhaYwCdUe5qTQDlq8JjEnhUbcvsYUBVT9aPfQgpRQ\n82tMkVD147p7qN78KnjE0z44LIa3FF72fbZf+/V9pjrBe2kfpwspUgkf+rmIy5/v+7yjdkZIAS+Q\n63QaQzYd53QHGANv+04o4ukcCCnU5HetJgF1GoQU6sCFzk+LFefQ7jue6O/OihXnUI0JctDuQSkK\nKAhCCur8xpMhobDpdfeAD8kV8Yz9+i6gSChQSUGP30PwzB/iCQkhBRWujnEhoZAhpFDeb0IRT1jg\nLAgojEUGtTjdsXukEgxH4eEJtHso5vuU5AQUVrEEASVYm1+5E1hBSOFw1lJB4XmEFA5krbHWmq/S\n40BNmJPCIaw1xnzH0/SHMgrPI6SwP2vN9H29Kaai8CraPeyMhMJnNg4p51wIId0NITjnZl/EiZBQ\n+NiWIeW9jzGmuyGEvu/ldrqBEyGhsIXN1rmnI13GcXTOGWPkb6mhpJ7y3qcHs+K8cSRUc0p9bDer\npKZpmn0DMUbJKWOMc24Yhq22Be1IKGxHy4UYZqizKkZCYVMcYIxNkVDY2o4h1XVduh1CyO8aLmnV\nLkmo0qNAO3ZcJ5XPQ8UY06y5mK7tNwwcx1prviShKKOwlR1DynvfdZ211lrbdV2aREebrE2nhSKh\nsKHd9ynKes75VjkzZ2OYijqBUksQip0+mFRqBwl1Dqc7fTAT540gobAzliDgAyQU9sdZEPAuEgqH\noN3DW7KEKj0UNI52D6+7TijKKOyKdg/vIKFwGEIKL2LRJo7FnBRewWQ5DsecFJ5GQqEE2j08h915\nKISQwhPYnYdymJPCIyQUimJOCquyfXmGhEIJtHt4jKkoFERI4T4aPShASOEOEgo6EFK4hYSCGuzd\nw8L1ZDlQFnv3cO0noSijoATtHjIkFPQhpPCDhIJKhBSMMTcSClCCkMLthKKMghKE1Ond2pdHQkEP\nQgrfmIqCTqyTOjcmy6Ee66ROjMly1IB276yYLEclCKlTIqFQD0LqfEgoVIWQOhnmoVAbQupM7pze\ngDIKmhFSZ0Sjh4oQUqfBVBTqdFBIOedCCMdsCyuYikJ1jggp732M8YAN4a5b1x+mjEIVdl9xPjv8\nBQVwCDFqtnslNU0TR8BowFQUKqXlAOMZcm0ztxo9oCIcYHw6lFGoi5ZKiszaxaKMIqFQHSqpdpFQ\naAKV1FmQUKjUQSG1zCBSaV/Ml6MVHBbTIho9NIR2r30kFKrGxHlzrssoEgq1o91rEFNRaEmxkLLX\nSg2jNdeH6VFGoQG0ew2h0UOLaPeaQkKhPYRUK7JGj4RCS1iC0A7my9Ek5qSaYK01XzR6aBLtXv1o\n9NA0QqoF0uiRUGgSc1KVyxo9oEnMSTWCMgqtot2r2U8ZRUKhYYRUtdIViUkoNI2QqhgLo3AGTJzX\niUYPp8HEecVIKJwB7V6FrLXmy0x/SCicASFVGxZG4WQIqSrR6OE8CKmq0OjhfGyRCWxry2y3etba\nyZBQKKLUx5ZKqh5SRgEnU6ySmn2FwuoxyigUVaqSYp1UJUgonBXtXg1o9HBiTJzXgDIKCjBxjjtI\nKJwbIaUdCYWTI6R04wL0OD1CSjEaPeDIkFqujQKAh44IqRACCfUyyijAGHPkEoR8/yVLEB7gQGLo\nU9kSBOdcCCHdDSE459IXQ2aDMQI4sXcOi/HexxjT3RBC3/dd1xlj+r6fpsk5t9X4zogyCsi8XL+l\n2aVxHCWM5G8pmqSe8t7ffCLt3mPWmmmy9l9CCtpU0+5N0zQbaIwxlU7OuWEY7j0xv2tXvTqqlnCh\nKiBXbJ3UtKrUqArjep/AQrGQonSa+2n0KKOA3Abnk5IpcxFCyO+uOG+59AhlFJDboJLK56FijDdn\nzZeopK5QRgF3bFBJee/TmvKu655cf0AlBeAZm+1TlPWcz26VJQi5rIyi14Na1Z/j/NUFnLMW77yZ\n9bNHD8BNXIhBBWajgHs4n1RR12UUvR6wVKySot37XbpJGQXcR7sHQDXavUIWZRS9HnAT7R4A1Wj3\nSmA2Cnga7Z4K9HrAPYTU4SijgFcQUgBUY+L8WLfKKHo9YAUT5wBUo907EMcSA68jpAqg1wOeR0gB\nUI2J86MwZQ68hYlzAKrR7h2CBZzAuwipYuj1gGcQUvujjAI+QEgBUI2Q2tmdBZz0esCTCKmD0OsB\n72Gd1J44Dgb4GOukjjAro+j1gOfR7gFQjZDazZ2VB5RRwEsIKQCqEVL7oIwCNkJIAVCNkNoBx8EA\n29k9pKy11lrv/d4b0o9eD3jDviHlnLtcLtM0DcMQQth1W1qwgBPYlN11UWUIwTkXQuj7Pt+Qtftu\ntySmzNGoUh/blyspCZ10V2IofTFk5MHe+77vx3HcasQATuW1w2K89zHGdFdKpK7rjDFSKznnZo8P\nITRbNC1RRgFbe6F+S4cEj+MoYSR/p6JJ6qabT8mfZRpu9wgptKuCdm+aptkQY4wpd5xzwzDcfIqY\nFVl21RvfiR4kFLAhzoIAQDXOJ7Uda800sYAT2NZHISVT5iKEkN99qIVUegK9HvChjxZz5vNQMcaX\nlpW3NAkFYD8fVVKywkAipuu62dT4uiYrKabMgc1tsE9R1nO+ttX2liDcmpAipNCSUh/bDSbOX00o\n0eDE+TUSCtgESxAAqMb5pDbD4gNgD8VCqqm9e4vTs9DrAVuh3QOgGu3e9iijgA1xWMzHOBoG2BPt\n3sYoo4Bt0e4BUI2Q2gC9HrAfQuoz14sP6PWAzTFxDkA1Js4BqEa79ykmpIBdEVIfYEIK2B8hBUA1\nJs4BqMbE+UfShBS9HrAT2r13LU7PAmAPhBQA1Qip99HrAQcgpACoRki9hQkp4CiEFADVWCf1Jiak\ngGMUq6Sma6WG8Q56PeBAtHsAVCOkPkKvB+yNkHoHp2cBDkNIvYgJKeBY+4ZUCEGuou6933VDRdDr\nAQewu+5Zc84557z31l5taHa3JtmlQAkpnEqpj+3ulZTkVNd1u27oSExIAUfaeDFnCCHdds6Zn2Jq\nGIZtN1SGTEjZf0uPAziRlysp51yeRKlWki+6jDHGWhtCYEIKwNteCynvfYwx3Q0h9H0vt9ON3DiO\nfd9ba5tp9+j1gIO9MBOWjrYbxzG1cuanxUtz5E++VJUT59Za88UhezinCibOlwfZxRglp4wxr048\n2VXPv85xshVSJBRwGC7EAEA1Vpw/J1seVXoowLl8VEnl0+EhhJdmx2s/nxSAY3xUSeXzUDHGl5Ya\n1HQ+KS6nDpTzUSXlvZej84wxXdelSfRnVFNJcUQxUNQG+xRlPedrW61lCUK+R4/zBePcSn1sC211\nschAaWYtyigSCqdVwTqpbVUwJ3VdRpUdC3BaLEG441ajB+B4XNLqljsJRa8HHI92b43OGkrpYUM/\nNA9P89iM+uGVQru38PNO0ZlQwNkQUtd+Gj0SClCCOanMzwF6xphlQjEhBRTBnNSVewUUCQWUorTd\nW59BfDi/+M7TsxPa/X7NWPmTJ9SHW39zeE/be+tnHt7ev/qHCnwuFCjW7ukyS6jJWGONMVRPQHHM\nSc1PCmyMMVbHHBkA5qR+RvMntXVUT4AqWg4wLka+ezXDATQ70VkQAOBJSvfuAYAgpACoViakZtdq\n1+OlS5weLF3RXufwEkUTjj+cc+mSjq+eRfYY3nudv1nv/eyamMd/cgssQZhdq10PuRy8XF1iGAZt\ns3V936fhhRB0przChDLGxBhfPQf/kZxzMcbL5aLwN5v/0EIIZT6507HSdsdxPHjT6y6XS/ppjOOo\nbYT58PLbqsjAFI5N4ZASebOl213XFR3Omq7rLpfL8ds9ut2TrR680Wc459JnLH2l0Fhu8N7Lu1n+\np33pEofHCCEMwzD7GWoghYl0fKp+p0KGF0KQXk9VGZWTMqpMQ3p8LkpIqapTcvIxK/I/xkMpAhT+\n9OSHltcFSqRfqM5CL/1/IzfUVlIFPxSE1BV5u+hMqETezaVHcUWGNI6jBIHOX+6kspHPf5sKI16U\nnWEgpH5p/n9MCgG5rfCtvKzQ9fx+x3HMB6NqbNPi869teKLUbJRgndQ3abblmsza9rCIdEV7hWNL\n76cUoHpmf0IIfd/LD01+y3rGZn6GpHZ4IsZYclRFotHo++9iORWteYTaxpYorPIm9T+6fG+DwuEV\n/51y7F5N0o6qwuOok6yGLT2Ku5QPryBCCoBqzEkBUI2QAqAaIQVANUIKgGqEFADVCCkAqhFSAFQj\npACoRkgBUI2QAqAaIQVANUIKgGqEFADVCCkAqhFSAFQjpACoRkgBUI2QAqAaIQVANUIKgGqEFADV\nCCkAqhFSAFT7P1QPpGWiq9j6AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "haeff_up  = TH1D(\"haeff_up\",\"effective area\", 24, 1,7 )\n",
    "f.roottree().Draw( \"log10(mc_trks[0].E)>>+haeff_up\",\"(mc_trks[0].dir.z > 0) * w[1]/mc_trks[0].E\" ,\"goff\")\n",
    "haeff_up.Scale ( 1 / ( ngen * log(10) * haeff.GetBinWidth(1) * 3600.0 * 24 * 365 * 2*pi )) # nb: 2 pi\n",
    "\n",
    "haeff_dn = TH1D(\"haeff_dn\",\"effective area\", 24, 1,7 )\n",
    "f.roottree().Draw( \"log10(mc_trks[0].E)>>+haeff_dn\",\"(mc_trks[0].dir.z < 0) * w[1]/mc_trks[0].E\" ,\"goff\")\n",
    "haeff_dn.Scale ( 1 / ( ngen * log(10) * haeff.GetBinWidth(1) * 3600.0 * 24 * 365 * 2*pi )) # nb: 2 pi\n",
    "\n",
    "c3 = ROOT.TCanvas(\"c3\",\"c3\",400,400)\n",
    "c3.SetLogy()\n",
    "haeff.Draw(\"Lhist\")\n",
    "haeff_up.SetLineColor(2)\n",
    "haeff_up.Draw(\"Lhist same\")\n",
    "haeff_dn.SetLineColor(3)\n",
    "haeff_dn.Draw(\"Lhist same\")\n",
    "\n",
    "c3.Draw()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Of course, we now need to check we did everyghing correctly by using the effective area to\n",
    "compute the number of detected events and comparing with the direct $w^2$-based calculation. It is a\n",
    "simple matter of figuring out how many neutrinos there are, during a year, in each bin and multiplying with\n",
    "the effective area."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning in <TFile::Append>: Replacing existing TH1: hcheck (Potential memory leak).\n",
      "Warning in <TCanvas::Constructor>: Deleting canvas with same name: c4\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAF0CAIAAAC/mplpAAAABmJLR0QAAAAAAAD5Q7t/AAAV2klE\nQVR4nO3dbZKbOLSAYenWXVcQM5Wt3KxiGphdJFtJuoHZGPfHmWhkIdMYI3QM71NdKQd/IAwcS0Lo\n2GmaDABo9T+lCwAASwhSAFQjSAFQjSAFQDWCFADVCFIAVCNIAVCNIAVANYIUTu7Xr1/WWmvtr1+/\nDltp27bOucNWd24EKZzc169fp2n6+fPn169fj1mjc67rus1v3xxVzxoZCVJ7GobBOTcMw6Nv9IfX\ntrcve3t7q6rq/f19/tQwDJtX17btwlPyrEsJX+bX3ratnJnOuYVPFvNvSb5A+f6jF//8+dMY888/\n/3z58iX5afIuWfuaVUermO/0YRiapln+nAXbouqTkVG1CZ+pqqppmjWvNMZUVdX3/UOf3/e9MUZW\n4R/syBjz5cuX5FNVVVVVteEz5SRMbqlsjl+1fCdeuHXzlzVNU1WVPF7eoujQ9W/xpZLPkdV9+fLl\n3jcgGyKr9pHl0w1PLgnXKBuy8DkLfv78OU3TX3/9da/MC2XbvFLNzhmk1oeJe68Ml0fnTN/39941\nP2/f398/XWl0Vj8fpMLPlybDr1+/otfIGZUMB59+e1Lge0EqDHwL4cZv9b3T/uPjI3xxuK7oLfJR\n8oL5pxlj5LS/V4boC5+XOdrM6C3JbVyOF+EH+scro+p0+23sFRk1O1uQ8uePP3RkL/oXhHvUv1KW\nyHv92Rv+PIY/1PcqodGnhZ8joSr65HmBoypVeCD6tfvTIxnOwlXI2/1/ozjV972vRCQLs1Ch8N9S\nMkiFyxeClNRc7m1I0zQ+SCWLFL4rCkxRwRa2KDo2RBgCwoMkigv3VuffON/wcO9I+ecfPq2Iqvc2\nhyD1GsxtxJEDzh8EyYX+EPdBapo1webxLnn6mdu2RhibpuBHb01NKhmkfLGTzbQwcPiV3qtJJbci\n/O+9+OI/OXlyhtsyLTb35O0fHx8mqDR93Jruf+cL+yIZ9dZsfh+Y7lfQ7j3+VLiDwtg0L/+9MJQ8\nPj2C1AvwHSViXgEJ40VVVf6VPnaY21pAFKTk8xdOgPCwe3t7k4XSY/3+/p783Z4eCVJT8GucXLsv\nm98WCVILBfafHEZwX8+aglN3uj1v79Ug1gepaZqiIBXWNcJ9FxVpmv20hCVJRvCVmx9Gh+RB4t8Y\nHRtrhAdA9HjlhyxUYE/shEEqFAaX6fbYivgg5T8teSD6VVSpDvLwBA57o3YMUmF1b772qK9kQ5CK\nhN+Vf0HYvT3/MY++w2RRfVtvXuzwjfeKFL5xXpd56JyfvzI8WuYHiX+juV+dXLO6e4/XfEh4bF/B\nCYcg+G3r+95fCDe/L5nLv9Fh4V+5TC42y+vHcazreuHF4zjKA6ksPMN/lDFGVjqO4/LQgc0DC8Lv\nJIra/tJ713VSpE+LcU/XdeGF/HvXzuU1YRQIw1ZVVV3XdV1X3fb1rOecm2+C/+/CQRIeV0cOTZIh\nI/4IPOWoqIQMga8kc9snFf5cm1mFqE/1Sc0/yqT6pJLVIjPrEppmfVLzMt+rSVVB70NUJUx+VHKl\n62tS0X/NZ7/VZl2fVNhi8s3G8DUSxKuqkkafj+m+uXevSD5gRcUwq/ukptu+oem2NZ08SKI3PlSd\neb4m9dAOOo2zBamwxRfuwvCAE+HPbxjX/AvMbU/WvPWx0NybbhsLUQ96JBmkoqarX3sYdufHaLhR\nUpKHgtS8X2ZB8htY+BLCbzsqeVTZlEZceDbeK9K9hQ/12kQVsXB3zA+S6F2PrujJIPXQDjqNc27n\n+kNn5SujFsfKd72/vyfHSa23oYt0ffEyfcJDPSyRcGzUtiI9dLktWsWGp4p4fhe/Fjvd+bkDthmG\noa7rUseVtbbv+6t01lzDCTvOUZZzrmmaNRcidjcMQ1VVRKiToSYFQDVqUgBUI0gBUI0gBUA1ghQA\n1QhSAFQjSAFQjSAFQDWCFADVCFIAVCNIAVCNIAVANYIUANX+t8harbVF1gvgGUXmIygTpEyhrQWw\nWam6Bc09AKoRpACoVqy5F1Udaf0BSKJPCoBqNPcAqLZzkHLOOefCSfg/Pj5k4TwHOgB8as8gJd1M\nzrmu6+Txx8fHn3/+KWkU//77b+IUgIftlcAvTMooeVb7vg/zRL69vZkgT+9e6wVwjFKn7W4d5/Nk\nZ845CVLy33Eco3zWAPC53cOeVKmapgkXSnjyOcePLxWAJ5U6MXcegiBdUWECW+mWkgj1xx9/hGFo\n31UDOKU9g5S1tqqqYRj8EolQb29vXdftuCIA17FbmvW2bbuuky5z4ZyTsCVd5kIqU9aS3v1CrP0x\nTf9XuhR4VqnTdrealFSg6rr2SyRgjeMozT1BbALwkJ2DVISQBOBJZepv84lpCGfnsX7WIXb6SynV\n3Ct27150lbFUMZDFNIV/1nyPliTDk7U/ji8p9OMGYwCqFZuqBYgahpMxxn5LvIyK9rVRk0Ihswbg\nylYhroYgBUA1ghQA1ZjjHIBqzHEOQDWaewBUYwgC9jcblvktOQqdu46xBkEK+4uij7XEI2xHcw+A\nalzdA6AaV/dQRup24kTXFUcJ6JNCGfNeqnTXVfJuPlxJ9gzGwzDIwuSseACwLG8G42EY/ITC4czC\nALDSbkFKak/DMLRtK7Oby2PJHzMMQ1VVYQ0LANbYLUg55yQtaLhkHEef2VhqWHutDsBF7JlmXeKR\n5LaKAtbcfJrzENf+AIjsGYzvIQwBWCNvBuOqqvxj6ZbacXUArmDnjvO2bYffzG0/1DiOdJwDeFTe\nDMYSs6QNWFWV70THpXB3MZ5xRLY/Gc95s9ZCWQZxBGu3JVCw9kdqxPnGT8PuSp22R9wWk6xAcYMx\ngDW4wRiAaswnhddDQvZLIUhBC/rXkUSQAqAaQQqAagQpAKqVGfgwv7uYi33nse/IpsUb0W9wCGVW\napxUsZrUdKtUMaDdNM3/rPkeL8R50dwDoBpBCoBqBClsx6BKHIAgBUA1MhgDUI0bjAGoRgZjrJMa\nrzSZ7BmG73R7xQnZ+cU7sf1HZ0kGYz+HlGTfM8aEC5n07vWkhmimp6k7oiyz9TI3Xn4nGczZtu04\njv6/UQZjMq0DeNTO2WKiJT6DsTwrU57vuEYAp7dnTWp+g4skMQ6Tx+y4OgBXsH8j01rb933Y/eSf\n8utaTl9suPanEH1Sl3eSPqmIc66qKqlhRSmtpkVZSwXgheQNUuM4+sAkTb+sqwNwPtnHSfmectKs\nv7rUkKV4vJJgtnLsKG+Q6vu+rmvfA0U77qXNQ4+1xCNkt3+QCiORc26aJi7tAdiMDMYAVOMGYwCq\nMZ8UANUIUng99NZfCkEKgGoEKQCqEaQAqMYc5wBUI4MxTouMW+dAcw+AagQpbMdQAByAIAVANYIU\nANW4ugdANW4wBqBamZnVSQ76epRnOvgstcd/NG+FbudJxOCci5LrOeckrfHu6wL+NU3zP2u+xwvx\ngvJmMDa/+56cc13XkRkUwKOOy2BsjBmGgUmEATxkzyAl7dUwVElUkoaeNPp2XB2AK8h7dW8cR2n9\nVVXVdd29zMZzdKvjIXdu04szbnFUvaLsQxB8c885V9e1jz6EIewoeYNOIuOW/XZQgbCf7CPOwwzG\nudcFbMBkCcrlDVJN03RdJ4+7riODMYBH5Q1ScnXPWis9UAxBAPCovBmMze/AxOADANscdO/ePEJx\ngzGANbjBGIBqzCcFQDWCFADVCFIAVCNIAVCNSe+wjvJJ71ZaPzeeYXq82HkmvVvJ3ipVDFzLyrnx\nCE+aMAQBF7JysgTDfAma0NxDgrU/ZvMHnKK5l5LYWHPm7d3scs09AFiDIAVANYIUANXIYAxANa7u\nXV5q/MdkmGkXWtDcQ2L0EEOHoMfOQWqevtgY07YtM94B2GbP5t48fbHw05wDwKN2C1L3bm3hlhcA\nz9gtSM3TFxtjJHdxmDMGCqVuFkndKZLKbQfklvHq3jAMXddN0yShKkIGYz3m0SeRVvO8rrOlLyrj\n1b26riV9sU8YEz47LcpXKgCvJe84qXEc67qWx2GOdQBYKWNNyleLmqYxtOAAbMJgTgCqMZ/U5aUm\nTkpPsXQpzCc1U+q05QZjAKpxgzEA1eiTAqAaQQqAanScXx49xEl8LTMkYgCABIIUANUYggBAtWI1\nKe4oxmu5k/0Y2dHcA6AaQQqAagQpAKoRpACoVmwwZ7SEvvNiGLWYtD6ByGW+vcsN5uTqHlQjYaoa\nNPcAqLZ/kIqSGA/D4JxzziVzxgDAsp1HnEdJjIdhqOu6aRpJb2V+Z+IDlFubi/CQwlzcnj1hvju8\n73vnnDFG/pWKVdu2kobPMAuCKnScr5OeUvlK394Zpg+eJzGWhp48HoahqqodVwfgCvLeYOwbd865\ncRz7vvdPkcEYwBrZZ0GQbikTtAEFYSi79SN9AMXyDkHwHefTNIURCgdhsA9eX96aVF3XVVWFgxII\nVQAekr25N46jNPcErTwAD9k/SIVhiJAE4EncFgNANeY4B6AaadaBVRLDzXGIYkEKB1h7AxqnHxQj\nSJ3ZPPpYS0jCi6HjHIBqBCkAqnF1D4BqXN0DoBrNPQCqEaQAqEaQAqBaseSg9Elld6Xpt4u50pd8\nueSgALAGQxAAqMYQBACqHZfBOFwIACvtHKSSGYzlcTiJMHBWqZkn8JQ9m3vzVHpt21ZVJXUo51zb\ntqRZB/CQPWtS0zRFPU3jOPr0MM65rut2XB2AK9BydS9CtzoAwdU9AKrlHcxZVZV/PAxD+F8AWCNv\nkAr7ocZxpNccJ2Rt+DeZb9GSf/+wVd4gJVf3rLXWWsm3nnV1wNGmKfqz5vt8YelSvra8GYyNMcMw\nyHjO3VcEFLc2H88hhTkrZkE4ryvdoK+HtT8S+XhOsS9KnbZahiAQswAkMQQBgGrMJwVANYIUANUI\nUgBUI0gBUI2rewBUY5zUeZ1ibM5JnGJfkC0GABIIUgBUI0gBUI0gBUA1ghQA1RiCAEA1bjAGoFr2\n5p6kL2biYADb5A1S0qaTmc6ZnBPABhmbe5K4WP41xpAZFMAGece5W2ubpnHO1XXt860bbos5xilu\nxTiJU+yLUqdt3rU658ZxlMd93/sW33L6YkO3+i5OcWKcxCn2xQnv3WvbdhzHvu+naWqapq7r8Nlp\nUb5SAXgtGUOj1JvCJp6vTNHcO8Ipfr1PIrUv0nllFDthTUraehKkZAgCF/gAPCrj1b22bYdh8K28\nvu/zrQvAWeUdce5HIVCHArDNETcYE6EAbFZs+uBoCf3o+6PjXI/Pxtz8R/EuO2HH+TLGHOBCpmn+\nZ833eCFSis2CAFyHtT9Si79FFSyiVBLZYs4gPeKG5p5uib2me5ddrrkHAGvQ3HtB88sOxhj7rURR\ngOyK1aTsrVLFeFWfdsHSEYuzYPpgAKrRJwVANfqkXtLsknZ8MVu81k32QBJB6iVF0cda4hFOi+Ye\nUAa/KysRpACoRgZjAKoxBAGAatmbe23bksEYwGZ5g5TPXUwG46zogsWJ5c1gPI6jNOuoTAHYJmNN\nyk9wLuHJ57YCgPWyB6m6riVnTNTcs4vylQrAa8l+dU+ae2Fuq3A5ACzLmxw0ekyLD8CjMgapsCuK\nDMbANnfmR7+QvM29pmnIYAzgGUfMrD7PYEwihqfonq4f26X2bDrLRglnTsRAKw/AZtxgDEA1bjAG\noBrJQV8QfVJntX4Y85VOWya9A9RI5SVL5Cu7GIIUANUIUgBU4+oeANW4ugdANZp7AFRjCMILYgjC\nxRU6ABiCAGC7E0+WQJACoBpBSrUT/zwCKzEEAYBqDEHQZHbr1mSMsd9KFAXQ4qDu+ui6AFf30mZX\nbfRMeAZFHkqntN+JVuq0PaImRYqq9WadUN+SXx6R68qs+Z5a/M3Mlk/mDNXw7EGKxMUPiaKPtcQj\nxJKHRPpQOUVfQd6re8MwdF3XNE3WtQA4sbxBqq7rpmmSc5yTwRjAGtmTgzrnJPVelBl0WpSvVABe\nS8bu+nmFqO97iVxc3Uvjpjxslb4QvOsRVeq0PSjvXl3XDEH4HEEK+zpFkOK2GACqMVWLJtSksC9q\nUgCQGzcYA1CNG4wBqEZzD4BqBCkAqhGkAKhGkAKgGlf3AKjG1T0AqtHcAxBTlaaIIAVciKrosxJB\nCoBqBCkAqhXrOAegxbqEj6UudRWbqiVawsU+Y5iqBXtbny5gRcLHy03VwqTmQHbTFP1Z832+UPlP\nY/Yg5ZxzzpF9D8A2eYOUNOucc13XkagKwAYZg5TUnoZhaNu273szy2oFQCdVebMzXt2b5wRNZgm9\npnQCIiCzVzzqjuiub9tWkq37nqlPm36n6kp/qJ17pg3HuZS6upd9nJTEozBCiVOFoU/NNpaaFLBS\n3iBlra2qiq6o1A1T35IVLCIXEMkYpKTq1LatD1KX7ZNKjYsjHgGrZAxSEpvquvZL+r6/bJwCsA0Z\njPPjZhecwuVuiwGANZjjHIBqzHEOQDWaewBUI0gBUI0gBUA1ghQA1bi6B0A1ru4BUI3m3s5eMfki\noBm3xTxhayoO4BVxW4wiD0zHvjIbR6nilaC5eJrLZtQXrxSCFADVaO4lrC1eoekNTvLtlaC5bIbi\n3cEQBACqMQQBgGpKMxgv9yB+2r/45Nuf+fDnX1C2eM9v3ZNv11y83Lv+U5rPi3yOy2BMpnUAG2Sf\n49xnYXj1OEUSKqCIjDWpMJPVa+VfmIwx1kZ/k/n2+8HvZwHkl/GaogQmH6qstf9li1F/hmsvH1DC\n2YYgSFeUPJZQ9V99Sv2lPe3lAy6DEecAVMs7hNQ38aKmHwCslLcm1TRNXdfW2nEco350nQFr26iu\nYwzDoLl4nsLhNs45+5vOazht2+rcs23b2lvHn7l5R5y3bdu2rZxd4cJxHLOudxtrbVVV0pXWdZ22\nMfF1XfviDcOgM8orjFDGmHEc5asrXZA059w4jk3TKNyz4Zc2DEOZM3c6ll9v3/cHr3pZ0zT+2+j7\nXlsJw+KFj1WRgiksm8IieXKw+cdVVRUtzpKqqpqmOX69R3ecy1oPXukazjl/jvklhcqS0LatHM3y\nS1tVVekSxYZh6Lou+g418FeWdbb1/JhnaeupqkaFpBpVpkF6fFyUIKWqnhKS06zIL8anfAhQ+O3J\nlxbWC5TwO1RnRc//3sgDtTWpgicFQeqGHC46I5QnR3PpUtyQIvV9L4FA586dVDbkw72pMMSLsj0M\nBKn/aP4dk4qAPFZ4KM9r6Hr2b9/3YWFUlW2anf/aiidK9UYJBnP+Sxrbci1S2xUWEQ3fV8UfTz6A\n6un9GYahrmv50mQv6ymb+V0ktcUT4ziWLFWR0Gj0/VzMu6I1l1Bb2TyFtbxJ/VcXXm1QWLzi+1T1\nnMqIxLdA4hHReD1tlBevIIIUANXokwKgGkEKgGoEKQCqEaQAqEaQAqAaQQqAagQpAKoRpACoRpAC\noBpBCoBqBCkAqhGkAKhGkAKgGkEKgGoEKQCqEaQAqEaQAqAaQQqAagQpAKoRpACoRpACoBpBCoBq\nBCkAqv0/EGYGyVBohY8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "hcheck = TH1D(\"hcheck\",\"event rate computed with Aeff\", 24, 1,7 )\n",
    "\n",
    "for b in range (1, hcheck.GetNbinsX()+1) :\n",
    "\n",
    "    E    = 10**hcheck.GetBinCenter(b)\n",
    "    flux = 1e-4 * E **-2\n",
    "    Aeff = haeff.GetBinContent(b)\n",
    "\n",
    "    # -- how many neutrinos in this bin in a year -- ?\n",
    "    \n",
    "    # How many GeV's in the bin?\n",
    "    Ngevs = ( 10**hcheck.GetXaxis().GetBinUpEdge(b) - 10**hcheck.GetXaxis().GetBinLowEdge(b) )\n",
    "    \n",
    "    # How many sterradians in the bin?\n",
    "    Omega_bin = 4*pi\n",
    "    \n",
    "    # The number of neutrinos in the bin for a diffuse flux \n",
    "    Nnus = flux * Ngevs * Omega_bin * 3600 * 24 * 365\n",
    "    \n",
    "    hcheck.SetBinContent(b, Nnus * Aeff )\n",
    "    \n",
    "c4 = ROOT.TCanvas(\"c4\",\"c4\",400,400)\n",
    "hevts.Draw()\n",
    "hcheck.SetLineColor(2)\n",
    "hcheck.Draw(\"same\")\n",
    "c4.Draw()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "The effecitve Area can also be used to compute the response of the detector to a point source of neutrinos. The difference with the foregoing is that the flux of a point source is expressed in $GeV^{-1} s^{-2} m^{-2}$ (note the absence of $sr^{-1}$).\n",
    "\n",
    "The following computes the response to a point source flux"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning in <TFile::Append>: Replacing existing TH1: hpnt (Potential memory leak).\n",
      "Warning in <TCanvas::Constructor>: Deleting canvas with same name: c4\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAF0CAIAAAC/mplpAAAABmJLR0QAAAAAAAD5Q7t/AAAVWklE\nQVR4nO3d7ZarqNaGYXjHPq4G0+dVQp9XFXhkvj/mXmxajMF8lDPmvsYaPewYI9H4BAgldp5nAwBa\n/d/RBQCALYQUANUIKQCqEVIAVCOkAKhGSAFQjZACoBohdVvO+Rc2+QWdpQohWGs7X9NaKy9bXrw8\nssF7b60NIfQ8Lef88/NjrfXer67tLOqzeO8XJanpPPXvbcYtxpiU0q5NnHPjOL6kNPcax7HzdPc/\nU54sC+Uo9Rwu+ew5524+TV7fOdc++Y7z8hQppY39HlWqEzttTSqE4L0v39X1V984jj8/P8aYnLN8\nK5ZvvxDC4kF5hcV3vjwthCCPj+Mom5SXnaYpxliqGIu91Mra9pFrpWrf2uom5QXlxeu30x4fY0xK\nafFIXZ6y1WK/8oTFUbr5luWgTdNUP14OY/0WpBo1TdM0TXJ467WLt1O2lSLVJa83vFYVWj2SpjlW\ni/Miazc+Las7vXaW6+X2Xdws2zkdnZIvYYyRuoz583VtjPn6+iprv7+/U0qy1jlnqlqAMWYcx/Jg\n/SKFrJUnfH19yYvLg+WVZe3qXgp58bK7eZ63S1XOWtmwLXN5EdlFWZa1ZZPF8SnL5a0tjmddkvpl\n26PUHsPF2SnVolKeeZ4Xh7Eu0vf3t7zm9/f34tDJ65S3ttj22omrn7l4m+3a9ljV5e//tJQibZzl\ncqzqV6sPbNnk2nk8pROGVN1aqS9ROYtyMcz/bkGsXjarn57F8+UFy8VTorA091b3UrS7qz9t5Y2U\np9Xpc63MzrnVkCqvtnp86hdpL2BZW4dvXbbFUVo9hqvvepEF7WEsm2839xaX9ziO11qs8i4WR7It\n2Pzv639xrPZ+WhbndPF22iNZluu9SxTWy6tlO6X/mNORynDd9Ztz/vr6+vvvv40x//zzj+SUtDXa\nHuKNPtFVMUbvvbx469pe2t3Nf76Zr7VENspQ1o7jGGPcLvDq8THGlDqU975tiMnLynNyzjHGumbX\nX1ppmMQY63JKO64+hosC3Hw7ZY/OudVWXilYCGHjxdsjWR8W0/SL93xaUkrDMFhrnXPy9qdpKkfP\ne3/zlC32WzcqF+dx76f3LZwwpOQaK1/79ZmTT0b5TEgnyN17kQX58H1/f18ul9Uw6tyLdC7cV5gi\n51wuqms2js/GJhIr9QG8r7SlQ0r+N8YYQvjrr7+MMaVZ96CNgg3DMI6jXN7DMGyX0znXmSA3yzPP\ns3QwDcNQjnzZ0d4XrDfZdR7f1dFVueeT0yY157pKvOhtafsd5rUq93yluVfX86VtItdY29xr91KU\nR0qZt0t1rblXXtZU7bLSUDL/buWtHp+6AbL6wZAHSy9bXci2udcew/p16kdK8RaHUZp+pru5VzeH\nFw2loj56q0/YPpLtser8tJimVXjtLG+8i3rZVI369nN+Pud8Y3VLpHxi6hARbVhf+9i1l9aix3fx\nUm3H9uqVX3+pLjJitVTXQmr7RUyVTXXxFu+r3aR2LW3bo7R6DOvz0h6BlFJdjWpffDuk6sMox/za\nRVvvwqx9Z7RlaI/VzZBafdc3z/Ki+dxGbV1Hvla2UzpnSIme0yZVg8dfqnT6loV6k+29tKs6SyXM\nnwrOYpPVV1iUarH25+enc6erHr9O6t/vXrT3jSNw7Uj2v/jGkzvP8s0dPVi2d2Rn7sz55qy1KaXT\n9kf8Io6kTqcdzPk5pH/36FKcAUdSJ2pSAFSjJgVANUIKgGqEFADVCCkAqhFSAFQjpACoRkgBUI2Q\nAqAaIQVANUIKgGqEFADVCCkAqu0OKZlUZ2MKncWUO9uzGwHAtn0hJffwlhs/r+bU4i73OedyJ+lh\nGMgpAHvtuFWLJI48P4QQY1xsW6YhKHcOqyfYkAkLyCkAu+yoSdX5snpvMLnXZ/2ITEwic7Fe2woA\ntvTfaXhxM3xz5d7vi8dX93XMWwXwmP64eKId8+7Vc5B11oy896WJt+g+J6qA97Ixx+1LvXYIwjRN\nJcjaeXEB4KYdISU/50lVKIRQZtzeGI5gqp6snsl1AWBh3zTr4ziWIQWlvRZjlKZc+/yU0jAMpZZI\nEw/AXvfMFrN31vm2A8taZqkB3sxRl+1BeyWkgHdz1GXL3+4BUG1fn9QTLX7OpGIFYNVhIUUqAehB\ncw+AaoQUANUIKQCqEVIAVDus4xyn1f9nqPx4gg6EFF6gJ30O+pN6vB3GSQFQjXFSOE5nZYqPymej\nuYeDdEYPrcKPx697AFQjpACoRkgBUI2QAqAaIQVANUIKgGoM5gSgGoM5AahGcw+AaoQUANV2h5TM\nA7oxa7H3vsxa3LkJAFyzL6Skt9t7H2NcDZ0QwjRN1zZZhBcA3LSj41wipgRNm1O2+VvQEIJzrmyy\nd+pjANgdUkJqRosnyA92dVRJKkmWSaPv/pIC+Eg75k2WiClRZa1NKbW5Uz9eAss5N01T+/g1DFB4\nY9Y+8w5Qz301POANpln33pf+JomqnpqRc26e55yzc24YhvL4vGnXewBwYi8fglCCjLYegHts12ja\nCk5KaZ5n55xUkVJK4ziuPmee53Ecyy7Mn1qVLO/aL97Jc08uHxU1jrps9/1ZzDiOpck2/2mUxRiv\ndYqHEHLOpQeKIQgA9rqnJ+yOkQSLTY7qgcNvoOP8pI66bA/aKyF1YoTUSb3Br3sA8Pu4nxQA1bif\nFADVaO4BUI2QAqAaIQVANUIKgGqHdZwDvW7dM+O/+CnmpAgp6NYZPZ1BhjfEOCkAqjFOCoBqdJwD\nUI2QAqAaIQVANUIKgGqEFADVGCeFPgxEwkEIKXRj1AiOwGBOAKoxmBOAanScA1CNkAKg2r6QkklA\nQwgbT2hnAA0hMMc6gPvsCCnp6vbexxhXcyqEME1T+3iMcfVxALipt+Nc6kelltTmlL0yjuba4wDQ\no7cmVTfiVttu8zy3P9hJkI3jeFfZAGBPSDnnZFlCqu17ajeJMV4bamA3dZYKwOn1hpT3vvQrSTzd\n7AsfhsE5l3NeNBXFvGnXewBwYq8dzDlN0zAMsjwMA+kDYK/empT0LkltKIQgTb+c88ZwhFItkj4p\nEgrAHXbUpMZxLNWikjgxRhk89fSSAYAxxu6t4OScH48ka3fvFwezVvVdEJQX7xSOumwP2ish9XaU\np4Dy4p3CUZctf7sHQDXuJwVANe4nBUA1mnsAVCOkAKhGSAFQjZACoBohBUA1QgqAaoyTwll03oaM\nT9q7YZwUTqHz48T9FN8QzT0AqhFSAFQjpACoRkgBUI2QAqAaIQVANUIKgGoM5gSgGoM5AahGcw+A\naoQUANV2h5RMBboxcbH3XiY6FjJP3/YmAHDNvpCS3m7vfYxxNXRCCNM0lf/NOQ/DIJOJXtsEALbM\n3VJK5fnjOLbbltdMKckjzjnnXLvJrv1ChXOcsnO8i4McddnuqEnVjbjVmdbLO6mfVmpPOWfnXP/u\nAMDsGoJQp4yElPQ3bWxSEsp7P02T1MWE3byzz8wABbxI/y2l+BDqsKMmJUEjy1Kr2k6o8kxrrSRU\n/fztCt6etwB0m+fef1DjtUMQpON8HMd5nnsSDQAW7K5qi7VWKkSluSfqn+3Kc2TZOVevLY9TXXoz\n1n5W/eLT3m+Hoy7bfX8WM47jMAyyXIobYyyx1ZqmqWxi6GwCsNM90Xizv/z2XqlJvZ1Pq1l82vvt\ncNRle9BeCam382kX7ae93w5HXbb87R4A1bif1MdjKjroxv2kwKhFqEZzD4BqhBQA1QgpAKoRUgBU\nI6QAqEZIAVCNcVIAVGOcFADVaO4BUI2QAqAaIQVANUIKgGqEFADVCCkAqhFSAFRjMCcA1RjMCUA1\nmnsAVCOkAKj2tJCS+UHryYrbJ+Scn7U7AB/iOSElveDe+xjjak6FEKZpesq+AHyUJ3ScS/2o1JLa\nnLJMmgTgXk+oSdWNuNXp1+d55rc8APd5Tkg552RZQqqn78luerxUAM7hCSHlvS/9TRJPq/WphXnT\n46UCHmVt1z+82GGDOQHVOr8pCanXe0JNSrrJpQ4VQpCmX855YzgCAHR6Tk1qHMdhGGS5NNZijDJ4\n6im7APCZ7BM7gHLOnZFk7TP3i4dY29u0QeuTjt5Rl+1BeyWk9Piky+z5PunoHXXZ8rd7AFTjflLn\nxQ9POAXuJ3VqHGS8P5p7AFQjpACoRkgBUI2QAqAaIQVANUIKgGqEFADVGMwJQDUGcwJQjeYeANUI\nKQCqEVIAVCOkAKhGSAFQjZACoBrjpACoxjgpAKrR3AOgGiEFQLXdISXzfV6bnbhdm1KSB1NK9xYS\nwOfaF1LS2+29jzG2OdWuTSldLhdZe7lcyCkAu83dJGJkeRzHxbara51zzjl50BhTL/fvF3fiIP+C\nTzrIR122O2pSOeey3E6nvrrWez9NU0pJIqxzEnYAKHbMmywRU8LIWiv9Tdtr6/FQZV/21ryV/aXC\nVZ80A/hhPukgv8E061ItkmUJo7pmtLrWe1+aeM65+vnbFbzH3hSA83jtEIRpmhZNv5fuDsD57Ku/\nlUZcadwJ+S2vXWutdc7V1S5ZPqre+Fk+qSVymFsdF//z/ufiqMt235/FjOM4DIMsl+LGGCWY2rU/\nPz+Xy6X0QBFMOJvOj3R/lqFxTzTmnDd+p2vXyk97Jb8MNanfQU1Kj1Oci6Mu24P2Skj9glNcGCdx\ninPxBr/uAcDv435SAFTjflIAVKO5B0A1QgqAaoQUANUIKQCqEVIAVCOkAKhGSAFQjcGcAFRjMCcA\n1WjuAVCNkAKgGiEFQDVCCoBqhBQA1QgpAKoxTgqAaoyTAqAazT0AqhFSAFTbHVIyD6hMWdy5dnsT\nANiwL6Skt9t7H2NsQ2d1bf2gzLEOAP12dJxLxJSgWeTU6toQgnOuPLg99TEAtHbUpOp6UJs1q2sl\nlSStaPEBuMOOeZMlekoYWWtTSiWPVtcOwyD/65ybpqk8fzFIqsUAhSc4xdTeJ3GKc/EG06x776dp\nkmUJo7o+dW2tc26e55yzc65kljFm3vTQewIUsrbrHxovH4KwqGoBn2ieu/5h1XaNpq3gpJTmeXbO\nSRUppTSO47W14ziWXZg/tSpZ3rVf3IOD/HZ0n7KjLtt9e5XQqdMtpVSyqV07z7Nzrn2QkPoNHOS3\no/uUHXXZ3tMTtj2SYHXt4sGjeuA+yyk6az+L7lN21GV70F4JqUf0d69ykN8LIbXmsLsg4CGKP8rA\nc3E/KQCqcT8pAKpxqxYAqhFSAFQjpACoRkgBUI2QAqAaIQVANUIKgGoM5gSgGoM5AahGcw+AaoQU\nANUIKQCqEVIAVCOkAKhGSAFQjXFSAFRjnBQA1WjuAVCNkAKg2r6Q8t5770MIu9aGEJhjHehibde/\nT7IjpKSr23sfY2yTaGNtjHGapgcLCpzfPHf9+zSdMx3LdOqyLNOpd65d3VH/frGCo/fhDvoAHHXZ\n9takcs5leXUW9dW1UqWS2AKAO+wIKeecLEsM1cG0ujalFGOcr9RO7ab73gyA8+kNKe996VeSeKpr\nTKtrL5eLcy6lJI9Ik7DYruA99J4AnMhrB3NO03S5XGT5crmQPgD2sv3BYa1NKck4A2NM/kM6ntq1\nZcMQwqLdZ+2O/WLJ2k/8iQfFQR+Aoy7bHTWpcRyHYZDlUtYYowTT6loAeNDuaMw5b4zM3F77v71S\nk3oENakP92E1qYP2Skg9gpD6cB8WUvztHgDVuJ+UJgwQAxrcT0oZDgt69H+fvf8n6rCQAnCn/tw5\nRd2cPikAqhFSAFQjpACoRkgBUI2QAqAaIQVANQZzAlCNwZwAVKO5B0A1QgqAaoQUANUIKQCqEVIA\nVCOkAKjGOCkAqjFOCoBqNPcAqEZIAVBtd0jJVKAya3HP2p+fnzJ76L2FBPC59oWU9HZ772OMbU61\na39+fv7++2/nnDHmn3/+IacA7DZ3SymV50vc3FzrnHPOyYNfX1/lCbv2+0E4LHiup36ijrpsd9Sk\ncs5luZ1LfXWtc06yyRgzTZNUqQDoYm3Pv6N+j98xb7JETwkja21KqeTRzbXTNH1/f18uF9MMkmr1\nl+pUmD8dz9X5iep72htMsy5BI8sSRnV96tran58fa22dUGK7gvfQewJQ9NSSdHvtEATpOP/6+prn\nuU4oAL9hnnv/Kbav/lYacaVxJ+S3vHattbbuljLGlOYe1aUVNPeg2GGX7a5u9noMgTwiP+qllFbX\nXtvd3v1+Cg4LFDvqsr0nGnPO7a97nWsFNal11KSg2FGX7UF7JaRWEVJQ7KjL9rC7IHwQ9b+eAJoR\nUr+C+hFwL256B0A1bnoHQDXuJwVANUIKgGqEFADVCCkAqhFSAFQjpACoxjgpAKoxTgqAajT3AKhG\nSAFQjZACoBohBUA1btXyAG4UBbweIbVixx0Iua9pQ3PxNJfNqC/eURgnBUA1xkkBUI2OcwCqEVIA\nVHt5SMmExjLFcT+7+cPZ9trHN3/kxR9/wrHFe/zdPbi55uK9+tTfpPm6eJ3XhpS8be99jHFvTgGA\neWnHec65/NcYQ04BuMMLa1IlnowxNydeV8faq//KWgCv98LBYxJMJaqstSkleVBt6xfAhrNNsy5d\nUbIsUVXqUwySAtCJIQgAVHvt3wqVJt6i6QcAnV5bkxrHcRgGa+00TYt+dJ2Bdd+ort+Rc9ZcvEJh\nh6P33v6h8zecEILOMxtCsP/2+1fua/92L4QQQpCrq35wmqaX7vc+1lrnnHSlxRi1dZwNw1CKl3PW\nmfIKE8oYM02THLqjC7LOez9N0ziOCs9sfdByzsdcufPvKvtNKf3yrreN41iORkpJWwnr4tXLqkjB\nFJZNYZEK+bCVZefcocXZ4pwbx/H39/vbHeey11/eaQ/vfbnGyiMHlWVFCEE+zfJN65w7ukRLOecY\n4+IYalB+WdbZ1itjnqWtp6oaVZNq1DEN0t/PRQkpVfWUmlxmh3xj3FQiQOHRk4NW1wuUKCdUZ0Wv\nfN/Igtqa1IEXBSH1L/Jx0ZlQhXyajy7Fv0iRUkoSBDpP7qyyIV+fTYURL47tYSCk/kfz95hUBGRZ\n4Ue5raHrOb8ppbowqso2N9e/tuKJo3qjBIM5/0sa2/JbpLZfWMRi+L4q5fNUAlRP70/OeRgGOWhy\nlvWUzfwpktriiWmajizVIdFo9H1dtF3RmkuorWyFwlrerP7Q1b82KCze4eeU2SneyeJPILHLYrye\nNsqLdyBCCoBq9EkBUI2QAqAaIQVANUIKgGqEFADVCCkAqhFSAFQjpACoRkgBUI2QAqAaIQVANUIK\ngGqEFADVCCkAqhFSAFQjpACoRkgBUI2QAqAaIQVANUIKgGqEFADVCCkAqhFSAFT7f1Ppa2JtPZv0\nAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "hpnt = TH1D(\"hpnt\",\"event rate computed with Aeff for a point source\", 24, 1,7 )\n",
    "\n",
    "for b in range (1, hpnt.GetNbinsX()+1) :\n",
    "\n",
    "    E    = 10**hcheck.GetBinCenter(b)\n",
    "    flux = 1e-5 * E **-2 # units are GeV-1 m-2 s-1\n",
    "    Aeff = haeff.GetBinContent(b)\n",
    "\n",
    "    # -- how many neutrinos in this bin in a year -- ?\n",
    "    \n",
    "    # How many GeV's in the bin?\n",
    "    Ngevs = ( 10**hcheck.GetXaxis().GetBinUpEdge(b) - 10**hcheck.GetXaxis().GetBinLowEdge(b) )\n",
    "    \n",
    "    # The number of neutrinos in the bin for a diffuse flux \n",
    "    Nnus = flux * Ngevs * 3600 * 24 * 365\n",
    "    \n",
    "    hpnt.SetBinContent(b, Nnus * Aeff )\n",
    "    \n",
    "c4 = ROOT.TCanvas(\"c4\",\"c4\",400,400)\n",
    "hpnt.SetLineColor(2)\n",
    "hpnt.Draw()\n",
    "c4.Draw()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "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.6"
  },
  "widgets": {
   "state": {},
   "version": "1.1.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}