The Truth Particles Collection

The Truth Particles Collection

The Truth particle collections are rich with information. Unfortunately, they are also a little tricky because what is actually available in an event depends on the sliming algorithms that are applied to the data. The code to access the data is always the same - but the fact that the data is there or not will depend on the file’s processing history. Protecting against that makes our code a bit more complex, as you’ll see below.

import matplotlib.pyplot as plt
from config import ds_jz2_exot15 as ds
import numpy as np
import pandas as pd
from particle import Particle
from hist import Hist

pd.options.display.max_rows = None

You’ll note above we are using ds_jz2_exot15. This is a Run 2 derivation, and not a DAOD_PHYS file. The reason is this has large TruthParticles and TruthVertices collections. DAOD_PHYS has TruthParticles collections, but only for select particles (TruthBottom, TruthTop, TruthBSM, for example). To demonstrate moving between particles and vertices we wanted something more complex to demonstrate.

We’ll start with particles.

all_particles = (ds
                 .SelectMany(lambda e: e.TruthParticles("TruthParticles"))
                 .Where(lambda tp: (tp.pt() / 1000) > 10)
                 .Select(lambda tp: {
                         'pt': tp.pt() / 1000.0,
                         'pdgId': abs(tp.pdgId()),
                 })
                 .AsAwkwardArray()
                 .value())
plt.hist(all_particles.pt, bins=100, range=(0, 30))
plt.xlabel('Truth Particle $p_T$ [GeV]')
plt.ylabel('Number of Particles')
_ = plt.title(r'Truth Particle $p_T$ distribution for $Z\rightarrow ee$ events')
../_images/truth_5_0.png

Lets take a very quick look at the pdgid’s that we have in the data we are looking at. Note that we already applied, above, abs to the pdgID to make this list more compact. We’ll use the amazing particle pypi package to get particle names. One thing is there are a few particles that aren’t in the database, so we need to protect that lookup during the lookup. Note that it is possible to create ones own named particles and add them - so this could be done in one’s own analysis if desired.

def get_pdg_name(pdgid: int) -> str:
    'Protected particle name from pdgid'
    if pdgid == 35:
        return 'LLP'
    try:
        return Particle.from_pdgid(pdgid).name
    except:
        # Protect against unknown particle id's
        return f'Unknown: {pdgid}'

