]>
Commit | Line | Data |
---|---|---|
1 | /* $Id: AliPhysicsSelection.cxx 35782 2009-10-22 11:54:31Z jgrosseo $ */ | |
2 | ||
3 | /************************************************************************** | |
4 | * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. * | |
5 | * * | |
6 | * Author: The ALICE Off-line Project. * | |
7 | * Contributors are mentioned in the code where appropriate. * | |
8 | * * | |
9 | * Permission to use, copy, modify and distribute this software and its * | |
10 | * documentation strictly for non-commercial purposes is hereby granted * | |
11 | * without fee, provided that the above copyright notice appears in all * | |
12 | * copies and that both the copyright notice and this permission notice * | |
13 | * appear in the supporting documentation. The authors make no claims * | |
14 | * about the suitability of this software for any purpose. It is * | |
15 | * provided "as is" without express or implied warranty. * | |
16 | **************************************************************************/ | |
17 | ||
18 | //------------------------------------------------------------------------- | |
19 | // Implementation of Class AliPhysicsSelection | |
20 | // This class selects collision candidates from data runs, applying selection cuts on triggers | |
21 | // and background rejection based on the content of the ESD | |
22 | // | |
23 | // Usage: | |
24 | // | |
25 | // Create the object and initialize it with the correct run number: | |
26 | // fPhysicsSelection = new AliPhysicsSelection; | |
27 | // fPhysicsSelection->Initialize(104160); | |
28 | // | |
29 | // To check if an event is a collision candidate, use: | |
30 | // fPhysicsSelection->IsCollisionCandidate(fESD) | |
31 | // | |
32 | // After processing save the resulting histograms to a file with (a folder physics_selection | |
33 | // will be created that contains the histograms): | |
34 | // fPhysicsSelection->SaveHistograms("physics_selection") | |
35 | // | |
36 | // To print statistics after processing use: | |
37 | // fPhysicsSelection->Print(); | |
38 | // | |
39 | // Origin: Jan Fiete Grosse-Oetringhaus, CERN | |
40 | //------------------------------------------------------------------------- | |
41 | ||
42 | #include <Riostream.h> | |
43 | #include <TH1F.h> | |
44 | #include <TH2F.h> | |
45 | #include <TList.h> | |
46 | #include <TIterator.h> | |
47 | #include <TDirectory.h> | |
48 | #include <TObjArray.h> | |
49 | ||
50 | #include <AliPhysicsSelection.h> | |
51 | ||
52 | #include <AliTriggerAnalysis.h> | |
53 | #include <AliLog.h> | |
54 | ||
55 | #include <AliESDEvent.h> | |
56 | ||
57 | ClassImp(AliPhysicsSelection) | |
58 | ||
59 | AliPhysicsSelection::AliPhysicsSelection() : | |
60 | AliAnalysisCuts("AliPhysicsSelection", "AliPhysicsSelection"), | |
61 | fCurrentRun(-1), | |
62 | fCollTrigClasses(), | |
63 | fBGTrigClasses(), | |
64 | fTriggerAnalysis(), | |
65 | fBackgroundIdentification(0), | |
66 | fHistStatistics(0), | |
67 | fHistBunchCrossing(0) | |
68 | { | |
69 | // constructor | |
70 | ||
71 | fCollTrigClasses.SetOwner(1); | |
72 | fBGTrigClasses.SetOwner(1); | |
73 | fTriggerAnalysis.SetOwner(1); | |
74 | ||
75 | AliLog::SetClassDebugLevel("AliPhysicsSelection", AliLog::kWarning); | |
76 | } | |
77 | ||
78 | AliPhysicsSelection::~AliPhysicsSelection() | |
79 | { | |
80 | // destructor | |
81 | ||
82 | fCollTrigClasses.Delete(); | |
83 | fBGTrigClasses.Delete(); | |
84 | fTriggerAnalysis.Delete(); | |
85 | ||
86 | if (fHistStatistics) | |
87 | { | |
88 | delete fHistStatistics; | |
89 | fHistStatistics = 0; | |
90 | } | |
91 | ||
92 | if (fHistBunchCrossing) | |
93 | { | |
94 | delete fHistBunchCrossing; | |
95 | fHistBunchCrossing = 0; | |
96 | } | |
97 | } | |
98 | ||
99 | Bool_t AliPhysicsSelection::CheckTriggerClass(const AliESDEvent* aEsd, const char* trigger) const | |
100 | { | |
101 | // checks if the given trigger class(es) are found for the current event | |
102 | // format of trigger: +TRIGGER1 -TRIGGER2 | |
103 | // requires TRIGGER1 and rejects TRIGGER2 | |
104 | ||
105 | TString str(trigger); | |
106 | TObjArray* tokens = str.Tokenize(" "); | |
107 | ||
108 | for (Int_t i=0; i < tokens->GetEntries(); i++) | |
109 | { | |
110 | TString str2(((TObjString*) tokens->At(i))->String()); | |
111 | ||
112 | if (str2[0] != '+' && str2[0] != '-') | |
113 | AliFatal(Form("Invalid trigger syntax: %s", trigger)); | |
114 | ||
115 | Bool_t flag = (str2[0] == '+'); | |
116 | ||
117 | str2.Remove(0, 1); | |
118 | ||
119 | if (flag && !aEsd->IsTriggerClassFired(str2)) | |
120 | { | |
121 | AliDebug(AliLog::kDebug, Form("Rejecting event because trigger class %s is not present", str2.Data())); | |
122 | delete tokens; | |
123 | return kFALSE; | |
124 | } | |
125 | if (!flag && aEsd->IsTriggerClassFired(str2)) | |
126 | { | |
127 | AliDebug(AliLog::kDebug, Form("Rejecting event because trigger class %s is present", str2.Data())); | |
128 | delete tokens; | |
129 | return kFALSE; | |
130 | } | |
131 | } | |
132 | ||
133 | delete tokens; | |
134 | return kTRUE; | |
135 | } | |
136 | ||
137 | Bool_t AliPhysicsSelection::IsCollisionCandidate(const AliESDEvent* aEsd) | |
138 | { | |
139 | // checks if the given event is a collision candidate | |
140 | ||
141 | const AliESDHeader* esdHeader = aEsd->GetHeader(); | |
142 | if (!esdHeader) | |
143 | { | |
144 | AliError("ESD Header could not be retrieved"); | |
145 | return kFALSE; | |
146 | } | |
147 | ||
148 | // check event type (should be PHYSICS = 7) | |
149 | if (esdHeader->GetEventType() != 7) | |
150 | return kFALSE; | |
151 | ||
152 | if (fCurrentRun != aEsd->GetRunNumber()) | |
153 | if (!Initialize(aEsd->GetRunNumber())) | |
154 | AliFatal(Form("Could not initialize for run %d", aEsd->GetRunNumber())); | |
155 | ||
156 | Bool_t accept = kFALSE; | |
157 | ||
158 | Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries(); | |
159 | for (Int_t i=0; i < count; i++) | |
160 | { | |
161 | const char* triggerClass = 0; | |
162 | if (i < fCollTrigClasses.GetEntries()) | |
163 | triggerClass = ((TObjString*) fCollTrigClasses.At(i))->String(); | |
164 | else | |
165 | triggerClass = ((TObjString*) fBGTrigClasses.At(i - fCollTrigClasses.GetEntries()))->String(); | |
166 | ||
167 | AliDebug(AliLog::kDebug, Form("Processing trigger class %s", triggerClass)); | |
168 | ||
169 | AliTriggerAnalysis* triggerAnalysis = static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i)); | |
170 | ||
171 | triggerAnalysis->FillTriggerClasses(aEsd); | |
172 | ||
173 | if (CheckTriggerClass(aEsd, triggerClass)) | |
174 | { | |
175 | triggerAnalysis->FillHistograms(aEsd); | |
176 | ||
177 | fHistStatistics->Fill(1, i); | |
178 | ||
179 | Int_t fastOR = triggerAnalysis->SPDFiredChips(aEsd, 0); | |
180 | Bool_t v0A = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0A); | |
181 | Bool_t v0C = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0C); | |
182 | Bool_t v0BG = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0ABG) || triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0CBG); | |
183 | ||
184 | if (fastOR > 0) | |
185 | fHistStatistics->Fill(2, i); | |
186 | if (fastOR > 1) | |
187 | fHistStatistics->Fill(3, i); | |
188 | ||
189 | if (v0A) | |
190 | fHistStatistics->Fill(4, i); | |
191 | if (v0C) | |
192 | fHistStatistics->Fill(5, i); | |
193 | if (v0BG) | |
194 | fHistStatistics->Fill(6, i); | |
195 | ||
196 | if ((fastOR > 0 || v0A || v0C) && !v0BG) | |
197 | fHistStatistics->Fill(7, i); | |
198 | ||
199 | if (fastOR > 0 && (v0A || v0C) && !v0BG) | |
200 | fHistStatistics->Fill(8, i); | |
201 | ||
202 | if (v0A && v0C && !v0BG) | |
203 | fHistStatistics->Fill(9, i); | |
204 | ||
205 | if (fastOR > 1 || (fastOR > 0 && (v0A || v0C)) || (v0A && v0C)) | |
206 | { | |
207 | if (!v0BG) | |
208 | { | |
209 | fHistStatistics->Fill(10, i); | |
210 | ||
211 | if (fBackgroundIdentification && !fBackgroundIdentification->IsSelected(const_cast<AliESDEvent*> (aEsd))) | |
212 | { | |
213 | AliDebug(AliLog::kDebug, "Rejecting event because of background identification"); | |
214 | fHistStatistics->Fill(11, i); | |
215 | } | |
216 | else | |
217 | { | |
218 | AliDebug(AliLog::kDebug, "Accepted event for histograms"); | |
219 | ||
220 | fHistStatistics->Fill(12, i); | |
221 | fHistBunchCrossing->Fill(aEsd->GetBunchCrossNumber(), i); | |
222 | if (i < fCollTrigClasses.GetEntries()) | |
223 | accept = kTRUE; | |
224 | } | |
225 | } | |
226 | else | |
227 | AliDebug(AliLog::kDebug, "Rejecting event because of V0 BG flag"); | |
228 | } | |
229 | else | |
230 | AliDebug(AliLog::kDebug, "Rejecting event because trigger condition is not fulfilled"); | |
231 | } | |
232 | } | |
233 | ||
234 | if (accept) | |
235 | AliDebug(AliLog::kDebug, "Accepted event as collision candidate"); | |
236 | ||
237 | return accept; | |
238 | } | |
239 | ||
240 | Bool_t AliPhysicsSelection::Initialize(UInt_t runNumber) | |
241 | { | |
242 | // initializes the object for the given run | |
243 | // TODO having the run number here and parameters hardcoded is clearly temporary, a way needs to be found to have a CDB-like configuration also for analysis | |
244 | ||
245 | if (fCurrentRun != -1) | |
246 | AliFatal("Processing several runs is not supported, yet"); | |
247 | ||
248 | AliInfo(Form("Initializing for run %d", runNumber)); | |
249 | fCurrentRun = runNumber; | |
250 | ||
251 | fTriggerAnalysis.Delete(); | |
252 | fCollTrigClasses.Delete(); | |
253 | fBGTrigClasses.Delete(); | |
254 | ||
255 | fCollTrigClasses.Add(new TObjString("+CINT1B-ABCE-NOPF-ALL")); | |
256 | fBGTrigClasses.Add(new TObjString("+CINT1A-ABCE-NOPF-ALL")); | |
257 | fBGTrigClasses.Add(new TObjString("+CINT1C-ABCE-NOPF-ALL")); | |
258 | fBGTrigClasses.Add(new TObjString("+CINT1-E-NOPF-ALL")); | |
259 | ||
260 | Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries(); | |
261 | ||
262 | for (Int_t i=0; i<count; i++) | |
263 | { | |
264 | AliTriggerAnalysis* triggerAnalysis = new AliTriggerAnalysis; | |
265 | triggerAnalysis->EnableHistograms(); | |
266 | triggerAnalysis->SetSPDGFOThreshhold(1); | |
267 | triggerAnalysis->SetV0TimeOffset(0); | |
268 | ||
269 | switch (runNumber) | |
270 | { | |
271 | case 104316: | |
272 | case 104320: | |
273 | case 104321: //OK | |
274 | case 104439: | |
275 | triggerAnalysis->SetV0TimeOffset(7.5); | |
276 | break; | |
277 | } | |
278 | ||
279 | fTriggerAnalysis.Add(triggerAnalysis); | |
280 | } | |
281 | ||
282 | if (fHistStatistics) | |
283 | delete fHistStatistics; | |
284 | ||
285 | fHistStatistics = new TH2F("fHistStatistics", "fHistStatistics;;", 12, 0.5, 12.5, count, -0.5, -0.5 + count); | |
286 | ||
287 | Int_t n = 1; | |
288 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "Correct trigger class(es)"); | |
289 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 1"); | |
290 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 2"); | |
291 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A"); | |
292 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0C"); | |
293 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0 BG"); | |
294 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "(FO >= 1 | V0A | VOC) & !V0 BG"); | |
295 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 1 & (V0A | VOC) & !V0 BG"); | |
296 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A & VOC & !V0 BG"); | |
297 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "(FO >= 2 | (FO >= 1 & (V0A | VOC)) | (V0A & VOC)) & !V0 BG"); | |
298 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "Background identification"); | |
299 | fHistStatistics->GetXaxis()->SetBinLabel(n++, "Accepted"); | |
300 | ||
301 | n = 1; | |
302 | for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++) | |
303 | fHistStatistics->GetYaxis()->SetBinLabel(n++, ((TObjString*) fCollTrigClasses.At(i))->String()); | |
304 | for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++) | |
305 | fHistStatistics->GetYaxis()->SetBinLabel(n++, ((TObjString*) fBGTrigClasses.At(i))->String()); | |
306 | ||
307 | if (fHistBunchCrossing) | |
308 | delete fHistBunchCrossing; | |
309 | ||
310 | fHistBunchCrossing = new TH2F("fHistBunchCrossing", "fHistBunchCrossing;bunch crossing number;", 4000, -0.5, 3999.5, count, -0.5, -0.5 + count); | |
311 | ||
312 | return kTRUE; | |
313 | } | |
314 | ||
315 | void AliPhysicsSelection::Print(Option_t* /* option */) const | |
316 | { | |
317 | // print the configuration | |
318 | ||
319 | Printf("Configuration initialized for run %d:", fCurrentRun); | |
320 | ||
321 | Printf("Collision trigger classes:"); | |
322 | for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++) | |
323 | Printf("%s", ((TObjString*) fCollTrigClasses.At(i))->String().Data()); | |
324 | ||
325 | Printf("Background trigger classes:"); | |
326 | for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++) | |
327 | Printf("%s", ((TObjString*) fBGTrigClasses.At(i))->String().Data()); | |
328 | ||
329 | AliTriggerAnalysis* triggerAnalysis = dynamic_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(0)); | |
330 | ||
331 | if (triggerAnalysis) | |
332 | { | |
333 | if (triggerAnalysis->GetV0TimeOffset() > 0) | |
334 | Printf("V0 time offset active: %.2f ns", triggerAnalysis->GetV0TimeOffset()); | |
335 | ||
336 | Printf("\nTotal available events:"); | |
337 | ||
338 | triggerAnalysis->PrintTriggerClasses(); | |
339 | } | |
340 | ||
341 | if (fHistStatistics) | |
342 | { | |
343 | Printf("\nSelection statistics for first collision trigger:"); | |
344 | ||
345 | Printf("Total events with correct trigger class: %d", (Int_t) fHistStatistics->GetBinContent(1, 1)); | |
346 | Printf("Selected collision candidates: %d", (Int_t) fHistStatistics->GetBinContent(12, 1)); | |
347 | } | |
348 | } | |
349 | ||
350 | Long64_t AliPhysicsSelection::Merge(TCollection* list) | |
351 | { | |
352 | // Merge a list of AliMultiplicityCorrection objects with this (needed for | |
353 | // PROOF). | |
354 | // Returns the number of merged objects (including this). | |
355 | ||
356 | if (!list) | |
357 | return 0; | |
358 | ||
359 | if (list->IsEmpty()) | |
360 | return 1; | |
361 | ||
362 | TIterator* iter = list->MakeIterator(); | |
363 | TObject* obj; | |
364 | ||
365 | // collections of all histograms | |
366 | const Int_t nHists = 9; | |
367 | TList collections[nHists]; | |
368 | ||
369 | Int_t count = 0; | |
370 | while ((obj = iter->Next())) { | |
371 | ||
372 | AliPhysicsSelection* entry = dynamic_cast<AliPhysicsSelection*> (obj); | |
373 | if (entry == 0) | |
374 | continue; | |
375 | ||
376 | Int_t n = 0; | |
377 | collections[n++].Add(&(entry->fTriggerAnalysis)); | |
378 | collections[n++].Add(entry->fHistStatistics); | |
379 | collections[n++].Add(entry->fHistBunchCrossing); | |
380 | if (entry->fBackgroundIdentification) | |
381 | collections[n++].Add(entry->fBackgroundIdentification); | |
382 | ||
383 | count++; | |
384 | } | |
385 | ||
386 | Int_t n = 0; | |
387 | fTriggerAnalysis.Merge(&collections[n++]); | |
388 | fHistStatistics->Merge(&collections[n++]); | |
389 | fHistBunchCrossing->Merge(&collections[n++]); | |
390 | if (fBackgroundIdentification) | |
391 | fBackgroundIdentification->Merge(&collections[n++]); | |
392 | ||
393 | delete iter; | |
394 | ||
395 | return count+1; | |
396 | } | |
397 | ||
398 | void AliPhysicsSelection::SaveHistograms(const char* folder) const | |
399 | { | |
400 | // write histograms to current directory | |
401 | ||
402 | if (!fHistStatistics) | |
403 | return; | |
404 | ||
405 | if (folder) | |
406 | { | |
407 | gDirectory->mkdir(folder); | |
408 | gDirectory->cd(folder); | |
409 | } | |
410 | ||
411 | fHistStatistics->Write(); | |
412 | fHistBunchCrossing->Write(); | |
413 | ||
414 | Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries(); | |
415 | for (Int_t i=0; i < count; i++) | |
416 | { | |
417 | TString triggerClass = "trigger_histograms_"; | |
418 | if (i < fCollTrigClasses.GetEntries()) | |
419 | triggerClass += ((TObjString*) fCollTrigClasses.At(i))->String(); | |
420 | else | |
421 | triggerClass += ((TObjString*) fBGTrigClasses.At(i - fCollTrigClasses.GetEntries()))->String(); | |
422 | ||
423 | gDirectory->mkdir(triggerClass); | |
424 | gDirectory->cd(triggerClass); | |
425 | ||
426 | static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i))->SaveHistograms(); | |
427 | ||
428 | gDirectory->cd(".."); | |
429 | } | |
430 | ||
431 | if (fBackgroundIdentification) | |
432 | { | |
433 | gDirectory->mkdir("background_identification"); | |
434 | gDirectory->cd("background_identification"); | |
435 | ||
436 | fBackgroundIdentification->GetOutput()->Write(); | |
437 | ||
438 | gDirectory->cd(".."); | |
439 | } | |
440 | ||
441 | if (folder) | |
442 | gDirectory->cd(".."); | |
443 | } |