ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitAna/TreeMod/src/AnaFwkMod.cc
Revision: 1.16
Committed: Fri Jul 22 18:21:02 2011 UTC (13 years, 9 months ago) by bendavid
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_025c_branch2, Mit_025c_branch1, Mit_025c_branch0, Mit_025d, Mit_025c, Mit_025b, Mit_025a, Mit_025, Mit_025pre2, Mit_024b, Mit_025pre1, Mit_024a, Mit_024
Branch point for: Mit_025c_branch
Changes since 1.15: +35 -3 lines
Log Message:
Add filling of pileup histograms

File Contents

# Content
1 // $Id: AnaFwkMod.cc,v 1.15 2011/03/11 04:03:54 bendavid Exp $
2
3 #include "MitAna/TreeMod/interface/AnaFwkMod.h"
4 #include "MitAna/DataUtil/interface/Debug.h"
5 #include "MitAna/DataTree/interface/Names.h"
6 #include "MitAna/DataTree/interface/PileupInfo.h"
7 #include <TFile.h>
8 #include <TH1D.h>
9 #include <TH3D.h>
10 #include <TStopwatch.h>
11 #include <TTree.h>
12
13 using namespace mithep;
14
15 ClassImp(mithep::AnaFwkMod)
16
17 //--------------------------------------------------------------------------------------------------
18 AnaFwkMod::AnaFwkMod(const char *name, const char *title) :
19 BaseMod(name,title),
20 fAllHeadTreeName(Names::gkAllEvtTreeName),
21 fAllHeadBrName(Names::gkAllEvtHeaderBrn),
22 fSkipNEvents(0),
23 fPrintScale(1),
24 fSWtotal(0),
25 fSWevent(0),
26 fAllHeaders(0,Names::gkSkimmedHeaders),
27 fAllHeadTree(0),
28 fAllEventHeader(0),
29 fReload(kFALSE),
30 fCurEnt(-2),
31 fNEventsSkimmed(0),
32 fNEventsSkipped(0),
33 fPileupInfoName("PileupInfo"),
34 fDoPUInfo(kFALSE),
35 hNPU(0),
36 hNPU50ns(0)
37 {
38 // Constructor.
39 }
40
41 //--------------------------------------------------------------------------------------------------
42 void AnaFwkMod::BeginRun()
43 {
44 // Get HLT tree and set branches if new file was opened. Read next entry in HLT key
45 // depending on entry in RunInfo.
46
47 if (fReload) {
48 // reset to be (re-)loaded variables
49 fReload = 0;
50 fAllHeadTree = 0;
51 fAllEventHeader = 0;
52 fCurEnt = 0;
53
54 // get current file
55 TFile *file = GetCurrentFile();
56 if (!file)
57 return;
58
59 // get all event header tree
60 fAllHeadTree = dynamic_cast<TTree*>(file->Get(fAllHeadTreeName));
61 if (!fAllHeadTree) {
62 SendError(kWarning, "BeginRun",
63 "Cannot find tree '%s' in file '%s'",
64 fAllHeadTreeName.Data(),file->GetName());
65 return;
66 }
67
68 // get all event header branch
69 if (fAllHeadTree->GetBranch(fAllHeadBrName)) {
70 fAllHeadTree->SetBranchAddress(fAllHeadBrName, &fAllEventHeader);
71 } else {
72 SendError(kWarning, "BeginRun",
73 "Cannot find branch '%s' in tree '%s'",
74 fAllHeadBrName.Data(), fAllHeadTreeName.Data());
75 return;
76 }
77 }
78 }
79
80 //--------------------------------------------------------------------------------------------------
81 void AnaFwkMod::CopyAllEventHeaders()
82 {
83 // Deal with the headers from all events if needed.
84
85 fAllHeaders.Reset();
86
87 const EventHeader *curev = GetEventHeader();
88 if (!curev) {
89 SendError(kWarning, "CopyAllEventHeaders", "Cannot obtain current event");
90 return;
91 }
92
93 if (fAllHeadTree) {
94 const Int_t nemax = fAllHeadTree->GetEntries();
95 if (fCurEnt == nemax) {
96 SendError(kAbortEvent, "CopyAllEventHeaders",
97 "End of all events tree reached (%d=%d)", fCurEnt, nemax);
98 return;
99 }
100 fAllHeadTree->GetEntry(fCurEnt++);
101 while(fCurEnt<=nemax && fAllEventHeader->Skimmed()) {
102 EventHeader *eh = fAllHeaders.AddNew();
103 eh->SetRunNum(fAllEventHeader->RunNum());
104 eh->SetEvtNum(fAllEventHeader->EvtNum());
105 eh->SetLumiSec(fAllEventHeader->LumiSec());
106 eh->SetRunEntry(fAllEventHeader->RunEntry());
107 eh->SetSkimmed(fAllEventHeader->Skimmed());
108 fAllHeadTree->GetEntry(fCurEnt++);
109 }
110 if ((fAllEventHeader->RunNum()!=curev->RunNum()) ||
111 (fAllEventHeader->EvtNum()!=curev->EvtNum()) ||
112 (fAllEventHeader->LumiSec()!=curev->LumiSec()) ||
113 (fAllEventHeader->RunEntry()!=curev->RunEntry())) {
114 SendError(kWarning, "CopyAllEventHeaders",
115 "Event header information for entry %d inconsistent: "
116 "%d==%d, %d==%d, %d==%d, %d==%d",
117 fCurEnt,
118 fAllEventHeader->RunNum(), curev->RunNum(),
119 fAllEventHeader->EvtNum(), curev->EvtNum(),
120 fAllEventHeader->LumiSec(), curev->LumiSec(),
121 fAllEventHeader->RunEntry(), curev->RunEntry());
122 return;
123 }
124
125 // read-ahead to check if more events are coming
126 if (fCurEnt<nemax) {
127 Int_t testEnt = fCurEnt;
128 fAllHeadTree->GetEntry(testEnt++);
129 while(testEnt<=nemax && fAllEventHeader->Skimmed())
130 fAllHeadTree->GetEntry(testEnt++);
131 if (testEnt==nemax+1) { // need to add remaining skimmed events
132 fAllHeadTree->GetEntry(fCurEnt++);
133 while(fCurEnt<=nemax) {
134 EventHeader *eh = fAllHeaders.AddNew();
135 eh->SetRunNum(fAllEventHeader->RunNum());
136 eh->SetEvtNum(fAllEventHeader->EvtNum());
137 eh->SetLumiSec(fAllEventHeader->LumiSec());
138 eh->SetRunEntry(fAllEventHeader->RunEntry());
139 eh->SetSkimmed(fAllEventHeader->Skimmed());
140 fAllHeadTree->GetEntry(fCurEnt++);
141 }
142 if (fCurEnt != nemax+1) {
143 SendError(kAbortEvent, "CopyAllEventHeaders",
144 "End of all events tree unexpectedly not reached (%d!=%d)", fCurEnt, nemax);
145 return;
146 }
147 }
148 }
149 }
150 }
151
152 //--------------------------------------------------------------------------------------------------
153 Bool_t AnaFwkMod::Notify()
154 {
155 // Make sure to get the new "AllEvents" tree when the file changes.
156
157 fReload = kTRUE;
158 return kTRUE;
159 }
160
161 //--------------------------------------------------------------------------------------------------
162 void AnaFwkMod::Process()
163 {
164 // Do event skipping and counting and print out timing information.
165
166 // counting processed events
167 IncNEventsProcessed();
168
169 // get skimmed event headers
170 CopyAllEventHeaders();
171 fNEventsSkimmed += fAllHeaders.GetEntries();
172
173 // check if events should be skipped
174 if (fNEventsSkipped<fSkipNEvents) {
175 ++fNEventsSkipped;
176 MDB(kAnalysis, 3) {
177 Info("Process", "Skipping (aborting) %d out of %lld first events.",
178 fNEventsSkipped, fSkipNEvents);
179 }
180 AbortEvent();
181 return;
182 }
183
184 // check if printout should be done
185 Bool_t doPrint = 0;
186 UInt_t nProcessed = fPrintScale;
187
188 MDB(kAnalysis, 4) {
189 if (GetNEventsProcessed() % (fPrintScale) == 0)
190 doPrint = 1;
191 } else {
192 MDB(kAnalysis, 3) {
193 if (GetNEventsProcessed() % (fPrintScale*10) == 0)
194 doPrint = 1;
195 nProcessed = fPrintScale*10;
196 } else {
197 MDB(kAnalysis, 2) {
198 if (GetNEventsProcessed() % (fPrintScale*100) == 0)
199 doPrint = 1;
200 nProcessed = fPrintScale*100;
201 } else {
202 MDB(kAnalysis, 1) {
203 if (GetNEventsProcessed() % (fPrintScale*1000) == 0)
204 doPrint = 1;
205 nProcessed = fPrintScale*1000;
206 }
207 }
208 }
209 }
210
211 if (GetEventHeader()->IsMC()) {
212 LoadBranch(fPileupInfoName);
213 Int_t npu[3] = {0,0,0};
214 for (UInt_t i=0; i<fPileupInfo->GetEntries(); ++i) {
215 const PileupInfo *puinfo = fPileupInfo->At(i);
216 if (puinfo->GetBunchCrossing()==0) npu[0]= puinfo->GetPU_NumInteractions();
217 else if (puinfo->GetBunchCrossing()==-1) npu[1] = puinfo->GetPU_NumInteractions();
218 else if (puinfo->GetBunchCrossing()==1) npu[2] = puinfo->GetPU_NumInteractions();
219 }
220
221 hNPU->Fill(npu[0]);
222 hNPU50ns->Fill(npu[0],npu[1],npu[2]);
223 }
224
225
226 if (doPrint) {
227 fSWevent->Stop();
228 Info("Process",
229 "Events %d -> %.2gs real, %.2gs cpu (%.2g real, %.2g cpu per event)",
230 GetNEventsProcessed(), fSWevent->RealTime(), fSWevent->CpuTime(),
231 fSWevent->RealTime()/nProcessed,
232 fSWevent->CpuTime()/nProcessed);
233 fSWevent->Start();
234 }
235
236
237 }
238
239 //--------------------------------------------------------------------------------------------------
240 void AnaFwkMod::SlaveBegin()
241 {
242 // Book our histogram and start the stop watches.
243
244 fSWtotal = new TStopwatch;
245 fSWevent = new TStopwatch;
246
247 if (!PublishObj(&fAllHeaders)) {
248 SendError(kAbortAnalysis, "SlaveBegin",
249 "Could not publish all event headers with name %s.", fAllHeaders.GetName());
250 return;
251 }
252
253 ReqBranch(fPileupInfoName, fPileupInfo);
254
255 hNPU = new TH1D("hNPU", "hNPU", 51, -0.5, 50.5);
256 AddOutput(hNPU);
257
258 hNPU50ns = new TH3D("hNPU50ns", "hNPU50ns", 51, -0.5, 50.5, 51, -0.5, 50.5, 51, -0.5, 50.5);
259 AddOutput(hNPU50ns);
260
261 }
262
263 //--------------------------------------------------------------------------------------------------
264 void AnaFwkMod::SlaveTerminate()
265 {
266 // Fill event histogram and printout timing information.
267
268 RetractObj(fAllHeaders.GetName());
269
270 SaveNEventsProcessed();
271 TH1D *hDAllEvents = new TH1D("hDAllEvents","Sum of processed and skimmed events",1,-0.5,0.5);
272 hDAllEvents->Fill(0.0,fNEventsSkimmed+GetNEventsProcessed());
273 hDAllEvents->SetEntries(fNEventsSkimmed+GetNEventsProcessed());
274 AddOutput(hDAllEvents);
275
276 TH1D *hDSkippedEvents = new TH1D("hDSkippedEvents","Number of skipped events",1,-0.5,0.5);
277 hDSkippedEvents->Fill(0.0,fNEventsSkipped);
278 hDSkippedEvents->SetEntries(fNEventsSkipped);
279 AddOutput(hDSkippedEvents);
280
281 fSWtotal->Stop();
282 fSWevent->Stop();
283
284 MDB(kAnalysis, 1)
285 Info("SlaveTerminate",
286 "Events %d -> %.2gs real, %.2gs cpu (%.2gs real, %.2gs cpu per event)",
287 GetNEventsProcessed(), fSWtotal->RealTime(), fSWtotal->CpuTime(),
288 fSWtotal->RealTime()/GetNEventsProcessed(),
289 fSWtotal->CpuTime()/GetNEventsProcessed());
290
291 delete fSWtotal;
292 delete fSWevent;
293 }