# Use a pandas DataFrame to make displaying a table of particles easy.
unique, counts = np.unique(all_particles.pdgId, return_counts=True)
df = pd.DataFrame({'PDGId' : unique, 'Count': counts, 'Name': [get_pdg_name(int(i)) for i in unique] })
df['PDGId'] = df['PDGId'].astype(int)
df[df.Count>100].sort_values('Name')
PDGId Count Name
50 533 8586 B(s)*0
49 531 18482 B(s)0
48 523 38349 B*+
46 513 37963 B*0
47 521 55403 B+
45 511 65420 B0
132 10411 263 D(0)*(2300)+
134 10421 370 D(0)*(2300)0
133 10413 730 D(1)(2420)+
135 10423 736 D(1)(2420)0
147 20423 410 D(1)(2430)0
35 415 708 D(2)*(2460)+
38 425 831 D(2)*(2460)0
40 433 17953 D(s)*+
39 431 33734 D(s)+
136 10431 236 D(s0)*(2317)+
148 20433 1260 D(s1)(2460)+
137 10433 252 D(s1)(2536)+
41 435 125 D(s2)*(2573)+
37 423 63360 D*(2007)0
34 413 64153 D*(2010)+
33 411 79878 D+
36 421 162787 D0
64 2214 70632 Delta(1232)+
65 2224 69776 Delta(1232)++
57 1114 67883 Delta(1232)-
61 2114 69400 Delta(1232)0
43 443 1456 J/psi(1S)
127 10311 797 K(0)*(1430)0
130 10323 1802 K(1)(1270)+
128 10313 1413 K(1)(1270)0
144 20313 3401 K(1)(1400)0
29 325 390 K(2)*(1430)+
26 315 489 K(2)*(1430)0
16 130 251501 K(L)0
23 310 256050 K(S)0
28 323 266493 K*(892)+
25 313 268251 K*(892)0
27 321 507919 K+
24 311 486569 K0
70 3122 161330 Lambda
109 5122 4817 Lambda(b)0
83 4122 13005 Lambda(c)+
97 4332 191 Omega(c)0
79 3334 674 Omega-
74 3224 15478 Sigma(1385)+
69 3114 14849 Sigma(1385)-
72 3214 15284 Sigma(1385)0
115 5224 407 Sigma(b)*+
108 5114 449 Sigma(b)*-
114 5222 207 Sigma(b)+
107 5112 229 Sigma(b)-
88 4212 749 Sigma(c)(2455)+
90 4222 694 Sigma(c)(2455)++
81 4112 790 Sigma(c)(2455)0
89 4214 949 Sigma(c)(2520)+
91 4224 841 Sigma(c)(2520)++
82 4114 1098 Sigma(c)(2520)0
73 3222 55237 Sigma+
68 3112 53338 Sigma-
71 3212 54348 Sigma0
146 20413 738 Unknown: 20413
112 5212 215 Unknown: 5212
113 5214 417 Unknown: 5214
76 3314 2966 Xi(1530)-
78 3324 3170 Xi(1530)0
110 5132 690 Xi(b)-
116 5232 675 Xi(b)0
95 4322 306 Xi(c)'+
93 4312 356 Xi(c)'0
96 4324 298 Xi(c)(2645)+
94 4314 317 Xi(c)(2645)0
92 4232 1761 Xi(c)+
85 4132 1818 Xi(c)0
75 3312 16656 Xi-
77 3322 17065 Xi0
142 20213 14951 a(1)(1260)+
4 5 1165840 b
3 4 1973433 c
138 10441 178 chi(c0)(1P)
149 20443 305 chi(c1)(1P)
44 445 140 chi(c2)(1P)
0 1 4295800 d
5 11 52546 e-
20 221 281030 eta
30 331 62947 eta'(958)
42 441 658 eta(c)(1S)
155 100441 140 eta(c)(2S)
162 9010221 297 f(0)(980)
11 21 52939440 g
12 22 830967 gamma
7 13 9119 mu-
60 2112 365179 n
6 12 9138 nu(e)
8 14 8453 nu(mu)
10 16 2043 nu(tau)
21 223 438327 omega(782)
63 2212 366263 p
31 333 40821 phi(1020)
17 211 1992975 pi+
13 111 1040796 pi0
156 100443 490 psi(2S)
153 30443 194 psi(3770)
18 213 920756 rho(770)+
14 113 475627 rho(770)0
2 3 2635153 s
9 15 4183 tau-
1 2 5474766 u

Truth Vertices

Lets explore the decay chain a bit. Lets look at what \(b\)-quarks hadrons decay into in this sample, and how far a way from their production the hadrons decay. The ATLAS data model has several ways to follow the decay chain. Because we want to understand decay distances, we will need to use the TruthVertex object.

The ATLAS truth data model is, in the world of infinite disk space, quite straight forward. Every TruthParticle has a production and a decay TruthVertex associated with it. As long as the particle isn’t stable (like an electron). However, in order to keep disk space under control in many cases the children or parent vertices may be missing. This leads us to a number of caveats that have to be encoded in our query:

  • In ATLAS, particles don’t have to have production or decay vertices. If this is the case, then the link to the production or decay vertex will be a null pointer (at least, in C++). And if you generate code that references it, you’ll get a hard segfault in ServiceX. So, always check to see if the vertex you are interested in exists!

  • \(b\)-quarks can radiate, and thus decay into another \(b\) quark and a \(\gamma\) or similar. We are not interested in those decays.

  • The \(b\)-quark hadronizes almost immediately, of course. So we are really interested in what the thing the \(b\)-quark hadronizes into is.

  • The hadronization isn’t straight forward, due to color strings, etc.

real_children = (
          ds
          .SelectMany(lambda e: e.TruthParticles("TruthParticles"))
          .Where(lambda tp: (abs(tp.pdgId()) == 5))
          .Where(lambda tp: tp.hasDecayVtx())
          .Where(lambda tp: tp.decayVtx().outgoingParticleLinks().Where(lambda dtp: dtp.isValid()).Where(lambda dtp: abs(dtp.pdgId()) == 5).Count() == 0)
          .SelectMany(lambda tp: tp.decayVtx().outgoingParticleLinks())
          .Where(lambda tp: tp.isValid())
          .Where(lambda tp: tp.hasBottom())
          .Where(lambda tp: tp.hasProdVtx() and tp.hasDecayVtx())
)

