AliPhysics  720d1f3 (720d1f3)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskHJetSpectra.cxx
Go to the documentation of this file.
1 #ifndef ALIANALYSISTASKSE_H
2 
3 #include <Riostream.h>
4 #include <TROOT.h>
5 #include <TFile.h>
6 #include <TChain.h>
7 #include <TTree.h>
8 #include <TKey.h>
9 #include <TProfile.h>
10 #include <TProfile2D.h>
11 #include <TH1.h>
12 #include <TH1F.h>
13 #include <TF1.h>
14 #include <TH2F.h>
15 #include <TH1D.h>
16 #include <TH1I.h>
17 #include <TArrayF.h>
18 #include <TArrayD.h>
19 #include <THnSparse.h>
20 #include <TCanvas.h>
21 #include <TList.h>
22 #include <TClonesArray.h>
23 #include <TObject.h>
24 #include <TMath.h>
25 #include <TSystem.h>
26 #include <TInterpreter.h>
27 #include "AliAnalysisTask.h"
28 #include "AliCentrality.h"
29 #include "AliStack.h"
30 #include "AliESDEvent.h"
31 #include "AliESDInputHandler.h"
32 #include "AliAODEvent.h"
33 #include "AliAODHandler.h"
34 #include "AliAnalysisManager.h"
35 #include "AliAnalysisTaskSE.h"
36 #include "AliAnalysisHelperJetTasks.h"
37 #include "AliParticleContainer.h"
38 #include "AliInputEventHandler.h"
39 #endif
40 
41 #include <time.h>
42 #include <TRandom3.h>
43 #include "AliGenEventHeader.h"
44 #include "AliGenPythiaEventHeader.h"
45 #include "AliGenHijingEventHeader.h"
46 #include "AliAODMCHeader.h"
47 #include "AliMCEvent.h"
48 #include "AliLog.h"
49 #include <AliEmcalJet.h>
50 #include <AliPicoTrack.h>
51 #include "AliVEventHandler.h"
52 #include "AliVParticle.h"
53 #include "AliAODMCParticle.h"
54 #include "AliAnalysisUtils.h"
55 #include "AliRhoParameter.h"
56 #include "TVector3.h"
57 #include "AliVVertex.h"
58 #include "AliExternalTrackParam.h"
59 
60 #include <stdio.h>
61 #include <stdlib.h>
62 
63 #include "AliGenDPMjetEventHeader.h"
64 #include "AliJetContainer.h"
65 #include "AliAnalysisTaskEmcal.h"
66 #include "AliAnalysisTaskEmcalJet.h"
67 #include "AliAnalysisTaskHJetSpectra.h"
68 #include "AliHeader.h"
69 #include "AliRunLoader.h"
70 #include "AliVVZERO.h"
71 #include "AliAODZDC.h"
72 #include "AliVZDC.h"
73 
74 using namespace std;
75 
76 // ANALYSIS OF HIGH PT HADRON TRIGGER ASSOCIATED SPECTRUM OF RECOIL JETS IN P+PB
77 // Author Filip Krizek (7.Oct. 2015)
78 
79 ClassImp(AliAnalysisTaskHJetSpectra)
80 //________________________________________________________________________________________
81 
82 AliAnalysisTaskHJetSpectra::AliAnalysisTaskHJetSpectra():
83 AliAnalysisTaskEmcalJet("AliAnalysisTaskHJetSpectra", kTRUE),
84  fCollisionSystem(0), fTypeOfData(0), fTypeOfAnal(0),
85  fUseDefaultVertexCut(1), fUsePileUpCut(1),
86  fRhoTaskName(), fRhoTaskNameMC(),
87  fSignalJetRadius(0.4), fSignalJetRadiusSquared(fSignalJetRadius*fSignalJetRadius),
88 fSignalJetEtaWindow(0.9 - fSignalJetRadius), fTrackEtaWindow(0.9), fMinTrackPt(0.150), fMinJetArea(0.0),
89 fCentralityType("V0A"), fMinFractionShared(0.5),
90 fCrossSection(0.0), fTrials(0.0), fImpParam(-1.0), fRandom(0), fHelperClass(0), fInitializedLocal(0),
91 fTTType(0), fDphiCut(TMath::Pi()-0.6), fUseDoubleBinPrecision(0),
92 fHistEvtSelection(0x0),fhKTAreaPt(0x0),
93  fhVertexZ(0x0), fhVertexXAccept(0x0), fhVertexYAccept(0x0), fhVertexZAccept(0x0),
94  fhVertexXAcceptTT(0x0), fhVertexYAcceptTT(0x0), fhVertexZAcceptTT(0x0),
95 fhVertexZMC(0x0), fhVertexZAcceptMC(0x0),
96  fhDphiTriggerJetAccept(0x0), fhJetPhiIncl(0x0), fhJetEtaIncl(0x0),
97 fhCentralityV0M(0x0), fhCentralityV0A(0x0), fhCentralityV0C(0x0), fhCentralityZNA(0x0),
98 fhDiffPtVsPtTrackTrue(0x0),
99 fhTrackPhiCG(0x0), fhTrackPhiTPCG(0x0),
100 fhDCAinXVsPt(0x0),
101 fhDCAinYVsPt(0x0),
102 fhDCAinXVsPtStrange(0x0),
103 fhDCAinYVsPtStrange(0x0),
104 fhDCAinXVsPtNonStrange(0x0),
105 fhDCAinYVsPtNonStrange(0x0),
106 fCentralityBins(kCAll),
107 fTrackPtRef(0x0),
108 fTrackPtSig(0x0),
109 fNofRandomCones(1),
110 fZVertexCut(10.0),fCutPhi(0.6),
111 fpyVtx(3)
112 {
113  //default constructor
114  for(Int_t it=0; it<kTT;it++){
115  fRhoRec[it].Set(kRho);
116  fRhoMC[it].Set(kRho);
117  fhDphiTTTT[it]=NULL;
118 
119  fhCentralityTT[it] = NULL;
120  }
121 
122  for(Int_t iw=0; iw<21; iw++){
123  fhPtJetPrimVsPtJetRec[iw] = NULL;
124  }
125 
126  for(Int_t ic =0; ic<kCAll; ic++){
127  for(Int_t it=0; it<kTT;it++){
128  fh1Ntriggers[ic][it]=NULL;
129  fh1TriggerMult[ic][it]=NULL;
130  fh1NtriggersGen[ic][it]=NULL;
131  fh1TriggerMultGen[ic][it]=NULL;
132 
133 
134  fhImpactParameterTT[ic][it]=NULL;
135 
136  for(Int_t ir=0; ir<kRho; ir++){
137  fhDphiTriggerJet[ic][it][ir]=NULL;
138  fhDphiTriggerJetGen[ic][it][ir]=NULL;
139  fHJetSpec[ic][it][ir]=NULL;
140  fHJetSpecGen[ic][it][ir]=NULL;
141 
142  fhDeltaPt[ic][it][ir]=NULL;
143  fhDeltaPtEmb[ic][it][ir]=NULL;
144  fhDeltaPtEmb2D[ic][it][ir]=NULL;
145  fhDeltaPtEmbPerp[ic][it][ir]=NULL;
146  fhDeltaPtEmbPerp2D[ic][it][ir]=NULL;
147  fhDeltaPtEmbBc2Bc[ic][it][ir]=NULL;
148  fhDeltaPtEmbBc2Bc2D[ic][it][ir]=NULL;
149 
150  fhRhoTT[ic][it][ir]=NULL;
151  fARhoTT[ic][it][ir]=NULL;
152 
153  }
154  fhVzeroATotMultTT[ic][it]=NULL;
155  fhZNAEnergyTT[ic][it]=NULL;
156  fhTrackMultiplicityTT[ic][it]=NULL;
157  fhZNAVzeroATrackTT[ic][it]=NULL;
158  }
159  fhJetPhi[ic]=NULL;
160  fhJetPhiGen[ic]=NULL;
161  fhTrackPhi[ic]=NULL;
162  fhJetEta[ic]=NULL;
163  fhJetEtaGen[ic]=NULL;
164  fhTrackEta[ic]=NULL;
165  fhTrackPt[ic]=NULL;
166  fhTrackPtGen[ic]=NULL;
167  fhCentrality[ic]=NULL;
168  fhVzeroATotMult[ic]=NULL;
169  fhZNAEnergy[ic]=NULL;
170  fhTrackMultiplicity[ic]=NULL;
171  fhZNAVzeroATrack[ic]=NULL;
172  fhImpactParameter[ic]=NULL;
173 
174  for(Int_t ir=0; ir<kRho; ir++){
175  fhJetPtGen[ic][ir]=NULL;
176  fhJetPtGenVsJetPtRec[ic][ir]=NULL;
177  fhJetPtResolutionVsPtGen[ic][ir]=NULL;
178 
179  }
180 
181  for(Int_t ir=0; ir<kRho-1; ir++){
182  fhRhoIncl[ic][ir]=NULL;
183  fhDeltaPtIncl[ic][ir]=NULL;
184  }
185 
186  fhPtTrkTruePrimRec[ic]=NULL;
187  fhPtTrkTruePrimGen[ic]=NULL;
188  fhPtTrkSecOrFakeRec[ic]=NULL;
189  }
190  //Centrality bin borders
191  fCentralityBins[0]=0.;
192  fCentralityBins[1]=20.;
193  fCentralityBins[2]=50.;
194  fCentralityBins[3]=100.;
195  fCentralityBins[4]=1e6; //centrality overflow bin
196 
197  ficb[0]=-1;
198  ficb[1]=-1;
199  ftmpArray[0]=0.;
200  ftmpArray[1]=0.;
201  ftmpArrayX[0]=0.;
202  ftmpArrayX[1]=0.;
203  ftmpArrayX[2]=0.;
204  fVtxArray[0]=0.;
205  fVtxArray[1]=0.;
206  fVtxArray[2]=0.;
207 
208 
209  for(Int_t i=0; i<999; i++){
210  frhovec[i] = 0.;
211 
212  for(Int_t it=0; it<kTT; it++){
213  fTrigTracksGen[it][i]=0x0;
214  fTrigTracks[it][i]=0x0;
215  }
216  }
217 
218  fTTlow[kRef] = 6.0;
219  fTThigh[kRef]= 7.0;
220  fTTlow[kSig] = 12.0;
221  fTThigh[kSig]= 50.0;
222 
223 
224  for(Int_t i=0; i<2; i++){
225  fhInvPtQVsPhi[i] = NULL;
226  fhInvPtQVsEta[i] = NULL;
227  fhInvPtQVsPhiASide[i] = NULL;
228  fhInvPtQVsPhiCSide[i] = NULL;
229  fhSigmaPtOverPtVsPt[i] = NULL;
230  }
231 }
232 
233 //________________________________________________________________________
234 AliAnalysisTaskHJetSpectra::AliAnalysisTaskHJetSpectra(const char *name) :
235 AliAnalysisTaskEmcalJet(name,kTRUE),
236 fCollisionSystem(0), fTypeOfData(0), fTypeOfAnal(0),
237  fUseDefaultVertexCut(1), fUsePileUpCut(1),
238  fRhoTaskName(), fRhoTaskNameMC(),
239 fSignalJetRadius(0.4), fSignalJetRadiusSquared(fSignalJetRadius*fSignalJetRadius),
240 fSignalJetEtaWindow(0.9 - fSignalJetRadius), fTrackEtaWindow(0.9), fMinTrackPt(0.150), fMinJetArea(0.0),
241 fCentralityType("V0A"), fMinFractionShared(0.5),
242 fCrossSection(0.0), fTrials(0.0), fImpParam(-1.0), fRandom(0), fHelperClass(0), fInitializedLocal(0),
243 fTTType(0), fDphiCut(TMath::Pi()-0.6), fUseDoubleBinPrecision(0),
244 fHistEvtSelection(0x0), fhKTAreaPt(0x0),
245  fhVertexZ(0x0), fhVertexXAccept(0x0), fhVertexYAccept(0x0), fhVertexZAccept(0x0),
246  fhVertexXAcceptTT(0x0), fhVertexYAcceptTT(0x0), fhVertexZAcceptTT(0x0),
247 fhVertexZMC(0x0), fhVertexZAcceptMC(0x0),
248 fhDphiTriggerJetAccept(0x0), fhJetPhiIncl(0x0), fhJetEtaIncl(0x0),
249 fhCentralityV0M(0x0), fhCentralityV0A(0x0), fhCentralityV0C(0x0), fhCentralityZNA(0x0),
250 /*fh1Xsec(0x0), fh1Trials(0x0), fh1PtHard(0x0),*/
251 fhDiffPtVsPtTrackTrue(0x0),
252 fhTrackPhiCG(0x0), fhTrackPhiTPCG(0x0),
253 fhDCAinXVsPt(0x0),
254 fhDCAinYVsPt(0x0),
255 fhDCAinXVsPtStrange(0x0),
256 fhDCAinYVsPtStrange(0x0),
257 fhDCAinXVsPtNonStrange(0x0),
258 fhDCAinYVsPtNonStrange(0x0),
259 fCentralityBins(kCAll),
260 fTrackPtRef(0x0),
261 fTrackPtSig(0x0),
262 fNofRandomCones(1),
263 fZVertexCut(10.0),fCutPhi(0.6),
264 fpyVtx(3)
265 {
266 //Constructor
267  for(Int_t it=0; it<kTT;it++){
268  fRhoRec[it].Set(kRho);
269  fRhoMC[it].Set(kRho);
270  fhDphiTTTT[it] = NULL;
271  fhCentralityTT[it] = NULL;
272  }
273 
274  for(Int_t iw=0; iw<21; iw++){
275  fhPtJetPrimVsPtJetRec[iw] = NULL;
276  }
277 
278  for(Int_t ic =0; ic<kCAll; ic++){
279  for(Int_t it=0; it<kTT;it++){
280  fh1Ntriggers[ic][it]=NULL;
281  fh1TriggerMult[ic][it]=NULL;
282  fh1NtriggersGen[ic][it]=NULL;
283  fh1TriggerMultGen[ic][it]=NULL;
284 
285  fhImpactParameterTT[ic][it]=NULL;
286 
287 
288  for(Int_t ir=0; ir<kRho; ir++){
289  fhDphiTriggerJet[ic][it][ir]=NULL;
290  fhDphiTriggerJetGen[ic][it][ir]=NULL;
291  fHJetSpec[ic][it][ir]=NULL;
292  fHJetSpecGen[ic][it][ir]=NULL;
293 
294  fhDeltaPt[ic][it][ir]=NULL;
295  fhDeltaPtEmb[ic][it][ir]=NULL;
296  fhDeltaPtEmb2D[ic][it][ir]=NULL;
297  fhDeltaPtEmbPerp[ic][it][ir]=NULL;
298  fhDeltaPtEmbPerp2D[ic][it][ir]=NULL;
299  fhDeltaPtEmbBc2Bc[ic][it][ir]=NULL;
300  fhDeltaPtEmbBc2Bc2D[ic][it][ir]=NULL;
301 
302  fhRhoTT[ic][it][ir]=NULL;
303  fARhoTT[ic][it][ir]=NULL;
304  }
305  fhVzeroATotMultTT[ic][it]=NULL;
306  fhZNAEnergyTT[ic][it]=NULL;
307  fhTrackMultiplicityTT[ic][it]=NULL;
308  fhZNAVzeroATrackTT[ic][it]=NULL;
309  }
310  fhJetPhi[ic]=NULL;
311  fhJetPhiGen[ic]=NULL;
312  fhTrackPhi[ic]=NULL;
313  fhJetEta[ic]=NULL;
314  fhJetEtaGen[ic]=NULL;
315  fhTrackEta[ic]=NULL;
316  fhTrackPt[ic]=NULL;
317  fhTrackPtGen[ic]=NULL;
318  fhCentrality[ic]=NULL;
319  fhVzeroATotMult[ic]=NULL;
320  fhZNAEnergy[ic]=NULL;
321  fhTrackMultiplicity[ic]=NULL;
322  fhZNAVzeroATrack[ic]=NULL;
323  fhImpactParameter[ic]=NULL;
324 
325  for(Int_t ir=0; ir<kRho; ir++){
326  fhJetPtGen[ic][ir]=NULL;
327  fhJetPtGenVsJetPtRec[ic][ir]=NULL;
328  fhJetPtResolutionVsPtGen[ic][ir]=NULL;
329 
330  }
331 
332  for(Int_t ir=0; ir<kRho-1; ir++){
333  fhRhoIncl[ic][ir]=NULL;
334  fhDeltaPtIncl[ic][ir]=NULL;
335  }
336 
337  fhPtTrkTruePrimRec[ic]=NULL;
338  fhPtTrkTruePrimGen[ic]=NULL;
339  fhPtTrkSecOrFakeRec[ic]=NULL;
340 
341  }
342 
343  //Centrality bin borders
344  fCentralityBins[0]=0.;
345  fCentralityBins[1]=20.;
346  fCentralityBins[2]=50.;
347  fCentralityBins[3]=100.;
348  fCentralityBins[4]=1e6; //centrality overflow bin
349 
350  ficb[0]=-1;
351  ficb[1]=-1;
352  ftmpArray[0]=0.;
353  ftmpArray[1]=0.;
354  ftmpArrayX[0]=0.;
355  ftmpArrayX[1]=0.;
356  ftmpArrayX[2]=0.;
357  fVtxArray[0]=0.;
358  fVtxArray[1]=0.;
359  fVtxArray[2]=0.;
360 
361 
362  for(Int_t i=0; i<999; i++){
363  frhovec[i] = 0.;
364  for(Int_t it=0; it<kTT; it++){
365  fTrigTracksGen[it][i]=0x0;
366  fTrigTracks[it][i]=0x0;
367  }
368  }
369 
370  fTTlow[kRef] = 6.0;
371  fTThigh[kRef]= 7.0;
372  fTTlow[kSig] = 12.0;
373  fTThigh[kSig]= 50.0;
374 
375  for(Int_t i=0; i<2; i++){
376  fhInvPtQVsPhi[i] = NULL;
377  fhInvPtQVsEta[i] = NULL;
378  fhInvPtQVsPhiASide[i] = NULL;
379  fhInvPtQVsPhiCSide[i] = NULL;
380  fhSigmaPtOverPtVsPt[i] = NULL;
381  }
382 
383  fTrackPtRef = new TF1("fTrackPtRef","[0]*exp(-x*[1])",fTTlow[kRef],fTThigh[kRef]);
384  fTrackPtRef->SetParameters(1.46886e+08,8.37087e-01);
385 
386  //inclusive pT spectrum times the boost function
387  fTrackPtSig = new TF1("fTrackPtSig","([0]*exp(-x*[1])+[2]*exp(-x*[3])+[4]*exp(-x*[5]))*([6]+x*[7]+x*x*[8])",fTTlow[kSig],fTThigh[kSig]);
388  fTrackPtSig->SetParameters(1.85619e+07,6.08943e-01,3.10336e+05,2.71288e-01,2.28454e+03,1.01523e-01, 8.43568e-01, 6.64789e-03, 2.13813e-04);
389 
390  DefineOutput(1, TList::Class());
391 }
392 //________________________________________________________________________
393 void AliAnalysisTaskHJetSpectra::SetTT(Double_t tlr, Double_t thr,Double_t tls, Double_t ths){
394  //set trigger track pT bins
395  fTTlow[kRef] = tlr;
396  fTThigh[kRef] = thr;
397  fTTlow[kSig] = tls;
398  fTThigh[kSig] = ths;
399 }
400 //________________________________________________________________________
401 Double_t AliAnalysisTaskHJetSpectra::GetConePt(Double_t eta, Double_t phi, Double_t radius, AliParticleContainer *recTrkCont, Bool_t isGen){
402  //sum up pt inside a cone
403  Double_t tmpConePt = 0.0;
404 
405  if(!recTrkCont) return 0.0;
406 
407  AliVParticle* tmpTrack=NULL;
408 
409  if(recTrkCont){
410  recTrkCont->ResetCurrentID();
411  while((tmpTrack = recTrkCont->GetNextAcceptParticle())){
412  if(!tmpTrack) continue;
413  if(IsTrackInAcceptance(tmpTrack, isGen)){
414  if(GetDeltaR(tmpTrack->Phi(), phi, tmpTrack->Eta(), eta) < radius){
415  tmpConePt = tmpConePt + tmpTrack->Pt();
416  }
417  }
418  }
419  }
420  return tmpConePt;
421 }
422 
423 //________________________________________________________________________
424 Double_t AliAnalysisTaskHJetSpectra::GetImpactParameter(){
425  //get impact parameter from hijing
426  AliGenHijingEventHeader* hijingHeader = dynamic_cast<AliGenHijingEventHeader*>(MCEvent()->GenEventHeader());
427  if(MCEvent()){
428  if(!hijingHeader){
429  // Check if AOD
430  AliAODMCHeader* aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
431 
432  if(aodMCH){
433  for(UInt_t i = 0;i<aodMCH->GetNCocktailHeaders();i++){
434  hijingHeader = dynamic_cast<AliGenHijingEventHeader*>(aodMCH->GetCocktailHeader(i));
435  if(hijingHeader) break;
436  }
437  }
438  }
439  }
440  if(hijingHeader){
441  return (Double_t) (hijingHeader->ImpactParameter());
442  }
443  AliWarning(Form("In task %s: GetImpactParameter() failed!", GetName()));
444  return -1.0;
445 }
446 //________________________________________________________________________
447 Double_t AliAnalysisTaskHJetSpectra::GetSimPrimaryVertex(){
448  //get generator level primary vertex
449 
450  AliGenEventHeader* mcHeader = NULL;
451  AliAODMCHeader* aodMCH = NULL;
452 
453  if(fTypeOfAnal == kKine){ //KINE
454  AliRunLoader *rl = AliRunLoader::Instance();
455  if(rl) mcHeader = dynamic_cast<AliGenPythiaEventHeader*>(rl->GetHeader()->GenEventHeader());
456 
457  } else {
458 
459  if(MCEvent()){
460  if(fTypeOfData == kPythia){
461  mcHeader = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
462  if(!mcHeader){
463  // Check if AOD
464  aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
465 
466  if(aodMCH){
467  for(UInt_t i = 0; i<aodMCH->GetNCocktailHeaders(); i++){
468  mcHeader = dynamic_cast<AliGenPythiaEventHeader*>(aodMCH->GetCocktailHeader(i));
469  if(mcHeader) break;
470  }
471  }
472  }
473  }else if(fTypeOfData == kHijing){
474  mcHeader = dynamic_cast<AliGenHijingEventHeader*>(MCEvent()->GenEventHeader());
475  if(!mcHeader){
476  // Check if AOD
477  aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
478 
479  if(aodMCH){
480  for(UInt_t i = 0; i<aodMCH->GetNCocktailHeaders(); i++){
481  mcHeader = dynamic_cast<AliGenHijingEventHeader*>(aodMCH->GetCocktailHeader(i));
482  if(mcHeader) break;
483  }
484  }
485  }
486  }else if(fTypeOfData == kDmpjet){
487  mcHeader = dynamic_cast<AliGenDPMjetEventHeader*>(MCEvent()->GenEventHeader());
488  if(!mcHeader){
489  // Check if AOD
490  aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
491 
492  if(aodMCH){
493  for(UInt_t i = 0; i<aodMCH->GetNCocktailHeaders(); i++){
494  mcHeader = dynamic_cast<AliGenDPMjetEventHeader*>(aodMCH->GetCocktailHeader(i));
495  if(mcHeader) break;
496  }
497  }
498  }
499  }
500  }
501  }
502 
503 
504  if(mcHeader){
505 
506  mcHeader->PrimaryVertex(fpyVtx);
507  return (Double_t) (fpyVtx[2]);
508  }
509  AliWarning(Form("In task %s: Pythia Vertex failed!", GetName()));
510  return 9999.0;
511 }
512 
513 //________________________________________________________________________
514 Bool_t AliAnalysisTaskHJetSpectra::IsMCEventInAcceptance(AliVEvent* event){
515  //EVENT SELECTION PURE MC
516 
517  if(!event) return kFALSE;
518 
519  //___________________________________________________
520 
521  if(fTypeOfAnal!= kKine && !MCEvent()) return kFALSE;
522 
523 
524  Double_t vtxMC = GetSimPrimaryVertex();
525  fhVertexZMC->Fill(vtxMC); //Fill BEFORE vertex cut
526 
527  if(TMath::Abs(vtxMC) > fZVertexCut){
528  fHistEvtSelection->Fill(3); //count events rejected by vertex cut
529  return kFALSE;
530  }
531  fhVertexZAcceptMC->Fill(vtxMC);//Fill AFTER vertex cut
532 
533  return kTRUE;
534 
535 }
536 //________________________________________________________________________
537 Bool_t AliAnalysisTaskHJetSpectra::IsEventInAcceptance(AliVEvent* event){
538  //EVENT SELECTION RECONSTRUCTED DATA
539 
540 
541  if(!event) return kFALSE;
542 
543 
544  //___________________________________________________
545  //TEST PILE UP
546  if(fUsePileUpCut){
547  if(!fHelperClass || fHelperClass->IsPileUpEvent(event)){
548  fHistEvtSelection->Fill(2); //count events rejected by pileup
549  return kFALSE;
550  }
551  }
552  //___________________________________________________
553  //BEFORE VERTEX CUT
554  fhVertexZ->Fill(event->GetPrimaryVertex()->GetZ());
555 
556  if(fUseDefaultVertexCut){
557  if(!fHelperClass || !fHelperClass->IsVertexSelected2013pA(event)){
558  fHistEvtSelection->Fill(3); //count events rejected by vertex cut
559  return kFALSE;
560  }
561  }else{
562  if(TMath::Abs(event->GetPrimaryVertex()->GetZ()) > fZVertexCut){
563  fHistEvtSelection->Fill(3); //count events rejected by vertex cut
564  return kFALSE;
565  }
566  if(event->GetPrimaryVertex()->GetNContributors()<1){
567  fHistEvtSelection->Fill(3); //count events rejected by vertex cut
568  return kFALSE;
569  }
570  }
571  //___________________________________________________
572  //AFTER VERTEX CUT
573  fhVertexXAccept->Fill(event->GetPrimaryVertex()->GetX());
574  fhVertexYAccept->Fill(event->GetPrimaryVertex()->GetY());
575  fhVertexZAccept->Fill(event->GetPrimaryVertex()->GetZ());
576  fVtxArray[0]=event->GetPrimaryVertex()->GetX();
577  fVtxArray[1]=event->GetPrimaryVertex()->GetY();
578  fVtxArray[2]=event->GetPrimaryVertex()->GetZ();
579 
580 
581  return kTRUE;
582 }
583 
584 //________________________________________________________________________
585 Bool_t AliAnalysisTaskHJetSpectra::IsTrackInAcceptance(AliVParticle* track, Bool_t isGen){
586  // Check if the track pt and eta range
587  if(!track) return kFALSE;
588 
589  if(isGen){ //pure MC select charged primary tracks
590  //Apply only for kine level or MC containers
591  if(!track->Charge()) return kFALSE;
592  if(fTypeOfAnal == kEff){
593  if(!(static_cast<AliAODMCParticle*>(track))->IsPhysicalPrimary()) return kFALSE;
594  }
595  }
596  if(TMath::Abs(track->Eta()) <= fTrackEtaWindow){ //APPLY TRACK ETA CUT
597  if(track->Pt() >= fMinTrackPt){ //APPLY TRACK CUT
598  return kTRUE;
599  }
600  }
601  return kFALSE;
602 }
603 
604 
605 //________________________________________________________________________
606 Bool_t AliAnalysisTaskHJetSpectra::IsSignalJetInAcceptance(AliEmcalJet *jet, Bool_t suppressGhost){
607  //select jets in acceptance
608  if(!jet) return kFALSE;
609  if(TMath::Abs(jet->Eta()) <= fSignalJetEtaWindow){
610  if(jet->Area() >= fMinJetArea){
611  if(suppressGhost){
612  if(jet->Pt() >= fMinTrackPt) return kTRUE;
613  }else{
614  return kTRUE;
615  }
616  }
617  }
618  return kFALSE;
619 }
620 
621 
622 //________________________________________________________________________
623 void AliAnalysisTaskHJetSpectra::ExecOnceLocal(){
624  // Initialization of jet containers done in AliAnalysisTaskEmcalJet::ExecOnce()
625  //Read arrays of jets and tracks
626  fInitializedLocal = kTRUE;
627 
628  // Initialize helper class (for vertex selection & pile up correction)
629  fHelperClass = new AliAnalysisUtils();
630  fHelperClass->SetCutOnZVertexSPD(kFALSE); // kFALSE: no cut; kTRUE: |zvtx-SPD - zvtx-TPC|<0.5cm
631 
632  return;
633 }
634 
635 
636 //________________________________________________________________________
637 void AliAnalysisTaskHJetSpectra::GetDeltaPt(Int_t nrho, TArrayD &rho, Double_t *dpt,
638  Double_t ttPhi, Double_t ttEta, AliParticleContainer *recTrkCont, Bool_t isGen){
639  //Double_t leadingJetExclusionProbability)
640 
641  //delta pt = pT in random cone - rho * Area of random cone
642  // processes real reconstructed data. Exclude region around jet with TT
643 
644  for(Int_t ir=0;ir<nrho;ir++){
645  dpt[ir] = -10000.0; // Set an invalid delta pt
646  }
647 
648  // Define random cone Eta+Phi
649  Bool_t coneValid = kTRUE;
650  Double_t tmpRandConeEta = -fSignalJetEtaWindow + fRandom->Rndm()*2*fSignalJetEtaWindow;
651  Double_t tmpRandConePhi = fRandom->Rndm()*TMath::TwoPi();
652 
653  if(ttEta > -2.0){ //make sure that one generates RC far away from the trigger track jet
654  while(GetDeltaR( tmpRandConePhi, ttPhi, tmpRandConeEta, ttEta)<fSignalJetRadius+0.1){
655  tmpRandConeEta = -fSignalJetEtaWindow + fRandom->Rndm()*2*fSignalJetEtaWindow;
656  tmpRandConePhi = fRandom->Rndm()*TMath::TwoPi();
657  }
658  }
659 
660  //if(fDebug>20) Printf("RND CONE RCphi=%f RCeta=%f TTphi=%f TTeta=%f",
661  // tmpRandConePhi,tmpRandConeEta,ttPhi,ttEta);
662 
663  //Based on Ncoll probability that Reject a jet based on the probability check for overlap if demanded
664  /*
665 
666  //WHAT TO DO HERE ??? SHOULD RANDOM CONES BE CORRELAtED WITH JETs and rejected based on 1./Ncoll ?
667 
668  if(leadingJetExclusionProbability>0){ //skips pp
669 
670  AliEmcalJet* tmpLeading = NULL;
671  Double_t lpt = -1.0;
672  // Get leading jet (regardless of pT)
673 
674  AliEmcalJet* tmpJet = NULL;
675  AliJetContainer *jetContRec = GetJetContainer(kContainerOne);
676  if(jetContRec){
677  jetContRec->ResetCurrentID();
678  while((tmpJet = jetContRec->GetNextAcceptJet())) {
679  if(!tmpJet) continue;
680  if((TMath::Abs(tmpJet->Eta()) <= fSignalJetEtaWindow) && (tmpJet->Area() >= fMinJetArea)){
681  if(tmpJet->Pt() > lpt){
682  tmpLeading = tmpJet;
683  lpt = tmpJet->Pt();
684  }
685  }
686  }
687 
688  if(tmpLeading){
689  Double_t tmpDeltaPhi = RelativePhi(tmpRandConePhi, tmpLeading->Phi());
690  Double_t tmpDeltaEta = tmpLeading->Eta()-tmpRandConeEta;
691 
692  // Check, if cone has overlap with jet
693  if(sqrt(tmpDeltaPhi*tmpDeltaPhi + tmpDeltaEta*tmpDeltaEta) <= fSignalJetRadius){
694  // probability to exclude the RC
695  // Only exclude cone with a given probability
696  if(fRandom->Rndm()<=leadingJetExclusionProbability) coneValid = kFALSE;
697  }
698  }
699  }
700  if(fRandom->Rndm()<=leadingJetExclusionProbability) coneValid = kFALSE;
701  }*/
702 
703  // Get the cones' pt and calculate delta pt
704  if(coneValid){
705  Double_t conePt = GetConePt(tmpRandConeEta,tmpRandConePhi, fSignalJetRadius, recTrkCont, isGen);
706  for(Int_t ir=0; ir < nrho; ir++){
707  dpt[ir] = conePt - (rho[ir]*fSignalJetRadiusSquared*TMath::Pi());
708  }
709  }
710 
711 }
712 
713 //________________________________________________________________________
714 Bool_t AliAnalysisTaskHJetSpectra::FillHistograms(){
715  // executed in each event
716  //called in AliAnalysisTaskEmcal::UserExec(Option_t *)
717  // Analyze the event and Fill histograms
718 
719  if(!InputEvent()){
720  AliError("??? Event pointer == 0 ???");
721  return kFALSE;
722  }
723 
724  //Execute only once: Get tracks, jets from arrays if not already given
725  if(!fInitializedLocal) ExecOnceLocal();
726 
727  //_________________________________________________________________
728  // fill MC information
729 
730  //if(fTypeOfAnal == kKine && fTypeOfData == kPythia){
731  // fh1PtHard->Fill(GetPtHard()); //Fills cross section
732  //}
733 
734  //_________________________________________________________________
735  // FILL EVENT STATISTICS
736  fHistEvtSelection->Fill(1); //Count input event
737 
738  if(fTypeOfData != kReal){ //Check MC event vertex
739  if(!IsMCEventInAcceptance(InputEvent())) return kFALSE; //post data is in UserExec
740  }
741 
742  //Check Reconstructed event vertex
743  if(fTypeOfAnal != kKine){ //Check MC event vertex
744  if(!IsEventInAcceptance(InputEvent())) return kFALSE; //post data is in UserExec
745  }
746 
747  //___________________
748  // Get centrality
749  Double_t centralityPercentile = -1.0; //KINE
750  Double_t centralityPercentileV0A = -1.0;
751  Double_t centralityPercentileV0C = -1.0;
752  Double_t centralityPercentileV0M = -1.0;
753  Double_t centralityPercentileZNA = -1.0;
754  Double_t multVzero = -1.0;
755  Double_t energyZdcNA = -1.0;
756  ficb[0] = 0; //0=min bias
757  ficb[1] = -1; //-1 assigned centrality bin
758 
759  if(InputEvent()->GetVZEROData()){
760  multVzero = InputEvent()->GetVZEROData()->GetMTotV0A();
761  }
762 
763  //if(InputEvent()->GetZDCData()){
764  // AliESDZDC.h : C --> 1 A--->2
765  //Double_t GetZNAEnergy() const {return (Double_t) fZDCN2Energy;}
766  //Double_t GetZNCEnergy() const {return (Double_t) fZDCN1Energy;}
767  // energyZdcNA = InputEvent()->GetZDCData()->GetZNAEnergy();
768  //energyZdcNA = InputEvent()->GetZDCN2Energy(); //should be ZNA
769  //}
770  AliVZDC *aodZDC = dynamic_cast<AliVZDC*> (InputEvent()->GetZDCData());
771  if(aodZDC){
772  const Double_t *ZNAtower = aodZDC->GetZNATowerEnergy();
773  energyZdcNA = ZNAtower[0]*4.*82./208./12.96; //ZNA energy in TeV
774  }
775 
776  if(fCollisionSystem != kpp){ //KINE Check MC event vertex
777  AliCentrality* tmpCentrality = InputEvent()->GetCentrality();
778  if(!tmpCentrality){
779  fHistEvtSelection->Fill(4);
780  return kFALSE; //post data is in UserExec
781  }else{
782  centralityPercentile = tmpCentrality->GetCentralityPercentile(fCentralityType.Data());
783  centralityPercentileV0A = tmpCentrality->GetCentralityPercentile("V0A");
784  centralityPercentileV0C = tmpCentrality->GetCentralityPercentile("V0C");
785  centralityPercentileV0M = tmpCentrality->GetCentralityPercentile("V0M");
786  centralityPercentileZNA = tmpCentrality->GetCentralityPercentile("ZNA");
787  }
788 
789  if((centralityPercentile < fCentralityBins[0]) || (centralityPercentile > fCentralityBins[4])){ //cut on centrality
790  AliWarning(Form("Centrality value not valid (c=%E)",centralityPercentile));
791  fHistEvtSelection->Fill(4);
792  return kFALSE;
793  }
794 
795  if(!(fTypeOfData == kPythia || fTypeOfData == kDmpjet || fTypeOfAnal == kKine)){
796  for(Int_t ic=0; ic<fCentralityBins.GetSize()-1;ic++){
797  if(fCentralityBins[ic] <= centralityPercentile &&
798  centralityPercentile < fCentralityBins[ic+1]){
799  ficb[1] = ic+1; //MB is ficb[0]=0; other centralities ficb[1]>0
800  }
801  }
802  }
803 
804  for(Int_t ic=0; ic<2; ic++){
805  if(ficb[ic]==-1) continue; //MB + CENT BIASED
806  fhCentrality[ficb[ic]]->Fill(centralityPercentile);
807  }
808 
809  //MB
810  fhCentralityV0M->Fill(centralityPercentileV0M);
811  fhCentralityV0A->Fill(centralityPercentileV0A);
812  fhCentralityV0C->Fill(centralityPercentileV0C);
813  fhCentralityZNA->Fill(centralityPercentileZNA);
814  }
815 
816  for(Int_t ic=0; ic<2; ic++){
817  if(ficb[ic]==-1) continue; //MB +CENT bias
818  fhVzeroATotMult[ficb[ic]]->Fill(multVzero);
819  fhZNAEnergy[ficb[ic]]->Fill(energyZdcNA);
820 
821  if(fTypeOfData == kHijing){
822  fImpParam = GetImpactParameter();
823  fhImpactParameter[ficb[ic]]->Fill(fImpParam);
824  }
825  }
826  //cout<<"energyZdcNA = "<<energyZdcNA<<endl;
827 
828  fHistEvtSelection->Fill(0); //Count Accepted input event
829 
830  // END EVENT SELECTION
831  //_________________________________________________________________
832  TClonesArray* arrayMC = 0; // array particles in the MC event
833  AliAODEvent* aodIn = 0;
834  Int_t iNTracksMC = 0; //number of particles in arrayMC
835  if(fTypeOfAnal == kEff){
836  aodIn = dynamic_cast<AliAODEvent*>(InputEvent()); // input AOD
837  if(!aodIn){
838  AliError("No input AOD found!");
839  return kFALSE;
840  }
841  arrayMC = (TClonesArray*) aodIn->FindListObject(AliAODMCParticle::StdBranchName());
842  if(!arrayMC){
843  AliError("No MC array found!");
844  return kFALSE;
845  }
846 
847  iNTracksMC = arrayMC->GetEntriesFast();
848  }
849 
850 
851  //_________________________________________________________________
852  // JET+TRACK CONTAINERS
853  AliJetContainer *jetContRec = NULL; //AKTjet container from reconstruced tracks
854  AliJetContainer *jetContGen = NULL; //AKT jet container from MC particles
855  AliJetContainer *jetContRecKT = NULL; //KT jet container from reconstruced tracks
856  AliJetContainer *jetContGenKT = NULL; //KT jet container from MC particles
857 
858  AliEmcalJet *jetGen = NULL; //jet pointer MC jet
859  AliEmcalJet *jetRec = NULL; //jet pointer real jet
860 
861  AliParticleContainer *trkContRec = NULL; //track contaier real reconstructed tracks
862  AliParticleContainer *parContGen = NULL; //particle container of MC particles
863 
864  AliVParticle *constTrackRec = NULL; //jet constituent
865  AliVParticle *constTrackGen = NULL; //jet constituent
866  //_________________________________________________________
867  //READ JET TRACK CONTAINERS
868  if(fTypeOfAnal == kRec){
869  //REAL DATA
870  jetContRec = GetJetContainer(kContainerOne); //AKT jet
871  trkContRec = GetParticleContainer(kContainerOne); //reconstructed particle container
872  jetContRecKT = GetJetContainer(kContainerTwo);//KT jet
873  }
874 
875  if(fTypeOfAnal == kEff || fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
876  //MC + EMB DATA
877  jetContRec = GetJetContainer(kContainerOne); //AKT
878  trkContRec = GetParticleContainer(kContainerOne); //reconstructed particle container
879  jetContGen = GetJetContainer(kContainerTwo); //AKT generator level jets
880  parContGen = GetParticleContainer(kContainerTwo); //true MC particle container
881 
882  jetContRecKT = GetJetContainer(kContainerThree); //KT reconstructed jets
883  jetContGenKT = GetJetContainer(kContainerFour); //KT generator level jets
884  }
885 
886  if(fTypeOfAnal == kKine){ // Kine written to the 0th container !!!!!!!!!!!!!!!!
887  //KINE
888  jetContGen = GetJetContainer(kContainerOne); //AKT jets
889  parContGen = GetParticleContainer(kContainerOne); //kine particle container
890 
891  jetContGenKT = GetJetContainer(kContainerTwo); //KT jets
892  }
893 
894  //if(fDebug>20) Printf("POINTER TO CONTAINERS JETrec=%p TRKrec=%p JETgen=%p TRKgen=%p",
895  // jetContRec,trkContRec,jetContGen, parContGen);
896 
897 
898  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
899  //LOOP OVER TRACKS SEARCH FOR TRIGGER CANDIDATES IN GEN TRACKS
900  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
901 
902  Int_t indexSingleRndTrigGen[kTT] = {-1,-1}; //index of single random trigger
903  AliVParticle* trackTTGen[kTT] = {NULL,NULL}; //pointer to the trigger hadron
904  Int_t ntriggersGen[kTT] = {0,0};
905 
906  Double_t sumFakeTrackPtInJetNotStrange = 0.;//sum of track pT from non strange secondaries
907  Double_t sumFakeTrackPtInJetStrange = 0.;//sum of track pT from strange secondaries
908  Double_t sumStrangeTrackPtInJet = 0.;//sum of track pT from strange secondaries
909  AliAODMCParticle* particleMC = NULL; //
910  AliAODMCParticle* particleMCMother = NULL;//
911  Int_t iIndexMother = 0; //
912  Int_t iPdgCode = 0;//
913  Int_t iPdgCodeMother = 0;//
914  Bool_t bStrange = 0;
915  Int_t lab = 0;
916  Double_t newJetPt=0., newSumFakePt=0., enhw =0.; //
917  TF1 *tmpFunc;
918  Double_t normprob, prob, probTT, dphiTTTT;
919 
920  if(fTypeOfAnal == kEff || fTypeOfAnal == kKine){
921 
922 
923  if(parContGen){
924  parContGen->ResetCurrentID();
925  while((constTrackGen = (AliVParticle*) parContGen->GetNextAcceptParticle())){
926  if(!constTrackGen) continue;
927 
928  if(IsTrackInAcceptance(constTrackGen, kTRUE)){
929  for(Int_t ic=0; ic<2; ic++){
930  if(ficb[ic]==-1) continue;
931  fhTrackPtGen[ficb[ic]]->Fill(constTrackGen->Pt()); //inclusive pT spectrum of tracks
932  }
933 
934  for(Int_t it=0; it<kTT; it++){
935  if((fTTlow[it] <= constTrackGen->Pt()) && (constTrackGen->Pt() < fTThigh[it])){
936  if(ntriggersGen[it]<999){
937  fTrigTracksGen[it][ntriggersGen[it]] = constTrackGen; //GEN trigger candidates
938  ntriggersGen[it]++;
939  }
940  }
941  }
942  }
943  }
944  }
945 
946 
947  for(Int_t it=0; it<kTT; it++){
948  for(Int_t ic=0; ic<2; ic++){
949  if(ficb[ic]==-1) continue;
950  fh1TriggerMultGen[ficb[ic]][it]->Fill(ntriggersGen[it]);
951  }
952 
953  if(ntriggersGen[it]==1){ //only one trigger candidate
954  indexSingleRndTrigGen[it] = 0;
955  trackTTGen[it] = (AliVParticle*) (fTrigTracksGen[it][indexSingleRndTrigGen[it]]);
956 
957  }else if(ntriggersGen[it]>1){ //more trigger candiates
958 
959  if(fTTType){ //TT selection weighted according to MB
960  tmpFunc = (it==kRef) ? fTrackPtRef : fTrackPtSig;
961  normprob = 0;
962  for(Int_t ik=0; ik<ntriggersGen[it]; ik++){
963  normprob += tmpFunc->Eval(((AliVParticle*) (fTrigTracksGen[it][ik]))->Pt());
964  }
965  prob = fRandom->Uniform(0,1);
966  probTT = 0;
967  if(normprob>0){
968  for(Int_t ik=0; ik<ntriggersGen[it]; ik++){
969  probTT += tmpFunc->Eval(((AliVParticle*) (fTrigTracksGen[it][ik]))->Pt())/normprob;
970  if(prob <= probTT){
971  indexSingleRndTrigGen[it] = ik;
972  break;
973  }
974  }
975  }else{
976  indexSingleRndTrigGen[it] = 0;
977  }
978  trackTTGen[it] = (AliVParticle*) (fTrigTracksGen[it][indexSingleRndTrigGen[it]]);
979  }else{ //TT selection without weighting
980  indexSingleRndTrigGen[it] = fRandom->Integer(ntriggersGen[it]); //Integer 0 ... ntriggers-1
981  trackTTGen[it] = (AliVParticle*) (fTrigTracksGen[it][indexSingleRndTrigGen[it]]);
982  }
983  }
984  }
985  }
986  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
987  //LOOP OVER TRACKS SEARCH FOR TRIGGER CANDIDATES IN REC TRACKS
988  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
989 
990  Int_t indexSingleRndTrig[kTT] = {-1,-1}; //index of single random trigger
991  AliVParticle* trackTT[kTT] = {NULL,NULL};
992  Int_t ntriggers[kTT] = {0, 0};
993 
994  if(fTypeOfAnal != kKine){
995 
996  if(trkContRec){
997  trkContRec->ResetCurrentID();
998  while((constTrackRec = (AliVParticle*) (trkContRec->GetNextAcceptParticle()))){
999  if(!constTrackRec) continue;
1000 
1001  if(fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
1002  //if(fDebug>99) Printf("TRACK LABEL %d", track->GetLabel());
1003  //in embed evts search for TT only among real tracks; do not consider embedded tracks as trigger
1004  if(TMath::Abs(constTrackRec->GetLabel()) == 99999) continue;//9999 set in AddTaskHJetSpectra.C
1005  }
1006 
1007 
1008  if(IsTrackInAcceptance(constTrackRec, kFALSE)){ //rec
1009  //Fill some inclusive spectra (Same analysis for gen level?)
1010  for(Int_t ic=0; ic<2; ic++){
1011  if(ficb[ic]==-1) continue;
1012  fhTrackPhi[ficb[ic]]->Fill(constTrackRec->Pt(), RelativePhi(constTrackRec->Phi(),0.0)); // phi = -pi,pi
1013  fhTrackEta[ficb[ic]]->Fill(constTrackRec->Pt(), constTrackRec->Eta());
1014  fhTrackPt[ficb[ic]]->Fill(constTrackRec->Pt());
1015  }
1016 
1017  for(Int_t it=0; it<kTT; it++){
1018  if(fTTlow[it] <= constTrackRec->Pt() && constTrackRec->Pt() < fTThigh[it]){
1019  if(ntriggers[it]<999){
1020  fTrigTracks[it][ntriggers[it]] = constTrackRec; //trigger candidates
1021  ntriggers[it]++;
1022  }
1023  }
1024  }
1025  }
1026  }
1027  }
1028 
1029  // SELECT SINGLE INCLUSIVE REC LEVEL TRIGGER
1030  for(Int_t it=0; it<kTT; it++){
1031  for(Int_t ic=0; ic<2; ic++){
1032  if(ficb[ic]==-1) continue;
1033  fh1TriggerMult[ficb[ic]][it]->Fill(ntriggers[it]);
1034  }
1035 
1036  // if(ntriggers[it]>0){
1037  // indexSingleRndTrig[it] = fRandom->Integer(ntriggers[it]); //Integer 0 ... ntriggers-1
1038  // trackTT[it] = (AliVParticle*) (fTrigTracks[it][indexSingleRndTrig[it]]);
1039  //if(fDebug>20) Printf("TT index = %d size =%d", indexSingleRndTrig, (int)fTrigTracks.size());
1040  // }
1041  if(ntriggers[it]==1){ //only one trigger candidate
1042  indexSingleRndTrig[it] = 0;
1043  trackTT[it] = (AliVParticle*) (fTrigTracks[it][indexSingleRndTrig[it]]);
1044 
1045  }else if(ntriggers[it]>1){ //more trigger candiates pickupt the trigger based on inclusive probabilities
1046 
1047  if(fTTType){ //TT selection weighted according to MB
1048  tmpFunc = (it==kRef) ? fTrackPtRef : fTrackPtSig;
1049  normprob = 0;
1050  for(Int_t ik=0; ik<ntriggers[it]; ik++){
1051  normprob += tmpFunc->Eval(((AliVParticle*) (fTrigTracks[it][ik]))->Pt());
1052 
1053  if(ik>0){ //azimuthal angle between triggers
1054  dphiTTTT = RelativePhi(((AliVParticle*) (fTrigTracks[it][0]))->Phi(),
1055  ((AliVParticle*) (fTrigTracks[it][ik]))->Phi());
1056  if(dphiTTTT<-0.5*TMath::Pi()) dphiTTTT += TMath::TwoPi();
1057  if(dphiTTTT> 1.5*TMath::Pi()) dphiTTTT -= TMath::TwoPi();
1058 
1059  fhDphiTTTT[it]->Fill(dphiTTTT);
1060  }
1061  }
1062  prob = fRandom->Uniform(0,1);
1063  probTT = 0;
1064  if(normprob>0){
1065  for(Int_t ik=0; ik<ntriggers[it]; ik++){
1066  probTT += tmpFunc->Eval(((AliVParticle*) (fTrigTracks[it][ik]))->Pt())/normprob;
1067  if(prob <= probTT){
1068  indexSingleRndTrig[it] = ik;
1069  break;
1070  }
1071  }
1072  }else{
1073  indexSingleRndTrig[it] = 0;
1074  }
1075  trackTT[it] = (AliVParticle*) (fTrigTracks[it][indexSingleRndTrig[it]]);
1076  }else{ //TT selection without weighting
1077  indexSingleRndTrig[it] = fRandom->Integer(ntriggers[it]);
1078  trackTT[it] = (AliVParticle*) (fTrigTracks[it][indexSingleRndTrig[it]]);
1079  }
1080  }
1081  }
1082 
1083  if(ntriggers[kRef]>0 || ntriggers[kSig]>0){
1084  fhVertexXAcceptTT->Fill(fVtxArray[0]);
1085  fhVertexYAcceptTT->Fill(fVtxArray[1]);
1086  fhVertexZAcceptTT->Fill(fVtxArray[2]);
1087  }
1088  }
1089 
1090 
1091  //___________________________________________________________
1092  // CALCULATE RHO
1093  for(Int_t it=0;it<kTT;it++){
1094  fRhoRec[it].Reset(0.);
1095  fRhoMC[it].Reset(0.);
1096  }
1097 
1098  for(Int_t it=0;it<kTT;it++){
1099  if(fTypeOfAnal == kRec || fTypeOfAnal == kEff || fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
1100  fRhoRec[it][kConeRho] = EstimateBgCone( jetContRec, trkContRec, trackTT[it], kFALSE); //container ID=0 reconstructed tracks
1101  fRhoRec[it][kCMSRho] = EstimateBgKTcms(jetContRecKT, trkContRec, trackTT[it]);
1102  fRhoRec[it][kKtRho] = EstimateBgKT( jetContRecKT, trkContRec, trackTT[it]);
1103  fRhoRec[it][kZeroRho] = 0.0;
1104  }
1105 
1106  if(fTypeOfAnal == kEff || fTypeOfAnal == kEmb){ //rho in MC events
1107  fRhoMC[it][kConeRho] = EstimateBgCone( jetContGen, parContGen, trackTTGen[it], kTRUE); //container ID=1 mc particles
1108  fRhoMC[it][kCMSRho] = EstimateBgKTcms(jetContGenKT, parContGen, trackTTGen[it]);
1109  fRhoMC[it][kKtRho] = EstimateBgKT( jetContGenKT, parContGen, trackTTGen[it]);
1110  fRhoMC[it][kZeroRho] = 0.0;
1111  }
1112 
1113  if(fTypeOfAnal == kEmbSingl){ //embedding single track
1114  fRhoMC[it][kConeRho] = 0.0;
1115  fRhoMC[it][kCMSRho] = 0.0;
1116  fRhoMC[it][kKtRho] = 0.0;
1117  fRhoMC[it][kZeroRho] = 0.0;
1118  }
1119 
1120  if(fTypeOfAnal == kKine){ //rho in KINE MC events
1121  fRhoMC[it][kConeRho] = EstimateBgCone( jetContGen, parContGen, trackTTGen[it], kTRUE); //container ID=1 mc particles
1122  fRhoMC[it][kCMSRho] = EstimateBgKTcms(jetContGenKT, parContGen, trackTTGen[it]);
1123  fRhoMC[it][kKtRho] = EstimateBgKT( jetContGenKT, parContGen, trackTTGen[it]);
1124  fRhoMC[it][kZeroRho] = 0.0;
1125  }
1126  }
1127 
1128 
1129  //_________________________________________________________
1130  //EVALUATE SINGLE PARTICLE EFFICIENCY + FILL RESPONSE MATRIX
1131  Bool_t bRecPrim = kFALSE; //tags the reconstructed primary particles
1132 
1133  if(fTypeOfAnal == kEff){
1134 
1135  //1) FILL HISTOS FOR SINGLE PARTICLE EFFICIENCY
1136  if(parContGen){
1137  parContGen->ResetCurrentID();
1138  while((constTrackGen = (AliVParticle*)(parContGen->GetNextAcceptParticle()))){
1139  if(!constTrackGen) continue;
1140  if(IsTrackInAcceptance(constTrackGen, kTRUE)){
1141  //pT spectrum of generator level particles
1142  for(Int_t ic=0; ic<2; ic++){
1143  if(ficb[ic]==-1) continue;
1144  fhPtTrkTruePrimGen[ficb[ic]]->Fill(constTrackGen->Pt(),constTrackGen->Eta());
1145  }
1146  }
1147  }
1148 
1149  //single particle efficiency and contamination
1150 
1151  if(trkContRec && parContGen){
1152  trkContRec->ResetCurrentID();
1153  while((constTrackRec =(AliVParticle*) (trkContRec->GetNextAcceptParticle()))){
1154  if(!constTrackRec) continue;
1155  if(!IsTrackInAcceptance(constTrackRec, kFALSE)) continue; //reconstructed level tracks
1156  bRecPrim = kFALSE; //not yet matched to generator level physical primary
1157 
1158  parContGen->ResetCurrentID();
1159  while((constTrackGen = (AliVParticle*)(parContGen->GetNextAcceptParticle()))){
1160  if(!constTrackGen) continue;
1161  if(!IsTrackInAcceptance(constTrackGen, kTRUE)) continue; //gen level physical primary
1162  if(TMath::Abs(constTrackRec->GetLabel()) == TMath::Abs(constTrackGen->GetLabel())){
1163  //has the same label as reconstr track
1164 
1165  bRecPrim = kTRUE;
1166  for(Int_t ic=0; ic<2; ic++){
1167  if(ficb[ic]==-1) continue;
1168  fhPtTrkTruePrimRec[ficb[ic]]->Fill(constTrackGen->Pt(),constTrackGen->Eta()); //this is well recontr phys primary
1169  }
1170  fhDiffPtVsPtTrackTrue->Fill(constTrackGen->Pt(), constTrackRec->Pt()-constTrackGen->Pt());
1171 
1172  break;
1173  }//same label with rec particle
1174  }//loop over gen tracks
1175  if(!bRecPrim){
1176  for(Int_t ic=0; ic<2; ic++){
1177  if(ficb[ic]==-1) continue;
1178  fhPtTrkSecOrFakeRec[ficb[ic]]->Fill(constTrackRec->Pt(),constTrackRec->Eta()); //matchnig to phys primary not found, this is fake or second.
1179  }
1180  }
1181  }//loop over rec tracks
1182  }//rec track array exists
1183  }//gen particle array exists
1184 
1185  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++
1186  //2) FILL JET RESPONSE MATRIX
1187  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++
1188  Double_t ptGenCorr; //GEN jet pt corrected for rho
1189  Double_t ptRecCorr; //REC jet pt corrected for rho
1190 
1191  //Response matrix normalization - spectrum of all generator level jets in acceptance
1192  if(jetContGen){
1193  jetContGen->ResetCurrentID();
1194  while((jetGen = jetContGen->GetNextAcceptJet())){
1195  if(!jetGen) continue;
1196  if(!IsSignalJetInAcceptance(jetGen,kTRUE)) continue; //cuts on eta, pT ,area
1197 
1198  //if(fDebug>20) Printf("GEN JET phi=%f eta=%f pt=%f", jetGen->Phi(), jetGen->Eta(), jetGen->Pt());
1199 
1200  for(Int_t ir=0; ir< kRho; ir++){
1201  ptGenCorr = jetGen->Pt() - jetGen->Area()*fRhoMC[kRef][ir]; // correct for rho
1202 
1203  for(Int_t ic=0; ic<2; ic++){
1204  if(ficb[ic]==-1) continue;
1205  fhJetPtGen[ficb[ic]][ir]->Fill(ptGenCorr);
1206  }
1207  }
1208  }
1209  }
1210 
1211  //Find closest gen level+rec level jets
1212  //Get momentum shift due to fake tracks
1213  if(jetContRec){
1214  jetContRec->ResetCurrentID();
1215  while((jetRec = jetContRec->GetNextAcceptJet())) {
1216  if(!jetRec) continue;
1217  if(!IsSignalJetInAcceptance(jetRec,kTRUE)) continue; //cuts on eta, pT ,area
1218 
1219 
1220  //Get momentum shift due to fake tracks
1221  sumFakeTrackPtInJetNotStrange = 0.;
1222  sumFakeTrackPtInJetStrange = 0.;
1223  sumStrangeTrackPtInJet = 0.;
1224  particleMC = NULL;
1225  particleMCMother = NULL;
1226  iIndexMother = 0;
1227  iPdgCode = 0;
1228  iPdgCodeMother = 0;
1229 
1230 
1231  for(Int_t iq=0; iq < jetRec->GetNumberOfTracks(); iq++) {
1232  constTrackRec = static_cast<AliVParticle*> (jetRec->TrackAt(iq,trkContRec->GetArray()));
1233  if(!constTrackRec) continue;
1234  bRecPrim = kFALSE; //not yet matched to generator level physical primary
1235 
1236  parContGen->ResetCurrentID();
1237  while((constTrackGen = (AliVParticle*)(parContGen->GetNextAcceptParticle()))){
1238  if(!constTrackGen) continue;
1239  if(!IsTrackInAcceptance(constTrackGen, kTRUE)) continue; //gen level physical primary
1240  if(TMath::Abs(constTrackRec->GetLabel()) == TMath::Abs(constTrackGen->GetLabel())){
1241  bRecPrim = kTRUE;
1242  break;
1243  }
1244  }
1245  if(!bRecPrim){ //this is a fake track
1246 
1247  lab = TMath::Abs(constTrackRec->GetLabel());
1248  bStrange = 0;
1249 
1250  if(lab < iNTracksMC){
1251  particleMC = (AliAODMCParticle*) arrayMC->At(lab);
1252  if(particleMC){
1253  iPdgCode = TMath::Abs(particleMC->GetPdgCode());
1254  iIndexMother = TMath::Abs(particleMC->GetMother());
1255  if(IsStrange(iPdgCode)){
1256  sumFakeTrackPtInJetStrange += constTrackRec->Pt();
1257  sumStrangeTrackPtInJet += constTrackRec->Pt();
1258  bStrange = 1;
1259  }
1260 
1261  if(iIndexMother < iNTracksMC && 0 <= iIndexMother && !bStrange){
1262  //check strangeness of mother of the secondary tracks
1263  particleMCMother = (AliAODMCParticle*) arrayMC->At(iIndexMother);
1264  if(particleMCMother){
1265  iPdgCodeMother = TMath::Abs(particleMCMother->GetPdgCode());
1266 
1267  if(IsStrange(iPdgCodeMother)){
1268  sumFakeTrackPtInJetStrange += constTrackRec->Pt();
1269  sumStrangeTrackPtInJet += constTrackRec->Pt();
1270  bStrange = 1;
1271  }
1272  }
1273  }
1274  }
1275  }
1276  if(!bStrange){
1277  sumFakeTrackPtInJetNotStrange += constTrackRec->Pt();
1278  }
1279  }else{
1280  //check strangeness of primary tracks
1281  lab = TMath::Abs(constTrackRec->GetLabel());
1282  if(lab < iNTracksMC){
1283  particleMC = (AliAODMCParticle*) arrayMC->At(lab);
1284  if(particleMC){
1285  iPdgCode = TMath::Abs(particleMC->GetPdgCode());
1286  if(IsStrange(iPdgCode)){
1287  sumStrangeTrackPtInJet += constTrackRec->Pt();
1288  }
1289  }
1290  }
1291  }
1292  }
1293 
1294  for(Int_t iw=0; iw<21; iw++){
1295  enhw = iw*0.05; //enhances weight of strangeness particles from 0 to 100% in steps of 5%
1296  newJetPt = jetRec->Pt() + enhw*sumStrangeTrackPtInJet; //jet contain 100 strg add only enhanc
1297  newSumFakePt = sumFakeTrackPtInJetNotStrange + (1.0+ enhw)*sumFakeTrackPtInJetStrange;
1298  fhPtJetPrimVsPtJetRec[iw]->Fill(newJetPt, newSumFakePt/newJetPt);
1299  }
1300 
1301  //if(fDebug>20) Printf("REC JET phi=%f eta=%f pt=%f",jetRec->Phi(), jetRec->Eta(), jetRec->Pt());
1302 
1303  //Find closest gen level+rec level jets
1304  jetGen = 0x0;
1305  jetGen = jetRec->ClosestJet();
1306 
1307  if(!jetGen){ //did not find matching generator level jet
1308  //if(fDebug>20) Printf("NO MATCH (NO SUCH GEN JET)");
1309 
1310  continue;
1311  }
1312  if(jetGen->Pt()<1e-3){
1313  //if(fDebug>20) Printf("SKIP MATCH WITH GHOST JET");
1314  continue;
1315  }
1316  //corresponding generator level jet found
1317  //if(fDebug>20) Printf("MATCHED WITH phi=%f eta=%f pt=%f",jetGen->Phi(), jetGen->Eta(), jetGen->Pt());
1318 
1319  if(!IsSignalJetInAcceptance(jetGen,kTRUE)) continue; //cuts on eta, pT ,area
1320 
1321  //check fraction of tracks from generator level jet in rec level jet
1322  //if(fDebug>20) Printf("FRACTIONH SHARED = %f ", GetFractionSharedPt(jetRec,jetContRec,jetGen, jetContGen));
1323 
1324  //FK// if(GetFractionSharedPt(jetRec, jetContRec, jetGen, jetContGen) < fMinFractionShared) continue;
1325 
1326  //if(fDebug>20) Printf("PASSED MIN FRACTION CRITERION ");
1327 
1328  for(Int_t ir=0; ir< kRho; ir++){
1329  ptGenCorr = jetGen->Pt() - jetGen->Area()*fRhoMC[kRef][ir]; //perp cone bg correction to pt
1330  ptRecCorr = jetRec->Pt() - jetRec->Area()*fRhoRec[kRef][ir]; //perp cone bg correction to pt
1331 
1332  for(Int_t ic=0; ic<2; ic++){
1333  if(ficb[ic]==-1) continue;
1334 
1335  fhJetPtGenVsJetPtRec[ficb[ic]][ir]->Fill(ptRecCorr, ptGenCorr); //response matrix
1336 
1337  if(ptGenCorr>0){
1338  fhJetPtResolutionVsPtGen[ficb[ic]][ir]->Fill(ptGenCorr,(ptRecCorr-ptGenCorr)/ptGenCorr); //jet pT resolution
1339  }
1340  }
1341  }
1342  }
1343  }//rec jet container exists
1344  }//analyze efficiency mode (response matrix + single particle efficiency)
1345 
1346 
1347  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1348  //if(bStop){
1349  // fHistEvtSelection->Fill(5); //Rejected events
1350  // return kTRUE; // SKIP EVENTS WHERE REF TT CLASS EVENT CONTAINS SIG TT TRACK
1351  // }
1352  fHistEvtSelection->Fill(6); //Accepted events
1353  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1354  Double_t areaJet, pTJet;
1355  Double_t dphi, dfi;
1356  Bool_t bFirstCycle = kTRUE;
1357 
1358 
1359  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1360  //H-JET CORRELATIONS IN MC TRUTH
1361  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1362  if(fTypeOfAnal == kEff || fTypeOfAnal == kKine){
1363  for(Int_t it=0; it< kTT; it++){
1364  if(parContGen && ntriggersGen[it] >0 && ((AliVParticle*)trackTTGen[it])){
1365  bFirstCycle = kTRUE;
1366 
1367  AliVParticle* triggerHadronGen = (AliVParticle*) trackTTGen[it]; //(fTrigTracksGen[it]);
1368  if(!triggerHadronGen) continue;
1369 
1370  for(Int_t ic=0; ic<2; ic++){
1371  if(ficb[ic]==-1) continue;
1372  fh1NtriggersGen[ficb[ic]][it]->Fill((Float_t) triggerHadronGen->Pt()); //trigger pT gen
1373 
1374  if( fTypeOfData == kHijing){ //impact parameter for triggered events
1375  fhImpactParameterTT[ficb[ic]][it]->Fill(fImpParam);
1376  }
1377  }
1378 
1379  //JET LOOP
1380  if(jetContGen){
1381  jetContGen->ResetCurrentID();
1382  while((jetGen = jetContGen->GetNextAcceptJet())) {
1383  if(!jetGen){
1384  AliError(Form("%s: Could not receive gen jet", GetName()));
1385  continue;
1386  }
1387  if(!IsSignalJetInAcceptance(jetGen,kTRUE)) continue;
1388 
1389  areaJet = jetGen->Area();
1390  pTJet = jetGen->Pt();
1391 
1392  if(bFirstCycle && it==kRef){
1393  for(Int_t ic=0; ic<2; ic++){
1394  if(ficb[ic]==-1) continue;
1395 
1396  fhJetPhiGen[ficb[ic]]->Fill( pTJet, RelativePhi(jetGen->Phi(),0.0));
1397  fhJetEtaGen[ficb[ic]]->Fill( pTJet, jetGen->Eta());
1398  }
1399  }
1400 
1401  dphi = RelativePhi(triggerHadronGen->Phi(), jetGen->Phi());
1402 
1403  dfi = dphi; //-0.5*pi to 1.5*Pi
1404  if(dfi < -0.5*TMath::Pi()) dfi += TMath::TwoPi();
1405  if(dfi > 1.5*TMath::Pi()) dfi -= TMath::TwoPi();
1406 
1407  for(Int_t ir=0; ir< kRho;ir++){
1408  for(Int_t ic=0; ic<2; ic++){
1409  if(ficb[ic]==-1) continue;
1410  fhDphiTriggerJetGen[ficb[ic]][it][ir]->Fill((Float_t) (pTJet - areaJet*fRhoMC[it][ir]), (Float_t) dfi); //Rongrong's analysis
1411  }
1412  }
1413  //-------------------------
1414 
1415  if(TMath::Abs(dphi) < fDphiCut) continue; //Dphi cut between trigger and assoc
1416 
1417  //Centrality, A, pTjet
1418  ftmpArray[0] = areaJet;
1419  for(Int_t ir=0; ir< kRho;ir++){
1420  ftmpArray[1] = pTJet - areaJet*fRhoMC[it][ir];
1421  for(Int_t ic=0; ic<2; ic++){
1422  if(ficb[ic]==-1) continue;
1423  fHJetSpecGen[ficb[ic]][it][ir]->Fill(ftmpArray[0],ftmpArray[1]);
1424  }
1425  }
1426  }//JET LOOP
1427  bFirstCycle = kFALSE;
1428  }//container exists
1429  }//TT exists
1430  }//TT loop
1431  }
1432 
1433  //++++++++++++++++++++++++++++++++++++++++++++++++++++
1434  if(fTypeOfAnal == kKine) return kTRUE;
1435  //++++++++++++++++++++++++++++++++++++++++++++++++++++
1436 
1437 
1438  for(Int_t it=0; it<kTT; it++){
1439  if((fTypeOfAnal==kEmb || fTypeOfAnal == kEmbSingl) && ((AliVParticle*) trackTT[it])){
1440  //delta pT using embedded pythia events
1441  //delta pT analyzed only in events with REAL EVENT TT present !!!!!!!!!!! (condition above)
1442  // PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDeltaPtJEmb.cxx
1443  //AliJetResponseMaker
1444  //get deltaPt from embedding
1445  //++++++++++++++++++++++++++++++++++++++++
1446  //ANALYZE DELTA PT FOR ALL EMB MC JETS
1447  //++++++++++++++++++++++++++++++++++++++++
1448 
1449  AliEmcalJet *jetEmb = NULL;
1450  Bool_t bEmbJetCloseToTT = kFALSE;
1451 
1452  if(jetContRec){
1453  jetContRec->ResetCurrentID();
1454  while((jetRec = jetContRec->GetNextAcceptJet())) { //loop over reconstructed jets
1455  if(!jetRec) continue;
1456  if(!IsSignalJetInAcceptance(jetRec,kTRUE)) continue; //apply cuts on eta, pT ,area
1457 
1458  //skip the jet that contains TT
1459  bEmbJetCloseToTT = kFALSE;
1460 
1461  if(jetRec->Pt() > trackTT[it]->Pt()*0.5){ // if jet contains TT it has to have at least pT of TT
1462 
1463  for(Int_t iq=0; iq < jetRec->GetNumberOfTracks(); iq++) {
1464  constTrackRec = static_cast<AliVParticle*> (jetRec->TrackAt(iq,trkContRec->GetArray())); //matched rec and emb tracks
1465  if(!constTrackRec) continue;
1466  if(constTrackRec != ((AliVParticle*) trackTT[it])) continue;
1467  //if(fDebug>21) Printf("EMB FIND TT COMPARE TRACK PT %f %f", constTrack->Pt(), triggerHadron->Pt());
1468  bEmbJetCloseToTT = kTRUE;
1469  break;
1470  }
1471  }
1472 
1473  if(bEmbJetCloseToTT) continue;//skip the jet that contains TT track
1474 
1475  jetEmb = NULL;
1476  jetEmb = jetRec->ClosestJet();
1477  if(!jetEmb){
1478  //if(fDebug>21) Printf("embedded jet does not exists, returning");
1479  continue;
1480  }
1481  if(jetEmb->Pt()<1e-3){
1482  //if(fDebug>21) Printf("SKIP MATCH WITH EMBEDDED GHOST JET");
1483  continue;
1484  }
1485  //if((fTypeOfAnal==kEmb) && (jetEmb->Pt()<6.0)) continue; // some hard cut on pT of emb jet ???
1486 
1487  if(!IsSignalJetInAcceptance(jetEmb,kTRUE)) continue; //apply cuts on eta, pT ,area on the embedded jet
1488 
1489  //Check fraction of tracks from generator level jet in rec level jet
1490  //At least 50% of embedded jet has to be in the reconstructed jet
1491  if(GetFractionSharedPt(jetRec, jetContRec, jetEmb, jetContGen) < fMinFractionShared) continue;
1492 
1493  for(Int_t ir=0; ir < kRho-1; ir++){
1494  for(Int_t ic=0; ic<2; ic++){
1495  if(ficb[ic]==-1) continue;
1496  //1 Dim distribution
1497  fhDeltaPtEmb[ficb[ic]][it][ir]->Fill(
1498  jetRec->Pt() - jetRec->Area() * fRhoRec[it][ir] - jetEmb->Pt());
1499 
1500  //2 Dim distribution
1501  fhDeltaPtEmb2D[ficb[ic]][it][ir]->Fill(jetEmb->Pt(),
1502  jetRec->Pt() - jetRec->Area() * fRhoRec[it][ir] - jetEmb->Pt());
1503 
1504 
1505  if(trackTT[it]){ //embedded track/jet is perp to TT
1506  dphi = TMath::Abs(RelativePhi(trackTT[it]->Phi(), jetEmb->Phi()));
1507  if(TMath::Pi()/4 < dphi && dphi < 3*TMath::Pi()/4){
1508  //1 Dim distribution
1509  fhDeltaPtEmbPerp[ficb[ic]][it][ir]->Fill(
1510  jetRec->Pt() - jetRec->Area() * fRhoRec[it][ir] - jetEmb->Pt());
1511 
1512  //2 Dim distribution
1513  fhDeltaPtEmbPerp2D[ficb[ic]][it][ir]->Fill(jetEmb->Pt(),
1514  jetRec->Pt() - jetRec->Area() * fRhoRec[it][ir] - jetEmb->Pt());
1515 
1516  }else if(dphi >= 3*TMath::Pi()/4 ){
1517  //embedded track back-to-back w.r.t the TT
1518  //1 Dim distribution
1519  fhDeltaPtEmbBc2Bc[ficb[ic]][it][ir]->Fill(
1520  jetRec->Pt() - jetRec->Area() * fRhoRec[it][ir] - jetEmb->Pt());
1521 
1522  //2 Dim distribution
1523  fhDeltaPtEmbBc2Bc2D[ficb[ic]][it][ir]->Fill(jetEmb->Pt(),
1524  jetRec->Pt() - jetRec->Area() * fRhoRec[it][ir] - jetEmb->Pt());
1525  }
1526  }
1527  }
1528  }
1529  }
1530  }
1531  }
1532  }//end of embedding
1533 
1534 
1535  // RECONSTRUCTED DATA ANALYSIS
1536 
1537  for(Int_t ir=0; ir < kRho-1; ir++){
1538  for(Int_t ic=0; ic<2; ic++){
1539  if(ficb[ic]==-1) continue;
1540 
1541  fhRhoIncl[ficb[ic]][ir]->Fill((Float_t) fRhoRec[kRef][ir]);
1542  }
1543  }
1544  for(Int_t it=0; it<kTT;it++){
1545  if(ntriggers[it]>0){
1546 
1547  for(Int_t ir=0; ir < kRho-1; ir++){
1548  //Estimate UE density in events with TT
1549  //Fill once per event
1550  for(Int_t ic=0; ic<2; ic++){
1551  if(ficb[ic]==-1) continue;
1552 
1553  fhRhoTT[ficb[ic]][it][ir]->Fill( (Float_t) fRhoRec[it][ir]);
1554  }
1555  }
1556  }
1557  }
1558  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1559  if(fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl) return kTRUE;
1560  //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1561 
1562 
1563  // CALCULATE DELTA PT IN RECONSTRUCTED DATA WITH RANDOM CONES
1564  Double_t deltapt[kRho-1], phiTT = 0., etaTT = -1000.;
1565  //Double_t ncoll = -1.0;
1566  //if(centralityPercentile>=0.) ncoll = GetNcoll(centralityPercentile);
1567  Bool_t bRecJetCloseToTT = kFALSE;
1568  //Exclude region around TT from delta pt calculation
1569 
1570  for(Int_t it=0; it<kTT;it++){
1571  bRecJetCloseToTT = kFALSE;
1572  phiTT = 0.;
1573  etaTT = -1000.;
1574 
1575  if(trackTT[it]){ //get phi and eta of the TT or the reconstructed jet that contains TT
1576  phiTT = trackTT[it]->Phi(); // TT is proxy for jet
1577  etaTT = trackTT[it]->Eta();
1578 
1579  if(jetContRec){ //find jet that contains TT if it exists
1580  jetContRec->ResetCurrentID();
1581  while((jetRec = jetContRec->GetNextAcceptJet())) { //loop over reconstructed jets
1582  if(!jetRec) continue;
1583  if(!IsSignalJetInAcceptance(jetRec,kTRUE)) continue;
1584 
1585  //skip the jet that contains TT
1586  bRecJetCloseToTT = kFALSE;
1587 
1588  if(jetRec->Pt() > trackTT[it]->Pt()*0.5){ // if jet contains TT it has to have at leat pT of TT
1589  for(Int_t iq=0; iq < jetRec->GetNumberOfTracks(); iq++){
1590  constTrackRec = (AliVParticle*) (jetRec->TrackAt(iq,trkContRec->GetArray())); //matched rec and emb tracks
1591  if(!constTrackRec) continue;
1592  if(constTrackRec != ((AliVParticle*) trackTT[it])) continue;
1593  phiTT = jetRec->Phi();
1594  etaTT = jetRec->Eta();
1595  bRecJetCloseToTT = kTRUE;
1596  break;
1597  }
1598  }
1599  if(bRecJetCloseToTT) break;
1600  }
1601  }
1602  }
1603 
1604  for(Int_t irc=0; irc<fNofRandomCones; irc++){
1605 
1606  //generate certain number of random cones per event
1607  GetDeltaPt(kRho-1, (fRhoRec[it]), &deltapt[0], phiTT, etaTT, trkContRec, kFALSE);// 1.0/ncoll);//FK//????? prob exlude RC ??? 1/Ncoll
1608 
1609 
1610  for(Int_t ir=0; ir < kRho-1; ir++){
1611  for(Int_t ic=0; ic<2; ic++){
1612  if(ficb[ic]==-1) continue;
1613  //fill delta pt histograms in inclusive events
1614  if(it==kRef) fhDeltaPtIncl[ficb[ic]][ir]->Fill(deltapt[ir]);
1615 
1616  if(ntriggers[it]>0){
1617  //fill delta pt histograms in events with TT (trigger track)
1618  fhDeltaPt[ficb[ic]][it][ir]->Fill( deltapt[ir]);
1619  }
1620  }
1621  }
1622  }
1623  }//trigger loop
1624  // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1625  // Track Multiplicity, Track momentum smearing
1626  Double_t xyz[50];
1627  Double_t pxpypz[50];
1628  Double_t cv[21];
1629  Int_t itrkq;
1630  AliAODTrack *atrk=NULL;
1631 
1632  if(trkContRec){
1633 
1634  Int_t mult = 0;
1635  trkContRec->ResetCurrentID();
1636  while((atrk = (AliAODTrack*) (trkContRec->GetNextAcceptParticle()))){
1637  if(!atrk) continue;
1638  if(!IsTrackInAcceptance((AliVParticle*) atrk, kFALSE)) continue; //reconstructed level tracks
1639  mult++;
1640 
1641  dphi = atrk->Phi();
1642  if(dphi < 0) dphi+= 2*TMath::Pi();
1643  if(dphi > 2*TMath::Pi()) dphi-= 2*TMath::Pi();
1644 
1645  if(atrk->IsGlobalConstrained()){
1646  fhTrackPhiCG->Fill(atrk->Pt(), atrk->Phi()); //global constrained
1647  }else{
1648  fhTrackPhiTPCG->Fill(atrk->Pt(), atrk->Phi()); //complementary
1649  }
1650 
1651  itrkq = (atrk->Charge()<0) ? 0 : 1;
1652 
1653  fhInvPtQVsPhi[itrkq]->Fill(atrk->Phi(), 1.0/atrk->Pt());
1654  fhInvPtQVsEta[itrkq]->Fill(atrk->Eta(), 1.0/atrk->Pt());
1655 
1656  if(atrk->Eta()>0){
1657  fhInvPtQVsPhiASide[itrkq]->Fill(atrk->Phi(), 1.0/atrk->Pt());
1658  }else{
1659  fhInvPtQVsPhiCSide[itrkq]->Fill(atrk->Phi(), 1.0/atrk->Pt());
1660  }
1661 
1662  //get sigma pT / pT
1663  //Taken from AliEMCalTriggerExtraCuts::CalculateTPCTrackLength
1664  memset(cv, 0, sizeof(Double_t) * 21); //cleanup arrays
1665  memset(pxpypz, 0, sizeof(Double_t) * 50);
1666  memset(xyz, 0, sizeof(Double_t) * 50);
1667  atrk->GetXYZ(xyz);
1668  atrk->GetPxPyPz(pxpypz);
1669  atrk->GetCovarianceXYZPxPyPz(cv);
1670 
1671  AliExternalTrackParam par(xyz, pxpypz, cv, atrk->Charge());
1672  fhSigmaPtOverPtVsPt[itrkq]->Fill(atrk->Pt(), TMath::Abs(sqrt(par.GetSigma1Pt2())/par.GetSigned1Pt()));
1673 
1674  //XXX DCA distributions
1675  fhDCAinXVsPt->Fill(atrk->Pt(), atrk->XAtDCA());
1676  fhDCAinYVsPt->Fill(atrk->Pt(), atrk->YAtDCA());
1677 
1678  if(fTypeOfAnal == kEff){
1679  lab = TMath::Abs(atrk->GetLabel());
1680 
1681  if(lab < iNTracksMC){
1682  particleMC = (AliAODMCParticle*) arrayMC->At(lab);
1683  if(particleMC){
1684  if(particleMC->IsPhysicalPrimary()){
1685  fhDCAinXVsPtNonStrange->Fill(atrk->Pt(), atrk->XAtDCA());
1686  fhDCAinYVsPtNonStrange->Fill(atrk->Pt(), atrk->YAtDCA());
1687  }else{
1688  fhDCAinXVsPtStrange->Fill(atrk->Pt(), atrk->XAtDCA());
1689  fhDCAinYVsPtStrange->Fill(atrk->Pt(), atrk->YAtDCA());
1690  }
1691  }
1692  }
1693  }
1694  }
1695 
1696  for(Int_t ic=0; ic<2; ic++){
1697  if(ficb[ic]==-1) continue;
1698  //MB or CENT biases
1699  fhTrackMultiplicity[ficb[ic]]->Fill(mult);
1700 
1701  ftmpArrayX[0] = energyZdcNA;
1702  ftmpArrayX[1] = multVzero;
1703  ftmpArrayX[2] = mult;
1704  fhZNAVzeroATrack[ficb[ic]]->Fill(ftmpArrayX);
1705 
1706  for(Int_t it=0; it<kTT;it++){
1707  if(trackTT[it]){//TT biased + TT&&CENT biased distributions
1708  fhTrackMultiplicityTT[ficb[ic]][it]->Fill(mult);
1709  fhVzeroATotMultTT[ficb[ic]][it]->Fill(multVzero);
1710  fhZNAEnergyTT[ficb[ic]][it]->Fill(energyZdcNA);
1711  fhZNAVzeroATrackTT[ficb[ic]][it]->Fill(ftmpArrayX);
1712  }
1713  }
1714  }
1715 
1716 
1717  for(Int_t it=0; it<kTT;it++){
1718  if(trackTT[it]){ //TT biased distributions
1719  fhCentralityTT[it]->Fill((Float_t) centralityPercentile);
1720  }
1721  }
1722  }
1723  // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1724  // INCLUSIVE JETS RECONSTRUCTED DATA
1725  if(jetContRec){
1726  jetContRec->ResetCurrentID();
1727  while((jetRec = jetContRec->GetNextAcceptJet())) {
1728  if(!jetRec){
1729  AliError(Form("%s: Could not receive jet", GetName()));
1730  continue;
1731  }
1732  if(!IsSignalJetInAcceptance(jetRec,kTRUE)) continue;
1733 
1734  pTJet = jetRec->Pt();
1735 
1736  fhJetPhiIncl->Fill( pTJet, RelativePhi(jetRec->Phi(),0.0));
1737  fhJetEtaIncl->Fill( pTJet, jetRec->Eta());
1738  }
1739  }
1740  // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1741  // H+JET IN RECONSTRUCTED DATA
1742 
1743  for(Int_t it=0; it<kTT;it++){
1744  if(ntriggers[it]>0 && ((AliVParticle*) trackTT[it])){
1745 
1746  bFirstCycle=kTRUE;
1747  if(trkContRec){
1748 
1749  AliVParticle* triggerHadron = (AliVParticle*) trackTT[it];// (fTrigTracks[it]);
1750  if(!triggerHadron) continue;
1751 
1752  for(Int_t ic=0; ic<2; ic++){
1753  if(ficb[ic]==-1) continue;
1754 
1755  fh1Ntriggers[ficb[ic]][it]->Fill((Float_t) triggerHadron->Pt()); //trigger p
1756  }
1757 
1758  //JET LOOP
1759  if(jetContRec){
1760  jetContRec->ResetCurrentID();
1761  while((jetRec = jetContRec->GetNextAcceptJet())) {
1762  if(!jetRec){
1763  AliError(Form("%s: Could not receive jet", GetName()));
1764  continue;
1765  }
1766  if(!IsSignalJetInAcceptance(jetRec,kTRUE)) continue;
1767 
1768  areaJet = jetRec->Area();
1769  pTJet = jetRec->Pt();
1770 
1771  if(bFirstCycle && it==kRef){
1772  for(Int_t ic=0; ic<2; ic++){
1773  if(ficb[ic]==-1) continue;
1774  fhJetPhi[ficb[ic]]->Fill( pTJet, RelativePhi(jetRec->Phi(),0.0));
1775  fhJetEta[ficb[ic]]->Fill( pTJet, jetRec->Eta());
1776  }
1777  }
1778 
1779  dphi = RelativePhi(triggerHadron->Phi(), jetRec->Phi());
1780 
1781  dfi = dphi; //-0.5*pi to 1.5*Pi
1782  if(dfi < -0.5*TMath::Pi()) dfi += TMath::TwoPi();
1783  if(dfi > 1.5*TMath::Pi()) dfi -= TMath::TwoPi();
1784 
1785  for(Int_t ic=0; ic<2; ic++){
1786  if(ficb[ic]==-1) continue;
1787 
1788  for(Int_t ir=0; ir< kRho;ir++){
1789  fhDphiTriggerJet[ficb[ic]][it][ir]->Fill((Float_t) (pTJet - areaJet*fRhoRec[it][ir]), (Float_t) dfi); //Rongrong's analysis
1790  }
1791  }
1792  //-------------------------
1793 
1794  if(TMath::Abs(dphi) < fDphiCut) continue; //Dphi cut between trigger and assoc
1795  if(it==kRef) fhDphiTriggerJetAccept->Fill(dfi); //Accepted
1796 
1797  //Centrality, A, pTjet
1798  ftmpArray[0] = areaJet;
1799  for(Int_t ir=0; ir< kRho;ir++){
1800  ftmpArray[1] = pTJet - areaJet*fRhoRec[it][ir];
1801 
1802  for(Int_t ic=0; ic<2; ic++){
1803  if(ficb[ic]==-1) continue;
1804 
1805  fHJetSpec[ficb[ic]][it][ir]->Fill(ftmpArray[0],ftmpArray[1]);
1806 
1807  if(ir<kRho-1){
1808  fARhoTT[ficb[ic]][it][ir]->Fill((Float_t) (areaJet*fRhoRec[it][ir]));
1809  }
1810  }
1811  }
1812  }//JET LOOP
1813 
1814  bFirstCycle = kFALSE;
1815  }//container exists
1816  }
1817  }
1818  } //TT loop
1819 
1820  return kTRUE;
1821 }
1822 
1823 //________________________________________________________________________
1824 /*Double_t AliAnalysisTaskHJetSpectra::GetNcoll(Double_t centr){
1825  //Get Ncoll for given centrality
1826  if(fCollisionSystem == kpPb){
1827  //convert centrality percentle to Ncoll (2014-Sep-26-analysis_note-AnalysisNoteDijetspPb.pdf table 1)
1828  // PLATI JEN PRO pA !!!!!!!!!! CO PP?
1829  if(centr < 0.0 || centr > 100.) return -1;
1830 
1831  if(centr < 5.0) return 14.7; //0-5%
1832  else if(centr < 10) return 13.0; //5-10%
1833  else if(centr < 20) return 11.7; //10-20%
1834  else if(centr < 40) return 9.38; //20-40%
1835  else if(centr < 60) return 6.49; //40-60%
1836  else if(centr < 80) return 3.96; //60-80%
1837  else return 1.52; //80-100%
1838  }
1839 
1840 
1841  return -1.; //pp, pbpb
1842 
1843 }*/
1844 //________________________________________________________________________
1845 void AliAnalysisTaskHJetSpectra::Terminate(Option_t *){
1846  //Treminate
1847  PostData(1, fOutput);
1848 
1849  // Mandatory
1850  fOutput = dynamic_cast<AliEmcalList*> (GetOutputData(1)); // '1' refers to the output slot
1851  if(!fOutput) {
1852  printf("ERROR: Output list not available\n");
1853  return;
1854  }
1855 }
1856 
1857 //________________________________________________________________________
1858 AliAnalysisTaskHJetSpectra::~AliAnalysisTaskHJetSpectra(){
1859  // Destructor. Clean-up the output list, but not the histograms that are put inside
1860  // (the list is owner and will clean-up these histograms). Protect in PROOF case.
1861  if(fOutput && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
1862  delete fOutput;
1863  }
1864  delete fRandom;
1865  delete fHelperClass;
1866 
1867 }
1868 
1869 //________________________________________________________________________
1870 void AliAnalysisTaskHJetSpectra::UserCreateOutputObjects(){
1871  // called once to create user defined output objects like histograms, plots etc.
1872  // and to put it on the output list.
1873  // Note: Saving to file with e.g. OpenFile(0) is must be before creating other objects.
1874  //fOutput TList defined in the mother class
1875  AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
1876 
1877  //Label of background
1878  TString bgtype[]={"Perp","CMS","KT","Zero"};
1879 
1880  fRandom = new TRandom3(0);
1881 
1882  Bool_t oldStatus = TH1::AddDirectoryStatus();
1883  TH1::AddDirectory(kFALSE);
1884  TString name;
1885 
1886  Bool_t bHistRec = (fTypeOfAnal != kEmb && fTypeOfAnal != kEmbSingl && fTypeOfAnal != kKine);
1887  Bool_t bNotKine = (fTypeOfAnal != kKine);
1888 
1889  Int_t icmax = (fTypeOfData == kPythia || fTypeOfData == kDmpjet ||fTypeOfAnal == kKine) ? 1 : kCAll;
1890  //__________________________________________________________
1891  // Event statistics
1892  fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 7, -0.5, 6.5);
1893  fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
1894  fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
1895  fHistEvtSelection->GetXaxis()->SetBinLabel(3,"pile up (rejected)");
1896  fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
1897  fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
1898  fHistEvtSelection->GetXaxis()->SetBinLabel(6,"sig TT (rejected)");
1899  fHistEvtSelection->GetXaxis()->SetBinLabel(7,"sig TT (accepted)");
1900 
1901  fOutput->Add(fHistEvtSelection);
1902  //___________________________________________________________
1903  // Hard trigger counter
1904  for(Int_t it =0; it<kTT; it++){
1905  for(Int_t ic =0; ic<icmax; ic++){
1906  name = (ic==0) ? Form("fh1NtriggersMB") :
1907  Form("fh1Ntriggers%d%d",TMath::Nint(fCentralityBins[ic-1]),TMath::Nint(fCentralityBins[ic]));
1908  name.Append(Form("TT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
1909  fh1Ntriggers[ic][it] = new TH1D(name.Data(),"# of triggers",50,0.0,50.0);
1910  if(bHistRec) fOutput->Add((TH1D*)fh1Ntriggers[ic][it]);
1911 
1912  name = (ic==0) ? Form("fh1TriggerMultMB") :
1913  Form("fh1TriggerMult%d%d",TMath::Nint(fCentralityBins[ic-1]),TMath::Nint(fCentralityBins[ic]));
1914  name.Append(Form("TT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
1915  fh1TriggerMult[ic][it] = new TH1D(name.Data(),"# of triggers",50,0.0,50.0);
1916  if(bHistRec) fOutput->Add((TH1D*)fh1TriggerMult[ic][it]);
1917  }
1918  }
1919  //___________________________________________________________
1920  // trigger associated jet spectra (jet pT not corrected for UE)
1921  Int_t bw = (fUseDoubleBinPrecision==0) ? 1 : 2; //make larger bin width
1922 
1923  //jet associated to given TT
1924  //A, pTjet
1925  const Int_t dimSpec = 2;
1926  const Int_t nBinsSpec[dimSpec] = { 50, bw*160};
1927  const Double_t lowBinSpec[dimSpec] = { 0.0, -20.0};
1928  const Double_t hiBinSpec[dimSpec] = { 2.0, 300.0};
1929 
1930  for(Int_t ic =0; ic<icmax; ic++){
1931  for(Int_t it =0; it<kTT; it++){
1932  for(Int_t ir=0; ir< kRho; ir++){
1933  name = (ic==0) ? "fHJetSpecMB" :
1934  Form("fHJetSpec%d%d",TMath::Nint(fCentralityBins[ic-1]),
1935  TMath::Nint(fCentralityBins[ic]));
1936  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
1937 
1938  fHJetSpec[ic][it][ir] = new TH2D(
1939  name.Data(),
1940  Form("Recoil jet spectrum TT%d%d [A,pTjet-A*rho%s]",
1941  TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()),
1942  nBinsSpec[0], lowBinSpec[0], hiBinSpec[0], nBinsSpec[1], lowBinSpec[1],hiBinSpec[1]);
1943  if(bHistRec) fOutput->Add((TH2D*) fHJetSpec[ic][it][ir]);
1944  }
1945  }
1946  }
1947  //____________________________________________________________________
1948  //UE from cell median [Centrality, rho, pTUe ]
1949 
1950  for(Int_t ic =0; ic<icmax; ic++){
1951  for(Int_t ir=0; ir< kRho-1; ir++){ //Skip Zero bg
1952  for(Int_t it=0; it<kTT;it++){
1953  name = (ic==0) ? "fhRhoMB" :
1954  Form("fhRho%d%d",TMath::Nint(fCentralityBins[ic-1]),
1955  TMath::Nint(fCentralityBins[ic]));
1956 
1957  name.Append(Form("TT%d%d%s",TMath::Nint(fTTlow[it]),TMath::Nint(fTThigh[it]),bgtype[ir].Data()));
1958 
1959  //rho in events with TT
1960  fhRhoTT[ic][it][ir] = new TH1F(name.Data(),
1961  Form("Rho%s",bgtype[ir].Data()),80, 0.0, 40.0);
1962  if(bNotKine) fOutput->Add((TH1F*) fhRhoTT[ic][it][ir]);
1963 
1964 
1965  // rho times area in events with TT
1966  name = (ic==0) ? "fARhoMB" :
1967  Form("fARho%d%d",TMath::Nint(fCentralityBins[ic-1]),
1968  TMath::Nint(fCentralityBins[ic]));
1969 
1970  name.Append(Form("TT%d%d%s",TMath::Nint(fTTlow[it]),TMath::Nint(fTThigh[it]),bgtype[ir].Data()));
1971  fARhoTT[ic][it][ir] = new TH1F(name.Data(),
1972  Form("Area times rho %s",bgtype[ir].Data()),80, 0.0, 40.0);
1973  if(bHistRec) fOutput->Add((TH1F*) fARhoTT[ic][it][ir]);
1974  }
1975  //rho in inclusive events
1976  name = (ic==0) ? Form("fhRhoInclMB%s",bgtype[ir].Data()) :
1977  Form("fhRhoIncl%d%d%s",TMath::Nint(fCentralityBins[ic-1]),
1978  TMath::Nint(fCentralityBins[ic]),bgtype[ir].Data());
1979 
1980  fhRhoIncl[ic][ir] = (TH1F*) fhRhoTT[ic][kRef][ir]->Clone(name.Data());
1981  if(bNotKine) fOutput->Add((TH1F*) fhRhoIncl[ic][ir]);
1982  }
1983  }
1984  //_______________________________________________________________________
1985  // Delta pt distributions
1986  for(Int_t it =0; it< kTT; it++){
1987  for(Int_t ic =0; ic<icmax; ic++){
1988  for(Int_t ir=0; ir< kRho-1; ir++){
1989  //events with TT tracks
1990  name = (ic==0) ? "fhDeltaPtMB" :
1991  Form("fhDeltaPt%d%d",TMath::Nint(fCentralityBins[ic-1]),
1992  TMath::Nint(fCentralityBins[ic]));
1993 
1994  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
1995  fhDeltaPt[ic][it][ir] = new TH1D(name.Data(),
1996  Form("DeltaPt%s",bgtype[ir].Data()), 150, -50, 100);
1997  if(bHistRec) fOutput->Add((TH1D*) fhDeltaPt[ic][it][ir]);
1998 
1999  if(it==kRef){
2000  //inclusive events
2001  name = (ic==0) ? "fhDeltaPtInclMB" :
2002  Form("fhDeltaPtIncl%d%d",TMath::Nint(fCentralityBins[ic-1]),
2003  TMath::Nint(fCentralityBins[ic]));
2004 
2005  name.Append(Form("%s", bgtype[ir].Data()));
2006  fhDeltaPtIncl[ic][ir] = (TH1D*) fhDeltaPt[ic][it][ir]->Clone(name.Data());
2007  if(bHistRec) fOutput->Add((TH1D*) fhDeltaPtIncl[ic][ir]);
2008  }
2009 
2010  if(fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
2011  //Embedded PYTHIA jets
2012  //1D
2013  name = (ic==0) ? "fhDeltaPtEmbMB" :
2014  Form("fhDeltaPtEmb%d%d",TMath::Nint(fCentralityBins[ic-1]),
2015  TMath::Nint(fCentralityBins[ic]));
2016 
2017  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
2018  fhDeltaPtEmb[ic][it][ir] = (TH1D*) fhDeltaPt[ic][it][ir]->Clone(name.Data());
2019  fOutput->Add((TH1D*) fhDeltaPtEmb[ic][it][ir]);
2020 
2021  //1D perp to TT
2022  name = (ic==0) ? "fhDeltaPtEmbPerpMB" :
2023  Form("fhDeltaPtEmbPerp%d%d",TMath::Nint(fCentralityBins[ic-1]),
2024  TMath::Nint(fCentralityBins[ic]));
2025 
2026  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
2027  fhDeltaPtEmbPerp[ic][it][ir] = (TH1D*) fhDeltaPt[ic][it][ir]->Clone(name.Data());
2028  fOutput->Add((TH1D*) fhDeltaPtEmbPerp[ic][it][ir]);
2029 
2030  //1D back-to back to TT
2031  name = (ic==0) ? "fhDeltaPtEmbBc2BcMB" :
2032  Form("fhDeltaPtEmbBc2Bc%d%d",TMath::Nint(fCentralityBins[ic-1]),
2033  TMath::Nint(fCentralityBins[ic]));
2034 
2035  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
2036  fhDeltaPtEmbBc2Bc[ic][it][ir]= (TH1D*) fhDeltaPt[ic][it][ir]->Clone(name.Data());
2037  fOutput->Add((TH1D*) fhDeltaPtEmbBc2Bc[ic][it][ir]);
2038 
2039  //2D
2040  name = (ic==0) ? "fhDeltaPtEmb2DMB" :
2041  Form("fhDeltaPtEmb2D%d%d",TMath::Nint(fCentralityBins[ic-1]),
2042  TMath::Nint(fCentralityBins[ic]));
2043 
2044  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
2045  fhDeltaPtEmb2D[ic][it][ir] = new TH2D(name.Data(),
2046  Form("fhDeltaPtEmb2D%s",bgtype[ir].Data()), 125,0, 250, 150, -50, 100);
2047  fOutput->Add((TH2D*)fhDeltaPtEmb2D[ic][it][ir]);
2048 
2049 
2050  //2D perp to TT
2051  name = (ic==0) ? "fhDeltaPtEmbPerp2DMB" :
2052  Form("fhDeltaPtEmbPerp2D%d%d",TMath::Nint(fCentralityBins[ic-1]),
2053  TMath::Nint(fCentralityBins[ic]));
2054 
2055  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
2056  fhDeltaPtEmbPerp2D[ic][it][ir] = (TH2D*) fhDeltaPtEmb2D[ic][it][ir]->Clone(name.Data());
2057  fOutput->Add((TH2D*)fhDeltaPtEmbPerp2D[ic][it][ir]);
2058 
2059 
2060  //2D back-to-back p to TT
2061  name = (ic==0) ? "fhDeltaPtEmbBc2Bc2DMB" :
2062  Form("fhDeltaPtEmbBc2Bc2D%d%d",TMath::Nint(fCentralityBins[ic-1]),
2063  TMath::Nint(fCentralityBins[ic]));
2064 
2065  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
2066  fhDeltaPtEmbBc2Bc2D[ic][it][ir]=(TH2D*) fhDeltaPtEmb2D[ic][it][ir]->Clone(name.Data());
2067  fOutput->Add((TH2D*)fhDeltaPtEmbBc2Bc2D[ic][it][ir]);
2068 
2069  }
2070  }
2071  }
2072  }
2073 
2074  //_______________________________________________________________________
2075  //KT jets area versus PT
2076  Double_t binsPt [] = {0, 0.05,0.15, 0.5, 1., 1.5, 2., 2.5, 3., 3.5, 4., 4.5, 5., 5.5, 6.,6.5, 7., 8., 9., 10.,15., 20., 50., 100.};
2077  Int_t nbinsPt = sizeof(binsPt)/sizeof(Double_t)-1;
2078 
2079  fhKTAreaPt = new TH2F("fhKTAreaPt","KT jet Area vs Pt",nbinsPt, binsPt, 50,0,2);
2080  if(fTypeOfAnal == kRec) fOutput->Add(fhKTAreaPt);
2081 
2082  //_______________________________________________________________________
2083  //inclusive azimuthal and pseudorapidity histograms
2084  for(Int_t ic =0; ic<icmax; ic++){
2085  name = (ic==0) ? Form("fhJetPhiMB") :
2086  Form("fhJetPhi%d%d",TMath::Nint(fCentralityBins[ic-1]),
2087  TMath::Nint(fCentralityBins[ic]));
2088 
2089  fhJetPhi[ic] = new TH2F(name.Data(),"Azim dist jets vs pTjet", 50, 0, 100, 50,-TMath::Pi(),TMath::Pi());
2090  if(bHistRec) fOutput->Add((TH2F*)fhJetPhi[ic]);
2091 
2092  //-------------------------
2093  name = (ic==0) ? Form("fhTrackPhiMB") :
2094  Form("fhTrackPhi%d%d",TMath::Nint(fCentralityBins[ic-1]),
2095  TMath::Nint(fCentralityBins[ic]));
2096 
2097  fhTrackPhi[ic] = new TH2F(name.Data(),"azim dist trig had vs pT,trk", 50, 0, 50, 50,-TMath::Pi(),TMath::Pi());
2098  if(bNotKine) fOutput->Add((TH2F*)fhTrackPhi[ic]);
2099  //-------------------------
2100  name = (ic==0) ? Form("fhJetEtaMB") :
2101  Form("fhJetEta%d%d",TMath::Nint(fCentralityBins[ic-1]),
2102  TMath::Nint(fCentralityBins[ic]));
2103 
2104  fhJetEta[ic] = new TH2F(name.Data(),"Eta dist jets vs pTjet", 50,0, 100, 40,-0.9,0.9);
2105  if(bHistRec) fOutput->Add((TH2F*)fhJetEta[ic]);
2106  //-------------------------
2107  name = (ic==0) ? Form("fhTrackEtaMB") :
2108  Form("fhTrackEta%d%d",TMath::Nint(fCentralityBins[ic-1]),
2109  TMath::Nint(fCentralityBins[ic]));
2110 
2111  fhTrackEta[ic] = new TH2F(name.Data(),"Eta dist trig had vs pT,trk", 50, 0, 50, 40,-0.9,0.9);
2112  if(bNotKine) fOutput->Add((TH2F*) fhTrackEta[ic]);
2113  //-------------------------
2114  name = (ic==0) ? Form("fhTrackPtMB") :
2115  Form("fhTrackPt%d%d",TMath::Nint(fCentralityBins[ic-1]),
2116  TMath::Nint(fCentralityBins[ic]));
2117 
2118  fhTrackPt[ic] = new TH1F(name.Data(),"pT,trk ", 50, 0, 50);
2119  if(bNotKine) fOutput->Add((TH1F*) fhTrackPt[ic]);
2120  }
2121  //-------------------------
2122  fhVertexZ = new TH1F("fhVertexZ","z vertex",40,-20,20);
2123  if(bNotKine) fOutput->Add(fhVertexZ);
2124  //-------------------------
2125  fhVertexXAccept = new TH1F("fhVertexXAccept","vertex after cut",600,-3,3);
2126  if(bNotKine) fOutput->Add(fhVertexXAccept);
2127 
2128  fhVertexYAccept = (TH1F*) fhVertexXAccept->Clone("fhVertexYAccept");
2129  if(bNotKine) fOutput->Add(fhVertexYAccept);
2130 
2131  fhVertexZAccept = new TH1F("fhVertexZAccept","z vertex after cut",40,-20,20);
2132  if(bNotKine) fOutput->Add(fhVertexZAccept);
2133 
2134  fhVertexXAcceptTT = (TH1F*) fhVertexXAccept->Clone("fhVertexXAcceptTT");
2135  if(bNotKine) fOutput->Add(fhVertexXAcceptTT);
2136 
2137  fhVertexYAcceptTT = (TH1F*) fhVertexXAccept->Clone("fhVertexYAcceptTT");
2138  if(bNotKine) fOutput->Add(fhVertexYAcceptTT);
2139 
2140  fhVertexZAcceptTT = (TH1F*) fhVertexZAccept->Clone("fhVertexZAcceptTT");
2141  if(bNotKine) fOutput->Add(fhVertexZAcceptTT);
2142 
2143  //-------------------------
2144 
2145  if(fTypeOfData!=kReal){
2146  fhVertexZMC = new TH1F("fhVertexZMC","z vertex",40,-20,20);
2147  fOutput->Add(fhVertexZMC);
2148  //-------------------------
2149  fhVertexZAcceptMC = new TH1F("fhVertexZAcceptMC","z vertex after cut",40,-20,20);
2150  fOutput->Add(fhVertexZAcceptMC);
2151  }
2152  //-------------------------
2153  for(Int_t ic =0; ic<icmax; ic++){
2154  for(Int_t it =0; it<kTT; it++){
2155  for(Int_t ir=0; ir < kRho; ir++){ //Rongrong's analysis
2156  name = (ic==0) ?
2157  Form("fhDphiTriggerJetMB") :
2158  Form("fhDphiTriggerJet%d%d",TMath::Nint(fCentralityBins[ic-1]), TMath::Nint(fCentralityBins[ic]));
2159 
2160  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]) ,bgtype[ir].Data()));
2161 
2162  fhDphiTriggerJet[ic][it][ir] = new TH2F(name.Data(),"Deltaphi trig-jet",75,-50,100, 100, -0.5*TMath::Pi(),1.5*TMath::Pi());
2163  if(bHistRec) fOutput->Add((TH2F*) fhDphiTriggerJet[ic][it][ir]);
2164  }
2165  }
2166  }
2167  //-------------------------
2168 
2169  fhDphiTriggerJetAccept = new TH1F("fhDphiTriggerJetAccept","Deltaphi trig-jet after cut",50, -0.5*TMath::Pi(),1.5*TMath::Pi());
2170  if(bHistRec) fOutput->Add(fhDphiTriggerJetAccept);
2171 
2172  fhDphiTTTT[kRef] = new TH1F("fhDphiTTTT67","Deltaphi between multiple TT",50, -0.5*TMath::Pi(),1.5*TMath::Pi());
2173  if(bHistRec) fOutput->Add((TH1F*) fhDphiTTTT[kRef]);
2174 
2175  fhDphiTTTT[kSig] = (TH1F*) fhDphiTTTT[kRef]->Clone("fhDphiTTTT1250");
2176  if(bHistRec) fOutput->Add((TH1F*) fhDphiTTTT[kSig]);
2177 
2178  fhJetPhiIncl = new TH2F("fhJetPhiIncl","Azim dist jets vs pTjet", 50, 0, 100, 50,-TMath::Pi(),TMath::Pi());
2179  if(bHistRec) fOutput->Add((TH2F*) fhJetPhiIncl);
2180 
2181  fhJetEtaIncl = new TH2F("fhJetEtaIncl","Eta dist inclusive jets vs pTjet", 50,0, 100, 40,-0.9,0.9);
2182  if(bHistRec) fOutput->Add((TH2F*) fhJetEtaIncl);
2183 
2184  fhTrackPhiCG = new TH2F("fhTrackPhiCG","azim dist trig had vs pT,trk Glob Const", 50, 0, 50, 50,0,2*TMath::Pi());
2185  if(bNotKine) fOutput->Add((TH2F*) fhTrackPhiCG);
2186 
2187  fhTrackPhiTPCG = new TH2F("fhTrackPhiTPCG","azim dist trig had vs pT,trk TPC Const", 50, 0, 50, 50,0,2*TMath::Pi());
2188  if(bNotKine) fOutput->Add((TH2F*) fhTrackPhiTPCG);
2189 
2190 
2191  for(Int_t i=0; i<2; i++){
2192  name = (i==0) ? "Neg" : "Pos";
2193  fhInvPtQVsPhi[i] = new TH2D(Form("fhInvPtVsPhiQ%s", name.Data()),
2194  Form("%s track 1/pt versus track phi", name.Data()), 36, 0, 2*TMath::Pi(), 40, 0, 0.4);
2195  if(bNotKine) fOutput->Add((TH2D*) fhInvPtQVsPhi[i]);
2196 
2197  fhInvPtQVsEta[i] = new TH2D(Form("fhInvPtVsEtaQ%s", name.Data()),
2198  Form("%s track 1/pt versus track eta", name.Data()), 20, -0.9, 0.9, 40, 0, 0.4);
2199  if(bNotKine) fOutput->Add((TH2D*) fhInvPtQVsEta[i]);
2200 
2201  fhInvPtQVsPhiASide[i] = new TH2D(Form("fhInvPtVsPhiASideQ%s", name.Data()),
2202  Form("%s track 1/pt versus track phi", name.Data()), 36, 0, 2*TMath::Pi(), 40, 0, 0.4);
2203  if(bNotKine) fOutput->Add((TH2D*) fhInvPtQVsPhiASide[i]);
2204 
2205  fhInvPtQVsPhiCSide[i] = new TH2D(Form("fhInvPtVsPhiCSideQ%s", name.Data()),
2206  Form("%s track 1/pt versus track phi", name.Data()), 36, 0, 2*TMath::Pi(), 40, 0, 0.4);
2207  if(bNotKine) fOutput->Add((TH2D*) fhInvPtQVsPhiCSide[i]);
2208 
2209  fhSigmaPtOverPtVsPt[i] = new TH2D(Form("fhSigmaPtOverPtVsPtQ%s", name.Data()),
2210  Form("%s track sigma(1/pt)/ 1/pt vs pt", name.Data()), 100, 0, 100, 250, 0, 1);
2211  if(bNotKine) fOutput->Add((TH2D*) fhSigmaPtOverPtVsPt[i]);
2212  }
2213  //-------------------------
2214  for(Int_t ic =0; ic<icmax; ic++){
2215  name = (ic==0) ? Form("fhCentralityMB") :
2216  Form("fhCentrality%d%d",TMath::Nint(fCentralityBins[ic-1]),
2217  TMath::Nint(fCentralityBins[ic]));
2218 
2219  fhCentrality[ic] = new TH1F(name.Data(),"Centrality",20,0,100);
2220  if(bNotKine) fOutput->Add((TH1F*) fhCentrality[ic]);
2221  }
2222  //-------------------------
2223  fhCentralityV0M = new TH1F("hCentralityV0M","hCentralityV0M",20,0,100);
2224  if(bNotKine)fOutput->Add(fhCentralityV0M);
2225  //-------------------------
2226  fhCentralityV0A = new TH1F("hCentralityV0A","hCentralityV0A",20,0,100);
2227  if(bNotKine) fOutput->Add(fhCentralityV0A);
2228  //-------------------------
2229  fhCentralityV0C = new TH1F("hCentralityV0C","hCentralityV0C",20,0,100);
2230  if(bNotKine) fOutput->Add(fhCentralityV0C);
2231  //-------------------------
2232  fhCentralityZNA = new TH1F("hCentralityZNA","hCentralityZNA",20,0,100);
2233  if(bNotKine) fOutput->Add(fhCentralityZNA);
2234  //-------------------------
2235  for(Int_t it =0; it<kTT;it++){
2236  name = Form("fhCentralityTT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]));
2237  fhCentralityTT[it] = (TH1F*) fhCentralityV0M->Clone(name.Data());
2238  if(bHistRec) fOutput->Add(fhCentralityTT[it]);
2239  }
2240  //-----------------------------------------------------
2241 
2242  for(Int_t ic =0; ic<icmax; ic++){
2243  name = (ic==0) ? Form("fhVzeroATotMultMB") :
2244  Form("fhVzeroATotMult%d%d",TMath::Nint(fCentralityBins[ic-1]),
2245  TMath::Nint(fCentralityBins[ic]));
2246 
2247  // vzero multiplicity
2248  fhVzeroATotMult[ic] = new TH1F(name.Data(),"hVzeroATotMult",1000,0,1000);
2249  if(bHistRec) fOutput->Add((TH1F*) fhVzeroATotMult[ic]);
2250 
2251  for(Int_t it=0; it<kTT; it++){
2252  name = (ic==0) ? "fhVzeroATotMultMB" :
2253  Form("fhVzeroATotMult%d%d",TMath::Nint(fCentralityBins[ic-1]),
2254  TMath::Nint(fCentralityBins[ic]));
2255 
2256  name.Append(Form("TT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
2257  fhVzeroATotMultTT[ic][it] = (TH1F*) fhVzeroATotMult[ic]->Clone(name.Data());
2258  if(bHistRec) fOutput->Add((TH1F*) fhVzeroATotMultTT[ic][it]);
2259  }
2260  }
2261 
2262 
2263  //-----------------------------------------------------
2264  // ZDC ZNA energy
2265 
2266  for(Int_t ic =0; ic<icmax; ic++){
2267  name = (ic==0) ? Form("fhZNAEnergyMB") :
2268  Form("fhZNAEnergy%d%d",TMath::Nint(fCentralityBins[ic-1]),
2269  TMath::Nint(fCentralityBins[ic]));
2270 
2271 
2272  fhZNAEnergy[ic] = new TH1F(name.Data(),"fhZNAEnergy",1000,0,200);
2273  if(bHistRec) fOutput->Add((TH1F*)fhZNAEnergy[ic]);
2274 
2275  for(Int_t it=0; it<kTT; it++){
2276  name = (ic==0) ? Form("fhZNAEnergyMB") :
2277  Form("fhZNAEnergy%d%d",TMath::Nint(fCentralityBins[ic-1]),
2278  TMath::Nint(fCentralityBins[ic]));
2279 
2280  name.Append(Form("TT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
2281  fhZNAEnergyTT[ic][it] = (TH1F*) fhZNAEnergy[ic]->Clone(name.Data());
2282  if(bHistRec) fOutput->Add((TH1F*) fhZNAEnergyTT[ic][it]);
2283  }
2284  }
2285 
2286  //-----------------------------------------------------
2287  // track multiplicity
2288 
2289  for(Int_t ic =0; ic<icmax; ic++){
2290  name = (ic==0) ? Form("fhTrackMultiplicityMB") :
2291  Form("fhTrackMultiplicity%d%d",TMath::Nint(fCentralityBins[ic-1]),
2292  TMath::Nint(fCentralityBins[ic]));
2293 
2294  fhTrackMultiplicity[ic] = new TH1D(name.Data(),"fhTrackMultiplicity",1000,0,1000);
2295  if(bNotKine) fOutput->Add((TH1D*)fhTrackMultiplicity[ic]);
2296 
2297  for(Int_t it=0; it<kTT; it++){
2298  name = (ic==0) ? Form("fhTrackMultiplicityMB") :
2299  Form("fhTrackMultiplicity%d%d",TMath::Nint(fCentralityBins[ic-1]),
2300  TMath::Nint(fCentralityBins[ic]));
2301 
2302  name.Append(Form("TT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
2303  fhTrackMultiplicityTT[ic][it] = (TH1D*)fhTrackMultiplicity[ic]->Clone(name.Data());
2304  if(bNotKine) fOutput->Add(fhTrackMultiplicityTT[ic][it]);
2305  }
2306  }
2307 
2308  //-----------------------------------------------------
2309  // ZNA energy versus Vzero mult. versus track mult. in all events
2310 
2311  const Int_t dimZ = 3;
2312  const Int_t nBinsZ[dimZ] = { 100, 100, 25 };
2313  const Double_t lowBinZ[dimZ] = { 0.0, 0.0, 0.0};
2314  const Double_t hiBinZ[dimZ] = { 200.0, 1000.0, 250};
2315 
2316  for(Int_t ic =0; ic<icmax; ic++){
2317  name = (ic==0) ? Form("fhZNAVzeroATrackMB") :
2318  Form("fhZNAVzeroATrack%d%d",TMath::Nint(fCentralityBins[ic-1]),
2319  TMath::Nint(fCentralityBins[ic]));
2320 
2321  fhZNAVzeroATrack[ic] = new THnSparseF(
2322  name.Data(),
2323  Form("ZNA, V0A mult, track mult"),
2324  dimZ, nBinsZ,lowBinZ,hiBinZ);
2325  if(fTypeOfAnal == kRec) fOutput->Add((THnSparseF*) fhZNAVzeroATrack[ic]);
2326 
2327 
2328  //-------------
2329  for(Int_t it=0; it<kTT; it++){
2330  // ZNA energy versus Vzero mult. versus track mult. in events with TT
2331  name = (ic==0) ? Form("fhZNAVzeroATrackMB") :
2332  Form("fhZNAVzeroATrack%d%d",TMath::Nint(fCentralityBins[ic-1]),
2333  TMath::Nint(fCentralityBins[ic]));
2334 
2335  name.Append(Form("TT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
2336 
2337  fhZNAVzeroATrackTT[ic][it] = (THnSparseF*) fhZNAVzeroATrack[ic]->Clone(name.Data());
2338  if(fTypeOfAnal == kRec) fOutput->Add((THnSparseF*) fhZNAVzeroATrackTT[ic][it]);
2339  }
2340  }
2341 
2342  //-----------------------------------------------------
2343 
2344  /*
2345  if(fTypeOfAnal == kKine){
2346  fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
2347  fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
2348  fOutput->Add(fh1Xsec);
2349 
2350  fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
2351  fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
2352  fOutput->Add(fh1Trials);
2353 
2354  //fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",1000,0,1000);
2355  //fOutput->Add(fh1PtHard);
2356  }
2357 */
2358  if(fTypeOfData == kHijing){
2359  for(Int_t ic =0; ic<kCAll; ic++){
2360  name = (ic==0) ? Form("fhImpactParameterMB") :
2361  Form("fhImpactParameter%d%d",TMath::Nint(fCentralityBins[ic-1]),
2362  TMath::Nint(fCentralityBins[ic]));
2363 
2364  fhImpactParameter[ic] = new TH1D(name.Data(),"impact parameter distribution from HIJING",50,0,10);
2365  fOutput->Add((TH1D*) fhImpactParameter[ic]);
2366 
2367 
2368  for(Int_t it=0; it<kTT;it++){
2369  name = (ic==0) ? Form("fhImpactParameterMB") :
2370  Form("fhImpactParameter%d%d",TMath::Nint(fCentralityBins[ic-1]),
2371  TMath::Nint(fCentralityBins[ic]));
2372 
2373  name.Append(Form("TT%d%d", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
2374  fhImpactParameterTT[ic][it] = new TH1D(name.Data(),"b versus TT",50,0,10);
2375  fOutput->Add((TH1D*) fhImpactParameterTT[ic][it]);
2376  }
2377  }
2378  }
2379 
2380  Double_t bins [] = {0, 0.2,0.4,0.6, 0.8, 1., 1.2, 1.4, 1.6, 1.8, 2., 2.5, 3., 3.5, 4., 5., 6., 8., 10., 20., 50.};
2381  Int_t nbins = sizeof(bins)/sizeof(Double_t)-1;
2382 
2383  if(fTypeOfAnal== kEff ){
2384  for(Int_t ic =0; ic<icmax; ic++){
2385  for(Int_t ir=0; ir < kRho; ir++){
2386  name = (ic==0) ? Form("fhJetPtGenMB%s",bgtype[ir].Data()) :
2387  Form("fhJetPtGen%d%d%s",TMath::Nint(fCentralityBins[ic-1]),
2388  TMath::Nint(fCentralityBins[ic]),bgtype[ir].Data());
2389 
2390  fhJetPtGen[ic][ir] = new TH1D(name.Data(),
2391  Form("Jet pT Gen %s",bgtype[ir].Data()),bw*160,-20,300);
2392  fOutput->Add((TH1D*) fhJetPtGen[ic][ir]);
2393 
2394 
2395  name = (ic==0) ? Form("fhJetPtGenVsJetPtRecMB%s",bgtype[ir].Data()) :
2396  Form("fhJetPtGenVsJetPtRec%d%d%s",TMath::Nint(fCentralityBins[ic-1]),
2397  TMath::Nint(fCentralityBins[ic]),bgtype[ir].Data());
2398 
2399  fhJetPtGenVsJetPtRec[ic][ir] = new TH2D(name.Data(),
2400  "", bw*160,-20,300, bw*160,-20,300);
2401  fOutput->Add((TH2D*) fhJetPtGenVsJetPtRec[ic][ir]);
2402 
2403 
2404  name = (ic==0) ? Form("fhJetPtResolutionVsPtGenMB%s",bgtype[ir].Data()) :
2405  Form("fhJetPtResolutionVsPtGen%d%d%s",TMath::Nint(fCentralityBins[ic-1]),
2406  TMath::Nint(fCentralityBins[ic]),bgtype[ir].Data());
2407 
2408  fhJetPtResolutionVsPtGen[ic][ir] = new TH2D(name.Data(),
2409  "Resolution", 20,0,100, 200,-1.,1.);
2410  fOutput->Add((TH2D*) fhJetPtResolutionVsPtGen[ic][ir]);
2411  }
2412  }
2413 
2414 
2415  Double_t binsDiff [] = {
2416  -50, -48., -46., -44, -42.,
2417  -40., -38., -36., -34., -32., -30., -28., -26., -25., -24., -23., -22., -21.,
2418  -20., -19.5, -19., -18.5, -18., -17.5, -17., -16.5, -16., -15.5,
2419  -15., -14.5, -14., -13.5, -13., -12.5, -12., -11.5, -11., -10.5,
2420  -10.0, -9.8, -9.6, -9.4, -9.2,
2421  -9.0, -8.8, -8.6, -8.4, -8.2,
2422  -8.0, -7.8, -7.6, -7.4, -7.2,
2423  -7.0, -6.8, -6.6, -6.4, -6.2,
2424  -6.0, -5.8, -5.6, -5.4, -5.2,
2425  -5.0, -4.8, -4.6, -4.4, -4.2,
2426  -4.0, -3.8, -3.6, -3.4, -3.2,
2427  -3.0, -2.8, -2.6, -2.4, -2.2,
2428  -2.0, -1.9, -1.8, -1.7, -1.6, -1.5, -1.4, -1.3, -1.2, -1.1,
2429  -1.0,-0.95, -0.9, -0.85, -0.8, -0.75, -0.7, -0.65, -0.6, -0.55,
2430  -0.5,-0.48, -0.46, -0.44, -0.42,
2431  -0.4, -0.38, -0.36, -0.34, -0.32,
2432  -0.3, -0.28, -0.26, -0.24, -0.22,
2433  -0.2, -0.18, -0.16, -0.14, -0.12,
2434  -0.1, -0.09, -0.08, -0.07, -0.06,
2435  -0.05, -0.045, -0.04,-0.035 ,-0.03, -0.025, -0.02, -0.015, -0.01, -0.005,
2436  0, 0.005, 0.01, 0.015, 0.02, 0.025, 0.03, 0.035, 0.04, 0.045,
2437  0.05, 0.06, 0.07, 0.08, 0.09,
2438  0.1, 0.12, 0.14, 0.16, 0.18,
2439  0.2, 0.22, 0.24, 0.26, 0.28,
2440  0.3, 0.32, 0.34, 0.36, 0.38,
2441  0.4, 0.42, 0.44, 0.46, 0.48,
2442  0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8,0.85, 0.9, 0.95,
2443  1.0,1.1, 1.2,1.3, 1.4,1.5, 1.6, 1.7, 1.8, 1.9,
2444  2.0, 2.2, 2.4, 2.6, 2.8,
2445  3.0, 3.2, 3.4, 3.6, 3.8,
2446  4.0, 4.2, 4.4, 4.6, 4.8,
2447  5.0, 5.2, 5.4, 5.6, 5.8,
2448  6.0, 6.2, 6.4, 6.6, 6.8,
2449  7.0, 7.2, 7.4, 7.6, 7.8,
2450  8.0, 8.2, 8.4, 8.6, 8.8,
2451  9.0, 9.2, 9.4, 9.6, 9.8,
2452  10., 10.5, 11., 11.5, 12., 12.5, 13., 13.5, 14., 14.5,
2453  15., 15.5, 16., 16.5, 17., 17.5, 18., 18.5, 19., 19.5,
2454  20., 21., 22., 23., 24, 25, 26., 28., 30., 32.,34., 36., 38.,
2455  40., 42., 44., 46., 48., 50.};
2456  Int_t nbinsDiff = sizeof(binsDiff)/sizeof(Double_t)-1;
2457 
2458  for(Int_t ic =0; ic<icmax; ic++){
2459  name = (ic==0) ? Form("fhPtTrkTruePrimRecMB") :
2460  Form("fhPtTrkTruePrimRec%d%d",TMath::Nint(fCentralityBins[ic-1]),
2461  TMath::Nint(fCentralityBins[ic]));
2462 
2463  fhPtTrkTruePrimRec[ic] = new TH2D(name.Data(),"",nbins, bins, 18,-0.9,0.9);
2464  fOutput->Add((TH2D*) fhPtTrkTruePrimRec[ic]);
2465 
2466 
2467  name = (ic==0) ? Form("fhPtTrkTruePrimGenMB") :
2468  Form("fhPtTrkTruePrimGen%d%d",TMath::Nint(fCentralityBins[ic-1]),
2469  TMath::Nint(fCentralityBins[ic]));
2470 
2471  fhPtTrkTruePrimGen[ic] = (TH2D*) fhPtTrkTruePrimRec[ic]->Clone(name.Data());
2472  fOutput->Add((TH2D*) fhPtTrkTruePrimGen[ic]);
2473 
2474  name = (ic==0) ? Form("fhPtTrkSecOrFakeRecMB") :
2475  Form("fhPtTrkSecOrFakeRec%d%d",TMath::Nint(fCentralityBins[ic-1]),
2476  TMath::Nint(fCentralityBins[ic]));
2477 
2478  fhPtTrkSecOrFakeRec[ic] = (TH2D*) fhPtTrkTruePrimRec[ic]->Clone(name.Data());
2479  fOutput->Add((TH2D*) fhPtTrkSecOrFakeRec[ic]);
2480  }
2481 
2482  for(Int_t iw=0; iw<21; iw++){
2483  name = Form("fhPtJetPrimVsPtJetRecMB%03d",TMath::Nint(100*iw*0.05));
2484  fhPtJetPrimVsPtJetRec[iw] = new TH2D(name.Data(), name.Data(),50,0,50,210,0,1.05);
2485  fOutput->Add((TH2D*) fhPtJetPrimVsPtJetRec[iw]);
2486  }
2487  name = Form("fhDiffPtVsPtTrackTrueMB");
2488  fhDiffPtVsPtTrackTrue = new TH2D(name.Data(), name.Data(),100,0,100, nbinsDiff,binsDiff);
2489  fOutput->Add((TH2D*) fhDiffPtVsPtTrackTrue);
2490  }
2491 
2492  //DCA DISTRIBUTIONS
2493  fhDCAinXVsPt = new TH2D("fhDCAinXVsPt","fhDCAinXVsPt",nbins, bins, 200, -10.,10);
2494  if(bNotKine) fOutput->Add((TH2D*) fhDCAinXVsPt);
2495  fhDCAinYVsPt = (TH2D*) fhDCAinXVsPt->Clone("fhDCAinYVsPt");
2496  if(bNotKine) fOutput->Add((TH2D*) fhDCAinYVsPt);
2497 
2498  fhDCAinXVsPtStrange = (TH2D*) fhDCAinXVsPt->Clone("fhDCAinXVsPtStrange");
2499  if(fTypeOfAnal == kEff) fOutput->Add((TH2D*) fhDCAinXVsPtStrange);
2500  fhDCAinYVsPtStrange = (TH2D*) fhDCAinXVsPt->Clone("fhDCAinYVsPtStrange");
2501  if(fTypeOfAnal == kEff) fOutput->Add((TH2D*) fhDCAinYVsPtStrange);
2502  fhDCAinXVsPtNonStrange = (TH2D*) fhDCAinXVsPt->Clone("fhDCAinXVsPtNonStrange");
2503  if(fTypeOfAnal == kEff) fOutput->Add((TH2D*) fhDCAinXVsPtNonStrange);
2504  fhDCAinYVsPtNonStrange = (TH2D*) fhDCAinXVsPt->Clone("fhDCAinYVsPtNonStrange");
2505  if(fTypeOfAnal == kEff) fOutput->Add((TH2D*) fhDCAinYVsPtNonStrange);
2506 
2507 
2508 
2509  if(fTypeOfAnal == kEff || fTypeOfAnal == kKine){
2510  for(Int_t ic =0; ic<icmax; ic++){
2511  for(Int_t it=0; it<kTT;it++){
2512  name = (ic==0) ? Form("fh1NtriggersGenMB") :
2513  Form("fh1NtriggersGen%d%d",TMath::Nint(fCentralityBins[ic-1]),
2514  TMath::Nint(fCentralityBins[ic]));
2515  name.Append(Form("TT%d%d",TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
2516 
2517  fh1NtriggersGen[ic][it] = (TH1D*) fh1Ntriggers[ic][it]->Clone(name.Data());
2518  fOutput->Add((TH1D*) fh1NtriggersGen[ic][it]);
2519 
2520 
2521  name = (ic==0) ? Form("fh1TriggerMultGenMB") :
2522  Form("fh1TriggerMultGen%d%d",TMath::Nint(fCentralityBins[ic-1]),
2523  TMath::Nint(fCentralityBins[ic]));
2524 
2525  name.Append(Form("TT%d%d",TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it])));
2526  fh1TriggerMultGen[ic][it] = (TH1D*) fh1TriggerMult[ic][it]->Clone();
2527  fOutput->Add((TH1D*) fh1TriggerMultGen[ic][it]);
2528 
2529  //-------------------------
2530  for(Int_t ir=0; ir< kRho; ir++){
2531  name = (ic==0) ? "fHJetSpecGenMB" :
2532  Form("fHJetSpecGen%d%d",TMath::Nint(fCentralityBins[ic-1]),
2533  TMath::Nint(fCentralityBins[ic]));
2534 
2535  name.Append(Form("TT%d%d%s", TMath::Nint(fTTlow[it]), TMath::Nint(fTThigh[it]), bgtype[ir].Data()));
2536  fHJetSpecGen[ic][it][ir] = (TH2D*)fHJetSpec[ic][it][ir]->Clone(name.Data());
2537  fOutput->Add((TH2D*) fHJetSpecGen[ic][it][ir]);
2538  }
2539  }
2540  //-------------------------
2541  name = (ic==0) ? Form("fhJetPhiGenMB") :
2542  Form("fhJetPhiGen%d%d",TMath::Nint(fCentralityBins[ic-1]),
2543  TMath::Nint(fCentralityBins[ic]));
2544 
2545  fhJetPhiGen[ic] = (TH2F*) fhJetPhi[ic]->Clone(name.Data());
2546  fOutput->Add((TH2F*) fhJetPhiGen[ic]);
2547  //-------------------------
2548  name = (ic==0) ? Form("fhJetEtaGenMB") :
2549  Form("fhJetEtaGen%d%d",TMath::Nint(fCentralityBins[ic-1]),
2550  TMath::Nint(fCentralityBins[ic]));
2551 
2552 
2553  fhJetEtaGen[ic] = (TH2F*) fhJetEta[ic]->Clone(name.Data());
2554  fOutput->Add((TH2F*) fhJetEtaGen[ic]);
2555  //-------------------------
2556  name = (ic==0) ? Form("fhTrackPtGenMB") :
2557  Form("fhTrackPtGen%d%d",TMath::Nint(fCentralityBins[ic-1]),
2558  TMath::Nint(fCentralityBins[ic]));
2559 
2560  fhTrackPtGen[ic] = (TH1F*) fhTrackPt[ic]->Clone(name.Data());
2561  fOutput->Add((TH1F*) fhTrackPtGen[ic]);
2562  //------------------------- Rongrong's analysis
2563  for(Int_t it=0; it< kTT; it++){
2564  for(Int_t ir=0; ir< kRho; ir++){
2565 
2566  name = Form("%sGen",fhDphiTriggerJet[ic][it][ir]->GetName());
2567  fhDphiTriggerJetGen[ic][it][ir] = (TH2F*) fhDphiTriggerJet[ic][it][ir]->Clone(name.Data());
2568  fOutput->Add((TH2F*) fhDphiTriggerJetGen[ic][it][ir]);
2569  }
2570  }
2571  }
2572  }
2573 
2574  // =========== Switch on Sumw2 for all histos ===========
2575  for(Int_t i=0; i<fOutput->GetEntries(); i++){
2576  TH1 *h1 = dynamic_cast<TH1*>(fOutput->At(i));
2577  if(h1){
2578  h1->Sumw2();
2579  continue;
2580  }
2581  THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i));
2582  if(hn){
2583  hn->Sumw2();
2584  }
2585  }
2586  TH1::AddDirectory(oldStatus);
2587 
2588 
2589  PostData(1, fOutput);
2590 }
2591 //________________________________________________________________________
2592 Bool_t AliAnalysisTaskHJetSpectra::RetrieveEventObjects() {
2593  //
2594  // retrieve event objects
2595  //
2596  if(!AliAnalysisTaskEmcalJet::RetrieveEventObjects()) return kFALSE;
2597 
2598  return kTRUE;
2599 }
2600 //________________________________________________________________________
2601 Bool_t AliAnalysisTaskHJetSpectra::Run()
2602 {
2603  // Run analysis code here, if needed. It will be executed before FillHistograms().
2604 
2605  return kTRUE;
2606 }
2607 
2608 //________________________________________________________________________
2609 
2610 Double_t AliAnalysisTaskHJetSpectra::RelativePhi(Double_t mphi,Double_t vphi){
2611  //Get relative azimuthal angle of two particles -pi to pi
2612  if (vphi < -TMath::Pi()) vphi += TMath::TwoPi();
2613  else if (vphi > TMath::Pi()) vphi -= TMath::TwoPi();
2614 
2615  if (mphi < -TMath::Pi()) mphi += TMath::TwoPi();
2616  else if (mphi > TMath::Pi()) mphi -= TMath::TwoPi();
2617 
2618  Double_t dphi = mphi - vphi;
2619  if (dphi < -TMath::Pi()) dphi += TMath::TwoPi();
2620  else if (dphi > TMath::Pi()) dphi -= TMath::TwoPi();
2621 
2622  return dphi;//dphi in [-Pi, Pi]
2623 }
2624 
2625 
2626 //________________________________________________________________________
2627 Double_t AliAnalysisTaskHJetSpectra::EstimateBgCone(AliJetContainer *jetCont, AliParticleContainer *trkCont, AliVParticle* triggerHadron, Bool_t isGen){
2628  //Estimate background rho by means of integrating track pT outside TT jet + recoil jet region
2629  //if TT exists find jet that containts TT and exclude range +- phiCut around the TT/TTjet in azimuth
2630 
2631  if(!trkCont) return 0.0;
2632 
2633  AliEmcalJet* jet = NULL;
2634  AliVParticle* track = NULL;
2635  Double_t phiTT = fRandom->Rndm()*TMath::TwoPi(); //in case of no TT make random dice
2636  Double_t etaTT = -fTrackEtaWindow + fRandom->Rndm()*2*fTrackEtaWindow;
2637  Bool_t bTTJetFound = kFALSE;
2638 
2639  if(triggerHadron){
2640 
2641  phiTT = triggerHadron->Phi();
2642  etaTT = triggerHadron->Eta();
2643 
2644  if(jetCont){
2645  //find ANY jet that contains TT if it exists
2646  jetCont->ResetCurrentID();
2647  while((jet = jetCont->GetNextAcceptJet())){ //loop over reconstructed jets
2648  if(!jet) continue;
2649  if(jet->Pt() < triggerHadron->Pt()*0.5) continue;
2650  //CUT ON JET ACCEPTANCE IS NOT NEEDED, WE SEARCH FOR ANY JET THAT CONTAINS TT
2651 
2652  for(Int_t iq=0; iq < jet->GetNumberOfTracks(); iq++){
2653  track = (AliVParticle*) (jet->TrackAt(iq,trkCont->GetArray())); //matched rec and emb tracks
2654  if(!track) continue;
2655  if(track != triggerHadron) continue;
2656 
2657  phiTT = jet->Phi(); // used phi,eta coordinates of the jet to exclude the TT jet
2658  etaTT = jet->Eta();
2659  bTTJetFound = kTRUE;
2660  break;
2661  }
2662  if(bTTJetFound) break; //skip the rest of jets when the jet with TT is found
2663  }
2664  }
2665  }
2666 
2667  phiTT = RelativePhi(phiTT,0.); //convert phi TT to (-pi,pi)
2668 
2669  if(TMath::Abs(etaTT) > fTrackEtaWindow){
2670  etaTT = (etaTT<0) ? -fTrackEtaWindow : fTrackEtaWindow;
2671  }
2672  //Sum pT outside TT+recoil jet region
2673  Double_t sumPt = 0.;
2674 
2675  trkCont->ResetCurrentID();
2676  while((track = trkCont->GetNextAcceptParticle())){
2677  if(!track) continue;
2678  if(!IsTrackInAcceptance(track, isGen)) continue;
2679 
2680  if(TMath::Abs(RelativePhi(phiTT+TMath::Pi(),track->Phi())) < fCutPhi) continue; //exclude recoil region of TT
2681  if(GetDeltaR(phiTT, track->Phi(), etaTT, track->Eta()) < fCutPhi) continue; //exclude region around TT
2682 
2683  if(fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
2684  if(TMath::Abs(track->GetLabel()) == 99999) continue;//reject embedded stuff 9999 set in AddTaskHJetSpectra.C
2685  }
2686 
2687  sumPt += track->Pt();
2688  }
2689  //Calculate area
2690  Double_t area = 2*fTrackEtaWindow*2*TMath::Pi();
2691  Double_t alpha;
2692 
2693  area -= 2*fTrackEtaWindow*2*fCutPhi; // subtract area of the recoil region
2694 
2695  if(TMath::Abs(etaTT) < fTrackEtaWindow - fCutPhi){ //TT cicle fully inside acceptance
2696  area -= fCutPhi*fCutPhi*TMath::Pi(); // subtract area of the trigger region
2697  }else{ //TT circle partly around acceptance
2698  alpha = TMath::ACos((fTrackEtaWindow - TMath::Abs(etaTT))/fCutPhi);
2699  area -= (fCutPhi*fCutPhi*(TMath::Pi()-alpha) + fCutPhi*TMath::Sin(alpha)*(fTrackEtaWindow - TMath::Abs(etaTT)));
2700  }
2701 
2702  return sumPt/area;
2703 
2704 }
2705 //________________________________________________________________________
2706 Double_t AliAnalysisTaskHJetSpectra::EstimateBgKT(AliJetContainer *jetCont, AliParticleContainer *trkCont, AliVParticle* triggerHadron){
2707  //Estimate rho from KT jet median. Ignore jet that contains TT
2708  Double_t rhoKT = 0.0;
2709 
2710  if(!jetCont) return rhoKT;
2711 
2712  AliEmcalJet* jet = NULL;
2713  AliVParticle* constTrack = NULL;
2714  Bool_t bKTJetCloseToTT = kFALSE;
2715  Int_t nJetAcc = 0;
2716  Double_t jetpt;
2717  Double_t sumEmbPt;
2718 
2719  jetCont->ResetCurrentID();
2720  while((jet = jetCont->GetNextAcceptJet())){ //loop over KT jets
2721  if(!jet) continue;
2722  if(!IsSignalJetInAcceptance(jet,kFALSE)) continue;
2723 
2724  bKTJetCloseToTT = kFALSE;
2725 
2726  if(triggerHadron){ //identify the KT jet which contains TT
2727  if(jet->Pt() > triggerHadron->Pt()*0.5){ //jet containing TT has pT larger than pT of TT
2728  for(Int_t iq=0; iq < jet->GetNumberOfTracks(); iq++) {
2729  constTrack = (AliVParticle*) (jet->TrackAt(iq,trkCont->GetArray())); //matched rec and emb tracks
2730  if(!constTrack) continue;
2731  if(constTrack == triggerHadron){
2732  bKTJetCloseToTT = kTRUE;
2733  break;
2734  }
2735  }
2736  }
2737  }
2738  if(bKTJetCloseToTT) continue; //skip the jet that contains TT
2739 
2740  sumEmbPt = 0.;//sum pt of embedded tracks in jet which is to be subtracted
2741  if(fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
2742  for(Int_t iq=0; iq < jet->GetNumberOfTracks(); iq++) {
2743  constTrack = (AliVParticle*) (jet->TrackAt(iq,trkCont->GetArray()));
2744  if(!constTrack) continue;
2745  if(TMath::Abs(constTrack->GetLabel()) == 99999){
2746  sumEmbPt += constTrack->Pt();
2747  }
2748  }
2749  }
2750 
2751  jetpt = jet->Pt()- sumEmbPt; //subtract embedded pt
2752  if(triggerHadron) fhKTAreaPt->Fill(jetpt,jet->Area());
2753 
2754  if(jetpt <0.005) jetpt = 0.; //set pt of ghost jets identical to zero
2755  frhovec[nJetAcc] = jetpt/jet->Area();
2756  nJetAcc++;
2757  }
2758 
2759  if(nJetAcc>0){
2760  rhoKT = TMath::Median(nJetAcc, frhovec);
2761  }
2762 
2763  return rhoKT;
2764 }
2765 //________________________________________________________________________
2766 Double_t AliAnalysisTaskHJetSpectra::EstimateBgKTcms(AliJetContainer *jetCont, AliParticleContainer *trkCont, AliVParticle* triggerHadron){
2767  //Estimate rho from KT jet median ala CMS. Ignore jet that contains TT
2768  Double_t rhoKTcms = 0.0;
2769 
2770  if(!jetCont) return rhoKTcms;
2771 
2772  AliEmcalJet* jet = NULL;
2773  AliVParticle* constTrack = NULL;
2774  Bool_t bKTJetCloseToTT = kFALSE;
2775  Int_t nJetAcc = 0;
2776  Double_t areaPhysJets = 0.0;
2777  Double_t areaAllJets = 0.0;
2778  Double_t jetpt;
2779  Double_t sumEmbPt = 0.;
2780 
2781  jetCont->ResetCurrentID();
2782  while((jet = jetCont->GetNextAcceptJet())){ //loop over KT jets
2783  if(!jet) continue;
2784  if(!IsSignalJetInAcceptance(jet,kFALSE)) continue;
2785 
2786  bKTJetCloseToTT = kFALSE;
2787 
2788  if(triggerHadron){ //identify the KT jet which contains TT
2789  if(jet->Pt() > triggerHadron->Pt()*0.5){ //jet containing TT has pT larger than pT of TT
2790  for(Int_t iq=0; iq < jet->GetNumberOfTracks(); iq++) {
2791  constTrack = (AliVParticle*) (jet->TrackAt(iq,trkCont->GetArray())); //matched rec and emb tracks
2792  if(!constTrack) continue;
2793  if(constTrack != triggerHadron) continue;
2794  bKTJetCloseToTT = kTRUE;
2795  break;
2796  }
2797  }
2798  }
2799  if(bKTJetCloseToTT) continue; //skip the jet that contains TT
2800 
2801  sumEmbPt = 0.;
2802  if(fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
2803  for(Int_t iq=0; iq < jet->GetNumberOfTracks(); iq++) {
2804  constTrack = (AliVParticle*) (jet->TrackAt(iq,trkCont->GetArray())); //matched rec and emb tracks
2805  if(!constTrack) continue;
2806  if(TMath::Abs(constTrack->GetLabel()) == 99999){
2807  sumEmbPt += constTrack->Pt();
2808  }
2809  }
2810  }
2811 
2812 
2813  areaAllJets += jet->Area();
2814 
2815  jetpt = jet->Pt()- sumEmbPt; //subtract pt of embedded tracks
2816 
2817  if(jetpt > 0.1){
2818  areaPhysJets += jet->Area();
2819  frhovec[nJetAcc] = jetpt/jet->Area();
2820  nJetAcc++;
2821  }
2822  }
2823 
2824  if(nJetAcc>0){
2825  rhoKTcms = TMath::Median(nJetAcc, frhovec)*(areaPhysJets/areaAllJets);
2826  }
2827 
2828  return rhoKTcms;
2829 }
2830 
2831 //________________________________________________________________________
2832 Double_t AliAnalysisTaskHJetSpectra::GetDeltaR(Double_t phi1, Double_t phi2, Double_t eta1, Double_t eta2){
2833  //angular distance between two jets
2834  Double_t dphi = RelativePhi(phi1,phi2);
2835  Double_t deta = eta1 - eta2;
2836  return sqrt(dphi*dphi + deta*deta);
2837 
2838 }
2839 
2840 //________________________________________________________________________
2841 Double_t AliAnalysisTaskHJetSpectra::GetFractionSharedPt(AliEmcalJet *jRec, AliJetContainer *jconRec, AliEmcalJet *jGen, AliJetContainer *jconGen){
2842 
2843  //get fraction of pT shared by reconstructed and generated level jet
2844  if(!jRec) return -1.0;
2845  if(!jconRec) return -1.0;
2846  if(!jGen) return -1.0;
2847  if(!jconGen) return -1.0;
2848 
2849  Double_t fraction = 0., sumPt = 0.;
2850  Double_t jetPt2 = jGen->Pt();
2851  //Int_t idxGen, idxRec;
2852  AliVParticle *pgen, *prec;
2853  if(jetPt2>0){
2854 
2855  for(Int_t ig=0; ig< jGen->GetNumberOfTracks(); ig++) {
2856  pgen = (AliVParticle*) (jGen->TrackAt(ig, jconGen->GetParticleContainer()->GetArray()));
2857  if(!pgen) continue;
2858 
2859  for(Int_t ir=0; ir< jRec->GetNumberOfTracks(); ir++){
2860  prec = (AliVParticle*) (jRec->TrackAt(ir, jconRec->GetParticleContainer()->GetArray()));
2861  if(!prec) continue;
2862 
2863  if(TMath::Abs(prec->GetLabel()) == TMath::Abs(pgen->GetLabel())){
2864 
2865  if(fTypeOfAnal == kEmb || fTypeOfAnal == kEmbSingl){
2866  //All embedded tracks have the same label check also spatial coordinates
2867  if(TMath::Abs(prec->Eta() - pgen->Eta()) > 1e-4) continue;
2868  if(TMath::Abs(RelativePhi(prec->Phi(), pgen->Phi())) > 1e-4) continue;
2869  if(TMath::Abs(prec->Pt() - pgen->Pt()) > 1e-4) continue;
2870  //if(fDebug>20){
2871  //Printf("fraction TRACK REC eta = %f; phi = %f; pt = %f", prec->Eta(), prec->Phi(), prec->Pt());
2872  //Printf("fraction TRACK GEN eta = %f; phi = %f; pt = %f", pgen->Eta(), pgen->Phi(), pgen->Pt());
2873  //}
2874 
2875  }
2876 
2877  sumPt += pgen->Pt();
2878  break;
2879  }
2880  }
2881  }
2882 
2883 
2884  fraction = sumPt/jetPt2;
2885  } else{
2886  fraction = -1;
2887  }
2888 
2889  //if(fDebug>20) Printf("fraction return = %f ",fraction);
2890 
2891  return fraction;
2892 
2893 }
2894 
2895 //________________________________________________________________________
2896 Bool_t AliAnalysisTaskHJetSpectra::IsStrange(Int_t ip){
2897  //Check whtether the particle is strange
2898  if(ip==321) return kTRUE; //K+
2899  if(ip==130) return kTRUE; //KOL
2900  if(ip==310) return kTRUE; //K0s
2901  if(ip==3112) return kTRUE; //Sigma -
2902  if(ip==3122) return kTRUE; //Lambda 0
2903  if(ip==3222) return kTRUE; //Sigma+
2904  if(ip==3312) return kTRUE; //Xi-
2905  if(ip==3322) return kTRUE; //Xi0
2906  if(ip==3334) return kTRUE; //Omega-
2907 
2908  return kFALSE;
2909 }
2910 
AliAnalysisTaskHJetSpectra::GetDeltaR
Double_t GetDeltaR(Double_t phi1, Double_t phi2, Double_t eta1, Double_t eta2)
Definition: AliAnalysisTaskHJetSpectra.cxx:2832
AliAnalysisTaskHJetSpectra::fhPtTrkSecOrFakeRec
TH2D * fhPtTrkSecOrFakeRec[kCAll]
pt spectrum of true generated primary track
Definition: AliAnalysisTaskHJetSpectra.h:298
AliAnalysisTaskHJetSpectra::fhVertexZAcceptTT
TH1F * fhVertexZAcceptTT
gc vertexZ accepted after vtx cut
Definition: AliAnalysisTaskHJetSpectra.h:256
AliAnalysisTaskHJetSpectra::SetTT
void SetTT(Double_t tlr, Double_t thr, Double_t tls, Double_t ths)
Definition: AliAnalysisTaskHJetSpectra.cxx:393
AliAnalysisTaskHJetSpectra::fhJetPhiGen
TH2F * fhJetPhiGen[kCAll]
gc jet phi vs jet pT
Definition: AliAnalysisTaskHJetSpectra.h:243
AliAnalysisTaskHJetSpectra::fhVertexZ
TH1F * fhVertexZ
gc X=centrality; Y= track pT
Definition: AliAnalysisTaskHJetSpectra.h:250
AliAnalysisTaskHJetSpectra::fTTlow
Double_t fTTlow[kTT]
gc trigger if tracks/jets are loaded initiates calling ExecOnce
Definition: AliAnalysisTaskHJetSpectra.h:209
AliEmcalJet::Area
Double_t Area() const
Definition: AliEmcalJet.h:123
AliParticleContainer::GetNextAcceptParticle
virtual AliVParticle * GetNextAcceptParticle()
Definition: AliParticleContainer.cxx:129
AliAnalysisTaskHJetSpectra::fhVzeroATotMultTT
TH1F * fhVzeroATotMultTT[kCAll][kTT]
V0A multiplicity for given V0A centrality selection.
Definition: AliAnalysisTaskHJetSpectra.h:275
AliPicoTrack.h
Double_t
double Double_t
Definition: External.C:58
AliAnalysisTaskHJetSpectra.h
AliAnalysisTaskHJetSpectra::fhDeltaPtEmbBc2Bc
TH1D * fhDeltaPtEmbBc2Bc[kCAll][kTT][kRho]
embedded delta pT versus pT of the embedded jet (emb track is perp to TT)
Definition: AliAnalysisTaskHJetSpectra.h:235
AliAnalysisTaskHJetSpectra::fTrigTracksGen
AliVParticle * fTrigTracksGen[kTT][999]
Definition: AliAnalysisTaskHJetSpectra.h:330
AliAnalysisTaskHJetSpectra::fhDCAinXVsPt
TH2D * fhDCAinXVsPt
hybrid TPC constrained track phi vs track pT
Definition: AliAnalysisTaskHJetSpectra.h:310
AliAnalysisTaskEmcalJet.h
TH2F
Definition: External.C:236
AliAnalysisTaskHJetSpectra::fhCentralityV0C
TH1F * fhCentralityV0C
centrality from V0A
Definition: AliAnalysisTaskHJetSpectra.h:271
AliAnalysisTaskHJetSpectra::fhVertexZAccept
TH1F * fhVertexZAccept
gc vertexZ accepted after vtx cut
Definition: AliAnalysisTaskHJetSpectra.h:253
AliEmcalJet::ClosestJet
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:219
AliAnalysisTaskHJetSpectra::fhTrackPhi
TH2F * fhTrackPhi[kCAll]
gc jet phi vs jet pT
Definition: AliAnalysisTaskHJetSpectra.h:244
AliAnalysisTaskEmcalJet::GetJetContainer
AliJetContainer * GetJetContainer(Int_t i=0) const
Definition: AliAnalysisTaskEmcalJet.cxx:361
AliAnalysisTaskHJetSpectra::EstimateBgKTcms
Double_t EstimateBgKTcms(AliJetContainer *jetCont, AliParticleContainer *trkArray, AliVParticle *triggerHadron)
Definition: AliAnalysisTaskHJetSpectra.cxx:2766
AliAnalysisTaskHJetSpectra::fhJetPtResolutionVsPtGen
TH2D * fhJetPtResolutionVsPtGen[kCAll][kRho]
pt jet gen level vs pT jet rec level
Definition: AliAnalysisTaskHJetSpectra.h:295
AliEmcalJet::Eta
Double_t Eta() const
Definition: AliEmcalJet.h:114
AliAnalysisTaskHJetSpectra::EstimateBgCone
Double_t EstimateBgCone(AliJetContainer *jetCont, AliParticleContainer *trkArray, AliVParticle *triggerHadron, Bool_t isGen=kFALSE)
Definition: AliAnalysisTaskHJetSpectra.cxx:2627
AliAnalysisTaskHJetSpectra::fhPtTrkTruePrimGen
TH2D * fhPtTrkTruePrimGen[kCAll]
pt spectrum of true reconstructed primary tracks
Definition: AliAnalysisTaskHJetSpectra.h:297
AliEmcalJet::Phi
Double_t Phi() const
Definition: AliEmcalJet.h:110
AliAnalysisTaskHJetSpectra::fhKTAreaPt
TH2F * fhKTAreaPt
delta pT from RndCone using rho from perp cone inclusive event
Definition: AliAnalysisTaskHJetSpectra.h:240
AliAnalysisTaskHJetSpectra::IsEventInAcceptance
Bool_t IsEventInAcceptance(AliVEvent *event)
Definition: AliAnalysisTaskHJetSpectra.cxx:537
AliAnalysisTaskHJetSpectra::fhDeltaPt
TH1D * fhDeltaPt[kCAll][kTT][kRho]
jet area times rho from perp cone
Definition: AliAnalysisTaskHJetSpectra.h:230
ClassImp
ClassImp(AliAnalysisTaskHJetSpectra) AliAnalysisTaskHJetSpectra
Definition: AliAnalysisTaskHJetSpectra.cxx:79
AliAnalysisTaskHJetSpectra::fhDeltaPtEmbPerp
TH1D * fhDeltaPtEmbPerp[kCAll][kTT][kRho]
embedded delta pT versus pT of the embedded jet
Definition: AliAnalysisTaskHJetSpectra.h:233
AliAnalysisTaskHJetSpectra::IsSignalJetInAcceptance
Bool_t IsSignalJetInAcceptance(AliEmcalJet *jet, Bool_t suppressGhost=1)
Definition: AliAnalysisTaskHJetSpectra.cxx:606
AliAnalysisTaskHJetSpectra::fhDeltaPtEmb2D
TH2D * fhDeltaPtEmb2D[kCAll][kTT][kRho]
embedded delta pT
Definition: AliAnalysisTaskHJetSpectra.h:232
AliAnalysisTaskHJetSpectra::fh1NtriggersGen
TH1D * fh1NtriggersGen[kCAll][kTT]
tirgger multiplicity in event
Definition: AliAnalysisTaskHJetSpectra.h:219
AliAnalysisTaskHJetSpectra::fhRhoIncl
TH1F * fhRhoIncl[kCAll][kRho]
gc X=rho from perp cone, Y=centrality
Definition: AliAnalysisTaskHJetSpectra.h:225
fCollisionSystem
Int_t fCollisionSystem
Definition: makeTFile4CutsLctoV0bachelor.C:26
AliAnalysisTaskHJetSpectra::fRhoRec
TArrayD fRhoRec[kTT]
Y DCA versus pT of non strange tracks.
Definition: AliAnalysisTaskHJetSpectra.h:318
AliAnalysisTaskHJetSpectra::fhVertexXAccept
TH1F * fhVertexXAccept
gc vertexZ inclusive
Definition: AliAnalysisTaskHJetSpectra.h:251
AliAnalysisTaskHJetSpectra::fhDeltaPtEmb
TH1D * fhDeltaPtEmb[kCAll][kTT][kRho]
delta pT
Definition: AliAnalysisTaskHJetSpectra.h:231
AliAnalysisTaskHJetSpectra::fh1TriggerMultGen
TH1D * fh1TriggerMultGen[kCAll][kTT]
trigger counter
Definition: AliAnalysisTaskHJetSpectra.h:220
AliAnalysisHelperJetTasks.h
AliAnalysisTaskHJetSpectra::fhJetEta
TH2F * fhJetEta[kCAll]
gc track phi vs track pT
Definition: AliAnalysisTaskHJetSpectra.h:245
AliAnalysisTaskHJetSpectra::fhJetPhi
TH2F * fhJetPhi[kCAll]
KT jets area versus PT.
Definition: AliAnalysisTaskHJetSpectra.h:242
AliAnalysisTaskHJetSpectra::fhTrackPt
TH1F * fhTrackPt[kCAll]
track eta vs track pT
Definition: AliAnalysisTaskHJetSpectra.h:248
AliAnalysisTaskHJetSpectra::fhDeltaPtIncl
TH1D * fhDeltaPtIncl[kCAll][kRho]
embedded delta pT versus pT of the embedded jet (emb track is backtoback in azimtuh w...
Definition: AliAnalysisTaskHJetSpectra.h:238
AliAnalysisTaskHJetSpectra::fhInvPtQVsPhi
TH2D * fhInvPtQVsPhi[2]
track Y= rec pt - true pt X= true track pT
Definition: AliAnalysisTaskHJetSpectra.h:302
AliAnalysisTaskHJetSpectra::fRandom
TRandom3 * fRandom
impact parameter from hijing
Definition: AliAnalysisTaskHJetSpectra.h:204
AliAnalysisTaskHJetSpectra::fh1TriggerMult
TH1D * fh1TriggerMult[kCAll][kTT]
trigger counter
Definition: AliAnalysisTaskHJetSpectra.h:218
AliAnalysisTaskHJetSpectra::fhSigmaPtOverPtVsPt
TH2D * fhSigmaPtOverPtVsPt[2]
q*1/pT versus eta
Definition: AliAnalysisTaskHJetSpectra.h:306
AliAnalysisTaskHJetSpectra::fhVertexXAcceptTT
TH1F * fhVertexXAcceptTT
gc vertexZ accepted after vtx cut
Definition: AliAnalysisTaskHJetSpectra.h:254
AliParticleContainer
Container for particles within the EMCAL framework.
Definition: AliParticleContainer.h:28
AliEmcalJet::GetNumberOfTracks
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:132
AliAnalysisTaskHJetSpectra::GetFractionSharedPt
Double_t GetFractionSharedPt(AliEmcalJet *jRec, AliJetContainer *jconRec, AliEmcalJet *jGen, AliJetContainer *jconGen)
Definition: AliAnalysisTaskHJetSpectra.cxx:2841
AliAnalysisTaskHJetSpectra::fhZNAVzeroATrackTT
THnSparse * fhZNAVzeroATrackTT[kCAll][kTT]
ZNA energy versus Vzero A mult versus track mult.
Definition: AliAnalysisTaskHJetSpectra.h:284
AliAnalysisTaskHJetSpectra::fhTrackMultiplicity
TH1D * fhTrackMultiplicity[kCAll]
ZDC A neutral energy.
Definition: AliAnalysisTaskHJetSpectra.h:280
AliAnalysisTaskEmcal::GetParticleContainer
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Definition: AliAnalysisTaskEmcal.cxx:1627
AliAnalysisTaskHJetSpectra::fhDphiTriggerJet
TH2F * fhDphiTriggerJet[kCAll][kTT][kRho]
gc vertexZ accepted after vtx cut in MC
Definition: AliAnalysisTaskHJetSpectra.h:259
AliAnalysisTaskHJetSpectra::fhCentrality
TH1F * fhCentrality[kCAll]
minimum bias eta inclusive
Definition: AliAnalysisTaskHJetSpectra.h:267
AliAnalysisTaskHJetSpectra::IsMCEventInAcceptance
Bool_t IsMCEventInAcceptance(AliVEvent *event)
Definition: AliAnalysisTaskHJetSpectra.cxx:514
AliJetContainer::GetParticleContainer
AliParticleContainer * GetParticleContainer() const
Definition: AliJetContainer.h:180
AliAnalysisTaskHJetSpectra::fhJetEtaGen
TH2F * fhJetEtaGen[kCAll]
jet eta vs jet pT
Definition: AliAnalysisTaskHJetSpectra.h:246
AliAnalysisTaskHJetSpectra::fhPtJetPrimVsPtJetRec
TH2D * fhPtJetPrimVsPtJetRec[21]
pt spectrum of reconstructed fake or secondary tracks
Definition: AliAnalysisTaskHJetSpectra.h:299
AliAnalysisTaskHJetSpectra::fhTrackPtGen
TH1F * fhTrackPtGen[kCAll]
gc X=centrality; Y= track pT
Definition: AliAnalysisTaskHJetSpectra.h:249
AliAnalysisTaskHJetSpectra::fhZNAVzeroATrack
THnSparse * fhZNAVzeroATrack[kCAll]
multiplicity of tracks in event with TT track
Definition: AliAnalysisTaskHJetSpectra.h:283
AliAnalysisTaskHJetSpectra::fh1Ntriggers
TH1D * fh1Ntriggers[kCAll][kTT]
gc event statistics
Definition: AliAnalysisTaskHJetSpectra.h:217
AliAnalysisTaskHJetSpectra::fhJetPtGen
TH1D * fhJetPtGen[kCAll][kRho]
impact parameter distribution hijing versus TT
Definition: AliAnalysisTaskHJetSpectra.h:293
AliAnalysisTaskHJetSpectra::fhDCAinYVsPtNonStrange
TH2D * fhDCAinYVsPtNonStrange
X DCA versus pT of non strange tracks.
Definition: AliAnalysisTaskHJetSpectra.h:315
Int_t
int Int_t
Definition: External.C:63
AliAnalysisTaskHJetSpectra::fhCentralityV0A
TH1F * fhCentralityV0A
centrality V0 multiplicity A+C
Definition: AliAnalysisTaskHJetSpectra.h:270
TH1I
Definition: External.C:204
UInt_t
unsigned int UInt_t
Definition: External.C:33
Float_t
float Float_t
Definition: External.C:68
AliAnalysisTaskHJetSpectra::IsTrackInAcceptance
Bool_t IsTrackInAcceptance(AliVParticle *track, Bool_t isGen=0)
Definition: AliAnalysisTaskHJetSpectra.cxx:585
AliAnalysisTaskHJetSpectra::fhZNAEnergy
TH1F * fhZNAEnergy[kCAll]
V0A multiplicity.
Definition: AliAnalysisTaskHJetSpectra.h:277
AliAnalysisTaskHJetSpectra::fhDeltaPtEmbPerp2D
TH2D * fhDeltaPtEmbPerp2D[kCAll][kTT][kRho]
embedded delta pT (emb track is perp to TT)
Definition: AliAnalysisTaskHJetSpectra.h:234
AliAnalysisTaskHJetSpectra::fHJetSpecGen
TH2D * fHJetSpecGen[kCAll][kTT][kRho]
TT associated spectrum of jets.
Definition: AliAnalysisTaskHJetSpectra.h:222
AliAnalysisTaskHJetSpectra::EstimateBgKT
Double_t EstimateBgKT(AliJetContainer *jetCont, AliParticleContainer *trkArray, AliVParticle *trackTT)
Definition: AliAnalysisTaskHJetSpectra.cxx:2706
AliAnalysisTaskHJetSpectra::fhVertexZAcceptMC
TH1F * fhVertexZAcceptMC
gc vertexZ inclusive in MC
Definition: AliAnalysisTaskHJetSpectra.h:258
AliAnalysisTaskHJetSpectra::fhInvPtQVsPhiCSide
TH2D * fhInvPtQVsPhiCSide[2]
q*1/pT versus eta
Definition: AliAnalysisTaskHJetSpectra.h:305
AliParticleContainer.h
AliAnalysisTaskHJetSpectra::fHelperClass
AliAnalysisUtils * fHelperClass
A random number.
Definition: AliAnalysisTaskHJetSpectra.h:205
TH2D
Definition: External.C:228
TH1D
Definition: External.C:212
AliJetContainer.h
AliAnalysisTaskHJetSpectra::GetConePt
Double_t GetConePt(Double_t eta, Double_t phi, Double_t radius, AliParticleContainer *trkArray, Bool_t isGen)
Definition: AliAnalysisTaskHJetSpectra.cxx:401
AliAnalysisTaskHJetSpectra::fhDphiTTTT
TH1F * fhDphiTTTT[kTT]
Dphi of accepted jets after dphi cut.
Definition: AliAnalysisTaskHJetSpectra.h:262
AliAnalysisTaskHJetSpectra::fhTrackMultiplicityTT
TH1D * fhTrackMultiplicityTT[kCAll][kTT]
multiplicity of tracks
Definition: AliAnalysisTaskHJetSpectra.h:281
AliRhoParameter.h
AliAnalysisTaskHJetSpectra::fhJetPtGenVsJetPtRec
TH2D * fhJetPtGenVsJetPtRec[kCAll][kRho]
pt distribution of generator level jets
Definition: AliAnalysisTaskHJetSpectra.h:294
AliAnalysisTaskHJetSpectra::fhDphiTriggerJetGen
TH2F * fhDphiTriggerJetGen[kCAll][kTT][kRho]
gc Delta phi versus jet pT
Definition: AliAnalysisTaskHJetSpectra.h:260
AliJetContainer::GetNextAcceptJet
AliEmcalJet * GetNextAcceptJet()
Definition: AliJetContainer.cxx:299
Data
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Definition: charmCutsOptimization.C:371
AliAnalysisTaskHJetSpectra::RelativePhi
Double_t RelativePhi(Double_t mphi, Double_t vphi)
Definition: AliAnalysisTaskHJetSpectra.cxx:2610
AliAnalysisTaskHJetSpectra::fhVertexYAccept
TH1F * fhVertexYAccept
gc vertexZ accepted after vtx cut
Definition: AliAnalysisTaskHJetSpectra.h:252
AliAnalysisTaskHJetSpectra::fhZNAEnergyTT
TH1F * fhZNAEnergyTT[kCAll][kTT]
ZDC A neutral energy for given V0A centrality selection.
Definition: AliAnalysisTaskHJetSpectra.h:278
AliEmcalJet::Pt
Double_t Pt() const
Definition: AliEmcalJet.h:102
AliAnalysisTaskHJetSpectra::fhDCAinXVsPtStrange
TH2D * fhDCAinXVsPtStrange
Y DCA versus pT.
Definition: AliAnalysisTaskHJetSpectra.h:312
AliEmcalList
Enhanced TList-derived class that implements correct merging for pt_hard binned production.
Definition: AliEmcalList.h:25
AliAnalysisTaskHJetSpectra::fhJetEtaIncl
TH2F * fhJetEtaIncl
minimum bias phi inclusive
Definition: AliAnalysisTaskHJetSpectra.h:265
AliAnalysisTaskHJetSpectra::fhRhoTT
TH1F * fhRhoTT[kCAll][kTT][kRho]
TT associated spectrum of jets.
Definition: AliAnalysisTaskHJetSpectra.h:224
AliAnalysisTaskHJetSpectra::fhDphiTriggerJetAccept
TH1F * fhDphiTriggerJetAccept
gc Delta phi versus jet pT
Definition: AliAnalysisTaskHJetSpectra.h:261
AliAnalysisTaskHJetSpectra::fhJetPhiIncl
TH2F * fhJetPhiIncl
Dphi between multiple trigger tracks.
Definition: AliAnalysisTaskHJetSpectra.h:264
AliAnalysisTaskHJetSpectra::fHJetSpec
TH2D * fHJetSpec[kCAll][kTT][kRho]
trigger multiplicity in event
Definition: AliAnalysisTaskHJetSpectra.h:221
AliAnalysisTaskHJetSpectra::fImpParam
Double_t fImpParam
gc value is filled, if pythia header is accessible
Definition: AliAnalysisTaskHJetSpectra.h:201
AliAnalysisTaskEmcal::fOutput
AliEmcalList * fOutput
!output list
Definition: AliAnalysisTaskEmcal.h:357
AliAnalysisTaskEmcal.h
AliAnalysisTaskHJetSpectra::fhDiffPtVsPtTrackTrue
TH2D * fhDiffPtVsPtTrackTrue
pt spectrum of reconstructed jets without fake track pT vs reconstructed jet pT
Definition: AliAnalysisTaskHJetSpectra.h:300
AliEmcalJet::TrackAt
Short_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:153
AliAnalysisTaskHJetSpectra::fARhoTT
TH1F * fARhoTT[kCAll][kTT][kRho]
gc X=rho from perp cone, Y=centrality
Definition: AliAnalysisTaskHJetSpectra.h:227
AliAnalysisTaskHJetSpectra::fhTrackEta
TH2F * fhTrackEta[kCAll]
jet eta vs jet pT
Definition: AliAnalysisTaskHJetSpectra.h:247
AliAnalysisTaskHJetSpectra::GetDeltaPt
void GetDeltaPt(Int_t nrho, TArrayD &rho, Double_t *dpt, Double_t ttPhi, Double_t ttEta, AliParticleContainer *trkArray, Bool_t isGen)
Definition: AliAnalysisTaskHJetSpectra.cxx:637
AliAnalysisTaskHJetSpectra::fhVertexYAcceptTT
TH1F * fhVertexYAcceptTT
gc vertexZ accepted after vtx cut
Definition: AliAnalysisTaskHJetSpectra.h:255
AliAnalysisTaskEmcalJet
Base task in the EMCAL jet framework.
Definition: AliAnalysisTaskEmcalJet.h:30
AliAnalysisTaskHJetSpectra::fhInvPtQVsPhiASide
TH2D * fhInvPtQVsPhiASide[2]
q*1/pT versus eta
Definition: AliAnalysisTaskHJetSpectra.h:304
AliEmcalJet
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:44
Option_t
const char Option_t
Definition: External.C:48
AliAnalysisTaskHJetSpectra::fhDCAinYVsPtStrange
TH2D * fhDCAinYVsPtStrange
X DCA versus pT of strange tracks.
Definition: AliAnalysisTaskHJetSpectra.h:313
nbins
const Int_t nbins
Definition: AverageDmesonRaa.C:63
AliAnalysisTaskHJetSpectra::fhVzeroATotMult
TH1F * fhVzeroATotMult[kCAll]
centrality from ZNA
Definition: AliAnalysisTaskHJetSpectra.h:274
Bool_t
bool Bool_t
Definition: External.C:53
AliAnalysisTaskHJetSpectra::fhVertexZMC
TH1F * fhVertexZMC
gc vertexZ accepted after vtx cut
Definition: AliAnalysisTaskHJetSpectra.h:257
AliAnalysisTaskHJetSpectra::fhInvPtQVsEta
TH2D * fhInvPtQVsEta[2]
q*1/pT versus phi
Definition: AliAnalysisTaskHJetSpectra.h:303
AliAnalysisTaskHJetSpectra::fTrigTracks
AliVParticle * fTrigTracks[kTT][999]
Definition: AliAnalysisTaskHJetSpectra.h:331
AliAnalysisTaskHJetSpectra::fhImpactParameter
TH1D * fhImpactParameter[kCAll]
ZNA energy versus Vzero mult. versus track mult. in events with TT.
Definition: AliAnalysisTaskHJetSpectra.h:290
AliAnalysisTaskHJetSpectra::fhTrackPhiTPCG
TH2F * fhTrackPhiTPCG
hybrid constrained global track phi vs track pT
Definition: AliAnalysisTaskHJetSpectra.h:308
AliAnalysisTaskHJetSpectra::fhDCAinXVsPtNonStrange
TH2D * fhDCAinXVsPtNonStrange
Y DCA versus pT of strange tracks.
Definition: AliAnalysisTaskHJetSpectra.h:314
AliAnalysisTaskHJetSpectra::fhDeltaPtEmbBc2Bc2D
TH2D * fhDeltaPtEmbBc2Bc2D[kCAll][kTT][kRho]
embedded delta pT (emb track is back-to-back in azimuth to TT)
Definition: AliAnalysisTaskHJetSpectra.h:236
AliAnalysisTaskHJetSpectra::fhCentralityZNA
TH1F * fhCentralityZNA
centrality from V0C
Definition: AliAnalysisTaskHJetSpectra.h:272
AliAnalysisTaskHJetSpectra::fhDCAinYVsPt
TH2D * fhDCAinYVsPt
X DCA versus pT.
Definition: AliAnalysisTaskHJetSpectra.h:311
AliAnalysisTaskHJetSpectra::fInitializedLocal
Bool_t fInitializedLocal
gc Vertex selection helper
Definition: AliAnalysisTaskHJetSpectra.h:206
AliEmcalJet.h
AliAnalysisTaskHJetSpectra::fhPtTrkTruePrimRec
TH2D * fhPtTrkTruePrimRec[kCAll]
pt jet resolution
Definition: AliAnalysisTaskHJetSpectra.h:296
AliJetContainer
Container for jet within the EMCAL jet framework.
Definition: AliJetContainer.h:37
TH1
Definition: External.C:196
AliAnalysisTaskHJetSpectra::fhImpactParameterTT
TH1D * fhImpactParameterTT[kCAll][kTT]
impact parameter distribution hijing
Definition: AliAnalysisTaskHJetSpectra.h:291
AliAnalysisTaskHJetSpectra::fhCentralityV0M
TH1F * fhCentralityV0M
centrality V0 multiplicity A+C when TT is present
Definition: AliAnalysisTaskHJetSpectra.h:269