Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
tree130_jets.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_tree
3///
4/// Usage of a Tree using the JetEvent class.
5///
6/// The JetEvent class has several collections (TClonesArray)
7/// and other collections (TRefArray) referencing objects
8/// in the TClonesArrays.
9/// The JetEvent class is in $ROOTSYS/tutorials/io/tree/JetEvent.h,cxx
10/// to execute the script, do
11/// ~~~
12/// .x tree130_jets.C
13/// ~~~
14///
15/// \macro_code
16///
17/// \author Rene Brun
18
19#ifdef JETS_SECOND_RUN
20
21#include "TFile.h"
22#include "TTree.h"
23#include "TRandom.h"
24#include "TROOT.h"
25#include "TSystem.h"
26#include "JetEvent.h"
27#include "Riostream.h"
28
29
30void write(Int_t nev=100)
31{
32 //write nev Jet events
33 TFile f("JetEvent.root","recreate");
34 auto T = new TTree("T", "Event example with Jets");
35 auto event = new JetEvent;
36 T->Branch("event", "JetEvent", &event, 8000, 2);
37
38 for (Int_t ev=0; ev<nev; ev++) {
39 event->Build();
40 T->Fill();
41 }
42
43 T->Print();
44 T->Write();
45}
46
47void read()
48{
49 //read the JetEvent file
50 TFile f("JetEvent.root");
51 auto T = f.Get<TTree>("T");
52 JetEvent *event = 0;
53 T->SetBranchAddress("event", &event);
54 Long64_t nentries = T->GetEntries();
55
56 for (Long64_t ev=0;ev<nentries;ev++) {
57 T->GetEntry(ev);
58 if (ev)
59 continue; //dump first event only
60 std::cout << " Event: "<< ev
61 << " Jets: " << event->GetNjet()
62 << " Tracks: " << event->GetNtrack()
63 << " Hits A: " << event->GetNhitA()
64 << " Hits B: " << event->GetNhitB() << std::endl;
65 }
66}
67
68void pileup(Int_t nev=200)
69{
70 //make nev pileup events, each build with LOOPMAX events selected
71 //randomly among the nentries
72 TFile f("JetEvent.root");
73 auto T = f.Get<TTree>("T");
74 // Long64_t nentries = T->GetEntries();
75
76 const Int_t LOOPMAX=10;
77 JetEvent *events[LOOPMAX];
78 Int_t loop;
79 for (loop=0; loop<LOOPMAX; loop++)
80 events[loop] = 0;
81 for (Long64_t ev=0; ev<nev; ev++) {
82 if (ev%10 == 0)
83 printf("building pileup: %lld\n", ev);
84 for (loop=0; loop<LOOPMAX; loop++) {
86 T->SetBranchAddress("event", &events[loop]);
87 T->GetEntry(rev);
88 }
89 }
90}
91
92void jets(Int_t nev=100, Int_t npileup=200, Bool_t secondrun = true)
93{
94 // Embedding these loads inside the first run of the script is not yet
95 // supported in v6
96 // gROOT->ProcessLine(".L $ROOTSYS/tutorials/io/tree/JetEvent.cxx+");
97 write(nev);
98 read();
100}
101
102#else
103
104//void jets(Int_t nev=100, Int_t npileup=200, Bool_t secondrun);
105void tree130_jets(Int_t nev = 100, Int_t npileup = 200)
106{
107 TString tutdir = gROOT->GetTutorialDir();
108 gROOT->ProcessLine(".L " + tutdir + "/io/tree/JetEvent.cxx+");
109 gROOT->ProcessLine("#define JETS_SECOND_RUN yes");
110 gROOT->ProcessLine("#include \"" __FILE__ "\"");
111 gROOT->ProcessLine("jets(100, 200, true)");
112}
113
114#endif
#define f(i)
Definition RSha256.hxx:104
bool Bool_t
Definition RtypesCore.h:63
int Int_t
Definition RtypesCore.h:45
long long Long64_t
Definition RtypesCore.h:69
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
int nentries
#define gROOT
Definition TROOT.h:414
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
void Build(Int_t jetm=3, Int_t trackm=10, Int_t hitam=100, Int_t hitbm=10)
Build one event.
Definition JetEvent.cxx:46
A ROOT file is an on-disk file, usually with extension .root, that stores objects in a file-system-li...
Definition TFile.h:131
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition TRandom.cxx:682
Basic string class.
Definition TString.h:139
A TTree represents a columnar dataset.
Definition TTree.h:84
void jets()
Definition jets.C:35