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