Note the protections in this code:

  • Use of hasDecayVtx and hasProdVtx - which tell us if the thing has any vertices at all.

  • Use of isValid, which is used against an ElementLink<X> object. If a particle has been removed during derivation, then ElementLink will point to nothing. isValid protects us from trying to follow that dead end link.

children = (real_children
    .Select(lambda tp: {
        'pdgId': abs(tp.pdgId()),
        'dx': abs(tp.prodVtx().x() - tp.decayVtx().x()),
        'dy': abs(tp.prodVtx().y() - tp.decayVtx().y()),
        'dz': abs(tp.prodVtx().z() - tp.decayVtx().z()),
    })
    .AsAwkwardArray()
    .value()
)

Dump out the names of all the hadrons. Partly, this is to make sure we got it right!

unique, counts = np.unique(children.pdgId, return_counts=True)
df = pd.DataFrame({'PDGId' : unique, 'Count': counts, 'Name': [get_pdg_name(int(i)) for i in unique] })
df['PDGId'] = df['PDGId'].astype(int)
df.sort_values('Count', ascending=False)
PDGId Count Name
3 523 65075 B*+
1 513 64754 B*0
2 521 29526 B+
0 511 29369 B0
5 533 14535 B(s)*0
4 531 6426 B(s)0
12 5122 5083 Lambda(b)0
18 5232 945 Xi(b)0
13 5132 941 Xi(b)-
11 5114 763 Sigma(b)*-
17 5224 722 Sigma(b)*+
15 5214 712 Unknown: 5214
10 5112 421 Sigma(b)-
14 5212 406 Unknown: 5212
16 5222 389 Sigma(b)+
20 5314 148 Unknown: 5314
22 5324 113 Unknown: 5324
7 543 84 Unknown: 543
19 5312 77 Unknown: 5312
21 5322 56 Unknown: 5322
6 541 38 B(c)+
24 5334 25 Unknown: 5334
23 5332 14 Omega(b)-
9 553 10 Upsilon(1S)
8 551 4 Unknown: 551

And, now we can calculate the decay distance:

decay_length_all = np.sqrt(children.dx**2 + children.dy**2 + children.dz**2)
pdg_name_all = np.array([get_pdg_name(int(i)) for i in children.pdgId])
good_decay = decay_length_all > 0.0
h = (
  Hist.new
  .Reg(50, 0 ,20, name="dl", label="Decay Length [mm]")
  .StrCategory([], name="pdg", label="Particle", growth=True)
  .Int64()
)
_ = h.fill(dl=decay_length_all[good_decay], pdg=pdg_name_all[good_decay])
plt.figure(figsize=(15, 7))
h.plot()
plt.yscale('log')
_ = plt.legend()
../_images/truth_17_0.png

I love this plot because it shows the different decay probabilities. There are differences in lifetime, but they are harder to see because we aren’t showing \(c\tau\) here (and correcting for boost, mass, etc.).

The Datamodel

The data model when this documentation was last built was:

from func_adl_servicex_xaodr21.xAOD.truthevent_v1 import TruthEvent_v1
help(TruthEvent_v1)
Help on class TruthEvent_v1 in module func_adl_servicex_xaodr21.xAOD.truthevent_v1:

