AliPhysics  e841483 (e841483)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalJetHadEPpid.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 //-------------------------------------------------------------------------
17 // 1) Analysis Task to perform Jet-Hadron Correlations
18 // 2) Event plane dependence task.
19 // 3) performs event plane resolution calculation
20 // 4) does PID of the associated pi/k/p hadrons
21 //
22 // Author: Joel Mazer (joel.mazer@cern.ch)
23 //-------------------------------------------------------------------------
24 
25 // task head include
26 #include "AliAnalysisTaskEmcalJetHadEPpid.h"
27 
28 // general ROOT includes
29 #include <TCanvas.h>
30 #include <TMath.h>
31 #include <TChain.h>
32 #include <TClonesArray.h>
33 #include <TH1F.h>
34 #include <TH2F.h>
35 #include <TH3F.h>
36 #include <TH1.h>
37 #include <TH2.h>
38 #include <TH3.h>
39 #include <TProfile.h>
40 #include <THnSparse.h>
41 #include <TList.h>
42 #include <TLorentzVector.h>
43 #include <TParameter.h>
44 #include <TParticle.h>
45 #include <TTree.h>
46 #include <TVector3.h>
47 #include <TObjArray.h>
48 
49 // AliROOT includes
50 #include "AliAODEvent.h"
51 #include "AliESDEvent.h"
52 #include "AliAnalysisManager.h"
53 #include "AliAnalysisTask.h"
54 #include "AliCentrality.h"
55 #include "AliEmcalJet.h"
56 #include "AliAODJet.h"
57 #include "AliVCluster.h"
58 #include "AliVTrack.h"
59 #include <AliVEvent.h>
60 #include <AliVParticle.h>
61 #include "AliRhoParameter.h"
62 #include "AliEmcalParticle.h"
63 
64 // Localized Rho includes
65 #include "AliLocalRhoParameter.h"
66 #include "AliAnalysisTaskLocalRho.h"
67 
68 // event handler (and pico's) includes
69 #include <AliInputEventHandler.h>
70 #include <AliVEventHandler.h>
71 #include "AliESDInputHandler.h"
72 #include "AliESDCaloCluster.h"
73 #include "AliPicoTrack.h"
74 #include "AliEventPoolManager.h"
75 #include "AliESDtrackCuts.h"
76 
77 // PID includes
78 #include "AliPIDResponse.h"
79 #include "AliTPCPIDResponse.h"
80 #include "AliESDpid.h"
81 
82 // magnetic field includes
83 #include "TGeoGlobalMagField.h"
84 #include "AliMagF.h"
85 
86 // container includes
87 #include "AliJetContainer.h"
88 #include "AliParticleContainer.h"
89 #include "AliClusterContainer.h"
90 #include "AliTrackContainer.h"
91 
92 using std::cout;
93 using std::endl;
94 
95 ClassImp(AliAnalysisTaskEmcalJetHadEPpid)
96 
97 //________________________________________________________________________
98 AliAnalysisTaskEmcalJetHadEPpid::AliAnalysisTaskEmcalJetHadEPpid() :
99  AliAnalysisTaskEmcalJet("correlations",kFALSE),
100  fPhimin(-10), fPhimax(10),
101  fEtamin(-0.9), fEtamax(0.9),
102  fAreacut(0.0), fTrkBias(5), fClusBias(5), fTrkEta(0.9),
103  fJetPtcut(15.0), fJetRad(0.4), fConstituentCut(3.0),
104  fesdTrackCuts(0),
105  fDetectorType(kVZEROComb),
106  fSoftTrackMinPt_ep(0.15), fSoftTrackMaxPt_ep(5.),
107  fNAcceptedTracks(0), fLeadingJet(0), fExcludeLeadingJetsFromFit(1.),
108  fCentralityClasses(0), fInCentralitySelection(-1),
109  fDoEventMixing(0), fMixingTracks(50000), fNMIXtracks(5000), fNMIXevents(5),
110  fCentBinSize(1), fReduceStatsCent(0),
111  fTriggerEventType(AliVEvent::kAny), fMixingEventType(AliVEvent::kAny),
112  fDoEffCorr(0), fcorrJetPt(0),
113  doPlotGlobalRho(0), doVariableBinning(0), dovarbinTHnSparse(0),
114  makeQAhistos(0), makeBIAShistos(0), makeextraCORRhistos(0), makeoldJEThadhistos(0),
115  allpidAXIS(0), fcutType("EMCAL"), doPID(0), doPIDtrackBIAS(0), doaltPIDbinning(0),
116  doEventPlaneRes(0),
117  doComments(0),
118  doFlavourJetAnalysis(0), fJetFlavTag(3),
119  douseOLDtrackFramework(kFALSE),
120  fBeam(kNA),
121  fLocalRhoVal(0),
122  fTracksName(""), fTracksNameME(""), fJetsName(""), fCaloClustersName(""),
123  event(0),
124  isPItpc(0), isKtpc(0), isPtpc(0), // pid TPC
125  isPIits(0), isKits(0), isPits(0), // pid ITS
126  isPItof(0), isKtof(0), isPtof(0), // pid TOF
127  fPoolMgr(0x0),
128  fPIDResponse(0x0), fTPCResponse(),
129  fESD(0), fAOD(0), fVevent(0),
130  fHistEventQA(0), fHistEventSelectionQA(0), fHistEventSelectionQAafterCuts(0),
131  fHistCentZvertGA(0), fHistCentZvertJE(0), fHistCentZvertMB(0), fHistCentZvertAny(0),
132  fHistTPCdEdX(0), fHistITSsignal(0), //fHistTOFsignal(0),
133  fHistRhovsCent(0), fHistNjetvsCent(0),
134  fHistNTrackPtNEW(0), fHistNTrackPhiNEW(0), fHistNTrackEtaNEW(0), fHistNTrackPhiEtaNEW(0),
135  fHistNTrackPt(0), fHistNTrackPhi(0), fHistNTrackEta(0), fHistNTrackPhiEta(0),
136  fHistNJetPt(0), fHistNJetPhi(0), fHistNJetEta(0), fHistNJetPhiEta(0),
137  fHistMult(0),
138  fHistJetPhi(0), fHistTrackPhi(0),
139  fHistLocalRhoJetpt(0),
140  fHistJetHaddPhiIN(0), fHistJetHaddPhiOUT(0), fHistJetHaddPhiMID(0),
141  fHistJetHaddPhiBias(0), fHistJetHaddPhiINBias(0), fHistJetHaddPhiOUTBias(0), fHistJetHaddPhiMIDBias(0),
142  fHistMEdPHI(0),
143  fHistTrackPtallcent(0),
144  fHistJetEtaPhi(0), fHistJetHEtaPhi(0),
145  fHistSEphieta(0), fHistMEphieta(0),
146  fHistJetHaddPHI(0),
147  fHistClusEtaPhiEnergy(0), fHistClusEnergy(0), fHistClusofJetEnergy(0),
148  fHistJetNClusterConstit(0), fHistJetNTrackConstit(0), fHistJetNConstit(0),
149  fHistPID(0),
150  fhnPID(0x0), fhnMixedEvents(0x0), fhnJH(0x0), fhnCorr(0x0),
151  fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0), fUseChiWeightForVZERO(kTRUE), // test
152  fTracksFromContainer(0),
153  fJetsCont(0), fTracksCont(0), fCaloClustersCont(0),
154  fContainerAllJets(0), fContainerPIDJets(1)//, fTrgJet(0)
155 {
156  // Default Constructor
157  for(Int_t ilab=0; ilab<4; ilab++){
158  for(Int_t ipta=0; ipta<7; ipta++){
159  //fHistTrackEtaPhi[ilab][ipta]=0; // keep out for now
160  } // end of pt-associated loop
161  } // end of lab loop
162 
163  for(Int_t cen = 0; cen<10; cen++){
164  fProfV2Resolution[cen]=0;
165  fProfV3Resolution[cen]=0;
166  fProfV4Resolution[cen]=0;
167  fProfV5Resolution[cen]=0;
168  }
169 
170  for(Int_t itrackpt=0; itrackpt<9; itrackpt++){
171  fHistJetHadbindPhi[itrackpt]=0;
172  fHistJetHadbindPhiIN[itrackpt]=0;
173  fHistJetHadbindPhiMID[itrackpt]=0;
174  fHistJetHadbindPhiOUT[itrackpt]=0;
175  } // end of trackpt bin loop
176 
177  for(Int_t icent = 0; icent<6; ++icent){
178  for(Int_t iptjet = 0; iptjet<5; ++iptjet){
179  for(Int_t ieta = 0; ieta<3; ++ieta){
180  fHistJetH[icent][iptjet][ieta]=0;
181  fHistJetHBias[icent][iptjet][ieta]=0;
182  fHistJetHTT[icent][iptjet][ieta]=0;
183  } // end of eta loop
184  } // end of pt-jet loop
185  } // end of centrality loop
186 
187  // centrality dependent histo's
188  for (Int_t i = 0;i<6;++i){
189  fHistJetPt[i] = 0;
190  fHistJetPtBias[i] = 0;
191  fHistJetPtTT[i] = 0;
192  fHistAreavsRawPt[i] = 0;
193  fHistJetPtvsTrackPt[i] = 0;
194  fHistTrackPt[i] = 0;
195  fHistEP0[i] = 0;
196  fHistEP0A[i] = 0;
197  fHistEP0C[i] = 0;
198  fHistEPAvsC[i] = 0;
199  fHistJetPtcorrGlRho[i] = 0;
200  fHistJetPtvsdEP[i] = 0;
201  fHistJetPtvsdEPBias[i] = 0;
202  fHistRhovsdEP[i] = 0;
203  fHistJetEtaPhiPt[i] = 0;
204  fHistJetEtaPhiPtBias[i] = 0;
205  fHistJetPtArea[i] = 0;
206  fHistJetPtAreaBias[i] = 0;
207  fHistJetPtNcon[i] = 0;
208  fHistJetPtNconBias[i] = 0;
209  fHistJetPtNconCh[i] = 0;
210  fHistJetPtNconBiasCh[i] = 0;
211  fHistJetPtNconEm[i] = 0;
212  fHistJetPtNconBiasEm[i] = 0;
213  fHistJetHaddPhiINcent[i] = 0;
214  fHistJetHaddPhiOUTcent[i] = 0;
215  fHistJetHaddPhiMIDcent[i] = 0;
216  }
217 
218  SetMakeGeneralHistograms(kTRUE);
219 
220 }
221 
222 //________________________________________________________________________
223 AliAnalysisTaskEmcalJetHadEPpid::AliAnalysisTaskEmcalJetHadEPpid(const char *name) :
224  AliAnalysisTaskEmcalJet(name,kTRUE),
225  fPhimin(-10), fPhimax(10),
226  fEtamin(-0.9), fEtamax(0.9),
227  fAreacut(0.0), fTrkBias(5), fClusBias(5), fTrkEta(0.9),
228  fJetPtcut(15.0), fJetRad(0.4), fConstituentCut(3.0),
229  fesdTrackCuts(0),
230  fDetectorType(kVZEROComb),
231  fSoftTrackMinPt_ep(0.15), fSoftTrackMaxPt_ep(5.),
232  fNAcceptedTracks(0), fLeadingJet(0), fExcludeLeadingJetsFromFit(1.),
233  fCentralityClasses(0), fInCentralitySelection(-1),
234  fDoEventMixing(0), fMixingTracks(50000), fNMIXtracks(5000), fNMIXevents(5),
235  fCentBinSize(1), fReduceStatsCent(0),
236  fTriggerEventType(AliVEvent::kAny), fMixingEventType(AliVEvent::kAny),
237  fDoEffCorr(0), fcorrJetPt(0),
238  doPlotGlobalRho(0), doVariableBinning(0), dovarbinTHnSparse(0),
239  makeQAhistos(0), makeBIAShistos(0), makeextraCORRhistos(0), makeoldJEThadhistos(0),
240  allpidAXIS(0), fcutType("EMCAL"), doPID(0), doPIDtrackBIAS(0), doaltPIDbinning(0),
241  doEventPlaneRes(0),
242  doComments(0),
243  doFlavourJetAnalysis(0), fJetFlavTag(3),
244  douseOLDtrackFramework(kFALSE),
245  fBeam(kNA),
246  fLocalRhoVal(0),
247  fTracksName(""), fTracksNameME(""), fJetsName(""), fCaloClustersName(""),
248  event(0),
249  isPItpc(0), isKtpc(0), isPtpc(0), // pid TPC
250  isPIits(0), isKits(0), isPits(0), // pid ITS
251  isPItof(0), isKtof(0), isPtof(0), // pid TOF
252  fPoolMgr(0x0),
253  fPIDResponse(0x0), fTPCResponse(),
254  fESD(0), fAOD(0), fVevent(0),
255  fHistEventQA(0), fHistEventSelectionQA(0), fHistEventSelectionQAafterCuts(0),
256  fHistCentZvertGA(0), fHistCentZvertJE(0), fHistCentZvertMB(0), fHistCentZvertAny(0),
257  fHistTPCdEdX(0), fHistITSsignal(0), //fHistTOFsignal(0),
258  fHistRhovsCent(0), fHistNjetvsCent(0),
259  fHistNTrackPtNEW(0), fHistNTrackPhiNEW(0), fHistNTrackEtaNEW(0), fHistNTrackPhiEtaNEW(0),
260  fHistNTrackPt(0), fHistNTrackPhi(0), fHistNTrackEta(0), fHistNTrackPhiEta(0),
261  fHistNJetPt(0), fHistNJetPhi(0), fHistNJetEta(0), fHistNJetPhiEta(0),
262  fHistMult(0),
263  fHistJetPhi(0), fHistTrackPhi(0),
264  fHistLocalRhoJetpt(0),
265  fHistJetHaddPhiIN(0), fHistJetHaddPhiOUT(0), fHistJetHaddPhiMID(0),
266  fHistJetHaddPhiBias(0), fHistJetHaddPhiINBias(0), fHistJetHaddPhiOUTBias(0), fHistJetHaddPhiMIDBias(0),
267  fHistMEdPHI(0),
268  fHistTrackPtallcent(0),
269  fHistJetEtaPhi(0), fHistJetHEtaPhi(0),
270  fHistSEphieta(0), fHistMEphieta(0),
271  fHistJetHaddPHI(0),
272  fHistClusEtaPhiEnergy(0), fHistClusEnergy(0), fHistClusofJetEnergy(0),
273  fHistJetNClusterConstit(0), fHistJetNTrackConstit(0), fHistJetNConstit(0),
274  fHistPID(0),
275  fhnPID(0x0), fhnMixedEvents(0x0), fhnJH(0x0), fhnCorr(0x0),
276  fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0), fUseChiWeightForVZERO(kTRUE), // test
277  fTracksFromContainer(0),
278  fJetsCont(0), fTracksCont(0), fCaloClustersCont(0),
279  fContainerAllJets(0), fContainerPIDJets(1)//, fTrgJet(0)
280 {
281  // Default Constructor
282  for(Int_t ilab=0; ilab<4; ilab++){
283  for(Int_t ipta=0; ipta<7; ipta++){
284  //fHistTrackEtaPhi[ilab][ipta]=0; //keep out for now
285  } // end of pt-associated loop
286  } // end of lab loop
287 
288  for(Int_t cen = 0; cen<10; cen++){
289  fProfV2Resolution[cen]=0;
290  fProfV3Resolution[cen]=0;
291  fProfV4Resolution[cen]=0;
292  fProfV5Resolution[cen]=0;
293  }
294 
295  for(Int_t itrackpt=0; itrackpt<9; itrackpt++){
296  fHistJetHadbindPhi[itrackpt]=0;
297  fHistJetHadbindPhiIN[itrackpt]=0;
298  fHistJetHadbindPhiMID[itrackpt]=0;
299  fHistJetHadbindPhiOUT[itrackpt]=0;
300  } // end of trackpt bin loop
301 
302  for(Int_t icent = 0; icent<6; ++icent){
303  for(Int_t iptjet = 0; iptjet<5; ++iptjet){
304  for(Int_t ieta = 0; ieta<3; ++ieta){
305  fHistJetH[icent][iptjet][ieta]=0;
306  fHistJetHBias[icent][iptjet][ieta]=0;
307  fHistJetHTT[icent][iptjet][ieta]=0;
308  } // end of eta loop
309  } // end of pt-jet loop
310  } // end of centrality loop
311 
312  // centrality dependent histo's
313  for (Int_t i = 0;i<6;++i){
314  fHistJetPt[i] = 0;
315  fHistJetPtBias[i] = 0;
316  fHistJetPtTT[i] = 0;
317  fHistAreavsRawPt[i] = 0;
318  fHistJetPtvsTrackPt[i] = 0;
319  fHistTrackPt[i] = 0;
320  fHistEP0[i] = 0;
321  fHistEP0A[i] = 0;
322  fHistEP0C[i] = 0;
323  fHistEPAvsC[i] = 0;
324  fHistJetPtcorrGlRho[i] = 0;
325  fHistJetPtvsdEP[i] = 0;
326  fHistJetPtvsdEPBias[i] = 0;
327  fHistRhovsdEP[i] = 0;
328  fHistJetEtaPhiPt[i] = 0;
329  fHistJetEtaPhiPtBias[i] = 0;
330  fHistJetPtArea[i] = 0;
331  fHistJetPtAreaBias[i] = 0;
332  fHistJetPtNcon[i] = 0;
333  fHistJetPtNconBias[i] = 0;
334  fHistJetPtNconCh[i] = 0;
335  fHistJetPtNconBiasCh[i] = 0;
336  fHistJetPtNconEm[i] = 0;
337  fHistJetPtNconBiasEm[i] = 0;
338  fHistJetHaddPhiINcent[i] = 0;
339  fHistJetHaddPhiOUTcent[i] = 0;
340  fHistJetHaddPhiMIDcent[i] = 0;
341  }
342 
343  SetMakeGeneralHistograms(kTRUE);
344 
345 }
346 
347 //_______________________________________________________________________
348 AliAnalysisTaskEmcalJetHadEPpid::~AliAnalysisTaskEmcalJetHadEPpid()
349 {
350  // destructor
351  if(fOutput) {delete fOutput; fOutput = 0;}
352  if(fCentralityClasses) {delete fCentralityClasses; fCentralityClasses = 0x0;} // test
353  if(fChi2A) {delete fChi2A; fChi2A = 0x0;}
354  if(fChi2C) {delete fChi2C; fChi2C = 0x0;}
355  if(fChi3A) {delete fChi3A; fChi3A = 0x0;}
356  if(fChi3C) {delete fChi3C; fChi3C = 0x0;}
357 }
358 
359 //________________________________________________________________________
360 void AliAnalysisTaskEmcalJetHadEPpid::UserCreateOutputObjects()
361 { // This is called ONCE!
362  if (!fCreateHisto) return;
363  AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
364  OpenFile(1); // do I need the 1?
365 
366  // char array for naming histograms
367  int nchar = 200;
368  char name[nchar];
369 
370  // strings for titles
371  TString name1;
372  TString title1;
373 
374  // Centrality and Zvertex distribution for various triggers - Event Mixing QA
375  fHistCentZvertGA = new TH2F("fHistCentZvertGA", "Centrality - Z-vertex distribution for GA trigger", 20, 0, 100, 10, -10, 10);
376  fOutput->Add(fHistCentZvertGA);
377  fHistCentZvertJE = new TH2F("fHistCentZvertJE", "Centrality - Z-vertex distribution for JE trigger", 20, 0, 100, 10, -10, 10);
378  fOutput->Add(fHistCentZvertJE);
379  fHistCentZvertMB = new TH2F("fHistCentZvertMB", "Centrality - Z-vertex distribution for MB trigger", 20, 0, 100, 10, -10, 10);
380  fOutput->Add(fHistCentZvertMB);
381  fHistCentZvertAny = new TH2F("fHistCentZvertAny", "Centrality - Z-vertex distribution for kAny trigger", 20, 0, 100, 10, -10, 10);
382  fOutput->Add(fHistCentZvertAny);
383 
384  // Event QA histo
385  fHistEventQA = new TH1F("fHistEventQA", "Event Counter at checkpoints in code", 20, 0.5, 20.5);
386  SetfHistQAcounterLabels(fHistEventQA);
387  fOutput->Add(fHistEventQA);
388 
389  // Event Selection QA histo
390  fHistEventSelectionQA = new TH1F("fHistEventSelectionQA", "Trigger Selection Counter", 20, 0.5, 20.5);
391  //SetfHistEvtSelQALabels(fHistEventSelectionQA);
392  fOutput->Add(fHistEventSelectionQA);
393 
394  // Event Selection QA histo after all cuts to return from task
395  fHistEventSelectionQAafterCuts = new TH1F("fHistEventSelectionQAafterCuts", "Trigger Selection Counter after Cuts", 20, 0.5, 20.5);
396  //SetfHistEvtSelQALabels(fHistEventSelectionQAafterCuts);
397  fOutput->Add(fHistEventSelectionQAafterCuts);
398 
399  // create histograms that arn't array
400  fHistNjetvsCent = new TH2F("NjetvsCent", "NjetvsCent", 100, 0.0, 100.0, 100, 0, 100);
401  fHistJetHaddPHI = new TH1F("fHistJetHaddPHI", "Jet-Hadron #Delta#varphi", 72,-0.5*TMath::Pi(),1.5*TMath::Pi());
402  fHistJetHaddPhiIN = new TH1F("fHistJetHaddPhiIN","Jet-Hadron #Delta#varphi IN PLANE", 72,-0.5*TMath::Pi(), 1.5*TMath::Pi());
403  fHistJetHaddPhiOUT = new TH1F("fHistJetHaddPhiOUT","Jet-Hadron #Delta#varphi OUT PLANE", 72,-0.5*TMath::Pi(), 1.5*TMath::Pi());
404  fHistJetHaddPhiMID = new TH1F("fHistJetHaddPhiMID","Jet-Hadron #Delta#varphi MIDDLE of PLANE", 72,-0.5*TMath::Pi(), 1.5*TMath::Pi());
405  fHistLocalRhoJetpt = new TH1F("fHistLocalRhoJetpt","Local Rho corrected Jet p_{T}", 50, -50, 200);
406 
407  // add to output lists
408  fOutput->Add(fHistNjetvsCent);
409  fOutput->Add(fHistJetHaddPHI);
410  fOutput->Add(fHistJetHaddPhiIN);
411  fOutput->Add(fHistJetHaddPhiOUT);
412  fOutput->Add(fHistJetHaddPhiMID);
413  fOutput->Add(fHistLocalRhoJetpt);
414 
415  // TPC signal (PID)
416  fHistTPCdEdX = new TH2F("TPCdEdX", "TPCdEdX", 400, 0.0, 20.0, 500, 0, 500);
417  fOutput->Add(fHistTPCdEdX);
418 
419  // TEST HISTOS - making sure container tracks are filtered properly
420  fHistNTrackPhiNEW = new TH1F("fHistNTrackPhiNEW", "NTrack vs #Psi (new framework)", 144, 0, 2*TMath::Pi());
421  fHistNTrackPtNEW = new TH1F("fHistNTrackPtNEW", "NTrack vs p_T (new framework)", 500, 0, 50);
422  fHistNTrackEtaNEW = new TH1F("fHistNTrackEtaNEW", "NTrack vs #eta (new framework)", 180, -0.9, 0.9);
423  fHistNTrackPhiEtaNEW = new TH2F("fHistNTrackPhiEtaNEW", "NTrack vs #phi vs #eta (new framework)", 72, 0, 2*TMath::Pi(), 90, -0.9, 0.9);
424  fOutput->Add(fHistNTrackPhiNEW);
425  fOutput->Add(fHistNTrackPtNEW);
426  fOutput->Add(fHistNTrackEtaNEW);
427  fOutput->Add(fHistNTrackPhiEtaNEW);
428 
429  // some default histos for QA - Added April 28, 2016
430  // tracks
431  fHistNTrackPt = new TH1F("fHistNTrackPt", "NTrack vs p_T", 500, 0, 50);
432  fHistNTrackPhi = new TH1F("fHistNTrackPhi", "NTrack vs #phi", 144, 0, 2*TMath::Pi());
433  fHistNTrackEta = new TH1F("fHistNTrackEta", "NTrack vs #eta", 180, -0.9, 0.9);
434  fHistNTrackPhiEta = new TH2F("fHistNTrackPhiEta", "NTrack vs #phi vs #eta", 72, 0, 2*TMath::Pi(), 90, -0.9, 0.9);
435  fOutput->Add(fHistNTrackPt);
436  fOutput->Add(fHistNTrackPhi);
437  fOutput->Add(fHistNTrackEta);
438  fOutput->Add(fHistNTrackPhiEta);
439 
440  // jets
441  fHistNJetPt = new TH1F("fHistNJetPt", "NJet vs p_T", 500, 0, 250);
442  fHistNJetPhi = new TH1F("fHistNJetPhi", "NJet vs #phi", 144, 0, 2*TMath::Pi());
443  fHistNJetEta = new TH1F("fHistNJetEta", "NJet vs #eta", 180, -0.9, 0.9);
444  fHistNJetPhiEta = new TH2F("fHistNJetPhiEta", "NJet vs #phi vs #eta", 72, 0, 2*TMath::Pi(), 90, -0.9, 0.9);
445  fOutput->Add(fHistNJetPt);
446  fOutput->Add(fHistNJetPhi);
447  fOutput->Add(fHistNJetEta);
448  fOutput->Add(fHistNJetPhiEta);
449 
450  if(doEventPlaneRes){
451  // Reaction Plane resolution as function of centrality - corrected for 2nd order event plane
452  for (Int_t i = 0; i<10; ++i){
453  fProfV2Resolution[i] = new TProfile(Form("fProfV2Resolution_%i", i), Form("fProfV2Resolution_%i", i), 11, -0.5, 10.5);
454  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(3, "<cos(2(#Psi_{VZEROA} - #Psi_{VZEROC}))>");
455  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(4, "<cos(2(#Psi_{VZEROC} - #Psi_{VZEROA}))>");
456  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(5, "<cos(2(#Psi_{VZEROA} - #Psi_{TPC}))>");
457  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(6, "<cos(2(#Psi_{TPC} - #Psi_{VZEROA}))>");
458  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(7, "<cos(2(#Psi_{VZEROC} - #Psi_{TPC}))>");
459  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(8, "<cos(2(#Psi_{TPC} - #Psi_{VZEROC}))>");
460  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(9, "<cos(2(#Psi_{VZERO} - #Psi_{TPC_A}))>");
461  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(10, "<cos(2(#Psi_{VZERO} - #Psi_{TPC_B}))>");
462  fProfV2Resolution[i]->GetXaxis()->SetBinLabel(11, "<cos(2(#Psi_{TPC_A} - #Psi_{TPC_B}))>");
463  fOutput->Add(fProfV2Resolution[i]);
464  fProfV3Resolution[i] = new TProfile(Form("fProfV3Resolution_%i", i), Form("fProfV3Resolution_%i", i), 11, -0.5, 10.5);
465  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(3, "<cos(3(#Psi_{VZEROA} - #Psi_{VZEROC}))>");
466  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(4, "<cos(3(#Psi_{VZEROC} - #Psi_{VZEROA}))>");
467  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(5, "<cos(3(#Psi_{VZEROA} - #Psi_{TPC}))>");
468  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(6, "<cos(3(#Psi_{TPC} - #Psi_{VZEROA}))>");
469  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(7, "<cos(3(#Psi_{VZEROC} - #Psi_{TPC}))>");
470  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(8, "<cos(3(#Psi_{TPC} - #Psi_{VZEROC}))>");
471  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(9, "<cos(3(#Psi_{VZERO} - #Psi_{TPC_A}))>");
472  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(10, "<cos(3(#Psi_{VZERO} - #Psi_{TPC_B}))>");
473  fProfV3Resolution[i]->GetXaxis()->SetBinLabel(11, "<cos(3(#Psi_{TPC_A} - #Psi_{TPC_B}))>");
474  fOutput->Add(fProfV3Resolution[i]);
475  fProfV4Resolution[i] = new TProfile(Form("fProfV4Resolution_%i", i), Form("fProfV4Resolution_%i", i), 11, -0.5, 10.5);
476  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(3, "<cos(4(#Psi_{VZEROA} - #Psi_{VZEROC}))>");
477  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(4, "<cos(4(#Psi_{VZEROC} - #Psi_{VZEROA}))>");
478  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(5, "<cos(4(#Psi_{VZEROA} - #Psi_{TPC}))>");
479  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(6, "<cos(4(#Psi_{TPC} - #Psi_{VZEROA}))>");
480  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(7, "<cos(4(#Psi_{VZEROC} - #Psi_{TPC}))>");
481  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(8, "<cos(4(#Psi_{TPC} - #Psi_{VZEROC}))>");
482  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(9, "<cos(4(#Psi_{VZERO} - #Psi_{TPC_A}))>");
483  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(10, "<cos(4(#Psi_{VZERO} - #Psi_{TPC_B}))>");
484  fProfV4Resolution[i]->GetXaxis()->SetBinLabel(11, "<cos(4(#Psi_{TPC_A} - #Psi_{TPC_B}))>");
485  fOutput->Add(fProfV4Resolution[i]);
486  fProfV5Resolution[i] = new TProfile(Form("fProfV5Resolution_%i", i), Form("fProfV5Resolution_%i", i), 11, -0.5, 10.5);
487  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(3, "<cos(5(#Psi_{VZEROA} - #Psi_{VZEROC}))>");
488  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(4, "<cos(5(#Psi_{VZEROC} - #Psi_{VZEROA}))>");
489  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(5, "<cos(5(#Psi_{VZEROA} - #Psi_{TPC}))>");
490  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(6, "<cos(5(#Psi_{TPC} - #Psi_{VZEROA}))>");
491  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(7, "<cos(5(#Psi_{VZEROC} - #Psi_{TPC}))>");
492  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(8, "<cos(5(#Psi_{TPC} - #Psi_{VZEROC}))>");
493  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(9, "<cos(5(#Psi_{VZERO} - #Psi_{TPC_A}))>");
494  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(10, "<cos(5(#Psi_{VZERO} - #Psi_{TPC_B}))>");
495  fProfV5Resolution[i]->GetXaxis()->SetBinLabel(11, "<cos(5(#Psi_{TPC_A} - #Psi_{TPC_B}))>");
496  fOutput->Add(fProfV5Resolution[i]);
497  }
498  }
499 
500  // create histo's used for general QA
501  if (makeQAhistos) {
502  //fHistTPCdEdX = new TH2F("TPCdEdX", "TPCdEdX", 2000, 0.0, 100.0, 500, 0, 500);
503  fHistITSsignal = new TH2F("ITSsignal", "ITSsignal", 2000, 0.0, 100.0, 500, 0, 500);
504  // fHistTOFsignal = new TH2F("TOFsignal", "TOFsignal", 2000, 0.0, 100.0, 500, 0, 500);
505  fHistJetPhi = new TH1F("fHistJetPhi", "Jet #phi Distribution", 72, 0., 1.0*TMath::Pi());
506  fHistTrackPhi = new TH1F("fHistTrackPhi", "Track #phi Distribution", 72, 0., 2.0*TMath::Pi());
507  fHistRhovsCent = new TH2F("RhovsCent", "RhovsCent", 100, 0.0, 100.0, 500, 0, 500);
508  fHistTrackPtallcent = new TH1F("fHistTrackPtallcent", "p_{T} distribution", 1000, 0.0, 100.0);
509  fHistJetEtaPhi = new TH2F("fHistJetEtaPhi", "Jet #eta-#phi",512,-1.8,1.8,512,-3.2,3.2);
510  fHistJetHEtaPhi = new TH2F("fHistJetHEtaPhi", "Jet-Hadron #Delta#eta-#Delta#phi", 72, -1.8, 1.8, 72, -1.6, 4.8);
511  fHistSEphieta = new TH2F("fHistSEphieta", "Single Event #phi-#eta distribution", 64,-0.5*TMath::Pi(), 1.5*TMath::Pi(), 64,-1.8,1.8); // was 64 bins
512  fHistMEphieta = new TH2F("fHistMEphieta", "Mixed Event #phi-#eta distribution", 64, -0.5*TMath::Pi(), 1.5*TMath::Pi(), 64,-1.8,1.8); // was 64 bins
513 
514  fHistClusEtaPhiEnergy = new TH3F("fHistClusEtaPhiEnergy", "Cluster eta-phi-energy", 900, -1.8, 1.8, 720, -3.2, 3.2, 20, 0, 20);
515  fHistClusEnergy = new TH1F("fHistClusEnergy", "Cluster Energy distribution", 200, 0, 20);
516  fHistClusofJetEnergy = new TH1F("fHistClusofJetEnergy", "Cluster of Jet Energy distribution", 200, 0, 20);
517  fHistJetNClusterConstit = new TH1F("fHistJetNClusterConstit", "Number of Cluster Constituents in Jet", 25, 0, 25);
518  fHistJetNTrackConstit = new TH1F("fHistJetNTrackConstit", "Number of Track Constituents in Jet", 25, 0, 25);
519  fHistJetNConstit = new TH1F("fHistJetNConstit", "Number of Total Constituents in Jet", 25, 0, 25);
520 
521  // add to output list
522  //fOutput->Add(fHistTPCdEdX);
523  fOutput->Add(fHistITSsignal);
524  //fOutput->Add(fHistTOFsignal);
525  fOutput->Add(fHistJetPhi);
526  fOutput->Add(fHistTrackPhi);
527  //fOutput->Add(fHistTrackEtaPhi);
528  fOutput->Add(fHistTrackPtallcent);
529  fOutput->Add(fHistJetEtaPhi);
530  fOutput->Add(fHistJetHEtaPhi);
531  fOutput->Add(fHistSEphieta);
532  fOutput->Add(fHistMEphieta);
533  fOutput->Add(fHistClusEtaPhiEnergy);
534  fOutput->Add(fHistClusEnergy);
535  fOutput->Add(fHistClusofJetEnergy);
536  fOutput->Add(fHistJetNClusterConstit);
537  fOutput->Add(fHistJetNTrackConstit);
538  fOutput->Add(fHistJetNConstit);
539 
540  //for(Int_t ipta=0; ipta<7; ++ipta){
541  // for(Int_t ilab=0; ilab<4; ++ilab){
542  // snprintf(name, nchar, "fHistTrackEtaPhi_%i_%i", ilab,ipta);
543  // fHistTrackEtaPhi[ilab][ipta] = new TH2F(name,name,400,-1,1,640,0.0,2.*TMath::Pi());
544  // fOutput->Add(fHistTrackEtaPhi[ilab][ipta]);
545  // } // end of lab loop
546  //} // end of pt-associated loop
547 
548  for (Int_t i = 0;i<6;++i){
549  name1 = TString(Form("EP0_%i",i));
550  title1 = TString(Form("EP V0 cent bin %i",i));
551  fHistEP0[i] = new TH1F(name1,title1,144,-TMath::Pi(),TMath::Pi());
552  fOutput->Add(fHistEP0[i]);
553 
554  name1 = TString(Form("EP0A_%i",i));
555  title1 = TString(Form("EP V0A cent bin %i",i));
556  fHistEP0A[i] = new TH1F(name1,title1,144,-TMath::Pi(),TMath::Pi());
557  fOutput->Add(fHistEP0A[i]);
558 
559  name1 = TString(Form("EP0C_%i",i));
560  title1 = TString(Form("EP V0C cent bin %i",i));
561  fHistEP0C[i] = new TH1F(name1,title1,144,-TMath::Pi(),TMath::Pi());
562  fOutput->Add(fHistEP0C[i]);
563 
564  name1 = TString(Form("EPAvsC_%i",i));
565  title1 = TString(Form("EP V0A vs V0C cent bin %i",i));
566  fHistEPAvsC[i] = new TH2F(name1,title1,144,-TMath::Pi(),TMath::Pi(),144,-TMath::Pi(),TMath::Pi());
567  fOutput->Add(fHistEPAvsC[i]);
568 
569  name1 = TString(Form("JetPtvsTrackPt_%i",i));
570  title1 = TString(Form("Jet p_{T} vs Leading Track p_{T} cent bin %i",i));
571  fHistJetPtvsTrackPt[i] = new TH2F(name1,title1,250,-50,200,250,0,50);
572  fOutput->Add(fHistJetPtvsTrackPt[i]);
573 
574  name1 = TString(Form("TrackPt_%i",i));
575  title1 = TString(Form("Track p_{T} cent bin %i",i));
576  fHistTrackPt[i] = new TH1F(name1,title1,1000,0,100); // up to 200?
577  fOutput->Add(fHistTrackPt[i]);
578 
579  name1 = TString(Form("JetPtcorrGLrho_%i",i));
580  title1 = TString(Form("Jet p_{T} corrected with Global #rho cent bin %i",i));
581  fHistJetPtcorrGlRho[i] = new TH1F(name1,title1,300,-100,200); // up to 200?
582  fOutput->Add(fHistJetPtcorrGlRho[i]);
583 
584  name1 = TString(Form("JetPtvsdEP_%i",i));
585  title1 = TString(Form("Jet p_{T} vs #DeltaEP cent bin %i",i));
586  fHistJetPtvsdEP[i] = new TH2F(name1,title1,250,-50,200,288,-2*TMath::Pi(),2*TMath::Pi());
587  fOutput->Add(fHistJetPtvsdEP[i]);
588 
589  name1 = TString(Form("RhovsdEP_%i",i));
590  title1 = TString(Form("#rho vs #DeltaEP cent bin %i",i));
591  fHistRhovsdEP[i] = new TH2F(name1,title1,500,0,500,288,-2*TMath::Pi(),2*TMath::Pi());
592  fOutput->Add(fHistRhovsdEP[i]);
593 
594  name1 = TString(Form("JetEtaPhiPt_%i",i));
595  title1 = TString(Form("Jet #eta-#phi p_{T} cent bin %i",i));
596  fHistJetEtaPhiPt[i] = new TH3F(name1,title1,250,-50,200,100,-1,1,64,-3.2,3.2);
597  fOutput->Add(fHistJetEtaPhiPt[i]);
598 
599  name1 = TString(Form("JetPtArea_%i",i));
600  title1 = TString(Form("Jet p_{T} Area cent bin %i",i));
601  fHistJetPtArea[i] = new TH2F(name1,title1,250,-50,200,100,0,1);
602  fOutput->Add(fHistJetPtArea[i]);
603 
604  snprintf(name, nchar, "fHistAreavsRawPt_%i",i);
605  fHistAreavsRawPt[i] = new TH2F(name,name,100,0,1,200,0,200);
606  fOutput->Add(fHistAreavsRawPt[i]);
607  } // loop over centrality
608 
609  } // QA histo switch
610 
611  if (makeBIAShistos) {
612  fHistJetHaddPhiBias = new TH1F("fHistJetHaddPhiBias","Jet-Hadron #Delta#varphi with bias",128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
613  fHistJetHaddPhiINBias = new TH1F("fHistJetHaddPhiINBias","Jet-Hadron #Delta#varphi IN PLANE with bias", 128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
614  fHistJetHaddPhiOUTBias = new TH1F("fHistJetHaddPhiOUTBias","Jet-Hadron #Delta#varphi OUT PLANE with bias",128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
615  fHistJetHaddPhiMIDBias = new TH1F("fHistJetHaddPhiMIDBias","Jet-Hadron #Delta#varphi MIDDLE of PLANE with bias",128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
616 
617  // add to output list
618  fOutput->Add(fHistJetHaddPhiBias);
619  fOutput->Add(fHistJetHaddPhiINBias);
620  fOutput->Add(fHistJetHaddPhiOUTBias);
621  fOutput->Add(fHistJetHaddPhiMIDBias);
622 
623  for (Int_t i = 0;i<6;++i){
624  name1 = TString(Form("JetPtvsdEPBias_%i",i));
625  title1 = TString(Form("Bias Jet p_{T} vs #DeltaEP cent bin %i",i));
626  fHistJetPtvsdEPBias[i] = new TH2F(name1,title1,250,-50,200,288,-2*TMath::Pi(),2*TMath::Pi());
627  fOutput->Add(fHistJetPtvsdEPBias[i]);
628 
629  name1 = TString(Form("JetEtaPhiPtBias_%i",i));
630  title1 = TString(Form("Jet #eta-#phi p_{T} Bias cent bin %i",i));
631  fHistJetEtaPhiPtBias[i] = new TH3F(name1,title1,250,-50,200,100,-1,1,64,-3.2,3.2);
632  fOutput->Add(fHistJetEtaPhiPtBias[i]);
633 
634  name1 = TString(Form("JetPtAreaBias_%i",i));
635  title1 = TString(Form("Jet p_{T} Area Bias cent bin %i",i));
636  fHistJetPtAreaBias[i] = new TH2F(name1,title1,250,-50,200,100,0,1);
637  fOutput->Add(fHistJetPtAreaBias[i]);
638  } // end of centrality loop
639  } // bias histos
640 
641  if (makeoldJEThadhistos) {
642  for(Int_t icent = 0; icent<6; ++icent){
643  snprintf(name, nchar, "fHistJetPtTT_%i",icent);
644  fHistJetPtTT[icent] = new TH1F(name,name,200,0,200);
645  fOutput->Add(fHistJetPtTT[icent]);
646 
647  snprintf(name, nchar, "fHistJetPt_%i",icent);
648  fHistJetPt[icent] = new TH1F(name,name,200,0,200);
649  fOutput->Add(fHistJetPt[icent]);
650 
651  snprintf(name, nchar, "fHistJetPtBias_%i",icent);
652  fHistJetPtBias[icent] = new TH1F(name,name,200,0,200);
653  fOutput->Add(fHistJetPtBias[icent]);
654 
655  for(Int_t iptjet = 0; iptjet<5; ++iptjet){
656  for(Int_t ieta = 0; ieta<3; ++ieta){
657  snprintf(name, nchar, "fHistJetH_%i_%i_%i",icent,iptjet,ieta);
658  fHistJetH[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
659  fOutput->Add(fHistJetH[icent][iptjet][ieta]);
660 
661  snprintf(name, nchar, "fHistJetHBias_%i_%i_%i",icent,iptjet,ieta);
662  fHistJetHBias[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
663  fOutput->Add(fHistJetHBias[icent][iptjet][ieta]);
664 
665  snprintf(name, nchar, "fHistJetHTT_%i_%i_%i",icent,iptjet,ieta);
666  fHistJetHTT[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
667  fOutput->Add(fHistJetHTT[icent][iptjet][ieta]);
668  } // end of eta loop
669  } // end of pt-jet loop
670  } // end of centrality loop
671  } // make JetHadhisto old
672 
673  if (makeextraCORRhistos) {
674  for(Int_t itrackpt=0; itrackpt<9; itrackpt++){
675  snprintf(name, nchar, "fHistJetHadbindPhi_%i",itrackpt);
676  fHistJetHadbindPhi[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
677  fOutput->Add(fHistJetHadbindPhi[itrackpt]);
678 
679  snprintf(name, nchar, "fHistJetHadbindPhiIN_%i",itrackpt);
680  fHistJetHadbindPhiIN[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
681  fOutput->Add(fHistJetHadbindPhiIN[itrackpt]);
682 
683  snprintf(name, nchar, "fHistJetHadbindPhiMID_%i",itrackpt);
684  fHistJetHadbindPhiMID[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
685  fOutput->Add(fHistJetHadbindPhiMID[itrackpt]);
686 
687  snprintf(name, nchar, "fHistJetHadbindPhiOUT_%i",itrackpt);
688  fHistJetHadbindPhiOUT[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
689  fOutput->Add(fHistJetHadbindPhiOUT[itrackpt]);
690  } // end of trackpt bin loop
691 
692  for (Int_t i = 0;i<6;++i){
693  name1 = TString(Form("JetHaddPhiINcent_%i",i));
694  title1 = TString(Form("Jet Hadron #Delta#varphi Distribution IN PLANE cent bin %i",i));
695  fHistJetHaddPhiINcent[i] = new TH1F(name1,title1,128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
696  fOutput->Add(fHistJetHaddPhiINcent[i]);
697 
698  name1 = TString(Form("JetHaddPhiOUTcent_%i",i));
699  title1 = TString(Form("Jet Hadron #Delta#varphi Distribution OUT PLANE cent bin %i",i));
700  fHistJetHaddPhiOUTcent[i] = new TH1F(name1,title1,128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
701  fOutput->Add(fHistJetHaddPhiOUTcent[i]);
702 
703  name1 = TString(Form("JetHaddPhiMIDcent_%i",i));
704  title1 = TString(Form("Jet Hadron #Delta#varphi Distribution MIDDLE of PLANE cent bin %i",i));
705  fHistJetHaddPhiMIDcent[i] = new TH1F(name1,title1,128,-0.5*TMath::Pi(), 1.5*TMath::Pi());
706  fOutput->Add(fHistJetHaddPhiMIDcent[i]);
707 
708  name1 = TString(Form("JetPtNcon_%i",i));
709  title1 = TString(Form("Jet p_{T} Ncon cent bin %i",i));
710  fHistJetPtNcon[i] = new TH2F(name1,title1,250,-50,200,100,0,2000);
711  fOutput->Add(fHistJetPtNcon[i]);
712 
713  name1 = TString(Form("JetPtNconBias_%i",i));
714  title1 = TString(Form("Jet p_{T} NconBias cent bin %i",i));
715  fHistJetPtNconBias[i] = new TH2F(name1,title1,250,-50,200,100,0,2000);
716  fOutput->Add(fHistJetPtNconBias[i]);
717 
718  name1 = TString(Form("JetPtNconCh_%i",i));
719  title1 = TString(Form("Jet p_{T} NconCh cent bin %i",i));
720  fHistJetPtNconCh[i] = new TH2F(name1,title1,250,-50,200,100,0,2000);
721  fOutput->Add(fHistJetPtNconCh[i]);
722 
723  name1 = TString(Form("JetPtNconBiasCh_%i",i));
724  title1 = TString(Form("Jet p_{T} NconBiasCh cent bin %i",i));
725  fHistJetPtNconBiasCh[i] = new TH2F(name1,title1,250,-50,200,100,0,2000);
726  fOutput->Add(fHistJetPtNconBiasCh[i]);
727 
728  name1 = TString(Form("JetPtNconEm_%i",i));
729  title1 = TString(Form("Jet p_{T} NconEm cent bin %i",i));
730  fHistJetPtNconEm[i] = new TH2F(name1,title1,250,-50,200,100,0,2000);
731  fOutput->Add(fHistJetPtNconEm[i]);
732 
733  name1 = TString(Form("JetPtNconBiasEm_%i",i));
734  title1 = TString(Form("Jet p_{T} NconBiasEm cent bin %i",i));
735  fHistJetPtNconBiasEm[i] = new TH2F(name1,title1,250,-50,200,100,0,2000);
736  fOutput->Add(fHistJetPtNconBiasEm[i]);
737  } // extra Correlations histos switch
738  }
739 
740  // variable binned pt for THnSparse's
741  Double_t xlowjetPT[] = {15, 16, 18, 20, 22, 24, 26, 28, 30, 32, 34, 36, 38, 40, 45, 50, 55, 60, 65, 70, 75, 80, 100, 150, 200, 300};
742  Double_t xlowtrPT[] = {0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0,4.2,4.4,4.6,4.8,5.0,5.5,6.0,6.5,7.0,7.5,8.0,9.0,10.0,12.0,14.0,16.0,18.0,20.0,25.0,30.0,40.0,50.0,75.0};
743 
744  // tracks: 51, jets: 26
745  // number of bins you tell histogram should be (# in array - 1) because the last bin
746  // is the right-most edge of the histogram
747  // i.e. this is for PT and there are 57 numbers (bins) thus we have 56 bins in our histo
748  Int_t nbinsjetPT = sizeof(xlowjetPT)/sizeof(Double_t) - 1;
749  Int_t nbinstrPT = sizeof(xlowtrPT)/sizeof(Double_t) - 1;
750 
751  // set up jet-hadron sparse
752  UInt_t bitcodeMESE = 0; // bit coded, see GetDimParams() below
753  UInt_t bitcodePID = 0; // bit coded, see GetDimParamsPID() below
754  UInt_t bitcodeCorr = 0; // bit coded, see GetDimparamsCorr() below
755  bitcodeMESE = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6; // | 1<<7 | 1<<8 | 1<<9;
756  if(fDoEventMixing) {
757  fhnJH = NewTHnSparseF("fhnJH", bitcodeMESE);
758 
759  if(dovarbinTHnSparse){
760  fhnJH->GetAxis(1)->Set(nbinsjetPT, xlowjetPT);
761  fhnJH->GetAxis(2)->Set(nbinstrPT, xlowtrPT);
762  }
763 
764  fOutput->Add(fhnJH);
765  }
766 
767  bitcodeCorr = 1<<0 | 1<<1 | 1<<2 | 1<<3; // | 1<<4 | 1<<5;
768  fhnCorr = NewTHnSparseFCorr("fhnCorr", bitcodeCorr);
769  if(dovarbinTHnSparse) fhnCorr->GetAxis(1)->Set(nbinsjetPT, xlowjetPT);
770  fOutput->Add(fhnCorr);
771 
772  /*
773  Double_t centralityBins[nCentralityBins+1];
774  for(Int_t ic=0; ic<nCentralityBins+1; ic++){
775  if(ic==nCentralityBins) centralityBins[ic]=500;
776  else centralityBins[ic]=10.0*ic;
777  }
778  */
779 
780  // set up centrality bins for mixed events
781  // for pp we need mult bins for event mixing. Create binning here, to also make a histogram from it
782  Int_t nCentralityBinspp = 8;
783  //Double_t centralityBinspp[nCentralityBinspp+1];
784  Double_t centralityBinspp[9] = {0.0, 4., 9, 15, 25, 35, 55, 100.0, 500.0};
785 
786  // Setup for Pb-Pb collisions
787  Int_t nCentralityBinsPbPb = 100;
788  Double_t mult = 1.0;
789  if(fCentBinSize==1) {
790  nCentralityBinsPbPb = 100;
791  mult = 1.0;
792  } else if(fCentBinSize==2){
793  nCentralityBinsPbPb = 50;
794  mult = 2.0;
795  } else if(fCentBinSize==5){
796  nCentralityBinsPbPb = 20;
797  mult = 5.0;
798  } else if(fCentBinSize==10){
799  nCentralityBinsPbPb = 10;
800  mult = 10.0;
801  }
802 
803  Double_t centralityBinsPbPb[nCentralityBinsPbPb]; // nCentralityBinsPbPb
804  for(Int_t ic=0; ic<nCentralityBinsPbPb; ic++){
805  centralityBinsPbPb[ic]=mult*ic;
806  }
807 
808 /*
809  Int_t nCentralityBinsPbPb = 10; //100;
810  Double_t centralityBinsPbPb[nCentralityBinsPbPb+1];
811  for(Int_t ic=0; ic<nCentralityBinsPbPb; ic++){
812  centralityBinsPbPb[ic]=10.0*ic; //1.0*ic;
813  }
814 */
815 
816  if(fBeam == 0) fHistMult = new TH1F("fHistMult","multiplicity",nCentralityBinspp,centralityBinspp);
817  if(fBeam == 1) fHistMult = new TH1F("fHistMult","multiplicity",nCentralityBinsPbPb,centralityBinsPbPb);
818 // fOutput->Add(fHistMult);
819 
820  // Event Mixing
821  Int_t trackDepth = fMixingTracks;
822  Int_t poolsize = 1000; // Maximum number of events, ignored in the present implemented of AliEventPoolManager
823  Int_t nZvtxBins = 5+1+5;
824  Double_t vertexBins[] = {-10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10};
825  Double_t* zvtxbin = vertexBins;
826  if(fBeam == 0) fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBinspp, centralityBinspp, nZvtxBins, zvtxbin);
827  if(fBeam == 1) fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBinsPbPb, centralityBinsPbPb, nZvtxBins, zvtxbin);
828 
829  // set up event mixing sparse
830  if(fDoEventMixing){
831  bitcodeMESE = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6; // | 1<<7 | 1<<8 | 1<<9;
832  fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", bitcodeMESE);
833 
834  // set up some variable binning of the sparse
835  if(dovarbinTHnSparse){
836  fhnMixedEvents->GetAxis(1)->Set(nbinsjetPT, xlowjetPT);
837  fhnMixedEvents->GetAxis(2)->Set(nbinstrPT, xlowtrPT);
838  }
839 
840  fOutput->Add(fhnMixedEvents);
841  } // end of do-eventmixing
842 
843  // set up PID sparse
844  if(doPID){
845  // ****************************** PID *****************************************************
846  // set up PID handler
847  AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
848  AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
849  if(!inputHandler) {
850  AliFatal("Input handler needed");
851  return;
852  }
853 
854  // PID response object
855  fPIDResponse = inputHandler->GetPIDResponse();
856  if (!fPIDResponse) {
857  AliError("PIDResponse object was not created");
858  return;
859  }
860  // *****************************************************************************************
861 
862  // PID counter
863  fHistPID = new TH1F("fHistPID","PID Counter", 15, 0.5, 25.5);
864  SetfHistPIDcounterLabels(fHistPID);
865  fOutput->Add(fHistPID);
866 
867  if(allpidAXIS) {
868  bitcodePID = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<7 | 1<<8 | 1<<9 |
869  1<<10 | 1<<11 | 1<<12 | 1<<13 | 1<<14 | 1<<15 | 1<<16 | 1<<17 | 1<<18; // | 1<<19 | 1<<20;
870  fhnPID = NewTHnSparseFPID("fhnPID", bitcodePID);
871  } else {
872  bitcodePID = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<7 | 1<<8 | 1<<9 |
873  1<<10 | 1<<11; // | 1<<12 | 1<<13;
874  fhnPID = NewTHnSparseFPID("fhnPID", bitcodePID);
875  }
876 
877  // set some variable binning of sparse
878  if(dovarbinTHnSparse){
879  fhnPID->GetAxis(1)->Set(nbinstrPT, xlowtrPT);
880  fhnPID->GetAxis(8)->Set(nbinsjetPT, xlowjetPT);
881  }
882 
883  fOutput->Add(fhnPID);
884  } // end of do-PID
885 
886  // =========== Switch on Sumw2 for all histos ===========
887  for (Int_t i=0; i<fOutput->GetEntries(); ++i) {
888  TH1 *h1 = dynamic_cast<TH1*>(fOutput->At(i));
889  if (h1){
890  h1->Sumw2();
891  continue;
892  }
893  TH2 *h2 = dynamic_cast<TH2*>(fOutput->At(i));
894  if (h2){
895  h2->Sumw2();
896  continue;
897  }
898  TH3 *h3 = dynamic_cast<TH3*>(fOutput->At(i));
899  if (h3){
900  h3->Sumw2();
901  continue;
902  }
903  THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i));
904  if(hn)hn->Sumw2();
905  }
906 
907  PostData(1, fOutput);
908 }
909 
910 //_________________________________________________________________________
911 void AliAnalysisTaskEmcalJetHadEPpid::ExecOnce()
912 {
913  AliAnalysisTaskEmcalJet::ExecOnce();
914 
915  // Added March14, 2016 - used for new framework, to filled ClonesArray with container contents
916  fTracksFromContainer = new TClonesArray("AliVParticle");
917  fTracksFromContainer->SetName("TracksFromContainer");
918  fTracksFromContainer->SetOwner(kTRUE);
919 
920  if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
921  if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
922  if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
923 }
924 
925 //_____________________________________________________________________________
926 Bool_t AliAnalysisTaskEmcalJetHadEPpid::Notify()
927 {
928  // determine the run number to see if the track and jet cuts should be refreshed for semi-good TPC runs
929  ReadVZEROCalibration2011h();
930  return kTRUE;
931 }
932 
933 //_________________________________________________________________________
934 Bool_t AliAnalysisTaskEmcalJetHadEPpid::Run()
935 { // Main loop called for each event
936  // TEST TEST TEST TEST for OBJECTS!
937 
938  fHistEventQA->Fill(1); // All Events that get entered
939 
940  // check and fill a Event Selection QA histogram for different trigger selections
941  UInt_t trig = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
942 
943  // fill Event Trigger QA
944  FillEventTriggerQA(fHistEventSelectionQA, trig);
945 
946  // see if we are running on PbPb and try and get LocalRho object
947  if(GetBeamType() == 1) {
948  if(!fLocalRho){
949  AliError(Form("Couldn't get fLocalRho object, try to get it from Event based on name\n"));
950  fLocalRho = GetLocalRhoFromEvent(fLocalRhoName);
951  if(!fLocalRho) return kTRUE;
952  }
953  } // check for LocalRho object if PbPb data
954 
955  if(!fTracks){
956  AliError(Form("No fTracks object!!\n"));
957  return kTRUE;
958  }
959  if(!fJets){
960  AliError(Form("No fJets object!!\n"));
961  //return kTRUE;
962  }
963 
964  fHistEventQA->Fill(2); // events after object check
965 
966  // what kind of event do we have: AOD or ESD?
967  Bool_t useAOD;
968  if (dynamic_cast<AliAODEvent*>(InputEvent())) useAOD = kTRUE;
969  else useAOD = kFALSE;
970 
971  // if we have ESD event, set up ESD object
972  if(!useAOD){
973  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
974  if (!fESD) {
975  AliError(Form("ERROR: fESD not available\n"));
976  return kTRUE;
977  }
978  }
979 
980  // if we have AOD event, set up AOD object
981  if(useAOD){
982  fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
983  if(!fAOD) {
984  AliError(Form("ERROR: fAOD not available\n"));
985  return kTRUE;
986  }
987  }
988 
989  fHistEventQA->Fill(3); // events after Aod/esd check
990 
991  // get centrality
992  Int_t centbin = GetCentBin(fCent);
993 
994  // to limit filling unused entries in sparse, only fill for certain centrality ranges
995  // ranges can be different than functional cent bin setter
996  Int_t cbin = -1;
997  if (fCent>=0 && fCent<10) cbin = 1;
998  else if (fCent>=10 && fCent<20) cbin = 2;
999  else if (fCent>=20 && fCent<30) cbin = 3;
1000  else if (fCent>=30 && fCent<50) cbin = 4;
1001  else if (fCent>=50 && fCent<90) cbin = 5;
1002  else cbin = -99;
1003 
1005 /*
1006  if(fCent <= fCentralityClasses->At(0) || fCent >= fCentralityClasses->At(fCentralityClasses->GetSize()-1) || TMath::Abs(fCent-InputEvent()->GetCentrality()->GetCentralityPercentile("TRK")) > 5.) return kFALSE;
1007  // determine centrality class
1008  fInCentralitySelection = -1;
1009  for(Int_t i(0); i < fCentralityClasses->GetSize()-1; i++) {
1010  if(fCent >= fCentralityClasses->At(i) && fCent <= fCentralityClasses->At(1+i)) {
1011  fInCentralitySelection = i;
1012  break;
1013  }
1014  }
1015  if(fInCentralitySelection<0) return kFALSE; // should be null op
1017 */
1018 
1019  // if we are on PbPb data do cut on centrality > 90%, else by default DON'T
1020  if (GetBeamType() == 1) {
1021  // apply cut to event on Centrality > 90%
1022  if(fCent>90) return kTRUE;
1023  }
1024 
1025  // BEAM TYPE enumerator: kNA = -1, kpp = 0, kAA = 1, kpA = 2
1026  // for pp analyses we will just use the first centrality bin
1027  if(GetBeamType() == 0) if (centbin == -1) centbin = 0;
1028  if(GetBeamType() == 1) if (centbin == -1) return kTRUE;
1029 
1030  fHistEventQA->Fill(4); // events after centrality check
1031 
1032  // get vertex information
1033  Double_t fvertex[3]={0,0,0};
1034  InputEvent()->GetPrimaryVertex()->GetXYZ(fvertex);
1035  Double_t zVtx=fvertex[2];
1036 
1037  // get z-vertex bin
1038  //Int_t zVbin = GetzVertexBin(zVtx);
1039 
1040  // apply zVtx cut
1041  if(fabs(zVtx)>10.0) return kTRUE;
1042 
1043  fHistEventQA->Fill(5); // events after zvertex check
1044 
1045  // create pointer to list of input event
1046  TList *list = InputEvent()->GetList();
1047  if(!list) {
1048  AliError(Form("ERROR: list not attached\n"));
1049  return kTRUE;
1050  }
1051 
1052  fHistEventQA->Fill(6); // events after list check
1053 
1054 //test
1055  GetTriggerList();
1056 
1057 // ==================================================================================================================================
1058 
1059  // initialize TClonesArray pointers to jets and tracks
1060  TClonesArray *jets = 0;
1061  TClonesArray *tracks = 0;
1062  TClonesArray *tracksME = 0;
1063  TClonesArray *clusters = 0;
1064 
1065  // TESTING HERE!!!
1066  AliTrackContainer* fTracksCont2 = GetTrackContainer("MyTrackContainer_JetHad");
1067  if(!fTracksCont2) {
1068  AliError(Form("Pointer to tracks: %s == 0", fTracksCont2->GetName()));
1069  } else {
1070  if(doComments) cout<<"#tracks = "<<fTracksCont2->GetNParticles()<<" #accepted tracks = "<<fTracksCont2->GetNAcceptedParticles()<<endl;
1071  }
1072 
1073  // method for filling collection output array of 'saved' tracks - Added March14, 2016 for new framework
1074  Int_t tacc = 0;
1075  fTracksFromContainer->Clear();
1076 
1077  // if we have track container loop over and add to TClonesArray
1078  if(fTracksCont2) {
1079  fTracksCont2->ResetCurrentID();
1080  AliVTrack *track = dynamic_cast<AliVTrack*>(fTracksCont2->GetNextAcceptParticle());
1081  while(track) {
1082  // add container tracks to clones array
1083  (*fTracksFromContainer)[tacc] = track;
1084  ++tacc;
1085 
1086  // apply track cuts
1087  if(TMath::Abs(track->Eta())>fTrkEta) continue;
1088  if (track->Pt()<0.15) continue;
1089 
1090  // calculate single particle tracking efficiency for correlations
1091  //Double_t efficiency = -999;
1092  //efficiency = EffCorrection(track->Eta(), track->Pt(), fDoEffCorr);
1093 
1094  // fill track distributions here (efficiency corrected)
1095  fHistNTrackPhiNEW->Fill(track->Phi());
1096  fHistNTrackPtNEW->Fill(track->Pt());
1097  fHistNTrackEtaNEW->Fill(track->Eta());
1098  fHistNTrackPhiEtaNEW->Fill(track->Phi(), track->Eta());
1099 
1100  // get next accepted track
1101  track = dynamic_cast<AliVTrack*>(fTracksCont2->GetNextAcceptParticle());
1102  } // accepted track loop
1103  }
1104 
1105  // get Tracks object - have switch for old/new framework version of tracks
1106  if(!douseOLDtrackFramework) { // NEW
1107  tracks = fTracksFromContainer; // added March14, 2016
1108  } else { // OLD
1109  tracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracks));
1110  }
1111  if (!tracks) {
1112  AliError(Form("Pointer to tracks: %s == 0", tracks->GetName()));
1113  return kTRUE;
1114  } // verify existence of tracks
1115 
1116  // get ME Tracks object - have switch for old/new framework version of tracks
1117  if(!douseOLDtrackFramework) { // NEW
1118  tracksME = fTracksFromContainer; // added March14, 2016
1119  } else { // OLD
1120  tracksME = dynamic_cast<TClonesArray*>(list->FindObject(fTracks));
1121  }
1122  if (!tracksME) {
1123  AliError(Form("Pointer to ME tracks: %s == 0", tracksME->GetName()));
1124  return kTRUE;
1125  } // verify existence of tracksME
1126 
1127 // ==================================================================================================================================
1128 
1129  // get Jets object
1130  Int_t Njets = 0;
1131  if(fJets) {
1132  jets = dynamic_cast<TClonesArray*>(list->FindObject(fJets));
1133  if(!jets){
1134  AliError(Form("Pointer to jets: %s == 0", fJets->GetName()));
1135  //return kTRUE;
1136  } // verify existence of jets
1137  Njets = jets->GetEntries();
1138  }
1139 
1140  fHistEventQA->Fill(7); // events after track/jet pointer check
1141 
1142  // get number of jets and tracks
1143 // Int_t Njets;
1144 // if(!fJets) { Njets = 0; }
1145 // else Njets = jets->GetEntries();
1146  //const Int_t Njets = jets->GetEntries();
1147  const Int_t Ntracks = tracks->GetEntries();
1150 
1151  fHistEventQA->Fill(8); // events after #track and jets < 1 check
1152 
1153  if (makeQAhistos) fHistMult->Fill(Ntracks); // fill multiplicity distribution
1154 
1155  // get cluster collection (mainly QA at this point)
1156  clusters = dynamic_cast<TClonesArray*>(list->FindObject(fCaloClustersName));
1157  if (!clusters) {
1158  AliError(Form("Pointer to clusters: %s == 0", fCaloClustersName.Data()));
1159  return kTRUE;
1160  } // verify existence of clusters
1161 
1162  // get clusters and loop over
1163  const Int_t Nclusters = clusters->GetEntries();
1164  for (Int_t iclus = 0; iclus < Nclusters; iclus++){
1165  AliVCluster* cluster = static_cast<AliVCluster*>(clusters->At(iclus));
1166  if(!cluster){
1167  AliError(Form("Couldn't get AliVCluster %d\n", iclus));
1168  continue;
1169  }
1170 
1171  // get some info on cluster and fill histo
1172  TLorentzVector nPart;
1173  cluster->GetMomentum(nPart, fvertex);
1174  if(makeQAhistos) fHistClusEtaPhiEnergy->Fill(nPart.Eta(), nPart.Phi(), nPart.E());
1175  if(makeQAhistos) fHistClusEnergy->Fill(nPart.E());
1176  }
1177 
1178  // get event object
1179  fVevent = dynamic_cast<AliVEvent*>(InputEvent());
1180  if (!fVevent) {
1181  AliError(Form("ERROR: AliVEvent (fVevent) not available! \n"));
1182  return kTRUE;
1183  }
1184 
1185  // fill event mixing QA
1186  if(trig & AliVEvent::kEMCEGA) fHistCentZvertGA->Fill(fCent, zVtx);
1187  if(trig & AliVEvent::kEMCEJE) fHistCentZvertJE->Fill(fCent, zVtx);
1188  if(trig & AliVEvent::kMB) fHistCentZvertMB->Fill(fCent, zVtx);
1189  if(trig & AliVEvent::kAny) fHistCentZvertAny->Fill(fCent, zVtx);
1190 
1191  // initialize track parameters
1192  Int_t iTT=-1;
1193  Double_t ptmax=-10;
1194  Int_t NtrackAcc = 0;
1195 
1196  // loop over tracks - to get hardest track (highest pt)
1197  for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++){
1198  AliVTrack *track = dynamic_cast<AliVTrack*>(tracks->At(iTracks));
1199  if(!track){
1200  AliError(Form("Couldn't get AliVTrack %d\n", iTracks));
1201  continue;
1202  }
1203 
1204  // apply track cuts
1205  if(TMath::Abs(track->Eta())>fTrkEta) continue;
1206  if (track->Pt()<0.15) continue;
1207 
1208  //iCount++;
1209  NtrackAcc++;
1210 
1211  if(track->Pt()>ptmax){
1212  ptmax=track->Pt(); // max pT track
1213  iTT=iTracks; // trigger tracks
1214  } // check if Pt>maxpt
1215 
1216  // calculate single particle tracking efficiency for correlations
1217  Double_t efficiency = -999;
1218  efficiency = EffCorrection(track->Eta(), track->Pt(), fDoEffCorr);
1219  //efficiency = 1.0;
1220 
1221  // fill track distributions here
1222  fHistNTrackPt->Fill(track->Pt());
1223  fHistNTrackPhi->Fill(track->Phi());
1224  fHistNTrackEta->Fill(track->Eta());
1225  fHistNTrackPhiEta->Fill(track->Phi(), track->Eta());
1226 
1227  if (makeQAhistos) fHistTrackPhi->Fill(track->Phi());
1228  if (makeQAhistos) fHistTrackPt[centbin]->Fill(track->Pt());
1229  if (makeQAhistos) fHistTrackPtallcent->Fill(track->Pt());
1230  } // end of loop over tracks
1231 
1232  // do event plane analysis for resolution parameter calculation
1233  if(doEventPlaneRes){
1234  // cut on event selection before calculating reaction plane and filling histo's
1235  if ((trig && fTriggerEventType)) {
1236  // cache the leading jet within acceptance
1237  fLeadingJet = GetLeadingJet();
1238 
1239  // set storage vectors for 2nd and 3rd order event plane values for different subevents
1240  Double_t vzero[2][2];
1241  /* for the combined vzero event plane
1242  * [0] psi2 [1] psi3
1243  */
1244  Double_t vzeroComb[2];
1245  // [0] psi2 [1] psi3
1246  Double_t tpc[2];
1247  // evaluate the actual event planes
1248  // grab the actual data
1249  CalculateEventPlaneVZERO(vzero);
1250  CalculateEventPlaneCombinedVZERO(vzeroComb);
1251  CalculateEventPlaneTPC(tpc);
1252  CalculateEventPlaneResolution(vzero, vzeroComb, tpc);
1253  }
1254  }
1255 
1256  // get rho from event and fill relative histo's
1257  fRho = GetRhoFromEvent(fRhoName);
1258  if(!fRho) {
1259  AliError(Form("Couldn't get fRho named %s\n", fRhoName.Data()));
1260  return kFALSE;
1261  }
1262  fRhoVal = fRho->GetVal();
1263 
1264  if (makeQAhistos) {
1265  fHistRhovsdEP[centbin]->Fill(fRhoVal,fEPV0); // Global Rho vs delta Event Plane angle
1266  fHistRhovsCent->Fill(fCent,fRhoVal); // Global Rho vs Centrality
1267  fHistEP0[centbin]->Fill(fEPV0);
1268  fHistEP0A[centbin]->Fill(fEPV0A);
1269  fHistEP0C[centbin]->Fill(fEPV0C);
1270  fHistEPAvsC[centbin]->Fill(fEPV0A,fEPV0C);
1271  }
1272 
1273  // initialize jet parameters
1274  Int_t ijethi = -1;
1275  Int_t passedTTcut = 0;
1276  Int_t NjetAcc = 0;
1277  Double_t highestjetpt = 0.0;
1278  Double_t leadhadronPT = 0.0;
1279  Double_t maxclusterpt = 0.0;
1280  Double_t maxtrackpt = 0.0;
1281 
1282  // loop over jets in an event - to find highest jet pT and apply some cuts
1283  for (Int_t ijet = 0; ijet < Njets; ijet++){
1284  // get our jets
1285  AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
1286  if (!jet) continue;
1287 
1288  // apply jet cuts
1289  if ((jet->Phi()<fPhimin)||(jet->Phi()>fPhimax)) continue;
1290  if ((jet->Eta()<fEtamin)||(jet->Eta()>fEtamax)) continue;
1291  if (makeQAhistos) fHistAreavsRawPt[centbin]->Fill(jet->Pt(),jet->Area());
1292  if(!AcceptMyJet(jet)) continue;
1293 
1294  NjetAcc++; // # of accepted jets
1295 
1296  // if FlavourJetAnalysis, get desired FlavTag and check against Jet
1297  if(doFlavourJetAnalysis) { if(!AcceptFlavourJet(jet, fJetFlavTag)) continue;}
1298 
1299  // fill track distributions here
1300  fHistNJetPt->Fill(jet->Pt());
1301  fHistNJetPhi->Fill(jet->Phi());
1302  fHistNJetEta->Fill(jet->Eta());
1303  fHistNJetPhiEta->Fill(jet->Phi(), jet->Eta());
1304 
1305  // use this to get total # of jets passing cuts in events!!!!!!!!
1306  if (makeQAhistos) fHistJetPhi->Fill(jet->Phi()); // Jet Phi histogram (filled)
1307 
1308  // get highest Pt jet in event
1309  if(highestjetpt<jet->Pt()){
1310  ijethi=ijet;
1311  highestjetpt=jet->Pt();
1312  }
1313  } // end of looping over jets
1314 
1315  // accepted jets
1316  fHistNjetvsCent->Fill(fCent,NjetAcc);
1317  Int_t NJETAcc = 0;
1318  fHistEventQA->Fill(9); // events after track/jet loop to get highest pt
1319 
1320  //cout<<"Event #: "<<event<<endl;
1321 
1322  // loop over jets in event and make appropriate cuts
1323  for (Int_t ijet = 0; ijet < Njets; ++ijet) {
1324  AliEmcalJet *jet = static_cast<AliEmcalJet*>(fJets->At(ijet));
1325  if (!jet) continue;
1326 
1327  // check our jet is in our selected trigger event
1328  if (!(trig && fTriggerEventType)) continue;
1329 
1330  // (should probably be higher..., but makes a cut on jet pT)
1331  if (jet->Pt()<0.1) continue;
1332  // do we accept jet? apply jet cuts
1333  if (!AcceptMyJet(jet)) continue;
1334 
1335  // if FlavourJetAnalysis, get desired FlavTag and check against Jet
1336  if(doFlavourJetAnalysis) { if(!AcceptFlavourJet(jet, fJetFlavTag)) continue;}
1337 
1338  fHistEventQA->Fill(10); // accepted jets
1339 
1340  // check on lead jet
1341  Double_t leadjet=0;
1342  if (ijet==ijethi) leadjet=1;
1343 
1344  // check on leading hadron pt
1345  if (ijet==ijethi) leadhadronPT = GetLeadingHadronPt(jet);
1346  if (ijet==ijethi) maxclusterpt = jet->MaxClusterPt();
1347  if (ijet==ijethi) maxtrackpt = jet->MaxTrackPt();
1348 
1349  // initialize and calculate various parameters: pt, eta, phi, rho, etc...
1350  NJETAcc++; // # accepted jets
1351  Double_t jetphi = jet->Phi(); // phi of jet
1352  Double_t jeteta = jet->Eta(); // ETA of jet
1353  Double_t jetPt = -500;
1354  Double_t jetPtGlobal = -500;
1355  Double_t jetPtLocal = -500; // initialize corr jet pt
1356  if(GetBeamType() == 1) {
1357  fLocalRhoVal = fLocalRho->GetLocalVal(jetphi, fJetRad); //GetJetRadius(0)); // get local rho value
1358  jetPtLocal = jet->Pt()-jet->Area()*fLocalRhoVal; // corrected pT of jet using Rho modulated for V2 and V3
1359  }
1360  jetPt = jet->Pt();
1361  jetPtGlobal = jet->Pt()-jet->Area()*fRhoVal; // corrected pT of jet from rho value
1362  Double_t dEP = -500; // initialize angle between jet and event plane
1363  dEP = RelativeEPJET(jetphi,fEPV0); // angle betweeen jet and event plane
1364 
1365  // make histo's
1366  if(makeQAhistos) fHistJetPtvsTrackPt[centbin]->Fill(jetPt,jet->MaxTrackPt());
1367  if(makeQAhistos) fHistJetPtcorrGlRho[centbin]->Fill(jetPtGlobal);
1368  if(makeQAhistos) fHistJetPtvsdEP[centbin]->Fill(jetPt, dEP);
1369  if(makeQAhistos) fHistJetEtaPhiPt[centbin]->Fill(jetPt,jet->Eta(),jet->Phi());
1370  if(makeQAhistos) fHistJetPtArea[centbin]->Fill(jetPt,jet->Area());
1371  if(makeQAhistos) fHistJetEtaPhi->Fill(jet->Eta(),jet->Phi()); // fill jet eta-phi distribution histo
1372  if(makeextraCORRhistos) fHistJetPtNcon[centbin]->Fill(jetPt,jet->GetNumberOfConstituents());
1373  if(makeextraCORRhistos) fHistJetPtNconCh[centbin]->Fill(jetPt,jet->GetNumberOfTracks());
1374  if(makeextraCORRhistos) fHistJetPtNconEm[centbin]->Fill(jetPt,jet->GetNumberOfClusters());
1375  if(makeoldJEThadhistos) fHistJetPt[centbin]->Fill(jet->Pt()); // fill #jets vs pT histo
1376  //fHistDeltaPtvsArea->Fill(jetPt,jet->Area());
1377 
1378  // make histo's with BIAS applied
1379  if (jet->MaxTrackPt()>fTrkBias){
1380  if(makeBIAShistos) fHistJetPtvsdEPBias[centbin]->Fill(jetPt, dEP);
1381  if(makeBIAShistos) fHistJetEtaPhiPtBias[centbin]->Fill(jetPt,jet->Eta(),jet->Phi());
1382  if(makeextraCORRhistos) fHistJetPtAreaBias[centbin]->Fill(jetPt,jet->Area());
1383  if(makeextraCORRhistos) fHistJetPtNconBias[centbin]->Fill(jetPt,jet->GetNumberOfConstituents());
1384  if(makeextraCORRhistos) fHistJetPtNconBiasCh[centbin]->Fill(jetPt,jet->GetNumberOfTracks());
1385  if(makeextraCORRhistos) fHistJetPtNconBiasEm[centbin]->Fill(jetPt,jet->GetNumberOfClusters());
1386  }
1387 
1388  //if(leadjet && centbin==0){
1389  // if(makeextraCORRhistos) fHistJetPt[centbin+1]->Fill(jet->Pt());
1390  //}
1391  if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
1392  if (makeoldJEThadhistos){
1393  fHistJetPtBias[centbin]->Fill(jet->Pt());
1394  //if(leadjet && centbin==0) fHistJetPtBias[centbin+1]->Fill(jet->Pt());
1395  }
1396  } // end of MaxTrackPt>ftrkBias or maxclusterPt>fclusBias
1397 
1398  // do we have trigger tracks
1399  if(iTT>0){
1400  AliVTrack* TT = static_cast<AliVTrack*>(tracks->At(iTT));
1401  if(TMath::Abs(jet->Phi()-TT->Phi()-TMath::Pi())<0.6) passedTTcut=1;
1402  else passedTTcut=0;
1403  } // end of check on iTT > 0
1404 
1405  if(passedTTcut) {
1406  if (makeoldJEThadhistos) fHistJetPtTT[centbin]->Fill(jet->Pt());
1407  }
1408 
1409  // cut on HIGHEST jet pt in event (15 GeV default)
1410  //if (highestjetpt>fJetPtcut) {
1411  if (jet->Pt() > fJetPtcut) {
1412  fHistEventQA->Fill(11); // jets meeting pt threshold
1413 
1414  // does our max track or cluster pass the bias?
1415  if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
1416  // set up and fill Jet-Hadron trigger jets THnSparse
1417  if(fcorrJetPt) {
1418  Double_t CorrEntries[4] = {fCent, jetPtLocal, dEP, zVtx};
1419  if(fReduceStatsCent > 0) {
1420  if(cbin == fReduceStatsCent) fhnCorr->Fill(CorrEntries); // fill Sparse Histo with trigger Jets entries
1421  } else fhnCorr->Fill(CorrEntries); // fill Sparse Histo with trigger Jets entries
1422  } else { // don't correct jet pt
1423  Double_t CorrEntries[4] = {fCent, jet->Pt(), dEP, zVtx};
1424  if(fReduceStatsCent > 0) {
1425  if(cbin == fReduceStatsCent) fhnCorr->Fill(CorrEntries); // fill Sparse Histo with trigger Jets entries
1426  } else fhnCorr->Fill(CorrEntries); // fill Sparse Histo with trigger Jets entries
1427  }
1428 
1429  if(GetBeamType() == 1) fHistLocalRhoJetpt->Fill(jetPtLocal);
1430  if(makeQAhistos) fHistJetNConstit->Fill(jet->GetNumberOfConstituents());
1431  if(makeQAhistos) fHistJetNTrackConstit->Fill(jet->GetNumberOfTracks());
1432  if(makeQAhistos) fHistJetNClusterConstit->Fill(jet->GetNumberOfClusters());
1433 
1434  // get clusters and loop over
1435  for (Int_t iclus = 0; iclus < Nclusters; iclus++){
1436  AliVCluster* clusterofJet = static_cast<AliVCluster*>(clusters->At(iclus));
1437  if(!clusterofJet){
1438  AliError(Form("Couldn't get AliVCluster %d\n", iclus));
1439  continue;
1440  }
1441 
1442  // get some info on cluster and fill histo
1443  TLorentzVector nPart;
1444  clusterofJet->GetMomentum(nPart, fvertex);
1445  //if(makeQAhistos) fHistClusEtaPhiEnergy->Fill(nPart.Eta(), nPart.Phi(), nPart.E());
1446  if(makeQAhistos) fHistClusofJetEnergy->Fill(nPart.E());
1447 
1448  } // loop over clusters for each jet
1449  } // check on max track and cluster pt/Et
1450 
1451  // loop over all track for an event containing a jet with a pT>fJetPtCut (15)GeV
1452  for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
1453  AliVTrack *track = dynamic_cast<AliVTrack*>(tracks->At(iTracks));
1454  if (!track) {
1455  AliError(Form("Couldn't get AliVtrack %d\n", iTracks));
1456  continue;
1457  }
1458 
1459  // apply track cuts
1460  if(TMath::Abs(track->Eta())>fTrkEta) continue;
1461  if (track->Pt()<0.15) continue;
1462 
1463  fHistEventQA->Fill(12); // accepted tracks in events from trigger jets
1464 
1465  // calculate and get some track parameters
1466  Double_t trCharge = -99;
1467  trCharge = track->Charge();
1468  Double_t tracketa=track->Eta(); // eta of track
1469  Double_t deta=tracketa-jeteta; // dETA between track and jet
1470  //Double_t dR=sqrt(deta*deta+dphijh*dphijh); // difference of R between jet and hadron track
1471 
1472  // calculate single particle tracking efficiency for correlations
1473  Double_t trefficiency = -999;
1474  trefficiency = EffCorrection(track->Eta(), track->Pt(), fDoEffCorr);
1475 
1476  Int_t ieta = -1; // initialize deta bin
1477  Int_t iptjet = -1; // initialize jet pT bin
1478  if (makeoldJEThadhistos) {
1479  ieta=GetEtaBin(deta); // bin of eta
1480  if(ieta<0) continue; // double check we don't have a negative array index
1481  iptjet=GetpTjetBin(jet->Pt()); // bin of jet pt
1482  if(iptjet<0) continue; // double check we don't have a negative array index
1483  }
1484 
1485  // dPHI between jet and hadron
1486  Double_t dphijh = RelativePhi(jet->Phi(), track->Phi()); // angle between jet and hadron
1487 
1488  // fill some jet-hadron histo's
1489  if (makeoldJEThadhistos) fHistJetH[centbin][iptjet][ieta]->Fill(dphijh,track->Pt()); // fill jet-hadron dPHI--track pT distribution
1490  if(makeQAhistos) fHistJetHEtaPhi->Fill(deta,dphijh); // fill jet-hadron eta--phi distribution
1491  fHistJetHaddPHI->Fill(dphijh);
1492  if(passedTTcut){
1493  if (makeoldJEThadhistos) fHistJetHTT[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
1494  }
1495 
1496  // does our max track or cluster pass the bias?
1497  if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
1498  // set up and fill Jet-Hadron THnSparse
1499  if(fcorrJetPt){
1500  Double_t triggerEntries[9] = {fCent, jetPtLocal, track->Pt(), deta, dphijh, dEP, zVtx, trCharge, leadjet};
1501  if(fDoEventMixing) {
1502  if(fReduceStatsCent > 0) {
1503  if(cbin == fReduceStatsCent) fhnJH->Fill(triggerEntries, 1.0/trefficiency); // fill Sparse Histo with trigger entries
1504  } else fhnJH->Fill(triggerEntries, 1.0/trefficiency); // fill Sparse Histo with trigger entries
1505  }
1506  } else { // don't correct jet pt
1507  Double_t triggerEntries[9] = {fCent, jet->Pt(), track->Pt(), deta, dphijh, dEP, zVtx, trCharge, leadjet};
1508  if(fDoEventMixing) {
1509  if(fReduceStatsCent > 0) {
1510  if(cbin == fReduceStatsCent) fhnJH->Fill(triggerEntries, 1.0/trefficiency); // fill Sparse Histo with trigger entries
1511  } else fhnJH->Fill(triggerEntries, 1.0/trefficiency); // fill Sparse Histo with trigger entries
1512  }
1513  }
1514 
1515  // fill histo's
1516  if(makeQAhistos) fHistSEphieta->Fill(dphijh, deta); // single event distribution
1517  if(makeoldJEThadhistos) fHistJetHBias[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
1518 
1519  if (makeBIAShistos) {
1520  fHistJetHaddPhiBias->Fill(dphijh);
1521 
1522  // in plane and out of plane histo's
1523  if( dEP>0 && dEP<=(TMath::Pi()/6) ){
1524  // we are IN plane
1525  fHistJetHaddPhiINBias->Fill(dphijh);
1526  }else if( dEP>(TMath::Pi()/3) && dEP<=(TMath::Pi()/2) ){
1527  // we are OUT of PLANE
1528  fHistJetHaddPhiOUTBias->Fill(dphijh);
1529  }else if( dEP>(TMath::Pi()/6) && dEP<=(TMath::Pi()/3) ){
1530  // we are in middle of plane
1531  fHistJetHaddPhiMIDBias->Fill(dphijh);
1532  }
1533  } // BIAS histos switch
1534 
1535  } // end of check maxtrackpt>ftrackbias or maxclusterpt>fclustbias
1536 
1537  // **************************************************************************************************************
1538  // *********************************** PID **********************************************************************
1539  // **************************************************************************************************************
1540  if(doPIDtrackBIAS){
1541  //if(ptmax < fTrkBias) continue; // force PID to happen when max track pt > 5.0 GeV
1542  if(leadhadronPT < fTrkBias) continue; // force PID to happen when lead hadron pt > 5.0 GeV
1543  }
1544 
1545  // some variables for PID
1546  Double_t pt = -999, dEdx = -999, ITSsig = -999, TOFsig = -999, charge = -999;
1547 
1548  // nSigma of particles in TPC, TOF, and ITS
1549  Double_t nSigmaPion_TPC, nSigmaProton_TPC, nSigmaKaon_TPC;
1550  Double_t nSigmaPion_TOF, nSigmaProton_TOF, nSigmaKaon_TOF;
1551  Double_t nSigmaPion_ITS, nSigmaProton_ITS, nSigmaKaon_ITS;
1552 
1553  if(fPIDResponse) fHistTPCdEdX->Fill(track->Pt(), track->GetTPCsignal());
1554 
1555  if(doPID && (fPIDResponse) && ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)) ){
1556  // get parameters of track
1557  charge = track->Charge(); // charge of track
1558  pt = track->Pt(); // pT of track
1559 
1560  //if (!fPIDResponse) return kTRUE; // just return, maybe put at beginning
1561 
1562  fHistEventQA->Fill(13); // check for AliVEvent and fPIDresponse objects
1563 
1564  // get detector signals
1565  dEdx = track->GetTPCsignal();
1566  ITSsig = track->GetITSsignal();
1567  TOFsig = track->GetTOFsignal();
1568 
1569  // TPC nSigma's
1570  nSigmaPion_TPC = fPIDResponse->NumberOfSigmasTPC(track,AliPID::kPion);
1571  nSigmaKaon_TPC = fPIDResponse->NumberOfSigmasTPC(track,AliPID::kKaon);
1572  nSigmaProton_TPC = fPIDResponse->NumberOfSigmasTPC(track,AliPID::kProton);
1573 
1574  // TOF nSigma's
1575  nSigmaPion_TOF = fPIDResponse->NumberOfSigmasTOF(track,AliPID::kPion);
1576  nSigmaKaon_TOF = fPIDResponse->NumberOfSigmasTOF(track,AliPID::kKaon);
1577  nSigmaProton_TOF = fPIDResponse->NumberOfSigmasTOF(track,AliPID::kProton);
1578 
1579  // ITS nSigma's
1580  nSigmaPion_ITS = fPIDResponse->NumberOfSigmasITS(track,AliPID::kPion);
1581  nSigmaKaon_ITS = fPIDResponse->NumberOfSigmasITS(track,AliPID::kKaon);
1582  nSigmaProton_ITS = fPIDResponse->NumberOfSigmasITS(track,AliPID::kProton);
1583 
1584  // fill detector signal histograms
1585  //if (makeQAhistos) fHistTPCdEdX->Fill(pt, dEdx); // temp out
1586  if (makeQAhistos) fHistITSsignal->Fill(pt, ITSsig);
1587  //if (makeQAhistos) fHistTOFsignal->Fill(pt, TOFsig);
1588 
1589  // Tests to PID pions, kaons, and protons, (default is undentified tracks)
1590  //Double_t nPIDtpc = 0, nPIDits = 0, nPIDtof = 0;
1591  Double_t nPID = -99;
1592 
1593  // check track has pT < 0.900 GeV - use TPC pid
1594  if (pt<0.900 && dEdx>0) {
1595  //nPIDtpc = 4;
1596  nPID = 1;
1597 
1598  // PION check - TPC
1599  if (TMath::Abs(nSigmaPion_TPC)<2 && TMath::Abs(nSigmaKaon_TPC)>2 && TMath::Abs(nSigmaProton_TPC)>2 ){
1600  isPItpc = kTRUE;
1601  //nPIDtpc = 1;
1602  nPID=2;
1603  }else isPItpc = kFALSE;
1604 
1605  // KAON check - TPC
1606  if (TMath::Abs(nSigmaKaon_TPC)<2 && TMath::Abs(nSigmaPion_TPC)>3 && TMath::Abs(nSigmaProton_TPC)>2 ){
1607  isKtpc = kTRUE;
1608  //nPIDtpc = 2;
1609  nPID=3;
1610  }else isKtpc = kFALSE;
1611 
1612  // PROTON check - TPC
1613  if (TMath::Abs(nSigmaProton_TPC)<2 && TMath::Abs(nSigmaPion_TPC)>3 && TMath::Abs(nSigmaKaon_TPC)>2 ){
1614  isPtpc = kTRUE;
1615  //nPIDtpc = 3;
1616  nPID=4;
1617  }else isPtpc = kFALSE;
1618  } // cut on track pT for TPC
1619 
1620  // check track has pT < 0.500 GeV - use ITS pid
1621  if (pt<0.500 && ITSsig>0) {
1622  //nPIDits = 4;
1623  nPID = 5;
1624 
1625  // PION check - ITS
1626  if (TMath::Abs(nSigmaPion_ITS)<2 && TMath::Abs(nSigmaKaon_ITS)>2 && TMath::Abs(nSigmaProton_ITS)>2 ){
1627  isPIits = kTRUE;
1628  //nPIDits = 1;
1629  nPID=6;
1630  }else isPIits = kFALSE;
1631 
1632  // KAON check - ITS
1633  if (TMath::Abs(nSigmaKaon_ITS)<2 && TMath::Abs(nSigmaPion_ITS)>3 && TMath::Abs(nSigmaProton_ITS)>2 ){
1634  isKits = kTRUE;
1635  //nPIDits = 2;
1636  nPID=7;
1637  }else isKits = kFALSE;
1638 
1639  // PROTON check - ITS
1640  if (TMath::Abs(nSigmaProton_ITS)<2 && TMath::Abs(nSigmaPion_ITS)>3 && TMath::Abs(nSigmaKaon_ITS)>2 ){
1641  isPits = kTRUE;
1642  //nPIDits = 3;
1643  nPID=8;
1644  }else isPits = kFALSE;
1645  } // cut on track pT for ITS
1646 
1647  // check track has 0.900 GeV < pT < 2.500 GeV - use TOF pid
1648  if (pt>0.900 && pt<2.500 && TOFsig>0) {
1649  //nPIDtof = 4;
1650  nPID = 9;
1651 
1652  // PION check - TOF
1653  if (TMath::Abs(nSigmaPion_TOF)<2 && TMath::Abs(nSigmaKaon_TOF)>2 && TMath::Abs(nSigmaProton_TOF)>2 ){
1654  isPItof = kTRUE;
1655  //nPIDtof = 1;
1656  nPID=10;
1657  }else isPItof = kFALSE;
1658 
1659  // KAON check - TOF
1660  if (TMath::Abs(nSigmaKaon_TOF)<2 && TMath::Abs(nSigmaPion_TOF)>3 && TMath::Abs(nSigmaProton_TOF)>2 ){
1661  isKtof = kTRUE;
1662  //nPIDtof = 2;
1663  nPID=11;
1664  }else isKtof = kFALSE;
1665 
1666  // PROTON check - TOF
1667  if (TMath::Abs(nSigmaProton_TOF)<2 && TMath::Abs(nSigmaPion_TOF)>3 && TMath::Abs(nSigmaKaon_TOF)>2 ){
1668  isPtof = kTRUE;
1669  //nPIDtof = 3;
1670  nPID=12;
1671  }else isPtof = kFALSE;
1672  } // cut on track pT for TOF
1673 
1674  // if not tagged - unidentifed particle
1675  if (nPID == -99) nPID = 14;
1676 
1677  // fill PID tagged histo
1678  fHistPID->Fill(nPID);
1679 
1680  // TOF PID cuts: (needs testing)
1681  // Pion ID
1682  if((pt > 0.4) && (pt <=1.2) && (nSigmaPion_TOF >= -5.0) && (nSigmaPion_TOF <= 5.0)) nPID = 21;
1683  if((pt > 1.2) && (pt <=1.6) && (nSigmaPion_TOF >= -4.0) && (nSigmaPion_TOF <= 4.0)) nPID = 21;
1684  if((pt > 1.6) && (pt <=2.0) && (nSigmaPion_TOF >= -4.0) && (nSigmaPion_TOF <= 2.0)) nPID = 21;
1685  if((pt > 2.0) && (pt <=3.6) && (nSigmaPion_TOF >= -4.0) && (nSigmaPion_TOF <= 0.0)) nPID = 21;
1686  // Proton ID
1687  if((pt > 0.4) && (pt <=0.9) && (nSigmaProton_TOF >= -3.0) && (nSigmaPion_TOF <= 3.0)) nPID = 22;
1688  if((pt > 0.9) && (pt <=1.4) && (nSigmaProton_TOF >= -4.0) && (nSigmaPion_TOF <= 4.0)) nPID = 22;
1689  if((pt > 1.4) && (pt <=2.2) && (nSigmaProton_TOF >= -5.0) && (nSigmaPion_TOF <= 5.0)) nPID = 22;
1690  if((pt > 2.2) && (pt <=2.4) && (nSigmaProton_TOF >= -4.0) && (nSigmaPion_TOF <= 5.0)) nPID = 22;
1691  if((pt > 2.4) && (pt <=3.0) && (nSigmaProton_TOF >= -2.0) && (nSigmaPion_TOF <= 5.0)) nPID = 22;
1692  if((pt > 3.0) && (pt <=3.6) && (nSigmaProton_TOF >= 0.0) && (nSigmaPion_TOF <= 5.0)) nPID = 22;
1693  // Kaon ID
1694  if((pt > 0.4) && (pt <=0.8) && (nSigmaKaon_TOF >= -5.0) && (nSigmaPion_TOF <= 5.0)) nPID = 23;
1695  if((pt > 0.8) && (pt <=1.2) && (nSigmaKaon_TOF >= -4.0) && (nSigmaPion_TOF <= 4.0)) nPID = 23;
1696  if((pt > 1.2) && (pt <=1.7) && (nSigmaKaon_TOF >= -3.0) && (nSigmaPion_TOF <= 4.0)) nPID = 23;
1697  if((pt > 1.7) && (pt <=2.3) && (nSigmaKaon_TOF >= 0.0) && (nSigmaPion_TOF <= 5.0)) nPID = 23;
1698  if((pt > 2.3) && (pt <=2.6) && (nSigmaKaon_TOF >= 0.0) && (nSigmaPion_TOF <= 3.0)) nPID = 23;
1699  if((pt > 2.6) && (pt <=3.6) && (nSigmaKaon_TOF >= 0.0) && (nSigmaPion_TOF <= 2.0)) nPID = 23;
1700 
1701 
1702  // PID sparse getting filled
1703  if(fcorrJetPt){ // subtract background density
1704  if (allpidAXIS) { // FILL ALL axis
1705  Double_t pid_EntriesALL[19] = {fCent,pt,charge,deta,dphijh,leadjet,zVtx,dEP,jetPtLocal,
1706  nSigmaPion_TPC, nSigmaPion_TOF, // pion nSig values in TPC/TOF
1707  nPID, //nPIDtpc, nPIDits, nPIDtof, // PID label for each detector
1708  nSigmaProton_TPC, nSigmaKaon_TPC, // nSig in TPC
1709  nSigmaPion_ITS, nSigmaProton_ITS, nSigmaKaon_ITS, // nSig in ITS
1710  nSigmaProton_TOF, nSigmaKaon_TOF, // nSig in TOF
1711  }; //array for PID sparse
1712  if(fReduceStatsCent > 0) {
1713  if(cbin == fReduceStatsCent) fhnPID->Fill(pid_EntriesALL, 1.0/trefficiency);
1714  } else fhnPID->Fill(pid_EntriesALL, 1.0/trefficiency);
1715  } else {
1716  // PID sparse getting filled
1717  Double_t pid_Entries[12] = {fCent,pt,charge,deta,dphijh,leadjet,zVtx,dEP,jetPtLocal,
1718  nSigmaPion_TPC, nSigmaPion_TOF, // pion nSig values in TPC/TOF
1719  nPID //nPIDtpc, nPIDits, nPIDtof // PID label for each detector
1720  }; //array for PID sparse
1721  if(fReduceStatsCent > 0) {
1722  if(cbin == fReduceStatsCent) fhnPID->Fill(pid_Entries, 1.0/trefficiency);
1723  } else fhnPID->Fill(pid_Entries, 1.0/trefficiency); // fill Sparse histo of PID tracks
1724  } // minimal pid sparse filling
1725  } else { // don't correct jet pt by background density
1726  if (allpidAXIS) { // FILL ALL axis
1727  Double_t pid_EntriesALL[19] = {fCent,pt,charge,deta,dphijh,leadjet,zVtx,dEP,jetPt,
1728  nSigmaPion_TPC, nSigmaPion_TOF, // pion nSig values in TPC/TOF
1729  nPID, //nPIDtpc, nPIDits, nPIDtof, // PID label for each detector
1730  nSigmaProton_TPC, nSigmaKaon_TPC, // nSig in TPC
1731  nSigmaPion_ITS, nSigmaProton_ITS, nSigmaKaon_ITS, // nSig in ITS
1732  nSigmaProton_TOF, nSigmaKaon_TOF, // nSig in TOF
1733  }; //array for PID sparse
1734  if(fReduceStatsCent > 0) {
1735  if(cbin == fReduceStatsCent) fhnPID->Fill(pid_EntriesALL, 1.0/trefficiency);
1736  } else fhnPID->Fill(pid_EntriesALL, 1.0/trefficiency);
1737  } else {
1738  // PID sparse getting filled
1739  Double_t pid_Entries[12] = {fCent,pt,charge,deta,dphijh,leadjet,zVtx,dEP,jetPt,
1740  nSigmaPion_TPC, nSigmaPion_TOF, // pion nSig values in TPC/TOF
1741  nPID //nPIDtpc, nPIDits, nPIDtof // PID label for each detector
1742  }; //array for PID sparse
1743  if(fReduceStatsCent > 0) {
1744  if(cbin == fReduceStatsCent) fhnPID->Fill(pid_Entries, 1.0/trefficiency);
1745  } else fhnPID->Fill(pid_Entries, 1.0/trefficiency);
1746  } // minimal pid sparse filling
1747  } // don't correct jet pet
1748  } // end of doPID check
1749 
1750  // get track pt bin
1751  Int_t itrackpt = -500; // initialize track pT bin
1752  itrackpt = GetpTtrackBin(track->Pt());
1753 
1754  // all tracks: jet hadron correlations in hadron pt bins
1755  if(makeextraCORRhistos) fHistJetHadbindPhi[itrackpt]->Fill(dphijh);
1756 
1757  // in plane and out of plane jet-hadron histo's
1758  if( dEP>0 && dEP<=(TMath::Pi()/6) ){
1759  // we are IN plane
1760  if(makeextraCORRhistos) fHistJetHaddPhiINcent[centbin]->Fill(dphijh);
1761  fHistJetHaddPhiIN->Fill(dphijh);
1762  if(makeextraCORRhistos) fHistJetHadbindPhiIN[itrackpt]->Fill(dphijh);
1763  }else if( dEP>(TMath::Pi()/3) && dEP<=(TMath::Pi()/2) ){
1764  // we are OUT of PLANE
1765  if(makeextraCORRhistos) fHistJetHaddPhiOUTcent[centbin]->Fill(dphijh);
1766  fHistJetHaddPhiOUT->Fill(dphijh);
1767  if(makeextraCORRhistos) fHistJetHadbindPhiOUT[itrackpt]->Fill(dphijh);
1768  }else if( dEP>(TMath::Pi()/6) && dEP<=(TMath::Pi()/3) ){
1769  // we are in the middle of plane
1770  if(makeextraCORRhistos) fHistJetHaddPhiMIDcent[centbin]->Fill(dphijh);
1771  fHistJetHaddPhiMID->Fill(dphijh);
1772  if(makeextraCORRhistos) fHistJetHadbindPhiMID[itrackpt]->Fill(dphijh);
1773  }
1774  } // loop over tracks found in event with highest JET pT > 10.0 GeV (change)
1775  } // jet pt cut
1776  } // jet loop
1777 
1778  fHistEventQA->Fill(14); // events right before event mixing
1779 
1780 // ***************************************************************************************************************
1781 // ******************************** Event MIXING *****************************************************************
1782 // ***************************************************************************************************************
1783 
1784  // initialize object array of cloned picotracks
1785  TObjArray* tracksClone = 0x0;
1786 
1787  // PbPb collisions - create cloned picotracks
1788  //if(GetBeamType() == 1) tracksClone = CloneAndReduceTrackList(tracks); // TEST
1789 
1790  //Prepare to do event mixing
1791  if(fDoEventMixing>0){
1792  // event mixing
1793 
1794  // 1. First get an event pool corresponding in mult (cent) and
1795  // zvertex to the current event. Once initialized, the pool
1796  // should contain nMix (reduced) events. This routine does not
1797  // pre-scan the chain. The first several events of every chain
1798  // will be skipped until the needed pools are filled to the
1799  // specified depth. If the pool categories are not too rare, this
1800  // should not be a problem. If they are rare, you could lose
1801  // statistics.
1802 
1803  // 2. Collect the whole pool's content of tracks into one TObjArray
1804  // (bgTracks), which is effectively a single background super-event.
1805 
1806  // 3. The reduced and bgTracks arrays must both be passed into
1807  // FillCorrelations(). Also nMix should be passed in, so a weight
1808  // of 1./nMix can be applied.
1809 
1810  // mix jets from triggered events with tracks from MB events
1811  // get the trigger bit, need to change trigger bits between different runs
1812  UInt_t trigger = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
1813  // if event was not selected (triggered) for any reseason (should never happen) then return
1814  if (trigger==0) return kTRUE;
1815 
1816  // initialize event pools
1817  AliEventPool* pool = 0x0;
1818  AliEventPool* poolpp = 0x0;
1819  Double_t Ntrks = -999;
1820 
1821  // pp collisions - get event pool
1822  if(GetBeamType() == 0) {
1823  Ntrks=(Double_t)Ntracks*1.0;
1824  //cout<<"Test.. Ntrks: "<<fPoolMgr->GetEventPool(Ntrks);
1825  poolpp = fPoolMgr->GetEventPool(Ntrks, zVtx); // for pp
1826  }
1827 
1828  // PbPb collisions - get event pool
1829  if(GetBeamType() == 1) pool = fPoolMgr->GetEventPool(fCent, zVtx); // for PbPb? fcent
1830 
1831  // if we don't have a pool, return
1832  if (!pool && !poolpp){
1833  if(GetBeamType() == 1) AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCent, zVtx));
1834  if(GetBeamType() == 0) AliFatal(Form("No pool found for multiplicity = %f, zVtx = %f", Ntrks, zVtx));
1835  return kTRUE;
1836  }
1837 
1838  fHistEventQA->Fill(15); // mixed events cases that have pool
1839 
1840  // initialize background tracks array
1841  TObjArray* bgTracks;
1842 
1843  // next line might not apply for PbPb collisions
1844  // use only jets from EMCal-triggered events (for lhc11a use AliVEvent::kEMC1)
1845  //check for a trigger jet
1846  // fmixingtrack/10 ??
1847  if(GetBeamType() == 1) if(trigger && fTriggerEventType) { //kEMCEJE)) {
1848  if (pool->IsReady() || pool->NTracksInPool() > fNMIXtracks || pool->GetCurrentNEvents() >= fNMIXevents) {
1849 
1850  // loop over jets (passing cuts?)
1851  for (Int_t ijet = 0; ijet < Njets; ijet++) {
1852  Double_t leadjet=0;
1853  if (ijet==ijethi) leadjet=1;
1854 
1855  // get jet object
1856  AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
1857  if (!jet) continue;
1858 
1859  // (should probably be higher..., but makes a cut on jet pT)
1860  if (jet->Pt()<0.1) continue;
1861  if (!AcceptMyJet(jet)) continue;
1862 
1863  Double_t jetPtLocalmix = -500; // initialize corr jet pt
1864  if(GetBeamType() == 1) {
1865  fLocalRhoVal = fLocalRho->GetLocalVal(jet->Phi(), fJetRad); //GetJetRadius(0)); // get local rho value
1866  jetPtLocalmix = jet->Pt()-jet->Area()*fLocalRhoVal; // corrected pT of jet using Rho modulated for V2 and V3
1867  }
1868 
1869  fHistEventQA->Fill(16); // event mixing jets
1870 
1871  // set cut to do event mixing only if we have a jet meeting our pt threshold (bias applied below)
1872  if (jet->Pt()<fJetPtcut) continue;
1873 
1874  // get number of current events in pool
1875  Int_t nMix = pool->GetCurrentNEvents(); // how many particles in pool to mix
1876 
1877  // Fill for biased jet triggers only
1878  if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)) { // && jet->Pt() > fJetPtcut) {
1879  // Fill mixed-event histos here
1880  for (Int_t jMix=0; jMix<nMix; jMix++) {
1881  fHistEventQA->Fill(17); // event mixing nMix
1882 
1883  // get jMix'th event
1884  bgTracks = pool->GetEvent(jMix);
1885  const Int_t Nbgtrks = bgTracks->GetEntries();
1886  for(Int_t ibg=0; ibg<Nbgtrks; ibg++) {
1887  AliPicoTrack *part = static_cast<AliPicoTrack*>(bgTracks->At(ibg));
1888  if(!part) continue;
1889  if(TMath::Abs(part->Eta())>0.9) continue;
1890  if(part->Pt()<0.15) continue;
1891 
1892  Double_t DEta = part->Eta()-jet->Eta(); // difference in eta
1893  Double_t DPhi = RelativePhi(jet->Phi(),part->Phi()); // difference in phi
1894  Double_t dEP = RelativeEPJET(jet->Phi(),fEPV0); // difference between jet and EP
1895  Double_t mixcharge = part->Charge();
1896  //Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta); // difference in R
1897 
1898  // calculate single particle tracking efficiency of mixed events for correlations
1899  Double_t mixefficiency = -999;
1900  mixefficiency = EffCorrection(part->Eta(), part->Pt(), fDoEffCorr);
1901 
1902  // create / fill mixed event sparse
1903  if(fcorrJetPt) {
1904  Double_t triggerEntries[9] = {fCent,jetPtLocalmix,part->Pt(),DEta,DPhi,dEP,zVtx, mixcharge, leadjet}; //array for ME sparse
1905  if(fReduceStatsCent > 0) {
1906  if(cbin == fReduceStatsCent) fhnMixedEvents->Fill(triggerEntries,1./(nMix*mixefficiency)); // fill Sparse histo of mixed events
1907  } else fhnMixedEvents->Fill(triggerEntries,1./(nMix*mixefficiency)); // fill Sparse histo of mixed events
1908  } else { // don't correct jet pt
1909  Double_t triggerEntries[9] = {fCent,jet->Pt(),part->Pt(),DEta,DPhi,dEP,zVtx, mixcharge, leadjet}; //array for ME sparse
1910  if(fReduceStatsCent > 0) {
1911  if(cbin == fReduceStatsCent) fhnMixedEvents->Fill(triggerEntries,1./(nMix*mixefficiency)); // fill Sparse histo of mixed events
1912  } else fhnMixedEvents->Fill(triggerEntries,1./(nMix*mixefficiency)); // fill Sparse histo of mixed events
1913  }
1914 
1915  fHistEventQA->Fill(18); // event mixing - nbgtracks
1916  if(makeQAhistos) fHistMEphieta->Fill(DPhi,DEta, 1./(nMix*mixefficiency));
1917  } // end of background track loop
1918  } // end of filling mixed-event histo's
1919  } // end of check for biased jet triggers
1920  } // end of jet loop
1921  } // end of check for pool being ready
1922  } // end EMC triggered loop
1923 
1924 //=============================================================================================================
1925 
1926  // use only jets from EMCal-triggered events (for lhc11a use AliVEvent::kEMC1)
1928  // pp collisions
1929  if(GetBeamType() == 0) if(trigger && fTriggerEventType) { //kEMC1)) {
1930  if (poolpp->IsReady() || poolpp->NTracksInPool() > fNMIXtracks || poolpp->GetCurrentNEvents() >= fNMIXevents) {
1931 
1932  // loop over jets (passing cuts?)
1933  for (Int_t ijet = 0; ijet < Njets; ijet++) {
1934  Double_t leadjet=0;
1935  if (ijet==ijethi) leadjet=1;
1936 
1937  // get jet object
1938  AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
1939  if (!jet) continue;
1940 
1941  // (should probably be higher..., but makes a cut on jet pT)
1942  if (jet->Pt()<0.1) continue;
1943  if (!AcceptMyJet(jet)) continue;
1944 
1945  fHistEventQA->Fill(16); // event mixing jets
1946 
1947  // set cut to do event mixing only if we have a jet meeting our pt threshold (bias applied below)
1948  if (jet->Pt()<fJetPtcut) continue;
1949 
1950  // get number of current events in pool
1951  Int_t nMix = poolpp->GetCurrentNEvents(); // how many particles in pool to mix
1952 
1953  // Fill for biased jet triggers only
1954  if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)) { // && jet->Pt() > fJetPtcut) {
1955  // Fill mixed-event histos here
1956  for (Int_t jMix=0; jMix<nMix; jMix++) {
1957  fHistEventQA->Fill(17); // event mixing nMix
1958 
1959  // get jMix'th event
1960  bgTracks = poolpp->GetEvent(jMix);
1961  const Int_t Nbgtrks = bgTracks->GetEntries();
1962  for(Int_t ibg=0; ibg<Nbgtrks; ibg++) {
1963  AliPicoTrack *part = static_cast<AliPicoTrack*>(bgTracks->At(ibg));
1964  if(!part) continue;
1965  if(TMath::Abs(part->Eta())>0.9) continue;
1966  if(part->Pt()<0.15) continue;
1967 
1968  Double_t DEta = part->Eta()-jet->Eta(); // difference in eta
1969  Double_t DPhi = RelativePhi(jet->Phi(),part->Phi()); // difference in phi
1970  Double_t dEP = RelativeEPJET(jet->Phi(),fEPV0); // difference between jet and EP
1971  Double_t mixcharge = part->Charge();
1972  //Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta); // difference in R
1973 
1974  // create / fill mixed event sparse
1975  Double_t triggerEntries[9] = {fCent,jet->Pt(),part->Pt(),DEta,DPhi,dEP,zVtx, mixcharge, leadjet}; //array for ME sparse
1976  fhnMixedEvents->Fill(triggerEntries,1./nMix); // fill Sparse histo of mixed events
1977 
1978  fHistEventQA->Fill(18); // event mixing - nbgtracks
1979  if(makeextraCORRhistos) fHistMEphieta->Fill(DPhi,DEta, 1./nMix);
1980  } // end of background track loop
1981  } // end of filling mixed-event histo's
1982  } // end of check for biased jet triggers
1983  } // end of jet loop
1984  } // end of check for pool being ready
1985  } //end EMC triggered loop
1986 
1987  // pp collisions
1988  if(GetBeamType() == 0) {
1989 
1990  // use only tracks from MB events (for lhc11a use AliVEvent::kMB) //kAnyINT as test
1991  if(trigger && fMixingEventType) { //kMB) {
1992 
1993  // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool
1994  tracksClone = CloneAndReduceTrackList(tracks);
1995 
1996  // update pool if jet in event or not
1997  poolpp->UpdatePool(tracksClone);
1998 
1999  } // check on track from MB events
2000  }
2001 
2002  // PbPb collisions
2003  if(GetBeamType() == 1) {
2004 
2005  // use only tracks from MB and Central and Semi-Central events
2006  if(trigger && fMixingEventType) { //kMB) {
2007 
2008  // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool
2009  tracksClone = CloneAndReduceTrackList(tracks);
2010 
2011  // update pool if jet in event or not
2012  pool->UpdatePool(tracksClone);
2013 
2014  } // MB and Central and Semi-Central events
2015  } // PbPb collisions
2016  } // end of event mixing
2017 
2018  // print some stats on the event
2019  event++;
2020  fHistEventQA->Fill(19); // events making it to end
2021 
2022  // fill Event Trigger QA (after cuts)
2023  FillEventTriggerQA(fHistEventSelectionQAafterCuts, trig);
2024 
2025  if (doComments) {
2026  cout<<"Event #: "<<event<<" Jet Radius: "<<fJetRad<<" Constituent Pt Cut: "<<fConstituentCut<<endl;
2027  cout<<"# of jets: "<<Njets<<" NjetAcc: "<<NjetAcc<<" Highest jet pt: "<<highestjetpt<<" leading hadron pt: "<<leadhadronPT<<endl;
2028  cout<<"# tracks: "<<Ntracks<<" NtrackAcc: "<<NtrackAcc<<" Highest track pt: "<<ptmax<<endl;
2029  cout<<"maxcluster pt = "<<maxclusterpt<<" maxtrack pt = "<<maxtrackpt<<endl;
2030  cout<<" =============================================== "<<endl;
2031  }
2032 
2033  return kTRUE; // used when the function is of type bool
2034 } // end of RUN
2035 
2036 //________________________________________________________________________
2037 Int_t AliAnalysisTaskEmcalJetHadEPpid::GetCentBin(Double_t cent) const
2038 { // Get centrality bin.
2039  Int_t centbin = -1;
2040  if (cent>=0 && cent<10) centbin = 0;
2041  else if (cent>=10 && cent<20) centbin = 1;
2042  else if (cent>=20 && cent<30) centbin = 2;
2043  else if (cent>=30 && cent<40) centbin = 3;
2044  else if (cent>=40 && cent<50) centbin = 4;
2045  else if (cent>=50 && cent<90) centbin = 5;
2046 
2047  return centbin;
2048 }
2049 
2050 //________________________________________________________________________
2051 Double_t AliAnalysisTaskEmcalJetHadEPpid::RelativePhi(Double_t mphi,Double_t vphi) const
2052 { // function to calculate relative PHI
2053  double dphi = mphi-vphi;
2054 
2055  // set dphi to operate on adjusted scale
2056  if(dphi<-0.5*TMath::Pi()) dphi+=2.*TMath::Pi();
2057  if(dphi>3./2.*TMath::Pi()) dphi-=2.*TMath::Pi();
2058 
2059  // test
2060  if( dphi < -1.*TMath::Pi()/2 || dphi > 3.*TMath::Pi()/2 )
2061  AliWarning(Form("%s: dPHI not in range [-0.5*Pi, 1.5*Pi]!", GetName()));
2062 
2063  return dphi; // dphi in [-0.5Pi, 1.5Pi]
2064 }
2065 
2066 //_________________________________________________________________________
2067 Double_t AliAnalysisTaskEmcalJetHadEPpid::RelativeEPJET(Double_t jetAng, Double_t EPAng) const
2068 { // function to calculate angle between jet and EP in the 1st quadrant (0,Pi/2)
2069  Double_t dphi = (EPAng - jetAng);
2070 
2071  // ran into trouble with a few dEP<-Pi so trying this...
2072  if( dphi<-1*TMath::Pi() ){
2073  dphi = dphi + 1*TMath::Pi();
2074  } // this assumes we are doing full jets currently
2075 
2076  if( (dphi>0) && (dphi<1*TMath::Pi()/2) ){
2077  // Do nothing! we are in quadrant 1
2078  }else if( (dphi>1*TMath::Pi()/2) && (dphi<1*TMath::Pi()) ){
2079  dphi = 1*TMath::Pi() - dphi;
2080  }else if( (dphi<0) && (dphi>-1*TMath::Pi()/2) ){
2081  dphi = fabs(dphi);
2082  }else if( (dphi<-1*TMath::Pi()/2) && (dphi>-1*TMath::Pi()) ){
2083  dphi = dphi + 1*TMath::Pi();
2084  }
2085 
2086  // test
2087  if( dphi < 0 || dphi > TMath::Pi()/2 )
2088  AliWarning(Form("%s: dPHI not in range [0, 0.5*Pi]!", GetName()));
2089 
2090  return dphi; // dphi in [0, Pi/2]
2091 }
2092 
2093 //Int_t ieta=GetEtaBin(deta);
2094 //________________________________________________________________________
2095 Int_t AliAnalysisTaskEmcalJetHadEPpid::GetEtaBin(Double_t eta) const
2096 {
2097  // Get eta bin for histos.
2098  Int_t etabin = -1;
2099  if (TMath::Abs(eta)<=0.4) etabin = 0;
2100  else if (TMath::Abs(eta)>0.4 && TMath::Abs(eta)<0.8) etabin = 1;
2101  else if (TMath::Abs(eta)>=0.8) etabin = 2;
2102 
2103  return etabin;
2104 } // end of get-eta-bin
2105 
2106 //________________________________________________________________________
2107 Int_t AliAnalysisTaskEmcalJetHadEPpid::GetpTjetBin(Double_t pt) const
2108 {
2109  // Get jet pt bin for histos.
2110  Int_t ptbin = -1;
2111  if (pt>=15 && pt<20) ptbin = 0;
2112  else if (pt>=20 && pt<25) ptbin = 1;
2113  else if (pt>=25 && pt<40) ptbin = 2;
2114  else if (pt>=40 && pt<60) ptbin = 3;
2115  else if (pt>=60) ptbin = 4;
2116 
2117  return ptbin;
2118 } // end of get-jet-pt-bin
2119 
2120 //________________________________________________________________________
2121 Int_t AliAnalysisTaskEmcalJetHadEPpid::GetpTtrackBin(Double_t pt) const
2122 {
2123  // May need to update bins for future runs... (testing locally)
2124 
2125  // Get track pt bin for histos.
2126  Int_t ptbin = -1;
2127  if (pt < 0.5) ptbin = 0;
2128  else if (pt>=0.5 && pt<1.0) ptbin = 1;
2129  else if (pt>=1.0 && pt<1.5) ptbin = 2;
2130  else if (pt>=1.5 && pt<2.0) ptbin = 3;
2131  else if (pt>=2.0 && pt<2.5) ptbin = 4;
2132  else if (pt>=2.5 && pt<3.0) ptbin = 5;
2133  else if (pt>=3.0 && pt<4.0) ptbin = 6;
2134  else if (pt>=4.0 && pt<5.0) ptbin = 7;
2135  else if (pt>=5.0) ptbin = 8;
2136 
2137  return ptbin;
2138 } // end of get-jet-pt-bin
2139 
2140 //___________________________________________________________________________
2141 Int_t AliAnalysisTaskEmcalJetHadEPpid::GetzVertexBin(Double_t zVtx) const
2142 {
2143  // get z-vertex bin for histo.
2144  int zVbin= -1;
2145  if (zVtx>=-10 && zVtx<-8) zVbin = 0;
2146  else if (zVtx>=-8 && zVtx<-6) zVbin = 1;
2147  else if (zVtx>=-6 && zVtx<-4) zVbin = 2;
2148  else if (zVtx>=-4 && zVtx<-2) zVbin = 3;
2149  else if (zVtx>=-2 && zVtx<0) zVbin = 4;
2150  else if (zVtx>=0 && zVtx<2) zVbin = 5;
2151  else if (zVtx>=2 && zVtx<4) zVbin = 6;
2152  else if (zVtx>=4 && zVtx<6) zVbin = 7;
2153  else if (zVtx>=6 && zVtx<8) zVbin = 8;
2154  else if (zVtx>=8 && zVtx<10) zVbin = 9;
2155  else zVbin = 10;
2156 
2157  return zVbin;
2158 } // end of get z-vertex bin
2159 
2160 //______________________________________________________________________
2161 THnSparse* AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseF(const char* name, UInt_t entries)
2162 {
2163  // generate new THnSparseF, axes are defined in GetDimParams()
2164  Int_t count = 0;
2165  UInt_t tmp = entries;
2166  while(tmp!=0){
2167  count++;
2168  tmp = tmp &~ -tmp; // clear lowest bit
2169  }
2170 
2171  TString hnTitle(name);
2172  const Int_t dim = count;
2173  Int_t nbins[dim];
2174  Double_t xmin[dim];
2175  Double_t xmax[dim];
2176 
2177  Int_t i=0;
2178  Int_t c=0;
2179  while(c<dim && i<32){
2180  if(entries&(1<<i)){
2181  TString label("");
2182  GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
2183  hnTitle += Form(";%s",label.Data());
2184  c++;
2185  }
2186 
2187  i++;
2188  }
2189  hnTitle += ";";
2190 
2191  return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
2192 } // end of NewTHnSparseF
2193 
2194 void AliAnalysisTaskEmcalJetHadEPpid::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
2195 {
2196  // stores label and binning of axis for THnSparse
2197  const Double_t pi = TMath::Pi();
2198 
2199  switch(iEntry){
2200 
2201  case 0:
2202  label = "V0 centrality (%)";
2203  nbins = 10;
2204  xmin = 0.;
2205  xmax = 100.;
2206  break;
2207 
2208  case 1:
2209  if(fcorrJetPt) { // correct jet pt
2210  label = "Jet Corrected p_{T}";
2211  nbins = 50;
2212  xmin = -50.;
2213  xmax = 200.;
2214  } else { // don't correct jet pt
2215  label = "Jet p_{T}";
2216  if(doaltPIDbinning) { // minimize bins (for PID)
2217  nbins = 20;
2218  xmin = 0.;
2219  xmax = 100.;
2220  } else { // don't minimize bins
2221  nbins = 50;
2222  xmin = 0.;
2223  xmax = 250.;
2224  }
2225  }
2226  break;
2227 
2228  case 2:
2229  label = "Track p_{T}";
2230  if(doaltPIDbinning) { // minimize bins (for PID)
2231  nbins = 100;
2232  xmin = 0.;
2233  xmax = 10.;
2234  } else { // don't minimize bins
2235  nbins = 80;
2236  xmin = 0.;
2237  xmax = 20.;
2238  }
2239  break;
2240 
2241  case 3:
2242  label = "Relative Eta";
2243  nbins = 56; // 48
2244  xmin = -1.4;
2245  xmax = 1.4;
2246  break;
2247 
2248  case 4:
2249  label = "Relative Phi";
2250  nbins = 72;
2251  xmin = -0.5*pi;
2252  xmax = 1.5*pi;
2253  break;
2254 
2255  case 5:
2256  label = "Relative angle of Jet and Reaction Plane";
2257  nbins = 3; // (12) 72
2258  xmin = 0;
2259  xmax = 0.5*pi;
2260  break;
2261 
2262  case 6:
2263  label = "z-vertex";
2264  nbins = 10;
2265  xmin = -10;
2266  xmax = 10;
2267  break;
2268 
2269  case 7:
2270  label = "track charge";
2271  nbins = 3;
2272  xmin = -1.5;
2273  xmax = 1.5;
2274  break;
2275 
2276  case 8:
2277  label = "leading jet";
2278  if(doaltPIDbinning) nbins = 1;
2279  else nbins = 3;
2280  xmin = -0.5;
2281  xmax = 2.5;
2282  break;
2283 
2284  case 9: // need to update
2285  label = "leading track";
2286  if(doaltPIDbinning) nbins = 1;
2287  else nbins = 10;
2288  xmin = 0;
2289  xmax = 50;
2290  break;
2291 
2292  } // end of switch
2293 } // end of getting dim-params
2294 
2295 //_________________________________________________
2296 // From CF event mixing code PhiCorrelations
2297 TObjArray* AliAnalysisTaskEmcalJetHadEPpid::CloneAndReduceTrackList(TObjArray* tracksME)
2298 {
2299  // clones a track list by using AliPicoTrack which uses much less memory (used for event mixing)
2300  TObjArray* tracksClone = new TObjArray;
2301  tracksClone->SetOwner(kTRUE);
2302 
2303 // ===============================
2304 
2305 // cout << "RM Hybrid track : " << i << " " << particle->Pt() << endl;
2306 
2307  //cout << "nEntries " << tracks->GetEntriesFast() <<endl;
2308  for (Int_t i=0; i<tracksME->GetEntriesFast(); i++) { // AOD/general case
2309  AliVParticle* particle = (AliVParticle*) tracksME->At(i); // AOD/general case
2310  if(TMath::Abs(particle->Eta())>fTrkEta) continue;
2311  if(particle->Pt()<0.15)continue;
2312 
2313 /*
2314 // DON'T USE
2315  Double_t trackpt=particle->Pt(); // track pT
2316 
2317  Int_t trklabel=-1;
2318  trklabel=particle->GetLabel();
2319  //cout << "TRACK_LABEL: " << particle->GetLabel()<<endl;
2320 
2321  Int_t hadbin=-1;
2322  if(trackpt<0.5) hadbin=0;
2323  else if(trackpt<1) hadbin=1;
2324  else if(trackpt<2) hadbin=2;
2325  else if(trackpt<3) hadbin=3;
2326  else if(trackpt<5) hadbin=4;
2327  else if(trackpt<8) hadbin=5;
2328  else if(trackpt<20) hadbin=6;
2329 // end of DON'T USE
2330 
2331 //feb10 comment out
2332  if(hadbin>-1 && trklabel>-1 && trklabel <3) fHistTrackEtaPhi[trklabel][hadbin]->Fill(particle->Eta(),particle->Phi());
2333  if(hadbin>-1) fHistTrackEtaPhi[3][hadbin]->Fill(particle->Eta(),particle->Phi());
2334 
2335  if(hadbin>-1) fHistTrackEtaPhi[hadbin]->Fill(particle->Eta(),particle->Phi()); // TEST
2336 */
2337 
2338  tracksClone->Add(new AliPicoTrack(particle->Pt(), particle->Eta(), particle->Phi(), particle->Charge(), 0, 0, 0, 0));
2339  } // end of looping through tracks
2340 
2341  return tracksClone;
2342 }
2343 
2344 //____________________________________________________________________________________________
2345 THnSparse* AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseFPID(const char* name, UInt_t entries)
2346 {
2347  // generate new THnSparseF PID, axes are defined in GetDimParams()
2348  Int_t count = 0;
2349  UInt_t tmp = entries;
2350  while(tmp!=0){
2351  count++;
2352  tmp = tmp &~ -tmp; // clear lowest bit
2353  }
2354 
2355  TString hnTitle(name);
2356  const Int_t dim = count;
2357  Int_t nbins[dim];
2358  Double_t xmin[dim];
2359  Double_t xmax[dim];
2360 
2361  Int_t i=0;
2362  Int_t c=0;
2363  while(c<dim && i<32){
2364  if(entries&(1<<i)){
2365  TString label("");
2366  GetDimParamsPID(i, label, nbins[c], xmin[c], xmax[c]);
2367  hnTitle += Form(";%s",label.Data());
2368  c++;
2369  }
2370 
2371  i++;
2372  }
2373  hnTitle += ";";
2374 
2375  return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
2376 } // end of NewTHnSparseF PID
2377 
2378 //________________________________________________________________________________
2379 void AliAnalysisTaskEmcalJetHadEPpid::GetDimParamsPID(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
2380 {
2381  // stores label and binning of axis for THnSparse
2382  const Double_t pi = TMath::Pi();
2383 
2384  switch(iEntry){
2385 
2386  case 0:
2387  label = "V0 centrality (%)";
2388  nbins = 10;
2389  xmin = 0.;
2390  xmax = 100.;
2391  break;
2392 
2393  case 1:
2394  label = "Track p_{T}";
2395  if(doaltPIDbinning) { // minimize bins
2396  nbins = 100;
2397  xmin = 0.;
2398  xmax = 10.;
2399  } else { // don't minmize
2400  nbins = 80;
2401  xmin = 0.;
2402  xmax = 20.;
2403  }
2404  break;
2405 
2406  case 2:
2407  label = "Charge of Track";
2408  nbins = 3;
2409  xmin = -1.5;
2410  xmax = 1.5;
2411  break;
2412 
2413  case 3:
2414  label = "Relative Eta of Track and Jet";
2415  nbins = 56; // 48
2416  xmin = -1.4;
2417  xmax = 1.4;
2418  break;
2419 
2420  case 4:
2421  label = "Relative Phi of Track and Jet";
2422  nbins = 72;
2423  xmin = -0.5*pi;
2424  xmax = 1.5*pi;
2425  break;
2426 
2427  case 5: // not used
2428  label = "leading jet";
2429  if(doaltPIDbinning) nbins = 1;
2430  else nbins = 3;
2431  xmin = -.5;
2432  xmax = 2.5;
2433  break;
2434 
2435  case 6:
2436  label = "Z-vertex";
2437  nbins = 10;
2438  xmin = -10.;
2439  xmax = 10.;
2440  break;
2441 
2442  case 7:
2443  label = "Relative angle: Jet and Reaction Plane";
2444  nbins = 3; // (12) 48
2445  xmin = 0.;
2446  xmax = 0.5*pi;
2447  break;
2448 
2449  case 8:
2450  if(fcorrJetPt) { // correct jet pt
2451  label = "Jet Corrected p_{T}";
2452  nbins = 50;
2453  xmin = -50.;
2454  xmax = 200.;
2455  } else { // don't correct jet pt
2456  label = "Jet p_{T}";
2457  if(doaltPIDbinning) { // minimize bins
2458  nbins = 20;
2459  xmin = 0.;
2460  xmax = 100.;
2461  } else { // don't minimize bins
2462  nbins = 50;
2463  xmin = 0.;
2464  xmax = 250.;
2465  }
2466  }
2467  break;
2468 
2469  case 9: // this may be temp, only using TOF
2470  label = "N-Sigma of pions in TPC";
2471  if(doaltPIDbinning) nbins = 1;
2472  else nbins = 200;
2473  xmin = -10.0;
2474  xmax = 10.0;
2475  break;
2476 
2477  case 10:
2478  label = "N-Sigma of pions in TOF";
2479  if(doaltPIDbinning) { // this may be temp, only using TOF
2480  nbins = 10;
2481  xmin = -5.;
2482  xmax = 5.;
2483  } else {
2484  nbins = 200;
2485  xmin = -10.;
2486  xmax = 10.;
2487  }
2488  break;
2489 
2490  case 11:
2491  label = "PID determination TPC 1-15, TOF 21-25";
2492  nbins = 25; //5;
2493  xmin = 0.5; //0.;
2494  xmax = 25.5; //5.;
2495  break;
2496 
2497 /*
2498  case 12:
2499  label = "ITS PID determination";
2500  nbins = 5;
2501  xmin = 0.;
2502  xmax = 5.;
2503  break;
2504 
2505  case 13:
2506  label = "TOF PID determination";
2507  nbins = 5;
2508  xmin = 0.;
2509  xmax = 5.;
2510  break;
2511 */
2512 
2513  case 12: // this may be temp, only using TOF
2514  label = "N-Sigma of protons in TPC";
2515  if(doaltPIDbinning) nbins = 1;
2516  else nbins = 200;
2517  xmin = -10.;
2518  xmax = 10.;
2519  break;
2520 
2521  case 13: // this may be temp, only using TOF
2522  label = "N-Sigma of kaons in TPC";
2523  if(doaltPIDbinning) nbins = 1;
2524  else nbins = 200;
2525  xmin = -10.;
2526  xmax = 10.;
2527  break;
2528 
2529  case 14: // this may be temp, only using TOF
2530  label = "N-Sigma of pions in ITS";
2531  if(doaltPIDbinning) nbins = 1;
2532  else nbins = 200;
2533  xmin = -10.0;
2534  xmax = 10.0;
2535  break;
2536 
2537  case 15: // this may be temp, only using TOF
2538  label = "N-Sigma of protons in ITS";
2539  if(doaltPIDbinning) nbins = 1;
2540  else nbins = 200;
2541  xmin = -10.;
2542  xmax = 10.;
2543  break;
2544 
2545  case 16: // this may be temp, only using TOF
2546  label = "N-Sigma of kaons in ITS";
2547  if(doaltPIDbinning) nbins = 1;
2548  else nbins = 200;
2549  xmin = -10.;
2550  xmax = 10.;
2551  break;
2552 
2553  case 17:
2554  label = "N-Sigma of protons in TOF";
2555  if(doaltPIDbinning) {
2556  nbins = 10;
2557  xmin = -5.;
2558  xmax = 5.;
2559  } else {
2560  nbins = 200;
2561  xmin = -10.;
2562  xmax = 10.;
2563  }
2564  break;
2565 
2566  case 18:
2567  label = "N-Sigma of kaons in TOF";
2568  if(doaltPIDbinning) {
2569  nbins = 10;
2570  xmin = -5.;
2571  xmax = 5.;
2572  } else {
2573  nbins = 200;
2574  xmin = -10.;
2575  xmax = 10.;
2576  }
2577  break;
2578 
2579  } // end of switch
2580 } // end of get dimension parameters PID
2581 
2582 void AliAnalysisTaskEmcalJetHadEPpid::Terminate(Option_t *) {
2583  cout<<"#########################"<<endl;
2584  cout<<"#### DONE RUNNING!!! ####"<<endl;
2585  cout<<"#########################"<<endl;
2586 } // end of terminate
2587 
2588 //________________________________________________________________________
2589 Int_t AliAnalysisTaskEmcalJetHadEPpid::AcceptMyJet(AliEmcalJet *jet) {
2590  //applies all jet cuts except pt
2591  if ((jet->Phi()<fPhimin)||(jet->Phi()>fPhimax)) return 0;
2592  if ((jet->Eta()<fEtamin)||(jet->Eta()>fEtamax)) return 0;
2593  if (jet->Area()<fAreacut) return 0;
2594  // prevents 0 area jets from sneaking by when area cut == 0
2595  if (jet->Area()==0) return 0;
2596  //exclude jets with extremely high pt tracks which are likely misreconstructed
2597  if(jet->MaxTrackPt()>100) return 0;
2598 
2599  //passed all above cuts
2600  return 1;
2601 }
2602 
2603 //void AliAnalysisTaskEmcalJetHadEPpid::FillAnalysisSummaryHistogram() const
2604 void AliAnalysisTaskEmcalJetHadEPpid::SetfHistPIDcounterLabels(TH1* h) const
2605 {
2606  // fill the analysis summary histrogram, saves all relevant analysis settigns
2607  h->GetXaxis()->SetBinLabel(1, "TPC: Unidentified");
2608  h->GetXaxis()->SetBinLabel(2, "TPC: Pion");
2609  h->GetXaxis()->SetBinLabel(3, "TPC: Kaon");
2610  h->GetXaxis()->SetBinLabel(4, "TPC: Proton");
2611  h->GetXaxis()->SetBinLabel(5, "ITS: Unidentified");
2612  h->GetXaxis()->SetBinLabel(6, "ITS: Pion");
2613  h->GetXaxis()->SetBinLabel(7, "ITS: Kaon");
2614  h->GetXaxis()->SetBinLabel(8, "ITS: Proton");
2615  h->GetXaxis()->SetBinLabel(9, "TOF: Unidentified");
2616  h->GetXaxis()->SetBinLabel(10, "TOF: Pion");
2617  h->GetXaxis()->SetBinLabel(11, "TOF: Kaon");
2618  h->GetXaxis()->SetBinLabel(12, "TOF: Proton");
2619  h->GetXaxis()->SetBinLabel(14, "Unidentified tracks");
2620 
2621  // set x-axis labels vertically
2622  h->LabelsOption("v");
2623 }
2624 
2625 //void AliAnalysisTaskEmcalJetHadEPpid::FillAnalysisSummaryHistogram() const
2626 void AliAnalysisTaskEmcalJetHadEPpid::SetfHistQAcounterLabels(TH1* h) const
2627 {
2628  // label bins of the analysis event summary
2629  h->GetXaxis()->SetBinLabel(1, "All events started");
2630  h->GetXaxis()->SetBinLabel(2, "object check");
2631  h->GetXaxis()->SetBinLabel(3, "aod/esd check");
2632  h->GetXaxis()->SetBinLabel(4, "centrality check");
2633  h->GetXaxis()->SetBinLabel(5, "zvertex check");
2634  h->GetXaxis()->SetBinLabel(6, "list check");
2635  h->GetXaxis()->SetBinLabel(7, "track/jet pointer check");
2636  h->GetXaxis()->SetBinLabel(8, "tracks & jets < than 1 check");
2637  h->GetXaxis()->SetBinLabel(9, "after track/jet loop to get highest pt");
2638  h->GetXaxis()->SetBinLabel(10, "accepted jets");
2639  h->GetXaxis()->SetBinLabel(11, "jets meeting pt threshold");
2640  h->GetXaxis()->SetBinLabel(12, "accepted tracks in events w/ trigger jet");
2641  h->GetXaxis()->SetBinLabel(13, "after AliVEvent & fPIDResponse");
2642  h->GetXaxis()->SetBinLabel(14, "events before event mixing");
2643  h->GetXaxis()->SetBinLabel(15, "mixed events w/ pool");
2644  h->GetXaxis()->SetBinLabel(16, "event mixing: jets");
2645  h->GetXaxis()->SetBinLabel(17, "event mixing: nMix");
2646  h->GetXaxis()->SetBinLabel(18, "event mixing: nbackground tracks");
2647  h->GetXaxis()->SetBinLabel(19, "event mixing: THE END");
2648 
2649  // set x-axis labels vertically
2650  h->LabelsOption("v");
2651 }
2652 
2653 void AliAnalysisTaskEmcalJetHadEPpid::SetfHistEvtSelQALabels(TH1* h) const
2654 {
2655  // label bins of the analysis trigger selection summary
2656  h->GetXaxis()->SetBinLabel(1, "no trigger");
2657  h->GetXaxis()->SetBinLabel(2, "kAny");
2658  h->GetXaxis()->SetBinLabel(3, "kAnyINT");
2659  h->GetXaxis()->SetBinLabel(4, "kMB");
2660  h->GetXaxis()->SetBinLabel(5, "kINT7");
2661  h->GetXaxis()->SetBinLabel(6, "kEMC1");
2662  h->GetXaxis()->SetBinLabel(7, "kEMC7");
2663  h->GetXaxis()->SetBinLabel(8, "kEMC8");
2664  h->GetXaxis()->SetBinLabel(9, "kEMCEJE");
2665  h->GetXaxis()->SetBinLabel(10, "kEMCEGA");
2666  h->GetXaxis()->SetBinLabel(11, "kCentral");
2667  h->GetXaxis()->SetBinLabel(12, "kSemiCentral");
2668  h->GetXaxis()->SetBinLabel(13, "kINT8");
2669  h->GetXaxis()->SetBinLabel(14, "kEMCEJE or kMB");
2670  h->GetXaxis()->SetBinLabel(15, "kEMCEGA or kMB");
2671  h->GetXaxis()->SetBinLabel(16, "kAnyINT or kMB");
2672  h->GetXaxis()->SetBinLabel(17, "kEMCEJE & (kMB or kCentral or kSemiCentral)");
2673  h->GetXaxis()->SetBinLabel(18, "kEMCEGA & (kMB or kCentral or kSemiCentral)");
2674  h->GetXaxis()->SetBinLabel(19, "kAnyINT & (kMB or kCentral or kSemiCentral)");
2675 
2676  // set x-axis labels vertically
2677  h->LabelsOption("v");
2678  //h->LabelsDeflate("X");
2679 }
2680 
2681 //______________________________________________________________________
2682 THnSparse* AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseFCorr(const char* name, UInt_t entries) {
2683  // generate new THnSparseD, axes are defined in GetDimParamsD()
2684  Int_t count = 0;
2685  UInt_t tmp = entries;
2686  while(tmp!=0){
2687  count++;
2688  tmp = tmp &~ -tmp; // clear lowest bit
2689  }
2690 
2691  TString hnTitle(name);
2692  const Int_t dim = count;
2693  Int_t nbins[dim];
2694  Double_t xmin[dim];
2695  Double_t xmax[dim];
2696 
2697  Int_t i=0;
2698  Int_t c=0;
2699  while(c<dim && i<32){
2700  if(entries&(1<<i)){
2701  TString label("");
2702  GetDimParamsCorr(i, label, nbins[c], xmin[c], xmax[c]);
2703  hnTitle += Form(";%s",label.Data());
2704  c++;
2705  }
2706 
2707  i++;
2708  }
2709  hnTitle += ";";
2710 
2711  return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
2712 } // end of NewTHnSparseF
2713 
2714 //______________________________________________________________________________________________
2715 void AliAnalysisTaskEmcalJetHadEPpid::GetDimParamsCorr(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
2716 {
2717  //stores label and binning of axis for THnSparse
2718  const Double_t pi = TMath::Pi();
2719 
2720  switch(iEntry){
2721 
2722  case 0:
2723  label = "V0 centrality (%)";
2724  nbins = 10;
2725  xmin = 0.;
2726  xmax = 100.;
2727  break;
2728 
2729  case 1:
2730  if(fcorrJetPt) { // correct jet pt
2731  label = "Jet Corrected p_{T}";
2732  nbins = 50;
2733  xmin = -50.;
2734  xmax = 200.;
2735  } else { // don't correct jet pt
2736  label = "Jet p_{T}";
2737  if(doaltPIDbinning) { // minimize bins
2738  nbins = 20;
2739  xmin = 0.;
2740  xmax = 100.;
2741  } else { // don't minimize bins
2742  nbins = 50;
2743  xmin = 0.;
2744  xmax = 250.;
2745  }
2746  }
2747  break;
2748 
2749  case 2:
2750  label = "Relative angle: Jet and Reaction Plane";
2751  nbins = 3; // (12) 48
2752  xmin = 0.;
2753  xmax = 0.5*pi;
2754  break;
2755 
2756  case 3:
2757  label = "Z-vertex";
2758  nbins = 10;
2759  xmin = -10.;
2760  xmax = 10.;
2761  break;
2762 
2763  case 4:
2764  label = "Jet p_{T} corrected with Local Rho";
2765  if(doaltPIDbinning) {
2766  nbins = 30; // 250
2767  xmin = -50.;
2768  xmax = 100.;
2769  break;
2770  } else {
2771  nbins = 50; // 250
2772  xmin = -50.;
2773  xmax = 200.;
2774  }
2775  break;
2776 
2777  case 5:
2778  label = "Jet p_{T} corrected with Global Rho";
2779  if(doaltPIDbinning) {
2780  nbins = 30; // 250
2781  xmin = -50.;
2782  xmax = 100.;
2783  break;
2784  } else {
2785  nbins = 50; // 250
2786  xmin = -50.;
2787  xmax = 200.;
2788  }
2789  break;
2790 
2791  }// end of switch
2792 } // end of Correction (ME) sparse
2793 
2794 //________________________________________________________________________
2795 //Int_t AliAnalysisTaskEmcalJetHadEPpid::AcceptFlavourJet(AliEmcalJet* fljet, Int_t NUM, Int_t NUM2, Int_t NUM3) {
2796 Int_t AliAnalysisTaskEmcalJetHadEPpid::AcceptFlavourJet(AliEmcalJet* fljet, Int_t NUM) {
2797  // Get jet if accepted for given flavour tag
2798  // If jet not accepted return 0
2799  if(!fljet) {
2800  AliError(Form("%s:Jet not found",GetName()));
2801  return 0;
2802  }
2803 
2804  Int_t flavNUM = -99;//, flavNUM2 = -99, flavNUM3 = -99; FIXME commented out to avoid compiler warning
2805  flavNUM = NUM;
2806  //flavNUM2 = NUM2;
2807  //flavNUM3 = NUM3;
2808 
2809 /*
2810  // from the AliEmcalJet class, the tagging enumerator
2811  enum EFlavourTag{
2812  kDStar = 1<<0, kD0 = 1<<1,
2813  kSig1 = 1<<2, kSig2 = 1<<3,
2814  kBckgrd1 = 1<<4, kBckgrd2 = 1<<5, kBckgrd3 = 1<<6
2815  };
2816  // bit 0 = no tag, bit 1 = Dstar, bit 2 = D0 and so forth...
2817 */
2818 
2819  // get flavour of jet (if any)
2820  Int_t flav = -999;
2821  flav = fljet->GetFlavour();
2822 
2823  // cases (for now..)
2824  // 3 = electron rich, 5 = hadron (bkgrd) rich
2825  // if flav < 1, its not tagged, so return kFALSE (0)
2826  if(flav < 1) return 0;
2827 
2828  // if flav is not equal to what we want then return kFALSE (0)
2829  //if(flav != flavNUM && flav != flavNUM2 && flav != flavNUM3) return 0;
2830  if(flav != flavNUM) return 0;
2831 
2832  // we have the flavour we want, so return kTRUE (1)
2833  //if(flav == flavNUM || flav == flavNUM2 || flav == flavNUM3) return 1;
2834  if(flav == flavNUM) return 1;
2835 
2836  // we by default have a flavour tagged jet
2837  return 1;
2838 }
2839 
2840 //________________________________________________________________________
2841 Double_t AliAnalysisTaskEmcalJetHadEPpid::EffCorrection(Double_t trackETA, Double_t trackPT, Int_t effSwitch) const {
2842  // default (current) parameters
2843 
2844 /*
2845  // min/max range of x & y (track Eta & track pt) used by Nat
2846  Double_t trETAmin = -0.9; Double_t trETAmax = 0.9;
2847  Double_t trPTmin = 0.2; Double_t trPTmax = 12.0; ??
2848 */
2849 
2850  // x-variable = track pt, y-variable = track eta
2851  Double_t x = trackPT;
2852  Double_t y = trackETA;
2853  //y = TMath::Abs(trackETA);
2854  Double_t TRefficiency = -999;
2855  Int_t runNUM = fCurrentRunNumber;
2856  Int_t runSwitchGood = -999;
2857  Int_t centbin = -99;
2858 
2859  Double_t etaaxis = 0;
2860  Double_t ptaxis = 0;
2861 
2862  if(effSwitch < 1) {
2863  // Semi-GOOD IROC C13 runlists
2864  // DO NOT USE THESE:
2865  // Runs: (169040, 169044, 169045, 169099, 169418, 169419, 169420, 169475, 169498, 169504, 169506, 169512, 169515, 169550, 169553, 169554, 169555, 169557, 169584, 169586, 169587, 169588, 169590, 169591)
2866 
2867  // Semi-GOOD OROC C08 runlists
2868  if ((runNUM == 169975 || runNUM == 169981 || runNUM == 170038 || runNUM == 170040 || runNUM == 170083 || runNUM == 170084 || runNUM == 170085 || runNUM == 170088 || runNUM == 170089 || runNUM == 170091 || runNUM == 170152 || runNUM == 170155 || runNUM == 170159 || runNUM == 170163 || runNUM == 170193 || runNUM == 170195 || runNUM == 170203 || runNUM == 170204 || runNUM == 170228 || runNUM == 170230 || runNUM == 170268 || runNUM == 170269 || runNUM == 170270 || runNUM == 170306 || runNUM == 170308 || runNUM == 170309)) runSwitchGood = 0;
2869 
2870  // GOOD runlists
2871  if ((runNUM == 167902 || runNUM == 167903 || runNUM == 167915 || runNUM == 167920 || runNUM == 167987 || runNUM == 167988 || runNUM == 168066 || runNUM == 168068 || runNUM == 168069 || runNUM == 168076 || runNUM == 168104 || runNUM == 168107 || runNUM == 168108 || runNUM == 168115 || runNUM == 168212 || runNUM == 168310 || runNUM == 168311 || runNUM == 168322 || runNUM == 168325 || runNUM == 168341 || runNUM == 168342 || runNUM == 168361 || runNUM == 168362 || runNUM == 168458 || runNUM == 168460 || runNUM == 168461 || runNUM == 168464 || runNUM == 168467 || runNUM == 168511 || runNUM == 168512 || runNUM == 168777 || runNUM == 168826 || runNUM == 168984 || runNUM == 168988 || runNUM == 168992 || runNUM == 169035 || runNUM == 169091 || runNUM == 169094 || runNUM == 169138 || runNUM == 169143 || runNUM == 169144 || runNUM == 169145 || runNUM == 169148 || runNUM == 169156 || runNUM == 169160 || runNUM == 169167 || runNUM == 169238 || runNUM == 169411 || runNUM == 169415 || runNUM == 169417 || runNUM == 169835 || runNUM == 169837 || runNUM == 169838 || runNUM == 169846 || runNUM == 169855 || runNUM == 169858 || runNUM == 169859 || runNUM == 169923 || runNUM == 169956 || runNUM == 170027 || runNUM == 170036 || runNUM == 170081)) runSwitchGood = 1;
2872 
2873  if(fCent>=0 && fCent<10) centbin = 0;
2874  else if (fCent>=10 && fCent<30) centbin = 1;
2875  else if (fCent>=30 && fCent<50) centbin = 2;
2876  else if (fCent>=50 && fCent<90) centbin = 3;
2877 
2878  if(runSwitchGood == 0 && centbin == 0) effSwitch = 2;
2879  if(runSwitchGood == 0 && centbin == 1) effSwitch = 3;
2880  if(runSwitchGood == 0 && centbin == 2) effSwitch = 4;
2881  if(runSwitchGood == 0 && centbin == 3) effSwitch = 5;
2882  if(runSwitchGood == 1 && centbin == 0) effSwitch = 6;
2883  if(runSwitchGood == 1 && centbin == 1) effSwitch = 7;
2884  if(runSwitchGood == 1 && centbin == 2) effSwitch = 8;
2885  if(runSwitchGood == 1 && centbin == 3) effSwitch = 9;
2886  }
2887 
2888  // 0-10% centrality: Semi-Good Runs
2889  Double_t p0_10SG[17] = {0.906767, 0.0754127, 1.11638, -0.0233078, 0.795454, 0.00935385, -0.000327857, 1.08903, 0.0107272, 0.443252, -0.143411, 0.965822, 0.359156, -0.581221, 1.0739, 0.00632828, 0.706356};
2890  // 10-30% centrality: Semi-Good Runs
2891  Double_t p10_30SG[17] = {0.908011, 0.0769254, 1.11912, -0.0249449, 0.741488, 0.0361252, -0.00367954, 1.10424, 0.011472, 0.452059, -0.133282, 0.980633, 0.358222, -0.620256, 1.06871, 0.00564449, 0.753168};
2892  // 30-50% centrality: Semi-Good Runs
2893  Double_t p30_50SG[17] = {0.958708, 0.0799197, 1.10817, -0.0357678, 0.75051, 0.0607808, -0.00929713, 0.998801, 0.00692244, 0.615452, -0.0480328, 0.968431, 0.321634, -0.619066, 1.03412, 0.00656201, 0.798666};
2894  // 50-90% centrality: Semi-Good Runs
2895  Double_t p50_90SG[17] = {0.944565, 0.0807258, 1.12709, -0.0324746, 0.666452, 0.0842476, -0.00963837, 1.02829, 0.00666852, 0.549625, -0.0603107, 0.981374, 0.309374, -0.619181, 1.05367, 0.005925, 0.744887};
2896 
2897  // 0-10% centrality: Good Runs
2898  Double_t p0_10G[17] = {0.971679, 0.0767571, 1.13355, -0.0274484, 0.856652, 0.00536795, 3.90795e-05, 1.06889, 0.011007, 0.447046, -0.146626, 0.919777, 0.192601, -0.268515, 1.00243, 0.00620849, 0.709477};
2899  // 10-30% centrality: Good Runs
2900  Double_t p10_30G[17] = {0.97929, 0.0776039, 1.12213, -0.0300645, 0.844722, 0.0134788, -0.0012333, 1.07955, 0.0116835, 0.456608, -0.132743, 0.930964, 0.174175, -0.267154, 0.993118, 0.00574892, 0.765256};
2901  // 30-50% centrality: Good Runs
2902  Double_t p30_50G[17] = {0.997696, 0.0816769, 1.14341, -0.0353734, 0.752151, 0.0744259, -0.0102926, 1.01561, 0.00713274, 0.57203, -0.0640248, 0.947747, 0.102007, -0.194698, 0.999164, 0.00568476, 0.7237};
2903  // 50-90% centrality: Good Runs
2904  Double_t p50_90G[17] = {0.97041, 0.0813559, 1.12151, -0.0368797, 0.709327, 0.0701501, -0.00784043, 1.06276, 0.00676173, 0.53607, -0.0703117, 0.982534, 0.0947881, -0.18073, 1.03229, 0.00580109, 0.737801};
2905 
2906  // set up a switch for different parameter values...
2907  switch(effSwitch) {
2908  case 1 :
2909  // first switch value - TRefficiency not used so = 1
2910  TRefficiency = 1.0;
2911  break;
2912 
2913  case 2 :
2914  // Parameter values for Semi-GOOD TPC (LHC11h) runs (0-10%):
2915  ptaxis = (x<2.9)*(p0_10SG[0]*exp(-pow(p0_10SG[1]/x,p0_10SG[2])) + p0_10SG[3]*x) + (x>=2.9)*(p0_10SG[4] + p0_10SG[5]*x + p0_10SG[6]*x*x);
2916  etaaxis = (y<-0.07)*(p0_10SG[7]*exp(-pow(p0_10SG[8]/TMath::Abs(y+0.91),p0_10SG[9])) + p0_10SG[10]*y) + (y>=-0.07 && y<=0.4)*(p0_10SG[11] + p0_10SG[12]*y + p0_10SG[13]*y*y) + (y>0.4)*(p0_10SG[14]*exp(-pow(p0_10SG[15]/TMath::Abs(-y+0.91),p0_10SG[16])));
2917  TRefficiency = ptaxis*etaaxis;
2918  break;
2919 
2920  case 3 :
2921  // Parameter values for Semi-GOOD TPC (LHC11h) runs (10-30%):
2922  ptaxis = (x<2.9)*(p10_30SG[0]*exp(-pow(p10_30SG[1]/x,p10_30SG[2])) + p10_30SG[3]*x) + (x>=2.9)*(p10_30SG[4] + p10_30SG[5]*x + p10_30SG[6]*x*x);
2923  etaaxis = (y<-0.07)*(p10_30SG[7]*exp(-pow(p10_30SG[8]/TMath::Abs(y+0.91),p10_30SG[9])) + p10_30SG[10]*y) + (y>=-0.07 && y<=0.4)*(p10_30SG[11] + p10_30SG[12]*y + p10_30SG[13]*y*y) + (y>0.4)*(p10_30SG[14]*exp(-pow(p10_30SG[15]/TMath::Abs(-y+0.91),p10_30SG[16])));
2924  TRefficiency = ptaxis*etaaxis;
2925  break;
2926 
2927  case 4 :
2928  // Parameter values for Semi-GOOD TPC (LHC11h) runs (30-50%):
2929  ptaxis = (x<2.9)*(p30_50SG[0]*exp(-pow(p30_50SG[1]/x,p30_50SG[2])) + p30_50SG[3]*x) + (x>=2.9)*(p30_50SG[4] + p30_50SG[5]*x + p30_50SG[6]*x*x);
2930  etaaxis = (y<-0.07)*(p30_50SG[7]*exp(-pow(p30_50SG[8]/TMath::Abs(y+0.91),p30_50SG[9])) + p30_50SG[10]*y) + (y>=-0.07 && y<=0.4)*(p30_50SG[11] + p30_50SG[12]*y + p30_50SG[13]*y*y) + (y>0.4)*(p30_50SG[14]*exp(-pow(p30_50SG[15]/TMath::Abs(-y+0.91),p30_50SG[16])));
2931  TRefficiency = ptaxis*etaaxis;
2932  break;
2933 
2934  case 5 :
2935  // Parameter values for Semi-GOOD TPC (LHC11h) runs (50-90%):
2936  ptaxis = (x<2.9)*(p50_90SG[0]*exp(-pow(p50_90SG[1]/x,p50_90SG[2])) + p50_90SG[3]*x) + (x>=2.9)*(p50_90SG[4] + p50_90SG[5]*x + p50_90SG[6]*x*x);
2937  etaaxis = (y<-0.07)*(p50_90SG[7]*exp(-pow(p50_90SG[8]/TMath::Abs(y+0.91),p50_90SG[9])) + p50_90SG[10]*y) + (y>=-0.07 && y<=0.4)*(p50_90SG[11] + p50_90SG[12]*y + p50_90SG[13]*y*y) + (y>0.4)*(p50_90SG[14]*exp(-pow(p50_90SG[15]/TMath::Abs(-y+0.91),p50_90SG[16])));
2938  TRefficiency = ptaxis*etaaxis;
2939  break;
2940 
2941  case 6 :
2942  // Parameter values for GOOD TPC (LHC11h) runs (0-10%):
2943  ptaxis = (x<2.9)*(p0_10G[0]*exp(-pow(p0_10G[1]/x,p0_10G[2])) + p0_10G[3]*x) + (x>=2.9)*(p0_10G[4] + p0_10G[5]*x + p0_10G[6]*x*x);
2944  etaaxis = (y<0.0)*(p0_10G[7]*exp(-pow(p0_10G[8]/TMath::Abs(y+0.91),p0_10G[9])) + p0_10G[10]*y) + (y>=0.0 && y<=0.4)*(p0_10G[11] + p0_10G[12]*y + p0_10G[13]*y*y) + (y>0.4)*(p0_10G[14]*exp(-pow(p0_10G[15]/TMath::Abs(-y+0.91),p0_10G[16])));
2945  TRefficiency = ptaxis*etaaxis;
2946  break;
2947 
2948  case 7 :
2949  // Parameter values for GOOD TPC (LHC11h) runs (10-30%):
2950  ptaxis = (x<2.9)*(p10_30G[0]*exp(-pow(p10_30G[1]/x,p10_30G[2])) + p10_30G[3]*x) + (x>=2.9)*(p10_30G[4] + p10_30G[5]*x + p10_30G[6]*x*x);
2951  etaaxis = (y<0.0)*(p10_30G[7]*exp(-pow(p10_30G[8]/TMath::Abs(y+0.91),p10_30G[9])) + p10_30G[10]*y) + (y>=0.0 && y<=0.4)*(p10_30G[11] + p10_30G[12]*y + p10_30G[13]*y*y) + (y>0.4)*(p10_30G[14]*exp(-pow(p10_30G[15]/TMath::Abs(-y+0.91),p10_30G[16])));
2952  TRefficiency = ptaxis*etaaxis;
2953  break;
2954 
2955  case 8 :
2956  // Parameter values for GOOD TPC (LHC11h) runs (30-50%):
2957  ptaxis = (x<2.9)*(p30_50G[0]*exp(-pow(p30_50G[1]/x,p30_50G[2])) + p30_50G[3]*x) + (x>=2.9)*(p30_50G[4] + p30_50G[5]*x + p30_50G[6]*x*x);
2958  etaaxis = (y<0.0)*(p30_50G[7]*exp(-pow(p30_50G[8]/TMath::Abs(y+0.91),p30_50G[9])) + p30_50G[10]*y) + (y>=0.0 && y<=0.4)*(p30_50G[11] + p30_50G[12]*y + p30_50G[13]*y*y) + (y>0.4)*(p30_50G[14]*exp(-pow(p30_50G[15]/TMath::Abs(-y+0.91),p30_50G[16])));
2959  TRefficiency = ptaxis*etaaxis;
2960  break;
2961 
2962  case 9 :
2963  // Parameter values for GOOD TPC (LHC11h) runs (50-90%):
2964  ptaxis = (x<2.9)*(p50_90G[0]*exp(-pow(p50_90G[1]/x,p50_90G[2])) + p50_90G[3]*x) + (x>=2.9)*(p50_90G[4] + p50_90G[5]*x + p50_90G[6]*x*x);
2965  etaaxis = (y<0.0)*(p50_90G[7]*exp(-pow(p50_90G[8]/TMath::Abs(y+0.91),p50_90G[9])) + p50_90G[10]*y) + (y>=0.0 && y<=0.4)*(p50_90G[11] + p50_90G[12]*y + p50_90G[13]*y*y) + (y>0.4)*(p50_90G[14]*exp(-pow(p50_90G[15]/TMath::Abs(-y+0.91),p50_90G[16])));
2966  TRefficiency = ptaxis*etaaxis;
2967  break;
2968 
2969  default :
2970  // no Efficiency Switch option selected.. therefore don't correct, and set eff = 1
2971  TRefficiency = 1.0;
2972  }
2973 
2974  //if(TRefficiency < 0.1) cout<<"pt: "<<x<<" eta: "<<y<<" cent: "<<fCent<<" Eff: "<<TRefficiency<<" Ptaxis: "<<ptaxis<<" etaaxis: "<<etaaxis<<endl;
2975  //if(TRefficiency > 0.90) cout<<"pt: "<<x<<" eta: "<<y<<" cent: "<<fCent<<" Eff: "<<TRefficiency<<" Ptaxis: "<<ptaxis<<" etaaxis: "<<etaaxis<<endl;
2976 
2977  return TRefficiency;
2978 }
2979 
2980 //_____________________________________________________________________________
2981 void AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneVZERO(Double_t vzero[2][2]) const
2982 {
2983  // by default use the ep from the event header (make sure EP selection task is enabled!)
2984  Double_t a(0), b(0), c(0), d(0), e(0), f(0), g(0), h(0);
2985  vzero[0][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 2, a, b);
2986  vzero[1][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 2, c, d);
2987  vzero[0][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 3, e, f);
2988  vzero[1][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 3, g, h);
2989 }
2990 //_____________________________________________________________________________
2991 void AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneCombinedVZERO(Double_t* comb) const
2992 {
2993  // define some placeholders
2994  Double_t Q2[] = {-999., -999.};
2995  Double_t Q3[] = {-999., -999.};
2996 
2997  // for all other types use calibrated event plane from the event header
2998  //
2999  // note that the code is a bit messy here, for 11h done all in this function
3000  // reason is that the procedure is much shorter as the calibration is done in another task
3001  //
3002  // define some pleaceholders to the values by reference
3003  Double_t qx2a(0), qy2a(0), qx2c(0), qy2c(0), qx3a(0), qy3a(0), qx3c(0), qy3c(0);
3004  // get the q-vectors by reference
3005  InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 2, qx2a, qy2a); // 2nd order VZERO Aside
3006  InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 2, qx2c, qy2c); // 2nd order VZERO Cside
3007  InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 3, qx3a, qy3a); // 3rd order VZERO Aside
3008  InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 3, qx3c, qy3c); // 3rd order VZERO Cside
3009 
3010  // get cache index and retrieve the chi weights for this centrality
3011 // Int_t VZEROcentralityBin(GetCentBin(fCent)); // won't work because I use different bins then calibration value bins
3012  Int_t VZEROcentralityBin(GetVZEROCentralityBin());
3013  Double_t chi2A(fChi2A->At(VZEROcentralityBin));
3014  Double_t chi2C(fChi2C->At(VZEROcentralityBin));
3015  Double_t chi3A(fChi3A->At(VZEROcentralityBin));
3016  Double_t chi3C(fChi3C->At(VZEROcentralityBin));
3017 
3018  // combine the vzeroa and vzeroc signal
3019  Q2[0] = chi2A*chi2A*qx2a+chi2C*chi2C*qx2c;
3020  Q2[1] = chi2A*chi2A*qy2a+chi2C*chi2C*qy2c;
3021  Q3[0] = chi3A*chi3A*qx3a+chi3C*chi3C*qx3c;
3022  Q3[1] = chi3A*chi3A*qy3a+chi3C*chi3C*qy3c;
3023 
3024  comb[0] = .5*TMath::ATan2(Q2[1], Q2[0]);
3025  comb[1] = (1./3.)*TMath::ATan2(Q3[1], Q3[0]);
3026 }
3027 
3028 //_____________________________________________________________________________
3029 void AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneTPC(Double_t* tpc)
3030 { // made change on Apr20, 2016 for updated framework
3031 
3032  // grab the TPC event plane
3033  fNAcceptedTracks = 0; // reset the track counter
3034  Double_t qx2(0), qy2(0); // for psi2
3035  Double_t qx3(0), qy3(0); // for psi3
3036 
3037  // first check for track object
3038  //if(tracks) {
3039  if(fTracksFromContainer) {
3040  Float_t excludeInEta = -999;
3041  if(fExcludeLeadingJetsFromFit > 0 ) { // remove the leading jet from ep estimate
3042  if(fLeadingJet) excludeInEta = fLeadingJet->Eta();
3043  }
3044 
3045  // loop over tracks in TPC
3046  //for(Int_t iTPC(0); iTPC < tracks->GetEntries(); iTPC++) {
3047  for(Int_t iTPC(0); iTPC < fTracksFromContainer->GetEntries(); iTPC++) {
3048  //AliVParticle *track = dynamic_cast<AliVParticle*>(tracks->At(iTPC));
3049  AliVParticle *track = dynamic_cast<AliVParticle*>(fTracksFromContainer->At(iTPC));
3050 
3051  // apply general track cuts
3052  if(TMath::Abs(track->Eta())>0.9) continue;
3053  if(track->Pt()<0.15) continue;
3054  if(track->Pt() < fSoftTrackMinPt_ep || track->Pt() > fSoftTrackMaxPt_ep) continue; // APPLY soft cut
3055  if(fExcludeLeadingJetsFromFit > 0 &&( (TMath::Abs(track->Eta() - excludeInEta) < fJetRad*fExcludeLeadingJetsFromFit ) || (TMath::Abs(track->Eta()) - fJetRad - fEtamax ) > 0 )) continue;
3056  fNAcceptedTracks++;
3057 
3058  // sum up q-vectors
3059  qx2+= TMath::Cos(2.*track->Phi());
3060  qy2+= TMath::Sin(2.*track->Phi());
3061  qx3+= TMath::Cos(3.*track->Phi());
3062  qy3+= TMath::Sin(3.*track->Phi());
3063  }
3064  }
3065  tpc[0] = .5*TMath::ATan2(qy2, qx2);
3066  tpc[1] = (1./3.)*TMath::ATan2(qy3, qx3);
3067 }
3068 
3069 //_____________________________________________________________________________
3070 void AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneResolution(Double_t vzero[2][2], Double_t* vzeroComb, Double_t* tpc)
3071 { // updated here Apr20, 2016 for updates to framework
3072  // fill the profiles for the resolution parameters
3073  Int_t centbin = GetCentBin(fCent);
3074  fInCentralitySelection = centbin;
3075 
3076  // R2 resolution for 2nd order event plane
3077  fProfV2Resolution[fInCentralitySelection]->Fill(2., TMath::Cos(2.*(vzero[0][0] - vzero[1][0])));
3078  fProfV2Resolution[fInCentralitySelection]->Fill(3., TMath::Cos(2.*(vzero[1][0] - vzero[0][0])));
3079  fProfV2Resolution[fInCentralitySelection]->Fill(4., TMath::Cos(2.*(vzero[0][0] - tpc[0])));
3080  fProfV2Resolution[fInCentralitySelection]->Fill(5., TMath::Cos(2.*(tpc[0] - vzero[0][0])));
3081  fProfV2Resolution[fInCentralitySelection]->Fill(6., TMath::Cos(2.*(vzero[1][0] - tpc[0])));
3082  fProfV2Resolution[fInCentralitySelection]->Fill(7., TMath::Cos(2.*(tpc[0] - vzero[1][0])));
3083  // R3 resolution for 2nd order event plane
3084  fProfV3Resolution[fInCentralitySelection]->Fill(2., TMath::Cos(3.*(vzero[0][0] - vzero[1][0])));
3085  fProfV3Resolution[fInCentralitySelection]->Fill(3., TMath::Cos(3.*(vzero[1][0] - vzero[0][0])));
3086  fProfV3Resolution[fInCentralitySelection]->Fill(4., TMath::Cos(3.*(vzero[0][0] - tpc[0])));
3087  fProfV3Resolution[fInCentralitySelection]->Fill(5., TMath::Cos(3.*(tpc[0] - vzero[0][0])));
3088  fProfV3Resolution[fInCentralitySelection]->Fill(6., TMath::Cos(3.*(vzero[1][0] - tpc[0])));
3089  fProfV3Resolution[fInCentralitySelection]->Fill(7., TMath::Cos(3.*(tpc[0] - vzero[1][0])));
3090  // R4 resolution for 2nd order event plane
3091  fProfV4Resolution[fInCentralitySelection]->Fill(2., TMath::Cos(4.*(vzero[0][0] - vzero[1][0])));
3092  fProfV4Resolution[fInCentralitySelection]->Fill(3., TMath::Cos(4.*(vzero[1][0] - vzero[0][0])));
3093  fProfV4Resolution[fInCentralitySelection]->Fill(4., TMath::Cos(4.*(vzero[0][0] - tpc[0])));
3094  fProfV4Resolution[fInCentralitySelection]->Fill(5., TMath::Cos(4.*(tpc[0] - vzero[0][0])));
3095  fProfV4Resolution[fInCentralitySelection]->Fill(6., TMath::Cos(4.*(vzero[1][0] - tpc[0])));
3096  fProfV4Resolution[fInCentralitySelection]->Fill(7., TMath::Cos(4.*(tpc[0] - vzero[1][0])));
3097  // R5 resolution for 2nd order event plane
3098  fProfV5Resolution[fInCentralitySelection]->Fill(2., TMath::Cos(5.*(vzero[0][0] - vzero[1][0])));
3099  fProfV5Resolution[fInCentralitySelection]->Fill(3., TMath::Cos(5.*(vzero[1][0] - vzero[0][0])));
3100  fProfV5Resolution[fInCentralitySelection]->Fill(4., TMath::Cos(5.*(vzero[0][0] - tpc[0])));
3101  fProfV5Resolution[fInCentralitySelection]->Fill(5., TMath::Cos(5.*(tpc[0] - vzero[0][0])));
3102  fProfV5Resolution[fInCentralitySelection]->Fill(6., TMath::Cos(5.*(vzero[1][0] - tpc[0])));
3103  fProfV5Resolution[fInCentralitySelection]->Fill(7., TMath::Cos(5.*(tpc[0] - vzero[1][0])));
3104 
3105  // for the resolution of the combined vzero event plane, use two tpc halves as uncorrelated subdetectors
3106  Double_t qx2a(0), qy2a(0); // for psi2a, negative eta
3107  Double_t qx3a(0), qy3a(0); // for psi3a, negative eta
3108  Double_t qx2b(0), qy2b(0); // for psi2c, positive eta
3109  Double_t qx3b(0), qy3b(0); // for psi3c, positive eta
3110  // not needed
3111  Double_t qx4a(0), qy4a(0); // for psi4a, negative eta
3112  Double_t qx5a(0), qy5a(0); // for psi4a, negative eta
3113  Double_t qx4b(0), qy4b(0); // for psi5c, positive eta
3114  Double_t qx5b(0), qy5b(0); // for psi5c, positive eta
3115 
3116  // check for track object and loop over tpc tracks
3117  //if(tracks) {
3118  if(fTracksFromContainer) {
3119  //Int_t iTracks(tracks->GetEntriesFast());
3120  Int_t iTracks(fTracksFromContainer->GetEntriesFast());
3121  for(Int_t iTPC(0); iTPC < iTracks; iTPC++) {
3122  //AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTPC));
3123  AliVTrack* track = static_cast<AliVTrack*>(fTracksFromContainer->At(iTPC));
3124 
3125  // apply track cuts
3126  if(TMath::Abs(track->Eta())>0.9) continue;
3127  if(track->Pt()<0.15) continue;
3128  if(track->Pt() < fSoftTrackMinPt_ep || track->Pt() > fSoftTrackMaxPt_ep) continue;
3129  if(track->Eta() < 0 ) {
3130  // negative Eta q-vect calculation
3131  qx2a+= TMath::Cos(2.*track->Phi());
3132  qy2a+= TMath::Sin(2.*track->Phi());
3133  qx3a+= TMath::Cos(3.*track->Phi());
3134  qy3a+= TMath::Sin(3.*track->Phi());
3135  qx4a+= TMath::Cos(4.*track->Phi());
3136  qy4a+= TMath::Sin(4.*track->Phi());
3137  qx5a+= TMath::Cos(5.*track->Phi());
3138  qy5a+= TMath::Sin(5.*track->Phi());
3139  } else if (track->Eta() > 0) {
3140  // positive Eta q-vect calculation
3141  qx2b+= TMath::Cos(2.*track->Phi());
3142  qy2b+= TMath::Sin(2.*track->Phi());
3143  qx3b+= TMath::Cos(3.*track->Phi());
3144  qy3b+= TMath::Sin(3.*track->Phi());
3145  qx4b+= TMath::Cos(4.*track->Phi());
3146  qy4b+= TMath::Sin(4.*track->Phi());
3147  qx5b+= TMath::Cos(5.*track->Phi());
3148  qy5b+= TMath::Sin(5.*track->Phi());
3149  }
3150  }
3151  }
3152 
3153  // combine x-y vectors for neg and pos Eta ranges
3154  Double_t tpca2(.5*TMath::ATan2(qy2a, qx2a));
3155  Double_t tpcb2(.5*TMath::ATan2(qy2b, qx2b));
3156 
3157  // not used, but can be for 3rd, 4th, and 5th order event planes
3158 /*
3159  Double_t tpca3((1./3.)*TMath::ATan2(qy3a, qx3a));
3160  Double_t tpcb3((1./3.)*TMath::ATan2(qy3b, qx3b));
3161  Double_t tpca4(.25*TMath::ATan2(qy4a, qx4a));
3162  Double_t tpcb4(.25*TMath::ATan2(qy4b, qx4b));
3163  Double_t tpca5((.2)*TMath::ATan2(qy5a, qx5a));
3164  Double_t tpcb5((.2)*TMath::ATan2(qy5b, qx5b));
3165 */
3166 
3167  // R2 resolution for 2nd order event plane
3168  fProfV2Resolution[fInCentralitySelection]->Fill(8., TMath::Cos(2.*(vzeroComb[0] - tpca2)));
3169  fProfV2Resolution[fInCentralitySelection]->Fill(9., TMath::Cos(2.*(vzeroComb[0] - tpcb2)));
3170  fProfV2Resolution[fInCentralitySelection]->Fill(10., TMath::Cos(2.*(tpca2 - tpcb2)));
3171  // R3 resolution for 2nd order event plane
3172  fProfV3Resolution[fInCentralitySelection]->Fill(8., TMath::Cos(3.*(vzeroComb[1] - tpca2)));
3173  fProfV3Resolution[fInCentralitySelection]->Fill(9., TMath::Cos(3.*(vzeroComb[1] - tpcb2)));
3174  fProfV3Resolution[fInCentralitySelection]->Fill(10., TMath::Cos(3.*(tpca2 - tpcb2)));
3175  // R4 resolution for 2nd order event plane
3176  fProfV4Resolution[fInCentralitySelection]->Fill(8., TMath::Cos(4.*(vzeroComb[0] - tpca2)));
3177  fProfV4Resolution[fInCentralitySelection]->Fill(9., TMath::Cos(4.*(vzeroComb[0] - tpcb2)));
3178  fProfV4Resolution[fInCentralitySelection]->Fill(10., TMath::Cos(4.*(tpca2 - tpcb2)));
3179  // R5 resolution for 2nd order event plane
3180  fProfV5Resolution[fInCentralitySelection]->Fill(8., TMath::Cos(5.*(vzeroComb[1] - tpca2)));
3181  fProfV5Resolution[fInCentralitySelection]->Fill(9., TMath::Cos(5.*(vzeroComb[1] - tpcb2)));
3182  fProfV5Resolution[fInCentralitySelection]->Fill(10., TMath::Cos(5.*(tpca2 - tpcb2)));
3183 }
3184 
3185 //_____________________________________________________________________________i
3186 void AliAnalysisTaskEmcalJetHadEPpid::ReadVZEROCalibration2011h()
3187 {
3188  // make sure calibration parameters are present for vzero ep 11h data
3189  //
3190  // chi values can be calculated using the static helper function
3191  // AliAnalysisTaskJetV2::CalculateEventPlaneChi(Double_t res) where res is the event plane
3192  // resolution in a given centrality bin
3193  //
3194  // the resolutions that were used for these defaults are
3195  // this might need a bit of updating as they were read 'by-eye' from a performance plot ..
3196  // Double_t R2VZEROA[] = {.35, .40, .48, .50, .48, .45, .38, .26, .16};
3197  // Double_t R2VZEROC[] = {.45, .60, .70, .73, .68, .60, .40, .36, .17};
3198  // Double_t R3VZEROA[] = {.22, .23, .22, .19, .15, .12, .08, .00, .00};
3199  // Double_t R3VZEROC[] = {.30, .30, .28, .25, .22, .17, .11, .00, .00};
3200 
3201  Double_t chiA2[] = {0.582214, 0.674622, 0.832214, 0.873962, 0.832214, 0.771423, 0.637146, 0.424255, 0.257385};
3202  Double_t chiC2[] = {0.771423, 1.10236, 1.38116, 1.48077, 1.31964, 1.10236, 0.674622, 0.600403, 0.273865};
3203  Double_t chiA3[] = {0.356628, 0.373474, 0.356628, 0.306702, 0.24115, 0.192322, 0.127869, 6.10352e-05, 6.10352e-05};
3204  Double_t chiC3[] = {0.493347, 0.493347, 0.458557, 0.407166, 0.356628, 0.273865, 0.176208, 6.10352e-05, 6.10352e-05};
3205 
3206  // when the user wants to, set the weights to 1 (effectively disabling them)
3207  if(!fUseChiWeightForVZERO) {
3208  for(Int_t i(0); i < 9; i++) {
3209  chiA2[i] = 1.;
3210  chiA3[i] = 1.;
3211  chiC2[i] = 1.;
3212  chiC3[i] = 1.;
3213  }
3214  }
3215 
3216  if(!fChi2A) fChi2A = new TArrayD(9, chiA2);
3217  if(!fChi2C) fChi2C = new TArrayD(9, chiC2);
3218  if(!fChi3A) fChi3A = new TArrayD(9, chiA3);
3219  if(!fChi3C) fChi3C = new TArrayD(9, chiC3);
3220 }
3221 
3222 //_____________________________________________________________________________
3223 Int_t AliAnalysisTaskEmcalJetHadEPpid::GetVZEROCentralityBin() const
3224 {
3225  // return cache index number corresponding to the event centrality
3226  Float_t v0Centr(InputEvent()->GetCentrality()->GetCentralityPercentile("V0M"));
3227  if(v0Centr < 5) return 0;
3228  else if(v0Centr < 10) return 1;
3229  else if(v0Centr < 20) return 2;
3230  else if(v0Centr < 30) return 3;
3231  else if(v0Centr < 40) return 4;
3232  else if(v0Centr < 50) return 5;
3233  else if(v0Centr < 60) return 6;
3234  else if(v0Centr < 70) return 7;
3235  else return 8;
3236 }
3237 
3238 //_____________________________________________________________________________
3239 AliEmcalJet* AliAnalysisTaskEmcalJetHadEPpid::GetLeadingJet(AliLocalRhoParameter* localRho) {
3240  // return pointer to the highest pt jet (before background subtraction) within acceptance
3241  // only rudimentary cuts are applied on this level, hence the implementation outside of
3242  // the framework
3243  if(fJets) {
3244  Int_t iJets(fJets->GetEntriesFast());
3245  Double_t pt(0);
3246  AliEmcalJet* leadingJet(0x0);
3247  if(!localRho) {
3248  for(Int_t i(0); i < iJets; i++) {
3249  AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(i));
3250  if(!AcceptMyJet(jet)) continue;
3251  if(jet->Pt() > pt) {
3252  leadingJet = jet;
3253  pt = leadingJet->Pt();
3254  }
3255  }
3256  return leadingJet;
3257  } else {
3258  // return leading jet after background subtraction
3259  Double_t rho(0);
3260  for(Int_t i(0); i < iJets; i++) {
3261  AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(i));
3262  if(!AcceptMyJet(jet)) continue;
3263  //rho = localRho->GetLocalVal(jet->Phi(), GetJetContainer()->GetJetRadius(), localRho->GetVal());
3264  rho = localRho->GetLocalVal(jet->Phi(), fJetRad);
3265  if((jet->Pt()-jet->Area()*rho) > pt) {
3266  leadingJet = jet;
3267  pt = (leadingJet->Pt()-jet->Area()*rho);
3268  }
3269  }
3270  return leadingJet;
3271 
3272  }
3273  }
3274 
3275  return 0x0;
3276 }
3277 
3278 //_____________________________________________________________________________
3279 Double_t AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneChi(Double_t res)
3280 {
3281  // return chi for given resolution to combine event plane estimates from two subevents
3282  // see Phys. Rev. C no. CS6346 (http://arxiv.org/abs/nucl-ex/9805001)
3283  Double_t chi(2.), delta(1.), con((TMath::Sqrt(TMath::Pi()))/(2.*TMath::Sqrt(2)));
3284  for (Int_t i(0); i < 15; i++) {
3285  chi = ((con*chi*TMath::Exp(-chi*chi/4.)*(TMath::BesselI0(chi*chi/4.)+TMath::BesselI1(chi*chi/4.))) < res) ? chi + delta : chi - delta;
3286  delta = delta / 2.;
3287  }
3288  return chi;
3289 }
3290 
3291 TH1* AliAnalysisTaskEmcalJetHadEPpid::FillEventTriggerQA(TH1* h, UInt_t trig) {
3292  // check and fill a Event Selection QA histogram for different trigger selections after cuts
3293  if(trig == 0) h->Fill(1);
3294  if(trig & AliVEvent::kAny) h->Fill(2);
3295  if(trig & AliVEvent::kAnyINT) h->Fill(3);
3296  if(trig & AliVEvent::kMB) h->Fill(4);
3297  if(trig & AliVEvent::kINT7) h->Fill(5);
3298  if(trig & AliVEvent::kEMC1) h->Fill(6);
3299  if(trig & AliVEvent::kEMC7) h->Fill(7);
3300  if(trig & AliVEvent::kEMC8) h->Fill(8);
3301  if(trig & AliVEvent::kEMCEJE) h->Fill(9);
3302  if(trig & AliVEvent::kEMCEGA) h->Fill(10);
3303  if(trig & AliVEvent::kCentral) h->Fill(11);
3304  if(trig & AliVEvent::kSemiCentral) h->Fill(12);
3305  if(trig & AliVEvent::kINT8) h->Fill(13);
3306 
3307  if(trig & (AliVEvent::kEMCEJE | AliVEvent::kMB)) h->Fill(14);
3308  if(trig & (AliVEvent::kEMCEGA | AliVEvent::kMB)) h->Fill(15);
3309  if(trig & (AliVEvent::kAnyINT | AliVEvent::kMB)) h->Fill(16);
3310 
3311  if(trig & (AliVEvent::kEMCEJE & (AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral))) h->Fill(17);
3312  if(trig & (AliVEvent::kEMCEGA & (AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral))) h->Fill(18);
3313  if(trig & (AliVEvent::kAnyINT & (AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral))) h->Fill(19);
3314 
3315  // label bins of the analysis trigger selection summary
3316  h->GetXaxis()->SetBinLabel(1, "no trigger");
3317  h->GetXaxis()->SetBinLabel(2, "kAny");
3318  h->GetXaxis()->SetBinLabel(3, "kAnyINT");
3319  h->GetXaxis()->SetBinLabel(4, "kMB");
3320  h->GetXaxis()->SetBinLabel(5, "kINT7");
3321  h->GetXaxis()->SetBinLabel(6, "kEMC1");
3322  h->GetXaxis()->SetBinLabel(7, "kEMC7");
3323  h->GetXaxis()->SetBinLabel(8, "kEMC8");
3324  h->GetXaxis()->SetBinLabel(9, "kEMCEJE");
3325  h->GetXaxis()->SetBinLabel(10, "kEMCEGA");
3326  h->GetXaxis()->SetBinLabel(11, "kCentral");
3327  h->GetXaxis()->SetBinLabel(12, "kSemiCentral");
3328  h->GetXaxis()->SetBinLabel(13, "kINT8");
3329  h->GetXaxis()->SetBinLabel(14, "kEMCEJE or kMB");
3330  h->GetXaxis()->SetBinLabel(15, "kEMCEGA or kMB");
3331  h->GetXaxis()->SetBinLabel(16, "kAnyINT or kMB");
3332  h->GetXaxis()->SetBinLabel(17, "kEMCEJE & (kMB or kCentral or kSemiCentral)");
3333  h->GetXaxis()->SetBinLabel(18, "kEMCEGA & (kMB or kCentral or kSemiCentral)");
3334  h->GetXaxis()->SetBinLabel(19, "kAnyINT & (kMB or kCentral or kSemiCentral)");
3335 
3336  // set x-axis labels vertically
3337  h->LabelsOption("v");
3338  //h->LabelsDeflate("X");
3339 
3340  return h;
3341 }
AliAnalysisTaskEmcalJetHadEPpid::CloneAndReduceTrackList
TObjArray * CloneAndReduceTrackList(TObjArray *tracks)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2297
charge
Int_t charge
Definition: MuonTrackingEfficiency.C:54
AliEmcalJet::Area
Double_t Area() const
Definition: AliEmcalJet.h:97
AliAnalysisTaskEmcalJetHadEPpid::GetpTtrackBin
virtual Int_t GetpTtrackBin(Double_t pt) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2121
AliAnalysisTaskEmcalJetHadEPpid::RelativeEPJET
Double_t RelativeEPJET(Double_t jetAng, Double_t EPAng) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2067
AliAnalysisTaskEmcalJetHadEPpid::fProfV5Resolution
TProfile * fProfV5Resolution[10]
resolution parameters for v4
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:458
AliAnalysisTaskEmcalJetHadEPpid::GetLeadingJet
AliEmcalJet * GetLeadingJet(AliLocalRhoParameter *localRho=0x0)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:3239
AliPicoTrack.h
AliAnalysisTaskEmcalJetHadEPpid::fhnJH
THnSparse * fhnJH
// mixed events matrix
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:451
AliAnalysisTaskEmcalJetHadEPpid::fLeadingJet
AliEmcalJet * fLeadingJet
number of accepted tracks
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:262
ClassImp
ClassImp(AliAnalysisTaskEmcalJetHadEPpid) AliAnalysisTaskEmcalJetHadEPpid
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:95
AliEmcalJet::Eta
Double_t Eta() const
Definition: AliEmcalJet.h:88
AliTrackContainer.h
AliEmcalJet::GetFlavour
Int_t GetFlavour() const
Definition: AliEmcalJet.h:216
AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneChi
static Double_t CalculateEventPlaneChi(Double_t res)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:3279
AliAnalysisTaskEmcalJetHadEPpid::fChi2A
TArrayD * fChi2A
resolution parameters for v5
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:459
AliEmcalJet::Phi
Double_t Phi() const
Definition: AliEmcalJet.h:84
AliParticleContainer::GetNParticles
Int_t GetNParticles() const
Definition: AliParticleContainer.h:69
AliTrackContainer
Container with name, TClonesArray and cuts for particles.
Definition: AliTrackContainer.h:29
AliLocalRhoParameter::GetLocalVal
Double_t GetLocalVal(Double_t phi, Double_t r, Double_t n) const
Definition: AliLocalRhoParameter.h:16
AliAnalysisTaskEmcal::fEPV0
Double_t fEPV0
!event plane V0
Definition: AliAnalysisTaskEmcal.h:301
list
TList * list
Definition: DrawAnaCaloTrackQA.C:18
AliAnalysisTaskEmcalJetHadEPpid::fHistTrackPtallcent
TH1 * fHistTrackPtallcent
// phi distrubtion of mixed events
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:415
AliEmcalParticle.h
AliAnalysisTaskEmcalJetHadEPpid::fJetsCont
AliJetContainer * fJetsCont
tracks from AliTrackContainer
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:469
AliAnalysisTaskEmcalJet::GetLeadingHadronPt
Double_t GetLeadingHadronPt(AliEmcalJet *jet, Int_t c=0)
Definition: AliAnalysisTaskEmcalJet.cxx:617
AliAnalysisTaskEmcalJetHadEPpid::GetEtaBin
virtual Int_t GetEtaBin(Double_t eta) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2095
AliAnalysisTaskEmcalJetHadEPpid::fProfV4Resolution
TProfile * fProfV4Resolution[10]
resolution parameters for v3
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:457
AliAnalysisTaskEmcalJet::GetLocalRhoFromEvent
AliLocalRhoParameter * GetLocalRhoFromEvent(const char *name)
Definition: AliAnalysisTaskEmcalJet.cxx:109
AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneResolution
void CalculateEventPlaneResolution(Double_t vzero[2][2], Double_t *vzeroComb, Double_t *tpc)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:3070
AliAnalysisTaskEmcalJetHadEPpid::GetDimParams
virtual void GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2194
AliAnalysisTaskEmcalJetHadEPpid::GetpTjetBin
virtual Int_t GetpTjetBin(Double_t pt) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2107
AliAnalysisTaskEmcal::fEPV0C
Double_t fEPV0C
!event plane V0C
Definition: AliAnalysisTaskEmcal.h:303
AliEmcalJet::GetNumberOfConstituents
UShort_t GetNumberOfConstituents() const
Definition: AliEmcalJet.h:107
AliEmcalJet::GetNumberOfTracks
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:106
AliAnalysisTaskEmcalJetHadEPpid::GetDimParamsPID
virtual void GetDimParamsPID(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2379
AliAnalysisTaskEmcalJet::fLocalRhoName
TString fLocalRhoName
name for local rho
Definition: AliAnalysisTaskEmcalJet.h:96
AliPicoTrack::Eta
Double_t Eta() const
Definition: AliPicoTrack.h:37
AliAnalysisTaskEmcalJet::GetRhoFromEvent
AliRhoParameter * GetRhoFromEvent(const char *name)
Definition: AliAnalysisTaskEmcalJet.cxx:89
AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseFCorr
virtual THnSparse * NewTHnSparseFCorr(const char *name, UInt_t entries)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2682
AliEmcalJet::GetNumberOfClusters
UShort_t GetNumberOfClusters() const
Definition: AliEmcalJet.h:105
AliAnalysisTaskEmcalJetHadEPpid::fHistJetHaddPHI
TH1F * fHistJetHaddPHI
// mixed events phi-eta distributions
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:435
AliAnalysisTaskEmcalJetHadEPpid::fProfV3Resolution
TProfile * fProfV3Resolution[10]
resolution parameters for v2
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:456
AliAnalysisTaskEmcalJetHadEPpid::GetDimParamsCorr
virtual void GetDimParamsCorr(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2715
AliAnalysisTaskEmcalJetHadEPpid::fHistJetPtvsTrackPt
TH2F * fHistJetPtvsTrackPt[6]
number of jets versus Centrality
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:365
ptmax
const Double_t ptmax
Definition: AddTaskCFDStar.C:39
AliAnalysisTaskEmcalJet::fRho
AliRhoParameter * fRho
! event rho
Definition: AliAnalysisTaskEmcalJet.h:100
AliParticleContainer.h
AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseF
virtual THnSparse * NewTHnSparseF(const char *name, UInt_t entries)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2161
AliJetContainer.h
AliEmcalJet::MaxTrackPt
Double_t MaxTrackPt() const
Definition: AliEmcalJet.h:122
AliPicoTrack::Phi
Double_t Phi() const
Definition: AliPicoTrack.h:33
AliAnalysisTaskEmcalJetHadEPpid::fPIDResponse
AliPIDResponse * fPIDResponse
// event pool Manager object
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:341
AliAnalysisTaskEmcal::fCent
Double_t fCent
!event centrality
Definition: AliAnalysisTaskEmcal.h:299
AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseFPID
virtual THnSparse * NewTHnSparseFPID(const char *name, UInt_t entries)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2345
AliAnalysisTaskEmcalJet::fLocalRho
AliLocalRhoParameter * fLocalRho
! local event rho
Definition: AliAnalysisTaskEmcalJet.h:101
AliRhoParameter.h
AliPicoTrack::Pt
Double_t Pt() const
Definition: AliPicoTrack.h:23
AliAnalysisTaskEmcalJetHadEPpid::fCaloClustersCont
AliClusterContainer * fCaloClustersCont
Tracks - Need this for applying track quality cuts.
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:472
AliAnalysisTaskEmcalJetHadEPpid::RelativePhi
Double_t RelativePhi(Double_t mphi, Double_t vphi) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2051
AliAnalysisTaskEmcalJetHadEPpid::GetzVertexBin
virtual Int_t GetzVertexBin(Double_t zVtx) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2141
AliAnalysisTaskEmcalJet::fJets
TClonesArray * fJets
! jets
Definition: AliAnalysisTaskEmcalJet.h:99
AliEmcalJet::Pt
Double_t Pt() const
Definition: AliEmcalJet.h:76
AliLocalRhoParameter.h
AliAnalysisTaskEmcalJetHadEPpid::EffCorrection
Double_t EffCorrection(Double_t trkETA, Double_t trkPT, Int_t effswitch) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2841
AliAnalysisTaskEmcal::fOutput
AliEmcalList * fOutput
!output list
Definition: AliAnalysisTaskEmcal.h:317
AliAnalysisTaskEmcalJetHadEPpid::AcceptFlavourJet
virtual Int_t AcceptFlavourJet(AliEmcalJet *jet, Int_t NUM)
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2796
AliAnalysisTaskEmcal::fTracks
TClonesArray * fTracks
!tracks
Definition: AliAnalysisTaskEmcal.h:294
AliAnalysisTaskEmcalJetHadEPpid::fProfV2Resolution
TProfile * fProfV2Resolution[10]
// sparse to get # jet triggers
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:455
AliAnalysisTaskEmcal::GetTrackContainer
AliTrackContainer * GetTrackContainer(Int_t i=0) const
Definition: AliAnalysisTaskEmcal.h:132
AliAnalysisTaskEmcal::fCreateHisto
Bool_t fCreateHisto
whether or not create histograms
Definition: AliAnalysisTaskEmcal.h:246
AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneCombinedVZERO
void CalculateEventPlaneCombinedVZERO(Double_t *comb) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2991
AliAnalysisTaskEmcal::SetMakeGeneralHistograms
void SetMakeGeneralHistograms(Bool_t g)
Definition: AliAnalysisTaskEmcal.h:148
AliAnalysisTaskEmcalJetHadEPpid.h
AliAnalysisTaskEmcalJet
Base task in the EMCAL jet framework.
Definition: AliAnalysisTaskEmcalJet.h:30
AliAnalysisTaskEmcal::fEPV0A
Double_t fEPV0A
!event plane V0A
Definition: AliAnalysisTaskEmcal.h:302
AliEmcalJet::MaxClusterPt
Double_t MaxClusterPt() const
Definition: AliEmcalJet.h:121
AliEmcalJet
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:44
pi
const Double_t pi
Definition: AliAnalysisTaskDijetHadron.cxx:72
nbins
const Int_t nbins
Definition: AverageDmesonRaa.C:63
AliTrackContainer::GetNextAcceptParticle
virtual AliVParticle * GetNextAcceptParticle()
Definition: AliTrackContainer.h:57
AliAnalysisTaskEmcalJetHadEPpid::fhnCorr
THnSparse * fhnCorr
// jet hadron events matrix
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:452
AliAnalysisTaskEmcalJetHadEPpid::fHistMEphieta
TH2F * fHistMEphieta
// single events phi-eta distributions
Definition: AliAnalysisTaskEmcalJetHadEPpid.h:434
AliParticleContainer::GetNAcceptedParticles
Int_t GetNAcceptedParticles() const
Definition: AliParticleContainer.cxx:414
AliAnalysisTaskEmcalJet::fRhoVal
Double_t fRhoVal
! event rho value, same for local rho
Definition: AliAnalysisTaskEmcalJet.h:102
AliAnalysisTaskEmcalJetHadEPpid::GetCentBin
virtual Int_t GetCentBin(Double_t cent) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2037
AliEmcalJet.h
AliAnalysisTaskEmcalJetHadEPpid::CalculateEventPlaneVZERO
void CalculateEventPlaneVZERO(Double_t vzero[2][2]) const
Definition: AliAnalysisTaskEmcalJetHadEPpid.cxx:2981
AliClusterContainer.h
AliPicoTrack::Charge
Short_t Charge() const
Definition: AliPicoTrack.h:39
AliAnalysisTaskLocalRho.h