PWGJE/AliAnalysisTaskJetCore.cxx

   1 // ******************************************
   2 // This task computes several jet observables like
   3 // the fraction of energy in inner and outer coronnas,
   4 // jet-track correlations,triggered jet shapes and
   5 // correlation strength distribution of particles inside jets.
   6 // Author: lcunquei@cern.ch
   7 // *******************************************
   8
   9
  10 /**************************************************************************
  11  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  12  *                                                                        *
  13  * Author: The ALICE Off-line Project.                                    *
  14  * Contributors are mentioned in the code where appropriate.              *
  15  *                                                                        *
  16  * Permission to use, copy, modify and distribute this software and its   *
  17  * documentation strictly for non-commercial purposes is hereby granted   *
  18  * without fee, provided that the above copyright notice appears in all   *
  19  * copies and that both the copyright notice and this permission notice   *
  20  * appear in the supporting documentation. The authors make no claims     *
  21  * about the suitability of this software for any purpose. It is          *
  22  * provided "as is" without express or implied warranty.                  *
  23  **************************************************************************/
  24
  25
  26 #include "TChain.h"
  27 #include "TTree.h"
  28 #include "TMath.h"
  29 #include "TH1F.h"
  30 #include "TH2F.h"
  31 #include "TH3F.h"
  32 #include "THnSparse.h"
  33 #include "TCanvas.h"
  34 #include "TRandom3.h"
  35 #include "AliLog.h"
  36
  37 #include "AliAnalysisTask.h"
  38 #include "AliAnalysisManager.h"
  39
  40 #include "AliVEvent.h"
  41 #include "AliESDEvent.h"
  42 #include "AliESDInputHandler.h"
  43 #include "AliCentrality.h"
  44 #include "AliAnalysisHelperJetTasks.h"
  45 #include "AliInputEventHandler.h"
  46 #include "AliAODJetEventBackground.h"
  47 #include "AliAODMCParticle.h"
  48 #include "AliAnalysisTaskFastEmbedding.h"
  49 #include "AliAODEvent.h"
  50 #include "AliAODHandler.h"
  51 #include "AliAODJet.h"
  52
  53 #include "AliAnalysisTaskJetCore.h"
  54
  55 using std::cout;
  56 using std::endl;
  57
  58 ClassImp(AliAnalysisTaskJetCore)
  59
  60 AliAnalysisTaskJetCore::AliAnalysisTaskJetCore() :
  61 AliAnalysisTaskSE(),
  62 fESD(0x0),
  63 fAODIn(0x0),
  64 fAODOut(0x0),
  65 fAODExtension(0x0),
  66 fBackgroundBranch(""),
  67 fNonStdFile(""),
  68 fIsPbPb(kTRUE),
  69 fOfflineTrgMask(AliVEvent::kAny),
  70 fMinContribVtx(1),
  71 fVtxZMin(-10.),
  72 fVtxZMax(10.),
  73 fEvtClassMin(0),
  74 fEvtClassMax(4),
  75 fFilterMask(0),
  76 fFilterMaskBestPt(0),
  77 fFilterType(0),
  78 fRadioFrac(0.2),
  79 fMinDist(0.1),
  80 fCentMin(0.),
  81 fCentMax(100.),
  82 fNInputTracksMin(0),
  83 fNInputTracksMax(-1),
  84 fRequireITSRefit(0),
  85 fApplySharedClusterCut(0),
  86 fAngStructCloseTracks(0),
  87 fCheckMethods(0),
  88 fDoEventMixing(0),
  89 fFlagPhiBkg(0),
  90 fFlagEtaBkg(0),
  91 fFlagJetHadron(0),
  92 fDodiHadron(0),
  93 fFrac(0.8),
  94 fTTLowRef(11),
  95 fTTUpRef(13),
  96 fTTLowSig(15),
  97 fTTUpSig(19),
  98 fHardest(0),
  99 fFlagRandom(0),
 100 fFlagOnlyRecoil(0),
 101 fFlagOnlyHardest(1),
 102 fTrackTypeRec(kTrackUndef),
 103 fRPAngle(0),
 104 fNRPBins(50),
 105 fSemigoodCorrect(0),
 106 fHolePos(4.71),
 107 fHoleWidth(0.2),
 108 fJetEtaMin(-.5),
 109 fJetEtaMax(.5),
 110 fNevents(0),
 111 fTindex(0),
 112 fTrigBufferIndex(0),
 113 fCountAgain(0),
 114 fJetPtMin(20.),
 115 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
 116 fJetPtFractionMin(0.5),
 117 fNMatchJets(4),
 118 fMatchMaxDist(0.8),
 119 fKeepJets(kFALSE),
 120 fRunAnaAzimuthalCorrelation(kFALSE),
 121 fkNbranches(2),
 122 fkEvtClasses(12),
 123 fOutputList(0x0),
 124 fbEvent(kTRUE),
 125 fHistEvtSelection(0x0),
 126 fhnDeltaR(0x0),
 127 fhnMixedEvents(0x0),
 128 fh2JetCoreMethod1C10(0x0),
 129 fh2JetCoreMethod2C10(0x0),
 130 fh2JetCoreMethod1C20(0x0),
 131 fh2JetCoreMethod2C20(0x0),
 132 fh2JetCoreMethod1C30(0x0),
 133 fh2JetCoreMethod2C30(0x0),
 134 fh2JetCoreMethod1C60(0x0),
 135 fh2JetCoreMethod2C60(0x0),
 136 fh3JetTrackC3060(0x0),
 137 fh3JetTrackC20(0x0),
 138 fh2AngStructpt1C10(0x0),
 139 fh2AngStructpt2C10(0x0),
 140 fh2AngStructpt3C10(0x0),
 141 fh2AngStructpt4C10(0x0),
 142 fh2AngStructpt1C20(0x0),
 143 fh2AngStructpt2C20(0x0),
 144 fh2AngStructpt3C20(0x0),
 145 fh2AngStructpt4C20(0x0),
 146 fh2AngStructpt1C30(0x0),
 147 fh2AngStructpt2C30(0x0),
 148 fh2AngStructpt3C30(0x0),
 149 fh2AngStructpt4C30(0x0),
 150 fh2AngStructpt1C60(0x0),
 151 fh2AngStructpt2C60(0x0),
 152 fh2AngStructpt3C60(0x0),
 153 fh2AngStructpt4C60(0x0),
 154 fh1TrigRef(0x0),
 155 fh1TrigSig(0x0),
 156 fh1TrackPhiDistance(0x0),
 157 fh1TrackRDistance(0x0),
 158 fh2Ntriggers(0x0),
 159 fh2Ntriggers2C10(0x0),
 160 fh2Ntriggers2C20(0x0),
 161 fh3JetDensity(0x0),
 162 fh3JetDensityA4(0x0),
 163 fh2RPJetsC10(0x0),
 164 fh2RPJetsC20(0x0),
 165 fh2RPTC10(0x0),
 166 fh2RPTC20(0x0),
 167 fHJetSpec(0x0),
 168 fhTTPt(0x0),
 169 fHJetPhiCorr(0x0),
 170 fRandom(0x0)
 171
 172 {
 173    // default Constructor
 174
 175
 176  // Trigger buffer.
 177    for(Int_t i=0; i<10; i++) {
 178                 for(Int_t j=0; j<6; j++) {
 179                         fTrigBuffer[i][j]=0;
 180                 }
 181         }
 182
 183
 184
 185
 186
 187    fJetBranchName[0] = "";
 188    fJetBranchName[1] = "";
 189
 190    fListJets[0] = new TList;
 191    fListJets[1] = new TList;
 192 }
 193
 194 AliAnalysisTaskJetCore::AliAnalysisTaskJetCore(const char *name) :
 195 AliAnalysisTaskSE(name),
 196 fESD(0x0),
 197 fAODIn(0x0),
 198 fAODOut(0x0),
 199 fAODExtension(0x0),
 200 fBackgroundBranch(""),
 201 fNonStdFile(""),
 202 fIsPbPb(kTRUE),
 203 fOfflineTrgMask(AliVEvent::kAny),
 204 fMinContribVtx(1),
 205 fVtxZMin(-10.),
 206 fVtxZMax(10.),
 207 fEvtClassMin(0),
 208 fEvtClassMax(4),
 209 fFilterMask(0),
 210 fFilterMaskBestPt(0),
 211 fFilterType(0),
 212 fRadioFrac(0.2),
 213 fMinDist(0.1),
 214 fCentMin(0.),
 215 fCentMax(100.),
 216 fNInputTracksMin(0),
 217 fNInputTracksMax(-1),
 218 fRequireITSRefit(0),
 219 fApplySharedClusterCut(0),
 220 fAngStructCloseTracks(0),
 221 fCheckMethods(0),
 222 fDoEventMixing(0),
 223 fFlagPhiBkg(0),
 224 fFlagEtaBkg(0),
 225 fFlagJetHadron(0),
 226 fDodiHadron(0),
 227 fFrac(0.8),
 228 fTTLowRef(11),
 229 fTTUpRef(13),
 230 fTTLowSig(15),
 231 fTTUpSig(19),
 232 fHardest(0),
 233 fFlagRandom(0),
 234 fFlagOnlyRecoil(0),
 235 fFlagOnlyHardest(1),
 236 fTrackTypeRec(kTrackUndef),
 237 fRPAngle(0),
 238 fNRPBins(50),
 239 fSemigoodCorrect(0),
 240 fHolePos(4.71),
 241 fHoleWidth(0.2),
 242 fJetEtaMin(-.5),
 243 fJetEtaMax(.5),
 244 fNevents(0),
 245 fTindex(0),
 246 fTrigBufferIndex(0),
 247 fCountAgain(0),
 248 fJetPtMin(20.),
 249 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
 250 fJetPtFractionMin(0.5),
 251 fNMatchJets(4),
 252 fMatchMaxDist(0.8),
 253 fKeepJets(kFALSE),
 254 fRunAnaAzimuthalCorrelation(kFALSE),
 255 fkNbranches(2),
 256 fkEvtClasses(12),
 257 fOutputList(0x0),
 258 fbEvent(kTRUE),
 259 fHistEvtSelection(0x0),
 260 fhnDeltaR(0x0),
 261 fhnMixedEvents(0x0),
 262 fh2JetCoreMethod1C10(0x0),
 263 fh2JetCoreMethod2C10(0x0),
 264 fh2JetCoreMethod1C20(0x0),
 265 fh2JetCoreMethod2C20(0x0),
 266 fh2JetCoreMethod1C30(0x0),
 267 fh2JetCoreMethod2C30(0x0),
 268 fh2JetCoreMethod1C60(0x0),
 269 fh2JetCoreMethod2C60(0x0),
 270 fh3JetTrackC3060(0x0),
 271 fh3JetTrackC20(0x0),
 272 fh2AngStructpt1C10(0x0),
 273 fh2AngStructpt2C10(0x0),
 274 fh2AngStructpt3C10(0x0),
 275 fh2AngStructpt4C10(0x0),
 276 fh2AngStructpt1C20(0x0),
 277 fh2AngStructpt2C20(0x0),
 278 fh2AngStructpt3C20(0x0),
 279 fh2AngStructpt4C20(0x0),
 280 fh2AngStructpt1C30(0x0),
 281 fh2AngStructpt2C30(0x0),
 282 fh2AngStructpt3C30(0x0),
 283 fh2AngStructpt4C30(0x0),
 284 fh2AngStructpt1C60(0x0),
 285 fh2AngStructpt2C60(0x0),
 286 fh2AngStructpt3C60(0x0),
 287 fh2AngStructpt4C60(0x0),
 288 fh1TrigRef(0x0),
 289 fh1TrigSig(0x0),
 290 fh1TrackPhiDistance(0x0),
 291 fh1TrackRDistance(0x0),
 292 fh2Ntriggers(0x0),
 293 fh2Ntriggers2C10(0x0),
 294 fh2Ntriggers2C20(0x0),
 295 fh3JetDensity(0x0),
 296 fh3JetDensityA4(0x0),
 297 fh2RPJetsC10(0x0),
 298 fh2RPJetsC20(0x0),
 299 fh2RPTC10(0x0),
 300 fh2RPTC20(0x0),
 301 fHJetSpec(0x0),
 302 fhTTPt(0x0),
 303 fHJetPhiCorr(0x0),
 304 fRandom(0x0)
 305
 306  {
 307    // Constructor
 308
 309
 310     for(Int_t i=0; i<10; i++) {
 311        for(Int_t j=0; j<6; j++) {
 312             fTrigBuffer[i][j]=0;
 313                 }
 314     }
 315
 316
 317
 318    fJetBranchName[0] = "";
 319    fJetBranchName[1] = "";
 320
 321    fListJets[0] = new TList;
 322    fListJets[1] = new TList;
 323
 324    DefineOutput(1, TList::Class());
 325 }
 326
 327 AliAnalysisTaskJetCore::~AliAnalysisTaskJetCore()
 328 {
 329    delete fListJets[0];
 330    delete fListJets[1];
 331 }
 332
 333 void AliAnalysisTaskJetCore::SetBranchNames(const TString &branch1, const TString &branch2)
 334 {
 335    fJetBranchName[0] = branch1;
 336    fJetBranchName[1] = branch2;
 337 }
 338
 339 void AliAnalysisTaskJetCore::Init()
 340 {
 341
 342    // check for jet branches
 343    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 344       AliError("Jet branch name not set.");
 345    }
 346
 347 }
 348
 349 void AliAnalysisTaskJetCore::UserCreateOutputObjects()
 350 {
 351    // Create histograms
 352    // Called once
 353    OpenFile(1);
 354    if(!fOutputList) fOutputList = new TList;
 355    fOutputList->SetOwner(kTRUE);
 356
 357    Bool_t oldStatus = TH1::AddDirectoryStatus();
 358    TH1::AddDirectory(kFALSE);
 359
 360
 361         // set seed
 362    fRandom = new TRandom3(0);
 363
 364
 365
 366
 367    fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
 368    fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
 369    fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
 370    fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
 371    fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
 372    fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
 373    fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
 374
 375      UInt_t entries = 0; // bit coded, see GetDimParams() below
 376      entries = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 |1<<7;
 377      fhnDeltaR = NewTHnSparseF("fhnDeltaR", entries);
 378
 379      //change binning in pTtrack
 380      Double_t *xPt3 = new Double_t[10];
 381      xPt3[0] = 0.;
 382      for(int i = 1; i<=9;i++){
 383       if(xPt3[i-1]<2)xPt3[i] = xPt3[i-1] + 0.4; // 1 - 5
 384       else if(xPt3[i-1]<11)xPt3[i] = xPt3[i-1] + 3; // 5 - 12
 385       else xPt3[i] = xPt3[i-1] + 150.; // 18
 386      }
 387     fhnDeltaR->SetBinEdges(2,xPt3);
 388     delete [] xPt3;
 389
 390     //change binning in HTI
 391      Double_t *xPt4 = new Double_t[14];
 392      xPt4[0] = 0.;
 393      for(int i = 1; i<=13;i++){
 394       if(xPt4[i-1]<10)xPt4[i] = xPt4[i-1] + 1; // 1 - 10
 395       else if(xPt4[i-1]<20)xPt4[i] = xPt4[i-1] + 5; // 10 - 12
 396       else xPt4[i] = xPt4[i-1] + 30.; // 13
 397      }
 398     fhnDeltaR->SetBinEdges(6,xPt4);
 399     delete [] xPt4;
 400
 401
 402
 403
 404
 405
 406
 407      if(fDoEventMixing){
 408      UInt_t cifras = 0; // bit coded, see GetDimParams() below
 409      cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7;
 410      fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras);}
 411
 412     if(fCheckMethods){
 413     fh2JetCoreMethod1C10 = new TH2F("JetCoreMethod1C10","",150, 0., 150.,100, 0., 1.5);
 414     fh2JetCoreMethod2C10 = new TH2F("JetCoreMethod2C10","",150, 0., 150.,100, 0., 1.5);
 415     fh2JetCoreMethod1C20 = new TH2F("JetCoreMethod1C20","",150, 0., 150.,100, 0., 1.5);
 416     fh2JetCoreMethod2C20 = new TH2F("JetCoreMethod2C20","",150, 0., 150.,100, 0., 1.5);
 417     fh2JetCoreMethod1C30 = new TH2F("JetCoreMethod1C30","",150, 0., 150.,100, 0., 1.5);
 418     fh2JetCoreMethod2C30 = new TH2F("JetCoreMethod2C30","",150, 0., 150.,100, 0., 1.5);
 419     fh2JetCoreMethod1C60 = new TH2F("JetCoreMethod1C60","",150, 0., 150.,100, 0., 1.5);
 420     fh2JetCoreMethod2C60 = new TH2F("JetCoreMethod2C60","",150, 0., 150.,100, 0., 1.5);}
 421
 422     fh3JetTrackC3060=new TH3F("JetTrackC3060","",50,0,50,150,0.,150.,35,0.,3.5);
 423     fh3JetTrackC20=new TH3F("JetTrackC20","",50,0,50,150,0.,150.,35,0.,3.5);
 424     if(fAngStructCloseTracks>0){
 425     fh2AngStructpt1C10 = new TH2F("Ang struct pt1 C10","",15,0.,1.5,150,0.,10.);
 426     fh2AngStructpt2C10 = new TH2F("Ang struct pt2 C10","",15,0.,1.5,150,0.,10.);
 427     fh2AngStructpt3C10 = new TH2F("Ang struct pt3 C10","",15,0.,1.5,150,0.,10.);
 428     fh2AngStructpt4C10 = new TH2F("Ang struct pt4 C10","",15,0.,1.5,150,0.,10.);
 429     fh2AngStructpt1C20 = new TH2F("Ang struct pt1 C20","",15,0.,1.5,150,0.,10.);
 430     fh2AngStructpt2C20 = new TH2F("Ang struct pt2 C20","",15,0.,1.5,150,0.,10.);
 431     fh2AngStructpt3C20 = new TH2F("Ang struct pt3 C20","",15,0.,1.5,150,0.,10.);
 432     fh2AngStructpt4C20 = new TH2F("Ang struct pt4 C20","",15,0.,1.5,150,0.,10.);
 433     fh2AngStructpt1C30 = new TH2F("Ang struct pt1 C30","",15,0.,1.5,150,0.,10.);
 434     fh2AngStructpt2C30 = new TH2F("Ang struct pt2 C30","",15,0.,1.5,150,0.,10.);
 435     fh2AngStructpt3C30 = new TH2F("Ang struct pt3 C30","",15,0.,1.5,150,0.,10.);
 436     fh2AngStructpt4C30 = new TH2F("Ang struct pt4 C30","",15,0.,1.5,150,0.,10.);
 437     fh2AngStructpt1C60 = new TH2F("Ang struct pt1 C60","",15,0.,1.5,150,0.,10.);
 438     fh2AngStructpt2C60 = new TH2F("Ang struct pt2 C60","",15,0.,1.5,150,0.,10.);
 439     fh2AngStructpt3C60 = new TH2F("Ang struct pt3 C60","",15,0.,1.5,150,0.,10.);
 440     fh2AngStructpt4C60 = new TH2F("Ang struct pt4 C60","",15,0.,1.5,150,0.,10.); }
 441
 442
 443     fh1TrigRef=new TH1D("Trig Ref","",10,0.,10);
 444     fh1TrigSig=new TH1D("Trig Sig","",10,0.,10);
 445     fh1TrackPhiDistance=new TH1D("PhiDistance","",35,0.,3.5);
 446     fh1TrackRDistance=new TH1D("RDistance","",10,0.,1);
 447     fh2Ntriggers=new TH2F("# of triggers","",100,0.,100.,50,0.,50.);
 448     fh2Ntriggers2C10=new TH2F("# of triggers2C10","",50,0.,50.,50,0.,50.);
 449     fh2Ntriggers2C20=new TH2F("# of triggers2C20","",50,0.,50.,50,0.,50.);
 450     fh3JetDensity=new TH3F("Jet density vs mutliplicity A>0.07","",100,0.,4000.,100,0.,5.,25,0.,50.);
 451     fh3JetDensityA4=new TH3F("Jet density vs multiplicity A>0.4","",100,0.,4000.,100,0.,5.,25,0.,50.);
 452     fh2RPJetsC10=new TH2F("RPJetC10","",35,0.,3.5,100,0.,100.);
 453     fh2RPJetsC20=new TH2F("RPJetC20","",35,0.,3.5,100,0.,100.);
 454     fh2RPTC10=new TH2F("RPTriggerC10","",35,0.,3.5,50,0.,50.);
 455     fh2RPTC20=new TH2F("RPTriggerC20","",35,0.,3.5,50,0.,50.);
 456
 457
 458
 459
 460    fOutputList->Add(fHistEvtSelection);
 461
 462    fOutputList->Add(fhnDeltaR);
 463
 464    fOutputList->Add(fhnMixedEvents);
 465
 466
 467
 468       if(fCheckMethods){
 469
 470       fOutputList->Add(fh2JetCoreMethod1C10);
 471       fOutputList->Add(fh2JetCoreMethod2C10);
 472       fOutputList->Add(fh2JetCoreMethod1C20);
 473       fOutputList->Add(fh2JetCoreMethod2C20);
 474       fOutputList->Add(fh2JetCoreMethod1C30);
 475       fOutputList->Add(fh2JetCoreMethod2C30);
 476       fOutputList->Add(fh2JetCoreMethod1C60);
 477       fOutputList->Add(fh2JetCoreMethod2C60);}
 478
 479       fOutputList->Add(fh3JetTrackC3060);
 480       fOutputList->Add(fh3JetTrackC20);
 481
 482
 483
 484
 485
 486         if(fAngStructCloseTracks>0){
 487        fOutputList->Add(fh2AngStructpt1C10);
 488        fOutputList->Add(fh2AngStructpt2C10);
 489        fOutputList->Add(fh2AngStructpt3C10);
 490        fOutputList->Add(fh2AngStructpt4C10);
 491        fOutputList->Add(fh2AngStructpt1C20);
 492        fOutputList->Add(fh2AngStructpt2C20);
 493        fOutputList->Add(fh2AngStructpt3C20);
 494        fOutputList->Add(fh2AngStructpt4C20);
 495        fOutputList->Add(fh2AngStructpt1C30);
 496        fOutputList->Add(fh2AngStructpt2C30);
 497        fOutputList->Add(fh2AngStructpt3C30);
 498        fOutputList->Add(fh2AngStructpt4C30);
 499        fOutputList->Add(fh2AngStructpt1C60);
 500        fOutputList->Add(fh2AngStructpt2C60);
 501        fOutputList->Add(fh2AngStructpt3C60);
 502        fOutputList->Add(fh2AngStructpt4C60);}
 503
 504
 505
 506
 507         fOutputList->Add(fh1TrigRef);
 508         fOutputList->Add(fh1TrigSig);
 509         fOutputList->Add(fh1TrackPhiDistance);
 510         fOutputList->Add(fh1TrackRDistance);
 511         fOutputList->Add(fh2Ntriggers);
 512         fOutputList->Add(fh2Ntriggers2C10);
 513         fOutputList->Add(fh2Ntriggers2C20);
 514         fOutputList->Add(fh3JetDensity);
 515         fOutputList->Add(fh3JetDensityA4);
 516         fOutputList->Add(fh2RPJetsC10);
 517         fOutputList->Add(fh2RPJetsC20);
 518          fOutputList->Add(fh2RPTC10);
 519         fOutputList->Add(fh2RPTC20);
 520
 521         const Int_t dimSpec = 5;
 522         const Int_t nBinsSpec[dimSpec]     = {100,6, 140, 50, fNRPBins};
 523         const Double_t lowBinSpec[dimSpec] = {0,0,-80, 0, 0};
 524         const Double_t hiBinSpec[dimSpec]  = {100,1, 200, 50,  static_cast<Double_t>(fNRPBins)};
 525         fHJetSpec = new THnSparseF("fHJetSpec","Recoil jet spectrum",dimSpec,nBinsSpec,lowBinSpec,hiBinSpec);
 526
 527              //change binning in jet area
 528      Double_t *xPt6 = new Double_t[7];
 529      xPt6[0] = 0.;
 530      xPt6[1]=0.07;
 531      xPt6[2]=0.2;
 532      xPt6[3]=0.4;
 533      xPt6[4]=0.6;
 534      xPt6[5]=0.8;
 535      xPt6[6]=1;
 536     fHJetSpec->SetBinEdges(1,xPt6);
 537     delete [] xPt6;
 538
 539
 540
 541
 542
 543         fOutputList->Add(fHJetSpec);
 544
 545
 546         if(fRunAnaAzimuthalCorrelation)
 547           {
 548             fhTTPt = new TH2F("fhTTPt","Trigger track p_{T} vs centrality",10,0,100,100,0,100);
 549             fOutputList->Add(fhTTPt);
 550
 551             const Int_t dimCor = 5;
 552             const Int_t nBinsCor[dimCor]     = {50, 200, 100,              8,   10};
 553             const Double_t lowBinCor[dimCor] = {0,  -50, -0.5*TMath::Pi(), 0,   0};
 554             const Double_t hiBinCor[dimCor]  = {50, 150, 1.5*TMath::Pi(),  0.8, 100};
 555             fHJetPhiCorr = new THnSparseF("fHJetPhiCorr","TT p_{T} vs jet p_{T} vs dPhi vs area vs centrality",dimCor,nBinsCor,lowBinCor,hiBinCor);
 556             fOutputList->Add(fHJetPhiCorr);
 557           }
 558
 559
 560
 561
 562    // =========== Switch on Sumw2 for all histos ===========
 563    for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
 564       TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
 565       if (h1){
 566          h1->Sumw2();
 567          continue;
 568       }
 569     THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
 570       if (hn){
 571          hn->Sumw2();
 572       }
 573    }
 574    TH1::AddDirectory(oldStatus);
 575
 576    PostData(1, fOutputList);
 577 }
 578
 579 void AliAnalysisTaskJetCore::UserExec(Option_t *)
 580 {
 581
 582
 583    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 584       AliError("Jet branch name not set.");
 585       return;
 586    }
 587
 588    fESD=dynamic_cast<AliESDEvent*>(InputEvent());
 589    if (!fESD) {
 590       AliError("ESD not available");
 591       fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
 592    }
 593       fAODOut = dynamic_cast<AliAODEvent*>(AODEvent());
 594
 595        static AliAODEvent* aod = 0;
 596        // take all other information from the aod we take the tracks from
 597        if(!aod){
 598        if(!fESD)aod = fAODIn;
 599        else aod = fAODOut;}
 600
 601
 602
 603     if(fNonStdFile.Length()!=0){
 604     // case that we have an AOD extension we need can fetch the jets from the extended output
 605     AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 606     fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
 607     if(!fAODExtension){
 608       if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
 609     }}
 610
 611
 612    // -- event selection --
 613    fHistEvtSelection->Fill(1); // number of events before event selection
 614
 615
 616        Bool_t selected=kTRUE;
 617        selected = AliAnalysisHelperJetTasks::Selected();
 618        if(!selected){
 619       // no selection by the service task, we continue
 620        PostData(1,fOutputList);
 621        return;}
 622
 623
 624
 625   // physics selection: this is now redundant, all should appear as accepted after service task selection
 626    AliInputEventHandler* inputHandler = (AliInputEventHandler*)
 627    ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
 628          std::cout<<inputHandler->IsEventSelected()<<" "<<fOfflineTrgMask<<std::endl;
 629    if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
 630       if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
 631       fHistEvtSelection->Fill(2);
 632       PostData(1, fOutputList);
 633       return;
 634    }
 635
 636
 637
 638    // vertex selection
 639    if(!aod){
 640      if(fDebug) Printf("%s:%d No AOD",(char*)__FILE__,__LINE__);
 641      fHistEvtSelection->Fill(3);
 642       PostData(1, fOutputList);
 643    }
 644    AliAODVertex* primVtx = aod->GetPrimaryVertex();
 645
 646    if(!primVtx){
 647      if(fDebug) Printf("%s:%d No primVtx",(char*)__FILE__,__LINE__);
 648      fHistEvtSelection->Fill(3);
 649      PostData(1, fOutputList);
 650      return;
 651    }
 652
 653    Int_t nTracksPrim = primVtx->GetNContributors();
 654    if ((nTracksPrim < fMinContribVtx) ||
 655          (primVtx->GetZ() < fVtxZMin) ||
 656          (primVtx->GetZ() > fVtxZMax) ){
 657       if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
 658       fHistEvtSelection->Fill(3);
 659       PostData(1, fOutputList);
 660       return;
 661    }
 662
 663    // event class selection (from jet helper task)
 664    Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
 665    if(fDebug) Printf("Event class %d", eventClass);
 666    if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
 667       fHistEvtSelection->Fill(4);
 668       PostData(1, fOutputList);
 669       return;
 670    }
 671
 672    // centrality selection
 673    AliCentrality *cent = 0x0;
 674    Double_t centValue = 0.;
 675    if(fIsPbPb){
 676    if(fESD) {cent = fESD->GetCentrality();
 677      if(cent) centValue = cent->GetCentralityPercentile("V0M");}
 678    else     centValue=((AliVAODHeader*)aod->GetHeader())->GetCentrality();
 679
 680    if(fDebug) printf("centrality: %f\n", centValue);
 681       if (centValue < fCentMin || centValue > fCentMax){
 682       fHistEvtSelection->Fill(4);
 683       PostData(1, fOutputList);
 684       return;
 685       }}
 686
 687
 688    fHistEvtSelection->Fill(0);
 689    // accepted events
 690    // -- end event selection --
 691
 692    // get background
 693    AliAODJetEventBackground* externalBackground = 0;
 694    if(fAODOut&&!externalBackground&&fBackgroundBranch.Length()){
 695       externalBackground =  (AliAODJetEventBackground*)(fAODOut->FindListObject(fBackgroundBranch.Data()));
 696       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 697    }
 698    if(fAODExtension&&!externalBackground&&fBackgroundBranch.Length()){
 699      externalBackground =  (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data()));
 700       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 701    }
 702
 703     if(fAODIn&&!externalBackground&&fBackgroundBranch.Length()){
 704       externalBackground =  (AliAODJetEventBackground*)(fAODIn->FindListObject(fBackgroundBranch.Data()));
 705       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 706     }
 707
 708    Float_t rho = 0;
 709
 710    if(fFlagRandom==0){
 711      if(externalBackground)rho = externalBackground->GetBackground(0);}
 712    if(fFlagRandom==2){
 713       if(externalBackground)rho = externalBackground->GetBackground(2);}
 714    if(fFlagRandom==3){
 715       if(externalBackground)rho = externalBackground->GetBackground(3);}
 716    // fetch jets
 717    TClonesArray *aodJets[2];
 718    aodJets[0]=0;
 719    if(fAODOut&&!aodJets[0]){
 720    aodJets[0] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[0].Data()));
 721    aodJets[1] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[1].Data()));  }
 722    if(fAODExtension && !aodJets[0]){
 723    aodJets[0] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[0].Data()));
 724    aodJets[1] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[1].Data()));  }
 725      if(fAODIn&&!aodJets[0]){
 726    aodJets[0] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[0].Data()));
 727    aodJets[1] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[1].Data()));  }
 728
 729
 730
 731    Int_t nT=0;
 732    TList ParticleList;
 733    Double_t minT=0;
 734    Double_t maxT=0;
 735    Int_t number=0;
 736    Double_t dice=fRandom->Uniform(0,1);
 737    if(dice>fFrac){ minT=fTTLowRef;
 738                    maxT=fTTUpRef;}
 739    if(dice<=fFrac){minT=fTTLowSig;
 740                    maxT=fTTUpSig;}
 741
 742
 743
 744    if(fHardest==1 || fHardest==2) nT = GetListOfTracks(&ParticleList);
 745    if(fHardest==0) nT=SelectTrigger(&ParticleList,minT,maxT,number);
 746    if(nT<0){
 747    PostData(1, fOutputList);
 748    return;}
 749
 750       if(dice>fFrac) fh1TrigRef->Fill(number);
 751       if(dice<=fFrac)fh1TrigSig->Fill(number);
 752
 753      for (Int_t iJetType = 0; iJetType < 2; iJetType++) {
 754       fListJets[iJetType]->Clear();
 755       if (!aodJets[iJetType]) continue;
 756       if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
 757       for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
 758          AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
 759          if (jet) fListJets[iJetType]->Add(jet);
 760          // if(iJetType==0){
 761          // ptsub[iJet]=jet->Pt()-rho*jet->EffectiveAreaCharged();}
 762       }}
 763
 764    Double_t etabig=0;
 765    Double_t ptbig=0;
 766    Double_t areabig=0;
 767    Double_t phibig=0.;
 768    Double_t etasmall=0;
 769    Double_t ptsmall=0;
 770    //   Double_t areasmall=0;
 771    Double_t phismall=0.;
 772
 773
 774
 775    Int_t iCount=0;
 776    Int_t trigJet=-1;
 777    Int_t trigBBTrack=-1;
 778           // Int_t trigInTrack=-1;
 779    fRPAngle = ((AliVAODHeader*)aod->GetHeader())->GetEventplane();
 780
 781    if(fHardest==0 || fHardest==1){
 782    AliVParticle *partback = (AliVParticle*)ParticleList.At(nT);
 783    if(!partback){
 784    PostData(1, fOutputList);
 785    return;}
 786
 787     if(fSemigoodCorrect){
 788    Double_t disthole=RelativePhi(partback->Phi(),fHolePos);
 789    if(TMath::Abs(disthole)+fHoleWidth>TMath::Pi()-0.6){
 790    PostData(1, fOutputList);
 791    return;}}
 792
 793
 794     }
 795
 796
 797    for(Int_t tt=0;tt<ParticleList.GetEntries();tt++){
 798      if(fHardest==0||fHardest==1){if(tt!=nT) continue;}
 799      AliVParticle *partback = (AliVParticle*)ParticleList.At(tt);
 800      if(!partback) continue;
 801      if(partback->Pt()<8) continue;
 802
 803    Double_t accep=2.*TMath::Pi()*1.8;
 804    Int_t injet4=0;
 805    Int_t injet=0;
 806
 807    if(fSemigoodCorrect){
 808    Double_t disthole=RelativePhi(partback->Phi(),fHolePos);
 809    if(TMath::Abs(disthole)+fHoleWidth>TMath::Pi()-0.6) continue;}
 810
 811
 812
 813    fh2Ntriggers->Fill(centValue,partback->Pt());
 814    Double_t phiBinT = RelativePhi(partback->Phi(),fRPAngle);
 815    if(centValue<20.) fh2RPTC20->Fill(TMath::Abs(phiBinT),partback->Pt());
 816    if(centValue<10.) fh2RPTC10->Fill(TMath::Abs(phiBinT),partback->Pt());
 817
 818
 819
 820    for(Int_t i=0; i<fListJets[0]->GetEntries(); ++i){
 821            AliAODJet* jetbig = (AliAODJet*)(fListJets[0]->At(i));
 822            etabig  = jetbig->Eta();
 823            phibig  = jetbig->Phi();
 824            ptbig   = jetbig->Pt();
 825            if(ptbig==0) continue;
 826            Double_t phiBin = RelativePhi(phibig,fRPAngle);
 827            areabig = jetbig->EffectiveAreaCharged();
 828            Double_t ptcorr=ptbig-rho*areabig;
 829            if((etabig<fJetEtaMin)||(etabig>fJetEtaMax)) continue;
 830            if(areabig>=0.07) injet=injet+1;
 831            if(areabig>=0.4) injet4=injet4+1;
 832            Double_t dphi=RelativePhi(partback->Phi(),phibig);
 833
 834            if(fFlagEtaBkg==1){
 835            Double_t etadif= partback->Eta()-etabig;
 836            if(TMath::Abs(etadif)<=0.5){
 837
 838            if(centValue<20.) fh3JetTrackC20->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));
 839            if(centValue>30. && centValue<60.) fh3JetTrackC3060->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));}}
 840            if(fFlagEtaBkg==0){
 841            if(centValue<20.) fh3JetTrackC20->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));
 842            if(centValue>30. && centValue<60.) fh3JetTrackC3060->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));}
 843
 844
 845            if(fFlagJetHadron==0){
 846            if(fFlagPhiBkg==1) if((TMath::Abs(dphi)<TMath::Pi()/2.-0.1)||(TMath::Abs(dphi)>TMath::Pi()/2.+0.1)) continue;
 847            if(fFlagPhiBkg==0) if(TMath::Abs(dphi)<TMath::Pi()-0.6) continue;
 848            if(fFlagPhiBkg==2) if(TMath::Abs(dphi)<TMath::Pi()-0.7) continue;
 849            if(fFlagPhiBkg==3) if(TMath::Abs(dphi)<TMath::Pi()-0.5) continue;}
 850
 851            if(fFlagJetHadron!=0) if(TMath::Abs(dphi)>0.4) continue;
 852
 853
 854            if(centValue<10.) fh2RPJetsC10->Fill(TMath::Abs(phiBin), ptcorr);
 855            if(centValue<20.) fh2RPJetsC20->Fill(TMath::Abs(phiBin), ptcorr);
 856                    Double_t dismin=100.;
 857                    Double_t ptmax=-10.;
 858                    Int_t index1=-1;
 859                    Int_t index2=-1;
 860
 861                    Float_t phitt=partback->Phi();
 862                    if(phitt<0)phitt+=TMath::Pi()*2.;
 863                    Int_t phiBintt = GetPhiBin(phitt-fRPAngle);
 864
 865                    Double_t fillspec[] = {centValue,jetbig->EffectiveAreaCharged(),ptcorr,partback->Pt(), static_cast<Double_t>(phiBintt)};
 866                   fHJetSpec->Fill(fillspec);
 867
 868
 869
 870                    if(ptcorr<=0) continue;
 871
 872                        AliAODTrack* leadtrack=0;
 873                        Int_t ippt=0;
 874                        Double_t ppt=-10;
 875                        if(fFlagJetHadron==0){
 876                          TRefArray *genTrackList = jetbig->GetRefTracks();
 877                          Int_t nTracksGenJet = genTrackList->GetEntriesFast();
 878                          AliAODTrack* genTrack;
 879                          for(Int_t ir=0; ir<nTracksGenJet; ++ir){
 880                            genTrack = (AliAODTrack*)(genTrackList->At(ir));
 881                            if(genTrack->Pt()>ppt){ppt=genTrack->Pt();
 882                              ippt=ir;}}
 883                          leadtrack=(AliAODTrack*)(genTrackList->At(ippt));
 884                          if(!leadtrack) continue;
 885                        }
 886
 887
 888
 889                        AliVParticle* leadtrackb=0;
 890                        if(fFlagJetHadron!=0){
 891                          Int_t nTb = GetHardestTrackBackToJet(jetbig);
 892                          leadtrackb = (AliVParticle*)ParticleList.At(nTb);
 893                          if(!leadtrackb) continue;
 894                        }
 895
 896
 897
 898
 899
 900                        //store one trigger info
 901                        if(iCount==0){
 902                          trigJet=i;
 903                          trigBBTrack=nT;
 904                          //                      trigInTrack=ippt;
 905                          iCount=iCount+1;}
 906
 907
 908                        if(fCheckMethods){
 909                          for(Int_t j=0; j<fListJets[1]->GetEntries(); ++j){
 910                            AliAODJet* jetsmall = (AliAODJet*)(fListJets[1]->At(j));
 911                            etasmall  = jetsmall->Eta();
 912                            phismall = jetsmall->Phi();
 913                            ptsmall   = jetsmall->Pt();
 914                            //                      areasmall = jetsmall->EffectiveAreaCharged();
 915                            Double_t tmpDeltaR=(phismall-phibig)*(phismall-phibig)+(etasmall-etabig)*(etasmall-etabig);
 916                            tmpDeltaR=TMath::Sqrt(tmpDeltaR);
 917                            //Fraction in the jet core
 918                            if((ptsmall>ptmax)&&(tmpDeltaR<=fRadioFrac)){ptmax=ptsmall;
 919                              index2=j;}
 920                     if(tmpDeltaR<=dismin){ dismin=tmpDeltaR;
 921                       index1=j;}} //en of loop over R=0.2 jets
 922                 //method1:most concentric jet=core
 923                 if(dismin<fMinDist){ AliAODJet* jetmethod1 = (AliAODJet*)(fListJets[1]->At(index1));
 924                   if(centValue<10) fh2JetCoreMethod1C10->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 925                   if((centValue>20)&&(centValue<40)) fh2JetCoreMethod1C20->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 926                   if((centValue>30)&&(centValue<60)) fh2JetCoreMethod1C30->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 927                   if(centValue>60) fh2JetCoreMethod1C60->Fill(ptcorr,jetmethod1->Pt()/ptbig); }
 928                 //method2:hardest contained jet=core
 929                 if(index2!=-1){
 930                   AliAODJet* jetmethod2 = (AliAODJet*)(fListJets[1]->At(index2));
 931                   if(centValue<10) fh2JetCoreMethod2C10->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 932                   if((centValue>20)&&(centValue<40)) fh2JetCoreMethod2C20->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 933                   if((centValue>30)&&(centValue<60)) fh2JetCoreMethod2C30->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 934                   if(centValue>60) fh2JetCoreMethod2C60->Fill(ptcorr,jetmethod2->Pt()/ptbig); }}
 935                   if(centValue<10&&leadtrack) fh2Ntriggers2C10->Fill(leadtrack->Pt(),partback->Pt());
 936                   if(centValue<20&&leadtrack) fh2Ntriggers2C20->Fill(leadtrack->Pt(),partback->Pt());
 937          if(fDoEventMixing==0 && fFlagOnlyRecoil==0){
 938          for(int it = 0;it<ParticleList.GetEntries();++it){
 939           AliVParticle *part = (AliVParticle*)ParticleList.At(it);
 940           Double_t deltaR = jetbig->DeltaR(part);
 941           Double_t deltaEta = etabig-part->Eta();
 942
 943           Double_t deltaPhi=phibig-part->Phi();
 944           if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
 945           if(deltaPhi>3./2.*TMath::Pi()) deltaPhi-=2.*TMath::Pi();
 946           Double_t pTcont=0;
 947           if(fFlagJetHadron==0) pTcont=leadtrack->Pt();
 948           if(fFlagJetHadron!=0) pTcont=leadtrackb->Pt();
 949            Double_t jetEntries[8] = {centValue,ptcorr,part->Pt(),deltaR,deltaEta,deltaPhi,pTcont,partback->Pt()};
 950            fhnDeltaR->Fill(jetEntries);}
 951
 952
 953           }
 954
 955
 956           //end of track loop, we only do it if EM is switched off
 957
 958
 959
 960
 961
 962
 963
 964
 965
 966    }
 967    if(injet>0) fh3JetDensity->Fill(ParticleList.GetEntries(),injet/accep,partback->Pt());
 968    if(injet4>0)fh3JetDensityA4->Fill(ParticleList.GetEntries(),injet4/accep,partback->Pt());
 969           //end of jet loop
 970
 971    //}
 972
 973
 974           if(fDoEventMixing>0){
 975             //check before if the trigger exists
 976             // fTrigBuffer[i][0] = zvtx
 977             // fTrigBuffer[i][1] = phi
 978             // fTrigBuffer[i][2] = eta
 979             // fTrigBuffer[i][3] = pt_jet
 980             // fTrigBuffer[i][4] = pt_trig
 981             // fTrigBuffer[i][5]= centrality
 982             if(fTindex==10) fTindex=0;
 983             if(fTrigBuffer[fTindex][3]>0){
 984             if (TMath::Abs(fTrigBuffer[fTindex][0]-primVtx->GetZ()<2.)){
 985             if (TMath::Abs(fTrigBuffer[fTindex][5]-centValue<5)){
 986
 987                         for(int it = 0;it<nT;++it){
 988                         AliVParticle *part = (AliVParticle*)ParticleList.At(it);
 989                         Double_t DPhi = fTrigBuffer[fTindex][1] - part->Phi();
 990                         Double_t DEta = fTrigBuffer[fTindex][2] - part->Eta();
 991                         Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta);
 992                         if(DPhi<-0.5*TMath::Pi()) DPhi+=2.*TMath::Pi();
 993                         if(DPhi>3./2.*TMath::Pi()) DPhi-=2.*TMath::Pi();
 994                         Double_t triggerEntries[7] = {centValue,fTrigBuffer[fTindex][3],part->Pt(),DR,DEta,DPhi,fTrigBuffer[fTindex][4]};
 995                         fhnMixedEvents->Fill(triggerEntries);
 996                         }
 997                         fNevents=fNevents+1;
 998                         if(fNevents==10) fTindex=fTindex+1;
 999             }}}
