AliPhysics  8bb951a (8bb951a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcal.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 #include <RVersion.h>
16 #include "AliAnalysisTaskEmcal.h"
17 
18 #include <TClonesArray.h>
19 #include <TList.h>
20 #include <TObject.h>
21 #include <TH1F.h>
22 #include <TProfile.h>
23 #include <TSystem.h>
24 #include <TFile.h>
25 #include <TChain.h>
26 #include <TKey.h>
27 
28 #include "AliStack.h"
29 #include "AliAODEvent.h"
30 #include "AliAnalysisManager.h"
31 #include "AliCentrality.h"
32 #include "AliEMCALGeometry.h"
33 #include "AliESDEvent.h"
34 #include "AliEmcalParticle.h"
35 #include "AliEventplane.h"
36 #include "AliInputEventHandler.h"
37 #include "AliLog.h"
38 #include "AliMCParticle.h"
39 #include "AliVCluster.h"
40 #include "AliVEventHandler.h"
41 #include "AliVParticle.h"
42 #include "AliAODTrack.h"
43 #include "AliVCaloTrigger.h"
44 #include "AliGenPythiaEventHeader.h"
45 #include "AliAODMCHeader.h"
46 #include "AliMCEvent.h"
47 #include "AliAnalysisUtils.h"
48 #include "AliEMCALTriggerPatchInfo.h"
49 #include "AliEmcalPythiaInfo.h"
50 
51 #include "AliMultSelection.h"
52 
54 
58 
63  AliAnalysisTaskSE("AliAnalysisTaskEmcal"),
64  fPythiaInfoName(""),
65  fForceBeamType(kNA),
66  fGeneralHistograms(kFALSE),
67  fInitialized(kFALSE),
68  fCreateHisto(kTRUE),
69  fCaloCellsName(),
70  fCaloTriggersName(),
71  fCaloTriggerPatchInfoName(),
72  fMinCent(-999),
73  fMaxCent(-999),
74  fMinVz(-999),
75  fMaxVz(999),
76  fTrackPtCut(0),
77  fMinNTrack(0),
78  fZvertexDiff(0.5),
79  fUseAliAnaUtils(kFALSE),
80  fRejectPileup(kFALSE),
81  fTklVsClusSPDCut(kFALSE),
82  fOffTrigger(AliVEvent::kAny),
83  fTrigClass(),
84  fTriggerTypeSel(kND),
85  fNbins(250),
86  fMinBinPt(0),
87  fMaxBinPt(250),
88  fMinPtTrackInEmcal(0),
89  fEventPlaneVsEmcal(-1),
90  fMinEventPlane(-1e6),
91  fMaxEventPlane(1e6),
92  fCentEst("V0M"),
93  fIsEmbedded(kFALSE),
94  fIsPythia(kFALSE),
95  fSelectPtHardBin(-999),
96  fMinMCLabel(0),
97  fMCLabelShift(0),
98  fNcentBins(4),
99  fNeedEmcalGeom(kTRUE),
100  fParticleCollArray(),
101  fClusterCollArray(),
102  fTriggers(0),
103  fEMCalTriggerMode(kOverlapWithLowThreshold),
104  fUseNewCentralityEstimation(kFALSE),
105  fGeneratePythiaInfoObject(kFALSE),
106  fAliAnalysisUtils(0x0),
107  fIsEsd(kFALSE),
108  fGeom(0),
109  fTracks(0),
110  fCaloClusters(0),
111  fCaloCells(0),
112  fCaloTriggers(0),
113  fTriggerPatchInfo(0),
114  fCent(0),
115  fCentBin(-1),
116  fEPV0(-1.0),
117  fEPV0A(-1.0),
118  fEPV0C(-1.0),
119  fNVertCont(0),
120  fNVertSPDCont(0),
121  fBeamType(kNA),
122  fPythiaHeader(0),
123  fPtHard(0),
124  fPtHardBin(0),
125  fNTrials(0),
126  fXsection(0),
127  fPythiaInfo(0),
128  fOutput(0),
129  fHistEventCount(0),
130  fHistTrialsAfterSel(0),
131  fHistEventsAfterSel(0),
132  fHistXsectionAfterSel(0),
133  fHistTrials(0),
134  fHistEvents(0),
135  fHistXsection(0),
136  fHistPtHard(0),
137  fHistCentrality(0),
138  fHistZVertex(0),
139  fHistEventPlane(0),
140  fHistEventRejection(0),
141  fHistTriggerClasses(0)
142 {
143  fVertex[0] = 0;
144  fVertex[1] = 0;
145  fVertex[2] = 0;
146  fVertexSPD[0] = 0;
147  fVertexSPD[1] = 0;
148  fVertexSPD[2] = 0;
149 
150  fParticleCollArray.SetOwner(kTRUE);
151  fClusterCollArray.SetOwner(kTRUE);
152 }
153 
164 AliAnalysisTaskEmcal::AliAnalysisTaskEmcal(const char *name, Bool_t histo) :
165  AliAnalysisTaskSE(name),
166  fPythiaInfoName(""),
167  fForceBeamType(kNA),
168  fGeneralHistograms(kFALSE),
169  fInitialized(kFALSE),
170  fCreateHisto(histo),
171  fCaloCellsName(),
172  fCaloTriggersName(),
173  fCaloTriggerPatchInfoName(),
174  fMinCent(-999),
175  fMaxCent(-999),
176  fMinVz(-999),
177  fMaxVz(999),
178  fTrackPtCut(0),
179  fMinNTrack(0),
180  fZvertexDiff(0.5),
181  fUseAliAnaUtils(kFALSE),
182  fRejectPileup(kFALSE),
183  fTklVsClusSPDCut(kFALSE),
184  fOffTrigger(AliVEvent::kAny),
185  fTrigClass(),
186  fTriggerTypeSel(kND),
187  fNbins(250),
188  fMinBinPt(0),
189  fMaxBinPt(250),
190  fMinPtTrackInEmcal(0),
191  fEventPlaneVsEmcal(-1),
192  fMinEventPlane(-1e6),
193  fMaxEventPlane(1e6),
194  fCentEst("V0M"),
195  fIsEmbedded(kFALSE),
196  fIsPythia(kFALSE),
197  fSelectPtHardBin(-999),
198  fMinMCLabel(0),
199  fMCLabelShift(0),
200  fNcentBins(4),
201  fNeedEmcalGeom(kTRUE),
202  fParticleCollArray(),
203  fClusterCollArray(),
204  fTriggers(0),
205  fEMCalTriggerMode(kOverlapWithLowThreshold),
206  fUseNewCentralityEstimation(kFALSE),
207  fGeneratePythiaInfoObject(kFALSE),
208  fAliAnalysisUtils(0x0),
209  fIsEsd(kFALSE),
210  fGeom(0),
211  fTracks(0),
212  fCaloClusters(0),
213  fCaloCells(0),
214  fCaloTriggers(0),
215  fTriggerPatchInfo(0),
216  fCent(0),
217  fCentBin(-1),
218  fEPV0(-1.0),
219  fEPV0A(-1.0),
220  fEPV0C(-1.0),
221  fNVertCont(0),
222  fNVertSPDCont(0),
223  fBeamType(kNA),
224  fPythiaHeader(0),
225  fPtHard(0),
226  fPtHardBin(0),
227  fNTrials(0),
228  fXsection(0),
229  fPythiaInfo(0),
230  fOutput(0),
231  fHistEventCount(0),
232  fHistTrialsAfterSel(0),
233  fHistEventsAfterSel(0),
234  fHistXsectionAfterSel(0),
235  fHistTrials(0),
236  fHistEvents(0),
237  fHistXsection(0),
238  fHistPtHard(0),
239  fHistCentrality(0),
240  fHistZVertex(0),
241  fHistEventPlane(0),
242  fHistEventRejection(0),
243  fHistTriggerClasses(0)
244 {
245  fVertex[0] = 0;
246  fVertex[1] = 0;
247  fVertex[2] = 0;
248  fVertexSPD[0] = 0;
249  fVertexSPD[1] = 0;
250  fVertexSPD[2] = 0;
251  fParticleCollArray.SetOwner(kTRUE);
252  fClusterCollArray.SetOwner(kTRUE);
253 
254  if (fCreateHisto) {
255  DefineOutput(1, TList::Class());
256  }
257 }
258 
263 {
264 }
265 
272 void AliAnalysisTaskEmcal::SetClusPtCut(Double_t cut, Int_t c)
273 {
275  if (cont) cont->SetClusPtCut(cut);
276  else AliError(Form("%s in SetClusPtCut(...): container %d not found",GetName(),c));
277 }
278 
286 void AliAnalysisTaskEmcal::SetClusTimeCut(Double_t min, Double_t max, Int_t c)
287 {
289  if (cont) cont->SetClusTimeCut(min,max);
290  else AliError(Form("%s in SetClusTimeCut(...): container %d not found",GetName(),c));
291 }
292 
299 void AliAnalysisTaskEmcal::SetTrackPtCut(Double_t cut, Int_t c)
300 {
302  if (cont) cont->SetParticlePtCut(cut);
303  else AliError(Form("%s in SetTrackPtCut(...): container %d not found",GetName(),c));
304 
305  fTrackPtCut = cut;
306 }
307 
315 void AliAnalysisTaskEmcal::SetTrackEtaLimits(Double_t min, Double_t max, Int_t c)
316 {
318  if (cont) cont->SetParticleEtaLimits(min,max);
319  else AliError(Form("%s in SetTrackPtCut(...): container %d not found",GetName(),c));
320 }
321 
329 void AliAnalysisTaskEmcal::SetTrackPhiLimits(Double_t min, Double_t max, Int_t c)
330 {
332  if (cont) cont->SetParticlePhiLimits(min,max);
333  else AliError(Form("%s in SetTrackPhiLimits(...): container %d not found",GetName(),c));
334 }
335 
356 {
357  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
358  if (mgr) {
359  AliVEventHandler *evhand = mgr->GetInputEventHandler();
360  if (evhand) {
361  if (evhand->InheritsFrom("AliESDInputHandler")) {
362  fIsEsd = kTRUE;
363  }
364  else {
365  fIsEsd = kFALSE;
366  }
367  }
368  else {
369  AliError("Event handler not found!");
370  }
371  }
372  else {
373  AliError("Analysis manager not found!");
374  }
375 
376  if (!fCreateHisto)
377  return;
378 
379  OpenFile(1);
380  fOutput = new TList();
381  fOutput->SetOwner();
382 
383  if (fForceBeamType == kpp)
384  fNcentBins = 1;
385 
386  if (!fGeneralHistograms)
387  return;
388 
389  if (fIsPythia) {
390  fHistTrialsAfterSel = new TH1F("fHistTrialsAfterSel", "fHistTrialsAfterSel", 11, 0, 11);
391  fHistTrialsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
392  fHistTrialsAfterSel->GetYaxis()->SetTitle("trials");
394 
395  fHistEventsAfterSel = new TH1F("fHistEventsAfterSel", "fHistEventsAfterSel", 11, 0, 11);
396  fHistEventsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
397  fHistEventsAfterSel->GetYaxis()->SetTitle("total events");
399 
400  fHistXsectionAfterSel = new TProfile("fHistXsectionAfterSel", "fHistXsectionAfterSel", 11, 0, 11);
401  fHistXsectionAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
402  fHistXsectionAfterSel->GetYaxis()->SetTitle("xsection");
404 
405  fHistTrials = new TH1F("fHistTrials", "fHistTrials", 11, 0, 11);
406  fHistTrials->GetXaxis()->SetTitle("p_{T} hard bin");
407  fHistTrials->GetYaxis()->SetTitle("trials");
408  fOutput->Add(fHistTrials);
409 
410  fHistEvents = new TH1F("fHistEvents", "fHistEvents", 11, 0, 11);
411  fHistEvents->GetXaxis()->SetTitle("p_{T} hard bin");
412  fHistEvents->GetYaxis()->SetTitle("total events");
413  fOutput->Add(fHistEvents);
414 
415  fHistXsection = new TProfile("fHistXsection", "fHistXsection", 11, 0, 11);
416  fHistXsection->GetXaxis()->SetTitle("p_{T} hard bin");
417  fHistXsection->GetYaxis()->SetTitle("xsection");
418  fOutput->Add(fHistXsection);
419 
420  const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
421  const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
422 
423  for (Int_t i = 1; i < 12; i++) {
424  fHistTrialsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
425  fHistEventsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
426 
427  fHistTrials->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
428  fHistXsection->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
429  fHistEvents->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
430  }
431 
432  fHistPtHard = new TH1F("fHistPtHard", "fHistPtHard", fNbins*2, fMinBinPt, fMaxBinPt*4);
433  fHistPtHard->GetXaxis()->SetTitle("p_{T,hard} (GeV/c)");
434  fHistPtHard->GetYaxis()->SetTitle("counts");
435  fOutput->Add(fHistPtHard);
436  }
437 
438  fHistZVertex = new TH1F("fHistZVertex","Z vertex position", 60, -30, 30);
439  fHistZVertex->GetXaxis()->SetTitle("z");
440  fHistZVertex->GetYaxis()->SetTitle("counts");
441  fOutput->Add(fHistZVertex);
442 
443  if (fForceBeamType != kpp) {
444  fHistCentrality = new TH1F("fHistCentrality","Event centrality distribution", 200, 0, 100);
445  fHistCentrality->GetXaxis()->SetTitle("Centrality (%)");
446  fHistCentrality->GetYaxis()->SetTitle("counts");
447  fOutput->Add(fHistCentrality);
448 
449  fHistEventPlane = new TH1F("fHistEventPlane","Event plane", 120, -TMath::Pi(), TMath::Pi());
450  fHistEventPlane->GetXaxis()->SetTitle("event plane");
451  fHistEventPlane->GetYaxis()->SetTitle("counts");
452  fOutput->Add(fHistEventPlane);
453  }
454 
455  fHistEventRejection = new TH1F("fHistEventRejection","Reasons to reject event",20,0,20);
456 #if ROOT_VERSION_CODE < ROOT_VERSION(6,4,2)
457  fHistEventRejection->SetBit(TH1::kCanRebin);
458 #else
459  fHistEventRejection->SetCanExtend(TH1::kAllAxes);
460 #endif
461  fHistEventRejection->GetXaxis()->SetBinLabel(1,"PhysSel");
462  fHistEventRejection->GetXaxis()->SetBinLabel(2,"trigger");
463  fHistEventRejection->GetXaxis()->SetBinLabel(3,"trigTypeSel");
464  fHistEventRejection->GetXaxis()->SetBinLabel(4,"Cent");
465  fHistEventRejection->GetXaxis()->SetBinLabel(5,"vertex contr.");
466  fHistEventRejection->GetXaxis()->SetBinLabel(6,"Vz");
467  fHistEventRejection->GetXaxis()->SetBinLabel(7,"VzSPD");
468  fHistEventRejection->GetXaxis()->SetBinLabel(8,"trackInEmcal");
469  fHistEventRejection->GetXaxis()->SetBinLabel(9,"minNTrack");
470  fHistEventRejection->GetXaxis()->SetBinLabel(10,"VtxSel2013pA");
471  fHistEventRejection->GetXaxis()->SetBinLabel(11,"PileUp");
472  fHistEventRejection->GetXaxis()->SetBinLabel(12,"EvtPlane");
473  fHistEventRejection->GetXaxis()->SetBinLabel(13,"SelPtHardBin");
474  fHistEventRejection->GetXaxis()->SetBinLabel(14,"Bkg evt");
475  fHistEventRejection->GetYaxis()->SetTitle("counts");
477 
478  fHistTriggerClasses = new TH1F("fHistTriggerClasses","fHistTriggerClasses",3,0,3);
479 #if ROOT_VERSION_CODE < ROOT_VERSION(6,4,2)
480  fHistTriggerClasses->SetBit(TH1::kCanRebin);
481 #else
482  fHistTriggerClasses->SetCanExtend(TH1::kAllAxes);
483 #endif
485 
486  fHistEventCount = new TH1F("fHistEventCount","fHistEventCount",2,0,2);
487  fHistEventCount->GetXaxis()->SetBinLabel(1,"Accepted");
488  fHistEventCount->GetXaxis()->SetBinLabel(2,"Rejected");
489  fHistEventCount->GetYaxis()->SetTitle("counts");
490  fOutput->Add(fHistEventCount);
491 
492  PostData(1, fOutput);
493 }
494 
509 {
510  if (fIsPythia) {
514  fHistPtHard->Fill(fPtHard);
515  }
516 
517  fHistZVertex->Fill(fVertex[2]);
518 
519  if (fForceBeamType != kpp) {
520  fHistCentrality->Fill(fCent);
521  fHistEventPlane->Fill(fEPV0);
522  }
523 
524  TObjArray* triggerClasses = InputEvent()->GetFiredTriggerClasses().Tokenize(" ");
525  TIter next(triggerClasses);
526  TObjString* triggerClass = 0;
527  while ((triggerClass = static_cast<TObjString*>(next()))) {
528  fHistTriggerClasses->Fill(triggerClass->GetString(), 1);
529  }
530  delete triggerClasses;
531  triggerClasses = 0;
532 
533  return kTRUE;
534 }
535 
555 void AliAnalysisTaskEmcal::UserExec(Option_t *option)
556 {
557  if (!fInitialized)
558  ExecOnce();
559 
560  if (!fInitialized)
561  return;
562 
563  if (!RetrieveEventObjects())
564  return;
565 
566  if (IsEventSelected()) {
567  if (fGeneralHistograms) fHistEventCount->Fill("Accepted",1);
568  }
569  else {
570  if (fGeneralHistograms) fHistEventCount->Fill("Rejected",1);
571  return;
572  }
573 
575  if (!FillGeneralHistograms())
576  return;
577  }
578 
579  if (!Run())
580  return;
581 
582  if (fCreateHisto) {
583  if (!FillHistograms())
584  return;
585  }
586 
587  if (fCreateHisto && fOutput) {
588  // information for this iteration of the UserExec in the container
589  PostData(1, fOutput);
590  }
591 }
592 
601 Bool_t AliAnalysisTaskEmcal::AcceptCluster(AliVCluster *clus, Int_t c) const
602 {
603  AliWarning("AliAnalysisTaskEmcal::AcceptCluster method is deprecated. Please use GetCusterContainer(c)->AcceptCluster(clus).");
604 
605  if (!clus) return kFALSE;
606 
608  if (!cont) {
609  AliError(Form("%s:Container %d not found",GetName(),c));
610  return 0;
611  }
612  UInt_t rejectionReason = 0;
613  return cont->AcceptCluster(clus, rejectionReason);
614 }
615 
624 Bool_t AliAnalysisTaskEmcal::AcceptTrack(AliVParticle *track, Int_t c) const
625 {
626  AliWarning("AliAnalysisTaskEmcal::AcceptTrack method is deprecated. Please use GetParticleContainer(c)->AcceptParticle(clus).");
627 
628  if (!track) return kFALSE;
629 
631  if (!cont) {
632  AliError(Form("%s:Container %d not found",GetName(),c));
633  return 0;
634  }
635 
636  UInt_t rejectionReason = 0;
637  return cont->AcceptParticle(track, rejectionReason);
638 }
639 
651 Bool_t AliAnalysisTaskEmcal::PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
652 {
653 
654  TString file(currFile);
655  fXsec = 0;
656  fTrials = 1;
657 
658  if (file.Contains(".zip#")) {
659  Ssiz_t pos1 = file.Index("root_archive",12,0,TString::kExact);
660  Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
661  Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
662  file.Replace(pos+1,pos2-pos1,"");
663  } else {
664  // not an archive take the basename....
665  file.ReplaceAll(gSystem->BaseName(file.Data()),"");
666  }
667  AliDebug(1,Form("File name: %s",file.Data()));
668 
669  // Get the pt hard bin
670  TString strPthard(file);
671 
672  strPthard.Remove(strPthard.Last('/'));
673  strPthard.Remove(strPthard.Last('/'));
674  if (strPthard.Contains("AOD")) strPthard.Remove(strPthard.Last('/'));
675  strPthard.Remove(0,strPthard.Last('/')+1);
676  if (strPthard.IsDec())
677  pthard = strPthard.Atoi();
678  else
679  AliWarning(Form("Could not extract file number from path %s", strPthard.Data()));
680 
681  // problem that we cannot really test the existance of a file in a archive so we have to live with open error message from root
682  TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root"));
683 
684  if (!fxsec) {
685  // next trial fetch the histgram file
686  fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
687  if (!fxsec) {
688  // not a severe condition but inciate that we have no information
689  return kFALSE;
690  } else {
691  // find the tlist we want to be independtent of the name so use the Tkey
692  TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
693  if (!key) {
694  fxsec->Close();
695  return kFALSE;
696  }
697  TList *list = dynamic_cast<TList*>(key->ReadObj());
698  if (!list) {
699  fxsec->Close();
700  return kFALSE;
701  }
702  fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
703  fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
704  fxsec->Close();
705  }
706  } else { // no tree pyxsec.root
707  TTree *xtree = (TTree*)fxsec->Get("Xsection");
708  if (!xtree) {
709  fxsec->Close();
710  return kFALSE;
711  }
712  UInt_t ntrials = 0;
713  Double_t xsection = 0;
714  xtree->SetBranchAddress("xsection",&xsection);
715  xtree->SetBranchAddress("ntrials",&ntrials);
716  xtree->GetEntry(0);
717  fTrials = ntrials;
718  fXsec = xsection;
719  fxsec->Close();
720  }
721  return kTRUE;
722 }
723 
738 {
740  return kTRUE;
741 
742  TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
743  if (!tree) {
744  AliError(Form("%s - UserNotify: No current tree!",GetName()));
745  return kFALSE;
746  }
747 
748  Float_t xsection = 0;
749  Float_t trials = 0;
750  Int_t pthardbin = 0;
751 
752  TFile *curfile = tree->GetCurrentFile();
753  if (!curfile) {
754  AliError(Form("%s - UserNotify: No current file!",GetName()));
755  return kFALSE;
756  }
757 
758  TChain *chain = dynamic_cast<TChain*>(tree);
759  if (chain) tree = chain->GetTree();
760 
761  Int_t nevents = tree->GetEntriesFast();
762 
763  PythiaInfoFromFile(curfile->GetName(), xsection, trials, pthardbin);
764 
765  // TODO: Workaround
766  if ((pthardbin < 0) || (pthardbin > 10)) pthardbin = 0;
767 
768  fHistTrials->Fill(pthardbin, trials);
769  fHistXsection->Fill(pthardbin, xsection);
770  fHistEvents->Fill(pthardbin, nevents);
771 
772  return kTRUE;
773 }
774 
780 {
781  if (!fPythiaInfoName.IsNull() && !fPythiaInfo) {
782  fPythiaInfo = dynamic_cast<AliEmcalPythiaInfo*>(event->FindListObject(fPythiaInfoName));
783  if (!fPythiaInfo) {
784  AliError(Form("%s: Could not retrieve parton infos! %s!", GetName(), fPythiaInfoName.Data()));
785  return;
786  }
787  }
788 }
789 
801 {
802  if (!InputEvent()) {
803  AliError(Form("%s: Could not retrieve event! Returning!", GetName()));
804  return;
805  }
806 
807  LoadPythiaInfo(InputEvent());
808 
809  if (fNeedEmcalGeom) {
810  fGeom = AliEMCALGeometry::GetInstanceFromRunNumber(InputEvent()->GetRunNumber());
811  if (!fGeom) {
812  AliFatal(Form("%s: Can not get EMCal geometry instance. If you do not need the EMCal geometry, disable it by setting task->SetNeedEmcalGeometry(kFALSE).", GetName()));
813  return;
814  }
815  }
816 
817  if (fEventPlaneVsEmcal >= 0) {
818  if (fGeom) {
819  Double_t ep = (fGeom->GetArm1PhiMax() + fGeom->GetArm1PhiMin()) / 2 * TMath::DegToRad() + fEventPlaneVsEmcal - TMath::Pi();
820  fMinEventPlane = ep - TMath::Pi() / 4;
821  fMaxEventPlane = ep + TMath::Pi() / 4;
822  }
823  else {
824  AliWarning("Could not set event plane limits because EMCal geometry was not loaded!");
825  }
826  }
827 
828  //Load all requested track branches - each container knows name already
829  for (Int_t i =0; i<fParticleCollArray.GetEntriesFast(); i++) {
830  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.At(i));
831  cont->SetArray(InputEvent());
832  }
833 
834  if (fParticleCollArray.GetEntriesFast()>0) {
836  if (!fTracks) {
837  AliError(Form("%s: Could not retrieve first track branch!", GetName()));
838  return;
839  }
840  }
841 
842  //Load all requested cluster branches - each container knows name already
843  for (Int_t i =0; i<fClusterCollArray.GetEntriesFast(); i++) {
844  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.At(i));
845  cont->SetArray(InputEvent());
846  }
847 
848  if (fClusterCollArray.GetEntriesFast()>0) {
850  if (!fCaloClusters) {
851  AliError(Form("%s: Could not retrieve first cluster branch!", GetName()));
852  return;
853  }
854  }
855 
856  if (!fCaloCellsName.IsNull() && !fCaloCells) {
857  fCaloCells = dynamic_cast<AliVCaloCells*>(InputEvent()->FindListObject(fCaloCellsName));
858  if (!fCaloCells) {
859  AliError(Form("%s: Could not retrieve cells %s!", GetName(), fCaloCellsName.Data()));
860  return;
861  }
862  }
863 
864  if (!fCaloTriggersName.IsNull() && !fCaloTriggers) {
865  fCaloTriggers = dynamic_cast<AliVCaloTrigger*>(InputEvent()->FindListObject(fCaloTriggersName));
866  if (!fCaloTriggers) {
867  AliError(Form("%s: Could not retrieve calo triggers %s!", GetName(), fCaloTriggersName.Data()));
868  return;
869  }
870  }
871 
872  if (!fCaloTriggerPatchInfoName.IsNull() && !fTriggerPatchInfo) {
873  fTriggerPatchInfo = GetArrayFromEvent(fCaloTriggerPatchInfoName.Data(),"AliEMCALTriggerPatchInfo");
874  if (!fTriggerPatchInfo) {
875  AliError(Form("%s: Could not retrieve calo trigger patch info %s!", GetName(), fCaloTriggerPatchInfoName.Data()));
876  return;
877  }
878 
879  }
880 
881  fInitialized = kTRUE;
882 }
883 
890 {
891 
892  if (fForceBeamType != kNA)
893  return fForceBeamType;
894 
895  AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
896  if (esd) {
897  const AliESDRun *run = esd->GetESDRun();
898  TString beamType = run->GetBeamType();
899  if (beamType == "p-p")
900  return kpp;
901  else if (beamType == "A-A")
902  return kAA;
903  else if (beamType == "p-A")
904  return kpA;
905  else
906  return kNA;
907  } else {
908  Int_t runNumber = InputEvent()->GetRunNumber();
909  if ((runNumber >= 136851 && runNumber <= 139517) || // LHC10h
910  (runNumber >= 166529 && runNumber <= 170593)) { // LHC11h
911  return kAA;
912  } else if ((runNumber>=188365 && runNumber <= 188366) || // LHC12g
913  (runNumber >= 195344 && runNumber <= 196608)) { // LHC13b-f
914  return kpA;
915  } else {
916  return kpp;
917  }
918  }
919 }
920 
928 {
929  if (!fTriggerPatchInfo)
930  return 0;
931 
932  //number of patches in event
933  Int_t nPatch = fTriggerPatchInfo->GetEntries();
934 
935  //loop over patches to define trigger type of event
936  Int_t nG1 = 0;
937  Int_t nG2 = 0;
938  Int_t nJ1 = 0;
939  Int_t nJ2 = 0;
940  Int_t nL0 = 0;
941  AliEMCALTriggerPatchInfo *patch;
942  for (Int_t iPatch = 0; iPatch < nPatch; iPatch++) {
943  patch = (AliEMCALTriggerPatchInfo*)fTriggerPatchInfo->At( iPatch );
944  if (patch->IsGammaHigh()) nG1++;
945  if (patch->IsGammaLow()) nG2++;
946  if (patch->IsJetHigh()) nJ1++;
947  if (patch->IsJetLow()) nJ2++;
948  if (patch->IsLevel0()) nL0++;
949  }
950 
951  AliDebug(2, "Patch summary: ");
952  AliDebug(2, Form("Number of patches: %d", nPatch));
953  AliDebug(2, Form("Jet: low[%d], high[%d]" ,nJ2, nJ1));
954  AliDebug(2, Form("Gamma: low[%d], high[%d]" ,nG2, nG1));
955 
956  ULong_t triggers(0);
957  if (nL0>0) SETBIT(triggers, kL0);
958  if (nG1>0) SETBIT(triggers, kG1);
959  if (nG2>0) SETBIT(triggers, kG2);
960  if (nJ1>0) SETBIT(triggers, kJ1);
961  if (nJ2>0) SETBIT(triggers, kJ2);
962  return triggers;
963 }
964 
972 {
973  //
974  if(trigger==kND) {
975  AliWarning(Form("%s: Requesting undefined trigger type!", GetName()));
976  return kFALSE;
977  }
978  //MV: removing this logic which as far as I can see doesn't make any sense
979  // if(trigger & kND){
980  // return fTriggers == 0;
981  // }
982  return TESTBIT(fTriggers, trigger);
983 }
984 
1007 {
1008  if (fOffTrigger != AliVEvent::kAny) {
1009  UInt_t res = 0;
1010  const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
1011  if (eev) {
1012  res = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
1013  } else {
1014  const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
1015  if (aev) {
1016  res = ((AliVAODHeader*)aev->GetHeader())->GetOfflineTrigger();
1017  }
1018  }
1019  if ((res & fOffTrigger) == 0) {
1020  if (fGeneralHistograms) fHistEventRejection->Fill("PhysSel",1);
1021  return kFALSE;
1022  }
1023  }
1024 
1025  if (!fTrigClass.IsNull()) {
1026  TString fired;
1027  const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
1028  if (eev) {
1029  fired = eev->GetFiredTriggerClasses();
1030  } else {
1031  const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
1032  if (aev) {
1033  fired = aev->GetFiredTriggerClasses();
1034  }
1035  }
1036  if (!fired.Contains("-B-")) {
1037  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
1038  return kFALSE;
1039  }
1040 
1041  TObjArray *arr = fTrigClass.Tokenize("|");
1042  if (!arr) {
1043  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
1044  return kFALSE;
1045  }
1046  Bool_t match = 0;
1047  for (Int_t i=0;i<arr->GetEntriesFast();++i) {
1048  TObject *obj = arr->At(i);
1049  if (!obj)
1050  continue;
1051 
1052  //Check if requested trigger was fired
1053  TString objStr = obj->GetName();
1055  (objStr.Contains("J1") || objStr.Contains("J2") || objStr.Contains("G1") || objStr.Contains("G2"))) {
1056  // This is relevant for EMCal triggers with 2 thresholds
1057  // If the kOverlapWithLowThreshold was requested than the overlap between the two triggers goes with the lower threshold trigger
1058  TString trigType1 = "J1";
1059  TString trigType2 = "J2";
1060  if(objStr.Contains("G")) {
1061  trigType1 = "G1";
1062  trigType2 = "G2";
1063  }
1064  if(objStr.Contains(trigType2) && fired.Contains(trigType2.Data())) { //requesting low threshold + overlap
1065  match = 1;
1066  break;
1067  }
1068  else if(objStr.Contains(trigType1) && fired.Contains(trigType1.Data()) && !fired.Contains(trigType2.Data())) { //high threshold only
1069  match = 1;
1070  break;
1071  }
1072  }
1073  else {
1074  // If this is not an EMCal trigger, or no particular treatment of EMCal triggers was requested,
1075  // simply check that the trigger was fired
1076  if (fired.Contains(obj->GetName())) {
1077  match = 1;
1078  break;
1079  }
1080  }
1081  }
1082  delete arr;
1083  if (!match) {
1084  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
1085  return kFALSE;
1086  }
1087  }
1088 
1089  if (fTriggerTypeSel != kND) {
1091  if (fGeneralHistograms) fHistEventRejection->Fill("trigTypeSel",1);
1092  return kFALSE;
1093  }
1094  }
1095 
1096  if ((fMinCent != -999) && (fMaxCent != -999)) {
1097  if (fCent<fMinCent || fCent>fMaxCent) {
1098  if (fGeneralHistograms) fHistEventRejection->Fill("Cent",1);
1099  return kFALSE;
1100  }
1101  }
1102 
1103  if (fUseAliAnaUtils) {
1104  if (!fAliAnalysisUtils)
1105  fAliAnalysisUtils = new AliAnalysisUtils();
1106  fAliAnalysisUtils->SetMinVtxContr(2);
1107  fAliAnalysisUtils->SetMaxVtxZ(999);
1108  if(fMinVz<-998.) fMinVz = -10.;
1109  if(fMaxVz>998.) fMaxVz = 10.;
1110 
1111  if (!fAliAnalysisUtils->IsVertexSelected2013pA(InputEvent())) {
1112  if (fGeneralHistograms) fHistEventRejection->Fill("VtxSel2013pA",1);
1113  return kFALSE;
1114  }
1115 
1116  if (fRejectPileup && fAliAnalysisUtils->IsPileUpEvent(InputEvent())) {
1117  if (fGeneralHistograms) fHistEventRejection->Fill("PileUp",1);
1118  return kFALSE;
1119  }
1120 
1121  if(fTklVsClusSPDCut && fAliAnalysisUtils->IsSPDClusterVsTrackletBG(InputEvent())) {
1122  if (fGeneralHistograms) fHistEventRejection->Fill("Bkg evt",1);
1123  return kFALSE;
1124  }
1125  }
1126 
1127  if ((fMinVz > -998.) && (fMaxVz < 998.)) {
1128  if (fNVertCont == 0 ) {
1129  if (fGeneralHistograms) fHistEventRejection->Fill("vertex contr.",1);
1130  return kFALSE;
1131  }
1132  Double_t vz = fVertex[2];
1133  if (vz < fMinVz || vz > fMaxVz) {
1134  if (fGeneralHistograms) fHistEventRejection->Fill("Vz",1);
1135  return kFALSE;
1136  }
1137 
1138  if (fNVertSPDCont > 0 && fZvertexDiff < 999) {
1139  Double_t vzSPD = fVertexSPD[2];
1140  Double_t dvertex = TMath::Abs(vz-vzSPD);
1141  //if difference larger than fZvertexDiff
1142  if (dvertex > fZvertexDiff) {
1143  if (fGeneralHistograms) fHistEventRejection->Fill("VzSPD",1);
1144  return kFALSE;
1145  }
1146  }
1147  }
1148 
1149  if (fMinPtTrackInEmcal > 0 && fGeom) {
1150  Bool_t trackInEmcalOk = kFALSE;
1151  Int_t ntracks = GetNParticles(0);
1152  for (Int_t i = 0; i < ntracks; i++) {
1153  AliVParticle *track = GetAcceptParticleFromArray(i,0);
1154  if (!track)
1155  continue;
1156 
1157  Double_t phiMin = fGeom->GetArm1PhiMin() * TMath::DegToRad();
1158  Double_t phiMax = fGeom->GetArm1PhiMax() * TMath::DegToRad();
1159  Int_t runNumber = InputEvent()->GetRunNumber();
1160  if (runNumber>=177295 && runNumber<=197470) { //small SM masked in 2012 and 2013
1161  phiMin = 1.4;
1162  phiMax = TMath::Pi();
1163  }
1164 
1165  if (track->Eta() < fGeom->GetArm1EtaMin() || track->Eta() > fGeom->GetArm1EtaMax() || track->Phi() < phiMin || track->Phi() > phiMax)
1166  continue;
1167  if (track->Pt() > fMinPtTrackInEmcal) {
1168  trackInEmcalOk = kTRUE;
1169  break;
1170  }
1171  }
1172  if (!trackInEmcalOk) {
1173  if (fGeneralHistograms) fHistEventRejection->Fill("trackInEmcal",1);
1174  return kFALSE;
1175  }
1176  }
1177 
1178  if (fMinNTrack > 0) {
1179  Int_t nTracksAcc = 0;
1180  Int_t ntracks = GetNParticles(0);
1181  for (Int_t i = 0; i < ntracks; i++) {
1182  AliVParticle *track = GetAcceptParticleFromArray(i,0);
1183  if (!track)
1184  continue;
1185  if (track->Pt() > fTrackPtCut) {
1186  nTracksAcc++;
1187  if (nTracksAcc>=fMinNTrack)
1188  break;
1189  }
1190  }
1191  if (nTracksAcc<fMinNTrack) {
1192  if (fGeneralHistograms) fHistEventRejection->Fill("minNTrack",1);
1193  return kFALSE;
1194  }
1195  }
1196 
1197  if (!(fEPV0 > fMinEventPlane && fEPV0 <= fMaxEventPlane) &&
1198  !(fEPV0 + TMath::Pi() > fMinEventPlane && fEPV0 + TMath::Pi() <= fMaxEventPlane) &&
1199  !(fEPV0 - TMath::Pi() > fMinEventPlane && fEPV0 - TMath::Pi() <= fMaxEventPlane))
1200  {
1201  if (fGeneralHistograms) fHistEventRejection->Fill("EvtPlane",1);
1202  return kFALSE;
1203  }
1204 
1205  if (fSelectPtHardBin != -999 && fSelectPtHardBin != fPtHardBin) {
1206  if (fGeneralHistograms) fHistEventRejection->Fill("SelPtHardBin",1);
1207  return kFALSE;
1208  }
1209 
1210  return kTRUE;
1211 }
1212 
1221 TClonesArray *AliAnalysisTaskEmcal::GetArrayFromEvent(const char *name, const char *clname)
1222 {
1223  TClonesArray *arr = 0;
1224  TString sname(name);
1225  if (!sname.IsNull()) {
1226  arr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(sname));
1227  if (!arr) {
1228  AliWarning(Form("%s: Could not retrieve array with name %s!", GetName(), name));
1229  return 0;
1230  }
1231  } else {
1232  return 0;
1233  }
1234 
1235  if (!clname)
1236  return arr;
1237 
1238  TString objname(arr->GetClass()->GetName());
1239  TClass cls(objname);
1240  if (!cls.InheritsFrom(clname)) {
1241  AliWarning(Form("%s: Objects of type %s in %s are not inherited from %s!",
1242  GetName(), cls.GetName(), name, clname));
1243  return 0;
1244  }
1245  return arr;
1246 }
1247 
1253 {
1254  fVertex[0] = 0;
1255  fVertex[1] = 0;
1256  fVertex[2] = 0;
1257  fNVertCont = 0;
1258 
1259  fVertexSPD[0] = 0;
1260  fVertexSPD[1] = 0;
1261  fVertexSPD[2] = 0;
1262  fNVertSPDCont = 0;
1263 
1264  if (fGeneratePythiaInfoObject && MCEvent()) {
1265  GeneratePythiaInfoObject(MCEvent());
1266  }
1267 
1268  const AliVVertex *vert = InputEvent()->GetPrimaryVertex();
1269  if (vert) {
1270  vert->GetXYZ(fVertex);
1271  fNVertCont = vert->GetNContributors();
1272  }
1273 
1274  const AliVVertex *vertSPD = InputEvent()->GetPrimaryVertexSPD();
1275  if (vertSPD) {
1276  vertSPD->GetXYZ(fVertexSPD);
1277  fNVertSPDCont = vertSPD->GetNContributors();
1278  }
1279 
1280  fBeamType = GetBeamType();
1281 
1282  if (fBeamType == kAA || fBeamType == kpA ) {
1284  AliMultSelection *MultSelection = static_cast<AliMultSelection*>(InputEvent()->FindListObject("MultSelection"));
1285  if (MultSelection) {
1286  fCent = MultSelection->GetMultiplicityPercentile(fCentEst.Data());
1287  }
1288  else {
1289  AliWarning(Form("%s: Could not retrieve centrality information! Assuming 99", GetName()));
1290  }
1291  }
1292  else { // old centrality estimation < 2015
1293  AliCentrality *aliCent = InputEvent()->GetCentrality();
1294  if (aliCent) {
1295  fCent = aliCent->GetCentralityPercentile(fCentEst.Data());
1296  }
1297  else {
1298  AliWarning(Form("%s: Could not retrieve centrality information! Assuming 99", GetName()));
1299  }
1300  }
1301 
1302  if (fNcentBins==4) {
1303  if (fCent >= 0 && fCent < 10) fCentBin = 0;
1304  else if (fCent >= 10 && fCent < 30) fCentBin = 1;
1305  else if (fCent >= 30 && fCent < 50) fCentBin = 2;
1306  else if (fCent >= 50 && fCent <= 100) fCentBin = 3;
1307  else {
1308  AliWarning(Form("%s: Negative centrality: %f. Assuming 99", GetName(), fCent));
1309  fCentBin = fNcentBins-1;
1310  }
1311  }
1312  else if (fNcentBins==5) { // for PbPb 2015
1313  if (fCent >= 0 && fCent < 10) fCentBin = 0;
1314  else if (fCent >= 10 && fCent < 30) fCentBin = 1;
1315  else if (fCent >= 30 && fCent < 50) fCentBin = 2;
1316  else if (fCent >= 50 && fCent <= 90) fCentBin = 3;
1317  else if (fCent > 90) {
1318  fCent = 99;
1319  fCentBin = 4;
1320  }
1321  else {
1322  AliWarning(Form("%s: Negative centrality: %f. Assuming 99", GetName(), fCent));
1323  fCentBin = fNcentBins-1;
1324  }
1325  }
1326  else {
1327  Double_t centWidth = (fMaxCent-fMinCent)/(Double_t)fNcentBins;
1328  if(centWidth>0.) {
1329  fCentBin = TMath::FloorNint(fCent/centWidth);
1330  }
1331  else {
1332  fCentBin = 0;
1333  }
1334  if (fCentBin>=fNcentBins) {
1335  AliWarning(Form("%s: fCentBin too large: cent = %f fCentBin = %d. Assuming 99", GetName(),fCent,fCentBin));
1336  fCentBin = fNcentBins-1;
1337  }
1338  }
1339 
1340  AliEventplane *aliEP = InputEvent()->GetEventplane();
1341  if (aliEP) {
1342  fEPV0 = aliEP->GetEventplane("V0" ,InputEvent());
1343  fEPV0A = aliEP->GetEventplane("V0A",InputEvent());
1344  fEPV0C = aliEP->GetEventplane("V0C",InputEvent());
1345  } else {
1346  AliWarning(Form("%s: Could not retrieve event plane information!", GetName()));
1347  }
1348  }
1349  else {
1350  fCent = 99;
1351  fCentBin = 0;
1352  }
1353 
1354  if (fIsPythia) {
1355 
1356  if (MCEvent()) {
1357  fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
1358  if (!fPythiaHeader) {
1359  // Check if AOD
1360  AliAODMCHeader* aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
1361 
1362  if (aodMCH) {
1363  for (UInt_t i = 0;i<aodMCH->GetNCocktailHeaders();i++) {
1364  fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(aodMCH->GetCocktailHeader(i));
1365  if (fPythiaHeader) break;
1366  }
1367  }
1368  }
1369  }
1370 
1371  if (fPythiaHeader) {
1372  fPtHard = fPythiaHeader->GetPtHard();
1373 
1374  const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
1375  const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
1376  for (fPtHardBin = 0; fPtHardBin < 11; fPtHardBin++) {
1377  if (fPtHard >= ptHardLo[fPtHardBin] && fPtHard < ptHardHi[fPtHardBin])
1378  break;
1379  }
1380 
1381  fXsection = fPythiaHeader->GetXsection();
1382  fNTrials = fPythiaHeader->Trials();
1383  }
1384  }
1385 
1387 
1388  AliEmcalContainer* cont = 0;
1389 
1390  TIter nextPartColl(&fParticleCollArray);
1391  while ((cont = static_cast<AliEmcalContainer*>(nextPartColl()))) cont->NextEvent();
1392 
1393  TIter nextClusColl(&fClusterCollArray);
1394  while ((cont = static_cast<AliParticleContainer*>(nextClusColl()))) cont->NextEvent();
1395 
1396  return kTRUE;
1397 }
1398 
1407 {
1408  if (TString(n).IsNull()) return 0;
1409 
1411 
1412  fParticleCollArray.Add(cont);
1413 
1414  return cont;
1415 }
1416 
1425 {
1426  if (TString(n).IsNull()) return 0;
1427 
1428  AliTrackContainer* cont = new AliTrackContainer(n);
1429 
1430  fParticleCollArray.Add(cont);
1431 
1432  return cont;
1433 }
1434 
1443 {
1444  if (TString(n).IsNull()) return 0;
1445 
1447 
1448  fParticleCollArray.Add(cont);
1449 
1450  return cont;
1451 }
1452 
1461 {
1462  if (TString(n).IsNull()) return 0;
1463 
1465 
1466  fClusterCollArray.Add(cont);
1467 
1468  return cont;
1469 }
1470 
1477 {
1478  if (i<0 || i>fParticleCollArray.GetEntriesFast()) return 0;
1479  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.At(i));
1480  return cont;
1481 }
1482 
1489 {
1490  if (i<0 || i>fClusterCollArray.GetEntriesFast()) return 0;
1491  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.At(i));
1492  return cont;
1493 }
1494 
1501 {
1502  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.FindObject(name));
1503  return cont;
1504 }
1505 
1512 {
1513  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.FindObject(name));
1514  return cont;
1515 }
1516 
1522 TClonesArray* AliAnalysisTaskEmcal::GetParticleArray(Int_t i) const
1523 {
1525  if (!cont) {
1526  AliError(Form("%s: Particle container %d not found",GetName(),i));
1527  return 0;
1528  }
1529  TString contName = cont->GetArrayName();
1530  return cont->GetArray();
1531 }
1532 
1538 TClonesArray* AliAnalysisTaskEmcal::GetClusterArray(Int_t i) const
1539 {
1541  if (!cont) {
1542  AliError(Form("%s:Cluster container %d not found",GetName(),i));
1543  return 0;
1544  }
1545  return cont->GetArray();
1546 }
1547 
1555 AliVParticle* AliAnalysisTaskEmcal::GetAcceptParticleFromArray(Int_t p, Int_t c) const
1556 {
1557 
1559  if (!cont) {
1560  AliError(Form("%s: Particle container %d not found",GetName(),c));
1561  return 0;
1562  }
1563  AliVParticle *vp = cont->GetAcceptParticle(p);
1564 
1565  return vp;
1566 }
1567 
1575 AliVCluster* AliAnalysisTaskEmcal::GetAcceptClusterFromArray(Int_t cl, Int_t c) const
1576 {
1578  if (!cont) {
1579  AliError(Form("%s: Cluster container %d not found",GetName(),c));
1580  return 0;
1581  }
1582  AliVCluster *vc = cont->GetAcceptCluster(cl);
1583 
1584  return vc;
1585 }
1586 
1593 {
1595  if (!cont) {
1596  AliError(Form("%s: Particle container %d not found",GetName(),i));
1597  return 0;
1598  }
1599  return cont->GetNEntries();
1600 }
1601 
1609 {
1611  if (!cont) {
1612  AliError(Form("%s: Cluster container %d not found",GetName(),i));
1613  return 0;
1614  }
1615  return cont->GetNEntries();
1616 }
1617 
1629 AliEMCALTriggerPatchInfo* AliAnalysisTaskEmcal::GetMainTriggerPatch(TriggerCategory trigger, Bool_t doSimpleOffline)
1630 {
1631 
1632  if (!fTriggerPatchInfo) {
1633  AliError(Form("%s: fTriggerPatchInfo not available",GetName()));
1634  return 0;
1635  }
1636 
1637  //number of patches in event
1638  Int_t nPatch = fTriggerPatchInfo->GetEntries();
1639 
1640  //extract main trigger patch(es)
1641  AliEMCALTriggerPatchInfo *patch(NULL), *selected(NULL);
1642  for (Int_t iPatch = 0; iPatch < nPatch; iPatch++) {
1643 
1644  patch = (AliEMCALTriggerPatchInfo*)fTriggerPatchInfo->At( iPatch );
1645  if (patch->IsMainTrigger()) {
1646  if(doSimpleOffline){
1647  if(patch->IsOfflineSimple()){
1648  switch(trigger){
1649  case kTriggerLevel0:
1650  // option not yet implemented in the trigger maker
1651  if(patch->IsLevel0()) selected = patch;
1652  break;
1653  case kTriggerLevel1Jet:
1654  if(patch->IsJetHighSimple() || patch->IsJetLowSimple()){
1655  if(!selected) selected = patch;
1656  else if(patch->GetADCOfflineAmp() > selected->GetADCOfflineAmp()) selected = patch;
1657  }
1658  break;
1659  case kTriggerLevel1Gamma:
1660  if(patch->IsGammaHighSimple() || patch->IsGammaLowSimple()){
1661  if(!selected) selected = patch;
1662  else if(patch->GetADCOfflineAmp() > selected->GetADCOfflineAmp()) selected = patch;
1663  }
1664  break;
1665  default: // Silence compiler warnings
1666  AliError("Untreated case: Main Patch is recalculated; should be in 'else' branch");
1667  };
1668  }
1669  } else { // Not OfflineSimple
1670  switch(trigger){
1671  case kTriggerLevel0:
1672  if(patch->IsLevel0()) selected = patch;
1673  break;
1674  case kTriggerLevel1Jet:
1675  if(patch->IsJetHigh() || patch->IsJetLow()){
1676  if(!selected) selected = patch;
1677  else if (patch->GetADCAmp() > selected->GetADCAmp())
1678  selected = patch;
1679  }
1680  break;
1681  case kTriggerLevel1Gamma:
1682  if(patch->IsGammaHigh() || patch->IsGammaLow()){
1683  if(!selected) selected = patch;
1684  else if (patch->GetADCAmp() > selected->GetADCAmp())
1685  selected = patch;
1686  }
1687  break;
1688  default:
1689  AliError("Untreated case: Main Patch is recalculated; should be in 'else' branch");
1690  };
1691  }
1692  }
1693  else if ((trigger == kTriggerRecalcJet && patch->IsRecalcJet()) ||
1694  (trigger == kTriggerRecalcGamma && patch->IsRecalcGamma())) { // recalculated patches
1695  if (doSimpleOffline && patch->IsOfflineSimple()) {
1696  if(!selected) selected = patch;
1697  else if (patch->GetADCOfflineAmp() > selected->GetADCOfflineAmp()) // this in fact should not be needed, but we have it in teh other branches as well, so keeping it for compleness
1698  selected = patch;
1699  }
1700  else if (!doSimpleOffline && !patch->IsOfflineSimple()) {
1701  if(!selected) selected = patch;
1702  else if (patch->GetADCAmp() > selected->GetADCAmp())
1703  selected = patch;
1704  }
1705  }
1706  }
1707  return selected;
1708 }
1709 
1715 void AliAnalysisTaskEmcal::AddObjectToEvent(TObject *obj, Bool_t attempt)
1716 {
1717  if (!(InputEvent()->FindListObject(obj->GetName()))) {
1718  InputEvent()->AddObject(obj);
1719  }
1720  else {
1721  if (!attempt) {
1722  AliFatal(Form("%s: Container with name %s already present. Aborting", GetName(), obj->GetName()));
1723  }
1724  }
1725 }
1726 
1734 Bool_t AliAnalysisTaskEmcal::IsTrackInEmcalAcceptance(AliVParticle* part, Double_t edges) const
1735 {
1736 
1737  if (!fGeom) {
1738  AliWarning(Form("%s - AliAnalysisTaskEmcal::IsTrackInEmcalAcceptance - Geometry is not available!", GetName()));
1739  return kFALSE;
1740  }
1741 
1742  Double_t minPhi = fGeom->GetArm1PhiMin() - edges;
1743  Double_t maxPhi = fGeom->GetArm1PhiMax() + edges;
1744 
1745  if (part->Phi() > minPhi && part->Phi() < maxPhi) {
1746  return kTRUE;
1747  }
1748  else {
1749  return kFALSE;
1750  }
1751 }
1752 
1754 {
1755  axis->SetBinLabel(1, "NullObject");
1756  axis->SetBinLabel(2, "Pt");
1757  axis->SetBinLabel(3, "Acceptance");
1758  axis->SetBinLabel(4, "MCLabel");
1759  axis->SetBinLabel(5, "BitMap");
1760  axis->SetBinLabel(6, "HF cut");
1761  axis->SetBinLabel(7, "Bit6");
1762  axis->SetBinLabel(8, "NotHybridTrack");
1763  axis->SetBinLabel(9, "MCFlag");
1764  axis->SetBinLabel(10, "MCGenerator");
1765  axis->SetBinLabel(11, "ChargeCut");
1766  axis->SetBinLabel(12, "MinDistanceTPCSectorEdge");
1767  axis->SetBinLabel(13, "Bit12");
1768  axis->SetBinLabel(14, "IsEMCal");
1769  axis->SetBinLabel(15, "Time");
1770  axis->SetBinLabel(16, "Energy");
1771  axis->SetBinLabel(17, "ExoticCut");
1772  axis->SetBinLabel(18, "Bit17");
1773  axis->SetBinLabel(19, "Area");
1774  axis->SetBinLabel(20, "AreaEmc");
1775  axis->SetBinLabel(21, "ZLeadingCh");
1776  axis->SetBinLabel(22, "ZLeadingEmc");
1777  axis->SetBinLabel(23, "NEF");
1778  axis->SetBinLabel(24, "MinLeadPt");
1779  axis->SetBinLabel(25, "MaxTrackPt");
1780  axis->SetBinLabel(26, "MaxClusterPt");
1781  axis->SetBinLabel(27, "Flavour");
1782  axis->SetBinLabel(28, "TagStatus");
1783  axis->SetBinLabel(29, "MinNConstituents");
1784  axis->SetBinLabel(30, "Bit29");
1785  axis->SetBinLabel(31, "Bit30");
1786  axis->SetBinLabel(32, "Bit31");
1787 }
1788 
1795 Double_t AliAnalysisTaskEmcal::GetParallelFraction(AliVParticle* part1, AliVParticle* part2)
1796 {
1797  TVector3 vect1(part1->Px(), part1->Py(), part1->Pz());
1798  TVector3 vect2(part2->Px(), part2->Py(), part2->Pz());
1799  Double_t z = (vect1 * vect2) / (vect2 * vect2);
1800  return z;
1801 }
1802 
1809 Double_t AliAnalysisTaskEmcal::GetParallelFraction(const TVector3& vect1, AliVParticle* part2)
1810 {
1811  TVector3 vect2(part2->Px(), part2->Py(), part2->Pz());
1812  Double_t z = (vect1 * vect2) / (vect2 * vect2);
1813  return z;
1814 }
1815 
1824 void AliAnalysisTaskEmcal::GetEtaPhiDiff(const AliVTrack *t, const AliVCluster *v, Double_t &phidiff, Double_t &etadiff)
1825 {
1826  phidiff = 999;
1827  etadiff = 999;
1828 
1829  if (!t||!v) return;
1830 
1831  Double_t veta = t->GetTrackEtaOnEMCal();
1832  Double_t vphi = t->GetTrackPhiOnEMCal();
1833 
1834  Float_t pos[3] = {0};
1835  v->GetPosition(pos);
1836  TVector3 cpos(pos);
1837  Double_t ceta = cpos.Eta();
1838  Double_t cphi = cpos.Phi();
1839  etadiff=veta-ceta;
1840  phidiff=TVector2::Phi_mpi_pi(vphi-cphi);
1841 }
1842 
1848 Byte_t AliAnalysisTaskEmcal::GetTrackType(const AliVTrack *t)
1849 {
1850  Byte_t ret = 0;
1851  if (t->TestBit(BIT(22)) && !t->TestBit(BIT(23)))
1852  ret = 1;
1853  else if (!t->TestBit(BIT(22)) && t->TestBit(BIT(23)))
1854  ret = 2;
1855  else if (t->TestBit(BIT(22)) && t->TestBit(BIT(23)))
1856  ret = 3;
1857  return ret;
1858 }
1859 
1869 Byte_t AliAnalysisTaskEmcal::GetTrackType(const AliAODTrack *aodTrack, UInt_t filterBit1, UInt_t filterBit2)
1870 {
1871 
1872  Int_t res = 0;
1873 
1874  if (aodTrack->TestFilterBit(filterBit1)) {
1875  res = 0;
1876  }
1877  else if (aodTrack->TestFilterBit(filterBit2)) {
1878  if ((aodTrack->GetStatus()&AliVTrack::kITSrefit)!=0) {
1879  res = 1;
1880  }
1881  else {
1882  res = 2;
1883  }
1884  }
1885  else {
1886  res = 3;
1887  }
1888 
1889  return res;
1890 }
1891 
1897 {
1898  if (!fPythiaInfo) {
1900  }
1901 
1902  AliStack* stack = mcEvent->Stack();
1903 
1904  const Int_t nprim = stack->GetNprimary();
1905  // reject if partons are missing from stack for some reason
1906  if (nprim < 8) return;
1907 
1908  TParticle *part6 = stack->Particle(6);
1909  TParticle *part7 = stack->Particle(7);
1910 
1911  fPythiaInfo->SetPartonFlag6(TMath::Abs(part6->GetPdgCode()));
1912  fPythiaInfo->SetParton6(part6->Pt(), part6->Eta(), part6->Phi(), part6->GetMass());
1913 
1914  fPythiaInfo->SetPartonFlag7(TMath::Abs(part7->GetPdgCode()));
1915  fPythiaInfo->SetParton7(part7->Pt(), part7->Eta(), part7->Phi(), part7->GetMass());
1916 
1917  AliGenPythiaEventHeader *pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(mcEvent->GenEventHeader());
1918  if(pythiaGenHeader){
1919  Float_t ptWeight=pythiaGenHeader->EventWeight();
1920  fPythiaInfo->SetPythiaEventWeight(ptWeight);}
1921 
1922 
1923 }
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Bool_t fGeneratePythiaInfoObject
Generate Pythia info object.
TObjArray fClusterCollArray
cluster collection array
Int_t fNVertSPDCont
!event SPD vertex number of contributors
void SetParticlePtCut(Double_t cut)
Bool_t fIsPythia
trigger, if it is a PYTHIA production
void SetParton7(Float_t pt, Float_t eta, Float_t phi, Float_t mass=0)
TH1 * fHistTrials
!trials from pyxsec.root
EMCAL Level1 gamma trigger, low threshold.
AliEmcalPythiaInfo * fPythiaInfo
!event parton info
Bool_t AcceptTrack(AliVParticle *track, Int_t c=0) const
EMCAL Level1 jet trigger, low threshold.
Bool_t HasTriggerType(TriggerType triggersel)
Int_t fNTrials
!event trials
UInt_t fOffTrigger
offline trigger for event selection
Double_t fVertexSPD[3]
!event Svertex
Double_t fMinCent
min centrality for event selection
Double_t fTrackPtCut
cut on track pt in event selection
Recalculated jet trigger patch; does not need to be above trigger threshold.
Base task in the EMCAL framework.
void SetPartonFlag7(Int_t flag7)
Container with name, TClonesArray and cuts for particles.
TSystem * gSystem
Double_t fPtHard
!event pt hard
void SetTrackPtCut(Double_t cut, Int_t c=0)
static void GetEtaPhiDiff(const AliVTrack *t, const AliVCluster *v, Double_t &phidiff, Double_t &etadiff)
Double_t fMinBinPt
min pt in histograms
Double_t fEPV0
!event plane V0
TList * list
Bool_t fGeneralHistograms
whether or not it should fill some general histograms
virtual void SetArray(AliVEvent *event)
Bool_t AcceptCluster(AliVCluster *clus, Int_t c=0) const
Int_t fCentBin
!event centrality bin
TH1 * fHistEventsAfterSel
!total number of events per pt hard bin after selection
Double_t fMinPtTrackInEmcal
min pt track in emcal
TH1 * fHistEventPlane
!event plane distribution
TList * fOutput
!output list
TH1 * fHistEvents
!total number of events per pt hard bin
void SetClusPtCut(Double_t cut, Int_t c=0)
AliClusterContainer * AddClusterContainer(const char *n)
Double_t fEPV0C
!event plane V0C
void SetParton6(Float_t pt, Float_t eta, Float_t phi, Float_t mass=0)
TH1 * fHistCentrality
!event centrality distribution
Container for particles within the EMCAL framework.
TObjArray fParticleCollArray
particle/track collection array
BeamType
Switch for the beam type.
void SetTrackEtaLimits(Double_t min, Double_t max, Int_t c=0)
AliParticleContainer * GetParticleContainer(Int_t i=0) const
TProfile * fHistXsectionAfterSel
!x section from pythia header
TriggerType
Switch for EMCAL trigger types.
EMCalTriggerMode_t fEMCalTriggerMode
EMCal trigger selection mode.
virtual Bool_t FillHistograms()
Int_t GetNParticles(Int_t i=0) const
TClonesArray * fCaloClusters
!clusters
Bool_t fUseNewCentralityEstimation
Use new centrality estimation (for 2015 data)
Bool_t IsTrackInEmcalAcceptance(AliVParticle *part, Double_t edges=0.9) const
TH1 * fHistTriggerClasses
!number of events in each trigger class
Double_t fMaxVz
max vertex for event selection
void GeneratePythiaInfoObject(AliMCEvent *mcEvent)
AliEMCALGeometry * fGeom
!emcal geometry
The overlap between low and high threshold trigger is assigned to the lower threshold only...
kRecalculated gamma trigger patch; does not need to be above trigger threshold
TString fCaloTriggerPatchInfoName
trigger patch info array name
TString fCaloTriggersName
name of calo triggers collection
AliGenPythiaEventHeader * fPythiaHeader
!event Pythia header
void SetTrackPhiLimits(Double_t min, Double_t max, Int_t c=0)
AliParticleContainer * AddParticleContainer(const char *n)
AliAnalysisUtils * fAliAnalysisUtils
!vertex selection (optional)
BeamType fForceBeamType
forced beam type
Int_t fNcentBins
how many centrality bins
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Base class for container structures within the EMCAL framework.
TriggerType fTriggerTypeSel
trigger type to select based on trigger patches
virtual Bool_t AcceptCluster(Int_t i, UInt_t &rejectionReason) const
virtual Bool_t FillGeneralHistograms()
TString fTrigClass
trigger class name for event selection
AliVCluster * GetAcceptCluster(Int_t i) const
TClonesArray * GetParticleArray(Int_t i=0) const
Double_t fMinVz
min vertex for event selection
virtual AliVParticle * GetAcceptParticle(Int_t i=-1) const
BeamType fBeamType
!event beam type
Double_t fCent
!event centrality
TClonesArray * GetArray() const
Double_t fMinEventPlane
minimum event plane value
TString fCaloCellsName
name of calo cell collection
Int_t GetNClusters(Int_t i=0) const
Int_t fNVertCont
!event vertex number of contributors
Double_t fZvertexDiff
upper limit for distance between primary and SPD vertex
virtual Bool_t AcceptParticle(const AliVParticle *vp, UInt_t &rejectionReason) const
EMCAL Level0 trigger.
EMCAL Level1 jet trigger, high threshold.
Int_t fSelectPtHardBin
select one pt hard bin for analysis
AliMCParticleContainer * AddMCParticleContainer(const char *n)
static Double_t GetParallelFraction(AliVParticle *part1, AliVParticle *part2)
virtual Bool_t RetrieveEventObjects()
Bool_t fRejectPileup
Reject pilup using function AliAnalysisUtils::IsPileUpEvent()
TProfile * fHistXsection
!x section from pyxsec.root
Bool_t PythiaInfoFromFile(const char *currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
void UserExec(Option_t *option)
void SetPartonFlag6(Int_t flag6)
AliVCaloCells * fCaloCells
!cells
const TString & GetArrayName() const
TClonesArray * GetArrayFromEvent(const char *name, const char *clname=0)
Double_t fEventPlaneVsEmcal
select events which have a certain event plane wrt the emcal
virtual Bool_t IsEventSelected()
TH1 * fHistPtHard
!pt hard distribution
void SetParticleEtaLimits(Double_t min, Double_t max)
Double_t fMaxBinPt
max pt in histograms
Int_t fPtHardBin
!event pt hard bin
TClonesArray * fTracks
!tracks
TH1 * fHistTrialsAfterSel
!total number of trials per pt hard bin after selection
void LoadPythiaInfo(AliVEvent *event)
Bool_t fIsEsd
!whether it's an ESD analysis
Double_t fVertex[3]
!event vertex
AliTrackContainer * AddTrackContainer(const char *n)
Bool_t fCreateHisto
whether or not create histograms
Store some informaion about a Pythia eventThis class is used to store some information about a Pythia...
TFile * file
TH1 * fHistEventRejection
!book keep reasons for rejecting event
Int_t nevents[nsamples]
TClonesArray * fTriggerPatchInfo
!trigger patch info array
TClonesArray * GetClusterArray(Int_t i=0) const
Double_t fEPV0A
!event plane V0A
TString fCentEst
name of V0 centrality estimator
TString fPythiaInfoName
name of pythia info object
Int_t GetRunNumber(TString)
Definition: PlotMuonQA.C:2235
Declaration of class AliEmcalPythiaInfo.
void SetClusPtCut(Double_t cut)
EMCAL Level1 gamma trigger, high threshold.
void AddObjectToEvent(TObject *obj, Bool_t attempt=kFALSE)
AliVCaloTrigger * fCaloTriggers
!calo triggers
void SetRejectionReasonLabels(TAxis *axis)
TH1 * fHistZVertex
!z vertex position
Int_t fMinNTrack
minimum nr of tracks in event with pT>fTrackPtCut
Int_t GetNEntries() const
static Byte_t GetTrackType(const AliVTrack *t)
Bool_t fUseAliAnaUtils
used for LHC13* data: z-vtx, Ncontributors, z-vtx resolution cuts
void SetClusTimeCut(Double_t min, Double_t max, Int_t c=0)
ULong_t fTriggers
list of fired triggers
Double_t fMaxEventPlane
maximum event plane value
void SetPythiaEventWeight(Float_t ptWeight)
Float_t fXsection
!x-section from pythia header
Bool_t fInitialized
whether or not the task has been already initialized
TH1 * fHistEventCount
!incoming and selected events
Double_t fMaxCent
max centrality for event selection
void SetClusTimeCut(Double_t min, Double_t max)
TriggerCategory
Online trigger categories.
void SetParticlePhiLimits(Double_t min, Double_t max)
AliVParticle * GetAcceptParticleFromArray(Int_t p, Int_t c=0) const
Container structure for EMCAL clusters.
virtual void NextEvent()
Container for MC-true particles within the EMCAL framework.
Bool_t fNeedEmcalGeom
whether or not the task needs the emcal geometry
AliVCluster * GetAcceptClusterFromArray(Int_t cl, Int_t c=0) const
Int_t fNbins
no. of pt bins
Bool_t fTklVsClusSPDCut
Apply tracklet-vs-cluster SPD cut to reject background events in pp.
AliEMCALTriggerPatchInfo * GetMainTriggerPatch(TriggerCategory triggersel=kTriggerLevel1Jet, Bool_t doOfflinSimple=kFALSE)
static Double_t fgkEMCalDCalPhiDivide
phi value used to distinguish between DCal and EMCal