class TruthEvent_v1(builtins.object)
 |  A class
 |  
 |  Methods defined here:
 |  
 |  beamParticle1Link(self) -> 'func_adl_servicex_xaodr21.elementlink_datavector_xaod_truthparticle_v1__.ElementLink_DataVector_xAOD_TruthParticle_v1__'
 |      A method
 |  
 |  beamParticle2Link(self) -> 'func_adl_servicex_xaodr21.elementlink_datavector_xaod_truthparticle_v1__.ElementLink_DataVector_xAOD_TruthParticle_v1__'
 |      A method
 |  
 |  clearDecorations(self) -> 'bool'
 |      A method
 |  
 |  crossSection(self) -> 'float'
 |      A method
 |  
 |  crossSectionError(self) -> 'float'
 |      A method
 |  
 |  hasNonConstStore(self) -> 'bool'
 |      A method
 |  
 |  hasStore(self) -> 'bool'
 |      A method
 |  
 |  index(self) -> 'int'
 |      A method
 |  
 |  nTruthParticles(self) -> 'int'
 |      A method
 |  
 |  nTruthVertices(self) -> 'int'
 |      A method
 |  
 |  pdfInfo(self) -> 'func_adl_servicex_xaodr21.xAOD.TruthEvent_v1.pdfinfo.PdfInfo'
 |      A method
 |  
 |  signalProcessVertex(self) -> 'func_adl_servicex_xaodr21.xAOD.truthvertex_v1.TruthVertex_v1'
 |      A method
 |  
 |  signalProcessVertexLink(self) -> 'func_adl_servicex_xaodr21.elementlink_datavector_xaod_truthvertex_v1__.ElementLink_DataVector_xAOD_TruthVertex_v1__'
 |      A method
 |  
 |  truthParticle(self, index: 'int') -> 'func_adl_servicex_xaodr21.xAOD.truthparticle_v1.TruthParticle_v1'
 |      A method
 |  
 |  truthParticleLink(self, index: 'int') -> 'func_adl_servicex_xaodr21.elementlink_datavector_xaod_truthparticle_v1__.ElementLink_DataVector_xAOD_TruthParticle_v1__'
 |      A method
 |  
 |  truthParticleLinks(self) -> 'func_adl_servicex_xaodr21.vector_elementlink_datavector_xaod_truthparticle_v1___.vector_ElementLink_DataVector_xAOD_TruthParticle_v1___'
 |      A method
 |  
 |  truthVertex(self, index: 'int') -> 'func_adl_servicex_xaodr21.xAOD.truthvertex_v1.TruthVertex_v1'
 |      A method
 |  
 |  truthVertexLink(self, index: 'int') -> 'func_adl_servicex_xaodr21.elementlink_datavector_xaod_truthvertex_v1__.ElementLink_DataVector_xAOD_TruthVertex_v1__'
 |      A method
 |  
 |  truthVertexLinks(self) -> 'func_adl_servicex_xaodr21.vector_elementlink_datavector_xaod_truthvertex_v1___.vector_ElementLink_DataVector_xAOD_TruthVertex_v1___'
 |      A method
 |  
 |  usingPrivateStore(self) -> 'bool'
 |      A method
 |  
 |  usingStandaloneStore(self) -> 'bool'
 |      A method
 |  
 |  weights(self) -> 'func_adl_servicex_xaodr21.vector_float_.vector_float_'
 |      A method
 |  
 |  ----------------------------------------------------------------------
 |  Readonly properties defined here:
 |  
 |  auxdataConst
 |      A method
 |  
 |  isAvailable
 |      A method
 |  
 |  ----------------------------------------------------------------------
 |  Data descriptors defined here:
 |  
 |  __dict__
 |      dictionary for instance variables (if defined)
 |  
 |  __weakref__
 |      list of weak references to the object (if defined)
from func_adl_servicex_xaodr21.xAOD.truthparticle_v1 import TruthParticle_v1
help(TruthParticle_v1)
Help on class TruthParticle_v1 in module func_adl_servicex_xaodr21.xAOD.truthparticle_v1:

