Fork me on GitHub

source: git/examples/DenseTracks.C

Last change on this file was affe5c7, checked in by Michele Selvaggi <michele.selvaggi@…>, 7 years ago

fix typo

  • Property mode set to 100644
File size: 3.9 KB
Line 
1/*
2This macro shows how to access the particle-level reference for reconstructed objects.
3It is also shown how to loop over the jet constituents.
4
5root -l examples/Example3.C'("delphes_output.root")'
6*/
7
8#ifdef __CLING__
9R__LOAD_LIBRARY(libDelphes)
10#include "classes/DelphesClasses.h"
11#include "external/ExRootAnalysis/ExRootTreeReader.h"
12#include "external/ExRootAnalysis/ExRootResult.h"
13#else
14class ExRootTreeReader;
15class ExRootResult;
16#endif
17
18//------------------------------------------------------------------------------
19
20struct TestPlots
21{
22 TH1 *fhistDen;
23 TH1 *fhistNum;
24
25};
26
27//------------------------------------------------------------------------------
28
29void BookHistograms(ExRootResult *result, TestPlots *plots)
30{
31 TLegend *legend;
32 TPaveText *comment;
33
34 plots->fhistDen = result->AddHist1D(
35 "den", "",
36 "", "",
37 200, -5.0, 0.0);
38
39 plots->fhistNum = result->AddHist1D(
40 "num", "",
41 "", "",
42 200, -5.0, 0.0);
43}
44
45//------------------------------------------------------------------------------
46
47void AnalyseEvents(ExRootTreeReader *treeReader, TestPlots *plots)
48{
49 TClonesArray *branchParticle = treeReader->UseBranch("Particle");
50 TClonesArray *branchEFlowTrack = treeReader->UseBranch("EFlowTrack");
51 TClonesArray *branchEFlowPhoton = treeReader->UseBranch("EFlowPhoton");
52 TClonesArray *branchEFlowNeutralHadron = treeReader->UseBranch("EFlowNeutralHadron");
53 TClonesArray *branchJet = treeReader->UseBranch("FatJet");
54
55 Long64_t allEntries = treeReader->GetEntries();
56
57 cout << "** Chain contains " << allEntries << " events" << endl;
58
59 GenParticle *particle, *gentrack;
60 Photon *photon;
61
62 Track *track;
63 Tower *tower;
64
65 Jet *jet;
66 TObject *object;
67
68 TLorentzVector momentum;
69
70 Long64_t entry;
71
72 Int_t i, j, k, pdgCode;
73
74 // Loop over all events
75 for(entry = 0; entry < allEntries; ++entry)
76 {
77 // Load selected branches with data from specified event
78 treeReader->ReadEntry(entry);
79
80 if(entry %100 == 0) cout<<entry<< endl;
81
82 jet = (Jet*) branchJet->At(0);
83
84 // Loop over all genparticles
85 for(k = 0; k < branchParticle->GetEntriesFast(); ++k)
86 {
87
88 particle = (GenParticle*) branchParticle->At(k);
89 if (particle->Status != 1) continue;
90 if (particle->Charge == 0) continue;
91 if (particle->PT < 1.0) continue;
92
93 plots->fhistDen->Fill(TMath::Log10(particle->P4().DeltaR(jet->P4())));
94
95 // Loop over all jet's constituents
96 for(j = 0; j < jet->Constituents.GetEntriesFast(); ++j)
97 {
98 object = jet->Constituents.At(j);
99
100 // Check if the constituent is accessible
101 if(object == 0) continue;
102
103 if(object->IsA() == Track::Class())
104 {
105 track = (Track*) object;
106 gentrack = (GenParticle*) track->Particle.GetObject();
107 if(gentrack->GetUniqueID() == particle->GetUniqueID())
108 {
109 //cout<<"found: "<<gentrack->GetUniqueID()<<endl;
110 plots->fhistNum->Fill(TMath::Log10(particle->P4().DeltaR(jet->P4())));
111 }
112 }
113 }
114 }
115 }
116 //plots->fhistNum->Divide(plots->fhistDen);
117}
118
119//------------------------------------------------------------------------------
120
121void PrintHistograms(ExRootResult *result, TestPlots *plots)
122{
123 result->Print("png");
124}
125
126//------------------------------------------------------------------------------
127
128void DenseTracks(const char *inputFile)
129{
130 gSystem->Load("libDelphes");
131
132 TChain *chain = new TChain("Delphes");
133 chain->Add(inputFile);
134
135 ExRootTreeReader *treeReader = new ExRootTreeReader(chain);
136 ExRootResult *result = new ExRootResult();
137
138 TestPlots *plots = new TestPlots;
139
140 BookHistograms(result, plots);
141
142 AnalyseEvents(treeReader, plots);
143
144 PrintHistograms(result, plots);
145
146 result->Write("results.root");
147
148 cout << "** Exiting..." << endl;
149
150 delete plots;
151 delete result;
152 delete treeReader;
153 delete chain;
154}
155
156//------------------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.