-
Notifications
You must be signed in to change notification settings - Fork 0
PythonCodeLoop
This method works really well if you want to explore the content of the ROOT file on an event by event basis, and then quickly make some histograms of the variables you discovered. It is also possible to write a full data analysis with the loop and this also works very well, but here the RDataFrame method will be a strong competitor, since the RDataFrames are likely to lead to more robust and faster code.
For more general details on how to use ROOT with Python, see the documentation for PyROOT:
Also, see Python notebook in Examples
In this example we open an HPS Monte Carlo output file, examine what is in it, explore the data in some of the events, and then make a histogram of some quantities in the file using a loop. It is a bit involved, but hopefully the comment explain sufficiently what is going on. I am committing the iPython prompts and code output responses so you can copy paste easily.
import ROOT as R
# Load the MiniDst library.
R.gSystem.Load("/data/HPS/lib/libMiniDst") # Note: the path to libMiniDst.so (or .dylib) will vary depending on where you installed it.
# Open a ROOT data file and get the TTree (you could also do a TChain with many files.)
root_file = R.TFile("pi0_455GeV_01_recon.root")
t = root_file.Get("MiniDST")
n_events = t.GetEntries() # Get the number of events.
# Initialize the MiniDst to be linked to the file.
mdst = R.MiniDst()
mdst.use_mc_particles=True # Turn *on* reading the MC Particles. This is not one of the branches that is always set to true.
mdst.DefineBranchMap() # This defines which branches will be read.
mdst.SetBranchAddressesOnTree(t) # Now link the TTree with the mdst object, so each branch is set on the mdst.
#
# Define example functions for printing the MC Particle tree
#
def print_daughters(i_part,indent):
print(" "*indent+f" {i_part:3d} pdg: {mdst.mc_part_pdg_id[i_part]:4d} E: {mdst.mc_part_energy[i_part]:9.6f} ")
if len(mdst.mc_part_daughters[i_part]) > 0:
print(" "*(indent+14) + "| ")
for i in range(len(mdst.mc_part_daughters[i_part])):
ii = mdst.mc_part_daughters[i_part][i] # Get the daughter reference
print_daughters(ii,indent+11) # Print by recursing
def print_mc_particle_tree():
for i in range(len(mdst.mc_part_parents)):
if len(mdst.mc_part_parents[i]) == 0: # top level particle
print_daughters(i,0)
#
# Now use the functions to print an event
#
t.GetEntry(1) # Get the 1 entry from the tree.
print_mc_particle_tree()
#
The output of the above code is shown here:
0 pdg: 2212 E: 1.084855
10 pdg: 111 E: 0.247149
|
7 pdg: 22 E: 0.109102
14 pdg: 22 E: 0.138047
|
5 pdg: -11 E: 0.029779
|
22 pdg: 22 E: 0.018695
8 pdg: 11 E: 0.108267
|
2 pdg: 22 E: 0.020494
4 pdg: 22 E: 0.014840
11 pdg: 22 E: 0.004657
13 pdg: 22 E: 0.016435
18 pdg: 22 E: 0.019782
12 pdg: 11 E: 4.155678
|
9 pdg: 22 E: 0.000000
|
1 pdg: -11 E: 0.000511
|
21 pdg: 22 E: 0.001406
|
15 pdg: 11 E: 0.001660
17 pdg: 22 E: 0.000000
|
3 pdg: 11 E: 0.000511
|
19 pdg: 22 E: 0.000000
|
20 pdg: 11 E: 0.000511
|
6 pdg: 22 E: 0.002974
|
16 pdg: 11 E: 0.002130