# 10: Scikit-HEP

The Scikit-HEP project is a well maintained collection of HEP packages in Python that are useful for analysis. It contains tools ranging from plotting helpers, PDG lookup, DecayLanguage converters over high-performance histogramming libraries and ROOT I/O up to likelihood fitting and statistics libraries.

This is a minimal overview over the packages that are available and have not yet been used in the other tutorials.

## formulate - converting expressions

The package formulate converts between different style of expressions, notably between ROOT and numexpr. The former was used in the previous notebooks to apply cuts and create new variables with Pandas DataFrames query and eval methods.

The standard usage is with a from_style and to_style, but it allows for more as it parses the string such as the variables used can be extracted.

:

import formulate

momentum = formulate.from_root('TMath::Sqrt(X_PX**2 + X_PY**2 + X_PZ**2)')
momentum

:

Expression<SQRT>(Expression<ADD>(Expression<POW>(Variable(X_PX), UnnamedConstant(2)), Expression<POW>(Variable(X_PY), UnnamedConstant(2)), Expression<POW>(Variable(X_PZ), UnnamedConstant(2))))

:

momentum.to_numexpr()  # as used in Pandas eval/query

:

'sqrt(((X_PX ** 2) + (X_PY ** 2) + (X_PZ ** 2)))'

:

momentum.to_root()  # as used in ROOT

:

'TMath::Sqrt(((X_PX ** 2) + (X_PY ** 2) + (X_PZ ** 2)))'


## Particle

particle provides an easy access to particle properties such as mass and decay width from PDG. One example is to get the actual particle in a Particle object which holds all the relevant information.

:

# Particle

from particle import Particle

piplus = Particle.from_pdgid(211)


Here we can access all the properties we wish to such as mass and decay width.

:

piplus.mass

:

139.57039

:

piplus.charge

:

1.0

:

piplus.width

:

2.5284e-14

:

piplus.name

:

'pi+'


It also serves as a data-base like structure and we can do searches. Let’s find all neutral beauty hadrons.

:

Particle.findall(lambda p: p.pdgid.has_bottom and p.charge==0)

:

[<Particle: name="B0", pdgid=511, mass=5279.66 ± 0.12 MeV>,
<Particle: name="B~0", pdgid=-511, mass=5279.66 ± 0.12 MeV>,
<Particle: name="B*0", pdgid=513, mass=5324.71 ± 0.21 MeV>,
<Particle: name="B*~0", pdgid=-513, mass=5324.71 ± 0.21 MeV>,
<Particle: name="B(2)*(5747)0", pdgid=515, mass=5739.5 ± 0.7 MeV>,
<Particle: name="B(2)*(5747)~0", pdgid=-515, mass=5739.5 ± 0.7 MeV>,
<Particle: name="B(s)0", pdgid=531, mass=5366.92 ± 0.10 MeV>,
<Particle: name="B(s)~0", pdgid=-531, mass=5366.92 ± 0.10 MeV>,
<Particle: name="B(s)*0", pdgid=533, mass=5415.4 + 1.8 - 1.5 MeV>,
<Particle: name="B(s)*~0", pdgid=-533, mass=5415.4 + 1.8 - 1.5 MeV>,
<Particle: name="B(s2)*(5840)0", pdgid=535, mass=5839.86 ± 0.12 MeV>,
<Particle: name="B(s2)*(5840)~0", pdgid=-535, mass=5839.86 ± 0.12 MeV>,
<Particle: name="eta(b)(1S)", pdgid=551, mass=9398.7 ± 2.0 MeV>,
<Particle: name="Upsilon(1S)", pdgid=553, mass=9460.40 ± 0.10 MeV>,
<Particle: name="chi(b2)(1P)", pdgid=555, mass=9912.2 ± 0.4 MeV>,
<Particle: name="Lambda(b)0", pdgid=5122, mass=5619.60 ± 0.17 MeV>,
<Particle: name="Lambda(b)~0", pdgid=-5122, mass=5619.60 ± 0.17 MeV>,
<Particle: name="Xi(b)0", pdgid=5232, mass=5791.9 ± 0.5 MeV>,
<Particle: name="Xi(b)~0", pdgid=-5232, mass=5791.9 ± 0.5 MeV>,
<Particle: name="chi(b0)(1P)", pdgid=10551, mass=9859.4 ± 0.5 MeV>,
<Particle: name="h(b)(1P)", pdgid=10553, mass=9899.3 ± 0.8 MeV>,
<Particle: name="chi(b1)(1P)", pdgid=20553, mass=9892.8 ± 0.4 MeV>,
<Particle: name="Upsilon(2)(1D)", pdgid=20555, mass=10163.7 ± 1.4 MeV>,
<Particle: name="Upsilon(2S)", pdgid=100553, mass=10023.4 ± 0.5 MeV>,
<Particle: name="chi(b2)(2P)", pdgid=100555, mass=10268.6 ± 0.5 MeV>,
<Particle: name="chi(b0)(2P)", pdgid=110551, mass=10232.5 ± 0.6 MeV>,
<Particle: name="h(b)(2P)", pdgid=110553, mass=10259.8 ± 1.2 MeV>,
<Particle: name="chi(b1)(2P)", pdgid=120553, mass=10255.5 ± 0.5 MeV>,
<Particle: name="Upsilon(3S)", pdgid=200553, mass=10355.1 ± 0.5 MeV>,
<Particle: name="chi(b2)(3P)", pdgid=200555, mass=10524.0 ± 0.8 MeV>,
<Particle: name="chi(b1)(3P)", pdgid=220553, mass=10513.4 ± 0.7 MeV>,
<Particle: name="Upsilon(4S)", pdgid=300553, mass=10579.4 ± 1.2 MeV>,
<Particle: name="Upsilon(10860)", pdgid=9000553, mass=10885.2 + 2.6 - 1.6 MeV>,
<Particle: name="Upsilon(11020)", pdgid=9010553, mass=11000 ± 4 MeV>]


