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