1 /* $Id: AliPhysicsSelection.cxx 35782 2009-10-22 11:54:31Z jgrosseo $ */
3 /**************************************************************************
4 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
6 * Author: The ALICE Off-line Project. *
7 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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
26 // fPhysicsSelection = new AliPhysicsSelection;
29 // fPhysicsSelection->SetAnalyzeMC()
31 // To check if an event is a collision candidate, use:
32 // fPhysicsSelection->IsCollisionCandidate(fESD)
34 // After processing save the resulting histograms to a file with (a folder physics_selection
35 // will be created that contains the histograms):
36 // fPhysicsSelection->SaveHistograms("physics_selection")
38 // To print statistics after processing use:
39 // fPhysicsSelection->Print();
41 // Usually the class selects the trigger scheme by itself depending on the run number.
42 // Nevertheless, you can do that manually by calling AddCollisionTriggerClass() and AddBGTriggerClass()
44 // To define the class CINT1B-ABCE-NOPF-ALL as collision trigger (those will be accepted as
45 // collision candidates when they pass the selection):
46 // AddCollisionTriggerClass("+CINT1B-ABCE-NOPF-ALL #769 #3119");
47 // To select on bunch crossing IDs in addition, use:
48 // AddCollisionTriggerClass("+CINT1B-ABCE-NOPF-ALL #769 #3119");
49 // To define the class CINT1A-ABCE-NOPF-ALL as a background trigger (those will only be counted
50 // for the control histograms):
51 // AddBGTriggerClass("+CINT1A-ABCE-NOPF-ALL");
52 // You can also specify more than one trigger class in a string or you can require that some are *not*
53 // present. The following line would require CSMBA-ABCE-NOPF-ALL, but CSMBB-ABCE-NOPF-ALL is not allowed
55 // AddBGTriggerClass("+CSMBA-ABCE-NOPF-ALL -CSMBB-ABCE-NOPF-ALL");
57 // Origin: Jan Fiete Grosse-Oetringhaus, CERN
58 //-------------------------------------------------------------------------
60 #include <Riostream.h>
64 #include <TIterator.h>
65 #include <TDirectory.h>
66 #include <TObjArray.h>
68 #include <AliPhysicsSelection.h>
70 #include <AliTriggerAnalysis.h>
73 #include <AliESDEvent.h>
75 ClassImp(AliPhysicsSelection)
77 AliPhysicsSelection::AliPhysicsSelection() :
78 AliAnalysisCuts("AliPhysicsSelection", "AliPhysicsSelection"),
84 fBackgroundIdentification(0),
86 fHistBunchCrossing(0),
87 fSkipTriggerClassSelection(0),
88 fUsingCustomClasses(0)
92 fCollTrigClasses.SetOwner(1);
93 fBGTrigClasses.SetOwner(1);
94 fTriggerAnalysis.SetOwner(1);
96 AliLog::SetClassDebugLevel("AliPhysicsSelection", AliLog::kWarning);
99 AliPhysicsSelection::~AliPhysicsSelection()
103 fCollTrigClasses.Delete();
104 fBGTrigClasses.Delete();
105 fTriggerAnalysis.Delete();
109 delete fHistStatistics;
113 if (fHistBunchCrossing)
115 delete fHistBunchCrossing;
116 fHistBunchCrossing = 0;
120 Bool_t AliPhysicsSelection::CheckTriggerClass(const AliESDEvent* aEsd, const char* trigger) const
122 // checks if the given trigger class(es) are found for the current event
123 // format of trigger: +TRIGGER1 -TRIGGER2
124 // requires TRIGGER1 and rejects TRIGGER2
126 Bool_t foundBCRequirement = kFALSE;
127 Bool_t foundCorrectBC = kFALSE;
129 TString str(trigger);
130 TObjArray* tokens = str.Tokenize(" ");
132 for (Int_t i=0; i < tokens->GetEntries(); i++)
134 TString str2(((TObjString*) tokens->At(i))->String());
136 if (str2[0] == '+' || str2[0] == '-')
138 Bool_t flag = (str2[0] == '+');
142 if (flag && !aEsd->IsTriggerClassFired(str2))
144 AliDebug(AliLog::kDebug, Form("Rejecting event because trigger class %s is not present", str2.Data()));
148 if (!flag && aEsd->IsTriggerClassFired(str2))
150 AliDebug(AliLog::kDebug, Form("Rejecting event because trigger class %s is present", str2.Data()));
155 else if (str2[0] == '#')
157 foundBCRequirement = kTRUE;
161 Int_t bcNumber = str2.Atoi();
162 AliDebug(AliLog::kDebug, Form("Checking for bunch crossing number %d", bcNumber));
164 if (aEsd->GetBunchCrossNumber() == bcNumber)
166 foundCorrectBC = kTRUE;
167 AliDebug(AliLog::kDebug, Form("Found correct bunch crossing %d", bcNumber));
171 AliFatal(Form("Invalid trigger syntax: %s", trigger));
176 if (foundBCRequirement && !foundCorrectBC)
182 Bool_t AliPhysicsSelection::IsCollisionCandidate(const AliESDEvent* aEsd)
184 // checks if the given event is a collision candidate
186 if (fCurrentRun != aEsd->GetRunNumber())
187 if (!Initialize(aEsd->GetRunNumber()))
188 AliFatal(Form("Could not initialize for run %d", aEsd->GetRunNumber()));
190 const AliESDHeader* esdHeader = aEsd->GetHeader();
193 AliError("ESD Header could not be retrieved");
197 // check event type; should be PHYSICS = 7 for data and 0 for MC
200 if (esdHeader->GetEventType() != 7)
205 if (esdHeader->GetEventType() != 0)
206 AliFatal(Form("Invalid event type for MC: %d", esdHeader->GetEventType()));
209 Bool_t accept = kFALSE;
211 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
212 for (Int_t i=0; i < count; i++)
214 const char* triggerClass = 0;
215 if (i < fCollTrigClasses.GetEntries())
216 triggerClass = ((TObjString*) fCollTrigClasses.At(i))->String();
218 triggerClass = ((TObjString*) fBGTrigClasses.At(i - fCollTrigClasses.GetEntries()))->String();
220 AliDebug(AliLog::kDebug, Form("Processing trigger class %s", triggerClass));
222 AliTriggerAnalysis* triggerAnalysis = static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i));
224 triggerAnalysis->FillTriggerClasses(aEsd);
226 if (CheckTriggerClass(aEsd, triggerClass))
228 triggerAnalysis->FillHistograms(aEsd);
230 fHistStatistics->Fill(1, i);
232 // hardware trigger (should only remove events for MC)
233 // replay CINT1B hardware trigger
234 // TODO this has to depend on the actual hardware trigger (and that depends on the run...)
235 Int_t fastORHW = triggerAnalysis->SPDFiredChips(aEsd, 1); // SPD number of chips from trigger bits (!)
236 Bool_t v0AOnline = (triggerAnalysis->V0Trigger(aEsd, AliTriggerAnalysis::kASide, kTRUE) == AliTriggerAnalysis::kV0BB);
237 Bool_t v0COnline = (triggerAnalysis->V0Trigger(aEsd, AliTriggerAnalysis::kCSide, kTRUE) == AliTriggerAnalysis::kV0BB);
239 if (fastORHW == 0 && !v0AOnline && !v0COnline)
241 AliDebug(AliLog::kDebug, "Rejecting event because hardware trigger is not fired");
245 fHistStatistics->Fill(2, i);
248 AliTriggerAnalysis::V0Decision v0ADecision = triggerAnalysis->V0Trigger(aEsd, AliTriggerAnalysis::kASide, kFALSE);
249 AliTriggerAnalysis::V0Decision v0CDecision = triggerAnalysis->V0Trigger(aEsd, AliTriggerAnalysis::kCSide, kFALSE);
251 Bool_t v0AOffline = (v0ADecision == AliTriggerAnalysis::kV0BB);
252 Bool_t v0COffline = (v0CDecision == AliTriggerAnalysis::kV0BB);
254 Int_t fastOROffline = triggerAnalysis->SPDFiredChips(aEsd, 0); // SPD number of chips from clusters (!)
255 Bool_t v0ABG = (v0ADecision == AliTriggerAnalysis::kV0BG);
256 Bool_t v0CBG = (v0CDecision == AliTriggerAnalysis::kV0BG);
257 Bool_t v0BG = v0ABG || v0CBG;
259 if (fastOROffline > 0)
260 fHistStatistics->Fill(3, i);
261 if (fastOROffline > 1)
262 fHistStatistics->Fill(4, i);
265 fHistStatistics->Fill(5, i);
267 fHistStatistics->Fill(6, i);
269 fHistStatistics->Fill(7, i);
271 fHistStatistics->Fill(8, i);
273 if (v0ADecision == AliTriggerAnalysis::kV0Fake)
274 fHistStatistics->Fill(9, i);
275 if (v0CDecision == AliTriggerAnalysis::kV0Fake)
276 fHistStatistics->Fill(10, i);
278 if (fastOROffline > 1 && !v0BG)
279 fHistStatistics->Fill(11, i);
281 if ((fastOROffline > 0 || v0AOffline || v0COffline) && !v0BG)
282 fHistStatistics->Fill(12, i);
284 if (fastOROffline > 0 && (v0AOffline || v0COffline) && !v0BG)
285 fHistStatistics->Fill(13, i);
287 if (v0AOffline && v0COffline && !v0BG)
288 fHistStatistics->Fill(14, i);
290 if (fastOROffline > 1 || (fastOROffline > 0 && (v0AOffline || v0COffline)) || (v0AOffline && v0COffline))
294 fHistStatistics->Fill(15, i);
296 if (fBackgroundIdentification && !fBackgroundIdentification->IsSelected(const_cast<AliESDEvent*> (aEsd)))
298 AliDebug(AliLog::kDebug, "Rejecting event because of background identification");
299 fHistStatistics->Fill(16, i);
303 AliDebug(AliLog::kDebug, "Accepted event for histograms");
305 fHistStatistics->Fill(17, i);
306 fHistBunchCrossing->Fill(aEsd->GetBunchCrossNumber(), i);
307 if (i < fCollTrigClasses.GetEntries() || fSkipTriggerClassSelection)
312 AliDebug(AliLog::kDebug, "Rejecting event because of V0 BG flag");
315 AliDebug(AliLog::kDebug, "Rejecting event because trigger condition is not fulfilled");
320 AliDebug(AliLog::kDebug, "Accepted event as collision candidate");
325 Int_t AliPhysicsSelection::GetTriggerScheme(UInt_t runNumber)
327 // returns the current trigger scheme (classes that are accepted/rejected)
332 // TODO dependent on run number
343 // default: CINT1 suite
347 Bool_t AliPhysicsSelection::Initialize(UInt_t runNumber)
349 // initializes the object for the given run
350 // 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
352 Bool_t oldStatus = TH1::AddDirectoryStatus();
353 TH1::AddDirectory(kFALSE);
355 Int_t triggerScheme = GetTriggerScheme(runNumber);
357 if (!fUsingCustomClasses && fCurrentRun != -1 && triggerScheme != GetTriggerScheme(fCurrentRun))
358 AliFatal("Processing several runs with different trigger schemes is not supported");
360 AliInfo(Form("Initializing for run %d", runNumber));
362 // initialize first time?
363 if (fCurrentRun == -1)
365 if (fUsingCustomClasses) {
366 AliInfo("Using user-provided trigger classes");
368 switch (triggerScheme)
371 fCollTrigClasses.Add(new TObjString(""));
375 fCollTrigClasses.Add(new TObjString("+CINT1B-ABCE-NOPF-ALL"));
376 fBGTrigClasses.Add(new TObjString("+CINT1A-ABCE-NOPF-ALL"));
377 fBGTrigClasses.Add(new TObjString("+CINT1C-ABCE-NOPF-ALL"));
378 fBGTrigClasses.Add(new TObjString("+CINT1-E-NOPF-ALL"));
382 fCollTrigClasses.Add(new TObjString("+CSMBB-ABCE-NOPF-ALL"));
383 fBGTrigClasses.Add(new TObjString("+CSMBA-ABCE-NOPF-ALL -CSMBB-ABCE-NOPF-ALL"));
384 fBGTrigClasses.Add(new TObjString("+CSMBC-ABCE-NOPF-ALL -CSMBB-ABCE-NOPF-ALL"));
388 AliFatal(Form("Unsupported trigger scheme %d", triggerScheme));
392 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
394 for (Int_t i=0; i<count; i++)
396 AliTriggerAnalysis* triggerAnalysis = new AliTriggerAnalysis;
397 triggerAnalysis->SetAnalyzeMC(fMC);
398 triggerAnalysis->EnableHistograms();
399 triggerAnalysis->SetSPDGFOThreshhold(1);
400 fTriggerAnalysis.Add(triggerAnalysis);
404 delete fHistStatistics;
406 fHistStatistics = new TH2F("fHistStatistics", "fHistStatistics;;", 17, 0.5, 17.5, count, -0.5, -0.5 + count);
409 fHistStatistics->GetXaxis()->SetBinLabel(n++, "Trigger class");
410 fHistStatistics->GetXaxis()->SetBinLabel(n++, "Hardware trigger");
411 fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 1");
412 fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 2");
413 fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A");
414 fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0C");
415 fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A BG");
416 fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0C BG");
417 fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A Fake");
418 fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0C Fake");
419 fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 2 &!V0 BG");
420 fHistStatistics->GetXaxis()->SetBinLabel(n++, "(FO >= 1 | V0A | V0C) & !V0 BG");
421 fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 1 & (V0A | V0C) & !V0 BG");
422 fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A & V0C & !V0 BG");
423 fHistStatistics->GetXaxis()->SetBinLabel(n++, "(FO >= 2 | (FO >= 1 & (V0A | V0C)) | (V0A & V0C)) & !V0 BG");
424 fHistStatistics->GetXaxis()->SetBinLabel(n++, "Background identification");
425 fHistStatistics->GetXaxis()->SetBinLabel(n++, "Accepted");
427 if (fHistBunchCrossing)
428 delete fHistBunchCrossing;
430 fHistBunchCrossing = new TH2F("fHistBunchCrossing", "fHistBunchCrossing;bunch crossing number;", 4000, -0.5, 3999.5, count, -0.5, -0.5 + count);
433 for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++)
435 fHistStatistics->GetYaxis()->SetBinLabel(n, ((TObjString*) fCollTrigClasses.At(i))->String());
436 fHistBunchCrossing->GetYaxis()->SetBinLabel(n, ((TObjString*) fCollTrigClasses.At(i))->String());
439 for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++)
441 fHistStatistics->GetYaxis()->SetBinLabel(n, ((TObjString*) fBGTrigClasses.At(i))->String());
442 fHistBunchCrossing->GetYaxis()->SetBinLabel(n, ((TObjString*) fBGTrigClasses.At(i))->String());
447 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
448 for (Int_t i=0; i<count; i++)
450 AliTriggerAnalysis* triggerAnalysis = static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i));
459 triggerAnalysis->SetV0TimeOffset(7.5);
462 triggerAnalysis->SetV0TimeOffset(0);
466 fCurrentRun = runNumber;
468 TH1::AddDirectory(oldStatus);
473 void AliPhysicsSelection::Print(Option_t* /* option */) const
475 // print the configuration
477 Printf("Configuration initialized for run %d (MC: %d):", fCurrentRun, fMC);
479 Printf("Collision trigger classes:");
480 for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++)
481 Printf("%s", ((TObjString*) fCollTrigClasses.At(i))->String().Data());
483 Printf("Background trigger classes:");
484 for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++)
485 Printf("%s", ((TObjString*) fBGTrigClasses.At(i))->String().Data());
487 AliTriggerAnalysis* triggerAnalysis = dynamic_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(0));
491 if (triggerAnalysis->GetV0TimeOffset() > 0)
492 Printf("V0 time offset active: %.2f ns", triggerAnalysis->GetV0TimeOffset());
494 Printf("\nTotal available events:");
496 triggerAnalysis->PrintTriggerClasses();
499 if (fHistStatistics && fCollTrigClasses.GetEntries() > 0)
501 Printf("\nSelection statistics for first collision trigger (%s):", ((TObjString*) fCollTrigClasses.First())->String().Data());
503 Printf("Total events with correct trigger class: %d", (Int_t) fHistStatistics->GetBinContent(1, 1));
504 Printf("Selected collision candidates: %d", (Int_t) fHistStatistics->GetBinContent(fHistStatistics->GetXaxis()->FindBin("Accepted"), 1));
507 if (fHistBunchCrossing)
509 Printf("\nBunch crossing statistics:");
511 for (Int_t i=1; i<=fHistBunchCrossing->GetNbinsY(); i++)
514 str.Form("Trigger %s has accepted events in the bunch crossings: ", fHistBunchCrossing->GetYaxis()->GetBinLabel(i));
516 for (Int_t j=1; j<=fHistBunchCrossing->GetNbinsX(); j++)
517 if (fHistBunchCrossing->GetBinContent(j, i) > 0)
518 str += Form("%d, ", (Int_t) fHistBunchCrossing->GetXaxis()->GetBinCenter(j));
520 Printf("%s", str.Data());
523 for (Int_t j=1; j<=fHistBunchCrossing->GetNbinsX(); j++)
526 for (Int_t i=1; i<=fHistBunchCrossing->GetNbinsY(); i++)
528 if (fHistBunchCrossing->GetBinContent(j, i) > 0)
532 Printf("WARNING: Bunch crossing %d has more than one trigger class active. Check BPTX functioning for this run!", (Int_t) fHistBunchCrossing->GetXaxis()->GetBinCenter(j));
536 if (fUsingCustomClasses)
537 Printf("WARNING: Using custom trigger classes!");
538 if (fSkipTriggerClassSelection)
539 Printf("WARNING: Skipping trigger class selection!");
542 Long64_t AliPhysicsSelection::Merge(TCollection* list)
544 // Merge a list of AliMultiplicityCorrection objects with this (needed for
546 // Returns the number of merged objects (including this).
554 TIterator* iter = list->MakeIterator();
557 // collections of all histograms
558 const Int_t nHists = 9;
559 TList collections[nHists];
562 while ((obj = iter->Next())) {
564 AliPhysicsSelection* entry = dynamic_cast<AliPhysicsSelection*> (obj);
568 collections[0].Add(&(entry->fTriggerAnalysis));
569 if (entry->fHistStatistics)
570 collections[1].Add(entry->fHistStatistics);
571 if (entry->fHistBunchCrossing)
572 collections[2].Add(entry->fHistBunchCrossing);
573 if (entry->fBackgroundIdentification)
574 collections[3].Add(entry->fBackgroundIdentification);
579 fTriggerAnalysis.Merge(&collections[0]);
581 fHistStatistics->Merge(&collections[1]);
582 if (fHistBunchCrossing)
583 fHistBunchCrossing->Merge(&collections[2]);
584 if (fBackgroundIdentification)
585 fBackgroundIdentification->Merge(&collections[3]);
592 void AliPhysicsSelection::SaveHistograms(const char* folder) const
594 // write histograms to current directory
596 if (!fHistStatistics)
601 gDirectory->mkdir(folder);
602 gDirectory->cd(folder);
605 fHistStatistics->Write();
606 fHistBunchCrossing->Write();
608 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
609 for (Int_t i=0; i < count; i++)
611 TString triggerClass = "trigger_histograms_";
612 if (i < fCollTrigClasses.GetEntries())
613 triggerClass += ((TObjString*) fCollTrigClasses.At(i))->String();
615 triggerClass += ((TObjString*) fBGTrigClasses.At(i - fCollTrigClasses.GetEntries()))->String();
617 gDirectory->mkdir(triggerClass);
618 gDirectory->cd(triggerClass);
620 static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i))->SaveHistograms();
622 gDirectory->cd("..");
625 if (fBackgroundIdentification)
627 gDirectory->mkdir("background_identification");
628 gDirectory->cd("background_identification");
630 fBackgroundIdentification->GetOutput()->Write();
632 gDirectory->cd("..");
636 gDirectory->cd("..");