## hepunits

Ever hardcoded a physics constant? Or needed to manually convert units? That’s what the neat hepunits package is for.

It contains constants as well as units that all get converted into the standard units (e.g. MeV) which allows for calculations without the need to explicitly convert by just adding the unit; furthermore it greatly increases readability

:

from hepunits.constants import c_light

c_light

:

299.792458

:

import hepunits as u  # u for units

:

150 * u.MeV + 1.1 * u.GeV  # result in MeV

:

1250.0


## Vector

Manipulating and calculation vector quantities can be done using the vector library. It provides functionality to load, combine and work with vectors.

:

import vector

:

# NumPy-like arguments (literally passed through to NumPy)
vector.array([
(1.1, 2.1), (1.2, 2.2), (1.3, 2.3), (1.4, 2.4), (1.5, 2.5)
], dtype=[("x", float), ("y", float)])

# Pandas-like arguments (dict from names to column arrays)
vector.array({"x": [1.1, 1.2, 1.3, 1.4, 1.5], "y": [2.1, 2.2, 2.3, 2.4, 2.5]})

# As with objects, the coordinate system and dimension is taken from the names of the fields.
vec1 = vector.array({
"x": [1.1, 1.2, 1.3, 1.4, 1.5],
"y": [2.1, 2.2, 2.3, 2.4, 2.5],
"z": [3.1, 3.2, 3.3, 3.4, 3.5],
"t": [4.1, 4.2, 4.3, 4.4, 4.5],
})

vector.array({
"pt": [1.1, 1.2, 1.3, 1.4, 1.5],
"phi": [2.1, 2.2, 2.3, 2.4, 2.5],
"eta": [3.1, 3.2, 3.3, 3.4, 3.5],
"M": [4.1, 4.2, 4.3, 4.4, 4.5],
})

:

MomentumNumpy4D([(1.1, 2.1, 3.1, 4.1), (1.2, 2.2, 3.2, 4.2), (1.3, 2.3, 3.3, 4.3),
(1.4, 2.4, 3.4, 4.4), (1.5, 2.5, 3.5, 4.5)],
dtype=[('rho', '<f8'), ('phi', '<f8'), ('eta', '<f8'), ('tau', '<f8')])


### Vector properties

Vectors can have properties and provide quantities as attributes

:

vector.obj(x=1, y=2, z=3).theta
vector.obj(x=1, y=2, z=3).eta

:

1.1035868415601453

:

vector.obj(x=3, y=4, z=-2, t=10)   # Cartesian 4D vector

:

VectorObject4D(x=3, y=4, z=-2, t=10)

[ ]: