diff --git a/README.md b/README.md index fa7ecdd7c07888a085abd779503082f977c8e379..f468da86f1b3f552231857bbd021083d3b8576f8 100644 --- a/README.md +++ b/README.md @@ -5,9 +5,9 @@ git clone https://eicweb.phy.anl.gov/jihee.kim/study.git cd study ``` -Input Root Files ----------------- - - from lAger event generator +Input File +---------- + - from `lAger` event generator Scripts ------- diff --git a/sample/input.root b/sample/input.root new file mode 100644 index 0000000000000000000000000000000000000000..1d939c5ccd4c9f2b40a1b355434f3de58de6421a Binary files /dev/null and b/sample/input.root differ diff --git a/sample/read_and_print.py b/sample/read_and_print.py new file mode 100644 index 0000000000000000000000000000000000000000..59151c0b3c7a982b96848ac79307531ca556e597 --- /dev/null +++ b/sample/read_and_print.py @@ -0,0 +1,56 @@ +############################# +# Read ROOT File +# Loop over event by event +# Print interested variables +############################# +import ROOT + +# ROOT file that you read +fname = "input.root" + +# Open file and get TTree +f = ROOT.TFile.Open(fname) +tree = f.Get("lAger") +print("=======") +print("Read in ROOT file: {}".format(fname)) + +# Total number of events in TTree +nevt = tree.GetEntries() +print ("total number of events: ",nevt) + +# Loop over event by event +print("Loop over event by event") +for ievt in range(0,nevt): + # Access info of an event + tree.GetEntry(ievt) + + if ievt % 5 == 0: # Printing + print("Processing event number: {}".format(ievt)) + + W = tree.W # invariant mass: sqrt(center of mass energy virtual photon + target (Helium-4)) + Q2 = tree.Q2 # virtual photon mass: four-momentum transfer of electron + t = tree.t # four-momentum transfer of Helium-4 + + # Total number of particles in an event + nptls = tree.particles.GetEntries() + print (nptls) + + ip = 0 # index of particles + # Loop over particles + for ptl in range(0,nptls): + # TODO: Access leaf and Save as in array + # infomration needed + # each particles have below info + # particles.fPdgCode / particles.fStatusCode / particles.fE / particles.fPx / particles.fPy / particles.fPz + + # Not Working... + #PID = tree.particles[ptl].fPdgCode + #PID[ip] = tree.GetLeaf("particles.fPdgCode").GetValue(ptl) + ip += 1 + #print (PID) + + # Just checking the first 10 events + if ievt > 9: + break + pass +print("=======")