Corrected cast
[u/mrichter/AliRoot.git] / PWG0 / AliTriggerAnalysis.cxx
1 /* $Id: AliTriggerAnalysis.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 AliTriggerAnalysis
20 //   This class provides function to check if events have been triggered based on the data in the ESD
21 //   The trigger bits, trigger class inputs and only the data (offline trigger) can be used
22 //   Origin: Jan Fiete Grosse-Oetringhaus, CERN
23 //-------------------------------------------------------------------------
24
25 #include <Riostream.h>
26 #include <TH1F.h>
27 #include <TH2F.h>
28 #include <TList.h>
29 #include <TIterator.h>
30 #include "TParameter.h"
31 #include <TMap.h>
32
33 #include <AliTriggerAnalysis.h>
34
35 #include <AliLog.h>
36
37 #include <AliESDEvent.h>
38
39 #include <AliMultiplicity.h>
40 #include <AliESDVZERO.h>
41 #include <AliESDZDC.h>
42 #include <AliESDFMD.h>
43
44 ClassImp(AliTriggerAnalysis)
45
46 AliTriggerAnalysis::AliTriggerAnalysis() :
47   fSPDGFOThreshold(2),
48   fV0AThreshold(1),
49   fV0CThreshold(1),
50   fFMDLowCut(0.2),
51   fFMDHitCut(0.5),
52   fHistBitsSPD(0),
53   fHistFiredBitsSPD(0),
54   fHistV0A(0),       
55   fHistV0C(0),    
56   fHistZDC(0),    
57   fHistFMDA(0),    
58   fHistFMDC(0),   
59   fHistFMDSingle(0),
60   fHistFMDSum(0),
61   fTriggerClasses(0)
62 {
63 }
64
65 void AliTriggerAnalysis::EnableHistograms()
66 {
67   // creates the monitoring histograms
68   
69   fHistBitsSPD = new TH2F("fHistBitsSPD", "SPD GFO;number of fired chips (offline);number of fired chips (hardware)", 1202, -1.5, 1200.5, 1202, -1.5, 1200.5);
70   fHistFiredBitsSPD = new TH1F("fHistFiredBitsSPD", "SPD GFO Hardware;chip number;events", 1200, -0.5, 1199.5);
71   fHistV0A = new TH1F("fHistV0A", "V0A;leading time (ns);events", 200, 0, 100);
72   fHistV0C = new TH1F("fHistV0C", "V0C;leading time (ns);events", 200, 0, 100);
73   fHistZDC = new TH1F("fHistZDC", "ZDC;trigger bits;events", 8, -1.5, 6.5);
74   
75   // TODO check limits
76   fHistFMDA = new TH1F("fHistFMDA", "FMDA;combinations above threshold;events", 102, -1.5, 100.5);
77   fHistFMDC = new TH1F("fHistFMDC", "FMDC;combinations above threshold;events", 102, -1.5, 100.5);
78   fHistFMDSingle = new TH1F("fHistFMDSingle", "FMD single;multiplicity value;counts", 1000, 0, 10);
79   fHistFMDSum = new TH1F("fHistFMDSum", "FMD sum;multiplicity value;counts", 1000, 0, 10);
80   
81   fTriggerClasses = new TMap;
82 }
83
84 //____________________________________________________________________
85 const char* AliTriggerAnalysis::GetTriggerName(Trigger trigger) 
86 {
87   // returns the name of the requested trigger
88   // the returned string will only be valid until the next call to this function [not thread-safe]
89   
90   static TString str;
91   
92   UInt_t triggerNoFlags = (UInt_t) trigger % (UInt_t) kStartOfFlags;
93   
94   switch (triggerNoFlags)
95   {
96     case kAcceptAll : str = "ACCEPT ALL (bypass!)"; break;
97     case kMB1 : str = "MB1"; break;
98     case kMB2 : str = "MB2"; break;
99     case kMB3 : str = "MB3"; break;
100     case kSPDGFO : str = "SPD GFO"; break;
101     case kSPDGFOBits : str = "SPD GFO Bits"; break;
102     case kV0A : str = "V0 A"; break;
103     case kV0C : str = "V0 C"; break;
104     case kZDC : str = "ZDC"; break;
105     case kZDCA : str = "ZDC A"; break;
106     case kZDCC : str = "ZDC C"; break;
107     case kFMDA : str = "FMD A"; break;
108     case kFMDC : str = "FMD C"; break;
109     case kFPANY : str = "SPD GFO | V0 | ZDC | FMD"; break;
110     default: str = ""; break;
111   }
112    
113   if (trigger & kOfflineFlag)
114     str += " OFFLINE";  
115   
116   return str;
117 }
118
119 Bool_t AliTriggerAnalysis::IsTriggerFired(const AliESDEvent* aEsd, Trigger trigger) const
120 {
121   // checks if an event has been triggered
122
123   if (trigger & kOfflineFlag)
124     return IsOfflineTriggerFired(aEsd, trigger);
125     
126   return IsTriggerBitFired(aEsd, trigger);
127 }
128
129 Bool_t AliTriggerAnalysis::IsTriggerBitFired(const AliESDEvent* aEsd, Trigger trigger) const
130 {
131   // checks if an event is fired using the trigger bits
132
133   return IsTriggerBitFired(aEsd->GetTriggerMask(), trigger);
134 }
135
136 Bool_t AliTriggerAnalysis::IsTriggerBitFired(ULong64_t triggerMask, Trigger trigger) const
137 {
138   // checks if an event is fired using the trigger bits
139   //
140   // this function needs the branch TriggerMask in the ESD
141   
142   // definitions from p-p.cfg
143   ULong64_t spdFO = (1 << 14);
144   ULong64_t v0left = (1 << 10);
145   ULong64_t v0right = (1 << 11);
146
147   switch (trigger)
148   {
149     case kAcceptAll:
150     {
151       return kTRUE;
152       break;
153     }
154     case kMB1:
155     {
156       if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
157         return kTRUE;
158       break;
159     }
160     case kMB2:
161     {
162       if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
163         return kTRUE;
164       break;
165     }
166     case kMB3:
167     {
168       if (triggerMask & spdFO && (triggerMask & v0left) && (triggerMask & v0right))
169         return kTRUE;
170       break;
171     }
172     case kSPDGFO:
173     {
174       if (triggerMask & spdFO)
175         return kTRUE;
176       break;
177     }
178     default:
179       Printf("IsEventTriggered: ERROR: Trigger type %d not implemented in this method", (Int_t) trigger);
180       break;
181   }
182
183   return kFALSE;
184 }
185
186 Bool_t AliTriggerAnalysis::IsTriggerBitFired(const AliESDEvent* aEsd, ULong64_t tclass) const
187 {
188   // Checks if corresponding bit in mask is on
189   
190   ULong64_t trigmask = aEsd->GetTriggerMask();
191   return (trigmask & (1ull << (tclass-1)));
192 }
193
194 Bool_t AliTriggerAnalysis::IsOfflineTriggerFired(const AliESDEvent* aEsd, Trigger trigger) const
195 {
196   // checks if an event has been triggered "offline"
197
198   UInt_t triggerNoFlags = (UInt_t) trigger % (UInt_t) kStartOfFlags;
199   
200   switch (triggerNoFlags)
201   {
202     case kAcceptAll:
203     {
204       return kTRUE;
205       break;
206     }
207     case kMB1:
208     {
209       if (SPDGFOTrigger(aEsd, 0) || V0Trigger(aEsd, kASide) == kV0BB || V0Trigger(aEsd, kCSide) == kV0BB)
210         return kTRUE;
211       break;
212     }
213     case kMB2:
214     {
215       if (SPDGFOTrigger(aEsd, 0) && (V0Trigger(aEsd, kASide) == kV0BB || V0Trigger(aEsd, kCSide) == kV0BB))
216         return kTRUE;
217       break;
218     }
219     case kMB3:
220     {
221       if (SPDGFOTrigger(aEsd, 0) && V0Trigger(aEsd, kASide) == kV0BB && V0Trigger(aEsd, kCSide) == kV0BB)
222         return kTRUE;
223       break;
224     }
225     case kSPDGFO:
226     {
227       if (SPDGFOTrigger(aEsd, 0))
228         return kTRUE;
229       break;
230     }
231     case kSPDGFOBits:
232     {
233       if (SPDGFOTrigger(aEsd, 1))
234         return kTRUE;
235       break;
236     }
237     case kV0A:
238     {
239       if (V0Trigger(aEsd, kASide) == kV0BB)
240         return kTRUE;
241       break;
242     }
243     case kV0C:
244     {
245       if (V0Trigger(aEsd, kCSide) == kV0BB)
246         return kTRUE;
247       break;
248     }
249     case kZDC:
250     {
251       if (ZDCTrigger(aEsd, kASide) || ZDCTrigger(aEsd, kCentralBarrel) || ZDCTrigger(aEsd, kCSide))
252         return kTRUE;
253       break;
254     }
255     case kZDCA:
256     {
257       if (ZDCTrigger(aEsd, kASide))
258         return kTRUE;
259       break;
260     }
261     case kZDCC:
262     {
263       if (ZDCTrigger(aEsd, kCSide))
264         return kTRUE;
265       break;
266     }
267     case kFMDA:
268     {
269       if (FMDTrigger(aEsd, kASide))
270         return kTRUE;
271       break;
272     }
273     case kFMDC:
274     {
275       if (FMDTrigger(aEsd, kCSide))
276         return kTRUE;
277       break;
278     }
279     case kFPANY:
280     {
281       if (SPDGFOTrigger(aEsd, 0) || V0Trigger(aEsd, kASide) == kV0BB || V0Trigger(aEsd, kCSide) == kV0BB || ZDCTrigger(aEsd, kASide) || ZDCTrigger(aEsd, kCentralBarrel) || ZDCTrigger(aEsd, kCSide) || FMDTrigger(aEsd, kASide) || FMDTrigger(aEsd, kCSide))
282         return kTRUE;
283       break;
284     }
285     default:
286     {
287       AliFatal(Form("Trigger type %d not implemented", triggerNoFlags));
288     }
289   }
290   
291   return kFALSE;
292 }
293
294
295 Bool_t AliTriggerAnalysis::IsTriggerClassFired(const AliESDEvent* aEsd, const Char_t* tclass) const 
296 {
297   // tclass is logical function of inputs, e.g. 01 && 02 || 03 && 11 && 21
298   // = L0 inp 1 && L0 inp 2 || L0 inp 3 && L1 inp 1 && L2 inp 1
299   // NO brackets in logical function !
300   // Spaces between operators and inputs.
301   // Not all logical functions are available in CTP= 
302   // =any function of first 4 inputs; 'AND' of other inputs, check not done
303   // This method will be replaced/complemened by similar one
304   // which works withh class and inputs names as in CTP cfg file
305   
306   TString TClass(tclass);
307   TObjArray* tcltokens = TClass.Tokenize(" ");
308   Char_t level=((TObjString*)tcltokens->At(0))->String()[0];
309   UInt_t input=atoi((((TObjString*)tcltokens->At(0))->String()).Remove(0));
310   Bool_t tcl = IsInputFired(aEsd,level,input);
311  
312   for (Int_t i=1;i<tcltokens->GetEntriesFast();i=i+2) {
313     level=((TObjString*)tcltokens->At(i+1))->String()[0];
314     input=atoi((((TObjString*)tcltokens->At(i+1))->String()).Remove(0));
315     Bool_t inpnext = IsInputFired(aEsd,level,input);
316     Char_t op =((TObjString*)tcltokens->At(i))->String()[0];
317     if (op == '&') tcl=tcl && inpnext;
318     else if (op == '|') tcl =tcl || inpnext;
319     else {
320        AliError(Form("Syntax error in %s", tclass));
321        tcltokens->Delete();
322        return kFALSE;
323     }
324   }
325   tcltokens->Delete();
326   return tcl;
327 }
328
329 Bool_t AliTriggerAnalysis::IsInputFired(const AliESDEvent* aEsd, Char_t level, UInt_t input) const
330 {
331   // Checks trigger input of any level
332   
333   switch (level)
334   {
335     case '0': return IsL0InputFired(aEsd,input);
336     case '1': return IsL1InputFired(aEsd,input);
337     case '2': return IsL2InputFired(aEsd,input);
338     default:
339       AliError(Form("Wrong level %i",level));
340       return kFALSE;
341   }
342 }
343
344 Bool_t AliTriggerAnalysis::IsL0InputFired(const AliESDEvent* aEsd, UInt_t input) const 
345 {
346   // Checks if corresponding bit in mask is on
347   
348   UInt_t inpmask = aEsd->GetHeader()->GetL0TriggerInputs();
349   return (inpmask & (1<<(input-1)));
350 }
351
352 Bool_t AliTriggerAnalysis::IsL1InputFired(const AliESDEvent* aEsd, UInt_t input) const
353 {
354   // Checks if corresponding bit in mask is on
355   
356   UInt_t inpmask = aEsd->GetHeader()->GetL1TriggerInputs();
357   return (inpmask & (1<<(input-1)));
358 }
359
360 Bool_t AliTriggerAnalysis::IsL2InputFired(const AliESDEvent* aEsd, UInt_t input) const 
361 {
362   // Checks if corresponding bit in mask is on
363   
364   UInt_t inpmask = aEsd->GetHeader()->GetL2TriggerInputs();
365   return (inpmask & (1<<(input-1)));
366 }
367
368 void AliTriggerAnalysis::FillHistograms(const AliESDEvent* aEsd) 
369 {
370   // fills the histograms with the info from the ESD
371   
372   fHistBitsSPD->Fill(SPDFiredChips(aEsd, 0), SPDFiredChips(aEsd, 1, kTRUE));
373   
374   V0Trigger(aEsd, kASide, kTRUE);
375   V0Trigger(aEsd, kCSide, kTRUE);
376   
377   AliESDZDC* zdcData = aEsd->GetESDZDC();
378   if (zdcData)
379   {
380     UInt_t quality = zdcData->GetESDQuality();
381     
382     // from Nora's presentation, general first physics meeting 16.10.09
383     static UInt_t zpc  = 0x20;
384     static UInt_t znc  = 0x10;
385     static UInt_t zem1 = 0x08;
386     static UInt_t zem2 = 0x04;
387     static UInt_t zpa  = 0x02;
388     static UInt_t zna  = 0x01;
389    
390     fHistZDC->Fill(1, quality & zna);
391     fHistZDC->Fill(2, quality & zpa);
392     fHistZDC->Fill(3, quality & zem2);
393     fHistZDC->Fill(4, quality & zem1);
394     fHistZDC->Fill(5, quality & znc);
395     fHistZDC->Fill(6, quality & zpc);
396   }
397   else
398   {
399     fHistZDC->Fill(-1);
400     AliError("AliESDZDC not available");
401   }
402   
403   fHistFMDA->Fill(FMDHitCombinations(aEsd, kASide, kTRUE));
404   fHistFMDC->Fill(FMDHitCombinations(aEsd, kCSide, kTRUE));
405 }
406   
407 void AliTriggerAnalysis::FillTriggerClasses(const AliESDEvent* aEsd)
408 {
409   // fills trigger classes map
410   
411   TParameter<Long64_t>* count = dynamic_cast<TParameter<Long64_t>*> (fTriggerClasses->GetValue(aEsd->GetFiredTriggerClasses().Data()));
412   if (!count)
413   {
414     count = new TParameter<Long64_t>(aEsd->GetFiredTriggerClasses(), 0);
415     fTriggerClasses->Add(new TObjString(aEsd->GetFiredTriggerClasses().Data()), count);
416   }
417   count->SetVal(count->GetVal() + 1);
418   
419   // TODO add first and last orbit number here
420 }
421
422 Int_t AliTriggerAnalysis::SPDFiredChips(const AliESDEvent* aEsd, Int_t origin, Bool_t fillHist) const
423 {
424   // returns the number of fired chips in the SPD
425   //
426   // origin = 0 --> aEsd->GetMultiplicity()->GetNumberOfFiredChips() (filled from clusters)
427   // origin = 1 --> aEsd->GetMultiplicity()->TestFastOrFiredChips() (from hardware bits)
428   
429   const AliMultiplicity* mult = aEsd->GetMultiplicity();
430   if (!mult)
431   {
432     AliError("AliMultiplicity not available");
433     return -1;
434   }
435   
436   if (origin == 0)
437     return mult->GetNumberOfFiredChips(0) + mult->GetNumberOfFiredChips(1);
438     
439   if (origin == 1)
440   {
441     Int_t nChips = 0;
442     for (Int_t i=0; i<1200; i++)
443       if (mult->TestFastOrFiredChips(i) == kTRUE)
444       {
445         nChips++;
446         if (fillHist)
447           fHistFiredBitsSPD->Fill(i);
448       }
449     return nChips;
450   }
451   
452   return -1;
453 }
454
455 Bool_t AliTriggerAnalysis::SPDGFOTrigger(const AliESDEvent* aEsd, Int_t origin) const
456 {
457   // Returns if the SPD gave a global Fast OR trigger
458   
459   Int_t firedChips = SPDFiredChips(aEsd, origin);
460   
461   if (firedChips >= fSPDGFOThreshold)
462     return kTRUE;
463   return kFALSE;
464 }
465
466 AliTriggerAnalysis::V0Decision AliTriggerAnalysis::V0Trigger(const AliESDEvent* aEsd, AliceSide side, Bool_t fillHists) const
467 {
468   // Returns the V0 trigger decision in V0A | V0C
469   //
470   // Based on algorithm by Cvetan Cheshkov
471   
472   AliESDVZERO* esdV0 = aEsd->GetVZEROData();
473   if (!esdV0)
474   {
475     AliError("AliESDVZERO not available");
476     return kV0Invalid;
477   }
478   
479   Int_t begin = -1;
480   Int_t end = -1;
481   
482   if (side == kASide)
483   {
484     begin = 32;
485     end = 64;
486   } 
487   else if (side == kCSide)
488   {
489     begin = 0;
490     end = 32;
491   }
492   else
493     return kV0Invalid;
494   
495   Float_t time = 0;
496   Int_t ntime = 0;
497   for (Int_t i = begin; i < end; ++i) {
498     if (esdV0->GetTime(i) > 1e-6) {
499       Float_t correctedTime = V0CorrectLeadingTime(i, esdV0->GetTime(i), esdV0->GetAdc(i));
500     
501       time += correctedTime;
502       ntime++;
503     }
504   }
505
506   if (ntime > 0) 
507     time /= ntime;
508
509   if (fillHists)
510   {
511     if (side == kASide && fHistV0A)
512       fHistV0A->Fill(time);
513     if (side == kCSide && fHistV0C)
514       fHistV0C->Fill(time);
515   }
516   
517   if (side == kASide)
518   {
519     if (time > 48 && time < 62)
520       return kV0BB;
521     if (time > 26 && time < 33) 
522       return kV0BG;
523   }
524   
525   if (side == kCSide)
526   {
527     if (time > 49 && time < 60)
528       return kV0BB;
529     if (time > 43 && time < 48.5)
530       return kV0BG; 
531   }
532   
533   return kV0Empty;
534 }
535
536 Float_t AliTriggerAnalysis::V0CorrectLeadingTime(Int_t i, Float_t time, Float_t adc) const
537 {
538   // Correct for slewing and align the channels
539   //
540   // Author: Cvetan Cheshkov
541
542   if (time == 0) return 0;
543
544   // Time offsets between channels
545   Float_t timeShift[64] = {30.2914 , 30.0019 , 30.7429 , 30.1997 , 30.1511 , 29.6437 , 30.0609 , 29.5452 , 30.1437 , 30.745 , 30.7537 , 30.446 , 30.2771 , 30.838 , 30.3748 , 30.0635 , 30.1786 , 30.282 , 31.0992 , 30.7491 , 30.624 , 30.9268 , 30.6585 , 30.4895 , 31.5815 , 31.3871 , 31.2032 , 31.5778 , 31.0838 , 31.2259 , 31.2122 , 31.5989 , 28.3792 , 28.8325 , 27.8719 , 28.3475 , 26.9925 , 27.9300 , 28.4223 , 28.4996 , 28.2934 , 28.1281 , 27.209 , 28.5327 , 28.1181 , 28.0888 , 29.5111 , 28.6601 , 29.7705 , 29.6531 , 30.3373 , 30.2345 , 30.5935 , 29.8164 , 30.2235 , 29.6505 , 30.1225 , 31.2045 , 30.8399 , 30.6789 , 30.2784 , 31.7028 , 31.4239 , 30.1814};
546   time -= timeShift[i];
547
548   // Slewing correction
549   if (adc == 0) return time;
550
551   time += 30.;
552   if (adc > 300.) adc = 300.;
553   if (adc > 70.) {
554     return (time -
555             2.93028e+01 +
556             adc*1.25188e-02 -
557             adc*adc*2.68348e-05);
558   }
559   else {
560     return (time -
561             3.52314e+01 +
562             adc*5.99289e-01 -
563             adc*adc*2.74668e-02 +
564             adc*adc*adc*6.61224e-04 -
565             adc*adc*adc*adc*7.77105e-06 +
566             adc*adc*adc*adc*adc*3.51229e-08);
567   }
568 }
569
570 Bool_t AliTriggerAnalysis::ZDCTrigger(const AliESDEvent* aEsd, AliceSide side) const
571 {
572   // Returns if ZDC triggered
573   
574   AliESDZDC* zdcData = aEsd->GetESDZDC();
575   if (!zdcData)
576   {
577     AliError("AliESDZDC not available");
578     return kFALSE;
579   }
580   
581   UInt_t quality = zdcData->GetESDQuality();
582   
583   // from Nora's presentation, general first physics meeting 16.10.09
584   static UInt_t zpc  = 0x20;
585   static UInt_t znc  = 0x10;
586   static UInt_t zem1 = 0x08;
587   static UInt_t zem2 = 0x04;
588   static UInt_t zpa  = 0x02;
589   static UInt_t zna  = 0x01;
590   
591   if (side == kASide && ((quality & zpa) || (quality & zna)))
592     return kTRUE;
593   if (side == kCentralBarrel && ((quality & zem1) || (quality & zem2)))
594     return kTRUE;
595   if (side == kCSide && ((quality & zpc) || (quality & znc)))
596     return kTRUE;
597   
598   return kFALSE;
599 }
600
601 Int_t AliTriggerAnalysis::FMDHitCombinations(const AliESDEvent* aEsd, AliceSide side, Bool_t fillHistograms) const
602 {
603   // returns number of hit combinations agove threshold
604   //
605   // Authors: FMD team, Hans Dalsgaard (code merged from FMD/AliFMDOfflineTrigger)
606
607   // Workaround for AliESDEvent::GetFMDData is not const!
608   const AliESDFMD* fmdData = (const_cast<AliESDEvent*>(aEsd))->GetFMDData();
609   if (!fmdData)
610   {
611     AliError("AliESDFMD not available");
612     return -1;
613   }
614
615   Int_t detFrom = (side == kASide) ? 1 : 3;
616   Int_t detTo   = (side == kASide) ? 2 : 3;
617
618   Int_t triggers = 0;
619   Float_t totalMult = 0;
620   for (UShort_t det=detFrom;det<=detTo;det++) {
621     Int_t nRings = (det == 1 ? 1 : 2);
622     for (UShort_t ir = 0; ir < nRings; ir++) {
623       Char_t   ring = (ir == 0 ? 'I' : 'O');
624       UShort_t nsec = (ir == 0 ? 20  : 40);
625       UShort_t nstr = (ir == 0 ? 512 : 256);
626       for (UShort_t sec =0; sec < nsec;  sec++) {
627         for (UShort_t strip = 0; strip < nstr; strip++) {
628           Float_t mult = fmdData->Multiplicity(det,ring,sec,strip);
629           if (mult == AliESDFMD::kInvalidMult) continue;
630           
631           if (fillHistograms)
632             fHistFMDSingle->Fill(mult);
633           
634           if (mult > fFMDLowCut)
635             totalMult = totalMult + mult;
636           else
637           {
638             if (totalMult > fFMDHitCut)
639               triggers++;
640               
641             if (fillHistograms)
642               fHistFMDSum->Fill(totalMult);
643               
644             totalMult = 0;
645           }
646         }
647       }
648     }
649   }
650   
651   return triggers;
652 }
653
654 Bool_t AliTriggerAnalysis::FMDTrigger(const AliESDEvent* aEsd, AliceSide side) const
655 {
656   // Returns if the FMD triggered
657   //
658   // Authors: FMD team, Hans Dalsgaard (code merged from FMD/AliFMDOfflineTrigger)
659
660   Int_t triggers = FMDHitCombinations(aEsd, side, kFALSE);
661   
662   if (triggers > 0)
663     return kTRUE;
664     
665   return kFALSE;
666 }
667
668 Long64_t AliTriggerAnalysis::Merge(TCollection* list)
669 {
670   // Merge a list of AliMultiplicityCorrection objects with this (needed for
671   // PROOF).
672   // Returns the number of merged objects (including this).
673
674   if (!list)
675     return 0;
676
677   if (list->IsEmpty())
678     return 1;
679
680   TIterator* iter = list->MakeIterator();
681   TObject* obj;
682
683   // collections of all histograms
684   const Int_t nHists = 9;
685   TList collections[nHists];
686
687   Int_t count = 0;
688   while ((obj = iter->Next())) {
689
690     AliTriggerAnalysis* entry = dynamic_cast<AliTriggerAnalysis*> (obj);
691     if (entry == 0) 
692       continue;
693
694     Int_t n = 0;
695     collections[n++].Add(entry->fHistV0A);
696     collections[n++].Add(entry->fHistV0C);
697     collections[n++].Add(entry->fHistZDC);
698     collections[n++].Add(entry->fHistFMDA);
699     collections[n++].Add(entry->fHistFMDC);
700     collections[n++].Add(entry->fHistFMDSingle);
701     collections[n++].Add(entry->fHistFMDSum);
702     collections[n++].Add(entry->fHistBitsSPD);
703     collections[n++].Add(entry->fHistFiredBitsSPD);
704
705     // merge fTriggerClasses
706     TIterator* iter2 = entry->fTriggerClasses->MakeIterator();
707     TObjString* obj = 0;
708     while ((obj = dynamic_cast<TObjString*> (iter2->Next())))
709     {
710       TParameter<Long64_t>* param2 = dynamic_cast<TParameter<Long64_t>*> (entry->fTriggerClasses->GetValue(obj));
711       
712       TParameter<Long64_t>* param1 = dynamic_cast<TParameter<Long64_t>*> (fTriggerClasses->GetValue(obj));
713       if (param1)
714       {
715         param1->SetVal(param1->GetVal() + param2->GetVal());
716       }
717       else
718       {
719         param1 = dynamic_cast<TParameter<Long64_t>*> (param2->Clone());
720         fTriggerClasses->Add(new TObjString(obj->String()), param1);
721       }
722     }
723     
724     delete iter2;
725   
726     count++;
727   }
728
729   Int_t n = 0;
730   fHistV0A->Merge(&collections[n++]);
731   fHistV0C->Merge(&collections[n++]);
732   fHistZDC->Merge(&collections[n++]);
733   fHistFMDA->Merge(&collections[n++]);
734   fHistFMDC->Merge(&collections[n++]);
735   fHistFMDSingle->Merge(&collections[n++]);
736   fHistFMDSum->Merge(&collections[n++]);
737   fHistBitsSPD->Merge(&collections[n++]);
738   fHistFiredBitsSPD->Merge(&collections[n++]);
739   
740   delete iter;
741
742   return count+1;
743 }
744
745 void AliTriggerAnalysis::WriteHistograms() const
746 {
747   // write histograms to current directory
748   
749   if (!fHistBitsSPD)
750     return;
751     
752   fHistBitsSPD->Write();
753   fHistBitsSPD->ProjectionX()->Write();
754   fHistBitsSPD->ProjectionY()->Write();
755   fHistFiredBitsSPD->Write();
756   fHistV0A->Write();
757   fHistV0C->Write();
758   fHistZDC->Write();
759   fHistFMDA->Write();
760   fHistFMDC->Write();
761   fHistFMDSingle->Write();
762   fHistFMDSum->Write();
763   
764   fTriggerClasses->Write("fTriggerClasses", TObject::kSingleKey);
765 }
766
767 void AliTriggerAnalysis::PrintTriggerClasses() const
768 {
769   // print trigger classes
770   
771   Printf("Trigger Classes:");
772   
773   TIterator* iter = fTriggerClasses->MakeIterator();
774   TObjString* obj = 0;
775   while ((obj = dynamic_cast<TObjString*> (iter->Next())))
776   {
777     TParameter<Long64_t>* param = dynamic_cast<TParameter<Long64_t>*> (fTriggerClasses->GetValue(obj));
778     
779     Printf("%s: %ld triggers", obj->String().Data(), param->GetVal());
780   }
781 }