class TruthParticle_v1(builtins.object)
 |  A class
 |  
 |  Methods defined here:
 |  
 |  absPdgId(self) -> 'int'
 |      A method
 |  
 |  abseta(self) -> 'float'
 |      A method
 |  
 |  absrapidity(self) -> 'float'
 |      A method
 |  
 |  barcode(self) -> 'int'
 |      A method
 |  
 |  charge(self) -> 'float'
 |      A method
 |  
 |  child(self, i: 'int') -> 'func_adl_servicex_xaodr21.xAOD.truthparticle_v1.TruthParticle_v1'
 |      A method
 |  
 |  clearDecorations(self) -> 'bool'
 |      A method
 |  
 |  decayVtx(self) -> 'func_adl_servicex_xaodr21.xAOD.truthvertex_v1.TruthVertex_v1'
 |      A method
 |  
 |  decayVtxLink(self) -> 'func_adl_servicex_xaodr21.elementlink_datavector_xaod_truthvertex_v1__.ElementLink_DataVector_xAOD_TruthVertex_v1__'
 |      A method
 |  
 |  e(self) -> 'float'
 |      A method
 |  
 |  eta(self) -> 'float'
 |      A method
 |  
 |  hasBottom(self) -> 'bool'
 |      A method
 |  
 |  hasCharm(self) -> 'bool'
 |      A method
 |  
 |  hasDecayVtx(self) -> 'bool'
 |      A method
 |  
 |  hasNonConstStore(self) -> 'bool'
 |      A method
 |  
 |  hasProdVtx(self) -> 'bool'
 |      A method
 |  
 |  hasStore(self) -> 'bool'
 |      A method
 |  
 |  hasStrange(self) -> 'bool'
 |      A method
 |  
 |  index(self) -> 'int'
 |      A method
 |  
 |  isBSM(self) -> 'bool'
 |      A method
 |  
 |  isBaryon(self) -> 'bool'
 |      A method
 |  
 |  isBottomBaryon(self) -> 'bool'
 |      A method
 |  
 |  isBottomHadron(self) -> 'bool'
 |      A method
 |  
 |  isBottomMeson(self) -> 'bool'
 |      A method
 |  
 |  isChLepton(self) -> 'bool'
 |      A method
 |  
 |  isCharged(self) -> 'bool'
 |      A method
 |  
 |  isCharmBaryon(self) -> 'bool'
 |      A method
 |  
 |  isCharmHadron(self) -> 'bool'
 |      A method
 |  
 |  isCharmMeson(self) -> 'bool'
 |      A method
 |  
 |  isElectron(self) -> 'bool'
 |      A method
 |  
 |  isGenSpecific(self) -> 'bool'
 |      A method
 |  
 |  isHadron(self) -> 'bool'
 |      A method
 |  
 |  isHeavyBaryon(self) -> 'bool'
 |      A method
 |  
 |  isHeavyHadron(self) -> 'bool'
 |      A method
 |  
 |  isHeavyMeson(self) -> 'bool'
 |      A method
 |  
 |  isHiggs(self) -> 'bool'
 |      A method
 |  
 |  isLepton(self) -> 'bool'
 |      A method
 |  
 |  isLightBaryon(self) -> 'bool'
 |      A method
 |  
 |  isLightHadron(self) -> 'bool'
 |      A method
 |  
 |  isLightMeson(self) -> 'bool'
 |      A method
 |  
 |  isMeson(self) -> 'bool'
 |      A method
 |  
 |  isMuon(self) -> 'bool'
 |      A method
 |  
 |  isNeutral(self) -> 'bool'
 |      A method
 |  
 |  isNeutrino(self) -> 'bool'
 |      A method
 |  
 |  isParton(self) -> 'bool'
 |      A method
 |  
 |  isPhoton(self) -> 'bool'
 |      A method
 |  
 |  isQuark(self) -> 'bool'
 |      A method
 |  
 |  isResonance(self) -> 'bool'
 |      A method
 |  
 |  isStrangeBaryon(self) -> 'bool'
 |      A method
 |  
 |  isStrangeHadron(self) -> 'bool'
 |      A method
 |  
 |  isStrangeMeson(self) -> 'bool'
 |      A method
 |  
 |  isTau(self) -> 'bool'
 |      A method
 |  
 |  isTop(self) -> 'bool'
 |      A method
 |  
 |  isW(self) -> 'bool'
 |      A method
 |  
 |  isZ(self) -> 'bool'
 |      A method
 |  
 |  m(self) -> 'float'
 |      A method
 |  
 |  nChildren(self) -> 'int'
 |      A method
 |  
 |  nParents(self) -> 'int'
 |      A method
 |  
 |  p4(self) -> 'func_adl_servicex_xaodr21.tlorentzvector.TLorentzVector'
 |      A method
 |  
 |  parent(self, i: 'int') -> 'func_adl_servicex_xaodr21.xAOD.truthparticle_v1.TruthParticle_v1'
 |      A method
 |  
 |  pdgId(self) -> 'int'
 |      A method
 |  
 |  phi(self) -> 'float'
 |      A method
 |  
 |  polarization(self) -> 'func_adl_servicex_xaodr21.xAOD.TruthParticle_v1.polarization.Polarization'
 |      A method
 |  
 |  prodVtx(self) -> 'func_adl_servicex_xaodr21.xAOD.truthvertex_v1.TruthVertex_v1'
 |      A method
 |  
 |  prodVtxLink(self) -> 'func_adl_servicex_xaodr21.elementlink_datavector_xaod_truthvertex_v1__.ElementLink_DataVector_xAOD_TruthVertex_v1__'
 |      A method
 |  
 |  pt(self) -> 'float'
 |      A method
 |  
 |  px(self) -> 'float'
 |      A method
 |  
 |  py(self) -> 'float'
 |      A method
 |  
 |  pz(self) -> 'float'
 |      A method
 |  
 |  rapidity(self) -> 'float'
 |      A method
 |  
 |  status(self) -> 'int'
 |      A method
 |  
 |  threeCharge(self) -> 'int'
 |      A method
 |  
 |  usingPrivateStore(self) -> 'bool'
 |      A method
 |  
 |  usingStandaloneStore(self) -> 'bool'
 |      A method
 |  
 |  ----------------------------------------------------------------------
 |  Readonly properties defined here:
 |  
 |  auxdataConst
 |      A method
 |  
 |  isAvailable
 |      A method
 |  
 |  ----------------------------------------------------------------------
 |  Data descriptors defined here:
 |  
 |  __dict__
 |      dictionary for instance variables (if defined)
 |  
 |  __weakref__
 |      list of weak references to the object (if defined)