1000
1001                if(fTindex==10&&fNevents==10) fCountAgain=0;
1002
1003                // Copy the triggers from the current event into the buffer.
1004                //again, only if the trigger exists:
1005                if(fCountAgain==0){
1006                 if(trigJet>-1){
1007                 AliAODJet* jetT = (AliAODJet*)(fListJets[0]->At(trigJet));                      AliVParticle *partT = (AliVParticle*)ParticleList.At(trigBBTrack);
1008                 fTrigBuffer[fTrigBufferIndex][0] = primVtx->GetZ();
1009                 fTrigBuffer[fTrigBufferIndex][1] = jetT->Phi();
1010                 fTrigBuffer[fTrigBufferIndex][2] = jetT->Eta();
1011                 fTrigBuffer[fTrigBufferIndex][3] = jetT->Pt()-rho*jetT->EffectiveAreaCharged();
1012                 fTrigBuffer[fTrigBufferIndex][4] = partT->Pt();
1013                 fTrigBuffer[fTrigBufferIndex][5] = centValue;
1014                 fTrigBufferIndex++;
1015                 if(fTrigBufferIndex==9) {fTrigBufferIndex=0;
1016                                          fCountAgain=1;}
1017                 }
1018                }
1019
1020           }
1021
1022           /////////////////////////////////////////////////////////////////////////////
1023           ////////////////////// Rongrong's analysis //////////////////////////////////
1024           if(fRunAnaAzimuthalCorrelation)
1025             {
1026               fhTTPt->Fill(centValue,partback->Pt());
1027               for(Int_t ij=0; ij<fListJets[0]->GetEntries(); ij++)
1028                 {
1029                   AliAODJet* jet = (AliAODJet*)(fListJets[0]->At(ij));
1030                   Double_t jetPt   = jet->Pt();
1031                   Double_t jetEta  = jet->Eta();
1032                   Double_t jetPhi  = jet->Phi();
1033                   if(jetPt==0) continue;
1034                   if((jetEta<fJetEtaMin)||(jetEta>fJetEtaMax)) continue;
1035                   Double_t jetArea = jet->EffectiveAreaCharged();
1036                   Double_t jetPtCorr=jetPt-rho*jetArea;
1037                   Double_t dPhi=jetPhi-partback->Phi();
1038                   if(dPhi>2*TMath::Pi()) dPhi -= 2*TMath::Pi();
1039                   if(dPhi<-2*TMath::Pi()) dPhi += 2*TMath::Pi();
1040                   if(dPhi<-0.5*TMath::Pi()) dPhi += 2*TMath::Pi();
1041                   if(dPhi>1.5*TMath::Pi()) dPhi -= 2*TMath::Pi();
1042
1043                   Double_t fill[] = {partback->Pt(),jetPtCorr,dPhi,jetArea,centValue};
1044                   fHJetPhiCorr->Fill(fill);
1045                 }
1046             }
1047           /////////////////////////////////////////////////////////////////////////////
1048           /////////////////////////////////////////////////////////////////////////////
1049
1050
1051      //////////////////ANGULAR STRUCTURE//////////////////////////////////////
1052
1053      //tracks up to R=0.8 distant from the jet axis
1054    //   if(fAngStructCloseTracks==1){
1055    //    TList CloseTrackList;
1056    //    Int_t nn=GetListOfTracksCloseToJet(&CloseTrackList,jetbig);
1057    //    Double_t difR=0.04;
1058    //    for(Int_t l=0;l<15;l++){
1059    //    Double_t rr=l*0.1+0.1;
1060    //     for(int it = 0;it<nn;++it){
1061    //         AliVParticle *part1 = (AliVParticle*)CloseTrackList.At(it);
1062    //         for(int itu=it+1;itu<CloseTrackList.GetEntries();itu++){
1063    //         AliVParticle *part2 = (AliVParticle*)CloseTrackList.At(itu);
1064    //         Double_t ptm=part1->Pt();
1065    //         Double_t ptn=part2->Pt();
1066    //         Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
1067    //                    Rnm=TMath::Sqrt(Rnm);
1068    //                    Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
1069    //                    Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
1070    //                    if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
1071    //                                                   down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
1072    //                    if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
1073    //                                                   down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
1074    //                    if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
1075    //                                                    down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
1076    //                    if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
1077    //                   down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}
1078    //   }
1079
1080    //   //only jet constituents
1081    //    if(fAngStructCloseTracks==2){
1082
1083    //    Double_t difR=0.04;
1084    //    for(Int_t l=0;l<15;l++){
1085    //    Double_t rr=l*0.1+0.1;
1086
1087
1088    //    AliAODTrack* part1;
1089    //    AliAODTrack* part2;
1090
1091    //        TRefArray *genTrackListb = jetbig->GetRefTracks();
1092    //        Int_t nTracksGenJetb = genTrackListb->GetEntriesFast();
1093
1094
1095
1096    //        for(Int_t it=0; it<nTracksGenJetb; ++it){
1097    //           part1 = (AliAODTrack*)(genTrackListb->At(it));
1098    //         for(Int_t itu=0; itu<nTracksGenJetb; ++itu){
1099    //           part2 = (AliAODTrack*)(genTrackListb->At(itu));
1100    //         Double_t ptm=part1->Pt();
1101    //         Double_t ptn=part2->Pt();
1102    //         Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
1103    //                    Rnm=TMath::Sqrt(Rnm);
1104    //                    Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
1105    //                    Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
1106    //                    if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
1107    //                                                   down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
1108    //                    if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
1109    //                                                   down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
1110    //                    if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
1111    //                                                    down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
1112    //                    if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
1113    //                   down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}}
1114    // }
1115      // //end loop over R=0.4 jets
1116      // if(fAngStructCloseTracks>0){
1117      // for(Int_t l=0;l<15;l++){
1118      // Double_t rr=l*0.1+0.1;
1119      //    if(down1[l]!=0){
1120      //         if(centValue<10.)fh2AngStructpt1C10->Fill(rr,rr*up1[l]/down1[l]);
1121      //    if(centValue>20. && centValue<40.) fh2AngStructpt1C20->Fill(rr,rr*up1[l]/down1[l]);
1122      //    if(centValue>30. && centValue<60.) fh2AngStructpt1C30->Fill(rr,rr*up1[l]/down1[l]);
1123      //    if(centValue>60.) fh2AngStructpt1C60->Fill(rr,rr*up1[l]/down1[l]);}
1124      //    if(down2[l]!=0){
1125      //         if(centValue<10.) fh2AngStructpt2C10->Fill(rr,rr*up2[l]/down2[l]);
1126      //    if(centValue>20. && centValue<40.) fh2AngStructpt2C20->Fill(rr,rr*up2[l]/down2[l]);
1127      //    if(centValue>30. && centValue<60.) fh2AngStructpt2C30->Fill(rr,rr*up2[l]/down2[l]);
1128      //    if(centValue>60.) fh2AngStructpt2C60->Fill(rr,rr*up2[l]/down2[l]);}
1129      //    if(down3[l]!=0){
1130      //         if(centValue<10.) fh2AngStructpt3C10->Fill(rr,rr*up3[l]/down3[l]);
1131      //    if(centValue>20. && centValue<40.) fh2AngStructpt3C20->Fill(rr,rr*up3[l]/down3[l]);
1132      //    if(centValue>30. && centValue<60.) fh2AngStructpt3C30->Fill(rr,rr*up3[l]/down3[l]);
1133      //    if(centValue>60.) fh2AngStructpt3C60->Fill(rr,rr*up3[l]/down3[l]);}
1134      //    if(down4[l]!=0){
1135      //         if(centValue<10.) fh2AngStructpt4C10->Fill(rr,rr*up4[l]/down4[l]);
1136      //    if(centValue>20. && centValue<40.) fh2AngStructpt4C20->Fill(rr,rr*up4[l]/down4[l]);
1137      //    if(centValue>30. && centValue<60.) fh2AngStructpt4C30->Fill(rr,rr*up4[l]/down4[l]);
1138      //    if(centValue>60.) fh2AngStructpt4C60->Fill(rr,rr*up4[l]/down4[l]);}}}
1139
1140
1141
1142    }
1143
1144
1145
1146    PostData(1, fOutputList);
1147 }
1148
1149 void AliAnalysisTaskJetCore::Terminate(const Option_t *)
1150 {
1151    // Draw result to the screen
1152    // Called once at the end of the query
1153
1154    if (!GetOutputData(1))
1155    return;
1156 }
1157
1158
1159
1160
1161
1162
1163 Int_t  AliAnalysisTaskJetCore::GetListOfTracks(TList *list){
1164
1165       Int_t iCount = 0;
1166       AliAODEvent *aod = 0;
1167
1168      if(!fESD)aod = fAODIn;
1169      else aod = fAODOut;
1170
1171      if(!aod)return 0;
1172
1173      Int_t index=-1;
1174      Double_t ptmax=-10;
1175
1176
1177
1178      for(int it = 0;it < aod->GetNumberOfTracks();++it){
1179       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1180       if(!tr) AliFatal("Not a standard AOD");
1181       Bool_t bGood = false;
1182       if(fFilterType == 0)bGood = true;
1183       else if(fFilterType == 1)bGood = tr->IsHybridTPCConstrainedGlobal();
1184       else if(fFilterType == 2)bGood = tr->IsHybridGlobalConstrainedGlobal();
1185       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1186       if(fRequireITSRefit==1){if((tr->GetStatus()&AliESDtrack::kITSrefit)==0)continue;}
1187       if(bGood==false) continue;
1188       if (fApplySharedClusterCut) {
1189            Double_t frac = Double_t(tr->GetTPCnclsS()) /Double_t(tr->GetTPCncls());
1190            if (frac > 0.4) continue;
1191       }
1192      if(TMath::Abs(tr->Eta())>0.9)continue;
1193       if(tr->Pt()<0.15)continue;
1194       list->Add(tr);
1195       iCount++;
1196       if(fFilterType==2 && fFilterMaskBestPt>0){// only set the trigger track index for good quality tracks
1197         if(tr->TestFilterBit(fFilterMaskBestPt)){
1198           if(tr->Pt()>ptmax){
1199             ptmax=tr->Pt();
1200             index=iCount-1;
1201           }
1202         }
1203       }
1204       else{
1205         if(tr->Pt()>ptmax){
1206           ptmax=tr->Pt();
1207           index=iCount-1;
1208         }
1209       }
1210      }
1211
1212
1213     // else if (type == kTrackAODMCCharged) {
1214     // TClonesArray *tca = dynamic_cast<TClonesArray*>(aod->FindListObject(AliAODMCParticle::StdBranchName()));
1215     // if(!tca)return iCount;
1216     // for(int it = 0;it < tca->GetEntriesFast();++it){
1217     //   AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
1218     //   if(!part)continue;
1219     //   if(part->Pt()<0.15)continue;
1220     //   if(!part->IsPhysicalPrimary())continue;
1221     //   if(part->Charge()==0)continue;
1222     //   if(TMath::Abs(part->Eta())>0.9)continue;
1223     //   list->Add(part);
1224     //   iCount++;
1225     //   if(part->Pt()>ptmax){ ptmax=part->Pt();
1226     //  index=iCount-1;}}}
1227       return index;
1228
1229 }
1230
1231
1232
1233 Int_t  AliAnalysisTaskJetCore::SelectTrigger(TList *list,Double_t minT,Double_t maxT,Int_t &number){
1234      Int_t iCount = 0;
1235      AliAODEvent *aod = 0;
1236      if(!fESD)aod = fAODIn;
1237      else aod = fAODOut;
1238      if(!aod)return 0;
1239      Int_t index=-1;
1240      Int_t triggers[100];
1241      Int_t triggers2[100];
1242      for(Int_t cr=0;cr<100;cr++){triggers[cr]=-1;
1243        triggers2[cr]=-1;}
1244      Int_t im=0;
1245      Int_t im2=0;
1246      for(int it = 0;it < aod->GetNumberOfTracks();++it){
1247       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1248       if(!tr) AliFatal("Not a standard AOD");
1249       Bool_t bGood = false;
1250       if(fFilterType == 0)bGood = true;
1251       else if(fFilterType == 1)bGood = tr->IsHybridTPCConstrainedGlobal();
1252       else if(fFilterType == 2)bGood = tr->IsHybridGlobalConstrainedGlobal();
1253       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1254        if(fRequireITSRefit==1){if((tr->GetStatus()&AliESDtrack::kITSrefit)==0)continue;}
1255       if(bGood==false) continue;
1256       if (fApplySharedClusterCut) {
1257            Double_t frac = Double_t(tr->GetTPCnclsS()) /Double_t(tr->GetTPCncls());
1258            if (frac > 0.4) continue;
1259       }
1260       if(TMath::Abs(tr->Eta())>0.9)continue;
1261       if(tr->Pt()<0.15)continue;
1262       list->Add(tr);
1263       iCount++;
1264
1265       if(tr->Pt()>=minT && tr->Pt()<maxT){
1266         if(tr->Pt()<20){
1267         triggers[im]=iCount-1;
1268         im=im+1;}
1269         if(tr->Pt()>=20.){triggers2[im2]=iCount-1;
1270           im2=im2+1;}
1271
1272       }}
1273
1274       number=im;
1275       Int_t rd=0;
1276       if(im2==0) rd=0;
1277       if(im2>0) rd=fRandom->Integer(im2);
1278       index=triggers2[rd];
1279       AliVParticle *tr1 = (AliVParticle*)list->At(index);
1280
1281
1282
1283       for(Int_t kk=0;kk<im;kk++){
1284         //if(kk==rd) continue;
1285         if(index==triggers[kk]) continue;
1286         Int_t lab=triggers[kk];
1287          AliVParticle *tr2 = (AliVParticle*)list->At(lab);
1288
1289        Double_t detat=tr1->Eta()-tr2->Eta();
1290        Double_t dphit=RelativePhi(tr1->Phi(),tr2->Phi());
1291        Double_t deltaRt=TMath::Sqrt(detat*detat+dphit*dphit);
1292        fh1TrackPhiDistance->Fill(TMath::Abs(dphit));
1293        fh1TrackRDistance->Fill(deltaRt);
1294
1295        if(fDodiHadron==1)  if(deltaRt>0.4) number=number-1;
1296        if(fDodiHadron==2) if((deltaRt>0.4) && (TMath::Abs(dphit)>TMath::Pi()-0.6)) number=number-1;}
1297
1298
1299
1300
1301
1302
1303       return index;
1304
1305 }
1306
1307
1308
1309
1310
1311
1312
1313
1314
1315
1316
1317
1318
1319
1320
1321
1322    Int_t  AliAnalysisTaskJetCore::GetHardestTrackBackToJet(AliAODJet *jetbig){
1323
1324     AliAODEvent *aod = 0;
1325     if(!fESD)aod = fAODIn;
1326     else aod = fAODOut;
1327     Int_t index=-1;
1328     Double_t ptmax=-10;
1329     Double_t dphi=0;
1330     // Double_t dif=0;
1331     Int_t iCount=0;
1332     for(int it = 0;it < aod->GetNumberOfTracks();++it){
1333       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1334       if(!tr) AliFatal("Not a standard AOD");
1335       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1336       if(TMath::Abs(tr->Eta())>0.9)continue;
1337       if(tr->Pt()<0.15)continue;
1338       iCount=iCount+1;
1339       dphi=RelativePhi(tr->Phi(),jetbig->Phi());
1340       if(TMath::Abs(dphi)<TMath::Pi()-0.6) continue;
1341       if(tr->Pt()>ptmax){ ptmax=tr->Pt();
1342       index=iCount-1;
1343       //  dif=dphi;
1344       }}
1345
1346       return index;
1347
1348    }
1349
1350
1351
1352
1353
1354
1355
1356
1357
1358  Int_t  AliAnalysisTaskJetCore::GetListOfTracksCloseToJet(TList *list,AliAODJet *jetbig){
1359
1360     Int_t iCount = 0;
1361       AliAODEvent *aod = 0;
1362      if(!fESD)aod = fAODIn;
1363      else aod = fAODOut;
1364
1365       for(int it = 0;it < aod->GetNumberOfTracks();++it){
1366       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1367       if(!tr) AliFatal("Not a standard AOD");
1368       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1369       if(TMath::Abs(tr->Eta())>0.9)continue;
1370       if(tr->Pt()<0.15)continue;
1371       Double_t disR=jetbig->DeltaR(tr);
1372       if(disR>0.8)  continue;
1373       list->Add(tr);
1374       //cout<<fAOD->GetNumberOfTracks()<<" "<<tr->Pt()<<endl;
1375       iCount++;
1376     }
1377
1378    list->Sort();
1379    return iCount;
1380
1381 }
1382
1383
1384
1385
1386
1387
1388
1389
1390
1391
1392
1393 Int_t AliAnalysisTaskJetCore::GetNInputTracks()
1394 {
1395
1396    Int_t nInputTracks = 0;
1397      AliAODEvent *aod = 0;
1398      if(!fESD)aod = fAODIn;
1399      else aod = fAODOut;
1400    TString jbname(fJetBranchName[1]);
1401    //needs complete event, use jets without background subtraction
1402    for(Int_t i=1; i<=3; ++i){
1403       if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
1404    }
1405    // use only HI event
1406    if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
1407    if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
1408
1409    if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
1410    TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(aod->FindListObject(jbname.Data()));
1411    if(!tmpAODjets){
1412       Printf("Jet branch %s not found", jbname.Data());
1413       Printf("AliAnalysisTaskJetCore::GetNInputTracks FAILED");
1414       return -1;
1415    }
1416
1417    for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
1418       AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
1419       if(!jet) continue;
1420       TRefArray *trackList = jet->GetRefTracks();
1421       Int_t nTracks = trackList->GetEntriesFast();
1422       nInputTracks += nTracks;
1423       if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
1424    }
1425    if(fDebug) Printf("---> input tracks: %d", nInputTracks);
1426
1427    return nInputTracks;
1428 }
1429
1430
1431
1432 Double_t AliAnalysisTaskJetCore::RelativePhi(Double_t mphi,Double_t vphi){
1433
1434   if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
1435   else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
1436   if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
1437   else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
1438   double dphi = mphi-vphi;
1439   if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
1440   else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
1441   return dphi;//dphi in [-Pi, Pi]
1442 }
1443
1444 Int_t AliAnalysisTaskJetCore::GetPhiBin(Double_t phi)
1445 {
1446     Int_t phibin=-1;
1447     if(!(TMath::Abs(phi)<=2*TMath::Pi())){AliError("phi w.r.t. RP out of defined range");return -1;}
1448     Double_t phiwrtrp=TMath::ACos(TMath::Abs(TMath::Cos(phi)));
1449     phibin=Int_t(fNRPBins*phiwrtrp/(0.5*TMath::Pi()));
1450     if(phibin<0||phibin>=fNRPBins){AliError("Phi Bin not defined");}
1451     return phibin;
1452 }
1453
1454
1455
1456
1457 THnSparse* AliAnalysisTaskJetCore::NewTHnSparseF(const char* name, UInt_t entries)
1458 {
1459    // generate new THnSparseF, axes are defined in GetDimParams()
1460
1461    Int_t count = 0;
1462    UInt_t tmp = entries;
1463    while(tmp!=0){
1464       count++;
1465       tmp = tmp &~ -tmp;  // clear lowest bit
1466    }
1467
1468    TString hnTitle(name);
1469    const Int_t dim = count;
1470    Int_t nbins[dim];
1471    Double_t xmin[dim];
1472    Double_t xmax[dim];
1473
1474    Int_t i=0;
1475    Int_t c=0;
1476    while(c<dim && i<32){
1477       if(entries&(1<<i)){
1478
1479          TString label("");
1480          GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
1481          hnTitle += Form(";%s",label.Data());
1482          c++;
1483       }
1484
1485       i++;
1486    }
1487    hnTitle += ";";
1488
1489    return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
1490 }
1491
1492 void AliAnalysisTaskJetCore::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
1493 {
1494    // stores label and binning of axis for THnSparse
1495
1496    const Double_t pi = TMath::Pi();
1497
1498    switch(iEntry){
1499
1500    case 0:
1501       label = "V0 centrality (%)";
1502
1503          nbins = 10;
1504          xmin = 0.;
1505          xmax = 100.;
1506          break;
1507
1508
1509    case 1:
1510       label = "corrected jet pt";
1511          nbins = 20;
1512          xmin = 0.;
1513          xmax = 200.;
1514           break;
1515
1516
1517    case 2:
1518       label = "track pT";
1519
1520          nbins = 9;
1521          xmin = 0.;
1522          xmax = 150;
1523          break;
1524
1525
1526     case 3:
1527       label = "deltaR";
1528       nbins = 15;
1529       xmin = 0.;
1530       xmax = 1.5;
1531       break;
1532
1533
1534
1535    case 4:
1536       label = "deltaEta";
1537       nbins = 8;
1538       xmin = -1.6;
1539       xmax = 1.6;
1540       break;
1541
1542
1543   case 5:
1544       label = "deltaPhi";
1545       nbins = 90;
1546       xmin = -0.5*pi;
1547       xmax = 1.5*pi;
1548       break;
1549
1550
1551
1552     case 6:
1553       label = "leading track";
1554       nbins = 13;
1555       xmin = 0;
1556       xmax = 50;
1557       break;
1558
1559      case 7:
1560
1561       label = "trigger track";
1562       nbins =10;
1563       xmin = 0;
1564       xmax = 50;
1565       break;
1566
1567
1568
1569
1570
1571
1572
1573
1574    }
1575
1576 }
1577