AliPhysics  6a0d37d (6a0d37d)
 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 <iostream>
17 #include <memory>
18 
19 #include <TClonesArray.h>
20 #include <AliEmcalList.h>
21 #include <TObject.h>
22 #include <TObjString.h>
23 #include <TH1F.h>
24 #include <TProfile.h>
25 #include <TSystem.h>
26 #include <TFile.h>
27 #include <TChain.h>
28 #include <TKey.h>
29 
30 #include "AliAnalysisTaskEmcal.h"
31 #include "AliAnalysisUtils.h"
32 #include "AliAODEvent.h"
33 #include "AliAODMCHeader.h"
34 #include "AliAODTrack.h"
35 #include "AliAnalysisManager.h"
36 #include "AliCentrality.h"
38 #include "AliEMCALGeometry.h"
39 #include "AliEmcalPythiaInfo.h"
40 #include "AliEMCALTriggerPatchInfo.h"
41 #include "AliESDEvent.h"
42 #include "AliAODInputHandler.h"
43 #include "AliESDInputHandler.h"
44 #include "AliEventplane.h"
45 #include "AliGenPythiaEventHeader.h"
46 #include "AliGenHerwigEventHeader.h"
47 #include "AliInputEventHandler.h"
48 #include "AliLog.h"
49 #include "AliMCEvent.h"
50 #include "AliMCParticle.h"
51 #include "AliMultiInputEventHandler.h"
52 #include "AliMultSelection.h"
53 #include "AliStack.h"
54 #include "AliVCaloTrigger.h"
55 #include "AliVCluster.h"
56 #include "AliVEventHandler.h"
57 #include "AliVParticle.h"
58 
60 
64 
69  AliAnalysisTaskSE("AliAnalysisTaskEmcal"),
70  fPythiaInfoName(""),
71  fForceBeamType(kNA),
72  fGeneralHistograms(kFALSE),
73  fLocalInitialized(kFALSE),
74  fCreateHisto(kTRUE),
75  fCaloCellsName(),
76  fCaloTriggersName(),
77  fCaloTriggerPatchInfoName(),
78  fMinCent(-999),
79  fMaxCent(-999),
80  fMinVz(-999),
81  fMaxVz(999),
82  fTrackPtCut(0),
83  fMinNTrack(0),
84  fZvertexDiff(0.5),
85  fUseAliAnaUtils(kFALSE),
86  fRejectPileup(kFALSE),
87  fTklVsClusSPDCut(kFALSE),
88  fOffTrigger(AliVEvent::kAny),
89  fTrigClass(),
90  fMinBiasRefTrigger("CINT7-B-NOPF-ALLNOTRD"),
91  fTriggerTypeSel(kND),
92  fNbins(250),
93  fMinBinPt(0),
94  fMaxBinPt(250),
95  fMinPtTrackInEmcal(0),
96  fEventPlaneVsEmcal(-1),
97  fMinEventPlane(-1e6),
98  fMaxEventPlane(1e6),
99  fCentEst("V0M"),
100  fIsEmbedded(kFALSE),
101  fIsPythia(kFALSE),
102  fIsHerwig(kFALSE),
103  fSelectPtHardBin(-999),
104  fMinMCLabel(0),
105  fMCLabelShift(0),
106  fNcentBins(4),
107  fNeedEmcalGeom(kTRUE),
108  fParticleCollArray(),
109  fClusterCollArray(),
110  fTriggers(0),
111  fEMCalTriggerMode(kOverlapWithLowThreshold),
112  fUseNewCentralityEstimation(kFALSE),
113  fGeneratePythiaInfoObject(kFALSE),
114  fUsePtHardBinScaling(kFALSE),
115  fUseXsecFromHeader(kFALSE),
116  fMCRejectFilter(kFALSE),
117  fCountDownscaleCorrectedEvents(kFALSE),
118  fPtHardAndJetPtFactor(0.),
119  fPtHardAndClusterPtFactor(0.),
120  fPtHardAndTrackPtFactor(0.),
121  fRunNumber(-1),
122  fAliAnalysisUtils(nullptr),
123  fIsEsd(kFALSE),
124  fGeom(nullptr),
125  fTracks(nullptr),
126  fCaloClusters(nullptr),
127  fCaloCells(nullptr),
128  fCaloTriggers(nullptr),
129  fTriggerPatchInfo(nullptr),
130  fCent(0),
131  fCentBin(-1),
132  fEPV0(-1.0),
133  fEPV0A(-1.0),
134  fEPV0C(-1.0),
135  fNVertCont(0),
136  fNVertSPDCont(0),
137  fBeamType(kNA),
138  fPythiaHeader(nullptr),
139  fHerwigHeader(nullptr),
140  fPtHard(0),
141  fPtHardBin(0),
142  fNPtHardBins(11),
143  fPtHardBinning(),
144  fNTrials(0),
145  fXsection(0),
146  fPythiaInfo(nullptr),
147  fOutput(nullptr),
148  fHistEventCount(nullptr),
149  fHistTrialsAfterSel(nullptr),
150  fHistEventsAfterSel(nullptr),
151  fHistXsectionAfterSel(nullptr),
152  fHistTrials(nullptr),
153  fHistEvents(nullptr),
154  fHistXsection(nullptr),
155  fHistPtHard(nullptr),
156  fHistCentrality(nullptr),
157  fHistZVertex(nullptr),
158  fHistEventPlane(nullptr),
159  fHistEventRejection(nullptr),
160  fHistTriggerClasses(nullptr),
161  fHistTriggerClassesCorr(nullptr)
162 {
163  fVertex[0] = 0;
164  fVertex[1] = 0;
165  fVertex[2] = 0;
166  fVertexSPD[0] = 0;
167  fVertexSPD[1] = 0;
168  fVertexSPD[2] = 0;
169 
170  fParticleCollArray.SetOwner(kTRUE);
171  fClusterCollArray.SetOwner(kTRUE);
172 }
173 
185  AliAnalysisTaskSE(name),
186  fPythiaInfoName(""),
187  fForceBeamType(kNA),
188  fGeneralHistograms(kFALSE),
189  fLocalInitialized(kFALSE),
190  fCreateHisto(histo),
191  fCaloCellsName(),
192  fCaloTriggersName(),
193  fCaloTriggerPatchInfoName(),
194  fMinCent(-999),
195  fMaxCent(-999),
196  fMinVz(-999),
197  fMaxVz(999),
198  fTrackPtCut(0),
199  fMinNTrack(0),
200  fZvertexDiff(0.5),
201  fUseAliAnaUtils(kFALSE),
202  fRejectPileup(kFALSE),
203  fTklVsClusSPDCut(kFALSE),
204  fOffTrigger(AliVEvent::kAny),
205  fTrigClass(),
206  fMinBiasRefTrigger("CINT7-B-NOPF-ALLNOTRD"),
207  fTriggerTypeSel(kND),
208  fNbins(250),
209  fMinBinPt(0),
210  fMaxBinPt(250),
211  fMinPtTrackInEmcal(0),
212  fEventPlaneVsEmcal(-1),
213  fMinEventPlane(-1e6),
214  fMaxEventPlane(1e6),
215  fCentEst("V0M"),
216  fIsEmbedded(kFALSE),
217  fIsPythia(kFALSE),
218  fIsHerwig(kFALSE),
219  fSelectPtHardBin(-999),
220  fMinMCLabel(0),
221  fMCLabelShift(0),
222  fNcentBins(4),
223  fNeedEmcalGeom(kTRUE),
224  fParticleCollArray(),
225  fClusterCollArray(),
226  fTriggers(0),
227  fEMCalTriggerMode(kOverlapWithLowThreshold),
228  fUseNewCentralityEstimation(kFALSE),
229  fGeneratePythiaInfoObject(kFALSE),
230  fUsePtHardBinScaling(kFALSE),
231  fUseXsecFromHeader(kFALSE),
232  fMCRejectFilter(kFALSE),
233  fCountDownscaleCorrectedEvents(kFALSE),
234  fPtHardAndJetPtFactor(0.),
235  fPtHardAndClusterPtFactor(0.),
236  fPtHardAndTrackPtFactor(0.),
237  fRunNumber(-1),
238  fAliAnalysisUtils(nullptr),
239  fIsEsd(kFALSE),
240  fGeom(nullptr),
241  fTracks(nullptr),
242  fCaloClusters(nullptr),
243  fCaloCells(nullptr),
244  fCaloTriggers(nullptr),
245  fTriggerPatchInfo(nullptr),
246  fCent(0),
247  fCentBin(-1),
248  fEPV0(-1.0),
249  fEPV0A(-1.0),
250  fEPV0C(-1.0),
251  fNVertCont(0),
252  fNVertSPDCont(0),
253  fBeamType(kNA),
254  fPythiaHeader(nullptr),
255  fHerwigHeader(nullptr),
256  fPtHard(0),
257  fPtHardBin(0),
258  fNPtHardBins(11),
259  fPtHardBinning(),
260  fNTrials(0),
261  fXsection(0),
262  fPythiaInfo(0),
263  fOutput(nullptr),
264  fHistEventCount(nullptr),
265  fHistTrialsAfterSel(nullptr),
266  fHistEventsAfterSel(nullptr),
267  fHistXsectionAfterSel(nullptr),
268  fHistTrials(nullptr),
269  fHistEvents(nullptr),
270  fHistXsection(nullptr),
271  fHistPtHard(nullptr),
272  fHistCentrality(nullptr),
273  fHistZVertex(nullptr),
274  fHistEventPlane(nullptr),
275  fHistEventRejection(nullptr),
276  fHistTriggerClasses(nullptr),
277  fHistTriggerClassesCorr(nullptr)
278 {
279  fVertex[0] = 0;
280  fVertex[1] = 0;
281  fVertex[2] = 0;
282  fVertexSPD[0] = 0;
283  fVertexSPD[1] = 0;
284  fVertexSPD[2] = 0;
285  fParticleCollArray.SetOwner(kTRUE);
286  fClusterCollArray.SetOwner(kTRUE);
287 
288  if (fCreateHisto) {
289  DefineOutput(1, AliEmcalList::Class());
290  }
291 }
292 
297 {
298 }
299 
307 {
309  if (cont) cont->SetClusPtCut(cut);
310  else AliError(Form("%s in SetClusPtCut(...): container %d not found",GetName(),c));
311 }
312 
321 {
323  if (cont) cont->SetClusTimeCut(min,max);
324  else AliError(Form("%s in SetClusTimeCut(...): container %d not found",GetName(),c));
325 }
326 
334 {
336  if (cont) cont->SetParticlePtCut(cut);
337  else AliError(Form("%s in SetTrackPtCut(...): container %d not found",GetName(),c));
338 
339  fTrackPtCut = cut;
340 }
341 
350 {
352  if (cont) cont->SetParticleEtaLimits(min,max);
353  else AliError(Form("%s in SetTrackPtCut(...): container %d not found",GetName(),c));
354 }
355 
364 {
366  if (cont) cont->SetParticlePhiLimits(min,max);
367  else AliError(Form("%s in SetTrackPhiLimits(...): container %d not found",GetName(),c));
368 }
369 
390 {
391  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
392  if (mgr) {
393  AliVEventHandler *evhand = mgr->GetInputEventHandler();
394  if (evhand) {
395  if (evhand->InheritsFrom("AliESDInputHandler")) {
396  fIsEsd = kTRUE;
397  }
398  else {
399  fIsEsd = kFALSE;
400  }
401  }
402  else {
403  AliError("Event handler not found!");
404  }
405  }
406  else {
407  AliError("Analysis manager not found!");
408  }
409 
410  if (!fCreateHisto)
411  return;
412 
413  OpenFile(1);
414  fOutput = new AliEmcalList();
416  fOutput->SetOwner();
417 
418  if (fForceBeamType == kpp)
419  fNcentBins = 1;
420 
421  if (!fGeneralHistograms)
422  return;
423 
424  if (fIsPythia || fIsHerwig) {
425  fHistTrialsAfterSel = new TH1F("fHistTrialsAfterSel", "fHistTrialsAfterSel", fNPtHardBins, 0, fNPtHardBins);
426  fHistTrialsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
427  fHistTrialsAfterSel->GetYaxis()->SetTitle("trials");
429 
430  fHistEventsAfterSel = new TH1F("fHistEventsAfterSel", "fHistEventsAfterSel", fNPtHardBins, 0, fNPtHardBins);
431  fHistEventsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
432  fHistEventsAfterSel->GetYaxis()->SetTitle("total events");
434 
435  fHistXsectionAfterSel = new TProfile("fHistXsectionAfterSel", "fHistXsectionAfterSel", fNPtHardBins, 0, fNPtHardBins);
436  fHistXsectionAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
437  fHistXsectionAfterSel->GetYaxis()->SetTitle("xsection");
439 
440  fHistTrials = new TH1F("fHistTrials", "fHistTrials", fNPtHardBins, 0, fNPtHardBins);
441  fHistTrials->GetXaxis()->SetTitle("p_{T} hard bin");
442  fHistTrials->GetYaxis()->SetTitle("trials");
443  fOutput->Add(fHistTrials);
444 
445  fHistEvents = new TH1F("fHistEvents", "fHistEvents", fNPtHardBins, 0, fNPtHardBins);
446  fHistEvents->GetXaxis()->SetTitle("p_{T} hard bin");
447  fHistEvents->GetYaxis()->SetTitle("total events");
448  fOutput->Add(fHistEvents);
449 
450  fHistXsection = new TProfile("fHistXsection", "fHistXsection", fNPtHardBins, 0, fNPtHardBins);
451  fHistXsection->GetXaxis()->SetTitle("p_{T} hard bin");
452  fHistXsection->GetYaxis()->SetTitle("xsection");
453  fOutput->Add(fHistXsection);
454 
455  // Set the bin labels
456  Bool_t binningAvailable = false;
457  if(fPtHardBinning.GetSize() > 0) {
458  AliInfoStream() << "Using custom pt-hard binning" << std::endl;
459  if(fPtHardBinning.GetSize() == fNPtHardBins + 1) binningAvailable = true;
460  else AliErrorStream() << "Pt-hard binning (" << fPtHardBinning.GetSize() -1 << ") does not match the amount of bins (" << fNPtHardBins << ")" << std::endl;
461  } else {
462  // Check if we fall back to the default binning
463  if(fNPtHardBins == 11) {
464  AliInfoStream() << "11 pt-hard bins - fall back to default binning for bin labels" << std::endl;
465  const Int_t kDefaultPtHardBinning[12] = {0,5,11,21,36,57, 84,117,152,191,234,1000000};
466  fPtHardBinning.Set(12);
467  for(Int_t ib = 0; ib < 12; ib++) fPtHardBinning[ib] = kDefaultPtHardBinning[ib];
468  binningAvailable = true;
469  } else {
470  AliErrorStream() << "No default binning available for " << fNPtHardBins << " pt-hard bins - bin labels will not be set." << std::endl;
471  }
472  }
473 
474  if(binningAvailable){
475  for (Int_t i = 0; i < fNPtHardBins; i++) {
476  fHistTrialsAfterSel->GetXaxis()->SetBinLabel(i+1, Form("%d-%d",fPtHardBinning[i],fPtHardBinning[i+1]));
477  fHistEventsAfterSel->GetXaxis()->SetBinLabel(i+1, Form("%d-%d",fPtHardBinning[i],fPtHardBinning[i+1]));
478  fHistXsectionAfterSel->GetXaxis()->SetBinLabel(i+1, Form("%d-%d",fPtHardBinning[i],fPtHardBinning[i+1]));
479 
480  fHistTrials->GetXaxis()->SetBinLabel(i+1, Form("%d-%d",fPtHardBinning[i],fPtHardBinning[i+1]));
481  fHistXsection->GetXaxis()->SetBinLabel(i+1, Form("%d-%d",fPtHardBinning[i],fPtHardBinning[i+1]));
482  fHistEvents->GetXaxis()->SetBinLabel(i+1, Form("%d-%d",fPtHardBinning[i],fPtHardBinning[i+1]));
483  }
484  } else {
485  AliErrorStream() << "No suitable binning available - skipping bin labels" << std::endl;
486  }
487 
488 
489  fHistPtHard = new TH1F("fHistPtHard", "fHistPtHard", fNbins*2, fMinBinPt, fMaxBinPt*4);
490  fHistPtHard->GetXaxis()->SetTitle("p_{T,hard} (GeV/c)");
491  fHistPtHard->GetYaxis()->SetTitle("counts");
492  fOutput->Add(fHistPtHard);
493  }
494 
495  fHistZVertex = new TH1F("fHistZVertex","Z vertex position", 60, -30, 30);
496  fHistZVertex->GetXaxis()->SetTitle("z");
497  fHistZVertex->GetYaxis()->SetTitle("counts");
498  fOutput->Add(fHistZVertex);
499 
500  if (fForceBeamType != kpp) {
501  fHistCentrality = new TH1F("fHistCentrality","Event centrality distribution", 200, 0, 100);
502  fHistCentrality->GetXaxis()->SetTitle("Centrality (%)");
503  fHistCentrality->GetYaxis()->SetTitle("counts");
504  fOutput->Add(fHistCentrality);
505 
506  fHistEventPlane = new TH1F("fHistEventPlane","Event plane", 120, -TMath::Pi(), TMath::Pi());
507  fHistEventPlane->GetXaxis()->SetTitle("event plane");
508  fHistEventPlane->GetYaxis()->SetTitle("counts");
509  fOutput->Add(fHistEventPlane);
510  }
511 
512  fHistEventRejection = new TH1F("fHistEventRejection","Reasons to reject event",20,0,20);
513 #if ROOT_VERSION_CODE < ROOT_VERSION(6,4,2)
514  fHistEventRejection->SetBit(TH1::kCanRebin);
515 #else
516  fHistEventRejection->SetCanExtend(TH1::kAllAxes);
517 #endif
518  fHistEventRejection->GetXaxis()->SetBinLabel(1,"PhysSel");
519  fHistEventRejection->GetXaxis()->SetBinLabel(2,"trigger");
520  fHistEventRejection->GetXaxis()->SetBinLabel(3,"trigTypeSel");
521  fHistEventRejection->GetXaxis()->SetBinLabel(4,"Cent");
522  fHistEventRejection->GetXaxis()->SetBinLabel(5,"vertex contr.");
523  fHistEventRejection->GetXaxis()->SetBinLabel(6,"Vz");
524  fHistEventRejection->GetXaxis()->SetBinLabel(7,"VzSPD");
525  fHistEventRejection->GetXaxis()->SetBinLabel(8,"trackInEmcal");
526  fHistEventRejection->GetXaxis()->SetBinLabel(9,"minNTrack");
527  fHistEventRejection->GetXaxis()->SetBinLabel(10,"VtxSel2013pA");
528  fHistEventRejection->GetXaxis()->SetBinLabel(11,"PileUp");
529  fHistEventRejection->GetXaxis()->SetBinLabel(12,"EvtPlane");
530  fHistEventRejection->GetXaxis()->SetBinLabel(13,"SelPtHardBin");
531  fHistEventRejection->GetXaxis()->SetBinLabel(14,"Bkg evt");
532  fHistEventRejection->GetYaxis()->SetTitle("counts");
534 
535  fHistTriggerClasses = new TH1F("fHistTriggerClasses","fHistTriggerClasses",3,0,3);
536 #if ROOT_VERSION_CODE < ROOT_VERSION(6,4,2)
537  fHistTriggerClasses->SetBit(TH1::kCanRebin);
538 #else
539  fHistTriggerClasses->SetCanExtend(TH1::kAllAxes);
540 #endif
542 
544  fHistTriggerClassesCorr = new TH1F("fHistTriggerClassesCorr","fHistTriggerClassesCorr",3,0,3);
545 #if ROOT_VERSION_CODE < ROOT_VERSION(6,4,2)
546  fHistTriggerClassesCorr->SetBit(TH1::kCanRebin);
547 #else
548  fHistTriggerClassesCorr->SetCanExtend(TH1::kAllAxes);
549 #endif
551  }
552 
553  fHistEventCount = new TH1F("fHistEventCount","fHistEventCount",2,0,2);
554  fHistEventCount->GetXaxis()->SetBinLabel(1,"Accepted");
555  fHistEventCount->GetXaxis()->SetBinLabel(2,"Rejected");
556  fHistEventCount->GetYaxis()->SetTitle("counts");
557  fOutput->Add(fHistEventCount);
558 
559  PostData(1, fOutput);
560 }
561 
576 {
577  if (fIsPythia || fIsHerwig) {
581  fHistPtHard->Fill(fPtHard);
582  }
583 
584 
585  fHistZVertex->Fill(fVertex[2]);
586 
587  if (fForceBeamType != kpp) {
588  fHistCentrality->Fill(fCent);
589  fHistEventPlane->Fill(fEPV0);
590  }
591 
592  std::unique_ptr<TObjArray> triggerClasses(InputEvent()->GetFiredTriggerClasses().Tokenize(" "));
593  TObjString* triggerClass(nullptr);
594  for(auto trg : *triggerClasses){
595  triggerClass = static_cast<TObjString*>(trg);
596  fHistTriggerClasses->Fill(triggerClass->GetString(), 1);
597  }
598 
600  // downscale-corrected number of events are calculated based on the min. bias reference
601  // Formula: N_corr = N_MB * d_Trg/d_{Min_Bias}
602  if(InputEvent()->GetFiredTriggerClasses().Contains(fMinBiasRefTrigger)){
604  Double_t downscaleref = downscalefactors->GetDownscaleFactorForTriggerClass(fMinBiasRefTrigger);
605  for(auto t : downscalefactors->GetTriggerClasses()){
606  Double_t downscaletrg = downscalefactors->GetDownscaleFactorForTriggerClass(t);
607  fHistTriggerClassesCorr->Fill(t, downscaletrg/downscaleref);
608  }
609  }
610  }
611 
612  return kTRUE;
613 }
614 
635 {
636  if (!fLocalInitialized){
637  ExecOnce();
638  UserExecOnce();
639  }
640 
641  if (!fLocalInitialized)
642  return;
643 
644  if (!RetrieveEventObjects())
645  return;
646 
647  if(InputEvent()->GetRunNumber() != fRunNumber){
648  fRunNumber = InputEvent()->GetRunNumber();
651  }
652 
653  // Apply fallback for pythia cross section if needed
655  AliDebugStream(1) << "Fallback to cross section from pythia header required" << std::endl;
656  // Get the pthard bin
657  Float_t pthard = fPythiaHeader->GetPtHard();
658  int pthardbin = 0;
659  if(fPtHardBinning.GetSize()){
660  for(int ib = 0; ib < fNPtHardBins; ib++){
661  if(pthard >= fPtHardBinning[ib] && pthard < fPtHardBinning[ib+1]) {
662  pthardbin = ib;
663  break;
664  }
665  }
666  }
667  fHistXsection->Fill(pthardbin, fPythiaHeader->GetXsection());
668  }
669 
670  if (IsEventSelected()) {
671  if (fGeneralHistograms) fHistEventCount->Fill("Accepted",1);
672  }
673  else {
674  if (fGeneralHistograms) fHistEventCount->Fill("Rejected",1);
675  return;
676  }
677 
679  if (!FillGeneralHistograms())
680  return;
681  }
682 
683  if (!Run())
684  return;
685 
686  if (fCreateHisto) {
687  if (!FillHistograms())
688  return;
689  }
690 
691  if (fCreateHisto && fOutput) {
692  // information for this iteration of the UserExec in the container
693  PostData(1, fOutput);
694  }
695 }
696 
706 {
707  AliWarning("AliAnalysisTaskEmcal::AcceptCluster method is deprecated. Please use GetCusterContainer(c)->AcceptCluster(clus).");
708 
709  if (!clus) return kFALSE;
710 
712  if (!cont) {
713  AliError(Form("%s:Container %d not found",GetName(),c));
714  return 0;
715  }
716  UInt_t rejectionReason = 0;
717  return cont->AcceptCluster(clus, rejectionReason);
718 }
719 
728 Bool_t AliAnalysisTaskEmcal::AcceptTrack(AliVParticle *track, Int_t c) const
729 {
730  AliWarning("AliAnalysisTaskEmcal::AcceptTrack method is deprecated. Please use GetParticleContainer(c)->AcceptParticle(clus).");
731 
732  if (!track) return kFALSE;
733 
735  if (!cont) {
736  AliError(Form("%s:Container %d not found",GetName(),c));
737  return 0;
738  }
739 
740  UInt_t rejectionReason = 0;
741  return cont->AcceptParticle(track, rejectionReason);
742 }
743 
755 Bool_t AliAnalysisTaskEmcal::PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
756 {
757  TString file(currFile);
758  fXsec = 0;
759  fTrials = 1;
760 
761  // Determine archive type
762  TString archivetype;
763  std::unique_ptr<TObjArray> walk(file.Tokenize("/"));
764  for(auto t : *walk){
765  TString &tok = static_cast<TObjString *>(t)->String();
766  if(tok.Contains(".zip")){
767  archivetype = tok;
768  Int_t pos = archivetype.Index(".zip");
769  archivetype.Replace(pos, archivetype.Length() - pos, "");
770  }
771  }
772  if(archivetype.Length()){
773  AliDebugStream(1) << "Auto-detected archive type " << archivetype << std::endl;
774  Ssiz_t pos1 = file.Index(archivetype,archivetype.Length(),0,TString::kExact);
775  Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
776  Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
777  file.Replace(pos+1,pos2-pos1,"");
778  } else {
779  // not an archive take the basename....
780  file.ReplaceAll(gSystem->BaseName(file.Data()),"");
781  }
782  AliDebugStream(1) << "File name: " << file << std::endl;
783 
784  // Build virtual file name
785  // Support for train tests
786  TString virtualFileName;
787  if(file.Contains("__alice")){
788  TString tmp(file);
789  Int_t pos = tmp.Index("__alice");
790  tmp.Replace(0, pos, "");
791  tmp.ReplaceAll("__", "/");
792  // cut out tag for archive and root file
793  // this needs a determin
794  std::unique_ptr<TObjArray> toks(tmp.Tokenize("/"));
795  TString tag = "_" + archivetype;
796  for(auto t : *toks){
797  TString &path = static_cast<TObjString *>(t)->String();
798  if(path.Contains(tag)){
799  Int_t posTag = path.Index(tag);
800  path.Replace(posTag, path.Length() - posTag, "");
801  }
802  virtualFileName += "/" + path;
803  }
804  } else {
805  virtualFileName = file;
806  }
807 
808  AliDebugStream(1) << "Physical file name " << file << ", virtual file name " << virtualFileName << std::endl;
809 
810  // Get the pt hard bin
811  TString strPthard(virtualFileName);
812 
813  /*
814  // Dead code - to be removed after testing phase
815  // Procedure will fail for everything else than the expected path name
816  strPthard.Remove(strPthard.Last('/'));
817  strPthard.Remove(strPthard.Last('/'));
818  if (strPthard.Contains("AOD")) strPthard.Remove(strPthard.Last('/'));
819  strPthard.Remove(0,strPthard.Last('/')+1);
820  if (strPthard.IsDec()) pthard = strPthard.Atoi();
821  else
822  AliWarningStream() << "Could not extract file number from path " << strPthard << std::endl;
823  */
824 
825  // New implementation : pattern matching
826  // Reason: Implementation valid only for old productions (new productions swap run number and pt-hard bin)
827  // Idea: Don't use the position in the string but the match differnet informations
828  // + Year clearly 2000+
829  // + Run number can be match to the one in the event
830  // + If we know it is not year or run number, it must be the pt-hard bin if we start from the beginning
831  // The procedure is only valid for the current implementations and unable to detect non-pt-hard bins
832  // It will also fail in case of arbitrary file names
833 
834  bool binfound = false;
835  std::unique_ptr<TObjArray> tokens(strPthard.Tokenize("/"));
836  for(auto t : *tokens) {
837  TString &tok = static_cast<TObjString *>(t)->String();
838  if(tok.IsDec()){
839  Int_t number = tok.Atoi();
840  if(number > 2000 && number < 3000){
841  // Year
842  continue;
843  } else if(number == fInputEvent->GetRunNumber()){
844  // Run number
845  } else {
846  if(!binfound){
847  // the first number that is not one of the two must be the pt-hard bin
848  binfound = true;
849  pthard = number;
850  break;
851  }
852  }
853  }
854  }
855  if(!binfound) {
856  AliErrorStream() << "Could not extract file number from path " << strPthard << std::endl;
857  } else {
858  AliInfoStream() << "Auto-detecting pt-hard bin " << pthard << std::endl;
859  }
860 
861  // 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
862  std::unique_ptr<TFile> fxsec(TFile::Open(Form("%s%s",file.Data(),"pyxsec.root")));
863 
864  if (!fxsec) {
865  // next trial fetch the histgram file
866  fxsec = std::unique_ptr<TFile>(TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root")));
867  if (!fxsec){
868  fUseXsecFromHeader = true;
869  return kFALSE; // not a severe condition but inciate that we have no information
870  }
871  else {
872  // find the tlist we want to be independtent of the name so use the Tkey
873  TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
874  if (!key) return kFALSE;
875  TList *list = dynamic_cast<TList*>(key->ReadObj());
876  if (!list) return kFALSE;
877  TProfile *xSecHist = static_cast<TProfile*>(list->FindObject("h1Xsec"));
878  // check for failure
879  if(!xSecHist->GetEntries()) {
880  // No cross seciton information available - fall back to raw
881  AliErrorStream() << "No cross section information available in file " << fxsec->GetName() <<" - fall back to cross section in PYTHIA header" << std::endl;
882  fUseXsecFromHeader = true;
883  } else {
884  // Cross section histogram filled - take it from there
885  fXsec = xSecHist->GetBinContent(1);
886  fUseXsecFromHeader = false;
887  }
888  fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
889  }
890  } else { // no tree pyxsec.root
891  TTree *xtree = (TTree*)fxsec->Get("Xsection");
892  if (!xtree) return kFALSE;
893  UInt_t ntrials = 0;
894  Double_t xsection = 0;
895  xtree->SetBranchAddress("xsection",&xsection);
896  xtree->SetBranchAddress("ntrials",&ntrials);
897  xtree->GetEntry(0);
898  fTrials = ntrials;
899  fXsec = xsection;
900  }
901  return kTRUE;
902 }
903 
918 {
920  return kTRUE;
921 
922  TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
923  if (!tree) {
924  AliError(Form("%s - UserNotify: No current tree!",GetName()));
925  return kFALSE;
926  }
927 
928  Float_t xsection = 0;
929  Float_t trials = 0;
930  Int_t pthardbin = 0;
931 
932  TFile *curfile = tree->GetCurrentFile();
933  if (!curfile) {
934  AliError(Form("%s - UserNotify: No current file!",GetName()));
935  return kFALSE;
936  }
937 
938  TChain *chain = dynamic_cast<TChain*>(tree);
939  if (chain) tree = chain->GetTree();
940 
941  Int_t nevents = tree->GetEntriesFast();
942 
943  PythiaInfoFromFile(curfile->GetName(), xsection, trials, pthardbin);
944 
945  if ((pthardbin < 0) || (pthardbin > fNPtHardBins-1)) pthardbin = 0;
946  fHistTrials->Fill(pthardbin, trials);
948  AliDebugStream(1) << "Using cross section from file pyxsec.root" << std::endl;
949  fHistXsection->Fill(pthardbin, xsection);
950  }
951  fHistEvents->Fill(pthardbin, nevents);
952 
953  return kTRUE;
954 }
955 
961 {
962  if (!fPythiaInfoName.IsNull() && !fPythiaInfo) {
963  fPythiaInfo = dynamic_cast<AliEmcalPythiaInfo*>(event->FindListObject(fPythiaInfoName));
964  if (!fPythiaInfo) {
965  AliError(Form("%s: Could not retrieve parton infos! %s!", GetName(), fPythiaInfoName.Data()));
966  return;
967  }
968  }
969 }
970 
982 {
983  if (!InputEvent()) {
984  AliError(Form("%s: Could not retrieve event! Returning!", GetName()));
985  return;
986  }
987 
988  LoadPythiaInfo(InputEvent());
989 
990  if (fNeedEmcalGeom) {
991  fGeom = AliEMCALGeometry::GetInstanceFromRunNumber(InputEvent()->GetRunNumber());
992  if (!fGeom) {
993  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()));
994  return;
995  }
996  }
997 
998  if (fEventPlaneVsEmcal >= 0) {
999  if (fGeom) {
1000  Double_t ep = (fGeom->GetArm1PhiMax() + fGeom->GetArm1PhiMin()) / 2 * TMath::DegToRad() + fEventPlaneVsEmcal - TMath::Pi();
1001  fMinEventPlane = ep - TMath::Pi() / 4;
1002  fMaxEventPlane = ep + TMath::Pi() / 4;
1003  }
1004  else {
1005  AliWarning("Could not set event plane limits because EMCal geometry was not loaded!");
1006  }
1007  }
1008 
1009  //Load all requested track branches - each container knows name already
1010  for (Int_t i =0; i<fParticleCollArray.GetEntriesFast(); i++) {
1011  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.At(i));
1012  cont->SetArray(InputEvent());
1013  }
1014 
1015  if (fParticleCollArray.GetEntriesFast()>0) {
1017  if (!fTracks) {
1018  AliError(Form("%s: Could not retrieve first track branch!", GetName()));
1019  return;
1020  }
1021  }
1022 
1023  //Load all requested cluster branches - each container knows name already
1024  for (Int_t i =0; i<fClusterCollArray.GetEntriesFast(); i++) {
1025  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.At(i));
1026  cont->SetArray(InputEvent());
1027  }
1028 
1029  if (fClusterCollArray.GetEntriesFast()>0) {
1031  if (!fCaloClusters) {
1032  AliError(Form("%s: Could not retrieve first cluster branch!", GetName()));
1033  return;
1034  }
1035  }
1036 
1037  if (!fCaloCellsName.IsNull() && !fCaloCells) {
1038  fCaloCells = dynamic_cast<AliVCaloCells*>(InputEvent()->FindListObject(fCaloCellsName));
1039  if (!fCaloCells) {
1040  AliError(Form("%s: Could not retrieve cells %s!", GetName(), fCaloCellsName.Data()));
1041  return;
1042  }
1043  }
1044 
1045  if (!fCaloTriggersName.IsNull() && !fCaloTriggers) {
1046  fCaloTriggers = dynamic_cast<AliVCaloTrigger*>(InputEvent()->FindListObject(fCaloTriggersName));
1047  if (!fCaloTriggers) {
1048  AliError(Form("%s: Could not retrieve calo triggers %s!", GetName(), fCaloTriggersName.Data()));
1049  return;
1050  }
1051  }
1052 
1053  if (!fCaloTriggerPatchInfoName.IsNull() && !fTriggerPatchInfo) {
1054  fTriggerPatchInfo = GetArrayFromEvent(fCaloTriggerPatchInfoName.Data(),"AliEMCALTriggerPatchInfo");
1055  if (!fTriggerPatchInfo) {
1056  AliError(Form("%s: Could not retrieve calo trigger patch info %s!", GetName(), fCaloTriggerPatchInfoName.Data()));
1057  return;
1058  }
1059 
1060  }
1061 
1062  fLocalInitialized = kTRUE;
1063 }
1064 
1071 {
1072  if (fForceBeamType != kNA)
1073  return fForceBeamType;
1074 
1075  AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
1076  if (esd) {
1077  const AliESDRun *run = esd->GetESDRun();
1078  TString beamType = run->GetBeamType();
1079  if (beamType == "p-p")
1080  return kpp;
1081  else if (beamType == "A-A")
1082  return kAA;
1083  else if (beamType == "p-A")
1084  return kpA;
1085  else
1086  return kNA;
1087  } else {
1088  Int_t runNumber = InputEvent()->GetRunNumber();
1089  // All run number ranges taken from the RCT
1090  if ((runNumber >= 136833 && runNumber <= 139517) || // LHC10h
1091  (runNumber >= 167693 && runNumber <= 170593) || // LHC11h
1092  (runNumber >= 244824 && runNumber <= 246994)) { // LHC15o
1093  return kAA;
1094  } else if ((runNumber >= 188356 && runNumber <= 188366) || // LHC12g
1095  (runNumber >= 195164 && runNumber <= 197388) || // LHC13b-f
1096  (runNumber >= 265015 && runNumber <= 267166)) { // LHC16q-t
1097  return kpA;
1098  } else {
1099  return kpp;
1100  }
1101  }
1102 }
1103 
1111 {
1112  if (!fTriggerPatchInfo)
1113  return 0;
1114 
1115  //number of patches in event
1116  Int_t nPatch = fTriggerPatchInfo->GetEntries();
1117 
1118  //loop over patches to define trigger type of event
1119  Int_t nG1 = 0;
1120  Int_t nG2 = 0;
1121  Int_t nJ1 = 0;
1122  Int_t nJ2 = 0;
1123  Int_t nL0 = 0;
1124  AliEMCALTriggerPatchInfo *patch;
1125  for (Int_t iPatch = 0; iPatch < nPatch; iPatch++) {
1126  patch = (AliEMCALTriggerPatchInfo*)fTriggerPatchInfo->At( iPatch );
1127  if (patch->IsGammaHigh()) nG1++;
1128  if (patch->IsGammaLow()) nG2++;
1129  if (patch->IsJetHigh()) nJ1++;
1130  if (patch->IsJetLow()) nJ2++;
1131  if (patch->IsLevel0()) nL0++;
1132  }
1133 
1134  AliDebug(2, "Patch summary: ");
1135  AliDebug(2, Form("Number of patches: %d", nPatch));
1136  AliDebug(2, Form("Jet: low[%d], high[%d]" ,nJ2, nJ1));
1137  AliDebug(2, Form("Gamma: low[%d], high[%d]" ,nG2, nG1));
1138 
1139  ULong_t triggers(0);
1140  if (nL0>0) SETBIT(triggers, kL0);
1141  if (nG1>0) SETBIT(triggers, kG1);
1142  if (nG2>0) SETBIT(triggers, kG2);
1143  if (nJ1>0) SETBIT(triggers, kJ1);
1144  if (nJ2>0) SETBIT(triggers, kJ2);
1145  return triggers;
1146 }
1147 
1155 {
1156  //
1157  if(trigger==kND) {
1158  AliWarning(Form("%s: Requesting undefined trigger type!", GetName()));
1159  return kFALSE;
1160  }
1161  //MV: removing this logic which as far as I can see doesn't make any sense
1162  // if(trigger & kND){
1163  // return fTriggers == 0;
1164  // }
1165  return TESTBIT(fTriggers, trigger);
1166 }
1167 
1190 {
1191  if (fOffTrigger != AliVEvent::kAny) {
1192  UInt_t res = 0;
1193  const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
1194  if (eev) {
1195  res = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
1196  } else {
1197  const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
1198  if (aev) {
1199  res = ((AliVAODHeader*)aev->GetHeader())->GetOfflineTrigger();
1200  }
1201  }
1202  if ((res & fOffTrigger) == 0) {
1203  if (fGeneralHistograms) fHistEventRejection->Fill("PhysSel",1);
1204  return kFALSE;
1205  }
1206  }
1207 
1208  if (!fTrigClass.IsNull()) {
1209  TString fired;
1210  const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
1211  if (eev) {
1212  fired = eev->GetFiredTriggerClasses();
1213  } else {
1214  const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
1215  if (aev) {
1216  fired = aev->GetFiredTriggerClasses();
1217  }
1218  }
1219  if (!fired.Contains("-B-")) {
1220  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
1221  return kFALSE;
1222  }
1223 
1224  std::unique_ptr<TObjArray> arr(fTrigClass.Tokenize("|"));
1225  if (!arr) {
1226  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
1227  return kFALSE;
1228  }
1229  Bool_t match = 0;
1230  for (Int_t i=0;i<arr->GetEntriesFast();++i) {
1231  TObject *obj = arr->At(i);
1232  if (!obj)
1233  continue;
1234 
1235  //Check if requested trigger was fired
1236  TString objStr = obj->GetName();
1238  (objStr.Contains("J1") || objStr.Contains("J2") || objStr.Contains("G1") || objStr.Contains("G2"))) {
1239  // This is relevant for EMCal triggers with 2 thresholds
1240  // If the kOverlapWithLowThreshold was requested than the overlap between the two triggers goes with the lower threshold trigger
1241  TString trigType1 = "J1";
1242  TString trigType2 = "J2";
1243  if(objStr.Contains("G")) {
1244  trigType1 = "G1";
1245  trigType2 = "G2";
1246  }
1247  if(objStr.Contains(trigType2) && fired.Contains(trigType2.Data())) { //requesting low threshold + overlap
1248  match = 1;
1249  break;
1250  }
1251  else if(objStr.Contains(trigType1) && fired.Contains(trigType1.Data()) && !fired.Contains(trigType2.Data())) { //high threshold only
1252  match = 1;
1253  break;
1254  }
1255  }
1256  else {
1257  // If this is not an EMCal trigger, or no particular treatment of EMCal triggers was requested,
1258  // simply check that the trigger was fired
1259  if (fired.Contains(obj->GetName())) {
1260  match = 1;
1261  break;
1262  }
1263  }
1264  }
1265  if (!match) {
1266  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
1267  return kFALSE;
1268  }
1269  }
1270 
1271  if (fTriggerTypeSel != kND) {
1273  if (fGeneralHistograms) fHistEventRejection->Fill("trigTypeSel",1);
1274  return kFALSE;
1275  }
1276  }
1277 
1278  if ((fMinCent != -999) && (fMaxCent != -999)) {
1279  if (fCent<fMinCent || fCent>fMaxCent) {
1280  if (fGeneralHistograms) fHistEventRejection->Fill("Cent",1);
1281  return kFALSE;
1282  }
1283  }
1284 
1285  if (fUseAliAnaUtils) {
1286  if (!fAliAnalysisUtils)
1287  fAliAnalysisUtils = new AliAnalysisUtils();
1288  fAliAnalysisUtils->SetMinVtxContr(2);
1289  fAliAnalysisUtils->SetMaxVtxZ(999);
1290  if(fMinVz<-998.) fMinVz = -10.;
1291  if(fMaxVz>998.) fMaxVz = 10.;
1292 
1293  if (!fAliAnalysisUtils->IsVertexSelected2013pA(InputEvent())) {
1294  if (fGeneralHistograms) fHistEventRejection->Fill("VtxSel2013pA",1);
1295  return kFALSE;
1296  }
1297 
1298  if (fRejectPileup && fAliAnalysisUtils->IsPileUpEvent(InputEvent())) {
1299  if (fGeneralHistograms) fHistEventRejection->Fill("PileUp",1);
1300  return kFALSE;
1301  }
1302 
1303  if(fTklVsClusSPDCut && fAliAnalysisUtils->IsSPDClusterVsTrackletBG(InputEvent())) {
1304  if (fGeneralHistograms) fHistEventRejection->Fill("Bkg evt",1);
1305  return kFALSE;
1306  }
1307  }
1308 
1309  if ((fMinVz > -998.) && (fMaxVz < 998.)) {
1310  if (fNVertCont == 0 ) {
1311  if (fGeneralHistograms) fHistEventRejection->Fill("vertex contr.",1);
1312  return kFALSE;
1313  }
1314  Double_t vz = fVertex[2];
1315  if (vz < fMinVz || vz > fMaxVz) {
1316  if (fGeneralHistograms) fHistEventRejection->Fill("Vz",1);
1317  return kFALSE;
1318  }
1319 
1320  if (fNVertSPDCont > 0 && fZvertexDiff < 999) {
1321  Double_t vzSPD = fVertexSPD[2];
1322  Double_t dvertex = TMath::Abs(vz-vzSPD);
1323  //if difference larger than fZvertexDiff
1324  if (dvertex > fZvertexDiff) {
1325  if (fGeneralHistograms) fHistEventRejection->Fill("VzSPD",1);
1326  return kFALSE;
1327  }
1328  }
1329  }
1330 
1331  if (fMinPtTrackInEmcal > 0 && fGeom) {
1332  Bool_t trackInEmcalOk = kFALSE;
1333  Int_t ntracks = GetNParticles(0);
1334  for (Int_t i = 0; i < ntracks; i++) {
1335  AliVParticle *track = GetAcceptParticleFromArray(i,0);
1336  if (!track)
1337  continue;
1338 
1339  Double_t phiMin = fGeom->GetArm1PhiMin() * TMath::DegToRad();
1340  Double_t phiMax = fGeom->GetArm1PhiMax() * TMath::DegToRad();
1341  Int_t runNumber = InputEvent()->GetRunNumber();
1342  if (runNumber>=177295 && runNumber<=197470) { //small SM masked in 2012 and 2013
1343  phiMin = 1.4;
1344  phiMax = TMath::Pi();
1345  }
1346 
1347  if (track->Eta() < fGeom->GetArm1EtaMin() || track->Eta() > fGeom->GetArm1EtaMax() || track->Phi() < phiMin || track->Phi() > phiMax)
1348  continue;
1349  if (track->Pt() > fMinPtTrackInEmcal) {
1350  trackInEmcalOk = kTRUE;
1351  break;
1352  }
1353  }
1354  if (!trackInEmcalOk) {
1355  if (fGeneralHistograms) fHistEventRejection->Fill("trackInEmcal",1);
1356  return kFALSE;
1357  }
1358  }
1359 
1360  if (fMinNTrack > 0) {
1361  Int_t nTracksAcc = 0;
1362  Int_t ntracks = GetNParticles(0);
1363  for (Int_t i = 0; i < ntracks; i++) {
1364  AliVParticle *track = GetAcceptParticleFromArray(i,0);
1365  if (!track)
1366  continue;
1367  if (track->Pt() > fTrackPtCut) {
1368  nTracksAcc++;
1369  if (nTracksAcc>=fMinNTrack)
1370  break;
1371  }
1372  }
1373  if (nTracksAcc<fMinNTrack) {
1374  if (fGeneralHistograms) fHistEventRejection->Fill("minNTrack",1);
1375  return kFALSE;
1376  }
1377  }
1378 
1379  if (!(fEPV0 > fMinEventPlane && fEPV0 <= fMaxEventPlane) &&
1380  !(fEPV0 + TMath::Pi() > fMinEventPlane && fEPV0 + TMath::Pi() <= fMaxEventPlane) &&
1381  !(fEPV0 - TMath::Pi() > fMinEventPlane && fEPV0 - TMath::Pi() <= fMaxEventPlane))
1382  {
1383  if (fGeneralHistograms) fHistEventRejection->Fill("EvtPlane",1);
1384  return kFALSE;
1385  }
1386 
1387  if (fSelectPtHardBin != -999 && fSelectPtHardBin != fPtHardBin) {
1388  if (fGeneralHistograms) fHistEventRejection->Fill("SelPtHardBin",1);
1389  return kFALSE;
1390  }
1391 
1392  // Reject filter for MC data
1393  if (!CheckMCOutliers()) return kFALSE;
1394 
1395  return kTRUE;
1396 }
1397 
1404 {
1405  if (!fPythiaHeader || !fMCRejectFilter) return kTRUE;
1406 
1407  // Condition 1: Pythia jet / pT-hard > factor
1408  if (fPtHardAndJetPtFactor > 0.) {
1409  AliTLorentzVector jet;
1410 
1411  Int_t nTriggerJets = fPythiaHeader->NTriggerJets();
1412 
1413  AliDebug(1,Form("Njets: %d, pT Hard %f",nTriggerJets, fPtHard));
1414 
1415  Float_t tmpjet[]={0,0,0,0};
1416  for (Int_t ijet = 0; ijet< nTriggerJets; ijet++) {
1417  fPythiaHeader->TriggerJet(ijet, tmpjet);
1418 
1419  jet.SetPxPyPzE(tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3]);
1420 
1421  AliDebug(1,Form("jet %d; pycell jet pT %f",ijet, jet.Pt()));
1422 
1423  //Compare jet pT and pt Hard
1424  if (jet.Pt() > fPtHardAndJetPtFactor * fPtHard) {
1425  AliInfo(Form("Reject jet event with : pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n", fPtHard, jet.Pt(), fPtHardAndJetPtFactor));
1426  return kFALSE;
1427  }
1428  }
1429  }
1430  // end condition 1
1431 
1432  // Condition 2 : Reconstructed EMCal cluster pT / pT-hard > factor
1433  if (fPtHardAndClusterPtFactor > 0.) {
1434  AliClusterContainer* mccluscont = GetClusterContainer(0);
1435  if ((Bool_t)mccluscont) {
1436  for (auto cluster : mccluscont->all()) {// Not cuts applied ; use accept for cuts
1437  Float_t ecluster = cluster->E();
1438 
1439  if (ecluster > (fPtHardAndClusterPtFactor * fPtHard)) {
1440  AliInfo(Form("Reject : ecluster %2.2f, calo %d, factor %2.2f, ptHard %f",ecluster,cluster->GetType(),fPtHardAndClusterPtFactor,fPtHard));
1441  return kFALSE;
1442  }
1443  }
1444  }
1445  }
1446  // end condition 2
1447 
1448  // condition 3 : Reconstructed track pT / pT-hard >factor
1449  if (fPtHardAndTrackPtFactor > 0.) {
1450  AliMCParticleContainer* mcpartcont = dynamic_cast<AliMCParticleContainer*>(GetParticleContainer(0));
1451  if ((Bool_t)mcpartcont) {
1452  for (auto mctrack : mcpartcont->all()) {// Not cuts applied ; use accept for cuts
1453  Float_t trackpt = mctrack->Pt();
1454  if (trackpt > (fPtHardAndTrackPtFactor * fPtHard) ) {
1455  AliInfo(Form("Reject : track %2.2f, factor %2.2f, ptHard %f", trackpt, fPtHardAndTrackPtFactor, fPtHard));
1456  return kFALSE;
1457  }
1458  }
1459  }
1460  }
1461  // end condition 3
1462 
1463  return kTRUE;
1464 }
1465 
1474 TClonesArray *AliAnalysisTaskEmcal::GetArrayFromEvent(const char *name, const char *clname)
1475 {
1476  TClonesArray *arr = 0;
1477  TString sname(name);
1478  if (!sname.IsNull()) {
1479  arr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(sname));
1480  if (!arr) {
1481  AliWarning(Form("%s: Could not retrieve array with name %s!", GetName(), name));
1482  return 0;
1483  }
1484  } else {
1485  return 0;
1486  }
1487 
1488  if (!clname)
1489  return arr;
1490 
1491  TString objname(arr->GetClass()->GetName());
1492  TClass cls(objname);
1493  if (!cls.InheritsFrom(clname)) {
1494  AliWarning(Form("%s: Objects of type %s in %s are not inherited from %s!",
1495  GetName(), cls.GetName(), name, clname));
1496  return 0;
1497  }
1498  return arr;
1499 }
1500 
1506 {
1507  fVertex[0] = 0;
1508  fVertex[1] = 0;
1509  fVertex[2] = 0;
1510  fNVertCont = 0;
1511 
1512  fVertexSPD[0] = 0;
1513  fVertexSPD[1] = 0;
1514  fVertexSPD[2] = 0;
1515  fNVertSPDCont = 0;
1516 
1517  if (fGeneratePythiaInfoObject && MCEvent()) {
1518  GeneratePythiaInfoObject(MCEvent());
1519  }
1520 
1521  const AliVVertex *vert = InputEvent()->GetPrimaryVertex();
1522  if (vert) {
1523  vert->GetXYZ(fVertex);
1524  fNVertCont = vert->GetNContributors();
1525  }
1526 
1527  const AliVVertex *vertSPD = InputEvent()->GetPrimaryVertexSPD();
1528  if (vertSPD) {
1529  vertSPD->GetXYZ(fVertexSPD);
1530  fNVertSPDCont = vertSPD->GetNContributors();
1531  }
1532 
1533  fBeamType = GetBeamType();
1534 
1535  if (fBeamType == kAA || fBeamType == kpA ) {
1537  AliMultSelection *MultSelection = static_cast<AliMultSelection*>(InputEvent()->FindListObject("MultSelection"));
1538  if (MultSelection) {
1539  fCent = MultSelection->GetMultiplicityPercentile(fCentEst.Data());
1540  }
1541  else {
1542  AliWarning(Form("%s: Could not retrieve centrality information! Assuming 99", GetName()));
1543  }
1544  }
1545  else { // old centrality estimation < 2015
1546  AliCentrality *aliCent = InputEvent()->GetCentrality();
1547  if (aliCent) {
1548  fCent = aliCent->GetCentralityPercentile(fCentEst.Data());
1549  }
1550  else {
1551  AliWarning(Form("%s: Could not retrieve centrality information! Assuming 99", GetName()));
1552  }
1553  }
1554 
1555  if (fNcentBins==4) {
1556  if (fCent >= 0 && fCent < 10) fCentBin = 0;
1557  else if (fCent >= 10 && fCent < 30) fCentBin = 1;
1558  else if (fCent >= 30 && fCent < 50) fCentBin = 2;
1559  else if (fCent >= 50 && fCent <= 100) fCentBin = 3;
1560  else {
1561  AliWarning(Form("%s: Negative centrality: %f. Assuming 99", GetName(), fCent));
1562  fCentBin = fNcentBins-1;
1563  }
1564  }
1565  else if (fNcentBins==5) { // for PbPb 2015
1566  if (fCent >= 0 && fCent < 10) fCentBin = 0;
1567  else if (fCent >= 10 && fCent < 30) fCentBin = 1;
1568  else if (fCent >= 30 && fCent < 50) fCentBin = 2;
1569  else if (fCent >= 50 && fCent <= 90) fCentBin = 3;
1570  else if (fCent > 90) {
1571  fCent = 99;
1572  fCentBin = 4;
1573  }
1574  else {
1575  AliWarning(Form("%s: Negative centrality: %f. Assuming 99", GetName(), fCent));
1576  fCentBin = fNcentBins-1;
1577  }
1578  }
1579  else {
1580  Double_t centWidth = (fMaxCent-fMinCent)/(Double_t)fNcentBins;
1581  if(centWidth>0.) {
1582  fCentBin = TMath::FloorNint(fCent/centWidth);
1583  }
1584  else {
1585  fCentBin = 0;
1586  }
1587  if (fCentBin>=fNcentBins) {
1588  AliWarning(Form("%s: fCentBin too large: cent = %f fCentBin = %d. Assuming 99", GetName(),fCent,fCentBin));
1589  fCentBin = fNcentBins-1;
1590  }
1591  }
1592 
1593  AliEventplane *aliEP = InputEvent()->GetEventplane();
1594  if (aliEP) {
1595  fEPV0 = aliEP->GetEventplane("V0" ,InputEvent());
1596  fEPV0A = aliEP->GetEventplane("V0A",InputEvent());
1597  fEPV0C = aliEP->GetEventplane("V0C",InputEvent());
1598  } else {
1599  AliWarning(Form("%s: Could not retrieve event plane information!", GetName()));
1600  }
1601  }
1602  else {
1603  fCent = 99;
1604  fCentBin = 0;
1605  }
1606 
1607  if (fIsPythia) {
1608  if (MCEvent()) {
1609  fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
1610  if (!fPythiaHeader) {
1611  // Check if AOD
1612  AliAODMCHeader* aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
1613 
1614  if (aodMCH) {
1615  for (UInt_t i = 0;i<aodMCH->GetNCocktailHeaders();i++) {
1616  fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(aodMCH->GetCocktailHeader(i));
1617  if (fPythiaHeader) break;
1618  }
1619  }
1620  }
1621  }
1622  }
1623 
1624  if (fPythiaHeader) {
1625  fPtHard = fPythiaHeader->GetPtHard();
1626 
1627  if(fPtHardBinning.GetSize()){
1628  // pt-hard binning defined for the corresponding dataset - automatically determine the bin
1629  for (fPtHardBin = 0; fPtHardBin < fNPtHardBins; fPtHardBin++) {
1631  break;
1632  }
1633  } else {
1634  // No pt-hard binning defined for the dataset - leaving the bin to 0
1635  fPtHardBin = 0;
1636  }
1637 
1638  fXsection = fPythiaHeader->GetXsection();
1639  fNTrials = fPythiaHeader->Trials();
1640  }
1641 
1642  if (fIsHerwig) {
1643  if (MCEvent()) {
1644  fHerwigHeader = dynamic_cast<AliGenHerwigEventHeader*>(MCEvent()->GenEventHeader());
1645 
1646  if (!fHerwigHeader) {
1647  // Check if AOD
1648  AliAODMCHeader* aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
1649 
1650  if (aodMCH) {
1651  for (UInt_t i = 0;i<aodMCH->GetNCocktailHeaders();i++) {
1652  fHerwigHeader = dynamic_cast<AliGenHerwigEventHeader*>(aodMCH->GetCocktailHeader(i));
1653  if (fHerwigHeader) break;
1654  }
1655  }
1656  }
1657  }
1658  }
1659 
1660  if (fHerwigHeader) {
1661  fPtHard = fHerwigHeader->GetPtHard();
1662 
1663  if(fPtHardBinning.GetSize()){
1664  // pt-hard binning defined for the corresponding dataset - automatically determine the bin
1665  for (fPtHardBin = 0; fPtHardBin < fNPtHardBins; fPtHardBin++) {
1667  break;
1668  }
1669  } else {
1670  // No pt-hard binning defined for the dataset - leaving the bin to 0
1671  fPtHardBin = 0;
1672  }
1673  fXsection = fHerwigHeader->Weight();
1674  fNTrials = fHerwigHeader->Trials();
1675  }
1676 
1677 
1679 
1680  AliEmcalContainer* cont = 0;
1681 
1682  TIter nextPartColl(&fParticleCollArray);
1683  while ((cont = static_cast<AliEmcalContainer*>(nextPartColl()))) cont->NextEvent();
1684 
1685  TIter nextClusColl(&fClusterCollArray);
1686  while ((cont = static_cast<AliParticleContainer*>(nextClusColl()))) cont->NextEvent();
1687 
1688  return kTRUE;
1689 }
1690 
1699 {
1700  if (TString(n).IsNull()) return 0;
1701 
1703 
1704  fParticleCollArray.Add(cont);
1705 
1706  return cont;
1707 }
1708 
1717 {
1718  if (TString(n).IsNull()) return 0;
1719 
1720  AliTrackContainer* cont = new AliTrackContainer(n);
1721 
1722  fParticleCollArray.Add(cont);
1723 
1724  return cont;
1725 }
1726 
1735 {
1736  if (TString(n).IsNull()) return 0;
1737 
1739 
1740  fParticleCollArray.Add(cont);
1741 
1742  return cont;
1743 }
1744 
1753 {
1754  if (TString(n).IsNull()) return 0;
1755 
1757 
1758  fClusterCollArray.Add(cont);
1759 
1760  return cont;
1761 }
1762 
1769 {
1770  if (i<0 || i>fParticleCollArray.GetEntriesFast()) return 0;
1771  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.At(i));
1772  return cont;
1773 }
1774 
1781 {
1782  if (i<0 || i>fClusterCollArray.GetEntriesFast()) return 0;
1783  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.At(i));
1784  return cont;
1785 }
1786 
1793 {
1794  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.FindObject(name));
1795  return cont;
1796 }
1797 
1804 {
1805  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.FindObject(name));
1806  return cont;
1807 }
1808 
1815 {
1817  if (!cont) {
1818  AliError(Form("%s: Particle container %d not found",GetName(),i));
1819  return 0;
1820  }
1821  TString contName = cont->GetArrayName();
1822  return cont->GetArray();
1823 }
1824 
1831 {
1833  if (!cont) {
1834  AliError(Form("%s:Cluster container %d not found",GetName(),i));
1835  return 0;
1836  }
1837  return cont->GetArray();
1838 }
1839 
1848 {
1849 
1851  if (!cont) {
1852  AliError(Form("%s: Particle container %d not found",GetName(),c));
1853  return 0;
1854  }
1855  AliVParticle *vp = cont->GetAcceptParticle(p);
1856 
1857  return vp;
1858 }
1859 
1868 {
1870  if (!cont) {
1871  AliError(Form("%s: Cluster container %d not found",GetName(),c));
1872  return 0;
1873  }
1874  AliVCluster *vc = cont->GetAcceptCluster(cl);
1875 
1876  return vc;
1877 }
1878 
1885 {
1887  if (!cont) {
1888  AliError(Form("%s: Particle container %d not found",GetName(),i));
1889  return 0;
1890  }
1891  return cont->GetNEntries();
1892 }
1893 
1901 {
1903  if (!cont) {
1904  AliError(Form("%s: Cluster container %d not found",GetName(),i));
1905  return 0;
1906  }
1907  return cont->GetNEntries();
1908 }
1909 
1921 AliEMCALTriggerPatchInfo* AliAnalysisTaskEmcal::GetMainTriggerPatch(TriggerCategory trigger, Bool_t doSimpleOffline)
1922 {
1923 
1924  if (!fTriggerPatchInfo) {
1925  AliError(Form("%s: fTriggerPatchInfo not available",GetName()));
1926  return 0;
1927  }
1928 
1929  //number of patches in event
1930  Int_t nPatch = fTriggerPatchInfo->GetEntries();
1931 
1932  //extract main trigger patch(es)
1933  AliEMCALTriggerPatchInfo *patch(NULL), *selected(NULL);
1934  for (Int_t iPatch = 0; iPatch < nPatch; iPatch++) {
1935 
1936  patch = (AliEMCALTriggerPatchInfo*)fTriggerPatchInfo->At( iPatch );
1937  if (patch->IsMainTrigger()) {
1938  if(doSimpleOffline){
1939  if(patch->IsOfflineSimple()){
1940  switch(trigger){
1941  case kTriggerLevel0:
1942  // option not yet implemented in the trigger maker
1943  if(patch->IsLevel0()) selected = patch;
1944  break;
1945  case kTriggerLevel1Jet:
1946  if(patch->IsJetHighSimple() || patch->IsJetLowSimple()){
1947  if(!selected) selected = patch;
1948  else if(patch->GetADCOfflineAmp() > selected->GetADCOfflineAmp()) selected = patch;
1949  }
1950  break;
1951  case kTriggerLevel1Gamma:
1952  if(patch->IsGammaHighSimple() || patch->IsGammaLowSimple()){
1953  if(!selected) selected = patch;
1954  else if(patch->GetADCOfflineAmp() > selected->GetADCOfflineAmp()) selected = patch;
1955  }
1956  break;
1957  default: // Silence compiler warnings
1958  AliError("Untreated case: Main Patch is recalculated; should be in 'else' branch");
1959  };
1960  }
1961  } else { // Not OfflineSimple
1962  switch(trigger){
1963  case kTriggerLevel0:
1964  if(patch->IsLevel0()) selected = patch;
1965  break;
1966  case kTriggerLevel1Jet:
1967  if(patch->IsJetHigh() || patch->IsJetLow()){
1968  if(!selected) selected = patch;
1969  else if (patch->GetADCAmp() > selected->GetADCAmp())
1970  selected = patch;
1971  }
1972  break;
1973  case kTriggerLevel1Gamma:
1974  if(patch->IsGammaHigh() || patch->IsGammaLow()){
1975  if(!selected) selected = patch;
1976  else if (patch->GetADCAmp() > selected->GetADCAmp())
1977  selected = patch;
1978  }
1979  break;
1980  default:
1981  AliError("Untreated case: Main Patch is recalculated; should be in 'else' branch");
1982  };
1983  }
1984  }
1985  else if ((trigger == kTriggerRecalcJet && patch->IsRecalcJet()) ||
1986  (trigger == kTriggerRecalcGamma && patch->IsRecalcGamma())) { // recalculated patches
1987  if (doSimpleOffline && patch->IsOfflineSimple()) {
1988  if(!selected) selected = patch;
1989  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
1990  selected = patch;
1991  }
1992  else if (!doSimpleOffline && !patch->IsOfflineSimple()) {
1993  if(!selected) selected = patch;
1994  else if (patch->GetADCAmp() > selected->GetADCAmp())
1995  selected = patch;
1996  }
1997  }
1998  }
1999  return selected;
2000 }
2001 
2008 {
2009  if (!(InputEvent()->FindListObject(obj->GetName()))) {
2010  InputEvent()->AddObject(obj);
2011  }
2012  else {
2013  if (!attempt) {
2014  AliFatal(Form("%s: Container with name %s already present. Aborting", GetName(), obj->GetName()));
2015  }
2016  }
2017 }
2018 
2027 {
2028 
2029  if (!fGeom) {
2030  AliWarning(Form("%s - AliAnalysisTaskEmcal::IsTrackInEmcalAcceptance - Geometry is not available!", GetName()));
2031  return kFALSE;
2032  }
2033 
2034  Double_t minPhi = fGeom->GetArm1PhiMin() - edges;
2035  Double_t maxPhi = fGeom->GetArm1PhiMax() + edges;
2036 
2037  if (part->Phi() > minPhi && part->Phi() < maxPhi) {
2038  return kTRUE;
2039  }
2040  else {
2041  return kFALSE;
2042  }
2043 }
2044 
2046 {
2047  axis->SetBinLabel(1, "NullObject");
2048  axis->SetBinLabel(2, "Pt");
2049  axis->SetBinLabel(3, "Acceptance");
2050  axis->SetBinLabel(4, "MCLabel");
2051  axis->SetBinLabel(5, "BitMap");
2052  axis->SetBinLabel(6, "HF cut");
2053  axis->SetBinLabel(7, "Bit6");
2054  axis->SetBinLabel(8, "NotHybridTrack");
2055  axis->SetBinLabel(9, "MCFlag");
2056  axis->SetBinLabel(10, "MCGenerator");
2057  axis->SetBinLabel(11, "ChargeCut");
2058  axis->SetBinLabel(12, "MinDistanceTPCSectorEdge");
2059  axis->SetBinLabel(13, "Bit12");
2060  axis->SetBinLabel(14, "IsEMCal");
2061  axis->SetBinLabel(15, "Time");
2062  axis->SetBinLabel(16, "Energy");
2063  axis->SetBinLabel(17, "ExoticCut");
2064  axis->SetBinLabel(18, "Bit17");
2065  axis->SetBinLabel(19, "Area");
2066  axis->SetBinLabel(20, "AreaEmc");
2067  axis->SetBinLabel(21, "ZLeadingCh");
2068  axis->SetBinLabel(22, "ZLeadingEmc");
2069  axis->SetBinLabel(23, "NEF");
2070  axis->SetBinLabel(24, "MinLeadPt");
2071  axis->SetBinLabel(25, "MaxTrackPt");
2072  axis->SetBinLabel(26, "MaxClusterPt");
2073  axis->SetBinLabel(27, "Flavour");
2074  axis->SetBinLabel(28, "TagStatus");
2075  axis->SetBinLabel(29, "MinNConstituents");
2076  axis->SetBinLabel(30, "Bit29");
2077  axis->SetBinLabel(31, "Bit30");
2078  axis->SetBinLabel(32, "Bit31");
2079 }
2080 
2087 Double_t AliAnalysisTaskEmcal::GetParallelFraction(AliVParticle* part1, AliVParticle* part2)
2088 {
2089  TVector3 vect1(part1->Px(), part1->Py(), part1->Pz());
2090  TVector3 vect2(part2->Px(), part2->Py(), part2->Pz());
2091  Double_t z = (vect1 * vect2) / (vect2 * vect2);
2092  return z;
2093 }
2094 
2101 Double_t AliAnalysisTaskEmcal::GetParallelFraction(const TVector3& vect1, AliVParticle* part2)
2102 {
2103  TVector3 vect2(part2->Px(), part2->Py(), part2->Pz());
2104  Double_t z = (vect1 * vect2) / (vect2 * vect2);
2105  return z;
2106 }
2107 
2116 void AliAnalysisTaskEmcal::GetEtaPhiDiff(const AliVTrack *t, const AliVCluster *v, Double_t &phidiff, Double_t &etadiff)
2117 {
2118  phidiff = 999;
2119  etadiff = 999;
2120 
2121  if (!t||!v) return;
2122 
2123  Double_t veta = t->GetTrackEtaOnEMCal();
2124  Double_t vphi = t->GetTrackPhiOnEMCal();
2125 
2126  Float_t pos[3] = {0};
2127  v->GetPosition(pos);
2128  TVector3 cpos(pos);
2129  Double_t ceta = cpos.Eta();
2130  Double_t cphi = cpos.Phi();
2131  etadiff=veta-ceta;
2132  phidiff=TVector2::Phi_mpi_pi(vphi-cphi);
2133 }
2134 
2140 Byte_t AliAnalysisTaskEmcal::GetTrackType(const AliVTrack *t)
2141 {
2142  Byte_t ret = 0;
2143  if (t->TestBit(BIT(22)) && !t->TestBit(BIT(23)))
2144  ret = 1;
2145  else if (!t->TestBit(BIT(22)) && t->TestBit(BIT(23)))
2146  ret = 2;
2147  else if (t->TestBit(BIT(22)) && t->TestBit(BIT(23)))
2148  ret = 3;
2149  return ret;
2150 }
2151 
2161 Byte_t AliAnalysisTaskEmcal::GetTrackType(const AliAODTrack *aodTrack, UInt_t filterBit1, UInt_t filterBit2)
2162 {
2163 
2164  Int_t res = 0;
2165 
2166  if (aodTrack->TestFilterBit(filterBit1)) {
2167  res = 0;
2168  }
2169  else if (aodTrack->TestFilterBit(filterBit2)) {
2170  if ((aodTrack->GetStatus()&AliVTrack::kITSrefit)!=0) {
2171  res = 1;
2172  }
2173  else {
2174  res = 2;
2175  }
2176  }
2177  else {
2178  res = 3;
2179  }
2180 
2181  return res;
2182 }
2183 
2189 {
2190  if (!fPythiaInfo) {
2192  }
2193 
2194  AliStack* stack = mcEvent->Stack();
2195 
2196  const Int_t nprim = stack->GetNprimary();
2197  // reject if partons are missing from stack for some reason
2198  if (nprim < 8) return;
2199 
2200  TParticle *part6 = stack->Particle(6);
2201  TParticle *part7 = stack->Particle(7);
2202 
2203  fPythiaInfo->SetPartonFlag6(TMath::Abs(part6->GetPdgCode()));
2204  fPythiaInfo->SetParton6(part6->Pt(), part6->Eta(), part6->Phi(), part6->GetMass());
2205 
2206  fPythiaInfo->SetPartonFlag7(TMath::Abs(part7->GetPdgCode()));
2207  fPythiaInfo->SetParton7(part7->Pt(), part7->Eta(), part7->Phi(), part7->GetMass());
2208 
2209  AliGenPythiaEventHeader *pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(mcEvent->GenEventHeader());
2210  if(pythiaGenHeader){
2211  Float_t ptWeight=pythiaGenHeader->EventWeight();
2212  fPythiaInfo->SetPythiaEventWeight(ptWeight);}
2213 }
2214 
2220 {
2221  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2222  if (!mgr) {
2223  ::Error("AddAODHandler", "No analysis manager to connect to.");
2224  return NULL;
2225  }
2226 
2227  AliAODInputHandler* aodHandler = new AliAODInputHandler();
2228 
2229  AliVEventHandler *inputHandler=mgr->GetInputEventHandler();
2230  if (inputHandler && (inputHandler->IsA() == AliMultiInputEventHandler::Class())) {
2231  AliMultiInputEventHandler *multiInputHandler=(AliMultiInputEventHandler*)inputHandler;
2232  multiInputHandler->AddInputEventHandler(aodHandler);
2233  }
2234  else {
2235  if (!inputHandler) {
2236  mgr->SetInputEventHandler(aodHandler);
2237  }
2238  else {
2239  ::Error("AddAODHandler", "inputHandler is NOT null. AOD handler was NOT added !!!");
2240  return NULL;
2241  }
2242  }
2243 
2244  return aodHandler;
2245 }
2246 
2252 {
2253  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2254  if (!mgr) {
2255  ::Error("AddESDHandler", "No analysis manager to connect to.");
2256  return NULL;
2257  }
2258 
2259  AliESDInputHandler *esdHandler = new AliESDInputHandler();
2260 
2261  AliVEventHandler *inputHandler=mgr->GetInputEventHandler();
2262  if (inputHandler && (inputHandler->IsA() == AliMultiInputEventHandler::Class())) {
2263  AliMultiInputEventHandler *multiInputHandler=(AliMultiInputEventHandler*)inputHandler;
2264  multiInputHandler->AddInputEventHandler(esdHandler);
2265  }
2266  else {
2267  if (!inputHandler) {
2268  mgr->SetInputEventHandler(esdHandler);
2269  }
2270  else {
2271  ::Error("AddESDHandler", "inputHandler is NOT null. ESD handler was NOT added !!!");
2272  return NULL;
2273  }
2274  }
2275 
2276  return esdHandler;
2277 }
2278 
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
double Double_t
Definition: External.C:58
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
static AliEmcalDownscaleFactorsOCDB * Instance()
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.
Bool_t fLocalInitialized
whether or not the task has been already initialized
Bool_t fUsePtHardBinScaling
Use -hard bin scaling in merging.
void SetPartonFlag7(Int_t flag7)
Container with name, TClonesArray and cuts for particles.
Bool_t fUseXsecFromHeader
! Use cross section from header instead of pyxsec.root (purely transient)
TSystem * gSystem
Double_t fPtHard
!event -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
Int_t fNPtHardBins
Number of -hard bins in the dataset.
Bool_t fGeneralHistograms
whether or not it should fill some general histograms
Bool_t AcceptCluster(AliVCluster *clus, Int_t c=0) const
TCanvas * c
Definition: TestFitELoss.C:172
Int_t fCentBin
!event centrality bin
TH1 * fHistEventsAfterSel
!total number of events per pt hard bin after selection
const AliMCParticleIterableContainer all() const
Float_t fPtHardAndClusterPtFactor
Factor between ptHard and cluster pT to reject/accept event.
Double_t fMinPtTrackInEmcal
min pt track in emcal
Double_t GetDownscaleFactorForTriggerClass(const TString &trigger) const
TH1 * fHistEventPlane
!event plane distribution
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
int Int_t
Definition: External.C:63
TH1 * fHistTriggerClasses
!number of events in each trigger class
Bool_t fCountDownscaleCorrectedEvents
Count event number corrected for downscaling.
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
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
TH1 * fHistTriggerClassesCorr
!corrected number of events in each trigger class
Bool_t fIsHerwig
trigger, if it is a HERWIG production
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
void SetUseScaling(Bool_t val)
Definition: AliEmcalList.h:31
Int_t fNcentBins
how many centrality bins
AliClusterContainer * GetClusterContainer(Int_t i=0) const
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
Float_t fPtHardAndJetPtFactor
Factor between ptHard and jet pT to reject/accept event.
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
std::vector< TString > GetTriggerClasses() const
Float_t fPtHardAndTrackPtFactor
Factor between ptHard and track pT to reject/accept event.
Double_t fCent
!event centrality
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
Bool_t fMCRejectFilter
enable the filtering of events by tail rejection
unsigned long ULong_t
Definition: External.C:38
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
TClonesArray * GetArrayFromEvent(const char *name, const char *clname=0)
Enhanced TList-derived class that implements correct merging for pt_hard binned production.
Definition: AliEmcalList.h:25
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)
AliEmcalList * fOutput
!output list
Double_t fMaxBinPt
max pt in histograms
Int_t fPtHardBin
!event -hard bin
virtual void RunChanged(Int_t)
TClonesArray * fTracks
!tracks
TH1 * fHistTrialsAfterSel
!total number of trials per pt hard bin after selection
AliGenHerwigEventHeader * fHerwigHeader
!event Herwig header
void LoadPythiaInfo(AliVEvent *event)
Bool_t fIsEsd
!whether it's an ESD analysis
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
static AliAODInputHandler * AddAODHandler()
Double_t fVertex[3]
!event vertex
AliTrackContainer * AddTrackContainer(const char *n)
Handler for downscale factors for various triggers obtained from the OCDB.
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.
const char Option_t
Definition: External.C:48
void SetClusPtCut(Double_t cut)
Int_t fRunNumber
!run number (triggering RunChanged()
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
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)
bool Bool_t
Definition: External.C:53
ULong_t fTriggers
list of fired triggers
static AliESDInputHandler * AddESDHandler()
Double_t fMaxEventPlane
maximum event plane value
void SetPythiaEventWeight(Float_t ptWeight)
Float_t fXsection
!x-section from pythia header
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.
TString fMinBiasRefTrigger
Name of the minmum bias reference trigger, used in the calculation of downscale-corrected event numbe...
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.
TArrayI fPtHardBinning
-hard binning
TList * OpenFile(const char *fname)
Definition: DrawAnaELoss.C:65
AliEMCALTriggerPatchInfo * GetMainTriggerPatch(TriggerCategory triggersel=kTriggerLevel1Jet, Bool_t doOfflinSimple=kFALSE)
static Double_t fgkEMCalDCalPhiDivide
phi value used to distinguish between DCal and EMCal