from func_adl_servicex_xaodr21.xAOD.truthvertex_v1 import TruthVertex_v1
help(TruthVertex_v1)
Help on class TruthVertex_v1 in module func_adl_servicex_xaodr21.xAOD.truthvertex_v1:

class TruthVertex_v1(builtins.object)
 |  A class
 |  
 |  Methods defined here:
 |  
 |  barcode(self) -> 'int'
 |      A method
 |  
 |  clearDecorations(self) -> 'bool'
 |      A method
 |  
 |  eta(self) -> 'float'
 |      A method
 |  
 |  hasNonConstStore(self) -> 'bool'
 |      A method
 |  
 |  hasStore(self) -> 'bool'
 |      A method
 |  
 |  id(self) -> 'int'
 |      A method
 |  
 |  incomingParticle(self, index: 'int') -> 'func_adl_servicex_xaodr21.xAOD.truthparticle_v1.TruthParticle_v1'
 |      A method
 |  
 |  incomingParticleLinks(self) -> 'func_adl_servicex_xaodr21.vector_elementlink_datavector_xaod_truthparticle_v1___.vector_ElementLink_DataVector_xAOD_TruthParticle_v1___'
 |      A method
 |  
 |  index(self) -> 'int'
 |      A method
 |  
 |  nIncomingParticles(self) -> 'int'
 |      A method
 |  
 |  nOutgoingParticles(self) -> 'int'
 |      A method
 |  
 |  outgoingParticle(self, index: 'int') -> 'func_adl_servicex_xaodr21.xAOD.truthparticle_v1.TruthParticle_v1'
 |      A method
 |  
 |  outgoingParticleLinks(self) -> 'func_adl_servicex_xaodr21.vector_elementlink_datavector_xaod_truthparticle_v1___.vector_ElementLink_DataVector_xAOD_TruthParticle_v1___'
 |      A method
 |  
 |  perp(self) -> 'float'
 |      A method
 |  
 |  phi(self) -> 'float'
 |      A method
 |  
 |  t(self) -> 'float'
 |      A method
 |  
 |  usingPrivateStore(self) -> 'bool'
 |      A method
 |  
 |  usingStandaloneStore(self) -> 'bool'
 |      A method
 |  
 |  v4(self) -> 'func_adl_servicex_xaodr21.tlorentzvector.TLorentzVector'
 |      A method
 |  
 |  x(self) -> 'float'
 |      A method
 |  
 |  y(self) -> 'float'
 |      A method
 |  
 |  z(self) -> 'float'
 |      A method
 |  
 |  ----------------------------------------------------------------------
 |  Readonly properties defined here:
 |  
 |  auxdataConst
 |      A method
 |  
 |  isAvailable
 |      A method
 |  
 |  ----------------------------------------------------------------------
 |  Data descriptors defined here:
 |  
 |  __dict__
 |      dictionary for instance variables (if defined)
 |  
 |  __weakref__
 |      list of weak references to the object (if defined)

Further Information