1 /* $Id: AliTriggerAnalysis.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 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 //-------------------------------------------------------------------------
25 #include <Riostream.h>
29 #include <TIterator.h>
30 #include "TParameter.h"
34 #include <AliTriggerAnalysis.h>
37 #include <AliAnalysisManager.h>
39 #include <AliESDEvent.h>
41 #include <AliMultiplicity.h>
42 #include <AliESDVZERO.h>
43 #include <AliESDZDC.h>
44 #include <AliESDFMD.h>
45 #include <AliESDVertex.h>
46 #include <AliESDtrackCuts.h>
48 #include <AliESDTrdTrack.h>
50 ClassImp(AliTriggerAnalysis)
52 AliTriggerAnalysis::AliTriggerAnalysis() :
60 fZDCCutRefSum(-568.5),
61 fZDCCutRefDelta(-2.1),
62 fZDCCutSigmaSum(3.25),
63 fZDCCutSigmaDelta(2.25),
64 fZDCCutRefSumCorr(-65.5),
65 fZDCCutRefDeltaCorr(-2.1),
66 fZDCCutSigmaSumCorr(6.0),
67 fZDCCutSigmaDeltaCorr(1.2),
68 fZDCCutZNATimeCorrMin(0.0),
69 fZDCCutZNATimeCorrMax(2.0),
70 fZDCCutZNCTimeCorrMin(0.0),
71 fZDCCutZNCTimeCorrMax(5.0),
109 AliTriggerAnalysis::~AliTriggerAnalysis()
119 if (fHistFiredBitsSPD)
121 delete fHistFiredBitsSPD;
122 fHistFiredBitsSPD = 0;
125 if (fHistSPDClsVsTrk)
127 delete fHistSPDClsVsTrk;
128 fHistSPDClsVsTrk = 0;
158 if (fHistTimeCorrZDC)
160 delete fHistTimeCorrZDC;
161 fHistTimeCorrZDC = 0;
178 delete fHistFMDSingle;
196 fTriggerClasses->DeleteAll();
197 delete fTriggerClasses;
202 delete fEsdTrackCuts;
207 void AliTriggerAnalysis::EnableHistograms(Bool_t isLowFlux)
209 // creates the monitoring histograms (dynamical range of histograms can be adapted for pp and pPb via isLowFlux flag)
211 // do not add this hists to the directory
212 Bool_t oldStatus = TH1::AddDirectoryStatus();
213 TH1::AddDirectory(kFALSE);
215 Int_t nBins = isLowFlux ? 600 : 1202;
216 fHistBitsSPD = new TH2F("fHistBitsSPD", "SPD GFO;number of fired chips (offline);number of fired chips (hardware)", nBins, -1.5, -1.5 + nBins, nBins, -1.5, -1.5+nBins);
218 fHistFiredBitsSPD = new TH1F("fHistFiredBitsSPD", "SPD GFO Hardware;chip number;events", 1200, -0.5, 1199.5);
220 Int_t nBinsX = isLowFlux ? 100 : 300;
221 Int_t nBinsY = isLowFlux ? 500 : 1000;
222 Float_t xMax = isLowFlux ? 400 : 2999.5;
223 Float_t yMax = isLowFlux ? 4000 : 9999.5;
224 fHistSPDClsVsTrk = new TH2F("fHistSPDClsVsTrk", "SPD Clusters vs Tracklets", nBinsX, -0.5, xMax, nBinsY, -0.5, yMax);
226 fHistV0A = new TH1F("fHistV0A", "V0A;leading time (ns);events", 400, -100, 100);
227 fHistV0C = new TH1F("fHistV0C", "V0C;leading time (ns);events", 400, -100, 100);
228 fHistZDC = new TH1F("fHistZDC", "ZDC;trigger bits;events", 8, -1.5, 6.5);
229 fHistTDCZDC = new TH1F("fHistTDCZDC", "ZDC;TDC bits;events", 32, -0.5, 32-0.5);
230 fHistTimeZDC = new TH2F("fHistTimeZDC", "ZDC;TDC timing A+C vs C-A; events", 120,-30,30,120,-600,-540);
231 fHistTimeCorrZDC = new TH2F("fHistTimeCorrZDC", "ZDC;Corrected TDC timing A+C vs C-A; events", 120,-30,30,260,-100,30);
234 fHistFMDA = new TH1F("fHistFMDA", "FMDA;combinations above threshold;events", 102, -1.5, 100.5);
235 fHistFMDC = new TH1F("fHistFMDC", "FMDC;combinations above threshold;events", 102, -1.5, 100.5);
236 fHistFMDSingle = new TH1F("fHistFMDSingle", "FMD single;multiplicity value;counts", 1000, 0, 10);
237 fHistFMDSum = new TH1F("fHistFMDSum", "FMD sum;multiplicity value;counts", 1000, 0, 10);
238 fHistT0 = new TH1F("fHistT0", "T0;time (ns);events", 100, -25, 25);
240 fTriggerClasses = new TMap;
241 fTriggerClasses->SetOwner();
243 TH1::AddDirectory(oldStatus);
246 //____________________________________________________________________
247 const char* AliTriggerAnalysis::GetTriggerName(Trigger trigger)
249 // returns the name of the requested trigger
250 // the returned string will only be valid until the next call to this function [not thread-safe]
254 UInt_t triggerNoFlags = (UInt_t) trigger % (UInt_t) kStartOfFlags;
256 switch (triggerNoFlags)
258 case kAcceptAll : str = "ACCEPT ALL (bypass!)"; break;
259 case kMB1 : str = "MB1"; break;
260 case kMB2 : str = "MB2"; break;
261 case kMB3 : str = "MB3"; break;
262 case kSPDGFO : str = "SPD GFO"; break;
263 case kSPDGFOBits : str = "SPD GFO Bits"; break;
264 case kSPDGFOL0 : str = "SPD GFO L0 (first layer)"; break;
265 case kSPDGFOL1 : str = "SPD GFO L1 (second layer)"; break;
266 case kSPDClsVsTrkBG : str = "Cluster vs Tracklets"; break;
267 case kV0A : str = "V0 A BB"; break;
268 case kV0C : str = "V0 C BB"; break;
269 case kV0OR : str = "V0 OR BB"; break;
270 case kV0AND : str = "V0 AND BB"; break;
271 case kV0ABG : str = "V0 A BG"; break;
272 case kV0CBG : str = "V0 C BG"; break;
273 case kZDC : str = "ZDC"; break;
274 case kZDCA : str = "ZDC A"; break;
275 case kZDCC : str = "ZDC C"; break;
276 case kZNA : str = "ZN A"; break;
277 case kZNC : str = "ZN C"; break;
278 case kZNABG : str = "ZN A BG"; break;
279 case kZNCBG : str = "ZN C BG"; break;
280 case kFMDA : str = "FMD A"; break;
281 case kFMDC : str = "FMD C"; break;
282 case kFPANY : str = "SPD GFO | V0 | ZDC | FMD"; break;
283 case kNSD1 : str = "NSD1"; break;
284 case kMB1Prime: str = "MB1prime"; break;
285 case kZDCTDCA : str = "ZDC TDC A"; break;
286 case kZDCTDCC : str = "ZDC TDC C"; break;
287 case kZDCTime : str = "ZDC Time Cut"; break;
288 case kCentral : str = "V0 Central"; break;
289 case kSemiCentral : str = "V0 Semi-central"; break;
290 case kEMCAL : str = "EMCAL"; break;
291 case kTRDHCO : str = "TRD cosmics"; break;
292 case kTRDHJT : str = "TRD Jet"; break;
293 case kTRDHSE : str = "TRD Single Electron"; break;
294 case kTRDHQU : str = "TRD Quarkonia"; break;
295 case kTRDHEE : str = "TRD Dielectron"; break;
296 default: str = ""; break;
299 if (trigger & kOfflineFlag)
302 if (trigger & kOneParticle)
303 str += " OneParticle";
305 if (trigger & kOneTrack)
311 Bool_t AliTriggerAnalysis::IsTriggerFired(const AliESDEvent* aEsd, Trigger trigger)
313 // checks if an event has been triggered
315 if (trigger & kOfflineFlag)
316 return IsOfflineTriggerFired(aEsd, trigger);
318 return IsTriggerBitFired(aEsd, trigger);
321 Bool_t AliTriggerAnalysis::IsTriggerBitFired(const AliESDEvent* aEsd, Trigger trigger) const
323 // checks if an event is fired using the trigger bits
325 return IsTriggerBitFired(aEsd->GetTriggerMask(), trigger);
328 Bool_t AliTriggerAnalysis::IsTriggerBitFired(ULong64_t triggerMask, Trigger trigger) const
330 // checks if an event is fired using the trigger bits
332 // this function needs the branch TriggerMask in the ESD
334 // definitions from p-p.cfg
335 ULong64_t spdFO = (1 << 14);
336 ULong64_t v0left = (1 << 10);
337 ULong64_t v0right = (1 << 11);
348 if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
354 if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
360 if (triggerMask & spdFO && (triggerMask & v0left) && (triggerMask & v0right))
366 if (triggerMask & spdFO)
371 Printf("IsEventTriggered: ERROR: Trigger type %d not implemented in this method", (Int_t) trigger);
378 Bool_t AliTriggerAnalysis::IsTriggerBitFired(const AliESDEvent* aEsd, ULong64_t tclass) const
380 // Checks if corresponding bit in mask is on
382 ULong64_t trigmask = aEsd->GetTriggerMask();
383 return (trigmask & (1ull << (tclass-1)));
386 Int_t AliTriggerAnalysis::EvaluateTrigger(const AliESDEvent* aEsd, Trigger trigger)
388 // evaluates a given trigger "offline"
389 // trigger combinations are not supported, for that see IsOfflineTriggerFired
391 UInt_t triggerNoFlags = (UInt_t) trigger % (UInt_t) kStartOfFlags;
392 Bool_t offline = kFALSE;
393 if (trigger & kOfflineFlag)
397 switch (triggerNoFlags)
401 decision = SPDFiredChips(aEsd, (offline) ? 0 : 1);
406 decision = SPDFiredChips(aEsd, (offline) ? 0 : 1, kFALSE, 1);
411 decision = SPDFiredChips(aEsd, (offline) ? 0 : 1, kFALSE, 2);
417 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
418 decision = IsSPDClusterVsTrackletBG(aEsd);
423 if (V0Trigger(aEsd, kASide, !offline) == kV0BB)
429 if (V0Trigger(aEsd, kCSide, !offline) == kV0BB)
435 if (V0Trigger(aEsd, kASide, !offline) == kV0BG)
441 if (V0Trigger(aEsd, kCSide, !offline) == kV0BG)
447 if (T0Trigger(aEsd, !offline) == kT0BB)
454 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
455 if (T0Trigger(aEsd, !offline) == kT0DecBG)
462 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
463 if (T0Trigger(aEsd, !offline) == kT0DecPileup)
470 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
471 if (ZDCTrigger(aEsd, kASide))
478 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
479 if (ZDCTrigger(aEsd, kCSide))
486 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
487 if (ZDCTDCTrigger(aEsd, kASide))
494 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
495 if (ZDCTDCTrigger(aEsd, kCSide))
502 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
503 if (ZDCTimeTrigger(aEsd))
510 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
511 if (ZDCTDCTrigger(aEsd,kASide,kTRUE,kFALSE,kFALSE))
518 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
519 if (ZDCTDCTrigger(aEsd,kCSide,kTRUE,kFALSE,kFALSE))
526 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
527 if (ZDCTimeBGTrigger(aEsd,kASide))
534 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
535 if (ZDCTimeBGTrigger(aEsd,kCSide))
542 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
543 if (FMDTrigger(aEsd, kASide))
550 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
551 if (FMDTrigger(aEsd, kCSide))
558 AliFatal(Form("Offline trigger not available for trigger %d", triggerNoFlags));
559 if (IsL0InputFired(aEsd, 2))
566 AliFatal(Form("Offline trigger not available for trigger %d", triggerNoFlags));
567 if (IsL0InputFired(aEsd, 3))
571 case kTPCLaserWarmUp:
574 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
575 return IsLaserWarmUpTPCEvent(aEsd);
580 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
581 return IsHVdipTPCEvent(aEsd);
586 AliFatal(Form("Offline trigger not available for trigger %d - use centrality selection", triggerNoFlags));
587 if (aEsd->GetVZEROData()) {
588 if (aEsd->GetVZEROData()->TestBit(AliESDVZERO::kTriggerChargeBitsFilled)) {
589 if (aEsd->GetVZEROData()->GetTriggerBits() & (1<<AliESDVZERO::kCTA2andCTC2))
593 AliWarning("V0 centrality trigger bits were not filled!");
600 AliFatal(Form("Offline trigger not available for trigger %d - use centrality selection", triggerNoFlags));
601 if (aEsd->GetVZEROData()) {
602 if (aEsd->GetVZEROData()->TestBit(AliESDVZERO::kTriggerChargeBitsFilled)) {
603 if (aEsd->GetVZEROData()->GetTriggerBits() & (1<<AliESDVZERO::kCTA1andCTC1))
607 AliWarning("V0 centrality trigger bits were not filled!");
614 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
615 if (EMCALCellsTrigger(aEsd))
622 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
623 if(TRDTrigger(aEsd,kTRDHCO))
630 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
631 if(TRDTrigger(aEsd,kTRDHJT))
638 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
639 if(TRDTrigger(aEsd,kTRDHSE))
646 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
647 if(TRDTrigger(aEsd,kTRDHQU))
654 AliFatal(Form("Online trigger not available for trigger %d", triggerNoFlags));
655 if(TRDTrigger(aEsd,kTRDHEE))
661 AliFatal(Form("Trigger type %d not implemented", triggerNoFlags));
668 Bool_t AliTriggerAnalysis::IsLaserWarmUpTPCEvent(const AliESDEvent* esd)
671 // This function flags noisy TPC events which can happen during laser warm-up.
674 Int_t trackCounter = 0;
675 for (Int_t i=0; i<esd->GetNumberOfTracks(); ++i)
677 AliESDtrack *track = esd->GetTrack(i);
681 if (track->GetTPCNcls() < 30) continue;
682 if (TMath::Abs(track->Eta()) > 0.005) continue;
683 if (track->Pt() < 4) continue;
684 if (track->GetKinkIndex(0) > 0) continue;
686 UInt_t status = track->GetStatus();
687 if ((status&AliESDtrack::kITSrefit)==AliESDtrack::kITSrefit) continue; // explicitly ask for tracks without ITS refit
688 if ((status&AliESDtrack::kTPCrefit)!=AliESDtrack::kTPCrefit) continue;
690 if (track->GetTPCsignal() > 10) continue; // explicitly ask for tracks without dE/dx
694 if (trackCounter > 15)
699 Bool_t AliTriggerAnalysis::IsHVdipTPCEvent(const AliESDEvent* esd) {
701 // This function flags events in which the TPC chamber HV is not at its nominal value
703 if (!esd->IsDetectorOn(AliDAQ::kTPC)) return kTRUE;
709 Bool_t AliTriggerAnalysis::IsOfflineTriggerFired(const AliESDEvent* aEsd, Trigger trigger)
711 // checks if an event has been triggered "offline"
713 UInt_t triggerNoFlags = (UInt_t) trigger % (UInt_t) kStartOfFlags;
715 Bool_t decision = kFALSE;
716 switch (triggerNoFlags)
725 if (SPDGFOTrigger(aEsd, 0) || V0Trigger(aEsd, kASide, kFALSE) == kV0BB || V0Trigger(aEsd, kCSide, kFALSE) == kV0BB)
731 if (SPDGFOTrigger(aEsd, 0) && (V0Trigger(aEsd, kASide, kFALSE) == kV0BB || V0Trigger(aEsd, kCSide, kFALSE) == kV0BB))
737 if (SPDGFOTrigger(aEsd, 0) && V0Trigger(aEsd, kASide, kFALSE) == kV0BB && V0Trigger(aEsd, kCSide, kFALSE) == kV0BB)
743 if (SPDGFOTrigger(aEsd, 0))
749 if (SPDGFOTrigger(aEsd, 1))
755 if (V0Trigger(aEsd, kASide, kFALSE) == kV0BB)
761 if (V0Trigger(aEsd, kCSide, kFALSE) == kV0BB)
767 if (V0Trigger(aEsd, kASide, kFALSE) == kV0BB || V0Trigger(aEsd, kCSide, kFALSE) == kV0BB)
773 if (V0Trigger(aEsd, kASide, kFALSE) == kV0BB && V0Trigger(aEsd, kCSide, kFALSE) == kV0BB)
779 if (V0Trigger(aEsd, kASide, kFALSE) == kV0BG)
785 if (V0Trigger(aEsd, kCSide, kFALSE) == kV0BG)
791 if (ZDCTrigger(aEsd, kASide) || ZDCTrigger(aEsd, kCentralBarrel) || ZDCTrigger(aEsd, kCSide))
797 if (ZDCTrigger(aEsd, kASide))
803 if (ZDCTrigger(aEsd, kCSide))
809 if (ZDCTDCTrigger(aEsd,kASide,kTRUE,kFALSE,kFALSE))
815 if (ZDCTDCTrigger(aEsd,kCSide,kTRUE,kFALSE,kFALSE))
821 if (ZDCTimeBGTrigger(aEsd,kASide))
827 if (ZDCTimeBGTrigger(aEsd,kCSide))
833 if (FMDTrigger(aEsd, kASide))
839 if (FMDTrigger(aEsd, kCSide))
845 if (SPDGFOTrigger(aEsd, 0) || V0Trigger(aEsd, kASide, kFALSE) == kV0BB || V0Trigger(aEsd, kCSide, kFALSE) == kV0BB || ZDCTrigger(aEsd, kASide) || ZDCTrigger(aEsd, kCentralBarrel) || ZDCTrigger(aEsd, kCSide) || FMDTrigger(aEsd, kASide) || FMDTrigger(aEsd, kCSide))
851 if (SPDFiredChips(aEsd, 0) >= 5 || (V0Trigger(aEsd, kASide, kFALSE) == kV0BB && V0Trigger(aEsd, kCSide, kFALSE) == kV0BB))
858 if (SPDGFOTrigger(aEsd, 0))
860 if (V0Trigger(aEsd, kASide, kFALSE) == kV0BB)
862 if (V0Trigger(aEsd, kCSide, kFALSE) == kV0BB)
872 if (EMCALCellsTrigger(aEsd))
878 if(TRDTrigger(aEsd,kTRDHCO))
884 if(TRDTrigger(aEsd,kTRDHJT))
890 if(TRDTrigger(aEsd,kTRDHSE))
896 if(TRDTrigger(aEsd,kTRDHQU))
902 if(TRDTrigger(aEsd,kTRDHEE))
908 AliFatal(Form("Trigger type %d not implemented", triggerNoFlags));
912 // hadron-level requirement
913 if (decision && (trigger & kOneParticle))
917 const AliESDVertex* vertex = aEsd->GetPrimaryVertexSPD();
918 const AliMultiplicity* mult = aEsd->GetMultiplicity();
920 if (mult && vertex && vertex->GetNContributors() > 0 && (!vertex->IsFromVertexerZ() || vertex->GetDispersion() < 0.02) && TMath::Abs(vertex->GetZv()) < 5.5)
922 for (Int_t i=0; i<mult->GetNumberOfTracklets(); ++i)
924 if (TMath::Abs(mult->GetEta(i)) < 1)
933 // hadron level definition for TPC tracks
935 if (decision && (trigger & kOneTrack))
938 const AliESDVertex* vertex =0x0;
939 vertex = aEsd->GetPrimaryVertexTracks();
940 if (!vertex || vertex->GetNContributors() <= 0)
942 vertex = aEsd->GetPrimaryVertexSPD();
944 Float_t ptmin, ptmax;
945 fEsdTrackCuts->GetPtRange(ptmin,ptmax);
946 AliDebug(3, Form("ptmin = %f, ptmax = %f\n",ptmin, ptmax));
948 if (vertex && vertex->GetNContributors() > 0 && (!vertex->IsFromVertexerZ() || vertex->GetDispersion() < 0.02) && TMath::Abs(vertex->GetZv()) < 10.) {
949 AliDebug(3,Form("Check on the vertex passed\n"));
950 for (Int_t i=0; i<aEsd->GetNumberOfTracks(); ++i){
951 if (fTPCOnly == kFALSE){
952 if (fEsdTrackCuts->AcceptTrack(aEsd->GetTrack(i))){
953 AliDebug(2, Form("pt of track = %f --> check passed\n",aEsd->GetTrack(i)->Pt()));
960 AliESDtrack *tpcTrack = fEsdTrackCuts->GetTPCOnlyTrack((AliESDEvent*)aEsd, i);
962 AliDebug(3,Form("track %d is NOT a TPC track",i));
966 AliDebug(3,Form("track %d IS a TPC track",i));
967 if (!(fEsdTrackCuts->AcceptTrack(tpcTrack))) {
968 AliDebug(2, Form("TPC track %d NOT ACCEPTED, pt = %f, eta = %f",i,tpcTrack->Pt(), tpcTrack->Eta()));
969 delete tpcTrack; tpcTrack = 0x0;
971 } // end if the TPC track is not accepted
973 AliDebug(2, Form("TPC track %d ACCEPTED, pt = %f, eta = %f",i,tpcTrack->Pt(), tpcTrack->Eta()));
975 delete tpcTrack; tpcTrack = 0x0;
977 } // end if the TPC track is accepted
978 } // end if it is a TPC track
979 } // end if you are looking at TPC only tracks
980 } // end loop on tracks
981 } // end check on vertex
983 AliDebug(4,Form("Check on the vertex not passed\n"));
984 for (Int_t i=0; i<aEsd->GetNumberOfTracks(); ++i){
985 if (fEsdTrackCuts->AcceptTrack(aEsd->GetTrack(i))){
986 AliDebug(4,Form("pt of track = %f --> check would be passed if the vertex was ok\n",aEsd->GetTrack(i)->Pt()));
991 if (!decision) AliDebug(3,("Check for kOneTrack NOT passed\n"));
997 Bool_t AliTriggerAnalysis::IsTriggerClassFired(const AliESDEvent* aEsd, const Char_t* tclass) const
999 // tclass is logical function of inputs, e.g. 01 && 02 || 03 && 11 && 21
1000 // = L0 inp 1 && L0 inp 2 || L0 inp 3 && L1 inp 1 && L2 inp 1
1001 // NO brackets in logical function !
1002 // Spaces between operators and inputs.
1003 // Not all logical functions are available in CTP=
1004 // =any function of first 4 inputs; 'AND' of other inputs, check not done
1005 // This method will be replaced/complemened by similar one
1006 // which works withh class and inputs names as in CTP cfg file
1008 TString TClass(tclass);
1009 TObjArray* tcltokens = TClass.Tokenize(" ");
1010 Char_t level=((TObjString*)tcltokens->At(0))->String()[0];
1011 UInt_t input=atoi((((TObjString*)tcltokens->At(0))->String()).Remove(0));
1012 Bool_t tcl = IsInputFired(aEsd,level,input);
1014 for (Int_t i=1;i<tcltokens->GetEntriesFast();i=i+2) {
1015 level=((TObjString*)tcltokens->At(i+1))->String()[0];
1016 input=atoi((((TObjString*)tcltokens->At(i+1))->String()).Remove(0));
1017 Bool_t inpnext = IsInputFired(aEsd,level,input);
1018 Char_t op =((TObjString*)tcltokens->At(i))->String()[0];
1019 if (op == '&') tcl=tcl && inpnext;
1020 else if (op == '|') tcl =tcl || inpnext;
1022 AliError(Form("Syntax error in %s", tclass));
1025 // tcltokens->Delete();
1031 // tcltokens->Delete();
1035 Bool_t AliTriggerAnalysis::IsInputFired(const AliESDEvent* aEsd, Char_t level, UInt_t input) const
1037 // Checks trigger input of any level
1041 case '0': return IsL0InputFired(aEsd,input);
1042 case '1': return IsL1InputFired(aEsd,input);
1043 case '2': return IsL2InputFired(aEsd,input);
1045 AliError(Form("Wrong level %i",level));
1050 Bool_t AliTriggerAnalysis::IsL0InputFired(const AliESDEvent* aEsd, UInt_t input) const
1052 // Checks if corresponding bit in mask is on
1054 UInt_t inpmask = aEsd->GetHeader()->GetL0TriggerInputs();
1055 return (inpmask & (1<<(input-1)));
1058 Bool_t AliTriggerAnalysis::IsL1InputFired(const AliESDEvent* aEsd, UInt_t input) const
1060 // Checks if corresponding bit in mask is on
1062 UInt_t inpmask = aEsd->GetHeader()->GetL1TriggerInputs();
1063 return (inpmask & (1<<(input-1)));
1066 Bool_t AliTriggerAnalysis::IsL2InputFired(const AliESDEvent* aEsd, UInt_t input) const
1068 // Checks if corresponding bit in mask is on
1070 UInt_t inpmask = aEsd->GetHeader()->GetL2TriggerInputs();
1071 return (inpmask & (1<<(input-1)));
1074 void AliTriggerAnalysis::FillHistograms(const AliESDEvent* aEsd)
1076 // fills the histograms with the info from the ESD
1078 fHistBitsSPD->Fill(SPDFiredChips(aEsd, 0), SPDFiredChips(aEsd, 1, kTRUE));
1080 V0Trigger(aEsd, kASide, kFALSE, kTRUE);
1081 V0Trigger(aEsd, kCSide, kFALSE, kTRUE);
1082 T0Trigger(aEsd, kFALSE, kTRUE);
1083 ZDCTDCTrigger(aEsd,kASide,kFALSE,kFALSE,kTRUE);
1084 ZDCTimeTrigger(aEsd,kTRUE);
1085 IsSPDClusterVsTrackletBG(aEsd, kTRUE);
1087 AliESDZDC* zdcData = aEsd->GetESDZDC();
1090 UInt_t quality = zdcData->GetESDQuality();
1092 // from Nora's presentation, general first physics meeting 16.10.09
1093 static UInt_t zpc = 0x20;
1094 static UInt_t znc = 0x10;
1095 static UInt_t zem1 = 0x08;
1096 static UInt_t zem2 = 0x04;
1097 static UInt_t zpa = 0x02;
1098 static UInt_t zna = 0x01;
1100 fHistZDC->Fill(1, (quality & zna) ? 1 : 0);
1101 fHistZDC->Fill(2, (quality & zpa) ? 1 : 0);
1102 fHistZDC->Fill(3, (quality & zem2) ? 1 : 0);
1103 fHistZDC->Fill(4, (quality & zem1) ? 1 : 0);
1104 fHistZDC->Fill(5, (quality & znc) ? 1 : 0);
1105 fHistZDC->Fill(6, (quality & zpc) ? 1 : 0);
1110 AliError("AliESDZDC not available");
1114 fHistFMDA->Fill(FMDHitCombinations(aEsd, kASide, kTRUE));
1115 fHistFMDC->Fill(FMDHitCombinations(aEsd, kCSide, kTRUE));
1119 void AliTriggerAnalysis::FillTriggerClasses(const AliESDEvent* aEsd)
1121 // fills trigger classes map
1123 TParameter<Long64_t>* count = dynamic_cast<TParameter<Long64_t>*> (fTriggerClasses->GetValue(aEsd->GetFiredTriggerClasses().Data()));
1126 count = new TParameter<Long64_t>(aEsd->GetFiredTriggerClasses(), 0);
1127 fTriggerClasses->Add(new TObjString(aEsd->GetFiredTriggerClasses().Data()), count);
1129 count->SetVal(count->GetVal() + 1);
1132 Int_t AliTriggerAnalysis::SSDClusters(const AliESDEvent* aEsd)
1134 // returns the number of clusters in the SSD
1135 const AliMultiplicity* mult = aEsd->GetMultiplicity();
1136 Int_t clusters = mult->GetNumberOfITSClusters(4)+mult->GetNumberOfITSClusters(5);
1141 Int_t AliTriggerAnalysis::SPDFiredChips(const AliESDEvent* aEsd, Int_t origin, Bool_t fillHists, Int_t layer)
1143 // returns the number of fired chips in the SPD
1145 // origin = 0 --> aEsd->GetMultiplicity()->GetNumberOfFiredChips() (filled from clusters)
1146 // origin = 1 --> aEsd->GetMultiplicity()->TestFastOrFiredChips() (from hardware bits)
1147 // layer = 0 --> both layers
1148 // layer = 1 --> inner
1149 // layer = 2 --> outer
1151 const AliMultiplicity* mult = aEsd->GetMultiplicity();
1154 AliError("AliMultiplicity not available");
1160 return mult->GetNumberOfFiredChips(0) + mult->GetNumberOfFiredChips(1);
1162 return mult->GetNumberOfFiredChips(layer-1);
1168 Int_t firstChip = 0;
1169 Int_t lastChip = 1200;
1175 for (Int_t i=firstChip; i<lastChip; i++)
1177 if (mult->TestFastOrFiredChips(i) == kTRUE)
1179 // efficiency simulation (if enabled)
1180 if (fSPDGFOEfficiency)
1182 if (gRandom->Uniform() > fSPDGFOEfficiency->GetBinContent(i+1))
1188 fHistFiredBitsSPD->Fill(i);
1197 Bool_t AliTriggerAnalysis::SPDGFOTrigger(const AliESDEvent* aEsd, Int_t origin)
1199 // Returns if the SPD gave a global Fast OR trigger
1201 Int_t firedChips = SPDFiredChips(aEsd, origin);
1203 if (firedChips >= fSPDGFOThreshold)
1208 Bool_t AliTriggerAnalysis::IsSPDClusterVsTrackletBG(const AliESDEvent* aEsd, Bool_t fillHists){
1209 //rejects BG based on the cluster vs tracklet correlation
1210 // returns true if the event is BG
1211 const AliMultiplicity* mult = aEsd->GetMultiplicity();
1213 AliFatal("No multiplicity object"); // TODO: Should this be fatal?
1215 Int_t ntracklet = mult->GetNumberOfTracklets();
1217 Int_t spdClusters = 0;
1218 for(Int_t ilayer = 0; ilayer < 2; ilayer++){
1219 spdClusters += mult->GetNumberOfITSClusters(ilayer);
1223 fHistSPDClsVsTrk->Fill(ntracklet,spdClusters);
1226 Bool_t isCvsTOk = kFALSE;
1227 Float_t limit = Float_t(fASPDCvsTCut) + Float_t(ntracklet) * fBSPDCvsTCut;
1228 if (spdClusters > limit) isCvsTOk = kTRUE;
1229 else isCvsTOk = kFALSE ;
1236 AliTriggerAnalysis::V0Decision AliTriggerAnalysis::V0Trigger(const AliESDEvent* aEsd, AliceSide side, Bool_t online, Bool_t fillHists)
1238 // Returns the V0 trigger decision in V0A | V0C
1240 // Returns kV0Fake if the calculated average time is in a window where neither BB nor BG is expected.
1241 // The rate of such triggers can be used to estimate the background. Note that the rate has to be
1242 // rescaled with the size of the windows (numerical values see below in the code)
1244 // argument 'online' is used as a switch between online and offline trigger algorithms
1246 // Based on an algorithm by Cvetan Cheshkov
1248 AliESDVZERO* esdV0 = aEsd->GetVZEROData();
1251 AliError("AliESDVZERO not available");
1254 AliDebug(2,Form("In V0Trigger: %f %f",esdV0->GetV0ATime(),esdV0->GetV0CTime()));
1264 else if (side == kCSide)
1272 if (esdV0->TestBit(AliESDVZERO::kDecisionFilled)) {
1274 if (esdV0->TestBit(AliESDVZERO::kOnlineBitsFilled)) {
1275 for (Int_t i = begin; i < end; ++i) {
1276 if (esdV0->GetBBFlag(i)) return kV0BB;
1278 for (Int_t i = begin; i < end; ++i) {
1279 if (esdV0->GetBGFlag(i)) return kV0BG;
1284 AliWarning("V0 online trigger analysis is not yet available!");
1291 if (side == kASide && fHistV0A)
1292 fHistV0A->Fill(esdV0->GetV0ATime());
1293 if (side == kCSide && fHistV0C)
1294 fHistV0C->Fill(esdV0->GetV0CTime());
1297 if (side == kASide) return (V0Decision)esdV0->GetV0ADecision();
1298 else if (side == kCSide) return (V0Decision)esdV0->GetV0CDecision();
1299 else return kV0Invalid;
1308 if (aEsd->GetRunNumber() <= 104803) runRange = 0;
1309 else if (aEsd->GetRunNumber() <= 104876) runRange = 1;
1312 Float_t factors[3][64] = {
1313 // runs: 104792-104803
1314 {4.6,5.9,6.3,6.0,4.7,5.9,4.9,5.4,4.8,4.1,4.9,4.6,4.5,5.5,5.1,5.8,4.3,4.0,4.0,3.3,3.1,2.9,3.0,5.6,3.3,4.9,3.9,5.3,4.1,4.4,3.9,5.5,5.7,9.5,5.1,5.3,6.6,7.1,8.9,4.4,4.1,5.9,9.0,4.5,4.1,6.0,4.7,7.1,4.2,4.7,3.9,6.3,5.9,4.8,4.7,4.5,4.7,5.4,5.8,5.0,5.1,5.9,5.3,3.6},
1315 // runs: 104841-104876
1316 {4.6,4.8,4.9,4.8,4.3,4.9,4.4,4.5,4.6,5.0,4.7,4.6,4.7,4.6,4.6,5.5,4.7,4.5,4.7,5.0,6.5,7.6,5.3,4.9,5.5,4.8,4.6,4.9,4.5,4.5,4.6,4.9,5.7,9.8,4.9,5.2,7.1,7.1,8.1,4.4,4.0,6.0,8.3,4.6,4.2,5.6,4.6,6.4,4.4,4.7,4.5,6.5,6.0,4.7,4.5,4.4,4.8,5.5,5.9,5.3,5.0,5.7,5.1,3.6},
1318 {4.7,5.2,4.8,5.0,4.4,5.0,4.4,4.6,4.6,4.5,4.4,4.6,4.5,4.6,4.8,5.5,4.8,4.5,4.4,4.3,5.4,7.7,5.6,5.0,5.4,4.3,4.5,4.8,4.5,4.5,4.6,5.3,5.7,9.6,4.9,5.4,6.1,7.2,8.6,4.4,4.0,5.4,8.8,4.4,4.2,5.8,4.7,6.7,4.3,4.7,4.0,6.1,6.0,4.9,4.8,4.6,4.7,5.2,5.7,5.0,5.0,5.8,5.3,3.6}
1320 Float_t dA = 77.4 - 11.0;
1321 Float_t dC = 77.4 - 2.9;
1322 // Time misalignment
1323 Float_t timeShift[64] = {0.477957 , 0.0889999 , 0.757669 , 0.205439 , 0.239666 , -0.183705 , 0.442873 , -0.281366 , 0.260976 , 0.788995 , 0.974758 , 0.548532 , 0.495023 , 0.868472 , 0.661167 , 0.358307 , 0.221243 , 0.530179 , 1.26696 , 1.33082 , 1.27086 , 1.77133 , 1.10253 , 0.634806 , 2.14838 , 1.50212 , 1.59253 , 1.66122 , 1.16957 , 1.52056 , 1.47791 , 1.81905 , -1.94123 , -1.29124 , -2.16045 , -1.78939 , -3.11111 , -1.87178 , -1.57671 , -1.70311 , -1.81208 , -1.94475 , -2.53058 , -1.7042 , -2.08109 , -1.84416 , -0.61073 , -1.77145 , 0.16999 , -0.0585339 , 0.00401133 , 0.397726 , 0.851111 , 0.264187 , 0.59573 , -0.158263 , 0.584362 , 1.20835 , 0.927573 , 1.13895 , 0.64648 , 2.18747 , 1.68909 , 0.451194};
1324 Float_t dA2 = 2.8, dC2 = 3.3;
1327 for (Int_t i = begin; i < end; ++i) {
1328 Float_t tempAdc = esdV0->GetAdc(i)/factors[runRange][i];
1329 Float_t tempTime = (i >= 32) ? esdV0->GetTime(i)+dA+timeShift[i]+dA2 : esdV0->GetTime(i)+dC+timeShift[i]+dC2;
1330 if (esdV0->GetTime(i) >= 1e-6 &&
1331 tempTime > fV0HwWinLow && tempTime < fV0HwWinHigh &&
1332 tempAdc > fV0HwAdcThr)
1338 for (Int_t i = begin; i < end; ++i) {
1339 Float_t tempAdc = esdV0->GetAdc(i)/factors[runRange][i];
1340 Float_t tempTime = (i >= 32) ? esdV0->GetTime(i)+dA : esdV0->GetTime(i)+dC;
1341 Float_t tempRawTime = (i >= 32) ? esdV0->GetTime(i)+dA+timeShift[i]+dA2 : esdV0->GetTime(i)+dC+timeShift[i]+dC2;
1342 if (esdV0->GetTime(i) >= 1e-6 &&
1343 tempRawTime < 125.0 &&
1344 tempAdc > fV0AdcThr) {
1353 for (Int_t i = begin; i < end; ++i) {
1354 if (esdV0->GetTime(i) >= 1e-6 &&
1355 esdV0->GetTime(i) > fV0HwWinLow && esdV0->GetTime(i) < fV0HwWinHigh &&
1356 esdV0->GetAdc(i) > fV0HwAdcThr)
1362 for (Int_t i = begin; i < end; ++i) {
1363 if (esdV0->GetTime(i) > 1e-6 && esdV0->GetAdc(i) > fV0AdcThr) {
1364 Float_t correctedTime = V0CorrectLeadingTime(i, esdV0->GetTime(i), esdV0->GetAdc(i),aEsd->GetRunNumber());
1365 Float_t timeWeight = V0LeadingTimeWeight(esdV0->GetAdc(i));
1366 time += correctedTime*timeWeight;
1368 weight += timeWeight;
1376 time += fV0TimeOffset;
1380 if (side == kASide && fHistV0A)
1381 fHistV0A->Fill(time);
1382 if (side == kCSide && fHistV0C)
1383 fHistV0C->Fill(time);
1388 if (time > 68 && time < 100)
1390 if (time > 54 && time < 57.5)
1392 if (time > 57.5 && time < 68)
1398 if (time > 75.5 && time < 100)
1400 if (time > 69.5 && time < 73)
1402 if (time > 55 && time < 69.5)
1409 Float_t AliTriggerAnalysis::V0CorrectLeadingTime(Int_t i, Float_t time, Float_t adc, Int_t runNumber) const
1411 // Correct for slewing and align the channels
1413 // Authors: Cvetan Cheshkov / Raphael Tieulent
1415 if (time == 0) return 0;
1418 Float_t timeShift[64] = {0.477957 , 0.0889999 , 0.757669 , 0.205439 , 0.239666 , -0.183705 , 0.442873 , -0.281366 , 0.260976 , 0.788995 , 0.974758 , 0.548532 , 0.495023 , 0.868472 , 0.661167 , 0.358307 , 0.221243 , 0.530179 , 1.26696 , 1.33082 , 1.27086 , 1.77133 , 1.10253 , 0.634806 , 2.14838 , 1.50212 , 1.59253 , 1.66122 , 1.16957 , 1.52056 , 1.47791 , 1.81905 , -1.94123 , -1.29124 , -2.16045 , -1.78939 , -3.11111 , -1.87178 , -1.57671 , -1.70311 , -1.81208 , -1.94475 , -2.53058 , -1.7042 , -2.08109 , -1.84416 , -0.61073 , -1.77145 , 0.16999 , -0.0585339 , 0.00401133 , 0.397726 , 0.851111 , 0.264187 , 0.59573 , -0.158263 , 0.584362 , 1.20835 , 0.927573 , 1.13895 , 0.64648 , 2.18747 , 1.68909 , 0.451194};
1420 if(runNumber < 106031)
1421 time -= timeShift[i];
1423 // Slewing correction
1424 if (adc == 0) return time;
1426 Float_t p1 = 1.57345e1;
1427 Float_t p2 =-4.25603e-1;
1429 if(runNumber >= 106031) adc *= (2.5/4.0);
1430 return (time - p1*TMath::Power(adc,p2));
1433 Float_t AliTriggerAnalysis::V0LeadingTimeWeight(Float_t adc) const
1435 if (adc < 1e-6) return 0;
1437 Float_t p1 = 40.211;
1438 Float_t p2 =-4.25603e-1;
1439 Float_t p3 = 0.5646;
1441 return 1./(p1*p1*TMath::Power(adc,2.*(p2-1.))+p3*p3);
1445 Bool_t AliTriggerAnalysis::ZDCTDCTrigger(const AliESDEvent* aEsd, AliceSide side, Bool_t useZN, Bool_t useZP, Bool_t fillHists) const
1447 // Returns if ZDC triggered, based on TDC information
1449 AliESDZDC *esdZDC = aEsd->GetESDZDC();
1451 Bool_t zdcNA = kFALSE;
1452 Bool_t zdcNC = kFALSE;
1453 Bool_t zdcPA = kFALSE;
1454 Bool_t zdcPC = kFALSE;
1457 // If it's MC, we use the energy
1458 Double_t minEnergy = 0;
1459 Double_t zNCEnergy = esdZDC->GetZDCN1Energy();
1460 Double_t zPCEnergy = esdZDC->GetZDCP1Energy();
1461 Double_t zNAEnergy = esdZDC->GetZDCN2Energy();
1462 Double_t zPAEnergy = esdZDC->GetZDCP2Energy();
1463 zdcNA = (zNAEnergy>minEnergy);
1464 zdcNC = (zNCEnergy>minEnergy);
1465 zdcPA = (zPAEnergy>minEnergy);
1466 zdcPC = (zPCEnergy>minEnergy);
1471 Bool_t tdc[32] = {kFALSE};
1472 for(Int_t itdc=0; itdc<32; itdc++){
1473 for(Int_t i=0; i<4; i++){
1474 if (esdZDC->GetZDCTDCData(itdc, i) != 0){
1478 if(fillHists && tdc[itdc]) {
1479 fHistTDCZDC->Fill(itdc);
1488 if (side == kASide) return ((useZP&&zdcPA) || (useZN&&zdcNA));
1489 if (side == kCSide) return ((useZP&&zdcPC) || (useZN&&zdcNC));
1493 Bool_t AliTriggerAnalysis::ZDCTimeTrigger(const AliESDEvent *aEsd, Bool_t fillHists) const
1495 // This method implements a selection
1496 // based on the timing in both sides of zdcN
1497 // It can be used in order to eliminate
1498 // parasitic collisions
1499 Bool_t zdcAccept = kFALSE;
1500 AliESDZDC *esdZDC = aEsd->GetESDZDC();
1503 UInt_t esdFlag = esdZDC->GetESDQuality();
1505 Bool_t znaFired=kFALSE, zpaFired=kFALSE;
1506 Bool_t zem1Fired=kFALSE, zem2Fired=kFALSE;
1507 Bool_t zncFired=kFALSE, zpcFired=kFALSE;
1509 // **** Trigger patterns
1510 if((esdFlag & 0x00000001) == 0x00000001) znaFired=kTRUE;
1511 if((esdFlag & 0x00000002) == 0x00000002) zpaFired=kTRUE;
1512 if((esdFlag & 0x00000004) == 0x00000004) zem1Fired=kTRUE;
1513 if((esdFlag & 0x00000008) == 0x00000008) zem2Fired=kTRUE;
1514 if((esdFlag & 0x00000010) == 0x00000010) zncFired=kTRUE;
1515 if((esdFlag & 0x00000020) == 0x00000020) zpcFired=kTRUE;
1516 zdcAccept = (znaFired | zncFired);
1519 for(Int_t i = 0; i < 4; ++i) {
1520 if (esdZDC->GetZDCTDCData(10,i) != 0) {
1521 Float_t tdcC = 0.025*(esdZDC->GetZDCTDCData(10,i)-esdZDC->GetZDCTDCData(14,i));
1522 Float_t tdcCcorr = esdZDC->GetZDCTDCCorrected(10,i);
1523 for(Int_t j = 0; j < 4; ++j) {
1524 if (esdZDC->GetZDCTDCData(12,j) != 0) {
1525 Float_t tdcA = 0.025*(esdZDC->GetZDCTDCData(12,j)-esdZDC->GetZDCTDCData(14,j));
1527 Float_t tdcAcorr = esdZDC->GetZDCTDCCorrected(12,j);
1529 fHistTimeZDC->Fill(tdcC-tdcA,tdcC+tdcA);
1530 fHistTimeCorrZDC->Fill(tdcCcorr-tdcAcorr,tdcCcorr+tdcAcorr);
1532 if (esdZDC->TestBit(AliESDZDC::kCorrectedTDCFilled)) {
1533 if (((tdcCcorr-tdcAcorr-fZDCCutRefDeltaCorr)*(tdcCcorr-tdcAcorr-fZDCCutRefDeltaCorr)/(fZDCCutSigmaDeltaCorr*fZDCCutSigmaDeltaCorr) +
1534 (tdcCcorr+tdcAcorr-fZDCCutRefSumCorr)*(tdcCcorr+tdcAcorr-fZDCCutRefSumCorr)/(fZDCCutSigmaSumCorr*fZDCCutSigmaSumCorr))< 1.0)
1538 if (((tdcC-tdcA-fZDCCutRefDelta)*(tdcC-tdcA-fZDCCutRefDelta)/(fZDCCutSigmaDelta*fZDCCutSigmaDelta) +
1539 (tdcC+tdcA-fZDCCutRefSum)*(tdcC+tdcA-fZDCCutRefSum)/(fZDCCutSigmaSum*fZDCCutSigmaSum))< 1.0)
1550 Bool_t AliTriggerAnalysis::ZDCTimeBGTrigger(const AliESDEvent *aEsd, AliceSide side) const
1552 // This method implements a selection
1553 // based on the timing in of zdcN
1554 // It can be used in order to flag background
1555 // ** So far only implemented for the 2012 pA run **
1557 if(fMC) return kFALSE;
1559 AliESDZDC *zdcData = aEsd->GetESDZDC();
1560 Bool_t zna = kFALSE;
1561 Bool_t znc = kFALSE;
1562 Bool_t znabadhit = kFALSE;
1563 Bool_t zncbadhit = kFALSE;
1565 Float_t tdcC=999, tdcCcorr=999, tdcA=999, tdcAcorr=999;
1566 for(Int_t i = 0; i < 4; ++i) {
1567 if (zdcData->GetZDCTDCData(10,i) != 0) {
1569 tdcC = 0.025*(zdcData->GetZDCTDCData(10,i)-zdcData->GetZDCTDCData(14,i));
1570 tdcCcorr = zdcData->GetZDCTDCCorrected(10,i);
1571 if((TMath::Abs(tdcCcorr)<fZDCCutZNCTimeCorrMax) && (TMath::Abs(tdcCcorr)>fZDCCutZNCTimeCorrMin)) zncbadhit = kTRUE;
1574 for(Int_t j = 0; j < 4; ++j) {
1575 if (zdcData->GetZDCTDCData(12,j) != 0) {
1577 tdcA = 0.025*(zdcData->GetZDCTDCData(12,j)-zdcData->GetZDCTDCData(14,j));
1578 tdcAcorr = zdcData->GetZDCTDCCorrected(12,j);
1579 if((TMath::Abs(tdcAcorr)<fZDCCutZNATimeCorrMax) && (TMath::Abs(tdcAcorr)>fZDCCutZNATimeCorrMin)) znabadhit = kTRUE;
1583 const Int_t runNumber = aEsd->GetRunNumber();
1584 if(runNumber<188124 || (runNumber>188374 && runNumber<194713)){ // FIXME: end of pA-run is not known
1585 AliDebug(3,Form(" ZN BG time cut not implemented for run %d",runNumber));
1589 Bool_t znabg = (zna && znabadhit);
1590 Bool_t zncbg = (znc && zncbadhit);
1592 // Printf("Checking ZN background (time) for run %d, A = %d, time=%2.2f, C = %d, time=%2.2f",runNumber,(Int_t)zna,tdcAcorr,(Int_t)znc,tdcCcorr);
1593 // Printf("Checking ZN background (time) for run %d, A-BG = %d, C-BG = %d",runNumber,(Int_t)znabg,(Int_t)zncbg);
1595 if (side == kASide) return znabg;
1596 if (side == kCSide) return zncbg;
1600 Bool_t AliTriggerAnalysis::ZDCTrigger(const AliESDEvent* aEsd, AliceSide side) const
1602 // Returns if ZDC triggered
1604 AliESDZDC* zdcData = aEsd->GetESDZDC();
1607 AliError("AliESDZDC not available");
1611 UInt_t quality = zdcData->GetESDQuality();
1613 // from Nora's presentation, general first physics meeting 16.10.09
1614 static UInt_t zpc = 0x20;
1615 static UInt_t znc = 0x10;
1616 static UInt_t zem1 = 0x08;
1617 static UInt_t zem2 = 0x04;
1618 static UInt_t zpa = 0x02;
1619 static UInt_t zna = 0x01;
1621 if (side == kASide && ((quality & zpa) || (quality & zna)))
1623 if (side == kCentralBarrel && ((quality & zem1) || (quality & zem2)))
1625 if (side == kCSide && ((quality & zpc) || (quality & znc)))
1631 Int_t AliTriggerAnalysis::FMDHitCombinations(const AliESDEvent* aEsd, AliceSide side, Bool_t fillHists)
1633 // returns number of hit combinations agove threshold
1635 // Authors: FMD team, Hans Dalsgaard (code merged from FMD/AliFMDOfflineTrigger)
1640 // Workaround for AliESDEvent::GetFMDData is not const!
1641 const AliESDFMD* fmdData = (const_cast<AliESDEvent*>(aEsd))->GetFMDData();
1644 AliError("AliESDFMD not available");
1648 Int_t detFrom = (side == kASide) ? 1 : 3;
1649 Int_t detTo = (side == kASide) ? 2 : 3;
1652 Float_t totalMult = 0;
1653 for (UShort_t det=detFrom;det<=detTo;det++) {
1654 Int_t nRings = (det == 1 ? 1 : 2);
1655 for (UShort_t ir = 0; ir < nRings; ir++) {
1656 Char_t ring = (ir == 0 ? 'I' : 'O');
1657 UShort_t nsec = (ir == 0 ? 20 : 40);
1658 UShort_t nstr = (ir == 0 ? 512 : 256);
1659 for (UShort_t sec =0; sec < nsec; sec++) {
1660 for (UShort_t strip = 0; strip < nstr; strip++) {
1661 Float_t mult = fmdData->Multiplicity(det,ring,sec,strip);
1662 if (mult == AliESDFMD::kInvalidMult) continue;
1665 fHistFMDSingle->Fill(mult);
1667 if (mult > fFMDLowCut)
1668 totalMult = totalMult + mult;
1671 if (totalMult > fFMDHitCut)
1675 fHistFMDSum->Fill(totalMult);
1687 Bool_t AliTriggerAnalysis::FMDTrigger(const AliESDEvent* aEsd, AliceSide side)
1689 // Returns if the FMD triggered
1691 // Authors: FMD team, Hans Dalsgaard (code merged from FMD/AliFMDOfflineTrigger)
1693 Int_t triggers = FMDHitCombinations(aEsd, side, kFALSE);
1701 Long64_t AliTriggerAnalysis::Merge(TCollection* list)
1703 // Merge a list of AliMultiplicityCorrection objects with this (needed for
1705 // Returns the number of merged objects (including this).
1710 if (list->IsEmpty())
1713 TIterator* iter = list->MakeIterator();
1716 // collections of all histograms
1717 const Int_t nHists = 14;
1718 TList collections[nHists];
1721 while ((obj = iter->Next())) {
1723 AliTriggerAnalysis* entry = dynamic_cast<AliTriggerAnalysis*> (obj);
1728 collections[n++].Add(entry->fHistV0A);
1729 collections[n++].Add(entry->fHistV0C);
1730 collections[n++].Add(entry->fHistZDC);
1731 collections[n++].Add(entry->fHistTDCZDC);
1732 collections[n++].Add(entry->fHistTimeZDC);
1733 collections[n++].Add(entry->fHistTimeCorrZDC);
1734 collections[n++].Add(entry->fHistFMDA);
1735 collections[n++].Add(entry->fHistFMDC);
1736 collections[n++].Add(entry->fHistFMDSingle);
1737 collections[n++].Add(entry->fHistFMDSum);
1738 collections[n++].Add(entry->fHistBitsSPD);
1739 collections[n++].Add(entry->fHistFiredBitsSPD);
1740 collections[n++].Add(entry->fHistSPDClsVsTrk);
1741 collections[n++].Add(entry->fHistT0);
1743 // merge fTriggerClasses
1744 TIterator* iter2 = entry->fTriggerClasses->MakeIterator();
1745 TObjString* obj2 = 0;
1746 while ((obj2 = dynamic_cast<TObjString*> (iter2->Next())))
1748 TParameter<Long64_t>* param2 = static_cast<TParameter<Long64_t>*> (entry->fTriggerClasses->GetValue(obj2));
1750 TParameter<Long64_t>* param1 = dynamic_cast<TParameter<Long64_t>*> (fTriggerClasses->GetValue(obj2));
1753 param1->SetVal(param1->GetVal() + param2->GetVal());
1757 param1 = dynamic_cast<TParameter<Long64_t>*> (param2->Clone());
1758 fTriggerClasses->Add(new TObjString(obj2->String()), param1);
1768 fHistV0A->Merge(&collections[n++]);
1769 fHistV0C->Merge(&collections[n++]);
1770 fHistZDC->Merge(&collections[n++]);
1771 fHistTDCZDC->Merge(&collections[n++]);
1773 fHistTimeZDC->Merge(&collections[n++]);
1776 if (fHistTimeCorrZDC)
1777 fHistTimeCorrZDC->Merge(&collections[n++]);
1780 fHistFMDA->Merge(&collections[n++]);
1781 fHistFMDC->Merge(&collections[n++]);
1782 fHistFMDSingle->Merge(&collections[n++]);
1783 fHistFMDSum->Merge(&collections[n++]);
1784 fHistBitsSPD->Merge(&collections[n++]);
1785 fHistFiredBitsSPD->Merge(&collections[n++]);
1786 fHistSPDClsVsTrk->Merge(&collections[n++]);
1787 fHistT0->Merge(&collections[n++]);
1793 void AliTriggerAnalysis::SaveHistograms() const
1795 // write histograms to current directory
1801 fHistBitsSPD->Write();
1802 //fHistBitsSPD->ProjectionX();
1803 //fHistBitsSPD->ProjectionY();
1805 else Printf("Cannot save fHistBitsSPD");
1806 if (fHistFiredBitsSPD) fHistFiredBitsSPD->Write();
1807 else Printf("Cannot save fHistFiredBitsSPD");
1808 if (fHistV0A) fHistV0A->Write();
1809 else Printf("Cannot save fHistV0A");
1810 if (fHistV0C) fHistV0C->Write();
1811 else Printf("Cannot save fHistV0C");
1812 if (fHistZDC) fHistZDC->Write();
1813 else Printf("Cannot save fHistZDC");
1814 if (fHistTDCZDC) fHistTDCZDC->Write();
1815 else Printf("Cannot save fHistTDCZDC");
1816 if (fHistTimeZDC) fHistTimeZDC->Write();
1817 else Printf("Cannot save fHistTimeZDC");
1818 if (fHistTimeCorrZDC) fHistTimeCorrZDC->Write();
1819 else Printf("Cannot save fHistTimeCorrZDC");
1820 if (fHistFMDA) fHistFMDA->Write();
1821 else Printf("Cannot save fHistFMDA");
1822 if (fHistFMDC) fHistFMDC->Write();
1823 else Printf("Cannot save fHistFMDC");
1824 if (fHistFMDSingle) fHistFMDSingle->Write();
1825 else Printf("Cannot save fHistFMDSingle");
1826 if (fHistFMDSum) fHistFMDSum->Write();
1827 else Printf("Cannot save fHistFMDSum");
1828 if (fSPDGFOEfficiency) fSPDGFOEfficiency->Write("fSPDGFOEfficiency");
1829 if (fHistSPDClsVsTrk) fHistSPDClsVsTrk->Write("fHistSPDClsVsTrk");
1830 if (fHistT0) fHistT0->Write("fHistT0");
1832 // else Printf("Cannot save fSPDGFOEfficiency");
1834 fTriggerClasses->Write("fTriggerClasses", TObject::kSingleKey);
1837 void AliTriggerAnalysis::PrintTriggerClasses() const
1839 // print trigger classes
1841 Printf("Trigger Classes:");
1843 Printf("Event count for trigger combinations:");
1846 singleTrigger.SetOwner();
1848 TIterator* iter = fTriggerClasses->MakeIterator();
1849 TObjString* obj = 0;
1850 while ((obj = dynamic_cast<TObjString*> (iter->Next())))
1852 TParameter<Long64_t>* param = static_cast<TParameter<Long64_t>*> (fTriggerClasses->GetValue(obj));
1854 Printf(" %s: %ld triggers", obj->String().Data(), (Long_t)param->GetVal());
1856 TObjArray* tokens = obj->String().Tokenize(" ");
1857 for (Int_t i=0; i<tokens->GetEntries(); i++)
1859 TParameter<Long64_t>* count = dynamic_cast<TParameter<Long64_t>*> (singleTrigger.GetValue(((TObjString*) tokens->At(i))->String().Data()));
1862 count = new TParameter<Long64_t>(((TObjString*) tokens->At(i))->String().Data(), 0);
1863 singleTrigger.Add(new TObjString(((TObjString*) tokens->At(i))->String().Data()), count);
1865 count->SetVal(count->GetVal() + param->GetVal());
1872 Printf("Event count for single trigger:");
1874 iter = singleTrigger.MakeIterator();
1875 while ((obj = dynamic_cast<TObjString*> (iter->Next())))
1877 TParameter<Long64_t>* param = static_cast<TParameter<Long64_t>*> (singleTrigger.GetValue(obj));
1879 Printf(" %s: %ld triggers", obj->String().Data(), (Long_t)param->GetVal());
1883 singleTrigger.DeleteAll();
1887 //----------------------------------------------------------------------------------------------------
1888 AliTriggerAnalysis::T0Decision AliTriggerAnalysis::T0Trigger(const AliESDEvent* aEsd, Bool_t online, Bool_t fillHists)
1890 // Returns the T0 TVDC trigger decision
1892 // argument 'online' is used as a switch between online and offline trigger algorithms
1893 // in online mode return 0TVX
1894 // in offline mode in addition check pile-up and background :
1895 // pile-up readed from ESD: check if TVDC (0TVX module name) has more 1 hit;
1896 // backgroud flag readed from ESD : check in given time interval OrA and OrC were correct but TVDC not
1898 // Based on an algorithm by Alla Maevskaya
1900 const AliESDTZERO* esdT0 = aEsd->GetESDTZERO();
1903 AliError("AliESDTZERO not available");
1906 //???? AliDebug(2,Form("In T0Trigger: %f %f",esdV0->GetV0ATime(),esdV0->GetV0CTime()));
1908 for (Int_t ii=0; ii<5; ii++)
1909 tvdc[ii] = esdT0->GetTVDC(ii);
1910 // Int_t trig=esdT0->GetT0Trig();
1911 // cout<<" T0 trig "<<trig<<endl;
1913 if(fillHists) fHistT0->Fill(tvdc[0]);
1916 if(aEsd->GetHeader()->GetFiredTriggerInputs().Contains("0TVX") ) return kT0BB;
1920 if (esdT0->GetPileupFlag()) return kT0DecPileup;
1921 if (esdT0->GetBackgroundFlag()) return kT0DecBG;
1922 if (tvdc[0]>-5 && tvdc[0]<5 && tvdc[0] != 0) return kT0BB;
1926 if( esdT0->GetT0zVertex()>-12.3 && esdT0->GetT0zVertex() < 10.3) return kT0BB;
1931 //----------------------------------------------------------------------------------------------------
1932 Bool_t AliTriggerAnalysis::EMCALCellsTrigger(const AliESDEvent *aEsd)
1935 // Returns the EMCAL trigger decision
1936 // so far only implemented for LHC11a data
1937 // see http://alisoft.cern.ch/viewvc/trunk/PWGGA/EMCALTasks/AliEmcalPhysicsSelection.cxx?view=markup&root=AliRoot Revision 56136
1940 Bool_t isFired = kTRUE;
1941 const Int_t runNumber = aEsd->GetRunNumber();
1944 // Load EMCAL branches from the manager
1945 AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
1946 am->LoadBranch("EMCALCells.");
1947 am->LoadBranch("CaloClusters");
1951 AliVCaloCells *cells = aEsd->GetEMCALCells();
1952 const Short_t nCells = cells->GetNumberOfCells();
1954 // count cells above threshold per sm
1955 Int_t nCellCount[10] = {0,0,0,0,0,0,0,0,0,0};
1956 for(Int_t iCell=0; iCell<nCells; ++iCell) {
1957 Short_t cellId = cells->GetCellNumber(iCell);
1958 Double_t cellE = cells->GetCellAmplitude(cellId);
1959 Int_t sm = cellId / (24*48);
1964 // Trigger decision for LHC11a
1965 Bool_t isLedEvent = kFALSE;
1966 if ((runNumber>=144871) && (runNumber<=146860)) {
1967 if (nCellCount[4] > 100)
1970 if ((runNumber>=146858) && (runNumber<=146860)) {
1971 if (nCellCount[3]>=35)
1984 //__________________________________________________________________________________________
1985 Bool_t AliTriggerAnalysis::TRDTrigger(const AliESDEvent *esd, Trigger trigger)
1987 // evaluate the TRD trigger conditions,
1988 // so far HCO, HSE, HQU, HJT, HEE
1990 Bool_t isFired = kFALSE;
1992 if(trigger!=kTRDHCO && trigger!=kTRDHJT && trigger!=kTRDHSE && trigger!=kTRDHQU && trigger!=kTRDHEE) {
1993 AliWarning("Beware you are erroneously trying to use this function (wrong trigger)");
1998 AliErrorClass("ESD event pointer is null");
2002 Int_t nTrdTracks = esd->GetNumberOfTrdTracks();
2004 if (nTrdTracks > 0 && (trigger==kTRDHCO) ) {
2009 if(trigger!=kTRDHJT) {
2010 for (Int_t iTrack = 0; iTrack < nTrdTracks; ++iTrack) {
2012 AliESDTrdTrack *trdTrack = esd->GetTrdTrack(iTrack);
2013 if (!trdTrack) continue;
2015 // for the electron triggers we only consider matched tracks
2016 if(trigger==kTRDHQU)
2017 if ( (TMath::Abs(trdTrack->Pt()) > fTRDptHQU) && (trdTrack->GetPID() > fTRDpidHQU) ) {
2022 if(trigger==kTRDHSE)
2023 if ( (TMath::Abs(trdTrack->Pt()) > fTRDptHSE) && (trdTrack->GetPID() > fTRDpidHSE) ) {
2028 if(trigger==kTRDHEE)
2029 if ( (trdTrack->GetSector() >= fTRDminSectorHEE) && (trdTrack->GetSector() <= fTRDmaxSectorHEE) &&
2030 (TMath::Abs(trdTrack->Pt()) > fTRDptHSE) && (trdTrack->GetPID() > fTRDpidHSE) ) {
2036 } else if(trigger==kTRDHJT) {
2038 Int_t nTracks[90] = { 0 }; // stack-wise counted number of tracks above pt threshold
2040 for (Int_t iTrack = 0; iTrack < nTrdTracks; ++iTrack) {
2042 AliESDTrdTrack *trdTrack = esd->GetTrdTrack(iTrack);
2043 if (!trdTrack) continue;
2045 Int_t globalStack = 5*trdTrack->GetSector() + trdTrack->GetStack();
2047 // stack-wise counting of tracks above pt threshold for jet trigger
2048 if (TMath::Abs(trdTrack->GetPt()) >= fTRDptHJT) {
2049 ++nTracks[globalStack];
2053 // check if HJT condition is fulfilled in any stack
2054 for (Int_t iStack = 0; iStack < 90; iStack++) {
2055 if (nTracks[iStack] >= fTRDnHJT) {