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