AliPhysics  c2a2810 (c2a2810)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliAnalysisTaskEmcal.cxx
Go to the documentation of this file.
1 //
2 // Emcal base analysis task.
3 //
4 // Author: S.Aiola, M. Verweij
5 
6 #include <RVersion.h>
7 #include "AliAnalysisTaskEmcal.h"
8 
9 #include <TClonesArray.h>
10 #include <TList.h>
11 #include <TObject.h>
12 #include <TH1F.h>
13 #include <TProfile.h>
14 #include <TSystem.h>
15 #include <TFile.h>
16 #include <TChain.h>
17 #include <TKey.h>
18 
19 #include "AliAODEvent.h"
20 #include "AliAnalysisManager.h"
21 #include "AliCentrality.h"
22 #include "AliEMCALGeometry.h"
23 #include "AliESDEvent.h"
24 #include "AliEmcalParticle.h"
25 #include "AliEventplane.h"
26 #include "AliInputEventHandler.h"
27 #include "AliLog.h"
28 #include "AliMCParticle.h"
29 #include "AliVCluster.h"
30 #include "AliVEventHandler.h"
31 #include "AliVParticle.h"
32 #include "AliVCaloTrigger.h"
33 #include "AliGenPythiaEventHeader.h"
34 #include "AliAODMCHeader.h"
35 #include "AliMCEvent.h"
36 #include "AliAnalysisUtils.h"
38 
39 #include "AliParticleContainer.h"
40 #include "AliClusterContainer.h"
41 
43 
44 //________________________________________________________________________
46  AliAnalysisTaskSE("AliAnalysisTaskEmcal"),
47  fForceBeamType(kNA),
48  fGeneralHistograms(kFALSE),
49  fInitialized(kFALSE),
50  fCreateHisto(kTRUE),
51  fCaloCellsName(),
52  fCaloTriggersName(),
53  fCaloTriggerPatchInfoName(),
54  fMinCent(-999),
55  fMaxCent(-999),
56  fMinVz(-999),
57  fMaxVz(-999),
58  fTrackPtCut(0),
59  fMinNTrack(0),
60  fUseAliAnaUtils(kFALSE),
61  fRejectPileup(kFALSE),
62  fTklVsClusSPDCut(kFALSE),
63  fOffTrigger(AliVEvent::kAny),
64  fTrigClass(),
65  fTriggerTypeSel(kND),
66  fNbins(500),
67  fMinBinPt(0),
68  fMaxBinPt(250),
69  fMinPtTrackInEmcal(0),
70  fEventPlaneVsEmcal(-1),
71  fMinEventPlane(-1e6),
72  fMaxEventPlane(1e6),
73  fCentEst("V0M"),
74  fIsEmbedded(kFALSE),
75  fIsPythia(kFALSE),
76  fSelectPtHardBin(-999),
77  fMinMCLabel(0),
78  fMCLabelShift(0),
79  fNcentBins(4),
80  fNeedEmcalGeom(kTRUE),
81  fParticleCollArray(),
82  fClusterCollArray(),
83  fTriggers(0),
84  fAliAnalysisUtils(0x0),
85  fIsEsd(kFALSE),
86  fGeom(0),
87  fTracks(0),
88  fCaloClusters(0),
89  fCaloCells(0),
90  fCaloTriggers(0),
91  fTriggerPatchInfo(0),
92  fCent(0),
93  fCentBin(-1),
94  fEPV0(-1.0),
95  fEPV0A(-1.0),
96  fEPV0C(-1.0),
97  fNVertCont(0),
98  fBeamType(kNA),
99  fPythiaHeader(0),
100  fPtHard(0),
101  fPtHardBin(0),
102  fNTrials(0),
103  fXsection(0),
104  fOutput(0),
105  fHistEventCount(0),
106  fHistTrialsAfterSel(0),
107  fHistEventsAfterSel(0),
108  fHistXsectionAfterSel(0),
109  fHistTrials(0),
110  fHistEvents(0),
111  fHistXsection(0),
112  fHistPtHard(0),
113  fHistCentrality(0),
114  fHistZVertex(0),
115  fHistEventPlane(0),
116  fHistEventRejection(0)
117 {
118  // Default constructor.
119 
120  fVertex[0] = 0;
121  fVertex[1] = 0;
122  fVertex[2] = 0;
123 
124  fParticleCollArray.SetOwner(kTRUE);
125  fClusterCollArray.SetOwner(kTRUE);
126 }
127 
128 //________________________________________________________________________
129 AliAnalysisTaskEmcal::AliAnalysisTaskEmcal(const char *name, Bool_t histo) :
130  AliAnalysisTaskSE(name),
131  fForceBeamType(kNA),
132  fGeneralHistograms(kFALSE),
133  fInitialized(kFALSE),
134  fCreateHisto(histo),
135  fCaloCellsName(),
136  fCaloTriggersName(),
137  fCaloTriggerPatchInfoName(),
138  fMinCent(-999),
139  fMaxCent(-999),
140  fMinVz(-999),
141  fMaxVz(-999),
142  fTrackPtCut(0),
143  fMinNTrack(0),
144  fUseAliAnaUtils(kFALSE),
145  fRejectPileup(kFALSE),
146  fTklVsClusSPDCut(kFALSE),
147  fOffTrigger(AliVEvent::kAny),
148  fTrigClass(),
149  fTriggerTypeSel(kND),
150  fNbins(500),
151  fMinBinPt(0),
152  fMaxBinPt(250),
153  fMinPtTrackInEmcal(0),
154  fEventPlaneVsEmcal(-1),
155  fMinEventPlane(-1e6),
156  fMaxEventPlane(1e6),
157  fCentEst("V0M"),
158  fIsEmbedded(kFALSE),
159  fIsPythia(kFALSE),
160  fSelectPtHardBin(-999),
161  fMinMCLabel(0),
162  fMCLabelShift(0),
163  fNcentBins(4),
164  fNeedEmcalGeom(kTRUE),
165  fParticleCollArray(),
166  fClusterCollArray(),
167  fTriggers(0),
168  fAliAnalysisUtils(0x0),
169  fIsEsd(kFALSE),
170  fGeom(0),
171  fTracks(0),
172  fCaloClusters(0),
173  fCaloCells(0),
174  fCaloTriggers(0),
175  fTriggerPatchInfo(0),
176  fCent(0),
177  fCentBin(-1),
178  fEPV0(-1.0),
179  fEPV0A(-1.0),
180  fEPV0C(-1.0),
181  fNVertCont(0),
182  fBeamType(kNA),
183  fPythiaHeader(0),
184  fPtHard(0),
185  fPtHardBin(0),
186  fNTrials(0),
187  fXsection(0),
188  fOutput(0),
189  fHistEventCount(0),
190  fHistTrialsAfterSel(0),
191  fHistEventsAfterSel(0),
192  fHistXsectionAfterSel(0),
193  fHistTrials(0),
194  fHistEvents(0),
195  fHistXsection(0),
196  fHistPtHard(0),
197  fHistCentrality(0),
198  fHistZVertex(0),
199  fHistEventPlane(0),
200  fHistEventRejection(0)
201 {
202  // Standard constructor.
203 
204  fVertex[0] = 0;
205  fVertex[1] = 0;
206  fVertex[2] = 0;
207 
208  fParticleCollArray.SetOwner(kTRUE);
209  fClusterCollArray.SetOwner(kTRUE);
210 
211  if (fCreateHisto) {
212  DefineOutput(1, TList::Class());
213  }
214 }
215 
216 //________________________________________________________________________
218 {
219  // Destructor
220 }
221 
222 //________________________________________________________________________
223 void AliAnalysisTaskEmcal::SetClusPtCut(Double_t cut, Int_t c)
224 {
226  if (cont) cont->SetClusPtCut(cut);
227  else AliError(Form("%s in SetClusPtCut(...): container %d not found",GetName(),c));
228 }
229 
230 //________________________________________________________________________
231 void AliAnalysisTaskEmcal::SetClusTimeCut(Double_t min, Double_t max, Int_t c)
232 {
234  if (cont) cont->SetClusTimeCut(min,max);
235  else AliError(Form("%s in SetClusTimeCut(...): container %d not found",GetName(),c));
236 }
237 
238 //________________________________________________________________________
239 void AliAnalysisTaskEmcal::SetTrackPtCut(Double_t cut, Int_t c)
240 {
242  if (cont) cont->SetParticlePtCut(cut);
243  else AliError(Form("%s in SetTrackPtCut(...): container %d not found",GetName(),c));
244 
245  fTrackPtCut = cut;
246 }
247 
248 //________________________________________________________________________
249 void AliAnalysisTaskEmcal::SetTrackEtaLimits(Double_t min, Double_t max, Int_t c)
250 {
252  if (cont) cont->SetParticleEtaLimits(min,max);
253  else AliError(Form("%s in SetTrackPtCut(...): container %d not found",GetName(),c));
254 }
255 
256 //________________________________________________________________________
257 void AliAnalysisTaskEmcal::SetTrackPhiLimits(Double_t min, Double_t max, Int_t c)
258 {
260  if (cont) cont->SetParticlePhiLimits(min,max);
261  else AliError(Form("%s in SetTrackPhiLimits(...): container %d not found",GetName(),c));
262 }
263 
264 //________________________________________________________________________
266 {
267  // Create user output.
268 
269  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
270  if (mgr) {
271  AliVEventHandler *evhand = mgr->GetInputEventHandler();
272  if (evhand) {
273  if (evhand->InheritsFrom("AliESDInputHandler")) {
274  fIsEsd = kTRUE;
275  }
276  else {
277  fIsEsd = kFALSE;
278  }
279  }
280  else {
281  AliError("Event handler not found!");
282  }
283  }
284  else {
285  AliError("Analysis manager not found!");
286  }
287 
288  if (!fCreateHisto)
289  return;
290 
291  OpenFile(1);
292  fOutput = new TList();
293  fOutput->SetOwner();
294 
295  if (fForceBeamType == kpp)
296  fNcentBins = 1;
297 
298  if (!fGeneralHistograms)
299  return;
300 
301  if (fIsPythia) {
302  fHistTrialsAfterSel = new TH1F("fHistTrialsAfterSel", "fHistTrialsAfterSel", 11, 0, 11);
303  fHistTrialsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
304  fHistTrialsAfterSel->GetYaxis()->SetTitle("trials");
306 
307  fHistEventsAfterSel = new TH1F("fHistEventsAfterSel", "fHistEventsAfterSel", 11, 0, 11);
308  fHistEventsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
309  fHistEventsAfterSel->GetYaxis()->SetTitle("total events");
311 
312  fHistXsectionAfterSel = new TProfile("fHistXsectionAfterSel", "fHistXsectionAfterSel", 11, 0, 11);
313  fHistXsectionAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
314  fHistXsectionAfterSel->GetYaxis()->SetTitle("xsection");
316 
317  fHistTrials = new TH1F("fHistTrials", "fHistTrials", 11, 0, 11);
318  fHistTrials->GetXaxis()->SetTitle("p_{T} hard bin");
319  fHistTrials->GetYaxis()->SetTitle("trials");
320  fOutput->Add(fHistTrials);
321 
322  fHistEvents = new TH1F("fHistEvents", "fHistEvents", 11, 0, 11);
323  fHistEvents->GetXaxis()->SetTitle("p_{T} hard bin");
324  fHistEvents->GetYaxis()->SetTitle("total events");
325  fOutput->Add(fHistEvents);
326 
327  fHistXsection = new TProfile("fHistXsection", "fHistXsection", 11, 0, 11);
328  fHistXsection->GetXaxis()->SetTitle("p_{T} hard bin");
329  fHistXsection->GetYaxis()->SetTitle("xsection");
330  fOutput->Add(fHistXsection);
331 
332  const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
333  const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
334 
335  for (Int_t i = 1; i < 12; i++) {
336  fHistTrialsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
337  fHistEventsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
338 
339  fHistTrials->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
340  fHistXsection->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
341  fHistEvents->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
342  }
343 
344  fHistPtHard = new TH1F("fHistPtHard", "fHistPtHard", fNbins*2, fMinBinPt, fMaxBinPt*4);
345  fHistPtHard->GetXaxis()->SetTitle("p_{T,hard} (GeV/c)");
346  fHistPtHard->GetYaxis()->SetTitle("counts");
347  fOutput->Add(fHistPtHard);
348  }
349 
350  fHistZVertex = new TH1F("fHistZVertex","Z vertex position", 60, -30, 30);
351  fHistZVertex->GetXaxis()->SetTitle("z");
352  fHistZVertex->GetYaxis()->SetTitle("counts");
353  fOutput->Add(fHistZVertex);
354 
355  if (fForceBeamType != kpp) {
356  fHistCentrality = new TH1F("fHistCentrality","Event centrality distribution", 200, 0, 100);
357  fHistCentrality->GetXaxis()->SetTitle("Centrality (%)");
358  fHistCentrality->GetYaxis()->SetTitle("counts");
359  fOutput->Add(fHistCentrality);
360 
361  fHistEventPlane = new TH1F("fHistEventPlane","Event plane", 120, -TMath::Pi(), TMath::Pi());
362  fHistEventPlane->GetXaxis()->SetTitle("event plane");
363  fHistEventPlane->GetYaxis()->SetTitle("counts");
364  fOutput->Add(fHistEventPlane);
365  }
366 
367  fHistEventRejection = new TH1F("fHistEventRejection","Reasons to reject event",20,0,20);
368 #if ROOT_VERSION_CODE < ROOT_VERSION(6,5,1)
369  fHistEventRejection->SetBit(TH1::kCanRebin);
370 #endif
371  fHistEventRejection->GetXaxis()->SetBinLabel(1,"PhysSel");
372  fHistEventRejection->GetXaxis()->SetBinLabel(2,"trigger");
373  fHistEventRejection->GetXaxis()->SetBinLabel(3,"trigTypeSel");
374  fHistEventRejection->GetXaxis()->SetBinLabel(4,"Cent");
375  fHistEventRejection->GetXaxis()->SetBinLabel(5,"vertex contr.");
376  fHistEventRejection->GetXaxis()->SetBinLabel(6,"Vz");
377  fHistEventRejection->GetXaxis()->SetBinLabel(7,"trackInEmcal");
378  fHistEventRejection->GetXaxis()->SetBinLabel(8,"minNTrack");
379  fHistEventRejection->GetXaxis()->SetBinLabel(9,"VtxSel2013pA");
380  fHistEventRejection->GetXaxis()->SetBinLabel(10,"PileUp");
381  fHistEventRejection->GetXaxis()->SetBinLabel(11,"EvtPlane");
382  fHistEventRejection->GetXaxis()->SetBinLabel(12,"SelPtHardBin");
383  fHistEventRejection->GetXaxis()->SetBinLabel(13,"Bkg evt");
384  fHistEventRejection->GetYaxis()->SetTitle("counts");
386 
387  fHistEventCount = new TH1F("fHistEventCount","fHistEventCount",2,0,2);
388  fHistEventCount->GetXaxis()->SetBinLabel(1,"Accepted");
389  fHistEventCount->GetXaxis()->SetBinLabel(2,"Rejected");
390  fHistEventCount->GetYaxis()->SetTitle("counts");
391  fOutput->Add(fHistEventCount);
392 
393  PostData(1, fOutput);
394 }
395 
396 //________________________________________________________________________
398 {
399  if (fIsPythia) {
403  fHistPtHard->Fill(fPtHard);
404  }
405 
406  fHistZVertex->Fill(fVertex[2]);
407 
408  if (fForceBeamType != kpp) {
409  fHistCentrality->Fill(fCent);
410  fHistEventPlane->Fill(fEPV0);
411  }
412 
413  return kTRUE;
414 }
415 
416 //________________________________________________________________________
418 {
419  // Main loop, called for each event.
420 
421  if (!fInitialized)
422  ExecOnce();
423 
424  if (!fInitialized)
425  return;
426 
427  if (!RetrieveEventObjects())
428  return;
429 
430  if (IsEventSelected()) {
431  if (fGeneralHistograms) fHistEventCount->Fill("Accepted",1);
432  }
433  else {
434  if (fGeneralHistograms) fHistEventCount->Fill("Rejected",1);
435  return;
436  }
437 
439  if (!FillGeneralHistograms())
440  return;
441  }
442 
443  if (!Run())
444  return;
445 
446  if (fCreateHisto) {
447  if (!FillHistograms())
448  return;
449  }
450 
451  if (fCreateHisto && fOutput) {
452  // information for this iteration of the UserExec in the container
453  PostData(1, fOutput);
454  }
455 }
456 
457 //________________________________________________________________________
458 Bool_t AliAnalysisTaskEmcal::AcceptCluster(AliVCluster *clus, Int_t c) const
459 {
460  // Return true if cluster is accepted.
461 
462  if (!clus)
463  return kFALSE;
464 
466  if (!cont) {
467  AliError(Form("%s:Container %d not found",GetName(),c));
468  return 0;
469  }
470 
471  return cont->AcceptCluster(clus);
472 }
473 
474 //________________________________________________________________________
475 Bool_t AliAnalysisTaskEmcal::AcceptTrack(AliVParticle *track, Int_t c) const
476 {
477  // Return true if track is accepted.
478 
479  if (!track)
480  return kFALSE;
481 
483  if (!cont) {
484  AliError(Form("%s:Container %d not found",GetName(),c));
485  return 0;
486  }
487 
488  return cont->AcceptParticle(track);
489 }
490 
491 //________________________________________________________________________
492 Bool_t AliAnalysisTaskEmcal::PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
493 {
494  //
495  // Get the cross section and the trails either from pyxsec.root or from pysec_hists.root
496  // Get the pt hard bin from the file path
497  // This is to called in Notify and should provide the path to the AOD/ESD file
498  // (Partially copied from AliAnalysisHelperJetTasks)
499 
500  TString file(currFile);
501  fXsec = 0;
502  fTrials = 1;
503 
504  if (file.Contains(".zip#")) {
505  Ssiz_t pos1 = file.Index("root_archive",12,0,TString::kExact);
506  Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
507  Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
508  file.Replace(pos+1,pos2-pos1,"");
509  } else {
510  // not an archive take the basename....
511  file.ReplaceAll(gSystem->BaseName(file.Data()),"");
512  }
513  AliDebug(1,Form("File name: %s",file.Data()));
514 
515  // Get the pt hard bin
516  TString strPthard(file);
517 
518  strPthard.Remove(strPthard.Last('/'));
519  strPthard.Remove(strPthard.Last('/'));
520  if (strPthard.Contains("AOD")) strPthard.Remove(strPthard.Last('/'));
521  strPthard.Remove(0,strPthard.Last('/')+1);
522  if (strPthard.IsDec())
523  pthard = strPthard.Atoi();
524  else
525  AliWarning(Form("Could not extract file number from path %s", strPthard.Data()));
526 
527  // 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
528  TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root"));
529 
530  if (!fxsec) {
531  // next trial fetch the histgram file
532  fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
533  if (!fxsec) {
534  // not a severe condition but inciate that we have no information
535  return kFALSE;
536  } else {
537  // find the tlist we want to be independtent of the name so use the Tkey
538  TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
539  if (!key) {
540  fxsec->Close();
541  return kFALSE;
542  }
543  TList *list = dynamic_cast<TList*>(key->ReadObj());
544  if (!list) {
545  fxsec->Close();
546  return kFALSE;
547  }
548  fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
549  fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
550  fxsec->Close();
551  }
552  } else { // no tree pyxsec.root
553  TTree *xtree = (TTree*)fxsec->Get("Xsection");
554  if (!xtree) {
555  fxsec->Close();
556  return kFALSE;
557  }
558  UInt_t ntrials = 0;
559  Double_t xsection = 0;
560  xtree->SetBranchAddress("xsection",&xsection);
561  xtree->SetBranchAddress("ntrials",&ntrials);
562  xtree->GetEntry(0);
563  fTrials = ntrials;
564  fXsec = xsection;
565  fxsec->Close();
566  }
567  return kTRUE;
568 }
569 
570 //________________________________________________________________________
572 {
573  // Called when file changes.
574 
576  return kTRUE;
577 
578  TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
579  if (!tree) {
580  AliError(Form("%s - UserNotify: No current tree!",GetName()));
581  return kFALSE;
582  }
583 
584  Float_t xsection = 0;
585  Float_t trials = 0;
586  Int_t pthardbin = 0;
587 
588  TFile *curfile = tree->GetCurrentFile();
589  if (!curfile) {
590  AliError(Form("%s - UserNotify: No current file!",GetName()));
591  return kFALSE;
592  }
593 
594  TChain *chain = dynamic_cast<TChain*>(tree);
595  if (chain) tree = chain->GetTree();
596 
597  Int_t nevents = tree->GetEntriesFast();
598 
599  PythiaInfoFromFile(curfile->GetName(), xsection, trials, pthardbin);
600 
601  // TODO: Workaround
602  if ((pthardbin < 0) || (pthardbin > 10)) pthardbin = 0;
603 
604  fHistTrials->Fill(pthardbin, trials);
605  fHistXsection->Fill(pthardbin, xsection);
606  fHistEvents->Fill(pthardbin, nevents);
607 
608  return kTRUE;
609 }
610 
611 //________________________________________________________________________
613 {
614  // Init the analysis.
615 
616  if (!InputEvent()) {
617  AliError(Form("%s: Could not retrieve event! Returning!", GetName()));
618  return;
619  }
620 
621  if (fNeedEmcalGeom) {
622  fGeom = AliEMCALGeometry::GetInstance();
623  if (!fGeom) {
624  AliError(Form("%s: Can not create geometry", GetName()));
625  return;
626  }
627  }
628 
629  if (fEventPlaneVsEmcal >= 0) {
630  if (fGeom) {
631  Double_t ep = (fGeom->GetArm1PhiMax() + fGeom->GetArm1PhiMin()) / 2 * TMath::DegToRad() + fEventPlaneVsEmcal - TMath::Pi();
632  fMinEventPlane = ep - TMath::Pi() / 4;
633  fMaxEventPlane = ep + TMath::Pi() / 4;
634  }
635  else {
636  AliWarning("Could not set event plane limits because EMCal geometry was not loaded!");
637  }
638  }
639 
640  //Load all requested track branches - each container knows name already
641  for (Int_t i =0; i<fParticleCollArray.GetEntriesFast(); i++) {
642  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.At(i));
643  cont->SetArray(InputEvent());
644  }
645 
646  if (fParticleCollArray.GetEntriesFast()>0) {
648  if (!fTracks) {
649  AliError(Form("%s: Could not retrieve first track branch!", GetName()));
650  return;
651  }
652  }
653 
654  //Load all requested cluster branches - each container knows name already
655  for (Int_t i =0; i<fClusterCollArray.GetEntriesFast(); i++) {
656  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.At(i));
657  cont->SetArray(InputEvent());
658  }
659 
660  if (fClusterCollArray.GetEntriesFast()>0) {
662  if (!fCaloClusters) {
663  AliError(Form("%s: Could not retrieve first cluster branch!", GetName()));
664  return;
665  }
666  }
667 
668  if (!fCaloCellsName.IsNull() && !fCaloCells) {
669  fCaloCells = dynamic_cast<AliVCaloCells*>(InputEvent()->FindListObject(fCaloCellsName));
670  if (!fCaloCells) {
671  AliError(Form("%s: Could not retrieve cells %s!", GetName(), fCaloCellsName.Data()));
672  return;
673  }
674  }
675 
676  if (!fCaloTriggersName.IsNull() && !fCaloTriggers) {
677  fCaloTriggers = dynamic_cast<AliVCaloTrigger*>(InputEvent()->FindListObject(fCaloTriggersName));
678  if (!fCaloTriggers) {
679  AliError(Form("%s: Could not retrieve calo triggers %s!", GetName(), fCaloTriggersName.Data()));
680  return;
681  }
682  }
683 
684  if (!fCaloTriggerPatchInfoName.IsNull() && !fTriggerPatchInfo) {
685  fTriggerPatchInfo = GetArrayFromEvent(fCaloTriggerPatchInfoName.Data(),"AliEmcalTriggerPatchInfo");
686  if (!fTriggerPatchInfo) {
687  AliError(Form("%s: Could not retrieve calo trigger patch info %s!", GetName(), fCaloTriggerPatchInfoName.Data()));
688  return;
689  }
690 
691  }
692 
693  fInitialized = kTRUE;
694 }
695 
696 //_____________________________________________________
698 {
699  // Get beam type : pp-AA-pA
700  // ESDs have it directly, AODs get it from hardcoded run number ranges
701 
702  if (fForceBeamType != kNA)
703  return fForceBeamType;
704 
705  AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
706  if (esd) {
707  const AliESDRun *run = esd->GetESDRun();
708  TString beamType = run->GetBeamType();
709  if (beamType == "p-p")
710  return kpp;
711  else if (beamType == "A-A")
712  return kAA;
713  else if (beamType == "p-A")
714  return kpA;
715  else
716  return kNA;
717  } else {
718  Int_t runNumber = InputEvent()->GetRunNumber();
719  if ((runNumber >= 136851 && runNumber <= 139517) || // LHC10h
720  (runNumber >= 166529 && runNumber <= 170593)) { // LHC11h
721  return kAA;
722  } else if ((runNumber>=188365 && runNumber <= 188366) || // LHC12g
723  (runNumber >= 195344 && runNumber <= 196608)) { // LHC13b-f
724  return kpA;
725  } else {
726  return kpp;
727  }
728  }
729 }
730 
731 //________________________________________________________________________
733 {
734  if (!fTriggerPatchInfo)
735  return 0;
736 
737  //number of patches in event
738  Int_t nPatch = fTriggerPatchInfo->GetEntries();
739 
740  //loop over patches to define trigger type of event
741  Int_t nG1 = 0;
742  Int_t nG2 = 0;
743  Int_t nJ1 = 0;
744  Int_t nJ2 = 0;
745  Int_t nL0 = 0;
747  for (Int_t iPatch = 0; iPatch < nPatch; iPatch++) {
748  patch = (AliEmcalTriggerPatchInfo*)fTriggerPatchInfo->At( iPatch );
749  if (patch->IsGammaHigh()) nG1++;
750  if (patch->IsGammaLow()) nG2++;
751  if (patch->IsJetHigh()) nJ1++;
752  if (patch->IsJetLow()) nJ2++;
753  if (patch->IsLevel0()) nL0++;
754  }
755 
756  AliDebug(2, "Patch summary: ");
757  AliDebug(2, Form("Number of patches: %d", nPatch));
758  AliDebug(2, Form("Jet: low[%d], high[%d]" ,nJ2, nJ1));
759  AliDebug(2, Form("Gamma: low[%d], high[%d]" ,nG2, nG1));
760 
761  ULong_t triggers(0);
762  if (nL0>0) SETBIT(triggers, kL0);
763  if (nG1>0) SETBIT(triggers, kG1);
764  if (nG2>0) SETBIT(triggers, kG2);
765  if (nJ1>0) SETBIT(triggers, kJ1);
766  if (nJ2>0) SETBIT(triggers, kJ2);
767  return triggers;
768 }
769 
770 //________________________________________________________________________
772 {
773  // Check if event has a given trigger type
774  if(trigger==kND) {
775  AliWarning(Form("%s: Requesting undefined trigger type!", GetName()));
776  return kFALSE;
777  }
778  //MV: removing this logic which as far as I can see doesn't make any sense
779  // if(trigger & kND){
780  // return fTriggers == 0;
781  // }
782  return TESTBIT(fTriggers, trigger);
783 }
784 
785 //________________________________________________________________________
787 {
788  // Check if event is selected
789 
790  if (fOffTrigger != AliVEvent::kAny) {
791  UInt_t res = 0;
792  const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
793  if (eev) {
794  res = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
795  } else {
796  const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
797  if (aev) {
798  res = ((AliVAODHeader*)aev->GetHeader())->GetOfflineTrigger();
799  }
800  }
801  if ((res & fOffTrigger) == 0) {
802  if (fGeneralHistograms) fHistEventRejection->Fill("PhysSel",1);
803  return kFALSE;
804  }
805  }
806 
807  if (!fTrigClass.IsNull()) {
808  TString fired;
809  const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
810  if (eev) {
811  fired = eev->GetFiredTriggerClasses();
812  } else {
813  const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
814  if (aev) {
815  fired = aev->GetFiredTriggerClasses();
816  }
817  }
818  if (!fired.Contains("-B-")) {
819  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
820  return kFALSE;
821  }
822 
823  TObjArray *arr = fTrigClass.Tokenize("|");
824  if (!arr) {
825  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
826  return kFALSE;
827  }
828  Bool_t match = 0;
829  for (Int_t i=0;i<arr->GetEntriesFast();++i) {
830  TObject *obj = arr->At(i);
831  if (!obj)
832  continue;
833 
834  //Check if requested trigger was fired, relevant for data sets with 2 emcal trigger thresholds
835  TString objStr = obj->GetName();
836  if(objStr.Contains("J1") || objStr.Contains("J2") || objStr.Contains("G1") || objStr.Contains("G2")) {
837  TString trigType1 = "J1";
838  TString trigType2 = "J2";
839  if(objStr.Contains("G")) {
840  trigType1 = "G1";
841  trigType2 = "G2";
842  }
843  if(objStr.Contains(trigType2) && fired.Contains(trigType2.Data())) { //requesting low threshold + overlap
844  match = 1;
845  break;
846  }
847  else if(objStr.Contains(trigType1) && fired.Contains(trigType1.Data()) && !fired.Contains(trigType2.Data())) { //high threshold only
848  match = 1;
849  break;
850  }
851  } else {
852  if (fired.Contains(obj->GetName())) {
853  match = 1;
854  break;
855  }
856  }
857  }
858  delete arr;
859  if (!match) {
860  if (fGeneralHistograms) fHistEventRejection->Fill("trigger",1);
861  return kFALSE;
862  }
863  }
864 
865  if (fTriggerTypeSel != kND) {
867  if (fGeneralHistograms) fHistEventRejection->Fill("trigTypeSel",1);
868  return kFALSE;
869  }
870  }
871 
872  if ((fMinCent != -999) && (fMaxCent != -999)) {
873  if (fCent<fMinCent || fCent>fMaxCent) {
874  if (fGeneralHistograms) fHistEventRejection->Fill("Cent",1);
875  return kFALSE;
876  }
877  }
878 
879  if (fUseAliAnaUtils) {
880  if (!fAliAnalysisUtils)
881  fAliAnalysisUtils = new AliAnalysisUtils();
882  fAliAnalysisUtils->SetMinVtxContr(2);
883  fAliAnalysisUtils->SetMaxVtxZ(999);
884  if(fMinVz<-10.) fMinVz = -10.;
885  if(fMinVz>10.) fMaxVz = 10.;
886 
887  if (!fAliAnalysisUtils->IsVertexSelected2013pA(InputEvent())) {
888  if (fGeneralHistograms) fHistEventRejection->Fill("VtxSel2013pA",1);
889  return kFALSE;
890  }
891 
892  if (fRejectPileup && fAliAnalysisUtils->IsPileUpEvent(InputEvent())) {
893  if (fGeneralHistograms) fHistEventRejection->Fill("PileUp",1);
894  return kFALSE;
895  }
896 
897  if(fTklVsClusSPDCut && fAliAnalysisUtils->IsSPDClusterVsTrackletBG(InputEvent())) {
898  if (fGeneralHistograms) fHistEventRejection->Fill("Bkg evt",1);
899  return kFALSE;
900  }
901  }
902 
903  if ((fMinVz != -999) && (fMaxVz != -999)) {
904  if (fNVertCont == 0 ) {
905  if (fGeneralHistograms) fHistEventRejection->Fill("vertex contr.",1);
906  return kFALSE;
907  }
908  Double_t vz = fVertex[2];
909  if (vz<fMinVz || vz>fMaxVz) {
910  if (fGeneralHistograms) fHistEventRejection->Fill("Vz",1);
911  return kFALSE;
912  }
913  }
914 
915  if (fMinPtTrackInEmcal > 0 && fGeom) {
916  Bool_t trackInEmcalOk = kFALSE;
917  Int_t ntracks = GetNParticles(0);
918  for (Int_t i = 0; i < ntracks; i++) {
919  AliVParticle *track = GetAcceptParticleFromArray(i,0);
920  if (!track)
921  continue;
922 
923  Double_t phiMin = fGeom->GetArm1PhiMin() * TMath::DegToRad();
924  Double_t phiMax = fGeom->GetArm1PhiMax() * TMath::DegToRad();
925  Int_t runNumber = InputEvent()->GetRunNumber();
926  if (runNumber>=177295 && runNumber<=197470) { //small SM masked in 2012 and 2013
927  phiMin = 1.4;
928  phiMax = TMath::Pi();
929  }
930 
931  if (track->Eta() < fGeom->GetArm1EtaMin() || track->Eta() > fGeom->GetArm1EtaMax() || track->Phi() < phiMin || track->Phi() > phiMax)
932  continue;
933  if (track->Pt() > fMinPtTrackInEmcal) {
934  trackInEmcalOk = kTRUE;
935  break;
936  }
937  }
938  if (!trackInEmcalOk) {
939  if (fGeneralHistograms) fHistEventRejection->Fill("trackInEmcal",1);
940  return kFALSE;
941  }
942  }
943 
944  if (fMinNTrack > 0) {
945  Int_t nTracksAcc = 0;
946  Int_t ntracks = GetNParticles(0);
947  for (Int_t i = 0; i < ntracks; i++) {
948  AliVParticle *track = GetAcceptParticleFromArray(i,0);
949  if (!track)
950  continue;
951  if (track->Pt() > fTrackPtCut) {
952  nTracksAcc++;
953  if (nTracksAcc>=fMinNTrack)
954  break;
955  }
956  }
957  if (nTracksAcc<fMinNTrack) {
958  if (fGeneralHistograms) fHistEventRejection->Fill("minNTrack",1);
959  return kFALSE;
960  }
961  }
962 
963  if (!(fEPV0 > fMinEventPlane && fEPV0 <= fMaxEventPlane) &&
964  !(fEPV0 + TMath::Pi() > fMinEventPlane && fEPV0 + TMath::Pi() <= fMaxEventPlane) &&
965  !(fEPV0 - TMath::Pi() > fMinEventPlane && fEPV0 - TMath::Pi() <= fMaxEventPlane))
966  {
967  if (fGeneralHistograms) fHistEventRejection->Fill("EvtPlane",1);
968  return kFALSE;
969  }
970 
971  if (fSelectPtHardBin != -999 && fSelectPtHardBin != fPtHardBin) {
972  if (fGeneralHistograms) fHistEventRejection->Fill("SelPtHardBin",1);
973  return kFALSE;
974  }
975 
976  return kTRUE;
977 }
978 
979 //________________________________________________________________________
980 TClonesArray *AliAnalysisTaskEmcal::GetArrayFromEvent(const char *name, const char *clname)
981 {
982  // Get array from event.
983 
984  TClonesArray *arr = 0;
985  TString sname(name);
986  if (!sname.IsNull()) {
987  arr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(sname));
988  if (!arr) {
989  AliWarning(Form("%s: Could not retrieve array with name %s!", GetName(), name));
990  return 0;
991  }
992  } else {
993  return 0;
994  }
995 
996  if (!clname)
997  return arr;
998 
999  TString objname(arr->GetClass()->GetName());
1000  TClass cls(objname);
1001  if (!cls.InheritsFrom(clname)) {
1002  AliWarning(Form("%s: Objects of type %s in %s are not inherited from %s!",
1003  GetName(), cls.GetName(), name, clname));
1004  return 0;
1005  }
1006  return arr;
1007 }
1008 
1009 //________________________________________________________________________
1011 {
1012  // Retrieve objects from event.
1013 
1014  fVertex[0] = 0;
1015  fVertex[1] = 0;
1016  fVertex[2] = 0;
1017  fNVertCont = 0;
1018 
1019  const AliVVertex *vert = InputEvent()->GetPrimaryVertex();
1020  if (vert) {
1021  vert->GetXYZ(fVertex);
1022  fNVertCont = vert->GetNContributors();
1023  }
1024 
1025  fBeamType = GetBeamType();
1026 
1027  if (fBeamType == kAA || fBeamType == kpA ) {
1028  AliCentrality *aliCent = InputEvent()->GetCentrality();
1029  if (aliCent) {
1030  fCent = aliCent->GetCentralityPercentile(fCentEst.Data());
1031  if (fNcentBins==4) {
1032  if (fCent >= 0 && fCent < 10) fCentBin = 0;
1033  else if (fCent >= 10 && fCent < 30) fCentBin = 1;
1034  else if (fCent >= 30 && fCent < 50) fCentBin = 2;
1035  else if (fCent >= 50 && fCent <= 100) fCentBin = 3;
1036  else {
1037  AliWarning(Form("%s: Negative centrality: %f. Assuming 99", GetName(), fCent));
1038  fCentBin = fNcentBins-1;
1039  }
1040  } else {
1041  Double_t centWidth = (fMaxCent-fMinCent)/(Double_t)fNcentBins;
1042  if(centWidth>0.)
1043  fCentBin = TMath::FloorNint(fCent/centWidth);
1044  else
1045  fCentBin = 0;
1046  if (fCentBin>=fNcentBins) {
1047  AliWarning(Form("%s: fCentBin too large: cent = %f fCentBin = %d. Assuming 99", GetName(),fCent,fCentBin));
1048  fCentBin = fNcentBins-1;
1049  }
1050  }
1051  } else {
1052  AliWarning(Form("%s: Could not retrieve centrality information! Assuming 99", GetName()));
1053  fCentBin = fNcentBins-1;
1054  }
1055  AliEventplane *aliEP = InputEvent()->GetEventplane();
1056  if (aliEP) {
1057  fEPV0 = aliEP->GetEventplane("V0" ,InputEvent());
1058  fEPV0A = aliEP->GetEventplane("V0A",InputEvent());
1059  fEPV0C = aliEP->GetEventplane("V0C",InputEvent());
1060  } else {
1061  AliWarning(Form("%s: Could not retrieve event plane information!", GetName()));
1062  }
1063  } else {
1064  fCent = 99;
1065  fCentBin = 0;
1066  }
1067 
1068  if (fIsPythia) {
1069 
1070  if (MCEvent()) {
1071  fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
1072  if (!fPythiaHeader) {
1073  // Check if AOD
1074  AliAODMCHeader* aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
1075 
1076  if (aodMCH) {
1077  for (UInt_t i = 0;i<aodMCH->GetNCocktailHeaders();i++) {
1078  fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(aodMCH->GetCocktailHeader(i));
1079  if (fPythiaHeader) break;
1080  }
1081  }
1082  }
1083  }
1084 
1085  if (fPythiaHeader) {
1086  fPtHard = fPythiaHeader->GetPtHard();
1087 
1088  const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
1089  const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
1090  for (fPtHardBin = 0; fPtHardBin < 11; fPtHardBin++) {
1091  if (fPtHard >= ptHardLo[fPtHardBin] && fPtHard < ptHardHi[fPtHardBin])
1092  break;
1093  }
1094 
1095  fXsection = fPythiaHeader->GetXsection();
1096  fNTrials = fPythiaHeader->Trials();
1097  }
1098  }
1099 
1101 
1102  return kTRUE;
1103 }
1104 
1105 //________________________________________________________________________
1107 {
1108  // Add particle container
1109  // will be called in AddTask macro
1110 
1111  TString tmp = TString(n);
1112  if (tmp.IsNull()) return 0;
1113 
1114  AliParticleContainer *cont = 0x0;
1115  cont = new AliParticleContainer();
1116  cont->SetArrayName(n);
1117  TString contName = cont->GetArrayName();
1118 
1119  fParticleCollArray.Add(cont);
1120 
1121  return cont;
1122 }
1123 
1124 //________________________________________________________________________
1126 {
1127  // Add cluster container
1128  // will be called in AddTask macro
1129 
1130  TString tmp = TString(n);
1131  if (tmp.IsNull()) return 0;
1132 
1133  AliClusterContainer *cont = 0x0;
1134  cont = new AliClusterContainer();
1135  cont->SetArrayName(n);
1136 
1137  fClusterCollArray.Add(cont);
1138 
1139  return cont;
1140 }
1141 
1142 //________________________________________________________________________
1144 {
1145  // Get i^th particle container
1146 
1147  if (i<0 || i>fParticleCollArray.GetEntriesFast()) return 0;
1148  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.At(i));
1149  return cont;
1150 }
1151 
1152 //________________________________________________________________________
1154 {
1155  // Get i^th cluster container
1156 
1157  if (i<0 || i>fClusterCollArray.GetEntriesFast()) return 0;
1158  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.At(i));
1159  return cont;
1160 }
1161 
1162 //________________________________________________________________________
1164 {
1165  // Get particle container with name
1166 
1167  AliParticleContainer *cont = static_cast<AliParticleContainer*>(fParticleCollArray.FindObject(name));
1168  return cont;
1169 }
1170 
1171 //________________________________________________________________________
1173 {
1174  // Get cluster container with name
1175 
1176  AliClusterContainer *cont = static_cast<AliClusterContainer*>(fClusterCollArray.FindObject(name));
1177  return cont;
1178 }
1179 
1180 //________________________________________________________________________
1181 TClonesArray* AliAnalysisTaskEmcal::GetParticleArray(Int_t i) const
1182 {
1183  // Get i^th TClonesArray with AliVParticle
1184 
1186  if (!cont) {
1187  AliError(Form("%s: Particle container %d not found",GetName(),i));
1188  return 0;
1189  }
1190  TString contName = cont->GetArrayName();
1191  return cont->GetArray();
1192 }
1193 
1194 //________________________________________________________________________
1195 TClonesArray* AliAnalysisTaskEmcal::GetClusterArray(Int_t i) const
1196 {
1197  // Get i^th TClonesArray with AliVCluster
1198 
1200  if (!cont) {
1201  AliError(Form("%s:Cluster container %d not found",GetName(),i));
1202  return 0;
1203  }
1204  return cont->GetArray();
1205 }
1206 
1207 //________________________________________________________________________
1208 AliVParticle* AliAnalysisTaskEmcal::GetAcceptParticleFromArray(Int_t p, Int_t c) const
1209 {
1210  // Get particle p if accepted from container c
1211  // If particle not accepted return 0
1212 
1214  if (!cont) {
1215  AliError(Form("%s: Particle container %d not found",GetName(),c));
1216  return 0;
1217  }
1218  AliVParticle *vp = cont->GetAcceptParticle(p);
1219 
1220  return vp;
1221 }
1222 
1223 //________________________________________________________________________
1224 AliVCluster* AliAnalysisTaskEmcal::GetAcceptClusterFromArray(Int_t cl, Int_t c) const
1225 {
1226  // Get particle p if accepted from container c
1227  // If particle not accepted return 0
1228 
1230  if (!cont) {
1231  AliError(Form("%s: Cluster container %d not found",GetName(),c));
1232  return 0;
1233  }
1234  AliVCluster *vc = cont->GetAcceptCluster(cl);
1235 
1236  return vc;
1237 }
1238 
1239 //________________________________________________________________________
1241 {
1242  // Get number of entries in particle array i
1243 
1245  if (!cont) {
1246  AliError(Form("%s: Particle container %d not found",GetName(),i));
1247  return 0;
1248  }
1249  return cont->GetNEntries();
1250 }
1251 
1252 //________________________________________________________________________
1254 {
1255  // Get number of entries in cluster array i
1256 
1258  if (!cont) {
1259  AliError(Form("%s: Cluster container %d not found",GetName(),i));
1260  return 0;
1261  }
1262  return cont->GetNEntries();
1263 }
1264 
1265 //________________________________________________________________________
1267 {
1268  // Get main trigger match
1269  //
1270  // For the selection of the main trigger patch, high and low threshold triggers of a given category are grouped
1271  // If there are more than 1 main patch of a given trigger category (i.e. different high and low threshold patches),
1272  // the highest one according to the ADC value is taken. In case doSimpleOffline is true, then only the patches from
1273  // the simple offline trigger are used.
1274 
1275  if (!fTriggerPatchInfo) {
1276  AliError(Form("%s: fTriggerPatchInfo not available",GetName()));
1277  return 0;
1278  }
1279 
1280  //number of patches in event
1281  Int_t nPatch = fTriggerPatchInfo->GetEntries();
1282 
1283  //extract main trigger patch(es)
1284  AliEmcalTriggerPatchInfo *patch(NULL), *selected(NULL);
1285  for (Int_t iPatch = 0; iPatch < nPatch; iPatch++) {
1286 
1287  patch = (AliEmcalTriggerPatchInfo*)fTriggerPatchInfo->At( iPatch );
1288  if (patch->IsMainTrigger()) {
1289  if(doSimpleOffline){
1290  if(patch->IsOfflineSimple()){
1291  switch(trigger){
1292  case kTriggerLevel0:
1293  // option not yet implemented in the trigger maker
1294  if(patch->IsLevel0()) selected = patch;
1295  break;
1296  case kTriggerLevel1Jet:
1297  if(patch->IsJetHighSimple() || patch->IsJetLowSimple()){
1298  if(!selected) selected = patch;
1299  else if(patch->GetADCOfflineAmp() > selected->GetADCOfflineAmp()) selected = patch;
1300  }
1301  break;
1302  case kTriggerLevel1Gamma:
1303  if(patch->IsGammaHighSimple() || patch->IsGammaLowSimple()){
1304  if(!selected) selected = patch;
1305  else if(patch->GetADCOfflineAmp() > selected->GetADCOfflineAmp()) selected = patch;
1306  }
1307  break;
1308  default: // Silence compiler warnings
1309  AliError("Untreated case: Main Patch is recalculated; should be in 'else' branch");
1310  };
1311  }
1312  } else { // Not OfflineSimple
1313  switch(trigger){
1314  case kTriggerLevel0:
1315  if(patch->IsLevel0()) selected = patch;
1316  break;
1317  case kTriggerLevel1Jet:
1318  if(patch->IsJetHigh() || patch->IsJetLow()){
1319  if(!selected) selected = patch;
1320  else if (patch->GetADCAmp() > selected->GetADCAmp())
1321  selected = patch;
1322  }
1323  break;
1324  case kTriggerLevel1Gamma:
1325  if(patch->IsGammaHigh() || patch->IsGammaLow()){
1326  if(!selected) selected = patch;
1327  else if (patch->GetADCAmp() > selected->GetADCAmp())
1328  selected = patch;
1329  }
1330  break;
1331  default:
1332  AliError("Untreated case: Main Patch is recalculated; should be in 'else' branch");
1333  };
1334  }
1335  }
1336  else if ((trigger == kTriggerRecalcJet && patch->IsRecalcJet()) ||
1337  (trigger == kTriggerRecalcGamma && patch->IsRecalcGamma())) { // recalculated patches
1338  if (doSimpleOffline && patch->IsOfflineSimple()) {
1339  if(!selected) selected = patch;
1340  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
1341  selected = patch;
1342  }
1343  else if (!doSimpleOffline && !patch->IsOfflineSimple()) {
1344  if(!selected) selected = patch;
1345  else if (patch->GetADCAmp() > selected->GetADCAmp())
1346  selected = patch;
1347  }
1348  }
1349  }
1350  return selected;
1351 }
1352 
1353 //________________________________________________________________________
1355 {
1356  // Add object to event
1357 
1358  if (!(InputEvent()->FindListObject(obj->GetName()))) {
1359  InputEvent()->AddObject(obj);
1360  } else {
1361  AliFatal(Form("%s: Container with name %s already present. Aborting", GetName(), obj->GetName()));
1362  }
1363 }
1364 
1365 //________________________________________________________________________
1367 {
1368  axis->SetBinLabel(1, "NullObject");
1369  axis->SetBinLabel(2, "Pt");
1370  axis->SetBinLabel(3, "Acceptance");
1371  axis->SetBinLabel(4, "BitMap");
1372  axis->SetBinLabel(5, "Bit4");
1373  axis->SetBinLabel(6, "Bit5");
1374  axis->SetBinLabel(7, "Bit6");
1375  axis->SetBinLabel(8, "NotHybridTrack");
1376  axis->SetBinLabel(9, "MCFlag");
1377  axis->SetBinLabel(10, "MCGenerator");
1378  axis->SetBinLabel(11, "ChargeCut");
1379  axis->SetBinLabel(12, "MinDistanceTPCSectorEdge");
1380  axis->SetBinLabel(13, "MinMCLabelAccept");
1381  axis->SetBinLabel(14, "IsEMCal");
1382  axis->SetBinLabel(15, "Time");
1383  axis->SetBinLabel(16, "Energy");
1384  axis->SetBinLabel(17, "Bit16");
1385  axis->SetBinLabel(18, "Bit17");
1386  axis->SetBinLabel(19, "Area");
1387  axis->SetBinLabel(20, "AreaEmc");
1388  axis->SetBinLabel(21, "ZLeadingCh");
1389  axis->SetBinLabel(22, "ZLeadingEmc");
1390  axis->SetBinLabel(23, "NEF");
1391  axis->SetBinLabel(24, "MinLeadPt");
1392  axis->SetBinLabel(25, "MaxTrackPt");
1393  axis->SetBinLabel(26, "MaxClusterPt");
1394  axis->SetBinLabel(27, "Flavour");
1395  axis->SetBinLabel(28, "TagStatus");
1396  axis->SetBinLabel(29, "MinNConstituents");
1397  axis->SetBinLabel(30, "Bit29");
1398  axis->SetBinLabel(31, "Bit30");
1399  axis->SetBinLabel(32, "Bit31");
1400 }
1401 
1402 //________________________________________________________________________
1403 Double_t AliAnalysisTaskEmcal::GetParallelFraction(AliVParticle* part1, AliVParticle* part2)
1404 {
1405  // Calculates the fraction of momentum of part 1 w.r.t. part 2 in the direction of part 2.
1406 
1407  TVector3 vect1(part1->Px(), part1->Py(), part1->Pz());
1408  TVector3 vect2(part2->Px(), part2->Py(), part2->Pz());
1409  Double_t z = (vect1 * vect2) / (vect2 * vect2);
1410  return z;
1411 }
1412 
1413 //________________________________________________________________________
1414 Double_t AliAnalysisTaskEmcal::GetParallelFraction(const TVector3& vect1, AliVParticle* part2)
1415 {
1416  // Calculates the fraction of momentum of vect 1 w.r.t. part 2 in the direction of part 2.
1417 
1418  TVector3 vect2(part2->Px(), part2->Py(), part2->Pz());
1419  Double_t z = (vect1 * vect2) / (vect2 * vect2);
1420  return z;
1421 }
Bool_t AcceptCluster(AliVCluster *vp)
void SetParticlePtCut(Double_t cut)
TH1 * fHistTrials
x section from pythia header
Bool_t AcceptTrack(AliVParticle *track, Int_t c=0) const
Bool_t HasTriggerType(TriggerType triggersel)
Int_t fNTrials
event pt hard bin
AliEmcalTriggerPatchInfo * GetMainTriggerPatch(TriggerCategory triggersel=kTriggerLevel1Jet, Bool_t doOfflinSimple=kFALSE)
Double_t fPtHard
event Pythia header
void SetTrackPtCut(Double_t cut, Int_t c=0)
Double_t fEPV0
event centrality bin
TList * list
virtual void SetArray(AliVEvent *event)
Bool_t AcceptCluster(AliVCluster *clus, Int_t c=0) const
void AddObjectToEvent(TObject *obj)
Int_t fCentBin
event centrality
TH1 * fHistEventsAfterSel
total number of trials per pt hard bin after selection
Bool_t AcceptParticle(AliVParticle *vp)
void SetArrayName(const char *n)
TH1 * fHistEventPlane
z vertex position
TList * fOutput
x-section from pythia header
TH1 * fHistEvents
trials from pyxsec.root
void SetClusPtCut(Double_t cut, Int_t c=0)
AliClusterContainer * AddClusterContainer(const char *n)
Double_t fEPV0C
event plane V0A
TH1 * fHistCentrality
pt hard distribution
void SetTrackEtaLimits(Double_t min, Double_t max, Int_t c=0)
AliParticleContainer * GetParticleContainer(Int_t i=0) const
TProfile * fHistXsectionAfterSel
total number of events per pt hard bin after selection
virtual Bool_t FillHistograms()
AliVCluster * GetAcceptCluster(Int_t i)
Int_t GetNParticles(Int_t i=0) const
TClonesArray * fCaloClusters
tracks
ClassImp(AliAnalysisTaskEmcal) AliAnalysisTaskEmcal
void SetArray(AliVEvent *event)
Main data structure storing all relevant information of EMCAL/DCAL trigger patches.
AliEMCALGeometry * fGeom
whether it's an ESD analysis
AliVParticle * GetAcceptParticle(Int_t i)
AliGenPythiaEventHeader * fPythiaHeader
event beam type
void SetTrackPhiLimits(Double_t min, Double_t max, Int_t c=0)
AliParticleContainer * AddParticleContainer(const char *n)
AliAnalysisUtils * fAliAnalysisUtils
AliClusterContainer * GetClusterContainer(Int_t i=0) const
virtual Bool_t FillGeneralHistograms()
TClonesArray * GetParticleArray(Int_t i=0) const
BeamType fBeamType
event vertex number of contributors
Double_t fCent
trigger patch info array
TClonesArray * GetArray() const
Int_t GetNClusters(Int_t i=0) const
Int_t fNVertCont
event vertex
static Double_t GetParallelFraction(AliVParticle *part1, AliVParticle *part2)
virtual Bool_t RetrieveEventObjects()
TProfile * fHistXsection
total number of events per pt hard bin
Bool_t PythiaInfoFromFile(const char *currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
void UserExec(Option_t *option)
AliVCaloCells * fCaloCells
clusters
const TString & GetArrayName() const
TClonesArray * GetArrayFromEvent(const char *name, const char *clname=0)
virtual Bool_t IsEventSelected()
TH1 * fHistPtHard
x section from pyxsec.root
void SetParticleEtaLimits(Double_t min, Double_t max)
Int_t fPtHardBin
event pt hard
TClonesArray * fTracks
emcal geometry
TH1 * fHistTrialsAfterSel
incoming and selected events
Bool_t fIsEsd
vertex selection (optional)
Double_t fVertex[3]
event plane V0C
TFile * file
TH1 * fHistEventRejection
event plane distribution
TClonesArray * fTriggerPatchInfo
calo triggers
TClonesArray * GetClusterArray(Int_t i=0) const
Double_t fEPV0A
event plane V0
Class to make array of trigger patch objects in AOD/ESD events.
void SetParticlePhiLimits(Double_t min, Double_t max, Double_t offset=0.)
void SetClusPtCut(Double_t cut)
AliVCaloTrigger * fCaloTriggers
cells
void SetRejectionReasonLabels(TAxis *axis)
TH1 * fHistZVertex
event centrality distribution
Int_t GetNEntries() const
void SetClusTimeCut(Double_t min, Double_t max, Int_t c=0)
Float_t fXsection
event trials
TH1 * fHistEventCount
output list
void SetClusTimeCut(Double_t min, Double_t max)
AliVParticle * GetAcceptParticleFromArray(Int_t p, Int_t c=0) const
AliVCluster * GetAcceptClusterFromArray(Int_t cl, Int_t c=0) const