Skip to content

Instantly share code, notes, and snippets.

@mickypaganini

mickypaganini/numpythia_pyjet.ipynb

Created December 19, 2017 22:03
Show Gist options
  • Save mickypaganini/05725807b4aff8b25d4e59c733f641b7 to your computer and use it in GitHub Desktop.
Save mickypaganini/05725807b4aff8b25d4e59c733f641b7 to your computer and use it in GitHub Desktop.
Download ZIP
numpythia and pyjet (aka: Pythia and FastJet)
Display the source blob
Display the rendered blob
Raw
numpythia_pyjet.ipynb
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Generate Events using `numpythia`"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": true,
"deletable": true,
"editable": true
},
"outputs": [],
"source": [
"from numpythia import Pythia"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false,
"deletable": true,
"editable": true
},
"outputs": [],
"source": [
"pythia = Pythia('/Users/mp744/Documents/CERN/external_dl_work/jet-generation/test.cmnd',\n",
" verbosity=1)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false,
"deletable": true,
"editable": true
},
"outputs": [],
"source": [
"import numpy as np\n",
"\n",
"events_particles = []\n",
"for event in pythia(events=10): # generate 10 events\n",
" events_particles.append(event.all())\n",
" \n",
"events_particles = np.array(events_particles)"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"(10,)"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"events_particles.shape # number of events (# of particles per event is variable)"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"array([ ( 6.50000000e+03, 0. , 0. , 6.49999993e+03, 9.38269998e-01, 0. , 0. , 0. , 0. , 2212, 4),\n",
" ( 1.52633227e+03, 1.06979794, 0.45483708, 1.52633183e+03, 3.73762473e-05, 0. , 0. , 0. , 0. , 1, 61),\n",
" ( 4.77244901e+02, 0.46232963, 3.50380663, 4.77231814e+02, -4.38275348e-05, 0. , 0. , 0. , 0. , 2, 61),\n",
" ...,\n",
" ( 7.47550827e-01, -0.04153067, 0.69428582, -2.73997519e-01, -1.05367121e-08, -15.72057326, 72.25061843, -13.65830637, 81.26422797, 22, 1),\n",
" ( 3.04572999e-01, -0.16297696, 0.24896064, 6.49755591e-02, 0.00000000e+00, -15.72057326, 72.25061843, -13.65830637, 81.26422797, 22, 1),\n",
" ( 2.83943705e-02, 0.01789206, 0.00421003, 2.16423190e-02, -4.65661287e-10, -15.72057326, 72.25061843, -13.65830637, 81.26422797, 22, 1)],\n",
" dtype=[('E', '<f8'), ('px', '<f8'), ('py', '<f8'), ('pz', '<f8'), ('mass', '<f8'), ('prodx', '<f8'), ('prody', '<f8'), ('prodz', '<f8'), ('prodt', '<f8'), ('pdgid', '<i4'), ('status', '<i4')])"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"events_particles[0] # all particles in event 0"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Cluster Jets in the 0th Event using `pyjet`"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"from pyjet import cluster"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"array([( 6.50000000e+03, 0. , 0. , 6.49999993e+03),\n",
" ( 1.52633227e+03, 1.06979794, 0.45483708, 1.52633183e+03),\n",
" ( 4.77244901e+02, 0.46232963, 3.50380663, 4.77231814e+02), ...,\n",
" ( 7.47550827e-01, -0.04153067, 0.69428582, -2.73997519e-01),\n",
" ( 3.04572999e-01, -0.16297696, 0.24896064, 6.49755591e-02),\n",
" ( 2.83943705e-02, 0.01789206, 0.00421003, 2.16423190e-02)],\n",
" dtype=[('E', '<f8'), ('px', '<f8'), ('py', '<f8'), ('pz', '<f8')])"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# how to extract the E, px, py, pz of all jets in the 0th event\n",
"events_particles[0][['E', 'px', 'py', 'pz']]"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"sequence = cluster(events_particles[0][['E', 'px', 'py', 'pz']],\n",
" algo='antikt',\n",
" ep=True,\n",
" R=1.0)"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"jets = sequence.inclusive_jets()"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"[PseudoJet(pt=9514.937, eta=-0.143, phi=2.402, mass=1084.248),\n",
" PseudoJet(pt=7942.452, eta=1.313, phi=-0.807, mass=2062.356),\n",
" PseudoJet(pt=622.976, eta=-1.129, phi=0.784, mass=330.919),\n",
" PseudoJet(pt=426.946, eta=1.328, phi=0.694, mass=204.502),\n",
" PseudoJet(pt=148.835, eta=-2.805, phi=1.531, mass=95.795),\n",
" PseudoJet(pt=139.375, eta=3.853, phi=-0.467, mass=93.964),\n",
" PseudoJet(pt=132.775, eta=4.957, phi=3.118, mass=86.112),\n",
" PseudoJet(pt=118.574, eta=0.016, phi=-1.573, mass=78.424),\n",
" PseudoJet(pt=98.295, eta=1.562, phi=-2.489, mass=54.040),\n",
" PseudoJet(pt=81.029, eta=-2.773, phi=-2.470, mass=56.464),\n",
" PseudoJet(pt=56.672, eta=2.057, phi=2.253, mass=36.112),\n",
" PseudoJet(pt=56.206, eta=2.506, phi=-1.192, mass=27.134),\n",
" PseudoJet(pt=47.054, eta=-1.616, phi=-1.754, mass=27.704),\n",
" PseudoJet(pt=46.337, eta=3.811, phi=1.327, mass=31.050),\n",
" PseudoJet(pt=32.825, eta=6.472, phi=-0.657, mass=27.099),\n",
" PseudoJet(pt=31.514, eta=-5.126, phi=-2.804, mass=20.859),\n",
" PseudoJet(pt=28.626, eta=-3.931, phi=2.712, mass=23.153),\n",
" PseudoJet(pt=26.634, eta=-1.484, phi=3.003, mass=21.662),\n",
" PseudoJet(pt=22.779, eta=0.132, phi=0.052, mass=13.168),\n",
" PseudoJet(pt=20.782, eta=2.555, phi=0.209, mass=9.057),\n",
" PseudoJet(pt=19.669, eta=0.195, phi=1.228, mass=10.364),\n",
" PseudoJet(pt=19.503, eta=-6.295, phi=-2.932, mass=7.497),\n",
" PseudoJet(pt=17.514, eta=-1.786, phi=-0.592, mass=9.562),\n",
" PseudoJet(pt=17.212, eta=5.665, phi=1.594, mass=8.971),\n",
" PseudoJet(pt=14.526, eta=-4.358, phi=1.472, mass=11.143),\n",
" PseudoJet(pt=12.627, eta=5.917, phi=-2.111, mass=6.795),\n",
" PseudoJet(pt=11.891, eta=-7.104, phi=2.663, mass=3.586),\n",
" PseudoJet(pt=11.349, eta=3.582, phi=3.019, mass=6.530),\n",
" PseudoJet(pt=10.842, eta=5.065, phi=0.362, mass=7.501),\n",
" PseudoJet(pt=8.518, eta=-6.508, phi=-0.466, mass=6.881),\n",
" PseudoJet(pt=8.408, eta=6.438, phi=3.051, mass=3.712),\n",
" PseudoJet(pt=7.669, eta=-4.011, phi=-1.179, mass=9.382),\n",
" PseudoJet(pt=6.176, eta=-4.783, phi=0.077, mass=4.470),\n",
" PseudoJet(pt=4.541, eta=-8.366, phi=-0.505, mass=3.769),\n",
" PseudoJet(pt=3.767, eta=-6.667, phi=1.190, mass=3.765),\n",
" PseudoJet(pt=3.486, eta=-7.686, phi=1.182, mass=1.600),\n",
" PseudoJet(pt=3.227, eta=4.536, phi=2.029, mass=1.138),\n",
" PseudoJet(pt=2.581, eta=4.794, phi=-1.750, mass=1.873),\n",
" PseudoJet(pt=2.355, eta=0.146, phi=-2.828, mass=0.169),\n",
" PseudoJet(pt=2.278, eta=6.119, phi=0.661, mass=0.656),\n",
" PseudoJet(pt=2.201, eta=7.765, phi=0.471, mass=1.150),\n",
" PseudoJet(pt=2.053, eta=8.018, phi=-2.185, mass=0.000),\n",
" PseudoJet(pt=1.842, eta=-5.225, phi=1.291, mass=0.456),\n",
" PseudoJet(pt=1.719, eta=-2.979, phi=0.135, mass=1.514),\n",
" PseudoJet(pt=1.582, eta=3.394, phi=-2.097, mass=1.298),\n",
" PseudoJet(pt=1.486, eta=7.143, phi=-2.384, mass=0.995),\n",
" PseudoJet(pt=1.088, eta=2.972, phi=2.104, mass=0.014),\n",
" PseudoJet(pt=1.081, eta=-2.347, phi=2.622, mass=0.293),\n",
" PseudoJet(pt=0.810, eta=-7.510, phi=-1.868, mass=1.272),\n",
" PseudoJet(pt=0.516, eta=7.326, phi=-0.953, mass=0.248),\n",
" PseudoJet(pt=0.406, eta=-1.269, phi=1.900, mass=0.000),\n",
" PseudoJet(pt=0.212, eta=-6.253, phi=-1.455, mass=0.938),\n",
" PseudoJet(pt=0.162, eta=2.565, phi=-3.072, mass=0.140),\n",
" PseudoJet(pt=0.144, eta=-4.792, phi=2.294, mass=0.000),\n",
" PseudoJet(pt=0.012, eta=2.395, phi=1.268, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=31.770, phi=0.335, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=34.615, phi=1.938, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=33.444, phi=-0.307, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-32.601, phi=0.945, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-32.702, phi=-1.425, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-35.221, phi=-2.805, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-35.800, phi=0.771, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=37.186, phi=-2.300, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-35.800, phi=3.087, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=2.834, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=37.218, phi=-1.719, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=34.522, phi=2.839, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=-0.278, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=38.171, phi=0.735, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=-2.316, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-35.800, phi=-0.043, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=36.630, phi=1.187, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=35.931, phi=-0.707, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=39.354, phi=2.547, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=-3.081, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=-0.241, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=34.849, phi=2.934, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=35.756, phi=0.059, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-35.800, phi=-1.951, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=35.219, phi=1.093, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-35.800, phi=-2.034, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=37.245, phi=2.578, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-32.480, phi=-0.895, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=2.099, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-32.854, phi=-0.136, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-34.857, phi=1.922, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=39.256, phi=1.373, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=2.046, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-35.800, phi=-0.232, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=39.846, phi=-0.936, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=40.118, phi=-0.934, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=37.767, phi=1.695, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=37.423, phi=-1.508, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-34.381, phi=-3.085, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=41.075, phi=-2.776, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=41.590, phi=0.078, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=38.829, phi=-2.356, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=42.100, phi=1.107, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=43.374, phi=3.142, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=42.311, phi=-2.266, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=-2.266, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=39.606, phi=0.588, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=-2.554, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=-37.332, phi=-1.571, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=45.314, phi=3.142, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.938),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.938),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=0.000),\n",
" PseudoJet(pt=0.000, eta=100000.000, phi=0.000, mass=-0.000)]"
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"jets"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"PseudoJet(pt=9514.937, eta=-0.143, phi=2.402, mass=1084.248)"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# leading jet\n",
"jets[0]"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"PseudoJet(pt=7942.452, eta=1.313, phi=-0.807, mass=2062.356)"
]
},
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# subleading jet\n",
"jets[1]"
]
},
{
"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.10"
},
"latex_envs": {
"bibliofile": "biblio.bib",
"cite_by": "apalike",
"current_citInitial": 1,
"eqLabelWithNumbers": true,
"eqNumInitial": 0
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment