added separate histogram entry for V0A and C BG flag
[u/mrichter/AliRoot.git] / PWG0 / AliPhysicsSelection.cxx
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 v0ABG = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0ABG);
183       Bool_t v0CBG = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0CBG);
184       Bool_t v0BG = v0ABG || v0CBG;
185     
186       if (fastOR > 0)
187         fHistStatistics->Fill(2, i);
188       if (fastOR > 1)
189         fHistStatistics->Fill(3, i);
190         
191       if (v0A)
192         fHistStatistics->Fill(4, i);
193       if (v0C)
194         fHistStatistics->Fill(5, i);
195       if (v0ABG)
196         fHistStatistics->Fill(6, i);
197       if (v0CBG)
198         fHistStatistics->Fill(7, i);
199         
200       if ((fastOR > 0 || v0A || v0C) && !v0BG)
201         fHistStatistics->Fill(8, i);
202     
203       if (fastOR > 0 && (v0A || v0C) && !v0BG)
204         fHistStatistics->Fill(9, i);
205   
206       if (v0A && v0C && !v0BG)
207         fHistStatistics->Fill(10, i);
208         
209       if (fastOR > 1 || (fastOR > 0 && (v0A || v0C)) || (v0A && v0C))
210       {
211         if (!v0BG)
212         {
213           fHistStatistics->Fill(11, i);
214       
215           if (fBackgroundIdentification && !fBackgroundIdentification->IsSelected(const_cast<AliESDEvent*> (aEsd)))
216           {
217             AliDebug(AliLog::kDebug, "Rejecting event because of background identification");
218             fHistStatistics->Fill(12, i);
219           }
220           else
221           {
222             AliDebug(AliLog::kDebug, "Accepted event for histograms");
223             
224             fHistStatistics->Fill(13, i);
225             fHistBunchCrossing->Fill(aEsd->GetBunchCrossNumber(), i);
226             if (i < fCollTrigClasses.GetEntries())
227               accept = kTRUE;
228           }
229         }
230         else
231           AliDebug(AliLog::kDebug, "Rejecting event because of V0 BG flag");
232       }
233       else
234         AliDebug(AliLog::kDebug, "Rejecting event because trigger condition is not fulfilled");
235     }
236   }
237  
238   if (accept)
239     AliDebug(AliLog::kDebug, "Accepted event as collision candidate");
240   
241   return accept;
242 }
243     
244 Bool_t AliPhysicsSelection::Initialize(UInt_t runNumber)
245 {
246   // initializes the object for the given run
247   // 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
248   
249   if (fCurrentRun != -1)
250     AliFatal("Processing several runs is not supported, yet");
251   
252   AliInfo(Form("Initializing for run %d", runNumber));
253   fCurrentRun = runNumber;
254   
255   fTriggerAnalysis.Delete();
256   fCollTrigClasses.Delete();
257   fBGTrigClasses.Delete();
258   
259   fCollTrigClasses.Add(new TObjString("+CINT1B-ABCE-NOPF-ALL"));
260   fBGTrigClasses.Add(new TObjString("+CINT1A-ABCE-NOPF-ALL"));
261   fBGTrigClasses.Add(new TObjString("+CINT1C-ABCE-NOPF-ALL"));
262   fBGTrigClasses.Add(new TObjString("+CINT1-E-NOPF-ALL"));
263   
264   Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
265   
266   for (Int_t i=0; i<count; i++)
267   {
268     AliTriggerAnalysis* triggerAnalysis = new AliTriggerAnalysis;
269     triggerAnalysis->EnableHistograms();
270     triggerAnalysis->SetSPDGFOThreshhold(1);
271     triggerAnalysis->SetV0TimeOffset(0);
272     
273     switch (runNumber)
274     {
275       case 104316:
276       case 104320:
277       case 104321: //OK
278       case 104439:
279         triggerAnalysis->SetV0TimeOffset(7.5);
280         break;
281     }
282   
283     fTriggerAnalysis.Add(triggerAnalysis);
284   }
285       
286   if (fHistStatistics)
287     delete fHistStatistics;
288
289   fHistStatistics = new TH2F("fHistStatistics", "fHistStatistics;;", 13, 0.5, 13.5, count, -0.5, -0.5 + count);
290     
291   Int_t n = 1;
292   fHistStatistics->GetXaxis()->SetBinLabel(n++, "Correct trigger class(es)");
293   fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 1");
294   fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 2");
295   fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A");
296   fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0C");
297   fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A BG");
298   fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0C BG");
299   fHistStatistics->GetXaxis()->SetBinLabel(n++, "(FO >= 1 | V0A | VOC) & !V0 BG");
300   fHistStatistics->GetXaxis()->SetBinLabel(n++, "FO >= 1 & (V0A | VOC) & !V0 BG");
301   fHistStatistics->GetXaxis()->SetBinLabel(n++, "V0A & VOC & !V0 BG");
302   fHistStatistics->GetXaxis()->SetBinLabel(n++, "(FO >= 2 | (FO >= 1 & (V0A | VOC)) | (V0A & VOC)) & !V0 BG");
303   fHistStatistics->GetXaxis()->SetBinLabel(n++, "Background identification");
304   fHistStatistics->GetXaxis()->SetBinLabel(n++, "Accepted");
305   
306   if (fHistBunchCrossing)
307     delete fHistBunchCrossing;
308
309   fHistBunchCrossing = new TH2F("fHistBunchCrossing", "fHistBunchCrossing;bunch crossing number;", 4000, -0.5, 3999.5,  count, -0.5, -0.5 + count);
310     
311   n = 1;
312   for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++)
313   {
314     fHistStatistics->GetYaxis()->SetBinLabel(n, ((TObjString*) fCollTrigClasses.At(i))->String());
315     fHistBunchCrossing->GetYaxis()->SetBinLabel(n, ((TObjString*) fCollTrigClasses.At(i))->String());
316     n++;
317   }
318   for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++)
319   {
320     fHistStatistics->GetYaxis()->SetBinLabel(n, ((TObjString*) fBGTrigClasses.At(i))->String());
321     fHistBunchCrossing->GetYaxis()->SetBinLabel(n, ((TObjString*) fBGTrigClasses.At(i))->String());
322     n++;
323   }
324     
325   return kTRUE;
326 }
327
328 void AliPhysicsSelection::Print(Option_t* /* option */) const
329 {
330   // print the configuration
331   
332   Printf("Configuration initialized for run %d:", fCurrentRun);
333   
334   Printf("Collision trigger classes:");
335   for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++)
336     Printf("%s", ((TObjString*) fCollTrigClasses.At(i))->String().Data());
337   
338   Printf("Background trigger classes:");
339   for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++)
340     Printf("%s", ((TObjString*) fBGTrigClasses.At(i))->String().Data());
341
342   AliTriggerAnalysis* triggerAnalysis = dynamic_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(0));
343   
344   if (triggerAnalysis)
345   {
346     if (triggerAnalysis->GetV0TimeOffset() > 0)
347       Printf("V0 time offset active: %.2f ns", triggerAnalysis->GetV0TimeOffset());
348     
349     Printf("\nTotal available events:");
350     
351     triggerAnalysis->PrintTriggerClasses();
352   }
353   
354   if (fHistStatistics)
355   {
356     Printf("\nSelection statistics for first collision trigger:");
357     
358     Printf("Total events with correct trigger class: %d", (Int_t) fHistStatistics->GetBinContent(1, 1));
359     Printf("Selected collision candidates: %d", (Int_t) fHistStatistics->GetBinContent(13, 1));
360   }
361 }
362
363 Long64_t AliPhysicsSelection::Merge(TCollection* list)
364 {
365   // Merge a list of AliMultiplicityCorrection objects with this (needed for
366   // PROOF).
367   // Returns the number of merged objects (including this).
368
369   if (!list)
370     return 0;
371
372   if (list->IsEmpty())
373     return 1;
374
375   TIterator* iter = list->MakeIterator();
376   TObject* obj;
377
378   // collections of all histograms
379   const Int_t nHists = 9;
380   TList collections[nHists];
381
382   Int_t count = 0;
383   while ((obj = iter->Next())) {
384
385     AliPhysicsSelection* entry = dynamic_cast<AliPhysicsSelection*> (obj);
386     if (entry == 0) 
387       continue;
388
389     Int_t n = 0;
390     collections[n++].Add(&(entry->fTriggerAnalysis));
391     collections[n++].Add(entry->fHistStatistics);
392     collections[n++].Add(entry->fHistBunchCrossing);
393     if (entry->fBackgroundIdentification)
394       collections[n++].Add(entry->fBackgroundIdentification);
395
396     count++;
397   }
398
399   Int_t n = 0;
400   fTriggerAnalysis.Merge(&collections[n++]);
401   fHistStatistics->Merge(&collections[n++]);
402   fHistBunchCrossing->Merge(&collections[n++]);
403   if (fBackgroundIdentification)
404     fBackgroundIdentification->Merge(&collections[n++]);
405   
406   delete iter;
407
408   return count+1;
409 }
410
411 void AliPhysicsSelection::SaveHistograms(const char* folder) const
412 {
413   // write histograms to current directory
414   
415   if (!fHistStatistics)
416     return;
417     
418   if (folder)
419   {
420     gDirectory->mkdir(folder);
421     gDirectory->cd(folder);
422   }
423   
424   fHistStatistics->Write();
425   fHistBunchCrossing->Write();
426   
427   Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
428   for (Int_t i=0; i < count; i++)
429   {
430     TString triggerClass = "trigger_histograms_";
431     if (i < fCollTrigClasses.GetEntries())
432       triggerClass += ((TObjString*) fCollTrigClasses.At(i))->String();
433     else
434       triggerClass += ((TObjString*) fBGTrigClasses.At(i - fCollTrigClasses.GetEntries()))->String();
435   
436     gDirectory->mkdir(triggerClass);
437     gDirectory->cd(triggerClass);
438   
439     static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i))->SaveHistograms();
440     
441     gDirectory->cd("..");
442   }
443  
444   if (fBackgroundIdentification)
445   {
446     gDirectory->mkdir("background_identification");
447     gDirectory->cd("background_identification");
448       
449     fBackgroundIdentification->GetOutput()->Write();
450       
451     gDirectory->cd("..");
452   }
453   
454   if (folder)
455     gDirectory->cd("..");
456 }