b7fbded40f455bf772040b98382eacd72d2e4f71
[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 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 }